Belle II Software  release-05-01-25
PXDInjectionDQMModule.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2019 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Bjoern Spruck *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include <pxd/modules/pxdDQM/PXDInjectionDQMModule.h>
12 #include "TDirectory.h"
13 
14 using namespace std;
15 using namespace Belle2;
16 using namespace Belle2::PXD;
17 using namespace Belle2::VXD;
18 
19 //-----------------------------------------------------------------
20 // Register the Module
21 //-----------------------------------------------------------------
22 REG_MODULE(PXDInjectionDQM)
23 
24 //-----------------------------------------------------------------
25 // Implementation
26 //-----------------------------------------------------------------
27 
28 PXDInjectionDQMModule::PXDInjectionDQMModule() : HistoModule() , m_vxdGeometry(VXD::GeoCache::getInstance())
29 {
30  //Set module properties
31  setDescription("Monitor Occupancy after Injection");
32  setPropertyFlags(c_ParallelProcessingCertified);
33  addParam("histogramDirectoryName", m_histogramDirectoryName, "Name of the directory where histograms will be placed",
34  std::string("PXDINJ"));
35  addParam("PXDRawHitsName", m_PXDRawHitsName, "Name of PXD raw hits", std::string(""));
36  addParam("PXDClusterName", m_PXDClustersName, "Name of PXD clusters", std::string(""));
37  addParam("eachModule", m_eachModule, "create for each module", false);
38  addParam("offlineStudy", m_offlineStudy, "use finest binning and larger range", false);
39  addParam("useClusters", m_useClusters, "use cluster instead of raw hits", false);
40  addParam("createMaxHist", m_createMaxHist, "create histo with max occupancy (not mp save!!!)", false);
41 
42 }
43 
44 void PXDInjectionDQMModule::defineHisto()
45 {
46  TDirectory* oldDir = gDirectory;
47  if (m_histogramDirectoryName != "") {
48  oldDir->mkdir(m_histogramDirectoryName.c_str());// do not rely on return value, might be ZERO
49  oldDir->cd(m_histogramDirectoryName.c_str());//changing to the right directory
50  }
51 
52  if (m_offlineStudy) {
53  hOccAfterInjLER = new TH1F("PXDOccInjLER", "PXDOccInjLER/Time;Time in #mus;Count/Time (0.5 #mus bins)", 100000, 0, 50000);
54  hOccAfterInjHER = new TH1F("PXDOccInjHER", "PXDOccInjHER/Time;Time in #mus;Count/Time (0.5 #mus bins)", 100000, 0, 50000);
55  hEOccAfterInjLER = new TH1I("PXDEOccInjLER", "PXDEOccInjLER/Time;Time in #mus;Triggers/Time (0.5 #mus bins)", 100000, 0, 50000);
56  hEOccAfterInjHER = new TH1I("PXDEOccInjHER", "PXDEOccInjHER/Time;Time in #mus;Triggers/Time (0.5 #mus bins)", 100000, 0, 50000);
57  if (m_createMaxHist) {
58  hMaxOccAfterInjLER = new TH1F("PXDMaxOccInjLER", "PXDMaxOccInjLER/Time;Time in #mus;Triggers/Time (0.5 #mus bins)", 100000, 0,
59  50000);
60  hMaxOccAfterInjHER = new TH1F("PXDMaxOccInjHER", "PXDMaxOccInjHER/Time;Time in #mus;Triggers/Time (0.5 #mus bins)", 100000, 0,
61  50000);
62  }
63 
64  } else {
65  hOccAfterInjLER = new TH1F("PXDOccInjLER", "PXDOccInjLER/Time;Time in #mus;Count/Time (5 #mus bins)", 4000, 0, 20000);
66  hOccAfterInjHER = new TH1F("PXDOccInjHER", "PXDOccInjHER/Time;Time in #mus;Count/Time (5 #mus bins)", 4000, 0, 20000);
67  hEOccAfterInjLER = new TH1I("PXDEOccInjLER", "PXDEOccInjLER/Time;Time in #mus;Triggers/Time (5 #mus bins)", 4000, 0, 20000);
68  hEOccAfterInjHER = new TH1I("PXDEOccInjHER", "PXDEOccInjHER/Time;Time in #mus;Triggers/Time (5 #mus bins)", 4000, 0, 20000);
69  if (m_createMaxHist) {
70  hMaxOccAfterInjLER = new TH1F("PXDMaxOccInjLER", "PXDMaxOccInjLER/Time;Time in #mus;Triggers/Time (5 #mus bins)", 4000, 0, 20000);
71  hMaxOccAfterInjHER = new TH1F("PXDMaxOccInjHER", "PXDMaxOccInjHER/Time;Time in #mus;Triggers/Time (5 #mus bins)", 4000, 0, 20000);
72  }
73 
74  }
75 
76  if (m_eachModule) {
77  std::vector<VxdID> sensors = m_vxdGeometry.getListOfSensors();
78  for (VxdID& avxdid : sensors) {
79  VXD::SensorInfoBase info = m_vxdGeometry.getSensorInfo(avxdid);
80  if (info.getType() != VXD::SensorInfoBase::PXD) continue;
81  // Only interested in PXD sensors
82 
83  TString buff = (std::string)avxdid;
84  TString bufful = buff;
85  bufful.ReplaceAll(".", "_");
86 
89  if (m_offlineStudy) {
90  hOccModAfterInjLER[avxdid] = new TH1F("PXDOccInjLER_" + bufful,
91  "PXDOccModInjLER " + buff + "/Time;Time in #mus;Count/Time (0.5 #mus bins)", 100000, 0, 50000);
92  hOccModAfterInjHER[avxdid] = new TH1F("PXDOccInjHER_" + bufful,
93  "PXDOccModInjHER " + buff + "/Time;Time in #mus;Count/Time (0.5 #mus bins)", 100000, 0, 50000);
94  if (m_createMaxHist) {
95  hMaxOccModAfterInjLER[avxdid] = new TH1F("PXDMaxOccInjLER_" + bufful,
96  "PXDMaxOccModInjLER " + buff + "/Time;Time in #mus;Count/Time (0.5 #mus bins)", 100000, 0, 50000);
97  hMaxOccModAfterInjHER[avxdid] = new TH1F("PXDMaxOccInjHER_" + bufful,
98  "PXDMaxOccModInjHER " + buff + "/Time;Time in #mus;Count/Time (0.5 #mus bins)", 100000, 0, 50000);
99  }
100  } else {
101  hOccModAfterInjLER[avxdid] = new TH1F("PXDOccInjLER_" + bufful,
102  "PXDOccModInjLER " + buff + "/Time;Time in #mus;Count/Time (5 #mus bins)", 4000, 0, 20000);
103  hOccModAfterInjHER[avxdid] = new TH1F("PXDOccInjHER_" + bufful,
104  "PXDOccModInjHER " + buff + "/Time;Time in #mus;Count/Time (5 #mus bins)", 4000, 0, 20000);
105  if (m_createMaxHist) {
106  hMaxOccModAfterInjLER[avxdid] = new TH1F("PXDMaxOccInjLER_" + bufful,
107  "PXDMaxOccModInjLER " + buff + "/Time;Time in #mus;Count/Time (5 #mus bins)", 4000, 0, 20000);
108  hMaxOccModAfterInjHER[avxdid] = new TH1F("PXDMaxOccInjHER_" + bufful,
109  "PXDMaxOccModInjHER " + buff + "/Time;Time in #mus;Count/Time (5 #mus bins)", 4000, 0, 20000);
110  }
111  }
112  }
113  }
114 
115 // hTrigAfterInjLER = new TH2F("TrigAfterInjLER",
116 // "Triggers for LER veto tuning;Time since last injection in #mus;Time within beam cycle in #mus", 500, 0, 30000, 100, 0,
117 // 5120 / 508.);
118 // hTrigAfterInjHER = new TH2F("TrigAfterInjHER",
119 // "Triggers for HER veto tuning;Time since last injection in #mus;Time within beam cycle in #mus", 500, 0, 30000, 100, 0,
120 // 5120 / 508.);
121 
122  hTriggersAfterTrigger = new TH1I("PXDTriggersAfterLast",
123  "PXD Trigger after Last Trigger;Time diff in #mus;Count/Time (0.5 #mus bins)", 100000, 0, 50000);
124  hTriggersPerBunch = new TH1I("PXDTriggerBunch", "PXD Trigger per Bunch;Bunch/4;Triggers", 1280, 0, 1280);
125 
126 
127  // cd back to root directory
128  oldDir->cd();
129 }
130 
131 void PXDInjectionDQMModule::initialize()
132 {
133  REG_HISTOGRAM
134  m_rawTTD.isOptional();
135  if (m_useClusters) {
136  m_storeClusters.isRequired(m_PXDClustersName);
137  } else {
138  m_storeRawHits.isRequired(m_PXDRawHitsName);
139  }
140 }
141 
142 void PXDInjectionDQMModule::beginRun()
143 {
144  // Do not assume that everthing is non-zero, e.g. Max might be nullptr
145  if (hOccAfterInjLER) hOccAfterInjLER->Reset();
146  if (hOccAfterInjHER) hOccAfterInjHER->Reset();
147  if (hEOccAfterInjLER) hEOccAfterInjLER->Reset();
148  if (hEOccAfterInjHER) hEOccAfterInjHER->Reset();
149  if (hMaxOccAfterInjLER) hMaxOccAfterInjLER->Reset();
150  if (hMaxOccAfterInjHER) hMaxOccAfterInjHER->Reset();
151  for (auto& a : hOccModAfterInjLER) if (a.second) a.second->Reset();
152  for (auto& a : hOccModAfterInjHER) if (a.second) a.second->Reset();
153  for (auto& a : hMaxOccModAfterInjLER) if (a.second) a.second->Reset();
154  for (auto& a : hMaxOccModAfterInjHER) if (a.second) a.second->Reset();
155  for (auto& a : hOccModAfterInjLERGate) if (a.second) a.second->Reset();
156  for (auto& a : hOccModAfterInjHERGate) if (a.second) a.second->Reset();
157 // hTrigAfterInjLER->Reset();
158 // hTrigAfterInjHER->Reset();
159  hTriggersAfterTrigger->Reset();
160  hTriggersPerBunch->Reset();
161 }
162 
163 void PXDInjectionDQMModule::event()
164 {
165 
166  for (auto& it : m_rawTTD) {
167  B2DEBUG(29, "TTD FTSW : " << hex << it.GetTTUtime(0) << " " << it.GetTTCtime(0) << " EvtNr " << it.GetEveNo(0) << " Type " <<
168  (it.GetTTCtimeTRGType(0) & 0xF) << " TimeSincePrev " << it.GetTimeSincePrevTrigger(0) << " TimeSinceInj " <<
169  it.GetTimeSinceLastInjection(0) << " IsHER " << it.GetIsHER(0) << " Bunch " << it.GetBunchNumber(0));
170 
171  // get last injection time
172  hTriggersAfterTrigger->Fill(it.GetTimeSincePrevTrigger(0) / 127.);
173  hTriggersPerBunch->Fill(it.GetBunchNumber(0));
174 
175  auto difference = it.GetTimeSinceLastInjection(0);
176  // check time overflow, too long ago
177  if (difference != 0x7FFFFFFF) {
178  // count raw pixel hits or clusters per module, only if necessary
179  unsigned int all = 0;
180  std::map <VxdID, int> freq;// count the number of RawHits per sensor
181  if (m_useClusters) {
182  for (auto& p : m_storeClusters) {
183  freq[p.getSensorID()]++;
184  all++;
185  }
186  } else {
187  for (auto& p : m_storeRawHits) {
188  freq[p.getSensorID()]++;
189  all++;
190  }
191  }
192  float diff2 = difference / 127.; // 127MHz clock ticks to us, inexact rounding
193  if (it.GetIsHER(0)) {
194  hOccAfterInjHER->Fill(diff2, all);
195  hEOccAfterInjHER->Fill(diff2);
196 // hTrigAfterInjHER->Fill(diff2, diff2 - int(diff2 / (5120 / 508.)) * (5120 / 508.));
197  if (m_createMaxHist) {
198  auto bin = hMaxOccAfterInjHER->FindBin(diff2);
199  auto value = hMaxOccAfterInjHER->GetBinContent(bin);
200  if (all > value) hMaxOccAfterInjHER->SetBinContent(bin, all);
201  }
202  for (auto& a : hOccModAfterInjHER) {
203  if (a.second) a.second->Fill(diff2, freq[a.first]);
204  }
205  if (m_createMaxHist) {
206  for (auto& a : hMaxOccModAfterInjHER) {
207  if (a.second) {
208  auto bin = a.second->FindBin(diff2);
209  auto value = a.second->GetBinContent(bin);
210  if (freq[a.first] > value) a.second->SetBinContent(bin, freq[a.first]);
211  }
212  }
213  }
214  } else {
215  hOccAfterInjLER->Fill(diff2, all);
216  hEOccAfterInjLER->Fill(diff2);
217 // hTrigAfterInjLER->Fill(diff2, diff2 - int(diff2 / (5120 / 508.)) * (5120 / 508.));
218  if (m_createMaxHist) {
219  auto bin = hMaxOccAfterInjLER->FindBin(diff2);
220  auto value = hMaxOccAfterInjLER->GetBinContent(bin);
221  if (all > value) hMaxOccAfterInjLER->SetBinContent(bin, all);
222  }
223  for (auto& a : hOccModAfterInjLER) {
224  if (a.second) a.second->Fill(diff2, freq[a.first]);
225  }
226  if (m_createMaxHist) {
227  for (auto& a : hMaxOccModAfterInjLER) {
228  if (a.second) {
229  auto bin = a.second->FindBin(diff2);
230  auto value = a.second->GetBinContent(bin);
231  if (freq[a.first] > value) a.second->SetBinContent(bin, freq[a.first]);
232  }
233  }
234 
235  }
236  }
237  }
238 
239  break;
240  }
241 }
Belle2::VxdID
Class to uniquely identify a any structure of the PXD and SVD.
Definition: VxdID.h:43
Belle2::PXD::PXDInjectionDQMModule
The PXD Occupancy after Injection DQM module.
Definition: PXDInjectionDQMModule.h:47
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::VXD::SensorInfoBase
Base class to provide Sensor Information for PXD and SVD.
Definition: SensorInfoBase.h:40
Belle2::VXD
Namespace to provide code needed by both Vertex Detectors, PXD and SVD, and also testbeam telescopes.
Definition: GeoCache.h:36
Belle2::PXD
Namespace to encapsulate code needed for simulation and reconstrucion of the PXD.
Definition: PXDCalibrationUtilities.h:28
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::VXD::GeoCache
Class to faciliate easy access to sensor information of the VXD like coordinate transformations or pi...
Definition: GeoCache.h:41
Belle2::HistoModule
HistoModule.h is supposed to be used instead of Module.h for the modules with histogram definitions t...
Definition: HistoModule.h:29