Belle II Software development
PXDDQMClustersModule.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/PXDDQMClustersModule.h"
10
11#include <framework/gearbox/Unit.h>
12#include <framework/datastore/StoreArray.h>
13#include <framework/datastore/RelationArray.h>
14
15#include <pxd/dataobjects/PXDDAQStatus.h>
16#include <pxd/dataobjects/PXDDigit.h>
17#include <pxd/dataobjects/PXDCluster.h>
18#include <pxd/geometry/SensorInfo.h>
19
20#include <vxd/geometry/GeoCache.h>
21#include <vxd/geometry/GeoTools.h>
22#include <pxd/unpacking/PXDMappingLookup.h>
23#include <pxd/reconstruction/PXDGainCalibrator.h>
24
25#include <boost/format.hpp>
26
27#include "TDirectory.h"
28
29using namespace std;
30using boost::format;
31using namespace Belle2;
32using namespace Belle2::PXD;
33
34//-----------------------------------------------------------------
35// Register the Module
36//-----------------------------------------------------------------
37REG_MODULE(PXDDQMClusters);
38
39
40//-----------------------------------------------------------------
41// Implementation
42//-----------------------------------------------------------------
43
45{
46 //Set module properties
47 setDescription("PXD DQM clusters module "
48 "Recommended Number of events for monitoring is 40 kEvents or more to fill all histograms "
49 );
50
51 setPropertyFlags(c_ParallelProcessingCertified); // specify this flag if you need parallel processing
52 addParam("CutMinCharge", m_CutMinCharge,
53 "cut on pixel charge for accepting good pixel, default >= 12", 12);
54 addParam("CutMinSeedCharge", m_CutMinSeedCharge,
55 "cut on cluster seed for accepting good cluster, default >= 12", 12);
56 addParam("CutMaxClusterSize", m_CutMaxClusterSize,
57 "cut on cluster size accepting good cluster, default <= 4", 4);
58 addParam("CutMinClusterCharge", m_CutMinClusterCharge,
59 "cut on cluster charge for accepting good cluster, default >= 12", 12);
60 addParam("histogramDirectoryName", m_histogramDirectoryName, "Name of the directory where histograms will be placed",
61 std::string("PXDDQMClusters"));
62}
63
64
65//------------------------------------------------------------------
66// Function to define histograms
67//-----------------------------------------------------------------
68
70{
72 if (gTools->getNumberOfLayers() == 0) {
73 B2FATAL("Missing geometry for VXD, check steering file.");
74 }
75 if (gTools->getNumberOfPXDLayers() == 0) {
76 B2WARNING("Missing geometry for PXD, PXD-DQM is skipped.");
77 return;
78 }
79
80 // Create a separate histogram directories and cd into it.
81 TDirectory* oldDir = gDirectory;
82 if (m_histogramDirectoryName != "") {
83 oldDir->mkdir(m_histogramDirectoryName.c_str());// do not use return value with ->cd(), its ZERO if dir already exists
84 oldDir->cd(m_histogramDirectoryName.c_str());
85 }
86
87 // basic constants presets:
88 int nPXDSensors = gTools->getNumberOfPXDSensors();
89 int nPXDChips = gTools->getTotalPXDChips();
90
91 // Create basic histograms:
92 m_hitMapCounts = new TH1D("DQM_PXD_PixelHitmapCounts", "PXD Integrated number of fired pixels per sensor",
93 nPXDSensors, 0, nPXDSensors);
94 m_hitMapCounts->GetXaxis()->SetTitle("Sensor ID");
95 m_hitMapCounts->GetYaxis()->SetTitle("counts");
96
97 m_hitMapFilterCounts = new TH1D("DQM_PXD_PixelHitmapFilterCounts", "PXD Integrated number of filtered pixels per sensor",
98 nPXDSensors, 0, nPXDSensors);
99 m_hitMapFilterCounts->GetXaxis()->SetTitle("Sensor ID");
100 m_hitMapFilterCounts->GetYaxis()->SetTitle("counts");
101
102 m_hitMapClCounts = new TH1D("DQM_PXD_ClusterHitmapCounts", "PXD Integrated number of clusters per sensor",
103 nPXDSensors, 0, nPXDSensors);
104 m_hitMapClCounts->GetXaxis()->SetTitle("Sensor ID");
105 m_hitMapClCounts->GetYaxis()->SetTitle("counts");
106
107 m_hitMapClFilterCounts = new TH1D("DQM_PXD_ClusterHitmapFilterCounts", "PXD Integrated number of filtered clusters per sensor",
108 nPXDSensors, 0, nPXDSensors);
109 m_hitMapClFilterCounts->GetXaxis()->SetTitle("Sensor ID");
110 m_hitMapClFilterCounts->GetYaxis()->SetTitle("counts");
111
112 // basic counters per chip:
113 m_hitMapCountsChip = new TH1D("DQM_PXD_PixelHitmapCountsChip", "PXD Integrated number of fired pixels per chip",
114 nPXDChips, 0, nPXDChips);
115 m_hitMapCountsChip->GetXaxis()->SetTitle("Chip ID");
116 m_hitMapCountsChip->GetYaxis()->SetTitle("counts");
117 m_hitMapClCountsChip = new TH1D("DQM_PXD_ClusterHitmapCountsChip", "PXD Integrated number of clusters per chip",
118 nPXDChips, 0, nPXDChips);
119 m_hitMapClCountsChip->GetXaxis()->SetTitle("Chip ID");
120 m_hitMapClCountsChip->GetYaxis()->SetTitle("counts");
121 for (int i = 0; i < nPXDChips; i++) {
122 VxdID id = gTools->getChipIDFromPXDIndex(i);
123 int iLayer = id.getLayerNumber();
124 int iLadder = id.getLadderNumber();
125 int iSensor = id.getSensorNumber();
126 int iChip = gTools->getPXDChipNumber(id);
127 int IsU = gTools->isPXDSideU(id);
128 TString AxisTicks = Form("%i_%i_%i_u%iDCD", iLayer, iLadder, iSensor, iChip);
129 if (!IsU)
130 AxisTicks = Form("%i_%i_%i_v%iSWB", iLayer, iLadder, iSensor, iChip);
131 m_hitMapCountsChip->GetXaxis()->SetBinLabel(i + 1, AxisTicks.Data());
132 m_hitMapClCountsChip->GetXaxis()->SetBinLabel(i + 1, AxisTicks.Data());
133 }
134
135 for (int i = 0; i < nPXDSensors; i++) {
136 VxdID id = gTools->getSensorIDFromPXDIndex(i);
137 int iLayer = id.getLayerNumber();
138 int iLadder = id.getLadderNumber();
139 int iSensor = id.getSensorNumber();
140 TString AxisTicks = Form("%i_%i_%i", iLayer, iLadder, iSensor);
141 m_hitMapCounts->GetXaxis()->SetBinLabel(i + 1, AxisTicks.Data());
142 m_hitMapClCounts->GetXaxis()->SetBinLabel(i + 1, AxisTicks.Data());
143 }
144
145 m_fired.resize(nPXDSensors);
146 m_goodfired.resize(nPXDSensors);
147 m_clusters.resize(nPXDSensors);
148 m_goodclusters.resize(nPXDSensors);
149 m_startRow.resize(nPXDSensors);
150 m_chargStartRow.resize(nPXDSensors);
151 m_startRowCount.resize(nPXDSensors);
152 m_clusterCharge.resize(nPXDSensors);
153 m_clusterEnergy.resize(nPXDSensors);
154 m_pixelSignal.resize(nPXDSensors);
155 m_clusterSizeU.resize(nPXDSensors);
156 m_clusterSizeV.resize(nPXDSensors);
157 m_clusterSizeUV.resize(nPXDSensors);
158
159 m_hitMapU.resize(nPXDSensors);
160 m_hitMapV.resize(nPXDSensors);
161 m_hitMap.resize(nPXDSensors);
162 m_hitMapUCl.resize(nPXDSensors);
163 m_hitMapVCl.resize(nPXDSensors);
164 m_hitMapCl.resize(nPXDSensors);
165 m_seed.resize(nPXDSensors);
166 for (int i = 0; i < nPXDSensors; i++) {
167 VxdID id = gTools->getSensorIDFromPXDIndex(i);
168 int iLayer = id.getLayerNumber();
169 int iLadder = id.getLadderNumber();
170 int iSensor = id.getSensorNumber();
171 VxdID sensorID(iLayer, iLadder, iSensor);
173 string sensorDescr = str(format("%1%_%2%_%3%") % iLayer % iLadder % iSensor);
174 auto nUPixels = SensorInfo.getUCells();
175 auto nVPixels = SensorInfo.getVCells();
176 //----------------------------------------------------------------
177 // Number of fired pixels per frame
178 //----------------------------------------------------------------
179 string name = str(format("DQM_PXD_%1%_Fired") % sensorDescr);
180 string title = str(format("PXD Sensor %1% Fired pixels") % sensorDescr);
181 m_fired[i] = new TH1D(name.c_str(), title.c_str(), 50, 0, 50);
182 m_fired[i]->GetXaxis()->SetTitle("# of fired pixels");
183 m_fired[i]->GetYaxis()->SetTitle("counts");
184 //----------------------------------------------------------------
185 // Number of good fired pixels per frame
186 //----------------------------------------------------------------
187 name = str(format("DQM_PXD_%1%_GoodFired") % sensorDescr);
188 title = str(format("PXD Sensor %1% Good Fired pixels") % sensorDescr);
189 m_goodfired[i] = new TH1D(name.c_str(), title.c_str(), 50, 0, 50);
190 m_goodfired[i]->GetXaxis()->SetTitle("# of fired pixels");
191 m_goodfired[i]->GetYaxis()->SetTitle("counts");
192 //----------------------------------------------------------------
193 // Number of clusters per frame
194 //----------------------------------------------------------------
195 name = str(format("DQM_PXD_%1%_Clusters") % sensorDescr);
196 title = str(format("PXD Sensor %1% Number of clusters") % sensorDescr);
197 m_clusters[i] = new TH1D(name.c_str(), title.c_str(), 20, 0, 20);
198 m_clusters[i]->GetXaxis()->SetTitle("# of clusters");
199 m_clusters[i]->GetYaxis()->SetTitle("counts");
200 //----------------------------------------------------------------
201 // Number of good clusters per frame
202 //----------------------------------------------------------------
203 name = str(format("DQM_PXD_%1%_GoodClusters") % sensorDescr);
204 title = str(format("PXD Sensor %1% Number of good clusters") % sensorDescr);
205 m_goodclusters[i] = new TH1D(name.c_str(), title.c_str(), 20, 0, 20);
206 m_goodclusters[i]->GetXaxis()->SetTitle("# of clusters");
207 m_goodclusters[i]->GetYaxis()->SetTitle("counts");
208 //----------------------------------------------------------------
209 // Start row distribution
210 //----------------------------------------------------------------
211 name = str(format("DQM_PXD_%1%_StartRow") % sensorDescr);
212 title = str(format("PXD Sensor %1% Start row distribution") % sensorDescr);
213
214 m_startRow[i] = new TH1D(name.c_str(), title.c_str(), nVPixels / 4, 0.0, nVPixels);
215 m_startRow[i]->GetXaxis()->SetTitle("start row [pitch units]");
216 m_startRow[i]->GetYaxis()->SetTitle("count");
217 //----------------------------------------------------------------
218 // Cluster seed charge by distance from the start row
219 //----------------------------------------------------------------
220 name = str(format("DQM_PXD_%1%_AverageSeedByStartRow") % sensorDescr);
221 title = str(format("PXD Sensor %1% Average seed charge by distance from the start row") % sensorDescr);
222 m_chargStartRow[i] = new TH1D(name.c_str(), title.c_str(), nVPixels / 4, 0.0, nVPixels);
223 m_chargStartRow[i]->GetXaxis()->SetTitle("distance from the start row [pitch units]");
224 m_chargStartRow[i]->GetYaxis()->SetTitle("average seed [ADU]");
225 name = str(format("DQM_PXD_%1%_SeedCountsByStartRow") % sensorDescr);
226 title = str(format("PXD Sensor %1% Seed charge count by distance from the start row") % sensorDescr);
227 m_startRowCount[i] = new TH1D(name.c_str(), title.c_str(), nVPixels / 4, 0.0, nVPixels);
228 m_startRowCount[i]->GetXaxis()->SetTitle("distance from the start row [pitch units]");
229 m_startRowCount[i]->GetYaxis()->SetTitle("count");
230 //----------------------------------------------------------------
231 // Cluster Charge
232 //----------------------------------------------------------------
233 name = str(format("DQM_PXD_%1%_ClusterCharge") % sensorDescr);
234 title = str(format("PXD Sensor %1% Cluster Charge") % sensorDescr);
235 m_clusterCharge[i] = new TH1D(name.c_str(), title.c_str(), 256, 0, 256);
236 m_clusterCharge[i]->GetXaxis()->SetTitle("charge of clusters [ADU]");
237 m_clusterCharge[i]->GetYaxis()->SetTitle("counts");
238 //----------------------------------------------------------------
239 // Cluster Energy
240 //----------------------------------------------------------------
241 name = str(format("DQM_PXD_%1%_ClusterEnergy") % sensorDescr);
242 title = str(format("PXD Sensor %1% Cluster Energy") % sensorDescr);
243 m_clusterEnergy[i] = new TH1D(name.c_str(), title.c_str(), 100, 0, 50);
244 m_clusterEnergy[i]->GetXaxis()->SetTitle("energy of clusters [keV]");
245 m_clusterEnergy[i]->GetYaxis()->SetTitle("counts");
246 //----------------------------------------------------------------
247 // Pixel Signal
248 //----------------------------------------------------------------
249 name = str(format("DQM_PXD_%1%_PixelSignal") % sensorDescr);
250 title = str(format("PXD Sensor %1% Pixel Signal") % sensorDescr);
251 m_pixelSignal[i] = new TH1D(name.c_str(), title.c_str(), 256, 0, 256);
252 m_pixelSignal[i]->GetXaxis()->SetTitle("signal of pixels [ADU]");
253 m_pixelSignal[i]->GetYaxis()->SetTitle("counts");
254 //----------------------------------------------------------------
255 // Cluster Size in U
256 //----------------------------------------------------------------
257 name = str(format("DQM_PXD_%1%_ClusterSizeU") % sensorDescr);
258 title = str(format("PXD Sensor %1% Cluster Size U") % sensorDescr);
259 m_clusterSizeU[i] = new TH1D(name.c_str(), title.c_str(), 10, 0, 10);
260 m_clusterSizeU[i]->GetXaxis()->SetTitle("size of u clusters");
261 m_clusterSizeU[i]->GetYaxis()->SetTitle("counts");
262 //----------------------------------------------------------------
263 // Cluster Size in V
264 //----------------------------------------------------------------
265 name = str(format("DQM_PXD_%1%_ClusterSizeV") % sensorDescr);
266 title = str(format("PXD Sensor %1% Cluster Size V") % sensorDescr);
267 m_clusterSizeV[i] = new TH1D(name.c_str(), title.c_str(), 10, 0, 10);
268 m_clusterSizeV[i]->GetXaxis()->SetTitle("size of v clusters");
269 m_clusterSizeV[i]->GetYaxis()->SetTitle("counts");
270 //----------------------------------------------------------------
271 // Cluster Size in U+V
272 //----------------------------------------------------------------
273 name = str(format("DQM_PXD_%1%_ClusterSizeUV") % sensorDescr);
274 title = str(format("PXD Sensor %1% Cluster Size U+V") % sensorDescr);
275 m_clusterSizeUV[i] = new TH1D(name.c_str(), title.c_str(), 10, 0, 10);
276 m_clusterSizeUV[i]->GetXaxis()->SetTitle("size of u+v clusters");
277 m_clusterSizeUV[i]->GetYaxis()->SetTitle("counts");
278
279 //----------------------------------------------------------------
280 // Hitmaps: Number of pixels by coordinate
281 //----------------------------------------------------------------
282 // Hitmaps in U
283 name = str(format("PXD_%1%_PixelHitmapU") % sensorDescr);
284 title = str(format("PXD Sensor %1% Pixel Hitmap in U") % sensorDescr);
285 m_hitMapU[i] = new TH1D(name.c_str(), title.c_str(), nUPixels, 0, nUPixels);
286 m_hitMapU[i]->GetXaxis()->SetTitle("u position [pitch units]");
287 m_hitMapU[i]->GetYaxis()->SetTitle("hits");
288 // Hitmaps in V
289 name = str(format("PXD_%1%_PixelHitmapV") % sensorDescr);
290 title = str(format("PXD Sensor %1% Pixel Hitmap in V") % sensorDescr);
291 m_hitMapV[i] = new TH1D(name.c_str(), title.c_str(), nVPixels, 0, nVPixels);
292 m_hitMapV[i]->GetXaxis()->SetTitle("v position [pitch units]");
293 m_hitMapV[i]->GetYaxis()->SetTitle("hits");
294 // Hitmaps in UV
295 name = str(format("PXD_%1%_PixelHitmap") % sensorDescr);
296 title = str(format("PXD Sensor %1% Pixel Hitmap") % sensorDescr);
297 m_hitMap[i] = new TH2D(name.c_str(), title.c_str(), nUPixels, 0, nUPixels, nVPixels, 0, nVPixels);
298 m_hitMap[i]->GetXaxis()->SetTitle("u position [pitch units]");
299 m_hitMap[i]->GetYaxis()->SetTitle("v position [pitch units]");
300 m_hitMap[i]->GetZaxis()->SetTitle("hits");
301
302 //----------------------------------------------------------------
303 // Hitmaps: Number of clusters by coordinate
304 //----------------------------------------------------------------
305 // Hitmaps in U
306 name = str(format("PXD_%1%_HitmapClstU") % sensorDescr);
307 title = str(format("PXD Sensor %1% Hitmap Clusters in U") % sensorDescr);
308 m_hitMapUCl[i] = new TH1D(name.c_str(), title.c_str(), nUPixels, 0, nUPixels);
309 m_hitMapUCl[i]->GetXaxis()->SetTitle("u position [pitch units]");
310 m_hitMapUCl[i]->GetYaxis()->SetTitle("hits");
311 // Hitmaps in V
312 name = str(format("PXD_%1%_HitmapClstV") % sensorDescr);
313 title = str(format("PXD Sensor %1% Hitmap Clusters in V") % sensorDescr);
314 m_hitMapVCl[i] = new TH1D(name.c_str(), title.c_str(), nVPixels, 0, nVPixels);
315 m_hitMapVCl[i]->GetXaxis()->SetTitle("v position [pitch units]");
316 m_hitMapVCl[i]->GetYaxis()->SetTitle("hits");
317 // Hitmaps in UV
318 name = str(format("PXD_%1%_HitmapClst") % sensorDescr);
319 title = str(format("PXD Sensor %1% Hitmap Clusters") % sensorDescr);
320 m_hitMapCl[i] = new TH2D(name.c_str(), title.c_str(), nUPixels, 0, nUPixels, nVPixels, 0, nVPixels);
321 m_hitMapCl[i]->GetXaxis()->SetTitle("u position [pitch units]");
322 m_hitMapCl[i]->GetYaxis()->SetTitle("v position [pitch units]");
323 m_hitMapCl[i]->GetZaxis()->SetTitle("hits");
324
325 //----------------------------------------------------------------
326 // Cluster seed charge distribution
327 //----------------------------------------------------------------
328 name = str(format("PXD_%1%_Seed") % sensorDescr);
329 title = str(format("PXD Sensor %1% Seed charge") % sensorDescr);
330 m_seed[i] = new TH1D(name.c_str(), title.c_str(), 256, 0, 256);
331 m_seed[i]->GetXaxis()->SetTitle("seed charge of clusters [ADU]");
332 m_seed[i]->GetYaxis()->SetTitle("count");
333
334 }
335 oldDir->cd();
336}
337
338
340{
341 // Register histograms (calls back defineHisto)
342 REG_HISTOGRAM
343
344 m_storeDAQEvtStats.isOptional();
345
346 auto gTools = VXD::GeoCache::getInstance().getGeoTools();
347 if (gTools->getNumberOfPXDLayers() != 0) {
348 //Register collections
351 RelationArray relPXDClusterDigits(storePXDClusters, storePXDDigits);
352 m_storePXDClustersName = storePXDClusters.getName();
353 m_relPXDClusterDigitName = relPXDClusterDigits.getName();
354
355 //Store names to speed up creation later
356 m_storePXDDigitsName = storePXDDigits.getName();
357 }
358}
359
361{
362 auto gTools = VXD::GeoCache::getInstance().getGeoTools();
363 if (gTools->getNumberOfPXDLayers() == 0) return;
364
365 if (m_hitMapCounts != nullptr) m_hitMapCounts->Reset();
366 if (m_hitMapFilterCounts != nullptr) m_hitMapFilterCounts->Reset();
367 if (m_hitMapClCounts != nullptr) m_hitMapClCounts->Reset();
368 if (m_hitMapClFilterCounts != nullptr) m_hitMapClFilterCounts->Reset();
369 if (m_hitMapCountsChip != nullptr) m_hitMapCountsChip->Reset();
370 if (m_hitMapClCountsChip != nullptr) m_hitMapClCountsChip->Reset();
371
372 for (int i = 0; i < gTools->getNumberOfPXDSensors(); i++) {
373 if (m_fired[i] != nullptr) m_fired[i]->Reset();
374 if (m_goodfired[i] != nullptr) m_goodfired[i]->Reset();
375 if (m_clusters[i] != nullptr) m_clusters[i]->Reset();
376 if (m_goodclusters[i] != nullptr) m_goodclusters[i]->Reset();
377 if (m_startRow[i] != nullptr) m_startRow[i]->Reset();
378 if (m_chargStartRow[i] != nullptr) m_chargStartRow[i]->Reset();
379 if (m_startRowCount[i] != nullptr) m_startRowCount[i]->Reset();
380 if (m_clusterCharge[i] != nullptr) m_clusterCharge[i]->Reset();
381 if (m_clusterEnergy[i] != nullptr) m_clusterEnergy[i]->Reset();
382 if (m_pixelSignal[i] != nullptr) m_pixelSignal[i]->Reset();
383 if (m_clusterSizeU[i] != nullptr) m_clusterSizeU[i]->Reset();
384 if (m_clusterSizeV[i] != nullptr) m_clusterSizeV[i]->Reset();
385 if (m_clusterSizeUV[i] != nullptr) m_clusterSizeUV[i]->Reset();
386
387 if (m_hitMapU[i] != nullptr) m_hitMapU[i]->Reset();
388 if (m_hitMapV[i] != nullptr) m_hitMapV[i]->Reset();
389 if (m_hitMap[i] != nullptr) m_hitMap[i]->Reset();
390 if (m_hitMapUCl[i] != nullptr) m_hitMapUCl[i]->Reset();
391 if (m_hitMapVCl[i] != nullptr) m_hitMapVCl[i]->Reset();
392 if (m_hitMapCl[i] != nullptr) m_hitMapCl[i]->Reset();
393 if (m_seed[i] != nullptr) m_seed[i]->Reset();
394
395 }
396}
397
398
400{
401 auto gTools = VXD::GeoCache::getInstance().getGeoTools();
402 if (gTools->getNumberOfPXDLayers() == 0) return;
403
404 const StoreArray<PXDDigit> storePXDDigits(m_storePXDDigitsName);
405 const StoreArray<PXDCluster> storePXDClusters(m_storePXDClustersName);
406 const RelationArray relPXDClusterDigits(storePXDClusters, storePXDDigits, m_relPXDClusterDigitName);
407
408 // If there are no digits, leave
409 if (!storePXDDigits || !storePXDDigits.getEntries()) return;
410
411 int firstPXDLayer = gTools->getFirstPXDLayer();
412 int lastPXDLayer = gTools->getLastPXDLayer();
413 int nPXDSensors = gTools->getNumberOfPXDSensors();
414
415 // PXD basic histograms:
416 // Fired pixels
417 vector< int > Pixels(nPXDSensors);
418 vector< int > GoodPixels(nPXDSensors);
419 for (const PXDDigit& digit : storePXDDigits) {
420 int iLayer = digit.getSensorID().getLayerNumber();
421 if ((iLayer < firstPXDLayer) || (iLayer > lastPXDLayer)) continue;
422 int iLadder = digit.getSensorID().getLadderNumber();
423 int iSensor = digit.getSensorID().getSensorNumber();
424 VxdID sensorID(iLayer, iLadder, iSensor);
425 int index = gTools->getPXDSensorIndex(sensorID);
427 Pixels[index]++;
428 int iChip = PXDMappingLookup::getDCDID(digit.getUCellID(), digit.getVCellID(), sensorID);
429 int indexChip = gTools->getPXDChipIndex(sensorID, kTRUE, iChip);
430 if (m_hitMapCountsChip != nullptr) m_hitMapCountsChip->Fill(indexChip);
431 iChip = PXDMappingLookup::getSWBID(digit.getVCellID());
432 indexChip = gTools->getPXDChipIndex(sensorID, kFALSE, iChip);
433 if (m_hitMapCountsChip != nullptr) m_hitMapCountsChip->Fill(indexChip);
434 if (m_pixelSignal[index] != nullptr) m_pixelSignal[index]->Fill(digit.getCharge());
435 if (m_hitMapCounts != nullptr) m_hitMapCounts->Fill(index);
436 // filter for pixel charge
437 if (digit.getCharge() >= m_CutMinCharge && digit.getCharge() < 255) {
438 GoodPixels[index]++;
439 if (m_hitMapFilterCounts != nullptr) m_hitMapFilterCounts->Fill(index);
440 if (m_hitMapU[index] != nullptr) m_hitMapU[index]->Fill(digit.getUCellID());
441 if (m_hitMapV[index] != nullptr) m_hitMapV[index]->Fill(digit.getVCellID());
442 if (m_hitMap[index] != nullptr) m_hitMap[index]->Fill(digit.getUCellID(), digit.getVCellID());
443 }
444 }
445 for (int i = 0; i < nPXDSensors; i++) {
446 if (m_fired[i] != nullptr && Pixels[i] > 0) m_fired[i]->Fill(Pixels[i]);
447 if (m_goodfired[i] != nullptr && GoodPixels[i] > 0) m_goodfired[i]->Fill(GoodPixels[i]);
448 }
449
450 // Hitmaps, Charge, Seed, Size, ...
451 vector< int > Clusters(nPXDSensors);
452 vector< int > GoodClusters(nPXDSensors);
453 for (const PXDCluster& cluster : storePXDClusters) {
454 int iLayer = cluster.getSensorID().getLayerNumber();
455 if ((iLayer < firstPXDLayer) || (iLayer > lastPXDLayer)) continue;
456 int iLadder = cluster.getSensorID().getLadderNumber();
457 int iSensor = cluster.getSensorID().getSensorNumber();
458 VxdID sensorID(iLayer, iLadder, iSensor);
459 int index = gTools->getPXDSensorIndex(sensorID);
461 Clusters[index]++;
462 int iChip = PXDMappingLookup::getDCDID(SensorInfo.getUCellID(cluster.getU()), SensorInfo.getVCellID(cluster.getV()), sensorID);
463 int indexChip = gTools->getPXDChipIndex(sensorID, kTRUE, iChip);
464 if (m_hitMapClCountsChip != nullptr) m_hitMapClCountsChip->Fill(indexChip);
465 iChip = PXDMappingLookup::getSWBID(SensorInfo.getVCellID(cluster.getV()));
466 indexChip = gTools->getPXDChipIndex(sensorID, kFALSE, iChip);
467 if (m_hitMapClCountsChip != nullptr) m_hitMapClCountsChip->Fill(indexChip);
468 if (m_hitMapClCounts != nullptr) m_hitMapClCounts->Fill(index);
469 if (m_clusterCharge[index] != nullptr) m_clusterCharge[index]->Fill(cluster.getCharge());
470 // FIXME: The cluster charge is stored in ADU. We have to extract the
471 // area dependent conversion factor ADU->eV. Assume this is the same for all pixels of the
472 // cluster.
473 auto cluster_uID = SensorInfo.getUCellID(cluster.getU());
474 auto cluster_vID = SensorInfo.getVCellID(cluster.getV());
475 auto ADUToEnergy = PXDGainCalibrator::getInstance().getADUToEnergy(sensorID, cluster_uID, cluster_vID);
476 if (m_clusterEnergy[index] != nullptr) m_clusterEnergy[index]->Fill(cluster.getCharge()* ADUToEnergy / Unit::keV);
477 if (m_clusterSizeU[index] != nullptr) m_clusterSizeU[index]->Fill(cluster.getUSize());
478 if (m_clusterSizeV[index] != nullptr) m_clusterSizeV[index]->Fill(cluster.getVSize());
479 if (m_clusterSizeUV[index] != nullptr) m_clusterSizeUV[index]->Fill(cluster.getSize());
480 if (m_seed[index] != nullptr) m_seed[index]->Fill(cluster.getSeedCharge());
481
482 // filter for cluster size and seed pixel
483 if (cluster.getSize() <= m_CutMaxClusterSize && cluster.getCharge() >= m_CutMinClusterCharge
484 && cluster.getSeedCharge() >= m_CutMinSeedCharge && cluster.getSeedCharge() < 255) {
485 GoodClusters[index]++;
486 if (m_hitMapClFilterCounts != nullptr) m_hitMapClFilterCounts->Fill(index);
487
488 if (m_hitMapUCl[index] != nullptr) m_hitMapUCl[index]->Fill(
489 SensorInfo.getUCellID(cluster.getU()));
490 if (m_hitMapVCl[index] != nullptr) m_hitMapVCl[index]->Fill(
491 SensorInfo.getVCellID(cluster.getV()));
492 if (m_hitMapCl[index] != nullptr) m_hitMapCl[index]->Fill(
493 SensorInfo.getUCellID(cluster.getU()),
494 SensorInfo.getVCellID(cluster.getV()));
495 }
496 }
497 for (int i = 0; i < nPXDSensors; i++) {
498 if (m_clusters[i] != nullptr && Clusters[i] > 0) m_clusters[i]->Fill(Clusters[i]);
499 if (m_goodclusters[i] != nullptr && GoodClusters[i] > 0) m_goodclusters[i]->Fill(GoodClusters[i]);
500 }
501
502 // Only fill start row (triggergate) histos when data is available
503 if (m_storeDAQEvtStats.isValid()) {
504 std::map<VxdID, unsigned short> startRows;
505 for (int index = 0; index < nPXDSensors; index++) {
506 VxdID id = gTools->getSensorIDFromPXDIndex(index);
507
508 const PXDDAQDHEStatus* dhe = (*m_storeDAQEvtStats).findDHE(id);
509 if (dhe != nullptr) {
510 auto startRow = dhe->getStartRow();
511 if (m_startRow[index] != nullptr) m_startRow[index]->Fill(startRow);
512 startRows.insert(std::make_pair(id, startRow));
513 } else {
514 B2WARNING("No PXDDAQDHEStatus for VXD Sensor " << id << " found.");
515 }
516 }
517
518 // Cluster seed charge by start row
519 for (auto& cluster : storePXDClusters) {
520 VxdID sensorID = cluster.getSensorID();
521 int index = gTools->getPXDSensorIndex(sensorID);
523
524 float fDistance = SensorInfo.getVCellID(cluster.getV()) - startRows[cluster.getSensorID()];
525 if (fDistance < 0) fDistance += SensorInfo.getVCells();
526 if (m_chargStartRow[index] != nullptr) m_chargStartRow[index]->Fill(fDistance, cluster.getSeedCharge());
527 if (m_startRowCount[index] != nullptr) m_startRowCount[index]->Fill(fDistance);
528 }
529 }
530
531}
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
The PXD Cluster class This class stores all information about reconstructed PXD clusters The position...
Definition PXDCluster.h:30
The PXD DAQ DHE Status class.
unsigned short getStartRow(void) const
get Trigger Start Row
std::vector< TH1D * > m_goodfired
Filtered fired pixels per event.
TH1D * m_hitMapClFilterCounts
Hitmaps of filtered Clusters.
void initialize() override final
Initialize.
std::string m_relPXDClusterDigitName
PXDClustersToPXDDigits RelationArray name.
std::vector< TH1D * > m_hitMapUCl
Hitmaps clusters for u.
std::vector< TH2D * > m_hitMapCl
Hitmaps clusters.
std::vector< TH1D * > m_fired
Fired pixels per event.
int m_CutMinSeedCharge
cut for accepting to filtered hitmap histogram, using cluster seed
std::vector< TH1D * > m_clusterSizeUV
Cluster size.
StoreObjPtr< PXDDAQStatus > m_storeDAQEvtStats
Input array for DAQ Status.
std::vector< TH2D * > m_hitMap
Hitmaps pixels.
std::vector< TH1D * > m_clusterEnergy
Deposited energy of clusters.
std::string m_storePXDClustersName
PXDClusters StoreArray name.
std::vector< TH1D * > m_hitMapVCl
Hitmaps clusters for v.
std::vector< TH1D * > m_chargStartRow
Cluster seed charge by distance from the start row.
std::vector< TH1D * > m_startRow
Start row distribution.
TH1D * m_hitMapFilterCounts
Hitmaps of filtered Digits.
void defineHisto() override final
Histogram definitions such as TH1(), TH2(), TNtuple(), TTree()....
int m_CutMinCharge
cut for accepting filtered pixel
std::vector< TH1D * > m_hitMapU
Hitmaps pixels for u.
void event() override final
Event.
std::vector< TH1D * > m_pixelSignal
Charge of pixels.
std::vector< TH1D * > m_clusterSizeV
v cluster size
TH1D * m_hitMapCounts
Hitmaps of Digits.
TH1D * m_hitMapCountsChip
Hitmaps of digits on chips.
TH1D * m_hitMapClCounts
Hitmaps of Clusters.
std::string m_histogramDirectoryName
Name of the histogram directory in ROOT file.
int m_CutMinClusterCharge
cut for accepting filtered cluster, using cluster charge
std::vector< TH1D * > m_startRowCount
counter for Cluster seed charge by distance from the start row
std::vector< TH1D * > m_hitMapV
Hitmaps pixels for v.
std::vector< TH1D * > m_goodclusters
filtered Clusters per event
std::vector< TH1D * > m_clusterSizeU
u cluster size
std::vector< TH1D * > m_clusters
Clusters per event.
void beginRun() override final
Begin run.
int m_CutMaxClusterSize
cut for accepting to filtered hitmap histogram, maximum cluster size
std::vector< TH1D * > m_seed
Seed.
TH1D * m_hitMapClCountsChip
Hitmaps of clusters on chips.
std::vector< TH1D * > m_clusterCharge
Charge of clusters.
std::string m_storePXDDigitsName
PXDDigits StoreArray name.
The PXD digit class.
Definition PXDDigit.h:27
float getADUToEnergy(VxdID id, unsigned int uid, unsigned int vid) const
Get conversion factor from ADU to energy.
static PXDGainCalibrator & getInstance()
Main (and only) way to access the PXDGainCalibrator.
static int getDCDID(const int u, const int v, const VxdID sensorID)
get ID of DCD for giving pixel, range: 1..4.
static int getSWBID(const int v)
get ID of SWB for giving pixel, range: 1..6.
Specific implementation of SensorInfo for PXD Sensors which provides additional pixel specific inform...
Definition SensorInfo.h:23
Low-level class to create/modify relations between StoreArrays.
const std::string & getName() const
Return name under which the object is saved in the DataStore.
Accessor to arrays stored in the data store.
Definition StoreArray.h:113
int getEntries() const
Get the number of objects in the array.
Definition StoreArray.h:216
static const double keV
[kiloelectronvolt]
Definition Unit.h:113
const SensorInfoBase & getSensorInfo(Belle2::VxdID id) const
Return a reference to the SensorInfo of a given SensorID.
Definition GeoCache.cc:67
static GeoCache & getInstance()
Return a reference to the singleton instance.
Definition GeoCache.cc:214
const GeoTools * getGeoTools()
Return a raw pointer to a GeoTools object.
Definition GeoCache.h:141
unsigned short getNumberOfPXDSensors() const
Get number of PXD sensors.
Definition GeoTools.h:133
int getVCells() const
Return number of pixel/strips in v direction.
int getUCells() const
Return number of pixel/strips in u direction.
int getVCellID(double v, bool clamp=false) const
Return the corresponding pixel/strip ID of a given v coordinate.
int getUCellID(double u, double v=0, bool clamp=false) const
Return the corresponding pixel/strip ID of a given u coordinate.
Class to uniquely identify a any structure of the PXD and SVD.
Definition VxdID.h:32
baseType getLayerNumber() const
Get the layer id.
Definition VxdID.h:95
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.
Abstract base class for different kinds of events.
STL namespace.