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

Diff of /DarthVader/TrackerLevel2/src/F77/functionspfa.f

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1.2 by pam-fi, Tue May 30 16:30:37 2006 UTC revision 1.4 by pam-fi, Wed Oct 11 06:53:02 2006 UTC
# Line 7  Line 7 
7    
8    
9  *** * * * *** * * * *** * * * *** * * * *** * * * *** * * * ***  *** * * * *** * * * *** * * * *** * * * *** * * * *** * * * ***
10        real function pfa_eta(ic,angle)        real function pfaeta(ic,angle)
11  *--------------------------------------------------------------  *--------------------------------------------------------------
12  *     this function returns the position (in strip units)  *     this function returns the position (in strip units)
13  *     it calls:  *     it calls:
14  *     - pfa_eta2(ic,angle)  *     - pfaeta2(ic,angle)
15  *     - pfa_eta3(ic,angle)  *     - pfaeta3(ic,angle)
16  *     - pfa_eta4(ic,angle)  *     - pfaeta4(ic,angle)
17  *     according to the angle  *     according to the angle
18  *--------------------------------------------------------------  *--------------------------------------------------------------
19        include 'commontracker.f'        include 'commontracker.f'
20  c      include 'calib.f'  c      include 'calib.f'
21        include 'level1.f'        include 'level1.f'
22                
23        pfa_eta = 0        pfaeta = 0
24    
25        if(mod(int(VIEW(ic)),2).eq.1)then !Y-view        if(mod(int(VIEW(ic)),2).eq.1)then !Y-view
26                
27           pfa_eta = pfa_eta2(ic,angle)           pfaeta = pfaeta2(ic,angle)
28        
29        else                      !X-view        else                      !X-view
30    
31           if(abs(angle).le.10.)then           if(abs(angle).le.10.)then
32              pfa_eta = pfa_eta2(ic,angle)              pfaeta = pfaeta2(ic,angle)
33           elseif(abs(angle).gt.10..and.abs(angle).le.15.)then           elseif(abs(angle).gt.10..and.abs(angle).le.15.)then
34              pfa_eta = pfa_eta3(ic,angle)              pfaeta = pfaeta3(ic,angle)
35           elseif(abs(angle).gt.15.)then           elseif(abs(angle).gt.15.)then
36              pfa_eta = pfa_eta4(ic,angle)              pfaeta = pfaeta4(ic,angle)
37           endif           endif
38                            
39        endif        endif
# Line 46  c      include 'calib.f' Line 46  c      include 'calib.f'
46        real function ris_eta(ic,angle)        real function ris_eta(ic,angle)
47  *--------------------------------------------------------------  *--------------------------------------------------------------
48  *     this function returns the average spatial resolution  *     this function returns the average spatial resolution
49  *     (in cm) for the ETA algorithm (function pfa_eta(ic,angle))  *     (in cm) for the ETA algorithm (function pfaeta(ic,angle))
50  *     it calls:  *     it calls:
51  *     - risx_eta2(angle)  *     - risx_eta2(angle)
52  *     - risy_eta2(angle)  *     - risy_eta2(angle)
# Line 129  c$$$     $     ,' -->',ris_eta Line 129  c$$$     $     ,' -->',ris_eta
129  c*****************************************************  c*****************************************************
130  cccccc 02/02/2006 modified by Elena Vannuccini --> (1)  cccccc 02/02/2006 modified by Elena Vannuccini --> (1)
131  c*****************************************************  c*****************************************************
132  c      real function pfa_eta2(cog2,view,lad,angle)  c      real function pfaeta2(cog2,view,lad,angle)
133        real function pfa_eta2(ic,angle) !(1)        real function pfaeta2(ic,angle) !(1)
134  *--------------------------------------------------------------  *--------------------------------------------------------------
135  *     this function returns  *     this function returns
136  *  *
# Line 155  c      common/dbg/DEBUG Line 155  c      common/dbg/DEBUG
155        iview = VIEW(ic)              !(1)        iview = VIEW(ic)              !(1)
156        lad = nld(MAXS(ic),VIEW(ic)) !(1)        lad = nld(MAXS(ic),VIEW(ic)) !(1)
157        cog2 = cog(2,ic)             !(1)        cog2 = cog(2,ic)             !(1)
158        pfa_eta2=cog2        pfaeta2=cog2
159    
160  *     find angular bin  *     find angular bin
161  *     (in futuro possiamo pensare di interpolare anche sull'angolo)  *     (in futuro possiamo pensare di interpolare anche sull'angolo)
# Line 167  c         print*,'~~~~~~~~~~~~ ',iang,an Line 167  c         print*,'~~~~~~~~~~~~ ',iang,an
167           endif           endif
168        enddo        enddo
169        if(DEBUG)        if(DEBUG)
170       $     print*,'pfa_eta2 *** warning *** angle out of range: ',angle       $     print*,'pfaeta2 *** warning *** angle out of range: ',angle
171        if(angle.lt.angL(1))iang=1        if(angle.lt.angL(1))iang=1
172        if(angle.gt.angR(nangbin))iang=nangbin        if(angle.gt.angR(nangbin))iang=nangbin
173   98   continue                  !jump here if ok   98   continue                  !jump here if ok
# Line 177  c$$$*     find extremes of interpolation Line 177  c$$$*     find extremes of interpolation
177  c$$$      iflag=0  c$$$      iflag=0
178  c$$$*     --------------------------------  c$$$*     --------------------------------
179  c$$$      if(cog2.lt.eta2(1,iang).or.cog2.gt.eta2(netaval,iang))then  c$$$      if(cog2.lt.eta2(1,iang).or.cog2.gt.eta2(netaval,iang))then
180  c$$$c         print*,'pfa_eta2 *** warning *** argument out of range: ',cog2  c$$$c         print*,'pfaeta2 *** warning *** argument out of range: ',cog2
181  c$$$*         goto 100  c$$$*         goto 100
182  c$$$*     ----------------------------------------------  c$$$*     ----------------------------------------------
183  c$$$*     non salto piu`, ma scalo di 1 o -1  c$$$*     non salto piu`, ma scalo di 1 o -1
# Line 230  c            print*,'-----',x1,x2,y1,y2 Line 230  c            print*,'-----',x1,x2,y1,y2
230        AA=(y2-y1)/(x2-x1)        AA=(y2-y1)/(x2-x1)
231        BB=y1-AA*x1        BB=y1-AA*x1
232    
233        pfa_eta2 = AA*cog2+BB        pfaeta2 = AA*cog2+BB
234        pfa_eta2 = pfa_eta2 - iadd        pfaeta2 = pfaeta2 - iadd
235    
236  c$$$      if(iflag.eq.1)then  c$$$      if(iflag.eq.1)then
237  c$$$         pfa_eta2=pfa_eta2-1.   !temp  c$$$         pfaeta2=pfaeta2-1.   !temp
238  c$$$         cog2=cog2-1.           !temp  c$$$         cog2=cog2-1.           !temp
239  c$$$      endif  c$$$      endif
240  c$$$      if(iflag.eq.-1)then  c$$$      if(iflag.eq.-1)then
241  c$$$         pfa_eta2=pfa_eta2+1.   !temp  c$$$         pfaeta2=pfaeta2+1.   !temp
242  c$$$         cog2=cog2+1.           !temp  c$$$         cog2=cog2+1.           !temp
243  c$$$      endif  c$$$      endif
244    
245        if(DEBUG)print*,'ETA2  (ic ',ic,' ang',angle,')'        if(DEBUG)print*,'ETA2  (ic ',ic,' ang',angle,')'
246       $     ,cog2-iadd,' -->',pfa_eta2       $     ,cog2-iadd,' -->',pfaeta2
247    
248    
249   100  return   100  return
# Line 253  c$$$      endif Line 253  c$$$      endif
253  c*****************************************************  c*****************************************************
254  cccccc 02/02/2006 modified by Elena Vannuccini --> (1)  cccccc 02/02/2006 modified by Elena Vannuccini --> (1)
255  c*****************************************************  c*****************************************************
256  c      real function pfa_eta3(cog3,view,lad,angle)  c      real function pfaeta3(cog3,view,lad,angle)
257        real function pfa_eta3(ic,angle) !(1)        real function pfaeta3(ic,angle) !(1)
258  *--------------------------------------------------------------  *--------------------------------------------------------------
259  *     this function returns  *     this function returns
260  *  *
# Line 280  c      common/dbg/DEBUG Line 280  c      common/dbg/DEBUG
280        iview = VIEW(ic)              !(1)        iview = VIEW(ic)              !(1)
281        lad = nld(MAXS(ic),VIEW(ic)) !(1)        lad = nld(MAXS(ic),VIEW(ic)) !(1)
282        cog3 = cog(3,ic)             !(1)        cog3 = cog(3,ic)             !(1)
283        pfa_eta3=cog3        pfaeta3=cog3
284    
285  *     find angular bin  *     find angular bin
286  *     (in futuro possiamo pensare di interpolare anche sull'angolo)  *     (in futuro possiamo pensare di interpolare anche sull'angolo)
# Line 292  c         print*,'~~~~~~~~~~~~ ',iang,an Line 292  c         print*,'~~~~~~~~~~~~ ',iang,an
292           endif           endif
293        enddo        enddo
294        if(DEBUG)        if(DEBUG)
295       $     print*,'pfa_eta3 *** warning *** angle out of range: ',angle       $     print*,'pfaeta3 *** warning *** angle out of range: ',angle
296        if(angle.lt.angL(1))iang=1        if(angle.lt.angL(1))iang=1
297        if(angle.gt.angR(nangbin))iang=nangbin        if(angle.gt.angR(nangbin))iang=nangbin
298   98   continue                  !jump here if ok   98   continue                  !jump here if ok
# Line 354  c            print*,'-----',x1,x2,y1,y2 Line 354  c            print*,'-----',x1,x2,y1,y2
354        AA=(y2-y1)/(x2-x1)        AA=(y2-y1)/(x2-x1)
355        BB=y1-AA*x1        BB=y1-AA*x1
356    
357        pfa_eta3 = AA*cog3+BB        pfaeta3 = AA*cog3+BB
358        pfa_eta3 = pfa_eta3 - iadd        pfaeta3 = pfaeta3 - iadd
359    
360  c$$$      if(iflag.eq.1)then  c$$$      if(iflag.eq.1)then
361  c$$$         pfa_eta2=pfa_eta2-1.   !temp  c$$$         pfaeta2=pfaeta2-1.   !temp
362  c$$$         cog2=cog2-1.           !temp  c$$$         cog2=cog2-1.           !temp
363  c$$$      endif  c$$$      endif
364  c$$$      if(iflag.eq.-1)then  c$$$      if(iflag.eq.-1)then
365  c$$$         pfa_eta2=pfa_eta2+1.   !temp  c$$$         pfaeta2=pfaeta2+1.   !temp
366  c$$$         cog2=cog2+1.           !temp  c$$$         cog2=cog2+1.           !temp
367  c$$$      endif  c$$$      endif
368    
369        if(DEBUG)print*,'ETA3  (ic ',ic,' ang',angle,')'        if(DEBUG)print*,'ETA3  (ic ',ic,' ang',angle,')'
370       $     ,cog3-iadd,' -->',pfa_eta3       $     ,cog3-iadd,' -->',pfaeta3
371    
372   100  return   100  return
373        end        end
# Line 376  c$$$      endif Line 376  c$$$      endif
376  c*****************************************************  c*****************************************************
377  cccccc 02/02/2006 modified by Elena Vannuccini --> (1)  cccccc 02/02/2006 modified by Elena Vannuccini --> (1)
378  c*****************************************************  c*****************************************************
379  c      real function pfa_eta4(cog4,view,lad,angle)  c      real function pfaeta4(cog4,view,lad,angle)
380        real function pfa_eta4(ic,angle) !(1)        real function pfaeta4(ic,angle) !(1)
381  *--------------------------------------------------------------  *--------------------------------------------------------------
382  *     this function returns  *     this function returns
383  *  *
# Line 399  c      real function pfa_eta4(cog4,view, Line 399  c      real function pfa_eta4(cog4,view,
399  c      logical DEBUG  c      logical DEBUG
400  c      common/dbg/DEBUG  c      common/dbg/DEBUG
401    
402        iview = VIEW(ic)              !(1)        iview = VIEW(ic)             !(1)
403        lad = nld(MAXS(ic),VIEW(ic)) !(1)        lad = nld(MAXS(ic),VIEW(ic)) !(1)
404        cog4=cog(4,ic)               !(1)        cog4=cog(4,ic)               !(1)
405        pfa_eta4=cog4        pfaeta4=cog4
406    
407  *     find angular bin  *     find angular bin
408  *     (in futuro possiamo pensare di interpolare anche sull'angolo)  *     (in futuro possiamo pensare di interpolare anche sull'angolo)
# Line 414  c         print*,'~~~~~~~~~~~~ ',iang,an Line 414  c         print*,'~~~~~~~~~~~~ ',iang,an
414           endif           endif
415        enddo        enddo
416        if(DEBUG)        if(DEBUG)
417       $     print*,'pfa_eta4 *** warning *** angle out of range: ',angle       $     print*,'pfaeta4 *** warning *** angle out of range: ',angle
418        if(angle.lt.angL(1))iang=1        if(angle.lt.angL(1))iang=1
419        if(angle.gt.angR(nangbin))iang=nangbin        if(angle.gt.angR(nangbin))iang=nangbin
420   98   continue                  !jump here if ok   98   continue                  !jump here if ok
# Line 476  c            print*,'-----',x1,x2,y1,y2 Line 476  c            print*,'-----',x1,x2,y1,y2
476        AA=(y2-y1)/(x2-x1)        AA=(y2-y1)/(x2-x1)
477        BB=y1-AA*x1        BB=y1-AA*x1
478    
479        pfa_eta4 = AA*cog4+BB        pfaeta4 = AA*cog4+BB
480        pfa_eta4 = pfa_eta4 - iadd        pfaeta4 = pfaeta4 - iadd
481    
482  c$$$      if(iflag.eq.1)then  c$$$      if(iflag.eq.1)then
483  c$$$         pfa_eta2=pfa_eta2-1.   !temp  c$$$         pfaeta2=pfaeta2-1.   !temp
484  c$$$         cog2=cog2-1.           !temp  c$$$         cog2=cog2-1.           !temp
485  c$$$      endif  c$$$      endif
486  c$$$      if(iflag.eq.-1)then  c$$$      if(iflag.eq.-1)then
487  c$$$         pfa_eta2=pfa_eta2+1.   !temp  c$$$         pfaeta2=pfaeta2+1.   !temp
488  c$$$         cog2=cog2+1.           !temp  c$$$         cog2=cog2+1.           !temp
489  c$$$      endif  c$$$      endif
490    
491        if(DEBUG)print*,'ETA4  (ic ',ic,' ang',angle,')'        if(DEBUG)print*,'ETA4  (ic ',ic,' ang',angle,')'
492       $     ,cog4-iadd,' -->',pfa_eta4       $     ,cog4-iadd,' -->',pfaeta4
493    
494   100  return   100  return
495        end        end
# Line 655  c      common/dbg/DEBUG Line 655  c      common/dbg/DEBUG
655                    
656           COG = 0.           COG = 0.
657                    
658           if(ncog.eq.2)then           if(ncog.eq.1)then
659                COG = 0.
660             elseif(ncog.eq.2)then
661              if(sl1.gt.sr1)then              if(sl1.gt.sr1)then
662                 COG = -sl1/(sl1+sc)                         COG = -sl1/(sl1+sc)        
663              elseif(sl1.le.sr1)then              elseif(sl1.le.sr1)then
# Line 686  ' Line 688  '
688           if(mod(iv,2).eq.1)incut=incuty           if(mod(iv,2).eq.1)incut=incuty
689           if(mod(iv,2).eq.0)incut=incutx           if(mod(iv,2).eq.0)incut=incutx
690    
691           istart=INDSTART(IC)           istart = INDSTART(IC)
692           istop=TOTCLLENGTH           istop  = TOTCLLENGTH
693           if(ic.lt.NCLSTR1)istop=INDSTART(IC+1)-1           if(ic.lt.NCLSTR1)istop=INDSTART(IC+1)-1
694           COG=0             COG = 0  
695           mu=0           mu  = 0
696           do i=istart,istop           do i = istart,istop
697              ipos=i-INDMAX(ic)              ipos = i-INDMAX(ic)
698              ivk=nvk(MAXS(ic)+ipos)  cc            ivk  = nvk(MAXS(ic)+ipos)
699              is=nst(MAXS(ic)+ipos)  cc            is   = nst(MAXS(ic)+ipos)
700  *            print*,'******************',istart,istop,ipos  *            print*,'******************',istart,istop,ipos
701  *     $           ,MAXS(ic)+ipos,iv,ivk,is  *     $           ,MAXS(ic)+ipos,iv,ivk,is
702              cut=incut*SIGMA(iv,ivk,is)  cc            cut  = incut*SIGMA(iv,ivk,is)
703                cut  = incut*CLSIGMA(i)
704    c            if(SIGMA(iv,ivk,is).ne.CLSIGMA(i))
705    c     $           print*,'cog(0,ic) --> hai fatto qualche cazzata'
706              if(CLSIGNAL(i).ge.cut)then              if(CLSIGNAL(i).ge.cut)then
707                 COG = COG + ipos*CLSIGNAL(i)                 COG = COG + ipos*CLSIGNAL(i)
708                 mu = mu + 1                 mu = mu + 1
709  c               print*,ipos,CLSIGNAL(i),incut,cut  c               print*,ipos,CLSIGNAL(i),incut,cut
710              endif              endif
711           enddo           enddo
712           COG=COG/DEDX(ic)           if(DEDX(ic).le.0)then
713                print*,'cog(0,ic) --> ic, dedx ',ic,DEDX(ic)
714                print*,(CLSIGNAL(i)/CLSIGMA(i),i=istart,istop)
715                print*,(CLSIGNAL(i),i=istart,istop)
716                print*,'cog(0,ic) --> NOT EVALUATED '
717             else
718                COG=COG/DEDX(ic)
719             endif
720  c         if(DEBUG)print*,'COG  (ic ',ic,' m',mu,')'  c         if(DEBUG)print*,'COG  (ic ',ic,' m',mu,')'
721  c     $        ,cog  c     $        ,cog
722                    
Legend:
Removed from v.1.2  
changed lines
  Added in v.1.4
  ViewVC Help
Powered by ViewVC 1.1.23