MaterialEffects/src/PropagatorWithMaterial.cc

0001 #include "TrackingTools/GeomPropagators/interface/PropagationDirectionFromPath.h"
0002 #include "TrackingTools/MaterialEffects/interface/PropagatorWithMaterial.h"
0003 #include "TrackingTools/GeomPropagators/interface/AnalyticalPropagator.h"
0004 #include "TrackingTools/MaterialEffects/interface/CombinedMaterialEffectsUpdator.h"
0005 #include "FWCore/Utilities/interface/Exception.h"
0006 #include <string>
0007
0008 using namespace std;
0009
0010 PropagatorWithMaterial::~PropagatorWithMaterial() {}
0011
0012 PropagatorWithMaterial::PropagatorWithMaterial(PropagationDirection dir,
0013                                                const float mass,
0014                                                const MagneticField* mf,
0015                                                const float maxDPhi,
0016                                                bool useRungeKutta,
0017                                                float ptMin,
0018                                                bool useOldAnalPropLogic)
0019     : Propagator(dir),
0020       rkProduct(mf, dir),
0021       theGeometricalPropagator(useRungeKutta ? rkProduct.propagator.clone()
0022                                              : new AnalyticalPropagator(mf, dir, maxDPhi, useOldAnalPropLogic)),
0023       theMEUpdator(new CombinedMaterialEffectsUpdator(mass, ptMin)),
0024       theMaterialLocation(atDestination),
0025       field(mf),
0026       useRungeKutta_(useRungeKutta) {}
0027
0028 pair<TrajectoryStateOnSurface, double> PropagatorWithMaterial::propagateWithPath(const FreeTrajectoryState& fts,
0029                                                                                  const Plane& plane) const {
0030   TsosWP newTsosWP = theGeometricalPropagator->propagateWithPath(fts, plane);
0031   if ((newTsosWP.first).isValid() && !materialAtSource()) {
0032     bool updateOk = theMEUpdator->updateStateInPlace(
0033         newTsosWP.first, PropagationDirectionFromPath()(newTsosWP.second, propagationDirection()));
0034     if UNLIKELY (!updateOk)
0035       newTsosWP.first = TrajectoryStateOnSurface();
0036   }
0037   return newTsosWP;
0038 }
0039
0040 pair<TrajectoryStateOnSurface, double> PropagatorWithMaterial::propagateWithPath(const FreeTrajectoryState& fts,
0041                                                                                  const Cylinder& cylinder) const {
0042   TsosWP newTsosWP = theGeometricalPropagator->propagateWithPath(fts, cylinder);
0043   if ((newTsosWP.first).isValid() && !materialAtSource()) {
0044     bool updateOk = theMEUpdator->updateStateInPlace(
0045         newTsosWP.first, PropagationDirectionFromPath()(newTsosWP.second, propagationDirection()));
0046     if UNLIKELY (!updateOk)
0047       newTsosWP.first = TrajectoryStateOnSurface();
0048   }
0049   return newTsosWP;
0050 }
0051
0052 pair<TrajectoryStateOnSurface, double> PropagatorWithMaterial::propagateWithPath(const TrajectoryStateOnSurface& tsos,
0053                                                                                  const Plane& plane) const {
0054   //
0055   // add material at starting surface, if requested
0056   //
0057   TsosWP newTsosWP(tsos, 0.);
0058   if (materialAtSource()) {
0059     bool updateOk = theMEUpdator->updateStateInPlace(newTsosWP.first, propagationDirection());
0060     if UNLIKELY (!updateOk)
0061       newTsosWP.first = TrajectoryStateOnSurface();
0062   }
0063   if UNLIKELY (!newTsosWP.first.isValid())
0064     return newTsosWP;
0065   //
0066   // geometrical propagation
0067   //
0068   newTsosWP = theGeometricalPropagator->propagateWithPath(newTsosWP.first, plane);
0069   if UNLIKELY (!newTsosWP.first.isValid() || materialAtSource())
0070     return newTsosWP;
0071   //
0072   // add material at destination surface, if requested
0073   //
0074   bool updateOk = theMEUpdator->updateStateInPlace(
0075       newTsosWP.first, PropagationDirectionFromPath()(newTsosWP.second, propagationDirection()));
0076   if UNLIKELY (!updateOk)
0077     newTsosWP.first = TrajectoryStateOnSurface();
0078   return newTsosWP;
0079 }
0080
0081 pair<TrajectoryStateOnSurface, double> PropagatorWithMaterial::propagateWithPath(const TrajectoryStateOnSurface& tsos,
0082                                                                                  const Cylinder& cylinder) const {
0083   //
0084   // add material at starting surface, if requested
0085   //
0086   TsosWP newTsosWP(tsos, 0.);
0087   if (materialAtSource()) {
0088     bool updateOk = theMEUpdator->updateStateInPlace(newTsosWP.first, propagationDirection());
0089     if UNLIKELY (!updateOk)
0090       newTsosWP.first = TrajectoryStateOnSurface();
0091   }
0092   if UNLIKELY (!newTsosWP.first.isValid())
0093     return newTsosWP;
0094   //
0095   // geometrical propagation
0096   //
0097   newTsosWP = theGeometricalPropagator->propagateWithPath(newTsosWP.first, cylinder);
0098   if UNLIKELY (!(newTsosWP.first).isValid() || materialAtSource())
0099     return newTsosWP;
0100   //
0101   // add material at destination surface, if requested
0102   //
0103   bool updateOk = theMEUpdator->updateStateInPlace(
0104       newTsosWP.first, PropagationDirectionFromPath()(newTsosWP.second, propagationDirection()));
0105   if UNLIKELY (!updateOk)
0106     newTsosWP.first = TrajectoryStateOnSurface();
0107   return newTsosWP;
0108 }
0109
0110 void PropagatorWithMaterial::setPropagationDirection(PropagationDirection dir) {
0111   theGeometricalPropagator->setPropagationDirection(dir);
0112   Propagator::setPropagationDirection(dir);
0113 }
0114
0115 bool PropagatorWithMaterial::materialAtSource() const {
0116   if UNLIKELY ((propagationDirection() == anyDirection) && (theMaterialLocation != atDestination))
0117     throw cms::Exception("TrackingTools/MaterialEffects",
0118                          "PropagatorWithMaterial: propagation direction = anyDirection is incompatible with adding of "
0119                          "material at source");
0120
0121   return theMaterialLocation == atSource ||
0122          (theMaterialLocation == fromDirection && propagationDirection() == alongMomentum);
0123 }