diff --git a/Makefile b/Makefile index 583996d3c4e14e461d1696e19877bb6e06e36908..14bc3941dfd27ea0975c78aa0ad6efdadba629e0 100644 --- a/Makefile +++ b/Makefile @@ -13,7 +13,8 @@ SRC = src/THcInterface.cxx src/THcParmList.cxx src/THcAnalyzer.cxx \ src/THcDCHit.cxx \ src/THcHitList.cxx src/THcDetectorMap.cxx src/THcHodoscope.cxx \ src/THcHallCSpectrometer.cxx src/THcDriftChamber.cxx \ - src/THcScintillatorPlane.cxx src/THcSignalHit.cxx + src/THcScintillatorPlane.cxx src/THcSignalHit.cxx \ + src/THcShower.cxx src/THcShowerHit.cxx # Name of your package. # The shared library that will be built will get the name lib$(PACKAGE).so diff --git a/examples/hodtest.C b/examples/hodtest.C index f655a8b15529bbbb72c9a27eec0c59613403b383..f730a8c5a076e1830498c48a4718cdde398292ea 100644 --- a/examples/hodtest.C +++ b/examples/hodtest.C @@ -17,6 +17,7 @@ // Add hodoscope HMS->AddDetector( new THcHodoscope("hod", "Hodoscope" )); + HMS->AddDetector( new THcShower("Cal", "Shower" )); // HMS->AddDetector( new THcDriftChamber("dc", "Drift Chambers" )); // Set up the analyzer - we use the standard one, diff --git a/src/HallC_LinkDef.h b/src/HallC_LinkDef.h index a8d80d8c6b54d291706a3a61a237af14e73c1e96..f904bb6cf104b8f5225f23974e0fb687b0fe3c7f 100644 --- a/src/HallC_LinkDef.h +++ b/src/HallC_LinkDef.h @@ -20,5 +20,7 @@ #pragma link C++ class THcHallCSpectrometer+; #pragma link C++ class THcScintillatorPlane+; #pragma link C++ class THcSignalHit+; +#pragma link C++ class THcShowerHit+; +#pragma link C++ class THcShower+; #endif diff --git a/src/THcShower.cxx b/src/THcShower.cxx new file mode 100644 index 0000000000000000000000000000000000000000..7d89c07841790e3623b53ffc03c46a10eed1a000 --- /dev/null +++ b/src/THcShower.cxx @@ -0,0 +1,339 @@ +/////////////////////////////////////////////////////////////////////////////// +// // +// THcShower // +// // +// Shower counter class, describing a generic segmented shower detector. // +// // +// // +// // +// // +/////////////////////////////////////////////////////////////////////////////// + +#include "THcShower.h" +#include "THaEvData.h" +#include "THaDetMap.h" +#include "THcDetectorMap.h" +#include "THcGlobals.h" +#include "THcParmList.h" +#include "VarDef.h" +#include "VarType.h" +#include "THaTrack.h" +#include "TClonesArray.h" +#include "TMath.h" + +#include "THaTrackProj.h" + +#include <cstring> +#include <cstdio> +#include <cstdlib> +#include <iostream> + +using namespace std; + +//_____________________________________________________________________________ +THcShower::THcShower( const char* name, const char* description, + THaApparatus* apparatus ) : + THaNonTrackingDetector(name,description,apparatus) +{ + // Constructor + +// fTrackProj = new TClonesArray( "THaTrackProj", 5 ); +} + +//_____________________________________________________________________________ +THcShower::THcShower( ) : + THaNonTrackingDetector() +{ + // Constructor +} + +//_____________________________________________________________________________ +THaAnalysisObject::EStatus THcShower::Init( const TDatime& date ) +{ + static const char* const here = "Init()"; + + if( THaNonTrackingDetector::Init( date ) ) + return fStatus; + + + // Should probably put this in ReadDatabase as we will know the + // maximum number of hits after setting up the detector map + + THcHitList::InitHitList(fDetMap, "THcShowerHit", 100); + + // 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 ) { + Error( Here(here), "Error filling detectormap for %s.", + "HCAL"); + return kInitError; + } + + return fStatus = kOK; +} + +//_____________________________________________________________________________ +Int_t THcShower::ReadDatabase( const TDatime& date ) +{ + // Read this detector's parameters from the database file 'fi'. + // This function is called by THaDetectorBase::Init() once at the + // beginning of the analysis. + // 'date' contains the date/time of the run being analyzed. + + // static const char* const here = "ReadDatabase()"; + + // Read data from database + // 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. + + // Will need to determine which spectrometer in order to construct + // the parameter names (e.g. hscin_1x_nr vs. sscin_1x_nr) + +cout << "THcShower::ReadDatabase called " << GetName() << endl; + + NLayers = 4; // Hardwire for now + + BlockThick = new Double_t [NLayers]; + + BlockThick[0] = *(Double_t *)gHcParms->Find("hcal_1pr_thick")->GetValuePointer(); + BlockThick[1] = *(Double_t *)gHcParms->Find("hcal_2ta_thick")->GetValuePointer(); + BlockThick[2] = *(Double_t *)gHcParms->Find("hcal_3ta_thick")->GetValuePointer(); + BlockThick[3] = *(Double_t *)gHcParms->Find("hcal_4ta_thick")->GetValuePointer(); + +cout << "Block thickness: " << BlockThick[2] << endl; + +fNBlocks = new Int_t [NLayers]; + + fNBlocks[0] = *(Int_t *)gHcParms->Find("hcal_1pr_nr")->GetValuePointer(); + fNBlocks[1] = *(Int_t *)gHcParms->Find("hcal_2ta_nr")->GetValuePointer(); + fNBlocks[2] = *(Int_t *)gHcParms->Find("hcal_3ta_nr")->GetValuePointer(); + fNBlocks[3] = *(Int_t *)gHcParms->Find("hcal_4ta_nr")->GetValuePointer(); + + cout << "Number of blocks per layer: " << fNBlocks[2] << endl; + + fNLayerZPos = new Double_t [NLayers]; + + fNLayerZPos[0] = *(Double_t *)gHcParms->Find("hcal_1pr_zpos")->GetValuePointer(); + fNLayerZPos[1] = *(Double_t *)gHcParms->Find("hcal_2ta_zpos")->GetValuePointer(); + fNLayerZPos[2] = *(Double_t *)gHcParms->Find("hcal_3ta_zpos")->GetValuePointer(); + fNLayerZPos[3] = *(Double_t *)gHcParms->Find("hcal_4ta_zpos")->GetValuePointer(); + +cout << "Z Position: " << fNLayerZPos[2] << endl; + + XPos = new Double_t [2*NLayers]; + + XPos[0] = *(Double_t *)gHcParms->Find("hcal_1pr_left")->GetValuePointer(); + XPos[1] = *(Double_t *)gHcParms->Find("hcal_1pr_right")->GetValuePointer(); + XPos[2] = *(Double_t *)gHcParms->Find("hcal_2ta_left")->GetValuePointer(); + XPos[3] = *(Double_t *)gHcParms->Find("hcal_2ta_right")->GetValuePointer(); + XPos[4] = *(Double_t *)gHcParms->Find("hcal_3ta_left")->GetValuePointer(); + XPos[5] = *(Double_t *)gHcParms->Find("hcal_3ta_right")->GetValuePointer(); + XPos[6] = *(Double_t *)gHcParms->Find("hcal_4ta_left")->GetValuePointer(); + XPos[7] = *(Double_t *)gHcParms->Find("hcal_4ta_right")->GetValuePointer(); + +cout << "X Positions: " << XPos[0] << ", " << XPos[1] << endl; + + YPos = new Double_t* [4]; + cout << "Y Positions:"; + + Double_t* p; + Int_t ilayer; + + ilayer = 0; + p = (Double_t *)gHcParms->Find("hcal_1pr_top")->GetValuePointer(); + YPos[ilayer] = new Double_t [fNBlocks[ilayer]]; + // Print out some parameters just to demonstrate that it works + + for(Int_t i=0;i<fNBlocks[ilayer];i++) { + YPos[ilayer][i] = p[i]; + cout << " " << YPos[ilayer][i]; + } + cout << endl; + + ilayer = 1; + p = (Double_t *)gHcParms->Find("hcal_2ta_top")->GetValuePointer(); + YPos[ilayer] = new Double_t [fNBlocks[ilayer]]; + + for(Int_t i=0;i<fNBlocks[ilayer];i++) { + YPos[ilayer][i] = p[i]; + } + + ilayer = 2; + p = (Double_t *)gHcParms->Find("hcal_3ta_top")->GetValuePointer(); + YPos[ilayer] = new Double_t [fNBlocks[ilayer]]; + + for(Int_t i=0;i<fNBlocks[ilayer];i++) { + YPos[ilayer][i] = p[i]; + } + + ilayer = 3; + p = (Double_t *)gHcParms->Find("hcal_4ta_top")->GetValuePointer(); + YPos[ilayer] = new Double_t [fNBlocks[ilayer]]; + + for(Int_t i=0;i<fNBlocks[ilayer];i++) { + YPos[ilayer][i] = p[i]; + } + + fIsInit = true; + + return kOK; +} + + + +//_____________________________________________________________________________ +Int_t THcShower::DefineVariables( EMode mode ) +{ + // Initialize global variables and lookup table for decoder + + if( mode == kDefine && fIsSetup ) return kOK; + fIsSetup = ( mode == kDefine ); + + // Register variables in global list + + // RVarDef vars[] = { + // { "nhit", "Number of hits", "fNhits" }, + // { "a", "Raw ADC amplitudes", "fA" }, + // { "a_p", "Ped-subtracted ADC amplitudes", "fA_p" }, + // { "a_c", "Calibrated ADC amplitudes", "fA_c" }, + // { "asum_p", "Sum of ped-subtracted ADCs", "fAsum_p" }, + // { "asum_c", "Sum of calibrated ADCs", "fAsum_c" }, + // { "nclust", "Number of clusters", "fNclust" }, + // { "e", "Energy (MeV) of largest cluster", "fE" }, + // { "x", "x-position (cm) of largest cluster", "fX" }, + // { "y", "y-position (cm) of largest cluster", "fY" }, + // { "mult", "Multiplicity of largest cluster", "fMult" }, + // { "nblk", "Numbers of blocks in main cluster", "fNblk" }, + // { "eblk", "Energies of blocks in main cluster", "fEblk" }, + // { "trx", "track x-position in det plane", "fTRX" }, + // { "try", "track y-position in det plane", "fTRY" }, + // { 0 } + // }; + // return DefineVarsFromList( vars, mode ); + return kOK; +} + +//_____________________________________________________________________________ +THcShower::~THcShower() +{ + // Destructor. Remove variables from global list. + + if( fIsSetup ) + RemoveVariables(); + if( fIsInit ) + DeleteArrays(); + if (fTrackProj) { + fTrackProj->Clear(); + delete fTrackProj; fTrackProj = 0; + } +} + +//_____________________________________________________________________________ +void THcShower::DeleteArrays() +{ + // Delete member arrays. Used by destructor. + + delete [] BlockThick; BlockThick = NULL; + delete [] fNBlocks; fNBlocks = NULL; + delete [] fNLayerZPos; fNLayerZPos = NULL; + delete [] XPos; XPos = NULL; + delete [] YPos; YPos = NULL; + //delete [] fSpacing; fSpacing = NULL; + //delete [] fCenter; fCenter = NULL; // This 2D. What is correct way to delete? +} + +//_____________________________________________________________________________ +inline +void THcShower::ClearEvent() +{ + // Reset per-event data. + + fTrackProj->Clear(); +} + +//_____________________________________________________________________________ +Int_t THcShower::Decode( const THaEvData& evdata ) +{ + + // Get the Hall C style hitlist (fRawHitList) for this event + Int_t nhits = THcHitList::DecodeToHitList(evdata); + +// fRawHitList is TClones array of THcShowerHit objects + for(Int_t ihit = 0; ihit < fNRawHits ; ihit++) { + THcShowerHit* hit = (THcShowerHit *) fRawHitList->At(ihit); + cout << ihit << " : " << hit->fPlane << ":" << hit->fCounter << " : " + << hit->fADC_pos << " " << hit->fADC_neg << " " << endl; + } + cout << endl; + + return nhits; +} + +//_____________________________________________________________________________ +Int_t THcShower::ApplyCorrections( void ) +{ + return(0); +} + +//_____________________________________________________________________________ +Double_t THcShower::TimeWalkCorrection(const Int_t& paddle, + const ESide side) +{ + return(0.0); +} + +//_____________________________________________________________________________ +Int_t THcShower::CoarseProcess( TClonesArray& ) //tracks +{ + // Calculation of coordinates of particle track cross point with shower + // plane in the detector coordinate system. For this, parameters of track + // reconstructed in THaVDC::CoarseTrack() are used. + // + // Apply corrections and reconstruct the complete hits. + // + // static const Double_t sqrt2 = TMath::Sqrt(2.); + + ApplyCorrections(); + + return 0; +} + +//_____________________________________________________________________________ +Int_t THcShower::FineProcess( TClonesArray& tracks ) +{ + // Reconstruct coordinates of particle track cross point with shower + // plane, and copy the data into the following local data structure: + // + // Units of measurements are meters. + + // Calculation of coordinates of particle track cross point with shower + // plane in the detector coordinate system. For this, parameters of track + // reconstructed in THaVDC::FineTrack() are used. + + return 0; +} + +ClassImp(THcShower) +//////////////////////////////////////////////////////////////////////////////// + diff --git a/src/THcShower.h b/src/THcShower.h new file mode 100644 index 0000000000000000000000000000000000000000..b29e44e0523ab66c5d8ec8f91e376c35f174294c --- /dev/null +++ b/src/THcShower.h @@ -0,0 +1,93 @@ +#ifndef ROOT_THcShower +#define ROOT_THcShower + +/////////////////////////////////////////////////////////////////////////////// +// // +// THcHodoscope // +// // +/////////////////////////////////////////////////////////////////////////////// + +#include "TClonesArray.h" +#include "THaNonTrackingDetector.h" +#include "THcHitList.h" +#include "THcShowerHit.h" +#include "THcScintillatorPlane.h" + +class THaScCalib; + +class THcShower : public THaNonTrackingDetector, public THcHitList { + +public: + THcShower( const char* name, const char* description = "", + THaApparatus* a = NULL ); + virtual ~THcShower(); + + 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 ); + + virtual Int_t ApplyCorrections( void ); + + // Int_t GetNHits() const { return fNhit; } + + Int_t GetNTracks() const { return fTrackProj->GetLast()+1; } + const TClonesArray* GetTrackHits() const { return fTrackProj; } + + friend class THaScCalib; + + THcShower(); // for ROOT I/O +protected: + + // Calibration + + // Per-event data + + + // Potential Hall C parameters. Mostly here for demonstration + char** LayerNames; + Int_t NLayers; + 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 + + TClonesArray* fTrackProj; // projection of track onto scintillator plane + // and estimated match to TOF paddle + // Useful derived quantities + // double tan_angle, sin_angle, cos_angle; + + // static const char NDEST = 2; + // struct DataDest { + // Int_t* nthit; + // Int_t* nahit; + // Double_t* tdc; + // Double_t* tdc_c; + // Double_t* adc; + // Double_t* adc_p; + // Double_t* adc_c; + // Double_t* offset; + // Double_t* ped; + // Double_t* gain; + // } fDataDest[NDEST]; // Lookup table for decoder + + void ClearEvent(); + void DeleteArrays(); + virtual Int_t ReadDatabase( const TDatime& date ); + virtual Int_t DefineVariables( EMode mode = kDefine ); + + enum ESide { kLeft = 0, kRight = 1 }; + + virtual Double_t TimeWalkCorrection(const Int_t& paddle, + const ESide side); + + ClassDef(THcShower,0) // Generic hodoscope class +}; + +//////////////////////////////////////////////////////////////////////////////// + +#endif + diff --git a/src/THcShowerHit.cxx b/src/THcShowerHit.cxx new file mode 100644 index 0000000000000000000000000000000000000000..c27357235534b2b8c67d13cdbcba0984c6fb8b09 --- /dev/null +++ b/src/THcShowerHit.cxx @@ -0,0 +1,75 @@ +/////////////////////////////////////////////////////////////////////////////// +// // +// THcShowerHit // +// // +// Class representing a single raw hit for a hodoscope paddle // +// // +// Contains plane, counter and pos/neg adc and tdc values // +// // +/////////////////////////////////////////////////////////////////////////////// + +#include "THcShowerHit.h" + +using namespace std; + + +void THcShowerHit::SetData(Int_t signal, Int_t data) { + if(signal==0) { + fADC_pos = data; + } else if (signal==1) { + fADC_neg = data; + } +} + +Int_t THcShowerHit::GetData(Int_t signal) { + if(signal==0) { + return(fADC_pos); + } else if (signal==1) { + return(fADC_neg); + } + + return(-1); // Actually should throw exception +} + +#if 0 +Int_t THcShowerHit::Compare(const TObject* obj) const +{ + // Compare to sort by plane and counter + + const THcShowerHit* hit = dynamic_cast<const THcShowerHit*>(obj); + + if(!hit) return -1; + Int_t p1 = fPlane; + Int_t p2 = hit->fPlane; + if(p1 < p2) return -1; + else if(p1 > p2) return 1; + else { + Int_t c1 = fCounter; + Int_t c2 = hit->fCounter; + if(c1 < c2) return -1; + else if (c1 == c2) return 0; + else return 1; + } +} +#endif +//_____________________________________________________________________________ +THcShowerHit& THcShowerHit::operator=( const THcShowerHit& rhs ) +{ + // Assignment operator. + + THcRawHit::operator=(rhs); + if ( this != &rhs ) { + fPlane = rhs.fPlane; + fCounter = rhs.fCounter; + fADC_pos = rhs.fADC_pos; + fADC_neg = rhs.fADC_neg; + + } + return *this; +} + + +////////////////////////////////////////////////////////////////////////// +ClassImp(THcShowerHit) + + diff --git a/src/THcShowerHit.h b/src/THcShowerHit.h new file mode 100644 index 0000000000000000000000000000000000000000..7d1a5f08ee650cf8eb3a96dd82ce4a9ff414171f --- /dev/null +++ b/src/THcShowerHit.h @@ -0,0 +1,36 @@ +#ifndef ROOT_THcShowerHit +#define ROOT_THcShowerHit + +#include "THcRawHit.h" + +class THcShowerHit : public THcRawHit { + + public: + + THcShowerHit(Int_t plane=0, Int_t counter=0) : THcRawHit(plane, counter), + fADC_pos(-1), fADC_neg(-1){ + } + THcShowerHit& operator=( const THcShowerHit& ); + virtual ~THcShowerHit() {} + + virtual void Clear( Option_t* opt="" ) + { fADC_pos = -1; fADC_neg = -1; } + + void SetData(Int_t signal, Int_t data); + Int_t GetData(Int_t signal); + + // virtual Bool_t IsSortable () const {return kTRUE; } + // virtual Int_t Compare(const TObject* obj) const; + + Int_t fADC_pos; + Int_t fADC_neg; + + protected: + + private: + + ClassDef(THcShowerHit, 0); // Shower hit class +}; + +#endif +