Belle II Software  release-05-01-25
BeamBkgHitRateMonitorModule.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2019 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Marko Staric *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 // Own include
12 #include <background/modules/BeamBkgHitRateMonitor/BeamBkgHitRateMonitorModule.h>
13 #include <background/modules/BeamBkgHitRateMonitor/PXDHitRateCounter.h>
14 #include <background/modules/BeamBkgHitRateMonitor/SVDHitRateCounter.h>
15 #include <background/modules/BeamBkgHitRateMonitor/CDCHitRateCounter.h>
16 #include <background/modules/BeamBkgHitRateMonitor/TOPHitRateCounter.h>
17 #include <background/modules/BeamBkgHitRateMonitor/ARICHHitRateCounter.h>
18 #include <background/modules/BeamBkgHitRateMonitor/ECLHitRateCounter.h>
19 #include <background/modules/BeamBkgHitRateMonitor/KLMHitRateCounter.h>
20 
21 // framework aux
22 #include <framework/logging/Logger.h>
23 
24 #include <framework/io/RootIOUtilities.h>
25 #include <framework/core/RandomNumbers.h>
26 #include <framework/core/Environment.h>
27 #include <framework/core/ModuleParam.templateDetails.h>
28 #include <framework/database/Database.h>
29 #include <framework/utilities/EnvironmentVariables.h>
30 
31 #include <boost/python.hpp>
32 #include <boost/optional.hpp>
33 #include <boost/filesystem/path.hpp>
34 #include <boost/filesystem/operations.hpp>
35 #include <boost/algorithm/string.hpp>
36 
37 using namespace std;
38 
39 namespace Belle2 {
45  //-----------------------------------------------------------------
46  // Register module
47  //-----------------------------------------------------------------
48 
49  REG_MODULE(BeamBkgHitRateMonitor)
50 
51  //-----------------------------------------------------------------
52  // Implementation
53  //-----------------------------------------------------------------
54 
56 
57  {
58  // set module description
59  setDescription("A module for off-line monitoring of beam background hit rates.");
60 
61  /* No multiprocessing allowed!
62  * setPropertyFlags(c_ParallelProcessingCertified);
63  */
64 
65  // Add parameters
66  addParam("outputFileName", m_outputFileName, "output file name",
67  string("beamBkgHitRates.root"));
68  addParam("treeName", m_treeName, "output tree name",
69  string("tree"));
70  m_trgTypes.push_back(TRGSummary::TTYP_DPHY);
71  m_trgTypes.push_back(TRGSummary::TTYP_RAND);
72  addParam("trgTypes", m_trgTypes,
73  "trigger types for event selection (see TRGSummary.h for definitions). "
74  "Empty list means all trigger types.",
75  m_trgTypes);
76  addParam("writeEmptyTimeStamps", m_writeEmptyTimeStamps,
77  "if true, write to ntuple also empty time stamps", false);
78  addParam("topTimeOffset", m_topTimeOffset,
79  "TOP: time offset of hits (to be subtracted) [ns]", 25.0);
80  addParam("topTimeWindow", m_topTimeWindow,
81  "TOP: time window in which to count hits [ns]", 100.0);
82  addParam("svdShaperDigitsName", m_svdShaperDigitsName,
83  "SVDShaperDigits collection name", string(""));
84  addParam("svdThrCharge", m_svdThrCharge,
85  "Energy cur on SVD Cluster charge in electrons", 15000.);
86  addParam("additionalDataDescription", m_additionalDataDescription,
87  "Additional dictionary of "
88  "name->value pairs to be added to the file metadata to describe the data",
89  m_additionalDataDescription);
90  addParam("cdcTimeWindowLowerEdgeSmallCell", m_cdcTimeWindowLowerEdgeSmallCell,
91  "CDC: lower edge of the time window for small cells [tdc count = ns]",
92  4550);
93  addParam("cdcTimeWindowUpperEdgeSmallCell", m_cdcTimeWindowUpperEdgeSmallCell,
94  "CDC: upper edge of the time window for small cells [tdc count = ns]",
95  5050);
96  addParam("cdcTimeWindowLowerEdgeNormalCell", m_cdcTimeWindowLowerEdgeNormalCell,
97  "CDC: lower edge of the time window for normal cells [tdc count = ns]",
98  4200);
99  addParam("cdcTimeWindowUpperEdgeNormalCell", m_cdcTimeWindowUpperEdgeNormalCell,
100  "CDC: upper edge of the time window for normal cells [tdc count = ns]",
101  5050);
102  addParam("cdcEnableBadWireTreatment", m_cdcEnableBadWireTreatment,
103  "CDC: flag to enable the bad wire treatment", true);
104  addParam("cdcEnableBackgroundHitFilter", m_cdcEnableBackgroundHitFilter,
105  "CDC: flag to enable the CDC background hit (crosstakl, noise) filter", true);
106  addParam("cdcEnableMarkBackgroundHit", m_cdcEnableMarkBackgroundHit,
107  "CDC: flag to enable to mark background flag on CDCHit (set 0x100 bit for CDCHit::m_status).", false);
108 
109  }
110 
111  BeamBkgHitRateMonitorModule::~BeamBkgHitRateMonitorModule()
112  {
113  for (auto& monitor : m_monitors) {
114  if (monitor) delete monitor;
115  }
116  }
117 
118  void BeamBkgHitRateMonitorModule::initialize()
119  {
120  // collections registration
121  m_eventMetaData.isRequired();
122  if (m_trgTypes.empty()) {
123  m_trgSummary.isOptional(); // enables to run the module when TRGSummary is absent
124  } else {
125  m_trgSummary.isRequired();
126  }
127  m_fileMetaData.isOptional(); // enables to run the module with simulation
128 
129  // create, set and append hit rate monitoring classes
130  auto* pxd = new Background::PXDHitRateCounter();
131  m_monitors.push_back(pxd);
132  auto* svd = new Background::SVDHitRateCounter(m_svdShaperDigitsName, m_svdThrCharge);
133  m_monitors.push_back(svd);
134  auto* cdc = new Background::CDCHitRateCounter(m_cdcTimeWindowLowerEdgeSmallCell, m_cdcTimeWindowUpperEdgeSmallCell,
135  m_cdcTimeWindowLowerEdgeNormalCell, m_cdcTimeWindowUpperEdgeNormalCell,
136  m_cdcEnableBadWireTreatment, m_cdcEnableBackgroundHitFilter,
137  m_cdcEnableMarkBackgroundHit);
138  m_monitors.push_back(cdc);
139  auto* top = new Background::TOPHitRateCounter(m_topTimeOffset, m_topTimeWindow);
140  m_monitors.push_back(top);
141  auto* arich = new Background::ARICHHitRateCounter();
142  m_monitors.push_back(arich);
143  auto* ecl = new Background::ECLHitRateCounter();
144  m_monitors.push_back(ecl);
145  auto* klm = new Background::KLMHitRateCounter();
146  m_monitors.push_back(klm);
147 
148  // open output root file
149  m_file = TFile::Open(m_outputFileName.c_str(), "RECREATE");
150  if (not m_file) {
151  B2FATAL("Cannot open output file '" << m_outputFileName << "' for writing");
152  }
153 
154  // create tree
155  m_tree = new TTree(m_treeName.c_str(), "hit rates of selected events");
156 
157  // create persistent tree to store fileMetaData
158  m_persistent = new TTree("persistent", "persistent data");
159  m_persistent->Branch("FileMetaData", &m_outputFileMetaData);
160 
161  // set tree branches
162  m_tree->Branch("run", &m_run, "run/I");
163  m_tree->Branch("numEvents", &m_numEvents, "numEvents/I");
164  m_tree->Branch("timeStamp", &m_timeStamp, "timeStamp/i");
165  m_tree->Branch("time", &m_time, "time/I");
166  for (auto& monitor : m_monitors) {
167  monitor->initialize(m_tree);
168  }
169 
170  // control histograms
171  m_trgAll = new TH1F("trgAll", "trigger types of all events", 16, -0.5, 15.5);
172  m_trgAll->SetXTitle("type of trigger timing source");
173  m_trgSel = new TH1F("trgSel", "trigger types of selected events", 16, -0.5, 15.5);
174  m_trgSel->SetXTitle("type of trigger timing source");
175 
176  }
177 
178  void BeamBkgHitRateMonitorModule::beginRun()
179  {
180  // clear buffers
181  for (auto& monitor : m_monitors) {
182  monitor->clear();
183  }
184  m_eventCounts.clear();
185 
186  // clear counters
187  m_numEventsSelected = 0;
188  m_trgTypesCount.clear();
189 
190  // set run number
191  m_run = m_eventMetaData->getRun();
192 
193  // set unix time of the first event in the run
194  unsigned utime = m_eventMetaData->getTime() / 1000000000;
195  m_utimeFirst = utime;
196  m_utimeMin = utime;
197  m_utimeMax = utime + 1;
198 
199  }
200 
201  void BeamBkgHitRateMonitorModule::event()
202  {
203  // get unix time of the event
204  unsigned utime = m_eventMetaData->getTime() / 1000000000;
205  m_utimeMin = std::min(m_utimeMin, utime);
206  m_utimeMax = std::max(m_utimeMax, utime + 1);
207 
208  // collect file meta data
209  collectFileMetaData();
210 
211  // event selection
212  if (not isEventSelected()) return;
213  m_numEventsSelected++;
214 
215  // accumulate
216  for (auto& monitor : m_monitors) {
217  monitor->accumulate(utime);
218  }
219  m_eventCounts[utime] += 1;
220 
221  }
222 
223  void BeamBkgHitRateMonitorModule::endRun()
224  {
225  // fill ntuple
226  for (unsigned utime = m_utimeMin; utime < m_utimeMax; utime++) {
227  if (not m_writeEmptyTimeStamps) {
228  if (m_eventCounts.find(utime) == m_eventCounts.end()) continue;
229  }
230  m_numEvents = m_eventCounts[utime];
231  m_timeStamp = utime;
232  m_time = utime - m_utimeMin;
233  for (auto& monitor : m_monitors) {
234  monitor->normalize(utime);
235  }
236  m_tree->Fill();
237  }
238 
239  // count selected events in all runs
240  m_allEventsSelected += m_numEventsSelected;
241 
242  // print a summary for this run
243  std::string trigs;
244  for (const auto& trgType : m_trgTypesCount) {
245  trigs += " trigger type " + std::to_string(trgType.first) + ": " +
246  std::to_string(trgType.second) + " events\n";
247  }
248  B2INFO("Run " << m_run << ": " << m_numEventsSelected
249  << " events selected for beam background hit rate monitoring.\n"
250  << trigs
251  << LogVar("first event utime ", m_utimeMin)
252  << LogVar("start utime ", m_utimeMin)
253  << LogVar("stop utime ", m_utimeMax)
254  << LogVar("duration [seconds]", m_utimeMax - m_utimeMin)
255  );
256  }
257 
258  void BeamBkgHitRateMonitorModule::terminate()
259  {
260  setFileMetaData();
261  m_persistent->Fill();
262 
263  // write to file and close
264  m_file->cd();
265  m_file->Write();
266  m_file->Close();
267 
268  B2INFO("Output file: " << m_outputFileName);
269  }
270 
271  bool BeamBkgHitRateMonitorModule::isEventSelected()
272  {
273  auto trgType = TRGSummary::TTYP_NONE;
274  if (m_trgSummary.isValid()) trgType = m_trgSummary->getTimType();
275  m_trgAll->Fill(trgType);
276 
277  if (m_trgTypes.empty()) {
278  m_trgTypesCount[trgType] += 1;
279  m_trgSel->Fill(trgType);
280  return true;
281  }
282  for (auto type : m_trgTypes) {
283  if (trgType == type) {
284  m_trgTypesCount[trgType] += 1;
285  m_trgSel->Fill(trgType);
286  return true;
287  }
288  }
289  return false;
290  }
291 
292 
293  void BeamBkgHitRateMonitorModule::collectFileMetaData()
294  {
295  // add file name to the list
296  if (m_fileMetaData.isValid()) {
297  std::string lfn = m_fileMetaData->getLfn();
298  if (not lfn.empty() and (m_parentLfns.empty() or (m_parentLfns.back() != lfn))) {
299  m_parentLfns.push_back(lfn);
300  }
301  }
302 
303  // low and high experiment, run and event numbers
304  unsigned long experiment = m_eventMetaData->getExperiment();
305  unsigned long run = m_eventMetaData->getRun();
306  unsigned long event = m_eventMetaData->getEvent();
307  if (m_experimentLow > m_experimentHigh) { //starting condition
308  m_experimentLow = m_experimentHigh = experiment;
309  m_runLow = m_runHigh = run;
310  m_eventLow = m_eventHigh = event;
311  } else {
312  if ((experiment < m_experimentLow) or ((experiment == m_experimentLow) and ((run < m_runLow) or ((run == m_runLow)
313  and (event < m_eventLow))))) {
314  m_experimentLow = experiment;
315  m_runLow = run;
316  m_eventLow = event;
317  }
318  if ((experiment > m_experimentHigh) or ((experiment == m_experimentHigh) and ((run > m_runHigh) or ((run == m_runHigh)
319  and (event > m_eventHigh))))) {
320  m_experimentHigh = experiment;
321  m_runHigh = run;
322  m_eventHigh = event;
323  }
324  }
325 
326  }
327 
328 
329  void BeamBkgHitRateMonitorModule::setFileMetaData()
330  {
331 
332  if (m_fileMetaData.isValid() and not m_fileMetaData->isMC()) {
333  m_outputFileMetaData.declareRealData();
334  }
335 
336  m_outputFileMetaData.setNEvents(m_allEventsSelected);
337 
338  if (m_experimentLow > m_experimentHigh) {
339  // starting condition so apparently no events at all
340  m_outputFileMetaData.setLow(-1, -1, 0);
341  m_outputFileMetaData.setHigh(-1, -1, 0);
342  } else {
343  m_outputFileMetaData.setLow(m_experimentLow, m_runLow, m_eventLow);
344  m_outputFileMetaData.setHigh(m_experimentHigh, m_runHigh, m_eventHigh);
345  }
346 
347  m_outputFileMetaData.setParents(m_parentLfns);
348  RootIOUtilities::setCreationData(m_outputFileMetaData);
349  m_outputFileMetaData.setRandomSeed(RandomNumbers::getSeed());
350  m_outputFileMetaData.setSteering(Environment::Instance().getSteering());
351  auto mcEvents = Environment::Instance().getNumberOfMCEvents();
352  m_outputFileMetaData.setMcEvents(mcEvents);
353  m_outputFileMetaData.setDatabaseGlobalTag(Database::Instance().getGlobalTags());
354 
355  for (const auto& item : m_additionalDataDescription) {
356  m_outputFileMetaData.setDataDescription(item.first, item.second);
357  }
358 
359  std::string lfn = m_file->GetName();
360  lfn = boost::filesystem::absolute(lfn, boost::filesystem::initial_path()).string();
361  std::string format = EnvironmentVariables::get("BELLE2_LFN_FORMATSTRING", "");
362  if (!format.empty()) {
363  auto format_filename = boost::python::import("B2Tools.format").attr("format_filename");
364  lfn = boost::python::extract<std::string>(format_filename(format, m_outputFileName, m_outputFileMetaData.getJsonStr()));
365  }
366  m_outputFileMetaData.setLfn(lfn);
367 
368  }
369 
370 
372 } // end Belle2 namespace
373 
Belle2::Background::ARICHHitRateCounter
Class for monitoring beam background hit rates of ARICH.
Definition: ARICHHitRateCounter.h:41
Belle2::Background::KLMHitRateCounter
Class for monitoring beam background hit rates of EKLM.
Definition: KLMHitRateCounter.h:47
Belle2::Background::SVDHitRateCounter
Class for monitoring beam background hit rates of SVD.
Definition: SVDHitRateCounter.h:43
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::Background::PXDHitRateCounter
Class for monitoring beam background hit rates of PXD.
Definition: PXDHitRateCounter.h:43
Belle2::Module
Base class for Modules.
Definition: Module.h:74
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
LogVar
Class to store variables with their name which were sent to the logging service.
Definition: LogVariableStream.h:24
Belle2::BeamBkgHitRateMonitorModule
A module to monitor detector hit rates of beam background Output is to a flat ntuple.
Definition: BeamBkgHitRateMonitorModule.h:45
Belle2::Background::ECLHitRateCounter
Class for monitoring beam background hit rates of ECL.
Definition: ECLHitRateCounter.h:38
Belle2::Background::CDCHitRateCounter
Class for monitoring beam background hit rates of CDC.
Definition: CDCHitRateCounter.h:38
Belle2::Background::TOPHitRateCounter
Class for monitoring beam background hit rates of TOP.
Definition: TOPHitRateCounter.h:43