/[PAMELA software]/DarthVader/ToFLevel2/src/tofl2com.for
ViewVC logotype

Contents of /DarthVader/ToFLevel2/src/tofl2com.for

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.15 - (show annotations) (download)
Thu Jan 16 15:29:35 2014 UTC (10 years, 10 months ago) by mocchiut
Branch: MAIN
CVS Tags: v10RED, v10REDr01, HEAD
Changes since 1.14: +4 -2 lines
Compilation warnings using GCC4.7 fixed

1
2 C******************************************************************************
3 C
4 C 08-12-06 WM: adc_c-bug : The raw ADc value was multiplied with cos(theta)
5 C and AFTER that there was an if statement "if tof32(right,i,iadc) < 4095"
6 C
7 C jan-07 GF: ADCflags(4,12) inserted to flag artificial ADC values
8 C jan-07 WM: artificial ADC values created using attenuation calibration
9 C jan-07 WM: modified xtofpos flag "101". xtofpos must be inside physical
10 C dimension of the paddle +/- 10 cm
11 C jan-07 WM: if xtofpos=101 then this paddle is not used for beta
12 C calculation
13 C jan-07 WM: the definition for a "hit" is changed: Now we must have a
14 C valid TDC signal on both sides
15 C jan-07 WM: flag for PMTs #10 and #35 added, TDC=819 due to bit-shift
16 C jan-07 WM: bug removed: in some cases tdc_tw was calculated due to a
17 C leftover "xhelp" value
18 C apr-07 WM: attenuation fit curve is now a double exponential fit
19 C conversion from raw ADC to pC using calibration function
20 C variables xtr_tof and ytr_tof inserted (filled with default)
21 C jan-08 WM: Major Update: Time Walk correction introduced
22 C Additionalyl we use the information from the "check_charge"
23 C function to fill artificial ADC values and make small corrections
24 C to the k1-parameter (for Z>2)
25 C feb-08 WM: Calculation of beta(13) changed: First a mean beta is calculated,
26 C then in a second step we check the residuals of the single
27 C measurements, reject if > 10 sigma, calculate chi2 and "quality"
28 C beta is taken as good if chi2<20 and quality>10
29 C mar-08 WM: Call to "newbeta" changed, now a flag tells the function if the
30 C call comes from "tofl2com" or form "toftrack"
31 C mar-08 WM: Bug found in dEdx if check_charge>1
32 C oct-08 WM: Calculation of zenith angle debugged, sometimes strange values
33 C were possible
34 C nov-09 WM: the dEdx part ("adctof_c") moved to the new dEdx routine from Napoli
35 C feb-10 WM: k1 values now for Z=1, Z=2, Z>2, k2 values are fix
36 C******************************************************************************
37
38 INTEGER FUNCTION TOFL2COM()
39 c
40 IMPLICIT NONE
41 C
42 include 'input_tof.txt'
43 include 'output_tof.txt'
44 include 'tofcomm.txt'
45
46 INTEGER icounter
47 DATA icounter / 0/
48
49 LOGICAL check
50 REAL secure
51
52 INTEGER j,hitvec(6)
53
54 REAL dx,dy,dr,ds
55 REAL yhelp,yhelp1,yhelp2,xhelp,xhelp1,xhelp2
56 REAL c1,c2
57 REAL dist
58
59 C REAL sw,sxw,w_i
60 C INTEGER icount
61 C REAL beta_mean
62
63 INTEGER tof11_j,tof21_j,tof31_j
64 INTEGER tof12_j,tof22_j,tof32_j
65
66 c value for status of each PM-data
67 c first index : 1 = left, 2 = right
68 c second index : 1... number of paddle
69 INTEGER tof11_event(2,8),tof12_event(2,6)
70 INTEGER tof21_event(2,2),tof22_event(2,2)
71 INTEGER tof31_event(2,3),tof32_event(2,3)
72
73
74 REAL y_coor_lin11c(8,2),x_coor_lin12c(6,2)
75 REAL x_coor_lin21c(2,2),y_coor_lin22c(2,2)
76 REAL y_coor_lin31c(3,2),x_coor_lin32c(3,2)
77
78 DATA y_coor_lin11c(1,1),y_coor_lin11c(1,2) /-20.66,-2.497/
79 DATA y_coor_lin11c(2,1),y_coor_lin11c(2,2) /-9.10, -2.52/
80 DATA y_coor_lin11c(3,1),y_coor_lin11c(3,2) /-24.07,-2.12/
81 DATA y_coor_lin11c(4,1),y_coor_lin11c(4,2) /-13.40,-2.47/
82 DATA y_coor_lin11c(5,1),y_coor_lin11c(5,2) /-31.07,-2.32/
83 DATA y_coor_lin11c(6,1),y_coor_lin11c(6,2) /-21.69,-2.63/
84 DATA y_coor_lin11c(7,1),y_coor_lin11c(7,2) /-12.37,-2.65/
85 DATA y_coor_lin11c(8,1),y_coor_lin11c(8,2) /-10.81,-3.15/
86
87 DATA x_coor_lin12c(1,1),x_coor_lin12c(1,2) /12.96, -2.65/
88 DATA x_coor_lin12c(2,1),x_coor_lin12c(2,2) /17.12,-2.44/
89 DATA x_coor_lin12c(3,1),x_coor_lin12c(3,2) /7.26, -1.98/
90 DATA x_coor_lin12c(4,1),x_coor_lin12c(4,2) /-22.52,-2.27/
91 DATA x_coor_lin12c(5,1),x_coor_lin12c(5,2) /-18.54,-2.28/
92 DATA x_coor_lin12c(6,1),x_coor_lin12c(6,2) /-7.67,-2.15/
93
94 DATA x_coor_lin21c(1,1),x_coor_lin21c(1,2) /22.56,-1.56/
95 DATA x_coor_lin21c(2,1),x_coor_lin21c(2,2) /13.94,-1.56/
96
97 DATA y_coor_lin22c(1,1),y_coor_lin22c(1,2) /-24.24,-2.23/
98 DATA y_coor_lin22c(2,1),y_coor_lin22c(2,2) /-45.99,-1.68/
99
100 DATA y_coor_lin31c(1,1),y_coor_lin31c(1,2) /-22.99,-3.54/
101 DATA y_coor_lin31c(2,1),y_coor_lin31c(2,2) /-42.28,-4.10/
102 DATA y_coor_lin31c(3,1),y_coor_lin31c(3,2) /-41.29,-3.69/
103
104 DATA x_coor_lin32c(1,1),x_coor_lin32c(1,2) /0.961, -3.22/
105 DATA x_coor_lin32c(2,1),x_coor_lin32c(2,2) /4.98,-3.48/
106 DATA x_coor_lin32c(3,1),x_coor_lin32c(3,2) /-22.08,-3.37/
107
108
109 REAL theta13
110
111 c DOUBLE PRECISION ZTOF(6)
112 REAL ZTOF(6) !EM GCC4.7
113 DATA ZTOF/53.74,53.04,23.94,23.44,-23.49,-24.34/ !Sergio 9.05.2006
114
115 REAL tofarm12
116 PARAMETER (tofarm12 = 29.70) ! from 53.39 to 23.69
117 REAL tofarm23
118 PARAMETER (tofarm23 = 47.61) ! from 23.69 to -23.92
119 REAL tofarm13
120 PARAMETER (tofarm13 = 77.31) ! from 53.39 to -23.92
121
122 REAL hepratio
123
124 INTEGER ihelp
125 REAL xkorr,btemp(12)
126
127 REAL atten,pc_adc,check_charge,newbeta
128
129 INTEGER IZ
130
131
132 INTEGER ifst
133 DATA ifst /0/
134
135 C---------------------------------------
136 C
137 C Begin !
138 C
139 TOFL2COM = 0
140 C
141 C CALCULATE COMMON VARIABLES
142 C
143 C-------------------------------------------------------------------
144
145 if (ifst.eq.0) then
146
147 ifst=1
148
149 C amplitude has to be 'secure' higher than pedestal for an adc event
150 secure = 2.
151
152 C ratio between helium and proton ca. 4
153 hepratio = 4. !
154 offset = 1
155 slope = 2
156 left = 1
157 right = 2
158 none_ev = 0
159 none_find = 0
160 tdc_ev = 1
161 adc_ev = 1
162 itdc = 1
163 iadc = 2
164
165 ENDIF
166 C---------------------------------------------------------------------
167
168 icounter = icounter + 1
169
170
171 do i=1,13
172 betatof_a(i) = 100. ! As in "troftrk.for"
173 enddo
174
175 do i=1,6
176 hitvec(i) = -1
177 enddo
178
179 do i=1,4
180 do j=1,12
181 adctof_c(i,j) = 1000.
182 enddo
183 enddo
184
185
186 do i=1,4
187 do j=1,12
188 tdc_c(i,j) = 4095.
189 enddo
190 enddo
191
192
193 do i=1,12
194 do j=1,4
195 tofmask(j,i) = 0
196 enddo
197 enddo
198
199
200 c gf adc falg:
201 do i=1,4
202 do j=1,12
203 adcflagtof(i,j) = 0
204 enddo
205 enddo
206
207 c gf tdc falg:
208 do i=1,4
209 do j=1,12
210 tdcflagtof(i,j) = 0
211 enddo
212 enddo
213
214
215 C--- Fill xtr_tof and ytr_tof: positions from tracker at ToF layers
216 C--- since this is standalone ToF fill with default values
217 do j=1,6
218 xtr_tof(j) = 101.
219 ytr_tof(j) = 101.
220 enddo
221
222 c the calibration files are read in the main program from xxx_tofcalib.rz
223
224 c-------------------------get ToF data --------------------------------
225
226 c put the adc and tdc values from ntuple into tofxx(i,j,k) variables
227 c adc valueas are then pC
228
229 do j=1,8
230 tof11(1,j,2) = pc_adc(adc(ch11a(j),hb11a(j)))
231 tof11(2,j,2) = pc_adc(adc(ch11b(j),hb11b(j)))
232 tof11(1,j,1) = (tdc(ch11a(j),hb11a(j)))
233 tof11(2,j,1) = (tdc(ch11b(j),hb11b(j)))
234 enddo
235
236
237 do j=1,6
238 tof12(1,j,2) = pc_adc(adc(ch12a(j),hb12a(j)))
239 tof12(2,j,2) = pc_adc(adc(ch12b(j),hb12b(j)))
240 tof12(1,j,1) = (tdc(ch12a(j),hb12a(j)))
241 tof12(2,j,1) = (tdc(ch12b(j),hb12b(j)))
242 enddo
243
244 do j=1,2
245 tof21(1,j,2) = pc_adc(adc(ch21a(j),hb21a(j)))
246 tof21(2,j,2) = pc_adc(adc(ch21b(j),hb21b(j)))
247 tof21(1,j,1) = (tdc(ch21a(j),hb21a(j)))
248 tof21(2,j,1) = (tdc(ch21b(j),hb21b(j)))
249 enddo
250
251 do j=1,2
252 tof22(1,j,2) = pc_adc(adc(ch22a(j),hb22a(j)))
253 tof22(2,j,2) = pc_adc(adc(ch22b(j),hb22b(j)))
254 tof22(1,j,1) = (tdc(ch22a(j),hb22a(j)))
255 tof22(2,j,1) = (tdc(ch22b(j),hb22b(j)))
256 enddo
257
258 do j=1,3
259 tof31(1,j,2) = pc_adc(adc(ch31a(j),hb31a(j)))
260 tof31(2,j,2) = pc_adc(adc(ch31b(j),hb31b(j)))
261 tof31(1,j,1) = (tdc(ch31a(j),hb31a(j)))
262 tof31(2,j,1) = (tdc(ch31b(j),hb31b(j)))
263 enddo
264
265 do j=1,3
266 tof32(1,j,2) = pc_adc(adc(ch32a(j),hb32a(j)))
267 tof32(2,j,2) = pc_adc(adc(ch32b(j),hb32b(j)))
268 tof32(1,j,1) = (tdc(ch32a(j),hb32a(j)))
269 tof32(2,j,1) = (tdc(ch32b(j),hb32b(j)))
270 enddo
271
272 C----------------------------------------------------------------------
273
274 DO i = 1,8
275 if (abs(tof11(1,i,itdc)).gt.10000.) tof11(1,i,itdc)= 10000.
276 if (abs(tof11(2,i,itdc)).gt.10000.) tof11(2,i,itdc)= 10000.
277 if (abs(tof11(1,i,iadc)).gt.10000.) tof11(1,i,iadc)= 10000.
278 if (abs(tof11(2,i,iadc)).gt.10000.) tof11(2,i,iadc)= 10000.
279 ENDDO
280
281 DO i = 1,6
282 if (abs(tof12(1,i,itdc)).gt.10000.) tof12(1,i,itdc)= 10000.
283 if (abs(tof12(2,i,itdc)).gt.10000.) tof12(2,i,itdc)= 10000.
284 if (abs(tof12(1,i,iadc)).gt.10000.) tof12(1,i,iadc)= 10000.
285 if (abs(tof12(2,i,iadc)).gt.10000.) tof12(2,i,iadc)= 10000.
286 ENDDO
287
288
289 DO i = 1,2
290 if (abs(tof21(1,i,itdc)).gt.10000.) tof21(1,i,itdc)= 10000.
291 if (abs(tof21(2,i,itdc)).gt.10000.) tof21(2,i,itdc)= 10000.
292 if (abs(tof21(1,i,iadc)).gt.10000.) tof21(1,i,iadc)= 10000.
293 if (abs(tof21(2,i,iadc)).gt.10000.) tof21(2,i,iadc)= 10000.
294 ENDDO
295
296 DO i = 1,2
297 if (abs(tof22(1,i,itdc)).gt.10000.) tof22(1,i,itdc)= 10000.
298 if (abs(tof22(2,i,itdc)).gt.10000.) tof22(2,i,itdc)= 10000.
299 if (abs(tof22(1,i,iadc)).gt.10000.) tof22(1,i,iadc)= 10000.
300 if (abs(tof22(2,i,iadc)).gt.10000.) tof22(2,i,iadc)= 10000.
301 ENDDO
302
303 DO i = 1,3
304 if (abs(tof31(1,i,itdc)).gt.10000.) tof31(1,i,itdc)= 10000.
305 if (abs(tof31(2,i,itdc)).gt.10000.) tof31(2,i,itdc)= 10000.
306 if (abs(tof31(1,i,iadc)).gt.10000.) tof31(1,i,iadc)= 10000.
307 if (abs(tof31(2,i,iadc)).gt.10000.) tof31(2,i,iadc)= 10000.
308 ENDDO
309
310 DO i = 1,3
311 if (abs(tof32(1,i,itdc)).gt.10000.) tof32(1,i,itdc)= 10000.
312 if (abs(tof32(2,i,itdc)).gt.10000.) tof32(2,i,itdc)= 10000.
313 if (abs(tof32(1,i,iadc)).gt.10000.) tof32(1,i,iadc)= 10000.
314 if (abs(tof32(2,i,iadc)).gt.10000.) tof32(2,i,iadc)= 10000.
315 ENDDO
316
317 C----------------------------------------------------------------------
318 C------------------ set ADC & TDC flag = 0 ------------------------
319 C----------------------------------------------------------------------
320
321 do j=1,8
322 if (adc(ch11a(j),hb11a(j)).LT.4096)adcflagtof(ch11a(j),hb11a(j))=0
323 if (adc(ch11b(j),hb11b(j)).LT.4096)adcflagtof(ch11b(j),hb11b(j))=0
324 if (tdc(ch11a(j),hb11a(j)).LT.4096)tdcflagtof(ch11a(j),hb11a(j))=0
325 if (tdc(ch11b(j),hb11b(j)).LT.4096)tdcflagtof(ch11b(j),hb11b(j))=0
326 enddo
327 do j=1,6
328 if (adc(ch12a(j),hb12a(j)).LT.4096)adcflagtof(ch12a(j),hb12a(j))=0
329 if (adc(ch12b(j),hb12b(j)).LT.4096)adcflagtof(ch12b(j),hb12b(j))=0
330 if (tdc(ch12a(j),hb12a(j)).LT.4096)tdcflagtof(ch12a(j),hb12a(j))=0
331 if (tdc(ch12b(j),hb12b(j)).LT.4096)tdcflagtof(ch12b(j),hb12b(j))=0
332 enddo
333 do j=1,2
334 if (adc(ch21a(j),hb21a(j)).LT.4096)adcflagtof(ch21a(j),hb21a(j))=0
335 if (adc(ch21b(j),hb21b(j)).LT.4096)adcflagtof(ch21b(j),hb21b(j))=0
336 if (tdc(ch21a(j),hb21a(j)).LT.4096)tdcflagtof(ch21a(j),hb21a(j))=0
337 if (tdc(ch21b(j),hb21b(j)).LT.4096)tdcflagtof(ch21b(j),hb21b(j))=0
338 enddo
339 do j=1,2
340 if (adc(ch22a(j),hb22a(j)).LT.4096)adcflagtof(ch22a(j),hb22a(j))=0
341 if (adc(ch22b(j),hb22b(j)).LT.4096)adcflagtof(ch22b(j),hb22b(j))=0
342 if (tdc(ch22a(j),hb22a(j)).LT.4096)tdcflagtof(ch22a(j),hb22a(j))=0
343 if (tdc(ch22b(j),hb22b(j)).LT.4096)tdcflagtof(ch22b(j),hb22b(j))=0
344 enddo
345 do j=1,3
346 if (adc(ch31a(j),hb31a(j)).LT.4096)adcflagtof(ch31a(j),hb31a(j))=0
347 if (adc(ch31b(j),hb31b(j)).LT.4096)adcflagtof(ch31b(j),hb31b(j))=0
348 if (tdc(ch31a(j),hb31a(j)).LT.4096)tdcflagtof(ch31a(j),hb31a(j))=0
349 if (tdc(ch31b(j),hb31b(j)).LT.4096)tdcflagtof(ch31b(j),hb31b(j))=0
350 enddo
351 do j=1,3
352 if (adc(ch32a(j),hb32a(j)).LT.4096)adcflagtof(ch32a(j),hb32a(j))=0
353 if (adc(ch32b(j),hb32b(j)).LT.4096)adcflagtof(ch32b(j),hb32b(j))=0
354 if (tdc(ch32a(j),hb32a(j)).LT.4096)tdcflagtof(ch32a(j),hb32a(j))=0
355 if (tdc(ch32b(j),hb32b(j)).LT.4096)tdcflagtof(ch32b(j),hb32b(j))=0
356 enddo
357
358 C----------------------------------------------------------------
359 C---------- Check PMTs 10 and 35 for strange TDC values----------
360 C----------------------------------------------------------------
361
362 C---- S116A TDC=819
363 if (tof11(1,6,1).EQ.819) then
364 tof11(1,6,1) = 4095
365 tdcflagtof(ch11a(6),hb11a(6))=2
366 endif
367
368 C---- S222B TDC=819
369 if (tof22(2,2,1).EQ.819) then
370 tof22(2,2,1) = 4095
371 tdcflagtof(ch22b(2),hb22b(2))=2
372 endif
373
374 C----------------------------------------------------------------
375 C------------ Check Paddles for hits -----------------------
376 C------ a "hit" means TDC values<4095 on both sides ------------
377 C----------------------------------------------------------------
378
379 C upper tof S11
380 DO i = 1,8
381
382 DO j = 1,2
383 tof11_event(j,i) = none_ev
384 IF ((tof11(j,i,itdc).LT.2000).AND.(tof11(j,i,itdc).GT.100))
385 + tof11_event(j,i) = tof11_event(j,i) + tdc_ev
386 ENDDO
387 ENDDO
388
389 c find single paddle in upper tof with tdc and adc signal
390 tof11_i = none_find
391 tof11_j = none_find
392 check = .TRUE.
393 DO i = 1, 8
394 IF ((tof11_event(left,i).GE.1).AND.(tof11_event(right,i).GE.1))
395 + THEN
396 c check if an other paddle has also an event - then set flag
397 tof11_j = tof11_j + 2**(i-1)
398 IF (check.EQV..TRUE.) THEN
399 IF (tof11_i.EQ.none_find) THEN
400 tof11_i = i
401 ELSE
402 tof11_i = -1
403 check = .FALSE.
404 ENDIF
405 ENDIF
406 ENDIF
407 ENDDO
408
409
410 C upper tof S12
411 DO i = 1,6
412 DO j = 1,2
413 tof12_event(j,i) = none_ev
414 IF ((tof12(j,i,itdc).LT.2000).AND.(tof12(j,i,itdc).GT.100))
415 + tof12_event(j,i) = tof12_event(j,i) + tdc_ev
416 ENDDO
417 ENDDO
418
419 c find single paddle in upper tof with tdc and adc signal
420 tof12_i = none_find
421 tof12_j = none_find
422 check = .TRUE.
423 DO i = 1, 6
424 IF ((tof12_event(left,i).GE.1).AND.(tof12_event(right,i).GE.1))
425 + THEN
426 c check if an other paddle has also an event - then set flag
427 tof12_j = tof12_j + 2**(i-1)
428 IF (check.EQV..TRUE.) THEN
429 IF (tof12_i.EQ.none_find) THEN
430 tof12_i = i
431 ELSE
432 tof12_i = -1
433 check = .FALSE.
434 ENDIF
435 ENDIF
436 ENDIF
437 ENDDO
438
439
440 C middle tof S21
441 DO i = 1,2
442 DO j = 1,2
443 tof21_event(j,i) = none_ev
444 IF ((tof21(j,i,itdc).LT.2000).AND.(tof21(j,i,itdc).GT.100))
445 + tof21_event(j,i) = tof21_event(j,i) + tdc_ev
446 ENDDO
447 ENDDO
448
449 c find single paddle in upper tof with tdc and adc signal
450 tof21_i = none_find
451 tof21_j = none_find
452 check = .TRUE.
453 DO i = 1, 2
454 IF ((tof21_event(left,i).GE.1).AND.(tof21_event(right,i).GE.1))
455 + THEN
456 c check if an other paddle has also an event - then set flag
457 tof21_j = tof21_j + 2**(i-1)
458 IF (check.EQV..TRUE.) THEN
459 IF (tof21_i.EQ.none_find) THEN
460 tof21_i = i
461 ELSE
462 tof21_i = -1
463 check = .FALSE.
464 ENDIF
465 ENDIF
466 ENDIF
467 ENDDO
468
469 C middle tof S22
470 DO i = 1,2
471 DO j = 1,2
472 tof22_event(j,i) = none_ev
473 IF ((tof22(j,i,itdc).LT.2000).AND.(tof22(j,i,itdc).GT.100))
474 + tof22_event(j,i) = tof22_event(j,i) + tdc_ev
475 ENDDO
476 ENDDO
477
478 c find single paddle in upper tof with tdc and adc signal
479 tof22_i = none_find
480 tof22_j = none_find
481 check = .TRUE.
482 DO i = 1, 2
483 IF ((tof22_event(left,i).GE.1).AND.(tof22_event(right,i).GE.1))
484 + THEN
485 c check if an other paddle has also an event - then set flag
486 tof22_j = tof22_j + 2**(i-1)
487 IF (check.EQV..TRUE.) THEN
488 IF (tof22_i.EQ.none_find) THEN
489 tof22_i = i
490 ELSE
491 tof22_i = -1
492 check = .FALSE.
493 ENDIF
494 ENDIF
495 ENDIF
496 ENDDO
497
498
499 C bottom tof S31
500 DO i = 1,3
501 DO j = 1,2
502 tof31_event(j,i) = none_ev
503 IF ((tof31(j,i,itdc).LT.2000).AND.(tof31(j,i,itdc).GT.100))
504 + tof31_event(j,i) = tof31_event(j,i) + tdc_ev
505 ENDDO
506 ENDDO
507
508 c find single paddle in upper tof with tdc and adc signal
509 tof31_i = none_find
510 tof31_j = none_find
511 check = .TRUE.
512 DO i = 1, 3
513 IF ((tof31_event(left,i).GE.1).AND.(tof31_event(right,i).GE.1))
514 + THEN
515 c check if an other paddle has also an event - then set flag
516 tof31_j = tof31_j + 2**(i-1)
517 IF (check.EQV..TRUE.) THEN
518 IF (tof31_i.EQ.none_find) THEN
519 tof31_i = i
520 ELSE
521 tof31_i = -1
522 check = .FALSE.
523 ENDIF
524 ENDIF
525 ENDIF
526 ENDDO
527
528 C bottom tof S32
529 DO i = 1,3
530 DO j = 1,2
531 tof32_event(j,i) = none_ev
532 IF ((tof32(j,i,itdc).LT.2000).AND.(tof32(j,i,itdc).GT.100))
533 + tof32_event(j,i) = tof32_event(j,i) + tdc_ev
534 ENDDO
535 ENDDO
536
537 c find single paddle in upper tof with tdc and adc signal
538 tof32_i = none_find
539 tof32_j = none_find
540 check = .TRUE.
541 DO i = 1, 3
542 IF ((tof32_event(left,i).GE.1).AND.(tof32_event(right,i).GE.1))
543 + THEN
544 c check if an other paddle has also an event - then set flag
545 tof32_j = tof32_j + 2**(i-1)
546 IF (check.EQV..TRUE.) THEN
547 IF (tof32_i.EQ.none_find) THEN
548 tof32_i = i
549 ELSE
550 tof32_i = -1
551 check = .FALSE.
552 ENDIF
553 ENDIF
554 ENDIF
555 ENDDO
556
557 do i=1,6
558 tof_i_flag(i)=0
559 tof_j_flag(i)=0
560 enddo
561
562 tof_i_flag(1)=tof11_i
563 tof_i_flag(2)=tof12_i
564 tof_i_flag(3)=tof21_i
565 tof_i_flag(4)=tof22_i
566 tof_i_flag(5)=tof31_i
567 tof_i_flag(6)=tof32_i
568
569 tof_j_flag(1)=tof11_j
570 tof_j_flag(2)=tof12_j
571 tof_j_flag(3)=tof21_j
572 tof_j_flag(4)=tof22_j
573 tof_j_flag(5)=tof31_j
574 tof_j_flag(6)=tof32_j
575
576 hitvec(1)=tof11_i
577 hitvec(2)=tof12_i
578 hitvec(3)=tof21_i
579 hitvec(4)=tof22_i
580 hitvec(5)=tof31_i
581 hitvec(6)=tof32_i
582
583
584 C------------------------------------------------------------------
585 C-- calculate track position in paddle using timing difference
586 C-- this calculation is preliminary and uses some standard
587 C-- calibration values, but we need to find a rough position to
588 C-- be able to calculate artificial ADC values (needed for the
589 C-- timewalk...
590 C------------------------------------------------------------------
591
592 do i=1,3
593 xtofpos(i)=100.
594 ytofpos(i)=100.
595 enddo
596
597 C-----------------------------S1 --------------------------------
598
599 IF (tof11_i.GT.none_find) THEN
600 ytofpos(1) = ((tof11(1,tof11_i,itdc)-tof11(2,tof11_i,itdc))/2.
601 + -y_coor_lin11c(tof11_i,offset))/y_coor_lin11c(tof11_i,slope)
602 endif
603
604 IF (tof12_i.GT.none_find) THEN
605 xtofpos(1) = ((tof12(1,tof12_i,itdc)-tof12(2,tof12_i,itdc))/2.
606 + -x_coor_lin12c(tof12_i,offset))/x_coor_lin12c(tof12_i,slope)
607 endif
608
609
610 C-----------------------------S2 --------------------------------
611
612 IF (tof21_i.GT.none_find) THEN
613 xtofpos(2) = ((tof21(1,tof21_i,itdc)-tof21(2,tof21_i,itdc))/2.
614 + -x_coor_lin21c(tof21_i,offset))/x_coor_lin21c(tof21_i,slope)
615 endif
616
617 IF (tof22_i.GT.none_find) THEN
618 ytofpos(2) = ((tof22(1,tof22_i,itdc)-tof22(2,tof22_i,itdc))/2.
619 + -y_coor_lin22c(tof22_i,offset))/y_coor_lin22c(tof22_i,slope)
620 endif
621
622
623 C-----------------------------S3 --------------------------------
624
625 IF (tof31_i.GT.none_find) THEN
626 ytofpos(3) = ((tof31(1,tof31_i,itdc)-tof31(2,tof31_i,itdc))/2.
627 + -y_coor_lin31c(tof31_i,offset))/y_coor_lin31c(tof31_i,slope)
628 endif
629
630 IF (tof32_i.GT.none_find) THEN
631 xtofpos(3) = ((tof32(1,tof32_i,itdc)-tof32(2,tof32_i,itdc))/2.
632 + -x_coor_lin32c(tof32_i,offset))/x_coor_lin32c(tof32_i,slope)
633 endif
634
635
636 C----------------------------------------------------------------------
637 C--------------------- zenith angle theta ---------------------------
638 C----------------------------------------------------------------------
639
640 xhelp1=0.
641 if (tof11_i.GT.none_find) xhelp1=tof11_x(tof11_i)
642 if (xtofpos(1).lt.100) xhelp1=xtofpos(1)
643
644 yhelp1=0.
645 if (tof12_i.GT.none_find) yhelp1=tof12_y(tof12_i)
646 if (ytofpos(1).lt.100) yhelp1=ytofpos(1)
647
648
649 yhelp2=0.
650 if (tof32_i.GT.none_find) yhelp2=tof32_y(tof32_i)
651 if (ytofpos(3).lt.100) yhelp2=ytofpos(3)
652
653 xhelp2=0.
654 if (tof31_i.GT.none_find) xhelp2=tof31_x(tof31_i)
655 if (xtofpos(3).lt.100) xhelp2=xtofpos(3)
656
657
658 dx=0.
659 dy=0.
660 dr=0.
661 theta13 = 0.
662
663 dx = xhelp1 - xhelp2
664 dy = yhelp1 - yhelp2
665 dr = sqrt(dx*dx+dy*dy)
666 theta13 = atan(dr/tofarm13)
667
668
669 C----------------------------------------------------------------------
670 C--- check charge:
671 C--- if Z=2 we should use the attenuation curve for helium to
672 C--- fill the artificail ADC values and NOT divide by "hepratio"
673 C--- if Z>2 we should do a correction to
674 C--- the k1 constants in the beta calculation
675 C----------------------------------------------------------------------
676
677 iz = int(check_charge(theta13,hitvec))
678 C write(*,*) 'charge in tofl2com',iz
679
680 C--------------------------------------------------------------------
681 C---- if TDCleft.and.TDCright and NO ADC insert artificial ADC
682 C---- values
683 C--------------------------------------------------------------------
684 c middle y (or x) position of the upper and middle ToF-Paddle
685 c DATA tof11_x/ -17.85,-12.75,-7.65,-2.55,2.55,7.65,12.75,17.85/
686 c DATA tof12_y/ -13.75,-8.25,-2.75,2.75,8.25,13.75/
687 c DATA tof21_y/ 3.75,-3.75/ ! paddles in different order
688 c DATA tof22_x/ -4.5,4.5/
689 c DATA tof31_x/ -6.0,0.,6.0/
690 c DATA tof32_y/ -5.0,0.0,5.0/
691
692
693 C---------------------------- S1 -------------------------------------
694
695 c yhelp=0.
696 yhelp=100. ! WM
697 if (tof12_i.GT.none_find) yhelp=tof12_y(tof12_i)
698 if (ytofpos(1).lt.100) yhelp=ytofpos(1)
699
700 IF (tof11_i.GT.none_find.AND.abs(yhelp).lt.100) THEN
701 i = tof11_i
702 if (adc(ch11a(i),hb11a(i)).eq.4095) then
703 xkorr = atten(left,11,i,yhelp)
704 if (iz.le.1) xkorr=xkorr/hepratio
705 tof11(left,i,iadc)=xkorr/cos(theta13)
706 adcflagtof(ch11a(i),hb11a(i)) = 1
707 endif
708 if (adc(ch11b(i),hb11b(i)).eq.4095) then
709 xkorr = atten(right,11,i,yhelp)
710 if (iz.le.1) xkorr=xkorr/hepratio
711 tof11(right,i,iadc)=xkorr/cos(theta13)
712 adcflagtof(ch11b(i),hb11b(i)) = 1
713 endif
714 ENDIF
715
716 c xhelp=0.
717 xhelp=100. ! WM
718 if (tof11_i.GT.none_find) xhelp=tof11_x(tof11_i)
719 if (xtofpos(1).lt.100) xhelp=xtofpos(1)
720
721 IF (tof12_i.GT.none_find.AND.abs(xhelp).lt.100) THEN
722 i = tof12_i
723 if (adc(ch12a(i),hb12a(i)).eq.4095) then
724 xkorr = atten(left,12,i,xhelp)
725 if (iz.le.1) xkorr=xkorr/hepratio
726 tof12(left,i,iadc) = xkorr/cos(theta13)
727 adcflagtof(ch12a(i),hb12a(i)) = 1
728 endif
729 if (adc(ch12b(i),hb12b(i)).eq.4095) then
730 xkorr = atten(right,12,i,xhelp)
731 if (iz.le.1) xkorr=xkorr/hepratio
732 tof12(right,i,iadc) = xkorr/cos(theta13)
733 adcflagtof(ch12b(i),hb12b(i)) = 1
734 endif
735 ENDIF
736
737 C-----------------------------S2 --------------------------------
738
739 c xhelp=0.
740 xhelp=100. ! WM
741 if (tof22_i.GT.none_find) xhelp=tof22_x(tof22_i)
742 if (xtofpos(2).lt.100) xhelp=xtofpos(2)
743
744 IF (tof21_i.GT.none_find.AND.abs(xhelp).lt.100) THEN
745 i = tof21_i
746 if (adc(ch21a(i),hb21a(i)).eq.4095) then
747 xkorr = atten(left,21,i,xhelp)
748 if (iz.le.1) xkorr=xkorr/hepratio
749 tof21(left,i,iadc) = xkorr/cos(theta13)
750 adcflagtof(ch21a(i),hb21a(i)) = 1
751 endif
752 if (adc(ch21b(i),hb21b(i)).eq.4095) then
753 xkorr = atten(right,21,i,xhelp)
754 if (iz.le.1) xkorr=xkorr/hepratio
755 tof21(right,i,iadc) = xkorr/cos(theta13)
756 adcflagtof(ch21b(i),hb21b(i)) = 1
757 endif
758 ENDIF
759
760
761 c yhelp=0.
762 yhelp=100. ! WM
763 if (tof21_i.GT.none_find) yhelp=tof21_y(tof21_i)
764 if (ytofpos(2).lt.100) yhelp=ytofpos(2)
765
766 IF (tof22_i.GT.none_find.AND.abs(yhelp).lt.100) THEN
767 i = tof22_i
768 if (adc(ch22a(i),hb22a(i)).eq.4095) then
769 xkorr = atten(left,22,i,yhelp)
770 if (iz.le.1) xkorr=xkorr/hepratio
771 tof22(left,i,iadc) = xkorr/cos(theta13)
772 adcflagtof(ch22a(i),hb22a(i)) = 1
773 endif
774 if (adc(ch22b(i),hb22b(i)).eq.4095) then
775 xkorr = atten(right,22,i,yhelp)
776 if (iz.le.1) xkorr=xkorr/hepratio
777 tof22(right,i,iadc) = xkorr/cos(theta13)
778 adcflagtof(ch22b(i),hb22b(i)) = 1
779 endif
780 ENDIF
781
782 C-----------------------------S3 --------------------------------
783
784 c yhelp=0.
785 yhelp=100. ! WM
786 if (tof32_i.GT.none_find) yhelp=tof32_y(tof32_i)
787 if (ytofpos(3).lt.100) yhelp=ytofpos(3)
788
789 IF (tof31_i.GT.none_find.AND.abs(yhelp).lt.100) THEN
790 i = tof31_i
791 if (adc(ch31a(i),hb31a(i)).eq.4095) then
792 xkorr = atten(left,31,i,yhelp)
793 if (iz.le.1) xkorr=xkorr/hepratio
794 tof31(left,i,iadc) = xkorr/cos(theta13)
795 adcflagtof(ch31a(i),hb31a(i)) = 1
796 endif
797 if (adc(ch31b(i),hb31b(i)).eq.4095) then
798 xkorr = atten(right,31,i,yhelp)
799 if (iz.le.1) xkorr=xkorr/hepratio
800 tof31(right,i,iadc) = xkorr/cos(theta13)
801 adcflagtof(ch31b(i),hb31b(i)) = 1
802 endif
803 ENDIF
804
805 c xhelp=0.
806 xhelp=100. ! WM
807 if (tof31_i.GT.none_find) xhelp=tof31_x(tof31_i)
808 if (xtofpos(3).lt.100) xhelp=xtofpos(3)
809
810 IF (tof32_i.GT.none_find.AND.abs(xhelp).lt.100) THEN
811 i = tof32_i
812 if (adc(ch32a(i),hb32a(i)).eq.4095) then
813 xkorr = atten(left,32,i,xhelp)
814 if (iz.le.1) xkorr=xkorr/hepratio
815 tof32(left,i,iadc) = xkorr/cos(theta13)
816 adcflagtof(ch32a(i),hb32a(i)) = 1
817 endif
818 if (adc(ch32b(i),hb32b(i)).eq.4095) then
819 xkorr = atten(right,32,i,xhelp)
820 if (iz.le.1) xkorr=xkorr/hepratio
821 tof32(right,i,iadc) = xkorr/cos(theta13)
822 adcflagtof(ch32b(i),hb32b(i)) = 1
823 endif
824 ENDIF
825
826
827 C-------------------------------------------------------------------
828 C--------------------Time walk correction -------------------------
829 C-------------------------------------------------------------------
830 C-------------------------------------------------------------------
831 C Now there is for each hitted paddle a TDC and ADC value, if the
832 C TDC was < 4095.
833 C There might be also TDC-ADC pairs in paddles not hitted
834
835 C-------------------------------------------------------------------
836 C If we have multiple paddles hit, so that no artificial ADC value
837 C is created, we set the raw TDC value as "tdc_c"
838 C-------------------------------------------------------------------
839 c
840 c do i=1,4
841 c do j=1,12
842 c tdc_c(i,j) = tdc(i,j)
843 c enddo
844 c enddo
845 c
846 C---- Let's correct the raw TDC value with the time walk ---------
847
848 DO i=1,8
849 if ((tdc(ch11a(i),hb11a(i)).lt.4095).and.
850 & (tof11(left,i,iadc).lt.3786)) THEN
851 xhelp = tw11(left,i)/(tof11(left,i,iadc)**0.5)
852 tof11(left,i,itdc) = tof11(left,i,itdc) + xhelp
853 tdc_c(ch11a(i),hb11a(i))=tof11(left,i,itdc)
854 ENDIF
855
856 if ((tdc(ch11b(i),hb11b(i)).lt.4095).and.
857 & (tof11(right,i,iadc).lt.3786)) THEN
858 xhelp = tw11(right,i)/(tof11(right,i,iadc)**0.5)
859 tof11(right,i,itdc) = tof11(right,i,itdc) + xhelp
860 tdc_c(ch11b(i),hb11b(i))=tof11(right,i,itdc)
861 ENDIF
862 ENDDO
863
864
865 DO i=1,6
866 if ((tdc(ch12a(i),hb12a(i)).lt.4095).and.
867 & (tof12(left,i,iadc).lt.3786)) THEN
868 xhelp = tw12(left,i)/(tof12(left,i,iadc)**0.5)
869 tof12(left,i,itdc) = tof12(left,i,itdc) + xhelp
870 tdc_c(ch12a(i),hb12a(i))=tof12(left,i,itdc)
871 ENDIF
872
873 if ((tdc(ch12b(i),hb12b(i)).lt.4095).and.
874 & (tof12(right,i,iadc).lt.3786)) THEN
875 xhelp = tw12(right,i)/(tof12(right,i,iadc)**0.5)
876 tof12(right,i,itdc) = tof12(right,i,itdc) + xhelp
877 tdc_c(ch12b(i),hb12b(i))=tof12(right,i,itdc)
878 ENDIF
879 ENDDO
880
881 C----
882 DO I=1,2
883 if ((tdc(ch21a(i),hb21a(i)).lt.4095).and.
884 & (tof21(left,i,iadc).lt.3786)) THEN
885 xhelp = tw21(left,i)/(tof21(left,i,iadc)**0.5)
886 tof21(left,i,itdc) = tof21(left,i,itdc) + xhelp
887 tdc_c(ch21a(i),hb21a(i))=tof21(left,i,itdc)
888 ENDIF
889
890 if ((tdc(ch21b(i),hb21b(i)).lt.4095).and.
891 & (tof21(right,i,iadc).lt.3786)) THEN
892 xhelp = tw21(right,i)/(tof21(right,i,iadc)**0.5)
893 tof21(right,i,itdc) = tof21(right,i,itdc) + xhelp
894 tdc_c(ch21b(i),hb21b(i))=tof21(right,i,itdc)
895 ENDIF
896 ENDDO
897
898 DO I=1,2
899 if ((tdc(ch22a(i),hb22a(i)).lt.4095).and.
900 & (tof22(left,i,iadc).lt.3786)) THEN
901 xhelp = tw22(left,i)/(tof22(left,i,iadc)**0.5)
902 tof22(left,i,itdc) = tof22(left,i,itdc) + xhelp
903 tdc_c(ch22a(i),hb22a(i))=tof22(left,i,itdc)
904 ENDIF
905
906 if ((tdc(ch22b(i),hb22b(i)).lt.4095).and.
907 & (tof22(right,i,iadc).lt.3786)) THEN
908 xhelp = tw22(right,i)/(tof22(right,i,iadc)**0.5)
909 tof22(right,i,itdc) = tof22(right,i,itdc) + xhelp
910 tdc_c(ch22b(i),hb22b(i))=tof22(right,i,itdc)
911 ENDIF
912 ENDDO
913
914 C----
915 DO I=1,3
916 if ((tdc(ch31a(i),hb31a(i)).lt.4095).and.
917 & (tof31(left,i,iadc).lt.3786)) THEN
918 xhelp = tw31(left,i)/(tof31(left,i,iadc)**0.5)
919 tof31(left,i,itdc) = tof31(left,i,itdc) + xhelp
920 tdc_c(ch31a(i),hb31a(i))=tof31(left,i,itdc)
921 ENDIF
922
923 if ((tdc(ch31b(i),hb31b(i)).lt.4095).and.
924 & (tof31(right,i,iadc).lt.3786)) THEN
925 xhelp = tw31(right,i)/(tof31(right,i,iadc)**0.5)
926 tof31(right,i,itdc) = tof31(right,i,itdc) + xhelp
927 tdc_c(ch31b(i),hb31b(i))=tof31(right,i,itdc)
928 ENDIF
929 ENDDO
930
931 DO I=1,3
932 if ((tdc(ch32a(i),hb32a(i)).lt.4095).and.
933 & (tof32(left,i,iadc).lt.3786)) THEN
934 xhelp = tw32(left,i)/(tof32(left,i,iadc)**0.5)
935 tof32(left,i,itdc) = tof32(left,i,itdc) + xhelp
936 tdc_c(ch32a(i),hb32a(i))=tof32(left,i,itdc)
937 ENDIF
938
939 if ((tdc(ch32b(i),hb32b(i)).lt.4095).and.
940 & (tof32(right,i,iadc).lt.3786)) THEN
941 xhelp = tw32(right,i)/(tof32(right,i,iadc)**0.5)
942 tof32(right,i,itdc) = tof32(right,i,itdc) + xhelp
943 tdc_c(ch32b(i),hb32b(i))=tof32(right,i,itdc)
944 ENDIF
945 ENDDO
946
947
948 C---------------------------------------------------------------
949 C--- calculate track position in paddle using timing difference
950 C--- now using the time-walk corrected TDC values
951 C---------------------------------------------------------------
952
953 do i=1,3
954 xtofpos(i)=100.
955 ytofpos(i)=100.
956 enddo
957
958 C-----------------------------S1 --------------------------------
959
960 IF (tof11_i.GT.none_find) THEN
961 ytofpos(1) = ((tof11(1,tof11_i,itdc)-tof11(2,tof11_i,itdc))/2.
962 + -y_coor_lin11(tof11_i,offset))/y_coor_lin11(tof11_i,slope)
963 i=tof11_i
964 endif
965
966 IF (tof12_i.GT.none_find) THEN
967 xtofpos(1) = ((tof12(1,tof12_i,itdc)-tof12(2,tof12_i,itdc))/2.
968 + -x_coor_lin12(tof12_i,offset))/x_coor_lin12(tof12_i,slope)
969 i=tof12_i
970 endif
971
972
973 C-----------------------------S2 --------------------------------
974
975 IF (tof21_i.GT.none_find) THEN
976 xtofpos(2) = ((tof21(1,tof21_i,itdc)-tof21(2,tof21_i,itdc))/2.
977 + -x_coor_lin21(tof21_i,offset))/x_coor_lin21(tof21_i,slope)
978 i=tof21_i
979 endif
980
981 IF (tof22_i.GT.none_find) THEN
982 ytofpos(2) = ((tof22(1,tof22_i,itdc)-tof22(2,tof22_i,itdc))/2.
983 + -y_coor_lin22(tof22_i,offset))/y_coor_lin22(tof22_i,slope)
984 i=tof22_i
985 endif
986
987
988 C-----------------------------S3 --------------------------------
989
990 IF (tof31_i.GT.none_find) THEN
991 ytofpos(3) = ((tof31(1,tof31_i,itdc)-tof31(2,tof31_i,itdc))/2.
992 + -y_coor_lin31(tof31_i,offset))/y_coor_lin31(tof31_i,slope)
993 i=tof31_i
994 endif
995
996 IF (tof32_i.GT.none_find) THEN
997 xtofpos(3) = ((tof32(1,tof32_i,itdc)-tof32(2,tof32_i,itdc))/2.
998 + -x_coor_lin32(tof32_i,offset))/x_coor_lin32(tof32_i,slope)
999 i=tof32_i
1000 endif
1001
1002
1003 c do i=1,3
1004 c if (abs(xtofpos(i)).gt.100.) then
1005 c xtofpos(i)=101.
1006 c endif
1007 c if (abs(ytofpos(i)).gt.100.) then
1008 c ytofpos(i)=101.
1009 c endif
1010 c enddo
1011
1012
1013 C-- restrict TDC measurements to physical paddle dimensions +/- 10 cm
1014 C-- this cut is now stronger than in the old versions
1015
1016 if (abs(xtofpos(1)).gt.31.) xtofpos(1)=101.
1017 if (abs(xtofpos(2)).gt.19.) xtofpos(2)=101.
1018 if (abs(xtofpos(3)).gt.19.) xtofpos(3)=101.
1019
1020 if (abs(ytofpos(1)).gt.26.) ytofpos(1)=101.
1021 if (abs(ytofpos(2)).gt.18.) ytofpos(2)=101.
1022 if (abs(ytofpos(3)).gt.18.) ytofpos(3)=101.
1023
1024 C----------------------------------------------------------------------
1025 C--------------------- zenith angle theta ---------------------------
1026 C----------------------------------------------------------------------
1027 C------------------- improved calculation ---------------------------
1028
1029 xhelp1=0.
1030 if (tof11_i.GT.none_find) xhelp1=tof11_x(tof11_i)
1031 if (xtofpos(1).lt.100) xhelp1=xtofpos(1)
1032
1033 yhelp1=0.
1034 if (tof12_i.GT.none_find) yhelp1=tof12_y(tof12_i)
1035 if (ytofpos(1).lt.100) yhelp1=ytofpos(1)
1036
1037 yhelp2=0.
1038 if (tof32_i.GT.none_find) yhelp2=tof32_y(tof32_i)
1039 if (ytofpos(3).lt.100) yhelp2=ytofpos(3)
1040
1041 xhelp2=0.
1042 if (tof31_i.GT.none_find) xhelp2=tof31_x(tof31_i)
1043 if (xtofpos(3).lt.100) xhelp2=xtofpos(3)
1044
1045
1046 dx=0.
1047 dy=0.
1048 dr=0.
1049 theta13 = 0.
1050
1051 dx = xhelp1 - xhelp2
1052 dy = yhelp1 - yhelp2
1053 dr = sqrt(dx*dx+dy*dy)
1054 theta13 = atan(dr/tofarm13)
1055
1056
1057 C------------------------------------------------------------------
1058 C------------------------------------------------------------------
1059 C-------angle and ADC(x) correction: moved to new dEdx routine
1060 C------------------------------------------------------------------
1061 C------------------------------------------------------------------
1062
1063 C--------------------------------------------------------------------
1064 C----------------------calculate Beta ------------------------------
1065 C--------------------------------------------------------------------
1066 C-------------------difference of sums -----------------------------
1067 C
1068 C DS = (t1+t2) - t3+t4)
1069 C DS = c1 + c2/beta*cos(theta)
1070 C c2 = 2d/c gives c2 = 2d/(c*TDCresolution) TDC=50ps/channel
1071 C => c2 = ca.60 for 0.45 m c2 = ca.109 for 0.81 m
1072 C since TDC resolution varies slightly c2 has to be calibrated
1073
1074 C S11 - S31
1075
1076 dist = ZTOF(1) - ZTOF(5)
1077 c2 = (2.*0.01*dist)/(3.E08*50.E-12 )
1078
1079 IF ((tof11_i.GT.none_find).AND.(tof31_i.GT.none_find).AND.
1080 & (ytofpos(1).NE.101.).AND.(ytofpos(3).NE.101.)) THEN
1081 xhelp1 = tof11(1,tof11_i,itdc)+tof11(2,tof11_i,itdc)
1082 xhelp2 = tof31(1,tof31_i,itdc)+tof31(2,tof31_i,itdc)
1083 ds = xhelp1-xhelp2
1084 ihelp=(tof11_i-1)*3+tof31_i
1085 if (iz.le.1) c1 = k_S11S31(1,ihelp)
1086 if (iz.eq.2) c1 = k_S11S31(2,ihelp)
1087 if (iz.gt.2) c1 = k_S11S31(3,ihelp)
1088 betatof_a(1) = c2/(cos(theta13)*(ds-c1))
1089
1090 C------- ToF Mask - S11 - S31
1091
1092 tofmask(ch11a(tof11_i),hb11a(tof11_i)) =
1093 $ tofmask(ch11a(tof11_i),hb11a(tof11_i)) + 1
1094 tofmask(ch11b(tof11_i),hb11b(tof11_i)) =
1095 $ tofmask(ch11b(tof11_i),hb11b(tof11_i)) + 1
1096
1097 tofmask(ch31a(tof31_i),hb31a(tof31_i)) =
1098 $ tofmask(ch31a(tof31_i),hb31a(tof31_i)) + 1
1099 tofmask(ch31b(tof31_i),hb31b(tof31_i)) =
1100 $ tofmask(ch31b(tof31_i),hb31b(tof31_i)) + 1
1101
1102 C-------
1103
1104 ENDIF
1105
1106 C S11 - S32
1107
1108 dist = ZTOF(1) - ZTOF(6)
1109 c2 = (2.*0.01*dist)/(3.E08*50.E-12 )
1110
1111 IF ((tof11_i.GT.none_find).AND.(tof32_i.GT.none_find).AND.
1112 & (ytofpos(1).NE.101.).AND.(xtofpos(3).NE.101.)) THEN
1113 xhelp1 = tof11(1,tof11_i,itdc)+tof11(2,tof11_i,itdc)
1114 xhelp2 = tof32(1,tof32_i,itdc)+tof32(2,tof32_i,itdc)
1115 ds = xhelp1-xhelp2
1116 ihelp=(tof11_i-1)*3+tof32_i
1117 if (iz.le.1) c1 = k_S11S32(1,ihelp)
1118 if (iz.eq.2) c1 = k_S11S32(2,ihelp)
1119 if (iz.gt.2) c1 = k_S11S32(3,ihelp)
1120 betatof_a(2) = c2/(cos(theta13)*(ds-c1))
1121
1122 C------- ToF Mask - S11 - S32
1123
1124 tofmask(ch11a(tof11_i),hb11a(tof11_i)) =
1125 $ tofmask(ch11a(tof11_i),hb11a(tof11_i)) + 1
1126 tofmask(ch11b(tof11_i),hb11b(tof11_i)) =
1127 $ tofmask(ch11b(tof11_i),hb11b(tof11_i)) + 1
1128
1129 tofmask(ch32a(tof32_i),hb32a(tof32_i)) =
1130 $ tofmask(ch32a(tof32_i),hb32a(tof32_i)) + 1
1131 tofmask(ch32b(tof32_i),hb32b(tof32_i)) =
1132 $ tofmask(ch32b(tof32_i),hb32b(tof32_i)) + 1
1133
1134 C-------
1135
1136 ENDIF
1137
1138 C S12 - S31
1139
1140 dist = ZTOF(2) - ZTOF(5)
1141 c2 = (2.*0.01*dist)/(3.E08*50.E-12 )
1142
1143 IF ((tof12_i.GT.none_find).AND.(tof31_i.GT.none_find).AND.
1144 & (xtofpos(1).NE.101.).AND.(ytofpos(3).NE.101.)) THEN
1145 xhelp1 = tof12(1,tof12_i,itdc)+tof12(2,tof12_i,itdc)
1146 xhelp2 = tof31(1,tof31_i,itdc)+tof31(2,tof31_i,itdc)
1147 ds = xhelp1-xhelp2
1148 ihelp=(tof12_i-1)*3+tof31_i
1149 if (iz.le.1) c1 = k_S12S31(1,ihelp)
1150 if (iz.eq.2) c1 = k_S12S31(2,ihelp)
1151 if (iz.gt.2) c1 = k_S12S31(3,ihelp)
1152 betatof_a(3) = c2/(cos(theta13)*(ds-c1))
1153
1154 C------- ToF Mask - S12 - S31
1155
1156 tofmask(ch12a(tof12_i),hb12a(tof12_i)) =
1157 $ tofmask(ch12a(tof12_i),hb12a(tof12_i)) + 1
1158 tofmask(ch12b(tof12_i),hb12b(tof12_i)) =
1159 $ tofmask(ch12b(tof12_i),hb12b(tof12_i)) + 1
1160
1161 tofmask(ch31a(tof31_i),hb31a(tof31_i)) =
1162 $ tofmask(ch31a(tof31_i),hb31a(tof31_i)) + 1
1163 tofmask(ch31b(tof31_i),hb31b(tof31_i)) =
1164 $ tofmask(ch31b(tof31_i),hb31b(tof31_i)) + 1
1165
1166 C-------
1167
1168 ENDIF
1169
1170 C S12 - S32
1171
1172 dist = ZTOF(2) - ZTOF(6)
1173 c2 = (2.*0.01*dist)/(3.E08*50.E-12 )
1174
1175 IF ((tof12_i.GT.none_find).AND.(tof32_i.GT.none_find).AND.
1176 & (xtofpos(1).NE.101.).AND.(xtofpos(3).NE.101.)) THEN
1177 xhelp1 = tof12(1,tof12_i,itdc)+tof12(2,tof12_i,itdc)
1178 xhelp2 = tof32(1,tof32_i,itdc)+tof32(2,tof32_i,itdc)
1179 ds = xhelp1-xhelp2
1180 ihelp=(tof12_i-1)*3+tof32_i
1181 if (iz.le.1) c1 = k_S12S32(1,ihelp)
1182 if (iz.eq.2) c1 = k_S12S32(2,ihelp)
1183 if (iz.gt.2) c1 = k_S12S32(3,ihelp)
1184 betatof_a(4) = c2/(cos(theta13)*(ds-c1))
1185
1186 C------- ToF Mask - S12 - S32
1187
1188 tofmask(ch12a(tof12_i),hb12a(tof12_i)) =
1189 $ tofmask(ch12a(tof12_i),hb12a(tof12_i)) + 1
1190 tofmask(ch12b(tof12_i),hb12b(tof12_i)) =
1191 $ tofmask(ch12b(tof12_i),hb12b(tof12_i)) + 1
1192
1193 tofmask(ch32a(tof32_i),hb32a(tof32_i)) =
1194 $ tofmask(ch32a(tof32_i),hb32a(tof32_i)) + 1
1195 tofmask(ch32b(tof32_i),hb32b(tof32_i)) =
1196 $ tofmask(ch32b(tof32_i),hb32b(tof32_i)) + 1
1197
1198 C-------
1199
1200 ENDIF
1201
1202 C S21 - S31
1203
1204 dist = ZTOF(3) - ZTOF(5)
1205 c2 = (2.*0.01*dist)/(3.E08*50.E-12 )
1206
1207 IF ((tof21_i.GT.none_find).AND.(tof31_i.GT.none_find).AND.
1208 & (xtofpos(2).NE.101.).AND.(ytofpos(3).NE.101.)) THEN
1209 xhelp1 = tof21(1,tof21_i,itdc)+tof21(2,tof21_i,itdc)
1210 xhelp2 = tof31(1,tof31_i,itdc)+tof31(2,tof31_i,itdc)
1211 ds = xhelp1-xhelp2
1212 ihelp=(tof21_i-1)*3+tof31_i
1213 if (iz.le.1) c1 = k_S21S31(1,ihelp)
1214 if (iz.eq.2) c1 = k_S21S31(2,ihelp)
1215 if (iz.gt.2) c1 = k_S21S31(3,ihelp)
1216 betatof_a(5) = c2/(cos(theta13)*(ds-c1))
1217
1218 C------- ToF Mask - S21 - S31
1219
1220 tofmask(ch21a(tof21_i),hb21a(tof21_i)) =
1221 $ tofmask(ch21a(tof21_i),hb21a(tof21_i)) + 1
1222 tofmask(ch21b(tof21_i),hb21b(tof21_i)) =
1223 $ tofmask(ch21b(tof21_i),hb21b(tof21_i)) + 1
1224
1225 tofmask(ch31a(tof31_i),hb31a(tof31_i)) =
1226 $ tofmask(ch31a(tof31_i),hb31a(tof31_i)) + 1
1227 tofmask(ch31b(tof31_i),hb31b(tof31_i)) =
1228 $ tofmask(ch31b(tof31_i),hb31b(tof31_i)) + 1
1229
1230 C-------
1231
1232 ENDIF
1233
1234 C S21 - S32
1235
1236 dist = ZTOF(3) - ZTOF(6)
1237 c2 = (2.*0.01*dist)/(3.E08*50.E-12 )
1238
1239
1240 IF ((tof21_i.GT.none_find).AND.(tof32_i.GT.none_find).AND.
1241 & (xtofpos(2).NE.101.).AND.(xtofpos(3).NE.101.)) THEN
1242 xhelp1 = tof21(1,tof21_i,itdc)+tof21(2,tof21_i,itdc)
1243 xhelp2 = tof32(1,tof32_i,itdc)+tof32(2,tof32_i,itdc)
1244 ds = xhelp1-xhelp2
1245 ihelp=(tof21_i-1)*3+tof32_i
1246 if (iz.le.1) c1 = k_S21S32(1,ihelp)
1247 if (iz.eq.2) c1 = k_S21S32(2,ihelp)
1248 if (iz.gt.2) c1 = k_S21S32(3,ihelp)
1249 betatof_a(6) = c2/(cos(theta13)*(ds-c1))
1250
1251 C------- ToF Mask - S21 - S32
1252
1253 tofmask(ch21a(tof21_i),hb21a(tof21_i)) =
1254 $ tofmask(ch21a(tof21_i),hb21a(tof21_i)) + 1
1255 tofmask(ch21b(tof21_i),hb21b(tof21_i)) =
1256 $ tofmask(ch21b(tof21_i),hb21b(tof21_i)) + 1
1257
1258 tofmask(ch32a(tof32_i),hb32a(tof32_i)) =
1259 $ tofmask(ch32a(tof32_i),hb32a(tof32_i)) + 1
1260 tofmask(ch32b(tof32_i),hb32b(tof32_i)) =
1261 $ tofmask(ch32b(tof32_i),hb32b(tof32_i)) + 1
1262
1263 C-------
1264
1265 ENDIF
1266
1267 C S22 - S31
1268
1269 dist = ZTOF(4) - ZTOF(5)
1270 c2 = (2.*0.01*dist)/(3.E08*50.E-12 )
1271
1272 IF ((tof22_i.GT.none_find).AND.(tof31_i.GT.none_find).AND.
1273 & (ytofpos(2).NE.101.).AND.(ytofpos(3).NE.101.)) THEN
1274 xhelp1 = tof22(1,tof22_i,itdc)+tof22(2,tof22_i,itdc)
1275 xhelp2 = tof31(1,tof31_i,itdc)+tof31(2,tof31_i,itdc)
1276 ds = xhelp1-xhelp2
1277 ihelp=(tof22_i-1)*3+tof31_i
1278 if (iz.le.1) c1 = k_S22S31(1,ihelp)
1279 if (iz.eq.2) c1 = k_S22S31(2,ihelp)
1280 if (iz.gt.2) c1 = k_S22S31(3,ihelp)
1281 betatof_a(7) = c2/(cos(theta13)*(ds-c1))
1282
1283 C------- ToF Mask - S22 - S31
1284
1285 tofmask(ch22a(tof22_i),hb22a(tof22_i)) =
1286 $ tofmask(ch22a(tof22_i),hb22a(tof22_i)) + 1
1287 tofmask(ch22b(tof22_i),hb22b(tof22_i)) =
1288 $ tofmask(ch22b(tof22_i),hb22b(tof22_i)) + 1
1289
1290 tofmask(ch31a(tof31_i),hb31a(tof31_i)) =
1291 $ tofmask(ch31a(tof31_i),hb31a(tof31_i)) + 1
1292 tofmask(ch31b(tof31_i),hb31b(tof31_i)) =
1293 $ tofmask(ch31b(tof31_i),hb31b(tof31_i)) + 1
1294
1295 C-------
1296
1297 ENDIF
1298
1299 C S22 - S32
1300
1301 dist = ZTOF(4) - ZTOF(6)
1302 c2 = (2.*0.01*dist)/(3.E08*50.E-12 )
1303
1304 IF ((tof22_i.GT.none_find).AND.(tof32_i.GT.none_find).AND.
1305 & (ytofpos(2).NE.101.).AND.(xtofpos(3).NE.101.)) THEN
1306 xhelp1 = tof22(1,tof22_i,itdc)+tof22(2,tof22_i,itdc)
1307 xhelp2 = tof32(1,tof32_i,itdc)+tof32(2,tof32_i,itdc)
1308 ds = xhelp1-xhelp2
1309 ihelp=(tof22_i-1)*3+tof32_i
1310 if (iz.le.1) c1 = k_S22S32(1,ihelp)
1311 if (iz.eq.2) c1 = k_S22S32(2,ihelp)
1312 if (iz.gt.2) c1 = k_S22S32(3,ihelp)
1313 betatof_a(8) = c2/(cos(theta13)*(ds-c1))
1314
1315 C------- ToF Mask - S22 - S32
1316
1317 tofmask(ch22a(tof22_i),hb22a(tof22_i)) =
1318 $ tofmask(ch22a(tof22_i),hb22a(tof22_i)) + 1
1319 tofmask(ch22b(tof22_i),hb22b(tof22_i)) =
1320 $ tofmask(ch22b(tof22_i),hb22b(tof22_i)) + 1
1321
1322 tofmask(ch32a(tof32_i),hb32a(tof32_i)) =
1323 $ tofmask(ch32a(tof32_i),hb32a(tof32_i)) + 1
1324 tofmask(ch32b(tof32_i),hb32b(tof32_i)) =
1325 $ tofmask(ch32b(tof32_i),hb32b(tof32_i)) + 1
1326
1327 C-------
1328
1329 ENDIF
1330
1331 C S11 - S21
1332
1333 dist = ZTOF(1) - ZTOF(3)
1334 c2 = (2.*0.01*dist)/(3.E08*50.E-12 )
1335
1336 IF ((tof11_i.GT.none_find).AND.(tof21_i.GT.none_find).AND.
1337 & (ytofpos(1).NE.101.).AND.(xtofpos(2).NE.101.)) THEN
1338 xhelp1 = tof11(1,tof11_i,itdc)+tof11(2,tof11_i,itdc)
1339 xhelp2 = tof21(1,tof21_i,itdc)+tof21(2,tof21_i,itdc)
1340 ds = xhelp1-xhelp2
1341 ihelp=(tof11_i-1)*2+tof21_i
1342 if (iz.le.1) c1 = k_S11S21(1,ihelp)
1343 if (iz.eq.2) c1 = k_S11S21(2,ihelp)
1344 if (iz.gt.2) c1 = k_S11S21(3,ihelp)
1345 betatof_a(9) = c2/(cos(theta13)*(ds-c1))
1346
1347 C------- ToF Mask - S11 - S21
1348
1349 tofmask(ch11a(tof11_i),hb11a(tof11_i)) =
1350 $ tofmask(ch11a(tof11_i),hb11a(tof11_i)) + 1
1351 tofmask(ch11b(tof11_i),hb11b(tof11_i)) =
1352 $ tofmask(ch11b(tof11_i),hb11b(tof11_i)) + 1
1353
1354 tofmask(ch21a(tof21_i),hb21a(tof21_i)) =
1355 $ tofmask(ch21a(tof21_i),hb21a(tof21_i)) + 1
1356 tofmask(ch21b(tof21_i),hb21b(tof21_i)) =
1357 $ tofmask(ch21b(tof21_i),hb21b(tof21_i)) + 1
1358
1359 C-------
1360
1361 ENDIF
1362
1363 C S11 - S22
1364
1365 dist = ZTOF(1) - ZTOF(4)
1366 c2 = (2.*0.01*dist)/(3.E08*50.E-12 )
1367
1368 IF ((tof11_i.GT.none_find).AND.(tof22_i.GT.none_find).AND.
1369 & (ytofpos(1).NE.101.).AND.(ytofpos(2).NE.101.)) THEN
1370 xhelp1 = tof11(1,tof11_i,itdc)+tof11(2,tof11_i,itdc)
1371 xhelp2 = tof22(1,tof22_i,itdc)+tof22(2,tof22_i,itdc)
1372 ds = xhelp1-xhelp2
1373 ihelp=(tof11_i-1)*2+tof22_i
1374 if (iz.le.1) c1 = k_S11S22(1,ihelp)
1375 if (iz.eq.2) c1 = k_S11S22(2,ihelp)
1376 if (iz.gt.2) c1 = k_S11S22(3,ihelp)
1377 betatof_a(10) = c2/(cos(theta13)*(ds-c1))
1378
1379 C------- ToF Mask - S11 - S22
1380
1381 tofmask(ch11a(tof11_i),hb11a(tof11_i)) =
1382 $ tofmask(ch11a(tof11_i),hb11a(tof11_i)) + 1
1383 tofmask(ch11b(tof11_i),hb11b(tof11_i)) =
1384 $ tofmask(ch11b(tof11_i),hb11b(tof11_i)) + 1
1385
1386 tofmask(ch22a(tof22_i),hb22a(tof22_i)) =
1387 $ tofmask(ch22a(tof22_i),hb22a(tof22_i)) + 1
1388 tofmask(ch22b(tof22_i),hb22b(tof22_i)) =
1389 $ tofmask(ch22b(tof22_i),hb22b(tof22_i)) + 1
1390
1391 C-------
1392
1393 ENDIF
1394
1395 C S12 - S21
1396
1397 dist = ZTOF(2) - ZTOF(3)
1398 c2 = (2.*0.01*dist)/(3.E08*50.E-12 )
1399
1400 IF ((tof12_i.GT.none_find).AND.(tof21_i.GT.none_find).AND.
1401 & (xtofpos(1).NE.101.).AND.(xtofpos(2).NE.101.)) THEN
1402 xhelp1 = tof12(1,tof12_i,itdc)+tof12(2,tof12_i,itdc)
1403 xhelp2 = tof21(1,tof21_i,itdc)+tof21(2,tof21_i,itdc)
1404 ds = xhelp1-xhelp2
1405 ihelp=(tof12_i-1)*2+tof21_i
1406 if (iz.le.1) c1 = k_S12S21(1,ihelp)
1407 if (iz.eq.2) c1 = k_S12S21(2,ihelp)
1408 if (iz.gt.2) c1 = k_S12S21(3,ihelp)
1409 betatof_a(11) = c2/(cos(theta13)*(ds-c1))
1410
1411 C------- ToF Mask - S12 - S21
1412
1413 tofmask(ch12a(tof12_i),hb12a(tof12_i)) =
1414 $ tofmask(ch12a(tof12_i),hb12a(tof12_i)) + 1
1415 tofmask(ch12b(tof12_i),hb12b(tof12_i)) =
1416 $ tofmask(ch12b(tof12_i),hb12b(tof12_i)) + 1
1417
1418 tofmask(ch21a(tof21_i),hb21a(tof21_i)) =
1419 $ tofmask(ch21a(tof21_i),hb21a(tof21_i)) + 1
1420 tofmask(ch21b(tof21_i),hb21b(tof21_i)) =
1421 $ tofmask(ch21b(tof21_i),hb21b(tof21_i)) + 1
1422
1423 C-------
1424
1425 ENDIF
1426
1427 C S12 - S22
1428
1429 dist = ZTOF(2) - ZTOF(4)
1430 c2 = (2.*0.01*dist)/(3.E08*50.E-12 )
1431
1432 IF ((tof12_i.GT.none_find).AND.(tof22_i.GT.none_find).AND.
1433 & (xtofpos(1).NE.101.).AND.(ytofpos(2).NE.101.)) THEN
1434 xhelp1 = tof12(1,tof12_i,itdc)+tof12(2,tof12_i,itdc)
1435 xhelp2 = tof22(1,tof22_i,itdc)+tof22(2,tof22_i,itdc)
1436 ds = xhelp1-xhelp2
1437 ihelp=(tof12_i-1)*2+tof22_i
1438 if (iz.le.1) c1 = k_S12S22(1,ihelp)
1439 if (iz.eq.2) c1 = k_S12S22(2,ihelp)
1440 if (iz.gt.2) c1 = k_S12S22(3,ihelp)
1441 betatof_a(12) = c2/(cos(theta13)*(ds-c1))
1442
1443 C------- ToF Mask - S12 - S22
1444
1445 tofmask(ch12a(tof12_i),hb12a(tof12_i)) =
1446 $ tofmask(ch12a(tof12_i),hb12a(tof12_i)) + 1
1447 tofmask(ch12b(tof12_i),hb12b(tof12_i)) =
1448 $ tofmask(ch12b(tof12_i),hb12b(tof12_i)) + 1
1449
1450 tofmask(ch22a(tof22_i),hb22a(tof22_i)) =
1451 $ tofmask(ch22a(tof22_i),hb22a(tof22_i)) + 1
1452 tofmask(ch22b(tof22_i),hb22b(tof22_i)) =
1453 $ tofmask(ch22b(tof22_i),hb22b(tof22_i)) + 1
1454
1455 C-------
1456
1457 ENDIF
1458
1459 C---------------------------------------------------------
1460 C
1461 C icount=0
1462 C sw=0.
1463 C sxw=0.
1464 C beta_mean=100.
1465 C
1466 C do i=1,12
1467 C if ((betatof_a(i).gt.-1.5).and.(betatof_a(i).lt.1.5)) then
1468 C icount= icount+1
1469 C if (i.le.4) w_i=1./(0.13**2.)
1470 C if ((i.ge.5).and.(i.le.8)) w_i=1./(0.16**2.)
1471 C if (i.ge.9) w_i=1./(0.25**2.) ! to be checked
1472 C sxw=sxw + betatof_a(i)*w_i
1473 C sw =sw + w_i
1474 C endif
1475 C enddo
1476 C
1477 C if (icount.gt.0) beta_mean=sxw/sw
1478 C betatof_a(13) = beta_mean
1479 C
1480
1481 C-------- New mean beta calculation -----------------------
1482
1483 do i=1,12
1484 btemp(i) = betatof_a(i)
1485 enddo
1486
1487 betatof_a(13)=newbeta(1,btemp,hitvec,10.,10.,20.)
1488
1489 C--------------------------------------------------------------
1490 C write(*,*) betatof_a
1491 c write(*,*) xtofpos
1492 c write(*,*) ytofpos
1493 c write(*,*)'tofl2com beta', betatof_a
1494 C write(*,*) adcflagtof
1495 c write(*,*) 'tofl2com'
1496 c write(*,*) xtofpos
1497 c write(*,*) ytofpos
1498 c write(*,*) xtr_tof
1499 c write(*,*) ytr_tof
1500
1501 c 100 continue
1502 continue
1503
1504 C
1505 RETURN
1506 END
1507
1508
1509 C------------------------------------------------------------------
1510 C------------------------------------------------------------------
1511
1512 function atten(is,ilay,ipad,x)
1513 include 'input_tof.txt'
1514 real atten
1515 real x
1516 real xmin,xmax
1517 integer ilay,ipad
1518
1519 * S11 8 paddles 33.0 x 5.1 cm
1520 * S12 6 paddles 40.8 x 5.5 cm
1521 * S21 2 paddles 18.0 x 7.5 cm
1522 * S22 2 paddles 15.0 x 9.0 cm
1523 * S31 3 paddles 15.0 x 6.0 cm
1524 * S32 3 paddles 18.0 x 5.0 cm
1525
1526
1527 c if (ilay.eq.11) write(*,*) 'start ',ipad,is,adcx11(is,ipad,1),
1528 c & adcx11(is,ipad,2),adcx11(is,ipad,3),adcx11(is,ipad,4)
1529 c if (ilay.eq.12) write(*,*) 'start ',ipad,is,adcx12(is,ipad,1),
1530 c & adcx12(is,ipad,2),adcx12(is,ipad,3),adcx12(is,ipad,4)
1531
1532
1533 if (ilay.eq.11) xmin=-33.0/2.
1534 if (ilay.eq.11) xmax= 33.0/2.
1535 if (ilay.eq.12) xmin=-40.8/2.
1536 if (ilay.eq.12) xmax= 40.8/2.
1537
1538 if (ilay.eq.21) xmin=-18.0/2.
1539 if (ilay.eq.21) xmax= 18.0/2.
1540 if (ilay.eq.22) xmin=-15.0/2.
1541 if (ilay.eq.22) xmax= 15.0/2.
1542
1543 if (ilay.eq.31) xmin=-15.0/2.
1544 if (ilay.eq.31) xmax= 15.0/2.
1545 if (ilay.eq.32) xmin=-18.0/2.
1546 if (ilay.eq.32) xmax= 18.0/2.
1547
1548 if (x .lt. xmin) x=xmin
1549 if (x .gt. xmax) x=xmax
1550
1551
1552 if (ilay.eq.11) atten=
1553 & adcx11(is,ipad,1)*exp(x*adcx11(is,ipad,2))
1554 & + adcx11(is,ipad,3)*exp(x*adcx11(is,ipad,4))
1555
1556 if (ilay.eq.12) atten=
1557 & adcx12(is,ipad,1)*exp(x*adcx12(is,ipad,2))
1558 & + adcx12(is,ipad,3)*exp(x*adcx12(is,ipad,4))
1559
1560 if (ilay.eq.21) atten=
1561 & adcx21(is,ipad,1)*exp(x*adcx21(is,ipad,2))
1562 & + adcx21(is,ipad,3)*exp(x*adcx21(is,ipad,4))
1563
1564 if (ilay.eq.22) atten=
1565 & adcx22(is,ipad,1)*exp(x*adcx22(is,ipad,2))
1566 & + adcx22(is,ipad,3)*exp(x*adcx22(is,ipad,4))
1567
1568 if (ilay.eq.31) atten=
1569 & adcx31(is,ipad,1)*exp(x*adcx31(is,ipad,2))
1570 & + adcx31(is,ipad,3)*exp(x*adcx31(is,ipad,4))
1571
1572 if (ilay.eq.32) atten=
1573 & adcx32(is,ipad,1)*exp(x*adcx32(is,ipad,2))
1574 & + adcx32(is,ipad,3)*exp(x*adcx32(is,ipad,4))
1575
1576 if (atten.gt.10000) atten=10000.
1577
1578 end
1579
1580 C------------------------------------------------------------------
1581 C------------------------------------------------------------------
1582
1583 function pc_adc(ix)
1584 include 'input_tof.txt'
1585 real pc_adc
1586 integer ix
1587
1588 pc_adc=28.0407 + 0.628929*ix
1589 & - 5.80901e-05*ix*ix + 3.14092e-08*ix*ix*ix
1590 c write(*,*) ix,pc_adc
1591 end
1592
1593 C------------------------------------------------------------------
1594 C------------------------------------------------------------------
1595
1596 function check_charge(theta,hitvec)
1597
1598 include 'input_tof.txt'
1599 include 'tofcomm.txt'
1600
1601 real check_charge
1602 integer hitvec(6)
1603 REAL CHARGE, theta
1604
1605 C upper and lower limits for the helium selection
1606 REAL A_l(24),A_h(24)
1607 DATA A_l /200,190,300,210,220,200,210,60,60,120,220,
1608 & 120,160,50,300,200,120,250,350,300,350,250,280,300/
1609 DATA A_h /550,490,800,600,650,600,600,260,200,380,
1610 & 620,380,550,200,850,560,400,750,900,800,880,800,750,800/
1611
1612 C The k1 constants for the beta calculation, only for S1-S3
1613 C k2 constant is taken to be the standard 2D/c
1614 REAL k1(84)
1615 DATA k1 /50,59.3296,28.4328,-26.0818,5.91253,-19.588,
1616 & -9.26316,24.7544,2.32465,-50.5058,-15.3195,-39.1443,
1617 & -91.2546,-58.6243,-84.5641,-63.1516,-32.2091,-58.3358,
1618 & 13.8084,45.5322,33.2416,-11.5313,51.3271,75,-14.1141,
1619 & 42.8466,15.1794,-63.6672,-6.07739,-32.164,-41.771,10.5274,
1620 & -9.46096,-81.7404,-28.783,-52.7167,-127.394,-69.6166,
1621 & -93.4655,-98.9543,-42.863,-67.8244,-19.3238,31.1221,8.7319,
1622 & -43.1627,5.55573,-14.4078,-83.4466,-47.4647,-77.8379,
1623 & -108.222,-75.986,-101.297,-96.0205,-63.1881,-90.1372,
1624 & -22.7347,8.31409,-19.6912,-7.49008,23.6979,-1.66677,
1625 & 1.81556,34.4668,6.23693,-100,-59.5861,-90.9159,-141.639,
1626 & -89.2521,-112.881,-130.199,-77.0357,-98.4632,-60.2086,
1627 & -4.82097,-29.3705,-43.6469,10.5884,-9.31304,-35.3329,
1628 & 25.2514,25.6/
1629
1630
1631
1632 REAL zin(6)
1633 DATA zin /53.74, 53.04, 23.94, 23.44, -23.49, -24.34/
1634
1635 REAL c1,c2,xhelp,xhelp1,xhelp2,ds,dist,F
1636 REAL sw,sxw,beta_mean_tof,w_i
1637 INTEGER ihelp
1638 INTEGER ipmt(4)
1639 REAL time(4),beta1(4)
1640
1641 REAL adca(48), tdca(48)
1642
1643 REAL a1,a2
1644 INTEGER jj
1645
1646 c get rid of warnings EMILIANO
1647 i = 0
1648 slope = 0
1649 offset = 0
1650 none_find = 0
1651 none_ev = 0
1652 adc_ev = 0
1653 tdc_ev = 0
1654 iadc = 0
1655 itdc = 0
1656 right = 0
1657 left = 0
1658 tof12_y(1) = tof12_y(1)
1659 tof11_x(1) = tof11_x(1)
1660 tof21_y(1) = tof21_y(1)
1661 tof22_x(1) = tof22_x(1)
1662 tof32_y(1) = tof32_y(1)
1663 tof31_x(1) = tof31_x(1)
1664 c get rid of warnings
1665
1666 C-----------------------------------------------------------
1667 C--- get data
1668 C-----------------------------------------------------------
1669
1670 do j=1,8
1671 ih = 1 + 0 +((j-1)*2)
1672 adca(ih) = adc(ch11a(j),hb11a(j))
1673 adca(ih+1) = adc(ch11b(j),hb11b(j))
1674 tdca(ih) = tdc(ch11a(j),hb11a(j))
1675 tdca(ih+1) = tdc(ch11b(j),hb11b(j))
1676 enddo
1677
1678 do j=1,6
1679 ih = 1 + 16+((j-1)*2)
1680 adca(ih) = adc(ch12a(j),hb12a(j))
1681 adca(ih+1) = adc(ch12b(j),hb12b(j))
1682 tdca(ih) = tdc(ch12a(j),hb12a(j))
1683 tdca(ih+1) = tdc(ch12b(j),hb12b(j))
1684 enddo
1685
1686 do j=1,2
1687 ih = 1 + 28+((j-1)*2)
1688 adca(ih) = adc(ch21a(j),hb21a(j))
1689 adca(ih+1) = adc(ch21b(j),hb21b(j))
1690 tdca(ih) = tdc(ch21a(j),hb21a(j))
1691 tdca(ih+1) = tdc(ch21b(j),hb21b(j))
1692 enddo
1693
1694 do j=1,2
1695 ih = 1 + 32+((j-1)*2)
1696 adca(ih) = adc(ch22a(j),hb22a(j))
1697 adca(ih+1) = adc(ch22b(j),hb22b(j))
1698 tdca(ih) = tdc(ch22a(j),hb22a(j))
1699 tdca(ih+1) = tdc(ch22b(j),hb22b(j))
1700 enddo
1701
1702 do j=1,3
1703 ih = 1 + 36+((j-1)*2)
1704 adca(ih) = adc(ch31a(j),hb31a(j))
1705 adca(ih+1) = adc(ch31b(j),hb31b(j))
1706 tdca(ih) = tdc(ch31a(j),hb31a(j))
1707 tdca(ih+1) = tdc(ch31b(j),hb31b(j))
1708 enddo
1709
1710 do j=1,3
1711 ih = 1 + 42+((j-1)*2)
1712 adca(ih) = adc(ch32a(j),hb32a(j))
1713 adca(ih+1) = adc(ch32b(j),hb32b(j))
1714 tdca(ih) = tdc(ch32a(j),hb32a(j))
1715 tdca(ih+1) = tdc(ch32b(j),hb32b(j))
1716 enddo
1717
1718
1719 c write(*,*) adca
1720 c write(*,*) tdca
1721
1722
1723 C============ calculate beta and select charge > Z=1 ===============
1724
1725 ICHARGE=1
1726
1727 C find hitted paddle by looking for ADC values on both sides
1728 C since we looking for Z>1 this gives decent results
1729
1730 tof11_i = hitvec(1)-1
1731 tof12_i = hitvec(2)-1
1732 tof21_i = hitvec(3)-1
1733 tof22_i = hitvec(4)-1
1734 tof31_i = hitvec(5)-1
1735 tof32_i = hitvec(6)-1
1736
1737 c write(*,*) ' in charge check'
1738 c write(*,*) theta,tof11_i,tof12_i,tof21_i,tof22_i,tof31_i,tof32_i
1739
1740 C----------------------------------------------------------------
1741
1742 beta_help=100.
1743 beta_mean_tof=100.
1744
1745 do jj=1,4
1746 beta1(jj) = 100.
1747 enddo
1748
1749 C----------------------------------------------------------------
1750 C--------- S1 - S3 ---------------------------------------------
1751 C----------------------------------------------------------------
1752
1753 C--------- S11 - S31 -------------------------------------------
1754
1755 if ((tof11_i.gt.-1).and.(tof31_i.gt.-1)) then
1756
1757 dist = zin(1) - zin(5)
1758 c2 = (2.*0.01*dist)/(3.E08*50.E-12)
1759 F = 1./cos(theta)
1760
1761 ipmt(1) = (tof11_i)*2+1
1762 ipmt(2) = (tof11_i)*2+2
1763 ipmt(3) = 36+(tof31_i)*2+1
1764 ipmt(4) = 36+(tof31_i)*2+2
1765
1766 c write(*,*) ipmt
1767
1768 do jj=1,4
1769 time(jj) = tdca(ipmt(jj))
1770 enddo
1771
1772 c write(*,*) time
1773
1774 if ((time(1).lt.4095).and.(time(2).lt.4095).and.
1775 & (time(3).lt.4095).and.(time(4).lt.4095)) then
1776 xhelp1 = time(1) + time(2)
1777 xhelp2 = time(3) + time(4)
1778 ds = xhelp1-xhelp2
1779 ihelp=0+(tof11_i)*3+tof31_i
1780 c1 = k1(ihelp+1)
1781 beta1(1) = c2*F/(ds-c1);
1782 endif
1783 c write(*,*) beta1(1)
1784 endif ! tof_....
1785
1786
1787 C--------- S11 - S32 -------------------------------------------
1788
1789 if ((tof11_i.gt.-1).and.(tof32_i.gt.-1)) then
1790
1791 dist = zin(1) - zin(6)
1792 c2 = (2.*0.01*dist)/(3.E08*50.E-12)
1793 F = 1./cos(theta)
1794
1795 ipmt(1) = (tof11_i)*2+1
1796 ipmt(2) = (tof11_i)*2+2
1797 ipmt(3) = 42+(tof32_i)*2+1
1798 ipmt(4) = 42+(tof32_i)*2+2
1799
1800 do jj=1,4
1801 time(jj) = tdca(ipmt(jj))
1802 enddo
1803
1804 if ((time(1).lt.4095).and.(time(2).lt.4095).and.
1805 & (time(3).lt.4095).and.(time(4).lt.4095)) then
1806 xhelp1 = time(1) + time(2)
1807 xhelp2 = time(3) + time(4)
1808 ds = xhelp1-xhelp2
1809 ihelp=24+(tof11_i)*3+tof32_i
1810 c1 = k1(ihelp+1)
1811 beta1(2) = c2*F/(ds-c1);
1812 endif
1813 endif ! tof_....
1814
1815
1816 C--------- S12 - S31 -------------------------------------------
1817
1818 if ((tof12_i.gt.-1).and.(tof31_i.gt.-1)) then
1819
1820 dist = zin(2) - zin(5)
1821 c2 = (2.*0.01*dist)/(3.E08*50.E-12)
1822 F = 1./cos(theta)
1823
1824 ipmt(1) = 16+(tof12_i)*2+1
1825 ipmt(2) = 16+(tof12_i)*2+2
1826 ipmt(3) = 36+(tof31_i)*2+1
1827 ipmt(4) = 36+(tof31_i)*2+2
1828
1829 do jj=1,4
1830 time(jj) = tdca(ipmt(jj))
1831 enddo
1832
1833 if ((time(1).lt.4095).and.(time(2).lt.4095).and.
1834 & (time(3).lt.4095).and.(time(4).lt.4095)) then
1835 xhelp1 = time(1) + time(2)
1836 xhelp2 = time(3) + time(4)
1837 ds = xhelp1-xhelp2
1838 ihelp=48+(tof12_i)*3+tof31_i
1839 c1 = k1(ihelp+1)
1840 beta1(3) = c2*F/(ds-c1);
1841 endif
1842 endif ! tof_....
1843
1844
1845 C--------- S12 - S32 -------------------------------------------
1846
1847 if ((tof12_i.gt.-1).and.(tof32_i.gt.-1)) then
1848
1849 dist = zin(2) - zin(6)
1850 c2 = (2.*0.01*dist)/(3.E08*50.E-12)
1851 F = 1./cos(theta)
1852
1853 ipmt(1) = 16+(tof12_i)*2+1
1854 ipmt(2) = 16+(tof12_i)*2+2
1855 ipmt(3) = 42+(tof32_i)*2+1
1856 ipmt(4) = 42+(tof32_i)*2+2
1857
1858 do jj=1,4
1859 time(jj) = tdca(ipmt(jj))
1860 enddo
1861
1862 if ((time(1).lt.4095).and.(time(2).lt.4095).and.
1863 & (time(3).lt.4095).and.(time(4).lt.4095)) then
1864 xhelp1 = time(1) + time(2)
1865 xhelp2 = time(3) + time(4)
1866 ds = xhelp1-xhelp2
1867 ihelp=56+(tof12_i)*3+tof32_i
1868 c1 = k1(ihelp+1)
1869 beta1(4) = c2*F/(ds-c1);
1870 endif
1871
1872 endif ! tof_....
1873
1874 c write(*,*) beta1
1875
1876 C---- calculate beta mean, only downward going particles are interesting ----
1877
1878 sw=0.
1879 sxw=0.
1880 beta_mean_tof=100.
1881
1882 do jj=1,4
1883 if ((beta1(jj).gt.0.1).and.(beta1(jj).lt.2.0)) then
1884 w_i=1./(0.13*0.13)
1885 sxw=sxw + beta1(jj)*w_i
1886 sw =sw + w_i ;
1887 endif
1888 enddo
1889
1890 if (sw.gt.0) beta_mean_tof=sxw/sw;
1891
1892 C write(*,*) 'beta_mean_tof ',beta_mean_tof
1893
1894 beta_help = beta_mean_tof ! pow(beta_mean_tof,1.0) gave best results
1895
1896 CCCCC endif ! if tof11_i > -1 && ...... beta calculation
1897
1898 C----------------------- Select charge --------------------------
1899
1900 charge=0
1901
1902 if ((beta_mean_tof.gt.0.2).and.(beta_mean_tof.lt.2.0)) then
1903
1904 icount1=0
1905 icount2=0
1906 icount3=0
1907
1908 do jj=0,23
1909 a1 = adca(2*jj+1)
1910 a2 = adca(2*jj+2)
1911 if ((a1.lt.4095).and.(a2.lt.4095)) then
1912 a1 = adca(2*jj+1)*cos(theta)
1913 a2 = adca(2*jj+2)*cos(theta)
1914 xhelp = 100000.
1915 xhelp1 = 100000.
1916 xhelp = sqrt(a1*a2) ! geometric mean
1917 xhelp1 = beta_help*xhelp
1918 C if geometric mean multiplied by beta_help is inside/outside helium
1919 C limits, increase counter
1920 if (xhelp1.lt.A_l(jj+1)) icount1=icount1+1
1921 if ((xhelp1.gt.A_l(jj+1)).and.(xhelp1.lt.A_h(jj+1)))
1922 & icount2=icount2+1
1923 if (xhelp1.gt.A_h(jj+1)) icount3=icount3+1
1924 endif
1925 enddo
1926
1927
1928 C if more than three paddles see the same...
1929
1930 if (icount1 .gt. 3) charge=1
1931 if (icount2 .gt. 3) charge=2
1932 if (icount3 .gt. 3) charge=3
1933
1934 endif ! 0.2<beta<2.0
1935
1936 C no beta found? Sum up geometric means of paddles and derive the mean...
1937
1938 if (beta_mean_tof.eq.100.) then
1939
1940 xhelp = 0.
1941 icount = 0
1942
1943 if (tof11_i.gt.-1) then
1944 jj=tof11_i
1945 a1 = adca(0+2*jj+1)
1946 a2 = adca(0+2*jj+2)
1947 if ((a1.lt.4095).and.(a2.lt.4095)) then
1948 a1 = a1*cos(theta)
1949 a2 = a2*cos(theta)
1950 xhelp = xhelp + sqrt(a1*a2)
1951 icount=icount+1
1952 endif
1953 endif
1954
1955 if (tof12_i.gt.-1) then
1956 jj=tof12_i
1957 a1 = adca(16+2*jj+1)
1958 a2 = adca(16+2*jj+2)
1959 if ((a1.lt.4095).and.(a2.lt.4095)) then
1960 a1 = a1*cos(theta)
1961 a2 = a2*cos(theta)
1962 xhelp = xhelp + sqrt(a1*a2)
1963 icount=icount+1
1964 endif
1965 endif
1966
1967 if (tof21_i.gt.-1) then
1968 jj=tof21_i
1969 a1 = adca(28+2*jj+1)
1970 a2 = adca(28+2*jj+2)
1971 if ((a1.lt.4095).and.(a2.lt.4095)) then
1972 a1 = a1*cos(theta)
1973 a2 = a2*cos(theta)
1974 xhelp = xhelp + sqrt(a1*a2)
1975 icount=icount+1
1976 endif
1977 endif
1978
1979 if (tof22_i.gt.-1) then
1980 jj=tof22_i
1981 a1 = adca(32+2*jj+1)
1982 a2 = adca(32+2*jj+2)
1983 if ((a1.lt.4095).and.(a2.lt.4095)) then
1984 a1 = a1*cos(theta)
1985 a2 = a2*cos(theta)
1986 xhelp = xhelp + sqrt(a1*a2)
1987 icount=icount+1
1988 endif
1989 endif
1990
1991 if (tof31_i.gt.-1) then
1992 jj=tof31_i
1993 a1 = adca(36+2*jj+1)
1994 a2 = adca(36+2*jj+2)
1995 if ((a1.lt.4095).and.(a2.lt.4095)) then
1996 a1 = a1*cos(theta)
1997 a2 = a2*cos(theta)
1998 xhelp = xhelp + sqrt(a1*a2)
1999 icount=icount+1
2000 endif
2001 endif
2002
2003 if (tof32_i.gt.-1) then
2004 jj=tof32_i
2005 a1 = adca(42+2*jj+1)
2006 a2 = adca(42+2*jj+2)
2007 if ((a1.lt.4095).and.(a2.lt.4095)) then
2008 a1 = a1*cos(theta)
2009 a2 = a2*cos(theta)
2010 xhelp = xhelp + sqrt(a1*a2)
2011 icount=icount+1
2012 endif
2013 endif
2014
2015
2016 if (icount.gt.0) xhelp=xhelp/icount
2017 if ((icount.gt.2).and.(xhelp.gt.1500.)) charge=3
2018
2019 endif ! beta_mean_tof.eq.100.
2020
2021 C write(*,*) 'in function charge: ',beta_mean_tof,charge
2022
2023 check_charge = charge
2024
2025
2026 END
2027
2028 C****************************************************************************
2029 C****************************************************************************
2030 C****************************************************************************
2031
2032 function newbeta(iflag,b,hitvec,resmax,qualitycut,chi2cut)
2033
2034 include 'input_tof.txt'
2035 include 'output_tof.txt'
2036 include 'tofcomm.txt'
2037
2038 REAL newbeta
2039 REAL resmax,qualitycut,chi2cut
2040 REAL w_i(12),w_il(6),quality,res,betachi,beta_mean_inv
2041 REAL sw,sxw,b(12),beta_mean,chi2,xhelp
2042 c REAL tdcfl(4,12)
2043 INTEGER tdcfl(4,12) !EM GCC4.7
2044
2045 INTEGER iflag,icount,hitvec(6)
2046
2047 INTEGER itop(12),ibot(12)
2048 DATA itop /1,1,2,2,3,3,4,4,1,1,2,2/
2049 DATA ibot /5,6,5,6,5,6,5,6,3,4,3,4/
2050
2051
2052 c get rid of warnings EMILIANO
2053 slope = 0
2054 offset = 0
2055 none_find = 0
2056 none_ev = 0
2057 adc_ev = 0
2058 tdc_ev = 0
2059 iadc = 0
2060 itdc = 0
2061 right = 0
2062 left = 0
2063 tof12_y(1) = tof12_y(1)
2064 tof11_x(1) = tof11_x(1)
2065 tof21_y(1) = tof21_y(1)
2066 tof22_x(1) = tof22_x(1)
2067 tof32_y(1) = tof32_y(1)
2068 tof31_x(1) = tof31_x(1)
2069 c get rid of warnings
2070
2071 C====================================================================
2072
2073 tof11_i = hitvec(1)
2074 tof12_i = hitvec(2)
2075 tof21_i = hitvec(3)
2076 tof22_i = hitvec(4)
2077 tof31_i = hitvec(5)
2078 tof32_i = hitvec(6)
2079
2080 if (iflag.eq.1) then ! call from tofl2com
2081 do i=1,4
2082 do j=1,12
2083 tdcfl(i,j) = tdcflagtof(i,j)
2084 enddo
2085 enddo
2086 endif
2087
2088 if (iflag.eq.2) then ! call from toftrk
2089 do i=1,4
2090 do j=1,12
2091 tdcfl(i,j) = tdcflag(i,j)
2092 enddo
2093 enddo
2094 endif
2095
2096
2097 C--- Find out ToF layers with artificial TDC values -------------
2098
2099 do jj=1,6
2100 w_il(jj) = 1000.
2101 enddo
2102
2103
2104 if (tof11_i.gt.0) then
2105 if ((tofmask(ch11a(tof11_i),hb11a(tof11_i)).gt.0).or.
2106 & (tofmask(ch11b(tof11_i),hb11b(tof11_i)).gt.0)) then
2107 w_il(1)=0
2108 i1=tdcfl(ch11a(tof11_i),hb11a(tof11_i))
2109 i2=tdcfl(ch11b(tof11_i),hb11b(tof11_i))
2110 if ((i1.eq.1).or.(i2.eq.1)) w_il(1) = 1 ! tdcflag
2111 endif
2112 endif
2113
2114 if (tof12_i.gt.0) then
2115 if ((tofmask(ch12a(tof12_i),hb12a(tof12_i)).gt.0).or.
2116 & (tofmask(ch12b(tof12_i),hb12b(tof12_i)).gt.0)) then
2117 w_il(2)=0
2118 i1=tdcfl(ch12a(tof12_i),hb12a(tof12_i))
2119 i2=tdcfl(ch12b(tof12_i),hb12b(tof12_i))
2120 if ((i1.eq.1).or.(i2.eq.1)) w_il(2) = 1 ! tdcflag
2121 endif
2122 endif
2123
2124 if (tof21_i.gt.0) then
2125 if ((tofmask(ch21a(tof21_i),hb21a(tof21_i)).gt.0).or.
2126 & (tofmask(ch21b(tof21_i),hb21b(tof21_i)).gt.0)) then
2127 w_il(3)=0
2128 i1=tdcfl(ch21a(tof21_i),hb21a(tof21_i))
2129 i2=tdcfl(ch21b(tof21_i),hb21b(tof21_i))
2130 if ((i1.eq.1).or.(i2.eq.1)) w_il(3) = 1 ! tdcflag
2131 endif
2132 endif
2133
2134 if (tof22_i.gt.0) then
2135 if ((tofmask(ch22a(tof22_i),hb22a(tof22_i)).gt.0).or.
2136 & (tofmask(ch22b(tof22_i),hb22b(tof22_i)).gt.0)) then
2137 w_il(4)=0
2138 i1=tdcfl(ch22a(tof22_i),hb22a(tof22_i))
2139 i2=tdcfl(ch22b(tof22_i),hb22b(tof22_i))
2140 if ((i1.eq.1).or.(i2.eq.1)) w_il(4) = 1 ! tdcflag
2141 endif
2142 endif
2143
2144 if (tof31_i.gt.0) then
2145 if ((tofmask(ch31a(tof31_i),hb11a(tof31_i)).gt.0).or.
2146 & (tofmask(ch31b(tof31_i),hb31b(tof31_i)).gt.0)) then
2147 w_il(5)=0
2148 i1=tdcfl(ch31a(tof31_i),hb31a(tof31_i))
2149 i2=tdcfl(ch31b(tof31_i),hb31b(tof31_i))
2150 if ((i1.eq.1).or.(i2.eq.1)) w_il(5) = 1 ! tdcflag
2151 endif
2152 endif
2153
2154 if (tof32_i.gt.0) then
2155 if ((tofmask(ch32a(tof32_i),hb32a(tof32_i)).gt.0).or.
2156 & (tofmask(ch32b(tof32_i),hb32b(tof32_i)).gt.0)) then
2157 w_il(6)=0
2158 i1=tdcfl(ch32a(tof32_i),hb32a(tof32_i))
2159 i2=tdcfl(ch32b(tof32_i),hb32b(tof32_i))
2160 if ((i1.eq.1).or.(i2.eq.1)) w_il(6) = 1 ! tdcflag
2161 endif
2162 endif
2163
2164 C------------------------------------------------------------------------
2165 C--- Set weights for the 12 measurements using information for top and bottom:
2166 C--- if no measurements: weight = set to very high value=> not used
2167 C--- top or bottom artificial: weight*sqrt(2)
2168 C--- top and bottom artificial: weight*sqrt(2)*sqrt(2)
2169
2170 DO jj=1,12
2171 if (jj.le.4) xhelp = 0.11 ! S1-S3
2172 if ((jj.gt.4).and.(jj.le.8)) xhelp = 0.18 ! S2-S3
2173 if (jj.gt.8) xhelp = 0.28 ! S1-S2
2174 if ((w_il(itop(jj)).eq.1000.).and.(w_il(ibot(jj)).eq.1000.))
2175 & xhelp = 1.E09
2176 if ((w_il(itop(jj)).eq.1).or.(w_il(ibot(jj)).eq.1.))
2177 & xhelp = xhelp*1.414
2178 if ((w_il(itop(jj)).eq.1).and.(w_il(ibot(jj)).eq.1.))
2179 & xhelp = xhelp*2.
2180 w_i(jj) = 1./xhelp
2181 ENDDO
2182
2183 C========================================================================
2184 C--- Calculate mean beta for the first time -----------------------------
2185 C--- We are using "1/beta" since its error is gaussian ------------------
2186
2187 icount=0
2188 sw=0.
2189 sxw=0.
2190 beta_mean=100.
2191
2192 DO jj=1,12
2193 IF ((abs(1./b(jj)).gt.0.1).and.(abs(1./b(jj)).lt.15.)) THEN
2194 icount = icount+1
2195 sxw = sxw + (1./b(jj))*w_i(jj)*w_i(jj)
2196 sw = sw + w_i(jj)*w_i(jj)
2197 ENDIF
2198 ENDDO
2199
2200 if (icount.gt.0) beta_mean=1./(sxw/sw)
2201 beta_mean_inv = 1./beta_mean
2202
2203
2204 C--- Calculate beta for the second time, use residuals of the single
2205 C--- measurements to get a chi2 value
2206
2207 icount = 0
2208 sw = 0.
2209 sxw = 0.
2210 betachi = 100.
2211 chi2 = 0.
2212 quality = 0.
2213
2214 DO jj=1,12
2215 IF ((abs(1./b(jj)).gt.0.1).and.(abs(1./b(jj)).lt.15.)
2216 & .and.(w_i(jj).GT.0.01)) THEN
2217 res = beta_mean_inv - (1./b(jj)) ;
2218 if (abs(res*w_i(jj)).lt.resmax) THEN
2219 chi2 = chi2 + (res*w_i(jj))**2.
2220 icount = icount+1
2221 sxw = sxw + (1./b(jj))*w_i(jj)*w_i(jj)
2222 sw = sw + w_i(jj)*w_i(jj)
2223 ENDIF
2224 ENDIF
2225 ENDDO
2226
2227 c quality = sw
2228 quality = sqrt(sw)
2229
2230 if (icount.eq.0) chi2 = 1000.
2231 if (icount.gt.0) chi2 = chi2/(icount)
2232
2233 if (icount.gt.0) betachi=1./(sxw/sw);
2234
2235 beta_mean=100.
2236 if ((chi2.lt.chi2cut).and.(quality.gt.qualitycut))
2237 & beta_mean = betachi
2238 newbeta = beta_mean
2239
2240 END
2241
2242 C****************************************************************************
2243

  ViewVC Help
Powered by ViewVC 1.1.23