/[PAMELA software]/DarthVader/TrackerLevel2/src/F77/functions.f
ViewVC logotype

Contents of /DarthVader/TrackerLevel2/src/F77/functions.f

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.3 - (show annotations) (download)
Fri Apr 27 12:11:03 2007 UTC (17 years, 7 months ago) by mocchiut
Branch: MAIN
CVS Tags: v3r04, v3r05, v3r06, v3r03
Changes since 1.2: +12 -12 lines
Comments commented

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

  ViewVC Help
Powered by ViewVC 1.1.23