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