Belle II Software  release-05-02-19
PXDEventPlotModule.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2018 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Bjoern Spruck *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include <pxd/modules/pxdHelper/PXDEventPlotModule.h>
12 #include <TBox.h>
13 #include <TCanvas.h>
14 #include <TStyle.h>
15 #include <TH2.h>
16 
17 
18 using namespace std;
19 using namespace Belle2;
20 using namespace Belle2::PXD;
21 
22 //-----------------------------------------------------------------
23 // Register the Module
24 //-----------------------------------------------------------------
25 REG_MODULE(PXDEventPlot)
26 
27 //-----------------------------------------------------------------
28 // Implementation
29 //-----------------------------------------------------------------
30 
31 PXDEventPlotModule::PXDEventPlotModule() : Module(), m_vxdGeometry(VXD::GeoCache::getInstance())
32 {
33  //Set module properties
34  setDescription("Plot Events on PXD Hit/Charge Maps and write pictures");
35  setPropertyFlags(c_ParallelProcessingCertified);
36  addParam("PXDRawHitsName", m_storeRawHitsName, "The name of the StoreArray of PXDRawHits to be processed", string(""));
37  addParam("GatedMode", m_gateModeFlag, "Extra histograms for gated mode", false);
38 }
39 
40 void PXDEventPlotModule::initialize()
41 {
42  m_eventMetaData.isRequired();
43  m_storeRawHits.isRequired(m_storeRawHitsName);
44  m_rawTTD.isRequired();
45  std::vector<VxdID> sensors = m_vxdGeometry.getListOfSensors();
46  for (VxdID& avxdid : sensors) {
47  VXD::SensorInfoBase info = m_vxdGeometry.getSensorInfo(avxdid);
48  if (info.getType() != VXD::SensorInfoBase::PXD) continue;
49 
50  m_histos[avxdid] = new TH2F("Chargemap for " + TString((std::string)avxdid),
51  "PXD Module Chargemap for " + TString((std::string)avxdid) + " ;VCell;UCell", 768, 0, 768, 250, 0, 256);
52  m_histos[avxdid]->SetContour(100);
53  if (m_gateModeFlag) {
54  m_histos_gm[avxdid] = new TH2F("Chargemap for " + TString((std::string)avxdid),
55  "PXD Module Chargemap for GM " + TString((std::string)avxdid) + " ;VCell-GM;UCell", 768, 0, 768, 250, 0, 256);
56  m_histos_gm[avxdid]->SetContour(100);
57  m_histos_gm2[avxdid] = new TH2F("Chargemap for " + TString((std::string)avxdid),
58  "PXD Module Chargemap for GM2 " + TString((std::string)avxdid) + " ;VCell+GM;UCell", 768, 0, 768, 250, 0, 256);
59  m_histos_gm2[avxdid]->SetContour(100);
60  }
61  }
62  m_c = new TCanvas("c1", "c1", 4000, 3000);
63  m_l1 = new TLine(0, 62, 768, 62);
64  m_l2 = new TLine(0, 125, 768, 125);
65  m_l3 = new TLine(0, 187, 768, 187);
66  m_l1->SetLineColor(kMagenta);
67  m_l2->SetLineColor(kMagenta);
68  m_l3->SetLineColor(kMagenta);
69 }
70 
71 void PXDEventPlotModule::event()
72 {
73  unsigned int evtNr = m_eventMetaData->getEvent();
74  unsigned int evtRun = m_eventMetaData->getRun();
75 
76  for (auto h : m_histos) {
77  if (h.second) h.second->Reset();
78  }
79 
80  std::map <VxdID, bool> cm_map;
81  bool cm_flag = false;
82  for (auto& pkt : *m_storeDAQEvtStats) {
83  for (auto& dhc : pkt) {
84  for (auto& dhe : dhc) {
85  for (auto cm = dhe.cm_begin(); cm < dhe.cm_end(); ++cm) {
86  // uint8_t, uint16_t, uint8_t ; tuple of Chip ID (2 bit), Row (10 bit), Common Mode (6 bit)
87  if (std::get<2>(*cm) == 63) {
88  m_histos[dhe.getSensorID()]->Fill(std::get<1>(*cm), 252 + std::get<0>(*cm), std::get<2>(*cm));
89  cm_flag = true;
90  cm_map[dhe.getSensorID()] = true;
91  }
92  }
93  }
94  }
95  }
96 
97 // if( !cm_flag ) return;
98 
99  int tinj = -1, rgate = -1;
100  for (auto& it : m_rawTTD) {
101  // get last injection time
102  auto difference = it.GetTimeSinceLastInjection(0);
103  // check time overflow, too long ago
104  if (difference != 0x7FFFFFFF) {
106  float diff2 = difference / (508.877 / 4.); // 127MHz clock ticks to us, inexact rounding
107  int bunch_trg = it.GetBunchNumber(0);
108  int time_inj = it.GetTimeSinceLastInjection(0);
109  int bunch_inj = (bunch_trg - time_inj) % 1280;
110  if (bunch_inj < 0) bunch_inj += 1280;
111  rgate = bunch_inj / (1280. / 96.); // 0-96 ?
112  tinj = diff2;
113  }
114  break;
115  }
116 
117  for (auto& pix : m_storeRawHits) {
118  if (m_histos[pix.getSensorID()]) {
119  m_histos[pix.getSensorID()]->Fill(pix.getRow(), pix.getColumn(), pix.getCharge());
120  }
121  if (m_gateModeFlag) {
122  if (m_histos_gm[pix.getSensorID()]) {
123  int v = int(pix.getVCellID()) - rgate * 4;
124  if (v < 0) v += 768; // work only for inner modules
125  m_histos_gm[pix.getSensorID()]->Fill(v, pix.getColumn(), pix.getCharge());
126  }
127  if (m_histos_gm2[pix.getSensorID()]) {
128  int v = int(pix.getVCellID()) + rgate * 4;
129  if (v >= 768) v -= 768; // work only for outer modules
130  m_histos_gm2[pix.getSensorID()]->Fill(v, pix.getColumn(), pix.getCharge());
131  }
132  }
133  }
134 
135  string canvasname;
136  canvasname = std::string(Form("Run_%d_Evt_%d", evtRun, evtNr));
137 
138  gStyle->SetPalette(55);
139  gStyle->SetOptStat(0);
140 
141  {
142  m_c->Clear();
143  m_c->Divide(4, 5);
144  m_c->cd(0);
145  if (cm_flag) {
146  m_c->Pad()->SetFrameLineColor(kRed);
147  m_c->Pad()->SetFrameBorderSize(4);
148  m_c->Pad()->SetFrameBorderMode(4);
149  m_c->SetFrameLineColor(kRed);
150  m_c->SetFrameBorderSize(4);
151  m_c->SetFrameBorderMode(4);
152  }
153  int i = 1;
154 
155  for (auto h : m_histos) {
156  m_c->cd(i++);
157  if (h.second) {
158  string abc = string(h.first);
159  auto dhe = (*m_storeDAQEvtStats).findDHE(h.first);
160  if (dhe) {
161 
162  auto tg = dhe->getTriggerGate();
163  auto fn = dhe->getFrameNr();
164  auto err = dhe->getEndErrorInfo();
165 
166  abc += Form(" ERR: $%X TG: %d FN: %d IN: %d GA: %d", err, tg, fn, tinj, rgate);
167  for (auto itdhp = dhe->cbegin(); itdhp != dhe->cend(); ++itdhp) {
168  abc += Form("(%d) %d ", itdhp->getChipID(), itdhp->getFrameNr());
169  }
170  if (err != 0) {
171  m_c->Pad()->SetFillColor(kYellow);
172  }
173  if (cm_map[h.first]) {
174  m_c->Pad()->SetFillColor(kRed);
175  m_c->Pad()->SetFrameFillColor(kWhite);
176  }
177  }
178  h.second->SetTitle(abc.data());
179  h.second->Draw("colz");
180  m_l1->Draw();
181  m_l2->Draw();
182  m_l3->Draw();
183  }
184  }
185  m_c->cd(0);
186  m_c->Print((canvasname + ".png").data());
187  m_c->Print((canvasname + ".pdf").data());
188  m_c->Print((canvasname + ".root").data());
189  }
190  if (m_gateModeFlag) {
191  m_c->Clear();
192  m_c->Divide(4, 5);
193  m_c->cd(0);
194  if (cm_flag) {
195  m_c->Pad()->SetFrameLineColor(kRed);
196  m_c->Pad()->SetFrameBorderSize(4);
197  m_c->Pad()->SetFrameBorderMode(4);
198  m_c->SetFrameLineColor(kRed);
199  m_c->SetFrameBorderSize(4);
200  m_c->SetFrameBorderMode(4);
201  }
202  int i = 1;
203 
204  for (auto h : m_histos_gm) {
205  m_c->cd(i++);
206  if (h.second) {
207  string abc = string(h.first);
208  auto dhe = (*m_storeDAQEvtStats).findDHE(h.first);
209  if (dhe) {
210 
211  auto tg = dhe->getTriggerGate();
212  auto fn = dhe->getFrameNr();
213  auto err = dhe->getEndErrorInfo();
214 
215  abc += Form(" ERR: $%X TG: %d FN: %d IN: %d GA: %d", err, tg, fn, tinj, rgate);
216  for (auto itdhp = dhe->cbegin(); itdhp != dhe->cend(); ++itdhp) {
217  abc += Form("(%d) %d ", itdhp->getChipID(), itdhp->getFrameNr());
218  }
219  if (err != 0) {
220  m_c->Pad()->SetFillColor(kYellow);
221  }
222  if (cm_map[h.first]) {
223  m_c->Pad()->SetFillColor(kRed);
224  m_c->Pad()->SetFrameFillColor(kWhite);
225  }
226  }
227  h.second->SetTitle(abc.data());
228  h.second->Draw("colz");
229  m_l1->Draw();
230  m_l2->Draw();
231  m_l3->Draw();
232  }
233  }
234  m_c->cd(0);
235  m_c->Print((canvasname + "_gm.png").data());
236  m_c->Print((canvasname + "_gm.pdf").data());
237  m_c->Print((canvasname + "_gm.root").data());
238  }
239  if (m_gateModeFlag) {
240  m_c->Clear();
241  m_c->Divide(4, 5);
242  m_c->cd(0);
243  if (cm_flag) {
244  m_c->Pad()->SetFrameLineColor(kRed);
245  m_c->Pad()->SetFrameBorderSize(4);
246  m_c->Pad()->SetFrameBorderMode(4);
247  m_c->SetFrameLineColor(kRed);
248  m_c->SetFrameBorderSize(4);
249  m_c->SetFrameBorderMode(4);
250  }
251  int i = 1;
252 
253  for (auto h : m_histos_gm2) {
254  m_c->cd(i++);
255  if (h.second) {
256  string abc = string(h.first);
257  auto dhe = (*m_storeDAQEvtStats).findDHE(h.first);
258  if (dhe) {
259 
260  auto tg = dhe->getTriggerGate();
261  auto fn = dhe->getFrameNr();
262  auto err = dhe->getEndErrorInfo();
263 
264  abc += Form(" ERR: $%X TG: %d FN: %d IN: %d GA: %d", err, tg, fn, tinj, rgate);
265  for (auto itdhp = dhe->cbegin(); itdhp != dhe->cend(); ++itdhp) {
266  abc += Form("(%d) %d ", itdhp->getChipID(), itdhp->getFrameNr());
267  }
268  if (err != 0) {
269  m_c->Pad()->SetFillColor(kYellow);
270  }
271  if (cm_map[h.first]) {
272  m_c->Pad()->SetFillColor(kRed);
273  m_c->Pad()->SetFrameFillColor(kWhite);
274  }
275  }
276  h.second->SetTitle(abc.data());
277  h.second->Draw("colz");
278  m_l1->Draw();
279  m_l2->Draw();
280  m_l3->Draw();
281  }
282  }
283  m_c->cd(0);
284  m_c->Print((canvasname + "_gm2.png").data());
285  m_c->Print((canvasname + "_gm2.pdf").data());
286  m_c->Print((canvasname + "_gm2.root").data());
287  }
288 }
Belle2::PXD::PXDEventPlotModule
Plot each event with ROI and Pixels.
Definition: PXDEventPlotModule.h:47
Belle2::VxdID
Class to uniquely identify a any structure of the PXD and SVD.
Definition: VxdID.h:43
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::Module
Base class for Modules.
Definition: Module.h:74
Belle2::PXD
Namespace to encapsulate code needed for simulation and reconstrucion of the PXD.
Definition: PXDCalibrationUtilities.h:28
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19