Merge pull request #647 from scipp/add-pulse-shaping-choppers

[ESSDIFFRACTION] feat: add workflow using lookup tables for pulse shaping modes

Author: jokasimr

packages/essdiffraction/docs/api-reference/index.md

@@ -123,8 +123,14 @@
 .. autosummary::
    :toctree: ../generated/functions
 
+   BeerMcStasWorkflowPulseShaping
+   BeerMcStasWorkflowPulseShapingAnalytical
    BeerModMcStasWorkflowKnownPeaks
    BeerModMcStasWorkflow
+   BeerPowderWorkflow
+   BeerPowderWorkflowAnalytical
+   BeerPowderMcStasWorkflow
+   BeerPowderMcStasWorkflowAnalytical
 ```
 
 ### Top-level functions
@@ -135,7 +141,8 @@
 .. autosummary::
    :toctree: ../generated/functions
 
-    load_beer_mcstas
+   dhkl_peaks_from_cif
+   load_beer_mcstas
 ```
 
 ### Submodules
@@ -150,6 +157,7 @@
    conversions
    data
    io
+   mcstas
    workflow
 ```
 

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

@@ -7,12 +7,17 @@
 
 import importlib.metadata
 
-from .io import load_beer_mcstas
+from .mcstas import load_beer_mcstas
 from .peakfinding import dhkl_peaks_from_cif
 from .workflow import (
     BeerMcStasWorkflowPulseShaping,
+    BeerMcStasWorkflowPulseShapingAnalytical,
     BeerModMcStasWorkflow,
     BeerModMcStasWorkflowKnownPeaks,
+    BeerPowderMcStasWorkflow,
+    BeerPowderMcStasWorkflowAnalytical,
+    BeerPowderWorkflow,
+    BeerPowderWorkflowAnalytical,
     default_parameters,
 )
 
@@ -25,8 +30,13 @@
 
 __all__ = [
     'BeerMcStasWorkflowPulseShaping',
+    'BeerMcStasWorkflowPulseShapingAnalytical',
     'BeerModMcStasWorkflow',
     'BeerModMcStasWorkflowKnownPeaks',
+    'BeerPowderMcStasWorkflow',
+    'BeerPowderMcStasWorkflowAnalytical',
+    'BeerPowderWorkflow',
+    'BeerPowderWorkflowAnalytical',
     '__version__',
     'default_parameters',
     'dhkl_peaks_from_cif',

packages/essdiffraction/src/ess/beer/clustering.py

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

packages/essdiffraction/src/ess/beer/conversions.py

@@ -19,7 +19,6 @@ def compute_wavelength_in_each_cluster(
     da: StreakClusteredData[RunType],
     chopper_delay: WavelengthDefinitionChopperDelay,
     mod_period: ModulationPeriod,
-    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.
@@ -38,7 +37,7 @@ def compute_wavelength_in_each_cluster(
     if isinstance(da, sc.DataGroup):
         return sc.DataGroup(
             {
-                k: compute_wavelength_in_each_cluster(v, mod_period)
+                k: compute_wavelength_in_each_cluster(v, chopper_delay, mod_period)
                 for k, v in da.items()
             }
         )
@@ -133,7 +132,7 @@ def _compute_d_given_list_of_peaks(
         / (scipp.constants.h / scipp.constants.m_n)
     ).to(unit=f'{time_of_arrival.unit}/angstrom')
     for dhkl in dhkl_list:
-        dt = sc.abs(t - dhkl * const)
+        dt = sc.abs(t - dhkl * const - pulse_length / 2)
         dt_in_range = dt < pulse_length / 2
         no_dt_found = sc.isnan(dtfound)
         dtfound = sc.where(dt_in_range, sc.where(no_dt_found, dt, dtfound), dtfound)
@@ -191,34 +190,16 @@ def _tof_from_dhkl(
     return out
 
 
-def t0_estimate(
-    wavelength_estimate: sc.Variable,
-    source_to_wavelength_definition_chopper_distance: sc.Variable,
-) -> sc.Variable:
-    """
-    Computes the time a neutron reaches a chopper at
-    ``source_to_wavelength_chopper_distance`` distance from the source
-    if it has wavelength ``wavelength_estimate``.
-    """
-    return (
-        sc.constants.m_n
-        / sc.constants.h
-        * wavelength_estimate
-        * source_to_wavelength_definition_chopper_distance.to(
-            unit=wavelength_estimate.unit
-        )
-    ).to(unit='s')
-
-
-def tof_from_t0_estimate_graph(
+def tof_from_nominal_time_at_chopper_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,
+        'tof': lambda time_of_arrival, nominal_time_at_chopper: (
+            time_of_arrival - nominal_time_at_chopper
+        ),
         'time_of_arrival': time_of_arrival,
     }
 
@@ -287,7 +268,7 @@ def wavelength_detector(
 )
 convert_pulse_shaping = (
     geometry_graph,
-    tof_from_t0_estimate_graph,
+    tof_from_nominal_time_at_chopper_graph,
     wavelength_detector,
 )
 providers = (compute_wavelength_in_each_cluster, geometry_graph)

packages/essdiffraction/src/ess/beer/data.py

@@ -31,7 +31,7 @@
         # Simulations from new model with corrected(?) moderator to detector
         # distance.
         # For correct reduction you need to use
-        # beer.io.mcstas_chopper_delay_from_mode_new_simulations
+        # beer.mcstas.mcstas_chopper_delay_from_mode_new_simulations
         # to obtain the correct WavelengthDefinitionChopperDelay for these files.
         "silicon-mode10-new-model.h5": "md5:98500830f27700fc719634e1acd49944",
         "silicon-mode16-new-model.h5": "md5:393f9287e7d3f97ceedbe64343918413",