Belle II Software  release-05-01-25
DQMHistAnalysisCDCMonObj.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2010 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Makoto Uchida *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 // Own include
12 #include <dqm/analysis/modules/DQMHistAnalysisCDCMonObj.h>
13 
14 //DQM
15 #include <dqm/analysis/modules/DQMHistAnalysis.h>
16 
17 // CDC geometry
18 #include <cdc/geometry/CDCGeometryPar.h>
19 
20 #include <TH1.h>
21 #include <TH1F.h>
22 #include <TH2F.h>
23 #include <TH2Poly.h>
24 #include <TEllipse.h>
25 #include <TF1.h>
26 #include <TCanvas.h>
27 #include <TLine.h>
28 #include <TClass.h>
29 #include <TROOT.h>
30 #include <TString.h>
31 #include <TFile.h>
32 #include <TStyle.h>
33 
34 #include <fstream>
35 #include <vector>
36 #include <algorithm>
37 #include <numeric>
38 
39 using namespace std;
40 using namespace Belle2;
41 
42 //-----------------------------------------------------------------
43 // Register module
44 //-----------------------------------------------------------------
45 
46 REG_MODULE(DQMHistAnalysisCDCMonObj);
47 
48 DQMHistAnalysisCDCMonObjModule::DQMHistAnalysisCDCMonObjModule()
50 {
51  // set module description (e.g. insert text)
52  setDescription("Modify and analyze the data quality histograms of CDCMonObj");
54  for (int i = 0; i < 300; i++) {
55  m_hADCs[i] = nullptr;
56  m_hADCTOTCuts[i] = nullptr;
57  m_hTDCs[i] = nullptr;
58  }
59  for (int i = 0; i < 56; i++) m_hHits[i] = nullptr;
60 }
61 
63 {
64 
65 }
66 
68 {
69 
71  if (!(*m_channelMapFromDB).isValid()) {
72  B2FATAL("Channel map is not valid");
73  }
74 
76  if (m_cdcGeo == nullptr) {
77  B2FATAL("CDCGeometryp is not valid");
78  }
79 
80 
81 
83 
84  gStyle->SetOptStat(0);
85  gStyle->SetPalette(kViridis);
86  gStyle->SetPadTopMargin(0.1);
87  gStyle->SetPadRightMargin(0.05);
88  gStyle->SetPadBottomMargin(0.1);
89  gStyle->SetPadLeftMargin(0.15);
90 
91  m_cMain = new TCanvas("cdc_main", "cdc_main", 1500, 1000);
93 
94  m_cADC = new TCanvas("cdc_adc", "cdc_adc", 2000, 10000);
96  m_cTDC = new TCanvas("cdc_tdc", "cdc_tdc", 2000, 10000);
98  m_cHit = new TCanvas("cdc_hit", "cdc_hit", 1500, 6000);
100 
101  B2DEBUG(20, "DQMHistAnalysisCDCMonObj: initialized.");
102 
103 }
104 
106 {
107  for (const auto& cm : (*m_channelMapFromDB)) {
108  const int isl = cm.getISuperLayer();
109  const int il = cm.getILayer();
110  const int iw = cm.getIWire();
111  const int iBoard = cm.getBoardID();
112  const int iCh = cm.getBoardChannel();
113  const WireID wireId(isl, il, iw);
114  m_chMap.insert(std::make_pair(wireId, std::make_pair(iBoard, iCh)));
115  }
116 
117 
118  const CDCGeometry& geom = **m_cdcGeo;
119 
120  for (const auto& sense : geom.getSenseLayers()) {
121  int i = sense.getId();
122  if (i < 0 || i > 55) {
123  B2FATAL("no such sense layer");
124  }
125  m_senseR[i] = sense.getR();
126  m_nSenseWires[i] = sense.getNWires();
127  m_offset[i] = sense.getOffset();
128  }
129 
130  for (const auto& field : geom.getFieldLayers()) {
131  int i = field.getId();
132  if (i < 0 || i > 54) {
133  B2FATAL("no such sense layer");
134  }
135  m_fieldR[i + 1] = field.getR();
136  }
137  m_fieldR[56] = geom.getOuterWall(0).getRmin();
138  m_fieldR[0] = geom.getInnerWall(0).getRmax();
139 
140 }
141 
143 {
144 }
145 
147 {
148  for (int ilayer = 0; ilayer < 56; ++ilayer) {
149  float dPhi = float(2 * M_PI / m_nSenseWires[ilayer]);
150  float r1 = m_fieldR[ilayer];
151  float r2 = m_fieldR[ilayer + 1];
152  for (int iwire = 0; iwire < m_nSenseWires[ilayer]; ++iwire) {
153  float phi = dPhi * (iwire + m_offset[ilayer]);
154  float phi1 = phi - dPhi * 0.5;
155  float phi2 = phi + dPhi * 0.5;
156  Double_t x_pos[] = {r1* (sin(phi)*tan(phi - phi1) + cos(phi)),
157  r2 * cos(phi1),
158  r2 * cos(phi2),
159  r1* (sin(phi)*tan(phi - phi2) + cos(phi))
160  };
161  Double_t y_pos[] = {r1* (-cos(phi)*tan(phi - phi1) + sin(phi)),
162  r2 * sin(phi1),
163  r2 * sin(phi2),
164  r1* (-cos(phi)*tan(phi - phi2) + sin(phi))
165  };
166  h->AddBin(4, x_pos, y_pos);
167  }
168  }
169 }
170 
171 
173 {
174  m_badChannels.clear();
175  for (int il = 0; il < 56; ++il) {
176  for (int iw = 0; iw < m_nSenseWires[il]; ++iw) {
177  const int y = m_hHits[il]->GetBinContent(iw + 1);
178  if (y == 0) {
179  m_badChannels.push_back(std::make_pair(il, iw));
180  }
181  }
182  }
183  B2DEBUG(20, "num bad wires " << m_badChannels.size());
184 }
185 
187 {
188  TH1D* hist = (TH1D*)h->Clone();
189  hist->SetBinContent(1, 0.0);
190  float m = hist->GetMean();
191  return m;
192 }
193 
194 
195 std::pair<int, int> DQMHistAnalysisCDCMonObjModule::getBoardChannel(unsigned short layer, unsigned short wire)
196 {
197  const WireID w(layer, wire);
198  decltype(m_chMap)::iterator it = m_chMap.find(w);
199  if (it != m_chMap.end()) {
200  return it->second;
201  } else {
202  B2ERROR("no corresponding board/channel found layer " << layer << " wire " << wire);
203  return std::make_pair(-1, -1);
204  }
205 }
206 
207 
209 {
210  B2DEBUG(20, "end run");
211  m_hfastTDC = (TH1F*)findHist("CDC/fast_tdc");
212  m_hADC = (TH2F*)findHist("CDC/hADC");
213  m_hADCTOTCut = (TH2F*)findHist("CDC/hADCTOTCut");
214  m_hTDC = (TH2F*)findHist("CDC/hTDC");
215  m_hHit = (TH2F*)findHist("CDC/hHit");
216  TF1* fitFunc[300] = {};
217  for (int i = 0; i < 300; ++i) {
218  fitFunc[i] = new TF1(Form("f%d", i), "[0]+[6]*x+[1]*(exp([2]*(x-[3]))/(1+exp(-([4]-x)/[5])))",
219  4921 - 100, 4921 + 100);
220  fitFunc[i]->SetParLimits(0, 0, 100);
221  fitFunc[i]->SetParLimits(6, 0, 0.1);
222  fitFunc[i]->SetParLimits(4, 4850., 5000.0);
223  fitFunc[i]->SetParLimits(5, 0, 50.0);
224  }
225 
226  int neve = m_hfastTDC->GetEntries();
227  if (neve == 0)neve = 1;
228 
229  B2DEBUG(20, "adc related");
230  int nDeadADC = -1; // bid 0 always empty
231  int nBadADC = 0;
232  TH1F* hADCMean = new TH1F("hADCMean", "ADC mean;board;adc mean", 300, 0, 300);
233  TH1F* hADC1000 = new TH1F("ADC1000", "ADC1000", 300, 0, 300);
234  TH1F* hADC0 = new TH1F("ADC0", "ADC0", 300, 0, 300);
235 
236  std::vector<float> means = {};
237  for (int i = 0; i < 300; ++i) {
238  m_hADCTOTCuts[i] = m_hADCTOTCut->ProjectionY(Form("hADCTOTCut%d", i), i + 1, i + 1, "");
239  m_hADCTOTCuts[i]->SetTitle(Form("hADCTOTCut%d", i));
240  m_hADCs[i] = m_hADC->ProjectionY(Form("hADC%d", i), i + 1, i + 1, "");
241  m_hADCs[i]->SetTitle(Form("hADC%d", i));
242  float n = static_cast<float>(m_hADCs[i]->GetEntries());
243  if (m_hADCs[i]->GetEntries() == 0) {
244  nDeadADC += 1;
245  hADC0->SetBinContent(i + 1, -0.1);
246  } else {
247  float n0 = static_cast<float>(m_hADCs[i]->GetBinContent(1));
248  if (n0 / n > 0.9) {
249  B2DEBUG(21, "bad adc bid " << i << " " << n0 << " " << n);
250  nBadADC += 1;
251  }
252  float bin1 = m_hADCs[i]->GetBinContent(1);
253  float m = getHistMean(m_hADCs[i]);
254  means.push_back(m);
255  hADCMean->SetBinContent(i + 1, m);
256  hADCMean->SetBinError(i + 1, 0);
257  double overflow = m_hADCs[i]->GetBinContent(1001);
258  hADC1000->SetBinContent(i + 1, overflow / (overflow + n));
259  hADC0->SetBinContent(i + 1, bin1 / (overflow + n));
261  }
262  }
263  // TDC related
264  B2DEBUG(20, "tdc related");
265  int nDeadTDC = 1; // bid 0 always empty
266  TH1F* hTDCEdge = new TH1F("hTDCEdge", "TDC edge;board;tdc edge [nsec]", 300, 0, 300);
267  TH1F* hTDCSlope = new TH1F("hTDCSlope", "TDC slope;board;tdc slope [nsec]", 300, 0, 300);
268  std::vector<float> tdcEdges = {};
269  std::vector<float> tdcSlopes = {};
270  for (int i = 0; i < 300; ++i) {
271  m_hTDCs[i] = m_hTDC->ProjectionY(Form("hTDC%d", i), i + 1, i + 1);
272  m_hTDCs[i]->SetTitle(Form("hTDC%d", i));
273  if (m_hTDCs[i]->GetEntries() == 0 || m_hTDCs[i] == nullptr) {
274  nDeadTDC += 1;
275  } else {
276  fitFunc[i]->SetParameters(0, 100, 0.01, 4700, 4900, 2, 0.01);
277  fitFunc[i]->SetParameter(0, 10);
278  fitFunc[i]->SetParameter(6, 0.02);
279  int xxx = -1;
280  if (i < 28) {
281  xxx = m_hTDCs[i]->Fit(fitFunc[i], "qM0", "", 4850, 5000);
282  } else {
283  xxx = m_hTDCs[i]->Fit(fitFunc[i], "qM0", "", 4800, 5000);
284  }
285  float p4 = fitFunc[i]->GetParameter(4);
286  float p5 = fitFunc[i]->GetParameter(5);
287  if (xxx != -1 && 4850 < p4 && p4 < 5000) {
288  hTDCEdge->SetBinContent(i + 1, p4);
289  hTDCEdge->SetBinError(i + 1, 0);
290  hTDCSlope->SetBinContent(i + 1, p5);
291  hTDCSlope->SetBinError(i + 1, 0);
292  tdcEdges.push_back(p4);
293  tdcSlopes.push_back(p5);
294  }
295  }
296 
297  }
298 
299 
300  // Hit related
301  B2DEBUG(20, "hit related");
302  TH1F* hHitPerLayer = new TH1F("hHitPerLayer", "hit/Layer;layer", 56, 0, 56);
303  int nHits = 0;
304  for (int i = 0; i < 56; ++i) {
305  m_hHits[i] = m_hHit->ProjectionY(Form("hHit%d", i), i + 1, i + 1);
306  m_hHits[i]->SetTitle(Form("hHit%d", i));
307  if (m_hHits[i]->GetEntries() > 0 && m_hHits[i] != nullptr) {
308  int nhitSumL = 0;
309  int nBins = m_nSenseWires[i];
310  for (int j = 0; j < nBins; ++j) {
311  nhitSumL += m_hHits[i]->GetBinContent(j + 1);
312  }
313  if (neve > 0) {
314  hHitPerLayer->SetBinContent(i + 1, static_cast<float>(nhitSumL / neve));
315  } else hHitPerLayer->SetBinContent(i + 1, static_cast<float>(nhitSumL));
316  hHitPerLayer->SetBinError(i + 1, 0);
317  nHits += nhitSumL;
318  }
319  }
320 
321  // Bad wires related
322  B2DEBUG(20, "bad wire related");
323  TH2F* hBadChannel = new TH2F("hbadch", "bad channel map;wire;layer", 400, 0, 400, 56, 0, 56);
324  for (int i = 0; i < 400; ++i) {
325  for (int j = 0; j < 56; ++j) {
326  hBadChannel->Fill(i, j, -1);
327  }
328  }
329 
330  TH2F* hBadChannelBC = new TH2F("hbadchBC", "bad channel map per board/channel;board;channel", 300, 0, 300, 48, 0, 48);
331  for (int i = 0; i < 300; ++i) {
332  for (int j = 0; j < 48; ++j) {
333  hBadChannelBC->Fill(i, j, -1);
334  }
335  }
336 
337  TH2Poly* h2p = new TH2Poly();
338  configureBins(h2p);
339  h2p->SetTitle("bad wires in xy view");
340  h2p->GetXaxis()->SetTitle("X [cm]");
341  h2p->GetYaxis()->SetTitle("Y [cm]");
343  for (const auto& lw : m_badChannels) {
344  const int l = lw.first;
345  const int w = lw.second;
346  B2DEBUG(21, "l " << l << " w " << w);
347  hBadChannel->Fill(w, l);
348  std::pair<int, int> bc = getBoardChannel(l, w);
349  hBadChannelBC->Fill(bc.first, bc.second);
350  float r = m_senseR[l];
351  float dPhi = static_cast<float>(2.0 * M_PI / m_nSenseWires[l]);
352  float phi = dPhi * (w + m_offset[l]);
353  float x = r * cos(phi);
354  float y = r * sin(phi);
355  h2p->Fill(x, y, 1.1);
356  }
357 
358  B2DEBUG(20, "writing");
359  m_cMain->Divide(3, 3);
360 
361  m_cMain->cd(1);
362  hADCMean->SetMinimum(0);
363  hADCMean->SetMaximum(300);
364  hADCMean->Draw();
365 
366  m_cMain->cd(2);
367  hTDCEdge->SetMinimum(4800);
368  hTDCEdge->SetMaximum(5000);
369  hTDCEdge->Draw();
370 
371  m_cMain->cd(3);
372  hTDCSlope->SetMinimum(0);
373  hTDCSlope->SetMaximum(50);
374  hTDCSlope->Draw();
375  m_cMain->cd(4);
376  hBadChannel->Draw("col");
377 
378  m_cMain->cd(5);
379  hBadChannelBC->Draw("col");
380  m_cMain->cd(9);
381  hHitPerLayer->Draw();
382  m_cMain->cd(7);
383  hADC1000->Draw();
384  m_cMain->cd(8);
385  hADC0->Draw();
386 
387  m_cHit->Divide(4, 14);
388  for (int i = 0; i < 56; i++) {
389  m_cHit->cd(i + 1);
390  m_hHits[i]->GetXaxis()->SetRangeUser(0, m_nSenseWires[i]);
391  m_hHits[i]->Draw("hist");
392  }
393 
394  m_cADC->Divide(6, 50, 0.0002, 0.0002);
395  m_cTDC->Divide(6, 50, 0.0002, 0.0002);
396 
397  for (int i = 0; i < 300; i++) {
398  m_cADC->cd(i + 1);
399  m_hADCTOTCuts[i]->SetFillColor(42);
400  Double_t max = m_hADCTOTCuts[i]->GetMaximum();
401  m_hADCs[i]->GetYaxis()->SetRangeUser(0, 3 * max);
402  m_hADCs[i]->Draw("hist");
403  m_hADCTOTCuts[i]->Draw("hist same");
404 
405  m_cTDC->cd(i + 1);
406  m_hTDCs[i]->Draw("hist");
407  fitFunc[i]->SetLineColor(kRed);
408  fitFunc[i]->Draw("same");
409  max = m_hTDCs[i]->GetMaximum();
410  TLine* l1 = new TLine(tdcEdges[i], 0, tdcEdges[i], max * 1.05);
411  l1->SetLineColor(kRed);
412  l1->Draw();
413  TLine* l0 = new TLine(4900, 0, 4900, max * 1.05);
414  l0->Draw();
415  }
416 
417  m_cMain->cd(6);
418  h2p->Draw("col");
419  float superLayerR[10] = {16.3, 24.3, 35.66, 46.63, 57.55, 68.47,
420  79.39, 90.31, 101.23, 112.05
421  };
422 
423  TEllipse* circs[10];
424  for (int i = 0; i < 10; ++i) {
425  circs[i] = new TEllipse(0, 0, superLayerR[i], superLayerR[i]);
426  circs[i]->SetFillStyle(4000);
427  circs[i]->SetLineStyle(kDashed);
428  circs[i]->SetLineColor(0);
429  circs[i]->Draw("same");
430  }
431 
432  m_monObj->setVariable("nEvents", neve);
433  m_monObj->setVariable("nHits", nHits / neve);
434  m_monObj->setVariable("nBadWires", m_badChannels.size());
435  m_monObj->setVariable("adcMean", std::accumulate(means.begin(), means.end(), 0.0) / means.size());
436  m_monObj->setVariable("nDeadADC", nDeadADC);
437  m_monObj->setVariable("nBadADC", nBadADC); //???? n_0/n_tot>0.9
438  m_monObj->setVariable("tdcEdge", std::accumulate(tdcEdges.begin(), tdcEdges.end(), 0.0) / tdcEdges.size());
439  m_monObj->setVariable("nDeadTDC", nDeadTDC);
440  m_monObj->setVariable("tdcSlope", std::accumulate(tdcSlopes.begin(), tdcSlopes.end(), 0.0) / tdcSlopes.size());
441 
442  delete hADCMean;
443  delete hADC1000;
444  delete hADC0;
445  delete hTDCEdge;
446  delete hTDCSlope;
447  delete hHitPerLayer;
448  delete hBadChannel;
449  delete hBadChannelBC;
450  delete h2p;
451 
452 }
453 
455 {
456 
457  B2DEBUG(20, "terminate called");
458 }
459 
Belle2::DQMHistAnalysisCDCMonObjModule::m_senseR
float m_senseR[56]
Radius of sense (+field) layer.
Definition: DQMHistAnalysisCDCMonObj.h:142
Belle2::DQMHistAnalysisCDCMonObjModule::m_cHit
TCanvas * m_cHit
main panel
Definition: DQMHistAnalysisCDCMonObj.h:125
Belle2::WireID
Class to identify a wire inside the CDC.
Definition: WireID.h:44
Belle2::DQMHistAnalysisCDCMonObjModule::m_hHits
TH1D * m_hHits[56]
hit histograms for each layer (0-55)
Definition: DQMHistAnalysisCDCMonObj.h:137
Belle2::DQMHistAnalysisCDCMonObjModule::getBoardChannel
std::pair< int, int > getBoardChannel(unsigned short layer, unsigned short wire)
Get board/channel from layer/wire.
Definition: DQMHistAnalysisCDCMonObj.cc:195
Belle2::Module::setDescription
void setDescription(const std::string &description)
Sets the description of the module.
Definition: Module.cc:216
Belle2::DQMHistAnalysisCDCMonObjModule::m_hADC
TH2F * m_hADC
Summary of ADC histograms
Definition: DQMHistAnalysisCDCMonObj.h:130
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::Module::c_ParallelProcessingCertified
@ 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:82
Belle2::DBArray
Class for accessing arrays of objects in the database.
Definition: DBArray.h:36
Belle2::DQMHistAnalysisCDCMonObjModule::m_hADCTOTCuts
TH1D * m_hADCTOTCuts[300]
ADC histograms with tot cut for each board (0-299)
Definition: DQMHistAnalysisCDCMonObj.h:135
Belle2::DQMHistAnalysisCDCMonObjModule::m_monObj
MonitoringObject * m_monObj
monitoring object
Definition: DQMHistAnalysisCDCMonObj.h:127
Belle2::DQMHistAnalysisCDCMonObjModule::m_fieldR
float m_fieldR[57]
Radius of field layer.
Definition: DQMHistAnalysisCDCMonObj.h:143
Belle2::DQMHistAnalysisCDCMonObjModule::m_hADCs
TH1D * m_hADCs[300]
ADC histograms for each board (0-299)
Definition: DQMHistAnalysisCDCMonObj.h:134
Belle2::DQMHistAnalysisModule::findHist
static TH1 * findHist(const std::string &histname)
Find histogram.
Definition: DQMHistAnalysis.cc:83
Belle2::DQMHistAnalysisCDCMonObjModule::m_hHit
TH2F * m_hHit
Summary of hit histograms.
Definition: DQMHistAnalysisCDCMonObj.h:133
Belle2::DQMHistAnalysisCDCMonObjModule::makeBadChannelList
void makeBadChannelList()
make bad channel list.
Definition: DQMHistAnalysisCDCMonObj.cc:172
Belle2::DQMHistAnalysisCDCMonObjModule::m_cADC
TCanvas * m_cADC
main panel
Definition: DQMHistAnalysisCDCMonObj.h:123
Belle2::DBObjPtr
Class for accessing objects in the database.
Definition: DBObjPtr.h:31
Belle2::DQMHistAnalysisCDCMonObjModule::m_nSenseWires
int m_nSenseWires[56]
number of wires for each layer.
Definition: DQMHistAnalysisCDCMonObj.h:145
Belle2::Module::setPropertyFlags
void setPropertyFlags(unsigned int propertyFlags)
Sets the flags for the module properties.
Definition: Module.cc:210
Belle2::CDCGeometry
The Class for CDC geometry.
Definition: CDCGeometry.h:37
Belle2::DQMHistAnalysisCDCMonObjModule::event
virtual void event() override
Event processor.
Definition: DQMHistAnalysisCDCMonObj.cc:142
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::DQMHistAnalysisCDCMonObjModule::m_cMain
TCanvas * m_cMain
main panel
Definition: DQMHistAnalysisCDCMonObj.h:122
Belle2::MonitoringObject::addCanvas
void addCanvas(TCanvas *canv)
Add Canvas to monitoring object.
Definition: MonitoringObject.h:58
Belle2::MonitoringObject::setVariable
void setVariable(const std::string &var, float val, float upErr=-1., float dwErr=-1)
set value to float variable (new variable is made if not yet existing)
Definition: MonitoringObject.h:134
Belle2::DQMHistAnalysisCDCMonObjModule::m_cdcGeo
DBObjPtr< CDCGeometry > * m_cdcGeo
Geometry of CDC.
Definition: DQMHistAnalysisCDCMonObj.h:146
Belle2::DQMHistAnalysisCDCMonObjModule::endRun
virtual void endRun() override
End-of-run action.
Definition: DQMHistAnalysisCDCMonObj.cc:208
Belle2::DQMHistAnalysisCDCMonObjModule::~DQMHistAnalysisCDCMonObjModule
virtual ~DQMHistAnalysisCDCMonObjModule()
Destructor.
Definition: DQMHistAnalysisCDCMonObj.cc:62
Belle2::DQMHistAnalysisCDCMonObjModule::m_hADCTOTCut
TH2F * m_hADCTOTCut
Summary of ADC histograms with tot cut.
Definition: DQMHistAnalysisCDCMonObj.h:131
Belle2::DQMHistAnalysisCDCMonObjModule::m_hfastTDC
TH1F * m_hfastTDC
Histogram of num.
Definition: DQMHistAnalysisCDCMonObj.h:129
Belle2::DQMHistAnalysisCDCMonObjModule::m_badChannels
std::vector< std::pair< int, int > > m_badChannels
bad wires list
Definition: DQMHistAnalysisCDCMonObj.h:139
Belle2::DQMHistAnalysisCDCMonObjModule::beginRun
virtual void beginRun() override
Called when entering a new run.
Definition: DQMHistAnalysisCDCMonObj.cc:105
Belle2::DQMHistAnalysisCDCMonObjModule::initialize
virtual void initialize() override
Initialize the Module.
Definition: DQMHistAnalysisCDCMonObj.cc:67
Belle2::DQMHistAnalysisCDCMonObjModule::m_offset
float m_offset[56]
Offset of sense layer
Definition: DQMHistAnalysisCDCMonObj.h:144
Belle2::DQMHistAnalysisCDCMonObjModule::m_hTDCs
TH1D * m_hTDCs[300]
TDC histograms for each board (0-299)
Definition: DQMHistAnalysisCDCMonObj.h:136
Belle2::DQMHistAnalysisCDCMonObjModule::m_cTDC
TCanvas * m_cTDC
bad wire panel
Definition: DQMHistAnalysisCDCMonObj.h:124
Belle2::DQMHistAnalysisCDCMonObjModule::m_channelMapFromDB
DBArray< CDCChannelMap > * m_channelMapFromDB
Channel map retrieved from DB.
Definition: DQMHistAnalysisCDCMonObj.h:141
Belle2::DQMHistAnalysisCDCMonObjModule::getHistMean
float getHistMean(TH1D *h)
Get mean of ADC histgram excluding 0-th bin.
Definition: DQMHistAnalysisCDCMonObj.cc:186
Belle2::DQMHistAnalysisCDCMonObjModule::m_chMap
std::map< WireID, std::pair< int, int > > m_chMap
Channel map retrieved
Definition: DQMHistAnalysisCDCMonObj.h:140
Belle2::DQMHistAnalysisCDCMonObjModule::configureBins
void configureBins(TH2Poly *h)
Configure bins of TH2Poly.
Definition: DQMHistAnalysisCDCMonObj.cc:146
Belle2::DQMHistAnalysisCDCMonObjModule::m_hTDC
TH2F * m_hTDC
Summary of TDC histograms.
Definition: DQMHistAnalysisCDCMonObj.h:132
Belle2::DQMHistAnalysisModule::getMonitoringObject
static MonitoringObject * getMonitoringObject(const std::string &histname)
Get MonitoringObject with given name (new object is created if non-existing)
Definition: DQMHistAnalysis.cc:55
Belle2::DQMHistAnalysisCDCMonObjModule::terminate
virtual void terminate() override
Termination action.
Definition: DQMHistAnalysisCDCMonObj.cc:454
Belle2::DQMHistAnalysisModule
The base class for the histogram analysis module.
Definition: DQMHistAnalysis.h:27