Belle II Software  release-05-02-19
V0findingPerformanceEvaluationModule.cc
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2010 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Giulia Casarosa *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #include <tracking/modules/trackingPerformanceEvaluation/V0findingPerformanceEvaluationModule.h>
12 
13 #include <tracking/dataobjects/MCParticleInfo.h>
14 #include <tracking/dataobjects/V0ValidationVertex.h>
15 
16 #include <framework/datastore/StoreArray.h>
17 #include <framework/datastore/RelationVector.h>
18 
19 #include <mdst/dataobjects/Track.h>
20 
21 #include <framework/geometry/BFieldManager.h>
22 
23 #include <root/TAxis.h>
24 #include <root/TObject.h>
25 
26 #include <boost/foreach.hpp>
27 #include <vector>
28 
29 using namespace Belle2;
30 
31 //-----------------------------------------------------------------
32 // Register the Module
33 //-----------------------------------------------------------------
34 REG_MODULE(V0findingPerformanceEvaluation)
35 
36 V0findingPerformanceEvaluationModule::V0findingPerformanceEvaluationModule() :
37  Module()
38 {
39 
40  setDescription("This module evaluates the V0 finding package performance");
41 
42  addParam("outputFileName", m_rootFileName, "Name of output root file.",
43  std::string("V0findingPerformanceEvaluation_output.root"));
44  addParam("V0sName", m_V0sName, "Name of V0 collection.", std::string("V0ValidationVertexs"));
45  addParam("MCParticlesName", m_MCParticlesName, "Name of MC Particle collection.", std::string(""));
46 
47 }
48 
49 V0findingPerformanceEvaluationModule::~V0findingPerformanceEvaluationModule()
50 {
51 
52 }
53 
55 {
56  StoreArray<MCParticle> mcParticles;
57  mcParticles.isRequired(m_MCParticlesName);
58 
59  StoreArray<V0ValidationVertex> v0ValidationVertices;
60  v0ValidationVertices.isRequired(m_V0sName);
61 
62  //create list of histograms to be saved in the rootfile
63  m_histoList = new TList;
64  m_histoList_multiplicity = new TList;
65  m_histoList_efficiency = new TList;
66  m_histoList_purity = new TList;
67  m_histoList_trkQuality = new TList;
68 
69  //set the ROOT File
70  m_rootFilePtr = new TFile(m_rootFileName.c_str(), "RECREATE");
71 
72  //now create histograms
73 
74  //multiplicity histograms
75  m_multiplicityV0s = createHistogram1D("h1nV0", "number of V0s per MC Particle", 8, -0.5, 7.5, "# V0s", m_histoList_multiplicity);
76 
77  m_multiplicityMCParticles = createHistogram1D("h1nMCPrtcl", "number of MCParticles per V0s", 5, -0.5, 4.5,
78  "# MCParticles", m_histoList_multiplicity);
79 
80  m_MCParticlePDGcode = createHistogram1D("h1PDGcode", "PDG code of MCParticles", 6244, -3122, 3122,
81  "PDG code", m_histoList_multiplicity);
82 
83 
84  //vertex and momentum parameters errors
85  m_h1_vtxX_err = createHistogram1D("h1vtxXerr", "vtxX error", 100, 0, 0.1, "#sigma_{vtxX} (cm)", m_histoList);
86  m_h1_vtxY_err = createHistogram1D("h1vtxYerr", "vtxY error", 100, 0, 0.1, "#sigma_{vtxY} (cm)", m_histoList);
87  m_h1_vtxZ_err = createHistogram1D("h1vtxZerr", "vtxZ error", 100, 0, 0.3, "#sigma_{vtxZ} (cm)", m_histoList);
88  m_h2_vtxTvsR_err = createHistogram2D("h2vtxTerrVsR", "vtxT error vs R", 100, 0, 100, "R (cm)", 100, 0, 0.3, "#sigma_{vtxT} (cm)",
89  m_histoList);
90  // m_h1_mom_err = createHistogram1D("h1momerr", "mom error", 100, 0, 0.1, "#sigma_{p} (GeV/c)", m_histoList);
91  // m_h1_mass_err = createHistogram1D("h1masserr", "mass error", 100, 0, 1, "#sigma_{m} (GeV/c2)", m_histoList);
92  //vertex and momentum parameters residuals
93  m_h1_vtxX_res = createHistogram1D("h1vtxXres", "vtxX resid", 100, -0.2, 0.2, "vtxX resid (cm)", m_histoList);
94  m_h1_vtxY_res = createHistogram1D("h1vtxYres", "vtxY resid", 100, -0.2, 0.2, "vtxY resid (cm)", m_histoList);
95  m_h1_vtxZ_res = createHistogram1D("h1vtxZres", "vtxZ resid", 100, -0.5, 0.5, "vtxZ resid (cm)", m_histoList);
96  m_h1_mom_res = createHistogram1D("h1momres", "mom resid", 1000, -0.5, 0.5, "mom resid (GeV/c)", m_histoList);
97  m_h1_mass_res = createHistogram1D("h1massres", "mass resid", 500, -0.3, 0.3, "mass resid (GeV/c)", m_histoList);
98 
99  //vertex and momentum parameters pulls
100  m_h1_vtxX_pll = createHistogram1D("h1vtxXpll", "vtxX pull", 100, -5, 5, "vtxX pull", m_histoList);
101  m_h1_vtxY_pll = createHistogram1D("h1vtxYpll", "vtxY pull", 100, -5, 5, "vtxY pull", m_histoList);
102  m_h1_vtxZ_pll = createHistogram1D("h1vtxZpll", "vtxZ pull", 100, -5, 5, "vtxZ pull", m_histoList);
103  // m_h1_mom_pll = createHistogram1D("h1mompll", "mom pull", 100, -5, 5, "momX pull", m_histoList);
104  // m_h1_mass_pll = createHistogram1D("h1masspll", "mass pull", 100, -5, 5, "momY pull", m_histoList);
105 
106 
107  m_h1_ChiSquare = createHistogram1D("h1Chi2", "Chi2 of the fit", 100, 0, 20, "Chi2", m_histoList_trkQuality);
108 
109  m_h1_nMatchedDau = createHistogram1D("h1nMatchedDau", "Number of Matched MCParticle Daughters", 3, -0.5, 2.5, "# matched dau",
110  m_histoList);
111 
112 
113  m_h2_mom = createHistogram2D("h2mom", "reco VS true momentum", 100, 0, 3, "V0 mom (GeV/c)", 100, 0, 3, "MC mom (GeV/c)",
114  m_histoList);
115  m_h2_mass = createHistogram2D("h2mass", "reco VS true mass", 100, 0, 1.5, "V0 mass (GeV/c2)", 100, 0, 1.5, "MC mass (GeV/c)",
116  m_histoList);
117 
118  //histograms to produce efficiency plots
119  Double_t bins_pt[9 + 1] = {0, 0.05, 0.1, 0.15, 0.2, 0.3, 0.5, 1, 2, 3.5}; //GeV/c
120  Double_t bins_theta[10 + 1] = {0, 0.25, 0.5, 0.75, 0.75 + 0.32, 0.75 + 2 * 0.32, 0.75 + 3 * 0.32, 0.75 + 4 * 0.32, 0.75 + 5 * 0.32, 2.65, TMath::Pi()};
121  Double_t bins_phi[14 + 1];
122  Double_t width_phi = 2 * TMath::Pi() / 14;
123  for (int bin = 0; bin < 14 + 1; bin++)
124  bins_phi[bin] = - TMath::Pi() + bin * width_phi;
125 
126  m_h1_MCParticle_R = createHistogram1D("h1nMCParticleVSr", "entry per MCParticles", 50, 0, 20, "transverse L", m_histoList);
127 
128  m_h1_V0sPerMCParticle_R = (TH1F*)duplicateHistogram("h1nV0perMCvsR", "entry per V0 related to a MCParticle", m_h1_MCParticle_R,
129  m_histoList);
130 
131 
132  m_h3_MCParticle = createHistogram3D("h3MCParticle", "entry per MCParticle",
133  9, bins_pt, "p_{t} (GeV/c)",
134  10, bins_theta, "#theta",
135  14, bins_phi, "#phi" /*, m_histoList*/);
136 
137  m_h3_V0sPerMCParticle = (TH3F*)duplicateHistogram("h3V0sPerMCParticle",
138  "entry per V0 connected to a MCParticle",
139  m_h3_MCParticle /*, m_histoList*/);
140 
141  m_h3_V0s = (TH3F*)duplicateHistogram("h3V0s", "entry per V0",
142  m_h3_MCParticle /*, m_histoList*/);
143 
144  //histograms to produce purity plots
145  m_h3_MCParticlesPerV0 = (TH3F*)duplicateHistogram("h3MCParticlesPerV0",
146  "entry per MCParticle connected to a V0",
147  m_h3_MCParticle /*, m_histoList*/);
148 }
149 
151 {
152 
153 }
154 
156 {
157 
159 
160  B2Vector3D magField = BFieldManager::getField(0, 0, 0) / Unit::T;
161 
162  B2DEBUG(99, "+++++ 1. loop on MCParticles");
163  BOOST_FOREACH(MCParticle & mcParticle, mcParticles) {
164 
165  if (! isV0(mcParticle))
166  continue;
167 
168  int nMatchedDau = nMatchedDaughters(mcParticle);
169  m_h1_nMatchedDau->Fill(nMatchedDau);
170 
171  //proceed only in case the MCParticle daughters have both one associated reconstructed track:
172  if (nMatchedDau != 2)
173  continue;
174 
175  int pdgCode = mcParticle.getPDG();
176  B2DEBUG(99, "MCParticle has PDG code " << pdgCode);
177  m_MCParticlePDGcode->Fill(mcParticle.getPDG());
178 
179  MCParticleInfo mcParticleInfo(mcParticle, magField);
180 
181  TVector3 MC_prodvtx = mcParticle.getVertex();
182  TVector3 MC_vtx = mcParticle.getDecayVertex();
183  float MC_mom = mcParticle.getMomentum().Mag();
184  float MC_mass = mcParticle.getMass();
185  TVector3 MC_FL = MC_vtx - MC_prodvtx;
186  float flightR = sqrt(MC_FL.X() * MC_FL.X() + MC_FL.Y() * MC_FL.Y());
187  m_h1_MCParticle_R->Fill(flightR);
188 
189  m_h3_MCParticle->Fill(mcParticleInfo.getPt(), mcParticleInfo.getPtheta(), mcParticleInfo.getPphi());
190 
191  //1. retrieve all the V0s related to the MCParticle
192 
193  //1.0 check if there is a V0
194  RelationVector<V0ValidationVertex> V0s_toMCParticle =
195  DataStore::getRelationsWithObj<V0ValidationVertex>(&mcParticle, m_V0sName);
196 
197  m_multiplicityV0s->Fill(V0s_toMCParticle.size());
198 
199  if (V0s_toMCParticle.size() > 0)
200  m_h3_V0sPerMCParticle->Fill(mcParticleInfo.getPt(), mcParticleInfo.getPtheta(), mcParticleInfo.getPphi());
201 
202  for (int v0 = 0; v0 < (int)V0s_toMCParticle.size(); v0++) {
203 
204  TVector3 V0_vtx = V0s_toMCParticle[v0]->getVertexPosition();
205  float V0_mom = V0s_toMCParticle[v0]->getFittedMomentum();
206  float V0_chi2 = V0s_toMCParticle[v0]->getVertexChi2();
207  float V0_mass = V0s_toMCParticle[v0]->getFittedInvariantMass();
208  TMatrixDSym V0_cov = V0s_toMCParticle[v0]->getVertexPositionCovariance();
209 
210  m_h1_vtxX_err->Fill(sqrt(V0_cov[0][0]));
211  m_h1_vtxY_err->Fill(sqrt(V0_cov[1][1]));
212  m_h1_vtxZ_err->Fill(sqrt(V0_cov[2][2]));
213  m_h2_vtxTvsR_err->Fill(flightR, sqrt(V0_cov[0][0] + V0_cov[1][1]));
214 
215  m_h1_V0sPerMCParticle_R->Fill(flightR);
216 
217  m_h1_vtxX_res->Fill(V0_vtx.X() - MC_vtx.X());
218  m_h1_vtxY_res->Fill(V0_vtx.Y() - MC_vtx.Y());
219  m_h1_vtxZ_res->Fill(V0_vtx.Z() - MC_vtx.Z());
220 
221  m_h1_mom_res->Fill(V0_mom - MC_mom);
222  m_h2_mom->Fill(V0_mom, MC_mom);
223  m_h1_mass_res->Fill(V0_mass - MC_mass);
224  m_h2_mass->Fill(V0_mass, MC_mass);
225 
226  m_h1_vtxX_pll->Fill((V0_vtx.X() - MC_vtx.X()) / sqrt(V0_cov[0][0]));
227  m_h1_vtxY_pll->Fill((V0_vtx.Y() - MC_vtx.Y()) / sqrt(V0_cov[1][1]));
228  m_h1_vtxZ_pll->Fill((V0_vtx.Z() - MC_vtx.Z()) / sqrt(V0_cov[2][2]));
229 
230  m_h1_ChiSquare->Fill(V0_chi2);
231 
232  }
233 
234 
235  }
236 
237 
238  B2DEBUG(99, "+++++ 2. loop on V0s");
239 
240 
242 
243 
244  BOOST_FOREACH(V0ValidationVertex & v0, V0s) {
245 
246  int nMCParticles = 0;
247 
248  // m_h3_V0s>Fill(mcParticleInfo.getPt(), mcParticleInfo.getPtheta(), mcParticleInfo.getPphi());
249  //2. retrieve all the MCParticles related to the V0s
250  RelationVector<MCParticle> MCParticles_fromV0 =
251  DataStore::getRelationsWithObj<MCParticle>(&v0, m_MCParticlesName);
252 
253  nMCParticles = MCParticles_fromV0.size();
254 
255  if (nMCParticles == 0)
256  continue;
257 
258  MCParticleInfo mcParticleInfo(* MCParticles_fromV0[0], magField);
259  m_h3_MCParticlesPerV0->Fill(mcParticleInfo.getPt(), mcParticleInfo.getPtheta(), mcParticleInfo.getPphi());
260  m_multiplicityMCParticles->Fill(nMCParticles);
261 
262  }
263 
264 }
265 
267 {
268 
269  double num = 0;
270  double den = 0;
271 
272  for (int bin = 1; bin < m_multiplicityV0s->GetNbinsX(); bin ++)
273  num += m_multiplicityV0s->GetBinContent(bin + 1);
274  den = m_multiplicityV0s->GetEntries();
275  double efficiency = num / den ;
276  double efficiencyErr = sqrt(efficiency * (1 - efficiency)) / sqrt(den);
277 
278  double nMCParticles = 0;
279  for (int bin = 1; bin < m_multiplicityMCParticles->GetNbinsX(); bin ++)
280  nMCParticles += m_multiplicityMCParticles->GetBinContent(bin + 1);
281  double purity = nMCParticles / m_multiplicityMCParticles->GetEntries();
282  double purityErr = sqrt(purity * (1 - purity)) / sqrt(m_multiplicityMCParticles->GetEntries());
283 
284  B2INFO("");
285  B2INFO("~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~");
286  B2INFO("~ V0 Finding Performance Evaluation ~ SHORT SUMMARY ~");
287  B2INFO("");
288  B2INFO(" + overall:");
289  B2INFO(" efficiency = (" << efficiency * 100 << " +/- " << efficiencyErr * 100 << ")% ");
290  B2INFO(" purity = (" << purity * 100 << " +/- " << purityErr * 100 << ")% ");
291  B2INFO("");
292  B2INFO("~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~");
293 }
294 
296 {
297 
298  TH1F* h_eff_R = (TH1F*)duplicateHistogram("h_eff_R", "efficiency vs R", m_h1_MCParticle_R, m_histoList_efficiency);
299 
300  for (int bin = 0; bin < h_eff_R->GetXaxis()->GetNbins(); bin++) {
301  float num = m_h1_V0sPerMCParticle_R->GetBinContent(bin + 1);
302  float den = m_h1_MCParticle_R->GetBinContent(bin + 1);
303  double eff = 0;
304  double err = 0;
305 
306  if (den > 0) {
307  eff = (double)num / den;
308  err = sqrt(eff * (1 - eff)) / sqrt(den);
309  }
310 
311  h_eff_R->SetBinContent(bin + 1, eff);
312  h_eff_R->SetBinError(bin + 1, err);
313  }
314 
316 
318 
319  // addPurityPlots(m_histoList_purity, m_h3_MCParticlesPerV0, m_h3_V0s);
320 
321  if (m_rootFilePtr != NULL) {
322  m_rootFilePtr->cd();
323 
324  TDirectory* oldDir = gDirectory;
325 
326  TDirectory* dir_multiplicity = oldDir->mkdir("multiplicity");
327  dir_multiplicity->cd();
328  TIter nextH_multiplicity(m_histoList_multiplicity);
329  TObject* obj;
330  while ((obj = nextH_multiplicity()))
331  obj->Write();
332 
333  TDirectory* dir_efficiency = oldDir->mkdir("efficiency");
334  dir_efficiency->cd();
335  TIter nextH_efficiency(m_histoList_efficiency);
336  while ((obj = nextH_efficiency()))
337  obj->Write();
338 
339  TDirectory* dir_purity = oldDir->mkdir("purity");
340  dir_purity->cd();
341  TIter nextH_purity(m_histoList_purity);
342  while ((obj = nextH_purity()))
343  obj->Write();
344 
345  TDirectory* dir_trkQuality = oldDir->mkdir("trkQuality");
346  dir_trkQuality->cd();
347  TIter nextH_trkQuality(m_histoList_trkQuality);
348  while ((obj = nextH_trkQuality()))
349  obj->Write();
350 
351  TIter nextH(m_histoList);
352  while ((obj = nextH()))
353  obj->Write();
354 
355 
356  m_rootFilePtr->Close();
357  }
358 
359 }
360 
361 
363 {
364 
365  bool isGamma = false;
366  if (abs(the_mcParticle.getPDG()) == 22)
367  isGamma = true;
368 
369  bool isK_S0 = false;
370  if (abs(the_mcParticle.getPDG()) == 310)
371  isK_S0 = true;
372 
373  bool isK_0 = false;
374  if (abs(the_mcParticle.getPDG()) == 311)
375  isK_0 = true;
376 
377  bool isLambda = false;
378  if (abs(the_mcParticle.getPDG()) == 3122)
379  isLambda = true;
380 
381  bool twoProngs = false;
382  bool twoChargedProngs = false;
383 
384  if (the_mcParticle.getDaughters().size() == 2)
385  twoProngs = true;
386 
387  if (twoProngs)
388  if (the_mcParticle.getDaughters()[0]->getCharge() * the_mcParticle.getDaughters()[1]->getCharge() < 0)
389  twoChargedProngs = true;
390 
391  return ((isGamma || isK_S0 || isK_0 || isLambda) && twoChargedProngs);
392 
393 }
394 
396 {
397 
398  int nMatchedDau = 0;
399 
400  std::vector< MCParticle* > MCPart_dau = the_mcParticle.getDaughters();
401 
402  bool first = false;
403  bool second = false;
404 
405  RelationVector<Track> Tracks_fromMCParticle_0 = DataStore::getRelationsWithObj<Track>(MCPart_dau[0]);
406  if (Tracks_fromMCParticle_0.size() > 0)
407  first = true;
408 
409  RelationVector<Track> Tracks_fromMCParticle_1 = DataStore::getRelationsWithObj<Track>(MCPart_dau[1]);
410  if (Tracks_fromMCParticle_1.size() > 0)
411  second = true;
412 
413 
414  if (first)
415  nMatchedDau++;
416 
417  if (second)
418  nMatchedDau++;
419 
420 
421  return nMatchedDau;
422 
423 }
Belle2::RelationVector::size
size_t size() const
Get number of relations.
Definition: RelationVector.h:98
Belle2::MCParticleInfo
This struct is used by the TrackingPerformanceEvaluation Module to save information of reconstructed ...
Definition: MCParticleInfo.h:34
Belle2::V0findingPerformanceEvaluationModule::m_h1_ChiSquare
TH1F * m_h1_ChiSquare
TH1F chi square.
Definition: V0findingPerformanceEvaluationModule.h:100
Belle2::Module::setDescription
void setDescription(const std::string &description)
Sets the description of the module.
Definition: Module.cc:216
Belle2::BFieldManager::getField
static void getField(const double *pos, double *field)
return the magnetic field at a given position.
Definition: BFieldManager.h:110
Belle2::V0findingPerformanceEvaluationModule::endRun
void endRun() override
This method is called if the current run ends.
Definition: V0findingPerformanceEvaluationModule.cc:266
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::V0findingPerformanceEvaluationModule::m_h3_MCParticlesPerV0
TH3F * m_h3_MCParticlesPerV0
V0-finding numerator.
Definition: V0findingPerformanceEvaluationModule.h:112
Belle2::V0findingPerformanceEvaluationModule::m_MCParticlePDGcode
TH1F * m_MCParticlePDGcode
MCParticle PDG code.
Definition: V0findingPerformanceEvaluationModule.h:76
Belle2::V0findingPerformanceEvaluationModule::m_h1_vtxX_pll
TH1F * m_h1_vtxX_pll
vtx pull
Definition: V0findingPerformanceEvaluationModule.h:94
Belle2::V0findingPerformanceEvaluationModule::m_h1_mom_res
TH1F * m_h1_mom_res
mom resid
Definition: V0findingPerformanceEvaluationModule.h:89
Belle2::PerformanceEvaluationBaseClass::createHistogram2D
TH2F * createHistogram2D(const char *name, const char *title, Int_t nbinsX, Double_t minX, Double_t maxX, const char *titleX, Int_t nbinsY, Double_t minY, Double_t maxY, const char *titleY, TList *histoList=NULL)
Create a 2D histogram and add it to the TList of 2D-histograms.
Definition: PerformanceEvaluationBaseClass.cc:58
Belle2::V0findingPerformanceEvaluationModule::nMatchedDaughters
int nMatchedDaughters(const MCParticle &the_mcParticle)
number of truth matched dauhters
Definition: V0findingPerformanceEvaluationModule.cc:395
Belle2::PerformanceEvaluationBaseClass::m_histoList_efficiency
TList * m_histoList_efficiency
List of efficiency histograms.
Definition: PerformanceEvaluationBaseClass.h:48
Belle2::V0ValidationVertex
Class which stores some additional information on V0 vertices.
Definition: V0ValidationVertex.h:18
Belle2::PerformanceEvaluationBaseClass::m_histoList_multiplicity
TList * m_histoList_multiplicity
List of multiplicity histograms.
Definition: PerformanceEvaluationBaseClass.h:42
Belle2::MCParticleInfo::getPtheta
double getPtheta()
Getter for theta of momentum vector.
Definition: MCParticleInfo.h:62
Belle2::PerformanceEvaluationBaseClass::addInefficiencyPlots
void addInefficiencyPlots(TList *graphList=NULL, TH3F *h3_xPerMCParticle=NULL, TH3F *h3_MCParticle=NULL)
Create pt-, theta- and phi-inefficiency 1D histograms and add them to the TList of 1D-histograms.
Definition: PerformanceEvaluationBaseClass.cc:271
Belle2::V0findingPerformanceEvaluationModule::m_h1_vtxX_res
TH1F * m_h1_vtxX_res
vtx resid
Definition: V0findingPerformanceEvaluationModule.h:86
Belle2::V0findingPerformanceEvaluationModule::beginRun
void beginRun() override
Called when entering a new run.
Definition: V0findingPerformanceEvaluationModule.cc:150
Belle2::V0findingPerformanceEvaluationModule::m_h1_vtxY_pll
TH1F * m_h1_vtxY_pll
vtx pull
Definition: V0findingPerformanceEvaluationModule.h:95
Belle2::V0findingPerformanceEvaluationModule::m_multiplicityV0s
TH1F * m_multiplicityV0s
number of V0s per MCParticles
Definition: V0findingPerformanceEvaluationModule.h:74
Belle2::V0findingPerformanceEvaluationModule::m_h1_mass_res
TH1F * m_h1_mass_res
mom resid
Definition: V0findingPerformanceEvaluationModule.h:91
Belle2::B2Vector3< double >
Belle2::PerformanceEvaluationBaseClass::addEfficiencyPlots
void addEfficiencyPlots(TList *graphList=NULL, TH3F *h3_xPerMCParticle=NULL, TH3F *h3_MCParticle=NULL)
Create pt-, theta- and phi-efficiency 1D histograms and add them to the TList of 1D-histograms.
Definition: PerformanceEvaluationBaseClass.cc:292
Belle2::MCParticle::getDecayVertex
TVector3 getDecayVertex() const
Return decay vertex.
Definition: MCParticle.h:230
Belle2::Module
Base class for Modules.
Definition: Module.h:74
Belle2::V0findingPerformanceEvaluationModule::m_h1_MCParticle_R
TH1F * m_h1_MCParticle_R
V0-finding denominator by radius.
Definition: V0findingPerformanceEvaluationModule.h:107
Belle2::Unit::T
static const double T
[tesla]
Definition: Unit.h:130
Belle2::MCParticleInfo::getPt
double getPt()
Getter for transverse momentum.
Definition: MCParticleInfo.h:48
Belle2::V0findingPerformanceEvaluationModule::m_h3_V0sPerMCParticle
TH3F * m_h3_V0sPerMCParticle
V0-finding numerator.
Definition: V0findingPerformanceEvaluationModule.h:106
Belle2::RelationVector
Class for type safe access to objects that are referred to in relations.
Definition: DataStore.h:38
Belle2::PerformanceEvaluationBaseClass::m_rootFileName
std::string m_rootFileName
root file name
Definition: PerformanceEvaluationBaseClass.h:139
Belle2::PerformanceEvaluationBaseClass::m_histoList_trkQuality
TList * m_histoList_trkQuality
List of track-quality histograms.
Definition: PerformanceEvaluationBaseClass.h:44
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::V0findingPerformanceEvaluationModule::initialize
void initialize() override
Initialize the Module.
Definition: V0findingPerformanceEvaluationModule.cc:54
Belle2::TrackingPerformanceEvaluationModule::m_MCParticlesName
std::string m_MCParticlesName
MCParticle StoreArray name.
Definition: TrackingPerformanceEvaluationModule.h:82
Belle2::MCParticle::getVertex
TVector3 getVertex() const
Return production vertex position, shorthand for getProductionVertex().
Definition: MCParticle.h:194
Belle2::MCParticleInfo::getPphi
double getPphi()
Getter for phi of momentum vector.
Definition: MCParticleInfo.h:64
Belle2::V0findingPerformanceEvaluationModule::terminate
void terminate() override
This method is called at the end of the event processing.
Definition: V0findingPerformanceEvaluationModule.cc:295
Belle2::V0findingPerformanceEvaluationModule::m_h1_vtxY_err
TH1F * m_h1_vtxY_err
vtx error
Definition: V0findingPerformanceEvaluationModule.h:80
Belle2::PerformanceEvaluationBaseClass::createHistogram3D
TH3F * createHistogram3D(const char *name, const char *title, Int_t nbinsX, Double_t minX, Double_t maxX, const char *titleX, Int_t nbinsY, Double_t minY, Double_t maxY, const char *titleY, Int_t nbinsZ, Double_t minZ, Double_t maxZ, const char *titleZ, TList *histoList=NULL)
Create a 3D histogram and add it to the TList of 3D-histograms.
Definition: PerformanceEvaluationBaseClass.cc:95
Belle2::MCParticle::getPDG
int getPDG() const
Return PDG code of particle.
Definition: MCParticle.h:123
Belle2::V0findingPerformanceEvaluationModule::m_h1_vtxX_err
TH1F * m_h1_vtxX_err
vtx error
Definition: V0findingPerformanceEvaluationModule.h:79
Belle2::V0findingPerformanceEvaluationModule::m_h1_vtxY_res
TH1F * m_h1_vtxY_res
vtx resid
Definition: V0findingPerformanceEvaluationModule.h:87
Belle2::V0findingPerformanceEvaluationModule::m_h2_mass
TH2F * m_h2_mass
mass reco VS true
Definition: V0findingPerformanceEvaluationModule.h:92
Belle2::MCParticle::getMass
float getMass() const
Return the particle mass in GeV.
Definition: MCParticle.h:146
Belle2::V0findingPerformanceEvaluationModule::m_h1_vtxZ_pll
TH1F * m_h1_vtxZ_pll
vtx pull
Definition: V0findingPerformanceEvaluationModule.h:96
Belle2::PerformanceEvaluationBaseClass::m_rootFilePtr
TFile * m_rootFilePtr
pointer at root file used for storing histograms
Definition: PerformanceEvaluationBaseClass.h:142
Belle2::PerformanceEvaluationBaseClass::m_histoList
TList * m_histoList
List of performance-evaluation histograms.
Definition: PerformanceEvaluationBaseClass.h:40
Belle2::V0findingPerformanceEvaluationModule::m_h3_V0s
TH3F * m_h3_V0s
V0-finding purity denominator.
Definition: V0findingPerformanceEvaluationModule.h:111
Belle2::V0findingPerformanceEvaluationModule::m_h1_V0sPerMCParticle_R
TH1F * m_h1_V0sPerMCParticle_R
V0-finding numerator by radius.
Definition: V0findingPerformanceEvaluationModule.h:108
Belle2::MCParticle::getDaughters
std::vector< Belle2::MCParticle * > getDaughters() const
Get vector of all daughter particles, empty vector if none.
Definition: MCParticle.cc:51
Belle2::V0findingPerformanceEvaluationModule::m_MCParticlesName
std::string m_MCParticlesName
MCParticle StoreArray name.
Definition: V0findingPerformanceEvaluationModule.h:69
Belle2::V0findingPerformanceEvaluationModule::m_h2_mom
TH2F * m_h2_mom
mom reco VS true
Definition: V0findingPerformanceEvaluationModule.h:90
Belle2::Module::addParam
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:562
Belle2::PerformanceEvaluationBaseClass::createHistogram1D
TH1F * createHistogram1D(const char *name, const char *title, Int_t nbins, Double_t min, Double_t max, const char *xtitle, TList *histoList=NULL)
Create a 1D histogram and add it to the TList of 1D-histograms.
Definition: PerformanceEvaluationBaseClass.cc:28
Belle2::MCParticle::getMomentum
TVector3 getMomentum() const
Return momentum.
Definition: MCParticle.h:209
Belle2::V0findingPerformanceEvaluationModule::m_h2_vtxTvsR_err
TH2F * m_h2_vtxTvsR_err
vtx error on transverse plane VS transverse flight length
Definition: V0findingPerformanceEvaluationModule.h:82
Belle2::V0findingPerformanceEvaluationModule::event
void event() override
This method is the core of the module.
Definition: V0findingPerformanceEvaluationModule.cc:155
Belle2::V0findingPerformanceEvaluationModule::m_h1_vtxZ_res
TH1F * m_h1_vtxZ_res
vtx resid
Definition: V0findingPerformanceEvaluationModule.h:88
Belle2::MCParticle
A Class to store the Monte Carlo particle information.
Definition: MCParticle.h:43
Belle2::StoreArray< MCParticle >
Belle2::PerformanceEvaluationBaseClass::duplicateHistogram
TH1 * duplicateHistogram(const char *newname, const char *newtitle, TH1 *h, TList *histoList=NULL)
Make a copy of a 1D histogram and add it to the TList of 1D-histograms.
Definition: PerformanceEvaluationBaseClass.cc:139
Belle2::V0findingPerformanceEvaluationModule::m_V0sName
std::string m_V0sName
MCTrackCand StoreArray name.
Definition: V0findingPerformanceEvaluationModule.h:70
Belle2::V0findingPerformanceEvaluationModule::m_multiplicityMCParticles
TH1F * m_multiplicityMCParticles
number of MCParticles per fitted Track
Definition: V0findingPerformanceEvaluationModule.h:75
Belle2::V0findingPerformanceEvaluationModule::m_h1_vtxZ_err
TH1F * m_h1_vtxZ_err
vtx error
Definition: V0findingPerformanceEvaluationModule.h:81
Belle2::V0findingPerformanceEvaluationModule::m_h3_MCParticle
TH3F * m_h3_MCParticle
V0-finding denominator.
Definition: V0findingPerformanceEvaluationModule.h:105
Belle2::V0findingPerformanceEvaluationModule::m_h1_nMatchedDau
TH1F * m_h1_nMatchedDau
TH1F n matched daughters.
Definition: V0findingPerformanceEvaluationModule.h:102
Belle2::V0findingPerformanceEvaluationModule::isV0
bool isV0(const MCParticle &the_mcParticle)
is V0
Definition: V0findingPerformanceEvaluationModule.cc:362
Belle2::PerformanceEvaluationBaseClass::m_histoList_purity
TList * m_histoList_purity
List of purity histograms.
Definition: PerformanceEvaluationBaseClass.h:49