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
//
//

#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 final : 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);

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

  ~SteppingHelixPropagator() override;

  const MagneticField* magneticField() const override { return field_; }

  using Propagator::propagate;
  using Propagator::propagateWithPath;

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

  //----------------------------------------------------------------------------------------------------------

  void propagate(const SteppingHelixStateInfo& ftsStart, const Surface& sDest, SteppingHelixStateInfo& out) const;
  void propagate(const SteppingHelixStateInfo& ftsStart, const Plane& pDest, SteppingHelixStateInfo& out) const;
  void propagate(const SteppingHelixStateInfo& ftsStart, const Cylinder& cDest, SteppingHelixStateInfo& out) const;
  void propagate(const SteppingHelixStateInfo& ftsStart, const GlobalPoint& pDest, SteppingHelixStateInfo& out) const;
  void propagate(const SteppingHelixStateInfo& ftsStart,
                 const GlobalPoint& pDest1,
                 const GlobalPoint& pDest2,
                 SteppingHelixStateInfo& out) 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;
  static const int MAX_POINTS = 7;
  typedef StateInfo StateArray[MAX_POINTS + 1];

  void initStateArraySHPSpecific(StateArray& svBuf, bool flagsOnly) const;

  void setIState(const SteppingHelixStateInfo& sStart, StateArray& svBuff, int& nPoints) const;

  Result propagate(StateArray& svBuff,
                   int& nPoints,
                   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;
  //mutable int nPoints_;
  //mutable StateInfo svBuf_[MAX_POINTS+1];

  StateInfo invalidState_;

  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