Core/Geometry/CylinderVolumeBoundsTests.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/unit_test.hpp>
0011
0012 #include "Acts/Definitions/Algebra.hpp"
0013 #include "Acts/Definitions/Direction.hpp"
0014 #include "Acts/Geometry/CylinderVolumeBounds.hpp"
0015 #include "Acts/Geometry/GeometryContext.hpp"
0016 #include "Acts/Geometry/VolumeBounds.hpp"
0017 #include "Acts/Surfaces/CylinderBounds.hpp"
0018 #include "Acts/Surfaces/RadialBounds.hpp"
0019 #include "Acts/Tests/CommonHelpers/FloatComparisons.hpp"
0020 #include "Acts/Utilities/BinningType.hpp"
0021 #include "Acts/Utilities/BoundingBox.hpp"
0022
0023 #include <algorithm>
0024 #include <array>
0025 #include <cmath>
0026 #include <memory>
0027 #include <stdexcept>
0028 #include <vector>
0029
0030 namespace bdata = boost::unit_test::data;
0031
0032 namespace Acts {
0033 namespace Test {
0034 BOOST_AUTO_TEST_SUITE(Geometry)
0035
0036 BOOST_AUTO_TEST_CASE(CylinderVolumeBoundsConstruction) {
0037   double rmin{10.}, rmax{20.}, halfz{30.}, halfphi{M_PI / 4}, avgphi{0.};
0038
0039   // Test different construction modes: solid
0040   CylinderVolumeBounds solidCylinder(0., rmax, halfz);
0041   BOOST_CHECK_EQUAL(solidCylinder.orientedSurfaces().size(), 3);
0042
0043   // Test different construction modes: sectoral solid
0044   CylinderVolumeBounds solidCylinderSector(0., rmax, halfz, halfphi);
0045   BOOST_CHECK_EQUAL(solidCylinderSector.orientedSurfaces().size(), 5);
0046
0047   // Test different construction modes: tube
0048   CylinderVolumeBounds tubeCylinder(rmin, rmax, halfz);
0049   BOOST_CHECK_EQUAL(tubeCylinder.orientedSurfaces().size(), 4);
0050
0051   // Test different construction modes: sectoral tube
0052   CylinderVolumeBounds tubeCylinderSector(rmin, rmax, halfz, halfphi);
0053   BOOST_CHECK_EQUAL(tubeCylinderSector.orientedSurfaces().size(), 6);
0054
0055   CylinderVolumeBounds original(rmin, rmax, halfz, halfphi, avgphi);
0056
0057   // Test construction from CylinderBounds and thickness
0058   double rmed = 0.5 * (rmin + rmax);
0059   double rthickness = (rmax - rmin);
0060   CylinderBounds cBounds(rmed, halfz, halfphi, avgphi);
0061   CylinderVolumeBounds fromCylinder(cBounds, rthickness);
0062   BOOST_CHECK_EQUAL(original, fromCylinder);
0063
0064   // Test construction from RadialBounds and thickness
0065   RadialBounds rBounds(rmin, rmax, halfphi, avgphi);
0066   CylinderVolumeBounds fromDisc(rBounds, 2 * halfz);
0067   BOOST_CHECK_EQUAL(original, fromDisc);
0068
0069   // Test the copy construction
0070   CylinderVolumeBounds copied(original);
0071   BOOST_CHECK_EQUAL(original, copied);
0072
0073   // Test the assignment
0074   CylinderVolumeBounds assigned = original;
0075   BOOST_CHECK_EQUAL(original, assigned);
0076 }
0077
0078 BOOST_AUTO_TEST_CASE(CylinderVolumeBoundsRecreation) {
0079   double rmin{10.}, rmax{20.}, halfz{30.}, halfphi{M_PI / 4}, avgphi{0.};
0080
0081   CylinderVolumeBounds original(rmin, rmax, halfz, halfphi, avgphi);
0082   std::array<double, CylinderVolumeBounds::eSize> values{};
0083   std::vector<double> valvector = original.values();
0084   std::copy_n(valvector.begin(), CylinderVolumeBounds::eSize, values.begin());
0085   CylinderVolumeBounds recreated(values);
0086   BOOST_CHECK_EQUAL(original, recreated);
0087 }
0088
0089 BOOST_AUTO_TEST_CASE(CylinderVolumeBoundsExceptions) {
0090   double rmin{10.}, rmax{20.}, halfz{30.}, halfphi{M_PI / 4}, avgphi{0.};
0091
0092   // Negative inner radius
0093   BOOST_CHECK_THROW(CylinderVolumeBounds(-rmin, rmax, halfz, halfphi, avgphi),
0094                     std::logic_error);
0095
0096   // Negative outer radius
0097   BOOST_CHECK_THROW(CylinderVolumeBounds(rmin, -rmax, halfz, halfphi, avgphi),
0098                     std::logic_error);
0099
0100   // Swapped radii
0101   BOOST_CHECK_THROW(CylinderVolumeBounds(rmax, rmin, halfz, halfphi, avgphi),
0102                     std::logic_error);
0103
0104   // Zero half length
0105   BOOST_CHECK_THROW(CylinderVolumeBounds(rmax, rmin, 0., halfphi, avgphi),
0106                     std::logic_error);
0107
0108   // Negative half length
0109   BOOST_CHECK_THROW(CylinderVolumeBounds(rmax, rmin, -halfz, halfphi, avgphi),
0110                     std::logic_error);
0111
0112   // Out of bounds half phi
0113   BOOST_CHECK_THROW(CylinderVolumeBounds(rmax, rmin, halfz, -4., avgphi),
0114                     std::logic_error);
0115
0116   // Wrong positioning phi
0117   BOOST_CHECK_THROW(CylinderVolumeBounds(rmax, rmin, halfz, halfphi, 4.),
0118                     std::logic_error);
0119
0120   // Test construction from CylinderBounds and thickness
0121   double rmed = 0.5 * (rmin + rmax);
0122   CylinderBounds cBounds(rmed, halfz, halfphi, avgphi);
0123   RadialBounds rBounds(rmin, rmax, halfphi, avgphi);
0124
0125   // Negative thickness
0126   BOOST_CHECK_THROW(CylinderVolumeBounds(cBounds, -1.), std::logic_error);
0127
0128   // Wrong thickness
0129   BOOST_CHECK_THROW(CylinderVolumeBounds(cBounds, 1000.), std::logic_error);
0130
0131   // Test construction from RadialBounds and thickness
0132   BOOST_CHECK_THROW(CylinderVolumeBounds(rBounds, -1), std::logic_error);
0133 }
0134
0135 BOOST_AUTO_TEST_CASE(CylinderVolumeBoundsAccess) {
0136   double rmin{10.}, rmax{20.}, halfz{30.}, halfphi{M_PI / 4}, avgphi{0.};
0137   CylinderVolumeBounds cvBounds(rmin, rmax, halfz, halfphi, avgphi);
0138
0139   // Test the accessors
0140   BOOST_CHECK_EQUAL(cvBounds.get(CylinderVolumeBounds::eMinR), rmin);
0141   BOOST_CHECK_EQUAL(cvBounds.get(CylinderVolumeBounds::eMaxR), rmax);
0142   BOOST_CHECK_EQUAL(cvBounds.get(CylinderVolumeBounds::eHalfLengthZ), halfz);
0143   BOOST_CHECK_EQUAL(cvBounds.get(CylinderVolumeBounds::eHalfPhiSector),
0144                     halfphi);
0145   BOOST_CHECK_EQUAL(cvBounds.get(CylinderVolumeBounds::eAveragePhi), avgphi);
0146 }
0147
0148 /// Unit test for testing the orientedSurfaces() function
0149 BOOST_DATA_TEST_CASE(
0150     CylinderVolumeBoundsOrientedSurfaces,
0151     bdata::random((bdata::engine = std::mt19937(), bdata::seed = 1,
0152                    bdata::distribution =
0153                        std::uniform_real_distribution<double>(-M_PI, M_PI))) ^
0154         bdata::random((bdata::engine = std::mt19937(), bdata::seed = 2,
0155                        bdata::distribution =
0156                            std::uniform_real_distribution<double>(-M_PI,
0157                                                                   M_PI))) ^
0158         bdata::random((bdata::engine = std::mt19937(), bdata::seed = 3,
0159                        bdata::distribution =
0160                            std::uniform_real_distribution<double>(-M_PI,
0161                                                                   M_PI))) ^
0162         bdata::random((bdata::engine = std::mt19937(), bdata::seed = 4,
0163                        bdata::distribution =
0164                            std::uniform_real_distribution<double>(-10., 10.))) ^
0165         bdata::random((bdata::engine = std::mt19937(), bdata::seed = 5,
0166                        bdata::distribution =
0167                            std::uniform_real_distribution<double>(-10., 10.))) ^
0168         bdata::random((bdata::engine = std::mt19937(), bdata::seed = 6,
0169                        bdata::distribution =
0170                            std::uniform_real_distribution<double>(-10., 10.))) ^
0171         bdata::xrange(100),
0172     alpha, beta, gamma, posX, posY, posZ, index) {
0173   (void)index;
0174
0175   // Create a test context
0176   GeometryContext tgContext = GeometryContext();
0177
0178   // position of volume
0179   const Vector3 pos(posX, posY, posZ);
0180   // rotation around x axis
0181   AngleAxis3 rotX(alpha, Vector3(1., 0., 0.));
0182   // rotation around y axis
0183   AngleAxis3 rotY(beta, Vector3(0., 1., 0.));
0184   // rotation around z axis
0185   AngleAxis3 rotZ(gamma, Vector3(0., 0., 1.));
0186
0187   // create the cylinder bounds
0188   double rmin = 1.;
0189   double rmax = 2.;
0190   double halfz = 3.;
0191   CylinderVolumeBounds cylBounds(rmin, rmax, halfz);
0192   // create the transformation matrix
0193   auto transform = Transform3(Translation3(pos));
0194   transform *= rotZ;
0195   transform *= rotY;
0196   transform *= rotX;
0197   // get the boundary surfaces
0198   auto boundarySurfaces = cylBounds.orientedSurfaces(transform);
0199   // Test
0200
0201   // check if difference is halfZ - sign and direction independent
0202   CHECK_CLOSE_REL(
0203       (pos - boundarySurfaces.at(0).surface->center(tgContext)).norm(),
0204       cylBounds.get(CylinderVolumeBounds::eHalfLengthZ), 1e-12);
0205   CHECK_CLOSE_REL(
0206       (pos - boundarySurfaces.at(1).surface->center(tgContext)).norm(),
0207       cylBounds.get(CylinderVolumeBounds::eHalfLengthZ), 1e-12);
0208   // transform to local
0209   double posDiscPosZ =
0210       (transform.inverse() * boundarySurfaces.at(1).surface->center(tgContext))
0211           .z();
0212   double centerPosZ = (transform.inverse() * pos).z();
0213   double negDiscPosZ =
0214       (transform.inverse() * boundarySurfaces.at(0).surface->center(tgContext))
0215           .z();
0216   // check if center of disc boundaries lies in the middle in z
0217   BOOST_CHECK_LT(centerPosZ, posDiscPosZ);
0218   BOOST_CHECK_GT(centerPosZ, negDiscPosZ);
0219   // check positions of disc boundarysurfaces
0220   // checks for zero value. double precision value is not exact.
0221   CHECK_CLOSE_ABS(
0222       negDiscPosZ + cylBounds.get(CylinderVolumeBounds::eHalfLengthZ),
0223       centerPosZ, 1e-12);
0224   CHECK_CLOSE_ABS(
0225       posDiscPosZ - cylBounds.get(CylinderVolumeBounds::eHalfLengthZ),
0226       centerPosZ, 1e-12);
0227   // orientation of disc surfaces
0228   // positive disc durface should point in positive direction in the frame of
0229   // the volume
0230   CHECK_CLOSE_REL(
0231       transform.rotation().col(2).dot(boundarySurfaces.at(1).surface->normal(
0232           tgContext, Acts::Vector2(0., 0.))),
0233       1., 1e-12);
0234   // negative disc durface should point in positive direction in the frame of
0235   // the volume
0236   CHECK_CLOSE_REL(
0237       transform.rotation().col(2).dot(boundarySurfaces.at(0).surface->normal(
0238           tgContext, Acts::Vector2(0., 0.))),
0239       1., 1e-12);
0240   // test in r
0241   CHECK_CLOSE_REL(boundarySurfaces.at(3).surface->center(tgContext), pos,
0242                   1e-12);
0243   CHECK_CLOSE_REL(boundarySurfaces.at(2).surface->center(tgContext), pos,
0244                   1e-12);
0245 }
0246
0247 BOOST_AUTO_TEST_CASE(CylinderVolumeBoundsBoundingBox) {
0248   GeometryContext tgContext = GeometryContext();
0249
0250   float tol = 1e-4;
0251
0252   CylinderVolumeBounds cvb(0., 5, 10);
0253   auto bb = cvb.boundingBox();
0254
0255   Transform3 rot;
0256   rot = AngleAxis3(M_PI / 2., Vector3::UnitX());
0257
0258   BOOST_CHECK_EQUAL(bb.entity(), nullptr);
0259   BOOST_CHECK_EQUAL(bb.max(), Vector3(5, 5, 10));
0260   BOOST_CHECK_EQUAL(bb.min(), Vector3(-5, -5, -10));
0261
0262   bb = cvb.boundingBox(&rot);
0263   BOOST_CHECK_EQUAL(bb.entity(), nullptr);
0264   CHECK_CLOSE_ABS(bb.max(), Vector3(5, 10, 5), tol);
0265   CHECK_CLOSE_ABS(bb.min(), Vector3(-5, -10, -5), tol);
0266
0267   cvb = CylinderVolumeBounds(5, 8, 12);
0268   bb = cvb.boundingBox();
0269   BOOST_CHECK_EQUAL(bb.entity(), nullptr);
0270   BOOST_CHECK_EQUAL(bb.max(), Vector3(8, 8, 12));
0271   BOOST_CHECK_EQUAL(bb.min(), Vector3(-8, -8, -12));
0272
0273   double angle = M_PI / 8.;
0274   cvb = CylinderVolumeBounds(5, 8, 13, angle);
0275   bb = cvb.boundingBox();
0276   BOOST_CHECK_EQUAL(bb.entity(), nullptr);
0277   CHECK_CLOSE_ABS(bb.max(), Vector3(8, 8 * std::sin(angle), 13), tol);
0278   CHECK_CLOSE_ABS(bb.min(),
0279                   Vector3(5 * std::cos(angle), -8 * std::sin(angle), -13), tol);
0280
0281   rot = AngleAxis3(M_PI / 2., Vector3::UnitZ());
0282   bb = cvb.boundingBox(&rot);
0283   BOOST_CHECK_EQUAL(bb.entity(), nullptr);
0284   CHECK_CLOSE_ABS(bb.max(), Vector3(8 * std::sin(angle), 8, 13), tol);
0285   CHECK_CLOSE_ABS(bb.min(),
0286                   Vector3(-8 * std::sin(angle), 5 * std::cos(angle), -13), tol);
0287
0288   rot = AngleAxis3(M_PI / 2., Vector3(-2, 4, 5).normalized());
0289   bb = cvb.boundingBox(&rot);
0290   BOOST_CHECK_EQUAL(bb.entity(), nullptr);
0291   CHECK_CLOSE_ABS(bb.max(), Vector3(8.40007, 15.2828, 3.88911), tol);
0292   CHECK_CLOSE_ABS(bb.min(), Vector3(-7.27834, -8.12028, -14.2182), tol);
0293 }
0294
0295 BOOST_AUTO_TEST_CASE(CylinderVolumeOrientedBoundaries) {
0296   GeometryContext tgContext = GeometryContext();
0297
0298   CylinderVolumeBounds cvb(5, 10, 20);
0299
0300   auto cvbOrientedSurfaces = cvb.orientedSurfaces(Transform3::Identity());
0301   BOOST_CHECK_EQUAL(cvbOrientedSurfaces.size(), 4);
0302
0303   auto geoCtx = GeometryContext();
0304   Vector3 xaxis(1., 0., 0.);
0305   Vector3 yaxis(0., 1., 0.);
0306   Vector3 zaxis(0., 0., 1.);
0307
0308   for (auto& os : cvbOrientedSurfaces) {
0309     auto onSurface = os.surface->binningPosition(geoCtx, binR);
0310     auto locPos =
0311         os.surface->globalToLocal(geoCtx, onSurface, Vector3::Zero()).value();
0312     auto osNormal = os.surface->normal(geoCtx, locPos);
0313     // Check if you step inside the volume with the oriented normal
0314     Vector3 insideCvb = onSurface + os.direction * osNormal;
0315     Vector3 outsideCvb = onSurface - os.direction * osNormal;
0316
0317     BOOST_CHECK(cvb.inside(insideCvb));
0318     BOOST_CHECK(!cvb.inside(outsideCvb));
0319
0320     // Test the orientation of the boundary surfaces
0321     auto rot = os.surface->transform(geoCtx).rotation();
0322     BOOST_CHECK(rot.col(0).isApprox(xaxis));
0323     BOOST_CHECK(rot.col(1).isApprox(yaxis));
0324     BOOST_CHECK(rot.col(2).isApprox(zaxis));
0325   }
0326 }
0327
0328 BOOST_AUTO_TEST_CASE(CylinderVolumeBoundsSetValues) {
0329   CylinderVolumeBounds cyl(100, 300, 200);
0330
0331   BOOST_CHECK_THROW(cyl.set(CylinderVolumeBounds::eMinR, 400),
0332                     std::invalid_argument);
0333   BOOST_CHECK_EQUAL(cyl.get(CylinderVolumeBounds::eMinR), 100);
0334
0335   cyl.set(CylinderVolumeBounds::eMinR, 200);
0336   BOOST_CHECK_EQUAL(cyl.get(CylinderVolumeBounds::eMinR), 200);
0337
0338   BOOST_CHECK_THROW(cyl.set(CylinderVolumeBounds::eMaxR, 50),
0339                     std::invalid_argument);
0340   BOOST_CHECK_EQUAL(cyl.get(CylinderVolumeBounds::eMaxR), 300);
0341
0342   cyl.set(CylinderVolumeBounds::eMaxR, 250);
0343   BOOST_CHECK_EQUAL(cyl.get(CylinderVolumeBounds::eMaxR), 250);
0344
0345   BOOST_CHECK_THROW(cyl.set(CylinderVolumeBounds::eHalfLengthZ, -200),
0346                     std::invalid_argument);
0347   BOOST_CHECK_EQUAL(cyl.get(CylinderVolumeBounds::eHalfLengthZ), 200);
0348
0349   cyl.set(CylinderVolumeBounds::eHalfLengthZ, 250);
0350   BOOST_CHECK_EQUAL(cyl.get(CylinderVolumeBounds::eHalfLengthZ), 250);
0351
0352   cyl.set(CylinderVolumeBounds::eHalfLengthZ, 150);
0353   BOOST_CHECK_EQUAL(cyl.get(CylinderVolumeBounds::eHalfLengthZ), 150);
0354
0355   BOOST_CHECK_THROW(cyl.set(CylinderVolumeBounds::eHalfPhiSector, -M_PI),
0356                     std::invalid_argument);
0357   BOOST_CHECK_EQUAL(cyl.get(CylinderVolumeBounds::eHalfPhiSector), M_PI);
0358
0359   BOOST_CHECK_THROW(cyl.set(CylinderVolumeBounds::eHalfPhiSector, 1.5 * M_PI),
0360                     std::invalid_argument);
0361   BOOST_CHECK_EQUAL(cyl.get(CylinderVolumeBounds::eHalfPhiSector), M_PI);
0362
0363   cyl.set(CylinderVolumeBounds::eHalfPhiSector, M_PI / 2);
0364   BOOST_CHECK_EQUAL(cyl.get(CylinderVolumeBounds::eHalfPhiSector), M_PI / 2);
0365
0366   for (auto bValue :
0367        {CylinderVolumeBounds::eAveragePhi, CylinderVolumeBounds::eBevelMaxZ,
0368         CylinderVolumeBounds::eBevelMinZ}) {
0369     BOOST_CHECK_THROW(cyl.set(bValue, -1.5 * M_PI), std::invalid_argument);
0370     BOOST_CHECK_EQUAL(cyl.get(bValue), 0);
0371
0372     BOOST_CHECK_THROW(cyl.set(bValue, 1.5 * M_PI), std::invalid_argument);
0373     BOOST_CHECK_EQUAL(cyl.get(bValue), 0);
0374
0375     cyl.set(bValue, 0.5 * M_PI);
0376     BOOST_CHECK_EQUAL(cyl.get(bValue), 0.5 * M_PI);
0377     cyl.set(bValue, -0.5 * M_PI);
0378     BOOST_CHECK_EQUAL(cyl.get(bValue), -0.5 * M_PI);
0379   }
0380
0381   cyl = CylinderVolumeBounds(100, 300, 200);
0382   auto previous = cyl.values();
0383
0384   BOOST_CHECK_THROW(cyl.set({
0385                         {CylinderVolumeBounds::eMinR, 50},
0386                         {CylinderVolumeBounds::eMaxR, 200},
0387                         {CylinderVolumeBounds::eHalfLengthZ, -1},
0388                     }),
0389                     std::logic_error);
0390   auto act = cyl.values();
0391   BOOST_CHECK_EQUAL_COLLECTIONS(previous.begin(), previous.end(), act.begin(),
0392                                 act.end());
0393
0394   cyl.set({
0395       {CylinderVolumeBounds::eMinR, 50},
0396       {CylinderVolumeBounds::eMaxR, 200},
0397       {CylinderVolumeBounds::eHalfLengthZ, 150},
0398   });
0399
0400   BOOST_CHECK_EQUAL(cyl.get(CylinderVolumeBounds::eMinR), 50);
0401   BOOST_CHECK_EQUAL(cyl.get(CylinderVolumeBounds::eMaxR), 200);
0402   BOOST_CHECK_EQUAL(cyl.get(CylinderVolumeBounds::eHalfLengthZ), 150);
0403   BOOST_CHECK_EQUAL(cyl.get(CylinderVolumeBounds::eHalfPhiSector), M_PI);
0404   BOOST_CHECK_EQUAL(cyl.get(CylinderVolumeBounds::eAveragePhi), 0);
0405   BOOST_CHECK_EQUAL(cyl.get(CylinderVolumeBounds::eBevelMinZ), 0);
0406   BOOST_CHECK_EQUAL(cyl.get(CylinderVolumeBounds::eBevelMaxZ), 0);
0407 }
0408
0409 BOOST_AUTO_TEST_SUITE_END()
0410 }  // namespace Test
0411 }  // namespace Acts