File indexing completed on 2024-04-06 12:31:31
0001 #include "TrackingTools/GsfTracking/interface/GsfTrajectoryFitter.h"
0002
0003 #include "TrackingTools/TrackFitters/interface/TrajectoryStateWithArbitraryError.h"
0004 #include "TrackingTools/TransientTrackingRecHit/interface/TransientTrackingRecHit.h"
0005 #include "DataFormats/GeometrySurface/interface/BoundPlane.h"
0006 #include "TrackingTools/GsfTracking/interface/MultiTrajectoryStateMerger.h"
0007 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0008
0009 #include "DataFormats/TrackerRecHit2D/interface/TkCloner.h"
0010 #include "DataFormats/TrackerRecHit2D/interface/BaseTrackerRecHit.h"
0011
0012 #include "TrackingTools/TrackFitters/interface/DebugHelpers.h"
0013
0014 GsfTrajectoryFitter::GsfTrajectoryFitter(const Propagator& aPropagator,
0015 const TrajectoryStateUpdator& aUpdator,
0016 const MeasurementEstimator& aEstimator,
0017 const MultiTrajectoryStateMerger& aMerger,
0018 const DetLayerGeometry* detLayerGeometry)
0019 : thePropagator(aPropagator.clone()),
0020 theUpdator(aUpdator.clone()),
0021 theEstimator(aEstimator.clone()),
0022 theMerger(aMerger.clone()),
0023 theGeometry(detLayerGeometry) {
0024 if (!theGeometry)
0025 theGeometry = &dummyGeometry;
0026 }
0027
0028 GsfTrajectoryFitter::~GsfTrajectoryFitter() {
0029 delete thePropagator;
0030 delete theUpdator;
0031 delete theEstimator;
0032 delete theMerger;
0033 }
0034
0035 Trajectory GsfTrajectoryFitter::fitOne(const Trajectory& aTraj, fitType type) const {
0036 if (aTraj.empty())
0037 return Trajectory();
0038
0039 TM const& firstTM = aTraj.firstMeasurement();
0040 TSOS firstTsos = TrajectoryStateWithArbitraryError()(firstTM.updatedState());
0041
0042 return fitOne(aTraj.seed(), aTraj.recHits(), firstTsos, type);
0043 }
0044
0045 Trajectory GsfTrajectoryFitter::fitOne(const TrajectorySeed& aSeed, const RecHitContainer& hits, fitType type) const {
0046 edm::LogError("GsfTrajectoryFitter") << "GsfTrajectoryFitter::fit(TrajectorySeed, vector<RecHit>) not implemented";
0047
0048 return Trajectory();
0049 }
0050
0051 Trajectory GsfTrajectoryFitter::fitOne(const TrajectorySeed& aSeed,
0052 const RecHitContainer& hits,
0053 const TrajectoryStateOnSurface& firstPredTsos,
0054 fitType) const {
0055 if (hits.empty())
0056 return Trajectory();
0057
0058 Trajectory myTraj(aSeed, propagator()->propagationDirection());
0059
0060 TSOS predTsos(firstPredTsos);
0061 if (!predTsos.isValid()) {
0062 edm::LogInfo("GsfTrackFitters") << "GsfTrajectoryFitter: predicted tsos of first measurement not valid!";
0063 return Trajectory();
0064 }
0065
0066 TSOS currTsos;
0067 if (hits.front()->isValid()) {
0068 auto const& ihit = hits.front();
0069
0070 assert((!(ihit)->canImproveWithTrack()) || (nullptr != theHitCloner));
0071 assert((!(ihit)->canImproveWithTrack()) || (nullptr != dynamic_cast<BaseTrackerRecHit const*>(ihit.get())));
0072 auto preciseHit = theHitCloner->makeShared(ihit, predTsos);
0073 dump(*preciseHit, 1, "GsfTrackFitters");
0074 { currTsos = updator()->update(predTsos, *preciseHit); }
0075 if (!predTsos.isValid() || !currTsos.isValid()) {
0076 edm::LogError("InvalidState") << "first hit";
0077 return Trajectory();
0078 }
0079 myTraj.push(TM(predTsos, currTsos, preciseHit, 0., theGeometry->idToLayer(preciseHit->geographicalId())),
0080 estimator()->estimate(predTsos, *preciseHit).second);
0081 } else {
0082 currTsos = predTsos;
0083 if (!predTsos.isValid()) {
0084 edm::LogError("InvalidState") << "first invalid hit";
0085 return Trajectory();
0086 }
0087 myTraj.push(TM(predTsos, hits.front(), 0., theGeometry->idToLayer((hits.front())->geographicalId())));
0088 }
0089
0090 int hitcounter = 1;
0091 for (RecHitContainer::const_iterator ihit = hits.begin() + 1; ihit != hits.end(); ihit++) {
0092 ++hitcounter;
0093
0094
0095
0096
0097 if ((**ihit).isValid() == false && (**ihit).det() == nullptr) {
0098 LogDebug("GsfTrackFitters") << " Error: invalid hit with no GeomDet attached .... skipping";
0099 continue;
0100 }
0101
0102 {
0103
0104 predTsos = propagator()->propagate(currTsos, (**ihit).det()->surface());
0105 }
0106 if (!predTsos.isValid()) {
0107 if (myTraj.foundHits() >= 3) {
0108 edm::LogInfo("GsfTrackFitters") << "GsfTrajectoryFitter: predicted tsos not valid! \n"
0109 << "Returning trajectory with " << myTraj.foundHits() << " found hits.";
0110 return myTraj;
0111 } else {
0112 edm::LogInfo("GsfTrackFitters") << "GsfTrajectoryFitter: predicted tsos not valid after " << myTraj.foundHits()
0113 << " hits, discarding candidate!";
0114 return Trajectory();
0115 }
0116 }
0117 if (merger())
0118 predTsos = merger()->merge(predTsos);
0119
0120 if ((**ihit).isValid()) {
0121
0122 assert((!(*ihit)->canImproveWithTrack()) || (nullptr != theHitCloner));
0123 assert((!(*ihit)->canImproveWithTrack()) || (nullptr != dynamic_cast<BaseTrackerRecHit const*>((*ihit).get())));
0124 if (!predTsos.isValid()) {
0125 return Trajectory();
0126 }
0127 auto preciseHit = theHitCloner->makeShared(*ihit, predTsos);
0128 dump(*preciseHit, hitcounter, "GsfTrackFitters");
0129 currTsos = updator()->update(predTsos, *preciseHit);
0130 if (!predTsos.isValid() || !currTsos.isValid()) {
0131 edm::LogError("InvalidState") << "inside hit";
0132 return Trajectory();
0133 }
0134 auto chi2 = estimator()->estimate(predTsos, *preciseHit).second;
0135 myTraj.push(TM(predTsos, currTsos, preciseHit, chi2, theGeometry->idToLayer(preciseHit->geographicalId())));
0136 LogDebug("GsfTrackFitters") << "added measurement with chi2 " << chi2;
0137 } else {
0138 currTsos = predTsos;
0139 if (!predTsos.isValid()) {
0140 edm::LogError("InvalidState") << "inside invalid hit";
0141 return Trajectory();
0142 }
0143 myTraj.push(TM(predTsos, *ihit, 0., theGeometry->idToLayer((*ihit)->geographicalId())));
0144 }
0145 dump(predTsos, "predTsos", "GsfTrackFitters");
0146 dump(currTsos, "currTsos", "GsfTrackFitters");
0147 }
0148 return myTraj;
0149 }
