Add dynamic pedestal subtraction to THcRawShowerHit class

  GetData accepts a pedestal range and an integration range.
  Samples in the integration range are added after subtracting the
  average of the pedestal range from each sample.
  Method also provided to return dynamic pedestal

src/THcRawShowerHit.cxx

@@ -6,10 +6,18 @@
 //                                                                           //
 // Contains plane, counter and pos/neg adc and tdc values                    //
 //                                                                           //
+// Enhanced to work with FADC250 data samples.  If fNPosSamples/fNNegSamples //
+// is greater than 1, assume that the data held in the hit is the sampled    //
+// waveform.  Signals 0,1 will return the integrated pulse with dynamic      //
+// pedestal subtraction (first four samples comprise the pedestal).  Signals //
+// 2,3 are reserved for time information.  Signals 4,5 are pedestals and     //
+// 6 and 7 are the straight sum of all the samples.                          //
+//                                                                           //
 ///////////////////////////////////////////////////////////////////////////////
 
 #include "THcRawShowerHit.h"
 #include <iostream>
+#include <cassert>
 
 using namespace std;
 
@@ -22,9 +30,56 @@ void THcRawShowerHit::SetData(Int_t signal, Int_t data) {
   } 
 }
 
+Int_t THcRawShowerHit::GetData(Int_t signal, Int_t isamplow, Int_t isamphigh,
+	      Int_t iintegrallow, Int_t iintegralhigh) {
+  Int_t adcsum=0;
+  Double_t pedestal=0.0;
+
+  if(signal==0 || signal == 1) {
+    pedestal = GetPedestal(signal, isamplow, isamphigh);
+
+#if 0
+    for(Int_t i=0;i<fNPosSamples;i++) {
+      cout << fCounter << " " << i;
+      if(i >= isamplow && i<=isamphigh) {
+	cout << "P ";
+      }	else if (i >= iintegrallow && i<=iintegralhigh) {
+	cout << "D ";
+      } else {
+	cout << "  ";
+      }
+      cout << fADC_pos[i] << " " << fADC_pos[i] - pedestal << endl;
+    }
+#endif
+  }
+  if(signal==4 || signal==5) {
+    pedestal = GetPedestal(signal-4, isamplow, isamphigh);
+    return(pedestal);
+  }
+  if(signal==0 || signal==6) {
+    assert(iintegralhigh<(Int_t) fNPosSamples);
+    for(UInt_t isample=iintegrallow;isample<=iintegralhigh;isample++) {
+      adcsum += fADC_pos[isample] - pedestal;
+    }
+    return(adcsum);
+  } else if (signal==1 || signal==7) {
+    assert(iintegralhigh<(Int_t) fNNegSamples);
+    for(UInt_t isample=iintegrallow;isample<=iintegralhigh;isample++) {
+      adcsum += fADC_neg[isample] - pedestal;
+    }
+    return(adcsum);
+  } 
+  return(-1); // Actually should throw exception
+}
+  
 // Return sum of samples
+// For Fastbus ADC, this will simply be the hardware ADC value as there
+// is just one sample.  For Flash ADCs, this should return the
+// integrated ADC value if the FADC provided that, otherwise the sum
+// of all the samples
 Int_t THcRawShowerHit::GetData(Int_t signal) {
   Int_t adcsum=0;
+
   if(signal==0) {
     for(UInt_t isample=0;isample<fNPosSamples;isample++) {
       adcsum += fADC_pos[isample];
@@ -40,7 +95,7 @@ Int_t THcRawShowerHit::GetData(Int_t signal) {
 }
 
 // Return a requested sample
-Int_t THcRawShowerHit::GetData(Int_t signal, UInt_t isample) {
+Int_t THcRawShowerHit::GetSample(Int_t signal, UInt_t isample) {
   if(signal==0) {
     if(isample >=0 && isample< fNPosSamples) {
       return(fADC_pos[isample]);
@@ -53,6 +108,27 @@ Int_t THcRawShowerHit::GetData(Int_t signal, UInt_t isample) {
   return(-1);
 }
 
+Double_t THcRawShowerHit::GetPedestal(Int_t signal, Int_t isamplow, Int_t isamphigh) {
+  // No error checking on pedestal range
+  Double_t pedestal=0.0;
+  if(signal==0 && fNPosSamples > 1) {
+    if(fNPosSamples > isamphigh) {
+      for(Int_t i = isamplow;i<=isamphigh;i++) {
+	pedestal += fADC_pos[i];
+      }
+      return(pedestal/(isamphigh-isamplow+1));
+    }
+  } else if (signal==1 && fNNegSamples > 1) {
+    if(fNNegSamples > isamphigh) {
+      for(Int_t i = isamplow;i<=isamphigh;i++) {
+	pedestal += fADC_neg[i];
+      }
+      return(pedestal/(isamphigh-isamplow+1));
+    }
+  }
+  return(pedestal);
+}
+      
 // Return the number of samples
 Int_t THcRawShowerHit::GetNSamples(Int_t signal) {
   if(signal==0) {

src/THcRawShowerHit.h

@@ -22,7 +22,10 @@ class THcRawShowerHit : public THcRawHit {
 
   void SetData(Int_t signal, Int_t data);
   Int_t GetData(Int_t signal);
-  Int_t GetData(Int_t signal, UInt_t isample);
+  Int_t GetData(Int_t signal, Int_t isamplow, Int_t isamphigh,
+		Int_t iintegrallow, Int_t iintegralhigh);
+  Int_t GetSample(Int_t signal, UInt_t isample);
+  Double_t GetPedestal(Int_t signal, Int_t isamplow, Int_t isamphigh);
   Int_t GetNSamples(Int_t signal);
 
   //  virtual Bool_t  IsSortable () const {return kTRUE; }

src/THcShowerArray.cxx

@@ -46,6 +46,7 @@ THcShowerArray::~THcShowerArray()
   delete fADCHits;
 
   delete [] fA;
+  delete [] fP;
 }
 
 //_____________________________________________________________________________
@@ -77,9 +78,19 @@ Int_t THcShowerArray::ReadDatabase( const TDatime& date )
   prefix[1]='\0';
 
   cout << "Parent name: " << GetParent()->GetPrefix() << endl;
+  fUsingFADC=0;
+  fPedSampLow=0;
+  fPedSampHigh=9;
+  fDataSampLow=23;
+  fDataSampHigh=49;
   DBRequest list[]={
     {"cal_nrows", &fNRows, kInt},
     {"cal_ncolumns", &fNColumns, kInt},
+    {"cal_using_fadc", &fUsingFADC, kInt, 0, 1},
+    {"cal_ped_sample_low", &fPedSampLow, kInt, 0, 1},
+    {"cal_ped_sample_high", &fPedSampHigh, kInt, 0, 1},
+    {"cal_data_sample_low", &fDataSampLow, kInt, 0, 1},
+    {"cal_data_sample_high", &fDataSampHigh, kInt, 0, 1},
     {0}
   };
   gHcParms->LoadParmValues((DBRequest*)&list, prefix);
@@ -89,7 +100,9 @@ Int_t THcShowerArray::ReadDatabase( const TDatime& date )
   // Here read the 2-D arras of pedestals, gains, etc.
 
 
+  // Event by event amplitude and pedestal
   fA = new Double_t[fNelem];
+  fP = new Double_t[fNelem];
 
   return kOK;
 }
@@ -106,6 +119,7 @@ Int_t THcShowerArray::DefineVariables( EMode mode )
   RVarDef vars[] = {
     {"adchits", "List of ADC hits", "fADCHits.THcSignalHit.GetPaddleNumber()"},
     {"a", "Raw ADC Amplitude", "fA"},
+    {"p", "Dynamic ADC Pedestal", "fP"},
     { 0 }
   };
 
@@ -173,7 +187,13 @@ Int_t THcShowerArray::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
     THcRawShowerHit* hit = (THcRawShowerHit *) rawhits->At(ihit);
 
     // Should probably check that counter # is in range
-    fA[hit->fCounter-1] = hit->GetData(0);
+    if(fUsingFADC) {
+      fA[hit->fCounter-1] = hit->GetData(0,fPedSampLow,fPedSampHigh,
+					 fDataSampLow,fDataSampHigh);
+      fP[hit->fCounter-1] = hit->GetPedestal(0,fPedSampLow,fPedSampHigh);
+    } else {
+          fA[hit->fCounter-1] = hit->GetData(0);
+    }
 
     // Do other stuff like comparison to thresholds, signal hits, energy sums
 

src/THcShowerArray.h

@@ -55,6 +55,7 @@ public:
 protected:
 
   Double_t* fA;                 // [fNelem] ADC amplitude of blocks
+  Double_t* fP;                 // [fNelem] Event by event pedestals
 
   TClonesArray* fADCHits;	// List of ADC hits
 
@@ -64,6 +65,11 @@ protected:
   // Parameters
   Int_t fNRows;
   Int_t fNColumns;
+  Int_t fUsingFADC;		// != 0 if using FADC in sample mode
+  Int_t fPedSampLow;		// Sample range for
+  Int_t fPedSampHigh;		// dynamic pedestal
+  Int_t fDataSampLow;		// Sample range for
+  Int_t fDataSampHigh;		// sample integration
 
   Int_t fLayerNum;		// 2 for SHMS
   // Accumulators for pedestals go here