Analysis Software
Documentation for sPHENIX simulation software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
CorrectedTransformFreeToBoundTests.cpp
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file CorrectedTransformFreeToBoundTests.cpp
1 // This file is part of the Acts project.
2 //
3 // Copyright (C) 2022 CERN for the benefit of the Acts project
4 //
5 // This Source Code Form is subject to the terms of the Mozilla Public
6 // License, v. 2.0. If a copy of the MPL was not distributed with this
7 // file, You can obtain one at http://mozilla.org/MPL/2.0/.
8 
9 #include <boost/test/data/test_case.hpp>
10 #include <boost/test/unit_test.hpp>
11 
22 
23 #include <algorithm>
24 #include <cmath>
25 #include <memory>
26 #include <optional>
27 #include <ostream>
28 #include <string>
29 #include <tuple>
30 
31 using namespace Acts;
32 using namespace Acts::UnitLiterals;
33 
34 namespace {
35 constexpr ActsScalar eps = 0.01;
36 }
37 
38 BOOST_AUTO_TEST_CASE(CorrectedFreeToBoundTrackParameters) {
40 
41  const auto loc0 = 0.0;
42  const auto loc1 = 0.0;
43  const auto phi = 0.0;
44  const auto theta = M_PI / 4;
45  const auto qOverP = 1 / 1_GeV;
46  const auto t = 1_ns;
47 
48  const auto resLoc0 = 0.0;
49  const auto resLoc1 = 0.0;
50  const auto resPhi = 0.25;
51  const auto resTheta = 0.25;
52  const auto resQOverP = 0.01 / 1_GeV;
53  const auto resTime = 0.01_ns;
54 
55  // construct two parallel plane surfaces with normal in x direction
56  ActsScalar distance = 10_mm;
57  auto sSurface =
58  Surface::makeShared<PlaneSurface>(Vector3(0, 0, 0), Vector3::UnitX());
59  auto eSurface = Surface::makeShared<PlaneSurface>(Vector3(distance, 0, 0),
60  Vector3::UnitX());
61 
62  // the bound parameters at the starting plane
63  BoundVector sBoundParams = BoundVector::Zero();
64  sBoundParams << loc0, loc1, phi, theta, qOverP, t;
65 
66  // the bound parameters covariance at the starting plane
67  BoundSquareMatrix sBoundCov = BoundSquareMatrix::Zero();
68  sBoundCov(eBoundLoc0, eBoundLoc0) = resLoc0 * resLoc0;
69  sBoundCov(eBoundLoc1, eBoundLoc1) = resLoc1 * resLoc1;
70  sBoundCov(eBoundPhi, eBoundPhi) = resPhi * resPhi;
71  sBoundCov(eBoundTheta, eBoundTheta) = resTheta * resTheta;
72  sBoundCov(eBoundQOverP, eBoundQOverP) = resQOverP * resQOverP;
73  sBoundCov(eBoundTime, eBoundTime) = resTime * resTime;
74 
75  Vector3 dir = makeDirectionFromPhiTheta(phi, theta);
76 
77  // the intersection of the track with the end surface
78  SurfaceIntersection intersection =
79  eSurface->intersect(geoCtx, Vector3(0, 0, 0), dir, true).closest();
80  Vector3 tpos = intersection.position();
81  auto s = intersection.pathLength();
82 
83  BOOST_CHECK_EQUAL(s, distance * std::sqrt(2));
84 
85  // construct the free parameters vector
86  FreeVector eFreeParams = FreeVector::Zero();
87  eFreeParams.segment<3>(0) = tpos;
88  eFreeParams[eFreeTime] = t;
89  eFreeParams.segment<3>(4) = dir;
90  eFreeParams[eFreeQOverP] = qOverP;
91 
92  // the jacobian from local to global at the starting position
93  BoundToFreeMatrix boundToFreeJac =
94  sSurface->boundToFreeJacobian(geoCtx, sBoundParams);
95 
96  // the transport jacobian without B field
97  FreeMatrix transportJac = FreeMatrix::Identity();
98  transportJac(eFreePos0, eFreeDir0) = s;
99  transportJac(eFreePos1, eFreeDir1) = s;
100  transportJac(eFreePos2, eFreeDir2) = s;
101 
102  // the free covariance at the start position
103  FreeSquareMatrix sFreeCov =
104  boundToFreeJac * sBoundCov * boundToFreeJac.transpose();
105  // the free covariance at the end position
106  FreeSquareMatrix eFreeCov =
107  transportJac * sFreeCov * transportJac.transpose();
108 
109  // convert free parameters to bound parameters with non-linear correction
110 
111  BOOST_TEST_INFO("Transform free parameters vector onto surface "
112  << eSurface->name());
113 
114  // the corrected transformation
115  auto freeToBoundCorrection = FreeToBoundCorrection(true);
116  BOOST_CHECK(freeToBoundCorrection);
117 
118  auto transformer =
119  detail::CorrectedFreeToBoundTransformer(freeToBoundCorrection);
120  auto correctedRes = transformer(eFreeParams, eFreeCov, *eSurface, geoCtx);
121 
122  BOOST_CHECK(correctedRes.has_value());
123  auto correctedValue = correctedRes.value();
124  BoundVector eCorrectedBoundParams = std::get<BoundVector>(correctedValue);
125  BoundSquareMatrix eCorrectedBoundCov =
126  std::get<BoundSquareMatrix>(correctedValue);
127 
128  // the loc0, phi are the same as that without correction
129  BOOST_CHECK_EQUAL(eCorrectedBoundParams[eBoundLoc0], loc0);
130  BOOST_CHECK_EQUAL(eCorrectedBoundParams[eBoundPhi], phi);
131  CHECK_CLOSE_REL(eCorrectedBoundParams[eBoundLoc1], 11.2563, eps);
132 
133  BOOST_TEST_INFO("Corrected Bound Params: \n" << eCorrectedBoundParams);
134  BOOST_TEST_INFO("Corrected Bound Covariance: \n" << eCorrectedBoundCov);
135  // the error for loc0 is the same as that without correction:
136  // loc0 at end position = distance * tan(theta) * sin(phi),
137  // dloc0/dphi = distance * tan(theta) * cos(phi) = distance,
138  // resolution of loc0 at end position = dloc0/dphi * resLoc0 = 2.5
139  CHECK_CLOSE_REL(eCorrectedBoundCov(eBoundLoc0, eBoundLoc0), pow(2.5, 2), eps);
140 }