CbmRoot
CbmKresFunctions.h
Go to the documentation of this file.
1 
9 #ifndef CBM_KRES_FUNCTIONS
10 #define CBM_KRES_FUNCTIONS
11 
12 #include "CbmGlobalTrack.h"
13 #include "CbmLmvmKinematicParams.h"
14 #include "CbmMCTrack.h"
15 #include "CbmStsKFTrackFitter.h"
16 #include "CbmVertex.h"
17 #include "TLorentzVector.h"
18 #include "TMath.h"
19 #include "TMatrixTSym.h"
20 #include "TVector3.h"
21 
22 #define M2E 2.6112004954086e-7
23 #define M2Pion 0.0194798351452
24 
26 public:
28  static TVector3
29  FitToVertex(CbmStsTrack* stsTrack, double x, double y, double z) {
30  CbmVertex* vtx = new CbmVertex();
31  TMatrixFSym* covMat = new TMatrixFSym(3);
32  vtx->SetVertex(x, y, z, 0, 0, 0, *covMat);
33  CbmStsKFTrackFitter fitter;
34  FairTrackParam neu_track;
35  fitter.FitToVertex(stsTrack, vtx, &neu_track);
36 
37  TVector3 momentum;
38 
39  neu_track.Momentum(momentum);
40 
41  delete vtx;
42  delete covMat;
43 
44  return momentum;
45  }
46 
48  static double
49  ChiToVertex(CbmStsTrack* stsTrack, double x, double y, double z) {
50  CbmVertex* vtx = new CbmVertex();
51  TMatrixFSym* covMat = new TMatrixFSym(3);
52  vtx->SetVertex(x, y, z, 0, 0, 0, *covMat);
53  CbmStsKFTrackFitter fitter;
54  FairTrackParam neu_track;
55  fitter.FitToVertex(stsTrack, vtx, &neu_track);
56 
57  double chi = fitter.GetChiToVertex(stsTrack, vtx);
58 
59  delete vtx;
60  delete covMat;
61 
62  return chi;
63  }
64 
65 
67  static TVector3 FitToVertexAndGetChi(CbmStsTrack* stsTrack,
68  double x,
69  double y,
70  double z,
71  double& chi) {
72  CbmVertex* vtx = new CbmVertex();
73  TMatrixFSym* covMat = new TMatrixFSym(3);
74  vtx->SetVertex(x, y, z, 0, 0, 0, *covMat);
75  CbmStsKFTrackFitter fitter;
76  FairTrackParam neu_track;
77  fitter.FitToVertex(stsTrack, vtx, &neu_track);
78 
79  chi = fitter.GetChiToVertex(stsTrack, vtx);
80 
81  TVector3 momentum;
82 
83  neu_track.Momentum(momentum);
84 
85  delete vtx;
86  delete covMat;
87 
88  return momentum;
89  }
90 
91  // calculation of invariant mass from 2 tracks using MCtrue momenta
92  static double Invmass_2particles_MC(const CbmMCTrack* mctrack1,
93  const CbmMCTrack* mctrack2) {
94  TLorentzVector lorVec1;
95  mctrack1->Get4Momentum(lorVec1);
96 
97  TLorentzVector lorVec2;
98  mctrack2->Get4Momentum(lorVec2);
99 
100  TLorentzVector sum;
101  sum = lorVec1 + lorVec2;
102 
103  return sum.Mag();
104  }
105 
106  // calculation of invariant mass from 2 tracks using reconstructed momenta
107  static double Invmass_2particles_RECO(const TVector3 part1,
108  const TVector3 part2) {
109  Double_t energy1 = TMath::Sqrt(part1.Mag2() + M2E);
110  TLorentzVector lorVec1(part1, energy1);
111 
112  Double_t energy2 = TMath::Sqrt(part2.Mag2() + M2E);
113  TLorentzVector lorVec2(part2, energy2);
114 
115  TLorentzVector sum;
116  sum = lorVec1 + lorVec2;
117 
118  return sum.Mag();
119  }
120 
121  // calculation of invariant mass from 4 tracks using MCtrue momenta
122  static double Invmass_4particles_MC(const CbmMCTrack* mctrack1,
123  const CbmMCTrack* mctrack2,
124  const CbmMCTrack* mctrack3,
125  const CbmMCTrack* mctrack4) {
126  TLorentzVector lorVec1;
127  mctrack1->Get4Momentum(lorVec1);
128 
129  TLorentzVector lorVec2;
130  mctrack2->Get4Momentum(lorVec2);
131 
132  TLorentzVector lorVec3;
133  mctrack3->Get4Momentum(lorVec3);
134 
135  TLorentzVector lorVec4;
136  mctrack4->Get4Momentum(lorVec4);
137 
138  TLorentzVector sum;
139  sum = lorVec1 + lorVec2 + lorVec3 + lorVec4;
140 
141  return sum.Mag();
142  }
143 
144  // calculation of invariant mass from 4 tracks using reconstructed momenta
145  static double Invmass_4particles_RECO(const TVector3 part1,
146  const TVector3 part2,
147  const TVector3 part3,
148  const TVector3 part4) {
149  Double_t energy1 = TMath::Sqrt(part1.Mag2() + M2E);
150  TLorentzVector lorVec1(part1, energy1);
151 
152  Double_t energy2 = TMath::Sqrt(part2.Mag2() + M2E);
153  TLorentzVector lorVec2(part2, energy2);
154 
155  Double_t energy3 = TMath::Sqrt(part3.Mag2() + M2E);
156  TLorentzVector lorVec3(part3, energy3);
157 
158  Double_t energy4 = TMath::Sqrt(part4.Mag2() + M2E);
159  TLorentzVector lorVec4(part4, energy4);
160 
161  TLorentzVector sum;
162  sum = lorVec1 + lorVec2 + lorVec3 + lorVec4;
163 
164  return sum.Mag();
165  }
166 
167  // gives amount of daughter particles belonging to current particle -> MCtrue info
168  static int NofDaughters(int motherId, vector<CbmMCTrack*> MC) {
169  int nofDaughters = 0;
170  for (size_t i = 0; i < MC.size(); i++) {
171  int motherId_temp = MC.at(i)->GetMotherId();
172  if (motherId == motherId_temp) nofDaughters++;
173  }
174  return nofDaughters;
175  }
176 
177 
178  // calculation of many interesting for analysis physics parameters of 2 tracks: Invariant mass, opening angle, rapidity, pt,
180  CalculateKinematicParamsReco(const TVector3 electron1,
181  const TVector3 electron2) {
182  CbmLmvmKinematicParams params;
183 
184  Double_t energyP = TMath::Sqrt(electron1.Mag2() + M2E);
185  TLorentzVector lorVecP(electron1, energyP);
186 
187  Double_t energyM = TMath::Sqrt(electron2.Mag2() + M2E);
188  TLorentzVector lorVecM(electron2, energyM);
189 
190  TVector3 momPair = electron1 + electron2;
191  Double_t energyPair = energyP + energyM;
192  Double_t ptPair = momPair.Perp();
193  Double_t pzPair = momPair.Pz();
194  Double_t yPair =
195  0.5 * TMath::Log((energyPair + pzPair) / (energyPair - pzPair));
196  Double_t anglePair = lorVecM.Angle(lorVecP.Vect());
197  Double_t theta = 180. * anglePair / TMath::Pi();
198  Double_t minv = 2. * TMath::Sin(anglePair / 2.)
199  * TMath::Sqrt(electron1.Mag() * electron2.Mag());
200 
201  params.fMomentumMag = momPair.Mag();
202  params.fPt = ptPair;
203  params.fRapidity = yPair;
204  params.fMinv = minv;
205  params.fAngle = theta;
206  return params;
207  }
208 
209 
210  // calculation of many interesting for analysis physics parameters of 4 tracks: Invariant mass, opening angle, rapidity, pt,
213  const TVector3 part2,
214  const TVector3 part3,
215  const TVector3 part4) {
216  CbmLmvmKinematicParams params;
217 
218  Double_t energy1 = TMath::Sqrt(part1.Mag2() + M2E);
219  TLorentzVector lorVec1(part1, energy1);
220 
221  Double_t energy2 = TMath::Sqrt(part2.Mag2() + M2E);
222  TLorentzVector lorVec2(part2, energy2);
223 
224  Double_t energy3 = TMath::Sqrt(part3.Mag2() + M2E);
225  TLorentzVector lorVec3(part3, energy3);
226 
227  Double_t energy4 = TMath::Sqrt(part4.Mag2() + M2E);
228  TLorentzVector lorVec4(part4, energy4);
229 
230  TLorentzVector sum;
231  sum = lorVec1 + lorVec2 + lorVec3 + lorVec4;
232 
233  TVector3 momPair = part1 + part2 + part3 + part4;
234  Double_t energyPair = energy1 + energy2 + energy3 + energy4;
235  Double_t ptPair = momPair.Perp();
236  Double_t pzPair = momPair.Pz();
237  Double_t yPair =
238  0.5 * TMath::Log((energyPair + pzPair) / (energyPair - pzPair));
239  Double_t minv = sum.Mag();
240 
241  params.fMomentumMag = momPair.Mag();
242  params.fPt = ptPair;
243  params.fRapidity = yPair;
244  params.fMinv = minv;
245  params.fAngle = 0;
246  return params;
247  }
248 
249 
250  // calculation of opening angle from 2 tracks using reconstructed momenta
251  static Double_t CalculateOpeningAngle_Reco(TVector3 electron1,
252  TVector3 electron2) {
253  Double_t energyP = TMath::Sqrt(electron1.Mag2() + M2E);
254  TLorentzVector lorVecP(electron1, energyP);
255 
256  Double_t energyM = TMath::Sqrt(electron2.Mag2() + M2E);
257  TLorentzVector lorVecM(electron2, energyM);
258 
259  Double_t anglePair = lorVecM.Angle(lorVecP.Vect());
260  Double_t theta = 180. * anglePair / TMath::Pi();
261 
262  return theta;
263  }
264 
265 
266  // calculation of opening angle from 2 tracks using MCtrue momenta
267  static Double_t CalculateOpeningAngle_MC(CbmMCTrack* mctrack1,
268  CbmMCTrack* mctrack2) {
269  TVector3 electron1;
270  mctrack1->GetMomentum(electron1);
271  Double_t energyP = TMath::Sqrt(electron1.Mag2() + M2E);
272  TLorentzVector lorVecP(electron1, energyP);
273 
274  TVector3 electron2;
275  mctrack2->GetMomentum(electron2);
276  Double_t energyM = TMath::Sqrt(electron2.Mag2() + M2E);
277  TLorentzVector lorVecM(electron2, energyM);
278 
279  Double_t anglePair = lorVecM.Angle(lorVecP.Vect());
280  Double_t theta = 180. * anglePair / TMath::Pi();
281 
282  return theta;
283  }
284 
285 
286  // calculation of opening angle from 2 photons using reconstructed momenta
287  static Double_t CalculateOpeningAngleBetweenGammas_Reco(TVector3 electron1,
288  TVector3 electron2,
289  TVector3 electron3,
290  TVector3 electron4) {
291  Double_t energy1 = TMath::Sqrt(electron1.Mag2() + M2E);
292  TLorentzVector lorVec1(electron1, energy1);
293 
294  Double_t energy2 = TMath::Sqrt(electron2.Mag2() + M2E);
295  TLorentzVector lorVec2(electron2, energy2);
296 
297  Double_t energy3 = TMath::Sqrt(electron3.Mag2() + M2E);
298  TLorentzVector lorVec3(electron3, energy3);
299 
300  Double_t energy4 = TMath::Sqrt(electron4.Mag2() + M2E);
301  TLorentzVector lorVec4(electron4, energy4);
302 
303  TLorentzVector lorPhoton1 = lorVec1 + lorVec2;
304  TLorentzVector lorPhoton2 = lorVec3 + lorVec4;
305 
306  Double_t angleBetweenPhotons = lorPhoton1.Angle(lorPhoton2.Vect());
307  Double_t theta = 180. * angleBetweenPhotons / TMath::Pi();
308 
309  return theta;
310  }
311 
312 
313  // calculation of invariant mass of 4 particles using reconstructed value of momenta
314  static double Invmass_2el_2pions_RECO(const TVector3 part1El,
315  const TVector3 part2El,
316  const TVector3 part3Pion,
317  const TVector3 part4Pion) {
318  Double_t energy1 = TMath::Sqrt(part1El.Mag2() + M2E);
319  TLorentzVector lorVec1(part1El, energy1);
320 
321  Double_t energy2 = TMath::Sqrt(part2El.Mag2() + M2E);
322  TLorentzVector lorVec2(part2El, energy2);
323 
324  Double_t energy3 = TMath::Sqrt(part3Pion.Mag2() + M2Pion);
325  TLorentzVector lorVec3(part3Pion, energy3);
326 
327  Double_t energy4 = TMath::Sqrt(part4Pion.Mag2() + M2Pion);
328  TLorentzVector lorVec4(part4Pion, energy4);
329 
330  TLorentzVector sum;
331  sum = lorVec1 + lorVec2 + lorVec3 + lorVec4;
332 
333  return sum.Mag();
334  }
335 
336 
337  // calculation of invariant mass of 6 particles using true value of momenta
338  static double Invmass_6particles_MC(const CbmMCTrack* mctrack1,
339  const CbmMCTrack* mctrack2,
340  const CbmMCTrack* mctrack3,
341  const CbmMCTrack* mctrack4,
342  const CbmMCTrack* mctrack5,
343  const CbmMCTrack* mctrack6) {
344  TLorentzVector lorVec1;
345  mctrack1->Get4Momentum(lorVec1);
346 
347  TLorentzVector lorVec2;
348  mctrack2->Get4Momentum(lorVec2);
349 
350  TLorentzVector lorVec3;
351  mctrack3->Get4Momentum(lorVec3);
352 
353  TLorentzVector lorVec4;
354  mctrack4->Get4Momentum(lorVec4);
355 
356  TLorentzVector lorVec5;
357  mctrack5->Get4Momentum(lorVec5);
358 
359  TLorentzVector lorVec6;
360  mctrack6->Get4Momentum(lorVec6);
361 
362  TLorentzVector sum;
363  sum = lorVec1 + lorVec2 + lorVec3 + lorVec4 + lorVec5 + lorVec6;
364 
365  return sum.Mag();
366  }
367 
368  // calculation of invariant mass of 6 particles using reconstructed value of momenta
369  static double Invmass_4el_2pions_RECO(const TVector3 part1El,
370  const TVector3 part2El,
371  const TVector3 part3El,
372  const TVector3 part4El,
373  const TVector3 part5Pion,
374  const TVector3 part6Pion) {
375  Double_t energy1 = TMath::Sqrt(part1El.Mag2() + M2E);
376  TLorentzVector lorVec1(part1El, energy1);
377 
378  Double_t energy2 = TMath::Sqrt(part2El.Mag2() + M2E);
379  TLorentzVector lorVec2(part2El, energy2);
380 
381  Double_t energy3 = TMath::Sqrt(part3El.Mag2() + M2E);
382  TLorentzVector lorVec3(part3El, energy3);
383 
384  Double_t energy4 = TMath::Sqrt(part4El.Mag2() + M2E);
385  TLorentzVector lorVec4(part4El, energy4);
386 
387  Double_t energy5 = TMath::Sqrt(part5Pion.Mag2() + M2Pion);
388  TLorentzVector lorVec5(part5Pion, energy5);
389 
390  Double_t energy6 = TMath::Sqrt(part6Pion.Mag2() + M2Pion);
391  TLorentzVector lorVec6(part6Pion, energy6);
392 
393  TLorentzVector sum;
394  sum = lorVec1 + lorVec2 + lorVec3 + lorVec4 + lorVec5 + lorVec6;
395 
396  return sum.Mag();
397  }
398 
399 
401  static Double_t CalculateOpeningAngleBetweenPions_Reco(TVector3 electron1,
402  TVector3 electron2) {
403  Double_t energyP = TMath::Sqrt(electron1.Mag2() + M2Pion);
404  TLorentzVector lorVecP(electron1, energyP);
405 
406  Double_t energyM = TMath::Sqrt(electron2.Mag2() + M2Pion);
407  TLorentzVector lorVecM(electron2, energyM);
408 
409  Double_t anglePair = lorVecM.Angle(lorVecP.Vect());
410  Double_t theta = 180. * anglePair / TMath::Pi();
411 
412  return theta;
413  }
414 
415 
418  CbmMCTrack* mctrack2) {
419  TVector3 electron1;
420  mctrack1->GetMomentum(electron1);
421  Double_t energyP = TMath::Sqrt(electron1.Mag2() + M2Pion);
422  TLorentzVector lorVecP(electron1, energyP);
423 
424  TVector3 electron2;
425  mctrack2->GetMomentum(electron2);
426  Double_t energyM = TMath::Sqrt(electron2.Mag2() + M2Pion);
427  TLorentzVector lorVecM(electron2, energyM);
428 
429  Double_t anglePair = lorVecM.Angle(lorVecP.Vect());
430  Double_t theta = 180. * anglePair / TMath::Pi();
431 
432  return theta;
433  }
434 
435 
436  // calculation of many interesting for analysis physics parameters of 2 leptons and 2 pions: Invariant mass, opening angle, rapidity, pt,
439  const TVector3 part2,
440  const TVector3 part3,
441  const TVector3 part4) {
442  CbmLmvmKinematicParams params;
443 
444  Double_t energy1 = TMath::Sqrt(part1.Mag2() + M2E);
445  TLorentzVector lorVec1(part1, energy1);
446 
447  Double_t energy2 = TMath::Sqrt(part2.Mag2() + M2E);
448  TLorentzVector lorVec2(part2, energy2);
449 
450  Double_t energy3 = TMath::Sqrt(part3.Mag2() + M2Pion);
451  TLorentzVector lorVec3(part3, energy3);
452 
453  Double_t energy4 = TMath::Sqrt(part4.Mag2() + M2Pion);
454  TLorentzVector lorVec4(part4, energy4);
455 
456  TLorentzVector sum;
457  sum = lorVec1 + lorVec2 + lorVec3 + lorVec4;
458 
459  TVector3 momPair = part1 + part2 + part3 + part4;
460  Double_t energyPair = energy1 + energy2 + energy3 + energy4;
461  Double_t ptPair = momPair.Perp();
462  Double_t pzPair = momPair.Pz();
463  Double_t yPair =
464  0.5 * TMath::Log((energyPair + pzPair) / (energyPair - pzPair));
465  Double_t minv = sum.Mag();
466 
467  params.fMomentumMag = momPair.Mag();
468  params.fPt = ptPair;
469  params.fRapidity = yPair;
470  params.fMinv = minv;
471  params.fAngle = 0;
472  return params;
473  }
474 };
475 
476 #endif
CbmVertex::SetVertex
void SetVertex(Double_t x, Double_t y, Double_t z, Double_t chi2, Int_t ndf, Int_t nTracks, const TMatrixFSym &covMat)
Definition: CbmVertex.cxx:129
CbmKresFunctions
Definition: CbmKresFunctions.h:25
CbmMCTrack::GetMomentum
void GetMomentum(TVector3 &momentum) const
Definition: CbmMCTrack.h:172
CbmVertex.h
CbmKresFunctions::Invmass_6particles_MC
static double Invmass_6particles_MC(const CbmMCTrack *mctrack1, const CbmMCTrack *mctrack2, const CbmMCTrack *mctrack3, const CbmMCTrack *mctrack4, const CbmMCTrack *mctrack5, const CbmMCTrack *mctrack6)
Definition: CbmKresFunctions.h:338
CbmKresFunctions::CalculateOpeningAngleBetweenPions_Reco
static Double_t CalculateOpeningAngleBetweenPions_Reco(TVector3 electron1, TVector3 electron2)
calculate opening angle between two pions using reconstructed momenta
Definition: CbmKresFunctions.h:401
CbmKresFunctions::Invmass_2el_2pions_RECO
static double Invmass_2el_2pions_RECO(const TVector3 part1El, const TVector3 part2El, const TVector3 part3Pion, const TVector3 part4Pion)
Definition: CbmKresFunctions.h:314
CbmKresFunctions::CalculateKinematicParamsReco
static CbmLmvmKinematicParams CalculateKinematicParamsReco(const TVector3 electron1, const TVector3 electron2)
Definition: CbmKresFunctions.h:180
CbmStsKFTrackFitter::GetChiToVertex
Double_t GetChiToVertex(CbmStsTrack *track, CbmVertex *vtx=0)
Definition: CbmStsKFTrackFitter.cxx:164
i
int i
Definition: L1/vectors/P4_F32vec4.h:25
CbmKresFunctions::Invmass_4el_2pions_RECO
static double Invmass_4el_2pions_RECO(const TVector3 part1El, const TVector3 part2El, const TVector3 part3El, const TVector3 part4El, const TVector3 part5Pion, const TVector3 part6Pion)
Definition: CbmKresFunctions.h:369
CbmGlobalTrack.h
M2Pion
#define M2Pion
Definition: CbmKresFunctions.h:23
CbmLmvmKinematicParams::fMinv
Double_t fMinv
Definition: CbmLmvmKinematicParams.h:22
CbmKresFunctions::ChiToVertex
static double ChiToVertex(CbmStsTrack *stsTrack, double x, double y, double z)
Definition: CbmKresFunctions.h:49
CbmKresFunctions::NofDaughters
static int NofDaughters(int motherId, vector< CbmMCTrack * > MC)
Definition: CbmKresFunctions.h:168
CbmStsKFTrackFitter
Definition: CbmStsKFTrackFitter.h:14
CbmKresFunctions::CalculateOpeningAngleBetweenPions_MC
static Double_t CalculateOpeningAngleBetweenPions_MC(CbmMCTrack *mctrack1, CbmMCTrack *mctrack2)
calculate opening angle between two pions using MCtrue momenta
Definition: CbmKresFunctions.h:417
CbmVertex
Definition: CbmVertex.h:26
CbmKresFunctions::FitToVertex
static TVector3 FitToVertex(CbmStsTrack *stsTrack, double x, double y, double z)
Definition: CbmKresFunctions.h:29
CbmKresFunctions::CalculateKinematicParams_4particles
static CbmLmvmKinematicParams CalculateKinematicParams_4particles(const TVector3 part1, const TVector3 part2, const TVector3 part3, const TVector3 part4)
Definition: CbmKresFunctions.h:212
CbmLmvmKinematicParams::fRapidity
Double_t fRapidity
Definition: CbmLmvmKinematicParams.h:21
xMath::Pi
double Pi()
Definition: xMath.h:5
CbmKresFunctions::Invmass_2particles_RECO
static double Invmass_2particles_RECO(const TVector3 part1, const TVector3 part2)
Definition: CbmKresFunctions.h:107
CbmKresFunctions::CalculateOpeningAngle_MC
static Double_t CalculateOpeningAngle_MC(CbmMCTrack *mctrack1, CbmMCTrack *mctrack2)
Definition: CbmKresFunctions.h:267
CbmMCTrack.h
CbmKresFunctions::Invmass_4particles_MC
static double Invmass_4particles_MC(const CbmMCTrack *mctrack1, const CbmMCTrack *mctrack2, const CbmMCTrack *mctrack3, const CbmMCTrack *mctrack4)
Definition: CbmKresFunctions.h:122
CbmKresFunctions::FitToVertexAndGetChi
static TVector3 FitToVertexAndGetChi(CbmStsTrack *stsTrack, double x, double y, double z, double &chi)
Definition: CbmKresFunctions.h:67
CbmMCTrack
Definition: CbmMCTrack.h:34
x
Double_t x
Definition: CbmMvdSensorDigiToHitTask.cxx:68
CbmStsKFTrackFitter::FitToVertex
Double_t FitToVertex(CbmStsTrack *track, CbmVertex *vtx, FairTrackParam *v_track)
Definition: CbmStsKFTrackFitter.cxx:200
y
Double_t y
Definition: CbmMvdSensorDigiToHitTask.cxx:68
CbmLmvmKinematicParams::fMomentumMag
Double_t fMomentumMag
Definition: CbmLmvmKinematicParams.h:19
CbmMCTrack::Get4Momentum
void Get4Momentum(TLorentzVector &momentum) const
Definition: CbmMCTrack.h:177
CbmStsTrack
Definition: CbmStsTrack.h:37
CbmKresFunctions::CalculateKinematicParams_2el_2pions
static CbmLmvmKinematicParams CalculateKinematicParams_2el_2pions(const TVector3 part1, const TVector3 part2, const TVector3 part3, const TVector3 part4)
Definition: CbmKresFunctions.h:438
M2E
#define M2E
Definition: CbmKresFunctions.h:22
CbmKresFunctions::CalculateOpeningAngle_Reco
static Double_t CalculateOpeningAngle_Reco(TVector3 electron1, TVector3 electron2)
Definition: CbmKresFunctions.h:251
CbmKresFunctions::Invmass_2particles_MC
static double Invmass_2particles_MC(const CbmMCTrack *mctrack1, const CbmMCTrack *mctrack2)
Definition: CbmKresFunctions.h:92
CbmLmvmKinematicParams
Definition: CbmLmvmKinematicParams.h:17
CbmLmvmKinematicParams::fPt
Double_t fPt
Definition: CbmLmvmKinematicParams.h:20
CbmStsKFTrackFitter.h
CbmKresFunctions::CalculateOpeningAngleBetweenGammas_Reco
static Double_t CalculateOpeningAngleBetweenGammas_Reco(TVector3 electron1, TVector3 electron2, TVector3 electron3, TVector3 electron4)
Definition: CbmKresFunctions.h:287
CbmKresFunctions::Invmass_4particles_RECO
static double Invmass_4particles_RECO(const TVector3 part1, const TVector3 part2, const TVector3 part3, const TVector3 part4)
Definition: CbmKresFunctions.h:145
CbmLmvmKinematicParams.h
CbmLmvmKinematicParams::fAngle
Double_t fAngle
Definition: CbmLmvmKinematicParams.h:23