From f56610d68372d86fcc13df8fa9383bfbd4117d8b Mon Sep 17 00:00:00 2001 From: Paul Gessinger Date: Fri, 11 Oct 2024 21:23:53 +0200 Subject: [PATCH] feat: Delayed Grid construction for Portals (#3718) This PR does three things: - `PortalLinkBase::merge` **no longer** does merging of grids or trivials into a combined grid. This has been observed to lead to accumulating floating point imprecision. - `CompositePortalLink` gets a method to make a grid **if** (and only if) it is composed of a set of `TrivialPortalLink`s. - `Portal::fuse` will attempt to convert `CompositePortalLink`s composed of only `TrivialPortalLink`s to a single `GridPortalLink` before fusing it with the other portal. In combination, this avoids the floating point issues and produces correctly sized grids. Part of #3502. --- This pull request introduces several enhancements and bug fixes to the `Core/include/Acts/Geometry` module, focusing on improving the `CompositePortalLink` and `GridPortalLink` classes. The most important changes include the addition of new constructors, the introduction of the `makeGrid` method, and various adjustments to ensure compatibility and correctness. ### Enhancements to `CompositePortalLink`: * Added new constructors to allow the creation of composite portals from multiple links and to handle nested composites. (`Core/include/Acts/Geometry/CompositePortalLink.hpp`, [[1]](diffhunk://#diff-248145d68399a17324b82d81d6e661a3ab739eb5b6f8d67f40145195ca465c36R55-R63) [[2]](diffhunk://#diff-5663ec0ea1d9723e610725aa0d0964e5c833bb90f431281c3802a9b9c5fa4314R101-R129) * Introduced the `makeGrid` method to convert composite portal links into grid portal links under specific conditions. (`Core/include/Acts/Geometry/CompositePortalLink.hpp`, [Core/src/Geometry/CompositePortalLink.cppR180-R297](diffhunk://#diff-5663ec0ea1d9723e610725aa0d0964e5c833bb90f431281c3802a9b9c5fa4314R180-R297)) ### Adjustments to `GridPortalLink`: * Changed the type of `atLocalBins` methods to return `const TrackingVolume*` instead of `TrackingVolume*`. (`Core/include/Acts/Geometry/GridPortalLink.hpp`, [[1]](diffhunk://#diff-5cc33f33e4da7753510a3c7bf2481d12c34dd9c1344bfd624c0b5db1d70f214fL384-R389) [[2]](diffhunk://#diff-5cc33f33e4da7753510a3c7bf2481d12c34dd9c1344bfd624c0b5db1d70f214fL402-R402) [[3]](diffhunk://#diff-5cc33f33e4da7753510a3c7bf2481d12c34dd9c1344bfd624c0b5db1d70f214fL548-R547) [[4]](diffhunk://#diff-5cc33f33e4da7753510a3c7bf2481d12c34dd9c1344bfd624c0b5db1d70f214fL560-R559) * Updated internal grid type to use `const TrackingVolume*`. (`Core/include/Acts/Geometry/GridPortalLink.hpp`, [Core/include/Acts/Geometry/GridPortalLink.hppL402-R402](diffhunk://#diff-5cc33f33e4da7753510a3c7bf2481d12c34dd9c1344bfd624c0b5db1d70f214fL402-R402)) ### Bug Fixes and Code Improvements: * Moved `mergedSurface` function to an anonymous namespace in the implementation file and refactored it for better error handling and type safety. (`Core/src/Geometry/CompositePortalLink.cpp`, [[1]](diffhunk://#diff-5663ec0ea1d9723e610725aa0d0964e5c833bb90f431281c3802a9b9c5fa4314R11-R78) [[2]](diffhunk://#diff-5663ec0ea1d9723e610725aa0d0964e5c833bb90f431281c3802a9b9c5fa4314L77-L106) * Improved the `isSameSurface` function to include detailed checks for surface bounds and transformations. (`Core/src/Geometry/Portal.cpp`, [Core/src/Geometry/Portal.cppL308-R357](diffhunk://#diff-e32791625fda93fd367fc971619ea03be19128e91bbca7e8a09b5af399beb461L308-R357)) * Enhanced logging and error messages for better debugging and clarity. (`Core/src/Geometry/Portal.cpp`, [[1]](diffhunk://#diff-e32791625fda93fd367fc971619ea03be19128e91bbca7e8a09b5af399beb461R274-R277) [[2]](diffhunk://#diff-e32791625fda93fd367fc971619ea03be19128e91bbca7e8a09b5af399beb461R293-R313) These changes collectively improve the functionality, safety, and maintainability of the `Acts` geometry module, particularly in handling complex portal link structures. --- .../Acts/Geometry/CompositePortalLink.hpp | 35 ++- Core/include/Acts/Geometry/GridPortalLink.hpp | 11 +- .../Acts/Geometry/TrivialPortalLink.hpp | 4 + Core/src/Geometry/CompositePortalLink.cpp | 232 +++++++++++++++--- Core/src/Geometry/GridPortalLink.cpp | 2 +- Core/src/Geometry/Portal.cpp | 55 ++++- Core/src/Geometry/PortalLinkBase.cpp | 18 +- Core/src/Geometry/TrivialPortalLink.cpp | 4 + .../Core/Geometry/PortalLinkTests.cpp | 207 ++++++++++++---- Tests/UnitTests/Core/Geometry/PortalTests.cpp | 67 ++++- 10 files changed, 521 insertions(+), 114 deletions(-) diff --git a/Core/include/Acts/Geometry/CompositePortalLink.hpp b/Core/include/Acts/Geometry/CompositePortalLink.hpp index 5ae6d76a1060..6586d875b794 100644 --- a/Core/include/Acts/Geometry/CompositePortalLink.hpp +++ b/Core/include/Acts/Geometry/CompositePortalLink.hpp @@ -17,6 +17,9 @@ namespace Acts { +class GridPortalLink; +class Surface; + /// Composite portal links can graft together other portal link instances, for /// example grids that could not be merged due to invalid binnings. /// @@ -49,6 +52,15 @@ class CompositePortalLink final : public PortalLinkBase { std::unique_ptr b, BinningValue direction, bool flatten = true); + /// Construct a composite portal from any number of arbitrary other portal + /// links. The only requirement is that the portal link surfaces are + /// mergeable. + /// @param links The portal links + /// @param direction The binning direction + /// @param flatten If true, the composite will flatten any nested composite + CompositePortalLink(std::vector> links, + BinningValue direction, bool flatten = true); + /// Print the composite portal link /// @param os The output stream void toStream(std::ostream& os) const override; @@ -81,19 +93,22 @@ class CompositePortalLink final : public PortalLinkBase { /// @return The number of children std::size_t size() const; - private: - /// Helper function to construct a merged surface from two portal links along - /// a given direction - /// @param a The first portal link - /// @param b The second portal link - /// @param direction The merging direction - /// @return The merged surface - static std::shared_ptr mergedSurface(const PortalLinkBase* a, - const PortalLinkBase* b, - BinningValue direction); + /// (Potentially) create a grid portal link that represents this composite + /// portal link. + /// @note This only works, if the composite is **flat** and only contains + /// **trivial portal links**. If these preconditions are not met, this + /// function returns a nullptr. + /// @param gctx The geometry context + /// @param logger The logger + /// @return The grid portal link + std::unique_ptr makeGrid(const GeometryContext& gctx, + const Logger& logger) const; + private: boost::container::small_vector, 4> m_children{}; + + BinningValue m_direction; }; } // namespace Acts diff --git a/Core/include/Acts/Geometry/GridPortalLink.hpp b/Core/include/Acts/Geometry/GridPortalLink.hpp index 0f84d496920a..411e6bb9962f 100644 --- a/Core/include/Acts/Geometry/GridPortalLink.hpp +++ b/Core/include/Acts/Geometry/GridPortalLink.hpp @@ -381,12 +381,12 @@ class GridPortalLink : public PortalLinkBase { /// Helper function to get grid bin content in type-eraased way. /// @param indices The bin indices /// @return The tracking volume at the bin - virtual TrackingVolume*& atLocalBins(IndexType indices) = 0; + virtual const TrackingVolume*& atLocalBins(IndexType indices) = 0; /// Helper function to get grid bin content in type-eraased way. /// @param indices The bin indices /// @return The tracking volume at the bin - virtual TrackingVolume* atLocalBins(IndexType indices) const = 0; + virtual const TrackingVolume* atLocalBins(IndexType indices) const = 0; private: BinningValue m_direction; @@ -399,7 +399,7 @@ template class GridPortalLinkT final : public GridPortalLink { public: /// The internal grid type - using GridType = Grid; + using GridType = Grid; /// The dimension of the grid static constexpr std::size_t DIM = sizeof...(Axes); @@ -530,7 +530,6 @@ class GridPortalLinkT final : public GridPortalLink { return m_grid.atPosition(m_projection(position)); } - protected: /// Type erased access to the number of bins /// @return The number of bins in each direction IndexType numLocalBins() const override { @@ -545,7 +544,7 @@ class GridPortalLinkT final : public GridPortalLink { /// Type erased local bin access /// @param indices The bin indices /// @return The tracking volume at the bin - TrackingVolume*& atLocalBins(IndexType indices) override { + const TrackingVolume*& atLocalBins(IndexType indices) override { throw_assert(indices.size() == DIM, "Invalid number of indices"); typename GridType::index_t idx; for (std::size_t i = 0; i < DIM; i++) { @@ -557,7 +556,7 @@ class GridPortalLinkT final : public GridPortalLink { /// Type erased local bin access /// @param indices The bin indices /// @return The tracking volume at the bin - TrackingVolume* atLocalBins(IndexType indices) const override { + const TrackingVolume* atLocalBins(IndexType indices) const override { throw_assert(indices.size() == DIM, "Invalid number of indices"); typename GridType::index_t idx; for (std::size_t i = 0; i < DIM; i++) { diff --git a/Core/include/Acts/Geometry/TrivialPortalLink.hpp b/Core/include/Acts/Geometry/TrivialPortalLink.hpp index 9f99881fdc7d..f68a7e326c68 100644 --- a/Core/include/Acts/Geometry/TrivialPortalLink.hpp +++ b/Core/include/Acts/Geometry/TrivialPortalLink.hpp @@ -58,6 +58,10 @@ class TrivialPortalLink final : public PortalLinkBase { const GeometryContext& gctx, const Vector3& position, double tolerance = s_onSurfaceTolerance) const override; + /// Get the single volume that this trivial portal link is associated with + /// @return The target volume + const TrackingVolume& volume() const; + private: TrackingVolume* m_volume; }; diff --git a/Core/src/Geometry/CompositePortalLink.cpp b/Core/src/Geometry/CompositePortalLink.cpp index 7d80e3b331e9..d8010bbbce24 100644 --- a/Core/src/Geometry/CompositePortalLink.cpp +++ b/Core/src/Geometry/CompositePortalLink.cpp @@ -8,21 +8,74 @@ #include "Acts/Geometry/CompositePortalLink.hpp" +#include "Acts/Geometry/GridPortalLink.hpp" #include "Acts/Geometry/PortalError.hpp" +#include "Acts/Geometry/TrivialPortalLink.hpp" #include "Acts/Surfaces/CylinderSurface.hpp" #include "Acts/Surfaces/DiscSurface.hpp" #include "Acts/Surfaces/PlaneSurface.hpp" +#include "Acts/Surfaces/RadialBounds.hpp" #include "Acts/Surfaces/RegularSurface.hpp" +#include "Acts/Utilities/Axis.hpp" +#include #include +#include #include +#include + namespace Acts { +namespace { +std::shared_ptr mergedSurface(const Surface& a, + const Surface& b, + BinningValue direction) { + if (a.type() != b.type()) { + throw std::invalid_argument{"Cannot merge surfaces of different types"}; + } + + if (const auto* cylinderA = dynamic_cast(&a); + cylinderA != nullptr) { + const auto& cylinderB = dynamic_cast(b); + + auto [merged, reversed] = cylinderA->mergedWith(cylinderB, direction, true); + return merged; + } else if (const auto* discA = dynamic_cast(&a); + discA != nullptr) { + const auto& discB = dynamic_cast(b); + auto [merged, reversed] = discA->mergedWith(discB, direction, true); + return merged; + } else if (const auto* planeA = dynamic_cast(&a); + planeA != nullptr) { + throw std::logic_error{"Plane surfaces not implemented yet"}; + } else { + throw std::invalid_argument{"Unsupported surface type"}; + } +} + +std::shared_ptr mergePortalLinks( + const std::vector>& links, + BinningValue direction) { + assert(std::ranges::all_of(links, + [](const auto& link) { return link != nullptr; })); + assert(!links.empty()); + + std::shared_ptr result = links.front()->surfacePtr(); + for (auto it = std::next(links.begin()); it != links.end(); ++it) { + assert(result != nullptr); + result = mergedSurface(*result, it->get()->surface(), direction); + } + + return result; +} +} // namespace + CompositePortalLink::CompositePortalLink(std::unique_ptr a, std::unique_ptr b, BinningValue direction, bool flatten) - : PortalLinkBase(mergedSurface(a.get(), b.get(), direction)) { + : PortalLinkBase(mergedSurface(a->surface(), b->surface(), direction)), + m_direction{direction} { if (!flatten) { m_children.push_back(std::move(a)); m_children.push_back(std::move(b)); @@ -45,6 +98,35 @@ CompositePortalLink::CompositePortalLink(std::unique_ptr a, } } +CompositePortalLink::CompositePortalLink( + std::vector> links, BinningValue direction, + bool flatten) + : PortalLinkBase(mergePortalLinks(links, direction)), + m_direction(direction) { + if (!flatten) { + for (auto& child : links) { + m_children.push_back(std::move(child)); + } + + } else { + auto handle = [&](std::unique_ptr link) { + if (auto* composite = dynamic_cast(link.get()); + composite != nullptr) { + m_children.insert( + m_children.end(), + std::make_move_iterator(composite->m_children.begin()), + std::make_move_iterator(composite->m_children.end())); + } else { + m_children.push_back(std::move(link)); + } + }; + + for (auto& child : links) { + handle(std::move(child)); + } + } +} + Result CompositePortalLink::resolveVolume( const GeometryContext& gctx, const Vector2& position, double tolerance) const { @@ -74,36 +156,6 @@ Result CompositePortalLink::resolveVolume( return PortalError::PositionNotOnAnyChildPortalLink; } -std::shared_ptr CompositePortalLink::mergedSurface( - const PortalLinkBase* a, const PortalLinkBase* b, BinningValue direction) { - assert(a != nullptr); - assert(b != nullptr); - if (a->surface().type() != b->surface().type()) { - throw std::invalid_argument{"Cannot merge surfaces of different types"}; - } - - if (const auto* cylinderA = - dynamic_cast(&a->surface()); - cylinderA != nullptr) { - const auto& cylinderB = dynamic_cast(b->surface()); - - auto [merged, reversed] = cylinderA->mergedWith(cylinderB, direction, true); - return merged; - } else if (const auto* discA = - dynamic_cast(&a->surface()); - discA != nullptr) { - const auto& discB = dynamic_cast(b->surface()); - auto [merged, reversed] = discA->mergedWith(discB, direction, true); - return merged; - } else if (const auto* planeA = - dynamic_cast(&a->surface()); - planeA != nullptr) { - throw std::logic_error{"Plane surfaces not implemented yet"}; - } else { - throw std::invalid_argument{"Unsupported surface type"}; - } -} - std::size_t CompositePortalLink::depth() const { std::size_t result = 1; @@ -125,4 +177,122 @@ void CompositePortalLink::toStream(std::ostream& os) const { os << "CompositePortalLink"; } +std::unique_ptr CompositePortalLink::makeGrid( + const GeometryContext& gctx, const Logger& logger) const { + ACTS_VERBOSE("Attempting to make a grid from a composite portal link"); + std::vector trivialLinks; + std::ranges::transform(m_children, std::back_inserter(trivialLinks), + [](const auto& child) { + return dynamic_cast(child.get()); + }); + + if (std::ranges::any_of(trivialLinks, + [](auto link) { return link == nullptr; })) { + ACTS_VERBOSE( + "Failed to make a grid from a composite portal link -> returning " + "nullptr"); + return nullptr; + } + + // we already know all children have surfaces that are compatible, we produced + // a merged surface for the overall dimensions. + + auto printEdges = [](const auto& edges) -> std::string { + std::vector strings; + std::ranges::transform(edges, std::back_inserter(strings), + [](const auto& v) { return std::to_string(v); }); + return boost::algorithm::join(strings, ", "); + }; + + if (surface().type() == Surface::SurfaceType::Cylinder) { + ACTS_VERBOSE("Combining composite into cylinder grid"); + + if (m_direction != BinningValue::binZ) { + ACTS_ERROR("Cylinder grid only supports binning in Z direction"); + throw std::runtime_error{"Unsupported binning direction"}; + } + + std::vector edges; + edges.reserve(m_children.size() + 1); + + const Transform3& groupTransform = m_surface->transform(gctx); + Transform3 itransform = groupTransform.inverse(); + + std::ranges::sort( + trivialLinks, [&itransform, &gctx](const auto& a, const auto& b) { + return (itransform * a->surface().transform(gctx)).translation()[eZ] < + (itransform * b->surface().transform(gctx)).translation()[eZ]; + }); + + for (const auto& [i, child] : enumerate(trivialLinks)) { + const auto& bounds = + dynamic_cast(child->surface().bounds()); + Transform3 ltransform = itransform * child->surface().transform(gctx); + ActsScalar hlZ = bounds.get(CylinderBounds::eHalfLengthZ); + ActsScalar minZ = ltransform.translation()[eZ] - hlZ; + ActsScalar maxZ = ltransform.translation()[eZ] + hlZ; + if (i == 0) { + edges.push_back(minZ); + } + edges.push_back(maxZ); + } + + ACTS_VERBOSE("~> Determined bin edges to be " << printEdges(edges)); + + Axis axis{AxisBound, edges}; + + auto grid = GridPortalLink::make(m_surface, m_direction, std::move(axis)); + for (const auto& [i, child] : enumerate(trivialLinks)) { + grid->atLocalBins({i + 1}) = &child->volume(); + } + + return grid; + + } else if (surface().type() == Surface::SurfaceType::Disc) { + ACTS_VERBOSE("Combining composite into disc grid"); + + if (m_direction != BinningValue::binR) { + ACTS_ERROR("Disc grid only supports binning in R direction"); + throw std::runtime_error{"Unsupported binning direction"}; + } + + std::vector edges; + edges.reserve(m_children.size() + 1); + + std::ranges::sort(trivialLinks, [](const auto& a, const auto& b) { + const auto& boundsA = + dynamic_cast(a->surface().bounds()); + const auto& boundsB = + dynamic_cast(b->surface().bounds()); + return boundsA.get(RadialBounds::eMinR) < + boundsB.get(RadialBounds::eMinR); + }); + + for (const auto& [i, child] : enumerate(trivialLinks)) { + const auto& bounds = + dynamic_cast(child->surface().bounds()); + + if (i == 0) { + edges.push_back(bounds.get(RadialBounds::eMinR)); + } + edges.push_back(bounds.get(RadialBounds::eMaxR)); + } + + ACTS_VERBOSE("~> Determined bin edges to be " << printEdges(edges)); + + Axis axis{AxisBound, edges}; + + auto grid = GridPortalLink::make(m_surface, m_direction, std::move(axis)); + for (const auto& [i, child] : enumerate(trivialLinks)) { + grid->atLocalBins({i + 1}) = &child->volume(); + } + + return grid; + } else if (surface().type() == Surface::SurfaceType::Plane) { + throw std::runtime_error{"Plane surfaces not implemented yet"}; + } else { + throw std::invalid_argument{"Unsupported surface type"}; + } +} + } // namespace Acts diff --git a/Core/src/Geometry/GridPortalLink.cpp b/Core/src/Geometry/GridPortalLink.cpp index 5e724b631209..57a3e1d31419 100644 --- a/Core/src/Geometry/GridPortalLink.cpp +++ b/Core/src/Geometry/GridPortalLink.cpp @@ -298,7 +298,7 @@ void GridPortalLink::fillGrid1dTo2d(FillDirection dir, assert(locDest.size() == 2); for (std::size_t i = 0; i <= locSource[0] + 1; ++i) { - TrackingVolume* source = grid1d.atLocalBins({i}); + const TrackingVolume* source = grid1d.atLocalBins({i}); if (dir == FillDirection::loc1) { for (std::size_t j = 0; j <= locDest[1] + 1; ++j) { diff --git a/Core/src/Geometry/Portal.cpp b/Core/src/Geometry/Portal.cpp index 2a9d09b5f68d..69cab13b8e00 100644 --- a/Core/src/Geometry/Portal.cpp +++ b/Core/src/Geometry/Portal.cpp @@ -9,13 +9,17 @@ #include "Acts/Geometry/Portal.hpp" #include "Acts/Definitions/Algebra.hpp" +#include "Acts/Definitions/Tolerance.hpp" +#include "Acts/Geometry/CompositePortalLink.hpp" #include "Acts/Geometry/GeometryContext.hpp" +#include "Acts/Geometry/GridPortalLink.hpp" #include "Acts/Geometry/PortalLinkBase.hpp" #include "Acts/Geometry/TrivialPortalLink.hpp" #include "Acts/Surfaces/RegularSurface.hpp" #include "Acts/Utilities/BinningType.hpp" +#include "Acts/Utilities/Zip.hpp" -#include +#include #include #include @@ -245,7 +249,7 @@ Portal Portal::fuse(const GeometryContext& gctx, Portal& aPortal, Portal& bPortal, const Logger& logger) { ACTS_VERBOSE("Fusing two portals"); if (&aPortal == &bPortal) { - ACTS_ERROR("Cannot merge a portal with itself"); + ACTS_ERROR("Cannot fuse a portal with itself"); throw PortalMergingException{}; } @@ -267,6 +271,10 @@ Portal Portal::fuse(const GeometryContext& gctx, Portal& aPortal, if (!isSameSurface(gctx, *aPortal.m_surface, *bPortal.m_surface)) { ACTS_ERROR("Portals have different surfaces"); + ACTS_ERROR("A: " << aPortal.m_surface->bounds()); + ACTS_ERROR("\n" << aPortal.m_surface->transform(gctx).matrix()); + ACTS_ERROR("B: " << bPortal.m_surface->bounds()); + ACTS_ERROR("\n" << bPortal.m_surface->transform(gctx).matrix()); throw PortalFusingException(); } @@ -282,16 +290,27 @@ Portal Portal::fuse(const GeometryContext& gctx, Portal& aPortal, throw PortalFusingException(); } + auto maybeConvertToGrid = [&](std::unique_ptr link) + -> std::unique_ptr { + auto* composite = dynamic_cast(link.get()); + if (composite == nullptr) { + return link; + } + + ACTS_VERBOSE("Converting composite to grid during portal fusing"); + return composite->makeGrid(gctx, logger); + }; + aPortal.m_surface.reset(); bPortal.m_surface.reset(); if (aHasAlongNormal) { ACTS_VERBOSE("Taking along normal from lhs, opposite normal from rhs"); - return Portal{gctx, std::move(aPortal.m_alongNormal), - std::move(bPortal.m_oppositeNormal)}; + return Portal{gctx, maybeConvertToGrid(std::move(aPortal.m_alongNormal)), + maybeConvertToGrid(std::move(bPortal.m_oppositeNormal))}; } else { ACTS_VERBOSE("Taking along normal from rhs, opposite normal from lhs"); - return Portal{gctx, std::move(bPortal.m_alongNormal), - std::move(aPortal.m_oppositeNormal)}; + return Portal{gctx, maybeConvertToGrid(std::move(bPortal.m_alongNormal)), + maybeConvertToGrid(std::move(aPortal.m_oppositeNormal))}; } } @@ -305,12 +324,30 @@ bool Portal::isSameSurface(const GeometryContext& gctx, const Surface& a, return false; } - if (a.bounds() != b.bounds()) { + std::vector aValues = a.bounds().values(); + std::vector bValues = b.bounds().values(); + bool different = false; + for (auto [aVal, bVal] : zip(aValues, bValues)) { + if (std::abs(aVal - bVal) > s_onSurfaceTolerance) { + different = true; + break; + } + } + + if (a.bounds().type() != b.bounds().type() || different) { return false; } - if (!a.transform(gctx).isApprox(b.transform(gctx), - s_transformEquivalentTolerance)) { + if (!a.transform(gctx).linear().isApprox(b.transform(gctx).linear(), + s_transformEquivalentTolerance)) { + return false; + } + + Vector3 delta = + (a.transform(gctx).translation() - b.transform(gctx).translation()) + .cwiseAbs(); + + if (delta.maxCoeff() > s_onSurfaceTolerance) { return false; } diff --git a/Core/src/Geometry/PortalLinkBase.cpp b/Core/src/Geometry/PortalLinkBase.cpp index 459c144272e2..bdc65f8a5520 100644 --- a/Core/src/Geometry/PortalLinkBase.cpp +++ b/Core/src/Geometry/PortalLinkBase.cpp @@ -109,8 +109,10 @@ std::unique_ptr PortalLinkBase::merge( } else if (const auto* bTrivial = dynamic_cast(b.get()); bTrivial != nullptr) { - ACTS_VERBOSE("Merging a grid portal with a trivial portal"); - return gridMerge(*aGrid, *bTrivial->makeGrid(direction)); + ACTS_VERBOSE( + "Merging a grid portal with a trivial portal (via composite)"); + return std::make_unique(std::move(a), std::move(b), + direction); } else if (dynamic_cast(b.get()) != nullptr) { ACTS_WARNING("Merging a grid portal with a composite portal"); @@ -126,15 +128,17 @@ std::unique_ptr PortalLinkBase::merge( aTrivial != nullptr) { if (const auto* bGrid = dynamic_cast(b.get()); bGrid) { - ACTS_VERBOSE("Merging a trivial portal with a grid portal"); - return gridMerge(*aTrivial->makeGrid(direction), *bGrid); + ACTS_VERBOSE( + "Merging a trivial portal with a grid portal (via composite)"); + return std::make_unique(std::move(a), std::move(b), + direction); } else if (const auto* bTrivial = dynamic_cast(b.get()); bTrivial != nullptr) { - ACTS_VERBOSE("Merging two trivial portals"); - return gridMerge(*aTrivial->makeGrid(direction), - *bTrivial->makeGrid(direction)); + ACTS_VERBOSE("Merging two trivial portals (via composite"); + return std::make_unique(std::move(a), std::move(b), + direction); } else if (dynamic_cast(b.get()) != nullptr) { ACTS_WARNING("Merging a trivial portal with a composite portal"); diff --git a/Core/src/Geometry/TrivialPortalLink.cpp b/Core/src/Geometry/TrivialPortalLink.cpp index d0e9fd70b8af..6976388dbc04 100644 --- a/Core/src/Geometry/TrivialPortalLink.cpp +++ b/Core/src/Geometry/TrivialPortalLink.cpp @@ -42,4 +42,8 @@ void TrivialPortalLink::toStream(std::ostream& os) const { os << "TrivialPortalLinkvolumeName() << ">"; } +const TrackingVolume& TrivialPortalLink::volume() const { + return *m_volume; +} + } // namespace Acts diff --git a/Tests/UnitTests/Core/Geometry/PortalLinkTests.cpp b/Tests/UnitTests/Core/Geometry/PortalLinkTests.cpp index dc032060a5bf..e5345268e848 100644 --- a/Tests/UnitTests/Core/Geometry/PortalLinkTests.cpp +++ b/Tests/UnitTests/Core/Geometry/PortalLinkTests.cpp @@ -6,6 +6,7 @@ // License, v. 2.0. If a copy of the MPL was not distributed with this // file, You can obtain one at https://mozilla.org/MPL/2.0/. +#include #include #include #include @@ -2135,7 +2136,161 @@ BOOST_AUTO_TEST_CASE(CompositeConstruction) { BOOST_AUTO_TEST_SUITE(PortalMerging) -BOOST_AUTO_TEST_CASE(TrivialTrivial) {} +BOOST_DATA_TEST_CASE(TrivialTrivial, + (boost::unit_test::data::make(0, -135, 180, 45) * + boost::unit_test::data::make(Vector3{0_mm, 0_mm, 0_mm}, + Vector3{20_mm, 0_mm, 0_mm}, + Vector3{0_mm, 20_mm, 0_mm}, + Vector3{20_mm, 20_mm, 0_mm}, + Vector3{0_mm, 0_mm, 20_mm})), + angle, offset) { + Transform3 base = + AngleAxis3(angle * 1_degree, Vector3::UnitX()) * Translation3(offset); + + BOOST_TEST_CONTEXT("RDirection") { + auto vol1 = std::make_shared( + base, std::make_shared(30_mm, 40_mm, 100_mm)); + + auto vol2 = std::make_shared( + base, std::make_shared(30_mm, 40_mm, 100_mm)); + + auto vol3 = std::make_shared( + base, std::make_shared(40_mm, 50_mm, 100_mm)); + + auto disc1 = + Surface::makeShared(Transform3::Identity(), 30_mm, 60_mm); + + auto disc2 = + Surface::makeShared(Transform3::Identity(), 60_mm, 90_mm); + + auto disc3 = + Surface::makeShared(Transform3::Identity(), 90_mm, 120_mm); + + auto trivial1 = std::make_unique(disc1, *vol1); + BOOST_REQUIRE(trivial1); + auto trivial2 = std::make_unique(disc2, *vol2); + BOOST_REQUIRE(trivial2); + auto trivial3 = std::make_unique(disc3, *vol3); + BOOST_REQUIRE(trivial3); + + auto grid1 = trivial1->makeGrid(BinningValue::binR); + auto compGridTrivial = PortalLinkBase::merge( + std::move(grid1), std::make_unique(*trivial2), + BinningValue::binR, *logger); + BOOST_REQUIRE(compGridTrivial); + BOOST_CHECK_EQUAL(dynamic_cast(*compGridTrivial) + .makeGrid(gctx, *logger), + nullptr); + + auto composite = PortalLinkBase::merge( + std::move(trivial1), std::move(trivial2), BinningValue::binR, *logger); + BOOST_REQUIRE(composite); + + auto grid12 = + dynamic_cast(*composite).makeGrid(gctx, *logger); + BOOST_REQUIRE(grid12); + + BOOST_CHECK_EQUAL( + grid12->resolveVolume(gctx, Vector2{40_mm, 0_degree}).value(), + vol1.get()); + + BOOST_CHECK_EQUAL( + grid12->resolveVolume(gctx, Vector2{70_mm, 0_degree}).value(), + vol2.get()); + + composite = PortalLinkBase::merge(std::move(composite), std::move(trivial3), + BinningValue::binR, *logger); + BOOST_REQUIRE(composite); + + auto grid123 = + dynamic_cast(*composite).makeGrid(gctx, *logger); + BOOST_REQUIRE(grid123); + + BOOST_CHECK_EQUAL( + grid123->resolveVolume(gctx, Vector2{40_mm, 0_degree}).value(), + vol1.get()); + + BOOST_CHECK_EQUAL( + grid123->resolveVolume(gctx, Vector2{70_mm, 0_degree}).value(), + vol2.get()); + + BOOST_CHECK_EQUAL( + grid123->resolveVolume(gctx, Vector2{100_mm, 0_degree}).value(), + vol3.get()); + } + + BOOST_TEST_CONTEXT("ZDirection") { + auto vol1 = std::make_shared( + base, std::make_shared(30_mm, 40_mm, 100_mm)); + + auto vol2 = std::make_shared( + base * Translation3{Vector3::UnitZ() * 200}, + std::make_shared(30_mm, 40_mm, 100_mm)); + + auto vol3 = std::make_shared( + base * Translation3{Vector3::UnitZ() * 400}, + std::make_shared(30_mm, 40_mm, 100_mm)); + + auto cyl1 = Surface::makeShared(base, 40_mm, 100_mm); + + auto cyl2 = Surface::makeShared( + base * Translation3{Vector3::UnitZ() * 200}, 40_mm, 100_mm); + + auto cyl3 = Surface::makeShared( + base * Translation3{Vector3::UnitZ() * 400}, 40_mm, 100_mm); + + auto trivial1 = std::make_unique(cyl1, *vol1); + BOOST_REQUIRE(trivial1); + auto trivial2 = std::make_unique(cyl2, *vol2); + BOOST_REQUIRE(trivial2); + auto trivial3 = std::make_unique(cyl3, *vol3); + BOOST_REQUIRE(trivial3); + + auto grid1 = trivial1->makeGrid(BinningValue::binZ); + auto compGridTrivial = PortalLinkBase::merge( + std::move(grid1), std::make_unique(*trivial2), + BinningValue::binZ, *logger); + BOOST_REQUIRE(compGridTrivial); + BOOST_CHECK_EQUAL(dynamic_cast(*compGridTrivial) + .makeGrid(gctx, *logger), + nullptr); + + auto composite = PortalLinkBase::merge( + std::move(trivial1), std::move(trivial2), BinningValue::binZ, *logger); + BOOST_REQUIRE(composite); + + auto grid12 = + dynamic_cast(*composite).makeGrid(gctx, *logger); + BOOST_REQUIRE(grid12); + + BOOST_CHECK_EQUAL( + grid12->resolveVolume(gctx, Vector2{40_mm, -40_mm}).value(), + vol1.get()); + + BOOST_CHECK_EQUAL( + grid12->resolveVolume(gctx, Vector2{40_mm, 40_mm}).value(), vol2.get()); + + composite = PortalLinkBase::merge(std::move(composite), std::move(trivial3), + BinningValue::binZ, *logger); + BOOST_REQUIRE(composite); + + auto grid123 = + dynamic_cast(*composite).makeGrid(gctx, *logger); + BOOST_REQUIRE(grid123); + + BOOST_CHECK_EQUAL( + grid123->resolveVolume(gctx, Vector2{40_mm, -110_mm}).value(), + vol1.get()); + + BOOST_CHECK_EQUAL( + grid123->resolveVolume(gctx, Vector2{40_mm, -10_mm}).value(), + vol2.get()); + + BOOST_CHECK_EQUAL( + grid123->resolveVolume(gctx, Vector2{40_mm, 190_mm}).value(), + vol3.get()); + } +} BOOST_AUTO_TEST_CASE(TrivialGridR) { auto vol1 = std::make_shared( @@ -2167,32 +2322,14 @@ BOOST_AUTO_TEST_CASE(TrivialGridR) { auto merged = PortalLinkBase::merge(copy(trivial), copy(gridR), BinningValue::binR, *logger); BOOST_REQUIRE(merged); - - auto* mergedGrid = dynamic_cast(merged.get()); - BOOST_REQUIRE(mergedGrid); - - mergedGrid->printContents(std::cout); - - BOOST_CHECK_EQUAL(mergedGrid->grid().axes().size(), 1); - Axis axisExpected{AxisBound, {30_mm, 45_mm, 60_mm, 90_mm}}; - BOOST_CHECK_EQUAL(*mergedGrid->grid().axes().front(), axisExpected); + BOOST_CHECK_NE(dynamic_cast(merged.get()), nullptr); } BOOST_TEST_CONTEXT("Orthogonal") { auto merged = PortalLinkBase::merge(copy(gridPhi), copy(trivial), BinningValue::binR, *logger); BOOST_REQUIRE(merged); - - auto* mergedGrid = dynamic_cast(merged.get()); - BOOST_REQUIRE(mergedGrid); - - mergedGrid->printContents(std::cout); - - BOOST_CHECK_EQUAL(mergedGrid->grid().axes().size(), 2); - Axis axis1Expected{AxisBound, 30_mm, 90_mm, 2}; - Axis axis2Expected{AxisClosed, -M_PI, M_PI, 2}; - BOOST_CHECK_EQUAL(*mergedGrid->grid().axes().front(), axis1Expected); - BOOST_CHECK_EQUAL(*mergedGrid->grid().axes().back(), axis2Expected); + BOOST_CHECK_NE(dynamic_cast(merged.get()), nullptr); } } @@ -2227,38 +2364,14 @@ BOOST_AUTO_TEST_CASE(TrivialGridPhi) { auto merged = PortalLinkBase::merge(copy(trivial), copy(gridPhi), BinningValue::binPhi, *logger); BOOST_REQUIRE(merged); - - auto* mergedGrid = dynamic_cast(merged.get()); - BOOST_REQUIRE(mergedGrid); - - mergedGrid->printContents(std::cout); - - BOOST_CHECK_EQUAL(mergedGrid->grid().axes().size(), 1); - const auto& axis = *mergedGrid->grid().axes().front(); - BOOST_CHECK_EQUAL(axis.getType(), AxisType::Variable); - BOOST_CHECK_EQUAL(axis.getBoundaryType(), AxisBoundaryType::Bound); - std::vector expectedBins{-90_degree, 30_degree, 60_degree, - 90_degree}; - CHECK_CLOSE_REL(axis.getBinEdges(), expectedBins, 1e-6); + BOOST_CHECK_NE(dynamic_cast(merged.get()), nullptr); } BOOST_TEST_CONTEXT("Orthogonal") { auto merged = PortalLinkBase::merge(copy(gridR), copy(trivial), BinningValue::binPhi, *logger); BOOST_REQUIRE(merged); - - auto* mergedGrid = dynamic_cast(merged.get()); - BOOST_REQUIRE(mergedGrid); - - mergedGrid->printContents(std::cout); - - BOOST_CHECK_EQUAL(mergedGrid->grid().axes().size(), 2); - const auto& axis1 = *mergedGrid->grid().axes().front(); - const auto& axis2 = *mergedGrid->grid().axes().back(); - Axis axis1Expected{AxisBound, 30_mm, 100_mm, 2}; - BOOST_CHECK_EQUAL(axis1, axis1Expected); - std::vector expectedBins{-90_degree, 30_degree, 90_degree}; - CHECK_CLOSE_REL(axis2.getBinEdges(), expectedBins, 1e-6); + BOOST_CHECK_NE(dynamic_cast(merged.get()), nullptr); } } diff --git a/Tests/UnitTests/Core/Geometry/PortalTests.cpp b/Tests/UnitTests/Core/Geometry/PortalTests.cpp index 490817486498..8997ccedda7a 100644 --- a/Tests/UnitTests/Core/Geometry/PortalTests.cpp +++ b/Tests/UnitTests/Core/Geometry/PortalTests.cpp @@ -12,6 +12,7 @@ #include #include "Acts/Definitions/Units.hpp" +#include "Acts/Geometry/CompositePortalLink.hpp" #include "Acts/Geometry/CylinderVolumeBounds.hpp" #include "Acts/Geometry/GeometryContext.hpp" #include "Acts/Geometry/GridPortalLink.hpp" @@ -171,10 +172,9 @@ BOOST_AUTO_TEST_CASE(Cylinder) { BOOST_CHECK_NE(merged12.getLink(Direction::AlongNormal), nullptr); BOOST_CHECK_EQUAL(merged12.getLink(Direction::OppositeNormal), nullptr); - auto grid12 = dynamic_cast( + auto composite12 = dynamic_cast( merged12.getLink(Direction::AlongNormal)); - BOOST_REQUIRE_NE(grid12, nullptr); - grid12->printContents(std::cout); + BOOST_REQUIRE_NE(composite12, nullptr); BOOST_CHECK_EQUAL( merged12 @@ -528,6 +528,67 @@ BOOST_AUTO_TEST_CASE(Material) { BOOST_CHECK_EQUAL(portal12.surface().surfaceMaterial(), material.get()); } +BOOST_AUTO_TEST_CASE(GridCreationOnFuse) { + Transform3 base = Transform3::Identity(); + + auto vol1 = std::make_shared( + base, std::make_shared(30_mm, 40_mm, 100_mm)); + + auto vol2 = std::make_shared( + base, std::make_shared(30_mm, 40_mm, 100_mm)); + + auto vol3 = std::make_shared( + base, std::make_shared(40_mm, 50_mm, 100_mm)); + + auto vol4 = std::make_shared( + base, std::make_shared(40_mm, 50_mm, 100_mm)); + + auto disc1 = + Surface::makeShared(Transform3::Identity(), 30_mm, 60_mm); + + auto disc2 = + Surface::makeShared(Transform3::Identity(), 60_mm, 90_mm); + + auto disc3 = + Surface::makeShared(Transform3::Identity(), 90_mm, 120_mm); + + auto trivial1 = std::make_unique(disc1, *vol1); + BOOST_REQUIRE(trivial1); + auto trivial2 = std::make_unique(disc2, *vol2); + BOOST_REQUIRE(trivial2); + auto trivial3 = std::make_unique(disc3, *vol3); + BOOST_REQUIRE(trivial3); + + std::vector> links; + links.push_back(std::move(trivial1)); + links.push_back(std::move(trivial2)); + links.push_back(std::move(trivial3)); + + auto composite = std::make_unique(std::move(links), + BinningValue::binR); + + auto discOpposite = + Surface::makeShared(Transform3::Identity(), 30_mm, 120_mm); + + auto trivialOpposite = + std::make_unique(discOpposite, *vol4); + + Portal aPortal{gctx, std::move(composite), nullptr}; + Portal bPortal{gctx, nullptr, std::move(trivialOpposite)}; + + Portal fused = Portal::fuse(gctx, aPortal, bPortal, *logger); + + BOOST_CHECK_NE(dynamic_cast( + fused.getLink(Direction::OppositeNormal)), + nullptr); + + const auto* grid = dynamic_cast( + fused.getLink(Direction::AlongNormal)); + BOOST_REQUIRE_NE(grid, nullptr); + + BOOST_CHECK_EQUAL(grid->grid().axes().front()->getNBins(), 3); +} + BOOST_AUTO_TEST_SUITE_END() // Fusing BOOST_AUTO_TEST_CASE(Construction) {