subroutine h_tof_fit(abort,errmsg,trk) *------------------------------------------------------------------- * author: John Arrington * created: 2/22/94 * * h_tof_fit fits the velocity of the paritcle from the corrected * times generated by h_tof. * * modifications: * $Log: h_tof_fit.f,v $ * Revision 1.10 1996/09/04 13:36:24 saw * (JRA) Include actual beta in calculation of focal plane time. * * Revision 1.9 1995/05/22 19:39:29 cdaq * (SAW) Split gen_data_data_structures into gen, hms, sos, and coin parts" * * Revision 1.8 1995/02/23 13:29:15 cdaq * (SAW) Cosmetic Changes * * Revision 1.7 1995/02/10 18:49:57 cdaq * (JRA) Add track index to hgood_scin_time * * Revision 1.6 1994/09/13 21:26:53 cdaq * (JRA) fix chisq calculation * * Revision 1.5 1994/07/13 15:05:08 cdaq * (SAW) Add abs around tmpdenom that I left out last update * * Revision 1.4 1994/07/11 18:34:35 cdaq * (JRA) Increase comparison of tmpdenom from 1e-15 to 1e-10 * * Revision 1.3 1994/07/08 19:42:31 cdaq * (JRA) Change fit from velocity to beta. Bad fits give beta=0 * * Revision 1.2 1994/06/14 04:53:41 cdaq * (DFG) Protect against divide by 0 in beta calc * * Revision 1.1 1994/04/13 16:29:15 cdaq * Initial revision * *------------------------------------------------------------------- implicit none include 'hms_data_structures.cmn' include 'hms_scin_parms.cmn' include 'hms_scin_tof.cmn' logical abort character*1024 errmsg character*20 here parameter (here = 'h_tof_fit') real*4 sumw, sumt, sumz, sumzz, sumtz real*4 scin_weight real*4 tmp, t0 ,tmpdenom real*4 pathnorm integer*4 hit, trk save sumw = 0. sumt = 0. sumz = 0. sumzz = 0. sumtz = 0. do hit = 1 , hscin_tot_hits if (hgood_scin_time(trk,hit)) then scin_weight = 1./hscin_sigma(hit)**2 sumw = sumw + scin_weight sumt = sumt + scin_weight * hscin_time(hit) sumz = sumz + scin_weight * hscin_zpos(hit) sumzz = sumzz + scin_weight * hscin_zpos(hit)**2 sumtz = sumtz + scin_weight * hscin_zpos(hit) * 1 hscin_time(hit) endif enddo * The formula for beta (and t0) come from taking chi-squared (as * defined below), and differentiating with respect to each * of the fit paramters (beta and t0 for fit to z=beta*(t-t0)). * Setting both of these derivatives to zero gives the minumum * chisquared (since they are quadratic in beta and t0), and * gives a solution for beta in terms of sums of z, t, and w. tmp = sumw*sumzz - sumz*sumz t0 = (sumt*sumzz - sumz*sumtz) / tmp tmpdenom = sumw*sumtz - sumz*sumt if(abs(tmpdenom) .gt. 1.e-10) then hbeta(trk) = tmp / tmpdenom !velocity in cm/ns. hbeta_chisq(trk) = 0. do hit = 1 , hscin_tot_hits if (hgood_scin_time(trk,hit)) then hbeta_chisq(trk) = hbeta_chisq(trk) + 1 (hscin_zpos(hit)/hbeta(trk) - 1 (hscin_time(hit) - t0))**2 / hscin_sigma(hit)**2 endif enddo pathnorm = sqrt(1 + hxp_fp(trk)**2 + hyp_fp(trk)**2) hbeta(trk) = hbeta(trk) * pathnorm !take angle into account hbeta(trk) = hbeta(trk) / 29.979 !velocity/c else hbeta(trk) = 0. ! set unphysical beta hbeta_chisq(trk) = -2 endif ! end if on denomimator = 0. return end