Geant4 Cross Reference |
1 // -*- C++ -*- 1 // -*- C++ -*- 2 // ------------------------------------------- 2 // --------------------------------------------------------------------------- 3 // 3 // 4 // This file is a part of the CLHEP - a Class 4 // This file is a part of the CLHEP - a Class Library for High Energy Physics. 5 // 5 // 6 // This is the implementation of methods of th 6 // This is the implementation of methods of the HepRotationX class which 7 // were introduced when ZOOM PhysicsVectors wa 7 // were introduced when ZOOM PhysicsVectors was merged in. 8 // 8 // 9 9 10 #include "CLHEP/Vector/RotationX.h" 10 #include "CLHEP/Vector/RotationX.h" 11 #include "CLHEP/Vector/AxisAngle.h" 11 #include "CLHEP/Vector/AxisAngle.h" 12 #include "CLHEP/Vector/EulerAngles.h" 12 #include "CLHEP/Vector/EulerAngles.h" 13 #include "CLHEP/Vector/LorentzRotation.h" 13 #include "CLHEP/Vector/LorentzRotation.h" 14 #include "CLHEP/Units/PhysicalConstants.h" 14 #include "CLHEP/Units/PhysicalConstants.h" 15 15 16 #include <cmath> 16 #include <cmath> 17 #include <stdlib.h> 17 #include <stdlib.h> 18 #include <iostream> 18 #include <iostream> 19 19 20 namespace CLHEP { 20 namespace CLHEP { 21 21 22 static inline double safe_acos (double x) { 22 static inline double safe_acos (double x) { 23 if (std::abs(x) <= 1.0) return std::acos(x); 23 if (std::abs(x) <= 1.0) return std::acos(x); 24 return ( (x>0) ? 0 : CLHEP::pi ); 24 return ( (x>0) ? 0 : CLHEP::pi ); 25 } 25 } 26 26 27 HepRotationX::HepRotationX(double ddelta) : 27 HepRotationX::HepRotationX(double ddelta) : 28 its_d(proper(ddelta)), its_s(std::sin(ddel 28 its_d(proper(ddelta)), its_s(std::sin(ddelta)), its_c(std::cos(ddelta)) 29 {} 29 {} 30 30 31 HepRotationX & HepRotationX::set ( double ddel 31 HepRotationX & HepRotationX::set ( double ddelta ) { 32 its_d = proper(ddelta); 32 its_d = proper(ddelta); 33 its_s = std::sin(its_d); 33 its_s = std::sin(its_d); 34 its_c = std::cos(its_d); 34 its_c = std::cos(its_d); 35 return *this; 35 return *this; 36 } 36 } 37 37 38 double HepRotationX::phi() const { 38 double HepRotationX::phi() const { 39 if ( (its_d > 0) && (its_d < CLHEP::pi) ) { 39 if ( (its_d > 0) && (its_d < CLHEP::pi) ) { 40 return CLHEP::pi; 40 return CLHEP::pi; 41 } else { 41 } else { 42 return 0.0; 42 return 0.0; 43 } 43 } 44 } // HepRotationX::phi() 44 } // HepRotationX::phi() 45 45 46 double HepRotationX::theta() const { 46 double HepRotationX::theta() const { 47 return std::fabs( its_d ); 47 return std::fabs( its_d ); 48 } // HepRotationX::theta() 48 } // HepRotationX::theta() 49 49 50 double HepRotationX::psi() const { 50 double HepRotationX::psi() const { 51 if ( (its_d > 0) && (its_d < CLHEP::pi) ) { 51 if ( (its_d > 0) && (its_d < CLHEP::pi) ) { 52 return CLHEP::pi; 52 return CLHEP::pi; 53 } else { 53 } else { 54 return 0.0; 54 return 0.0; 55 } 55 } 56 } // HepRotationX::psi() 56 } // HepRotationX::psi() 57 57 58 HepEulerAngles HepRotationX::eulerAngles() con 58 HepEulerAngles HepRotationX::eulerAngles() const { 59 return HepEulerAngles( phi(), theta(), psi 59 return HepEulerAngles( phi(), theta(), psi() ); 60 } // HepRotationX::eulerAngles() 60 } // HepRotationX::eulerAngles() 61 61 62 62 63 // From the defining code in the implementatio 63 // From the defining code in the implementation of CLHEP (in Rotation.cc) 64 // it is clear that thetaX, phiX form the pola 64 // it is clear that thetaX, phiX form the polar angles in the original 65 // coordinate system of the new X axis (and si 65 // coordinate system of the new X axis (and similarly for phiY and phiZ). 66 // 66 // 67 // This code is taken directly from the origin 67 // This code is taken directly from the original CLHEP. However, there are as 68 // shown opportunities for significant speed i 68 // shown opportunities for significant speed improvement. 69 69 70 double HepRotationX::phiX() const { 70 double HepRotationX::phiX() const { 71 return (yx() == 0.0 && xx() == 0.0) ? 0.0 : 71 return (yx() == 0.0 && xx() == 0.0) ? 0.0 : std::atan2(yx(),xx()); 72 // or ---- return 0; 72 // or ---- return 0; 73 } 73 } 74 74 75 double HepRotationX::phiY() const { 75 double HepRotationX::phiY() const { 76 return (yy() == 0.0 && xy() == 0.0) ? 0.0 : 76 return (yy() == 0.0 && xy() == 0.0) ? 0.0 : std::atan2(yy(),xy()); 77 // or ---- return (yy() == 0.0) ? 0.0 : s 77 // or ---- return (yy() == 0.0) ? 0.0 : std::atan2(yy(),xy()); 78 } 78 } 79 79 80 double HepRotationX::phiZ() const { 80 double HepRotationX::phiZ() const { 81 return (yz() == 0.0 && xz() == 0.0) ? 0.0 : 81 return (yz() == 0.0 && xz() == 0.0) ? 0.0 : std::atan2(yz(),xz()); 82 // or ---- return (yz() == 0.0) ? 0.0 : s 82 // or ---- return (yz() == 0.0) ? 0.0 : std::atan2(yz(),xz()); 83 } 83 } 84 84 85 double HepRotationX::thetaX() const { 85 double HepRotationX::thetaX() const { 86 return safe_acos(zx()); 86 return safe_acos(zx()); 87 // or ---- return CLHEP::halfpi; 87 // or ---- return CLHEP::halfpi; 88 } 88 } 89 89 90 double HepRotationX::thetaY() const { 90 double HepRotationX::thetaY() const { 91 return safe_acos(zy()); 91 return safe_acos(zy()); 92 } 92 } 93 93 94 double HepRotationX::thetaZ() const { 94 double HepRotationX::thetaZ() const { 95 return safe_acos(zz()); 95 return safe_acos(zz()); 96 // or ---- return d; 96 // or ---- return d; 97 } 97 } 98 98 99 void HepRotationX::setDelta ( double ddelta ) 99 void HepRotationX::setDelta ( double ddelta ) { 100 set(ddelta); 100 set(ddelta); 101 } 101 } 102 102 103 void HepRotationX::decompose 103 void HepRotationX::decompose 104 (HepAxisAngle & rotation, Hep3Vector & boost 104 (HepAxisAngle & rotation, Hep3Vector & boost) const { 105 boost.set(0,0,0); 105 boost.set(0,0,0); 106 rotation = axisAngle(); 106 rotation = axisAngle(); 107 } 107 } 108 108 109 void HepRotationX::decompose 109 void HepRotationX::decompose 110 (Hep3Vector & boost, HepAxisAngle & rotation 110 (Hep3Vector & boost, HepAxisAngle & rotation) const { 111 boost.set(0,0,0); 111 boost.set(0,0,0); 112 rotation = axisAngle(); 112 rotation = axisAngle(); 113 } 113 } 114 114 115 void HepRotationX::decompose 115 void HepRotationX::decompose 116 (HepRotation & rotation, HepBoost & bo 116 (HepRotation & rotation, HepBoost & boost) const { 117 boost.set(0,0,0); 117 boost.set(0,0,0); 118 rotation = HepRotation(*this); 118 rotation = HepRotation(*this); 119 } 119 } 120 120 121 void HepRotationX::decompose 121 void HepRotationX::decompose 122 (HepBoost & boost, HepRotation & rotat 122 (HepBoost & boost, HepRotation & rotation) const { 123 boost.set(0,0,0); 123 boost.set(0,0,0); 124 rotation = HepRotation(*this); 124 rotation = HepRotation(*this); 125 } 125 } 126 126 127 double HepRotationX::distance2( const HepRotat 127 double HepRotationX::distance2( const HepRotationX & r ) const { 128 double answer = 2.0 * ( 1.0 - ( its_s * r.it 128 double answer = 2.0 * ( 1.0 - ( its_s * r.its_s + its_c * r.its_c ) ) ; 129 return (answer >= 0) ? answer : 0; 129 return (answer >= 0) ? answer : 0; 130 } 130 } 131 131 132 double HepRotationX::distance2( const HepRotat 132 double HepRotationX::distance2( const HepRotation & r ) const { 133 double sum = r.xx() + 133 double sum = r.xx() + 134 yy() * r.yy() + yz() * r 134 yy() * r.yy() + yz() * r.yz() 135 + zy() * r.zy() + zz() * r 135 + zy() * r.zy() + zz() * r.zz(); 136 double answer = 3.0 - sum; 136 double answer = 3.0 - sum; 137 return (answer >= 0 ) ? answer : 0; 137 return (answer >= 0 ) ? answer : 0; 138 } 138 } 139 139 140 double HepRotationX::distance2( const HepLoren 140 double HepRotationX::distance2( const HepLorentzRotation & lt ) const { 141 HepAxisAngle a; 141 HepAxisAngle a; 142 Hep3Vector b; 142 Hep3Vector b; 143 lt.decompose(b, a); 143 lt.decompose(b, a); 144 double bet = b.beta(); 144 double bet = b.beta(); 145 double bet2 = bet*bet; 145 double bet2 = bet*bet; 146 HepRotation r(a); 146 HepRotation r(a); 147 return bet2/(1-bet2) + distance2(r); 147 return bet2/(1-bet2) + distance2(r); 148 } 148 } 149 149 150 double HepRotationX::distance2( const HepBoost 150 double HepRotationX::distance2( const HepBoost & lt ) const { 151 return distance2( HepLorentzRotation(lt)); 151 return distance2( HepLorentzRotation(lt)); 152 } 152 } 153 153 154 double HepRotationX::howNear( const HepRotatio 154 double HepRotationX::howNear( const HepRotationX & r ) const { 155 return std::sqrt(distance2(r)); 155 return std::sqrt(distance2(r)); 156 } 156 } 157 double HepRotationX::howNear( const HepRotatio 157 double HepRotationX::howNear( const HepRotation & r ) const { 158 return std::sqrt(distance2(r)); 158 return std::sqrt(distance2(r)); 159 } 159 } 160 double HepRotationX::howNear( const HepBoost & 160 double HepRotationX::howNear( const HepBoost & b ) const { 161 return std::sqrt(distance2(b)); 161 return std::sqrt(distance2(b)); 162 } 162 } 163 double HepRotationX::howNear( const HepLorentz 163 double HepRotationX::howNear( const HepLorentzRotation & lt ) const { 164 return std::sqrt(distance2(lt)); 164 return std::sqrt(distance2(lt)); 165 } 165 } 166 bool HepRotationX::isNear(const HepRotationX & 166 bool HepRotationX::isNear(const HepRotationX & r,double epsilon)const{ 167 return (distance2(r) <= epsilon*epsilon); 167 return (distance2(r) <= epsilon*epsilon); 168 } 168 } 169 bool HepRotationX::isNear(const HepRotation & 169 bool HepRotationX::isNear(const HepRotation & r,double epsilon) const{ 170 return (distance2(r) <= epsilon*epsilon); 170 return (distance2(r) <= epsilon*epsilon); 171 } 171 } 172 bool HepRotationX::isNear( const HepBoost & lt 172 bool HepRotationX::isNear( const HepBoost & lt,double epsilon) const { 173 return (distance2(lt) <= epsilon*epsilon); 173 return (distance2(lt) <= epsilon*epsilon); 174 } 174 } 175 175 176 bool HepRotationX::isNear( const HepLorentzRot 176 bool HepRotationX::isNear( const HepLorentzRotation & lt, 177 double ep 177 double epsilon ) const { 178 return (distance2(lt) <= epsilon*epsilon); 178 return (distance2(lt) <= epsilon*epsilon); 179 } 179 } 180 180 181 double HepRotationX::norm2() const { 181 double HepRotationX::norm2() const { 182 return 2.0 - 2.0 * its_c; 182 return 2.0 - 2.0 * its_c; 183 } 183 } 184 184 185 std::ostream & HepRotationX::print( std::ostre 185 std::ostream & HepRotationX::print( std::ostream & os ) const { 186 os << "\nRotation about X (" << its_d << 186 os << "\nRotation about X (" << its_d << 187 ") [cos d = " << its_c << " sin d = " << i 187 ") [cos d = " << its_c << " sin d = " << its_s << "]\n"; 188 return os; 188 return os; 189 } 189 } 190 190 191 } // namespace CLHEP 191 } // namespace CLHEP 192 192 193 193