Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2025-10-19 07:58:54

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/data/test_case.hpp>
0010 #include <boost/test/tools/old/interface.hpp>
0011 #include <boost/test/unit_test.hpp>
0012 #include <boost/test/unit_test_suite.hpp>
0013 
0014 #include "Acts/Definitions/Units.hpp"
0015 #include "Acts/Geometry/Blueprint.hpp"
0016 #include "Acts/Geometry/BlueprintNode.hpp"
0017 #include "Acts/Geometry/ContainerBlueprintNode.hpp"
0018 #include "Acts/Geometry/CylinderVolumeBounds.hpp"
0019 #include "Acts/Geometry/CylinderVolumeStack.hpp"
0020 #include "Acts/Geometry/GeometryContext.hpp"
0021 #include "Acts/Geometry/LayerBlueprintNode.hpp"
0022 #include "Acts/Geometry/MaterialDesignatorBlueprintNode.hpp"
0023 #include "Acts/Geometry/StaticBlueprintNode.hpp"
0024 #include "Acts/Geometry/TrackingVolume.hpp"
0025 #include "Acts/Geometry/TrapezoidVolumeBounds.hpp"
0026 #include "Acts/Geometry/VolumeAttachmentStrategy.hpp"
0027 #include "Acts/Material/BinnedSurfaceMaterial.hpp"
0028 #include "Acts/Material/HomogeneousSurfaceMaterial.hpp"
0029 #include "Acts/Material/Material.hpp"
0030 #include "Acts/Material/MaterialSlab.hpp"
0031 #include "Acts/Material/ProtoSurfaceMaterial.hpp"
0032 #include "Acts/Surfaces/RectangleBounds.hpp"
0033 #include "Acts/Utilities/BinningType.hpp"
0034 #include "Acts/Utilities/Logger.hpp"
0035 #include "Acts/Utilities/ProtoAxis.hpp"
0036 #include "ActsTests/CommonHelpers/DetectorElementStub.hpp"
0037 
0038 #include <fstream>
0039 #include <memory>
0040 #include <stdexcept>
0041 #include <vector>
0042 
0043 using namespace Acts;
0044 using namespace Acts::UnitLiterals;
0045 using Experimental::Blueprint;
0046 using Experimental::BlueprintNode;
0047 using Experimental::BlueprintOptions;
0048 using Experimental::LayerBlueprintNode;
0049 using Experimental::MaterialDesignatorBlueprintNode;
0050 using Experimental::StaticBlueprintNode;
0051 
0052 namespace ActsTests {
0053 
0054 auto logger = getDefaultLogger("UnitTests", Logging::VERBOSE);
0055 
0056 GeometryContext gctx;
0057 
0058 auto nameLookup(const TrackingGeometry& geo) {
0059   return [&](const std::string& name) -> const TrackingVolume& {
0060     const TrackingVolume* volume = nullptr;
0061 
0062     geo.visitVolumes([&](const TrackingVolume* v) {
0063       if (v->volumeName() == name) {
0064         volume = v;
0065       }
0066     });
0067 
0068     if (volume == nullptr) {
0069       throw std::runtime_error("Volume not found: " + name);
0070     }
0071     return *volume;
0072   };
0073 }
0074 
0075 std::size_t countVolumes(const TrackingGeometry& geo) {
0076   std::size_t nVolumes = 0;
0077   geo.visitVolumes([&](const TrackingVolume* /*volume*/) { nVolumes++; });
0078   return nVolumes;
0079 }
0080 
0081 BOOST_AUTO_TEST_SUITE(GeometrySuite);
0082 
0083 BOOST_AUTO_TEST_CASE(InvalidRoot) {
0084   Logging::ScopedFailureThreshold threshold{Logging::Level::FATAL};
0085 
0086   Blueprint::Config cfg;
0087   Blueprint root{cfg};
0088   BOOST_CHECK_THROW(root.construct({}, gctx, *logger), std::logic_error);
0089 
0090   // More than one child is also invalid
0091   auto cylBounds = std::make_shared<CylinderVolumeBounds>(10_mm, 20_mm, 100_mm);
0092   root.addChild(
0093       std::make_unique<StaticBlueprintNode>(std::make_unique<TrackingVolume>(
0094           Transform3::Identity(), cylBounds, "child1")));
0095   root.addChild(
0096       std::make_unique<StaticBlueprintNode>(std::make_unique<TrackingVolume>(
0097           Transform3::Identity(), cylBounds, "child2")));
0098 
0099   BOOST_CHECK_THROW(root.construct({}, gctx, *logger), std::logic_error);
0100 }
0101 
0102 class DummyNode : public BlueprintNode {
0103  public:
0104   explicit DummyNode(const std::string& name) : m_name(name) {}
0105 
0106   const std::string& name() const override { return m_name; }
0107 
0108   Volume& build(const BlueprintOptions& /*options*/,
0109                 const GeometryContext& /*gctx*/,
0110                 const Logger& /*logger*/) override {
0111     throw std::logic_error("Not implemented");
0112   }
0113 
0114   PortalShellBase& connect(const BlueprintOptions& /*options*/,
0115                            const GeometryContext& /*gctx*/,
0116                            const Logger& /*logger */) override {
0117     throw std::logic_error("Not implemented");
0118   }
0119 
0120   void finalize(const BlueprintOptions& /*options*/,
0121                 const GeometryContext& /*gctx*/, TrackingVolume& /*parent*/,
0122                 const Logger& /*logger*/) override {
0123     throw std::logic_error("Not implemented");
0124   }
0125 
0126  private:
0127   std::string m_name;
0128 };
0129 
0130 BOOST_AUTO_TEST_CASE(RootCannotBeChild) {
0131   auto node = std::make_shared<DummyNode>("node");
0132   Blueprint::Config cfg;
0133   auto root = std::make_shared<Blueprint>(cfg);
0134 
0135   BOOST_CHECK_THROW(node->addChild(root), std::invalid_argument);
0136 }
0137 
0138 BOOST_AUTO_TEST_CASE(AddChildInvalid) {
0139   auto node = std::make_shared<DummyNode>("node");
0140 
0141   // Add self
0142   BOOST_CHECK_THROW(node->addChild(node), std::invalid_argument);
0143 
0144   // Add nullptr
0145   BOOST_CHECK_THROW(node->addChild(nullptr), std::invalid_argument);
0146 
0147   auto nodeB = std::make_shared<DummyNode>("nodeB");
0148   auto nodeC = std::make_shared<DummyNode>("nodeC");
0149 
0150   node->addChild(nodeB);
0151   nodeB->addChild(nodeC);
0152   BOOST_CHECK_THROW(nodeC->addChild(node), std::invalid_argument);
0153 
0154   // already has parent, can't be added as a child anywhere else
0155   BOOST_CHECK_THROW(node->addChild(nodeC), std::invalid_argument);
0156 }
0157 
0158 BOOST_AUTO_TEST_CASE(Depth) {
0159   auto node1 = std::make_shared<DummyNode>("node1");
0160   auto node2 = std::make_shared<DummyNode>("node2");
0161   auto node3 = std::make_shared<DummyNode>("node3");
0162 
0163   BOOST_CHECK_EQUAL(node1->depth(), 0);
0164   BOOST_CHECK_EQUAL(node2->depth(), 0);
0165   BOOST_CHECK_EQUAL(node3->depth(), 0);
0166 
0167   node2->addChild(node3);
0168   BOOST_CHECK_EQUAL(node2->depth(), 0);
0169   BOOST_CHECK_EQUAL(node3->depth(), 1);
0170 
0171   node1->addChild(node2);
0172   BOOST_CHECK_EQUAL(node1->depth(), 0);
0173   BOOST_CHECK_EQUAL(node2->depth(), 1);
0174   BOOST_CHECK_EQUAL(node3->depth(), 2);
0175 }
0176 
0177 BOOST_AUTO_TEST_CASE(Static) {
0178   Blueprint::Config cfg;
0179   cfg.envelope[AxisDirection::AxisZ] = {20_mm, 20_mm};
0180   cfg.envelope[AxisDirection::AxisR] = {1_mm, 2_mm};
0181   Blueprint root{cfg};
0182 
0183   double hlZ = 30_mm;
0184   auto cylBounds = std::make_shared<CylinderVolumeBounds>(10_mm, 20_mm, hlZ);
0185   auto cyl = std::make_unique<TrackingVolume>(Transform3::Identity(), cylBounds,
0186                                               "child");
0187 
0188   root.addStaticVolume(std::move(cyl));
0189 
0190   BOOST_CHECK_EQUAL(root.children().size(), 1);
0191 
0192   auto tGeometry = root.construct({}, gctx, *logger);
0193 
0194   BOOST_REQUIRE(tGeometry);
0195 
0196   BOOST_CHECK(tGeometry->geometryVersion() ==
0197               TrackingGeometry::GeometryVersion::Gen3);
0198 
0199   BOOST_CHECK_EQUAL(tGeometry->highestTrackingVolume()->volumes().size(), 1);
0200 
0201   BOOST_CHECK_EQUAL(countVolumes(*tGeometry), 2);
0202 
0203   auto lookup = nameLookup(*tGeometry);
0204   BOOST_CHECK_EQUAL(lookup("child").volumeBounds().type(),
0205                     VolumeBounds::BoundsType::eCylinder);
0206   auto actCyl =
0207       dynamic_cast<const CylinderVolumeBounds&>(lookup("child").volumeBounds());
0208   // Size as given
0209   BOOST_CHECK_EQUAL(actCyl.get(CylinderVolumeBounds::eMinR), 10_mm);
0210   BOOST_CHECK_EQUAL(actCyl.get(CylinderVolumeBounds::eMaxR), 20_mm);
0211   BOOST_CHECK_EQUAL(actCyl.get(CylinderVolumeBounds::eHalfLengthZ), hlZ);
0212 
0213   BOOST_CHECK_EQUAL(lookup("World").volumeBounds().type(),
0214                     VolumeBounds::BoundsType::eCylinder);
0215   auto worldCyl =
0216       dynamic_cast<const CylinderVolumeBounds&>(lookup("World").volumeBounds());
0217   BOOST_CHECK_EQUAL(worldCyl.get(CylinderVolumeBounds::eMinR), 9_mm);
0218   BOOST_CHECK_EQUAL(worldCyl.get(CylinderVolumeBounds::eMaxR), 22_mm);
0219   BOOST_CHECK_EQUAL(worldCyl.get(CylinderVolumeBounds::eHalfLengthZ),
0220                     hlZ + 20_mm);
0221 
0222   BOOST_CHECK_EQUAL(lookup("World").portals().size(), 8);
0223 }
0224 
0225 BOOST_AUTO_TEST_CASE(CylinderContainer) {
0226   Blueprint::Config cfg;
0227   cfg.envelope[AxisDirection::AxisZ] = {20_mm, 20_mm};
0228   cfg.envelope[AxisDirection::AxisR] = {2_mm, 20_mm};
0229   auto root = std::make_unique<Blueprint>(cfg);
0230 
0231   auto& cyl = root->addCylinderContainer("Container", AxisDirection::AxisZ);
0232   cyl.setAttachmentStrategy(VolumeAttachmentStrategy::Gap);
0233 
0234   double z0 = -200_mm;
0235   double hlZ = 30_mm;
0236   auto cylBounds = std::make_shared<CylinderVolumeBounds>(10_mm, 20_mm, hlZ);
0237   for (std::size_t i = 0; i < 3; i++) {
0238     auto childCyl = std::make_unique<TrackingVolume>(
0239         Transform3::Identity() *
0240             Translation3{Vector3{0, 0, z0 + i * 2 * hlZ * 1.2}},
0241         cylBounds, "child" + std::to_string(i));
0242     cyl.addStaticVolume(std::move(childCyl));
0243   }
0244 
0245   auto tGeometry = root->construct({}, gctx, *logger);
0246   BOOST_CHECK(tGeometry->geometryVersion() ==
0247               TrackingGeometry::GeometryVersion::Gen3);
0248 
0249   // 4 real volumes + 2 gaps
0250   BOOST_CHECK_EQUAL(countVolumes(*tGeometry), 6);
0251 
0252   auto lookup = nameLookup(*tGeometry);
0253   BOOST_CHECK_EQUAL(lookup("World").volumeBounds().type(),
0254                     VolumeBounds::BoundsType::eCylinder);
0255   auto worldCyl =
0256       dynamic_cast<const CylinderVolumeBounds&>(lookup("World").volumeBounds());
0257   BOOST_CHECK_EQUAL(worldCyl.get(CylinderVolumeBounds::eMinR), 8_mm);
0258   BOOST_CHECK_EQUAL(worldCyl.get(CylinderVolumeBounds::eMaxR), 40_mm);
0259   BOOST_CHECK_EQUAL(worldCyl.get(CylinderVolumeBounds::eHalfLengthZ), 122_mm);
0260 
0261   BOOST_CHECK_EQUAL(lookup("World").portals().size(), 8);
0262 
0263   for (std::size_t i = 0; i < 3; i++) {
0264     const auto& vol{lookup("child" + std::to_string(i))};
0265     BOOST_CHECK_EQUAL(vol.volumeBounds().type(),
0266                       VolumeBounds::BoundsType::eCylinder);
0267     auto actCyl = dynamic_cast<const CylinderVolumeBounds&>(vol.volumeBounds());
0268     BOOST_CHECK_EQUAL(actCyl.get(CylinderVolumeBounds::eMinR), 10_mm);
0269     BOOST_CHECK_EQUAL(actCyl.get(CylinderVolumeBounds::eMaxR), 20_mm);
0270     BOOST_CHECK_EQUAL(actCyl.get(CylinderVolumeBounds::eHalfLengthZ), hlZ);
0271   }
0272 
0273   for (std::size_t i = 0; i < 2; i++) {
0274     const auto& vol{lookup("Container::Gap" + std::to_string(i + 1))};
0275     BOOST_CHECK_EQUAL(vol.volumeBounds().type(),
0276                       VolumeBounds::BoundsType::eCylinder);
0277     auto gapCyl = dynamic_cast<const CylinderVolumeBounds&>(vol.volumeBounds());
0278     BOOST_CHECK_EQUAL(gapCyl.get(CylinderVolumeBounds::eMinR), 10_mm);
0279     BOOST_CHECK_EQUAL(gapCyl.get(CylinderVolumeBounds::eMaxR), 20_mm);
0280     BOOST_CHECK_EQUAL(gapCyl.get(CylinderVolumeBounds::eHalfLengthZ), 6_mm);
0281   }
0282 }
0283 
0284 BOOST_AUTO_TEST_CASE(Confined) {
0285   Transform3 base{Transform3::Identity()};
0286 
0287   Blueprint::Config cfg;
0288   cfg.envelope[AxisDirection::AxisZ] = {20_mm, 20_mm};
0289   cfg.envelope[AxisDirection::AxisR] = {2_mm, 20_mm};
0290   auto root = std::make_unique<Blueprint>(cfg);
0291 
0292   root->addStaticVolume(
0293       base, std::make_shared<CylinderVolumeBounds>(50_mm, 400_mm, 1000_mm),
0294       "PixelWrapper", [&](auto& wrap) {
0295         double rMin = 100_mm;
0296         double rMax = 350_mm;
0297         double hlZ = 100_mm;
0298 
0299         wrap.addStaticVolume(
0300             base * Translation3{Vector3{0, 0, -600_mm}},
0301             std::make_shared<CylinderVolumeBounds>(rMin, rMax, hlZ),
0302             "PixelNeg1");
0303 
0304         wrap.addStaticVolume(
0305             base * Translation3{Vector3{0, 0, -200_mm}},
0306             std::make_shared<CylinderVolumeBounds>(rMin, rMax, hlZ),
0307             "PixelNeg2");
0308 
0309         wrap.addStaticVolume(
0310             base * Translation3{Vector3{0, 0, 200_mm}},
0311             std::make_shared<CylinderVolumeBounds>(rMin, rMax, hlZ),
0312             "PixelPos1");
0313 
0314         wrap.addStaticVolume(
0315             base * Translation3{Vector3{0, 0, 600_mm}},
0316             std::make_shared<CylinderVolumeBounds>(rMin, rMax, hlZ),
0317             "PixelPos2");
0318       });
0319 
0320   auto trackingGeometry = root->construct({}, gctx, *logger);
0321 
0322   // overall dimensions are the wrapper volume + envelope
0323   auto lookup = nameLookup(*trackingGeometry);
0324 
0325   BOOST_CHECK_EQUAL(lookup("World").volumeBounds().type(),
0326                     VolumeBounds::BoundsType::eCylinder);
0327   auto worldCyl =
0328       dynamic_cast<const CylinderVolumeBounds&>(lookup("World").volumeBounds());
0329   BOOST_CHECK_EQUAL(worldCyl.get(CylinderVolumeBounds::eMinR), 48_mm);
0330   BOOST_CHECK_EQUAL(worldCyl.get(CylinderVolumeBounds::eMaxR), 420_mm);
0331   BOOST_CHECK_EQUAL(worldCyl.get(CylinderVolumeBounds::eHalfLengthZ), 1020_mm);
0332 
0333   // 4 outer portals and 4 inner
0334   BOOST_CHECK_EQUAL(lookup("World").portals().size(), 8);
0335   BOOST_CHECK_EQUAL(lookup("World").volumes().size(), 1);
0336 
0337   BOOST_CHECK_EQUAL(lookup("PixelWrapper").volumeBounds().type(),
0338                     VolumeBounds::BoundsType::eCylinder);
0339 
0340   auto wrapperCyl = dynamic_cast<const CylinderVolumeBounds&>(
0341       lookup("PixelWrapper").volumeBounds());
0342   BOOST_CHECK_EQUAL(wrapperCyl.get(CylinderVolumeBounds::eMinR), 50_mm);
0343   BOOST_CHECK_EQUAL(wrapperCyl.get(CylinderVolumeBounds::eMaxR), 400_mm);
0344   BOOST_CHECK_EQUAL(wrapperCyl.get(CylinderVolumeBounds::eHalfLengthZ),
0345                     1000_mm);
0346   BOOST_CHECK_EQUAL(lookup("PixelWrapper").portals().size(), 4 + 4 * 4);
0347   BOOST_CHECK_EQUAL(lookup("PixelWrapper").volumes().size(), 4);
0348 
0349   for (const auto& name :
0350        {"PixelNeg1", "PixelNeg2", "PixelPos1", "PixelPos2"}) {
0351     BOOST_CHECK_EQUAL(lookup(name).volumeBounds().type(),
0352                       VolumeBounds::BoundsType::eCylinder);
0353     auto actCyl =
0354         dynamic_cast<const CylinderVolumeBounds&>(lookup(name).volumeBounds());
0355     BOOST_CHECK_EQUAL(actCyl.get(CylinderVolumeBounds::eMinR), 100_mm);
0356     BOOST_CHECK_EQUAL(actCyl.get(CylinderVolumeBounds::eMaxR), 350_mm);
0357     BOOST_CHECK_EQUAL(actCyl.get(CylinderVolumeBounds::eHalfLengthZ), 100_mm);
0358     BOOST_CHECK_EQUAL(lookup(name).portals().size(), 4);
0359   }
0360 }
0361 
0362 BOOST_AUTO_TEST_CASE(ConfinedWithShared) {
0363   Transform3 base{Transform3::Identity()};
0364 
0365   constexpr double rMin = 100_mm;
0366   constexpr double rMax = 350_mm;
0367   constexpr double hlZ = 100_mm;
0368 
0369   auto sharedBounds = std::make_shared<CylinderVolumeBounds>(rMin, rMax, hlZ);
0370   Blueprint::Config cfg;
0371   cfg.envelope[AxisDirection::AxisZ] = {20_mm, 20_mm};
0372   cfg.envelope[AxisDirection::AxisR] = {2_mm, 20_mm};
0373   auto root = std::make_unique<Blueprint>(cfg);
0374 
0375   root->addCylinderContainer(
0376       "PixelWrapper", AxisDirection::AxisZ, [&](auto& wrap) {
0377         wrap.addStaticVolume(base * Translation3{Vector3{0, 0, -750_mm}},
0378                              sharedBounds, "PixelNeg1");
0379 
0380         wrap.addStaticVolume(base * Translation3{Vector3{0, 0, -200_mm}},
0381                              sharedBounds, "PixelNeg2");
0382 
0383         wrap.addStaticVolume(base * Translation3{Vector3{0, 0, 200_mm}},
0384                              sharedBounds, "PixelPos1");
0385 
0386         wrap.addStaticVolume(base * Translation3{Vector3{0, 0, 975_mm}},
0387                              sharedBounds, "PixelPos2");
0388       });
0389   auto trackingGeometry = root->construct({}, gctx, *logger);
0390   // overall dimensions are the wrapper volume + envelope
0391   auto lookup = nameLookup(*trackingGeometry);
0392   BOOST_CHECK_EQUAL(lookup("World").volumeBounds().type(),
0393                     VolumeBounds::BoundsType::eCylinder);
0394   auto worldCyl =
0395       dynamic_cast<const CylinderVolumeBounds&>(lookup("World").volumeBounds());
0396   BOOST_CHECK_EQUAL(worldCyl.get(CylinderVolumeBounds::eMinR), 98_mm);
0397   BOOST_CHECK_EQUAL(worldCyl.get(CylinderVolumeBounds::eMaxR), 370_mm);
0398   BOOST_CHECK_EQUAL(worldCyl.get(CylinderVolumeBounds::eHalfLengthZ), 982.5_mm);
0399   // 4 outer portals and 4 inner
0400   BOOST_CHECK_EQUAL(lookup("World").portals().size(), 8);
0401   BOOST_CHECK_EQUAL(lookup("World").volumes().size(), 4);
0402 
0403   constexpr std::array<double, 4> expHalfL{187.5_mm, 237.5_mm, 293.75_mm,
0404                                            243.75_mm};
0405   const std::array<std::string, 4> volNames{"PixelNeg1", "PixelNeg2",
0406                                             "PixelPos1", "PixelPos2"};
0407   for (std::size_t v = 0; v < 4; ++v) {
0408     const auto& testMe{lookup(volNames[v])};
0409     BOOST_CHECK_EQUAL(testMe.volumeBounds().type(),
0410                       VolumeBounds::BoundsType::eCylinder);
0411     BOOST_CHECK_EQUAL(testMe.volumeBoundsPtr() != sharedBounds, true);
0412 
0413     auto actCyl =
0414         dynamic_cast<const CylinderVolumeBounds&>(testMe.volumeBounds());
0415     BOOST_CHECK_EQUAL(actCyl.get(CylinderVolumeBounds::eMinR), 100_mm);
0416     BOOST_CHECK_EQUAL(actCyl.get(CylinderVolumeBounds::eMaxR), 350_mm);
0417     BOOST_CHECK_EQUAL(actCyl.get(CylinderVolumeBounds::eHalfLengthZ),
0418                       expHalfL[v]);
0419     BOOST_CHECK_EQUAL(testMe.portals().size(), 4);
0420     if (v + 1 == 4) {
0421       break;
0422     }
0423     const auto& nextVol = lookup(volNames[(v + 1)]);
0424     const Vector3 outside =
0425         testMe.transform().translation() +
0426         Vector3{150_mm, 0.,
0427                 actCyl.get(CylinderVolumeBounds::eHalfLengthZ) - 0.5_mm};
0428     BOOST_CHECK_EQUAL(nextVol.inside(outside), false);
0429     const Vector3 inside =
0430         testMe.transform().translation() +
0431         Vector3{150_mm, 0.,
0432                 actCyl.get(CylinderVolumeBounds::eHalfLengthZ) + 0.5_mm};
0433     BOOST_CHECK_EQUAL(nextVol.inside(inside), true);
0434   }
0435 }
0436 
0437 BOOST_AUTO_TEST_CASE(DiscLayer) {
0438   double yrot = 45_degree;
0439   Transform3 base = Transform3::Identity() * AngleAxis3{yrot, Vector3::UnitY()};
0440 
0441   std::vector<std::shared_ptr<Surface>> surfaces;
0442   std::vector<std::unique_ptr<DetectorElementBase>> elements;
0443   double r = 300_mm;
0444   std::size_t nSensors = 8;
0445   double thickness = 2.5_mm;
0446   auto recBounds = std::make_shared<RectangleBounds>(40_mm, 60_mm);
0447 
0448   double deltaPhi = 2 * std::numbers::pi / nSensors;
0449   for (std::size_t i = 0; i < nSensors; i++) {
0450     // Create a fan of sensors
0451 
0452     Transform3 trf = base * AngleAxis3{deltaPhi * i, Vector3::UnitZ()} *
0453                      Translation3(Vector3::UnitX() * r);
0454 
0455     if (i % 2 == 0) {
0456       trf = trf * Translation3{Vector3::UnitZ() * 5_mm};
0457     }
0458 
0459     auto& element = elements.emplace_back(
0460         std::make_unique<DetectorElementStub>(trf, recBounds, thickness));
0461 
0462     element->surface().assignDetectorElement(*element);
0463 
0464     surfaces.push_back(element->surface().getSharedPtr());
0465   }
0466 
0467   std::function<void(LayerBlueprintNode&)> withSurfaces =
0468       [&surfaces, &base](LayerBlueprintNode& layer) {
0469         layer.setSurfaces(surfaces)
0470             .setLayerType(LayerBlueprintNode::LayerType::Disc)
0471             .setEnvelope(ExtentEnvelope{{
0472                 .z = {0.1_mm, 0.1_mm},
0473                 .r = {1_mm, 1_mm},
0474             }})
0475             .setTransform(base);
0476       };
0477 
0478   std::function<void(LayerBlueprintNode&)> withProtoLayer =
0479       [&surfaces, &base](LayerBlueprintNode& layer) {
0480         MutableProtoLayer protoLayer{gctx, surfaces, base.inverse()};
0481         layer.setProtoLayer(protoLayer)
0482             .setLayerType(LayerBlueprintNode::LayerType::Disc)
0483             .setEnvelope(ExtentEnvelope{{
0484                 .z = {0.1_mm, 0.1_mm},
0485                 .r = {1_mm, 1_mm},
0486             }});
0487       };
0488 
0489   for (const auto& func : {withSurfaces, withProtoLayer}) {
0490     Blueprint root{{.envelope = ExtentEnvelope{{
0491                         .z = {2_mm, 2_mm},
0492                         .r = {3_mm, 5_mm},
0493                     }}}};
0494 
0495     root.addLayer("Layer0", [&](auto& layer) { func(layer); });
0496 
0497     auto trackingGeometry = root.construct({}, gctx, *logger);
0498 
0499     std::size_t nSurfaces = 0;
0500 
0501     trackingGeometry->visitSurfaces([&](const Surface* surface) {
0502       if (surface->associatedDetectorElement() != nullptr) {
0503         nSurfaces++;
0504       }
0505     });
0506 
0507     BOOST_CHECK_EQUAL(nSurfaces, surfaces.size());
0508     BOOST_CHECK_EQUAL(countVolumes(*trackingGeometry), 2);
0509     auto lookup = nameLookup(*trackingGeometry);
0510 
0511     BOOST_CHECK_EQUAL(lookup("Layer0").volumeBounds().type(),
0512                       VolumeBounds::BoundsType::eCylinder);
0513     auto layerCyl = dynamic_cast<const CylinderVolumeBounds&>(
0514         lookup("Layer0").volumeBounds());
0515     BOOST_CHECK_CLOSE(layerCyl.get(CylinderVolumeBounds::eMinR), 258.9999999_mm,
0516                       1e-6);
0517     BOOST_CHECK_CLOSE(layerCyl.get(CylinderVolumeBounds::eMaxR),
0518                       346.25353003_mm, 1e-6);
0519     BOOST_CHECK_CLOSE(layerCyl.get(CylinderVolumeBounds::eHalfLengthZ), 3.85_mm,
0520                       1e-6);
0521   }
0522 }
0523 
0524 BOOST_AUTO_TEST_CASE(CylinderLayer) {
0525   double yrot = 0_degree;
0526   Transform3 base = Transform3::Identity() * AngleAxis3{yrot, Vector3::UnitY()};
0527 
0528   std::vector<std::shared_ptr<Surface>> surfaces;
0529   std::vector<std::unique_ptr<DetectorElementBase>> elements;
0530 
0531   double r = 300_mm;
0532   std::size_t nStaves = 10;
0533   int nSensorsPerStave = 8;
0534   double thickness = 0;
0535   double hlPhi = 40_mm;
0536   double hlZ = 60_mm;
0537   auto recBounds = std::make_shared<RectangleBounds>(hlPhi, hlZ);
0538 
0539   double deltaPhi = 2 * std::numbers::pi / nStaves;
0540 
0541   for (std::size_t istave = 0; istave < nStaves; istave++) {
0542     for (int isensor = -nSensorsPerStave; isensor <= nSensorsPerStave;
0543          isensor++) {
0544       double z = isensor * (2 * hlZ + 5_mm);
0545 
0546       Transform3 trf = base * Translation3(Vector3::UnitZ() * z) *
0547                        AngleAxis3{deltaPhi * istave, Vector3::UnitZ()} *
0548                        Translation3(Vector3::UnitX() * r) *
0549                        AngleAxis3{10_degree, Vector3::UnitZ()} *
0550                        AngleAxis3{90_degree, Vector3::UnitY()} *
0551                        AngleAxis3{90_degree, Vector3::UnitZ()};
0552       auto& element = elements.emplace_back(
0553           std::make_unique<DetectorElementStub>(trf, recBounds, thickness));
0554       element->surface().assignDetectorElement(*element);
0555       surfaces.push_back(element->surface().getSharedPtr());
0556     }
0557   }
0558 
0559   std::function<void(LayerBlueprintNode&)> withSurfaces =
0560       [&surfaces, &base](LayerBlueprintNode& layer) {
0561         layer.setSurfaces(surfaces)
0562             .setLayerType(LayerBlueprintNode::LayerType::Cylinder)
0563             .setEnvelope(ExtentEnvelope{{
0564                 .z = {10_mm, 10_mm},
0565                 .r = {20_mm, 10_mm},
0566             }})
0567             .setTransform(base);
0568       };
0569 
0570   std::function<void(LayerBlueprintNode&)> withProtoLayer =
0571       [&surfaces, &base](LayerBlueprintNode& layer) {
0572         MutableProtoLayer protoLayer{gctx, surfaces, base.inverse()};
0573         layer.setProtoLayer(protoLayer)
0574             .setLayerType(LayerBlueprintNode::LayerType::Cylinder)
0575             .setEnvelope(ExtentEnvelope{{
0576                 .z = {10_mm, 10_mm},
0577                 .r = {20_mm, 10_mm},
0578             }});
0579       };
0580 
0581   for (const auto& func : {withSurfaces, withProtoLayer}) {
0582     Blueprint root{{.envelope = ExtentEnvelope{{
0583                         .z = {2_mm, 2_mm},
0584                         .r = {3_mm, 5_mm},
0585                     }}}};
0586 
0587     root.addLayer("Layer0", [&](auto& layer) { func(layer); });
0588 
0589     auto trackingGeometry = root.construct({}, gctx, *logger);
0590 
0591     std::size_t nSurfaces = 0;
0592 
0593     trackingGeometry->visitSurfaces([&](const Surface* surface) {
0594       if (surface->associatedDetectorElement() != nullptr) {
0595         nSurfaces++;
0596       }
0597     });
0598 
0599     BOOST_CHECK_EQUAL(nSurfaces, surfaces.size());
0600     BOOST_CHECK_EQUAL(countVolumes(*trackingGeometry), 2);
0601     auto lookup = nameLookup(*trackingGeometry);
0602 
0603     BOOST_CHECK_EQUAL(lookup("Layer0").volumeBounds().type(),
0604                       VolumeBounds::BoundsType::eCylinder);
0605     auto layerCyl = dynamic_cast<const CylinderVolumeBounds&>(
0606         lookup("Layer0").volumeBounds());
0607     BOOST_CHECK_EQUAL(lookup("Layer0").portals().size(), 4);
0608     BOOST_CHECK_CLOSE(layerCyl.get(CylinderVolumeBounds::eMinR), 275.6897761_mm,
0609                       1e-6);
0610     BOOST_CHECK_CLOSE(layerCyl.get(CylinderVolumeBounds::eMaxR), 319.4633358_mm,
0611                       1e-6);
0612     BOOST_CHECK_CLOSE(layerCyl.get(CylinderVolumeBounds::eHalfLengthZ), 1070_mm,
0613                       1e-6);
0614   }
0615 }
0616 
0617 BOOST_AUTO_TEST_CASE(MaterialTesting) {
0618   Blueprint::Config cfg;
0619   cfg.envelope[AxisDirection::AxisZ] = {20_mm, 20_mm};
0620   cfg.envelope[AxisDirection::AxisR] = {1_mm, 2_mm};
0621   Blueprint root{cfg};
0622 
0623   double hlZ = 30_mm;
0624   auto cylBounds = std::make_shared<CylinderVolumeBounds>(10_mm, 20_mm, hlZ);
0625   auto cyl = std::make_unique<TrackingVolume>(Transform3::Identity(), cylBounds,
0626                                               "child");
0627 
0628   using enum AxisDirection;
0629   using enum CylinderVolumeBounds::Face;
0630   using enum AxisBoundaryType;
0631 
0632   root.addMaterial("Material", [&](auto& mat) {
0633     mat.configureFace(NegativeDisc, {AxisR, Bound, 5}, {AxisPhi, Bound, 10});
0634     mat.configureFace(PositiveDisc, {AxisR, Bound, 15}, {AxisPhi, Bound, 20});
0635     mat.configureFace(OuterCylinder, {AxisRPhi, Bound, 25}, {AxisZ, Bound, 30});
0636 
0637     mat.addStaticVolume(std::move(cyl));
0638   });
0639 
0640   auto trackingGeometry = root.construct({}, gctx, *logger);
0641 
0642   BOOST_CHECK_EQUAL(countVolumes(*trackingGeometry), 2);
0643   auto lookup = nameLookup(*trackingGeometry);
0644   auto& child = lookup("child");
0645 
0646   // Check negative disc material
0647   const auto* negDisc = child.portals()
0648                             .at(static_cast<std::size_t>(NegativeDisc))
0649                             .surface()
0650                             .surfaceMaterial();
0651   BOOST_REQUIRE_NE(negDisc, nullptr);
0652   const auto& negDiscMat =
0653       dynamic_cast<const ProtoGridSurfaceMaterial&>(*negDisc);
0654   // Check positive disc material
0655   const auto* posDisc = child.portals()
0656                             .at(static_cast<std::size_t>(PositiveDisc))
0657                             .surface()
0658                             .surfaceMaterial();
0659   BOOST_REQUIRE_NE(posDisc, nullptr);
0660   const auto& posDiscMat =
0661       dynamic_cast<const ProtoGridSurfaceMaterial&>(*posDisc);
0662 
0663   BOOST_CHECK_EQUAL(negDiscMat.binning().at(0).getAxis().getNBins(), 5);
0664   BOOST_CHECK_EQUAL(negDiscMat.binning().at(1).getAxis().getNBins(), 10);
0665   BOOST_CHECK_EQUAL(posDiscMat.binning().at(0).getAxis().getNBins(), 15);
0666   BOOST_CHECK_EQUAL(posDiscMat.binning().at(1).getAxis().getNBins(), 20);
0667 
0668   // Check outer cylinder material
0669   const auto* outerCyl = child.portals()
0670                              .at(static_cast<std::size_t>(OuterCylinder))
0671                              .surface()
0672                              .surfaceMaterial();
0673   BOOST_REQUIRE_NE(outerCyl, nullptr);
0674   const auto& outerCylMat =
0675       dynamic_cast<const ProtoGridSurfaceMaterial&>(*outerCyl);
0676   BOOST_CHECK_EQUAL(outerCylMat.binning().at(0).getAxis().getNBins(), 25);
0677   BOOST_CHECK_EQUAL(outerCylMat.binning().at(1).getAxis().getNBins(), 30);
0678 
0679   // Check that other faces have no material
0680   for (std::size_t i = 0; i < child.portals().size(); i++) {
0681     if (i != static_cast<std::size_t>(NegativeDisc) &&
0682         i != static_cast<std::size_t>(PositiveDisc) &&
0683         i != static_cast<std::size_t>(OuterCylinder)) {
0684       BOOST_CHECK_EQUAL(child.portals().at(i).surface().surfaceMaterial(),
0685                         nullptr);
0686     }
0687   }
0688 }
0689 
0690 BOOST_AUTO_TEST_CASE(MaterialInvalidAxisDirections) {
0691   Blueprint::Config cfg;
0692   cfg.envelope[AxisDirection::AxisZ] = {20_mm, 20_mm};
0693   cfg.envelope[AxisDirection::AxisR] = {1_mm, 2_mm};
0694   Blueprint root{cfg};
0695 
0696   using enum AxisDirection;
0697   using enum AxisBoundaryType;
0698 
0699   // Test invalid axis direction combinations for cylinder faces
0700   BOOST_CHECK_THROW(
0701       root.addMaterial("Material",
0702                        [&](auto& mat) {
0703                          mat.configureFace(
0704                              CylinderVolumeBounds::Face::NegativeDisc,
0705                              {AxisZ, Bound, 5}, {AxisPhi, Bound, 10});
0706                        }),
0707       std::invalid_argument);
0708 
0709   BOOST_CHECK_THROW(
0710       root.addMaterial("Material",
0711                        [&](auto& mat) {
0712                          mat.configureFace(
0713                              CylinderVolumeBounds::Face::OuterCylinder,
0714                              {AxisR, Bound, 5}, {AxisR, Bound, 10});
0715                        }),
0716       std::invalid_argument);
0717 
0718   // Test invalid axis direction combinations for cuboid faces
0719   BOOST_CHECK_THROW(
0720       root.addMaterial("Material",
0721                        [&](auto& mat) {
0722                          mat.configureFace(
0723                              CuboidVolumeBounds::Face::NegativeXFace,
0724                              {AxisX, Bound, 5}, {AxisZ, Bound, 10});
0725                        }),
0726       std::invalid_argument);
0727 
0728   BOOST_CHECK_THROW(
0729       root.addMaterial("Material",
0730                        [&](auto& mat) {
0731                          mat.configureFace(
0732                              CuboidVolumeBounds::Face::PositiveYFace,
0733                              {AxisY, Bound, 5}, {AxisX, Bound, 10});
0734                        }),
0735       std::invalid_argument);
0736 
0737   BOOST_CHECK_THROW(
0738       root.addMaterial("Material",
0739                        [&](auto& mat) {
0740                          mat.configureFace(
0741                              CuboidVolumeBounds::Face::NegativeZFace,
0742                              {AxisZ, Bound, 5}, {AxisY, Bound, 10});
0743                        }),
0744       std::invalid_argument);
0745 }
0746 
0747 BOOST_AUTO_TEST_CASE(MaterialMixedVolumeTypes) {
0748   Blueprint::Config cfg;
0749   cfg.envelope[AxisDirection::AxisZ] = {20_mm, 20_mm};
0750   cfg.envelope[AxisDirection::AxisR] = {1_mm, 2_mm};
0751   Blueprint root{cfg};
0752 
0753   using enum AxisDirection;
0754   using enum AxisBoundaryType;
0755 
0756   // Configure for cylinder first, then try to add cuboid - should throw
0757   BOOST_CHECK_THROW(
0758       root.addMaterial(
0759           "Material",
0760           [&](auto& mat) {
0761             mat.configureFace(CylinderVolumeBounds::Face::NegativeDisc,
0762                               {AxisR, Bound, 5}, {AxisPhi, Bound, 10});
0763             mat.configureFace(CuboidVolumeBounds::Face::NegativeXFace,
0764                               {AxisX, Bound, 5}, {AxisY, Bound, 10});
0765           }),
0766       std::runtime_error);
0767 
0768   // Configure for cuboid first, then try to add cylinder - should throw
0769   BOOST_CHECK_THROW(
0770       root.addMaterial(
0771           "Material",
0772           [&](auto& mat) {
0773             mat.configureFace(CuboidVolumeBounds::Face::NegativeXFace,
0774                               {AxisX, Bound, 5}, {AxisY, Bound, 10});
0775             mat.configureFace(CylinderVolumeBounds::Face::NegativeDisc,
0776                               {AxisR, Bound, 5}, {AxisPhi, Bound, 10});
0777           }),
0778       std::runtime_error);
0779 }
0780 
0781 BOOST_AUTO_TEST_CASE(MaterialCuboid) {
0782   Blueprint::Config cfg;
0783   cfg.envelope[AxisDirection::AxisZ] = {20_mm, 20_mm};
0784   cfg.envelope[AxisDirection::AxisR] = {1_mm, 2_mm};
0785   Blueprint root{cfg};
0786 
0787   using enum AxisDirection;
0788   using enum AxisBoundaryType;
0789   using enum CuboidVolumeBounds::Face;
0790 
0791   double hlX = 30_mm;
0792   double hlY = 40_mm;
0793   double hlZ = 50_mm;
0794   auto cuboidBounds = std::make_shared<CuboidVolumeBounds>(hlX, hlY, hlZ);
0795   auto cuboid = std::make_unique<TrackingVolume>(Transform3::Identity(),
0796                                                  cuboidBounds, "child");
0797 
0798   auto mat = std::make_shared<MaterialDesignatorBlueprintNode>("Material");
0799 
0800   // Configure material for different faces with different binning
0801   mat->configureFace(NegativeXFace, {AxisX, Bound, 5}, {AxisY, Bound, 10});
0802   mat->configureFace(PositiveXFace, {AxisX, Bound, 15}, {AxisY, Bound, 20});
0803   mat->configureFace(NegativeYFace, {AxisX, Bound, 25}, {AxisY, Bound, 30});
0804   mat->configureFace(PositiveYFace, {AxisX, Bound, 35}, {AxisY, Bound, 40});
0805   mat->configureFace(NegativeZFace, {AxisX, Bound, 45}, {AxisY, Bound, 50});
0806   mat->configureFace(PositiveZFace, {AxisX, Bound, 55}, {AxisY, Bound, 60});
0807 
0808   mat->addChild(std::make_shared<StaticBlueprintNode>(std::move(cuboid)));
0809 
0810   root.addChild(mat);
0811 
0812   auto trackingGeometry =
0813       root.construct({}, gctx, *logger->clone(std::nullopt, Logging::VERBOSE));
0814 
0815   BOOST_REQUIRE(trackingGeometry);
0816 
0817   auto lookup = nameLookup(*trackingGeometry);
0818   auto& child = lookup("child");
0819 
0820   // Check that material is attached to all faces
0821   for (std::size_t i = 0; i < child.portals().size(); i++) {
0822     const auto* material = child.portals().at(i).surface().surfaceMaterial();
0823     BOOST_REQUIRE_NE(material, nullptr);
0824 
0825     const auto& gridMaterial =
0826         dynamic_cast<const ProtoGridSurfaceMaterial&>(*material);
0827 
0828     // Check binning based on face
0829     CuboidVolumeBounds::Face face = static_cast<CuboidVolumeBounds::Face>(i);
0830     switch (face) {
0831       case NegativeXFace:
0832         BOOST_CHECK_EQUAL(gridMaterial.binning().at(0).getAxis().getNBins(), 5);
0833         BOOST_CHECK_EQUAL(gridMaterial.binning().at(1).getAxis().getNBins(),
0834                           10);
0835         break;
0836       case PositiveXFace:
0837         BOOST_CHECK_EQUAL(gridMaterial.binning().at(0).getAxis().getNBins(),
0838                           15);
0839         BOOST_CHECK_EQUAL(gridMaterial.binning().at(1).getAxis().getNBins(),
0840                           20);
0841         break;
0842       case NegativeYFace:
0843         BOOST_CHECK_EQUAL(gridMaterial.binning().at(0).getAxis().getNBins(),
0844                           25);
0845         BOOST_CHECK_EQUAL(gridMaterial.binning().at(1).getAxis().getNBins(),
0846                           30);
0847         break;
0848       case PositiveYFace:
0849         BOOST_CHECK_EQUAL(gridMaterial.binning().at(0).getAxis().getNBins(),
0850                           35);
0851         BOOST_CHECK_EQUAL(gridMaterial.binning().at(1).getAxis().getNBins(),
0852                           40);
0853         break;
0854       case NegativeZFace:
0855         BOOST_CHECK_EQUAL(gridMaterial.binning().at(0).getAxis().getNBins(),
0856                           45);
0857         BOOST_CHECK_EQUAL(gridMaterial.binning().at(1).getAxis().getNBins(),
0858                           50);
0859         break;
0860       case PositiveZFace:
0861         BOOST_CHECK_EQUAL(gridMaterial.binning().at(0).getAxis().getNBins(),
0862                           55);
0863         BOOST_CHECK_EQUAL(gridMaterial.binning().at(1).getAxis().getNBins(),
0864                           60);
0865         break;
0866     }
0867   }
0868 }
0869 
0870 BOOST_AUTO_TEST_CASE(HomogeneousMaterialCylinder) {
0871   Blueprint::Config cfg;
0872   cfg.envelope[AxisDirection::AxisZ] = {20_mm, 20_mm};
0873   cfg.envelope[AxisDirection::AxisR] = {1_mm, 2_mm};
0874   Blueprint root{cfg};
0875 
0876   double hlZ = 30_mm;
0877   auto cylBounds = std::make_shared<CylinderVolumeBounds>(10_mm, 20_mm, hlZ);
0878   auto cyl = std::make_unique<TrackingVolume>(Transform3::Identity(), cylBounds,
0879                                               "child");
0880 
0881   using enum CylinderVolumeBounds::Face;
0882 
0883   // Create some homogeneous materials with different properties
0884   auto testMaterial = Material::fromMolarDensity(
0885       9.370_cm, 46.52_cm, 28.0855, 14, (2.329 / 28.0855) * 1_mol / 1_cm3);
0886 
0887   auto negDiscMat = std::make_shared<HomogeneousSurfaceMaterial>(
0888       MaterialSlab(testMaterial, 0.1_mm));
0889   auto posDiscMat = std::make_shared<HomogeneousSurfaceMaterial>(
0890       MaterialSlab(testMaterial, 0.2_mm));
0891   auto outerCylMat = std::make_shared<HomogeneousSurfaceMaterial>(
0892       MaterialSlab(testMaterial, 0.3_mm));
0893 
0894   root.addMaterial("Material", [&](auto& mat) {
0895     mat.configureFace(NegativeDisc, negDiscMat);
0896     mat.configureFace(PositiveDisc, posDiscMat);
0897     mat.configureFace(OuterCylinder, outerCylMat);
0898 
0899     mat.addStaticVolume(std::move(cyl));
0900   });
0901 
0902   auto trackingGeometry = root.construct({}, gctx, *logger);
0903 
0904   BOOST_CHECK_EQUAL(countVolumes(*trackingGeometry), 2);
0905   auto lookup = nameLookup(*trackingGeometry);
0906   auto& child = lookup("child");
0907 
0908   // Check negative disc material
0909   const auto* negDisc = child.portals()
0910                             .at(static_cast<std::size_t>(NegativeDisc))
0911                             .surface()
0912                             .surfaceMaterial();
0913   BOOST_REQUIRE_NE(negDisc, nullptr);
0914   BOOST_CHECK_EQUAL(negDisc, negDiscMat.get());
0915 
0916   // Check positive disc material
0917   const auto* posDisc = child.portals()
0918                             .at(static_cast<std::size_t>(PositiveDisc))
0919                             .surface()
0920                             .surfaceMaterial();
0921   BOOST_REQUIRE_NE(posDisc, nullptr);
0922   BOOST_CHECK_EQUAL(posDisc, posDiscMat.get());
0923 
0924   // Check outer cylinder material
0925   const auto* outerCyl = child.portals()
0926                              .at(static_cast<std::size_t>(OuterCylinder))
0927                              .surface()
0928                              .surfaceMaterial();
0929   BOOST_REQUIRE_NE(outerCyl, nullptr);
0930   BOOST_CHECK_EQUAL(outerCyl, outerCylMat.get());
0931 
0932   // Check that other faces have no material
0933   for (std::size_t i = 0; i < child.portals().size(); i++) {
0934     if (i != static_cast<std::size_t>(NegativeDisc) &&
0935         i != static_cast<std::size_t>(PositiveDisc) &&
0936         i != static_cast<std::size_t>(OuterCylinder)) {
0937       BOOST_CHECK_EQUAL(child.portals().at(i).surface().surfaceMaterial(),
0938                         nullptr);
0939     }
0940   }
0941 }
0942 
0943 BOOST_AUTO_TEST_CASE(HomogeneousMaterialCuboid) {
0944   Blueprint::Config cfg;
0945   cfg.envelope[AxisDirection::AxisZ] = {20_mm, 20_mm};
0946   cfg.envelope[AxisDirection::AxisR] = {1_mm, 2_mm};
0947   Blueprint root{cfg};
0948 
0949   using enum CuboidVolumeBounds::Face;
0950 
0951   double hlX = 30_mm;
0952   double hlY = 40_mm;
0953   double hlZ = 50_mm;
0954   auto cuboidBounds = std::make_shared<CuboidVolumeBounds>(hlX, hlY, hlZ);
0955   auto cuboid = std::make_unique<TrackingVolume>(Transform3::Identity(),
0956                                                  cuboidBounds, "child");
0957 
0958   // Create different homogeneous materials for each face
0959   auto testMaterial = Material::fromMolarDensity(
0960       9.370_cm, 46.52_cm, 28.0855, 14, (2.329 / 28.0855) * 1_mol / 1_cm3);
0961 
0962   auto negXMat = std::make_shared<HomogeneousSurfaceMaterial>(
0963       MaterialSlab(testMaterial, 0.1_mm));
0964   auto posXMat = std::make_shared<HomogeneousSurfaceMaterial>(
0965       MaterialSlab(testMaterial, 0.2_mm));
0966   auto negYMat = std::make_shared<HomogeneousSurfaceMaterial>(
0967       MaterialSlab(testMaterial, 0.3_mm));
0968   auto posYMat = std::make_shared<HomogeneousSurfaceMaterial>(
0969       MaterialSlab(testMaterial, 0.4_mm));
0970   auto negZMat = std::make_shared<HomogeneousSurfaceMaterial>(
0971       MaterialSlab(testMaterial, 0.5_mm));
0972   auto posZMat = std::make_shared<HomogeneousSurfaceMaterial>(
0973       MaterialSlab(testMaterial, 0.6_mm));
0974 
0975   root.addMaterial("Material", [&](auto& mat) {
0976     mat.configureFace(NegativeXFace, negXMat);
0977     mat.configureFace(PositiveXFace, posXMat);
0978     mat.configureFace(NegativeYFace, negYMat);
0979     mat.configureFace(PositiveYFace, posYMat);
0980     mat.configureFace(NegativeZFace, negZMat);
0981     mat.configureFace(PositiveZFace, posZMat);
0982 
0983     mat.addStaticVolume(std::move(cuboid));
0984   });
0985 
0986   auto trackingGeometry = root.construct({}, gctx, *logger);
0987 
0988   BOOST_REQUIRE(trackingGeometry);
0989 
0990   auto lookup = nameLookup(*trackingGeometry);
0991   auto& child = lookup("child");
0992 
0993   // Check that material is attached to all faces with correct properties
0994   for (std::size_t i = 0; i < child.portals().size(); i++) {
0995     const auto* material = child.portals().at(i).surface().surfaceMaterial();
0996     BOOST_REQUIRE_NE(material, nullptr);
0997 
0998     const auto* homMaterial =
0999         dynamic_cast<const HomogeneousSurfaceMaterial*>(material);
1000     BOOST_REQUIRE_NE(homMaterial, nullptr);
1001 
1002     // Check thickness based on face
1003     CuboidVolumeBounds::Face face = static_cast<CuboidVolumeBounds::Face>(i);
1004     switch (face) {
1005       case NegativeXFace:
1006         BOOST_CHECK_EQUAL(homMaterial, negXMat.get());
1007         break;
1008       case PositiveXFace:
1009         BOOST_CHECK_EQUAL(homMaterial, posXMat.get());
1010         break;
1011       case NegativeYFace:
1012         BOOST_CHECK_EQUAL(homMaterial, negYMat.get());
1013         break;
1014       case PositiveYFace:
1015         BOOST_CHECK_EQUAL(homMaterial, posYMat.get());
1016         break;
1017       case NegativeZFace:
1018         BOOST_CHECK_EQUAL(homMaterial, negZMat.get());
1019         break;
1020       case PositiveZFace:
1021         BOOST_CHECK_EQUAL(homMaterial, posZMat.get());
1022         break;
1023     }
1024   }
1025 }
1026 
1027 BOOST_AUTO_TEST_CASE(HomogeneousMaterialMixedVolumeTypes) {
1028   Blueprint::Config cfg;
1029   cfg.envelope[AxisDirection::AxisZ] = {20_mm, 20_mm};
1030   cfg.envelope[AxisDirection::AxisR] = {1_mm, 2_mm};
1031   Blueprint root{cfg};
1032 
1033   auto testMaterial = Material::fromMolarDensity(
1034       9.370_cm, 46.52_cm, 28.0855, 14, (2.329 / 28.0855) * 1_mol / 1_cm3);
1035 
1036   auto material = std::make_shared<HomogeneousSurfaceMaterial>(
1037       MaterialSlab(testMaterial, 0.1_mm));
1038 
1039   // Configure for cylinder first, then try to add cuboid - should throw
1040   BOOST_CHECK_THROW(
1041       root.addMaterial("Material",
1042                        [&](auto& mat) {
1043                          mat.configureFace(
1044                              CylinderVolumeBounds::Face::NegativeDisc,
1045                              material);
1046                          mat.configureFace(
1047                              CuboidVolumeBounds::Face::NegativeXFace, material);
1048                        }),
1049       std::runtime_error);
1050 
1051   // Configure for cuboid first, then try to add cylinder - should throw
1052   BOOST_CHECK_THROW(
1053       root.addMaterial("Material",
1054                        [&](auto& mat) {
1055                          mat.configureFace(
1056                              CuboidVolumeBounds::Face::NegativeXFace, material);
1057                          mat.configureFace(
1058                              CylinderVolumeBounds::Face::NegativeDisc,
1059                              material);
1060                        }),
1061       std::runtime_error);
1062 }
1063 
1064 BOOST_AUTO_TEST_CASE(LayerCenterOfGravity) {
1065   // Test disc layer with center of gravity disabled
1066   {
1067     double yrot = 45_degree;
1068     Transform3 base =
1069         Transform3::Identity() * AngleAxis3{yrot, Vector3::UnitY()};
1070 
1071     std::vector<std::shared_ptr<Surface>> surfaces;
1072     std::vector<std::unique_ptr<DetectorElementBase>> elements;
1073     double r = 300_mm;
1074     std::size_t nSensors = 8;
1075     double thickness = 2.5_mm;
1076     auto recBounds = std::make_shared<RectangleBounds>(40_mm, 60_mm);
1077 
1078     double deltaPhi = 2 * std::numbers::pi / nSensors;
1079     for (std::size_t i = 0; i < nSensors; i++) {
1080       Transform3 trf = base * AngleAxis3{deltaPhi * i, Vector3::UnitZ()} *
1081                        Translation3(Vector3::UnitX() * r);
1082 
1083       if (i % 2 == 0) {
1084         trf = trf * Translation3{Vector3::UnitZ() * 5_mm};
1085       }
1086 
1087       auto& element = elements.emplace_back(
1088           std::make_unique<DetectorElementStub>(trf, recBounds, thickness));
1089 
1090       element->surface().assignDetectorElement(*element);
1091       surfaces.push_back(element->surface().getSharedPtr());
1092     }
1093 
1094     Blueprint root{{.envelope = ExtentEnvelope{{
1095                         .z = {2_mm, 2_mm},
1096                         .r = {3_mm, 5_mm},
1097                     }}}};
1098 
1099     root.addLayer("Layer0", [&](auto& layer) {
1100       layer.setSurfaces(surfaces)
1101           .setLayerType(LayerBlueprintNode::LayerType::Disc)
1102           .setEnvelope(ExtentEnvelope{{
1103               .z = {0.1_mm, 0.1_mm},
1104               .r = {1_mm, 1_mm},
1105           }})
1106           .setTransform(base)
1107           .setUseCenterOfGravity(false, false, false);  // Disable all axes
1108     });
1109 
1110     auto trackingGeometry = root.construct({}, gctx, *logger);
1111     auto lookup = nameLookup(*trackingGeometry);
1112 
1113     BOOST_CHECK_EQUAL(lookup("Layer0").volumeBounds().type(),
1114                       VolumeBounds::BoundsType::eCylinder);
1115 
1116     auto layerCyl = dynamic_cast<const CylinderVolumeBounds&>(
1117         lookup("Layer0").volumeBounds());
1118 
1119     // With center of gravity disabled, the layer should be at the origin
1120     BOOST_CHECK_CLOSE(layerCyl.get(CylinderVolumeBounds::eMinR), 258.9999999_mm,
1121                       1e-6);
1122     BOOST_CHECK_CLOSE(layerCyl.get(CylinderVolumeBounds::eMaxR),
1123                       346.25353003_mm, 1e-6);
1124     BOOST_CHECK_CLOSE(layerCyl.get(CylinderVolumeBounds::eHalfLengthZ), 3.85_mm,
1125                       1e-6);
1126   }
1127 
1128   // Test cylinder layer with center of gravity disabled
1129   {
1130     double yrot = 0_degree;
1131     Transform3 base =
1132         Transform3::Identity() * AngleAxis3{yrot, Vector3::UnitY()};
1133 
1134     std::vector<std::shared_ptr<Surface>> surfaces;
1135     std::vector<std::unique_ptr<DetectorElementBase>> elements;
1136 
1137     double r = 300_mm;
1138     std::size_t nStaves = 10;
1139     int nSensorsPerStave = 8;
1140     double thickness = 0;
1141     double hlPhi = 40_mm;
1142     double hlZ = 60_mm;
1143     auto recBounds = std::make_shared<RectangleBounds>(hlPhi, hlZ);
1144 
1145     double deltaPhi = 2 * std::numbers::pi / nStaves;
1146 
1147     for (std::size_t istave = 0; istave < nStaves; istave++) {
1148       for (int isensor = -nSensorsPerStave; isensor <= nSensorsPerStave;
1149            isensor++) {
1150         double z = isensor * (2 * hlZ + 5_mm);
1151 
1152         Transform3 trf = base * Translation3(Vector3::UnitZ() * z) *
1153                          AngleAxis3{deltaPhi * istave, Vector3::UnitZ()} *
1154                          Translation3(Vector3::UnitX() * r) *
1155                          AngleAxis3{10_degree, Vector3::UnitZ()} *
1156                          AngleAxis3{90_degree, Vector3::UnitY()} *
1157                          AngleAxis3{90_degree, Vector3::UnitZ()};
1158         auto& element = elements.emplace_back(
1159             std::make_unique<DetectorElementStub>(trf, recBounds, thickness));
1160         element->surface().assignDetectorElement(*element);
1161         surfaces.push_back(element->surface().getSharedPtr());
1162       }
1163     }
1164 
1165     Blueprint root{{.envelope = ExtentEnvelope{{
1166                         .z = {2_mm, 2_mm},
1167                         .r = {3_mm, 5_mm},
1168                     }}}};
1169 
1170     root.addLayer("Layer0", [&](auto& layer) {
1171       layer.setSurfaces(surfaces)
1172           .setLayerType(LayerBlueprintNode::LayerType::Cylinder)
1173           .setEnvelope(ExtentEnvelope{{
1174               .z = {10_mm, 10_mm},
1175               .r = {20_mm, 10_mm},
1176           }})
1177           .setTransform(base)
1178           .setUseCenterOfGravity(false, false, false);  // Disable all axes
1179     });
1180 
1181     auto trackingGeometry = root.construct({}, gctx, *logger);
1182     auto lookup = nameLookup(*trackingGeometry);
1183     BOOST_CHECK_EQUAL(lookup("Layer0").volumeBounds().type(),
1184                       VolumeBounds::BoundsType::eCylinder);
1185 
1186     auto layerCyl = dynamic_cast<const CylinderVolumeBounds&>(
1187         lookup("Layer0").volumeBounds());
1188 
1189     // With center of gravity disabled, the layer should be at the origin
1190     BOOST_CHECK_EQUAL(lookup("Layer0").portals().size(), 4);
1191     BOOST_CHECK_CLOSE(layerCyl.get(CylinderVolumeBounds::eMinR), 275.6897761_mm,
1192                       1e-6);
1193     BOOST_CHECK_CLOSE(layerCyl.get(CylinderVolumeBounds::eMaxR), 319.4633358_mm,
1194                       1e-6);
1195     BOOST_CHECK_CLOSE(layerCyl.get(CylinderVolumeBounds::eHalfLengthZ), 1070_mm,
1196                       1e-6);
1197   }
1198 }
1199 
1200 BOOST_AUTO_TEST_CASE(GeometryIdnetifiersForPortals) {
1201   Blueprint::Config cfg;
1202   cfg.envelope[AxisDirection::AxisZ] = {20_mm, 20_mm};
1203   cfg.envelope[AxisDirection::AxisR] = {1_mm, 2_mm};
1204   Blueprint root{cfg};
1205 
1206   auto& cubcontainer =
1207       root.addCuboidContainer("CuboidContainer", AxisDirection::AxisX);
1208   auto parentBounds = std::make_shared<CuboidVolumeBounds>(1_m, 20_mm, 20_mm);
1209   auto parentVol = std::make_unique<TrackingVolume>(Transform3::Identity(),
1210                                                     parentBounds, "parent");
1211   parentVol->assignGeometryId(GeometryIdentifier{}.withVolume(1));
1212   auto parentNode = std::make_shared<StaticBlueprintNode>(std::move(parentVol));
1213   std::size_t nChambers = 50;
1214   // start from the edge of the parent volume
1215   double startX = -1000. + 3. + 0.5;
1216   Transform3 trf = Transform3(Translation3(startX, 0, 0));
1217   auto tbounds =
1218       std::make_shared<TrapezoidVolumeBounds>(3_mm, 3_mm, 10_mm, 15_mm);
1219 
1220   for (std::size_t i = 0; i < nChambers; i++) {
1221     // move the chambers position
1222     trf.translation() += Vector3::UnitX() * i * 7_mm;
1223 
1224     auto childVol = std::make_unique<TrackingVolume>(
1225         trf, tbounds, "child" + std::to_string(i));
1226     childVol->assignGeometryId(
1227         GeometryIdentifier{}.withVolume(1).withLayer(i + 1));
1228     auto childNode = std::make_shared<StaticBlueprintNode>(std::move(childVol));
1229     parentNode->addChild(std::move(childNode));
1230   }
1231 
1232   cubcontainer.addChild(std::move(parentNode));
1233 
1234   auto trackingGeometry = root.construct({}, gctx, *logger);
1235 }
1236 
1237 BOOST_AUTO_TEST_SUITE_END();
1238 
1239 }  // namespace ActsTests