SteppingHelixPropagator.h

Go to the documentation of this file.

#ifndef SteppingHelixPropagator_SteppingHelixPropagator_h
#define SteppingHelixPropagator_SteppingHelixPropagator_h



//
// Original Author:  Vyacheslav Krutelyov
//         Created:  Fri Mar  3 16:01:24 CST 2006
// $Id: SteppingHelixPropagator.h,v 1.30 2010/03/10 22:35:20 slava77 Exp $
//
//



#include "DataFormats/GeometryVector/interface/GlobalVector.h"
#include "DataFormats/GeometryVector/interface/GlobalPoint.h"
#include "TrackingTools/GeomPropagators/interface/Propagator.h"
#include "TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h"
#include "DataFormats/TrajectorySeed/interface/PropagationDirection.h"
#include "DataFormats/BeamSpot/interface/BeamSpot.h"

#include "CLHEP/Matrix/SymMatrix.h"
#include "CLHEP/Matrix/Matrix.h"
#include "CLHEP/Vector/ThreeVector.h"


#include "TrackPropagation/SteppingHelixPropagator/interface/SteppingHelixStateInfo.h"

class MagneticField;
class VolumeBasedMagneticField;
class MagVolume;

class SteppingHelixPropagator : public Propagator {
 public:
  typedef CLHEP::Hep3Vector Vector;
  typedef CLHEP::Hep3Vector  Point;

  typedef SteppingHelixStateInfo StateInfo;
  typedef SteppingHelixStateInfo::Result Result;

  struct Basis {
    Vector lX;
    Vector lY;
    Vector lZ;
  };
  
  enum Pars {
    RADIUS_P=0,
    Z_P = 0,
    PATHL_P = 0
  };

  enum DestType {
    RADIUS_DT=0,
    Z_DT,
    PLANE_DT,
    CONE_DT,
    CYLINDER_DT,
    PATHL_DT,
    POINT_PCA_DT,
    LINE_PCA_DT,
    UNDEFINED_DT
  };

  enum Fancy {
    HEL_AS_F=0, //simple analytical helix, eloss at end of step
    HEL_ALL_F,  //analytical helix with linear eloss
    POL_1_F, //1st order approximation, straight line
    POL_2_F,//2nd order
    POL_M_F //highest available
  };

  SteppingHelixPropagator();
  SteppingHelixPropagator(const MagneticField* field, PropagationDirection dir = alongMomentum);

  virtual SteppingHelixPropagator* clone() const {return new SteppingHelixPropagator(*this);}

  ~SteppingHelixPropagator() {}
  
  virtual const MagneticField* magneticField() const { return field_;}
  
  virtual TrajectoryStateOnSurface 
    propagate(const FreeTrajectoryState& ftsStart, const Plane& pDest) const;
  virtual TrajectoryStateOnSurface 
    propagate(const FreeTrajectoryState& ftsStart, const Cylinder& cDest) const;
  virtual FreeTrajectoryState 
    propagate(const FreeTrajectoryState& ftsStart, const GlobalPoint& pDest) const;
  virtual FreeTrajectoryState 
    propagate(const FreeTrajectoryState& ftsStart, 
          const GlobalPoint& pDest1, const GlobalPoint& pDest2) const;
  virtual FreeTrajectoryState 
    propagate(const FreeTrajectoryState& ftsStart, 
          const reco::BeamSpot& beamSpot) const;

  virtual std::pair<TrajectoryStateOnSurface, double> 
    propagateWithPath(const FreeTrajectoryState& ftsStart, const Plane& pDest) const;
  virtual std::pair<TrajectoryStateOnSurface, double> 
    propagateWithPath(const FreeTrajectoryState& ftsStart, const Cylinder& cDest) const;
  virtual std::pair<FreeTrajectoryState, double> 
    propagateWithPath(const FreeTrajectoryState& ftsStart, const GlobalPoint& pDest) const;
  virtual std::pair<FreeTrajectoryState, double> 
    propagateWithPath(const FreeTrajectoryState& ftsStart, 
              const GlobalPoint& pDest1, const GlobalPoint& pDest2) const;
  virtual std::pair<FreeTrajectoryState, double> 
    propagateWithPath(const FreeTrajectoryState& ftsStart, 
              const reco::BeamSpot& beamSpot) const;
    
    
  const SteppingHelixStateInfo& 
    propagate(const SteppingHelixStateInfo& ftsStart, const Surface& sDest) const;
  const SteppingHelixStateInfo& 
    propagate(const SteppingHelixStateInfo& ftsStart, const Plane& pDest) const;
  const SteppingHelixStateInfo& 
    propagate(const SteppingHelixStateInfo& ftsStart, const Cylinder& cDest) const;
  const SteppingHelixStateInfo& 
    propagate(const SteppingHelixStateInfo& ftsStart, const GlobalPoint& pDest) const;
  const SteppingHelixStateInfo& 
    propagate(const SteppingHelixStateInfo& ftsStart, 
          const GlobalPoint& pDest1, const GlobalPoint& pDest2) const;

  
  void setDebug(bool debug){ debug_ = debug;}

  void setMaterialMode(bool noMaterial) { noMaterialMode_ = noMaterial;}

  void setNoErrorPropagation(bool noErrorPropagation) { noErrorPropagation_ = noErrorPropagation;}

  void applyRadX0Correction(bool applyRadX0Correction) { applyRadX0Correction_ = applyRadX0Correction;}

  void setUseMagVolumes(bool val){ useMagVolumes_ = val;}

  void setUseMatVolumes(bool val){ useMatVolumes_ = val;}

  void setReturnTangentPlane(bool val){ returnTangentPlane_ = val;}

  void setSendLogWarning(bool val){ sendLogWarning_ = val;}

  void setUseIsYokeFlag(bool val){ useIsYokeFlag_ = val;}

  void setUseTuningForL2Speed(bool val){ useTuningForL2Speed_ = val;}

  void setVBFPointer(const VolumeBasedMagneticField* val) { vbField_ = val;}

  void setUseInTeslaFromMagField(bool val) { useInTeslaFromMagField_ = val;}

  void setEndcapShiftsInZPosNeg(double valPos, double valNeg){
    ecShiftPos_ = valPos; ecShiftNeg_ = valNeg;
  }

 protected:
  typedef SteppingHelixStateInfo::VolumeBounds MatBounds;
  void setIState(const SteppingHelixStateInfo& sStart) const;

  Result propagate(SteppingHelixPropagator::DestType type, 
           const double pars[6], double epsilon = 1e-3) const;

  void loadState(SteppingHelixPropagator::StateInfo& svCurrent, 
         const SteppingHelixPropagator::Vector& p3, 
         const SteppingHelixPropagator::Point& r3, int charge,
         PropagationDirection dir,
         const AlgebraicSymMatrix55& covCurv) const;

  void getNextState(const SteppingHelixPropagator::StateInfo& svPrevious,
            SteppingHelixPropagator::StateInfo& svNext,           
            double dP, 
            const SteppingHelixPropagator::Vector& tau, const SteppingHelixPropagator::Vector& drVec, 
            double dS, double dX0,
            const AlgebraicMatrix55& dCovCurv) const;
  
  void setRep(SteppingHelixPropagator::Basis& rep,
          const SteppingHelixPropagator::Vector& tau) const;

  bool makeAtomStep(SteppingHelixPropagator::StateInfo& svCurrent, 
            SteppingHelixPropagator::StateInfo& svNext, 
            double dS, PropagationDirection dir, 
            SteppingHelixPropagator::Fancy fancy) const;

  double getDeDx(const SteppingHelixPropagator::StateInfo& sv, 
         double& dEdXPrime, double& radX0, MatBounds& rzLims) const;

  int cIndex_(int ind) const;

  Result refToDest(DestType dest, const SteppingHelixPropagator::StateInfo& sv,
           const double pars[6], 
           double& dist, double& tanDist, PropagationDirection& refDirection,
           double fastSkipDist = 1e12) const;

  Result refToMagVolume(const SteppingHelixPropagator::StateInfo& sv,
            PropagationDirection dir,
            double& dist, double& tanDist,
            double fastSkipDist = 1e12, bool expectNewMagVolume = false, double maxStep = 1e12) const;

  Result refToMatVolume(const SteppingHelixPropagator::StateInfo& sv,
            PropagationDirection dir,
            double& dist, double& tanDist,
            double fastSkipDist = 1e12) const;

  bool isYokeVolume(const MagVolume* vol) const;


 private:
  typedef std::pair<TrajectoryStateOnSurface, double> TsosPP;
  typedef std::pair<FreeTrajectoryState, double> FtsPP;
  static const int MAX_STEPS = 10000;
  static const int MAX_POINTS = 7;
  mutable int nPoints_;
  mutable StateInfo svBuf_[MAX_POINTS+1];

  StateInfo invalidState_;

  mutable AlgebraicMatrix55 covCurvRot_;
  mutable AlgebraicMatrix55 dCCurvTransform_;

  const MagneticField* field_;
  const VolumeBasedMagneticField* vbField_;
  const AlgebraicSymMatrix55 unit55_;
  bool debug_;
  bool noMaterialMode_;
  bool noErrorPropagation_;
  bool applyRadX0Correction_;
  bool useMagVolumes_;
  bool useIsYokeFlag_;
  bool useMatVolumes_;
  bool useInTeslaFromMagField_;
  bool returnTangentPlane_;
  bool sendLogWarning_;
  bool useTuningForL2Speed_;

  double defaultStep_;

  double ecShiftPos_;
  double ecShiftNeg_;
};

#endif