Core/EventData/CorrectedTransformFreeToBoundTests.cpp

0001 // This file is part of the Acts project.
0002 //
0003 // Copyright (C) 2022 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 #include <boost/test/data/test_case.hpp>
0010 #include <boost/test/unit_test.hpp>
0011
0012 #include "Acts/Definitions/Algebra.hpp"
0013 #include "Acts/Definitions/TrackParametrization.hpp"
0014 #include "Acts/Definitions/Units.hpp"
0015 #include "Acts/EventData/detail/CorrectedTransformationFreeToBound.hpp"
0016 #include "Acts/Geometry/GeometryContext.hpp"
0017 #include "Acts/Surfaces/BoundaryCheck.hpp"
0018 #include "Acts/Surfaces/PlaneSurface.hpp"
0019 #include "Acts/Surfaces/Surface.hpp"
0020 #include "Acts/Tests/CommonHelpers/FloatComparisons.hpp"
0021 #include "Acts/Utilities/Intersection.hpp"
0022 #include "Acts/Utilities/UnitVectors.hpp"
0023
0024 #include <algorithm>
0025 #include <cmath>
0026 #include <memory>
0027 #include <optional>
0028 #include <ostream>
0029 #include <string>
0030 #include <tuple>
0031
0032 using namespace Acts;
0033 using namespace Acts::UnitLiterals;
0034
0035 namespace {
0036 constexpr ActsScalar eps = 0.01;
0037 }
0038
0039 BOOST_AUTO_TEST_CASE(CorrectedFreeToBoundTrackParameters) {
0040   GeometryContext geoCtx;
0041
0042   const auto loc0 = 0.0;
0043   const auto loc1 = 0.0;
0044   const auto phi = 0.0;
0045   const auto theta = M_PI / 4;
0046   const auto qOverP = 1 / 1_GeV;
0047   const auto t = 1_ns;
0048
0049   const auto resLoc0 = 0.0;
0050   const auto resLoc1 = 0.0;
0051   const auto resPhi = 0.25;
0052   const auto resTheta = 0.25;
0053   const auto resQOverP = 0.01 / 1_GeV;
0054   const auto resTime = 0.01_ns;
0055
0056   // construct two parallel plane surfaces with normal in x direction
0057   ActsScalar distance = 10_mm;
0058   auto eSurface = Surface::makeShared<PlaneSurface>(Vector3(distance, 0, 0),
0059                                                     Vector3::UnitX());
0060
0061   // the bound parameters at the starting plane
0062   BoundVector sBoundParams = BoundVector::Zero();
0063   sBoundParams << loc0, loc1, phi, theta, qOverP, t;
0064
0065   // the bound parameters covariance at the starting  plane
0066   BoundSquareMatrix sBoundCov = BoundSquareMatrix::Zero();
0067   sBoundCov(eBoundLoc0, eBoundLoc0) = resLoc0 * resLoc0;
0068   sBoundCov(eBoundLoc1, eBoundLoc1) = resLoc1 * resLoc1;
0069   sBoundCov(eBoundPhi, eBoundPhi) = resPhi * resPhi;
0070   sBoundCov(eBoundTheta, eBoundTheta) = resTheta * resTheta;
0071   sBoundCov(eBoundQOverP, eBoundQOverP) = resQOverP * resQOverP;
0072   sBoundCov(eBoundTime, eBoundTime) = resTime * resTime;
0073
0074   Vector3 dir = makeDirectionFromPhiTheta(phi, theta);
0075
0076   // the intersection of the track with the end surface
0077   SurfaceIntersection intersection =
0078       eSurface->intersect(geoCtx, Vector3(0, 0, 0), dir, BoundaryCheck(true))
0079           .closest();
0080   Vector3 tpos = intersection.position();
0081   auto s = intersection.pathLength();
0082
0083   BOOST_CHECK_EQUAL(s, distance * std::sqrt(2));
0084
0085   // construct the free parameters vector
0086   FreeVector eFreeParams = FreeVector::Zero();
0087   eFreeParams.segment<3>(eFreePos0) = tpos;
0088   eFreeParams[eFreeTime] = t;
0089   eFreeParams.segment<3>(eFreeDir0) = dir;
0090   eFreeParams[eFreeQOverP] = qOverP;
0091
0092   // the jacobian from local to global at the starting position
0093   BoundToFreeMatrix boundToFreeJac =
0094       eSurface->boundToFreeJacobian(geoCtx, tpos, dir);
0095
0096   // the transport jacobian without B field
0097   FreeMatrix transportJac = FreeMatrix::Identity();
0098   transportJac(eFreePos0, eFreeDir0) = s;
0099   transportJac(eFreePos1, eFreeDir1) = s;
0100   transportJac(eFreePos2, eFreeDir2) = s;
0101
0102   // the free covariance at the start position
0103   FreeSquareMatrix sFreeCov =
0104       boundToFreeJac * sBoundCov * boundToFreeJac.transpose();
0105   // the free covariance at the end position
0106   FreeSquareMatrix eFreeCov =
0107       transportJac * sFreeCov * transportJac.transpose();
0108
0109   // convert free parameters to bound parameters with non-linear correction
0110
0111   BOOST_TEST_INFO("Transform free parameters vector onto surface "
0112                   << eSurface->name());
0113
0114   // the corrected transformation
0115   auto freeToBoundCorrection = FreeToBoundCorrection(true);
0116   BOOST_CHECK(freeToBoundCorrection);
0117
0118   auto transformer =
0119       detail::CorrectedFreeToBoundTransformer(freeToBoundCorrection);
0120   auto correctedRes = transformer(eFreeParams, eFreeCov, *eSurface, geoCtx);
0121
0122   BOOST_CHECK(correctedRes.has_value());
0123   auto correctedValue = correctedRes.value();
0124   BoundVector eCorrectedBoundParams = std::get<BoundVector>(correctedValue);
0125   BoundSquareMatrix eCorrectedBoundCov =
0126       std::get<BoundSquareMatrix>(correctedValue);
0127
0128   // the loc0, phi are the same as that without correction
0129   BOOST_CHECK_EQUAL(eCorrectedBoundParams[eBoundLoc0], loc0);
0130   BOOST_CHECK_EQUAL(eCorrectedBoundParams[eBoundPhi], phi);
0131   CHECK_CLOSE_REL(eCorrectedBoundParams[eBoundLoc1], 11.2563, eps);
0132
0133   BOOST_TEST_INFO("Corrected Bound Params: \n" << eCorrectedBoundParams);
0134   BOOST_TEST_INFO("Corrected Bound Covariance: \n" << eCorrectedBoundCov);
0135   // the error for loc0 is the same as that without correction:
0136   // loc0 at end position = distance * tan(theta) * sin(phi),
0137   // dloc0/dphi = distance * tan(theta) * cos(phi) = distance,
0138   // resolution of loc0 at end position = dloc0/dphi * resLoc0 = 2.5
0139   CHECK_CLOSE_REL(eCorrectedBoundCov(eBoundLoc0, eBoundLoc0), pow(2.5, 2), eps);
0140 }