From e7807c3546b98f774078c372f14113a5a8066342 Mon Sep 17 00:00:00 2001 From: Steven Doran <78985334+S81D@users.noreply.github.com> Date: Mon, 6 Jan 2025 13:17:12 -0600 Subject: [PATCH 1/5] Adjust EventSelector logic to handle MC waveforms Another downstream tool modification to adjust how the tool uses MC hits that are coming from the waveform generator. --- UserTools/EventSelector/EventSelector.cpp | 131 ++++++++++++++++++---- UserTools/EventSelector/EventSelector.h | 1 + 2 files changed, 108 insertions(+), 24 deletions(-) diff --git a/UserTools/EventSelector/EventSelector.cpp b/UserTools/EventSelector/EventSelector.cpp index 97efba7a6..d65196f72 100644 --- a/UserTools/EventSelector/EventSelector.cpp +++ b/UserTools/EventSelector/EventSelector.cpp @@ -14,6 +14,7 @@ bool EventSelector::Initialise(std::string configfile, DataModel &data){ fPMTMRDOffset = false; fIsMC = true; + fMCWaveform = false; fPMTMRDOffset = 755; fRecoPDG = -1; @@ -50,6 +51,7 @@ bool EventSelector::Initialise(std::string configfile, DataModel &data){ } m_variables.Get("SaveStatusToStore", fSaveStatusToStore); m_variables.Get("IsMC",fIsMC); + m_variables.Get("PMTWaveformSim",fMCWaveform); m_variables.Get("RecoPDG",fRecoPDG); m_variables.Get("TriggerExtendedWindow",fTriggerExtended); m_variables.Get("BeamOK",fBeamOK); @@ -632,8 +634,13 @@ bool EventSelector::EventSelectionByMCProjectedMRDHit() { bool EventSelector::EventSelectionByPMTMRDCoinc() { if (fIsMC){ - bool has_clustered_pmt = m_data->CStore.Get("ClusterMapMC",m_all_clusters_MC); - if (not has_clustered_pmt) { Log("EventSelector Tool: Error retrieving ClusterMapMC from CStore, did you run ClusterFinder beforehand?",v_error,verbosity); return false; } + if (fMCWaveform){ + bool has_clustered_pmt = m_data->CStore.Get("ClusterMap",m_all_clusters); + if (not has_clustered_pmt) { Log("EventSelector Tool: MCWaveform - Error retrieving ClusterMap from CStore, did you run ClusterFinder beforehand?",v_error,verbosity); return false; } + } else { + bool has_clustered_pmt = m_data->CStore.Get("ClusterMapMC",m_all_clusters_MC); + if (not has_clustered_pmt) { Log("EventSelector Tool: Error retrieving ClusterMapMC from CStore, did you run ClusterFinder beforehand?",v_error,verbosity); return false; } + } } else { bool has_clustered_pmt = m_data->CStore.Get("ClusterMap",m_all_clusters); if (not has_clustered_pmt) { Log("EventSelector Tool: Error retrieving ClusterMap from CStore, did you run ClusterFinder beforehand?",v_error,verbosity); return false; } @@ -649,8 +656,15 @@ bool EventSelector::EventSelectionByPMTMRDCoinc() { } int pmt_cluster_size; - if (fIsMC) pmt_cluster_size = (int) m_all_clusters_MC->size(); - else pmt_cluster_size = (int) m_all_clusters->size(); + if (fIsMC) { + if (fMCWaveform) { + pmt_cluster_size = (int) m_all_clusters->size(); + } else { + pmt_cluster_size = (int) m_all_clusters_MC->size(); + } + } else { + pmt_cluster_size = (int) m_all_clusters->size(); + } m_data->Stores["RecoEvent"]->Set("NumPMTClusters",pmt_cluster_size); vec_pmtclusters_charge->clear(); vec_pmtclusters_time->clear(); @@ -668,17 +682,27 @@ bool EventSelector::EventSelectionByPMTMRDCoinc() { pmt_time = -1; if (fIsMC){ - if (m_all_clusters_MC->size()){ + if (fMCWaveform) { + if (m_all_clusters->size()){ double cluster_time; - for(std::pair>&& apair : *m_all_clusters_MC){ - std::vector&MCHits = apair.second; + for(std::pair>&& apair : *m_all_clusters){ + std::vector&Hits = apair.second; double time_temp = 0; double charge_temp = 0; - for (unsigned int i_hit = 0; i_hit < MCHits.size(); i_hit++){ - time_temp+=MCHits.at(i_hit).GetTime(); - charge_temp+=MCHits.at(i_hit).GetCharge(); + for (unsigned int i_hit = 0; i_hit < Hits.size(); i_hit++){ + time_temp+=Hits.at(i_hit).GetTime(); + int tube = Hits.at(i_hit).GetTubeId(); + // check if PMT is present in the map before accessing it + auto it = ChannelNumToTankPMTSPEChargeMap->find(tube); + if (it != ChannelNumToTankPMTSPEChargeMap->end()) { + double charge_pe = Hits.at(i_hit).GetCharge() / it->second; + charge_temp += charge_pe; + } else { + std::cerr << "PMT channel with hit not found in ChannelNumToTankPMTSPEChargeMap. Skipping this hit." << std::endl; + continue; + } } - if (MCHits.size()>0) time_temp/=MCHits.size(); + if (Hits.size()>0) time_temp/=Hits.size(); vec_pmtclusters_charge->push_back(charge_temp); vec_pmtclusters_time->push_back(time_temp); if (time_temp > 2000.) continue; //not a prompt event @@ -686,7 +710,31 @@ bool EventSelector::EventSelectionByPMTMRDCoinc() { max_charge = charge_temp; prompt_cluster = true; pmt_time = time_temp; - n_hits = int(MCHits.size()); + n_hits = int(Hits.size()); + } + } + } + } else { + if (m_all_clusters_MC->size()){ + double cluster_time; + for(std::pair>&& apair : *m_all_clusters_MC){ + std::vector&MCHits = apair.second; + double time_temp = 0; + double charge_temp = 0; + for (unsigned int i_hit = 0; i_hit < MCHits.size(); i_hit++){ + time_temp+=MCHits.at(i_hit).GetTime(); + charge_temp+=MCHits.at(i_hit).GetCharge(); + } + if (MCHits.size()>0) time_temp/=MCHits.size(); + vec_pmtclusters_charge->push_back(charge_temp); + vec_pmtclusters_time->push_back(time_temp); + if (time_temp > 2000.) continue; //not a prompt event + if (charge_temp > max_charge){ + max_charge = charge_temp; + prompt_cluster = true; + pmt_time = time_temp; + n_hits = int(MCHits.size()); + } } } } @@ -762,10 +810,14 @@ bool EventSelector::EventSelectionByPMTMRDCoinc() { } m_data->Stores["RecoEvent"]->Set("MRDClustersTime",vec_mrdclusters_time, true); - if (fIsMC){ - if (MrdTimeClusters.size() == 0 || m_all_clusters_MC->size() == 0) return false; + if (fIsMC) { + if (MrdTimeClusters.size() == 0 || (fMCWaveform ? m_all_clusters->size() == 0 : m_all_clusters_MC->size() == 0)) { + return false; + } } else { - if (MrdTimeClusters.size() == 0 || m_all_clusters->size() == 0) return false; + if (MrdTimeClusters.size() == 0 || m_all_clusters->size() == 0) { + return false; + } } pmtmrd_coinc_min = fPMTMRDOffset - 50; @@ -1027,8 +1079,13 @@ bool EventSelector::FindPaddleChankey(double x, double y, int layer, unsigned lo bool EventSelector::EventSelectionByRecoPDG(int recoPDG, std::vector & cluster_reco_pdg){ if (fIsMC){ - bool has_clustered_pmt = m_data->CStore.Get("ClusterMapMC",m_all_clusters_MC); - if (not has_clustered_pmt) { Log("EventSelector Tool: Error retrieving ClusterMapMC from CStore, did you run ClusterFinder beforehand?",v_error,verbosity); return false; } + if (fMCWaveform) { + bool has_clustered_pmt = m_data->CStore.Get("ClusterMap",m_all_clusters); + if (not has_clustered_pmt) { Log("EventSelector Tool: MCWaveform --> Error retrieving ClusterMap from CStore, did you run ClusterFinder beforehand?",v_error,verbosity); return false; } + } else { + bool has_clustered_pmt = m_data->CStore.Get("ClusterMapMC",m_all_clusters_MC); + if (not has_clustered_pmt) { Log("EventSelector Tool: Error retrieving ClusterMapMC from CStore, did you run ClusterFinder beforehand?",v_error,verbosity); return false; } + } } else { bool has_clustered_pmt = m_data->CStore.Get("ClusterMap",m_all_clusters); if (not has_clustered_pmt) { Log("EventSelector Tool: Error retrieving ClusterMap from CStore, did you run ClusterFinder beforehand?",v_error,verbosity); return false; } @@ -1046,26 +1103,52 @@ bool EventSelector::EventSelectionByRecoPDG(int recoPDG, std::vector & c if (fabs(recoPDG)==2112){ if (fIsMC){ - if (m_all_clusters_MC->size()){ - for(std::pair>&& apair : *m_all_clusters_MC){ + if (fMCWaveform) { + if (m_all_clusters->size()){ + double cluster_time; + for(std::pair>&& apair : *m_all_clusters){ double cluster_time = apair.first; double charge_balance = ClusterChargeBalances.at(cluster_time); - std::vector&MCHits = apair.second; + std::vector&Hits = apair.second; double time_temp = 0; double charge_temp = 0; - for (unsigned int i_hit = 0; i_hit < MCHits.size(); i_hit++){ - time_temp+=MCHits.at(i_hit).GetTime(); - charge_temp+=MCHits.at(i_hit).GetCharge(); + for (unsigned int i_hit = 0; i_hit < Hits.size(); i_hit++){ + time_temp+=Hits.at(i_hit).GetTime(); + int tube = Hits.at(i_hit).GetTubeId(); + double charge_pe = Hits.at(i_hit).GetCharge()/ChannelNumToTankPMTSPEChargeMap->at(tube); + charge_temp+=charge_pe; } if (cluster_time > 10000 && charge_balance < 0.4 && charge_temp < 120) { cluster_reco_pdg.push_back(cluster_time); found_pdg = true; - std::cout <<"Found neutron candidate for cluster at time = "<size()){ + for(std::pair>&& apair : *m_all_clusters_MC){ + double cluster_time = apair.first; + double charge_balance = ClusterChargeBalances.at(cluster_time); + std::vector&MCHits = apair.second; + double time_temp = 0; + double charge_temp = 0; + for (unsigned int i_hit = 0; i_hit < MCHits.size(); i_hit++){ + time_temp+=MCHits.at(i_hit).GetTime(); + charge_temp+=MCHits.at(i_hit).GetCharge(); + } + if (cluster_time > 10000 && charge_balance < 0.4 && charge_temp < 120) { + cluster_reco_pdg.push_back(cluster_time); + found_pdg = true; + std::cout <<"Found neutron candidate for cluster at time = "<size()){ double cluster_time; diff --git a/UserTools/EventSelector/EventSelector.h b/UserTools/EventSelector/EventSelector.h index 82a18d3cf..d6d1f6455 100644 --- a/UserTools/EventSelector/EventSelector.h +++ b/UserTools/EventSelector/EventSelector.h @@ -249,6 +249,7 @@ class EventSelector: public Tool { bool fThroughGoing = false; bool fEventCutStatus; bool fIsMC; + bool fMCWaveform; int fTriggerWord; int fRecoPDG; bool fTriggerExtended = false; From a3c1d8c8f2c8d5b596f966ffa0d33ada4b582dcf Mon Sep 17 00:00:00 2001 From: Steven Doran <78985334+S81D@users.noreply.github.com> Date: Mon, 6 Jan 2025 13:18:57 -0600 Subject: [PATCH 2/5] Update README.md Added config variable that tells EventSelector to use MC waveform hits rather than the default MC hits. --- UserTools/EventSelector/README.md | 1 + 1 file changed, 1 insertion(+) diff --git a/UserTools/EventSelector/README.md b/UserTools/EventSelector/README.md index efbea06ad..69513e7cf 100644 --- a/UserTools/EventSelector/README.md +++ b/UserTools/EventSelector/README.md @@ -106,5 +106,6 @@ ThroughGoing TriggerWord RecoPDG IsMC +PMTWaveformSim SaveStatusToStore ``` From a7773faefc955943e59cb1f535c47cc29002ddbf Mon Sep 17 00:00:00 2001 From: Steven Doran <78985334+S81D@users.noreply.github.com> Date: Mon, 8 Sep 2025 14:54:06 -0500 Subject: [PATCH 3/5] Update EventSelector.cpp Added skip event for PMTWaveformSim --- UserTools/EventSelector/EventSelector.cpp | 9 +++++++++ 1 file changed, 9 insertions(+) diff --git a/UserTools/EventSelector/EventSelector.cpp b/UserTools/EventSelector/EventSelector.cpp index d65196f72..ba23b85aa 100644 --- a/UserTools/EventSelector/EventSelector.cpp +++ b/UserTools/EventSelector/EventSelector.cpp @@ -97,6 +97,15 @@ bool EventSelector::Execute(){ return false; }; + // an upstream tool may have told us to skip this event (PMTWaveformSim) + bool skip = false; + bool goodSkipStatus = m_data->Stores.at("ANNIEEvent")->Get("SkipExecute", skip); + if (goodSkipStatus && skip) { + logmessage = "EventSelector: An upstream tool told me to skip this event."; + Log(logmessage, v_warning, verbosity); + return true; + } + // ANNIE Event number m_data->Stores.at("ANNIEEvent")->Get("EventNumber",fEventNumber); From 3b1b36ab5665894e71e308b21d250e4b3f4db44e Mon Sep 17 00:00:00 2001 From: Steven Doran <78985334+S81D@users.noreply.github.com> Date: Wed, 22 Oct 2025 10:29:41 -0500 Subject: [PATCH 4/5] Update EventSelector.cpp Cleaned up logic --- UserTools/EventSelector/EventSelector.cpp | 127 ++++++++-------------- 1 file changed, 44 insertions(+), 83 deletions(-) diff --git a/UserTools/EventSelector/EventSelector.cpp b/UserTools/EventSelector/EventSelector.cpp index ba23b85aa..8b279fa4a 100644 --- a/UserTools/EventSelector/EventSelector.cpp +++ b/UserTools/EventSelector/EventSelector.cpp @@ -642,17 +642,12 @@ bool EventSelector::EventSelectionByMCProjectedMRDHit() { bool EventSelector::EventSelectionByPMTMRDCoinc() { - if (fIsMC){ - if (fMCWaveform){ - bool has_clustered_pmt = m_data->CStore.Get("ClusterMap",m_all_clusters); + if (!fIsMC || fMCWaveform) { + bool has_clustered_pmt = m_data->CStore.Get("ClusterMap",m_all_clusters); if (not has_clustered_pmt) { Log("EventSelector Tool: MCWaveform - Error retrieving ClusterMap from CStore, did you run ClusterFinder beforehand?",v_error,verbosity); return false; } - } else { - bool has_clustered_pmt = m_data->CStore.Get("ClusterMapMC",m_all_clusters_MC); - if (not has_clustered_pmt) { Log("EventSelector Tool: Error retrieving ClusterMapMC from CStore, did you run ClusterFinder beforehand?",v_error,verbosity); return false; } - } } else { - bool has_clustered_pmt = m_data->CStore.Get("ClusterMap",m_all_clusters); - if (not has_clustered_pmt) { Log("EventSelector Tool: Error retrieving ClusterMap from CStore, did you run ClusterFinder beforehand?",v_error,verbosity); return false; } + bool has_clustered_pmt = m_data->CStore.Get("ClusterMapMC",m_all_clusters_MC); + if (not has_clustered_pmt) { Log("EventSelector Tool: Error retrieving ClusterMapMC from CStore, did you run ClusterFinder beforehand?",v_error,verbosity); return false; } } bool has_clustered_mrd = m_data->CStore.Get("MrdTimeClusters",MrdTimeClusters); @@ -665,15 +660,12 @@ bool EventSelector::EventSelectionByPMTMRDCoinc() { } int pmt_cluster_size; - if (fIsMC) { - if (fMCWaveform) { - pmt_cluster_size = (int) m_all_clusters->size(); - } else { - pmt_cluster_size = (int) m_all_clusters_MC->size(); - } + if (!fIsMC || fMCWaveform) { + pmt_cluster_size = (int) m_all_clusters->size(); } else { - pmt_cluster_size = (int) m_all_clusters->size(); + pmt_cluster_size = (int) m_all_clusters_MC->size(); } + m_data->Stores["RecoEvent"]->Set("NumPMTClusters",pmt_cluster_size); vec_pmtclusters_charge->clear(); vec_pmtclusters_time->clear(); @@ -690,84 +682,53 @@ bool EventSelector::EventSelectionByPMTMRDCoinc() { pmt_time = -1; - if (fIsMC){ - if (fMCWaveform) { - if (m_all_clusters->size()){ - double cluster_time; - for(std::pair>&& apair : *m_all_clusters){ - std::vector&Hits = apair.second; - double time_temp = 0; - double charge_temp = 0; - for (unsigned int i_hit = 0; i_hit < Hits.size(); i_hit++){ - time_temp+=Hits.at(i_hit).GetTime(); - int tube = Hits.at(i_hit).GetTubeId(); - // check if PMT is present in the map before accessing it + // MC Waveform or Data + if (!fIsMC || fMCWaveform) { + if (m_all_clusters->size()){ + double cluster_time; + for(std::pair>&& apair : *m_all_clusters){ + std::vector&Hits = apair.second; + double time_temp = 0; + double charge_temp = 0; + for (unsigned int i_hit = 0; i_hit < Hits.size(); i_hit++){ + time_temp+=Hits.at(i_hit).GetTime(); + int tube = Hits.at(i_hit).GetTubeId(); + // check if PMT is present in the map before accessing it auto it = ChannelNumToTankPMTSPEChargeMap->find(tube); if (it != ChannelNumToTankPMTSPEChargeMap->end()) { - double charge_pe = Hits.at(i_hit).GetCharge() / it->second; - charge_temp += charge_pe; + double charge_pe = Hits.at(i_hit).GetCharge() / it->second; + charge_temp += charge_pe; } else { - std::cerr << "PMT channel with hit not found in ChannelNumToTankPMTSPEChargeMap. Skipping this hit." << std::endl; - continue; + std::cerr << "PMT channel with hit not found in ChannelNumToTankPMTSPEChargeMap. Skipping this hit." << std::endl; + continue; } - } - if (Hits.size()>0) time_temp/=Hits.size(); - vec_pmtclusters_charge->push_back(charge_temp); - vec_pmtclusters_time->push_back(time_temp); - if (time_temp > 2000.) continue; //not a prompt event - if (charge_temp > max_charge){ - max_charge = charge_temp; - prompt_cluster = true; - pmt_time = time_temp; - n_hits = int(Hits.size()); - } } - } - } else { - if (m_all_clusters_MC->size()){ - double cluster_time; - for(std::pair>&& apair : *m_all_clusters_MC){ - std::vector&MCHits = apair.second; - double time_temp = 0; - double charge_temp = 0; - for (unsigned int i_hit = 0; i_hit < MCHits.size(); i_hit++){ - time_temp+=MCHits.at(i_hit).GetTime(); - charge_temp+=MCHits.at(i_hit).GetCharge(); - } - if (MCHits.size()>0) time_temp/=MCHits.size(); - vec_pmtclusters_charge->push_back(charge_temp); - vec_pmtclusters_time->push_back(time_temp); - if (time_temp > 2000.) continue; //not a prompt event - if (charge_temp > max_charge){ - max_charge = charge_temp; - prompt_cluster = true; - pmt_time = time_temp; - n_hits = int(MCHits.size()); - } - } + if (Hits.size()>0) time_temp/=Hits.size(); + vec_pmtclusters_charge->push_back(charge_temp); + vec_pmtclusters_time->push_back(time_temp); + if (time_temp > 2000.) continue; //not a prompt event + if (charge_temp > max_charge){ + max_charge = charge_temp; + prompt_cluster = true; + pmt_time = time_temp; + n_hits = int(Hits.size()); } } + } + + // MC (parametric) } else { - if (m_all_clusters->size()){ + if (m_all_clusters_MC->size()){ double cluster_time; - for(std::pair>&& apair : *m_all_clusters){ - std::vector&Hits = apair.second; + for(std::pair>&& apair : *m_all_clusters_MC){ + std::vector&MCHits = apair.second; double time_temp = 0; double charge_temp = 0; - for (unsigned int i_hit = 0; i_hit < Hits.size(); i_hit++){ - time_temp+=Hits.at(i_hit).GetTime(); - int tube = Hits.at(i_hit).GetTubeId(); - // check if PMT is present in the map before accessing it - auto it = ChannelNumToTankPMTSPEChargeMap->find(tube); - if (it != ChannelNumToTankPMTSPEChargeMap->end()) { - double charge_pe = Hits.at(i_hit).GetCharge() / it->second; - charge_temp += charge_pe; - } else { - std::cerr << "PMT channel with hit not found in ChannelNumToTankPMTSPEChargeMap. Skipping this hit." << std::endl; - continue; - } + for (unsigned int i_hit = 0; i_hit < MCHits.size(); i_hit++){ + time_temp+=MCHits.at(i_hit).GetTime(); + charge_temp+=MCHits.at(i_hit).GetCharge(); } - if (Hits.size()>0) time_temp/=Hits.size(); + if (MCHits.size()>0) time_temp/=MCHits.size(); vec_pmtclusters_charge->push_back(charge_temp); vec_pmtclusters_time->push_back(time_temp); if (time_temp > 2000.) continue; //not a prompt event @@ -775,7 +736,7 @@ bool EventSelector::EventSelectionByPMTMRDCoinc() { max_charge = charge_temp; prompt_cluster = true; pmt_time = time_temp; - n_hits = int(Hits.size()); + n_hits = int(MCHits.size()); } } } From eab30b8980adccbbda4002c419adb9913e4c6031 Mon Sep 17 00:00:00 2001 From: Steven Doran <78985334+S81D@users.noreply.github.com> Date: Tue, 4 Nov 2025 09:42:09 -0600 Subject: [PATCH 5/5] Update EventSelector.cpp - code condensed by merging data and MC logic together - removed the event skipping check --- UserTools/EventSelector/EventSelector.cpp | 103 +++++++--------------- 1 file changed, 31 insertions(+), 72 deletions(-) diff --git a/UserTools/EventSelector/EventSelector.cpp b/UserTools/EventSelector/EventSelector.cpp index 8b279fa4a..e06fb99c5 100644 --- a/UserTools/EventSelector/EventSelector.cpp +++ b/UserTools/EventSelector/EventSelector.cpp @@ -96,15 +96,6 @@ bool EventSelector::Execute(){ Log("EventSelector Tool: No ANNIEEvent store!",v_error,verbosity); return false; }; - - // an upstream tool may have told us to skip this event (PMTWaveformSim) - bool skip = false; - bool goodSkipStatus = m_data->Stores.at("ANNIEEvent")->Get("SkipExecute", skip); - if (goodSkipStatus && skip) { - logmessage = "EventSelector: An upstream tool told me to skip this event."; - Log(logmessage, v_warning, verbosity); - return true; - } // ANNIE Event number m_data->Stores.at("ANNIEEvent")->Get("EventNumber",fEventNumber); @@ -1048,17 +1039,12 @@ bool EventSelector::FindPaddleChankey(double x, double y, int layer, unsigned lo bool EventSelector::EventSelectionByRecoPDG(int recoPDG, std::vector & cluster_reco_pdg){ - if (fIsMC){ - if (fMCWaveform) { - bool has_clustered_pmt = m_data->CStore.Get("ClusterMap",m_all_clusters); - if (not has_clustered_pmt) { Log("EventSelector Tool: MCWaveform --> Error retrieving ClusterMap from CStore, did you run ClusterFinder beforehand?",v_error,verbosity); return false; } - } else { - bool has_clustered_pmt = m_data->CStore.Get("ClusterMapMC",m_all_clusters_MC); - if (not has_clustered_pmt) { Log("EventSelector Tool: Error retrieving ClusterMapMC from CStore, did you run ClusterFinder beforehand?",v_error,verbosity); return false; } - } + if (!fIsMC || fMCWaveform) { + bool has_clustered_pmt = m_data->CStore.Get("ClusterMap",m_all_clusters); + if (not has_clustered_pmt) { Log("EventSelector Tool: Error retrieving ClusterMap from CStore, did you run ClusterFinder beforehand?",v_error,verbosity); return false; } } else { - bool has_clustered_pmt = m_data->CStore.Get("ClusterMap",m_all_clusters); - if (not has_clustered_pmt) { Log("EventSelector Tool: Error retrieving ClusterMap from CStore, did you run ClusterFinder beforehand?",v_error,verbosity); return false; } + bool has_clustered_pmt = m_data->CStore.Get("ClusterMapMC",m_all_clusters_MC); + if (not has_clustered_pmt) { Log("EventSelector Tool: Error retrieving ClusterMapMC from CStore, did you run ClusterFinder beforehand?",v_error,verbosity); return false; } } std::map ClusterMaxPEs; @@ -1072,33 +1058,31 @@ bool EventSelector::EventSelectionByRecoPDG(int recoPDG, std::vector & c bool found_pdg = false; if (fabs(recoPDG)==2112){ - if (fIsMC){ - if (fMCWaveform) { - if (m_all_clusters->size()){ - double cluster_time; - for(std::pair>&& apair : *m_all_clusters){ - double cluster_time = apair.first; - double charge_balance = ClusterChargeBalances.at(cluster_time); - std::vector&Hits = apair.second; - double time_temp = 0; - double charge_temp = 0; - for (unsigned int i_hit = 0; i_hit < Hits.size(); i_hit++){ - time_temp+=Hits.at(i_hit).GetTime(); - int tube = Hits.at(i_hit).GetTubeId(); - double charge_pe = Hits.at(i_hit).GetCharge()/ChannelNumToTankPMTSPEChargeMap->at(tube); - charge_temp+=charge_pe; - } - if (cluster_time > 10000 && charge_balance < 0.4 && charge_temp < 120) { - cluster_reco_pdg.push_back(cluster_time); - found_pdg = true; - std::cout <<"Found neutron candidate for cluster at time = "<size()){ + if (!fIsMC || fMCWaveform) { // data-like + if (m_all_clusters->size()){ + for(std::pair>&& apair : *m_all_clusters){ + double cluster_time = apair.first; + double charge_balance = ClusterChargeBalances.at(cluster_time); + std::vector&Hits = apair.second; + double time_temp = 0; + double charge_temp = 0; + for (unsigned int i_hit = 0; i_hit < Hits.size(); i_hit++){ + time_temp+=Hits.at(i_hit).GetTime(); + int tube = Hits.at(i_hit).GetTubeId(); + double charge_pe = Hits.at(i_hit).GetCharge()/ChannelNumToTankPMTSPEChargeMap->at(tube); + charge_temp+=charge_pe; + } + if (cluster_time > 10000 && charge_balance < 0.4 && charge_temp < 120) { + cluster_reco_pdg.push_back(cluster_time); + found_pdg = true; + std::cout <<"Found neutron candidate for cluster at time = "<size()){ for(std::pair>&& apair : *m_all_clusters_MC){ double cluster_time = apair.first; double charge_balance = ClusterChargeBalances.at(cluster_time); @@ -1118,32 +1102,7 @@ bool EventSelector::EventSelectionByRecoPDG(int recoPDG, std::vector & c } } } - } - } else { - if (m_all_clusters->size()){ - double cluster_time; - for(std::pair>&& apair : *m_all_clusters){ - double cluster_time = apair.first; - double charge_balance = ClusterChargeBalances.at(cluster_time); - std::vector&Hits = apair.second; - double time_temp = 0; - double charge_temp = 0; - for (unsigned int i_hit = 0; i_hit < Hits.size(); i_hit++){ - time_temp+=Hits.at(i_hit).GetTime(); - int tube = Hits.at(i_hit).GetTubeId(); - double charge_pe = Hits.at(i_hit).GetCharge()/ChannelNumToTankPMTSPEChargeMap->at(tube); - charge_temp+=charge_pe; - } - if (cluster_time > 10000 && charge_balance < 0.4 && charge_temp < 120) { - cluster_reco_pdg.push_back(cluster_time); - found_pdg = true; - std::cout <<"Found neutron candidate for cluster at time = "<