Belle II Software  release-05-01-25
DQMHistOutputToEPICS.cc
1 //+
2 // File : DQMHistOutputToEPICS.cc
3 // Description : Write Histogram Content to EPICS Arrays
4 //
5 // Author : Bjoern Spruck, University Mainz
6 // Date : 2019
7 //-
8 
9 
10 #include <dqm/analysis/modules/DQMHistOutputToEPICS.h>
11 
12 using namespace std;
13 using namespace Belle2;
14 
15 //-----------------------------------------------------------------
16 // Register the Module
17 //-----------------------------------------------------------------
18 REG_MODULE(DQMHistOutputToEPICS)
19 
20 //-----------------------------------------------------------------
21 // Implementation
22 //-----------------------------------------------------------------
23 
26 {
27  // This module CAN NOT be run in parallel!
28 
29  //Parameter definition
30  addParam("HistoList", m_histlist, "histname, pvname");
31  B2DEBUG(99, "DQMHistOutputToEPICS: Constructor done.");
32 }
33 
34 DQMHistOutputToEPICSModule::~DQMHistOutputToEPICSModule()
35 {
36 #ifdef _BELLE2_EPICS
37  if (ca_current_context()) ca_context_destroy();
38 #endif
39 }
40 
41 void DQMHistOutputToEPICSModule::initialize()
42 {
43 #ifdef _BELLE2_EPICS
44  if (!ca_current_context()) SEVCHK(ca_context_create(ca_disable_preemptive_callback), "ca_context_create");
45 #endif
46  for (auto& it : m_histlist) {
47  if (it.size() < 2) {
48  B2WARNING("Histolist with wrong nr of parameters " << it.size());
49  continue;
50  }
51 #ifdef _BELLE2_EPICS
52  auto n = new MYNODE;
53  n->histoname = it.at(0);
54  SEVCHK(ca_create_channel(it.at(1).c_str(), NULL, NULL, 10, &n->mychid), "ca_create_channel failure");
55  if (it.size() >= 3) {
56  SEVCHK(ca_create_channel(it.at(2).c_str(), NULL, NULL, 10, &n->mychid_last), "ca_create_channel failure");
57  } else {
58  n->mychid_last = 0;
59  }
60  pmynode.push_back(n);
61 #endif
62  }
63 
64 #ifdef _BELLE2_EPICS
65  SEVCHK(ca_pend_io(5.0), "ca_pend_io failure");
66 #endif
67  B2DEBUG(99, "DQMHistOutputToEPICS: initialized.");
68 }
69 
70 void DQMHistOutputToEPICSModule::cleanPVs(void)
71 {
72 #ifdef _BELLE2_EPICS
73  for (auto* it : pmynode) {
74  if (!it->mychid) continue;
75  int length = int(ca_element_count(it->mychid));
76  if (length > 0) {
77  it->data.resize(length, 0.0);
78  SEVCHK(ca_array_put(DBR_DOUBLE, length, it->mychid, (void*)(it->data.data())), "ca_put failure");
79  if (it->mychid_last) {
80  if (length == int(ca_element_count(it->mychid_last))) {
81  SEVCHK(ca_array_put(DBR_DOUBLE, length, it->mychid_last, (void*)(it->data.data())), "ca_put failure");
82  }
83  }
84  }
85  }
86 #endif
87 }
88 
89 void DQMHistOutputToEPICSModule::beginRun()
90 {
91  B2DEBUG(99, "DQMHistOutputToEPICS: beginRun called.");
92  cleanPVs();
93  m_dirty = true;
94 }
95 
96 void DQMHistOutputToEPICSModule::event()
97 {
98 #ifdef _BELLE2_EPICS
99  for (auto& it : pmynode) {
100  if (!it->mychid) continue;
101  int length = it->data.size();
102  TH1* hh1 = findHist(it->histoname);
103  if (hh1 && hh1->GetNcells() > 2 && length > 0 && length == int(ca_element_count(it->mychid))) {
104  // If bin count doesnt match, we loose bins but otherwise ca_array_put will complain
105  // We fill up the array with ZEROs otherwise
106  if (hh1->GetDimension() == 1) {
107  int i = 0;
108  int nx = hh1->GetNbinsX() + 1;
109  for (int x = 1; x < nx && i < length ; x++) {
110  it->data[i++] = hh1->GetBinContent(x);
111  }
112 
113  } else if (hh1->GetDimension() == 2) {
114  int i = 0;
115  int nx = hh1->GetNbinsX() + 1;
116  int ny = hh1->GetNbinsY() + 1;
117  for (int y = 1; y < ny && i < length; y++) {
118  for (int x = 1; x < nx && i < length ; x++) {
119  it->data[i++] = hh1->GetBinContent(x, y);
120  }
121  }
122  }
123 
124  SEVCHK(ca_array_put(DBR_DOUBLE, length, it->mychid, (void*)it->data.data()), "ca_set failure");
125  }
126  }
127  SEVCHK(ca_pend_io(5.0), "ca_pend_io failure");
128 #endif
129 }
130 
131 void DQMHistOutputToEPICSModule::copyToLast(void)
132 {
133  if (!m_dirty) return;
134 #ifdef _BELLE2_EPICS
135  for (auto* it : pmynode) {
136  if (it->mychid && it->mychid_last) {
137  // copy PVs to last-run-PV if existing
138  int length = it->data.size();
139  if (length > 0 && length == int(ca_element_count(it->mychid_last))) {
140  SEVCHK(ca_array_put(DBR_DOUBLE, length, it->mychid_last, (void*)it->data.data()), "ca_put failure");
141  }
142  }
143  }
144  SEVCHK(ca_pend_io(5.0), "ca_pend_io failure");
145 #endif
146 
147  m_dirty = false;
148 }
149 
150 void DQMHistOutputToEPICSModule::endRun()
151 {
152  B2DEBUG(99, "DQMHistOutputToEPICS: endRun called");
153  copyToLast();
154 }
155 
156 void DQMHistOutputToEPICSModule::terminate()
157 {
158  B2DEBUG(99, "DQMHistOutputToEPICS: terminate called");
159  copyToLast();
160  // the following belongs to terminate
161 #ifdef _BELLE2_EPICS
162  for (auto* it : pmynode) {
163  if (it->mychid) SEVCHK(ca_clear_channel(it->mychid), "ca_clear_channel failure");
164  if (it->mychid_last) SEVCHK(ca_clear_channel(it->mychid_last), "ca_clear_channel failure");
165  }
166  SEVCHK(ca_pend_io(5.0), "ca_pend_io failure");
167 #endif
168 }
169 
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::DQMHistOutputToEPICSModule
Write DQM Histogram Content to EPICS Arrays.
Definition: DQMHistOutputToEPICS.h:28
Belle2::DQMHistAnalysisModule
The base class for the histogram analysis module.
Definition: DQMHistAnalysis.h:27