From 4fe33ab73a55f56f77d0db876e21af780ff20e69 Mon Sep 17 00:00:00 2001
From: Carlos Yero <cyero002@fiu.edu>
Date: Mon, 8 May 2017 08:15:52 -0400
Subject: [PATCH] Updates to THcShowerPlane   Add "Good" variables to
 calorimeter class      occupancies and multiplicities      required to pass
 time and threshold cut   Change some variables from arrays to vectors  
 Modify how thresholds are calculated   Remove obsolete variables: posadchits,
 negadchits   Add flash adc debug flag   Correctly delete pointers in
 destructor

---
 src/THcShowerArray.cxx | 202 ++++++++++++++-------
 src/THcShowerArray.h   |  38 ++--
 src/THcShowerPlane.cxx | 399 +++++++++++++++++++++++++++--------------
 src/THcShowerPlane.h   |  68 +++++--
 4 files changed, 470 insertions(+), 237 deletions(-)

diff --git a/src/THcShowerArray.cxx b/src/THcShowerArray.cxx
index 57a231f..39944a4 100644
--- a/src/THcShowerArray.cxx
+++ b/src/THcShowerArray.cxx
@@ -56,11 +56,23 @@ THcShowerArray::THcShowerArray( const char* name,
 THcShowerArray::~THcShowerArray()
 {
   // Destructor
-  delete fXPos;
-  delete fYPos;
-  delete fZPos;
+ 
+  for (UInt_t i=0; i<fNRows; i++) {
+    delete [] fXPos[i];
+    delete [] fYPos[i];
+    delete [] fZPos[i];
+  }
+
+  delete [] fPedLimit;
+  delete [] fGain;
+  delete [] fPedSum;
+  delete [] fPedSum2;
+  delete [] fPedCount;
+  delete [] fSig;
+  delete [] fPed;
+  delete [] fThresh;
 
-  delete fADCHits;
+  delete fADCHits; fADCHits = NULL;
 
   delete frAdcPedRaw; frAdcPedRaw = NULL;
   delete frAdcErrorFlag; frAdcErrorFlag = NULL;
@@ -72,11 +84,11 @@ THcShowerArray::~THcShowerArray()
   delete frAdcPulseInt; frAdcPulseInt = NULL;
   delete frAdcPulseAmp; frAdcPulseAmp = NULL;
 
-  delete [] fA;
-  delete [] fP;
-  delete [] fA_p;
+  //  delete [] fA;
+  //delete [] fP;
+  // delete [] fA_p;
 
-  delete [] fE;
+  //delete [] fE;
   delete [] fBlock_ClusterID;
 }
 
@@ -135,9 +147,13 @@ Int_t THcShowerArray::ReadDatabase( const TDatime& date )
     {"cal_ped_sample_high", &fPedSampHigh, kInt, 0, 1},
     {"cal_data_sample_low", &fDataSampLow, kInt, 0, 1},
     {"cal_data_sample_high", &fDataSampHigh, kInt, 0, 1},
+    {"cal_debug_adc", &fDebugAdc, kInt, 0, 1},
     {0}
   };
-  gHcParms->LoadParmValues((DBRequest*)&list, prefix);
+  
+  fDebugAdc = 0;  // Set ADC debug parameter to false unless set in parameter file
+
+  gHcParms->LoadParmValues((DBRequest*)&list, prefix);  
   fADCMode=kADCDynamicPedestal;
   fAdcTimeWindowMin=0;
   fAdcTimeWindowMax=10000;
@@ -307,14 +323,24 @@ Int_t THcShowerArray::ReadDatabase( const TDatime& date )
   InitializePedestals();
 
   // Event by event amplitude and pedestal
-  fA = new Double_t[fNelem];
-  fP = new Double_t[fNelem];
-  fA_p = new Double_t[fNelem];
+  //fA = new Double_t[fNelem];
+  //fP = new Double_t[fNelem];
+  //fA_p = new Double_t[fNelem];
+  
+  fE                       = vector<Double_t> (fNelem, 0.0);
+  fNumGoodAdcHits          = vector<Int_t>    (fNelem, 0.0);
+  fGoodAdcPulseIntRaw      = vector<Double_t> (fNelem, 0.0);
+  fGoodAdcPed              = vector<Double_t> (fNelem, 0.0);
+  fGoodAdcPulseInt         = vector<Double_t> (fNelem, 0.0);
+  fGoodAdcPulseAmp         = vector<Double_t> (fNelem, 0.0);
+  fGoodAdcPulseTime        = vector<Double_t> (fNelem, 0.0); 
+
+
   fBlock_ClusterID = new Int_t[fNelem];
 
   // Energy depositions per block.
 
-  fE = new Double_t[fNelem];
+  //fE = new Double_t[fNelem];
 
 #ifdef HITPIC
   hitpic = new char*[fNRows];
@@ -345,35 +371,49 @@ Int_t THcShowerArray::DefineVariables( EMode mode )
   fIsSetup = ( mode == kDefine );
 
   // Register variables in global list
-  RVarDef vars[] = {
-    {"adchits", "List of ADC hits", "fADCHits.THcSignalHit.GetPaddleNumber()"},
-    {"a", "Raw ADC Amplitude", "fA"},
-    {"p", "Dynamic ADC Pedestal", "fP"},
-    {"a_p", "Sparsified, ped-subtracted ADC Amplitudes", "fA_p"},
-    { "nhits", "Number of hits", "fNhits" },
-    { "nghits", "Number of good hits ( pass threshold on raw ADC)", "fNgoodhits" },
-    { "nclust", "Number of clusters", "fNclust" },
-    {"e", "Energy Depositions per block", "fE"},
-    {"block_clusterID", "Cluster ID number", "fBlock_ClusterID"},
-    {"earray", "Energy Deposition in array", "fEarray"},
-    { "ntracks", "Number of shower tracks", "fNtracks" },
-
-    {"adcCounter",      "List of ADC counter numbers.",      "frAdcPulseIntRaw.THcSignalHit.GetPaddleNumber()"},
-
-    {"adcPedRaw",       "List of raw ADC pedestals",         "frAdcPedRaw.THcSignalHit.GetData()"},
-    {"adcErrorFlag",       "List of raw ADC pedestals",      "frAdcErrorFlag.THcSignalHit.GetData()"},
-    {"adcPulseIntRaw",  "List of raw ADC pulse integrals.",  "frAdcPulseIntRaw.THcSignalHit.GetData()"},
-    {"adcPulseAmpRaw",  "List of raw ADC pulse amplitudes.", "frAdcPulseAmpRaw.THcSignalHit.GetData()"},
-    {"adcPulseTimeRaw", "List of raw ADC pulse times.",      "frAdcPulseTimeRaw.THcSignalHit.GetData()"},
-
-    {"adcPed",          "List of ADC pedestals",             "frAdcPed.THcSignalHit.GetData()"},
-    {"adcPulseInt",     "List of ADC pulse integrals.",      "frAdcPulseInt.THcSignalHit.GetData()"},
-    {"adcPulseAmp",     "List of ADC pulse amplitudes.",     "frAdcPulseAmp.THcSignalHit.GetData()"},
-
-    { 0 }
-  };
+  vector<RVarDef> vars;
+
+    //{"adchits", "List of ADC hits", "fADCHits.THcSignalHit.GetPaddleNumber()"}, // appears an empty histogram in the root file
+    
+    vars.push_back(RVarDef{"adcErrorFlag",       "Error Flag When FPGA Fails",      "frAdcErrorFlag.THcSignalHit.GetData()"});
+    
+    vars.push_back(RVarDef{"adcCounter",      "List of ADC counter numbers.",      "frAdcPulseIntRaw.THcSignalHit.GetPaddleNumber()"});  //raw occupancy
+    vars.push_back(RVarDef{"numGoodAdcHits", "Number of Good ADC Hits per PMT", "fNumGoodAdcHits" });                                   //good occupancy
+
+    vars.push_back(RVarDef{"totNumAdcHits", "Total Number of ADC Hits", "fTotNumAdcHits" });                                            // raw multiplicity
+    vars.push_back(RVarDef{"totNumGoodAdcHits", "Total Number of Good ADC Hits", "fTotNumGoodAdcHits" });                               // good multiplicity
+
+
+    vars.push_back(RVarDef{"goodAdcPulseIntRaw", "Good Raw ADC Pulse Integrals", "fGoodAdcPulseIntRaw"});    //this is defined as pulseIntRaw, NOT ADC Amplitude in FillADC_DynamicPedestal() method
+    vars.push_back(RVarDef{"goodAdcPed", "Good ADC Pedestals", "fGoodAdcPed"});   
+    vars.push_back(RVarDef{"goodAdcPulseInt", "Good ADC Pulse Integrals", "fGoodAdcPulseInt"});     //this is defined as pulseInt, which is the pedestal subtracted pulse integrals, and is defined if threshold is passed
+    vars.push_back(RVarDef{"goodAdcPulseAmp", "Good ADC Pulse Amplitudes", "fGoodAdcPulseAmp"});   
+    vars.push_back(RVarDef{"goodAdcPulseTime", "Good ADC Pulse Times", "fGoodAdcPulseTime"});     //this is defined as pulseInt, which is the pedestal subtracted pulse integrals, and is defined if threshold is passed
+
+    
+    vars.push_back(RVarDef{"e", "Energy Depositions per block", "fE"});       //defined as fE = fA_p*fGain = pulseInt * Gain
+    vars.push_back(RVarDef{"earray", "Energy Deposition in Shower Array", "fEarray"});   //defined as a Double_t and represents a sum of the total deposited energy in the shower counter
+
+    vars.push_back(RVarDef{"nclust", "Number of clusters", "fNclust" });       //what is the difference between nclust defined here and that in THcShower.cxx ?
+    vars.push_back(RVarDef{"block_clusterID", "Cluster ID number", "fBlock_ClusterID"}); // im NOT very clear about this. it is histogrammed at wither -1 or 0.
+    vars.push_back(RVarDef{"ntracks", "Number of shower tracks", "fNtracks" }); //number of cluster-to-track associations
+ 
+    if (fDebugAdc) {
+    vars.push_back(RVarDef{"adcPedRaw",       "List of raw ADC pedestals",         "frAdcPedRaw.THcSignalHit.GetData()"});  
+    vars.push_back(RVarDef{"adcPulseIntRaw",  "List of raw ADC pulse integrals.",  "frAdcPulseIntRaw.THcSignalHit.GetData()"});
+    vars.push_back(RVarDef{"adcPulseAmpRaw",  "List of raw ADC pulse amplitudes.", "frAdcPulseAmpRaw.THcSignalHit.GetData()"});
+    vars.push_back(RVarDef{"adcPulseTimeRaw", "List of raw ADC pulse times.",      "frAdcPulseTimeRaw.THcSignalHit.GetData()"});
+
+    vars.push_back(RVarDef{"adcPed",          "List of ADC pedestals",             "frAdcPed.THcSignalHit.GetData()"});
+    vars.push_back(RVarDef{"adcPulseInt",     "List of ADC pulse integrals.",      "frAdcPulseInt.THcSignalHit.GetData()"});
+    vars.push_back(RVarDef{"adcPulseAmp",     "List of ADC pulse amplitudes.",     "frAdcPulseAmp.THcSignalHit.GetData()"});
+
+    }
 
-  return DefineVarsFromList( vars, mode );
+    RVarDef end {0};
+    vars.push_back(end);
+
+    return DefineVarsFromList(vars.data(), mode );
 }
 
 //_____________________________________________________________________________
@@ -382,8 +422,8 @@ void THcShowerArray::Clear( Option_t* )
   // Clears the hit lists
   fADCHits->Clear();
 
-  fNhits = 0;
-  fNgoodhits = 0;
+  fTotNumAdcHits = 0;
+  fTotNumGoodAdcHits = 0;
   fNclust = 0;
   fClustSize = 0;
   fNtracks = 0;
@@ -408,6 +448,16 @@ void THcShowerArray::Clear( Option_t* )
   frAdcPulseInt->Clear();
   frAdcPulseAmp->Clear();
 
+  for (UInt_t ielem = 0; ielem < fGoodAdcPed.size(); ielem++) {
+    fGoodAdcPulseIntRaw.at(ielem)      = 0.0;   
+    fGoodAdcPed.at(ielem)              = 0.0;
+    fGoodAdcPulseInt.at(ielem)         = 0.0;
+    fGoodAdcPulseAmp.at(ielem)         = 0.0;   
+    fGoodAdcPulseTime.at(ielem)        = 0.0;
+    fNumGoodAdcHits.at(ielem)          = 0.0;
+    fE.at(ielem)                       = 0.0;
+  }
+
 }
 
 //_____________________________________________________________________________
@@ -424,7 +474,6 @@ Int_t THcShowerArray::CoarseProcess( TClonesArray& tracks )
 
   // Fill set of unclustered shower array hits.
   // Reuse hit class pertained to the HMS/SOS type calorimeters.
-  // Save Y coordinate of the hit in Z parameter of the class.
   // Save energy deposition in the module as hit mean energy, do not use
   // positive and negative side energies.
 
@@ -434,7 +483,7 @@ Int_t THcShowerArray::CoarseProcess( TClonesArray& tracks )
   for(UInt_t j=0; j < fNColumns; j++) {
     for (UInt_t i=0; i<fNRows; i++) {
 
-      if (fA_p[k] > 0) {    //hit
+      if (fGoodAdcPulseInt.at(k) > 0) {    //hit
 
 	THcShowerHit* hit =
 	  new THcShowerHit(i, j, fXPos[i][j], fYPos[i][j], fZPos[i][j], fE[k], 0., 0.);
@@ -454,13 +503,22 @@ Int_t THcShowerArray::CoarseProcess( TClonesArray& tracks )
     cout << "Debug output from THcShowerArray::CoarseProcess for " << GetName()
 	 << endl;
 
-    cout << "  List of unclustered hits. Total hits:     " << fNhits << endl;
+    cout << "  List of unclustered hits. Total hits:     " << fTotNumAdcHits << endl;
     THcShowerHitIt it = HitSet.begin();    //<set> version
-    for (Int_t i=0; i!=fNgoodhits; i++) {
+    for (Int_t i=0; i!=fTotNumGoodAdcHits; i++) {
       cout << "  hit " << i << ": ";
       (*(it++))->show();
     }
   }
+    
+  ////Sanity check. (Vardan)
+  
+  // if ((int)HitSet.size() != fTotNumGoodAdcHits) {
+  //	cout << "***" << endl;
+  //	cout << "*** THcShowerArray::CoarseProcess: HitSet.size = " << HitSet.size()
+  //	     << " != fTotNumGoodAdcHits = " << fTotNumGoodAdcHits << endl;
+  //	cout << "***" << endl;
+  //    }
 
   // Cluster hits and fill list of clusters.
 
@@ -718,7 +776,7 @@ Int_t THcShowerArray::CoarseProcessHits()
     Int_t k=0;
     for(UInt_t j=0; j < fNColumns; j++) {
     for (UInt_t i=0; i<fNRows; i++) {
-     if(fA[k] > fThresh[k]) {
+      if(fGoodAdcPulseIntRaw.at(k) > fThresh[k]) {
 	cout << "  counter =  " << k
 	     << "  E = " << fE[k]
 	     << endl;
@@ -763,19 +821,30 @@ void THcShowerArray::FillADC_DynamicPedestal()
 {
   for (Int_t ielem=0;ielem<frAdcPulseInt->GetEntries();ielem++) {
     Int_t npad = ((THcSignalHit*) frAdcPulseInt->ConstructedAt(ielem))->GetPaddleNumber() - 1;
-    Double_t pulseInt = ((THcSignalHit*) frAdcPulseInt->ConstructedAt(ielem))->GetData();
     Double_t pulseIntRaw = ((THcSignalHit*) frAdcPulseIntRaw->ConstructedAt(ielem))->GetData();
+    Double_t pulsePed     = ((THcSignalHit*) frAdcPed->ConstructedAt(ielem))->GetData(); 
+    Double_t pulseInt = ((THcSignalHit*) frAdcPulseInt->ConstructedAt(ielem))->GetData();
+    Double_t pulseAmp     = ((THcSignalHit*) frAdcPulseAmp->ConstructedAt(ielem))->GetData();
     Double_t pulseTime = ((THcSignalHit*) frAdcPulseTimeRaw->ConstructedAt(ielem))->GetData();
-    Double_t errorflag = ((THcSignalHit*) frAdcErrorFlag->ConstructedAt(ielem))->GetData();
+    Bool_t errorflag = ((THcSignalHit*) frAdcErrorFlag->ConstructedAt(ielem))->GetData();
     Bool_t pulseTimeCut = (pulseTime > fAdcTimeWindowMin) &&  (pulseTime < fAdcTimeWindowMax);
-    if (errorflag==0 && pulseTimeCut) {
-      fNhits++;
-      fA[npad] =pulseIntRaw;
-      if(fA[npad] >  fThresh[npad]) {
-       fNgoodhits++;
-       fA_p[npad] =pulseInt ;
-       fE[npad] = fA_p[npad]*fGain[npad];
-       fEarray += fE[npad];
+    if (!errorflag && pulseTimeCut) {
+      fTotNumAdcHits++;
+      fGoodAdcPulseIntRaw.at(npad) = pulseIntRaw;
+      
+      if(fGoodAdcPulseIntRaw.at(npad) >  fThresh[npad]) {
+       fTotNumGoodAdcHits++;
+       fGoodAdcPulseInt.at(npad) = pulseInt;
+       fE.at(npad) = fGoodAdcPulseInt.at(npad)*fGain[npad];
+       fEarray += fE.at(npad);
+     
+       fGoodAdcPed.at(npad) = pulsePed; 
+       fGoodAdcPulseAmp.at(npad) = pulseAmp; 
+       fGoodAdcPulseTime.at(npad) = pulseTime;  
+
+       fNumGoodAdcHits.at(npad) = npad + 1;
+
+
       }
      }        
    }
@@ -804,9 +873,10 @@ Int_t THcShowerArray::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
   frAdcPulseAmp->Clear();
 
   for(Int_t i=0;i<fNelem;i++) {
-    fA[i] = 0;
-    fA_p[i] = 0;
-    fE[i] = 0;
+    //fA[i] = 0;
+    //fA_p[i] = 0;
+    //fE[i] = 0;
+
     fBlock_ClusterID[i] = -1;
   }
 
@@ -854,7 +924,7 @@ Int_t THcShowerArray::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
   }
 
 #if 0
-  if(fNgoodhits > 0) {
+  if(fTotNumGoodAdcHits > 0) {
     cout << "+";
     for(Int_t column=0;column<fNColumns;column++) {
       cout << "-";
@@ -864,7 +934,7 @@ Int_t THcShowerArray::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
       cout << "|";
       for(Int_t column=0;column<fNColumns;column++) {
 	Int_t counter = column*fNRows + row;
-	if(fA[counter]>threshold) {
+	if(fGoodAdcPulseIntRaw.at(counter) > threshold) {
 	  cout << "X";
 	} else {
 	  cout << " ";
@@ -875,14 +945,14 @@ Int_t THcShowerArray::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
   }
 #endif
 #ifdef HITPIC
-  if(fNgoodhits > 0) {
+  if(fTotNumGoodAdcHits > 0) {
     for(Int_t row=0;row<fNRows;row++) {
       if(piccolumn==0) {
 	hitpic[row][0] = '|';
       }
       for(Int_t column=0;column<fNColumns;column++) {
 	Int_t counter = column*fNRows+row;
-	if(fA[counter] > threshold) {
+	if(fGoodAdcPulseIntRaw.at(counter) > threshold) {
 	  hitpic[row][piccolumn*(fNColumns+1)+column+1] = 'X';
 	} else {
 	  hitpic[row][piccolumn*(fNColumns+1)+column+1] = ' ';
@@ -898,7 +968,7 @@ Int_t THcShowerArray::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
 	  cout << "-";
 	}
 	cout << "+";
-      }
+0      }
       cout << endl;
       for(Int_t row=0;row<fNRows;row++) {
 	hitpic[row][(piccolumn+1)*(fNColumns+1)+1] = '\0';
diff --git a/src/THcShowerArray.h b/src/THcShowerArray.h
index a7e2a0e..8e45937 100644
--- a/src/THcShowerArray.h
+++ b/src/THcShowerArray.h
@@ -95,11 +95,27 @@ protected:
   char **hitpic;
   Int_t piccolumn;
 #endif
-  Double_t* fA;               // [fNelem] ADC amplitudes of blocks
-  Double_t* fP;               // [fNelem] Event by event (FADC) pedestals
-  Double_t* fA_p;	      // [fNelem] sparsified, pedestal subtracted
-                              // (FASTBUS) ADC amplitudes
 
+
+ 
+
+
+  //counting variables
+  Int_t fTotNumAdcHits;              // Total number of ADC hits
+  Int_t fTotNumGoodAdcHits;          // Total number of good ADC hits (pass threshold)
+
+  vector<Int_t>         fNumGoodAdcHits;          // shower good occupancy
+  vector<Double_t>      fGoodAdcPulseIntRaw;      // [fNelem] Good Raw ADC pulse Integrals of blocks
+  
+  vector<Double_t>      fGoodAdcPed;             // [fNelem] Event by event (FADC) good pulse pedestals
+  vector<Double_t>      fGoodAdcPulseInt;       // [fNelem] good pedestal subtracted pulse integrals
+  vector<Double_t>      fGoodAdcPulseAmp;
+  vector<Double_t>      fGoodAdcPulseTime;
+
+  vector<Double_t>      fE;                    //[fNelem] energy deposition in shower blocks
+
+  Int_t* fBlock_ClusterID;              // [fNelem] Cluster ID of the block -1 then not in a cluster
+  Double_t  fEarray;                          // Total Energy deposition in the array.
   TClonesArray* fADCHits;	// List of ADC hits
 
   // Parameters
@@ -128,7 +144,9 @@ protected:
   Double_t fAdcTimeWindowMin ;
   Double_t fAdcTimeWindowMax ;
   Double_t fAdcThreshold ;
-Int_t fPedSampLow;		// Sample range for
+
+  Int_t fDebugAdc;
+  Int_t fPedSampLow;		// Sample range for
   Int_t fPedSampHigh;		// dynamic pedestal
   Int_t fDataSampLow;		// Sample range for
   Int_t fDataSampHigh;		// sample integration
@@ -152,15 +170,7 @@ Int_t fPedSampLow;		// Sample range for
   Float_t *fThresh;          // [fNelem] ADC thresholds
 
   Double_t* fGain;           // [fNelem] Gain constants from calibration
-
-  //Energy depositions.
-
-  Double_t* fE;              // [fNelem] energy depositions in the blocks.
-  Int_t* fBlock_ClusterID;              // [fNelem] Cluster ID of the block -1 then not in a cluster
-  Double_t  fEarray;         // Total Energy deposition in the array.
-
-  Int_t fNhits;              // Total number of hits
-  Int_t fNgoodhits;              // Total number of good hits (pass threshold)
+  
   Int_t fNclust;             // Number of hit clusters
   Int_t fNtracks;            // Number of shower tracks, i.e. number of
                              // cluster-to-track associations
diff --git a/src/THcShowerPlane.cxx b/src/THcShowerPlane.cxx
index 53b1c82..866b981 100644
--- a/src/THcShowerPlane.cxx
+++ b/src/THcShowerPlane.cxx
@@ -24,6 +24,7 @@ One plane of shower blocks with side readout
 #include <iostream>
 
 #include <fstream>
+
 using namespace std;
 
 ClassImp(THcShowerPlane)
@@ -41,6 +42,7 @@ THcShowerPlane::THcShowerPlane( const char* name,
 
   frPosAdcErrorFlag = new TClonesArray("THcSignalHit", 16);
   frPosAdcPedRaw = new TClonesArray("THcSignalHit", 16);
+  frPosAdcThreshold = new TClonesArray("THcSignalHit", 16);
   frPosAdcPulseIntRaw = new TClonesArray("THcSignalHit", 16);
   frPosAdcPulseAmpRaw = new TClonesArray("THcSignalHit", 16);
   frPosAdcPulseTimeRaw = new TClonesArray("THcSignalHit", 16);
@@ -51,6 +53,7 @@ THcShowerPlane::THcShowerPlane( const char* name,
 
   frNegAdcErrorFlag = new TClonesArray("THcSignalHit", 16);
   frNegAdcPedRaw = new TClonesArray("THcSignalHit", 16);
+  frNegAdcThreshold = new TClonesArray("THcSignalHit", 16);
   frNegAdcPulseIntRaw = new TClonesArray("THcSignalHit", 16);
   frNegAdcPulseAmpRaw = new TClonesArray("THcSignalHit", 16);
   frNegAdcPulseTimeRaw = new TClonesArray("THcSignalHit", 16);
@@ -71,37 +74,52 @@ THcShowerPlane::THcShowerPlane( const char* name,
 THcShowerPlane::~THcShowerPlane()
 {
   // Destructor
-  delete fPosADCHits;
-  delete fNegADCHits;
-
-  frPosAdcErrorFlag = NULL;
-  frPosAdcPedRaw = NULL;
-  frPosAdcPulseIntRaw = NULL;
-  frPosAdcPulseAmpRaw = NULL;
-  frPosAdcPulseTimeRaw = NULL;
-
-  frPosAdcPed = NULL;
-  frPosAdcPulseInt = NULL;
-  frPosAdcPulseAmp = NULL;
-
-  frNegAdcErrorFlag = NULL;
-  frNegAdcPedRaw = NULL;
-  frNegAdcPulseIntRaw = NULL;
-  frNegAdcPulseAmpRaw = NULL;
-  frNegAdcPulseTimeRaw = NULL;
-
-  frNegAdcPed = NULL;
-  frNegAdcPulseInt = NULL;
-  frNegAdcPulseAmp = NULL;
-
-  delete [] fA_Pos;
-  delete [] fA_Neg;
-  delete [] fA_Pos_p;
-  delete [] fA_Neg_p;
-
-  delete [] fEpos;
-  delete [] fEneg;
-  delete [] fEmean;
+  delete fPosADCHits; fPosADCHits = NULL;
+  delete fNegADCHits; fNegADCHits = NULL;
+
+  delete  frPosAdcErrorFlag; frPosAdcErrorFlag = NULL;
+  delete  frPosAdcPedRaw; frPosAdcPedRaw = NULL;
+  delete  frPosAdcThreshold; frPosAdcThreshold = NULL;
+  delete  frPosAdcPulseIntRaw; frPosAdcPulseIntRaw = NULL;
+  delete  frPosAdcPulseAmpRaw; frPosAdcPulseAmpRaw = NULL;
+  delete  frPosAdcPulseTimeRaw; frPosAdcPulseTimeRaw = NULL;
+
+  delete  frPosAdcPed; frPosAdcPed = NULL;
+  delete  frPosAdcPulseInt; frPosAdcPulseInt = NULL;
+  delete  frPosAdcPulseAmp; frPosAdcPulseAmp = NULL;
+
+  delete  frNegAdcErrorFlag; frNegAdcErrorFlag = NULL;
+  delete  frNegAdcPedRaw; frNegAdcPedRaw = NULL;
+  delete  frNegAdcThreshold; frNegAdcThreshold = NULL;
+  delete  frNegAdcPulseIntRaw; frNegAdcPulseIntRaw = NULL;
+  delete  frNegAdcPulseAmpRaw; frNegAdcPulseAmpRaw = NULL;
+  delete  frNegAdcPulseTimeRaw; frNegAdcPulseTimeRaw = NULL;
+
+  delete  frNegAdcPed; frNegAdcPed = NULL;
+  delete  frNegAdcPulseInt; frNegAdcPulseInt = NULL;
+  delete  frNegAdcPulseAmp; frNegAdcPulseAmp = NULL;
+
+  delete [] fPosPedSum;
+  delete [] fPosPedSum2;
+  delete [] fPosPedLimit;
+  delete [] fPosPedCount;
+
+  delete [] fNegPedSum;
+  delete [] fNegPedSum2;
+  delete [] fNegPedLimit;
+  delete [] fNegPedCount;
+
+  delete [] fPosPed;
+  delete [] fPosSig;
+  delete [] fPosThresh;
+
+  delete [] fNegPed;
+  delete [] fNegSig;
+  delete [] fNegThresh;
+
+//  delete [] fEpos;
+  // delete [] fEneg;
+  // delete [] fEmean;
 }
 
 //_____________________________________________________________________________
@@ -138,8 +156,8 @@ Int_t THcShowerPlane::ReadDatabase( const TDatime& date )
   fPedSampHigh=9;
   fDataSampLow=23;
   fDataSampHigh=49;
-    fAdcNegThreshold=0.;
-    fAdcPosThreshold=0.;
+  fAdcNegThreshold=0.;
+  fAdcPosThreshold=0.;
   DBRequest list[]={
     {"cal_AdcNegThreshold", &fAdcNegThreshold, kDouble, 0, 1},
     {"cal_AdcPosThreshold", &fAdcPosThreshold, kDouble, 0, 1},
@@ -147,8 +165,12 @@ Int_t THcShowerPlane::ReadDatabase( const TDatime& date )
     {"cal_ped_sample_high", &fPedSampHigh, kInt, 0, 1},
     {"cal_data_sample_low", &fDataSampLow, kInt, 0, 1},
     {"cal_data_sample_high", &fDataSampHigh, kInt, 0, 1},
+    {"cal_debug_adc", &fDebugAdc, kInt, 0, 1},
     {0}
   };
+
+  fDebugAdc   = 0; // Set ADC debug parameter to false unless set in parameter file
+
   gHcParms->LoadParmValues((DBRequest*)&list, prefix);
 
   // Retrieve more parameters we need from parent class
@@ -196,18 +218,30 @@ Int_t THcShowerPlane::ReadDatabase( const TDatime& date )
 
   InitializePedestals();
 
-  // ADC amplitudes per channel.
+  fNumGoodPosAdcHits          = vector<Int_t>    (fNelem, 0.0);
+  fNumGoodNegAdcHits          = vector<Int_t>    (fNelem, 0.0);
 
-  fA_Pos = new Double_t[fNelem];
-  fA_Neg = new Double_t[fNelem];
-  fA_Pos_p = new Double_t[fNelem];
-  fA_Neg_p = new Double_t[fNelem];
+  fGoodPosAdcPulseIntRaw      = vector<Double_t> (fNelem, 0.0);
+  fGoodNegAdcPulseIntRaw      = vector<Double_t> (fNelem, 0.0);
+
+  fGoodPosAdcPed              = vector<Double_t> (fNelem, 0.0);
+  fGoodPosAdcPulseInt         = vector<Double_t> (fNelem, 0.0);
+  fGoodPosAdcPulseAmp         = vector<Double_t> (fNelem, 0.0);
+  fGoodPosAdcPulseTime        = vector<Double_t> (fNelem, 0.0);
+
+  fGoodNegAdcPed              = vector<Double_t> (fNelem, 0.0);
+  fGoodNegAdcPulseInt         = vector<Double_t> (fNelem, 0.0);
+  fGoodNegAdcPulseAmp         = vector<Double_t> (fNelem, 0.0);
+  fGoodNegAdcPulseTime        = vector<Double_t> (fNelem, 0.0);
 
   // Energy depositions per block (not corrected for track coordinate)
 
-  fEpos = new Double_t[fNelem];
-  fEneg = new Double_t[fNelem];
-  fEmean= new Double_t[fNelem];
+  fEpos                       = vector<Double_t> (fNelem, 0.0);
+  fEneg                       = vector<Double_t> (fNelem, 0.0);
+  fEmean                      = vector<Double_t> (fNelem, 0.0);
+  //  fEpos = new Double_t[fNelem];
+  // fEneg = new Double_t[fNelem];
+  // fEmean= new Double_t[fNelem];
 
   // Debug output.
 
@@ -245,58 +279,83 @@ Int_t THcShowerPlane::DefineVariables( EMode mode )
   fIsSetup = ( mode == kDefine );
 
   // Register variables in global list
-  RVarDef vars[] = {
-    {"posadchits", "List of Positive ADC hits","fPosADCHits.THcSignalHit.GetPaddleNumber()"},
-    {"negadchits", "List of Negative ADC hits","fNegADCHits.THcSignalHit.GetPaddleNumber()"},
-    {"apos",       "Raw Positive ADC Amplitudes",                   "fA_Pos"},
-    {"aneg",       "Raw Negative ADC Amplitudes",                   "fA_Neg"},
-    {"apos_p",     "Ped-subtracted Positive ADC Amplitudes",        "fA_Pos_p"},
-    {"aneg_p",     "Ped-subtracted Negative ADC Amplitudes",        "fA_Neg_p"},
-    {"epos",       "Energy Depositions from Positive Side PMTs",    "fEpos"},
-    {"eneg",       "Energy Depositions from Negative Side PMTs",    "fEneg"},
-    {"emean",      "Mean Energy Depositions",                       "fEmean"},
-    {"eplane",     "Energy Deposition per plane",                   "fEplane"},
-    {"eplane_pos", "Energy Deposition per plane from pos. PMTs","fEplane_pos"},
-    {"eplane_neg", "Energy Deposition per plane from neg. PMTs","fEplane_neg"},
-
-    {"posAdcCounter",      "List of positive ADC counter numbers.",      "frPosAdcPulseIntRaw.THcSignalHit.GetPaddleNumber()"},
-    {"negAdcCounter",      "List of negative ADC counter numbers.",      "frNegAdcPulseIntRaw.THcSignalHit.GetPaddleNumber()"},
-
-    {"posAdcErrorFlag",       "List of positive raw ADC Error Flags",         "frPosAdcErrorFlag.THcSignalHit.GetData()"},
-    {"posAdcPedRaw",       "List of positive raw ADC pedestals",         "frPosAdcPedRaw.THcSignalHit.GetData()"},
-    {"posAdcPulseIntRaw",  "List of positive raw ADC pulse integrals.",  "frPosAdcPulseIntRaw.THcSignalHit.GetData()"},
-    {"posAdcPulseAmpRaw",  "List of positive raw ADC pulse amplitudes.", "frPosAdcPulseAmpRaw.THcSignalHit.GetData()"},
-    {"posAdcPulseTimeRaw", "List of positive raw ADC pulse times.",      "frPosAdcPulseTimeRaw.THcSignalHit.GetData()"},
-
-    {"posAdcPed",          "List of positive ADC pedestals",             "frPosAdcPed.THcSignalHit.GetData()"},
-    {"posAdcPulseInt",     "List of positive ADC pulse integrals.",      "frPosAdcPulseInt.THcSignalHit.GetData()"},
-    {"posAdcPulseAmp",     "List of positive ADC pulse amplitudes.",     "frPosAdcPulseAmp.THcSignalHit.GetData()"},
-
-    {"negAdcErrorFlag",       "List of negative raw ADC Error Flag ",         "frNegAdcErrorFlag.THcSignalHit.GetData()"},
-    {"negAdcPedRaw",       "List of negative raw ADC pedestals",         "frNegAdcPedRaw.THcSignalHit.GetData()"},
-    {"negAdcPulseIntRaw",  "List of negative raw ADC pulse integrals.",  "frNegAdcPulseIntRaw.THcSignalHit.GetData()"},
-    {"negAdcPulseAmpRaw",  "List of negative raw ADC pulse amplitudes.", "frNegAdcPulseAmpRaw.THcSignalHit.GetData()"},
-    {"negAdcPulseTimeRaw", "List of negative raw ADC pulse times.",      "frNegAdcPulseTimeRaw.THcSignalHit.GetData()"},
-
-    {"negAdcPed",          "List of negative ADC pedestals",             "frNegAdcPed.THcSignalHit.GetData()"},
-    {"negAdcPulseInt",     "List of negative ADC pulse integrals.",      "frNegAdcPulseInt.THcSignalHit.GetData()"},
-    {"negAdcPulseAmp",     "List of negative ADC pulse amplitudes.",     "frNegAdcPulseAmp.THcSignalHit.GetData()"},
-
-    { 0 }
-  };
+  vector<RVarDef> vars;
+
+  vars.push_back(RVarDef{"posAdcErrorFlag",       "List of positive raw ADC Error Flags",         "frPosAdcErrorFlag.THcSignalHit.GetData()"});
+  vars.push_back(RVarDef{"negAdcErrorFlag",       "List of negative raw ADC Error Flags ",        "frNegAdcErrorFlag.THcSignalHit.GetData()"});
+
+  vars.push_back(RVarDef{"posAdcCounter",      "List of positive ADC counter numbers.",      "frPosAdcPulseIntRaw.THcSignalHit.GetPaddleNumber()"}); //PreSh+ raw occupancy
+  vars.push_back(RVarDef{"negAdcCounter",      "List of negative ADC counter numbers.",      "frNegAdcPulseIntRaw.THcSignalHit.GetPaddleNumber()"}); //PreSh- raw occupancy
+
+  vars.push_back(RVarDef{"totNumPosAdcHits", "Total Number of Positive ADC Hits",   "fTotNumPosAdcHits"}); // PreSh+ raw multiplicity
+  vars.push_back(RVarDef{"totNumNegAdcHits", "Total Number of Negative ADC Hits",   "fTotNumNegAdcHits"}); // PreSh+ raw multiplicity
+  vars.push_back(RVarDef{"totnumAdcHits",   "Total Number of ADC Hits Per PMT",      "fTotNumAdcHits"});    // PreSh raw multiplicity
+
+  vars.push_back(RVarDef{"numGoodPosAdcHits",    "Number of Good Positive ADC Hits Per PMT", "fNumGoodPosAdcHits"});    // PreSh occupancy
+  vars.push_back(RVarDef{"numGoodNegAdcHits",    "Number of Good Negative ADC Hits Per PMT", "fNumGoodNegAdcHits"});    // PreSh occupancy
+  vars.push_back(RVarDef{"totNumGoodPosAdcHits", "Total Number of Good Positive ADC Hits",   "fTotNumGoodPosAdcHits"}); // PreSh multiplicity
+  vars.push_back(RVarDef{"totNumGoodNegAdcHits", "Total Number of Good Negative ADC Hits",   "fTotNumGoodNegAdcHits"}); // PreSh multiplicity
+  vars.push_back(RVarDef{"totnumGoodAdcHits",   "TotalNumber of Good ADC Hits Per PMT",      "fTotNumGoodAdcHits"});    // PreSh multiplicity
+
+  vars.push_back(RVarDef{"goodPosAdcPulseIntRaw",       "Good Positive Raw ADC Pulse Integrals",                   "fGoodPosAdcPulseIntRaw"});
+  vars.push_back(RVarDef{"goodNegAdcPulseIntRaw",       "Good Negative Raw ADC Pulse Integrals",                   "fGoodNegAdcPulseIntRaw"});
+
+  vars.push_back(RVarDef{"goodPosAdcPed",         "Good Positive ADC pedestals",           "fGoodPosAdcPed"});
+  vars.push_back(RVarDef{"goodPosAdcPulseInt",         "Good Positive ADC integrals",           "fGoodPosAdcPulseInt"});
+  vars.push_back(RVarDef{"goodPosAdcPulseAmp",         "Good Positive ADC amplitudes",           "fGoodPosAdcPulseAmp"});
+  vars.push_back(RVarDef{"goodPosAdcPulseTime",         "Good Positive ADC times",           "fGoodPosAdcPulseTime"});
 
-  return DefineVarsFromList( vars, mode );
+  vars.push_back(RVarDef{"goodNegAdcPed",         "Good Negative ADC pedestals",           "fGoodNegAdcPed"});
+  vars.push_back(RVarDef{"goodNegAdcPulseInt",         "Good Negative ADC integrals",           "fGoodNegAdcPulseInt"});
+  vars.push_back(RVarDef{"goodNegAdcPulseAmp",         "Good Negative ADC amplitudes",           "fGoodNegAdcPulseAmp"});
+  vars.push_back(RVarDef{"goodNegAdcPulseTime",         "Good Negative ADC times",           "fGoodNegAdcPulseTime"});
+
+  vars.push_back(RVarDef{"epos",       "Energy Depositions from Positive Side PMTs",    "fEpos"});
+  vars.push_back(RVarDef{"eneg",       "Energy Depositions from Negative Side PMTs",    "fEneg"});
+  vars.push_back(RVarDef{"emean",      "Mean Energy Depositions",                       "fEmean"});
+  vars.push_back(RVarDef{"eplane",     "Energy Deposition per plane",                   "fEplane"});
+  vars.push_back(RVarDef{"eplane_pos", "Energy Deposition per plane from pos. PMTs","fEplane_pos"});
+  vars.push_back(RVarDef{"eplane_neg", "Energy Deposition per plane from neg. PMTs","fEplane_neg"});
+
+
+  if (fDebugAdc) {
+  vars.push_back(RVarDef{"posAdcPedRaw",       "List of positive raw ADC pedestals",         "frPosAdcPedRaw.THcSignalHit.GetData()"});
+  vars.push_back(RVarDef{"posAdcPulseIntRaw",  "List of positive raw ADC pulse integrals.",  "frPosAdcPulseIntRaw.THcSignalHit.GetData()"});
+  vars.push_back(RVarDef{"posAdcPulseAmpRaw",  "List of positive raw ADC pulse amplitudes.", "frPosAdcPulseAmpRaw.THcSignalHit.GetData()"});
+  vars.push_back(RVarDef{"posAdcPulseTimeRaw", "List of positive raw ADC pulse times.",      "frPosAdcPulseTimeRaw.THcSignalHit.GetData()"});
+
+  vars.push_back(RVarDef{"posAdcPed",          "List of positive ADC pedestals",             "frPosAdcPed.THcSignalHit.GetData()"});
+  vars.push_back(RVarDef{"posAdcPulseInt",     "List of positive ADC pulse integrals.",      "frPosAdcPulseInt.THcSignalHit.GetData()"});
+  vars.push_back(RVarDef{"posAdcPulseAmp",     "List of positive ADC pulse amplitudes.",     "frPosAdcPulseAmp.THcSignalHit.GetData()"});
+
+  vars.push_back(RVarDef{"negAdcPedRaw",       "List of negative raw ADC pedestals",         "frNegAdcPedRaw.THcSignalHit.GetData()"});
+  vars.push_back(RVarDef{"negAdcPulseIntRaw",  "List of negative raw ADC pulse integrals.",  "frNegAdcPulseIntRaw.THcSignalHit.GetData()"});
+  vars.push_back(RVarDef{"negAdcPulseAmpRaw",  "List of negative raw ADC pulse amplitudes.", "frNegAdcPulseAmpRaw.THcSignalHit.GetData()"});
+  vars.push_back(RVarDef{"negAdcPulseTimeRaw", "List of negative raw ADC pulse times.",      "frNegAdcPulseTimeRaw.THcSignalHit.GetData()"});
+
+  vars.push_back(RVarDef{"negAdcPed",          "List of negative ADC pedestals",             "frNegAdcPed.THcSignalHit.GetData()"});
+  vars.push_back(RVarDef{"negAdcPulseInt",     "List of negative ADC pulse integrals.",      "frNegAdcPulseInt.THcSignalHit.GetData()"});
+  vars.push_back(RVarDef{"negAdcPulseAmp",     "List of negative ADC pulse amplitudes.",     "frNegAdcPulseAmp.THcSignalHit.GetData()"});
+
+   }
+
+  RVarDef end {0};
+  vars.push_back(end);
+
+  return DefineVarsFromList(vars.data(), mode);
 }
 
 //_____________________________________________________________________________
 void THcShowerPlane::Clear( Option_t* )
 {
   // Clears the hit lists
+
   fPosADCHits->Clear();
   fNegADCHits->Clear();
 
   frPosAdcErrorFlag->Clear();
   frPosAdcPedRaw->Clear();
+  frPosAdcThreshold->Clear();
   frPosAdcPulseIntRaw->Clear();
   frPosAdcPulseAmpRaw->Clear();
   frPosAdcPulseTimeRaw->Clear();
@@ -307,6 +366,7 @@ void THcShowerPlane::Clear( Option_t* )
 
   frNegAdcErrorFlag->Clear();
   frNegAdcPedRaw->Clear();
+  frNegAdcThreshold->Clear();
   frNegAdcPulseIntRaw->Clear();
   frNegAdcPulseAmpRaw->Clear();
   frNegAdcPulseTimeRaw->Clear();
@@ -315,7 +375,40 @@ void THcShowerPlane::Clear( Option_t* )
   frNegAdcPulseInt->Clear();
   frNegAdcPulseAmp->Clear();
 
-  // Debug output.
+  for (UInt_t ielem = 0; ielem < fGoodPosAdcPed.size(); ielem++) {
+    fGoodPosAdcPed.at(ielem)              = 0.0;
+    fGoodPosAdcPulseIntRaw.at(ielem)      = 0.0;
+    fGoodPosAdcPulseInt.at(ielem)         = 0.0;
+    fGoodPosAdcPulseAmp.at(ielem)         = 0.0;
+    fGoodPosAdcPulseTime.at(ielem)        = 0.0;
+    fEpos.at(ielem)                       = 0.0;
+    fNumGoodPosAdcHits.at(ielem)          = 0.0;
+  }
+
+  for (UInt_t ielem = 0; ielem < fGoodNegAdcPed.size(); ielem++) {
+    fGoodNegAdcPed.at(ielem)              = 0.0;
+    fGoodNegAdcPulseIntRaw.at(ielem)      = 0.0;
+    fGoodNegAdcPulseInt.at(ielem)         = 0.0;
+    fGoodNegAdcPulseAmp.at(ielem)         = 0.0;
+    fGoodNegAdcPulseTime.at(ielem)        = 0.0;
+    fEneg.at(ielem)                       = 0.0;
+    fNumGoodNegAdcHits.at(ielem)          = 0.0;
+  }
+
+  for (UInt_t ielem = 0; ielem < fEmean.size(); ielem++) {
+    fEmean.at(ielem)                                = 0.0;
+  }
+
+  fTotNumAdcHits       = 0;
+  fTotNumPosAdcHits    = 0;
+  fTotNumNegAdcHits    = 0;
+
+  fTotNumGoodAdcHits       = 0;
+  fTotNumGoodPosAdcHits    = 0;
+  fTotNumGoodNegAdcHits    = 0;
+
+
+ // Debug output.
   if ( ((THcShower*) GetParent())->fdbg_decoded_cal ) {
     cout << "---------------------------------------------------------------\n";
     cout << "Debug output from THcShowerPlane::Clear for "
@@ -374,6 +467,7 @@ Int_t THcShowerPlane::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
 
   frPosAdcErrorFlag->Clear();
   frPosAdcPedRaw->Clear();
+  frPosAdcThreshold->Clear();
   frPosAdcPulseIntRaw->Clear();
   frPosAdcPulseAmpRaw->Clear();
   frPosAdcPulseTimeRaw->Clear();
@@ -384,6 +478,7 @@ Int_t THcShowerPlane::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
 
   frNegAdcErrorFlag->Clear();
   frNegAdcPedRaw->Clear();
+  frNegAdcThreshold->Clear();
   frNegAdcPulseIntRaw->Clear();
   frNegAdcPulseAmpRaw->Clear();
   frNegAdcPulseTimeRaw->Clear();
@@ -392,15 +487,14 @@ Int_t THcShowerPlane::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
   frNegAdcPulseInt->Clear();
   frNegAdcPulseAmp->Clear();
 
-  for(Int_t i=0;i<fNelem;i++) {
-    fA_Pos[i] = 0;
-    fA_Neg[i] = 0;
-    fA_Pos_p[i] = 0;
-    fA_Neg_p[i] = 0;
+  /*
+    for(Int_t i=0;i<fNelem;i++) {
+
     fEpos[i] = 0;
     fEneg[i] = 0;
     fEmean[i] = 0;
   }
+  */
 
   fEplane = 0;
   fEplane_pos = 0;
@@ -430,8 +524,7 @@ Int_t THcShowerPlane::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
     THcRawAdcHit& rawPosAdcHit = hit->GetRawAdcHitPos();
     for (UInt_t thit=0; thit<rawPosAdcHit.GetNPulses(); ++thit) {
       ((THcSignalHit*) frPosAdcPedRaw->ConstructedAt(nrPosAdcHits))->Set(padnum, rawPosAdcHit.GetPedRaw());
-       fPosThresh[hit->fCounter -1]=rawPosAdcHit.GetPedRaw()*rawPosAdcHit.GetF250_PeakPedestalRatio()+fAdcPosThreshold;
-       
+      ((THcSignalHit*) frPosAdcThreshold->ConstructedAt(nrPosAdcHits))->Set(padnum,rawPosAdcHit.GetPedRaw()*rawPosAdcHit.GetF250_PeakPedestalRatio()+fAdcPosThreshold);
       ((THcSignalHit*) frPosAdcPed->ConstructedAt(nrPosAdcHits))->Set(padnum, rawPosAdcHit.GetPed());
 
       ((THcSignalHit*) frPosAdcPulseIntRaw->ConstructedAt(nrPosAdcHits))->Set(padnum, rawPosAdcHit.GetPulseIntRaw(thit));
@@ -442,16 +535,19 @@ Int_t THcShowerPlane::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
 
       ((THcSignalHit*) frPosAdcPulseTimeRaw->ConstructedAt(nrPosAdcHits))->Set(padnum, rawPosAdcHit.GetPulseTimeRaw(thit));
       if (rawPosAdcHit.GetPulseAmp(thit)>0&&rawPosAdcHit.GetPulseIntRaw(thit)>0) {
-	((THcSignalHit*) frPosAdcErrorFlag->ConstructedAt(nrPosAdcHits))->Set(padnum,0);        
+	((THcSignalHit*) frPosAdcErrorFlag->ConstructedAt(nrPosAdcHits))->Set(padnum,0);
       } else {
-	((THcSignalHit*) frPosAdcErrorFlag->ConstructedAt(nrPosAdcHits))->Set(padnum,1);        
+	((THcSignalHit*) frPosAdcErrorFlag->ConstructedAt(nrPosAdcHits))->Set(padnum,1);
       }
       ++nrPosAdcHits;
+      fTotNumAdcHits++;
+      fTotNumPosAdcHits++;
+
     }
     THcRawAdcHit& rawNegAdcHit = hit->GetRawAdcHitNeg();
     for (UInt_t thit=0; thit<rawNegAdcHit.GetNPulses(); ++thit) {
       ((THcSignalHit*) frNegAdcPedRaw->ConstructedAt(nrNegAdcHits))->Set(padnum, rawNegAdcHit.GetPedRaw());
-       fNegThresh[hit->fCounter -1]=rawNegAdcHit.GetPedRaw()*rawNegAdcHit.GetF250_PeakPedestalRatio()+fAdcNegThreshold;
+      ((THcSignalHit*) frNegAdcThreshold->ConstructedAt(nrNegAdcHits))->Set(padnum,rawNegAdcHit.GetPedRaw()*rawNegAdcHit.GetF250_PeakPedestalRatio()+fAdcNegThreshold);
       ((THcSignalHit*) frNegAdcPed->ConstructedAt(nrNegAdcHits))->Set(padnum, rawNegAdcHit.GetPed());
 
       ((THcSignalHit*) frNegAdcPulseIntRaw->ConstructedAt(nrNegAdcHits))->Set(padnum, rawNegAdcHit.GetPulseIntRaw(thit));
@@ -467,12 +563,11 @@ Int_t THcShowerPlane::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
 	((THcSignalHit*) frNegAdcErrorFlag->ConstructedAt(nrNegAdcHits))->Set(padnum,1);
       }
       ++nrNegAdcHits;
+      fTotNumAdcHits++;
+      fTotNumNegAdcHits++;
     }
-
     ihit++;
   }
-
-
   return(ihit);
 }
 //_____________________________________________________________________________
@@ -487,7 +582,7 @@ Int_t THcShowerPlane::CoarseProcessHits()
       FillADC_SampleIntegral();
     } else if (ADCMode == kADCSampIntDynPed) {
       FillADC_SampIntDynPed();
-        } else {
+    } else {
       FillADC_Standard();
     }
     //
@@ -552,26 +647,26 @@ void THcShowerPlane::FillADC_Standard()
   for (Int_t ielem=0;ielem<frNegAdcPulseIntRaw->GetEntries();ielem++) {
     Int_t npad = ((THcSignalHit*) frNegAdcPulseIntRaw->ConstructedAt(ielem))->GetPaddleNumber() - 1;
     Double_t pulseIntRaw = ((THcSignalHit*) frNegAdcPulseIntRaw->ConstructedAt(ielem))->GetData();
-      fA_Neg[npad] =pulseIntRaw;
-      if(fA_Neg[npad] >  fNegThresh[npad]) {
-       fA_Neg_p[npad] = pulseIntRaw-fNegPed[npad];
-       fEneg[npad] = fA_Neg_p[npad]*fParent->GetGain(npad,fLayerNum-1,1);
-       fEmean[npad] += fEneg[npad];
-       fEplane_neg += fEneg[npad];
+    fGoodNegAdcPulseIntRaw.at(npad) = pulseIntRaw;
+      if(fGoodNegAdcPulseIntRaw.at(npad) >  fNegThresh[npad]) {
+	fGoodNegAdcPulseInt.at(npad) = pulseIntRaw-fNegPed[npad];
+	fEneg.at(npad) = fGoodNegAdcPulseInt.at(npad)*fParent->GetGain(npad,fLayerNum-1,1);
+	fEmean.at(npad) += fEneg.at(npad);
+	fEplane_neg += fEneg.at(npad);
       }
   }
   for (Int_t ielem=0;ielem<frPosAdcPulseIntRaw->GetEntries();ielem++) {
     Int_t npad = ((THcSignalHit*) frPosAdcPulseIntRaw->ConstructedAt(ielem))->GetPaddleNumber() - 1;
     Double_t pulseIntRaw = ((THcSignalHit*) frPosAdcPulseIntRaw->ConstructedAt(ielem))->GetData();
-      fA_Pos[npad] =pulseIntRaw;
-      if(fA_Pos[npad] > fPosThresh[npad]) {
-       fA_Pos_p[npad] =pulseIntRaw-fPosPed[npad] ;
-       fEpos[npad] = fA_Pos_p[npad]*fParent->GetGain(npad,fLayerNum-1,0);
-       fEmean[npad] += fEpos[npad];
-       fEplane_pos += fEpos[npad];
-      }
+    fGoodPosAdcPulseIntRaw.at(npad) =pulseIntRaw;
+    if(fGoodPosAdcPulseIntRaw.at(npad) > fPosThresh[npad]) {
+      fGoodPosAdcPulseInt.at(npad) =pulseIntRaw-fPosPed[npad] ;
+      fEpos.at(npad) =fGoodPosAdcPulseInt.at(npad)*fParent->GetGain(npad,fLayerNum-1,0);
+      fEmean.at(npad) += fEpos.at(npad);
+      fEplane_pos += fEpos.at(npad);
+    }
   }
-    fEplane= fEplane_neg+fEplane_pos;
+  fEplane= fEplane_neg+fEplane_pos;
 }
 //_____________________________________________________________________________
 void THcShowerPlane::FillADC_DynamicPedestal()
@@ -581,40 +676,68 @@ void THcShowerPlane::FillADC_DynamicPedestal()
   Double_t AdcTimeWindowMin=fParent->GetAdcTimeWindowMin();
   Double_t AdcTimeWindowMax=fParent->GetAdcTimeWindowMax();
   for (Int_t ielem=0;ielem<frNegAdcPulseInt->GetEntries();ielem++) {
-    Int_t npad = ((THcSignalHit*) frNegAdcPulseInt->ConstructedAt(ielem))->GetPaddleNumber() - 1;
-    Double_t pulseInt = ((THcSignalHit*) frNegAdcPulseInt->ConstructedAt(ielem))->GetData();
-    Double_t pulseIntRaw = ((THcSignalHit*) frNegAdcPulseIntRaw->ConstructedAt(ielem))->GetData();
-    Double_t pulseTime = ((THcSignalHit*) frNegAdcPulseTimeRaw->ConstructedAt(ielem))->GetData();
-    Double_t errorflag = ((THcSignalHit*) frNegAdcErrorFlag->ConstructedAt(ielem))->GetData();
-    Bool_t pulseTimeCut = (pulseTime > AdcTimeWindowMin) &&  (pulseTime < AdcTimeWindowMax);
-    if (errorflag==0 && pulseTimeCut) {
-      fA_Neg[npad] =pulseIntRaw;
-      if(fA_Neg[npad] >  fNegThresh[npad]) {
-       fA_Neg_p[npad] =pulseInt ;
-       fEneg[npad] = fA_Neg_p[npad]*fParent->GetGain(npad,fLayerNum-1,1);
-       fEmean[npad] += fEneg[npad];
-       fEplane_neg += fEneg[npad];
+    Int_t    npad         = ((THcSignalHit*) frNegAdcPulseInt->ConstructedAt(ielem))->GetPaddleNumber() - 1;
+    Double_t pulseInt     = ((THcSignalHit*) frNegAdcPulseInt->ConstructedAt(ielem))->GetData();
+    Double_t pulsePed     = ((THcSignalHit*) frNegAdcPed->ConstructedAt(ielem))->GetData();
+    Double_t pulseAmp     = ((THcSignalHit*) frNegAdcPulseAmp->ConstructedAt(ielem))->GetData();
+    Double_t pulseIntRaw  = ((THcSignalHit*) frNegAdcPulseIntRaw->ConstructedAt(ielem))->GetData();
+    Double_t pulseTime    = ((THcSignalHit*) frNegAdcPulseTimeRaw->ConstructedAt(ielem))->GetData();
+    Double_t threshold    = ((THcSignalHit*) frNegAdcThreshold->ConstructedAt(ielem))->GetData();
+    Bool_t   errorflag    = ((THcSignalHit*) frNegAdcErrorFlag->ConstructedAt(ielem))->GetData();
+    Bool_t   pulseTimeCut = (pulseTime > AdcTimeWindowMin) &&  (pulseTime < AdcTimeWindowMax);
+    if (!errorflag && pulseTimeCut) {
+      fGoodNegAdcPulseIntRaw.at(npad) =pulseIntRaw;
+
+      if(fGoodNegAdcPulseIntRaw.at(npad) >  threshold && fGoodNegAdcPulseInt.at(npad)==0) {
+        fGoodNegAdcPulseInt.at(npad) =pulseInt ;
+	fEneg.at(npad) =  fGoodNegAdcPulseInt.at(npad)*fParent->GetGain(npad,fLayerNum-1,1);
+	fEmean.at(npad) += fEneg.at(npad);
+	fEplane_neg += fEneg.at(npad);
+
+	fGoodNegAdcPed.at(npad) = pulsePed;
+      fGoodNegAdcPulseAmp.at(npad) = pulseAmp;
+      fGoodNegAdcPulseTime.at(npad) = pulseTime;
+
+      fTotNumGoodAdcHits++;
+      fTotNumGoodNegAdcHits++;
+      fNumGoodNegAdcHits.at(npad)++;
+
       }
-     }        
+
     }
+  }
   //
   for (Int_t ielem=0;ielem<frPosAdcPulseInt->GetEntries();ielem++) {
-    Int_t npad = ((THcSignalHit*) frPosAdcPulseInt->ConstructedAt(ielem))->GetPaddleNumber() - 1;
-    Double_t pulseInt = ((THcSignalHit*) frPosAdcPulseInt->ConstructedAt(ielem))->GetData();
-    Double_t pulseIntRaw = ((THcSignalHit*) frPosAdcPulseIntRaw->ConstructedAt(ielem))->GetData();
-    Double_t pulseTime = ((THcSignalHit*) frPosAdcPulseTimeRaw->ConstructedAt(ielem))->GetData();
-    Double_t errorflag = ((THcSignalHit*) frPosAdcErrorFlag->ConstructedAt(ielem))->GetData();
-    Bool_t pulseTimeCut = pulseTime > AdcTimeWindowMin &&  pulseTime < AdcTimeWindowMax;
-    if (errorflag==0 && pulseTimeCut) {
-      fA_Pos[npad] =pulseIntRaw;
-      if(fA_Pos[npad] >  fPosThresh[npad]) {
-       fA_Pos_p[npad] =pulseInt ;
-       fEpos[npad] = fA_Pos_p[npad]*fParent->GetGain(npad,fLayerNum-1,0);
-       fEmean[npad] += fEpos[npad];
-       fEplane_pos += fEpos[npad];
+    Int_t    npad         = ((THcSignalHit*) frPosAdcPulseInt->ConstructedAt(ielem))->GetPaddleNumber() - 1;
+    Double_t pulsePed     = ((THcSignalHit*) frPosAdcPed->ConstructedAt(ielem))->GetData();
+    Double_t threshold    = ((THcSignalHit*) frPosAdcThreshold->ConstructedAt(ielem))->GetData();
+    Double_t pulseAmp     = ((THcSignalHit*) frPosAdcPulseAmp->ConstructedAt(ielem))->GetData();
+    Double_t pulseInt     = ((THcSignalHit*) frPosAdcPulseInt->ConstructedAt(ielem))->GetData();
+    Double_t pulseIntRaw  = ((THcSignalHit*) frPosAdcPulseIntRaw->ConstructedAt(ielem))->GetData();
+    Double_t pulseTime    = ((THcSignalHit*) frPosAdcPulseTimeRaw->ConstructedAt(ielem))->GetData();
+    Bool_t   errorflag    = ((THcSignalHit*) frPosAdcErrorFlag->ConstructedAt(ielem))->GetData();
+    Bool_t   pulseTimeCut = (pulseTime > AdcTimeWindowMin) &&  (pulseTime < AdcTimeWindowMax);
+    if (!errorflag && pulseTimeCut) {
+      fGoodPosAdcPulseIntRaw.at(npad) = pulseIntRaw;
+
+      if(fGoodPosAdcPulseIntRaw.at(npad) >  threshold && fGoodPosAdcPulseInt.at(npad)==0) {
+
+       	fGoodPosAdcPulseInt.at(npad) =pulseInt ;
+	fEpos.at(npad) = fGoodPosAdcPulseInt.at(npad)*fParent->GetGain(npad,fLayerNum-1,0);
+	fEmean.at(npad) += fEpos[npad];
+	fEplane_pos += fEpos.at(npad);
+
+	fGoodPosAdcPed.at(npad) = pulsePed;
+	fGoodPosAdcPulseAmp.at(npad) = pulseAmp;
+	fGoodPosAdcPulseTime.at(npad) = pulseTime;
+
+	fTotNumGoodAdcHits++;
+	fTotNumGoodPosAdcHits++;
+	fNumGoodPosAdcHits.at(npad)++;
+
       }
-     }        
     }
+  }
   //
     fEplane= fEplane_neg+fEplane_pos;
 
diff --git a/src/THcShowerPlane.h b/src/THcShowerPlane.h
index 1c84d1a..6a029e1 100644
--- a/src/THcShowerPlane.h
+++ b/src/THcShowerPlane.h
@@ -16,9 +16,11 @@
 #include "TClonesArray.h"
 
 #include <iostream>
-
+#include <vector>
 #include <fstream>
 
+using namespace std;
+
 class THaEvData;
 class THaSignalHit;
 
@@ -76,19 +78,19 @@ public:
   };
 
   Double_t GetAposP(Int_t i) {
-    return fA_Pos_p[i];
+    return fGoodPosAdcPulseInt[i];
   };
 
   Double_t GetAnegP(Int_t i) {
-    return fA_Neg_p[i];
+    return fGoodNegAdcPulseInt[i];
   };
 
   Double_t GetApos(Int_t i) {
-    return fA_Pos[i];
+    return fGoodPosAdcPulseIntRaw[i];
   };
 
   Double_t GetAneg(Int_t i) {
-    return fA_Neg[i];
+    return fGoodNegAdcPulseIntRaw[i];
   };
 
   Double_t GetPosThr(Int_t i) {
@@ -114,28 +116,54 @@ protected:
    //  1 == Use the pulse int - pulse ped
     //  2 == Use the sample integral - known ped
     //  3 == Use the sample integral - sample ped
- static const Int_t kADCStandard=0;
+  static const Int_t kADCStandard=0;
   static const Int_t kADCDynamicPedestal=1;
   static const Int_t kADCSampleIntegral=2;
   static const Int_t kADCSampIntDynPed=3;
-   Int_t fPedSampLow;		// Sample range for
+
+  Int_t fDebugAdc;              // fADC debug flag
+  Int_t fPedSampLow;		// Sample range for
   Int_t fPedSampHigh;		// dynamic pedestal
   Int_t fDataSampLow;		// Sample range for
   Int_t fDataSampHigh;		// sample integration
-  Double_t fAdcNegThreshold;		// 
-  Double_t fAdcPosThreshold;		// 
-
-  Double_t*   fA_Pos;         // [fNelem] ADC amplitudes of blocks
-  Double_t*   fA_Neg;         // [fNelem] ADC amplitudes of blocks
-  Double_t*   fA_Pos_p;	      // [fNelem] pedestal subtracted ADC amplitudes
-  Double_t*   fA_Neg_p;	      // [fNelem] pedestal subtracted ADC amplitudes
-
-  Double_t* fEpos;     // [fNelem] energy depositions seen by positive PMTs
-  Double_t* fEneg;     // [fNelem] energy depositions seen by negative PMTs
-  Double_t* fEmean;    // [fNelem] mean energy depositions (pos + neg)
-  Double_t  fEplane;   // Energy deposition in the plane
+  Double_t fAdcNegThreshold;		//
+  Double_t fAdcPosThreshold;		//
+
+  //counting variables
+  Int_t     fTotNumPosAdcHits;
+  Int_t     fTotNumNegAdcHits;
+  Int_t     fTotNumAdcHits;
+
+  Int_t     fTotNumGoodPosAdcHits;
+  Int_t     fTotNumGoodNegAdcHits;
+  Int_t     fTotNumGoodAdcHits;
+
+   //individual pmt data objects
+  vector<Int_t>    fNumGoodPosAdcHits;
+  vector<Int_t>    fNumGoodNegAdcHits;
+
+  vector<Double_t>      fGoodPosAdcPed;
+  vector<Double_t>      fGoodPosAdcPulseInt;
+  vector<Double_t>      fGoodPosAdcPulseAmp;
+  vector<Double_t>      fGoodPosAdcPulseTime;
+
+  vector<Double_t>      fGoodNegAdcPed;
+  vector<Double_t>      fGoodNegAdcPulseInt;
+  vector<Double_t>      fGoodNegAdcPulseAmp;
+  vector<Double_t>      fGoodNegAdcPulseTime;
+
+  vector<Double_t>      fGoodPosAdcPulseIntRaw;
+  vector<Double_t>      fGoodNegAdcPulseIntRaw;
+
+
+  vector<Double_t>      fEpos;        // [fNelem] energy depositions seen by positive PMTs
+  vector<Double_t>      fEneg;        // [fNelem] energy depositions seen by negative PMTs
+  vector<Double_t>      fEmean;        // [fNelem] mean energy depositions (pos + neg)
+
+
   Double_t  fEplane_pos;   // Energy deposition in the plane from positive PMTs
   Double_t  fEplane_neg;   // Energy deposition in the plane from negative PMTs
+  Double_t  fEplane;
 
   // These lists are not used actively for now.
   TClonesArray* fPosADCHits;    // List of positive ADC hits
@@ -163,6 +191,7 @@ protected:
 
   TClonesArray* frPosAdcErrorFlag;
   TClonesArray* frPosAdcPedRaw;
+  TClonesArray* frPosAdcThreshold;
   TClonesArray* frPosAdcPulseIntRaw;
   TClonesArray* frPosAdcPulseAmpRaw;
   TClonesArray* frPosAdcPulseTimeRaw;
@@ -173,6 +202,7 @@ protected:
 
   TClonesArray* frNegAdcErrorFlag;
   TClonesArray* frNegAdcPedRaw;
+  TClonesArray* frNegAdcThreshold;
   TClonesArray* frNegAdcPulseIntRaw;
   TClonesArray* frNegAdcPulseAmpRaw;
   TClonesArray* frNegAdcPulseTimeRaw;
-- 
GitLab