Logo ROOT   6.18/05
Reference Guide
EulerAngles.h
Go to the documentation of this file.
1// @(#)root/mathcore:$Id$
2// Authors: W. Brown, M. Fischler, L. Moneta 2005
3
4 /**********************************************************************
5 * *
6 * Copyright (c) 2005 , LCG ROOT MathLib Team *
7 * *
8 * *
9 **********************************************************************/
10
11// Header file for class EulerAngles
12//
13// Created by: Lorenzo Moneta at Tue May 10 17:55:10 2005
14//
15// Last update: Tue May 10 17:55:10 2005
16//
17#ifndef ROOT_Math_GenVector_EulerAngles
18#define ROOT_Math_GenVector_EulerAngles 1
19
25#include <algorithm>
26#include <cassert>
27
28namespace ROOT {
29namespace Math {
30
31
32//__________________________________________________________________________________________
33 /**
34 EulerAngles class describing rotation as three angles (Euler Angles).
35 The Euler angles definition matches that of Classical Mechanics (Goldstein).
36 It is also the same convention defined in
37 <A HREF="http://mathworld.wolfram.com/EulerAngles.html">mathworld</A>
38 and used in Mathematica and CLHEP. Note that the ROOT class TRotation defines
39 a slightly different convention.
40
41 @ingroup GenVector
42 */
44
45public:
46
47 typedef double Scalar;
48
49 /**
50 Default constructor
51 */
52 EulerAngles() : fPhi(0.0), fTheta(0.0), fPsi(0.0) { }
53
54 /**
55 Constructor from phi, theta and psi
56 */
57 EulerAngles( Scalar phi, Scalar theta, Scalar psi ) :
58 fPhi(phi), fTheta(theta), fPsi(psi)
59 {Rectify();} // Added 27 Jan. 06 JMM
60
61 /**
62 Construct given a pair of pointers or iterators defining the
63 beginning and end of an array of three Scalars, to be treated as
64 the angles phi, theta and psi.
65 */
66 template<class IT>
67 EulerAngles(IT begin, IT end) { SetComponents(begin,end); }
68
69 // The compiler-generated copy ctor, copy assignment, and dtor are OK.
70
71 /**
72 Re-adjust components place angles in canonical ranges
73 */
74 void Rectify();
75
76
77 // ======== Construction and assignement from any other rotation ==================
78
79 /**
80 Create from any other supported rotation (see gv_detail::convert )
81 */
82 template <class OtherRotation>
83 explicit EulerAngles(const OtherRotation & r) {gv_detail::convert(r,*this);}
84
85 /**
86 Assign from any other rotation (see gv_detail::convert )
87 */
88 template <class OtherRotation>
89 EulerAngles & operator=( OtherRotation const & r ) {
90 gv_detail::convert(r,*this);
91 return *this;
92 }
93
94#ifdef OLD
95 explicit EulerAngles(const Rotation3D & r) {gv_detail::convert(r,*this);}
96
97 /**
98 Construct from a rotation matrix
99 */
100 explicit EulerAngles(const Rotation3D & r) {gv_detail::convert(r,*this);}
101
102 /**
103 Construct from a rotation represented by a Quaternion
104 */
105 explicit EulerAngles(const Quaternion & q) {gv_detail::convert(q,*this);}
106
107 /**
108 Construct from an AxisAngle
109 */
110 explicit EulerAngles(const AxisAngle & a ) { gv_detail::convert(a, *this); }
111
112 /**
113 Construct from an axial rotation
114 */
115 explicit EulerAngles( RotationZ const & r ) { gv_detail::convert(r, *this); }
116 explicit EulerAngles( RotationY const & r ) { gv_detail::convert(r, *this); }
117 explicit EulerAngles( RotationX const & r ) { gv_detail::convert(r, *this); }
118
119
120 /**
121 Assign from an AxisAngle
122 */
124 operator=( AxisAngle const & a ) { return operator=(EulerAngles(a)); }
125
126 /**
127 Assign from a Quaternion
128 */
130 operator=( Quaternion const & q ) {return operator=(EulerAngles(q)); }
131
132 /**
133 Assign from an axial rotation
134 */
136 operator=( RotationZ const & r ) { return operator=(EulerAngles(r)); }
138 operator=( RotationY const & r ) { return operator=(EulerAngles(r)); }
140 operator=( RotationX const & r ) { return operator=(EulerAngles(r)); }
141
142#endif
143
144 // ======== Components ==============
145
146 /**
147 Set the three Euler angles given a pair of pointers or iterators
148 defining the beginning and end of an array of three Scalars.
149 */
150 template<class IT>
151#ifndef NDEBUG
152 void SetComponents(IT begin, IT end) {
153#else
154 void SetComponents(IT begin, IT ) {
155#endif
156 fPhi = *begin++;
157 fTheta = *begin++;
158 fPsi = *begin++;
159 assert(begin == end);
160 Rectify(); // Added 27 Jan. 06 JMM
161 }
162
163 /**
164 Get the axis and then the angle into data specified by an iterator begin
165 and another to the end of the desired data (4 past start).
166 */
167 template<class IT>
168#ifndef NDEBUG
169 void GetComponents(IT begin, IT end) const {
170#else
171 void GetComponents(IT begin, IT ) const {
172#endif
173 *begin++ = fPhi;
174 *begin++ = fTheta;
175 *begin++ = fPsi;
176 assert(begin == end);
177 }
178
179 /**
180 Get the axis and then the angle into data specified by an iterator begin
181 */
182 template<class IT>
183 void GetComponents(IT begin) const {
184 *begin++ = fPhi;
185 *begin++ = fTheta;
186 *begin = fPsi;
187 }
188
189 /**
190 Set the components phi, theta, psi based on three Scalars.
191 */
192 void SetComponents(Scalar phi, Scalar theta, Scalar psi) {
193 fPhi=phi; fTheta=theta; fPsi=psi;
194 Rectify(); // Added 27 Jan. 06 JMM
195 }
196
197 /**
198 Get the components phi, theta, psi into three Scalars.
199 */
200 void GetComponents(Scalar & phi, Scalar & theta, Scalar & psi) const {
201 phi=fPhi; theta=fTheta; psi=fPsi;
202 }
203
204 /**
205 Set Phi Euler angle // JMM 30 Jan. 2006
206 */
207 void SetPhi(Scalar phi) { fPhi=phi; Rectify(); }
208
209 /**
210 Return Phi Euler angle
211 */
212 Scalar Phi() const { return fPhi; }
213
214 /**
215 Set Theta Euler angle // JMM 30 Jan. 2006
216 */
217 void SetTheta(Scalar theta) { fTheta=theta; Rectify(); }
218
219 /**
220 Return Theta Euler angle
221 */
222 Scalar Theta() const { return fTheta; }
223
224 /**
225 Set Psi Euler angle // JMM 30 Jan. 2006
226 */
227 void SetPsi(Scalar psi) { fPsi=psi; Rectify(); }
228
229 /**
230 Return Psi Euler angle
231 */
232 Scalar Psi() const { return fPsi; }
233
234 // =========== operations ==============
235
236
237 /**
238 Rotation operation on a displacement vector in any coordinate system and tag
239 */
240 template <class CoordSystem, class U>
243 return Rotation3D(*this) ( v );
244 }
245
246 /**
247 Rotation operation on a position vector in any coordinate system
248 */
249 template <class CoordSystem, class U>
254 return PositionVector3D<CoordSystem,U> ( rxyz );
255 }
256
257 /**
258 Rotation operation on a Lorentz vector in any 4D coordinate system
259 */
260 template <class CoordSystem>
264 xyz = operator()(xyz);
265 LorentzVector< PxPyPzE4D<double> > xyzt (xyz.X(), xyz.Y(), xyz.Z(), v.E());
266 return LorentzVector<CoordSystem> ( xyzt );
267 }
268
269 /**
270 Rotation operation on an arbitrary vector v.
271 Preconditions: v must implement methods x(), y(), and z()
272 and the arbitrary vector type must have a constructor taking (x,y,z)
273 */
274 template <class ForeignVector>
275 ForeignVector
276 operator() (const ForeignVector & v) const {
279 return ForeignVector ( rxyz.X(), rxyz.Y(), rxyz.Z() );
280 }
281
282 /**
283 Overload operator * for rotation on a vector
284 */
285 template <class AVector>
286 inline
287 AVector operator* (const AVector & v) const
288 {
289 return operator()(v);
290 }
291
292 /**
293 Invert a rotation in place
294 */
295 // theta stays the same and negative rotation in Theta is done via a rotation
296 // of + PI in phi and Psi
297 void Invert() {
298 Scalar tmp = -fPhi;
299 fPhi = -fPsi + Pi();
300 fPsi=tmp + Pi();
301 }
302
303 /**
304 Return inverse of a rotation
305 */
306 EulerAngles Inverse() const { return EulerAngles(-fPsi + Pi(), fTheta, -fPhi + Pi()); }
307
308 // ========= Multi-Rotation Operations ===============
309
310 /**
311 Multiply (combine) two rotations
312 */
313 EulerAngles operator * (const Rotation3D & r) const;
314 EulerAngles operator * (const AxisAngle & a) const;
315 EulerAngles operator * (const EulerAngles & e) const;
316 EulerAngles operator * (const Quaternion & q) const;
317 EulerAngles operator * (const RotationX & rx) const;
318 EulerAngles operator * (const RotationY & ry) const;
319 EulerAngles operator * (const RotationZ & rz) const;
320
321 /**
322 Post-Multiply (on right) by another rotation : T = T*R
323 */
324 template <class R>
325 EulerAngles & operator *= (const R & r) { return *this = (*this)*r; }
326
327 /**
328 Distance between two rotations
329 */
330 template <class R>
331 Scalar Distance ( const R & r ) const {return gv_detail::dist(*this,r);}
332
333 /**
334 Equality/inequality operators
335 */
336 bool operator == (const EulerAngles & rhs) const {
337 if( fPhi != rhs.fPhi ) return false;
338 if( fTheta != rhs.fTheta ) return false;
339 if( fPsi != rhs.fPsi ) return false;
340 return true;
341 }
342 bool operator != (const EulerAngles & rhs) const {
343 return ! operator==(rhs);
344 }
345
346private:
347
348 double fPhi; // Z rotation angle (first) defined in [-PI,PI]
349 double fTheta; // X rotation angle (second) defined only [0,PI]
350 double fPsi; // Z rotation angle (third) defined in [-PI,PI]
351
352 static double Pi() { return M_PI; }
353
354}; // EulerAngles
355
356/**
357 Distance between two rotations
358 */
359template <class R>
360inline
361typename EulerAngles::Scalar
362Distance ( const EulerAngles& r1, const R & r2) {return gv_detail::dist(r1,r2);}
363
364/**
365 Multiplication of an axial rotation by an AxisAngle
366 */
367EulerAngles operator* (RotationX const & r1, EulerAngles const & r2);
368EulerAngles operator* (RotationY const & r1, EulerAngles const & r2);
369EulerAngles operator* (RotationZ const & r1, EulerAngles const & r2);
370
371/**
372 Stream Output and Input
373 */
374 // TODO - I/O should be put in the manipulator form
375
376std::ostream & operator<< (std::ostream & os, const EulerAngles & e);
377
378} // namespace Math
379} // namespace ROOT
380
381
382#endif /* ROOT_Math_GenVector_EulerAngles */
SVector< double, 2 > v
Definition: Dict.h:5
ROOT::R::TRInterface & r
Definition: Object.C:4
#define R(a, b, c, d, e, f, g, h, i)
Definition: RSha256.hxx:110
#define e(i)
Definition: RSha256.hxx:103
#define M_PI
Definition: Rotated.cxx:105
float * q
Definition: THbookFile.cxx:87
AxisAngle class describing rotation represented with direction axis (3D Vector) and an angle of rotat...
Definition: AxisAngle.h:41
Class describing a generic displacement vector in 3 dimensions.
Scalar X() const
Cartesian X, converting if necessary from internal coordinate system.
Scalar Y() const
Cartesian Y, converting if necessary from internal coordinate system.
Scalar Z() const
Cartesian Z, converting if necessary from internal coordinate system.
EulerAngles class describing rotation as three angles (Euler Angles).
Definition: EulerAngles.h:43
Scalar Psi() const
Return Psi Euler angle.
Definition: EulerAngles.h:232
bool operator!=(const EulerAngles &rhs) const
Definition: EulerAngles.h:342
Scalar Theta() const
Return Theta Euler angle.
Definition: EulerAngles.h:222
void SetPsi(Scalar psi)
Set Psi Euler angle // JMM 30 Jan.
Definition: EulerAngles.h:227
void SetTheta(Scalar theta)
Set Theta Euler angle // JMM 30 Jan.
Definition: EulerAngles.h:217
void GetComponents(IT begin, IT end) const
Get the axis and then the angle into data specified by an iterator begin and another to the end of th...
Definition: EulerAngles.h:169
void GetComponents(Scalar &phi, Scalar &theta, Scalar &psi) const
Get the components phi, theta, psi into three Scalars.
Definition: EulerAngles.h:200
EulerAngles()
Default constructor.
Definition: EulerAngles.h:52
void SetComponents(Scalar phi, Scalar theta, Scalar psi)
Set the components phi, theta, psi based on three Scalars.
Definition: EulerAngles.h:192
void GetComponents(IT begin) const
Get the axis and then the angle into data specified by an iterator begin.
Definition: EulerAngles.h:183
void SetPhi(Scalar phi)
Set Phi Euler angle // JMM 30 Jan.
Definition: EulerAngles.h:207
EulerAngles & operator*=(const R &r)
Post-Multiply (on right) by another rotation : T = T*R.
Definition: EulerAngles.h:325
EulerAngles & operator=(OtherRotation const &r)
Assign from any other rotation (see gv_detail::convert )
Definition: EulerAngles.h:89
AVector operator*(const AVector &v) const
Overload operator * for rotation on a vector.
Definition: EulerAngles.h:287
EulerAngles(const OtherRotation &r)
Create from any other supported rotation (see gv_detail::convert )
Definition: EulerAngles.h:83
void Invert()
Invert a rotation in place.
Definition: EulerAngles.h:297
EulerAngles(IT begin, IT end)
Construct given a pair of pointers or iterators defining the beginning and end of an array of three S...
Definition: EulerAngles.h:67
bool operator==(const EulerAngles &rhs) const
Equality/inequality operators.
Definition: EulerAngles.h:336
EulerAngles Inverse() const
Return inverse of a rotation.
Definition: EulerAngles.h:306
static double Pi()
Definition: EulerAngles.h:352
void Rectify()
Re-adjust components place angles in canonical ranges.
Definition: EulerAngles.cxx:38
EulerAngles(Scalar phi, Scalar theta, Scalar psi)
Constructor from phi, theta and psi.
Definition: EulerAngles.h:57
DisplacementVector3D< CoordSystem, U > operator()(const DisplacementVector3D< CoordSystem, U > &v) const
Rotation operation on a displacement vector in any coordinate system and tag.
Definition: EulerAngles.h:242
Scalar Distance(const R &r) const
Distance between two rotations.
Definition: EulerAngles.h:331
Scalar Phi() const
Return Phi Euler angle.
Definition: EulerAngles.h:212
void SetComponents(IT begin, IT end)
Set the three Euler angles given a pair of pointers or iterators defining the beginning and end of an...
Definition: EulerAngles.h:152
Class describing a generic position vector (point) in 3 dimensions.
Rotation class with the (3D) rotation represented by a unit quaternion (u, i, j, k).
Definition: Quaternion.h:47
Rotation class with the (3D) rotation represented by a 3x3 orthogonal matrix.
Definition: Rotation3D.h:65
Rotation class representing a 3D rotation about the X axis by the angle of rotation.
Definition: RotationX.h:43
Rotation class representing a 3D rotation about the Y axis by the angle of rotation.
Definition: RotationY.h:43
Rotation class representing a 3D rotation about the Z axis by the angle of rotation.
Definition: RotationZ.h:43
Namespace for new Math classes and functions.
double dist(Rotation3D const &r1, Rotation3D const &r2)
Definition: 3DDistances.cxx:48
void convert(R1 const &, R2 const)
Definition: 3DConversions.h:41
std::ostream & operator<<(std::ostream &os, const AxisAngle &a)
Stream Output and Input.
Definition: AxisAngle.cxx:91
Namespace for new ROOT classes and functions.
Definition: StringConv.hxx:21
auto * a
Definition: textangle.C:12