Belle II Software  release-05-01-25
DQMHistComparitor.cc
1 //+
2 // File : DQMHistComparitor.cc
3 // Description : DQM Histogram analysis module, compares histo with references
4 //
5 // Author : Bjoern Spruck, Uni Mainz
6 // Date : yesterday
7 //-
8 
9 
10 #include <framework/core/ModuleParam.templateDetails.h>
11 #include <dqm/analysis/modules/DQMHistComparitor.h>
12 #include <daq/slc/base/StringUtil.h>
13 #include <TROOT.h>
14 #include <TStyle.h>
15 #include <TClass.h>
16 #include <TDirectory.h>
17 #include <TH1F.h>
18 #include <TH2F.h>
19 #include <TKey.h>
20 
21 using namespace std;
22 using namespace Belle2;
23 
24 //-----------------------------------------------------------------
25 // Register the Module
26 //-----------------------------------------------------------------
27 REG_MODULE(DQMHistComparitor)
28 
29 //-----------------------------------------------------------------
30 // Implementation
31 //-----------------------------------------------------------------
32 
35 {
36  //Parameter definition
37  addParam("HistoList", m_histlist, "['histname', 'refhistname', 'canvas', 'minentry', 'warnlvl', 'errlvl', 'pvname'],[...],...");
38  addParam("RefHistoFile", m_refFileName, "Reference histrogram file name", std::string("refHisto.root"));
39  addParam("ColorAlert", m_color, "Whether to show the color alert", true);
40  B2DEBUG(1, "DQMHistComparitor: Constructor done.");
41 }
42 
43 DQMHistComparitorModule::~DQMHistComparitorModule()
44 {
45 #ifdef _BELLE2_EPICS
46  if (ca_current_context()) ca_context_destroy();
47 #endif
48 }
49 
50 TH1* DQMHistComparitorModule::find_histo_in_canvas(TString histo_name)
51 {
52  StringList s = StringUtil::split(histo_name.Data(), '/');
53  std::string dirname = s[0];
54  std::string hname = s[1];
55  std::string canvas_name = dirname + "/c_" + hname;
56 
57  TIter nextckey(gROOT->GetListOfCanvases());
58  TObject* cobj = NULL;
59 
60  while ((cobj = (TObject*)nextckey())) {
61  if (cobj->IsA()->InheritsFrom("TCanvas")) {
62  if (cobj->GetName() == canvas_name)
63  break;
64  }
65  }
66  if (cobj == NULL) return NULL;
67 
68  TIter nextkey(((TCanvas*)cobj)->GetListOfPrimitives());
69  TObject* obj = NULL;
70 
71  while ((obj = (TObject*)nextkey())) {
72  if (obj->IsA()->InheritsFrom("TH1")) {
73  if (obj->GetName() == histo_name)
74  return (TH1*)obj;
75  }
76  }
77  return NULL;
78 }
79 
80 TH1* DQMHistComparitorModule::GetHisto(TString histoname)
81 {
82  TH1* hh1;
83  gROOT->cd();
84  hh1 = findHist(histoname.Data());
85  if (hh1 == NULL) {
86  B2DEBUG(20, "findHisto failed " << histoname << " not in memfile");
87  // the following code sux ... is there no root function for that?
88 
89 
90  // first search reference root file ... if ther is one
91  if (m_refFile && m_refFile->IsOpen()) {
92  TDirectory* d = m_refFile;
93  TString myl = histoname;
94  TString tok;
95  Ssiz_t from = 0;
96  B2DEBUG(20, myl);
97  while (myl.Tokenize(tok, from, "/")) {
98  TString dummy;
99  Ssiz_t f;
100  f = from;
101  if (myl.Tokenize(dummy, f, "/")) { // check if its the last one
102  auto e = d->GetDirectory(tok);
103  if (e) {
104  B2DEBUG(20, "Cd Dir " << tok << " from " << d->GetPath());
105  d = e;
106  } else {
107  B2DEBUG(20, "cd failed " << tok << " from " << d->GetPath());
108  }
109  } else {
110  break;
111  }
112  }
113  TObject* obj = d->FindObject(tok);
114  if (obj != NULL) {
115  if (obj->IsA()->InheritsFrom("TH1")) {
116  B2DEBUG(20, "Histo " << histoname << " found in ref file");
117  hh1 = (TH1*)obj;
118  } else {
119  B2DEBUG(20, "Histo " << histoname << " found in ref file but wrong type");
120  }
121  } else {
122  // seems find will only find objects, not keys, thus get the object on first access
123  TIter next(d->GetListOfKeys());
124  TKey* key;
125  while ((key = (TKey*)next())) {
126  TObject* obj2 = key->ReadObj() ;
127  if (obj2->InheritsFrom("TH1")) {
128  if (obj2->GetName() == tok) {
129  hh1 = (TH1*)obj2;
130  B2DEBUG(20, "Histo " << histoname << " found as key -> readobj");
131  break;
132  }
133  }
134  }
135  if (hh1 == NULL) B2DEBUG(20, "Histo " << histoname << " NOT found in ref file " << tok);
136  }
137  }
138 
139  if (hh1 == NULL) {
140  B2DEBUG(20, "Histo " << histoname << " not in memfile or ref file");
141  // the following code sux ... is there no root function for that?
142 
143  TDirectory* d = gROOT;
144  TString myl = histoname;
145  TString tok;
146  Ssiz_t from = 0;
147  while (myl.Tokenize(tok, from, "/")) {
148  TString dummy;
149  Ssiz_t f;
150  f = from;
151  if (myl.Tokenize(dummy, f, "/")) { // check if its the last one
152  auto e = d->GetDirectory(tok);
153  if (e) {
154  B2DEBUG(20, "Cd Dir " << tok);
155  d = e;
156  } else B2DEBUG(20, "cd failed " << tok);
157  d->cd();
158  } else {
159  break;
160  }
161  }
162  TObject* obj = d->FindObject(tok);
163  if (obj != NULL) {
164  if (obj->IsA()->InheritsFrom("TH1")) {
165  B2DEBUG(20, "Histo " << histoname << " found in mem");
166  hh1 = (TH1*)obj;
167  }
168  } else {
169  B2DEBUG(20, "Histo " << histoname << " NOT found in mem");
170  }
171  }
172  }
173 
174  if (hh1 == NULL) {
175  B2DEBUG(20, "Histo " << histoname << " not found");
176  }
177 
178  return hh1;
179 }
180 
181 void DQMHistComparitorModule::initialize()
182 {
183  m_refFile = NULL;
184  if (m_refFileName != "") {
185  m_refFile = new TFile(m_refFileName.data());
186  }
187 
188  gStyle->SetOptStat(0);
189  gStyle->SetStatStyle(1);
190  gStyle->SetOptDate(22);// Date and Time in Bottom Right, does no work
191 
192 #ifdef _BELLE2_EPICS
193  if (!ca_current_context()) SEVCHK(ca_context_create(ca_disable_preemptive_callback), "ca_context_create");
194 #endif
195  for (auto& it : m_histlist) {
196  if (it.size() != 7) {
197  B2WARNING("Histolist with wrong nr of parameters (" << it.size() << "), hist1=" << it.at(0));
198  continue;
199  }
200 
201  // Do not check for histogram existance here, as it might not be
202  // in memfile when analysis task is started
203  TString a = it.at(2).c_str();
204 
205  auto n = new CMPNODE;
206  n->histo1 = it.at(0);
207  n->histo2 = it.at(1);
208  TCanvas* c = new TCanvas(a);
209  n->canvas = c;
210 
211  n->min_entries = atoi(it.at(3).c_str());
212  n->warning = atof(it.at(4).c_str());
213  n->error = atof(it.at(5).c_str());
214  n->epicsflag = false;
215 
216 #ifdef _BELLE2_EPICS
217  if (it.at(6) != "") {
218  SEVCHK(ca_create_channel(it.at(6).c_str(), NULL, NULL, 10, &n->mychid), "ca_create_channel failure");
219  n->epicsflag = true;
220  }
221 #endif
222  m_pnode.push_back(n);
223  }
224 #ifdef _BELLE2_EPICS
225  SEVCHK(ca_pend_io(5.0), "ca_pend_io failure");
226 #endif
227 
228  B2DEBUG(20, "DQMHistComparitor: initialized.");
229 }
230 
231 
232 void DQMHistComparitorModule::beginRun()
233 {
234  B2DEBUG(20, "DQMHistComparitor: beginRun called.");
235 }
236 
237 void DQMHistComparitorModule::event()
238 {
239  for (auto& it : m_pnode) {
240 
241  TH1* hist1, *hist2;
242 
243  B2DEBUG(20, "== Search for " << it->histo1 << " with ref " << it->histo2 << "==");
244  hist1 = find_histo_in_canvas(it->histo1);
245  if (!hist1) B2DEBUG(20, "NOT Found " << it->histo1);
246  hist2 = GetHisto(it->histo2);
247  if (!hist2) B2DEBUG(20, "NOT Found " << it->histo2);
248  // Dont compare if any of hists is missing ... TODO We should clean CANVAS, raise some error if we cannot compare
249  if (!hist1) continue;
250  if (!hist2) continue;
251 
252  B2DEBUG(20, "Compare " << it->histo1 << " with ref " << it->histo2);
253  // if compare normalized ... does not work!
254  // hist2->Scale(hist1->GetEntries()/hist2->GetEntries());
255 
256  double data = hist1->KolmogorovTest(hist2, ""); // returns p value (0 bad, 1 good), N - do not compare normalized
257  // data = hist1->Chi2Test(hist2);// return p value (0 bad, 1 good), ignores normalization
258  // data= BinByBinTest(hits1,hist2);// user function (like Peters test)
259  // printf(" %.2f %.2f %.2f\n",(float)data,it->warning,it->error);
260 #ifdef _BELLE2_EPICS
261  if (it->epicsflag) SEVCHK(ca_put(DBR_DOUBLE, it->mychid, (void*)&data), "ca_set failure");
262 #endif
263  it->canvas->cd();
264  hist2->SetLineStyle(3);// 2 or 3
265  hist2->SetLineColor(3);
266 
267  TIter nextkey(it->canvas->GetListOfPrimitives());
268  TObject* obj = NULL;
269  while ((obj = (TObject*)nextkey())) {
270  if (obj->IsA()->InheritsFrom("TH1")) {
271  if (string(obj->GetName()) == string(hist2->GetName())) {
272  delete obj;
273  }
274  }
275  }
276 
277  // if draw normalized
278  TH1* h;
279  if (1) {
280  h = (TH1*)hist2->Clone(); // Anoying ... Maybe an memory leak? TODO
281  if (abs(hist2->Integral()) > 0)
282  h->Scale(hist1->Integral() / hist2->Integral());
283  } else {
284  hist2->Draw("hist");
285  }
286 
287  h->SetFillColor(0);
288  h->SetStats(kFALSE);
289  hist1->SetFillColor(0);
290  if (h->GetMaximum() > hist1->GetMaximum())
291  hist1->SetMaximum(1.1 * h->GetMaximum());
292  hist1->Draw("hist");
293  h->Draw("hist,same");
294 
295  it->canvas->Pad()->SetFrameFillColor(10);
296  if (m_color) {
297  if (hist1->GetEntries() < it->min_entries) {
298  // not enough Entries
299  it->canvas->Pad()->SetFillColor(6);// Magenta
300  } else {
301  if (data < it->error) {
302  it->canvas->Pad()->SetFillColor(2);// Red
303  } else if (data < it->warning) {
304  it->canvas->Pad()->SetFillColor(5);// Yellow
305  } else {
306  it->canvas->Pad()->SetFillColor(0);// White
307  }
308  }
309  } else {
310  it->canvas->Pad()->SetFillColor(0);// White
311  }
312  it->canvas->Modified();
313  it->canvas->Update();
314  }
315 #ifdef _BELLE2_EPICS
316  SEVCHK(ca_pend_io(5.0), "ca_pend_io failure");
317 #endif
318 }
319 
320 void DQMHistComparitorModule::endRun()
321 {
322  B2DEBUG(20, "DQMHistComparitor: endRun called");
323 }
324 
325 
326 void DQMHistComparitorModule::terminate()
327 {
328  B2DEBUG(20, "DQMHistComparitor: terminate called");
329  if (m_refFile) delete m_refFile;
330 }
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::DQMHistComparitorModule::CMPNODE
The struct for reference histogram comparison.
Definition: DQMHistComparitor.h:37
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::DQMHistComparitorModule
Class definition for the reference histogram display.
Definition: DQMHistComparitor.h:32
Belle2::DQMHistAnalysisModule
The base class for the histogram analysis module.
Definition: DQMHistAnalysis.h:27