Belle II Software development
PXDROIDQMModule.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/pxdDQM/PXDROIDQMModule.h>
10#include <pxd/dataobjects/PXDRawROIs.h>
11#include "TDirectory.h"
12
13using namespace std;
14using namespace Belle2;
15using namespace Belle2::PXD;
16using namespace Belle2::VXD;
17
18//-----------------------------------------------------------------
19// Register the Module
20//-----------------------------------------------------------------
21REG_MODULE(PXDROIDQM);
22
23//-----------------------------------------------------------------
24// Implementation
25//-----------------------------------------------------------------
26
28{
29 //Set module properties
30 setDescription("Monitor ROIs");
32 addParam("histogramDirectoryName", m_histogramDirectoryName, "Name of the directory where histograms will be placed",
33 std::string("pxdrawroi"));
34 addParam("PXDRawROIsName", m_PXDRawROIsName, "The name of the StoreArray of PXDRawROIs to be processed", std::string(""));
35 addParam("eachModule", m_eachModule, "create for each module", false);
36 addParam("offlineDQM", m_offlineDQM, "offline DQM, use 2D plots", false);
37}
38
40{
41 TDirectory* oldDir = gDirectory;
42 if (m_histogramDirectoryName != "") {
43 oldDir->mkdir(m_histogramDirectoryName.c_str());// do not rely on return value, might be ZERO
44 oldDir->cd(m_histogramDirectoryName.c_str());//changing to the right directory
45 }
46
47 hrawROIcount = new TH1F("hrawROIcount", "ROI count;Nr per Event", 250, 0, 250);
48 hrawROItype = new TH1F("hrawROItype", "ROI type;Nr per Event", 2, 0, 2);
49
50 hrawROIHLT_DHHID = new TH1F("hrawROIHLT_DHHID", "HLT ROI per Module;DHH ID;ROIs per Module", 64, 0, 64);
51 hrawROIDC_DHHID = new TH1F("hrawROIDC_DHHID", "DATCON ROI per Module;DHH ID;ROIs per Module", 64, 0, 64);
52
53 if (m_offlineDQM) {
54 hrawROIHLTmap = new TH2F("hrawROIHLTmap", "HLT ROI Middle Map;Ucell;Vcell", 250 / 5, 0, 250, 768 / 4, 0, 768);
55 hrawROIHLTsize = new TH2F("hrawROIHLTsize", "HLT ROI Size Map;;Ucell;Vcell", 50, 0, 200, 50, 0, 200);
56 hrawROIHLTminV = new TH1F("hrawROIHLTminV", "HLT ROI minV;V", 768, 0, 768);
57 hrawROIHLTmaxV = new TH1F("hrawROIHLTmaxV", "HLT ROI maxV;V", 768, 0, 768);
58 hrawROIHLTminU = new TH1F("hrawROIHLTminU", "HLT ROI minU;U", 250, 0, 250);
59 hrawROIHLTmaxU = new TH1F("hrawROIHLTmaxU", "HLT ROI maxU;U", 250, 0, 250);
60 hrawROIHLTsizeV = new TH1F("hrawROIHLTsizeV", "HLT ROI size;V", 768, 0, 768);
61 hrawROIHLTsizeU = new TH1F("hrawROIHLTsizeU", "HLT ROI size;U", 250, 0, 250);
62
63 hrawROIDCmap = new TH2F("hrawROIDCmap", "DATCON ROI Middle Map ;U;V", 250 / 5, 0, 250, 768 / 4, 0, 768);
64 hrawROIDCsize = new TH2F("hrawROIDCsize", "DATCON ROI Size Map ;U;V", 50, 0, 200, 50, 0, 200);
65 hrawROIDCminV = new TH1F("hrawROIDCminV", "DATCON ROI minV;V", 768, 0, 768);
66 hrawROIDCmaxV = new TH1F("hrawROIDCmaxV", "DATCON ROI maxV;V", 768, 0, 768);
67 hrawROIDCminU = new TH1F("hrawROIDCminU", "DATCON ROI minU;U", 250, 0, 250);
68 hrawROIDCmaxU = new TH1F("hrawROIDCmaxU", "DATCON ROI maxU;U", 250, 0, 250);
69 hrawROIDCsizeV = new TH1F("hrawROIDCsizeV", "DATCON ROI size;V", 768, 0, 768);
70 hrawROIDCsizeU = new TH1F("hrawROIDCsizeU", "DATCON ROI size;U", 250, 0, 250);
71
72 hrawROINrDCvsNrHLT = new TH2F("hrawROINrDCvsNrHLT", "Nr DATCON ROI vs Nr HLT ROI; Nr HLT ROI;Nr DATCON ROI",
73 100, 0, 100, 100, 0, 100);
74 }
76 std::vector<VxdID> sensors = m_vxdGeometry.getListOfSensors();
77 for (VxdID& avxdid : sensors) {
78 VXD::SensorInfoBase info = m_vxdGeometry.getSensorInfo(avxdid);
79 if (info.getType() != VXD::SensorInfoBase::PXD) continue;
80 // Only interested in PXD sensors
81
82 TString buff = (std::string)avxdid;
83 TString buffus = buff;
84 buffus.ReplaceAll(".", "_");
85
86 int dhh_id = ((avxdid.getLayerNumber() - 1) << 5) | ((avxdid.getLadderNumber()) << 1) | (avxdid.getSensorNumber() - 1);
87
88 hrawROIHLTmapModule[dhh_id] = new TH2F("hrawROIHLTmap_" + buffus,
89 "HLT ROI Middle Map " + buff + ";Ucell;Vcell", 250, 0, 250, 768, 0, 768);
90 hrawROIHLTsizeModule[dhh_id] = new TH2F("hrawROIHLTsize" + buffus,
91 "HLT ROI Size Map " + buff + ";Ucell;Vcell", 50, 0, 200, 50, 0, 200);
92 hrawROIDCmapModule[dhh_id] = new TH2F("hrawROIDCmap_" + buffus,
93 "DC ROI Middle Map " + buff + ";Ucell;Vcell", 250, 0, 250, 768, 0, 768);
94 hrawROIDCsizeModule[dhh_id] = new TH2F("hrawROIDCsize" + buffus,
95 "DC ROI Size Map " + buff + ";Ucell;Vcell", 50, 0, 200, 50, 0, 200);
96 }
97 }
98
99 // cd back to root directory
100 oldDir->cd();
101}
102
104{
105 REG_HISTOGRAM
106 m_storeROIs.isRequired(m_PXDRawROIsName);
107}
108
110{
111 // Just to make sure, reset all the histograms.
112 if (hrawROIcount) hrawROIcount->Reset();
113 if (hrawROItype) hrawROItype->Reset();
115 if (hrawROIDC_DHHID) hrawROIDC_DHHID->Reset();
116
117 for (auto& a : hrawROIHLTmapModule) if (a.second) a.second->Reset();
118 for (auto& a : hrawROIHLTsizeModule) if (a.second) a.second->Reset();
119 if (hrawROIHLTmap) hrawROIHLTmap->Reset();
120 if (hrawROIHLTsize) hrawROIHLTsize->Reset();
121 if (hrawROIHLTminV) hrawROIHLTminV->Reset();
122 if (hrawROIHLTmaxV) hrawROIHLTmaxV->Reset();
123 if (hrawROIHLTminU) hrawROIHLTminU->Reset();
124 if (hrawROIHLTmaxU) hrawROIHLTmaxU->Reset();
125 if (hrawROIHLTsizeV) hrawROIHLTsizeV->Reset();
126 if (hrawROIHLTsizeU) hrawROIHLTsizeU->Reset();
127
128 for (auto& a : hrawROIDCmapModule) if (a.second) a.second->Reset();
129 for (auto& a : hrawROIDCsizeModule) if (a.second) a.second->Reset();
130 if (hrawROIDCmap) hrawROIDCmap->Reset();
131 if (hrawROIDCsize) hrawROIDCsize->Reset();
132 if (hrawROIDCminV) hrawROIDCminV->Reset();
133 if (hrawROIDCmaxV) hrawROIDCmaxV->Reset();
134 if (hrawROIDCminU) hrawROIDCminU->Reset();
135 if (hrawROIDCmaxU) hrawROIDCmaxU->Reset();
136 if (hrawROIDCsizeV) hrawROIDCsizeV->Reset();
137 if (hrawROIDCsizeU) hrawROIDCsizeU->Reset();
138
140}
141
143{
144 int nr_HLT = 0;
145 int nr_DC = 0;
146 if (hrawROIcount) hrawROIcount->Fill(-1); // misuse underflow for event count
147 for (auto& it : m_storeROIs) {
148 int nr;
149 nr = it.getNrROIs();
150 if (hrawROIcount) hrawROIcount->Fill(nr);
151 if (hrawROIDC_DHHID) hrawROIDC_DHHID->Fill(-1); // misuse underflow for roi raw packet"event" count
152 if (hrawROIHLT_DHHID) hrawROIHLT_DHHID->Fill(-1); // misuse underflow for roi raw packet"event" count
153 for (auto j = 0; j < nr; j++) {
154 if (hrawROItype) hrawROItype->Fill(it.getType(j));
155 int Vmin, Vmax, Umin, Umax, Vmean, Umean, Vsize, Usize;
156 Vmin = it.getMinVid(j);
157 Vmax = it.getMaxVid(j);
158 Vsize = Vmax - Vmin;
159 Vmean = (Vmin + Vmax) / 2;
160 Umin = it.getMinUid(j);
161 Umax = it.getMaxUid(j);
162 Usize = Umax - Umin;
163 Umean = (Umin + Umax) / 2;
164 auto id = it.getDHHID(j);
165 if (it.getType(j)) {
166 nr_DC++;
167 if (hrawROIDC_DHHID) hrawROIDC_DHHID->Fill(id);
168 if (hrawROIDCmapModule[id]) hrawROIDCmapModule[id]->Fill(Umean, Vmean);
169 if (hrawROIDCsizeModule[id]) hrawROIDCsizeModule[id]->Fill(Usize, Vsize);
170 if (hrawROIDCmap) hrawROIDCmap->Fill(Umean, Vmean);
171 if (hrawROIDCsize) hrawROIDCsize->Fill(Usize, Vsize);
172 if (hrawROIDCminV) hrawROIDCminV->Fill(Vmin);
173 if (hrawROIDCmaxV) hrawROIDCmaxV->Fill(Vmax);
174 if (hrawROIDCminU) hrawROIDCminU->Fill(Umin);
175 if (hrawROIDCmaxU) hrawROIDCmaxU->Fill(Umax);
176 if (hrawROIDCsizeV) hrawROIDCsizeV->Fill(Vsize);
177 if (hrawROIDCsizeU) hrawROIDCsizeU->Fill(Usize);
178 } else {
179 nr_HLT++;
180 if (hrawROIHLT_DHHID) hrawROIHLT_DHHID->Fill(id);
181 if (hrawROIHLTmapModule[id]) hrawROIHLTmapModule[id]->Fill(Umean, Vmean);
182 if (hrawROIHLTsizeModule[id]) hrawROIHLTsizeModule[id]->Fill(Usize, Vsize);
183 if (hrawROIHLTmap) hrawROIHLTmap->Fill(Umean, Vmean);
184 if (hrawROIHLTsize) hrawROIHLTsize->Fill(Usize, Vsize);
185 if (hrawROIHLTminV) hrawROIHLTminV->Fill(Vmin);
186 if (hrawROIHLTmaxV) hrawROIHLTmaxV->Fill(Vmax);
187 if (hrawROIHLTminU) hrawROIHLTminU->Fill(Umin);
188 if (hrawROIHLTmaxU) hrawROIHLTmaxU->Fill(Umax);
189 if (hrawROIHLTsizeV) hrawROIHLTsizeV->Fill(Vsize);
190 if (hrawROIHLTsizeU) hrawROIHLTsizeU->Fill(Usize);
191 }
192 }
193 }
194 if (hrawROINrDCvsNrHLT) hrawROINrDCvsNrHLT->Fill(nr_HLT, nr_DC);
195}
HistoModule()
Constructor.
Definition HistoModule.h:32
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
@ 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
void initialize() override final
Initialize.
StoreArray< PXDRawROIs > m_storeROIs
Storearray for ROIs.
void defineHisto() override final
Define histograms.
TH2F * hrawROINrDCvsNrHLT
Histogram.
bool m_eachModule
create a histo per module
bool m_offlineDQM
create 2d plots for offline checks
void event() override final
Event.
std::string m_histogramDirectoryName
Name of the histogram directory in ROOT file.
std::map< int, TH2F * > hrawROIHLTmapModule
Histogram.
void beginRun() override final
Begin run.
VXD::GeoCache & m_vxdGeometry
the VXD geometry
TH1F * hrawROIHLT_DHHID
Histogram.
std::map< int, TH2F * > hrawROIHLTsizeModule
Histogram.
TH1F * hrawROIcount
Histogram 2d hitmap.
PXDROIDQMModule()
Constructor defining the parameters.
std::map< int, TH2F * > hrawROIDCmapModule
Histogram.
std::map< int, TH2F * > hrawROIDCsizeModule
Histogram.
std::string m_PXDRawROIsName
RawROI StoreArray name.
Class to facilitate easy access to sensor information of the VXD like coordinate transformations or p...
Definition GeoCache.h:38
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.