Core/Surfaces/DiscSurfaceTests.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/data/test_case.hpp>
0010 #include <boost/test/tools/output_test_stream.hpp>
0011 #include <boost/test/unit_test.hpp>
0012
0013 #include "Acts/Definitions/Algebra.hpp"
0014 #include "Acts/Definitions/Alignment.hpp"
0015 #include "Acts/Definitions/Tolerance.hpp"
0016 #include "Acts/Definitions/Units.hpp"
0017 #include "Acts/Geometry/Extent.hpp"
0018 #include "Acts/Geometry/GeometryContext.hpp"
0019 #include "Acts/Geometry/Polyhedron.hpp"
0020 #include "Acts/Surfaces/AnnulusBounds.hpp"
0021 #include "Acts/Surfaces/DiscSurface.hpp"
0022 #include "Acts/Surfaces/RadialBounds.hpp"
0023 #include "Acts/Surfaces/Surface.hpp"
0024 #include "Acts/Surfaces/SurfaceBounds.hpp"
0025 #include "Acts/Surfaces/SurfaceMergingException.hpp"
0026 #include "Acts/Tests/CommonHelpers/DetectorElementStub.hpp"
0027 #include "Acts/Tests/CommonHelpers/FloatComparisons.hpp"
0028 #include "Acts/Utilities/BinningType.hpp"
0029 #include "Acts/Utilities/Intersection.hpp"
0030 #include "Acts/Utilities/Result.hpp"
0031 #include "Acts/Utilities/ThrowAssert.hpp"
0032 #include "Acts/Utilities/detail/periodic.hpp"
0033
0034 #include <cmath>
0035 #include <memory>
0036 #include <numbers>
0037 #include <ostream>
0038 #include <string>
0039
0040 using namespace Acts::UnitLiterals;
0041
0042 namespace Acts::Test {
0043 // Create a test context
0044 GeometryContext tgContext = GeometryContext();
0045 auto logger = Acts::getDefaultLogger("UnitTests", Acts::Logging::VERBOSE);
0046
0047 BOOST_AUTO_TEST_SUITE(Surfaces)
0048 /// Unit tests for creating DiscSurface object
0049 BOOST_AUTO_TEST_CASE(DiscSurfaceConstruction) {
0050   /// Test default construction
0051   // default construction is deleted
0052
0053   const double rMin = 1.;
0054   const double rMax = 5.;
0055   const double halfPhiSector = std::numbers::pi / 8.;
0056
0057   /// Test DiscSurface constructor with default halfPhiSector
0058   BOOST_CHECK_NO_THROW(
0059       Surface::makeShared<DiscSurface>(Transform3::Identity(), rMin, rMax));
0060
0061   /// Test DiscSurface constructor with a transform specified
0062   Translation3 translation{0., 1., 2.};
0063   auto pTransform = Transform3(translation);
0064   BOOST_CHECK_NO_THROW(
0065       Surface::makeShared<DiscSurface>(pTransform, rMin, rMax, halfPhiSector));
0066
0067   /// Copy constructed DiscSurface
0068   auto anotherDiscSurface =
0069       Surface::makeShared<DiscSurface>(pTransform, rMin, rMax, halfPhiSector);
0070   // N.B. Just using
0071   // BOOST_CHECK_NO_THROW(Surface::makeShared<DiscSurface>(anotherDiscSurface))
0072   // tries to call the (deleted) default constructor.
0073   auto copiedSurface = Surface::makeShared<DiscSurface>(*anotherDiscSurface);
0074   BOOST_TEST_MESSAGE("Copy constructed DiscSurface ok");
0075
0076   /// Copied and transformed DiscSurface
0077   BOOST_CHECK_NO_THROW(Surface::makeShared<DiscSurface>(
0078       tgContext, *anotherDiscSurface, pTransform));
0079
0080   /// Construct with nullptr bounds
0081   DetectorElementStub detElem;
0082   BOOST_CHECK_THROW(
0083       auto nullBounds = Surface::makeShared<DiscSurface>(nullptr, detElem),
0084       AssertionFailureException);
0085 }
0086
0087 /// Unit tests of all named methods
0088 BOOST_AUTO_TEST_CASE(DiscSurfaceProperties) {
0089   const double rMin = 1.;
0090   const double rMax = 5.;
0091   const double halfPhiSector = std::numbers::pi / 8.;
0092
0093   const Vector3 origin3D{0, 0, 0};
0094
0095   auto discSurfaceObject = Surface::makeShared<DiscSurface>(
0096       Transform3::Identity(), rMin, rMax, halfPhiSector);
0097
0098   /// Test type
0099   BOOST_CHECK_EQUAL(discSurfaceObject->type(), Surface::Disc);
0100
0101   /// Test normal, no local position specified
0102   Vector3 zAxis{0, 0, 1};
0103   BOOST_CHECK_EQUAL(discSurfaceObject->normal(tgContext), zAxis);
0104
0105   /// Test normal, local position specified
0106   Vector2 lpos(2., 0.05);
0107   BOOST_CHECK_EQUAL(discSurfaceObject->normal(tgContext, lpos), zAxis);
0108
0109   /// Test referencePosition
0110   BOOST_CHECK_EQUAL(
0111       discSurfaceObject->referencePosition(tgContext, AxisDirection::AxisRPhi),
0112       origin3D);
0113
0114   /// Test bounds
0115   BOOST_CHECK_EQUAL(discSurfaceObject->bounds().type(), SurfaceBounds::eDisc);
0116
0117   Vector3 ignoredMomentum{0., 0., 0.};
0118   /// Test isOnSurface()
0119   Vector3 point3DNotInSector{0., 1.2, 0};
0120   Vector3 point3DOnSurface{1.2, 0., 0};
0121   BOOST_CHECK(!discSurfaceObject->isOnSurface(tgContext, point3DNotInSector,
0122                                               ignoredMomentum,
0123                                               BoundaryTolerance::None()));
0124   BOOST_CHECK(!discSurfaceObject->isOnSurface(tgContext, point3DNotInSector,
0125                                               BoundaryTolerance::None()));
0126   BOOST_CHECK(discSurfaceObject->isOnSurface(
0127       tgContext, point3DOnSurface, ignoredMomentum, BoundaryTolerance::None()));
0128   BOOST_CHECK(discSurfaceObject->isOnSurface(tgContext, point3DOnSurface,
0129                                              BoundaryTolerance::None()));
0130
0131   /// Test localToGlobal
0132   Vector3 returnedPosition{10.9, 8.7, 6.5};
0133   Vector3 expectedPosition{1.2, 0, 0};
0134   Vector2 rPhiOnDisc{1.2, 0.};
0135   Vector2 rPhiNotInSector{
0136       1.2, std::numbers::pi};  // outside sector at Phi=0, +/- pi/8
0137   returnedPosition =
0138       discSurfaceObject->localToGlobal(tgContext, rPhiOnDisc, ignoredMomentum);
0139   CHECK_CLOSE_ABS(returnedPosition, expectedPosition, 1e-6);
0140
0141   returnedPosition = discSurfaceObject->localToGlobal(
0142       tgContext, rPhiNotInSector, ignoredMomentum);
0143   Vector3 expectedNonPosition{-1.2, 0, 0};
0144   CHECK_CLOSE_ABS(returnedPosition, expectedNonPosition, 1e-6);
0145
0146   /// Test globalToLocal
0147   Vector2 returnedLocalPosition{33., 44.};
0148   Vector2 expectedLocalPosition{1.2, 0.};
0149   returnedLocalPosition =
0150       discSurfaceObject
0151           ->globalToLocal(tgContext, point3DOnSurface, ignoredMomentum)
0152           .value();
0153   CHECK_CLOSE_ABS(returnedLocalPosition, expectedLocalPosition, 1e-6);
0154
0155   // Global to local does not check inside bounds
0156   returnedLocalPosition =
0157       discSurfaceObject
0158           ->globalToLocal(tgContext, point3DNotInSector, ignoredMomentum)
0159           .value();
0160
0161   Vector3 pointOutsideR{0., 100., 0};
0162   returnedLocalPosition =
0163       discSurfaceObject
0164           ->globalToLocal(tgContext, pointOutsideR, ignoredMomentum)
0165           .value();
0166
0167   /// Test localPolarToCartesian
0168   Vector2 rPhi1_1{std::numbers::sqrt2, std::numbers::pi / 4.};
0169   Vector2 cartesian1_1{1., 1.};
0170   CHECK_CLOSE_REL(discSurfaceObject->localPolarToCartesian(rPhi1_1),
0171                   cartesian1_1, 1e-6);
0172
0173   /// Test localCartesianToPolar
0174   CHECK_CLOSE_REL(discSurfaceObject->localCartesianToPolar(cartesian1_1),
0175                   rPhi1_1, 1e-6);
0176
0177   /// Test localPolarToLocalCartesian
0178   CHECK_CLOSE_REL(discSurfaceObject->localPolarToLocalCartesian(rPhi1_1),
0179                   cartesian1_1, 1e-6);
0180
0181   /// Test localCartesianToGlobal
0182   Vector3 cartesian3D1_1{1., 1., 0.};
0183   CHECK_CLOSE_ABS(
0184       discSurfaceObject->localCartesianToGlobal(tgContext, cartesian1_1),
0185       cartesian3D1_1, 1e-6);
0186
0187   /// Test globalToLocalCartesian
0188   CHECK_CLOSE_REL(
0189       discSurfaceObject->globalToLocalCartesian(tgContext, cartesian3D1_1),
0190       cartesian1_1, 1e-6);
0191
0192   /// Test pathCorrection
0193   double projected3DMomentum = std::numbers::sqrt3 * 1.e6;
0194   Vector3 momentum{projected3DMomentum, projected3DMomentum,
0195                    projected3DMomentum};
0196   Vector3 ignoredPosition = discSurfaceObject->center(tgContext);
0197   CHECK_CLOSE_REL(discSurfaceObject->pathCorrection(tgContext, ignoredPosition,
0198                                                     momentum.normalized()),
0199                   std::numbers::sqrt3, 0.01);
0200
0201   /// intersection test
0202   Vector3 globalPosition{1.2, 0., -10.};
0203   Vector3 direction{0., 0., 1.};  // must be normalised
0204   Vector3 expected{1.2, 0., 0.};
0205
0206   // intersect is a struct of (Vector3) position, pathLength, distance and
0207   // (bool) valid, it's contained in a Surface intersection
0208   auto sfIntersection = discSurfaceObject
0209                             ->intersect(tgContext, globalPosition, direction,
0210                                         BoundaryTolerance::Infinite())
0211                             .closest();
0212   Intersection3D expectedIntersect{Vector3{1.2, 0., 0.}, 10.,
0213                                    IntersectionStatus::reachable};
0214   BOOST_CHECK(sfIntersection.isValid());
0215   CHECK_CLOSE_ABS(sfIntersection.position(), expectedIntersect.position(),
0216                   1e-9);
0217   CHECK_CLOSE_ABS(sfIntersection.pathLength(), expectedIntersect.pathLength(),
0218                   1e-9);
0219   BOOST_CHECK_EQUAL(sfIntersection.object(), discSurfaceObject.get());
0220
0221   /// Test name
0222   boost::test_tools::output_test_stream nameOuput;
0223   nameOuput << discSurfaceObject->name();
0224   BOOST_CHECK(nameOuput.is_equal("Acts::DiscSurface"));
0225 }
0226
0227 /// Unit test for testing DiscSurface assignment and equality
0228 BOOST_AUTO_TEST_CASE(DiscSurfaceAssignment) {
0229   const double rMin = 1.;
0230   const double rMax = 5.;
0231   const double halfPhiSector = std::numbers::pi / 8.;
0232
0233   auto discSurfaceObject = Surface::makeShared<DiscSurface>(
0234       Transform3::Identity(), rMin, rMax, halfPhiSector);
0235   auto assignedDisc =
0236       Surface::makeShared<DiscSurface>(Transform3::Identity(), 2.2, 4.4, 0.07);
0237
0238   BOOST_CHECK_NO_THROW(*assignedDisc = *discSurfaceObject);
0239   BOOST_CHECK((*assignedDisc) == (*discSurfaceObject));
0240 }
0241
0242 /// Unit test for testing DiscSurface assignment and equality
0243 BOOST_AUTO_TEST_CASE(DiscSurfaceExtent) {
0244   const double rMin = 1.;
0245   const double rMax = 5.;
0246
0247   auto pDisc =
0248       Surface::makeShared<DiscSurface>(Transform3::Identity(), 0., rMax);
0249   auto pDiscExtent = pDisc->polyhedronRepresentation(tgContext, 1).extent();
0250
0251   CHECK_CLOSE_ABS(0., pDiscExtent.min(AxisDirection::AxisZ),
0252                   s_onSurfaceTolerance);
0253   CHECK_CLOSE_ABS(0., pDiscExtent.max(AxisDirection::AxisZ),
0254                   s_onSurfaceTolerance);
0255   CHECK_CLOSE_ABS(0., pDiscExtent.min(AxisDirection::AxisR),
0256                   s_onSurfaceTolerance);
0257   CHECK_CLOSE_ABS(rMax, pDiscExtent.max(AxisDirection::AxisR),
0258                   s_onSurfaceTolerance);
0259   CHECK_CLOSE_ABS(-rMax, pDiscExtent.min(AxisDirection::AxisX),
0260                   s_onSurfaceTolerance);
0261   CHECK_CLOSE_ABS(rMax, pDiscExtent.max(AxisDirection::AxisX),
0262                   s_onSurfaceTolerance);
0263   CHECK_CLOSE_ABS(-rMax, pDiscExtent.min(AxisDirection::AxisY),
0264                   s_onSurfaceTolerance);
0265   CHECK_CLOSE_ABS(rMax, pDiscExtent.max(AxisDirection::AxisY),
0266                   s_onSurfaceTolerance);
0267   CHECK_CLOSE_ABS(-std::numbers::pi, pDiscExtent.min(AxisDirection::AxisPhi),
0268                   s_onSurfaceTolerance);
0269   CHECK_CLOSE_ABS(std::numbers::pi, pDiscExtent.max(AxisDirection::AxisPhi),
0270                   s_onSurfaceTolerance);
0271
0272   auto pRing =
0273       Surface::makeShared<DiscSurface>(Transform3::Identity(), rMin, rMax);
0274   auto pRingExtent = pRing->polyhedronRepresentation(tgContext, 1).extent();
0275
0276   CHECK_CLOSE_ABS(0., pRingExtent.min(AxisDirection::AxisZ),
0277                   s_onSurfaceTolerance);
0278   CHECK_CLOSE_ABS(0., pRingExtent.max(AxisDirection::AxisZ),
0279                   s_onSurfaceTolerance);
0280   CHECK_CLOSE_ABS(rMin, pRingExtent.min(AxisDirection::AxisR),
0281                   s_onSurfaceTolerance);
0282   CHECK_CLOSE_ABS(rMax, pRingExtent.max(AxisDirection::AxisR),
0283                   s_onSurfaceTolerance);
0284   CHECK_CLOSE_ABS(-rMax, pRingExtent.min(AxisDirection::AxisX),
0285                   s_onSurfaceTolerance);
0286   CHECK_CLOSE_ABS(rMax, pRingExtent.max(AxisDirection::AxisX),
0287                   s_onSurfaceTolerance);
0288   CHECK_CLOSE_ABS(-rMax, pRingExtent.min(AxisDirection::AxisY),
0289                   s_onSurfaceTolerance);
0290   CHECK_CLOSE_ABS(rMax, pRingExtent.max(AxisDirection::AxisY),
0291                   s_onSurfaceTolerance);
0292 }
0293
0294 /// Unit test for testing DiscSurface alignment derivatives
0295 BOOST_AUTO_TEST_CASE(DiscSurfaceAlignment) {
0296   Translation3 translation{0., 1., 2.};
0297   Transform3 transform(translation);
0298   const double rMin = 1.;
0299   const double rMax = 5.;
0300   const double halfPhiSector = std::numbers::pi / 8.;
0301
0302   auto discSurfaceObject =
0303       Surface::makeShared<DiscSurface>(transform, rMin, rMax, halfPhiSector);
0304
0305   const auto& rotation = transform.rotation();
0306   // The local frame z axis
0307   const Vector3 localZAxis = rotation.col(2);
0308   // Check the local z axis is aligned to global z axis
0309   CHECK_CLOSE_ABS(localZAxis, Vector3(0., 0., 1.), 1e-15);
0310
0311   // Define the track (global) position and direction
0312   Vector3 globalPosition{0, 4, 2};
0313   Vector3 momentum{0, 0, 1};
0314   Vector3 direction = momentum.normalized();
0315
0316   // (a) Test the derivative of path length w.r.t. alignment parameters
0317   const AlignmentToPathMatrix& alignToPath =
0318       discSurfaceObject->alignmentToPathDerivative(tgContext, globalPosition,
0319                                                    direction);
0320   // The expected results
0321   AlignmentToPathMatrix expAlignToPath = AlignmentToPathMatrix::Zero();
0322   expAlignToPath << 0, 0, 1, 3, 0, 0;
0323   // Check if the calculated derivative is as expected
0324   CHECK_CLOSE_ABS(alignToPath, expAlignToPath, 1e-10);
0325
0326   // (b) Test the derivative of bound track parameters local position w.r.t.
0327   // position in local 3D Cartesian coordinates
0328   const auto& loc3DToLocBound =
0329       discSurfaceObject->localCartesianToBoundLocalDerivative(tgContext,
0330                                                               globalPosition);
0331   // Check if the result is as expected
0332   ActsMatrix<2, 3> expLoc3DToLocBound = ActsMatrix<2, 3>::Zero();
0333   expLoc3DToLocBound << 0, 1, 0, -1. / 3, 0, 0;
0334   CHECK_CLOSE_ABS(loc3DToLocBound, expLoc3DToLocBound, 1e-10);
0335 }
0336
0337 BOOST_AUTO_TEST_CASE(DiscSurfaceBinningPosition) {
0338   using namespace Acts::UnitLiterals;
0339   Vector3 s{5_mm, 7_mm, 10_cm};
0340   Transform3 trf;
0341   trf = Translation3(s) * AngleAxis3{0.5, Vector3::UnitZ()};
0342
0343   double minR = 300;
0344   double maxR = 330;
0345
0346   {
0347     // Radial Bounds
0348     auto bounds =
0349         std::make_shared<RadialBounds>(minR, maxR, std::numbers::pi / 8, 0.1);
0350     auto disc = Acts::Surface::makeShared<Acts::DiscSurface>(trf, bounds);
0351
0352     Vector3 bp = disc->referencePosition(tgContext, AxisDirection::AxisR);
0353     double r = (bounds->rMax() + bounds->rMin()) / 2.0;
0354     double phi = bounds->get(RadialBounds::eAveragePhi);
0355     Vector3 exp = Vector3{r * std::cos(phi), r * std::sin(phi), 0};
0356     exp = trf * exp;
0357
0358     BOOST_CHECK_EQUAL(bp, exp);
0359     BOOST_CHECK_EQUAL(
0360         disc->referencePositionValue(tgContext, AxisDirection::AxisR),
0361         VectorHelpers::perp(exp));
0362
0363     bp = disc->referencePosition(tgContext, AxisDirection::AxisPhi);
0364     BOOST_CHECK_EQUAL(bp, exp);
0365     BOOST_CHECK_EQUAL(
0366         disc->referencePositionValue(tgContext, AxisDirection::AxisPhi),
0367         VectorHelpers::phi(exp));
0368
0369     for (auto b :
0370          {AxisDirection::AxisX, AxisDirection::AxisY, AxisDirection::AxisZ,
0371           AxisDirection::AxisEta, AxisDirection::AxisRPhi,
0372           AxisDirection::AxisTheta, AxisDirection::AxisMag}) {
0373       BOOST_TEST_CONTEXT("binValue: " << b) {
0374         BOOST_CHECK_EQUAL(disc->referencePosition(tgContext, b),
0375                           disc->center(tgContext));
0376       }
0377     }
0378   }
0379
0380   {
0381     // Annulus Bounds
0382     double minPhiRel = -0.3;
0383     double maxPhiRel = 0.2;
0384     Vector2 origin{5_mm, 5_mm};
0385     auto bounds = std::make_shared<AnnulusBounds>(minR, maxR, minPhiRel,
0386                                                   maxPhiRel, origin);
0387
0388     auto disc = Acts::Surface::makeShared<Acts::DiscSurface>(trf, bounds);
0389
0390     Vector3 bp = disc->referencePosition(tgContext, AxisDirection::AxisR);
0391     double r = (bounds->rMax() + bounds->rMin()) / 2.0;
0392     double phi = bounds->get(AnnulusBounds::eAveragePhi);
0393     Vector3 exp = Vector3{r * std::cos(phi), r * std::sin(phi), 0};
0394     exp = trf * exp;
0395
0396     BOOST_CHECK_EQUAL(bp, exp);
0397
0398     bp = disc->referencePosition(tgContext, AxisDirection::AxisPhi);
0399     BOOST_CHECK_EQUAL(bp, exp);
0400
0401     for (auto b :
0402          {AxisDirection::AxisX, AxisDirection::AxisY, AxisDirection::AxisZ,
0403           AxisDirection::AxisEta, AxisDirection::AxisRPhi,
0404           AxisDirection::AxisTheta, AxisDirection::AxisMag}) {
0405       BOOST_TEST_CONTEXT("binValue: " << b) {
0406         BOOST_CHECK_EQUAL(disc->referencePosition(tgContext, b),
0407                           disc->center(tgContext));
0408       }
0409     }
0410   }
0411 }
0412
0413 BOOST_AUTO_TEST_SUITE(DiscSurfaceMerging)
0414
0415 namespace {
0416 std::shared_ptr<DiscSurface> makeDisc(const Transform3& transform, double rMin,
0417                                       double rMax,
0418                                       double halfPhi = std::numbers::pi,
0419                                       double avgPhi = 0) {
0420   return Surface::makeShared<DiscSurface>(
0421       transform, std::make_shared<RadialBounds>(rMin, rMax, halfPhi, avgPhi));
0422 }
0423
0424 }  // namespace
0425
0426 BOOST_AUTO_TEST_CASE(IncompatibleBounds) {
0427   Logging::ScopedFailureThreshold ft{Logging::FATAL};
0428   Transform3 base = Transform3::Identity();
0429   auto discRadial = makeDisc(base, 30_mm, 100_mm);
0430   auto discTrap =
0431       Surface::makeShared<DiscSurface>(base, 20_mm, 40_mm, 100_mm, 150_mm);
0432   auto discTrap2 =
0433       Surface::makeShared<DiscSurface>(base, 20_mm, 40_mm, 30_mm, 100_mm);
0434
0435   BOOST_CHECK_THROW(
0436       discRadial->mergedWith(*discTrap, AxisDirection::AxisR, false, *logger),
0437
0438       SurfaceMergingException);
0439
0440   BOOST_CHECK_THROW(
0441       discTrap2->mergedWith(*discTrap, AxisDirection::AxisR, false, *logger),
0442       SurfaceMergingException);
0443 }
0444
0445 BOOST_AUTO_TEST_CASE(InvalidDetectorElement) {
0446   DetectorElementStub detElem;
0447
0448   auto bounds1 = std::make_shared<RadialBounds>(30_mm, 100_mm);
0449   auto disc1 = Surface::makeShared<DiscSurface>(bounds1, detElem);
0450
0451   auto bounds2 = std::make_shared<RadialBounds>(100_mm, 150_mm);
0452   auto disc2 = Surface::makeShared<DiscSurface>(bounds2, detElem);
0453
0454   BOOST_CHECK_THROW(
0455       disc1->mergedWith(*disc2, AxisDirection::AxisR, false, *logger),
0456       SurfaceMergingException);
0457 }
0458
0459 BOOST_DATA_TEST_CASE(IncompatibleRDirection,
0460                      (boost::unit_test::data::xrange(-135, 180, 45) *
0461                       boost::unit_test::data::make(Vector3{0_mm, 0_mm, 0_mm},
0462                                                    Vector3{20_mm, 0_mm, 0_mm},
0463                                                    Vector3{0_mm, 20_mm, 0_mm},
0464                                                    Vector3{20_mm, 20_mm, 0_mm},
0465                                                    Vector3{0_mm, 0_mm, 20_mm})),
0466                      angle, offset) {
0467   Logging::ScopedFailureThreshold ft{Logging::FATAL};
0468
0469   Transform3 base =
0470       AngleAxis3(angle * 1_degree, Vector3::UnitX()) * Translation3(offset);
0471
0472   auto disc = makeDisc(base, 30_mm, 100_mm);
0473
0474   // Disc with overlap in r
0475   auto discOverlap = makeDisc(base, 90_mm, 150_mm);
0476   BOOST_CHECK_THROW(disc->mergedWith(*discOverlap, Acts::AxisDirection::AxisR,
0477                                      false, *logger),
0478                     SurfaceMergingException);
0479
0480   // Disc with gap in r
0481   auto discGap = makeDisc(base, 110_mm, 150_mm);
0482   BOOST_CHECK_THROW(
0483       disc->mergedWith(*discGap, Acts::AxisDirection::AxisR, false, *logger),
0484       SurfaceMergingException);
0485
0486   auto discShiftedZ = Surface::makeShared<DiscSurface>(
0487       base * Translation3{Vector3::UnitZ() * 10_mm}, 100_mm, 150_mm);
0488   BOOST_CHECK_THROW(disc->mergedWith(*discShiftedZ, Acts::AxisDirection::AxisR,
0489                                      false, *logger),
0490                     SurfaceMergingException);
0491
0492   auto discShiftedXy = makeDisc(
0493       base * Translation3{Vector3{1_mm, 2_mm, 200_mm}}, 100_mm, 150_mm);
0494   BOOST_CHECK_THROW(disc->mergedWith(*discShiftedXy, Acts::AxisDirection::AxisZ,
0495                                      false, *logger),
0496                     SurfaceMergingException);
0497
0498   auto discRotatedZ = makeDisc(base * AngleAxis3{10_degree, Vector3::UnitZ()},
0499                                100_mm, 150_mm, 60_degree, 0_degree);
0500   BOOST_CHECK_THROW(disc->mergedWith(*discRotatedZ, Acts::AxisDirection::AxisR,
0501                                      false, *logger),
0502                     SurfaceMergingException);
0503
0504   auto discRotatedX =
0505       makeDisc(base * AngleAxis3{10_degree, Vector3::UnitX()}, 100_mm, 150_mm);
0506   BOOST_CHECK_THROW(disc->mergedWith(*discRotatedX, Acts::AxisDirection::AxisR,
0507                                      false, *logger),
0508                     SurfaceMergingException);
0509
0510   // Test not same phi sector
0511   auto discPhi1 = makeDisc(base, 30_mm, 100_mm, 10_degree, 40_degree);
0512   auto discPhi2 = makeDisc(base, 100_mm, 160_mm, 20_degree, 40_degree);
0513   auto discPhi3 = makeDisc(base, 100_mm, 160_mm, 10_degree, 50_degree);
0514   BOOST_CHECK_THROW(
0515       discPhi1->mergedWith(*discPhi2, AxisDirection::AxisR, false, *logger),
0516       SurfaceMergingException);
0517
0518   BOOST_CHECK_THROW(
0519       discPhi1->mergedWith(*discPhi3, AxisDirection::AxisR, false, *logger),
0520       SurfaceMergingException);
0521 }
0522
0523 BOOST_DATA_TEST_CASE(RDirection,
0524                      (boost::unit_test::data::xrange(-135, 180, 45) *
0525                       boost::unit_test::data::make(Vector3{0_mm, 0_mm, 0_mm},
0526                                                    Vector3{20_mm, 0_mm, 0_mm},
0527                                                    Vector3{0_mm, 20_mm, 0_mm},
0528                                                    Vector3{20_mm, 20_mm, 0_mm},
0529                                                    Vector3{0_mm, 0_mm, 20_mm})),
0530                      angle, offset) {
0531   Transform3 base =
0532       AngleAxis3(angle * 1_degree, Vector3::UnitX()) * Translation3(offset);
0533
0534   auto disc = makeDisc(base, 30_mm, 100_mm);
0535
0536   auto disc2 =
0537       makeDisc(base * AngleAxis3(14_degree, Vector3::UnitZ()), 100_mm, 150_mm);
0538
0539   auto [disc3, reversed] =
0540       disc->mergedWith(*disc2, Acts::AxisDirection::AxisR, false, *logger);
0541   BOOST_REQUIRE_NE(disc3, nullptr);
0542   BOOST_CHECK(!reversed);
0543
0544   auto [disc3Reversed, reversed2] =
0545       disc2->mergedWith(*disc, Acts::AxisDirection::AxisR, false, *logger);
0546   BOOST_REQUIRE_NE(disc3Reversed, nullptr);
0547   BOOST_CHECK(disc3->bounds() == disc3Reversed->bounds());
0548   BOOST_CHECK(reversed2);
0549
0550   const auto* bounds = dynamic_cast<const RadialBounds*>(&disc3->bounds());
0551   BOOST_REQUIRE_NE(bounds, nullptr);
0552
0553   BOOST_CHECK_EQUAL(bounds->get(RadialBounds::eMinR), 30_mm);
0554   BOOST_CHECK_EQUAL(bounds->get(RadialBounds::eMaxR), 150_mm);
0555
0556   // Disc did not move
0557   BOOST_CHECK_EQUAL(base.matrix(), disc3->transform(tgContext).matrix());
0558
0559   // Rotation in z depends on the ordering, the left side "wins"
0560   Transform3 expected12 = base;
0561   BOOST_CHECK_EQUAL(expected12.matrix(), disc3->transform(tgContext).matrix());
0562
0563   Transform3 expected21 = base * AngleAxis3(14_degree, Vector3::UnitZ());
0564   CHECK_CLOSE_OR_SMALL(disc3Reversed->transform(tgContext).matrix(),
0565                        expected21.matrix(), 1e-6, 1e-10);
0566
0567   // Test r merging with phi sectors (matching)
0568   auto discPhi1 = makeDisc(base, 30_mm, 100_mm, 10_degree, 40_degree);
0569   auto discPhi2 = makeDisc(base, 100_mm, 160_mm, 10_degree, 40_degree);
0570   auto [discPhi12, reversedPhi12] =
0571       discPhi1->mergedWith(*discPhi2, AxisDirection::AxisR, false, *logger);
0572   BOOST_REQUIRE_NE(discPhi12, nullptr);
0573
0574   const auto* boundsPhi12 =
0575       dynamic_cast<const RadialBounds*>(&discPhi12->bounds());
0576   BOOST_REQUIRE_NE(boundsPhi12, nullptr);
0577
0578   BOOST_CHECK_EQUAL(boundsPhi12->get(RadialBounds::eMinR), 30_mm);
0579   BOOST_CHECK_EQUAL(boundsPhi12->get(RadialBounds::eMaxR), 160_mm);
0580   BOOST_CHECK_EQUAL(boundsPhi12->get(RadialBounds::eHalfPhiSector), 10_degree);
0581 }
0582
0583 BOOST_DATA_TEST_CASE(IncompatiblePhiDirection,
0584                      (boost::unit_test::data::xrange(-135, 180, 45) *
0585                       boost::unit_test::data::make(Vector3{0_mm, 0_mm, 0_mm},
0586                                                    Vector3{20_mm, 0_mm, 0_mm},
0587                                                    Vector3{0_mm, 20_mm, 0_mm},
0588                                                    Vector3{20_mm, 20_mm, 0_mm},
0589                                                    Vector3{0_mm, 0_mm, 20_mm}) *
0590                       boost::unit_test::data::xrange(-1300, 1300, 104)),
0591                      angle, offset, phiShift) {
0592   Logging::ScopedFailureThreshold ft{Logging::FATAL};
0593   Transform3 base =
0594       AngleAxis3(angle * 1_degree, Vector3::UnitX()) * Translation3(offset);
0595
0596   auto a = [phiShift](double v) {
0597     return detail::radian_sym(v + phiShift * 1_degree);
0598   };
0599
0600   auto discPhi = makeDisc(base, 30_mm, 100_mm, 10_degree, a(40_degree));
0601
0602   // Disc with overlap in phi
0603   auto discPhi2 = makeDisc(base, 30_mm, 100_mm, 45_degree, a(85_degree));
0604   BOOST_CHECK_THROW(discPhi->mergedWith(*discPhi2, Acts::AxisDirection::AxisPhi,
0605                                         false, *logger),
0606                     SurfaceMergingException);
0607
0608   // Disc with gap in phi
0609   auto discPhi3 = makeDisc(base, 30_mm, 100_mm, 45_degree, a(105_degree));
0610   BOOST_CHECK_THROW(discPhi->mergedWith(*discPhi3, Acts::AxisDirection::AxisPhi,
0611                                         false, *logger),
0612                     SurfaceMergingException);
0613
0614   // Disc with a z shift
0615   auto discPhi4 = makeDisc(base * Translation3{Vector3::UnitZ() * 20_mm}, 30_mm,
0616                            100_mm, 45_degree, a(95_degree));
0617   BOOST_CHECK_THROW(discPhi->mergedWith(*discPhi4, Acts::AxisDirection::AxisPhi,
0618                                         false, *logger),
0619                     SurfaceMergingException);
0620
0621   // Disc with different r bounds: could be merged in r but not in phi
0622   auto discPhi5 = makeDisc(base, 100_mm, 150_mm, 45_degree, a(95_degree));
0623   BOOST_CHECK_THROW(discPhi->mergedWith(*discPhi5, Acts::AxisDirection::AxisPhi,
0624                                         false, *logger),
0625                     SurfaceMergingException);
0626 }
0627
0628 BOOST_DATA_TEST_CASE(PhiDirection,
0629                      (boost::unit_test::data::xrange(-135, 180, 45) *
0630                       boost::unit_test::data::make(Vector3{0_mm, 0_mm, 0_mm},
0631                                                    Vector3{20_mm, 0_mm, 0_mm},
0632                                                    Vector3{0_mm, 20_mm, 0_mm},
0633                                                    Vector3{20_mm, 20_mm, 0_mm},
0634                                                    Vector3{0_mm, 0_mm, 20_mm}) *
0635                       boost::unit_test::data::xrange(-1300, 1300, 104)),
0636                      angle, offset, phiShift) {
0637   Transform3 base =
0638       AngleAxis3(angle * 1_degree, Vector3::UnitX()) * Translation3(offset);
0639
0640   auto a = [phiShift](double v) {
0641     return detail::radian_sym(v + phiShift * 1_degree);
0642   };
0643
0644   BOOST_TEST_CONTEXT("Internal rotation") {
0645     auto disc = makeDisc(base, 30_mm, 100_mm, 10_degree, a(40_degree));
0646     auto disc2 = makeDisc(base, 30_mm, 100_mm, 45_degree, a(95_degree));
0647
0648     auto [disc3, reversed] =
0649         disc->mergedWith(*disc2, Acts::AxisDirection::AxisPhi, false, *logger);
0650     BOOST_REQUIRE_NE(disc3, nullptr);
0651     BOOST_CHECK_EQUAL(base.matrix(), disc3->transform(tgContext).matrix());
0652     BOOST_CHECK(reversed);
0653
0654     auto [disc3Reversed, reversed2] =
0655         disc2->mergedWith(*disc, Acts::AxisDirection::AxisPhi, false, *logger);
0656     BOOST_REQUIRE_NE(disc3Reversed, nullptr);
0657     BOOST_CHECK(*disc3 == *disc3Reversed);
0658     BOOST_CHECK(!reversed2);
0659
0660     const auto* bounds = dynamic_cast<const RadialBounds*>(&disc3->bounds());
0661     BOOST_REQUIRE_NE(bounds, nullptr);
0662
0663     BOOST_CHECK_SMALL(
0664         detail::difference_periodic(bounds->get(RadialBounds::eAveragePhi),
0665                                     a(85_degree), 2 * std::numbers::pi),
0666         1e-6);
0667     BOOST_CHECK_CLOSE(bounds->get(RadialBounds::eHalfPhiSector), 55_degree,
0668                       1e-6);
0669
0670     auto disc4 = makeDisc(base, 30_mm, 100_mm, 20_degree, a(170_degree));
0671     auto disc5 = makeDisc(base, 30_mm, 100_mm, 10_degree, a(-160_degree));
0672     auto [disc45, reversed45] =
0673         disc4->mergedWith(*disc5, Acts::AxisDirection::AxisPhi, false, *logger);
0674     BOOST_REQUIRE_NE(disc45, nullptr);
0675     BOOST_CHECK_EQUAL(base.matrix(), disc45->transform(tgContext).matrix());
0676     BOOST_CHECK(reversed45);
0677
0678     auto [disc54, reversed54] =
0679         disc5->mergedWith(*disc4, Acts::AxisDirection::AxisPhi, false, *logger);
0680     BOOST_REQUIRE_NE(disc54, nullptr);
0681     BOOST_CHECK(!reversed54);
0682
0683     BOOST_CHECK(*disc54 == *disc45);
0684
0685     const auto* bounds45 = dynamic_cast<const RadialBounds*>(&disc45->bounds());
0686     BOOST_REQUIRE_NE(bounds, nullptr);
0687
0688     BOOST_CHECK_SMALL(
0689         detail::difference_periodic(bounds45->get(RadialBounds::eAveragePhi),
0690                                     a(180_degree), 2 * std::numbers::pi),
0691         1e-6);
0692     BOOST_CHECK_CLOSE(bounds45->get(RadialBounds::eHalfPhiSector), 30_degree,
0693                       1e-6);
0694
0695     auto disc6 = makeDisc(base, 30_mm, 100_mm, 90_degree, a(0_degree));
0696     auto disc7 = makeDisc(base, 30_mm, 100_mm, 90_degree, a(180_degree));
0697
0698     auto [disc67, reversed67] =
0699         disc6->mergedWith(*disc7, Acts::AxisDirection::AxisPhi, false, *logger);
0700     BOOST_REQUIRE_NE(disc67, nullptr);
0701     CHECK_CLOSE_OR_SMALL(disc67->transform(tgContext).matrix(), base.matrix(),
0702                          1e-6, 1e-10);
0703     BOOST_CHECK(!reversed67);
0704
0705     auto [disc76, reversed76] =
0706         disc7->mergedWith(*disc6, Acts::AxisDirection::AxisPhi, false, *logger);
0707     BOOST_REQUIRE_NE(disc76, nullptr);
0708     // surfaces are not equal because bounds are not equal
0709     BOOST_CHECK(*disc76 != *disc67);
0710     // bounds are different because of avg phi
0711     BOOST_CHECK_NE(disc76->bounds(), disc67->bounds());
0712     // transforms should be the same
0713     BOOST_CHECK_EQUAL(disc76->transform(tgContext).matrix(),
0714                       disc67->transform(tgContext).matrix());
0715     // not reversed either because you get the ordering you put in
0716     BOOST_CHECK(!reversed76);
0717
0718     const auto* bounds67 = dynamic_cast<const RadialBounds*>(&disc67->bounds());
0719     BOOST_REQUIRE_NE(bounds67, nullptr);
0720     BOOST_CHECK_SMALL(
0721         detail::difference_periodic(bounds67->get(RadialBounds::eAveragePhi),
0722                                     a(90_degree), 2 * std::numbers::pi),
0723         1e-6);
0724     BOOST_CHECK_CLOSE(bounds67->get(RadialBounds::eHalfPhiSector), 180_degree,
0725                       1e-6);
0726   }
0727
0728   BOOST_TEST_CONTEXT("External rotation") {
0729     Transform3 trf1 = base * AngleAxis3(a(40_degree), Vector3::UnitZ());
0730     auto disc = makeDisc(trf1, 30_mm, 100_mm, 10_degree, 0_degree);
0731     Transform3 trf2 = base * AngleAxis3(a(95_degree), Vector3::UnitZ());
0732     auto disc2 = makeDisc(trf2, 30_mm, 100_mm, 45_degree, 0_degree);
0733
0734     auto [disc3, reversed] =
0735         disc->mergedWith(*disc2, Acts::AxisDirection::AxisPhi, true, *logger);
0736     BOOST_REQUIRE_NE(disc3, nullptr);
0737     Transform3 trfExpected12 =
0738         base * AngleAxis3(a(85_degree), Vector3::UnitZ());
0739     CHECK_CLOSE_OR_SMALL(disc3->transform(tgContext).matrix(),
0740                          trfExpected12.matrix(), 1e-6, 1e-10);
0741     BOOST_CHECK(reversed);
0742
0743     auto [disc3Reversed, reversed2] =
0744         disc2->mergedWith(*disc, Acts::AxisDirection::AxisPhi, true, *logger);
0745     BOOST_REQUIRE_NE(disc3Reversed, nullptr);
0746     BOOST_CHECK(*disc3 == *disc3Reversed);
0747     BOOST_CHECK(!reversed2);
0748
0749     const auto* bounds = dynamic_cast<const RadialBounds*>(&disc3->bounds());
0750     BOOST_REQUIRE_NE(bounds, nullptr);
0751
0752     BOOST_CHECK_EQUAL(bounds->get(RadialBounds::eAveragePhi), 0);
0753     BOOST_CHECK_CLOSE(bounds->get(RadialBounds::eHalfPhiSector), 55_degree,
0754                       1e-6);
0755
0756     Transform3 trf4 = base * AngleAxis3(a(170_degree), Vector3::UnitZ());
0757     auto disc4 = makeDisc(trf4, 30_mm, 100_mm, 20_degree, 0_degree);
0758     Transform3 trf5 = base * AngleAxis3(a(-160_degree), Vector3::UnitZ());
0759     auto disc5 = makeDisc(trf5, 30_mm, 100_mm, 10_degree, 0_degree);
0760     auto [disc45, reversed45] =
0761         disc4->mergedWith(*disc5, Acts::AxisDirection::AxisPhi, true, *logger);
0762     BOOST_REQUIRE_NE(disc45, nullptr);
0763     Transform3 trfExpected45 =
0764         base * AngleAxis3(a(180_degree), Vector3::UnitZ());
0765     CHECK_CLOSE_OR_SMALL(disc45->transform(tgContext).matrix(),
0766                          trfExpected45.matrix(), 1e-6, 1e-10);
0767     BOOST_CHECK(reversed45);
0768
0769     auto [disc54, reversed54] =
0770         disc5->mergedWith(*disc4, Acts::AxisDirection::AxisPhi, true, *logger);
0771     BOOST_REQUIRE_NE(disc54, nullptr);
0772     BOOST_CHECK(!reversed54);
0773
0774     BOOST_CHECK(*disc54 == *disc45);
0775
0776     const auto* bounds45 = dynamic_cast<const RadialBounds*>(&disc45->bounds());
0777     BOOST_REQUIRE_NE(bounds, nullptr);
0778
0779     BOOST_CHECK_EQUAL(bounds45->get(RadialBounds::eAveragePhi), 0);
0780     BOOST_CHECK_CLOSE(bounds45->get(RadialBounds::eHalfPhiSector), 30_degree,
0781                       1e-6);
0782
0783     Transform3 trf6 = base * AngleAxis3(a(0_degree), Vector3::UnitZ());
0784     auto disc6 = makeDisc(trf6, 30_mm, 100_mm, 90_degree, 0_degree);
0785     Transform3 trf7 = base * AngleAxis3(a(180_degree), Vector3::UnitZ());
0786     auto disc7 = makeDisc(trf7, 30_mm, 100_mm, 90_degree, 0_degree);
0787     auto [disc67, reversed67] =
0788         disc6->mergedWith(*disc7, Acts::AxisDirection::AxisPhi, true, *logger);
0789     BOOST_REQUIRE_NE(disc67, nullptr);
0790     Transform3 trfExpected67 =
0791         base * AngleAxis3(a(90_degree), Vector3::UnitZ());
0792     CHECK_CLOSE_OR_SMALL(disc67->transform(tgContext).matrix(),
0793                          trfExpected67.matrix(), 1e-6, 1e-10);
0794     BOOST_CHECK(!reversed67);
0795
0796     auto [disc76, reversed76] =
0797         disc7->mergedWith(*disc6, Acts::AxisDirection::AxisPhi, true, *logger);
0798     BOOST_REQUIRE_NE(disc76, nullptr);
0799     // surfaces are not equal due to different transforms
0800     BOOST_CHECK(*disc76 != *disc67);
0801     BOOST_CHECK_NE(disc76->transform(tgContext).matrix(),
0802                    disc67->transform(tgContext).matrix());
0803     // bounds should be equal however
0804     BOOST_CHECK_EQUAL(disc76->bounds(), disc67->bounds());
0805
0806     // not reversed either because you get the ordering you put in
0807     BOOST_CHECK(!reversed76);
0808
0809     const auto* bounds67 = dynamic_cast<const RadialBounds*>(&disc67->bounds());
0810     BOOST_REQUIRE_NE(bounds67, nullptr);
0811     BOOST_CHECK_EQUAL(bounds67->get(RadialBounds::eAveragePhi), 0);
0812     BOOST_CHECK_CLOSE(bounds67->get(RadialBounds::eHalfPhiSector), 180_degree,
0813                       1e-6);
0814   }
0815 }
0816
0817 BOOST_AUTO_TEST_SUITE_END()
0818
0819 BOOST_AUTO_TEST_SUITE_END()
0820
0821 }  // namespace Acts::Test