Core/Propagator/JacobianEngineTests.cpp

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 <boost/test/unit_test.hpp>
0010
0011 #include "Acts/Definitions/Algebra.hpp"
0012 #include "Acts/Definitions/TrackParametrization.hpp"
0013 #include "Acts/Geometry/GeometryContext.hpp"
0014 #include "Acts/Propagator/detail/JacobianEngine.hpp"
0015 #include "Acts/Surfaces/CurvilinearSurface.hpp"
0016 #include "Acts/Surfaces/PlaneSurface.hpp"
0017 #include "Acts/Surfaces/Surface.hpp"
0018
0019 namespace Acts::Test {
0020
0021 /// These tests do not test for a correct covariance transport but only for the
0022 /// correct conservation or modification of certain variables. A test suite for
0023 /// the numerical correctness is performed in the integration tests.
0024 ///
0025 BOOST_AUTO_TEST_CASE(jacobian_engine_to_bound) {
0026   // Create a test context
0027   GeometryContext tgContext = GeometryContext();
0028
0029   // Build a start vector
0030   Vector3 position{1., 2., 3.};
0031   double time = 4.;
0032   Vector3 direction = Vector3(5., 2., 7.).normalized();
0033   double qop = 0.125;
0034
0035   // Build a surface
0036   std::shared_ptr<PlaneSurface> pSurface =
0037       CurvilinearSurface(position, direction).planeSurface();
0038
0039   // Other rotated surface
0040   Vector3 odirection = Vector3(6., 2., 8.).normalized();
0041   std::shared_ptr<PlaneSurface> oSurface =
0042       CurvilinearSurface(position, odirection).planeSurface();
0043
0044   // The free parameter vector
0045   FreeVector freeParameters;
0046   freeParameters << position[0], position[1], position[2], time, direction[0],
0047       direction[1], direction[2], qop;
0048   // And its associated bound vector
0049   BoundVector boundParameters;
0050   boundParameters << 0., 0., VectorHelpers::phi(direction),
0051       VectorHelpers::theta(direction), qop, time;
0052
0053   // Build covariance matrices for bound and free case
0054   BoundSquareMatrix boundCovariance = 0.025 * BoundSquareMatrix::Identity();
0055   FreeSquareMatrix freeCovariance = 0.025 * FreeSquareMatrix::Identity();
0056
0057   FreeMatrix noTransportJacobian = FreeMatrix::Identity();
0058   FreeMatrix realTransportJacobian = 2 * FreeMatrix::Identity();
0059
0060   FreeToPathMatrix freeToPathDerivatives =
0061       pSurface->freeToPathDerivative(tgContext, position, direction);
0062   BoundToFreeMatrix boundToFreeJacobian =
0063       pSurface->boundToFreeJacobian(tgContext, position, direction);
0064
0065   // (1) curvilinear/bound to bound transport jacobian
0066   // a) test without actual transport into the same surface
0067   BoundMatrix b2bTransportJacobian;
0068   FreeToBoundMatrix freeToBoundJacobian;
0069   detail::boundToBoundTransportJacobian(
0070       tgContext, *pSurface, freeParameters, boundToFreeJacobian,
0071       noTransportJacobian, freeToBoundJacobian, freeToPathDerivatives,
0072       b2bTransportJacobian);
0073   BoundSquareMatrix newBoundCovariance =
0074       b2bTransportJacobian * boundCovariance * b2bTransportJacobian.transpose();
0075   BOOST_CHECK(boundCovariance.isApprox(newBoundCovariance));
0076   // b) test without actual transport but to a new surface
0077   detail::boundToBoundTransportJacobian(
0078       tgContext, *oSurface, freeParameters, boundToFreeJacobian,
0079       noTransportJacobian, freeToBoundJacobian, freeToPathDerivatives,
0080       b2bTransportJacobian);
0081   newBoundCovariance =
0082       b2bTransportJacobian * boundCovariance * b2bTransportJacobian.transpose();
0083   BOOST_CHECK(!boundCovariance.isApprox(newBoundCovariance));
0084   // c) test to "the same" surface with transport
0085   // (not really senseful, but should give a different result)
0086   detail::boundToBoundTransportJacobian(
0087       tgContext, *pSurface, freeParameters, boundToFreeJacobian,
0088       realTransportJacobian, freeToBoundJacobian, freeToPathDerivatives,
0089       b2bTransportJacobian);
0090   newBoundCovariance =
0091       b2bTransportJacobian * boundCovariance * b2bTransportJacobian.transpose();
0092   BOOST_CHECK(!boundCovariance.isApprox(newBoundCovariance));
0093
0094   FreeToBoundMatrix f2bTransportJacobian;
0095   detail::freeToBoundTransportJacobian(
0096       tgContext, *pSurface, freeParameters, noTransportJacobian,
0097       freeToPathDerivatives, f2bTransportJacobian);
0098
0099   newBoundCovariance =
0100       f2bTransportJacobian * freeCovariance * f2bTransportJacobian.transpose();
0101   BOOST_CHECK(!boundCovariance.isApprox(newBoundCovariance));
0102 }
0103
0104 /// These tests do not test for a correct covariance transport but only for the
0105 /// correct conservation or modification of certain variables. A test suite for
0106 /// the numerical correctness is performed in the integration tests.
0107 ///
0108 BOOST_AUTO_TEST_CASE(jacobian_engine_to_curvilinear) {
0109   // Create a test context
0110   GeometryContext tgContext = GeometryContext();
0111
0112   // Build a start vector
0113   Vector3 position{1., 2., 3.};
0114   Vector3 direction = Vector3(5., 2., 7.).normalized();
0115
0116   // Build a surface, starting surface for curvilinear
0117   std::shared_ptr<PlaneSurface> pSurface =
0118       CurvilinearSurface(position, direction).planeSurface();
0119
0120   // Build covariance matrices for bound and free case
0121   BoundSquareMatrix boundCovariance = 0.025 * BoundSquareMatrix::Identity();
0122   FreeSquareMatrix freeCovariance = 0.025 * FreeSquareMatrix::Identity();
0123
0124   FreeMatrix noTransportJacobian = FreeMatrix::Identity();
0125   FreeToBoundMatrix freeToBoundJacobian;
0126
0127   FreeToPathMatrix freeToPathDerivatives =
0128       pSurface->freeToPathDerivative(tgContext, position, direction);
0129   BoundToFreeMatrix boundToFreeJacobian =
0130       CurvilinearSurface(direction).boundToFreeJacobian();
0131
0132   // (1) curvilinear/bound to curvilinear transport jacobian
0133   // a) test without actual transport into the same surface
0134   BoundMatrix b2cTransportJacobian;
0135   detail::boundToCurvilinearTransportJacobian(
0136       direction, boundToFreeJacobian, noTransportJacobian, freeToBoundJacobian,
0137       freeToPathDerivatives, b2cTransportJacobian);
0138   BoundSquareMatrix newBoundCovariance =
0139       b2cTransportJacobian * boundCovariance * b2cTransportJacobian.transpose();
0140   BOOST_CHECK(boundCovariance.isApprox(newBoundCovariance));
0141   // b) test to another curvilinear frame at the same point (no transport)
0142   detail::boundToCurvilinearTransportJacobian(
0143       Vector3(4., 5., 6.).normalized(), boundToFreeJacobian,
0144       noTransportJacobian, freeToBoundJacobian, freeToPathDerivatives,
0145       b2cTransportJacobian);
0146   newBoundCovariance =
0147       b2cTransportJacobian * boundCovariance * b2cTransportJacobian.transpose();
0148   BOOST_CHECK(!boundCovariance.isApprox(newBoundCovariance));
0149
0150   FreeToBoundMatrix f2cTransportJacobian =
0151       detail::freeToCurvilinearTransportJacobian(direction, noTransportJacobian,
0152                                                  freeToPathDerivatives);
0153
0154   newBoundCovariance =
0155       f2cTransportJacobian * freeCovariance * f2cTransportJacobian.transpose();
0156   BOOST_CHECK(!boundCovariance.isApprox(newBoundCovariance));
0157 }
0158
0159 /// These tests do not test for a correct covariance transport but only for the
0160 /// correct conservation or modification of certain variables. A test suite for
0161 /// the numerical correctness is performed in the integration tests.
0162 ///
0163 BOOST_AUTO_TEST_CASE(jacobian_engine_to_free) {
0164   // Create a test context
0165   GeometryContext tgContext = GeometryContext();
0166
0167   // Build a start vector
0168   Vector3 position{1., 2., 3.};
0169   Vector3 direction = Vector3(5., 2., 7.).normalized();
0170
0171   // Build a surface, starting surface for curvilinear
0172   std::shared_ptr<PlaneSurface> pSurface =
0173       CurvilinearSurface(position, direction).planeSurface();
0174
0175   // Build covariance matrices for bound and free case
0176   BoundSquareMatrix boundCovariance = 0.025 * BoundSquareMatrix::Identity();
0177   FreeSquareMatrix freeCovariance = 0.025 * FreeSquareMatrix::Identity();
0178
0179   FreeMatrix noTransportJacobian = FreeMatrix::Identity();
0180
0181   BoundToFreeMatrix boundToFreeJacobian =
0182       pSurface->boundToFreeJacobian(tgContext, position, direction);
0183
0184   // (1) bound to free
0185   BoundToFreeMatrix b2fTransportJacobian = detail::boundToFreeTransportJacobian(
0186       boundToFreeJacobian, noTransportJacobian);
0187
0188   FreeMatrix newFreeCovariance1 =
0189       b2fTransportJacobian * boundCovariance * b2fTransportJacobian.transpose();
0190   BOOST_CHECK(!newFreeCovariance1.isApprox(freeCovariance));
0191
0192   // (2) curvilinear to free
0193   boundToFreeJacobian = CurvilinearSurface(direction).boundToFreeJacobian();
0194   BoundToFreeMatrix c2fTransportJacobian = detail::boundToFreeTransportJacobian(
0195       boundToFreeJacobian, noTransportJacobian);
0196
0197   FreeMatrix newFreeCovariance2 =
0198       c2fTransportJacobian * boundCovariance * c2fTransportJacobian.transpose();
0199   BOOST_CHECK(!newFreeCovariance2.isApprox(freeCovariance));
0200   // But those should be similar/equal
0201   BOOST_CHECK(newFreeCovariance1.isApprox(newFreeCovariance2));
0202 }
0203
0204 }  // namespace Acts::Test