ReferenceTrajectories/plugins/DualBzeroTrajectoryFactory.cc

0001 #include "FWCore/Framework/interface/ESHandle.h"
0002 #include "FWCore/Framework/interface/EventSetup.h"
0003 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0004 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0005 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0006 #include "Alignment/ReferenceTrajectories/interface/TrajectoryFactoryPlugin.h"
0007 #include "TrackingTools/GeomPropagators/interface/AnalyticalPropagator.h"
0008
0009 #include <algorithm>
0010
0011 #include "Alignment/ReferenceTrajectories/interface/DualBzeroReferenceTrajectory.h"
0012
0013 #include "Alignment/ReferenceTrajectories/interface/TrajectoryFactoryBase.h"
0014
0015 /// A factory that produces instances of class ReferenceTrajectory from a given TrajTrackPairCollection.
0016
0017 class DualBzeroTrajectoryFactory : public TrajectoryFactoryBase {
0018 public:
0019   DualBzeroTrajectoryFactory(const edm::ParameterSet &config, edm::ConsumesCollector &iC);
0020   ~DualBzeroTrajectoryFactory() override;
0021   const edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> m_MagFieldToken;
0022
0023   /// Produce the reference trajectories.
0024   const ReferenceTrajectoryCollection trajectories(const edm::EventSetup &setup,
0025                                                    const ConstTrajTrackPairCollection &tracks,
0026                                                    const reco::BeamSpot &beamSpot) const override;
0027
0028   const ReferenceTrajectoryCollection trajectories(const edm::EventSetup &setup,
0029                                                    const ConstTrajTrackPairCollection &tracks,
0030                                                    const ExternalPredictionCollection &external,
0031                                                    const reco::BeamSpot &beamSpot) const override;
0032
0033   DualBzeroTrajectoryFactory *clone() const override { return new DualBzeroTrajectoryFactory(*this); }
0034
0035 protected:
0036   struct DualBzeroTrajectoryInput {
0037     TrajectoryStateOnSurface refTsos;
0038     TransientTrackingRecHit::ConstRecHitContainer fwdRecHits;
0039     TransientTrackingRecHit::ConstRecHitContainer bwdRecHits;
0040   };
0041
0042   const DualBzeroTrajectoryInput referenceStateAndRecHits(const ConstTrajTrackPair &track) const;
0043
0044   const TrajectoryStateOnSurface propagateExternal(const TrajectoryStateOnSurface &external,
0045                                                    const Surface &surface,
0046                                                    const MagneticField *magField) const;
0047
0048   double theMass;
0049   double theMomentumEstimate;
0050 };
0051
0052 ////////////////////////////////////////////////////////////////////////////////
0053 ////////////////////////////////////////////////////////////////////////////////
0054 ////////////////////////////////////////////////////////////////////////////////
0055
0056 DualBzeroTrajectoryFactory::DualBzeroTrajectoryFactory(const edm::ParameterSet &config, edm::ConsumesCollector &iC)
0057     : TrajectoryFactoryBase(config, iC), m_MagFieldToken(iC.esConsumes()) {
0058   theMass = config.getParameter<double>("ParticleMass");
0059   theMomentumEstimate = config.getParameter<double>("MomentumEstimate");
0060 }
0061
0062 DualBzeroTrajectoryFactory::~DualBzeroTrajectoryFactory(void) {}
0063
0064 const DualBzeroTrajectoryFactory::ReferenceTrajectoryCollection DualBzeroTrajectoryFactory::trajectories(
0065     const edm::EventSetup &setup, const ConstTrajTrackPairCollection &tracks, const reco::BeamSpot &beamSpot) const {
0066   ReferenceTrajectoryCollection trajectories;
0067
0068   const MagneticField *magneticField = &setup.getData(m_MagFieldToken);
0069
0070   ConstTrajTrackPairCollection::const_iterator itTracks = tracks.begin();
0071
0072   while (itTracks != tracks.end()) {
0073     const DualBzeroTrajectoryInput input = this->referenceStateAndRecHits(*itTracks);
0074     // Check input: If all hits were rejected, the TSOS is initialized as invalid.
0075     if (input.refTsos.isValid()) {
0076       ReferenceTrajectoryBase::Config config(materialEffects(), propagationDirection(), theMass, theMomentumEstimate);
0077       config.useBeamSpot = useBeamSpot_;
0078       config.includeAPEs = includeAPEs_;
0079       config.allowZeroMaterial = allowZeroMaterial_;
0080       ReferenceTrajectoryPtr ptr(new DualBzeroReferenceTrajectory(
0081           input.refTsos, input.fwdRecHits, input.bwdRecHits, magneticField, beamSpot, config));
0082       trajectories.push_back(ptr);
0083     }
0084
0085     ++itTracks;
0086   }
0087
0088   return trajectories;
0089 }
0090
0091 const DualBzeroTrajectoryFactory::ReferenceTrajectoryCollection DualBzeroTrajectoryFactory::trajectories(
0092     const edm::EventSetup &setup,
0093     const ConstTrajTrackPairCollection &tracks,
0094     const ExternalPredictionCollection &external,
0095     const reco::BeamSpot &beamSpot) const {
0096   ReferenceTrajectoryCollection trajectories;
0097
0098   if (tracks.size() != external.size()) {
0099     edm::LogInfo("ReferenceTrajectories")
0100         << "@SUB=DualBzeroTrajectoryFactory::trajectories"
0101         << "Inconsistent input:\n"
0102         << "\tnumber of tracks = " << tracks.size() << "\tnumber of external predictions = " << external.size();
0103     return trajectories;
0104   }
0105
0106   const MagneticField *magneticField = &setup.getData(m_MagFieldToken);
0107
0108   ConstTrajTrackPairCollection::const_iterator itTracks = tracks.begin();
0109   ExternalPredictionCollection::const_iterator itExternal = external.begin();
0110
0111   while (itTracks != tracks.end()) {
0112     const DualBzeroTrajectoryInput input = referenceStateAndRecHits(*itTracks);
0113     // Check input: If all hits were rejected, the TSOS is initialized as invalid.
0114     if (input.refTsos.isValid()) {
0115       if ((*itExternal).isValid()) {
0116         TrajectoryStateOnSurface propExternal = propagateExternal(*itExternal, input.refTsos.surface(), magneticField);
0117
0118         if (!propExternal.isValid())
0119           continue;
0120
0121         ReferenceTrajectoryBase::Config config(materialEffects(), propagationDirection(), theMass, theMomentumEstimate);
0122         config.useBeamSpot = useBeamSpot_;
0123         config.includeAPEs = includeAPEs_;
0124         config.allowZeroMaterial = allowZeroMaterial_;
0125         ReferenceTrajectoryPtr ptr(new DualBzeroReferenceTrajectory(
0126             propExternal, input.fwdRecHits, input.bwdRecHits, magneticField, beamSpot, config));
0127
0128         AlgebraicSymMatrix externalParamErrors(asHepMatrix<5>(propExternal.localError().matrix()));
0129         ptr->setParameterErrors(externalParamErrors.sub(2, 5));
0130         trajectories.push_back(ptr);
0131       } else {
0132         ReferenceTrajectoryBase::Config config(materialEffects(), propagationDirection(), theMass, theMomentumEstimate);
0133         config.useBeamSpot = useBeamSpot_;
0134         config.includeAPEs = includeAPEs_;
0135         config.allowZeroMaterial = allowZeroMaterial_;
0136         ReferenceTrajectoryPtr ptr(new DualBzeroReferenceTrajectory(
0137             input.refTsos, input.fwdRecHits, input.bwdRecHits, magneticField, beamSpot, config));
0138
0139         trajectories.push_back(ptr);
0140       }
0141     }
0142
0143     ++itTracks;
0144     ++itExternal;
0145   }
0146
0147   return trajectories;
0148 }
0149
0150 const DualBzeroTrajectoryFactory::DualBzeroTrajectoryInput DualBzeroTrajectoryFactory::referenceStateAndRecHits(
0151     const ConstTrajTrackPair &track) const {
0152   DualBzeroTrajectoryInput input;
0153
0154   // get the trajectory measurements in the correct order, i.e. reverse if needed
0155   Trajectory::DataContainer allTrajMeas = this->orderedTrajectoryMeasurements(*track.first);
0156   Trajectory::DataContainer usedTrajMeas;
0157   Trajectory::DataContainer::iterator itM;
0158   // get all relevant trajectory measurements
0159   for (itM = allTrajMeas.begin(); itM != allTrajMeas.end(); itM++) {
0160     if (useRecHit((*itM).recHit()))
0161       usedTrajMeas.push_back(*itM);
0162   }
0163
0164   unsigned int iMeas = 0;
0165   unsigned int nMeas = usedTrajMeas.size();
0166   unsigned int nRefStateMeas = nMeas / 2;
0167   // get the valid RecHits
0168   for (itM = usedTrajMeas.begin(); itM != usedTrajMeas.end(); itM++, iMeas++) {
0169     TransientTrackingRecHit::ConstRecHitPointer aRecHit = (*itM).recHit();
0170
0171     if (iMeas < nRefStateMeas) {
0172       input.bwdRecHits.push_back(aRecHit);
0173     } else if (iMeas > nRefStateMeas) {
0174       input.fwdRecHits.push_back(aRecHit);
0175     } else {  // iMeas == nRefStateMeas
0176       if ((*itM).updatedState().isValid()) {
0177         input.refTsos = (*itM).updatedState();
0178         input.bwdRecHits.push_back(aRecHit);
0179         input.fwdRecHits.push_back(aRecHit);
0180       } else {
0181         // if the tsos of the middle hit is not valid, try the next one ...
0182         nRefStateMeas++;
0183         input.bwdRecHits.push_back(aRecHit);
0184       }
0185     }
0186   }
0187
0188   // bring input.fwdRecHits into correct order
0189   std::reverse(input.bwdRecHits.begin(), input.bwdRecHits.end());
0190
0191   return input;
0192 }
0193
0194 const TrajectoryStateOnSurface DualBzeroTrajectoryFactory::propagateExternal(const TrajectoryStateOnSurface &external,
0195                                                                              const Surface &surface,
0196                                                                              const MagneticField *magField) const {
0197   AnalyticalPropagator propagator(magField, anyDirection);
0198   const std::pair<TrajectoryStateOnSurface, double> tsosWithPath = propagator.propagateWithPath(external, surface);
0199   return tsosWithPath.first;
0200 }
0201
0202 DEFINE_EDM_PLUGIN(TrajectoryFactoryPlugin, DualBzeroTrajectoryFactory, "DualBzeroTrajectoryFactory");