Added hits plotting

Makefile

@@ -14,7 +14,7 @@ SRC  =  src/THcInterface.cxx src/THcParmList.cxx src/THcAnalyzer.cxx \
 	src/THcHitList.cxx src/THcDetectorMap.cxx src/THcHodoscope.cxx \
 	src/THcHallCSpectrometer.cxx src/THcDriftChamber.cxx \
 	src/THcScintillatorPlane.cxx src/THcSignalHit.cxx \
-	src/THcShower.cxx src/THcShowerHit.cxx 
+	src/THcShower.cxx src/THcShowerHit.cxx src/THcShowerPlane.cxx  
 
 # Name of your package. 
 # The shared library that will be built will get the name lib$(PACKAGE).so

examples/calhit.C

+{
+  TFile* f = new TFile("hodtest.root");
+ 
+  TCanvas *c1 = new TCanvas("c1", "Shower Hit Maps", 800, 800); 
+  c1->Divide(1, 4);
+
+  TH1F* h[6];
+
+  h[0] = chposadc1;
+  h[1] = chnegadc1;
+  h[2] = chposadc2;
+  h[3] = chnegadc2;
+  h[4] = chposadc3;
+  h[5] = chposadc4;
+
+  
+    c1->cd(1);
+    h[0]->SetFillColor(kGreen);
+    h[0]->SetFillStyle(3345);
+    h[0]->Draw();
+    h[1]->SetFillColor(kBlue);
+    h[1]->SetFillStyle(3354);
+    h[1]->Draw("same");
+  
+    c1->cd(2);
+    h[2]->SetFillColor(kGreen);
+    h[2]->SetFillStyle(3345);
+    h[2]->Draw();
+    h[3]->SetFillColor(kBlue);
+    h[3]->SetFillStyle(3354);
+    h[3]->Draw("same");
+
+    c1->cd(3);
+    h[4]->SetFillColor(kGreen);
+    h[4]->SetFillStyle(3345);
+    h[4]->Draw();
+
+    c1->cd(4);
+    h[5]->SetFillColor(kGreen);
+    h[5]->SetFillStyle(3345);
+    h[5]->Draw();
+}
+ 

examples/output.def

@@ -2,6 +2,7 @@
 #
 block H.hod.*
 
+
 # TDC hits per paddle
 TH1F hpostdc1 'HMS s1x+ TDC hits' H.hod.1x.postdchits 16 0.5 16.5
 TH1F hnegtdc1 'HMS s1x- TDC hits' H.hod.1x.negtdchits 16 0.5 16.5
@@ -21,3 +22,12 @@ TH1F hposadc3 'HMS s2x+ ADC hits' H.hod.2x.posadchits 16 0.5 16.5
 TH1F hnegadc3 'HMS s2x- ADC hits' H.hod.2x.negadchits 16 0.5 16.5
 TH1F hposadc4 'HMS s2y+ ADC hits' H.hod.2y.posadchits 10 0.5 10.5
 TH1F hnegadc4 'HMS s2y- ADC hits' H.hod.2y.negadchits 10 0.5 10.5
+
+block H.Cal.*
+# ADC hits per Calorimeter layer
+TH1F chposadc1 'HMS Cal 1z+ ADC hits' H.Cal.1z.posadchits 13 0.5 13.5
+TH1F chnegadc1 'HMS Cal 1z- ADC hits' H.Cal.1z.negadchits 13 0.5 13.5
+TH1F chposadc2 'HMS Cal 2z+ ADC hits' H.Cal.2z.posadchits 13 0.5 13.5
+TH1F chnegadc2 'HMS Cal 2z- ADC hits' H.Cal.2z.negadchits 13 0.5 13.5
+TH1F chposadc3 'HMS Cal 3z+ ADC hits' H.Cal.3z.posadchits 13 0.5 13.5
+TH1F chposadc4 'HMS Cal 4z+ ADC hits' H.Cal.4z.posadchits 13 0.5 13.5

src/HallC_LinkDef.h

@@ -19,6 +19,7 @@
 #pragma link C++ class THcDetectorMap+;
 #pragma link C++ class THcHallCSpectrometer+;
 #pragma link C++ class THcScintillatorPlane+;
+#pragma link C++ class THcShowerPlane+;
 #pragma link C++ class THcSignalHit+;
 #pragma link C++ class THcShowerHit+;
 #pragma link C++ class THcShower+;

src/THcShower.cxx

@@ -36,7 +36,7 @@ THcShower::THcShower( const char* name, const char* description,
   THaNonTrackingDetector(name,description,apparatus)
 {
   // Constructor
-
+  Setup(name, description);
 //  fTrackProj = new TClonesArray( "THaTrackProj", 5 );
 }
 
@@ -47,11 +47,51 @@ THcShower::THcShower( ) :
   // Constructor
 }
 
+void THcShower::Setup(const char* name, const char* description)
+{
+
+  static const char* const here = "Setup()";
+  static const char* const message = 
+    "Must construct %s detector with valid name! Object construction failed.";
+
+  // Base class constructor failed?
+  if( IsZombie()) return;
+
+  fNLayers = 4;	// Eventually get # layers and layer names from a DB
+  fLayerNames = new char* [fNLayers];
+  for(Int_t i=0;i<fNLayers;i++) {fLayerNames[i] = new char[3];}
+  strcpy(fLayerNames[0],"1z");
+  strcpy(fLayerNames[1],"2z");
+  strcpy(fLayerNames[2],"3z");
+  strcpy(fLayerNames[3],"4z");
+
+  size_t nlen = strlen(name);
+  size_t slen = 0;
+  for(Int_t i=0;i < fNLayers;i++)
+    {slen = TMath::Max(slen,strlen(fLayerNames[i]));}
+  size_t len = nlen+slen+1;
+
+  // Probably shouldn't assume that description is defined
+  char* desc = new char[strlen(description)+50+slen];
+  fPlanes = new THcShowerPlane* [fNLayers];
+  for(Int_t i=0;i < fNLayers;i++) {
+    strcpy(desc, description);
+    strcat(desc, " Plane ");
+    strcat(desc, fLayerNames[i]);
+
+    fPlanes[i] = new THcShowerPlane(fLayerNames[i], desc, i+1, this); 
+    cout << "Created Shower Plane " << fLayerNames[i] << ", " << desc << endl;
+  }
+}
+
+
 //_____________________________________________________________________________
 THaAnalysisObject::EStatus THcShower::Init( const TDatime& date )
 {
   static const char* const here = "Init()";
 
+  cout << "THcHodoscope::Init " << GetName() << endl;
+
   if( THaNonTrackingDetector::Init( date ) )
     return fStatus;
 
@@ -61,6 +101,15 @@ THaAnalysisObject::EStatus THcShower::Init( const TDatime& date )
 
   THcHitList::InitHitList(fDetMap, "THcShowerHit", 100);
 
+  EStatus status;
+  if( (status = THaNonTrackingDetector::Init( date )) )
+    return fStatus=status;
+
+  for(Int_t ip=0;ip<fNLayers;ip++) {
+    if((status = fPlanes[ip]->Init( date ))) {
+      return fStatus=status;
+    }
+  }
   // Will need to determine which apparatus it belongs to and use the
   // appropriate detector ID in the FillMap call
   if( gHcDetectorMap->FillMap(fDetMap, "HCAL") < 0 ) {
@@ -86,24 +135,6 @@ Int_t THcShower::ReadDatabase( const TDatime& date )
   // Pull values from the THcParmList instead of reading a database
   // file like Hall A does.
 
-  //  DBRequest list[] = {
-  //    { "TDC_offsetsL", fLOff, kDouble, fNelem },
-  //    { "TDC_offsetsR", fROff, kDouble, fNelem },
-  //    { "ADC_pedsL", fLPed, kDouble, fNelem },
-  //    { "ADC_pedsR", fRPed, kDouble, fNelem },
-  //    { "ADC_coefL", fLGain, kDouble, fNelem },
-  //    { "ADC_coefR", fRGain, kDouble, fNelem },
-  //    { "TDC_res",   &fTdc2T },
-  //    { "TransSpd",  &fCn },
-  //    { "AdcMIP",    &fAdcMIP },
-  //    { "NTWalk",    &fNTWalkPar, kInt },
-  //    { "Timewalk",  fTWalkPar, kDouble, 2*fNelem },
-  //    { "ReTimeOff", fTrigOff, kDouble, fNelem },
-  //    { "AvgRes",    &fResolution },
-  //    { "Atten",     &fAttenuation },
-  //    { 0 }
-  //  };
-
   // We will probably want to add some kind of method to gHcParms to allow
   // bulk retrieval of parameters of interest.
 
@@ -112,9 +143,9 @@ Int_t THcShower::ReadDatabase( const TDatime& date )
 
 cout << "THcShower::ReadDatabase called " << GetName() << endl;
 
-  NLayers = 4;			// Hardwire for now
+  fNLayers = 4;			// Hardwire for now
 
-  BlockThick = new Double_t [NLayers];
+  BlockThick = new Double_t [fNLayers];
 
   BlockThick[0] = *(Double_t *)gHcParms->Find("hcal_1pr_thick")->GetValuePointer();
   BlockThick[1] = *(Double_t *)gHcParms->Find("hcal_2ta_thick")->GetValuePointer();
@@ -123,7 +154,7 @@ cout << "THcShower::ReadDatabase called " << GetName() << endl;
 
 cout << "Block thickness: " << BlockThick[2]  << endl;
 
-fNBlocks = new Int_t [NLayers];
+fNBlocks = new Int_t [fNLayers];
 
   fNBlocks[0] = *(Int_t *)gHcParms->Find("hcal_1pr_nr")->GetValuePointer();
   fNBlocks[1] = *(Int_t *)gHcParms->Find("hcal_2ta_nr")->GetValuePointer();
@@ -132,7 +163,7 @@ fNBlocks = new Int_t [NLayers];
 
   cout << "Number of blocks per layer: " << fNBlocks[2]  << endl;
 
-  fNLayerZPos = new Double_t [NLayers];
+  fNLayerZPos = new Double_t [fNLayers];
 
   fNLayerZPos[0] = *(Double_t *)gHcParms->Find("hcal_1pr_zpos")->GetValuePointer();
   fNLayerZPos[1] = *(Double_t *)gHcParms->Find("hcal_2ta_zpos")->GetValuePointer();
@@ -141,7 +172,7 @@ fNBlocks = new Int_t [NLayers];
 
 cout << "Z Position: " << fNLayerZPos[2]  << endl;
 
-  XPos = new Double_t [2*NLayers];
+  XPos = new Double_t [2*fNLayers];
 
   XPos[0] = *(Double_t *)gHcParms->Find("hcal_1pr_left")->GetValuePointer();
   XPos[1] = *(Double_t *)gHcParms->Find("hcal_1pr_right")->GetValuePointer();
@@ -210,6 +241,8 @@ Int_t THcShower::DefineVariables( EMode mode )
   if( mode == kDefine && fIsSetup ) return kOK;
   fIsSetup = ( mode == kDefine );
 
+  cout << "THcShower::DefineVariables called " << GetName() << endl;
+
   // Register variables in global list
 
  //  RVarDef vars[] = {
@@ -278,7 +311,11 @@ Int_t THcShower::Decode( const THaEvData& evdata )
 
   // Get the Hall C style hitlist (fRawHitList) for this event
   Int_t nhits = THcHitList::DecodeToHitList(evdata);
-
+  Int_t nexthit = 0;
+  for(Int_t ip=0;ip<fNLayers;ip++) {
+    nexthit = fPlanes[ip]->ProcessHits(fRawHitList, nexthit);
+  }
+/*
 //   fRawHitList is TClones array of THcShowerHit objects
   for(Int_t ihit = 0; ihit < fNRawHits ; ihit++) {
     THcShowerHit* hit = (THcShowerHit *) fRawHitList->At(ihit);
@@ -286,7 +323,7 @@ Int_t THcShower::Decode( const THaEvData& evdata )
 	 << hit->fADC_pos << " " << hit->fADC_neg << " "  << endl;
   }
   cout << endl;
-
+*/
   return nhits;
 }
 

src/THcShower.h

@@ -3,7 +3,7 @@
 
 ///////////////////////////////////////////////////////////////////////////////
 //                                                                           //
-// THcHodoscope                                                              //
+// THcShower                                                                 //
 //                                                                           //
 ///////////////////////////////////////////////////////////////////////////////
 
@@ -11,7 +11,7 @@
 #include "THaNonTrackingDetector.h"
 #include "THcHitList.h"
 #include "THcShowerHit.h"
-#include "THcScintillatorPlane.h"
+#include "THcShowerPlane.h"
 
 class THaScCalib;
 
@@ -45,15 +45,15 @@ protected:
 
 
   // Potential Hall C parameters.  Mostly here for demonstration
-  char** LayerNames;
-  Int_t NLayers;
+  char** fLayerNames;
+  Int_t fNLayers;
   Double_t* fNLayerZPos;		// Z position of front of shower counter layers
   Double_t* BlockThick;		// Thickness of shower counter blocks, blocks
   Int_t* fNBlocks;           // Number of shower counter blocks per layer
   Double_t** YPos;		//X,Y positions of shower counter blocks
   Double_t* XPos;
 
-  THcScintillatorPlane** fPlane; // List of plane objects
+  THcShowerPlane** fPlanes; // List of plane objects
 
   TClonesArray*  fTrackProj;  // projection of track onto scintillator plane
                               // and estimated match to TOF paddle
@@ -84,6 +84,8 @@ protected:
   virtual  Double_t TimeWalkCorrection(const Int_t& paddle,
 					   const ESide side);
 
+void Setup(const char* name, const char* description);
+
   ClassDef(THcShower,0)   // Generic hodoscope class
 };
 

src/THcShowerPlane.cxx

+//*-- Author :
+
+//////////////////////////////////////////////////////////////////////////
+//
+// THcShowerPlane
+//
+//////////////////////////////////////////////////////////////////////////
+
+#include "THcShowerPlane.h"
+#include "TClonesArray.h"
+#include "THcSignalHit.h"
+#include "THcGlobals.h"
+#include "THcParmList.h"
+#include "THcHitList.h"
+#include "THcShower.h"
+#include "TClass.h"
+
+#include <cstring>
+#include <cstdio>
+#include <cstdlib>
+#include <iostream>
+
+using namespace std;
+
+ClassImp(THcShowerPlane)
+
+//______________________________________________________________________________
+THcShowerPlane::THcShowerPlane( const char* name, 
+					    const char* description,
+					    const Int_t layernum,
+					    THaDetectorBase* parent )
+  : THaSubDetector(name,description,parent)
+{
+  // Normal constructor with name and description
+  fPosADCHits = new TClonesArray("THcSignalHit",16);
+  fNegADCHits = new TClonesArray("THcSignalHit",16);
+  fPosADCHitsClass = fPosADCHits->GetClass();
+  fNegADCHitsClass = fNegADCHits->GetClass();
+  fLayerNum = layernum;
+}
+
+//______________________________________________________________________________
+THcShowerPlane::~THcShowerPlane()
+{
+  // Destructor
+  delete fPosADCHits;
+  delete fNegADCHits;
+
+}
+THaAnalysisObject::EStatus THcShowerPlane::Init( const TDatime& date )
+{
+  // Extra initialization for shower layer: set up DataDest map
+
+  cout << "THcShowerPlane::Init called " << GetName() << endl;
+
+  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 THcShowerPlane::ReadDatabase( const TDatime& date )
+{
+
+  // See what file it looks for
+  
+  static const char* const here = "ReadDatabase()";
+  char prefix[2];
+  char parname[100];
+  
+  prefix[0]=tolower(GetParent()->GetPrefix()[0]);
+  prefix[1]='\0';
+
+  strcpy(parname,prefix);
+  strcat(parname,"cal_");
+  strcat(parname,GetName());
+  Int_t plen=strlen(parname);
+
+  strcat(parname,"_nr");
+  cout << " Getting value of SHOWER!!!" << parname << endl;
+//   fNelem = *(Int_t *)gHcParms->Find(parname)->GetValuePointer();
+// 
+//   parname[plen]='\0';
+//   strcat(parname,"_spacing");
+// 
+//   fSpacing =  gHcParms->Find(parname)->GetValue(0);
+  
+  // First letter of GetParent()->GetPrefix() tells us what prefix to
+  // use on parameter names.  
+
+
+  //  Find the number of elements
+  
+  // Create arrays to hold results here
+
+
+  return kOK;
+}
+//_____________________________________________________________________________
+Int_t THcShowerPlane::DefineVariables( EMode mode )
+{
+  // Initialize global variables and lookup table for decoder
+
+  cout << "THcShowerPlane::DefineVariables called " << GetName() << endl;
+
+  if( mode == kDefine && fIsSetup ) return kOK;
+  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()"},
+    { 0 }
+  };
+
+  return DefineVarsFromList( vars, mode );
+}
+
+//_____________________________________________________________________________
+void THcShowerPlane::Clear( Option_t* )
+{
+  //cout << " Calling THcShowerPlane::Clear " << GetName() << endl;
+  // Clears the hit lists
+  fPosADCHits->Clear();
+  fNegADCHits->Clear();
+}
+
+//_____________________________________________________________________________
+Int_t THcShowerPlane::Decode( const THaEvData& evdata )
+{
+  // Doesn't actually get called.  Use Fill method instead
+  cout << " Calling THcShowerPlane::Decode " << GetName() << endl;
+
+  return 0;
+}
+//_____________________________________________________________________________
+Int_t THcShowerPlane::CoarseProcess( TClonesArray& tracks )
+{
+ 
+  //  HitCount();
+
+ return 0;
+}
+
+//_____________________________________________________________________________
+Int_t THcShowerPlane::FineProcess( TClonesArray& tracks )
+{
+  return 0;
+}
+Int_t THcShowerPlane::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.
+
+  Int_t nPosADCHits=0;
+  Int_t nNegADCHits=0;
+  fPosADCHits->Clear();
+  fNegADCHits->Clear();
+
+  Int_t nrawhits = rawhits->GetLast()+1;
+
+  Int_t ihit = nexthit;
+  while(ihit < nrawhits) {
+    THcShowerHit* hit = (THcShowerHit *) rawhits->At(ihit);
+    if(hit->fPlane > fLayerNum) {
+      break;
+    }
+
+
+if(hit->fADC_pos >  0) {
+#if ROOT_VERSION_CODE >= ROOT_VERSION(5,32,0)
+	THcSignalHit *sighit = (THcSignalHit*) fPosADCHits->ConstructedAt(nPosADCHits++);
+	sighit->Set(hit->fCounter, hit->fADC_pos);
+#else
+	TObject* obj = (*fPosTDCHits)[nPosADCHits++];
+	R__ASSERT( obj );
+if(!obj->TestBit (TObject::kNotDeleted))
+	fPosADCHitsClass->New(obj);
+	THcSignalHit *sighit = (THcSignalHit*)obj;
+#endif
+   sighit->Set(hit->fCounter, hit->fADC_pos);
+}
+
+if(hit->fADC_neg >  0) {
+#if ROOT_VERSION_CODE >= ROOT_VERSION(5,32,0)
+	THcSignalHit *sighit = (THcSignalHit*) fNegADCHits->ConstructedAt(nNegADCHits++);
+	sighit->Set(hit->fCounter, hit->fADC_neg);
+#else
+	TObject* obj = (*fPosTDCHits)[nNegADCHits++];
+	R__ASSERT( obj );
+if(!obj->TestBit (TObject::kNotDeleted))
+	fNegADCHitsClass->New(obj);
+	THcSignalHit *sighit = (THcSignalHit*)obj;
+#endif
+   sighit->Set(hit->fCounter, hit->fADC_neg);
+}
+
+    ihit++;
+  }
+  return(ihit);
+}
+
+    
+  
+  
+

src/THcShowerPlane.h

+#ifndef ROOT_THcShowerPlane
+#define ROOT_THcShowerPlane
+
+//////////////////////////////////////////////////////////////////////////////
+//                         
+// THcShowerPlane
+//
+// A Hall C scintillator plane
+//
+// May want to later inherit from a THcPlane class if there are similarities
+// in what a plane is shared with other detector types (shower, etc.)
+// 
+//////////////////////////////////////////////////////////////////////////////
+
+#include "THaSubDetector.h"
+#include "TClonesArray.h"
+
+class THaEvData;
+class THaSignalHit;
+
+class THcShowerPlane : public THaSubDetector {
+  
+ public:
+  THcShowerPlane( const char* name, const char* description,
+			Int_t planenum, THaDetectorBase* parent = NULL);
+  virtual ~THcShowerPlane();
+
+  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);
+
+  Double_t fSpacing;
+
+  TClonesArray* fParentHitList;
+
+ protected:
+
+  TClonesArray* fPosADCHits;
+  TClonesArray* fNegADCHits;
+
+  TClass* fPosADCHitsClass;
+  TClass* fNegADCHitsClass;
+
+  Int_t fLayerNum;
+
+  virtual Int_t  ReadDatabase( const TDatime& date );
+  virtual Int_t  DefineVariables( EMode mode = kDefine );
+
+  ClassDef(THcShowerPlane,0)
+};
+#endif
+
+
+