10 #include <framework/core/ModuleParam.templateDetails.h>
11 #include <dqm/analysis/modules/DQMHistComparitor.h>
12 #include <daq/slc/base/StringUtil.h>
16 #include <TDirectory.h>
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.");
43 DQMHistComparitorModule::~DQMHistComparitorModule()
46 if (ca_current_context()) ca_context_destroy();
50 TH1* DQMHistComparitorModule::find_histo_in_canvas(TString histo_name)
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;
57 TIter nextckey(gROOT->GetListOfCanvases());
60 while ((cobj = (TObject*)nextckey())) {
61 if (cobj->IsA()->InheritsFrom(
"TCanvas")) {
62 if (cobj->GetName() == canvas_name)
66 if (cobj == NULL)
return NULL;
68 TIter nextkey(((TCanvas*)cobj)->GetListOfPrimitives());
71 while ((obj = (TObject*)nextkey())) {
72 if (obj->IsA()->InheritsFrom(
"TH1")) {
73 if (obj->GetName() == histo_name)
80 TH1* DQMHistComparitorModule::GetHisto(TString histoname)
84 hh1 = findHist(histoname.Data());
86 B2DEBUG(20,
"findHisto failed " << histoname <<
" not in memfile");
91 if (m_refFile && m_refFile->IsOpen()) {
92 TDirectory* d = m_refFile;
93 TString myl = histoname;
97 while (myl.Tokenize(tok, from,
"/")) {
101 if (myl.Tokenize(dummy, f,
"/")) {
102 auto e = d->GetDirectory(tok);
104 B2DEBUG(20,
"Cd Dir " << tok <<
" from " << d->GetPath());
107 B2DEBUG(20,
"cd failed " << tok <<
" from " << d->GetPath());
113 TObject* obj = d->FindObject(tok);
115 if (obj->IsA()->InheritsFrom(
"TH1")) {
116 B2DEBUG(20,
"Histo " << histoname <<
" found in ref file");
119 B2DEBUG(20,
"Histo " << histoname <<
" found in ref file but wrong type");
123 TIter next(d->GetListOfKeys());
125 while ((key = (TKey*)next())) {
126 TObject* obj2 = key->ReadObj() ;
127 if (obj2->InheritsFrom(
"TH1")) {
128 if (obj2->GetName() == tok) {
130 B2DEBUG(20,
"Histo " << histoname <<
" found as key -> readobj");
135 if (hh1 == NULL) B2DEBUG(20,
"Histo " << histoname <<
" NOT found in ref file " << tok);
140 B2DEBUG(20,
"Histo " << histoname <<
" not in memfile or ref file");
143 TDirectory* d = gROOT;
144 TString myl = histoname;
147 while (myl.Tokenize(tok, from,
"/")) {
151 if (myl.Tokenize(dummy, f,
"/")) {
152 auto e = d->GetDirectory(tok);
154 B2DEBUG(20,
"Cd Dir " << tok);
156 }
else B2DEBUG(20,
"cd failed " << tok);
162 TObject* obj = d->FindObject(tok);
164 if (obj->IsA()->InheritsFrom(
"TH1")) {
165 B2DEBUG(20,
"Histo " << histoname <<
" found in mem");
169 B2DEBUG(20,
"Histo " << histoname <<
" NOT found in mem");
175 B2DEBUG(20,
"Histo " << histoname <<
" not found");
181 void DQMHistComparitorModule::initialize()
184 if (m_refFileName !=
"") {
185 m_refFile =
new TFile(m_refFileName.data());
188 gStyle->SetOptStat(0);
189 gStyle->SetStatStyle(1);
190 gStyle->SetOptDate(22);
193 if (!ca_current_context()) SEVCHK(ca_context_create(ca_disable_preemptive_callback),
"ca_context_create");
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));
203 TString a = it.at(2).c_str();
206 n->histo1 = it.at(0);
207 n->histo2 = it.at(1);
208 TCanvas* c =
new TCanvas(a);
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;
217 if (it.at(6) !=
"") {
218 SEVCHK(ca_create_channel(it.at(6).c_str(), NULL, NULL, 10, &n->mychid),
"ca_create_channel failure");
222 m_pnode.push_back(n);
225 SEVCHK(ca_pend_io(5.0),
"ca_pend_io failure");
228 B2DEBUG(20,
"DQMHistComparitor: initialized.");
232 void DQMHistComparitorModule::beginRun()
234 B2DEBUG(20,
"DQMHistComparitor: beginRun called.");
237 void DQMHistComparitorModule::event()
239 for (
auto& it : m_pnode) {
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);
249 if (!hist1)
continue;
250 if (!hist2)
continue;
252 B2DEBUG(20,
"Compare " << it->histo1 <<
" with ref " << it->histo2);
256 double data = hist1->KolmogorovTest(hist2,
"");
261 if (it->epicsflag) SEVCHK(ca_put(DBR_DOUBLE, it->mychid, (
void*)&data),
"ca_set failure");
264 hist2->SetLineStyle(3);
265 hist2->SetLineColor(3);
267 TIter nextkey(it->canvas->GetListOfPrimitives());
269 while ((obj = (TObject*)nextkey())) {
270 if (obj->IsA()->InheritsFrom(
"TH1")) {
271 if (
string(obj->GetName()) ==
string(hist2->GetName())) {
280 h = (TH1*)hist2->Clone();
281 if (abs(hist2->Integral()) > 0)
282 h->Scale(hist1->Integral() / hist2->Integral());
289 hist1->SetFillColor(0);
290 if (h->GetMaximum() > hist1->GetMaximum())
291 hist1->SetMaximum(1.1 * h->GetMaximum());
293 h->Draw(
"hist,same");
295 it->canvas->Pad()->SetFrameFillColor(10);
297 if (hist1->GetEntries() < it->min_entries) {
299 it->canvas->Pad()->SetFillColor(6);
301 if (data < it->error) {
302 it->canvas->Pad()->SetFillColor(2);
303 }
else if (data < it->warning) {
304 it->canvas->Pad()->SetFillColor(5);
306 it->canvas->Pad()->SetFillColor(0);
310 it->canvas->Pad()->SetFillColor(0);
312 it->canvas->Modified();
313 it->canvas->Update();
316 SEVCHK(ca_pend_io(5.0),
"ca_pend_io failure");
320 void DQMHistComparitorModule::endRun()
322 B2DEBUG(20,
"DQMHistComparitor: endRun called");
326 void DQMHistComparitorModule::terminate()
328 B2DEBUG(20,
"DQMHistComparitor: terminate called");
329 if (m_refFile)
delete m_refFile;