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