Go to the documentation of this file.00001 #include "TrackPropagation/RungeKutta/interface/RKPropagatorInR.h"
00002 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h"
00003
00004 #include "RKCylindricalDistance.h"
00005 #include "CylindricalLorentzForce.h"
00006 #include "RKLocalFieldProvider.h"
00007 #include "DataFormats/GeometrySurface/interface/Cylinder.h"
00008 #include "RKAdaptiveSolver.h"
00009 #include "RKOne4OrderStep.h"
00010 #include "RKOneCashKarpStep.h"
00011 #include "CylindricalState.h"
00012 #include "MagneticField/VolumeGeometry/interface/MagVolume.h"
00013
00014 TrajectoryStateOnSurface
00015 RKPropagatorInR::propagate (const FreeTrajectoryState& ts, const Cylinder& cyl) const
00016 {
00017
00018 typedef RKAdaptiveSolver<double,RKOne4OrderStep, 5> Solver;
00019
00020 typedef double Scalar;
00021 typedef Solver::Vector RKVector;
00022
00023 GlobalPoint pos( ts.position());
00024 GlobalVector mom( ts.momentum());
00025
00026 LocalPoint startpos = cyl.toLocal(pos);
00027 LocalVector startmom = cyl.toLocal(mom);
00028
00029 CylindricalState startState( startpos, startmom, ts.charge());
00030 RKVector start = startState.parameters();
00031
00032 RKLocalFieldProvider localField( *theVolume, cyl);
00033 CylindricalLorentzForce<double,5> deriv(localField);
00034 RKCylindricalDistance<double,5> dist;
00035 double eps = 1.e-5;
00036 Solver solver;
00037 try {
00038 Scalar step = cyl.radius() - startState.rho();
00039 RKVector rkresult = solver( startState.rho(), start, step, deriv, dist, eps);
00040 CylindricalState endState( cyl.radius(), rkresult, startState.prSign());
00041 return TrajectoryStateOnSurface( GlobalTrajectoryParameters( cyl.toGlobal( endState.position()),
00042 cyl.toGlobal( endState.momentum()),
00043 TrackCharge(endState.charge()),
00044 theVolume),
00045 cyl);
00046 }
00047 catch (CylindricalLorentzForceException& e) {
00048
00049
00050
00051 return TrajectoryStateOnSurface();
00052 }
00053 }
00054
00055 TrajectoryStateOnSurface
00056 RKPropagatorInR::propagate (const FreeTrajectoryState&, const Plane&) const
00057 {
00058 return TrajectoryStateOnSurface();
00059 }
00060
00061 std::pair< TrajectoryStateOnSurface, double>
00062 RKPropagatorInR::propagateWithPath (const FreeTrajectoryState&, const Plane&) const
00063 {
00064 return std::pair< TrajectoryStateOnSurface, double>();
00065 }
00066
00067 std::pair< TrajectoryStateOnSurface, double>
00068 RKPropagatorInR::propagateWithPath (const FreeTrajectoryState&, const Cylinder&) const
00069 {
00070 return std::pair< TrajectoryStateOnSurface, double>();
00071 }
00072
00073 Propagator * RKPropagatorInR::clone() const
00074 {
00075 return new RKPropagatorInR(*this);
00076 }