interface/private/Basic3DVectorLD.h

0001 #ifndef GeometryVector_Basic3DVectorLD_h
0002 #define GeometryVector_Basic3DVectorLD_h
0003
0004 #ifdef __clang__
0005 #pragma clang diagnostic push
0006 #pragma clang diagnostic ignored "-Wunused-private-field"
0007 #endif
0008
0009 #include "extBasic3DVector.h"
0010 // long double specialization
0011 template <>
0012 class Basic3DVector<long double> {
0013 public:
0014   typedef long double T;
0015   typedef T ScalarType;
0016   typedef Geom::Cylindrical2Cartesian<T> Cylindrical;
0017   typedef Geom::Spherical2Cartesian<T> Spherical;
0018   typedef Spherical Polar;  // synonym
0019
0020   typedef Basic3DVector<T> MathVector;
0021
0022   /** default constructor uses default constructor of T to initialize the
0023    *  components. For built-in floating-point types this means initialization
0024    * to zero??? (force init to 0)
0025    */
0026   Basic3DVector() : theX(0), theY(0), theZ(0), theW(0) {}
0027
0028   /// Copy constructor from same type. Should not be needed but for gcc bug 12685
0029   Basic3DVector(const Basic3DVector& p) : theX(p.x()), theY(p.y()), theZ(p.z()), theW(0) {}
0030
0031   /// Copy constructor and implicit conversion from Basic3DVector of different precision
0032   template <class U>
0033   Basic3DVector(const Basic3DVector<U>& p) : theX(p.x()), theY(p.y()), theZ(p.z()), theW(0) {}
0034
0035   /// constructor from 2D vector (X and Y from 2D vector, z set to zero)
0036   Basic3DVector(const Basic2DVector<T>& p) : theX(p.x()), theY(p.y()), theZ(0), theW(0) {}
0037
0038   /// Assignment operator
0039   Basic3DVector& operator=(const Basic3DVector&) = default;
0040
0041   /** Explicit constructor from other (possibly unrelated) vector classes
0042    *  The only constraint on the argument type is that it has methods
0043    *  x(), y() and z(), and that these methods return a type convertible to T.
0044    *  Examples of use are
0045    *   <BR> construction from a Basic3DVector with different precision
0046    *   <BR> construction from a Hep3Vector
0047    *   <BR> construction from a coordinate system converter
0048    */
0049   template <class OtherPoint>
0050   explicit Basic3DVector(const OtherPoint& p) : theX(p.x()), theY(p.y()), theZ(p.z()), theW(0) {}
0051
0052 #ifdef USE_SSEVECT
0053   // constructor from Vec4
0054   template <typename U>
0055   Basic3DVector(mathSSE::Vec4<U> const& iv) : theX(iv.arr[0]), theY(iv.arr[1]), theZ(iv.arr[2]), theW(0) {}
0056 #endif
0057 #ifdef USE_EXTVECT
0058   // constructor from Vec4
0059   template <typename U>
0060   Basic3DVector(Vec4<U> const& iv) : theX(iv[0]), theY(iv[1]), theZ(iv[2]), theW(0) {}
0061 #endif
0062
0063   /// construct from cartesian coordinates
0064   Basic3DVector(const T& x, const T& y, const T& z) : theX(x), theY(y), theZ(z), theW(0) {}
0065
0066   /** Deprecated construct from polar coordinates, use
0067    *  <BR> Basic3DVector<T>( Basic3DVector<T>::Polar( theta, phi, r))
0068    *  instead.
0069    */
0070   template <typename U>
0071   Basic3DVector(const Geom::Theta<U>& theta, const Geom::Phi<U>& phi, const T& r) {
0072     Polar p(theta.value(), phi.value(), r);
0073     theX = p.x();
0074     theY = p.y();
0075     theZ = p.z();
0076   }
0077
0078   /// Cartesian x coordinate
0079   T x() const { return theX; }
0080
0081   /// Cartesian y coordinate
0082   T y() const { return theY; }
0083
0084   /// Cartesian z coordinate
0085   T z() const { return theZ; }
0086
0087   Basic2DVector<T> xy() const { return Basic2DVector<T>(theX, theY); }
0088
0089   // equality
0090   bool operator==(const Basic3DVector& rh) const { return x() == rh.x() && y() == rh.y() && z() == rh.z(); }
0091
0092   /// The vector magnitude squared. Equivalent to vec.dot(vec)
0093   T mag2() const { return x() * x() + y() * y() + z() * z(); }
0094
0095   /// The vector magnitude. Equivalent to sqrt(vec.mag2())
0096   T mag() const { return std::sqrt(mag2()); }
0097
0098   /// Squared magnitude of transverse component
0099   T perp2() const { return x() * x() + y() * y(); }
0100
0101   /// Magnitude of transverse component
0102   T perp() const { return std::sqrt(perp2()); }
0103
0104   /// Another name for perp()
0105   T transverse() const { return perp(); }
0106
0107   /** Azimuthal angle. The value is returned in radians, in the range (-pi,pi].
0108    *  Same precision as the system atan2(x,y) function.
0109    *  The return type is Geom::Phi<T>, see it's documentation.
0110    */
0111   T barePhi() const { return std::atan2(y(), x()); }
0112   Geom::Phi<T> phi() const { return Geom::Phi<T>(barePhi()); }
0113
0114   /** Polar angle. The value is returned in radians, in the range [0,pi]
0115    *  Same precision as the system atan2(x,y) function.
0116    *  The return type is Geom::Phi<T>, see it's documentation.
0117    */
0118   T bareTheta() const { return std::atan2(perp(), z()); }
0119   Geom::Theta<T> theta() const { return Geom::Theta<T>(std::atan2(perp(), z())); }
0120
0121   /** Pseudorapidity.
0122    *  Does not check for zero transverse component; in this case the behavior
0123    *  is as for divide-by zero, i.e. system-dependent.
0124    */
0125   // T eta() const { return -log( tan( theta()/2.));}
0126   T eta() const { return detailsBasic3DVector::eta(x(), y(), z()); }  // correct
0127
0128   /** Unit vector parallel to this.
0129    *  If mag() is zero, a zero vector is returned.
0130    */
0131   Basic3DVector unit() const {
0132     T my_mag = mag2();
0133     if (my_mag == 0)
0134       return *this;
0135     my_mag = T(1) / std::sqrt(my_mag);
0136     Basic3DVector ret(*this);
0137     return ret *= my_mag;
0138   }
0139
0140   /** Operator += with a Basic3DVector of possibly different precision.
0141    */
0142   template <class U>
0143   Basic3DVector& operator+=(const Basic3DVector<U>& p) {
0144     theX += p.x();
0145     theY += p.y();
0146     theZ += p.z();
0147     return *this;
0148   }
0149
0150   /** Operator -= with a Basic3DVector of possibly different precision.
0151    */
0152   template <class U>
0153   Basic3DVector& operator-=(const Basic3DVector<U>& p) {
0154     theX -= p.x();
0155     theY -= p.y();
0156     theZ -= p.z();
0157     return *this;
0158   }
0159
0160   /// Unary minus, returns a vector with components (-x(),-y(),-z())
0161   Basic3DVector operator-() const { return Basic3DVector(-x(), -y(), -z()); }
0162
0163   /// Scaling by a scalar value (multiplication)
0164   Basic3DVector& operator*=(T t) {
0165     theX *= t;
0166     theY *= t;
0167     theZ *= t;
0168     return *this;
0169   }
0170
0171   /// Scaling by a scalar value (division)
0172   Basic3DVector& operator/=(T t) {
0173     t = T(1) / t;
0174     theX *= t;
0175     theY *= t;
0176     theZ *= t;
0177     return *this;
0178   }
0179
0180   /// Scalar product, or "dot" product, with a vector of same type.
0181   T dot(const Basic3DVector& v) const { return x() * v.x() + y() * v.y() + z() * v.z(); }
0182
0183   /** Scalar (or dot) product with a vector of different precision.
0184    *  The product is computed without loss of precision. The type
0185    *  of the returned scalar is the more precise of the scalar types
0186    *  of the two vectors.
0187    */
0188   template <class U>
0189   typename PreciseFloatType<T, U>::Type dot(const Basic3DVector<U>& v) const {
0190     return x() * v.x() + y() * v.y() + z() * v.z();
0191   }
0192
0193   /// Vector product, or "cross" product, with a vector of same type.
0194   Basic3DVector cross(const Basic3DVector& v) const {
0195     return Basic3DVector(y() * v.z() - v.y() * z(), z() * v.x() - v.z() * x(), x() * v.y() - v.x() * y());
0196   }
0197
0198   /** Vector (or cross) product with a vector of different precision.
0199    *  The product is computed without loss of precision. The type
0200    *  of the returned vector is the more precise of the types
0201    *  of the two vectors.
0202    */
0203   template <class U>
0204   Basic3DVector<typename PreciseFloatType<T, U>::Type> cross(const Basic3DVector<U>& v) const {
0205     return Basic3DVector<typename PreciseFloatType<T, U>::Type>(
0206         y() * v.z() - v.y() * z(), z() * v.x() - v.z() * x(), x() * v.y() - v.x() * y());
0207   }
0208
0209 private:
0210   T theX;
0211   T theY;
0212   T theZ;
0213   T theW;
0214 }
0215 #ifndef __CINT__
0216 __attribute__((aligned(16)))
0217 #endif
0218 ;
0219
0220 /// vector sum and subtraction of vectors of possibly different precision
0221 inline Basic3DVector<long double> operator+(const Basic3DVector<long double>& a, const Basic3DVector<long double>& b) {
0222   typedef Basic3DVector<long double> RT;
0223   return RT(a.x() + b.x(), a.y() + b.y(), a.z() + b.z());
0224 }
0225 inline Basic3DVector<long double> operator-(const Basic3DVector<long double>& a, const Basic3DVector<long double>& b) {
0226   typedef Basic3DVector<long double> RT;
0227   return RT(a.x() - b.x(), a.y() - b.y(), a.z() - b.z());
0228 }
0229
0230 template <class U>
0231 inline Basic3DVector<typename PreciseFloatType<long double, U>::Type> operator+(const Basic3DVector<long double>& a,
0232                                                                                 const Basic3DVector<U>& b) {
0233   typedef Basic3DVector<typename PreciseFloatType<long double, U>::Type> RT;
0234   return RT(a.x() + b.x(), a.y() + b.y(), a.z() + b.z());
0235 }
0236
0237 template <class U>
0238 inline Basic3DVector<typename PreciseFloatType<long double, U>::Type> operator+(const Basic3DVector<U>& a,
0239                                                                                 const Basic3DVector<long double>& b) {
0240   typedef Basic3DVector<typename PreciseFloatType<long double, U>::Type> RT;
0241   return RT(a.x() + b.x(), a.y() + b.y(), a.z() + b.z());
0242 }
0243
0244 template <class U>
0245 inline Basic3DVector<typename PreciseFloatType<long double, U>::Type> operator-(const Basic3DVector<long double>& a,
0246                                                                                 const Basic3DVector<U>& b) {
0247   typedef Basic3DVector<typename PreciseFloatType<long double, U>::Type> RT;
0248   return RT(a.x() - b.x(), a.y() - b.y(), a.z() - b.z());
0249 }
0250
0251 template <class U>
0252 inline Basic3DVector<typename PreciseFloatType<long double, U>::Type> operator-(const Basic3DVector<U>& a,
0253                                                                                 const Basic3DVector<long double>& b) {
0254   typedef Basic3DVector<typename PreciseFloatType<long double, U>::Type> RT;
0255   return RT(a.x() - b.x(), a.y() - b.y(), a.z() - b.z());
0256 }
0257
0258 /// scalar product of vectors of same precision
0259 // template <>
0260 inline long double operator*(const Basic3DVector<long double>& v1, const Basic3DVector<long double>& v2) {
0261   return v1.dot(v2);
0262 }
0263
0264 /// scalar product of vectors of different precision
0265 template <class U>
0266 inline typename PreciseFloatType<long double, U>::Type operator*(const Basic3DVector<long double>& v1,
0267                                                                  const Basic3DVector<U>& v2) {
0268   return v1.x() * v2.x() + v1.y() * v2.y() + v1.z() * v2.z();
0269 }
0270
0271 template <class U>
0272 inline typename PreciseFloatType<long double, U>::Type operator*(const Basic3DVector<U>& v1,
0273                                                                  const Basic3DVector<long double>& v2) {
0274   return v1.x() * v2.x() + v1.y() * v2.y() + v1.z() * v2.z();
0275 }
0276
0277 /** Multiplication by scalar, does not change the precision of the vector.
0278  *  The return type is the same as the type of the vector argument.
0279  */
0280 //template <>
0281 inline Basic3DVector<long double> operator*(const Basic3DVector<long double>& v, long double t) {
0282   return Basic3DVector<long double>(v.x() * t, v.y() * t, v.z() * t);
0283 }
0284
0285 /// Same as operator*( Vector, Scalar)
0286 // template <>
0287 inline Basic3DVector<long double> operator*(long double t, const Basic3DVector<long double>& v) {
0288   return Basic3DVector<long double>(v.x() * t, v.y() * t, v.z() * t);
0289 }
0290
0291 template <typename S>
0292 inline Basic3DVector<long double> operator*(S t, const Basic3DVector<long double>& v) {
0293   return static_cast<long double>(t) * v;
0294 }
0295
0296 template <typename S>
0297 inline Basic3DVector<long double> operator*(const Basic3DVector<long double>& v, S t) {
0298   return static_cast<long double>(t) * v;
0299 }
0300
0301 /** Division by scalar, does not change the precision of the vector.
0302  *  The return type is the same as the type of the vector argument.
0303  */
0304 template <typename S>
0305 inline Basic3DVector<long double> operator/(const Basic3DVector<long double>& v, S s) {
0306   long double t = 1 / s;
0307   return v * t;
0308 }
0309
0310 typedef Basic3DVector<long double> Basic3DVectorLD;
0311
0312 #ifdef __clang__
0313 #pragma clang diagnostic pop
0314 #endif
0315
0316 #endif  // GeometryVector_Basic3DVectorLD_h