Belle II Software  release-05-02-19
ROIDQMModule.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2013 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Giulia Casarosa, Eugenio Paoloni, Bjoern Spruck *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include <tracking/modules/pxdDataReduction/ROIDQMModule.h>
12 #include <vxd/geometry/GeoCache.h>
13 
14 #include <TDirectory.h>
15 #include <TH2F.h>
16 
17 using namespace std;
18 using namespace Belle2;
19 
20 //-----------------------------------------------------------------
21 // Register the Module
22 //-----------------------------------------------------------------
23 REG_MODULE(ROIDQM)
24 
25 //-----------------------------------------------------------------
26 // Implementation
27 //-----------------------------------------------------------------
28 
29 ROIDQMModule::ROIDQMModule()
30  : HistoModule()
31  , m_InterDir(NULL)
32  , m_ROIDir(NULL)
33  , m_hInterDictionary(40, [](const Belle2::VxdID & vxdid) {return (size_t)vxdid.getID(); })
34 , m_hROIDictionary(40, [](const Belle2::VxdID& vxdid) {return (size_t)vxdid.getID(); })
35 , m_hROIDictionaryEvt(40, [](const Belle2::VxdID& vxdid) {return (size_t)vxdid.getID(); })
36 , m_numModules(0)
37 , m_hnROIs(NULL)
38 , m_hnInter(NULL)
39 , m_harea(NULL)
40 , m_hredFactor(NULL)
41 {
42  //Set module properties
43  setDescription("Monitor of the ROI creation on HLT");
44  setPropertyFlags(c_ParallelProcessingCertified);
45 
46  addParam("PXDDigitsName", m_PXDDigitsName,
47  "name of the list of PXDDigits", std::string(""));
48 
49  addParam("InterceptsName", m_InterceptsName,
50  "name of the list of Interceptions", std::string(""));
51 
52  addParam("ROIsName", m_ROIsName,
53  "name of the list of ROIs", std::string(""));
54 
55 }
56 
58 {
59 
60  // Create a separate histogram directory and cd into it.
61  TDirectory* oldDir = gDirectory;
62  oldDir->mkdir("intercept");
63  oldDir->cd("intercept");
64  m_InterDir = gDirectory;
65  oldDir->mkdir("roi");
66  oldDir->cd("roi");
67  m_ROIDir = gDirectory;
68 
69  m_InterDir->cd();
70  m_hnInter = new TH1F("hnInter", "number of intercepts", 100, 0, 100);
71 
72  m_ROIDir->cd();
73  m_hnROIs = new TH1F("hnROIs", "number of ROIs", 100, 0, 100);
74  m_harea = new TH1F("harea", "ROI area", 100, 0, 100000);
75  m_hredFactor = new TH1F("hredFactor", "ROI reduction factor", 1000, 0, 1);
76 
78 
79  oldDir->cd();
80 
81 }
82 
84 {
85  REG_HISTOGRAM
86 
87  m_pxdDigits.isOptional();
88  m_roiIDs.isRequired(m_ROIsName);
89  m_pxdIntercept.isRequired(m_InterceptsName);
90 
91 }
92 
94 {
95 
96  m_hnInter->Fill(m_pxdIntercept.getEntries());
97 
98  for (auto& it : m_pxdIntercept)
100 
101 
102  for (auto it = m_hROIDictionaryEvt.begin(); it != m_hROIDictionaryEvt.end(); ++it)
103  (it->second).value = 0;
104 
105  int ROIarea = 0;
106  double redFactor = 0;
107 
108  for (auto& it : m_roiIDs) {
109  fillSensorROIHistos(&it);
110 
111  const VXD::SensorInfoBase& aSensorInfo = m_aGeometry.getSensorInfo(it.getSensorID());
112  const int nPixelsU = aSensorInfo.getUCells();
113  const int nPixelsV = aSensorInfo.getVCells();
114 
115  const int minU = it.getMinUid();
116  const int minV = it.getMinVid();
117  const int maxU = it.getMaxUid();
118  const int maxV = it.getMaxVid();
119 
120  int tmpROIarea = (maxU - minU) * (maxV - minV);
121  ROIarea += tmpROIarea;
122  redFactor += (double)tmpROIarea / (nPixelsU * nPixelsV * m_numModules);
123 
124  }
125 
126  m_hnROIs->Fill(m_roiIDs.getEntries());
127 
128  m_harea->Fill((double)ROIarea);
129 
130  m_hredFactor->Fill((double)redFactor);
131 
132 
133  for (auto it = m_hROIDictionaryEvt.begin(); it != m_hROIDictionaryEvt.end(); ++it) {
134  ROIHistoAccumulateAndFill aROIHistoAccumulateAndFill = it->second;
135  aROIHistoAccumulateAndFill.fill(aROIHistoAccumulateAndFill.hPtr, aROIHistoAccumulateAndFill.value);
136  }
137 
138 }
139 
141 {
142 
143  // VXD::GeoCache& aGeometry = VXD::GeoCache::getInstance();
144 
145  string name; //name of the histogram
146  string title; //title of the histogram
147  TH2F* tmp2D; //temporary 2D histo used to set axis title
148  TH1F* tmp1D; //temporary 1D histo used to set axis title
149 
150  m_numModules = 0;
151 
152  std::set<Belle2::VxdID> pxdLayers = m_aGeometry.getLayers(VXD::SensorInfoBase::PXD);
153  std::set<Belle2::VxdID>::iterator itPxdLayers = pxdLayers.begin();
154 
155  while (itPxdLayers != pxdLayers.end()) {
156 
157  std::set<Belle2::VxdID> pxdLadders = m_aGeometry.getLadders(*itPxdLayers);
158  std::set<Belle2::VxdID>::iterator itPxdLadders = pxdLadders.begin();
159 
160  while (itPxdLadders != pxdLadders.end()) {
161 
162  std::set<Belle2::VxdID> pxdSensors = m_aGeometry.getSensors(*itPxdLadders);
163  std::set<Belle2::VxdID>::iterator itPxdSensors = pxdSensors.begin();
164 
165  while (itPxdSensors != pxdSensors.end()) {
166 
167  m_numModules++; //counting the total number of modules
168 
169  const VXD::SensorInfoBase& wSensorInfo = m_aGeometry.getSensorInfo(*itPxdSensors);
170 
171  const int nPixelsU = wSensorInfo.getUCells();
172  const int nPixelsV = wSensorInfo.getVCells();
173  string sensorid = std::to_string(itPxdSensors->getLayerNumber()) + "_" + std::to_string(itPxdSensors->getLadderNumber()) + "_" +
174  std::to_string(itPxdSensors->getSensorNumber());
175 
176 
177  // ------ HISTOGRAMS WITH AN ACCUMULATE PER ROI AND A FILL PER EVENT -------
178  m_ROIDir->cd();
179 
180  name = "hNROIs_" + sensorid;
181  title = "number of m_roiIDs for sensor " + sensorid;
182  double value = 0;
184  new TH1F(name.c_str(), title.c_str(), 25, 0, 25),
185  [](const ROIid*, double & val) {val++;},
186  [](TH1 * hPtr, double & val) { hPtr->Fill(val); },
187  value
188  };
189  m_hROIDictionaryEvt.insert(pair< Belle2::VxdID, ROIHistoAccumulateAndFill& > ((Belle2::VxdID)*itPxdSensors, *aHAAF));
190 
191 
192 
193 
194  // ------ HISTOGRAMS WITH A FILL PER INTERCEPT -------
195  m_InterDir->cd();
196 
197  // coor U and V
198  name = "hCoorU_" + sensorid;
199  title = "U coordinate of the extrapolation in U for sensor " + sensorid;
200  m_hInterDictionary.insert(pair< Belle2::VxdID, InterHistoAndFill >
201  (
202  (Belle2::VxdID)*itPxdSensors,
204  new TH1F(name.c_str(), title.c_str(), 100, -5, 5),
205  [](TH1 * hPtr, const PXDIntercept * inter) { hPtr->Fill(inter->getCoorU()); }
206  )
207  )
208  );
209 
210  name = "hCoorV_" + sensorid;
211  title = "V coordinate of the extrapolation in V for sensor " + sensorid;
212  m_hInterDictionary.insert(pair< Belle2::VxdID, InterHistoAndFill >
213  (
214  (Belle2::VxdID)*itPxdSensors,
216  new TH1F(name.c_str(), title.c_str(), 100, -5, 5),
217  [](TH1 * hPtr, const PXDIntercept * inter) { hPtr->Fill(inter->getCoorV()); }
218  )
219  )
220  );
221 
222  // Intercept U vs V coordinate
223  name = "hCoorU_vs_CoorV_" + sensorid;
224  title = "U vs V intercept (cm) " + sensorid;
225  tmp2D = new TH2F(name.c_str(), title.c_str(), 100, -5, 5, 100, -5, 5);
226  tmp2D->GetXaxis()->SetTitle("intercept U coor (cm)");
227  tmp2D->GetYaxis()->SetTitle("intercept V coor (cm)");
228  m_hInterDictionary.insert(pair< Belle2::VxdID, InterHistoAndFill >
229  (
230  (Belle2::VxdID)*itPxdSensors,
232  tmp2D,
233  [](TH1 * hPtr, const PXDIntercept * inter) { hPtr->Fill(inter->getCoorU(), inter->getCoorV()); }
234  )
235  )
236  );
237 
238 
239  // sigma U and V
240  name = "hStatErrU_" + sensorid;
241  title = "stat error of the extrapolation in U for sensor " + sensorid;
242  m_hInterDictionary.insert(pair< Belle2::VxdID, InterHistoAndFill >
243  (
244  (Belle2::VxdID)*itPxdSensors,
246  new TH1F(name.c_str(), title.c_str(), 100, 0, 0.35),
247  [](TH1 * hPtr, const PXDIntercept * inter) { hPtr->Fill(inter->getSigmaU()); }
248  )
249  )
250  );
251  name = "hStatErrV_" + sensorid;
252  title = "stat error of the extrapolation in V for sensor " + sensorid;
253  m_hInterDictionary.insert(pair< Belle2::VxdID, InterHistoAndFill >
254  (
255  (Belle2::VxdID)*itPxdSensors,
257  new TH1F(name.c_str(), title.c_str(), 100, 0, 0.35),
258  [](TH1 * hPtr, const PXDIntercept * inter) { hPtr->Fill(inter->getSigmaV()); }
259  )
260  )
261  );
262 
263  //1D residuals
264  name = "hResidU_" + sensorid;
265  title = "U residuals = intercept - digit, for sensor " + sensorid;
266  tmp1D = new TH1F(name.c_str(), title.c_str(), 1000, -5, 5);
267  m_hInterDictionary.insert(pair< Belle2::VxdID, InterHistoAndFill >
268  (
269  (Belle2::VxdID)*itPxdSensors,
271  tmp1D,
272  [this](TH1 * hPtr, const PXDIntercept * inter) {
273  for (auto& it : m_pxdDigits)
274  if ((int)it.getSensorID() == (int)inter->getSensorID()) {
275  const VXD::SensorInfoBase& aSensorInfo = m_aGeometry.getSensorInfo(it.getSensorID());
276  hPtr->Fill(inter->getCoorU() - aSensorInfo.getUCellPosition(it.getUCellID(), it.getVCellID()));
277  }
278  }
279  )
280  )
281  );
282 
283  name = "hResidV_" + sensorid;
284  title = "V residuals = intercept - digit, for sensor " + sensorid;
285  tmp1D = new TH1F(name.c_str(), title.c_str(), 1000, -5, 5);
286  m_hInterDictionary.insert(pair< Belle2::VxdID, InterHistoAndFill >
287  (
288  (Belle2::VxdID)*itPxdSensors,
290  tmp1D,
291  [this](TH1 * hPtr, const PXDIntercept * inter) {
292  for (auto& it : m_pxdDigits)
293  if ((int)it.getSensorID() == (int)inter->getSensorID()) {
294  const VXD::SensorInfoBase& aSensorInfo = m_aGeometry.getSensorInfo(it.getSensorID());
295  hPtr->Fill(inter->getCoorV() - aSensorInfo.getVCellPosition(it.getVCellID()));
296  }
297  }
298  )
299  )
300  );
301 
302  name = "hResidV_vs_ResidU_" + sensorid;
303  title = "V vs U residuals = intercept - digit, for sensor " + sensorid;
304  tmp2D = new TH2F(name.c_str(), title.c_str(), 1000, -5, 5, 1000, -5, 5);
305  tmp2D->GetXaxis()->SetTitle("U resid (cm)");
306  tmp2D->GetYaxis()->SetTitle("V resid (cm)");
307  m_hInterDictionary.insert(pair< Belle2::VxdID, InterHistoAndFill >
308  (
309  (Belle2::VxdID)*itPxdSensors,
311  tmp2D,
312  [this](TH1 * hPtr, const PXDIntercept * inter) {
313  for (auto& it : m_pxdDigits)
314  if ((int)it.getSensorID() == (int)inter->getSensorID()) {
315  const VXD::SensorInfoBase& aSensorInfo = m_aGeometry.getSensorInfo(it.getSensorID());
316  double residU = inter->getCoorU() - aSensorInfo.getUCellPosition(it.getUCellID(), it.getVCellID());
317  double residV = inter->getCoorV() - aSensorInfo.getVCellPosition(it.getVCellID());
318  hPtr->Fill(residU, residV);
319  }
320  }
321  )
322  )
323  );
324 
325  name = "hResidVm_vs_ResidU_" + sensorid;
326  title = "V vs U residuals = intercept - digit, for sensor " + sensorid;
327  tmp2D = new TH2F(name.c_str(), title.c_str(), 1000, -5, 5, 1000, -5, 5);
328  tmp2D->GetXaxis()->SetTitle("U resid (cm)");
329  tmp2D->GetYaxis()->SetTitle("V* resid (cm)");
330  m_hInterDictionary.insert(pair< Belle2::VxdID, InterHistoAndFill >
331  (
332  (Belle2::VxdID)*itPxdSensors,
334  tmp2D,
335  [this](TH1 * hPtr, const PXDIntercept * inter) {
336  for (auto& it : m_pxdDigits)
337  if ((int)it.getSensorID() == (int)inter->getSensorID()) {
338  const VXD::SensorInfoBase& aSensorInfo = m_aGeometry.getSensorInfo(it.getSensorID());
339  double residU = inter->getCoorU() - aSensorInfo.getUCellPosition(it.getUCellID(), it.getVCellID());
340  double residV = inter->getCoorV() + aSensorInfo.getVCellPosition(it.getVCellID());
341  hPtr->Fill(residU, residV);
342  }
343  }
344  )
345  )
346  );
347 
348  name = "hResidV_vs_ResidUm_" + sensorid;
349  title = "V vs U residuals = intercept - digit, for sensor " + sensorid;
350  tmp2D = new TH2F(name.c_str(), title.c_str(), 1000, -5, 5, 1000, -5, 5);
351  tmp2D->GetXaxis()->SetTitle("U* resid (cm)");
352  tmp2D->GetYaxis()->SetTitle("V resid (cm)");
353  m_hInterDictionary.insert(pair< Belle2::VxdID, InterHistoAndFill >
354  (
355  (Belle2::VxdID)*itPxdSensors,
357  tmp2D,
358  [this](TH1 * hPtr, const PXDIntercept * inter) {
359  for (auto& it : m_pxdDigits)
360  if ((int)it.getSensorID() == (int)inter->getSensorID()) {
361  const VXD::SensorInfoBase& aSensorInfo = m_aGeometry.getSensorInfo(it.getSensorID());
362  double residU = inter->getCoorU() + aSensorInfo.getUCellPosition(it.getUCellID(), it.getVCellID());
363  double residV = inter->getCoorV() - aSensorInfo.getVCellPosition(it.getVCellID());
364  hPtr->Fill(residU, residV);
365  }
366  }
367  )
368  )
369  );
370 
371  name = "hResidVm_vs_ResidUm_" + sensorid;
372  title = "V vs U residuals = intercept - digit, for sensor " + sensorid;
373  tmp2D = new TH2F(name.c_str(), title.c_str(), 1000, -5, 5, 1000, -5, 5);
374  tmp2D->GetXaxis()->SetTitle("U* resid (cm)");
375  tmp2D->GetYaxis()->SetTitle("V* resid (cm)");
376  m_hInterDictionary.insert(pair< Belle2::VxdID, InterHistoAndFill >
377  (
378  (Belle2::VxdID)*itPxdSensors,
380  tmp2D,
381  [this](TH1 * hPtr, const PXDIntercept * inter) {
382  for (auto& it : m_pxdDigits)
383  if ((int)it.getSensorID() == (int)inter->getSensorID()) {
384  const VXD::SensorInfoBase& aSensorInfo = m_aGeometry.getSensorInfo(it.getSensorID());
385  double residU = inter->getCoorU() + aSensorInfo.getUCellPosition(it.getUCellID(), it.getVCellID());
386  double residV = inter->getCoorV() + aSensorInfo.getVCellPosition(it.getVCellID());
387  hPtr->Fill(residU, residV);
388  }
389  }
390  )
391  )
392  );
393 
394  //residual U,V vs coordinate U,V
395  name = "hResidU_vs_CoorU_" + sensorid;
396  title = "U residual (cm) vs coor U (cm) " + sensorid;
397  tmp2D = new TH2F(name.c_str(), title.c_str(), 1000, -5, 5, 1000, -5, 5);
398  tmp2D->GetYaxis()->SetTitle("U resid (cm)");
399  tmp2D->GetXaxis()->SetTitle("U coor (cm)");
400  m_hInterDictionary.insert(pair< Belle2::VxdID, InterHistoAndFill >
401  (
402  (Belle2::VxdID)*itPxdSensors,
404  tmp2D,
405  [this](TH1 * hPtr, const PXDIntercept * inter) {
406  for (auto& it : m_pxdDigits)
407  if ((int)it.getSensorID() == (int)inter->getSensorID()) {
408  const VXD::SensorInfoBase& aSensorInfo = m_aGeometry.getSensorInfo(it.getSensorID());
409  double resid = inter->getCoorU() - aSensorInfo.getUCellPosition(it.getUCellID(), it.getVCellID());
410  hPtr->Fill(inter->getCoorU(), resid);
411  }
412  }
413  )
414  )
415  );
416 
417  name = "hResidV_vs_CoorV_" + sensorid;
418  title = "V residual (cm) vs coor V (cm) " + sensorid;
419  tmp2D = new TH2F(name.c_str(), title.c_str(), 1000, -5, 5, 1000, -5, 5);
420  tmp2D->GetYaxis()->SetTitle("V resid (cm)");
421  tmp2D->GetXaxis()->SetTitle("V coor (cm)");
422  m_hInterDictionary.insert(pair< Belle2::VxdID, InterHistoAndFill >
423  (
424  (Belle2::VxdID)*itPxdSensors,
426  tmp2D,
427  [this](TH1 * hPtr, const PXDIntercept * inter) {
428  for (auto& it : m_pxdDigits)
429  if ((int)it.getSensorID() == (int)inter->getSensorID()) {
430  const VXD::SensorInfoBase& aSensorInfo = m_aGeometry.getSensorInfo(it.getSensorID());
431  double resid = inter->getCoorV() - aSensorInfo.getVCellPosition(it.getVCellID());
432  hPtr->Fill(inter->getCoorV(), resid);
433  }
434  }
435  )
436  )
437  );
438 
439  //residual U,V vs coordinate V,U
440  name = "hResidU_vs_CoorV_" + sensorid;
441  title = "U residual (cm) vs coor V (cm) " + sensorid;
442  tmp2D = new TH2F(name.c_str(), title.c_str(), 1000, -5, 5, 1000, -5, 5);
443  tmp2D->GetYaxis()->SetTitle("U resid (cm)");
444  tmp2D->GetXaxis()->SetTitle("V coor (cm)");
445  m_hInterDictionary.insert(pair< Belle2::VxdID, InterHistoAndFill >
446  (
447  (Belle2::VxdID)*itPxdSensors,
449  tmp2D,
450  [this](TH1 * hPtr, const PXDIntercept * inter) {
451  for (auto& it : m_pxdDigits)
452  if ((int)it.getSensorID() == (int)inter->getSensorID()) {
453  const VXD::SensorInfoBase& aSensorInfo = m_aGeometry.getSensorInfo(it.getSensorID());
454  double resid = inter->getCoorU() - aSensorInfo.getUCellPosition(it.getUCellID(), it.getVCellID());
455  hPtr->Fill(inter->getCoorV(), resid);
456  }
457  }
458  )
459  )
460  );
461 
462  name = "hResidV_vs_CoorU_" + sensorid;
463  title = "V residual (cm) vs coor U (cm) " + sensorid;
464  tmp2D = new TH2F(name.c_str(), title.c_str(), 1000, -5, 5, 1000, -5, 5);
465  tmp2D->GetYaxis()->SetTitle("V resid (cm)");
466  tmp2D->GetXaxis()->SetTitle("U coor (cm)");
467  m_hInterDictionary.insert(pair< Belle2::VxdID, InterHistoAndFill >
468  (
469  (Belle2::VxdID)*itPxdSensors,
471  tmp2D,
472  [this](TH1 * hPtr, const PXDIntercept * inter) {
473  for (auto& it : m_pxdDigits)
474  if ((int)it.getSensorID() == (int)inter->getSensorID()) {
475  const VXD::SensorInfoBase& aSensorInfo = m_aGeometry.getSensorInfo(it.getSensorID());
476  double resid = inter->getCoorV() - aSensorInfo.getVCellPosition(it.getVCellID());
477  hPtr->Fill(inter->getCoorU(), resid);
478  }
479  }
480  )
481  )
482  );
483 
484 
485 
486  //residual vs charge
487  name = "hResidU_vs_charge_" + sensorid;
488  title = "U residual (cm) vs charge " + sensorid;
489  tmp2D = new TH2F(name.c_str(), title.c_str(), 250, 0, 250, 100, -5, 5);
490  tmp2D->GetYaxis()->SetTitle("U resid (cm)");
491  tmp2D->GetXaxis()->SetTitle("charge");
492  m_hInterDictionary.insert(pair< Belle2::VxdID, InterHistoAndFill >
493  (
494  (Belle2::VxdID)*itPxdSensors,
496  tmp2D,
497  [this](TH1 * hPtr, const PXDIntercept * inter) {
498  for (auto& it : m_pxdDigits)
499  if ((int)it.getSensorID() == (int)inter->getSensorID()) {
500  const VXD::SensorInfoBase& aSensorInfo = m_aGeometry.getSensorInfo(it.getSensorID());
501  double resid = inter->getCoorU() - aSensorInfo.getUCellPosition(it.getUCellID(), it.getVCellID());
502  hPtr->Fill(it.getCharge(), resid);
503  }
504  }
505  )
506  )
507  );
508 
509  name = "hResidV_vs_charge_" + sensorid;
510  title = "V residual (cm) vs charge " + sensorid;
511  tmp2D = new TH2F(name.c_str(), title.c_str(), 250, 0, 250, 100, -5, 5);
512  tmp2D->GetYaxis()->SetTitle("V resid (cm)");
513  tmp2D->GetXaxis()->SetTitle("charge");
514  m_hInterDictionary.insert(pair< Belle2::VxdID, InterHistoAndFill >
515  (
516  (Belle2::VxdID)*itPxdSensors,
518  tmp2D,
519  [this](TH1 * hPtr, const PXDIntercept * inter) {
520  for (auto& it : m_pxdDigits)
521  if ((int)it.getSensorID() == (int)inter->getSensorID()) {
522  const VXD::SensorInfoBase& aSensorInfo = m_aGeometry.getSensorInfo(it.getSensorID());
523  double resid = inter->getCoorV() - aSensorInfo.getVCellPosition(it.getVCellID());
524  hPtr->Fill(it.getCharge(), resid);
525  }
526  }
527  )
528  )
529  );
530 
531 
532  // scatter plot: U,V intercept in cm VS U,V cell position
533  name = "hCoorU_vs_UDigit_" + sensorid;
534  title = "U intercept (cm) vs U Digit (ID) " + sensorid;
535  tmp2D = new TH2F(name.c_str(), title.c_str(), 1000, -5, 5, 1000, -5, 5);
536  tmp2D->GetXaxis()->SetTitle("intercept U coor (cm)");
537  tmp2D->GetYaxis()->SetTitle("digit U coor (cm)");
538  m_hInterDictionary.insert(pair< Belle2::VxdID, InterHistoAndFill >
539  (
540  (Belle2::VxdID)*itPxdSensors,
542  tmp2D,
543  [this](TH1 * hPtr, const PXDIntercept * inter) {
544  for (auto& it : m_pxdDigits)
545  if ((int)it.getSensorID() == (int)inter->getSensorID()) {
546  const VXD::SensorInfoBase& aSensorInfo = m_aGeometry.getSensorInfo(it.getSensorID());
547  hPtr->Fill(inter->getCoorU(), aSensorInfo.getUCellPosition(it.getUCellID(), it.getVCellID()));
548  // hPtr->Fill( inter->getCoorU(), it.getVCellID()*75e-4 );
549  }
550  }
551  )
552  )
553  );
554 
555  name = "hCoorV_vs_VDigit_" + sensorid;
556  title = "V intercept (cm) vs V Digit (ID) " + sensorid;
557  tmp2D = new TH2F(name.c_str(), title.c_str(), 1000, -5, 5, 1000, -5, 5);
558  tmp2D->GetXaxis()->SetTitle("intercept V coor (cm)");
559  tmp2D->GetYaxis()->SetTitle("digi V coor (cm)");
560  m_hInterDictionary.insert(pair< Belle2::VxdID, InterHistoAndFill >
561  (
562  (Belle2::VxdID)*itPxdSensors,
564  tmp2D,
565  [this](TH1 * hPtr, const PXDIntercept * inter) {
566  for (auto& it : m_pxdDigits) {
567  if ((int)it.getSensorID() == (int)inter->getSensorID()) {
568  const VXD::SensorInfoBase& aSensorInfo = m_aGeometry.getSensorInfo(it.getSensorID());
569  hPtr->Fill(inter->getCoorV(), aSensorInfo.getVCellPosition(it.getVCellID()));
570  // hPtr->Fill( inter->getCoorV(), it.getUCellID()*50e-4 );
571  }
572  }
573  }
574  )
575  )
576  );
577 
578 
579 
580  // ------ HISTOGRAMS WITH A FILL PER ROI -------
581  m_ROIDir->cd();
582 
583  // MIN in U and V
584  name = "hminU_" + sensorid;
585  title = "ROI min in U for sensor " + sensorid;
586  m_hROIDictionary.insert(pair< Belle2::VxdID, ROIHistoAndFill >
587  (
588  (Belle2::VxdID)*itPxdSensors,
590  new TH1F(name.c_str(), title.c_str(), nPixelsU, 0, nPixelsU),
591  [](TH1 * hPtr, const ROIid * roi) { hPtr->Fill(roi->getMinUid()); }
592  )
593  )
594  );
595  name = "hminV_" + sensorid;
596  title = "ROI min in V for sensor " + sensorid;
597  m_hROIDictionary.insert(pair< Belle2::VxdID, ROIHistoAndFill >
598  (
599  (Belle2::VxdID)*itPxdSensors,
601  new TH1F(name.c_str(), title.c_str(), nPixelsV, 0, nPixelsV),
602  [](TH1 * hPtr, const ROIid * roi) { hPtr->Fill(roi->getMinVid()); }
603  )
604  )
605  );
606  //--------------------------
607  // MAX in U and V
608  name = "hmaxU_" + sensorid;
609  title = "ROI max in U for sensor " + sensorid;
610  m_hROIDictionary.insert(pair< Belle2::VxdID, ROIHistoAndFill >
611  (
612  (Belle2::VxdID)*itPxdSensors,
614  new TH1F(name.c_str(), title.c_str(), nPixelsU, 0, nPixelsU),
615  [](TH1 * hPtr, const ROIid * roi) { hPtr->Fill(roi->getMaxUid()); }
616  )
617  )
618  );
619  name = "hmaxV_" + sensorid;
620  title = "ROI max in V for sensor " + sensorid;
621  m_hROIDictionary.insert(pair< Belle2::VxdID, ROIHistoAndFill >
622  (
623  (Belle2::VxdID)*itPxdSensors,
625  new TH1F(name.c_str(), title.c_str(), nPixelsV, 0, nPixelsV),
626  [](TH1 * hPtr, const ROIid * roi) { hPtr->Fill(roi->getMaxVid()); }
627  )
628  )
629  );
630  //--------------------------
631 
632  // WIDTH in U and V
633  name = "hwidthU_" + sensorid;
634  title = "ROI width in U for sensor " + sensorid;
635  m_hROIDictionary.insert(pair< Belle2::VxdID, ROIHistoAndFill >
636  (
637  (Belle2::VxdID)*itPxdSensors,
639  new TH1F(name.c_str(), title.c_str(), nPixelsU, 0, nPixelsU),
640  [](TH1 * hPtr, const ROIid * roi) { hPtr->Fill(roi->getMaxUid() - roi->getMinUid()); }
641  )
642  )
643  );
644  name = "hwidthV_" + sensorid;
645  title = "ROI width in V for sensor " + sensorid;
646  m_hROIDictionary.insert(pair< Belle2::VxdID, ROIHistoAndFill >
647  (
648  (Belle2::VxdID)*itPxdSensors,
650  new TH1F(name.c_str(), title.c_str(), nPixelsV, 0, nPixelsV),
651  [](TH1 * hPtr, const ROIid * roi) { hPtr->Fill(roi->getMaxVid() - roi->getMinVid()); }
652  )
653  )
654  );
655 
656  // ROI center
657  name = "hROIcenter_" + sensorid;
658  title = "ROI center " + sensorid;
659  tmp2D = new TH2F(name.c_str(), title.c_str(), nPixelsU, 0, nPixelsU, nPixelsV, 0, nPixelsV);
660  tmp2D->GetXaxis()->SetTitle(" U (ID)");
661  tmp2D->GetYaxis()->SetTitle(" V (ID)");
662  m_hROIDictionary.insert(pair< Belle2::VxdID, ROIHistoAndFill >
663  (
664  (Belle2::VxdID)*itPxdSensors,
666  tmp2D,
667  [](TH1 * hPtr, const ROIid * roi) { hPtr->Fill((roi->getMaxUid() + roi->getMinUid()) / 2, (roi->getMaxVid() + roi->getMinVid()) / 2); }
668  )
669  )
670  );
671 
672  //--------------------------
673 
674  ++itPxdSensors;
675  }
676  ++itPxdLadders;
677  }
678  ++itPxdLayers;
679  }
680 
681 }
682 
684 {
685 
686  auto its = m_hInterDictionary.equal_range(inter->getSensorID());
687 
688  for (auto it = its.first; it != its.second; ++it) {
689  InterHistoAndFill aInterHistoAndFill = it->second;
690  aInterHistoAndFill.second(aInterHistoAndFill.first, inter);
691  }
692 
693 }
694 
696 {
697 
698  auto its = m_hROIDictionary.equal_range(roi->getSensorID());
699 
700  for (auto it = its.first; it != its.second; ++it) {
701  ROIHistoAndFill aROIHistoAndFill = it->second;
702  aROIHistoAndFill.second(aROIHistoAndFill.first, roi);
703  }
704 
705  auto itsEvt = m_hROIDictionaryEvt.equal_range(roi->getSensorID());
706  for (auto it = itsEvt.first; it != itsEvt.second; ++it)
707  (it->second).accumulate(roi, (it->second).value);
708 }
709 
711 {
712  for (auto it = m_hROIDictionaryEvt.begin(); it != m_hROIDictionaryEvt.end(); ++it)
713  delete &(it->second);
714 }
Belle2::VXD::SensorInfoBase::getUCells
int getUCells() const
Return number of pixel/strips in u direction.
Definition: SensorInfoBase.h:223
Belle2::VxdID
Class to uniquely identify a any structure of the PXD and SVD.
Definition: VxdID.h:43
Belle2::ROIDQMModule::ROIHistoAccumulateAndFill::value
double value
value used to fill
Definition: ROIDQMModule.h:85
Belle2::ROIDQMModule::ROIHistoAccumulateAndFill::hPtr
TH1 * hPtr
histogram pointer
Definition: ROIDQMModule.h:82
Belle2::ROIDQMModule::m_hnROIs
TH1F * m_hnROIs
number of ROIs
Definition: ROIDQMModule.h:97
Belle2::ROIDQMModule::m_hROIDictionaryEvt
std::unordered_multimap< Belle2::VxdID, ROIHistoAccumulateAndFill &, std::function< size_t(const Belle2::VxdID &) > > m_hROIDictionaryEvt
map of histograms to be filled once per event
Definition: ROIDQMModule.h:89
Belle2::ROIDQMModule::m_numModules
int m_numModules
number of modules
Definition: ROIDQMModule.h:95
Belle2::ROIDQMModule::m_roiIDs
StoreArray< ROIid > m_roiIDs
the ROIids dataobjects collection
Definition: ROIDQMModule.h:57
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::ROIDQMModule::terminate
void terminate(void) override final
Function to terminate module.
Definition: ROIDQMModule.cc:710
Belle2::ROIDQMModule::ROIHistoAccumulateAndFill::fill
std::function< void(TH1 *, double &) > fill
fill function
Definition: ROIDQMModule.h:84
Belle2::ROIDQMModule::m_hInterDictionary
std::unordered_multimap< Belle2::VxdID, InterHistoAndFill, std::function< size_t(const Belle2::VxdID &)> > m_hInterDictionary
map of histograms to be filled once per intercept
Definition: ROIDQMModule.h:73
Belle2::ROIDQMModule::ROIHistoAccumulateAndFill
struct: histograms to be filled once per event + filling fucntion + accumulate function
Definition: ROIDQMModule.h:81
Belle2::ROIDQMModule::fillSensorInterHistos
void fillSensorInterHistos(const PXDIntercept *inter)
fill histograms per sensor, filled once per intercept
Definition: ROIDQMModule.cc:683
Belle2::VXD::SensorInfoBase
Base class to provide Sensor Information for PXD and SVD.
Definition: SensorInfoBase.h:40
Belle2::ROIDQMModule::m_hROIDictionary
std::unordered_multimap< Belle2::VxdID, ROIHistoAndFill, std::function< size_t(const Belle2::VxdID &)> > m_hROIDictionary
map of histograms to be filled once per roi
Definition: ROIDQMModule.h:78
Belle2::VXD::GeoCache::getLayers
const std::set< Belle2::VxdID > getLayers(SensorInfoBase::SensorType sensortype=SensorInfoBase::VXD)
Return a set of all known Layers.
Definition: GeoCache.cc:177
Belle2::ROIDQMModule::m_aGeometry
VXD::GeoCache & m_aGeometry
the geometry
Definition: ROIDQMModule.h:62
Belle2::ROIDQMModule::m_ROIsName
std::string m_ROIsName
Name of the ROIid StoreArray.
Definition: ROIDQMModule.h:64
Belle2::ROIDQMModule::m_pxdIntercept
StoreArray< PXDIntercept > m_pxdIntercept
the PXDIntercepts dataobjects collection
Definition: ROIDQMModule.h:58
Belle2::ROIDQMModule::m_hnInter
TH1F * m_hnInter
number of intercpets
Definition: ROIDQMModule.h:98
Belle2::VXD::GeoCache::getSensors
const std::set< Belle2::VxdID > & getSensors(Belle2::VxdID ladder) const
Return a set of all sensor IDs belonging to a given ladder.
Definition: GeoCache.cc:205
Belle2::ROIDQMModule::m_ROIDir
TDirectory * m_ROIDir
ROI directory in the root file.
Definition: ROIDQMModule.h:68
Belle2::ROIDQMModule::createHistosDictionaries
void createHistosDictionaries()
create the dictionary
Definition: ROIDQMModule.cc:140
Belle2::ROIid
ROIid stores the U and V ids and the sensor id of the Region Of Interest.
Definition: ROIid.h:35
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::ROIDQMModule::m_InterDir
TDirectory * m_InterDir
intercepts directory in the root file
Definition: ROIDQMModule.h:67
Belle2::ROIDQMModule::m_hredFactor
TH1F * m_hredFactor
reduction factor
Definition: ROIDQMModule.h:100
Belle2::ROIDQMModule::fillSensorROIHistos
void fillSensorROIHistos(const ROIid *roi)
fill histograms per sensor, filled once per ROI
Definition: ROIDQMModule.cc:695
Belle2::ROIDQMModule::m_InterceptsName
std::string m_InterceptsName
Name of the PXDIntercept StoreArray.
Definition: ROIDQMModule.h:65
Belle2::ROIDQMModule::initialize
void initialize(void) override final
Function for dynamic initialization of module.
Definition: ROIDQMModule.cc:83
Belle2::VXD::SensorInfoBase::PXD
@ PXD
PXD Sensor.
Definition: SensorInfoBase.h:44
Belle2::VXD::GeoCache::getSensorInfo
const SensorInfoBase & getSensorInfo(Belle2::VxdID id) const
Return a referecne to the SensorInfo of a given SensorID.
Definition: GeoCache.cc:68
Belle2::ROIDQMModule::defineHisto
void defineHisto() override final
define histograms
Definition: ROIDQMModule.cc:57
Belle2::ROIDQMModule::event
void event(void) override final
Function to process event record.
Definition: ROIDQMModule.cc:93
Belle2::ROIDQMModule::m_pxdDigits
StoreArray< PXDDigit > m_pxdDigits
the PXDDigits dataobjects collection
Definition: ROIDQMModule.h:56
Belle2::PXDIntercept
PXDIntercept stores the U,V coordinates and uncertainties of the intersection of a track with an PXD ...
Definition: PXDIntercept.h:32
Belle2::VXD::GeoCache::getLadders
const std::set< Belle2::VxdID > & getLadders(Belle2::VxdID layer) const
Return a set of all ladder IDs belonging to a given layer.
Definition: GeoCache.cc:194
Belle2::ROIDQMModule::InterHistoAndFill
std::pair< TH1 *, std::function< void(TH1 *, const PXDIntercept *) > > InterHistoAndFill
typedef: histograms to be filled once per intercept + filling function
Definition: ROIDQMModule.h:71
Belle2::ROIDQMModule::ROIHistoAndFill
std::pair< TH1 *, std::function< void(TH1 *, const ROIid *) > > ROIHistoAndFill
typedef: histograms to be filled once per roi + filling function
Definition: ROIDQMModule.h:76
Belle2::ROIDQMModule::m_harea
TH1F * m_harea
ROis area.
Definition: ROIDQMModule.h:99
Belle2::HistoModule
HistoModule.h is supposed to be used instead of Module.h for the modules with histogram definitions t...
Definition: HistoModule.h:29