Belle II Software development
EventT0DQM.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/EventT0DQM/EventT0DQM.h>
10#include <framework/core/HistoModule.h>
11
12using namespace Belle2;
13
14REG_MODULE(EventT0DQM);
15
16//---------------------------------
18{
19 setPropertyFlags(c_ParallelProcessingCertified); // parallel processing
20 setDescription("Make data quality monitoring plots for EventT0 for bhabha, mu mu, and hadron samples for different trigger (time) sources.");
21}
22
23//---------------------------------
25
26
27//---------------------------------
29{
30
31 TDirectory* oldDir = gDirectory;
32 oldDir->mkdir("EventT0")->cd();
33
34 int nBins = 400 ;
35 double minT0 = -100 ;
36 double maxT0 = 100 ;
37
38 m_histEventT0_ECL_bhabha_L1_ECLTRG = new TH1F("m_histEventT0_ECL_bhabha_L1_ECLTRG",
39 "ECL EventT0, L1TRG from ECL, HLT Bhabha;EventT0 [ns];events / 0.5 ns",
40 nBins, minT0, maxT0);
41 m_histEventT0_CDC_bhabha_L1_ECLTRG = new TH1F("m_histEventT0_CDC_bhabha_L1_ECLTRG",
42 "CDC EventT0, L1TRG from ECL, HLT Bhabha;EventT0 [ns];events / 0.5 ns",
43 nBins, minT0, maxT0);
44 m_histEventT0_TOP_bhabha_L1_ECLTRG = new TH1F("m_histEventT0_TOP_bhabha_L1_ECLTRG",
45 "TOP EventT0, L1TRG from ECL, HLT Bhabha;EventT0 [ns];events / 0.5 ns",
46 nBins, minT0, maxT0);
47 m_histEventT0_SVD_bhabha_L1_ECLTRG = new TH1F("m_histEventT0_SVD_bhabha_L1_ECLTRG",
48 "SVD EventT0, L1TRG from ECL, HLT Bhabha;EventT0 [ns];events / 0.5 ns",
49 nBins, minT0, maxT0);
50
51 m_histEventT0_ECL_hadron_L1_ECLTRG = new TH1F("m_histEventT0_ECL_hadron_L1_ECLTRG",
52 "ECL EventT0, L1TRG from ECL, HLT Hadron;EventT0 [ns];events / 0.5 ns",
53 nBins, minT0, maxT0);
54 m_histEventT0_CDC_hadron_L1_ECLTRG = new TH1F("m_histEventT0_CDC_hadron_L1_ECLTRG",
55 "CDC EventT0, L1TRG from ECL, HLT Hadron;EventT0 [ns];events / 0.5 ns",
56 nBins, minT0, maxT0);
57 m_histEventT0_TOP_hadron_L1_ECLTRG = new TH1F("m_histEventT0_TOP_hadron_L1_ECLTRG",
58 "TOP EventT0, L1TRG from ECL, HLT Hadron;EventT0 [ns];events / 0.5 ns",
59 nBins, minT0, maxT0);
60 m_histEventT0_SVD_hadron_L1_ECLTRG = new TH1F("m_histEventT0_SVD_hadron_L1_ECLTRG",
61 "SVD EventT0, L1TRG from ECL, HLT Hadron;EventT0 [ns];events / 0.5 ns",
62 nBins, minT0, maxT0);
63
64 m_histEventT0_ECL_mumu_L1_ECLTRG = new TH1F("m_histEventT0_ECL_mumu_L1_ECLTRG",
65 "ECL EventT0, L1TRG from ECL, HLT mumu;EventT0 [ns];events / 0.5 ns",
66 nBins, minT0, maxT0);
67 m_histEventT0_CDC_mumu_L1_ECLTRG = new TH1F("m_histEventT0_CDC_mumu_L1_ECLTRG",
68 "CDC EventT0, L1TRG from ECL, HLT mumu;EventT0 [ns];events / 0.5 ns",
69 nBins, minT0, maxT0);
70 m_histEventT0_TOP_mumu_L1_ECLTRG = new TH1F("m_histEventT0_TOP_mumu_L1_ECLTRG",
71 "TOP EventT0, L1TRG from ECL, HLT mumu;EventT0 [ns];events / 0.5 ns",
72 nBins, minT0, maxT0);
73 m_histEventT0_SVD_mumu_L1_ECLTRG = new TH1F("m_histEventT0_SVD_mumu_L1_ECLTRG",
74 "SVD EventT0, L1TRG from ECL, HLT mumu;EventT0 [ns];events / 0.5 ns",
75 nBins, minT0, maxT0);
76
77
78 m_histEventT0_ECL_bhabha_L1_CDCTRG = new TH1F("m_histEventT0_ECL_bhabha_L1_CDCTRG",
79 "ECL EventT0, L1TRG from CDC, HLT Bhabha;EventT0 [ns];events / 0.5 ns",
80 nBins, minT0, maxT0);
81 m_histEventT0_CDC_bhabha_L1_CDCTRG = new TH1F("m_histEventT0_CDC_bhabha_L1_CDCTRG",
82 "CDC EventT0, L1TRG from CDC, HLT Bhabha;EventT0 [ns];events / 0.5 ns",
83 nBins, minT0, maxT0);
84 m_histEventT0_TOP_bhabha_L1_CDCTRG = new TH1F("m_histEventT0_TOP_bhabha_L1_CDCTRG",
85 "TOP EventT0, L1TRG from CDC, HLT Bhabha;EventT0 [ns];events / 0.5 ns",
86 nBins, minT0, maxT0);
87 m_histEventT0_SVD_bhabha_L1_CDCTRG = new TH1F("m_histEventT0_SVD_bhabha_L1_CDCTRG",
88 "SVD EventT0, L1TRG from CDC, HLT Bhabha;EventT0 [ns];events / 0.5 ns",
89 nBins, minT0, maxT0);
90
91 m_histEventT0_ECL_hadron_L1_CDCTRG = new TH1F("m_histEventT0_ECL_hadron_L1_CDCTRG",
92 "ECL EventT0, L1TRG from CDC, HLT Hadron;EventT0 [ns];events / 0.5 ns",
93 nBins, minT0, maxT0);
94 m_histEventT0_CDC_hadron_L1_CDCTRG = new TH1F("m_histEventT0_CDC_hadron_L1_CDCTRG",
95 "CDC EventT0, L1TRG from CDC, HLT Hadron;EventT0 [ns];events / 0.5 ns",
96 nBins, minT0, maxT0);
97 m_histEventT0_TOP_hadron_L1_CDCTRG = new TH1F("m_histEventT0_TOP_hadron_L1_CDCTRG",
98 "TOP EventT0, L1TRG from CDC, HLT Hadron;EventT0 [ns];events / 0.5 ns",
99 nBins, minT0, maxT0);
100 m_histEventT0_SVD_hadron_L1_CDCTRG = new TH1F("m_histEventT0_SVD_hadron_L1_CDCTRG",
101 "SVD EventT0, L1TRG from CDC, HLT Hadron;EventT0 [ns];events / 0.5 ns",
102 nBins, minT0, maxT0);
103
104 m_histEventT0_ECL_mumu_L1_CDCTRG = new TH1F("m_histEventT0_ECL_mumu_L1_CDCTRG",
105 "ECL EventT0, L1TRG from CDC, HLT mumu;EventT0 [ns];events / 0.5 ns",
106 nBins, minT0, maxT0);
107 m_histEventT0_CDC_mumu_L1_CDCTRG = new TH1F("m_histEventT0_CDC_mumu_L1_CDCTRG",
108 "CDC EventT0, L1TRG from CDC, HLT mumu;EventT0 [ns];events / 0.5 ns",
109 nBins, minT0, maxT0);
110 m_histEventT0_TOP_mumu_L1_CDCTRG = new TH1F("m_histEventT0_TOP_mumu_L1_CDCTRG",
111 "TOP EventT0, L1TRG from CDC, HLT mumu;EventT0 [ns];events / 0.5 ns",
112 nBins, minT0, maxT0);
113 m_histEventT0_SVD_mumu_L1_CDCTRG = new TH1F("m_histEventT0_SVD_mumu_L1_CDCTRG",
114 "SVD EventT0, L1TRG from CDC, HLT mumu;EventT0 [ns];events / 0.5 ns",
115 nBins, minT0, maxT0);
116
117
118 m_histEventT0_ECL_bhabha_L1_TOPTRG = new TH1F("m_histEventT0_ECL_bhabha_L1_TOPTRG",
119 "ECL EventT0, L1TRG from TOP, HLT Bhabha;EventT0 [ns];events / 0.5 ns",
120 nBins, minT0, maxT0);
121 m_histEventT0_CDC_bhabha_L1_TOPTRG = new TH1F("m_histEventT0_CDC_bhabha_L1_TOPTRG",
122 "CDC EventT0, L1TRG from TOP, HLT Bhabha;EventT0 [ns];events / 0.5 ns",
123 nBins, minT0, maxT0);
124 m_histEventT0_TOP_bhabha_L1_TOPTRG = new TH1F("m_histEventT0_TOP_bhabha_L1_TOPTRG",
125 "TOP EventT0, L1TRG from TOP, HLT Bhabha;EventT0 [ns];events / 0.5 ns",
126 nBins, minT0, maxT0);
127 m_histEventT0_SVD_bhabha_L1_TOPTRG = new TH1F("m_histEventT0_SVD_bhabha_L1_TOPTRG",
128 "SVD EventT0, L1TRG from TOP, HLT Bhabha;EventT0 [ns];events / 0.5 ns",
129 nBins, minT0, maxT0);
130
131 m_histEventT0_ECL_hadron_L1_TOPTRG = new TH1F("m_histEventT0_ECL_hadron_L1_TOPTRG",
132 "ECL EventT0, L1TRG from TOP, HLT Hadron;EventT0 [ns];events / 0.5 ns",
133 nBins, minT0, maxT0);
134 m_histEventT0_CDC_hadron_L1_TOPTRG = new TH1F("m_histEventT0_CDC_hadron_L1_TOPTRG",
135 "CDC EventT0, L1TRG from TOP, HLT Hadron;EventT0 [ns];events / 0.5 ns",
136 nBins, minT0, maxT0);
137 m_histEventT0_TOP_hadron_L1_TOPTRG = new TH1F("m_histEventT0_TOP_hadron_L1_TOPTRG",
138 "TOP EventT0, L1TRG from TOP, HLT Hadron;EventT0 [ns];events / 0.5 ns",
139 nBins, minT0, maxT0);
140 m_histEventT0_SVD_hadron_L1_TOPTRG = new TH1F("m_histEventT0_SVD_hadron_L1_TOPTRG",
141 "SVD EventT0, L1TRG from TOP, HLT Hadron;EventT0 [ns];events / 0.5 ns",
142 nBins, minT0, maxT0);
143
144 m_histEventT0_ECL_mumu_L1_TOPTRG = new TH1F("m_histEventT0_ECL_mumu_L1_TOPTRG",
145 "ECL EventT0, L1TRG from TOP, HLT mumu;EventT0 [ns];events / 0.5 ns",
146 nBins, minT0, maxT0);
147 m_histEventT0_CDC_mumu_L1_TOPTRG = new TH1F("m_histEventT0_CDC_mumu_L1_TOPTRG",
148 "CDC EventT0, L1TRG from TOP, HLT mumu;EventT0 [ns];events / 0.5 ns",
149 nBins, minT0, maxT0);
150 m_histEventT0_TOP_mumu_L1_TOPTRG = new TH1F("m_histEventT0_TOP_mumu_L1_TOPTRG",
151 "TOP EventT0, L1TRG from TOP, HLT mumu;EventT0 [ns];events / 0.5 ns",
152 nBins, minT0, maxT0);
153 m_histEventT0_SVD_mumu_L1_TOPTRG = new TH1F("m_histEventT0_SVD_mumu_L1_TOPTRG",
154 "SVD EventT0, L1TRG from TOP, HLT mumu;EventT0 [ns];events / 0.5 ns",
155 nBins, minT0, maxT0);
156
158 new TH1D("AlgorithmSourceFractionsHadronL1ECLTRG",
159 "Fraction of events with EventT0 from each algorithm for hadronic events triggerd by ECL;Algorithm;Fraction",
160 6, 0, 6);
162 new TH1D("AlgorithmSourceFractionsHadronL1CDCTRG",
163 "Fraction of events with EventT0 from each algorithm for hadronic events triggerd by CDC;Algorithm;Fraction",
164 6, 0, 6);
166 new TH1D("AlgorithmSourceFractionsHadronL1TOPTRG",
167 "Fraction of events with EventT0 from each algorithm for hadronic events triggerd by TOP;Algorithm;Fraction",
168 6, 0, 6);
170 new TH1D("AlgorithmSourceFractionsBhaBhaL1ECLTRG",
171 "Fraction of events with EventT0 from each algorithm for Bhabha events triggerd by ECL;Algorithm;Fraction",
172 6, 0, 6);
174 new TH1D("AlgorithmSourceFractionsBhaBhaL1CDCTRG",
175 "Fraction of events with EventT0 from each algorithm for Bhabha events triggerd by CDC;Algorithm;Fraction",
176 6, 0, 6);
178 new TH1D("AlgorithmSourceFractionsBhaBhaL1TOPTRG",
179 "Fraction of events with EventT0 from each algorithm for Bhabha events triggerd by TOP;Algorithm;Fraction",
180 6, 0, 6);
182 new TH1D("AlgorithmSourceFractionsMuMuL1ECLTRG",
183 "Fraction of events with EventT0 from each algorithm for #mu#mu events triggerd by ECL;Algorithm;Fraction",
184 6, 0, 6);
186 new TH1D("AlgorithmSourceFractionsMuMuL1CDCTRG",
187 "Fraction of events with EventT0 from each algorithm for #mu#mu events triggerd by CDC;Algorithm;Fraction",
188 6, 0, 6);
190 new TH1D("AlgorithmSourceFractionsMuMuL1TOPTRG",
191 "Fraction of events with EventT0 from each algorithm for #mu#mu events triggerd by TOP;Algorithm;Fraction",
192 6, 0, 6);
193
194 for (uint i = 0; i < 6; i++) {
195 m_histAlgorithmSourceFractionsHadronL1ECLTRG->GetXaxis()->SetBinLabel(i + 1, c_eventT0Algorithms[i]);
196 m_histAlgorithmSourceFractionsHadronL1CDCTRG->GetXaxis()->SetBinLabel(i + 1, c_eventT0Algorithms[i]);
197 m_histAlgorithmSourceFractionsHadronL1TOPTRG->GetXaxis()->SetBinLabel(i + 1, c_eventT0Algorithms[i]);
198 m_histAlgorithmSourceFractionsBhaBhaL1ECLTRG->GetXaxis()->SetBinLabel(i + 1, c_eventT0Algorithms[i]);
199 m_histAlgorithmSourceFractionsBhaBhaL1CDCTRG->GetXaxis()->SetBinLabel(i + 1, c_eventT0Algorithms[i]);
200 m_histAlgorithmSourceFractionsBhaBhaL1TOPTRG->GetXaxis()->SetBinLabel(i + 1, c_eventT0Algorithms[i]);
201 m_histAlgorithmSourceFractionsMuMuL1ECLTRG->GetXaxis()->SetBinLabel(i + 1, c_eventT0Algorithms[i]);
202 m_histAlgorithmSourceFractionsMuMuL1CDCTRG->GetXaxis()->SetBinLabel(i + 1, c_eventT0Algorithms[i]);
203 m_histAlgorithmSourceFractionsMuMuL1TOPTRG->GetXaxis()->SetBinLabel(i + 1, c_eventT0Algorithms[i]);
204 }
205
206 oldDir->cd();
207
208}
209
210
211//---------------------------------
213{
214
215 m_TrgResult.isOptional();
216 m_eventT0.isOptional();
217
218 REG_HISTOGRAM
219
220}
221
222
223
224//---------------------------------
226{
227 if (!m_eventT0.isValid()) {
228 B2WARNING("Missing EventT0, EventT0DQM is skipped.");
229 return;
230 }
231
232 if (!m_TrgResult.isValid()) {
233 B2WARNING("Missing SoftwareTriggerResult, EventT0DQM is skipped.");
234 return;
235 }
236
241
246
251
256
261
266
271
276
281
291
292}
293
294
295//---------------------------------
297{
298 if (!m_objTrgSummary.isValid()) {
299 B2WARNING("TRGSummary object not available but required to indicate which detector provided the L1 trigger time");
300 return;
301 }
302
303 // Skip this event if there is no event t0, to avoid crashing other DQM
304 if (!m_eventT0.isOptional()) {
305 B2WARNING("Missing EventT0, EventT0DQM is skipped.");
306 return;
307 }
308
309 // Determine if there is a valid event t0 to use and then extract the information about it
310 if (!m_eventT0.isValid()) {
311 return ;
312 }
313
314 if (!m_TrgResult.isValid()) {
315 B2WARNING("SoftwareTriggerResult object not available but require to select bhabha/mumu/hadron events for this module");
316 return;
317 }
318
319 const std::map<std::string, int>& fresults = m_TrgResult->getResults();
320 if ((fresults.find("software_trigger_cut&skim&accept_bhabha") == fresults.end()) ||
321 (fresults.find("software_trigger_cut&skim&accept_mumu_2trk") == fresults.end()) ||
322 (fresults.find("software_trigger_cut&skim&accept_hadron") == fresults.end())) {
323 B2WARNING("EventT0DQMModule: Can't find required bhabha or mumu or hadron trigger identifier");
324 return;
325 }
326
327 // for L1 timing source is "ecl trigger" -> TRGSummary::ETimingType is 0
328 const bool IsECLL1TriggerSource = (m_objTrgSummary->getTimType() == TRGSummary::ETimingType::TTYP_ECL);
329 // for L1 timing source is "top trigger" -> TRGSummary::ETimingType is 1
330 const bool IsTOPL1TriggerSource = (m_objTrgSummary->getTimType() == TRGSummary::ETimingType::TTYP_TOP);
331 // for L1 timing source is "cdc trigger" -> TRGSummary::ETimingType is 3
332 const bool IsCDCL1TriggerSource = (m_objTrgSummary->getTimType() == TRGSummary::ETimingType::TTYP_CDC);
333
334 B2DEBUG(20, "IsECLL1TriggerSource = " << IsECLL1TriggerSource);
335 B2DEBUG(20, "IsTOPL1TriggerSource = " << IsTOPL1TriggerSource);
336 B2DEBUG(20, "IsCDCL1TriggerSource = " << IsCDCL1TriggerSource);
337
338 // determine if the event was part of the hadron skim or bhabha skim or mumu skim
339 const bool IsEvtAcceptedBhabha = (m_TrgResult->getResult("software_trigger_cut&skim&accept_bhabha") ==
341 const bool IsEvtAcceptedHadron = (m_TrgResult->getResult("software_trigger_cut&skim&accept_hadron") ==
343 const bool IsEvtAcceptedMumu = (m_TrgResult->getResult("software_trigger_cut&skim&accept_mumu_2trk") ==
345
346
347 B2DEBUG(20, "bhabha trigger result = " << static_cast<std::underlying_type<SoftwareTriggerCutResult>::type>
348 (m_TrgResult->getResult("software_trigger_cut&skim&accept_bhabha"))) ;
349 B2DEBUG(20, "hadron trigger result = " << static_cast<std::underlying_type<SoftwareTriggerCutResult>::type>
350 (m_TrgResult->getResult("software_trigger_cut&skim&accept_hadron"))) ;
351 B2DEBUG(20, "mu mu trigger result = " << static_cast<std::underlying_type<SoftwareTriggerCutResult>::type>
352 (m_TrgResult->getResult("software_trigger_cut&skim&accept_mumu_2trk"))) ;
353 B2DEBUG(20, "bhabha trigger comparison bool = " << IsEvtAcceptedBhabha) ;
354 B2DEBUG(20, "hadron trigger comparison bool = " << IsEvtAcceptedHadron) ;
355 B2DEBUG(20, "mumu trigger comparison bool = " << IsEvtAcceptedMumu) ;
356
357
358 // Set the different EventT0 values, default is -1000 in case there are no information based on a given detector
359 const double eventT0ECL =
360 m_eventT0->hasTemporaryEventT0(Const::EDetector::ECL) ? m_eventT0->getBestECLTemporaryEventT0()->eventT0 : -1000;
361 const double eventT0CDC =
362 m_eventT0->hasTemporaryEventT0(Const::EDetector::CDC) ? m_eventT0->getBestCDCTemporaryEventT0()->eventT0 : -1000;
363 const double eventT0TOP =
364 m_eventT0->hasTemporaryEventT0(Const::EDetector::TOP) ? m_eventT0->getBestTOPTemporaryEventT0()->eventT0 : -1000;
365 const double eventT0SVD =
366 m_eventT0->hasTemporaryEventT0(Const::EDetector::SVD) ? m_eventT0->getBestSVDTemporaryEventT0()->eventT0 : -1000;
367
368 const auto checkForCDCAlgorithm = [cdcEventT0s = m_eventT0->getTemporaryEventT0s(Const::EDetector::CDC)](
369 const std::string & algorithm) {
370 for (const auto& evtt0 : cdcEventT0s) {
371 if (evtt0.algorithm == algorithm) {
372 return true;
373 }
374 }
375 return false;
376 };
377
378 const bool hasCDCHitBasedEventT0 = checkForCDCAlgorithm("hit based");
379 const bool hasCDCFullGridEventT0 = checkForCDCAlgorithm("chi2");
380 // We are interested if an EventT0 is set, not whether temporary EventT0s exist that might not be used
381 const bool hasAnyEventT0 = m_eventT0->hasEventT0();
382 const bool hasECLEventT0 = m_eventT0->hasTemporaryEventT0(Const::EDetector::ECL);
383 const bool hasSVDEventT0 = m_eventT0->hasTemporaryEventT0(Const::EDetector::SVD);
384 const bool hasTOPEventT0 = m_eventT0->hasTemporaryEventT0(Const::EDetector::TOP);
385
386 // Fill the plots that used the ECL trigger as the L1 timing source
387 if (IsECLL1TriggerSource) {
388 // Fill the histograms with the event t0 values
389 if (IsEvtAcceptedBhabha) { // fill the bha bha skim event t0s
390 m_histEventT0_ECL_bhabha_L1_ECLTRG->Fill(eventT0ECL);
391 m_histEventT0_CDC_bhabha_L1_ECLTRG->Fill(eventT0CDC);
392 m_histEventT0_TOP_bhabha_L1_ECLTRG->Fill(eventT0TOP);
393 m_histEventT0_SVD_bhabha_L1_ECLTRG->Fill(eventT0SVD);
394 fillHistogram(m_histAlgorithmSourceFractionsBhaBhaL1ECLTRG, hasAnyEventT0, hasECLEventT0, hasSVDEventT0,
395 hasCDCHitBasedEventT0, hasCDCFullGridEventT0, hasTOPEventT0);
396 }
397
398 if (IsEvtAcceptedHadron) { // fill the hadron skim event t0s
399 m_histEventT0_ECL_hadron_L1_ECLTRG->Fill(eventT0ECL);
400 m_histEventT0_CDC_hadron_L1_ECLTRG->Fill(eventT0CDC);
401 m_histEventT0_TOP_hadron_L1_ECLTRG->Fill(eventT0TOP);
402 m_histEventT0_SVD_hadron_L1_ECLTRG->Fill(eventT0SVD);
403 fillHistogram(m_histAlgorithmSourceFractionsHadronL1ECLTRG, hasAnyEventT0, hasECLEventT0, hasSVDEventT0,
404 hasCDCHitBasedEventT0, hasCDCFullGridEventT0, hasTOPEventT0);
405 }
406
407 if (IsEvtAcceptedMumu) { // fill the mumu skim event t0s
408 m_histEventT0_ECL_mumu_L1_ECLTRG->Fill(eventT0ECL);
409 m_histEventT0_CDC_mumu_L1_ECLTRG->Fill(eventT0CDC);
410 m_histEventT0_TOP_mumu_L1_ECLTRG->Fill(eventT0TOP);
411 m_histEventT0_SVD_mumu_L1_ECLTRG->Fill(eventT0SVD);
412 fillHistogram(m_histAlgorithmSourceFractionsMuMuL1ECLTRG, hasAnyEventT0, hasECLEventT0, hasSVDEventT0,
413 hasCDCHitBasedEventT0, hasCDCFullGridEventT0, hasTOPEventT0);
414 }
415 }
416 // Fill the plots that used the TOP trigger as the L1 timing source
417 else if (IsTOPL1TriggerSource) {
418 // Fill the histograms with the event t0 values
419 if (IsEvtAcceptedBhabha) { // fill the bha bha skim event t0s
420 m_histEventT0_ECL_bhabha_L1_TOPTRG->Fill(eventT0ECL);
421 m_histEventT0_CDC_bhabha_L1_TOPTRG->Fill(eventT0CDC);
422 m_histEventT0_TOP_bhabha_L1_TOPTRG->Fill(eventT0TOP);
423 m_histEventT0_SVD_bhabha_L1_TOPTRG->Fill(eventT0SVD);
424 fillHistogram(m_histAlgorithmSourceFractionsBhaBhaL1TOPTRG, hasAnyEventT0, hasECLEventT0, hasSVDEventT0,
425 hasCDCHitBasedEventT0, hasCDCFullGridEventT0, hasTOPEventT0);
426 }
427
428 if (IsEvtAcceptedHadron) { // fill the hadron skim event t0s
429 m_histEventT0_ECL_hadron_L1_TOPTRG->Fill(eventT0ECL);
430 m_histEventT0_CDC_hadron_L1_TOPTRG->Fill(eventT0CDC);
431 m_histEventT0_TOP_hadron_L1_TOPTRG->Fill(eventT0TOP);
432 m_histEventT0_SVD_hadron_L1_TOPTRG->Fill(eventT0SVD);
433 fillHistogram(m_histAlgorithmSourceFractionsHadronL1TOPTRG, hasAnyEventT0, hasECLEventT0, hasSVDEventT0,
434 hasCDCHitBasedEventT0, hasCDCFullGridEventT0, hasTOPEventT0);
435 }
436
437 if (IsEvtAcceptedMumu) { // fill the mumu skim event t0s
438 m_histEventT0_ECL_mumu_L1_TOPTRG->Fill(eventT0ECL);
439 m_histEventT0_CDC_mumu_L1_TOPTRG->Fill(eventT0CDC);
440 m_histEventT0_TOP_mumu_L1_TOPTRG->Fill(eventT0TOP);
441 m_histEventT0_SVD_mumu_L1_TOPTRG->Fill(eventT0SVD);
442 fillHistogram(m_histAlgorithmSourceFractionsMuMuL1TOPTRG, hasAnyEventT0, hasECLEventT0, hasSVDEventT0,
443 hasCDCHitBasedEventT0, hasCDCFullGridEventT0, hasTOPEventT0);
444 }
445 }
446 // Fill the plots that used the CDC trigger as the L1 timing source
447 else if (IsCDCL1TriggerSource) {
448 // Fill the histograms with the event t0 values
449 if (IsEvtAcceptedBhabha) { // fill the bha bha skim event t0s
450 m_histEventT0_ECL_bhabha_L1_CDCTRG->Fill(eventT0ECL);
451 m_histEventT0_CDC_bhabha_L1_CDCTRG->Fill(eventT0CDC);
452 m_histEventT0_TOP_bhabha_L1_CDCTRG->Fill(eventT0TOP);
453 m_histEventT0_SVD_bhabha_L1_CDCTRG->Fill(eventT0SVD);
454 fillHistogram(m_histAlgorithmSourceFractionsBhaBhaL1CDCTRG, hasAnyEventT0, hasECLEventT0, hasSVDEventT0,
455 hasCDCHitBasedEventT0, hasCDCFullGridEventT0, hasTOPEventT0);
456 }
457
458 if (IsEvtAcceptedHadron) { // fill the hadron skim event t0s
459 m_histEventT0_ECL_hadron_L1_CDCTRG->Fill(eventT0ECL);
460 m_histEventT0_CDC_hadron_L1_CDCTRG->Fill(eventT0CDC);
461 m_histEventT0_TOP_hadron_L1_CDCTRG->Fill(eventT0TOP);
462 m_histEventT0_SVD_hadron_L1_CDCTRG->Fill(eventT0SVD);
463 fillHistogram(m_histAlgorithmSourceFractionsHadronL1CDCTRG, hasAnyEventT0, hasECLEventT0, hasSVDEventT0,
464 hasCDCHitBasedEventT0, hasCDCFullGridEventT0, hasTOPEventT0);
465 }
466
467 if (IsEvtAcceptedMumu) { // fill the mumu skim event t0s
468 m_histEventT0_ECL_mumu_L1_CDCTRG->Fill(eventT0ECL);
469 m_histEventT0_CDC_mumu_L1_CDCTRG->Fill(eventT0CDC);
470 m_histEventT0_TOP_mumu_L1_CDCTRG->Fill(eventT0TOP);
471 m_histEventT0_SVD_mumu_L1_CDCTRG->Fill(eventT0SVD);
472 fillHistogram(m_histAlgorithmSourceFractionsMuMuL1CDCTRG, hasAnyEventT0, hasECLEventT0, hasSVDEventT0,
473 hasCDCHitBasedEventT0, hasCDCFullGridEventT0, hasTOPEventT0);
474 }
475 }
476
477 B2DEBUG(20, "eventT0ECL = " << eventT0ECL << " ns") ;
478 B2DEBUG(20, "eventT0CDC = " << eventT0CDC << " ns") ;
479 B2DEBUG(20, "eventT0TOP = " << eventT0TOP << " ns") ;
480 B2DEBUG(20, "eventT0SVD = " << eventT0SVD << " ns") ;
481}
482
483void EventT0DQMModule::fillHistogram(TH1D* hist, const bool hasAnyT0, const bool hasECLT0, const bool hasSVDT0,
484 const bool hasCDCHitT0, const bool hasCDCGridT0, const bool hasTOPT0)
485{
486 hist->Fill(-1); // counting events for normalisation
487 hist->Fill(c_eventT0Algorithms[0], hasECLT0);
488 hist->Fill(c_eventT0Algorithms[1], hasSVDT0);
489 hist->Fill(c_eventT0Algorithms[2], hasCDCHitT0);
490 hist->Fill(c_eventT0Algorithms[3], hasCDCGridT0);
491 hist->Fill(c_eventT0Algorithms[4], hasTOPT0);
492 hist->Fill(c_eventT0Algorithms[5], hasAnyT0);
493}
TH1F * m_histEventT0_CDC_bhabha_L1_TOPTRG
event t0 histogram for CDC, HLT bha bha events, L1 time by TOP trigger
Definition: EventT0DQM.h:98
TH1D * m_histAlgorithmSourceFractionsBhaBhaL1CDCTRG
Fraction of events with EventT0 from a given algorithm, HLT bhabha events, L1 time by CDC trigger.
Definition: EventT0DQM.h:123
TH1F * m_histEventT0_CDC_bhabha_L1_CDCTRG
event t0 histogram for CDC, HLT bha bha events, L1 time by CDC trigger
Definition: EventT0DQM.h:82
virtual ~EventT0DQMModule()
Destructor.
Definition: EventT0DQM.cc:24
TH1F * m_histEventT0_TOP_bhabha_L1_TOPTRG
event t0 histogram for TOP, HLT bha bha events, L1 time by TOP trigger
Definition: EventT0DQM.h:99
TH1F * m_histEventT0_ECL_hadron_L1_ECLTRG
event t0 histogram for ECL, HLT hadronic events, L1 time by ECL trigger
Definition: EventT0DQM.h:70
TH1F * m_histEventT0_ECL_bhabha_L1_CDCTRG
event t0 histogram for ECL, HLT bha bha events, L1 time by CDC trigger
Definition: EventT0DQM.h:81
TH1F * m_histEventT0_SVD_mumu_L1_ECLTRG
event t0 histogram for SVD, HLT mu mu events, L1 time by ECL trigger
Definition: EventT0DQM.h:78
StoreObjPtr< EventT0 > m_eventT0
Store array for event t0.
Definition: EventT0DQM.h:62
const char * c_eventT0Algorithms[6]
EventT0 algorithms for which to calculate fractions of abundance.
Definition: EventT0DQM.h:113
TH1F * m_histEventT0_ECL_hadron_L1_TOPTRG
event t0 histogram for ECL, HLT hadronic events, L1 time by TOP trigger
Definition: EventT0DQM.h:102
TH1F * m_histEventT0_CDC_mumu_L1_TOPTRG
event t0 histogram for CDC, HLT mu mu events, L1 time by TOP trigger
Definition: EventT0DQM.h:108
TH1F * m_histEventT0_CDC_hadron_L1_ECLTRG
event t0 histogram for CDC, HLT hadronic events, L1 time by ECL trigger
Definition: EventT0DQM.h:71
virtual void initialize() override
Initialize the module.
Definition: EventT0DQM.cc:212
TH1F * m_histEventT0_SVD_bhabha_L1_CDCTRG
event t0 histogram for SVD, HLT bha bha events, L1 time by CDC trigger
Definition: EventT0DQM.h:84
TH1F * m_histEventT0_CDC_mumu_L1_CDCTRG
event t0 histogram for CDC, HLT mu mu events, L1 time by CDC trigger
Definition: EventT0DQM.h:92
virtual void event() override
This method is called for each event.
Definition: EventT0DQM.cc:296
TH1D * m_histAlgorithmSourceFractionsHadronL1CDCTRG
Fraction of events with EventT0 from a given algorithm, HLT hadronic events, L1 time by CDC trigger.
Definition: EventT0DQM.h:117
TH1F * m_histEventT0_SVD_bhabha_L1_TOPTRG
event t0 histogram for SVD, HLT bha bha events, L1 time by TOP trigger
Definition: EventT0DQM.h:100
TH1F * m_histEventT0_SVD_hadron_L1_ECLTRG
event t0 histogram for SVD, HLT hadronic events, L1 time by ECL trigger
Definition: EventT0DQM.h:73
TH1F * m_histEventT0_CDC_mumu_L1_ECLTRG
event t0 histogram for CDC, HLT mu mu events, L1 time by ECL trigger
Definition: EventT0DQM.h:76
TH1F * m_histEventT0_SVD_bhabha_L1_ECLTRG
event t0 histogram for SVD, HLT bha bha events, L1 time by ECL trigger
Definition: EventT0DQM.h:68
TH1F * m_histEventT0_SVD_mumu_L1_CDCTRG
event t0 histogram for SVD, HLT mu mu events, L1 time by CDC trigger
Definition: EventT0DQM.h:94
TH1F * m_histEventT0_TOP_bhabha_L1_ECLTRG
event t0 histogram for TOP, HLT bha bha events, L1 time by ECL trigger
Definition: EventT0DQM.h:67
TH1F * m_histEventT0_ECL_mumu_L1_TOPTRG
event t0 histogram for ECL, HLT mu mu events, L1 time by TOP trigger
Definition: EventT0DQM.h:107
StoreObjPtr< TRGSummary > m_objTrgSummary
Trigger Summary data object.
Definition: EventT0DQM.h:59
TH1D * m_histAlgorithmSourceFractionsHadronL1ECLTRG
Fraction of events with EventT0 from a given algorithm, HLT hadronic events, L1 time by ECL trigger.
Definition: EventT0DQM.h:115
virtual void beginRun() override
This method is called for each run.
Definition: EventT0DQM.cc:225
void fillHistogram(TH1D *hist, const bool hasAnyT0, const bool hasECLT0, const bool hasSVDT0, const bool hasCDCHitT0, const bool hasCDCGridT0, const bool hasTOPT0)
Fill fraction histogram with values.
Definition: EventT0DQM.cc:483
TH1F * m_histEventT0_TOP_mumu_L1_ECLTRG
event t0 histogram for TOP, HLT mu mu events, L1 time by ECL trigger
Definition: EventT0DQM.h:77
TH1F * m_histEventT0_TOP_hadron_L1_CDCTRG
event t0 histogram for TOP, HLT hadronic events, L1 time by CDC trigger
Definition: EventT0DQM.h:88
TH1D * m_histAlgorithmSourceFractionsMuMuL1CDCTRG
Fraction of events with EventT0 from a given algorithm, HLT mumu events, L1 time by CDC trigger.
Definition: EventT0DQM.h:129
TH1F * m_histEventT0_SVD_hadron_L1_CDCTRG
event t0 histogram for SVD, HLT hadronic events, L1 time by CDC trigger
Definition: EventT0DQM.h:89
TH1F * m_histEventT0_TOP_hadron_L1_TOPTRG
event t0 histogram for TOP, HLT hadronic events, L1 time by TOP trigger
Definition: EventT0DQM.h:104
TH1D * m_histAlgorithmSourceFractionsBhaBhaL1ECLTRG
Fraction of events with EventT0 from a given algorithm, HLT bhabha events, L1 time by ECL trigger.
Definition: EventT0DQM.h:121
TH1D * m_histAlgorithmSourceFractionsBhaBhaL1TOPTRG
Fraction of events with EventT0 from a given algorithm, HLT bhabha events, L1 time by TOP trigger.
Definition: EventT0DQM.h:125
TH1F * m_histEventT0_SVD_mumu_L1_TOPTRG
event t0 histogram for SVD, HLT mu mu events, L1 time by TOP trigger
Definition: EventT0DQM.h:110
TH1F * m_histEventT0_SVD_hadron_L1_TOPTRG
event t0 histogram for SVD, HLT hadronic events, L1 time by TOP trigger
Definition: EventT0DQM.h:105
TH1F * m_histEventT0_ECL_mumu_L1_CDCTRG
event t0 histogram for ECL, HLT mu mu events, L1 time by CDC trigger
Definition: EventT0DQM.h:91
TH1F * m_histEventT0_TOP_mumu_L1_TOPTRG
event t0 histogram for TOP, HLT mu mu events, L1 time by TOP trigger
Definition: EventT0DQM.h:109
TH1F * m_histEventT0_ECL_mumu_L1_ECLTRG
event t0 histogram for ECL, HLT mu mu events, L1 time by ECL trigger
Definition: EventT0DQM.h:75
EventT0DQMModule()
Default constructor.
Definition: EventT0DQM.cc:17
TH1F * m_histEventT0_TOP_mumu_L1_CDCTRG
event t0 histogram for TOP, HLT mu mu events, L1 time by CDC trigger
Definition: EventT0DQM.h:93
TH1F * m_histEventT0_TOP_bhabha_L1_CDCTRG
event t0 histogram for TOP, HLT bha bha events, L1 time by CDC trigger
Definition: EventT0DQM.h:83
TH1F * m_histEventT0_TOP_hadron_L1_ECLTRG
event t0 histogram for TOP, HLT hadronic events, L1 time by ECL trigger
Definition: EventT0DQM.h:72
TH1F * m_histEventT0_ECL_bhabha_L1_TOPTRG
event t0 histogram for ECL, HLT bha bha events, L1 time by TOP trigger
Definition: EventT0DQM.h:97
TH1D * m_histAlgorithmSourceFractionsMuMuL1TOPTRG
Fraction of events with EventT0 from a given algorithm, HLT mumu events, L1 time by TOP trigger.
Definition: EventT0DQM.h:131
TH1F * m_histEventT0_CDC_hadron_L1_TOPTRG
event t0 histogram for CDC, HLT hadronic events, L1 time by TOP trigger
Definition: EventT0DQM.h:103
TH1F * m_histEventT0_CDC_hadron_L1_CDCTRG
event t0 histogram for CDC, HLT hadronic events, L1 time by CDC trigger
Definition: EventT0DQM.h:87
StoreObjPtr< SoftwareTriggerResult > m_TrgResult
Store array for Trigger selection.
Definition: EventT0DQM.h:61
TH1D * m_histAlgorithmSourceFractionsMuMuL1ECLTRG
Fraction of events with EventT0 from a given algorithm, HLT mumu events, L1 time by ECL trigger.
Definition: EventT0DQM.h:127
TH1F * m_histEventT0_CDC_bhabha_L1_ECLTRG
event t0 histogram for CDC, HLT bha bha events, L1 time by ECL trigger
Definition: EventT0DQM.h:66
TH1F * m_histEventT0_ECL_hadron_L1_CDCTRG
event t0 histogram for ECL, HLT hadronic events, L1 time by CDC trigger
Definition: EventT0DQM.h:86
virtual void defineHisto() override
Defination of histograms.
Definition: EventT0DQM.cc:28
TH1F * m_histEventT0_ECL_bhabha_L1_ECLTRG
event t0 histogram for ECL, HLT bha bha events, L1 time by ECL trigger
Definition: EventT0DQM.h:65
TH1D * m_histAlgorithmSourceFractionsHadronL1TOPTRG
Fraction of events with EventT0 from a given algorithm, HLT hadronic events, L1 time by TOP trigger.
Definition: EventT0DQM.h:119
HistoModule.h is supposed to be used instead of Module.h for the modules with histogram definitions t...
Definition: HistoModule.h:29
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
@ 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
@ TTYP_TOP
events triggered by TOP timing
Definition: TRGSummary.h:61
@ TTYP_CDC
events triggered by CDC timing
Definition: TRGSummary.h:63
@ TTYP_ECL
events triggered by ECL timing
Definition: TRGSummary.h:45
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:650
@ c_accept
Accept this event.
Abstract base class for different kinds of events.