1 // -*- C++ -*-
2 // CLASSDOC OFF
3 // ---------------------------------------------------------------------------
4 // CLASSDOC ON
5 //
6 // This file is a part of the CLHEP - a Class Library for High Energy Physics.
7 //
8 // This is the definition of the HepRotationZ class for performing rotations
9 // around the X axis on objects of the Hep3Vector (and HepLorentzVector) class.
10 //
11 // HepRotationZ is a concrete implementation of Hep3RotationInterface.
12 //
13 // .SS See Also
14 // RotationInterfaces.h
15 // ThreeVector.h, LorentzVector.h, LorentzRotation.h
16 //
17 // .SS Author
18 // Mark Fischler
19 
20 #ifndef HEP_ROTATIONZ_H
21 #define HEP_ROTATIONZ_H
22 
23 #ifdef GNUPRAGMA
24 #pragma interface
25 #endif
26 
27 #include "CLHEP/Vector/defs.h"
28 #include "CLHEP/Vector/RotationInterfaces.h"
29 
30 namespace CLHEP {
31 
32 class HepRotationZ;
33 class HepRotation;
34 class HepBoost;
35 
36 inline HepRotationZ inverseOf(const HepRotationZ & r);
37 // Returns the inverse of a RotationZ.
38 
39 /**
40  * @author
41  * @ingroup vector
42  */
43 class HepRotationZ {
44 
45 public:
46 
47   // ----------  Constructors and Assignment:
48 
49   inline HepRotationZ();
50   // Default constructor. Gives an identity rotation.
51 
52   HepRotationZ(double delta);
53   // supply angle of rotation
54 
55   inline HepRotationZ(const HepRotationZ & orig);
56   // Copy constructor.
57 
58   inline HepRotationZ & operator = (const HepRotationZ & r);
59   // Assignment from a Rotation, which must be RotationZ
60 
61   HepRotationZ & set ( double delta );
62   // set angle of rotation
63 
64   inline ~HepRotationZ();
65   // Trivial destructor.
66 
67   // ----------  Accessors:
68 
69   inline Hep3Vector colX() const;
70   inline Hep3Vector colY() const;
71   inline Hep3Vector colZ() const;
72   // orthogonal unit-length column vectors
73 
74   inline Hep3Vector rowX() const;
75   inline Hep3Vector rowY() const;
76   inline Hep3Vector rowZ() const;
77   // orthogonal unit-length row vectors
78 
79   inline double xx() const;
80   inline double xy() const;
81   inline double xz() const;
82   inline double yx() const;
83   inline double yy() const;
84   inline double yz() const;
85   inline double zx() const;
86   inline double zy() const;
87   inline double zz() const;
88   // Elements of the rotation matrix (Geant4).
89 
90   inline HepRep3x3 rep3x3() const;
91   //   3x3 representation:
92 
93   // ------------  Euler angles:
94   inline  double getPhi  () const;
95   inline  double getTheta() const;
96   inline  double getPsi  () const;
97   double    phi  () const;
98   double    theta() const;
99   double    psi  () const;
100   HepEulerAngles eulerAngles() const;
101 
102   // ------------  axis & angle of rotation:
103   inline  double  getDelta() const;
104   inline  Hep3Vector getAxis () const;
105   inline  double     delta() const;
106   inline  Hep3Vector    axis () const;
107   inline  HepAxisAngle  axisAngle() const;
108   inline  void getAngleAxis(double & delta, Hep3Vector & axis) const;
109   // Returns the rotation angle and rotation axis (Geant4).
110 
111   // ------------- Angles of rotated axes
112   double phiX() const;
113   double phiY() const;
114   double phiZ() const;
115   double thetaX() const;
116   double thetaY() const;
117   double thetaZ() const;
118   // Return angles (RADS) made by rotated axes against original axes (Geant4).
119 
120   // ----------  Other accessors treating pure rotation as a 4-rotation
121 
122   inline HepLorentzVector col1() const;
123   inline HepLorentzVector col2() const;
124   inline HepLorentzVector col3() const;
125   //  orthosymplectic 4-vector columns - T component will be zero
126 
127   inline HepLorentzVector col4() const;
128   // Will be (0,0,0,1) for this pure Rotation.
129 
130   inline HepLorentzVector row1() const;
131   inline HepLorentzVector row2() const;
132   inline HepLorentzVector row3() const;
133   //  orthosymplectic 4-vector rows - T component will be zero
134 
135   inline HepLorentzVector row4() const;
136   // Will be (0,0,0,1) for this pure Rotation.
137 
138   inline double xt() const;
139   inline double yt() const;
140   inline double zt() const;
141   inline double tx() const;
142   inline double ty() const;
143   inline double tz() const;
144   // Will be zero for this pure Rotation
145 
146   inline double tt() const;
147   // Will be one for this pure Rotation
148 
149   inline HepRep4x4 rep4x4() const;
150   //   4x4 representation.
151 
152   // ---------   Mutators
153 
154   void setDelta (double delta);
155   // change angle of rotation, leaving rotation axis unchanged.
156 
157   // ----------  Decomposition:
158 
159   void decompose (HepAxisAngle & rotation, Hep3Vector & boost) const;
160   void decompose (Hep3Vector & boost, HepAxisAngle & rotation) const;
161   void decompose (HepRotation & rotation, HepBoost & boost) const;
162   void decompose (HepBoost & boost, HepRotation & rotation) const;
163    // These are trivial, as the boost vector is 0.
164 
165   // ----------  Comparisons:
166 
167   inline bool isIdentity() const;
168   // Returns true if the identity matrix (Geant4).
169 
170   inline int compare( const HepRotationZ & r  ) const;
171   // Dictionary-order comparison, in order of delta
172   // Used in operator<, >, <=, >=
173 
174   inline bool operator== ( const HepRotationZ & r ) const;
175   inline bool operator!= ( const HepRotationZ & r ) const;
176   inline bool operator<  ( const HepRotationZ & r ) const;
177   inline bool operator>  ( const HepRotationZ & r ) const;
178   inline bool operator<= ( const HepRotationZ & r ) const;
179   inline bool operator>= ( const HepRotationZ & r ) const;
180 
181   double distance2( const HepRotationZ & r  ) const;
182   // 3 - Tr ( this/r )
183 
184   double distance2( const HepRotation &  r  ) const;
185   // 3 - Tr ( this/r ) -- This works with RotationY or Z also
186 
187   double howNear( const HepRotationZ & r ) const;
188   double howNear( const HepRotation  & r ) const;
189   bool isNear( const HepRotationZ & r,
190              double epsilon=Hep4RotationInterface::tolerance) const;
191   bool isNear( const HepRotation  & r,
192              double epsilon=Hep4RotationInterface::tolerance) const;
193 
194   double distance2( const HepBoost           & lt  ) const;
195   // 3 - Tr ( this ) + |b|^2 / (1-|b|^2)
196   double distance2( const HepLorentzRotation & lt  ) const;
197   // 3 - Tr ( this/r ) + |b|^2 / (1-|b|^2) where b is the boost vector of lt
198 
199   double howNear( const HepBoost           & lt ) const;
200   double howNear( const HepLorentzRotation & lt ) const;
201   bool isNear( const HepBoost           & lt,
202              double epsilon=Hep4RotationInterface::tolerance) const;
203   bool isNear( const HepLorentzRotation & lt,
204              double epsilon=Hep4RotationInterface::tolerance) const;
205 
206   // ----------  Properties:
207 
208   double norm2() const;
209   // distance2 (IDENTITY), which is 3 - Tr ( *this )
210 
211   inline void rectify();
212   // non-const but logically moot correction for accumulated roundoff errors
213 
214   // ---------- Application:
215 
216   inline Hep3Vector operator() (const Hep3Vector & p) const;
217   // Rotate a Hep3Vector.
218 
219   inline  Hep3Vector operator * (const Hep3Vector & p) const;
220   // Multiplication with a Hep3Vector.
221 
222   inline HepLorentzVector operator()( const HepLorentzVector & w ) const;
223   // Rotate (the space part of) a HepLorentzVector.
224 
225   inline  HepLorentzVector operator* ( const HepLorentzVector & w ) const;
226   // Multiplication with a HepLorentzVector.
227 
228  // ---------- Operations in the group of Rotations
229 
230   inline HepRotationZ operator * (const HepRotationZ & rz) const;
231   // Product of two Z rotations:  (this) * rz is known to be RotationZ.
232 
233   // Product of two rotations (this) * b - matrix multiplication
234 
235   inline  HepRotationZ & operator *= (const HepRotationZ & r);
236   inline  HepRotationZ & transform   (const HepRotationZ & r);
237   // Matrix multiplication.
238   // Note a *= b; <=> a = a * b; while a.transform(b); <=> a = b * a;
239   // However, in this special case, they commute:  Both just add deltas.
240 
241   inline HepRotationZ inverse() const;
242   // Returns the inverse.
243 
244   friend HepRotationZ inverseOf(const HepRotationZ & r);
245   // Returns the inverse of a RotationZ.
246 
247   inline HepRotationZ & invert();
248   // Inverts the Rotation matrix (be negating delta).
249 
250   // ---------- I/O:
251 
252   std::ostream & print( std::ostream & os ) const;
253   // Output, identifying type of rotation and delta.
254 
255   // ---------- Tolerance
256 
257   static inline double getTolerance();
258   static inline double setTolerance(double tol);
259 
260 protected:
261 
262   double its_d;
263   // The angle of rotation.
264 
265   double its_s;
266   double its_c;
267   // Cache the trig functions, for rapid operations.
268 
269   inline HepRotationZ ( double dd, double ss, double cc );
270   // Unchecked load-the-data-members
271 
272   static inline double proper (double delta);
273   // Put an angle into the range of (-PI, PI].  Useful helper method.
274 
275 };  // HepRotationZ
276 
277 inline
278 std::ostream & operator <<
279 	( std::ostream & os, const HepRotationZ & r ) {return r.print(os);}
280 
281 // ---------- Free-function operations in the group of Rotations
282 
283 }  // namespace CLHEP
284 
285 #include "CLHEP/Vector/RotationZ.icc"
286 
287 #ifdef ENABLE_BACKWARDS_COMPATIBILITY
288 //  backwards compatibility will be enabled ONLY in CLHEP 1.9
289 using namespace CLHEP;
290 #endif
291 
292 #endif /* HEP_ROTATIONZ_H */
293 
294