Core/Surfaces/PlaneSurfaceTests.cpp

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