/[PAMELA software]/tof/ground/toftrack.f
ViewVC logotype

Contents of /tof/ground/toftrack.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.1 - (show annotations) (download)
Thu Mar 9 12:31:47 2006 UTC (18 years, 8 months ago) by pam-de
Branch point for: v_tag, MAIN
Initial revision

1 C------------------------------------------------
2 PROGRAM TOFTRACK
3 C------------------------------------------------
4 C W. Menn
5 C
6 C Based on the "template" progran from Elena Vannuccini
7 C
8 C Version 1.00 August 2005
9 C Version 1.01 07-oct-2005: Some bugs found
10 C - new rz file must be named "***_toftrack.rz'
11 C - booking of the EVENT and CPU blocks changed according
12 C to Elena's "book_level2" routine
13 C
14 C------------------------------------------------
15
16 include 'trk_level2.f'
17 include 'common_tofroutine.f'
18 include 'common_tof.f'
19
20 c----- HBOOK
21 INTEGER HMEM
22 parameter (NWPAWC=8500000)
23 common/PAWC/HMEM(NWPAWC)
24
25 parameter (ntp_level2=22)
26 parameter (ntp_new=2)
27
28 c------------------------------------------------------------------------
29 c local variables
30 c------------------------------------------------------------------------
31 character*74 data_file !data file name
32 character*74 tof_file !tofdata file name
33 character*74 data_dir !data directory
34 character*74 data_file_level2
35 character*74 data_file_tof
36 character*74 data_file_tofcalib
37 character*74 data_file_new
38
39 parameter (lun_data_level2=72) !data file id number
40
41 c------------------------------------------------------------------
42 COMMON/QUEST/IQUEST(100)
43
44 C------------------------------------------------------------------
45 c =======================================
46 c variables for tracking routine
47 c =======================================
48 parameter(NPOINT_MAX=100)
49 DOUBLE PRECISION ZIN(NPOINT_MAX)
50 DOUBLE PRECISION XOUT1(NPOINT_MAX),YOUT1(NPOINT_MAX)
51 REAL XOUT(NPOINT_MAX),YOUT(NPOINT_MAX)
52 DOUBLE PRECISION AL_P(5)
53 REAL XFIT(6),YFIT(6),ZFIT(6),ALPHA(5),XGOOD1(6),YGOOD1(6)
54 c =======================================
55
56 c define TOF Z-coordinates
57 parameter (NPTOF=3)
58 REAL ZTOF(NPTOF)
59 c DATA ZTOF/55.,25.,-23/ !NB!!! from TOP to BOTTOM
60 DATA ZTOF/53.4,23.9,-23.4/ !Stockholm 29.9.2005
61
62 c------------------------------------------------------------------------
63 C
64 CHARACTER*100 name
65 C
66
67 INTEGER lrec, openntuple
68 INTEGER istat, ierr, icycle
69 PARAMETER (lrec=4096)
70
71
72 integer patterntrig1(6)
73 real chi
74
75 character*35 block1,block2
76
77 C-------------------------------------------------------------------------
78 c
79 c HBOOK initialization
80 c
81 c------------------------------------------------------------------------
82
83 call HLIMIT(NWPAWC)
84
85 C HERE A PATH where you can write!
86 C
87 c print *
88 c print *,'File saved in: ',name
89 C
90 C CHANGE here (hbname and hbnt if you want another ntuple id number
91 C
92
93 111 format(a)
94 print*,'data dir:'
95 read(*,111)data_dir
96 c read(*,*)data_dir
97 print*,data_dir
98 print*,'tof file?'
99 print*,'(without estention: output_YYMMDD_XXX )'
100 read(*,*)tof_file
101
102 print*,tof_file
103 print*,' '
104
105 print*,'first event:'
106 read(*,*) minevent
107 print*,minevent
108 c print*,' '
109
110 print*,'number of events to be analysed'
111 read(*,*) ntotev
112 print*,ntotev
113
114
115 c------------------------------------------------------------------------
116 c--------------- open ToF rz file
117 c------------------------------------------------------------------------
118
119 print*,
120 $ data_dir(1:LNBLNK(data_dir))//tof_file(1:LNBLNK(tof_file))
121
122
123 505 format(a,'DW_',a,'_tof.rz')
124 write(data_file_tof,505)
125 $ data_dir(1:LNBLNK(data_dir))
126 $ ,tof_file(1:LNBLNK(tof_file))
127 print*,'__________ opening TOF rz file __________'
128 print*,data_file_tof
129
130
131 call HROPEN(59,
132 $ 'EVENT',data_file_tof,'QP',4096,istat) !opens rz
133
134
135 if(istat.ne.0) stop
136 print*,' '
137 print*,'reading TOF n-tuple...'
138 call HRIN(ntp_tof,9999,0)
139 CALL HBNAME(ntp_tof,' ',0,'$CLEAR')
140 CALL HBNAME(ntp_tof,'EVENT',good,'$SET')
141 CALL HBNAME(ntp_tof,'TRIGGER',trig_evcount,'$SET')
142 CALL HBNAME(ntp_tof,'TOF',tdcid,'$SET')
143 call HPRNTU(ntp_tof)
144 call HNOENT(ntp_tof,iemax0)
145 write(*,*) 'Number of events ToF ',iemax0
146
147
148 c------------------------------------------------------------------------
149 c--------------- open ToF calib rz file
150 c------------------------------------------------------------------------
151
152 print*,
153 $ data_dir(1:LNBLNK(data_dir))//tof_file(1:LNBLNK(tof_file))
154
155 506 format(a,'DW_',a,'_tofcalib.rz')
156 write(data_file_tofcalib,506)
157 $ data_dir(1:LNBLNK(data_dir))
158 $ ,tof_file(1:LNBLNK(tof_file))
159 print*,'__________ opening TOF calib rz file __________'
160 print*,data_file_tofcalib
161
162 CALL HROPEN(59,'TOF K1A',data_file_tofcalib,'QP',4096,istat)
163
164 if (istat.ne.0) then ! check if HROPEN was OK
165 write(*,*) 'Can''t open correct ToF calibration File !!!'
166 write(*,*) 'Will use standard file tofcalib.rz !!!'
167 507 format(a,'tofcalib.rz')
168 write(data_file_tofcalib,507)
169 $ data_dir(1:LNBLNK(data_dir))
170 print*,'__________ opening TOF calib rz file __________'
171 print*,data_file_tofcalib
172
173 CALL HROPEN(59,'TOF K1A',data_file_tofcalib,'QP',4096,istat)
174 endif
175
176 print*,' reading TOF CALIB n-tuple...'
177
178 call HRIN(ntp_tofcalib,9999,0)
179
180 call HBNAME(ntp_tofcalib,' ',0,'$CLEAR')
181
182 call HBNAME(ntp_tofcalib,'TOFK1A',k1_s11s31,'$SET')
183 call HBNAME(ntp_tofcalib,'TOFK1B',k1_s12s32,'$SET')
184 call HBNAME(ntp_tofcalib,'TOFK1C',k1_s21s31,'$SET')
185 call HBNAME(ntp_tofcalib,'TOFK1D',k1_s22s32,'$SET')
186
187 call HBNAME(ntp_tofcalib,'TOFLIN11',y_coor_lin11,'$SET')
188 call HBNAME(ntp_tofcalib,'TOFLIN12',x_coor_lin12,'$SET')
189 call HBNAME(ntp_tofcalib,'TOFLIN21',x_coor_lin21,'$SET')
190 call HBNAME(ntp_tofcalib,'TOFLIN22',y_coor_lin22,'$SET')
191 call HBNAME(ntp_tofcalib,'TOFLIN31',y_coor_lin31,'$SET')
192 call HBNAME(ntp_tofcalib,'TOFLIN32',x_coor_lin32,'$SET')
193
194 call HBNAME(ntp_tofcalib,'TOFTW11',tw11,'$SET')
195 call HBNAME(ntp_tofcalib,'TOFTW12',tw12,'$SET')
196 call HBNAME(ntp_tofcalib,'TOFTW21',tw21,'$SET')
197 call HBNAME(ntp_tofcalib,'TOFTW22',tw22,'$SET')
198 call HBNAME(ntp_tofcalib,'TOFTW31',tw31,'$SET')
199 call HBNAME(ntp_tofcalib,'TOFTW32',tw32,'$SET')
200
201 call HBNAME(ntp_tofcalib,'TOFADC11',adcx11,'$SET')
202 call HBNAME(ntp_tofcalib,'TOFADC12',adcx12,'$SET')
203 call HBNAME(ntp_tofcalib,'TOFADC21',adcx21,'$SET')
204 call HBNAME(ntp_tofcalib,'TOFADC22',adcx22,'$SET')
205 call HBNAME(ntp_tofcalib,'TOFADC31',adcx31,'$SET')
206 call HBNAME(ntp_tofcalib,'TOFADC32',adcx32,'$SET')
207
208 call HPRNTU(ntp_tofcalib)
209
210 call HNOENT(ntp_tofcalib,iemax_cal)
211 c write(*,*) 'Number of Events CALIB ',iemax_cal
212
213 do iev=1,iemax_cal
214 call HGNT(ntp_tofcalib,iev,ierr) !reads an event
215 enddo
216
217 call hrout(ntp_tofcalib,icycle,' ')
218 call hrend('TOF K1A')
219 close(0)
220
221
222
223 c------------------------------------------------------------------------
224 c--------------- open Track Level2 rz file
225 c------------------------------------------------------------------------
226
227 data_file=tof_file
228 print*,'data dir:'
229 c read(*,111)data_dir
230 print*,data_dir
231 c print*,'data file?'
232 c print*,'(without estention: output_YYMMDD_XXX )'
233 c read(*,*)data_file
234
235 print*,'data_file:'
236 print*,data_file
237 print*,' '
238
239 print*,
240 $ data_dir(1:LNBLNK(data_dir))//data_file(1:LNBLNK(data_file))
241
242 c------------------------------------------------------------------------
243 504 format(a,'DW_',a,'_level2.rz')
244 write(data_file_level2,504)
245 $ data_dir(1:LNBLNK(data_dir))
246 $ ,data_file(1:LNBLNK(data_file))
247 print*,'__________ opening LEVEL2 rz file __________'
248 print*,data_file_level2
249 print*,'__________ opening LEVEL2 rz file __________'
250 print*,data_file_level2
251 IQUEST(10)=65000
252 call HROPEN(lun_data_level2,
253 $ 'LEVEL2',data_file_level2,'QP',4096,istat) !opens rz
254 if(istat.ne.0) goto 19
255 print*,'reading LEVEL2 n-tuple...'
256 call HRIN(ntp_level2,9999,0)
257
258 * -----------------------------------------------
259 CALL HBNAME(ntp_level2,' ',0,'$CLEAR')
260 CALL HBNAME(ntp_level2,'EVENT',GOOD2,'$SET')
261 CALL HBNAME(ntp_level2,'CPU',pkt_type,'$SET')
262 CALL HBNAME(ntp_level2,'TRACKS',ntrk,'$SET')
263 CALL HBNAME(ntp_level2,'SINGLETS',nclsx,'$SET')
264 * -----------------------------------------------
265
266 call HPRNTU(ntp_level2)
267 call HNOENT(ntp_level2,iemax2)
268
269 write(*,*) 'Number of events LEVEL2 ',iemax2
270 print*,'ok'
271 print*,' '
272
273
274
275 c------------------------------------------------------------------------
276 c read magnetic field map
277 c------------------------------------------------------------------------
278 c
279 c =======================
280
281
282 c read magnetic field map
283 c =======================
284 c
285 c------------------------------------------------------------------------
286 print*,'- read magnetic field map'
287 print*,' '
288 call read_B
289 print*,' '
290
291 c print*,' '
292 c print*,'- read magnetic field map'
293 c print*,' '
294 c call read_B_2maps
295 c
296
297
298 c------------------------------------------------------------------------
299 c open new RZ file
300 c------------------------------------------------------------------------
301
302 503 format(a,'DW_',a,'_toftrack.rz')
303 write(data_file_new,503)
304 $ data_dir(1:LNBLNK(data_dir))
305 $ ,data_file(1:LNBLNK(data_file))
306 print*,'__________ opening NEW rz file __________'
307 print*,data_file_new
308
309
310 call HROPEN(10,'TEST',data_file_new,'NP',4096,istat) !opens rz
311
312 call HBNT(ntp_new,'TOF',' ')
313
314
315 c call HBNAME(ntp_new,'EVENT',good,'GOOD:L,NEV_TRK:I')
316
317 call HBNAME(ntp_new,'EVENT',good2,'GOOD2:L,NEV2:I')
318
319
320 c call HBNAME(ntp_new,'CPU',pkt_type
321 c $ ,'PKT_TYPE:I
322 c $ ,PKT_NUM:I
323 c $ ,OBT:I
324 c $ ,WHICH_CALIB:I')
325
326 call HBNAME(ntp_new,'CPU',pkt_type
327 $ ,'PKT_TYPE:I::[0,50]
328 $ ,PKT_NUM:I
329 $ ,OBT:I
330 $ ,WHICH_CALIB:I::[0,50]')
331
332
333
334 call HBNAME(ntp_new,'TOF',tdcid
335 $ ,'TDCID(12):I
336 $ ,EVCOUNT(12):I
337 $ ,TDCMASK(12):I
338 $ ,ADC(4,12):I
339 $ ,TDC(4,12):I
340 $ ,TEMP1(12):I
341 $ ,TEMP2(12):I')
342
343 call HBNAME(ntp_new,'TOF',beta_a,'BETA(5)')
344 call HBNAME(ntp_new,'TOF',xtofpos,'XTOF(3)')
345 call HBNAME(ntp_new,'TOF',ytofpos,'YTOF(3)')
346 call HBNAME(ntp_new,'TOF',adc_c,'ADC_C(4,12)')
347 call HBNAME(ntp_new,'TOF',tof_i_flag,'IFLAG(6)')
348 call HBNAME(ntp_new,'TOF',tof_j_flag,'JFLAG(6)')
349
350 call HBNAME(ntp_new,'TOF',xout,'XOUT(3)')
351 call HBNAME(ntp_new,'TOF',yout,'YOUT(3)')
352
353
354 call HBNAME(ntp_new,'TRIGGER',trig_evcount
355 $ ,'TRIG_EVCOUNT:I
356 $ ,PMTPL(3):I
357 $ ,TRIGRATE(6):I
358 $ ,DLTIME(2):I
359 $ ,S4CALCOUNT(2):I
360 $ ,PMTCOUNT1(24):I
361 $ ,PMTCOUNT2(24):I
362 $ ,PATTERNBUSY(3):I
363 $ ,PATTERNTRIG(6):I
364 $ ,TRIGCONF:I')
365
366
367 417 format('NTRK:I::[0,',I4,']')
368 418 format(',IMAGE(NTRK):I::[0,',I4,']')
369 write(block1,417)NTRKMAX
370 write(block2,418)NTRKMAX
371 call HBNAME(ntp_new,'TRACKS',NTRK,
372 $ block1//
373 $ block2//'
374 $ ,XM(6,NTRK):R
375 $ ,YM(6,NTRK):R
376 $ ,ZM(6,NTRK):R
377 $ ,RESX(6,NTRK):R
378 $ ,RESY(6,NTRK):R
379 $ ,AL(5,NTRK):R
380 $ ,COVAL(5,5,NTRK):R
381 $ ,CHI2(NTRK):R
382 $ ,XGOOD(6,NTRK):I::[0,1]
383 $ ,YGOOD(6,NTRK):I::[0,1]
384 $ ,XV(6,NTRK):R
385 $ ,YV(6,NTRK):R
386 $ ,ZV(6,NTRK):R
387 $ ,AXV(6,NTRK):R
388 $ ,AYV(6,NTRK):R
389 $ ,DEDXP(6,NTRK):R
390 $ ')
391 call HBNAME(ntp_new,'SINGLETS',nclsx,
392 $ 'NCLSX(6):I,NCLSY(6):I')
393
394
395 call HPRNTU(ntp_new)
396
397 c------------------------------------------------------------------------
398
399 call HCDIR('//LEVEL2',' ')
400
401 c------------------------------------------------------------------------
402 c start loop on events
403 c------------------------------------------------------------------------
404
405 maxevent=minevent+ntotev
406 igoodevent=0
407
408 c do iev = 1,iemax0
409 do iev = minevent,MIN(iemax0,maxevent) !loop on events
410 Call Hcdir('//EVENT',' ')
411 Call Hgnt(ntp_tof,iev,Ierr)
412
413 do i=1,6
414 patterntrig1(i) = patterntrig(i)
415 enddo
416
417 * ----------------------------------------------
418 Call Hcdir('//LEVEL2',' ')
419 call HGNT(ntp_level2,iev,ierr) !reads an event
420
421 if(ierr.ne.0) goto 21
422 * ----------------------------------------------
423
424 do i=1,nptof
425 xout(i)=1000.
426 yout(i)=1000.
427 enddo
428 do i=1,5
429 ALPHA(i) = 1000.
430 enddo
431 do i=1,6
432 XFIT(i) = 1000.
433 YFIT(i) = 1000.
434 ZFIT(i) = 1000.
435 XGOOD1(i) = 1000.
436 YGOOD1(i) = 1000.
437 enddo
438
439 chi = 1000.
440
441 *=======> INSTERT HERE YOUR CODE
442
443 * ---- example ----
444 * write(*,*) 'GOOD2 ',good2
445 * if(GOOD2)then
446 if(NTRK.ne.0)then
447 print*,'Event ',iev,' # tracks ',ntrk
448 if(ntrk.eq.1)then
449
450 chi = CHI2(1)
451 do i=1,5
452 ALPHA(i) = AL(i,1)
453 enddo
454 do i=1,6
455 XFIT(i) = XV(i,1)
456 YFIT(i) = YV(i,1)
457 ZFIT(i) = ZV(i,1)
458 XGOOD1(i)=XGOOD(i,1)
459 YGOOD1(i)=YGOOD(i,1)
460 enddo
461
462 igoodevent = igoodevent+1
463 * assigned input parameters for track routine
464 * 1) Z-coordinates where the trajectory is evaluated
465 do itof=1,NPTOF
466 ZIN(itof) = ZTOF(itof)
467 enddo
468 * 2) track status vector
469 do i=1,5
470 AL_P(i) = AL(i,1)
471 enddo
472 * -------- *** tracking routine *** --------
473 call track(NPTOF,ZIN,XOUT1,YOUT1,AL_P,IFAIL)
474 * ------------------------------------------
475 do itof=1,NPTOF
476 XOUT(itof)=XOUT1(itof)
477 YOUT(itof)=YOUT1(itof)
478 enddo
479
480 do itof=1,NPTOF
481 c print*,' ',itof,ZIN(itof),XOUT(itof),YOUT(itof)
482 enddo
483 endif
484 else
485 print*,'Event ',iev,' *** BAD ***'
486 endif
487
488
489 call tofroutine(xout,yout,alpha)
490
491
492 c--------------------- fill new rz file -----------------------------
493
494 call HFNT(2)
495
496 enddo !end loop on events
497 GOTO 9000 !got to end
498
499 c------------------------------------------------------------------------
500 c
501 c data file opening error
502 c
503 c------------------------------------------------------------------------
504 19 continue
505
506 print*,' '
507 print*,'ERROR OPENING DATA FILE: ',data_file
508 print*,' '
509 print*,' '
510
511 goto 9000 !the end
512
513 c------------------------------------------------------------------------
514 c
515 c level2 ntuple event reading error
516 c
517 c------------------------------------------------------------------------
518
519 21 continue
520
521 print*,' '
522 print*,'ERROR WHILE READING LEVEL2 NTUPLE, AT EVENT
523 $ : ',iev
524 print*,' '
525 print*,' '
526
527 goto 9000 !the end
528
529 c------------------------------------------------------------------------
530 c
531 c closes files and exits
532 c
533 c------------------------------------------------------------------------
534
535 9000 continue
536
537 write(*,*) 'Good Events ',igoodevent
538
539
540 c--------------------------------------------------------------
541
542 call HCDIR('//EVENT',' ')
543 CALL HROUT(0,icycle,' ')
544 call HREND('EVENT')
545 close(0)
546
547 call HCDIR('//LEVEL2',' ')
548 call HREND('level2')
549 close(lun_data_level2)
550
551
552 call HPRNTU(2)
553 call HCDIR('//TEST',' ')
554 CALL HROUT(0,icycle,' ')
555 call HREND('TEST')
556 close(0)
557
558 STOP
559 END
560
561
562 include 'tofroutine.f'
563

  ViewVC Help
Powered by ViewVC 1.1.23