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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.2 - (hide annotations) (download)
Fri Apr 27 10:39:58 2007 UTC (17 years, 7 months ago) by pam-fi
Branch: MAIN
Changes since 1.1: +179 -162 lines
v3r00: new hough parameters, new variables, and other things...

1 mocchiut 1.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     if((is.le.3).or.(is.ge.1022)) then !X has 1018 strips...
167     print*,'functions: WARNING: false X strip: strip ',is
168     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     if((is.le.3).or.(is.ge.1022)) then !X has 1018 strips...
238     print*,'functions: WARNING: false X strip: strip ',is
239     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 pam-fi 1.2 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 mocchiut 1.1
396 pam-fi 1.2 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 mocchiut 1.1
448 pam-fi 1.2 integer nbad
449     nbad = 0
450 mocchiut 1.1
451 pam-fi 1.2 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     if(ncog.ge.5)then
513     print*,'function CLBAD(NCOG,IC) ==> WARNING!! NCOG=',NCOG
514     $ ,' not implemented'
515     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     print*,'function CLBAD(NCOG,IC) ==> WARNING!! NCOG=',NCOG
558     $ ,' not implemented'
559    
560    
561     endif
562    
563     nbadstrips = nbad
564    
565     return
566     end

  ViewVC Help
Powered by ViewVC 1.1.23