CLHEP VERSION Reference Documentation
CLHEP Home Page CLHEP Documentation CLHEP Bug Reports |
00001 // -*- C++ -*- 00002 // $Id: testThreeVector.cc,v 1.3 2003/08/08 13:47:09 garren Exp $ 00003 // --------------------------------------------------------------------------- 00004 // 00005 // This file is a part of the CLHEP - a Class Library for High Energy Physics. 00006 // 00007 // This is a small program for testing the Hep3Vector class 00008 // and the interaction with the HepRotation class. 00009 // 00010 00011 #include "CLHEP/Units/GlobalSystemOfUnits.h" // to see shadowing problems 00012 #include "CLHEP/Units/GlobalPhysicalConstants.h" 00013 #include "CLHEP/Vector/ThreeVector.h" 00014 #include "CLHEP/Vector/TwoVector.h" 00015 #include "CLHEP/Vector/Rotation.h" 00016 00017 #include <cmath> 00018 #include <iostream> 00019 #include <stdlib.h> // for exit 00020 00021 using namespace CLHEP; 00022 00023 00024 #define DEPS 1.0e-14 00025 #define FEPS 1.0e-6 00026 00027 bool approx(double a, double b, double eps) { 00028 return bool( std::abs(a-b) < eps ); 00029 } 00030 00031 bool 00032 test(const Hep3Vector & p, double x, double y, double z, 00033 double eps) { 00034 return bool( approx(p.x(), x, eps) && approx(p.y(), y, eps) && 00035 approx(p.z(), z, eps) ); 00036 } 00037 00038 bool 00039 test2(const Hep2Vector & p, double x, double y, double eps) { 00040 return bool( approx(p.x(), x, eps) && approx(p.y(), y, eps) ); 00041 } 00042 00043 int main () { 00044 00045 // test constructors: 00046 00047 Hep3Vector d0; if ( !test(d0, 0.0, 0.0, 0.0, DEPS) ) exit(1); 00048 Hep3Vector f0; if ( !test(f0, 0.0, 0.0, 0.0, FEPS) ) exit(1); 00049 Hep3Vector d1(1.0); if ( !test(d1, 1.0, 0.0, 0.0, DEPS) ) exit(1); 00050 Hep3Vector f1(1.0); if ( !test(f1, 1.0, 0.0, 0.0, FEPS) ) exit(1); 00051 Hep3Vector d2(1.0, 1.0); if ( !test(d2, 1.0, 1.0, 0.0, DEPS) ) exit(1); 00052 Hep3Vector f2(1.0, 1.0); if ( !test(f2, 1.0, 1.0, 0.0, FEPS) ) exit(1); 00053 Hep3Vector d3(1.0, 1.0, 1.0); if ( !test(d3, 1.0, 1.0, 1.0, DEPS) ) exit(1); 00054 Hep3Vector f3(1.0, 1.0, 1.0); if ( !test(f3, 1.0, 1.0, 1.0, FEPS) ) exit(1); 00055 Hep3Vector d4(f3); if ( !test(d4, 1.0, 1.0, 1.0, DEPS) ) exit(1); 00056 Hep3Vector f4(d3); if ( !test(f4, 1.0, 1.0, 1.0, FEPS) ) exit(1); 00057 00058 Hep2Vector t3(d4); if ( !test2(t3, 1.0, 1.0, DEPS) ) exit(1); 00059 00060 // test input/output from a stream 00061 00062 std::cin >> d0; if ( !test(d0, 1.1, 2.2, 3.3, DEPS) ) exit(1); 00063 std::cin >> f0; if ( !test(f0, 3.0, 2.0, 1.0, FEPS) ) exit(1); 00064 std::cout << d0 << std::endl; 00065 std::cout << f0 << std::endl; 00066 00067 // test assignment: 00068 00069 d4 = d1; if ( !test(d4, 1.0, 0.0, 0.0, DEPS) ) exit(2); 00070 f4 = f1; if ( !test(f4, 1.0, 0.0, 0.0, FEPS) ) exit(2); 00071 d4 = f1; if ( !test(d4, 1.0, 0.0, 0.0, FEPS) ) exit(2); 00072 f4 = d1; if ( !test(f4, 1.0, 0.0, 0.0, FEPS) ) exit(2); 00073 00074 // test addition: 00075 00076 d4 = d1 + d2; if ( !test(d4, 2.0, 1.0, 0.0, DEPS) ) exit(3); 00077 d4 = f1 + d2; if ( !test(d4, 2.0, 1.0, 0.0, FEPS) ) exit(3); 00078 d4 = d1 + f2; if ( !test(d4, 2.0, 1.0, 0.0, FEPS) ) exit(3); 00079 d4 = f1 + f2; if ( !test(d4, 2.0, 1.0, 0.0, FEPS) ) exit(3); 00080 d4 += d3; if ( !test(d4, 3.0, 2.0, 1.0, FEPS) ) exit(3); 00081 d4 += f3; if ( !test(d4, 4.0, 3.0, 2.0, FEPS) ) exit(3); 00082 f4 = d1 + d2; if ( !test(f4, 2.0, 1.0, 0.0, FEPS) ) exit(3); 00083 f4 = f1 + d2; if ( !test(f4, 2.0, 1.0, 0.0, FEPS) ) exit(3); 00084 f4 = d1 + f2; if ( !test(f4, 2.0, 1.0, 0.0, FEPS) ) exit(3); 00085 f4 = f1 + f2; if ( !test(f4, 2.0, 1.0, 0.0, FEPS) ) exit(3); 00086 f4 += d3; if ( !test(f4, 3.0, 2.0, 1.0, FEPS) ) exit(3); 00087 f4 += f3; if ( !test(f4, 4.0, 3.0, 2.0, FEPS) ) exit(3); 00088 00089 // test subtraction 00090 00091 d4 -= d3; if ( !test(d4, 3.0, 2.0, 1.0, FEPS) ) exit(4); 00092 d4 -= f3; if ( !test(d4, 2.0, 1.0, 0.0, FEPS) ) exit(4); 00093 f4 -= d3; if ( !test(f4, 3.0, 2.0, 1.0, FEPS) ) exit(4); 00094 f4 -= f3; if ( !test(f4, 2.0, 1.0, 0.0, FEPS) ) exit(4); 00095 d4 = d1 - d2; if ( !test(d4, 0.0, -1.0, 0.0, DEPS) ) exit(4); 00096 d4 = f1 - d2; if ( !test(d4, 0.0, -1.0, 0.0, FEPS) ) exit(4); 00097 d4 = d1 - f2; if ( !test(d4, 0.0, -1.0, 0.0, FEPS) ) exit(4); 00098 d4 = f1 - f2; if ( !test(d4, 0.0, -1.0, 0.0, FEPS) ) exit(4); 00099 f4 = d1 - d2; if ( !test(f4, 0.0, -1.0, 0.0, FEPS) ) exit(4); 00100 f4 = f1 - d2; if ( !test(f4, 0.0, -1.0, 0.0, FEPS) ) exit(4); 00101 f4 = d1 - f2; if ( !test(f4, 0.0, -1.0, 0.0, FEPS) ) exit(4); 00102 f4 = f1 - f2; if ( !test(f4, 0.0, -1.0, 0.0, FEPS) ) exit(4); 00103 00104 // test unary minus: 00105 00106 if ( !test(-d3, -1.0, -1.0, -1.0, DEPS) ) exit(5); 00107 if ( !test(-f3, -1.0, -1.0, -1.0, FEPS) ) exit(5); 00108 if ( !test(-d1, -1.0, 0.0, 0.0, DEPS) ) exit(5); 00109 if ( !test(-f1, -1.0, 0.0, 0.0, FEPS) ) exit(5); 00110 00111 // test scaling: 00112 00113 if ( !test(d3*2.0, 2.0, 2.0, 2.0, DEPS) ) exit(6); 00114 if ( !test(2.0*d3, 2.0, 2.0, 2.0, DEPS) ) exit(6); 00115 if ( !test(d1*2.0, 2.0, 0.0, 0.0, DEPS) ) exit(6); 00116 if ( !test(2.0*d1, 2.0, 0.0, 0.0, DEPS) ) exit(6); 00117 if ( !test(f3*2.0f, 2.0, 2.0, 2.0, FEPS) ) exit(6); 00118 if ( !test(2.0f*f3, 2.0, 2.0, 2.0, FEPS) ) exit(6); 00119 if ( !test(f1*2.0f, 2.0, 0.0, 0.0, FEPS) ) exit(6); 00120 if ( !test(2.0f*f1, 2.0, 0.0, 0.0, FEPS) ) exit(6); 00121 if ( !test(d4*=2.0, 0.0, -2.0, 0.0, FEPS) ) exit(6); 00122 if ( !test(f4*=2.0, 0.0, -2.0, 0.0, FEPS) ) exit(6); 00123 00124 // testing scalar and vector product: 00125 00126 if ( !approx(d4*d1, 0.0, DEPS) ) exit(7); 00127 if ( !approx(d4*f1, 0.0, FEPS) ) exit(7); 00128 if ( !approx(f4*d1, 0.0, FEPS) ) exit(7); 00129 if ( !approx(f4*f1, 0.0, FEPS) ) exit(7); 00130 if ( !approx(d4.dot(d1), 0.0, DEPS) ) exit(7); 00131 if ( !approx(d4.dot(f1), 0.0, FEPS) ) exit(7); 00132 if ( !approx(f4.dot(d1), 0.0, FEPS) ) exit(7); 00133 if ( !approx(f4.dot(f1), 0.0, FEPS) ) exit(7); 00134 if ( !approx(d4*d2, -2.0, DEPS) ) exit(7); 00135 if ( !approx(d4*f2, -2.0, FEPS) ) exit(7); 00136 if ( !approx(f4*d2, -2.0, FEPS) ) exit(7); 00137 if ( !approx(f4*f2, -2.0, FEPS) ) exit(7); 00138 if ( !approx(d4.dot(d2), -2.0, DEPS) ) exit(7); 00139 if ( !approx(d4.dot(f2), -2.0, FEPS) ) exit(7); 00140 if ( !approx(f4.dot(d2), -2.0, FEPS) ) exit(7); 00141 if ( !approx(f4.dot(f2), -2.0, FEPS) ) exit(7); 00142 d4 = d1.cross(d2); if ( !test(d4, 0.0, 0.0, 1.0, DEPS) ) exit(7); 00143 d4 = d2.cross(d1); if ( !test(d4, 0.0, 0.0, -1.0, DEPS) ) exit(7); 00144 f4 = f1.cross(d2); if ( !test(f4, 0.0, 0.0, 1.0, FEPS) ) exit(7); 00145 f4 = d2.cross(f1); if ( !test(f4, 0.0, 0.0, -1.0, FEPS) ) exit(7); 00146 00147 // testing ptot and pt: 00148 00149 d4 = d1 + f2 + d3; 00150 f4 = d1 + f2 + d3; 00151 if ( !approx(d4.mag2(), 14.0, FEPS) ) exit(8); 00152 if ( !approx(d4.mag(), std::sqrt(14.0), FEPS) ) exit(8); 00153 if ( !approx(d4.perp2(), 13.0, FEPS) ) exit(8); 00154 if ( !approx(d4.perp(), std::sqrt(13.0), FEPS) ) exit(8); 00155 if ( !approx(f4.mag2(), 14.0, FEPS) ) exit(8); 00156 if ( !approx(f4.mag(), std::sqrt(14.0), FEPS) ) exit(8); 00157 if ( !approx(f4.perp2(), 13.0, FEPS) ) exit(8); 00158 if ( !approx(f4.perp(), std::sqrt(13.0), FEPS) ) exit(8); 00159 00160 // testing angles: 00161 00162 d4 = d2 - 2.0 * d1; 00163 f4 = d2 - 2.0f * f1; 00164 if ( !approx(d1.phi(), 0.0, DEPS) ) exit(9); 00165 if ( !approx(d1.theta(), CLHEP::halfpi, DEPS) ) exit(9); 00166 if ( !approx(d1.cosTheta(), 0.0, DEPS) ) exit(9); 00167 if ( !approx(d2.phi(), CLHEP::halfpi*0.5, DEPS) ) exit(9); 00168 if ( !approx(d2.theta(), CLHEP::halfpi, DEPS) ) exit(9); 00169 if ( !approx(d2.cosTheta(), 0.0, DEPS) ) exit(9); 00170 if ( !approx((-d2).phi(), -3.0*CLHEP::halfpi*0.5, DEPS) ) exit(9); 00171 if ( !approx(d4.phi(), 3.0*CLHEP::halfpi*0.5, DEPS) ) exit(9); 00172 00173 if ( !approx(f1.phi(), 0.0, FEPS) ) exit(9); 00174 if ( !approx(f1.theta(), CLHEP::halfpi, FEPS) ) exit(9); 00175 if ( !approx(f1.cosTheta(), 0.0, FEPS) ) exit(9); 00176 if ( !approx(f2.phi(), CLHEP::halfpi*0.5, FEPS) ) exit(9); 00177 if ( !approx(f2.theta(), CLHEP::halfpi, FEPS) ) exit(9); 00178 if ( !approx(f2.cosTheta(), 0.0, FEPS) ) exit(9); 00179 if ( !approx((-f2).phi(), -3.0*CLHEP::halfpi*0.5, FEPS) ) exit(9); 00180 if ( !approx(f4.phi(), 3.0*CLHEP::halfpi*0.5, FEPS) ) exit(9); 00181 00182 d4 = d3 - d1; if ( !approx(d4.theta(), CLHEP::halfpi*0.5, DEPS) ) exit(9); 00183 if ( !approx((-d4).theta(), 3.0*CLHEP::halfpi*0.5, DEPS) ) exit(9); 00184 if ( !approx((-d4).cosTheta(), -std::sqrt(0.5), DEPS) ) exit(9); 00185 d4 = d3 - d2; if ( !approx(d4.theta(), 0.0, DEPS) ) exit(9); 00186 if ( !approx(d4.cosTheta(), 1.0, DEPS) ) exit(9); 00187 if ( !approx((-d4).theta(), CLHEP::pi, DEPS) ) exit(9); 00188 if ( !approx((-d4).cosTheta(), -1.0, DEPS) ) exit(9); 00189 f4 = d3 - d1; if ( !approx(f4.theta(), CLHEP::halfpi*0.5, FEPS) ) exit(9); 00190 if ( !approx((-f4).theta(), 3.0*CLHEP::halfpi*0.5, FEPS) ) exit(9); 00191 if ( !approx((-f4).cosTheta(), -std::sqrt(0.5), FEPS) ) exit(9); 00192 f4 = d3 - d2; if ( !approx(f4.theta(), 0.0, FEPS) ) exit(9); 00193 if ( !approx(f4.cosTheta(), 1.0, FEPS) ) exit(9); 00194 if ( !approx((-f4).theta(), CLHEP::pi, FEPS) ) exit(9); 00195 if ( !approx((-f4).cosTheta(), -1.0, FEPS) ) exit(9); 00196 00197 d4 = d2 - 2.0*d1; if ( !approx(d4.angle(d2), CLHEP::halfpi, DEPS) ) exit(9); 00198 f4 = d2 - 2.0*d1; if ( !approx(f4.angle(f2), CLHEP::halfpi, FEPS) ) exit(9); 00199 00200 // testing rotations 00201 00202 d4 = d1; 00203 d4.rotateZ(CLHEP::halfpi); if ( !test(d4, 0.0, 1.0, 0.0, DEPS) ) exit(10); 00204 d4.rotateY(25.3); if ( !test(d4, 0.0, 1.0, 0.0, DEPS) ) exit(10); 00205 d4.rotateZ(CLHEP::halfpi); if ( !test(d4, -1.0, 0.0, 0.0, DEPS) ) exit(10); 00206 d4.rotateY(CLHEP::halfpi); if ( !test(d4, 0.0, 0.0, 1.0, DEPS) ) exit(10); 00207 d4.rotateZ(2.6); if ( !test(d4, 0.0, 0.0, 1.0, DEPS) ) exit(10); 00208 d4.rotateY(CLHEP::pi*0.25); 00209 if ( !test(d4, std::sqrt(0.5), 0.0, std::sqrt(0.5), DEPS) ) exit(10); 00210 f4 = f1; 00211 f4.rotateZ(CLHEP::halfpi); if ( !test(f4, 0.0, 1.0, 0.0, FEPS) ) exit(10); 00212 f4.rotateY(25.3); if ( !test(f4, 0.0, 1.0, 0.0, FEPS) ) exit(10); 00213 f4.rotateZ(CLHEP::halfpi); if ( !test(f4, -1.0, 0.0, 0.0, FEPS) ) exit(10); 00214 f4.rotateY(CLHEP::halfpi); if ( !test(f4, 0.0, 0.0, 1.0, FEPS) ) exit(10); 00215 f4.rotateZ(2.6); if ( !test(f4, 0.0, 0.0, 1.0, FEPS) ) exit(10); 00216 f4.rotateY(CLHEP::pi*0.25); 00217 if ( !test(f4, std::sqrt(0.5), 0.0, std::sqrt(0.5), FEPS) ) exit(10); 00218 00219 d4 = d1; 00220 d4.rotate(d4.angle(d3), d4.cross(d3)); 00221 d4 *= d3.mag(); 00222 if ( !test(d4, 1.0, 1.0, 1.0, DEPS) ) exit(10); 00223 d4 = d1; 00224 d4.rotate(0.23, d4.cross(d3)); 00225 if ( !approx(d4.angle(d1), 0.23, DEPS) ) exit(10); 00226 f4 = d1; 00227 f4.rotate(f4.angle(d3), f4.cross(d3)); 00228 f4 *= f3.mag(); 00229 if ( !test(f4, 1.0, 1.0, 1.0, FEPS) ) exit(10); 00230 f4 = f1; 00231 f4.rotate(0.23, f4.cross(d3)); 00232 if ( !approx(f4.angle(f1), 0.23, FEPS) ) exit(10); 00233 if ( !approx(f4.angle(d3), f1.angle(d3) - 0.23, FEPS) ) exit(10); 00234 00235 // test rotation maticies: 00236 00237 d4 = d1; 00238 00239 HepRotation r0, r1, r2, r3, r4, r5; 00240 r1.rotateZ(CLHEP::halfpi); 00241 r2.rotateY(CLHEP::halfpi); 00242 r4.rotate(d4.angle(d3), d4.cross(d3)); 00243 r5.rotate(0.23, d4.cross(d3)); 00244 d4 = r4.inverse() * d3; 00245 if ( !test(d4, d3.mag(), 0.0, 0.0, DEPS) ) exit(11); 00246 d4 = r5 * d3; 00247 if ( !approx(d1.angle(d4), d1.angle(d3)+0.23, DEPS) ) exit(11); 00248 f4 = r4.inverse() * f3; 00249 if ( !test(f4, f3.mag(), 0.0, 0.0, FEPS) ) exit(11); 00250 f4 = r5 * d3; 00251 if ( !approx(d1.angle(f4), f1.angle(f3)+0.23, FEPS) ) exit(11); 00252 r5 = r2 * r1 * r3.inverse() * r0 * r0.inverse(); 00253 d4 = d3; 00254 d4 *= r3.inverse(); 00255 d4 *= r1; 00256 d4 *= r2; 00257 if ( !test(d4, 1.0, 1.0, 1.0, DEPS) ) exit(11); 00258 r5.invert(); 00259 d4 = r5 * d4; 00260 if ( !test(d4, 1.0, 1.0, 1.0, DEPS) ) exit(11); 00261 d1 = d2 = Hep3Vector(1.0, -0.5, 2.1); 00262 d3 = Hep3Vector(-0.3, 1.1, 1.5); 00263 d4 = d3.unit(); 00264 d4 *= d3.mag(); 00265 if ( !test(d4, d3.x(), d3.y(), d3.z(), DEPS) ) exit(11); 00266 r0.rotate(0.10, d1.cross(d3)); 00267 d1 *= r0; 00268 if ( !approx(d1.angle(d3), d2.angle(d3)-0.1, DEPS) ) exit(12); 00269 if ( !approx(d1.angle(d2), 0.1, DEPS) ) exit(12); 00270 00271 return 0; 00272 00273 }