text: #ifndef _COMMONRECO_ANALYTICALPROPAGATOR_H_
  #define _COMMONRECO_ANALYTICALPROPAGATOR_H_
  
  #include "TrackingTools/GeomPropagators/interface/Propagator.h"
  #include "TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h"
  #include "FWCore/Utilities/interface/Visibility.h"
  
  #include <cmath>
  #include <cfloat>
  
  class Surface;
  class Cylinder;
  class Plane;
  class HelixPlaneCrossing;
  class MagneticField; 
  
  
  /** (Mostly) analytical helix propagation to cylindrical or planar surfaces.
   *  Based on GtfGeometricalPropagator with successive replacement of components
   *  (currently: propagation to arbitrary plane).
   */
  
  class AnalyticalPropagator final : public Propagator {
  
  public:
    
    AnalyticalPropagator( const MagneticField* field,
  		        PropagationDirection dir = alongMomentum,
  			float maxDPhi = 1.6,bool isOld=true) :
      Propagator(dir),
      theMaxDPhi2(maxDPhi*maxDPhi),
      theMaxDBzRatio(0.5),
      theField(field),
      isOldPropagationType(isOld) {}
  
    ~AnalyticalPropagator() override {}
  
    //
    // use base class methods where necessary:
    // - propagation from TrajectoryStateOnSurface 
    //     (will use propagation from FreeTrajectoryState)
    // - propagation to general Surface
    //     (will use specialised methods for planes or cylinders)
    //
    using Propagator::propagate;
    using Propagator::propagateWithPath;
    
   private:
    /// propagation to plane with path length  
    std::pair<TrajectoryStateOnSurface,double> 
      propagateWithPath(const FreeTrajectoryState& fts, 
  		    const Plane& plane) const override; 
    
  
    /// propagation to cylinder with path length
    std::pair<TrajectoryStateOnSurface,double> 
    propagateWithPath(const FreeTrajectoryState& fts, 
  		    const Cylinder& cylinder) const override;
  
  
   public:
    /** limitation of change in transverse direction
     *  (to avoid loops).
     */
    bool setMaxDirectionChange( float phiMax) override { 
      theMaxDPhi2 = phiMax*phiMax;
      return true;
    }
    
    AnalyticalPropagator * clone() const override
    {
      return new AnalyticalPropagator(*this);
    }
    
    /** Set the maximum relative change in Bz (Bz_at_end-Bz_at_start)/Bz_at_start
     * for a single propagation. The default is no limit.
     * NB: this propagator assumes constant, non-zero magnetic field parallel to the z-axis!
     **/
    void setMaxRelativeChangeInBz (const float maxDBz) {
      theMaxDBzRatio = maxDBz;
    }
  
  
  private:
    /// propagation of errors (if needed) and generation of a new TSOS
    std::pair<TrajectoryStateOnSurface,double> 
    propagatedStateWithPath (const FreeTrajectoryState& fts, 
  			   const Surface& surface, 
  			   const GlobalTrajectoryParameters& gtp, 
  			   const double& s) const dso_internal;
  
    /// parameter propagation to cylinder (returns position, momentum and path length)
    bool propagateParametersOnCylinder(const FreeTrajectoryState& fts, 
  				     const Cylinder& cylinder, 
  				     GlobalPoint& x, 
  				     GlobalVector& p, 
  				     double& s) const dso_internal;
  
    /// parameter propagation to plane (returns position, momentum and path length)
    bool propagateParametersOnPlane(const FreeTrajectoryState& fts, 
  				  const Plane& plane, 
  				  GlobalPoint& x, 
  				  GlobalVector& p, 
  				  double& s) const dso_internal;
    
    /// straight line parameter propagation to a plane
    bool propagateWithLineCrossing(const GlobalPoint&, const GlobalVector&, 
  				 const Plane&, GlobalPoint&, double&) const dso_internal;
    /// straight line parameter propagation to a cylinder
    bool propagateWithLineCrossing(const GlobalPoint&, const GlobalVector&, 
  				 const Cylinder&, GlobalPoint&, double&) const dso_internal;
    /// helix parameter propagation to a plane using HelixPlaneCrossing
    bool propagateWithHelixCrossing(HelixPlaneCrossing&, const Plane&, const float,
  				  GlobalPoint&, GlobalVector&, double& s) const dso_internal;
  
    const MagneticField* magneticField() const override {return theField;}
  
  private:
    typedef std::pair<TrajectoryStateOnSurface,double> TsosWP;
    float theMaxDPhi2;
    float theMaxDBzRatio;
    const MagneticField* theField;
    bool isOldPropagationType;
  };
  
  #endif
  