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