[ESSDIFFRACTION] feat: add provider for dspacing range of instrument

packages/essdiffraction/docs/user-guide/beer/beer_modulation_mcstas.ipynb

@@ -713,7 +713,7 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython3",
-   "version": "3.11.13"
+   "version": "3.12.13"
   }
  },
  "nbformat": 4,

packages/essdiffraction/docs/user-guide/beer/beer_powder_mcstas_analytical.ipynb

@@ -80,7 +80,12 @@
    "outputs": [],
    "source": [
     "workflow[DetectorBank] = DetectorBank.north\n",
-    "workflow[DspacingBins] = sc.linspace(\"dspacing\", 0.6, 2.0, 701, unit=\"angstrom\")"
+    "\n",
+    "# Select the number of dspacing bins.\n",
+    "# By default the range of `DspacingBins` (the binning in dspacing) is\n",
+    "# selected to be the dspacing range covered by the instrument.\n",
+    "workflow[DspacingNBins] = 1000\n",
+    "workflow.compute(DspacingBins)"
    ]
   },
   {
@@ -198,8 +203,6 @@
     "    title=\"Intensity over dspacing\",\n",
     "    xlabel=r\"$d$-spacing [$\\AA$]\",\n",
     "    ylabel=f\"Intensity [{histogram.unit}]\",\n",
-    "    xmin=sc.scalar(0.6, unit=\"angstrom\"),\n",
-    "    xmax=sc.scalar(2.0, unit=\"angstrom\"),\n",
     ")\n",
     "\n",
     "ymin, ymax = fig.ax.get_ylim()\n",
@@ -308,8 +311,6 @@
     "    },\n",
     "    xlabel=r\"$d$-spacing [$\\AA$]\",\n",
     "    ylabel=f\"Intensity [{histogram.unit}]\",\n",
-    "    xmin=sc.scalar(0.6, unit=\"angstrom\"),\n",
-    "    xmax=sc.scalar(2.0, unit=\"angstrom\"),\n",
     ")\n",
     "\n",
     "ymin, ymax = fig.ax.get_ylim()\n",
@@ -337,7 +338,7 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython3",
-   "version": "3.11.13"
+   "version": "3.12.13"
   }
  },
  "nbformat": 4,

packages/essdiffraction/src/ess/beer/mcstas/__init__.py

@@ -11,6 +11,8 @@
 from .load import (
     chopper_mode_from_mcstas_mode,
     load_beer_mcstas,
+    load_beer_mcstas_geometry,
+    load_beer_mcstas_geometry_provider,
     load_beer_mcstas_monitor,
     load_beer_mcstas_monitor_provider,
     load_beer_mcstas_provider,
@@ -26,6 +28,8 @@
     'PulseShapingMode',
     'chopper_mode_from_mcstas_mode',
     'load_beer_mcstas',
+    'load_beer_mcstas_geometry',
+    'load_beer_mcstas_geometry_provider',
     'load_beer_mcstas_monitor',
     'load_beer_mcstas_monitor_provider',
     'load_beer_mcstas_provider',

packages/essdiffraction/src/ess/beer/mcstas/load.py

@@ -7,7 +7,7 @@
 import mcstastox
 import scipp as sc
 import scippnexus as snx
-from ess.powder.types import CaveMonitor, RunType, WavelengthMonitor
+from ess.powder.types import CaveMonitor, EmptyDetector, RunType, WavelengthMonitor
 
 from ess.reduce.nexus.types import DiskChoppers, Position
 from ess.reduce.unwrap.types import DetectorLtotal
@@ -24,6 +24,8 @@
 __all__ = [
     'chopper_mode_from_mcstas_mode',
     'load_beer_mcstas',
+    'load_beer_mcstas_geometry',
+    'load_beer_mcstas_geometry_provider',
     'load_beer_mcstas_monitor',
     'load_beer_mcstas_monitor_provider',
     'load_beer_mcstas_provider',
@@ -62,6 +64,25 @@ def _detector_components(data: mcstastox.Read, bank: DetectorBank) -> list[str]:
     )
 
 
+def _load_position_table(data: mcstastox.Read, components: list[str]) -> sc.DataArray:
+    """Load pixel positions for selected detector components."""
+    position_tables = []
+    for component in components:
+        positions = data.get_component_global(component)
+        begin, end = data.pixel_range[component]
+        position_tables.append(
+            sc.DataArray(
+                sc.vectors(dims=['pixel_id'], values=positions, unit='m'),
+                coords={
+                    'pixel_id': sc.arange(
+                        'pixel_id', int(begin), int(end) + 1, dtype='int32'
+                    )
+                },
+            )
+        )
+    return sc.sort(sc.concat(position_tables, dim='pixel_id'), key='pixel_id')
+
+
 def _load_events(data: mcstastox.Read, components: list[str]) -> sc.DataArray:
     """Load weighted events and pixel positions for selected components only.
 
@@ -86,30 +107,22 @@ def _load_events(data: mcstastox.Read, components: list[str]) -> sc.DataArray:
             ),
             't': sc.array(dims=['events'], values=event_data['t'], unit='s'),
         },
-    ).group('pixel_id')
-
-    position_tables = []
-    for component in components:
-        positions = data.get_component_global(component)
-        begin, end = data.pixel_range[component]
-        position_tables.append(
-            sc.DataArray(
-                sc.vectors(dims=['pixel_id'], values=positions, unit='m'),
-                coords={
-                    'pixel_id': sc.arange(
-                        'pixel_id', int(begin), int(end) + 1, dtype='int32'
-                    )
-                },
-            )
-        )
+    )
 
-    position_table = sc.sort(sc.concat(position_tables, dim='pixel_id'), key='pixel_id')
-    events.coords['position'] = sc.lookup(
-        position_table, dim='pixel_id', mode='previous'
-    )[events.coords['pixel_id']]
+    position_table = _load_position_table(data, components)
+    # Group by the complete list of detector pixels, rather than by the IDs that
+    # occur in the event table. This retains pixels with no events as empty bins and
+    # keeps detector geometry independent of the measured event population.
+    events = events.group(position_table.coords['pixel_id'])
+    events.coords['position'] = position_table.data
     return events
 
 
+def _load_mode(data: mcstastox.Read) -> str:
+    mode = data.file['entry1/simulation/Param/mode'][0]
+    return mode.decode() if isinstance(mode, bytes) else str(mode)
+
+
 def load_beer_mcstas(
     filename: str | Path,
     bank: DetectorBank,
@@ -132,8 +145,7 @@ def load_beer_mcstas(
 
     filename = Path(filename)
     with mcstastox.Read(filename.parent, filename.name) as data:
-        mode = data.file['entry1/simulation/Param/mode'][0]
-        mode = mode.decode() if isinstance(mode, bytes) else str(mode)
+        mode = _load_mode(data)
         events = _load_events(data, _detector_components(data, bank))
 
     events.coords['mode'] = sc.scalar(mode)
@@ -143,6 +155,40 @@ def load_beer_mcstas(
     return events
 
 
+def load_beer_mcstas_geometry(
+    filename: str | Path,
+    bank: DetectorBank,
+) -> sc.DataArray:
+    """Load detector geometry from a BEER McStas file without loading events."""
+    if not isinstance(bank, DetectorBank):
+        raise ValueError(
+            'bank must be ``DetectorBank.north``, '
+            '``DetectorBank.south``, or ``DetectorBank.both``'
+        )
+
+    if bank == DetectorBank.both:
+        return sc.concat(
+            [
+                load_beer_mcstas_geometry(filename, bank)
+                for bank in (DetectorBank.south, DetectorBank.north)
+            ],
+            dim='pixel_id',
+        )
+
+    filename = Path(filename)
+    with mcstastox.Read(filename.parent, filename.name) as data:
+        positions = _load_position_table(data, _detector_components(data, bank))
+
+    pixel_id = positions.coords['pixel_id']
+    return sc.DataArray(
+        pixel_id,
+        coords={
+            'pixel_id': pixel_id,
+            'position': positions.data,
+        },
+    )
+
+
 def _to_edges(centers: sc.Variable) -> sc.Variable:
     interior_edges = sc.midpoints(centers)
     return sc.concat(
@@ -201,6 +247,13 @@ def load_beer_mcstas_provider(
     return load_beer_mcstas(fname, bank)
 
 
+def load_beer_mcstas_geometry_provider(
+    fname: Filename[RunType], bank: DetectorBank
+) -> EmptyDetector[RunType]:
+    """Sciline provider for loading BEER McStas detector geometry."""
+    return EmptyDetector[RunType](load_beer_mcstas_geometry(fname, bank))
+
+
 def load_beer_mcstas_monitor_provider(
     fname: Filename[RunType],
 ) -> WavelengthMonitor[RunType, CaveMonitor]:
@@ -287,22 +340,34 @@ def mcstas_choppers(
     return simulation_choppers(mode, source_position)
 
 
+def _mcstas_choppers_from_file(
+    fname: Filename[RunType],
+    source_position: Position[snx.NXsource, RunType],
+) -> DiskChoppers[RunType]:
+    """Return BEER choppers without loading detector events."""
+    fname = Path(fname)
+    with mcstastox.Read(fname.parent, fname.name) as data:
+        mode = chopper_mode_from_mcstas_mode(_load_mode(data))
+    return simulation_choppers(mode, source_position)
+
+
 def mcstas_detector_ltotal(
-    da: RawDetector[RunType],
+    detector: RawDetector[RunType],
     source_position: Position[snx.NXsource, RunType],
     sample_position: Position[snx.NXsample, RunType],
 ) -> DetectorLtotal[RunType]:
     """Return moderator-to-detector flight path lengths for BEER McStas data."""
     source_to_sample = sc.norm(sample_position - source_position)
-    sample_to_detector = sc.norm(da.coords['position'] - sample_position)
+    sample_to_detector = sc.norm(detector.coords['position'] - sample_position)
     return source_to_sample + sample_to_detector
 
 
 mcstas_providers = (
     load_beer_mcstas_provider,
+    load_beer_mcstas_geometry_provider,
     load_beer_mcstas_monitor_provider,
     mcstas_source_position,
     mcstas_sample_position,
-    mcstas_choppers,
+    _mcstas_choppers_from_file,
     mcstas_detector_ltotal,
 )

packages/essdiffraction/src/ess/powder/__init__.py

@@ -7,6 +7,7 @@
 import importlib.metadata
 
 from . import (
+    binning,
     calibration,
     conversion,
     correction,
@@ -27,6 +28,7 @@
 del importlib
 
 providers = (
+    *binning.providers,
     *calibration.providers,
     *conversion.providers,
     *correction.providers,
@@ -41,6 +43,7 @@
 __all__ = [
     "RunNormalization",
     "__version__",
+    "binning",
     "calibration",
     "conversion",
     "correction",

packages/essdiffraction/src/ess/powder/binning.py

@@ -0,0 +1,61 @@
+# SPDX-License-Identifier: BSD-3-Clause
+# Copyright (c) 2026 Scipp contributors (https://github.com/scipp)
+"""Automatic bin-edge generation for powder diffraction."""
+
+import scipp as sc
+import scippneutron as scn
+
+from ess.reduce.unwrap import ChopperFrameSequence
+
+from .types import (
+    DetectorTwoTheta,
+    DspacingBins,
+    DspacingNBins,
+    RunType,
+    SampleRun,
+    WavelengthRange,
+)
+
+
+def wavelength_range_from_chopper_frames(
+    frames: ChopperFrameSequence[RunType],
+) -> WavelengthRange[RunType]:
+    """Return the wavelength envelope transmitted by the chopper cascade."""
+    return WavelengthRange[RunType](frames[-1].bounds()['wavelength'])
+
+
+def dspacing_bins_from_wavelength_and_two_theta(
+    wavelength_range: WavelengthRange[SampleRun],
+    two_theta: DetectorTwoTheta[SampleRun],
+    nbins: DspacingNBins,
+) -> DspacingBins:
+    """Make d-spacing bin edges from wavelength and detector geometry.
+
+    Wavelength and two-theta are independent geometry and beamline quantities.
+    Pairing their opposite extrema and transforming them to d-spacing gives an
+    envelope that covers every detector pixel.
+    """
+    if nbins < 1:
+        raise ValueError(f'DspacingNBins must be positive, got {nbins}.')
+    dspacing_range = scn.conversion.tof.dspacing_from_wavelength(
+        wavelength=sc.concat(
+            [wavelength_range.nanmin(), wavelength_range.nanmax()], dim='bound'
+        ),
+        two_theta=sc.concat([two_theta.nanmax(), two_theta.nanmin()], dim='bound'),
+    )
+    return DspacingBins(
+        sc.linspace(
+            dim='dspacing',
+            start=dspacing_range[0].value,
+            stop=dspacing_range[-1].value,
+            num=nbins + 1,
+            unit=dspacing_range.unit,
+        )
+    )
+
+
+providers = (
+    wavelength_range_from_chopper_frames,
+    dspacing_bins_from_wavelength_and_two_theta,
+)
+"""Sciline providers for automatic powder-diffraction binning."""

packages/essdiffraction/src/ess/powder/types.py

@@ -62,6 +62,9 @@
 DspacingBins = NewType("DspacingBins", sc.Variable)
 """Bin edges for d-spacing."""
 
+DspacingNBins = NewType("DspacingNBins", int)
+"""Number of d-spacing bins."""
+
 OutFilename = NewType("OutFilename", str)
 """Filename of the output."""
 
@@ -115,6 +118,10 @@ class DetectorTwoTheta(sciline.Scope[RunType, sc.Variable], sc.Variable):
     """Scattering angle (two-theta) for each detector pixel."""
 
 
+class WavelengthRange(sciline.Scope[RunType, sc.Variable], sc.Variable):
+    """Minimum and maximum wavelength transmitted by the chopper cascade."""
+
+
 class ElasticCoordTransformGraph(sciline.Scope[RunType, dict], dict):
     """Graph for transforming coordinates in elastic scattering."""