From 084c08890e1b5a7f01a6e6a3524277b8e2c2f7d7 Mon Sep 17 00:00:00 2001 From: Matei Climescu Date: Sat, 11 Jul 2026 12:11:39 +0200 Subject: [PATCH 1/3] feat(ms)!: Add configurable iron-block muon shield Replace the GDML-derived station stub with a configurable block list driven by MS.toml (MuonShieldConfig / readMuonShieldConfig, toml++), mirroring the calorimeter's calo.toml pattern. Each block is positioned by its upstream face and may be rotated or tapered (GeoBox, or GeoTrd when tapered). The shipped MS.toml reproduces the FairShip TRY_2026 magnets (FairShip PR #1334) as a solid-block approximation: each magnet is one iron block sized to its FairShip outer envelope, dropping the aperture, yoke shape and field (field is handled elsewhere). The FairShip absorber is the magnetised hadron stopper and stays in the Target subsystem. Place the container at the envelope centre reported by the factory instead of the hard-coded GDML z, and nest the neutrino detector inside the shield container as an embedded daughter while keeping it an independent subsystem. BREAKING CHANGE: MuonShieldFactory's public API and the muon-shield geometry and placement have changed; the old station/piece structs are removed. --- src/SHiPGeometry.cpp | 24 +- subsystems/MuonShield/CMakeLists.txt | 39 +++ subsystems/MuonShield/MS.toml | 58 ++++ subsystems/MuonShield/README.md | 132 ++++++-- .../include/MuonShield/MuonShieldConfig.h | 90 ++++++ .../include/MuonShield/MuonShieldFactory.h | 114 +++++-- .../MuonShield/src/MuonShieldConfig.cpp | 201 ++++++++++++ .../MuonShield/src/MuonShieldFactory.cpp | 289 ++++++++---------- subsystems/MuonShield/test_muonshield.cpp | 169 +++++++++- tests/test_consistency.cpp | 19 +- 10 files changed, 893 insertions(+), 242 deletions(-) create mode 100644 subsystems/MuonShield/MS.toml create mode 100644 subsystems/MuonShield/include/MuonShield/MuonShieldConfig.h create mode 100644 subsystems/MuonShield/src/MuonShieldConfig.cpp diff --git a/src/SHiPGeometry.cpp b/src/SHiPGeometry.cpp index 7c81059..4c594a0 100644 --- a/src/SHiPGeometry.cpp +++ b/src/SHiPGeometry.cpp @@ -48,22 +48,24 @@ GeoPhysVol* SHiPGeometryBuilder::build() { // Note: These are relative to the cave origin placeChild(world, target, "/SHiP/target", 1, GeoTrf::Translate3D(0.0, -14.45 * cm, 43.25 * cm)); - // Build and place MuonShieldArea - // GDML z range: 204–3148.66 cm → centre: 1676.33 cm = 16763.3 mm from world origin + // Build the muon shield, with the neutrino detector embedded inside it. + // + // The SND is an independent subsystem, but in volume terms it is a daughter + // of the muon-shield container: it sits within the shield region (WARM SND + // slot 26.40–31.50 m → centre 28.95 m). We build it first, then nest it in + // the shield container. The muon-shield block list must leave this slot + // free of iron. MuonShieldFactory muonShieldFactory(materials); - GeoPhysVol* muonShield = muonShieldFactory.build(); - placeChild(world, muonShield, "/SHiP/muon_shield", 2, - GeoTrf::Translate3D(0.0, 0.0, 16763.3 * mm)); - // Build and place the Scattering and Neutrino Detector (SND). - // Z: 26.40 to 31.50 m (WARM muon-shield configuration) → centre 28.95 m. - // The SND sits within the downstream end of the muon-shield region, so its - // envelope overlaps the muon-shield container by design (see test_consistency). NeutrinoDetectorFactory neutrinoDetectorFactory(materials); GeoPhysVol* neutrinoDetector = neutrinoDetectorFactory.build(); - placeChild(world, neutrinoDetector, "/SHiP/neutrino_detector", 9, - GeoTrf::Translate3D(0.0, 0.0, 28.95 * m)); + muonShieldFactory.embedDaughter(neutrinoDetector, 28.95 * 1000.0, "/SHiP/neutrino_detector"); + // The container is built centred on its own origin, so it is placed at the + // envelope centre reported by the factory after build(). + GeoPhysVol* muonShield = muonShieldFactory.build(); + placeChild(world, muonShield, "/SHiP/muon_shield", 2, + GeoTrf::Translate3D(0.0, 0.0, muonShieldFactory.centreZ_mm())); // Build and place UpstreamTagger (sensitive scintillator slab) // Z: 32.52 to 32.92 m → centre: 32.72 m SHiPUBTManager ubtManager; diff --git a/subsystems/MuonShield/CMakeLists.txt b/subsystems/MuonShield/CMakeLists.txt index 16222ae..8507094 100644 --- a/subsystems/MuonShield/CMakeLists.txt +++ b/subsystems/MuonShield/CMakeLists.txt @@ -1,4 +1,43 @@ # SPDX-License-Identifier: LGPL-3.0-or-later # Copyright (C) CERN for the benefit of the SHiP Collaboration +# Shared recipe: library (globs src/*.cpp), include dirs, kernel link, Catch test. ship_add_subsystem(MuonShield) + +# --- Extras: MuonShield parses MS.toml via MuonShieldConfig (toml++), same as +# the calorimeter. Everything below mirrors subsystems/Calorimeter/CMakeLists.txt. + +find_package(tomlplusplus REQUIRED) + +# toml++ is a build-time-only dependency (used only inside MuonShieldConfig.cpp, +# never in a public header). BUILD_INTERFACE keeps it out of the install export. +target_link_libraries( + MuonShield + PRIVATE $ +) + +# Stage MS.toml next to the tests and into the top-level build dir (for +# test_builder / OverlapCheck), so it is found when running from either. +configure_file( + ${CMAKE_CURRENT_SOURCE_DIR}/MS.toml + ${CMAKE_CURRENT_BINARY_DIR}/MS.toml + COPYONLY +) +configure_file( + ${CMAKE_CURRENT_SOURCE_DIR}/MS.toml + ${CMAKE_BINARY_DIR}/MS.toml + COPYONLY +) + +# Absolute fallbacks so the factory finds MS.toml when the CWD does not have it. +target_compile_definitions( + MuonShield + PRIVATE + MS_TOML_DEFAULT_PATH="${CMAKE_CURRENT_SOURCE_DIR}/MS.toml" + MS_TOML_INSTALL_PATH="${CMAKE_INSTALL_PREFIX}/${CMAKE_INSTALL_DATADIR}/SHiPGeometry/MS.toml" +) + +install( + FILES ${CMAKE_CURRENT_SOURCE_DIR}/MS.toml + DESTINATION ${CMAKE_INSTALL_DATADIR}/SHiPGeometry +) diff --git a/subsystems/MuonShield/MS.toml b/subsystems/MuonShield/MS.toml new file mode 100644 index 0000000..73d8bb9 --- /dev/null +++ b/subsystems/MuonShield/MS.toml @@ -0,0 +1,58 @@ +# SPDX-License-Identifier: LGPL-3.0-or-later +# Copyright (C) CERN for the benefit of the SHiP Collaboration +# +# Muon shield geometry — solid-block approximation of the FairShip TRY_2026 +# magnet layout (FairShip PR #1334). Each magnet is one solid iron block from +# its outer envelope (dX*(1+ratio)+midGap+ceil(max(100/dY,gap)) x dY+dY_yoke); +# the aperture, yoke shape and field are dropped (field handled elsewhere). +# +# The muon shield is defined INDEPENDENTLY of the neutrino detector. The SND +# reservation is declared in NeutrinoDetector/SD.toml and subtracted from the +# iron at build time (A - B), so this file needs no knowledge of the SND. +# +# The FairShip absorber (row 0) is the magnetised hadron stopper and lives in +# the Target subsystem (2.14-4.44 m); the shield front is at 2.14 m, so magnet 1 +# begins at 4.59 m. + +block_material = "Iron" + +envelope_half_x_mm = 1760.0 +envelope_half_y_mm = 1320.0 +envelope_z_start_m = 4.5400 +envelope_z_end_m = 32.0800 + +# magnet 1 +[[block]] +start = [0.0, 0.0, 4590.0] +size = [2720.0, 1600.0, 3000.0] + +# magnet 2 +[[block]] +start = [0.0, 0.0, 7740.0] +size = [2980.0, 1780.0, 4500.0] + +# magnet 3 +[[block]] +start = [0.0, 0.0, 12400.0] +size = [3104.0, 1700.0, 4500.0] + +# magnet 4 +[[block]] +start = [0.0, 0.0, 17050.0] +size = [3184.0, 1200.0, 4500.0] + +# magnet 5 +[[block]] +start = [0.0, 0.0, 21710.0] +size = [2680.0, 1040.0, 2620.0] + +# magnet 6 (widens in Y toward downstream) +[[block]] +start = [0.0, 0.0, 24480.0] +size = [3000.0, 1060.0, 3360.0] +taper = [0.0000, 1.7047] + +# magnet 7 +[[block]] +start = [0.0, 0.0, 28030.0] +size = [3420.0, 2540.0, 4000.0] diff --git a/subsystems/MuonShield/README.md b/subsystems/MuonShield/README.md index 2f5fbbf..2c4b096 100644 --- a/subsystems/MuonShield/README.md +++ b/subsystems/MuonShield/README.md @@ -1,47 +1,125 @@ # MuonShield -Active muon shield magnets for background suppression. +Passive iron muon shield for background suppression. ## Description -The MuonShield subsystem implements 6 magnet stations that deflect muons away from the detector acceptance. The current implementation uses bounding-box approximations of the GDML arb8 (trapezoidal) shapes — each station contains 8 iron pieces (left/right middle magnets, left/right return yokes, and 4 top/bottom corner pieces). +The MuonShield subsystem reserves an Air **envelope** and places an explicit +list of iron **blocks** inside it. Each block is a box, or a symmetric +trapezoid/frustum (`GeoTrd`) when tapered, anchored by its upstream face and +optionally rotated. The block list is the sole source of iron: holes and +apertures are represented by leaving a region free of blocks, or by arranging +several blocks around the gap. + +The default `MS.toml` is the **FairShip TRY_2026** magnet layout (FairShip +PR #1334) as a solid-block approximation: the 7 shield magnets, each reproduced +as one iron block from its outer envelope (aperture, yoke and field dropped). +It is generated by `fairship_to_mstoml.py`, which converts the FairShip +`ShipMuonShield` `params` table into blocks. The FairShip absorber (row 0) is +the magnetised hadron stopper and lives in the Target subsystem, so it is not +included here. ## Geometry Tree +The default `MS.toml` places the 7 TRY_2026 magnets as solid blocks. The SND +reservation (`NeutrinoDetector/SD.toml`) is subtracted from the iron at build +time, and the SND is embedded as a daughter that sits in the resulting cavity: + +```text +muon_shield (Air envelope, 3520 × 2640 × 27540 mm) centre z = 18.31 m + ├─ block_0 magnet 1 (Iron box) + ├─ … magnets 2–5 (Iron boxes; magnet 6 is a GeoTrd, widens in Y) + ├─ block_5 magnet 6 (Iron − SND box) ┐ carved by A − B + ├─ block_6 magnet 7 (Iron − SND box) ┘ (magnets the reservation intersects) + └─ neutrino_detector (Air box, 26.40 → 31.50 m) embedded, sits in the cavity ``` -MuonShieldArea (Air, 3620×3400×29448 mm) - ├─ MagnAbsorb_container (Air) z = -13568.3 mm - │ └─ 8 × Iron boxes (MiddleMagL/R, MagRetL/R, MagTop/BotLeft/Right) - ├─ Magn1_container (Air) z = -7263.3 mm - │ └─ 8 × Iron boxes - ├─ Magn2_container (Air) z = +591.5 mm - │ └─ 8 × Iron boxes - ├─ Magn3_container (Air) z = +5821.6 mm - │ └─ 8 × Iron boxes - ├─ Magn4_container (Air) z = +9096.9 mm - │ └─ 8 × Iron boxes - └─ Magn5_container (Air) z = +12385.1 mm - └─ 8 × Iron boxes + +The envelope values in `MS.toml` are generated to contain the magnets +(3520 × 2640 mm, 4.54–32.08 m); the container is placed in the world at its Z +centre (18.31 m) by `SHiPGeometryBuilder`. The SND cavity is a real Boolean +subtraction (like FairShip's `SetSNDSpace`): the reservation box is subtracted +from whichever magnets it intersects — so which magnets get carved follows from +the box position, not a hard-coded list. See the **Reserved space** section +below. + +## Configuration (`MS.toml`) + +| Key | Fallback | Meaning | +|----------------------|----------|------------------------------------------------| +| `block_material` | "Iron" | Absorber material (must exist in SHiPMaterials) | +| `envelope_half_x_mm` | 2100 | Container half-X (mm) | +| `envelope_half_y_mm` | 2300 | Container half-Y (mm) | +| `envelope_z_start_m` | 4.74 | Envelope start along Z (m, world coords) | +| `envelope_z_end_m` | 32.22 | Envelope end along Z (m, world coords) | +| `[[block]]` | — | One table per iron block (see below) | + +The **Fallback** column is the parser default used only when a key is *omitted*. +The shipped `MS.toml` sets the envelope explicitly (half-sizes 1760 × 1320 mm, +z = 4.54–32.08 m — the values in the geometry tree above), generated by +`fairship_to_mstoml.py` to contain the blocks. The **envelope** is the Air +container and the subsystem's world placement. The parser rejects an +inverted/zero envelope and any block with a non-positive size, a taper that +collapses the downstream face, or a face outside the envelope in Z. + +### Blocks + +Each `[[block]]` is a table with these fields (all lengths in mm, angles in +degrees, world/beamline coordinates): + +```toml +[[block]] +start = [0.0, 0.0, 4740.0] # centre of the UPSTREAM (-z) face +size = [3000.0, 2000.0, 21660.0] # near-face full x, y, and z length +rotation = [0.0, 2.0, 0.0] # optional: extrinsic X→Y→Z about `start` +taper = [5.0, 0.0] # optional: symmetric half-opening angles ``` -Position in world: z = 16763.3 mm (centre of MuonShieldArea). +- `start` is the centre of the block's upstream (−z) face; the block extends + downstream by `size[2]` before rotation. +- `rotation` (default `[0,0,0]`) is applied about `start`, extrinsically about + the world X then Y then Z axes. +- `taper` (default `[0,0]`) gives symmetric half-opening angles: the section + widens (or, if negative, narrows) toward the downstream face, so the far + half-x is `size[0]/2 + size[2]·tan(taper[0])` and similarly for Y. `[0,0]` + yields a `GeoBox`; otherwise a `GeoTrd`. + +### Reserved space (SND cavity) + +The shield iron is carved by Boolean subtraction (A − B) so external detectors +can define their space independently. `MuonShieldFactory::reserveSpace(centre, +size, rotation)` registers a box (world coords, mm); at build time every magnet +the box intersects has it subtracted (`GeoShapeSubtraction`), while other +magnets are untouched. The **neutrino detector** reservation is declared in +`NeutrinoDetector/SD.toml` and wired in by `SHiPGeometryBuilder`, which reads the +envelope, calls `reserveSpace(...)`, then `embedDaughter(...)` to place the SND +in the cavity. Because the carve follows from the box position, the muon shield +never names which magnets hold the SND — moving the box in `SD.toml` moves the +cavity. + +### Embedded daughters + +Other subsystems can be nested inside the shield container while remaining +independent subsystems. `MuonShieldFactory::embedDaughter(vol, worldCentreZ_mm, +name)` registers a pre-built volume, placed at build time at the given world-Z +inside the container. The SND keeps its own factory, config, and +`/SHiP/neutrino_detector` naming — only its position in the volume tree changes +(it becomes a daughter of `/SHiP/muon_shield` rather than a direct child of the +world), sitting in the cavity carved by its `reserveSpace` reservation. ## Materials -| Material | Density | Usage | -|----------|-----------|------------------| -| Air | 1.29 mg/cm³ | Container volumes | -| Iron | 7.87 g/cm³ | Magnet pieces | +| Material | Density | Usage | +|----------|-------------|-------------------| +| Air | 1.29 mg/cm³ | Container volume | +| Iron | 7.87 g/cm³ | Shield blocks | ## Status -- [x] C++ implementation (box approximations) -- [ ] Replace boxes with proper GeoTrap/arb8 shapes -- [ ] Verify field map integration points -- [ ] Verification against GDML +- [x] Explicit block list driven by `MS.toml` (position, rotation, taper) +- [ ] Magnetic field (handled elsewhere; blocks are passive iron) +- [ ] Populate the block list from the engineering reference layout ## TODO -- Replace bounding-box iron pieces with proper trapezoidal (GeoTrap or GeoGenericTrap) shapes matching the GDML arb8 vertices -- Add magnetic field regions (currently geometry only, no field) -- Verify station z-positions and piece dimensions against GDML reference +- Replace the placeholder default blocks with the real muon-shield magnet + layout supplied by the subsystem coordinator. diff --git a/subsystems/MuonShield/include/MuonShield/MuonShieldConfig.h b/subsystems/MuonShield/include/MuonShield/MuonShieldConfig.h new file mode 100644 index 0000000..62921b4 --- /dev/null +++ b/subsystems/MuonShield/include/MuonShield/MuonShieldConfig.h @@ -0,0 +1,90 @@ +// SPDX-License-Identifier: LGPL-3.0-or-later +// Copyright (C) CERN for the benefit of the SHiP Collaboration + +#pragma once + +#include +#include +#include + +namespace SHiPGeometry { + +/// Degrees-to-radians conversion, shared by the config parser and the factory. +inline constexpr double kDegToRad = 3.14159265358979323846 / 180.0; + +/** + * @brief A single muon-shield iron block. + * + * A block is a box (or, when tapered, a symmetric trapezoid/frustum) placed by + * its upstream (−z) face. All lengths are in millimetres and all angles in + * degrees, in world/beamline coordinates (target-front-face origin). + * + * Geometry conventions: + * - @c start is the centre of the block's upstream face; the block extends + * downstream (+z) by @c size[2] before rotation. + * - @c rotation is applied about @c start, extrinsically about the world X, Y + * then Z axes. + * - @c size = {sx, sy, sz}: the full X and Y of the *upstream* face and the Z + * length. + * - @c taper = {ax, ay}: symmetric half-opening angles. The cross-section + * grows (or, if negative, shrinks) toward the downstream face, so the far + * half-width is (sx/2 + sz·tan(ax)) and similarly for Y. {0, 0} → a box. + */ +struct MuonShieldBlock { + std::array start = {0.0, 0.0, 0.0}; ///< upstream-face centre (mm) + std::array size = {0.0, 0.0, 0.0}; ///< near-face sx, sy and length sz (mm) + std::array rotation_deg = {0.0, 0.0, 0.0}; ///< rotation about x, y, z (deg) + std::array taper_deg = {0.0, 0.0}; ///< half-opening angles ax, ay (deg) +}; + +/** + * @brief Configuration for the SHiP muon shield geometry. + * + * Populated by readMuonShieldConfig() from an MS.toml file. + * + * The shield is an explicit list of iron @c blocks placed inside an Air + * @c envelope container. The envelope sizes the container and fixes the + * subsystem's world placement (its Z centre); the blocks are the only iron. + * Holes/apertures are represented simply by the absence of a block, or by + * arranging several blocks around the gap. + */ +struct MuonShieldConfig { + // Absorber material for all blocks (must exist in SHiPMaterials). + std::string block_material = "Iron"; + + // Container (envelope) transverse half-sizes (mm). + double envelope_half_x_mm = 2100.0; + double envelope_half_y_mm = 2300.0; + + // Envelope Z extent in world coordinates (target-front-face origin, m). + // This is the Air container and the world placement of the subsystem. + double envelope_z_start_m = 4.74; + double envelope_z_end_m = 32.22; + + // The iron blocks. Empty = no iron. + std::vector blocks; + + // ── Derived helpers ───────────────────────────────────────────────── + /// Full Z length of the envelope (mm). + double envelopeLengthZ_mm() const { return (envelope_z_end_m - envelope_z_start_m) * 1000.0; } + /// World-Z centre of the envelope (mm). + double envelopeCentreZ_mm() const { + return 0.5 * (envelope_z_start_m + envelope_z_end_m) * 1000.0; + } +}; + +/** + * @brief Parse an MS.toml file and return a MuonShieldConfig. + * + * Uses toml++ for parsing. Unknown top-level keys are reported via stderr + * (helpful when a stale or mistyped key would otherwise be silently ignored), + * but do not cause the parse to fail. + * + * @throws std::runtime_error if the file cannot be opened, contains malformed + * TOML, has a non-positive or inverted Z envelope, or a block with a + * non-positive size, a taper that collapses the far face, or an + * upstream face outside the envelope in Z. + */ +MuonShieldConfig readMuonShieldConfig(const std::string& path); + +} // namespace SHiPGeometry diff --git a/subsystems/MuonShield/include/MuonShield/MuonShieldFactory.h b/subsystems/MuonShield/include/MuonShield/MuonShieldFactory.h index 366aad8..7e0f650 100644 --- a/subsystems/MuonShield/include/MuonShield/MuonShieldFactory.h +++ b/subsystems/MuonShield/include/MuonShield/MuonShieldFactory.h @@ -3,60 +3,110 @@ #pragma once -#include +#include +#include +#include class GeoPhysVol; namespace SHiPGeometry { class SHiPMaterials; +struct MuonShieldConfig; /** - * @brief Factory for the MuonShield (muon shield magnets) geometry + * @brief Factory for the MuonShield geometry. * - * Creates a MuonShieldArea (Air) container holding 6 magnet stations. - * Each station contains 8 Iron pieces — bounding-box approximations of the - * GDML arb8 shapes. Station z-positions (from GDML MuonShieldArea origin, cm): - * MagnAbsorb: 319.5, Magn1: 950, Magn2: 1735.48, - * Magn3: 2258.49, Magn4: 2586.02, Magn5: 2914.84 + * Builds an Air container spanning the muon-shield Z envelope and places an + * explicit list of iron blocks inside it (from MS.toml). Each block is a box, + * or a symmetric trapezoid/frustum (GeoTrd) when tapered, anchored by its + * upstream face and optionally rotated. + * + * Other subsystems can be embedded as daughter volumes of the shield container + * via embedDaughter() — used, for example, to nest the neutrino detector inside + * the muon shield while keeping it an independent subsystem. + * + * The container is built centred on its own origin; SHiPGeometryBuilder places + * it in the world at centreZ_mm() (the envelope centre from MS.toml). + * + * The config file is resolved at build() time: + * 1. "MS.toml" relative to the current working directory (works when running + * from the build directory, where CMake stages the file). + * 2. The absolute source-tree path baked in via MS_TOML_DEFAULT_PATH. + * 3. The installed data-dir path via MS_TOML_INSTALL_PATH. */ class MuonShieldFactory { public: - explicit MuonShieldFactory(SHiPMaterials& materials); + explicit MuonShieldFactory(SHiPMaterials& materials, std::string configPath = "MS.toml"); ~MuonShieldFactory() = default; - [[nodiscard]] GeoPhysVol* build(); + /** + * @brief Embed an external volume as a daughter of the shield container. + * + * Registers @p daughter to be placed inside the shield at world-Z + * @p worldCentreZ_mm. The daughter remains an independent subsystem; this + * only nests it in the volume tree. Since the iron is defined by an explicit + * block list, the caller is responsible for not listing blocks over the + * daughter's footprint. Call before build(). + * + * @param daughter Pre-built volume (centred on its own origin). + * @param worldCentreZ_mm - private: - SHiPMaterials& m_materials; - // Unit shorthand (GeoModel's native length unit is mm) - static constexpr double mm = GeoModelKernelUnits::mm; + centre where the daughter is placed (mm). + * @param name Name tag for the placement. + */ + void embedDaughter(GeoPhysVol* daughter, double worldCentreZ_mm, const std::string& name); - struct PieceData { - double halfX, halfY, halfZ; // bounding-box half-sizes - double centX, centY; // centre offset in station XY frame - const char* name; // piece name suffix - }; + /** + * @brief Reserve (carve) a box out of the shield iron via Boolean A - B. + * + * Registers a box to subtract from every magnet it intersects, so external + * detectors (e.g. the neutrino detector, from SD.toml) can define their + * space independently and have it removed from the iron at build() time. + * Positions are in world/beamline coordinates. Call before build(). + * + * @param worldCentre_mm Box centre in world coordinates (mm). + * @param size_mm Full box dimensions {x, y, z} (mm). + * @param rotation_deg Optional extrinsic X->Y->Z rotation about the centre. + */ + void reserveSpace(const std::array& worldCentre_mm, + const std::array& size_mm, + const std::array& rotation_deg = {0.0, 0.0, 0.0}); - struct StationData { - const char* name; // station name - double stationZ; // z in MuonShieldArea frame - double containerHalfX; - double containerHalfY; - double containerHalfZ; - PieceData pieces[8]; + /** Build and return the muon-shield container volume (centred on origin). */ + GeoPhysVol* build(); + + /** + * @brief World-Z centre (mm) at which the container should be placed. + * + * Valid only after build() has been called (it is read from MS.toml). + */ + double centreZ_mm() const { return m_centreZ_mm; } + + /** Return the config path that will actually be opened (after resolution). */ + std::string resolvedConfigPath() const; + + private: + struct EmbeddedDaughter { + GeoPhysVol* volume; + double worldCentreZ_mm; + std::string name; }; - // GDML-derived station data for all 6 stations × 8 pieces - static const StationData k_stations[6]; + struct ReservedBox { + std::array centre_mm; + std::array size_mm; + std::array rotation_deg; + }; - GeoPhysVol* buildStation(const StationData& station); + SHiPMaterials& m_materials; + std::string m_configPath; + std::vector m_daughters; + std::vector m_reservations; - // MuonShieldArea container dimensions - static constexpr double s_areaHalfX = 1810.0 * mm; - static constexpr double s_areaHalfY = 1700.0 * mm; - static constexpr double s_areaHalfZ = 14724.0 * mm; + // World-Z centre of the container, populated by build() from the config. + double m_centreZ_mm = 0.0; }; } // namespace SHiPGeometry diff --git a/subsystems/MuonShield/src/MuonShieldConfig.cpp b/subsystems/MuonShield/src/MuonShieldConfig.cpp new file mode 100644 index 0000000..67ee191 --- /dev/null +++ b/subsystems/MuonShield/src/MuonShieldConfig.cpp @@ -0,0 +1,201 @@ +// SPDX-License-Identifier: LGPL-3.0-or-later +// Copyright (C) CERN for the benefit of the SHiP Collaboration +// +// Muon shield configuration parser. Reads an MS.toml file and populates a +// MuonShieldConfig struct. Built on toml++ +// (https://github.com/marzer/tomlplusplus, single-header MIT, provided by the +// tomlplusplus package), so the file format is standard TOML. +// +// Behaviour notes +// =============== +// * Unknown top-level keys are NOT silently ignored. Each one is reported on +// stderr at parse time. This catches typos and stale keys. +// +// * Numeric fields accept both TOML integers and floats. + +#include "MuonShield/MuonShieldConfig.h" + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +namespace SHiPGeometry { + +namespace { + +using namespace std::string_view_literals; + +// Recognised top-level keys (sorted for binary search). Anything outside this +// set triggers a warning. "block" is the [[block]] array of tables. +static constexpr std::array kKnownKeys = { + "block"sv, + "block_material"sv, + "envelope_half_x_mm"sv, + "envelope_half_y_mm"sv, + "envelope_z_end_m"sv, + "envelope_z_start_m"sv, +}; + +// Mapping from TOML key name to MuonShieldConfig double member pointer. +struct NumericField { + const char* key; + double MuonShieldConfig::* member; +}; + +static constexpr NumericField kNumericFields[] = { + {"envelope_half_x_mm", &MuonShieldConfig::envelope_half_x_mm}, + {"envelope_half_y_mm", &MuonShieldConfig::envelope_half_y_mm}, + {"envelope_z_start_m", &MuonShieldConfig::envelope_z_start_m}, + {"envelope_z_end_m", &MuonShieldConfig::envelope_z_end_m}, +}; + +// Extract a TOML numeric (integer or float) as a double, if the node holds one. +std::optional asDouble(const toml::node* node) { + if (node) { + if (auto d = node->value()) + return *d; + if (auto i = node->value()) + return static_cast(*i); + } + return std::nullopt; +} + +// Read a double or integer as a double. +double readNumeric(const toml::node_view& node, const std::string& key) { + if (auto v = asDouble(node.node())) + return *v; + throw std::runtime_error("MuonShieldConfig: '" + key + "' must be a number"); +} + +// Read a single toml node as a double (int or float). +double nodeToNumber(const toml::node* node, const std::string& context, const std::string& path) { + if (auto v = asDouble(node)) + return *v; + throw std::runtime_error("MuonShieldConfig: " + context + " values must be numbers in " + path); +} + +// Read a fixed-length numeric array (e.g. size = [x, y, z]) from a block table. +template +std::array readFixedArray(const toml::table& block, const char* key, + const std::string& path, bool required, + const std::array& fallback) { + auto node = block[key]; + if (!node) { + if (required) + throw std::runtime_error("MuonShieldConfig: block is missing required '" + + std::string(key) + "' in " + path); + return fallback; + } + const toml::array* arr = node.as_array(); + if (!arr || arr->size() != N) + throw std::runtime_error("MuonShieldConfig: block '" + std::string(key) + "' must be an " + + std::to_string(N) + "-element array in " + path); + std::array out{}; + for (std::size_t i = 0; i < N; ++i) + out[i] = nodeToNumber(arr->get(i), "block '" + std::string(key) + "'", path); + return out; +} + +} // namespace + +MuonShieldConfig readMuonShieldConfig(const std::string& path) { + MuonShieldConfig cfg; + + toml::table table; + try { + table = toml::parse_file(path); + } catch (const toml::parse_error& e) { + throw std::runtime_error("MuonShieldConfig: failed to parse " + path + ": " + + std::string(e.description())); + } + + // First pass: warn about unknown keys. + for (const auto& [k, _] : table) { + if (!std::ranges::binary_search(kKnownKeys, std::string_view{k})) { + std::cerr << "MuonShieldConfig: warning: unknown key '" << k << "' in " << path + << " (typo? stale field? — value will be ignored)\n"; + } + } + + // Numeric (double) envelope fields. + for (const auto& [key, member] : kNumericFields) + if (auto n = table[key]; n) + cfg.*member = readNumeric(n, key); + + // String field. + if (auto n = table["block_material"]; n) { + if (auto s = n.value(); s) + cfg.block_material = *s; + else + throw std::runtime_error("MuonShieldConfig: 'block_material' must be a string"); + } + + // Blocks: an array of tables ([[block]]). + if (auto blocksNode = table["block"]; blocksNode) { + const toml::array* arr = blocksNode.as_array(); + if (!arr) + throw std::runtime_error("MuonShieldConfig: 'block' must be an array of tables in " + + path); + for (const auto& elem : *arr) { + const toml::table* bt = elem.as_table(); + if (!bt) + throw std::runtime_error( + "MuonShieldConfig: each [[block]] entry must be a table in " + path); + MuonShieldBlock block; + block.start = readFixedArray<3>(*bt, "start", path, true, {0.0, 0.0, 0.0}); + block.size = readFixedArray<3>(*bt, "size", path, true, {0.0, 0.0, 0.0}); + block.rotation_deg = readFixedArray<3>(*bt, "rotation", path, false, {0.0, 0.0, 0.0}); + block.taper_deg = readFixedArray<2>(*bt, "taper", path, false, {0.0, 0.0}); + cfg.blocks.push_back(block); + } + } + + // ── Validation ────────────────────────────────────────────────────── + if (cfg.envelope_z_end_m <= cfg.envelope_z_start_m) + throw std::runtime_error( + "MuonShieldConfig: envelope_z_end_m must be greater than envelope_z_start_m in " + + path); + if (cfg.envelope_half_x_mm <= 0.0 || cfg.envelope_half_y_mm <= 0.0) + throw std::runtime_error( + "MuonShieldConfig: envelope_half_x_mm and envelope_half_y_mm must be positive in " + + path); + + const double envStartMm = cfg.envelope_z_start_m * 1000.0; + const double envEndMm = cfg.envelope_z_end_m * 1000.0; + constexpr double kEps = 1e-6; // mm + + for (std::size_t i = 0; i < cfg.blocks.size(); ++i) { + const auto& b = cfg.blocks[i]; + const std::string where = "block " + std::to_string(i) + " in " + path; + if (b.size[0] <= 0.0 || b.size[1] <= 0.0 || b.size[2] <= 0.0) + throw std::runtime_error("MuonShieldConfig: " + where + " has a non-positive size"); + + // Taper must not collapse (or invert) the downstream face. + const double farHalfX = 0.5 * b.size[0] + b.size[2] * std::tan(b.taper_deg[0] * kDegToRad); + const double farHalfY = 0.5 * b.size[1] + b.size[2] * std::tan(b.taper_deg[1] * kDegToRad); + if (farHalfX <= 0.0 || farHalfY <= 0.0) + throw std::runtime_error("MuonShieldConfig: " + where + + " has a taper that collapses its downstream face"); + + // The block must sit within the envelope in Z (light sanity check; + // rotations are not accounted for here). + if (b.start[2] < envStartMm - kEps || b.start[2] > envEndMm + kEps) + throw std::runtime_error("MuonShieldConfig: " + where + " has its upstream face (z = " + + std::to_string(b.start[2]) + " mm) outside the envelope"); + if (b.start[2] + b.size[2] > envEndMm + kEps) + throw std::runtime_error( + "MuonShieldConfig: " + where + " has its downstream face (z = " + + std::to_string(b.start[2] + b.size[2]) + " mm) beyond the envelope end"); + } + + return cfg; +} + +} // namespace SHiPGeometry diff --git a/subsystems/MuonShield/src/MuonShieldFactory.cpp b/subsystems/MuonShield/src/MuonShieldFactory.cpp index f77bd59..d1cf46e 100644 --- a/subsystems/MuonShield/src/MuonShieldFactory.cpp +++ b/subsystems/MuonShield/src/MuonShieldFactory.cpp @@ -3,182 +3,159 @@ #include "MuonShield/MuonShieldFactory.h" +#include "MuonShield/MuonShieldConfig.h" +#include "SHiPGeometry/ConfigPath.h" #include "SHiPGeometry/SHiPMaterials.h" #include #include #include #include +#include #include #include +#include +#include +#include #include +#include +#include +#include #include +// Absolute fallback path baked in by CMake so out-of-source builds always find +// MS.toml even when the CWD doesn't contain a copy of it. +#ifndef MS_TOML_DEFAULT_PATH +#define MS_TOML_DEFAULT_PATH "MS.toml" +#endif +// Install-time data directory path, set by CMake during install configuration. +#ifndef MS_TOML_INSTALL_PATH +#define MS_TOML_INSTALL_PATH "" +#endif + namespace SHiPGeometry { -// --------------------------------------------------------------------------- -// GDML-derived station data (all dimensions in mm) -// Bounding-box centres and half-sizes are computed from the arb8 vertices in -// ship_geometry.gdml. Station z-positions are relative to the MuonShieldArea -// container centre (GDML station_z_cm - 1676.33 cm) × 10 mm/cm. -// --------------------------------------------------------------------------- -const MuonShieldFactory::StationData MuonShieldFactory::k_stations[6] = { - // ── MagnAbsorb (GDML z = 319.5 cm, dz = 115.5 cm) ────────────────── - {"magn_absorb", - -13568.3 * mm, - 1020.0 * mm, - 1691.0 * mm, - 1155.0 * mm, - { - {250.0 * mm, 1690.0 * mm, 1155.0 * mm, 250.0 * mm, -1.0 * mm, "middle_mag_l"}, - {250.0 * mm, 1690.0 * mm, 1155.0 * mm, -250.0 * mm, 1.0 * mm, "middle_mag_r"}, - {250.0 * mm, 1690.0 * mm, 1155.0 * mm, 770.0 * mm, 0.0, "mag_ret_l"}, - {250.0 * mm, 1690.0 * mm, 1155.0 * mm, -770.0 * mm, 0.0, "mag_ret_r"}, - {510.0 * mm, 250.0 * mm, 1155.0 * mm, 510.0 * mm, 1440.0 * mm, "mag_top_left"}, - {510.0 * mm, 250.0 * mm, 1155.0 * mm, -510.0 * mm, 1440.0 * mm, "mag_top_right"}, - {510.0 * mm, 250.0 * mm, 1155.0 * mm, 510.0 * mm, -1440.0 * mm, "mag_bot_left"}, - {510.0 * mm, 250.0 * mm, 1155.0 * mm, -510.0 * mm, -1440.0 * mm, "mag_bot_right"}, - }}, - - // ── Magn1 (GDML z = 950 cm, dz = 495 cm) ──────────────────────────── - {"magn_1", - -7263.3 * mm, - 1697.0 * mm, - 1230.0 * mm, - 4950.0 * mm, - { - {399.6 * mm, 1228.2 * mm, 4950.0 * mm, 399.6 * mm, 0.0, "middle_mag_l"}, - {399.6 * mm, 1228.2 * mm, 4950.0 * mm, -399.6 * mm, 0.0, "middle_mag_r"}, - {487.7 * mm, 1229.2 * mm, 4950.0 * mm, 1208.7 * mm, 0.0, "mag_ret_l"}, - {487.7 * mm, 1229.2 * mm, 4950.0 * mm, -1208.7 * mm, 0.0, "mag_ret_r"}, - {848.2 * mm, 479.6 * mm, 4950.0 * mm, 848.2 * mm, 749.6 * mm, "mag_top_left"}, - {848.2 * mm, 479.6 * mm, 4950.0 * mm, -848.2 * mm, 749.6 * mm, "mag_top_right"}, - {848.2 * mm, 479.6 * mm, 4950.0 * mm, 848.2 * mm, -749.6 * mm, "mag_bot_left"}, - {848.2 * mm, 479.6 * mm, 4950.0 * mm, -848.2 * mm, -749.6 * mm, "mag_bot_right"}, - }}, - - // ── Magn2 (GDML z = 1735.48 cm, dz = 280.48 cm) ───────────────────── - {"magn_2", - 591.5 * mm, - 1736.0 * mm, - 1056.0 * mm, - 2804.8 * mm, - { - {265.6 * mm, 1054.6 * mm, 2804.8 * mm, 265.6 * mm, 0.0, "middle_mag_l"}, - {265.6 * mm, 1054.6 * mm, 2804.8 * mm, -265.6 * mm, 0.0, "middle_mag_r"}, - {594.7 * mm, 1055.6 * mm, 2804.8 * mm, 1140.3 * mm, 0.0, "mag_ret_l"}, - {594.7 * mm, 1055.6 * mm, 2804.8 * mm, -1140.3 * mm, 0.0, "mag_ret_r"}, - {867.5 * mm, 312.8 * mm, 2804.8 * mm, 867.5 * mm, 742.8 * mm, "mag_top_left"}, - {867.5 * mm, 312.8 * mm, 2804.8 * mm, -867.5 * mm, 742.8 * mm, "mag_top_right"}, - {867.5 * mm, 312.8 * mm, 2804.8 * mm, 867.5 * mm, -742.8 * mm, "mag_bot_left"}, - {867.5 * mm, 312.8 * mm, 2804.8 * mm, -867.5 * mm, -742.8 * mm, "mag_bot_right"}, - }}, - - // ── Magn3 (GDML z = 2258.49 cm, dz = 232.53 cm) ───────────────────── - {"magn_3", - 5821.6 * mm, - 1781.0 * mm, - 597.0 * mm, - 2325.3 * mm, - { - {18.4 * mm, 595.8 * mm, 2325.3 * mm, 23.4 * mm, 0.0, "middle_mag_l"}, - {18.4 * mm, 595.8 * mm, 2325.3 * mm, -23.4 * mm, 0.0, "middle_mag_r"}, - {849.3 * mm, 596.8 * mm, 2325.3 * mm, 931.6 * mm, 0.0, "mag_ret_l"}, - {849.3 * mm, 596.8 * mm, 2325.3 * mm, -931.6 * mm, 0.0, "mag_ret_r"}, - {888.0 * mm, 18.4 * mm, 2325.3 * mm, 893.0 * mm, 578.4 * mm, "mag_top_left"}, - {888.0 * mm, 18.4 * mm, 2325.3 * mm, -893.0 * mm, 578.4 * mm, "mag_top_right"}, - {888.0 * mm, 18.4 * mm, 2325.3 * mm, 893.0 * mm, -578.4 * mm, "mag_bot_left"}, - {888.0 * mm, 18.4 * mm, 2325.3 * mm, -893.0 * mm, -578.4 * mm, "mag_bot_right"}, - }}, - - // ── Magn4 (GDML z = 2586.02 cm, dz = 85 cm) ───────────────────────── - {"magn_4", - 9096.9 * mm, - 1797.0 * mm, - 1332.0 * mm, - 850.0 * mm, - { - {535.6 * mm, 1330.2 * mm, 850.0 * mm, 535.6 * mm, 0.0, "middle_mag_l"}, - {535.6 * mm, 1330.2 * mm, 850.0 * mm, -535.6 * mm, 0.0, "middle_mag_r"}, - {713.0 * mm, 1331.2 * mm, 850.0 * mm, 1083.0 * mm, 0.0, "mag_ret_l"}, - {713.0 * mm, 1331.2 * mm, 850.0 * mm, -1083.0 * mm, 0.0, "mag_ret_r"}, - {898.0 * mm, 385.6 * mm, 850.0 * mm, 898.0 * mm, 945.6 * mm, "mag_top_left"}, - {898.0 * mm, 385.6 * mm, 850.0 * mm, -898.0 * mm, 945.6 * mm, "mag_top_right"}, - {898.0 * mm, 385.6 * mm, 850.0 * mm, 898.0 * mm, -945.6 * mm, "mag_bot_left"}, - {898.0 * mm, 385.6 * mm, 850.0 * mm, -898.0 * mm, -945.6 * mm, "mag_bot_right"}, - }}, - - // ── Magn5 (GDML z = 2914.84 cm, dz = 233.82 cm) ───────────────────── - {"magn_5", - 12385.1 * mm, - 1808.0 * mm, - 960.0 * mm, - 2338.2 * mm, - { - {200.0 * mm, 959.0 * mm, 2338.2 * mm, 200.0 * mm, 0.0, "middle_mag_l"}, - {200.0 * mm, 959.0 * mm, 2338.2 * mm, -200.0 * mm, 0.0, "middle_mag_r"}, - {728.9 * mm, 960.0 * mm, 2338.2 * mm, 1079.2 * mm, 0.0, "mag_ret_l"}, - {728.9 * mm, 960.0 * mm, 2338.2 * mm, -1079.2 * mm, 0.0, "mag_ret_r"}, - {904.0 * mm, 200.0 * mm, 2338.2 * mm, 904.0 * mm, 760.0 * mm, "mag_top_left"}, - {904.0 * mm, 200.0 * mm, 2338.2 * mm, -904.0 * mm, 760.0 * mm, "mag_top_right"}, - {904.0 * mm, 200.0 * mm, 2338.2 * mm, 904.0 * mm, -760.0 * mm, "mag_bot_left"}, - {904.0 * mm, 200.0 * mm, 2338.2 * mm, -904.0 * mm, -760.0 * mm, "mag_bot_right"}, - }}, -}; - -// --------------------------------------------------------------------------- - -MuonShieldFactory::MuonShieldFactory(SHiPMaterials& materials) : m_materials(materials) {} - -GeoPhysVol* MuonShieldFactory::buildStation(const StationData& station) { - const GeoMaterial* air = m_materials.requireMaterial("Air"); - const GeoMaterial* iron = m_materials.requireMaterial("Iron"); - - // Air container that spans all 8 pieces of this station - auto* stationBox = - new GeoBox(station.containerHalfX, station.containerHalfY, station.containerHalfZ); - std::string containerName = "/SHiP/muon_shield/" + std::string(station.name); - auto* stationLog = new GeoLogVol(containerName, stationBox, air); - auto* stationPhys = new GeoPhysVol(stationLog); - - // Place 8 Iron bounding-box approximations - for (const PieceData& piece : station.pieces) { - auto* pieceBox = new GeoBox(piece.halfX, piece.halfY, piece.halfZ); - std::string pieceName = "/SHiP/muon_shield/" + std::string(station.name) + "/" + piece.name; - auto* pieceLog = new GeoLogVol(pieceName, pieceBox, iron); - auto* piecePhys = new GeoPhysVol(pieceLog); - - GeoTrf::Transform3D trf = GeoTrf::Translate3D(piece.centX, piece.centY, 0.0); - stationPhys->add(new GeoNameTag(pieceName)); - stationPhys->add(new GeoIdentifierTag(static_cast(&piece - &station.pieces[0]))); - stationPhys->add(new GeoTransform(trf)); - stationPhys->add(piecePhys); - } +// ── constructor / accessors ────────────────────────────────────────────────── + +MuonShieldFactory::MuonShieldFactory(SHiPMaterials& materials, std::string configPath) + : m_materials(materials), m_configPath(std::move(configPath)) {} + +void MuonShieldFactory::embedDaughter(GeoPhysVol* daughter, double worldCentreZ_mm, + const std::string& name) { + m_daughters.push_back({daughter, worldCentreZ_mm, name}); +} - return stationPhys; +void MuonShieldFactory::reserveSpace(const std::array& worldCentre_mm, + const std::array& size_mm, + const std::array& rotation_deg) { + m_reservations.push_back({worldCentre_mm, size_mm, rotation_deg}); } +std::string MuonShieldFactory::resolvedConfigPath() const { + return resolveConfigPath(m_configPath, MS_TOML_DEFAULT_PATH, MS_TOML_INSTALL_PATH); +} + +// ── build ──────────────────────────────────────────────────────────────────── + GeoPhysVol* MuonShieldFactory::build() { - const GeoMaterial* air = m_materials.requireMaterial("Air"); - - // Overall MuonShieldArea container (Air) - auto* areaBox = new GeoBox(s_areaHalfX, s_areaHalfY, s_areaHalfZ); - auto* areaLog = new GeoLogVol("/SHiP/muon_shield", areaBox, air); - auto* areaPhys = new GeoPhysVol(areaLog); - - // Build and place 6 stations - for (const StationData& station : k_stations) { - GeoPhysVol* stationPhys = buildStation(station); - std::string stationName = "/SHiP/muon_shield/" + std::string(station.name); - GeoTrf::Transform3D trf = GeoTrf::Translate3D(0.0, 0.0, station.stationZ); - areaPhys->add(new GeoNameTag(stationName)); - areaPhys->add(new GeoIdentifierTag(static_cast(&station - &k_stations[0]))); - areaPhys->add(new GeoTransform(trf)); - areaPhys->add(stationPhys); + const MuonShieldConfig cfg = readMuonShieldConfig(resolvedConfigPath()); + + GeoMaterial* air = m_materials.requireMaterial("Air"); + GeoMaterial* absorber = m_materials.requireMaterial(cfg.block_material); + + // Envelope-derived dimensions (mm). + const double halfZ = 0.5 * cfg.envelopeLengthZ_mm(); + m_centreZ_mm = cfg.envelopeCentreZ_mm(); + + // Air container spanning the full muon-shield envelope, centred on its own + // origin (SHiPGeometryBuilder places it at centreZ_mm()). + auto* containerBox = new GeoBox(cfg.envelope_half_x_mm, cfg.envelope_half_y_mm, halfZ); + auto* containerLog = new GeoLogVol("/SHiP/muon_shield", containerBox, air); + auto* containerPhys = new GeoPhysVol(containerLog); + + // Place each iron block. A block is anchored by its upstream (−z) face at + // its `start`, extends downstream by size[2], is rotated about `start` + // (extrinsic X→Y→Z), and may taper symmetrically into a GeoTrd. + int childId = 0; + for (const auto& b : cfg.blocks) { + const double halfX = 0.5 * b.size[0]; + const double halfY = 0.5 * b.size[1]; + const double halfLenZ = 0.5 * b.size[2]; + + const double farHalfX = halfX + b.size[2] * std::tan(b.taper_deg[0] * kDegToRad); + const double farHalfY = halfY + b.size[2] * std::tan(b.taper_deg[1] * kDegToRad); + + // Placement transform: shift the solid (centred at its origin) so its + // upstream face sits at `start`, rotate about `start` (extrinsic + // X→Y→Z), then translate into place. Container is at m_centreZ_mm in Z. + const double anchorX = b.start[0]; + const double anchorY = b.start[1]; + const double anchorZ = b.start[2] - m_centreZ_mm; + const GeoTrf::Transform3D rotation = GeoTrf::RotateZ3D(b.rotation_deg[2] * kDegToRad) * + GeoTrf::RotateY3D(b.rotation_deg[1] * kDegToRad) * + GeoTrf::RotateX3D(b.rotation_deg[0] * kDegToRad); + const GeoTrf::Transform3D trf = GeoTrf::Translate3D(anchorX, anchorY, anchorZ) * rotation * + GeoTrf::Translate3D(0.0, 0.0, halfLenZ); + + GeoShape* shape = nullptr; + if (b.taper_deg[0] == 0.0 && b.taper_deg[1] == 0.0) + shape = new GeoBox(halfX, halfY, halfLenZ); + else + shape = new GeoTrd(halfX, farHalfX, halfY, farHalfY, halfLenZ); + + // Carve out any reserved box that intersects this magnet (A - B). The + // reserved box, given in world coords, is expressed in the block's own + // (shape) frame before subtracting. + const GeoTrf::Vector3D blockCentre = trf.translation(); + const double blockMaxHX = std::max(halfX, farHalfX); + const double blockMaxHY = std::max(halfY, farHalfY); + for (const auto& r : m_reservations) { + const double rcx = r.centre_mm[0]; + const double rcy = r.centre_mm[1]; + const double rcz = r.centre_mm[2] - m_centreZ_mm; + const double rhx = 0.5 * r.size_mm[0]; + const double rhy = 0.5 * r.size_mm[1]; + const double rhz = 0.5 * r.size_mm[2]; + const bool intersects = std::abs(blockCentre.x() - rcx) < blockMaxHX + rhx && + std::abs(blockCentre.y() - rcy) < blockMaxHY + rhy && + std::abs(blockCentre.z() - rcz) < halfLenZ + rhz; + if (!intersects) + continue; + const GeoTrf::Transform3D resRot = GeoTrf::RotateZ3D(r.rotation_deg[2] * kDegToRad) * + GeoTrf::RotateY3D(r.rotation_deg[1] * kDegToRad) * + GeoTrf::RotateX3D(r.rotation_deg[0] * kDegToRad); + const GeoTrf::Transform3D resInBlock = + trf.inverse() * GeoTrf::Translate3D(rcx, rcy, rcz) * resRot; + auto* voidBox = new GeoBox(rhx, rhy, rhz); + shape = new GeoShapeSubtraction(shape, new GeoShapeShift(voidBox, resInBlock)); + } + + const std::string blockName = "/SHiP/muon_shield/block_" + std::to_string(childId); + auto* blockLog = new GeoLogVol(blockName, shape, absorber); + auto* blockPhys = new GeoPhysVol(blockLog); + + containerPhys->add(new GeoNameTag(blockName)); + containerPhys->add(new GeoIdentifierTag(childId++)); + containerPhys->add(new GeoTransform(trf)); + containerPhys->add(blockPhys); + } + + // Place the embedded daughters inside the container, keeping them as + // independent, named subsystems. (Callers must not list iron blocks over a + // daughter's footprint — the block list is the sole source of iron.) + for (const auto& d : m_daughters) { + const double localZ = d.worldCentreZ_mm - m_centreZ_mm; + if (std::abs(localZ) > halfZ) + throw std::runtime_error("MuonShieldFactory: embedded daughter '" + d.name + + "' lies outside the shield envelope in Z"); + containerPhys->add(new GeoNameTag(d.name)); + containerPhys->add(new GeoIdentifierTag(childId++)); + containerPhys->add(new GeoTransform(GeoTrf::Translate3D(0.0, 0.0, localZ))); + containerPhys->add(d.volume); } - return areaPhys; + return containerPhys; } } // namespace SHiPGeometry diff --git a/subsystems/MuonShield/test_muonshield.cpp b/subsystems/MuonShield/test_muonshield.cpp index 71dbdee..0345234 100644 --- a/subsystems/MuonShield/test_muonshield.cpp +++ b/subsystems/MuonShield/test_muonshield.cpp @@ -1,26 +1,183 @@ // SPDX-License-Identifier: LGPL-3.0-or-later // Copyright (C) CERN for the benefit of the SHiP Collaboration +#include "MuonShield/MuonShieldConfig.h" #include "MuonShield/MuonShieldFactory.h" #include "SHiPGeometry/SHiPMaterials.h" #include #include #include +#include +#include #include +#include +#include +#include +#include +#include +using SHiPGeometry::MuonShieldConfig; +using SHiPGeometry::MuonShieldFactory; +using SHiPGeometry::readMuonShieldConfig; using SHiPGeometry::SHiPMaterials; -// MuonShieldArea container halfX ≤ 2100 (CSV WARM max half-width), -// halfY ≤ 2300 (CSV WARM max half-height) -TEST_CASE("MuonShieldWithinEnvelope", "[muonshield]") { +namespace { +// Write a temp MS.toml with the given body and return its path. +std::string writeTempToml(const std::string& name, const std::string& body) { + std::ofstream out(name); + out << body; + out.close(); + return name; +} +} // namespace + +// Default MS.toml: the 7 FairShip TRY_2026 magnets (solid-block approximation) +// inside an auto-sized envelope (4.54–32.08 m, 1760 × 1320 mm half-sizes). +TEST_CASE("MuonShieldBuilds", "[muonshield]") { SHiPMaterials materials; - SHiPGeometry::MuonShieldFactory factory(materials); + MuonShieldFactory factory(materials); GeoPhysVol* ms = factory.build(); REQUIRE(ms != nullptr); + auto* box = dynamic_cast(ms->getLogVol()->getShape()); REQUIRE(box != nullptr); - CHECK(box->getXHalfLength() <= 2100.0); - CHECK(box->getYHalfLength() <= 2300.0); + CHECK_THAT(box->getXHalfLength(), Catch::Matchers::WithinAbs(1760.0, 1e-6)); + CHECK_THAT(box->getYHalfLength(), Catch::Matchers::WithinAbs(1320.0, 1e-6)); + CHECK_THAT(box->getZHalfLength(), Catch::Matchers::WithinAbs(13770.0, 1e-6)); + + // Envelope centre = (4.54 + 32.08)/2 m = 18.31 m. + CHECK_THAT(factory.centreZ_mm(), Catch::Matchers::WithinAbs(18310.0, 1e-6)); +} + +TEST_CASE("MuonShieldDefaultLayout", "[muonshield]") { + SHiPMaterials materials; + MuonShieldFactory factory(materials); + GeoPhysVol* ms = factory.build(); + REQUIRE(ms != nullptr); + // 7 solid magnets (the SND cavity is carved by reserveSpace, not here). + REQUIRE(ms->getNChildVols() == 7u); // NOLINT(readability/check) + + // Magnet 1: straight box, upstream face at z = 4.59 m, 2720 × 1600 × 3000 mm. + auto* block0 = dynamic_cast(ms->getChildVol(0)->getLogVol()->getShape()); + REQUIRE(block0 != nullptr); + CHECK_THAT(block0->getXHalfLength(), Catch::Matchers::WithinAbs(1360.0, 1e-6)); + CHECK_THAT(block0->getYHalfLength(), Catch::Matchers::WithinAbs(800.0, 1e-6)); + CHECK_THAT(block0->getZHalfLength(), Catch::Matchers::WithinAbs(1500.0, 1e-6)); + + // Anchor (upstream face) at world 4590 → local (4590-18310); centre + halfLen. + const double blockCentreLocalZ = (4590.0 - 18310.0) + 1500.0; // = -12220 + CHECK_THAT(ms->getXToChildVol(0).translation().z(), + Catch::Matchers::WithinAbs(blockCentreLocalZ, 1e-3)); + + // Magnet 6 (index 5) widens in Y → GeoTrd. + CHECK(dynamic_cast(ms->getChildVol(5)->getLogVol()->getShape()) != nullptr); +} + +TEST_CASE("MuonShieldReservationCarvesIron", "[muonshield]") { + // A reserved box (the SND envelope: 800 × 800 × 5100 mm at z = 28.95 m) is + // subtracted (A - B) from every magnet it intersects, leaving upstream + // magnets untouched. + SHiPMaterials materials; + MuonShieldFactory factory(materials); + factory.reserveSpace({0.0, 0.0, 28950.0}, {800.0, 800.0, 5100.0}); + GeoPhysVol* ms = factory.build(); + REQUIRE(ms != nullptr); + REQUIRE(ms->getNChildVols() == 7u); + + // Magnet 1 (upstream, far from the SND) is untouched → still a plain box. + CHECK(dynamic_cast(ms->getChildVol(0)->getLogVol()->getShape()) != nullptr); + // The last magnet intersects the reservation → its shape is a subtraction. + CHECK(dynamic_cast(ms->getChildVol(6)->getLogVol()->getShape()) != + nullptr); +} + +TEST_CASE("MuonShieldTaperMakesTrd", "[muonshield]") { + // A block with an X taper becomes a GeoTrd that widens downstream. + const std::string path = writeTempToml( + "MS_taper.toml", + "[[block]]\nstart = [0,0,12000]\nsize = [3000,2000,2000]\ntaper = [5.0, 0.0]\n"); + SHiPMaterials materials; + MuonShieldFactory factory(materials, path); + GeoPhysVol* ms = factory.build(); + REQUIRE(ms != nullptr); + REQUIRE(ms->getNChildVols() == 1u); + + auto* trd = dynamic_cast(ms->getChildVol(0)->getLogVol()->getShape()); + REQUIRE(trd != nullptr); + CHECK_THAT(trd->getXHalfLength1(), Catch::Matchers::WithinAbs(1500.0, 1e-6)); + const double farHalfX = 1500.0 + 2000.0 * std::tan(5.0 * 3.14159265358979323846 / 180.0); + CHECK_THAT(trd->getXHalfLength2(), Catch::Matchers::WithinAbs(farHalfX, 1e-3)); + CHECK_THAT(trd->getYHalfLength1(), Catch::Matchers::WithinAbs(1000.0, 1e-6)); + CHECK_THAT(trd->getYHalfLength2(), Catch::Matchers::WithinAbs(1000.0, 1e-6)); +} + +TEST_CASE("MuonShieldParsesRotation", "[muonshield]") { + // The optional `rotation` field is parsed into rotation_deg (degrees, about + // x, y, z). The factory applies it as an extrinsic X→Y→Z rotation about the + // block's start; asserting the built orientation would require reading the + // child transform back, which is left to OverlapCheck / a future test. + const std::string path = writeTempToml( + "MS_rot.toml", + "[[block]]\nstart = [0,0,12000]\nsize = [200,200,400]\nrotation = [10.0, 20.0, 30.0]\n"); + const MuonShieldConfig cfg = readMuonShieldConfig(path); + REQUIRE(cfg.blocks.size() == 1u); + CHECK_THAT(cfg.blocks[0].rotation_deg[0], Catch::Matchers::WithinAbs(10.0, 1e-9)); + CHECK_THAT(cfg.blocks[0].rotation_deg[1], Catch::Matchers::WithinAbs(20.0, 1e-9)); + CHECK_THAT(cfg.blocks[0].rotation_deg[2], Catch::Matchers::WithinAbs(30.0, 1e-9)); + + // A rotated block still builds and yields a valid shape. + SHiPMaterials materials; + MuonShieldFactory factory(materials, path); + GeoPhysVol* ms = factory.build(); + REQUIRE(ms != nullptr); + REQUIRE(ms->getNChildVols() == 1u); + CHECK(dynamic_cast(ms->getChildVol(0)->getLogVol()->getShape()) != nullptr); +} + +TEST_CASE("MuonShieldEmbedsDaughter", "[muonshield]") { + SHiPMaterials materials; + auto* dBox = new GeoBox(100.0, 100.0, 500.0); + auto* dLog = new GeoLogVol("/SHiP/dummy", dBox, materials.requireMaterial("Air")); + auto* dPhys = new GeoPhysVol(dLog); + + MuonShieldFactory factory(materials); // default 7 solid magnets + factory.embedDaughter(dPhys, 28.95 * 1000.0, "/SHiP/dummy"); + GeoPhysVol* ms = factory.build(); + REQUIRE(ms != nullptr); + // 7 iron magnets + the embedded daughter. + CHECK(ms->getNChildVols() == 8u); // NOLINT(readability/check) + + bool found = false; + for (unsigned i = 0; i < ms->getNChildVols(); ++i) + if (ms->getChildVol(i)->getLogVol()->getName() == "/SHiP/dummy") + found = true; + CHECK(found); +} + +TEST_CASE("MuonShieldRejectsNonPositiveSize", "[muonshield]") { + const std::string path = writeTempToml( + "MS_badsize.toml", "[[block]]\nstart = [0,0,12000]\nsize = [-3000,2000,2000]\n"); + CHECK_THROWS_AS(readMuonShieldConfig(path), std::runtime_error); +} + +TEST_CASE("MuonShieldRejectsCollapsingTaper", "[muonshield]") { + const std::string path = writeTempToml( + "MS_badtaper.toml", + "[[block]]\nstart = [0,0,12000]\nsize = [3000,2000,2000]\ntaper = [-45.0, 0.0]\n"); + CHECK_THROWS_AS(readMuonShieldConfig(path), std::runtime_error); +} + +TEST_CASE("MuonShieldRejectsBlockOutsideEnvelope", "[muonshield]") { + const std::string path = writeTempToml( + "MS_outside.toml", "[[block]]\nstart = [0,0,40000]\nsize = [3000,2000,1000]\n"); + CHECK_THROWS_AS(readMuonShieldConfig(path), std::runtime_error); +} + +TEST_CASE("MuonShieldEmptyBlockList", "[muonshield]") { + // No blocks → container with no iron (valid; e.g. a placeholder). + const std::string path = writeTempToml("MS_empty.toml", "block_material = \"Iron\"\n"); + MuonShieldConfig cfg = readMuonShieldConfig(path); + CHECK(cfg.blocks.empty()); } diff --git a/tests/test_consistency.cpp b/tests/test_consistency.cpp index ad8130c..4497ddd 100644 --- a/tests/test_consistency.cpp +++ b/tests/test_consistency.cpp @@ -76,9 +76,10 @@ TEST_CASE("ConsistencyTest.ExpectedSubsystemCount", "[consistency]") { REQUIRE(world != nullptr); auto subsystems = collectSubsystems(world); - // 9 subsystems: target, muon_shield, neutrino_detector, upstream_tagger, - // decay_volume, trackers, magnet, timing_detector, calorimeter - CHECK(subsystems.size() == 9u); // NOLINT(readability/check) + // 8 subsystems: target, muon_shield, upstream_tagger, decay_volume, + // trackers, magnet, timing_detector, calorimeter. The neutrino detector is + // nested inside the muon-shield container, so it is not a direct world child. + CHECK(subsystems.size() == 8u); // NOLINT(readability/check) } TEST_CASE("ConsistencyTest.SubsystemsGenerallyInZOrder", "[consistency]") { @@ -110,13 +111,9 @@ TEST_CASE("ConsistencyTest.NoUnexpectedZOverlaps", "[consistency]") { // The trackers container intentionally spans across the magnet // (stations 1-2 before, stations 3-4 after), so that pair is allowed to overlap. - // The SND sits inside the downstream end of the muon-shield region, so the - // muon_shield / neutrino_detector pair is also an intentional overlap. auto isAllowedOverlap = [](const std::string& a, const std::string& b) { return (a == "/SHiP/trackers" && b == "/SHiP/magnet") || - (a == "/SHiP/magnet" && b == "/SHiP/trackers") || - (a == "/SHiP/muon_shield" && b == "/SHiP/neutrino_detector") || - (a == "/SHiP/neutrino_detector" && b == "/SHiP/muon_shield"); + (a == "/SHiP/magnet" && b == "/SHiP/trackers"); }; // Sort by Z centre @@ -156,8 +153,10 @@ TEST_CASE("ConsistencyTest.PositionsSanity", "[consistency]") { // Centres as placed in SHiPGeometryBuilder::build() std::vector expected = { {"/SHiP/target", 432.5, 500.0}, - {"/SHiP/muon_shield", 16763.3, 500.0}, - {"/SHiP/neutrino_detector", 28950.0, 500.0}, + {"/SHiP/muon_shield", 18310.0, 500.0}, + // The neutrino detector is nested inside the muon shield (see + // MuonShieldFactory::embedDaughter), not a direct world child, so it is + // not among collectSubsystems(world). {"/SHiP/upstream_tagger", 32720.0, 500.0}, {"/SHiP/decay_volume", 58120.0, 500.0}, {"/SHiP/trackers", 89570.0, 500.0}, From 42b93a98d5c83c0e653eafe5520f7301dec4078c Mon Sep 17 00:00:00 2001 From: Matei Climescu Date: Sat, 11 Jul 2026 12:11:45 +0200 Subject: [PATCH 2/3] feat(ms): Reserve SND space via Boolean subtraction Define the muon shield and the neutrino detector independently and carve the SND cavity out of the iron at build time. MuonShieldFactory::reserveSpace() registers a box; build() subtracts it (GeoShapeSubtraction) from every magnet it intersects and leaves the rest untouched, so which magnets are carved follows from the reservation rather than a hard-coded list. The reservation envelope is declared in NeutrinoDetector/SD.toml, parsed by SNDEnvelope via the shared resolveConfigPath helper, and wired in by SHiPGeometryBuilder before the shield is built. --- subsystems/NeutrinoDetector/CMakeLists.txt | 41 ++++++++ subsystems/NeutrinoDetector/README.md | 16 ++-- subsystems/NeutrinoDetector/SD.toml | 16 ++++ .../include/NeutrinoDetector/SNDEnvelope.h | 39 ++++++++ .../NeutrinoDetector/src/SNDEnvelope.cpp | 95 +++++++++++++++++++ 5 files changed, 201 insertions(+), 6 deletions(-) create mode 100644 subsystems/NeutrinoDetector/SD.toml create mode 100644 subsystems/NeutrinoDetector/include/NeutrinoDetector/SNDEnvelope.h create mode 100644 subsystems/NeutrinoDetector/src/SNDEnvelope.cpp diff --git a/subsystems/NeutrinoDetector/CMakeLists.txt b/subsystems/NeutrinoDetector/CMakeLists.txt index 95ce231..a1b4c8a 100644 --- a/subsystems/NeutrinoDetector/CMakeLists.txt +++ b/subsystems/NeutrinoDetector/CMakeLists.txt @@ -1,4 +1,45 @@ # SPDX-License-Identifier: LGPL-3.0-or-later # Copyright (C) CERN for the benefit of the SHiP Collaboration +# Shared recipe: library (globs src/*.cpp, incl. SNDEnvelope.cpp), include dirs, +# kernel link, Catch test. ship_add_subsystem(NeutrinoDetector) + +# --- Extras: NeutrinoDetector parses SD.toml (the SND reservation envelope) via +# SNDEnvelope (toml++), the same way MuonShield/Calorimeter parse their configs. + +find_package(tomlplusplus REQUIRED) + +# toml++ is a build-time-only dependency (used only inside SNDEnvelope.cpp, never +# in a public header). BUILD_INTERFACE keeps it out of the install export set. +target_link_libraries( + NeutrinoDetector + PRIVATE $ +) + +# Stage SD.toml next to the tests and into the top-level build dir (for +# test_builder / OverlapCheck), so it is found when running from either. +configure_file( + ${CMAKE_CURRENT_SOURCE_DIR}/SD.toml + ${CMAKE_CURRENT_BINARY_DIR}/SD.toml + COPYONLY +) +configure_file( + ${CMAKE_CURRENT_SOURCE_DIR}/SD.toml + ${CMAKE_BINARY_DIR}/SD.toml + COPYONLY +) + +install( + FILES ${CMAKE_CURRENT_SOURCE_DIR}/SD.toml + DESTINATION ${CMAKE_INSTALL_DATADIR}/SHiPGeometry +) + +# Absolute fallbacks for resolveConfigPath(), so the parser finds SD.toml when +# the CWD does not have it (mirrors the calorimeter's CALO_TOML_* defines). +target_compile_definitions( + NeutrinoDetector + PRIVATE + SD_TOML_DEFAULT_PATH="${CMAKE_CURRENT_SOURCE_DIR}/SD.toml" + SD_TOML_INSTALL_PATH="${CMAKE_INSTALL_PREFIX}/${CMAKE_INSTALL_DATADIR}/SHiPGeometry/SD.toml" +) diff --git a/subsystems/NeutrinoDetector/README.md b/subsystems/NeutrinoDetector/README.md index 99ca33c..1d1d0ef 100644 --- a/subsystems/NeutrinoDetector/README.md +++ b/subsystems/NeutrinoDetector/README.md @@ -35,7 +35,7 @@ bounded for navigation. ## Geometry tree ```text -/SHiP/neutrino_detector (Air, 640 x 640 x 4100 mm box) +/SHiP/neutrino_detector (Air, 800 x 800 x 5100 mm box) |- veto/P{0..2}_B{0..6} 21 PVT bars (3 planes x 7) |- target/L{0..119}/{W,Si_X,Si_Y} 120 W + 120 Si-X + 120 Si-Y \- hcal/S{0..2}_L{0..13}/ @@ -48,11 +48,15 @@ bounded for navigation. Direct container children: 21 + 360 + 4438 = 4819. Total volumes including the individual fibres: 340819 (~336k fibres). -Position in world: z = 28 950 mm (centre of the 26.40-31.50 m WARM SND slot). -The SND is part of the muon shield, so its box container sits within the -downstream end of the muon-shield region; the resulting envelope overlap is -declared in `tests/test_consistency.cpp`, analogous to the trackers/magnet -overlap. +Position: z = 28 950 mm (centre of the 26.40-31.50 m WARM SND slot). The SND is +part of the muon shield: in volume terms its container is a **daughter of the +muon-shield container**, nested by `SHiPGeometryBuilder` via +`MuonShieldFactory::embedDaughter`, which only places the pre-built daughter in +the shield's volume tree. The matching cavity in the iron is carved separately +by `MuonShieldFactory::reserveSpace` (a Boolean subtraction of the SND +reservation envelope from `SD.toml`), so the two do not overlap. It remains an +independent subsystem with its own factory, config, and +`/SHiP/neutrino_detector` naming — only its place in the volume tree changes. ## Materials diff --git a/subsystems/NeutrinoDetector/SD.toml b/subsystems/NeutrinoDetector/SD.toml new file mode 100644 index 0000000..ae6ba9a --- /dev/null +++ b/subsystems/NeutrinoDetector/SD.toml @@ -0,0 +1,16 @@ +# SPDX-License-Identifier: LGPL-3.0-or-later +# Copyright (C) CERN for the benefit of the SHiP Collaboration +# +# Scattering & Neutrino Detector (SND) reservation envelope. +# +# This declares the box that is subtracted from the muon-shield iron (A - B) to +# make room for the SND, so the two subsystems can be defined independently and +# merged at build time. For now this is only the reservation envelope; the SND's +# internal geometry will move into this file in a later, dedicated PR. +# +# Values are in world/beamline coordinates (target-front-face origin), mm. +# Defaults match the current SND container (800 x 800 x 5100 mm at z = 28.95 m). + +centre = [0.0, 0.0, 28950.0] # world centre of the reservation box +size = [800.0, 800.0, 5100.0] # full x, y, z +# rotation = [0.0, 0.0, 0.0] # optional, degrees, extrinsic X -> Y -> Z about centre diff --git a/subsystems/NeutrinoDetector/include/NeutrinoDetector/SNDEnvelope.h b/subsystems/NeutrinoDetector/include/NeutrinoDetector/SNDEnvelope.h new file mode 100644 index 0000000..cb1f26b --- /dev/null +++ b/subsystems/NeutrinoDetector/include/NeutrinoDetector/SNDEnvelope.h @@ -0,0 +1,39 @@ +// SPDX-License-Identifier: LGPL-3.0-or-later +// Copyright (C) CERN for the benefit of the SHiP Collaboration + +#pragma once + +#include +#include + +namespace SHiPGeometry { + +/** + * @brief The SND reservation envelope: a box subtracted from the muon-shield + * iron to make room for the neutrino detector. + * + * All values are in world/beamline coordinates (target-front-face origin). + * Populated by readSNDEnvelope() from SD.toml. + */ +struct SNDEnvelope { + std::array centre_mm = {0.0, 0.0, 28950.0}; ///< world centre (mm) + std::array size_mm = {800.0, 800.0, 5100.0}; ///< full x, y, z (mm) + std::array rotation_deg = {0.0, 0.0, 0.0}; ///< extrinsic X->Y->Z (deg) +}; + +/** + * @brief Parse an SD.toml file into an SNDEnvelope. + * + * The path is resolved with the shared resolveConfigPath() fallback chain + * (CWD, then the source tree, then the installed data dir), so callers can + * simply use the default. + * + * @throws std::runtime_error on malformed TOML, a wrong-length array, or a + * non-positive size. + */ +SNDEnvelope readSNDEnvelope(const std::string& path = "SD.toml"); + +/// Resolve the SD.toml path that readSNDEnvelope() would actually open. +std::string resolvedSNDConfigPath(const std::string& path = "SD.toml"); + +} // namespace SHiPGeometry diff --git a/subsystems/NeutrinoDetector/src/SNDEnvelope.cpp b/subsystems/NeutrinoDetector/src/SNDEnvelope.cpp new file mode 100644 index 0000000..f2c58d9 --- /dev/null +++ b/subsystems/NeutrinoDetector/src/SNDEnvelope.cpp @@ -0,0 +1,95 @@ +// SPDX-License-Identifier: LGPL-3.0-or-later +// Copyright (C) CERN for the benefit of the SHiP Collaboration + +#include "NeutrinoDetector/SNDEnvelope.h" + +#include "SHiPGeometry/ConfigPath.h" + +#include +#include +#include +#include +#include +#include +#include + +// Absolute fallback paths baked in by CMake so out-of-source builds always find +// SD.toml even when the CWD doesn't contain a copy of it. +#ifndef SD_TOML_DEFAULT_PATH +#define SD_TOML_DEFAULT_PATH "SD.toml" +#endif +#ifndef SD_TOML_INSTALL_PATH +#define SD_TOML_INSTALL_PATH "" +#endif + +namespace SHiPGeometry { + +namespace { + +double toNumber(const toml::node* node, const std::string& key, const std::string& path) { + if (node) { + if (auto d = node->value()) + return *d; + if (auto i = node->value()) + return static_cast(*i); + } + throw std::runtime_error("SNDEnvelope: '" + key + "' values must be numbers in " + path); +} + +std::array readVec3(const toml::table& table, const char* key, const std::string& path, + bool required, const std::array& fallback) { + auto node = table[key]; + if (!node) { + if (required) + throw std::runtime_error("SNDEnvelope: missing required '" + std::string(key) + + "' in " + path); + return fallback; + } + const toml::array* arr = node.as_array(); + if (!arr || arr->size() != 3) + throw std::runtime_error("SNDEnvelope: '" + std::string(key) + + "' must be a 3-element array in " + path); + std::array out{}; + for (std::size_t i = 0; i < 3; ++i) + out[i] = toNumber(arr->get(i), key, path); + return out; +} + +} // namespace + +std::string resolvedSNDConfigPath(const std::string& path) { + return resolveConfigPath(path, SD_TOML_DEFAULT_PATH, SD_TOML_INSTALL_PATH); +} + +SNDEnvelope readSNDEnvelope(const std::string& rawPath) { + const std::string path = resolvedSNDConfigPath(rawPath); + toml::table table; + try { + table = toml::parse_file(path); + } catch (const toml::parse_error& e) { + throw std::runtime_error("SNDEnvelope: failed to parse " + path + ": " + + std::string(e.description())); + } + + SNDEnvelope env; + + // Warn about unknown top-level keys (typos / stale fields), like + // readMuonShieldConfig. Only centre, size and rotation are recognised. + for (const auto& [k, _] : table) { + const std::string_view key{k}; + if (key != "centre" && key != "size" && key != "rotation") + std::cerr << "SNDEnvelope: warning: unknown key '" << k << "' in " << path + << " (ignored)\n"; + } + + env.centre_mm = readVec3(table, "centre", path, true, env.centre_mm); + env.size_mm = readVec3(table, "size", path, true, env.size_mm); + env.rotation_deg = readVec3(table, "rotation", path, false, env.rotation_deg); + + if (env.size_mm[0] <= 0.0 || env.size_mm[1] <= 0.0 || env.size_mm[2] <= 0.0) + throw std::runtime_error("SNDEnvelope: size must be positive in " + path); + + return env; +} + +} // namespace SHiPGeometry From 54490ea9dda1ea4e30ee060832badadcebdc45e1 Mon Sep 17 00:00:00 2001 From: Matei Climescu Date: Fri, 17 Jul 2026 14:17:49 +0200 Subject: [PATCH 3/3] fix: updated MS to fit fairship update --- subsystems/MuonShield/MS.toml | 6 +++--- subsystems/MuonShield/test_muonshield.cpp | 2 +- 2 files changed, 4 insertions(+), 4 deletions(-) diff --git a/subsystems/MuonShield/MS.toml b/subsystems/MuonShield/MS.toml index 73d8bb9..d853921 100644 --- a/subsystems/MuonShield/MS.toml +++ b/subsystems/MuonShield/MS.toml @@ -46,11 +46,11 @@ size = [3184.0, 1200.0, 4500.0] start = [0.0, 0.0, 21710.0] size = [2680.0, 1040.0, 2620.0] -# magnet 6 (widens in Y toward downstream) +# magnet 6 [[block]] start = [0.0, 0.0, 24480.0] -size = [3000.0, 1060.0, 3360.0] -taper = [0.0000, 1.7047] +size = [3000.0, 1260.0, 3360.0] + # magnet 7 [[block]] diff --git a/subsystems/MuonShield/test_muonshield.cpp b/subsystems/MuonShield/test_muonshield.cpp index 0345234..215786a 100644 --- a/subsystems/MuonShield/test_muonshield.cpp +++ b/subsystems/MuonShield/test_muonshield.cpp @@ -72,7 +72,7 @@ TEST_CASE("MuonShieldDefaultLayout", "[muonshield]") { Catch::Matchers::WithinAbs(blockCentreLocalZ, 1e-3)); // Magnet 6 (index 5) widens in Y → GeoTrd. - CHECK(dynamic_cast(ms->getChildVol(5)->getLogVol()->getShape()) != nullptr); + CHECK(dynamic_cast(ms->getChildVol(5)->getLogVol()->getShape()) != nullptr); } TEST_CASE("MuonShieldReservationCarvesIron", "[muonshield]") {