Belle II Software  release-05-01-25
DQMHistAnalysisPXDCM.cc
1 //+
2 // File : DQMHistAnalysisPXDCM.cc
3 // Description : Analysis of PXD Common Modes
4 //
5 // Author : Bjoern Spruck, University Mainz
6 // Date : 2018
7 //-
8 
9 
10 #include <dqm/analysis/modules/DQMHistAnalysisPXDCM.h>
11 #include <TROOT.h>
12 #include <TLatex.h>
13 #include <vxd/geometry/GeoCache.h>
14 
15 using namespace std;
16 using namespace Belle2;
17 
18 //-----------------------------------------------------------------
19 // Register the Module
20 //-----------------------------------------------------------------
21 REG_MODULE(DQMHistAnalysisPXDCM)
22 
23 //-----------------------------------------------------------------
24 // Implementation
25 //-----------------------------------------------------------------
26 
29 {
30  // This module CAN NOT be run in parallel!
31 
32  //Parameter definition
33  addParam("histogramDirectoryName", m_histogramDirectoryName, "Name of Histogram dir", std::string("PXDCM"));
34  addParam("PVPrefix", m_pvPrefix, "PV Prefix", std::string("DQM:PXD:CommonMode:"));
35  addParam("useEpics", m_useEpics, "useEpics", true);
36  addParam("minEntries", m_minEntries, "minimum number of new entries for last time slot", 10000);
37  addParam("warnMeanAdhoc", m_warnMeanAdhoc, "warn level for peak position", 1.0);
38  addParam("errorMeanAdhoc", m_errorMeanAdhoc, "error level for peak position", 2.0);
39  B2DEBUG(99, "DQMHistAnalysisPXDCM: Constructor done.");
40 }
41 
42 DQMHistAnalysisPXDCMModule::~DQMHistAnalysisPXDCMModule()
43 {
44 #ifdef _BELLE2_EPICS
45  if (m_useEpics) {
46  if (ca_current_context()) ca_context_destroy();
47  }
48 #endif
49 }
50 
51 void DQMHistAnalysisPXDCMModule::initialize()
52 {
53  m_monObj = getMonitoringObject("pxd");
54  const VXD::GeoCache& geo = VXD::GeoCache::getInstance();
55 
56  // collect the list of all PXD Modules in the geometry here
57  std::vector<VxdID> sensors = geo.getListOfSensors();
58  for (VxdID& aVxdID : sensors) {
59  VXD::SensorInfoBase info = geo.getSensorInfo(aVxdID);
60  if (info.getType() != VXD::SensorInfoBase::PXD) continue;
61  m_PXDModules.push_back(aVxdID); // reorder, sort would be better
62  }
63  std::sort(m_PXDModules.begin(), m_PXDModules.end()); // back to natural order
64 
65  gROOT->cd(); // this seems to be important, or strange things happen
66 
67  m_cCommonMode = new TCanvas((m_histogramDirectoryName + "/c_CommonMode").data());
68  m_cCommonModeDelta = new TCanvas((m_histogramDirectoryName + "/c_CommonModeDelta").data());
69 
70  m_hCommonMode = new TH2D("CommonMode", "CommonMode; Module; CommonMode", m_PXDModules.size(), 0, m_PXDModules.size(), 63, 0, 63);
71  m_hCommonMode->SetDirectory(0);// dont mess with it, this is MY histogram
72  m_hCommonMode->SetStats(false);
73  m_hCommonModeDelta = new TH2D("CommonModeAdhoc", "CommonMode Adhoc; Module; CommonMode", m_PXDModules.size(), 0,
74  m_PXDModules.size(), 63, 0, 63);
75  m_hCommonModeDelta->SetDirectory(0);// dont mess with it, this is MY histogram
76  m_hCommonModeDelta->SetStats(false);
77  m_hCommonModeOld = new TH2D("CommonModeOld", "CommonMode Old; Module; CommonMode", m_PXDModules.size(), 0, m_PXDModules.size(),
78  63, 0, 63);
79  m_hCommonModeOld->SetDirectory(0);// dont mess with it, this is MY histogram
80  m_hCommonModeOld->SetStats(false);
81 
82  for (unsigned int i = 0; i < m_PXDModules.size(); i++) {
83  TString ModuleName = (std::string)m_PXDModules[i];
84  m_hCommonMode->GetXaxis()->SetBinLabel(i + 1, ModuleName);
85  m_hCommonModeDelta->GetXaxis()->SetBinLabel(i + 1, ModuleName);
86  m_hCommonModeOld->GetXaxis()->SetBinLabel(i + 1, ModuleName);
87  }
88  //Unfortunately this only changes the labels, but can't fill the bins by the VxdIDs
89  m_hCommonMode->Draw("colz");
90  m_hCommonModeDelta->Draw("colz");
91  m_hCommonModeOld->Draw("colz");
92 
93  m_monObj->addCanvas(m_cCommonMode);
94 
96  m_line1 = new TLine(0, 10, m_PXDModules.size(), 10);
97  m_line2 = new TLine(0, 16, m_PXDModules.size(), 16);
98 // m_line3 = new TLine(0, 3, m_PXDModules.size(), 3);
99  m_line1->SetHorizontal(true);
100  m_line1->SetLineColor(3);// Green
101  m_line1->SetLineWidth(3);
102  m_line2->SetHorizontal(true);
103  m_line2->SetLineColor(1);// Black
104  m_line2->SetLineWidth(3);
105 // m_line3->SetHorizontal(true);
106 // m_line3->SetLineColor(1);
107 // m_line3->SetLineWidth(3);
108 
109 
110 #ifdef _BELLE2_EPICS
111  if (m_useEpics) {
112  if (!ca_current_context()) SEVCHK(ca_context_create(ca_disable_preemptive_callback), "ca_context_create");
113  mychid.resize(4);
114  SEVCHK(ca_create_channel((m_pvPrefix + "Outside").data(), NULL, NULL, 10, &mychid[0]), "ca_create_channel failure");
115  SEVCHK(ca_create_channel((m_pvPrefix + "Status").data(), NULL, NULL, 10, &mychid[1]), "ca_create_channel failure");
116  SEVCHK(ca_create_channel((m_pvPrefix + "CM63").data(), NULL, NULL, 10, &mychid[2]), "ca_create_channel failure");
117  SEVCHK(ca_create_channel((m_pvPrefix + "Status_Adhoc").data(), NULL, NULL, 10, &mychid[3]), "ca_create_channel failure");
118 
119  for (VxdID& aPXDModule : m_PXDModules) {
120  TString buff = (std::string)aPXDModule;
121  buff.ReplaceAll(".", "_");
122  auto& my = mychid_mean[aPXDModule];
123  SEVCHK(ca_create_channel((m_pvPrefix + "Mean:" + buff).Data(), NULL, NULL, 10, &my), "ca_create_channel failure");
124  }
125  SEVCHK(ca_pend_io(5.0), "ca_pend_io failure");
126  }
127 #endif
128  B2DEBUG(99, "DQMHistAnalysisPXDCM: initialized.");
129 }
130 
131 void DQMHistAnalysisPXDCMModule::beginRun()
132 {
133  B2DEBUG(99, "DQMHistAnalysisPXDCM: beginRun called.");
134 
135  m_cCommonMode->Clear();
136  m_cCommonModeDelta->Clear();
137  m_cCommonMode->SetLogz();
138  m_cCommonModeDelta->SetLogz();
139 
140  // this is needed at least for the "Old" and "Delta" one or update doesnt work
141  m_hCommonMode->Reset();
142  m_hCommonModeDelta->Reset();
143  m_hCommonModeOld->Reset();
144 }
145 
146 void DQMHistAnalysisPXDCMModule::event()
147 {
148  double all_outside = 0.0, all = 0.0;
149  double all_cm = 0.0;
150  bool error_flag = false;
151  bool warn_flag = false;
152  bool error_adhoc_flag = false;
153  bool warn_adhoc_flag = false;
154  bool anyupdate = false;
155  if (!m_cCommonMode) return;
156  m_hCommonMode->Reset(); // dont sum up!!!
157 
158  for (unsigned int i = 0; i < m_PXDModules.size(); i++) {
159  std::string name = "PXDDAQCM2_" + (std::string)m_PXDModules[i ];
160  // std::replace( name.begin(), name.end(), '.', '_');
161 
162  TH1* hh1 = findHist(name);
163  if (hh1 == NULL) {
164  hh1 = findHist(m_histogramDirectoryName, name);
165  }
166  if (hh1) {
167  double current_full = 0.0;
168  double outside_full = 0.0;
169 
170  auto nevent = hh1->GetBinContent(0); // misuse underflow as event counter
171  bool update = nevent - m_hCommonModeOld->GetBinContent(i + 1, 0) > m_minEntries ;
172  anyupdate |= update;
173  if (update) m_hCommonModeOld->SetBinContent(i + 1, 0, nevent);
174  for (int bin = 1; bin <= 63; bin++) { // we ignore CM63!!!
175  double v;
176  v = hh1->GetBinContent(bin);
177  m_hCommonMode->SetBinContent(i + 1, bin, v); // attention, mixing bin nr and index
178  current_full += v;
179  if (nevent < m_minEntries) {
180  m_hCommonModeDelta->SetBinContent(i + 1, bin, v); // attention, mixing bin nr and index
181  } else if (update) {
182  auto old = m_hCommonModeOld->GetBinContent(i + 1, bin); // attention, mixing bin nr and index
183  m_hCommonModeDelta->SetBinContent(i + 1, bin, v - old); // attention, mixing bin nr and index
184  m_hCommonModeOld->SetBinContent(i + 1, bin, v); // attention, mixing bin nr and index
185  }
186  }
187 
189  // Attention, Integral uses the bin nr, not the value!
190  outside_full += hh1->Integral(16, 63);
191  // FIXME currently we have to much noise below the line ... thus excluding this to avoid false alarms
192  // outside_full += hh1->Integral(1 /*0*/, 5); /// FIXME we exclude bin 0 as we use it for debugging/timing pixels
193  all_outside += outside_full;
194  all += current_full;
195  double dhpc = hh1->GetBinContent(64);
196  all_cm += dhpc;
197  if (current_full > 1) {
198  error_flag |= (outside_full / current_full > 1e-5);
199  warn_flag |= (outside_full / current_full > 1e-6);
200 // error_flag |= (dhpc / current_full > 1e-5); // DHP Fifo overflow ... might be critical/unrecoverable
201 // warn_flag |= (dhpc / current_full > 1e-6); // DHP Fifo overflow ... might be critical/unrecoverable
202  }
203 
204  if (update) {
205  Double_t mean_adhoc = 0.;
206  Double_t entries_adhoc = 0.;
207  Double_t outside_adhoc = 0.;
208  // Attention, Bins
209  for (int cm_y = 0; cm_y < 16; cm_y++) {
210  auto v = m_hCommonModeDelta->GetBinContent(m_hCommonModeDelta->GetBin(i + 1, cm_y + 1));
211  entries_adhoc += v;
212  mean_adhoc += v * (cm_y + 1);
213  }
214  // Attention, Bins
215  for (int cm_y = 16; cm_y < 64; cm_y++) {
216  auto v = m_hCommonModeDelta->GetBinContent(m_hCommonModeDelta->GetBin(i + 1, cm_y + 1));
217  entries_adhoc += v;
218  outside_adhoc += v;
219  }
220  if (entries_adhoc > 0) { // ignore 1.3.2
221  mean_adhoc /= entries_adhoc; // calculate mean
222  warn_adhoc_flag |= entries_adhoc > 1000 && fabs(10.0 - mean_adhoc) > m_warnMeanAdhoc;
223  error_adhoc_flag |= entries_adhoc > 1000 && fabs(10.0 - mean_adhoc) > m_errorMeanAdhoc;
224  m_monObj->setVariable(("cm_" + (std::string)m_PXDModules[i]).c_str(), mean_adhoc);
225 #ifdef _BELLE2_EPICS
226  if (m_useEpics) {
227  auto my = mychid_mean[m_PXDModules[i]];
228  // B2ERROR("Mean "<< name << " " << mean_adhoc << " " << outside_adhoc << " " << entries_adhoc << " " <<warn_adhoc_flag <<error_adhoc_flag );
229  if (my) SEVCHK(ca_put(DBR_DOUBLE, my, (void*)&mean_adhoc), "ca_set failure");
230  }
231 #endif
232  }
233  }
234  }
235  }
236 
237  int status = 0;
238  {
239  m_cCommonMode->cd();
240  // not enough Entries
241  if (all < 100.) {
242  m_cCommonMode->Pad()->SetFillColor(kGray);// Magenta or Gray
243  status = 0; // default
244  } else {
246  if (all_outside / all > 1e-5 || /*all_cm / all > 1e-5 ||*/ error_flag) {
247  m_cCommonMode->Pad()->SetFillColor(kRed);// Red
248  status = 4;
249  } else if (all_outside / all > 1e-6 || /*all_cm / all > 1e-6 ||*/ warn_flag) {
250  m_cCommonMode->Pad()->SetFillColor(kYellow);// Yellow
251  status = 3;
252  } else if (all_outside == 0. /*&& all_cm == 0.*/) {
253  m_cCommonMode->Pad()->SetFillColor(kGreen);// Green
254  status = 2;
255  } else { // between 0 and 50 ...
256  m_cCommonMode->Pad()->SetFillColor(kWhite);// White
257  status = 1;
258  }
259  }
260 
261  if (m_hCommonMode) {
262  m_hCommonMode->Draw("colz");
263  m_line1->Draw();
264  m_line2->Draw();
265 // m_line3->Draw();
266  }
267 
268  auto tt = new TLatex(5.5, 3, "1.3.2 Module is excluded, please ignore");
269  tt->SetTextAngle(90);// Rotated
270  tt->SetTextAlign(12);// Centered
271  tt->Draw();
272 
273  m_cCommonMode->Modified();
274  m_cCommonMode->Update();
275  }
276 
277  {
278  int status_adhoc = 0;
279  m_cCommonModeDelta->cd();
280  // not enough Entries
281 
282  if (all < 100.) { // delta cannot be more than all
283  m_cCommonModeDelta->Pad()->SetFillColor(kGray);// Magenta or Gray
284  status_adhoc = 0; // default
285  } else {
287  if (error_adhoc_flag) {
288  m_cCommonModeDelta->Pad()->SetFillColor(kRed);// Red
289  status_adhoc = 4;
290  } else if (warn_adhoc_flag) {
291  m_cCommonModeDelta->Pad()->SetFillColor(kYellow);// Yellow
292  status_adhoc = 3;
293  } else {
294  m_cCommonModeDelta->Pad()->SetFillColor(kGreen);// Green
295  status_adhoc = 2;
296  /* } else { // between 0 and 50 ...
297  m_cCommonModeDelta->Pad()->SetFillColor(kWhite);// White
298  status_adhoc = 1;*/
299  }
300  }
301 #ifdef _BELLE2_EPICS
302  if (m_useEpics && anyupdate) SEVCHK(ca_put(DBR_INT, mychid[3], (void*)&status_adhoc), "ca_set failure");
303 #endif
304  if (m_hCommonModeDelta) {
305  m_hCommonModeDelta->Draw("colz");
306  m_line1->Draw();
307  m_line2->Draw();
308 // m_line3->Draw();
309  }
310 
311  auto tt = new TLatex(5.5, 3, "1.3.2 Module is excluded, please ignore");
312  tt->SetTextAngle(90);// Rotated
313  tt->SetTextAlign(12);// Centered
314  tt->Draw();
315 
316  m_cCommonModeDelta->Modified();
317  m_cCommonModeDelta->Update();
318  }
319 
320 #ifdef _BELLE2_EPICS
321  if (m_useEpics) {
322  double data = all > 0 ? (all_outside / all) : 0;
323  double data2 = all > 0 ? (all_cm / all) : 0;
324  SEVCHK(ca_put(DBR_DOUBLE, mychid[0], (void*)&data), "ca_set failure");
325  SEVCHK(ca_put(DBR_INT, mychid[1], (void*)&status), "ca_set failure");
326  SEVCHK(ca_put(DBR_DOUBLE, mychid[2], (void*)&data2), "ca_set failure");
327  SEVCHK(ca_pend_io(5.0), "ca_pend_io failure");
328  }
329 #endif
330 }
331 
332 void DQMHistAnalysisPXDCMModule::terminate()
333 {
334  B2DEBUG(99, "DQMHistAnalysisPXDCM: terminate called");
335  // should delete canvas here, maybe hist, too? Who owns it?
336 #ifdef _BELLE2_EPICS
337  if (m_useEpics) {
338  for (auto m : mychid) SEVCHK(ca_clear_channel(m), "ca_clear_channel failure");
339  for (auto& m : mychid_mean) SEVCHK(ca_clear_channel(m.second), "ca_clear_channel failure");
340  SEVCHK(ca_pend_io(5.0), "ca_pend_io failure");
341  }
342 #endif
343 }
344 
Belle2::VxdID
Class to uniquely identify a any structure of the PXD and SVD.
Definition: VxdID.h:43
Belle2::DQMHistAnalysisPXDCMModule
DQM Histogram Analysis for PXD Common Modes.
Definition: DQMHistAnalysisPXDCM.h:32
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::VXD::SensorInfoBase
Base class to provide Sensor Information for PXD and SVD.
Definition: SensorInfoBase.h:40
Belle2::VXD::GeoCache::getListOfSensors
const std::vector< VxdID > getListOfSensors() const
Get list of all sensors.
Definition: GeoCache.cc:60
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::VXD::GeoCache
Class to faciliate easy access to sensor information of the VXD like coordinate transformations or pi...
Definition: GeoCache.h:41
Belle2::VXD::GeoCache::getSensorInfo
const SensorInfoBase & getSensorInfo(Belle2::VxdID id) const
Return a referecne to the SensorInfo of a given SensorID.
Definition: GeoCache.cc:68
Belle2::DQMHistAnalysisModule
The base class for the histogram analysis module.
Definition: DQMHistAnalysis.h:27