Core/Geometry/AlignmentContextTests.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/Definitions/Units.hpp"
0013 #include "Acts/Geometry/Blueprint.hpp"
0014 #include "Acts/Geometry/BlueprintNode.hpp"
0015 #include "Acts/Geometry/ContainerBlueprintNode.hpp"
0016 #include "Acts/Geometry/CuboidVolumeBounds.hpp"
0017 #include "Acts/Geometry/DiamondVolumeBounds.hpp"
0018 #include "Acts/Geometry/GeometryContext.hpp"
0019 #include "Acts/Geometry/StaticBlueprintNode.hpp"
0020 #include "Acts/Geometry/TrackingVolume.hpp"
0021 #include "Acts/Geometry/detail/TrackingGeometryPrintVisitor.hpp"
0022 #include "Acts/Surfaces/PlaneSurface.hpp"
0023 #include "Acts/Surfaces/RectangleBounds.hpp"
0024 #include "Acts/Surfaces/Surface.hpp"
0025 #include "Acts/Utilities/Result.hpp"
0026 #include "Acts/Utilities/StringHelpers.hpp"
0027
0028 #include <array>
0029 #include <memory>
0030 #include <utility>
0031
0032 namespace {
0033 /// Checks whether the two transforms are the same
0034 /// @param a: First transform to check
0035 /// @param b: Second transform to check
0036 inline bool isSame(const Acts::Transform3& a, const Acts::Transform3& b) {
0037   const Acts::Transform3 c = a * b.inverse();
0038   if (c.translation().norm() > Acts::s_onSurfaceTolerance) {
0039     return false;
0040   }
0041   for (std::size_t d = 0; d < 3; ++d) {
0042     const Acts::Vector3 e = Acts::Vector3::Unit(d);
0043     if (std::abs(e.dot(c * e) - 1.) > Acts::s_onSurfaceTolerance) {
0044       return false;
0045     }
0046   }
0047   return true;
0048 }
0049
0050 }  // namespace
0051
0052 const Acts::Logger& logger() {
0053   static const auto logObj =
0054       Acts::getDefaultLogger("UnitTests", Acts::Logging::VERBOSE);
0055   return *logObj;
0056 }
0057
0058 namespace Acts {
0059 class PlanarBounds;
0060 }  // namespace Acts
0061
0062 using namespace Acts;
0063 using namespace UnitLiterals;
0064
0065 namespace ActsTests {
0066
0067 /// @class AlignmentContext
0068 struct AlignmentContext {
0069   /// Return the geometry context
0070   Acts::GeometryContext getContext() const {
0071     return Acts::GeometryContext{this};
0072   }
0073   /// We have 2 different transforms
0074   using StoreType_t = std::array<Transform3, 2>;
0075   std::shared_ptr<const StoreType_t> detElementAlignment = nullptr;
0076   // Two transforms for the local -> global trf of the volume
0077   std::shared_ptr<StoreType_t> volumeLocToGlobAlign = nullptr;
0078   // Two transforms of the global -> local trf of the volume
0079   std::shared_ptr<StoreType_t> volumeGlobToLocAlign = nullptr;
0080   // List of portal transforms
0081   std::vector<StoreType_t> portalAlignments{};
0082
0083   /// Context index
0084   unsigned int alignmentIndex{0};
0085
0086   /// Default constructor
0087   AlignmentContext() = default;
0088
0089   /// Constructor with Store and context index
0090   explicit AlignmentContext(
0091       std::shared_ptr<const std::array<Transform3, 2>> aStore,
0092       unsigned int aIndex = 0)
0093       : detElementAlignment(std::move(aStore)), alignmentIndex(aIndex) {}
0094 };
0095
0096 /// @class AlignableDetectorElement
0097 ///
0098 /// This is a lightweight type of detector element,
0099 /// it simply implements the base class.
0100 class AlignableDetectorElement : public SurfacePlacementBase {
0101  public:
0102   // Deleted default constructor
0103   AlignableDetectorElement() = delete;
0104
0105   /// Constructor for single sided detector element
0106   /// - bound to a Plane Surface
0107   ///
0108   /// @param transform is the transform that element the layer in 3D frame
0109   /// @param pBounds is the planar bounds for the planar detector element
0110   /// @param thickness is the module thickness
0111   AlignableDetectorElement(std::shared_ptr<const Transform3> transform,
0112                            const std::shared_ptr<const PlanarBounds>& pBounds,
0113                            double thickness)
0114       : m_elementTransform(std::move(transform)) {
0115     m_elementSurface = Surface::makeShared<PlaneSurface>(pBounds, *this);
0116     m_elementSurface->assignThickness(thickness);
0117   }
0118
0119   ///  Destructor
0120   ~AlignableDetectorElement() override = default;
0121
0122   /// Return local to global transform associated with this identifier
0123   ///
0124   /// @param gctx The current geometry context object, e.g. alignment
0125   ///
0126   /// @note this is called from the surface().transform() in the PROXY mode
0127   const Transform3& localToGlobalTransform(
0128       const GeometryContext& gctx) const override;
0129
0130   /// Return surface associated with this detector element
0131   const Surface& surface() const override;
0132
0133   /// Non-const access to the surface associated with this detector element
0134   Surface& surface() override;
0135
0136   /// Is the detector element a sensitive element
0137   bool isSensitive() const override { return true; }
0138
0139  private:
0140   /// the transform for positioning in 3D space
0141   std::shared_ptr<const Transform3> m_elementTransform;
0142   /// the surface represented by it
0143   std::shared_ptr<Surface> m_elementSurface{nullptr};
0144 };
0145
0146 inline const Transform3& AlignableDetectorElement::localToGlobalTransform(
0147     const GeometryContext& gctx) const {
0148   const auto* alignContext = gctx.get<const AlignmentContext*>();
0149   if (alignContext != nullptr && alignContext->detElementAlignment != nullptr &&
0150       alignContext->alignmentIndex < 2) {
0151     return alignContext->detElementAlignment->at(alignContext->alignmentIndex);
0152   }
0153   return (*m_elementTransform);
0154 }
0155
0156 inline const Surface& AlignableDetectorElement::surface() const {
0157   return *m_elementSurface;
0158 }
0159
0160 inline Surface& AlignableDetectorElement::surface() {
0161   return *m_elementSurface;
0162 }
0163
0164 class AlignableVolumePlacement : public VolumePlacementBase {
0165  public:
0166   ///
0167   explicit AlignableVolumePlacement(const Transform3& volTrf)
0168       : m_locToGlob{volTrf} {}
0169
0170   const Transform3& localToGlobalTransform(
0171       const GeometryContext& gctx) const override {
0172     const auto* alignContext = gctx.get<const AlignmentContext*>();
0173     if (alignContext != nullptr &&
0174         alignContext->volumeLocToGlobAlign != nullptr &&
0175         alignContext->alignmentIndex < 2) {
0176       return alignContext->volumeLocToGlobAlign->at(
0177           alignContext->alignmentIndex);
0178     }
0179     return m_locToGlob;
0180   }
0181
0182   const Transform3& globalToLocalTransform(
0183       const GeometryContext& gctx) const override {
0184     const auto* alignContext = gctx.get<const AlignmentContext*>();
0185     if (alignContext != nullptr &&
0186         alignContext->volumeGlobToLocAlign != nullptr &&
0187         alignContext->alignmentIndex < 2) {
0188       return alignContext->volumeGlobToLocAlign->at(
0189           alignContext->alignmentIndex);
0190     }
0191     return m_globToLoc;
0192   }
0193
0194   const Transform3& portalLocalToGlobal(
0195       const GeometryContext& gctx, const std::size_t portalIdx) const override {
0196     const auto* alignContext = gctx.get<const AlignmentContext*>();
0197     if (alignContext != nullptr && alignContext->alignmentIndex < 2 &&
0198         alignContext->portalAlignments.size() > portalIdx) {
0199       return alignContext
0200           ->portalAlignments[portalIdx][alignContext->alignmentIndex];
0201     }
0202     assert(portalIdx < m_portalTrfs.size());
0203     return m_portalTrfs[portalIdx];
0204   }
0205
0206   void setAlignmentDelta(AlignmentContext& context, const Transform3& delta,
0207                          const std::size_t idx) {
0208     using StoreType_t = AlignmentContext::StoreType_t;
0209     assert(idx < 2);
0210
0211     if (context.volumeLocToGlobAlign == nullptr) {
0212       context.volumeLocToGlobAlign = std::make_shared<StoreType_t>();
0213       context.volumeGlobToLocAlign = std::make_shared<StoreType_t>();
0214     }
0215     context.volumeLocToGlobAlign->at(idx) = m_locToGlob * delta;
0216     context.volumeGlobToLocAlign->at(idx) =
0217         context.volumeLocToGlobAlign->at(idx).inverse();
0218
0219     context.portalAlignments.resize(nPortalPlacements());
0220     for (std::size_t portal = 0ul; portal < nPortalPlacements(); ++portal) {
0221       context.portalAlignments[portal][idx] =
0222           alignPortal(context.getContext(), portal);
0223     }
0224   }
0225
0226   void makePortalsAlignable(const GeometryContext& gctx,
0227                             const std::vector<std::shared_ptr<RegularSurface>>&
0228                                 portalsToAlign) override {
0229     VolumePlacementBase::makePortalsAlignable(gctx, portalsToAlign);
0230     for (std::size_t portal = 0ul; portal < portalsToAlign.size(); ++portal) {
0231       m_portalTrfs.push_back(alignPortal(gctx, portal));
0232     }
0233   }
0234
0235  private:
0236   Transform3 m_locToGlob{Transform3::Identity()};
0237   Transform3 m_globToLoc{m_locToGlob.inverse()};
0238   std::vector<Transform3> m_portalTrfs{};
0239 };
0240
0241 }  // namespace ActsTests
0242
0243 using namespace ActsTests;
0244
0245 BOOST_AUTO_TEST_SUITE(GeometrySuite);
0246
0247 /// Unit test for creating compliant/non-compliant Surface object
0248 BOOST_AUTO_TEST_CASE(AlignmentContextTests) {
0249   // The nominal and alignments
0250   Vector3 nominalCenter(0., 0., 0.);
0251   Vector3 negativeCenter(0., 0., -1.);
0252   Vector3 positiveCenter(0., 0., 1.);
0253
0254   // Checkpoints
0255   Vector3 onNominal(3., 3., 0.);
0256   Vector3 onNegative(3., 3., -1.);
0257   Vector3 onPositive(3., 3., 1.);
0258
0259   // Local position
0260   Vector2 localPosition(3., 3.);
0261
0262   // A position placeholder and dummy momentum
0263   Vector3 globalPosition(100_cm, 100_cm, 100_cm);
0264   Vector3 dummyMomentum(4., 4., 4.);
0265
0266   Transform3 negativeTransform = Transform3::Identity();
0267   negativeTransform.translation() = negativeCenter;
0268
0269   Transform3 positiveTransform = Transform3::Identity();
0270   positiveTransform.translation() = positiveCenter;
0271
0272   std::array<Transform3, 2> alignmentArray = {negativeTransform,
0273                                               positiveTransform};
0274
0275   auto detElementAlignment =
0276       std::make_shared<const std::array<Transform3, 2>>(alignmentArray);
0277
0278   // The detector element at nominal position
0279   AlignableDetectorElement alignedElement(
0280       std::make_shared<const Transform3>(Transform3::Identity()),
0281       std::make_shared<const RectangleBounds>(100_cm, 100_cm), 1_mm);
0282
0283   const auto& alignedSurface = alignedElement.surface();
0284
0285   // The alignment contexts
0286   AlignmentContext alignDefault{};
0287   AlignmentContext alignNegative{detElementAlignment, 0};
0288   AlignmentContext alignPositive{detElementAlignment, 1};
0289
0290   GeometryContext defaultContext{alignDefault.getContext()};
0291   GeometryContext negativeContext{alignNegative.getContext()};
0292   GeometryContext positiveContext{alignPositive.getContext()};
0293
0294   // Test the transforms
0295   BOOST_CHECK(isSame(alignedSurface.localToGlobalTransform(defaultContext),
0296                      Transform3::Identity()));
0297   BOOST_CHECK(isSame(alignedSurface.localToGlobalTransform(negativeContext),
0298                      negativeTransform));
0299   BOOST_CHECK(isSame(alignedSurface.localToGlobalTransform(positiveContext),
0300                      positiveTransform));
0301
0302   // Test the centers
0303   BOOST_CHECK_EQUAL(alignedSurface.center(defaultContext), nominalCenter);
0304   BOOST_CHECK_EQUAL(alignedSurface.center(negativeContext), negativeCenter);
0305   BOOST_CHECK_EQUAL(alignedSurface.center(positiveContext), positiveCenter);
0306
0307   // Test OnSurface
0308   BOOST_CHECK(
0309       alignedSurface.isOnSurface(defaultContext, onNominal, dummyMomentum));
0310   BOOST_CHECK(
0311       alignedSurface.isOnSurface(negativeContext, onNegative, dummyMomentum));
0312   BOOST_CHECK(
0313       alignedSurface.isOnSurface(positiveContext, onPositive, dummyMomentum));
0314
0315   // Test local to Global and vice versa
0316   globalPosition = alignedSurface.localToGlobal(defaultContext, localPosition,
0317                                                 dummyMomentum);
0318   BOOST_CHECK_EQUAL(globalPosition, onNominal);
0319   localPosition =
0320       alignedSurface.globalToLocal(defaultContext, onNominal, dummyMomentum)
0321           .value();
0322   BOOST_CHECK_EQUAL(localPosition, Vector2(3., 3.));
0323
0324   globalPosition = alignedSurface.localToGlobal(negativeContext, localPosition,
0325                                                 dummyMomentum);
0326   BOOST_CHECK_EQUAL(globalPosition, onNegative);
0327   localPosition =
0328       alignedSurface.globalToLocal(negativeContext, onNegative, dummyMomentum)
0329           .value();
0330   BOOST_CHECK_EQUAL(localPosition, Vector2(3., 3.));
0331
0332   globalPosition = alignedSurface.localToGlobal(positiveContext, localPosition,
0333                                                 dummyMomentum);
0334   BOOST_CHECK_EQUAL(globalPosition, onPositive);
0335   localPosition =
0336       alignedSurface.globalToLocal(positiveContext, onPositive, dummyMomentum)
0337           .value();
0338   BOOST_CHECK_EQUAL(localPosition, Vector2(3., 3.));
0339 }
0340
0341 BOOST_AUTO_TEST_CASE(AlignVolumeTests) {
0342   Transform3 volTrf1{Transform3::Identity()};
0343   volTrf1.translation() = Vector3{100._mm, 200._mm, 300._mm};
0344   AlignableVolumePlacement volumePlacement{volTrf1};
0345
0346   /// define the diamond volume bounds
0347   constexpr double halfX1 = 10.0_cm;
0348   constexpr double halfX2 = 12.0_cm;
0349   constexpr double halfX3 = 12.0_cm;
0350   constexpr double halfY1 = 5.0_cm;
0351   constexpr double halfY2 = 10.0_cm;
0352   constexpr double halfZ = 2.0_cm;
0353
0354   /// Create an alignable volume
0355   Volume alignedVol1{volumePlacement,
0356                      std::make_shared<DiamondVolumeBounds>(
0357                          halfX1, halfX2, halfX3, halfY1, halfY2, halfZ)};
0358
0359   BOOST_CHECK_EQUAL(alignedVol1.isAlignable(), true);
0360   BOOST_CHECK_THROW(alignedVol1.setTransform(volTrf1), std::runtime_error);
0361   const AlignmentContext defaultContext{};
0362   // Check that the transform of the alignable volume is defined by the
0363   // placement
0364
0365   BOOST_CHECK_EQUAL(
0366       &alignedVol1.localToGlobalTransform(defaultContext.getContext()),
0367       &volumePlacement.localToGlobalTransform(defaultContext.getContext()));
0368
0369   BOOST_CHECK_EQUAL(
0370       &alignedVol1.globalToLocalTransform(defaultContext.getContext()),
0371       &volumePlacement.globalToLocalTransform(defaultContext.getContext()));
0372
0373   BOOST_CHECK(isSame(
0374       volumePlacement.localToGlobalTransform(defaultContext.getContext()),
0375       volTrf1));
0376   // Then fetch the oriented surfaces
0377   std::vector<OrientedSurface> orientedSurfaces =
0378       alignedVol1.volumeBounds().orientedSurfaces(
0379           alignedVol1.localToGlobalTransform(defaultContext.getContext()));
0380
0381   for (const OrientedSurface& surface : orientedSurfaces) {
0382     BOOST_CHECK_EQUAL(surface.surface->isAlignable(), false);
0383     BOOST_CHECK_EQUAL(surface.surface->isSensitive(), false);
0384   }
0385
0386   std::vector<std::shared_ptr<RegularSurface>> portalSurfaces{};
0387   std::ranges::transform(
0388       orientedSurfaces, std::back_inserter(portalSurfaces),
0389       [](const OrientedSurface& surface) { return surface.surface; });
0390
0391   volumePlacement.makePortalsAlignable(defaultContext.getContext(),
0392                                        portalSurfaces);
0393
0394   // Reset the surfaces
0395   orientedSurfaces = alignedVol1.volumeBounds().orientedSurfaces(
0396       alignedVol1.localToGlobalTransform(defaultContext.getContext()));
0397
0398   for (std::size_t portal = 0ul; portal < portalSurfaces.size(); ++portal) {
0399     // The portal mut be alignable
0400     BOOST_CHECK_EQUAL(portalSurfaces[portal]->isAlignable(), true);
0401     // But not sensitive
0402     BOOST_CHECK_EQUAL(portalSurfaces[portal]->isSensitive(), false);
0403     // Then check that the surface are placed at the same spot
0404     BOOST_CHECK(isSame(orientedSurfaces[portal].surface->localToGlobalTransform(
0405                            defaultContext.getContext()),
0406                        portalSurfaces[portal]->localToGlobalTransform(
0407                            defaultContext.getContext())));
0408     // Also check that their bounds are the same
0409     BOOST_CHECK(orientedSurfaces[portal].surface->bounds() ==
0410                 portalSurfaces[portal]->bounds());
0411   }
0412
0413   AlignmentContext alignedContext{};
0414   const Transform3 rotationDelta{Translation3{0., -200._mm, -300._mm} *
0415                                  AngleAxis3{25._degree, Vector3{1., 0., 0.}}};
0416   ACTS_INFO("Rotation correction: " << toString(rotationDelta));
0417   volumePlacement.setAlignmentDelta(alignedContext, rotationDelta, 0);
0418
0419   const Transform3 alignedTrf =
0420       volumePlacement.localToGlobalTransform(alignedContext.getContext());
0421   const Transform3 assembledTrf = volTrf1 * rotationDelta;
0422   ACTS_INFO("Test whether the aligned volume transform is what's expected \n"
0423             << " ---- " << toString(alignedTrf) << "\n ---- "
0424             << toString(assembledTrf));
0425   BOOST_CHECK(isSame(alignedTrf, assembledTrf));
0426
0427   orientedSurfaces = alignedVol1.volumeBounds().orientedSurfaces(
0428       alignedVol1.localToGlobalTransform(alignedContext.getContext()));
0429
0430   for (std::size_t portal = 0ul; portal < portalSurfaces.size(); ++portal) {
0431     BOOST_CHECK(isSame(orientedSurfaces[portal].surface->localToGlobalTransform(
0432                            alignedContext.getContext()),
0433                        portalSurfaces[portal]->localToGlobalTransform(
0434                            alignedContext.getContext())));
0435   }
0436
0437   // Ensure that the bound values can no longer be updated
0438   BOOST_CHECK_THROW(
0439       alignedVol1.assignVolumeBounds(std::make_shared<DiamondVolumeBounds>(
0440           halfX1, 2. * halfX2, halfX3, halfY1, halfY2, halfZ)),
0441       std::runtime_error);
0442   // But equivalent bounds can be pushed
0443   BOOST_CHECK_NO_THROW(
0444       alignedVol1.assignVolumeBounds(std::make_shared<DiamondVolumeBounds>(
0445           halfX1, halfX2, halfX3, halfY1, halfY2, halfZ)));
0446 }
0447
0448 BOOST_AUTO_TEST_CASE(ConfinedVolumes) {
0449   AlignableVolumePlacement outsidePlacement{Transform3::Identity()};
0450   std::vector<std::unique_ptr<AlignableVolumePlacement>> innerPlacements{};
0451
0452   const AlignmentContext gctx{};
0453   using namespace Acts::Experimental;
0454
0455   constexpr double hX = 10._cm;
0456   constexpr double hY = 20._cm;
0457   constexpr double hZ = 30._cm;
0458
0459   auto innerBounds = std::make_shared<CuboidVolumeBounds>(hX, hY, hZ);
0460   auto outerBounds =
0461       std::make_shared<CuboidVolumeBounds>(10. * hX, 2. * hY, 2. * hZ);
0462
0463   Blueprint::Config cfg;
0464   cfg.envelope[AxisDirection::AxisZ] = {20_m, 20_m};
0465   cfg.envelope[AxisDirection::AxisR] = {0, 20_m};
0466   auto root = std::make_unique<Blueprint>(cfg);
0467
0468   // auto& cubcontainer = root->addCuboidContainer("CuboidContainer",
0469   // AxisDirection::AxisX);
0470   auto parentVol =
0471       std::make_unique<TrackingVolume>(outsidePlacement, outerBounds, "parent");
0472   parentVol->assignGeometryId(GeometryIdentifier{}.withVolume(5));
0473   auto parentNode = std::make_shared<StaticBlueprintNode>(std::move(parentVol));
0474   // start from the edge of the parent volume
0475   const double xShift =
0476       5._mm + innerBounds->get(CuboidVolumeBounds::eHalfLengthX);
0477   double startX = -outerBounds->get(CuboidVolumeBounds::eHalfLengthX) + xShift;
0478
0479   unsigned i = 0;
0480   while (startX < outerBounds->get(CuboidVolumeBounds::eHalfLengthX)) {
0481     const Transform3 trf{Translation3{startX, 0., 0.}};
0482     startX += 2. * xShift;
0483     std::unique_ptr<TrackingVolume> childVol{};
0484     auto& placement = innerPlacements.emplace_back(
0485         std::make_unique<AlignableVolumePlacement>(trf));
0486     childVol = std::make_unique<TrackingVolume>(
0487         *placement, innerBounds, std::format("alignable_child_{:}", i));
0488
0489     ++i;
0490     childVol->assignGeometryId(
0491         GeometryIdentifier{}.withVolume(10).withLayer(i + 1));
0492     auto childNode = std::make_shared<StaticBlueprintNode>(std::move(childVol));
0493     parentNode->addChild(std::move(childNode));
0494   }
0495
0496   root->addChild(std::move(parentNode));
0497   auto trackingGeometry = root->construct({}, gctx.getContext(), logger());
0498   const auto* rootVol =
0499       trackingGeometry->findVolume(GeometryIdentifier{}.withVolume(5));
0500   BOOST_CHECK_NE(rootVol, nullptr);
0501   // Inspection time
0502   {
0503     Acts::detail::TrackingGeometryPrintVisitor printVisitor{gctx.getContext()};
0504     rootVol->apply(printVisitor);
0505     ACTS_INFO("Constructed tracking geometry: \n"
0506               << printVisitor.stream().str());
0507   }
0508   // Ensure that the volume is alignable
0509   BOOST_CHECK_EQUAL(rootVol->isAlignable(), true);
0510   BOOST_CHECK_EQUAL(rootVol->volumePlacement(), &outsidePlacement);
0511   BOOST_CHECK_EQUAL(outsidePlacement.nPortalPlacements(), 6);
0512   for (const auto& [child, placement] :
0513        zip(rootVol->volumes(), innerPlacements)) {
0514     BOOST_CHECK_EQUAL(child.isAlignable(), true);
0515     BOOST_CHECK_EQUAL(child.volumePlacement(), placement.get());
0516     BOOST_CHECK_EQUAL(placement->nPortalPlacements(), 6);
0517   }
0518 }
0519
0520 BOOST_AUTO_TEST_SUITE_END();