Belle II Software  release-05-01-25
DQMHistAnalysisPXDTrackCharge.cc
1 //+
2 // File : DQMHistAnalysisPXDTrackCharge.cc
3 // Description : Analysis of PXD Cluster Charge
4 //
5 // Author : Bjoern Spruck, University Mainz
6 // Date : 2019
7 //-
8 
9 
10 #include <dqm/analysis/modules/DQMHistAnalysisPXDTrackCharge.h>
11 #include <TROOT.h>
12 #include <TStyle.h>
13 #include <TClass.h>
14 #include <TLatex.h>
15 #include <vxd/geometry/GeoCache.h>
16 #include <TKey.h>
17 
18 #include <RooDataHist.h>
19 #include <RooAbsPdf.h>
20 #include <RooPlot.h>
21 
22 
23 using namespace std;
24 using namespace Belle2;
25 
26 //-----------------------------------------------------------------
27 // Register the Module
28 //-----------------------------------------------------------------
29 REG_MODULE(DQMHistAnalysisPXDTrackCharge)
30 
31 //-----------------------------------------------------------------
32 // Implementation
33 //-----------------------------------------------------------------
34 
37 {
38  // This module CAN NOT be run in parallel!
39 
40  //Parameter definition
41  addParam("histogramDirectoryName", m_histogramDirectoryName, "Name of Histogram dir", std::string("PXDER"));
42  addParam("RangeLow", m_rangeLow, "Lower border for fit", 20.);
43  addParam("RangeHigh", m_rangeHigh, "High border for fit", 80.);
44 // addParam("PeakBefore", m_peakBefore, "Range for fit before peak (positive)", 5.);
45 // addParam("PeakAfter", m_peakAfter, "Range for after peak", 40.);
46  addParam("PVPrefix", m_pvPrefix, "PV Prefix", std::string("DQM:PXD:TrackCharge:"));
47  addParam("useEpics", m_useEpics, "useEpics", true);
48  addParam("RefHistoFile", m_refFileName, "Reference histrogram file name", std::string("refHisto.root"));
49  addParam("ColorAlert", m_color, "Whether to show the color alert", true);
50 
51  B2DEBUG(99, "DQMHistAnalysisPXDTrackCharge: Constructor done.");
52 }
53 
54 DQMHistAnalysisPXDTrackChargeModule::~DQMHistAnalysisPXDTrackChargeModule()
55 {
56 #ifdef _BELLE2_EPICS
57  if (m_useEpics) {
58  if (ca_current_context()) ca_context_destroy();
59  }
60 #endif
61 }
62 
63 void DQMHistAnalysisPXDTrackChargeModule::initialize()
64 {
65  B2DEBUG(99, "DQMHistAnalysisPXDTrackCharge: initialized.");
66 
67  m_refFile = NULL;
68  if (m_refFileName != "") {
69  m_refFile = new TFile(m_refFileName.data());
70  }
71 
72  m_monObj = getMonitoringObject("pxd");
73  const VXD::GeoCache& geo = VXD::GeoCache::getInstance();
74 
75  m_rfws = new RooWorkspace("w");
76  m_rfws->factory("Landau::landau(x[0,100],ml[20,10,50],sl[5,1,30])");
77  m_rfws->factory("Gaussian::gauss(x,mg[0],sg[2,0.1,10])");
78  m_rfws->factory("FCONV::lxg(x,landau,gauss)");
79 
80  m_x = m_rfws->var("x");
81  m_x->setRange("signal", m_rangeLow, m_rangeHigh);
82 
83  // collect the list of all PXD Modules in the geometry here
84  std::vector<VxdID> sensors = geo.getListOfSensors();
85  for (VxdID& aVxdID : sensors) {
86  VXD::SensorInfoBase info = geo.getSensorInfo(aVxdID);
87  if (info.getType() != VXD::SensorInfoBase::PXD) continue;
88  m_PXDModules.push_back(aVxdID); // reorder, sort would be better
89 
90  std::string name = "PXD_Track_Cluster_Charge_" + (std::string)aVxdID;
91  std::replace(name.begin(), name.end(), '.', '_');
92  m_cChargeMod[aVxdID] = new TCanvas((m_histogramDirectoryName + "/c_Fit_" + name).data());
93  if (aVxdID == VxdID("1.5.1")) {
94  for (int s = 0; s < 6; s++) {
95  for (int d = 0; d < 4; d++) {
96  m_cChargeModASIC[aVxdID][s][d] = new TCanvas((m_histogramDirectoryName + "/c_Fit_" + name + Form("_s%d_d%d", s + 1, d + 1)).data());
97  }
98  }
99  m_hChargeModASIC2d[aVxdID] = new TH2F(("PXD_TCChargeMPV_" + name).data(), ("PXD TCCharge MPV " + name + ";Switcher;DCD;MPV").data(),
100  6, 0.5, 6.5, 4, 0.5, 4.5);
101  m_cChargeModASIC2d[aVxdID] = new TCanvas((m_histogramDirectoryName + "/c_TCCharge_MPV_" + name).data());
102  }
103  }
104  std::sort(m_PXDModules.begin(), m_PXDModules.end()); // back to natural order
105 
106  gROOT->cd(); // this seems to be important, or strange things happen
107 
108  m_cTrackedClusters = new TCanvas((m_histogramDirectoryName + "/c_TrackedClusters").data());
109  m_hTrackedClusters = new TH1F("Tracked_Clusters", "Tracked Clusters/Event;Module", 40, 0, 40);
110  m_hTrackedClusters->Draw();
111  auto ax = m_hTrackedClusters->GetXaxis();
112  if (ax) {
113  ax->Set(m_PXDModules.size(), 0, m_PXDModules.size());
114  for (unsigned int i = 0; i < m_PXDModules.size(); i++) {
115  TString ModuleName = (std::string)m_PXDModules[i];
116  ax->SetBinLabel(i + 1, ModuleName);
117  }
118  } else B2ERROR("no axis");
119 
120  m_cCharge = new TCanvas((m_histogramDirectoryName + "/c_TrackCharge").data());
121  m_monObj->addCanvas(m_cCharge);
122 
123  m_gCharge = new TGraphErrors();
124  m_gCharge->SetName("Track_Cluster_Charge");
125  m_gCharge->SetTitle("Track Cluster Charge");
126 
128  m_line_up = new TLine(0, 10, m_PXDModules.size(), 10);
129  m_line_mean = new TLine(0, 16, m_PXDModules.size(), 16);
130  m_line_low = new TLine(0, 3, m_PXDModules.size(), 3);
131  m_line_up->SetHorizontal(true);
132  m_line_up->SetLineColor(kMagenta);// Green
133  m_line_up->SetLineWidth(3);
134  m_line_up->SetLineStyle(7);
135  m_line_mean->SetHorizontal(true);
136  m_line_mean->SetLineColor(kGreen);// Black
137  m_line_mean->SetLineWidth(3);
138  m_line_mean->SetLineStyle(4);
139  m_line_low->SetHorizontal(true);
140  m_line_low->SetLineColor(kMagenta);
141  m_line_low->SetLineWidth(3);
142  m_line_low->SetLineStyle(7);
143 
144  m_fMean = new TF1("f_Mean", "pol0", 0, m_PXDModules.size());
145  m_fMean->SetParameter(0, 50);
146  m_fMean->SetLineColor(kYellow);
147  m_fMean->SetLineWidth(3);
148  m_fMean->SetLineStyle(7);
149  m_fMean->SetNpx(m_PXDModules.size());
150  m_fMean->SetNumberFitPoints(m_PXDModules.size());
151 
152 #ifdef _BELLE2_EPICS
153  if (m_useEpics) {
154  if (!ca_current_context()) SEVCHK(ca_context_create(ca_disable_preemptive_callback), "ca_context_create");
155  mychid.resize(3);
156  SEVCHK(ca_create_channel((m_pvPrefix + "Mean").data(), NULL, NULL, 10, &mychid[0]), "ca_create_channel failure");
157  SEVCHK(ca_create_channel((m_pvPrefix + "Diff").data(), NULL, NULL, 10, &mychid[1]), "ca_create_channel failure");
158  SEVCHK(ca_create_channel((m_pvPrefix + "Status").data(), NULL, NULL, 10, &mychid[2]), "ca_create_channel failure");
159  SEVCHK(ca_pend_io(5.0), "ca_pend_io failure");
160  }
161 #endif
162 }
163 
164 
165 void DQMHistAnalysisPXDTrackChargeModule::beginRun()
166 {
167  B2DEBUG(99, "DQMHistAnalysisPXDTrackCharge: beginRun called.");
168 
169  m_cCharge->Clear();
170 }
171 
172 void DQMHistAnalysisPXDTrackChargeModule::event()
173 {
174  if (!m_cCharge) return;
175 
176  gStyle->SetOptStat(0);
177  gStyle->SetStatStyle(1);
178  gStyle->SetOptDate(22);// Date and Time in Bottom Right, does no work
179 
180  bool enough = false;
181 
182 // auto landau = m_rfws->pdf("landau");
183 // auto gauss = m_rfws->pdf("gauss");
184  auto model = m_rfws->pdf("lxg");
185 
186  auto ml = m_rfws->var("ml");
187 // auto sl = m_rfws->var("sl");
188 // auto mg = m_rfws->var("mg");
189 // auto sg = m_rfws->var("sg");
190 
191  {
192  m_cTrackedClusters->Clear();
193  m_cTrackedClusters->cd();
194  m_hTrackedClusters->Reset();
195 
196  std::string name = "Tracked_Clusters"; // new name
197  TH1* hh2 = findHist(m_histogramDirectoryName, "PXD_Tracked_Clusters");
198  if (hh2) {
199  auto scale = hh2->GetBinContent(0);// overflow misused as event counter!
200  if (scale > 0) {
201  int j = 1;
202  for (int i = 0; i < 64; i++) {
203  auto layer = (((i >> 5) & 0x1) + 1);
204  auto ladder = ((i >> 1) & 0xF);
205  auto sensor = ((i & 0x1) + 1);
206 
207  auto id = Belle2::VxdID(layer, ladder, sensor);
208  // Check if sensor exist
209  if (Belle2::VXD::GeoCache::getInstance().validSensorID(id)) {
210  m_hTrackedClusters->SetBinContent(j, hh2->GetBinContent(i + 1) / scale);
211  j++;
212  }
213  }
214  }
215  m_hTrackedClusters->SetName(name.data());
216  m_hTrackedClusters->SetTitle("Tracked Clusters/Event");
217  m_hTrackedClusters->SetFillColor(kWhite);
218  m_hTrackedClusters->SetStats(kFALSE);
219  m_hTrackedClusters->SetLineStyle(1);// 2 or 3
220  m_hTrackedClusters->SetLineColor(kBlue);
221  m_hTrackedClusters->Draw("hist");
222 
223  TH1* href2 = GetHisto("ref/" + m_histogramDirectoryName + "/" + name);
224 
225  if (href2) {
226  href2->SetLineStyle(3);// 2 or 3
227  href2->SetLineColor(kBlack);
228  href2->Draw("same,hist");
229  }
230 
231  auto tt = new TLatex(5.5, 0.1, "1.3.2 Module is excluded, please ignore");
232  tt->SetTextAngle(90);// Rotated
233  tt->SetTextAlign(12);// Centered
234  tt->Draw();
235  }
236  }
237 
238  m_gCharge->Set(0);
239 
240  for (unsigned int i = 0; i < m_PXDModules.size(); i++) {
241  TCanvas* canvas = m_cChargeMod[m_PXDModules[i]];
242  if (canvas == nullptr) continue;
243 
244  std::string name = "PXD_Track_Cluster_Charge_" + (std::string)m_PXDModules[i];
245  std::replace(name.begin(), name.end(), '.', '_');
246 
247  canvas->cd();
248  canvas->Clear();
249 
250  TH1* hh1 = findHist(m_histogramDirectoryName, name);
251  if (hh1) {
252 
253  if (hh1->GetEntries() > 50) {
254 
255  auto hdata = new RooDataHist(hh1->GetName(), hh1->GetTitle(), *m_x, (const TH1*) hh1);
256  auto plot = m_x->frame(RooFit::Title(hh1->GetTitle()));
257  /*auto r =*/ model->fitTo(*hdata, RooFit::Range("signal"));
258 
259  model->paramOn(plot, RooFit::Format("NELU", RooFit::AutoPrecision(2)), RooFit::Layout(0.6, 0.9, 0.9));
260  hdata->plotOn(plot, RooFit::LineColor(kBlue)/*, RooFit::Range("plot"), RooFit::NormRange("signal")*/);
261  model->plotOn(plot, RooFit::LineColor(kRed), RooFit::Range("signal"), RooFit::NormRange("signal"));
262 
263  plot->Draw("");
264 
265 // model->Print("");
266 // ml->Print("");
267 // sl->Print("");
268 // mg->Print("");
269 // sg->Print("");
270 // cout << "ZZZ , " << Form("%d%02d%d ,", std::get<0>(t), std::get<1>(t), std::get<2>(t)) << ml->getValV() << "," << ml->getError() << "," << sl->getValV() << "," << sl->getError() << "," << sg->getValV() << "," << sg->getError() << endl;
271 
272 
273  int p = m_gCharge->GetN();
274  m_gCharge->SetPoint(p, i + 0.49, ml->getValV());
275  m_gCharge->SetPointError(p, 0.1, ml->getError()); // error in x is useless
276  m_monObj->setVariable(("trackcharge_" + (std::string)m_PXDModules[i]).c_str(), ml->getValV(), ml->getError());
277  }
278 
279  TH1* hist2 = GetHisto("ref/" + m_histogramDirectoryName + "/" + name);
280 
281  if (hist2) {
282 // B2INFO("Draw Normalized " << hist2->GetName());
283  hist2->SetLineStyle(3);// 2 or 3
284  hist2->SetLineColor(kBlack);
285 
286 // TIter nextkey(canvas->GetListOfPrimitives());
287 // TObject* obj = NULL;
288 // while ((obj = (TObject*)nextkey())) {
289 // if (obj->IsA()->InheritsFrom("TH1")) {
290 // if (string(obj->GetName()) == string(hist2->GetName())) {
291 // delete obj;
292 // }
293 // }
294 // }
295 
296  canvas->cd();
297 
298  // if draw normalized
299  TH1* h = (TH1*)hist2->Clone(); // Annoying ... Maybe an memory leak? TODO
300  if (abs(hist2->GetEntries()) > 0) h->Scale(hh1->GetEntries() / hist2->GetEntries());
301 
302  h->SetStats(kFALSE);
303  h->Draw("same,hist");
304  }
305 
306  // add coloring, cuts? based on fit, compare with ref?
307  canvas->Pad()->SetFrameFillColor(10);
308  if (m_color) {
309  if (hh1->GetEntries() < 100) {
310  // not enough Entries
311  canvas->Pad()->SetFillColor(kGray);
312  } else {
313  canvas->Pad()->SetFillColor(kGreen);
314  }
315  } else {
316  canvas->Pad()->SetFillColor(kWhite);// White
317  }
318 
319  canvas->Modified();
320  canvas->Update();
321 
322  // means if ANY plot is > 100 entries, all plots are assumed to be o.k.
323  if (hh1->GetEntries() >= 1000) enough = true;
324  }
325  }
326 
327  // now loop per module over asics pairs (1.5.1)
328  for (unsigned int i = 0; i < m_PXDModules.size(); i++) {
329 // TCanvas* canvas = m_cChargeMod[m_PXDModules[i]];
330  VxdID& aVxdID = m_PXDModules[i];
331 
332  if (m_hChargeModASIC2d[aVxdID]) m_hChargeModASIC2d[aVxdID]->Reset();
333  if (m_cChargeModASIC2d[aVxdID]) m_cChargeModASIC2d[aVxdID]->Clear();
334 
335  for (int s = 1; s <= 6; s++) {
336  for (int d = 1; d <= 4; d++) {
337  std::string name = "PXD_Track_Cluster_Charge_" + (std::string)m_PXDModules[i] + Form("_sw%d_dcd%d", s, d);
338  std::replace(name.begin(), name.end(), '.', '_');
339 
340  if (m_cChargeModASIC[aVxdID][s - 1][d - 1]) {
341  m_cChargeModASIC[aVxdID][s - 1][d - 1]->Clear();
342  m_cChargeModASIC[aVxdID][s - 1][d - 1]->cd();
343  }
344 
345  TH1* hh1 = findHist(m_histogramDirectoryName, name);
346  if (hh1) {
347  double mpv = 0.0;
348  if (hh1->GetEntries() > 50) {
349  auto hdata = new RooDataHist(hh1->GetName(), hh1->GetTitle(), *m_x, (const TH1*) hh1);
350  auto plot = m_x->frame(RooFit::Title(hh1->GetTitle()));
351  /*auto r =*/ model->fitTo(*hdata, RooFit::Range("signal"));
352 
353  if (m_cChargeModASIC[aVxdID][s - 1][d - 1]) {
354  model->paramOn(plot, RooFit::Format("NELU", RooFit::AutoPrecision(2)), RooFit::Layout(0.6, 0.9, 0.9));
355  hdata->plotOn(plot, RooFit::LineColor(kBlue)/*, RooFit::Range("plot"), RooFit::NormRange("signal")*/);
356  model->plotOn(plot, RooFit::LineColor(kRed), RooFit::Range("signal"), RooFit::NormRange("signal"));
357  }
358  plot->Draw("");
359 
360  mpv = ml->getValV();
361  }
362 
363  if (m_hChargeModASIC2d[aVxdID]) {
364  if (mpv > 0.0) m_hChargeModASIC2d[aVxdID]->Fill(s, d, mpv); // TODO check what is s, d
365  }
366  }
367  }
368  }
369 
370  // Overview map of ASCI combinations
371  if (m_hChargeModASIC2d[aVxdID] && m_cChargeModASIC2d[aVxdID]) {
372  m_cChargeModASIC2d[aVxdID]->cd();
373  m_hChargeModASIC2d[aVxdID]->Draw("colz");
374  }
375  }
376 
377  m_cCharge->cd();
378  m_cCharge->Clear();
379  m_gCharge->SetMinimum(0);
380  m_gCharge->SetMaximum(70);
381  auto ax = m_gCharge->GetXaxis();
382  if (ax) {
383  ax->Set(m_PXDModules.size(), 0, m_PXDModules.size());
384  for (unsigned int i = 0; i < m_PXDModules.size(); i++) {
385  TString ModuleName = (std::string)m_PXDModules[i];
386  ax->SetBinLabel(i + 1, ModuleName);
387  }
388  } else B2ERROR("no axis");
389 
390  m_gCharge->SetLineColor(4);
391  m_gCharge->SetLineWidth(2);
392  m_gCharge->SetMarkerStyle(8);
393  m_gCharge->Draw("AP");
394  {
395  auto tt = new TLatex(5.5, 0.1, "1.3.2 Module is excluded, please ignore");
396  tt->SetTextAngle(90);// Rotated
397  tt->SetTextAlign(12);// Centered
398  tt->Draw();
399  }
400  m_cCharge->cd(0);
401  m_cCharge->Modified();
402  m_cCharge->Update();
403 
404  double data = 0;
405  double diff = 0;
406  if (m_gCharge && enough) {
407 // double currentMin, currentMax;
408  m_gCharge->Fit(m_fMean, "R");
409  double mean = m_gCharge->GetMean(2);
410  double maxi = mean + 15;
411  double mini = mean - 15;
412  m_line_up->SetY1(maxi);
413  m_line_up->SetY2(maxi);
414  m_line_mean->SetY1(mean);
415  m_line_mean->SetY2(mean);
416  m_line_low->SetY1(mini);
417  m_line_low->SetY2(mini);
418  data = mean; // m_fMean->GetParameter(0); // we are more interessted in the maximum deviation from mean
419  // m_gCharge->GetMinimumAndMaximum(currentMin, currentMax);
420  diff = m_gCharge->GetRMS(2);// RMS of Y
421  // better, max deviation as fabs(data - currentMin) > fabs(currentMax - data) ? fabs(data - currentMin) : fabs(currentMax - data);
422  m_line_up->Draw();
423  m_line_mean->Draw();
424  m_line_low->Draw();
425 
426  m_monObj->setVariable("trackcharge", mean, diff);
427  }
428 
429 #ifdef _BELLE2_EPICS
430  if (m_useEpics) {
431  SEVCHK(ca_put(DBR_DOUBLE, mychid[0], (void*)&data), "ca_set failure");
432  SEVCHK(ca_put(DBR_DOUBLE, mychid[1], (void*)&diff), "ca_set failure");
433  }
434 #endif
435 
436  int status = 0;
437 
438  if (!enough) {
439  // not enough Entries
440  m_cCharge->Pad()->SetFillColor(kGray);// Magenta or Gray
441  status = 0; // default
442  } else {
444  if (fabs(data - 30.) > 20. || diff > 12) {
445  m_cCharge->Pad()->SetFillColor(kRed);// Red
446  status = 4;
447  } else if (fabs(data - 30) > 15. || diff > 8) {
448  m_cCharge->Pad()->SetFillColor(kYellow);// Yellow
449  status = 3;
450  } else {
451  m_cCharge->Pad()->SetFillColor(kGreen);// Green
452  status = 2;
453  }
454 
455  // FIXME overwrite for now
456  // m_cCharge->Pad()->SetFillColor(kGreen);// Green
457 
458  }
459 
460 #ifdef _BELLE2_EPICS
461  if (m_useEpics) {
462  SEVCHK(ca_put(DBR_INT, mychid[2], (void*)&status), "ca_set failure");
463  SEVCHK(ca_pend_io(5.0), "ca_pend_io failure");
464  }
465 #endif
466 
467 }
468 
469 void DQMHistAnalysisPXDTrackChargeModule::endRun()
470 {
471  B2DEBUG(99, "DQMHistAnalysisPXDTrackCharge : endRun called");
472 }
473 
474 
475 void DQMHistAnalysisPXDTrackChargeModule::terminate()
476 {
477  B2DEBUG(99, "DQMHistAnalysisPXDTrackCharge: terminate called");
478  // should delete canvas here, maybe hist, too? Who owns it?
479 #ifdef _BELLE2_EPICS
480  if (m_useEpics) {
481  for (auto m : mychid) SEVCHK(ca_clear_channel(m), "ca_clear_channel failure");
482  SEVCHK(ca_pend_io(5.0), "ca_pend_io failure");
483  }
484 #endif
485  if (m_refFile) delete m_refFile;
486 
487 }
488 
489 TH1* DQMHistAnalysisPXDTrackChargeModule::GetHisto(TString histoname)
490 {
491  TH1* hh1 = nullptr;
492  gROOT->cd();
493 // hh1 = findHist(histoname.Data());
494  // cppcheck-suppress knownConditionTrueFalse
495  if (hh1 == NULL) {
496  B2DEBUG(20, "findHisto failed " << histoname << " not in memfile");
497  // the following code sux ... is there no root function for that?
498 
499 
500  // first search reference root file ... if ther is one
501  if (m_refFile && m_refFile->IsOpen()) {
502  TDirectory* d = m_refFile;
503  TString myl = histoname;
504  TString tok;
505  Ssiz_t from = 0;
506  B2DEBUG(20, myl);
507  while (myl.Tokenize(tok, from, "/")) {
508  TString dummy;
509  Ssiz_t f;
510  f = from;
511  if (myl.Tokenize(dummy, f, "/")) { // check if its the last one
512  auto e = d->GetDirectory(tok);
513  if (e) {
514  B2DEBUG(20, "Cd Dir " << tok << " from " << d->GetPath());
515  d = e;
516  } else {
517  B2DEBUG(20, "cd failed " << tok << " from " << d->GetPath());
518  }
519  } else {
520  break;
521  }
522  }
523  TObject* obj = d->FindObject(tok);
524  if (obj != NULL) {
525  if (obj->IsA()->InheritsFrom("TH1")) {
526  B2DEBUG(20, "Histo " << histoname << " found in ref file");
527  hh1 = (TH1*)obj;
528  } else {
529  B2DEBUG(20, "Histo " << histoname << " found in ref file but wrong type");
530  }
531  } else {
532  // seems find will only find objects, not keys, thus get the object on first access
533  TIter next(d->GetListOfKeys());
534  TKey* key;
535  while ((key = (TKey*)next())) {
536  TObject* obj2 = key->ReadObj() ;
537  if (obj2->InheritsFrom("TH1")) {
538  if (obj2->GetName() == tok) {
539  hh1 = (TH1*)obj2;
540  B2DEBUG(20, "Histo " << histoname << " found as key -> readobj");
541  break;
542  }
543  }
544  }
545  if (hh1 == NULL) B2DEBUG(20, "Histo " << histoname << " NOT found in ref file " << tok);
546  }
547  }
548 
549  if (hh1 == NULL) {
550  B2DEBUG(20, "Histo " << histoname << " not in memfile or ref file");
551  // the following code sux ... is there no root function for that?
552 
553  TDirectory* d = gROOT;
554  TString myl = histoname;
555  TString tok;
556  Ssiz_t from = 0;
557  while (myl.Tokenize(tok, from, "/")) {
558  TString dummy;
559  Ssiz_t f;
560  f = from;
561  if (myl.Tokenize(dummy, f, "/")) { // check if its the last one
562  auto e = d->GetDirectory(tok);
563  if (e) {
564  B2DEBUG(20, "Cd Dir " << tok);
565  d = e;
566  } else B2DEBUG(20, "cd failed " << tok);
567  d->cd();
568  } else {
569  break;
570  }
571  }
572  TObject* obj = d->FindObject(tok);
573  if (obj != NULL) {
574  if (obj->IsA()->InheritsFrom("TH1")) {
575  B2DEBUG(20, "Histo " << histoname << " found in mem");
576  hh1 = (TH1*)obj;
577  }
578  } else {
579  B2DEBUG(20, "Histo " << histoname << " NOT found in mem");
580  }
581  }
582  }
583 
584  if (hh1 == NULL) {
585  B2DEBUG(20, "Histo " << histoname << " not found");
586  }
587 
588  return hh1;
589 }
Belle2::VxdID
Class to uniquely identify a any structure of the PXD and SVD.
Definition: VxdID.h:43
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::GeoCache::getListOfSensors
const std::vector< VxdID > getListOfSensors() const
Get list of all sensors.
Definition: GeoCache.cc:60
Belle2::VXD::GeoCache::getInstance
static GeoCache & getInstance()
Return a reference to the singleton instance.
Definition: GeoCache.cc:215
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::DQMHistAnalysisPXDTrackChargeModule
DQM Histogram Analysis for PXD Cluster Charge.
Definition: DQMHistAnalysisPXDTrackCharge.h:39
Belle2::VXD::GeoCache
Class to faciliate easy access to sensor information of the VXD like coordinate transformations or pi...
Definition: GeoCache.h:41
Belle2::VXD::GeoCache::getSensorInfo
const SensorInfoBase & getSensorInfo(Belle2::VxdID id) const
Return a referecne to the SensorInfo of a given SensorID.
Definition: GeoCache.cc:68
Belle2::DQMHistAnalysisModule
The base class for the histogram analysis module.
Definition: DQMHistAnalysis.h:27