PWG2/FEMTOSCOPY/AliFemto/Infrastructure/AliFmHelix.h

   1 /**
   2  * \class AliFmHelix
   3  * \author Thomas Ullrich, Sep 26 1997
   4  *
   5  * Parametrization of a helix. Can also cope with straight tracks, i.e.
   6  * with zero curvature. This represents only the mathematical model of
   7  * a helix. See the SCL user guide for more.
   8  */
   9 /***************************************************************************
  10  *
  11  * $Id$
  12  *
  13  * Author: Thomas Ullrich, Sep 1997
  14  ***************************************************************************
  15  *
  16  * Description: Parametrization of a helix
  17  *
  18  ***************************************************************************
  19  *
  20  * $Log$
  21  * Revision 1.1.1.1  2007/04/25 15:38:41  panos
  22  * Importing the HBT code dir
  23  *
  24  * Revision 1.1.1.1  2007/03/07 10:14:49  mchojnacki
  25  * First version on CVS
  26  *
  27  * Revision 1.11  2005/10/13 22:23:27  genevb
  28  * NoSolution is public
  29  *
  30  * Revision 1.10  2005/07/06 18:49:56  fisyak
  31  * Replace AliFmHelixD, AliFmLorentzVectorD,AliFmLorentzVectorF,AliFmMatrixD,AliFmMatrixF,AliFmPhysicalHelixD,AliFmThreeVectorD,AliFmThreeVectorF by templated version
  32  *
  33  * Revision 1.9  2004/12/02 02:51:16  ullrich
  34  * Added option to pathLenghth() and distance() to search for
  35  * DCA only within one period. Default stays as it was.
  36  *
  37  * Revision 1.8  2003/10/30 20:06:46  perev
  38  * Check of quality added
  39  *
  40  * Revision 1.7  2002/06/21 17:49:25  genevb
  41  * Some minor speed improvements
  42  *
  43  * Revision 1.6  2002/04/24 02:41:55  ullrich
  44  * Restored old format.
  45  *
  46  **************************************************************************/
  47
  48 #ifndef ST_HELIX_HH
  49 #define ST_HELIX_HH
  50
  51 #include <math.h>
  52 #include <utility>
  53 #include <algorithm>
  54 #include "AliFmThreeVector.h"
  55 #if !defined(ST_NO_NAMESPACES)
  56 using std::pair;
  57 using std::swap;
  58 using std::max;
  59 #endif
  60
  61 #ifdef WIN32
  62 #include "gcc2vs.h"
  63 #endif
  64
  65 class AliFmHelix {
  66 public:
  67     /// curvature, dip angle, phase, origin, h
  68     AliFmHelix(double c, double dip, double phase,
  69                const AliFmThreeVector<double>& o, int h=-1);
  70
  71     virtual ~AliFmHelix();
  72     // AliFmHelix(const AliFmHelix&);                   // use default
  73     // AliFmHelix& operator=(const AliFmHelix&);        // use default
  74
  75     double       dipAngle()   const;
  76     double       curvature()  const;    /// 1/R in xy-plane
  77     double       phase()      const;    /// aziumth in xy-plane measured from ring center
  78     double       xcenter()    const;    /// x-center of circle in xy-plane
  79     double       ycenter()    const;    /// y-center of circle in xy-plane
  80     int          h()          const;    /// -sign(q*B);
  81
  82     const AliFmThreeVector<double>& origin() const;     /// starting point
  83
  84     void setParameters(double c, double dip, double phase, const AliFmThreeVector<double>& o, int h);
  85
  86     double       x(double s)  const;
  87     double       y(double s)  const;
  88     double       z(double s)  const;
  89
  90     AliFmThreeVector<double>  at(double s) const;
  91
  92     /// returns period length of helix
  93     double       period()       const;
  94
  95     /// path length at given r (cylindrical r)
  96     pair<double, double> pathLength(double r)   const;
  97
  98     /// path length at given r (cylindrical r, cylinder axis at x,y)
  99     pair<double, double> pathLength(double r, double x, double y, bool scanPeriods = true);
 100
 101     /// path length at distance of closest approach to a given point
 102     double       pathLength(const AliFmThreeVector<double>& p, bool scanPeriods = true) const;
 103
 104     /// path length at intersection with plane
 105     double       pathLength(const AliFmThreeVector<double>& r,
 106                             const AliFmThreeVector<double>& n) const;
 107
 108     /// path length at distance of closest approach in the xy-plane to a given point
 109     double       pathLength(double x, double y) const;
 110
 111     /// path lengths at dca between two helices
 112     pair<double, double> pathLengths(const AliFmHelix&, bool scanPeriods = true) const;
 113
 114     /// minimal distance between point and helix
 115     double       distance(const AliFmThreeVector<double>& p, bool scanPeriods = true) const;
 116
 117     /// checks for valid parametrization
 118     bool         valid(double world = 1.e+5) const {return !bad(world);}
 119     int            bad(double world = 1.e+5) const;
 120
 121     /// move the origin along the helix to s which becomes then s=0
 122     virtual void moveOrigin(double s);
 123
 124     static const double NoSolution;
 125
 126 protected:
 127     AliFmHelix();
 128
 129     void setCurvature(double);  /// performs also various checks
 130     void setPhase(double);
 131     void setDipAngle(double);
 132
 133     /// value of S where distance in x-y plane is minimal
 134     double fudgePathLength(const AliFmThreeVector<double>&) const;
 135
 136 protected:
 137     bool                   mSingularity;        // true for straight line case (B=0)
 138     AliFmThreeVector<double>  mOrigin;
 139     double                 mDipAngle;
 140     double                 mCurvature;
 141     double                 mPhase;
 142     int                    mH;                  // -sign(q*B);
 143
 144     double                 mCosDipAngle;
 145     double                 mSinDipAngle;
 146     double                 mCosPhase;
 147     double                 mSinPhase;
 148 #ifdef __ROOT__
 149   ClassDef(AliFmHelix,1)
 150 #endif
 151 };
 152
 153 //
 154 //     Non-member functions
 155 //
 156 int operator== (const AliFmHelix&, const AliFmHelix&);
 157 int operator!= (const AliFmHelix&, const AliFmHelix&);
 158 ostream& operator<<(ostream&, const AliFmHelix&);
 159
 160 //
 161 //     Inline functions
 162 //
 163 inline int AliFmHelix::h() const {return mH;}
 164
 165 inline double AliFmHelix::dipAngle() const {return mDipAngle;}
 166
 167 inline double AliFmHelix::curvature() const {return mCurvature;}
 168
 169 inline double AliFmHelix::phase() const {return mPhase;}
 170
 171 inline double AliFmHelix::x(double s) const
 172 {
 173     if (mSingularity)
 174         return mOrigin.x() - s*mCosDipAngle*mSinPhase;
 175     else
 176         return mOrigin.x() + (cos(mPhase + s*mH*mCurvature*mCosDipAngle)-mCosPhase)/mCurvature;
 177 }
 178
 179 inline double AliFmHelix::y(double s) const
 180 {
 181     if (mSingularity)
 182         return mOrigin.y() + s*mCosDipAngle*mCosPhase;
 183     else
 184         return mOrigin.y() + (sin(mPhase + s*mH*mCurvature*mCosDipAngle)-mSinPhase)/mCurvature;
 185 }
 186
 187 inline double AliFmHelix::z(double s) const
 188 {
 189     return mOrigin.z() + s*mSinDipAngle;
 190 }
 191
 192 inline const AliFmThreeVector<double>& AliFmHelix::origin() const {return mOrigin;}
 193
 194 inline AliFmThreeVector<double> AliFmHelix::at(double s) const
 195 {
 196     return AliFmThreeVector<double>(x(s), y(s), z(s));
 197 }
 198
 199 inline double AliFmHelix::pathLength(double x, double y) const
 200 {
 201     return fudgePathLength(AliFmThreeVector<double>(x, y, 0));
 202 }
 203 inline int AliFmHelix::bad(double WorldSize) const
 204 {
 205
 206     int ierr;
 207     if (!::finite(mDipAngle    ))       return   11;
 208     if (!::finite(mCurvature   ))       return   12;
 209
 210     ierr = mOrigin.bad(WorldSize);
 211     if (ierr)                           return    3+ierr*100;
 212
 213     if (::fabs(mDipAngle)  >1.58)       return   21;
 214     double qwe = ::fabs(::fabs(mDipAngle)-M_PI/2);
 215     if (qwe < 1./WorldSize      )       return   31;
 216
 217     if (::fabs(mCurvature) > WorldSize) return   22;
 218     if (mCurvature < 0          )       return   32;
 219
 220     if (abs(mH) != 1            )       return   24;
 221
 222     return 0;
 223 }
 224
 225 #endif