Skip to content

Commit

Permalink
ndf in container
Browse files Browse the repository at this point in the history
  • Loading branch information
AJPfleger committed Oct 31, 2024
1 parent 8d7ec5a commit d1906ed
Show file tree
Hide file tree
Showing 2 changed files with 32 additions and 39 deletions.
58 changes: 26 additions & 32 deletions Core/include/Acts/TrackFitting/GlobalChiSquareFitter.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -317,6 +317,24 @@ struct Gx2fSystem {
// Accessor for a modifiable reference to the vector.
Eigen::VectorXd& bVector() { return m_bVector; }

// It automatically deduces if we want to fit e.g. q/p and adjusts itself
// later. We have only 3 cases, because we always have l0, l1, phi, theta:
// - 4: no magnetic field -> q/p is empty
// - 5: no time measurement -> time not fittable
// - 6: full fit
std::size_t findNdf() {
std::size_t ndfSystem;
if (m_aMatrix(4, 4) == 0) {
ndfSystem = 4;
} else if (m_aMatrix(5, 5) == 0) {
ndfSystem = 5;
} else {
ndfSystem = 6;
}

return ndfSystem;
}

private:
/// Number of dimensions of the (extended) system
std::size_t m_nDims;
Expand Down Expand Up @@ -537,9 +555,7 @@ void addMaterialToGx2fSums(
/// @tparam track_proxy_t The type of the track proxy
///
/// @param track A mutable track proxy to operate on
/// @param aMatrixExtended The aMatrix, summing over second derivatives for the fitting system
/// @param bVectorExtended The bVector, summing over first derivatives for the fitting system
/// @param chi2sum Accumulated chi2 value of the system
/// @param extendedSystem All parameters of the current equation system
/// @param countNdf The number of degrees of freedom counted so far
/// @param multipleScattering Flag to consider multiple scattering in the calculation
/// @param scatteringMap Map of geometry identifiers to scattering properties, containing all scattering angles and covariances
Expand Down Expand Up @@ -632,13 +648,11 @@ void fillGx2fSystem(
/// no qop/time fit)
///
/// @param fullCovariancePredicted The covariance matrix to update
/// @param aMatrixExtended The matrix containing the coefficients of the linear system.
/// @param ndfSystem The number of degrees of freedom, determining the size of meaning full block
/// @param extendedSystem All parameters of the current equation system
///
/// @return deltaParams The calculated delta parameters.
void updateGx2fCovarianceParams(BoundMatrix& fullCovariancePredicted,
Eigen::MatrixXd& aMatrixExtended,
const std::size_t ndfSystem);
Gx2fSystem& extendedSystem);

/// Global Chi Square fitter (GX2F) implementation.
///
Expand Down Expand Up @@ -1183,10 +1197,6 @@ class Gx2Fitter {
// and used for the final track
std::size_t tipIndex = Acts::MultiTrajectoryTraits::kInvalid;

// Here we will store, the ndf of the system. It automatically deduces if we
// want to fit e.g. q/p and adjusts itself later.
std::size_t ndfSystem = std::numeric_limits<std::size_t>::max();

// The scatteringMap stores for each visited surface their scattering
// properties
std::unordered_map<GeometryIdentifier, ScatteringProperties> scatteringMap;
Expand Down Expand Up @@ -1331,19 +1341,6 @@ class Gx2Fitter {

chi2sum = extendedSystem.chi2();

// Get required number of degrees of freedom ndfSystem.
// We have only 3 cases, because we always have l0, l1, phi, theta
// 4: no magnetic field -> q/p is empty
// 5: no time measurement -> time not fittable
// 6: full fit
if (extendedSystem.aMatrix()(4, 4) == 0) {
ndfSystem = 4;
} else if (extendedSystem.aMatrix()(5, 5) == 0) {
ndfSystem = 5;
} else {
ndfSystem = 6;
}

// This check takes into account the evaluated dimensions of the
// measurements. To fit, we need at least NDF+1 measurements. However, we
// count n-dimensional measurements for n measurements, reducing the
Expand All @@ -1353,9 +1350,9 @@ class Gx2Fitter {
// We skip the check during the first iteration, since we cannot guarantee
// to hit all/enough measurement surfaces with the initial parameter
// guess.
if ((nUpdate > 0) && (ndfSystem + 1 > countNdf)) {
if ((nUpdate > 0) && (extendedSystem.findNdf() + 1 > countNdf)) {
ACTS_INFO("Not enough measurements. Require "
<< ndfSystem + 1 << ", but only " << countNdf
<< extendedSystem.findNdf() + 1 << ", but only " << countNdf
<< " could be used.");
return Experimental::GlobalChiSquareFitterError::NotEnoughMeasurements;
}
Expand Down Expand Up @@ -1384,16 +1381,14 @@ class Gx2Fitter {
ACTS_INFO("Abort with relChi2changeCutOff after "
<< nUpdate + 1 << "/" << gx2fOptions.nUpdateMax
<< " iterations.");
updateGx2fCovarianceParams(fullCovariancePredicted,
extendedSystem.aMatrix(), ndfSystem);
updateGx2fCovarianceParams(fullCovariancePredicted, extendedSystem);
break;
}

if (extendedSystem.chi2() > oldChi2sum + 1e-5) {
ACTS_DEBUG("chi2 not converging monotonically");

updateGx2fCovarianceParams(fullCovariancePredicted,
extendedSystem.aMatrix(), ndfSystem);
updateGx2fCovarianceParams(fullCovariancePredicted, extendedSystem);
break;
}

Expand All @@ -1411,8 +1406,7 @@ class Gx2Fitter {
return Experimental::GlobalChiSquareFitterError::DidNotConverge;
}

updateGx2fCovarianceParams(fullCovariancePredicted,
extendedSystem.aMatrix(), ndfSystem);
updateGx2fCovarianceParams(fullCovariancePredicted, extendedSystem);
break;
}

Expand Down
13 changes: 6 additions & 7 deletions Core/src/TrackFitting/GlobalChiSquareFitter.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -11,18 +11,17 @@
#include "Acts/Definitions/TrackParametrization.hpp"

void Acts::Experimental::updateGx2fCovarianceParams(
BoundMatrix& fullCovariancePredicted, Eigen::MatrixXd& aMatrixExtended,
const std::size_t ndfSystem) {
BoundMatrix& fullCovariancePredicted, Gx2fSystem& extendedSystem) {
// make invertible
for (int i = 0; i < aMatrixExtended.rows(); ++i) {
if (aMatrixExtended(i, i) == 0.) {
aMatrixExtended(i, i) = 1.;
for (std::size_t i = 0; i < extendedSystem.nDims(); ++i) {
if (extendedSystem.aMatrix()(i, i) == 0.) {
extendedSystem.aMatrix()(i, i) = 1.;
}
}

visit_measurement(ndfSystem, [&](auto N) {
visit_measurement(extendedSystem.findNdf(), [&](auto N) {
fullCovariancePredicted.topLeftCorner<N, N>() =
aMatrixExtended.inverse().topLeftCorner<N, N>();
extendedSystem.aMatrix().inverse().topLeftCorner<N, N>();
});

return;
Expand Down

0 comments on commit d1906ed

Please sign in to comment.