CLHEP 2.0.4.7 Reference Documentation
CLHEP Home Page CLHEP Documentation CLHEP Bug Reports |
00001 // -*- C++ -*- 00002 // $Id: testLorentzVector.cc,v 1.2 2003/08/13 20:00:14 garren Exp $ 00003 // --------------------------------------------------------------------------- 00004 // 00005 // This file is a part of what might become CLHEP - 00006 // a Class Library for High Energy Physics. 00007 // 00008 // This is a small program for testing the HepLorentzVector class 00009 // and the interaction with the HepLorentzRotation class. 00010 // 00011 00012 #include "CLHEP/Vector/defs.h" 00013 #include "CLHEP/Vector/LorentzVector.h" 00014 #include "CLHEP/Vector/LorentzRotation.h" 00015 #include "CLHEP/Vector/Sqr.h" 00016 #include <iostream> 00017 #include <cmath> 00018 #include <stdlib.h> 00019 #include <assert.h> 00020 00021 using namespace CLHEP; 00022 00023 #define DEPS 1.0e-14 00024 #define FEPS 1.0e-6 00025 00026 bool approx(double a, double b, double eps) { 00027 return bool( std::abs(a-b) < eps ); 00028 } 00029 00030 bool 00031 test(const HepLorentzVector & p, double x, double y, double z, double e, 00032 double eps) { 00033 bool t = bool( approx(p.x(), x, eps) && approx(p.y(), y, eps) && 00034 approx(p.z(), z, eps) && approx(p.t(), e, eps)); 00035 if ( !t ) std::cerr << p << std::endl 00036 << x << '\t' << y << '\t' << z << '\t' << e 00037 << std::endl; 00038 return t; 00039 } 00040 00041 void conversion_test(Hep3Vector & v3, HepLorentzVector & v4) { 00042 v3 = Hep3Vector(3.,2.,1.); 00043 assert (v3.x() == v4.x() && v3.y() == v4.y() && v3.z() == v4.z()); 00044 } 00045 00046 void conversion_test(const Hep3Vector & v3, const HepLorentzVector & v4) { 00047 assert (v3.x() == v4.x() && v3.y() == v4.y() && v3.z() == v4.z()); 00048 } 00049 00050 bool 00051 test(const HepLorentzVector & p, const HepLorentzVector & q, double eps) { 00052 bool t = bool( approx(p.x(), q.x(), eps) && 00053 approx(p.y(), q.y(), eps) && 00054 approx(p.z(), q.z(), eps) && 00055 approx(p.t(), q.t(), eps)); 00056 if ( !t ) std::cerr << p << std::endl 00057 << q << std::endl; 00058 return t; 00059 } 00060 00061 int main () { 00062 HepLorentzVector v4(1.,2.,3.,4.); 00063 const HepLorentzVector v4const(1.,2.,3.,4); 00064 conversion_test(v4,v4); 00065 conversion_test(v4const, v4const); 00066 00067 Hep3Vector f3x(1.0), f3y(0.0, 1.0), f3z(0.0, 0.0, 1.0); 00068 Hep3Vector d30, d3x(1.0), d3y(0.0, 1.0), d3z(0.0, 0.0, 1.0); 00069 00070 // test constructors: 00071 00072 HepLorentzVector d0; 00073 if ( !test(d0, 0.0, 0.0, 0.0, 0.0, DEPS) ) exit(1); 00074 HepLorentzVector d1(d3x, 1.0); 00075 if ( !test(d1, 1.0, 0.0, 0.0, 1.0, DEPS) ) exit(1); 00076 HepLorentzVector d2(d3x + d3y, sqrt(2.0)); 00077 if ( !test(d2, 1.0, 1.0, 0.0, sqrt(2.0), DEPS) ) exit(1); 00078 HepLorentzVector d3(d3z + d2, sqrt(3.0)); 00079 if ( !test(d3, 1.0, 1.0, 1.0, sqrt(3.0), DEPS) ) exit(1); 00080 HepLorentzVector d4(0.0, 0.0, 0.0, 1.0); 00081 if ( !test(d4,0.0, 0.0, 0.0, 1.0, DEPS) ) exit(1); 00082 HepLorentzVector d5(f3x, f3x.mag()); if ( !test(d5, d1, FEPS) ) exit(1); 00083 HepLorentzVector d6(d3x+f3y, (d3x+f3y).mag()); 00084 if ( !test(d6, d2, FEPS) ) exit(1); 00085 HepLorentzVector d7(f3x+f3y+f3z, (f3x+f3y+f3z).mag()); 00086 if ( !test(d7, d3, FEPS) ) exit(1); 00087 00088 HepLorentzVector f0; if ( !test(f0, 0.0, 0.0, 0.0, 0.0, FEPS) ) exit(1); 00089 HepLorentzVector f1(f3x, 1.0); 00090 if ( !test(f1, 1.0, 0.0, 0.0, 1.0, FEPS) ) exit(1); 00091 HepLorentzVector f2(f3x + f3y, sqrt(2.0)); 00092 if ( !test(f2, 1.0, 1.0, 0.0, sqrt(2.0), FEPS) ) exit(1); 00093 HepLorentzVector f3(f3z + f2, sqrt(3.0)); 00094 if ( !test(f3, 1.0, 1.0, 1.0, sqrt(3.0), FEPS) ) exit(1); 00095 HepLorentzVector f4(0.0, 0.0, 0.0, 1.0); 00096 if ( !test(f4,0.0, 0.0, 0.0, 1.0, FEPS) ) exit(1); 00097 HepLorentzVector f5(d3x, d3x.mag()); if ( !test(f5, f1, FEPS) ) exit(1); 00098 HepLorentzVector f6(f3x+d3y, (f3x+d3y).mag()); 00099 if ( !test(f6, f2, FEPS) ) exit(1); 00100 HepLorentzVector f7(d3x+d3y+d3z, (d3x+d3y+d3z).mag()); 00101 if ( !test(f7, f3, FEPS) ) exit(1); 00102 00103 HepLorentzVector d8(f7); if ( !test(d8, d7, FEPS) ) exit(1); 00104 HepLorentzVector d9(d7); if ( !test(d9, d7, DEPS) ) exit(1); 00105 HepLorentzVector f8(f7); if ( !test(f8, d7, FEPS) ) exit(1); 00106 HepLorentzVector f9(d7); if ( !test(f9, d7, FEPS) ) exit(1); 00107 00108 HepLorentzVector d10(1.0, 1.0, 1.0, sqrt(3.0)); 00109 if ( !test(d10, d7, FEPS) ) exit(1); 00110 HepLorentzVector f10(1.0, 1.0, 1.0, sqrt(3.0)); 00111 if ( !test(f10, f7, FEPS) ) exit(1); 00112 00113 HepLorentzVector d11(d3x+d3y+d3z, 1.0); 00114 if ( !test(d11, 1.0, 1.0, 1.0, 1.0, DEPS) ) exit(1); 00115 HepLorentzVector f11(d3x+d3y+d3z, 1.0); 00116 if ( !test(f11, 1.0, 1.0, 1.0, 1.0, FEPS) ) exit(1); 00117 00118 // test input/output from a stream 00119 00120 std::cin >> d0; if ( !test(d0, 1.1, 2.2, 3.3, 4.4, DEPS) ) exit(1); 00121 std::cin >> f0; if ( !test(f0, 4.0, 3.0, 2.0, 1.0, FEPS) ) exit(1); 00122 std::cout << d0 << std::endl; 00123 std::cout << f0 << std::endl; 00124 00125 // testing assignment 00126 00127 d6 = d7; if ( !test(d6, d7, DEPS) ) exit(2); 00128 d6 = f7; if ( !test(d6, d7, FEPS) ) exit(2); 00129 f6 = d7; if ( !test(f6, f7, FEPS) ) exit(2); 00130 f6 = f7; if ( !test(f6, f7, FEPS) ) exit(2); 00131 00132 //testing addition and subtraction: 00133 00134 d11 = d3 + d7 + f3; 00135 if ( !test(d11, 3.0, 3.0, 3.0, sqrt(27.0), FEPS) ) exit(4); 00136 f11 = d3 + d7 + f3; 00137 if ( !test(f11, 3.0, 3.0, 3.0, sqrt(27.0), FEPS) ) exit(4); 00138 d11 += d3; 00139 if ( !test(d11, 4.0, 4.0, 4.0, sqrt(48.0), FEPS) ) exit(4); 00140 f11 += f3; 00141 if ( !test(f11, 4.0, 4.0, 4.0, sqrt(48.0), FEPS) ) exit(4); 00142 d11 = d3 + d7 - f3; 00143 if ( !test(d11, 1.0, 1.0, 1.0, sqrt(3.0), FEPS) ) exit(4); 00144 if ( !test(-d11, -1.0, -1.0, -1.0, -sqrt(3.0), FEPS) ) exit(4); 00145 f11 = d3 + f7 - d3; 00146 if ( !test(f11, 1.0, 1.0, 1.0, sqrt(3.0), FEPS) ) exit(4); 00147 if ( !test(-f11, -1.0, -1.0, -1.0, -sqrt(3.0), FEPS) ) exit(4); 00148 d11 -= d3; 00149 if ( !test(d11, 0.0, 0.0, 0.0, 0.0, FEPS) ) exit(4); 00150 f11 -= f3; 00151 if ( !test(f11, 0.0, 0.0, 0.0, 0.0, FEPS) ) exit(4); 00152 00153 d11 = HepLorentzVector(1.0, 2.0, 3.0, 4.0); 00154 d11 *= 2.; 00155 if ( !test(d11, 2.0, 4.0, 6.0, 8.0, DEPS) ) exit(4); 00156 d11 = 2.*HepLorentzVector(1.0, 2.0, 3.0, 4.0); 00157 if ( !test(d11, 2.0, 4.0, 6.0, 8.0, DEPS) ) exit(4); 00158 d11 = HepLorentzVector(1.0, 2.0, 3.0, 4.0)*2.; 00159 if ( !test(d11, 2.0, 4.0, 6.0, 8.0, DEPS) ) exit(4); 00160 00161 // testing scalar products: 00162 00163 if ( !approx(d1 * d2, sqrt(2.0)-1.0, DEPS) ) exit(5); 00164 if ( !approx(d3.dot(d7), 0.0, FEPS) ) exit(5); 00165 if ( !approx(d2 * f1, sqrt(2.0)-1.0, FEPS) ) exit(5); 00166 if ( !approx(f3.dot(d7), 0.0, FEPS) ) exit(5); 00167 00168 // testing components: 00169 00170 d11 = HepLorentzVector(1.0, 1.0, 1.0, sqrt(7.0)); 00171 if ( !approx(d11.mag2(), 4.0, DEPS) ) exit(6); 00172 if ( !approx(d11.mag(), 2.0, DEPS) ) exit(6); 00173 if ( !approx(Hep3Vector(d11).mag2(), 3.0, DEPS) ) exit(6); 00174 if ( !approx(Hep3Vector(d11).mag(), sqrt(3.0), DEPS) ) exit(6); 00175 if ( !approx(d11.perp2(), 2.0, DEPS) ) exit(6); 00176 if ( !approx(d11.perp(), sqrt(2.0), DEPS) ) exit(6); 00177 f11 = HepLorentzVector(1.0, 1.0, 1.0, sqrt(7.0)); 00178 if ( !approx(f11.mag2(), 4.0, FEPS) ) exit(6); 00179 if ( !approx(f11.mag(), 2.0, FEPS) ) exit(6); 00180 if ( !approx(f11.vect().mag2(), 3.0, FEPS) ) exit(6); 00181 if ( !approx(f11.vect().mag(), sqrt(3.0), FEPS) ) exit(6); 00182 if ( !approx(f11.perp2(), 2.0, FEPS) ) exit(6); 00183 if ( !approx(f11.perp(), sqrt(2.0), FEPS) ) exit(6); 00184 00185 // testing boosts: 00186 00187 d5 = d3 = d1 = HepLorentzVector(1.0, 2.0, -1.0, 3.0); 00188 d6 = d4 = d2 = HepLorentzVector(-1.0, 1.0, 2.0, 4.0); 00189 double M = (d1 + d2).mag(); 00190 double m1 = d1.mag(); 00191 double m2 = d2.mag(); 00192 double p2 = (sqr(M)-sqr(m1+m2))*(sqr(M)-sqr(m1-m2))/(4.0*sqr(M)); 00193 d30 = -(d1 + d2).boostVector(); 00194 d1.boost(d30); 00195 double phi = d1.phi(); 00196 double theta = d1.theta(); 00197 d1.rotateZ(-phi); 00198 d1.rotateY(-theta); 00199 HepRotation r; 00200 r.rotateZ(-phi); 00201 HepLorentzRotation r1(d30), r2(r), r3, r4, r5; 00202 r3.rotateY(-theta); 00203 r4 = r3 * r2 * r1; 00204 d2 *= r4; 00205 if ( !test(d1, 0.0, 0.0, sqrt(p2), sqrt(p2 + sqr(m1)), DEPS) ) exit(7); 00206 if ( !test(d2, 0.0, 0.0, -sqrt(p2), sqrt(p2 + sqr(m2)), DEPS) ) exit(7); 00207 d1.transform(r4.inverse()); 00208 if ( !test(d1, d3, DEPS) ) exit(7); 00209 r5 *= r3; 00210 r5 *= r; 00211 r5 *= r1; 00212 r5.invert(); 00213 d2 *= r5; 00214 if ( !test(d2, d4, DEPS) ) exit(7); 00215 r4 = r1; 00216 r4.rotateZ(-phi); 00217 r4.rotateY(-theta); 00218 d3 *= r4; 00219 d4 = r4 * d6; 00220 if ( !test(d3, 0.0, 0.0, sqrt(p2), sqrt(p2 + sqr(m1)), DEPS) ) exit(7); 00221 if ( !test(d4, 0.0, 0.0, -sqrt(p2), sqrt(p2 + sqr(m2)), DEPS) ) exit(7); 00222 r5 = r1.inverse(); 00223 r5 *= r.inverse(); 00224 r5 *= r3.inverse(); 00225 d4.transform(r5); 00226 d3.transform(r5); 00227 00228 if ( !test(d4, d6, DEPS) ) exit(7); 00229 if ( !test(d3, d5, DEPS) ) exit(7); 00230 00231 r5 = r1; 00232 r5.transform(r); 00233 r5.transform(r3); 00234 d4.transform(r5); 00235 d3.transform(r5); 00236 if ( !test(d3, 0.0, 0.0, sqrt(p2), sqrt(p2 + sqr(m1)), DEPS) ) exit(7); 00237 if ( !test(d4, 0.0, 0.0, -sqrt(p2), sqrt(p2 + sqr(m2)), DEPS) ) exit(7); 00238 00239 return 0; 00240 }