Belle II Software  release-05-01-25
EventT0DQM.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2012 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Ewan Hill (ehill@mail.ubc.ca) *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include <reconstruction/modules/EventT0DQM/EventT0DQM.h>
12 #include <framework/core/HistoModule.h>
13 
14 using namespace Belle2;
15 
16 REG_MODULE(EventT0DQM)
17 
18 //---------------------------------
20 {
21  setPropertyFlags(c_ParallelProcessingCertified); // parallel processing
22  setDescription("Make data quality monitoring plots for event t0 for bhabha, mu mu, and hadron samples seeded by different trigger times.");
23 }
24 
25 //---------------------------------
27 
28 
29 //---------------------------------
31 {
32 
33  TDirectory* oldDir = gDirectory;
34  oldDir->mkdir("EventT0DQMdir")->cd();
35 
36  int nBins = 400 ;
37  double minT0 = -100 ;
38  double maxT0 = 100 ;
39 
40  m_histEventT0_ECL_bhabha_L1_ECLTRG = new TH1F("m_histEventT0_ECL_bhabha_L1_ECLTRG",
41  "ECL event t0 - bhabhas - ECLTRG time;event t0 [ns];events / 0.5 ns",
42  nBins, minT0, maxT0);
43  m_histEventT0_CDC_bhabha_L1_ECLTRG = new TH1F("m_histEventT0_CDC_bhabha_L1_ECLTRG",
44  "CDC event t0 - bhabhas - ECLTRG time;event t0 [ns];events / 0.5 ns",
45  nBins, minT0, maxT0);
46  m_histEventT0_TOP_bhabha_L1_ECLTRG = new TH1F("m_histEventT0_TOP_bhabha_L1_ECLTRG",
47  "TOP event t0 - bhabhas - ECLTRG time;event t0 [ns];events / 0.5 ns",
48  nBins, minT0, maxT0);
49  m_histEventT0_ECL_hadron_L1_ECLTRG = new TH1F("m_histEventT0_ECL_hadron_L1_ECLTRG",
50  "ECL event t0 - hadrons - ECLTRG time;event t0 [ns];events / 0.5 ns",
51  nBins, minT0, maxT0);
52  m_histEventT0_CDC_hadron_L1_ECLTRG = new TH1F("m_histEventT0_CDC_hadron_L1_ECLTRG",
53  "CDC event t0 - hadrons - ECLTRG time;event t0 [ns];events / 0.5 ns",
54  nBins, minT0, maxT0);
55  m_histEventT0_TOP_hadron_L1_ECLTRG = new TH1F("m_histEventT0_TOP_hadron_L1_ECLTRG",
56  "TOP event t0 - hadrons - ECLTRG time;event t0 [ns];events / 0.5 ns",
57  nBins, minT0, maxT0);
58  m_histEventT0_ECL_mumu_L1_ECLTRG = new TH1F("m_histEventT0_ECL_mumu_L1_ECLTRG",
59  "ECL event t0 - mu mu - ECLTRG time;event t0 [ns];events / 0.5 ns",
60  nBins, minT0, maxT0);
61  m_histEventT0_CDC_mumu_L1_ECLTRG = new TH1F("m_histEventT0_CDC_mumu_L1_ECLTRG",
62  "CDC event t0 - mu mu - ECLTRG time;event t0 [ns];events / 0.5 ns",
63  nBins, minT0, maxT0);
64  m_histEventT0_TOP_mumu_L1_ECLTRG = new TH1F("m_histEventT0_TOP_mumu_L1_ECLTRG",
65  "TOP event t0 - mu mu - ECLTRG time;event t0 [ns];events / 0.5 ns",
66  nBins, minT0, maxT0);
67 
68 
69 
70 
71 
72  m_histEventT0_ECL_bhabha_L1_CDCTRG = new TH1F("m_histEventT0_ECL_bhabha_L1_CDCTRG",
73  "ECL event t0 - bhabhas - CDCTRG time;event t0 [ns];events / 0.5 ns",
74  nBins, minT0, maxT0);
75  m_histEventT0_CDC_bhabha_L1_CDCTRG = new TH1F("m_histEventT0_CDC_bhabha_L1_CDCTRG",
76  "CDC event t0 - bhabhas - CDCTRG time;event t0 [ns];events / 0.5 ns",
77  nBins, minT0, maxT0);
78  m_histEventT0_TOP_bhabha_L1_CDCTRG = new TH1F("m_histEventT0_TOP_bhabha_L1_CDCTRG",
79  "TOP event t0 - bhabhas - CDCTRG time;event t0 [ns];events / 0.5 ns",
80  nBins, minT0, maxT0);
81  m_histEventT0_ECL_hadron_L1_CDCTRG = new TH1F("m_histEventT0_ECL_hadron_L1_CDCTRG",
82  "ECL event t0 - hadrons - CDCTRG time;event t0 [ns];events / 0.5 ns",
83  nBins, minT0, maxT0);
84  m_histEventT0_CDC_hadron_L1_CDCTRG = new TH1F("m_histEventT0_CDC_hadron_L1_CDCTRG",
85  "CDC event t0 - hadrons - CDCTRG time;event t0 [ns];events / 0.5 ns",
86  nBins, minT0, maxT0);
87  m_histEventT0_TOP_hadron_L1_CDCTRG = new TH1F("m_histEventT0_TOP_hadron_L1_CDCTRG",
88  "TOP event t0 - hadrons - CDCTRG time;event t0 [ns];events / 0.5 ns",
89  nBins, minT0, maxT0);
90  m_histEventT0_ECL_mumu_L1_CDCTRG = new TH1F("m_histEventT0_ECL_mumu_L1_CDCTRG",
91  "ECL event t0 - mu mu - CDCTRG time;event t0 [ns];events / 0.5 ns",
92  nBins, minT0, maxT0);
93  m_histEventT0_CDC_mumu_L1_CDCTRG = new TH1F("m_histEventT0_CDC_mumu_L1_CDCTRG",
94  "CDC event t0 - mu mu - CDCTRG time;event t0 [ns];events / 0.5 ns",
95  nBins, minT0, maxT0);
96  m_histEventT0_TOP_mumu_L1_CDCTRG = new TH1F("m_histEventT0_TOP_mumu_L1_CDCTRG",
97  "TOP event t0 - mu mu - CDCTRG time;event t0 [ns];events / 0.5 ns",
98  nBins, minT0, maxT0);
99 
100  oldDir->cd();
101 
102 }
103 
104 
105 //---------------------------------
107 {
108 
109  m_TrgResult.isOptional();
110  m_eventT0.isOptional();
111 
112  REG_HISTOGRAM
113 
114 }
115 
116 
117 
118 //---------------------------------
120 {
121  if (!m_eventT0.isValid()) {
122  B2WARNING("Missing event t0, EventT0DQM is skipped.");
123  return;
124  }
125 
126  if (!m_TrgResult.isValid()) {
127  B2WARNING("Missing SoftwareTriggerResult, EventT0DQM is skipped.");
128  return;
129  }
130 
140 
150 
151 }
152 
153 
154 //---------------------------------
156 {
157  if (!m_objTrgSummary.isValid()) {
158  B2WARNING("TRGSummary object not available but required to indicate which detector provided the L1 trigger time");
159  return;
160  } else {
161  m_L1TimingSrc = m_objTrgSummary->getTimType();
162  }
163 
164  bool Is_ECL_L1TriggerSource = false ;
165  bool Is_CDC_L1TriggerSource = false ;
166  if (m_L1TimingSrc == 0) { // for L1 timing source is "ecl trigger"
167  Is_ECL_L1TriggerSource = true ;
168  } else if (m_L1TimingSrc == 3) { // for L1 timing source is "cdc trigger"
169  Is_CDC_L1TriggerSource = true ;
170  }
171  // else if(m_L1TimingSrc==5){ // for L1 timing source is "delayed Bhabha" }
172  B2DEBUG(20, "Is_ECL_L1TriggerSource = " << Is_ECL_L1TriggerSource) ;
173  B2DEBUG(20, "Is_CDC_L1TriggerSource= " << Is_CDC_L1TriggerSource) ;
174 
175 
176  if (!m_TrgResult.isValid()) {
177  B2WARNING("SoftwareTriggerResult object not available but require to select bhabha/mumu/hadron events for this module");
178  return;
179  }
180 
181  const std::map<std::string, int>& fresults = m_TrgResult->getResults();
182  if ((fresults.find("software_trigger_cut&skim&accept_bhabha") == fresults.end()) ||
183  (fresults.find("software_trigger_cut&skim&accept_mumu_2trk") == fresults.end()) ||
184  (fresults.find("software_trigger_cut&skim&accept_hadron") == fresults.end())) {
185  B2WARNING("EventT0DQMModule: Can't find required bhabha or mumu or hadron trigger identifier");
186  return;
187  }
188 
189 
190 
191  // Skip this event if there is no event t0, to avoid crashing other DQM
192  if (!m_eventT0.isOptional()) {
193  B2WARNING("Missing event t0, EventT0DQM is skipped.");
194  return;
195  }
196 
197  // Determine if there is a valid event t0 to use and then extract the information about it
198  if (!m_eventT0.isValid()) {
199  return ;
200  }
201 
202  // determine if the event was part of the hadron skim or bhabha skim or mumu skim
203  const bool IsEvtAcceptedBhabha = (m_TrgResult->getResult("software_trigger_cut&skim&accept_bhabha") ==
205  const bool IsEvtAcceptedHadron = (m_TrgResult->getResult("software_trigger_cut&skim&accept_hadron") ==
207  const bool IsEvtAcceptedMumu = (m_TrgResult->getResult("software_trigger_cut&skim&accept_mumu_2trk") ==
209 
210 
211  B2DEBUG(20, "bhabha trigger result = " << static_cast<std::underlying_type<SoftwareTriggerCutResult>::type>
212  (m_TrgResult->getResult("software_trigger_cut&skim&accept_bhabha"))) ;
213  B2DEBUG(20, "hadron trigger result = " << static_cast<std::underlying_type<SoftwareTriggerCutResult>::type>
214  (m_TrgResult->getResult("software_trigger_cut&skim&accept_hadron"))) ;
215  B2DEBUG(20, "mu mu trigger result = " << static_cast<std::underlying_type<SoftwareTriggerCutResult>::type>
216  (m_TrgResult->getResult("software_trigger_cut&skim&accept_mumu_2trk"))) ;
217  B2DEBUG(20, "bhabha trigger comparison bool = " << IsEvtAcceptedBhabha) ;
218  B2DEBUG(20, "hadron trigger comparison bool = " << IsEvtAcceptedHadron) ;
219  B2DEBUG(20, "mumu trigger comparison bool = " << IsEvtAcceptedMumu) ;
220 
221 
222  // default values of the event t0 given that there may not be a value for every event depending on the detector measuring it.
223  double eventT0_ECL = -1000 ;
224  double eventT0_CDC = -1000 ;
225  double eventT0_TOP = -1000 ;
226 
227  // Set the CDC event t0 value if it exists
228  if (m_eventT0->hasTemporaryEventT0(Const::EDetector::CDC)) {
229  auto evtT0List_CDC = m_eventT0->getTemporaryEventT0s(Const::EDetector::CDC) ;
230 
231  // set the CDC event t0 value for filling into the histogram
232  // The most accurate CDC event t0 value is the last one in the list.
233  eventT0_CDC = evtT0List_CDC.back().eventT0 ;
234  }
235 
236  // Set the ECL event t0 value if it exists
237  if (m_eventT0->hasTemporaryEventT0(Const::EDetector::ECL)) {
238  // Get the list of ECL event t0 values. There are several event t0 values, not just one.
239  auto evtT0List_ECL = m_eventT0->getTemporaryEventT0s(Const::EDetector::ECL) ;
240 
241  // Select the event t0 value from the ECL as the one with the smallest chi squared value (defined as ".quality")
242  double smallest_ECL_t0_minChi2 = evtT0List_ECL[0].quality ;
243  int smallest_ECL_t0_minChi2_idx = 0 ;
244  for (long unsigned int ECLi = 0; ECLi < evtT0List_ECL.size(); ECLi++) {
245  if (evtT0List_ECL[ECLi].quality < smallest_ECL_t0_minChi2) {
246  smallest_ECL_t0_minChi2 = evtT0List_ECL[ECLi].quality ;
247  smallest_ECL_t0_minChi2_idx = ECLi ;
248  }
249  }
250  // set the ECL event t0 value for filling into the histogram
251  // It is the value found to have the small chi square
252  eventT0_ECL = evtT0List_ECL[smallest_ECL_t0_minChi2_idx].eventT0 ;
253  }
254 
255  // Set the TOP event t0 value if it exists
256  if (m_eventT0->hasTemporaryEventT0(Const::EDetector::TOP)) {
257  auto evtT0List_TOP = m_eventT0->getTemporaryEventT0s(Const::EDetector::TOP) ;
258 
259  // set the TOP event t0 value for filling into the histogram
260  // There should only be at most one value in the list per event
261  eventT0_TOP = evtT0List_TOP.back().eventT0 ;
262  }
263 
264 
265 
266 
267  // Fill the plots that used the ECL trigger as the L1 timing source
268  if (Is_ECL_L1TriggerSource) {
269  // Fill the histograms with the event t0 values
270  if (IsEvtAcceptedBhabha) { // fill the bha bha skim event t0s
271  m_histEventT0_ECL_bhabha_L1_ECLTRG->Fill(eventT0_ECL);
272  m_histEventT0_CDC_bhabha_L1_ECLTRG->Fill(eventT0_CDC);
273  m_histEventT0_TOP_bhabha_L1_ECLTRG->Fill(eventT0_TOP);
274  }
275 
276  if (IsEvtAcceptedHadron) { // fill the hadron skim event t0s
277  m_histEventT0_ECL_hadron_L1_ECLTRG->Fill(eventT0_ECL);
278  m_histEventT0_CDC_hadron_L1_ECLTRG->Fill(eventT0_CDC);
279  m_histEventT0_TOP_hadron_L1_ECLTRG->Fill(eventT0_TOP);
280  }
281 
282  if (IsEvtAcceptedMumu) { // fill the mumu skim event t0s
283  m_histEventT0_ECL_mumu_L1_ECLTRG->Fill(eventT0_ECL);
284  m_histEventT0_CDC_mumu_L1_ECLTRG->Fill(eventT0_CDC);
285  m_histEventT0_TOP_mumu_L1_ECLTRG->Fill(eventT0_TOP);
286  }
287  }
288  // Fill the plots that used the CDC trigger as the L1 timing source
289  else if (Is_CDC_L1TriggerSource) {
290  // Fill the histograms with the event t0 values
291  if (IsEvtAcceptedBhabha) { // fill the bha bha skim event t0s
292  m_histEventT0_ECL_bhabha_L1_CDCTRG->Fill(eventT0_ECL);
293  m_histEventT0_CDC_bhabha_L1_CDCTRG->Fill(eventT0_CDC);
294  m_histEventT0_TOP_bhabha_L1_CDCTRG->Fill(eventT0_TOP);
295  }
296 
297  if (IsEvtAcceptedHadron) { // fill the hadron skim event t0s
298  m_histEventT0_ECL_hadron_L1_CDCTRG->Fill(eventT0_ECL);
299  m_histEventT0_CDC_hadron_L1_CDCTRG->Fill(eventT0_CDC);
300  m_histEventT0_TOP_hadron_L1_CDCTRG->Fill(eventT0_TOP);
301  }
302 
303  if (IsEvtAcceptedMumu) { // fill the mumu skim event t0s
304  m_histEventT0_ECL_mumu_L1_CDCTRG->Fill(eventT0_ECL);
305  m_histEventT0_CDC_mumu_L1_CDCTRG->Fill(eventT0_CDC);
306  m_histEventT0_TOP_mumu_L1_CDCTRG->Fill(eventT0_TOP);
307  }
308  }
309 
310  B2DEBUG(20, "eventT0_ECL = " << eventT0_ECL << " ns") ;
311  B2DEBUG(20, "eventT0_CDC = " << eventT0_CDC << " ns") ;
312  B2DEBUG(20, "eventT0_TOP = " << eventT0_TOP << " ns") ;
313 }
314 
315 
316 //---------------------------------
318 {
319 
320 }
321 
322 
323 //---------------------------------
325 {
326 
327 }
Belle2::EventT0DQMModule::m_histEventT0_TOP_bhabha_L1_ECLTRG
TH1F * m_histEventT0_TOP_bhabha_L1_ECLTRG
event t0 histogram for TOP and bha bha events wrt the ECL trigger time
Definition: EventT0DQM.h:90
Belle2::EventT0DQMModule::terminate
virtual void terminate() override
End of the event processing.
Definition: EventT0DQM.cc:324
Belle2::EventT0DQMModule::~EventT0DQMModule
virtual ~EventT0DQMModule()
Destructor.
Definition: EventT0DQM.cc:26
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::EventT0DQMModule::m_histEventT0_ECL_mumu_L1_ECLTRG
TH1F * m_histEventT0_ECL_mumu_L1_ECLTRG
event t0 histogram for ECL and mu mu events wrt the ECL trigger time
Definition: EventT0DQM.h:94
Belle2::EventT0DQMModule::m_eventT0
StoreObjPtr< EventT0 > m_eventT0
Store array for event t0.
Definition: EventT0DQM.h:86
Belle2::EventT0DQMModule::beginRun
virtual void beginRun() override
This method is called for each run.
Definition: EventT0DQM.cc:119
Belle2::EventT0DQMModule::m_histEventT0_ECL_bhabha_L1_ECLTRG
TH1F * m_histEventT0_ECL_bhabha_L1_ECLTRG
event t0 histogram for ECL and bha bha events wrt the ECL trigger time
Definition: EventT0DQM.h:88
Belle2::EventT0DQMModule::m_L1TimingSrc
int m_L1TimingSrc
L1 timing source from getTimeType() in TRGSummary See ETimingTYpe in mdst/dataobjects/include/TRGSumm...
Definition: EventT0DQM.h:83
Belle2::EventT0DQMModule::m_histEventT0_CDC_bhabha_L1_CDCTRG
TH1F * m_histEventT0_CDC_bhabha_L1_CDCTRG
event t0 histogram for CDC and bha bha events wrt the CDC trigger time
Definition: EventT0DQM.h:99
Belle2::EventT0DQMModule::m_histEventT0_TOP_hadron_L1_CDCTRG
TH1F * m_histEventT0_TOP_hadron_L1_CDCTRG
event t0 histogram for TOP and hadronic events wrt the CDC trigger time
Definition: EventT0DQM.h:103
Belle2::SoftwareTriggerCutResult::c_accept
@ c_accept
Accept this event.
Belle2::EventT0DQMModule::event
virtual void event() override
This method is called for each event.
Definition: EventT0DQM.cc:155
Belle2::EventT0DQMModule::m_histEventT0_CDC_mumu_L1_ECLTRG
TH1F * m_histEventT0_CDC_mumu_L1_ECLTRG
event t0 histogram for CDC and mu mu events wrt the ECL trigger time
Definition: EventT0DQM.h:95
Belle2::EventT0DQMModule::m_histEventT0_ECL_hadron_L1_ECLTRG
TH1F * m_histEventT0_ECL_hadron_L1_ECLTRG
event t0 histogram for ECL and hadronic events wrt the ECL trigger time
Definition: EventT0DQM.h:91
Belle2::EventT0DQMModule::m_histEventT0_CDC_hadron_L1_CDCTRG
TH1F * m_histEventT0_CDC_hadron_L1_CDCTRG
event t0 histogram for CDC and hadronic events wrt the CDC trigger time
Definition: EventT0DQM.h:102
Belle2::EventT0DQMModule::m_histEventT0_ECL_hadron_L1_CDCTRG
TH1F * m_histEventT0_ECL_hadron_L1_CDCTRG
event t0 histogram for ECL and hadronic events wrt the CDC trigger time
Definition: EventT0DQM.h:101
Belle2::EventT0DQMModule::m_histEventT0_CDC_mumu_L1_CDCTRG
TH1F * m_histEventT0_CDC_mumu_L1_CDCTRG
event t0 histogram for CDC and mu mu events wrt the CDC trigger time
Definition: EventT0DQM.h:105
Belle2::EventT0DQMModule::m_histEventT0_TOP_bhabha_L1_CDCTRG
TH1F * m_histEventT0_TOP_bhabha_L1_CDCTRG
event t0 histogram for TOP and bha bha events wrt the CDC trigger time
Definition: EventT0DQM.h:100
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::EventT0DQMModule::endRun
virtual void endRun() override
This method is called at the end of each run.
Definition: EventT0DQM.cc:317
Belle2::EventT0DQMModule::m_histEventT0_TOP_hadron_L1_ECLTRG
TH1F * m_histEventT0_TOP_hadron_L1_ECLTRG
event t0 histogram for TOP and hadronic events wrt the ECL trigger time
Definition: EventT0DQM.h:93
Belle2::EventT0DQMModule::defineHisto
virtual void defineHisto() override
Defination of histograms.
Definition: EventT0DQM.cc:30
Belle2::EventT0DQMModule::m_histEventT0_CDC_bhabha_L1_ECLTRG
TH1F * m_histEventT0_CDC_bhabha_L1_ECLTRG
event t0 histogram for CDC and bha bha events wrt the ECL trigger time
Definition: EventT0DQM.h:89
Belle2::EventT0DQMModule::initialize
virtual void initialize() override
Initialize the module.
Definition: EventT0DQM.cc:106
Belle2::EventT0DQMModule::m_histEventT0_CDC_hadron_L1_ECLTRG
TH1F * m_histEventT0_CDC_hadron_L1_ECLTRG
event t0 histogram for CDC and hadronic events wrt the ECL trigger time
Definition: EventT0DQM.h:92
Belle2::EventT0DQMModule::m_histEventT0_TOP_mumu_L1_CDCTRG
TH1F * m_histEventT0_TOP_mumu_L1_CDCTRG
event t0 histogram for TOP and mu mu events wrt the CDC trigger time
Definition: EventT0DQM.h:106
Belle2::EventT0DQMModule::m_TrgResult
StoreObjPtr< SoftwareTriggerResult > m_TrgResult
Store array for Trigger selection.
Definition: EventT0DQM.h:85
Belle2::EventT0DQMModule::m_objTrgSummary
StoreObjPtr< TRGSummary > m_objTrgSummary
Trigger Summary data object.
Definition: EventT0DQM.h:76
Belle2::EventT0DQMModule::m_histEventT0_ECL_mumu_L1_CDCTRG
TH1F * m_histEventT0_ECL_mumu_L1_CDCTRG
event t0 histogram for ECL and mu mu events wrt the CDC trigger time
Definition: EventT0DQM.h:104
Belle2::EventT0DQMModule::m_histEventT0_TOP_mumu_L1_ECLTRG
TH1F * m_histEventT0_TOP_mumu_L1_ECLTRG
event t0 histogram for TOP and mu mu events wrt the ECL trigger time
Definition: EventT0DQM.h:96
Belle2::EventT0DQMModule::m_histEventT0_ECL_bhabha_L1_CDCTRG
TH1F * m_histEventT0_ECL_bhabha_L1_CDCTRG
event t0 histogram for ECL and bha bha events wrt the CDC trigger time
Definition: EventT0DQM.h:98
Belle2::HistoModule
HistoModule.h is supposed to be used instead of Module.h for the modules with histogram definitions t...
Definition: HistoModule.h:29
Belle2::EventT0DQMModule
This module to design collect the event t0 values base on different detectors and physics processes.
Definition: EventT0DQM.h:42