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