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