Skip to content

Commit

Permalink
refactor!: Clean track EDM projector
Browse files Browse the repository at this point in the history
  • Loading branch information
andiwand committed Sep 11, 2024
1 parent ae0c21d commit bd418ee
Show file tree
Hide file tree
Showing 18 changed files with 72 additions and 199 deletions.
130 changes: 17 additions & 113 deletions Core/include/Acts/EventData/TrackStateProxy.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -133,11 +133,6 @@ struct TrackStateTraits {
typename detail_lt::DynamicSizeTypes<ReadOnly>::CoefficientsMap;
using EffectiveCalibratedCovariance =
typename detail_lt::DynamicSizeTypes<ReadOnly>::CovarianceMap;

constexpr static auto ProjectorFlags = Eigen::RowMajor | Eigen::AutoAlign;
using Projector = Eigen::Matrix<Scalar, M, eBoundSize, ProjectorFlags>;
using EffectiveProjector = Eigen::Matrix<Scalar, Eigen::Dynamic, eBoundSize,
ProjectorFlags, M, eBoundSize>;
};

/// Proxy object to access a single point on the trajectory.
Expand Down Expand Up @@ -213,17 +208,6 @@ class TrackStateProxy {
/// Sentinel value that indicates an invalid index
static constexpr IndexType kInvalid = kTrackIndexInvalid;

/// Matrix representing the projector (measurement mapping function) for a
/// measurement. This is not a map type, but an actual matrix. This matrix
/// is always \f$M \times M\f$, even if the local measurement dimension is lower.
/// The actual \f$N\times M\f$ projector is given by the top \f$N\f$ rows.
using Projector = typename TrackStateTraits<M, ReadOnly>::Projector;

/// Dynamic variant of the projector matrix
/// @warning Using this type is discouraged, as it has a runtime overhead
using EffectiveProjector =
typename TrackStateTraits<M, ReadOnly>::EffectiveProjector;

/// The track state container backend given as a template parameter
using Trajectory = trajectory_t;

Expand Down Expand Up @@ -621,145 +605,65 @@ class TrackStateProxy {
///
/// @{

/// Returns the projector (measurement mapping function) for this track
/// state. It is derived from the uncalibrated measurement
/// @note This function returns the overallocated projector. This means it
/// is of dimension MxM, where M is the maximum number of measurement
/// dimensions. The NxM submatrix, where N is the actual dimension of the
/// measurement, is located in the top left corner, everything else is zero.
/// @return The overallocated projector
Projector projector() const;

/// Returns whether a projector is set
/// @return Whether it is set
bool hasProjector() const { return has<hashString("projector")>(); }

/// Returns the projector (measurement mapping function) for this track
/// state. It is derived from the uncalibrated measurement
/// @warning This function returns the effective projector. This means it
/// is of dimension \f$N\times M\f$, where \f$N\f$ is the actual dimension of the
/// measurement.
/// @return The effective projector
EffectiveProjector effectiveProjector() const {
return projector().topLeftCorner(calibratedSize(), M);
}

/// Set the projector on this track state
/// This will convert the projector to a more compact bitset representation
/// and store it.
/// @param projector The projector in the form of a dense matrix
/// @note @p projector is assumed to only have 0s or 1s as components.
template <typename Derived>
[[deprecated("use setProjector(span) instead")]] void setProjector(
const Eigen::MatrixBase<Derived>& projector)
requires(!ReadOnly)
{
constexpr int rows = Eigen::MatrixBase<Derived>::RowsAtCompileTime;
constexpr int cols = Eigen::MatrixBase<Derived>::ColsAtCompileTime;

static_assert(rows != -1 && cols != -1,
"Assignment of dynamic matrices is currently not supported.");

BoundSubspaceIndices projectorSubspaceIndices() const {
assert(has<hashString("projector")>());

static_assert(rows <= M, "Given projector has too many rows");
static_assert(cols <= eBoundSize, "Given projector has too many columns");

// set up full size projector with only zeros
typename TrackStateProxy::Projector fullProjector =
decltype(fullProjector)::Zero();

// assign (potentially) smaller actual projector to matrix, preserving
// zeroes outside of smaller matrix block.
fullProjector.template topLeftCorner<rows, cols>() = projector;

// convert to bitset before storing
ProjectorBitset projectorBitset = matrixToBitset(fullProjector).to_ulong();
setProjectorBitset(projectorBitset);
}

/// Directly get the projector bitset, a compressed form of a projection
/// matrix
/// @note This is mainly to copy explicitly a projector from one state
/// to another. Use the `projector` or `effectiveProjector` method if
/// you want to access the matrix.
/// @return The projector bitset
[[deprecated("use projector() instead")]] ProjectorBitset projectorBitset()
const {
return variableBoundSubspaceHelper().projectorBitset();
return deserializeSubspaceIndices<eBoundSize>(
component<SerializedSubspaceIndices, hashString("projector")>());
}

/// Set the projector bitset, a compressed form of a projection matrix
/// @param proj The projector bitset
///
/// @note This is mainly to copy explicitly a projector from one state
/// to another. If you have a projection matrix, set it with
/// `setProjector`.
[[deprecated("use setProjector(span) instead")]] void setProjectorBitset(
ProjectorBitset proj)
void setProjectorSubspaceIndices(BoundSubspaceIndices boundSubspace)
requires(!ReadOnly)
{
BoundMatrix projMatrix = bitsetToMatrix<BoundMatrix>(proj);
BoundSubspaceIndices boundSubspace =
projectorToSubspaceIndices<eBoundSize>(projMatrix);
setBoundSubspaceIndices(boundSubspace);
}

BoundSubspaceIndices boundSubspaceIndices() const {
assert(has<hashString("projector")>());
return deserializeSubspaceIndices<eBoundSize>(
component<SerializedSubspaceIndices, hashString("projector")>());
component<SerializedSubspaceIndices, hashString("projector")>() =
serializeSubspaceIndices(boundSubspace);
}

template <std::size_t measdim>
SubspaceIndices<measdim> subspaceIndices() const {
SubspaceIndices<measdim> projectorSubspaceIndices() const {
BoundSubspaceIndices boundSubspace = BoundSubspaceIndices();
SubspaceIndices<measdim> subspace;
std::copy(boundSubspace.begin(), boundSubspace.begin() + measdim,
subspace.begin());
return subspace;
}

void setBoundSubspaceIndices(BoundSubspaceIndices boundSubspace)
requires(!ReadOnly)
{
assert(has<hashString("projector")>());
component<SerializedSubspaceIndices, hashString("projector")>() =
serializeSubspaceIndices(boundSubspace);
}

template <std::size_t measdim>
void setSubspaceIndices(SubspaceIndices<measdim> subspace)
void setProjectorSubspaceIndices(SubspaceIndices<measdim> subspace)
requires(!ReadOnly && measdim <= eBoundSize)
{
assert(has<hashString("projector")>());
BoundSubspaceIndices boundSubspace{};
std::copy(subspace.begin(), subspace.end(), boundSubspace.begin());
setBoundSubspaceIndices(boundSubspace);
setProjectorSubspaceIndices(boundSubspace);
}

template <std::size_t measdim, typename index_t>
void setSubspaceIndices(std::array<index_t, measdim> subspaceIndices)
void setProjectorSubspaceIndices(std::array<index_t, measdim> subspaceIndices)
requires(!ReadOnly && measdim <= eBoundSize)
{
assert(has<hashString("projector")>());
BoundSubspaceIndices boundSubspace{};
std::transform(subspaceIndices.begin(), subspaceIndices.end(),
boundSubspace.begin(),
[](index_t i) { return static_cast<std::uint8_t>(i); });
setBoundSubspaceIndices(boundSubspace);
setProjectorSubspaceIndices(boundSubspace);
}

VariableBoundSubspaceHelper variableBoundSubspaceHelper() const {
BoundSubspaceIndices boundSubspace = boundSubspaceIndices();
VariableBoundSubspaceHelper projectorSubspaceHelper() const {
BoundSubspaceIndices boundSubspace = projectorSubspaceIndices();
std::span<std::uint8_t> validSubspaceIndices(
boundSubspace.begin(), boundSubspace.begin() + calibratedSize());
return VariableBoundSubspaceHelper(validSubspaceIndices);
}

template <std::size_t measdim>
FixedBoundSubspaceHelper<measdim> fixedBoundSubspaceHelper() const {
SubspaceIndices<measdim> subspace = subspaceIndices<measdim>();
FixedBoundSubspaceHelper<measdim> projectorSubspaceHelper() const {
SubspaceIndices<measdim> subspace = projectorSubspaceIndices<measdim>();
return FixedBoundSubspaceHelper<measdim>(subspace);
}

Expand Down Expand Up @@ -1025,7 +929,7 @@ class TrackStateProxy {
other.template calibratedCovariance<measdim>();
});

setBoundSubspaceIndices(other.boundSubspaceIndices());
setProjectorSubspaceIndices(other.projectorSubspaceIndices());
}
} else {
if (ACTS_CHECK_BIT(mask, PM::Predicted) &&
Expand Down Expand Up @@ -1072,7 +976,7 @@ class TrackStateProxy {
other.template calibratedCovariance<measdim>();
});

setBoundSubspaceIndices(other.boundSubspaceIndices());
setProjectorSubspaceIndices(other.projectorSubspaceIndices());
}
}

Expand Down
5 changes: 0 additions & 5 deletions Core/include/Acts/EventData/TrackStateProxy.ipp
Original file line number Diff line number Diff line change
Expand Up @@ -61,11 +61,6 @@ inline auto TrackStateProxy<D, M, ReadOnly>::covariance() const
}
}

template <typename D, std::size_t M, bool ReadOnly>
inline auto TrackStateProxy<D, M, ReadOnly>::projector() const -> Projector {
return variableBoundSubspaceHelper().fullProjector();
}

template <typename D, std::size_t M, bool ReadOnly>
inline auto TrackStateProxy<D, M, ReadOnly>::getUncalibratedSourceLink() const
-> SourceLink {
Expand Down
6 changes: 0 additions & 6 deletions Core/include/Acts/EventData/TrackStateProxyConcept.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -111,12 +111,6 @@ concept TrackStateProxyConcept =
{ cv.hasProjector() } -> std::same_as<bool>;
{ v.hasProjector() } -> std::same_as<bool>;

{ cv.effectiveProjector() } -> std::same_as<detail::EffectiveProjector>;
{ v.effectiveProjector() } -> std::same_as<detail::EffectiveProjector>;

{ cv.projectorBitset() } -> std::same_as<ProjectorBitset>;
{ v.projectorBitset() } -> std::same_as<ProjectorBitset>;

{ cv.getUncalibratedSourceLink() } -> std::same_as<SourceLink>;
{ v.getUncalibratedSourceLink() } -> std::same_as<SourceLink>;

Expand Down
31 changes: 12 additions & 19 deletions Core/include/Acts/EventData/detail/MultiTrajectoryTestsCommon.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -444,8 +444,8 @@ class MultiTrajectoryTestsCommon {
BOOST_CHECK_EQUAL(
ts.getUncalibratedSourceLink().template get<TestSourceLink>().sourceId,
pc.sourceLink.sourceId);
testSourceLinkCalibratorReturn<trajectory_t>(
gctx, cctx, SourceLink{ttsb.sourceLink}, ts);
testSourceLinkCalibrator<trajectory_t>(gctx, cctx,
SourceLink{ttsb.sourceLink}, ts);
BOOST_CHECK_EQUAL(
ts.getUncalibratedSourceLink().template get<TestSourceLink>().sourceId,
ttsb.sourceLink.sourceId);
Expand Down Expand Up @@ -524,19 +524,8 @@ class MultiTrajectoryTestsCommon {
}

BOOST_CHECK(ts.hasProjector());
ActsMatrix<MultiTrajectoryTraits::MeasurementSizeMax, eBoundSize> fullProj;
fullProj.setZero();
{
Acts::GeometryContext gctx;
Acts::CalibrationContext cctx;
// create a temporary measurement to extract the projector matrix
testSourceLinkCalibratorReturn<trajectory_t>(
gctx, cctx, SourceLink{pc.sourceLink}, ts);
fullProj = ts.projector();
}
BOOST_CHECK_EQUAL(ts.effectiveProjector(),
fullProj.topLeftCorner(nMeasurements, eBoundSize));
BOOST_CHECK_EQUAL(ts.projector(), fullProj);
BOOST_CHECK_EQUAL(ts.projectorSubspaceIndices()[0], eBoundLoc0);
BOOST_CHECK_EQUAL(ts.projectorSubspaceIndices()[1], eBoundLoc1);
}

void testTrackStateProxyAllocations(std::default_random_engine& rng) {
Expand Down Expand Up @@ -755,7 +744,8 @@ class MultiTrajectoryTestsCommon {
});

BOOST_CHECK_NE(ts1.calibratedSize(), ts2.calibratedSize());
BOOST_CHECK_NE(ts1.projector(), ts2.projector());
BOOST_CHECK_NE(ts1.projectorSubspaceIndices(),
ts2.projectorSubspaceIndices());

BOOST_CHECK_NE(ts1.jacobian(), ts2.jacobian());
BOOST_CHECK_NE(ts1.chi2(), ts2.chi2());
Expand Down Expand Up @@ -784,7 +774,8 @@ class MultiTrajectoryTestsCommon {
});

BOOST_CHECK_EQUAL(ts1.calibratedSize(), ts2.calibratedSize());
BOOST_CHECK_EQUAL(ts1.projector(), ts2.projector());
BOOST_CHECK_EQUAL(ts1.projectorSubspaceIndices(),
ts2.projectorSubspaceIndices());

BOOST_CHECK_EQUAL(ts1.jacobian(), ts2.jacobian());
BOOST_CHECK_EQUAL(ts1.chi2(), ts2.chi2());
Expand All @@ -809,7 +800,8 @@ class MultiTrajectoryTestsCommon {
});

BOOST_CHECK_NE(ts1.calibratedSize(), ts2.calibratedSize());
BOOST_CHECK_NE(ts1.projector(), ts2.projector());
BOOST_CHECK_NE(ts1.projectorSubspaceIndices(),
ts2.projectorSubspaceIndices());

BOOST_CHECK_NE(ts1.jacobian(), ts2.jacobian());
BOOST_CHECK_NE(ts1.chi2(), ts2.chi2());
Expand All @@ -831,7 +823,8 @@ class MultiTrajectoryTestsCommon {
});

BOOST_CHECK_EQUAL(ts1.calibratedSize(), ts2.calibratedSize());
BOOST_CHECK_EQUAL(ts1.projector(), ts2.projector());
BOOST_CHECK_EQUAL(ts1.projectorSubspaceIndices(),
ts2.projectorSubspaceIndices());

BOOST_CHECK_EQUAL(ts1.jacobian(), ts2.jacobian());
BOOST_CHECK_EQUAL(ts1.chi2(), ts2.chi2()); // always copied
Expand Down
19 changes: 3 additions & 16 deletions Core/include/Acts/EventData/detail/TestSourceLink.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -125,7 +125,7 @@ inline std::ostream& operator<<(std::ostream& os,
/// @param gctx Unused
/// @param trackState TrackState to calibrated
template <typename trajectory_t>
void testSourceLinkCalibratorReturn(
void testSourceLinkCalibrator(
const GeometryContext& /*gctx*/, const CalibrationContext& /*cctx*/,
const SourceLink& sourceLink,
typename trajectory_t::TrackStateProxy trackState) {
Expand All @@ -138,31 +138,18 @@ void testSourceLinkCalibratorReturn(
trackState.allocateCalibrated(2);
trackState.template calibrated<2>() = sl.parameters;
trackState.template calibratedCovariance<2>() = sl.covariance;
trackState.template setSubspaceIndices(
trackState.template setProjectorSubspaceIndices(
std::array{sl.indices[0], sl.indices[1]});
} else if (sl.indices[0] != Acts::eBoundSize) {
trackState.allocateCalibrated(1);
trackState.template calibrated<1>() = sl.parameters.head<1>();
trackState.template calibratedCovariance<1>() =
sl.covariance.topLeftCorner<1, 1>();
trackState.template setSubspaceIndices(std::array{sl.indices[0]});
trackState.template setProjectorSubspaceIndices(std::array{sl.indices[0]});
} else {
throw std::runtime_error(
"Tried to extract measurement from invalid TestSourceLink");
}
}

/// Extract the measurement from a TestSourceLink.
///
/// @param gctx Unused
/// @param trackState TrackState to calibrated
template <typename trajectory_t>
void testSourceLinkCalibrator(
const GeometryContext& gctx, const CalibrationContext& cctx,
const SourceLink& sourceLink,
typename trajectory_t::TrackStateProxy trackState) {
testSourceLinkCalibratorReturn<trajectory_t>(gctx, cctx, sourceLink,
trackState);
}

} // namespace Acts::detail::Test
2 changes: 1 addition & 1 deletion Core/include/Acts/TrackFinding/MeasurementSelector.ipp
Original file line number Diff line number Diff line change
Expand Up @@ -63,7 +63,7 @@ MeasurementSelector::select(
trackState.effectiveCalibrated().data(),
trackState.effectiveCalibratedCovariance().data(),
trackState.predicted(), trackState.predictedCovariance(),
trackState.boundSubspaceIndices(), trackState.calibratedSize());
trackState.projectorSubspaceIndices(), trackState.calibratedSize());
trackState.chi2() = chi2;

if (chi2 < minChi2) {
Expand Down
2 changes: 1 addition & 1 deletion Core/include/Acts/TrackFitting/GainMatrixUpdater.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -80,7 +80,7 @@ class GainMatrixUpdater {
// shape later
trackState.effectiveCalibrated().data(),
trackState.effectiveCalibratedCovariance().data(),
trackState.boundSubspaceIndices(),
trackState.projectorSubspaceIndices(),
trackState.predicted(),
trackState.predictedCovariance(),
trackState.filtered(),
Expand Down
2 changes: 1 addition & 1 deletion Core/include/Acts/TrackFitting/GlobalChiSquareFitter.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -350,7 +350,7 @@ void addMeasurementToGx2fSums(Eigen::MatrixXd& aMatrixExtended,
trackState.template calibrated<kMeasDim>();

const ActsMatrix<kMeasDim, eBoundSize> projector =
trackState.projector().template topLeftCorner<kMeasDim, eBoundSize>();
trackState.template projectorSubspaceHelper<kMeasDim>().projector();

const Eigen::MatrixXd projJacobian = projector * extendedJacobian;

Expand Down
Loading

0 comments on commit bd418ee

Please sign in to comment.