feat(ms)!: Add configurable iron-block muon shield (FairShip TRY_2026)

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;

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 $<BUILD_INTERFACE:tomlplusplus::tomlplusplus>
+)
+
+# 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
+)

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
+[[block]]
+start = [0.0, 0.0, 24480.0]
+size  = [3000.0, 1260.0, 3360.0]
+
+
+# magnet 7
+[[block]]
+start = [0.0, 0.0, 28030.0]
+size  = [3420.0, 2540.0, 4000.0]

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.

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 <array>
+#include <string>
+#include <vector>
+
+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<double, 3> start = {0.0, 0.0, 0.0};         ///< upstream-face centre (mm)
+    std::array<double, 3> size = {0.0, 0.0, 0.0};          ///< near-face sx, sy and length sz (mm)
+    std::array<double, 3> rotation_deg = {0.0, 0.0, 0.0};  ///< rotation about x, y, z (deg)
+    std::array<double, 2> 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<MuonShieldBlock> 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