Core/Utilities/VectorHelpersTests.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/Utilities/AxisDefinitions.hpp"
0013 #include "Acts/Utilities/VectorHelpers.hpp"
0014 #include "ActsTests/CommonHelpers/FloatComparisons.hpp"
0015
0016 #include <cmath>
0017 #include <limits>
0018 #include <numbers>
0019
0020 using namespace Acts;
0021 using namespace Acts::VectorHelpers;
0022
0023 namespace ActsTests {
0024
0025 BOOST_AUTO_TEST_SUITE(UtilitiesSuite)
0026
0027 BOOST_AUTO_TEST_CASE(PhiFromVector) {
0028   // Test phi calculation for Eigen vectors
0029   CHECK_CLOSE_ABS(0.0, phi(Vector2{1, 0}), 1e-6);
0030   CHECK_CLOSE_ABS(std::numbers::pi / 2, phi(Vector2{0, 1}), 1e-6);
0031   CHECK_CLOSE_ABS(std::numbers::pi, phi(Vector2{-1, 0}), 1e-6);
0032   CHECK_CLOSE_ABS(-std::numbers::pi / 2, phi(Vector2{0, -1}), 1e-6);
0033   CHECK_CLOSE_ABS(std::numbers::pi / 4, phi(Vector2{1, 1}), 1e-6);
0034
0035   // Test with 3D vectors
0036   CHECK_CLOSE_ABS(0.0, phi(Vector3{1, 0, 5}), 1e-6);
0037   CHECK_CLOSE_ABS(std::numbers::pi / 2, phi(Vector3{0, 1, -2}), 1e-6);
0038
0039   // Test with dynamic vectors
0040   ActsDynamicVector dynVec2(2);
0041   dynVec2 << 1, 1;
0042   CHECK_CLOSE_ABS(std::numbers::pi / 4, phi(dynVec2), 1e-6);
0043
0044   ActsDynamicVector dynVec3(3);
0045   dynVec3 << 0, 1, 5;
0046   CHECK_CLOSE_ABS(std::numbers::pi / 2, phi(dynVec3), 1e-6);
0047 }
0048
0049 BOOST_AUTO_TEST_CASE(PhiFromObject) {
0050   // Test phi calculation for objects with phi() method
0051   struct MockObject {
0052     double phi() const { return 1.5; }
0053   };
0054
0055   MockObject obj;
0056   CHECK_CLOSE_ABS(1.5, phi(obj), 1e-6);
0057 }
0058
0059 BOOST_AUTO_TEST_CASE(PerpFromVector) {
0060   // Test transverse radius calculation
0061   CHECK_CLOSE_ABS(1.0, perp(Vector2{1, 0}), 1e-6);
0062   CHECK_CLOSE_ABS(1.0, perp(Vector2{0, 1}), 1e-6);
0063   CHECK_CLOSE_ABS(std::sqrt(2.0), perp(Vector2{1, 1}), 1e-6);
0064   CHECK_CLOSE_ABS(5.0, perp(Vector2{3, 4}), 1e-6);
0065
0066   // Test with 3D vectors (should ignore z component)
0067   CHECK_CLOSE_ABS(1.0, perp(Vector3{1, 0, 10}), 1e-6);
0068   CHECK_CLOSE_ABS(5.0, perp(Vector3{3, 4, 100}), 1e-6);
0069
0070   // Test with dynamic vectors
0071   ActsDynamicVector dynVec(3);
0072   dynVec << 3, 4, 5;
0073   CHECK_CLOSE_ABS(5.0, perp(dynVec), 1e-6);
0074 }
0075
0076 BOOST_AUTO_TEST_CASE(ThetaFromVector) {
0077   // Test theta calculation
0078   CHECK_CLOSE_ABS(std::numbers::pi / 2, theta(Vector3{1, 0, 0}), 1e-6);
0079   CHECK_CLOSE_ABS(std::numbers::pi / 2, theta(Vector3{0, 1, 0}), 1e-6);
0080   CHECK_CLOSE_ABS(0.0, theta(Vector3{0, 0, 1}), 1e-6);
0081   CHECK_CLOSE_ABS(std::numbers::pi, theta(Vector3{0, 0, -1}), 1e-6);
0082   CHECK_CLOSE_ABS(std::numbers::pi / 4, theta(Vector3{1, 0, 1}), 1e-6);
0083   CHECK_CLOSE_ABS(3 * std::numbers::pi / 4, theta(Vector3{1, 0, -1}), 1e-6);
0084
0085   // Test with dynamic vectors
0086   ActsDynamicVector dynVec(3);
0087   dynVec << 1, 0, 1;
0088   CHECK_CLOSE_ABS(std::numbers::pi / 4, theta(dynVec), 1e-6);
0089
0090   dynVec << 0, 0, 1;
0091   CHECK_CLOSE_ABS(0.0, theta(dynVec), 1e-6);
0092 }
0093
0094 BOOST_AUTO_TEST_CASE(EtaFromVector) {
0095   // Test pseudorapidity calculation
0096   CHECK_CLOSE_ABS(0.0, eta(Vector3{1, 0, 0}), 1e-6);
0097   CHECK_CLOSE_ABS(0.0, eta(Vector3{0, 1, 0}), 1e-6);
0098
0099   // Test special case: vector parallel to z-axis
0100   BOOST_CHECK_EQUAL(eta(Vector3{0, 0, 1}),
0101                     std::numeric_limits<double>::infinity());
0102   BOOST_CHECK_EQUAL(eta(Vector3{0, 0, -1}),
0103                     -std::numeric_limits<double>::infinity());
0104
0105   // Test finite values
0106   CHECK_CLOSE_ABS(std::asinh(1.0), eta(Vector3{1, 0, 1}), 1e-6);
0107   CHECK_CLOSE_ABS(std::asinh(-1.0), eta(Vector3{1, 0, -1}), 1e-6);
0108   CHECK_CLOSE_ABS(std::asinh(2.0), eta(Vector3{1, 1, 2 * std::sqrt(2)}), 1e-6);
0109
0110   // Test with dynamic vectors
0111   ActsDynamicVector dynVec(3);
0112   dynVec << 1, 0, 1;
0113   CHECK_CLOSE_ABS(std::asinh(1.0), eta(dynVec), 1e-6);
0114
0115   dynVec << 0, 0, 1;
0116   BOOST_CHECK_EQUAL(eta(dynVec), std::numeric_limits<double>::infinity());
0117
0118   dynVec << 0, 0, -1;
0119   BOOST_CHECK_EQUAL(eta(dynVec), -std::numeric_limits<double>::infinity());
0120 }
0121
0122 BOOST_AUTO_TEST_CASE(EvaluateTrigonomics) {
0123   // Test trigonometric evaluation
0124   Vector3 dir{1, 0, 0};  // pointing in x direction
0125   auto trig = evaluateTrigonomics(dir);
0126
0127   CHECK_CLOSE_ABS(1.0, trig[0], 1e-6);  // cos(phi)
0128   CHECK_CLOSE_ABS(0.0, trig[1], 1e-6);  // sin(phi)
0129   CHECK_CLOSE_ABS(0.0, trig[2], 1e-6);  // cos(theta)
0130   CHECK_CLOSE_ABS(1.0, trig[3], 1e-6);  // sin(theta)
0131
0132   // Test with y direction
0133   dir = Vector3{0, 1, 0};
0134   trig = evaluateTrigonomics(dir);
0135
0136   CHECK_CLOSE_ABS(0.0, trig[0], 1e-6);  // cos(phi)
0137   CHECK_CLOSE_ABS(1.0, trig[1], 1e-6);  // sin(phi)
0138   CHECK_CLOSE_ABS(0.0, trig[2], 1e-6);  // cos(theta)
0139   CHECK_CLOSE_ABS(1.0, trig[3], 1e-6);  // sin(theta)
0140
0141   // Test with diagonal direction (avoid z-axis singularity)
0142   dir = Vector3{1, 1, 1};
0143   dir.normalize();
0144   trig = evaluateTrigonomics(dir);
0145
0146   // For normalized {1,1,1}, cos(theta) = 1/sqrt(3), sin(theta) = sqrt(2/3)
0147   // cos(phi) = sin(phi) = 1/sqrt(2) (45 degree angle in xy plane)
0148   CHECK_CLOSE_ABS(1.0 / std::sqrt(2), trig[0], 1e-6);    // cos(phi)
0149   CHECK_CLOSE_ABS(1.0 / std::sqrt(2), trig[1], 1e-6);    // sin(phi)
0150   CHECK_CLOSE_ABS(1.0 / std::sqrt(3), trig[2], 1e-6);    // cos(theta)
0151   CHECK_CLOSE_ABS(std::sqrt(2.0 / 3.0), trig[3], 1e-6);  // sin(theta)
0152 }
0153
0154 BOOST_AUTO_TEST_CASE(CastAxisDirection) {
0155   using enum AxisDirection;
0156   Vector3 pos{3, 4, 5};
0157
0158   // Test all axis directions
0159   CHECK_CLOSE_ABS(3.0, cast(pos, AxisX), 1e-6);
0160   CHECK_CLOSE_ABS(4.0, cast(pos, AxisY), 1e-6);
0161   CHECK_CLOSE_ABS(5.0, cast(pos, AxisZ), 1e-6);
0162   CHECK_CLOSE_ABS(5.0, cast(pos, AxisR), 1e-6);  // sqrt(3²+4²)
0163   CHECK_CLOSE_ABS(std::atan2(4, 3), cast(pos, AxisPhi), 1e-6);
0164   CHECK_CLOSE_ABS(5.0 * std::atan2(4, 3), cast(pos, AxisRPhi), 1e-6);
0165   CHECK_CLOSE_ABS(std::atan2(5, 5), cast(pos, AxisTheta), 1e-6);
0166   CHECK_CLOSE_ABS(std::asinh(1.0), cast(pos, AxisEta), 1e-6);
0167   CHECK_CLOSE_ABS(std::sqrt(50), cast(pos, AxisMag), 1e-6);
0168 }
0169
0170 BOOST_AUTO_TEST_CASE(CrossProduct) {
0171   SquareMatrix3 m;
0172   m << 1, 0, 0, 0, 1, 0, 0, 0, 1;
0173
0174   Vector3 v{1, 0, 0};
0175
0176   SquareMatrix3 result = cross(m, v);
0177
0178   // Check that each column is the cross product of the corresponding column of
0179   // m with v
0180   Vector3 expected_col0 = m.col(0).cross(v);  // [1,0,0] x [1,0,0] = [0,0,0]
0181   Vector3 expected_col1 = m.col(1).cross(v);  // [0,1,0] x [1,0,0] = [0,0,-1]
0182   Vector3 expected_col2 = m.col(2).cross(v);  // [0,0,1] x [1,0,0] = [0,1,0]
0183
0184   CHECK_CLOSE_ABS(expected_col0[0], result.col(0)[0], 1e-6);
0185   CHECK_CLOSE_ABS(expected_col0[1], result.col(0)[1], 1e-6);
0186   CHECK_CLOSE_ABS(expected_col0[2], result.col(0)[2], 1e-6);
0187
0188   CHECK_CLOSE_ABS(expected_col1[0], result.col(1)[0], 1e-6);
0189   CHECK_CLOSE_ABS(expected_col1[1], result.col(1)[1], 1e-6);
0190   CHECK_CLOSE_ABS(expected_col1[2], result.col(1)[2], 1e-6);
0191
0192   CHECK_CLOSE_ABS(expected_col2[0], result.col(2)[0], 1e-6);
0193   CHECK_CLOSE_ABS(expected_col2[1], result.col(2)[1], 1e-6);
0194   CHECK_CLOSE_ABS(expected_col2[2], result.col(2)[2], 1e-6);
0195 }
0196
0197 BOOST_AUTO_TEST_CASE(PositionFromVector4) {
0198   Vector4 pos4{1, 2, 3, 4};
0199   auto pos3 = position(pos4);
0200
0201   CHECK_CLOSE_ABS(1.0, pos3[0], 1e-6);
0202   CHECK_CLOSE_ABS(2.0, pos3[1], 1e-6);
0203   CHECK_CLOSE_ABS(3.0, pos3[2], 1e-6);
0204 }
0205
0206 BOOST_AUTO_TEST_CASE(PositionFromFreeVector) {
0207   FreeVector params = FreeVector::Zero();
0208   params[eFreePos0] = 1;
0209   params[eFreePos1] = 2;
0210   params[eFreePos2] = 3;
0211
0212   auto pos3 = position(params);
0213
0214   CHECK_CLOSE_ABS(1.0, pos3[0], 1e-6);
0215   CHECK_CLOSE_ABS(2.0, pos3[1], 1e-6);
0216   CHECK_CLOSE_ABS(3.0, pos3[2], 1e-6);
0217 }
0218
0219 BOOST_AUTO_TEST_CASE(MakeVector4) {
0220   Vector3 vec3{1, 2, 3};
0221   double w = 4.0;
0222
0223   auto vec4 = makeVector4(vec3, w);
0224
0225   CHECK_CLOSE_ABS(1.0, vec4[ePos0], 1e-6);
0226   CHECK_CLOSE_ABS(2.0, vec4[ePos1], 1e-6);
0227   CHECK_CLOSE_ABS(3.0, vec4[ePos2], 1e-6);
0228   CHECK_CLOSE_ABS(4.0, vec4[eTime], 1e-6);
0229 }
0230
0231 BOOST_AUTO_TEST_CASE(IncidentAngles) {
0232   Vector3 direction{1, 1, 1};
0233   direction.normalize();
0234
0235   // Identity rotation (global == local)
0236   RotationMatrix3 globalToLocal = RotationMatrix3::Identity();
0237
0238   auto angles = incidentAngles(direction, globalToLocal);
0239
0240   // In local frame, direction is still {1,1,1}/sqrt(3)
0241   // phi = atan2(z, x) = atan2(1/sqrt(3), 1/sqrt(3)) = pi/4
0242   // theta = atan2(z, y) = atan2(1/sqrt(3), 1/sqrt(3)) = pi/4
0243   CHECK_CLOSE_ABS(std::numbers::pi / 4, angles.first, 1e-6);
0244   CHECK_CLOSE_ABS(std::numbers::pi / 4, angles.second, 1e-6);
0245
0246   // Test with rotation around z-axis by 90 degrees
0247   RotationMatrix3 rotZ;
0248   rotZ << 0, 1, 0, -1, 0, 0, 0, 0, 1;
0249
0250   angles = incidentAngles(Vector3{1, 0, 0}, rotZ);
0251   // After rotation: {1,0,0} -> {0,-1,0}
0252   // phi = atan2(0, 0) = 0, theta = atan2(0, -1) = pi (atan2 of 0,-1)
0253   CHECK_CLOSE_ABS(0.0, angles.first, 1e-6);
0254   CHECK_CLOSE_ABS(std::numbers::pi, angles.second, 1e-6);
0255 }
0256
0257 BOOST_AUTO_TEST_SUITE_END()
0258
0259 }  // namespace ActsTests