Belle II Software  release-05-01-25
PXDRawDQMCorrModule.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2016 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Bjoern Spruck *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include "pxd/modules/pxdDQM/PXDRawDQMCorrModule.h"
12 
13 #include <framework/core/HistoModule.h>
14 
15 #include <framework/datastore/StoreArray.h>
16 
17 #include <vector>
18 #include <boost/format.hpp>
19 
20 #include "TH1F.h"
21 #include "TH2F.h"
22 #include "TDirectory.h"
23 
24 using namespace std;
25 using boost::format;
26 using namespace Belle2;
27 
28 //-----------------------------------------------------------------
29 // Register the Module
30 //-----------------------------------------------------------------
31 REG_MODULE(PXDRawDQMCorr)
32 
33 
34 //-----------------------------------------------------------------
35 // Implementation
36 //-----------------------------------------------------------------
37 
39 {
40  //Set module properties
41  setDescription("PXD DQM Correlation module");
42  setPropertyFlags(c_ParallelProcessingCertified); // specify this flag if you need parallel processing
43  addParam("PXDRawHitsName", m_storeRawHitsName, "The name of the StoreArray of PXDRawHits to be processed", string(""));
44  addParam("histogramDirectoryName", m_histogramDirectoryName, "Name of the directory where histograms will be placed",
45  std::string("pxdrawcorr"));
46 }
47 
48 
49 //------------------------------------------------------------------
50 // Function to define histograms
51 //-----------------------------------------------------------------
52 
53 void PXDRawDQMCorrModule::defineHisto()
54 {
55  // Create a separate histogram directory and cd into it.
56  TDirectory* oldDir = gDirectory;
57  if (m_histogramDirectoryName != "") {
58  oldDir->mkdir(m_histogramDirectoryName.c_str());
59  oldDir->cd(m_histogramDirectoryName.c_str());
60  }
61  //----------------------------------------------------------------
62 
63  int nPixelsU1 = 0;
64  int nPixelsV1 = 0;
65  int nPixelsU2 = 250;
66  int nPixelsV2 = 768;
67  m_CorrelationU = new TH2F("RawCorrelationU", "Raw Correlation of U;U1/cells;U2/cells", nPixelsU1 + nPixelsU2, nPixelsU1, nPixelsU2,
68  nPixelsU1 + nPixelsU2, nPixelsU1, nPixelsU2);
69  m_CorrelationV = new TH2F("RawCorrelationV", "Raw Correlation of V;V1/cells;V2/cells", nPixelsV1 + nPixelsV2, nPixelsV1, nPixelsV2,
70  nPixelsV1 + nPixelsV2, nPixelsV1, nPixelsV2);
71  m_DeltaU = new TH1F("RawDeltaU", "Raw Correlation of U2-U1;Udiff/cells", 2 * nPixelsU1 + 2 * nPixelsU2, -nPixelsU2 + nPixelsU1,
72  nPixelsU2 - nPixelsU1);
73  m_DeltaV = new TH1F("RawDeltaV", "Raw Correlation of V2-V1;Vdiff/cells", 2 * nPixelsV1 + 2 * nPixelsV2, -nPixelsV2 + nPixelsV1,
74  nPixelsV2 - nPixelsV1);
75 
76  m_In1CorrelationU = new TH2F("RawIn1CorrelationU", "Raw In1 Correlation of U;U1/cells;U2/cells", nPixelsU1 + nPixelsU2, nPixelsU1,
77  nPixelsU2, nPixelsU1 + nPixelsU2, nPixelsU1, nPixelsU2);
78  m_In1CorrelationV = new TH2F("RawIn1CorrelationV", "Raw In1 Correlation of V;V1/cells;V2/cells", nPixelsV1 + nPixelsV2, nPixelsV1,
79  nPixelsV2, nPixelsV1 + nPixelsV2, nPixelsV1, nPixelsV2);
80  m_In1DeltaU = new TH1F("RawIn1DeltaU", "Raw In1 Correlation of U2-U1;Udiff/cells", 2 * nPixelsU1 + 2 * nPixelsU2,
81  -nPixelsU2 + nPixelsU1, nPixelsU2 - nPixelsU1);
82  m_In1DeltaV = new TH1F("RawIn1DeltaV", "Raw In1 Correlation of V2-V1;Vdiff/cells", 2 * nPixelsV1 + 2 * nPixelsV2,
83  -nPixelsV2 + nPixelsV1, nPixelsV2 - nPixelsV1);
84 
85  m_In2CorrelationU = new TH2F("RawIn2CorrelationU", "Raw In2 Correlation of U;U1/cells;U2/cells", nPixelsU1 + nPixelsU2, nPixelsU1,
86  nPixelsU2, nPixelsU1 + nPixelsU2, nPixelsU1, nPixelsU2);
87  m_In2CorrelationV = new TH2F("RawIn2CorrelationV", "Raw In2 Correlation of V;V1/cells;V2/cells", nPixelsV1 + nPixelsV2, nPixelsV1,
88  nPixelsV2, nPixelsV1 + nPixelsV2, nPixelsV1, nPixelsV2);
89  m_In2DeltaU = new TH1F("RawIn2DeltaU", "Raw In2 Correlation of U2-U1;Udiff/cells", 2 * nPixelsU1 + 2 * nPixelsU2,
90  -nPixelsU2 + nPixelsU1, nPixelsU2 - nPixelsU1);
91  m_In2DeltaV = new TH1F("RawIn2DeltaV", "Raw In2 Correlation of V2-V1;Vdiff/cells", 2 * nPixelsV1 + 2 * nPixelsV2,
92  -nPixelsV2 + nPixelsV1, nPixelsV2 - nPixelsV1);
93 
94  // cd back to root directory
95  oldDir->cd();
96 }
97 
98 
99 void PXDRawDQMCorrModule::initialize()
100 {
101  // Register histograms (calls back defineHisto)
102  REG_HISTOGRAM
103 
104  m_storeRawHits.isRequired(m_storeRawHitsName);
105 }
106 
107 void PXDRawDQMCorrModule::beginRun()
108 {
109  // Just to make sure, reset all the histograms.
110  m_CorrelationU->Reset();
111  m_CorrelationV->Reset();
112  m_DeltaU->Reset();
113  m_DeltaV->Reset();
114 
115  m_In1CorrelationU->Reset();
116  m_In1CorrelationV->Reset();
117  m_In1DeltaU->Reset();
118  m_In1DeltaV->Reset();
119 
120  m_In2CorrelationU->Reset();
121  m_In2CorrelationV->Reset();
122  m_In2DeltaU->Reset();
123  m_In2DeltaV->Reset();
124 
125 }
126 
127 
128 void PXDRawDQMCorrModule::event()
129 {
130  for (auto& hit1 : m_storeRawHits) {
131  int iPlane1 = hit1.getSensorID().getLayerNumber();
132  if (iPlane1 == 0) {
133  for (auto& hit2 : m_storeRawHits) {
134  int iPlane2 = hit2.getSensorID().getLayerNumber();
135  if (iPlane2 == 1) {
136  m_CorrelationU->Fill(hit1.getColumn(), hit2.getColumn());
137  m_CorrelationV->Fill(hit1.getRow(), hit2.getRow());
138  m_DeltaU->Fill(hit2.getColumn() - hit1.getColumn());
139  m_DeltaV->Fill(hit2.getRow() - hit1.getRow());
140  } else { // iPlane2=0
141  if (hit1.getColumn() != hit2.getColumn()) {
142  m_In1CorrelationU->Fill(hit1.getColumn(), hit2.getColumn());
143  m_In1DeltaU->Fill(hit2.getColumn() - hit1.getColumn());
144  }
145  if (hit1.getRow() != hit2.getRow()) {
146  m_In1CorrelationV->Fill(hit1.getRow(), hit2.getRow());
147  m_In1DeltaV->Fill(hit2.getRow() - hit1.getRow());
148  }
149  }
150  }
151  } else { // index 1=1
152  for (auto& hit2 : m_storeRawHits) {
153  int iPlane2 = hit2.getSensorID().getLayerNumber();
154  if (iPlane2 == 1) {
155  if (hit1.getColumn() != hit2.getColumn()) {
156  m_In2CorrelationU->Fill(hit1.getColumn(), hit2.getColumn());
157  m_In2DeltaU->Fill(hit2.getColumn() - hit1.getColumn());
158  }
159  if (hit1.getRow() != hit2.getRow()) {
160  m_In2CorrelationV->Fill(hit1.getRow(), hit2.getRow());
161  m_In2DeltaV->Fill(hit2.getRow() - hit1.getRow());
162  }
163  }
164  }
165  }
166  }
167 }
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::PXDRawDQMCorrModule
PXD DQM Corr Module.
Definition: PXDRawDQMCorrModule.h:40
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::HistoModule
HistoModule.h is supposed to be used instead of Module.h for the modules with histogram definitions t...
Definition: HistoModule.h:29