TrackKinematicStatePropagator Class Reference

#include <TrackKinematicStatePropagator.h>

Inheritance diagram for TrackKinematicStatePropagator:

Public Member Functions
KinematicStatePropagator *	clone () const override
KinematicState	propagateToTheTransversePCA (const KinematicState &state, const GlobalPoint &referencePoint) const override
	TrackKinematicStatePropagator ()
bool	willPropagateToTheTransversePCA (const KinematicState &state, const GlobalPoint &point) const override
	~TrackKinematicStatePropagator () override
Public Member Functions inherited from KinematicStatePropagator
	KinematicStatePropagator ()
virtual	~KinematicStatePropagator ()

Private Types
typedef Point3DBase< double, GlobalTag >	GlobalPointDouble
typedef Vector3DBase< double, GlobalTag >	GlobalVectorDouble

Private Member Functions
virtual KinematicState	propagateToTheTransversePCACharged (const KinematicState &state, const GlobalPoint &referencePoint) const
virtual KinematicState	propagateToTheTransversePCANeutral (const KinematicState &state, const GlobalPoint &referencePoint) const

Detailed Description

Propagator for TransientTrack based KinematicStates. Does not include the material.

Definition at line 17 of file TrackKinematicStatePropagator.h.

Member Typedef Documentation

typedef Point3DBase<double, GlobalTag> TrackKinematicStatePropagator::GlobalPointDouble

private

Definition at line 49 of file TrackKinematicStatePropagator.h.

typedef Vector3DBase<double, GlobalTag> TrackKinematicStatePropagator::GlobalVectorDouble

private

Definition at line 50 of file TrackKinematicStatePropagator.h.

Constructor & Destructor Documentation

TrackKinematicStatePropagator::TrackKinematicStatePropagator ( )

inline

Definition at line 19 of file TrackKinematicStatePropagator.h.

Referenced by clone().

{}

TrackKinematicStatePropagator::~TrackKinematicStatePropagator ( )

inlineoverride

Definition at line 21 of file TrackKinematicStatePropagator.h.

References point, propagateToTheTransversePCA(), and willPropagateToTheTransversePCA().

{}

Member Function Documentation

KinematicStatePropagator* TrackKinematicStatePropagator::clone ( void  ) const

inlineoverridevirtual

Clone method reimplemented from abstract class

Implements KinematicStatePropagator.

Definition at line 36 of file TrackKinematicStatePropagator.h.

References propagateToTheTransversePCACharged(), propagateToTheTransversePCANeutral(), and TrackKinematicStatePropagator().

{ return new TrackKinematicStatePropagator(*this); }

KinematicState TrackKinematicStatePropagator::propagateToTheTransversePCA ( const KinematicState &  state, const GlobalPoint &  referencePoint  ) const

overridevirtual

Propagation to the point of closest approach in transverse plane to the given point

Implements KinematicStatePropagator.

Definition at line 11 of file TrackKinematicStatePropagator.cc.

References anyDirection, Plane::build(), FreeTrajectoryState::charge(), kappa, GlobalTrajectoryParameters::magneticFieldInInverseGeV(), FreeTrajectoryState::momentum(), FreeTrajectoryState::parameters(), KinematicState::particleCharge(), point, FreeTrajectoryState::position(), makeMuonMisalignmentScenario::rot, FreeTrajectoryState::transverseCurvature(), Vector3DBase< T, FrameTag >::unit(), X, PV3DBase< T, PVType, FrameType >::x(), DOFs::Y, PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by ~TrackKinematicStatePropagator().

                                                                                                                   {
  if (state.particleCharge() == 0.) {
    return propagateToTheTransversePCANeutral(state, referencePoint);
  } else {
    return propagateToTheTransversePCACharged(state, referencePoint);
  }
}

KinematicState TrackKinematicStatePropagator::propagateToTheTransversePCACharged ( const KinematicState &  state, const GlobalPoint &  referencePoint  ) const

privatevirtual

Internal private methods, distinguishing between the propagation of neutrals and propagation of cahrged tracks.

Definition at line 68 of file TrackKinematicStatePropagator.cc.

References FreeTrajectoryState::curvilinearError(), HelixBarrelPlaneCrossingByCircle::direction(), KinematicState::freeTrajectoryState(), hcaldqm::quantity::fState, KinematicState::globalMomentum(), mps_fire::i, JacobianCartesianToCurvilinear::jacobian(), JacobianCurvilinearToCartesian::jacobian(), AnalyticalCurvilinearJacobian::jacobian(), KinematicState::kinematicParametersError(), LogDebug, PV3DBase< T, PVType, FrameType >::mag(), KinematicState::magneticField(), EgHLTOffHistBins_cfi::mass, KinematicState::mass(), KinematicParametersError::matrix(), CurvilinearTrajectoryError::matrix(), KinematicState::particleCharge(), HelixBarrelPlaneCrossingByCircle::pathLength(), HelixBarrelPlaneCrossingByCircle::position(), alignCSCRings::s, mathSSE::sqrt(), KinematicState::trajectoryParameters(), PV3DBase< T, PVType, FrameType >::x(), Basic3DVector< T >::x(), Basic3DVector< T >::y(), and Basic3DVector< T >::z().

Referenced by clone().

                                                                          {
  //first use the existing FTS propagator to obtain parameters at PCA
  //in transverse plane to the given point

  //final parameters and covariance

  //initial parameters as class and vectors:
  //making a free trajectory state and looking
  //for helix barrel plane crossing
  FreeTrajectoryState const& fState = state.freeTrajectoryState();
  const GlobalPoint& iP = referencePoint;
  std::pair<HelixBarrelPlaneCrossingByCircle, BoundPlane::BoundPlanePointer> cros = planeCrossing(fState, iP);

  HelixBarrelPlaneCrossingByCircle planeCrossing = cros.first;
  BoundPlane::BoundPlanePointer plane = cros.second;
  std::pair<bool, double> propResult = planeCrossing.pathLength(*plane);
  if (!propResult.first) {
    LogDebug("RecoVertex/TrackKinematicStatePropagator") << "Propagation failed! State is invalid\n";
    return KinematicState();
  }
  double s = propResult.second;

  GlobalTrajectoryParameters const& inPar = state.trajectoryParameters();
  ParticleMass mass = state.mass();
  GlobalVector inMom = state.globalMomentum();

  HelixPlaneCrossing::PositionType xGen = planeCrossing.position(s);
  GlobalPoint nPosition(xGen.x(), xGen.y(), xGen.z());
  HelixPlaneCrossing::DirectionType pGen = planeCrossing.direction(s);
  pGen *= inMom.mag() / pGen.mag();
  GlobalVector nMomentum(pGen.x(), pGen.y(), pGen.z());
  AlgebraicVector7 par;
  AlgebraicSymMatrix77 cov;
  par(0) = nPosition.x();
  par(1) = nPosition.y();
  par(2) = nPosition.z();
  par(3) = nMomentum.x();
  par(4) = nMomentum.y();
  par(5) = nMomentum.z();
  par(6) = mass;

  //covariance matrix business
  //elements of 7x7 covariance matrix responcible for the mass and
  //mass - momentum projections corellations do change under such a transformation:
  //special Jacobian needed
  GlobalTrajectoryParameters fPar(nPosition, nMomentum, state.particleCharge(), state.magneticField());

  // check if correlation are present between mass and others
  bool thereIsNoCorr = true;

  for (auto i = 0; i < 6; ++i)
    thereIsNoCorr &= (0 == state.kinematicParametersError().matrix()(i, 6));

  if (thereIsNoCorr) {
    //  easy life
    AnalyticalCurvilinearJacobian prop(inPar, nPosition, nMomentum, s);
    AlgebraicSymMatrix55 cov2 = ROOT::Math::Similarity(prop.jacobian(), fState.curvilinearError().matrix());
    FreeTrajectoryState fts(fPar, CurvilinearTrajectoryError(cov2));

    return KinematicState(fts, state.mass(), std::sqrt(state.kinematicParametersError().matrix()(6, 6)));

    //KinematicState kRes(fts, state.mass(), std::sqrt(state.kinematicParametersError().matrix()(6,6)));
    //std::cout << "\n\ncart from final Kstate\n" << kRes.kinematicParametersError().matrix() << std::endl;
    // std::cout << "curv from final K\n" << kRes.freeTrajectoryState().curvilinearError().matrix() << std::endl;

  } else {
    JacobianCartesianToCurvilinear cart2curv(inPar);
    JacobianCurvilinearToCartesian curv2cart(fPar);

    AlgebraicMatrix67 ca2cu;
    AlgebraicMatrix76 cu2ca;
    ca2cu.Place_at(cart2curv.jacobian(), 0, 0);
    cu2ca.Place_at(curv2cart.jacobian(), 0, 0);
    ca2cu(5, 6) = 1;
    cu2ca(6, 5) = 1;

    //now both transformation jacobians: cartesian to curvilinear and back are done
    //We transform matrix to curv frame, then propagate it and translate it back to
    //cartesian frame.
    AlgebraicSymMatrix66 cov1 = ROOT::Math::Similarity(ca2cu, state.kinematicParametersError().matrix());

    /*
      std::cout << "\n\ncurv from Kstate\n" << cov1 << std::endl;
      std::cout << "curv from fts\n" << fState.curvilinearError().matrix() << std::endl;
    */

    //propagation jacobian
    AnalyticalCurvilinearJacobian prop(inPar, nPosition, nMomentum, s);
    AlgebraicMatrix66 pr;
    pr(5, 5) = 1;
    pr.Place_at(prop.jacobian(), 0, 0);

    //transportation
    AlgebraicSymMatrix66 cov2 = ROOT::Math::Similarity(pr, cov1);

    //now geting back to 7-parametrization from curvilinear
    cov = ROOT::Math::Similarity(cu2ca, cov2);

    /*
      std::cout << "curv prop \n" << cov2 << std::endl;
   std::cout << "cart prop\n" << cov << std::endl;
    */

    //return parameters as a kiematic state
    KinematicParameters resPar(par);
    KinematicParametersError resEr(cov);

    return KinematicState(resPar, resEr, state.particleCharge(), state.magneticField());

    /*
    KinematicState resK(resPar,resEr,state.particleCharge(), state.magneticField());
    std::cout << "\ncurv from K prop\n" << resK.freeTrajectoryState().curvilinearError().matrix() << std::endl;
    return resK;
  */
  }
}

KinematicState TrackKinematicStatePropagator::propagateToTheTransversePCANeutral ( const KinematicState &  state, const GlobalPoint &  referencePoint  ) const

privatevirtual

Definition at line 186 of file TrackKinematicStatePropagator.cc.

References funct::cos(), KinematicState::freeTrajectoryState(), hcaldqm::quantity::fState, JacobianCurvilinearToCartesian::jacobian(), KinematicState::kinematicParametersError(), KinematicState::magneticField(), KinematicState::mass(), KinematicParametersError::matrix(), FreeTrajectoryState::momentum(), KinematicState::particleCharge(), PV3DBase< T, PVType, FrameType >::phi(), FreeTrajectoryState::position(), alignCSCRings::s, funct::sin(), funct::tan(), PV3DBase< T, PVType, FrameType >::theta(), theta(), KinematicState::trajectoryParameters(), PV3DBase< T, PVType, FrameType >::x(), PV3DBase< T, PVType, FrameType >::y(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by clone().

                                                                          {
  //new parameters vector and covariance:
  AlgebraicVector7 par;
  AlgebraicSymMatrix77 cov;

  //AlgebraicVector7 inStatePar = state.kinematicParameters().vector();
  GlobalTrajectoryParameters const& inPar = state.trajectoryParameters();

  //first making a free trajectory state and propagating it according
  //to the algorithm provided by Thomas Speer and Wolfgang Adam
  FreeTrajectoryState const& fState = state.freeTrajectoryState();

  GlobalPoint xvecOrig = fState.position();
  double phi = fState.momentum().phi();
  double theta = fState.momentum().theta();
  double xOrig = xvecOrig.x();
  double yOrig = xvecOrig.y();
  double zOrig = xvecOrig.z();
  double xR = referencePoint.x();
  double yR = referencePoint.y();

  double s2D = (xR - xOrig) * cos(phi) + (yR - yOrig) * sin(phi);
  double s = s2D / sin(theta);
  double xGen = xOrig + s2D * cos(phi);
  double yGen = yOrig + s2D * sin(phi);
  double zGen = zOrig + s2D / tan(theta);
  GlobalPoint xPerigee = GlobalPoint(xGen, yGen, zGen);

  //new parameters
  GlobalVector pPerigee = fState.momentum();
  par(0) = xPerigee.x();
  par(1) = xPerigee.y();
  par(2) = xPerigee.z();
  par(3) = pPerigee.x();
  par(4) = pPerigee.y();
  par(5) = pPerigee.z();
  // par(6) = inStatePar(7);
  par(6) = state.mass();

  //covariance matrix business:
  //everything lake it was before: jacobains are smart enouhg to
  //distinguish between neutral and charged states themselves

  GlobalTrajectoryParameters fPar(xPerigee, pPerigee, state.particleCharge(), state.magneticField());

  JacobianCartesianToCurvilinear cart2curv(inPar);
  JacobianCurvilinearToCartesian curv2cart(fPar);

  AlgebraicMatrix67 ca2cu;
  AlgebraicMatrix76 cu2ca;
  ca2cu.Place_at(cart2curv.jacobian(), 0, 0);
  cu2ca.Place_at(curv2cart.jacobian(), 0, 0);
  ca2cu(5, 6) = 1;
  cu2ca(6, 5) = 1;

  //now both transformation jacobians: cartesian to curvilinear and back are done
  //We transform matrix to curv frame, then propagate it and translate it back to
  //cartesian frame.
  AlgebraicSymMatrix66 cov1 = ROOT::Math::Similarity(ca2cu, state.kinematicParametersError().matrix());

  //propagation jacobian
  AnalyticalCurvilinearJacobian prop(inPar, xPerigee, pPerigee, s);
  AlgebraicMatrix66 pr;
  pr(5, 5) = 1;
  pr.Place_at(prop.jacobian(), 0, 0);

  //transportation
  AlgebraicSymMatrix66 cov2 = ROOT::Math::Similarity(pr, cov1);

  //now geting back to 7-parametrization from curvilinear
  cov = ROOT::Math::Similarity(cu2ca, cov2);

  // FreeTrajectoryState fts(fPar);

  //return parameters as a kiematic state
  KinematicParameters resPar(par);
  KinematicParametersError resEr(cov);
  return KinematicState(resPar, resEr, state.particleCharge(), state.magneticField());

  //return  KinematicState(fts,state.mass(), cov(6,6));
}

bool TrackKinematicStatePropagator::willPropagateToTheTransversePCA ( const KinematicState &  state, const GlobalPoint &  point  ) const

overridevirtual

Reimplemented from KinematicStatePropagator.

Definition at line 53 of file TrackKinematicStatePropagator.cc.

References KinematicState::freeTrajectoryState(), hcaldqm::quantity::fState, KinematicState::particleCharge(), and HelixBarrelPlaneCrossingByCircle::pathLength().

Referenced by KinematicParticleVertexFitter::fit(), and ~TrackKinematicStatePropagator().

                                                                                                    {
  if (state.particleCharge() == 0.)
    return true;

  // copied from below...
  FreeTrajectoryState const& fState = state.freeTrajectoryState();
  std::pair<HelixBarrelPlaneCrossingByCircle, BoundPlane::BoundPlanePointer> cros = planeCrossing(fState, point);

  HelixBarrelPlaneCrossingByCircle planeCrossing = cros.first;
  BoundPlane::BoundPlanePointer plane = cros.second;
  std::pair<bool, double> propResult = planeCrossing.pathLength(*plane);
  return propResult.first;
}