From 57ad820f21be4a26595790ca343182276c680a9c Mon Sep 17 00:00:00 2001 From: AJPfleger Date: Wed, 30 Nov 2022 10:54:57 +0100 Subject: [PATCH 1/9] fix: update --- Core/include/Acts/TrackFitting/Chi2Fitter.hpp | 79 +++++++++++-------- 1 file changed, 48 insertions(+), 31 deletions(-) diff --git a/Core/include/Acts/TrackFitting/Chi2Fitter.hpp b/Core/include/Acts/TrackFitting/Chi2Fitter.hpp index bd624a2c835..e91c1f38876 100644 --- a/Core/include/Acts/TrackFitting/Chi2Fitter.hpp +++ b/Core/include/Acts/TrackFitting/Chi2Fitter.hpp @@ -245,6 +245,7 @@ class Chi2Fitter { /// Whether to consider energy loss. bool energyLoss = false; // TODO: add later + int updateNumber; /// Whether to include non-linear correction during global to local /// transformation FreeToBoundCorrection freeToBoundCorrection; @@ -321,10 +322,6 @@ class Chi2Fitter { // We need the full jacobianFromStart, so we'll need to calculate it no // matter if we have a measurement or not. - // Transport the covariance to the surface - stepper.transportCovarianceToBound(state.stepping, *surface, - freeToBoundCorrection); - // Update state and stepper with pre material effects materialInteractor(surface, state, stepper, MaterialUpdateStage::PreUpdate); @@ -333,7 +330,7 @@ class Chi2Fitter { // is called with fullUpdate *after* retrieving the boundState. // Bind the transported state to the current surface - auto res = stepper.boundState(state.stepping, *surface, false); + auto res = stepper.boundState(state.stepping, *surface, true); if (!res.ok()) { return res.error(); } @@ -352,13 +349,27 @@ class Chi2Fitter { // add a full TrackState entry multi trajectory // (this allocates storage for all components, we will set them later) - result.lastTrackIndex = result.fittedStates->addTrackState( - ~(TrackStatePropMask::Smoothed | TrackStatePropMask::Filtered), - result.lastTrackIndex); + + size_t index; + if (updateNumber == 0) { + result.lastTrackIndex = + result.fittedStates + ->addTrackState( // which lastTrackIndex is this? + ~(TrackStatePropMask::Smoothed | + TrackStatePropMask::Filtered), + result.lastTrackIndex); + index = result.lastTrackIndex; + } else { + result.fittedStates->visitBackwards( + result.lastTrackIndex, [&](auto proxy) { + if (&proxy.referenceSurface() == surface) { + index = proxy.index(); + } + }); + } // now get track state proxy back - auto trackStateProxy = - result.fittedStates->getTrackState(result.lastTrackIndex); + auto trackStateProxy = result.fittedStates->getTrackState(index); trackStateProxy.setReferenceSurface(surface->getSharedPtr()); @@ -390,10 +401,9 @@ class Chi2Fitter { trackStateProxy .template calibrated(); // 2x1 or 1x1 - const auto& covariance = - trackStateProxy.template calibratedCovariance< - kMeasurementSize>(); // 2x2 or 1x1. Should - // be diagonal. + const auto covariance = trackStateProxy.template calibratedCovariance< + kMeasurementSize>(); // 2x2 or 1x1. Should + // be diagonal. const auto covInv = covariance.inverse(); auto residuals = @@ -597,11 +607,12 @@ class Chi2Fitter { /// the fit. /// /// @return the output as an output track - template Result> fit( source_link_iterator_t it, source_link_iterator_t end, - const start_parameters_t& sParameters, + const BoundTrackParameters& sParameters, + // const start_parameters_t& sParameters, const Chi2FitterOptions& chi2FitterOptions, std::shared_ptr trajectory = {}) const { const auto& logger = chi2FitterOptions.logger; @@ -632,6 +643,18 @@ class Chi2Fitter { // the result object which will be returned. Overridden every iteration. Chi2Result c2r; + if (not trajectory) { + trajectory = std::make_shared(); + }; + + using paramType = typename std::conditional< + std::is_same::value, + std::variant, + std::variant>::type; + + paramType vParams = sParameters; + auto updatedStartParameters = sParameters; + for (int i = 0; i <= chi2FitterOptions.nUpdates; ++i) { // Create relevant options for the propagation options PropagatorOptions propOptions( @@ -648,32 +671,24 @@ class Chi2Fitter { chi2Actor.energyLoss = chi2FitterOptions.energyLoss; chi2Actor.freeToBoundCorrection = chi2FitterOptions.freeToBoundCorrection; chi2Actor.extensions = chi2FitterOptions.extensions; + chi2Actor.updateNumber = i; typename propagator_t::template action_list_t_result_t< CurvilinearTrackParameters, Actors> inputResult; auto& r = inputResult.template get>(); - - if (trajectory) { - r.fittedStates = trajectory; - } else { - r.fittedStates = std::make_shared(); + r.fittedStates = trajectory; + if (i > 0) { + r.lastTrackIndex = c2r.lastTrackIndex; } - using paramType = typename std::conditional< - std::is_same::value, - std::variant, - std::variant>::type; - paramType vParams = sParameters; // start_parameters_t and parameter_t can be the same - - auto result = m_propagator.template propagate(sParameters, propOptions, - std::move(inputResult)); - + auto result = m_propagator.template propagate( + updatedStartParameters, propOptions, std::move(inputResult)); if (!result.ok()) { ACTS_ERROR("chi2 | it=" << i - << " | propapation failed: " << result.error()); + << " | propagation failed: " << result.error()); return result.error(); } @@ -706,6 +721,7 @@ class Chi2Fitter { c2rCurrent.covariance.inverse() * c2rCurrent.residuals; ACTS_VERBOSE("chi2 | it=" << i << " | χ² = " << c2rCurrent.chisquare); + // chisquare is here the global chisquare // copy over data from previous runs (namely chisquares vector) c2rCurrent.chisquares.reserve(c2r.chisquares.size() + 1); @@ -751,6 +767,7 @@ class Chi2Fitter { vParams); vParams = c2r.fittedParameters.value(); // passed to next iteration + updatedStartParameters = c2r.fittedParameters.value(); } // Return the converted track From 2f793fbc9f18c09abe34417814ab2f31f01f8574 Mon Sep 17 00:00:00 2001 From: AJPfleger Date: Wed, 30 Nov 2022 11:04:29 +0100 Subject: [PATCH 2/9] improve variable name currentTrackIndex --- Core/include/Acts/TrackFitting/Chi2Fitter.hpp | 9 +++++---- 1 file changed, 5 insertions(+), 4 deletions(-) diff --git a/Core/include/Acts/TrackFitting/Chi2Fitter.hpp b/Core/include/Acts/TrackFitting/Chi2Fitter.hpp index e91c1f38876..5b08e68d01a 100644 --- a/Core/include/Acts/TrackFitting/Chi2Fitter.hpp +++ b/Core/include/Acts/TrackFitting/Chi2Fitter.hpp @@ -350,7 +350,7 @@ class Chi2Fitter { // add a full TrackState entry multi trajectory // (this allocates storage for all components, we will set them later) - size_t index; + size_t currentTrackIndex; if (updateNumber == 0) { result.lastTrackIndex = result.fittedStates @@ -358,18 +358,19 @@ class Chi2Fitter { ~(TrackStatePropMask::Smoothed | TrackStatePropMask::Filtered), result.lastTrackIndex); - index = result.lastTrackIndex; + currentTrackIndex = result.lastTrackIndex; } else { result.fittedStates->visitBackwards( result.lastTrackIndex, [&](auto proxy) { if (&proxy.referenceSurface() == surface) { - index = proxy.index(); + currentTrackIndex = proxy.index(); } }); } // now get track state proxy back - auto trackStateProxy = result.fittedStates->getTrackState(index); + auto trackStateProxy = + result.fittedStates->getTrackState(currentTrackIndex); trackStateProxy.setReferenceSurface(surface->getSharedPtr()); From 86656fa8e5ac001e24683a5078f7709d23a1c860 Mon Sep 17 00:00:00 2001 From: AJPfleger Date: Wed, 30 Nov 2022 13:08:47 +0100 Subject: [PATCH 3/9] clang-tidy --- Core/include/Acts/TrackFitting/Chi2Fitter.hpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/Core/include/Acts/TrackFitting/Chi2Fitter.hpp b/Core/include/Acts/TrackFitting/Chi2Fitter.hpp index 5b08e68d01a..00316078267 100644 --- a/Core/include/Acts/TrackFitting/Chi2Fitter.hpp +++ b/Core/include/Acts/TrackFitting/Chi2Fitter.hpp @@ -245,7 +245,7 @@ class Chi2Fitter { /// Whether to consider energy loss. bool energyLoss = false; // TODO: add later - int updateNumber; + int updateNumber = SIZE_MAX; /// Whether to include non-linear correction during global to local /// transformation FreeToBoundCorrection freeToBoundCorrection; @@ -350,7 +350,7 @@ class Chi2Fitter { // add a full TrackState entry multi trajectory // (this allocates storage for all components, we will set them later) - size_t currentTrackIndex; + size_t currentTrackIndex = SIZE_MAX; if (updateNumber == 0) { result.lastTrackIndex = result.fittedStates From 4c6d622f126db3e90c0cdd4336ac01aef01d3f21 Mon Sep 17 00:00:00 2001 From: AJPfleger Date: Wed, 30 Nov 2022 13:32:02 +0100 Subject: [PATCH 4/9] clang-tidy --- Core/include/Acts/TrackFitting/Chi2Fitter.hpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Core/include/Acts/TrackFitting/Chi2Fitter.hpp b/Core/include/Acts/TrackFitting/Chi2Fitter.hpp index 00316078267..bf93e7e429c 100644 --- a/Core/include/Acts/TrackFitting/Chi2Fitter.hpp +++ b/Core/include/Acts/TrackFitting/Chi2Fitter.hpp @@ -245,7 +245,7 @@ class Chi2Fitter { /// Whether to consider energy loss. bool energyLoss = false; // TODO: add later - int updateNumber = SIZE_MAX; + int updateNumber = -1; /// Whether to include non-linear correction during global to local /// transformation FreeToBoundCorrection freeToBoundCorrection; From 8325eae1e8b1bf85b730e4d8f9191e7be9954aad Mon Sep 17 00:00:00 2001 From: AJPfleger Date: Thu, 1 Dec 2022 14:47:04 +0100 Subject: [PATCH 5/9] implement review benjaminhuth --- Core/include/Acts/TrackFitting/Chi2Fitter.hpp | 45 +++++++++++-------- 1 file changed, 26 insertions(+), 19 deletions(-) diff --git a/Core/include/Acts/TrackFitting/Chi2Fitter.hpp b/Core/include/Acts/TrackFitting/Chi2Fitter.hpp index bf93e7e429c..892d796c3da 100644 --- a/Core/include/Acts/TrackFitting/Chi2Fitter.hpp +++ b/Core/include/Acts/TrackFitting/Chi2Fitter.hpp @@ -165,7 +165,7 @@ struct Chi2FitterResult { // This correspond to the last state in the multitrajectory. // Since this GX2F only stores one trajectory, it is unambiguous. // SIZE_MAX is the start of a trajectory. - size_t lastTrackIndex = SIZE_MAX; + size_t lastTrackIndex = Acts::MultiTrajectoryTraits::kInvalid; // The optional Parameters at the provided surface std::optional fittedParameters; @@ -350,22 +350,33 @@ class Chi2Fitter { // add a full TrackState entry multi trajectory // (this allocates storage for all components, we will set them later) - size_t currentTrackIndex = SIZE_MAX; + size_t currentTrackIndex = Acts::MultiTrajectoryTraits::kInvalid; + bool foundExistingSurface = + false; // Checks if during the update an existing surface is found. + // If not, there will be a new index generated afterwards + if (updateNumber == 0) { - result.lastTrackIndex = - result.fittedStates - ->addTrackState( // which lastTrackIndex is this? - ~(TrackStatePropMask::Smoothed | - TrackStatePropMask::Filtered), - result.lastTrackIndex); + result.lastTrackIndex = result.fittedStates->addTrackState( + ~(TrackStatePropMask::Smoothed | TrackStatePropMask::Filtered), + result.lastTrackIndex); currentTrackIndex = result.lastTrackIndex; } else { result.fittedStates->visitBackwards( result.lastTrackIndex, [&](auto proxy) { if (&proxy.referenceSurface() == surface) { currentTrackIndex = proxy.index(); + foundExistingSurface = true; } }); + + if (!foundExistingSurface) { + ACTS_VERBOSE( + "chi2 | processSurface: Found new surface during update."); + result.lastTrackIndex = result.fittedStates->addTrackState( + ~(TrackStatePropMask::Smoothed | TrackStatePropMask::Filtered), + result.lastTrackIndex); + currentTrackIndex = result.lastTrackIndex; + } } // now get track state proxy back @@ -608,12 +619,10 @@ class Chi2Fitter { /// the fit. /// /// @return the output as an output track - template + template Result> fit( source_link_iterator_t it, source_link_iterator_t end, const BoundTrackParameters& sParameters, - // const start_parameters_t& sParameters, const Chi2FitterOptions& chi2FitterOptions, std::shared_ptr trajectory = {}) const { const auto& logger = chi2FitterOptions.logger; @@ -634,8 +643,8 @@ class Chi2Fitter { // propagation. Use dynamic Matrix instead? Performance? // Create the ActionList and AbortList - using Chi2Aborter = Aborter; - using Chi2Actor = Actor; + using Chi2Aborter = Aborter; + using Chi2Actor = Actor; using Chi2Result = typename Chi2Actor::result_type; using Actors = ActionList; @@ -644,14 +653,12 @@ class Chi2Fitter { // the result object which will be returned. Overridden every iteration. Chi2Result c2r; - if (not trajectory) { - trajectory = std::make_shared(); - }; + trajectory = std::make_shared(); using paramType = typename std::conditional< - std::is_same::value, + std::is_same::value, std::variant, - std::variant>::type; + std::variant>::type; paramType vParams = sParameters; auto updatedStartParameters = sParameters; @@ -718,11 +725,11 @@ class Chi2Fitter { c2rCurrent.collectorCovariance.size()); c2rCurrent.covariance = variance.asDiagonal(); + // calculate the global chisquare c2rCurrent.chisquare = c2rCurrent.residuals.transpose() * c2rCurrent.covariance.inverse() * c2rCurrent.residuals; ACTS_VERBOSE("chi2 | it=" << i << " | χ² = " << c2rCurrent.chisquare); - // chisquare is here the global chisquare // copy over data from previous runs (namely chisquares vector) c2rCurrent.chisquares.reserve(c2r.chisquares.size() + 1); From f7b82a9bc7d0845005149d1e0f3fbdd49c92c152 Mon Sep 17 00:00:00 2001 From: AJPfleger Date: Mon, 5 Dec 2022 10:59:18 +0100 Subject: [PATCH 6/9] changes review benjaminhuth --- Core/include/Acts/TrackFitting/Chi2Fitter.hpp | 7 +------ 1 file changed, 1 insertion(+), 6 deletions(-) diff --git a/Core/include/Acts/TrackFitting/Chi2Fitter.hpp b/Core/include/Acts/TrackFitting/Chi2Fitter.hpp index 892d796c3da..cf38f5a1d76 100644 --- a/Core/include/Acts/TrackFitting/Chi2Fitter.hpp +++ b/Core/include/Acts/TrackFitting/Chi2Fitter.hpp @@ -655,12 +655,7 @@ class Chi2Fitter { trajectory = std::make_shared(); - using paramType = typename std::conditional< - std::is_same::value, - std::variant, - std::variant>::type; - - paramType vParams = sParameters; + std::variant vParams = sParameters; auto updatedStartParameters = sParameters; for (int i = 0; i <= chi2FitterOptions.nUpdates; ++i) { From 2319d56d141f42633cd33a3286f2e6b5e24c16c8 Mon Sep 17 00:00:00 2001 From: AJPfleger Date: Mon, 5 Dec 2022 14:57:23 +0100 Subject: [PATCH 7/9] remove variant/visit --- Core/include/Acts/TrackFitting/Chi2Fitter.hpp | 30 +++++++------------ 1 file changed, 10 insertions(+), 20 deletions(-) diff --git a/Core/include/Acts/TrackFitting/Chi2Fitter.hpp b/Core/include/Acts/TrackFitting/Chi2Fitter.hpp index cf38f5a1d76..cd6ca35fc30 100644 --- a/Core/include/Acts/TrackFitting/Chi2Fitter.hpp +++ b/Core/include/Acts/TrackFitting/Chi2Fitter.hpp @@ -655,7 +655,7 @@ class Chi2Fitter { trajectory = std::make_shared(); - std::variant vParams = sParameters; + BoundTrackParameters vParams = sParameters; auto updatedStartParameters = sParameters; for (int i = 0; i <= chi2FitterOptions.nUpdates; ++i) { @@ -738,13 +738,9 @@ class Chi2Fitter { if (i == chi2FitterOptions.nUpdates) { // don't update parameters in last iteration - c2r.fittedParameters = std::visit( - [](auto&& prevParams) { - return BoundTrackParameters( - prevParams.referenceSurface().getSharedPtr(), - prevParams.parameters(), prevParams.covariance()); - }, - vParams); + c2r.fittedParameters = + BoundTrackParameters(vParams.referenceSurface().getSharedPtr(), + vParams.parameters(), vParams.covariance()); break; // TODO: verify if another step would be useful, e.g. by comparing the @@ -756,18 +752,12 @@ class Chi2Fitter { c2r.collectorDerive2Chi2Sum.colPivHouseholderQr().solve( c2r.collectorDerive1Chi2Sum); - c2r.fittedParameters = std::visit( - [delta_start_parameters, logger, i](auto&& prevParams) { - BoundVector newParamsVec = - prevParams.parameters() - delta_start_parameters; - ACTS_VERBOSE("chi2 | it=" << i << " | updated parameters = " - << newParamsVec.transpose()); - - return BoundTrackParameters( - prevParams.referenceSurface().getSharedPtr(), newParamsVec, - prevParams.covariance()); - }, - vParams); + BoundVector newParamsVec = vParams.parameters() - delta_start_parameters; + ACTS_VERBOSE("chi2 | it=" << i << " | updated parameters = " + << newParamsVec.transpose()); + c2r.fittedParameters = + BoundTrackParameters(vParams.referenceSurface().getSharedPtr(), + newParamsVec, vParams.covariance()); vParams = c2r.fittedParameters.value(); // passed to next iteration updatedStartParameters = c2r.fittedParameters.value(); From 73b69c9da0069c8b17d52680ff24c8969106120b Mon Sep 17 00:00:00 2001 From: AJPfleger Date: Mon, 5 Dec 2022 15:49:58 +0100 Subject: [PATCH 8/9] improve unittest --- Tests/UnitTests/Core/TrackFitting/Chi2FitterTests.cpp | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/Tests/UnitTests/Core/TrackFitting/Chi2FitterTests.cpp b/Tests/UnitTests/Core/TrackFitting/Chi2FitterTests.cpp index c67eae76141..4cf21c482d5 100644 --- a/Tests/UnitTests/Core/TrackFitting/Chi2FitterTests.cpp +++ b/Tests/UnitTests/Core/TrackFitting/Chi2FitterTests.cpp @@ -186,7 +186,8 @@ BOOST_AUTO_TEST_CASE(ZeroFieldNoSurfaceForward) { LoggerWrapper{*chi2Logger}, PropagatorPlainOptions()); - // chi2Options.nUpdates = 2; // χ² = 17.9695 -> 11.0035 -> 11.0035 ... + chi2Options.nUpdates = 5; + // χ² = 9.14513 -> 6.34088 -> 6.34088 ... // BOOST_TEST_INFO("Test Case ZeroFieldNoSurfaceForward: running .fit()..."); From 68b0a150d48fb268a7d29de74d0fa219c76b1e25 Mon Sep 17 00:00:00 2001 From: AJPfleger Date: Thu, 8 Dec 2022 16:06:47 +0100 Subject: [PATCH 9/9] add: write correct chi2 to root --- Core/include/Acts/TrackFitting/Chi2Fitter.hpp | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/Core/include/Acts/TrackFitting/Chi2Fitter.hpp b/Core/include/Acts/TrackFitting/Chi2Fitter.hpp index cd6ca35fc30..8fd14f290fa 100644 --- a/Core/include/Acts/TrackFitting/Chi2Fitter.hpp +++ b/Core/include/Acts/TrackFitting/Chi2Fitter.hpp @@ -428,6 +428,10 @@ class Chi2Fitter { result.collectorDerive1Chi2Sum += derive1Chi2; result.collectorDerive2Chi2Sum += derive2Chi2; + double localChi2 = + (residuals.transpose() * covInv * residuals).eval()(0); + trackStateProxy.chi2() = localChi2; + for (int i = 0; i < localMeasurements.rows(); ++i) { result.collectorMeasurements.push_back(localMeasurements(i)); result.collectorResiduals.push_back(residuals(i));