[ESSREFLECTOMETRY] Freia followup

packages/essreflectometry/docs/user-guide/freia/freia-detector-diagnostics.ipynb

@@ -0,0 +1,132 @@
+{
+ "cells": [
+  {
+   "cell_type": "markdown",
+   "id": "0",
+   "metadata": {},
+   "source": [
+    "# FREIA detector diagnostics\n",
+    "\n",
+    "Use an interactive inspector to compare wavelength spectra from selected regions of the FREIA detector. This can help identify where features in a spectrum originate on the detector."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "1",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "%matplotlib widget\n",
+    "import plopp as pp\n",
+    "import scipp as sc\n",
+    "\n",
+    "from ess import freia\n",
+    "from ess.freia import data\n",
+    "from ess.reduce.nexus.types import Filename\n",
+    "from ess.reduce.unwrap import TimeResolution, WavelengthDetector\n",
+    "from ess.reflectometry.types import SampleRun"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "2",
+   "metadata": {},
+   "source": [
+    "## Load detector events\n",
+    "\n",
+    "Use the simulated direct-beam run and reconstruct the wavelength of each event from its arrival time and the WFM chopper settings. To inspect another run, replace the example data with its file path."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "3",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "workflow = freia.FreiaMcStasWorkflow(wavelength_from='analytical')\n",
+    "workflow[Filename[SampleRun]] = data.freia_mcstas_sample_run()\n",
+    "workflow[TimeResolution] = sc.scalar(20.0, unit='us')\n",
+    "\n",
+    "events = workflow.compute(WavelengthDetector[SampleRun])"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "4",
+   "metadata": {},
+   "source": [
+    "## Histogram by detector position and wavelength\n",
+    "\n",
+    "The inspector expects three-dimensional data. Histogram the events by `longitude` and `height` in the detector's local frame, and by reconstructed `wavelength`."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "5",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "wavelength_bins = sc.linspace(\n",
+    "    'wavelength', start=1.0, stop=12.0, num=221, unit='angstrom'\n",
+    ")\n",
+    "detector_data = events.hist(\n",
+    "    height=64,\n",
+    "    longitude=120,\n",
+    "    wavelength=wavelength_bins,\n",
+    "    dim=events.dims,\n",
+    ")"
+   ]
+  },
+  {
+   "cell_type": "markdown",
+   "id": "6",
+   "metadata": {},
+   "source": [
+    "## Inspect detector regions\n",
+    "\n",
+    "The left panel shows the detector image summed over the wavelength range selected with the slider. Activate the rectangle tool in its toolbar, then left-click to add rectangles around regions of interest. The right panel shows the wavelength spectrum summed within each rectangle. Drag a rectangle with the right mouse button, resize it by dragging its vertices, or delete it with the middle mouse button."
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "7",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "pp.inspector(\n",
+    "    detector_data,\n",
+    "    dim='wavelength',\n",
+    "    mode='rectangle',\n",
+    "    logc=True,\n",
+    "    cmin=1.0,\n",
+    "    title='FREIA detector',\n",
+    ")"
+   ]
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "Python 3 (ipykernel)",
+   "language": "python",
+   "name": "python3"
+  },
+  "language_info": {
+   "codemirror_mode": {
+    "name": "ipython",
+    "version": 3
+   },
+   "file_extension": ".py",
+   "mimetype": "text/x-python",
+   "name": "python",
+   "nbconvert_exporter": "python",
+   "pygments_lexer": "ipython3",
+   "version": "3.12.13"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}

packages/essreflectometry/docs/user-guide/freia/freia-reflectivity.ipynb

@@ -42,7 +42,6 @@
     "    Sample,\n",
     "    SampleRun,\n",
     "    SampleSize,\n",
-    "    WavelengthBins,\n",
     ")"
    ]
   },
@@ -68,7 +67,6 @@
     "workflow[Filename[ReferenceRun]] = data.freia_mcstas_reference_run()\n",
     "workflow[NeXusName[IncidentMonitor]] = 'SampleLambda'\n",
     "\n",
-    "workflow[WavelengthBins] = sc.linspace('wavelength', 2.0, 10.0, 81, unit='angstrom')\n",
     "workflow[QBins] = sc.geomspace('Q', 0.003, 0.5, 151, unit='1/angstrom')"
    ]
   },

packages/essreflectometry/docs/user-guide/freia/index.md

@@ -8,7 +8,8 @@ The examples use local McStas files. Set the input paths in each notebook.
 ```{toctree}
 :maxdepth: 1
 
-freia-mcstas-visualization
+freia-visualizations
+freia-detector-diagnostics
 freia-wavelength-lookup-table
 freia-reflectivity
 ```

packages/essreflectometry/src/ess/freia/corrections.py

@@ -15,23 +15,21 @@
     Sample,
     SampleRun,
     SampleSize,
-    WavelengthBins,
     WavelengthDetector,
 )
 from .conversions import add_coords
-from .maskings import add_masks
+from .maskings import add_roi_masks
 from .types import DetectorRegionOfInterest, WavelengthMonitor
 
 
 def add_coords_and_masks(
     da: WavelengthDetector[RunType],
     graph: CoordTransformationGraph[RunType],
     roi: DetectorRegionOfInterest[RunType],
-    wavelength_bins: WavelengthBins,
 ) -> CorrectedDetector[RunType]:
     """Transform coordinates and mask events before run normalization."""
     da = add_coords(da, graph)
-    return CorrectedDetector[RunType](add_masks(da, roi, wavelength_bins))
+    return CorrectedDetector[RunType](add_roi_masks(da, roi))
 
 
 def normalize_by_monitor_histogram(
@@ -78,7 +76,7 @@ def insert_run_normalization(
     )
 
 
-def prepare_sample(
+def correct_sample(
     sample: ReducibleData[SampleRun],
     beam_size: BeamSize[SampleRun],
     sample_size: SampleSize[SampleRun],
@@ -102,4 +100,4 @@ def prepare_sample(
     return Sample(corrected.bins.assign_masks(footprint=invalid))
 
 
-providers = (add_coords_and_masks, prepare_sample)
+providers = (add_coords_and_masks, correct_sample)

packages/essreflectometry/src/ess/freia/maskings.py

@@ -3,12 +3,14 @@
 import scipp as sc
 
 
-def add_masks(
+def add_roi_masks(
     da: sc.DataArray,
     roi: dict,
-    wavelength_bins: sc.Variable,
 ) -> sc.DataArray:
-    """Mask events outside the ROI and wavelength range."""
+    """Mask events outside the ROI.
+    When a coordinate with the same name exists both on the bins
+    and on the events the method masks the event coordinate.
+    """
     masks = {}
     event_masks = {}
     for name, (low, high) in roi.items():
@@ -20,12 +22,4 @@ def add_masks(
         (event_masks if is_event_coord else masks)[f'roi_{name}'] = ~(
             (coord >= low) & (coord <= high)
         )
-    wavelength = da.bins.coords['wavelength']
-    return da.assign_masks(masks).bins.assign_masks(
-        event_masks,
-        wavelength=~(
-            sc.isfinite(wavelength)
-            & (wavelength >= wavelength_bins[0].to(unit=wavelength.unit))
-            & (wavelength < wavelength_bins[-1].to(unit=wavelength.unit))
-        ),
-    )
+    return da.assign_masks(masks).bins.assign_masks(event_masks)

packages/essreflectometry/src/ess/freia/normalization.py

@@ -21,6 +21,9 @@ def evaluate_direct_beam(
     normal = direct_beam.coords['sample_surface_normal']
     normal = normal / sc.norm(normal)
     outgoing = direct_beam.bins.coords['outgoing_direction']
+    # The reflected direction is the direction the direct beam
+    # would have gone if it had been reflected.
+    # It is the direct beam reflected by the sample surface.
     reflected = outgoing - 2 * sc.dot(outgoing, normal) * normal
     reflection_angle = theta(
         outgoing_direction=reflected,
@@ -37,11 +40,7 @@ def reduce_sample_over_q(
     reference: Reference,
     qbins: QBins,
 ) -> ReflectivityOverQ:
-    """Divide ROI intensities on a common Q grid, propagating both variances.
-
-    Histogram before division so changing Q bin widths does not rescale R.
-    Empty, masked, or nonfinite direct-beam bins provide no normalization.
-    """
+    """Divide ROI intensities on a common Q grid, propagating both variances."""
     numerator = sample.hist(Q=qbins, dim=sample.dims)
     denominator = reference.hist(Q=qbins, dim=reference.dims)
     valid = sc.isfinite(denominator.data) & (