TrackFitting/src/KalmanFitterFunction.cpp

0001 // This file is part of the ACTS project.
0002 //
0003 // Copyright (C) 2016 CERN for the benefit of the ACTS project
0004 //
0005 // This Source Code Form is subject to the terms of the Mozilla Public
0006 // License, v. 2.0. If a copy of the MPL was not distributed with this
0007 // file, You can obtain one at https://mozilla.org/MPL/2.0/.
0008
0009 #include "Acts/Definitions/TrackParametrization.hpp"
0010 #include "Acts/EventData/MultiTrajectory.hpp"
0011 #include "Acts/EventData/TrackContainer.hpp"
0012 #include "Acts/EventData/VectorMultiTrajectory.hpp"
0013 #include "Acts/EventData/VectorTrackContainer.hpp"
0014 #include "Acts/EventData/detail/CorrectedTransformationFreeToBound.hpp"
0015 #include "Acts/Geometry/GeometryIdentifier.hpp"
0016 #include "Acts/Propagator/DirectNavigator.hpp"
0017 #include "Acts/Propagator/Navigator.hpp"
0018 #include "Acts/Propagator/Propagator.hpp"
0019 #include "Acts/Propagator/SympyStepper.hpp"
0020 #include "Acts/TrackFitting/GainMatrixUpdater.hpp"
0021 #include "Acts/TrackFitting/KalmanFitter.hpp"
0022 #include "Acts/TrackFitting/MbfSmoother.hpp"
0023 #include "Acts/Utilities/Delegate.hpp"
0024 #include "Acts/Utilities/Logger.hpp"
0025 #include "ActsExamples/EventData/IndexSourceLink.hpp"
0026 #include "ActsExamples/EventData/MeasurementCalibration.hpp"
0027 #include "ActsExamples/EventData/Track.hpp"
0028 #include "ActsExamples/TrackFitting/RefittingCalibrator.hpp"
0029 #include "ActsExamples/TrackFitting/TrackFitterFunction.hpp"
0030
0031 #include <cmath>
0032 #include <functional>
0033 #include <memory>
0034 #include <utility>
0035 #include <vector>
0036
0037 namespace {
0038
0039 using Stepper = Acts::SympyStepper;
0040 using Propagator = Acts::Propagator<Stepper, Acts::Navigator>;
0041 using Fitter = Acts::KalmanFitter<Propagator, Acts::VectorMultiTrajectory>;
0042 using DirectPropagator = Acts::Propagator<Stepper, Acts::DirectNavigator>;
0043 using DirectFitter =
0044     Acts::KalmanFitter<DirectPropagator, Acts::VectorMultiTrajectory>;
0045
0046 using TrackContainer =
0047     Acts::TrackContainer<Acts::VectorTrackContainer,
0048                          Acts::VectorMultiTrajectory, std::shared_ptr>;
0049
0050 struct SimpleReverseFilteringLogic {
0051   double momentumThreshold = 0;
0052
0053   bool doBackwardFiltering(
0054       Acts::VectorMultiTrajectory::ConstTrackStateProxy trackState) const {
0055     auto momentum = std::abs(1 / trackState.filtered()[Acts::eBoundQOverP]);
0056     return (momentum <= momentumThreshold);
0057   }
0058 };
0059
0060 struct SimpleOutlierFinder {
0061   double chi2Cut = std::numeric_limits<double>::infinity();
0062
0063   bool isOutlier(
0064       Acts::VectorMultiTrajectory::ConstTrackStateProxy trackState) const {
0065     double chi2 = Acts::calculatePredictedChi2(trackState);
0066     return chi2 > chi2Cut;
0067   }
0068 };
0069
0070 using namespace ActsExamples;
0071
0072 struct KalmanFitterFunctionImpl final : public TrackFitterFunction {
0073   Fitter fitter;
0074   DirectFitter directFitter;
0075
0076   Acts::GainMatrixUpdater kfUpdater;
0077   Acts::MbfSmoother kfSmoother;
0078   SimpleReverseFilteringLogic reverseFilteringLogic;
0079   double reverseFilteringCovarianceScaling = 1.0;
0080   SimpleOutlierFinder outlierFinder;
0081
0082   bool multipleScattering = false;
0083   bool energyLoss = false;
0084   Acts::FreeToBoundCorrection freeToBoundCorrection;
0085
0086   IndexSourceLink::SurfaceAccessor slSurfaceAccessor;
0087
0088   KalmanFitterFunctionImpl(Fitter&& f, DirectFitter&& df,
0089                            const Acts::TrackingGeometry& trkGeo)
0090       : fitter(std::move(f)),
0091         directFitter(std::move(df)),
0092         slSurfaceAccessor{trkGeo} {}
0093
0094   template <typename calibrator_t>
0095   auto makeKfOptions(const GeneralFitterOptions& options,
0096                      const calibrator_t& calibrator) const {
0097     Acts::KalmanFitterExtensions<Acts::VectorMultiTrajectory> extensions;
0098     extensions.updater.connect<
0099         &Acts::GainMatrixUpdater::operator()<Acts::VectorMultiTrajectory>>(
0100         &kfUpdater);
0101     extensions.smoother
0102         .connect<&Acts::MbfSmoother::operator()<Acts::VectorMultiTrajectory>>(
0103             &kfSmoother);
0104     extensions.reverseFilteringLogic
0105         .connect<&SimpleReverseFilteringLogic::doBackwardFiltering>(
0106             &reverseFilteringLogic);
0107     extensions.outlierFinder.connect<&SimpleOutlierFinder::isOutlier>(
0108         &outlierFinder);
0109
0110     Acts::KalmanFitterOptions<Acts::VectorMultiTrajectory> kfOptions(
0111         options.geoContext, options.magFieldContext, options.calibrationContext,
0112         extensions, options.propOptions, &(*options.referenceSurface));
0113
0114     kfOptions.referenceSurfaceStrategy =
0115         Acts::KalmanFitterTargetSurfaceStrategy::first;
0116     kfOptions.multipleScattering = multipleScattering;
0117     kfOptions.energyLoss = energyLoss;
0118     kfOptions.freeToBoundCorrection = freeToBoundCorrection;
0119     kfOptions.extensions.calibrator.connect<&calibrator_t::calibrate>(
0120         &calibrator);
0121     kfOptions.reversedFilteringCovarianceScaling =
0122         reverseFilteringCovarianceScaling;
0123
0124     if (options.doRefit) {
0125       kfOptions.extensions.surfaceAccessor
0126           .connect<&RefittingCalibrator::accessSurface>();
0127     } else {
0128       kfOptions.extensions.surfaceAccessor
0129           .connect<&IndexSourceLink::SurfaceAccessor::operator()>(
0130               &slSurfaceAccessor);
0131     }
0132
0133     return kfOptions;
0134   }
0135
0136   TrackFitterResult operator()(const std::vector<Acts::SourceLink>& sourceLinks,
0137                                const TrackParameters& initialParameters,
0138                                const GeneralFitterOptions& options,
0139                                const MeasurementCalibratorAdapter& calibrator,
0140                                TrackContainer& tracks) const override {
0141     const auto kfOptions = makeKfOptions(options, calibrator);
0142     return fitter.fit(sourceLinks.begin(), sourceLinks.end(), initialParameters,
0143                       kfOptions, tracks);
0144   }
0145
0146   TrackFitterResult operator()(
0147       const std::vector<Acts::SourceLink>& sourceLinks,
0148       const TrackParameters& initialParameters,
0149       const GeneralFitterOptions& options,
0150       const RefittingCalibrator& calibrator,
0151       const std::vector<const Acts::Surface*>& surfaceSequence,
0152       TrackContainer& tracks) const override {
0153     const auto kfOptions = makeKfOptions(options, calibrator);
0154     return directFitter.fit(sourceLinks.begin(), sourceLinks.end(),
0155                             initialParameters, kfOptions, surfaceSequence,
0156                             tracks);
0157   }
0158 };
0159
0160 }  // namespace
0161
0162 std::shared_ptr<ActsExamples::TrackFitterFunction>
0163 ActsExamples::makeKalmanFitterFunction(
0164     std::shared_ptr<const Acts::TrackingGeometry> trackingGeometry,
0165     std::shared_ptr<const Acts::MagneticFieldProvider> magneticField,
0166     bool multipleScattering, bool energyLoss,
0167     double reverseFilteringMomThreshold,
0168     double reverseFilteringCovarianceScaling,
0169     Acts::FreeToBoundCorrection freeToBoundCorrection, double chi2Cut,
0170     const Acts::Logger& logger) {
0171   // Stepper should be copied into the fitters
0172   const Stepper stepper(std::move(magneticField));
0173
0174   // Standard fitter
0175   const auto& geo = *trackingGeometry;
0176   Acts::Navigator::Config cfg{std::move(trackingGeometry)};
0177   cfg.resolvePassive = false;
0178   cfg.resolveMaterial = true;
0179   cfg.resolveSensitive = true;
0180   Acts::Navigator navigator(cfg, logger.cloneWithSuffix("Navigator"));
0181   Propagator propagator(stepper, std::move(navigator),
0182                         logger.cloneWithSuffix("Propagator"));
0183   Fitter trackFitter(std::move(propagator), logger.cloneWithSuffix("Fitter"));
0184
0185   // Direct fitter
0186   Acts::DirectNavigator directNavigator{
0187       logger.cloneWithSuffix("DirectNavigator")};
0188   DirectPropagator directPropagator(stepper, std::move(directNavigator),
0189                                     logger.cloneWithSuffix("DirectPropagator"));
0190   DirectFitter directTrackFitter(std::move(directPropagator),
0191                                  logger.cloneWithSuffix("DirectFitter"));
0192
0193   // build the fitter function. owns the fitter object.
0194   auto fitterFunction = std::make_shared<KalmanFitterFunctionImpl>(
0195       std::move(trackFitter), std::move(directTrackFitter), geo);
0196   fitterFunction->multipleScattering = multipleScattering;
0197   fitterFunction->energyLoss = energyLoss;
0198   fitterFunction->reverseFilteringLogic.momentumThreshold =
0199       reverseFilteringMomThreshold;
0200   fitterFunction->freeToBoundCorrection = freeToBoundCorrection;
0201   fitterFunction->reverseFilteringCovarianceScaling =
0202       reverseFilteringCovarianceScaling;
0203   fitterFunction->outlierFinder.chi2Cut = chi2Cut;
0204
0205   return fitterFunction;
0206 }