diff --git a/src/THcShowerArray.cxx b/src/THcShowerArray.cxx index 57a231f2a078c8b7604317d903a3e08bfbae08ed..39944a40f4d78b3e4a67b51fd8479d22e8f6aa89 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 a7e2a0e4bcca8aa16ea80a17f5283737a05a74e9..8e45937cb1e1e15c7c38c007cc4a7f3ea289638c 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 53b1c824a05403dc7728917d71be62443a861be5..866b981d4c8b721d7ffb22ddf7e13a6c7478c48e 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 1c84d1a4cb5a48764e83dda78deea3c67731bba8..6a029e13ee808848d1b85d63eb58db573929c8a1 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;