Belle II Software  release-05-02-19
DQMHistAnalysisSVDDose.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2021 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Ludovico Massaccesi *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include <dqm/analysis/modules/DQMHistAnalysisSVDDose.h>
12 #include <framework/utilities/Utils.h>
13 #include <TROOT.h>
14 #include <TText.h>
15 
16 using namespace std;
17 using namespace Belle2;
18 
19 // Utility function
20 inline double getClockSeconds() { return Utils::getClock() / 1e9; }
21 
22 REG_MODULE(DQMHistAnalysisSVDDose)
23 
24 DQMHistAnalysisSVDDoseModule::DQMHistAnalysisSVDDoseModule()
25 {
26  setDescription("Monitoring of SVD Dose with events from Poisson trigger w/o inj. veto. See also SVDDQMDoseModule.");
27  // THIS MODULE CAN NOT BE RUN IN PARALLEL
28  addParam("pvPrefix", m_pvPrefix, "Prefix for EPICS PVs.", std::string("DQM:SVD:"));
29  addParam("useEpics", m_useEpics, "Whether to update EPICS PVs.", true);
30  addParam("epicsUpdateSeconds", m_epicsUpdateSeconds,
31  "Minimum interval between two successive PV updates (in seconds).", 1000.0);
32  addParam("pvSuffix", m_pvSuffix, "Suffix for EPICS PVs.", std::string(":Occ:Pois:Avg"));
33  addParam("deltaTPVSuffix", m_deltaTPVSuffix, "Suffix for the PV that monitors the update interval of the PVs.",
34  std::string("Occ:Pois:UpdateInterval"));
35  addParam("statePVSuffix", m_statePVSuffix, "Suffix for the PV with the state of the monitoring.",
36  std::string("Occ:Pois:State"));
37 }
38 
39 DQMHistAnalysisSVDDoseModule::~DQMHistAnalysisSVDDoseModule()
40 {
41 #ifdef _BELLE2_EPICS
42  if (m_useEpics && ca_current_context()) ca_context_destroy();
43 #endif
44 }
45 
46 void DQMHistAnalysisSVDDoseModule::initialize()
47 {
48  B2DEBUG(18, "DQMHistAnalysisSVDDose: initialize");
49 
50  gROOT->cd(); // Don't know why I need this, but DQMHistAnalysisSVDOnMiraBelle uses it
51 
52  m_monObj = getMonitoringObject("svd"); // To write to MiraBelle
53 
54  m_c_instOccu.reserve(c_sensorGroups.size());
55  m_c_occuLER.reserve(c_sensorGroups.size());
56  m_c_occuHER.reserve(c_sensorGroups.size());
57  m_c_occuLER1.reserve(c_sensorGroups.size());
58  m_c_occuHER1.reserve(c_sensorGroups.size());
59  m_c_instOccuAll.reserve(c_sensorGroups.size());
60  m_c_occuLERAll.reserve(c_sensorGroups.size());
61  m_c_occuHERAll.reserve(c_sensorGroups.size());
62  m_c_occuLER1All.reserve(c_sensorGroups.size());
63  m_c_occuHER1All.reserve(c_sensorGroups.size());
64  for (const auto& group : c_sensorGroups) {
65  TCanvas* c = new TCanvas("svd_instOccupancy_" + group.nameSuffix + "_pois",
66  "Instantaneous occupancy (Pois. trig.) " + group.titleSuffix,
67  0, 0, 800, 600);
68  m_c_instOccu.push_back(c);
69  m_monObj->addCanvas(c);
70  c = new TCanvas("svd_occuLER_" + group.nameSuffix + "_pois",
71  "Occupancy vs time since LER inj. (Pois. trig.) " + group.titleSuffix,
72  0, 0, 800, 600);
73  m_c_occuLER.push_back(c);
74  // m_monObj->addCanvas(c);
75  c = new TCanvas("svd_occuHER_" + group.nameSuffix + "_pois",
76  "Occupancy vs time since HER inj. (Pois. trig.) " + group.titleSuffix,
77  0, 0, 800, 600);
78  m_c_occuHER.push_back(c);
79  // m_monObj->addCanvas(c);
80  c = new TCanvas("svd_1DoccuLER_" + group.nameSuffix + "_pois",
81  "Occupancy vs time since LER inj. (Pois. trig.) " + group.titleSuffix,
82  0, 0, 800, 600);
83  m_c_occuLER1.push_back(c);
84  m_monObj->addCanvas(c);
85  c = new TCanvas("svd_1DoccuHER_" + group.nameSuffix + "_pois",
86  "Occupancy vs time since HER inj. (Pois. trig.) " + group.titleSuffix,
87  0, 0, 800, 600);
88  m_c_occuHER1.push_back(c);
89  m_monObj->addCanvas(c);
90 
91  c = new TCanvas("svd_instOccupancy_" + group.nameSuffix + "_all",
92  "Instantaneous occupancy (all events) " + group.titleSuffix,
93  0, 0, 800, 600);
94  m_c_instOccuAll.push_back(c);
95  c = new TCanvas("svd_occuLER_" + group.nameSuffix + "_all",
96  "Occupancy vs time since LER inj. (all events) " + group.titleSuffix,
97  0, 0, 800, 600);
98  m_c_occuLERAll.push_back(c);
99  c = new TCanvas("svd_occuHER_" + group.nameSuffix + "_all",
100  "Occupancy vs time since HER inj. (all events) " + group.titleSuffix,
101  0, 0, 800, 600);
102  m_c_occuHERAll.push_back(c);
103  c = new TCanvas("svd_1DoccuLER_" + group.nameSuffix + "_all",
104  "Occupancy vs time since LER inj. (all events) " + group.titleSuffix,
105  0, 0, 800, 600);
106  m_c_occuLER1All.push_back(c);
107  c = new TCanvas("svd_1DoccuHER_" + group.nameSuffix + "_all",
108  "Occupancy vs time since HER inj. (all events) " + group.titleSuffix,
109  0, 0, 800, 600);
110  m_c_occuHER1All.push_back(c);
111  }
112  m_h_occuLER.resize(c_sensorGroups.size(), nullptr);
113  m_h_occuHER.resize(c_sensorGroups.size(), nullptr);
114  m_h_occuLER1.resize(c_sensorGroups.size(), nullptr);
115  m_h_occuHER1.resize(c_sensorGroups.size(), nullptr);
116  m_h_occuLERAll.resize(c_sensorGroups.size(), nullptr);
117  m_h_occuHERAll.resize(c_sensorGroups.size(), nullptr);
118  m_h_occuLER1All.resize(c_sensorGroups.size(), nullptr);
119  m_h_occuHER1All.resize(c_sensorGroups.size(), nullptr);
120 
121  // The legend need to be memory-leaked, so we make it once and use it evey time
122  m_legend = new TPaveText(0.53, 0.73, 0.68, 0.88, "brNDC");
123  m_legend->AddText("LER inj."); ((TText*)m_legend->GetListOfLines()->Last())->SetTextColor(kRed);
124  m_legend->AddText("HER inj."); ((TText*)m_legend->GetListOfLines()->Last())->SetTextColor(kAzure);
125  m_legend->AddText("No inj."); ((TText*)m_legend->GetListOfLines()->Last())->SetTextColor(kBlack);
126 
127 #ifdef _BELLE2_EPICS
128  if (m_useEpics) {
129  if (!ca_current_context())
130  SEVCHK(ca_context_create(ca_disable_preemptive_callback), "ca_context_create");
131  // Channels for the occupancies
132  m_myPVs.resize(c_sensorGroups.size());
133  for (unsigned int g = 0; g < c_sensorGroups.size(); g++)
134  SEVCHK(ca_create_channel((m_pvPrefix + c_sensorGroups[g].pvMiddle + m_pvSuffix).data(),
135  NULL, NULL, 10, &m_myPVs[g].mychid), "ca_create_channel");
136  // Channels for update interval and state
137  SEVCHK(ca_create_channel((m_pvPrefix + m_deltaTPVSuffix).data(),
138  NULL, NULL, 10, &m_timeSinceLastPVUpdateChan), "ca_create_channel");
139  SEVCHK(ca_create_channel((m_pvPrefix + m_statePVSuffix).data(),
140  NULL, NULL, 10, &m_stateChan), "ca_create_channel");
141  // Actually do create the channels, communicating with the IOC
142  SEVCHK(ca_pend_io(5.0), "ca_pend_io");
143  // Get the state enum
144  SEVCHK(ca_get(DBR_CTRL_ENUM, m_stateChan, &m_stateCtrl), "ca_get");
145  SEVCHK(ca_pend_io(5.0), "ca_pend_io");
146  B2DEBUG(19, "State PV initialized (ca_get)" << LogVar("value", m_stateCtrl.value));
147  // First update should happen m_epicsUpdateSeconds from now
148  m_lastPVUpdate = getClockSeconds();
149  }
150 #endif
151 }
152 
153 void DQMHistAnalysisSVDDoseModule::beginRun()
154 {
155  // Set status PV to running, reset last update time
156 #ifdef _BELLE2_EPICS
157  if (m_useEpics) {
158  B2DEBUG(19, "beginRun: setting state PV to RUNNING");
159  m_stateCtrl.value = 1;
160  SEVCHK(ca_put(DBR_ENUM, m_stateChan, &m_stateCtrl.value), "ca_put");
161  SEVCHK(ca_pend_io(5.0), "ca_pend_io");
162  // First update should happen m_epicsUpdateSeconds from now
163  m_lastPVUpdate = getClockSeconds();
164  }
165 #endif
166 }
167 
168 void DQMHistAnalysisSVDDoseModule::event()
169 {
170  // Update PVs ("write" to EPICS)
171 #ifdef _BELLE2_EPICS
172  double timeSinceLastPVUpdate = getClockSeconds() - m_lastPVUpdate;
173  if (m_useEpics && timeSinceLastPVUpdate >= m_epicsUpdateSeconds) {
174  for (unsigned int g = 0; g < c_sensorGroups.size(); g++) {
175  const auto& group = c_sensorGroups[g];
176  double nHits = 0.0, nEvts = 0.0;
177  for (TString dir : {"SVDDoseLERInjPois", "SVDDoseHERInjPois", "SVDDoseNoInjPois"}) {
178  auto hHits = findHistT<TH2F>(dir + "/SVDHitsVsTime_" + group.nameSuffix);
179  auto hEvts = findHistT<TH2F>(dir + "/SVDEvtsVsTime");
180  if (!hHits || !hEvts) {
181  B2WARNING("Histograms needed for Average Poisson Occupancy U-side not found.");
182  nEvts = 0.0;
183  break;
184  }
185  nHits += hHits->GetEntries();
186  nEvts += hEvts->GetEntries();
187  }
188 
189  B2DEBUG(19, "DQMHistAnalysisSVDDose: PV write"
190  << LogVar("group", group.nameSuffix.Data())
191  << LogVar("nEvts", nEvts) << LogVar("nHits", nHits));
192 
193  auto& pv = m_myPVs[g];
194  double delta_nHits = nHits - pv.lastNHits;
195  double delta_nEvts = nEvts - pv.lastNEvts;
196  double occ = delta_nEvts > 0.0 ? (delta_nHits / delta_nEvts * 100.0 / group.nStrips) : -1.0;
197  if (pv.mychid)
198  SEVCHK(ca_put(DBR_DOUBLE, pv.mychid, (void*)&occ), "ca_put");
199 
200  pv.lastNEvts = nEvts;
201  pv.lastNHits = nHits;
202  }
203  if (m_timeSinceLastPVUpdateChan)
204  SEVCHK(ca_put(DBR_DOUBLE, m_timeSinceLastPVUpdateChan, (void*)&timeSinceLastPVUpdate), "ca_put");
205  SEVCHK(ca_pend_io(5.0), "ca_pend_io");
206  m_lastPVUpdate = getClockSeconds();
207  }
208 #endif
209 
210  updateCanvases();
211 }
212 
213 void DQMHistAnalysisSVDDoseModule::endRun()
214 {
215  B2DEBUG(18, "DQMHistAnalysisSVDDose: endRun");
216 
217  // EPICS: reset the counters used for the delta computation, set state to NOT RUNNING
218 #ifdef _BELLE2_EPICS
219  if (m_useEpics) {
220  B2DEBUG(19, "endRun: setting state PV to NOT RUNNING");
221  m_stateCtrl.value = 0;
222  SEVCHK(ca_put(DBR_ENUM, m_stateChan, &m_stateCtrl.value), "ca_put");
223  SEVCHK(ca_pend_io(5.0), "ca_pend_io");
224  // Reset events and hits counters
225  for (auto& pv : m_myPVs)
226  pv.lastNEvts = pv.lastNHits = 0.0;
227  }
228 #endif
229 
230  // Write to MiraBelle
231  for (unsigned int g = 0; g < c_sensorGroups.size(); g++) {
232  const auto& group = c_sensorGroups[g];
233  double nHits = 0.0, nEvts = 0.0;
234  for (TString dir : {"SVDDoseLERInjPois", "SVDDoseHERInjPois", "SVDDoseNoInjPois"}) {
235  auto hHits = findHistT<TH2F>(dir + "/SVDHitsVsTime_" + group.nameSuffix);
236  auto hEvts = findHistT<TH2F>(dir + "/SVDEvtsVsTime");
237  if (!hHits || !hEvts) {
238  B2WARNING("Histograms needed for Average Poisson Occupancy U-side not found.");
239  nEvts = 0.0;
240  break;
241  }
242  nHits += hHits->GetEntries();
243  nEvts += hEvts->GetEntries();
244  }
245 
246  B2DEBUG(19, "DQMHistAnalysisSVDDose: MonObj write"
247  << LogVar("group", group.nameSuffix.Data())
248  << LogVar("nEvts", nEvts) << LogVar("nHits", nHits));
249 
250  double occ = nEvts ? (nHits / nEvts * 100.0 / group.nStrips) : -1.0;
251  TString vName = group.nameSuffix + "OccPoisAvg"; // e.g. L3XXUOccPoisAvg
252  m_monObj->setVariable(vName.Data(), occ);
253  }
254 
255  updateCanvases();
256 }
257 
258 void DQMHistAnalysisSVDDoseModule::updateCanvases()
259 {
260  B2DEBUG(18, "DQMHistAnalysisSVDDose: updating canvases");
261 
262  for (unsigned int g = 0; g < c_sensorGroups.size(); g++) {
263  const auto& group = c_sensorGroups[g];
264 
265  auto c = m_c_instOccu[g];
266  auto hLER = findHistT<TH1F>("SVDDoseLERInjPois/SVDInstOccu_" + group.nameSuffix);
267  auto hHER = findHistT<TH1F>("SVDDoseHERInjPois/SVDInstOccu_" + group.nameSuffix);
268  auto hNo = findHistT<TH1F>("SVDDoseNoInjPois/SVDInstOccu_" + group.nameSuffix);
269  if (hLER && hHER && hNo) {
270  hLER->SetLineColor(kRed);
271  hHER->SetLineColor(kAzure);
272  hNo->SetLineColor(kBlack);
273  carryOverflowOver(hLER);
274  carryOverflowOver(hHER);
275  carryOverflowOver(hNo);
276  c->Clear();
277  c->cd(0);
278  hNo->SetTitle("SVD instantaneous occu. " + group.titleSuffix + " U-side Pois. trig.");
279  hNo->Draw("hist"); // hNo usually has the larger maximum by far
280  hLER->Draw("hist same");
281  hHER->Draw("hist same");
282  c->SetLogy();
283  m_legend->Draw();
284  }
285 
286  c = m_c_occuLER[g];
287  auto hHits = findHistT<TH2F>("SVDDoseLERInjPois/SVDHitsVsTime_" + group.nameSuffix);
288  auto hEvts = findHistT<TH2F>("SVDDoseLERInjPois/SVDEvtsVsTime");
289  if (hHits && hEvts) {
290  if (m_h_occuLER[g]) delete m_h_occuLER[g];
291  auto hOccu = m_h_occuLER[g] = divide(hHits, hEvts, 100.0f / group.nStrips);
292  hOccu->SetTitle("SVD Occupancy " + group.titleSuffix + " - LER inj. Pois. trig."
293  ";Time since last injection [#mus];Time in beam cycle [#mus]"
294  ";Occupancy [%]");
295  hOccu->SetMinimum(1e-3);
296  hOccu->SetMaximum(10);
297  c->Clear();
298  c->cd(0);
299  c->SetRightMargin(0.16); // For the colorbar
300  hOccu->Draw("COLZ");
301  c->SetLogz();
302  }
303 
304  c = m_c_occuLER1[g];
305  auto hpEvts = findHistT<TH1F>("SVDDoseLERInjPois/SVDEvtsVsTime1");
306  auto hpHits = findHistT<TH1F>("SVDDoseLERInjPois/SVDHitsVsTime1_" + group.nameSuffix);
307  if (hpHits && hpEvts) {
308  if (m_h_occuLER1[g]) delete m_h_occuLER1[g];
309  auto hpOccu = m_h_occuLER1[g] = divide(hpHits, hpEvts, 100.0f / group.nStrips);
310  hpOccu->SetTitle("SVD Occupancy " + group.titleSuffix + " - LER inj. Pois. trig."
311  ";Time since last injection [#mus];Occupancy [%]");
312  hpOccu->SetMinimum(1e-3);
313  hpOccu->SetMaximum(10);
314  c->Clear();
315  c->cd(0);
316  hpOccu->SetMarkerStyle(7);
317  hpOccu->Draw("hist P");
318  c->SetLogy();
319  }
320 
321  c = m_c_occuHER[g];
322  hHits = findHistT<TH2F>("SVDDoseHERInjPois/SVDHitsVsTime_" + group.nameSuffix);
323  hEvts = findHistT<TH2F>("SVDDoseHERInjPois/SVDEvtsVsTime");
324  if (hHits && hEvts) {
325  if (m_h_occuHER[g]) delete m_h_occuHER[g];
326  auto hOccu = m_h_occuHER[g] = divide(hHits, hEvts, 100.0f / group.nStrips);
327  hOccu->SetTitle("SVD Occupancy " + group.titleSuffix + " - HER inj. Pois. trig."
328  ";Time since last injection [#mus];Time in beam cycle [#mus]"
329  ";Occupancy [%]");
330  hOccu->SetMinimum(1e-3);
331  hOccu->SetMaximum(10);
332  c->Clear();
333  c->cd(0);
334  c->SetRightMargin(0.16); // For the colorbar
335  hOccu->Draw("COLZ");
336  c->SetLogz();
337  }
338 
339  c = m_c_occuHER1[g];
340  hpEvts = findHistT<TH1F>("SVDDoseHERInjPois/SVDEvtsVsTime1");
341  hpHits = findHistT<TH1F>("SVDDoseHERInjPois/SVDHitsVsTime1_" + group.nameSuffix);
342  if (hpHits && hpEvts) {
343  if (m_h_occuHER1[g]) delete m_h_occuHER1[g];
344  auto hpOccu = m_h_occuHER1[g] = divide(hpHits, hpEvts, 100.0f / group.nStrips);
345  hpOccu->SetTitle("SVD Occupancy " + group.titleSuffix + " - HER inj. Pois. trig."
346  ";Time since last injection [#mus];Occupancy [%]");
347  hpOccu->SetMinimum(1e-3);
348  hpOccu->SetMaximum(10);
349  c->Clear();
350  c->cd(0);
351  hpOccu->SetMarkerStyle(7);
352  hpOccu->Draw("hist P");
353  c->SetLogy();
354  }
355 
356  // ========== All events =============
357  c = m_c_instOccuAll[g];
358  hLER = findHistT<TH1F>("SVDDoseLERInjAll/SVDInstOccu_" + group.nameSuffix);
359  hHER = findHistT<TH1F>("SVDDoseHERInjAll/SVDInstOccu_" + group.nameSuffix);
360  hNo = findHistT<TH1F>("SVDDoseNoInjAll/SVDInstOccu_" + group.nameSuffix);
361  if (hLER && hHER && hNo) {
362  hLER->SetLineColor(kRed);
363  hHER->SetLineColor(kAzure);
364  hNo->SetLineColor(kBlack);
365  carryOverflowOver(hLER);
366  carryOverflowOver(hHER);
367  carryOverflowOver(hNo);
368  c->Clear();
369  c->cd(0);
370  hNo->SetTitle("SVD instantaneous occu. " + group.titleSuffix + " U-side all events");
371  hNo->Draw("hist"); // hNo usually has the larger maximum by far
372  hLER->Draw("hist same");
373  hHER->Draw("hist same");
374  c->SetLogy();
375  m_legend->Draw();
376  }
377 
378  c = m_c_occuLERAll[g];
379  hHits = findHistT<TH2F>("SVDDoseLERInjAll/SVDHitsVsTime_" + group.nameSuffix);
380  hEvts = findHistT<TH2F>("SVDDoseLERInjAll/SVDEvtsVsTime");
381  if (hHits && hEvts) {
382  if (m_h_occuLERAll[g]) delete m_h_occuLERAll[g];
383  auto hOccu = m_h_occuLERAll[g] = divide(hHits, hEvts, 100.0f / group.nStrips);
384  hOccu->SetTitle("SVD Occupancy " + group.titleSuffix + " - LER inj. all events"
385  ";Time since last injection [#mus];Time in beam cycle [#mus]"
386  ";Occupancy [%]");
387  hOccu->SetMinimum(1e-3);
388  hOccu->SetMaximum(10);
389  c->Clear();
390  c->cd(0);
391  c->SetRightMargin(0.16); // For the colorbar
392  hOccu->Draw("COLZ");
393  c->SetLogz();
394  }
395 
396  c = m_c_occuLER1All[g];
397  hpEvts = findHistT<TH1F>("SVDDoseLERInjAll/SVDEvtsVsTime1");
398  hpHits = findHistT<TH1F>("SVDDoseLERInjAll/SVDHitsVsTime1_" + group.nameSuffix);
399  if (hpHits && hpEvts) {
400  if (m_h_occuLER1All[g]) delete m_h_occuLER1All[g];
401  auto hpOccu = m_h_occuLER1All[g] = divide(hpHits, hpEvts, 100.0f / group.nStrips);
402  hpOccu->SetTitle("SVD Occupancy " + group.titleSuffix + " - LER inj. all events"
403  ";Time since last injection [#mus];Occupancy [%]");
404  hpOccu->SetMinimum(1e-3);
405  hpOccu->SetMaximum(10);
406  c->Clear();
407  c->cd(0);
408  hpOccu->SetMarkerStyle(7);
409  hpOccu->Draw("hist P");
410  c->SetLogy();
411  }
412 
413  c = m_c_occuHERAll[g];
414  hHits = findHistT<TH2F>("SVDDoseHERInjAll/SVDHitsVsTime_" + group.nameSuffix);
415  hEvts = findHistT<TH2F>("SVDDoseHERInjAll/SVDEvtsVsTime");
416  if (hHits && hEvts) {
417  if (m_h_occuHERAll[g]) delete m_h_occuHERAll[g];
418  auto hOccu = m_h_occuHERAll[g] = divide(hHits, hEvts, 100.0f / group.nStrips);
419  hOccu->SetTitle("SVD Occupancy " + group.titleSuffix + " - HER inj. all events"
420  ";Time since last injection [#mus];Time in beam cycle [#mus]"
421  ";Occupancy [%]");
422  hOccu->SetMinimum(1e-3);
423  hOccu->SetMaximum(10);
424  c->Clear();
425  c->cd(0);
426  c->SetRightMargin(0.16); // For the colorbar
427  hOccu->Draw("COLZ");
428  c->SetLogz();
429  }
430 
431  c = m_c_occuHER1All[g];
432  hpEvts = findHistT<TH1F>("SVDDoseHERInjAll/SVDEvtsVsTime1");
433  hpHits = findHistT<TH1F>("SVDDoseHERInjAll/SVDHitsVsTime1_" + group.nameSuffix);
434  if (hpHits && hpEvts) {
435  if (m_h_occuHER1All[g]) delete m_h_occuHER1All[g];
436  auto hpOccu = m_h_occuHER1All[g] = divide(hpHits, hpEvts, 100.0f / group.nStrips);
437  hpOccu->SetTitle("SVD Occupancy " + group.titleSuffix + " - HER inj. all events"
438  ";Time since last injection [#mus];Occupancy [%]");
439  hpOccu->SetMinimum(1e-3);
440  hpOccu->SetMaximum(10);
441  c->Clear();
442  c->cd(0);
443  hpOccu->SetMarkerStyle(7);
444  hpOccu->Draw("hist P");
445  c->SetLogy();
446  }
447  }
448 }
449 
450 void DQMHistAnalysisSVDDoseModule::carryOverflowOver(TH1F* h)
451 {
452  int i = h->GetNbinsX();
453  float t = h->GetBinContent(i) + h->GetBinContent(i + 1);
454  h->SetBinContent(i, t);
455  h->SetBinContent(i + 1, 0);
456 }
457 
458 const vector<DQMHistAnalysisSVDDoseModule::SensorGroup> DQMHistAnalysisSVDDoseModule::c_sensorGroups = {
459  {"L31XU", "L3.1", "L3:1", 768 * 2},
460  {"L32XU", "L3.2", "L3:2", 768 * 2},
461  {"L3XXU", "L3 avg.", "L3", 768 * 14},
462  {"L4XXU", "L4 avg.", "L4", 768 * 30},
463  {"L5XXU", "L5 avg.", "L5", 768 * 48},
464  {"L6XXU", "L6 avg.", "L6", 768 * 80}
465 };
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
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