CbmRoot
CbmAnaConversionReco.cxx
Go to the documentation of this file.
1 
11 #include "CbmAnaConversionReco.h"
12 
13 // standard includes
14 #include <iostream>
15 
16 // includes from ROOT
17 #include "TRandom3.h"
18 
19 // included from CbmRoot
20 #include "CbmGlobalTrack.h"
21 #include "CbmMCTrack.h"
22 #include "CbmRichElectronIdAnn.h"
23 #include "CbmRichHit.h"
24 #include "CbmRichRing.h"
25 #include "CbmRichRingLight.h"
26 #include "CbmRichUtil.h"
27 #include "FairRootManager.h"
28 
29 #include "CbmAnaConversionCutSettings.h"
30 #include "CbmAnaConversionGlobalFunctions.h"
31 
32 #define M2E 2.6112004954086e-7
33 using namespace std;
34 
35 
36 CbmAnaConversionReco::CbmAnaConversionReco()
37  : fMcTracks(NULL)
38  , fGlobalTracks(NULL)
39  , fRichRings(NULL)
40  , fRichHits(NULL)
41  , fMCTracklist_all()
42  , fRecoTracklistEPEM()
43  , fRecoTracklistEPEM_ids()
44  , fRecoTracklistEPEM_chi()
45  , fRecoTracklistEPEM_gtid()
46  , fRecoMomentum()
47  , fRecoRefittedMomentum()
48  , fHistoList_MC()
49  , fHistoList_reco()
50  , fHistoList_reco_mom()
51  , fHistoList_gg()
52  , fHistoList_gee()
53  , fHistoList_eeee()
54  , fHistoList_all()
55  , fHistoList_eta()
56  , fhInvariantMass_MC_all(NULL)
57  , fhInvariantMass_MC_pi0(NULL)
58  , fhInvariantMass_MC_pi0_epem(NULL)
59  , fhInvariantMass_MC_pi0_gepem(NULL)
60  , fhInvariantMass_MC_pi0_gg(NULL)
61  , fhInvariantMass_MC_eta(NULL)
62  , fhInvariantMass_MC_etaPrime(NULL)
63  , fhMC_electrons_theta(NULL)
64  , fhMC_electrons_p(NULL)
65  , fhMC_electrons_theta_vs_p(NULL)
66  , fhEta_openingAngleGG(NULL)
67  , fhMC_grandmotherPDGs(NULL)
68  , fhInvariantMassReco_pi0(NULL)
69  , fhMCtest(NULL)
70  , fhEPEM_invmass_gg_mc(NULL)
71  , fhEPEM_invmass_gg_refitted(NULL)
72  , fhEPEM_invmass_gee_mc(NULL)
73  , fhEPEM_invmass_gee_refitted(NULL)
74  , fhEPEM_invmass_eeee_mc(NULL)
75  , fhEPEM_invmass_eeee_refitted(NULL)
76  , fhEPEM_invmass_all_mc(NULL)
77  , fhEPEM_invmass_all_refitted(NULL)
78  , fhEPEM_openingAngle_gg_mc(NULL)
79  , fhEPEM_openingAngle_gg_refitted(NULL)
80  , fhEPEM_openingAngle_gee_mc(NULL)
81  , fhEPEM_openingAngle_gee_refitted(NULL)
82  , fhEPEM_openingAngle_gee_mc_dalitz(NULL)
83  , fhEPEM_openingAngle_gee_refitted_dalitz(NULL)
84  , fhEPEM_openingAngle_vs_pt_gg_mc(NULL)
85  , fhEPEM_openingAngle_vs_pt_gg_reco(NULL)
86  , fhEPEM_openingAngle_betweenGammas_mc(NULL)
87  , fhEPEM_openingAngle_betweenGammas_reco(NULL)
88  , fhPi0_pt_vs_rap_gg(NULL)
89  , fhPi0_pt_vs_rap_gee(NULL)
90  , fhPi0_pt_vs_rap_all(NULL)
91  , fhPi0_pt_gg(NULL)
92  , fhPi0_pt_gee(NULL)
93  , fhPi0_pt_all(NULL)
94  , fhEPEM_efficiencyCuts(NULL)
95  , fhEPEM_efficiencyCuts2(NULL)
96  , fhEPEM_pi0_nofLeptons_ann(NULL)
97  , fhEPEM_pi0_ANNvalues_noCuts(NULL)
98  , fhEPEM_pi0_ANNvalues_angleCut(NULL)
99  , fhEPEM_pi0_ANNefficiencies(NULL)
100  , fhEPEM_rap_vs_chi(NULL)
101  , fhEPEM_rap_vs_invmass(NULL)
102  , fhInvMass_EPEM_mc(NULL)
103  , fhInvMass_EPEM_stsMomVec(NULL)
104  , fhInvMass_EPEM_refitted(NULL)
105  , fhInvMass_EPEM_error_stsMomVec(NULL)
106  , fhInvMass_EPEM_error_refitted(NULL)
107  , fhInvMass_EPEM_openingAngleRef(NULL)
108  , fhUsedMomenta_stsMomVec(NULL)
109  , fhUsedMomenta_mc(NULL)
110  , fhUsedMomenta_error_stsMomVec(NULL)
111  , fhUsedMomenta_error_refitted(NULL)
112  , fhUsedMomenta_errorX_stsMomVec(NULL)
113  , fhUsedMomenta_vsX_stsMomVec(NULL)
114  , fhUsedMomenta_errorY_stsMomVec(NULL)
115  , fhUsedMomenta_vsY_stsMomVec(NULL)
116  , fhUsedMomenta_errorZ_stsMomVec(NULL)
117  , fhUsedMomenta_vsZ_stsMomVec(NULL)
118  , fhUsedMomenta_errorX_refitted(NULL)
119  , fhUsedMomenta_vsX_refitted(NULL)
120  , fhUsedMomenta_errorY_refitted(NULL)
121  , fhUsedMomenta_vsY_refitted(NULL)
122  , fhUsedMomenta_errorZ_refitted(NULL)
123  , fhUsedMomenta_vsZ_refitted(NULL)
124  , fhInvariantMass_pi0epem(NULL)
125  , fhPi0_startvertex(NULL)
126  , fhPi0_startvertexElectrons_all(NULL)
127  , fhPi0_startvertexElectrons_gg(NULL)
128  , fhPi0_startvertexElectrons_gee(NULL)
129  , fhPi0_startvertex_vs_chi(NULL)
130  , fhPi0_startvertex_vs_momentum(NULL)
131  , fhInvMassWithFullRecoCuts(NULL)
132  , fhEPEM_InDetector_invmass_gg_mc(NULL)
133  , fhEPEM_InDetector_invmass_gg_refitted(NULL)
134  , fhEPEM_InDetector_invmass_gee_mc(NULL)
135  , fhEPEM_InDetector_invmass_gee_refitted(NULL)
136  , fhEPEM_InDetector_invmass_all_mc(NULL)
137  , fhEPEM_InDetector_invmass_all_refitted(NULL)
138  , fhEPEM_pt_vs_p_all_mc(NULL)
139  , fhEPEM_pt_vs_p_all_refitted(NULL)
140  , fhEPEM_missingLepton_nofRingHits(NULL)
141  , fhEPEM_missingLepton_ringMid(NULL)
142  , fhEPEM_missingLepton_ringRadius(NULL)
143  , fhEPEM_missingLepton_distance(NULL)
144  , fhEPEM_missingLepton_selectionNN(NULL)
145  , fhEPEM_missingLepton_rings(NULL)
146  , fhEPEM_missingLepton_radius_vs_p(NULL)
147  , fhEPEM_missingLepton_ANNvalue(NULL)
148  , fhEPEM_identifiedLepton_nofRingHits(NULL)
149  , fhEPEM_identifiedLepton_ringMid(NULL)
150  , fhEPEM_identifiedLepton_ringRadius(NULL)
151  , fhEPEM_identifiedLepton_distance(NULL)
152  , fhEPEM_identifiedLepton_selectionNN(NULL)
153  , fhEPEM_identifiedLepton_rings(NULL)
154  , fhEPEM_identifiedLepton_radius_vs_p(NULL)
155  , fhEPEM_identifiedLepton_ANNvalue(NULL)
156  , fhEPEM_invmass_eta_mc(NULL)
157  , fhEPEM_invmass_eta_refitted(NULL)
158  , fhEPEM_efficiencyCuts_eta(NULL)
159  , timer()
160  , fTime(0.) {}
161 
162 CbmAnaConversionReco::~CbmAnaConversionReco() {}
163 
164 
165 void CbmAnaConversionReco::Init() {
166  FairRootManager* ioman = FairRootManager::Instance();
167  if (NULL == ioman) {
168  Fatal("CbmAnaConversion::Init", "RootManager not instantised!");
169  }
170 
171  fMcTracks = (TClonesArray*) ioman->GetObject("MCTrack");
172  if (NULL == fMcTracks) {
173  Fatal("CbmAnaConversion::Init", "No MCTrack array!");
174  }
175 
176  fGlobalTracks = (TClonesArray*) ioman->GetObject("GlobalTrack");
177  if (NULL == fGlobalTracks) {
178  Fatal("CbmAnaConversion::Init", "No GlobalTrack array!");
179  }
180 
181  fRichRings = (TClonesArray*) ioman->GetObject("RichRing");
182  if (NULL == fRichRings) {
183  Fatal("CbmAnaConversion::Init", "No RichRing array!");
184  }
185 
186  fRichHits = (TClonesArray*) ioman->GetObject("RichHit");
187  if (NULL == fRichHits) {
188  Fatal("CbmAnaConversion::Init", "No RichHit array!");
189  }
190 
191 
192  InitHistos();
193 }
194 
195 
196 void CbmAnaConversionReco::InitHistos() {
197  fHistoList_MC.clear();
198  fHistoList_reco.clear();
199  fHistoList_reco_mom.clear();
200 
201  fHistoList_gg.clear();
202  fHistoList_gee.clear();
203  fHistoList_eeee.clear();
204  fHistoList_all.clear();
205 
206  fHistoList_eta.clear();
207 
208 
209  Double_t invmassSpectra_nof = 800;
210  Double_t invmassSpectra_start = -0.00125;
211  Double_t invmassSpectra_end = 1.99875;
212 
213 
214  fhInvariantMass_MC_all =
215  new TH1D("fhInvariantMass_MC_all",
216  "fhInvariantMass_MC_all;invariant mass in GeV/c^{2};#",
217  2001,
218  -0.0005,
219  2.0005);
220  fhInvariantMass_MC_pi0 =
221  new TH1D("fhInvariantMass_MC_pi0",
222  "fhInvariantMass_MC_pi0;invariant mass in GeV/c^{2};#",
223  2001,
224  -0.0005,
225  2.0005);
226  fhInvariantMass_MC_pi0_epem =
227  new TH1D("fhInvariantMass_MC_pi0_epem",
228  "fhInvariantMass_MC_pi0_epem;invariant mass in GeV/c^{2};#",
229  2001,
230  -0.0005,
231  2.0005);
232  fhInvariantMass_MC_pi0_gepem =
233  new TH1D("fhInvariantMass_MC_pi0_gepem",
234  "fhInvariantMass_MC_pi0_gepem;invariant mass in GeV/c^{2};#",
235  2001,
236  -0.0005,
237  2.0005);
238  fhInvariantMass_MC_pi0_gg =
239  new TH1D("fhInvariantMass_MC_pi0_gg",
240  "fhInvariantMass_MC_pi0_gg;invariant mass in GeV/c^{2};#",
241  2001,
242  -0.0005,
243  2.0005);
244  fhInvariantMass_MC_eta =
245  new TH1D("fhInvariantMass_MC_eta",
246  "fhInvariantMass_MC_eta;invariant mass in GeV/c^{2};#",
247  2001,
248  -0.0005,
249  2.0005);
250  fhInvariantMass_MC_etaPrime =
251  new TH1D("fhInvariantMass_MC_etaPrime",
252  "fhInvariantMass_MC_etaPrime;invariant mass in GeV/c^{2};#",
253  2001,
254  -0.0005,
255  2.0005);
256  fHistoList_MC.push_back(fhInvariantMass_MC_all);
257  fHistoList_MC.push_back(fhInvariantMass_MC_pi0);
258  fHistoList_MC.push_back(fhInvariantMass_MC_pi0_epem);
259  fHistoList_MC.push_back(fhInvariantMass_MC_pi0_gepem);
260  fHistoList_MC.push_back(fhInvariantMass_MC_pi0_gg);
261  fHistoList_MC.push_back(fhInvariantMass_MC_eta);
262  fHistoList_MC.push_back(fhInvariantMass_MC_etaPrime);
263 
264 
265  fhMC_electrons_theta = new TH1D(
266  "fhMC_electrons_theta", "fhMC_electrons_theta;theta in deg;#", 90, 0., 90.);
267  fhMC_electrons_p = new TH1D(
268  "fhMC_electrons_p", "fhMC_electrons_p;momentum p in GeV/c;#", 100, 0., 10.);
269  fhMC_electrons_theta_vs_p =
270  new TH2D("fhMC_electrons_theta_vs_p",
271  "fhMC_electrons_theta_vs_p;theta in deg;momentum p in GeV/c",
272  90,
273  0.,
274  90.,
275  100,
276  0.,
277  10.);
278  fHistoList_MC.push_back(fhMC_electrons_theta);
279  fHistoList_MC.push_back(fhMC_electrons_p);
280  fHistoList_MC.push_back(fhMC_electrons_theta_vs_p);
281 
282  fhEta_openingAngleGG = new TH1D("fhEta_openingAngleGG",
283  "fhEta_openingAngleGG;theta in deg;#",
284  900,
285  0.,
286  90.);
287  fHistoList_MC.push_back(fhEta_openingAngleGG);
288 
289  fhMC_grandmotherPDGs = new TH1D("fhMC_grandmotherPDGs",
290  "fhMC_grandmotherPDGs;particle pdg;#",
291  1000,
292  0.,
293  1000.);
294  fHistoList_MC.push_back(fhMC_grandmotherPDGs);
295 
296  fhMCtest = new TH1D(
297  "fhMCtest", "fhMCtest;invariant mass in GeV/c^{2};#", 2000, 0., 2.);
298  fHistoList_MC.push_back(fhMCtest);
299 
300 
301  fhEPEM_invmass_gg_mc =
302  new TH1D("fhEPEM_invmass_gg_mc",
303  "fhEPEM_invmass_gg_mc;invariant mass in GeV/c^{2};#",
304  invmassSpectra_nof,
305  invmassSpectra_start,
306  invmassSpectra_end);
307  fhEPEM_invmass_gg_refitted =
308  new TH1D("fhEPEM_invmass_gg_refitted",
309  "fhEPEM_invmass_gg_refitted;invariant mass in GeV/c^{2};#",
310  invmassSpectra_nof,
311  invmassSpectra_start,
312  invmassSpectra_end);
313  fhEPEM_invmass_gee_mc =
314  new TH1D("fhEPEM_invmass_gee_mc",
315  "fhEPEM_invmass_gee_mc;invariant mass in GeV/c^{2};#",
316  invmassSpectra_nof,
317  invmassSpectra_start,
318  invmassSpectra_end);
319  fhEPEM_invmass_gee_refitted =
320  new TH1D("fhEPEM_invmass_gee_refitted",
321  "fhEPEM_invmass_gee_refitted;invariant mass in GeV/c^{2};#",
322  invmassSpectra_nof,
323  invmassSpectra_start,
324  invmassSpectra_end);
325  fhEPEM_invmass_eeee_mc =
326  new TH1D("fhEPEM_invmass_eeee_mc",
327  "fhEPEM_invmass_eeee_mc;invariant mass in GeV/c^{2};#",
328  invmassSpectra_nof,
329  invmassSpectra_start,
330  invmassSpectra_end);
331  fhEPEM_invmass_eeee_refitted =
332  new TH1D("fhEPEM_invmass_eeee_refitted",
333  "fhEPEM_invmass_eeee_refitted;invariant mass in GeV/c^{2};#",
334  invmassSpectra_nof,
335  invmassSpectra_start,
336  invmassSpectra_end);
337  fhEPEM_invmass_all_mc =
338  new TH1D("fhEPEM_invmass_all_mc",
339  "fhEPEM_invmass_all_mc;invariant mass in GeV/c^{2};#",
340  invmassSpectra_nof,
341  invmassSpectra_start,
342  invmassSpectra_end);
343  fhEPEM_invmass_all_refitted =
344  new TH1D("fhEPEM_invmass_all_refitted",
345  "fhEPEM_invmass_all_refitted;invariant mass in GeV/c^{2};#",
346  invmassSpectra_nof,
347  invmassSpectra_start,
348  invmassSpectra_end);
349  fHistoList_gg.push_back(fhEPEM_invmass_gg_mc);
350  fHistoList_gg.push_back(fhEPEM_invmass_gg_refitted);
351  fHistoList_gee.push_back(fhEPEM_invmass_gee_mc);
352  fHistoList_gee.push_back(fhEPEM_invmass_gee_refitted);
353  fHistoList_eeee.push_back(fhEPEM_invmass_eeee_mc);
354  fHistoList_eeee.push_back(fhEPEM_invmass_eeee_refitted);
355  fHistoList_all.push_back(fhEPEM_invmass_all_mc);
356  fHistoList_all.push_back(fhEPEM_invmass_all_refitted);
357 
358  fhEPEM_openingAngle_gg_mc =
359  new TH1D("fhEPEM_openingAngle_gg_mc",
360  "fhEPEM_openingAngle_gg_mc (between e+e- from g);angle [deg];#",
361  1010,
362  -0.1,
363  100.9);
364  fhEPEM_openingAngle_gg_refitted = new TH1D(
365  "fhEPEM_openingAngle_gg_refitted",
366  "fhEPEM_openingAngle_gg_refitted (between e+e- from g);angle [deg];#",
367  1010,
368  -0.1,
369  100.9);
370  fhEPEM_openingAngle_gee_mc =
371  new TH1D("fhEPEM_openingAngle_gee_mc",
372  "fhEPEM_openingAngle_gee_mc (between e+e- from g);angle [deg];#",
373  1010,
374  -0.1,
375  100.9);
376  fhEPEM_openingAngle_gee_refitted = new TH1D(
377  "fhEPEM_openingAngle_gee_refitted",
378  "fhEPEM_openingAngle_gee_refitted (between e+e- from g);angle [deg];#",
379  1010,
380  -0.1,
381  100.9);
382  fHistoList_gg.push_back(fhEPEM_openingAngle_gg_mc);
383  fHistoList_gg.push_back(fhEPEM_openingAngle_gg_refitted);
384  fHistoList_gee.push_back(fhEPEM_openingAngle_gee_mc);
385  fHistoList_gee.push_back(fhEPEM_openingAngle_gee_refitted);
386  fhEPEM_openingAngle_gee_mc_dalitz = new TH1D(
387  "fhEPEM_openingAngle_gee_mc_dalitz",
388  "fhEPEM_openingAngle_gee_mc_dalitz (between e+e- from pi0);angle [deg];#",
389  1010,
390  -0.1,
391  100.9);
392  fhEPEM_openingAngle_gee_refitted_dalitz =
393  new TH1D("fhEPEM_openingAngle_gee_refitted_dalitz",
394  "fhEPEM_openingAngle_gee_refitted_dalitz (between e+e- from "
395  "pi0);angle [deg];#",
396  1010,
397  -0.1,
398  100.9);
399  fHistoList_gee.push_back(fhEPEM_openingAngle_gee_mc_dalitz);
400  fHistoList_gee.push_back(fhEPEM_openingAngle_gee_refitted_dalitz);
401 
402  fhEPEM_openingAngle_vs_pt_gg_mc =
403  new TH2D("fhEPEM_openingAngle_vs_pt_gg_mc",
404  "fhEPEM_openingAngle_vs_pt_gg_mc;pt [GeV]; opening angle [deg]",
405  220,
406  -1.,
407  10.,
408  100,
409  0.,
410  10.);
411  fHistoList_gg.push_back(fhEPEM_openingAngle_vs_pt_gg_mc);
412  fhEPEM_openingAngle_vs_pt_gg_reco =
413  new TH2D("fhEPEM_openingAngle_vs_pt_gg_reco",
414  "fhEPEM_openingAngle_vs_pt_gg_reco;pt [GeV]; opening angle [deg]",
415  220,
416  -1.,
417  10.,
418  100,
419  0.,
420  10.);
421  fHistoList_gg.push_back(fhEPEM_openingAngle_vs_pt_gg_reco);
422 
423  fhEPEM_openingAngle_betweenGammas_mc =
424  new TH1D("fhEPEM_openingAngle_betweenGammas_mc",
425  "fhEPEM_openingAngle_betweenGammas_mc;angle [deg];#",
426  1010,
427  -0.1,
428  100.9);
429  fhEPEM_openingAngle_betweenGammas_reco =
430  new TH1D("fhEPEM_openingAngle_betweenGammas_reco",
431  "fhEPEM_openingAngle_betweenGammas_reco;angle [deg];#",
432  1010,
433  -0.1,
434  100.9);
435  fHistoList_gg.push_back(fhEPEM_openingAngle_betweenGammas_mc);
436  fHistoList_gg.push_back(fhEPEM_openingAngle_betweenGammas_reco);
437 
438 
439  fhPi0_pt_vs_rap_gg = new TH2D("fhPi0_pt_vs_rap_gg",
440  "fhPi0_pt_vs_rap_gg;pt [GeV]; rap [GeV]",
441  240,
442  -2.,
443  10.,
444  270,
445  -2.,
446  7.);
447  fhPi0_pt_vs_rap_gee = new TH2D("fhPi0_pt_vs_rap_gee",
448  "fhPi0_pt_vs_rap_gee;pt [GeV]; rap [GeV]",
449  240,
450  -2.,
451  10.,
452  270,
453  -2.,
454  7.);
455  fhPi0_pt_vs_rap_all = new TH2D("fhPi0_pt_vs_rap_all",
456  "fhPi0_pt_vs_rap_all;pt [GeV]; rap [GeV]",
457  240,
458  -2.,
459  10.,
460  270,
461  -2.,
462  7.);
463  fHistoList_gg.push_back(fhPi0_pt_vs_rap_gg);
464  fHistoList_gee.push_back(fhPi0_pt_vs_rap_gee);
465  fHistoList_all.push_back(fhPi0_pt_vs_rap_all);
466 
467  fhPi0_pt_gg = new TH1D("fhPi0_pt_gg", "fhPi0_pt_gg;pt [GeV];#", 200, 0., 10.);
468  fhPi0_pt_gee =
469  new TH1D("fhPi0_pt_gee", "fhPi0_pt_gee;pt [GeV];#", 200, 0., 10.);
470  fhPi0_pt_all =
471  new TH1D("fhPi0_pt_all", "fhPi0_pt_all;pt [GeV];#", 200, 0., 10.);
472  fHistoList_gg.push_back(fhPi0_pt_gg);
473  fHistoList_gee.push_back(fhPi0_pt_gee);
474  fHistoList_all.push_back(fhPi0_pt_all);
475 
476  fhEPEM_efficiencyCuts =
477  new TH1D("fhEPEM_efficiencyCuts", "fhEPEM_efficiencyCuts;;#", 13, 0., 13.);
478  fHistoList_all.push_back(fhEPEM_efficiencyCuts);
479  fhEPEM_efficiencyCuts->GetXaxis()->SetBinLabel(1, "no cuts");
480  fhEPEM_efficiencyCuts->GetXaxis()->SetBinLabel(2, "ANN: 4 rich electrons");
481  fhEPEM_efficiencyCuts->GetXaxis()->SetBinLabel(3, "ANN: #chi^{2}-cut");
482  fhEPEM_efficiencyCuts->GetXaxis()->SetBinLabel(
483  4, "ANN: #theta of e^{+}e^{-} pairs");
484  fhEPEM_efficiencyCuts->GetXaxis()->SetBinLabel(
485  5, "ANN: m_{inv} of e^{+}e^{-} pairs");
486  fhEPEM_efficiencyCuts->GetXaxis()->SetBinLabel(6, "Normal: 4 rich electrons");
487  fhEPEM_efficiencyCuts->GetXaxis()->SetBinLabel(7, "Normal: #chi^{2}-cut");
488  fhEPEM_efficiencyCuts->GetXaxis()->SetBinLabel(
489  8, "Normal: #theta of e^{+}e^{-} pairs");
490  fhEPEM_efficiencyCuts->GetXaxis()->SetBinLabel(
491  9, "Normal: m_{inv} of e^{+}e^{-} pairs");
492  fhEPEM_efficiencyCuts->GetXaxis()->SetBinLabel(10, "MC: 4 rich electrons");
493  fhEPEM_efficiencyCuts->GetXaxis()->SetBinLabel(11, "MC: #chi^{2}-cut");
494  fhEPEM_efficiencyCuts->GetXaxis()->SetBinLabel(
495  12, "MC: #theta of e^{+}e^{-} pairs");
496  fhEPEM_efficiencyCuts->GetXaxis()->SetBinLabel(
497  13, "MC: m_{inv} of e^{+}e^{-} pairs");
498 
499  fhEPEM_efficiencyCuts2 = new TH1D(
500  "fhEPEM_efficiencyCuts2", "fhEPEM_efficiencyCuts2;;#", 10, 0., 10.);
501  fHistoList_all.push_back(fhEPEM_efficiencyCuts2);
502 
503  fhEPEM_pi0_nofLeptons_ann = new TH1D(
504  "fhEPEM_pi0_nofLeptons_ann", "fhEPEM_pi0_nofLeptons_ann;;#", 5, -0.5, 4.5);
505  fhEPEM_pi0_ANNvalues_noCuts = new TH1D("fhEPEM_pi0_ANNvalues_noCuts",
506  "fhEPEM_pi0_ANNvalues_noCuts;;#",
507  400,
508  -2,
509  2);
510  fhEPEM_pi0_ANNvalues_angleCut = new TH1D("fhEPEM_pi0_ANNvalues_angleCut",
511  "fhEPEM_pi0_ANNvalues_angleCut;;#",
512  400,
513  -2,
514  2);
515  fhEPEM_pi0_ANNefficiencies = new TH1D(
516  "fhEPEM_pi0_ANNefficiencies", "fhEPEM_pi0_ANNefficiencies;;#", 15, 0, 15);
517  fHistoList_all.push_back(fhEPEM_pi0_ANNvalues_noCuts);
518  fHistoList_all.push_back(fhEPEM_pi0_ANNvalues_angleCut);
519  fHistoList_all.push_back(fhEPEM_pi0_nofLeptons_ann);
520  fHistoList_all.push_back(fhEPEM_pi0_ANNefficiencies);
521 
522  fhEPEM_pi0_ANNefficiencies->GetXaxis()->SetBinLabel(1, "no cuts");
523  fhEPEM_pi0_ANNefficiencies->GetXaxis()->SetBinLabel(2, "-1.0");
524  fhEPEM_pi0_ANNefficiencies->GetXaxis()->SetBinLabel(3, "-0.9");
525  fhEPEM_pi0_ANNefficiencies->GetXaxis()->SetBinLabel(4, "-0.8");
526  fhEPEM_pi0_ANNefficiencies->GetXaxis()->SetBinLabel(5, "-0.7");
527  fhEPEM_pi0_ANNefficiencies->GetXaxis()->SetBinLabel(6, "-0.6");
528  fhEPEM_pi0_ANNefficiencies->GetXaxis()->SetBinLabel(7, "-0.5");
529  fhEPEM_pi0_ANNefficiencies->GetXaxis()->SetBinLabel(8, "-0.4");
530  fhEPEM_pi0_ANNefficiencies->GetXaxis()->SetBinLabel(9, "-0.3");
531  fhEPEM_pi0_ANNefficiencies->GetXaxis()->SetBinLabel(10, "-0.2");
532  fhEPEM_pi0_ANNefficiencies->GetXaxis()->SetBinLabel(11, "-0.1");
533  fhEPEM_pi0_ANNefficiencies->GetXaxis()->SetBinLabel(12, "0.0");
534 
535  fhEPEM_rap_vs_chi = new TH2D("fhEPEM_rap_vs_chi",
536  "fhEPEM_rap_vs_chi; rap [GeV]; chi of electrons",
537  300,
538  0.,
539  10.,
540  100,
541  0.,
542  100.);
543  fHistoList_all.push_back(fhEPEM_rap_vs_chi);
544  fhEPEM_rap_vs_invmass = new TH2D("fhEPEM_rap_vs_invmass",
545  "fhEPEM_rap_vs_invmass; rap [GeV]; invmass",
546  300,
547  0.,
548  10.,
549  100,
550  0.,
551  10.);
552  fHistoList_all.push_back(fhEPEM_rap_vs_invmass);
553 
554 
555  fhInvMass_EPEM_mc =
556  new TH1D("fhInvMass_EPEM_mc",
557  "fhInvariantMass_recoMomentum1 (mc);mass [GeV/c^2];#",
558  invmassSpectra_nof,
559  invmassSpectra_start,
560  invmassSpectra_end);
561  fhInvMass_EPEM_stsMomVec =
562  new TH1D("fhInvMass_EPEM_stsMomVec",
563  "fhInvariantMass_recoMomentum2 (stsMomentumVec);mass [GeV/c^2];#",
564  invmassSpectra_nof,
565  invmassSpectra_start,
566  invmassSpectra_end);
567  fhInvMass_EPEM_refitted = new TH1D(
568  "fhInvMass_EPEM_refitted",
569  "fhInvariantMass_recoMomentum3 (refitted at primary);mass [GeV/c^2];#",
570  invmassSpectra_nof,
571  invmassSpectra_start,
572  invmassSpectra_end);
573  fhInvMass_EPEM_error_stsMomVec = new TH1D(
574  "fhInvMass_EPEM_error_stsMomVec",
575  "fhInvariantMass_recoMomentum4 (error, stsMomentumVec);(mc-reco)/mc;#",
576  500,
577  -0.005,
578  4.995);
579  fhInvMass_EPEM_error_refitted =
580  new TH1D("fhInvMass_EPEM_error_refitted",
581  "fhInvariantMass_recoMomentum5 (error, refitted);(mc-reco)/mc;#",
582  500,
583  -0.005,
584  4.995);
585  fhInvMass_EPEM_openingAngleRef =
586  new TH1D("fhInvMass_EPEM_openingAngleRef",
587  "fhInvariantMass_openingAngleRef;angle [deg];#",
588  1010,
589  -0.1,
590  100.9);
591  fhUsedMomenta_stsMomVec = new TH1D(
592  "fhUsedMomenta_stsMomVec", "fhMomentumtest1;momentum;#", 100, 0., 2.);
593  fhUsedMomenta_mc =
594  new TH1D("fhUsedMomenta_mc", "fhMomentumtest2;momentum;#", 100, 0., 2.);
595  fhUsedMomenta_error_stsMomVec =
596  new TH1D("fhUsedMomenta_error_stsMomVec",
597  "fhMomentumtest3 (error);(mc-reco)/mc;#",
598  400,
599  -2.005,
600  1.995);
601  fhUsedMomenta_error_refitted =
602  new TH1D("fhUsedMomenta_error_refitted",
603  "fhMomentumtest4 (error);(mc-reco_refitted)/mc;#",
604  400,
605  -2.005,
606  1.995);
607 
608  fhUsedMomenta_errorX_stsMomVec =
609  new TH1D("fhUsedMomenta_errorX_stsMomVec",
610  "fhMomentumtest5 (error of x-momentum);(mc-reco_reco)/mc;#",
611  401,
612  -4.01,
613  4.01);
614  fhUsedMomenta_vsX_stsMomVec =
615  new TH2D("fhUsedMomenta_vsX_stsMomVec",
616  "fhMomentumtest5vs (error of x-momentum);mc;reco",
617  101,
618  -1.01,
619  1.01,
620  101,
621  -1.01,
622  1.01);
623  fhUsedMomenta_errorY_stsMomVec =
624  new TH1D("fhUsedMomenta_errorY_stsMomVec",
625  "fhMomentumtest6 (error of y-momentum);(mc-reco_reco)/mc;#",
626  401,
627  -4.01,
628  4.01);
629  fhUsedMomenta_vsY_stsMomVec =
630  new TH2D("fhUsedMomenta_vsY_stsMomVec",
631  "fhMomentumtest6vs (error of y-momentum);mc;reco",
632  101,
633  -1.01,
634  1.01,
635  101,
636  -1.01,
637  1.01);
638  fhUsedMomenta_errorZ_stsMomVec =
639  new TH1D("fhUsedMomenta_errorZ_stsMomVec",
640  "fhMomentumtest7 (error of z-momentum);(mc-reco_reco)/mc;#",
641  401,
642  -4.01,
643  4.01);
644  fhUsedMomenta_vsZ_stsMomVec =
645  new TH2D("fhUsedMomenta_vsZ_stsMomVec",
646  "fhMomentumtest7vs (error of z-momentum);mc;reco",
647  201,
648  -0.01,
649  4.01,
650  201,
651  -0.01,
652  4.01);
653  fhUsedMomenta_errorX_refitted =
654  new TH1D("fhUsedMomenta_errorX_refitted",
655  "fhMomentumtest5 (error of x-momentum);(mc-reco_reco)/mc;#",
656  401,
657  -4.01,
658  4.01);
659  fhUsedMomenta_vsX_refitted =
660  new TH2D("fhUsedMomenta_vsX_refitted",
661  "fhMomentumtest5vs (error of x-momentum);mc;reco",
662  101,
663  -1.01,
664  1.01,
665  101,
666  -1.01,
667  1.01);
668  fhUsedMomenta_errorY_refitted =
669  new TH1D("fhUsedMomenta_errorY_refitted",
670  "fhMomentumtest6 (error of y-momentum);(mc-reco_reco)/mc;#",
671  401,
672  -4.01,
673  4.01);
674  fhUsedMomenta_vsY_refitted =
675  new TH2D("fhUsedMomenta_vsY_refitted",
676  "fhMomentumtest6vs (error of y-momentum);mc;reco",
677  101,
678  -1.01,
679  1.01,
680  101,
681  -1.01,
682  1.01);
683  fhUsedMomenta_errorZ_refitted =
684  new TH1D("fhUsedMomenta_errorZ_refitted",
685  "fhMomentumtest7 (error of z-momentum);(mc-reco_reco)/mc;#",
686  401,
687  -4.01,
688  4.01);
689  fhUsedMomenta_vsZ_refitted =
690  new TH2D("fhUsedMomenta_vsZ_refitted",
691  "fhMomentumtest7vs (error of z-momentum);mc;reco",
692  201,
693  -0.01,
694  4.01,
695  201,
696  -0.01,
697  4.01);
698  fHistoList_reco.push_back(fhInvMass_EPEM_mc);
699  fHistoList_reco.push_back(fhInvMass_EPEM_stsMomVec);
700  fHistoList_reco.push_back(fhInvMass_EPEM_refitted);
701  fHistoList_reco.push_back(fhInvMass_EPEM_error_stsMomVec);
702  fHistoList_reco.push_back(fhInvMass_EPEM_error_refitted);
703  fHistoList_reco.push_back(fhInvMass_EPEM_openingAngleRef);
704  fHistoList_reco_mom.push_back(fhUsedMomenta_stsMomVec);
705  fHistoList_reco_mom.push_back(fhUsedMomenta_mc);
706  fHistoList_reco_mom.push_back(fhUsedMomenta_error_stsMomVec);
707  fHistoList_reco_mom.push_back(fhUsedMomenta_error_refitted);
708  fHistoList_reco_mom.push_back(
709  fhUsedMomenta_errorX_stsMomVec); // error of x-component of reconstructed momentum
710  fHistoList_reco_mom.push_back(
711  fhUsedMomenta_errorY_stsMomVec); // error of y-component of reconstructed momentum
712  fHistoList_reco_mom.push_back(
713  fhUsedMomenta_errorZ_stsMomVec); // error of z-component of reconstructed momentum
714  fHistoList_reco_mom.push_back(
715  fhUsedMomenta_vsX_stsMomVec); // x-component of reconstructed momentum vs mc-momentum
716  fHistoList_reco_mom.push_back(
717  fhUsedMomenta_vsY_stsMomVec); // y-component of reconstructed momentum vs mc-momentum
718  fHistoList_reco_mom.push_back(
719  fhUsedMomenta_vsZ_stsMomVec); // z-component of reconstructed momentum vs mc-momentum
720  fHistoList_reco_mom.push_back(
721  fhUsedMomenta_errorX_refitted); // error of x-component of reconstructed momentum
722  fHistoList_reco_mom.push_back(
723  fhUsedMomenta_errorY_refitted); // error of y-component of reconstructed momentum
724  fHistoList_reco_mom.push_back(
725  fhUsedMomenta_errorZ_refitted); // error of z-component of reconstructed momentum
726  fHistoList_reco_mom.push_back(
727  fhUsedMomenta_vsX_refitted); // x-component of reconstructed momentum vs mc-momentum
728  fHistoList_reco_mom.push_back(
729  fhUsedMomenta_vsY_refitted); // y-component of reconstructed momentum vs mc-momentum
730  fHistoList_reco_mom.push_back(
731  fhUsedMomenta_vsZ_refitted); // z-component of reconstructed momentum vs mc-momentum
732 
733  fhInvariantMass_pi0epem = new TH1D("fhInvariantMass_pi0epem",
734  "fhInvariantMass_pi0epem;mass [GeV/c^2];#",
735  400,
736  0.,
737  2.);
738  fHistoList_reco.push_back(fhInvariantMass_pi0epem);
739 
740  fhPi0_startvertex = new TH1D(
741  "fhPi0_startvertex", "fhPi0_startvertex;z in cm;#", 210, -5., 100.);
742  fHistoList_reco.push_back(fhPi0_startvertex);
743 
744  fhPi0_startvertexElectrons_all =
745  new TH1D("fhPi0_startvertexElectrons_all",
746  "fhPi0_startvertexElectrons_all;z in cm;#",
747  411,
748  -5.25,
749  200.25);
750  fHistoList_reco.push_back(fhPi0_startvertexElectrons_all);
751 
752  fhPi0_startvertexElectrons_gg =
753  new TH1D("fhPi0_startvertexElectrons_gg",
754  "fhPi0_startvertexElectrons_gg;z in cm;#",
755  411,
756  -5.25,
757  200.25);
758  fHistoList_reco.push_back(fhPi0_startvertexElectrons_gg);
759 
760  fhPi0_startvertexElectrons_gee =
761  new TH1D("fhPi0_startvertexElectrons_gee",
762  "fhPi0_startvertexElectrons_gee;z in cm;#",
763  411,
764  -5.25,
765  200.25);
766  fHistoList_reco.push_back(fhPi0_startvertexElectrons_gee);
767 
768  fhPi0_startvertex_vs_chi = new TH2D("fhPi0_startvertex_vs_chi",
769  "fhPi0_startvertex_vs_chi;z[cm];chi",
770  210,
771  -5.,
772  100.,
773  1000,
774  0.,
775  100.);
776  fHistoList_reco.push_back(fhPi0_startvertex_vs_chi);
777 
778  fhPi0_startvertex_vs_momentum =
779  new TH2D("fhPi0_startvertex_vs_momentum",
780  "fhPi0_startvertex_vs_momentum;z in cm;momentum (MC-true)",
781  210,
782  -5.,
783  100.,
784  1000,
785  0.,
786  100.);
787  fHistoList_reco.push_back(fhPi0_startvertex_vs_momentum);
788 
789  fhInvMassWithFullRecoCuts =
790  new TH1D("fhInvMassWithFullRecoCuts",
791  "fhInvMassWithFullRecoCuts;mass [GeV/c^2];#",
792  800,
793  0.,
794  2.);
795  fHistoList_reco.push_back(fhInvMassWithFullRecoCuts);
796 
797 
798  fhEPEM_InDetector_invmass_gg_mc =
799  new TH1D("fhEPEM_InDetector_invmass_gg_mc",
800  "fhEPEM_InDetector_invmass_gg_mc;mass [GeV/c^2];#",
801  invmassSpectra_nof,
802  invmassSpectra_start,
803  invmassSpectra_end);
804  fhEPEM_InDetector_invmass_gg_refitted =
805  new TH1D("fhEPEM_InDetector_invmass_gg_refitted",
806  "fhEPEM_InDetector_invmass_gg_refitted;mass [GeV/c^2];#",
807  invmassSpectra_nof,
808  invmassSpectra_start,
809  invmassSpectra_end);
810  fhEPEM_InDetector_invmass_gee_mc =
811  new TH1D("fhEPEM_InDetector_invmass_gee_mc",
812  "fhEPEM_InDetector_invmass_gee_mc;mass [GeV/c^2];#",
813  invmassSpectra_nof,
814  invmassSpectra_start,
815  invmassSpectra_end);
816  fhEPEM_InDetector_invmass_gee_refitted =
817  new TH1D("fhEPEM_InDetector_invmass_gee_refitted",
818  "fhEPEM_InDetector_invmass_gee_refitted;mass [GeV/c^2];#",
819  invmassSpectra_nof,
820  invmassSpectra_start,
821  invmassSpectra_end);
822  fhEPEM_InDetector_invmass_all_mc =
823  new TH1D("fhEPEM_InDetector_invmass_all_mc",
824  "fhEPEM_InDetector_invmass_all_mc;mass [GeV/c^2];#",
825  invmassSpectra_nof,
826  invmassSpectra_start,
827  invmassSpectra_end);
828  fhEPEM_InDetector_invmass_all_refitted =
829  new TH1D("fhEPEM_InDetector_invmass_all_refitted",
830  "fhEPEM_InDetector_invmass_all_refitted;mass [GeV/c^2];#",
831  invmassSpectra_nof,
832  invmassSpectra_start,
833  invmassSpectra_end);
834  fHistoList_gg.push_back(fhEPEM_InDetector_invmass_gg_mc);
835  fHistoList_gg.push_back(fhEPEM_InDetector_invmass_gg_refitted);
836  fHistoList_gee.push_back(fhEPEM_InDetector_invmass_gee_mc);
837  fHistoList_gee.push_back(fhEPEM_InDetector_invmass_gee_refitted);
838  fHistoList_all.push_back(fhEPEM_InDetector_invmass_all_mc);
839  fHistoList_all.push_back(fhEPEM_InDetector_invmass_all_refitted);
840 
841 
842  fhEPEM_pt_vs_p_all_mc =
843  new TH2D("fhEPEM_pt_vs_p_all_mc",
844  "fhEPEM_pt_vs_p_all_mc;p_{t} in GeV/c; p in GeV/c",
845  240,
846  -2.,
847  10.,
848  360,
849  -2.,
850  16.);
851  fHistoList_all.push_back(fhEPEM_pt_vs_p_all_mc);
852  fhEPEM_pt_vs_p_all_refitted =
853  new TH2D("fhEPEM_pt_vs_p_all_refitted",
854  "fhTest2_electrons_pt_vs_p;p_{t} in GeV/c; p in GeV/c",
855  240,
856  -2.,
857  10.,
858  360,
859  -2.,
860  16.);
861  fHistoList_all.push_back(fhEPEM_pt_vs_p_all_refitted);
862 
863  fhEPEM_missingLepton_nofRingHits =
864  new TH1D("fhEPEM_missingLepton_nofRingHits",
865  "fhEPEM_missingLepton_nofRingHits;nofringhits;#",
866  41,
867  -0.5,
868  40.5);
869  fhEPEM_missingLepton_ringMid = new TH2D("fhEPEM_missingLepton_ringMid",
870  "fhEPEM_missingLepton_ringMid;X;Y",
871  100,
872  -150,
873  150,
874  150,
875  -250,
876  250);
877  fhEPEM_missingLepton_ringRadius =
878  new TH1D("fhEPEM_missingLepton_ringRadius",
879  "fhEPEM_missingLepton_ringRadius;ringRadius;#",
880  100,
881  0,
882  10);
883  fhEPEM_missingLepton_distance =
884  new TH1D("fhEPEM_missingLepton_distance",
885  "fhEPEM_missingLepton_distance;distance;#",
886  100,
887  0,
888  20);
889  fhEPEM_missingLepton_selectionNN =
890  new TH1D("fhEPEM_missingLepton_selectionNN",
891  "fhEPEM_missingLepton_selectionNN;selectionNN;#",
892  40,
893  -2,
894  2);
895  fhEPEM_missingLepton_rings =
896  new TH2D("fhEPEM_missingLepton_rings",
897  "fhEPEM_missingLepton_rings;selectionNN;#",
898  400,
899  -200,
900  200,
901  600,
902  -300,
903  300);
904  fhEPEM_missingLepton_radius_vs_p = new TH2D(
905  "fhEPEM_missingLepton_radius_vs_p",
906  "fhEPEM_missingLepton_radius_vs_p;momentum p in GeV/c;ring radius in cm",
907  120,
908  0,
909  12,
910  100,
911  0,
912  10);
913  fhEPEM_missingLepton_ANNvalue =
914  new TH1D("fhEPEM_missingLepton_ANNvalue",
915  "fhEPEM_missingLepton_ANNvalue;ANNvalue;#",
916  40,
917  -2,
918  2);
919  fHistoList_all.push_back(fhEPEM_missingLepton_nofRingHits);
920  fHistoList_all.push_back(fhEPEM_missingLepton_ringMid);
921  fHistoList_all.push_back(fhEPEM_missingLepton_ringRadius);
922  fHistoList_all.push_back(fhEPEM_missingLepton_distance);
923  fHistoList_all.push_back(fhEPEM_missingLepton_selectionNN);
924  fHistoList_all.push_back(fhEPEM_missingLepton_rings);
925  fHistoList_all.push_back(fhEPEM_missingLepton_radius_vs_p);
926  fHistoList_all.push_back(fhEPEM_missingLepton_ANNvalue);
927 
928 
929  fhEPEM_identifiedLepton_nofRingHits =
930  new TH1D("fhEPEM_identifiedLepton_nofRingHits",
931  "fhEPEM_identifiedLepton_nofRingHits;nofringhits;#",
932  41,
933  -0.5,
934  40.5);
935  fhEPEM_identifiedLepton_ringMid =
936  new TH2D("fhEPEM_identifiedLepton_ringMid",
937  "fhEPEM_identifiedLepton_ringMid;X;Y",
938  100,
939  -150,
940  150,
941  150,
942  -250,
943  250);
944  fhEPEM_identifiedLepton_ringRadius =
945  new TH1D("fhEPEM_identifiedLepton_ringRadius",
946  "fhEPEM_identifiedLepton_ringRadius;ringRadius;#",
947  100,
948  0,
949  10);
950  fhEPEM_identifiedLepton_distance =
951  new TH1D("fhEPEM_identifiedLepton_distance",
952  "fhEPEM_identifiedLepton_distance;distance;#",
953  100,
954  0,
955  20);
956  fhEPEM_identifiedLepton_selectionNN =
957  new TH1D("fhEPEM_identifiedLepton_selectionNN",
958  "fhEPEM_identifiedLepton_selectionNN;selectionNN;#",
959  40,
960  -2,
961  2);
962  fhEPEM_identifiedLepton_rings =
963  new TH2D("fhEPEM_identifiedLepton_rings",
964  "fhEPEM_identifiedLepton_rings;selectionNN;#",
965  400,
966  -200,
967  200,
968  600,
969  -300,
970  300);
971  fhEPEM_identifiedLepton_radius_vs_p = new TH2D(
972  "fhEPEM_identifiedLepton_radius_vs_p",
973  "fhEPEM_identifiedLepton_radius_vs_p;momentum p in GeV/c;ring radius in cm",
974  120,
975  0,
976  12,
977  100,
978  0,
979  10);
980  fhEPEM_identifiedLepton_ANNvalue =
981  new TH1D("fhEPEM_identifiedLepton_ANNvalue",
982  "fhEPEM_identifiedLepton_ANNvalue;ANNvalue;#",
983  40,
984  -2,
985  2);
986  fHistoList_all.push_back(fhEPEM_identifiedLepton_nofRingHits);
987  fHistoList_all.push_back(fhEPEM_identifiedLepton_ringMid);
988  fHistoList_all.push_back(fhEPEM_identifiedLepton_ringRadius);
989  fHistoList_all.push_back(fhEPEM_identifiedLepton_distance);
990  fHistoList_all.push_back(fhEPEM_identifiedLepton_selectionNN);
991  fHistoList_all.push_back(fhEPEM_identifiedLepton_rings);
992  fHistoList_all.push_back(fhEPEM_identifiedLepton_radius_vs_p);
993  fHistoList_all.push_back(fhEPEM_identifiedLepton_ANNvalue);
994 
995 
996  // histograms for eta analysis
997  fhEPEM_invmass_eta_mc = new TH1D("fhEPEM_invmass_eta_mc",
998  "fhEPEM_Iinvmass_eta_mc;mass [GeV/c^2];#",
999  invmassSpectra_nof,
1000  invmassSpectra_start,
1001  invmassSpectra_end);
1002  fhEPEM_invmass_eta_refitted =
1003  new TH1D("fhEPEM_invmass_eta_refitted",
1004  "fhEPEM_Iinvmass_eta_refitted;mass [GeV/c^2];#",
1005  invmassSpectra_nof,
1006  invmassSpectra_start,
1007  invmassSpectra_end);
1008  fHistoList_eta.push_back(fhEPEM_invmass_eta_mc);
1009  fHistoList_eta.push_back(fhEPEM_invmass_eta_refitted);
1010 
1011  fhEPEM_efficiencyCuts_eta = new TH1D(
1012  "fhEPEM_efficiencyCuts_eta", "fhEPEM_efficiencyCuts_eta;;#", 13, 0., 13.);
1013  fHistoList_eta.push_back(fhEPEM_efficiencyCuts_eta);
1014  fhEPEM_efficiencyCuts_eta->GetXaxis()->SetBinLabel(1, "no cuts");
1015  fhEPEM_efficiencyCuts_eta->GetXaxis()->SetBinLabel(2,
1016  "ANN: 4 rich electrons");
1017  fhEPEM_efficiencyCuts_eta->GetXaxis()->SetBinLabel(3, "ANN: #chi^{2}-cut");
1018  fhEPEM_efficiencyCuts_eta->GetXaxis()->SetBinLabel(
1019  4, "ANN: #theta of e^{+}e^{-} pairs");
1020  fhEPEM_efficiencyCuts_eta->GetXaxis()->SetBinLabel(
1021  5, "ANN: m_{inv} of e^{+}e^{-} pairs");
1022  fhEPEM_efficiencyCuts_eta->GetXaxis()->SetBinLabel(
1023  6, "Normal: 4 rich electrons");
1024  fhEPEM_efficiencyCuts_eta->GetXaxis()->SetBinLabel(7, "Normal: #chi^{2}-cut");
1025  fhEPEM_efficiencyCuts_eta->GetXaxis()->SetBinLabel(
1026  8, "Normal: #theta of e^{+}e^{-} pairs");
1027  fhEPEM_efficiencyCuts_eta->GetXaxis()->SetBinLabel(
1028  9, "Normal: m_{inv} of e^{+}e^{-} pairs");
1029  fhEPEM_efficiencyCuts_eta->GetXaxis()->SetBinLabel(10,
1030  "MC: 4 rich electrons");
1031  fhEPEM_efficiencyCuts_eta->GetXaxis()->SetBinLabel(11, "MC: #chi^{2}-cut");
1032  fhEPEM_efficiencyCuts_eta->GetXaxis()->SetBinLabel(
1033  12, "MC: #theta of e^{+}e^{-} pairs");
1034  fhEPEM_efficiencyCuts_eta->GetXaxis()->SetBinLabel(
1035  13, "MC: m_{inv} of e^{+}e^{-} pairs");
1036 }
1037 
1038 
1039 void CbmAnaConversionReco::Finish() {
1040  //gDirectory->cd("analysis-conversion");
1041  gDirectory->mkdir("MCreco");
1042  gDirectory->cd("MCreco");
1043  for (UInt_t i = 0; i < fHistoList_MC.size(); i++) {
1044  fHistoList_MC[i]->Write();
1045  }
1046  gDirectory->cd("..");
1047 
1048  gDirectory->mkdir("Reconstruction2");
1049  gDirectory->cd("Reconstruction2");
1050 
1051  gDirectory->mkdir("pi0 -> gg");
1052  gDirectory->cd("pi0 -> gg");
1053  for (UInt_t i = 0; i < fHistoList_gg.size(); i++) {
1054  fHistoList_gg[i]->Write();
1055  }
1056  gDirectory->cd("..");
1057 
1058  gDirectory->mkdir("pi0 -> ge+e-");
1059  gDirectory->cd("pi0 -> ge+e-");
1060  for (UInt_t i = 0; i < fHistoList_gee.size(); i++) {
1061  fHistoList_gee[i]->Write();
1062  }
1063  gDirectory->cd("..");
1064 
1065  gDirectory->mkdir("pi0 -> e+e-e+e-");
1066  gDirectory->cd("pi0 -> e+e-e+e-");
1067  for (UInt_t i = 0; i < fHistoList_eeee.size(); i++) {
1068  fHistoList_eeee[i]->Write();
1069  }
1070  gDirectory->cd("..");
1071 
1072  gDirectory->mkdir("pi0 -> all");
1073  gDirectory->cd("pi0 -> all");
1074  for (UInt_t i = 0; i < fHistoList_all.size(); i++) {
1075  fHistoList_all[i]->Write();
1076  }
1077  gDirectory->cd("..");
1078 
1079  gDirectory->mkdir("eta");
1080  gDirectory->cd("eta");
1081  for (UInt_t i = 0; i < fHistoList_eta.size(); i++) {
1082  fHistoList_eta[i]->Write();
1083  }
1084  gDirectory->cd("..");
1085 
1086  for (UInt_t i = 0; i < fHistoList_reco.size(); i++) {
1087  fHistoList_reco[i]->Write();
1088  }
1089  gDirectory->mkdir("Momenta2");
1090  gDirectory->cd("Momenta2");
1091  for (UInt_t i = 0; i < fHistoList_reco_mom.size(); i++) {
1092  fHistoList_reco_mom[i]->Write();
1093  }
1094  gDirectory->cd("../..");
1095 
1096  cout << "CbmAnaConversionReco: Realtime - " << fTime << endl;
1097  //timer.Print();
1098 }
1099 
1100 
1101 void CbmAnaConversionReco::SetTracklistMC(vector<CbmMCTrack*> MCTracklist) {
1102  fMCTracklist_all = MCTracklist;
1103 }
1104 
1105 
1106 void CbmAnaConversionReco::SetTracklistReco(vector<CbmMCTrack*> MCTracklist,
1107  vector<TVector3> RecoTracklist1,
1108  vector<TVector3> RecoTracklist2,
1109  vector<int> ids,
1110  vector<Double_t> chi,
1111  vector<Int_t> GlobalTrackId) {
1112  fRecoTracklistEPEM = MCTracklist;
1113  fRecoTracklistEPEM_ids = ids;
1114  fRecoTracklistEPEM_chi = chi;
1115  fRecoTracklistEPEM_gtid = GlobalTrackId;
1116  fRecoMomentum = RecoTracklist1;
1117  fRecoRefittedMomentum = RecoTracklist2;
1118 }
1119 
1120 
1121 void CbmAnaConversionReco::InvariantMassMC_all()
1122 // calculation of invariant mass via X -> gamma gamma -> e+ e- e+ e-, only MC data with several cuts (see FillMCTrackslists())
1123 {
1124  timer.Start();
1125 
1126 
1127  cout << "CbmAnaConversionReco: InvariantMassTestMC - Start..." << endl;
1128  cout << "CbmAnaConversionReco: InvariantMassTestMC - Size of "
1129  "fTracklistMC_all:\t "
1130  << fMCTracklist_all.size() << endl;
1131  if (fMCTracklist_all.size() >= 4) {
1132  for (unsigned int i = 0; i < fMCTracklist_all.size() - 3; i++) {
1133  if (i % 10 == 0)
1134  cout << "CbmAnaConversionReco: InvariantMassTestMC - iteration i = "
1135  << i << endl;
1136  for (unsigned int j = i + 1; j < fMCTracklist_all.size() - 2; j++) {
1137  for (unsigned int k = j + 1; k < fMCTracklist_all.size() - 1; k++) {
1138  for (unsigned int l = k + 1; l < fMCTracklist_all.size(); l++) {
1139 
1140  if (fMCTracklist_all[i]->GetPdgCode()
1141  + fMCTracklist_all[j]->GetPdgCode()
1142  + fMCTracklist_all[k]->GetPdgCode()
1143  + fMCTracklist_all[l]->GetPdgCode()
1144  != 0)
1145  continue;
1146 
1147  if (i == j || i == k || i == l || j == k || j == l || k == l)
1148  continue;
1149 
1150  int motherId1 = fMCTracklist_all[i]->GetMotherId();
1151  int motherId2 = fMCTracklist_all[j]->GetMotherId();
1152  int motherId3 = fMCTracklist_all[k]->GetMotherId();
1153  int motherId4 = fMCTracklist_all[l]->GetMotherId();
1154 
1155  TVector3 momentum1, momentum2, momentum3, momentum4;
1156  fMCTracklist_all[i]->GetMomentum(momentum1);
1157  fMCTracklist_all[j]->GetMomentum(momentum2);
1158  fMCTracklist_all[k]->GetMomentum(momentum3);
1159  fMCTracklist_all[l]->GetMomentum(momentum4);
1160 
1161 
1162  // decay pi0 -> e+ e- e+ e-
1163  if (motherId1 == motherId2 && motherId1 == motherId3
1164  && motherId1 == motherId4) {
1165  cout << "testxyz" << endl;
1166  Double_t invmass = Invmass_4particles(fMCTracklist_all[i],
1167  fMCTracklist_all[j],
1168  fMCTracklist_all[k],
1169  fMCTracklist_all[l]);
1170  fhMCtest->Fill(invmass);
1171 
1172  if (motherId1 != -1) {
1173  int mcMotherPdg1 = -1;
1174  CbmMCTrack* mother1 = (CbmMCTrack*) fMcTracks->At(motherId1);
1175  if (NULL != mother1) mcMotherPdg1 = mother1->GetPdgCode();
1176  if (mcMotherPdg1 == 111) {
1177  Double_t invmass2 = 0;
1178  invmass2 = Invmass_4particles(fMCTracklist_all[i],
1179  fMCTracklist_all[j],
1180  fMCTracklist_all[k],
1181  fMCTracklist_all[l]);
1182  fhMCtest->Fill(invmass2);
1183  fhInvariantMass_MC_pi0->Fill(invmass2);
1184  fhInvariantMass_MC_pi0_epem->Fill(invmass2);
1185  cout << "####################################### Decay pi0 "
1186  "-> e+e-e+e- detected!\t\t"
1187  << invmass2 << endl;
1188  }
1189  } else {
1190  continue;
1191  }
1192  }
1193 
1194 
1195  int grandmotherId1 = -1;
1196  int grandmotherId2 = -1;
1197  int grandmotherId3 = -1;
1198  int grandmotherId4 = -1;
1199 
1200  int mcMotherPdg1 = -1;
1201  int mcMotherPdg2 = -1;
1202  int mcMotherPdg3 = -1;
1203  // int mcMotherPdg4 = -1;
1204  int mcGrandmotherPdg1 = -1;
1205  // int mcGrandmotherPdg2 = -1;
1206  // int mcGrandmotherPdg3 = -1;
1207  // int mcGrandmotherPdg4 = -1;
1208 
1209 
1210  if (motherId1 != -1) {
1211  CbmMCTrack* mother1 = (CbmMCTrack*) fMcTracks->At(motherId1);
1212  if (NULL != mother1) mcMotherPdg1 = mother1->GetPdgCode();
1213  grandmotherId1 = mother1->GetMotherId();
1214  if (grandmotherId1 != -1) {
1215  CbmMCTrack* grandmother1 =
1216  (CbmMCTrack*) fMcTracks->At(grandmotherId1);
1217  if (NULL != grandmother1)
1218  mcGrandmotherPdg1 = grandmother1->GetPdgCode();
1219  }
1220  }
1221  if (motherId2 != -1) {
1222  CbmMCTrack* mother2 = (CbmMCTrack*) fMcTracks->At(motherId2);
1223  if (NULL != mother2) mcMotherPdg2 = mother2->GetPdgCode();
1224  grandmotherId2 = mother2->GetMotherId();
1225  if (grandmotherId2 != -1) {
1226  // CbmMCTrack* grandmother2 = (CbmMCTrack*) fMcTracks->At(grandmotherId2);
1227  // if (NULL != grandmother2) mcGrandmotherPdg2 = grandmother2->GetPdgCode();
1228  }
1229  }
1230  if (motherId3 != -1) {
1231  CbmMCTrack* mother3 = (CbmMCTrack*) fMcTracks->At(motherId3);
1232  if (NULL != mother3) mcMotherPdg3 = mother3->GetPdgCode();
1233  grandmotherId3 = mother3->GetMotherId();
1234  if (grandmotherId3 != -1) {
1235  // CbmMCTrack* grandmother3 = (CbmMCTrack*) fMcTracks->At(grandmotherId3);
1236  // if (NULL != grandmother3) mcGrandmotherPdg3 = grandmother3->GetPdgCode();
1237  }
1238  }
1239  if (motherId4 != -1) {
1240  CbmMCTrack* mother4 = (CbmMCTrack*) fMcTracks->At(motherId4);
1241  // if (NULL != mother4) mcMotherPdg4 = mother4->GetPdgCode();
1242  grandmotherId4 = mother4->GetMotherId();
1243  if (grandmotherId4 != -1) {
1244  // CbmMCTrack* grandmother4 = (CbmMCTrack*) fMcTracks->At(grandmotherId4);
1245  // if (NULL != grandmother4) mcGrandmotherPdg4 = grandmother4->GetPdgCode();
1246  }
1247  }
1248 
1249 
1250  if (motherId1 == motherId2 && motherId1 == motherId3) {}
1251  if (motherId1 == motherId2 && motherId1 == motherId4) {}
1252  if (motherId1 == motherId3 && motherId1 == motherId4) {}
1253  if (motherId2 == motherId3 && motherId2 == motherId4) {}
1254 
1255  // if( ((motherId1 == motherId2 && motherId3 == motherId4) && (mcMotherPdg1 == 22 || mcMotherPdg1 == 111) && (mcMotherPdg3 == 22 || mcMotherPdg3 == 111))
1256  // || ((motherId1 == motherId3 && motherId2 == motherId4) && (mcMotherPdg1 == 22 || mcMotherPdg1 == 111) && (mcMotherPdg2 == 22 || mcMotherPdg2 == 111))
1257  // || ((motherId1 == motherId4 && motherId2 == motherId3) && (mcMotherPdg1 == 22 || mcMotherPdg1 == 111) && (mcMotherPdg2 == 22 || mcMotherPdg2 == 111))) {
1258 
1259  // decay X -> gg -> e+ e- e+ e-
1260  if (((motherId1 == motherId2 && motherId3 == motherId4)
1261  && (mcMotherPdg1 == 22) && (mcMotherPdg3 == 22)
1262  && grandmotherId1 == grandmotherId3)
1263  || ((motherId1 == motherId3 && motherId2 == motherId4)
1264  && (mcMotherPdg1 == 22) && (mcMotherPdg2 == 22)
1265  && grandmotherId1 == grandmotherId2)
1266  || ((motherId1 == motherId4 && motherId2 == motherId3)
1267  && (mcMotherPdg1 == 22) && (mcMotherPdg2 == 22)
1268  && grandmotherId1 == grandmotherId2)) {
1269  Double_t invmass = Invmass_4particles(fMCTracklist_all[i],
1270  fMCTracklist_all[j],
1271  fMCTracklist_all[k],
1272  fMCTracklist_all[l]);
1273  //fhMCtest->Fill(invmass);
1274  fhInvariantMass_MC_all->Fill(invmass);
1275  if (mcGrandmotherPdg1 == 111) {
1276  fhInvariantMass_MC_pi0->Fill(invmass);
1277  fhInvariantMass_MC_pi0_gg->Fill(invmass);
1278  fhMC_electrons_theta->Fill(momentum1.Theta() * 180.
1279  / TMath::Pi());
1280  fhMC_electrons_theta->Fill(momentum2.Theta() * 180.
1281  / TMath::Pi());
1282  fhMC_electrons_theta->Fill(momentum3.Theta() * 180.
1283  / TMath::Pi());
1284  fhMC_electrons_theta->Fill(momentum4.Theta() * 180.
1285  / TMath::Pi());
1286  fhMC_electrons_p->Fill(momentum1.Mag());
1287  fhMC_electrons_p->Fill(momentum2.Mag());
1288  fhMC_electrons_p->Fill(momentum3.Mag());
1289  fhMC_electrons_p->Fill(momentum4.Mag());
1290  fhMC_electrons_theta_vs_p->Fill(
1291  momentum1.Theta() * 180. / TMath::Pi(), momentum1.Mag());
1292  fhMC_electrons_theta_vs_p->Fill(
1293  momentum2.Theta() * 180. / TMath::Pi(), momentum2.Mag());
1294  fhMC_electrons_theta_vs_p->Fill(
1295  momentum3.Theta() * 180. / TMath::Pi(), momentum3.Mag());
1296  fhMC_electrons_theta_vs_p->Fill(
1297  momentum4.Theta() * 180. / TMath::Pi(), momentum4.Mag());
1298  }
1299  if (mcGrandmotherPdg1 == 221) {
1300  fhInvariantMass_MC_eta->Fill(invmass);
1301 
1302  Double_t opening_angle_gg = 0;
1303  if (motherId1 == motherId2) {
1304  opening_angle_gg =
1305  CbmAnaConversionGlobalFunctions::OpeningAngleBetweenGamma(
1306  momentum1, momentum2, momentum3, momentum4);
1307  }
1308  if (motherId1 == motherId3) {
1309  opening_angle_gg =
1310  CbmAnaConversionGlobalFunctions::OpeningAngleBetweenGamma(
1311  momentum1, momentum3, momentum2, momentum4);
1312  }
1313  if (motherId1 == motherId4) {
1314  opening_angle_gg =
1315  CbmAnaConversionGlobalFunctions::OpeningAngleBetweenGamma(
1316  momentum1, momentum4, momentum2, momentum3);
1317  }
1318 
1319  fhEta_openingAngleGG->Fill(opening_angle_gg);
1320  }
1321  if (TMath::Abs(mcGrandmotherPdg1) == 331) { // eta prime (958)
1322  fhInvariantMass_MC_etaPrime->Fill(invmass);
1323  }
1324  fhMC_grandmotherPDGs->Fill(TMath::Abs(mcGrandmotherPdg1));
1325  }
1326 
1327 
1328  // decay pi0 -> g e+ e- -> e+ e- e+ e-
1329  if (((motherId1 == motherId2 && motherId3 == motherId4)
1330  && (mcMotherPdg1 == 22) && (mcMotherPdg3 == 111)
1331  && grandmotherId1 == motherId3)
1332  || ((motherId1 == motherId2 && motherId3 == motherId4)
1333  && (mcMotherPdg1 == 111) && (mcMotherPdg3 == 22)
1334  && grandmotherId3 == motherId1)
1335  || ((motherId1 == motherId3 && motherId2 == motherId4)
1336  && (mcMotherPdg1 == 22) && (mcMotherPdg2 == 111)
1337  && grandmotherId1 == motherId2)
1338  || ((motherId1 == motherId3 && motherId2 == motherId4)
1339  && (mcMotherPdg1 == 111) && (mcMotherPdg2 == 22)
1340  && grandmotherId2 == motherId1)
1341  || ((motherId1 == motherId4 && motherId2 == motherId3)
1342  && (mcMotherPdg1 == 22) && (mcMotherPdg2 == 111)
1343  && grandmotherId1 == motherId2)
1344  || ((motherId1 == motherId4 && motherId2 == motherId3)
1345  && (mcMotherPdg1 == 111) && (mcMotherPdg2 == 22)
1346  && grandmotherId2 == motherId1)) {
1347  Double_t invmass = Invmass_4particles(fMCTracklist_all[i],
1348  fMCTracklist_all[j],
1349  fMCTracklist_all[k],
1350  fMCTracklist_all[l]);
1351  fhInvariantMass_MC_pi0_gepem->Fill(invmass);
1352  fhInvariantMass_MC_pi0->Fill(invmass);
1353 
1354  fhMC_electrons_theta->Fill(momentum1.Theta() * 180.
1355  / TMath::Pi());
1356  fhMC_electrons_theta->Fill(momentum2.Theta() * 180.
1357  / TMath::Pi());
1358  fhMC_electrons_theta->Fill(momentum3.Theta() * 180.
1359  / TMath::Pi());
1360  fhMC_electrons_theta->Fill(momentum4.Theta() * 180.
1361  / TMath::Pi());
1362  fhMC_electrons_p->Fill(momentum1.Mag());
1363  fhMC_electrons_p->Fill(momentum2.Mag());
1364  fhMC_electrons_p->Fill(momentum3.Mag());
1365  fhMC_electrons_p->Fill(momentum4.Mag());
1366  fhMC_electrons_theta_vs_p->Fill(
1367  momentum1.Theta() * 180. / TMath::Pi(), momentum1.Mag());
1368  fhMC_electrons_theta_vs_p->Fill(
1369  momentum2.Theta() * 180. / TMath::Pi(), momentum2.Mag());
1370  fhMC_electrons_theta_vs_p->Fill(
1371  momentum3.Theta() * 180. / TMath::Pi(), momentum3.Mag());
1372  fhMC_electrons_theta_vs_p->Fill(
1373  momentum4.Theta() * 180. / TMath::Pi(), momentum4.Mag());
1374  }
1375  // decay eta -> g e+ e- -> e+ e- e+ e-
1376  if (((motherId1 == motherId2 && motherId3 == motherId4)
1377  && (mcMotherPdg1 == 22) && (mcMotherPdg3 == 221)
1378  && grandmotherId1 == motherId3)
1379  || ((motherId1 == motherId2 && motherId3 == motherId4)
1380  && (mcMotherPdg1 == 221) && (mcMotherPdg3 == 22)
1381  && grandmotherId3 == motherId1)
1382  || ((motherId1 == motherId3 && motherId2 == motherId4)
1383  && (mcMotherPdg1 == 22) && (mcMotherPdg2 == 221)
1384  && grandmotherId1 == motherId2)
1385  || ((motherId1 == motherId3 && motherId2 == motherId4)
1386  && (mcMotherPdg1 == 221) && (mcMotherPdg2 == 22)
1387  && grandmotherId2 == motherId1)
1388  || ((motherId1 == motherId4 && motherId2 == motherId3)
1389  && (mcMotherPdg1 == 22) && (mcMotherPdg2 == 221)
1390  && grandmotherId1 == motherId2)
1391  || ((motherId1 == motherId4 && motherId2 == motherId3)
1392  && (mcMotherPdg1 == 221) && (mcMotherPdg2 == 22)
1393  && grandmotherId2 == motherId1)) {
1394  Double_t invmass = Invmass_4particles(fMCTracklist_all[i],
1395  fMCTracklist_all[j],
1396  fMCTracklist_all[k],
1397  fMCTracklist_all[l]);
1398  fhInvariantMass_MC_eta->Fill(invmass);
1399  }
1400  }
1401  }
1402  }
1403  }
1404  }
1405  cout << "CbmAnaConversionReco: InvariantMassTestMC - End!" << endl;
1406 
1407  timer.Stop();
1408  fTime += timer.RealTime();
1409 }
1410 
1411 
1412 Double_t CbmAnaConversionReco::Invmass_4particles(const CbmMCTrack* mctrack1,
1413  const CbmMCTrack* mctrack2,
1414  const CbmMCTrack* mctrack3,
1415  const CbmMCTrack* mctrack4)
1416 // calculation of invariant mass from four electrons/positrons
1417 {
1418  /* TVector3 mom1;
1419  mctrack1->GetMomentum(mom1);
1420  TVector3 tempmom1;
1421  tempmom1.SetX(SmearValue(mom1.X()));
1422  tempmom1.SetY(SmearValue(mom1.Y()));
1423  tempmom1.SetZ(SmearValue(mom1.Z()));
1424  Double_t energy1 = TMath::Sqrt(tempmom1.Mag2() + M2E);
1425  TLorentzVector lorVec1(tempmom1, energy1);
1426 
1427  TVector3 mom2;
1428  mctrack2->GetMomentum(mom2);
1429  TVector3 tempmom2;
1430  tempmom2.SetX(SmearValue(mom2.X()));
1431  tempmom2.SetY(SmearValue(mom2.Y()));
1432  tempmom2.SetZ(SmearValue(mom2.Z()));
1433  Double_t energy2 = TMath::Sqrt(tempmom2.Mag2() + M2E);
1434  TLorentzVector lorVec2(tempmom2, energy2);
1435 
1436  TVector3 mom3;
1437  mctrack3->GetMomentum(mom3);
1438  TVector3 tempmom3;
1439  tempmom3.SetX(SmearValue(mom3.X()));
1440  tempmom3.SetY(SmearValue(mom3.Y()));
1441  tempmom3.SetZ(SmearValue(mom3.Z()));
1442  Double_t energy3 = TMath::Sqrt(tempmom3.Mag2() + M2E);
1443  TLorentzVector lorVec3(tempmom3, energy3);
1444 
1445  TVector3 mom4;
1446  mctrack4->GetMomentum(mom4);
1447  TVector3 tempmom4;
1448  tempmom4.SetX(SmearValue(mom4.X()));
1449  tempmom4.SetY(SmearValue(mom4.Y()));
1450  tempmom4.SetZ(SmearValue(mom4.Z()));
1451  Double_t energy4 = TMath::Sqrt(tempmom4.Mag2() + M2E);
1452  TLorentzVector lorVec4(tempmom4, energy4);
1453 */
1454  TLorentzVector lorVec1;
1455  mctrack1->Get4Momentum(lorVec1);
1456 
1457  TLorentzVector lorVec2;
1458  mctrack2->Get4Momentum(lorVec2);
1459 
1460  TLorentzVector lorVec3;
1461  mctrack3->Get4Momentum(lorVec3);
1462 
1463  TLorentzVector lorVec4;
1464  mctrack4->Get4Momentum(lorVec4);
1465 
1466 
1467  TLorentzVector sum;
1468  sum = lorVec1 + lorVec2 + lorVec3 + lorVec4;
1469  //cout << "mc: \t" << sum.Px() << " / " << sum.Py() << " / " << sum.Pz() << " / " << sum.E() << "\t => mag = " << sum.Mag() << endl;
1470 
1471 
1472  return sum.Mag();
1473 }
1474 
1475 
1476 Double_t CbmAnaConversionReco::SmearValue(Double_t value) {
1477  TRandom3 generator(0);
1478  Double_t result = 0;
1479  // Double_t smear = 0;
1480  Int_t plusminus = 0;
1481  while (plusminus == 0) { // should be either 1 or -1, not 0
1482  plusminus = generator.Uniform(-2, 2);
1483  }
1484  // Double_t gaus = generator.Gaus(1,1);
1485  // smear = gaus * plusminus;
1486  // result = value * (1. + 1.0*smear/100); //smearing as wished
1487 
1488  result = value; // -> no smearing
1489 
1490  return result;
1491 }
1492 
1493 
1494 Double_t CbmAnaConversionReco::Invmass_4particlesRECO(const TVector3 part1,
1495  const TVector3 part2,
1496  const TVector3 part3,
1497  const TVector3 part4)
1498 // calculation of invariant mass from four electrons/positrons
1499 {
1500  Double_t energy1 = TMath::Sqrt(part1.Mag2() + M2E);
1501  TLorentzVector lorVec1(part1, energy1);
1502 
1503  Double_t energy2 = TMath::Sqrt(part2.Mag2() + M2E);
1504  TLorentzVector lorVec2(part2, energy2);
1505 
1506  Double_t energy3 = TMath::Sqrt(part3.Mag2() + M2E);
1507  TLorentzVector lorVec3(part3, energy3);
1508 
1509  Double_t energy4 = TMath::Sqrt(part4.Mag2() + M2E);
1510  TLorentzVector lorVec4(part4, energy4);
1511 
1512  TLorentzVector sum;
1513  sum = lorVec1 + lorVec2 + lorVec3 + lorVec4;
1514  //cout << "reco: \t" << sum.Px() << " / " << sum.Py() << " / " << sum.Pz() << " / " << sum.E() << "\t => mag = " << sum.Mag() << endl;
1515 
1516 
1517  return sum.Mag();
1518 }
1519 
1520 
1521 void CbmAnaConversionReco::InvariantMassTest_4epem()
1522 // Calculating invariant mass of 4 ep/em, using MC data AND reconstructed momentum
1523 {
1524  timer.Start();
1525 
1526 
1527  cout << "CbmAnaConversionReco: InvariantMassTest_4epem - Start..." << endl;
1528  cout << "CbmAnaConversionReco: InvariantMassTest_4epem - "
1529  << fRecoTracklistEPEM.size() << "\t" << fRecoMomentum.size() << endl;
1530  int fill = 0;
1531  if (fRecoTracklistEPEM.size() < 4) return;
1532  for (unsigned int i = 0; i < fRecoTracklistEPEM.size(); i++) {
1533  if (i % 10 == 0)
1534  cout << "CbmAnaConversionReco: InvariantMassTest_4epem - iteration i = "
1535  << i << endl;
1536  for (unsigned int j = i + 1; j < fRecoTracklistEPEM.size(); j++) {
1537  for (unsigned int k = j + 1; k < fRecoTracklistEPEM.size(); k++) {
1538  for (unsigned int l = k + 1; l < fRecoTracklistEPEM.size(); l++) {
1539  if (fRecoTracklistEPEM[i]->GetPdgCode()
1540  + fRecoTracklistEPEM[j]->GetPdgCode()
1541  + fRecoTracklistEPEM[k]->GetPdgCode()
1542  + fRecoTracklistEPEM[l]->GetPdgCode()
1543  != 0)
1544  continue;
1545 
1546  if (fRecoTracklistEPEM.size() != fRecoMomentum.size()
1547  || fRecoTracklistEPEM.size() != fRecoRefittedMomentum.size()) {
1548  cout << "CbmAnaConversionReco: InvariantMassTest_4epem - not "
1549  "matching number of entries!"
1550  << endl;
1551  continue;
1552  }
1553 
1554 
1555  // starting points of each electron (-> i.e. conversion points of gamma OR decay points of pi0, depending on decay channel)
1556  TVector3 pi0start_i;
1557  fRecoTracklistEPEM[i]->GetStartVertex(pi0start_i);
1558  TVector3 pi0start_j;
1559  fRecoTracklistEPEM[j]->GetStartVertex(pi0start_j);
1560  TVector3 pi0start_k;
1561  fRecoTracklistEPEM[k]->GetStartVertex(pi0start_k);
1562  TVector3 pi0start_l;
1563  fRecoTracklistEPEM[l]->GetStartVertex(pi0start_l);
1564 
1565 
1566  int motherId1 = fRecoTracklistEPEM[i]->GetMotherId();
1567  int motherId2 = fRecoTracklistEPEM[j]->GetMotherId();
1568  int motherId3 = fRecoTracklistEPEM[k]->GetMotherId();
1569  int motherId4 = fRecoTracklistEPEM[l]->GetMotherId();
1570 
1571  // decay pi0 -> e+ e- e+ e-
1572  if (motherId1 == motherId2 && motherId1 == motherId3
1573  && motherId1 == motherId4) {
1574  if (motherId1 != -1) {
1575  int mcMotherPdg1 = -1;
1576  CbmMCTrack* mother1 = (CbmMCTrack*) fMcTracks->At(motherId1);
1577  if (NULL != mother1) mcMotherPdg1 = mother1->GetPdgCode();
1578  if ((mcMotherPdg1 == 111
1579  || mcMotherPdg1
1580  == 221)) { // && NofDaughters(motherId1) == 4) {
1581  Double_t invmass1 = 0;
1582  // Double_t invmass2 = 0; // momenta from stsMomentumVec
1583  Double_t invmass3 = 0;
1584  // invmass2 = Invmass_4particlesRECO(fRecoMomentum[i], fRecoMomentum[j], fRecoMomentum[k], fRecoMomentum[l]);
1585  invmass1 =
1586  Invmass_4particles(fRecoTracklistEPEM[i],
1587  fRecoTracklistEPEM[j],
1588  fRecoTracklistEPEM[k],
1589  fRecoTracklistEPEM[l]); // true MC values
1590  invmass3 = Invmass_4particlesRECO(fRecoRefittedMomentum[i],
1591  fRecoRefittedMomentum[j],
1592  fRecoRefittedMomentum[k],
1593  fRecoRefittedMomentum[l]);
1594  // fhInvariantMass_pi0epem->Fill(invmass2);
1595  cout << "Decay pi0 -> e+e-e+e- detected!\t\t mc mass: "
1596  << invmass1 << "\t, reco mass: " << invmass3 << endl;
1597  cout << "motherids: " << motherId1 << "/" << motherId2 << "/"
1598  << motherId3 << "/" << motherId4
1599  << "\t motherpdg: " << mcMotherPdg1
1600  << "\t mctrack mass: " << mother1->GetMass()
1601  << "\t nofdaughters: " << NofDaughters(motherId1) << endl;
1602  cout << "pdgs " << fRecoTracklistEPEM[i]->GetPdgCode() << "/"
1603  << fRecoTracklistEPEM[j]->GetPdgCode() << "/"
1604  << fRecoTracklistEPEM[k]->GetPdgCode() << "/"
1605  << fRecoTracklistEPEM[l]->GetPdgCode() << endl;
1606  TVector3 start1;
1607  fRecoTracklistEPEM[i]->GetStartVertex(start1);
1608  TVector3 start2;
1609  fRecoTracklistEPEM[j]->GetStartVertex(start2);
1610  TVector3 start3;
1611  fRecoTracklistEPEM[k]->GetStartVertex(start3);
1612  TVector3 start4;
1613  fRecoTracklistEPEM[l]->GetStartVertex(start4);
1614  cout << "start: " << start1.Z() << "/" << start2.Z() << "/"
1615  << start3.Z() << "/" << start4.Z() << endl;
1616 
1617  fhEPEM_invmass_eeee_mc->Fill(invmass1);
1618  fhEPEM_invmass_eeee_refitted->Fill(invmass3);
1619  //fhEPEM_invmass_all_mc->Fill(invmass1);
1620  //fhEPEM_invmass_all_refitted->Fill(invmass3);
1621  }
1622  } else { // all 4 particles come directly from the vertex
1623  continue;
1624  }
1625  }
1626 
1627  if ((motherId1 == motherId2 && motherId3 == motherId4)
1628  || (motherId1 == motherId3 && motherId2 == motherId4)
1629  || (motherId1 == motherId4 && motherId2 == motherId3)) {
1630 
1631 
1632  int grandmotherId1 = -1;
1633  int grandmotherId2 = -1;
1634  int grandmotherId3 = -1;
1635  int grandmotherId4 = -1;
1636 
1637  int mcMotherPdg1 = -1;
1638  int mcMotherPdg2 = -1;
1639  int mcMotherPdg3 = -1;
1640  int mcMotherPdg4 = -1;
1641  int mcGrandmotherPdg1 = -1;
1642  // int mcGrandmotherPdg2 = -1;
1643  // int mcGrandmotherPdg3 = -1;
1644  // int mcGrandmotherPdg4 = -1;
1645 
1646  CbmMCTrack* grandmother1 = NULL;
1647 
1648  if (motherId1 != -1) {
1649  CbmMCTrack* mother1 = (CbmMCTrack*) fMcTracks->At(motherId1);
1650  if (NULL != mother1) mcMotherPdg1 = mother1->GetPdgCode();
1651  grandmotherId1 = mother1->GetMotherId();
1652  if (grandmotherId1 != -1) {
1653  grandmother1 = (CbmMCTrack*) fMcTracks->At(grandmotherId1);
1654  if (NULL != grandmother1)
1655  mcGrandmotherPdg1 = grandmother1->GetPdgCode();
1656  }
1657  }
1658  if (motherId2 != -1) {
1659  CbmMCTrack* mother2 = (CbmMCTrack*) fMcTracks->At(motherId2);
1660  if (NULL != mother2) mcMotherPdg2 = mother2->GetPdgCode();
1661  grandmotherId2 = mother2->GetMotherId();
1662  if (grandmotherId2 != -1) {
1663  // CbmMCTrack* grandmother2 = (CbmMCTrack*) fMcTracks->At(grandmotherId2);
1664  // if (NULL != grandmother2) mcGrandmotherPdg2 = grandmother2->GetPdgCode();
1665  }
1666  }
1667  if (motherId3 != -1) {
1668  CbmMCTrack* mother3 = (CbmMCTrack*) fMcTracks->At(motherId3);
1669  if (NULL != mother3) mcMotherPdg3 = mother3->GetPdgCode();
1670  grandmotherId3 = mother3->GetMotherId();
1671  if (grandmotherId3 != -1) {
1672  // CbmMCTrack* grandmother3 = (CbmMCTrack*) fMcTracks->At(grandmotherId3);
1673  // if (NULL != grandmother3) mcGrandmotherPdg3 = grandmother3->GetPdgCode();
1674  }
1675  }
1676  if (motherId4 != -1) {
1677  CbmMCTrack* mother4 = (CbmMCTrack*) fMcTracks->At(motherId4);
1678  if (NULL != mother4) mcMotherPdg4 = mother4->GetPdgCode();
1679  grandmotherId4 = mother4->GetMotherId();
1680  if (grandmotherId4 != -1) {
1681  // CbmMCTrack* grandmother4 = (CbmMCTrack*) fMcTracks->At(grandmotherId4);
1682  // if (NULL != grandmother4) mcGrandmotherPdg4 = grandmother4->GetPdgCode();
1683  }
1684  }
1685 
1686  //if(grandmotherId1 == -1) continue;
1687 
1688  if (motherId1 == motherId2 && motherId3 == motherId4) {
1689  if (NofDaughters(motherId1) != 2 || NofDaughters(motherId3) != 2)
1690  continue;
1691  if ((grandmotherId1 == motherId3 && mcMotherPdg3 == 111)
1692  || (motherId1 == grandmotherId3 && mcMotherPdg1 == 111)) {
1693 
1694  fhPi0_startvertexElectrons_gee->Fill(pi0start_i.Z());
1695  fhPi0_startvertexElectrons_gee->Fill(pi0start_j.Z());
1696  fhPi0_startvertexElectrons_gee->Fill(pi0start_k.Z());
1697  fhPi0_startvertexElectrons_gee->Fill(pi0start_l.Z());
1698  fhPi0_startvertexElectrons_all->Fill(pi0start_i.Z());
1699  fhPi0_startvertexElectrons_all->Fill(pi0start_j.Z());
1700  fhPi0_startvertexElectrons_all->Fill(pi0start_k.Z());
1701  fhPi0_startvertexElectrons_all->Fill(pi0start_l.Z());
1702 
1703  Double_t invmass1 = 0;
1704  Double_t invmass3 = 0;
1705  invmass1 =
1706  Invmass_4particles(fRecoTracklistEPEM[i],
1707  fRecoTracklistEPEM[j],
1708  fRecoTracklistEPEM[k],
1709  fRecoTracklistEPEM[l]); // true MC values
1710  invmass3 = Invmass_4particlesRECO(fRecoRefittedMomentum[i],
1711  fRecoRefittedMomentum[j],
1712  fRecoRefittedMomentum[k],
1713  fRecoRefittedMomentum[l]);
1714 
1715  // consider only electrons from the target (only then the momenta are correctly refitted/reconstructed)
1716  if (pi0start_i.Z() > 1 || pi0start_j.Z() > 1
1717  || pi0start_k.Z() > 1 || pi0start_l.Z() > 1) {
1718  fhEPEM_InDetector_invmass_gee_mc->Fill(invmass1);
1719  fhEPEM_InDetector_invmass_gee_refitted->Fill(invmass3);
1720  fhEPEM_InDetector_invmass_all_mc->Fill(invmass1);
1721  fhEPEM_InDetector_invmass_all_refitted->Fill(invmass3);
1722  continue;
1723  }
1724 
1725 
1726  cout << "HURRAY! .-.-.-.-.-.-..-.-.-.-.-.-.-.-.-.-.-.-."
1727  << endl;
1728  CutEfficiencyStudies(
1729  i, j, k, l, motherId1, motherId2, motherId3, motherId4);
1730  fhInvariantMass_pi0epem->Fill(invmass1);
1731 
1732  fhEPEM_invmass_gee_mc->Fill(invmass1);
1733  fhEPEM_invmass_gee_refitted->Fill(invmass3);
1734  fhEPEM_invmass_all_mc->Fill(invmass1);
1735  fhEPEM_invmass_all_refitted->Fill(invmass3);
1736 
1737  fhEPEM_pt_vs_p_all_mc->Fill(fRecoTracklistEPEM[i]->GetPt(),
1738  fRecoTracklistEPEM[i]->GetP());
1739  fhEPEM_pt_vs_p_all_mc->Fill(fRecoTracklistEPEM[j]->GetPt(),
1740  fRecoTracklistEPEM[j]->GetP());
1741  fhEPEM_pt_vs_p_all_mc->Fill(fRecoTracklistEPEM[k]->GetPt(),
1742  fRecoTracklistEPEM[k]->GetP());
1743  fhEPEM_pt_vs_p_all_mc->Fill(fRecoTracklistEPEM[l]->GetPt(),
1744  fRecoTracklistEPEM[l]->GetP());
1745  fhEPEM_pt_vs_p_all_refitted->Fill(
1746  fRecoRefittedMomentum[i].Perp(),
1747  fRecoRefittedMomentum[i].Mag());
1748  fhEPEM_pt_vs_p_all_refitted->Fill(
1749  fRecoRefittedMomentum[j].Perp(),
1750  fRecoRefittedMomentum[j].Mag());
1751  fhEPEM_pt_vs_p_all_refitted->Fill(
1752  fRecoRefittedMomentum[k].Perp(),
1753  fRecoRefittedMomentum[k].Mag());
1754  fhEPEM_pt_vs_p_all_refitted->Fill(
1755  fRecoRefittedMomentum[l].Perp(),
1756  fRecoRefittedMomentum[l].Mag());
1757 
1758  CbmLmvmKinematicParams params1 =
1759  CalculateKinematicParams_4particles(fRecoRefittedMomentum[i],
1760  fRecoRefittedMomentum[j],
1761  fRecoRefittedMomentum[k],
1762  fRecoRefittedMomentum[l]);
1763  fhPi0_pt_vs_rap_gee->Fill(params1.fPt, params1.fRapidity);
1764  fhPi0_pt_vs_rap_all->Fill(params1.fPt, params1.fRapidity);
1765  fhPi0_pt_gee->Fill(params1.fPt);
1766  fhPi0_pt_all->Fill(params1.fPt);
1767 
1768  if (
1769  mcGrandmotherPdg1
1770  == 111) { // case: i,j = electrons from gamma, k,l = electrons from pi0
1771  Double_t opening_angle1_mc = 0;
1772  Double_t opening_angle1_refitted = 0;
1773  opening_angle1_mc = CalculateOpeningAngleMC(
1774  fRecoTracklistEPEM[i], fRecoTracklistEPEM[j]);
1775  opening_angle1_refitted = CalculateOpeningAngleReco(
1776  fRecoRefittedMomentum[i], fRecoRefittedMomentum[j]);
1777  fhEPEM_openingAngle_gee_mc->Fill(opening_angle1_mc);
1778  fhEPEM_openingAngle_gee_refitted->Fill(
1779  opening_angle1_refitted);
1780 
1781  Double_t opening_angle1_mc_dalitz = 0;
1782  Double_t opening_angle1_refitted_dalitz = 0;
1783  opening_angle1_mc_dalitz = CalculateOpeningAngleMC(
1784  fRecoTracklistEPEM[k], fRecoTracklistEPEM[l]);
1785  opening_angle1_refitted_dalitz = CalculateOpeningAngleReco(
1786  fRecoRefittedMomentum[k], fRecoRefittedMomentum[l]);
1787  fhEPEM_openingAngle_gee_mc_dalitz->Fill(
1788  opening_angle1_mc_dalitz);
1789  fhEPEM_openingAngle_gee_refitted_dalitz->Fill(
1790  opening_angle1_refitted_dalitz);
1791  }
1792 
1793  if (
1794  mcMotherPdg1
1795  == 111) { // case: i,j = electrons from pi0, k,l = electrons from gamma
1796  Double_t opening_angle1_mc = 0;
1797  Double_t opening_angle1_refitted = 0;
1798  opening_angle1_mc = CalculateOpeningAngleMC(
1799  fRecoTracklistEPEM[k], fRecoTracklistEPEM[l]);
1800  opening_angle1_refitted = CalculateOpeningAngleReco(
1801  fRecoRefittedMomentum[k], fRecoRefittedMomentum[l]);
1802  fhEPEM_openingAngle_gee_mc->Fill(opening_angle1_mc);
1803  fhEPEM_openingAngle_gee_refitted->Fill(
1804  opening_angle1_refitted);
1805 
1806  Double_t opening_angle1_mc_dalitz = 0;
1807  Double_t opening_angle1_refitted_dalitz = 0;
1808  opening_angle1_mc_dalitz = CalculateOpeningAngleMC(
1809  fRecoTracklistEPEM[i], fRecoTracklistEPEM[j]);
1810  opening_angle1_refitted_dalitz = CalculateOpeningAngleReco(
1811  fRecoRefittedMomentum[i], fRecoRefittedMomentum[j]);
1812  fhEPEM_openingAngle_gee_mc_dalitz->Fill(
1813  opening_angle1_mc_dalitz);
1814  fhEPEM_openingAngle_gee_refitted_dalitz->Fill(
1815  opening_angle1_refitted_dalitz);
1816  }
1817  }
1818  }
1819  if (motherId1 == motherId3 && motherId2 == motherId4) {
1820  if (NofDaughters(motherId1) != 2 || NofDaughters(motherId2) != 2)
1821  continue;
1822  if ((grandmotherId1 == motherId2 && mcMotherPdg2 == 111)
1823  || (motherId1 == grandmotherId2 && mcMotherPdg1 == 111)) {
1824 
1825  fhPi0_startvertexElectrons_gee->Fill(pi0start_i.Z());
1826  fhPi0_startvertexElectrons_gee->Fill(pi0start_j.Z());
1827  fhPi0_startvertexElectrons_gee->Fill(pi0start_k.Z());
1828  fhPi0_startvertexElectrons_gee->Fill(pi0start_l.Z());
1829  fhPi0_startvertexElectrons_all->Fill(pi0start_i.Z());
1830  fhPi0_startvertexElectrons_all->Fill(pi0start_j.Z());
1831  fhPi0_startvertexElectrons_all->Fill(pi0start_k.Z());
1832  fhPi0_startvertexElectrons_all->Fill(pi0start_l.Z());
1833 
1834  Double_t invmass1 = 0;
1835  Double_t invmass3 = 0;
1836  invmass1 =
1837  Invmass_4particles(fRecoTracklistEPEM[i],
1838  fRecoTracklistEPEM[j],
1839  fRecoTracklistEPEM[k],
1840  fRecoTracklistEPEM[l]); // true MC values
1841  invmass3 = Invmass_4particlesRECO(fRecoRefittedMomentum[i],
1842  fRecoRefittedMomentum[j],
1843  fRecoRefittedMomentum[k],
1844  fRecoRefittedMomentum[l]);
1845 
1846  // consider only electrons from the target (only then the momenta are correctly refitted/reconstructed)
1847  if (pi0start_i.Z() > 1 || pi0start_j.Z() > 1
1848  || pi0start_k.Z() > 1 || pi0start_l.Z() > 1) {
1849  fhEPEM_InDetector_invmass_gee_mc->Fill(invmass1);
1850  fhEPEM_InDetector_invmass_gee_refitted->Fill(invmass3);
1851  fhEPEM_InDetector_invmass_all_mc->Fill(invmass1);
1852  fhEPEM_InDetector_invmass_all_refitted->Fill(invmass3);
1853  continue;
1854  }
1855 
1856 
1857  cout << "HURRAY! .-.-.-.-.-.-..-.-.-.-.-.-.-.-.-.-.-.-."
1858  << endl;
1859  CutEfficiencyStudies(
1860  i, j, k, l, motherId1, motherId2, motherId3, motherId4);
1861  fhInvariantMass_pi0epem->Fill(invmass1);
1862 
1863  fhEPEM_invmass_gee_mc->Fill(invmass1);
1864  fhEPEM_invmass_gee_refitted->Fill(invmass3);
1865  fhEPEM_invmass_all_mc->Fill(invmass1);
1866  fhEPEM_invmass_all_refitted->Fill(invmass3);
1867 
1868  fhEPEM_pt_vs_p_all_mc->Fill(fRecoTracklistEPEM[i]->GetPt(),
1869  fRecoTracklistEPEM[i]->GetP());
1870  fhEPEM_pt_vs_p_all_mc->Fill(fRecoTracklistEPEM[j]->GetPt(),
1871  fRecoTracklistEPEM[j]->GetP());
1872  fhEPEM_pt_vs_p_all_mc->Fill(fRecoTracklistEPEM[k]->GetPt(),
1873  fRecoTracklistEPEM[k]->GetP());
1874  fhEPEM_pt_vs_p_all_mc->Fill(fRecoTracklistEPEM[l]->GetPt(),
1875  fRecoTracklistEPEM[l]->GetP());
1876  fhEPEM_pt_vs_p_all_refitted->Fill(
1877  fRecoRefittedMomentum[i].Perp(),
1878  fRecoRefittedMomentum[i].Mag());
1879  fhEPEM_pt_vs_p_all_refitted->Fill(
1880  fRecoRefittedMomentum[j].Perp(),
1881  fRecoRefittedMomentum[j].Mag());
1882  fhEPEM_pt_vs_p_all_refitted->Fill(
1883  fRecoRefittedMomentum[k].Perp(),
1884  fRecoRefittedMomentum[k].Mag());
1885  fhEPEM_pt_vs_p_all_refitted->Fill(
1886  fRecoRefittedMomentum[l].Perp(),
1887  fRecoRefittedMomentum[l].Mag());
1888 
1889  CbmLmvmKinematicParams params1 =
1890  CalculateKinematicParams_4particles(fRecoRefittedMomentum[i],
1891  fRecoRefittedMomentum[j],
1892  fRecoRefittedMomentum[k],
1893  fRecoRefittedMomentum[l]);
1894  fhPi0_pt_vs_rap_gee->Fill(params1.fPt, params1.fRapidity);
1895  fhPi0_pt_vs_rap_all->Fill(params1.fPt, params1.fRapidity);
1896  fhPi0_pt_gee->Fill(params1.fPt);
1897  fhPi0_pt_all->Fill(params1.fPt);
1898 
1899  if (mcGrandmotherPdg1 == 111) {
1900  Double_t opening_angle1_mc = 0;
1901  Double_t opening_angle1_refitted = 0;
1902  opening_angle1_mc = CalculateOpeningAngleMC(
1903  fRecoTracklistEPEM[i], fRecoTracklistEPEM[k]);
1904  opening_angle1_refitted = CalculateOpeningAngleReco(
1905  fRecoRefittedMomentum[i], fRecoRefittedMomentum[k]);
1906  fhEPEM_openingAngle_gee_mc->Fill(opening_angle1_mc);
1907  fhEPEM_openingAngle_gee_refitted->Fill(
1908  opening_angle1_refitted);
1909 
1910  Double_t opening_angle1_mc_dalitz = 0;
1911  Double_t opening_angle1_refitted_dalitz = 0;
1912  opening_angle1_mc_dalitz = CalculateOpeningAngleMC(
1913  fRecoTracklistEPEM[j], fRecoTracklistEPEM[l]);
1914  opening_angle1_refitted_dalitz = CalculateOpeningAngleReco(
1915  fRecoRefittedMomentum[j], fRecoRefittedMomentum[l]);
1916  fhEPEM_openingAngle_gee_mc_dalitz->Fill(
1917  opening_angle1_mc_dalitz);
1918  fhEPEM_openingAngle_gee_refitted_dalitz->Fill(
1919  opening_angle1_refitted_dalitz);
1920  }
1921 
1922  if (mcMotherPdg1 == 111) {
1923  Double_t opening_angle1_mc = 0;
1924  Double_t opening_angle1_refitted = 0;
1925  opening_angle1_mc = CalculateOpeningAngleMC(
1926  fRecoTracklistEPEM[j], fRecoTracklistEPEM[l]);
1927  opening_angle1_refitted = CalculateOpeningAngleReco(
1928  fRecoRefittedMomentum[j], fRecoRefittedMomentum[l]);
1929  fhEPEM_openingAngle_gee_mc->Fill(opening_angle1_mc);
1930  fhEPEM_openingAngle_gee_refitted->Fill(
1931  opening_angle1_refitted);
1932 
1933  Double_t opening_angle1_mc_dalitz = 0;
1934  Double_t opening_angle1_refitted_dalitz = 0;
1935  opening_angle1_mc_dalitz = CalculateOpeningAngleMC(
1936  fRecoTracklistEPEM[i], fRecoTracklistEPEM[k]);
1937  opening_angle1_refitted_dalitz = CalculateOpeningAngleReco(
1938  fRecoRefittedMomentum[i], fRecoRefittedMomentum[k]);
1939  fhEPEM_openingAngle_gee_mc_dalitz->Fill(
1940  opening_angle1_mc_dalitz);
1941  fhEPEM_openingAngle_gee_refitted_dalitz->Fill(
1942  opening_angle1_refitted_dalitz);
1943  }
1944  }
1945  }
1946  if (motherId1 == motherId4 && motherId2 == motherId3) {
1947  if (NofDaughters(motherId1) != 2 || NofDaughters(motherId2) != 2)
1948  continue;
1949  if ((grandmotherId1 == motherId2 && mcMotherPdg2 == 111)
1950  || (motherId1 == grandmotherId2 && mcMotherPdg1 == 111)) {
1951 
1952  fhPi0_startvertexElectrons_gee->Fill(pi0start_i.Z());
1953  fhPi0_startvertexElectrons_gee->Fill(pi0start_j.Z());
1954  fhPi0_startvertexElectrons_gee->Fill(pi0start_k.Z());
1955  fhPi0_startvertexElectrons_gee->Fill(pi0start_l.Z());
1956  fhPi0_startvertexElectrons_all->Fill(pi0start_i.Z());
1957  fhPi0_startvertexElectrons_all->Fill(pi0start_j.Z());
1958  fhPi0_startvertexElectrons_all->Fill(pi0start_k.Z());
1959  fhPi0_startvertexElectrons_all->Fill(pi0start_l.Z());
1960 
1961  Double_t invmass1 = 0;
1962  Double_t invmass3 = 0;
1963  invmass1 =
1964  Invmass_4particles(fRecoTracklistEPEM[i],
1965  fRecoTracklistEPEM[j],
1966  fRecoTracklistEPEM[k],
1967  fRecoTracklistEPEM[l]); // true MC values
1968  invmass3 = Invmass_4particlesRECO(fRecoRefittedMomentum[i],
1969  fRecoRefittedMomentum[j],
1970  fRecoRefittedMomentum[k],
1971  fRecoRefittedMomentum[l]);
1972 
1973  // consider only electrons from the target (only then the momenta are correctly refitted/reconstructed)
1974  if (pi0start_i.Z() > 1 || pi0start_j.Z() > 1
1975  || pi0start_k.Z() > 1 || pi0start_l.Z() > 1) {
1976  fhEPEM_InDetector_invmass_gee_mc->Fill(invmass1);
1977  fhEPEM_InDetector_invmass_gee_refitted->Fill(invmass3);
1978  fhEPEM_InDetector_invmass_all_mc->Fill(invmass1);
1979  fhEPEM_InDetector_invmass_all_refitted->Fill(invmass3);
1980  continue;
1981  }
1982 
1983 
1984  cout << "HURRAY! .-.-.-.-.-.-..-.-.-.-.-.-.-.-.-.-.-.-."
1985  << endl;
1986  CutEfficiencyStudies(
1987  i, j, k, l, motherId1, motherId2, motherId3, motherId4);
1988  fhInvariantMass_pi0epem->Fill(invmass1);
1989 
1990  fhEPEM_invmass_gee_mc->Fill(invmass1);
1991  fhEPEM_invmass_gee_refitted->Fill(invmass3);
1992  fhEPEM_invmass_all_mc->Fill(invmass1);
1993  fhEPEM_invmass_all_refitted->Fill(invmass3);
1994 
1995  fhEPEM_pt_vs_p_all_mc->Fill(fRecoTracklistEPEM[i]->GetPt(),
1996  fRecoTracklistEPEM[i]->GetP());
1997  fhEPEM_pt_vs_p_all_mc->Fill(fRecoTracklistEPEM[j]->GetPt(),
1998  fRecoTracklistEPEM[j]->GetP());
1999  fhEPEM_pt_vs_p_all_mc->Fill(fRecoTracklistEPEM[k]->GetPt(),
2000  fRecoTracklistEPEM[k]->GetP());
2001  fhEPEM_pt_vs_p_all_mc->Fill(fRecoTracklistEPEM[l]->GetPt(),
2002  fRecoTracklistEPEM[l]->GetP());
2003  fhEPEM_pt_vs_p_all_refitted->Fill(
2004  fRecoRefittedMomentum[i].Perp(),
2005  fRecoRefittedMomentum[i].Mag());
2006  fhEPEM_pt_vs_p_all_refitted->Fill(
2007  fRecoRefittedMomentum[j].Perp(),
2008  fRecoRefittedMomentum[j].Mag());
2009  fhEPEM_pt_vs_p_all_refitted->Fill(
2010  fRecoRefittedMomentum[k].Perp(),
2011  fRecoRefittedMomentum[k].Mag());
2012  fhEPEM_pt_vs_p_all_refitted->Fill(
2013  fRecoRefittedMomentum[l].Perp(),
2014  fRecoRefittedMomentum[l].Mag());
2015 
2016  CbmLmvmKinematicParams params1 =
2017  CalculateKinematicParams_4particles(fRecoRefittedMomentum[i],
2018  fRecoRefittedMomentum[j],
2019  fRecoRefittedMomentum[k],
2020  fRecoRefittedMomentum[l]);
2021  fhPi0_pt_vs_rap_gee->Fill(params1.fPt, params1.fRapidity);
2022  fhPi0_pt_vs_rap_all->Fill(params1.fPt, params1.fRapidity);
2023  fhPi0_pt_gee->Fill(params1.fPt);
2024  fhPi0_pt_all->Fill(params1.fPt);
2025 
2026  if (mcGrandmotherPdg1 == 111) {
2027  Double_t opening_angle1_mc = 0;
2028  Double_t opening_angle1_refitted = 0;
2029  opening_angle1_mc = CalculateOpeningAngleMC(
2030  fRecoTracklistEPEM[i], fRecoTracklistEPEM[l]);
2031  opening_angle1_refitted = CalculateOpeningAngleReco(
2032  fRecoRefittedMomentum[i], fRecoRefittedMomentum[l]);
2033  fhEPEM_openingAngle_gee_mc->Fill(opening_angle1_mc);
2034  fhEPEM_openingAngle_gee_refitted->Fill(
2035  opening_angle1_refitted);
2036 
2037  Double_t opening_angle1_mc_dalitz = 0;
2038  Double_t opening_angle1_refitted_dalitz = 0;
2039  opening_angle1_mc_dalitz = CalculateOpeningAngleMC(
2040  fRecoTracklistEPEM[j], fRecoTracklistEPEM[k]);
2041  opening_angle1_refitted_dalitz = CalculateOpeningAngleReco(
2042  fRecoRefittedMomentum[j], fRecoRefittedMomentum[k]);
2043  fhEPEM_openingAngle_gee_mc_dalitz->Fill(
2044  opening_angle1_mc_dalitz);
2045  fhEPEM_openingAngle_gee_refitted_dalitz->Fill(
2046  opening_angle1_refitted_dalitz);
2047  }
2048 
2049  if (mcMotherPdg1 == 111) {
2050  Double_t opening_angle1_mc = 0;
2051  Double_t opening_angle1_refitted = 0;
2052  opening_angle1_mc = CalculateOpeningAngleMC(
2053  fRecoTracklistEPEM[j], fRecoTracklistEPEM[k]);
2054  opening_angle1_refitted = CalculateOpeningAngleReco(
2055  fRecoRefittedMomentum[j], fRecoRefittedMomentum[k]);
2056  fhEPEM_openingAngle_gee_mc->Fill(opening_angle1_mc);
2057  fhEPEM_openingAngle_gee_refitted->Fill(
2058  opening_angle1_refitted);
2059 
2060  Double_t opening_angle1_mc_dalitz = 0;
2061  Double_t opening_angle1_refitted_dalitz = 0;
2062  opening_angle1_mc_dalitz = CalculateOpeningAngleMC(
2063  fRecoTracklistEPEM[i], fRecoTracklistEPEM[l]);
2064  opening_angle1_refitted_dalitz = CalculateOpeningAngleReco(
2065  fRecoRefittedMomentum[i], fRecoRefittedMomentum[l]);
2066  fhEPEM_openingAngle_gee_mc_dalitz->Fill(
2067  opening_angle1_mc_dalitz);
2068  fhEPEM_openingAngle_gee_refitted_dalitz->Fill(
2069  opening_angle1_refitted_dalitz);
2070  }
2071  }
2072  }
2073 
2074 
2075  // ===================================================================================================
2076  // HERE DECAY pi0 -> gamma gamma -> e+e- e+e-
2077  if (grandmotherId1 == grandmotherId2
2078  && grandmotherId1 == grandmotherId3
2079  && grandmotherId1 == grandmotherId4) {
2080  // if(mcGrandmotherPdg1 != 111 && mcGrandmotherPdg1 != 221) continue; // 111 = pi0, 221 = eta
2081 
2082 
2083  if (pi0start_i.Z() <= 1 && pi0start_j.Z() <= 1
2084  && pi0start_k.Z() <= 1 && pi0start_l.Z() <= 1
2085  && mcGrandmotherPdg1 == 221) {
2086  Double_t invmass_eta_mc =
2087  Invmass_4particles(fRecoTracklistEPEM[i],
2088  fRecoTracklistEPEM[j],
2089  fRecoTracklistEPEM[k],
2090  fRecoTracklistEPEM[l]);
2091  Double_t invmass_eta_reco =
2092  Invmass_4particlesRECO(fRecoRefittedMomentum[i],
2093  fRecoRefittedMomentum[j],
2094  fRecoRefittedMomentum[k],
2095  fRecoRefittedMomentum[l]);
2096 
2097  fhEPEM_invmass_eta_mc->Fill(invmass_eta_mc);
2098  fhEPEM_invmass_eta_refitted->Fill(invmass_eta_reco);
2099 
2100  CutEfficiencyStudies(
2101  i, j, k, l, motherId1, motherId2, motherId3, motherId4, 1);
2102  }
2103 
2104 
2105  if (mcGrandmotherPdg1 != 111) continue; // 111 = pi0, 221 = eta
2106 
2107  TVector3 pi0start;
2108  grandmother1->GetStartVertex(pi0start);
2109  fhPi0_startvertex->Fill(pi0start.Z());
2110 
2111  fhPi0_startvertexElectrons_gg->Fill(pi0start_i.Z());
2112  fhPi0_startvertexElectrons_gg->Fill(pi0start_j.Z());
2113  fhPi0_startvertexElectrons_gg->Fill(pi0start_k.Z());
2114  fhPi0_startvertexElectrons_gg->Fill(pi0start_l.Z());
2115  fhPi0_startvertexElectrons_all->Fill(pi0start_i.Z());
2116  fhPi0_startvertexElectrons_all->Fill(pi0start_j.Z());
2117  fhPi0_startvertexElectrons_all->Fill(pi0start_k.Z());
2118  fhPi0_startvertexElectrons_all->Fill(pi0start_l.Z());
2119 
2120  fhPi0_startvertex_vs_chi->Fill(pi0start_i.Z(),
2121  fRecoTracklistEPEM_chi[i]);
2122  fhPi0_startvertex_vs_chi->Fill(pi0start_j.Z(),
2123  fRecoTracklistEPEM_chi[j]);
2124  fhPi0_startvertex_vs_chi->Fill(pi0start_k.Z(),
2125  fRecoTracklistEPEM_chi[l]);
2126  fhPi0_startvertex_vs_chi->Fill(pi0start_l.Z(),
2127  fRecoTracklistEPEM_chi[l]);
2128 
2129  fhPi0_startvertex_vs_momentum->Fill(
2130  pi0start_i.Z(), fRecoTracklistEPEM[i]->GetP());
2131  fhPi0_startvertex_vs_momentum->Fill(
2132  pi0start_j.Z(), fRecoTracklistEPEM[j]->GetP());
2133  fhPi0_startvertex_vs_momentum->Fill(
2134  pi0start_k.Z(), fRecoTracklistEPEM[k]->GetP());
2135  fhPi0_startvertex_vs_momentum->Fill(
2136  pi0start_l.Z(), fRecoTracklistEPEM[l]->GetP());
2137 
2138  Double_t invmass1 = 0; // true MC values
2139  invmass1 = Invmass_4particles(fRecoTracklistEPEM[i],
2140  fRecoTracklistEPEM[j],
2141  fRecoTracklistEPEM[k],
2142  fRecoTracklistEPEM[l]);
2143  Double_t invmass2 = 0; // momenta from stsMomentumVec
2144  invmass2 = Invmass_4particlesRECO(fRecoMomentum[i],
2145  fRecoMomentum[j],
2146  fRecoMomentum[k],
2147  fRecoMomentum[l]);
2148  Double_t invmass3 = 0; // momenta from refitted at primary vertex
2149  invmass3 = Invmass_4particlesRECO(fRecoRefittedMomentum[i],
2150  fRecoRefittedMomentum[j],
2151  fRecoRefittedMomentum[k],
2152  fRecoRefittedMomentum[l]);
2153 
2154  // consider only electrons from the target (only then the momenta are correctly refitted/reconstructed)
2155  if (pi0start_i.Z() > 1 || pi0start_j.Z() > 1 || pi0start_k.Z() > 1
2156  || pi0start_l.Z() > 1) {
2157  fhEPEM_InDetector_invmass_gg_mc->Fill(invmass1);
2158  fhEPEM_InDetector_invmass_gg_refitted->Fill(invmass3);
2159  fhEPEM_InDetector_invmass_all_mc->Fill(invmass1);
2160  fhEPEM_InDetector_invmass_all_refitted->Fill(invmass3);
2161  continue;
2162  }
2163 
2164  fhEPEM_pt_vs_p_all_mc->Fill(fRecoTracklistEPEM[i]->GetPt(),
2165  fRecoTracklistEPEM[i]->GetP());
2166  fhEPEM_pt_vs_p_all_mc->Fill(fRecoTracklistEPEM[j]->GetPt(),
2167  fRecoTracklistEPEM[j]->GetP());
2168  fhEPEM_pt_vs_p_all_mc->Fill(fRecoTracklistEPEM[k]->GetPt(),
2169  fRecoTracklistEPEM[k]->GetP());
2170  fhEPEM_pt_vs_p_all_mc->Fill(fRecoTracklistEPEM[l]->GetPt(),
2171  fRecoTracklistEPEM[l]->GetP());
2172  fhEPEM_pt_vs_p_all_refitted->Fill(fRecoRefittedMomentum[i].Perp(),
2173  fRecoRefittedMomentum[i].Mag());
2174  fhEPEM_pt_vs_p_all_refitted->Fill(fRecoRefittedMomentum[j].Perp(),
2175  fRecoRefittedMomentum[j].Mag());
2176  fhEPEM_pt_vs_p_all_refitted->Fill(fRecoRefittedMomentum[k].Perp(),
2177  fRecoRefittedMomentum[k].Mag());
2178  fhEPEM_pt_vs_p_all_refitted->Fill(fRecoRefittedMomentum[l].Perp(),
2179  fRecoRefittedMomentum[l].Mag());
2180 
2181  cout << "########################################################"
2182  "##############"
2183  << endl;
2184  cout << fRecoMomentum[i].X() << "\t"
2185  << fRecoRefittedMomentum[i].X() << endl;
2186  cout << fRecoMomentum[i].Y() << "\t"
2187  << fRecoRefittedMomentum[i].Y() << endl;
2188  cout << fRecoMomentum[i].Z() << "\t"
2189  << fRecoRefittedMomentum[i].Z() << endl;
2190 
2191  if (fRecoTracklistEPEM_ids.size() != fRecoTracklistEPEM.size()) {
2192  cout << "size mismatch! " << fRecoTracklistEPEM_ids.size()
2193  << "/" << fRecoTracklistEPEM.size() << endl;
2194  }
2195 
2196  cout << "########################################################"
2197  "##############"
2198  << endl;
2199  cout << "index: " << i << "\t" << j << "\t" << k << "\t" << l
2200  << endl;
2201  cout << "mc id: " << fRecoTracklistEPEM_ids[i] << "\t"
2202  << fRecoTracklistEPEM_ids[j] << "\t"
2203  << fRecoTracklistEPEM_ids[k] << "\t"
2204  << fRecoTracklistEPEM_ids[l] << endl;
2205  cout << "motherid: " << motherId1 << "\t" << motherId2 << "\t"
2206  << motherId3 << "\t" << motherId4 << endl;
2207  cout << "motherpdg: " << mcMotherPdg1 << "\t" << mcMotherPdg2
2208  << "\t" << mcMotherPdg3 << "\t" << mcMotherPdg4 << endl;
2209  cout << "grandmotherid: " << grandmotherId1 << "\t"
2210  << grandmotherId2 << "\t" << grandmotherId3 << "\t"
2211  << grandmotherId4 << endl;
2212  cout << "pdg: " << fRecoTracklistEPEM[i]->GetPdgCode() << "\t"
2213  << fRecoTracklistEPEM[j]->GetPdgCode() << "\t"
2214  << fRecoTracklistEPEM[k]->GetPdgCode() << "\t"
2215  << fRecoTracklistEPEM[l]->GetPdgCode() << endl;
2216  cout << "invmass reco: " << invmass2
2217  << "\t invmass mc: " << invmass1 << endl;
2218 
2219  fhInvMass_EPEM_mc->Fill(invmass1);
2220  fhInvMass_EPEM_stsMomVec->Fill(invmass2);
2221  fhInvMass_EPEM_refitted->Fill(invmass3);
2222  fhInvMass_EPEM_error_stsMomVec->Fill(
2223  1.0 * TMath::Abs(invmass1 - invmass2) / invmass1);
2224  fhInvMass_EPEM_error_refitted->Fill(
2225  1.0 * TMath::Abs(invmass1 - invmass3) / invmass1);
2226 
2227  fhUsedMomenta_stsMomVec->Fill(fRecoMomentum[i].Mag());
2228  fhUsedMomenta_stsMomVec->Fill(fRecoMomentum[j].Mag());
2229  fhUsedMomenta_stsMomVec->Fill(fRecoMomentum[k].Mag());
2230  fhUsedMomenta_stsMomVec->Fill(fRecoMomentum[l].Mag());
2231 
2232  fhUsedMomenta_mc->Fill(fRecoTracklistEPEM[i]->GetP());
2233  fhUsedMomenta_mc->Fill(fRecoTracklistEPEM[j]->GetP());
2234  fhUsedMomenta_mc->Fill(fRecoTracklistEPEM[k]->GetP());
2235  fhUsedMomenta_mc->Fill(fRecoTracklistEPEM[l]->GetP());
2236 
2237  fhUsedMomenta_error_stsMomVec->Fill(
2238  TMath::Abs(fRecoTracklistEPEM[i]->GetP()
2239  - fRecoMomentum[i].Mag())
2240  / fRecoTracklistEPEM[i]->GetP());
2241  fhUsedMomenta_error_stsMomVec->Fill(
2242  TMath::Abs(fRecoTracklistEPEM[j]->GetP()
2243  - fRecoMomentum[j].Mag())
2244  / fRecoTracklistEPEM[j]->GetP());
2245  fhUsedMomenta_error_stsMomVec->Fill(
2246  TMath::Abs(fRecoTracklistEPEM[k]->GetP()
2247  - fRecoMomentum[k].Mag())
2248  / fRecoTracklistEPEM[k]->GetP());
2249  fhUsedMomenta_error_stsMomVec->Fill(
2250  TMath::Abs(fRecoTracklistEPEM[l]->GetP()
2251  - fRecoMomentum[l].Mag())
2252  / fRecoTracklistEPEM[l]->GetP());
2253 
2254  // fhUsedMomenta_error_refitted->Fill(TMath::Abs(fRecoRefittedMomentum[i].Mag() - fRecoMomentum[i].Mag())/fRecoRefittedMomentum[i].Mag());
2255  // fhUsedMomenta_error_refitted->Fill(TMath::Abs(fRecoRefittedMomentum[j].Mag() - fRecoMomentum[j].Mag())/fRecoRefittedMomentum[j].Mag());
2256  // fhUsedMomenta_error_refitted->Fill(TMath::Abs(fRecoRefittedMomentum[k].Mag() - fRecoMomentum[k].Mag())/fRecoRefittedMomentum[k].Mag());
2257  // fhUsedMomenta_error_refitted->Fill(TMath::Abs(fRecoRefittedMomentum[l].Mag() - fRecoMomentum[l].Mag())/fRecoRefittedMomentum[l].Mag());
2258 
2259  fhUsedMomenta_error_refitted->Fill(
2260  TMath::Abs(fRecoTracklistEPEM[i]->GetP()
2261  - fRecoRefittedMomentum[i].Mag())
2262  / fRecoTracklistEPEM[i]->GetP());
2263  fhUsedMomenta_error_refitted->Fill(
2264  TMath::Abs(fRecoTracklistEPEM[j]->GetP()
2265  - fRecoRefittedMomentum[j].Mag())
2266  / fRecoTracklistEPEM[j]->GetP());
2267  fhUsedMomenta_error_refitted->Fill(
2268  TMath::Abs(fRecoTracklistEPEM[k]->GetP()
2269  - fRecoRefittedMomentum[k].Mag())
2270  / fRecoTracklistEPEM[k]->GetP());
2271  fhUsedMomenta_error_refitted->Fill(
2272  TMath::Abs(fRecoTracklistEPEM[l]->GetP()
2273  - fRecoRefittedMomentum[l].Mag())
2274  / fRecoTracklistEPEM[l]->GetP());
2275 
2276 
2277  TVector3 momentumtest5a;
2278  fRecoTracklistEPEM[i]->GetMomentum(momentumtest5a);
2279  fhUsedMomenta_errorX_stsMomVec->Fill(
2280  TMath::Abs(momentumtest5a.X() - fRecoMomentum[i].X())
2281  / momentumtest5a.X());
2282  fhUsedMomenta_errorY_stsMomVec->Fill(
2283  TMath::Abs(momentumtest5a.Y() - fRecoMomentum[i].Y())
2284  / momentumtest5a.Y());
2285  fhUsedMomenta_errorZ_stsMomVec->Fill(
2286  TMath::Abs(momentumtest5a.Z() - fRecoMomentum[i].Z())
2287  / momentumtest5a.Z());
2288  fhUsedMomenta_vsX_stsMomVec->Fill(momentumtest5a.X(),
2289  fRecoMomentum[i].X());
2290  fhUsedMomenta_vsY_stsMomVec->Fill(momentumtest5a.Y(),
2291  fRecoMomentum[i].Y());
2292  fhUsedMomenta_vsZ_stsMomVec->Fill(momentumtest5a.Z(),
2293  fRecoMomentum[i].Z());
2294  fhUsedMomenta_errorX_refitted->Fill(
2295  TMath::Abs(momentumtest5a.X() - fRecoRefittedMomentum[i].X())
2296  / momentumtest5a.X());
2297  fhUsedMomenta_errorY_refitted->Fill(
2298  TMath::Abs(momentumtest5a.Y() - fRecoRefittedMomentum[i].Y())
2299  / momentumtest5a.Y());
2300  fhUsedMomenta_errorZ_refitted->Fill(
2301  TMath::Abs(momentumtest5a.Z() - fRecoRefittedMomentum[i].Z())
2302  / momentumtest5a.Z());
2303  fhUsedMomenta_vsX_refitted->Fill(momentumtest5a.X(),
2304  fRecoRefittedMomentum[i].X());
2305  fhUsedMomenta_vsY_refitted->Fill(momentumtest5a.Y(),
2306  fRecoRefittedMomentum[i].Y());
2307  fhUsedMomenta_vsZ_refitted->Fill(momentumtest5a.Z(),
2308  fRecoRefittedMomentum[i].Z());
2309 
2310  TVector3 momentumtest5b;
2311  fRecoTracklistEPEM[j]->GetMomentum(momentumtest5b);
2312  fhUsedMomenta_errorX_stsMomVec->Fill(
2313  TMath::Abs(momentumtest5b.X() - fRecoMomentum[j].X())
2314  / momentumtest5b.X());
2315  fhUsedMomenta_errorY_stsMomVec->Fill(
2316  TMath::Abs(momentumtest5b.Y() - fRecoMomentum[j].Y())
2317  / momentumtest5b.Y());
2318  fhUsedMomenta_errorZ_stsMomVec->Fill(
2319  TMath::Abs(momentumtest5b.Z() - fRecoMomentum[j].Z())
2320  / momentumtest5b.Z());
2321  fhUsedMomenta_vsX_stsMomVec->Fill(momentumtest5b.X(),
2322  fRecoMomentum[j].X());
2323  fhUsedMomenta_vsY_stsMomVec->Fill(momentumtest5b.Y(),
2324  fRecoMomentum[j].Y());
2325  fhUsedMomenta_vsZ_stsMomVec->Fill(momentumtest5b.Z(),
2326  fRecoMomentum[j].Z());
2327  fhUsedMomenta_errorX_refitted->Fill(
2328  TMath::Abs(momentumtest5b.X() - fRecoRefittedMomentum[j].X())
2329  / momentumtest5b.X());
2330  fhUsedMomenta_errorY_refitted->Fill(
2331  TMath::Abs(momentumtest5b.Y() - fRecoRefittedMomentum[j].Y())
2332  / momentumtest5b.Y());
2333  fhUsedMomenta_errorZ_refitted->Fill(
2334  TMath::Abs(momentumtest5b.Z() - fRecoRefittedMomentum[j].Z())
2335  / momentumtest5b.Z());
2336  fhUsedMomenta_vsX_refitted->Fill(momentumtest5b.X(),
2337  fRecoRefittedMomentum[j].X());
2338  fhUsedMomenta_vsY_refitted->Fill(momentumtest5b.Y(),
2339  fRecoRefittedMomentum[j].Y());
2340  fhUsedMomenta_vsZ_refitted->Fill(momentumtest5b.Z(),
2341  fRecoRefittedMomentum[j].Z());
2342 
2343  TVector3 momentumtest5c;
2344  fRecoTracklistEPEM[k]->GetMomentum(momentumtest5c);
2345  fhUsedMomenta_errorX_stsMomVec->Fill(
2346  TMath::Abs(momentumtest5c.X() - fRecoMomentum[k].X())
2347  / momentumtest5c.X());
2348  fhUsedMomenta_errorY_stsMomVec->Fill(
2349  TMath::Abs(momentumtest5c.Y() - fRecoMomentum[k].Y())
2350  / momentumtest5c.Y());
2351  fhUsedMomenta_errorZ_stsMomVec->Fill(
2352  TMath::Abs(momentumtest5c.Z() - fRecoMomentum[k].Z())
2353  / momentumtest5c.Z());
2354  fhUsedMomenta_vsX_stsMomVec->Fill(momentumtest5c.X(),
2355  fRecoMomentum[k].X());
2356  fhUsedMomenta_vsY_stsMomVec->Fill(momentumtest5c.Y(),
2357  fRecoMomentum[k].Y());
2358  fhUsedMomenta_vsZ_stsMomVec->Fill(momentumtest5c.Z(),
2359  fRecoMomentum[k].Z());
2360  fhUsedMomenta_errorX_refitted->Fill(
2361  TMath::Abs(momentumtest5c.X() - fRecoRefittedMomentum[k].X())
2362  / momentumtest5c.X());
2363  fhUsedMomenta_errorY_refitted->Fill(
2364  TMath::Abs(momentumtest5c.Y() - fRecoRefittedMomentum[k].Y())
2365  / momentumtest5c.Y());
2366  fhUsedMomenta_errorZ_refitted->Fill(
2367  TMath::Abs(momentumtest5c.Z() - fRecoRefittedMomentum[k].Z())
2368  / momentumtest5c.Z());
2369  fhUsedMomenta_vsX_refitted->Fill(momentumtest5c.X(),
2370  fRecoRefittedMomentum[k].X());
2371  fhUsedMomenta_vsY_refitted->Fill(momentumtest5c.Y(),
2372  fRecoRefittedMomentum[k].Y());
2373  fhUsedMomenta_vsZ_refitted->Fill(momentumtest5c.Z(),
2374  fRecoRefittedMomentum[k].Z());
2375 
2376  TVector3 momentumtest5d;
2377  fRecoTracklistEPEM[l]->GetMomentum(momentumtest5d);
2378  fhUsedMomenta_errorX_stsMomVec->Fill(
2379  TMath::Abs(momentumtest5d.X() - fRecoMomentum[l].X())
2380  / momentumtest5d.X());
2381  fhUsedMomenta_errorY_stsMomVec->Fill(
2382  TMath::Abs(momentumtest5d.Y() - fRecoMomentum[l].Y())
2383  / momentumtest5d.Y());
2384  fhUsedMomenta_errorZ_stsMomVec->Fill(
2385  TMath::Abs(momentumtest5d.Z() - fRecoMomentum[l].Z())
2386  / momentumtest5d.Z());
2387  fhUsedMomenta_vsX_stsMomVec->Fill(momentumtest5d.X(),
2388  fRecoMomentum[l].X());
2389  fhUsedMomenta_vsY_stsMomVec->Fill(momentumtest5d.Y(),
2390  fRecoMomentum[l].Y());
2391  fhUsedMomenta_vsZ_stsMomVec->Fill(momentumtest5d.Z(),
2392  fRecoMomentum[l].Z());
2393  fhUsedMomenta_errorX_refitted->Fill(
2394  TMath::Abs(momentumtest5d.X() - fRecoRefittedMomentum[l].X())
2395  / momentumtest5d.X());
2396  fhUsedMomenta_errorY_refitted->Fill(
2397  TMath::Abs(momentumtest5d.Y() - fRecoRefittedMomentum[l].Y())
2398  / momentumtest5d.Y());
2399  fhUsedMomenta_errorZ_refitted->Fill(
2400  TMath::Abs(momentumtest5d.Z() - fRecoRefittedMomentum[l].Z())
2401  / momentumtest5d.Z());
2402  fhUsedMomenta_vsX_refitted->Fill(momentumtest5d.X(),
2403  fRecoRefittedMomentum[l].X());
2404  fhUsedMomenta_vsY_refitted->Fill(momentumtest5d.Y(),
2405  fRecoRefittedMomentum[l].Y());
2406  fhUsedMomenta_vsZ_refitted->Fill(momentumtest5d.Z(),
2407  fRecoRefittedMomentum[l].Z());
2408 
2409 
2410  Double_t opening_angle1_mc = 0;
2411  Double_t opening_angle2_mc = 0;
2412  Double_t opening_angle1_refitted = 0;
2413  Double_t opening_angle2_refitted = 0;
2414 
2415  if (motherId1 == motherId2) {
2416  opening_angle1_mc = CalculateOpeningAngleMC(
2417  fRecoTracklistEPEM[i], fRecoTracklistEPEM[j]);
2418  opening_angle2_mc = CalculateOpeningAngleMC(
2419  fRecoTracklistEPEM[k], fRecoTracklistEPEM[l]);
2420  opening_angle1_refitted = CalculateOpeningAngleReco(
2421  fRecoRefittedMomentum[i], fRecoRefittedMomentum[j]);
2422  opening_angle2_refitted = CalculateOpeningAngleReco(
2423  fRecoRefittedMomentum[k], fRecoRefittedMomentum[l]);
2424  }
2425  if (motherId1 == motherId3) {
2426  opening_angle1_mc = CalculateOpeningAngleMC(
2427  fRecoTracklistEPEM[i], fRecoTracklistEPEM[k]);
2428  opening_angle2_mc = CalculateOpeningAngleMC(
2429  fRecoTracklistEPEM[j], fRecoTracklistEPEM[l]);
2430  opening_angle1_refitted = CalculateOpeningAngleReco(
2431  fRecoRefittedMomentum[i], fRecoRefittedMomentum[k]);
2432  opening_angle2_refitted = CalculateOpeningAngleReco(
2433  fRecoRefittedMomentum[j], fRecoRefittedMomentum[l]);
2434  }
2435  if (motherId1 == motherId4) {
2436  opening_angle1_mc = CalculateOpeningAngleMC(
2437  fRecoTracklistEPEM[i], fRecoTracklistEPEM[l]);
2438  opening_angle2_mc = CalculateOpeningAngleMC(
2439  fRecoTracklistEPEM[j], fRecoTracklistEPEM[k]);
2440  opening_angle1_refitted = CalculateOpeningAngleReco(
2441  fRecoRefittedMomentum[i], fRecoRefittedMomentum[l]);
2442  opening_angle2_refitted = CalculateOpeningAngleReco(
2443  fRecoRefittedMomentum[j], fRecoRefittedMomentum[k]);
2444  }
2445 
2446 
2447  fhInvMass_EPEM_openingAngleRef->Fill(opening_angle1_refitted);
2448  fhInvMass_EPEM_openingAngleRef->Fill(opening_angle2_refitted);
2449 
2450 
2451  fhEPEM_invmass_gg_mc->Fill(invmass1);
2452  fhEPEM_invmass_gg_refitted->Fill(invmass3);
2453  fhEPEM_invmass_all_mc->Fill(invmass1);
2454  fhEPEM_invmass_all_refitted->Fill(invmass3);
2455  fhEPEM_openingAngle_gg_mc->Fill(opening_angle1_mc);
2456  fhEPEM_openingAngle_gg_mc->Fill(opening_angle2_mc);
2457  fhEPEM_openingAngle_gg_refitted->Fill(opening_angle1_refitted);
2458  fhEPEM_openingAngle_gg_refitted->Fill(opening_angle2_refitted);
2459 
2460  Double_t openingAngleBetweenGammas =
2461  CalculateOpeningAngleBetweenGammasMC(fRecoTracklistEPEM[i],
2462  fRecoTracklistEPEM[j],
2463  fRecoTracklistEPEM[k],
2464  fRecoTracklistEPEM[l]);
2465  fhEPEM_openingAngle_betweenGammas_mc->Fill(
2466  openingAngleBetweenGammas);
2467  Double_t openingAngleBetweenGammasReco =
2468  CalculateOpeningAngleBetweenGammasReco(
2469  fRecoRefittedMomentum[i],
2470  fRecoRefittedMomentum[j],
2471  fRecoRefittedMomentum[k],
2472  fRecoRefittedMomentum[l]);
2473  fhEPEM_openingAngle_betweenGammas_reco->Fill(
2474  openingAngleBetweenGammasReco);
2475 
2476 
2477  CbmLmvmKinematicParams params1 =
2478  CalculateKinematicParams_4particles(fRecoRefittedMomentum[i],
2479  fRecoRefittedMomentum[j],
2480  fRecoRefittedMomentum[k],
2481  fRecoRefittedMomentum[l]);
2482  fhPi0_pt_vs_rap_gg->Fill(params1.fPt, params1.fRapidity);
2483  fhPi0_pt_vs_rap_all->Fill(params1.fPt, params1.fRapidity);
2484  fhPi0_pt_gg->Fill(params1.fPt);
2485  fhPi0_pt_all->Fill(params1.fPt);
2486 
2487 
2488  fhEPEM_openingAngle_vs_pt_gg_reco->Fill(params1.fPt,
2489  opening_angle1_refitted);
2490  fhEPEM_openingAngle_vs_pt_gg_reco->Fill(params1.fPt,
2491  opening_angle2_refitted);
2492 
2493  fhEPEM_rap_vs_chi->Fill(params1.fRapidity,
2494  fRecoTracklistEPEM_chi[i]);
2495  fhEPEM_rap_vs_chi->Fill(params1.fRapidity,
2496  fRecoTracklistEPEM_chi[j]);
2497  fhEPEM_rap_vs_chi->Fill(params1.fRapidity,
2498  fRecoTracklistEPEM_chi[k]);
2499  fhEPEM_rap_vs_chi->Fill(params1.fRapidity,
2500  fRecoTracklistEPEM_chi[l]);
2501 
2502  fhEPEM_rap_vs_invmass->Fill(params1.fRapidity, invmass3);
2503 
2504 
2505  // ####################################
2506  // STUDIES: efficiency of cuts
2507  // ####################################
2508  //fRichElIdAnn = new CbmRichElectronIdAnn();
2509  //fRichElIdAnn->Init();
2510  //Double_t ann = fRichElIdAnn->DoSelect(ring, momentum);
2511  //if (ann > fCuts.fRichAnnCut) return true; // cut = 0.0
2512 
2513 
2514  CutEfficiencyStudies(
2515  i, j, k, l, motherId1, motherId2, motherId3, motherId4);
2516 
2517  //Bool_t IsRichElectron1 = CbmLitGlobalElectronId::GetInstance().IsRichElectron(fRecoTracklistEPEM_gtid[i], fRecoRefittedMomentum[i].Mag());
2518  //Bool_t IsRichElectron2 = CbmLitGlobalElectronId::GetInstance().IsRichElectron(fRecoTracklistEPEM_gtid[j], fRecoRefittedMomentum[j].Mag());
2519  //Bool_t IsRichElectron3 = CbmLitGlobalElectronId::GetInstance().IsRichElectron(fRecoTracklistEPEM_gtid[k], fRecoRefittedMomentum[k].Mag());
2520  //Bool_t IsRichElectron4 = CbmLitGlobalElectronId::GetInstance().IsRichElectron(fRecoTracklistEPEM_gtid[l], fRecoRefittedMomentum[l].Mag());
2521 
2522  Bool_t IsRichElectron1ann = IsRichElectronANN(
2523  fRecoTracklistEPEM_gtid[i], fRecoRefittedMomentum[i].Mag());
2524  Bool_t IsRichElectron2ann = IsRichElectronANN(
2525  fRecoTracklistEPEM_gtid[j], fRecoRefittedMomentum[j].Mag());
2526  Bool_t IsRichElectron3ann = IsRichElectronANN(
2527  fRecoTracklistEPEM_gtid[k], fRecoRefittedMomentum[k].Mag());
2528  Bool_t IsRichElectron4ann = IsRichElectronANN(
2529  fRecoTracklistEPEM_gtid[l], fRecoRefittedMomentum[l].Mag());
2530 
2531  Bool_t IsRichElectron1normal = IsRichElectronNormal(
2532  fRecoTracklistEPEM_gtid[i], fRecoRefittedMomentum[i].Mag());
2533  Bool_t IsRichElectron2normal = IsRichElectronNormal(
2534  fRecoTracklistEPEM_gtid[j], fRecoRefittedMomentum[j].Mag());
2535  Bool_t IsRichElectron3normal = IsRichElectronNormal(
2536  fRecoTracklistEPEM_gtid[k], fRecoRefittedMomentum[k].Mag());
2537  Bool_t IsRichElectron4normal = IsRichElectronNormal(
2538  fRecoTracklistEPEM_gtid[l], fRecoRefittedMomentum[l].Mag());
2539 
2540  Bool_t IsRichElectron1MC =
2541  (TMath::Abs(fRecoTracklistEPEM[i]->GetPdgCode()) == 11);
2542  Bool_t IsRichElectron2MC =
2543  (TMath::Abs(fRecoTracklistEPEM[j]->GetPdgCode()) == 11);
2544  Bool_t IsRichElectron3MC =
2545  (TMath::Abs(fRecoTracklistEPEM[k]->GetPdgCode()) == 11);
2546  Bool_t IsRichElectron4MC =
2547  (TMath::Abs(fRecoTracklistEPEM[l]->GetPdgCode()) == 11);
2548 
2549  CbmLmvmKinematicParams paramsCut1;
2550  CbmLmvmKinematicParams paramsCut2;
2551 
2552  if (motherId1 == motherId2) {
2553  paramsCut1 = CalculateKinematicParamsReco(
2554  fRecoRefittedMomentum[i], fRecoRefittedMomentum[j]);
2555  paramsCut2 = CalculateKinematicParamsReco(
2556  fRecoRefittedMomentum[k], fRecoRefittedMomentum[l]);
2557  }
2558  if (motherId1 == motherId3) {
2559  paramsCut1 = CalculateKinematicParamsReco(
2560  fRecoRefittedMomentum[i], fRecoRefittedMomentum[k]);
2561  paramsCut2 = CalculateKinematicParamsReco(
2562  fRecoRefittedMomentum[j], fRecoRefittedMomentum[l]);
2563  }
2564  if (motherId1 == motherId4) {
2565  paramsCut1 = CalculateKinematicParamsReco(
2566  fRecoRefittedMomentum[i], fRecoRefittedMomentum[l]);
2567  paramsCut2 = CalculateKinematicParamsReco(
2568  fRecoRefittedMomentum[j], fRecoRefittedMomentum[k]);
2569  }
2570 
2571 
2572  Double_t Value_invariantMassCut = 0.03;
2573  Double_t Value_openingAngleCut1 = 1.8 - 0.6 * paramsCut1.fPt;
2574  Double_t Value_openingAngleCut2 = 1.8 - 0.6 * paramsCut2.fPt;
2575 
2576  Bool_t OpeningAngleCut1 =
2577  (opening_angle1_refitted < Value_openingAngleCut1);
2578  Bool_t OpeningAngleCut2 =
2579  (opening_angle2_refitted < Value_openingAngleCut2);
2580  Bool_t InvariantMassCut1 =
2581  (paramsCut1.fMinv < Value_invariantMassCut);
2582  Bool_t InvariantMassCut2 =
2583  (paramsCut2.fMinv < Value_invariantMassCut);
2584 
2585 
2586  fhEPEM_efficiencyCuts2->Fill(0); // no further cuts applied
2587  // first ANN usage for electron identification
2588  if (
2589  IsRichElectron1ann && IsRichElectron2ann && IsRichElectron3ann
2590  && IsRichElectron4ann) { // all 4 electrons correctly identified with the RICH via ANN
2591  fhEPEM_efficiencyCuts2->Fill(1);
2592  if (
2593  OpeningAngleCut1
2594  && OpeningAngleCut2) { // opening angle of e+e- pairs below x
2595  fhEPEM_efficiencyCuts2->Fill(2);
2596  if (InvariantMassCut1 && InvariantMassCut2) {
2597  fhEPEM_efficiencyCuts2->Fill(3);
2598  }
2599  }
2600  }
2601  // then standard method for electron identification
2602  if (
2603  IsRichElectron1normal && IsRichElectron2normal
2604  && IsRichElectron3normal
2605  && IsRichElectron4normal) { // all 4 electrons correctly identified with the RICH via "normal way"
2606  fhEPEM_efficiencyCuts2->Fill(4);
2607  if (
2608  OpeningAngleCut1
2609  && OpeningAngleCut2) { // opening angle of e+e- pairs below x
2610  fhEPEM_efficiencyCuts2->Fill(5);
2611  if (InvariantMassCut1 && InvariantMassCut2) {
2612  fhEPEM_efficiencyCuts2->Fill(6);
2613  }
2614  }
2615  }
2616  // MC-true data for electron identification
2617  if (
2618  IsRichElectron1MC && IsRichElectron2MC && IsRichElectron3MC
2619  && IsRichElectron4MC) { // all 4 electrons correctly identified with the RICH via MC-true data
2620  fhEPEM_efficiencyCuts2->Fill(7);
2621  if (
2622  OpeningAngleCut1
2623  && OpeningAngleCut2) { // opening angle of e+e- pairs below x
2624  fhEPEM_efficiencyCuts2->Fill(8);
2625  if (InvariantMassCut1 && InvariantMassCut2) {
2626  fhEPEM_efficiencyCuts2->Fill(9);
2627  }
2628  }
2629  }
2630 
2631 
2632  cout << "reco/mc: " << fRecoMomentum[i].Mag() << " / "
2633  << fRecoTracklistEPEM[i]->GetP() << " ### "
2634  << fRecoMomentum[j].Mag() << " / "
2635  << fRecoTracklistEPEM[j]->GetP() << " ### "
2636  << fRecoMomentum[k].Mag() << " / "
2637  << fRecoTracklistEPEM[k]->GetP() << " ### "
2638  << fRecoMomentum[l].Mag() << " / "
2639  << fRecoTracklistEPEM[l]->GetP() << endl;
2640 
2641  fill++;
2642  }
2643  }
2644  }
2645  }
2646  }
2647  }
2648  cout << "CbmAnaConversionReco: InvariantMassTest_4epem - Filled events: "
2649  << fill << endl;
2650  cout << "CbmAnaConversionReco: InvariantMassTest_4epem - End!" << endl;
2651 
2652 
2653  timer.Stop();
2654  fTime += timer.RealTime();
2655 }
2656 
2657 
2658 void CbmAnaConversionReco::CutEfficiencyStudies(int e1,
2659  int e2,
2660  int e3,
2661  int e4,
2662  int motherE1,
2663  int motherE2,
2664  int motherE3,
2665  int motherE4,
2666  int IsEta) {
2667  // ####################################
2668  // STUDIES: efficiency of cuts
2669  // ####################################
2670  //fRichElIdAnn = new CbmRichElectronIdAnn();
2671  //fRichElIdAnn->Init();
2672  //Double_t ann = fRichElIdAnn->DoSelect(ring, momentum);
2673  //if (ann > fCuts.fRichAnnCut) return true; // cut = 0.0
2674 
2675  Double_t opening_angle1_refitted = 0;
2676  Double_t opening_angle2_refitted = 0;
2677 
2678  if (motherE1 == motherE2) {
2679  opening_angle1_refitted = CalculateOpeningAngleReco(
2680  fRecoRefittedMomentum[e1], fRecoRefittedMomentum[e2]);
2681  opening_angle2_refitted = CalculateOpeningAngleReco(
2682  fRecoRefittedMomentum[e3], fRecoRefittedMomentum[e4]);
2683  }
2684  if (motherE1 == motherE3) {
2685  opening_angle1_refitted = CalculateOpeningAngleReco(
2686  fRecoRefittedMomentum[e1], fRecoRefittedMomentum[e3]);
2687  opening_angle2_refitted = CalculateOpeningAngleReco(
2688  fRecoRefittedMomentum[e2], fRecoRefittedMomentum[e4]);
2689  }
2690  if (motherE1 == motherE4) {
2691  opening_angle1_refitted = CalculateOpeningAngleReco(
2692  fRecoRefittedMomentum[e1], fRecoRefittedMomentum[e4]);
2693  opening_angle2_refitted = CalculateOpeningAngleReco(
2694  fRecoRefittedMomentum[e2], fRecoRefittedMomentum[e3]);
2695  }
2696 
2697 
2698  Bool_t IsWithinChiCut1 =
2699  (fRecoTracklistEPEM_chi[e1] < CbmAnaConversionCutSettings::CalcChiCut(
2700  fRecoRefittedMomentum[e1].Perp()));
2701  Bool_t IsWithinChiCut2 =
2702  (fRecoTracklistEPEM_chi[e2] < CbmAnaConversionCutSettings::CalcChiCut(
2703  fRecoRefittedMomentum[e2].Perp()));
2704  Bool_t IsWithinChiCut3 =
2705  (fRecoTracklistEPEM_chi[e3] < CbmAnaConversionCutSettings::CalcChiCut(
2706  fRecoRefittedMomentum[e3].Perp()));
2707  Bool_t IsWithinChiCut4 =
2708  (fRecoTracklistEPEM_chi[e4] < CbmAnaConversionCutSettings::CalcChiCut(
2709  fRecoRefittedMomentum[e4].Perp()));
2710 
2711  Bool_t AllWithinChiCut =
2712  (IsWithinChiCut1 && IsWithinChiCut2 && IsWithinChiCut3 && IsWithinChiCut4);
2713 
2714 
2715  //Bool_t IsRichElectron1 = CbmLitGlobalElectronId::GetInstance().IsRichElectron(fRecoTracklistEPEM_gtid[i], fRecoRefittedMomentum[i].Mag());
2716  //Bool_t IsRichElectron2 = CbmLitGlobalElectronId::GetInstance().IsRichElectron(fRecoTracklistEPEM_gtid[j], fRecoRefittedMomentum[j].Mag());
2717  //Bool_t IsRichElectron3 = CbmLitGlobalElectronId::GetInstance().IsRichElectron(fRecoTracklistEPEM_gtid[k], fRecoRefittedMomentum[k].Mag());
2718  //Bool_t IsRichElectron4 = CbmLitGlobalElectronId::GetInstance().IsRichElectron(fRecoTracklistEPEM_gtid[l], fRecoRefittedMomentum[l].Mag());
2719 
2720  Bool_t IsRichElectron1ann = IsRichElectronANN(
2721  fRecoTracklistEPEM_gtid[e1], fRecoRefittedMomentum[e1].Mag());
2722  Bool_t IsRichElectron2ann = IsRichElectronANN(
2723  fRecoTracklistEPEM_gtid[e2], fRecoRefittedMomentum[e2].Mag());
2724  Bool_t IsRichElectron3ann = IsRichElectronANN(
2725  fRecoTracklistEPEM_gtid[e3], fRecoRefittedMomentum[e3].Mag());
2726  Bool_t IsRichElectron4ann = IsRichElectronANN(
2727  fRecoTracklistEPEM_gtid[e4], fRecoRefittedMomentum[e4].Mag());
2728 
2729  Bool_t IsRichElectron1normal = IsRichElectronNormal(
2730  fRecoTracklistEPEM_gtid[e1], fRecoRefittedMomentum[e1].Mag());
2731  Bool_t IsRichElectron2normal = IsRichElectronNormal(
2732  fRecoTracklistEPEM_gtid[e2], fRecoRefittedMomentum[e2].Mag());
2733  Bool_t IsRichElectron3normal = IsRichElectronNormal(
2734  fRecoTracklistEPEM_gtid[e3], fRecoRefittedMomentum[e3].Mag());
2735  Bool_t IsRichElectron4normal = IsRichElectronNormal(
2736  fRecoTracklistEPEM_gtid[e4], fRecoRefittedMomentum[e4].Mag());
2737 
2738  Bool_t IsRichElectron1MC =
2739  (TMath::Abs(fRecoTracklistEPEM[e1]->GetPdgCode()) == 11);
2740  Bool_t IsRichElectron2MC =
2741  (TMath::Abs(fRecoTracklistEPEM[e2]->GetPdgCode()) == 11);
2742  Bool_t IsRichElectron3MC =
2743  (TMath::Abs(fRecoTracklistEPEM[e3]->GetPdgCode()) == 11);
2744  Bool_t IsRichElectron4MC =
2745  (TMath::Abs(fRecoTracklistEPEM[e4]->GetPdgCode()) == 11);
2746 
2747  Double_t ANNvalueE1 = ElectronANNvalue(fRecoTracklistEPEM_gtid[e1],
2748  fRecoRefittedMomentum[e1].Mag());
2749  Double_t ANNvalueE2 = ElectronANNvalue(fRecoTracklistEPEM_gtid[e2],
2750  fRecoRefittedMomentum[e2].Mag());
2751  Double_t ANNvalueE3 = ElectronANNvalue(fRecoTracklistEPEM_gtid[e3],
2752  fRecoRefittedMomentum[e3].Mag());
2753  Double_t ANNvalueE4 = ElectronANNvalue(fRecoTracklistEPEM_gtid[e4],
2754  fRecoRefittedMomentum[e4].Mag());
2755 
2756  CbmLmvmKinematicParams paramsCut1;
2757  CbmLmvmKinematicParams paramsCut2;
2758 
2759  if (motherE1 == motherE2) {
2760  paramsCut1 = CalculateKinematicParamsReco(fRecoRefittedMomentum[e1],
2761  fRecoRefittedMomentum[e2]);
2762  paramsCut2 = CalculateKinematicParamsReco(fRecoRefittedMomentum[e3],
2763  fRecoRefittedMomentum[e4]);
2764  }
2765  if (motherE1 == motherE3) {
2766  paramsCut1 = CalculateKinematicParamsReco(fRecoRefittedMomentum[e1],
2767  fRecoRefittedMomentum[e3]);
2768  paramsCut2 = CalculateKinematicParamsReco(fRecoRefittedMomentum[e2],
2769  fRecoRefittedMomentum[e4]);
2770  }
2771  if (motherE1 == motherE4) {
2772  paramsCut1 = CalculateKinematicParamsReco(fRecoRefittedMomentum[e1],
2773  fRecoRefittedMomentum[e4]);
2774  paramsCut2 = CalculateKinematicParamsReco(fRecoRefittedMomentum[e2],
2775  fRecoRefittedMomentum[e3]);
2776  }
2777 
2778 
2779  Double_t Value_invariantMassCut = 0.03;
2780  Double_t Value_openingAngleCut1 = 1.8 - 0.6 * paramsCut1.fPt;
2781  Double_t Value_openingAngleCut2 = 1.8 - 0.6 * paramsCut2.fPt;
2782 
2783  Bool_t OpeningAngleCut1 = (opening_angle1_refitted < Value_openingAngleCut1);
2784  Bool_t OpeningAngleCut2 = (opening_angle2_refitted < Value_openingAngleCut2);
2785  Bool_t InvariantMassCut1 = (paramsCut1.fMinv < Value_invariantMassCut);
2786  Bool_t InvariantMassCut2 = (paramsCut2.fMinv < Value_invariantMassCut);
2787 
2788 
2789  if (IsEta == 0) {
2790  fhEPEM_pi0_nofLeptons_ann->Fill(IsRichElectron1ann + IsRichElectron2ann
2791  + IsRichElectron3ann + IsRichElectron4ann);
2792 
2793  if (IsRichElectron1ann + IsRichElectron2ann + IsRichElectron3ann
2794  + IsRichElectron4ann
2795  == 3) {
2796  CbmGlobalTrack* gTrack = NULL;
2797  if (IsRichElectron1ann == 0)
2798  gTrack =
2799  (CbmGlobalTrack*) fGlobalTracks->At(fRecoTracklistEPEM_gtid[e1]);
2800  if (IsRichElectron2ann == 0)
2801  gTrack =
2802  (CbmGlobalTrack*) fGlobalTracks->At(fRecoTracklistEPEM_gtid[e2]);
2803  if (IsRichElectron3ann == 0)
2804  gTrack =
2805  (CbmGlobalTrack*) fGlobalTracks->At(fRecoTracklistEPEM_gtid[e3]);
2806  if (IsRichElectron4ann == 0)
2807  gTrack =
2808  (CbmGlobalTrack*) fGlobalTracks->At(fRecoTracklistEPEM_gtid[e4]);
2809  int richInd = gTrack->GetRichRingIndex();
2810  CbmRichRing* richRing = (CbmRichRing*) fRichRings->At(richInd);
2811  Int_t nofringhits = richRing->GetNofHits();
2812  fhEPEM_missingLepton_nofRingHits->Fill(nofringhits);
2813  fhEPEM_missingLepton_ringMid->Fill(richRing->GetCenterX(),
2814  richRing->GetCenterY());
2815  fhEPEM_missingLepton_ringRadius->Fill(richRing->GetRadius());
2816  // CbmRichRing::GetDistance() method is no longer supported
2817  // If you wan to use cuts update code using CbmRichUtil::GetRingTrackDistance()
2818  fhEPEM_missingLepton_distance->Fill(1. /*richRing->GetDistance()*/);
2819  fhEPEM_missingLepton_selectionNN->Fill(richRing->GetSelectionNN());
2820  if (IsRichElectron1ann == 0) {
2821  fhEPEM_missingLepton_radius_vs_p->Fill(fRecoRefittedMomentum[e1].Mag(),
2822  richRing->GetRadius());
2823  fhEPEM_missingLepton_ANNvalue->Fill(ANNvalueE1);
2824  }
2825  if (IsRichElectron2ann == 0) {
2826  fhEPEM_missingLepton_radius_vs_p->Fill(fRecoRefittedMomentum[e2].Mag(),
2827  richRing->GetRadius());
2828  fhEPEM_missingLepton_ANNvalue->Fill(ANNvalueE2);
2829  }
2830  if (IsRichElectron3ann == 0) {
2831  fhEPEM_missingLepton_radius_vs_p->Fill(fRecoRefittedMomentum[e3].Mag(),
2832  richRing->GetRadius());
2833  fhEPEM_missingLepton_ANNvalue->Fill(ANNvalueE3);
2834  }
2835  if (IsRichElectron4ann == 0) {
2836  fhEPEM_missingLepton_radius_vs_p->Fill(fRecoRefittedMomentum[e4].Mag(),
2837  richRing->GetRadius());
2838  fhEPEM_missingLepton_ANNvalue->Fill(ANNvalueE4);
2839  }
2840 
2841  for (int i = 0; i < nofringhits; i++) {
2842  Int_t hitInd = richRing->GetHit(i);
2843  CbmRichHit* hit = (CbmRichHit*) fRichHits->At(hitInd);
2844  if (NULL == hit) continue;
2845  //CbmRichHitLight hl(hit->GetX(), hit->GetY());
2846  fhEPEM_missingLepton_rings->Fill(hit->GetX(), hit->GetY());
2847  }
2848  }
2849 
2850 
2851  if (IsRichElectron1ann > 0) {
2852  CbmGlobalTrack* gTrack =
2853  (CbmGlobalTrack*) fGlobalTracks->At(fRecoTracklistEPEM_gtid[e1]);
2854  int richInd = gTrack->GetRichRingIndex();
2855  CbmRichRing* richRing = (CbmRichRing*) fRichRings->At(richInd);
2856  Int_t nofringhits = richRing->GetNofHits();
2857  fhEPEM_identifiedLepton_nofRingHits->Fill(nofringhits);
2858  fhEPEM_identifiedLepton_ringMid->Fill(richRing->GetCenterX(),
2859  richRing->GetCenterY());
2860  fhEPEM_identifiedLepton_ringRadius->Fill(richRing->GetRadius());
2861  // CbmRichRing::GetDistance() method is no longer supported
2862  // If you wan to use cuts update code using CbmRichUtil::GetRingTrackDistance()
2863  fhEPEM_identifiedLepton_distance->Fill(1. /*richRing->GetDistance()*/);
2864  fhEPEM_identifiedLepton_selectionNN->Fill(richRing->GetSelectionNN());
2865  fhEPEM_identifiedLepton_radius_vs_p->Fill(fRecoRefittedMomentum[e1].Mag(),
2866  richRing->GetRadius());
2867  fhEPEM_identifiedLepton_ANNvalue->Fill(ANNvalueE1);
2868  }
2869  if (IsRichElectron2ann > 0) {
2870  CbmGlobalTrack* gTrack =
2871  (CbmGlobalTrack*) fGlobalTracks->At(fRecoTracklistEPEM_gtid[e2]);
2872  int richInd = gTrack->GetRichRingIndex();
2873  CbmRichRing* richRing = (CbmRichRing*) fRichRings->At(richInd);
2874  Int_t nofringhits = richRing->GetNofHits();
2875  fhEPEM_identifiedLepton_nofRingHits->Fill(nofringhits);
2876  fhEPEM_identifiedLepton_ringMid->Fill(richRing->GetCenterX(),
2877  richRing->GetCenterY());
2878  fhEPEM_identifiedLepton_ringRadius->Fill(richRing->GetRadius());
2879  // CbmRichRing::GetDistance() method is no longer supported
2880  // If you wan to use cuts update code using CbmRichUtil::GetRingTrackDistance()
2881  fhEPEM_identifiedLepton_distance->Fill(1. /*richRing->GetDistance()*/);
2882  fhEPEM_identifiedLepton_selectionNN->Fill(richRing->GetSelectionNN());
2883  fhEPEM_identifiedLepton_radius_vs_p->Fill(fRecoRefittedMomentum[e2].Mag(),
2884  richRing->GetRadius());
2885  fhEPEM_identifiedLepton_ANNvalue->Fill(ANNvalueE2);
2886  }
2887  if (IsRichElectron3ann > 0) {
2888  CbmGlobalTrack* gTrack =
2889  (CbmGlobalTrack*) fGlobalTracks->At(fRecoTracklistEPEM_gtid[e3]);
2890  int richInd = gTrack->GetRichRingIndex();
2891  CbmRichRing* richRing = (CbmRichRing*) fRichRings->At(richInd);
2892  Int_t nofringhits = richRing->GetNofHits();
2893  fhEPEM_identifiedLepton_nofRingHits->Fill(nofringhits);
2894  fhEPEM_identifiedLepton_ringMid->Fill(richRing->GetCenterX(),
2895  richRing->GetCenterY());
2896  fhEPEM_identifiedLepton_ringRadius->Fill(richRing->GetRadius());
2897  // CbmRichRing::GetDistance() method is no longer supported
2898  // If you wan to use cuts update code using CbmRichUtil::GetRingTrackDistance()
2899  fhEPEM_identifiedLepton_distance->Fill(1. /*richRing->GetDistance() */);
2900  fhEPEM_identifiedLepton_selectionNN->Fill(richRing->GetSelectionNN());
2901  fhEPEM_identifiedLepton_radius_vs_p->Fill(fRecoRefittedMomentum[e3].Mag(),
2902  richRing->GetRadius());
2903  fhEPEM_identifiedLepton_ANNvalue->Fill(ANNvalueE3);
2904  }
2905  if (IsRichElectron4ann > 0) {
2906  CbmGlobalTrack* gTrack =
2907  (CbmGlobalTrack*) fGlobalTracks->At(fRecoTracklistEPEM_gtid[e4]);
2908  int richInd = gTrack->GetRichRingIndex();
2909  CbmRichRing* richRing = (CbmRichRing*) fRichRings->At(richInd);
2910  Int_t nofringhits = richRing->GetNofHits();
2911  fhEPEM_identifiedLepton_nofRingHits->Fill(nofringhits);
2912  fhEPEM_identifiedLepton_ringMid->Fill(richRing->GetCenterX(),
2913  richRing->GetCenterY());
2914  fhEPEM_identifiedLepton_ringRadius->Fill(richRing->GetRadius());
2915  // CbmRichRing::GetDistance() method is no longer supported
2916  // If you wan to use cuts update code using CbmRichUtil::GetRingTrackDistance()
2917  fhEPEM_identifiedLepton_distance->Fill(1. /*richRing->GetDistance()*/);
2918  fhEPEM_identifiedLepton_selectionNN->Fill(richRing->GetSelectionNN());
2919  fhEPEM_identifiedLepton_radius_vs_p->Fill(fRecoRefittedMomentum[e4].Mag(),
2920  richRing->GetRadius());
2921  fhEPEM_identifiedLepton_ANNvalue->Fill(ANNvalueE4);
2922  }
2923 
2924 
2925  fhEPEM_pi0_ANNvalues_noCuts->Fill(ANNvalueE1);
2926  fhEPEM_pi0_ANNvalues_noCuts->Fill(ANNvalueE2);
2927  fhEPEM_pi0_ANNvalues_noCuts->Fill(ANNvalueE3);
2928  fhEPEM_pi0_ANNvalues_noCuts->Fill(ANNvalueE4);
2929  if (OpeningAngleCut1 && OpeningAngleCut2) {
2930  fhEPEM_pi0_ANNvalues_angleCut->Fill(ANNvalueE1);
2931  fhEPEM_pi0_ANNvalues_angleCut->Fill(ANNvalueE2);
2932  fhEPEM_pi0_ANNvalues_angleCut->Fill(ANNvalueE3);
2933  fhEPEM_pi0_ANNvalues_angleCut->Fill(ANNvalueE4);
2934  }
2935 
2936  fhEPEM_pi0_ANNefficiencies->Fill(0);
2937  if (ANNvalueE1 > -1 && ANNvalueE2 > -1 && ANNvalueE3 > -1
2938  && ANNvalueE4 > -1) {
2939  fhEPEM_pi0_ANNefficiencies->Fill(1);
2940  }
2941  if (ANNvalueE1 > -0.9 && ANNvalueE2 > -0.9 && ANNvalueE3 > -0.9
2942  && ANNvalueE4 > -0.9) {
2943  fhEPEM_pi0_ANNefficiencies->Fill(2);
2944  }
2945  if (ANNvalueE1 > -0.8 && ANNvalueE2 > -0.8 && ANNvalueE3 > -0.8
2946  && ANNvalueE4 > -0.8) {
2947  fhEPEM_pi0_ANNefficiencies->Fill(3);
2948  }
2949  if (ANNvalueE1 > -0.7 && ANNvalueE2 > -0.7 && ANNvalueE3 > -0.7
2950  && ANNvalueE4 > -0.7) {
2951  fhEPEM_pi0_ANNefficiencies->Fill(4);
2952  }
2953  if (ANNvalueE1 > -0.6 && ANNvalueE2 > -0.6 && ANNvalueE3 > -0.6
2954  && ANNvalueE4 > -0.6) {
2955  fhEPEM_pi0_ANNefficiencies->Fill(5);
2956  }
2957  if (ANNvalueE1 > -0.5 && ANNvalueE2 > -0.5 && ANNvalueE3 > -0.5
2958  && ANNvalueE4 > -0.5) {
2959  fhEPEM_pi0_ANNefficiencies->Fill(6);
2960  }
2961  if (ANNvalueE1 > -0.4 && ANNvalueE2 > -0.4 && ANNvalueE3 > -0.4
2962  && ANNvalueE4 > -0.4) {
2963  fhEPEM_pi0_ANNefficiencies->Fill(7);
2964  }
2965  if (ANNvalueE1 > -0.3 && ANNvalueE2 > -0.3 && ANNvalueE3 > -0.3
2966  && ANNvalueE4 > -0.3) {
2967  fhEPEM_pi0_ANNefficiencies->Fill(8);
2968  }
2969  if (ANNvalueE1 > -0.2 && ANNvalueE2 > -0.2 && ANNvalueE3 > -0.2
2970  && ANNvalueE4 > -0.2) {
2971  fhEPEM_pi0_ANNefficiencies->Fill(9);
2972  }
2973  if (ANNvalueE1 > -0.1 && ANNvalueE2 > -0.1 && ANNvalueE3 > -0.1
2974  && ANNvalueE4 > -0.1) {
2975  fhEPEM_pi0_ANNefficiencies->Fill(10);
2976  }
2977  if (ANNvalueE1 > -0.0 && ANNvalueE2 > -0.0 && ANNvalueE3 > -0.0
2978  && ANNvalueE4 > -0.0) {
2979  fhEPEM_pi0_ANNefficiencies->Fill(11);
2980  }
2981 
2982  fhEPEM_efficiencyCuts->Fill(0); // no further cuts applied
2983  // first ANN usage for electron identification
2984  if (
2985  IsRichElectron1ann && IsRichElectron2ann && IsRichElectron3ann
2986  && IsRichElectron4ann) { // all 4 electrons correctly identified with the RICH via ANN
2987  fhEPEM_efficiencyCuts->Fill(1);
2988  if (AllWithinChiCut) { // refitted momenta are within chi cut
2989  fhEPEM_efficiencyCuts->Fill(2);
2990  if (OpeningAngleCut1
2991  && OpeningAngleCut2) { // opening angle of e+e- pairs below x
2992  fhEPEM_efficiencyCuts->Fill(3);
2993  if (InvariantMassCut1 && InvariantMassCut2) {
2994  fhEPEM_efficiencyCuts->Fill(4);
2995  }
2996  }
2997  }
2998  }
2999  // then standard method for electron identification
3000  if (
3001  IsRichElectron1normal && IsRichElectron2normal && IsRichElectron3normal
3002  && IsRichElectron4normal) { // all 4 electrons correctly identified with the RICH via "normal way"
3003  fhEPEM_efficiencyCuts->Fill(5);
3004  if (AllWithinChiCut) { // refitted momenta are within chi cut
3005  fhEPEM_efficiencyCuts->Fill(6);
3006  if (OpeningAngleCut1
3007  && OpeningAngleCut2) { // opening angle of e+e- pairs below x
3008  fhEPEM_efficiencyCuts->Fill(7);
3009  if (InvariantMassCut1 && InvariantMassCut2) {
3010  fhEPEM_efficiencyCuts->Fill(8);
3011  }
3012  }
3013  }
3014  }
3015  // MC-true data for electron identification
3016  if (
3017  IsRichElectron1MC && IsRichElectron2MC && IsRichElectron3MC
3018  && IsRichElectron4MC) { // all 4 electrons correctly identified with the RICH via MC-true data
3019  fhEPEM_efficiencyCuts->Fill(9);
3020  if (AllWithinChiCut) { // refitted momenta are within chi cut
3021  fhEPEM_efficiencyCuts->Fill(10);
3022  if (OpeningAngleCut1
3023  && OpeningAngleCut2) { // opening angle of e+e- pairs below x
3024  fhEPEM_efficiencyCuts->Fill(11);
3025  if (InvariantMassCut1 && InvariantMassCut2) {
3026  fhEPEM_efficiencyCuts->Fill(12);
3027  }
3028  }
3029  }
3030  }
3031  }
3032 
3033  if (IsEta == 1) {
3034  fhEPEM_efficiencyCuts_eta->Fill(0); // no further cuts applied
3035  // first ANN usage for electron identification
3036  if (
3037  IsRichElectron1ann && IsRichElectron2ann && IsRichElectron3ann
3038  && IsRichElectron4ann) { // all 4 electrons correctly identified with the RICH via ANN
3039  fhEPEM_efficiencyCuts_eta->Fill(1);
3040  if (AllWithinChiCut) { // refitted momenta are within chi cut
3041  fhEPEM_efficiencyCuts_eta->Fill(2);
3042  if (OpeningAngleCut1
3043  && OpeningAngleCut2) { // opening angle of e+e- pairs below x
3044  fhEPEM_efficiencyCuts_eta->Fill(3);
3045  if (InvariantMassCut1 && InvariantMassCut2) {
3046  fhEPEM_efficiencyCuts_eta->Fill(4);
3047  }
3048  }
3049  }
3050  }
3051  // then standard method for electron identification
3052  if (
3053  IsRichElectron1normal && IsRichElectron2normal && IsRichElectron3normal
3054  && IsRichElectron4normal) { // all 4 electrons correctly identified with the RICH via "normal way"
3055  fhEPEM_efficiencyCuts_eta->Fill(5);
3056  if (AllWithinChiCut) { // refitted momenta are within chi cut
3057  fhEPEM_efficiencyCuts_eta->Fill(6);
3058  if (OpeningAngleCut1
3059  && OpeningAngleCut2) { // opening angle of e+e- pairs below x
3060  fhEPEM_efficiencyCuts_eta->Fill(7);
3061  if (InvariantMassCut1 && InvariantMassCut2) {
3062  fhEPEM_efficiencyCuts_eta->Fill(8);
3063  }
3064  }
3065  }
3066  }
3067  // MC-true data for electron identification
3068  if (
3069  IsRichElectron1MC && IsRichElectron2MC && IsRichElectron3MC
3070  && IsRichElectron4MC) { // all 4 electrons correctly identified with the RICH via MC-true data
3071  fhEPEM_efficiencyCuts_eta->Fill(9);
3072  if (AllWithinChiCut) { // refitted momenta are within chi cut
3073  fhEPEM_efficiencyCuts_eta->Fill(10);
3074  if (OpeningAngleCut1
3075  && OpeningAngleCut2) { // opening angle of e+e- pairs below x
3076  fhEPEM_efficiencyCuts_eta->Fill(11);
3077  if (InvariantMassCut1 && InvariantMassCut2) {
3078  fhEPEM_efficiencyCuts_eta->Fill(12);
3079  }
3080  }
3081  }
3082  }
3083  }
3084 }
3085 
3086 
3087 Int_t CbmAnaConversionReco::NofDaughters(Int_t motherId) {
3088  Int_t nofDaughters = 0;
3089  for (unsigned int i = 0; i < fRecoTracklistEPEM.size(); i++) {
3090  Int_t motherId_temp = fRecoTracklistEPEM[i]->GetMotherId();
3091  if (motherId == motherId_temp) nofDaughters++;
3092  }
3093 
3094 
3095  return nofDaughters;
3096 }
3097 
3098 
3099 Double_t CbmAnaConversionReco::CalculateOpeningAngleReco(TVector3 electron1,
3100  TVector3 electron2) {
3101  Double_t energyP = TMath::Sqrt(electron1.Mag2() + M2E);
3102  TLorentzVector lorVecP(electron1, energyP);
3103 
3104  Double_t energyM = TMath::Sqrt(electron2.Mag2() + M2E);
3105  TLorentzVector lorVecM(electron2, energyM);
3106 
3107  Double_t anglePair = lorVecM.Angle(lorVecP.Vect());
3108  Double_t theta = 180. * anglePair / TMath::Pi();
3109 
3110  return theta;
3111 }
3112 
3113 
3114 Double_t CbmAnaConversionReco::CalculateOpeningAngleMC(CbmMCTrack* mctrack1,
3115  CbmMCTrack* mctrack2) {
3116  TVector3 electron1;
3117  mctrack1->GetMomentum(electron1);
3118  Double_t energyP = TMath::Sqrt(electron1.Mag2() + M2E);
3119  TLorentzVector lorVecP(electron1, energyP);
3120 
3121  TVector3 electron2;
3122  mctrack2->GetMomentum(electron2);
3123  Double_t energyM = TMath::Sqrt(electron2.Mag2() + M2E);
3124  TLorentzVector lorVecM(electron2, energyM);
3125 
3126  Double_t anglePair = lorVecM.Angle(lorVecP.Vect());
3127  Double_t theta = 180. * anglePair / TMath::Pi();
3128 
3129  return theta;
3130 }
3131 
3132 
3133 Double_t CbmAnaConversionReco::CalculateOpeningAngleBetweenGammasMC(
3134  CbmMCTrack* mctrack1,
3135  CbmMCTrack* mctrack2,
3136  CbmMCTrack* mctrack3,
3137  CbmMCTrack* mctrack4) {
3138  Double_t openingAngle;
3139  TLorentzVector gamma1;
3140  TLorentzVector gamma2;
3141 
3142  if (mctrack1->GetMotherId() == mctrack2->GetMotherId()
3143  && mctrack3->GetMotherId() == mctrack4->GetMotherId()) {
3144  CbmMCTrack* mother1 = (CbmMCTrack*) fMcTracks->At(mctrack1->GetMotherId());
3145  mother1->Get4Momentum(gamma1);
3146  CbmMCTrack* mother2 = (CbmMCTrack*) fMcTracks->At(mctrack3->GetMotherId());
3147  mother2->Get4Momentum(gamma2);
3148  }
3149  if (mctrack1->GetMotherId() == mctrack3->GetMotherId()
3150  && mctrack2->GetMotherId() == mctrack4->GetMotherId()) {
3151  CbmMCTrack* mother1 = (CbmMCTrack*) fMcTracks->At(mctrack1->GetMotherId());
3152  mother1->Get4Momentum(gamma1);
3153  CbmMCTrack* mother2 = (CbmMCTrack*) fMcTracks->At(mctrack2->GetMotherId());
3154  mother2->Get4Momentum(gamma2);
3155  }
3156  if (mctrack1->GetMotherId() == mctrack4->GetMotherId()
3157  && mctrack2->GetMotherId() == mctrack3->GetMotherId()) {
3158  CbmMCTrack* mother1 = (CbmMCTrack*) fMcTracks->At(mctrack1->GetMotherId());
3159  mother1->Get4Momentum(gamma1);
3160  CbmMCTrack* mother2 = (CbmMCTrack*) fMcTracks->At(mctrack2->GetMotherId());
3161  mother2->Get4Momentum(gamma2);
3162  }
3163 
3164  Double_t angle = gamma1.Angle(gamma2.Vect());
3165  openingAngle = 180. * angle / TMath::Pi();
3166 
3167  return openingAngle;
3168 }
3169 
3170 
3171 Double_t CbmAnaConversionReco::CalculateOpeningAngleBetweenGammasReco(
3172  TVector3 electron1,
3173  TVector3 electron2,
3174  TVector3 electron3,
3175  TVector3 electron4) {
3176  Double_t energy1 = TMath::Sqrt(electron1.Mag2() + M2E);
3177  TLorentzVector lorVec1(electron1, energy1);
3178 
3179  Double_t energy2 = TMath::Sqrt(electron2.Mag2() + M2E);
3180  TLorentzVector lorVec2(electron2, energy2);
3181 
3182  Double_t energy3 = TMath::Sqrt(electron3.Mag2() + M2E);
3183  TLorentzVector lorVec3(electron3, energy3);
3184 
3185  Double_t energy4 = TMath::Sqrt(electron4.Mag2() + M2E);
3186  TLorentzVector lorVec4(electron4, energy4);
3187 
3188  TLorentzVector lorPhoton1 = lorVec1 + lorVec2;
3189  TLorentzVector lorPhoton2 = lorVec3 + lorVec4;
3190 
3191  Double_t angleBetweenPhotons = lorPhoton1.Angle(lorPhoton2.Vect());
3192  Double_t theta = 180. * angleBetweenPhotons / TMath::Pi();
3193 
3194  return theta;
3195 }
3196 
3197 
3198 void CbmAnaConversionReco::CalculateInvMassWithFullRecoCuts() {
3199  Int_t nof = fRecoTracklistEPEM.size();
3200  cout
3201  << "CbmAnaConversionReco::CalculateInvMassWithFullRecoCuts: nof entries - "
3202  << nof << endl;
3203  if (nof >= 4) {
3204  for (int a = 0; a < nof - 3; a++) {
3205  for (int b = a + 1; b < nof - 2; b++) {
3206  for (int c = b + 1; c < nof - 1; c++) {
3207  for (int d = c + 1; d < nof; d++) {
3208  Int_t check1 = (fRecoTracklistEPEM[a]->GetPdgCode() > 0);
3209  Int_t check2 = (fRecoTracklistEPEM[b]->GetPdgCode() > 0);
3210  Int_t check3 = (fRecoTracklistEPEM[c]->GetPdgCode() > 0);
3211  Int_t check4 = (fRecoTracklistEPEM[d]->GetPdgCode() > 0);
3212  Int_t test = check1 + check2 + check3 + check4;
3213  if (test != 2) continue; // need two electrons and two positrons
3214 
3215 
3216  Double_t invmass = Invmass_4particlesRECO(fRecoRefittedMomentum[a],
3217  fRecoRefittedMomentum[b],
3218  fRecoRefittedMomentum[c],
3219  fRecoRefittedMomentum[d]);
3220  //fhElectrons_invmass->Fill(invmass);
3221 
3222 
3223  CbmLmvmKinematicParams params1 = CalculateKinematicParamsReco(
3224  fRecoRefittedMomentum[a], fRecoRefittedMomentum[b]);
3225  CbmLmvmKinematicParams params2 = CalculateKinematicParamsReco(
3226  fRecoRefittedMomentum[a], fRecoRefittedMomentum[c]);
3227  CbmLmvmKinematicParams params3 = CalculateKinematicParamsReco(
3228  fRecoRefittedMomentum[a], fRecoRefittedMomentum[d]);
3229  CbmLmvmKinematicParams params4 = CalculateKinematicParamsReco(
3230  fRecoRefittedMomentum[b], fRecoRefittedMomentum[c]);
3231  CbmLmvmKinematicParams params5 = CalculateKinematicParamsReco(
3232  fRecoRefittedMomentum[b], fRecoRefittedMomentum[d]);
3233  CbmLmvmKinematicParams params6 = CalculateKinematicParamsReco(
3234  fRecoRefittedMomentum[c], fRecoRefittedMomentum[d]);
3235 
3236  Double_t openingAngleCut = 1;
3237  Int_t IsPhoton_openingAngle1 = (params1.fAngle < openingAngleCut);
3238  Int_t IsPhoton_openingAngle2 = (params2.fAngle < openingAngleCut);
3239  Int_t IsPhoton_openingAngle3 = (params3.fAngle < openingAngleCut);
3240  Int_t IsPhoton_openingAngle4 = (params4.fAngle < openingAngleCut);
3241  Int_t IsPhoton_openingAngle5 = (params5.fAngle < openingAngleCut);
3242  Int_t IsPhoton_openingAngle6 = (params6.fAngle < openingAngleCut);
3243 
3244  Double_t invMassCut = 0.03;
3245  Int_t IsPhoton_invMass1 = (params1.fMinv < invMassCut);
3246  Int_t IsPhoton_invMass2 = (params2.fMinv < invMassCut);
3247  Int_t IsPhoton_invMass3 = (params3.fMinv < invMassCut);
3248  Int_t IsPhoton_invMass4 = (params4.fMinv < invMassCut);
3249  Int_t IsPhoton_invMass5 = (params5.fMinv < invMassCut);
3250  Int_t IsPhoton_invMass6 = (params6.fMinv < invMassCut);
3251 
3252  if (IsPhoton_openingAngle1 && IsPhoton_openingAngle6
3253  && IsPhoton_invMass1 && IsPhoton_invMass6
3254  && (check1 + check2 == 1) && (check3 + check4 == 1)) {
3255  fhInvMassWithFullRecoCuts->Fill(invmass);
3256  }
3257  if (IsPhoton_openingAngle2 && IsPhoton_openingAngle5
3258  && IsPhoton_invMass2 && IsPhoton_invMass5
3259  && (check1 + check3 == 1) && (check2 + check4 == 1)) {
3260  fhInvMassWithFullRecoCuts->Fill(invmass);
3261  }
3262  if (IsPhoton_openingAngle3 && IsPhoton_openingAngle4
3263  && IsPhoton_invMass3 && IsPhoton_invMass4
3264  && (check1 + check4 == 1) && (check2 + check3 == 1)) {
3265  fhInvMassWithFullRecoCuts->Fill(invmass);
3266  }
3267  }
3268  }
3269  }
3270  }
3271  }
3272 }
3273 
3274 
3275 CbmLmvmKinematicParams
3276 CbmAnaConversionReco::CalculateKinematicParamsReco(const TVector3 electron1,
3277  const TVector3 electron2) {
3278  CbmLmvmKinematicParams params;
3279 
3280  Double_t energyP = TMath::Sqrt(electron1.Mag2() + M2E);
3281  TLorentzVector lorVecP(electron1, energyP);
3282 
3283  Double_t energyM = TMath::Sqrt(electron2.Mag2() + M2E);
3284  TLorentzVector lorVecM(electron2, energyM);
3285 
3286  TVector3 momPair = electron1 + electron2;
3287  Double_t energyPair = energyP + energyM;
3288  Double_t ptPair = momPair.Perp();
3289  Double_t pzPair = momPair.Pz();
3290  Double_t yPair =
3291  0.5 * TMath::Log((energyPair + pzPair) / (energyPair - pzPair));
3292  Double_t anglePair = lorVecM.Angle(lorVecP.Vect());
3293  Double_t theta = 180. * anglePair / TMath::Pi();
3294  Double_t minv = 2. * TMath::Sin(anglePair / 2.)
3295  * TMath::Sqrt(electron1.Mag() * electron2.Mag());
3296 
3297  params.fMomentumMag = momPair.Mag();
3298  params.fPt = ptPair;
3299  params.fRapidity = yPair;
3300  params.fMinv = minv;
3301  params.fAngle = theta;
3302  return params;
3303 }
3304 
3305 
3306 CbmLmvmKinematicParams
3307 CbmAnaConversionReco::CalculateKinematicParams_4particles(
3308  const TVector3 part1,
3309  const TVector3 part2,
3310  const TVector3 part3,
3311  const TVector3 part4) {
3312  CbmLmvmKinematicParams params;
3313 
3314  Double_t energy1 = TMath::Sqrt(part1.Mag2() + M2E);
3315  TLorentzVector lorVec1(part1, energy1);
3316 
3317  Double_t energy2 = TMath::Sqrt(part2.Mag2() + M2E);
3318  TLorentzVector lorVec2(part2, energy2);
3319 
3320  Double_t energy3 = TMath::Sqrt(part3.Mag2() + M2E);
3321  TLorentzVector lorVec3(part3, energy3);
3322 
3323  Double_t energy4 = TMath::Sqrt(part4.Mag2() + M2E);
3324  TLorentzVector lorVec4(part4, energy4);
3325 
3326  TLorentzVector sum;
3327  sum = lorVec1 + lorVec2 + lorVec3 + lorVec4;
3328 
3329  TVector3 momPair = part1 + part2 + part3 + part4;
3330  Double_t energyPair = energy1 + energy2 + energy3 + energy4;
3331  Double_t ptPair = momPair.Perp();
3332  Double_t pzPair = momPair.Pz();
3333  Double_t yPair =
3334  0.5 * TMath::Log((energyPair + pzPair) / (energyPair - pzPair));
3335  // Double_t anglePair = 0;
3336  // Double_t theta = 180.*anglePair/TMath::Pi();
3337  Double_t minv = sum.Mag();
3338 
3339  params.fMomentumMag = momPair.Mag();
3340  params.fPt = ptPair;
3341  params.fRapidity = yPair;
3342  params.fMinv = minv;
3343  params.fAngle = 0;
3344  return params;
3345 }
3346 
3347 
3348 Bool_t CbmAnaConversionReco::IsRichElectronANN(Int_t globalTrackIndex,
3349  Double_t momentum) {
3350  if (NULL == fGlobalTracks || NULL == fRichRings) return false;
3351  //CbmGlobalTrack* globalTrack = (CbmGlobalTrack*) fGlobalTracks->At(globalTrackIndex);
3352  const CbmGlobalTrack* globalTrack =
3353  static_cast<const CbmGlobalTrack*>(fGlobalTracks->At(globalTrackIndex));
3354  Int_t richId = globalTrack->GetRichRingIndex();
3355  if (richId < 0) return false;
3356  CbmRichRing* ring = static_cast<CbmRichRing*>(fRichRings->At(richId));
3357  if (NULL == ring) return false;
3358 
3359  Double_t fRichAnnCut = -0.8;
3360  Double_t ann = CbmRichElectronIdAnn::GetInstance().CalculateAnnValue(
3361  globalTrackIndex, momentum);
3362  if (ann > fRichAnnCut)
3363  return true;
3364  else
3365  return false;
3366 }
3367 
3368 
3369 Double_t CbmAnaConversionReco::ElectronANNvalue(Int_t globalTrackIndex,
3370  Double_t momentum) {
3371  if (NULL == fGlobalTracks || NULL == fRichRings) return -2;
3372  //CbmGlobalTrack* globalTrack = (CbmGlobalTrack*) fGlobalTracks->At(globalTrackIndex);
3373  const CbmGlobalTrack* globalTrack =
3374  static_cast<const CbmGlobalTrack*>(fGlobalTracks->At(globalTrackIndex));
3375  Int_t richId = globalTrack->GetRichRingIndex();
3376  if (richId < 0) return -2;
3377  CbmRichRing* ring = static_cast<CbmRichRing*>(fRichRings->At(richId));
3378  if (NULL == ring) return -2;
3379 
3380  Double_t ann = CbmRichElectronIdAnn::GetInstance().CalculateAnnValue(
3381  globalTrackIndex, momentum);
3382  return ann;
3383 }
3384 
3385 
3386 Bool_t CbmAnaConversionReco::IsRichElectronNormal(Int_t globalTrackIndex,
3387  Double_t momentum) {
3388  if (NULL == fGlobalTracks || NULL == fRichRings) return false;
3389  CbmGlobalTrack* globalTrack =
3390  (CbmGlobalTrack*) fGlobalTracks->At(globalTrackIndex);
3391  Int_t richId = globalTrack->GetRichRingIndex();
3392  if (richId < 0) return false;
3393  CbmRichRing* ring = static_cast<CbmRichRing*>(fRichRings->At(richId));
3394  if (NULL == ring) return false;
3395 
3396  Double_t axisA = ring->GetAaxis();
3397  Double_t axisB = ring->GetBaxis();
3398  Double_t dist = CbmRichUtil::GetRingTrackDistance(globalTrackIndex);
3399 
3400  Double_t fMeanA = 4.95;
3401  Double_t fMeanB = 4.54;
3402  Double_t fRmsA = 0.30;
3403  Double_t fRmsB = 0.22;
3404  Double_t fRmsCoeff = 4.5;
3405  Double_t fDistCut = 1.;
3406 
3407  Bool_t isElectronRICH = 0;
3408 
3409  if (momentum < 5.5) {
3410  if (fabs(axisA - fMeanA) < fRmsCoeff * fRmsA
3411  && fabs(axisB - fMeanB) < fRmsCoeff * fRmsB && dist < fDistCut)
3412  isElectronRICH = 1;
3413  } else {
3415  // Double_t polAaxis = 5.80008 - 4.10118 / (momentum - 3.67402);
3416  // Double_t polBaxis = 5.58839 - 4.75980 / (momentum - 3.31648);
3417  // Double_t polRaxis = 5.87252 - 7.64641/(momentum - 1.62255);
3419  Double_t polAaxis = 5.64791 - 4.24077 / (momentum - 3.65494);
3420  Double_t polBaxis = 5.41106 - 4.49902 / (momentum - 3.52450);
3421  //Double_t polRaxis = 5.66516 - 6.62229/(momentum - 2.25304);
3422  if (axisA < (fMeanA + fRmsCoeff * fRmsA) && axisA > polAaxis
3423  && axisB < (fMeanB + fRmsCoeff * fRmsB) && axisB > polBaxis
3424  && dist < fDistCut)
3425  isElectronRICH = 1;
3426  }
3427  return isElectronRICH;
3428 }
CbmMCTrack::GetMotherId
Int_t GetMotherId() const
Definition: CbmMCTrack.h:71
CbmMCTrack::GetMomentum
void GetMomentum(TVector3 &momentum) const
Definition: CbmMCTrack.h:172
CbmMCTrack::GetMass
Double_t GetMass() const
Mass of the associated particle.
Definition: CbmMCTrack.cxx:114
CbmAnaConversionReco::~CbmAnaConversionReco
virtual ~CbmAnaConversionReco()
Definition: CbmAnaConversionReco.cxx:162
CbmAnaConversionReco::fhEPEM_pi0_ANNvalues_angleCut
TH1D * fhEPEM_pi0_ANNvalues_angleCut
Definition: CbmAnaConversionReco.h:176
CbmAnaConversionReco::fhEPEM_pt_vs_p_all_mc
TH2D * fhEPEM_pt_vs_p_all_mc
Definition: CbmAnaConversionReco.h:223
CbmAnaConversionReco::fMCTracklist_all
std::vector< CbmMCTrack * > fMCTracklist_all
Definition: CbmAnaConversionReco.h:98
CbmAnaConversionReco::fhInvariantMass_MC_all
TH1D * fhInvariantMass_MC_all
Definition: CbmAnaConversionReco.h:120
CbmAnaConversionReco::fhEPEM_invmass_all_mc
TH1D * fhEPEM_invmass_all_mc
Definition: CbmAnaConversionReco.h:148
CbmPixelHit::GetX
Double_t GetX() const
Definition: CbmPixelHit.h:83
CbmAnaConversionReco::fhEPEM_invmass_eeee_mc
TH1D * fhEPEM_invmass_eeee_mc
Definition: CbmAnaConversionReco.h:146
CbmAnaConversionReco::fHistoList_MC
std::vector< TH1 * > fHistoList_MC
Definition: CbmAnaConversionReco.h:108
CbmAnaConversionReco::fRecoTracklistEPEM_ids
std::vector< int > fRecoTracklistEPEM_ids
Definition: CbmAnaConversionReco.h:100
CbmAnaConversionReco::CalculateOpeningAngleReco
Double_t CalculateOpeningAngleReco(TVector3 electron1, TVector3 electron2)
Definition: CbmAnaConversionReco.cxx:3099
CbmAnaConversionReco::fhEPEM_invmass_eta_mc
TH1D * fhEPEM_invmass_eta_mc
Definition: CbmAnaConversionReco.h:248
CbmAnaConversionReco::fhEPEM_missingLepton_ringMid
TH2D * fhEPEM_missingLepton_ringMid
Definition: CbmAnaConversionReco.h:229
CbmAnaConversionReco::fhEPEM_invmass_gg_refitted
TH1D * fhEPEM_invmass_gg_refitted
Definition: CbmAnaConversionReco.h:143
CbmAnaConversionReco::fhMC_electrons_theta_vs_p
TH2D * fhMC_electrons_theta_vs_p
Definition: CbmAnaConversionReco.h:130
CbmPixelHit::GetY
Double_t GetY() const
Definition: CbmPixelHit.h:84
CbmAnaConversionReco::fhEPEM_identifiedLepton_distance
TH1D * fhEPEM_identifiedLepton_distance
Definition: CbmAnaConversionReco.h:240
CbmAnaConversionReco::fhUsedMomenta_vsZ_refitted
TH2D * fhUsedMomenta_vsZ_refitted
Definition: CbmAnaConversionReco.h:204
CbmAnaConversionReco::CbmAnaConversionReco
CbmAnaConversionReco()
Definition: CbmAnaConversionReco.cxx:36
CbmAnaConversionReco::fhInvMass_EPEM_error_stsMomVec
TH1D * fhInvMass_EPEM_error_stsMomVec
Definition: CbmAnaConversionReco.h:186
CbmAnaConversionReco::IsRichElectronANN
Bool_t IsRichElectronANN(Int_t globalTrackIndex, Double_t momentum)
Definition: CbmAnaConversionReco.cxx:3348
CbmAnaConversionReco::InvariantMassMC_all
void InvariantMassMC_all()
Definition: CbmAnaConversionReco.cxx:1121
CbmAnaConversionReco::fGlobalTracks
TClonesArray * fGlobalTracks
Definition: CbmAnaConversionReco.h:93
CbmAnaConversionReco::fhInvMassWithFullRecoCuts
TH1D * fhInvMassWithFullRecoCuts
Definition: CbmAnaConversionReco.h:214
CbmAnaConversionReco::fhInvMass_EPEM_openingAngleRef
TH1D * fhInvMass_EPEM_openingAngleRef
Definition: CbmAnaConversionReco.h:188
CbmGlobalTrack::GetRichRingIndex
Int_t GetRichRingIndex() const
Definition: CbmGlobalTrack.h:41
CbmAnaConversionReco::fHistoList_all
std::vector< TH1 * > fHistoList_all
Definition: CbmAnaConversionReco.h:116
CbmAnaConversionReco::fHistoList_gee
std::vector< TH1 * > fHistoList_gee
Definition: CbmAnaConversionReco.h:114
CbmAnaConversionReco::fhInvariantMass_MC_etaPrime
TH1D * fhInvariantMass_MC_etaPrime
Definition: CbmAnaConversionReco.h:126
CbmAnaConversionReco::fhEPEM_efficiencyCuts_eta
TH1D * fhEPEM_efficiencyCuts_eta
Definition: CbmAnaConversionReco.h:250
CbmAnaConversionReco::fhEPEM_InDetector_invmass_all_mc
TH1D * fhEPEM_InDetector_invmass_all_mc
Definition: CbmAnaConversionReco.h:220
CbmAnaConversionReco::fhPi0_startvertexElectrons_all
TH1D * fhPi0_startvertexElectrons_all
Definition: CbmAnaConversionReco.h:208
CbmMCTrack::GetPdgCode
Int_t GetPdgCode() const
Definition: CbmMCTrack.h:70
CbmAnaConversionReco::fhEPEM_missingLepton_rings
TH2D * fhEPEM_missingLepton_rings
Definition: CbmAnaConversionReco.h:233
CbmAnaConversionReco::fhUsedMomenta_vsY_refitted
TH2D * fhUsedMomenta_vsY_refitted
Definition: CbmAnaConversionReco.h:202
i
int i
Definition: L1/vectors/P4_F32vec4.h:25
CbmAnaConversionReco::fhInvMass_EPEM_refitted
TH1D * fhInvMass_EPEM_refitted
Definition: CbmAnaConversionReco.h:185
CbmAnaConversionReco::CalculateKinematicParamsReco
CbmLmvmKinematicParams CalculateKinematicParamsReco(const TVector3 electron1, const TVector3 electron2)
Definition: CbmAnaConversionReco.cxx:3276
CbmAnaConversionGlobalFunctions::OpeningAngleBetweenGamma
static Double_t OpeningAngleBetweenGamma(const TVector3 part11, const TVector3 part12, const TVector3 part21, const TVector3 part22)
Definition: CbmAnaConversionGlobalFunctions.h:83
CbmAnaConversionReco::fhUsedMomenta_vsX_refitted
TH2D * fhUsedMomenta_vsX_refitted
Definition: CbmAnaConversionReco.h:200
CbmGlobalTrack.h
CbmAnaConversionReco::NofDaughters
Int_t NofDaughters(Int_t motherId)
Definition: CbmAnaConversionReco.cxx:3087
CbmAnaConversionReco::fhEPEM_invmass_gg_mc
TH1D * fhEPEM_invmass_gg_mc
Definition: CbmAnaConversionReco.h:142
CbmAnaConversionReco::fhInvariantMass_pi0epem
TH1D * fhInvariantMass_pi0epem
Definition: CbmAnaConversionReco.h:206
CbmAnaConversionReco::fhPi0_pt_all
TH1D * fhPi0_pt_all
Definition: CbmAnaConversionReco.h:170
CbmAnaConversionReco::fhUsedMomenta_error_refitted
TH1D * fhUsedMomenta_error_refitted
Definition: CbmAnaConversionReco.h:192
CbmAnaConversionReco::fhPi0_pt_vs_rap_gg
TH2D * fhPi0_pt_vs_rap_gg
Definition: CbmAnaConversionReco.h:164
CbmAnaConversionReco::CutEfficiencyStudies
void CutEfficiencyStudies(int e1, int e2, int e3, int e4, int motherE1, int motherE2, int motherE3, int motherE4, int IsEta=0)
Definition: CbmAnaConversionReco.cxx:2658
CbmAnaConversionReco::CalculateOpeningAngleBetweenGammasMC
Double_t CalculateOpeningAngleBetweenGammasMC(CbmMCTrack *mctrack1, CbmMCTrack *mctrack2, CbmMCTrack *mctrack3, CbmMCTrack *mctrack4)
Definition: CbmAnaConversionReco.cxx:3133
CbmRichRing::GetNofHits
Int_t GetNofHits() const
Definition: CbmRichRing.h:40
CbmRichRing
Definition: CbmRichRing.h:17
CbmAnaConversionReco::fHistoList_reco
std::vector< TH1 * > fHistoList_reco
Definition: CbmAnaConversionReco.h:110
CbmAnaConversionReco::fhEPEM_missingLepton_selectionNN
TH1D * fhEPEM_missingLepton_selectionNN
Definition: CbmAnaConversionReco.h:232
CbmAnaConversionReco::fhEPEM_openingAngle_gg_refitted
TH1D * fhEPEM_openingAngle_gg_refitted
Definition: CbmAnaConversionReco.h:152
CbmRichRing::GetSelectionNN
Double_t GetSelectionNN() const
Definition: CbmRichRing.h:94
CbmRichRing.h
CbmAnaConversionReco::fhPi0_pt_vs_rap_gee
TH2D * fhPi0_pt_vs_rap_gee
Definition: CbmAnaConversionReco.h:165
CbmAnaConversionReco::InitHistos
void InitHistos()
Definition: CbmAnaConversionReco.cxx:196
CbmAnaConversionReco::fhUsedMomenta_mc
TH1D * fhUsedMomenta_mc
Definition: CbmAnaConversionReco.h:190
CbmAnaConversionReco::fRecoTracklistEPEM
std::vector< CbmMCTrack * > fRecoTracklistEPEM
Definition: CbmAnaConversionReco.h:99
CbmLmvmKinematicParams::fMinv
Double_t fMinv
Definition: CbmLmvmKinematicParams.h:22
CbmAnaConversionReco::fhEPEM_InDetector_invmass_gee_refitted
TH1D * fhEPEM_InDetector_invmass_gee_refitted
Definition: CbmAnaConversionReco.h:219
CbmAnaConversionReco::fhUsedMomenta_errorX_stsMomVec
TH1D * fhUsedMomenta_errorX_stsMomVec
Definition: CbmAnaConversionReco.h:193
CbmAnaConversionReco::fTime
Double_t fTime
Definition: CbmAnaConversionReco.h:254
CbmAnaConversionReco::fhEPEM_missingLepton_ringRadius
TH1D * fhEPEM_missingLepton_ringRadius
Definition: CbmAnaConversionReco.h:230
CbmAnaConversionReco::fhEPEM_rap_vs_chi
TH2D * fhEPEM_rap_vs_chi
Definition: CbmAnaConversionReco.h:179
CbmAnaConversionReco::fhEPEM_openingAngle_betweenGammas_mc
TH1D * fhEPEM_openingAngle_betweenGammas_mc
Definition: CbmAnaConversionReco.h:161
CbmAnaConversionReco::fhEPEM_InDetector_invmass_gee_mc
TH1D * fhEPEM_InDetector_invmass_gee_mc
Definition: CbmAnaConversionReco.h:218
CbmAnaConversionReco::fhEPEM_pi0_ANNefficiencies
TH1D * fhEPEM_pi0_ANNefficiencies
Definition: CbmAnaConversionReco.h:177
CbmAnaConversionCutSettings.h
CbmAnaConversionReco::fhEPEM_invmass_gee_refitted
TH1D * fhEPEM_invmass_gee_refitted
Definition: CbmAnaConversionReco.h:145
CbmAnaConversionReco::fhEPEM_identifiedLepton_ringMid
TH2D * fhEPEM_identifiedLepton_ringMid
Definition: CbmAnaConversionReco.h:238
CbmRichUtil::GetRingTrackDistance
static Double_t GetRingTrackDistance(Int_t globalTrackId)
Definition: alignment/CbmRichUtil.h:20
CbmAnaConversionReco::fhInvariantMass_MC_eta
TH1D * fhInvariantMass_MC_eta
Definition: CbmAnaConversionReco.h:125
CbmAnaConversionReco::Init
void Init()
Definition: CbmAnaConversionReco.cxx:165
CbmAnaConversionReco::fhEPEM_openingAngle_gee_refitted
TH1D * fhEPEM_openingAngle_gee_refitted
Definition: CbmAnaConversionReco.h:154
CbmAnaConversionReco::fhInvMass_EPEM_mc
TH1D * fhInvMass_EPEM_mc
Definition: CbmAnaConversionReco.h:183
CbmAnaConversionReco::fhUsedMomenta_stsMomVec
TH1D * fhUsedMomenta_stsMomVec
Definition: CbmAnaConversionReco.h:189
CbmRichRing::GetHit
UInt_t GetHit(Int_t i) const
Definition: CbmRichRing.h:42
CbmAnaConversionReco::fhInvariantMass_MC_pi0_epem
TH1D * fhInvariantMass_MC_pi0_epem
Definition: CbmAnaConversionReco.h:122
CbmAnaConversionReco::fhEPEM_pt_vs_p_all_refitted
TH2D * fhEPEM_pt_vs_p_all_refitted
Definition: CbmAnaConversionReco.h:224
CbmAnaConversionReco::Finish
void Finish()
Definition: CbmAnaConversionReco.cxx:1039
CbmAnaConversionReco::CalculateKinematicParams_4particles
CbmLmvmKinematicParams CalculateKinematicParams_4particles(const TVector3 part1, const TVector3 part2, const TVector3 part3, const TVector3 part4)
Definition: CbmAnaConversionReco.cxx:3307
CbmAnaConversionReco::Invmass_4particles
Double_t Invmass_4particles(const CbmMCTrack *mctrack1, const CbmMCTrack *mctrack2, const CbmMCTrack *mctrack3, const CbmMCTrack *mctrack4)
Definition: CbmAnaConversionReco.cxx:1412
CbmAnaConversionReco::fRecoMomentum
std::vector< TVector3 > fRecoMomentum
Definition: CbmAnaConversionReco.h:103
CbmAnaConversionReco::fhUsedMomenta_errorZ_refitted
TH1D * fhUsedMomenta_errorZ_refitted
Definition: CbmAnaConversionReco.h:203
CbmAnaConversionReco::fhEPEM_identifiedLepton_ANNvalue
TH1D * fhEPEM_identifiedLepton_ANNvalue
Definition: CbmAnaConversionReco.h:244
d
double d
Definition: P4_F64vec2.h:24
CbmAnaConversionReco::fhEPEM_efficiencyCuts
TH1D * fhEPEM_efficiencyCuts
Definition: CbmAnaConversionReco.h:172
CbmAnaConversionReco::fhEPEM_identifiedLepton_ringRadius
TH1D * fhEPEM_identifiedLepton_ringRadius
Definition: CbmAnaConversionReco.h:239
CbmAnaConversionReco::ElectronANNvalue
Double_t ElectronANNvalue(Int_t globalTrackIndex, Double_t momentum)
Definition: CbmAnaConversionReco.cxx:3369
CbmAnaConversionReco::fhPi0_startvertex_vs_momentum
TH2D * fhPi0_startvertex_vs_momentum
Definition: CbmAnaConversionReco.h:212
CbmAnaConversionReco::fhEPEM_missingLepton_radius_vs_p
TH2D * fhEPEM_missingLepton_radius_vs_p
Definition: CbmAnaConversionReco.h:234
CbmAnaConversionReco::fhUsedMomenta_errorY_refitted
TH1D * fhUsedMomenta_errorY_refitted
Definition: CbmAnaConversionReco.h:201
CbmAnaConversionReco::fhEPEM_openingAngle_gg_mc
TH1D * fhEPEM_openingAngle_gg_mc
Definition: CbmAnaConversionReco.h:151
CbmAnaConversionReco::timer
TStopwatch timer
Definition: CbmAnaConversionReco.h:253
CbmAnaConversionReco::fhPi0_startvertex_vs_chi
TH2D * fhPi0_startvertex_vs_chi
Definition: CbmAnaConversionReco.h:211
CbmRichRingLight.h
CbmAnaConversionReco::fhPi0_pt_gee
TH1D * fhPi0_pt_gee
Definition: CbmAnaConversionReco.h:169
CbmAnaConversionReco::CalculateInvMassWithFullRecoCuts
void CalculateInvMassWithFullRecoCuts()
Definition: CbmAnaConversionReco.cxx:3198
CbmAnaConversionReco::fhMC_electrons_theta
TH1D * fhMC_electrons_theta
Definition: CbmAnaConversionReco.h:128
CbmAnaConversionReco::fHistoList_reco_mom
std::vector< TH1 * > fHistoList_reco_mom
Definition: CbmAnaConversionReco.h:112
CbmAnaConversionReco::fRecoTracklistEPEM_chi
std::vector< Double_t > fRecoTracklistEPEM_chi
Definition: CbmAnaConversionReco.h:101
CbmAnaConversionReco::fhEta_openingAngleGG
TH1D * fhEta_openingAngleGG
Definition: CbmAnaConversionReco.h:132
CbmAnaConversionReco::fhEPEM_identifiedLepton_radius_vs_p
TH2D * fhEPEM_identifiedLepton_radius_vs_p
Definition: CbmAnaConversionReco.h:243
CbmAnaConversionCutSettings::CalcChiCut
static Double_t CalcChiCut(Double_t pt)
Definition: CbmAnaConversionCutSettings.h:69
M2E
#define M2E
Definition: CbmAnaConversionReco.cxx:32
CbmAnaConversionReco::fMcTracks
TClonesArray * fMcTracks
Definition: CbmAnaConversionReco.h:92
CbmAnaConversionReco::fhMC_electrons_p
TH1D * fhMC_electrons_p
Definition: CbmAnaConversionReco.h:129
CbmAnaConversionReco::fhEPEM_identifiedLepton_rings
TH2D * fhEPEM_identifiedLepton_rings
Definition: CbmAnaConversionReco.h:242
CbmAnaConversionReco::fhInvariantMass_MC_pi0_gg
TH1D * fhInvariantMass_MC_pi0_gg
Definition: CbmAnaConversionReco.h:124
CbmAnaConversionReco::fRecoTracklistEPEM_gtid
std::vector< Int_t > fRecoTracklistEPEM_gtid
Definition: CbmAnaConversionReco.h:102
CbmLmvmKinematicParams::fRapidity
Double_t fRapidity
Definition: CbmLmvmKinematicParams.h:21
CbmAnaConversionReco::fHistoList_eeee
std::vector< TH1 * > fHistoList_eeee
Definition: CbmAnaConversionReco.h:115
CbmAnaConversionReco.h
CbmMCTrack::GetStartVertex
void GetStartVertex(TVector3 &vertex) const
Definition: CbmMCTrack.h:182
xMath::Pi
double Pi()
Definition: xMath.h:5
CbmAnaConversionReco::fhInvariantMass_MC_pi0_gepem
TH1D * fhInvariantMass_MC_pi0_gepem
Definition: CbmAnaConversionReco.h:123
CbmAnaConversionReco::fhEPEM_rap_vs_invmass
TH2D * fhEPEM_rap_vs_invmass
Definition: CbmAnaConversionReco.h:180
CbmRichRing::GetAaxis
Double_t GetAaxis() const
Definition: CbmRichRing.h:83
CbmAnaConversionReco::fhEPEM_openingAngle_gee_refitted_dalitz
TH1D * fhEPEM_openingAngle_gee_refitted_dalitz
Definition: CbmAnaConversionReco.h:156
CbmAnaConversionReco::fhInvMass_EPEM_error_refitted
TH1D * fhInvMass_EPEM_error_refitted
Definition: CbmAnaConversionReco.h:187
CbmAnaConversionReco::fhEPEM_openingAngle_gee_mc
TH1D * fhEPEM_openingAngle_gee_mc
Definition: CbmAnaConversionReco.h:153
CbmAnaConversionReco::fhEPEM_invmass_eeee_refitted
TH1D * fhEPEM_invmass_eeee_refitted
Definition: CbmAnaConversionReco.h:147
CbmAnaConversionReco::fhMCtest
TH1D * fhMCtest
Definition: CbmAnaConversionReco.h:138
CbmAnaConversionReco::fhMC_grandmotherPDGs
TH1D * fhMC_grandmotherPDGs
Definition: CbmAnaConversionReco.h:134
CbmAnaConversionReco::fhEPEM_invmass_all_refitted
TH1D * fhEPEM_invmass_all_refitted
Definition: CbmAnaConversionReco.h:149
CbmAnaConversionReco::fRichHits
TClonesArray * fRichHits
Definition: CbmAnaConversionReco.h:95
CbmAnaConversionReco::IsRichElectronNormal
Bool_t IsRichElectronNormal(Int_t globalTrackIndex, Double_t momentum)
Definition: CbmAnaConversionReco.cxx:3386
CbmAnaConversionReco::fhEPEM_openingAngle_vs_pt_gg_mc
TH2D * fhEPEM_openingAngle_vs_pt_gg_mc
Definition: CbmAnaConversionReco.h:158
CbmRichRing::GetBaxis
Double_t GetBaxis() const
Definition: CbmRichRing.h:84
CbmAnaConversionReco::fhUsedMomenta_errorZ_stsMomVec
TH1D * fhUsedMomenta_errorZ_stsMomVec
Definition: CbmAnaConversionReco.h:197
CbmAnaConversionReco::InvariantMassTest_4epem
void InvariantMassTest_4epem()
Definition: CbmAnaConversionReco.cxx:1521
CbmAnaConversionReco::fhPi0_startvertexElectrons_gg
TH1D * fhPi0_startvertexElectrons_gg
Definition: CbmAnaConversionReco.h:209
CbmGlobalTrack
Definition: CbmGlobalTrack.h:26
CbmRichElectronIdAnn::GetInstance
static CbmRichElectronIdAnn & GetInstance()
Definition: CbmRichElectronIdAnn.h:43
CbmAnaConversionReco::fhEPEM_missingLepton_nofRingHits
TH1D * fhEPEM_missingLepton_nofRingHits
Definition: CbmAnaConversionReco.h:228
CbmAnaConversionReco::fhEPEM_missingLepton_ANNvalue
TH1D * fhEPEM_missingLepton_ANNvalue
Definition: CbmAnaConversionReco.h:235
CbmMCTrack.h
CbmAnaConversionReco::fhEPEM_invmass_gee_mc
TH1D * fhEPEM_invmass_gee_mc
Definition: CbmAnaConversionReco.h:144
CbmAnaConversionReco::fhInvariantMass_MC_pi0
TH1D * fhInvariantMass_MC_pi0
Definition: CbmAnaConversionReco.h:121
CbmMCTrack
Definition: CbmMCTrack.h:34
fabs
friend F32vec4 fabs(const F32vec4 &a)
Definition: L1/vectors/P4_F32vec4.h:60
CbmAnaConversionReco::fhEPEM_pi0_nofLeptons_ann
TH1D * fhEPEM_pi0_nofLeptons_ann
Definition: CbmAnaConversionReco.h:174
CbmAnaConversionReco::fhEPEM_openingAngle_betweenGammas_reco
TH1D * fhEPEM_openingAngle_betweenGammas_reco
Definition: CbmAnaConversionReco.h:162
CbmAnaConversionReco::fhEPEM_efficiencyCuts2
TH1D * fhEPEM_efficiencyCuts2
Definition: CbmAnaConversionReco.h:173
CbmAnaConversionReco::SetTracklistReco
void SetTracklistReco(std::vector< CbmMCTrack * > MCTracklist, std::vector< TVector3 > RecoTracklist1, std::vector< TVector3 > RecoTracklist2, std::vector< int > ids, std::vector< Double_t > chi, std::vector< Int_t > GlobalTrackId)
Definition: CbmAnaConversionReco.cxx:1106
CbmRichElectronIdAnn::CalculateAnnValue
double CalculateAnnValue(int globalTrackIndex, double momentum)
Calculate output value of the ANN.
Definition: CbmRichElectronIdAnn.cxx:84
CbmAnaConversionReco::fhEPEM_missingLepton_distance
TH1D * fhEPEM_missingLepton_distance
Definition: CbmAnaConversionReco.h:231
CbmAnaConversionReco::fHistoList_eta
std::vector< TH1 * > fHistoList_eta
Definition: CbmAnaConversionReco.h:117
CbmAnaConversionReco::fhEPEM_invmass_eta_refitted
TH1D * fhEPEM_invmass_eta_refitted
Definition: CbmAnaConversionReco.h:249
CbmRichRing::GetRadius
Float_t GetRadius() const
Definition: CbmRichRing.h:82
CbmAnaConversionGlobalFunctions.h
CbmAnaConversionReco::fhPi0_pt_gg
TH1D * fhPi0_pt_gg
Definition: CbmAnaConversionReco.h:168
CbmRichRing::GetCenterY
Float_t GetCenterY() const
Definition: CbmRichRing.h:81
CbmLmvmKinematicParams::fMomentumMag
Double_t fMomentumMag
Definition: CbmLmvmKinematicParams.h:19
CbmAnaConversionReco::fhPi0_pt_vs_rap_all
TH2D * fhPi0_pt_vs_rap_all
Definition: CbmAnaConversionReco.h:166
CbmAnaConversionReco::fRichRings
TClonesArray * fRichRings
Definition: CbmAnaConversionReco.h:94
CbmRichElectronIdAnn.h
Implementation of the electron identification algorithm in the RICH detector using Artificial Neural ...
CbmAnaConversionReco::fhEPEM_InDetector_invmass_all_refitted
TH1D * fhEPEM_InDetector_invmass_all_refitted
Definition: CbmAnaConversionReco.h:221
CbmMCTrack::Get4Momentum
void Get4Momentum(TLorentzVector &momentum) const
Definition: CbmMCTrack.h:177
CbmAnaConversionReco::SmearValue
Double_t SmearValue(Double_t value)
Definition: CbmAnaConversionReco.cxx:1476
CbmAnaConversionReco::fhUsedMomenta_errorY_stsMomVec
TH1D * fhUsedMomenta_errorY_stsMomVec
Definition: CbmAnaConversionReco.h:195
CbmRichHit.h
CbmAnaConversionReco::fhEPEM_openingAngle_vs_pt_gg_reco
TH2D * fhEPEM_openingAngle_vs_pt_gg_reco
Definition: CbmAnaConversionReco.h:159
CbmAnaConversionReco::fhEPEM_identifiedLepton_selectionNN
TH1D * fhEPEM_identifiedLepton_selectionNN
Definition: CbmAnaConversionReco.h:241
CbmAnaConversionReco::fhEPEM_InDetector_invmass_gg_mc
TH1D * fhEPEM_InDetector_invmass_gg_mc
Definition: CbmAnaConversionReco.h:216
CbmAnaConversionReco::CalculateOpeningAngleBetweenGammasReco
Double_t CalculateOpeningAngleBetweenGammasReco(TVector3 electron1, TVector3 electron2, TVector3 electron3, TVector3 electron4)
Definition: CbmAnaConversionReco.cxx:3171
CbmAnaConversionReco::fhUsedMomenta_errorX_refitted
TH1D * fhUsedMomenta_errorX_refitted
Definition: CbmAnaConversionReco.h:199
CbmAnaConversionReco::SetTracklistMC
void SetTracklistMC(std::vector< CbmMCTrack * > MCTracklist)
Definition: CbmAnaConversionReco.cxx:1101
CbmAnaConversionReco::fhUsedMomenta_error_stsMomVec
TH1D * fhUsedMomenta_error_stsMomVec
Definition: CbmAnaConversionReco.h:191
CbmAnaConversionReco::fhEPEM_openingAngle_gee_mc_dalitz
TH1D * fhEPEM_openingAngle_gee_mc_dalitz
Definition: CbmAnaConversionReco.h:155
CbmAnaConversionReco::fhUsedMomenta_vsY_stsMomVec
TH2D * fhUsedMomenta_vsY_stsMomVec
Definition: CbmAnaConversionReco.h:196
CbmAnaConversionReco::Invmass_4particlesRECO
Double_t Invmass_4particlesRECO(const TVector3 part1, const TVector3 part2, const TVector3 part3, const TVector3 part4)
Definition: CbmAnaConversionReco.cxx:1494
CbmAnaConversionReco::fHistoList_gg
std::vector< TH1 * > fHistoList_gg
Definition: CbmAnaConversionReco.h:113
CbmAnaConversionReco::fhEPEM_InDetector_invmass_gg_refitted
TH1D * fhEPEM_InDetector_invmass_gg_refitted
Definition: CbmAnaConversionReco.h:217
CbmAnaConversionReco::fhPi0_startvertex
TH1D * fhPi0_startvertex
Definition: CbmAnaConversionReco.h:207
CbmAnaConversionReco::fhEPEM_pi0_ANNvalues_noCuts
TH1D * fhEPEM_pi0_ANNvalues_noCuts
Definition: CbmAnaConversionReco.h:175
CbmAnaConversionReco::CalculateOpeningAngleMC
Double_t CalculateOpeningAngleMC(CbmMCTrack *mctrack1, CbmMCTrack *mctrack2)
Definition: CbmAnaConversionReco.cxx:3114
CbmRichRing::GetCenterX
Float_t GetCenterX() const
Definition: CbmRichRing.h:80
CbmLmvmKinematicParams
Definition: CbmLmvmKinematicParams.h:17
CbmLmvmKinematicParams::fPt
Double_t fPt
Definition: CbmLmvmKinematicParams.h:20
CbmAnaConversionReco::fhUsedMomenta_vsX_stsMomVec
TH2D * fhUsedMomenta_vsX_stsMomVec
Definition: CbmAnaConversionReco.h:194
CbmRichHit
Definition: CbmRichHit.h:19
CbmAnaConversionReco::fhInvMass_EPEM_stsMomVec
TH1D * fhInvMass_EPEM_stsMomVec
Definition: CbmAnaConversionReco.h:184
CbmAnaConversionReco::fRecoRefittedMomentum
std::vector< TVector3 > fRecoRefittedMomentum
Definition: CbmAnaConversionReco.h:104
CbmLmvmKinematicParams::fAngle
Double_t fAngle
Definition: CbmLmvmKinematicParams.h:23
CbmAnaConversionReco::fhPi0_startvertexElectrons_gee
TH1D * fhPi0_startvertexElectrons_gee
Definition: CbmAnaConversionReco.h:210
CbmAnaConversionReco::fhEPEM_identifiedLepton_nofRingHits
TH1D * fhEPEM_identifiedLepton_nofRingHits
Definition: CbmAnaConversionReco.h:237
CbmAnaConversionReco::fhUsedMomenta_vsZ_stsMomVec
TH2D * fhUsedMomenta_vsZ_stsMomVec
Definition: CbmAnaConversionReco.h:198