diff --git a/examples/hodtest.C b/examples/hodtest.C
index 4ae408b4713b8626ab5a9d5a2fd9d7f0ca60bea8..977b1e4b13880a83ded7f7946b4c546288b7c661 100644
--- a/examples/hodtest.C
+++ b/examples/hodtest.C
@@ -64,7 +64,7 @@
analyzer->SetEvent( event );
analyzer->SetOutFile( "hodtest.root" );
analyzer->SetOdefFile("output.def");
- // analyzer->SetCutFile("cuts_example.def"); // optional
+ analyzer->SetCutFile("hodtest_cuts.def"); // optional
// File to record cuts accounting information
// analyzer->SetSummaryFile("summary_example.log"); // optional
diff --git a/examples/hodtest_cuts.def b/examples/hodtest_cuts.def
new file mode 100644
index 0000000000000000000000000000000000000000..f67e2b34b7ba962a9018c3feb4f68919e385d2a0
--- /dev/null
+++ b/examples/hodtest_cuts.def
@@ -0,0 +1,11 @@
+# Demo cuts for hodtest
+#
+
+Block: RawDecode
+
+Pedestal_event g.evtyp==4
+RawDecode_master 1
+
+Block: Decode
+Decode_master !xPedestal_event
+
diff --git a/src/THcAnalyzer.h b/src/THcAnalyzer.h
index 172d12560cea9e381cae9465e9209313b38c8b39..3506ac50ef33db55e40ad69cb2a264b7c6c43b0b 100644
--- a/src/THcAnalyzer.h
+++ b/src/THcAnalyzer.h
@@ -16,7 +16,11 @@ public:
THcAnalyzer();
virtual ~THcAnalyzer();
+ void SetPedestalEvtype( Int_t evtype ) { fPedestalEvtype = evtype; }
+
protected:
+
+ Int_t fPedestalEvtype;
private:
// THcAnalyzer( const THcAnalyzer& );
diff --git a/src/THcHodoscope.cxx b/src/THcHodoscope.cxx
index 769d92d0ac346a888dfa00ffbd52840024cb0195..d1931a9b847efd2eb0bc7348a9995a580d0305fb 100644
--- a/src/THcHodoscope.cxx
+++ b/src/THcHodoscope.cxx
@@ -13,6 +13,8 @@
#include "THaEvData.h"
#include "THaDetMap.h"
#include "THcDetectorMap.h"
+#include "THaGlobals.h"
+#include "THaCutList.h"
#include "THcGlobals.h"
#include "THcParmList.h"
#include "VarDef.h"
@@ -330,6 +332,21 @@ Int_t THcHodoscope::Decode( const THaEvData& evdata )
// Get the Hall C style hitlist (fRawHitList) for this event
Int_t nhits = THcHitList::DecodeToHitList(evdata);
+ if(gHaCuts->Result("Pedestal_event")) {
+ Int_t nexthit = 0;
+ for(Int_t ip=0;ip<fNPlanes;ip++) {
+ nexthit = fPlanes[ip]->AccumulatePedestals(fRawHitList, nexthit);
+ }
+ fAnalyzePedestals = 1; // Analyze pedestals first normal events
+ return(0);
+ }
+ if(fAnalyzePedestals) {
+ for(Int_t ip=0;ip<fNPlanes;ip++) {
+ fPlanes[ip]->CalculatePedestals();
+ }
+ fAnalyzePedestals = 0; // Don't analyze pedestals next event
+ }
+
// Let each plane get its hits
Int_t nexthit = 0;
for(Int_t ip=0;ip<fNPlanes;ip++) {
diff --git a/src/THcHodoscope.h b/src/THcHodoscope.h
index e0937b89af863f4abdb727306d8315574f980726..23540d1b573b9960527a960f6c098eeb8eb4f879 100644
--- a/src/THcHodoscope.h
+++ b/src/THcHodoscope.h
@@ -41,6 +41,8 @@ public:
THcHodoscope(); // for ROOT I/O
protected:
+ Int_t fAnalyzePedestals;
+
// Calibration
// Per-event data
diff --git a/src/THcScintillatorPlane.cxx b/src/THcScintillatorPlane.cxx
index 7f966c20a54b61f1168a702bdf432ce5e18fb6c5..350c3de69bc038de337c9de0ca40e142f84935f7 100644
--- a/src/THcScintillatorPlane.cxx
+++ b/src/THcScintillatorPlane.cxx
@@ -87,7 +87,7 @@ Int_t THcScintillatorPlane::ReadDatabase( const TDatime& date )
// See what file it looks for
- static const char* const here = "ReadDatabase()";
+ // static const char* const here = "ReadDatabase()";
char prefix[2];
char parname[100];
@@ -124,7 +124,7 @@ Int_t THcScintillatorPlane::ReadDatabase( const TDatime& date )
// Create arrays to hold results here
-
+ InitializePedestals();
return kOK;
}
@@ -187,6 +187,8 @@ Int_t THcScintillatorPlane::FineProcess( TClonesArray& tracks )
{
return 0;
}
+
+//_____________________________________________________________________________
Int_t THcScintillatorPlane::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
{
// Extract the data for this plane from hit list
@@ -271,6 +273,98 @@ Int_t THcScintillatorPlane::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
return(ihit);
}
+//_____________________________________________________________________________
+Int_t THcScintillatorPlane::AccumulatePedestals(TClonesArray* rawhits, Int_t nexthit)
+{
+ // Extract the data for this plane from hit list, accumulating into
+ // arrays for calculating pedestals.
+
+ Int_t nrawhits = rawhits->GetLast()+1;
+ // cout << "THcScintillatorPlane::AcculatePedestals " << fPlaneNum << " " << nexthit << "/" << nrawhits << endl;
+
+ Int_t ihit = nexthit;
+ while(ihit < nrawhits) {
+ THcHodoscopeHit* hit = (THcHodoscopeHit *) rawhits->At(ihit);
+ if(hit->fPlane > fPlaneNum) {
+ break;
+ }
+ Int_t element = hit->fCounter - 1; // Should check if in range
+ Int_t adcpos = hit->fADC_pos;
+ Int_t adcneg = hit->fADC_neg;
+
+ if(adcpos <= fPosPedLimit[element]) {
+ fPosPedSum[element] += adcpos;
+ fPosPedSum2[element] += adcpos*adcpos;
+ fPosPedCount[element]++;
+ if(fPosPedCount[element] == fMinPeds/5) {
+ fPosPedLimit[element] = 100 + fPosPedSum[element]/fPosPedCount[element];
+ }
+ }
+ if(adcneg <= fNegPedLimit[element]) {
+ fNegPedSum[element] += adcneg;
+ fNegPedSum2[element] += adcneg*adcneg;
+ fNegPedCount[element]++;
+ if(fNegPedCount[element] == fMinPeds/5) {
+ fNegPedLimit[element] = 100 + fNegPedSum[element]/fNegPedCount[element];
+ }
+ }
+ ihit++;
+ }
+
+ fNPedestalEvents++;
+
+ return(ihit);
+}
+
+//_____________________________________________________________________________
+void THcScintillatorPlane::CalculatePedestals( )
+{
+ // Use the accumulated pedestal data to calculate pedestals
+ // Later add check to see if pedestals have drifted ("Danger Will Robinson!")
+ // cout << "Plane: " << fPlaneNum << endl;
+ for(Int_t i=0; i<fNelem;i++) {
+ // Positive tubes
+ fPosPed[i] = ((Double_t) fPosPedSum[i]) / TMath::Max(1, fPosPedCount[i]);
+ fPosThresh[i] = fPosPed[i] + 15;
+
+ // Negative tubes
+ fNegPed[i] = ((Double_t) fNegPedSum[i]) / TMath::Max(1, fNegPedCount[i]);
+ fNegThresh[i] = fNegPed[i] + 15;
+
+ // cout << i+1 << " " << fPosPed[i] << " " << fNegPed[i] << endl;
+ }
+ // cout << " " << endl;
-
+}
+
+//_____________________________________________________________________________
+void THcScintillatorPlane::InitializePedestals( )
+{
+ fNPedestalEvents = 0;
+ fMinPeds = 500; // In engine, this is set in parameter file
+ fPosPedSum = new Int_t [fNelem];
+ fPosPedSum2 = new Int_t [fNelem];
+ fPosPedLimit = new Int_t [fNelem];
+ fPosPedCount = new Int_t [fNelem];
+ fNegPedSum = new Int_t [fNelem];
+ fNegPedSum2 = new Int_t [fNelem];
+ fNegPedLimit = new Int_t [fNelem];
+ fNegPedCount = new Int_t [fNelem];
+
+ fPosPed = new Double_t [fNelem];
+ fNegPed = new Double_t [fNelem];
+ fPosThresh = new Double_t [fNelem];
+ fNegThresh = new Double_t [fNelem];
+ for(Int_t i=0;i<fNelem;i++) {
+ fPosPedSum[i] = 0;
+ fPosPedSum2[i] = 0;
+ fPosPedLimit[i] = 1000; // In engine, this are set in parameter file
+ fPosPedCount[i] = 0;
+ fNegPedSum[i] = 0;
+ fNegPedSum2[i] = 0;
+ fNegPedLimit[i] = 1000; // In engine, this are set in parameter file
+ fNegPedCount[i] = 0;
+ }
+}
+
diff --git a/src/THcScintillatorPlane.h b/src/THcScintillatorPlane.h
index dda5931fed5480320f08a38a0ce9eb620da187d8..431fbd85a21f4638e6e760aa18d0670ddda02d0c 100644
--- a/src/THcScintillatorPlane.h
+++ b/src/THcScintillatorPlane.h
@@ -35,6 +35,8 @@ class THcScintillatorPlane : public THaSubDetector {
virtual Bool_t IsPid() { return kFALSE; }
virtual Int_t ProcessHits(TClonesArray* rawhits, Int_t nexthit);
+ virtual Int_t AccumulatePedestals(TClonesArray* rawhits, Int_t nexthit);
+ virtual void CalculatePedestals( );
Double_t fSpacing;
@@ -51,11 +53,31 @@ class THcScintillatorPlane : public THaSubDetector {
TClass* fNegTDCHitsClass;
TClass* fPosADCHitsClass;
TClass* fNegADCHitsClass;
-
- Int_t fPlaneNum;
-
+
+ Int_t fPlaneNum; /* Which plane am I 1-4 */
+ Int_t fNelem; /* Need since we don't inherit from
+ detector base class */
+ Int_t fNPedestalEvents; /* Number of pedestal events */
+ Int_t fMinPeds; /* Only analyze/update if num events > */
+ Int_t *fPosPedSum; /* Accumulators for pedestals */
+ Int_t *fPosPedSum2;
+ Int_t *fPosPedLimit;
+ Int_t *fPosPedCount;
+ Int_t *fNegPedSum;
+ Int_t *fNegPedSum2;
+ Int_t *fNegPedLimit;
+ Int_t *fNegPedCount;
+
+ Double_t *fPosPed;
+ Double_t *fPosSig;
+ Double_t *fPosThresh;
+ Double_t *fNegPed;
+ Double_t *fNegSig;
+ Double_t *fNegThresh;
+
virtual Int_t ReadDatabase( const TDatime& date );
virtual Int_t DefineVariables( EMode mode = kDefine );
+ virtual void InitializePedestals( );
ClassDef(THcScintillatorPlane,0)
};