Belle II Software development
PXDRawDQMCorrModule.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/PXDRawDQMCorrModule.h"
10#include <pxd/dataobjects/PXDRawHit.h>
11
12#include <vector>
13#include <boost/format.hpp>
14
15#include "TH1F.h"
16#include "TH2F.h"
17#include "TDirectory.h"
18
19using namespace std;
20using boost::format;
21using namespace Belle2;
22
23//-----------------------------------------------------------------
24// Register the Module
25//-----------------------------------------------------------------
26REG_MODULE(PXDRawDQMCorr);
27
28
29//-----------------------------------------------------------------
30// Implementation
31//-----------------------------------------------------------------
32
34{
35 //Set module properties
36 setDescription("PXD DQM Correlation module");
37 setPropertyFlags(c_ParallelProcessingCertified); // specify this flag if you need parallel processing
38 addParam("PXDRawHitsName", m_storeRawHitsName, "The name of the StoreArray of PXDRawHits to be processed", string(""));
39 addParam("histogramDirectoryName", m_histogramDirectoryName, "Name of the directory where histograms will be placed",
40 std::string("pxdrawcorr"));
41}
42
43
44//------------------------------------------------------------------
45// Function to define histograms
46//-----------------------------------------------------------------
47
49{
50 // Create a separate histogram directory and cd into it.
51 TDirectory* oldDir = gDirectory;
52 if (m_histogramDirectoryName != "") {
53 oldDir->mkdir(m_histogramDirectoryName.c_str());
54 oldDir->cd(m_histogramDirectoryName.c_str());
55 }
56 //----------------------------------------------------------------
57
58 int nPixelsU1 = 0;
59 int nPixelsV1 = 0;
60 int nPixelsU2 = 250;
61 int nPixelsV2 = 768;
62 m_CorrelationU = new TH2F("RawCorrelationU", "Raw Correlation of U;U1/cells;U2/cells", nPixelsU1 + nPixelsU2, nPixelsU1, nPixelsU2,
63 nPixelsU1 + nPixelsU2, nPixelsU1, nPixelsU2);
64 m_CorrelationV = new TH2F("RawCorrelationV", "Raw Correlation of V;V1/cells;V2/cells", nPixelsV1 + nPixelsV2, nPixelsV1, nPixelsV2,
65 nPixelsV1 + nPixelsV2, nPixelsV1, nPixelsV2);
66 m_DeltaU = new TH1F("RawDeltaU", "Raw Correlation of U2-U1;Udiff/cells", 2 * nPixelsU1 + 2 * nPixelsU2, -nPixelsU2 + nPixelsU1,
67 nPixelsU2 - nPixelsU1);
68 m_DeltaV = new TH1F("RawDeltaV", "Raw Correlation of V2-V1;Vdiff/cells", 2 * nPixelsV1 + 2 * nPixelsV2, -nPixelsV2 + nPixelsV1,
69 nPixelsV2 - nPixelsV1);
70
71 m_In1CorrelationU = new TH2F("RawIn1CorrelationU", "Raw In1 Correlation of U;U1/cells;U2/cells", nPixelsU1 + nPixelsU2, nPixelsU1,
72 nPixelsU2, nPixelsU1 + nPixelsU2, nPixelsU1, nPixelsU2);
73 m_In1CorrelationV = new TH2F("RawIn1CorrelationV", "Raw In1 Correlation of V;V1/cells;V2/cells", nPixelsV1 + nPixelsV2, nPixelsV1,
74 nPixelsV2, nPixelsV1 + nPixelsV2, nPixelsV1, nPixelsV2);
75 m_In1DeltaU = new TH1F("RawIn1DeltaU", "Raw In1 Correlation of U2-U1;Udiff/cells", 2 * nPixelsU1 + 2 * nPixelsU2,
76 -nPixelsU2 + nPixelsU1, nPixelsU2 - nPixelsU1);
77 m_In1DeltaV = new TH1F("RawIn1DeltaV", "Raw In1 Correlation of V2-V1;Vdiff/cells", 2 * nPixelsV1 + 2 * nPixelsV2,
78 -nPixelsV2 + nPixelsV1, nPixelsV2 - nPixelsV1);
79
80 m_In2CorrelationU = new TH2F("RawIn2CorrelationU", "Raw In2 Correlation of U;U1/cells;U2/cells", nPixelsU1 + nPixelsU2, nPixelsU1,
81 nPixelsU2, nPixelsU1 + nPixelsU2, nPixelsU1, nPixelsU2);
82 m_In2CorrelationV = new TH2F("RawIn2CorrelationV", "Raw In2 Correlation of V;V1/cells;V2/cells", nPixelsV1 + nPixelsV2, nPixelsV1,
83 nPixelsV2, nPixelsV1 + nPixelsV2, nPixelsV1, nPixelsV2);
84 m_In2DeltaU = new TH1F("RawIn2DeltaU", "Raw In2 Correlation of U2-U1;Udiff/cells", 2 * nPixelsU1 + 2 * nPixelsU2,
85 -nPixelsU2 + nPixelsU1, nPixelsU2 - nPixelsU1);
86 m_In2DeltaV = new TH1F("RawIn2DeltaV", "Raw In2 Correlation of V2-V1;Vdiff/cells", 2 * nPixelsV1 + 2 * nPixelsV2,
87 -nPixelsV2 + nPixelsV1, nPixelsV2 - nPixelsV1);
88
89 // cd back to root directory
90 oldDir->cd();
91}
92
93
95{
96 // Register histograms (calls back defineHisto)
97 REG_HISTOGRAM
98
100}
101
103{
104 // Just to make sure, reset all the histograms.
105 m_CorrelationU->Reset();
106 m_CorrelationV->Reset();
107 m_DeltaU->Reset();
108 m_DeltaV->Reset();
109
110 m_In1CorrelationU->Reset();
111 m_In1CorrelationV->Reset();
112 m_In1DeltaU->Reset();
113 m_In1DeltaV->Reset();
114
115 m_In2CorrelationU->Reset();
116 m_In2CorrelationV->Reset();
117 m_In2DeltaU->Reset();
118 m_In2DeltaV->Reset();
119
120}
121
122
124{
125 for (auto& hit1 : m_storeRawHits) {
126 int iPlane1 = hit1.getSensorID().getLayerNumber();
127 if (iPlane1 == 0) {
128 for (auto& hit2 : m_storeRawHits) {
129 int iPlane2 = hit2.getSensorID().getLayerNumber();
130 if (iPlane2 == 1) {
131 m_CorrelationU->Fill(hit1.getColumn(), hit2.getColumn());
132 m_CorrelationV->Fill(hit1.getRow(), hit2.getRow());
133 m_DeltaU->Fill(hit2.getColumn() - hit1.getColumn());
134 m_DeltaV->Fill(hit2.getRow() - hit1.getRow());
135 } else { // iPlane2=0
136 if (hit1.getColumn() != hit2.getColumn()) {
137 m_In1CorrelationU->Fill(hit1.getColumn(), hit2.getColumn());
138 m_In1DeltaU->Fill(hit2.getColumn() - hit1.getColumn());
139 }
140 if (hit1.getRow() != hit2.getRow()) {
141 m_In1CorrelationV->Fill(hit1.getRow(), hit2.getRow());
142 m_In1DeltaV->Fill(hit2.getRow() - hit1.getRow());
143 }
144 }
145 }
146 } else { // index 1=1
147 for (auto& hit2 : m_storeRawHits) {
148 int iPlane2 = hit2.getSensorID().getLayerNumber();
149 if (iPlane2 == 1) {
150 if (hit1.getColumn() != hit2.getColumn()) {
151 m_In2CorrelationU->Fill(hit1.getColumn(), hit2.getColumn());
152 m_In2DeltaU->Fill(hit2.getColumn() - hit1.getColumn());
153 }
154 if (hit1.getRow() != hit2.getRow()) {
155 m_In2CorrelationV->Fill(hit1.getRow(), hit2.getRow());
156 m_In2DeltaV->Fill(hit2.getRow() - hit1.getRow());
157 }
158 }
159 }
160 }
161 }
162}
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
TH2F * m_In1CorrelationU
Correlation Sensor 1 vs 2.
TH2F * m_CorrelationV
Correlation Sensor 1 vs 2.
void initialize() override final
Initialize.
TH1F * m_In1DeltaV
Correlation Sensor 1 vs 2.
TH2F * m_In2CorrelationV
Correlation Sensor 1 vs 2.
TH1F * m_In2DeltaV
Correlation Sensor 1 vs 2.
void defineHisto() override final
Histogram definitions such as TH1(), TH2()TNtuple(), TTree()....
std::string m_storeRawHitsName
PXDRawHits StoreArray name.
void event() override final
Event.
TH2F * m_In1CorrelationV
Correlation Sensor 1 vs 2.
std::string m_histogramDirectoryName
Name of the histogram directory in ROOT file.
TH1F * m_In1DeltaU
Correlation Sensor 1 vs 2.
TH2F * m_CorrelationU
Correlation Sensor 1 vs 2.
TH1F * m_DeltaV
Correlation Sensor 1 vs 2.
void beginRun() override final
Begin run.
TH1F * m_DeltaU
Correlation Sensor 1 vs 2.
StoreArray< PXDRawHit > m_storeRawHits
Storearray for raw pixels.
TH2F * m_In2CorrelationU
Correlation Sensor 1 vs 2.
TH1F * m_In2DeltaU
Correlation Sensor 1 vs 2.
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
Abstract base class for different kinds of events.
STL namespace.