Beer powdiff workflow

docs/user-guide/beer/beer_modulation_mcstas.ipynb

@@ -22,6 +22,7 @@
     "from ess.beer import BeerModMcStasWorkflow, BeerModMcStasWorkflowKnownPeaks\n",
     "from ess.beer.data import mcstas_silicon_medium_resolution, mcstas_duplex, duplex_peaks_array, silicon_peaks_array\n",
     "from ess.beer.types import *\n",
+    "from ess.powder.types import *\n",
     "\n",
     "# Default bin edges for our d_hkl histograms\n",
     "dspacing = sc.linspace('dspacing', 0.8, 2.2, 2000, unit='angstrom')\n",
@@ -163,7 +164,8 @@
     "wf = BeerModMcStasWorkflowKnownPeaks()\n",
     "wf[DetectorBank] = DetectorBank.north\n",
     "wf[Filename[SampleRun]] = mcstas_silicon_medium_resolution()\n",
-    "da = wf.compute(RawDetector[SampleRun])\n",
+    "wf[DHKLList] = silicon_peaks_array()\n",
+    "da = wf.compute(WavelengthDetector[SampleRun])\n",
     "da.masks.clear()\n",
     "da.hist(two_theta=400, event_time_offset=1000).plot(norm='log', cmin=1.0e-3)"
    ]
@@ -183,13 +185,12 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "wf[DHKLList] = silicon_peaks_array()\n",
     "wf[Filename[SampleRun]] = mcstas_silicon_medium_resolution()\n",
     "\n",
     "results = {}\n",
     "for bank in DetectorBank:\n",
     "    wf[DetectorBank] = bank\n",
-    "    da = wf.compute(TofDetector[SampleRun])\n",
+    "    da = wf.compute(WavelengthDetector[SampleRun])\n",
     "    results[bank] = (\n",
     "        da\n",
     "        .transform_coords(('dspacing',), graph=scn.conversion.graph.tof.elastic('tof'),)\n",
@@ -220,7 +221,7 @@
     "results = {}\n",
     "for bank in DetectorBank:\n",
     "    wf[DetectorBank] = bank\n",
-    "    da = wf.compute(TofDetector[SampleRun])\n",
+    "    da = wf.compute(WavelengthDetector[SampleRun])\n",
     "    results[bank] = (\n",
     "        da\n",
     "        .transform_coords(('dspacing',), graph=scn.conversion.graph.tof.elastic('tof'),)\n",
@@ -279,7 +280,8 @@
     "wf = BeerModMcStasWorkflowKnownPeaks()\n",
     "wf[DetectorBank] = DetectorBank.south\n",
     "wf[Filename[SampleRun]] = mcstas_duplex(8)\n",
-    "wf.compute(RawDetector[SampleRun]).hist(two_theta=400, event_time_offset=1000).plot(norm='log', cmin=1.0e-2)"
+    "wf[DHKLList] = duplex_peaks_array()\n",
+    "wf.compute(WavelengthDetector[SampleRun]).hist(two_theta=400, event_time_offset=1000).plot(norm='log', cmin=1.0e-2)"
    ]
   },
   {
@@ -297,12 +299,10 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "wf[DHKLList] = duplex_peaks_array()\n",
-    "\n",
     "results = {}\n",
     "for bank in DetectorBank:\n",
     "    wf[DetectorBank] = bank\n",
-    "    da = wf.compute(TofDetector[SampleRun])\n",
+    "    da = wf.compute(WavelengthDetector[SampleRun])\n",
     "    results[bank] = (\n",
     "        da\n",
     "        .transform_coords(('dspacing',), graph=scn.conversion.graph.tof.elastic('tof'),)\n",
@@ -333,7 +333,7 @@
     "results = {}\n",
     "for bank in DetectorBank:\n",
     "    wf[DetectorBank] = bank\n",
-    "    da = wf.compute(TofDetector[SampleRun])\n",
+    "    da = wf.compute(WavelengthDetector[SampleRun])\n",
     "    results[bank] = (\n",
     "        da\n",
     "        .transform_coords(('dspacing',), graph=scn.conversion.graph.tof.elastic('tof'),)\n",
@@ -392,7 +392,8 @@
     "wf = BeerModMcStasWorkflowKnownPeaks()\n",
     "wf[DetectorBank] = DetectorBank.south\n",
     "wf[Filename[SampleRun]] = mcstas_duplex(9)\n",
-    "wf.compute(RawDetector[SampleRun]).hist(two_theta=400, event_time_offset=1000).plot(norm='log', cmin=1.0e-3)"
+    "wf[DHKLList] = duplex_peaks_array()\n",
+    "wf.compute(WavelengthDetector[SampleRun]).hist(two_theta=400, event_time_offset=1000).plot(norm='log', cmin=1.0e-3)"
    ]
   },
   {
@@ -410,12 +411,10 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "wf[DHKLList] = duplex_peaks_array()\n",
-    "\n",
     "results = {}\n",
     "for bank in DetectorBank:\n",
     "    wf[DetectorBank] = bank\n",
-    "    da = wf.compute(TofDetector[SampleRun])\n",
+    "    da = wf.compute(WavelengthDetector[SampleRun])\n",
     "    results[bank] = (\n",
     "        da\n",
     "        .transform_coords(('dspacing',), graph=scn.conversion.graph.tof.elastic('tof'),)\n",
@@ -446,7 +445,7 @@
     "results = {}\n",
     "for bank in DetectorBank:\n",
     "    wf[DetectorBank] = bank\n",
-    "    da = wf.compute(TofDetector[SampleRun])\n",
+    "    da = wf.compute(WavelengthDetector[SampleRun])\n",
     "    results[bank] = (\n",
     "        da\n",
     "        .transform_coords(('dspacing',), graph=scn.conversion.graph.tof.elastic('tof'),)\n",
@@ -505,7 +504,8 @@
     "wf = BeerModMcStasWorkflowKnownPeaks()\n",
     "wf[DetectorBank] = DetectorBank.south\n",
     "wf[Filename[SampleRun]] = mcstas_duplex(10)\n",
-    "wf.compute(RawDetector[SampleRun]).hist(two_theta=400, event_time_offset=1000).plot(norm='log', cmin=1.0e-3)"
+    "wf[DHKLList] = duplex_peaks_array()\n",
+    "wf.compute(WavelengthDetector[SampleRun]).hist(two_theta=400, event_time_offset=1000).plot(norm='log', cmin=1.0e-3)"
    ]
   },
   {
@@ -523,12 +523,10 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "wf[DHKLList] = duplex_peaks_array()\n",
-    "\n",
     "results = {}\n",
     "for bank in DetectorBank:\n",
     "    wf[DetectorBank] = bank\n",
-    "    da = wf.compute(TofDetector[SampleRun])\n",
+    "    da = wf.compute(WavelengthDetector[SampleRun])\n",
     "    results[bank] = (\n",
     "        da\n",
     "        .transform_coords(('dspacing',), graph=scn.conversion.graph.tof.elastic('tof'),)\n",
@@ -559,7 +557,7 @@
     "results = {}\n",
     "for bank in DetectorBank:\n",
     "    wf[DetectorBank] = bank\n",
-    "    da = wf.compute(TofDetector[SampleRun])\n",
+    "    da = wf.compute(WavelengthDetector[SampleRun])\n",
     "    results[bank] = (\n",
     "        da\n",
     "        .transform_coords(('dspacing',), graph=scn.conversion.graph.tof.elastic('tof'),)\n",
@@ -618,7 +616,8 @@
     "wf = BeerModMcStasWorkflowKnownPeaks()\n",
     "wf[DetectorBank] = DetectorBank.south\n",
     "wf[Filename[SampleRun]] = mcstas_duplex(16)\n",
-    "wf.compute(RawDetector[SampleRun]).hist(two_theta=400, event_time_offset=1000).plot(norm='log', cmin=1.0e-3)"
+    "wf[DHKLList] = duplex_peaks_array()\n",
+    "wf.compute(WavelengthDetector[SampleRun]).hist(two_theta=400, event_time_offset=1000).plot(norm='log', cmin=1.0e-3)"
    ]
   },
   {
@@ -636,12 +635,10 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "wf[DHKLList] = duplex_peaks_array()\n",
-    "\n",
     "results = {}\n",
     "for bank in DetectorBank:\n",
     "    wf[DetectorBank] = bank\n",
-    "    da = wf.compute(TofDetector[SampleRun])\n",
+    "    da = wf.compute(WavelengthDetector[SampleRun])\n",
     "    results[bank] = (\n",
     "        da\n",
     "        .transform_coords(('dspacing',), graph=scn.conversion.graph.tof.elastic('tof'),)\n",
@@ -672,7 +669,7 @@
     "results = {}\n",
     "for bank in DetectorBank:\n",
     "    wf[DetectorBank] = bank\n",
-    "    da = wf.compute(TofDetector[SampleRun])\n",
+    "    da = wf.compute(WavelengthDetector[SampleRun])\n",
     "    results[bank] = (\n",
     "        da\n",
     "        .transform_coords(('dspacing',), graph=scn.conversion.graph.tof.elastic('tof'),)\n",

src/ess/beer/clustering.py

@@ -1,19 +1,20 @@
 import scipp as sc
 from scipy.signal import find_peaks, medfilt
 
+from ess.powder.types import RunType
+
 from .conversions import tof_from_t0_estimate_graph
 from .types import (
     GeometryCoordTransformGraph,
     RawDetector,
-    RunType,
     StreakClusteredData,
 )
 
 
 def cluster_events_by_streak(
     da: RawDetector[RunType], gg: GeometryCoordTransformGraph
 ) -> StreakClusteredData[RunType]:
-    graph = tof_from_t0_estimate_graph(gg)
+    graph = tof_from_t0_estimate_graph(da, gg)
 
     da = da.transform_coords(['dspacing'], graph=graph)
     da.bins.coords['coarse_d'] = da.bins.coords.pop('dspacing').to(unit='angstrom')

src/ess/beer/conversions.py

@@ -2,25 +2,26 @@
 import scipp.constants
 from scippneutron.conversion import graph
 
+from ess.powder.types import ElasticCoordTransformGraph, RunType
+
 from .types import (
     DHKLList,
     GeometryCoordTransformGraph,
     ModulationPeriod,
     PulseLength,
     RawDetector,
-    RunType,
     StreakClusteredData,
-    TofCoordTransformGraph,
-    TofDetector,
     WavelengthDefinitionChopperDelay,
+    WavelengthDetector,
 )
 
 
-def compute_tof_in_each_cluster(
+def compute_wavelength_in_each_cluster(
     da: StreakClusteredData[RunType],
     chopper_delay: WavelengthDefinitionChopperDelay,
     mod_period: ModulationPeriod,
-) -> TofDetector[RunType]:
+    graph: GeometryCoordTransformGraph,
+) -> WavelengthDetector[RunType]:
     """Fits a line through each cluster, the intercept of the line is t0.
     The line is fitted using linear regression with an outlier removal procedure.
 
@@ -37,7 +38,10 @@ def compute_tof_in_each_cluster(
     """
     if isinstance(da, sc.DataGroup):
         return sc.DataGroup(
-            {k: compute_tof_in_each_cluster(v, mod_period) for k, v in da.items()}
+            {
+                k: compute_wavelength_in_each_cluster(v, mod_period)
+                for k, v in da.items()
+            }
         )
 
     max_distance_from_streak_line = mod_period / 3
@@ -185,17 +189,48 @@ def _tof_from_dhkl(
     return out
 
 
+def t0_estimate(
+    wavelength_estimate: sc.Variable,
+    L0: sc.Variable,
+    Ltotal: sc.Variable,
+) -> sc.Variable:
+    """Estimates the time-at-chopper by assuming the wavelength."""
+    return (
+        sc.constants.m_n
+        / sc.constants.h
+        * wavelength_estimate
+        * (L0 - Ltotal).to(unit=wavelength_estimate.unit)
+    ).to(unit='s')
+
+
+def tof_from_t0_estimate_graph(
+    da: RawDetector[RunType],
+    gg: GeometryCoordTransformGraph,
+) -> ElasticCoordTransformGraph[RunType]:
+    """Graph for computing ``wavelength`` in pulse shaping chopper modes."""
+    return {
+        **gg,
+        't0': t0_estimate,
+        'tof': lambda time_of_arrival, t0: time_of_arrival - t0,
+        'time_of_arrival': time_of_arrival,
+    }
+
+
 def geometry_graph() -> GeometryCoordTransformGraph:
-    return graph.beamline.beamline(scatter=True)
+    return {
+        **graph.beamline.beamline(scatter=True),
+        **graph.tof.elastic("tof"),
+    }
 
 
 def tof_from_known_dhkl_graph(
+    da: RawDetector[RunType],
     mod_period: ModulationPeriod,
     pulse_length: PulseLength,
     chopper_delay: WavelengthDefinitionChopperDelay,
     dhkl_list: DHKLList,
     gg: GeometryCoordTransformGraph,
-) -> TofCoordTransformGraph:
+) -> ElasticCoordTransformGraph[RunType]:
     """Graph computing ``tof`` in modulation chopper modes using
     list of peak positions."""
 
@@ -221,7 +256,6 @@ def _compute_coarse_dspacing(
 
     return {
         **gg,
-        **graph.tof.elastic("tof"),
         'pulse_length': lambda: pulse_length,
         'mod_period': lambda: mod_period,
         'chopper_delay': lambda: chopper_delay,
@@ -232,56 +266,21 @@ def _compute_coarse_dspacing(
     }
 
 
-def t0_estimate(
-    wavelength_estimate: sc.Variable,
-    L0: sc.Variable,
-    Ltotal: sc.Variable,
-) -> sc.Variable:
-    """Estimates the time-at-chopper by assuming the wavelength."""
-    return (
-        sc.constants.m_n
-        / sc.constants.h
-        * wavelength_estimate
-        * (L0 - Ltotal).to(unit=wavelength_estimate.unit)
-    ).to(unit='s')
-
-
-def _tof_from_t0(
-    time_of_arrival: sc.Variable,
-    t0: sc.Variable,
-) -> sc.Variable:
-    """Computes time-of-flight by subtracting a start time."""
-    return time_of_arrival - t0
-
-
-def tof_from_t0_estimate_graph(
-    gg: GeometryCoordTransformGraph,
-) -> TofCoordTransformGraph:
-    """Graph for computing ``tof`` in pulse shaping chopper modes."""
-    return {
-        **gg,
-        **graph.tof.elastic("tof"),
-        't0': t0_estimate,
-        'tof': _tof_from_t0,
-        'time_of_arrival': time_of_arrival,
-    }
-
-
-def compute_tof(
-    da: RawDetector[RunType], graph: TofCoordTransformGraph
-) -> TofDetector[RunType]:
-    """Uses the transformation graph to compute ``tof``."""
-    return da.transform_coords(('tof',), graph=graph)
+def wavelength_detector(
+    da: RawDetector[RunType], graph: ElasticCoordTransformGraph[RunType]
+) -> WavelengthDetector[RunType]:
+    """Applies the transformation graph to compute ``wavelength``."""
+    return da.transform_coords(('wavelength',), graph=graph)
 
 
 convert_from_known_peaks_providers = (
     geometry_graph,
     tof_from_known_dhkl_graph,
-    compute_tof,
+    wavelength_detector,
 )
 convert_pulse_shaping = (
     geometry_graph,
     tof_from_t0_estimate_graph,
-    compute_tof,
+    wavelength_detector,
 )
-providers = (compute_tof_in_each_cluster, geometry_graph)
+providers = (compute_wavelength_in_each_cluster, geometry_graph)