CLHEP 2.0.4.7 Reference Documentation
CLHEP Home Page CLHEP Documentation CLHEP Bug Reports |
00001 // -*- C++ -*- 00002 // CLASSDOC OFF 00003 // $Id: ThreeVector.h,v 1.2.4.1.2.2 2009/11/11 23:18:26 garren Exp $ 00004 // --------------------------------------------------------------------------- 00005 // CLASSDOC ON 00006 // 00007 // This file is a part of the CLHEP - a Class Library for High Energy Physics. 00008 // 00009 // Hep3Vector is a general 3-vector class defining vectors in three 00010 // dimension using double components. Rotations of these vectors are 00011 // performed by multiplying with an object of the HepRotation class. 00012 // 00013 // .SS See Also 00014 // LorentzVector.h, Rotation.h, LorentzRotation.h 00015 // 00016 // .SS Authors 00017 // Leif Lonnblad and Anders Nilsson; Modified by Evgueni Tcherniaev; 00018 // ZOOM additions by Mark Fischler 00019 // 00020 00021 #ifndef HEP_THREEVECTOR_H 00022 #define HEP_THREEVECTOR_H 00023 00024 #ifdef GNUPRAGMA 00025 #pragma interface 00026 #endif 00027 00028 #include <iostream> 00029 #include "CLHEP/Vector/defs.h" 00030 00031 namespace CLHEP { 00032 00033 class HepRotation; 00034 class HepEulerAngles; 00035 class HepAxisAngle; 00036 00041 class Hep3Vector { 00042 00043 public: 00044 00045 // Basic properties and operations on 3-vectors: 00046 00047 enum { X=0, Y=1, Z=2, NUM_COORDINATES=3, SIZE=NUM_COORDINATES }; 00048 // Safe indexing of the coordinates when using with matrices, arrays, etc. 00049 // (BaBar) 00050 00051 inline Hep3Vector(double x = 0.0, double y = 0.0, double z = 0.0); 00052 // The constructor. 00053 00054 inline Hep3Vector(const Hep3Vector &); 00055 // The copy constructor. 00056 00057 inline ~Hep3Vector(); 00058 // The destructor. Not virtual - inheritance from this class is dangerous. 00059 00060 double operator () (int) const; 00061 // Get components by index -- 0-based (Geant4) 00062 00063 inline double operator [] (int) const; 00064 // Get components by index -- 0-based (Geant4) 00065 00066 double & operator () (int); 00067 // Set components by index. 0-based. 00068 00069 inline double & operator [] (int); 00070 // Set components by index. 0-based. 00071 00072 inline double x() const; 00073 inline double y() const; 00074 inline double z() const; 00075 // The components in cartesian coordinate system. Same as getX() etc. 00076 00077 inline void setX(double); 00078 inline void setY(double); 00079 inline void setZ(double); 00080 // Set the components in cartesian coordinate system. 00081 00082 inline void set( double x, double y, double z); 00083 // Set all three components in cartesian coordinate system. 00084 00085 inline double phi() const; 00086 // The azimuth angle. 00087 00088 inline double theta() const; 00089 // The polar angle. 00090 00091 inline double cosTheta() const; 00092 // Cosine of the polar angle. 00093 00094 inline double cos2Theta() const; 00095 // Cosine squared of the polar angle - faster than cosTheta(). (ZOOM) 00096 00097 inline double mag2() const; 00098 // The magnitude squared (r^2 in spherical coordinate system). 00099 00100 inline double mag() const; 00101 // The magnitude (r in spherical coordinate system). 00102 00103 inline void setPhi(double); 00104 // Set phi keeping mag and theta constant (BaBar). 00105 00106 inline void setTheta(double); 00107 // Set theta keeping mag and phi constant (BaBar). 00108 00109 void setMag(double); 00110 // Set magnitude keeping theta and phi constant (BaBar). 00111 00112 inline double perp2() const; 00113 // The transverse component squared (rho^2 in cylindrical coordinate system). 00114 00115 inline double perp() const; 00116 // The transverse component (rho in cylindrical coordinate system). 00117 00118 inline void setPerp(double); 00119 // Set the transverse component keeping phi and z constant. 00120 00121 void setCylTheta(double); 00122 // Set theta while keeping transvers component and phi fixed 00123 00124 inline double perp2(const Hep3Vector &) const; 00125 // The transverse component w.r.t. given axis squared. 00126 00127 inline double perp(const Hep3Vector &) const; 00128 // The transverse component w.r.t. given axis. 00129 00130 inline Hep3Vector & operator = (const Hep3Vector &); 00131 // Assignment. 00132 00133 inline bool operator == (const Hep3Vector &) const; 00134 inline bool operator != (const Hep3Vector &) const; 00135 // Comparisons (Geant4). 00136 00137 bool isNear (const Hep3Vector &, double epsilon=tolerance) const; 00138 // Check for equality within RELATIVE tolerance (default 2.2E-14). (ZOOM) 00139 // |v1 - v2|**2 <= epsilon**2 * |v1.dot(v2)| 00140 00141 double howNear(const Hep3Vector & v ) const; 00142 // sqrt ( |v1-v2|**2 / v1.dot(v2) ) with a maximum of 1. 00143 // If v1.dot(v2) is negative, will return 1. 00144 00145 double deltaR(const Hep3Vector & v) const; 00146 // sqrt( pseudorapity_difference**2 + deltaPhi **2 ) 00147 00148 inline Hep3Vector & operator += (const Hep3Vector &); 00149 // Addition. 00150 00151 inline Hep3Vector & operator -= (const Hep3Vector &); 00152 // Subtraction. 00153 00154 inline Hep3Vector operator - () const; 00155 // Unary minus. 00156 00157 inline Hep3Vector & operator *= (double); 00158 // Scaling with real numbers. 00159 00160 Hep3Vector & operator /= (double); 00161 // Division by (non-zero) real number. 00162 00163 inline Hep3Vector unit() const; 00164 // Vector parallel to this, but of length 1. 00165 00166 inline Hep3Vector orthogonal() const; 00167 // Vector orthogonal to this (Geant4). 00168 00169 inline double dot(const Hep3Vector &) const; 00170 // double product. 00171 00172 inline Hep3Vector cross(const Hep3Vector &) const; 00173 // Cross product. 00174 00175 double angle(const Hep3Vector &) const; 00176 // The angle w.r.t. another 3-vector. 00177 00178 double pseudoRapidity() const; 00179 // Returns the pseudo-rapidity, i.e. -ln(tan(theta/2)) 00180 00181 void setEta ( double p ); 00182 // Set pseudo-rapidity, keeping magnitude and phi fixed. (ZOOM) 00183 00184 void setCylEta ( double p ); 00185 // Set pseudo-rapidity, keeping transverse component and phi fixed. (ZOOM) 00186 00187 Hep3Vector & rotateX(double); 00188 // Rotates the Hep3Vector around the x-axis. 00189 00190 Hep3Vector & rotateY(double); 00191 // Rotates the Hep3Vector around the y-axis. 00192 00193 Hep3Vector & rotateZ(double); 00194 // Rotates the Hep3Vector around the z-axis. 00195 00196 Hep3Vector & rotateUz(const Hep3Vector&); 00197 // Rotates reference frame from Uz to newUz (unit vector) (Geant4). 00198 00199 Hep3Vector & rotate(double, const Hep3Vector &); 00200 // Rotates around the axis specified by another Hep3Vector. 00201 // (Uses methods of HepRotation, forcing linking in of Rotation.cc.) 00202 00203 Hep3Vector & operator *= (const HepRotation &); 00204 Hep3Vector & transform(const HepRotation &); 00205 // Transformation with a Rotation matrix. 00206 00207 00208 // = = = = = = = = = = = = = = = = = = = = = = = = 00209 // 00210 // Esoteric properties and operations on 3-vectors: 00211 // 00212 // 1 - Set vectors in various coordinate systems 00213 // 2 - Synonyms for accessing coordinates and properties 00214 // 3 - Comparisions (dictionary, near-ness, and geometric) 00215 // 4 - Intrinsic properties 00216 // 5 - Properties releative to z axis and arbitrary directions 00217 // 6 - Polar and azimuthal angle decomposition and deltaPhi 00218 // 7 - Rotations 00219 // 00220 // = = = = = = = = = = = = = = = = = = = = = = = = 00221 00222 // 1 - Set vectors in various coordinate systems 00223 00224 inline void setRThetaPhi (double r, double theta, double phi); 00225 // Set in spherical coordinates: Angles are measured in RADIANS 00226 00227 inline void setREtaPhi ( double r, double eta, double phi ); 00228 // Set in spherical coordinates, but specify peudorapidiy to determine theta. 00229 00230 inline void setRhoPhiZ (double rho, double phi, double z); 00231 // Set in cylindrical coordinates: Phi angle is measured in RADIANS 00232 00233 void setRhoPhiTheta ( double rho, double phi, double theta); 00234 // Set in cylindrical coordinates, but specify theta to determine z. 00235 00236 void setRhoPhiEta ( double rho, double phi, double eta); 00237 // Set in cylindrical coordinates, but specify pseudorapidity to determine z. 00238 00239 // 2 - Synonyms for accessing coordinates and properties 00240 00241 inline double getX() const; 00242 inline double getY() const; 00243 inline double getZ() const; 00244 // x(), y(), and z() 00245 00246 inline double getR () const; 00247 inline double getTheta() const; 00248 inline double getPhi () const; 00249 // mag(), theta(), and phi() 00250 00251 inline double r () const; 00252 // mag() 00253 00254 inline double rho () const; 00255 inline double getRho () const; 00256 // perp() 00257 00258 double eta () const; 00259 double getEta () const; 00260 // pseudoRapidity() 00261 00262 inline void setR ( double s ); 00263 // setMag() 00264 00265 inline void setRho ( double s ); 00266 // setPerp() 00267 00268 // 3 - Comparisions (dictionary, near-ness, and geometric) 00269 00270 int compare (const Hep3Vector & v) const; 00271 bool operator > (const Hep3Vector & v) const; 00272 bool operator < (const Hep3Vector & v) const; 00273 bool operator>= (const Hep3Vector & v) const; 00274 bool operator<= (const Hep3Vector & v) const; 00275 // dictionary ordering according to z, then y, then x component 00276 00277 inline double diff2 (const Hep3Vector & v) const; 00278 // |v1-v2|**2 00279 00280 static double setTolerance (double tol); 00281 static inline double getTolerance (); 00282 // Set the tolerance used in isNear() for Hep3Vectors 00283 00284 bool isParallel (const Hep3Vector & v, double epsilon=tolerance) const; 00285 // Are the vectors parallel, within the given tolerance? 00286 00287 bool isOrthogonal (const Hep3Vector & v, double epsilon=tolerance) const; 00288 // Are the vectors orthogonal, within the given tolerance? 00289 00290 double howParallel (const Hep3Vector & v) const; 00291 // | v1.cross(v2) / v1.dot(v2) |, to a maximum of 1. 00292 00293 double howOrthogonal (const Hep3Vector & v) const; 00294 // | v1.dot(v2) / v1.cross(v2) |, to a maximum of 1. 00295 00296 enum { ToleranceTicks = 100 }; 00297 00298 // 4 - Intrinsic properties 00299 00300 double beta () const; 00301 // relativistic beta (considering v as a velocity vector with c=1) 00302 // Same as mag() but will object if >= 1 00303 00304 double gamma() const; 00305 // relativistic gamma (considering v as a velocity vector with c=1) 00306 00307 double coLinearRapidity() const; 00308 // inverse tanh (beta) 00309 00310 // 5 - Properties relative to Z axis and to an arbitrary direction 00311 00312 // Note that the non-esoteric CLHEP provides 00313 // theta(), cosTheta(), cos2Theta, and angle(const Hep3Vector&) 00314 00315 inline double angle() const; 00316 // angle against the Z axis -- synonym for theta() 00317 00318 inline double theta(const Hep3Vector & v2) const; 00319 // synonym for angle(v2) 00320 00321 double cosTheta (const Hep3Vector & v2) const; 00322 double cos2Theta(const Hep3Vector & v2) const; 00323 // cos and cos^2 of the angle between two vectors 00324 00325 inline Hep3Vector project () const; 00326 Hep3Vector project (const Hep3Vector & v2) const; 00327 // projection of a vector along a direction. 00328 00329 inline Hep3Vector perpPart() const; 00330 inline Hep3Vector perpPart (const Hep3Vector & v2) const; 00331 // vector minus its projection along a direction. 00332 00333 double rapidity () const; 00334 // inverse tanh(v.z()) 00335 00336 double rapidity (const Hep3Vector & v2) const; 00337 // rapidity with respect to specified direction: 00338 // inverse tanh (v.dot(u)) where u is a unit in the direction of v2 00339 00340 double eta(const Hep3Vector & v2) const; 00341 // - ln tan of the angle beween the vector and the ref direction. 00342 00343 // 6 - Polar and azimuthal angle decomposition and deltaPhi 00344 00345 // Decomposition of an angle within reference defined by a direction: 00346 00347 double polarAngle (const Hep3Vector & v2) const; 00348 // The reference direction is Z: the polarAngle is abs(v.theta()-v2.theta()). 00349 00350 double deltaPhi (const Hep3Vector & v2) const; 00351 // v.phi()-v2.phi(), brought into the range (-PI,PI] 00352 00353 double azimAngle (const Hep3Vector & v2) const; 00354 // The reference direction is Z: the azimAngle is the same as deltaPhi 00355 00356 double polarAngle (const Hep3Vector & v2, 00357 const Hep3Vector & ref) const; 00358 // For arbitrary reference direction, 00359 // polarAngle is abs(v.angle(ref) - v2.angle(ref)). 00360 00361 double azimAngle (const Hep3Vector & v2, 00362 const Hep3Vector & ref) const; 00363 // To compute azimangle, project v and v2 into the plane normal to 00364 // the reference direction. Then in that plane take the angle going 00365 // clockwise around the direction from projection of v to that of v2. 00366 00367 // 7 - Rotations 00368 00369 // These mehtods **DO NOT** use anything in the HepRotation class. 00370 // Thus, use of v.rotate(axis,delta) does not force linking in Rotation.cc. 00371 00372 Hep3Vector & rotate (const Hep3Vector & axis, double delta); 00373 // Synonym for rotate (delta, axis) 00374 00375 Hep3Vector & rotate (const HepAxisAngle & ax); 00376 // HepAxisAngle is a struct holding an axis direction and an angle. 00377 00378 Hep3Vector & rotate (const HepEulerAngles & e); 00379 Hep3Vector & rotate (double phi, 00380 double theta, 00381 double psi); 00382 // Rotate via Euler Angles. Our Euler Angles conventions are 00383 // those of Goldstein Classical Mechanics page 107. 00384 00385 protected: 00386 void setSpherical (double r, double theta, double phi); 00387 void setCylindrical (double r, double phi, double z); 00388 double negativeInfinity() const; 00389 00390 protected: 00391 00392 double dx; 00393 double dy; 00394 double dz; 00395 // The components. 00396 00397 static double tolerance; 00398 // default tolerance criterion for isNear() to return true. 00399 }; // Hep3Vector 00400 00401 // Global Methods 00402 00403 Hep3Vector rotationXOf (const Hep3Vector & vec, double delta); 00404 Hep3Vector rotationYOf (const Hep3Vector & vec, double delta); 00405 Hep3Vector rotationZOf (const Hep3Vector & vec, double delta); 00406 00407 Hep3Vector rotationOf (const Hep3Vector & vec, 00408 const Hep3Vector & axis, double delta); 00409 Hep3Vector rotationOf (const Hep3Vector & vec, const HepAxisAngle & ax); 00410 00411 Hep3Vector rotationOf (const Hep3Vector & vec, 00412 double phi, double theta, double psi); 00413 Hep3Vector rotationOf (const Hep3Vector & vec, const HepEulerAngles & e); 00414 // Return a new vector based on a rotation of the supplied vector 00415 00416 std::ostream & operator << (std::ostream &, const Hep3Vector &); 00417 // Output to a stream. 00418 00419 std::istream & operator >> (std::istream &, Hep3Vector &); 00420 // Input from a stream. 00421 00422 extern const Hep3Vector HepXHat, HepYHat, HepZHat; 00423 00424 typedef Hep3Vector HepThreeVectorD; 00425 typedef Hep3Vector HepThreeVectorF; 00426 00427 Hep3Vector operator / (const Hep3Vector &, double a); 00428 // Division of 3-vectors by non-zero real number 00429 00430 inline Hep3Vector operator + (const Hep3Vector &, const Hep3Vector &); 00431 // Addition of 3-vectors. 00432 00433 inline Hep3Vector operator - (const Hep3Vector &, const Hep3Vector &); 00434 // Subtraction of 3-vectors. 00435 00436 inline double operator * (const Hep3Vector &, const Hep3Vector &); 00437 // double product of 3-vectors. 00438 00439 inline Hep3Vector operator * (const Hep3Vector &, double a); 00440 inline Hep3Vector operator * (double a, const Hep3Vector &); 00441 // Scaling of 3-vectors with a real number 00442 00443 } // namespace CLHEP 00444 00445 #include "CLHEP/Vector/ThreeVector.icc" 00446 00447 #ifdef ENABLE_BACKWARDS_COMPATIBILITY 00448 // backwards compatibility will be enabled ONLY in CLHEP 1.9 00449 using namespace CLHEP; 00450 #endif 00451 00452 #endif /* HEP_THREEVECTOR_H */