Belle II Software  release-05-01-25
IPDQMExpressRecoModule.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2017 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Gaetano de Marino *
7  * *
8  * *
9  * *
10  * This software is provided "as is" without any warranty. *
11  **************************************************************************/
12 
13 #include <dqm/modules/PhysicsObjectsDQM/IPDQMExpressRecoModule.h>
14 #include <analysis/dataobjects/ParticleList.h>
15 #include <analysis/utility/ReferenceFrame.h>
16 #include <TLorentzVector.h>
17 #include <TDirectory.h>
18 
19 using namespace Belle2;
20 
21 //-----------------------------------------------------------------
22 // Register the Module
23 //-----------------------------------------------------------------
24 REG_MODULE(IPDQMExpressReco)
25 
26 //-----------------------------------------------------------------
27 // Implementation
28 //-----------------------------------------------------------------
29 
30 IPDQMExpressRecoModule::IPDQMExpressRecoModule() : HistoModule()
31 {
32  //Set module properties
33 
34  setDescription("Monitor the position of the beamspot using mu+mu- events");
35  setPropertyFlags(c_ParallelProcessingCertified);
36 
37  addParam("Y4SPListName", m_Y4SPListName, "Name of the Y4S particle list", std::string("Upsilon(4S):IPDQM"));
38 }
39 
40 void IPDQMExpressRecoModule::defineHisto()
41 {
42  TDirectory* oldDir = gDirectory;
43  oldDir->mkdir("IPMonitoring")->cd();
44 
45  m_h_x = new TH1F("Y4S_Vertex.X", "IP position - coord. X", 1000, -0.5, 0.5);
46  m_h_x->SetXTitle("IP_coord. X [cm]");
47  m_h_y = new TH1F("Y4S_Vertex.Y", "IP position - coord. Y", 1000, -0.5, 0.5);
48  m_h_y->SetXTitle("IP_coord. Y [cm]");
49  m_h_z = new TH1F("Y4S_Vertex.Z", "IP position - coord. Z", 2000, -2, 2);
50  m_h_z->SetXTitle("IP_coord. Z [cm]");
51  m_h_px = new TH1F("Y4S_Vertex.pX", "Total momentum in lab. frame - coord. X", 100, -2, 2);
52  m_h_px->SetXTitle("pX [GeV/c]");
53  m_h_py = new TH1F("Y4S_Vertex.pY", "Total momentum in lab. frame - coord. Y", 100, -2, 2);
54  m_h_py->SetXTitle("pY [GeV/c]");
55  m_h_pz = new TH1F("Y4S_Vertex.pZ", "Total momentum in lab. frame - coord. Z", 100, 1, 5);
56  m_h_pz->SetXTitle("pZ [GeV/c]");
57  m_h_E = new TH1F("Y4S_Vertex.E", "Energy in lab. frame", 100, 8, 13);
58  m_h_E->SetXTitle("E [GeV]");
59  m_h_pull = new TH1F("Y4S_pull.Y", "Pull - coord. Y", 100, -10, 10);
60  m_h_pull->SetXTitle("(y_{reco}- y_{med})/#sigma_{y}^{reco}");
61  m_h_pull->GetXaxis()->SetTitleOffset(1.2);
62  m_h_pull->GetXaxis()->SetTitleSize(.04);
63  m_h_y_risol = new TH1F("Y4S_Res.Y", "Resolution - coord. Y", 1000, -1, 1);
64  m_h_y_risol->GetXaxis()->SetTitleOffset(1.2);
65  m_h_y_risol->SetXTitle("(y_{reco}- y_{med}) [cm]");
66  m_h_temp = new TH1F("Y4S_temp.Y", "temp - coord. Y", 2000, -2, 2);
67  m_h_temp->SetXTitle("IP_coord. Y [cm]");
68  m_h_xx = new TH1F("Y4S_Prod.XX", "IP position - prod. XX", 1000, 0, 0.5);
69  m_h_xx->SetXTitle("IP_prod. XX [cm^{2} ]");
70  m_h_yy = new TH1F("Y4S_Prod.YY", "IP position - prod. YY", 1000, 0, 0.5);
71  m_h_yy->SetXTitle("IP_prod. YY [cm^{2} ]");
72  m_h_zz = new TH1F("Y4S_Prod.ZZ", "IP position - prod. ZZ", 2000, 0, 4);
73  m_h_zz->SetXTitle("IP_prod. ZZ [cm^{2} ]");
74  m_h_xy = new TH1F("Y4S_Prod.XY", "IP position - prod. XY", 1000, -1, 1);
75  m_h_xy->SetXTitle("IP_prod. XY [cm^{2} ]");
76  m_h_yz = new TH1F("Y4S_Prod.YZ", "IP position - prod. YZ", 1000, -1, 1);
77  m_h_yz->SetXTitle("IP_prod. YZ [cm^{2} ]");
78  m_h_xz = new TH1F("Y4S_Prod.XZ", "IP position - prod. XZ", 1000, -1, 1);
79  m_h_xz->SetXTitle("IP_prod. XZ [cm^{2} ]");
80  m_h_cov_x_x = new TH1F("Var.X", "X Variance", 500, 0., 0.005);
81  m_h_cov_x_x->SetXTitle("Var. X [cm^{2} ]");
82  m_h_cov_y_y = new TH1F("Var.Y", "Y Variance", 500, 0., 0.005);
83  m_h_cov_y_y->SetXTitle("Var. Y [cm^{2} ]");
84  m_h_cov_z_z = new TH1F("Var.Z", "Z Variance", 500, 0., 0.005);
85  m_h_cov_z_z->SetXTitle("Var. Z [cm^{2} ]");
86  m_h_cov_x_y = new TH1F("Covar.XY", "XY Covariance", 1000, -0.005, 0.005);
87  m_h_cov_x_y->SetXTitle("Var. XY [cm^{2} ]");
88  m_h_cov_x_z = new TH1F("Covar.XZ", "XZ Covariance", 1000, -0.005, 0.005);
89  m_h_cov_x_z->SetXTitle("Covar. XZ [cm^{2} ]");
90  m_h_cov_y_z = new TH1F("Covar.YZ", "YZ Covariance", 1000, -0.005, 0.005);
91  m_h_cov_y_z->SetXTitle("Covar. YZ [cm^{2} ]");
92 
93  oldDir->cd();
94 }
95 
96 
97 void IPDQMExpressRecoModule::initialize()
98 {
99  REG_HISTOGRAM
100 }
101 
102 
103 void IPDQMExpressRecoModule::beginRun()
104 {
105  m_h_x->Reset();
106  m_h_y->Reset();
107  m_h_z->Reset();
108  m_h_pull->Reset();
109  m_h_temp->Reset();
110  m_h_y_risol->Reset();
111  m_h_cov_x_x->Reset();
112  m_h_cov_y_y->Reset();
113  m_h_cov_z_z->Reset();
114  m_h_cov_x_y->Reset();
115  m_h_cov_x_z->Reset();
116  m_h_cov_y_z->Reset();
117  m_h_xx->Reset();
118  m_h_yy->Reset();
119  m_h_zz->Reset();
120  m_h_xy->Reset();
121  m_h_xz->Reset();
122  m_h_yz->Reset();
123  m_h_px->Reset();
124  m_h_py->Reset();
125  m_h_pz->Reset();
126  m_h_E->Reset();
127  m_v_y.clear();
128  m_err_y.clear();
129  m_r = 0;
130 }
131 
132 
133 void IPDQMExpressRecoModule::endRun()
134 {
135  m_h_x->GetXaxis()->SetRangeUser(m_h_x->GetMean(1) - 5 * m_h_x->GetRMS(1), m_h_x->GetMean(1) + 5 * m_h_x->GetRMS(1));
136  m_h_y->GetXaxis()->SetRangeUser(m_h_y->GetMean(1) - 5 * m_h_y->GetRMS(1), m_h_y->GetMean(1) + 5 * m_h_y->GetRMS(1));
137  m_h_z->GetXaxis()->SetRangeUser(m_h_z->GetMean(1) - 5 * m_h_z->GetRMS(1), m_h_z->GetMean(1) + 5 * m_h_z->GetRMS(1));
138  m_h_xy->GetXaxis()->SetRangeUser(m_h_xy->GetMean(1) - 5 * m_h_xy->GetRMS(1), m_h_xy->GetMean(1) + 5 * m_h_xy->GetRMS(1));
139  m_h_xz->GetXaxis()->SetRangeUser(m_h_xz->GetMean(1) - 5 * m_h_xz->GetRMS(1), m_h_xz->GetMean(1) + 5 * m_h_xz->GetRMS(1));
140  m_h_yz->GetXaxis()->SetRangeUser(m_h_yz->GetMean(1) - 5 * m_h_yz->GetRMS(1), m_h_yz->GetMean(1) + 5 * m_h_yz->GetRMS(1));
141  m_h_xx->GetXaxis()->SetRangeUser(m_h_xx->GetMean(1) - 5 * m_h_xx->GetRMS(1), m_h_xx->GetMean(1) + 5 * m_h_xx->GetRMS(1));
142  m_h_yy->GetXaxis()->SetRangeUser(m_h_yy->GetMean(1) - 5 * m_h_yy->GetRMS(1), m_h_yy->GetMean(1) + 5 * m_h_yy->GetRMS(1));
143  m_h_zz->GetXaxis()->SetRangeUser(m_h_zz->GetMean(1) - 5 * m_h_zz->GetRMS(1), m_h_zz->GetMean(1) + 5 * m_h_zz->GetRMS(1));
144  m_h_pull->GetXaxis()->SetRangeUser(m_h_pull->GetMean(1) - 5 * m_h_pull->GetRMS(1), m_h_pull->GetMean(1) + 5 * m_h_pull->GetRMS(1));
145  m_h_y_risol->GetXaxis()->SetRangeUser(m_h_y_risol->GetMean(1) - 5 * m_h_y_risol->GetRMS(1),
146  m_h_y_risol->GetMean(1) + 5 * m_h_y_risol->GetRMS(1));
147  m_h_cov_x_x->GetXaxis()->SetRangeUser(m_h_cov_x_x->GetMean(1) - 5 * m_h_cov_x_x->GetRMS(1),
148  m_h_cov_x_x->GetMean(1) + 5 * m_h_cov_x_x->GetRMS(1));
149  m_h_cov_y_y->GetXaxis()->SetRangeUser(m_h_cov_y_y->GetMean(1) - 5 * m_h_cov_y_y->GetRMS(1),
150  m_h_cov_y_y->GetMean(1) + 5 * m_h_cov_y_y->GetRMS(1));
151  m_h_cov_z_z->GetXaxis()->SetRangeUser(m_h_cov_z_z->GetMean(1) - 5 * m_h_cov_z_z->GetRMS(1),
152  m_h_cov_z_z->GetMean(1) + 5 * m_h_cov_z_z->GetRMS(1));
153  m_h_cov_x_z->GetXaxis()->SetRangeUser(m_h_cov_x_z->GetMean(1) - 5 * m_h_cov_x_z->GetRMS(1),
154  m_h_cov_x_z->GetMean(1) + 5 * m_h_cov_x_z->GetRMS(1));
155  m_h_cov_x_y->GetXaxis()->SetRangeUser(m_h_cov_x_y->GetMean(1) - 5 * m_h_cov_x_y->GetRMS(1),
156  m_h_cov_x_y->GetMean(1) + 5 * m_h_cov_x_y->GetRMS(1));
157  m_h_cov_y_z->GetXaxis()->SetRangeUser(m_h_cov_y_z->GetMean(1) - 5 * m_h_cov_y_z->GetRMS(1),
158  m_h_cov_y_z->GetMean(1) + 5 * m_h_cov_y_z->GetRMS(1));
159 }
160 
161 void IPDQMExpressRecoModule::terminate()
162 {
163 }
164 
165 void IPDQMExpressRecoModule::event()
166 {
167 
168  StoreObjPtr<ParticleList> Y4SParticles(m_Y4SPListName);
169  const auto& frame = ReferenceFrame::GetCurrent();
170 
171 
172  if (Y4SParticles.isValid() && abs(Y4SParticles->getPDGCode()) == 300553) {
173 
174  for (unsigned int i = 0; i < Y4SParticles->getListSize(); i++) {
175 
176  Particle* Y4S = Y4SParticles->getParticle(i);
177  TVector3 IPVertex = frame.getVertex(Y4S);
178  const auto& errMatrix = Y4S->getVertexErrorMatrix();
179  m_h_x->Fill(IPVertex.X());
180  m_h_y->Fill(IPVertex.Y());
181  m_h_z->Fill(IPVertex.Z());
182  m_h_xy->Fill(IPVertex.X()*IPVertex.Y());
183  m_h_xz->Fill(IPVertex.X()*IPVertex.Z());
184  m_h_yz->Fill(IPVertex.Y()*IPVertex.Z());
185  m_h_xx->Fill(IPVertex.X()*IPVertex.X());
186  m_h_yy->Fill(IPVertex.Y()*IPVertex.Y());
187  m_h_zz->Fill(IPVertex.Z()*IPVertex.Z());
188  m_h_temp->Fill(IPVertex.Y());
189  m_h_cov_x_x->Fill(errMatrix(0, 0));
190  m_h_cov_y_y->Fill(errMatrix(1, 1));
191  m_h_cov_z_z->Fill(errMatrix(2, 2));
192  m_h_cov_x_y->Fill(errMatrix(0, 1));
193  m_h_cov_x_z->Fill(errMatrix(0, 2));
194  m_h_cov_y_z->Fill(errMatrix(1, 2));
195  m_h_px->Fill(frame.getMomentum(Y4S).Px());
196  m_h_py->Fill(frame.getMomentum(Y4S).Py());
197  m_h_pz->Fill(frame.getMomentum(Y4S).Pz());
198  m_h_E->Fill(frame.getMomentum(Y4S).E());
199  m_err_y.push_back(std::sqrt(errMatrix(1, 1)));
200  m_v_y.push_back(IPVertex.Y());
201 
202  if (m_r == m_size_per_unit) {
203  m_h_temp->GetQuantiles(1, &m_median, &m_quantile);
204  for (unsigned int u = 0; u < m_v_y.size(); u++) {
205  m_h_y_risol->Fill(m_v_y.at(u) - m_median);
206  m_h_pull->Fill((m_v_y.at(u) - m_median) / m_err_y.at(u));
207  }
208  m_r = 0;
209  m_v_y.clear();
210  m_err_y.clear();
211  }
212  m_r += 1;
213  }
214 
215  }
216 
217 }
Belle2::IPDQMExpressRecoModule::m_h_pz
TH1F * m_h_pz
z coord momentum in LAB frame
Definition: IPDQMExpressRecoModule.h:94
Belle2::IPDQMExpressRecoModule::m_h_cov_y_y
TH1F * m_h_cov_y_y
Var y.
Definition: IPDQMExpressRecoModule.h:119
Belle2::IPDQMExpressRecoModule::m_err_y
std::vector< float > m_err_y
store the y errors for the pull
Definition: IPDQMExpressRecoModule.h:131
Belle2::IPDQMExpressRecoModule::m_Y4SPListName
std::string m_Y4SPListName
Name of the Y4S particle list.
Definition: IPDQMExpressRecoModule.h:140
Belle2::IPDQMExpressRecoModule::m_h_E
TH1F * m_h_E
Energy in LAB frame.
Definition: IPDQMExpressRecoModule.h:96
Belle2::IPDQMExpressRecoModule::m_h_z
TH1F * m_h_z
z coord
Definition: IPDQMExpressRecoModule.h:88
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::IPDQMExpressRecoModule::m_h_yy
TH1F * m_h_yy
yy coord
Definition: IPDQMExpressRecoModule.h:106
Belle2::IPDQMExpressRecoModule::defineHisto
void defineHisto() override
Defining the histograms.
Definition: IPDQMExpressRecoModule.cc:40
Belle2::IPDQMExpressRecoModule::m_h_y
TH1F * m_h_y
y coord
Definition: IPDQMExpressRecoModule.h:86
Belle2::IPDQMExpressRecoModule::m_h_zz
TH1F * m_h_zz
zz coord
Definition: IPDQMExpressRecoModule.h:108
Belle2::IPDQMExpressRecoModule::m_h_cov_x_z
TH1F * m_h_cov_x_z
Cov xz.
Definition: IPDQMExpressRecoModule.h:123
Belle2::IPDQMExpressRecoModule::m_quantile
Double_t m_quantile
The 0.5 quantile for the median.
Definition: IPDQMExpressRecoModule.h:134
Belle2::IPDQMExpressRecoModule::m_h_cov_x_y
TH1F * m_h_cov_x_y
Cov xy.
Definition: IPDQMExpressRecoModule.h:127
Belle2::IPDQMExpressRecoModule::event
void event() override
Event processor The main analysis happens here.
Definition: IPDQMExpressRecoModule.cc:165
Belle2::IPDQMExpressRecoModule::m_h_y_risol
TH1F * m_h_y_risol
y resolution
Definition: IPDQMExpressRecoModule.h:100
Belle2::IPDQMExpressRecoModule::m_h_cov_z_z
TH1F * m_h_cov_z_z
Var z.
Definition: IPDQMExpressRecoModule.h:121
Belle2::IPDQMExpressRecoModule::m_h_px
TH1F * m_h_px
x coord momentum in LAB frame
Definition: IPDQMExpressRecoModule.h:90
Belle2::IPDQMExpressRecoModule::m_h_pull
TH1F * m_h_pull
y pull
Definition: IPDQMExpressRecoModule.h:98
Belle2::IPDQMExpressRecoModule::endRun
void endRun() override
End-of-run action.
Definition: IPDQMExpressRecoModule.cc:133
Belle2::IPDQMExpressRecoModule::m_h_xx
TH1F * m_h_xx
xx coord
Definition: IPDQMExpressRecoModule.h:104
Belle2::IPDQMExpressRecoModule::initialize
void initialize() override
Initialize the module.
Definition: IPDQMExpressRecoModule.cc:97
Belle2::IPDQMExpressRecoModule::m_h_py
TH1F * m_h_py
y coord momentum in LAB frame
Definition: IPDQMExpressRecoModule.h:92
Belle2::IPDQMExpressRecoModule::m_h_cov_x_x
TH1F * m_h_cov_x_x
Var x.
Definition: IPDQMExpressRecoModule.h:117
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::IPDQMExpressRecoModule
This Module, made for ExpressReco, monitors the position and the dimension of the beamspot using mu+m...
Definition: IPDQMExpressRecoModule.h:38
Belle2::StoreObjPtr
Type-safe access to single objects in the data store.
Definition: ParticleList.h:33
Belle2::IPDQMExpressRecoModule::m_median
Double_t m_median
The median of y coord.
Definition: IPDQMExpressRecoModule.h:133
Belle2::IPDQMExpressRecoModule::m_h_yz
TH1F * m_h_yz
yz coord
Definition: IPDQMExpressRecoModule.h:112
Belle2::Particle
Class to store reconstructed particles.
Definition: Particle.h:77
Belle2::IPDQMExpressRecoModule::m_v_y
std::vector< float > m_v_y
store the y coordinates for the pull
Definition: IPDQMExpressRecoModule.h:129
Belle2::Particle::getVertexErrorMatrix
TMatrixFSym getVertexErrorMatrix() const
Returns the 3x3 position error sub-matrix.
Definition: Particle.cc:413
Belle2::ReferenceFrame::GetCurrent
static const ReferenceFrame & GetCurrent()
Get current rest frame.
Definition: ReferenceFrame.cc:28
Belle2::IPDQMExpressRecoModule::terminate
void terminate() override
Termination action.
Definition: IPDQMExpressRecoModule.cc:161
Belle2::IPDQMExpressRecoModule::m_h_cov_y_z
TH1F * m_h_cov_y_z
Cov yz.
Definition: IPDQMExpressRecoModule.h:125
Belle2::IPDQMExpressRecoModule::m_h_x
TH1F * m_h_x
x coord
Definition: IPDQMExpressRecoModule.h:84
Belle2::IPDQMExpressRecoModule::m_h_xz
TH1F * m_h_xz
xz coord
Definition: IPDQMExpressRecoModule.h:110
Belle2::IPDQMExpressRecoModule::m_size_per_unit
Int_t m_size_per_unit
Size for sampling per each unit.
Definition: IPDQMExpressRecoModule.h:138
Belle2::IPDQMExpressRecoModule::beginRun
void beginRun() override
Called when entering a new run Reset the histograms and counter m_r and clear the vectors.
Definition: IPDQMExpressRecoModule.cc:103
Belle2::IPDQMExpressRecoModule::m_h_xy
TH1F * m_h_xy
xy coord
Definition: IPDQMExpressRecoModule.h:114
Belle2::IPDQMExpressRecoModule::m_r
Int_t m_r
Counter for sampling.
Definition: IPDQMExpressRecoModule.h:135
Belle2::HistoModule
HistoModule.h is supposed to be used instead of Module.h for the modules with histogram definitions t...
Definition: HistoModule.h:29
Belle2::IPDQMExpressRecoModule::m_h_temp
TH1F * m_h_temp
initial histogram for median calculation
Definition: IPDQMExpressRecoModule.h:102
Belle2::StoreObjPtr::isValid
bool isValid() const
Check whether the object was created.
Definition: StoreObjPtr.h:120