/[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.10 by pam-fi, Thu Nov 2 11:19:47 2006 UTC revision 1.13 by pam-fi, Wed Nov 8 16:42:28 2006 UTC
# Line 183  c      iflag=0 Line 183  c      iflag=0
183  c      include 'momanhough_init.f'  c      include 'momanhough_init.f'
184                
185        logical FIMAGE            !        logical FIMAGE            !
186          real*8 AL_GUESS(5)
187    
188  *-------------------------------------------------------------------------------  *-------------------------------------------------------------------------------
189  *     STEP 4   (ITERATED until any other physical track isn't found)  *     STEP 4   (ITERATED until any other physical track isn't found)
# Line 281  c$$$         if(ibest.eq.0)goto 880 !>> Line 282  c$$$         if(ibest.eq.0)goto 880 !>>
282  *     **********************************************************  *     **********************************************************
283  *     ************************** FIT *** FIT *** FIT *** FIT ***  *     ************************** FIT *** FIT *** FIT *** FIT ***
284  *     **********************************************************  *     **********************************************************
285             call guess()
286             do i=1,5
287                AL_GUESS(i)=AL(i)
288             enddo
289    
290           do i=1,5           do i=1,5
291              AL(i)=dble(AL_STORE(i,icand))                          AL(i)=dble(AL_STORE(i,icand))            
292           enddo           enddo
293            
294           IDCAND = icand         !fitted track-candidate           IDCAND = icand         !fitted track-candidate
295           ifail=0                !error flag in chi2 computation           ifail=0                !error flag in chi2 computation
296           jstep=0                !# minimization steps           jstep=0                !# minimization steps
297    
298           iprint=0           iprint=0
299           if(DEBUG)iprint=1  c         if(DEBUG)iprint=1
300             if(VERBOSE)iprint=1
301           call mini2(jstep,ifail,iprint)           call mini2(jstep,ifail,iprint)
302           if(ifail.ne.0) then           if(ifail.ne.0) then
303              if(DEBUG)then              if(.true.)then
304                 print *,                 print *,
305       $              '*** MINIMIZATION FAILURE *** (mini2) '       $              '*** MINIMIZATION FAILURE *** (after refinement) '
306       $              ,iev       $              ,iev
307    
308    c$$$               print*,'guess:   ',(al_guess(i),i=1,5)
309    c$$$               print*,'previous: ',(al_store(i,icand),i=1,5)
310    c$$$               print*,'result:   ',(al(i),i=1,5)
311    c$$$               print*,'xgood ',xgood
312    c$$$               print*,'ygood ',ygood
313    c$$$               print*,'----------------------------------------------'
314              endif              endif
315              chi2=-chi2  c            chi2=-chi2
316           endif           endif
317                    
318           if(DEBUG)then           if(DEBUG)then
# Line 660  c      include 'level1.f' Line 675  c      include 'level1.f'
675        include 'common_align.f'        include 'common_align.f'
676        include 'common_mech.f'        include 'common_mech.f'
677        include 'common_xyzPAM.f'        include 'common_xyzPAM.f'
678        include 'common_resxy.f'  c      include 'common_resxy.f'
679    
680  c      logical DEBUG  c      logical DEBUG
681  c      common/dbg/DEBUG  c      common/dbg/DEBUG
# Line 1756  c      include 'level1.f' Line 1771  c      include 'level1.f'
1771        do iv=1,nviews        do iv=1,nviews
1772           if( ncl_view(iv).gt. nclustermax)then           if( ncl_view(iv).gt. nclustermax)then
1773              mask_view(iv) = 1              mask_view(iv) = 1
1774              if(VERBOSE)print*,' * WARNING * cl_to_couple: n.clusters > '              if(DEBUG)print*,' * WARNING * cl_to_couple: n.clusters > '
1775       $           ,nclustermax,' on view ', iv,' --> masked!'       $           ,nclustermax,' on view ', iv,' --> masked!'
1776           endif           endif
1777        enddo        enddo
# Line 2912  c      include 'momanhough_init.f' Line 2927  c      include 'momanhough_init.f'
2927  *     -----------------------------------------------------------  *     -----------------------------------------------------------
2928  *     variables for track fitting  *     variables for track fitting
2929        double precision AL_INI(5)        double precision AL_INI(5)
2930        double precision tath  c      double precision tath
2931  *     -----------------------------------------------------------  *     -----------------------------------------------------------
2932  c      real fitz(nplanes)        !z coordinates of the planes in cm  c      real fitz(nplanes)        !z coordinates of the planes in cm
2933    
# Line 3005  c$$$  print*,'6 -- ',(cly(6,i),i=1,ncp_p Line 3020  c$$$  print*,'6 -- ',(cly(6,i),i=1,ncp_p
3020  c$$$  print*,'~~~~~~~~~~~~~~~~~~~~~~~~~'  c$$$  print*,'~~~~~~~~~~~~~~~~~~~~~~~~~'
3021                            
3022  *     -------> INITIAL GUESS <-------  *     -------> INITIAL GUESS <-------
3023              AL_INI(1) = dreal(alfaxz1_av(ixz))  cccc       SBAGLIATO
3024              AL_INI(2) = dreal(alfayz1_av(iyz))  c$$$            AL_INI(1) = dreal(alfaxz1_av(ixz))
3025              AL_INI(4) = PIGR + datan(dreal(alfayz2_av(iyz))  c$$$            AL_INI(2) = dreal(alfayz1_av(iyz))
3026       $           /dreal(alfaxz2_av(ixz)))  c$$$            AL_INI(4) = PIGR + datan(dreal(alfayz2_av(iyz))
3027              tath      = -dreal(alfaxz2_av(ixz))/dcos(AL_INI(4))  c$$$     $           /dreal(alfaxz2_av(ixz)))
3028              AL_INI(3) = tath/sqrt(1+tath**2)  c$$$            tath      = -dreal(alfaxz2_av(ixz))/dcos(AL_INI(4))
3029              AL_INI(5) = (1.e2*alfaxz3_av(ixz))/(0.3*0.43) !0.  c$$$            AL_INI(3) = tath/sqrt(1+tath**2)
3030                c$$$            AL_INI(5) = (1.e2*alfaxz3_av(ixz))/(0.3*0.43) !0.
3031  c     print*,'*******',AL_INI(5)  cccc       GIUSTO (ma si sua guess())
3032              if(AL_INI(5).gt.defmax)goto 888 !next cloud  c$$$            AL_INI(1) = dreal(alfaxz1_av(ixz))
3033                c$$$            AL_INI(2) = dreal(alfayz1_av(iyz))
3034  c     print*,'alfaxz2, alfayz2 '  c$$$            tath      = -dreal(alfaxz2_av(ixz))/dcos(AL_INI(4))
3035  c     $              ,alfaxz2_av(ixz),alfayz2_av(iyz)  c$$$            AL_INI(3) = tath/sqrt(1+tath**2)
3036                c$$$            IF(alfaxz2_av(ixz).NE.0)THEN
3037  *     -------> INITIAL GUESS <-------  c$$$            AL_INI(4) = PIGR + datan(dreal(alfayz2_av(iyz))
3038  c     print*,'AL_INI ',(al_ini(i),i=1,5)  c$$$     $           /dreal(alfaxz2_av(ixz)))
3039                c$$$            ELSE
3040    c$$$               AL_INI(4) = acos(-1.)/2
3041    c$$$               IF(alfayz2_av(iyz).LT.0)AL_INI(4) = AL_INI(4)+acos(-1.)
3042    c$$$            ENDIF
3043    c$$$            IF(alfaxz2_av(ixz).LT.0)AL_INI(4)= acos(-1.)+ AL_INI(4)
3044    c$$$            AL_INI(4) = -acos(-1.) + AL_INI(4) !from incidence direction to tracking rs
3045    c$$$            
3046    c$$$            AL_INI(5) = (1.e2*alfaxz3_av(ixz))/(0.3*0.43) !0.
3047    c$$$            
3048    c$$$            if(AL_INI(5).gt.defmax)goto 888 !next cloud
3049                            
3050              if(DEBUG)then              if(DEBUG)then
3051                 print*,'1 >>> ',(cp_match(6,i),i=1,ncp_match(6))                 print*,'1 >>> ',(cp_match(6,i),i=1,ncp_match(6))
3052                 print*,'2 >>> ',(cp_match(5,i),i=1,ncp_match(5))                 print*,'2 >>> ',(cp_match(5,i),i=1,ncp_match(5))
# Line 3089  c     $                                 Line 3114  c     $                                
3114  *     **********************************************************  *     **********************************************************
3115  *     ************************** FIT *** FIT *** FIT *** FIT ***  *     ************************** FIT *** FIT *** FIT *** FIT ***
3116  *     **********************************************************  *     **********************************************************
3117    cccc  scommentare se si usa al_ini della nuvola
3118    c$$$                              do i=1,5
3119    c$$$                                 AL(i)=AL_INI(i)
3120    c$$$                              enddo
3121                                  call guess()
3122                                do i=1,5                                do i=1,5
3123                                   AL(i)=AL_INI(i)                                   AL_INI(i)=AL(i)
3124                                enddo                                enddo
3125                                ifail=0 !error flag in chi^2 computation                                ifail=0 !error flag in chi^2 computation
3126                                jstep=0 !number of  minimization steps                                jstep=0 !number of  minimization steps
3127                                iprint=0                                iprint=0
3128    c                              if(DEBUG)iprint=1
3129                                if(DEBUG)iprint=1                                if(DEBUG)iprint=1
3130                                call mini2(jstep,ifail,iprint)                                call mini2(jstep,ifail,iprint)
3131                                if(ifail.ne.0) then                                if(ifail.ne.0) then
3132                                   if(DEBUG)then                                   if(DEBUG)then
3133                                      print *,                                      print *,
3134       $                              '*** MINIMIZATION FAILURE *** '       $                              '*** MINIMIZATION FAILURE *** '
3135       $                              //'(mini2 in clouds_to_ctrack)'       $                              //'(clouds_to_ctrack)'
3136                                        print*,'initial guess: '
3137    
3138                                        print*,'AL_INI(1) = ',AL_INI(1)
3139                                        print*,'AL_INI(2) = ',AL_INI(2)
3140                                        print*,'AL_INI(3) = ',AL_INI(3)
3141                                        print*,'AL_INI(4) = ',AL_INI(4)
3142                                        print*,'AL_INI(5) = ',AL_INI(5)
3143                                   endif                                   endif
3144                                   chi2=-chi2  c                                 chi2=-chi2
3145                                endif                                endif
3146  *     **********************************************************  *     **********************************************************
3147  *     ************************** FIT *** FIT *** FIT *** FIT ***  *     ************************** FIT *** FIT *** FIT *** FIT ***
# Line 3245  c      include 'level1.f' Line 3283  c      include 'level1.f'
3283        call track_init        call track_init
3284        do ip=1,nplanes           !loop on planes        do ip=1,nplanes           !loop on planes
3285    
3286    *     |||||||||||||||||||||||||||||||||||||||||||||||||
3287  *     -------------------------------------------------  *     -------------------------------------------------
3288  *     If the plane has been already included, it just  *     If the plane has been already included, it just
3289  *     computes again the coordinates of the x-y couple  *     computes again the coordinates of the x-y couple
3290  *     using improved PFAs  *     using improved PFAs
3291  *     -------------------------------------------------  *     -------------------------------------------------
3292    *     |||||||||||||||||||||||||||||||||||||||||||||||||
3293           if(XGOOD_STORE(nplanes-ip+1,ibest).eq.1..and.           if(XGOOD_STORE(nplanes-ip+1,ibest).eq.1..and.
3294       $        YGOOD_STORE(nplanes-ip+1,ibest).eq.1. )then       $        YGOOD_STORE(nplanes-ip+1,ibest).eq.1. )then
3295                            
# Line 3284  c            dedxtrk(nplanes-ip+1) = (de Line 3324  c            dedxtrk(nplanes-ip+1) = (de
3324              dedxtrk_x(nplanes-ip+1)=dedx(icx)/mip(VIEW(icx),LADDER(icx)) !(1)(2)              dedxtrk_x(nplanes-ip+1)=dedx(icx)/mip(VIEW(icx),LADDER(icx)) !(1)(2)
3325              dedxtrk_y(nplanes-ip+1)=dedx(icy)/mip(VIEW(icy),LADDER(icy)) !(1)(2)              dedxtrk_y(nplanes-ip+1)=dedx(icy)/mip(VIEW(icy),LADDER(icy)) !(1)(2)
3326                            
3327    *     |||||||||||||||||||||||||||||||||||||||||||||||||
3328  *     -------------------------------------------------  *     -------------------------------------------------
3329  *     If the plane has NOT  been already included,  *     If the plane has NOT  been already included,
3330  *     it tries to include a COUPLE or a single cluster  *     it tries to include a COUPLE or a single cluster
3331  *     -------------------------------------------------  *     -------------------------------------------------
3332    *     |||||||||||||||||||||||||||||||||||||||||||||||||
3333           else                             else                  
3334                                
3335              xgood(nplanes-ip+1)=0              xgood(nplanes-ip+1)=0
# Line 3342  c     $              'ETA2','ETA2', Line 3384  c     $              'ETA2','ETA2',
3384       $              AYV_STORE(nplanes-ip+1,ibest))       $              AYV_STORE(nplanes-ip+1,ibest))
3385                                
3386                 distance = distance_to(XP,YP)                 distance = distance_to(XP,YP)
3387                   distance = distance / RCHI2_STORE(ibest)!<<< MS
3388                 id=id_cp(ip,icp,ist)                 id=id_cp(ip,icp,ist)
3389                 if(DEBUG)print*,'( couple ',id                 if(DEBUG)print*,'( couple ',id
3390       $              ,' ) normalized distance ',distance       $              ,' ) normalized distance ',distance
# Line 3412  c     $              'ETA2','ETA2', Line 3455  c     $              'ETA2','ETA2',
3455       $              PFA,PFA,       $              PFA,PFA,
3456       $              AXV_STORE(nplanes-ip+1,ibest),0.)                     $              AXV_STORE(nplanes-ip+1,ibest),0.)              
3457                 distance = distance_to(XP,YP)                 distance = distance_to(XP,YP)
3458  c     if(DEBUG)print*,'normalized distance ',distance                 distance = distance / RCHI2_STORE(ibest)!<<< MS
3459                 if(DEBUG)print*,'( cl-X ',icx                 if(DEBUG)print*,'( cl-X ',icx
3460       $              ,' in cp ',id,' ) normalized distance ',distance       $              ,' in cp ',id,' ) normalized distance ',distance
3461                 if(distance.lt.distmin)then                 if(distance.lt.distmin)then
# Line 3440  c     $              'ETA2','ETA2', Line 3483  c     $              'ETA2','ETA2',
3483       $              PFA,PFA,       $              PFA,PFA,
3484       $              0.,AYV_STORE(nplanes-ip+1,ibest))       $              0.,AYV_STORE(nplanes-ip+1,ibest))
3485                 distance = distance_to(XP,YP)                 distance = distance_to(XP,YP)
3486                   distance = distance / RCHI2_STORE(ibest)!<<< MS
3487                 if(DEBUG)print*,'( cl-Y ',icy                 if(DEBUG)print*,'( cl-Y ',icy
3488       $              ,' in cp ',id,' ) normalized distance ',distance       $              ,' in cp ',id,' ) normalized distance ',distance
3489                 if(distance.lt.distmin)then                 if(distance.lt.distmin)then
# Line 3479  c     $                 'ETA2','ETA2', Line 3523  c     $                 'ETA2','ETA2',
3523                 endif                 endif
3524    
3525                 distance = distance_to(XP,YP)                 distance = distance_to(XP,YP)
3526                   distance = distance / RCHI2_STORE(ibest)!<<< MS
3527                 if(DEBUG)print*,'( cl-s ',icl                 if(DEBUG)print*,'( cl-s ',icl
3528       $              ,' ) normalized distance ',distance       $              ,' ) normalized distance ',distance
3529                 if(distance.lt.distmin)then                 if(distance.lt.distmin)then
# Line 3508  c                  dedxmm = dedx(icl)   Line 3553  c                  dedxmm = dedx(icl)  
3553                                
3554                 CLS_STORE(nplanes-ip+1,ibest)=iclm !<<<<                     CLS_STORE(nplanes-ip+1,ibest)=iclm !<<<<    
3555  *              ----------------------------  *              ----------------------------
3556    c               print*,'~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~'
3557                 if(mod(VIEW(iclm),2).eq.0)then                 if(mod(VIEW(iclm),2).eq.0)then
3558                    XGOOD(nplanes-ip+1)=1.                    XGOOD(nplanes-ip+1)=1.
3559                    resx(nplanes-ip+1)=rxmm                    resx(nplanes-ip+1)=rxmm
3560                    if(DEBUG)print*,'%%%% included X-cl ',iclm                    if(DEBUG)print*,'%%%% included X-cl ',iclm
3561       $                 ,' ( norm.dist.= ',distmin,', cut ',clinc,' )'  c                  if(.true.)print*,'%%%% included X-cl ',iclm
3562         $                 ,'( chi^2, ',RCHI2_STORE(ibest)
3563         $                 ,', norm.dist.= ',distmin
3564         $                 ,', cut ',clinc,' )'
3565                 else                 else
3566                    YGOOD(nplanes-ip+1)=1.                    YGOOD(nplanes-ip+1)=1.
3567                    resy(nplanes-ip+1)=rymm                    resy(nplanes-ip+1)=rymm
3568                    if(DEBUG)print*,'%%%% included Y-cl ',iclm                    if(DEBUG)print*,'%%%% included Y-cl ',iclm
3569       $                 ,' ( norm.dist.= ',distmin,', cut ',clinc,' )'  c                  if(.true.)print*,'%%%% included Y-cl ',iclm
3570         $                 ,'( chi^2, ',RCHI2_STORE(ibest)
3571         $                 ,', norm.dist.= ', distmin
3572         $                 ,', cut ',clinc,' )'
3573                 endif                 endif
3574    c               print*,'~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~'
3575  *              ----------------------------  *              ----------------------------
3576                 xm_A(nplanes-ip+1) = xmm_A                 xm_A(nplanes-ip+1) = xmm_A
3577                 ym_A(nplanes-ip+1) = ymm_A                 ym_A(nplanes-ip+1) = ymm_A
Legend:
Removed from v.1.10  
changed lines
  Added in v.1.13
  ViewVC Help
Powered by ViewVC 1.1.23