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

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

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

revision 1.1 by mocchiut, Fri May 19 13:15:54 2006 UTC revision 1.6 by pam-fi, Mon Oct 2 13:12:48 2006 UTC
# Line 22  Line 22 
22    
23        include 'momanhough_init.f'        include 'momanhough_init.f'
24                
25        logical DEBUG  c      logical DEBUG
26        common/dbg/DEBUG  c      common/dbg/DEBUG
27    
28  *-------------------------------------------------------------------------------  *-------------------------------------------------------------------------------
29  *     STEP 1  *     STEP 1
# Line 154  c      iflag=0 Line 154  c      iflag=0
154    
155        include 'momanhough_init.f'        include 'momanhough_init.f'
156                
157        logical DEBUG  c      logical DEBUG
158        common/dbg/DEBUG  c      common/dbg/DEBUG
159    
160        logical FIMAGE            !        logical FIMAGE            !
161    
# Line 237  c         iflag=0 Line 237  c         iflag=0
237  *     ************************** FIT *** FIT *** FIT *** FIT ***  *     ************************** FIT *** FIT *** FIT *** FIT ***
238  *     **********************************************************  *     **********************************************************
239           do i=1,5           do i=1,5
240              AL(i)=dble(AL_STORE(i,icand))              AL(i)=dble(AL_STORE(i,icand))            
241           enddo           enddo
242             IDCAND = icand         !fitted track-candidate
243           ifail=0                !error flag in chi2 computation           ifail=0                !error flag in chi2 computation
244           jstep=0                !# minimization steps           jstep=0                !# minimization steps
245    
# Line 310  c         print*,'++++++++++ iimage,fima Line 311  c         print*,'++++++++++ iimage,fima
311  c     $        ,iimage,fimage,ntrk,image(ntrk)  c     $        ,iimage,fimage,ntrk,image(ntrk)
312    
313           if(ntrk.eq.NTRKMAX)then           if(ntrk.eq.NTRKMAX)then
314              if(DEBUG)              if(verbose)
315       $           print*,       $           print*,
316       $           '** warning ** number of identified '//       $           '** warning ** number of identified '//
317       $           'tracks exceeds vector dimension '       $           'tracks exceeds vector dimension '
# Line 613  c*************************************** Line 614  c***************************************
614        include 'common_xyzPAM.f'        include 'common_xyzPAM.f'
615        include 'common_resxy.f'        include 'common_resxy.f'
616    
617        logical DEBUG  c      logical DEBUG
618        common/dbg/DEBUG  c      common/dbg/DEBUG
619    
620        integer icx,icy           !X-Y cluster ID        integer icx,icy           !X-Y cluster ID
621        integer sensor        integer sensor
# Line 701  c            resxPAM = resxPAM*fbad_cog( Line 702  c            resxPAM = resxPAM*fbad_cog(
702           endif           endif
703    
704        endif        endif
705    c      if(icy.eq.0.and.icx.ne.0)
706    c     $     print*,PFAx,icx,angx,stripx,resxPAM,'***'
707                
708  *     -----------------  *     -----------------
709  *     CLUSTER Y  *     CLUSTER Y
# Line 772  C======================================= Line 775  C=======================================
775  c------------------------------------------------------------------------  c------------------------------------------------------------------------
776  c     (xi,yi,zi) = mechanical coordinates in the silicon sensor frame  c     (xi,yi,zi) = mechanical coordinates in the silicon sensor frame
777  c------------------------------------------------------------------------  c------------------------------------------------------------------------
778             if(((mod(int(stripx+0.5)-1,1024)+1).le.3)
779         $        .or.((mod(int(stripx+0.5)-1,1024)+1).ge.1022)) then !X has 1018 strips...
780                print*,'xyz_PAM (couple):',
781         $          ' WARNING: false X strip: strip ',stripx
782             endif
783           xi = acoordsi(stripx,viewx)           xi = acoordsi(stripx,viewx)
784           yi = acoordsi(stripy,viewy)           yi = acoordsi(stripy,viewy)
785           zi = 0.           zi = 0.
# Line 858  C======================================= Line 866  C=======================================
866              nldy = nldx              nldy = nldx
867              viewy = viewx - 1              viewy = viewx - 1
868    
869    c            print*,'X-singlet ',icx,nplx,nldx,viewx,stripx
870    c            if((stripx.le.3).or.(stripx.ge.1022)) then !X has 1018 strips...
871                if(((mod(int(stripx+0.5)-1,1024)+1).le.3)
872         $           .or.((mod(int(stripx+0.5)-1,1024)+1).ge.1022)) then !X has 1018 strips...
873                   print*,'xyz_PAM (X-singlet):',
874         $             ' WARNING: false X strip: strip ',stripx
875                endif
876              xi   = acoordsi(stripx,viewx)              xi   = acoordsi(stripx,viewx)
877    
878              xi_A = xi              xi_A = xi
# Line 1136  c$$$         print*,' resolution ',resxP Line 1151  c$$$         print*,' resolution ',resxP
1151  c------------------------------------------------------------------------  c------------------------------------------------------------------------
1152  c     (xi,yi,zi) = mechanical coordinates in the silicon sensor frame  c     (xi,yi,zi) = mechanical coordinates in the silicon sensor frame
1153  c------------------------------------------------------------------------  c------------------------------------------------------------------------
1154                   if(((mod(int(stripx+0.5)-1,1024)+1).le.3)
1155         $              .or.((mod(int(stripx+0.5)-1,1024)+1).ge.1022)) then !X has 1018 strips...
1156    c     if((stripx.le.3).or.(stripx.ge.1022)) then !X has 1018 strips...
1157                      print*,'whichsensor: ',
1158         $                ' WARNING: false X strip: strip ',stripx
1159                   endif
1160                 xi = acoordsi(stripx,viewx)                 xi = acoordsi(stripx,viewx)
1161                 yi = acoordsi(stripy,viewy)                 yi = acoordsi(stripy,viewy)
1162                 zi = 0.                 zi = 0.
# Line 1644  c$$$      end Line 1665  c$$$      end
1665        include 'calib.f'        include 'calib.f'
1666        include 'level1.f'        include 'level1.f'
1667    
1668        logical DEBUG  c      logical DEBUG
1669        common/dbg/DEBUG  c      common/dbg/DEBUG
1670    
1671  *     output flag  *     output flag
1672  *     --------------  *     --------------
# Line 1654  c$$$      end Line 1675  c$$$      end
1675  *     --------------  *     --------------
1676        integer iflag        integer iflag
1677    
1678        integer badseed,badcl        integer badseed,badclx,badcly
1679    
1680  *     init variables  *     init variables
1681        ncp_tot=0        ncp_tot=0
# Line 1681  c$$$      end Line 1702  c$$$      end
1702                    
1703  *     ----------------------------------------------------  *     ----------------------------------------------------
1704  *     cut on charge (X VIEW)  *     cut on charge (X VIEW)
1705    *     ----------------------------------------------------
1706           if(dedx(icx).lt.dedx_x_min)then           if(dedx(icx).lt.dedx_x_min)then
1707              cl_single(icx)=0              cl_single(icx)=0
1708              goto 10              goto 10
1709           endif           endif
1710    *     ----------------------------------------------------
1711  *     cut BAD (X VIEW)              *     cut BAD (X VIEW)            
1712    *     ----------------------------------------------------
1713           badseed=BAD(VIEW(icx),nvk(MAXS(icx)),nst(MAXS(icx)))           badseed=BAD(VIEW(icx),nvk(MAXS(icx)),nst(MAXS(icx)))
1714           ifirst=INDSTART(icx)           ifirst=INDSTART(icx)
1715           if(icx.ne.nclstr1) then           if(icx.ne.nclstr1) then
# Line 1693  c$$$      end Line 1717  c$$$      end
1717           else           else
1718              ilast=TOTCLLENGTH              ilast=TOTCLLENGTH
1719           endif           endif
1720           badcl=badseed           badclx=badseed
1721           do igood=-ngoodstr,ngoodstr           do igood=-ngoodstr,ngoodstr
1722              ibad=1              ibad=1
1723              if((INDMAX(icx)+igood).gt.ifirst.and.              if((INDMAX(icx)+igood).gt.ifirst.and.
# Line 1703  c$$$      end Line 1727  c$$$      end
1727       $              nvk(MAXS(icx)+igood),       $              nvk(MAXS(icx)+igood),
1728       $              nst(MAXS(icx)+igood))       $              nst(MAXS(icx)+igood))
1729              endif              endif
1730              badcl=badcl*ibad              badclx=badclx*ibad
1731           enddo           enddo
1732    *     ----------------------------------------------------
1733    *     >>> eliminato il taglio sulle BAD <<<
1734    *     ----------------------------------------------------
1735  c     if(badcl.eq.0)then  c     if(badcl.eq.0)then
1736  c     cl_single(icx)=0  c     cl_single(icx)=0
1737  c     goto 10  c     goto 10
# Line 1720  c     endif Line 1747  c     endif
1747                            
1748  *     ----------------------------------------------------  *     ----------------------------------------------------
1749  *     cut on charge (Y VIEW)  *     cut on charge (Y VIEW)
1750    *     ----------------------------------------------------
1751              if(dedx(icy).lt.dedx_y_min)then              if(dedx(icy).lt.dedx_y_min)then
1752                 cl_single(icy)=0                 cl_single(icy)=0
1753                 goto 20                 goto 20
1754              endif              endif
1755    *     ----------------------------------------------------
1756  *     cut BAD (Y VIEW)              *     cut BAD (Y VIEW)            
1757    *     ----------------------------------------------------
1758              badseed=BAD(VIEW(icy),nvk(MAXS(icy)),nst(MAXS(icy)))              badseed=BAD(VIEW(icy),nvk(MAXS(icy)),nst(MAXS(icy)))
1759              ifirst=INDSTART(icy)              ifirst=INDSTART(icy)
1760              if(icy.ne.nclstr1) then              if(icy.ne.nclstr1) then
# Line 1732  c     endif Line 1762  c     endif
1762              else              else
1763                 ilast=TOTCLLENGTH                 ilast=TOTCLLENGTH
1764              endif              endif
1765              badcl=badseed              badcly=badseed
1766              do igood=-ngoodstr,ngoodstr              do igood=-ngoodstr,ngoodstr
1767                 ibad=1                 ibad=1
1768                 if((INDMAX(icy)+igood).gt.ifirst.and.                 if((INDMAX(icy)+igood).gt.ifirst.and.
# Line 1741  c     endif Line 1771  c     endif
1771       $              ibad=BAD(VIEW(icy),       $              ibad=BAD(VIEW(icy),
1772       $              nvk(MAXS(icy)+igood),       $              nvk(MAXS(icy)+igood),
1773       $              nst(MAXS(icy)+igood))       $              nst(MAXS(icy)+igood))
1774                 badcl=badcl*ibad                 badcly=badcly*ibad
1775              enddo              enddo
1776    *     ----------------------------------------------------
1777    *     >>> eliminato il taglio sulle BAD <<<
1778    *     ----------------------------------------------------
1779  c     if(badcl.eq.0)then  c     if(badcl.eq.0)then
1780  c     cl_single(icy)=0  c     cl_single(icy)=0
1781  c     goto 20  c     goto 20
1782  c     endif  c     endif
1783  *     ----------------------------------------------------  *     ----------------------------------------------------
1784                            
               
1785              cl_good(icy)=1                                cl_good(icy)=1                  
1786              nply=npl(VIEW(icy))              nply=npl(VIEW(icy))
1787              nldy=nld(MAXS(icy),VIEW(icy))              nldy=nld(MAXS(icy),VIEW(icy))
# Line 1760  c     endif Line 1792  c     endif
1792  *     geometrical consistency (same plane and ladder)  *     geometrical consistency (same plane and ladder)
1793              if(nply.eq.nplx.and.nldy.eq.nldx)then              if(nply.eq.nplx.and.nldy.eq.nldx)then
1794  *     charge correlation  *     charge correlation
1795                 ddd=(dedx(icy)  *     (modified to be applied only below saturation... obviously)
1796       $              -kch(nplx,nldx)*dedx(icx)-cch(nplx,nldx))  
1797                 ddd=ddd/sqrt(kch(nplx,nldx)**2+1)  *     -------------------------------------------------------------
1798                 cut=chcut*sch(nplx,nldx)  *     >>> eliminata (TEMPORANEAMENTE) la correlazione di carica <<<
1799                 if(abs(ddd).gt.cut)goto 20 !charge not consistent  *     -------------------------------------------------------------
1800                                 if(  .not.(dedx(icy).gt.chsaty.and.dedx(icx).gt.chsatx)
1801         $              .and.
1802         $              .not.(dedx(icy).lt.chmipy.and.dedx(icx).lt.chmipx)
1803         $              .and.
1804         $              (badclx.eq.1.and.badcly.eq.1)
1805         $              .and.
1806         $              .true.)then
1807    
1808                      ddd=(dedx(icy)
1809         $                 -kch(nplx,nldx)*dedx(icx)-cch(nplx,nldx))
1810                      ddd=ddd/sqrt(kch(nplx,nldx)**2+1)
1811    
1812    c                  cut = chcut * sch(nplx,nldx)
1813    
1814                      sss=(kch(nplx,nldx)*dedx(icy)+dedx(icx)
1815         $                 -kch(nplx,nldx)*cch(nplx,nldx))
1816                      sss=sss/sqrt(kch(nplx,nldx)**2+1)
1817                      cut = chcut * (16 + sss/50.)
1818    
1819                      if(abs(ddd).gt.cut)then
1820                         goto 20    !charge not consistent
1821                      endif
1822                   endif
1823                                
1824  *     ------------------> COUPLE <------------------  *     ------------------> COUPLE <------------------
1825  *     check to do not overflow vector dimentions  *     check to do not overflow vector dimentions
1826                 if(ncp_plane(nplx).gt.ncouplemax)then  c$$$               if(ncp_plane(nplx).gt.ncouplemax)then
1827                    if(DEBUG)print*,  c$$$                  if(DEBUG)print*,
1828       $                    ' ** warning ** number of identified'//  c$$$     $                    ' ** warning ** number of identified'//
1829       $                    ' couples on plane ',nplx,  c$$$     $                    ' couples on plane ',nplx,
1830       $                    ' exceeds vector dimention'//  c$$$     $                    ' exceeds vector dimention'//
1831       $                    ' ( ',ncouplemax,' )'  c$$$     $                    ' ( ',ncouplemax,' )'
1832  c     good2=.false.  c$$$c     good2=.false.
1833  c     goto 880   !fill ntp and go to next event  c$$$c     goto 880   !fill ntp and go to next event
1834                    iflag=1  c$$$                  iflag=1
1835                    return  c$$$                  return
1836                 endif  c$$$               endif
1837                                
1838                 if(ncp_plane(nplx).eq.ncouplemax)then                 if(ncp_plane(nplx).eq.ncouplemax)then
1839                    if(DEBUG)print*,                    if(verbose)print*,
1840       $                 '** warning ** number of identified '//       $                 '** warning ** number of identified '//
1841       $                 'couples on plane ',nplx,       $                 'couples on plane ',nplx,
1842       $                 'exceeds vector dimention '       $                 'exceeds vector dimention '
# Line 1830  c     goto 880   !fill ntp and go to nex Line 1884  c     goto 880   !fill ntp and go to nex
1884  c     if(ncp_tot.gt.ncp_max)goto 100!next event (TEMPORANEO!!!)  c     if(ncp_tot.gt.ncp_max)goto 100!next event (TEMPORANEO!!!)
1885                
1886        if(ncp_tot.gt.ncp_max)then        if(ncp_tot.gt.ncp_max)then
1887           if(DEBUG)print*,           if(verbose)print*,
1888       $           '** warning ** number of identified '//       $           '** warning ** number of identified '//
1889       $           'couples exceeds upper limit for Hough tr. '       $           'couples exceeds upper limit for Hough tr. '
1890       $           ,'( ',ncp_max,' )'                   $           ,'( ',ncp_max,' )'            
# Line 1859  c     goto 880       !fill ntp and go to Line 1913  c     goto 880       !fill ntp and go to
1913        include 'calib.f'        include 'calib.f'
1914        include 'level1.f'        include 'level1.f'
1915    
1916        logical DEBUG  c      logical DEBUG
1917        common/dbg/DEBUG  c      common/dbg/DEBUG
1918    
1919  *     output flag  *     output flag
1920  *     --------------  *     --------------
# Line 1989  c$$$               if(abs(ddd).gt.cut)go Line 2043  c$$$               if(abs(ddd).gt.cut)go
2043                                
2044  *     ------------------> COUPLE <------------------  *     ------------------> COUPLE <------------------
2045  *     check to do not overflow vector dimentions  *     check to do not overflow vector dimentions
2046                 if(ncp_plane(nplx).gt.ncouplemax)then  c$$$               if(ncp_plane(nplx).gt.ncouplemax)then
2047                    if(DEBUG)print*,  c$$$                  if(DEBUG)print*,
2048       $                    ' ** warning ** number of identified'//  c$$$     $                    ' ** warning ** number of identified'//
2049       $                    ' couples on plane ',nplx,  c$$$     $                    ' couples on plane ',nplx,
2050       $                    ' exceeds vector dimention'//  c$$$     $                    ' exceeds vector dimention'//
2051       $                    ' ( ',ncouplemax,' )'  c$$$     $                    ' ( ',ncouplemax,' )'
2052  c     good2=.false.  c$$$c     good2=.false.
2053  c     goto 880   !fill ntp and go to next event  c$$$c     goto 880   !fill ntp and go to next event
2054                    iflag=1  c$$$                  iflag=1
2055                    return  c$$$                  return
2056                 endif  c$$$               endif
2057                                
2058                 if(ncp_plane(nplx).eq.ncouplemax)then                 if(ncp_plane(nplx).eq.ncouplemax)then
2059                    if(DEBUG)print*,                    if(verbose)print*,
2060       $                 '** warning ** number of identified '//       $                 '** warning ** number of identified '//
2061       $                 'couples on plane ',nplx,       $                 'couples on plane ',nplx,
2062       $                 'exceeds vector dimention '       $                 'exceeds vector dimention '
# Line 2050  c     goto 880   !fill ntp and go to nex Line 2104  c     goto 880   !fill ntp and go to nex
2104  c     if(ncp_tot.gt.ncp_max)goto 100!next event (TEMPORANEO!!!)  c     if(ncp_tot.gt.ncp_max)goto 100!next event (TEMPORANEO!!!)
2105                
2106        if(ncp_tot.gt.ncp_max)then        if(ncp_tot.gt.ncp_max)then
2107           if(DEBUG)print*,           if(verbose)print*,
2108       $           '** warning ** number of identified '//       $           '** warning ** number of identified '//
2109       $           'couples exceeds upper limit for Hough tr. '       $           'couples exceeds upper limit for Hough tr. '
2110       $           ,'( ',ncp_max,' )'                   $           ,'( ',ncp_max,' )'            
# Line 2063  c     goto 880       !fill ntp and go to Line 2117  c     goto 880       !fill ntp and go to
2117        return        return
2118        end        end
2119    
 c$$$      subroutine cl_to_couples_2(iflag)  
 c$$$  
 c$$$      include 'commontracker.f'  
 c$$$      include 'common_momanhough.f'  
 c$$$      include 'momanhough_init.f'  
 c$$$      include 'calib.f'  
 c$$$      include 'level1.f'  
 c$$$  
 c$$$      logical DEBUG  
 c$$$      common/dbg/DEBUG  
 c$$$  
 c$$$*     output flag  
 c$$$*     --------------  
 c$$$*     0 = good event  
 c$$$*     1 = bad event  
 c$$$*     --------------  
 c$$$      integer iflag  
 c$$$  
 c$$$      integer badseed,badcl  
 c$$$  
 c$$$*     init variables  
 c$$$      ncp_tot=0  
 c$$$      do ip=1,nplanes  
 c$$$         do ico=1,ncouplemax  
 c$$$            clx(ip,ico)=0  
 c$$$            cly(ip,ico)=0  
 c$$$         enddo  
 c$$$         ncp_plane(ip)=0  
 c$$$         do icl=1,nclstrmax_level2  
 c$$$            cls(ip,icl)=1  
 c$$$         enddo  
 c$$$         ncls(ip)=0  
 c$$$      enddo  
 c$$$      do icl=1,nclstrmax_level2  
 c$$$         cl_single(icl)=1  
 c$$$         cl_good(icl)=0  
 c$$$      enddo  
 c$$$        
 c$$$*     start association  
 c$$$      ncouples=0  
 c$$$      do icx=1,nclstr1          !loop on cluster (X)  
 c$$$         if(mod(VIEW(icx),2).eq.1)goto 10  
 c$$$          
 c$$$*     ----------------------------------------------------  
 c$$$*     cut on charge (X VIEW)  
 c$$$         if(dedx(icx).lt.dedx_x_min)then  
 c$$$            cl_single(icx)=0  
 c$$$            goto 10  
 c$$$         endif  
 c$$$*     cut BAD (X VIEW)              
 c$$$         badseed=BAD(VIEW(icx),nvk(MAXS(icx)),nst(MAXS(icx)))  
 c$$$         ifirst=INDSTART(icx)  
 c$$$         if(icx.ne.nclstr1) then  
 c$$$            ilast=INDSTART(icx+1)-1  
 c$$$         else  
 c$$$            ilast=TOTCLLENGTH  
 c$$$         endif  
 c$$$         badcl=badseed  
 c$$$         do igood=-ngoodstr,ngoodstr  
 c$$$            ibad=1  
 c$$$            if((INDMAX(icx)+igood).gt.ifirst.and.  
 c$$$     $           (INDMAX(icx)+igood).lt.ilast.and.  
 c$$$     $           .true.)then  
 c$$$               ibad=BAD(VIEW(icx),  
 c$$$     $              nvk(MAXS(icx)+igood),  
 c$$$     $              nst(MAXS(icx)+igood))  
 c$$$            endif  
 c$$$            badcl=badcl*ibad  
 c$$$         enddo  
 c$$$*         print*,'icx ',icx,badcl  
 c$$$         if(badcl.eq.0)then  
 c$$$            cl_single(icx)=0  
 c$$$            goto 10  
 c$$$         endif  
 c$$$*     ----------------------------------------------------  
 c$$$          
 c$$$         cl_good(icx)=1  
 c$$$         nplx=npl(VIEW(icx))  
 c$$$         nldx=nld(MAXS(icx),VIEW(icx))  
 c$$$          
 c$$$         do icy=1,nclstr1       !loop on cluster (Y)  
 c$$$            if(mod(VIEW(icy),2).eq.0)goto 20  
 c$$$              
 c$$$*     ----------------------------------------------------  
 c$$$*     cut on charge (Y VIEW)  
 c$$$            if(dedx(icy).lt.dedx_y_min)then  
 c$$$               cl_single(icy)=0  
 c$$$               goto 20  
 c$$$            endif  
 c$$$*     cut BAD (Y VIEW)              
 c$$$            badseed=BAD(VIEW(icy),nvk(MAXS(icy)),nst(MAXS(icy)))  
 c$$$            ifirst=INDSTART(icy)  
 c$$$            if(icy.ne.nclstr1) then  
 c$$$               ilast=INDSTART(icy+1)-1  
 c$$$            else  
 c$$$               ilast=TOTCLLENGTH  
 c$$$            endif  
 c$$$            badcl=badseed  
 c$$$            do igood=-ngoodstr,ngoodstr  
 c$$$               ibad=1  
 c$$$               if((INDMAX(icy)+igood).gt.ifirst.and.  
 c$$$     $              (INDMAX(icy)+igood).lt.ilast.and.  
 c$$$     $              .true.)  
 c$$$     $              ibad=BAD(VIEW(icy),  
 c$$$     $              nvk(MAXS(icy)+igood),  
 c$$$     $              nst(MAXS(icy)+igood))  
 c$$$               badcl=badcl*ibad  
 c$$$            enddo  
 c$$$*            print*,'icy ',icy,badcl  
 c$$$            if(badcl.eq.0)then  
 c$$$               cl_single(icy)=0  
 c$$$               goto 20  
 c$$$            endif  
 c$$$*     ----------------------------------------------------  
 c$$$              
 c$$$              
 c$$$            cl_good(icy)=1                    
 c$$$            nply=npl(VIEW(icy))  
 c$$$            nldy=nld(MAXS(icy),VIEW(icy))  
 c$$$              
 c$$$*     ----------------------------------------------  
 c$$$*     CONDITION TO FORM A COUPLE  
 c$$$*     ----------------------------------------------  
 c$$$*     geometrical consistency (same plane and ladder)  
 c$$$            if(nply.eq.nplx.and.nldy.eq.nldx)then  
 c$$$  
 c$$$c$$$*     charge correlation  
 c$$$c$$$               ddd=(dedx(icy)  
 c$$$c$$$     $              -kch(nplx,nldx)*dedx(icx)-cch(nplx,nldx))  
 c$$$c$$$               ddd=ddd/sqrt(kch(nplx,nldx)**2+1)  
 c$$$c$$$               cut=chcut*sch(nplx,nldx)  
 c$$$c$$$               if(abs(ddd).gt.cut)goto 20 !charge not consistent  
 c$$$                
 c$$$*     ------------------> COUPLE <------------------  
 c$$$*     check to do not overflow vector dimentions  
 c$$$               if(ncp_plane(nplx).gt.ncouplemax)then  
 c$$$                  if(DEBUG)print*,  
 c$$$     $                    ' ** warning ** number of identified'//  
 c$$$     $                    ' couples on plane ',nplx,  
 c$$$     $                    ' exceeds vector dimention'//  
 c$$$     $                    ' ( ',ncouplemax,' )'  
 c$$$c     good2=.false.  
 c$$$c     goto 880   !fill ntp and go to next event  
 c$$$                  iflag=1  
 c$$$                  return  
 c$$$               endif  
 c$$$                
 c$$$               if(ncp_plane(nplx).eq.ncouplemax)then  
 c$$$                  if(DEBUG)print*,  
 c$$$     $                 '** warning ** number of identified '//  
 c$$$     $                 'couples on plane ',nplx,  
 c$$$     $                 'exceeds vector dimention '  
 c$$$     $                 ,'( ',ncouplemax,' )'  
 c$$$c     good2=.false.  
 c$$$c     goto 880   !fill ntp and go to next event                      
 c$$$                  iflag=1  
 c$$$                  return  
 c$$$               endif  
 c$$$                
 c$$$               ncp_plane(nplx) = ncp_plane(nplx) + 1  
 c$$$               clx(nplx,ncp_plane(nplx))=icx  
 c$$$               cly(nply,ncp_plane(nplx))=icy  
 c$$$               cl_single(icx)=0  
 c$$$               cl_single(icy)=0  
 c$$$c               print*,'couple ',nplx,ncp_plane(nplx),' --- ',icx,icy  
 c$$$            endif                                
 c$$$*     ----------------------------------------------  
 c$$$  
 c$$$ 20         continue  
 c$$$         enddo                  !end loop on clusters(Y)  
 c$$$          
 c$$$ 10      continue  
 c$$$      enddo                     !end loop on clusters(X)  
 c$$$        
 c$$$        
 c$$$      do icl=1,nclstr1  
 c$$$         if(cl_single(icl).eq.1)then  
 c$$$            ip=npl(VIEW(icl))  
 c$$$            ncls(ip)=ncls(ip)+1  
 c$$$            cls(ip,ncls(ip))=icl  
 c$$$         endif  
 c$$$      enddo  
 c$$$        
 c$$$        
 c$$$      if(DEBUG)then  
 c$$$         print*,'clusters  ',nclstr1  
 c$$$         print*,'good    ',(cl_good(i),i=1,nclstr1)  
 c$$$         print*,'singles ',(cl_single(i),i=1,nclstr1)  
 c$$$         print*,'couples per plane: ',(ncp_plane(ip),ip=1,nplanes)  
 c$$$      endif  
 c$$$        
 c$$$      do ip=1,6  
 c$$$         ncp_tot=ncp_tot+ncp_plane(ip)  
 c$$$      enddo  
 c$$$c     if(ncp_tot.gt.ncp_max)goto 100!next event (TEMPORANEO!!!)  
 c$$$        
 c$$$      if(ncp_tot.gt.ncp_max)then  
 c$$$         if(DEBUG)print*,  
 c$$$     $           '** warning ** number of identified '//  
 c$$$     $           'couples exceeds upper limit for Hough tr. '  
 c$$$     $           ,'( ',ncp_max,' )'              
 c$$$c            good2=.false.  
 c$$$c     goto 880       !fill ntp and go to next event  
 c$$$         iflag=1  
 c$$$         return  
 c$$$      endif  
 c$$$        
 c$$$      return  
 c$$$      end  
2120                
2121  ***************************************************  ***************************************************
2122  *                                                 *  *                                                 *
# Line 2295  c*************************************** Line 2140  c***************************************
2140        include 'calib.f'        include 'calib.f'
2141        include 'level1.f'        include 'level1.f'
2142    
2143        logical DEBUG  c      logical DEBUG
2144        common/dbg/DEBUG  c      common/dbg/DEBUG
2145    
2146  *     output flag  *     output flag
2147  *     --------------  *     --------------
# Line 2364  c     $                       (icx2,icy2 Line 2209  c     $                       (icx2,icy2
2209  *     (2 couples needed)  *     (2 couples needed)
2210  *     - - - - - - - - - - - - - - - - - - - - - - - - - - - -  *     - - - - - - - - - - - - - - - - - - - - - - - - - - - -
2211                          if(ndblt.eq.ndblt_max)then                          if(ndblt.eq.ndblt_max)then
2212                             if(DEBUG)print*,                             if(verbose)print*,
2213       $                          '** warning ** number of identified '//       $                          '** warning ** number of identified '//
2214       $                          'doublets exceeds vector dimention '       $                          'doublets exceeds vector dimention '
2215       $                          ,'( ',ndblt_max,' )'       $                          ,'( ',ndblt_max,' )'
# Line 2434  c     $                                 Line 2279  c     $                                
2279  *     (3 couples needed)  *     (3 couples needed)
2280  *     - - - - - - - - - - - - - - - - - - - - - - - - - - - -  *     - - - - - - - - - - - - - - - - - - - - - - - - - - - -
2281                                   if(ntrpt.eq.ntrpt_max)then                                   if(ntrpt.eq.ntrpt_max)then
2282                                      if(DEBUG)print*,                                      if(verbose)print*,
2283       $                     '** warning ** number of identified '//       $                     '** warning ** number of identified '//
2284       $                     'triplets exceeds vector dimention '       $                     'triplets exceeds vector dimention '
2285       $                    ,'( ',ntrpt_max,' )'       $                    ,'( ',ntrpt_max,' )'
# Line 2517  c     goto 880               !ntp fill Line 2362  c     goto 880               !ntp fill
2362        include 'common_momanhough.f'        include 'common_momanhough.f'
2363        include 'momanhough_init.f'        include 'momanhough_init.f'
2364    
2365        logical DEBUG  c      logical DEBUG
2366        common/dbg/DEBUG  c      common/dbg/DEBUG
2367    
2368  *     output flag  *     output flag
2369  *     --------------  *     --------------
# Line 2661  c     print*,'>>>> ',ncpused,npt,nplused Line 2506  c     print*,'>>>> ',ncpused,npt,nplused
2506  *     >>> NEW CLOUD <<<  *     >>> NEW CLOUD <<<
2507    
2508           if(nclouds_yz.ge.ncloyz_max)then           if(nclouds_yz.ge.ncloyz_max)then
2509              if(DEBUG)print*,              if(verbose)print*,
2510       $           '** warning ** number of identified '//       $           '** warning ** number of identified '//
2511       $           'YZ clouds exceeds vector dimention '       $           'YZ clouds exceeds vector dimention '
2512       $           ,'( ',ncloyz_max,' )'       $           ,'( ',ncloyz_max,' )'
# Line 2733  c$$$     $           ,(db_cloud(iii),iii Line 2578  c$$$     $           ,(db_cloud(iii),iii
2578        include 'common_momanhough.f'        include 'common_momanhough.f'
2579        include 'momanhough_init.f'        include 'momanhough_init.f'
2580    
2581        logical DEBUG  c      logical DEBUG
2582        common/dbg/DEBUG  c      common/dbg/DEBUG
2583    
2584  *     output flag  *     output flag
2585  *     --------------  *     --------------
# Line 2873  c     print*,'check cp_used' Line 2718  c     print*,'check cp_used'
2718  *     ~~~~~~~~~~~~~~~~~  *     ~~~~~~~~~~~~~~~~~
2719  *     >>> NEW CLOUD <<<  *     >>> NEW CLOUD <<<
2720           if(nclouds_xz.ge.ncloxz_max)then           if(nclouds_xz.ge.ncloxz_max)then
2721              if(DEBUG)print*,              if(verbose)print*,
2722       $           '** warning ** number of identified '//       $           '** warning ** number of identified '//
2723       $           'XZ clouds exceeds vector dimention '       $           'XZ clouds exceeds vector dimention '
2724       $           ,'( ',ncloxz_max,' )'       $           ,'( ',ncloxz_max,' )'
# Line 2947  c*************************************** Line 2792  c***************************************
2792        include 'common_mech.f'        include 'common_mech.f'
2793        include 'momanhough_init.f'        include 'momanhough_init.f'
2794    
2795        logical DEBUG  c      logical DEBUG
2796        common/dbg/DEBUG  c      common/dbg/DEBUG
2797    
2798  *     output flag  *     output flag
2799  *     --------------  *     --------------
# Line 3173  c     $                                 Line 3018  c     $                                
3018  *     --------------------------  *     --------------------------
3019                                if(ntracks.eq.NTRACKSMAX)then                                if(ntracks.eq.NTRACKSMAX)then
3020                                                                    
3021                                   if(DEBUG)print*,                                   if(verbose)print*,
3022       $                 '** warning ** number of candidate tracks '//       $                 '** warning ** number of candidate tracks '//
3023       $                 ' exceeds vector dimension '       $                 ' exceeds vector dimension '
3024       $                ,'( ',NTRACKSMAX,' )'       $                ,'( ',NTRACKSMAX,' )'
# Line 3273  c$$$  rchi2=chi2/dble(ndof) Line 3118  c$$$  rchi2=chi2/dble(ndof)
3118  c******************************************************  c******************************************************
3119  cccccc 06/10/2005 modified by elena vannuccini ---> (1)  cccccc 06/10/2005 modified by elena vannuccini ---> (1)
3120  cccccc 31/01/2006 modified by elena vannuccini ---> (2)  cccccc 31/01/2006 modified by elena vannuccini ---> (2)
3121    cccccc 12/08/2006 modified by elena vannucicni ---> (3)
3122  c******************************************************  c******************************************************
3123    
3124        include 'commontracker.f'        include 'commontracker.f'
# Line 3284  c*************************************** Line 3130  c***************************************
3130        include 'level1.f'        include 'level1.f'
3131        include 'calib.f'        include 'calib.f'
3132    
3133        logical DEBUG  c      logical DEBUG
3134        common/dbg/DEBUG  c      common/dbg/DEBUG
3135    
3136  *     flag to chose PFA  *     flag to chose PFA
3137        character*10 PFA        character*10 PFA
# Line 3383  c            if(DEBUG)print*,'>>>> try t Line 3229  c            if(DEBUG)print*,'>>>> try t
3229                 icx=clx(ip,icp)                 icx=clx(ip,icp)
3230                 icy=cly(ip,icp)                 icy=cly(ip,icp)
3231                 if(LADDER(icx).ne.nldt.or. !If the ladder number does not match                 if(LADDER(icx).ne.nldt.or. !If the ladder number does not match
3232       $              cl_used(icx).eq.1.or. !or the X cluster is already used  c     $              cl_used(icx).eq.1.or. !or the X cluster is already used
3233       $              cl_used(icy).eq.1.or. !or the Y cluster is already used  c     $              cl_used(icy).eq.1.or. !or the Y cluster is already used
3234         $              cl_used(icx).ne.0.or. !or the X cluster is already used !(3)
3235         $              cl_used(icy).ne.0.or. !or the Y cluster is already used !(3)
3236       $              .false.)goto 1188 !then jump to next couple.       $              .false.)goto 1188 !then jump to next couple.
3237  *            *          
3238                 call xyz_PAM(icx,icy,ist,                 call xyz_PAM(icx,icy,ist,
# Line 3456  c            if(DEBUG)print*,'>>>> try t Line 3304  c            if(DEBUG)print*,'>>>> try t
3304                 if(LADDER(icx).ne.nldt)goto 11882 !if the ladder number does not match                 if(LADDER(icx).ne.nldt)goto 11882 !if the ladder number does not match
3305  *                                                !jump to the next couple  *                                                !jump to the next couple
3306  *----- try cluster x -----------------------------------------------  *----- try cluster x -----------------------------------------------
3307                 if(cl_used(icx).eq.1)goto 11881 !if the X cluster is already used  c               if(cl_used(icx).eq.1)goto 11881 !if the X cluster is already used
3308                   if(cl_used(icx).ne.0)goto 11881 !if the X cluster is already used  !(3)
3309  *                                              !jump to the Y cluster  *                                              !jump to the Y cluster
3310                 call xyz_PAM(icx,0,ist,                 call xyz_PAM(icx,0,ist,
3311  c     $              'ETA2','ETA2',  c     $              'ETA2','ETA2',
# Line 3483  c                  dedxmm = dedx(icx) !( Line 3332  c                  dedxmm = dedx(icx) !(
3332                 endif                                   endif                  
3333  11881          continue  11881          continue
3334  *----- try cluster y -----------------------------------------------  *----- try cluster y -----------------------------------------------
3335                 if(cl_used(icy).eq.1)goto 11882 !if the Y cluster is already used  c               if(cl_used(icy).eq.1)goto 11882 !if the Y cluster is already used
3336                   if(cl_used(icy).ne.0)goto 11882 !if the Y cluster is already used !(3)
3337  *                                              !jump to the next couple  *                                              !jump to the next couple
3338                 call xyz_PAM(0,icy,ist,                 call xyz_PAM(0,icy,ist,
3339  c     $              'ETA2','ETA2',  c     $              'ETA2','ETA2',
# Line 3512  c                 dedxmm = dedx(icy)  !( Line 3362  c                 dedxmm = dedx(icy)  !(
3362  *----- single clusters -----------------------------------------------      *----- single clusters -----------------------------------------------    
3363              do ic=1,ncls(ip)    !loop on single clusters              do ic=1,ncls(ip)    !loop on single clusters
3364                 icl=cls(ip,ic)                 icl=cls(ip,ic)
3365                 if(cl_used(icl).eq.1.or.     !if the cluster is already used  c               if(cl_used(icl).eq.1.or.     !if the cluster is already used
3366                   if(cl_used(icl).ne.0.or.     !if the cluster is already used !(3)
3367       $              LADDER(icl).ne.nldt.or. !or the ladder number does not match       $              LADDER(icl).ne.nldt.or. !or the ladder number does not match
3368       $              .false.)goto 18882      !jump to the next singlet       $              .false.)goto 18882      !jump to the next singlet
3369                 if(mod(VIEW(icl),2).eq.0)then!<---- X view                 if(mod(VIEW(icl),2).eq.0)then!<---- X view
# Line 3596  c              dedxtrk(nplanes-ip+1) = d Line 3447  c              dedxtrk(nplanes-ip+1) = d
3447  *                                                 *  *                                                 *
3448  *                                                 *  *                                                 *
3449  **************************************************  **************************************************
3450    cccccc 12/08/2006 modified by elena ---> (1)
3451    *
3452        subroutine clean_XYclouds(ibest,iflag)        subroutine clean_XYclouds(ibest,iflag)
3453    
3454        include 'commontracker.f'        include 'commontracker.f'
3455        include 'common_momanhough.f'        include 'common_momanhough.f'
3456        include 'momanhough_init.f'        include 'momanhough_init.f'
3457          include 'level2.f'        !(1)
3458  c      include 'calib.f'  c      include 'calib.f'
3459  c      include 'level1.f'  c      include 'level1.f'
3460    
3461        logical DEBUG  c      logical DEBUG
3462        common/dbg/DEBUG  c      common/dbg/DEBUG
3463    
3464    
3465        do ip=1,nplanes           !loop on planes        do ip=1,nplanes           !loop on planes
# Line 3617  c      include 'level1.f' Line 3470  c      include 'level1.f'
3470              if(id.ne.0)then              if(id.ne.0)then
3471                 iclx=clx(ip,icp_cp(id))                 iclx=clx(ip,icp_cp(id))
3472                 icly=cly(ip,icp_cp(id))                 icly=cly(ip,icp_cp(id))
3473                 cl_used(iclx)=1  !tag used clusters  c               cl_used(iclx)=1  !tag used clusters
3474                 cl_used(icly)=1  !tag used clusters  c               cl_used(icly)=1  !tag used clusters
3475                   cl_used(iclx)=ntrk  !tag used clusters !(1)
3476                   cl_used(icly)=ntrk  !tag used clusters !(1)
3477              elseif(icl.ne.0)then              elseif(icl.ne.0)then
3478                 cl_used(icl)=1   !tag used clusters  c               cl_used(icl)=1   !tag used clusters
3479                   cl_used(icl)=ntrk   !tag used clusters !1)
3480              endif              endif
3481                            
3482  c               if(DEBUG)then  c               if(DEBUG)then
# Line 3788  c*************************************** Line 3644  c***************************************
3644        include 'level2.f'        include 'level2.f'
3645        include 'level1.f'        include 'level1.f'
3646    
3647          do i=1,nviews
3648             good2(i)=good1(i)
3649          enddo
3650    
3651    c      good2 = 0!.false.
       good2 = 0!.false.  
3652  c$$$      nev2 = nev1  c$$$      nev2 = nev1
3653    
3654  c$$$# ifndef TEST2003  c$$$# ifndef TEST2003
# Line 3817  c*************************************** Line 3675  c***************************************
3675        do it=1,NTRKMAX!NTRACKSMAX        do it=1,NTRKMAX!NTRACKSMAX
3676           IMAGE(IT)=0           IMAGE(IT)=0
3677           CHI2_nt(IT) = -100000.           CHI2_nt(IT) = -100000.
3678           BdL(IT) = 0.  c         BdL(IT) = 0.
3679           do ip=1,nplanes           do ip=1,nplanes
3680              XM_nt(IP,IT) = 0              XM_nt(IP,IT) = 0
3681              YM_nt(IP,IT) = 0              YM_nt(IP,IT) = 0
# Line 3831  cccccc 11/9/2005 modified by david fedel Line 3689  cccccc 11/9/2005 modified by david fedel
3689              DEDX_X(IP,IT) = 0              DEDX_X(IP,IT) = 0
3690              DEDX_Y(IP,IT) = 0              DEDX_Y(IP,IT) = 0
3691  c******************************************************  c******************************************************
3692    cccccc 17/8/2006 modified by elena
3693                CLTRX(IP,IT) = 0
3694                CLTRY(IP,IT) = 0
3695           enddo           enddo
3696           do ipa=1,5           do ipa=1,5
3697              AL_nt(IPA,IT) = 0              AL_nt(IPA,IT) = 0
# Line 3882  c*************************************** Line 3743  c***************************************
3743    
3744            
3745        include 'commontracker.f'        include 'commontracker.f'
3746          include 'level1.f'
3747        include 'level2.f'        include 'level2.f'
3748        include 'common_mini_2.f'        include 'common_mini_2.f'
3749          include 'common_momanhough.f'
3750        real sinth,phi,pig        !(4)        real sinth,phi,pig        !(4)
3751        pig=acos(-1.)        pig=acos(-1.)
3752    
3753        good2=1!.true.  c      good2=1!.true.
3754        chi2_nt(ntr)        = sngl(chi2)        chi2_nt(ntr)        = sngl(chi2)
3755          nstep_nt(ntr)       = 0!nstep
3756    
3757        phi   = al(4)             !(4)        phi   = al(4)             !(4)
3758        sinth = al(3)             !(4)        sinth = al(3)             !(4)
# Line 3926  c     print*,al_nt(i,ntr) Line 3790  c     print*,al_nt(i,ntr)
3790           ayv_nt(ip,ntr)   = sngl(ayv(ip))           ayv_nt(ip,ntr)   = sngl(ayv(ip))
3791  c        dedxp(ip,ntr)    = sngl(dedxtrk(ip))   !(1)  c        dedxp(ip,ntr)    = sngl(dedxtrk(ip))   !(1)
3792           dedx_x(ip,ntr)   = sngl(dedxtrk_x(ip)) !(2)           dedx_x(ip,ntr)   = sngl(dedxtrk_x(ip)) !(2)
3793           dedx_y(ip,ntr)   = sngl(dedxtrk_y(ip)) !(2)                     dedx_y(ip,ntr)   = sngl(dedxtrk_y(ip)) !(2)  
3794      
3795             id  = CP_STORE(ip,IDCAND)
3796             icl = CLS_STORE(ip,IDCAND)
3797             if(id.ne.0)then
3798                cltrx(ip,ntr)   = clx(nplanes-ip+1,icp_cp(id))
3799                cltry(ip,ntr)   = cly(nplanes-ip+1,icp_cp(id))
3800    c            print*,ip,' ',cltrx(ip,ntr),cltry(ip,ntr)
3801             elseif(icl.ne.0)then
3802                if(mod(VIEW(icl),2).eq.0)cltrx(ip,ntr)=icl
3803                if(mod(VIEW(icl),2).eq.1)cltry(ip,ntr)=icl
3804    c            print*,ip,' ',cltrx(ip,ntr),cltry(ip,ntr)
3805             endif          
3806    
3807        enddo        enddo
3808  c      call CalcBdL(100,xxxx,IFAIL)  c      call CalcBdL(100,xxxx,IFAIL)
3809  c      if(ifps(xxxx).eq.1)BdL(ntr) = xxxx  c      if(ifps(xxxx).eq.1)BdL(ntr) = xxxx
# Line 3959  c*************************************** Line 3836  c***************************************
3836        include 'common_xyzPAM.f'        include 'common_xyzPAM.f'
3837    
3838  *     count #cluster per plane not associated to any track  *     count #cluster per plane not associated to any track
3839        good2=1!.true.  c      good2=1!.true.
3840        nclsx = 0        nclsx = 0
3841        nclsy = 0        nclsy = 0
3842    
# Line 3970  c*************************************** Line 3847  c***************************************
3847                 nclsx = nclsx + 1                 nclsx = nclsx + 1
3848                 planex(nclsx) = ip                 planex(nclsx) = ip
3849                 sgnlxs(nclsx) = dedx(icl)/mip(VIEW(icl),LADDER(icl))!(2)                 sgnlxs(nclsx) = dedx(icl)/mip(VIEW(icl),LADDER(icl))!(2)
3850                   clsx(nclsx)   = icl
3851                 do is=1,2                 do is=1,2
3852  c                  call xyz_PAM(icl,0,is,'COG1',' ',0.,0.)  c                  call xyz_PAM(icl,0,is,'COG1',' ',0.,0.)
3853                    call xyz_PAM(icl,0,is,PFAdef,' ',0.,0.)                    call xyz_PAM(icl,0,is,PFAdef,' ',0.,0.)
# Line 3984  c$$$               print*,'xs(2,nclsx)   Line 3862  c$$$               print*,'xs(2,nclsx)  
3862                 nclsy = nclsy + 1                 nclsy = nclsy + 1
3863                 planey(nclsy) = ip                 planey(nclsy) = ip
3864                 sgnlys(nclsy) = dedx(icl)/mip(VIEW(icl),LADDER(icl))!(2)                 sgnlys(nclsy) = dedx(icl)/mip(VIEW(icl),LADDER(icl))!(2)
3865                   clsy(nclsy)   = icl
3866                 do is=1,2                 do is=1,2
3867  c                  call xyz_PAM(0,icl,is,' ','COG1',0.,0.)  c                  call xyz_PAM(0,icl,is,' ','COG1',0.,0.)
3868                    call xyz_PAM(0,icl,is,' ',PFAdef,0.,0.)                    call xyz_PAM(0,icl,is,' ',PFAdef,0.,0.)
# Line 3997  c$$$               print*,'ys(2,nclsy)   Line 3876  c$$$               print*,'ys(2,nclsy)  
3876              endif              endif
3877           endif           endif
3878  c      print*,icl,cl_used(icl),cl_good(icl),ip,VIEW(icl)!nclsx(ip),nclsy(ip)  c      print*,icl,cl_used(icl),cl_good(icl),ip,VIEW(icl)!nclsx(ip),nclsy(ip)
3879    
3880    ***** LO METTO QUI PERCHE` NON SO DOVE METTERLO
3881             whichtrack(icl) = cl_used(icl)
3882    
3883        enddo        enddo
3884        end        end
3885    
# Line 4004  c      print*,icl,cl_used(icl),cl_good(i Line 3887  c      print*,icl,cl_used(icl),cl_good(i
3887    
3888    
3889    
3890    
Legend:
Removed from v.1.1  
changed lines
  Added in v.1.6
  ViewVC Help
Powered by ViewVC 1.1.23