00001 #ifndef GeometryVector_oldBasic3DVector_h
00002 #define GeometryVector_oldBasic3DVector_h
00003 #if ( defined(IN_DICTBUILD) || defined(__CINT__) ) && !defined(__REFLEX__)
00004 #define __REFLEX__
00005 #endif
00006 #include "DataFormats/GeometryVector/interface/Basic2DVector.h"
00007 #include "DataFormats/GeometryVector/interface/Theta.h"
00008 #include "DataFormats/GeometryVector/interface/Phi.h"
00009 #include "DataFormats/GeometryVector/interface/PreciseFloatType.h"
00010 #include "DataFormats/GeometryVector/interface/CoordinateSets.h"
00011 #ifndef __REFLEX__
00012 #include "DataFormats/Math/interface/SIMDVec.h"
00013 #endif
00014 #include <iosfwd>
00015 #include <cmath>
00016
00017
00018 namespace detailsBasic3DVector {
00019 inline float __attribute__((always_inline)) __attribute__ ((pure))
00020 eta(float x, float y, float z) { float t(z/std::sqrt(x*x+y*y)); return ::asinhf(t);}
00021 inline double __attribute__((always_inline)) __attribute__ ((pure))
00022 eta(double x, double y, double z) { double t(z/std::sqrt(x*x+y*y)); return ::asinh(t);}
00023 inline long double __attribute__((always_inline)) __attribute__ ((pure))
00024 eta(long double x, long double y, long double z) { long double t(z/std::sqrt(x*x+y*y)); return ::asinhl(t);}
00025 }
00026
00027
00028 template < typename T>
00029 class Basic3DVector {
00030 public:
00031 typedef Basic3DVector<T> MathVector;
00032
00033
00034 typedef T ScalarType;
00035 typedef Geom::Cylindrical2Cartesian<T> Cylindrical;
00036 typedef Geom::Spherical2Cartesian<T> Spherical;
00037 typedef Spherical Polar;
00038
00043 Basic3DVector() : theX(0), theY(0), theZ(0), theW(0) {}
00044
00046 Basic3DVector( const Basic3DVector & p) :
00047 theX(p.x()), theY(p.y()), theZ(p.z()), theW(p.w()) {}
00048
00050 template <class U>
00051 Basic3DVector( const Basic3DVector<U> & p) :
00052 theX(p.x()), theY(p.y()), theZ(p.z()), theW(p.w()) {}
00053
00055 Basic3DVector( const Basic2DVector<T> & p) :
00056 theX(p.x()), theY(p.y()), theZ(0), theW(0) {}
00057
00066 template <class OtherPoint>
00067 explicit Basic3DVector( const OtherPoint& p) :
00068 theX(p.x()), theY(p.y()), theZ(p.z()), theW(0) {}
00069
00070
00071 #if defined(USE_EXTVECT)
00072 template<typename U>
00073 Basic3DVector(Vec4<U> const& iv) :
00074 theX(iv[0]), theY(iv[1]), theZ(iv[2]), theW(0) {}
00075 #elif defined(USE_SSEVECT)
00076
00077 template<typename U>
00078 Basic3DVector(mathSSE::Vec4<U> const& iv) :
00079 theX(iv.arr[0]), theY(iv.arr[1]), theZ(iv.arr[2]), theW(0) {}
00080 #endif
00081
00082 #ifndef __REFLEX__
00083
00084 Basic3DVector( const T& x, const T& y, const T& z, const T& w=0) :
00085 theX(x), theY(y), theZ(z), theW(w) {}
00086 #else
00087
00088 Basic3DVector( const T& x, const T& y, const T& z) :
00089 theX(x), theY(y), theZ(z), theW(0) {}
00090 Basic3DVector( const T& x, const T& y, const T& z, const T& w) :
00091 theX(x), theY(y), theZ(z), theW(w) {}
00092 #endif
00093
00094
00095
00100 template <typename U>
00101 Basic3DVector( const Geom::Theta<U>& theta,
00102 const Geom::Phi<U>& phi, const T& r) {
00103 Polar p( theta.value(), phi.value(), r);
00104 theX = p.x(); theY = p.y(); theZ = p.z();
00105 }
00106
00107
00108 T operator[](int i) const { return *((&theX)+i) ;}
00109 T & operator[](int i) { return *((&theX)+i);}
00110
00111
00112
00114 T x() const { return theX;}
00115
00117 T y() const { return theY;}
00118
00120 T z() const { return theZ;}
00121
00122 T w() const { return theW;}
00123
00124
00125 Basic2DVector<T> xy() const { return Basic2DVector<T>(theX,theY);}
00126
00127
00128
00129 bool operator==(const Basic3DVector& rh) const {
00130 return x()==rh.x() && y()==rh.y() && z()==rh.z();
00131 }
00132
00134 T mag2() const { return x()*x() + y()*y()+z()*z();}
00135
00137 T mag() const { return std::sqrt( mag2());}
00138
00140 T perp2() const { return x()*x() + y()*y();}
00141
00143 T perp() const { return std::sqrt( perp2());}
00144
00146 T transverse() const { return perp();}
00147
00152 T barePhi() const {return std::atan2(y(),x());}
00153 Geom::Phi<T> phi() const {return Geom::Phi<T>(barePhi());}
00154
00159 T bareTheta() const {return std::atan2(perp(),z());}
00160 Geom::Theta<T> theta() const {return Geom::Theta<T>(std::atan2(perp(),z()));}
00161
00166
00167 T eta() const { return detailsBasic3DVector::eta(x(),y(),z());}
00168
00172 Basic3DVector unit() const {
00173 T my_mag = mag2();
00174 if (my_mag==0) return *this;
00175 my_mag = T(1)/std::sqrt(my_mag);
00176 return *this * my_mag;
00177 }
00178
00181 template <class U>
00182 Basic3DVector& operator+= ( const Basic3DVector<U>& p) {
00183 theX += p.x();
00184 theY += p.y();
00185 theZ += p.z();
00186 theW += p.w();
00187 return *this;
00188 }
00189
00192 template <class U>
00193 Basic3DVector& operator-= ( const Basic3DVector<U>& p) {
00194 theX -= p.x();
00195 theY -= p.y();
00196 theZ -= p.z();
00197 theW -= p.w();
00198 return *this;
00199 }
00200
00202 Basic3DVector operator-() const { return Basic3DVector(-x(),-y(),-z());}
00203
00205 Basic3DVector& operator*= ( T t) {
00206 theX *= t;
00207 theY *= t;
00208 theZ *= t;
00209 theW *= t;;
00210 return *this;
00211 }
00212
00214 Basic3DVector& operator/= ( T t) {
00215 t = T(1)/t;
00216 theX *= t;
00217 theY *= t;
00218 theZ *= t;
00219 theW *= t;;
00220 return *this;
00221 }
00222
00224 T dot( const Basic3DVector& v) const {
00225 return x()*v.x() + y()*v.y() + z()*v.z();
00226 }
00227
00233 template <class U>
00234 typename PreciseFloatType<T,U>::Type dot( const Basic3DVector<U>& v) const {
00235 return x()*v.x() + y()*v.y() + z()*v.z();
00236 }
00237
00239 Basic3DVector cross( const Basic3DVector& v) const {
00240 return Basic3DVector( y()*v.z() - v.y()*z(),
00241 z()*v.x() - v.z()*x(),
00242 x()*v.y() - v.x()*y());
00243 }
00244
00245
00251 template <class U>
00252 Basic3DVector<typename PreciseFloatType<T,U>::Type>
00253 cross( const Basic3DVector<U>& v) const {
00254 return Basic3DVector<typename PreciseFloatType<T,U>::Type>( y()*v.z() - v.y()*z(),
00255 z()*v.x() - v.z()*x(),
00256 x()*v.y() - v.x()*y());
00257 }
00258
00259 private:
00260 T theX;
00261 T theY;
00262 T theZ;
00263 T theW;
00264 }
00265 #ifndef __CINT__
00266 __attribute__ ((aligned (16)))
00267 #endif
00268 ;
00269
00270
00271 namespace geometryDetails {
00272 std::ostream & print3D(std::ostream& s, double x, double y, double z);
00273 }
00274
00276 template <class T>
00277 inline std::ostream & operator<<( std::ostream& s, const Basic3DVector<T>& v) {
00278 return geometryDetails::print3D(s, v.x(),v.y(), v.z());
00279 }
00280
00281
00283 template <class T, class U>
00284 inline Basic3DVector<typename PreciseFloatType<T,U>::Type>
00285 operator+( const Basic3DVector<T>& a, const Basic3DVector<U>& b) {
00286 typedef Basic3DVector<typename PreciseFloatType<T,U>::Type> RT;
00287 return RT(a.x()+b.x(), a.y()+b.y(), a.z()+b.z(), a.w()+b.w());
00288 }
00289
00290 template <class T, class U>
00291 inline Basic3DVector<typename PreciseFloatType<T,U>::Type>
00292 operator-( const Basic3DVector<T>& a, const Basic3DVector<U>& b) {
00293 typedef Basic3DVector<typename PreciseFloatType<T,U>::Type> RT;
00294 return RT(a.x()-b.x(), a.y()-b.y(), a.z()-b.z(), a.w()-b.w());
00295 }
00296
00298 template <class T>
00299 inline T operator*( const Basic3DVector<T>& v1, const Basic3DVector<T>& v2) {
00300 return v1.dot(v2);
00301 }
00302
00304 template <class T, class U>
00305 inline typename PreciseFloatType<T,U>::Type operator*( const Basic3DVector<T>& v1,
00306 const Basic3DVector<U>& v2) {
00307 return v1.x()*v2.x() + v1.y()*v2.y() + v1.z()*v2.z();
00308 }
00309
00313 template <class T>
00314 inline Basic3DVector<T> operator*( const Basic3DVector<T>& v, T t) {
00315 return Basic3DVector<T>(v.x()*t, v.y()*t, v.z()*t, v.w()*t);
00316 }
00317
00319 template <class T>
00320 inline Basic3DVector<T> operator*(T t, const Basic3DVector<T>& v) {
00321 return Basic3DVector<T>(v.x()*t, v.y()*t, v.z()*t, v.w()*t);
00322 }
00323
00324 template <class T, typename S>
00325 inline Basic3DVector<T> operator*(S t, const Basic3DVector<T>& v) {
00326 return static_cast<T>(t)*v;
00327 }
00328
00329 template <class T, typename S>
00330 inline Basic3DVector<T> operator*(const Basic3DVector<T>& v, S t) {
00331 return static_cast<T>(t)*v;
00332 }
00333
00334
00338 template <class T, typename S>
00339 inline Basic3DVector<T> operator/( const Basic3DVector<T>& v, S s) {
00340 T t = T(1)/s;
00341 return v*t;
00342 }
00343
00344
00345 typedef Basic3DVector<float> Basic3DVectorF;
00346 typedef Basic3DVector<double> Basic3DVectorD;
00347 typedef Basic3DVector<long double> Basic3DVectorLD;
00348
00349
00350 #endif // GeometryVector_Basic3DVector_h
00351
00352