Skip to content

Commit

Permalink
feat: Validate track parameters (acts-project#3756)
Browse files Browse the repository at this point in the history
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
- acts-project#3777
  • Loading branch information
andiwand authored Oct 29, 2024
1 parent d155d23 commit 9c0bc47
Show file tree
Hide file tree
Showing 13 changed files with 169 additions and 30 deletions.
7 changes: 5 additions & 2 deletions Core/include/Acts/EventData/GenericBoundTrackParameters.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -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"
Expand All @@ -18,7 +19,6 @@
#include <cassert>
#include <cmath>
#include <memory>
#include <type_traits>

namespace Acts {

Expand Down Expand Up @@ -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();
}

Expand Down
11 changes: 8 additions & 3 deletions Core/include/Acts/EventData/GenericFreeTrackParameters.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -11,17 +11,16 @@
#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"
#include "Acts/Utilities/MathHelpers.hpp"
#include "Acts/Utilities/UnitVectors.hpp"
#include "Acts/Utilities/VectorHelpers.hpp"

#include <cassert>
#include <cmath>
#include <optional>
#include <type_traits>

namespace Acts {

Expand Down Expand Up @@ -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.
///
Expand All @@ -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.
Expand All @@ -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.
Expand Down
30 changes: 30 additions & 0 deletions Core/include/Acts/EventData/TrackParameterHelpers.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -11,8 +11,38 @@
#include "Acts/Definitions/TrackParametrization.hpp"
#include "Acts/Utilities/detail/periodic.hpp"

#include <limits>

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<double>::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<double>::infinity());

/// Normalize the bound parameter angles
///
/// @param boundParams The bound parameters to normalize
Expand Down
3 changes: 1 addition & 2 deletions Core/include/Acts/EventData/TrackParametersConcept.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -12,9 +12,7 @@
#include "Acts/Definitions/TrackParametrization.hpp"
#include "Acts/Geometry/GeometryContext.hpp"

#include <concepts>
#include <optional>
#include <type_traits>

namespace Acts {

Expand Down Expand Up @@ -68,4 +66,5 @@ concept BoundTrackParametersConcept =
{ p.position(c) } -> std::same_as<Vector3>;
};
};

} // namespace Acts
8 changes: 8 additions & 0 deletions Core/include/Acts/EventData/detail/GenerateParameters.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -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 <typename generator_t>
inline BoundVector generateBoundParameters(
generator_t& rng, const GenerateBoundParametersOptions& options) {
Expand Down Expand Up @@ -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 <typename generator_t>
inline FreeVector generateFreeParameters(
generator_t& rng, const GenerateFreeParametersOptions& options) {
Expand Down
6 changes: 3 additions & 3 deletions Core/include/Acts/EventData/detail/TestTrackState.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -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<double>(measdim)(rng)),
Expand Down
1 change: 1 addition & 0 deletions Core/src/EventData/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -8,4 +8,5 @@ target_sources(
TrackStatePropMask.cpp
VectorMultiTrajectory.cpp
VectorTrackContainer.cpp
TrackParameterHelpers.cpp
)
64 changes: 64 additions & 0 deletions Core/src/EventData/TrackParameterHelpers.cpp
Original file line number Diff line number Diff line change
@@ -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 <numbers>

bool Acts::isBoundVectorValid(const BoundVector& v, bool validateAngleRange,
double epsilon, double maxAbsEta) {
constexpr auto pi = std::numbers::pi_v<double>;

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;
}
10 changes: 10 additions & 0 deletions Tests/UnitTests/Core/EventData/TrackParameterHelpersTests.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -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,
Expand Down
45 changes: 29 additions & 16 deletions Tests/UnitTests/Core/Propagator/MultiStepperTests.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -39,6 +39,7 @@
#include <cstddef>
#include <memory>
#include <optional>
#include <random>
#include <stdexcept>
#include <tuple>
#include <type_traits>
Expand Down Expand Up @@ -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.})
Expand Down Expand Up @@ -430,7 +452,8 @@ void test_multi_stepper_surface_status_update() {
std::vector<std::tuple<double, BoundVector, std::optional<BoundSquareMatrix>>>
cmps(2, {0.5, BoundVector::Zero(), std::nullopt});
std::get<BoundVector>(cmps[0])[eBoundTheta] = M_PI_2;
std::get<BoundVector>(cmps[1])[eBoundTheta] = -M_PI_2;
std::get<BoundVector>(cmps[1])[eBoundPhi] = M_PI;
std::get<BoundVector>(cmps[1])[eBoundTheta] = M_PI_2;
std::get<BoundVector>(cmps[0])[eBoundQOverP] = 1.0;
std::get<BoundVector>(cmps[1])[eBoundQOverP] = 1.0;

Expand Down Expand Up @@ -541,7 +564,8 @@ void test_component_bound_state() {
std::vector<std::tuple<double, BoundVector, std::optional<BoundSquareMatrix>>>
cmps(2, {0.5, BoundVector::Zero(), std::nullopt});
std::get<BoundVector>(cmps[0])[eBoundTheta] = M_PI_2;
std::get<BoundVector>(cmps[1])[eBoundTheta] = -M_PI_2;
std::get<BoundVector>(cmps[1])[eBoundPhi] = M_PI;
std::get<BoundVector>(cmps[1])[eBoundTheta] = M_PI_2;
std::get<BoundVector>(cmps[0])[eBoundQOverP] = 1.0;
std::get<BoundVector>(cmps[1])[eBoundQOverP] = 1.0;

Expand Down Expand Up @@ -703,18 +727,7 @@ void test_single_component_interface_function() {
using MultiState = typename multi_stepper_t::State;
using MultiStepper = multi_stepper_t;

std::vector<std::tuple<double, BoundVector, std::optional<BoundSquareMatrix>>>
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);
Expand Down
6 changes: 4 additions & 2 deletions Tests/UnitTests/Core/Surfaces/LineSurfaceTests.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -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"
Expand Down Expand Up @@ -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;
Expand Down
4 changes: 4 additions & 0 deletions Tests/UnitTests/Core/Vertexing/ImpactPointEstimatorTests.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -43,6 +43,7 @@
#include <cmath>
#include <limits>
#include <memory>
#include <numbers>
#include <optional>
#include <tuple>
#include <utility>
Expand Down Expand Up @@ -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;
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -8,12 +8,12 @@

#include <boost/test/unit_test.hpp>

#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 <fstream>
#include <memory>

#include <nlohmann/json.hpp>
Expand Down

0 comments on commit 9c0bc47

Please sign in to comment.