diff --git a/src/THcShower.cxx b/src/THcShower.cxx
index 6c7d1204714b23d30f5402db3c0e8b91fdbd50a4..60fc502245cc358e27f4cbc122b1442e4eb0aa1a 100644
--- a/src/THcShower.cxx
+++ b/src/THcShower.cxx
@@ -209,9 +209,12 @@ Int_t THcShower::ReadDatabase( const TDatime& date )
{
DBRequest list[]={
{"cal_num_neg_columns", &fNegCols, kInt, 0, 1},
- {"cal_slop", &fSlop, kDouble},
+ {"cal_slop", &fSlop, kDouble,0,1},
{"cal_fv_test", &fvTest, kInt,0,1},
- {"cal_fv_delta", &fvDelta, kDouble},
+ {"cal_fv_delta", &fvDelta, kDouble,0,1},
+ {"cal_ADCMode", &fADCMode, kInt, 0, 1},
+ {"cal_AdcTimeWindowMin", &fAdcTimeWindowMin, kDouble, 0, 1},
+ {"cal_AdcTimeWindowMax", &fAdcTimeWindowMax, kDouble, 0, 1},
{"dbg_raw_cal", &fdbg_raw_cal, kInt, 0, 1},
{"dbg_decoded_cal", &fdbg_decoded_cal, kInt, 0, 1},
{"dbg_sparsified_cal", &fdbg_sparsified_cal, kInt, 0, 1},
@@ -227,7 +230,9 @@ Int_t THcShower::ReadDatabase( const TDatime& date )
fdbg_clusters_cal = 0;
fdbg_tracks_cal = 0;
fdbg_init_cal = 0;
-
+ fAdcTimeWindowMin=0;
+ fAdcTimeWindowMax=10000;
+ fADCMode=kADCDynamicPedestal;
gHcParms->LoadParmValues((DBRequest*)&list, prefix);
}
@@ -343,6 +348,10 @@ Int_t THcShower::ReadDatabase( const TDatime& date )
//Pedestal limits from hcal.param.
fShPosPedLimit = new Int_t [fNTotBlocks];
fShNegPedLimit = new Int_t [fNTotBlocks];
+ for (UInt_t i;i<fNTotBlocks;i++) {
+ fShPosPedLimit[i]=0.;
+ fShNegPedLimit[i]=0.;
+ }
//Calibration constants
fPosGain = new Double_t [fNTotBlocks];
@@ -350,31 +359,22 @@ Int_t THcShower::ReadDatabase( const TDatime& date )
//Read in parameters from hcal.param
Double_t hcal_pos_cal_const[fNTotBlocks];
- // Double_t hcal_pos_gain_ini[fNTotBlocks]; not used
- // Double_t hcal_pos_gain_cur[fNTotBlocks]; not used
- // Int_t hcal_pos_ped_limit[fNTotBlocks]; not used
Double_t hcal_pos_gain_cor[fNTotBlocks];
Double_t hcal_neg_cal_const[fNTotBlocks];
- // Double_t hcal_neg_gain_ini[fNTotBlocks]; not used
- // Double_t hcal_neg_gain_cur[fNTotBlocks]; not used
- // Int_t hcal_neg_ped_limit[fNTotBlocks]; not used
Double_t hcal_neg_gain_cor[fNTotBlocks];
DBRequest list[]={
{"cal_pos_cal_const", hcal_pos_cal_const, kDouble, fNTotBlocks},
- // {"cal_pos_gain_ini", hcal_pos_gain_ini, kDouble, fNTotBlocks},
- // {"cal_pos_gain_cur", hcal_pos_gain_cur, kDouble, fNTotBlocks},
- {"cal_pos_ped_limit", fShPosPedLimit, kInt, fNTotBlocks},
+ {"cal_pos_ped_limit", fShPosPedLimit, kInt, fNTotBlocks,1},
{"cal_pos_gain_cor", hcal_pos_gain_cor, kDouble, fNTotBlocks},
{"cal_neg_cal_const", hcal_neg_cal_const, kDouble, fNTotBlocks},
- // {"cal_neg_gain_ini", hcal_neg_gain_ini, kDouble, fNTotBlocks},
- // {"cal_neg_gain_cur", hcal_neg_gain_cur, kDouble, fNTotBlocks},
- {"cal_neg_ped_limit", fShNegPedLimit, kInt, fNTotBlocks},
+ {"cal_neg_ped_limit", fShNegPedLimit, kInt, fNTotBlocks,1},
{"cal_neg_gain_cor", hcal_neg_gain_cor, kDouble, fNTotBlocks},
- {"cal_min_peds", &fShMinPeds, kInt},
+ {"cal_min_peds", &fShMinPeds, kInt,0,1},
{0}
};
+ fShMinPeds=0.;
gHcParms->LoadParmValues((DBRequest*)&list, prefix);
// Debug output.
@@ -507,9 +507,13 @@ Int_t THcShower::DefineVariables( EMode mode )
RVarDef vars[] = {
{ "nhits", "Number of hits", "fNhits" },
{ "nclust", "Number of clusters", "fNclust" },
+ { "nclusttrack", "Number of cluster which matched best track","fNclustTrack" },
+ { "xclusttrack", "X pos of cluster which matched best track","fXclustTrack" },
+ { "xtrack", "X pos of track which matched cluster", "fXTrack" },
{ "etot", "Total energy", "fEtot" },
{ "etotnorm", "Total energy divided by Central Momentum", "fEtotNorm" },
{ "etrack", "Track energy", "fEtrack" },
+ { "etracknorm", "Total energy divided by track momentum", "fEtrackNorm" },
{ "ntracks", "Number of shower tracks", "fNtracks" },
{ 0 }
};
@@ -560,10 +564,14 @@ void THcShower::Clear(Option_t* opt)
fNhits = 0;
fNclust = 0;
+ fNclustTrack = -2;
+ fXclustTrack = -1000;
+ fXTrack = -1000;
fNtracks = 0;
fEtot = 0.;
fEtotNorm = 0.;
fEtrack = 0.;
+ fEtrackNorm = 0.;
// Purge cluster list
@@ -584,7 +592,7 @@ Int_t THcShower::Decode( const THaEvData& evdata )
// Get the Hall C style hitlist (fRawHitList) for this event
Int_t nhits = DecodeToHitList(evdata);
-
+
fEvent = evdata.GetEvNum();
if(gHaCuts->Result("Pedestal_event")) {
@@ -612,14 +620,10 @@ Int_t THcShower::Decode( const THaEvData& evdata )
Int_t nexthit = 0;
for(UInt_t ip=0;ip<fNLayers;ip++) {
nexthit = fPlanes[ip]->ProcessHits(fRawHitList, nexthit);
- fEtot += fPlanes[ip]->GetEplane();
}
if(fHasArray) {
nexthit = fArray->ProcessHits(fRawHitList, nexthit);
- fEtot += fArray->GetEarray();
}
- THcHallCSpectrometer *app = static_cast<THcHallCSpectrometer*>(GetApparatus());
- fEtotNorm=fEtot/(app->GetPcentral());
return nhits;
}
@@ -637,27 +641,25 @@ Int_t THcShower::CoarseProcess( TClonesArray& tracks)
//
// Fill set of unclustered hits.
-
+ for(UInt_t ip=0;ip<fNLayers;ip++) {
+ fPlanes[ip]->CoarseProcessHits();
+ fEtot += fPlanes[ip]->GetEplane();
+ }
+ if(fHasArray) {
+ fArray->CoarseProcessHits();
+ fEtot += fArray->GetEarray();
+ }
+ THcHallCSpectrometer *app = static_cast<THcHallCSpectrometer*>(GetApparatus());
+ fEtotNorm=fEtot/(app->GetPcentral());
+ //
THcShowerHitSet HitSet;
for(UInt_t j=0; j < fNLayers; j++) {
for (UInt_t i=0; i<fNBlocks[j]; i++) {
- //May be should be done this way.
- //
- // Double_t Edep = fPlanes[j]->GetEmean(i);
- // if (Edep > 0.) { //hit
- // Double_t x = fYPos[j][i] + BlockThick[j]/2.; //top + thick/2
- // Double_t z = fLayerZPos[j] + BlockThick[j]/2.;//front + thick/2
- // THcShowerHit* hit = new THcShowerHit(i,j,x,z,Edep);
-
- //ENGINE way.
//
- if (fPlanes[j]->GetApos(i) - fPlanes[j]->GetPosPed(i) >
- fPlanes[j]->GetPosThr(i) - fPlanes[j]->GetPosPed(i) ||
- fPlanes[j]->GetAneg(i) - fPlanes[j]->GetNegPed(i) >
- fPlanes[j]->GetNegThr(i) - fPlanes[j]->GetNegPed(i)) { //hit
+ if (fPlanes[j]->GetAposP(i) > 0 || fPlanes[j]->GetAnegP(i) >0) { //hit
Double_t Edep = fPlanes[j]->GetEmean(i);
Double_t Epos = fPlanes[j]->GetEpos(i);
@@ -682,6 +684,7 @@ Int_t THcShower::CoarseProcess( TClonesArray& tracks)
cout << "Debug output from THcShower::CoarseProcess for "
<< GetApparatus()->GetName() << endl;
+ cout << " event = " << fEvent << endl;
cout << " List of unclustered hits. Total hits: " << fNhits << endl;
THcShowerHitIt it = HitSet.begin(); //<set> version
for (Int_t i=0; i!=fNhits; i++) {
@@ -700,7 +703,7 @@ Int_t THcShower::CoarseProcess( TClonesArray& tracks)
if (fdbg_clusters_cal) {
- cout << " Clustered hits. Number of clusters: " << fNclust << endl;
+ cout << " event = " << fEvent << " Clustered hits. Number of clusters: " << fNclust << endl;
UInt_t i = 0;
for (THcShowerClusterListIt ppcl = (*fClusterList).begin();
@@ -729,7 +732,18 @@ Int_t THcShower::CoarseProcess( TClonesArray& tracks)
// Do same for Array.
if(fHasArray) fArray->CoarseProcess(tracks);
+ //
+ Int_t Ntracks = tracks.GetLast()+1; // Number of reconstructed tracks
+ Double_t save_energy=0;
+ for (Int_t itrk=0; itrk<Ntracks; itrk++) {
+ THaTrack* theTrack = static_cast<THaTrack*>( tracks[itrk] );
+ save_energy = GetShEnergy(theTrack);
+ if (fHasArray) save_energy += fArray->GetShEnergy(theTrack);
+ theTrack->SetEnergy(save_energy);
+ } //over tracks
+
+ //
return 0;
}
@@ -761,7 +775,6 @@ void THcShower::ClusterHits(THcShowerHitSet& HitSet,
k++) {
if ((**i).isNeighbour(*k)) {
-
(*cluster).insert(*i); //If the hit #i is neighbouring a hit
HitSet.erase(i); //in the cluster, then move it
//into the cluster.
@@ -945,7 +958,7 @@ Int_t THcShower::MatchCluster(THaTrack* Track,
Double_t dx = TMath::Abs( clX(cluster) - XTrFront );
if (dx <= (0.5*BlockThick[0] + fSlop)) {
- fNtracks++; // number of shower tracks (Consistent with engine)
+ fNtracks++; // lumber of shower tracks (Consistent with engine)
if (dx < deltaX) {
mclust = i;
deltaX = dx;
@@ -960,6 +973,7 @@ Int_t THcShower::MatchCluster(THaTrack* Track,
cout << "---------------------------------------------------------------\n";
cout << "Debug output from THcShower::MatchCluster for "
<< GetApparatus()->GetName() << endl;
+ cout << " event = " << fEvent << endl;
cout << " Track at DC:"
<< " X = " << Track->GetX()
@@ -1026,6 +1040,7 @@ Float_t THcShower::GetShEnergy(THaTrack* Track) {
corneg = 0.;
}
+ // cout << ip << " clust energy pos = " << clEplane(cluster,ip,0)<< " clust energy pos = " << clEplane(cluster,ip,1) << endl;
Etrk += clEplane(cluster,ip,0) * corpos;
Etrk += clEplane(cluster,ip,1) * corneg;
@@ -1039,6 +1054,7 @@ Float_t THcShower::GetShEnergy(THaTrack* Track) {
cout << "---------------------------------------------------------------\n";
cout << "Debug output from THcShower::GetShEnergy for "
<< GetApparatus()->GetName() << endl;
+ cout << " event = " << fEvent << endl;
cout << " Track at the calorimeter: X = " << Xtr << " Y = " << Ytr;
if (mclust >= 0)
@@ -1061,15 +1077,23 @@ Int_t THcShower::FineProcess( TClonesArray& tracks )
//
Int_t Ntracks = tracks.GetLast()+1; // Number of reconstructed tracks
-
for (Int_t itrk=0; itrk<Ntracks; itrk++) {
-
THaTrack* theTrack = static_cast<THaTrack*>( tracks[itrk] );
-
- fEtrack = GetShEnergy(theTrack);
- if (fHasArray) fEtrack += fArray->GetShEnergy(theTrack);
- theTrack->SetEnergy(fEtrack);
-
+ if (theTrack->GetIndex()==0) {
+ fEtrack=theTrack->GetEnergy();
+ fEtrackNorm=fEtrack/theTrack->GetP();
+
+ Double_t Xtr = -100.;
+ Double_t Ytr = -100.;
+ fNclustTrack = MatchCluster(theTrack, Xtr, Ytr);
+ fXTrack=Xtr;
+ if (fNclustTrack>=0) {
+ THcShowerCluster* cluster = *(fClusterList->begin()+fNclustTrack);
+ Double_t dx = TMath::Abs( clX(cluster) - Xtr );
+ fXclustTrack=clX(cluster);
+ //cout << fNclustTrack << " " << Xtr << " " << clX(cluster) << " " << dx << " " << Ytr << " " << fEtrack<< endl;
+ }
+ }
} //over tracks
//Debug output.
@@ -1078,7 +1102,7 @@ Int_t THcShower::FineProcess( TClonesArray& tracks )
cout << "---------------------------------------------------------------\n";
cout << "Debug output from THcShower::FineProcess for "
<< GetApparatus()->GetName() << endl;
-
+ cout << " event = " << fEvent << endl;
cout << " Number of tracks = " << Ntracks << endl;
for (Int_t itrk=0; itrk<Ntracks; itrk++) {
@@ -1088,7 +1112,8 @@ Int_t THcShower::FineProcess( TClonesArray& tracks )
<< " Y = " << theTrack->GetY()
<< " Theta = " << theTrack->GetTheta()
<< " Phi = " << theTrack->GetPhi()
- << " Energy = " << theTrack->GetEnergy() << endl;
+ << " Energy = " << theTrack->GetEnergy()
+ << " Energy/Ptrack = " << fEtrackNorm << endl;
}
cout << "---------------------------------------------------------------\n";
diff --git a/src/THcShower.h b/src/THcShower.h
index b4ad6aa74a1d423f7546efa18731d8250a8d10fe..5770cb497761bfcc5b9854ec689b2e1420703559 100644
--- a/src/THcShower.h
+++ b/src/THcShower.h
@@ -71,6 +71,15 @@ public:
return ( Side == 0 ? fPosGain[nelem] : fNegGain[nelem]);
}
+ Int_t GetADCMode() {
+ return fADCMode;
+ }
+ Double_t GetAdcTimeWindowMin() {
+ return fAdcTimeWindowMin;
+ }
+ Double_t GetAdcTimeWindowMax() {
+ return fAdcTimeWindowMax;
+ }
Int_t GetMinPeds() {
return fShMinPeds;
}
@@ -78,6 +87,9 @@ public:
Int_t GetNLayers() {
return fNLayers;
}
+ Int_t GetNBlocks(Int_t layer) {
+ return fNBlocks[layer];
+ }
//Coordinate correction for single PMT modules.
//PMT attached at right (positive) side.
@@ -109,6 +121,16 @@ public:
protected:
Int_t fEvent;
+ Int_t fADCMode; // != 0 if using FADC
+ // 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 kADCDynamicPedestal=1;
+ static const Int_t kADCSampleIntegral=2;
+ static const Int_t kADCSampIntDynPed=3;
+ Double_t fAdcTimeWindowMin;
+ Double_t fAdcTimeWindowMax;
Int_t fAnalyzePedestals; // Flag for pedestal analysis.
@@ -124,11 +146,15 @@ protected:
Int_t fNhits; // Total number of hits
Int_t fNclust; // Number of clusters
+ Int_t fNclustTrack; // NUmber of cluster that match best track
+ Double_t fXclustTrack; // X pos of cluster that match best track
+ Double_t fXTrack; // X pos of best track that match cluster
Int_t fNtracks; // Number of shower tracks, i.e. number of
// cluster-to-track association
Double_t fEtot; // Total energy
Double_t fEtotNorm; // Total energy divided by spec central momentum
- Double_t fEtrack; // Cluster energy associated to the last track
+ Double_t fEtrack; // Cluster energy associated to the best track
+ Double_t fEtrackNorm; // Cluster energy associated to the best track
THcShowerClusterList* fClusterList; // List of hit clusters
diff --git a/src/THcShowerArray.cxx b/src/THcShowerArray.cxx
index da4257281a27931d590bfcb4e033bacb37a8d460..63e53574f0b190c6a34ed29b0adbf4cba353c486 100644
--- a/src/THcShowerArray.cxx
+++ b/src/THcShowerArray.cxx
@@ -39,7 +39,8 @@ THcShowerArray::THcShowerArray( const char* name,
fADCHits = new TClonesArray("THcSignalHit",100);
fLayerNum = layernum;
- frAdcPedRaw = new TClonesArray("THcSignalHit", 16);
+ frAdcPedRaw = new TClonesArray("THcSignalHit", 16);
+ frAdcErrorFlag = new TClonesArray("THcSignalHit", 16);
frAdcPulseIntRaw = new TClonesArray("THcSignalHit", 16);
frAdcPulseAmpRaw = new TClonesArray("THcSignalHit", 16);
frAdcPulseTimeRaw = new TClonesArray("THcSignalHit", 16);
@@ -61,6 +62,7 @@ THcShowerArray::~THcShowerArray()
delete fADCHits;
delete frAdcPedRaw; frAdcPedRaw = NULL;
+ delete frAdcErrorFlag; frAdcErrorFlag = NULL;
delete frAdcPulseIntRaw; frAdcPulseIntRaw = NULL;
delete frAdcPulseAmpRaw; frAdcPulseAmpRaw = NULL;
delete frAdcPulseTimeRaw; frAdcPulseTimeRaw = NULL;
@@ -123,6 +125,9 @@ Int_t THcShowerArray::ReadDatabase( const TDatime& date )
{"cal_arr_ystep", &fYStep, kDouble},
{"cal_arr_zsize", &fZSize, kDouble},
{"cal_using_fadc", &fUsingFADC, kInt, 0, 1},
+ {"cal_arr_ADCMode", &fADCMode, kInt, 0, 1},
+ {"cal_arr_AdcTimeWindowMin", &fAdcTimeWindowMin, kDouble, 0, 1},
+ {"cal_arr_AdcTimeWindowMax", &fAdcTimeWindowMax, kDouble, 0, 1},
{"cal_ped_sample_low", &fPedSampLow, kInt, 0, 1},
{"cal_ped_sample_high", &fPedSampHigh, kInt, 0, 1},
{"cal_data_sample_low", &fDataSampLow, kInt, 0, 1},
@@ -130,7 +135,9 @@ Int_t THcShowerArray::ReadDatabase( const TDatime& date )
{0}
};
gHcParms->LoadParmValues((DBRequest*)&list, prefix);
-
+ fADCMode=kADCDynamicPedestal;
+ fAdcTimeWindowMin=0;
+ fAdcTimeWindowMax=10000;
fNelem = fNRows*fNColumns;
fXPos = new Double_t* [fNRows];
@@ -216,10 +223,9 @@ Int_t THcShowerArray::ReadDatabase( const TDatime& date )
Double_t cal_arr_gain_cor[fNelem];
DBRequest list1[]={
- {"cal_arr_ped_limit", fPedLimit, kInt, static_cast<UInt_t>(fNelem)},
+ {"cal_arr_ped_limit", fPedLimit, kInt, static_cast<UInt_t>(fNelem),1},
{"cal_arr_cal_const", cal_arr_cal_const, kDouble, static_cast<UInt_t>(fNelem)},
{"cal_arr_gain_cor", cal_arr_gain_cor, kDouble, static_cast<UInt_t>(fNelem)},
- // {"cal_min_peds", &fShMinPeds, kInt},
{0}
};
gHcParms->LoadParmValues((DBRequest*)&list1, prefix);
@@ -261,7 +267,7 @@ Int_t THcShowerArray::ReadDatabase( const TDatime& date )
// Calibration constants (GeV / ADC channel).
fGain = new Double_t [fNelem];
- for (UInt_t i=0; i<fNelem; i++) {
+ for (Int_t i=0; i<fNelem; i++) {
fGain[i] = cal_arr_cal_const[i] * cal_arr_gain_cor[i];
}
@@ -328,6 +334,7 @@ Int_t THcShowerArray::DefineVariables( EMode mode )
{"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"},
{"earray", "Energy Deposition in array", "fEarray"},
@@ -336,6 +343,7 @@ Int_t THcShowerArray::DefineVariables( EMode mode )
{"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()"},
@@ -357,6 +365,7 @@ void THcShowerArray::Clear( Option_t* )
fADCHits->Clear();
fNhits = 0;
+ fNgoodhits = 0;
fNclust = 0;
fNtracks = 0;
@@ -368,6 +377,7 @@ void THcShowerArray::Clear( Option_t* )
fClusterList->clear();
frAdcPedRaw->Clear();
+ frAdcErrorFlag->Clear();
frAdcPulseIntRaw->Clear();
frAdcPulseAmpRaw->Clear();
frAdcPulseTimeRaw->Clear();
@@ -413,9 +423,6 @@ Int_t THcShowerArray::CoarseProcess( TClonesArray& tracks )
k++;
}
}
-
- fNhits = HitSet.size();
-
//Debug output, print out hits before clustering.
THcShower* fParent = (THcShower*) GetParent();
@@ -427,7 +434,7 @@ Int_t THcShowerArray::CoarseProcess( TClonesArray& tracks )
cout << " List of unclustered hits. Total hits: " << fNhits << endl;
THcShowerHitIt it = HitSet.begin(); //<set> version
- for (Int_t i=0; i!=fNhits; i++) {
+ for (Int_t i=0; i!=fNgoodhits; i++) {
cout << " hit " << i << ": ";
(*(it++))->show();
}
@@ -643,6 +650,96 @@ Int_t THcShowerArray::FineProcess( TClonesArray& tracks )
return 0;
}
+//_____________________________________________________________________________
+Int_t THcShowerArray::CoarseProcessHits()
+{
+ THcShower* fParent;
+ fParent = (THcShower*) GetParent();
+ Int_t ADCMode=fParent->GetADCMode();
+ if(ADCMode == kADCDynamicPedestal) {
+ FillADC_DynamicPedestal();
+ } else if (ADCMode == kADCSampleIntegral) {
+ FillADC_SampleIntegral();
+ } else if (ADCMode == kADCSampIntDynPed) {
+ FillADC_SampIntDynPed();
+ } else {
+ FillADC_Standard();
+ }
+ //
+ if (fParent->fdbg_decoded_cal) {
+
+ cout << "---------------------------------------------------------------\n";
+ cout << "Debug output from THcShowerArray::ProcessHits for "
+ << fParent->GetPrefix() << ":" << endl;
+
+ cout << " Sparsified hits for shower array, plane #" << fLayerNum
+ << ", " << GetName() << ":" << endl;
+
+ Int_t nspar = 0;
+ 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]) {
+ cout << " counter = " << k
+ << " E = " << fE[k]
+ << endl;
+ nspar++;
+ }
+ k++;
+ }
+ }
+
+ if (nspar == 0) cout << " No hits\n";
+
+ cout << " Earray = " << fEarray << endl;
+ cout << "---------------------------------------------------------------\n";
+ }
+ //
+ return 1;
+}
+//_____________________________________________________________________________
+void THcShowerArray::FillADC_SampIntDynPed()
+{
+ // adc_pos = hit->GetRawAdcHitPos().GetSampleInt();
+ // adc_neg = hit->GetRawAdcHitNeg().GetSampleInt();
+ // adc_pos_pedsub = hit->GetRawAdcHitPos().GetSampleIntRaw();
+ // adc_neg_pedsub = hit->GetRawAdcHitNeg().GetSampleIntRaw();
+ // Need to create
+}
+//_____________________________________________________________________________
+void THcShowerArray::FillADC_SampleIntegral()
+{
+ /// adc_pos_pedsub = hit->GetRawAdcHitPos().GetSampleIntRaw() - fPosPed[hit->fCounter -1];
+ // adc_neg_pedsub = hit->GetRawAdcHitNeg().GetSampleIntRaw() - fNegPed[hit->fCounter -1];
+ // adc_pos = hit->GetRawAdcHitPos().GetSampleIntRaw();
+ // adc_neg = hit->GetRawAdcHitNeg().GetSampleIntRaw();
+ // need to create
+}
+//_____________________________________________________________________________
+void THcShowerArray::FillADC_Standard()
+{
+}
+//_____________________________________________________________________________
+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 pulseTime = ((THcSignalHit*) frAdcPulseTimeRaw->ConstructedAt(ielem))->GetData();
+ Double_t errorflag = ((THcSignalHit*) frAdcErrorFlag->ConstructedAt(ielem))->GetData();
+ Bool_t pulseTimeCut = (pulseTime > fAdcTimeWindowMin) && (pulseTime < fAdcTimeWindowMax);
+ if (errorflag==0 && pulseTimeCut) {
+ fA[npad] =pulseIntRaw;
+ if(fA[npad] > fThresh[npad]) {
+ fA_p[npad] =pulseInt ;
+ fE[npad] = fA_p[npad]*fGain[npad];
+ fEarray += fE[npad];
+ }
+ }
+ }
+ //
+}
//_____________________________________________________________________________
Int_t THcShowerArray::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
{
@@ -656,6 +753,7 @@ Int_t THcShowerArray::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
fADCHits->Clear();
frAdcPedRaw->Clear();
+ frAdcErrorFlag->Clear();
frAdcPulseIntRaw->Clear();
frAdcPulseAmpRaw->Clear();
frAdcPulseTimeRaw->Clear();
@@ -679,9 +777,6 @@ Int_t THcShowerArray::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
Int_t ihit = nexthit;
- Int_t ngood = 0;
- Int_t threshold = 100;
-
UInt_t nrAdcHits = 0;
while(ihit < nrawhits) {
@@ -692,11 +787,12 @@ Int_t THcShowerArray::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
}
Int_t padnum = hit->fCounter;
-
THcRawAdcHit& rawAdcHit = hit->GetRawAdcHitPos();
+ //
for (UInt_t thit=0; thit<rawAdcHit.GetNPulses(); ++thit) {
((THcSignalHit*) frAdcPedRaw->ConstructedAt(nrAdcHits))->Set(padnum, rawAdcHit.GetPedRaw());
- ((THcSignalHit*) frAdcPed->ConstructedAt(nrAdcHits))->Set(padnum, rawAdcHit.GetPed());
+ fThresh[padnum]=rawAdcHit.GetPedRaw()*rawAdcHit.GetF250_PeakPedestalRatio()+250.;
+ ((THcSignalHit*) frAdcPed->ConstructedAt(nrAdcHits))->Set(padnum, rawAdcHit.GetPed());
((THcSignalHit*) frAdcPulseIntRaw->ConstructedAt(nrAdcHits))->Set(padnum, rawAdcHit.GetPulseIntRaw(thit));
((THcSignalHit*) frAdcPulseInt->ConstructedAt(nrAdcHits))->Set(padnum, rawAdcHit.GetPulseInt(thit));
@@ -705,48 +801,18 @@ Int_t THcShowerArray::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
((THcSignalHit*) frAdcPulseAmp->ConstructedAt(nrAdcHits))->Set(padnum, rawAdcHit.GetPulseAmp(thit));
((THcSignalHit*) frAdcPulseTimeRaw->ConstructedAt(nrAdcHits))->Set(padnum, rawAdcHit.GetPulseTimeRaw(thit));
-
+ if (rawAdcHit.GetPulseAmp(thit)>0&&rawAdcHit.GetPulseIntRaw(thit)>0) {
+ ((THcSignalHit*) frAdcErrorFlag->ConstructedAt(nrAdcHits))->Set(padnum,0);
+ } else {
+ ((THcSignalHit*) frAdcErrorFlag->ConstructedAt(nrAdcHits))->Set(padnum,1);
+ }
++nrAdcHits;
}
- fADCMode=kADCDynamicPedestal;
- Double_t adc;
- Double_t adc_pedsub;
- if(fADCMode == kADCDynamicPedestal) {
- adc_pedsub = hit->GetRawAdcHitPos().GetPulseInt();
- adc= hit->GetRawAdcHitPos().GetPulseIntRaw();
- } else if (fADCMode == kADCSampleIntegral) {
- adc_pedsub = hit->GetRawAdcHitPos().GetSampleIntRaw() - fPed[hit->fCounter -1];
- adc = hit->GetRawAdcHitPos().GetSampleIntRaw();
- } else if (fADCMode == kADCSampIntDynPed) {
- adc = hit->GetRawAdcHitPos().GetSampleInt();
- adc_pedsub = hit->GetRawAdcHitPos().GetSampleIntRaw();
- } else {
- adc_pedsub = hit->GetRawAdcHitPos().GetPulseIntRaw()-fPed[hit->fCounter -1];
- adc = hit->GetRawAdcHitPos().GetPulseIntRaw();
- }
-
- fA[hit->fCounter-1] = adc;
- threshold=hit->GetRawAdcHitPos().GetPedRaw()+100;
- if(fA[hit->fCounter-1] > threshold) {
- ngood++;
- }
-
- // Sparsify hits, fill the hit list, compute the energy depostion.
- fThresh[hit->fCounter -1] = hit->GetRawAdcHitPos().GetPedRaw()+100;
- if(fA[hit->fCounter-1] > fThresh[hit->fCounter -1]) {
- fA_p[hit->fCounter-1] = adc_pedsub;
- fE[hit->fCounter-1] += fA_p[hit->fCounter-1] * fGain[hit->fCounter-1];
- }
-
- // Accumulate energies in the plane.
-
- fEarray += fE[hit->fCounter-1];
-
ihit++;
}
#if 0
- if(ngood > 0) {
+ if(fNgoodhits > 0) {
cout << "+";
for(Int_t column=0;column<fNColumns;column++) {
cout << "-";
@@ -767,7 +833,7 @@ Int_t THcShowerArray::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
}
#endif
#ifdef HITPIC
- if(ngood > 0) {
+ if(fNgoodhits > 0) {
for(Int_t row=0;row<fNRows;row++) {
if(piccolumn==0) {
hitpic[row][0] = '|';
@@ -803,37 +869,6 @@ Int_t THcShowerArray::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
//Debug output.
- if (fParent->fdbg_decoded_cal) {
-
- cout << "---------------------------------------------------------------\n";
- cout << "Debug output from THcShowerArray::ProcessHits for "
- << fParent->GetPrefix() << ":" << endl;
-
- cout << " nrawhits = " << nrawhits << " nexthit = " << nexthit << endl;
- cout << " Sparsified hits for shower array, plane #" << fLayerNum
- << ", " << GetName() << ":" << endl;
-
- Int_t nspar = 0;
- for (Int_t jhit = nexthit; jhit < nrawhits; jhit++) {
-
- THcRawShowerHit* hit = (THcRawShowerHit *) rawhits->At(jhit);
- if(hit->fPlane != fLayerNum) {
- break;
- }
-
- if(fA[hit->fCounter-1] > fThresh[hit->fCounter -1]) {
- cout << " counter = " << hit->fCounter
- << " E = " << fE[hit->fCounter-1]
- << endl;
- nspar++;
- }
- }
-
- if (nspar == 0) cout << " No hits\n";
-
- cout << " Earray = " << fEarray << endl;
- cout << "---------------------------------------------------------------\n";
- }
return(ihit);
}
diff --git a/src/THcShowerArray.h b/src/THcShowerArray.h
index ad0167030189296b445014c8880208e9f00f461a..4b178cbcf9c43ff304b32b7bce40ce3d286ee9e2 100644
--- a/src/THcShowerArray.h
+++ b/src/THcShowerArray.h
@@ -45,10 +45,15 @@ public:
virtual Bool_t IsPid() { return kFALSE; }
virtual Int_t ProcessHits(TClonesArray* rawhits, Int_t nexthit);
+ virtual Int_t CoarseProcessHits();
virtual Int_t AccumulatePedestals(TClonesArray* rawhits, Int_t nexthit);
virtual void CalculatePedestals( );
virtual void InitializePedestals( );
-
+ virtual void FillADC_DynamicPedestal( );
+ virtual void FillADC_SampleIntegral( );
+ virtual void FillADC_SampIntDynPed( );
+ virtual void FillADC_Standard( );
+
// Cluster to track association method.
Int_t MatchCluster(THaTrack*, Double_t&, Double_t&);
@@ -106,7 +111,9 @@ protected:
static const Int_t kADCDynamicPedestal=1;
static const Int_t kADCSampleIntegral=2;
static const Int_t kADCSampIntDynPed=3;
- Int_t fPedSampLow; // Sample range for
+ Double_t fAdcTimeWindowMin ;
+ Double_t fAdcTimeWindowMax ;
+Int_t fPedSampLow; // Sample range for
Int_t fPedSampHigh; // dynamic pedestal
Int_t fDataSampLow; // Sample range for
Int_t fDataSampHigh; // sample integration
@@ -137,6 +144,7 @@ protected:
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
@@ -144,6 +152,7 @@ protected:
THcShowerClusterList* fClusterList; // List of hit clusters
TClonesArray* frAdcPedRaw;
+ TClonesArray* frAdcErrorFlag;
TClonesArray* frAdcPulseIntRaw;
TClonesArray* frAdcPulseAmpRaw;
TClonesArray* frAdcPulseTimeRaw;
diff --git a/src/THcShowerHit.cxx b/src/THcShowerHit.cxx
index 574d05646a1c4b567ec69ca4bac4f14298dfdaaa..e7bc7bf88428e0cafacfe7ce954ec8189cf76773 100644
--- a/src/THcShowerHit.cxx
+++ b/src/THcShowerHit.cxx
@@ -44,7 +44,7 @@ bool THcShowerHit::isNeighbour(THcShowerHit* hit1) {
//Print out hit information
//
void THcShowerHit::show() {
- cout << "row=" << fRow << " column=" << fCol
+ cout << "row=" << fRow+1 << " column=" << fCol+1
<< " x=" << fX << " z=" << fZ
<< " E=" << fE << " Epos=" << fEpos << " Eneg=" << fEneg << endl;
}
diff --git a/src/THcShowerPlane.cxx b/src/THcShowerPlane.cxx
index 44eba1a16fea009b9fb5a524d01e0277ad622d4c..0a9452d50ab57599e3a35bcd45242ecb2c193fae 100644
--- a/src/THcShowerPlane.cxx
+++ b/src/THcShowerPlane.cxx
@@ -39,6 +39,7 @@ THcShowerPlane::THcShowerPlane( const char* name,
fPosADCHits = new TClonesArray("THcSignalHit",fNelem);
fNegADCHits = new TClonesArray("THcSignalHit",fNelem);
+ frPosAdcErrorFlag = new TClonesArray("THcSignalHit", 16);
frPosAdcPedRaw = new TClonesArray("THcSignalHit", 16);
frPosAdcPulseIntRaw = new TClonesArray("THcSignalHit", 16);
frPosAdcPulseAmpRaw = new TClonesArray("THcSignalHit", 16);
@@ -48,6 +49,7 @@ THcShowerPlane::THcShowerPlane( const char* name,
frPosAdcPulseInt = new TClonesArray("THcSignalHit", 16);
frPosAdcPulseAmp = new TClonesArray("THcSignalHit", 16);
+ frNegAdcErrorFlag = new TClonesArray("THcSignalHit", 16);
frNegAdcPedRaw = new TClonesArray("THcSignalHit", 16);
frNegAdcPulseIntRaw = new TClonesArray("THcSignalHit", 16);
frNegAdcPulseAmpRaw = new TClonesArray("THcSignalHit", 16);
@@ -72,6 +74,7 @@ THcShowerPlane::~THcShowerPlane()
delete fPosADCHits;
delete fNegADCHits;
+ frPosAdcErrorFlag = NULL;
frPosAdcPedRaw = NULL;
frPosAdcPulseIntRaw = NULL;
frPosAdcPulseAmpRaw = NULL;
@@ -81,6 +84,7 @@ THcShowerPlane::~THcShowerPlane()
frPosAdcPulseInt = NULL;
frPosAdcPulseAmp = NULL;
+ frNegAdcErrorFlag = NULL;
frNegAdcPedRaw = NULL;
frNegAdcPulseIntRaw = NULL;
frNegAdcPulseAmpRaw = NULL;
@@ -256,6 +260,7 @@ Int_t THcShowerPlane::DefineVariables( EMode mode )
{"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()"},
@@ -265,6 +270,7 @@ Int_t THcShowerPlane::DefineVariables( EMode mode )
{"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()"},
@@ -287,6 +293,7 @@ void THcShowerPlane::Clear( Option_t* )
fPosADCHits->Clear();
fNegADCHits->Clear();
+ frPosAdcErrorFlag->Clear();
frPosAdcPedRaw->Clear();
frPosAdcPulseIntRaw->Clear();
frPosAdcPulseAmpRaw->Clear();
@@ -296,6 +303,7 @@ void THcShowerPlane::Clear( Option_t* )
frPosAdcPulseInt->Clear();
frPosAdcPulseAmp->Clear();
+ frNegAdcErrorFlag->Clear();
frNegAdcPedRaw->Clear();
frNegAdcPulseIntRaw->Clear();
frNegAdcPulseAmpRaw->Clear();
@@ -365,6 +373,7 @@ Int_t THcShowerPlane::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
fPosADCHits->Clear();
fNegADCHits->Clear();
+ frPosAdcErrorFlag->Clear();
frPosAdcPedRaw->Clear();
frPosAdcPulseIntRaw->Clear();
frPosAdcPulseAmpRaw->Clear();
@@ -374,6 +383,7 @@ Int_t THcShowerPlane::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
frPosAdcPulseInt->Clear();
frPosAdcPulseAmp->Clear();
+ frNegAdcErrorFlag->Clear();
frNegAdcPedRaw->Clear();
frNegAdcPulseIntRaw->Clear();
frNegAdcPulseAmpRaw->Clear();
@@ -421,6 +431,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()+250.;
((THcSignalHit*) frPosAdcPed->ConstructedAt(nrPosAdcHits))->Set(padnum, rawPosAdcHit.GetPed());
((THcSignalHit*) frPosAdcPulseIntRaw->ConstructedAt(nrPosAdcHits))->Set(padnum, rawPosAdcHit.GetPulseIntRaw(thit));
@@ -430,12 +441,17 @@ Int_t THcShowerPlane::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
((THcSignalHit*) frPosAdcPulseAmp->ConstructedAt(nrPosAdcHits))->Set(padnum, rawPosAdcHit.GetPulseAmp(thit));
((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);
+ } else {
+ ((THcSignalHit*) frPosAdcErrorFlag->ConstructedAt(nrPosAdcHits))->Set(padnum,1);
+ }
++nrPosAdcHits;
}
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()+250.;
((THcSignalHit*) frNegAdcPed->ConstructedAt(nrNegAdcHits))->Set(padnum, rawNegAdcHit.GetPed());
((THcSignalHit*) frNegAdcPulseIntRaw->ConstructedAt(nrNegAdcHits))->Set(padnum, rawNegAdcHit.GetPulseIntRaw(thit));
@@ -445,107 +461,54 @@ Int_t THcShowerPlane::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
((THcSignalHit*) frNegAdcPulseAmp->ConstructedAt(nrNegAdcHits))->Set(padnum, rawNegAdcHit.GetPulseAmp(thit));
((THcSignalHit*) frNegAdcPulseTimeRaw->ConstructedAt(nrNegAdcHits))->Set(padnum, rawNegAdcHit.GetPulseTimeRaw(thit));
-
+ if (rawNegAdcHit.GetPulseAmp(thit)>0&&rawNegAdcHit.GetPulseIntRaw(thit)>0) {
+ ((THcSignalHit*) frNegAdcErrorFlag->ConstructedAt(nrNegAdcHits))->Set(padnum,0);
+ } else {
+ ((THcSignalHit*) frNegAdcErrorFlag->ConstructedAt(nrNegAdcHits))->Set(padnum,1);
+ }
++nrNegAdcHits;
}
- fADCMode=kADCDynamicPedestal;
- Double_t adc_pos;
- Double_t adc_neg;
- Double_t adc_pos_pedsub;
- Double_t adc_neg_pedsub;
- if(fADCMode == kADCDynamicPedestal) {
- adc_pos_pedsub = hit->GetRawAdcHitPos().GetPulseInt();
- adc_neg_pedsub = hit->GetRawAdcHitNeg().GetPulseInt();
- adc_pos = hit->GetRawAdcHitPos().GetPulseIntRaw();
- adc_neg = hit->GetRawAdcHitNeg().GetPulseIntRaw();
- } else if (fADCMode == kADCSampleIntegral) {
- adc_pos_pedsub = hit->GetRawAdcHitPos().GetSampleIntRaw() - fPosPed[hit->fCounter -1];
- adc_neg_pedsub = hit->GetRawAdcHitNeg().GetSampleIntRaw() - fNegPed[hit->fCounter -1];
- adc_pos = hit->GetRawAdcHitPos().GetSampleIntRaw();
- adc_neg = hit->GetRawAdcHitNeg().GetSampleIntRaw();
- } else if (fADCMode == kADCSampIntDynPed) {
- adc_pos = hit->GetRawAdcHitPos().GetSampleInt();
- adc_neg = hit->GetRawAdcHitNeg().GetSampleInt();
- adc_pos_pedsub = hit->GetRawAdcHitPos().GetSampleIntRaw();
- adc_neg_pedsub = hit->GetRawAdcHitNeg().GetSampleIntRaw();
- } else {
- adc_pos_pedsub = hit->GetRawAdcHitPos().GetPulseIntRaw()-fPosPed[hit->fCounter -1];
- adc_neg_pedsub = hit->GetRawAdcHitNeg().GetPulseIntRaw()-fNegPed[hit->fCounter -1];
- adc_pos = hit->GetRawAdcHitPos().GetPulseIntRaw();
- adc_neg = hit->GetRawAdcHitNeg().GetPulseIntRaw();
- }
- //
- // Sparsify positive side hits, fill the hit list, compute the
- // energy depostion from positive side for the counter.
- fA_Pos[hit->fCounter-1] =adc_pos;
- fA_Neg[hit->fCounter-1] =adc_neg;
-
- Double_t thresh_pos = fPosThresh[hit->fCounter -1];
- thresh_pos =0;
- if(fA_Pos[hit->fCounter-1] > thresh_pos) {
-
- fA_Pos_p[hit->fCounter-1] =adc_pos_pedsub ;
-
- // cout << " pos " << hit->fCounter << " " << fA_Pos_p[hit->fCounter-1] << " " << fA_Pos[hit->fCounter-1] << " " << fPosPed[hit->fCounter -1] << "" << frNegAdcPulseAmp->ConstructedAt(nrNegAdcHits-1))->Get(padnum)<< endl;
- fEpos[hit->fCounter-1] += fA_Pos_p[hit->fCounter-1]*
- fParent->GetGain(hit->fCounter-1,fLayerNum-1,0);
- }
-
- // Sparsify negative side hits, fill the hit list, compute the
- // energy depostion from negative side for the counter.
-
- Double_t thresh_neg = fNegThresh[hit->fCounter -1];
- thresh_neg=0;
- if(fA_Neg[hit->fCounter-1] > thresh_neg) {
-
-
- fA_Neg_p[hit->fCounter-1] =adc_neg_pedsub ;
-
- fEneg[hit->fCounter-1] += fA_Neg_p[hit->fCounter-1]*
- fParent->GetGain(hit->fCounter-1,fLayerNum-1,1);
- }
-
- // Mean energy in the counter.
-
- fEmean[hit->fCounter-1] += (fEpos[hit->fCounter-1] + fEneg[hit->fCounter-1]);
-
- // Accumulate energies in the plane.
-
- fEplane += fEmean[hit->fCounter-1];
- fEplane_pos += fEpos[hit->fCounter-1];
- fEplane_neg += fEneg[hit->fCounter-1];
-
ihit++;
}
- //Debug output.
+ return(ihit);
+}
+//_____________________________________________________________________________
+Int_t THcShowerPlane::CoarseProcessHits()
+{
+ THcShower* fParent;
+ fParent = (THcShower*) GetParent();
+ Int_t ADCMode=fParent->GetADCMode();
+ if(ADCMode == kADCDynamicPedestal) {
+ FillADC_DynamicPedestal();
+ } else if (ADCMode == kADCSampleIntegral) {
+ FillADC_SampleIntegral();
+ } else if (ADCMode == kADCSampIntDynPed) {
+ FillADC_SampIntDynPed();
+ } else {
+ FillADC_Standard();
+ }
+ //
if (fParent->fdbg_decoded_cal) {
cout << "---------------------------------------------------------------\n";
cout << "Debug output from THcShowerPlane::ProcessHits for "
<< fParent->GetPrefix() << ":" << endl;
- cout << " nrawhits = " << nrawhits << " nexthit = " << nexthit << endl;
cout << " Sparsified hits for HMS calorimeter plane #" << fLayerNum
<< ", " << GetName() << ":" << endl;
Int_t nspar = 0;
- for (Int_t jhit = nexthit; jhit < nrawhits; jhit++) {
- THcRawShowerHit* hit = (THcRawShowerHit *) rawhits->At(jhit);
- if(hit->fPlane > fLayerNum) {
- break;
- }
+ for (UInt_t i=0; i<fParent->GetNBlocks(fLayerNum-1); i++) {
- if(fA_Pos[hit->fCounter-1] > fPosThresh[hit->fCounter -1] ||
- fA_Neg[hit->fCounter-1] > fNegThresh[hit->fCounter -1]) {
- cout << " plane = " << hit->fPlane
- << " counter = " << hit->fCounter
- << " Emean = " << fEmean[hit->fCounter-1]
- << " Epos = " << fEpos[hit->fCounter-1]
- << " Eneg = " << fEneg[hit->fCounter-1]
+ if (GetAposP(i) > 0 || GetAnegP(i) >0) { //hit
+ cout << " counter = " << i
+ << " Emean = " << fEmean[i]
+ << " Epos = " << fEpos[i]
+ << " Eneg = " << fEneg[i]
<< endl;
nspar++;
}
@@ -559,10 +522,103 @@ Int_t THcShowerPlane::ProcessHits(TClonesArray* rawhits, Int_t nexthit)
<< endl;
cout << "---------------------------------------------------------------\n";
}
-
- return(ihit);
+ //
+ return 1;
+ //
+}
+//_____________________________________________________________________________
+void THcShowerPlane::FillADC_SampIntDynPed()
+{
+ // adc_pos = hit->GetRawAdcHitPos().GetSampleInt();
+ // adc_neg = hit->GetRawAdcHitNeg().GetSampleInt();
+ // adc_pos_pedsub = hit->GetRawAdcHitPos().GetSampleIntRaw();
+ // adc_neg_pedsub = hit->GetRawAdcHitNeg().GetSampleIntRaw();
+ // Need to create
}
+//_____________________________________________________________________________
+void THcShowerPlane::FillADC_SampleIntegral()
+{
+ /// adc_pos_pedsub = hit->GetRawAdcHitPos().GetSampleIntRaw() - fPosPed[hit->fCounter -1];
+ // adc_neg_pedsub = hit->GetRawAdcHitNeg().GetSampleIntRaw() - fNegPed[hit->fCounter -1];
+ // adc_pos = hit->GetRawAdcHitPos().GetSampleIntRaw();
+ // adc_neg = hit->GetRawAdcHitNeg().GetSampleIntRaw();
+ // need to create
+}
+//_____________________________________________________________________________
+void THcShowerPlane::FillADC_Standard()
+{
+ THcShower* fParent;
+ fParent = (THcShower*) GetParent();
+ 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];
+ }
+ }
+ 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];
+ }
+ }
+ fEplane= fEplane_neg+fEplane_pos;
+}
+//_____________________________________________________________________________
+void THcShowerPlane::FillADC_DynamicPedestal()
+{
+ THcShower* fParent;
+ fParent = (THcShower*) GetParent();
+ 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];
+ }
+ }
+ }
+ //
+ 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];
+ }
+ }
+ }
+ //
+ fEplane= fEplane_neg+fEplane_pos;
+}
//_____________________________________________________________________________
Int_t THcShowerPlane::AccumulatePedestals(TClonesArray* rawhits, Int_t nexthit)
{
diff --git a/src/THcShowerPlane.h b/src/THcShowerPlane.h
index 8cb4346aa21ab3dd609d91ed79556275d0fb4332..fe00766d7e3c1fa6dfe81a52bf934a40b400f147 100644
--- a/src/THcShowerPlane.h
+++ b/src/THcShowerPlane.h
@@ -39,6 +39,7 @@ public:
virtual Bool_t IsPid() { return kFALSE; }
virtual Int_t ProcessHits(TClonesArray* rawhits, Int_t nexthit);
+ virtual Int_t CoarseProcessHits();
virtual Int_t AccumulatePedestals(TClonesArray* rawhits, Int_t nexthit);
virtual void CalculatePedestals( );
@@ -110,7 +111,6 @@ protected:
// Flash ADC parameters
Int_t fUsingFADC; // != 0 if using FADC in sample mode
- Int_t fADCMode; //
// 1 == Use the pulse int - pulse ped
// 2 == Use the sample integral - known ped
// 3 == Use the sample integral - sample ped
@@ -159,6 +159,7 @@ protected:
Float_t *fNegSig;
Float_t *fNegThresh;
+ TClonesArray* frPosAdcErrorFlag;
TClonesArray* frPosAdcPedRaw;
TClonesArray* frPosAdcPulseIntRaw;
TClonesArray* frPosAdcPulseAmpRaw;
@@ -168,6 +169,7 @@ protected:
TClonesArray* frPosAdcPulseInt;
TClonesArray* frPosAdcPulseAmp;
+ TClonesArray* frNegAdcErrorFlag;
TClonesArray* frNegAdcPedRaw;
TClonesArray* frNegAdcPulseIntRaw;
TClonesArray* frNegAdcPulseAmpRaw;
@@ -180,6 +182,10 @@ protected:
virtual Int_t ReadDatabase( const TDatime& date );
virtual Int_t DefineVariables( EMode mode = kDefine );
virtual void InitializePedestals( );
+ virtual void FillADC_DynamicPedestal( );
+ virtual void FillADC_SampleIntegral( );
+ virtual void FillADC_SampIntDynPed( );
+ virtual void FillADC_Standard( );
ClassDef(THcShowerPlane,0); // Calorimeter bars in a plane
};
#endif