Logo ROOT   6.08/07
Reference Guide
TVector3.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 #ifndef ROOT_TVector3
12 #define ROOT_TVector3
13 
14 #ifndef ROOT_TError
15 #include "TError.h"
16 #endif
17 #ifndef ROOT_TVector2
18 #include "TVector2.h"
19 #endif
20 #ifndef ROOT_TMatrix
21 #include "TMatrix.h"
22 #endif
23 #ifndef ROOT_TMath
24 #include "TMath.h"
25 #endif
26 
27 class TRotation;
28 
29 
30 class TVector3 : public TObject {
31 
32 public:
33 
34  typedef Double_t Scalar; // to be able to use it with the ROOT::Math::VectorUtil functions
35 
36  TVector3();
37 
39  // The constructor.
40 
41  TVector3(const Double_t *);
42  TVector3(const Float_t *);
43  // Constructors from an array
44 
45  TVector3(const TVector3 &);
46  // The copy constructor.
47 
48  virtual ~TVector3() {};
49  // Destructor
50 
51  Double_t operator () (int) const;
52  inline Double_t operator [] (int) const;
53  // Get components by index (Geant4).
54 
55  Double_t & operator () (int);
56  inline Double_t & operator [] (int);
57  // Set components by index.
58 
59  inline Double_t x() const;
60  inline Double_t y() const;
61  inline Double_t z() const;
62  inline Double_t X() const;
63  inline Double_t Y() const;
64  inline Double_t Z() const;
65  inline Double_t Px() const;
66  inline Double_t Py() const;
67  inline Double_t Pz() const;
68  // The components in cartesian coordinate system.
69 
70  inline void SetX(Double_t);
71  inline void SetY(Double_t);
72  inline void SetZ(Double_t);
73  inline void SetXYZ(Double_t x, Double_t y, Double_t z);
74  void SetPtEtaPhi(Double_t pt, Double_t eta, Double_t phi);
75  void SetPtThetaPhi(Double_t pt, Double_t theta, Double_t phi);
76 
77  inline void GetXYZ(Double_t *carray) const;
78  inline void GetXYZ(Float_t *carray) const;
79  // Get the components into an array
80  // not checked!
81 
82  Double_t Phi() const;
83  // The azimuth angle. returns phi from -pi to pi
84 
85  Double_t Theta() const;
86  // The polar angle.
87 
88  inline Double_t CosTheta() const;
89  // Cosine of the polar angle.
90 
91  inline Double_t Mag2() const;
92  // The magnitude squared (rho^2 in spherical coordinate system).
93 
94  Double_t Mag() const { return TMath::Sqrt(Mag2()); }
95  // The magnitude (rho in spherical coordinate system).
96 
97  void SetPhi(Double_t);
98  // Set phi keeping mag and theta constant (BaBar).
99 
100  void SetTheta(Double_t);
101  // Set theta keeping mag and phi constant (BaBar).
102 
103  inline void SetMag(Double_t);
104  // Set magnitude keeping theta and phi constant (BaBar).
105 
106  inline Double_t Perp2() const;
107  // The transverse component squared (R^2 in cylindrical coordinate system).
108 
109  inline Double_t Pt() const;
110  Double_t Perp() const;
111  // The transverse component (R in cylindrical coordinate system).
112 
113  inline void SetPerp(Double_t);
114  // Set the transverse component keeping phi and z constant.
115 
116  inline Double_t Perp2(const TVector3 &) const;
117  // The transverse component w.r.t. given axis squared.
118 
119  inline Double_t Pt(const TVector3 &) const;
120  Double_t Perp(const TVector3 &) const;
121  // The transverse component w.r.t. given axis.
122 
123  inline Double_t DeltaPhi(const TVector3 &) const;
124  Double_t DeltaR(const TVector3 &) const;
125  inline Double_t DrEtaPhi(const TVector3 &) const;
126  inline TVector2 EtaPhiVector() const;
127  void SetMagThetaPhi(Double_t mag, Double_t theta, Double_t phi);
128 
129  inline TVector3 & operator = (const TVector3 &);
130  // Assignment.
131 
132  inline Bool_t operator == (const TVector3 &) const;
133  inline Bool_t operator != (const TVector3 &) const;
134  // Comparisons (Geant4).
135 
136  inline TVector3 & operator += (const TVector3 &);
137  // Addition.
138 
139  inline TVector3 & operator -= (const TVector3 &);
140  // Subtraction.
141 
142  inline TVector3 operator - () const;
143  // Unary minus.
144 
145  inline TVector3 & operator *= (Double_t);
146  // Scaling with real numbers.
147 
148  TVector3 Unit() const;
149  // Unit vector parallel to this.
150 
151  inline TVector3 Orthogonal() const;
152  // Vector orthogonal to this (Geant4).
153 
154  inline Double_t Dot(const TVector3 &) const;
155  // Scalar product.
156 
157  inline TVector3 Cross(const TVector3 &) const;
158  // Cross product.
159 
160  Double_t Angle(const TVector3 &) const;
161  // The angle w.r.t. another 3-vector.
162 
163  Double_t PseudoRapidity() const;
164  // Returns the pseudo-rapidity, i.e. -ln(tan(theta/2))
165 
166  inline Double_t Eta() const;
167 
168  void RotateX(Double_t);
169  // Rotates the Hep3Vector around the x-axis.
170 
171  void RotateY(Double_t);
172  // Rotates the Hep3Vector around the y-axis.
173 
174  void RotateZ(Double_t);
175  // Rotates the Hep3Vector around the z-axis.
176 
177  void RotateUz(const TVector3&);
178  // Rotates reference frame from Uz to newUz (unit vector) (Geant4).
179 
180  void Rotate(Double_t, const TVector3 &);
181  // Rotates around the axis specified by another Hep3Vector.
182 
183  TVector3 & operator *= (const TRotation &);
184  TVector3 & Transform(const TRotation &);
185  // Transformation with a Rotation matrix.
186 
187  inline TVector2 XYvector() const;
188 
189  void Print(Option_t* option="") const;
190 
191 private:
192 
194  // The components.
195 
196  ClassDef(TVector3,3) // A 3D physics vector
197 
198  // make TLorentzVector a friend class
199  friend class TLorentzVector;
200 };
201 
202 TVector3 operator + (const TVector3 &, const TVector3 &);
203 // Addition of 3-vectors.
204 
205 TVector3 operator - (const TVector3 &, const TVector3 &);
206 // Subtraction of 3-vectors.
207 
208 Double_t operator * (const TVector3 &, const TVector3 &);
209 // Scalar product of 3-vectors.
210 
211 TVector3 operator * (const TVector3 &, Double_t a);
212 TVector3 operator * (Double_t a, const TVector3 &);
213 // Scaling of 3-vectors with a real number
214 
215 TVector3 operator * (const TMatrix &, const TVector3 &);
216 
217 
218 inline Double_t & TVector3::operator[] (int i) { return operator()(i); }
219 inline Double_t TVector3::operator[] (int i) const { return operator()(i); }
220 
221 inline Double_t TVector3::x() const { return fX; }
222 inline Double_t TVector3::y() const { return fY; }
223 inline Double_t TVector3::z() const { return fZ; }
224 inline Double_t TVector3::X() const { return fX; }
225 inline Double_t TVector3::Y() const { return fY; }
226 inline Double_t TVector3::Z() const { return fZ; }
227 inline Double_t TVector3::Px() const { return fX; }
228 inline Double_t TVector3::Py() const { return fY; }
229 inline Double_t TVector3::Pz() const { return fZ; }
230 
231 inline void TVector3::SetX(Double_t xx) { fX = xx; }
232 inline void TVector3::SetY(Double_t yy) { fY = yy; }
233 inline void TVector3::SetZ(Double_t zz) { fZ = zz; }
234 
235 inline void TVector3::SetXYZ(Double_t xx, Double_t yy, Double_t zz) {
236  fX = xx;
237  fY = yy;
238  fZ = zz;
239 }
240 
241 inline void TVector3::GetXYZ(Double_t *carray) const {
242  carray[0] = fX;
243  carray[1] = fY;
244  carray[2] = fZ;
245 }
246 
247 inline void TVector3::GetXYZ(Float_t *carray) const {
248  carray[0] = fX;
249  carray[1] = fY;
250  carray[2] = fZ;
251 }
252 
253 ////////////////////////////////////////////////////////////////////////////////
254 /// Constructors
256 : fX(0.0), fY(0.0), fZ(0.0) {}
257 
258 inline TVector3::TVector3(const TVector3 & p) : TObject(p),
259  fX(p.fX), fY(p.fY), fZ(p.fZ) {}
260 
262 : fX(xx), fY(yy), fZ(zz) {}
263 
264 inline TVector3::TVector3(const Double_t * x0)
265 : fX(x0[0]), fY(x0[1]), fZ(x0[2]) {}
266 
267 inline TVector3::TVector3(const Float_t * x0)
268 : fX(x0[0]), fY(x0[1]), fZ(x0[2]) {}
269 
270 
271 inline Double_t TVector3::operator () (int i) const {
272  switch(i) {
273  case 0:
274  return fX;
275  case 1:
276  return fY;
277  case 2:
278  return fZ;
279  default:
280  Error("operator()(i)", "bad index (%d) returning 0",i);
281  }
282  return 0.;
283 }
284 
286  switch(i) {
287  case 0:
288  return fX;
289  case 1:
290  return fY;
291  case 2:
292  return fZ;
293  default:
294  Error("operator()(i)", "bad index (%d) returning &fX",i);
295  }
296  return fX;
297 }
298 
300  fX = p.fX;
301  fY = p.fY;
302  fZ = p.fZ;
303  return *this;
304 }
305 
306 inline Bool_t TVector3::operator == (const TVector3& v) const {
307  return (v.fX==fX && v.fY==fY && v.fZ==fZ) ? kTRUE : kFALSE;
308 }
309 
310 inline Bool_t TVector3::operator != (const TVector3& v) const {
311  return (v.fX!=fX || v.fY!=fY || v.fZ!=fZ) ? kTRUE : kFALSE;
312 }
313 
315  fX += p.fX;
316  fY += p.fY;
317  fZ += p.fZ;
318  return *this;
319 }
320 
322  fX -= p.fX;
323  fY -= p.fY;
324  fZ -= p.fZ;
325  return *this;
326 }
327 
329  return TVector3(-fX, -fY, -fZ);
330 }
331 
333  fX *= a;
334  fY *= a;
335  fZ *= a;
336  return *this;
337 }
338 
339 inline Double_t TVector3::Dot(const TVector3 & p) const {
340  return fX*p.fX + fY*p.fY + fZ*p.fZ;
341 }
342 
343 inline TVector3 TVector3::Cross(const TVector3 & p) const {
344  return TVector3(fY*p.fZ-p.fY*fZ, fZ*p.fX-p.fZ*fX, fX*p.fY-p.fX*fY);
345 }
346 
347 inline Double_t TVector3::Mag2() const { return fX*fX + fY*fY + fZ*fZ; }
348 
349 
351  Double_t xx = fX < 0.0 ? -fX : fX;
352  Double_t yy = fY < 0.0 ? -fY : fY;
353  Double_t zz = fZ < 0.0 ? -fZ : fZ;
354  if (xx < yy) {
355  return xx < zz ? TVector3(0,fZ,-fY) : TVector3(fY,-fX,0);
356  } else {
357  return yy < zz ? TVector3(-fZ,0,fX) : TVector3(fY,-fX,0);
358  }
359 }
360 
361 inline Double_t TVector3::Perp2() const { return fX*fX + fY*fY; }
362 
363 
364 inline Double_t TVector3::Pt() const { return Perp(); }
365 
366 inline Double_t TVector3::Perp2(const TVector3 & p) const {
367  Double_t tot = p.Mag2();
368  Double_t ss = Dot(p);
369  Double_t per = Mag2();
370  if (tot > 0.0) per -= ss*ss/tot;
371  if (per < 0) per = 0;
372  return per;
373 }
374 
375 inline Double_t TVector3::Pt(const TVector3 & p) const {
376  return Perp(p);
377 }
378 
379 inline Double_t TVector3::CosTheta() const {
380  Double_t ptot = Mag();
381  return ptot == 0.0 ? 1.0 : fZ/ptot;
382 }
383 
384 inline void TVector3::SetMag(Double_t ma) {
385  Double_t factor = Mag();
386  if (factor == 0) {
387  Warning("SetMag","zero vector can't be stretched");
388  } else {
389  factor = ma/factor;
390  SetX(fX*factor);
391  SetY(fY*factor);
392  SetZ(fZ*factor);
393  }
394 }
395 
397  Double_t p = Perp();
398  if (p != 0.0) {
399  fX *= r/p;
400  fY *= r/p;
401  }
402 }
403 
404 inline Double_t TVector3::DeltaPhi(const TVector3 & v) const {
405  return TVector2::Phi_mpi_pi(Phi()-v.Phi());
406 }
407 
408 inline Double_t TVector3::Eta() const {
409  return PseudoRapidity();
410 }
411 
412 inline Double_t TVector3::DrEtaPhi(const TVector3 & v) const{
413  return DeltaR(v);
414 }
415 
416 
418  return TVector2 (Eta(),Phi());
419 }
420 
421 inline TVector2 TVector3::XYvector() const {
422  return TVector2(fX,fY);
423 }
424 
425 
426 #endif
TVector3 & operator=(const TVector3 &)
Definition: TVector3.h:299
Double_t X() const
Definition: TVector3.h:224
Double_t Dot(const TVector3 &) const
Definition: TVector3.h:339
void RotateUz(const TVector3 &)
NewUzVector must be normalized !
Definition: TVector3.cxx:299
virtual ~TVector3()
Definition: TVector3.h:48
Double_t DeltaPhi(const TVector3 &) const
Definition: TVector3.h:404
void SetY(Double_t)
Definition: TVector3.h:232
void SetPerp(Double_t)
Definition: TVector3.h:396
float Float_t
Definition: RtypesCore.h:53
Bool_t operator==(const TVector3 &) const
Definition: TVector3.h:306
const char Option_t
Definition: RtypesCore.h:62
Double_t operator()(int) const
Definition: TVector3.h:271
Double_t CosTheta() const
Definition: TVector3.h:379
Double_t Pt() const
Definition: TVector3.h:364
Double_t Eta() const
Definition: TVector3.h:408
Double_t Mag2() const
Definition: TVector3.h:347
bool Bool_t
Definition: RtypesCore.h:59
TArc * a
Definition: textangle.C:12
const Bool_t kFALSE
Definition: Rtypes.h:92
Double_t fY
Definition: TVector3.h:193
Double_t y() const
Definition: TVector3.h:222
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
Double_t DrEtaPhi(const TVector3 &) const
Definition: TVector3.h:412
TVector2 is a general two vector class, which can be used for the description of different vectors in...
Definition: TVector2.h:18
Double_t Mag() const
Definition: TVector3.h:94
TMatrixT.
Definition: TMatrixDfwd.h:24
TVector3 Unit() const
Return unit vector parallel to this.
Definition: TVector3.cxx:246
TVector3 operator-() const
Definition: TVector3.h:328
#define ClassDef(name, id)
Definition: Rtypes.h:254
Double_t fZ
Definition: TVector3.h:193
TVector3 & operator*=(Double_t)
Definition: TVector3.h:332
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
TVector2 EtaPhiVector() const
Definition: TVector3.h:417
TVector2 XYvector() const
Definition: TVector3.h:421
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
TVector3 Cross(const TVector3 &) const
Definition: TVector3.h:343
Double_t Z() const
Definition: TVector3.h:226
The TRotation class describes a rotation of objects of the TVector3 class.
Definition: TRotation.h:22
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...
SVector< double, 2 > v
Definition: Dict.h:5
TVector3 & operator+=(const TVector3 &)
Definition: TVector3.h:314
Double_t Theta() const
Return the polar angle.
Definition: TVector3.cxx:238
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:925
void Print(Option_t *option="") const
Print vector parameters.
Definition: TVector3.cxx:460
Double_t Pz() const
Definition: TVector3.h:229
TVector3 & operator-=(const TVector3 &)
Definition: TVector3.h:321
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 Py() const
Definition: TVector3.h:228
Double_t z() const
Definition: TVector3.h:223
Bool_t operator!=(const TVector3 &) const
Definition: TVector3.h:310
Double_t fX
Definition: TVector3.h:193
double Double_t
Definition: RtypesCore.h:55
void SetZ(Double_t)
Definition: TVector3.h:233
TVector3 Orthogonal() const
Definition: TVector3.h:350
TVector3 & Transform(const TRotation &)
Transform this vector with a TRotation.
Definition: TVector3.cxx:190
void RotateZ(Double_t)
Rotate vector around Z.
Definition: TVector3.cxx:279
Mother of all ROOT objects.
Definition: TObject.h:37
void SetPhi(Double_t)
Set phi keeping mag and theta constant (BaBar).
Definition: TVector3.cxx:362
Double_t DeltaR(const TVector3 &) const
Return deltaR with respect to v.
Definition: TVector3.cxx:372
Double_t operator[](int) const
Definition: TVector3.h:219
Double_t Phi() const
Return the azimuth angle. Returns phi from -pi to pi.
Definition: TVector3.cxx:230
Double_t x() const
Definition: TVector3.h:221
TVector3()
Constructors.
Definition: TVector3.h:255
Double_t Px() const
Definition: TVector3.h:227
Double_t Perp2() const
Definition: TVector3.h:361
void SetMag(Double_t)
Definition: TVector3.h:384
void SetX(Double_t)
Definition: TVector3.h:231
Double_t Sqrt(Double_t x)
Definition: TMath.h:464
void SetPtThetaPhi(Double_t pt, Double_t theta, Double_t phi)
Set Pt, Theta and Phi.
Definition: TVector3.cxx:340
const Bool_t kTRUE
Definition: Rtypes.h:91
void SetMagThetaPhi(Double_t mag, Double_t theta, Double_t phi)
Setter with mag, theta, phi.
Definition: TVector3.cxx:382
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition: TObject.cxx:911
Double_t Scalar
Definition: TVector3.h:34
void SetPtEtaPhi(Double_t pt, Double_t eta, Double_t phi)
Set Pt, Eta and Phi.
Definition: TVector3.cxx:332