Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-01-19 09:23:18

0001 // This file is part of the Acts project.
0002 //
0003 // Copyright (C) 2019-2024 CERN for the benefit of the Acts project
0004 //
0005 // This Source Code Form is subject to the terms of the Mozilla Public
0006 // License, v. 2.0. If a copy of the MPL was not distributed with this
0007 // file, You can obtain one at http://mozilla.org/MPL/2.0/.
0008 
0009 #pragma once
0010 
0011 #include "Acts/Definitions/TrackParametrization.hpp"
0012 
0013 #include <type_traits>
0014 
0015 namespace Acts::detail {
0016 
0017 /// @brief check and correct covariance matrix
0018 ///
0019 /// @tparam CovMatrix_t The type of covariance matrix
0020 /// @tparam NumIter The number of iterations to run the correction
0021 ///
0022 /// Invocation:
0023 ///   - covariance_helper<CovMatrix_t, numIter>::validate(covariance)
0024 ///    The 'covariance' is checked against semi-positivedefiniteness
0025 ///    and limited number of iterations to replace it with the
0026 ///    closest semi-positivedefinite are made if it's not
0027 ///
0028 /// @return The (corrected) covariance is semi-positivedefinite or not
0029 template <typename CovMatrix_t, signed int NumIter = 1>
0030 struct CovarianceHelper {
0031   /// check if the covariance is semi-positive and correction is attempted
0032   static bool validate(CovMatrix_t& covariance) {
0033     if (covariance.hasNaN()) {
0034       return false;
0035     }
0036     std::size_t nIteration = 0;
0037     while (nIteration < NumIter) {
0038       if (isSemiPositive(covariance)) {
0039         return true;
0040       } else {
0041         Eigen::JacobiSVD<CovMatrix_t> svdCov(
0042             covariance, Eigen::ComputeFullU | Eigen::ComputeFullV);
0043         CovMatrix_t S = svdCov.singularValues().asDiagonal();
0044         CovMatrix_t V = svdCov.matrixV();
0045         CovMatrix_t H = V * S * V.transpose();
0046         covariance = (covariance + H) / 2;
0047         nIteration++;
0048       }
0049     }
0050     /// check again after the iterations
0051     return isSemiPositive(covariance);
0052   }
0053 
0054   /// check if the covariance is semi-positive
0055   static bool isSemiPositive(const CovMatrix_t& covariance) {
0056     if (covariance.hasNaN()) {
0057       return false;
0058     }
0059     Eigen::LDLT<CovMatrix_t> ldltCov(covariance);
0060     return ldltCov.isPositive();
0061   }
0062 
0063   /// check if the covariance is positive
0064   static bool isPositive(const CovMatrix_t& covariance) {
0065     if (covariance.hasNaN()) {
0066       return false;
0067     }
0068     Eigen::LLT<CovMatrix_t> lltCov(covariance);
0069     return lltCov.info() == Eigen::Success ? true : false;
0070   }
0071 };
0072 
0073 }  // namespace Acts::detail