From 42bc7a0d3487bd6b7c83dcb02f662eee72ab22cf Mon Sep 17 00:00:00 2001
From: Mark Jones <jones@jlab.org>
Date: Mon, 15 Apr 2019 16:28:05 -0400
Subject: [PATCH] Add THcScintPlaneClust, Modify THcHodoscope and
 THcScintillatorPlane

THcScintPlaneCluster is a new class to hold information
about scintillator plane clusters of paddle hits.
The clusters of paddles for each plane are formed in
THcHodoscope::TrackEffTest which is used in the tracking
efficiency.

Modify THcScintillatorPlane so that the cluster information
(number of clusters, size,pos, flag if cluster used)
is in the tree.

Modify THcHodoscope::TrackEffTest so position
difference  are taken for all clusters in X1 (or Y1)
with all clusters in X2 (Y2).
---
 src/THcHodoscope.cxx         | 136 +++++++++++++++++++++++++----------
 src/THcHodoscope.h           |   4 +-
 src/THcScintPlaneCluster.cxx |   8 +++
 src/THcScintPlaneCluster.h   |  43 +++++++++++
 src/THcScintillatorPlane.cxx |  23 ++++--
 src/THcScintillatorPlane.h   |  11 +++
 6 files changed, 182 insertions(+), 43 deletions(-)
 create mode 100644 src/THcScintPlaneCluster.cxx
 create mode 100644 src/THcScintPlaneCluster.h

diff --git a/src/THcHodoscope.cxx b/src/THcHodoscope.cxx
index e28702d..3603d8e 100644
--- a/src/THcHodoscope.cxx
+++ b/src/THcHodoscope.cxx
@@ -10,6 +10,7 @@ hodoscope array, not just one plane.
 */
 
 #include "THcSignalHit.h"
+#include "THcScintPlaneCluster.h"
 #include "THcShower.h"
 #include "THcCherenkov.h"
 #include "THcHallCSpectrometer.h"
@@ -229,6 +230,8 @@ Int_t THcHodoscope::ReadDatabase( const TDatime& date )
   prefix[1]='\0';
   strcpy(parname,prefix);
   strcat(parname,"scin_");
+  //
+  //
   //  Int_t plen=strlen(parname);
   // cout << " readdatabse hodo fnplanes = " << fNPlanes << endl;
 
@@ -335,6 +338,8 @@ Int_t THcHodoscope::ReadDatabase( const TDatime& date )
     {"yloscin",                          &fyLoScin[0],            kInt,     (UInt_t) fNHodoscopes},
     {"yhiscin",                          &fyHiScin[0],            kInt,     (UInt_t) fNHodoscopes},
     {"track_eff_test_num_scin_planes",   &fTrackEffTestNScinPlanes,                 kInt},
+    {"trackeff_scint_ydiff_max",     &trackeff_scint_ydiff_max,        kDouble,         0,  1},
+    {"trackeff_scint_xdiff_max",     &trackeff_scint_xdiff_max,        kDouble,         0,  1},
     {"cer_npe",                          &fNCerNPE,               kDouble,         0,  1},
     {"normalized_energy_tot",            &fNormETot,              kDouble,         0,  1},
     {"hodo_slop",                        fHodoSlop,               kDouble,  (UInt_t) fNPlanes},
@@ -356,7 +361,8 @@ Int_t THcHodoscope::ReadDatabase( const TDatime& date )
   };
 
   // Defaults if not defined in parameter file
-
+  trackeff_scint_ydiff_max=20.;
+  trackeff_scint_xdiff_max=20.;
   for(UInt_t ip=0;ip<fMaxHodoScin;ip++) {
     fHodoPosAdcTimeWindowMin[ip] = -1000.;
     fHodoPosAdcTimeWindowMax[ip] = 1000.;
@@ -520,6 +526,15 @@ Int_t THcHodoscope::ReadDatabase( const TDatime& date )
     fTofTolerance= 3.0;
     cout << "*** USING DEFAULT 3 NSEC WINDOW FOR FP NO_TRACK CALCULATIONS!! ***\n";
   }
+  //
+  fRatio_xpfp_to_xfp=0.00;
+  TString SHMS="p";
+  TString HMS="h";
+  TString test=prefix[0];
+  if (test==SHMS ) fRatio_xpfp_to_xfp=0.0018; // SHMS 
+  if (test == HMS ) fRatio_xpfp_to_xfp=0.0011; // HMS 
+  cout << " fRatio_xpfp_to_xfp= " << fRatio_xpfp_to_xfp << endl;
+  //
   fIsInit = true;
   return kOK;
 }
@@ -1093,12 +1108,13 @@ Int_t THcHodoscope::CoarseProcess( TClonesArray& tracks )
 
 	  // Index to access the 2d arrays of paddle/scintillator properties
 	  Int_t fPIndex = GetScinIndex(ip,paddle);
+          Double_t betatrack = theTrack->GetP()/TMath::Sqrt(theTrack->GetP()*theTrack->GetP()+fPartMass*fPartMass);
           
 	  if ( TMath::Abs( scinCenter - scinTrnsCoord ) <
 	       ( fPlanes[ip]->GetSize() * 0.5 + fPlanes[ip]->GetHodoSlop() ) ){ // Line 293
 
 	    fTOFPInfo[ihhit].onTrack = kTRUE;
-	    Double_t zcor = zposition/(29.979*fBetaNominal)*
+	    Double_t zcor = zposition/(29.979*betatrack)*
 	      TMath::Sqrt(1. + theTrack->GetTheta()*theTrack->GetTheta()
 			  + theTrack->GetPhi()*theTrack->GetPhi());
 	    fTOFPInfo[ihhit].zcor = zcor;
@@ -1480,6 +1496,8 @@ void THcHodoscope::TrackEffTest(void)
   PadHigh[1]=fyHiScin[0];
   PadHigh[3]=fyHiScin[1];
   //
+  Bool_t efftest_debug = kFALSE;
+  if (efftest_debug) cout << " spec = " << GetApparatus()->GetName()[0] << endl;
   Double_t PadPosLo[4];
   Double_t PadPosHi[4];
   for (Int_t ip = 0; ip < fNumPlanesBetaCalc; ip++ ){
@@ -1511,70 +1529,112 @@ void THcHodoscope::TrackEffTest(void)
       if ( hit->GetTwoGoodTimes() ) {
 	if ( padnum==prev_padnum+1 ) {
 	  fClustSize[ip][fNClust[ip]-1]=fClustSize[ip][fNClust[ip]-1]+1;
-	  fClustPos[ip][fNClust[ip]-1]=fClustPos[ip][fNClust[ip]-1]+fPlanes[ip]->GetPosCenter(padnum-1)+ fPlanes[ip]->GetPosOffset();
-	  //  cout << "Add to cluster  pl = " << ip+1 << " hit = " << iphit << " pad = " << padnum << " clus =  " << fNClust[ip] << " cl size = " << fClustSize[ip][fNClust[ip]-1] << " pos " << fPlanes[ip]->GetPosCenter(padnum-1)+ fPlanes[ip]->GetPosOffset() << endl;
+	  fClustPos[ip][fNClust[ip]-1]+=fPlanes[ip]->GetPosCenter(padnum-1)+ fPlanes[ip]->GetPosOffset();
+	   if (efftest_debug) cout << "Add to cluster  pl = " << ip+1 << " hit = " << iphit << " pad = " << padnum << " clus =  " << fNClust[ip] << " cl size = " << fClustSize[ip][fNClust[ip]-1] << " pos " << fPlanes[ip]->GetPosCenter(padnum-1)+ fPlanes[ip]->GetPosOffset() << endl;
 	} else {
 	  if (fNClust[ip]<MaxNClus) fNClust[ip]++;
 	  fClustSize[ip][fNClust[ip]-1]=1;
 	  fClustPos[ip][fNClust[ip]-1]=fPlanes[ip]->GetPosCenter(padnum-1)+ fPlanes[ip]->GetPosOffset();
-	  //  cout << " New clus pl = " << ip+1 << " hit = " << iphit << " pad = " << padnum << " clus = " << fNClust[ip] << " cl size = " << fClustSize[ip][fNClust[ip]] << " pos " << fPlanes[ip]->GetPosCenter(padnum-1)+ fPlanes[ip]->GetPosOffset() << endl;
+	   if (efftest_debug) cout << " New clus pl = " << ip+1 << " hit = " << iphit << " pad = " << padnum << " clus = " << fNClust[ip] << " cl size = " << fClustSize[ip][fNClust[ip]-1] << " pos " << fPlanes[ip]->GetPosCenter(padnum-1)+ fPlanes[ip]->GetPosOffset() << endl;
 	}
 	prev_padnum=padnum;
       }
+      if (!(hit->GetTwoGoodTimes()) && efftest_debug)  cout << "no two good times  plane = " << ip+1 << " hit = " << iphit << endl;
     }
   }
   //
-  Bool_t inside_bound[4]={kFALSE,kFALSE,kFALSE,kFALSE};
+  Bool_t inside_bound[4][MaxNClus];
   for(Int_t ip = 0; ip < fNumPlanesBetaCalc; ip++ ) {	 
+    fPlanes[ip]->SetNumberClusters(fNClust[ip]);
     for(Int_t ic = 0; ic <fNClust[ip] ; ic++ ) {
       fClustPos[ip][ic]=fClustPos[ip][ic]/fClustSize[ip][ic];
-      inside_bound[ip] = fClustPos[ip][ic]>=PadPosLo[ip] &&  fClustPos[ip][ic]<=PadPosHi[ip];
-      //cout << "plane = " << ip+1 << " Cluster = " << ic+1 << " size = " << fClustSize[ip][ic]<< " pos = " << fClustPos[ip][ic] << " inside = " << inside_bound[ip] << " lo = " << PadPosLo[ip]<< " hi = " << PadPosHi[ip]<< endl;
+      fPlanes[ip]->SetCluster(ic,fClustPos[ip][ic]);
+      fPlanes[ip]->SetClusterSize(ic,fClustSize[ip][ic]);
+     inside_bound[ip][ic] = fClustPos[ip][ic]>=PadPosLo[ip] &&  fClustPos[ip][ic]<=PadPosHi[ip];
+      if (efftest_debug) cout << "plane = " << ip+1 << " Cluster = " << ic+1 << " size = " << fClustSize[ip][ic]<< " pos = " << fClustPos[ip][ic] << " inside = " << inside_bound[ip][ic] << " lo = " << PadPosLo[ip]<< " hi = " << PadPosHi[ip]<< endl;
     }
   }
   //
-  Int_t good_for_track_test[4]={0,0,0,0};
-  Int_t sum_good_track_test=0;
+  Int_t MaxClusterSize=3;
+  Int_t good_for_track_test[4][MaxNClus];
+  Int_t sum_good_track_test[4]={0,0,0,0};
+  Int_t num_good_plane_hit=0;
   for(Int_t ip = 0; ip < fNumPlanesBetaCalc; ip++ ) {
-    if (fNClust[ip]==1 && inside_bound[ip] && fClustSize[ip][0]<=2) good_for_track_test[ip]=1;
-    //cout << " good for track = " << good_for_track_test[ip] << endl;
-    sum_good_track_test+=good_for_track_test[ip];
+    for(Int_t ic = 0; ic <fNClust[ip] ; ic++ ) {
+      if (inside_bound[ip][ic] && fClustSize[ip][ic]<=MaxClusterSize) {
+         fPlanes[ip]->SetClusterFlag(ic,1.);
+         good_for_track_test[ip][ic]=1;
+	  sum_good_track_test[ip]++;
+	  if (sum_good_track_test[ip]==1) num_good_plane_hit++;
+      } else {
+           good_for_track_test[ip][ic]=0;
+     }
+      if (efftest_debug) cout << " ip " << ip+1 << " clus = " << ic << " good for track = " << good_for_track_test[ip][ic] << endl;
+    }
+    if (efftest_debug) cout << " ip = " << ip+1 << "  sum_good_track_test = " << sum_good_track_test[ip] << endl;
   }	 
+  if (efftest_debug) cout << " number of planes hits = " << num_good_plane_hit << endl;
   //
-  Double_t trackeff_scint_ydiff_max= 10. ;
-  Double_t trackeff_scint_xdiff_max= 10. ;
   Bool_t xdiffTest=kFALSE;
   Bool_t ydiffTest=kFALSE;
   fGoodScinHits = 0;
-  if (fTrackEffTestNScinPlanes == 4) {
-    if (fTrackEffTestNScinPlanes==sum_good_track_test) {
-      xdiffTest= TMath::Abs(fClustPos[0][0]-fClustPos[2][0])<trackeff_scint_xdiff_max;
-      ydiffTest= TMath::Abs(fClustPos[1][0]-fClustPos[3][0])<trackeff_scint_ydiff_max;
-      if (xdiffTest && ydiffTest) fGoodScinHits = 1;
+  if (efftest_debug) cout << " fTrackEffTestNScinPlanes = " << fTrackEffTestNScinPlanes << endl;
+  if ( (fTrackEffTestNScinPlanes == 4 || fTrackEffTestNScinPlanes == 3) && num_good_plane_hit==4) {
+    
+    // check for matching clusters in the X planes assumed to be planes 0 and 2
+    for(Int_t ic0 = 0; ic0 <fNClust[0] ; ic0++ ) {
+    for(Int_t ic2 = 0; ic2 <fNClust[2] ; ic2++ ) {
+      if (good_for_track_test[0][ic0] && good_for_track_test[2][ic2]) {
+           Double_t x1_proj = fClustPos[0][ic0]*(1+fRatio_xpfp_to_xfp*(fPlanes[2]->GetZpos()-fPlanes[0]->GetZpos())); // project X1 to X2 Z position
+           xdiffTest= TMath::Abs(x1_proj-fClustPos[2][ic2])<trackeff_scint_xdiff_max;
+          if (xdiffTest) fPlanes[0]->SetClusterUsedFlag(ic0,1.);
+          if (xdiffTest) fPlanes[2]->SetClusterUsedFlag(ic2,1.);
+      }
     }
+    }
+    // check for matching clusters in the Y planes assumed to be planes 1 and 3
+    for(Int_t ic1 = 0; ic1 <fNClust[1] ; ic1++ ) {
+    for(Int_t ic3 = 0; ic3 <fNClust[3] ; ic3++ ) {
+       if (good_for_track_test[1][ic1] && good_for_track_test[3][ic3]) {
+           ydiffTest= TMath::Abs(fClustPos[1][ic1]-fClustPos[3][ic3])<trackeff_scint_ydiff_max;
+          if (ydiffTest) fPlanes[1]->SetClusterUsedFlag(ic1,1.);
+          if (ydiffTest) fPlanes[3]->SetClusterUsedFlag(ic3,1.);
+       }
+    }
+    }
+    if (xdiffTest && ydiffTest) fGoodScinHits = 1;
+    if (efftest_debug) cout << " 4 good planes  xdiff = " << xdiffTest << " ydiff = " <<  ydiffTest << endl;
   }
   //
-  if (fTrackEffTestNScinPlanes == 3) {
-    if (fTrackEffTestNScinPlanes==sum_good_track_test) {
-      if(good_for_track_test[0]==1&&good_for_track_test[2]==1) {
-	xdiffTest= TMath::Abs(fClustPos[0][0]-fClustPos[2][0])<trackeff_scint_xdiff_max;
-	ydiffTest=kTRUE;
-      }
-      if (good_for_track_test[1]==1&&good_for_track_test[3]==1) {
-	xdiffTest=kTRUE;
-	ydiffTest= TMath::Abs(fClustPos[1][0]-fClustPos[3][0])<trackeff_scint_ydiff_max;
-      }
-      if (xdiffTest && ydiffTest) fGoodScinHits = 1;
+  if (fTrackEffTestNScinPlanes == 3 && num_good_plane_hit==3) {
+    xdiffTest=kFALSE;
+    ydiffTest=kFALSE;
+    // Check if two X planes hit
+    if (sum_good_track_test[0]>0&&sum_good_track_test[2]>0) {
+       for(Int_t ic0 = 0; ic0 <fNClust[0] ; ic0++ ) {
+       for(Int_t ic2 = 0; ic2 <fNClust[2] ; ic2++ ) {
+         if (good_for_track_test[0][ic0] && good_for_track_test[2][ic2]) {
+          xdiffTest= TMath::Abs(fClustPos[0][ic0]-fClustPos[2][ic2])<trackeff_scint_xdiff_max;
+         }
+       }
+       }
+      ydiffTest = kTRUE;
+    }
+    // Check if two Y planes hit
+    if ((sum_good_track_test[1]>0||sum_good_track_test[3]>0)) {
+    for(Int_t ic1 = 0; ic1 <fNClust[1] ; ic1++ ) {
+    for(Int_t ic3 = 0; ic3 <fNClust[3] ; ic3++ ) {
+       if (good_for_track_test[1][ic1] && good_for_track_test[3][ic3]) {
+           ydiffTest= TMath::Abs(fClustPos[1][ic1]-fClustPos[3][ic3])<trackeff_scint_ydiff_max;
+       }
     }
-    if (sum_good_track_test==4) {
-      xdiffTest= TMath::Abs(fClustPos[0][0]-fClustPos[2][0])<trackeff_scint_xdiff_max;
-      ydiffTest= TMath::Abs(fClustPos[1][0]-fClustPos[3][0])<trackeff_scint_ydiff_max;
-      if (xdiffTest && ydiffTest) fGoodScinHits = 1;
+      xdiffTest = kTRUE;
     }
+    }  
+     if (xdiffTest && ydiffTest) fGoodScinHits = 1;
+    if (efftest_debug) cout << " 3 good planes  xdiff = " << xdiffTest << " ydiff = " <<  ydiffTest << endl;
   }
-  //       
-  //	cout << " good scin = " << fGoodScinHits << " " << sum_good_track_test << " " << xdiffTest  << " " << ydiffTest<< endl;
-  //cout << " ************" << endl;
+  if (efftest_debug) cout << " ************" << endl;
   //
 }
 //
diff --git a/src/THcHodoscope.h b/src/THcHodoscope.h
index 25635c6..3712855 100644
--- a/src/THcHodoscope.h
+++ b/src/THcHodoscope.h
@@ -148,7 +148,9 @@ protected:
 
   TH1F *hTime;
   // Calibration
-
+  Double_t fRatio_xpfp_to_xfp;
+  Double_t trackeff_scint_ydiff_max ;
+  Double_t trackeff_scint_xdiff_max ;
   // Per-event data
   Bool_t fSHMS;
   Bool_t fGoodStartTime;
diff --git a/src/THcScintPlaneCluster.cxx b/src/THcScintPlaneCluster.cxx
new file mode 100644
index 0000000..f5d5681
--- /dev/null
+++ b/src/THcScintPlaneCluster.cxx
@@ -0,0 +1,8 @@
+/** \class THcScintPlaneCluster
+    \ingroup DetSupport
+
+*/
+
+#include "THcScintPlaneCluster.h"
+
+ClassImp(THcScintPlaneCluster)
diff --git a/src/THcScintPlaneCluster.h b/src/THcScintPlaneCluster.h
new file mode 100644
index 0000000..d6c31aa
--- /dev/null
+++ b/src/THcScintPlaneCluster.h
@@ -0,0 +1,43 @@
+#ifndef ROOT_THcScintPlaneCluster
+#define ROOT_THcScintPlaneCluster
+
+/////////////////////////////////////////////////////////////////////////////
+//                                                                         //
+// THcScintPlaneCluster                                                             //
+//                                                                         //
+/////////////////////////////////////////////////////////////////////////////
+
+#include "TObject.h"
+#include <cstdio>
+
+class THcScintPlaneCluster : public TObject {
+
+ public:
+ THcScintPlaneCluster(Int_t cluster=0, Double_t pos=0.0) :
+  fClusterNumber(cluster), fClusterPosition(pos) {}
+  virtual ~THcScintPlaneCluster() {}
+
+  Int_t GetClusterNumber() {return fClusterNumber;}
+  Double_t GetClusterPosition() {return fClusterPosition;}
+  Double_t GetClusterSize() {return fClusterSize;}
+  Double_t GetClusterFlag() {return fClusterFlag;}
+  Double_t GetClusterUsedFlag() {return fClusterUsedFlag;}
+
+  virtual void Set(Int_t cluster, Double_t pos)
+  { fClusterNumber=cluster; fClusterPosition=pos; fClusterFlag=0.; fClusterUsedFlag=0.;}
+
+  void SetSize(Double_t size) {fClusterSize=size;}
+  void SetFlag(Double_t flag) {fClusterFlag=flag;}
+  void SetUsedFlag(Double_t flag) {fClusterUsedFlag=flag;}
+
+ private:
+  Int_t fClusterNumber;
+  Double_t fClusterPosition;
+  Double_t fClusterSize;
+  Double_t fClusterFlag;
+  Double_t fClusterUsedFlag;
+
+  ClassDef(THcScintPlaneCluster,0); // Single signal value and wire/counter number
+};
+/////////////////////////////////////////////////////////////////
+#endif
diff --git a/src/THcScintillatorPlane.cxx b/src/THcScintillatorPlane.cxx
index 175ab79..6de7f7a 100644
--- a/src/THcScintillatorPlane.cxx
+++ b/src/THcScintillatorPlane.cxx
@@ -8,6 +8,7 @@ The THcHodoscope class instatiates one of these objects per plane.
 */
 #include "TMath.h"
 #include "THcScintillatorPlane.h"
+#include "THcScintPlaneCluster.h"
 #include "TClonesArray.h"
 #include "THcSignalHit.h"
 #include "THcHodoHit.h"
@@ -34,7 +35,7 @@ THcScintillatorPlane::THcScintillatorPlane( const char* name,
 					    const Int_t planenum,
 					    THaDetectorBase* parent )
 : THaSubDetector(name,description,parent),
-  fParentHitList(0), frPosAdcErrorFlag(0), frNegAdcErrorFlag(0),
+  fParentHitList(0), fCluster(0),frPosAdcErrorFlag(0), frNegAdcErrorFlag(0),
   frPosTDCHits(0), frNegTDCHits(0), frPosADCHits(0), frNegADCHits(0),
   frPosADCSums(0), frNegADCSums(0), frPosADCPeds(0), frNegADCPeds(0),
   fHodoHits(0), frPosTdcTimeRaw(0), frPosAdcPedRaw(0),
@@ -61,6 +62,8 @@ THcScintillatorPlane::THcScintillatorPlane( const char* name,
   // Normal constructor with name and description
   fHodoHits = new TClonesArray("THcHodoHit",16);
 
+  fCluster = new TClonesArray("THcScintPlaneCluster", 10);
+
   frPosAdcErrorFlag = new TClonesArray("THcSignalHit", 16);
   frNegAdcErrorFlag = new TClonesArray("THcSignalHit", 16);
 
@@ -97,6 +100,7 @@ THcScintillatorPlane::THcScintillatorPlane( const char* name,
   frNegAdcPulseAmp  = new TClonesArray("THcSignalHit", 16);
   frNegAdcPulseTime = new TClonesArray("THcSignalHit", 16);
 
+
   fPlaneNum = planenum;
   fTotPlanes = planenum;
   fNScinHits = 0;
@@ -111,6 +115,8 @@ THcScintillatorPlane::~THcScintillatorPlane()
   delete  frPosAdcErrorFlag; frPosAdcErrorFlag = NULL;
   delete  frNegAdcErrorFlag; frNegAdcErrorFlag = NULL;
 
+  delete  fCluster; fCluster = NULL;
+
   delete fHodoHits;
   delete frPosTDCHits;
   delete frNegTDCHits;
@@ -384,6 +390,7 @@ Int_t THcScintillatorPlane::ReadDatabase( const TDatime& date )
   // Create arrays to hold results here
   InitializePedestals();
 
+
   fNumGoodPosAdcHits     = vector<Int_t> (fNelem, 0.0);
   fNumGoodNegAdcHits     = vector<Int_t> (fNelem, 0.0);
   fNumGoodPosTdcHits     = vector<Int_t> (fNelem, 0.0);
@@ -526,7 +533,12 @@ Int_t THcScintillatorPlane::DefineVariables( EMode mode )
     {"DiffDisTrackCorr",   "TW Corrected Dist Difference between track and scintillator position (cm)", "fGoodDiffDistTrack"},
     {"TrackXPos",   "Track X position at plane (cm)", "fTrackXPosition"},
     {"TrackYPos",   "Track Y position at plane (cm)", "fTrackYPosition"},
-    //{"ngoodhits", "Number of paddle hits (passed tof tolerance and used to determine the focal plane time )",           "GetNGoodHits() "},
+    {"NumClus",   "Number of clusters", "fNumberClusters"},
+    {"Clus.Pos",   "Position of each paddle clusters", "fCluster.THcScintPlaneCluster.GetClusterPosition()"},
+    {"Clus.Size",   "Size of each paddle clusters", "fCluster.THcScintPlaneCluster.GetClusterSize()"},
+    {"Clus.Flag",   "Flag of each paddle clusters", "fCluster.THcScintPlaneCluster.GetClusterFlag()"},
+    {"Clus.UsedFlag",   "USed Flag of each paddle clusters", "fCluster.THcScintPlaneCluster.GetClusterUsedFlag()"},
+   //{"ngoodhits", "Number of paddle hits (passed tof tolerance and used to determine the focal plane time )",           "GetNGoodHits() "},
     { 0 }
   };
 
@@ -542,6 +554,8 @@ void THcScintillatorPlane::Clear( Option_t* )
    */
   //cout << " Calling THcScintillatorPlane::Clear " << GetName() << endl;
   // Clears the hit lists
+  fCluster->Clear();
+
   frPosAdcErrorFlag->Clear();
   frNegAdcErrorFlag->Clear();
 
@@ -575,6 +589,7 @@ void THcScintillatorPlane::Clear( Option_t* )
   frNegAdcPulseAmp->Clear();
   frNegAdcPulseTime->Clear();
 
+
   //Clear occupancies
   for (UInt_t ielem = 0; ielem < fNumGoodPosAdcHits.size(); ielem++)
     fNumGoodPosAdcHits.at(ielem) = 0;
@@ -629,6 +644,7 @@ void THcScintillatorPlane::Clear( Option_t* )
   fHitDistance = kBig;
   fTrackXPosition = kBig;
   fTrackYPosition = kBig;
+  fNumberClusters=0;
 }
 
 //_____________________________________________________________________________
@@ -724,7 +740,6 @@ Int_t THcScintillatorPlane::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
   frNegAdcPulseInt->Clear();
   frNegAdcPulseAmp->Clear();
   frNegAdcPulseTime->Clear();
-
   //stripped
   fNScinHits=0;
 
@@ -1030,7 +1045,7 @@ Int_t THcScintillatorPlane::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
       }
 
       // Do corrections if valid TDC on both ends of bar
-      if(btdcraw_pos && btdcraw_neg && badcraw_pos && badcraw_neg) {
+      if( (btdcraw_pos && btdcraw_neg) && (badcraw_pos && badcraw_neg) ) {
 	// Do the pulse height correction to the time.  (Position dependent corrections later)
 	Double_t timec_pos, timec_neg;
 	if(fTofUsingInvAdc) {
diff --git a/src/THcScintillatorPlane.h b/src/THcScintillatorPlane.h
index e39d514..7dc7e30 100644
--- a/src/THcScintillatorPlane.h
+++ b/src/THcScintillatorPlane.h
@@ -14,6 +14,7 @@
 
 #include "THaSubDetector.h"
 #include "TClonesArray.h"
+#include "THcScintPlaneCluster.h"
 
 using namespace std;
 
@@ -57,22 +58,31 @@ class THcScintillatorPlane : public THaSubDetector {
   Double_t GetPosOffset() {return fPosOffset;};
   Double_t GetPosCenter(Int_t PaddleNo) {return fPosCenter[PaddleNo];}; // counting from zero!
   Double_t GetFpTime() {return fFptime;};
+  Double_t GetNumberClusters() {return fNumberClusters;}; // number of paddle clusters
 
   void SetFpTime(Double_t f) {fFptime=f;};
   void SetNGoodHits(Int_t ng) {fNGoodHits=ng;};
   void SetHitDistance(Double_t f) {fHitDistance=f;}; // Distance between track and hit paddle
   void SetTrackXPosition(Double_t f) {fTrackXPosition=f;}; // Distance track X position at plane
   void SetTrackYPosition(Double_t f) {fTrackYPosition=f;}; // Distance track Y position at plane
+  void SetNumberClusters(Int_t nclus) {fNumberClusters=nclus;}; // number of paddle 
+  void SetCluster(Int_t ic,Double_t pos) {((THcScintPlaneCluster*) fCluster->ConstructedAt(ic))->Set(ic,pos);}
+  void SetClusterSize(Int_t ic,Double_t size) {((THcScintPlaneCluster*) fCluster->ConstructedAt(ic))->SetSize(size);}
+  void SetClusterFlag(Int_t ic,Double_t flag) {((THcScintPlaneCluster*) fCluster->ConstructedAt(ic))->SetFlag(flag);}
+  void SetClusterUsedFlag(Int_t ic,Double_t flag) {((THcScintPlaneCluster*) fCluster->ConstructedAt(ic))->SetUsedFlag(flag);}
 
   TClonesArray* fParentHitList;
 
   TClonesArray* GetHits() { return fHodoHits;};
 
+  TClonesArray* fCluster;
+
  protected:
 
   TClonesArray* frPosAdcErrorFlag;
   TClonesArray* frNegAdcErrorFlag;
 
+
   TClonesArray* frPosTDCHits;
   TClonesArray* frNegTDCHits;
   TClonesArray* frPosADCHits;
@@ -111,6 +121,7 @@ class THcScintillatorPlane : public THaSubDetector {
   Int_t fTotNumPosAdcHits;
   Int_t fTotNumNegAdcHits;
   Int_t fTotNumAdcHits;
+  Int_t fNumberClusters;
 
   Int_t fTotNumPosTdcHits;
   Int_t fTotNumNegTdcHits;
-- 
GitLab