From 201bcda27107be6190d13a57ed00682536fb2859 Mon Sep 17 00:00:00 2001
From: Cdaq Account <cdaq@cdaql1.jlab.org>
Date: Mon, 6 Mar 2017 08:38:30 -0500
Subject: [PATCH] Modified hodoscope and scintillatorplane code for hodo
 calibration Temporarily set THcHodoscope.cxx so that only the first three
 planes are used to calculate the focal plane time and beta. In THcHodoscope ,
 Modified format for writing out the dump file used in hodo calibration. Add
 cosmicflag to THcScintillatorPLane so the corrected time is done right for
 cosmics. Modified hc_hodo_calib/tofcal.f so that it uses detector 10 as the
 reference and that it writes out better defualts for parameters in there are
 not enough events for a PMT.

---
 hc_hodo_calib/tofcal.f       | 16 +++----
 src/THcHitList.cxx           |  6 +--
 src/THcHodoscope.cxx         | 84 ++++++++++++++++++++----------------
 src/THcHodoscope.h           |  4 +-
 src/THcScintillatorPlane.cxx | 24 ++++++++---
 src/THcScintillatorPlane.h   |  1 +
 6 files changed, 80 insertions(+), 55 deletions(-)

diff --git a/hc_hodo_calib/tofcal.f b/hc_hodo_calib/tofcal.f
index 39c7917..5fe197d 100644
--- a/hc_hodo_calib/tofcal.f
+++ b/hc_hodo_calib/tofcal.f
@@ -50,7 +50,7 @@
       open(unit=9,file=filename,err=998)
        if(filename(1:1).ne.'h'.and.filename(1:1).ne.'s') then
         write(6,'(1x,''error, input file name '',
-     >   '' should start with h or s'')')
+     >   '' sould start with h or s'')')
         goto 999
        endif
        write(*,*) ' filling adc histograms'
@@ -102,7 +102,7 @@ c         write(13,'(1x,''adcmin='',f6.1)') adcmin(i)
 
 ! Do everything twice, 2nd time with tighter tolerance (ttol) on
 ! time differences, based on first pass
-      ttol=5.0
+      ttol=20.0
       do iloop=1,2
          write(*,*) ' loop = ',iloop
 ! Initialize the fitting arrays
@@ -128,7 +128,6 @@ c         write(13,'(1x,''adcmin='',f6.1)') adcmin(i)
           n=n+1
           read(string,'(1x,i1,2i3,5f10.3)') ipn,ipl,idt, 
      >      tr(n),p(n),zc(n),tc1(n),adc(n)
-
 ! linearize the detector numbers
           idet(n) = 100*(ipn-1) + 20*(ipl-1) + idt
           if(idet(n).lt.1.or.idet(n).gt.200) write(6,'(1x,
@@ -149,7 +148,6 @@ c          adc(n) = min(500, max(0., adc(n)))
 
 ! Loop over all pairs, if at least 6
  11     if(n.ge.6) then
-
 ! see if this is first time a detector is used
           do j=1,n
             nhit(idet(j))=nhit(idet(j))+1
@@ -158,7 +156,7 @@ c          adc(n) = min(500, max(0., adc(n)))
 ! Note that detector had has a fixed time offset (ip1) of zero
 ! since all times are relative)
             if(nhit(idet(j)).eq.1) then
-              if(idet(j).eq.4) then
+              if(idet(j).eq.10) then
                 ip1(idet(j))=0
               else
 
@@ -334,17 +332,17 @@ c          adc(n) = min(500, max(0., adc(n)))
       enddo
 
       write(10+iloop,'(/a,''hodo_pos_invadc_linear ='',3(f8.2,'',''),
-     >  f8.2)')filename(1:1),( -1./min(-0.02,vel(i)),i=1,80,20)
+     >  f8.2)')filename(1:1),( -1./min(-1./15.,vel(i)),i=1,80,20)
       do j=2,16
        write(10+iloop,'(1x,''                        '',3(f8.2,'',''),
-     >  f8.2)')(-1./min(-0.02,vel(i)),i=j,79+j,20)
+     >  f8.2)')(-1./min(-1./15.,vel(i)),i=j,79+j,20)
       enddo
 
       write(10+iloop,'(/a,''hodo_neg_invadc_linear ='',3(f8.2,'',''),
-     >  f8.2)')filename(1:1),( -1./min(-0.02,vel(i)),i=101,180,20)
+     >  f8.2)')filename(1:1),( -1./min(-1./15.,vel(i)),i=101,180,20)
       do j=2,16
        write(10+iloop,'(1x,''                        '',3(f8.2,'',''),
-     >  f8.2)')(-1./min(-0.02,vel(i)),i=100+j,179+j,20)
+     >  f8.2)')(-1./min(-1./15.,vel(i)),i=100+j,179+j,20)
       enddo
 
       write(10+iloop,'(/a,''hodo_pos_invadc_adc='',3(f9.2,'',''),
diff --git a/src/THcHitList.cxx b/src/THcHitList.cxx
index eafad16..41f7815 100644
--- a/src/THcHitList.cxx
+++ b/src/THcHitList.cxx
@@ -203,7 +203,7 @@ Int_t THcHitList::DecodeToHitList( const THaEvData& evdata ) {
 	    Int_t reftime = evdata.GetData(d->crate, d->slot, d->refchan, 0);
 	    rawhit->SetReference(signal, reftime);
 	  } else {
-	    cout << "HitList: refchan " << d->refchan <<
+	    cout << "HitList(event=" << evdata.GetEvNum() << "): refchan " << d->refchan <<
 	      " missing for (" << d->crate << ", " << d->slot <<
 	      ", " << chan << ")" << endl;
 	  }
@@ -212,7 +212,7 @@ Int_t THcHitList::DecodeToHitList( const THaEvData& evdata ) {
 	    if(fRefIndexMaps[d->refindex].hashit) {
 	      rawhit->SetReference(signal, fRefIndexMaps[d->refindex].reftime);
 	    } else {
-	      cout << "HitList: refindex " << d->refindex <<
+	      cout << "HitList(event=" << evdata.GetEvNum() << "): refindex " << d->refindex <<
           " (" << fRefIndexMaps[d->refindex].crate <<
           ", " << fRefIndexMaps[d->refindex].slot <<
           ", " << fRefIndexMaps[d->refindex].channel << ")" <<
@@ -223,7 +223,7 @@ Int_t THcHitList::DecodeToHitList( const THaEvData& evdata ) {
 	}
       } else {			// This is a Flash ADC
 
-        if (fPSE125 && fPSE125->IsPresent(d->crate)) {  // Set F250 parameters.
+        if (fPSE125) {  // Set F250 parameters.
           rawhit->SetF250Params(
             fPSE125->GetNSA(d->crate),
             fPSE125->GetNSB(d->crate),
diff --git a/src/THcHodoscope.cxx b/src/THcHodoscope.cxx
index cb4f531..d82373e 100644
--- a/src/THcHodoscope.cxx
+++ b/src/THcHodoscope.cxx
@@ -40,6 +40,7 @@ hodoscope array, not just one plane.
 #include <cstdio>
 #include <cstdlib>
 #include <iostream>
+#include <iomanip>
 #include <fstream>
 
 using namespace std;
@@ -85,6 +86,10 @@ void THcHodoscope::Setup(const char* name, const char* description)
   prefix[0]=tolower(GetApparatus()->GetName()[0]);
   prefix[1]='\0';
 
+  TString temp(prefix[0]);
+  fSHMS=kFALSE;
+  if (temp == "p" ) fSHMS=kTRUE;
+  cout << " fSHMS = " << fSHMS << endl;
   string planenamelist;
   DBRequest listextra[]={
     {"hodo_num_planes", &fNPlanes, kInt},
@@ -360,9 +365,9 @@ Int_t THcHodoscope::ReadDatabase( const TDatime& date )
   if (fCosmicFlag==1) cout << " setup for cosmics in TOF"<< endl;
   cout << " cosmic flag = " << fCosmicFlag << endl;
   if(fDumpTOF) {
-    fDumpOut.open(fTOFDumpFile.c_str());
-    if(fDumpOut.is_open()) {
-      fDumpOut << "Hodoscope Time of Flight calibration data" << endl;
+      fDumpOut.open(fTOFDumpFile.c_str());
+   if(fDumpOut.is_open()) {
+     //fDumpOut << "Hodoscope Time of Flight calibration data" << endl;
     } else {
       fDumpTOF = 0;
       cout << "WARNING: Unable to open TOF Dump file " << fTOFDumpFile << endl;
@@ -561,7 +566,7 @@ void THcHodoscope::ClearEvent()
     fFPTime[ip]=0.;
     fPlaneCenter[ip]=0.;
     fPlaneSpacing[ip]=0.;
-    for(Int_t iPaddle=0;iPaddle<fNPaddle[ip]; ++iPaddle) {
+    for(UInt_t iPaddle=0;iPaddle<fNPaddle[ip]; ++iPaddle) {
       fScinHitPaddle[ip][iPaddle]=0;
     }
   }
@@ -706,7 +711,9 @@ void THcHodoscope::EstimateFocalPlaneTime( void )
 
   Bool_t goodplanetime[fNPlanes];
   Bool_t twogoodtimes[nscinhits];
-  for(Int_t ip=0;ip<fNPlanes;ip++) {
+  Int_t temp_planes=fNPlanes;
+  if (fSHMS) temp_planes=3;
+  for(Int_t ip=0;ip<temp_planes;ip++) {
     goodplanetime[ip] = kFALSE;
     Int_t nphits=fPlanes[ip]->GetNScinHits();
     TClonesArray* hodoHits = fPlanes[ip]->GetHits();
@@ -722,10 +729,8 @@ void THcHodoscope::EstimateFocalPlaneTime( void )
       if ((postime>tmin) && (postime<tmin+fTofTolerance) &&
 	  (negtime>tmin) && (negtime<tmin+fTofTolerance)) {
 	hit->SetTwoGoodTimes(kTRUE);
-	twogoodtimes[ihit] = kTRUE;
-	// Both tubes fired
-	Int_t index=hit->GetPaddleNumber()-1;
-	// Need to put this in a multihit histo
+	twogoodtimes[ihit] = kTRUE;	// Both tubes fired
+        Int_t index=hit->GetPaddleNumber()-1;	 //
         Double_t fptime;
         if(fCosmicFlag==1) {
 	   fptime = hit->GetScinCorrectedTime()
@@ -772,7 +777,9 @@ void THcHodoscope::EstimateFocalPlaneTime( void )
     Double_t sumTZ = 0.;
     Int_t ihhit = 0;
 
-    for(Int_t ip=0;ip<fNPlanes;ip++) {
+    temp_planes=fNPlanes;
+  if (fSHMS) temp_planes=3;
+    for(Int_t ip=0;ip<temp_planes;ip++) {
       Int_t nphits=fPlanes[ip]->GetNScinHits();
       TClonesArray* hodoHits = fPlanes[ip]->GetHits();
 
@@ -785,8 +792,7 @@ void THcHodoscope::EstimateFocalPlaneTime( void )
 						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;
+ 	  //	  cout << "hit = " << ihhit + 1 << "   zpos = " << zPosition << "   sigma = " << sigma << endl;
 
 	  sumW  += scinWeight;
 	  sumT  += scinWeight * ((THcHodoHit*)hodoHits->At(i))->GetScinCorrectedTime();
@@ -809,7 +815,7 @@ void THcHodoscope::EstimateFocalPlaneTime( void )
       fBetaNoTrkChiSq = 0.;
       ihhit = 0;
 
-      for (Int_t ip = 0; ip < fNPlanes; ip++ ){                           // Loop over planes
+      for (Int_t ip = 0; ip < temp_planes; ip++ ){                           // Loop over planes
 	Int_t nphits=fPlanes[ip]->GetNScinHits();
 	TClonesArray* hodoHits = fPlanes[ip]->GetHits();
 
@@ -842,6 +848,9 @@ void THcHodoscope::EstimateFocalPlaneTime( void )
       fBetaNoTrkChiSq = -2.;
     } // else condition for fTmpDenom
   }
+   fGoodEventTOFCalib=kFALSE;
+  if (!fSHMS&&goodplanetime[0]&&goodplanetime[1]&&goodplanetime[2]&&goodplanetime[3]&&fPlanes[0]->GetNGoodHits()==1&&fPlanes[1]->GetNGoodHits()==1&&fPlanes[2]->GetNGoodHits()==1&&fPlanes[3]->GetNGoodHits()==1) fGoodEventTOFCalib=kTRUE;
+  if (fSHMS&&goodplanetime[0]&&goodplanetime[1]&&goodplanetime[2]&&fPlanes[0]->GetNGoodHits()==1&&fPlanes[1]->GetNGoodHits()==1&&fPlanes[2]->GetNGoodHits()==1) fGoodEventTOFCalib=kTRUE;
 
 }
 //_____________________________________________________________________________
@@ -874,6 +883,8 @@ Int_t THcHodoscope::FineProcess( TClonesArray& tracks )
   // -------------------------------------------------
 
   //  fDumpOut << " ntrack =  " << ntracks  << endl;
+  Int_t temp_planes=fNPlanes;
+  if (fSHMS) temp_planes=3;
 
   if (tracks.GetLast()+1 > 0 ) {
 
@@ -887,7 +898,7 @@ Int_t THcHodoscope::FineProcess( TClonesArray& tracks )
       THaTrack* theTrack = dynamic_cast<THaTrack*>( tracks.At(itrack) );
       if (!theTrack) return -1;
 
-      for (Int_t ip = 0; ip < fNPlanes; ip++ ){
+      for (Int_t ip = 0; ip < temp_planes; ip++ ){
 	fGoodPlaneTime[ip] = kFALSE;
 	fNScinHits[ip] = 0;
 	fNPlaneTime[ip] = 0;
@@ -928,7 +939,7 @@ Int_t THcHodoscope::FineProcess( TClonesArray& tracks )
       fTOFPInfo.clear();  // SAW - combine these two?
       Int_t ihhit = 0;		// Hit # overall
 
-      for(Int_t ip = 0; ip < fNPlanes; ip++ ) {
+     for(Int_t ip = 0; ip < temp_planes; ip++ ) {
 
 	std::vector<GoodFlags> goodflagstmp2;
 	fGoodFlags[itrack].push_back(goodflagstmp2);
@@ -989,12 +1000,16 @@ Int_t THcHodoscope::FineProcess( TClonesArray& tracks )
 	       ( fPlanes[ip]->GetSize() * 0.5 + fPlanes[ip]->GetHodoSlop() ) ){ // Line 293
 
 	    fTOFPInfo[ihhit].onTrack = kTRUE;
-
 	    Double_t zcor = zposition/(29.979*fBetaP)*
 		TMath::Sqrt(1. + theTrack->GetTheta()*theTrack->GetTheta()
 			    + theTrack->GetPhi()*theTrack->GetPhi());
 	    fTOFPInfo[ihhit].zcor = zcor;
-
+	    if (fCosmicFlag) {
+	    Double_t zcor = -zposition/(29.979*1.0)*
+		TMath::Sqrt(1. + theTrack->GetTheta()*theTrack->GetTheta()
+			    + theTrack->GetPhi()*theTrack->GetPhi());
+	    fTOFPInfo[ihhit].zcor = zcor;
+	    }
 	    Double_t tdc_pos = hit->GetPosTDC();
 	    if(tdc_pos >=fScinTdcMin && tdc_pos <= fScinTdcMax ) {
 	      Double_t adc_pos = hit->GetPosADC();
@@ -1013,7 +1028,7 @@ Int_t THcHodoscope::FineProcess( TClonesArray& tracks )
 		  + fHodoPosTimeOffset[fPIndex];
 	      }
 	      fTOFPInfo[ihhit].scin_pos_time = timep;
-	      timep -= zcor;
+ 	      timep -= zcor;
 	      fTOFPInfo[ihhit].time_pos = timep;
 
 	      for ( Int_t k = 0; k < 200; k++ ){ // Line 211
@@ -1122,21 +1137,17 @@ Int_t THcHodoscope::FineProcess( TClonesArray& tracks )
 	  if ( fTOFPInfo[ih].keep_pos ) { // 301
 	    fTOFCalc[ih].good_tdc_pos = kTRUE;
 	    fGoodFlags[itrack][ip][iphit].goodTdcPos = kTRUE;
-	    if(fDumpTOF && ntracks==1) {
-	      fDumpOut << "1 " << ip+1 << " " << paddle+1 << " " <<
-		hit->GetPosTDC()*fScinTdcToTime << " " << fTOFPInfo[ih].pathp <<
-		" "  << fTOFPInfo[ih].zcor << " " << fTOFPInfo[ih].time_pos <<
-		" " << hit->GetPosADC() << endl;
+	    if(fDumpTOF && ntracks==1 && fGoodEventTOFCalib) {
+              fDumpOut << fixed << setprecision(2);
+             fDumpOut  << showpoint << " 1" << setw(3) << ip+1 << setw(3) << hit->GetPaddleNumber()  << setw(10) << hit->GetPosTDC()*fScinTdcToTime  << setw(10) << fTOFPInfo[ih].pathp << setw(10) << fTOFPInfo[ih].zcor  << setw(10) << fTOFPInfo[ih].time_pos << setw(10) << hit->GetPosADC() << endl;
 	    }
 	  }
 	  if ( fTOFPInfo[ih].keep_neg ) { //
 	    fTOFCalc[ih].good_tdc_neg = kTRUE;
 	    fGoodFlags[itrack][ip][iphit].goodTdcNeg = kTRUE;
-	    if(fDumpTOF && ntracks==1) {
-	      fDumpOut << "2 " << ip+1 << " " << paddle+1 << " " <<
-		hit->GetNegTDC()*fScinTdcToTime << " " << fTOFPInfo[ih].pathn <<
-		" "  << fTOFPInfo[ih].zcor << " " << fTOFPInfo[ih].time_neg <<
-		" " << hit->GetNegADC() << endl;
+	    if(fDumpTOF && ntracks==1 && fGoodEventTOFCalib) {
+             fDumpOut << fixed << setprecision(2);
+             fDumpOut  << showpoint << " 2" << setw(3) << ip+1 << setw(3) << hit->GetPaddleNumber() << setw(10) << hit->GetNegTDC()*fScinTdcToTime  << setw(10) << fTOFPInfo[ih].pathn << setw(10) << fTOFPInfo[ih].zcor  << setw(10) << fTOFPInfo[ih].time_neg << setw(10) << hit->GetNegADC() << endl;
 	    }
 	  }
 
@@ -1225,11 +1236,10 @@ Int_t THcHodoscope::FineProcess( TClonesArray& tracks )
 	}
 
       } // Second loop over hits of a scintillator plane ends here
-
-      theTrack->SetGoodPlane3( fGoodPlaneTime[2] ? 1 : 0 );
+       theTrack->SetGoodPlane3( fGoodPlaneTime[2] ? 1 : 0 );
       theTrack->SetGoodPlane4( fGoodPlaneTime[3] ? 1 : 0 );
-
-      //------------------------------------------------------------------------------
+     //
+       //------------------------------------------------------------------------------
       //------------------------------------------------------------------------------
       //------------------------------------------------------------------------------
       //------------------------------------------------------------------------------
@@ -1317,7 +1327,7 @@ Int_t THcHodoscope::FineProcess( TClonesArray& tracks )
 
       Double_t FPTimeSum=0.0;
       Int_t nFPTimeSum=0;
-      for (Int_t ip = 0; ip < fNPlanes; ip++ ){
+      for (Int_t ip = 0; ip < temp_planes; ip++ ){
 	if ( fNPlaneTime[ip] != 0 ){
 	  fFPTime[ip] = ( fSumPlaneTime[ip] / fNPlaneTime[ip] );
 	  FPTimeSum += fSumPlaneTime[ip];
@@ -1348,8 +1358,8 @@ Int_t THcHodoscope::FineProcess( TClonesArray& tracks )
 
   } // If condition for at least one track
 
-   if(fDumpTOF && ntracks==1 ) {
-	      fDumpOut << "0 "  << endl;
+   if(fDumpTOF && ntracks==1 && fGoodEventTOFCalib) {
+         fDumpOut << "0 "  << endl;
 	    }
 
   //-----------------------------------------------------------------------
@@ -1362,7 +1372,7 @@ Int_t THcHodoscope::FineProcess( TClonesArray& tracks )
   //  *second, we move the scintillators.  here we use scintillator cuts to see
   //  *if a track should have been found.
 
-  for(Int_t ip = 0; ip < fNPlanes; ip++ ) {
+  for(Int_t ip = 0; ip < temp_planes; ip++ ) {
     if (!fPlanes[ip])
       return -1;
 
@@ -1383,7 +1393,7 @@ Int_t THcHodoscope::FineProcess( TClonesArray& tracks )
   //  *Wwe count the number of three adjacent scintillators too.  (A signle track
   //  *shouldn't fire three adjacent scintillators.
 
-  for(Int_t ip = 0; ip < fNPlanes; ip++ ) {
+  for(Int_t ip = 0; ip < temp_planes; ip++ ) {
     // Planes ip = 0 = 1X
     // Planes ip = 2 = 2X
     if (!fPlanes[ip]) return -1;
diff --git a/src/THcHodoscope.h b/src/THcHodoscope.h
index 27a9d6b..4e3cd4e 100644
--- a/src/THcHodoscope.h
+++ b/src/THcHodoscope.h
@@ -126,6 +126,7 @@ protected:
   // Calibration
 
   // Per-event data
+  Bool_t fSHMS;
   Bool_t fGoodStartTime;
   Double_t fStartTime;
   Double_t fFPTimeAll;
@@ -212,8 +213,9 @@ protected:
   Int_t        fNHodoscopes;
 
   Int_t        fDumpTOF;
-  ofstream      fDumpOut;
+  ofstream    fDumpOut;
   string       fTOFDumpFile;
+  Bool_t      fGoodEventTOFCalib;
 
 
   Int_t fHitSweet1X;
diff --git a/src/THcScintillatorPlane.cxx b/src/THcScintillatorPlane.cxx
index ffa690b..29fec53 100644
--- a/src/THcScintillatorPlane.cxx
+++ b/src/THcScintillatorPlane.cxx
@@ -214,15 +214,19 @@ Int_t THcScintillatorPlane::ReadDatabase( const TDatime& date )
     {"hodo_adc_mode", &fADCMode, kInt, 0, 1},
     {"hodo_pedestal_scale", &fADCPedScaleFactor, kDouble, 0, 1},
     {"hodo_adc_diag_cut", &fADCDiagCut, kInt, 0, 1},
-    {0}
+   {"cosmicflag",                       &fCosmicFlag,            kInt,            0,  1},
+     {0}
   };
 
   fTofUsingInvAdc = 1;
   fADCMode = kADCStandard;
   fADCPedScaleFactor = 1.0;
   fADCDiagCut = 50.0;
-  gHcParms->LoadParmValues((DBRequest*)&list,prefix);
-  // fetch the parameter from the temporary list
+  fCosmicFlag=0;
+   gHcParms->LoadParmValues((DBRequest*)&list,prefix);
+   if (fCosmicFlag==1) cout << " setup for cosmics in scint plane"<< endl;
+  cout << " cosmic flag = " << fCosmicFlag << endl;
+ // fetch the parameter from the temporary list
 
   // Retrieve parameters we need from parent class
   // Common for all planes
@@ -675,14 +679,24 @@ Int_t THcScintillatorPlane::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
 	  timec_neg -= (hit_position-fPosRight)/
 	    fHodoNegInvAdcLinear[index];
 	  scin_corrected_time = 0.5*(timec_pos+timec_neg);
+	  if (fCosmicFlag) {
+	  postime = timec_pos + (fZpos+(index%2)*fDzpos)/(29.979*fBetaNominal);
+	  negtime = timec_neg + (fZpos+(index%2)*fDzpos)/(29.979*fBetaNominal);
+	  } else {
 	  postime = timec_pos - (fZpos+(index%2)*fDzpos)/(29.979*fBetaNominal);
 	  negtime = timec_neg - (fZpos+(index%2)*fDzpos)/(29.979*fBetaNominal);
+	  }
 	} else {
 	  postime=timec_pos-(fPosLeft-hit_position)/fHodoVelLight[index];
 	  negtime=timec_neg-(hit_position-fPosRight)/fHodoVelLight[index];
 	  scin_corrected_time = 0.5*(postime+negtime);
-	  postime = postime-(fZpos+(index%2)*fDzpos)/(29.979*fBetaNominal);
-	  negtime = negtime-(fZpos+(index%2)*fDzpos)/(29.979*fBetaNominal);
+	  if (fCosmicFlag) {
+	  postime = timec_pos + (fZpos+(index%2)*fDzpos)/(29.979*fBetaNominal);
+	  negtime = timec_neg + (fZpos+(index%2)*fDzpos)/(29.979*fBetaNominal);
+	  } else {
+	  postime = timec_pos - (fZpos+(index%2)*fDzpos)/(29.979*fBetaNominal);
+	  negtime = timec_neg - (fZpos+(index%2)*fDzpos)/(29.979*fBetaNominal);
+	  }
 	}
 	//        cout << fNScinHits<< " " << timec_pos << " " << timec_neg << endl;
         ((THcHodoHit*) fHodoHits->At(fNScinHits))->SetPaddleCenter(fPosCenter[index]);
diff --git a/src/THcScintillatorPlane.h b/src/THcScintillatorPlane.h
index 123eaf5..e14fa39 100644
--- a/src/THcScintillatorPlane.h
+++ b/src/THcScintillatorPlane.h
@@ -94,6 +94,7 @@ class THcScintillatorPlane : public THaSubDetector {
   TClonesArray* frNegAdcPulseInt;
   TClonesArray* frNegAdcPulseAmp;
 
+  Int_t fCosmicFlag; //
   Int_t fPlaneNum;		/* Which plane am I 1-4 */
   UInt_t fTotPlanes;            /* so we can read variables that are not indexed by plane id */
   UInt_t fNelem;		/* Need since we don't inherit from
-- 
GitLab