Geant4  9.6.p02
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
ThreeVector.h
Go to the documentation of this file.
1 // -*- C++ -*-
2 // CLASSDOC OFF
3 // $Id:$
4 // ---------------------------------------------------------------------------
5 // CLASSDOC ON
6 //
7 // This file is a part of the CLHEP - a Class Library for High Energy Physics.
8 //
9 // Hep3Vector is a general 3-vector class defining vectors in three
10 // dimension using double components. Rotations of these vectors are
11 // performed by multiplying with an object of the HepRotation class.
12 //
13 // .SS See Also
14 // LorentzVector.h, Rotation.h, LorentzRotation.h
15 //
16 // .SS Authors
17 // Leif Lonnblad and Anders Nilsson; Modified by Evgueni Tcherniaev;
18 // ZOOM additions by Mark Fischler
19 //
20 
21 #ifndef HEP_THREEVECTOR_H
22 #define HEP_THREEVECTOR_H
23 
24 #ifdef GNUPRAGMA
25 #pragma interface
26 #endif
27 
28 #include <iostream>
29 #include "CLHEP/Utility/defs.h"
30 
31 namespace CLHEP {
32 
33 class HepRotation;
34 class HepEulerAngles;
35 class HepAxisAngle;
36 
41 class Hep3Vector {
42 
43 public:
44 
45 // Basic properties and operations on 3-vectors:
46 
47  enum { X=0, Y=1, Z=2, NUM_COORDINATES=3, SIZE=NUM_COORDINATES };
48  // Safe indexing of the coordinates when using with matrices, arrays, etc.
49  // (BaBar)
50 
51  Hep3Vector();
52  explicit Hep3Vector(double x);
53  Hep3Vector(double x, double y);
54  Hep3Vector(double x, double y, double z);
55  // The constructor.
56 
57  inline Hep3Vector(const Hep3Vector &);
58  // The copy constructor.
59 
60  inline ~Hep3Vector();
61  // The destructor. Not virtual - inheritance from this class is dangerous.
62 
63  double operator () (int) const;
64  // Get components by index -- 0-based (Geant4)
65 
66  inline double operator [] (int) const;
67  // Get components by index -- 0-based (Geant4)
68 
69  double & operator () (int);
70  // Set components by index. 0-based.
71 
72  inline double & operator [] (int);
73  // Set components by index. 0-based.
74 
75  inline double x() const;
76  inline double y() const;
77  inline double z() const;
78  // The components in cartesian coordinate system. Same as getX() etc.
79 
80  inline void setX(double);
81  inline void setY(double);
82  inline void setZ(double);
83  // Set the components in cartesian coordinate system.
84 
85  inline void set( double x, double y, double z);
86  // Set all three components in cartesian coordinate system.
87 
88  inline double phi() const;
89  // The azimuth angle.
90 
91  inline double theta() const;
92  // The polar angle.
93 
94  inline double cosTheta() const;
95  // Cosine of the polar angle.
96 
97  inline double cos2Theta() const;
98  // Cosine squared of the polar angle - faster than cosTheta(). (ZOOM)
99 
100  inline double mag2() const;
101  // The magnitude squared (r^2 in spherical coordinate system).
102 
103  inline double mag() const;
104  // The magnitude (r in spherical coordinate system).
105 
106  inline void setPhi(double);
107  // Set phi keeping mag and theta constant (BaBar).
108 
109  inline void setTheta(double);
110  // Set theta keeping mag and phi constant (BaBar).
111 
112  void setMag(double);
113  // Set magnitude keeping theta and phi constant (BaBar).
114 
115  inline double perp2() const;
116  // The transverse component squared (rho^2 in cylindrical coordinate system).
117 
118  inline double perp() const;
119  // The transverse component (rho in cylindrical coordinate system).
120 
121  inline void setPerp(double);
122  // Set the transverse component keeping phi and z constant.
123 
124  void setCylTheta(double);
125  // Set theta while keeping transvers component and phi fixed
126 
127  inline double perp2(const Hep3Vector &) const;
128  // The transverse component w.r.t. given axis squared.
129 
130  inline double perp(const Hep3Vector &) const;
131  // The transverse component w.r.t. given axis.
132 
133  inline Hep3Vector & operator = (const Hep3Vector &);
134  // Assignment.
135 
136  inline bool operator == (const Hep3Vector &) const;
137  inline bool operator != (const Hep3Vector &) const;
138  // Comparisons (Geant4).
139 
140  bool isNear (const Hep3Vector &, double epsilon=tolerance) const;
141  // Check for equality within RELATIVE tolerance (default 2.2E-14). (ZOOM)
142  // |v1 - v2|**2 <= epsilon**2 * |v1.dot(v2)|
143 
144  double howNear(const Hep3Vector & v ) const;
145  // sqrt ( |v1-v2|**2 / v1.dot(v2) ) with a maximum of 1.
146  // If v1.dot(v2) is negative, will return 1.
147 
148  double deltaR(const Hep3Vector & v) const;
149  // sqrt( pseudorapity_difference**2 + deltaPhi **2 )
150 
151  inline Hep3Vector & operator += (const Hep3Vector &);
152  // Addition.
153 
154  inline Hep3Vector & operator -= (const Hep3Vector &);
155  // Subtraction.
156 
157  inline Hep3Vector operator - () const;
158  // Unary minus.
159 
160  inline Hep3Vector & operator *= (double);
161  // Scaling with real numbers.
162 
163  Hep3Vector & operator /= (double);
164  // Division by (non-zero) real number.
165 
166  inline Hep3Vector unit() const;
167  // Vector parallel to this, but of length 1.
168 
169  inline Hep3Vector orthogonal() const;
170  // Vector orthogonal to this (Geant4).
171 
172  inline double dot(const Hep3Vector &) const;
173  // double product.
174 
175  inline Hep3Vector cross(const Hep3Vector &) const;
176  // Cross product.
177 
178  double angle(const Hep3Vector &) const;
179  // The angle w.r.t. another 3-vector.
180 
181  double pseudoRapidity() const;
182  // Returns the pseudo-rapidity, i.e. -ln(tan(theta/2))
183 
184  void setEta ( double p );
185  // Set pseudo-rapidity, keeping magnitude and phi fixed. (ZOOM)
186 
187  void setCylEta ( double p );
188  // Set pseudo-rapidity, keeping transverse component and phi fixed. (ZOOM)
189 
190  Hep3Vector & rotateX(double);
191  // Rotates the Hep3Vector around the x-axis.
192 
193  Hep3Vector & rotateY(double);
194  // Rotates the Hep3Vector around the y-axis.
195 
196  Hep3Vector & rotateZ(double);
197  // Rotates the Hep3Vector around the z-axis.
198 
199  Hep3Vector & rotateUz(const Hep3Vector&);
200  // Rotates reference frame from Uz to newUz (unit vector) (Geant4).
201 
202  Hep3Vector & rotate(double, const Hep3Vector &);
203  // Rotates around the axis specified by another Hep3Vector.
204  // (Uses methods of HepRotation, forcing linking in of Rotation.cc.)
205 
207  Hep3Vector & transform(const HepRotation &);
208  // Transformation with a Rotation matrix.
209 
210 // = = = = = = = = = = = = = = = = = = = = = = = =
211 //
212 // Esoteric properties and operations on 3-vectors:
213 //
214 // 1 - Set vectors in various coordinate systems
215 // 2 - Synonyms for accessing coordinates and properties
216 // 3 - Comparisions (dictionary, near-ness, and geometric)
217 // 4 - Intrinsic properties
218 // 5 - Properties releative to z axis and arbitrary directions
219 // 6 - Polar and azimuthal angle decomposition and deltaPhi
220 // 7 - Rotations
221 //
222 // = = = = = = = = = = = = = = = = = = = = = = = =
223 
224 // 1 - Set vectors in various coordinate systems
225 
226  inline void setRThetaPhi (double r, double theta, double phi);
227  // Set in spherical coordinates: Angles are measured in RADIANS
228 
229  inline void setREtaPhi ( double r, double eta, double phi );
230  // Set in spherical coordinates, but specify peudorapidiy to determine theta.
231 
232  inline void setRhoPhiZ (double rho, double phi, double z);
233  // Set in cylindrical coordinates: Phi angle is measured in RADIANS
234 
235  void setRhoPhiTheta ( double rho, double phi, double theta);
236  // Set in cylindrical coordinates, but specify theta to determine z.
237 
238  void setRhoPhiEta ( double rho, double phi, double eta);
239  // Set in cylindrical coordinates, but specify pseudorapidity to determine z.
240 
241 // 2 - Synonyms for accessing coordinates and properties
242 
243  inline double getX() const;
244  inline double getY() const;
245  inline double getZ() const;
246  // x(), y(), and z()
247 
248  inline double getR () const;
249  inline double getTheta() const;
250  inline double getPhi () const;
251  // mag(), theta(), and phi()
252 
253  inline double r () const;
254  // mag()
255 
256  inline double rho () const;
257  inline double getRho () const;
258  // perp()
259 
260  double eta () const;
261  double getEta () const;
262  // pseudoRapidity()
263 
264  inline void setR ( double s );
265  // setMag()
266 
267  inline void setRho ( double s );
268  // setPerp()
269 
270 // 3 - Comparisions (dictionary, near-ness, and geometric)
271 
272  int compare (const Hep3Vector & v) const;
273  bool operator > (const Hep3Vector & v) const;
274  bool operator < (const Hep3Vector & v) const;
275  bool operator>= (const Hep3Vector & v) const;
276  bool operator<= (const Hep3Vector & v) const;
277  // dictionary ordering according to z, then y, then x component
278 
279  inline double diff2 (const Hep3Vector & v) const;
280  // |v1-v2|**2
281 
282  static double setTolerance (double tol);
283  static inline double getTolerance ();
284  // Set the tolerance used in isNear() for Hep3Vectors
285 
286  bool isParallel (const Hep3Vector & v, double epsilon=tolerance) const;
287  // Are the vectors parallel, within the given tolerance?
288 
289  bool isOrthogonal (const Hep3Vector & v, double epsilon=tolerance) const;
290  // Are the vectors orthogonal, within the given tolerance?
291 
292  double howParallel (const Hep3Vector & v) const;
293  // | v1.cross(v2) / v1.dot(v2) |, to a maximum of 1.
294 
295  double howOrthogonal (const Hep3Vector & v) const;
296  // | v1.dot(v2) / v1.cross(v2) |, to a maximum of 1.
297 
298  enum { ToleranceTicks = 100 };
299 
300 // 4 - Intrinsic properties
301 
302  double beta () const;
303  // relativistic beta (considering v as a velocity vector with c=1)
304  // Same as mag() but will object if >= 1
305 
306  double gamma() const;
307  // relativistic gamma (considering v as a velocity vector with c=1)
308 
309  double coLinearRapidity() const;
310  // inverse tanh (beta)
311 
312 // 5 - Properties relative to Z axis and to an arbitrary direction
313 
314  // Note that the non-esoteric CLHEP provides
315  // theta(), cosTheta(), cos2Theta, and angle(const Hep3Vector&)
316 
317  inline double angle() const;
318  // angle against the Z axis -- synonym for theta()
319 
320  inline double theta(const Hep3Vector & v2) const;
321  // synonym for angle(v2)
322 
323  double cosTheta (const Hep3Vector & v2) const;
324  double cos2Theta(const Hep3Vector & v2) const;
325  // cos and cos^2 of the angle between two vectors
326 
327  inline Hep3Vector project () const;
328  Hep3Vector project (const Hep3Vector & v2) const;
329  // projection of a vector along a direction.
330 
331  inline Hep3Vector perpPart() const;
332  inline Hep3Vector perpPart (const Hep3Vector & v2) const;
333  // vector minus its projection along a direction.
334 
335  double rapidity () const;
336  // inverse tanh(v.z())
337 
338  double rapidity (const Hep3Vector & v2) const;
339  // rapidity with respect to specified direction:
340  // inverse tanh (v.dot(u)) where u is a unit in the direction of v2
341 
342  double eta(const Hep3Vector & v2) const;
343  // - ln tan of the angle beween the vector and the ref direction.
344 
345 // 6 - Polar and azimuthal angle decomposition and deltaPhi
346 
347  // Decomposition of an angle within reference defined by a direction:
348 
349  double polarAngle (const Hep3Vector & v2) const;
350  // The reference direction is Z: the polarAngle is abs(v.theta()-v2.theta()).
351 
352  double deltaPhi (const Hep3Vector & v2) const;
353  // v.phi()-v2.phi(), brought into the range (-PI,PI]
354 
355  double azimAngle (const Hep3Vector & v2) const;
356  // The reference direction is Z: the azimAngle is the same as deltaPhi
357 
358  double polarAngle (const Hep3Vector & v2,
359  const Hep3Vector & ref) const;
360  // For arbitrary reference direction,
361  // polarAngle is abs(v.angle(ref) - v2.angle(ref)).
362 
363  double azimAngle (const Hep3Vector & v2,
364  const Hep3Vector & ref) const;
365  // To compute azimangle, project v and v2 into the plane normal to
366  // the reference direction. Then in that plane take the angle going
367  // clockwise around the direction from projection of v to that of v2.
368 
369 // 7 - Rotations
370 
371 // These mehtods **DO NOT** use anything in the HepRotation class.
372 // Thus, use of v.rotate(axis,delta) does not force linking in Rotation.cc.
373 
374  Hep3Vector & rotate (const Hep3Vector & axis, double delta);
375  // Synonym for rotate (delta, axis)
376 
377  Hep3Vector & rotate (const HepAxisAngle & ax);
378  // HepAxisAngle is a struct holding an axis direction and an angle.
379 
380  Hep3Vector & rotate (const HepEulerAngles & e);
381  Hep3Vector & rotate (double phi,
382  double theta,
383  double psi);
384  // Rotate via Euler Angles. Our Euler Angles conventions are
385  // those of Goldstein Classical Mechanics page 107.
386 
387 protected:
388  void setSpherical (double r, double theta, double phi);
389  void setCylindrical (double r, double phi, double z);
390  double negativeInfinity() const;
391 
392 protected:
393 
394  double dx;
395  double dy;
396  double dz;
397  // The components.
398 
399  DLL_API static double tolerance;
400  // default tolerance criterion for isNear() to return true.
401 }; // Hep3Vector
402 
403 // Global Methods
404 
405 Hep3Vector rotationXOf (const Hep3Vector & vec, double delta);
406 Hep3Vector rotationYOf (const Hep3Vector & vec, double delta);
407 Hep3Vector rotationZOf (const Hep3Vector & vec, double delta);
408 
409 Hep3Vector rotationOf (const Hep3Vector & vec,
410  const Hep3Vector & axis, double delta);
411 Hep3Vector rotationOf (const Hep3Vector & vec, const HepAxisAngle & ax);
412 
413 Hep3Vector rotationOf (const Hep3Vector & vec,
414  double phi, double theta, double psi);
415 Hep3Vector rotationOf (const Hep3Vector & vec, const HepEulerAngles & e);
416 // Return a new vector based on a rotation of the supplied vector
417 
418 std::ostream & operator << (std::ostream &, const Hep3Vector &);
419 // Output to a stream.
420 
421 std::istream & operator >> (std::istream &, Hep3Vector &);
422 // Input from a stream.
423 
424 extern DLL_API const Hep3Vector HepXHat, HepYHat, HepZHat;
425 
428 
429 Hep3Vector operator / (const Hep3Vector &, double a);
430 // Division of 3-vectors by non-zero real number
431 
432 inline Hep3Vector operator + (const Hep3Vector &, const Hep3Vector &);
433 // Addition of 3-vectors.
434 
435 inline Hep3Vector operator - (const Hep3Vector &, const Hep3Vector &);
436 // Subtraction of 3-vectors.
437 
438 inline double operator * (const Hep3Vector &, const Hep3Vector &);
439 // double product of 3-vectors.
440 
441 inline Hep3Vector operator * (const Hep3Vector &, double a);
442 inline Hep3Vector operator * (double a, const Hep3Vector &);
443 // Scaling of 3-vectors with a real number
444 
445 } // namespace CLHEP
446 
447 #include "CLHEP/Vector/ThreeVector.icc"
448 
449 #endif /* HEP_THREEVECTOR_H */