#include <IterativeHelixExtrapolatorToLine.h>

Inheritance diagram for IterativeHelixExtrapolatorToLine:

Public Member Functions
DirectionType	direction (double s) const override

DirectionTypeDouble	directionInDouble (double s) const

	IterativeHelixExtrapolatorToLine (const PositionType &point, const DirectionType &direction, const float curvature, const PropagationDirection propDir=anyDirection)

std::pair< bool, double >	pathLength (const GlobalPoint &point) const override

std::pair< bool, double >	pathLength (const Line &line) const override

PositionType	position (double s) const override

PositionTypeDouble	positionInDouble (double s) const

	~IterativeHelixExtrapolatorToLine () override

Public Member Functions inherited from HelixLineExtrapolation
virtual	~HelixLineExtrapolation ()=default

Private Member Functions
template<class T >
std::pair< bool, double >	genericPathLength (const T &object) const
	common functionality for extrapolation to line or point More...

Private Attributes
double	theCachedCDPhi

double	theCachedDPhi

double	theCachedS

double	theCachedSDPhi

double	theCosPhi0

double	theCosTheta

const PropagationDirection	thePropDir

HelixExtrapolatorToLine2Order	theQuadraticSolutionFromStart

const double	theRho

double	theSinPhi0

double	theSinTheta

const double	theX0

const double	theY0

const double	theZ0

Additional Inherited Members
Public Types inherited from HelixLineExtrapolation
typedef Basic3DVector< float >	DirectionType

typedef Basic3DVector< double >	DirectionTypeDouble

typedef Basic3DVector< float >	PositionType

typedef Basic3DVector< double >	PositionTypeDouble

Detailed Description

Calculates closest approach of a helix to a line or a point by iterative use of a 2nd order expansion of the helix.

Definition at line 13 of file IterativeHelixExtrapolatorToLine.h.

Constructor & Destructor Documentation

IterativeHelixExtrapolatorToLine::IterativeHelixExtrapolatorToLine	(	const PositionType &	point,
		const DirectionType &	direction,
		const float	curvature,
		const PropagationDirection	propDir = `anyDirection`
	)

Constructor using point, direction and (transverse!) curvature.

Definition at line 4 of file IterativeHelixExtrapolatorToLine.cc.

References AlCaHLTBitMon_ParallelJobs::p, DiDispStaMuonMonitor_cfi::pt, mathSSE::sqrt(), theCosPhi0, theCosTheta, theSinPhi0, theSinTheta, Basic3DVector< T >::x(), Basic3DVector< T >::y(), and Basic3DVector< T >::z().

     : theX0(point.x()),
       theY0(point.y()),
       theZ0(point.z()),
       theRho(curvature),
       theQuadraticSolutionFromStart(point, direction, curvature, propDir),
       thePropDir(propDir),
       theCachedS(0),
       theCachedDPhi(0.),
       theCachedSDPhi(0.),
       theCachedCDPhi(1.) {
   //
   // Components of direction vector (with correct normalisation)
   //
   double px = direction.x();
   double py = direction.y();
   double pz = direction.z();
   double pt = px * px + py * py;
   double p = sqrt(pt + pz * pz);
   pt = sqrt(pt);
   theCosPhi0 = px / pt;
   theSinPhi0 = py / pt;
   theCosTheta = pz / p;
   theSinTheta = pt / p;
 }

IterativeHelixExtrapolatorToLine::~IterativeHelixExtrapolatorToLine ( )

inlineoverride

Definition at line 22 of file IterativeHelixExtrapolatorToLine.h.

22 {}

Member Function Documentation

HelixLineExtrapolation::DirectionType IterativeHelixExtrapolatorToLine::direction ( double s ) const

overridevirtual

Direction at pathlength s from the starting point.

Implements HelixLineExtrapolation.

Definition at line 174 of file IterativeHelixExtrapolatorToLine.cc.

References directionInDouble().

Referenced by AnalyticalImpactPointExtrapolator::propagateWithHelix(), and AnalyticalTrajectoryExtrapolatorToLine::propagateWithHelix().

                                                                                               {
   // use result in double precision
   //   DirectionTypeDouble dir = directionInDouble(s);
   //   return DirectionType(dir.x(),dir.y(),dir.z());
   return DirectionType(directionInDouble(s));
 }

HelixLineExtrapolation::DirectionTypeDouble IterativeHelixExtrapolatorToLine::directionInDouble ( double s ) const

Direction at pathlength s from the starting point.

Definition at line 184 of file IterativeHelixExtrapolatorToLine.cc.

References funct::cos(), HelixExtrapolatorToLine2Order::directionInDouble(), alignCSCRings::e, alignCSCRings::s, funct::sin(), theCachedCDPhi, theCachedDPhi, theCachedS, theCachedSDPhi, theCosPhi0, theCosTheta, theQuadraticSolutionFromStart, theRho, theSinPhi0, and theSinTheta.

Referenced by direction(), and genericPathLength().

                                                                                                             {
   //
   // Calculate delta phi (if not already available)
   //
   if (s != theCachedS) {
     theCachedS = s;
     theCachedDPhi = theCachedS * theRho * theSinTheta;
     theCachedSDPhi = sin(theCachedDPhi);
     theCachedCDPhi = cos(theCachedDPhi);
   }
 
   if (fabs(theCachedDPhi) > 1.e-4) {
     // full helix formula
     return DirectionTypeDouble(theCosPhi0 * theCachedCDPhi - theSinPhi0 * theCachedSDPhi,
                                theSinPhi0 * theCachedCDPhi + theCosPhi0 * theCachedSDPhi,
                                theCosTheta / theSinTheta);
   } else {
     // 1st order
     return theQuadraticSolutionFromStart.directionInDouble(theCachedS);
   }
 }

template<class T >

std::pair< bool, double > IterativeHelixExtrapolatorToLine::genericPathLength ( const T & object ) const

private

common functionality for extrapolation to line or point

Definition at line 53 of file IterativeHelixExtrapolatorToLine.cc.

References alongMomentum, anyDirection, directionInDouble(), first, align_cfg::iteration, oppositeToMomentum, HelixExtrapolatorToLine2Order::pathLength(), positionInDouble(), theCosPhi0, theCosTheta, thePropDir, theQuadraticSolutionFromStart, theRho, theSinPhi0, theSinTheta, theX0, theY0, theZ0, Basic3DVector< T >::x(), Basic3DVector< T >::y(), and Basic3DVector< T >::z().

Referenced by pathLength().

                                                                                                {
   //
   // Constants used for control of convergence
   //
   const int maxIterations(100);
   //
   // Prepare internal value of the propagation direction and position / direction vectors for iteration
   //
   PropagationDirection propDir = thePropDir;
   PositionTypeDouble xnew(theX0, theY0, theZ0);
   DirectionTypeDouble pnew(theCosPhi0, theSinPhi0, theCosTheta / theSinTheta);
   //
   // Prepare iterations: count and total pathlength
   //
   unsigned int iteration(maxIterations + 1);
   double dSTotal(0.);
   //
   // Convergence criterion: maximal lateral displacement in a step < 1um
   //
   double maxDeltaS2 = 2 * 1.e-4 / theSinTheta / theSinTheta / fabs(theRho);
   //
   bool first(true);
   while (true) {
     //
     // return empty solution vector if no convergence after maxIterations iterations
     //
     if (--iteration == 0) {
       return std::pair<bool, double>(false, 0);
     }
     //
     // Use existing second order object at first pass, create temporary object
     // for subsequent passes.
     //
     std::pair<bool, double> deltaS1;
     if (first) {
       first = false;
       deltaS1 = theQuadraticSolutionFromStart.pathLength(object);
     } else {
       HelixExtrapolatorToLine2Order linearCrossing(
           xnew.x(), xnew.y(), xnew.z(), pnew.x(), pnew.y(), theCosTheta, theSinTheta, theRho, anyDirection);
       deltaS1 = linearCrossing.pathLength(object);
     }
     if (!deltaS1.first)
       return deltaS1;
     //
     // Calculate total pathlength
     //
     dSTotal += deltaS1.second;
     PropagationDirection newDir = dSTotal >= 0 ? alongMomentum : oppositeToMomentum;
     if (propDir == anyDirection) {
       propDir = newDir;
     } else {
       if (newDir != propDir)
         return std::pair<bool, double>(false, 0);
     }
     //
     // Check convergence
     //
     if (deltaS1.second * deltaS1.second < maxDeltaS2)
       break;
     //
     // Step forward by dSTotal.
     //
     xnew = positionInDouble(dSTotal);
     pnew = directionInDouble(dSTotal);
   }
   //
   // Return result
   //
   return std::pair<bool, double>(true, dSTotal);
 }

std::pair< bool, double > IterativeHelixExtrapolatorToLine::pathLength ( const GlobalPoint & point ) const

overridevirtual

Propagation status (true if valid) and (signed) path length along the helix from the starting point to the closest approach. to the point. The starting point is given in the constructor.

Implements HelixLineExtrapolation.

Definition at line 37 of file IterativeHelixExtrapolatorToLine.cc.

References genericPathLength().

Referenced by AnalyticalImpactPointExtrapolator::propagateWithHelix(), and AnalyticalTrajectoryExtrapolatorToLine::propagateWithHelix().

                                                                                                  {
   return genericPathLength(point);
 }

std::pair< bool, double > IterativeHelixExtrapolatorToLine::pathLength ( const Line & line ) const

overridevirtual

Propagation status (true if valid) and (signed) path length along the helix from the starting point to the closest approach to the line. The starting point is given in the constructor.

Implements HelixLineExtrapolation.

Definition at line 45 of file IterativeHelixExtrapolatorToLine.cc.

References genericPathLength().

                                                                                          {
   return genericPathLength(line);
 }

HelixLineExtrapolation::PositionType IterativeHelixExtrapolatorToLine::position ( double s ) const

overridevirtual

Position at pathlength s from the starting point.

Implements HelixLineExtrapolation.

Definition at line 128 of file IterativeHelixExtrapolatorToLine.cc.

References positionInDouble().

Referenced by AnalyticalImpactPointExtrapolator::propagateWithHelix(), and AnalyticalTrajectoryExtrapolatorToLine::propagateWithHelix().

                                                                                             {
   // use result in double precision
   return PositionType(positionInDouble(s));
 }

HelixLineExtrapolation::PositionTypeDouble IterativeHelixExtrapolatorToLine::positionInDouble ( double s ) const

Position at pathlength s from the starting point.

Definition at line 136 of file IterativeHelixExtrapolatorToLine.cc.

References funct::cos(), alignCSCRings::e, HelixExtrapolatorToLine2Order::positionInDouble(), alignCSCRings::s, funct::sin(), theCachedCDPhi, theCachedDPhi, theCachedS, theCachedSDPhi, theCosPhi0, theCosTheta, theQuadraticSolutionFromStart, theRho, theSinPhi0, theSinTheta, theX0, theY0, and theZ0.

Referenced by genericPathLength(), and position().

                                                                                                           {
   //
   // Calculate delta phi (if not already available)
   //
   if (s != theCachedS) {
     theCachedS = s;
     theCachedDPhi = theCachedS * theRho * theSinTheta;
     theCachedSDPhi = sin(theCachedDPhi);
     theCachedCDPhi = cos(theCachedDPhi);
   }
   //
   // Calculate with appropriate formulation of full helix formula or with
   //   1st order approximation.
   //
   //    if ( fabs(theCachedDPhi)>1.e-1 ) {
   if (fabs(theCachedDPhi) > 1.e-4) {
     // "standard" helix formula
     return PositionTypeDouble(theX0 + (-theSinPhi0 * (1. - theCachedCDPhi) + theCosPhi0 * theCachedSDPhi) / theRho,
                               theY0 + (theCosPhi0 * (1. - theCachedCDPhi) + theSinPhi0 * theCachedSDPhi) / theRho,
                               theZ0 + theCachedS * theCosTheta);
   }
   //    else if ( fabs(theCachedDPhi)>theNumericalPrecision ) {
   //      // full helix formula, but avoiding (1-cos(deltaPhi)) for small angles
   //      return PositionTypeDouble(theX0+(-theSinPhi0*theCachedSDPhi*theCachedSDPhi/(1.+theCachedCDPhi)+
   //                     theCosPhi0*theCachedSDPhi)/theRho,
   //                  theY0+(theCosPhi0*theCachedSDPhi*theCachedSDPhi/(1.+theCachedCDPhi)+
   //                     theSinPhi0*theCachedSDPhi)/theRho,
   //                  theZ0+theCachedS*theCosTheta);
   //    }
   else {
     // Use 1st order.
     return theQuadraticSolutionFromStart.positionInDouble(theCachedS);
   }
 }