Math/GenVector/CylindricalEta3D.h

0001 // @(#)root/mathcore:$Id$
0002 // Authors: W. Brown, M. Fischler, L. Moneta    2005
0003
0004  /**********************************************************************
0005   *                                                                    *
0006   * Copyright (c) 2005 , LCG ROOT MathLib Team  and                    *
0007   *                      FNAL LCG ROOT MathLib Team                    *
0008   *                                                                    *
0009   *                                                                    *
0010   **********************************************************************/
0011
0012 // Header file for class CylindricalEta3D
0013 //
0014 // Created by: Lorenzo Moneta  at Mon May 30 11:58:46 2005
0015 // Major revamp:  M. Fischler  at Fri Jun 10 2005
0016 //
0017 // Last update: $Id$
0018
0019 //
0020 #ifndef ROOT_Math_GenVector_CylindricalEta3D
0021 #define ROOT_Math_GenVector_CylindricalEta3D  1
0022
0023 #include "Math/Math.h"
0024
0025 #include "Math/GenVector/etaMax.h"
0026
0027 #include "TMath.h"
0028 #include <limits>
0029 #include <cmath>
0030
0031
0032 namespace ROOT {
0033
0034 namespace Math {
0035
0036 //__________________________________________________________________________________________
0037   /**
0038       Class describing a cylindrical coordinate system based on eta (pseudorapidity) instead of z.
0039       The base coordinates are rho (transverse component) , eta and phi
0040       Phi is restricted to be in the range [-PI,PI)
0041
0042       @ingroup GenVector
0043
0044       @see GenVector
0045   */
0046
0047 template <class T>
0048 class CylindricalEta3D {
0049
0050 public :
0051
0052   typedef T Scalar;
0053   static constexpr unsigned int Dimension = 3U;
0054
0055   /**
0056      Default constructor with rho=eta=phi=0
0057    */
0058   constexpr CylindricalEta3D() noexcept = default;
0059
0060   /**
0061      Construct from rho eta and phi values
0062    */
0063   constexpr CylindricalEta3D(Scalar rho, Scalar eta, Scalar phi) noexcept : fRho(rho), fEta(eta), fPhi(phi)
0064   {
0065      Restrict();
0066   }
0067
0068   /**
0069      Construct from any Vector or coordinate system implementing
0070      Rho(), Eta() and Phi()
0071     */
0072   template <class CoordSystem >
0073   explicit CylindricalEta3D( const CoordSystem & v ) :
0074      fRho(v.Rho() ),  fEta(v.Eta() ),  fPhi(v.Phi() )
0075   {
0076      using std::log;
0077      static Scalar bigEta = Scalar(-0.3) * log(std::numeric_limits<Scalar>::epsilon());
0078      if (std::fabs(fEta) > bigEta) {
0079         // This gives a small absolute adjustment in rho,
0080         // which, for large eta, results in a significant
0081         // improvement in the faithfullness of reproducing z.
0082         fRho *= v.Z() / Z();
0083     }
0084   }
0085
0086    /**
0087       Set internal data based on an array of 3 Scalar numbers
0088    */
0089    void SetCoordinates( const Scalar src[] )
0090    { fRho=src[0]; fEta=src[1]; fPhi=src[2]; Restrict(); }
0091
0092    /**
0093       get internal data into an array of 3 Scalar numbers
0094    */
0095    void GetCoordinates( Scalar dest[] ) const
0096    { dest[0] = fRho; dest[1] = fEta; dest[2] = fPhi; }
0097
0098    /**
0099       Set internal data based on 3 Scalar numbers
0100    */
0101    void SetCoordinates(Scalar  rho, Scalar  eta, Scalar  phi)
0102    { fRho=rho; fEta=eta; fPhi=phi; Restrict(); }
0103
0104    /**
0105       get internal data into 3 Scalar numbers
0106    */
0107    void GetCoordinates(Scalar& rho, Scalar& eta, Scalar& phi) const
0108    {rho=fRho; eta=fEta; phi=fPhi;}
0109
0110 private:
0111    inline static Scalar pi() { return TMath::Pi(); }
0112    inline void Restrict() {
0113       using std::floor;
0114       if (fPhi <= -pi() || fPhi > pi()) fPhi = fPhi - floor(fPhi / (2 * pi()) + .5) * 2 * pi();
0115       return;
0116    }
0117 public:
0118
0119    // accessors
0120
0121    T Rho()   const { return fRho; }
0122    T Eta()   const { return fEta; }
0123    T Phi()   const { return fPhi; }
0124    T X() const { using std::cos; return fRho * cos(fPhi); }
0125    T Y() const { using std::sin; return fRho * sin(fPhi); }
0126    T Z() const
0127    {
0128       using std::sinh;
0129       return fRho > 0 ? fRho * sinh(fEta) : fEta == 0 ? 0 : fEta > 0 ? fEta - etaMax<T>() : fEta + etaMax<T>();
0130    }
0131    T R() const
0132    {
0133       using std::cosh;
0134       return fRho > 0 ? fRho * cosh(fEta)
0135                       : fEta > etaMax<T>() ? fEta - etaMax<T>() : fEta < -etaMax<T>() ? -fEta - etaMax<T>() : 0;
0136    }
0137    T Mag2() const
0138    {
0139       const Scalar r = R();
0140       return r * r;
0141    }
0142    T Perp2() const { return fRho*fRho;            }
0143    T Theta() const { using std::atan; return fRho > 0 ? 2 * atan(exp(-fEta)) : (fEta >= 0 ? 0 : pi()); }
0144
0145    // setters (only for data members)
0146
0147
0148    /**
0149        set the rho coordinate value keeping eta and phi constant
0150    */
0151    void SetRho(T rho) {
0152       fRho = rho;
0153    }
0154
0155    /**
0156        set the eta coordinate value keeping rho and phi constant
0157    */
0158    void SetEta(T eta) {
0159       fEta = eta;
0160    }
0161
0162    /**
0163        set the phi coordinate value keeping rho and eta constant
0164    */
0165    void SetPhi(T phi) {
0166       fPhi = phi;
0167       Restrict();
0168    }
0169
0170    /**
0171        set all values using cartesian coordinates
0172    */
0173    void SetXYZ(Scalar x, Scalar y, Scalar z);
0174
0175
0176    /**
0177       scale by a scalar quantity a --
0178       for cylindrical eta coords, as long as a >= 0, only rho changes!
0179    */
0180    void Scale (T a) {
0181       if (a < 0) {
0182          Negate();
0183          a = -a;
0184       }
0185       // angles do not change when scaling by a positive quantity
0186       if (fRho > 0) {
0187          fRho *= a;
0188       } else if ( fEta >  etaMax<T>() ) {
0189          fEta =  ( fEta-etaMax<T>())*a + etaMax<T>();
0190       } else if ( fEta < -etaMax<T>() ) {
0191          fEta =  ( fEta+etaMax<T>())*a - etaMax<T>();
0192       } // when rho==0 and eta is not above etaMax, vector represents 0
0193       // and remains unchanged
0194    }
0195
0196    /**
0197       negate the vector
0198    */
0199    void Negate ( ) {
0200       fPhi = ( fPhi > 0 ? fPhi - pi() : fPhi + pi()  );
0201       fEta = -fEta;
0202    }
0203
0204    // assignment operators
0205    /**
0206       generic assignment operator from any coordinate system
0207    */
0208    template <class CoordSystem >
0209    CylindricalEta3D & operator= ( const CoordSystem & c ) {
0210       fRho = c.Rho();
0211       fEta = c.Eta();
0212       fPhi = c.Phi();
0213       return *this;
0214    }
0215
0216    /**
0217       Exact component-by-component equality
0218       Note: Peculiar representations of the zero vector such as (0,1,0) will
0219       not test as equal to one another.
0220    */
0221    bool operator==(const CylindricalEta3D & rhs) const {
0222       return fRho == rhs.fRho && fEta == rhs.fEta && fPhi == rhs.fPhi;
0223    }
0224    bool operator!= (const CylindricalEta3D & rhs) const
0225    {return !(operator==(rhs));}
0226
0227
0228    // ============= Compatibility section ==================
0229
0230    // The following make this coordinate system look enough like a CLHEP
0231    // vector that an assignment member template can work with either
0232    T x() const { return X();}
0233    T y() const { return Y();}
0234    T z() const { return Z(); }
0235
0236    // ============= Specializations for improved speed ==================
0237
0238    // (none)
0239
0240 #if defined(__MAKECINT__) || defined(G__DICTIONARY)
0241
0242    // ====== Set member functions for coordinates in other systems =======
0243
0244    void SetX(Scalar x);
0245
0246    void SetY(Scalar y);
0247
0248    void SetZ(Scalar z);
0249
0250    void SetR(Scalar r);
0251
0252    void SetTheta(Scalar theta);
0253
0254
0255 #endif
0256
0257   private:
0258      T fRho = 0;
0259      T fEta = 0;
0260      T fPhi = 0;
0261 };
0262
0263   } // end namespace Math
0264
0265 } // end namespace ROOT
0266
0267
0268 // move implementations here to avoid circle dependencies
0269
0270 #include "Math/GenVector/Cartesian3D.h"
0271
0272 #if defined(__MAKECINT__) || defined(G__DICTIONARY)
0273 #include "Math/GenVector/GenVector_exception.h"
0274 #include "Math/GenVector/Polar3D.h"
0275 #endif
0276
0277 namespace ROOT {
0278
0279   namespace Math {
0280
0281 template <class T>
0282 void CylindricalEta3D<T>::SetXYZ(Scalar xx, Scalar yy, Scalar zz) {
0283    *this = Cartesian3D<Scalar>(xx, yy, zz);
0284 }
0285
0286 #if defined(__MAKECINT__) || defined(G__DICTIONARY)
0287
0288
0289      // ====== Set member functions for coordinates in other systems =======
0290
0291
0292 template <class T>
0293 void CylindricalEta3D<T>::SetX(Scalar xx) {
0294    GenVector_exception e("CylindricalEta3D::SetX() is not supposed to be called");
0295    throw e;
0296    Cartesian3D<Scalar> v(*this); v.SetX(xx);
0297    *this = CylindricalEta3D<Scalar>(v);
0298 }
0299 template <class T>
0300 void CylindricalEta3D<T>::SetY(Scalar yy) {
0301    GenVector_exception e("CylindricalEta3D::SetY() is not supposed to be called");
0302    throw e;
0303    Cartesian3D<Scalar> v(*this); v.SetY(yy);
0304    *this = CylindricalEta3D<Scalar>(v);
0305 }
0306 template <class T>
0307 void CylindricalEta3D<T>::SetZ(Scalar zz) {
0308    GenVector_exception e("CylindricalEta3D::SetZ() is not supposed to be called");
0309    throw e;
0310    Cartesian3D<Scalar> v(*this); v.SetZ(zz);
0311    *this = CylindricalEta3D<Scalar>(v);
0312 }
0313 template <class T>
0314 void CylindricalEta3D<T>::SetR(Scalar r) {
0315    GenVector_exception e("CylindricalEta3D::SetR() is not supposed to be called");
0316    throw e;
0317    Polar3D<Scalar> v(*this); v.SetR(r);
0318    *this = CylindricalEta3D<Scalar>(v);
0319 }
0320 template <class T>
0321 void CylindricalEta3D<T>::SetTheta(Scalar theta) {
0322    GenVector_exception e("CylindricalEta3D::SetTheta() is not supposed to be called");
0323    throw e;
0324    Polar3D<Scalar> v(*this); v.SetTheta(theta);
0325    *this = CylindricalEta3D<Scalar>(v);
0326 }
0327
0328 #endif
0329
0330
0331   } // end namespace Math
0332
0333 } // end namespace ROOT
0334
0335
0336
0337 #endif /* ROOT_Math_GenVector_CylindricalEta3D  */