Belle II Software development
EventT0Validation.cc
1/**************************************************************************
2 * basf2 (Belle II Analysis Software Framework) *
3 * Author: The Belle II Collaboration *
4 * *
5 * See git log for contributors and copyright holders. *
6 * This file is licensed under LGPL-3.0, see LICENSE.md. *
7 **************************************************************************/
8
9#include <reconstruction/modules/EventT0Validation/EventT0Validation.h>
10
11#include <framework/dataobjects/EventT0.h>
12
13#include <TEfficiency.h>
14
15using namespace Belle2;
16
17REG_MODULE(EventT0Validation);
18
19//---------------------------------
21{
22 setPropertyFlags(c_ParallelProcessingCertified); // parallel processing
23 setDescription("Make data quality monitoring plots for EventT0 for bhabha, mu mu, and hadron samples for different trigger (time) sources.");
24 addParam("RootFileName", m_RootFileName, "Name of the ROOT output file.", m_RootFileName);
25}
26
27//---------------------------------
29
30
31//---------------------------------
33{
34 m_outputFile = new TFile(m_RootFileName.c_str(), "RECREATE");
35 if (m_outputFile != nullptr) {
36 m_outputFile->cd();
37
38 int nBins = 400 ;
39 double minT0 = -100 ;
40 double maxT0 = 100 ;
41
42 m_histECLEventT0 = new TH1F("histECLEventT0", "ECL EventT0;EventT0 [ns];events / 0.5 ns", nBins, minT0, maxT0);
43 setPlotMetaData(m_histECLEventT0, "Distribution of ECL EventT0s.", "Should be centered around 0.", m_contact);
44 m_histSVDEventT0 = new TH1F("histSVDEventT0", "SVD EventT0;EventT0 [ns];events / 0.5 ns", nBins, minT0, maxT0);
45 setPlotMetaData(m_histSVDEventT0, "Distribution of SVD EventT0s.", "Should be centered around 0.", m_contact);
46 m_histTOPEventT0 = new TH1F("histTOPEventT0", "TOP EventT0;EventT0 [ns];events / 0.5 ns", nBins, minT0, maxT0);
47 setPlotMetaData(m_histTOPEventT0, "Distribution of TOP EventT0s.", "Should be centered around 0.", m_contact);
48 m_histCDCEventT0 = new TH1F("histCDCEventT0", "CDC EventT0;EventT0 [ns];events / 0.5 ns", nBins, minT0, maxT0);
49 setPlotMetaData(m_histCDCEventT0, "Distribution of CDC EventT0s.", "Should be centered around 0.", m_contact);
50 m_histCDCHitBasedEventT0 = new TH1F("histCDCHitBasedEventT0", "CDC hit based EventT0;EventT0 [ns];events / 0.5 ns",
51 nBins, minT0, maxT0);
52 setPlotMetaData(m_histCDCHitBasedEventT0, "Distribution of CDC hit based EventT0s.", "Should be centered around 0.", m_contact);
53 m_histCDCChi2EventT0 = new TH1F("histCDCChi2EventT0", "CDC FullGrid #chi^{2} EventT0;EventT0 [ns];events / 0.5 ns",
54 nBins, minT0, maxT0);
55 setPlotMetaData(m_histCDCChi2EventT0, "Distribution of CDC FullGrid #chi^{2} EventT0s.", "Should be centered around 0.",
56 m_contact, "");
57 m_histCDCGridEventT0 = new TH1F("histCDCGridEventT0", "CDC Grid search EventT0;EventT0 [ns];events / 0.5 ns", nBins, minT0,
58 maxT0);
59 setPlotMetaData(m_histCDCGridEventT0, "Distribution of CDC Grid search EventT0s.", "Should be centered around 0.", m_contact, "");
60
62 new TH1D("histAlgorithmSourceCounts",
63 "Number of events with EventT0 from each algorithm;Algorithm;Count",
64 11, 0, 11);
66 "Number of events in which an EventT0 was found by each algorithm. "\
67 "Some of the CDC algorithms are only executed if no SVD EventT0 is found.",
68 "Values should be around 1000 (1000 events in validation).", m_contact);
70 new TH1D("histAlgorithmSourceCountsActive",
71 "Number of events with EventT0 from each algorithm where it was active;Algorithm;Count",
72 11, 0, 11);
74 "Number of events in which an EventT0 was found by each algorithm if they are executed. "\
75 "Some of the CDC algorithms are only executed if no SVD EventT0 is found.",
76 "Values should be around 1000 (1000 events in validation).", m_contact);
77
78 for (uint i = 0; i < 11; i++) {
79 m_histAlgorithmSourceCounts->GetXaxis()->SetBinLabel(i + 1, c_eventT0Algorithms[i]);
80 m_histAlgorithmSourceCountsActive->GetXaxis()->SetBinLabel(i + 1, c_eventT0Algorithms[i]);
81 }
82 m_histAlgorithmSourceCounts->GetXaxis()->CenterTitle(kTRUE);
83 m_histAlgorithmSourceCounts->GetXaxis()->SetTitleOffset(1.2);
84 m_histAlgorithmSourceCountsActive->GetXaxis()->CenterTitle(kTRUE);
85 m_histAlgorithmSourceCountsActive->GetXaxis()->SetTitleOffset(1.2);
86 }
87
88 m_eventT0.isRequired();
89}
90
91
92
93//---------------------------------
95{
96 if (!m_eventT0.isValid()) {
97 B2WARNING("Missing EventT0, EventT0Validation is skipped.");
98 return;
99 }
100}
101
102
103//---------------------------------
105{
106 // Determine if there is a valid event t0 to use and then extract the information about it
107 if (!m_eventT0.isValid()) {
108 B2WARNING("Missing EventT0, EventT0Validation is skipped.");
109 return ;
110 }
111
112 // Set the different EventT0 values, default is -1000 in case there are no information based on a given detector
113 const double eventT0ECL =
114 m_eventT0->hasTemporaryEventT0(Const::EDetector::ECL) ? m_eventT0->getBestECLTemporaryEventT0()->eventT0 : -1000;
115 const double eventT0CDC =
116 m_eventT0->hasTemporaryEventT0(Const::EDetector::CDC) ? m_eventT0->getBestCDCTemporaryEventT0()->eventT0 : -1000;
117 const double eventT0TOP =
118 m_eventT0->hasTemporaryEventT0(Const::EDetector::TOP) ? m_eventT0->getBestTOPTemporaryEventT0()->eventT0 : -1000;
119 const double eventT0SVD =
120 m_eventT0->hasTemporaryEventT0(Const::EDetector::SVD) ? m_eventT0->getBestSVDTemporaryEventT0()->eventT0 : -1000;
121
122 const auto getCDCEventT0sForAlgorithm = [cdcEventT0s = m_eventT0->getTemporaryEventT0s(Const::EDetector::CDC)](
123 const std::string & algorithm) {
124 std::vector<EventT0::EventT0Component> temporaries;
125 temporaries.reserve(cdcEventT0s.size());
126 for (const auto& evtt0 : cdcEventT0s) {
127 if (evtt0.algorithm == algorithm) {
128 temporaries.push_back(evtt0);
129 }
130 }
131 return temporaries;
132 };
133
134 m_histECLEventT0->Fill(eventT0ECL);
135 m_histSVDEventT0->Fill(eventT0SVD);
136 m_histTOPEventT0->Fill(eventT0TOP);
137 m_histCDCEventT0->Fill(eventT0CDC);
138
139 const auto hitBasedCDCT0 = getCDCEventT0sForAlgorithm("hit based");
140 const auto chi2CDCT0 = getCDCEventT0sForAlgorithm("chi2");
141 const auto gridCDCT0 = getCDCEventT0sForAlgorithm("grid");
142 m_histCDCHitBasedEventT0->Fill(hitBasedCDCT0.empty() ? -1000 : hitBasedCDCT0.back().eventT0);
143 m_histCDCChi2EventT0->Fill(chi2CDCT0.empty() ? -1000 : chi2CDCT0.back().eventT0);
144 m_histCDCGridEventT0->Fill(gridCDCT0.empty() ? -1000 : gridCDCT0.back().eventT0);
145
146 B2DEBUG(20, "eventT0ECL = " << eventT0ECL << " ns") ;
147 B2DEBUG(20, "eventT0CDC = " << eventT0CDC << " ns") ;
148 B2DEBUG(20, "eventT0TOP = " << eventT0TOP << " ns") ;
149 B2DEBUG(20, "eventT0SVD = " << eventT0SVD << " ns") ;
150
151 const bool hasECLEventT0 = m_eventT0->hasTemporaryEventT0(Const::EDetector::ECL);
152 const bool hasSVDEventT0 = m_eventT0->hasTemporaryEventT0(Const::EDetector::SVD);
153 const bool hasTOPEventT0 = m_eventT0->hasTemporaryEventT0(Const::EDetector::TOP);
154 const bool hasCDCEventT0 = m_eventT0->hasTemporaryEventT0(Const::EDetector::CDC);
155 const bool hasCDCHitBasedEventT0 = not hitBasedCDCT0.empty();
156 const bool hasCDCFullGridChi2EventT0 = not chi2CDCT0.empty();
157 const bool hasCDCGridEventT0 = not gridCDCT0.empty();
158
161 m_histAlgorithmSourceCounts->Fill(c_eventT0Algorithms[1], hasECLEventT0);
162 m_histAlgorithmSourceCounts->Fill(c_eventT0Algorithms[2], hasSVDEventT0);
163 m_histAlgorithmSourceCounts->Fill(c_eventT0Algorithms[3], hasTOPEventT0);
164 m_histAlgorithmSourceCounts->Fill(c_eventT0Algorithms[4], hasCDCEventT0);
165 m_histAlgorithmSourceCounts->Fill(c_eventT0Algorithms[5], hasCDCHitBasedEventT0);
166 m_histAlgorithmSourceCounts->Fill(c_eventT0Algorithms[6], hasCDCHitBasedEventT0);
167 m_histAlgorithmSourceCounts->Fill(c_eventT0Algorithms[7], hasCDCFullGridChi2EventT0);
168 m_histAlgorithmSourceCounts->Fill(c_eventT0Algorithms[8], hasCDCFullGridChi2EventT0);
169 m_histAlgorithmSourceCounts->Fill(c_eventT0Algorithms[9], hasCDCGridEventT0);
170 m_histAlgorithmSourceCounts->Fill(c_eventT0Algorithms[10], hasCDCGridEventT0);
171
178 m_histAlgorithmSourceCountsActive->Fill(c_eventT0Algorithms[5], 1); // We always execute hit based search
179 m_histAlgorithmSourceCountsActive->Fill(c_eventT0Algorithms[6], 1); // We always execute hit based search
180 // We only execute the chi2 algorithm if no SVD value is found, but this is the "all" column
182 // We only execute the chi2 algorithm if no SVD value is found
183 m_histAlgorithmSourceCountsActive->Fill(c_eventT0Algorithms[8], not hasSVDEventT0);
184 // We only execute the chi2 algorithm if no SVD value is found, but this is the "all" column
186 // We only execute the grid algorithm if no SVD value is found
187 m_histAlgorithmSourceCountsActive->Fill(c_eventT0Algorithms[10], not hasSVDEventT0);
188}
189
190
192{
193 if (m_outputFile != nullptr) {
194 m_outputFile->cd();
195
196 TEfficiency* algorithmEffi = new TEfficiency(*m_histAlgorithmSourceCounts, *m_histAlgorithmSourceCountsActive);
197 algorithmEffi->SetTitle("Efficiency of finding an EventT0 per Algorithm");
198 algorithmEffi->GetListOfFunctions()->Add(new TNamed("Contact", m_contact.c_str()));
199 algorithmEffi->GetListOfFunctions()->Add(new TNamed("Description", "Efficiencies of the various EventT0 algorithms."));
200 algorithmEffi->GetListOfFunctions()->Add(new TNamed("Check", "Efficiency should be close to 1 for all active algorithms."));
201 algorithmEffi->GetListOfFunctions()->Add(new TNamed("MetaOptions", "shifter"));
202 algorithmEffi->Write("EventT0AlgorithmEfficiency");
203
204 m_outputFile->Write();
205 m_outputFile->Close();
206 }
207}
208
209
210void EventT0ValidationModule::setPlotMetaData(TH1* hist, const std::string& description, const std::string& check,
211 const std::string& contact, const std::string& shifter)
212{
213 hist->GetListOfFunctions()->Add(new TNamed("Contact", contact.c_str()));
214 hist->GetListOfFunctions()->Add(new TNamed("Description", description.c_str()));
215 hist->GetListOfFunctions()->Add(new TNamed("Check", check.c_str()));
216 hist->GetListOfFunctions()->Add(new TNamed("MetaOptions", shifter.c_str()));
217}
TH1F * m_histTOPEventT0
TOP based EventT0 histograms.
TFile * m_outputFile
Output ROOT file.
virtual ~EventT0ValidationModule()
Destructor.
StoreObjPtr< EventT0 > m_eventT0
Store array for event t0.
TH1F * m_histCDCEventT0
CDC based EventT0 histograms.
std::string m_contact
For simplicity, just set the contact once.
virtual void initialize() override
Initialize the module.
TH1D * m_histAlgorithmSourceCounts
Count of events with EventT0 from a given algorithm, numerator for efficiency calculation.
virtual void event() override
This method is called for each event.
TH1F * m_histCDCHitBasedEventT0
CDC hit based EventT0 histograms.
TH1F * m_histCDCChi2EventT0
CDC chi2 based EventT0 histograms.
virtual void endRun() override
This method is called for each run.
TH1F * m_histCDCGridEventT0
CDC grid EventT0 histograms.
std::string m_RootFileName
Name of the output ROOT file.
TH1F * m_histECLEventT0
ECL based EventT0 histograms.
virtual void beginRun() override
This method is called for each run.
const char * c_eventT0Algorithms[11]
EventT0 algorithms for which to calculate fractions of abundance.
EventT0ValidationModule()
Default constructor.
TH1D * m_histAlgorithmSourceCountsActive
Counts of events with EventT0 from a given algorithm, denominator for efficiency calculation.
void setPlotMetaData(TH1 *hist, const std::string &description, const std::string &check, const std::string &contact, const std::string &shifter="shifter")
Set additional plot info: contact, description, shifter.
TH1F * m_histSVDEventT0
SVD based EventT0 histograms.
void setDescription(const std::string &description)
Sets the description of the module.
Definition Module.cc:214
void setPropertyFlags(unsigned int propertyFlags)
Sets the flags for the module properties.
Definition Module.cc:208
Module()
Constructor.
Definition Module.cc:30
@ c_ParallelProcessingCertified
This module can be run in parallel processing mode safely (All I/O must be done through the data stor...
Definition Module.h:80
void addParam(const std::string &name, T &paramVariable, const std::string &description, const T &defaultValue)
Adds a new parameter to the module.
Definition Module.h:559
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition Module.h:649
Abstract base class for different kinds of events.