Belle II Software  release-05-01-25
DQMHistAnalysisInputPVSrv.cc
1 //+
2 // File : DQMHistAnalysisInputPVSrv.cc
3 // Description : DQM input module, convert epics PVs to histograms for analysis
4 //
5 // Author : B. Spruck
6 // Date : 25 - Mar - 2017
7 //-
8 
9 
10 #include <framework/core/ModuleParam.templateDetails.h>
11 #include <dqm/analysis/modules/DQMHistAnalysisInputPVSrv.h>
12 #include <TSystem.h>
13 #include <TDirectory.h>
14 #include <TH1F.h>
15 #include <TH2F.h>
16 
17 using namespace std;
18 using namespace Belle2;
19 
20 //-----------------------------------------------------------------
21 // Register the Module
22 //-----------------------------------------------------------------
23 REG_MODULE(DQMHistAnalysisInputPVSrv)
24 
25 //-----------------------------------------------------------------
26 // Implementation
27 //-----------------------------------------------------------------
28 
29 #ifdef _BELLE2_EPICS
30 static void printChidInfo(chid ichid, const char* message)
31 {
32  B2DEBUG(20, message);
33  B2DEBUG(20, "pv: " << ca_name(ichid) << " type(" << ca_field_type(ichid) << ") nelements(" << ca_element_count(
34  ichid) << ") host(" << ca_host_name(ichid)
35  << ") read(" << ca_read_access(ichid) << ") write(" << ca_write_access(ichid) << ") state(" << ca_state(ichid) << ")");
36 }
37 
38 static void exceptionCallback(struct exception_handler_args args)
39 {
40  chid ichid = args.chid;
41  long stat = args.stat; /* Channel access status code*/
42  const char* channel;
43  const char* noname = "unknown";
44 
45  channel = (ichid ? ca_name(ichid) : noname);
46 
47 
48  if (ichid) printChidInfo(ichid, "exceptionCallback");
49  printf("exceptionCallback stat %s channel %s\n", ca_message(stat), channel);
50 }
51 
52 static void connectionCallback(struct connection_handler_args args)
53 {
54  chid ichid = args.chid;
55 
56  printChidInfo(ichid, "connectionCallback");
57 }
58 
59 static void accessRightsCallback(struct access_rights_handler_args args)
60 {
61  chid ichid = args.chid;
62 
63  printChidInfo(ichid, "accessRightsCallback");
64 }
65 static void eventCallback(struct event_handler_args eha)
66 {
67  chid ichid = eha.chid;
68  MYNODE* n = (MYNODE*)eha.usr;
69 
70  if (eha.status != ECA_NORMAL) {
71  printChidInfo(ichid, "eventCallback");
72  } else {
73 // char *pdata = (char *)eha.dbr;
74 // printf("Event Callback: %s = %s (%d,%d)\n",ca_name(eha.chid),pdata,(int)eha.type,(int)eha.count);
75 // printf("Event Callback: %s (%ld,%s,%ld)\n",ca_name(n->mychid),ca_field_type(n->mychid),dbr_type_to_text(ca_field_type(n->mychid)),ca_element_count(n->mychid));
76  n->changed = true;
77  }
78 }
79 
80 #endif
81 
82 DQMHistAnalysisInputPVSrvModule::DQMHistAnalysisInputPVSrvModule()
84 {
85  //Parameter definition
86  addParam("RefreshInterval", m_interval, "Refresh interval of histograms in ms", 2000);
87  addParam("HistoList", m_histlist, "pvname, histname, histtitle, (bins,min,max[,bins,min,max])");
88  addParam("Callback", m_callback, "Using EPICS callback for changes", true);
89  addParam("Server", m_server, "Start http server on port 8082", false);
90  B2DEBUG(20, "DQMHistAnalysisInputPVSrv: Constructor done.");
91 }
92 
93 
94 DQMHistAnalysisInputPVSrvModule::~DQMHistAnalysisInputPVSrvModule()
95 {
96 #ifdef _BELLE2_EPICS
97  if (ca_current_context()) ca_context_destroy();
98 #endif
99 }
100 
102 {
103  m_eventMetaDataPtr.registerInDataStore();
104  //if (m_server) m_serv = new THttpServer("http:8082");
105 
106 #ifdef _BELLE2_EPICS
107  if (!ca_current_context()) SEVCHK(ca_context_create(ca_disable_preemptive_callback), "ca_context_create");
108  SEVCHK(ca_add_exception_event(exceptionCallback, NULL), "ca_add_exception_event");
109  for (auto& it : m_histlist) {
110  if (it.size() != 4 && it.size() != 5) {
111  B2WARNING("Histolist with wrong nr of parameters " << it.size());
112  continue;
113  }
114  auto n = (MYNODE*) callocMustSucceed(1, sizeof(MYNODE), "caMonitor");
115  pmynode.push_back(n);
116 
117  {
118  // the following code sux ... is there no root function for that?
119  TDirectory* oldDir = gDirectory;
120  TDirectory* d = oldDir;
121  TString myl = it.at(1).c_str();
122  TString tok;
123  Ssiz_t from = 0;
124  while (myl.Tokenize(tok, from, "/")) {
125  TString dummy;
126  Ssiz_t f;
127  f = from;
128  if (myl.Tokenize(dummy, f, "/")) { // check if its the last one
129  TDirectory* e;
130  e = d->GetDirectory(tok);
131  if (e) {
132  B2DEBUG(20, "Cd Dir " << tok);
133  d = e;
134  } else {
135  B2DEBUG(20, "Create Dir " << tok);
136  d = d->mkdir(tok);
137  }
138  d->cd();
139  } else {
140  break;
141  }
142  }
143 
144  B2DEBUG(20, "Create Histo " << tok);
145 
146  Int_t x;
147  Double_t xmin, xmax;
148  strncpy(n->name, it.at(0).c_str(), MAX_PV_NAME_LEN - 1);
149  istringstream is(it.at(3));
150  is >> x;
151  is >> xmin;
152  is >> xmax;
153  if (it.size() == 4) {
154  n->histo = (TH1*)new TH1F(tok, it.at(2).c_str(), x, xmin, xmax);
155  n->binx = x;
156  n->biny = 0;
157  n->binmax = x;
158  } else {
159  Int_t y;
160  Double_t ymin, ymax;
161  istringstream iss(it.at(4));
162  iss >> y;
163  iss >> ymin;
164  iss >> ymax;
165  n->histo = (TH1*)new TH2F(tok, it.at(2).c_str(), x, xmin, xmax, y, ymin, ymax);
166  n->binx = x;
167  n->biny = y;
168  n->binmax = x * y;
169  }
170 
171  // cd back to root directory
172  oldDir->cd();
173  }
174 
175  }
176 
177  for (auto n : pmynode) {
178  SEVCHK(ca_create_channel(n->name, connectionCallback, n, 20, &n->mychid), "ca_create_channel");
179  SEVCHK(ca_replace_access_rights_event(n->mychid, accessRightsCallback), "ca_replace_access_rights_event");
180  if (m_callback) {
181  SEVCHK(ca_create_subscription(DBR_STRING, 1, n->mychid, DBE_VALUE, eventCallback, n, &n->myevid), "ca_create_subscription");
182  }
183  }
184 
185 #endif
186  B2DEBUG(20, "DQMHistAnalysisInputPVSrv: initialized.");
187 }
188 
189 
191 {
192  B2DEBUG(20, "DQMHistAnalysisInputPVSrv: beginRun called.");
193 }
194 
196 {
197  m_count++;
198  m_eventMetaDataPtr.create();
199  m_eventMetaDataPtr->setExperiment(m_expno);
200  m_eventMetaDataPtr->setRun(m_runno);
201  m_eventMetaDataPtr->setEvent(m_count);
202 
203  TTimer t(m_interval, kFALSE);// in ms
204 
205 #ifdef _BELLE2_EPICS
206  SEVCHK(ca_pend_event(0.0001), "ca_pend_event");
207 
208  for (auto n : pmynode) {
209  if (m_callback && !n->changed) continue;
210  n->changed = false;
211 
212  auto bufferorg = new char[dbr_size_n(ca_field_type(n->mychid), ca_element_count(n->mychid))];
213  void* buffer = (void*) bufferorg;
214  int status;
215 
216  if (ca_field_type(n->mychid) != DBF_LONG && ca_field_type(n->mychid) != DBF_FLOAT) continue;
217  status = ca_array_get(ca_field_type(n->mychid), ca_element_count(n->mychid), n->mychid, buffer);
218  SEVCHK(status, "ca_array_get()");
219  status = ca_pend_io(15.0);
220  if (status != ECA_NORMAL) {
221  B2WARNING("EPICS ca_array_get " << ca_name(n->mychid) << " didn't return a value.");
222  } else {
223  if (!n->histo) {
224  // this should NEVER happen
225  continue;
226 // B2DEBUG(20, "Create Histo " << tok);
227 // n->histo=new TH1F(ca_name(n->mychid),ca_name(n->mychid),ca_element_count(n->mychid),0,ca_element_count(n->mychid));
228  }
229  unsigned int bins;
230  bins = ca_element_count(n->mychid) < n->binmax ? ca_element_count(n->mychid) : n->binmax;
231  TH1* histo = n->histo;
232  for (unsigned int j = ca_element_count(n->mychid); j < n->binmax; j++) histo->SetBinContent(j + 1, 0); // zero out undefined bins
233  switch (ca_field_type(n->mychid)) {
234  case DBF_CHAR: {
235  dbr_char_t* b = (dbr_char_t*)buffer;
236  for (unsigned int j = 0; j < bins; j++) {
237  histo->SetBinContent(j + 1, b[j]);
238  }
239  }; break;
240 // case DBF_INT:
241  case DBF_SHORT: { // same as INT
242  dbr_short_t* b = (dbr_short_t*)buffer;
243  for (unsigned int j = 0; j < bins; j++) {
244  histo->SetBinContent(j + 1, b[j]);
245  }
246  }; break;
247  case DBF_LONG: {
248  dbr_long_t* b = (dbr_long_t*)buffer;
249  for (unsigned int j = 0; j < bins; j++) {
250  histo->SetBinContent(j + 1, b[j]);
251  }
252  }; break;
253  case DBF_FLOAT: {
254  dbr_float_t* b = (dbr_float_t*)buffer;
255  for (unsigned int j = 0; j < bins; j++) {
256  histo->SetBinContent(j + 1, b[j]);
257  }
258  }; break;
259  case DBF_DOUBLE: {
260  dbr_double_t* b = (dbr_double_t*)buffer;
261  for (unsigned int j = 0; j < bins; j++) {
262  histo->SetBinContent(j + 1, b[j]);
263  }
264  }; break;
265  default:
266  // type not supported
267  break;
268  }
269  }
270  delete[] bufferorg;
271  }
272 #endif
273  do { // call at least once!
274  //if (m_serv) m_serv->ProcessRequests();
275  gSystem->Sleep(10); // 10 ms sleep
276  } while (!t.CheckTimer(gSystem->Now()));
277 
278 }
279 
281 {
282  B2DEBUG(20, "DQMHistAnalysisInputPVSrv: endRun called");
283 }
284 
285 
287 {
288  B2DEBUG(20, "DQMHistAnalysisInputPVSrv: terminate called");
289 }
290 
291 
292 
Belle2::DQMHistAnalysisInputPVSrvModule::m_expno
unsigned int m_expno
Exp number.
Definition: DQMHistAnalysisInputPVSrv.h:88
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::DQMHistAnalysisInputPVSrvModule::event
virtual void event() override
This method is the core of the module.
Definition: DQMHistAnalysisInputPVSrv.cc:195
Belle2::DQMHistAnalysisInputPVSrvModule::beginRun
virtual void beginRun() override
Module functions to be called from event process.
Definition: DQMHistAnalysisInputPVSrv.cc:190
Belle2::DQMHistAnalysisInputPVSrvModule::m_server
bool m_server
Whether to start http server on port 8082.
Definition: DQMHistAnalysisInputPVSrv.h:73
Belle2::DQMHistAnalysisInputPVSrvModule::m_runno
unsigned int m_runno
Run number.
Definition: DQMHistAnalysisInputPVSrv.h:90
Belle2::DQMHistAnalysisInputPVSrvModule::m_callback
bool m_callback
Whether to use EPICS callback for changes.
Definition: DQMHistAnalysisInputPVSrv.h:71
Belle2::Module::addParam
void addParam(const std::string &name, T &paramVariable, const std::string &description, const T &defaultValue)
Adds a new parameter to the module.
Definition: Module.h:562
Belle2::DQMHistAnalysisInputPVSrvModule::m_histlist
std::vector< std::vector< std::string > > m_histlist
Parameter list for histograms.
Definition: DQMHistAnalysisInputPVSrv.h:78
Belle2::DQMHistAnalysisInputPVSrvModule::m_interval
int m_interval
The refresh interval in ms.
Definition: DQMHistAnalysisInputPVSrv.h:69
Belle2::DQMHistAnalysisInputPVSrvModule::m_eventMetaDataPtr
StoreObjPtr< EventMetaData > m_eventMetaDataPtr
The metadata for each event.
Definition: DQMHistAnalysisInputPVSrv.h:85
Belle2::DQMHistAnalysisInputPVSrvModule::initialize
virtual void initialize() override
Module functions to be called from main process.
Definition: DQMHistAnalysisInputPVSrv.cc:101
Belle2::DQMHistAnalysisInputPVSrvModule::terminate
virtual void terminate() override
This method is called at the end of the event processing.
Definition: DQMHistAnalysisInputPVSrv.cc:286
Belle2::DQMHistAnalysisInputPVSrvModule::endRun
virtual void endRun() override
This method is called if the current run ends.
Definition: DQMHistAnalysisInputPVSrv.cc:280
Belle2::DQMHistAnalysisInputPVSrvModule::m_count
unsigned int m_count
Event number.
Definition: DQMHistAnalysisInputPVSrv.h:92
Belle2::DQMHistAnalysisModule
The base class for the histogram analysis module.
Definition: DQMHistAnalysis.h:27