Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-07-11 07:50:14

0001 // This file is part of the ACTS project.
0002 //
0003 // Copyright (C) 2016 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 https://mozilla.org/MPL/2.0/.
0008 
0009 #include "Acts/Propagator/detail/CovarianceEngine.hpp"
0010 
0011 #include "Acts/Definitions/Common.hpp"
0012 #include "Acts/Definitions/TrackParametrization.hpp"
0013 #include "Acts/EventData/GenericBoundTrackParameters.hpp"
0014 #include "Acts/EventData/TrackParameters.hpp"
0015 #include "Acts/EventData/TransformationHelpers.hpp"
0016 #include "Acts/EventData/detail/CorrectedTransformationFreeToBound.hpp"
0017 #include "Acts/Propagator/detail/JacobianEngine.hpp"
0018 #include "Acts/Utilities/Result.hpp"
0019 
0020 #include <optional>
0021 #include <system_error>
0022 #include <utility>
0023 
0024 namespace Acts {
0025 
0026 /// Some type defs
0027 using Jacobian = BoundMatrix;
0028 using BoundState = std::tuple<BoundTrackParameters, Jacobian, double>;
0029 
0030 Result<BoundState> detail::boundState(
0031     const GeometryContext& geoContext, const Surface& surface,
0032     BoundSquareMatrix& boundCovariance, BoundMatrix& fullTransportJacobian,
0033     FreeMatrix& freeTransportJacobian, FreeVector& freeToPathDerivatives,
0034     BoundToFreeMatrix& boundToFreeJacobian,
0035     const std::optional<FreeMatrix>& additionalFreeCovariance,
0036     FreeVector& freeParameters, const ParticleHypothesis& particleHypothesis,
0037     bool covTransport, double accumulatedPath,
0038     const FreeToBoundCorrection& freeToBoundCorrection) {
0039   // Create the bound parameters
0040   Result<BoundVector> bv =
0041       transformFreeToBoundParameters(freeParameters, surface, geoContext);
0042   if (!bv.ok()) {
0043     return bv.error();
0044   }
0045 
0046   // Covariance transport
0047   std::optional<BoundSquareMatrix> cov = std::nullopt;
0048   if (covTransport) {
0049     // Calculate the jacobian and transport the covarianceMatrix to final local.
0050     // Then reinitialize the transportJacobian, derivatives and the
0051     // boundToFreeJacobian
0052     transportCovarianceToBound(
0053         geoContext, surface, boundCovariance, fullTransportJacobian,
0054         freeTransportJacobian, freeToPathDerivatives, boundToFreeJacobian,
0055         additionalFreeCovariance, freeParameters, freeToBoundCorrection);
0056     cov = boundCovariance;
0057   }
0058 
0059   // Create the bound state
0060   return std::make_tuple(
0061       BoundTrackParameters(surface.getSharedPtr(), *bv, std::move(cov),
0062                            particleHypothesis),
0063       fullTransportJacobian, accumulatedPath);
0064 }
0065 
0066 BoundState detail::curvilinearState(
0067     BoundSquareMatrix& boundCovariance, BoundMatrix& fullTransportJacobian,
0068     FreeMatrix& freeTransportJacobian, FreeVector& freeToPathDerivatives,
0069     BoundToFreeMatrix& boundToFreeJacobian,
0070     const std::optional<FreeMatrix>& additionalFreeCovariance,
0071     const FreeVector& freeParameters,
0072     const ParticleHypothesis& particleHypothesis, bool covTransport,
0073     double accumulatedPath) {
0074   const Vector3& direction = freeParameters.segment<3>(eFreeDir0);
0075 
0076   // Covariance transport
0077   std::optional<BoundSquareMatrix> cov = std::nullopt;
0078   if (covTransport) {
0079     // Calculate the jacobian and transport the covarianceMatrix to final local.
0080     // Then reinitialize the transportJacobian, derivatives and the
0081     // boundToFreeJacobian
0082     transportCovarianceToCurvilinear(boundCovariance, fullTransportJacobian,
0083                                      freeTransportJacobian,
0084                                      freeToPathDerivatives, boundToFreeJacobian,
0085                                      additionalFreeCovariance, direction);
0086     cov = boundCovariance;
0087   }
0088 
0089   // Create the curvilinear parameters
0090   Vector4 pos4 = Vector4::Zero();
0091   pos4[ePos0] = freeParameters[eFreePos0];
0092   pos4[ePos1] = freeParameters[eFreePos1];
0093   pos4[ePos2] = freeParameters[eFreePos2];
0094   pos4[eTime] = freeParameters[eFreeTime];
0095   BoundTrackParameters curvilinearParams =
0096       BoundTrackParameters::createCurvilinear(
0097           pos4, direction, freeParameters[eFreeQOverP], std::move(cov),
0098           particleHypothesis);
0099   // Create the curvilinear state
0100   return {std::move(curvilinearParams), fullTransportJacobian, accumulatedPath};
0101 }
0102 
0103 void detail::transportCovarianceToBound(
0104     const GeometryContext& geoContext, const Surface& surface,
0105     BoundSquareMatrix& boundCovariance, BoundMatrix& fullTransportJacobian,
0106     FreeMatrix& freeTransportJacobian, FreeVector& freeToPathDerivatives,
0107     BoundToFreeMatrix& boundToFreeJacobian,
0108     const std::optional<FreeMatrix>& additionalFreeCovariance,
0109     FreeVector& freeParameters,
0110     const FreeToBoundCorrection& freeToBoundCorrection) {
0111   FreeToBoundMatrix freeToBoundJacobian;
0112 
0113   // Calculate the full jacobian from local parameters at the start surface to
0114   // current bound parameters
0115   boundToBoundTransportJacobian(geoContext, surface, freeParameters,
0116                                 boundToFreeJacobian, freeTransportJacobian,
0117                                 freeToBoundJacobian, freeToPathDerivatives,
0118                                 fullTransportJacobian);
0119 
0120   bool correction = false;
0121   if (freeToBoundCorrection) {
0122     BoundToFreeMatrix startBoundToFinalFreeJacobian =
0123         freeTransportJacobian * boundToFreeJacobian;
0124     FreeSquareMatrix freeCovariance = startBoundToFinalFreeJacobian *
0125                                       boundCovariance *
0126                                       startBoundToFinalFreeJacobian.transpose();
0127 
0128     auto transformer =
0129         detail::CorrectedFreeToBoundTransformer(freeToBoundCorrection);
0130     auto correctedRes =
0131         transformer(freeParameters, freeCovariance, surface, geoContext);
0132 
0133     if (correctedRes.has_value()) {
0134       auto correctedValue = correctedRes.value();
0135       BoundVector boundParams = std::get<BoundVector>(correctedValue);
0136       // 1. Update the free parameters with the corrected bound parameters
0137       freeParameters =
0138           transformBoundToFreeParameters(surface, geoContext, boundParams);
0139 
0140       // 2. Update the bound covariance
0141       boundCovariance = std::get<BoundSquareMatrix>(correctedValue);
0142 
0143       correction = true;
0144     }
0145   }
0146 
0147   if (!correction) {
0148     // Apply the actual covariance transport to get covariance of the current
0149     // bound parameters
0150     boundCovariance = fullTransportJacobian * boundCovariance *
0151                       fullTransportJacobian.transpose();
0152   }
0153 
0154   if (additionalFreeCovariance) {
0155     boundCovariance += freeToBoundJacobian * (*additionalFreeCovariance) *
0156                        freeToBoundJacobian.transpose();
0157   }
0158 
0159   // Reinitialize jacobian components:
0160   // ->The transportJacobian is reinitialized to Identity
0161   // ->The derivatives is reinitialized to Zero
0162   // ->The boundToFreeJacobian is initialized to that at the current surface
0163   reinitializeJacobians(geoContext, surface, freeTransportJacobian,
0164                         freeToPathDerivatives, boundToFreeJacobian,
0165                         freeParameters);
0166 }
0167 
0168 void detail::transportCovarianceToCurvilinear(
0169     BoundSquareMatrix& boundCovariance, BoundMatrix& fullTransportJacobian,
0170     FreeMatrix& freeTransportJacobian, FreeVector& freeToPathDerivatives,
0171     BoundToFreeMatrix& boundToFreeJacobian,
0172     const std::optional<FreeMatrix>& additionalFreeCovariance,
0173     const Vector3& direction) {
0174   FreeToBoundMatrix freeToBoundJacobian;
0175 
0176   // Calculate the full jacobian from local parameters at the start surface to
0177   // current curvilinear parameters
0178   boundToCurvilinearTransportJacobian(
0179       direction, boundToFreeJacobian, freeTransportJacobian,
0180       freeToBoundJacobian, freeToPathDerivatives, fullTransportJacobian);
0181 
0182   // Apply the actual covariance transport to get covariance of the current
0183   // curvilinear parameters
0184   boundCovariance = fullTransportJacobian * boundCovariance *
0185                     fullTransportJacobian.transpose();
0186 
0187   if (additionalFreeCovariance) {
0188     boundCovariance += freeToBoundJacobian * (*additionalFreeCovariance) *
0189                        freeToBoundJacobian.transpose();
0190   }
0191 
0192   // Reinitialize jacobian components:
0193   // ->The free transportJacobian is reinitialized to Identity
0194   // ->The path derivatives is reinitialized to Zero
0195   // ->The boundToFreeJacobian is reinitialized to that at the current
0196   // curvilinear surface
0197   reinitializeJacobians(freeTransportJacobian, freeToPathDerivatives,
0198                         boundToFreeJacobian, direction);
0199 }
0200 
0201 Result<BoundTrackParameters> detail::boundToBoundConversion(
0202     const GeometryContext& gctx, const BoundTrackParameters& boundParameters,
0203     const Surface& targetSurface, const Vector3& bField) {
0204   const auto& sourceSurface = boundParameters.referenceSurface();
0205 
0206   FreeVector freePars = transformBoundToFreeParameters(
0207       sourceSurface, gctx, boundParameters.parameters());
0208 
0209   auto res = transformFreeToBoundParameters(freePars, targetSurface, gctx);
0210 
0211   if (!res.ok()) {
0212     return res.error();
0213   }
0214   BoundVector parOut = *res;
0215 
0216   std::optional<BoundMatrix> covOut = std::nullopt;
0217 
0218   if (boundParameters.covariance().has_value()) {
0219     BoundToFreeMatrix boundToFreeJacobian = sourceSurface.boundToFreeJacobian(
0220         gctx, freePars.segment<3>(eFreePos0), freePars.segment<3>(eFreeDir0));
0221 
0222     FreeMatrix freeTransportJacobian = FreeMatrix::Identity();
0223 
0224     FreeVector freeToPathDerivatives = FreeVector::Zero();
0225     freeToPathDerivatives.head<3>() = freePars.segment<3>(eFreeDir0);
0226 
0227     freeToPathDerivatives.segment<3>(eFreeDir0) =
0228         bField.cross(freePars.segment<3>(eFreeDir0));
0229 
0230     BoundMatrix boundToBoundJac;
0231     FreeToBoundMatrix freeToBoundJacobian;
0232     detail::boundToBoundTransportJacobian(
0233         gctx, targetSurface, freePars, boundToFreeJacobian,
0234         freeTransportJacobian, freeToBoundJacobian, freeToPathDerivatives,
0235         boundToBoundJac);
0236 
0237     covOut = boundToBoundJac * (*boundParameters.covariance()) *
0238              boundToBoundJac.transpose();
0239   }
0240 
0241   return BoundTrackParameters{targetSurface.getSharedPtr(), parOut, covOut,
0242                               boundParameters.particleHypothesis()};
0243 }
0244 
0245 }  // namespace Acts