1 |
************************************************************************* |
2 |
* |
3 |
* functions.f |
4 |
* |
5 |
* - !??? |
6 |
* |
7 |
* needs: |
8 |
* - !??? |
9 |
* |
10 |
* output variables: |
11 |
* - !??? |
12 |
* |
13 |
* to be called inside !??? |
14 |
* |
15 |
* |
16 |
* MODIFIED in order to have in input a |
17 |
* REAL-defined strip number instead of INTEGER |
18 |
* |
19 |
************************************************************************* |
20 |
|
21 |
|
22 |
function pitch(view) !it gives the strip pitch, knowing the view number |
23 |
|
24 |
real pitch |
25 |
integer view |
26 |
|
27 |
include 'commontracker.f' |
28 |
|
29 |
if(mod(view,2).eq.0) then !X |
30 |
pitch=pitchX |
31 |
else !Y |
32 |
pitch=pitchY |
33 |
endif |
34 |
|
35 |
end |
36 |
|
37 |
|
38 |
|
39 |
c------------------------------------------------------------------------ |
40 |
|
41 |
|
42 |
|
43 |
function npl(view) !it gives the plane number, knowing the view number. |
44 |
! plane 1 = views 11+12, calorimeter side |
45 |
! ... |
46 |
! plane 6 = views 1+2, TRD side |
47 |
integer npl,view |
48 |
|
49 |
npl=7-(INT((view-1)/2)+1) |
50 |
|
51 |
end |
52 |
|
53 |
|
54 |
|
55 |
c------------------------------------------------------------------------ |
56 |
|
57 |
|
58 |
|
59 |
function nld(istrip,view) |
60 |
* it gives the number of the ladder, knowing the |
61 |
* strip number (1..3072) and the view number. |
62 |
* the first strip belongs to ladder 1 |
63 |
|
64 |
integer istrip,view,nld |
65 |
|
66 |
include 'commontracker.f' |
67 |
|
68 |
|
69 |
nld=INT((istrip-1)/nstrips_ladder)+1 |
70 |
|
71 |
end |
72 |
|
73 |
|
74 |
c------------------------------------------------------------------------ |
75 |
|
76 |
|
77 |
function nviewx(iplane) !it gives the view number of a X plane |
78 |
|
79 |
integer nviewx,iplane |
80 |
|
81 |
nviewx=2*(7-iplane) |
82 |
|
83 |
end |
84 |
|
85 |
|
86 |
c------------------------------------------------------------------------ |
87 |
|
88 |
function nviewy(iplane) !it gives the view number of a Y plane |
89 |
|
90 |
integer nviewy,iplane |
91 |
|
92 |
nviewy=2*(7-iplane)-1 |
93 |
|
94 |
end |
95 |
|
96 |
c------------------------------------------------------------------------ |
97 |
|
98 |
|
99 |
|
100 |
|
101 |
function nvk(istrip) |
102 |
|
103 |
* it gives the number of the VA1, knowing the strip |
104 |
* number (1..3072). |
105 |
* the first strip belongs to VA1 1 |
106 |
integer istrip,nvk |
107 |
|
108 |
include 'commontracker.f' |
109 |
|
110 |
nvk=INT((istrip-1)/nstrips_va1)+1 |
111 |
|
112 |
end |
113 |
|
114 |
|
115 |
|
116 |
c------------------------------------------------------------------------ |
117 |
|
118 |
|
119 |
|
120 |
function nst(istrip) |
121 |
|
122 |
* it gives the VA1 strip, knowing the strip number |
123 |
* (1..3072). |
124 |
* the first strip belongs to VA1 1 |
125 |
|
126 |
integer istrip,nst |
127 |
|
128 |
include 'commontracker.f' |
129 |
|
130 |
nst=INT(mod((istrip-1),nstrips_va1))+1 |
131 |
|
132 |
|
133 |
end |
134 |
|
135 |
|
136 |
|
137 |
function coordsi(istrip,view) |
138 |
* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * |
139 |
* it gives the strip coordinate in micrometers, |
140 |
* knowing the strip number (1..3072) and the view |
141 |
* number. the origin of the coordinate is on the |
142 |
* centre of the sensor the strip belongs to. |
143 |
* the axes directions are the same as in the PAMELA |
144 |
* reference frame (i.e.: the 11th view coordinate |
145 |
* direction has to be inverted here) |
146 |
* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * |
147 |
|
148 |
c integer is,view,istrip |
149 |
|
150 |
integer view,is,istrip |
151 |
real coordsi |
152 |
|
153 |
include 'commontracker.f' |
154 |
|
155 |
c NB mettere il 1024 nel commontracker...!??? |
156 |
|
157 |
|
158 |
|
159 |
is=istrip !it stores istrip number |
160 |
is=mod(is-1,1024)+1 !it puts all clusters on a single ladder |
161 |
|
162 |
coordsi=0. |
163 |
|
164 |
if(mod(view,2).eq.0) then !X view |
165 |
|
166 |
c if((is.le.3).or.(is.ge.1022)) then !X has 1018 strips... |
167 |
c print*,'functions: WARNING: false X strip: strip ',is |
168 |
c endif |
169 |
|
170 |
is=is-3 !4 =< is =< 1021 --> 1 =< is =< 1018 |
171 |
|
172 |
edge=edgeX |
173 |
dim=SiDimX |
174 |
|
175 |
elseif(mod(view,2).eq.1) then !Y view |
176 |
|
177 |
edge=edgeY |
178 |
dim=SiDimY |
179 |
|
180 |
c$$$ if(view.eq.11) then !INVERSIONE!??? |
181 |
c$$$ is=1025-is |
182 |
c$$$ endif |
183 |
|
184 |
endif |
185 |
|
186 |
p=pitch(view) |
187 |
|
188 |
coord1=(is-1)*p !referred to 1st sensor strip |
189 |
coord1=coord1+edge !referred to sensor edge |
190 |
|
191 |
coordsi=coord1-dim/2 !referred to the centre of the sensor |
192 |
|
193 |
if(view.eq.11) then !INVERSION: it puts y axis in the same direction for all views |
194 |
coordsi=-coordsi |
195 |
endif |
196 |
|
197 |
end |
198 |
|
199 |
|
200 |
c------------------------------------------------------------------------ |
201 |
|
202 |
|
203 |
function acoordsi(strip,view) |
204 |
* |
205 |
* same as COORDSI, but accept a real value of strip!!! |
206 |
* |
207 |
* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * |
208 |
* it gives the strip coordinate in micrometers, |
209 |
* knowing the strip number (1..3072) and the view |
210 |
* number. the origin of the coordinate is on the |
211 |
* centre of the sensor the strip belongs to. |
212 |
* the axes directions are the same as in the PAMELA |
213 |
* reference frame (i.e.: the 11th view coordinate |
214 |
* direction has to be inverted here) |
215 |
* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * |
216 |
|
217 |
c integer is,view,istrip |
218 |
|
219 |
integer view,is,istrip |
220 |
real coordsi,acoordsi |
221 |
real strip,stripladder |
222 |
|
223 |
|
224 |
include 'commontracker.f' |
225 |
|
226 |
c NB mettere il 1024 nel commontracker...!??? |
227 |
|
228 |
istrip = int(strip+0.5) !istrip stores the closest integer to strip |
229 |
|
230 |
is=istrip !it stores istrip number |
231 |
is=mod(is-1,1024)+1 !it puts all clusters on a single ladder |
232 |
|
233 |
coordsi=0. |
234 |
|
235 |
if(mod(view,2).eq.0) then !X view |
236 |
|
237 |
c if((is.le.3).or.(is.ge.1022)) then !X has 1018 strips... |
238 |
c print*,'functions: WARNING: false X strip: strip ',is |
239 |
c endif |
240 |
|
241 |
is=is-3 !4 =< is =< 1021 --> 1 =< is =< 1018 |
242 |
|
243 |
edge=edgeX |
244 |
dim=SiDimX |
245 |
|
246 |
elseif(mod(view,2).eq.1) then !Y view |
247 |
|
248 |
edge=edgeY |
249 |
dim=SiDimY |
250 |
|
251 |
c$$$ if(view.eq.11) then !INVERSIONE!??? |
252 |
c$$$ is=1025-is |
253 |
c$$$ endif |
254 |
|
255 |
endif |
256 |
|
257 |
|
258 |
stripladder = float(is)+(strip-float(istrip))!cluster position relative to ladder |
259 |
p=pitch(view) |
260 |
|
261 |
ccccc coord1=(is-1)*p !referred to 1st sensor strip |
262 |
coord1=(stripladder-1)*p !referred to 1st sensor strip |
263 |
coord1=coord1+edge !referred to sensor edge |
264 |
acoordsi=coord1-dim/2 !referred to the centre of the sensor |
265 |
|
266 |
if(view.eq.11) then !INVERSION: it puts y axis in the same direction for all views |
267 |
acoordsi=-acoordsi |
268 |
endif |
269 |
|
270 |
end |
271 |
|
272 |
|
273 |
|
274 |
c------------------------------------------------------------------------ |
275 |
|
276 |
|
277 |
function coord(coordsi,view,ladder,sen) |
278 |
* it gives the coordinate in |
279 |
* micrometers, knowing the coordinate in the sensor |
280 |
* frame, the view, the ladder and the sensor numbers. |
281 |
* the origin is in the centre of the magnet (PAMELA |
282 |
* reference frame) |
283 |
|
284 |
include 'commontracker.f' |
285 |
include 'common_tracks.f' |
286 |
|
287 |
integer view,ladder,sen |
288 |
integer sx,sy,sz |
289 |
|
290 |
real coord,coordsi,trasl |
291 |
|
292 |
c$$$c parameter (offset=4365.) !??? ! in um !CONTROLLARE SE HA SENSO: |
293 |
c$$$ ! dalle misure sul piano dovrebbe essere 4970, |
294 |
c$$$ ! dallo shift dei residui viene 4365 |
295 |
c$$$ ! va messo .ne.0. se in mech_sensor assegno ai |
296 |
c$$$ ! sensori del sesto piano coordinate Y uguali |
297 |
c$$$ ! a quelle degli altri sensori |
298 |
c$$$ parameter (offset=0.) !??? altrimenti se il sesto piano ha coordinate |
299 |
c$$$ ! Y diverse offset dovrebbe essere .eq.0. |
300 |
c$$$ ! CONTROLLARE CON I GRAFICI DEI RESIDUI!!! |
301 |
|
302 |
|
303 |
coord=0. |
304 |
|
305 |
sx=ladder |
306 |
sy=sen |
307 |
sz=npl(view) |
308 |
|
309 |
if(mod(view,2).eq.0) then !X view |
310 |
|
311 |
trasl=x_mech_sensor(sz,sx,sy) !in mm |
312 |
|
313 |
elseif(mod(view,2).eq.1) then !Y view |
314 |
|
315 |
trasl=y_mech_sensor(sz,sx,sy) !in mm |
316 |
|
317 |
c$$$ if(view.eq.11) then !INVERSIONE!???INUTILE, ne e' gia' tenuto conto |
318 |
c$$$ coordsi=coordsi+offset ! in y_mech_pos... |
319 |
c$$$ endif |
320 |
|
321 |
endif |
322 |
|
323 |
coord=coordsi+trasl*1000. |
324 |
|
325 |
end |
326 |
|
327 |
|
328 |
c------------------------------------------------------------------------ |
329 |
c------------------------------------------------------------------------ |
330 |
|
331 |
c double precision version of the above subroutine |
332 |
|
333 |
double precision function dcoord(coordsi,view,ladder,sen) !it gives the coordinate in |
334 |
! micrometers, knowing the coordinate in the sensor |
335 |
! frame, the view, the ladder and the sensor numbers. |
336 |
! the origin is in the centre of the magnet (PAMELA |
337 |
! reference frame) |
338 |
|
339 |
include 'commontracker.f' |
340 |
include 'common_tracks.f' |
341 |
|
342 |
integer view,ladder,sen |
343 |
integer sx,sy,sz |
344 |
|
345 |
c double precision dcoord |
346 |
double precision coordsi,trasl |
347 |
|
348 |
c$$$c parameter (offset=4365.) !??? ! in um !CONTROLLARE SE HA SENSO: |
349 |
c$$$ ! dalle misure sul piano dovrebbe essere 4970, |
350 |
c$$$ ! dallo shift dei residui viene 4365 |
351 |
c$$$ ! va messo .ne.0. se in mech_sensor assegno ai |
352 |
c$$$ ! sensori del sesto piano coordinate Y uguali |
353 |
c$$$ ! a quelle degli altri sensori |
354 |
c$$$ parameter (offset=0.) !??? altrimenti se il sesto piano ha coordinate |
355 |
c$$$ ! Y diverse offset dovrebbe essere .eq.0. |
356 |
c$$$ ! CONTROLLARE CON I GRAFICI DEI RESIDUI!!! |
357 |
|
358 |
|
359 |
dcoord=0. |
360 |
|
361 |
sx=ladder |
362 |
sy=sen |
363 |
sz=npl(view) |
364 |
|
365 |
if(mod(view,2).eq.0) then !X view |
366 |
|
367 |
trasl=x_mech_sensor(sz,sx,sy) !in mm |
368 |
|
369 |
elseif(mod(view,2).eq.1) then !Y view |
370 |
|
371 |
trasl=y_mech_sensor(sz,sx,sy) !in mm |
372 |
|
373 |
c$$$ if(view.eq.11) then !INVERSIONE!???INUTILE, ne e' gia' tenuto conto |
374 |
c$$$ dcoordsi=dcoordsi+offset ! in y_mech_pos... |
375 |
c$$$ endif |
376 |
|
377 |
endif |
378 |
|
379 |
dcoord=coordsi+trasl*1000. |
380 |
|
381 |
end |
382 |
|
383 |
|
384 |
c------------------------------------------------------------------------ |
385 |
integer function nsatstrips(ic) |
386 |
*-------------------------------------------------------------- |
387 |
* this function returns the number of saturated strips |
388 |
* inside a cluster |
389 |
*-------------------------------------------------------------- |
390 |
include 'commontracker.f' |
391 |
include 'level1.f' |
392 |
include 'calib.f' |
393 |
|
394 |
|
395 |
|
396 |
integer nsat |
397 |
nsat = 0 |
398 |
iv=VIEW(ic) |
399 |
if(mod(iv,2).eq.1)incut=incuty |
400 |
if(mod(iv,2).eq.0)incut=incutx |
401 |
istart = INDSTART(IC) |
402 |
istop = TOTCLLENGTH |
403 |
if(ic.lt.NCLSTR1)istop=INDSTART(IC+1)-1 |
404 |
do i = INDMAX(IC),istart,-1 |
405 |
c cut = incut*CLSIGMA(i) |
406 |
c if(CLSIGNAL(i).ge.cut)then |
407 |
if( (mod(iv,2).eq.1.and.CLADC(i).lt.ADCsatx) |
408 |
$ .or. |
409 |
$ (mod(iv,2).eq.0.and.CLADC(i).gt.ADCsaty) )then |
410 |
nsat = nsat +1 |
411 |
else |
412 |
goto 10 |
413 |
endif |
414 |
enddo |
415 |
10 continue |
416 |
do i = INDMAX(IC)+1,istop |
417 |
c cut = incut*CLSIGMA(i) |
418 |
c if(CLSIGNAL(i).ge.cut)then |
419 |
if( (mod(iv,2).eq.1.and.CLADC(i).lt.ADCsatx) |
420 |
$ .or. |
421 |
$ (mod(iv,2).eq.0.and.CLADC(i).gt.ADCsaty) )then |
422 |
nsat = nsat +1 |
423 |
else |
424 |
goto 20 |
425 |
endif |
426 |
enddo |
427 |
20 continue |
428 |
|
429 |
nsatstrips = nsat |
430 |
return |
431 |
end |
432 |
|
433 |
c------------------------------------------------------------------------ |
434 |
integer function nbadstrips(ncog,ic) |
435 |
*-------------------------------------------------------------- |
436 |
* this function returns the number of BAD strips |
437 |
* inside a cluster: |
438 |
* - if NCOG=0, the number BAD strips inside the whole cluster |
439 |
* are given, according to the cluster multiplicity |
440 |
* |
441 |
* - if NCOG>0, the number BAD strips is evaluated using NCOG |
442 |
* strips, even if they have a negative signal (according to Landi) |
443 |
*-------------------------------------------------------------- |
444 |
include 'commontracker.f' |
445 |
include 'level1.f' |
446 |
include 'calib.f' |
447 |
|
448 |
integer nbad |
449 |
nbad = 0 |
450 |
|
451 |
if (ncog.gt.0) then |
452 |
|
453 |
* --> signal of the central strip |
454 |
sc = CLSIGNAL(INDMAX(ic)) !center |
455 |
* signal of adjacent strips |
456 |
sl1 = -100000 !left 1 |
457 |
if( |
458 |
$ (INDMAX(ic)-1).ge.INDSTART(ic) |
459 |
$ ) |
460 |
$ sl1 = CLSIGNAL(INDMAX(ic)-1) |
461 |
|
462 |
sl2 = -100000 !left 2 |
463 |
if( |
464 |
$ (INDMAX(ic)-2).ge.INDSTART(ic) |
465 |
$ ) |
466 |
$ sl2 = CLSIGNAL(INDMAX(ic)-2) |
467 |
|
468 |
sr1 = -100000 !right 1 |
469 |
if( |
470 |
$ (ic.ne.NCLSTR1.and.(INDMAX(ic)+1).lt.INDSTART(ic+1)) |
471 |
$ .or. |
472 |
$ (ic.eq.NCLSTR1.and.(INDMAX(ic)+1).le.TOTCLLENGTH) |
473 |
$ ) |
474 |
$ sr1 = CLSIGNAL(INDMAX(ic)+1) |
475 |
|
476 |
sr2 = -100000 !right 2 |
477 |
if( |
478 |
$ (ic.ne.NCLSTR1.and.(INDMAX(ic)+2).lt.INDSTART(ic+1)) |
479 |
$ .or. |
480 |
$ (ic.eq.NCLSTR1.and.(INDMAX(ic)+2).le.TOTCLLENGTH) |
481 |
$ ) |
482 |
$ sr2 = CLSIGNAL(INDMAX(ic)+2) |
483 |
|
484 |
if(ncog.ge.1)nbad = nbad+1-CLBAD(INDMAX(ic)) |
485 |
|
486 |
if(ncog.ge.2)then |
487 |
if(sl1.gt.sr1.and.(INDMAX(ic)-1).ge.1)then |
488 |
nbad=nbad+1-CLBAD(INDMAX(ic)-1) |
489 |
elseif(sl1.le.sr1.and.(INDMAX(ic)+1).le.NCLSTR1)then |
490 |
nbad=nbad+1-CLBAD(INDMAX(ic)+1) |
491 |
endif |
492 |
endif |
493 |
|
494 |
if(ncog.ge.3)then |
495 |
if(sl1.gt.sr1.and.(INDMAX(ic)+1).le.NCLSTR1)then |
496 |
nbad=nbad+1-CLBAD(INDMAX(ic)+1) |
497 |
elseif(sl1.le.sr1.and.(INDMAX(ic)-1).ge.1)then |
498 |
c if(INDMAX(ic)-1.eq.0) |
499 |
c $ print*,' ======= ',sl2,sl1,sc,sr1,sr2 |
500 |
nbad=nbad+1-CLBAD(INDMAX(ic)-1) |
501 |
endif |
502 |
endif |
503 |
|
504 |
if(ncog.ge.4)then |
505 |
if(sl2.gt.sr2.and.(INDMAX(ic)-2).ge.1)then |
506 |
nbad=nbad+1-CLBAD(INDMAX(ic)-2) |
507 |
elseif(sl2.le.sr2.and.(INDMAX(ic)+2).le.NCLSTR1)then |
508 |
nbad=nbad+1-CLBAD(INDMAX(ic)+2) |
509 |
endif |
510 |
endif |
511 |
|
512 |
c if(ncog.ge.5)then |
513 |
c print*,'function CLBAD(NCOG,IC) ==> WARNING!! NCOG=',NCOG |
514 |
c $ ,' not implemented' |
515 |
c endif |
516 |
|
517 |
elseif(ncog.eq.0)then |
518 |
* ========================= |
519 |
* COG computation |
520 |
* ========================= |
521 |
|
522 |
|
523 |
|
524 |
iv=VIEW(ic) |
525 |
if(mod(iv,2).eq.1)incut=incuty |
526 |
if(mod(iv,2).eq.0)incut=incutx |
527 |
|
528 |
istart = INDSTART(IC) |
529 |
istop = TOTCLLENGTH |
530 |
if(ic.lt.NCLSTR1)istop=INDSTART(IC+1)-1 |
531 |
nbad = 0 |
532 |
c$$$ do i=istart,istop |
533 |
c$$$ cut = incut*CLSIGMA(i) |
534 |
c$$$ if(CLSIGNAL(i).ge.cut)nbad = nbad +1 -CLBAD(i) |
535 |
c$$$ enddo |
536 |
do i = INDMAX(IC),istart,-1 |
537 |
cut = incut*CLSIGMA(i) |
538 |
if(CLSIGNAL(i).ge.cut)then |
539 |
nbad = nbad +1 -CLBAD(i) |
540 |
else |
541 |
goto 10 |
542 |
endif |
543 |
enddo |
544 |
10 continue |
545 |
do i = INDMAX(IC)+1,istop |
546 |
cut = incut*CLSIGMA(i) |
547 |
if(CLSIGNAL(i).ge.cut)then |
548 |
nbad = nbad +1 -CLBAD(i) |
549 |
else |
550 |
goto 20 |
551 |
endif |
552 |
enddo |
553 |
20 continue |
554 |
|
555 |
else |
556 |
|
557 |
c print*,'function CLBAD(NCOG,IC) ==> WARNING!! NCOG=',NCOG |
558 |
c $ ,' not implemented' |
559 |
|
560 |
|
561 |
endif |
562 |
|
563 |
nbadstrips = nbad |
564 |
|
565 |
return |
566 |
end |