Very basic support for the SHMS shower counter

  Added as an optional sub detector THcShowerArray to the THcShower class
  Get hit data with GetData method instead of hit member variables
    so as to be compatible with FADC compatible hit list

Makefile

@@ -21,6 +21,7 @@ SRC  =  src/THcInterface.cxx src/THcParmList.cxx src/THcAnalyzer.cxx \
 	src/THcDCLookupTTDConv.cxx src/THcDCTimeToDistConv.cxx \
 	src/THcSpacePoint.cxx src/THcDCTrack.cxx \
 	src/THcShower.cxx src/THcShowerPlane.cxx \
+	src/THcShowerArray.cxx \
 	src/THcRawShowerHit.cxx \
 	src/THcAerogel.cxx src/THcAerogelHit.cxx \
 	src/THcCherenkov.cxx src/THcCherenkovHit.cxx \

SConscript.py

@@ -17,6 +17,7 @@ hcheaders = Split("""
         src/THcDC.h src/THcDriftChamberPlane.h 
 	src/THcDriftChamber.h src/THcRawDCHit.h src/THcDCHit.h src/THcDCWire.h src/THcSpacePoint.h 
 	src/THcDCLookupTTDConv.h src/THcDCTimeToDistConv.h src/THcShower.h src/THcShowerPlane.h 
+        src/THcShowerArray.h
 	src/THcRawShowerHit.h src/THcAerogel.h src/THcAerogelHit.h src/THcCherenkov.h src/THcCherenkovHit.h
         src/THcGlobals.h src/THcDCTrack.h src/THcFormula.h
         src/THcRaster.h src/THcRasteredBeam.h src/THcRasterRawHit.h src/THcScalerEvtHandler.h

src/HallC_LinkDef.h

@@ -43,6 +43,7 @@
 #pragma link C++ class THcDCTrack+;
 #pragma link C++ class THcShower+;
 #pragma link C++ class THcShowerPlane+;
+#pragma link C++ class THcShowerArray+;
 #pragma link C++ class THcRawShowerHit+;
 #pragma link C++ class THcAerogel+;
 #pragma link C++ class THcAerogelHit+;

src/THcRawShowerHit.h

@@ -7,6 +7,7 @@ class THcRawShowerHit : public THcRawHit {
 
  public:
   friend class THcShowerPlane;
+  friend class THcShowerArray;
 
   THcRawShowerHit(Int_t plane=0, Int_t counter=0) : THcRawHit(plane, counter), 
     fADC_pos(-1), fADC_neg(-1){

src/THcShower.cxx

@@ -39,6 +39,8 @@ THcShower::THcShower( const char* name, const char* description,
 {
   // Constructor
   fNLayers = 0;			// No layers until we make them
+  fNTotLayers = 0;
+  fHasArray = 0;
 
   fClusterList = new THcShowerClusterList;
 }
@@ -58,26 +60,29 @@ void THcShower::Setup(const char* name, const char* description)
   prefix[0] = tolower(GetApparatus()->GetName()[0]);
   prefix[1] = '\0';
 
+  
   string layernamelist;
+  fHasArray = 0;		// Flag for presence of fly's eye array
   DBRequest list[]={
     {"cal_num_layers", &fNLayers, kInt},
     {"cal_layer_names", &layernamelist, kString},
+    {"cal_array",&fHasArray, kInt,0, 1}, 
     {0}
   };
 
   gHcParms->LoadParmValues((DBRequest*)&list,prefix);
-
+  fNTotLayers = (fNLayers+(fHasArray!=0?1:0));
   vector<string> layer_names = vsplit(layernamelist);
 
-  if(layer_names.size() != (UInt_t) fNLayers) {
-    cout << "THcShower::Setup ERROR: Number of layers " << fNLayers 
+  if(layer_names.size() != fNTotLayers) {
+    cout << "THcShower::Setup ERROR: Number of layers " << fNTotLayers 
 	 << " doesn't agree with number of layer names "
 	 << layer_names.size() << endl;
     // Should quit.  Is there an official way to quit?
   }
 
-  fLayerNames = new char* [fNLayers];
-  for(UInt_t i=0;i<fNLayers;i++) {
+  fLayerNames = new char* [fNTotLayers];
+  for(UInt_t i=0;i<fNTotLayers;i++) {
     fLayerNames[i] = new char[layer_names[i].length()+1];
     strcpy(fLayerNames[i], layer_names[i].c_str());
   }
@@ -93,6 +98,16 @@ void THcShower::Setup(const char* name, const char* description)
     fPlanes[i] = new THcShowerPlane(fLayerNames[i], desc, i+1, this); 
 
   }
+  if(fHasArray) {
+    strcpy(desc, description);
+    strcat(desc, " Array ");
+    strcat(desc, fLayerNames[fNTotLayers-1]);
+
+    fArray = new THcShowerArray(fLayerNames[fNTotLayers-1], desc, fNTotLayers, this);
+  } else {
+    fArray = 0;
+  }
+
   delete [] desc;
 
   cout << "---------------------------------------------------------------\n";
@@ -127,6 +142,11 @@ THaAnalysisObject::EStatus THcShower::Init( const TDatime& date )
       return fStatus=status;
     }
   }
+  if(fHasArray) {
+    if((status = fArray->Init( date ))) {
+      return fStatus = status;
+    }
+  }
 
   char EngineDID[] = " CAL";
   EngineDID[0] = toupper(GetApparatus()->GetName()[0]);
@@ -171,9 +191,11 @@ Int_t THcShower::ReadDatabase( const TDatime& date )
   prefix[0]=tolower(GetApparatus()->GetName()[0]);
   prefix[1]='\0';
 
+  fNegCols = fNLayers;		// If not defined assume tube on each end
+                                // for every layer
   {
     DBRequest list[]={
-      {"cal_num_neg_columns", &fNegCols, kInt},
+      {"cal_num_neg_columns", &fNegCols, kInt, 0, 1},
       {"cal_slop", &fSlop, kDouble},
       {"cal_fv_test", &fvTest, kInt,0,1},
       {"cal_fv_delta", &fvDelta, kDouble},
@@ -518,6 +540,9 @@ void THcShower::Clear(Option_t* opt)
   for(UInt_t ip=0;ip<fNLayers;ip++) {
     fPlanes[ip]->Clear(opt);
   }
+  if(fHasArray) {
+    fArray->Clear(opt);
+  }
 
   fNhits = 0;
   fNclust = 0;
@@ -552,6 +577,9 @@ Int_t THcShower::Decode( const THaEvData& evdata )
     for(UInt_t ip=0;ip<fNLayers;ip++) {
       nexthit = fPlanes[ip]->AccumulatePedestals(fRawHitList, nexthit);
     }
+    if(fHasArray) {
+      nexthit = fArray->AccumulatePedestals(fRawHitList, nexthit);
+    }
     fAnalyzePedestals = 1;	// Analyze pedestals first normal events
     return(0);
   }
@@ -560,6 +588,9 @@ Int_t THcShower::Decode( const THaEvData& evdata )
     for(UInt_t ip=0;ip<fNLayers;ip++) {
       fPlanes[ip]->CalculatePedestals();
     }
+    if(fHasArray) {
+      fArray->CalculatePedestals();
+    }
     fAnalyzePedestals = 0;	// Don't analyze pedestals next event
   }
 
@@ -568,6 +599,10 @@ Int_t THcShower::Decode( const THaEvData& evdata )
     nexthit = fPlanes[ip]->ProcessHits(fRawHitList, nexthit);
     fEtot += fPlanes[ip]->GetEplane();
   }
+  if(fHasArray) {
+    nexthit = fArray->ProcessHits(fRawHitList, nexthit);
+    fEtot += fArray->GetEplane();
+  }
   THcHallCSpectrometer *app = static_cast<THcHallCSpectrometer*>(GetApparatus());
   fEtotNorm=fEtot/(app->GetPcentral());
  
@@ -588,6 +623,8 @@ Int_t THcShower::CoarseProcess( TClonesArray& tracks)
 
   // Fill set of unclustered hits.
 
+  // Ignore shower array (SHMS) for now
+
   THcShowerHitSet HitSet;
 
   for(UInt_t j=0; j < fNLayers; j++) {

src/THcShower.h

@@ -11,6 +11,7 @@
 #include "THaNonTrackingDetector.h"
 #include "THcHitList.h"
 #include "THcShowerPlane.h"
+#include "THcShowerArray.h"
 #include "TMath.h"
 
 
@@ -256,7 +257,10 @@ protected:
 
   char** fLayerNames;
   UInt_t fNLayers;	        // Number of layers in the calorimeter
+  UInt_t fNTotLayers;	        // Number of layers including array
+  UInt_t fHasArray;		// If !=0 fly's eye array behind preshower
   Double_t* fNLayerZPos;	// Z positions of fronts of layers
+  // Following apply to just sideways readout layers
   Double_t* BlockThick;		// Thickness of blocks
   UInt_t* fNBlocks;              // [fNLayers] number of blocks per layer
   UInt_t fNtotBlocks;            // Total number of shower counter blocks
@@ -287,6 +291,7 @@ protected:
   Double_t fDcor[2];
 
   THcShowerPlane** fPlanes;     // [fNLayers] Shower Plane objects
+  THcShowerArray* fArray;
 
   TClonesArray*  fTrackProj;    // projection of track onto plane
 

src/THcShowerArray.cxx

+//*-- Author :
+
+//////////////////////////////////////////////////////////////////////////
+//
+// THcShowerArray
+//
+//////////////////////////////////////////////////////////////////////////
+
+#include "THcShowerArray.h"
+#include "TClonesArray.h"
+#include "THcSignalHit.h"
+#include "THcGlobals.h"
+#include "THcParmList.h"
+#include "THcHitList.h"
+#include "THcShower.h"
+#include "THcRawShowerHit.h"
+#include "TClass.h"
+#include "math.h"
+#include "THaTrack.h"
+#include "THaTrackProj.h"
+
+#include <cstring>
+#include <cstdio>
+#include <cstdlib>
+#include <iostream>
+
+#include <fstream>
+using namespace std;
+
+ClassImp(THcShowerArray)
+
+//______________________________________________________________________________
+THcShowerArray::THcShowerArray( const char* name, 
+				const char* description,
+				const Int_t layernum,
+				THaDetectorBase* parent )
+  : THaSubDetector(name,description,parent)
+{
+  fADCHits = new TClonesArray("THcSignalHit",100);
+  fLayerNum = layernum;
+}
+
+THcShowerArray::~THcShowerArray()
+{
+  // Destructor
+  delete fADCHits;
+
+  delete [] fA;
+}
+
+//_____________________________________________________________________________
+THaAnalysisObject::EStatus THcShowerArray::Init( const TDatime& date )
+{
+  // Extra initialization for shower layer: set up DataDest map
+
+  if( IsZombie())
+    return fStatus = kInitError;
+
+  // How to get information for parent
+  //  if( GetParent() )
+  //    fOrigin = GetParent()->GetOrigin();
+
+  EStatus status;
+  if( (status=THaSubDetector::Init( date )) )
+    return fStatus = status;
+
+  return fStatus = kOK;
+
+}
+
+//_____________________________________________________________________________
+Int_t THcShowerArray::ReadDatabase( const TDatime& date )
+{
+
+  char prefix[2];
+  prefix[0]=tolower(GetParent()->GetPrefix()[0]);
+  prefix[1]='\0';
+
+  cout << "Parent name: " << GetParent()->GetPrefix() << endl;
+  DBRequest list[]={
+    {"cal_nrows", &fNRows, kInt},
+    {"cal_ncolumns", &fNColumns, kInt},
+    {0}
+  };
+  gHcParms->LoadParmValues((DBRequest*)&list, prefix);
+
+  fNelem = fNRows*fNColumns;
+
+  // Here read the 2-D arras of pedestals, gains, etc.
+
+
+  fA = new Double_t[fNelem];
+
+  return kOK;
+}
+
+//_____________________________________________________________________________
+Int_t THcShowerArray::DefineVariables( EMode mode )
+{
+  // Initialize global variables
+
+  if( mode == kDefine && fIsSetup ) return kOK;
+  fIsSetup = ( mode == kDefine );
+
+  // Register variables in global list
+  RVarDef vars[] = {
+    {"adchits", "List of ADC hits", "fADCHits.THcSignalHit.GetPaddleNumber()"},
+    {"a", "Raw ADC Amplitude", "fA"},
+    { 0 }
+  };
+
+  return DefineVarsFromList( vars, mode );
+}
+
+//_____________________________________________________________________________
+void THcShowerArray::Clear( Option_t* )
+{
+  // Clears the hit lists
+  fADCHits->Clear();
+
+}
+
+//_____________________________________________________________________________
+Int_t THcShowerArray::Decode( const THaEvData& evdata )
+{
+  // Doesn't actually get called.  Use Fill method instead
+
+  return 0;
+}
+
+//_____________________________________________________________________________
+Int_t THcShowerArray::CoarseProcess( TClonesArray& tracks )
+{
+
+  // Nothing is done here. See ProcessHits method instead.
+  //  
+
+ return 0;
+}
+
+//_____________________________________________________________________________
+Int_t THcShowerArray::FineProcess( TClonesArray& tracks )
+{
+  return 0;
+}
+
+//_____________________________________________________________________________
+Int_t THcShowerArray::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
+{
+  // Extract the data for this layer from hit list
+  // Assumes that the hit list is sorted by layer, so we stop when the
+  // plane doesn't agree and return the index for the next hit.
+
+  THcShower* fParent;
+  fParent = (THcShower*) GetParent();
+
+  // Initialize variables.
+
+  fADCHits->Clear();
+
+  for(Int_t i=0;i<fNelem;i++) {
+    fA[i] = 0;
+  }
+
+  // Process raw hits. Get ADC hits for the plane, assign variables for each
+  // channel.
+
+  Int_t nrawhits = rawhits->GetLast()+1;
+
+  Int_t ihit = nexthit;
+
+  while(ihit < nrawhits) {
+    THcRawShowerHit* hit = (THcRawShowerHit *) rawhits->At(ihit);
+
+    // Should probably check that counter # is in range
+    fA[hit->fCounter-1] = hit->fADC_pos;
+
+    // Do other stuff like comparison to thresholds, signal hits, energy sums
+
+    ihit++;
+  }
+
+  return(ihit);
+}
+//_____________________________________________________________________________
+Int_t THcShowerArray::AccumulatePedestals(TClonesArray* rawhits, Int_t nexthit)
+{
+  // Doesn't do anything yet except skip over hits
+
+  Int_t nrawhits = rawhits->GetLast()+1;
+
+  Int_t ihit = nexthit;
+
+  while(ihit < nrawhits) {
+    THcRawShowerHit* hit = (THcRawShowerHit *) rawhits->At(ihit);
+
+    // OK for hit list sorted by layer.
+    if(hit->fPlane > fLayerNum) {
+      break;
+    }
+    ihit++;
+  }
+  fNPedestalEvents++;
+
+  return(ihit);
+}
+//_____________________________________________________________________________
+void THcShowerArray::CalculatePedestals( )
+{
+  // Doesn't do anything yet
+
+}
+//_____________________________________________________________________________
+void THcShowerArray::InitializePedestals( )
+{
+  // Doesn't do anything yet
+  fNPedestalEvents = 0;
+} 

src/THcShowerArray.h

+#ifndef ROOT_THcShowerArray
+#define ROOT_THcShowerArray 
+
+//////////////////////////////////////////////////////////////////////////////
+//                         
+// THcShowerArray
+//
+// A Hall C Fly's Eye Shower Array
+//
+// Subdetector for the fly's eye part of the SHMS shower counter.
+// 
+//////////////////////////////////////////////////////////////////////////////
+
+#include "THaSubDetector.h"
+#include "TClonesArray.h"
+
+#include <iostream>
+
+#include <fstream>
+
+class THaEvData;
+class THaSignalHit;
+
+class THcShowerArray : public THaSubDetector {
+
+public:
+  THcShowerArray( const char* name, const char* description,
+		  Int_t planenum, THaDetectorBase* parent = NULL);
+  virtual ~THcShowerArray();
+
+  virtual void    Clear( Option_t* opt="" );
+  virtual Int_t Decode( const THaEvData& );
+  virtual EStatus Init( const TDatime& run_time );
+
+  virtual Int_t CoarseProcess( TClonesArray& tracks );
+  virtual Int_t FineProcess( TClonesArray& tracks );
+  Bool_t   IsTracking() { return kFALSE; }
+  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( );
+  virtual void  InitializePedestals( );
+
+  //  Double_t fSpacing;   not used
+
+  TClonesArray* fParentHitList;
+
+  // Return zero for now
+  Double_t GetEplane() {
+    return 0.0;
+  };
+
+
+protected:
+
+  Double_t* fA;                 // [fNelem] ADC amplitude of blocks
+
+  TClonesArray* fADCHits;	// List of ADC hits
+
+  Int_t fNPedestalEvents;	/* Pedestal event counter */
+  Int_t fMinPeds;		/* Only analyze/update if num events > */
+
+  // Parameters
+  Int_t fNRows;
+  Int_t fNColumns;
+
+  Int_t fLayerNum;		// 2 for SHMS
+  // Accumulators for pedestals go here
+  // 2D arrays
+
+  virtual Int_t  ReadDatabase( const TDatime& date );
+  virtual Int_t  DefineVariables( EMode mode = kDefine );
+  //virtual void  InitializePedestals( );
+  ClassDef(THcShowerArray,0); // Fly;s Eye calorimeter array
+};
+#endif