Belle II Software  release-05-01-25
DQMHistAnalysisIP.cc
1 //+
2 // File : DQMHistAnalysisIP.cc
3 // Description : Mean for IP position with delta histogramming
4 //
5 // Author : Bjoern Spruck, Mainz Univerisity
6 // Date : 2019
7 //-
8 
9 
10 #include <dqm/analysis/modules/DQMHistAnalysisIP.h>
11 #include <TROOT.h>
12 #include <TClass.h>
13 
14 using namespace std;
15 using namespace Belle2;
16 
17 //-----------------------------------------------------------------
18 // Register the Module
19 //-----------------------------------------------------------------
20 REG_MODULE(DQMHistAnalysisIP)
21 
22 //-----------------------------------------------------------------
23 // Implementation
24 //-----------------------------------------------------------------
25 
28 {
29  // This module CAN NOT be run in parallel!
30 
31  //Parameter definition
32  addParam("HistoName", m_histoname, "Name of Histogram (incl dir)", std::string(""));
33  addParam("PVName", m_pvPrefix, "PV Prefix", std::string("DQM:TEST:hist:"));
34  addParam("MonitorPrefix", m_monPrefix, "Monitor Prefix");// force to be set!
35  addParam("useEpics", m_useEpics, "useEpics", true);
36  addParam("minEntries", m_minEntries, "minimum number of new Entries for a fit", 1000);
37  B2DEBUG(20, "DQMHistAnalysisIP: Constructor done.");
38 }
39 
40 
41 DQMHistAnalysisIPModule::~DQMHistAnalysisIPModule()
42 {
43 #ifdef _BELLE2_EPICS
44  if (m_useEpics) {
45  if (ca_current_context()) ca_context_destroy();
46  }
47 #endif
48 }
49 
50 void DQMHistAnalysisIPModule::initialize()
51 {
52  B2DEBUG(20, "DQMHistAnalysisIP: initialized.");
53 
54  m_monObj = getMonitoringObject("ip");
55 
56  TString a;
57  a = m_histoname;
58  m_c1 = new TCanvas(a + "_fit"); // prefer to change canvas name to monitorPrefix, but then chnages on the web gui are needed :-(
59 
60  m_line = new TLine(0, 10, 0, 0);
61  m_line->SetVertical(true);
62  m_line->SetLineColor(8);
63  m_line->SetLineWidth(3);
64 
65  m_monObj->addCanvas(m_c1);
66 
67  // need the function to get parameter names
68 #ifdef _BELLE2_EPICS
69  if (m_useEpics) {
70  if (!ca_current_context()) SEVCHK(ca_context_create(ca_disable_preemptive_callback), "ca_context_create");
71  std::string aa;
72  aa = m_pvPrefix + "Mean";
73  SEVCHK(ca_create_channel(aa.c_str(), NULL, NULL, 10, &mychid[0]), "ca_create_channel failure");
74  aa = m_pvPrefix + "RMS";
75  SEVCHK(ca_create_channel(aa.c_str(), NULL, NULL, 10, &mychid[1]), "ca_create_channel failure");
76  // Read LO and HI limits from EPICS, seems this needs additional channels?
77  // SEVCHK(ca_get(DBR_DOUBLE,mychid[i],(void*)&data),"ca_get failure"); // data is only valid after ca_pend_io!!
78  SEVCHK(ca_pend_io(5.0), "ca_pend_io failure");
79  }
80 #endif
81 }
82 
83 
84 void DQMHistAnalysisIPModule::beginRun()
85 {
86  //m_serv->SetTimer(100, kFALSE);
87  B2DEBUG(20, "DQMHistAnalysisIP: beginRun called.");
88  m_c1->Clear();
89 
90  TH1* hh1;
91  hh1 = findHist(m_histoname.c_str());
92 
93  if (hh1 == NULL) {
94  B2DEBUG(20, "Histo " << m_histoname << " not in memfile");
95  // the following code sux ... is there no root function for that?
96  TDirectory* d = gROOT;
97  TString myl = m_histoname;
98  TString tok;
99  Ssiz_t from = 0;
100  while (myl.Tokenize(tok, from, "/")) {
101  TString dummy;
102  Ssiz_t f;
103  f = from;
104  if (myl.Tokenize(dummy, f, "/")) { // check if its the last one
105  auto e = d->GetDirectory(tok);
106  if (e) {
107  B2DEBUG(20, "Cd Dir " << tok);
108  d = e;
109  }
110  d->cd();
111  } else {
112  break;
113  }
114  }
115  TObject* obj = d->FindObject(tok);
116  if (obj != NULL) {
117  if (obj->IsA()->InheritsFrom("TH1")) {
118  B2DEBUG(20, "Histo " << m_histoname << " found in mem");
119  hh1 = (TH1*)obj;
120  }
121  } else {
122  B2DEBUG(20, "Histo " << m_histoname << " NOT found in mem");
123  }
124  }
125 
126  if (hh1 != NULL) {
127  m_c1->cd();
128  hh1->Draw();
129  m_line->Draw();
130  } else {
131  B2DEBUG(20, "Histo " << m_histoname << " not found");
132  }
133 
134  if (m_h_last) m_h_last->Reset();
135 }
136 
137 void DQMHistAnalysisIPModule::event()
138 {
139  TH1* hh1;
140  bool flag = false;
141 
142  hh1 = findHist(m_histoname.c_str());
143  if (hh1 == NULL) {
144  B2DEBUG(20, "Histo " << m_histoname << " not in memfile");
145  // the following code sux ... is there no root function for that?
146  TDirectory* d = gROOT;
147  TString myl = m_histoname;
148  TString tok;
149  Ssiz_t from = 0;
150  while (myl.Tokenize(tok, from, "/")) {
151  TString dummy;
152  Ssiz_t f;
153  f = from;
154  if (myl.Tokenize(dummy, f, "/")) { // check if its the last one
155  auto e = d->GetDirectory(tok);
156  if (e) {
157  B2DEBUG(20, "Cd Dir " << tok);
158  d = e;
159  }
160  d->cd();
161  } else {
162  break;
163  }
164  }
165  TObject* obj = d->FindObject(tok);
166  if (obj != NULL) {
167  if (obj->IsA()->InheritsFrom("TH1")) {
168  B2DEBUG(20, "Histo " << m_histoname << " found in mem");
169  hh1 = (TH1*)obj;
170  flag = true;
171  }
172  } else {
173  B2DEBUG(20, "Histo " << m_histoname << " NOT found in mem");
174  }
175  }
176  if (hh1 != NULL) {
177 
178  if (!m_h_last) {
179  B2DEBUG(20, "Initial Clone");
180  gROOT->cd();
181  m_h_last = (TH1*)hh1->Clone();
182  m_h_last->SetName(TString(hh1->GetName()) + "_last");
183  m_h_last->Reset();
184  }
185 
186  double last = m_h_last->GetEntries();
187  double current = hh1->GetEntries();
188  B2DEBUG(20, "Entries: " << last << "," << current);
189  if (current - last >= m_minEntries) {
190  B2DEBUG(20, "Update Delta & Fit");
191  gROOT->cd();
192  TH1* delta = (TH1*)hh1->Clone();
193  delta->Add(m_h_last, -1.);
194 
195  m_c1->Clear();
196  m_c1->cd();// necessary!
197  delta->Draw("hist");
198  m_h_last->Reset();
199  m_h_last->Add(hh1);
200 
201  double x = delta->GetMean();// must be double bc of EPICS below
202  double w = delta->GetRMS();// must be double bc of EPICS below
203  double y1 = delta->GetMaximum();
204  double y2 = delta->GetMinimum();
205  B2DEBUG(20, "Fit " << x << "," << w << "," << y1 << "," << y2);
206  m_line->SetY1(y1 + (y1 - y2) * 0.05);
207  m_line->SetX1(x);
208  m_line->SetX2(x);
209  delta->GetXaxis()->SetRangeUser(x - 3 * w, x + 3 * w);
210  if (!flag) {
211  // dont add another line...
212  m_line->Draw();
213  }
214  m_c1->Modified();
215  m_c1->Update();
216 
217  m_monObj->setVariable(m_monPrefix + "_mean", x);
218  m_monObj->setVariable(m_monPrefix + "_width", w);
219 
220 #ifdef _BELLE2_EPICS
221  if (m_useEpics) {
222  B2INFO("Update EPICS");
223  if (mychid[0]) SEVCHK(ca_put(DBR_DOUBLE, mychid[0], (void*)&x), "ca_set failure");
224  if (mychid[1]) SEVCHK(ca_put(DBR_DOUBLE, mychid[1], (void*)&w), "ca_set failure");
225  SEVCHK(ca_pend_io(5.0), "ca_pend_io failure");
226  }
227 #endif
228  }
229  } else {
230  B2DEBUG(20, "Histo " << m_histoname << " not found");
231  }
232 
233 }
234 
235 void DQMHistAnalysisIPModule::terminate()
236 {
237 #ifdef _BELLE2_EPICS
238  if (m_useEpics) {
239  // cppcheck-suppress knownConditionTrueFalse
240  if (m_parameters > 0) {
241  for (auto i = 0; i < m_parameters; i++) {
242  if (mychid[i]) SEVCHK(ca_clear_channel(mychid[i]), "ca_clear_channel failure");
243  }
244  SEVCHK(ca_pend_io(5.0), "ca_pend_io failure");
245  }
246  }
247 #endif
248  B2DEBUG(20, "DQMHistAnalysisIP: terminate called");
249 }
250 
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
Belle2::DQMHistAnalysisIPModule
Class definition for the output module of Sequential ROOT I/O.
Definition: DQMHistAnalysisIP.h:30
Belle2::DQMHistAnalysisModule
The base class for the histogram analysis module.
Definition: DQMHistAnalysis.h:27