From 69a8cb44dd96d4d1d05e9bc34e6dd6480b423e84 Mon Sep 17 00:00:00 2001 From: Andreas Stefl Date: Tue, 29 Oct 2024 10:15:07 +0100 Subject: [PATCH] feat: Validate track parameters (#3756) Adds two free functions which allow checking `BoundVector`s and `FreeVectors`s for valid track parameters. This is then used as `assert`s in the track parameter EDM. blocked by - https://github.com/acts-project/acts/pull/3777 --- .../EventData/GenericBoundTrackParameters.hpp | 7 +- .../EventData/GenericFreeTrackParameters.hpp | 11 +++- .../Acts/EventData/TrackParameterHelpers.hpp | 30 +++++++++ .../Acts/EventData/TrackParametersConcept.hpp | 3 +- .../EventData/detail/GenerateParameters.hpp | 8 +++ .../Acts/EventData/detail/TestTrackState.hpp | 6 +- Core/src/EventData/CMakeLists.txt | 1 + Core/src/EventData/TrackParameterHelpers.cpp | 64 +++++++++++++++++++ .../EventData/TrackParameterHelpersTests.cpp | 10 +++ .../Core/Propagator/MultiStepperTests.cpp | 45 ++++++++----- .../Core/Surfaces/LineSurfaceTests.cpp | 6 +- .../Vertexing/ImpactPointEstimatorTests.cpp | 4 ++ .../TrackParametersJsonConverterTests.cpp | 4 +- 13 files changed, 169 insertions(+), 30 deletions(-) create mode 100644 Core/src/EventData/TrackParameterHelpers.cpp diff --git a/Core/include/Acts/EventData/GenericBoundTrackParameters.hpp b/Core/include/Acts/EventData/GenericBoundTrackParameters.hpp index f624cd88aee1..38f297764b1b 100644 --- a/Core/include/Acts/EventData/GenericBoundTrackParameters.hpp +++ b/Core/include/Acts/EventData/GenericBoundTrackParameters.hpp @@ -9,6 +9,7 @@ #pragma once #include "Acts/Definitions/Tolerance.hpp" +#include "Acts/EventData/TrackParameterHelpers.hpp" #include "Acts/EventData/TransformationHelpers.hpp" #include "Acts/EventData/detail/PrintParameters.hpp" #include "Acts/Surfaces/Surface.hpp" @@ -18,7 +19,6 @@ #include #include #include -#include namespace Acts { @@ -61,7 +61,10 @@ class GenericBoundTrackParameters { m_cov(std::move(cov)), m_surface(std::move(surface)), m_particleHypothesis(std::move(particleHypothesis)) { - assert(m_surface); + // TODO set `validateAngleRange` to `true` after fixing caller code + assert(isBoundVectorValid(m_params, false) && + "Invalid bound parameters vector"); + assert(m_surface != nullptr && "Reference surface must not be null"); normalizePhiTheta(); } diff --git a/Core/include/Acts/EventData/GenericFreeTrackParameters.hpp b/Core/include/Acts/EventData/GenericFreeTrackParameters.hpp index 214aa2e1551e..7a2b4b46522e 100644 --- a/Core/include/Acts/EventData/GenericFreeTrackParameters.hpp +++ b/Core/include/Acts/EventData/GenericFreeTrackParameters.hpp @@ -11,6 +11,7 @@ #include "Acts/Definitions/Algebra.hpp" #include "Acts/Definitions/Common.hpp" #include "Acts/Definitions/TrackParametrization.hpp" +#include "Acts/EventData/TrackParameterHelpers.hpp" #include "Acts/EventData/TrackParametersConcept.hpp" #include "Acts/EventData/TransformationHelpers.hpp" #include "Acts/EventData/detail/PrintParameters.hpp" @@ -18,10 +19,8 @@ #include "Acts/Utilities/UnitVectors.hpp" #include "Acts/Utilities/VectorHelpers.hpp" -#include #include #include -#include namespace Acts { @@ -55,7 +54,9 @@ class GenericFreeTrackParameters { ParticleHypothesis particleHypothesis) : m_params(params), m_cov(std::move(cov)), - m_particleHypothesis(std::move(particleHypothesis)) {} + m_particleHypothesis(std::move(particleHypothesis)) { + assert(isFreeVectorValid(m_params) && "Invalid free parameters vector"); + } /// Construct from four-position, direction, absolute momentum, and charge. /// @@ -78,6 +79,8 @@ class GenericFreeTrackParameters { m_params[eFreeDir1] = dir[eMom1]; m_params[eFreeDir2] = dir[eMom2]; m_params[eFreeQOverP] = qOverP; + + assert(isFreeVectorValid(m_params) && "Invalid free parameters vector"); } /// Construct from four-position, angles, absolute momentum, and charge. @@ -103,6 +106,8 @@ class GenericFreeTrackParameters { m_params[eFreeDir1] = dir[eMom1]; m_params[eFreeDir2] = dir[eMom2]; m_params[eFreeQOverP] = qOverP; + + assert(isFreeVectorValid(m_params) && "Invalid free parameters vector"); } /// Converts a free track parameter with a different hypothesis. diff --git a/Core/include/Acts/EventData/TrackParameterHelpers.hpp b/Core/include/Acts/EventData/TrackParameterHelpers.hpp index c9089ff9ab9f..df40c934824d 100644 --- a/Core/include/Acts/EventData/TrackParameterHelpers.hpp +++ b/Core/include/Acts/EventData/TrackParameterHelpers.hpp @@ -11,8 +11,38 @@ #include "Acts/Definitions/TrackParametrization.hpp" #include "Acts/Utilities/detail/periodic.hpp" +#include + namespace Acts { +/// Check if a bound vector is valid. This checks the following: +/// - All values are finite +/// - (optionally) The phi value is in the range [-pi, pi) +/// - (optionally) The theta value is in the range [0, pi] +/// +/// @param v The bound vector to check +/// @param validateAngleRange If true, the phi and theta values are range checked +/// @param epsilon The epsilon to use for the checks +/// @param maxAbsEta The maximum allowed eta value +/// +/// @return True if the bound vector is valid +bool isBoundVectorValid( + const BoundVector& v, bool validateAngleRange, double epsilon = 1e-6, + double maxAbsEta = std::numeric_limits::infinity()); + +/// Check if a free vector is valid. This checks the following: +/// - All values are finite +/// - Direction is normalized +/// +/// @param v The free vector to check +/// @param epsilon The epsilon to use for the checks +/// @param maxAbsEta The maximum allowed eta value +/// +/// @return True if the free vector is valid +bool isFreeVectorValid( + const FreeVector& v, double epsilon = 1e-6, + double maxAbsEta = std::numeric_limits::infinity()); + /// Normalize the bound parameter angles /// /// @param boundParams The bound parameters to normalize diff --git a/Core/include/Acts/EventData/TrackParametersConcept.hpp b/Core/include/Acts/EventData/TrackParametersConcept.hpp index b39255e5ef70..5e4982545a40 100644 --- a/Core/include/Acts/EventData/TrackParametersConcept.hpp +++ b/Core/include/Acts/EventData/TrackParametersConcept.hpp @@ -12,9 +12,7 @@ #include "Acts/Definitions/TrackParametrization.hpp" #include "Acts/Geometry/GeometryContext.hpp" -#include #include -#include namespace Acts { @@ -68,4 +66,5 @@ concept BoundTrackParametersConcept = { p.position(c) } -> std::same_as; }; }; + } // namespace Acts diff --git a/Core/include/Acts/EventData/detail/GenerateParameters.hpp b/Core/include/Acts/EventData/detail/GenerateParameters.hpp index 023314e2527f..1b591dbcd7ae 100644 --- a/Core/include/Acts/EventData/detail/GenerateParameters.hpp +++ b/Core/include/Acts/EventData/detail/GenerateParameters.hpp @@ -195,6 +195,10 @@ struct GenerateBoundParametersOptions { GenerateQoverPOptions qOverP; }; +inline BoundVector someBoundParametersA() { + return {0.0, 0.0, 0.0, std::numbers::pi / 2, 1.0, 0.0}; +} + template inline BoundVector generateBoundParameters( generator_t& rng, const GenerateBoundParametersOptions& options) { @@ -242,6 +246,10 @@ struct GenerateFreeParametersOptions { GenerateQoverPOptions qOverP; }; +inline FreeVector someFreeParametersA() { + return {0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 1.0}; +} + template inline FreeVector generateFreeParameters( generator_t& rng, const GenerateFreeParametersOptions& options) { diff --git a/Core/include/Acts/EventData/detail/TestTrackState.hpp b/Core/include/Acts/EventData/detail/TestTrackState.hpp index c1471566c5bf..514ae073c2f9 100644 --- a/Core/include/Acts/EventData/detail/TestTrackState.hpp +++ b/Core/include/Acts/EventData/detail/TestTrackState.hpp @@ -38,11 +38,11 @@ struct TestTrackState { : surface( CurvilinearSurface(Vector3::Zero(), Vector3::UnitZ()).surface()), // set bogus parameters first since they are not default-constructible - predicted(surface, BoundVector::Zero(), std::nullopt, + predicted(surface, someBoundParametersA(), std::nullopt, ParticleHypothesis::pion()), - filtered(surface, BoundVector::Zero(), std::nullopt, + filtered(surface, someBoundParametersA(), std::nullopt, ParticleHypothesis::pion()), - smoothed(surface, BoundVector::Zero(), std::nullopt, + smoothed(surface, someBoundParametersA(), std::nullopt, ParticleHypothesis::pion()), jacobian(BoundMatrix::Identity()), chi2(std::chi_squared_distribution(measdim)(rng)), diff --git a/Core/src/EventData/CMakeLists.txt b/Core/src/EventData/CMakeLists.txt index 08a697113784..c425d533649c 100644 --- a/Core/src/EventData/CMakeLists.txt +++ b/Core/src/EventData/CMakeLists.txt @@ -8,4 +8,5 @@ target_sources( TrackStatePropMask.cpp VectorMultiTrajectory.cpp VectorTrackContainer.cpp + TrackParameterHelpers.cpp ) diff --git a/Core/src/EventData/TrackParameterHelpers.cpp b/Core/src/EventData/TrackParameterHelpers.cpp new file mode 100644 index 000000000000..9e7ddb40abf6 --- /dev/null +++ b/Core/src/EventData/TrackParameterHelpers.cpp @@ -0,0 +1,64 @@ +// This file is part of the ACTS project. +// +// Copyright (C) 2016 CERN for the benefit of the ACTS project +// +// This Source Code Form is subject to the terms of the Mozilla Public +// License, v. 2.0. If a copy of the MPL was not distributed with this +// file, You can obtain one at https://mozilla.org/MPL/2.0/. + +#include "Acts/EventData/TrackParameterHelpers.hpp" + +#include "Acts/Utilities/AngleHelpers.hpp" +#include "Acts/Utilities/VectorHelpers.hpp" + +#include + +bool Acts::isBoundVectorValid(const BoundVector& v, bool validateAngleRange, + double epsilon, double maxAbsEta) { + constexpr auto pi = std::numbers::pi_v; + + bool loc0Valid = std::isfinite(v[eBoundLoc0]); + bool loc1Valid = std::isfinite(v[eBoundLoc1]); + bool phiValid = std::isfinite(v[eBoundPhi]); + bool thetaValid = std::isfinite(v[eBoundTheta]); + bool qOverPValid = std::isfinite(v[eBoundQOverP]); + bool timeValid = std::isfinite(v[eBoundTime]); + + if (validateAngleRange) { + phiValid = phiValid && (v[eBoundPhi] + epsilon >= -pi) && + (v[eBoundPhi] - epsilon < pi); + thetaValid = thetaValid && (v[eBoundTheta] + epsilon >= 0) && + (v[eBoundTheta] - epsilon <= pi); + } + + bool etaValid = true; + if (std::isfinite(maxAbsEta)) { + double eta = AngleHelpers::etaFromTheta(v[eBoundTheta]); + etaValid = std::isfinite(eta) && (std::abs(eta) - epsilon <= maxAbsEta); + } + + return loc0Valid && loc1Valid && phiValid && thetaValid && qOverPValid && + timeValid && etaValid; +} + +bool Acts::isFreeVectorValid(const FreeVector& v, double epsilon, + double maxAbsEta) { + bool pos0Valid = std::isfinite(v[eFreePos0]); + bool pos1Valid = std::isfinite(v[eFreePos1]); + bool pos2Valid = std::isfinite(v[eFreePos2]); + bool dir0Valid = std::isfinite(v[eFreeDir0]); + bool dir1Valid = std::isfinite(v[eFreeDir1]); + bool dir2Valid = std::isfinite(v[eFreeDir2]); + bool dirValid = (std::abs(v.segment<3>(eFreeDir0).norm() - 1) - epsilon <= 0); + bool qOverPValid = std::isfinite(v[eFreeQOverP]); + bool timeValid = std::isfinite(v[eFreeTime]); + + bool etaValid = true; + if (std::isfinite(maxAbsEta)) { + double eta = VectorHelpers::eta(v.segment<3>(eFreeDir0)); + etaValid = std::isfinite(eta) && (std::abs(eta) - epsilon <= maxAbsEta); + } + + return pos0Valid && pos1Valid && pos2Valid && dir0Valid && dir1Valid && + dir2Valid && dirValid && qOverPValid && timeValid && etaValid; +} diff --git a/Tests/UnitTests/Core/EventData/TrackParameterHelpersTests.cpp b/Tests/UnitTests/Core/EventData/TrackParameterHelpersTests.cpp index 0a6f30645948..bba77ff8ef2d 100644 --- a/Tests/UnitTests/Core/EventData/TrackParameterHelpersTests.cpp +++ b/Tests/UnitTests/Core/EventData/TrackParameterHelpersTests.cpp @@ -14,6 +14,16 @@ BOOST_AUTO_TEST_SUITE(TrackParameterHelpers) +BOOST_AUTO_TEST_CASE(isBoundVectorValid) { + BOOST_CHECK(!Acts::isBoundVectorValid({1, 2, 3, 4, 5, 6}, true)); + BOOST_CHECK(Acts::isBoundVectorValid({1, 2, 1, 1, 5, 6}, true)); +} + +BOOST_AUTO_TEST_CASE(isFreeVectorValid) { + BOOST_CHECK(!Acts::isFreeVectorValid({1, 2, 3, 4, 5, 6, 7, 8})); + BOOST_CHECK(Acts::isFreeVectorValid({1, 2, 3, 4, 1, 0, 0, 8})); +} + BOOST_AUTO_TEST_CASE(normalizeBoundParameters) { CHECK_CLOSE_OR_SMALL(Acts::normalizeBoundParameters({1, 2, 3, 4, 5, 6}), Acts::BoundVector(1, 2, -0.141593, 2.28319, 5, 6), 1e-3, diff --git a/Tests/UnitTests/Core/Propagator/MultiStepperTests.cpp b/Tests/UnitTests/Core/Propagator/MultiStepperTests.cpp index b9842b55420e..e7edefb1f03e 100644 --- a/Tests/UnitTests/Core/Propagator/MultiStepperTests.cpp +++ b/Tests/UnitTests/Core/Propagator/MultiStepperTests.cpp @@ -39,6 +39,7 @@ #include #include #include +#include #include #include #include @@ -124,9 +125,30 @@ auto makeDefaultBoundPars(bool cov = true, std::size_t n = 4, return c; }; + // note that we are using the default random device + std::mt19937 gen; + std::uniform_real_distribution<> locDis(-10.0, 10.0); + std::uniform_real_distribution<> phiDis(-M_PI, M_PI); + std::uniform_real_distribution<> thetaDis(0, M_PI); + std::uniform_real_distribution<> qOverPDis(-10.0, 10.0); + std::uniform_real_distribution<> timeDis(0.0, 100.0); + for (auto i = 0ul; i < n; ++i) { - cmps.push_back({1. / n, ext_pars ? *ext_pars : BoundVector::Random(), - cov ? Opt{make_random_sym_matrix()} : Opt{}}); + BoundVector params = BoundVector::Zero(); + + if (ext_pars) { + params = *ext_pars; + } else { + params[eBoundLoc0] = locDis(gen); + params[eBoundLoc1] = locDis(gen); + params[eBoundPhi] = phiDis(gen); + params[eBoundTheta] = thetaDis(gen); + params[eBoundQOverP] = qOverPDis(gen); + params[eBoundTime] = timeDis(gen); + } + + cmps.push_back( + {1. / n, params, cov ? Opt{make_random_sym_matrix()} : Opt{}}); } auto surface = Acts::CurvilinearSurface(Vector3::Zero(), Vector3{1., 0., 0.}) @@ -430,7 +452,8 @@ void test_multi_stepper_surface_status_update() { std::vector>> cmps(2, {0.5, BoundVector::Zero(), std::nullopt}); std::get(cmps[0])[eBoundTheta] = M_PI_2; - std::get(cmps[1])[eBoundTheta] = -M_PI_2; + std::get(cmps[1])[eBoundPhi] = M_PI; + std::get(cmps[1])[eBoundTheta] = M_PI_2; std::get(cmps[0])[eBoundQOverP] = 1.0; std::get(cmps[1])[eBoundQOverP] = 1.0; @@ -541,7 +564,8 @@ void test_component_bound_state() { std::vector>> cmps(2, {0.5, BoundVector::Zero(), std::nullopt}); std::get(cmps[0])[eBoundTheta] = M_PI_2; - std::get(cmps[1])[eBoundTheta] = -M_PI_2; + std::get(cmps[1])[eBoundPhi] = M_PI; + std::get(cmps[1])[eBoundTheta] = M_PI_2; std::get(cmps[0])[eBoundQOverP] = 1.0; std::get(cmps[1])[eBoundQOverP] = 1.0; @@ -703,18 +727,7 @@ void test_single_component_interface_function() { using MultiState = typename multi_stepper_t::State; using MultiStepper = multi_stepper_t; - std::vector>> - cmps; - for (int i = 0; i < 4; ++i) { - cmps.push_back({0.25, BoundVector::Random(), BoundSquareMatrix::Random()}); - } - - auto surface = - Acts::CurvilinearSurface(Vector3::Zero(), Vector3::Ones().normalized()) - .planeSurface(); - - MultiComponentBoundTrackParameters multi_pars(surface, cmps, - particleHypothesis); + MultiComponentBoundTrackParameters multi_pars = makeDefaultBoundPars(true, 4); MultiState multi_state(geoCtx, magCtx, defaultBField, multi_pars, defaultStepSize); diff --git a/Tests/UnitTests/Core/Surfaces/LineSurfaceTests.cpp b/Tests/UnitTests/Core/Surfaces/LineSurfaceTests.cpp index 0175578a423a..15d930929af2 100644 --- a/Tests/UnitTests/Core/Surfaces/LineSurfaceTests.cpp +++ b/Tests/UnitTests/Core/Surfaces/LineSurfaceTests.cpp @@ -14,6 +14,7 @@ #include "Acts/Definitions/TrackParametrization.hpp" #include "Acts/EventData/ParticleHypothesis.hpp" #include "Acts/EventData/TrackParameters.hpp" +#include "Acts/EventData/detail/GenerateParameters.hpp" #include "Acts/Geometry/GeometryContext.hpp" #include "Acts/Material/HomogeneousSurfaceMaterial.hpp" #include "Acts/Propagator/Propagator.hpp" @@ -343,8 +344,9 @@ BOOST_AUTO_TEST_CASE(LineSurfaceIntersection) { .closest(); CHECK_CLOSE_ABS(intersection.pathLength(), pathLimit, eps); - BoundTrackParameters endParameters{surface, BoundVector::Zero(), std::nullopt, - ParticleHypothesis::pion()}; + BoundTrackParameters endParameters{surface, + detail::Test::someBoundParametersA(), + std::nullopt, ParticleHypothesis::pion()}; { PropagatorOptions options(tgContext, {}); options.direction = Acts::Direction::Forward; diff --git a/Tests/UnitTests/Core/Vertexing/ImpactPointEstimatorTests.cpp b/Tests/UnitTests/Core/Vertexing/ImpactPointEstimatorTests.cpp index c9743e90d7f9..279f2c1a739d 100644 --- a/Tests/UnitTests/Core/Vertexing/ImpactPointEstimatorTests.cpp +++ b/Tests/UnitTests/Core/Vertexing/ImpactPointEstimatorTests.cpp @@ -43,6 +43,7 @@ #include #include #include +#include #include #include #include @@ -371,6 +372,9 @@ BOOST_DATA_TEST_CASE(VertexCompatibility4D, IPs* vertices, d0, l0, vx0, vy0, BoundVector paramVecClose = BoundVector::Zero(); paramVecClose[eBoundLoc0] = d0; paramVecClose[eBoundLoc1] = l0; + paramVecClose[eBoundPhi] = 0; + paramVecClose[eBoundTheta] = std::numbers::pi / 2; + paramVecClose[eBoundQOverP] = 0; paramVecClose[eBoundTime] = vt0 + sgnClose * timeDiffClose; BoundVector paramVecFar = paramVecClose; diff --git a/Tests/UnitTests/Plugins/Json/TrackParametersJsonConverterTests.cpp b/Tests/UnitTests/Plugins/Json/TrackParametersJsonConverterTests.cpp index 566464a18892..5c0882cb12c4 100644 --- a/Tests/UnitTests/Plugins/Json/TrackParametersJsonConverterTests.cpp +++ b/Tests/UnitTests/Plugins/Json/TrackParametersJsonConverterTests.cpp @@ -8,12 +8,12 @@ #include +#include "Acts/EventData/ParticleHypothesis.hpp" +#include "Acts/EventData/TrackParameters.hpp" #include "Acts/Plugins/Json/TrackParametersJsonConverter.hpp" #include "Acts/Surfaces/PlaneSurface.hpp" #include "Acts/Surfaces/RectangleBounds.hpp" -#include "Acts/Tests/CommonHelpers/FloatComparisons.hpp" -#include #include #include