Core/Surfaces/PlaneSurfaceTests.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/tools/output_test_stream.hpp>
0010 #include <boost/test/unit_test.hpp>
0011
0012 #include "Acts/Definitions/Algebra.hpp"
0013 #include "Acts/Definitions/Alignment.hpp"
0014 #include "Acts/Definitions/Tolerance.hpp"
0015 #include "Acts/Definitions/TrackParametrization.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/PlaneSurface.hpp"
0021 #include "Acts/Surfaces/RectangleBounds.hpp"
0022 #include "Acts/Surfaces/Surface.hpp"
0023 #include "Acts/Surfaces/SurfaceBounds.hpp"
0024 #include "Acts/Surfaces/TrapezoidBounds.hpp"
0025 #include "Acts/Tests/CommonHelpers/DetectorElementStub.hpp"
0026 #include "Acts/Tests/CommonHelpers/FloatComparisons.hpp"
0027 #include "Acts/Utilities/BinningType.hpp"
0028 #include "Acts/Utilities/Intersection.hpp"
0029 #include "Acts/Utilities/Result.hpp"
0030 #include "Acts/Utilities/ThrowAssert.hpp"
0031
0032 #include <cmath>
0033 #include <memory>
0034 #include <numbers>
0035 #include <string>
0036
0037 using namespace Acts::UnitLiterals;
0038
0039 namespace Acts::Test {
0040
0041 // Create a test context
0042 GeometryContext tgContext = GeometryContext();
0043
0044 BOOST_AUTO_TEST_SUITE(PlaneSurfaces)
0045 /// Unit test for creating compliant/non-compliant PlaneSurface object
0046 BOOST_AUTO_TEST_CASE(PlaneSurfaceConstruction) {
0047   /// Test default construction
0048   // default construction is deleted
0049
0050   // bounds object, rectangle type
0051   auto rBounds = std::make_shared<const RectangleBounds>(3., 4.);
0052   /// Constructor with transform and bounds
0053   Translation3 translation{0., 1., 2.};
0054   auto pTransform = Transform3(translation);
0055
0056   /// Constructor with transform
0057   BOOST_CHECK_EQUAL(
0058       Surface::makeShared<PlaneSurface>(pTransform, rBounds)->type(),
0059       Surface::Plane);
0060
0061   /// Copy constructor
0062   auto planeSurfaceObject =
0063       Surface::makeShared<PlaneSurface>(pTransform, rBounds);
0064   auto copiedPlaneSurface =
0065       Surface::makeShared<PlaneSurface>(*planeSurfaceObject);
0066   BOOST_CHECK_EQUAL(copiedPlaneSurface->type(), Surface::Plane);
0067   BOOST_CHECK(*copiedPlaneSurface == *planeSurfaceObject);
0068
0069   /// Copied and transformed
0070   auto copiedTransformedPlaneSurface = Surface::makeShared<PlaneSurface>(
0071       tgContext, *planeSurfaceObject, pTransform);
0072   BOOST_CHECK_EQUAL(copiedTransformedPlaneSurface->type(), Surface::Plane);
0073
0074   /// Construct with nullptr bounds
0075   DetectorElementStub detElem;
0076   BOOST_CHECK_THROW(
0077       auto nullBounds = Surface::makeShared<PlaneSurface>(nullptr, detElem),
0078       AssertionFailureException);
0079 }
0080
0081 /// Unit test for testing PlaneSurface properties
0082 BOOST_AUTO_TEST_CASE(PlaneSurfaceProperties) {
0083   // bounds object, rectangle type
0084   auto rBounds = std::make_shared<const RectangleBounds>(3., 4.);
0085
0086   /// Test clone method
0087   Translation3 translation{0., 1., 2.};
0088   auto pTransform = Transform3(translation);
0089   auto planeSurfaceObject =
0090       Surface::makeShared<PlaneSurface>(pTransform, rBounds);
0091   // Is it in the right place?
0092   Translation3 translation2{0., 2., 4.};
0093   auto pTransform2 = Transform3(translation2);
0094   auto planeSurfaceObject2 =
0095       Surface::makeShared<PlaneSurface>(pTransform2, rBounds);
0096
0097   /// Test type (redundant)
0098   BOOST_CHECK_EQUAL(planeSurfaceObject->type(), Surface::Plane);
0099
0100   /// Test referencePosition
0101   Vector3 referencePosition{0., 1., 2.};
0102   BOOST_CHECK_EQUAL(
0103       planeSurfaceObject->referencePosition(tgContext, AxisDirection::AxisX),
0104       referencePosition);
0105
0106   /// Test referenceFrame
0107   Vector3 arbitraryGlobalPosition{2., 2., 2.};
0108   Vector3 momentum{1.e6, 1.e6, 1.e6};
0109   RotationMatrix3 expectedFrame;
0110   expectedFrame << 1., 0., 0., 0., 1., 0., 0., 0., 1.;
0111
0112   CHECK_CLOSE_OR_SMALL(planeSurfaceObject->referenceFrame(
0113                            tgContext, arbitraryGlobalPosition, momentum),
0114                        expectedFrame, 1e-6, 1e-9);
0115
0116   /// Test normal, given 3D position
0117   Vector3 normal3D(0., 0., 1.);
0118   BOOST_CHECK_EQUAL(planeSurfaceObject->normal(tgContext), normal3D);
0119
0120   /// Test bounds
0121   BOOST_CHECK_EQUAL(planeSurfaceObject->bounds().type(),
0122                     SurfaceBounds::eRectangle);
0123
0124   /// Test localToGlobal
0125   Vector2 localPosition{1.5, 1.7};
0126   Vector3 globalPosition =
0127       planeSurfaceObject->localToGlobal(tgContext, localPosition, momentum);
0128   // expected position is the translated one
0129   Vector3 expectedPosition{1.5 + translation.x(), 1.7 + translation.y(),
0130                            translation.z()};
0131
0132   CHECK_CLOSE_REL(globalPosition, expectedPosition, 1e-2);
0133
0134   /// Testing globalToLocal
0135   localPosition =
0136       planeSurfaceObject->globalToLocal(tgContext, globalPosition, momentum)
0137           .value();
0138   Vector2 expectedLocalPosition{1.5, 1.7};
0139
0140   CHECK_CLOSE_REL(localPosition, expectedLocalPosition, 1e-2);
0141
0142   Vector3 globalPositionOff =
0143       globalPosition +
0144       planeSurfaceObject->normal(tgContext, localPosition) * 0.1;
0145
0146   BOOST_CHECK(
0147       planeSurfaceObject->globalToLocal(tgContext, globalPositionOff, momentum)
0148           .error());
0149   BOOST_CHECK(planeSurfaceObject
0150                   ->globalToLocal(tgContext, globalPositionOff, momentum, 0.05)
0151                   .error());
0152   BOOST_CHECK(planeSurfaceObject
0153                   ->globalToLocal(tgContext, globalPositionOff, momentum, 0.2)
0154                   .ok());
0155
0156   /// Test isOnSurface
0157   Vector3 offSurface{0, 1, -2.};
0158   BOOST_CHECK(planeSurfaceObject->isOnSurface(
0159       tgContext, globalPosition, momentum, BoundaryTolerance::None()));
0160   BOOST_CHECK(planeSurfaceObject->isOnSurface(tgContext, globalPosition,
0161                                               BoundaryTolerance::None()));
0162   BOOST_CHECK(!planeSurfaceObject->isOnSurface(tgContext, offSurface, momentum,
0163                                                BoundaryTolerance::None()));
0164   BOOST_CHECK(!planeSurfaceObject->isOnSurface(tgContext, offSurface,
0165                                                BoundaryTolerance::None()));
0166
0167   /// Test intersection
0168   Vector3 direction{0., 0., 1.};
0169   auto sfIntersection = planeSurfaceObject
0170                             ->intersect(tgContext, offSurface, direction,
0171                                         BoundaryTolerance::None())
0172                             .closest();
0173   Intersection3D expectedIntersect{Vector3{0, 1, 2}, 4.,
0174                                    IntersectionStatus::reachable};
0175   BOOST_CHECK(sfIntersection.isValid());
0176   BOOST_CHECK_EQUAL(sfIntersection.position(), expectedIntersect.position());
0177   BOOST_CHECK_EQUAL(sfIntersection.pathLength(),
0178                     expectedIntersect.pathLength());
0179   BOOST_CHECK_EQUAL(sfIntersection.object(), planeSurfaceObject.get());
0180
0181   /// Test pathCorrection
0182   CHECK_CLOSE_REL(planeSurfaceObject->pathCorrection(tgContext, offSurface,
0183                                                      momentum.normalized()),
0184                   std::numbers::sqrt3, 0.01);
0185
0186   /// Test name
0187   BOOST_CHECK_EQUAL(planeSurfaceObject->name(),
0188                     std::string("Acts::PlaneSurface"));
0189
0190   /// Test dump
0191   boost::test_tools::output_test_stream dumpOutput;
0192   dumpOutput << planeSurfaceObject->toStream(tgContext);
0193   BOOST_CHECK(dumpOutput.is_equal(
0194       "Acts::PlaneSurface\n"
0195       "     Center position  (x, y, z) = (0.0000, 1.0000, 2.0000)\n"
0196       "     Rotation:             colX = (1.000000, 0.000000, 0.000000)\n"
0197       "                           colY = (0.000000, 1.000000, 0.000000)\n"
0198       "                           colZ = (0.000000, 0.000000, 1.000000)\n"
0199       "     Bounds  : Acts::RectangleBounds:  (hlX, hlY) = (3.0000000, "
0200       "4.0000000)\n"
0201       "(lower left, upper right):\n"
0202       "-3.0000000 -4.0000000\n"
0203       "3.0000000 4.0000000"));
0204 }
0205
0206 BOOST_AUTO_TEST_CASE(PlaneSurfaceEqualityOperators) {
0207   // rectangle bounds
0208   auto rBounds = std::make_shared<const RectangleBounds>(3., 4.);
0209   Translation3 translation{0., 1., 2.};
0210   auto pTransform = Transform3(translation);
0211   auto planeSurfaceObject =
0212       Surface::makeShared<PlaneSurface>(pTransform, rBounds);
0213   auto planeSurfaceObject2 =
0214       Surface::makeShared<PlaneSurface>(pTransform, rBounds);
0215
0216   /// Test equality operator
0217   BOOST_CHECK(*planeSurfaceObject == *planeSurfaceObject2);
0218
0219   BOOST_TEST_CHECKPOINT(
0220       "Create and then assign a PlaneSurface object to the existing one");
0221
0222   /// Test assignment
0223   auto assignedPlaneSurface =
0224       Surface::makeShared<PlaneSurface>(Transform3::Identity(), nullptr);
0225   *assignedPlaneSurface = *planeSurfaceObject;
0226
0227   /// Test equality of assigned to original
0228   BOOST_CHECK(*assignedPlaneSurface == *planeSurfaceObject);
0229 }
0230
0231 /// Unit test for testing PlaneSurface extent via Polyhedron representation
0232 BOOST_AUTO_TEST_CASE(PlaneSurfaceExtent) {
0233   // First test - non-rotated
0234   static const Transform3 planeZX =
0235       AngleAxis3(-std::numbers::pi / 2., Vector3::UnitX()) *
0236       AngleAxis3(-std::numbers::pi / 2., Vector3::UnitZ()) *
0237       Transform3::Identity();
0238
0239   double rHx = 2.;
0240   double rHy = 4.;
0241   double yPs = 3.;
0242   auto rBounds = std::make_shared<RectangleBounds>(rHx, rHy);
0243
0244   auto plane = Surface::makeShared<PlaneSurface>(
0245       Transform3(Translation3(Vector3(0., yPs, 0.)) * planeZX), rBounds);
0246
0247   auto planeExtent = plane->polyhedronRepresentation(tgContext, 1).extent();
0248
0249   CHECK_CLOSE_ABS(planeExtent.min(AxisDirection::AxisZ), -rHx,
0250                   s_onSurfaceTolerance);
0251   CHECK_CLOSE_ABS(planeExtent.max(AxisDirection::AxisZ), rHx,
0252                   s_onSurfaceTolerance);
0253   CHECK_CLOSE_ABS(planeExtent.min(AxisDirection::AxisX), -rHy,
0254                   s_onSurfaceTolerance);
0255   CHECK_CLOSE_ABS(planeExtent.max(AxisDirection::AxisX), rHy,
0256                   s_onSurfaceTolerance);
0257   CHECK_CLOSE_ABS(planeExtent.min(AxisDirection::AxisY), yPs,
0258                   s_onSurfaceTolerance);
0259   CHECK_CLOSE_ABS(planeExtent.max(AxisDirection::AxisY), yPs,
0260                   s_onSurfaceTolerance);
0261   CHECK_CLOSE_ABS(planeExtent.min(AxisDirection::AxisR), yPs,
0262                   s_onSurfaceTolerance);
0263   CHECK_CLOSE_ABS(planeExtent.max(AxisDirection::AxisR), std::hypot(yPs, rHy),
0264                   s_onSurfaceTolerance);
0265
0266   // Now rotate
0267   double alpha = 0.123;
0268   auto planeRot = Surface::makeShared<PlaneSurface>(
0269       Transform3(Translation3(Vector3(0., yPs, 0.)) *
0270                  AngleAxis3(alpha, Vector3(0., 0., 1.)) * planeZX),
0271       rBounds);
0272
0273   auto planeExtentRot =
0274       planeRot->polyhedronRepresentation(tgContext, 1).extent();
0275   CHECK_CLOSE_ABS(planeExtentRot.min(AxisDirection::AxisZ), -rHx,
0276                   s_onSurfaceTolerance);
0277   CHECK_CLOSE_ABS(planeExtentRot.max(AxisDirection::AxisZ), rHx,
0278                   s_onSurfaceTolerance);
0279   CHECK_CLOSE_ABS(planeExtentRot.min(AxisDirection::AxisX),
0280                   -rHy * std::cos(alpha), s_onSurfaceTolerance);
0281   CHECK_CLOSE_ABS(planeExtentRot.max(AxisDirection::AxisX),
0282                   rHy * std::cos(alpha), s_onSurfaceTolerance);
0283   CHECK_CLOSE_ABS(planeExtentRot.min(AxisDirection::AxisY),
0284                   yPs - rHy * std::sin(alpha), s_onSurfaceTolerance);
0285   CHECK_CLOSE_ABS(planeExtentRot.max(AxisDirection::AxisY),
0286                   yPs + rHy * std::sin(alpha), s_onSurfaceTolerance);
0287   CHECK_CLOSE_ABS(planeExtentRot.min(AxisDirection::AxisR),
0288                   yPs * std::cos(alpha), s_onSurfaceTolerance);
0289 }
0290
0291 BOOST_AUTO_TEST_CASE(RotatedTrapezoid) {
0292   const double shortHalfX = 100.;
0293   const double longHalfX = 200.;
0294   const double halfY = 300.;
0295   const double rotAngle = 45._degree;
0296
0297   Vector2 edgePoint{longHalfX - 10., halfY};
0298
0299   std::shared_ptr<TrapezoidBounds> bounds =
0300       std::make_shared<TrapezoidBounds>(shortHalfX, longHalfX, halfY);
0301
0302   BOOST_CHECK(bounds->inside(edgePoint, BoundaryTolerance::None()));
0303   BOOST_CHECK(!bounds->inside(Eigen::Rotation2D(-rotAngle) * edgePoint,
0304                               BoundaryTolerance::None()));
0305
0306   std::shared_ptr<TrapezoidBounds> rotatedBounds =
0307       std::make_shared<TrapezoidBounds>(shortHalfX, longHalfX, halfY, rotAngle);
0308
0309   BOOST_CHECK(!rotatedBounds->inside(edgePoint, BoundaryTolerance::None()));
0310   BOOST_CHECK(rotatedBounds->inside(Eigen::Rotation2D(-rotAngle) * edgePoint,
0311                                     BoundaryTolerance::None()));
0312 }
0313
0314 /// Unit test for testing PlaneSurface alignment derivatives
0315 BOOST_AUTO_TEST_CASE(PlaneSurfaceAlignment) {
0316   // bounds object, rectangle type
0317   auto rBounds = std::make_shared<const RectangleBounds>(3., 4.);
0318   // Test clone method
0319   Translation3 translation{0., 1., 2.};
0320   const double rotationAngle = std::numbers::pi / 2.;
0321   AngleAxis3 rotation(rotationAngle, Vector3::UnitY());
0322   RotationMatrix3 rotationMat = rotation.toRotationMatrix();
0323
0324   auto pTransform = Transform3{translation * rotationMat};
0325   auto planeSurfaceObject =
0326       Surface::makeShared<PlaneSurface>(pTransform, rBounds);
0327
0328   // The local frame z axis
0329   const Vector3 localZAxis = rotationMat.col(2);
0330   // Check the local z axis is aligned to global x axis
0331   CHECK_CLOSE_ABS(localZAxis, Vector3(1., 0., 0.), 1e-15);
0332
0333   // Define the track (local) position and direction
0334   Vector2 localPosition{1, 2};
0335   Vector3 momentum{1, 0, 0};
0336   Vector3 direction = momentum.normalized();
0337   // Get the global position
0338   Vector3 globalPosition =
0339       planeSurfaceObject->localToGlobal(tgContext, localPosition, momentum);
0340
0341   // (a) Test the derivative of path length w.r.t. alignment parameters
0342   const AlignmentToPathMatrix& alignToPath =
0343       planeSurfaceObject->alignmentToPathDerivative(tgContext, globalPosition,
0344                                                     direction);
0345   // The expected results
0346   AlignmentToPathMatrix expAlignToPath = AlignmentToPathMatrix::Zero();
0347   expAlignToPath << 1, 0, 0, 2, -1, 0;
0348
0349   // Check if the calculated derivative is as expected
0350   CHECK_CLOSE_ABS(alignToPath, expAlignToPath, 1e-10);
0351
0352   // (b) Test the derivative of bound track parameters local position w.r.t.
0353   // position in local 3D Cartesian coordinates
0354   const auto& loc3DToLocBound =
0355       planeSurfaceObject->localCartesianToBoundLocalDerivative(tgContext,
0356                                                                globalPosition);
0357   // For plane surface, this should be identity matrix
0358   CHECK_CLOSE_ABS(loc3DToLocBound, (ActsMatrix<2, 3>::Identity()), 1e-10);
0359
0360   // (c) Test the derivative of bound parameters (only test loc0, loc1 here)
0361   // w.r.t. alignment parameters
0362   FreeVector derivatives = FreeVector::Zero();
0363   derivatives.head<3>() = direction;
0364   const AlignmentToBoundMatrix& alignToBound =
0365       planeSurfaceObject->alignmentToBoundDerivative(tgContext, globalPosition,
0366                                                      direction, derivatives);
0367   const AlignmentToPathMatrix alignToloc0 =
0368       alignToBound.block<1, 6>(eBoundLoc0, eAlignmentCenter0);
0369   const AlignmentToPathMatrix alignToloc1 =
0370       alignToBound.block<1, 6>(eBoundLoc1, eAlignmentCenter0);
0371   // The expected results
0372   AlignmentToPathMatrix expAlignToloc0;
0373   expAlignToloc0 << 0, 0, 1, 0, 0, 2;
0374   AlignmentToPathMatrix expAlignToloc1;
0375   expAlignToloc1 << 0, -1, 0, 0, 0, -1;
0376   // Check if the calculated derivatives are as expected
0377   CHECK_CLOSE_ABS(alignToloc0, expAlignToloc0, 1e-10);
0378   CHECK_CLOSE_ABS(alignToloc1, expAlignToloc1, 1e-10);
0379 }
0380
0381 BOOST_AUTO_TEST_SUITE_END()
0382
0383 }  // namespace Acts::Test