Belle II Software  release-05-02-19
TrackingPerformanceEvaluationModule.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/TrackingPerformanceEvaluationModule.h>
12 
13 #include <framework/datastore/StoreArray.h>
14 #include <framework/datastore/RelationVector.h>
15 
16 #include <framework/geometry/BFieldManager.h>
17 
18 #include <vxd/geometry/GeoCache.h>
19 
20 #include <mdst/dataobjects/MCParticle.h>
21 #include <mdst/dataobjects/Track.h>
22 #include <mdst/dataobjects/HitPatternCDC.h>
23 #include <mdst/dataobjects/HitPatternVXD.h>
24 
25 #include <pxd/reconstruction/PXDRecoHit.h>
26 #include <svd/reconstruction/SVDRecoHit.h>
27 #include <svd/reconstruction/SVDRecoHit2D.h>
28 #include <cdc/dataobjects/CDCRecoHit.h>
29 
30 #include <tracking/dataobjects/RecoTrack.h>
31 
32 #include <pxd/dataobjects/PXDTrueHit.h>
33 #include <pxd/dataobjects/PXDCluster.h>
34 #include <svd/dataobjects/SVDCluster.h>
35 #include <cdc/dataobjects/CDCHit.h>
36 
37 #include <genfit/KalmanFitterInfo.h>
38 
39 #include <root/TObject.h>
40 
41 #include <boost/foreach.hpp>
42 
43 using namespace Belle2;
44 
45 //-----------------------------------------------------------------
46 // Register the Module
47 //-----------------------------------------------------------------
48 REG_MODULE(TrackingPerformanceEvaluation)
49 
50 TrackingPerformanceEvaluationModule::TrackingPerformanceEvaluationModule() :
51  Module()
52 {
53 
54  setDescription("This module evaluates the tracking package performance");
55 
56  addParam("outputFileName", m_rootFileName, "Name of output root file.",
57  std::string("TrackingPerformanceEvaluation_output.root"));
58  addParam("MCParticlesName", m_MCParticlesName, "Name of MC Particle collection.", std::string(""));
59  addParam("TracksName", m_TracksName, "Name of Track collection.", std::string(""));
60  addParam("RecoTracksName", m_RecoTracksName, "Name of RecoTrack collection.", std::string("RecoTracks"));
61  addParam("MCRecoTracksName", m_MCRecoTracksName, "Name of MCRecoTrack collection.", std::string("MCRecoTracks"));
62  addParam("ParticleHypothesis", m_ParticleHypothesis, "Particle Hypothesis used in the track fit.", int(211));
63 
64 }
65 
66 TrackingPerformanceEvaluationModule::~TrackingPerformanceEvaluationModule()
67 {
68 
69 }
70 
72 {
73  // MCParticles, Tracks, RecoTracks, MCRecoTracks needed for this module
74  StoreArray<MCParticle> mcParticles;
75  mcParticles.isRequired(m_MCParticlesName);
76 
77  StoreArray<RecoTrack> recoTracks;
78  recoTracks.isRequired(m_RecoTracksName);
79 
80  StoreArray<RecoTrack> mcRecoTracks;
81  mcRecoTracks.isRequired(m_MCRecoTracksName);
82 
83  StoreArray<Track> tracks;
84  tracks.isRequired(m_TracksName);
85 
86  //create list of histograms to be saved in the rootfile
87  m_histoList = new TList;
88  m_histoList_multiplicity = new TList;
90  m_histoList_trkQuality = new TList;
91  m_histoList_firstHit = new TList;
92  m_histoList_pr = new TList;
93  m_histoList_fit = new TList;
94  m_histoList_efficiency = new TList;
95  m_histoList_purity = new TList;
96  m_histoList_others = new TList;
97 
98  //set the ROOT File
99  m_rootFilePtr = new TFile(m_rootFileName.c_str(), "RECREATE");
100 
101  //now create histograms
102 
103  //multiplicity histograms
104  m_multiplicityTracks = createHistogram1D("h1nTrk", "number of tracks per MC Particle", 8, -0.5, 7.5, "# tracks",
106 
107  m_multiplicityRecoTracks = createHistogram1D("h1nRcTrk", "number of recoTracks per MC Particle", 8, -0.5, 7.5, "# tracks",
109 
110  m_multiplicityMCRecoTracks = createHistogram1D("h1nMCRcTrk", "number of MC recoTracks per MC Particle", 8, -0.5, 7.5, "# tracks",
112 
113  m_multiplicityFittedTracks = createHistogram1D("h1nFitTrk", "number of fitted tracks per MC Particle", 5, -0.5, 4.5,
114  "# fitted tracks", m_histoList_multiplicity);
115  m_multiplicityFittedTracksPerMCRT = createHistogram1D("h1nFitTrkMCRT", "number of fitted tracks per MCRecoTrack", 5, -0.5, 4.5,
116  "# fitted tracks", m_histoList_multiplicity);
117  m_multiplicityMCParticlesPerTrack = createHistogram1D("h1nMCPrtcl", "number of MCParticles per fitted tracks", 5, -0.5, 4.5,
118  "# MCParticles", m_histoList_multiplicity);
119  m_multiplicityRecoTracksPerMCRT = createHistogram1D("h1nRecoTrack", "number of RecoTrack per MCRecoTrack", 5, -0.5, 4.5,
120  "# RecoTrack",
122  m_multiplicityMCRecoTracksPerRT = createHistogram1D("h1nMCRecoTrack", "number of MCRecoTrack per RecoTrack", 5, -0.5, 4.5,
123  "# MCRecoTrack", m_histoList_multiplicity);
124 
125  //tracks pValue
126  m_h1_pValue = createHistogram1D("h1pValue", "pValue of the fit", 100, 0, 1, "pValue", m_histoList_trkQuality);
127 
128 
129  //track parameters errors
130  m_h1_d0_err = createHistogram1D("h1d0err", "d0 error", 100, 0, 0.1, "#sigma_{d0} (cm)", m_histoList_trkQuality);
131  m_h1_phi_err = createHistogram1D("h1phierr", "#phi error", 100, 0, 0.02, "#sigma_{#phi} (rad)", m_histoList_trkQuality);
132  m_h1_omega_err = createHistogram1D("h1omegaerr", "#omega error", 100, 0, 0.002, "#sigma_{#omega} (cm^{-1})",
134  m_h1_z0_err = createHistogram1D("h1z0err", "z0 error", 100, 0, 0.1, "#sigma_{z0} (cm)", m_histoList_trkQuality);
135  m_h1_cotTheta_err = createHistogram1D("h1cotThetaerr", "cot#theta error", 100, 0, 0.03, "#sigma_{cot#theta}",
137 
138  //track parameters residuals
139  m_h1_d0_res = createHistogram1D("h1d0res", "d0 residuals", 100, -0.1, 0.1, "d0 resid (cm)", m_histoList_trkQuality);
140  m_h1_phi_res = createHistogram1D("h1phires", "#phi residuals", 100, -0.1, 0.1, "#phi resid (rad)", m_histoList_trkQuality);
141  m_h1_omega_res = createHistogram1D("h1omegares", "#omega residuals", 100, -0.0005, 0.0005, "#omega resid (cm^{-1})",
143  m_h1_z0_res = createHistogram1D("h1z0res", "z0 residuals", 100, -0.1, 0.1, "z0 resid (cm)", m_histoList_trkQuality);
144  m_h1_cotTheta_res = createHistogram1D("h1cotThetares", "cot#theta residuals", 100, -0.1, 0.1, "cot#theta resid",
146 
147  //track parameters residuals - momentum
148  m_h1_px_res = createHistogram1D("h1pxres", "px absolute residuals", 200, -0.05, 0.05, "px_{reco} - px_{MC} (GeV/c)",
150  m_h1_py_res = createHistogram1D("h1pyres", "py absolute residuals", 200, -0.05, 0.05, "py_{reco} - py_{MC} (GeV/c)",
152  m_h1_pz_res = createHistogram1D("h1pzres", "pz absolute residuals", 200, -0.05, 0.05, "pz_{reco} - pz_{MC} (GeV/c)",
154  m_h1_p_res = createHistogram1D("h1pres", "p relative residuals", 200, -0.05, 0.05, "p_{reco} - p_{MC} / p_{MC}",
156  m_h1_pt_res = createHistogram1D("h1ptres", "pt relative residuals", 200, -0.05, 0.05, "pt_{reco} - pt_{MC} / pt_{MC}",
158  //track parameters residuals - position
159  m_h1_x_res = createHistogram1D("h1xres", " residuals", 200, -0.05, 0.05, "x_{reco} - x_{MC} (cm)", m_histoList_trkQuality);
160  m_h1_y_res = createHistogram1D("h1yres", " residuals", 200, -0.05, 0.05, "y_{reco} - y_{MC} (cm)", m_histoList_trkQuality);
161  m_h1_z_res = createHistogram1D("h1zres", " residuals", 200, -0.05, 0.05, "z_{reco} - z_{MC} (cm)", m_histoList_trkQuality);
162  m_h1_r_res = createHistogram1D("h1rres", " residuals", 200, -0.05, 0.05, "r_{reco} - r_{MC} (cm)", m_histoList_trkQuality);
163  m_h1_rtot_res = createHistogram1D("h1rtotres", " residuals", 200, -0.05, 0.05, "rtot_{reco} - rtot_{MC} (cm)",
165 
166  m_h2_chargeVSchargeMC = createHistogram2D("h2chargecheck", "chargeVSchargeMC", 3, -1.5, 1.5, "charge MC", 3, -1.5, 1.5,
167  "charge reco", m_histoList_trkQuality);
168 
169  //track parameters pulls
170  m_h1_d0_pll = createHistogram1D("h1d0pll", "d0 pulls", 100, -5, 5, "d0 pull", m_histoList_trkQuality);
171  m_h1_phi_pll = createHistogram1D("h1phipll", "#phi pulls", 100, -5, 5, "#phi pull", m_histoList_trkQuality);
172  m_h1_omega_pll = createHistogram1D("h1omegapll", "#omega pulls", 100, -5, 5, "#omega pull", m_histoList_trkQuality);
173  m_h1_z0_pll = createHistogram1D("h1z0pll", "z0 pulls", 100, -5, 5, "z0 pull", m_histoList_trkQuality);
174  m_h1_cotTheta_pll = createHistogram1D("h1cotThetapll", "cot#theta pulls", 100, -5, 5, "cot#theta pull", m_histoList_trkQuality);
175 
176 
177  //first hit position using track parameters errors
178  m_h2_d0errphi0err_xy = createHistogram2D("h2d0errphierrXY", "#sigma_{d0}/#sigma_{#phi} projected on x,y",
179  2000, -10, 10, "x (cm)",
180  2000, -10, 10, "y (cm)",
182 
183  m_h2_d0errphi0err_rz = createHistogram2D("h2d0errphierrRZ", "#sigma_{d0}/#sigma_{#phi} projected on z and r_{t}=#sqrt{x^{2}+y^{2}}",
184  2000, -30, 40, "z (cm)",
185  2000, 0, 15, "r_{t} (cm)",
187 
188  m_h2_z0errcotThetaerr_xy = (TH2F*)duplicateHistogram("h2z0errcotThetaerrXY", "#sigma_{z0}/#sigma_{cot#theta} projected on x,y",
191 
192  m_h2_OmegaerrOmegaVSpt = createHistogram2D("h2OmegaerrOmegaVSpt", "#sigma_{#omega}/#omega VS p_{t}",
193  100, 0, 3, "p_{t} (GeV/c)",
194  1000, 0, 0.2, "#sigma_{#omega}/#omega",
196 
197 
198  m_h2_z0errVSpt = createHistogram2D("h2z0errVSpt", "#sigma_{z0} VS p_{t}",
199  100, 0, 3, "p_{t} (GeV/c)",
200  100, 0, 0.1, "#sigma_{z0} (cm)",
202 
203  m_h2_z0errVSpt_wtpxd = (TH2F*)duplicateHistogram("h2z0errVSpt_wTruePXD", "#sigma_{z0} VS p_{t}, with True PXD hits", m_h2_z0errVSpt,
205  m_h2_z0errVSpt_wfpxd = (TH2F*)duplicateHistogram("h2z0errVSpt_wFalsePXD", "#sigma_{z0} VS p_{t}, with False PXD hits",
208  m_h2_z0errVSpt_wpxd = (TH2F*)duplicateHistogram("h2z0errVSpt_wPXD", "#sigma_{z0} VS p_{t}, with PXD hits", m_h2_z0errVSpt,
210 
211  m_h2_z0errVSpt_wopxd = (TH2F*)duplicateHistogram("h2z0errVSpt_woPXD", "#sigma_{z0} VS p_{t}, no PXD hits", m_h2_z0errVSpt,
213 
214  m_h2_d0errVSpt = createHistogram2D("h2d0errVSpt", "#sigma_{d0} VS p_{t}",
215  100, 0, 3, "p_{t} (GeV/c)",
216  100, 0, 0.1, "#sigma_{d0} (cm)",
218  m_h2_d0errVSpt_wtpxd = (TH2F*)duplicateHistogram("h2d0errVSpt_wTruePXD", "#sigma_{d0} VS p_{t}, with True PXD hits", m_h2_d0errVSpt,
220  m_h2_d0errVSpt_wfpxd = (TH2F*)duplicateHistogram("h2d0errVSpt_wFalsePXD", "#sigma_{d0} VS p_{t}, with False PXD hits",
223  m_h2_d0errVSpt_wpxd = (TH2F*)duplicateHistogram("h2d0errVSpt_wPXD", "#sigma_{d0} VS p_{t}, with PXD hits", m_h2_d0errVSpt,
225 
226  m_h2_d0errVSpt_wopxd = (TH2F*)duplicateHistogram("h2d0errVSpt_woPXD", "#sigma_{d0} VS p_{t}, no PXD hits", m_h2_d0errVSpt,
228  m_h2_d0errMSVSpt = createHistogram2D("h2d0errMSVSpt", "#sigma_{d0} * #betapsin^{3/2}#theta VS p_{t}",
229  50, 0, 2.5, "p_{t} (GeV/c)",
230  500, 0, 1, "cm",
232 
233  //hits used in the fit
234  m_h2_TrackPointFitWeightVXD = createHistogram2D("h2TPFitWeightVXD", "VXD TrackPoint Fit Weight", 6, 0.5, 6.5, "VXD layer", 20, 0, 1,
235  "weight", m_histoList);
236  m_h2_TrackPointFitWeightCDC = createHistogram2D("h2TPFitWeightCDC", "CDC TrackPoint Fit Weight", 56, -0.5, 55.5, "CDC layer", 20, 0,
237  1, "weight", m_histoList);
238 
239  m_h1_nHitDetID = createHistogram1D("h1nHitDetID", "detector ID per hit", 4, -0.5, 3.5, "0=PXD, 1=SVD2D, 2=SVD,3=CDC", m_histoList);
240  m_h1_nCDChitsPR = createHistogram1D("h1nCDCHitsPR", "number of CDC hits from PR per Layer", 56, -0.5, 55.5, "CDC Layer",
241  m_histoList);
242  m_h1_nCDChitsWeighted = (TH1F*)duplicateHistogram("h1nCDCHitsWeighted", "CDC hits used in the fit per Layer, weighted",
244  m_h1_nCDChitsUsed = (TH1F*)duplicateHistogram("h1nCDCHitsUsed",
245  "approximated number of CDC hits used in the fit per Layer, weighted", m_h1_nCDChitsPR, m_histoList);
246  m_h1_nVXDhitsPR = createHistogram1D("h1nVXDHitsPR", "number of VXD hits from PR per Layer", 6, 0.5, 6.5, "VXD Layer", m_histoList);
247  m_h1_nVXDhitsWeighted = (TH1F*)duplicateHistogram("h1nVXDHitsWeighted", "number of VXD hits used in the fit per Layer, weighted",
249  m_h1_nVXDhitsUsed = (TH1F*)duplicateHistogram("h1nVXDHitsUsed",
250  "approximate number of VXD hits used in the fit per Layer, weighted", m_h1_nVXDhitsPR, m_histoList);
251  m_h2_VXDhitsPR_xy = createHistogram2D("h2hitsPRXY", "Pattern Recognition hits, transverse plane",
252  2000, -15, 15, "x (cm)",
253  2000, -15, 15, "y (cm)",
254  m_histoList);
255 
256  m_h2_VXDhitsPR_rz = createHistogram2D("h2hitsPRRZ", "Pattern Recognition Hits, r_{t} z",
257  2000, -30, 40, "z (cm)",
258  2000, 0, 15, "r_{t} (cm)",
259  m_histoList);
260 
261 
262 
263  //histograms to produce efficiency plots
264  Double_t bins_pt[10 + 1] = {0, 0.05, 0.1, 0.15, 0.2, 0.3, 0.5, 1, 1.5, 2, 3.5}; //GeV/c
265  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()};
266  Double_t bins_phi[14 + 1];
267  Double_t width_phi = 2 * TMath::Pi() / 14;
268  for (int bin = 0; bin < 14 + 1; bin++)
269  bins_phi[bin] = - TMath::Pi() + bin * width_phi;
270 
271 
272  m_h3_MCParticle = createHistogram3D("h3MCParticle", "entry per MCParticle",
273  10, bins_pt, "p_{t} (GeV/c)",
274  10, bins_theta, "#theta",
275  14, bins_phi, "#phi" /*, m_histoList*/);
276 
277  m_h3_TracksPerMCParticle = (TH3F*)duplicateHistogram("h3TracksPerMCParticle",
278  "entry per Track connected to a MCParticle",
279  m_h3_MCParticle /*, m_histoList*/);
280 
281  m_h3_TrackswPXDHitsPerMCParticle = (TH3F*)duplicateHistogram("h3TrackswPXDHitsPerMCParticle",
282  "entry per Track with PXD hits connected to a MCParticle",
283  m_h3_MCParticle /*, m_histoList*/);
284 
285  m_h3_RecoTrackswPXDHitsPerMCParticle = (TH3F*)duplicateHistogram("h3RecoTrackswPXDHitsPerMCParticle",
286  "entry per RecoTrack with PXD hits connected to a MCParticle",
287  m_h3_MCParticle /*, m_histoList*/);
288 
289  m_h3_RecoTrackswPXDHitsPerMCParticlewPXDHits = (TH3F*)duplicateHistogram("h3RecoTrackswPXDHitsPerMCParticlewPXDHits",
290  "entry per RecoTrack with PXD hits connected to a MCParticle with PXD hits",
291  m_h3_MCParticle /*, m_histoList*/);
292 
293  m_h3_MCParticleswPXDHits = (TH3F*)duplicateHistogram("h3MCParticleswPXDHitsPerMCParticle",
294  "entry per MCParticle with PXD hits",
295  m_h3_MCParticle /*, m_histoList*/);
296 
297  m_h3_MCRecoTrack = (TH3F*)duplicateHistogram("h3MCRecoTrack",
298  "entry per MCRecoTrack connected to the MCParticle",
299  m_h3_MCParticle /*, m_histoList*/);
300 
301  m_h3_TracksPerMCRecoTrack = (TH3F*)duplicateHistogram("h3TracksPerMCRecoTrack",
302  "entry per Track connected to an MCRecoTrack",
303  m_h3_MCParticle /*, m_histoList*/);
304  //plus
305  m_h3_MCParticle_plus = (TH3F*)duplicateHistogram("h3MCParticle_plus", "entry per positive MCParticle",
306  m_h3_MCParticle /*, m_histoList*/);
307 
308  m_h3_TracksPerMCParticle_plus = (TH3F*)duplicateHistogram("h3TracksPerMCParticle_plus",
309  "entry per Track connected to a positive MCParticle",
310  m_h3_MCParticle /*, m_histoList*/);
311 
312  m_h3_MCRecoTrack_plus = (TH3F*)duplicateHistogram("h3MCRecoTrack_plus",
313  "entry per MCRecoTrack connected to the positive MCParticle",
314  m_h3_MCParticle /*, m_histoList*/);
315 
316  m_h3_TracksPerMCRecoTrack_plus = (TH3F*)duplicateHistogram("h3TracksPerMCRecoTrack_plus",
317  "entry per Track connected to a positive MCRecoTrack",
318  m_h3_MCParticle /*, m_histoList*/);
319 
320 
321  //minus
322  m_h3_MCParticle_minus = (TH3F*)duplicateHistogram("h3MCParticlee_minus", "entry per negative MCParticle",
323  m_h3_MCParticle /*, m_histoList*/);
324 
325  m_h3_TracksPerMCParticle_minus = (TH3F*)duplicateHistogram("h3TracksPerMCParticle_minus",
326  "entry per Track connected to a negative MCParticle",
327  m_h3_MCParticle /*, m_histoList*/);
328 
329  m_h3_MCRecoTrack_minus = (TH3F*)duplicateHistogram("h3MCRecoTrack_minus",
330  "entry per MCRecoTrack connected to the negative MCParticle",
331  m_h3_MCParticle /*, m_histoList*/);
332 
333  m_h3_TracksPerMCRecoTrack_minus = (TH3F*)duplicateHistogram("h3TracksPerMCRecoTrack_minus",
334  "entry per Track connected to a negative MCRecoTrack",
335  m_h3_MCParticle /*, m_histoList*/);
336 
337  //histograms to produce efficiency plots
338  m_h1_HitsRecoTrackPerMCRecoTrack = createHistogram1D("h1hitsTCperMCRT", "RecoTrack per MCRecoTrack Hit in VXD layers", 6, 0.5, 6.5,
339  "# VXD layer" /*, m_histoList*/);
340 
341  m_h1_HitsMCRecoTrack = (TH1F*) duplicateHistogram("h1hitsMCRT", " MCRecoTrack Hit in VXD layers",
342  m_h1_HitsRecoTrackPerMCRecoTrack /*, m_histoList*/);
343 
344 
345  //histograms to produce purity plots
346  m_h3_Tracks = (TH3F*)duplicateHistogram("h3Tracks", "entry per Track",
347  m_h3_MCParticle /*, m_histoList*/);
348 
349  m_h3_MCParticlesPerTrack = (TH3F*)duplicateHistogram("h3MCParticlesPerTrack",
350  "entry per MCParticle connected to a Track",
351  m_h3_MCParticle /*, m_histoList*/);
352 }
353 
355 {
356 
357 }
358 
360 {
361 
363 
364  B2Vector3D magField = BFieldManager::getField(0, 0, 0) / Unit::T;
365 
366  bool hasTrack = false;
367  B2DEBUG(99, "+++++ 1. loop on MCParticles");
368  BOOST_FOREACH(MCParticle & mcParticle, mcParticles) {
369 
370  if (! isTraceable(mcParticle))
371  continue;
372 
373  int pdgCode = mcParticle.getPDG();
374  B2DEBUG(99, "MCParticle has PDG code " << pdgCode);
375 
376  int nFittedTracksMCRT = 0;
377  int nFittedTracks = 0;
378  int nFittedTrackswPXDHits = 0;
379 
380  MCParticleInfo mcParticleInfo(mcParticle, magField);
381 
382  hasTrack = false;
383 
384  m_h3_MCParticle->Fill(mcParticleInfo.getPt(), mcParticleInfo.getPtheta(), mcParticleInfo.getPphi());
385 
386  if (mcParticleInfo.getCharge() > 0)
387  m_h3_MCParticle_plus->Fill(mcParticleInfo.getPt(), mcParticleInfo.getPtheta(), mcParticleInfo.getPphi());
388  else if (mcParticleInfo.getCharge() < 0)
389  m_h3_MCParticle_minus->Fill(mcParticleInfo.getPt(), mcParticleInfo.getPtheta(), mcParticleInfo.getPphi());
390  else
391  continue;
392 
393  if (mcParticle.hasSeenInDetector(Const::PXD))
394  m_h3_MCParticleswPXDHits->Fill(mcParticleInfo.getPt(), mcParticleInfo.getPtheta(), mcParticleInfo.getPphi());
395 
396  //1. retrieve all the Tracks related to the MCParticle
397 
398  //1.0 check if there is a RecoTrack
399  RelationVector<RecoTrack> MCRecoTracks_fromMCParticle =
400  DataStore::getRelationsWithObj<RecoTrack>(&mcParticle, m_MCRecoTracksName);
401 
402  if (MCRecoTracks_fromMCParticle.size() > 0)
403  if (MCRecoTracks_fromMCParticle[0]->hasPXDHits()) {
404  m_h3_RecoTrackswPXDHitsPerMCParticle->Fill(mcParticleInfo.getPt(), mcParticleInfo.getPtheta(), mcParticleInfo.getPphi());
405  if (mcParticle.hasSeenInDetector(Const::PXD))
406  m_h3_RecoTrackswPXDHitsPerMCParticlewPXDHits->Fill(mcParticleInfo.getPt(), mcParticleInfo.getPtheta(), mcParticleInfo.getPphi());
407  }
408  m_multiplicityMCRecoTracks->Fill(MCRecoTracks_fromMCParticle.size());
409 
410  RelationVector<RecoTrack> RecoTracks_fromMCParticle =
411  DataStore::getRelationsWithObj<RecoTrack>(&mcParticle, m_RecoTracksName);
412 
413  m_multiplicityRecoTracks->Fill(RecoTracks_fromMCParticle.size());
414 
415  if (MCRecoTracks_fromMCParticle.size() > 0) {
416  m_h3_MCRecoTrack->Fill(mcParticleInfo.getPt(), mcParticleInfo.getPtheta(), mcParticleInfo.getPphi());
417 
418  if (mcParticleInfo.getCharge() > 0)
419  m_h3_MCRecoTrack_plus->Fill(mcParticleInfo.getPt(), mcParticleInfo.getPtheta(), mcParticleInfo.getPphi());
420  else if (mcParticleInfo.getCharge() < 0)
421  m_h3_MCRecoTrack_minus->Fill(mcParticleInfo.getPt(), mcParticleInfo.getPtheta(), mcParticleInfo.getPphi());
422  else
423  continue;
424  }
425 
426  //1.a retrieve all Tracks related to the MCParticle
427  RelationVector<Track> Tracks_fromMCParticle = DataStore::getRelationsWithObj<Track>(&mcParticle);
428  m_multiplicityTracks->Fill(Tracks_fromMCParticle.size());
429 
430  B2DEBUG(99, Tracks_fromMCParticle.size() << " Tracks related to this MCParticle");
431 
432  for (int trk = 0; trk < (int)Tracks_fromMCParticle.size(); trk++) {
433 
434  const TrackFitResult* fitResult = Tracks_fromMCParticle[trk]->getTrackFitResult(Const::ChargedStable(m_ParticleHypothesis));
435 
436  if ((fitResult == NULL) || (fitResult->getParticleType() != Const::ChargedStable(m_ParticleHypothesis)))
437  B2WARNING(" the TrackFitResult is not found!");
438 
439  else { // valid TrackFitResult found
440 
441  if (!hasTrack) {
442 
443  hasTrack = true;
444 
445  nFittedTracks++;
446 
447  if (fitResult->getHitPatternVXD().getNPXDHits() > 0) {
448  m_h3_TrackswPXDHitsPerMCParticle->Fill(mcParticleInfo.getPt(), mcParticleInfo.getPtheta(), mcParticleInfo.getPphi());
449  nFittedTrackswPXDHits++;
450  }
451 
452  m_h3_TracksPerMCParticle->Fill(mcParticleInfo.getPt(), mcParticleInfo.getPtheta(), mcParticleInfo.getPphi());
453  if (mcParticleInfo.getCharge() > 0)
454  m_h3_TracksPerMCParticle_plus->Fill(mcParticleInfo.getPt(), mcParticleInfo.getPtheta(), mcParticleInfo.getPphi());
455  else if (mcParticleInfo.getCharge() < 0)
456  m_h3_TracksPerMCParticle_minus->Fill(mcParticleInfo.getPt(), mcParticleInfo.getPtheta(), mcParticleInfo.getPphi());
457  else
458  continue;
459 
460  if (MCRecoTracks_fromMCParticle.size() > 0) {
461  nFittedTracksMCRT++;
462  m_h3_TracksPerMCRecoTrack->Fill(mcParticleInfo.getPt(), mcParticleInfo.getPtheta(), mcParticleInfo.getPphi());
463  if (mcParticleInfo.getCharge() > 0)
464  m_h3_TracksPerMCRecoTrack_plus->Fill(mcParticleInfo.getPt(), mcParticleInfo.getPtheta(), mcParticleInfo.getPphi());
465  else if (mcParticleInfo.getCharge() < 0)
466  m_h3_TracksPerMCRecoTrack_minus->Fill(mcParticleInfo.getPt(), mcParticleInfo.getPtheta(), mcParticleInfo.getPphi());
467  else
468  continue;
469  }
470 
471 
472  }
473 
474 
475  fillTrackParams1DHistograms(fitResult, mcParticleInfo);
476 
477  }
478  }
479 
480  m_multiplicityFittedTracks->Fill(nFittedTracks);
481  if (MCRecoTracks_fromMCParticle.size() > 0)
482  m_multiplicityFittedTracksPerMCRT->Fill(nFittedTracksMCRT);
483 
484  }
485 
486 
487  B2DEBUG(99, "+++++ 2. loop on Tracks");
488 
489  //2. retrieve all the MCParticles related to the Tracks
491 
492  BOOST_FOREACH(Track & track, tracks) {
493 
494  int nMCParticles = 0;
495 
496  //check if the track has been fitted
497  const TrackFitResult* fitResult = track.getTrackFitResult(Const::ChargedStable(m_ParticleHypothesis));
498 
499  if ((fitResult == NULL) || (fitResult->getParticleType() != Const::ChargedStable(m_ParticleHypothesis)))
500  continue;
501 
502  m_h1_pValue->Fill(fitResult->getPValue());
503 
504  TVector3 momentum = fitResult->getMomentum();
505  m_h3_Tracks->Fill(momentum.Pt(), momentum.Theta(), momentum.Phi());
506 
508 
510 
511  for (int layer = 0; layer < 56; layer++) {
512  if (fitResult->getHitPatternCDC().hasLayer(layer))
513  m_h1_nCDChitsUsed->Fill(layer);
514  }
515  for (int layer = 1; layer <= 2; layer++) {
516  for (int i = 0; i < fitResult->getHitPatternVXD().getPXDLayer(layer); i++)
517  m_h1_nVXDhitsUsed->Fill(layer);
518  }
519  for (int layer = 3; layer <= 6; layer++) {
520  int n1 = fitResult->getHitPatternVXD().getSVDLayer(layer).first;
521  int n2 = fitResult->getHitPatternVXD().getSVDLayer(layer).second;
522  int N = n1 + n2;
523 
524  for (int i = 0; i < N; i++)
525  m_h1_nVXDhitsUsed->Fill(layer);
526  }
527 
528 
529  RelationVector<MCParticle> MCParticles_fromTrack = DataStore::getRelationsWithObj<MCParticle>(&track);
530 
531  for (int mcp = 0; mcp < (int)MCParticles_fromTrack.size(); mcp++)
532  if (isTraceable(*MCParticles_fromTrack[mcp])) {
533  nMCParticles ++;
534  m_h3_MCParticlesPerTrack->Fill(momentum.Pt(), momentum.Theta(), momentum.Phi());
535  }
536  // }
537 
538  m_multiplicityMCParticlesPerTrack->Fill(nMCParticles);
539 
540  }
541 
542 
543  B2DEBUG(99, "+++++ 3. loop on MCRecoTracks");
544 
545  //3. retrieve all MCRecoTracks
547  StoreArray<PXDCluster> pxdClusters;
548  StoreArray<SVDCluster> svdClusters;
549  StoreArray<CDCHit> cdcHit;
550 
551 
552  BOOST_FOREACH(RecoTrack & mcRecoTrack, mcRecoTracks) {
553 
554  int nRecoTrack = 0;
555  bool hasRecoTrack = false;
556 
557  //3.a retrieve the RecoTrack
558  RelationVector<RecoTrack> RecoTracks_fromMCRecoTrack = DataStore::getRelationsWithObj<RecoTrack>(&mcRecoTrack);
559  B2DEBUG(99, "~ " << RecoTracks_fromMCRecoTrack.size() << " RecoTracks related to this MCRecoTrack");
560  m_multiplicityRecoTracksPerMCRT->Fill(RecoTracks_fromMCRecoTrack.size());
561 
562  //3.a retrieve the MCParticle
563  RelationVector<MCParticle> MCParticles_fromMCRecoTrack = DataStore::getRelationsWithObj<MCParticle>(&mcRecoTrack);
564 
565  B2DEBUG(99, "~~~ " << MCParticles_fromMCRecoTrack.size() << " MCParticles related to this MCRecoTrack");
566  for (int mcp = 0; mcp < (int)MCParticles_fromMCRecoTrack.size(); mcp++) {
567 
568  //3.b retrieve all RecoTracks related to the MCRecoTrack
569  RelationVector<RecoTrack> RecoTracks_fromMCParticle = DataStore::getRelationsWithObj<RecoTrack>
570  (MCParticles_fromMCRecoTrack[mcp]);
571 
572  B2DEBUG(99, "~~~~~ " << RecoTracks_fromMCParticle.size() << " RecoTracks related to this MCParticle");
573  for (int tc = 0; tc < (int)RecoTracks_fromMCParticle.size(); tc++)
574  if (!hasRecoTrack) {
575  hasRecoTrack = true;
576  nRecoTrack++;
577  }
578 
579  }
580 
581  }
582 
583 
584  B2DEBUG(99, "+++++ 4. loop on RecoTracks");
585 
586  //4. retrieve all RecoTracks
587  StoreArray<RecoTrack> RecoTracks;
588 
589  BOOST_FOREACH(RecoTrack & recoTrack, RecoTracks) {
590 
591  // int nMCRecoTrack = 0;
592 
593  // retrieve the MCRecoTrack
594  RelationVector<RecoTrack> MCRecoTracks_fromRecoTrack = DataStore::getRelationsWithObj<RecoTrack>(&recoTrack, m_MCRecoTracksName);
595  m_multiplicityMCRecoTracksPerRT->Fill(MCRecoTracks_fromRecoTrack.size());
596 
597 
598  /*
599  //4.a retrieve the MCParticle
600  RelationVector<MCParticle> MCParticles_fromRecoTrack = DataStore::getRelationsWithObj<MCParticle>(&recoTrack);
601 
602  B2DEBUG(99, "~~~ " << MCParticles_fromRecoTrack.size() << " MCParticles related to this RecoTrack");
603  for (int mcp = 0; mcp < (int)MCParticles_fromRecoTrack.size(); mcp++) {
604 
605  //4.b retrieve all MCRecoTracks related to the RecoTrack
606  RelationVector<RecoTrack> mcRecoTracks_fromMCParticle = DataStore::getRelationsWithObj<RecoTrack>
607  (MCParticles_fromRecoTrack[mcp], m_MCRecoTracksName);
608 
609  B2DEBUG(99, "~~~~~ " << mcRecoTracks_fromMCParticle.size() << " MCRecoTracks related to this MCParticle");
610  for (int mctc = 0; mctc < (int)mcRecoTracks_fromMCParticle.size(); mctc++) {
611  nMCRecoTrack++;
612 
613  }
614  }
615 
616  // m_multiplicityMCRecoTracksPerRT->Fill(nMCRecoTrack);
617  */
618  }
619 
620 }
621 
623 {
624 
625  double num = 0;
626  double den = 0;
627 
628  for (int bin = 1; bin < m_multiplicityFittedTracks->GetNbinsX(); bin ++)
629  num += m_multiplicityFittedTracks->GetBinContent(bin + 1);
630  den = m_multiplicityFittedTracks->GetEntries();
631  double efficiency = num / den ;
632  double efficiencyErr = sqrt(efficiency * (1 - efficiency)) / sqrt(den);
633 
634  double nFittedTracksMCRT = 0;
635  for (int bin = 1; bin < m_multiplicityFittedTracksPerMCRT->GetNbinsX(); bin ++)
636  nFittedTracksMCRT += m_multiplicityFittedTracksPerMCRT->GetBinContent(bin + 1);
637  double efficiencyMCRT = nFittedTracksMCRT / m_multiplicityFittedTracksPerMCRT->GetEntries();
638  double efficiencyMCRTErr = sqrt(efficiencyMCRT * (1 - efficiencyMCRT)) / sqrt(m_multiplicityFittedTracksPerMCRT->GetEntries());
639 
640  double nRecoTrack = 0;
641  for (int bin = 1; bin < m_multiplicityRecoTracksPerMCRT->GetNbinsX(); bin ++)
642  nRecoTrack += m_multiplicityRecoTracksPerMCRT->GetBinContent(bin + 1);
643  double efficiencyPR = nRecoTrack / m_multiplicityRecoTracksPerMCRT->GetEntries();
644  double efficiencyPRErr = sqrt(efficiencyPR * (1 - efficiencyPR)) / sqrt(m_multiplicityRecoTracksPerMCRT->GetEntries());
645 
646  double nMCRecoTrack = 0;
647  for (int bin = 1; bin < m_multiplicityMCRecoTracksPerRT->GetNbinsX(); bin ++)
648  nMCRecoTrack += m_multiplicityMCRecoTracksPerRT->GetBinContent(bin + 1);
649  double purityPR = nMCRecoTrack / m_multiplicityMCRecoTracksPerRT->GetEntries();
650  double purityPRErr = sqrt(purityPR * (1 - purityPR)) / sqrt(m_multiplicityMCRecoTracksPerRT->GetEntries());
651 
652  double nMCParticles = 0;
653  for (int bin = 1; bin < m_multiplicityMCParticlesPerTrack->GetNbinsX(); bin ++)
654  nMCParticles += m_multiplicityMCParticlesPerTrack->GetBinContent(bin + 1);
655  double purity = nMCParticles / m_multiplicityMCParticlesPerTrack->GetEntries();
656  double purityErr = sqrt(purity * (1 - purity)) / sqrt(m_multiplicityMCParticlesPerTrack->GetEntries());
657 
658  B2INFO("");
659  B2INFO("~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~");
660  B2INFO("~ Tracking Performance Evaluation ~ SHORT SUMMARY ~");
661  B2INFO("");
662  B2INFO(" + overall:");
663  B2INFO(" efficiency = (" << efficiency * 100 << " +/- " << efficiencyErr * 100 << ")% ");
664  B2INFO(" purity = " << purity * 100 << " +/- " << purityErr * 100 << ")% ");
665  B2INFO("");
666  B2INFO(" + factorizing geometrical acceptance:");
667  B2INFO(" efficiency = " << efficiencyMCRT * 100 << " +/- " << efficiencyMCRTErr * 100 << ")% ");
668  B2INFO("");
669  B2INFO(" + pattern recognition:");
670  B2INFO(" efficiency = " << efficiencyPR * 100 << " +/- " << efficiencyPRErr * 100 << ")% ");
671  B2INFO(" purity = " << purityPR * 100 << " +/- " << purityPRErr * 100 << ")% ");
672  B2INFO("");
673  B2INFO("~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~");
674 }
675 
677 {
678 
679 
682 
684 
685  if (m_rootFilePtr != NULL) {
686  m_rootFilePtr->cd();
687 
688  TDirectory* oldDir = gDirectory;
689 
690  TDirectory* dir_multiplicity = oldDir->mkdir("multiplicity");
691  dir_multiplicity->cd();
692  TIter nextH_multiplicity(m_histoList_multiplicity);
693  TObject* obj;
694  while ((obj = nextH_multiplicity()))
695  obj->Write();
696 
697  TDirectory* dir_efficiency = oldDir->mkdir("efficiency");
698  dir_efficiency->cd();
699  TIter nextH_efficiency(m_histoList_efficiency);
700  while ((obj = nextH_efficiency()))
701  obj->Write();
702 
703  TDirectory* dir_trkQuality = oldDir->mkdir("trkQuality");
704  dir_trkQuality->cd();
705  TIter nextH_trkQuality(m_histoList_trkQuality);
706  while ((obj = nextH_trkQuality()))
707  obj->Write();
708 
709  TDirectory* dir_firstHit = oldDir->mkdir("firstHit");
710  dir_firstHit->cd();
711  TIter nextH_firstHit(m_histoList_firstHit);
712  while ((obj = nextH_firstHit()))
713  obj->Write();
714 
715 
716 
717 
718  m_rootFilePtr->Close();
719  }
720 
721 }
722 
723 
725  MCParticleInfo mcParticleInfo)
726 {
727 
728  //track parameters errors
729  double d0_err = sqrt((fitResult->getCovariance5())[0][0]);
730  double phi_err = sqrt((fitResult->getCovariance5())[1][1]);
731  double omega_err = sqrt((fitResult->getCovariance5())[2][2]);
732  double z0_err = sqrt((fitResult->getCovariance5())[3][3]);
733  double cotTheta_err = sqrt((fitResult->getCovariance5())[4][4]);
734 
735  //track parameters residuals:
736  double d0_res = fitResult->getD0() - mcParticleInfo.getD0();
737  double phi_res = TMath::ASin(TMath::Sin(fitResult->getPhi() - mcParticleInfo.getPhi())); //giulia
738  double omega_res = fitResult->getOmega() - mcParticleInfo.getOmega();
739  double z0_res = fitResult->getZ0() - mcParticleInfo.getZ0();
740  double cotTheta_res = fitResult->getCotTheta() - mcParticleInfo.getCotTheta();
741 
742  //track parameters residuals in momentum:
743  double px_res = fitResult->getMomentum().X() - mcParticleInfo.getPx();
744  double py_res = fitResult->getMomentum().Y() - mcParticleInfo.getPy();
745  double pz_res = fitResult->getMomentum().Z() - mcParticleInfo.getPz();
746  double p_res = (fitResult->getMomentum().Mag() - mcParticleInfo.getP()) / mcParticleInfo.getP();
747  double pt_res = (fitResult->getMomentum().Pt() - mcParticleInfo.getPt()) / mcParticleInfo.getPt();
748 
749  //track parameters residuals in position:
750  double x_res = fitResult->getPosition().X() - mcParticleInfo.getX();
751  double y_res = fitResult->getPosition().Y() - mcParticleInfo.getY();
752  double z_res = fitResult->getPosition().Z() - mcParticleInfo.getZ();
753  double r_res = fitResult->getPosition().Perp() - sqrt(mcParticleInfo.getX() * mcParticleInfo.getX() + mcParticleInfo.getY() *
754  mcParticleInfo.getY());
755  double rtot_res = fitResult->getPosition().Mag() - sqrt(mcParticleInfo.getX() * mcParticleInfo.getX() + mcParticleInfo.getY() *
756  mcParticleInfo.getY() + mcParticleInfo.getZ() * mcParticleInfo.getZ());
757 
758  m_h1_d0_err->Fill(d0_err);
759  m_h1_phi_err->Fill(phi_err);
760  m_h1_omega_err->Fill(omega_err);
761  m_h1_z0_err->Fill(z0_err);
762  m_h1_cotTheta_err->Fill(cotTheta_err);
763 
764  m_h1_d0_res->Fill(d0_res);
765  m_h1_phi_res->Fill(phi_res);
766  m_h1_omega_res->Fill(omega_res);
767  m_h1_z0_res->Fill(z0_res);
768  m_h1_cotTheta_res->Fill(cotTheta_res);
769 
770  m_h1_px_res->Fill(px_res);
771  m_h1_py_res->Fill(py_res);
772  m_h1_pz_res->Fill(pz_res);
773  m_h1_p_res->Fill(p_res);
774  m_h1_pt_res->Fill(pt_res);
775 
776  m_h1_x_res->Fill(x_res);
777  m_h1_y_res->Fill(y_res);
778  m_h1_z_res->Fill(z_res);
779  m_h1_r_res->Fill(r_res);
780  m_h1_rtot_res->Fill(rtot_res);
781 
782  m_h2_chargeVSchargeMC->Fill(mcParticleInfo.getCharge(), fitResult->getChargeSign());
783 
784  m_h1_d0_pll->Fill(d0_res / d0_err);
785  m_h1_phi_pll->Fill(phi_res / phi_err);
786  m_h1_omega_pll->Fill(omega_res / omega_err);
787  m_h1_z0_pll->Fill(z0_res / z0_err);
788  m_h1_cotTheta_pll->Fill(cotTheta_res / cotTheta_err);
789 
790 
791  m_h2_OmegaerrOmegaVSpt->Fill(fitResult->getMomentum().Pt(), omega_err / mcParticleInfo.getOmega());
792 
793 
794 }
795 
797 {
798 
799 
800  double d0_err = sqrt((fitResult->getCovariance5())[0][0]);
801  double phi_err = sqrt((fitResult->getCovariance5())[1][1]);
802  double z0_err = sqrt((fitResult->getCovariance5())[3][3]);
803  double cotTheta_err = sqrt((fitResult->getCovariance5())[4][4]);
804 
805  TVector3 momentum = fitResult->getMomentum();
806 
807  double px = momentum.Px();
808  double py = momentum.Py();
809  double pz = momentum.Pz();
810  double pt = momentum.Pt();
811  double p = momentum.Mag();
812  double mass = fitResult->getParticleType().getMass();
813  double beta = p / sqrt(p * p + mass * mass);
814  double sinTheta = TMath::Sin(momentum.Theta());
815 
816  m_h2_d0errphi0err_xy->Fill(d0_err / phi_err * px / pt,
817  d0_err / phi_err * py / pt);
818  m_h2_z0errcotThetaerr_xy->Fill(z0_err / cotTheta_err * px / pt,
819  z0_err / cotTheta_err * py / pt);
820  m_h2_d0errphi0err_rz->Fill(d0_err / phi_err * pz / pt,
821  d0_err / phi_err);
822 
823  m_h2_d0errVSpt->Fill(pt, d0_err);
824 
825  m_h2_z0errVSpt->Fill(pt, z0_err);
826 
827  m_h2_d0errMSVSpt->Fill(pt, d0_err * beta * p * pow(sinTheta, 3 / 2) / 0.0136);
828 
829 }
830 
832 {
833 
834  //hits used in the fit
835 
837 
839 
840  bool hasPXDhit = false;
841  bool isTrueHit = false;
842 
843  double d0_err = -999;
844  double z0_err = -999;
845  double pt = -999;
846 
847 
848  if ((fitResult != NULL)
849  || (fitResult->getParticleType() != Const::ChargedStable(m_ParticleHypothesis))) { // valid TrackFitResult found
850  d0_err = sqrt((fitResult->getCovariance5())[0][0]);
851  z0_err = sqrt((fitResult->getCovariance5())[3][3]);
852  pt = fitResult->getMomentum().Pt();
853  }
854 
855  bool hasCDChit[56] = { false };
856 
857  RelationVector<RecoTrack> RecoTracks_fromTrack = DataStore::getRelationsWithObj<RecoTrack>(&theTrack);
858 
859  for (int tc = 0; tc < (int)RecoTracks_fromTrack.size(); tc++) {
860 
861  const std::vector< genfit::TrackPoint* >& tp_vector = RecoTracks_fromTrack[tc]->getHitPointsWithMeasurement();
862  for (int i = 0; i < (int) tp_vector.size(); i++) {
863  genfit::TrackPoint* tp = tp_vector[i];
864 
865  int nMea = tp->getNumRawMeasurements();
866  for (int mea = 0; mea < nMea; mea++) {
867 
868  genfit::AbsMeasurement* absMeas = tp->getRawMeasurement(mea);
869  double weight = 0;
870 
871  std::vector<double> weights;
872  genfit::KalmanFitterInfo* kalmanInfo = tp->getKalmanFitterInfo();
873  if (kalmanInfo)
874  weights = kalmanInfo->getWeights();
875  else //no kalman fitter info, fill the weights vector with 0 (VXD), or 0,0 (CDC)
876  B2WARNING(" No KalmanFitterInfo associated to the TrackPoint!");
877 
878  double detId(-999);
879  TVector3 globalHit(-999, -999, -999);
880 
881  PXDRecoHit* pxdHit = dynamic_cast<PXDRecoHit*>(absMeas);
882  SVDRecoHit2D* svdHit2D = dynamic_cast<SVDRecoHit2D*>(absMeas);
883  SVDRecoHit* svdHit = dynamic_cast<SVDRecoHit*>(absMeas);
884  CDCRecoHit* cdcHit = dynamic_cast<CDCRecoHit*>(absMeas);
885 
886  if (pxdHit) {
887  hasPXDhit = true;
888  isTrueHit = false;
889 
890  if (kalmanInfo)
891  weight = weights.at(mea);
892 
893  detId = 0;
894  double uCoor = pxdHit->getU();
895  double vCoor = pxdHit->getV();
896  VxdID sensor = pxdHit->getSensorID();
897 
898  m_h1_nVXDhitsPR->Fill(sensor.getLayerNumber());
899 
900  m_h1_nVXDhitsWeighted->Fill(sensor.getLayerNumber() , weight);
901 
902  m_h2_TrackPointFitWeightVXD->Fill(sensor.getLayerNumber(), weight);
903  const VXD::SensorInfoBase& aSensorInfo = aGeometry.getSensorInfo(sensor);
904  globalHit = aSensorInfo.pointToGlobal(TVector3(uCoor, vCoor, 0), true);
905 
906 
907  const PXDCluster* pxdcl = pxdHit->getCluster();
908  RelationVector<PXDTrueHit> pxdth_fromcl = DataStore::getRelationsWithObj<PXDTrueHit>(pxdcl);
909 
910  if ((int)pxdth_fromcl.size() != 0) {
911  const PXDTrueHit* trueHit = pxdth_fromcl[0];
912 
913  if (trueHit) {
914  int trueHitIndex = trueHit->getArrayIndex();
915  RelationVector<MCParticle> MCParticles_fromTrack = DataStore::getRelationsWithObj<MCParticle>(&theTrack);
916  for (int mcp = 0; mcp < (int)MCParticles_fromTrack.size(); mcp++) {
917  RelationVector<PXDTrueHit> trueHit_fromMCParticles = DataStore::getRelationsWithObj<PXDTrueHit>(MCParticles_fromTrack[mcp]);
918  for (int th = 0; th < (int)trueHit_fromMCParticles.size(); th++) {
919  if (trueHit_fromMCParticles[th]->getArrayIndex() == trueHitIndex)
920  isTrueHit = true;
921  }
922  }
923  }
924  }
925 
926  } else if (svdHit2D) {
927 
928  if (kalmanInfo)
929  weight = weights.at(mea);
930 
931  detId = 1;
932  double uCoor = svdHit2D->getU();
933  double vCoor = svdHit2D->getV();
934  VxdID sensor = svdHit2D->getSensorID();
935 
936  m_h1_nVXDhitsPR->Fill(sensor.getLayerNumber());
937 
938  m_h1_nVXDhitsWeighted->Fill(sensor.getLayerNumber() , weight);
939 
940  m_h2_TrackPointFitWeightVXD->Fill(sensor.getLayerNumber(), weight);
941 
942  const VXD::SensorInfoBase& aSensorInfo = aGeometry.getSensorInfo(sensor);
943  globalHit = aSensorInfo.pointToGlobal(TVector3(uCoor, vCoor, 0), true);
944 
945  } else if (svdHit) {
946 
947  if (kalmanInfo)
948  weight = weights.at(mea);
949 
950  detId = 2;
951  double uCoor = 0;
952  double vCoor = 0;
953  if (svdHit->isU())
954  uCoor = svdHit->getPosition();
955  else
956  vCoor = svdHit->getPosition();
957 
958  VxdID sensor = svdHit->getSensorID();
959  m_h1_nVXDhitsPR->Fill(sensor.getLayerNumber());
960 
961  m_h1_nVXDhitsWeighted->Fill(sensor.getLayerNumber() , weight);
962 
963  m_h2_TrackPointFitWeightVXD->Fill(sensor.getLayerNumber(), weight);
964  const VXD::SensorInfoBase& aSensorInfo = aGeometry.getSensorInfo(sensor);
965  globalHit = aSensorInfo.pointToGlobal(TVector3(uCoor, vCoor, 0), true);
966  } else if (cdcHit) {
967 
968  if (kalmanInfo)
969  weight = weights.at(mea);
970 
971  WireID wire = cdcHit->getWireID();
972  if (! hasCDChit[wire.getICLayer()]) { //needed to validate the HitPatternCDC filling
973  m_h1_nCDChitsPR->Fill(wire.getICLayer());
974 
975  m_h1_nCDChitsWeighted->Fill(wire.getICLayer() , weight);
976  // hasCDChit[wire.getICLayer()] = true; //to validate the HitPatternCDC filling: uncomment this
977  }
978  m_h2_TrackPointFitWeightCDC->Fill(wire.getICLayer(), weight);
979  detId = 3;
980  }
981 
982 
983  m_h1_nHitDetID ->Fill(detId);
984 
985  m_h2_VXDhitsPR_xy->Fill(globalHit.X(), globalHit.Y());
986 
987  m_h2_VXDhitsPR_rz->Fill(globalHit.Z(), globalHit.Perp());
988 
989  }
990 
991  }
992  }
993 
994  if ((fitResult != NULL) && (fitResult->getParticleType() == Const::ChargedStable(m_ParticleHypothesis))) {
995  if (hasPXDhit) {
996  m_h2_d0errVSpt_wpxd->Fill(pt, d0_err);
997  m_h2_z0errVSpt_wpxd->Fill(pt, z0_err);
998  if (isTrueHit) {
999  m_h2_d0errVSpt_wtpxd->Fill(pt, d0_err);
1000  m_h2_z0errVSpt_wtpxd->Fill(pt, z0_err);
1001  } else {
1002  m_h2_d0errVSpt_wfpxd->Fill(pt, d0_err);
1003  m_h2_z0errVSpt_wfpxd->Fill(pt, z0_err);
1004  }
1005  } else {
1006  m_h2_d0errVSpt_wopxd->Fill(pt, d0_err);
1007  m_h2_z0errVSpt_wopxd->Fill(pt, z0_err);
1008  }
1009  }
1010 
1011 }
1012 
1014 {
1015 
1016  //normalized to MCParticles
1017  TH1F* h_ineff_pt = createHistogramsRatio("hineffpt", "inefficiency VS pt, normalized to MCParticles", m_h3_TracksPerMCParticle,
1018  m_h3_MCParticle, false, 0);
1019  histoList->Add(h_ineff_pt);
1020 
1021  TH1F* h_ineff_theta = createHistogramsRatio("hinefftheta", "inefficiency VS #theta, normalized to MCParticles",
1023  histoList->Add(h_ineff_theta);
1024 
1025  TH1F* h_ineff_phi = createHistogramsRatio("hineffphi", "inefficiency VS #phi, normalized to MCParticles", m_h3_TracksPerMCParticle,
1026  m_h3_MCParticle, false, 2);
1027  histoList->Add(h_ineff_phi);
1028 
1029  //normalized to MCRecoTracks
1030  TH1F* h_ineffMCRT_pt = createHistogramsRatio("hineffMCRTpt", "inefficiency VS pt, normalized to MCRecoTrack",
1032  histoList->Add(h_ineffMCRT_pt);
1033 
1034  TH1F* h_ineffMCRT_theta = createHistogramsRatio("hineffMCRTtheta", "inefficiency VS #theta, normalized to MCRecoTrack",
1036  histoList->Add(h_ineffMCRT_theta);
1037 
1038  TH1F* h_ineffMCRT_phi = createHistogramsRatio("hineffMCRTphi", "inefficiency VS #phi, normalized to MCRecoTrack",
1040  histoList->Add(h_ineffMCRT_phi);
1041 
1042 }
1043 
1045 {
1046 
1047 
1048  TH1F* h_MCPwPXDhits_pt = createHistogramsRatio("hMCPwPXDhits", "fraction of MCParticles with PXD hits VS pt",
1050  m_h3_MCParticle, true, 0);
1051  histoList->Add(h_MCPwPXDhits_pt);
1052 
1053  TH1F* h_RTwPXDhitsMCPwPXDHits_pt = createHistogramsRatio("hRecoTrkswPXDhitsMCPwPXDHits",
1054  "fraction of MCParticles with PXD Hits with RecoTracks with PXD hits VS pt",
1056  m_h3_MCParticleswPXDHits, true, 0);
1057  histoList->Add(h_RTwPXDhitsMCPwPXDHits_pt);
1058 
1059  TH1F* h_wPXDhits_pt = createHistogramsRatio("hTrkswPXDhits", "fraction of tracks with PXD hits VS pt",
1061  m_h3_TracksPerMCParticle, true, 0);
1062  histoList->Add(h_wPXDhits_pt);
1063 
1064  //normalized to MCParticles
1065  TH1F* h_eff_pt = createHistogramsRatio("heffpt", "efficiency VS pt, normalized to MCParticles", m_h3_TracksPerMCParticle,
1066  m_h3_MCParticle, true, 0);
1067  histoList->Add(h_eff_pt);
1068  // B2INFO(" efficiency in pt, NUM = "<<m_nFittedTracks<<", DEN = "<<m_nMCParticles<<", eff integrata = "<<(double)m_nFittedTracks/m_nMCParticles);
1069 
1070  TH1F* h_eff_theta = createHistogramsRatio("hefftheta", "efficiency VS #theta, normalized to MCParticles", m_h3_TracksPerMCParticle,
1071  m_h3_MCParticle, true, 1);
1072  histoList->Add(h_eff_theta);
1073 
1074  TH1F* h_eff_phi = createHistogramsRatio("heffphi", "efficiency VS #phi, normalized to MCParticles", m_h3_TracksPerMCParticle,
1075  m_h3_MCParticle, true, 2);
1076  histoList->Add(h_eff_phi);
1077 
1078  //normalized to MCRecoTracks
1079  TH1F* h_effMCRT_pt = createHistogramsRatio("heffMCRTpt", "efficiency VS pt, normalized to MCRecoTrack", m_h3_TracksPerMCRecoTrack,
1080  m_h3_MCRecoTrack, true, 0);
1081  histoList->Add(h_effMCRT_pt);
1082 
1083  TH1F* h_effMCRT_theta = createHistogramsRatio("heffMCRTtheta", "efficiency VS #theta, normalized to MCRecoTrack",
1085  histoList->Add(h_effMCRT_theta);
1086 
1087  TH1F* h_effMCRT_phi = createHistogramsRatio("heffMCRTphi", "efficiency VS #phi, normalized to MCRecoTrack",
1089  histoList->Add(h_effMCRT_phi);
1090 
1091  // plus
1092 
1093  //normalized to MCParticles
1094  TH1F* h_eff_pt_plus = createHistogramsRatio("heffpt_plus", "efficiency VS pt, normalized to positive MCParticles",
1096  histoList->Add(h_eff_pt_plus);
1097  // B2INFO(" efficiency in pt, NUM = "<<m_nFittedTracks<<", DEN = "<<m_nMCParticles<<", eff integrata = "<<(double)m_nFittedTracks/m_nMCParticles);
1098 
1099  TH1F* h_eff_theta_plus = createHistogramsRatio("hefftheta_plus", "efficiency VS #theta, normalized to positive MCParticles",
1101  histoList->Add(h_eff_theta_plus);
1102 
1103  TH1F* h_eff_phi_plus = createHistogramsRatio("heffphi_plus", "efficiency VS #phi, normalized to positive MCParticles",
1105  histoList->Add(h_eff_phi_plus);
1106 
1107  //normalized to MCRecoTracks
1108  TH1F* h_effMCRT_pt_plus = createHistogramsRatio("heffMCRTpt_plus", "efficiency VS pt, normalized to positive MCRecoTrack",
1110  histoList->Add(h_effMCRT_pt_plus);
1111 
1112  TH1F* h_effMCRT_theta_plus = createHistogramsRatio("heffMCRTtheta_plus", "efficiency VS #theta, normalized to positive MCRecoTrack",
1114  histoList->Add(h_effMCRT_theta_plus);
1115 
1116  TH1F* h_effMCRT_phi_plus = createHistogramsRatio("heffMCRTphi_plus", "efficiency VS #phi, normalized to positive MCRecoTrack",
1118  histoList->Add(h_effMCRT_phi_plus);
1119 
1120  // minus
1121 
1122  //normalized to MCParticles
1123  TH1F* h_eff_pt_minus = createHistogramsRatio("heffpt_minus", "efficiency VS pt, normalized to positive MCParticles",
1125  histoList->Add(h_eff_pt_minus);
1126  // B2INFO(" efficiency in pt, NUM = "<<m_nFittedTracks<<", DEN = "<<m_nMCParticles<<", eff integrata = "<<(double)m_nFittedTracks/m_nMCParticles);
1127 
1128  TH1F* h_eff_theta_minus = createHistogramsRatio("hefftheta_minus", "efficiency VS #theta, normalized to positive MCParticles",
1130  histoList->Add(h_eff_theta_minus);
1131 
1132  TH1F* h_eff_phi_minus = createHistogramsRatio("heffphi_minus", "efficiency VS #phi, normalized to positive MCParticles",
1134  histoList->Add(h_eff_phi_minus);
1135 
1136  //normalized to MCRecoTracks
1137  TH1F* h_effMCRT_pt_minus = createHistogramsRatio("heffMCRTpt_minus", "efficiency VS pt, normalized to positive MCRecoTrack",
1139  histoList->Add(h_effMCRT_pt_minus);
1140 
1141  TH1F* h_effMCRT_theta_minus = createHistogramsRatio("heffMCRTtheta_minus",
1142  "efficiency VS #theta, normalized to positive MCRecoTrack", m_h3_TracksPerMCRecoTrack_minus, m_h3_MCRecoTrack_minus, true, 1);
1143  histoList->Add(h_effMCRT_theta_minus);
1144 
1145  TH1F* h_effMCRT_phi_minus = createHistogramsRatio("heffMCRTphi_minus", "efficiency VS #phi, normalized to positive MCRecoTrack",
1147  histoList->Add(h_effMCRT_phi_minus);
1148 
1149  //pattern recognition efficiency
1150  TH1F* h_effPR = createHistogramsRatio("heffPR", "PR efficiency VS VXD Layer, normalized to MCRecoTrack",
1152  histoList->Add(h_effPR);
1153 
1154  //tracks used in the fit
1155  TH1F* h_effVXDHitFit = createHistogramsRatio("heffVXDHitFit",
1156  "weighted hits used in the fit VS VXD Layer, normalized to hits form PR", m_h1_nVXDhitsWeighted, m_h1_nVXDhitsPR, true, 0);
1157  histoList->Add(h_effVXDHitFit);
1158 
1159  TH1F* h_effCDCHitFit = createHistogramsRatio("heffCDCHitFit",
1160  "weighted hits used in the fit VS CDC Layer, normalized to hits form PR", m_h1_nCDChitsWeighted, m_h1_nCDChitsPR, true, 0);
1161  histoList->Add(h_effCDCHitFit);
1162 
1163 }
1164 
1165 
1166 
1167 
1169 {
1170 
1171  bool isChargedStable = Const::chargedStableSet.find(abs(the_mcParticle.getPDG())) != Const::invalidParticle;
1172 
1173  bool isPrimary = the_mcParticle.hasStatus(MCParticle::c_PrimaryParticle);
1174 
1175  return (isPrimary && isChargedStable);
1176 
1177 }
Belle2::RelationVector::size
size_t size() const
Get number of relations.
Definition: RelationVector.h:98
Belle2::TrackingPerformanceEvaluationModule::m_h2_d0errphi0err_xy
TH2F * m_h2_d0errphi0err_xy
error
Definition: TrackingPerformanceEvaluationModule.h:132
Belle2::TrackingPerformanceEvaluationModule::m_h1_cotTheta_res
TH1F * m_h1_cotTheta_res
error
Definition: TrackingPerformanceEvaluationModule.h:110
Belle2::TrackingPerformanceEvaluationModule::beginRun
void beginRun() override
Called when entering a new run.
Definition: TrackingPerformanceEvaluationModule.cc:354
Belle2::TrackingPerformanceEvaluationModule::m_h1_phi_pll
TH1F * m_h1_phi_pll
pull distribution phi
Definition: TrackingPerformanceEvaluationModule.h:125
genfit::TrackPoint
Object containing AbsMeasurement and AbsFitterInfo objects.
Definition: TrackPoint.h:46
Belle2::TrackingPerformanceEvaluationModule::m_h3_MCRecoTrack_minus
TH3F * m_h3_MCRecoTrack_minus
efficiency
Definition: TrackingPerformanceEvaluationModule.h:185
Belle2::TrackingPerformanceEvaluationModule::m_h2_d0errVSpt_wfpxd
TH2F * m_h2_d0errVSpt_wfpxd
error
Definition: TrackingPerformanceEvaluationModule.h:158
Belle2::MCParticleInfo::getZ
double getZ()
Getter for z component of vertex.
Definition: MCParticleInfo.h:59
Belle2::WireID
Class to identify a wire inside the CDC.
Definition: WireID.h:44
Belle2::SVDRecoHit
SVDRecoHit - an extended form of SVDHit containing geometry information.
Definition: SVDRecoHit.h:57
Belle2::MCParticleInfo
This struct is used by the TrackingPerformanceEvaluation Module to save information of reconstructed ...
Definition: MCParticleInfo.h:34
Belle2::TrackingPerformanceEvaluationModule::m_h1_cotTheta_err
TH1F * m_h1_cotTheta_err
error
Definition: TrackingPerformanceEvaluationModule.h:104
Belle2::VxdID
Class to uniquely identify a any structure of the PXD and SVD.
Definition: VxdID.h:43
Belle2::TrackingPerformanceEvaluationModule::m_h1_r_res
TH1F * m_h1_r_res
R residual (in cylindrical coordinates)
Definition: TrackingPerformanceEvaluationModule.h:121
Belle2::TrackingPerformanceEvaluationModule::m_h1_pValue
TH1F * m_h1_pValue
p val
Definition: TrackingPerformanceEvaluationModule.h:148
Belle2::TrackingPerformanceEvaluationModule::m_h3_MCParticlesPerTrack
TH3F * m_h3_MCParticlesPerTrack
purity
Definition: TrackingPerformanceEvaluationModule.h:189
Belle2::SVDRecoHit2D::getU
float getU() const
Get u coordinate.
Definition: SVDRecoHit2D.h:119
Belle2::HitPatternVXD::getNPXDHits
unsigned short getNPXDHits() const
Get total number of hits in the PXD.
Definition: HitPatternVXD.cc:147
Belle2::TrackingPerformanceEvaluationModule::m_TracksName
std::string m_TracksName
Track StoreArray name.
Definition: TrackingPerformanceEvaluationModule.h:85
Belle2::PXDTrueHit
Class PXDTrueHit - Records of tracks that either enter or leave the sensitive volume.
Definition: PXDTrueHit.h:35
Belle2::TrackingPerformanceEvaluationModule::m_h2_TrackPointFitWeightCDC
TH2F * m_h2_TrackPointFitWeightCDC
TP.
Definition: TrackingPerformanceEvaluationModule.h:146
Belle2::BFieldManager::getField
static void getField(const double *pos, double *field)
return the magnetic field at a given position.
Definition: BFieldManager.h:110
Belle2::TrackingPerformanceEvaluationModule::m_multiplicityFittedTracks
TH1F * m_multiplicityFittedTracks
number of fitted tracks per MCParticles
Definition: TrackingPerformanceEvaluationModule.h:93
Belle2::TrackingPerformanceEvaluationModule::m_h1_nCDChitsPR
TH1F * m_h1_nCDChitsPR
PR.
Definition: TrackingPerformanceEvaluationModule.h:141
Belle2::TrackingPerformanceEvaluationModule::m_h2_TrackPointFitWeightVXD
TH2F * m_h2_TrackPointFitWeightVXD
TP.
Definition: TrackingPerformanceEvaluationModule.h:145
Belle2::Track::getTrackFitResult
const TrackFitResult * getTrackFitResult(const Const::ChargedStable &chargedStable) const
Access to TrackFitResults.
Definition: Track.cc:19
Belle2::TrackingPerformanceEvaluationModule::m_h3_MCParticleswPXDHits
TH3F * m_h3_MCParticleswPXDHits
efficiency
Definition: TrackingPerformanceEvaluationModule.h:171
Belle2::Const::invalidParticle
static const ParticleType invalidParticle
Invalid particle, used internally.
Definition: Const.h:554
Belle2::TrackFitResult::getMomentum
TVector3 getMomentum() const
Getter for vector of momentum at closest approach of track in r/phi projection.
Definition: TrackFitResult.h:116
Belle2::TrackingPerformanceEvaluationModule::m_multiplicityTracks
TH1F * m_multiplicityTracks
number of tracks per MCParticles
Definition: TrackingPerformanceEvaluationModule.h:90
Belle2::TrackingPerformanceEvaluationModule::m_RecoTracksName
std::string m_RecoTracksName
RecoTrack StoreArray name.
Definition: TrackingPerformanceEvaluationModule.h:84
Belle2::TrackingPerformanceEvaluationModule::m_h1_y_res
TH1F * m_h1_y_res
y residual
Definition: TrackingPerformanceEvaluationModule.h:119
Belle2::TrackingPerformanceEvaluationModule::m_h3_MCParticle
TH3F * m_h3_MCParticle
efficiency
Definition: TrackingPerformanceEvaluationModule.h:170
REG_MODULE
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition: Module.h:652
Belle2::TrackingPerformanceEvaluationModule::m_h3_MCRecoTrack_plus
TH3F * m_h3_MCRecoTrack_plus
efficiency
Definition: TrackingPerformanceEvaluationModule.h:181
Belle2::TrackFitResult::getOmega
double getOmega() const
Getter for omega.
Definition: TrackFitResult.h:194
Belle2::PerformanceEvaluationBaseClass::m_histoList_others
TList * m_histoList_others
List of other performance-evaluation histograms.
Definition: PerformanceEvaluationBaseClass.h:50
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::TrackingPerformanceEvaluationModule::m_h1_nVXDhitsUsed
TH1F * m_h1_nVXDhitsUsed
hits used
Definition: TrackingPerformanceEvaluationModule.h:140
Belle2::TrackFitResult::getPValue
double getPValue() const
Getter for Chi2 Probability of the track fit.
Definition: TrackFitResult.h:163
Belle2::PerformanceEvaluationBaseClass::m_histoList_efficiency
TList * m_histoList_efficiency
List of efficiency histograms.
Definition: PerformanceEvaluationBaseClass.h:48
Belle2::MCParticleInfo::getP
double getP()
Getter for magnitut of momentum.
Definition: MCParticleInfo.h:50
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::createHistogramsRatio
TH1F * createHistogramsRatio(const char *name, const char *title, TH1 *hNum, TH1 *hDen, bool isEffPlot, int axisRef)
Make a new 1D histogram from the ratio of two others and add it to the TList of 1D-histograms.
Definition: PerformanceEvaluationBaseClass.cc:166
Belle2::TrackingPerformanceEvaluationModule::m_h1_z0_err
TH1F * m_h1_z0_err
error
Definition: TrackingPerformanceEvaluationModule.h:103
Belle2::TrackingPerformanceEvaluationModule::m_h1_phi_err
TH1F * m_h1_phi_err
error
Definition: TrackingPerformanceEvaluationModule.h:101
genfit::TrackPoint::getKalmanFitterInfo
KalmanFitterInfo * getKalmanFitterInfo(const AbsTrackRep *rep=nullptr) const
Helper to avoid casting.
Definition: TrackPoint.cc:180
Belle2::TrackingPerformanceEvaluationModule::m_h3_TracksPerMCRecoTrack
TH3F * m_h3_TracksPerMCRecoTrack
efficiency
Definition: TrackingPerformanceEvaluationModule.h:177
Belle2::TrackingPerformanceEvaluationModule::fillTrackParams1DHistograms
void fillTrackParams1DHistograms(const TrackFitResult *fitResult, MCParticleInfo mcParticleInfo)
fills err, resid and pull TH1F for each of the 5 track parameters
Definition: TrackingPerformanceEvaluationModule.cc:724
Belle2::MCParticleInfo::getPy
double getPy()
Getter for y component of momentum.
Definition: MCParticleInfo.h:44
Belle2::TrackingPerformanceEvaluationModule::m_h1_rtot_res
TH1F * m_h1_rtot_res
r residual (3D distance)
Definition: TrackingPerformanceEvaluationModule.h:122
Belle2::Const::chargedStableSet
static const ParticleSet chargedStableSet
set of charged stable particles
Definition: Const.h:494
Belle2::VXD::SensorInfoBase
Base class to provide Sensor Information for PXD and SVD.
Definition: SensorInfoBase.h:40
Belle2::MCParticle::hasSeenInDetector
bool hasSeenInDetector(Const::DetectorSet set) const
Return if the seen-in flag for a specific subdetector is set or not.
Definition: MCParticle.h:318
Belle2::TrackingPerformanceEvaluationModule::m_h2_OmegaerrOmegaVSpt
TH2F * m_h2_OmegaerrOmegaVSpt
error
Definition: TrackingPerformanceEvaluationModule.h:150
Belle2::TrackingPerformanceEvaluationModule::m_h1_phi_res
TH1F * m_h1_phi_res
error
Definition: TrackingPerformanceEvaluationModule.h:107
Belle2::TrackingPerformanceEvaluationModule::m_h1_py_res
TH1F * m_h1_py_res
py residual
Definition: TrackingPerformanceEvaluationModule.h:113
Belle2::PXDRecoHit
PXDRecoHit - an extended form of PXDCluster containing geometry information.
Definition: PXDRecoHit.h:61
Belle2::TrackingPerformanceEvaluationModule::m_h1_z_res
TH1F * m_h1_z_res
z residual
Definition: TrackingPerformanceEvaluationModule.h:120
Belle2::TrackingPerformanceEvaluationModule::m_multiplicityMCRecoTracks
TH1F * m_multiplicityMCRecoTracks
number of MCRecoTracks per MCParticles
Definition: TrackingPerformanceEvaluationModule.h:92
Belle2::TrackingPerformanceEvaluationModule::m_h2_d0errVSpt_wtpxd
TH2F * m_h2_d0errVSpt_wtpxd
error
Definition: TrackingPerformanceEvaluationModule.h:157
Belle2::TrackFitResult
Values of the result of a track fit with a given particle hypothesis.
Definition: TrackFitResult.h:59
Belle2::TrackingPerformanceEvaluationModule::endRun
void endRun() override
This method is called if the current run ends.
Definition: TrackingPerformanceEvaluationModule.cc:622
Belle2::TrackFitResult::getPosition
TVector3 getPosition() const
Getter for vector of position at closest approach of track in r/phi projection.
Definition: TrackFitResult.h:109
Belle2::MCParticleInfo::getCharge
double getCharge()
Getter for electric charge of particle.
Definition: MCParticleInfo.h:67
Belle2::MCParticleInfo::getY
double getY()
Getter for y component of vertex.
Definition: MCParticleInfo.h:57
Belle2::MCParticleInfo::getZ0
double getZ0()
Getter for Z0.
Definition: MCParticleInfo.cc:76
Belle2::PerformanceEvaluationBaseClass::m_histoList_fit
TList * m_histoList_fit
List of track-fit histograms.
Definition: PerformanceEvaluationBaseClass.h:47
Belle2::TrackingPerformanceEvaluationModule::m_h1_z0_res
TH1F * m_h1_z0_res
error
Definition: TrackingPerformanceEvaluationModule.h:109
Belle2::TrackFitResult::getHitPatternCDC
HitPatternCDC getHitPatternCDC() const
Getter for the hit pattern in the CDC;.
Definition: TrackFitResult.cc:120
Belle2::TrackFitResult::getZ0
double getZ0() const
Getter for z0.
Definition: TrackFitResult.h:200
Belle2::TrackingPerformanceEvaluationModule::m_h2_d0errVSpt
TH2F * m_h2_d0errVSpt
error
Definition: TrackingPerformanceEvaluationModule.h:161
Belle2::TrackingPerformanceEvaluationModule::m_h2_d0errphi0err_rz
TH2F * m_h2_d0errphi0err_rz
error
Definition: TrackingPerformanceEvaluationModule.h:133
Belle2::MCParticleInfo::getPhi
double getPhi()
Getter for Phi.
Definition: MCParticleInfo.cc:47
Belle2::CDCRecoHit
This class is used to transfer CDC information to the track fit.
Definition: CDCRecoHit.h:43
Belle2::B2Vector3< double >
Belle2::TrackingPerformanceEvaluationModule::m_h1_HitsRecoTrackPerMCRecoTrack
TH1F * m_h1_HitsRecoTrackPerMCRecoTrack
hits
Definition: TrackingPerformanceEvaluationModule.h:167
Belle2::Module
Base class for Modules.
Definition: Module.h:74
Belle2::TrackingPerformanceEvaluationModule::m_h1_d0_res
TH1F * m_h1_d0_res
error
Definition: TrackingPerformanceEvaluationModule.h:106
Belle2::MCParticle::hasStatus
bool hasStatus(unsigned short int bitmask) const
Return if specific status bit is set.
Definition: MCParticle.h:140
Belle2::TrackingPerformanceEvaluationModule::m_multiplicityRecoTracksPerMCRT
TH1F * m_multiplicityRecoTracksPerMCRT
number of RecoTracks per MCRecoTracks
Definition: TrackingPerformanceEvaluationModule.h:96
genfit::AbsMeasurement
Contains the measurement and covariance in raw detector coordinates.
Definition: AbsMeasurement.h:42
Belle2::TrackFitResult::getCotTheta
double getCotTheta() const
Getter for tanLambda with CDF naming convention.
Definition: TrackFitResult.h:209
Belle2::RecoTrack
This is the Reconstruction Event-Data Model Track.
Definition: RecoTrack.h:78
Belle2::Unit::T
static const double T
[tesla]
Definition: Unit.h:130
Belle2::HitPatternVXD::getSVDLayer
std::pair< const unsigned short, const unsigned short > getSVDLayer(const unsigned short layerId) const
Get the number of hits in a specific layer of the SVD.
Definition: HitPatternVXD.cc:65
Belle2::MCParticleInfo::getPt
double getPt()
Getter for transverse momentum.
Definition: MCParticleInfo.h:48
Belle2::CDCRecoHit::getWireID
WireID getWireID() const
Getter for WireID object.
Definition: CDCRecoHit.h:60
Belle2::TrackingPerformanceEvaluationModule::m_h1_d0_err
TH1F * m_h1_d0_err
error
Definition: TrackingPerformanceEvaluationModule.h:100
Belle2::RelationVector
Class for type safe access to objects that are referred to in relations.
Definition: DataStore.h:38
genfit::KalmanFitterInfo
Collects information needed and produced by a AbsKalmanFitter implementations and is specific to one ...
Definition: KalmanFitterInfo.h:44
Belle2::SVDRecoHit2D::getSensorID
VxdID getSensorID() const
Get the compact ID.
Definition: SVDRecoHit2D.h:109
Belle2::TrackingPerformanceEvaluationModule::fillTrackErrParams2DHistograms
void fillTrackErrParams2DHistograms(const TrackFitResult *fitResult)
fills TH2F
Definition: TrackingPerformanceEvaluationModule.cc:796
Belle2::TrackingPerformanceEvaluationModule::m_h1_d0_pll
TH1F * m_h1_d0_pll
pull distribution d0
Definition: TrackingPerformanceEvaluationModule.h:124
Belle2::PerformanceEvaluationBaseClass::m_rootFileName
std::string m_rootFileName
root file name
Definition: PerformanceEvaluationBaseClass.h:139
Belle2::TrackingPerformanceEvaluationModule::addMoreInefficiencyPlots
void addMoreInefficiencyPlots(TList *histoList)
add inefficiency plots
Definition: TrackingPerformanceEvaluationModule.cc:1013
Belle2::TrackFitResult::getParticleType
Const::ParticleType getParticleType() const
Getter for ParticleType of the mass hypothesis of the track fit.
Definition: TrackFitResult.h:154
Belle2::PerformanceEvaluationBaseClass::m_histoList_trkQuality
TList * m_histoList_trkQuality
List of track-quality histograms.
Definition: PerformanceEvaluationBaseClass.h:44
Belle2::Const::ParticleSet::find
const ParticleType & find(int pdg) const
Returns particle in set with given PDG code, or invalidParticle if not found.
Definition: Const.h:447
Belle2::TrackingPerformanceEvaluationModule::m_h1_cotTheta_pll
TH1F * m_h1_cotTheta_pll
pull distribution cotTheta
Definition: TrackingPerformanceEvaluationModule.h:128
Belle2::TrackingPerformanceEvaluationModule::m_h1_z0_pll
TH1F * m_h1_z0_pll
pull distribution z0
Definition: TrackingPerformanceEvaluationModule.h:127
Belle2::VXD::GeoCache::getInstance
static GeoCache & getInstance()
Return a reference to the singleton instance.
Definition: GeoCache.cc:215
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::TrackingPerformanceEvaluationModule::m_ParticleHypothesis
int m_ParticleHypothesis
Particle Hypothesis for the track fit (default: 211)
Definition: TrackingPerformanceEvaluationModule.h:86
Belle2::SVDRecoHit2D::getV
float getV() const
Get v coordinate.
Definition: SVDRecoHit2D.h:121
Belle2::PerformanceEvaluationBaseClass::addPurityPlots
void addPurityPlots(TList *graphList=NULL, TH3F *h3_xPerMCParticle=NULL, TH3F *h3_MCParticle=NULL)
Create pt-, theta- and phi-purity 1D histograms and add them to the TList of 1D-histograms.
Definition: PerformanceEvaluationBaseClass.cc:314
Belle2::TrackingPerformanceEvaluationModule::m_MCParticlesName
std::string m_MCParticlesName
MCParticle StoreArray name.
Definition: TrackingPerformanceEvaluationModule.h:82
Belle2::TrackingPerformanceEvaluationModule::m_h1_p_res
TH1F * m_h1_p_res
p residual
Definition: TrackingPerformanceEvaluationModule.h:115
Belle2::MCParticleInfo::getPphi
double getPphi()
Getter for phi of momentum vector.
Definition: MCParticleInfo.h:64
Belle2::TrackingPerformanceEvaluationModule::m_h2_z0errVSpt_wtpxd
TH2F * m_h2_z0errVSpt_wtpxd
error
Definition: TrackingPerformanceEvaluationModule.h:152
Belle2::TrackingPerformanceEvaluationModule::m_h1_omega_err
TH1F * m_h1_omega_err
error
Definition: TrackingPerformanceEvaluationModule.h:102
Belle2::SVDRecoHit::isU
bool isU() const
Is the coordinate u or v?
Definition: SVDRecoHit.h:100
Belle2::TrackingPerformanceEvaluationModule::m_h3_Tracks
TH3F * m_h3_Tracks
purity
Definition: TrackingPerformanceEvaluationModule.h:190
Belle2::MCParticleInfo::getOmega
double getOmega()
Getter for Omega.
Definition: MCParticleInfo.cc:69
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::TrackingPerformanceEvaluationModule::m_h3_TracksPerMCParticle_minus
TH3F * m_h3_TracksPerMCParticle_minus
efficiency
Definition: TrackingPerformanceEvaluationModule.h:184
Belle2::TrackingPerformanceEvaluationModule::m_h3_RecoTrackswPXDHitsPerMCParticlewPXDHits
TH3F * m_h3_RecoTrackswPXDHitsPerMCParticlewPXDHits
efficiency
Definition: TrackingPerformanceEvaluationModule.h:175
Belle2::TrackingPerformanceEvaluationModule::m_multiplicityFittedTracksPerMCRT
TH1F * m_multiplicityFittedTracksPerMCRT
number of fitted tracks per MCRecoTrack
Definition: TrackingPerformanceEvaluationModule.h:94
Belle2::TrackingPerformanceEvaluationModule::m_h1_pz_res
TH1F * m_h1_pz_res
pz residual
Definition: TrackingPerformanceEvaluationModule.h:114
Belle2::TrackingPerformanceEvaluationModule::m_h3_MCParticle_minus
TH3F * m_h3_MCParticle_minus
efficiency
Definition: TrackingPerformanceEvaluationModule.h:183
Belle2::TrackingPerformanceEvaluationModule::m_h2_d0errVSpt_wpxd
TH2F * m_h2_d0errVSpt_wpxd
error
Definition: TrackingPerformanceEvaluationModule.h:159
Belle2::PXDRecoHit::getU
float getU() const
Get u coordinate.
Definition: PXDRecoHit.h:119
Belle2::PerformanceEvaluationBaseClass::m_rootFilePtr
TFile * m_rootFilePtr
pointer at root file used for storing histograms
Definition: PerformanceEvaluationBaseClass.h:142
Belle2::TrackingPerformanceEvaluationModule::m_multiplicityMCParticlesPerTrack
TH1F * m_multiplicityMCParticlesPerTrack
number of MCParticles per fitted Track
Definition: TrackingPerformanceEvaluationModule.h:95
Belle2::TrackingPerformanceEvaluationModule::m_h2_z0errVSpt_wfpxd
TH2F * m_h2_z0errVSpt_wfpxd
error
Definition: TrackingPerformanceEvaluationModule.h:153
Belle2::PXDRecoHit::getSensorID
VxdID getSensorID() const
Get the compact ID.
Definition: PXDRecoHit.h:111
Belle2::PXDCluster
The PXD Cluster class This class stores all information about reconstructed PXD clusters The position...
Definition: PXDCluster.h:41
Belle2::TrackingPerformanceEvaluationModule::m_h3_TracksPerMCRecoTrack_plus
TH3F * m_h3_TracksPerMCRecoTrack_plus
efficiency
Definition: TrackingPerformanceEvaluationModule.h:182
Belle2::TrackingPerformanceEvaluationModule::event
void event() override
This method is the core of the module.
Definition: TrackingPerformanceEvaluationModule.cc:359
Belle2::TrackingPerformanceEvaluationModule::m_h2_VXDhitsPR_xy
TH2F * m_h2_VXDhitsPR_xy
PR.
Definition: TrackingPerformanceEvaluationModule.h:136
Belle2::TrackingPerformanceEvaluationModule::m_h2_VXDhitsPR_rz
TH2F * m_h2_VXDhitsPR_rz
PR.
Definition: TrackingPerformanceEvaluationModule.h:137
Belle2::RelationsInterface::getArrayIndex
int getArrayIndex() const
Returns this object's array index (in StoreArray), or -1 if not found.
Definition: RelationsObject.h:387
Belle2::MCParticleInfo::getX
double getX()
Getter for x component of vertex.
Definition: MCParticleInfo.h:55
Belle2::PerformanceEvaluationBaseClass::m_histoList
TList * m_histoList
List of performance-evaluation histograms.
Definition: PerformanceEvaluationBaseClass.h:40
Belle2::TrackingPerformanceEvaluationModule::m_h1_nVXDhitsPR
TH1F * m_h1_nVXDhitsPR
PR.
Definition: TrackingPerformanceEvaluationModule.h:138
Belle2::VXD::SensorInfoBase::pointToGlobal
TVector3 pointToGlobal(const TVector3 &local, bool reco=false) const
Convert a point from local to global coordinates.
Definition: SensorInfoBase.h:373
Belle2::TrackingPerformanceEvaluationModule::fillHitsUsedInTrackFitHistograms
void fillHitsUsedInTrackFitHistograms(const Track &track)
fill TH2F
Definition: TrackingPerformanceEvaluationModule.cc:831
Belle2::PXDRecoHit::getCluster
const PXDCluster * getCluster() const
Get pointer to the Cluster used when creating this RecoHit, can be NULL if created from something els...
Definition: PXDRecoHit.h:116
Belle2::TrackingPerformanceEvaluationModule::m_h3_TracksPerMCParticle_plus
TH3F * m_h3_TracksPerMCParticle_plus
efficiency
Definition: TrackingPerformanceEvaluationModule.h:180
Belle2::TrackFitResult::getHitPatternVXD
HitPatternVXD getHitPatternVXD() const
Getter for the hit pattern in the VXD;.
Definition: TrackFitResult.cc:125
Belle2::MCParticleInfo::getPx
double getPx()
Getter for x component of momentum.
Definition: MCParticleInfo.h:42
Belle2::SVDRecoHit2D
SVDRecoHit - an extended form of SVDHit containing geometry information.
Definition: SVDRecoHit2D.h:56
Belle2::TrackingPerformanceEvaluationModule::m_h2_chargeVSchargeMC
TH2F * m_h2_chargeVSchargeMC
charge comparison
Definition: TrackingPerformanceEvaluationModule.h:164
Belle2::TrackingPerformanceEvaluationModule::initialize
void initialize() override
Initialize the Module.
Definition: TrackingPerformanceEvaluationModule.cc:71
Belle2::TrackingPerformanceEvaluationModule::m_h2_z0errVSpt
TH2F * m_h2_z0errVSpt
error
Definition: TrackingPerformanceEvaluationModule.h:156
Belle2::PerformanceEvaluationBaseClass::m_histoList_pr
TList * m_histoList_pr
List of pattern-recognition histograms.
Definition: PerformanceEvaluationBaseClass.h:46
Belle2::TrackingPerformanceEvaluationModule::m_h1_nCDChitsWeighted
TH1F * m_h1_nCDChitsWeighted
weighted
Definition: TrackingPerformanceEvaluationModule.h:142
Belle2::TrackingPerformanceEvaluationModule::m_h1_omega_res
TH1F * m_h1_omega_res
error
Definition: TrackingPerformanceEvaluationModule.h:108
Belle2::TrackingPerformanceEvaluationModule::m_h3_TracksPerMCRecoTrack_minus
TH3F * m_h3_TracksPerMCRecoTrack_minus
efficiency
Definition: TrackingPerformanceEvaluationModule.h:186
Belle2::TrackFitResult::getCovariance5
TMatrixDSym getCovariance5() const
Getter for covariance matrix of perigee parameters in matrix form.
Definition: TrackFitResult.cc:106
Belle2::VXD::GeoCache
Class to faciliate easy access to sensor information of the VXD like coordinate transformations or pi...
Definition: GeoCache.h:41
Belle2::TrackingPerformanceEvaluationModule::m_multiplicityMCRecoTracksPerRT
TH1F * m_multiplicityMCRecoTracksPerRT
number of MCRecoTracks per RecoTracks
Definition: TrackingPerformanceEvaluationModule.h:97
Belle2::TrackingPerformanceEvaluationModule::m_h3_RecoTrackswPXDHitsPerMCParticle
TH3F * m_h3_RecoTrackswPXDHitsPerMCParticle
efficiency
Definition: TrackingPerformanceEvaluationModule.h:174
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::Const::ChargedStable
Provides a type-safe way to pass members of the chargedStableSet set.
Definition: Const.h:465
Belle2::TrackingPerformanceEvaluationModule::m_multiplicityRecoTracks
TH1F * m_multiplicityRecoTracks
number of recoTracks per MCParticles
Definition: TrackingPerformanceEvaluationModule.h:91
Belle2::TrackingPerformanceEvaluationModule::addMoreEfficiencyPlots
void addMoreEfficiencyPlots(TList *histoList)
add efficiency plots
Definition: TrackingPerformanceEvaluationModule.cc:1044
Belle2::HitPatternCDC::hasLayer
bool hasLayer(const unsigned short layer) const
Getter for single layer.
Definition: HitPatternCDC.cc:92
Belle2::TrackingPerformanceEvaluationModule::m_h2_d0errVSpt_wopxd
TH2F * m_h2_d0errVSpt_wopxd
error
Definition: TrackingPerformanceEvaluationModule.h:160
Belle2::TrackingPerformanceEvaluationModule::m_h2_z0errVSpt_wpxd
TH2F * m_h2_z0errVSpt_wpxd
error
Definition: TrackingPerformanceEvaluationModule.h:154
Belle2::TrackingPerformanceEvaluationModule::m_h1_pt_res
TH1F * m_h1_pt_res
pt residual
Definition: TrackingPerformanceEvaluationModule.h:116
Belle2::TrackingPerformanceEvaluationModule::m_h3_TracksPerMCParticle
TH3F * m_h3_TracksPerMCParticle
efficiency
Definition: TrackingPerformanceEvaluationModule.h:172
Belle2::TrackingPerformanceEvaluationModule::m_h1_px_res
TH1F * m_h1_px_res
px residual
Definition: TrackingPerformanceEvaluationModule.h:112
Belle2::PXDRecoHit::getV
float getV() const
Get v coordinate.
Definition: PXDRecoHit.h:121
Belle2::TrackingPerformanceEvaluationModule::m_h2_z0errVSpt_wopxd
TH2F * m_h2_z0errVSpt_wopxd
error
Definition: TrackingPerformanceEvaluationModule.h:155
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::TrackingPerformanceEvaluationModule::m_h3_TrackswPXDHitsPerMCParticle
TH3F * m_h3_TrackswPXDHitsPerMCParticle
efficiency
Definition: TrackingPerformanceEvaluationModule.h:173
Belle2::Track
Class that bundles various TrackFitResults.
Definition: Track.h:35
Belle2::TrackingPerformanceEvaluationModule::m_h2_z0errcotThetaerr_xy
TH2F * m_h2_z0errcotThetaerr_xy
error
Definition: TrackingPerformanceEvaluationModule.h:134
Belle2::MCParticle
A Class to store the Monte Carlo particle information.
Definition: MCParticle.h:43
Belle2::PerformanceEvaluationBaseClass::m_histoList_evtCharacterization
TList * m_histoList_evtCharacterization
List of event-characterization histograms.
Definition: PerformanceEvaluationBaseClass.h:43
genfit::KalmanFitterInfo::getWeights
std::vector< double > getWeights() const
Get weights of measurements.
Definition: KalmanFitterInfo.cc:168
Belle2::StoreArray< MCParticle >
Belle2::TrackingPerformanceEvaluationModule::m_h2_d0errMSVSpt
TH2F * m_h2_d0errMSVSpt
error
Definition: TrackingPerformanceEvaluationModule.h:162
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::TrackingPerformanceEvaluationModule::isTraceable
bool isTraceable(const MCParticle &the_mcParticle)
is traceable
Definition: TrackingPerformanceEvaluationModule.cc:1168
Belle2::TrackingPerformanceEvaluationModule::m_h1_nHitDetID
TH1F * m_h1_nHitDetID
det ID
Definition: TrackingPerformanceEvaluationModule.h:144
Belle2::TrackingPerformanceEvaluationModule::m_h1_HitsMCRecoTrack
TH1F * m_h1_HitsMCRecoTrack
hits
Definition: TrackingPerformanceEvaluationModule.h:168
Belle2::SVDRecoHit::getPosition
float getPosition() const
Get coordinate.
Definition: SVDRecoHit.h:103
Belle2::MCParticleInfo::getPz
double getPz()
Getter for z component of momentum.
Definition: MCParticleInfo.h:46
Belle2::MCParticleInfo::getD0
double getD0()
Getter for D0.
Definition: MCParticleInfo.cc:27
Belle2::TrackingPerformanceEvaluationModule::m_h1_omega_pll
TH1F * m_h1_omega_pll
pull distribution omega
Definition: TrackingPerformanceEvaluationModule.h:126
Belle2::TrackingPerformanceEvaluationModule::m_MCRecoTracksName
std::string m_MCRecoTracksName
MCRecoTrack StoreArray name.
Definition: TrackingPerformanceEvaluationModule.h:83
Belle2::WireID::getICLayer
unsigned short getICLayer() const
Getter for continuous layer numbering.
Definition: WireID.cc:26
Belle2::TrackingPerformanceEvaluationModule::m_h1_nCDChitsUsed
TH1F * m_h1_nCDChitsUsed
used
Definition: TrackingPerformanceEvaluationModule.h:143
Belle2::TrackingPerformanceEvaluationModule::m_h3_MCParticle_plus
TH3F * m_h3_MCParticle_plus
efficiency
Definition: TrackingPerformanceEvaluationModule.h:179
Belle2::SVDRecoHit::getSensorID
VxdID getSensorID() const
Get the compact ID.
Definition: SVDRecoHit.h:91
Belle2::PerformanceEvaluationBaseClass::m_histoList_firstHit
TList * m_histoList_firstHit
List of first-hit-position histograms.
Definition: PerformanceEvaluationBaseClass.h:45
Belle2::TrackingPerformanceEvaluationModule::m_h1_x_res
TH1F * m_h1_x_res
x residual
Definition: TrackingPerformanceEvaluationModule.h:118
Belle2::TrackFitResult::getChargeSign
short getChargeSign() const
Return track charge (1 or -1).
Definition: TrackFitResult.h:160
Belle2::MCParticle::c_PrimaryParticle
@ c_PrimaryParticle
bit 0: Particle is primary particle.
Definition: MCParticle.h:58
Belle2::TrackingPerformanceEvaluationModule::terminate
void terminate() override
This method is called at the end of the event processing.
Definition: TrackingPerformanceEvaluationModule.cc:676
Belle2::TrackFitResult::getPhi
double getPhi() const
Getter for phi0 with CDF naming convention.
Definition: TrackFitResult.h:187
Belle2::PerformanceEvaluationBaseClass::m_histoList_purity
TList * m_histoList_purity
List of purity histograms.
Definition: PerformanceEvaluationBaseClass.h:49
Belle2::TrackingPerformanceEvaluationModule::m_h1_nVXDhitsWeighted
TH1F * m_h1_nVXDhitsWeighted
weighted
Definition: TrackingPerformanceEvaluationModule.h:139
Belle2::TrackFitResult::getD0
double getD0() const
Getter for d0.
Definition: TrackFitResult.h:178
Belle2::MCParticleInfo::getCotTheta
double getCotTheta()
Getter for Theta.
Definition: MCParticleInfo.cc:82
Belle2::Const::ParticleType::getMass
double getMass() const
Particle mass.
Definition: UnitConst.cc:310
Belle2::TrackingPerformanceEvaluationModule::m_h3_MCRecoTrack
TH3F * m_h3_MCRecoTrack
efficiency
Definition: TrackingPerformanceEvaluationModule.h:176