/[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.5 - (hide annotations) (download)
Tue Sep 4 15:44:50 2007 UTC (17 years, 3 months ago) by bongi
Branch: MAIN
CVS Tags: v4r00
Changes since 1.4: +80 -0 lines
real to double precision change

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 bongi 1.5 c------------------------------------------------------------------------
137 mocchiut 1.1
138     function coordsi(istrip,view)
139     * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
140     * it gives the strip coordinate in micrometers,
141     * knowing the strip number (1..3072) and the view
142     * number. the origin of the coordinate is on the
143     * centre of the sensor the strip belongs to.
144     * the axes directions are the same as in the PAMELA
145     * reference frame (i.e.: the 11th view coordinate
146     * direction has to be inverted here)
147     * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
148    
149     c integer is,view,istrip
150    
151     integer view,is,istrip
152     real coordsi
153    
154     include 'commontracker.f'
155    
156     c NB mettere il 1024 nel commontracker...!???
157    
158    
159    
160     is=istrip !it stores istrip number
161     is=mod(is-1,1024)+1 !it puts all clusters on a single ladder
162    
163     coordsi=0.
164    
165     if(mod(view,2).eq.0) then !X view
166    
167 mocchiut 1.3 c if((is.le.3).or.(is.ge.1022)) then !X has 1018 strips...
168     c print*,'functions: WARNING: false X strip: strip ',is
169     c endif
170 mocchiut 1.1
171     is=is-3 !4 =< is =< 1021 --> 1 =< is =< 1018
172    
173     edge=edgeX
174     dim=SiDimX
175    
176     elseif(mod(view,2).eq.1) then !Y view
177    
178     edge=edgeY
179     dim=SiDimY
180    
181     c$$$ if(view.eq.11) then !INVERSIONE!???
182     c$$$ is=1025-is
183     c$$$ endif
184    
185     endif
186    
187     p=pitch(view)
188    
189     coord1=(is-1)*p !referred to 1st sensor strip
190     coord1=coord1+edge !referred to sensor edge
191    
192     coordsi=coord1-dim/2 !referred to the centre of the sensor
193    
194     if(view.eq.11) then !INVERSION: it puts y axis in the same direction for all views
195     coordsi=-coordsi
196     endif
197    
198     end
199    
200    
201     c------------------------------------------------------------------------
202    
203    
204     function acoordsi(strip,view)
205     *
206     * same as COORDSI, but accept a real value of strip!!!
207     *
208     * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
209     * it gives the strip coordinate in micrometers,
210     * knowing the strip number (1..3072) and the view
211     * number. the origin of the coordinate is on the
212     * centre of the sensor the strip belongs to.
213     * the axes directions are the same as in the PAMELA
214     * reference frame (i.e.: the 11th view coordinate
215     * direction has to be inverted here)
216     * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
217    
218     c integer is,view,istrip
219    
220     integer view,is,istrip
221     real coordsi,acoordsi
222     real strip,stripladder
223    
224    
225     include 'commontracker.f'
226    
227     c NB mettere il 1024 nel commontracker...!???
228    
229     istrip = int(strip+0.5) !istrip stores the closest integer to strip
230    
231     is=istrip !it stores istrip number
232     is=mod(is-1,1024)+1 !it puts all clusters on a single ladder
233    
234     coordsi=0.
235    
236     if(mod(view,2).eq.0) then !X view
237    
238 mocchiut 1.3 c if((is.le.3).or.(is.ge.1022)) then !X has 1018 strips...
239     c print*,'functions: WARNING: false X strip: strip ',is
240     c endif
241 mocchiut 1.1
242     is=is-3 !4 =< is =< 1021 --> 1 =< is =< 1018
243    
244     edge=edgeX
245     dim=SiDimX
246    
247     elseif(mod(view,2).eq.1) then !Y view
248    
249     edge=edgeY
250     dim=SiDimY
251    
252     c$$$ if(view.eq.11) then !INVERSIONE!???
253     c$$$ is=1025-is
254     c$$$ endif
255    
256     endif
257    
258    
259     stripladder = float(is)+(strip-float(istrip))!cluster position relative to ladder
260     p=pitch(view)
261    
262     ccccc coord1=(is-1)*p !referred to 1st sensor strip
263     coord1=(stripladder-1)*p !referred to 1st sensor strip
264     coord1=coord1+edge !referred to sensor edge
265     acoordsi=coord1-dim/2 !referred to the centre of the sensor
266    
267     if(view.eq.11) then !INVERSION: it puts y axis in the same direction for all views
268     acoordsi=-acoordsi
269 bongi 1.5 endif
270    
271     end
272    
273    
274    
275     c------------------------------------------------------------------------
276    
277    
278     double precision function dcoordsi(rstrip,view)
279     *
280     * same as COORDSI, but accept a real value of strip!!!
281     * and gives a double precision output
282     *
283     * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
284     * it gives the strip coordinate in micrometers,
285     * knowing the strip number (1..3072) and the view
286     * number. the origin of the coordinate is on the
287     * centre of the sensor the strip belongs to.
288     * the axes directions are the same as in the PAMELA
289     * reference frame (i.e.: the 11th view coordinate
290     * direction has to be inverted here)
291     * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *
292    
293     c integer is,view,istrip
294    
295     integer view,is,istrip
296     real rstrip
297     double precision strip,stripladder,p
298     double precision edge,dim
299     double precision coord1
300    
301    
302     include 'commontracker.f'
303    
304     c NB mettere il 1024 nel commontracker...!???
305    
306     strip=DBLE(rstrip)
307    
308     istrip = int(strip+0.5) !istrip stores the closest integer to strip
309    
310     is=istrip !it stores istrip number
311     is=mod(is-1,1024)+1 !it puts all clusters on a single ladder
312    
313     dcoordsi=0.
314    
315     if(mod(view,2).eq.0) then !X view
316    
317     c if((is.le.3).or.(is.ge.1022)) then !X has 1018 strips...
318     c print*,'functions: WARNING: false X strip: strip ',is
319     c endif
320    
321     is=is-3 !4 =< is =< 1021 --> 1 =< is =< 1018
322    
323     edge=edgeX
324     dim=SiDimX
325    
326     elseif(mod(view,2).eq.1) then !Y view
327    
328     edge=edgeY
329     dim=SiDimY
330    
331     c$$$ if(view.eq.11) then !INVERSIONE!???
332     c$$$ is=1025-is
333     c$$$ endif
334    
335     endif
336    
337    
338     stripladder = DBLE(is)+(strip-DBLE(istrip))!cluster position relative to ladder
339     p=pitch(view)
340    
341     ccccc coord1=(is-1)*p !referred to 1st sensor strip
342     coord1=(stripladder-1)*p !referred to 1st sensor strip
343     coord1=coord1+edge !referred to sensor edge
344     dcoordsi=coord1-dim/2 !referred to the centre of the sensor
345    
346     if(view.eq.11) then !INVERSION: it puts y axis in the same direction for all views
347     dcoordsi=-dcoordsi
348 mocchiut 1.1 endif
349    
350     end
351    
352    
353    
354     c------------------------------------------------------------------------
355    
356    
357     function coord(coordsi,view,ladder,sen)
358     * it gives the coordinate in
359     * micrometers, knowing the coordinate in the sensor
360     * frame, the view, the ladder and the sensor numbers.
361     * the origin is in the centre of the magnet (PAMELA
362     * reference frame)
363    
364     include 'commontracker.f'
365     include 'common_tracks.f'
366    
367     integer view,ladder,sen
368     integer sx,sy,sz
369    
370     real coord,coordsi,trasl
371    
372     c$$$c parameter (offset=4365.) !??? ! in um !CONTROLLARE SE HA SENSO:
373     c$$$ ! dalle misure sul piano dovrebbe essere 4970,
374     c$$$ ! dallo shift dei residui viene 4365
375     c$$$ ! va messo .ne.0. se in mech_sensor assegno ai
376     c$$$ ! sensori del sesto piano coordinate Y uguali
377     c$$$ ! a quelle degli altri sensori
378     c$$$ parameter (offset=0.) !??? altrimenti se il sesto piano ha coordinate
379     c$$$ ! Y diverse offset dovrebbe essere .eq.0.
380     c$$$ ! CONTROLLARE CON I GRAFICI DEI RESIDUI!!!
381    
382    
383     coord=0.
384    
385     sx=ladder
386     sy=sen
387     sz=npl(view)
388    
389     if(mod(view,2).eq.0) then !X view
390    
391     trasl=x_mech_sensor(sz,sx,sy) !in mm
392    
393     elseif(mod(view,2).eq.1) then !Y view
394    
395     trasl=y_mech_sensor(sz,sx,sy) !in mm
396    
397     c$$$ if(view.eq.11) then !INVERSIONE!???INUTILE, ne e' gia' tenuto conto
398     c$$$ coordsi=coordsi+offset ! in y_mech_pos...
399     c$$$ endif
400    
401     endif
402    
403     coord=coordsi+trasl*1000.
404    
405     end
406    
407    
408     c------------------------------------------------------------------------
409     c------------------------------------------------------------------------
410    
411     c double precision version of the above subroutine
412    
413     double precision function dcoord(coordsi,view,ladder,sen) !it gives the coordinate in
414     ! micrometers, knowing the coordinate in the sensor
415     ! frame, the view, the ladder and the sensor numbers.
416     ! the origin is in the centre of the magnet (PAMELA
417     ! reference frame)
418    
419     include 'commontracker.f'
420     include 'common_tracks.f'
421    
422     integer view,ladder,sen
423     integer sx,sy,sz
424    
425     c double precision dcoord
426     double precision coordsi,trasl
427    
428     c$$$c parameter (offset=4365.) !??? ! in um !CONTROLLARE SE HA SENSO:
429     c$$$ ! dalle misure sul piano dovrebbe essere 4970,
430     c$$$ ! dallo shift dei residui viene 4365
431     c$$$ ! va messo .ne.0. se in mech_sensor assegno ai
432     c$$$ ! sensori del sesto piano coordinate Y uguali
433     c$$$ ! a quelle degli altri sensori
434     c$$$ parameter (offset=0.) !??? altrimenti se il sesto piano ha coordinate
435     c$$$ ! Y diverse offset dovrebbe essere .eq.0.
436     c$$$ ! CONTROLLARE CON I GRAFICI DEI RESIDUI!!!
437    
438    
439     dcoord=0.
440    
441 pam-fi 1.4 if(
442     $ ladder.lt.1.or.
443     $ ladder.gt.nladders_view.or.
444     $ sen.lt.1.or.
445     $ sen.gt.2.or.
446     $ view.lt.1.or.
447     $ view.gt.nviews.or.
448     $ .false.)then
449     print*,'dcoord ---> wrong input: ladder ',ladder
450     $ ,' sensor ',sen
451     $ ,' view ',view
452     return
453     endif
454    
455 mocchiut 1.1 sx=ladder
456     sy=sen
457     sz=npl(view)
458    
459 pam-fi 1.4
460 mocchiut 1.1 if(mod(view,2).eq.0) then !X view
461    
462     trasl=x_mech_sensor(sz,sx,sy) !in mm
463    
464     elseif(mod(view,2).eq.1) then !Y view
465    
466     trasl=y_mech_sensor(sz,sx,sy) !in mm
467    
468     c$$$ if(view.eq.11) then !INVERSIONE!???INUTILE, ne e' gia' tenuto conto
469     c$$$ dcoordsi=dcoordsi+offset ! in y_mech_pos...
470     c$$$ endif
471    
472     endif
473    
474     dcoord=coordsi+trasl*1000.
475    
476     end
477    
478    
479     c------------------------------------------------------------------------
480 pam-fi 1.2 integer function nsatstrips(ic)
481     *--------------------------------------------------------------
482     * this function returns the number of saturated strips
483     * inside a cluster
484     *--------------------------------------------------------------
485     include 'commontracker.f'
486     include 'level1.f'
487     include 'calib.f'
488    
489    
490 mocchiut 1.1
491 pam-fi 1.2 integer nsat
492     nsat = 0
493     iv=VIEW(ic)
494     if(mod(iv,2).eq.1)incut=incuty
495     if(mod(iv,2).eq.0)incut=incutx
496     istart = INDSTART(IC)
497     istop = TOTCLLENGTH
498     if(ic.lt.NCLSTR1)istop=INDSTART(IC+1)-1
499     do i = INDMAX(IC),istart,-1
500     c cut = incut*CLSIGMA(i)
501     c if(CLSIGNAL(i).ge.cut)then
502     if( (mod(iv,2).eq.1.and.CLADC(i).lt.ADCsatx)
503     $ .or.
504     $ (mod(iv,2).eq.0.and.CLADC(i).gt.ADCsaty) )then
505     nsat = nsat +1
506     else
507     goto 10
508     endif
509     enddo
510     10 continue
511     do i = INDMAX(IC)+1,istop
512     c cut = incut*CLSIGMA(i)
513     c if(CLSIGNAL(i).ge.cut)then
514     if( (mod(iv,2).eq.1.and.CLADC(i).lt.ADCsatx)
515     $ .or.
516     $ (mod(iv,2).eq.0.and.CLADC(i).gt.ADCsaty) )then
517     nsat = nsat +1
518     else
519     goto 20
520     endif
521     enddo
522     20 continue
523    
524     nsatstrips = nsat
525     return
526     end
527    
528     c------------------------------------------------------------------------
529     integer function nbadstrips(ncog,ic)
530     *--------------------------------------------------------------
531     * this function returns the number of BAD strips
532     * inside a cluster:
533     * - if NCOG=0, the number BAD strips inside the whole cluster
534     * are given, according to the cluster multiplicity
535     *
536     * - if NCOG>0, the number BAD strips is evaluated using NCOG
537     * strips, even if they have a negative signal (according to Landi)
538     *--------------------------------------------------------------
539     include 'commontracker.f'
540     include 'level1.f'
541     include 'calib.f'
542 mocchiut 1.1
543 pam-fi 1.2 integer nbad
544     nbad = 0
545 mocchiut 1.1
546 pam-fi 1.2 if (ncog.gt.0) then
547    
548     * --> signal of the central strip
549     sc = CLSIGNAL(INDMAX(ic)) !center
550     * signal of adjacent strips
551     sl1 = -100000 !left 1
552     if(
553     $ (INDMAX(ic)-1).ge.INDSTART(ic)
554     $ )
555     $ sl1 = CLSIGNAL(INDMAX(ic)-1)
556    
557     sl2 = -100000 !left 2
558     if(
559     $ (INDMAX(ic)-2).ge.INDSTART(ic)
560     $ )
561     $ sl2 = CLSIGNAL(INDMAX(ic)-2)
562    
563     sr1 = -100000 !right 1
564     if(
565     $ (ic.ne.NCLSTR1.and.(INDMAX(ic)+1).lt.INDSTART(ic+1))
566     $ .or.
567     $ (ic.eq.NCLSTR1.and.(INDMAX(ic)+1).le.TOTCLLENGTH)
568     $ )
569     $ sr1 = CLSIGNAL(INDMAX(ic)+1)
570    
571     sr2 = -100000 !right 2
572     if(
573     $ (ic.ne.NCLSTR1.and.(INDMAX(ic)+2).lt.INDSTART(ic+1))
574     $ .or.
575     $ (ic.eq.NCLSTR1.and.(INDMAX(ic)+2).le.TOTCLLENGTH)
576     $ )
577     $ sr2 = CLSIGNAL(INDMAX(ic)+2)
578    
579     if(ncog.ge.1)nbad = nbad+1-CLBAD(INDMAX(ic))
580    
581     if(ncog.ge.2)then
582     if(sl1.gt.sr1.and.(INDMAX(ic)-1).ge.1)then
583     nbad=nbad+1-CLBAD(INDMAX(ic)-1)
584     elseif(sl1.le.sr1.and.(INDMAX(ic)+1).le.NCLSTR1)then
585     nbad=nbad+1-CLBAD(INDMAX(ic)+1)
586     endif
587     endif
588    
589     if(ncog.ge.3)then
590     if(sl1.gt.sr1.and.(INDMAX(ic)+1).le.NCLSTR1)then
591     nbad=nbad+1-CLBAD(INDMAX(ic)+1)
592     elseif(sl1.le.sr1.and.(INDMAX(ic)-1).ge.1)then
593     c if(INDMAX(ic)-1.eq.0)
594     c $ print*,' ======= ',sl2,sl1,sc,sr1,sr2
595     nbad=nbad+1-CLBAD(INDMAX(ic)-1)
596     endif
597     endif
598    
599     if(ncog.ge.4)then
600     if(sl2.gt.sr2.and.(INDMAX(ic)-2).ge.1)then
601     nbad=nbad+1-CLBAD(INDMAX(ic)-2)
602     elseif(sl2.le.sr2.and.(INDMAX(ic)+2).le.NCLSTR1)then
603     nbad=nbad+1-CLBAD(INDMAX(ic)+2)
604     endif
605     endif
606    
607 mocchiut 1.3 c if(ncog.ge.5)then
608     c print*,'function CLBAD(NCOG,IC) ==> WARNING!! NCOG=',NCOG
609     c $ ,' not implemented'
610     c endif
611 pam-fi 1.2
612     elseif(ncog.eq.0)then
613     * =========================
614     * COG computation
615     * =========================
616    
617    
618    
619     iv=VIEW(ic)
620     if(mod(iv,2).eq.1)incut=incuty
621     if(mod(iv,2).eq.0)incut=incutx
622    
623     istart = INDSTART(IC)
624     istop = TOTCLLENGTH
625     if(ic.lt.NCLSTR1)istop=INDSTART(IC+1)-1
626     nbad = 0
627     c$$$ do i=istart,istop
628     c$$$ cut = incut*CLSIGMA(i)
629     c$$$ if(CLSIGNAL(i).ge.cut)nbad = nbad +1 -CLBAD(i)
630     c$$$ enddo
631     do i = INDMAX(IC),istart,-1
632     cut = incut*CLSIGMA(i)
633     if(CLSIGNAL(i).ge.cut)then
634     nbad = nbad +1 -CLBAD(i)
635     else
636     goto 10
637     endif
638     enddo
639     10 continue
640     do i = INDMAX(IC)+1,istop
641     cut = incut*CLSIGMA(i)
642     if(CLSIGNAL(i).ge.cut)then
643     nbad = nbad +1 -CLBAD(i)
644     else
645     goto 20
646     endif
647     enddo
648     20 continue
649    
650     else
651    
652 mocchiut 1.3 c print*,'function CLBAD(NCOG,IC) ==> WARNING!! NCOG=',NCOG
653     c $ ,' not implemented'
654 pam-fi 1.2
655    
656     endif
657    
658     nbadstrips = nbad
659    
660     return
661     end

  ViewVC Help
Powered by ViewVC 1.1.23