diff --git a/Core/include/Acts/Vertexing/AdaptiveMultiVertexFinder.hpp b/Core/include/Acts/Vertexing/AdaptiveMultiVertexFinder.hpp index 8c2ca780850..23f6759085c 100644 --- a/Core/include/Acts/Vertexing/AdaptiveMultiVertexFinder.hpp +++ b/Core/include/Acts/Vertexing/AdaptiveMultiVertexFinder.hpp @@ -348,18 +348,19 @@ class AdaptiveMultiVertexFinder final : public IVertexFinder { /// @param fitterState The vertex fitter state /// /// @return Keep new vertex - bool keepNewVertex(Vertex& vtx, const std::vector& allVertices, - VertexFitterState& fitterState) const; + Result keepNewVertex(Vertex& vtx, + const std::vector& allVertices, + VertexFitterState& fitterState) const; /// @brief Method that evaluates if the new vertex candidate is /// merged with one of the previously found vertices /// /// @param vtx The vertex candidate - /// @param allVertices All so far found vertices + /// @param allVertices All vertices that were found so far /// - /// @return Vertex is merged - bool isMergedVertex(const Vertex& vtx, - const std::vector& allVertices) const; + /// @return Bool indicating whether the vertex is merged + Result isMergedVertex(const Vertex& vtx, + const std::vector& allVertices) const; /// @brief Method that deletes last vertex from list of all vertices /// and refits all vertices afterwards diff --git a/Core/include/Acts/Vertexing/VertexingError.hpp b/Core/include/Acts/Vertexing/VertexingError.hpp index 3ef90ed463f..3e4931a12ad 100644 --- a/Core/include/Acts/Vertexing/VertexingError.hpp +++ b/Core/include/Acts/Vertexing/VertexingError.hpp @@ -21,6 +21,9 @@ enum class VertexingError { NotConverged, ElementNotFound, NoCovariance, + SingularMatrix, + NonPositiveVariance, + MatrixNotPositiveDefinite, InvalidInput, }; diff --git a/Core/src/Vertexing/AdaptiveMultiVertexFinder.cpp b/Core/src/Vertexing/AdaptiveMultiVertexFinder.cpp index f2d3835a989..87c8d8eed84 100644 --- a/Core/src/Vertexing/AdaptiveMultiVertexFinder.cpp +++ b/Core/src/Vertexing/AdaptiveMultiVertexFinder.cpp @@ -125,10 +125,12 @@ Acts::Result> AdaptiveMultiVertexFinder::find( break; } } - // MARK: fpeMaskBegin(FLTUND, 1, #2590) - bool keepVertex = isGoodVertex && - keepNewVertex(vtxCandidate, allVerticesPtr, fitterState); - // MARK: fpeMaskEnd(FLTUND) + auto keepNewVertexResult = + keepNewVertex(vtxCandidate, allVerticesPtr, fitterState); + if (!keepNewVertexResult.ok()) { + return keepNewVertexResult.error(); + } + bool keepVertex = isGoodVertex && *keepNewVertexResult; ACTS_DEBUG("New vertex will be saved: " << keepVertex); // Delete vertex from allVertices list again if it's not kept @@ -464,7 +466,7 @@ bool AdaptiveMultiVertexFinder::removeTrackIfIncompatible( return true; } -bool AdaptiveMultiVertexFinder::keepNewVertex( +Result AdaptiveMultiVertexFinder::keepNewVertex( Vertex& vtx, const std::vector& allVertices, VertexFitterState& fitterState) const { double contamination = 0.; @@ -475,23 +477,30 @@ bool AdaptiveMultiVertexFinder::keepNewVertex( fitterState.tracksAtVerticesMap.at(std::make_pair(trk, &vtx)); double trackWeight = trkAtVtx.trackWeight; contaminationNum += trackWeight * (1. - trackWeight); + // MARK: fpeMaskBegin(FLTUND, 1, #2590) contaminationDeNom += trackWeight * trackWeight; + // MARK: fpeMaskEnd(FLTUND) } if (contaminationDeNom != 0) { contamination = contaminationNum / contaminationDeNom; } if (contamination > m_cfg.maximumVertexContamination) { - return false; + return Result::success(false); } - if (isMergedVertex(vtx, allVertices)) { - return false; + auto isMergedVertexResult = isMergedVertex(vtx, allVertices); + if (!isMergedVertexResult.ok()) { + return Result::failure(isMergedVertexResult.error()); } - return true; + if (*isMergedVertexResult) { + return Result::success(false); + } + + return Result::success(true); } -bool AdaptiveMultiVertexFinder::isMergedVertex( +Result AdaptiveMultiVertexFinder::isMergedVertex( const Vertex& vtx, const std::vector& allVertices) const { const Vector4& candidatePos = vtx.fullPosition(); const SquareMatrix4& candidateCov = vtx.fullCovariance(); @@ -505,14 +514,26 @@ bool AdaptiveMultiVertexFinder::isMergedVertex( double significance = 0; if (!m_cfg.doFullSplitting) { - const double deltaZPos = otherPos[eZ] - candidatePos[eZ]; - const double sumVarZ = otherCov(eZ, eZ) + candidateCov(eZ, eZ); - if (sumVarZ <= 0) { - // TODO FIXME this should never happen - continue; + if (m_cfg.useTime) { + const Vector2 deltaZT = otherPos.tail<2>() - candidatePos.tail<2>(); + SquareMatrix2 sumCovZT = candidateCov.bottomRightCorner<2, 2>() + + otherCov.bottomRightCorner<2, 2>(); + auto sumCovZTInverse = safeInverse(sumCovZT); + if (!sumCovZTInverse) { + ACTS_ERROR("Vertex z-t covariance matrix is singular."); + return Result::failure(VertexingError::SingularMatrix); + } + significance = std::sqrt(deltaZT.dot(*sumCovZTInverse * deltaZT)); + } else { + const double deltaZPos = otherPos[eZ] - candidatePos[eZ]; + const double sumVarZ = otherCov(eZ, eZ) + candidateCov(eZ, eZ); + if (sumVarZ <= 0) { + ACTS_ERROR("Variance of the vertex's z-coordinate is not positive."); + return Result::failure(VertexingError::SingularMatrix); + } + // Use only z significance + significance = std::abs(deltaZPos) / std::sqrt(sumVarZ); } - // Use only z significance - significance = std::abs(deltaZPos) / std::sqrt(sumVarZ); } else { if (m_cfg.useTime) { // Use 4D information for significance @@ -520,8 +541,8 @@ bool AdaptiveMultiVertexFinder::isMergedVertex( SquareMatrix4 sumCov = candidateCov + otherCov; auto sumCovInverse = safeInverse(sumCov); if (!sumCovInverse) { - // TODO FIXME this should never happen - continue; + ACTS_ERROR("Vertex 4D covariance matrix is singular."); + return Result::failure(VertexingError::SingularMatrix); } significance = std::sqrt(deltaPos.dot(*sumCovInverse * deltaPos)); } else { @@ -531,17 +552,25 @@ bool AdaptiveMultiVertexFinder::isMergedVertex( candidateCov.topLeftCorner<3, 3>() + otherCov.topLeftCorner<3, 3>(); auto sumCovInverse = safeInverse(sumCov); if (!sumCovInverse) { - // TODO FIXME this should never happen - continue; + ACTS_ERROR("Vertex 3D covariance matrix is singular."); + return Result::failure(VertexingError::SingularMatrix); } significance = std::sqrt(deltaPos.dot(*sumCovInverse * deltaPos)); } } + if (significance < 0.) { + ACTS_ERROR( + "Found a negative significance (i.e., a negative chi2) when checking " + "if vertices are merged. This should never happen since the vertex " + "covariance should be positive definite, and thus its inverse " + "should be positive definite as well."); + return Result::failure(VertexingError::MatrixNotPositiveDefinite); + } if (significance < m_cfg.maxMergeVertexSignificance) { - return true; + return Result::success(true); } } - return false; + return Result::success(false); } Acts::Result AdaptiveMultiVertexFinder::deleteLastVertex( diff --git a/Core/src/Vertexing/VertexingError.cpp b/Core/src/Vertexing/VertexingError.cpp index fbc8bc01aa3..922b8a83dc9 100644 --- a/Core/src/Vertexing/VertexingError.cpp +++ b/Core/src/Vertexing/VertexingError.cpp @@ -34,6 +34,12 @@ class VertexingErrorCategory : public std::error_category { return "Unable to find element."; case VertexingError::NoCovariance: return "No covariance provided."; + case VertexingError::SingularMatrix: + return "Encountered non-invertible matrix."; + case VertexingError::NonPositiveVariance: + return "Encountered negative or zero variance."; + case VertexingError::MatrixNotPositiveDefinite: + return "Encountered a matrix that is not positive definite."; case VertexingError::InvalidInput: return "Invalid input provided."; default: diff --git a/Examples/Algorithms/Vertexing/src/AdaptiveMultiVertexFinderAlgorithm.cpp b/Examples/Algorithms/Vertexing/src/AdaptiveMultiVertexFinderAlgorithm.cpp index bba33aa436c..5e236f9fa40 100644 --- a/Examples/Algorithms/Vertexing/src/AdaptiveMultiVertexFinderAlgorithm.cpp +++ b/Examples/Algorithms/Vertexing/src/AdaptiveMultiVertexFinderAlgorithm.cpp @@ -138,7 +138,6 @@ auto ActsExamples::AdaptiveMultiVertexFinderAlgorithm::makeVertexFinder() const // coordinate. We thus need to increase tracksMaxSignificance (i.e., the // maximum chi2 that a track can have to be associated with a vertex). finderConfig.tracksMaxSignificance = 7.5; - // Check if vertices are merged in space and time finderConfig.doFullSplitting = true; // Reset the maximum significance that two vertices can have before they // are considered as merged. The default value 3 is tuned for comparing