Belle II Software  release-05-01-25
DQMHistAnalysisPXDDAQ.cc
1 //+
2 // File : DQMHistAnalysisPXDDAQ.cc
3 // Description : Analysis of PXD Reduction
4 //
5 // Author : Bjoern Spruck, Univerisity Mainz
6 // Date : 2020
7 //-
8 
9 
10 #include <dqm/analysis/modules/DQMHistAnalysisPXDDAQ.h>
11 #include <TROOT.h>
12 
13 using namespace std;
14 using namespace Belle2;
15 
16 //-----------------------------------------------------------------
17 // Register the Module
18 //-----------------------------------------------------------------
19 REG_MODULE(DQMHistAnalysisPXDDAQ)
20 
21 //-----------------------------------------------------------------
22 // Implementation
23 //-----------------------------------------------------------------
24 
27 {
28  // This module CAN NOT be run in parallel!
29 
30  //Parameter definition
31  addParam("histogramDirectoryName", m_histogramDirectoryName, "Name of Histogram dir", std::string("PXDDAQ"));
32  addParam("PVPrefix", m_pvPrefix, "PV Prefix", std::string("DQM:PXD:DAQ:"));
33  addParam("useEpics", m_useEpics, "useEpics", true);
34  B2DEBUG(1, "DQMHistAnalysisPXDDAQ: Constructor done.");
35 
36 }
37 
38 DQMHistAnalysisPXDDAQModule::~DQMHistAnalysisPXDDAQModule()
39 {
40 #ifdef _BELLE2_EPICS
41  if (m_useEpics) {
42  if (ca_current_context()) ca_context_destroy();
43  }
44 #endif
45 }
46 
47 void DQMHistAnalysisPXDDAQModule::initialize()
48 {
49  B2DEBUG(1, "DQMHistAnalysisPXDDAQ: initialized.");
50 
51  m_monObj = getMonitoringObject("pxd");
52 
53  gROOT->cd(); // this seems to be important, or strange things happen
54 
55  m_cMissingDHC = new TCanvas((m_histogramDirectoryName + "/c_MissingDHC").data());
56  m_cMissingDHE = new TCanvas((m_histogramDirectoryName + "/c_MissingDHE").data());
57  m_cMissingDHP = new TCanvas((m_histogramDirectoryName + "/c_MissingDHP").data());
58  m_cStatistic = new TCanvas((m_histogramDirectoryName + "/c_Statistic").data());
59 
60  m_monObj->addCanvas(m_cMissingDHC);
61  m_monObj->addCanvas(m_cMissingDHE);
62  m_monObj->addCanvas(m_cMissingDHP);
63  m_monObj->addCanvas(m_cStatistic);
64 
65  m_hMissingDHC = new TH2F("hPXDMissingDHC", "hPXDMissingDHC", 16, 0, 16, 2, 0, 2);
66  m_hMissingDHE = new TH2F("hPXDMissingDHE", "hPXDMissingDHE", 64, 0, 64, 2, 0, 2);
67 
68 #ifdef _BELLE2_EPICS
69  mychid.resize(20);
70  if (m_useEpics) {
71  if (!ca_current_context()) SEVCHK(ca_context_create(ca_disable_preemptive_callback), "ca_context_create");
72  SEVCHK(ca_create_channel((m_pvPrefix + "HLTRej").data(), NULL, NULL, 10, &mychid[0]), "ca_create_channel failure");
73  SEVCHK(ca_create_channel((m_pvPrefix + "Trunc").data(), NULL, NULL, 10, &mychid[1]), "ca_create_channel failure");
74  SEVCHK(ca_create_channel((m_pvPrefix + "HER_Trunc").data(), NULL, NULL, 10, &mychid[2]), "ca_create_channel failure");
75  SEVCHK(ca_create_channel((m_pvPrefix + "LER_Trunc").data(), NULL, NULL, 10, &mychid[3]), "ca_create_channel failure");
76  SEVCHK(ca_create_channel((m_pvPrefix + "CM63").data(), NULL, NULL, 10, &mychid[4]), "ca_create_channel failure");
77  SEVCHK(ca_create_channel((m_pvPrefix + "HER_CM63").data(), NULL, NULL, 10, &mychid[5]), "ca_create_channel failure");
78  SEVCHK(ca_create_channel((m_pvPrefix + "LER_CM63").data(), NULL, NULL, 10, &mychid[6]), "ca_create_channel failure");
79  SEVCHK(ca_create_channel((m_pvPrefix + "HER_CM63_1ms").data(), NULL, NULL, 10, &mychid[7]), "ca_create_channel failure");
80  SEVCHK(ca_create_channel((m_pvPrefix + "LER_CM63_1ms").data(), NULL, NULL, 10, &mychid[8]), "ca_create_channel failure");
81  SEVCHK(ca_create_channel((m_pvPrefix + "HER_Trunc_1ms").data(), NULL, NULL, 10, &mychid[9]), "ca_create_channel failure");
82  SEVCHK(ca_create_channel((m_pvPrefix + "LER_Trunc_1ms").data(), NULL, NULL, 10, &mychid[10]), "ca_create_channel failure");
83  SEVCHK(ca_create_channel((m_pvPrefix + "MissFrame").data(), NULL, NULL, 10, &mychid[11]), "ca_create_channel failure");
84  SEVCHK(ca_create_channel((m_pvPrefix + "Timeout").data(), NULL, NULL, 10, &mychid[12]), "ca_create_channel failure");
85  SEVCHK(ca_create_channel((m_pvPrefix + "LinkDown").data(), NULL, NULL, 10, &mychid[13]), "ca_create_channel failure");
86  SEVCHK(ca_create_channel((m_pvPrefix + "Mismatch").data(), NULL, NULL, 10, &mychid[14]), "ca_create_channel failure");
87  SEVCHK(ca_create_channel((m_pvPrefix + "HER_Miss").data(), NULL, NULL, 10, &mychid[15]), "ca_create_channel failure");
88  SEVCHK(ca_create_channel((m_pvPrefix + "LER_Miss").data(), NULL, NULL, 10, &mychid[16]), "ca_create_channel failure");
89  SEVCHK(ca_create_channel((m_pvPrefix + "HER_Miss_1ms").data(), NULL, NULL, 10, &mychid[17]), "ca_create_channel failure");
90  SEVCHK(ca_create_channel((m_pvPrefix + "LER_Miss_1ms").data(), NULL, NULL, 10, &mychid[18]), "ca_create_channel failure");
91  SEVCHK(ca_create_channel((m_pvPrefix + "unused").data(), NULL, NULL, 10, &mychid[19]), "ca_create_channel failure");
92  SEVCHK(ca_pend_io(5.0), "ca_pend_io failure");
93  }
94 #endif
95 }
96 
97 void DQMHistAnalysisPXDDAQModule::beginRun()
98 {
99  B2DEBUG(1, "DQMHistAnalysisPXDDAQ: beginRun called.");
100 
101  m_cMissingDHP->Clear();
102 }
103 
104 void DQMHistAnalysisPXDDAQModule::event()
105 {
106 // double data = 0.0;
107  if (m_cMissingDHP == nullptr || m_cMissingDHE == nullptr || m_cMissingDHC == nullptr
108  || m_cStatistic == nullptr) return; // we could assume this
109 
110  m_cMissingDHP->Clear();
111  m_cMissingDHE->Clear();
112  m_cMissingDHC->Clear();
113  m_cStatistic->Clear();
114 
115 
116  {
117  // DHC histogram
118  std::string name = "PXDDAQDHCError";
119 
120  TH1* hh1 = findHist(name);
121  if (hh1 == NULL) {
122  hh1 = findHist(m_histogramDirectoryName, name);
123  }
124  m_cMissingDHC->cd();
125  if (hh1) {
126  auto events = hh1->GetBinContent(hh1->GetBin(-1, -1));
127  // first, we have to relate the per-DHC overflow (DHC object count) to the overall overflow (event count)
128  // second, we have to relate the "fake data" DHC bin to the per-DHC overflow (DHC object count)
129  m_hMissingDHC->Reset();
130  for (int i = 0; i < 16; i++) {
131  auto dhecount = hh1->GetBinContent(hh1->GetBin(i, -1));
132  if (events > 0) m_hMissingDHC->Fill((double)i, 0.0, 1.0 - dhecount / events);
133  // c_FAKE_NO_DATA_TRIG = 1ull << 29,
134  if (dhecount > 0) m_hMissingDHC->Fill((double)i, 1.0, hh1->GetBinContent(hh1->GetBin(i, 29) / dhecount));
135  }
136  m_hMissingDHC->Draw("text");
137  }
138  }
139 
140  {
141  // DHE histogram
142  std::string name = "PXDDAQDHEError";
143 
144  TH1* hh1 = findHist(name);
145  if (hh1 == NULL) {
146  hh1 = findHist(m_histogramDirectoryName, name);
147  }
148  m_cMissingDHE->cd();
149  if (hh1) {
150  auto events = hh1->GetBinContent(hh1->GetBin(-1, -1));
151  // first, we have to relate the per-DHE overflow (DHE object count) to the overall overflow (event count)
152  // second, we have to relate the "fake data" DHE bin to the per-DHE overflow (DHE object count)
153  m_hMissingDHE->Reset();
154  for (int i = 0; i < 64; i++) {
155  auto dhecount = hh1->GetBinContent(hh1->GetBin(i, -1));
156  if (events > 0) m_hMissingDHE->Fill((double)i, 0.0, 1.0 - dhecount / events);
157  // c_FAKE_NO_DATA_TRIG = 1ull << 29,
158  if (dhecount > 0) m_hMissingDHE->Fill((double)i, 1.0, hh1->GetBinContent(hh1->GetBin(i, 29) / dhecount));
159  }
160  m_hMissingDHE->Draw("text");
161  }
162  }
163 
164  {
165  // DHP histogram
166  if (m_hMissingDHP) { delete m_hMissingDHP; m_hMissingDHP = nullptr;}
167 
168  std::string name = "PXDDAQDHPDataMissing";
169 
170  TH1* hh1 = findHist(name);
171  if (hh1 == NULL) {
172  hh1 = findHist(m_histogramDirectoryName, name);
173  }
174  m_cMissingDHP->cd();
175  if (hh1) {
176  m_hMissingDHP = (TH1F*)hh1->DrawClone("text");
177  if (m_hMissingDHP->GetBinContent(0)) {
178  m_hMissingDHP->Scale(1.0 / m_hMissingDHP->GetBinContent(0));
179  m_hMissingDHP->Draw("text");
180  }
181  }
182  // double data = m_hMissingDHP->Max???;
183  //
184  // m_monObj->setVariable("missingDHPFraction", data);
185  //
186  m_cMissingDHP->Modified();
187  m_cMissingDHP->Update();
188  }
189 
190  double data_HLTRej = 0.0;
191  double data_Trunc = 0.0;
192  double data_HER_Trunc = 0.0;
193  double data_LER_Trunc = 0.0;
194  double data_CM63 = 0.0;
195  double data_HER_CM63 = 0.0;
196  double data_LER_CM63 = 0.0;
197  double data_HER_CM63_1ms = 0.0;
198  double data_LER_CM63_1ms = 0.0;
199  double data_HER_Trunc_1ms = 0.0;
200  double data_LER_Trunc_1ms = 0.0;
201  double data_MissFrame = 0.0;
202  double data_Timeout = 0.0;
203  double data_LinkDown = 0.0;
204  double data_Mismatch = 0.0;
205  double data_HER_Miss = 0.0;
206  double data_LER_Miss = 0.0;
207  double data_HER_Miss_1ms = 0.0;
208  double data_LER_Miss_1ms = 0.0;
209  double data_unused = 0.0;
210 
211  // Stat histogram
212  if (m_hStatistic) { delete m_hStatistic; m_hStatistic = nullptr;}
213 
214  std::string name = "PXDDAQStat";
215 
216  TH1* hh1 = findHist(name);
217  if (hh1 == NULL) {
218  hh1 = findHist(m_histogramDirectoryName, name);
219  }
220  m_cStatistic->cd();
221  if (hh1) {
222  m_hStatistic = (TH1F*)hh1->DrawClone("text");
223  if (m_hStatistic->GetBinContent(0)) {
224  m_hStatistic->Scale(1.0 / m_hStatistic->GetBinContent(0));
225  m_hStatistic->Draw("text");
226  }
227  double scale = hh1->GetBinContent(0);// underflow is event counter
228  if (scale != 0.0) scale = 1.0 / scale; // just avoid dive by zero
229  data_HLTRej = hh1->GetBinContent(1 + 0) * scale;
230  data_Trunc = hh1->GetBinContent(1 + 1) * scale;
231  data_HER_Trunc = hh1->GetBinContent(1 + 2) * scale;
232  data_LER_Trunc = hh1->GetBinContent(1 + 3) * scale;
233  data_CM63 = hh1->GetBinContent(1 + 4) * scale;
234  data_HER_CM63 = hh1->GetBinContent(1 + 5) * scale;
235  data_LER_CM63 = hh1->GetBinContent(1 + 6) * scale;
236  data_HER_CM63_1ms = hh1->GetBinContent(1 + 7) * scale;
237  data_LER_CM63_1ms = hh1->GetBinContent(1 + 8) * scale;
238  data_HER_Trunc_1ms = hh1->GetBinContent(1 + 9) * scale;
239  data_LER_Trunc_1ms = hh1->GetBinContent(1 + 10) * scale;
240  data_MissFrame = hh1->GetBinContent(1 + 11) * scale;
241  data_Timeout = hh1->GetBinContent(1 + 12) * scale;
242  data_LinkDown = hh1->GetBinContent(1 + 13) * scale;
243  data_Mismatch = hh1->GetBinContent(1 + 14) * scale;
244  data_HER_Miss = hh1->GetBinContent(1 + 15) * scale;
245  data_LER_Miss = hh1->GetBinContent(1 + 16) * scale;
246  data_HER_Miss_1ms = hh1->GetBinContent(1 + 17) * scale;
247  data_LER_Miss_1ms = hh1->GetBinContent(1 + 18) * scale;
248  data_unused = hh1->GetBinContent(1 + 19) * scale;
249  }
250 
251  m_monObj->setVariable("HLTReject", data_HLTRej);
252  m_monObj->setVariable("Trunc", data_Trunc);
253  m_monObj->setVariable("HER_Trunc", data_HER_Trunc);
254  m_monObj->setVariable("LER_Trunc", data_LER_Trunc);
255  m_monObj->setVariable("CM63", data_CM63);
256  m_monObj->setVariable("HER_CM63", data_HER_CM63);
257  m_monObj->setVariable("LER_CM63", data_LER_CM63);
258  m_monObj->setVariable("HER_CM63_1ms", data_HER_CM63_1ms);
259  m_monObj->setVariable("LER_CM63_1ms", data_LER_CM63_1ms);
260  m_monObj->setVariable("HER_Trunc_1ms", data_HER_Trunc_1ms);
261  m_monObj->setVariable("LER_Trunc_1ms", data_LER_Trunc_1ms);
262  m_monObj->setVariable("MissFrame", data_MissFrame);
263  m_monObj->setVariable("Timeout", data_Timeout);
264  m_monObj->setVariable("LinkDown", data_LinkDown);
265  m_monObj->setVariable("Mismatch", data_Mismatch);
266  m_monObj->setVariable("HER_Miss", data_HER_Miss);
267  m_monObj->setVariable("LER_Miss", data_LER_Miss);
268  m_monObj->setVariable("HER_Miss_1ms", data_HER_Miss_1ms);
269  m_monObj->setVariable("LER_Miss_1ms", data_LER_Miss_1ms);
270 
271  m_cStatistic->Modified();
272  m_cStatistic->Update();
273 
274 #ifdef _BELLE2_EPICS
275  if (m_useEpics) {
276  SEVCHK(ca_put(DBR_DOUBLE, mychid[0], (void*)&data_HLTRej), "ca_set failure");
277  SEVCHK(ca_put(DBR_DOUBLE, mychid[1], (void*)&data_Trunc), "ca_set failure");
278  SEVCHK(ca_put(DBR_DOUBLE, mychid[2], (void*)&data_HER_Trunc), "ca_set failure");
279  SEVCHK(ca_put(DBR_DOUBLE, mychid[3], (void*)&data_LER_Trunc), "ca_set failure");
280  SEVCHK(ca_put(DBR_DOUBLE, mychid[4], (void*)&data_CM63), "ca_set failure");
281  SEVCHK(ca_put(DBR_DOUBLE, mychid[5], (void*)&data_HER_CM63), "ca_set failure");
282  SEVCHK(ca_put(DBR_DOUBLE, mychid[6], (void*)&data_LER_CM63), "ca_set failure");
283  SEVCHK(ca_put(DBR_DOUBLE, mychid[7], (void*)&data_HER_CM63_1ms), "ca_set failure");
284  SEVCHK(ca_put(DBR_DOUBLE, mychid[8], (void*)&data_LER_CM63_1ms), "ca_set failure");
285  SEVCHK(ca_put(DBR_DOUBLE, mychid[9], (void*)&data_HER_Trunc_1ms), "ca_set failure");
286  SEVCHK(ca_put(DBR_DOUBLE, mychid[10], (void*)&data_LER_Trunc_1ms), "ca_set failure");
287  SEVCHK(ca_put(DBR_DOUBLE, mychid[11], (void*)&data_MissFrame), "ca_set failure");
288  SEVCHK(ca_put(DBR_DOUBLE, mychid[12], (void*)&data_Timeout), "ca_set failure");
289  SEVCHK(ca_put(DBR_DOUBLE, mychid[13], (void*)&data_LinkDown), "ca_set failure");
290  SEVCHK(ca_put(DBR_DOUBLE, mychid[14], (void*)&data_Mismatch), "ca_set failure");
291  SEVCHK(ca_put(DBR_DOUBLE, mychid[15], (void*)&data_HER_Miss), "ca_set failure");
292  SEVCHK(ca_put(DBR_DOUBLE, mychid[16], (void*)&data_LER_Miss), "ca_set failure");
293  SEVCHK(ca_put(DBR_DOUBLE, mychid[17], (void*)&data_HER_Miss_1ms), "ca_set failure");
294  SEVCHK(ca_put(DBR_DOUBLE, mychid[18], (void*)&data_LER_Miss_1ms), "ca_set failure");
295  SEVCHK(ca_put(DBR_DOUBLE, mychid[19], (void*)&data_unused), "ca_set failure");
296  // write out
297  SEVCHK(ca_pend_io(5.0), "ca_pend_io failure");
298  }
299 #endif
300 }
301 
302 void DQMHistAnalysisPXDDAQModule::terminate()
303 {
304  B2DEBUG(1, "DQMHistAnalysisPXDDAQ: terminate called");
305  // should delete canvas here, maybe hist, too? Who owns it?
306 }
307 
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::DQMHistAnalysisPXDDAQModule
DQM Histogram Analysis for PXD DAQ.
Definition: DQMHistAnalysisPXDDAQ.h:29
Belle2::DQMHistAnalysisModule
The base class for the histogram analysis module.
Definition: DQMHistAnalysis.h:27