Beta no track added.

src/THcHodoscope.cxx

@@ -259,8 +259,10 @@ Int_t THcHodoscope::ReadDatabase( const TDatime& date )
   //  Int_t plen=strlen(parname);
   cout << " readdatabse hodo fnplanes = " << fNPlanes << endl;
 
+  fBetaNoTrk = 0.;
+  fBetaNoTrkChiSq = 0.;
+  
   fNPaddle = new UInt_t [fNPlanes];
-
   fFPTime = new Double_t [fNPlanes];
   fPlaneCenter = new Double_t[fNPlanes];
   fPlaneSpacing = new Double_t[fNPlanes];
@@ -430,6 +432,8 @@ Int_t THcHodoscope::DefineVariables( EMode mode )
 
   RVarDef vars[] = {
     // Move these into THcHallCSpectrometer using track fTracks
+    {"betanotrack",       "Beta from scintillator hits",                "fBetaNoTrk"},
+    {"betachisqnotrack",  "Chi square of beta from scintillator hits",  "fBetaNoTrkChiSq"},
     {"fpHitsTime",      "Time at focal plane from all hits",         "fFPTime"},
     {"starttime",       "Hodoscope Start Time",                      "fStartTime"},
     {"goodstarttime",   "Hodoscope Good Start Time",                 "fGoodStartTime"},
@@ -504,12 +508,9 @@ void THcHodoscope::DeleteArrays()
 inline 
 void THcHodoscope::ClearEvent()
 {
-  // Reset per-event data.
 
-  //  for ( Int_t imaxhit = 0; imaxhit < MAXHODHITS; imaxhit++ ){
-  //    fBeta[imaxhit] = 0.;
-  //    fBetaChisq[imaxhit] = 0.;
-  //  }
+  fBetaNoTrk = 0.0;
+  fBetaNoTrkChiSq = 0.0;
 
   for(Int_t ip=0;ip<fNPlanes;ip++) {
     fPlanes[ip]->Clear();
@@ -590,6 +591,7 @@ Int_t THcHodoscope::Decode( const THaEvData& evdata )
 //_____________________________________________________________________________
 void THcHodoscope::EstimateFocalPlaneTime( void )
 {
+
   Int_t timehist[200];
 
   for (Int_t i=0;i<200;i++) {
@@ -629,17 +631,27 @@ void THcHodoscope::EstimateFocalPlaneTime( void )
 
   ihit = 0;
   Double_t fpTimeSum = 0.0;
+  Int_t jhit = 0;
   fNfptimes=0;
     
+  fNoTrkPlaneInfo.clear();
+  fNoTrkHitInfo.clear();
   for(Int_t ip=0;ip<fNPlanes;ip++) {
+    fNoTrkPlaneInfo.push_back(NoTrkPlaneInfo());
+    fNoTrkPlaneInfo[ip].goodplanetime = kFALSE;
     Int_t nphits=fPlanes[ip]->GetNScinHits();
     TClonesArray* hodoHits = fPlanes[ip]->GetHits();
     for(Int_t i=0;i<nphits;i++) {
+      fNoTrkHitInfo.push_back(NoTrkHitInfo());
+      fNoTrkHitInfo[jhit].goodtwotimes = kFALSE;
+      fNoTrkHitInfo[jhit].goodscintime = kFALSE;
       Double_t tmin = 0.5*binmax;
       Double_t postime=((THcHodoHit*) hodoHits->At(i))->GetPosTOFCorrectedTime();
       Double_t negtime=((THcHodoHit*) hodoHits->At(i))->GetNegTOFCorrectedTime();
       if ((postime>tmin) && (postime<tmin+fTofTolerance) &&
 	  (negtime>tmin) && (negtime<tmin+fTofTolerance)) {
+	fNoTrkHitInfo[jhit].goodtwotimes = kTRUE;
+	fNoTrkHitInfo[jhit].goodscintime = kTRUE;
 	// Both tubes fired
 	Int_t index=((THcHodoHit*)hodoHits->At(i))->GetPaddleNumber()-1;
 	Double_t fptime = ((THcHodoHit*)hodoHits->At(i))->GetScinCorrectedTime() 
@@ -649,11 +661,14 @@ void THcHodoscope::EstimateFocalPlaneTime( void )
 	  // Should also fill the all FP times histogram
 	  fpTimeSum += fptime;
 	  fNfptimes++;
+	  fNoTrkPlaneInfo[ip].goodplanetime = kTRUE;
 	}
       }
+      jhit++;
     }
     ihit++;
   }
+
   if(fNfptimes>0) {
     fStartTime = fpTimeSum/fNfptimes;
     fGoodStartTime=kTRUE;
@@ -661,6 +676,89 @@ void THcHodoscope::EstimateFocalPlaneTime( void )
     fStartTime = fStartTimeCenter;
     fGoodStartTime=kFALSE;
   }
+
+
+  if ( ( fNoTrkPlaneInfo[0].goodplanetime || fNoTrkPlaneInfo[1].goodplanetime ) &&
+       ( fNoTrkPlaneInfo[2].goodplanetime || fNoTrkPlaneInfo[3].goodplanetime ) ){
+
+    Double_t sumW = 0.;
+    Double_t sumT = 0.;
+    Double_t sumZ = 0.;
+    Double_t sumZZ = 0.;
+    Double_t sumTZ = 0.;    
+    Int_t ihhit = 0;  
+
+    for(Int_t ip=0;ip<fNPlanes;ip++) {
+      Int_t nphits=fPlanes[ip]->GetNScinHits();
+      TClonesArray* hodoHits = fPlanes[ip]->GetHits();
+
+      for(Int_t i=0;i<nphits;i++) {	
+	Int_t index=((THcHodoHit*)hodoHits->At(i))->GetPaddleNumber()-1;
+	    
+	if ( fNoTrkHitInfo[ihhit].goodscintime ) {
+	  
+	  Double_t sigma = 0.5 * ( TMath::Sqrt( TMath::Power( fHodoPosSigma[GetScinIndex(ip,index)],2) + 
+						TMath::Power( fHodoNegSigma[GetScinIndex(ip,index)],2) ) );
+	  Double_t scinWeight = 1 / TMath::Power(sigma,2);
+	  Double_t zPosition = fPlanes[ip]->GetZpos() + (index%2)*fPlanes[ip]->GetDzpos();
+	  
+	  //	  cout << "hit = " << ihhit + 1 << "   zpos = " << zPosition << "   sigma = " << sigma << endl;
+
+	  sumW  += scinWeight;
+	  sumT  += scinWeight * ((THcHodoHit*)hodoHits->At(i))->GetScinCorrectedTime();
+	  sumZ  += scinWeight * zPosition;
+	  sumZZ += scinWeight * ( zPosition * zPosition );
+	  sumTZ += scinWeight * zPosition * ((THcHodoHit*)hodoHits->At(i))->GetScinCorrectedTime();
+	  
+	} // condition of good scin time
+	ihhit ++;
+      } // loop over hits of plane
+    } // loop over planes
+
+    Double_t tmp = sumW * sumZZ - sumZ * sumZ ;
+    Double_t t0 = ( sumT * sumZZ - sumZ * sumTZ ) / tmp ;
+    Double_t tmpDenom = sumW * sumTZ - sumZ * sumT;
+    
+    if ( TMath::Abs( tmpDenom ) > ( 1 / 10000000000.0 ) ) {
+      
+      fBetaNoTrk = tmp / tmpDenom;
+      fBetaNoTrkChiSq = 0.;	  
+      ihhit = 0;
+      
+      for (Int_t ip = 0; ip < fNPlanes; ip++ ){                           // Loop over planes
+	Int_t nphits=fPlanes[ip]->GetNScinHits();
+	TClonesArray* hodoHits = fPlanes[ip]->GetHits();
+
+	for(Int_t i=0;i<nphits;i++) {
+	  Int_t index=((THcHodoHit*)hodoHits->At(i))->GetPaddleNumber()-1;
+	  
+	  if ( fNoTrkHitInfo[ihhit].goodscintime ) {
+	    
+	    Double_t zPosition = fPlanes[ip]->GetZpos() + (index%2)*fPlanes[ip]->GetDzpos();
+	    Double_t timeDif = ( ((THcHodoHit*)hodoHits->At(i))->GetScinCorrectedTime() - t0 );		
+	    Double_t sigma = 0.5 * ( TMath::Sqrt( TMath::Power( fHodoPosSigma[GetScinIndex(ip,index)],2) + 
+						  TMath::Power( fHodoNegSigma[GetScinIndex(ip,index)],2) ) );
+	    fBetaNoTrkChiSq += ( ( zPosition / fBetaNoTrk - timeDif ) *  
+				 ( zPosition / fBetaNoTrk - timeDif ) ) / ( sigma * sigma );
+	    
+	    
+	  } // condition for good scin time
+	  ihhit++;
+	} // loop over hits of a plane
+      } // loop over planes
+	  
+      Double_t pathNorm = 1.0;
+      
+      fBetaNoTrk = fBetaNoTrk * pathNorm;
+      fBetaNoTrk = fBetaNoTrk / 29.979;    // velocity / c	  
+      
+    }  // condition for fTmpDenom	
+    else {
+      fBetaNoTrk = 0.;
+      fBetaNoTrkChiSq = -2.;
+    } // else condition for fTmpDenom
+  }  
+
 }
 //_____________________________________________________________________________
 Int_t THcHodoscope::ApplyCorrections( void )

src/THcHodoscope.h

@@ -69,7 +69,7 @@ public:
   Double_t GetHodoVelLight(Int_t iii) const {return fHodoVelLight[iii];}
   Double_t GetStartTimeCenter() const {return fStartTimeCenter;}
   Double_t GetStartTimeSlop() const {return fStartTimeSlop;}
-  Double_t GetBetaNotrk() const {return fBetaNotrk;}
+  Double_t GetBetaNotrk() const {return fBetaNoTrk;}
 
   Int_t GetGoodRawPad(Int_t iii){return fTOFCalc[iii].good_raw_pad;}
   Int_t GetGoodRawPlane(Int_t iii){return fTOFCalc[iii].pindex;}
@@ -107,7 +107,8 @@ protected:
   Double_t fStartTime; 
   Int_t fNfptimes;
 
-  Double_t fBetaNotrk;
+  Double_t     fBetaNoTrk;
+  Double_t     fBetaNoTrkChiSq;
   // Per-event data
 
   // Potential Hall C parameters.  Mostly here for demonstration
@@ -223,6 +224,22 @@ protected:
   //    Double_t*  gain;
   //  } fDataDest[NDEST];     // Lookup table for decoder
 
+  // Inforamtion for each plane
+  struct NoTrkPlaneInfo {
+    Bool_t goodplanetime;
+    NoTrkPlaneInfo () : goodplanetime(kFALSE) {}
+  };
+  std::vector<NoTrkPlaneInfo> fNoTrkPlaneInfo;
+  
+  // Inforamtion for each plane
+  struct NoTrkHitInfo {
+    Bool_t goodtwotimes;
+    Bool_t goodscintime;
+    NoTrkHitInfo () : goodtwotimes(kFALSE) {}
+  };
+  std::vector<NoTrkHitInfo> fNoTrkHitInfo;
+
+
   // Used in TOF calculation (FineProcess) to hold information about hits
   // within a given plane
   struct TOFPInfo {