CLHEP 2.4.6.4
C++ Class Library for High Energy Physics
Loading...
Searching...
No Matches
Functions | Variables
eulerTest.cc File Reference
#include "CLHEP/Vector/Rotation.h"
#include "CLHEP/Vector/EulerAngles.h"
#include "CLHEP/Units/PhysicalConstants.h"
#include <iostream>
#include <math.h>

Go to the source code of this file.

Functions

void perturb (int i, const HepRotation &r, HepRotation &rp, double &del)
 
bool compareR (const HepRotation &r1, const HepRotation &r2, double tol)
 
bool test (double phi, double theta, double psi)
 
int main ()
 

Variables

const int Nperturbs = 24
 

Function Documentation

◆ compareR()

bool compareR ( const HepRotation &  r1,
const HepRotation &  r2,
double  tol 
)

Definition at line 144 of file eulerTest.cc.

144 {
145 HepRep3x3 m1 = r1.rep3x3();
146 HepRep3x3 m2 = r2.rep3x3();
147 double flaw = 0;
148 flaw = max (flaw, (m1.xx_ - m2.xx_));
149 flaw = max (flaw, (m1.xy_ - m2.xy_));
150 flaw = max (flaw, (m1.xz_ - m2.xz_));
151 flaw = max (flaw, (m1.yx_ - m2.yx_));
152 flaw = max (flaw, (m1.yy_ - m2.yy_));
153 flaw = max (flaw, (m1.yz_ - m2.yz_));
154 flaw = max (flaw, (m1.zx_ - m2.zx_));
155 flaw = max (flaw, (m1.zy_ - m2.zy_));
156 flaw = max (flaw, (m1.zz_ - m2.zz_));
157 if (flaw > 20*std::sqrt(tol)) {
158 cout << "???????? comparison flaw at level of " << flaw << "\n"
159 << r1 << r2;
160 }
161 cout << "flaw size is " << flaw << " (" << flaw/std::sqrt(tol) << ")\n\n";
162 return (flaw <= tol);
163}
CLHEP::HepRotation::rep3x3
HepRep3x3 rep3x3() const

Referenced by test().

◆ main()

int main ( )

Definition at line 73 of file eulerTest.cc.

73 {
74
75 bool res = true;
76 double PI = CLHEP::pi;
77
78 // Some cases not in the potentially unstable region:
79 res &= test ( .05, PI/5, .1 );
80 res &= test ( .4, PI/7, -.35 );
81 res &= test ( PI/3, PI/6, -PI/3 );
82 res &= test ( PI/5, PI/2, -PI/2.5 );
83 res &= test ( -PI/5, PI/2, PI/3 );
84 res &= test ( -4*PI/5, PI/2, PI/2.5 );
85 res &= test ( 4*PI/5, PI/2, 2*PI/3 );
86 res &= test ( -3*PI/4, PI/2, -PI/3 );
87 res &= test ( 5*PI/6, PI/2, -4*PI/5 );
88 res &= test ( 5*PI/6, PI/2, -PI/3 );
89
90 // Specialized cases
91 res &= test ( .05, 0, .1 );
92 res &= test ( .2, 0, 1.1 );
93 res &= test ( -.4, 0, .4 );
94 res &= test ( -2.4, 0, 2.0 );
95 res &= test ( -2.4, 0, -2.0 );
96 res &= test ( -2.2, 0, 1.8 );
97 res &= test ( -2.2, 0, -1.8 );
98 res &= test ( .05, PI, .1 );
99 res &= test ( .2, PI, 1.1 );
100 res &= test ( -.4, PI, .4 );
101 res &= test ( -2.4, PI, 2.0 );
102 res &= test ( -2.4, PI, -2.0 );
103 res &= test ( -2.2, PI, 1.8 );
104 res &= test ( -2.2, PI, -1.8 );
105
106 // Cases near zero
107 res &= test ( .1, .0000000004, .5 );
108 res &= test ( -1.2, .0000000004, .5 );
109 res &= test ( .7, .0000000004, -.6 );
110 res &= test ( 1.5, .0000000004, -1.1 );
111 res &= test ( 1.4, .0000000004, -1.5 );
112 res &= test ( -.1, .0000000000028, .5 );
113 res &= test ( -1.2, .0000000000028, -.5 );
114 res &= test ( .7, .0000000000028, -.6 );
115 res &= test ( -1.5, .0000000000028, -1.1 );
116 res &= test ( 1.4, .0000000000028, 1.5 );
117
118 // Cases near PI
119 double nearPI = PI - .00000002;
120 res &= test ( .1, nearPI, .5 );
121 res &= test ( -1.2, nearPI, .5 );
122 res &= test ( .7, nearPI, -.6 );
123 res &= test ( 1.5, nearPI, -1.1 );
124 res &= test ( 1.4, nearPI, -1.5 );
125 res &= test ( 2.4, nearPI, -1.6 );
126 res &= test ( 2.3, nearPI, 1.9 );
127 res &= test ( -2.8, nearPI, .6 );
128 res &= test ( -.4, nearPI, -3.1 );
129 nearPI = PI - .000000000009;
130 res &= test ( .1, nearPI, -.5 );
131 res &= test ( 1.2, nearPI, .5 );
132 res &= test ( .7, nearPI, -.6 );
133 res &= test ( 1.5, nearPI, 1.1 );
134 res &= test ( -1.4, nearPI, -1.5 );
135 res &= test ( 2.1, nearPI, -1.2 );
136 res &= test ( 2.9, nearPI, .9 );
137 res &= test ( -2.8, nearPI, 1.6 );
138 res &= test ( -.4, nearPI, -3.0 );
139
140 if (!res) return -1;
141 return 0;
142}
test
bool test(double phi, double theta, double psi)
Definition: eulerTest.cc:25

◆ perturb()

void perturb ( int  i,
const HepRotation &  r,
HepRotation &  rp,
double &  del 
)

Definition at line 165 of file eulerTest.cc.

165 {
166
167 HepRep3x3 p0 ( 1, 3, -2,
168 -1, -2, 4,
169 2, -1, -1 );
170
171 HepRep3x3 p1 ( 1, -1, -2,
172 1, 3, -1,
173 2, -1, -3 );
174
175 HepRep3x3 p2 ( 5, 1, -5,
176 -3, -2, 3,
177 -1, 4, -1 );
178
179 HepRep3x3 p3 ( -2, -2, 1,
180 -1, -2, -4,
181 4, 2, -2 );
182
183 HepRep3x3 p[4];
184 p[0] = p0;
185 p[1] = p1;
186 p[2] = p2;
187 p[3] = p3;
188
189 int cycle = i/4;
190 double q;
191 switch (cycle){
192 case 0: q = 1.0e-14; break;
193 case 1: q = 1.0e-12; break;
194 case 2: q = 1.0e-10; break;
195 case 3: q = 1.0e-8; break;
196 case 4: q = 1.0e-6; break;
197 case 5: q = 1.0e-4; break;
198 }
199 HepRep3x3 d = p[i%4];
200 HepRep3x3 m = r.rep3x3();
201 if ((m.zz_ + q*d.zz_) < -1) {
202 d.zz_ = -d.zz_;
203 }
204 cout << "i = " << i << " q is " << q << "\n";
205 rp.set (HepRep3x3 ( m.xx_ + q*d.xx_ , m.xy_ + q*d.xy_ , m.xz_ + q*d.xz_ ,
206 m.yx_ + q*d.yx_ , m.yy_ + q*d.yy_ , m.yz_ + q*d.yz_ ,
207 m.zx_ + q*d.zx_ , m.zy_ + q*d.zy_ , m.zz_ + q*d.zz_ ) );
208 del = q;
209}
CLHEP::HepRotation::set
HepRotation & set(const Hep3Vector &axis, double delta)
Definition: RotationA.cc:24

Referenced by test().

◆ test()

bool test ( double  phi,
double  theta,
double  psi 
)

Definition at line 25 of file eulerTest.cc.

25 {
26 HepRotation r (phi, theta, psi);
27 HepRotation rp;
28 HepRotation rpe;
29 HepEulerAngles e;
30 bool retval = true;
31 double del;
32 cout << "\n\n -------------------------------- \n\n";
33 for ( int i = 0; i < Nperturbs; ++i ) {
34 perturb ( i, r, rp, del );
35 e = rp.eulerAngles();
36 cout << "(" << phi <<", " << theta << ", " << psi << ") --> "<< e << "\n";
37 cout << " e.phi() = " << e.phi()
38 << " rp.phi() = " << rp.phi()
39 << " difference is " << (e.phi() - rp.phi())/del << " * del\n";
40 if ( std::fabs (e.phi() - rp.phi()) > 200*del ) {
41 cout << "phi discrepancy in " << rp << "\n";
42 cout << " e.phi() = " << e.phi()
43 << " rp.phi() = " << rp.phi()
44 << " difference is " << e.phi() - rp.phi() << "\n";
45 if (del < 1.0e-4) cout << "??????????\n";
46 retval = false;
47 }
48 if ( std::fabs (e.theta() - rp.theta()) > 200*del ) {
49 cout << "theta discrepancy in " << rp << "\n";
50 cout << " e.theta() = " << e.theta()
51 << " rp.theta() = " << rp.theta()
52 << " difference is " << e.theta() - rp.theta() << "\n";
53 if (del < 1.0e-4) cout << "??????????\n";
54 retval = false;
55 }
56 cout << " e.psi() = " << e.psi()
57 << " rp.psi() = " << rp.psi()
58 << " difference is " << (e.psi() - rp.psi())/del << " * del\n";
59 if ( std::fabs (e.psi() - rp.psi()) > 200*del ) {
60 cout << "psi discrepancy in " << rp << "\n";
61 cout << " e.psi() = " << e.psi()
62 << " rp.psi() = " << rp.psi()
63 << " difference is " << e.psi() - rp.psi() << "\n";
64 if (del < 1.0e-4) cout << "??????????\n";
65 retval = false;
66 }
67 rpe.set(e);
68 retval &= compareR (rpe, rp, del);
69 }
70 return retval;
71}
CLHEP::HepEulerAngles::phi
double phi() const
CLHEP::HepRotation::eulerAngles
HepEulerAngles eulerAngles() const
Definition: RotationE.cc:203
CLHEP::HepRotation::phi
double phi() const
Definition: RotationE.cc:70
perturb
void perturb(int i, const HepRotation &r, HepRotation &rp, double &del)
Definition: eulerTest.cc:165
Nperturbs
const int Nperturbs
Definition: eulerTest.cc:20
compareR
bool compareR(const HepRotation &r1, const HepRotation &r2, double tol)
Definition: eulerTest.cc:144

Referenced by main().

Variable Documentation

◆ Nperturbs

const int Nperturbs = 24

Definition at line 20 of file eulerTest.cc.

Referenced by test().