#include <KFTrajectorySmoother.h>

Inheritance diagram for KFTrajectorySmoother:

Public Member Functions
const Propagator *	alongPropagator () const

KFTrajectorySmoother *	clone () const override

const MeasurementEstimator *	estimator () const

	KFTrajectorySmoother (const Propagator &aPropagator, const TrajectoryStateUpdator &aUpdator, const MeasurementEstimator &aEstimator, float errorRescaling=100.f, int minHits=3)

	KFTrajectorySmoother (const Propagator aPropagator, const TrajectoryStateUpdator aUpdator, const MeasurementEstimator aEstimator, float errorRescaling=100.f, int minHits=3, const DetLayerGeometry detLayerGeometry=nullptr, TkCloner const *hc=nullptr)

const Propagator *	oppositePropagator () const

void	setHitCloner (TkCloner const *hc) override

Trajectory	trajectory (const Trajectory &aTraj) const override

const TrajectoryStateUpdator *	updator () const

	~KFTrajectorySmoother () override

Public Member Functions inherited from TrajectorySmoother
virtual TrajectoryContainer	trajectories (const Trajectory &traj) const

virtual	~TrajectorySmoother ()

Private Types
typedef FreeTrajectoryState	FTS

typedef TrajectoryMeasurement	TM

typedef TrajectoryStateOnSurface	TSOS

Private Attributes
const DetLayerGeometry	dummyGeometry

int	minHits_

const Propagator *	theAlongPropagator

float	theErrorRescaling

const MeasurementEstimator *	theEstimator

const DetLayerGeometry *	theGeometry

TkCloner const *	theHitCloner = nullptr

const Propagator *	theOppositePropagator

const TrajectoryStateUpdator *	theUpdator

Additional Inherited Members
Public Types inherited from TrajectorySmoother
typedef std::vector< Trajectory >	TrajectoryContainer

typedef TrajectoryContainer::iterator	TrajectoryIterator

Detailed Description

A Standard Kalman smoother. The forward fit is not redone, only the backward smoothing. Ported from ORCA

Author: todorov, cerati

Definition at line 20 of file KFTrajectorySmoother.h.

Member Typedef Documentation

typedef FreeTrajectoryState KFTrajectorySmoother::FTS

private

Definition at line 23 of file KFTrajectorySmoother.h.

typedef TrajectoryMeasurement KFTrajectorySmoother::TM

private

Definition at line 24 of file KFTrajectorySmoother.h.

typedef TrajectoryStateOnSurface KFTrajectorySmoother::TSOS

private

Definition at line 22 of file KFTrajectorySmoother.h.

Constructor & Destructor Documentation

KFTrajectorySmoother::KFTrajectorySmoother	(	const Propagator &	aPropagator,
		const TrajectoryStateUpdator &	aUpdator,
		const MeasurementEstimator &	aEstimator,
		float	errorRescaling = `100.f`,
		int	minHits = `3`
	)

inline

Definition at line 27 of file KFTrajectorySmoother.h.

References alongMomentum, Propagator::clone(), dummyGeometry, oppositeToMomentum, AlCaHLTBitMon_ParallelJobs::p, theAlongPropagator, theGeometry, and theOppositePropagator.

Referenced by clone().

       : theAlongPropagator(nullptr),
         theOppositePropagator(nullptr),
         theUpdator(aUpdator.clone()),
         theEstimator(aEstimator.clone()),
         theErrorRescaling(errorRescaling),
         minHits_(minHits),
         // to be fixed. Why this first constructor is needed? who is using it? Can it be removed?
         theGeometry(nullptr) {
     if (!theGeometry)
       theGeometry = &dummyGeometry;
     auto p = aPropagator.clone();
     p->setPropagationDirection(alongMomentum);
     theAlongPropagator = p;
     p = aPropagator.clone();
     p->setPropagationDirection(oppositeToMomentum);
     theOppositePropagator = p;
   }

KFTrajectorySmoother::KFTrajectorySmoother	(	const Propagator *	aPropagator,
		const TrajectoryStateUpdator *	aUpdator,
		const MeasurementEstimator *	aEstimator,
		float	errorRescaling = `100.f`,
		int	minHits = `3`,
		const DetLayerGeometry *	detLayerGeometry = `nullptr`,
		TkCloner const *	hc = `nullptr`
	)

inline

Definition at line 50 of file KFTrajectorySmoother.h.

References alongMomentum, Propagator::clone(), dummyGeometry, oppositeToMomentum, AlCaHLTBitMon_ParallelJobs::p, theAlongPropagator, theGeometry, and theOppositePropagator.

       : theAlongPropagator(nullptr),
         theOppositePropagator(nullptr),
         theUpdator(aUpdator->clone()),
         theEstimator(aEstimator->clone()),
         theHitCloner(hc),
         theErrorRescaling(errorRescaling),
         minHits_(minHits),
         theGeometry(detLayerGeometry) {
     if (!theGeometry)
       theGeometry = &dummyGeometry;
     auto p = aPropagator->clone();
     p->setPropagationDirection(alongMomentum);
     theAlongPropagator = p;
     p = aPropagator->clone();
     p->setPropagationDirection(oppositeToMomentum);
     theOppositePropagator = p;
   }

KFTrajectorySmoother::~KFTrajectorySmoother ( )

override

Definition at line 12 of file KFTrajectorySmoother.cc.

References theAlongPropagator, theEstimator, theOppositePropagator, and theUpdator.

                                             {
   delete theAlongPropagator;
   delete theOppositePropagator;
   delete theUpdator;
   delete theEstimator;
 }

Member Function Documentation

const Propagator* KFTrajectorySmoother::alongPropagator ( ) const

inline

Definition at line 79 of file KFTrajectorySmoother.h.

References theAlongPropagator.

79 { return theAlongPropagator; }

KFTrajectorySmoother::theAlongPropagator

const Propagator * theAlongPropagator

Definition: KFTrajectorySmoother.h:95

KFTrajectorySmoother* KFTrajectorySmoother::clone ( void ) const

inlineoverridevirtual

Implements TrajectorySmoother.

Definition at line 85 of file KFTrajectorySmoother.h.

References KFTrajectorySmoother(), minHits_, theAlongPropagator, theErrorRescaling, theEstimator, theGeometry, theHitCloner, and theUpdator.

                                                {
     return new KFTrajectorySmoother(
         theAlongPropagator, theUpdator, theEstimator, theErrorRescaling, minHits_, theGeometry, theHitCloner);
   }

const MeasurementEstimator* KFTrajectorySmoother::estimator ( void ) const

inline

Definition at line 83 of file KFTrajectorySmoother.h.

References theEstimator.

Referenced by trajectory().

83 { return theEstimator; }

KFTrajectorySmoother::theEstimator

const MeasurementEstimator * theEstimator

Definition: KFTrajectorySmoother.h:98

const Propagator* KFTrajectorySmoother::oppositePropagator ( ) const

inline

Definition at line 80 of file KFTrajectorySmoother.h.

References theOppositePropagator.

80 { return theOppositePropagator; }

KFTrajectorySmoother::theOppositePropagator

const Propagator * theOppositePropagator

Definition: KFTrajectorySmoother.h:96

void KFTrajectorySmoother::setHitCloner ( TkCloner const * hc )

inlineoverridevirtual

Implements TrajectorySmoother.

Definition at line 91 of file KFTrajectorySmoother.h.

References theHitCloner.

Referenced by CosmicTrajectoryBuilder::init().

91 { theHitCloner = hc; }

KFTrajectorySmoother::theHitCloner

TkCloner const * theHitCloner

Definition: KFTrajectorySmoother.h:99

Trajectory KFTrajectorySmoother::trajectory ( const Trajectory & aTraj ) const

overridevirtual

Implements TrajectorySmoother.

Definition at line 19 of file KFTrajectorySmoother.cc.

References alongMomentum, cms::cuda::assert(), heavyFlavorValidationHarvestingSequence_cff::combiner, Trajectory::direction(), remoteMonitoring_LASER_era2018_cfg::dump, Trajectory::empty(), MeasurementEstimator::estimate(), estimator(), TrajectoryStateOnSurface::globalMomentum(), TrajectoryStateOnSurface::globalPosition(), DetLayerGeometry::idToLayer(), TrajectoryStateOnSurface::isValid(), dqmiolumiharvest::j, TrajectoryStateOnSurface::localPosition(), LogDebug, LogTrace, TkCloner::makeShared(), Trajectory::measurements(), minHits_, PV3DBase< T, PVType, FrameType >::perp(), Propagator::propagate(), Propagator::propagationDirection(), TrajectoryStateOnSurface::rescaleError(), runTheMatrix::ret, Trajectory::seed(), command_line::start, theAlongPropagator, theErrorRescaling, theGeometry, theHitCloner, theOppositePropagator, UNLIKELY, TrajectoryStateUpdator::update(), and updator().

                                                                          {
   if (aTraj.empty())
     return Trajectory();
 
   const Propagator* usePropagator = theAlongPropagator;
   if (aTraj.direction() == alongMomentum) {
     usePropagator = theOppositePropagator;
   }
 
   const std::vector<TM>& avtm = aTraj.measurements();
 
 #ifdef EDM_ML_DEBUG
   LogDebug("TrackFitters") << "KFTrajectorySmoother::trajectories starting with " << avtm.size() << " HITS\n";
   for (unsigned int j = 0; j < avtm.size(); j++) {
     if (avtm[j].recHit()->det())
       LogTrace("TrackFitters") << "hit #:" << j + 1 << " rawId=" << avtm[j].recHit()->det()->geographicalId().rawId()
                                << " validity=" << avtm[j].recHit()->isValid();
     else
       LogTrace("TrackFitters") << "hit #:" << j + 1 << " Hit with no Det information";
   }
 #endif  // EDM_ML_DEBUG
 
   TrajectoryStateCombiner combiner;
   bool retry = false;
   auto start = avtm.rbegin();
 
   do {
     auto hitSize = avtm.rend() - start;
     if UNLIKELY (hitSize < minHits_) {
       LogDebug("TrackFitters") << " killing trajectory"
                                << "\n";
       return Trajectory();
     }
     Trajectory ret(aTraj.seed(), usePropagator->propagationDirection());
     Trajectory& myTraj = ret;
     myTraj.reserve(hitSize);
     retry = false;
 
     TSOS predTsos = (*start).forwardPredictedState();
     predTsos.rescaleError(theErrorRescaling);
     TSOS currTsos;
 
     auto hitCounter = hitSize;
     for (std::vector<TM>::const_reverse_iterator itm = start; itm != (avtm.rend()); ++itm, --hitCounter) {
       TransientTrackingRecHit::ConstRecHitPointer hit = itm->recHit();
 
       //check surface just for safety: should never be ==0 because they are skipped in the fitter
       // if UNLIKELY(hit->det() == nullptr) continue;
       if UNLIKELY (hit->surface() == nullptr) {
         LogDebug("TrackFitters") << " Error: invalid hit with no GeomDet attached .... skipping";
         continue;
       }
 
       if (itm != start)  //no propagation needed for first smoothed (==last fitted) hit
         predTsos = usePropagator->propagate(currTsos, *(hit->surface()));
 
       if UNLIKELY (!predTsos.isValid()) {
         LogDebug("TrackFitters") << "KFTrajectorySmoother: predicted tsos not valid!";
         LogDebug("TrackFitters") << " retry with last hit removed"
                                  << "\n";
         LogDebug("TrackFitters")
             // std::cout
             << "tsos not valid " << currTsos.globalMomentum().perp() << ' ' << hitSize << ' ' << hitCounter << ' '
             << int(hit->geographicalId()) << ' ' << hit->surface()->position().perp() << ' ' << hit->surface()->eta()
             << ' ' << hit->surface()->phi() << std::endl;
         start++;
         retry = true;
         break;
       }
 
       if (hit->isValid()) {
         TSOS combTsos, smooTsos;
 
         //3 different possibilities to calculate smoothed state:
         //1: update combined predictions with hit
         //2: combine fwd-prediction with bwd-filter
         //3: combine bwd-prediction with fwd-filter
 
         //combTsos is the predicted state with N-1 hits information. this means:
         //forward predicted state for first smoothed (last fitted) hit
         //backward predicted state for last smoothed (first fitted) hit
         //combination of forward and backward predictions for other hits
         if (itm == start)
           combTsos = itm->forwardPredictedState();
         else if (hitCounter == 1)
           combTsos = predTsos;
         else
           combTsos = combiner(predTsos, itm->forwardPredictedState());
 
         if UNLIKELY (!combTsos.isValid()) {
           LogDebug("TrackFitters") << "KFTrajectorySmoother: combined tsos not valid!\n"
                                    << "pred Tsos pos: " << predTsos.globalPosition() << "\n"
                                    << "pred Tsos mom: " << predTsos.globalMomentum() << "\n"
                                    << "TrackingRecHit: " << hit->surface()->toGlobal(hit->localPosition()) << "\n";
           start++;
           retry = true;
           break;
         }
 
         assert((hit->geographicalId() != 0U) | (!hit->canImproveWithTrack()));
         assert(hit->surface() != nullptr);
         assert((!(hit)->canImproveWithTrack()) | (nullptr != theHitCloner));
         assert((!(hit)->canImproveWithTrack()) | (nullptr != dynamic_cast<BaseTrackerRecHit const*>(hit.get())));
         auto preciseHit = theHitCloner->makeShared(hit, combTsos);
         assert(preciseHit->isValid());
         assert((preciseHit->geographicalId() != 0U) | (!preciseHit->canImproveWithTrack()));
         assert(preciseHit->surface() != nullptr);
 
         dump(*hit, hitCounter, "TrackFitters");
 
         if UNLIKELY (!preciseHit->isValid()) {
           LogTrace("TrackFitters") << "THE Precise HIT IS NOT VALID: using currTsos = predTsos"
                                    << "\n";
           currTsos = predTsos;
           myTraj.push(TM(predTsos, hit, 0, theGeometry->idToLayer(hit->geographicalId())));
         } else {
           LogTrace("TrackFitters") << "THE Precise HIT IS VALID: updating currTsos"
                                    << "\n";
 
           //update backward predicted tsos with the hit
           currTsos = updator()->update(predTsos, *preciseHit);
           if UNLIKELY (!currTsos.isValid()) {
             currTsos = predTsos;
             edm::LogWarning("KFSmoother_UpdateFailed")
                 << "Failed updating state with hit. Rolling back to non-updated state.\n"
                 << "State: " << predTsos << "Hit local pos:  " << hit->localPosition() << "\n"
                 << "Hit local err:  " << hit->localPositionError() << "\n"
                 << "Hit global pos: " << hit->globalPosition() << "\n"
                 << "Hit global err: " << hit->globalPositionError().matrix() << "\n";
           }
 
           //smooTsos updates the N-1 hits prediction with the hit
           if (itm == start)
             smooTsos = itm->updatedState();
           else if (hitCounter == 1)
             smooTsos = currTsos;
           else
             smooTsos = combiner(itm->forwardPredictedState(), currTsos);
 
           if UNLIKELY (!smooTsos.isValid()) {
             LogDebug("TrackFitters") << "KFTrajectorySmoother: smoothed tsos not valid!";
             start++;
             retry = true;
             break;
           }
 
           double estimate;
           if (itm != start)
             estimate = estimator()->estimate(combTsos, *preciseHit).second;  //correct?
           else
             estimate = itm->estimate();
 
           LogTrace("TrackFitters") << "predTsos !"
                                    << "\n"
                                    << predTsos << " with local position " << predTsos.localPosition() << "\n\n"
                                    << "currTsos !"
                                    << "\n"
                                    << currTsos << "\n"
                                    << " with local position " << currTsos.localPosition() << "\n\n"
                                    << "smooTsos !"
                                    << "\n"
                                    << smooTsos << " with local position " << smooTsos.localPosition() << "\n\n"
                                    << "smoothing estimate (with combTSOS)=" << estimate << "\n"
                                    << "filtering estimate=" << itm->estimate() << "\n";
 
           //check for valid hits with no det (refitter with constraints)
           if (preciseHit->det())
             myTraj.push(TM(itm->forwardPredictedState(),
                            predTsos,
                            smooTsos,
                            preciseHit,
                            estimate,
                            theGeometry->idToLayer(preciseHit->geographicalId())),
                         estimator()->estimate(predTsos, *preciseHit).second);
           else
             myTraj.push(TM(itm->forwardPredictedState(), predTsos, smooTsos, preciseHit, estimate),
                         estimator()->estimate(predTsos, *preciseHit).second);
           //itm->estimate());
         }
       } else {
         LogDebug("TrackFitters") << "----------------- HIT #" << hitCounter << " (INVALID)-----------------------";
 
         //no update
         currTsos = predTsos;
         TSOS combTsos;
         if (itm == start)
           combTsos = itm->forwardPredictedState();
         else if (hitCounter == 1)
           combTsos = predTsos;
         else
           combTsos = combiner(predTsos, itm->forwardPredictedState());
 
         if UNLIKELY (!combTsos.isValid()) {
           LogDebug("TrackFitters") << "KFTrajectorySmoother: combined tsos not valid!";
           return Trajectory();
         }
         assert((hit->det() == nullptr) || hit->geographicalId() != 0U);
         if (hit->det())
           myTraj.push(TM(
               itm->forwardPredictedState(), predTsos, combTsos, hit, 0, theGeometry->idToLayer(hit->geographicalId())));
         else
           myTraj.push(TM(itm->forwardPredictedState(), predTsos, combTsos, hit, 0));
       }
     }  // for loop
 
     if (!retry)
       return ret;
   } while (true);
 
   return Trajectory();
 }