CbmRoot
CbmAnaConversionKinematicParams.h
Go to the documentation of this file.
1 
7 #ifndef CBM_ANA_CONVERSION_KINEMATIC_PARAMS
8 #define CBM_ANA_CONVERSION_KINEMATIC_PARAMS
9 
10 #include "CbmMCTrack.h"
11 #include "TLorentzVector.h"
12 #include "TMath.h"
13 
14 #define M2E 2.6112004954086e-7
15 
17 public:
18  Double_t fMomentumMag = 0; // Absolute value of momentum
19  Double_t fPt = 0; // Transverse momentum
20  Double_t fRapidity = 0; // Rapidity
21  Double_t fMinv = 0; // Invariant mass
22  Double_t fAngle = 0; // Opening angle
23  Double_t fRap2 = 0;
24  Double_t fPseudoRap2 = 0;
25 
26  /*
27  * Calculate kinematic parameters for MC tracks.
28  */
31  const CbmMCTrack* mctrack2,
32  const CbmMCTrack* mctrack3,
33  const CbmMCTrack* mctrack4) {
35 
36  TLorentzVector lorVec1;
37  mctrack1->Get4Momentum(lorVec1);
38  TVector3 part1 = lorVec1.Vect();
39  Double_t energy1 = lorVec1.Energy();
40 
41  TLorentzVector lorVec2;
42  mctrack2->Get4Momentum(lorVec2);
43  TVector3 part2 = lorVec2.Vect();
44  Double_t energy2 = lorVec2.Energy();
45 
46  TLorentzVector lorVec3;
47  mctrack3->Get4Momentum(lorVec3);
48  TVector3 part3 = lorVec3.Vect();
49  Double_t energy3 = lorVec3.Energy();
50 
51  TLorentzVector lorVec4;
52  mctrack4->Get4Momentum(lorVec4);
53  TVector3 part4 = lorVec4.Vect();
54  Double_t energy4 = lorVec4.Energy();
55 
56  TLorentzVector sum;
57  sum = lorVec1 + lorVec2 + lorVec3 + lorVec4;
58 
59  TVector3 momPair = part1 + part2 + part3 + part4;
60  Double_t energyPair = energy1 + energy2 + energy3 + energy4;
61  Double_t pzPair = momPair.Pz();
62  Double_t yPair =
63  0.5 * TMath::Log((energyPair + pzPair) / (energyPair - pzPair));
64 
65  Double_t invmass = sum.Mag();
66 
67  Double_t perp = sum.Perp();
68  //Double_t pt = TMath::Sqrt(sum.X() * sum.X() + sum.Y() * sum.Y() );
69 
70 
71  params.fMomentumMag = momPair.Mag();
72  params.fPt = perp;
73  params.fRapidity = yPair;
74  params.fMinv = invmass;
75  params.fRap2 = sum.Rapidity();
76  params.fPseudoRap2 = sum.PseudoRapidity();
77  return params;
78  }
79 
80  /*
81  * Calculate kinematic parameters for reconstructed momenta
82  */
84  KinematicParams_4particles_Reco(const TVector3 part1,
85  const TVector3 part2,
86  const TVector3 part3,
87  const TVector3 part4) {
89 
90  Double_t energy1 = TMath::Sqrt(part1.Mag2() + M2E);
91  TLorentzVector lorVec1(part1, energy1);
92 
93  Double_t energy2 = TMath::Sqrt(part2.Mag2() + M2E);
94  TLorentzVector lorVec2(part2, energy2);
95 
96  Double_t energy3 = TMath::Sqrt(part3.Mag2() + M2E);
97  TLorentzVector lorVec3(part3, energy3);
98 
99  Double_t energy4 = TMath::Sqrt(part4.Mag2() + M2E);
100  TLorentzVector lorVec4(part4, energy4);
101 
102  TLorentzVector sum;
103  sum = lorVec1 + lorVec2 + lorVec3 + lorVec4;
104 
105  TVector3 momPair = part1 + part2 + part3 + part4;
106  Double_t energyPair = energy1 + energy2 + energy3 + energy4;
107  Double_t pzPair = momPair.Pz();
108  Double_t yPair =
109  0.5 * TMath::Log((energyPair + pzPair) / (energyPair - pzPair));
110 
111  Double_t invmass = sum.Mag();
112 
113  Double_t perp = sum.Perp();
114  //Double_t pt = TMath::Sqrt(sum.X() * sum.X() + sum.Y() * sum.Y() );
115 
116 
117  params.fMomentumMag = momPair.Mag();
118  params.fPt = perp;
119  params.fRapidity = yPair;
120  params.fMinv = invmass;
121  params.fRap2 = sum.Rapidity();
122  params.fPseudoRap2 = sum.PseudoRapidity();
123  return params;
124  }
125 
126 
128  KinematicParams_2particles_Reco(const TVector3 electron1,
129  const TVector3 electron2) {
131 
132  Double_t energyP = TMath::Sqrt(electron1.Mag2() + M2E);
133  TLorentzVector lorVecP(electron1, energyP);
134 
135  Double_t energyM = TMath::Sqrt(electron2.Mag2() + M2E);
136  TLorentzVector lorVecM(electron2, energyM);
137 
138  TVector3 momPair = electron1 + electron2;
139  Double_t energyPair = energyP + energyM;
140  Double_t ptPair = momPair.Perp();
141  Double_t pzPair = momPair.Pz();
142  Double_t yPair =
143  0.5 * TMath::Log((energyPair + pzPair) / (energyPair - pzPair));
144  Double_t anglePair = lorVecM.Angle(lorVecP.Vect());
145  Double_t theta = 180. * anglePair / TMath::Pi();
146  Double_t minv = 2. * TMath::Sin(anglePair / 2.)
147  * TMath::Sqrt(electron1.Mag() * electron2.Mag());
148 
149  params.fMomentumMag = momPair.Mag();
150  params.fPt = ptPair;
151  params.fRapidity = yPair;
152  params.fMinv = minv;
153  params.fAngle = theta;
154  return params;
155  }
156 };
157 
158 #endif
CbmAnaConversionKinematicParams::fPt
Double_t fPt
Definition: CbmAnaConversionKinematicParams.h:19
CbmAnaConversionKinematicParams::fAngle
Double_t fAngle
Definition: CbmAnaConversionKinematicParams.h:22
CbmAnaConversionKinematicParams::fRap2
Double_t fRap2
Definition: CbmAnaConversionKinematicParams.h:23
CbmAnaConversionKinematicParams::KinematicParams_2particles_Reco
static CbmAnaConversionKinematicParams KinematicParams_2particles_Reco(const TVector3 electron1, const TVector3 electron2)
Definition: CbmAnaConversionKinematicParams.h:128
CbmAnaConversionKinematicParams::fMinv
Double_t fMinv
Definition: CbmAnaConversionKinematicParams.h:21
CbmAnaConversionKinematicParams::fPseudoRap2
Double_t fPseudoRap2
Definition: CbmAnaConversionKinematicParams.h:24
CbmAnaConversionKinematicParams::KinematicParams_4particles_Reco
static CbmAnaConversionKinematicParams KinematicParams_4particles_Reco(const TVector3 part1, const TVector3 part2, const TVector3 part3, const TVector3 part4)
Definition: CbmAnaConversionKinematicParams.h:84
CbmAnaConversionKinematicParams::fMomentumMag
Double_t fMomentumMag
Definition: CbmAnaConversionKinematicParams.h:18
CbmAnaConversionKinematicParams::KinematicParams_4particles_MC
static CbmAnaConversionKinematicParams KinematicParams_4particles_MC(const CbmMCTrack *mctrack1, const CbmMCTrack *mctrack2, const CbmMCTrack *mctrack3, const CbmMCTrack *mctrack4)
Definition: CbmAnaConversionKinematicParams.h:30
xMath::Pi
double Pi()
Definition: xMath.h:5
CbmAnaConversionKinematicParams
Definition: CbmAnaConversionKinematicParams.h:16
CbmMCTrack.h
CbmMCTrack
Definition: CbmMCTrack.h:34
CbmMCTrack::Get4Momentum
void Get4Momentum(TLorentzVector &momentum) const
Definition: CbmMCTrack.h:177
CbmAnaConversionKinematicParams::fRapidity
Double_t fRapidity
Definition: CbmAnaConversionKinematicParams.h:20
M2E
#define M2E
Definition: CbmAnaConversionKinematicParams.h:14