CLHEP 2.4.6.4
C++ Class Library for High Energy Physics
Loading...
Searching...
No Matches
CLHEP::HepRotation Class Reference

#include <Rotation.h>

+ Inheritance diagram for CLHEP::HepRotation:

Classes

class  HepRotation_row
 

Public Member Functions

 HepRotation ()
 
 HepRotation (const HepRotation &m)
 
 HepRotation (HepRotation &&m)=default
 
 HepRotation (const HepRotationX &m)
 
 HepRotation (const HepRotationY &m)
 
 HepRotation (const HepRotationZ &m)
 
HepRotationset (const Hep3Vector &axis, double delta)
 
 HepRotation (const Hep3Vector &axis, double delta)
 
HepRotationset (const HepAxisAngle &ax)
 
 HepRotation (const HepAxisAngle &ax)
 
HepRotationset (double phi, double theta, double psi)
 
 HepRotation (double phi, double theta, double psi)
 
HepRotationset (const HepEulerAngles &e)
 
 HepRotation (const HepEulerAngles &e)
 
 HepRotation (const Hep3Vector &colX, const Hep3Vector &colY, const Hep3Vector &colZ)
 
HepRotationset (const Hep3Vector &colX, const Hep3Vector &colY, const Hep3Vector &colZ)
 
HepRotationsetRows (const Hep3Vector &rowX, const Hep3Vector &rowY, const Hep3Vector &rowZ)
 
HepRotationset (const HepRotationX &r)
 
HepRotationset (const HepRotationY &r)
 
HepRotationset (const HepRotationZ &r)
 
HepRotationoperator= (const HepRotation &r)
 
HepRotationoperator= (HepRotation &&r)=default
 
HepRotationoperator= (const HepRotationX &r)
 
HepRotationoperator= (const HepRotationY &r)
 
HepRotationoperator= (const HepRotationZ &r)
 
HepRotationset (const HepRep3x3 &m)
 
 HepRotation (const HepRep3x3 &m)
 
 ~HepRotation ()
 
Hep3Vector colX () const
 
Hep3Vector colY () const
 
Hep3Vector colZ () const
 
Hep3Vector rowX () const
 
Hep3Vector rowY () const
 
Hep3Vector rowZ () const
 
double xx () const
 
double xy () const
 
double xz () const
 
double yx () const
 
double yy () const
 
double yz () const
 
double zx () const
 
double zy () const
 
double zz () const
 
HepRep3x3 rep3x3 () const
 
const HepRotation_row operator[] (int) const
 
double operator() (int, int) const
 
double getPhi () const
 
double getTheta () const
 
double getPsi () const
 
double phi () const
 
double theta () const
 
double psi () const
 
HepEulerAngles eulerAngles () const
 
double getDelta () const
 
Hep3Vector getAxis () const
 
double delta () const
 
Hep3Vector axis () const
 
HepAxisAngle axisAngle () const
 
void getAngleAxis (double &delta, Hep3Vector &axis) const
 
double phiX () const
 
double phiY () const
 
double phiZ () const
 
double thetaX () const
 
double thetaY () const
 
double thetaZ () const
 
HepLorentzVector col1 () const
 
HepLorentzVector col2 () const
 
HepLorentzVector col3 () const
 
HepLorentzVector col4 () const
 
HepLorentzVector row1 () const
 
HepLorentzVector row2 () const
 
HepLorentzVector row3 () const
 
HepLorentzVector row4 () const
 
double xt () const
 
double yt () const
 
double zt () const
 
double tx () const
 
double ty () const
 
double tz () const
 
double tt () const
 
HepRep4x4 rep4x4 () const
 
void setPhi (double phi)
 
void setTheta (double theta)
 
void setPsi (double psi)
 
void setAxis (const Hep3Vector &axis)
 
void setDelta (double delta)
 
void decompose (HepAxisAngle &rotation, Hep3Vector &boost) const
 
void decompose (Hep3Vector &boost, HepAxisAngle &rotation) const
 
bool isIdentity () const
 
int compare (const HepRotation &r) const
 
bool operator== (const HepRotation &r) const
 
bool operator!= (const HepRotation &r) const
 
bool operator< (const HepRotation &r) const
 
bool operator> (const HepRotation &r) const
 
bool operator<= (const HepRotation &r) const
 
bool operator>= (const HepRotation &r) const
 
double distance2 (const HepRotation &r) const
 
double howNear (const HepRotation &r) const
 
bool isNear (const HepRotation &r, double epsilon=Hep4RotationInterface::tolerance) const
 
double distance2 (const HepBoost &lt) const
 
double distance2 (const HepLorentzRotation &lt) const
 
double howNear (const HepBoost &lt) const
 
double howNear (const HepLorentzRotation &lt) const
 
bool isNear (const HepBoost &lt, double epsilon=Hep4RotationInterface::tolerance) const
 
bool isNear (const HepLorentzRotation &lt, double epsilon=Hep4RotationInterface::tolerance) const
 
double norm2 () const
 
void rectify ()
 
Hep3Vector operator() (const Hep3Vector &p) const
 
Hep3Vector operator* (const Hep3Vector &p) const
 
HepLorentzVector operator() (const HepLorentzVector &w) const
 
HepLorentzVector operator* (const HepLorentzVector &w) const
 
HepRotation operator* (const HepRotation &r) const
 
HepRotation operator* (const HepRotationX &rx) const
 
HepRotation operator* (const HepRotationY &ry) const
 
HepRotation operator* (const HepRotationZ &rz) const
 
HepRotationoperator*= (const HepRotation &r)
 
HepRotationtransform (const HepRotation &r)
 
HepRotationoperator*= (const HepRotationX &r)
 
HepRotationoperator*= (const HepRotationY &r)
 
HepRotationoperator*= (const HepRotationZ &r)
 
HepRotationtransform (const HepRotationX &r)
 
HepRotationtransform (const HepRotationY &r)
 
HepRotationtransform (const HepRotationZ &r)
 
HepRotationrotateX (double delta)
 
HepRotationrotateY (double delta)
 
HepRotationrotateZ (double delta)
 
HepRotationrotate (double delta, const Hep3Vector &axis)
 
HepRotationrotate (double delta, const Hep3Vector *axis)
 
HepRotationrotateAxes (const Hep3Vector &newX, const Hep3Vector &newY, const Hep3Vector &newZ)
 
HepRotation inverse () const
 
HepRotationinvert ()
 
std::ostream & print (std::ostream &os) const
 

Static Public Member Functions

static double getTolerance ()
 
static double setTolerance (double tol)
 

Static Public Attributes

static const HepRotation IDENTITY
 

Protected Member Functions

 HepRotation (double mxx, double mxy, double mxz, double myx, double myy, double myz, double mzx, double mzy, double mzz)
 

Protected Attributes

double rxx
 
double rxy
 
double rxz
 
double ryx
 
double ryy
 
double ryz
 
double rzx
 
double rzy
 
double rzz
 

Friends

HepRotation operator* (const HepRotationX &rx, const HepRotation &r)
 
HepRotation operator* (const HepRotationY &ry, const HepRotation &r)
 
HepRotation operator* (const HepRotationZ &rz, const HepRotation &r)
 

Detailed Description

Author

Definition at line 44 of file Rotation.h.

Constructor & Destructor Documentation

◆ HepRotation() [1/13]

CLHEP::HepRotation::HepRotation ( )
inline

Referenced by rotateAxes().

◆ HepRotation() [2/13]

CLHEP::HepRotation::HepRotation ( const HepRotation m)
inline

◆ HepRotation() [3/13]

CLHEP::HepRotation::HepRotation ( HepRotation &&  m)
inlinedefault

◆ HepRotation() [4/13]

CLHEP::HepRotation::HepRotation ( const HepRotationX m)
inline

◆ HepRotation() [5/13]

CLHEP::HepRotation::HepRotation ( const HepRotationY m)
inline

◆ HepRotation() [6/13]

CLHEP::HepRotation::HepRotation ( const HepRotationZ m)
inline

◆ HepRotation() [7/13]

CLHEP::HepRotation::HepRotation ( const Hep3Vector axis,
double  delta 
)

Definition at line 51 of file RotationA.cc.

52{
53 set( aaxis, ddelta );
54}
HepRotation & set(const Hep3Vector &axis, double delta)
Definition: RotationA.cc:24

◆ HepRotation() [8/13]

CLHEP::HepRotation::HepRotation ( const HepAxisAngle ax)

Definition at line 58 of file RotationA.cc.

59{
60 set ( ax.axis(), ax.delta() );
61}

◆ HepRotation() [9/13]

CLHEP::HepRotation::HepRotation ( double  phi,
double  theta,
double  psi 
)

Definition at line 56 of file RotationE.cc.

57{
58 set (phi1, theta1, psi1);
59}

◆ HepRotation() [10/13]

CLHEP::HepRotation::HepRotation ( const HepEulerAngles e)

Definition at line 63 of file RotationE.cc.

64{
65 set(e.phi(), e.theta(), e.psi());
66}

◆ HepRotation() [11/13]

CLHEP::HepRotation::HepRotation ( const Hep3Vector colX,
const Hep3Vector colY,
const Hep3Vector colZ 
)

Definition at line 131 of file RotationC.cc.

134{
135 set (ccolX, ccolY, ccolZ);
136}

◆ HepRotation() [12/13]

CLHEP::HepRotation::HepRotation ( const HepRep3x3 m)
inline

◆ ~HepRotation()

CLHEP::HepRotation::~HepRotation ( )
inline

◆ HepRotation() [13/13]

CLHEP::HepRotation::HepRotation ( double  mxx,
double  mxy,
double  mxz,
double  myx,
double  myy,
double  myz,
double  mzx,
double  mzy,
double  mzz 
)
inlineprotected

Member Function Documentation

◆ axis()

Hep3Vector CLHEP::HepRotation::axis ( ) const

Definition at line 76 of file RotationA.cc.

76 {
77
78 const double eps = 1e-15;
79
80 double Ux = rzy - ryz;
81 double Uy = rxz - rzx;
82 double Uz = ryx - rxy;
83 if (std::abs(Ux) < eps && std::abs(Uy) < eps && std::abs(Uz) < eps) {
84
85 double cosdelta = (rxx + ryy + rzz - 1.0) / 2.0;
86 if (cosdelta > 0.0) return Hep3Vector(0,0,1); // angle = 0, any axis is good
87
88 double xxt = (rxx + 1)/2;
89 double yyt = (ryy + 1)/2;
90 double zzt = (rzz + 1)/2;
91 double xyt = (rxy + ryx)/4;
92 double xzt = (rxz + rzx)/4;
93 double yzt = (ryz + rzy)/4;
94 double x, y, z;
95
96 if (xxt > ryy && xxt > rzz) {
97 x = std::sqrt(xxt);
98 if (rzy - ryz < 0) x = -x;
99 y = xyt/x;
100 z = xzt/x;
101 return Hep3Vector( x, y, z ).unit();
102 } else if (yyt > zzt) {
103 y = std::sqrt(yyt);
104 if (rxz - rzx < 0) y = -y;
105 x = xyt/y;
106 z = yzt/y;
107 return Hep3Vector( x, y, z ).unit();
108 } else {
109 z = std::sqrt(zzt);
110 if (ryx - rxy < 0) z = -z;
111 x = xzt/z;
112 y = yzt/z;
113 return Hep3Vector( x, y, z ).unit();
114 }
115 } else {
116 return Hep3Vector( Ux, Uy, Uz ).unit();
117 }
118
119} // axis()

Referenced by axisAngle(), main(), rectify(), and setDelta().

◆ axisAngle()

HepAxisAngle CLHEP::HepRotation::axisAngle ( ) const

Definition at line 121 of file RotationA.cc.

121 {
122
123 return HepAxisAngle (axis(), delta());
124
125} // axisAngle()
Hep3Vector axis() const
Definition: RotationA.cc:76
double delta() const
Definition: RotationA.cc:63

Referenced by CLHEP::HepLorentzRotation::decompose(), decompose(), and XF::Pow::operator()().

◆ col1()

HepLorentzVector CLHEP::HepRotation::col1 ( ) const
inline

◆ col2()

HepLorentzVector CLHEP::HepRotation::col2 ( ) const
inline

◆ col3()

HepLorentzVector CLHEP::HepRotation::col3 ( ) const
inline

◆ col4()

HepLorentzVector CLHEP::HepRotation::col4 ( ) const
inline

◆ colX()

Hep3Vector CLHEP::HepRotation::colX ( ) const
inline

◆ colY()

Hep3Vector CLHEP::HepRotation::colY ( ) const
inline

◆ colZ()

Hep3Vector CLHEP::HepRotation::colZ ( ) const
inline

◆ compare()

int CLHEP::HepRotation::compare ( const HepRotation r) const

Definition at line 174 of file Rotation.cc.

174 {
175 if (rzz<r.rzz) return -1; else if (rzz>r.rzz) return 1;
176 else if (rzy<r.rzy) return -1; else if (rzy>r.rzy) return 1;
177 else if (rzx<r.rzx) return -1; else if (rzx>r.rzx) return 1;
178 else if (ryz<r.ryz) return -1; else if (ryz>r.ryz) return 1;
179 else if (ryy<r.ryy) return -1; else if (ryy>r.ryy) return 1;
180 else if (ryx<r.ryx) return -1; else if (ryx>r.ryx) return 1;
181 else if (rxz<r.rxz) return -1; else if (rxz>r.rxz) return 1;
182 else if (rxy<r.rxy) return -1; else if (rxy>r.rxy) return 1;
183 else if (rxx<r.rxx) return -1; else if (rxx>r.rxx) return 1;
184 else return 0;
185}

◆ decompose() [1/2]

void CLHEP::HepRotation::decompose ( Hep3Vector boost,
HepAxisAngle rotation 
) const

Definition at line 26 of file RotationP.cc.

26 {
27 boost.set(0,0,0);
28 rotation = axisAngle();
29}
HepAxisAngle axisAngle() const
Definition: RotationA.cc:121

◆ decompose() [2/2]

void CLHEP::HepRotation::decompose ( HepAxisAngle rotation,
Hep3Vector boost 
) const

Definition at line 21 of file RotationP.cc.

21 {
22 boost.set(0,0,0);
23 rotation = axisAngle();
24}

◆ delta()

double CLHEP::HepRotation::delta ( ) const

Definition at line 63 of file RotationA.cc.

63 {
64
65 double cosdelta = (rxx + ryy + rzz - 1.0) / 2.0;
66 if (cosdelta > 1.0) {
67 return 0;
68 } else if (cosdelta < -1.0) {
69 return CLHEP::pi;
70 } else {
71 return std::acos( cosdelta ); // Already safe due to the cosdelta > 1 check
72 }
73
74} // delta()

Referenced by axisAngle(), main(), rectify(), and setAxis().

◆ distance2() [1/3]

double CLHEP::HepRotation::distance2 ( const HepBoost lt) const

Definition at line 36 of file RotationL.cc.

36 {
37 return distance2( HepLorentzRotation(lt));
38}
double distance2(const HepRotation &r) const
Definition: RotationP.cc:31

◆ distance2() [2/3]

double CLHEP::HepRotation::distance2 ( const HepLorentzRotation lt) const

Definition at line 26 of file RotationL.cc.

26 {
27 HepAxisAngle a;
28 Hep3Vector b;
29 lt.decompose(b, a);
30 double bet = b.beta();
31 double bet2 = bet*bet;
32 HepRotation r(a);
33 return bet2/(1-bet2) + distance2(r);
34}

◆ distance2() [3/3]

double CLHEP::HepRotation::distance2 ( const HepRotation r) const

Definition at line 31 of file RotationP.cc.

31 {
32 double sum = rxx * r.rxx + rxy * r.rxy + rxz * r.rxz
33 + ryx * r.ryx + ryy * r.ryy + ryz * r.ryz
34 + rzx * r.rzx + rzy * r.rzy + rzz * r.rzz;
35 double answer = 3.0 - sum;
36 return (answer >= 0 ) ? answer : 0;
37}

Referenced by distance2(), CLHEP::HepLorentzRotation::distance2(), howNear(), isNear(), and CLHEP::HepLorentzRotation::isNear().

◆ eulerAngles()

HepEulerAngles CLHEP::HepRotation::eulerAngles ( ) const

Definition at line 203 of file RotationE.cc.

203 {
204
205 // Please see the mathematical justification in eulerAngleComputations.ps
206
207 double phi1, theta1, psi1;
208 double psiPlusPhi, psiMinusPhi;
209
210 theta1 = safe_acos( rzz );
211
212 if (rzz > 1 || rzz < -1) {
213 ZMthrowC ( ZMxpvImproperRotation (
214 "HepRotation::eulerAngles() finds | rzz | > 1 "));
215 }
216
217 double cosTheta = rzz;
218 if (cosTheta > 1) cosTheta = 1;
219 if (cosTheta < -1) cosTheta = -1;
220
221 if (cosTheta == 1) {
222 psiPlusPhi = std::atan2 ( rxy - ryx, rxx + ryy );
223 psiMinusPhi = 0;
224
225 } else if (cosTheta >= 0) {
226
227 // In this realm, the atan2 expression for psi + phi is numerically stable
228 psiPlusPhi = std::atan2 ( rxy - ryx, rxx + ryy );
229
230 // psi - phi is potentially more subtle, but when unstable it is moot
231 double s1 = -rxy - ryx; // std::sin (psi-phi) * (1 - cos theta)
232 double c = rxx - ryy; // std::cos (psi-phi) * (1 - cos theta)
233 psiMinusPhi = std::atan2 ( s1, c );
234
235 } else if (cosTheta > -1) {
236
237 // In this realm, the atan2 expression for psi - phi is numerically stable
238 psiMinusPhi = std::atan2 ( -rxy - ryx, rxx - ryy );
239
240 // psi + phi is potentially more subtle, but when unstable it is moot
241 double s1 = rxy - ryx; // std::sin (psi+phi) * (1 + cos theta)
242 double c = rxx + ryy; // std::cos (psi+phi) * (1 + cos theta)
243 psiPlusPhi = std::atan2 ( s1, c );
244
245 } else { // cosTheta == -1
246
247 psiMinusPhi = std::atan2 ( -rxy - ryx, rxx - ryy );
248 psiPlusPhi = 0;
249
250 }
251
252 psi1 = .5 * (psiPlusPhi + psiMinusPhi);
253 phi1 = .5 * (psiPlusPhi - psiMinusPhi);
254
255 // Now correct by pi if we have managed to get a value of psiPlusPhi
256 // or psiMinusPhi that was off by 2 pi:
257 correctPsiPhi ( rxz, rzx, ryz, rzy, psi1, phi1 );
258
259 return HepEulerAngles( phi1, theta1, psi1 );
260
261} // eulerAngles()
#define ZMthrowC(A)
Definition: ZMxpv.h:133

Referenced by phi(), psi(), and test().

◆ getAngleAxis()

void CLHEP::HepRotation::getAngleAxis ( double delta,
Hep3Vector axis 
) const

Definition at line 149 of file Rotation.cc.

149 {
150 double cosa = 0.5*(xx()+yy()+zz()-1);
151 double cosa1 = 1-cosa;
152 if (cosa1 <= 0) {
153 angle = 0;
154 aaxis = Hep3Vector(0,0,1);
155 }else{
156 double x=0, y=0, z=0;
157 if (xx() > cosa) x = std::sqrt((xx()-cosa)/cosa1);
158 if (yy() > cosa) y = std::sqrt((yy()-cosa)/cosa1);
159 if (zz() > cosa) z = std::sqrt((zz()-cosa)/cosa1);
160 if (zy() < yz()) x = -x;
161 if (xz() < zx()) y = -y;
162 if (yx() < xy()) z = -z;
163 angle = (cosa < -1.) ? std::acos(-1.) : std::acos(cosa);
164 aaxis = Hep3Vector(x,y,z);
165 }
166}
double zz() const
double yz() const
double zx() const
double yx() const
double zy() const
double xx() const
double yy() const
double xz() const
double xy() const

◆ getAxis()

Hep3Vector CLHEP::HepRotation::getAxis ( ) const
inline

◆ getDelta()

double CLHEP::HepRotation::getDelta ( ) const
inline

◆ getPhi()

double CLHEP::HepRotation::getPhi ( ) const
inline

◆ getPsi()

double CLHEP::HepRotation::getPsi ( ) const
inline

◆ getTheta()

double CLHEP::HepRotation::getTheta ( ) const
inline

◆ getTolerance()

static double CLHEP::HepRotation::getTolerance ( )
inlinestatic

◆ howNear() [1/3]

double CLHEP::HepRotation::howNear ( const HepBoost lt) const

Definition at line 44 of file RotationL.cc.

44 {
45 return std::sqrt( distance2( lt ) );
46}

◆ howNear() [2/3]

double CLHEP::HepRotation::howNear ( const HepLorentzRotation lt) const

Definition at line 40 of file RotationL.cc.

40 {
41 return std::sqrt( distance2( lt ) );
42}

◆ howNear() [3/3]

double CLHEP::HepRotation::howNear ( const HepRotation r) const

Definition at line 39 of file RotationP.cc.

39 {
40 return std::sqrt( distance2( r ) );
41}

◆ inverse()

HepRotation CLHEP::HepRotation::inverse ( ) const
inline

Referenced by main().

◆ invert()

HepRotation & CLHEP::HepRotation::invert ( )
inline

Referenced by main(), and setRows().

◆ isIdentity()

bool CLHEP::HepRotation::isIdentity ( ) const

Definition at line 168 of file Rotation.cc.

168 {
169 return (rxx == 1.0 && rxy == 0.0 && rxz == 0.0 &&
170 ryx == 0.0 && ryy == 1.0 && ryz == 0.0 &&
171 rzx == 0.0 && rzy == 0.0 && rzz == 1.0) ? true : false;
172}

◆ isNear() [1/3]

bool CLHEP::HepRotation::isNear ( const HepBoost lt,
double  epsilon = Hep4RotationInterface::tolerance 
) const

Definition at line 53 of file RotationL.cc.

54 {
55 return distance2( lt ) <= epsilon*epsilon;
56}

◆ isNear() [2/3]

bool CLHEP::HepRotation::isNear ( const HepLorentzRotation lt,
double  epsilon = Hep4RotationInterface::tolerance 
) const

Definition at line 48 of file RotationL.cc.

49 {
50 return distance2( lt ) <= epsilon*epsilon;
51}

◆ isNear() [3/3]

bool CLHEP::HepRotation::isNear ( const HepRotation r,
double  epsilon = Hep4RotationInterface::tolerance 
) const

Definition at line 43 of file RotationP.cc.

44 {
45 return distance2( r ) <= epsilon*epsilon;
46}

◆ norm2()

◆ operator!=()

bool CLHEP::HepRotation::operator!= ( const HepRotation r) const
inline

◆ operator()() [1/3]

Hep3Vector CLHEP::HepRotation::operator() ( const Hep3Vector p) const
inline

◆ operator()() [2/3]

HepLorentzVector CLHEP::HepRotation::operator() ( const HepLorentzVector w) const
inline

◆ operator()() [3/3]

double CLHEP::HepRotation::operator() ( int  i,
int  j 
) const

Definition at line 25 of file Rotation.cc.

25 {
26 if (i == 0) {
27 if (j == 0) { return xx(); }
28 if (j == 1) { return xy(); }
29 if (j == 2) { return xz(); }
30 } else if (i == 1) {
31 if (j == 0) { return yx(); }
32 if (j == 1) { return yy(); }
33 if (j == 2) { return yz(); }
34 } else if (i == 2) {
35 if (j == 0) { return zx(); }
36 if (j == 1) { return zy(); }
37 if (j == 2) { return zz(); }
38 }
39 std::cerr << "HepRotation subscripting: bad indices "
40 << "(" << i << "," << j << ")" << std::endl;
41 return 0.0;
42}

◆ operator*() [1/6]

Hep3Vector CLHEP::HepRotation::operator* ( const Hep3Vector p) const
inline

◆ operator*() [2/6]

HepLorentzVector CLHEP::HepRotation::operator* ( const HepLorentzVector w) const
inline

◆ operator*() [3/6]

HepRotation CLHEP::HepRotation::operator* ( const HepRotation r) const
inline

◆ operator*() [4/6]

HepRotation CLHEP::HepRotation::operator* ( const HepRotationX rx) const
inline

◆ operator*() [5/6]

HepRotation CLHEP::HepRotation::operator* ( const HepRotationY ry) const
inline

◆ operator*() [6/6]

HepRotation CLHEP::HepRotation::operator* ( const HepRotationZ rz) const
inline

◆ operator*=() [1/4]

HepRotation & CLHEP::HepRotation::operator*= ( const HepRotation r)
inline

◆ operator*=() [2/4]

HepRotation & CLHEP::HepRotation::operator*= ( const HepRotationX r)
inline

◆ operator*=() [3/4]

HepRotation & CLHEP::HepRotation::operator*= ( const HepRotationY r)
inline

◆ operator*=() [4/4]

HepRotation & CLHEP::HepRotation::operator*= ( const HepRotationZ r)
inline

◆ operator<()

bool CLHEP::HepRotation::operator< ( const HepRotation r) const
inline

◆ operator<=()

bool CLHEP::HepRotation::operator<= ( const HepRotation r) const
inline

◆ operator=() [1/5]

HepRotation & CLHEP::HepRotation::operator= ( const HepRotation r)
inline

◆ operator=() [2/5]

HepRotation & CLHEP::HepRotation::operator= ( const HepRotationX r)
inline

◆ operator=() [3/5]

HepRotation & CLHEP::HepRotation::operator= ( const HepRotationY r)
inline

◆ operator=() [4/5]

HepRotation & CLHEP::HepRotation::operator= ( const HepRotationZ r)
inline

◆ operator=() [5/5]

HepRotation & CLHEP::HepRotation::operator= ( HepRotation &&  r)
inlinedefault

◆ operator==()

bool CLHEP::HepRotation::operator== ( const HepRotation r) const
inline

◆ operator>()

bool CLHEP::HepRotation::operator> ( const HepRotation r) const
inline

◆ operator>=()

bool CLHEP::HepRotation::operator>= ( const HepRotation r) const
inline

◆ operator[]()

const HepRotation_row CLHEP::HepRotation::operator[] ( int  ) const
inline

◆ phi()

double CLHEP::HepRotation::phi ( ) const

Definition at line 70 of file RotationE.cc.

70 {
71
72 double s2 = 1.0 - rzz*rzz;
73 if (s2 < 0) {
74 ZMthrowC ( ZMxpvImproperRotation (
75 "HepRotation::phi() finds | rzz | > 1 "));
76 s2 = 0;
77 }
78 const double sinTheta = std::sqrt( s2 );
79
80 if (sinTheta < .01) { // For theta close to 0 or PI, use the more stable
81 // algorithm to get all three Euler angles
82 HepEulerAngles ea = eulerAngles();
83 return ea.phi();
84 }
85
86 const double cscTheta = 1/sinTheta;
87 double cosabsphi = - rzy * cscTheta;
88 if ( std::fabs(cosabsphi) > 1 ) { // NaN-proofing
89 ZMthrowC ( ZMxpvImproperRotation (
90 "HepRotation::phi() finds | cos phi | > 1 "));
91 cosabsphi = 1;
92 }
93 const double absPhi = std::acos ( cosabsphi );
94 if (rzx > 0) {
95 return absPhi;
96 } else if (rzx < 0) {
97 return -absPhi;
98 } else {
99 return (rzy < 0) ? 0 : CLHEP::pi;
100 }
101
102} // phi()
HepEulerAngles eulerAngles() const
Definition: RotationE.cc:203

Referenced by main(), setPsi(), setTheta(), and test().

◆ phiX()

double CLHEP::HepRotation::phiX ( ) const

Definition at line 125 of file Rotation.cc.

125 {
126 return (yx() == 0.0 && xx() == 0.0) ? 0.0 : std::atan2(yx(),xx());
127}

◆ phiY()

double CLHEP::HepRotation::phiY ( ) const

Definition at line 129 of file Rotation.cc.

129 {
130 return (yy() == 0.0 && xy() == 0.0) ? 0.0 : std::atan2(yy(),xy());
131}

◆ phiZ()

double CLHEP::HepRotation::phiZ ( ) const

Definition at line 133 of file Rotation.cc.

133 {
134 return (yz() == 0.0 && xz() == 0.0) ? 0.0 : std::atan2(yz(),xz());
135}

◆ print()

std::ostream & CLHEP::HepRotation::print ( std::ostream &  os) const

Definition at line 18 of file RotationIO.cc.

18 {
19 os << "\n [ ( " <<
20 std::setw(11) << std::setprecision(6) << xx() << " " <<
21 std::setw(11) << std::setprecision(6) << xy() << " " <<
22 std::setw(11) << std::setprecision(6) << xz() << ")\n"
23 << " ( " <<
24 std::setw(11) << std::setprecision(6) << yx() << " " <<
25 std::setw(11) << std::setprecision(6) << yy() << " " <<
26 std::setw(11) << std::setprecision(6) << yz() << ")\n"
27 << " ( " <<
28 std::setw(11) << std::setprecision(6) << zx() << " " <<
29 std::setw(11) << std::setprecision(6) << zy() << " " <<
30 std::setw(11) << std::setprecision(6) << zz() << ") ]\n";
31 return os;
32}

Referenced by main().

◆ psi()

double CLHEP::HepRotation::psi ( ) const

Definition at line 110 of file RotationE.cc.

110 {
111
112 double sinTheta;
113 if ( std::fabs(rzz) > 1 ) { // NaN-proofing
114 ZMthrowC ( ZMxpvImproperRotation (
115 "HepRotation::psi() finds | rzz | > 1"));
116 sinTheta = 0;
117 } else {
118 sinTheta = std::sqrt( 1.0 - rzz*rzz );
119 }
120
121 if (sinTheta < .01) { // For theta close to 0 or PI, use the more stable
122 // algorithm to get all three Euler angles
123 HepEulerAngles ea = eulerAngles();
124 return ea.psi();
125 }
126
127 const double cscTheta = 1/sinTheta;
128 double cosabspsi = ryz * cscTheta;
129 if ( std::fabs(cosabspsi) > 1 ) { // NaN-proofing
130 ZMthrowC ( ZMxpvImproperRotation (
131 "HepRotation::psi() finds | cos psi | > 1"));
132 cosabspsi = 1;
133 }
134 const double absPsi = std::acos ( cosabspsi );
135 if (rxz > 0) {
136 return absPsi;
137 } else if (rxz < 0) {
138 return -absPsi;
139 } else {
140 return (ryz > 0) ? 0 : CLHEP::pi;
141 }
142
143} // psi()

Referenced by main(), setPhi(), setTheta(), and test().

◆ rectify()

void CLHEP::HepRotation::rectify ( )

Definition at line 149 of file RotationC.cc.

149 {
150 // Assuming the representation of this is close to a true Rotation,
151 // but may have drifted due to round-off error from many operations,
152 // this forms an "exact" orthonormal matrix for the rotation again.
153
154 // The first step is to average with the transposed inverse. This
155 // will correct for small errors such as those occuring when decomposing
156 // a LorentzTransformation. Then we take the bull by the horns and
157 // formally extract the axis and delta (assuming the Rotation were true)
158 // and re-setting the rotation according to those.
159
160 double det = rxx * ryy * rzz +
161 rxy * ryz * rzx +
162 rxz * ryx * rzy -
163 rxx * ryz * rzy -
164 rxy * ryx * rzz -
165 rxz * ryy * rzx ;
166 if (det <= 0) {
167 ZMthrowA(ZMxpvImproperRotation(
168 "Attempt to rectify a Rotation with determinant <= 0\n"));
169 return;
170 }
171 double di = 1.0 / det;
172
173 // xx, xy, ... are components of inverse matrix:
174 double xx1 = (ryy * rzz - ryz * rzy) * di;
175 double xy1 = (rzy * rxz - rzz * rxy) * di;
176 double xz1 = (rxy * ryz - rxz * ryy) * di;
177 double yx1 = (ryz * rzx - ryx * rzz) * di;
178 double yy1 = (rzz * rxx - rzx * rxz) * di;
179 double yz1 = (rxz * ryx - rxx * ryz) * di;
180 double zx1 = (ryx * rzy - ryy * rzx) * di;
181 double zy1 = (rzx * rxy - rzy * rxx) * di;
182 double zz1 = (rxx * ryy - rxy * ryx) * di;
183
184 // Now average with the TRANSPOSE of that:
185 rxx = .5*(rxx + xx1);
186 rxy = .5*(rxy + yx1);
187 rxz = .5*(rxz + zx1);
188 ryx = .5*(ryx + xy1);
189 ryy = .5*(ryy + yy1);
190 ryz = .5*(ryz + zy1);
191 rzx = .5*(rzx + xz1);
192 rzy = .5*(rzy + yz1);
193 rzz = .5*(rzz + zz1);
194
195 // Now force feed this improved rotation
196 double del = delta();
197 Hep3Vector u = axis();
198 u = u.unit(); // Because if the rotation is inexact, then the
199 // axis() returned will not have length 1!
200 set(u, del);
201
202} // rectify()
#define ZMthrowA(A)
Definition: ZMxpv.h:128

Referenced by CLHEP::HepLorentzRotation::decompose(), and CLHEP::HepLorentzRotation::rectify().

◆ rep3x3()

HepRep3x3 CLHEP::HepRotation::rep3x3 ( ) const
inline

Referenced by compareR(), and perturb().

◆ rep4x4()

HepRep4x4 CLHEP::HepRotation::rep4x4 ( ) const
inline

◆ rotate() [1/2]

HepRotation & CLHEP::HepRotation::rotate ( double  delta,
const Hep3Vector axis 
)

Definition at line 44 of file Rotation.cc.

44 {
45 if (a != 0.0) {
46 double ll = aaxis.mag();
47 if (ll == 0.0) {
48 ZMthrowC (ZMxpvZeroVector("HepRotation: zero axis"));
49 }else{
50 double sa = std::sin(a), ca = std::cos(a);
51 double dx = aaxis.x()/ll, dy = aaxis.y()/ll, dz = aaxis.z()/ll;
52 HepRotation m1(
53 ca+(1-ca)*dx*dx, (1-ca)*dx*dy-sa*dz, (1-ca)*dx*dz+sa*dy,
54 (1-ca)*dy*dx+sa*dz, ca+(1-ca)*dy*dy, (1-ca)*dy*dz-sa*dx,
55 (1-ca)*dz*dx-sa*dy, (1-ca)*dz*dy+sa*dx, ca+(1-ca)*dz*dz );
56 transform(m1);
57 }
58 }
59 return *this;
60}
HepRotation & transform(const HepRotation &r)

Referenced by main(), and CLHEP::Hep3Vector::rotate().

◆ rotate() [2/2]

HepRotation & CLHEP::HepRotation::rotate ( double  delta,
const Hep3Vector axis 
)
inline

◆ rotateAxes()

HepRotation & CLHEP::HepRotation::rotateAxes ( const Hep3Vector newX,
const Hep3Vector newY,
const Hep3Vector newZ 
)

Definition at line 101 of file Rotation.cc.

103 {
104 double del = 0.001;
105 Hep3Vector w = newX.cross(newY);
106
107 if (std::abs(newZ.x()-w.x()) > del ||
108 std::abs(newZ.y()-w.y()) > del ||
109 std::abs(newZ.z()-w.z()) > del ||
110 std::abs(newX.mag2()-1.) > del ||
111 std::abs(newY.mag2()-1.) > del ||
112 std::abs(newZ.mag2()-1.) > del ||
113 std::abs(newX.dot(newY)) > del ||
114 std::abs(newY.dot(newZ)) > del ||
115 std::abs(newZ.dot(newX)) > del) {
116 std::cerr << "HepRotation::rotateAxes: bad axis vectors" << std::endl;
117 return *this;
118 }else{
119 return transform(HepRotation(newX.x(), newY.x(), newZ.x(),
120 newX.y(), newY.y(), newZ.y(),
121 newX.z(), newY.z(), newZ.z()));
122 }
123}

◆ rotateX()

HepRotation & CLHEP::HepRotation::rotateX ( double  delta)

Definition at line 62 of file Rotation.cc.

62 {
63 double c1 = std::cos(a);
64 double s1 = std::sin(a);
65 double x1 = ryx, y1 = ryy, z1 = ryz;
66 ryx = c1*x1 - s1*rzx;
67 ryy = c1*y1 - s1*rzy;
68 ryz = c1*z1 - s1*rzz;
69 rzx = s1*x1 + c1*rzx;
70 rzy = s1*y1 + c1*rzy;
71 rzz = s1*z1 + c1*rzz;
72 return *this;
73}

◆ rotateY()

HepRotation & CLHEP::HepRotation::rotateY ( double  delta)

Definition at line 75 of file Rotation.cc.

75 {
76 double c1 = std::cos(a);
77 double s1 = std::sin(a);
78 double x1 = rzx, y1 = rzy, z1 = rzz;
79 rzx = c1*x1 - s1*rxx;
80 rzy = c1*y1 - s1*rxy;
81 rzz = c1*z1 - s1*rxz;
82 rxx = s1*x1 + c1*rxx;
83 rxy = s1*y1 + c1*rxy;
84 rxz = s1*z1 + c1*rxz;
85 return *this;
86}

Referenced by main().

◆ rotateZ()

HepRotation & CLHEP::HepRotation::rotateZ ( double  delta)

Definition at line 88 of file Rotation.cc.

88 {
89 double c1 = std::cos(a);
90 double s1 = std::sin(a);
91 double x1 = rxx, y1 = rxy, z1 = rxz;
92 rxx = c1*x1 - s1*ryx;
93 rxy = c1*y1 - s1*ryy;
94 rxz = c1*z1 - s1*ryz;
95 ryx = s1*x1 + c1*ryx;
96 ryy = s1*y1 + c1*ryy;
97 ryz = s1*z1 + c1*ryz;
98 return *this;
99}

Referenced by main().

◆ row1()

HepLorentzVector CLHEP::HepRotation::row1 ( ) const
inline

◆ row2()

HepLorentzVector CLHEP::HepRotation::row2 ( ) const
inline

◆ row3()

HepLorentzVector CLHEP::HepRotation::row3 ( ) const
inline

◆ row4()

HepLorentzVector CLHEP::HepRotation::row4 ( ) const
inline

◆ rowX()

Hep3Vector CLHEP::HepRotation::rowX ( ) const
inline

◆ rowY()

Hep3Vector CLHEP::HepRotation::rowY ( ) const
inline

◆ rowZ()

Hep3Vector CLHEP::HepRotation::rowZ ( ) const
inline

◆ set() [1/9]

HepRotation & CLHEP::HepRotation::set ( const Hep3Vector axis,
double  delta 
)

Definition at line 24 of file RotationA.cc.

24 {
25
26 double sinDelta = std::sin(ddelta), cosDelta = std::cos(ddelta);
27 double oneMinusCosDelta = 1.0 - cosDelta;
28
29 Hep3Vector u = aaxis.unit();
30
31 double uX = u.getX();
32 double uY = u.getY();
33 double uZ = u.getZ();
34
35 rxx = oneMinusCosDelta * uX * uX + cosDelta;
36 rxy = oneMinusCosDelta * uX * uY - sinDelta * uZ;
37 rxz = oneMinusCosDelta * uX * uZ + sinDelta * uY;
38
39 ryx = oneMinusCosDelta * uY * uX + sinDelta * uZ;
40 ryy = oneMinusCosDelta * uY * uY + cosDelta;
41 ryz = oneMinusCosDelta * uY * uZ - sinDelta * uX;
42
43 rzx = oneMinusCosDelta * uZ * uX - sinDelta * uY;
44 rzy = oneMinusCosDelta * uZ * uY + sinDelta * uX;
45 rzz = oneMinusCosDelta * uZ * uZ + cosDelta;
46
47 return *this;
48
49} // HepRotation::set(axis, delta)

Referenced by CLHEP::HepLorentzRotation::decompose(), HepRotation(), main(), perturb(), rectify(), set(), setAxis(), setDelta(), setPhi(), setPsi(), setRows(), setTheta(), and test().

◆ set() [2/9]

HepRotation & CLHEP::HepRotation::set ( const Hep3Vector colX,
const Hep3Vector colY,
const Hep3Vector colZ 
)

Definition at line 72 of file RotationC.cc.

74 {
75 Hep3Vector ucolX = ccolX.unit();
76 Hep3Vector ucolY = ccolY.unit();
77 Hep3Vector ucolZ = ccolZ.unit();
78
79 double u1u2 = ucolX.dot(ucolY);
80 double f12 = std::fabs(u1u2);
82 ZMthrowC (ZMxpvNotOrthogonal(
83 "col's X and Y supplied for Rotation are not close to orthogonal"));
84 }
85 double u1u3 = ucolX.dot(ucolZ);
86 double f13 = std::fabs(u1u3);
88 ZMthrowC (ZMxpvNotOrthogonal(
89 "col's X and Z supplied for Rotation are not close to orthogonal"));
90 }
91 double u2u3 = ucolY.dot(ucolZ);
92 double f23 = std::fabs(u2u3);
94 ZMthrowC (ZMxpvNotOrthogonal(
95 "col's Y and Z supplied for Rotation are not close to orthogonal"));
96 }
97
98 Hep3Vector v1, v2, v3;
99 bool isRotation;
100 if ( (f12 <= f13) && (f12 <= f23) ) {
101 isRotation = setCols ( ucolX, ucolY, ucolZ, u1u2, v1, v2, v3 );
102 if ( !isRotation ) {
103 ZMthrowC (ZMxpvImproperRotation(
104 "col's X Y and Z supplied form closer to a reflection than a Rotation "
105 "\n col Z is set to col X cross col Y"));
106 }
107 } else if ( f13 <= f23 ) {
108 isRotation = setCols ( ucolZ, ucolX, ucolY, u1u3, v3, v1, v2 );
109 if ( !isRotation ) {
110 ZMthrowC (ZMxpvImproperRotation(
111 "col's X Y and Z supplied form closer to a reflection than a Rotation "
112 "\n col Y is set to col Z cross col X"));
113 }
114 } else {
115 isRotation = setCols ( ucolY, ucolZ, ucolX, u2u3, v2, v3, v1 );
116 if ( !isRotation ) {
117 ZMthrowC (ZMxpvImproperRotation(
118 "col's X Y and Z supplied form closer to a reflection than a Rotation "
119 "\n col X is set to col Y cross col Z"));
120 }
121 }
122
123 rxx = v1.x(); ryx = v1.y(); rzx = v1.z();
124 rxy = v2.x(); ryy = v2.y(); rzy = v2.z();
125 rxz = v3.x(); ryz = v3.y(); rzz = v3.z();
126
127 return *this;
128
129} // HepRotation::set(colX, colY, colZ)

◆ set() [3/9]

HepRotation & CLHEP::HepRotation::set ( const HepAxisAngle ax)

Definition at line 55 of file RotationA.cc.

55 {
56 return set ( ax.axis(), ax.delta() );
57}

◆ set() [4/9]

HepRotation & CLHEP::HepRotation::set ( const HepEulerAngles e)

Definition at line 60 of file RotationE.cc.

60 {
61 return set(e.phi(), e.theta(), e.psi());
62}

◆ set() [5/9]

HepRotation & CLHEP::HepRotation::set ( const HepRep3x3 m)
inline

◆ set() [6/9]

HepRotation & CLHEP::HepRotation::set ( const HepRotationX r)
inline

◆ set() [7/9]

HepRotation & CLHEP::HepRotation::set ( const HepRotationY r)
inline

◆ set() [8/9]

HepRotation & CLHEP::HepRotation::set ( const HepRotationZ r)
inline

◆ set() [9/9]

HepRotation & CLHEP::HepRotation::set ( double  phi,
double  theta,
double  psi 
)

Definition at line 34 of file RotationE.cc.

34 {
35
36 double sinPhi = std::sin( phi1 ), cosPhi = std::cos( phi1 );
37 double sinTheta = std::sin( theta1 ), cosTheta = std::cos( theta1 );
38 double sinPsi = std::sin( psi1 ), cosPsi = std::cos( psi1 );
39
40 rxx = cosPsi * cosPhi - cosTheta * sinPhi * sinPsi;
41 rxy = cosPsi * sinPhi + cosTheta * cosPhi * sinPsi;
42 rxz = sinPsi * sinTheta;
43
44 ryx = - sinPsi * cosPhi - cosTheta * sinPhi * cosPsi;
45 ryy = - sinPsi * sinPhi + cosTheta * cosPhi * cosPsi;
46 ryz = cosPsi * sinTheta;
47
48 rzx = sinTheta * sinPhi;
49 rzy = - sinTheta * cosPhi;
50 rzz = cosTheta;
51
52 return *this;
53
54} // Rotation::set(phi, theta, psi)

◆ setAxis()

void CLHEP::HepRotation::setAxis ( const Hep3Vector axis)

Definition at line 128 of file RotationA.cc.

128 {
129 set ( aaxis, delta() );
130}

◆ setDelta()

void CLHEP::HepRotation::setDelta ( double  delta)

Definition at line 132 of file RotationA.cc.

132 {
133 set ( axis(), ddelta );
134}

◆ setPhi()

void CLHEP::HepRotation::setPhi ( double  phi)

Definition at line 265 of file RotationE.cc.

265 {
266 set ( phi1, theta(), psi() );
267}
double psi() const
Definition: RotationE.cc:110
double theta() const
Definition: RotationE.cc:104

◆ setPsi()

void CLHEP::HepRotation::setPsi ( double  psi)

Definition at line 273 of file RotationE.cc.

273 {
274 set ( phi(), theta(), psi1 );
275}
double phi() const
Definition: RotationE.cc:70

◆ setRows()

HepRotation & CLHEP::HepRotation::setRows ( const Hep3Vector rowX,
const Hep3Vector rowY,
const Hep3Vector rowZ 
)

Definition at line 138 of file RotationC.cc.

140 {
141 set (rrowX, rrowY, rrowZ);
142 invert();
143 return *this;
144}
HepRotation & invert()

◆ setTheta()

void CLHEP::HepRotation::setTheta ( double  theta)

Definition at line 269 of file RotationE.cc.

269 {
270 set ( phi(), theta1, psi() );
271}

◆ setTolerance()

static double CLHEP::HepRotation::setTolerance ( double  tol)
inlinestatic

◆ theta()

double CLHEP::HepRotation::theta ( ) const

Definition at line 104 of file RotationE.cc.

104 {
105
106 return safe_acos( rzz );
107
108} // theta()

Referenced by main(), setPhi(), setPsi(), and test().

◆ thetaX()

double CLHEP::HepRotation::thetaX ( ) const

Definition at line 137 of file Rotation.cc.

137 {
138 return safe_acos(zx());
139}

◆ thetaY()

double CLHEP::HepRotation::thetaY ( ) const

Definition at line 141 of file Rotation.cc.

141 {
142 return safe_acos(zy());
143}

◆ thetaZ()

double CLHEP::HepRotation::thetaZ ( ) const

Definition at line 145 of file Rotation.cc.

145 {
146 return safe_acos(zz());
147}

◆ transform() [1/4]

HepRotation & CLHEP::HepRotation::transform ( const HepRotation r)
inline

Referenced by rotate(), and rotateAxes().

◆ transform() [2/4]

HepRotation & CLHEP::HepRotation::transform ( const HepRotationX r)
inline

◆ transform() [3/4]

HepRotation & CLHEP::HepRotation::transform ( const HepRotationY r)
inline

◆ transform() [4/4]

HepRotation & CLHEP::HepRotation::transform ( const HepRotationZ r)
inline

◆ tt()

double CLHEP::HepRotation::tt ( ) const
inline

Referenced by CLHEP::operator*().

◆ tx()

double CLHEP::HepRotation::tx ( ) const
inline

Referenced by CLHEP::operator*().

◆ ty()

double CLHEP::HepRotation::ty ( ) const
inline

Referenced by CLHEP::operator*().

◆ tz()

double CLHEP::HepRotation::tz ( ) const
inline

Referenced by CLHEP::operator*().

◆ xt()

double CLHEP::HepRotation::xt ( ) const
inline

Referenced by CLHEP::operator*().

◆ xx()

◆ xy()

◆ xz()

◆ yt()

double CLHEP::HepRotation::yt ( ) const
inline

Referenced by CLHEP::operator*().

◆ yx()

◆ yy()

◆ yz()

◆ zt()

double CLHEP::HepRotation::zt ( ) const
inline

Referenced by CLHEP::operator*().

◆ zx()

◆ zy()

◆ zz()

Friends And Related Function Documentation

◆ operator* [1/3]

HepRotation operator* ( const HepRotationX rx,
const HepRotation r 
)
friend

◆ operator* [2/3]

HepRotation operator* ( const HepRotationY ry,
const HepRotation r 
)
friend

◆ operator* [3/3]

HepRotation operator* ( const HepRotationZ rz,
const HepRotation r 
)
friend

Member Data Documentation

◆ IDENTITY

const HepRotation CLHEP::HepRotation::IDENTITY
static

Definition at line 368 of file Rotation.h.

◆ rxx

double CLHEP::HepRotation::rxx
protected

◆ rxy

double CLHEP::HepRotation::rxy
protected

◆ rxz

double CLHEP::HepRotation::rxz
protected

◆ ryx

double CLHEP::HepRotation::ryx
protected

◆ ryy

double CLHEP::HepRotation::ryy
protected

◆ ryz

double CLHEP::HepRotation::ryz
protected

◆ rzx

double CLHEP::HepRotation::rzx
protected

◆ rzy

double CLHEP::HepRotation::rzy
protected

◆ rzz

double CLHEP::HepRotation::rzz
protected

The documentation for this class was generated from the following files: