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

testThreeVector.cc
Go to the documentation of this file.
1 // -*- C++ -*-
2 // $Id: testThreeVector.cc,v 1.3 2003/08/08 13:47:09 garren Exp $
3 // ---------------------------------------------------------------------------
4 //
5 // This file is a part of the CLHEP - a Class Library for High Energy Physics.
6 //
7 // This is a small program for testing the Hep3Vector class
8 // and the interaction with the HepRotation class.
9 //
10 
11 #include "CLHEP/Units/GlobalSystemOfUnits.h" // to see shadowing problems
12 #include "CLHEP/Units/GlobalPhysicalConstants.h"
13 #include "CLHEP/Vector/ThreeVector.h"
14 #include "CLHEP/Vector/TwoVector.h"
15 #include "CLHEP/Vector/Rotation.h"
16 
17 #include <cmath>
18 #include <iostream>
19 #include <stdlib.h> // for exit
20 
21 using namespace CLHEP;
22 
23 
24 #define DEPS 1.0e-14
25 #define FEPS 1.0e-6
26 
27 bool approx(double a, double b, double eps) {
28  return bool( std::abs(a-b) < eps );
29 }
30 
31 bool
32 test(const Hep3Vector & p, double x, double y, double z,
33  double eps) {
34  return bool( approx(p.x(), x, eps) && approx(p.y(), y, eps) &&
35  approx(p.z(), z, eps) );
36 }
37 
38 bool
39 test2(const Hep2Vector & p, double x, double y, double eps) {
40  return bool( approx(p.x(), x, eps) && approx(p.y(), y, eps) );
41 }
42 
43 int main () {
44 
45 // test constructors:
46 
47  Hep3Vector d0; if ( !test(d0, 0.0, 0.0, 0.0, DEPS) ) exit(1);
48  Hep3Vector f0; if ( !test(f0, 0.0, 0.0, 0.0, FEPS) ) exit(1);
49  Hep3Vector d1(1.0); if ( !test(d1, 1.0, 0.0, 0.0, DEPS) ) exit(1);
50  Hep3Vector f1(1.0); if ( !test(f1, 1.0, 0.0, 0.0, FEPS) ) exit(1);
51  Hep3Vector d2(1.0, 1.0); if ( !test(d2, 1.0, 1.0, 0.0, DEPS) ) exit(1);
52  Hep3Vector f2(1.0, 1.0); if ( !test(f2, 1.0, 1.0, 0.0, FEPS) ) exit(1);
53  Hep3Vector d3(1.0, 1.0, 1.0); if ( !test(d3, 1.0, 1.0, 1.0, DEPS) ) exit(1);
54  Hep3Vector f3(1.0, 1.0, 1.0); if ( !test(f3, 1.0, 1.0, 1.0, FEPS) ) exit(1);
55  Hep3Vector d4(f3); if ( !test(d4, 1.0, 1.0, 1.0, DEPS) ) exit(1);
56  Hep3Vector f4(d3); if ( !test(f4, 1.0, 1.0, 1.0, FEPS) ) exit(1);
57 
58  Hep2Vector t3(d4); if ( !test2(t3, 1.0, 1.0, DEPS) ) exit(1);
59 
60 // test input/output from a stream
61 
62  std::cin >> d0; if ( !test(d0, 1.1, 2.2, 3.3, DEPS) ) exit(1);
63  std::cin >> f0; if ( !test(f0, 3.0, 2.0, 1.0, FEPS) ) exit(1);
64  std::cout << d0 << std::endl;
65  std::cout << f0 << std::endl;
66 
67 // test assignment:
68 
69  d4 = d1; if ( !test(d4, 1.0, 0.0, 0.0, DEPS) ) exit(2);
70  f4 = f1; if ( !test(f4, 1.0, 0.0, 0.0, FEPS) ) exit(2);
71  d4 = f1; if ( !test(d4, 1.0, 0.0, 0.0, FEPS) ) exit(2);
72  f4 = d1; if ( !test(f4, 1.0, 0.0, 0.0, FEPS) ) exit(2);
73 
74 // test addition:
75 
76  d4 = d1 + d2; if ( !test(d4, 2.0, 1.0, 0.0, DEPS) ) exit(3);
77  d4 = f1 + d2; if ( !test(d4, 2.0, 1.0, 0.0, FEPS) ) exit(3);
78  d4 = d1 + f2; if ( !test(d4, 2.0, 1.0, 0.0, FEPS) ) exit(3);
79  d4 = f1 + f2; if ( !test(d4, 2.0, 1.0, 0.0, FEPS) ) exit(3);
80  d4 += d3; if ( !test(d4, 3.0, 2.0, 1.0, FEPS) ) exit(3);
81  d4 += f3; if ( !test(d4, 4.0, 3.0, 2.0, FEPS) ) exit(3);
82  f4 = d1 + d2; if ( !test(f4, 2.0, 1.0, 0.0, FEPS) ) exit(3);
83  f4 = f1 + d2; if ( !test(f4, 2.0, 1.0, 0.0, FEPS) ) exit(3);
84  f4 = d1 + f2; if ( !test(f4, 2.0, 1.0, 0.0, FEPS) ) exit(3);
85  f4 = f1 + f2; if ( !test(f4, 2.0, 1.0, 0.0, FEPS) ) exit(3);
86  f4 += d3; if ( !test(f4, 3.0, 2.0, 1.0, FEPS) ) exit(3);
87  f4 += f3; if ( !test(f4, 4.0, 3.0, 2.0, FEPS) ) exit(3);
88 
89 // test subtraction
90 
91  d4 -= d3; if ( !test(d4, 3.0, 2.0, 1.0, FEPS) ) exit(4);
92  d4 -= f3; if ( !test(d4, 2.0, 1.0, 0.0, FEPS) ) exit(4);
93  f4 -= d3; if ( !test(f4, 3.0, 2.0, 1.0, FEPS) ) exit(4);
94  f4 -= f3; if ( !test(f4, 2.0, 1.0, 0.0, FEPS) ) exit(4);
95  d4 = d1 - d2; if ( !test(d4, 0.0, -1.0, 0.0, DEPS) ) exit(4);
96  d4 = f1 - d2; if ( !test(d4, 0.0, -1.0, 0.0, FEPS) ) exit(4);
97  d4 = d1 - f2; if ( !test(d4, 0.0, -1.0, 0.0, FEPS) ) exit(4);
98  d4 = f1 - f2; if ( !test(d4, 0.0, -1.0, 0.0, FEPS) ) exit(4);
99  f4 = d1 - d2; if ( !test(f4, 0.0, -1.0, 0.0, FEPS) ) exit(4);
100  f4 = f1 - d2; if ( !test(f4, 0.0, -1.0, 0.0, FEPS) ) exit(4);
101  f4 = d1 - f2; if ( !test(f4, 0.0, -1.0, 0.0, FEPS) ) exit(4);
102  f4 = f1 - f2; if ( !test(f4, 0.0, -1.0, 0.0, FEPS) ) exit(4);
103 
104 // test unary minus:
105 
106  if ( !test(-d3, -1.0, -1.0, -1.0, DEPS) ) exit(5);
107  if ( !test(-f3, -1.0, -1.0, -1.0, FEPS) ) exit(5);
108  if ( !test(-d1, -1.0, 0.0, 0.0, DEPS) ) exit(5);
109  if ( !test(-f1, -1.0, 0.0, 0.0, FEPS) ) exit(5);
110 
111 // test scaling:
112 
113  if ( !test(d3*2.0, 2.0, 2.0, 2.0, DEPS) ) exit(6);
114  if ( !test(2.0*d3, 2.0, 2.0, 2.0, DEPS) ) exit(6);
115  if ( !test(d1*2.0, 2.0, 0.0, 0.0, DEPS) ) exit(6);
116  if ( !test(2.0*d1, 2.0, 0.0, 0.0, DEPS) ) exit(6);
117  if ( !test(f3*2.0f, 2.0, 2.0, 2.0, FEPS) ) exit(6);
118  if ( !test(2.0f*f3, 2.0, 2.0, 2.0, FEPS) ) exit(6);
119  if ( !test(f1*2.0f, 2.0, 0.0, 0.0, FEPS) ) exit(6);
120  if ( !test(2.0f*f1, 2.0, 0.0, 0.0, FEPS) ) exit(6);
121  if ( !test(d4*=2.0, 0.0, -2.0, 0.0, FEPS) ) exit(6);
122  if ( !test(f4*=2.0, 0.0, -2.0, 0.0, FEPS) ) exit(6);
123 
124 // testing scalar and vector product:
125 
126  if ( !approx(d4*d1, 0.0, DEPS) ) exit(7);
127  if ( !approx(d4*f1, 0.0, FEPS) ) exit(7);
128  if ( !approx(f4*d1, 0.0, FEPS) ) exit(7);
129  if ( !approx(f4*f1, 0.0, FEPS) ) exit(7);
130  if ( !approx(d4.dot(d1), 0.0, DEPS) ) exit(7);
131  if ( !approx(d4.dot(f1), 0.0, FEPS) ) exit(7);
132  if ( !approx(f4.dot(d1), 0.0, FEPS) ) exit(7);
133  if ( !approx(f4.dot(f1), 0.0, FEPS) ) exit(7);
134  if ( !approx(d4*d2, -2.0, DEPS) ) exit(7);
135  if ( !approx(d4*f2, -2.0, FEPS) ) exit(7);
136  if ( !approx(f4*d2, -2.0, FEPS) ) exit(7);
137  if ( !approx(f4*f2, -2.0, FEPS) ) exit(7);
138  if ( !approx(d4.dot(d2), -2.0, DEPS) ) exit(7);
139  if ( !approx(d4.dot(f2), -2.0, FEPS) ) exit(7);
140  if ( !approx(f4.dot(d2), -2.0, FEPS) ) exit(7);
141  if ( !approx(f4.dot(f2), -2.0, FEPS) ) exit(7);
142  d4 = d1.cross(d2); if ( !test(d4, 0.0, 0.0, 1.0, DEPS) ) exit(7);
143  d4 = d2.cross(d1); if ( !test(d4, 0.0, 0.0, -1.0, DEPS) ) exit(7);
144  f4 = f1.cross(d2); if ( !test(f4, 0.0, 0.0, 1.0, FEPS) ) exit(7);
145  f4 = d2.cross(f1); if ( !test(f4, 0.0, 0.0, -1.0, FEPS) ) exit(7);
146 
147 // testing ptot and pt:
148 
149  d4 = d1 + f2 + d3;
150  f4 = d1 + f2 + d3;
151  if ( !approx(d4.mag2(), 14.0, FEPS) ) exit(8);
152  if ( !approx(d4.mag(), std::sqrt(14.0), FEPS) ) exit(8);
153  if ( !approx(d4.perp2(), 13.0, FEPS) ) exit(8);
154  if ( !approx(d4.perp(), std::sqrt(13.0), FEPS) ) exit(8);
155  if ( !approx(f4.mag2(), 14.0, FEPS) ) exit(8);
156  if ( !approx(f4.mag(), std::sqrt(14.0), FEPS) ) exit(8);
157  if ( !approx(f4.perp2(), 13.0, FEPS) ) exit(8);
158  if ( !approx(f4.perp(), std::sqrt(13.0), FEPS) ) exit(8);
159 
160 // testing angles:
161 
162  d4 = d2 - 2.0 * d1;
163  f4 = d2 - 2.0f * f1;
164  if ( !approx(d1.phi(), 0.0, DEPS) ) exit(9);
165  if ( !approx(d1.theta(), CLHEP::halfpi, DEPS) ) exit(9);
166  if ( !approx(d1.cosTheta(), 0.0, DEPS) ) exit(9);
167  if ( !approx(d2.phi(), CLHEP::halfpi*0.5, DEPS) ) exit(9);
168  if ( !approx(d2.theta(), CLHEP::halfpi, DEPS) ) exit(9);
169  if ( !approx(d2.cosTheta(), 0.0, DEPS) ) exit(9);
170  if ( !approx((-d2).phi(), -3.0*CLHEP::halfpi*0.5, DEPS) ) exit(9);
171  if ( !approx(d4.phi(), 3.0*CLHEP::halfpi*0.5, DEPS) ) exit(9);
172 
173  if ( !approx(f1.phi(), 0.0, FEPS) ) exit(9);
174  if ( !approx(f1.theta(), CLHEP::halfpi, FEPS) ) exit(9);
175  if ( !approx(f1.cosTheta(), 0.0, FEPS) ) exit(9);
176  if ( !approx(f2.phi(), CLHEP::halfpi*0.5, FEPS) ) exit(9);
177  if ( !approx(f2.theta(), CLHEP::halfpi, FEPS) ) exit(9);
178  if ( !approx(f2.cosTheta(), 0.0, FEPS) ) exit(9);
179  if ( !approx((-f2).phi(), -3.0*CLHEP::halfpi*0.5, FEPS) ) exit(9);
180  if ( !approx(f4.phi(), 3.0*CLHEP::halfpi*0.5, FEPS) ) exit(9);
181 
182  d4 = d3 - d1; if ( !approx(d4.theta(), CLHEP::halfpi*0.5, DEPS) ) exit(9);
183  if ( !approx((-d4).theta(), 3.0*CLHEP::halfpi*0.5, DEPS) ) exit(9);
184  if ( !approx((-d4).cosTheta(), -std::sqrt(0.5), DEPS) ) exit(9);
185  d4 = d3 - d2; if ( !approx(d4.theta(), 0.0, DEPS) ) exit(9);
186  if ( !approx(d4.cosTheta(), 1.0, DEPS) ) exit(9);
187  if ( !approx((-d4).theta(), CLHEP::pi, DEPS) ) exit(9);
188  if ( !approx((-d4).cosTheta(), -1.0, DEPS) ) exit(9);
189  f4 = d3 - d1; if ( !approx(f4.theta(), CLHEP::halfpi*0.5, FEPS) ) exit(9);
190  if ( !approx((-f4).theta(), 3.0*CLHEP::halfpi*0.5, FEPS) ) exit(9);
191  if ( !approx((-f4).cosTheta(), -std::sqrt(0.5), FEPS) ) exit(9);
192  f4 = d3 - d2; if ( !approx(f4.theta(), 0.0, FEPS) ) exit(9);
193  if ( !approx(f4.cosTheta(), 1.0, FEPS) ) exit(9);
194  if ( !approx((-f4).theta(), CLHEP::pi, FEPS) ) exit(9);
195  if ( !approx((-f4).cosTheta(), -1.0, FEPS) ) exit(9);
196 
197  d4 = d2 - 2.0*d1; if ( !approx(d4.angle(d2), CLHEP::halfpi, DEPS) ) exit(9);
198  f4 = d2 - 2.0*d1; if ( !approx(f4.angle(f2), CLHEP::halfpi, FEPS) ) exit(9);
199 
200 // testing rotations
201 
202  d4 = d1;
203  d4.rotateZ(CLHEP::halfpi); if ( !test(d4, 0.0, 1.0, 0.0, DEPS) ) exit(10);
204  d4.rotateY(25.3); if ( !test(d4, 0.0, 1.0, 0.0, DEPS) ) exit(10);
205  d4.rotateZ(CLHEP::halfpi); if ( !test(d4, -1.0, 0.0, 0.0, DEPS) ) exit(10);
206  d4.rotateY(CLHEP::halfpi); if ( !test(d4, 0.0, 0.0, 1.0, DEPS) ) exit(10);
207  d4.rotateZ(2.6); if ( !test(d4, 0.0, 0.0, 1.0, DEPS) ) exit(10);
208  d4.rotateY(CLHEP::pi*0.25);
209  if ( !test(d4, std::sqrt(0.5), 0.0, std::sqrt(0.5), DEPS) ) exit(10);
210  f4 = f1;
211  f4.rotateZ(CLHEP::halfpi); if ( !test(f4, 0.0, 1.0, 0.0, FEPS) ) exit(10);
212  f4.rotateY(25.3); if ( !test(f4, 0.0, 1.0, 0.0, FEPS) ) exit(10);
213  f4.rotateZ(CLHEP::halfpi); if ( !test(f4, -1.0, 0.0, 0.0, FEPS) ) exit(10);
214  f4.rotateY(CLHEP::halfpi); if ( !test(f4, 0.0, 0.0, 1.0, FEPS) ) exit(10);
215  f4.rotateZ(2.6); if ( !test(f4, 0.0, 0.0, 1.0, FEPS) ) exit(10);
216  f4.rotateY(CLHEP::pi*0.25);
217  if ( !test(f4, std::sqrt(0.5), 0.0, std::sqrt(0.5), FEPS) ) exit(10);
218 
219  d4 = d1;
220  d4.rotate(d4.angle(d3), d4.cross(d3));
221  d4 *= d3.mag();
222  if ( !test(d4, 1.0, 1.0, 1.0, DEPS) ) exit(10);
223  d4 = d1;
224  d4.rotate(0.23, d4.cross(d3));
225  if ( !approx(d4.angle(d1), 0.23, DEPS) ) exit(10);
226  f4 = d1;
227  f4.rotate(f4.angle(d3), f4.cross(d3));
228  f4 *= f3.mag();
229  if ( !test(f4, 1.0, 1.0, 1.0, FEPS) ) exit(10);
230  f4 = f1;
231  f4.rotate(0.23, f4.cross(d3));
232  if ( !approx(f4.angle(f1), 0.23, FEPS) ) exit(10);
233  if ( !approx(f4.angle(d3), f1.angle(d3) - 0.23, FEPS) ) exit(10);
234 
235 // test rotation maticies:
236 
237  d4 = d1;
238 
239  HepRotation r0, r1, r2, r3, r4, r5;
240  r1.rotateZ(CLHEP::halfpi);
241  r2.rotateY(CLHEP::halfpi);
242  r4.rotate(d4.angle(d3), d4.cross(d3));
243  r5.rotate(0.23, d4.cross(d3));
244  d4 = r4.inverse() * d3;
245  if ( !test(d4, d3.mag(), 0.0, 0.0, DEPS) ) exit(11);
246  d4 = r5 * d3;
247  if ( !approx(d1.angle(d4), d1.angle(d3)+0.23, DEPS) ) exit(11);
248  f4 = r4.inverse() * f3;
249  if ( !test(f4, f3.mag(), 0.0, 0.0, FEPS) ) exit(11);
250  f4 = r5 * d3;
251  if ( !approx(d1.angle(f4), f1.angle(f3)+0.23, FEPS) ) exit(11);
252  r5 = r2 * r1 * r3.inverse() * r0 * r0.inverse();
253  d4 = d3;
254  d4 *= r3.inverse();
255  d4 *= r1;
256  d4 *= r2;
257  if ( !test(d4, 1.0, 1.0, 1.0, DEPS) ) exit(11);
258  r5.invert();
259  d4 = r5 * d4;
260  if ( !test(d4, 1.0, 1.0, 1.0, DEPS) ) exit(11);
261  d1 = d2 = Hep3Vector(1.0, -0.5, 2.1);
262  d3 = Hep3Vector(-0.3, 1.1, 1.5);
263  d4 = d3.unit();
264  d4 *= d3.mag();
265  if ( !test(d4, d3.x(), d3.y(), d3.z(), DEPS) ) exit(11);
266  r0.rotate(0.10, d1.cross(d3));
267  d1 *= r0;
268  if ( !approx(d1.angle(d3), d2.angle(d3)-0.1, DEPS) ) exit(12);
269  if ( !approx(d1.angle(d2), 0.1, DEPS) ) exit(12);
270 
271  return 0;
272 
273 }
double x() const
double mag2() const
bool test2(const Hep2Vector &p, double x, double y, double eps)
#define FEPS
Hep3Vector cross(const Hep3Vector &) const
void(* f1)()
double perp() const
double y() const
double cosTheta() const
double x() const
#define exit(x)
int(* f2)(int)
double z() const
double phi() const
#define DEPS
Hep3Vector & rotate(double, const Hep3Vector &)
Definition: ThreeVectorR.cc:29
Hep3Vector & rotateZ(double)
Definition: ThreeVector.cc:144
int(* f3)(int, bool)
HepRotation & rotateY(double delta)
Definition: Rotation.cc:79
void f(void g())
Definition: excDblThrow.cc:38
bool approx(double a, double b, double eps)
double angle(const Hep3Vector &) const
HepRotation & rotateZ(double delta)
Definition: Rotation.cc:92
HepRotation & rotate(double delta, const Hep3Vector &axis)
Definition: Rotation.cc:48
double theta() const
HepRotation inverse() const
double mag() const
double y() const
double dot(const Hep3Vector &) const
Hep3Vector & rotateY(double)
Definition: ThreeVector.cc:134
bool test(const Hep3Vector &p, double x, double y, double z, double eps)
double perp2() const
HepRotation & invert()
int main()