8 #include "Acts/EventData/Types.hpp"
9 #include "GaudiKernel/TypeNameString.h"
14 #include "GaudiKernel/EventContext.h"
20 #include "Acts/EventData/Types.hpp"
21 #include "Acts/Definitions/TrackParametrization.hpp"
22 #include "Acts/Definitions/Units.hpp"
23 #include "Acts/Propagator/SympyStepper.hpp"
24 #include "Acts/Propagator/Navigator.hpp"
25 #include "Acts/Propagator/Propagator.hpp"
26 #include "Acts/Surfaces/PerigeeSurface.hpp"
27 #include "Acts/Surfaces/Surface.hpp"
28 #include "Acts/TrackFitting/KalmanFitter.hpp"
29 #include "Acts/Utilities/Helpers.hpp"
30 #include "Acts/Utilities/Logger.hpp"
31 #include "Acts/Utilities/CalibrationContext.hpp"
32 #include "Acts/EventData/VectorTrackContainer.hpp"
48 #include "Acts/Propagator/DirectNavigator.hpp"
60 template <
typename trajectory_t>
62 const Acts::CalibrationContext& ,
63 const Acts::SourceLink& sl,
64 typename trajectory_t::TrackStateProxy trackState)
const {
67 trackState.setUncalibratedSourceLink(Acts::SourceLink{sl});
69 const Acts::BoundTrackParameters actsParam(trackState.referenceSurface().getSharedPtr(),
70 trackState.predicted(),
71 trackState.predictedCovariance(),
91 if (slsurf !=
nullptr) {
105 }
else if (plsurf !=
nullptr) {
106 if ((trkParam.get())->covariance() !=
nullptr) {
128 int dim = (*rot).localParameters().dimension();
129 trackState.allocateCalibrated(
dim);
134 throw std::runtime_error(
"Cannot create dim 0 measurement");
135 }
else if (
dim == 1) {
136 trackState.template calibrated<1>() = (*rot).localParameters().template head<1>();
137 trackState.template calibratedCovariance<1>() = (*rot).localCovariance().template topLeftCorner<1, 1>();
138 Acts::BoundSubspaceIndices subspaceIndices;
140 subspaceIndices = {Acts::eBoundLoc1};
142 subspaceIndices = {Acts::eBoundLoc0};
144 trackState.setProjectorSubspaceIndices(subspaceIndices);
148 trackState.template calibrated<2>() = (*rot).localParameters().template head<2>();
149 trackState.template calibratedCovariance<2>() = (*rot).localCovariance().template topLeftCorner<2, 2>();
150 Acts::BoundSubspaceIndices subspaceIndices = {Acts::eBoundLoc0, Acts::eBoundLoc1};
151 trackState.setProjectorSubspaceIndices(subspaceIndices);
156 " currently not supported.");
170 const std::string&
n,
171 const IInterface*
p) :
189 auto field = std::make_shared<ATLASMagneticFieldWrapper>();
192 Acts::SympyStepper stepper(
field);
194 logger().cloneWithSuffix(
"Navigator"));
195 Acts::Propagator<Acts::SympyStepper, Acts::Navigator> propagator(stepper,
196 std::move(navigator),
197 logger().cloneWithSuffix(
"Prop"));
199 m_fitter = std::make_unique<Fitter>(std::move(propagator),
200 logger().cloneWithSuffix(
"KalmanFitter"));
203 Acts::DirectNavigator directNavigator(
logger().cloneWithSuffix(
"DirectNavigator") );
204 Acts::Propagator<Acts::SympyStepper, Acts::DirectNavigator> directPropagator(std::move(stepper),
205 std::move(directNavigator),
206 logger().cloneWithSuffix(
"DirectPropagator"));
208 m_directFitter = std::make_unique<DirectFitter>(std::move(directPropagator),
209 logger().cloneWithSuffix(
"DirectKalmanFitter"));
218 m_kfExtensions.updater.connect<&ActsTrk::detail::FitterHelperFunctions::gainMatrixUpdate<ActsTrk::MutableTrackStateBackend>>();
219 m_kfExtensions.smoother.connect<&ActsTrk::detail::FitterHelperFunctions::mbfSmoother<ActsTrk::MutableTrackStateBackend>>();
220 m_kfExtensions.calibrator.connect<&ActsTrk::detail::TrkMeasurementCalibrator::calibrate<ActsTrk::MutableTrackStateBackend>>(
m_calibrator.get());
222 return StatusCode::SUCCESS;
227 std::unique_ptr<Trk::Track>
233 std::unique_ptr<Trk::Track>
track =
nullptr;
234 ATH_MSG_VERBOSE (
"--> enter KalmanFitter::fit(Track,,) with Track from author = "
239 ATH_MSG_DEBUG(
"called to refit empty track or track with too little information, reject fit");
245 ATH_MSG_DEBUG(
"input fails to provide track parameters for seeding the KF, reject fit");
250 auto pSurface = Acts::Surface::makeShared<Acts::PerigeeSurface>(
251 Acts::Vector3{0., 0., 0.});
256 Acts::CalibrationContext calContext = Acts::CalibrationContext();
258 Acts::KalmanFitterExtensions<ActsTrk::MutableTrackStateBackend> kfExtensions =
m_kfExtensions;
261 kfExtensions.surfaceAccessor.connect<&ATLASSourceLinkSurfaceAccessor::operator()>(&surfaceAccessor);
263 Acts::PropagatorPlainOptions propagationOption(tgContext, mfContext);
266 Acts::KalmanFitterOptions
267 kfOptions(tgContext, mfContext, calContext,
272 std::vector<Acts::SourceLink> trackSourceLinks =
m_ATLASConverterTool->trkTrackToSourceLinks(tgContext,inputTrack);
274 if (trackSourceLinks.empty()) {
275 ATH_MSG_DEBUG(
"input contain measurement but no source link created, probable issue with the converter, reject fit ");
283 Acts::BoundSquareMatrix scaledCov = Acts::BoundSquareMatrix::Identity();
284 for (
int i=0;
i<6; ++
i) {
286 (scaledCov)(
i,
i) =
scale * initialParams.covariance().value()(
i,
i);
292 const Acts::BoundTrackParameters scaledInitialParams(initialParams.referenceSurface().getSharedPtr(),
293 initialParams.parameters(),
300 auto result =
m_fitter->fit(trackSourceLinks.begin(), trackSourceLinks.end(),
301 scaledInitialParams, kfOptions, tracks);
310 std::unique_ptr<Trk::Track>
317 std::unique_ptr<Trk::Track>
track =
nullptr;
320 if (inputMeasSet.size() < 2) {
321 ATH_MSG_DEBUG(
"called to refit empty measurement set or a measurement set with too little information, reject fit");
326 auto pSurface = Acts::Surface::makeShared<Acts::PerigeeSurface>(
327 Acts::Vector3{0., 0., 0.});
332 Acts::CalibrationContext calContext = Acts::CalibrationContext();
334 Acts::KalmanFitterExtensions<ActsTrk::MutableTrackStateBackend> kfExtensions =
m_kfExtensions;
337 kfExtensions.surfaceAccessor.connect<&ATLASSourceLinkSurfaceAccessor::operator()>(&surfaceAccessor);
339 Acts::PropagatorPlainOptions propagationOption(tgContext, mfContext);
342 Acts::KalmanFitterOptions
343 kfOptions(tgContext, mfContext, calContext,
348 std::vector<Acts::SourceLink> trackSourceLinks;
349 trackSourceLinks.reserve(inputMeasSet.size());
351 for (
auto it = inputMeasSet.begin();
it != inputMeasSet.end(); ++
it){
355 if (trackSourceLinks.empty()) {
356 ATH_MSG_DEBUG(
"input contain measurement but no source link created, probable issue with the converter, reject fit ");
360 const auto& initialParams =
m_ATLASConverterTool->trkTrackParametersToActsParameters(estimatedStartParameters, tgContext);
365 auto result =
m_fitter->fit(trackSourceLinks.begin(), trackSourceLinks.end(),
366 initialParams, kfOptions, tracks);
375 std::unique_ptr<Trk::Track>
385 std::unique_ptr<Trk::Track>
track =
nullptr;
388 auto pSurface = Acts::Surface::makeShared<Acts::PerigeeSurface>(
389 Acts::Vector3{0., 0., 0.});
394 Acts::CalibrationContext calContext = Acts::CalibrationContext();
396 Acts::KalmanFitterExtensions<ActsTrk::MutableTrackStateBackend> kfExtensions =
m_kfExtensions;
402 kfExtensions.calibrator.connect<&PRDSourceLinkCalibrator::calibrate<ActsTrk::MutableTrackStateBackend>>(&calibrator);
405 kfExtensions.surfaceAccessor.connect<&PRDSourceLinkSurfaceAccessor::operator()>(&surfaceAccessor);
407 Acts::PropagatorPlainOptions propagationOption(tgContext, mfContext);
410 Acts::KalmanFitterOptions
411 kfOptions(tgContext, mfContext, calContext,
417 std::vector<Acts::SourceLink> trackSourceLinks;
418 trackSourceLinks.reserve(inputPRDColl.size());
421 trackSourceLinks.push_back(Acts::SourceLink{
PRDSourceLink{prd}});
424 if (trackSourceLinks.empty()) {
425 ATH_MSG_WARNING(
"input contain measurement but no source link created, probable issue with the converter, reject fit ");
430 const auto& initialParams =
m_ATLASConverterTool->trkTrackParametersToActsParameters(estimatedStartParameters, tgContext);
435 auto result =
m_fitter->fit(trackSourceLinks.begin(), trackSourceLinks.end(),
436 initialParams, kfOptions, tracks);
445 std::unique_ptr< ActsTrk::MutableTrackContainer >
447 const std::vector< ActsTrk::ATLASUncalibSourceLink> & clusterList,
448 const Acts::BoundTrackParameters& initialParams,
449 const Acts::GeometryContext& tgContext,
450 const Acts::MagneticFieldContext& mfContext,
451 const Acts::CalibrationContext& calContext,
453 const Acts::Surface* targetSurface)
const{
454 ATH_MSG_DEBUG(
"--> entering KalmanFitter::fit(xAODMeasure...things,TP,)");
456 std::vector<Acts::SourceLink> sourceLinks;
457 sourceLinks.reserve(clusterList.size());
459 std::vector<const Acts::Surface*> surfaces;
460 surfaces.reserve(clusterList.size());
462 const Acts::TrackingGeometry *
464 if (!actsTrackingGeometry) {
465 throw std::runtime_error(
"No Acts tracking geometry.");
469 sourceLinks.emplace_back(
el );
473 Acts::KalmanFitterExtensions<ActsTrk::MutableTrackStateBackend> kfExtensions =
m_kfExtensions;
476 kfExtensions.surfaceAccessor.connect<&ActsTrk::ATLASUncalibSourceLinkSurfaceAccessor::operator()>(&surfaceAccessor);
479 ::NoCalibration(*actsTrackingGeometry, detectorElementToGeometryIdMap);
483 Acts::PropagatorPlainOptions propagationOption(tgContext, mfContext);
487 std::shared_ptr<Acts::Surface> pSurface{
nullptr};
489 pSurface = Acts::Surface::makeShared<Acts::PerigeeSurface>(Acts::Vector3{0., 0., 0.});
490 targetSurface = pSurface.get();
494 Acts::KalmanFitterOptions<ActsTrk::MutableTrackStateBackend>
495 kfOptions(tgContext, mfContext, calContext,
500 std::unique_ptr< ActsTrk::MutableTrackContainer > tracks = std::make_unique< ActsTrk::MutableTrackContainer >();
521 std::unique_ptr<Trk::Track>
532 if (addMeasColl.empty()) {
533 ATH_MSG_DEBUG(
"client tries to add an empty MeasurementSet to the track fit." );
534 return fit(ctx,inputTrack);
539 ATH_MSG_DEBUG(
"called to refit empty track or track with too little information, reject fit");
545 ATH_MSG_DEBUG(
"input fails to provide track parameters for seeding the KF, reject fit");
549 std::unique_ptr<Trk::Track>
track =
nullptr;
552 auto pSurface = Acts::Surface::makeShared<Acts::PerigeeSurface>(
553 Acts::Vector3{0., 0., 0.});
558 Acts::CalibrationContext calContext = Acts::CalibrationContext();
560 Acts::KalmanFitterExtensions<ActsTrk::MutableTrackStateBackend> kfExtensions =
m_kfExtensions;
562 Acts::PropagatorPlainOptions propagationOption(tgContext, mfContext);
565 Acts::KalmanFitterOptions
566 kfOptions(tgContext, mfContext, calContext,
571 std::vector<Acts::SourceLink> trackSourceLinks =
m_ATLASConverterTool->trkTrackToSourceLinks(tgContext, inputTrack);
574 for (
auto it = addMeasColl.begin();
it != addMeasColl.end(); ++
it)
579 if (trackSourceLinks.empty()) {
580 ATH_MSG_DEBUG(
"input contain measurement but no source link created, probable issue with the converter, reject fit ");
586 auto result =
m_fitter->fit(trackSourceLinks.begin(), trackSourceLinks.end(),
587 initialParams, kfOptions, tracks);
596 std::unique_ptr<Trk::Track>
603 ATH_MSG_DEBUG(
"Fit of Track with additional PrepRawDataSet not yet implemented");
609 std::unique_ptr<Trk::Track>
622 ATH_MSG_DEBUG(
"input #2 is empty try to fit track 1 alone" );
623 return fit(ctx,intrk1);
628 ATH_MSG_DEBUG(
"input #1 is empty try to fit track 2 alone" );
629 return fit(ctx,intrk2);
634 ATH_MSG_DEBUG(
"input #1 fails to provide track parameters for seeding the KF, reject fit");
638 std::unique_ptr<Trk::Track>
track =
nullptr;
641 auto pSurface = Acts::Surface::makeShared<Acts::PerigeeSurface>(
642 Acts::Vector3{0., 0., 0.});
647 Acts::CalibrationContext calContext = Acts::CalibrationContext();
649 Acts::KalmanFitterExtensions<ActsTrk::MutableTrackStateBackend> kfExtensions =
m_kfExtensions;
651 Acts::PropagatorPlainOptions propagationOption(tgContext, mfContext);
654 Acts::KalmanFitterOptions
655 kfOptions(tgContext, mfContext, calContext,
660 std::vector<Acts::SourceLink> trackSourceLinks =
m_ATLASConverterTool->trkTrackToSourceLinks(tgContext, intrk1);
661 std::vector<Acts::SourceLink> trackSourceLinks2 =
m_ATLASConverterTool->trkTrackToSourceLinks(tgContext, intrk2);
662 trackSourceLinks.insert(trackSourceLinks.end(), trackSourceLinks2.begin(), trackSourceLinks2.end());
664 if (trackSourceLinks.empty()) {
665 ATH_MSG_DEBUG(
"input contain measurement but no source link created, probable issue with the converter, reject fit ");
673 Acts::BoundSquareMatrix scaledCov = Acts::BoundSquareMatrix::Identity();
674 for (
int i=0;
i<6; ++
i) {
676 (scaledCov)(
i,
i) =
scale * initialParams.covariance().value()(
i,
i);
679 const Acts::BoundTrackParameters scaledInitialParams(initialParams.referenceSurface().getSharedPtr(),
680 initialParams.parameters(),
687 auto result =
m_fitter->fit(trackSourceLinks.begin(), trackSourceLinks.end(),
688 scaledInitialParams, kfOptions, tracks);
695 std::unique_ptr<Trk::Track>
697 Acts::GeometryContext& tgContext,
699 Acts::Result<ActsTrk::MutableTrackContainer::TrackProxy, std::error_code>& fitResult,
700 bool SourceLinkType)
const {
702 if (not fitResult.ok())
705 std::unique_ptr<Trk::Track> newtrack =
nullptr;
707 const auto& acts_track = fitResult.value();
708 auto finalTrajectory = std::make_unique<Trk::TrackStates>();
710 int numberOfDeadPixel = 0;
711 int numberOfDeadSCT = 0;
713 std::vector<std::unique_ptr<const Acts::BoundTrackParameters>> actsSmoothedParam;
715 tracks.trackStateContainer().visitBackwards(acts_track.tipIndex(),
716 [&] (
const auto &state) ->
void
719 auto flag = state.typeFlags();
720 const auto* associatedDetEl = state.referenceSurface().associatedDetectorElement();
721 if (not associatedDetEl)
724 const auto* actsElement = dynamic_cast<const ActsDetectorElement*>(associatedDetEl);
728 const auto* upstreamDetEl = actsElement->upstreamDetectorElement();
729 if (not upstreamDetEl)
732 ATH_MSG_VERBOSE(
"Try casting to TRT for if");
733 if (dynamic_cast<const InDetDD::TRT_BaseElement*>(upstreamDetEl))
736 const auto* trkDetElem = dynamic_cast<const Trk::TrkDetElementBase*>(upstreamDetEl);
740 ATH_MSG_VERBOSE(
"trkDetElem type: " << static_cast<std::underlying_type_t<Trk::DetectorElemType>>(trkDetElem->detectorType()));
742 ATH_MSG_VERBOSE(
"Try casting to SiDetectorElement");
743 const auto* detElem = dynamic_cast<const InDetDD::SiDetectorElement*>(upstreamDetEl);
746 ATH_MSG_VERBOSE(
"detElem = " << detElem);
749 std::bitset<Trk::TrackStateOnSurface::NumberOfTrackStateOnSurfaceTypes> typePattern;
750 std::unique_ptr<Trk::TrackParameters> parm;
753 if (flag.test(Acts::TrackStateFlag::HoleFlag)){
754 ATH_MSG_VERBOSE(
"State is a hole (no associated measurement), use predicted parameters");
755 const Acts::BoundTrackParameters actsParam(state.referenceSurface().getSharedPtr(),
757 state.predictedCovariance(),
758 acts_track.particleHypothesis());
759 parm = m_ATLASConverterTool->actsTrackParametersToTrkParameters(actsParam, tgContext);
760 auto boundaryCheck = m_boundaryCheckTool->boundaryCheck(*parm);
763 ATH_MSG_VERBOSE(
"Check if this is a hole, a dead sensors or a state outside the sensor boundary");
764 if(boundaryCheck == Trk::BoundaryCheckResult::DeadElement){
765 if (detElem->isPixel()) {
768 else if (detElem->isSCT()) {
773 } else if (boundaryCheck != Trk::BoundaryCheckResult::Candidate){
777 typePattern.set(Trk::TrackStateOnSurface::Hole);
780 else if (
flag.test(Acts::TrackStateFlag::OutlierFlag) or !state.hasSmoothed()) {
781 ATH_MSG_VERBOSE(
"The state was tagged as an outlier or was missed in the reverse filtering, use filtered parameters");
782 const Acts::BoundTrackParameters actsParam(state.referenceSurface().getSharedPtr(),
784 state.filteredCovariance(),
785 acts_track.particleHypothesis());
786 parm = m_ATLASConverterTool->actsTrackParametersToTrkParameters(actsParam, tgContext);
787 typePattern.set(Trk::TrackStateOnSurface::Outlier);
791 ATH_MSG_VERBOSE(
"The state is a measurement state, use smoothed parameters");
793 const Acts::BoundTrackParameters actsParam(state.referenceSurface().getSharedPtr(),
795 state.smoothedCovariance(),
796 acts_track.particleHypothesis());
798 actsSmoothedParam.push_back(std::make_unique<const Acts::BoundTrackParameters>(Acts::BoundTrackParameters(actsParam)));
799 parm = m_ATLASConverterTool->actsTrackParametersToTrkParameters(actsParam, tgContext);
800 typePattern.set(Trk::TrackStateOnSurface::Measurement);
802 std::unique_ptr<Trk::MeasurementBase> measState;
803 if (state.hasUncalibratedSourceLink() && !SourceLinkType){
804 auto sl = state.getUncalibratedSourceLink().template get<ATLASSourceLink>();
806 measState = sl->uniqueClone();
808 else if (state.hasUncalibratedSourceLink() && SourceLinkType){
809 auto sl = state.getUncalibratedSourceLink().template get<PRDSourceLink>().prd;
812 const IdentifierHash idHash = sl->detectorElement()->identifyHash();
813 int dim = state.calibratedSize();
814 std::unique_ptr<Trk::RIO_OnTrack> rot;
816 const InDet::SCT_Cluster* sct_Cluster = dynamic_cast<const InDet::SCT_Cluster*>(sl);
818 ATH_MSG_ERROR(
"ERROR could not cast PRD to SCT_Cluster");
821 rot = std::make_unique<InDet::SCT_ClusterOnTrack>(sct_Cluster,Trk::LocalParameters(Trk::DefinedParameter(state.template calibrated<1>()[0], Trk::loc1)), state.template calibratedCovariance<1>(),idHash);
827 ATH_MSG_VERBOSE(
"Dimension is 2 but we need SCT_Cluster for this measurment");
832 rot = std::make_unique<InDet::PixelClusterOnTrack>(
pixelCluster,
Trk::LocalParameters(state.template calibrated<2>()),state.template calibratedCovariance<2>(),idHash);
836 throw std::domain_error(
"Cannot handle measurement dim>2");
838 measState = rot->uniqueClone();
840 double nDoF = state.calibratedSize();
845 ATH_MSG_VERBOSE(
"State succesfully creates, adding it to the trajectory");
846 finalTrajectory->insert(finalTrajectory->begin(), perState);
850 const Acts::BoundTrackParameters actsPer(acts_track.referenceSurface().getSharedPtr(),
851 acts_track.parameters(),
852 acts_track.covariance(),
853 acts_track.particleHypothesis());
854 std::unique_ptr<Trk::TrackParameters> per = m_ATLASConverterTool->actsTrackParametersToTrkParameters(actsPer, tgContext);
855 std::bitset<Trk::TrackStateOnSurface::NumberOfTrackStateOnSurfaceTypes> typePattern;
858 if (perState) finalTrajectory->insert(finalTrajectory->begin(), perState);
863 newtrack = std::make_unique<Trk::Track>(newInfo, std::move(finalTrajectory),
nullptr);
866 if (!newtrack->trackSummary()) {
867 newtrack->setTrackSummary(std::make_unique<Trk::TrackSummary>());
874 m_trkSummaryTool->updateTrackSummary(ctx, *newtrack,
true);
879 std::unique_ptr< ActsTrk::MutableTrackContainer >
882 const Acts::BoundTrackParameters& initialParams,
883 const Acts::GeometryContext& tgContext,
884 const Acts::MagneticFieldContext& mfContext,
885 const Acts::CalibrationContext& calContext,
888 const Acts::TrackingGeometry *
889 actsTrackingGeometry = m_trackingGeometryTool->trackingGeometry().get();
890 if (!actsTrackingGeometry) {
891 throw std::runtime_error(
"No Acts tracking geometry.");
894 std::vector<ActsTrk::ATLASUncalibSourceLink> sourceLinks;
895 sourceLinks.reserve(6);
897 std::vector<const Acts::Surface*> surfaces;
900 const auto& sps = seed.sp();
902 const auto& measurements = sp->measurements();
905 sourceLinks.emplace_back(
el );
909 return fit(ctx, sourceLinks, initialParams, tgContext, mfContext, calContext, detectorElementToGeometryIdMap, surfaces.front());