CLHEP VERSION Reference Documentation
   
CLHEP Home Page     CLHEP Documentation     CLHEP Bug Reports

Vector/Vector/Rotation.h
Go to the documentation of this file.
1// -*- C++ -*-
2// CLASSDOC OFF
3// $Id: Rotation.h,v 1.3 2003/10/23 21:29:52 garren Exp $
4// ---------------------------------------------------------------------------
5// CLASSDOC ON
6//
7// This file is a part of the CLHEP - a Class Library for High Energy Physics.
8//
9// This is the definition of the HepRotation class for performing rotations
10// on objects of the Hep3Vector (and HepLorentzVector) class.
11//
12// HepRotation is a concrete implementation of Hep3RotationInterface.
13//
14// .SS See Also
15// RotationInterfaces.h
16// ThreeVector.h, LorentzVector.h, LorentzRotation.h
17//
18// .SS Author
19// Leif Lonnblad, Mark Fischler
20
21#ifndef HEP_ROTATION_H
22#define HEP_ROTATION_H
23
24#ifdef GNUPRAGMA
25#pragma interface
26#endif
27
28#include "CLHEP/Vector/defs.h"
29#include "CLHEP/Vector/RotationInterfaces.h"
30#include "CLHEP/Vector/RotationX.h"
31#include "CLHEP/Vector/RotationY.h"
32#include "CLHEP/Vector/RotationZ.h"
33#include "CLHEP/Vector/LorentzVector.h"
34
35namespace CLHEP {
36
37// Declarations of classes and global methods
38class HepRotation;
39inline HepRotation inverseOf ( const HepRotation & r );
40inline HepRotation operator * (const HepRotationX & rx, const HepRotation & r);
41inline HepRotation operator * (const HepRotationY & ry, const HepRotation & r);
42inline HepRotation operator * (const HepRotationZ & rz, const HepRotation & r);
43
48class HepRotation {
49
50public:
51
52 // ---------- Constructors and Assignment:
53
54 inline HepRotation();
55 // Default constructor. Gives a unit matrix.
56
57 inline HepRotation(const HepRotation & m);
58 // Copy constructor.
59
60 inline HepRotation(const HepRotationX & m);
61 inline HepRotation(const HepRotationY & m);
62 inline HepRotation(const HepRotationZ & m);
63 // Construct from specialized rotation.
64
65 HepRotation & set( const Hep3Vector & axis, double delta );
66 HepRotation ( const Hep3Vector & axis, double delta );
67 // Construct from axis and angle.
68
69 HepRotation & set( const HepAxisAngle & ax );
70 HepRotation ( const HepAxisAngle & ax );
71 // Construct from AxisAngle structure.
72
73 HepRotation & set( double phi, double theta, double psi );
74 HepRotation ( double phi, double theta, double psi );
75 // Construct from three Euler angles (in radians).
76
79 // Construct from EulerAngles structure.
80
82 const Hep3Vector & colY,
83 const Hep3Vector & colZ );
84 // Construct from three *orthogonal* unit vector columns.
85 // NOTE:
86 // This constructor, and the two set methods below,
87 // will check that the columns (or rows) form an orthonormal
88 // matrix, and will adjust values so that this relation is
89 // as exact as possible.
90
92 const Hep3Vector & colY,
93 const Hep3Vector & colZ );
94 // supply three *orthogonal* unit vectors for the columns.
95
97 const Hep3Vector & rowY,
98 const Hep3Vector & rowZ );
99 // supply three *orthogonal* unit vectors for the rows.
100
101 inline HepRotation & set(const HepRotationX & r);
102 inline HepRotation & set(const HepRotationY & r);
103 inline HepRotation & set(const HepRotationZ & r);
104 // set from specialized rotation.
105
107 // Assignment.
108
112 // Assignment from specialized rotation.
113
114 inline HepRotation &set( const HepRep3x3 & m );
115 inline HepRotation ( const HepRep3x3 & m );
116 // WARNING - NO CHECKING IS DONE!
117 // Constructon directly from from a 3x3 representation,
118 // which is required to be an orthogonal matrix.
119
120 inline ~HepRotation();
121 // Trivial destructor.
122
123 // ---------- Accessors:
124
125 inline Hep3Vector colX() const;
126 inline Hep3Vector colY() const;
127 inline Hep3Vector colZ() const;
128 // orthogonal unit-length column vectors
129
130 inline Hep3Vector rowX() const;
131 inline Hep3Vector rowY() const;
132 inline Hep3Vector rowZ() const;
133 // orthogonal unit-length row vectors
134
135 inline double xx() const;
136 inline double xy() const;
137 inline double xz() const;
138 inline double yx() const;
139 inline double yy() const;
140 inline double yz() const;
141 inline double zx() const;
142 inline double zy() const;
143 inline double zz() const;
144 // Elements of the rotation matrix (Geant4).
145
146 inline HepRep3x3 rep3x3() const;
147 // 3x3 representation:
148
149 // ------------ Subscripting:
150
151 class HepRotation_row {
152 public:
153 inline HepRotation_row(const HepRotation &, int);
154 inline double operator [] (int) const;
155 private:
156 const HepRotation & rr;
157 int ii;
158 };
159 // Helper class for implemention of C-style subscripting r[i][j]
160
161 inline const HepRotation_row operator [] (int) const;
162 // Returns object of the helper class for C-style subscripting r[i][j]
163 // i and j range from 0 to 2.
164
165 double operator () (int, int) const;
166 // Fortran-style subscripting: returns (i,j) element of the rotation matrix.
167 // Note: i and j still range from 0 to 2. [Rotation.cc]
168
169 // ------------ Euler angles:
170 inline double getPhi () const;
171 inline double getTheta() const;
172 inline double getPsi () const;
173 double phi () const;
174 double theta() const;
175 double psi () const;
177
178 // ------------ axis & angle of rotation:
179 inline double getDelta() const;
180 inline Hep3Vector getAxis () const;
181 double delta() const;
182 Hep3Vector axis () const;
184 void getAngleAxis(double & delta, Hep3Vector & axis) const;
185 // Returns the rotation angle and rotation axis (Geant4). [Rotation.cc]
186
187 // ------------- Angles of rotated axes
188 double phiX() const;
189 double phiY() const;
190 double phiZ() const;
191 double thetaX() const;
192 double thetaY() const;
193 double thetaZ() const;
194 // Return angles (RADS) made by rotated axes against original axes (Geant4).
195 // [Rotation.cc]
196
197 // ---------- Other accessors treating pure rotation as a 4-rotation
198
199 inline HepLorentzVector col1() const;
200 inline HepLorentzVector col2() const;
201 inline HepLorentzVector col3() const;
202 // orthosymplectic 4-vector columns - T component will be zero
203
204 inline HepLorentzVector col4() const;
205 // Will be (0,0,0,1) for this pure Rotation.
206
207 inline HepLorentzVector row1() const;
208 inline HepLorentzVector row2() const;
209 inline HepLorentzVector row3() const;
210 // orthosymplectic 4-vector rows - T component will be zero
211
212 inline HepLorentzVector row4() const;
213 // Will be (0,0,0,1) for this pure Rotation.
214
215 inline double xt() const;
216 inline double yt() const;
217 inline double zt() const;
218 inline double tx() const;
219 inline double ty() const;
220 inline double tz() const;
221 // Will be zero for this pure Rotation
222
223 inline double tt() const;
224 // Will be one for this pure Rotation
225
226 inline HepRep4x4 rep4x4() const;
227 // 4x4 representation.
228
229 // --------- Mutators
230
231 void setPhi (double phi);
232 // change Euler angle phi, leaving theta and psi unchanged.
233
234 void setTheta (double theta);
235 // change Euler angle theta, leaving phi and psi unchanged.
236
237 void setPsi (double psi);
238 // change Euler angle psi, leaving theta and phi unchanged.
239
240 void setAxis (const Hep3Vector & axis);
241 // change rotation axis, leaving delta unchanged.
242
243 void setDelta (double delta);
244 // change angle of rotation, leaving rotation axis unchanged.
245
246 // ---------- Decomposition:
247
248 void decompose (HepAxisAngle & rotation, Hep3Vector & boost) const;
249 void decompose (Hep3Vector & boost, HepAxisAngle & rotation) const;
250 // These are trivial, as the boost vector is 0. [RotationP.cc]
251
252 // ---------- Comparisons:
253
254 bool isIdentity() const;
255 // Returns true if the identity matrix (Geant4). [Rotation.cc]
256
257 int compare( const HepRotation & r ) const;
258 // Dictionary-order comparison, in order zz, zy, zx, yz, ... xx
259 // Used in operator<, >, <=, >=
260
261 inline bool operator== ( const HepRotation & r ) const;
262 inline bool operator!= ( const HepRotation & r ) const;
263 inline bool operator< ( const HepRotation & r ) const;
264 inline bool operator> ( const HepRotation & r ) const;
265 inline bool operator<= ( const HepRotation & r ) const;
266 inline bool operator>= ( const HepRotation & r ) const;
267
268 double distance2( const HepRotation & r ) const;
269 // 3 - Tr ( this/r ) -- This works with RotationX, Y or Z also
270
271 double howNear( const HepRotation & r ) const;
272 bool isNear( const HepRotation & r,
273 double epsilon=Hep4RotationInterface::tolerance) const;
274
275 double distance2( const HepBoost & lt ) const;
276 // 3 - Tr ( this ) + |b|^2 / (1-|b|^2)
277 double distance2( const HepLorentzRotation & lt ) const;
278 // 3 - Tr ( this/r ) + |b|^2 / (1-|b|^2) where b is the boost vector of lt
279
280 double howNear( const HepBoost & lt ) const;
281 double howNear( const HepLorentzRotation & lt ) const;
282 bool isNear( const HepBoost & lt,
283 double epsilon=Hep4RotationInterface::tolerance) const;
284 bool isNear( const HepLorentzRotation & lt,
285 double epsilon=Hep4RotationInterface::tolerance) const;
286
287 // ---------- Properties:
288
289 double norm2() const;
290 // distance2 (IDENTITY), which is 3 - Tr ( *this )
291
292 void rectify();
293 // non-const but logically moot correction for accumulated roundoff errors
294 // rectify averages the matrix with the transpose of its actual
295 // inverse (absent accumulated roundoff errors, the transpose IS
296 // the inverse)); this removes to first order those errors.
297 // Then it formally extracts axis and delta, and forms a true
298 // HepRotation with those values of axis and delta.
299
300 // ---------- Application:
301
302 inline Hep3Vector operator() (const Hep3Vector & p) const;
303 // Rotate a Hep3Vector.
304
305 inline Hep3Vector operator * (const Hep3Vector & p) const;
306 // Multiplication with a Hep3Vector.
307
309 // Rotate (the space part of) a HepLorentzVector.
310
312 // Multiplication with a HepLorentzVector.
313
314 // ---------- Operations in the group of Rotations
315
316 inline HepRotation operator * (const HepRotation & r) const;
317 // Product of two rotations (this) * r - matrix multiplication
318
319 inline HepRotation operator * (const HepRotationX & rx) const;
320 inline HepRotation operator * (const HepRotationY & ry) const;
321 inline HepRotation operator * (const HepRotationZ & rz) const;
322 // Product of two rotations (this) * r - faster when specialized type
323
325 inline HepRotation & transform (const HepRotation & r);
326 // Matrix multiplication.
327 // Note a *= b; <=> a = a * b; while a.transform(b); <=> a = b * a;
328
332 inline HepRotation & transform (const HepRotationX & r);
333 inline HepRotation & transform (const HepRotationY & r);
334 inline HepRotation & transform (const HepRotationZ & r);
335 // Matrix multiplication by specialized matrices
336
338 // Rotation around the x-axis; equivalent to R = RotationX(delta) * R
339
341 // Rotation around the y-axis; equivalent to R = RotationY(delta) * R
342
344 // Rotation around the z-axis; equivalent to R = RotationZ(delta) * R
345
347 inline HepRotation & rotate(double delta, const Hep3Vector * axis);
348 // Rotation around a specified vector.
349 // r.rotate(d,a) is equivalent to r = Rotation(d,a) * r
350
352 const Hep3Vector & newY,
353 const Hep3Vector & newZ);
354 // Rotation of local axes defined by 3 orthonormal vectors (Geant4).
355 // Equivalent to r = Rotation (newX, newY, newZ) * r
356
357 inline HepRotation inverse() const;
358 // Returns the inverse.
359
360 inline HepRotation & invert();
361 // Inverts the Rotation matrix.
362
363 // ---------- I/O:
364
365 std::ostream & print( std::ostream & os ) const;
366 // Aligned six-digit-accurate output of the rotation matrix. [RotationIO.cc]
367
368 // ---------- Identity Rotation:
369
370 static const HepRotation IDENTITY;
371
372 // ---------- Tolerance
373
374 static inline double getTolerance();
375 static inline double setTolerance(double tol);
376
377protected:
378
379 inline HepRotation(double mxx, double mxy, double mxz,
380 double myx, double myy, double myz,
381 double mzx, double mzy, double mzz);
382 // Protected constructor.
383 // DOES NOT CHECK FOR VALIDITY AS A ROTATION.
384
385 friend HepRotation operator* (const HepRotationX & rx, const HepRotation & r);
386 friend HepRotation operator* (const HepRotationY & ry, const HepRotation & r);
387 friend HepRotation operator* (const HepRotationZ & rz, const HepRotation & r);
388
389 double rxx, rxy, rxz,
390 ryx, ryy, ryz,
391 rzx, rzy, rzz;
392 // The matrix elements.
393
394private:
395 bool
396 setCols ( const Hep3Vector & u1, // Vectors assume to be of unit length
397 const Hep3Vector & u2,
398 const Hep3Vector & u3,
399 double u1u2,
400 Hep3Vector & v1, // Returned vectors
401 Hep3Vector & v2,
402 Hep3Vector & v3 ) const;
403 void setArbitrarily (const Hep3Vector & colX, // assumed to be of unit length
404 Hep3Vector & v1,
405 Hep3Vector & v2,
406 Hep3Vector & v3) const;
407}; // HepRotation
408
409inline
410std::ostream & operator <<
411 ( std::ostream & os, const HepRotation & r ) {return r.print(os);}
412
413} // namespace CLHEP
414
415#include "CLHEP/Vector/Rotation.icc"
416
417#ifdef ENABLE_BACKWARDS_COMPATIBILITY
418// backwards compatibility will be enabled ONLY in CLHEP 1.9
419using namespace CLHEP;
420#endif
421
422#endif /* HEP_ROTATION_H */
423
HepRotation_row(const HepRotation &, int)
HepAxisAngle axisAngle() const
double operator()(int, int) const
Definition Rotation.cc:29
double getPsi() const
double getPhi() const
static double getTolerance()
double zz() const
HepEulerAngles eulerAngles() const
HepRotation(double phi, double theta, double psi)
double yz() const
Hep3Vector colX() const
HepRotation & rotate(double delta, const Hep3Vector *axis)
static const HepRotation IDENTITY
bool operator==(const HepRotation &r) const
bool operator!=(const HepRotation &r) const
HepLorentzVector col3() const
HepRotation & set(double phi, double theta, double psi)
HepRotation & rotateAxes(const Hep3Vector &newX, const Hep3Vector &newY, const Hep3Vector &newZ)
double zx() const
void setPsi(double psi)
std::ostream & print(std::ostream &os) const
HepRotation(double mxx, double mxy, double mxz, double myx, double myy, double myz, double mzx, double mzy, double mzz)
Hep3Vector axis() const
Definition RotationA.cc:82
HepRotation & setRows(const Hep3Vector &rowX, const Hep3Vector &rowY, const Hep3Vector &rowZ)
double phi() const
Definition RotationE.cc:73
HepRotation & set(const HepRotationY &r)
HepRotation(const Hep3Vector &axis, double delta)
double thetaY() const
HepRotation inverse() const
double tz() const
HepRotation & set(const HepEulerAngles &e)
double distance2(const HepLorentzRotation &lt) const
HepRotation(const HepRotationY &m)
HepRotation & rotate(double delta, const Hep3Vector &axis)
double howNear(const HepBoost &lt) const
Hep3Vector rowY() const
Hep3Vector colY() const
double howNear(const HepLorentzRotation &lt) const
bool isNear(const HepLorentzRotation &lt, double epsilon=Hep4RotationInterface::tolerance) const
void decompose(Hep3Vector &boost, HepAxisAngle &rotation) const
Hep3Vector getAxis() const
const HepRotation_row operator[](int) const
bool isIdentity() const
HepRotation & operator*=(const HepRotation &r)
double ty() const
double distance2(const HepBoost &lt) const
double distance2(const HepRotation &r) const
HepRotation & transform(const HepRotationZ &r)
HepLorentzVector col1() const
double delta() const
Definition RotationA.cc:69
HepRotation & set(const Hep3Vector &colX, const Hep3Vector &colY, const Hep3Vector &colZ)
HepLorentzVector col2() const
HepLorentzVector col4() const
HepRotation & set(const HepAxisAngle &ax)
double norm2() const
static double setTolerance(double tol)
HepRotation(const Hep3Vector &colX, const Hep3Vector &colY, const Hep3Vector &colZ)
HepRotation(const HepRotationZ &m)
double tx() const
HepRotation(const HepAxisAngle &ax)
HepRotation & transform(const HepRotation &r)
double phiY() const
double yx() const
HepRep3x3 rep3x3() const
HepRotation(const HepRotationX &m)
HepRotation & set(const HepRotationX &r)
double zy() const
HepRep4x4 rep4x4() const
HepLorentzVector row1() const
HepRotation(const HepRep3x3 &m)
bool operator>(const HepRotation &r) const
void decompose(HepAxisAngle &rotation, Hep3Vector &boost) const
double thetaX() const
double zt() const
HepRotation & set(const Hep3Vector &axis, double delta)
void getAngleAxis(double &delta, Hep3Vector &axis) const
HepRotation & operator=(const HepRotation &r)
bool isNear(const HepRotation &r, double epsilon=Hep4RotationInterface::tolerance) const
void setDelta(double delta)
Hep3Vector colZ() const
double xx() const
HepLorentzVector row2() const
double howNear(const HepRotation &r) const
HepRotation & rotateX(double delta)
bool isNear(const HepBoost &lt, double epsilon=Hep4RotationInterface::tolerance) const
double psi() const
Definition RotationE.cc:113
HepRotation & set(const HepRep3x3 &m)
HepRotation & transform(const HepRotationX &r)
Hep3Vector rowX() const
double theta() const
Definition RotationE.cc:107
double phiX() const
double tt() const
HepRotation(const HepEulerAngles &e)
HepRotation & rotateZ(double delta)
double yy() const
void setAxis(const Hep3Vector &axis)
double xz() const
void setPhi(double phi)
HepLorentzVector row3() const
double yt() const
HepRotation & rotateY(double delta)
Hep3Vector rowZ() const
double getTheta() const
bool operator<(const HepRotation &r) const
HepLorentzVector row4() const
HepLorentzVector operator()(const HepLorentzVector &w) const
double thetaZ() const
HepRotation & set(const HepRotationZ &r)
HepRotation & transform(const HepRotationY &r)
int compare(const HepRotation &r) const
HepRotation & invert()
HepRotation(const HepRotation &m)
double xy() const
friend HepRotation operator*(const HepRotationX &rx, const HepRotation &r)
bool operator>=(const HepRotation &r) const
bool operator<=(const HepRotation &r) const
double getDelta() const
void setTheta(double theta)
double phiZ() const
double xt() const
HepLorentzRotation operator*(const HepRotation &r, const HepLorentzRotation &lt)
HepBoost inverseOf(const HepBoost &lt)