Logo ROOT   6.13/01
Reference Guide
TLorentzVector.cxx
Go to the documentation of this file.
1 // @(#)root/physics:$Id$
2 // Author: Pasha Murat , Peter Malzacher 12/02/99
3 // Oct 8 1999: changed Warning to Error and
4 // return fX in Double_t & operator()
5 // Oct 20 1999: dito in Double_t operator()
6 // Jan 25 2000: implemented as (fP,fE) instead of (fX,fY,fZ,fE)
7 
8 
9 /** \class TLorentzVector
10  \ingroup Physics
11 TLorentzVector is a general four-vector class, which can be used
12 either for the description of position and time (x,y,z,t) or momentum and
13 energy (px,py,pz,E).
14 
15 ### Declaration
16 TLorentzVector has been implemented as a set a TVector3 and a Double_t variable.
17 By default all components are initialized by zero.
18 
19 ~~~ {.cpp}
20  TLorentzVector v1; // initialized by (0., 0., 0., 0.)
21  TLorentzVector v2(1., 1., 1., 1.);
22  TLorentzVector v3(v1);
23  TLorentzVector v4(TVector3(1., 2., 3.),4.);
24 ~~~
25 
26 For backward compatibility there are two constructors from an Double_t
27 and Float_t C array.
28 
29 
30 ### Access to the components
31 There are two sets of access functions to the components of a LorentzVector:
32 X(), Y(), Z(), T() and Px(),
33 Py(), Pz() and E(). Both sets return the same values
34 but the first set is more relevant for use where TLorentzVector
35 describes a combination of position and time and the second set is more
36 relevant where TLorentzVector describes momentum and energy:
37 
38 ~~~ {.cpp}
39  Double_t xx =v.X();
40  ...
41  Double_t tt = v.T();
42 
43  Double_t px = v.Px();
44  ...
45  Double_t ee = v.E();
46 ~~~
47 
48 The components of TLorentzVector can also accessed by index:
49 
50 ~~~ {.cpp}
51  xx = v(0); or xx = v[0];
52  yy = v(1); yy = v[1];
53  zz = v(2); zz = v[2];
54  tt = v(3); tt = v[3];
55 ~~~
56 
57 You can use the Vect() member function to get the vector component
58 of TLorentzVector:
59 
60 ~~~ {.cpp}
61  TVector3 p = v.Vect();
62 ~~~
63 
64 For setting components also two sets of member functions can be used:
65 
66 ~~~ {.cpp}
67  v.SetX(1.); or v.SetPx(1.);
68  ... ...
69  v.SetT(1.); v.SetE(1.);
70 ~~~
71 
72 To set more the one component by one call you can use the SetVect()
73 function for the TVector3 part or SetXYZT(), SetPxPyPzE(). For convenience there is
74 
75 also a SetXYZM():
76 
77 ~~~ {.cpp}
78  v.SetVect(TVector3(1,2,3));
79  v.SetXYZT(x,y,z,t);
80  v.SetPxPyPzE(px,py,pz,e);
81  v.SetXYZM(x,y,z,m); // -> v=(x,y,z,e=Sqrt(x*x+y*y+z*z+m*m))
82 ~~~
83 
84 ### Vector components in non-cartesian coordinate systems
85 There are a couple of member functions to get and set the TVector3
86 part of the parameters in
87 spherical coordinate systems:
88 
89 ~~~ {.cpp}
90  Double_t m, theta, cost, phi, pp, pp2, ppv2, pp2v2;
91  m = v.Rho();
92  t = v.Theta();
93  cost = v.CosTheta();
94  phi = v.Phi();
95 
96  v.SetRho(10.);
97  v.SetTheta(TMath::Pi()*.3);
98  v.SetPhi(TMath::Pi());
99 ~~~
100 
101 or get information about the r-coordinate in cylindrical systems:
102 
103 ~~~ {.cpp}
104  Double_t pp, pp2, ppv2, pp2v2;
105  pp = v.Perp(); // get transvers component
106  pp2 = v.Perp2(); // get transverse component squared
107  ppv2 = v.Perp(v1); // get transvers component with
108  // respect to another vector
109  pp2v2 = v.Perp(v1);
110 ~~~
111 
112 for convenience there are two more set functions SetPtEtaPhiE(pt,eta,phi,e);
113 and SetPtEtaPhiM(pt,eta,phi,m);
114 
115 ### Arithmetic and comparison operators
116 The TLorentzVector class provides operators to add, subtract or
117 compare four-vectors:
118 
119 ~~~ {.cpp}
120  v3 = -v1;
121  v1 = v2+v3;
122  v1+= v3;
123  v1 = v2 + v3;
124  v1-= v3;
125 
126  if (v1 == v2) {...}
127  if(v1 != v3) {...}
128 ~~~
129 
130 ### Magnitude/Invariant mass, beta, gamma, scalar product
131 The scalar product of two four-vectors is calculated with the (-,-,-,+)
132 metric,
133 
134  i.e. `s = v1*v2 = t1*t2-x1*x2-y1*y2-z1*z2`
135 The magnitude squared mag2 of a four-vector is therefore:
136 
137 ~~~ {.cpp}
138  mag2 = v*v = t*t-x*x-y*y-z*z
139 ~~~
140 It mag2 is negative mag = -Sqrt(-mag*mag). The member
141 functions are:
142 
143 ~~~ {.cpp}
144  Double_t s, s2;
145  s = v1.Dot(v2); // scalar product
146  s = v1*v2; // scalar product
147  s2 = v.Mag2(); or s2 = v.M2();
148  s = v.Mag(); s = v.M();
149 ~~~
150 
151 Since in case of momentum and energy the magnitude has the meaning of
152 invariant mass TLorentzVector provides the more meaningful aliases
153 M2() and M();
154 The member functions Beta() and Gamma() returns
155 beta and gamma = 1/Sqrt(1-beta*beta).
156 ### Lorentz boost
157 A boost in a general direction can be parameterised with three parameters
158 which can be taken as the components of a three vector b = (bx,by,bz).
159 With x = (x,y,z) and gamma = 1/Sqrt(1-beta*beta) (beta being the module of vector b),
160 an arbitrary active Lorentz boost transformation (from the rod frame
161 to the original frame) can be written as:
162 
163 ~~~ {.cpp}
164  x = x' + (gamma-1)/(beta*beta) * (b*x') * b + gamma * t' * b
165  t = gamma (t'+ b*x').
166 ~~~
167 
168 The member function Boost() performs a boost transformation
169 from the rod frame to the original frame. BoostVector() returns
170 a TVector3 of the spatial components divided by the time component:
171 
172 ~~~ {.cpp}
173  TVector3 b;
174  v.Boost(bx,by,bz);
175  v.Boost(b);
176  b = v.BoostVector(); // b=(x/t,y/t,z/t)
177 ~~~
178 
179 ### Rotations
180 There are four sets of functions to rotate the TVector3 component
181 of a TLorentzVector:
182 
183 #### rotation around axes
184 
185 ~~~ {.cpp}
186  v.RotateX(TMath::Pi()/2.);
187  v.RotateY(.5);
188  v.RotateZ(.99);
189 ~~~
190 
191 #### rotation around an arbitrary axis
192  v.Rotate(TMath::Pi()/4., v1); // rotation around v1
193 
194 #### transformation from rotated frame
195 
196 ~~~ {.cpp}
197  v.RotateUz(direction); // direction must be a unit TVector3
198 ~~~
199 
200 #### by TRotation (see TRotation)
201 
202 ~~~ {.cpp}
203  TRotation r;
204  v.Transform(r); or v *= r; // Attention v=M*v
205 ~~~
206 
207 ### Misc
208 
209 #### Angle between two vectors
210 
211 ~~~ {.cpp}
212  Double_t a = v1.Angle(v2.Vect()); // get angle between v1 and v2
213 ~~~
214 
215 #### Light-cone components
216 Member functions Plus() and Minus() return the positive
217 and negative light-cone components:
218 
219 ~~~ {.cpp}
220  Double_t pcone = v.Plus();
221  Double_t mcone = v.Minus();
222 ~~~
223 
224 CAVEAT: The values returned are T{+,-}Z. It is known that some authors
225 find it easier to define these components as (T{+,-}Z)/sqrt(2). Thus
226 check what definition is used in the physics you're working in and adapt
227 your code accordingly.
228 
229 #### Transformation by TLorentzRotation
230 A general Lorentz transformation see class TLorentzRotation can
231 be used by the Transform() member function, the *= or
232 * operator of the TLorentzRotation class:
233 
234 ~~~ {.cpp}
235  TLorentzRotation l;
236  v.Transform(l);
237  v = l*v; or v *= l; // Attention v = l*v
238 ~~~
239 */
240 
241 #include "TLorentzVector.h"
242 
243 #include "TBuffer.h"
244 #include "TString.h"
245 #include "TLorentzRotation.h"
246 
247 ClassImp(TLorentzVector);
248 
249 
250 void TLorentzVector::Boost(Double_t bx, Double_t by, Double_t bz)
251 {
252  //Boost this Lorentz vector
253  Double_t b2 = bx*bx + by*by + bz*bz;
254  Double_t gamma = 1.0 / TMath::Sqrt(1.0 - b2);
255  Double_t bp = bx*X() + by*Y() + bz*Z();
256  Double_t gamma2 = b2 > 0 ? (gamma - 1.0)/b2 : 0.0;
257 
258  SetX(X() + gamma2*bp*bx + gamma*bx*T());
259  SetY(Y() + gamma2*bp*by + gamma*by*T());
260  SetZ(Z() + gamma2*bp*bz + gamma*bz*T());
261  SetT(gamma*(T() + bp));
262 }
263 
264 Double_t TLorentzVector::Rapidity() const
265 {
266  //return rapidity
267  return 0.5*log( (E()+Pz()) / (E()-Pz()) );
268 }
269 
271 {
272  //multiply this Lorentzvector by m
273  return *this = m.VectorMultiplication(*this);
274 }
275 
277 {
278  //Transform this Lorentzvector
279  return *this = m.VectorMultiplication(*this);
280 }
281 
282 void TLorentzVector::Streamer(TBuffer &R__b)
283 {
284  // Stream an object of class TLorentzVector.
285  Double_t x, y, z;
286  UInt_t R__s, R__c;
287  if (R__b.IsReading()) {
288  Version_t R__v = R__b.ReadVersion(&R__s, &R__c);
289  if (R__v > 3) {
290  R__b.ReadClassBuffer(TLorentzVector::Class(), this, R__v, R__s, R__c);
291  return;
292  }
293  //====process old versions before automatic schema evolution
294  if (R__v != 2) TObject::Streamer(R__b);
295  R__b >> x;
296  R__b >> y;
297  R__b >> z;
298  fP.SetXYZ(x,y,z);
299  R__b >> fE;
300  R__b.CheckByteCount(R__s, R__c, TLorentzVector::IsA());
301  } else {
302  R__b.WriteClassBuffer(TLorentzVector::Class(),this);
303  }
304 }
305 
306 
307 ////////////////////////////////////////////////////////////////////////////////
308 /// Print the TLorentz vector components as (x,y,z,t) and (P,eta,phi,E) representations
309 
310 void TLorentzVector::Print(Option_t *) const
311 {
312  Printf("(x,y,z,t)=(%f,%f,%f,%f) (P,eta,phi,E)=(%f,%f,%f,%f)",
313  fP.x(),fP.y(),fP.z(),fE,
314  P(),Eta(),Phi(),fE);
315 }
The TLorentzRotation class describes Lorentz transformations including Lorentz boosts and rotations (...
void Boost(Double_t, Double_t, Double_t)
Double_t Z() const
virtual void Print(Option_t *option="") const
Print the TLorentz vector components as (x,y,z,t) and (P,eta,phi,E) representations.
Double_t y() const
Definition: TVector3.h:214
you should not use this method at all Int_t y
Definition: TRolke.cxx:630
TLorentzVector & operator*=(Double_t a)
Double_t Pz() const
void SetXYZ(Double_t x, Double_t y, Double_t z)
Definition: TVector3.h:227
you should not use this method at all Int_t Int_t Double_t Double_t Double_t Int_t Double_t Double_t Double_t Double_t Int_t m
Definition: TRolke.cxx:637
double gamma(double x)
TLorentzVector is a general four-vector class, which can be used either for the description of positi...
Double_t Eta() const
void SetX(Double_t a)
Double_t P() const
* x
Deprecated and error prone model selection interface.
Definition: TRolke.cxx:630
Double_t z() const
Definition: TVector3.h:215
void SetT(Double_t a)
Double_t Y() const
void SetY(Double_t a)
Double_t Phi() const
Double_t T() const
void SetZ(Double_t a)
you should not use this method at all Int_t Int_t z
Definition: TRolke.cxx:630
TLorentzVector VectorMultiplication(const TLorentzVector &) const
Double_t x() const
Definition: TVector3.h:213
Double_t X() const
TLorentzVector & Transform(const TRotation &)
Double_t E() const
Double_t Sqrt(Double_t x)
Definition: TMath.h:690
Double_t Rapidity() const
double log(double)