Logo ROOT   6.08/07
Reference Guide
TLorentzVector.h
Go to the documentation of this file.
1 // @(#)root/physics:$Id$
2 // Author: Pasha Murat , Peter Malzacher 12/02/99
3 
4 /*************************************************************************
5  * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *
6  * All rights reserved. *
7  * *
8  * For the licensing terms see $ROOTSYS/LICENSE. *
9  * For the list of contributors see $ROOTSYS/README/CREDITS. *
10  *************************************************************************/
11 
12 #ifndef ROOT_TLorentzVector
13 #define ROOT_TLorentzVector
14 
15 
16 //////////////////////////////////////////////////////////////////////////
17 // //
18 // TLorentzVector //
19 // //
20 // Place holder for real lorentz vector class. //
21 // //
22 //////////////////////////////////////////////////////////////////////////
23 
24 #ifndef ROOT_TMath
25 #include "TMath.h"
26 #endif
27 #ifndef ROOT_TVector3
28 #include "TVector3.h"
29 #endif
30 #ifndef ROOT_TRotation
31 #include "TRotation.h"
32 #endif
33 
34 
35 class TLorentzRotation;
36 
37 
38 class TLorentzVector : public TObject {
39 
40 private:
41 
42  TVector3 fP; // 3 vector component
43  Double_t fE; // time or energy of (x,y,z,t) or (px,py,pz,e)
44 
45 public:
46 
47  typedef Double_t Scalar; // to be able to use it with the ROOT::Math::VectorUtil functions
48 
49  enum { kX=0, kY=1, kZ=2, kT=3, kNUM_COORDINATES=4, kSIZE=kNUM_COORDINATES };
50  // Safe indexing of the coordinates when using with matrices, arrays, etc.
51 
53 
55  // Constructor giving the components x, y, z, t.
56 
57  TLorentzVector(const Double_t * carray);
58  TLorentzVector(const Float_t * carray);
59  // Constructor from an array, not checked!
60 
61  TLorentzVector(const TVector3 & vector3, Double_t t);
62  // Constructor giving a 3-Vector and a time component.
63 
64  TLorentzVector(const TLorentzVector & lorentzvector);
65  // Copy constructor.
66 
67  virtual ~TLorentzVector(){};
68  // Destructor
69 
70  // inline operator TVector3 () const;
71  // inline operator TVector3 & ();
72  // Conversion (cast) to TVector3.
73 
74  inline Double_t X() const;
75  inline Double_t Y() const;
76  inline Double_t Z() const;
77  inline Double_t T() const;
78  // Get position and time.
79 
80  inline void SetX(Double_t a);
81  inline void SetY(Double_t a);
82  inline void SetZ(Double_t a);
83  inline void SetT(Double_t a);
84  // Set position and time.
85 
86  inline Double_t Px() const;
87  inline Double_t Py() const;
88  inline Double_t Pz() const;
89  inline Double_t P() const;
90  inline Double_t E() const;
91  inline Double_t Energy() const;
92  // Get momentum and energy.
93 
94  inline void SetPx(Double_t a);
95  inline void SetPy(Double_t a);
96  inline void SetPz(Double_t a);
97  inline void SetE(Double_t a);
98  // Set momentum and energy.
99 
100  inline TVector3 Vect() const ;
101  // Get spatial component.
102 
103  inline void SetVect(const TVector3 & vect3);
104  // Set spatial component.
105 
106  inline Double_t Theta() const;
107  inline Double_t CosTheta() const;
108  inline Double_t Phi() const; //returns phi from -pi to pi
109  inline Double_t Rho() const;
110  // Get spatial vector components in spherical coordinate system.
111 
112  inline void SetTheta(Double_t theta);
113  inline void SetPhi(Double_t phi);
114  inline void SetRho(Double_t rho);
115  // Set spatial vector components in spherical coordinate system.
116 
117  inline void SetPxPyPzE(Double_t px, Double_t py, Double_t pz, Double_t e);
118  inline void SetXYZT(Double_t x, Double_t y, Double_t z, Double_t t);
119  inline void SetXYZM(Double_t x, Double_t y, Double_t z, Double_t m);
120  inline void SetPtEtaPhiM(Double_t pt, Double_t eta, Double_t phi, Double_t m);
121  inline void SetPtEtaPhiE(Double_t pt, Double_t eta, Double_t phi, Double_t e);
122  // Setters to provide the functionality (but a more meanigful name) of
123  // the previous version eg SetV4... PsetV4...
124 
125  inline void GetXYZT(Double_t *carray) const;
126  inline void GetXYZT(Float_t *carray) const;
127  // Getters into an arry
128  // no checking!
129 
130  Double_t operator () (int i) const;
131  inline Double_t operator [] (int i) const;
132  // Get components by index.
133 
134  Double_t & operator () (int i);
135  inline Double_t & operator [] (int i);
136  // Set components by index.
137 
138  inline TLorentzVector & operator = (const TLorentzVector &);
139  // Assignment.
140 
141  inline TLorentzVector operator + (const TLorentzVector &) const;
142  inline TLorentzVector & operator += (const TLorentzVector &);
143  // Additions.
144 
145  inline TLorentzVector operator - (const TLorentzVector &) const;
146  inline TLorentzVector & operator -= (const TLorentzVector &);
147  // Subtractions.
148 
149  inline TLorentzVector operator - () const;
150  // Unary minus.
151 
152  inline TLorentzVector operator * (Double_t a) const;
154  // Scaling with real numbers.
155 
156  inline Bool_t operator == (const TLorentzVector &) const;
157  inline Bool_t operator != (const TLorentzVector &) const;
158  // Comparisons.
159 
160  inline Double_t Perp2() const;
161  // Transverse component of the spatial vector squared.
162 
163  inline Double_t Pt() const;
164  inline Double_t Perp() const;
165  // Transverse component of the spatial vector (R in cylindrical system).
166 
167  inline void SetPerp(Double_t);
168  // Set the transverse component of the spatial vector.
169 
170  inline Double_t Perp2(const TVector3 & v) const;
171  // Transverse component of the spatial vector w.r.t. given axis squared.
172 
173  inline Double_t Pt(const TVector3 & v) const;
174  inline Double_t Perp(const TVector3 & v) const;
175  // Transverse component of the spatial vector w.r.t. given axis.
176 
177  inline Double_t Et2() const;
178  // Transverse energy squared.
179 
180  inline Double_t Et() const;
181  // Transverse energy.
182 
183  inline Double_t Et2(const TVector3 &) const;
184  // Transverse energy w.r.t. given axis squared.
185 
186  inline Double_t Et(const TVector3 &) const;
187  // Transverse energy w.r.t. given axis.
188 
189  inline Double_t DeltaPhi(const TLorentzVector &) const;
190  inline Double_t DeltaR(const TLorentzVector &) const;
191  inline Double_t DrEtaPhi(const TLorentzVector &) const;
192  inline TVector2 EtaPhiVector();
193 
194  inline Double_t Angle(const TVector3 & v) const;
195  // Angle wrt. another vector.
196 
197  inline Double_t Mag2() const;
198  inline Double_t M2() const;
199  // Invariant mass squared.
200 
201  inline Double_t Mag() const;
202  inline Double_t M() const;
203  // Invariant mass. If mag2() is negative then -sqrt(-mag2()) is returned.
204 
205  inline Double_t Mt2() const;
206  // Transverse mass squared.
207 
208  inline Double_t Mt() const;
209  // Transverse mass.
210 
211  inline Double_t Beta() const;
212  inline Double_t Gamma() const;
213 
214  inline Double_t Dot(const TLorentzVector &) const;
215  inline Double_t operator * (const TLorentzVector &) const;
216  // Scalar product.
217 
218  inline void SetVectMag(const TVector3 & spatial, Double_t magnitude);
219  inline void SetVectM(const TVector3 & spatial, Double_t mass);
220  // Copy spatial coordinates, and set energy = sqrt(mass^2 + spatial^2)
221 
222  inline Double_t Plus() const;
223  inline Double_t Minus() const;
224  // Returns t +/- z.
225  // Related to the positive/negative light-cone component,
226  // which some define this way and others define as (t +/- z)/sqrt(2)
227 
228  inline TVector3 BoostVector() const ;
229  // Returns the spatial components divided by the time component.
230 
232  inline void Boost(const TVector3 &);
233  // Lorentz boost.
234 
235  Double_t Rapidity() const;
236  // Returns the rapidity, i.e. 0.5*ln((E+pz)/(E-pz))
237 
238  inline Double_t Eta() const;
239  inline Double_t PseudoRapidity() const;
240  // Returns the pseudo-rapidity, i.e. -ln(tan(theta/2))
241 
242  inline void RotateX(Double_t angle);
243  // Rotate the spatial component around the x-axis.
244 
245  inline void RotateY(Double_t angle);
246  // Rotate the spatial component around the y-axis.
247 
248  inline void RotateZ(Double_t angle);
249  // Rotate the spatial component around the z-axis.
250 
251  inline void RotateUz(TVector3 & newUzVector);
252  // Rotates the reference frame from Uz to newUz (unit vector).
253 
254  inline void Rotate(Double_t, const TVector3 &);
255  // Rotate the spatial component around specified axis.
256 
257  inline TLorentzVector & operator *= (const TRotation &);
258  inline TLorentzVector & Transform(const TRotation &);
259  // Transformation with HepRotation.
260 
263  // Transformation with HepLorenzRotation.
264 
265  virtual void Print(Option_t *option="") const;
266 
267  ClassDef(TLorentzVector,4) // A four vector with (-,-,-,+) metric
268 };
269 
270 
271 //inline TLorentzVector operator * (const TLorentzVector &, Double_t a);
272 // moved to TLorentzVector::operator * (Double_t a)
274 // Scaling LorentzVector with a real number
275 
276 
277 inline Double_t TLorentzVector::X() const { return fP.X(); }
278 inline Double_t TLorentzVector::Y() const { return fP.Y(); }
279 inline Double_t TLorentzVector::Z() const { return fP.Z(); }
280 inline Double_t TLorentzVector::T() const { return fE; }
281 
282 inline void TLorentzVector::SetX(Double_t a) { fP.SetX(a); }
283 inline void TLorentzVector::SetY(Double_t a) { fP.SetY(a); }
284 inline void TLorentzVector::SetZ(Double_t a) { fP.SetZ(a); }
285 inline void TLorentzVector::SetT(Double_t a) { fE = a; }
286 
287 inline Double_t TLorentzVector::Px() const { return X(); }
288 inline Double_t TLorentzVector::Py() const { return Y(); }
289 inline Double_t TLorentzVector::Pz() const { return Z(); }
290 inline Double_t TLorentzVector::P() const { return fP.Mag(); }
291 inline Double_t TLorentzVector::E() const { return T(); }
292 inline Double_t TLorentzVector::Energy() const { return T(); }
293 
294 inline void TLorentzVector::SetPx(Double_t a) { SetX(a); }
295 inline void TLorentzVector::SetPy(Double_t a) { SetY(a); }
296 inline void TLorentzVector::SetPz(Double_t a) { SetZ(a); }
297 inline void TLorentzVector::SetE(Double_t a) { SetT(a); }
298 
299 inline TVector3 TLorentzVector::Vect() const { return fP; }
300 
301 inline void TLorentzVector::SetVect(const TVector3 &p) { fP = p; }
302 
304  return fP.Phi();
305 }
306 
308  return fP.Theta();
309 }
310 
312  return fP.CosTheta();
313 }
314 
315 
317  return fP.Mag();
318 }
319 
321  fP.SetTheta(th);
322 }
323 
325  fP.SetPhi(phi);
326 }
327 
329  fP.SetMag(rho);
330 }
331 
333  fP.SetXYZ(x, y, z);
334  SetT(t);
335 }
336 
338  SetXYZT(px, py, pz, e);
339 }
340 
342  if ( m >= 0 )
343  SetXYZT( x, y, z, TMath::Sqrt(x*x+y*y+z*z+m*m) );
344  else
345  SetXYZT( x, y, z, TMath::Sqrt( TMath::Max((x*x+y*y+z*z-m*m), 0. ) ) );
346 }
347 
349  pt = TMath::Abs(pt);
350  SetXYZM(pt*TMath::Cos(phi), pt*TMath::Sin(phi), pt*sinh(eta) ,m);
351 }
352 
354  pt = TMath::Abs(pt);
355  SetXYZT(pt*TMath::Cos(phi), pt*TMath::Sin(phi), pt*sinh(eta) ,e);
356 }
357 
358 inline void TLorentzVector::GetXYZT(Double_t *carray) const {
359  fP.GetXYZ(carray);
360  carray[3] = fE;
361 }
362 
363 inline void TLorentzVector::GetXYZT(Float_t *carray) const{
364  fP.GetXYZ(carray);
365  carray[3] = fE;
366 }
367 
368 inline Double_t & TLorentzVector::operator [] (int i) { return (*this)(i); }
369 inline Double_t TLorentzVector::operator [] (int i) const { return (*this)(i); }
370 
372  fP = q.Vect();
373  fE = q.T();
374  return *this;
375 }
376 
378  return TLorentzVector(fP+q.Vect(), fE+q.T());
379 }
380 
382  fP += q.Vect();
383  fE += q.T();
384  return *this;
385 }
386 
388  return TLorentzVector(fP-q.Vect(), fE-q.T());
389 }
390 
392  fP -= q.Vect();
393  fE -= q.T();
394  return *this;
395 }
396 
398  return TLorentzVector(-X(), -Y(), -Z(), -T());
399 }
400 
402  fP *= a;
403  fE *= a;
404  return *this;
405 }
406 
408  return TLorentzVector(a*X(), a*Y(), a*Z(), a*T());
409 }
410 
412  return (Vect() == q.Vect() && T() == q.T());
413 }
414 
416  return (Vect() != q.Vect() || T() != q.T());
417 }
418 
419 inline Double_t TLorentzVector::Perp2() const { return fP.Perp2(); }
420 
421 inline Double_t TLorentzVector::Perp() const { return fP.Perp(); }
422 
423 inline Double_t TLorentzVector::Pt() const { return Perp(); }
424 
426  fP.SetPerp(r);
427 }
428 
429 inline Double_t TLorentzVector::Perp2(const TVector3 &v) const {
430  return fP.Perp2(v);
431 }
432 
433 inline Double_t TLorentzVector::Perp(const TVector3 &v) const {
434  return fP.Perp(v);
435 }
436 
437 inline Double_t TLorentzVector::Pt(const TVector3 &v) const {
438  return Perp(v);
439 }
440 
442  Double_t pt2 = fP.Perp2();
443  return pt2 == 0 ? 0 : E()*E() * pt2/(pt2+Z()*Z());
444 }
445 
446 inline Double_t TLorentzVector::Et() const {
447  Double_t etet = Et2();
448  return E() < 0.0 ? -sqrt(etet) : sqrt(etet);
449 }
450 
451 inline Double_t TLorentzVector::Et2(const TVector3 & v) const {
452  Double_t pt2 = fP.Perp2(v);
453  Double_t pv = fP.Dot(v.Unit());
454  return pt2 == 0 ? 0 : E()*E() * pt2/(pt2+pv*pv);
455 }
456 
457 inline Double_t TLorentzVector::Et(const TVector3 & v) const {
458  Double_t etet = Et2(v);
459  return E() < 0.0 ? -sqrt(etet) : sqrt(etet);
460 }
461 
463  return TVector2::Phi_mpi_pi(Phi()-v.Phi());
464 }
465 
467  return PseudoRapidity();
468 }
470  Double_t deta = Eta()-v.Eta();
471  Double_t dphi = TVector2::Phi_mpi_pi(Phi()-v.Phi());
472  return TMath::Sqrt( deta*deta+dphi*dphi );
473 }
474 
476  return DeltaR(v);
477 }
478 
480  return TVector2 (Eta(),Phi());
481 }
482 
483 
484 inline Double_t TLorentzVector::Angle(const TVector3 &v) const {
485  return fP.Angle(v);
486 }
487 
489  return T()*T() - fP.Mag2();
490 }
491 
493  Double_t mm = Mag2();
494  return mm < 0.0 ? -TMath::Sqrt(-mm) : TMath::Sqrt(mm);
495 }
496 
497 inline Double_t TLorentzVector::M2() const { return Mag2(); }
498 inline Double_t TLorentzVector::M() const { return Mag(); }
499 
501  return E()*E() - Z()*Z();
502 }
503 
504 inline Double_t TLorentzVector::Mt() const {
505  Double_t mm = Mt2();
506  return mm < 0.0 ? -TMath::Sqrt(-mm) : TMath::Sqrt(mm);
507 }
508 
510  return fP.Mag() / fE;
511 }
512 
514  Double_t b = Beta();
515  return 1.0/TMath::Sqrt(1- b*b);
516 }
517 
518 inline void TLorentzVector::SetVectMag(const TVector3 & spatial, Double_t magnitude) {
519  SetXYZM(spatial.X(), spatial.Y(), spatial.Z(), magnitude);
520 }
521 
522 inline void TLorentzVector::SetVectM(const TVector3 & spatial, Double_t mass) {
523  SetVectMag(spatial, mass);
524 }
525 
527  return T()*q.T() - Z()*q.Z() - Y()*q.Y() - X()*q.X();
528 }
529 
531  return Dot(q);
532 }
533 
534 //Member functions Plus() and Minus() return the positive and negative
535 //light-cone components:
536 //
537 // Double_t pcone = v.Plus();
538 // Double_t mcone = v.Minus();
539 //
540 //CAVEAT: The values returned are T{+,-}Z. It is known that some authors
541 //find it easier to define these components as (T{+,-}Z)/sqrt(2). Thus
542 //check what definition is used in the physics you're working in and adapt
543 //your code accordingly.
544 
546  return T() + Z();
547 }
548 
550  return T() - Z();
551 }
552 
554  return TVector3(X()/T(), Y()/T(), Z()/T());
555 }
556 
557 inline void TLorentzVector::Boost(const TVector3 & b) {
558  Boost(b.X(), b.Y(), b.Z());
559 }
560 
562  return fP.PseudoRapidity();
563 }
564 
565 inline void TLorentzVector::RotateX(Double_t angle) {
566  fP.RotateX(angle);
567 }
568 
569 inline void TLorentzVector::RotateY(Double_t angle) {
570  fP.RotateY(angle);
571 }
572 
573 inline void TLorentzVector::RotateZ(Double_t angle) {
574  fP.RotateZ(angle);
575 }
576 
577 inline void TLorentzVector::RotateUz(TVector3 &newUzVector) {
578  fP.RotateUz(newUzVector);
579 }
580 
582  fP.Rotate(a,v);
583 }
584 
586  fP *= m;
587  return *this;
588 }
589 
591  fP.Transform(m);
592  return *this;
593 }
594 
596  return TLorentzVector(a*p.X(), a*p.Y(), a*p.Z(), a*p.T());
597 }
598 
600  : fP(), fE(0.0) {}
601 
603  : fP(x,y,z), fE(t) {}
604 
606  : fP(x0), fE(x0[3]) {}
607 
609  : fP(x0), fE(x0[3]) {}
610 
612  : fP(p), fE(e) {}
613 
615  , fP(p.Vect()), fE(p.T()) {}
616 
617 
618 
620 {
621  //dereferencing operator const
622  switch(i) {
623  case kX:
624  return fP.X();
625  case kY:
626  return fP.Y();
627  case kZ:
628  return fP.Z();
629  case kT:
630  return fE;
631  default:
632  Error("operator()()", "bad index (%d) returning 0",i);
633  }
634  return 0.;
635 }
636 
638 {
639  //dereferencing operator
640  switch(i) {
641  case kX:
642  return fP.fX;
643  case kY:
644  return fP.fY;
645  case kZ:
646  return fP.fZ;
647  case kT:
648  return fE;
649  default:
650  Error("operator()()", "bad index (%d) returning &fE",i);
651  }
652  return fE;
653 }
654 
655 #endif
Double_t X() const
Definition: TVector3.h:224
Double_t Dot(const TVector3 &) const
Definition: TVector3.h:339
The TLorentzRotation class describes Lorentz transformations including Lorentz boosts and rotations (...
void SetPxPyPzE(Double_t px, Double_t py, Double_t pz, Double_t e)
TLorentzVector operator+(const TLorentzVector &) const
void RotateUz(const TVector3 &)
NewUzVector must be normalized !
Definition: TVector3.cxx:299
TLorentzVector & operator+=(const TLorentzVector &)
TLorentzVector operator-() const
void Boost(Double_t, Double_t, Double_t)
Double_t Z() const
Double_t Plus() const
Double_t Rho() const
void SetPx(Double_t a)
void SetY(Double_t)
Definition: TVector3.h:232
void SetPerp(Double_t)
Definition: TVector3.h:396
TVector2 EtaPhiVector()
float Float_t
Definition: RtypesCore.h:53
Double_t Theta() const
const char Option_t
Definition: RtypesCore.h:62
Double_t CosTheta() const
Definition: TVector3.h:379
Double_t Gamma() const
TVector3 Vect() const
Double_t M2() const
Double_t Perp2() const
Double_t Mag2() const
Definition: TVector3.h:347
Double_t Angle(const TVector3 &v) const
Double_t Pt() const
bool Bool_t
Definition: RtypesCore.h:59
TArc * a
Definition: textangle.C:12
void RotateZ(Double_t angle)
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 fY
Definition: TVector3.h:193
void RotateY(Double_t angle)
Double_t Perp() const
void RotateX(Double_t)
Rotate vector around X.
Definition: TVector3.cxx:257
Double_t Angle(const TVector3 &) const
Return the angle w.r.t. another 3-vector.
Definition: TVector3.cxx:197
Short_t Abs(Short_t d)
Definition: TMathBase.h:110
Double_t CosTheta() const
TVector2 is a general two vector class, which can be used for the description of different vectors in...
Definition: TVector2.h:18
void SetVectMag(const TVector3 &spatial, Double_t magnitude)
TLorentzVector & operator*=(Double_t a)
Double_t Mag() const
Definition: TVector3.h:94
Double_t DeltaPhi(const TLorentzVector &) const
TVector3 Unit() const
Return unit vector parallel to this.
Definition: TVector3.cxx:246
Double_t DeltaR(const TLorentzVector &) const
double sqrt(double)
Double_t x[n]
Definition: legend1.C:17
Double_t Beta() const
Double_t Pz() const
#define ClassDef(name, id)
Definition: Rtypes.h:254
Double_t fZ
Definition: TVector3.h:193
void RotateUz(TVector3 &newUzVector)
Double_t Et() const
double sinh(double)
Double_t Y() const
Definition: TVector3.h:225
void Rotate(Double_t, const TVector3 &)
Rotate vector.
Definition: TVector3.cxx:290
static Double_t Phi_mpi_pi(Double_t x)
Returns phi angle in the interval [-PI,PI)
Definition: TVector2.cxx:102
void SetXYZ(Double_t x, Double_t y, Double_t z)
Definition: TVector3.h:235
void GetXYZ(Double_t *carray) const
Definition: TVector3.h:241
Double_t PseudoRapidity() const
Double_t Dot(const TLorentzVector &) const
void SetRho(Double_t rho)
void SetPhi(Double_t phi)
Double_t Minus() const
void SetPy(Double_t a)
Double_t Perp() const
Return the transverse component (R in cylindrical coordinate system)
Definition: TVector3.cxx:213
void RotateY(Double_t)
Rotate vector around Y.
Definition: TVector3.cxx:268
TVector3 is a general three vector class, which can be used for the description of different vectors ...
Definition: TVector3.h:30
void SetPerp(Double_t)
TLorentzVector & operator=(const TLorentzVector &)
Double_t Z() const
Definition: TVector3.h:226
The TRotation class describes a rotation of objects of the TVector3 class.
Definition: TRotation.h:22
TVector3 BoostVector() const
TPaveText * pt
void SetTheta(Double_t)
Set theta keeping mag and phi constant (BaBar).
Definition: TVector3.cxx:350
TRandom2 r(17)
TLorentzVector is a general four-vector class, which can be used either for the description of positi...
Double_t Eta() const
Double_t operator()(int i) const
void SetXYZT(Double_t x, Double_t y, Double_t z, Double_t t)
SVector< double, 2 > v
Definition: Dict.h:5
Double_t Mt() const
TLorentzVector & operator-=(const TLorentzVector &)
Bool_t operator!=(const TLorentzVector &) const
virtual ~TLorentzVector()
Double_t Theta() const
Return the polar angle.
Definition: TVector3.cxx:238
TMarker * m
Definition: textangle.C:8
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:925
TLorentzVector operator*(Double_t a) const
void SetX(Double_t a)
Double_t Energy() const
Double_t Mag2() const
void SetPtEtaPhiE(Double_t pt, Double_t eta, Double_t phi, Double_t e)
Double_t Mt2() const
Double_t Et2() const
Double_t PseudoRapidity() const
Double_t m = Mag(); return 0.5*log( (m+fZ)/(m-fZ) ); guard against Pt=0.
Definition: TVector3.cxx:320
Double_t P() const
void SetPz(Double_t a)
Double_t Cos(Double_t)
Definition: TMath.h:424
void SetTheta(Double_t theta)
Double_t fX
Definition: TVector3.h:193
double Double_t
Definition: RtypesCore.h:55
void SetE(Double_t a)
void SetZ(Double_t)
Definition: TVector3.h:233
void Rotate(Double_t, const TVector3 &)
Double_t Py() const
void SetXYZM(Double_t x, Double_t y, Double_t z, Double_t m)
void SetT(Double_t a)
Double_t Y() const
Double_t M() const
TVector3 & Transform(const TRotation &)
Transform this vector with a TRotation.
Definition: TVector3.cxx:190
Double_t y[n]
Definition: legend1.C:17
void SetY(Double_t a)
you should not use this method at all Int_t Int_t Double_t Double_t Double_t e
Definition: TRolke.cxx:630
Double_t Px() const
Double_t Phi() const
Double_t T() const
void RotateZ(Double_t)
Rotate vector around Z.
Definition: TVector3.cxx:279
Mother of all ROOT objects.
Definition: TObject.h:37
void SetZ(Double_t a)
you should not use this method at all Int_t Int_t z
Definition: TRolke.cxx:630
void SetPhi(Double_t)
Set phi keeping mag and theta constant (BaBar).
Definition: TVector3.cxx:362
void GetXYZT(Double_t *carray) const
Double_t Phi() const
Return the azimuth angle. Returns phi from -pi to pi.
Definition: TVector3.cxx:230
Double_t Perp2() const
Definition: TVector3.h:361
Short_t Max(Short_t a, Short_t b)
Definition: TMathBase.h:202
void SetMag(Double_t)
Definition: TVector3.h:384
Double_t Sin(Double_t)
Definition: TMath.h:421
Bool_t operator==(const TLorentzVector &) const
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 b
Definition: TRolke.cxx:630
void SetX(Double_t)
Definition: TVector3.h:231
Double_t Mag() const
Double_t X() const
TLorentzVector & Transform(const TRotation &)
void SetVectM(const TVector3 &spatial, Double_t mass)
Double_t E() const
void SetPtEtaPhiM(Double_t pt, Double_t eta, Double_t phi, Double_t m)
Double_t Sqrt(Double_t x)
Definition: TMath.h:464
Double_t Rapidity() const
float * q
Definition: THbookFile.cxx:87
Double_t DrEtaPhi(const TLorentzVector &) const
Double_t operator[](int i) const
void SetVect(const TVector3 &vect3)
void RotateX(Double_t angle)