(file) Return to h_track_fit.f CVS log (file) (dir) Up to [HallC] / Analyzer / HTRACKING

Diff for /Analyzer/HTRACKING/h_track_fit.f between version 1.1 and 1.7

version 1.1, 1994/02/19 06:20:53 version 1.7, 1995/01/27 20:26:50
Line 9 
Line 9 
 *                              remove minuit. Make fit linear *                              remove minuit. Make fit linear
 *                              still does not do errors properly *                              still does not do errors properly
 * $Log$ * $Log$
   * Revision 1.7  1995/01/27 20:26:50  cdaq
   * (JRA) Remove Mack's personal focalplane diceamatic (z slicer) code
   *
   * Revision 1.6  1994/12/06  15:45:27  cdaq
   * (DJM) Take slices in Z to look for best focus
   *
   * Revision 1.5  1994/10/12  18:52:06  cdaq
   * (DJM) Initialize some variables
   * (SAW) Prettify indentation
   *
   * Revision 1.4  1994/09/01  12:29:07  cdaq
   * (DJM) Make registered versions of residuals
   *
   * Revision 1.3  1994/08/18  02:45:53  cdaq
   * (DM) Add calculation of residuals
   *
   * Revision 1.2  1994/02/22  05:25:00  cdaq
   * (SAW) Remove dfloat calls with floating args
   *
 * Revision 1.1  1994/02/19 06:20:53  cdaq * Revision 1.1  1994/02/19 06:20:53  cdaq
 * Initial revision * Initial revision
 * *
Line 17 
Line 36 
       include "gen_data_structures.cmn"       include "gen_data_structures.cmn"
       include "hms_tracking.cmn"       include "hms_tracking.cmn"
       include "hms_geometry.cmn"       include "hms_geometry.cmn"
         external H_DPSIFUN
 * *
 * *
 *     local variables *     local variables
Line 28 
Line 48 
       integer*4 itrack                        ! track loop index       integer*4 itrack                        ! track loop index
       integer*4 ihit,ierr       integer*4 ihit,ierr
       integer*4 hit,plane       integer*4 hit,plane
       integer*4 i,j                             ! loop index        integer*4 i,j,k                             ! loop index
   *      real*4 z_slice
   
         real*8   H_DPSIFUN
         real*8   pos
       real*4   ray(hnum_fpray_param)       real*4   ray(hnum_fpray_param)
         real*8   ray1(4)
         real*8   ray2(4)
       real*8 TT(hnum_fpray_param)       real*8 TT(hnum_fpray_param)
       real*8 AA(hnum_fpray_param,hnum_fpray_param)       real*8 AA(hnum_fpray_param,hnum_fpray_param)
       real*8 dray(hnum_fpray_param)       real*8 dray(hnum_fpray_param)
Line 43 
Line 69 
 * *
       ABORT= .FALSE.       ABORT= .FALSE.
       ierr=0       ierr=0
   *  initailize residuals
   
         do itrack=1,HNTRACKS_MAX
           do plane=1,HMAX_NUM_DC_PLANES
             hdc_residual(itrack,plane)=1000
             hdc_sing_res(itrack,plane)=1000
             hdc1_sing_res(plane)=1000
             hdc2_sing_res(plane)=1000
             hdc1_dbl_res(plane)=1000
             hdc2_dbl_res(plane)=1000
   
           enddo
         enddo
   
   
 *     test for no tracks *     test for no tracks
       if(HNTRACKS_FP.ge.1) then       if(HNTRACKS_FP.ge.1) then
         do itrack=1,HNTRACKS_FP         do itrack=1,HNTRACKS_FP
           chi2= dummychi2           chi2= dummychi2
           htrack_fit_num=itrack           htrack_fit_num=itrack
   
 *     are there enough degrees of freedom *     are there enough degrees of freedom
           HNFREE_FP(itrack)=HNTRACK_HITS(itrack,1)-hnum_fpray_param           HNFREE_FP(itrack)=HNTRACK_HITS(itrack,1)-hnum_fpray_param
           if(HNFREE_FP(itrack).gt.0) then           if(HNFREE_FP(itrack).gt.0) then
   
 *     initialize parameters *     initialize parameters
            do i=1,hnum_fpray_param            do i=1,hnum_fpray_param
             TT(i)=0.             TT(i)=0.
             do ihit=2,HNTRACK_HITS(itrack,1)+1             do ihit=2,HNTRACK_HITS(itrack,1)+1
               hit=HNTRACK_HITS(itrack,ihit)               hit=HNTRACK_HITS(itrack,ihit)
               plane=HDC_PLANE_NUM(hit)               plane=HDC_PLANE_NUM(hit)
               TT(i)=TT(i)+DFLOAT((HDC_WIRE_COORD(hit)*                  TT(i)=TT(i)+((HDC_WIRE_COORD(hit)*
      &          hplane_coeff(remap(i),plane))      &          hplane_coeff(remap(i),plane))
      &          /(hdc_sigma(plane)*hdc_sigma(plane)))      &          /(hdc_sigma(plane)*hdc_sigma(plane)))
             enddo             enddo
Line 65 
Line 108 
           do i=1,hnum_fpray_param           do i=1,hnum_fpray_param
            do j=1,hnum_fpray_param            do j=1,hnum_fpray_param
            AA(i,j)=0.            AA(i,j)=0.
                   if(j.lt.i)then
                     AA(i,j)=AA(j,i)
                   else
             do ihit=2,HNTRACK_HITS(itrack,1)+1             do ihit=2,HNTRACK_HITS(itrack,1)+1
               hit=HNTRACK_HITS(itrack,ihit)               hit=HNTRACK_HITS(itrack,ihit)
               plane=HDC_PLANE_NUM(hit)               plane=HDC_PLANE_NUM(hit)
               AA(i,j)=AA(i,j) + DFLOAT(                      AA(i,j)=AA(i,j) + (
      &           hplane_coeff(remap(i),plane)*hplane_coeff(remap(j),plane)       &                   hplane_coeff(remap(i),plane)*hplane_coeff(remap(j)
      &          /(hdc_sigma(plane)*hdc_sigma(plane)))       $                   ,plane)/(hdc_sigma(plane)*hdc_sigma(plane)))
             enddo               ! end loop on ihit             enddo               ! end loop on ihit
                   endif                   ! end test on j .lt. i
            enddo                ! end loop on j            enddo                ! end loop on j
           enddo                 ! end loop on i           enddo                 ! end loop on i
 * *
Line 86 
Line 133 
       else       else
 *      calculate chi2 *      calculate chi2
         chi2=0.         chi2=0.
   
         ray(1)=dray(1)         ray(1)=dray(1)
         ray(2)=dray(2)         ray(2)=dray(2)
         ray(3)=dray(3)         ray(3)=dray(3)
Line 93 
Line 141 
          do ihit=2,HNTRACK_HITS(itrack,1)+1          do ihit=2,HNTRACK_HITS(itrack,1)+1
               hit=HNTRACK_HITS(itrack,ihit)               hit=HNTRACK_HITS(itrack,ihit)
               plane=HDC_PLANE_NUM(hit)               plane=HDC_PLANE_NUM(hit)
   
 * note chi2 is single precision * note chi2 is single precision
            chi2=chi2+                  hdc_sing_res(itrack,plane)=
      &         (HDC_WIRE_COORD(hit)      &         (HDC_WIRE_COORD(hit)
      &         -hplane_coeff(remap(1),plane)*ray(1)      &         -hplane_coeff(remap(1),plane)*ray(1)
      &         -hplane_coeff(remap(2),plane)*ray(2)      &         -hplane_coeff(remap(2),plane)*ray(2)
      &         -hplane_coeff(remap(3),plane)*ray(3)      &         -hplane_coeff(remap(3),plane)*ray(3)
      &         -hplane_coeff(remap(4),plane)*ray(4))**2       &               -hplane_coeff(remap(4),plane)*ray(4))
      &         /(hdc_sigma(plane)*hdc_sigma(plane))                  chi2=chi2+((hdc_sing_res(itrack,plane))**2
        &               )/(hdc_sigma(plane)*hdc_sigma(plane))
   *djm
                   if(itrack.eq.1 .and. plane.ge.1 .and. plane.le.6)
        &               hdc1_sing_res(plane) = hdc_sing_res(itrack,plane)
                   if(itrack.eq.2 .and. plane.ge.7 .and. plane.le.12)
        &               hdc2_sing_res(plane) = hdc_sing_res(itrack,plane)
         enddo         enddo
       endif       endif
  
Line 113 
Line 168 
           HCHI2_FP(itrack)=chi2           HCHI2_FP(itrack)=chi2
         enddo                                 ! end loop over tracks         enddo                                 ! end loop over tracks
       endif       endif
   
   **     A reasonable selection of slices is presently -80,-60,-40,-20,0,20,40
   **     ,60,80.   Zero is the nominal midplane between the chambers, -80
   **     corresponds closely to the exit flange position.  This slice pattern
   **     is created with hz_wild=-80., hdelta_z_wild=20., and hnum_zslice=9
   *      if(hz_slice_enable.ne.0)then
   *        do k=1,hnum_zslice
   *          z_slice = hz_wild + (k-1)*hdelta_z_wild
   *          hx_fp_wild(k) = hx_fp(1) + hxp_fp(1)*z_slice
   *          hy_fp_wild(k) = hy_fp(1) + hyp_fp(1)*z_slice
   *        enddo
   *      endif
   
   * calculate residuals for each chamber if in single stub mode
   * and there were 2 tracks found one in first chanber and one in the second
   
         if (hsingle_stub.ne.0) then
           if (HNTRACKS_FP.eq.2) then
             itrack=1
             ihit=2
             hit=HNTRACK_HITS(itrack,ihit)
             plane=HDC_PLANE_NUM(hit)
             if (plane.le.6) then
               itrack=2
               hit=HNTRACK_HITS(itrack,ihit)
               plane=HDC_PLANE_NUM(hit)
               if (plane.ge.7) then
   
   * condition of above met calculating residuals
   * assigning rays to tracks in each chamber
   * ray1 is ray from first chamber fit
   * ray2 is ray from second chamber fit
   
                 ray1(1)=dble(HX_FP(1))
                 ray1(2)=dble(HY_FP(1))
                 ray1(3)=dble(HXP_FP(1))
                 ray1(4)=dble(HYP_FP(1))
                 ray2(1)=dble(HX_FP(2))
                 ray2(2)=dble(HY_FP(2))
                 ray2(3)=dble(HXP_FP(2))
                 ray2(4)=dble(HYP_FP(2))
   
                 itrack=1
   * loop over hits in second chamber
                 do ihit=1,HNTRACK_HITS(itrack+1,1)
   
   * calculate residual in second chamber from first chamber track
                   hit=HNTRACK_HITS(itrack+1,ihit+1)
                   plane=HDC_PLANE_NUM(hit)
                   pos=H_DPSIFUN(ray1,plane)
                   hdc_residual(itrack,plane)=HDC_WIRE_COORD(hit)-pos
   * djm 8/31/94 stuff this variable into 1d array we can register
                   hdc2_dbl_res(plane) = hdc_residual(1,plane)
   
                 enddo
   
                 itrack=2
   * loop over hits in first chamber
                 do ihit=1,HNTRACK_HITS(itrack-1,1)
   
   * calculate residual in first chamber from second chamber track
                   hit=HNTRACK_HITS(itrack-1,ihit+1)
                   plane=HDC_PLANE_NUM(hit)
                   pos=H_DPSIFUN(ray2,plane)
                   hdc_residual(itrack,plane)=HDC_WIRE_COORD(hit)-pos
   * djm 8/31/94 stuff this variable into 1d array we can register
                   hdc1_dbl_res(plane) = hdc_residual(2,plane)
   
                 enddo
               endif                       ! end plane ge 7
             endif                         ! end plane le 6
           endif                           ! end HNTRACKS_FP eq 2
         endif                             ! end hsignle_stub .ne. 0
   
 *     test if we want to dump out trackfit results *     test if we want to dump out trackfit results
       if(hdebugtrackprint.ne.0) then       if(hdebugtrackprint.ne.0) then
          call h_print_tracks          call h_print_tracks
Line 121 
Line 250 
       end       end
  
 * *
   

Legend:
Removed from v.1.1  
changed lines
  Added in v.1.7
Analyzer/Replay: Mark Jones, Documents: Stephen Wood
 Powered by
ViewCVS 0.9.2-cvsgraph-1.4.0