/[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.21 by mocchiut, Fri Apr 27 12:11:02 2007 UTC revision 1.25 by pam-fi, Thu May 24 13:29:09 2007 UTC
# Line 578  c     $        rchi2best.lt.15..and. Line 578  c     $        rchi2best.lt.15..and.
578                
579    
580        parameter (ndivx=30)        parameter (ndivx=30)
581    
582    
583    c$$$      print*,icx,icy,sensor,PFAx,PFAy,ax,ay,bfx,bfy
584                
585        resxPAM = 0        resxPAM = 0
586        resyPAM = 0        resyPAM = 0
# Line 609  c      print*,'## xyz_PAM: ',icx,icy,sen Line 612  c      print*,'## xyz_PAM: ',icx,icy,sen
612  *        --------------------------  *        --------------------------
613           angtemp  = ax           angtemp  = ax
614           bfytemp  = bfy           bfytemp  = bfy
615           if(nplx.eq.6) angtemp = -1. * ax  *        /////////////////////////////////
616           if(nplx.eq.6) bfytemp = -1. * bfy  *        AAAAHHHHHHHH!!!!!!!!!!!!!!!!!!!!!
617    *        *grvzkkjsdgjhhhgngbn###>:(
618    *        /////////////////////////////////
619    c         if(nplx.eq.6) angtemp = -1. * ax
620    c         if(nplx.eq.6) bfytemp = -1. * bfy
621             if(viewx.eq.12) angtemp = -1. * ax
622             if(viewx.eq.12) bfytemp = -1. * bfy
623           tgtemp   = tan(angtemp*acos(-1.)/180.) + pmuH_h*bfytemp*0.00001           tgtemp   = tan(angtemp*acos(-1.)/180.) + pmuH_h*bfytemp*0.00001
624           angx     = 180.*atan(tgtemp)/acos(-1.)           angx     = 180.*atan(tgtemp)/acos(-1.)
625           stripx   = stripx - 0.5*pmuH_h*bfytemp*0.00001*SiDimZ/pitchX           stripx   = stripx - 0.5*pmuH_h*bfytemp*0.00001*SiDimZ/pitchX
626  c$$$         print*,nplx,ax,bfy/10.  c$$$         print*,nplx,ax,bfy/10.
627  c$$$         print*,angx,0.5*pmuH_h*bfytemp*0.00001*SiDimZ/pitchX  c$$$         print*,angx,0.5*pmuH_h*bfytemp*0.00001*SiDimZ/pitchX
628  c$$$         print*,'========================'  c$$$         print*,'========================'
629    c$$$         if(bfy.ne.0.)print*,viewx,'-x- '
630    c$$$     $        ,bfy,-1*0.5*pmuH_h*bfytemp*0.00001*SiDimZ
631  *        --------------------------  *        --------------------------
632    
633  c$$$         print*,'--- x-cl ---'  c$$$         print*,'--- x-cl ---'
# Line 656  c$$$         print*,fbad_cog(4,icx) Line 667  c$$$         print*,fbad_cog(4,icx)
667           elseif(PFAx.eq.'ETA2')then           elseif(PFAx.eq.'ETA2')then
668    
669              stripx  = stripx + pfaeta2(icx,angx)                        stripx  = stripx + pfaeta2(icx,angx)          
670              resxPAM = risx_eta2(abs(angx))              resxPAM = risxeta2(abs(angx))
671              resxPAM = resxPAM*fbad_cog(2,icx)              resxPAM = resxPAM*fbad_cog(2,icx)
672              if(DEBUG.and.fbad_cog(2,icx).ne.1)              if(DEBUG.and.fbad_cog(2,icx).ne.1)
673       $           print*,'BAD icx >>> ',viewx,fbad_cog(2,icx)       $           print*,'BAD icx >>> ',viewx,fbad_cog(2,icx)
# Line 664  c$$$         print*,fbad_cog(4,icx) Line 675  c$$$         print*,fbad_cog(4,icx)
675           elseif(PFAx.eq.'ETA3')then                                   elseif(PFAx.eq.'ETA3')then                        
676    
677              stripx  = stripx + pfaeta3(icx,angx)                        stripx  = stripx + pfaeta3(icx,angx)          
678              resxPAM = risx_eta3(abs(angx))                                    resxPAM = risxeta3(abs(angx))                      
679              resxPAM = resxPAM*fbad_cog(3,icx)                            resxPAM = resxPAM*fbad_cog(3,icx)              
680              if(DEBUG.and.fbad_cog(3,icx).ne.1)                          if(DEBUG.and.fbad_cog(3,icx).ne.1)            
681       $           print*,'BAD icx >>> ',viewx,fbad_cog(3,icx)       $           print*,'BAD icx >>> ',viewx,fbad_cog(3,icx)
# Line 672  c$$$         print*,fbad_cog(4,icx) Line 683  c$$$         print*,fbad_cog(4,icx)
683           elseif(PFAx.eq.'ETA4')then                                   elseif(PFAx.eq.'ETA4')then                        
684    
685              stripx  = stripx + pfaeta4(icx,angx)                          stripx  = stripx + pfaeta4(icx,angx)            
686              resxPAM = risx_eta4(abs(angx))                                    resxPAM = risxeta4(abs(angx))                      
687              resxPAM = resxPAM*fbad_cog(4,icx)                            resxPAM = resxPAM*fbad_cog(4,icx)              
688              if(DEBUG.and.fbad_cog(4,icx).ne.1)                            if(DEBUG.and.fbad_cog(4,icx).ne.1)              
689       $           print*,'BAD icx >>> ',viewx,fbad_cog(4,icx)       $           print*,'BAD icx >>> ',viewx,fbad_cog(4,icx)
# Line 680  c$$$         print*,fbad_cog(4,icx) Line 691  c$$$         print*,fbad_cog(4,icx)
691           elseif(PFAx.eq.'ETA')then             elseif(PFAx.eq.'ETA')then  
692    
693              stripx  = stripx + pfaeta(icx,angx)                          stripx  = stripx + pfaeta(icx,angx)            
694              resxPAM = ris_eta(icx,angx)                      c            resxPAM = riseta(icx,angx)                    
695                resxPAM = riseta(viewx,angx)                    
696              resxPAM = resxPAM*fbad_eta(icx,angx)                          resxPAM = resxPAM*fbad_eta(icx,angx)            
697              if(DEBUG.and.fbad_cog(2,icx).ne.1)                            if(DEBUG.and.fbad_cog(2,icx).ne.1)              
698       $           print*,'BAD icx >>> ',viewx,fbad_cog(2,icx)       $           print*,'BAD icx >>> ',viewx,fbad_cog(2,icx)
# Line 734  c$$$            print*,icx,' *** ',resxP Line 746  c$$$            print*,icx,' *** ',resxP
746           tgtemp = tan(ay*acos(-1.)/180.)+pmuH_e*bfx*0.00001                   tgtemp = tan(ay*acos(-1.)/180.)+pmuH_e*bfx*0.00001        
747           angy    = 180.*atan(tgtemp)/acos(-1.)           angy    = 180.*atan(tgtemp)/acos(-1.)
748           stripy = stripy + 0.5*pmuH_e*bfx*0.00001*SiDimZ/pitchY           stripy = stripy + 0.5*pmuH_e*bfx*0.00001*SiDimZ/pitchY
749    c$$$         if(bfx.ne.0.)print*,viewy,'-y- '
750    c$$$     $        ,bfx,0.5*pmuH_e*bfx*0.00001*SiDimZ
751  *        --------------------------  *        --------------------------
752                    
753  c$$$         print*,'--- y-cl ---'  c$$$         print*,'--- y-cl ---'
# Line 772  c$$$         print*,fbad_cog(4,icy) Line 786  c$$$         print*,fbad_cog(4,icy)
786           elseif(PFAy.eq.'ETA2')then           elseif(PFAy.eq.'ETA2')then
787    
788              stripy  = stripy + pfaeta2(icy,angy)                        stripy  = stripy + pfaeta2(icy,angy)          
789              resyPAM = risy_eta2(abs(angy))                            resyPAM = risyeta2(abs(angy))              
790              resyPAM = resyPAM*fbad_cog(2,icy)              resyPAM = resyPAM*fbad_cog(2,icy)
791              if(DEBUG.and.fbad_cog(2,icy).ne.1)              if(DEBUG.and.fbad_cog(2,icy).ne.1)
792       $           print*,'BAD icy >>> ',viewy,fbad_cog(2,icy)       $           print*,'BAD icy >>> ',viewy,fbad_cog(2,icy)
# Line 794  c$$$         print*,fbad_cog(4,icy) Line 808  c$$$         print*,fbad_cog(4,icy)
808           elseif(PFAy.eq.'ETA')then           elseif(PFAy.eq.'ETA')then
809    
810              stripy  = stripy + pfaeta(icy,angy)              stripy  = stripy + pfaeta(icy,angy)
811              resyPAM = ris_eta(icy,angy)    c            resyPAM = riseta(icy,angy)  
812                resyPAM = riseta(viewy,angy)  
813              resyPAM = resyPAM*fbad_eta(icy,angy)              resyPAM = resyPAM*fbad_eta(icy,angy)
814              if(DEBUG.and.fbad_cog(2,icy).ne.1)              if(DEBUG.and.fbad_cog(2,icy).ne.1)
815       $           print*,'BAD icy >>> ',viewy,fbad_cog(2,icy)       $           print*,'BAD icy >>> ',viewy,fbad_cog(2,icy)
# Line 824  c$$$            print*,icy,' *** ',resyP Line 839  c$$$            print*,icy,' *** ',resyP
839    
840        endif        endif
841    
842  c      print*,'## stripx,stripy ',stripx,stripy  c$$$      print*,'## stripx,stripy ',stripx,stripy
843    
844  c===========================================================  c===========================================================
845  C     COUPLE  C     COUPLE
# Line 1040  c      print*,'## xPAM_B,yPAM_B,zPAM_B ' Line 1055  c      print*,'## xPAM_B,yPAM_B,zPAM_B '
1055   100  continue   100  continue
1056        end        end
1057    
1058    ************************************************************************
1059    *     Call xyz_PAM subroutine with default PFA and fill the mini2 common.
1060    *     (done to be called from c/c++)
1061    ************************************************************************
1062    
1063          subroutine xyzpam(ip,icx,icy,lad,sensor,ax,ay,bfx,bfy)
1064    
1065          include 'commontracker.f'
1066          include 'level1.f'
1067          include 'common_mini_2.f'
1068          include 'common_xyzPAM.f'
1069          include 'common_mech.f'
1070          include 'calib.f'
1071          
1072    *     flag to chose PFA
1073    c$$$      character*10 PFA
1074    c$$$      common/FINALPFA/PFA
1075    
1076          integer icx,icy           !X-Y cluster ID
1077          integer sensor
1078          character*4 PFAx,PFAy     !PFA to be used
1079          real ax,ay                !X-Y geometric angle
1080          real bfx,bfy              !X-Y b-field components
1081    
1082          ipx=0
1083          ipy=0      
1084          
1085    c$$$      PFAx = 'COG4'!PFA
1086    c$$$      PFAy = 'COG4'!PFA
1087          
1088          call idtoc(pfaid,PFAx)
1089          call idtoc(pfaid,PFAy)
1090    
1091          call xyz_PAM(icx,icy,sensor,PFAx,PFAy,ax,ay,bfx,bfy)
1092    
1093    c$$$      print*,icx,icy,sensor,PFAx,PFAy,ax,ay,bfx,bfy
1094          
1095          if(icx.ne.0.and.icy.ne.0)then
1096    
1097             ipx=npl(VIEW(icx))
1098             ipy=npl(VIEW(icy))
1099             if( (nplanes-ipx+1).ne.ip.or.(nplanes-ipy+1).ne.ip )
1100         $        print*,'xyzpam: ***WARNING*** clusters ',icx,icy
1101         $        ,' does not belong to the correct plane: ',ip,ipx,ipy
1102    
1103             xgood(ip) = 1.
1104             ygood(ip) = 1.
1105             resx(ip)  = resxPAM
1106             resy(ip)  = resyPAM
1107    
1108             xm(ip) = xPAM
1109             ym(ip) = yPAM
1110             zm(ip) = zPAM
1111             xm_A(ip) = 0.
1112             ym_A(ip) = 0.
1113             xm_B(ip) = 0.
1114             ym_B(ip) = 0.
1115    
1116    c         zv(ip) = zPAM
1117    
1118          elseif(icx.eq.0.and.icy.ne.0)then
1119    
1120             ipy=npl(VIEW(icy))
1121             if((nplanes-ipy+1).ne.ip)
1122         $        print*,'xyzpam: ***WARNING*** clusters ',icx,icy
1123         $        ,' does not belong to the correct plane: ',ip,ipx,ipy
1124    
1125             xgood(ip) = 0.
1126             ygood(ip) = 1.
1127             resx(ip)  = 1000.
1128             resy(ip)  = resyPAM
1129    
1130             xm(ip) = -100.
1131             ym(ip) = -100.
1132             zm(ip) = (zPAM_A+zPAM_B)/2.
1133             xm_A(ip) = xPAM_A
1134             ym_A(ip) = yPAM_A
1135             xm_B(ip) = xPAM_B
1136             ym_B(ip) = yPAM_B
1137    
1138    c         zv(ip) = (zPAM_A+zPAM_B)/2.
1139            
1140          elseif(icx.ne.0.and.icy.eq.0)then
1141    
1142             ipx=npl(VIEW(icx))
1143             if((nplanes-ipx+1).ne.ip)
1144         $        print*,'xyzpam: ***WARNING*** clusters ',icx,icy
1145         $        ,' does not belong to the correct plane: ',ip,ipx,ipy
1146    
1147             xgood(ip) = 1.
1148             ygood(ip) = 0.
1149             resx(ip)  = resxPAM
1150             resy(ip)  = 1000.
1151    
1152             xm(ip) = -100.
1153             ym(ip) = -100.
1154             zm(ip) = (zPAM_A+zPAM_B)/2.
1155             xm_A(ip) = xPAM_A
1156             ym_A(ip) = yPAM_A
1157             xm_B(ip) = xPAM_B
1158             ym_B(ip) = yPAM_B
1159            
1160    c         zv(ip) = (zPAM_A+zPAM_B)/2.
1161    
1162          else
1163    
1164             il = 2
1165             if(lad.ne.0)il=lad
1166             is = 1
1167             if(sensor.ne.0)is=sensor
1168    c         print*,nplanes-ip+1,il,is
1169    
1170             xgood(ip) = 0.
1171             ygood(ip) = 0.
1172             resx(ip)  = 1000.
1173             resy(ip)  = 1000.
1174    
1175             xm(ip) = -100.
1176             ym(ip) = -100.          
1177             zm(ip) = z_mech_sensor(nplanes-ip+1,il,is)*1000./1.d4
1178             xm_A(ip) = 0.
1179             ym_A(ip) = 0.
1180             xm_B(ip) = 0.
1181             ym_B(ip) = 0.
1182    
1183    c         zv(ip) = z_mech_sensor(nplanes-ip+1,il,is)*1000./1.d4
1184    
1185          endif
1186    
1187          if(DEBUG)then
1188    c         print*,'----------------------------- track coord'
1189    22222    format(i2,' * ',3f10.4,' --- ',4f10.4,' --- ',2f4.0,2f10.5)
1190             write(*,22222)ip,zm(ip),xm(ip),ym(ip)
1191         $        ,xm_A(ip),ym_A(ip),xm_B(ip),ym_B(ip)
1192         $        ,xgood(ip),ygood(ip),resx(ip),resy(ip)
1193    c$$$         print*,'-----------------------------'
1194          endif
1195          end
1196    
1197  ********************************************************************************  ********************************************************************************
1198  ********************************************************************************  ********************************************************************************
Legend:
Removed from v.1.21  
changed lines
  Added in v.1.25
  ViewVC Help
Powered by ViewVC 1.1.23