diff --git a/docs/user-guide/common/beam-center-finder.ipynb b/docs/user-guide/common/beam-center-finder.ipynb index 379f9c64..609e71b6 100644 --- a/docs/user-guide/common/beam-center-finder.ipynb +++ b/docs/user-guide/common/beam-center-finder.ipynb @@ -76,7 +76,7 @@ "metadata": {}, "outputs": [], "source": [ - "raw = workflow.compute(RawData[SampleRun])['spectrum', :61440]\n", + "raw = workflow.compute(RawDetector[SampleRun])['spectrum', :61440]\n", "\n", "p = isis.plot_flat_detector_xy(raw.hist(), norm='log')\n", "p.ax.plot(0, 0, '+', color='k', ms=10)\n", @@ -128,6 +128,7 @@ "metadata": {}, "outputs": [], "source": [ + "workflow[BeamCenter] = sc.vector(value=[0, 0, 0], unit='m')\n", "masked = workflow.compute(MaskedData[SampleRun])['spectrum', :61440].copy()\n", "masked.masks['low_counts'] = masked.hist().data < sc.scalar(80.0, unit='counts')\n", "\n", @@ -156,9 +157,16 @@ "metadata": {}, "outputs": [], "source": [ - "# Insert the provider that computes the center-of-mass of the pixels\n", - "workflow.insert(sans.beam_center_finder.beam_center_from_center_of_mass)\n", - "workflow.visualize(BeamCenter, graph_attr={'rankdir': 'LR'})" + "# The center-of-mass approach is based on the MaskedData\n", + "workflow.visualize(MaskedData[SampleRun])" + ] + }, + { + "cell_type": "markdown", + "id": "d7aac9d0", + "metadata": {}, + "source": [ + "We use the workflow to fill in the required arguments and call the function:" ] }, { @@ -168,7 +176,7 @@ "metadata": {}, "outputs": [], "source": [ - "com = workflow.compute(BeamCenter)\n", + "com = sans.beam_center_finder.beam_center_from_center_of_mass(workflow)\n", "com" ] }, @@ -245,11 +253,7 @@ "workflow = isis.sans2d.Sans2dTutorialWorkflow()\n", "# For real data use:\n", "# workflow = isis.sans2d.Sans2dWorkflow()\n", - "\n", - "workflow.insert(isis.data.transmission_from_sample_run)\n", - "# Insert the provider that computes the center from I(Q) curves\n", - "workflow.insert(sans.beam_center_finder.beam_center_from_iofq)\n", - "workflow.visualize(BeamCenter, graph_attr={'rankdir': 'LR'})" + "workflow.insert(isis.data.transmission_from_sample_run)" ] }, { @@ -280,12 +284,8 @@ "\n", "workflow[CorrectForGravity] = True\n", "workflow[UncertaintyBroadcastMode] = UncertaintyBroadcastMode.upper_bound\n", - "workflow[sans.beam_center_finder.BeamCenterFinderQBins] = sc.linspace(\n", - " 'Q', 0.02, 0.25, 71, unit='1/angstrom'\n", - ")\n", - "workflow[sans.beam_center_finder.BeamCenterFinderMinimizer] = None\n", - "workflow[sans.beam_center_finder.BeamCenterFinderTolerance] = None\n", - "workflow[DirectBeam] = None" + "workflow[DirectBeam] = None\n", + "workflow[isis.sans2d.LowCountThreshold] = sc.scalar(-1, unit=\"counts\")" ] }, { @@ -293,7 +293,7 @@ "id": "641be8e9", "metadata": {}, "source": [ - "Finally, we set the data to be used, including overriding with data that contains the new mask defined earlier:" + "Finally, we set the data to be used, including overriding with the new mask defined earlier:" ] }, { @@ -304,7 +304,10 @@ "outputs": [], "source": [ "workflow[Filename[SampleRun]] = isis.data.sans2d_tutorial_sample_run()\n", - "workflow[MaskedData[SampleRun]] = masked" + "detector = workflow.compute(RawDetector[SampleRun])['spectrum', :61440].assign_masks(\n", + " masked.masks\n", + ")\n", + "workflow[RawDetector[SampleRun]] = detector" ] }, { @@ -322,11 +325,15 @@ { "cell_type": "code", "execution_count": null, - "id": "c087798a-fc33-4292-924c-eed8e3baf36a", + "id": "329a0cd9", "metadata": {}, "outputs": [], "source": [ - "iofq_center = workflow.compute(BeamCenter)\n", + "q_bins = sc.linspace('Q', 0.02, 0.25, 71, unit='1/angstrom')\n", + "workflow[BeamCenter] = sc.vector([0, 0, 0], unit='m')\n", + "iofq_center = sans.beam_center_finder.beam_center_from_iofq(\n", + " workflow=workflow, q_bins=q_bins\n", + ")\n", "iofq_center" ] }, @@ -433,14 +440,16 @@ "metadata": {}, "outputs": [], "source": [ + "workflow = workflow.copy()\n", + "workflow[QBins] = q_bins\n", + "workflow[ReturnEvents] = False\n", + "workflow[DimsToKeep] = ()\n", + "workflow[WavelengthMask] = None\n", + "workflow[WavelengthBands] = None\n", "kwargs = dict( # noqa: C408\n", - " data=masked,\n", + " workflow=workflow,\n", + " detector=detector,\n", " norm=workflow.compute(NormWavelengthTerm[SampleRun]),\n", - " graph=workflow.compute(sans.conversions.ElasticCoordTransformGraph),\n", - " q_bins=workflow.compute(sans.beam_center_finder.BeamCenterFinderQBins),\n", - " wavelength_bins=workflow.compute(WavelengthBins),\n", - " transform=workflow.compute(LabFrameTransform[SampleRun]),\n", - " pixel_shape=workflow.compute(DetectorPixelShape[SampleRun]),\n", ")" ] }, @@ -628,7 +637,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.10.12" + "version": "3.10.14" } }, "nbformat": 4, diff --git a/docs/user-guide/isis/sans2d.ipynb b/docs/user-guide/isis/sans2d.ipynb index 580a8f72..c67717b3 100644 --- a/docs/user-guide/isis/sans2d.ipynb +++ b/docs/user-guide/isis/sans2d.ipynb @@ -72,8 +72,7 @@ "outputs": [], "source": [ "workflow.insert(isis.data.transmission_from_background_run)\n", - "workflow.insert(isis.data.transmission_from_sample_run)\n", - "workflow.insert(sans.beam_center_from_center_of_mass)" + "workflow.insert(isis.data.transmission_from_sample_run)" ] }, { @@ -167,6 +166,44 @@ "workflow[Filename[EmptyBeamRun]] = isis.data.sans2d_tutorial_empty_beam_run()" ] }, + { + "cell_type": "markdown", + "id": "981f1f28", + "metadata": {}, + "source": [ + "The beam center is not set yet.\n", + "Unless we have a previously computed value, we can do so now:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "064f87a4", + "metadata": {}, + "outputs": [], + "source": [ + "center = sans.beam_center_from_center_of_mass(workflow)\n", + "center" + ] + }, + { + "cell_type": "markdown", + "id": "7d7ce61a", + "metadata": {}, + "source": [ + "Remember to update the workflow with the new center:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "0a4c43f9", + "metadata": {}, + "outputs": [], + "source": [ + "workflow[BeamCenter] = center" + ] + }, { "cell_type": "markdown", "id": "c19eeaf0", @@ -391,7 +428,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.10.12" + "version": "3.10.14" } }, "nbformat": 4, diff --git a/docs/user-guide/isis/zoom.ipynb b/docs/user-guide/isis/zoom.ipynb index 03e347e1..30d6a81c 100644 --- a/docs/user-guide/isis/zoom.ipynb +++ b/docs/user-guide/isis/zoom.ipynb @@ -71,8 +71,7 @@ "metadata": {}, "source": [ "We can insert steps for configuring the workflow.\n", - "In this case, we would like to use the transmission monitor from the regular background and sample runs since there was no separate transmission run.\n", - "We also want to compute the beam center using a simple center-of-mass estimation:" + "In this case, we would like to use the transmission monitor from the regular background and sample runs since there was no separate transmission run." ] }, { @@ -83,8 +82,7 @@ "outputs": [], "source": [ "workflow.insert(isis.data.transmission_from_background_run)\n", - "workflow.insert(isis.data.transmission_from_sample_run)\n", - "workflow.insert(sans.beam_center_from_center_of_mass)" + "workflow.insert(isis.data.transmission_from_sample_run)" ] }, { @@ -176,6 +174,29 @@ "source": [ "## Use the workflow\n", "\n", + "### Set or compute the beam center\n", + "\n", + "The beam center is not set by default.\n", + "We can either set it to a known value, or compute it from the data:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "568cd4be", + "metadata": {}, + "outputs": [], + "source": [ + "workflow[BeamCenter] = sans.beam_center_from_center_of_mass(workflow)" + ] + }, + { + "cell_type": "markdown", + "id": "13ddd330", + "metadata": {}, + "source": [ + "\n", + "\n", "### Compute final result\n", "\n", "We can now compute $I(Q)$:" diff --git a/docs/user-guide/loki/loki-direct-beam.ipynb b/docs/user-guide/loki/loki-direct-beam.ipynb index 8c95902a..9b18aa35 100644 --- a/docs/user-guide/loki/loki-direct-beam.ipynb +++ b/docs/user-guide/loki/loki-direct-beam.ipynb @@ -190,36 +190,8 @@ "and inserting the provider that will compute the beam center.\n", "\n", "For now, we compute the beam center only for the rear detector (named 'larmor_detector') but apply it to all banks (currently there is only one bank).\n", - "The beam center may need to be computed or applied differently to each bank, see [scipp/esssans#28](https://github.com/scipp/esssans/issues/28)." - ] - }, - { - "cell_type": "code", - "execution_count": null, - "id": "bdc15ab4-6ec6-4a29-81dc-1d00a8e890dc", - "metadata": {}, - "outputs": [], - "source": [ - "bc_workflow = workflow.copy()\n", - "bc_workflow.insert(sans.beam_center_from_center_of_mass)" - ] - }, - { - "cell_type": "code", - "execution_count": null, - "id": "7cc2cf3b-973e-4425-93ea-deb30b12c956", - "metadata": {}, - "outputs": [], - "source": [ - "bc_workflow.visualize(BeamCenter, compact=True, graph_attr={'rankdir': 'LR'})" - ] - }, - { - "cell_type": "markdown", - "id": "b6bf5342-6768-4b65-882b-aefc9b583724", - "metadata": {}, - "source": [ - "We can now compute the value for the center" + "The beam center may need to be computed or applied differently to each bank, see [scipp/esssans#28](https://github.com/scipp/esssans/issues/28).\n", + "We use a center-of-mass approach to find the beam center:" ] }, { @@ -229,7 +201,7 @@ "metadata": {}, "outputs": [], "source": [ - "center = bc_workflow.compute(BeamCenter)\n", + "center = sans.beam_center_from_center_of_mass(workflow)\n", "center" ] }, @@ -238,7 +210,7 @@ "id": "0814805a-9611-4951-a015-c12ae254b099", "metadata": {}, "source": [ - "and set that value onto our original workflow" + "and set that value in our workflow" ] }, { @@ -639,7 +611,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.10.12" + "version": "3.10.14" } }, "nbformat": 4, diff --git a/src/ess/isissans/components.py b/src/ess/isissans/components.py index 9f5b2b14..8141aac5 100644 --- a/src/ess/isissans/components.py +++ b/src/ess/isissans/components.py @@ -6,9 +6,10 @@ import scipp as sc from ..sans.types import ( + BeamCenter, + CalibratedDetector, MonitorType, RawDetector, - RawDetectorData, RawMonitor, RawMonitorData, RunType, @@ -24,11 +25,9 @@ class MonitorOffset(sciline.Scope[MonitorType, sc.Variable], sc.Variable): DetectorBankOffset = NewType('DetectorBankOffset', sc.Variable) -def apply_component_user_offsets_to_raw_data( - data: RawDetector[ScatteringRunType], - sample_offset: SampleOffset, - detector_bank_offset: DetectorBankOffset, -) -> RawDetectorData[ScatteringRunType]: +def apply_beam_center( + data: RawDetector[ScatteringRunType], beam_center: BeamCenter +) -> CalibratedDetector[ScatteringRunType]: """Apply user offsets to raw data. Parameters @@ -40,14 +39,9 @@ def apply_component_user_offsets_to_raw_data( detector_bank_offset: Detector bank offset. """ - data = data.copy(deep=False) - sample_pos = data.coords['sample_position'] - data.coords['sample_position'] = sample_pos + sample_offset.to( - unit=sample_pos.unit, copy=False + return CalibratedDetector[ScatteringRunType]( + data.assign_coords(position=data.coords['position'] - beam_center) ) - pos = data.coords['position'] - data.coords['position'] = pos + detector_bank_offset.to(unit=pos.unit, copy=False) - return RawDetectorData[ScatteringRunType](data) def apply_component_user_offsets_to_raw_monitor( @@ -69,6 +63,6 @@ def apply_component_user_offsets_to_raw_monitor( providers = ( - apply_component_user_offsets_to_raw_data, + apply_beam_center, apply_component_user_offsets_to_raw_monitor, ) diff --git a/src/ess/isissans/general.py b/src/ess/isissans/general.py index f937f6d4..5d16247b 100644 --- a/src/ess/isissans/general.py +++ b/src/ess/isissans/general.py @@ -7,7 +7,9 @@ import scipp as sc from ..sans.types import ( + CalibratedDetector, CorrectForGravity, + DetectorData, DetectorIDs, DetectorPixelShape, DimsToKeep, @@ -17,7 +19,6 @@ NeXusMonitorName, NonBackgroundWavelengthRange, RawDetector, - RawDetectorData, RawMonitor, RawMonitorData, RunNumber, @@ -53,8 +54,35 @@ def default_parameters() -> dict: def get_detector_data( dg: LoadedFileContents[RunType], + sample_offset: SampleOffset, + detector_bank_offset: DetectorBankOffset, ) -> RawDetector[RunType]: - return RawDetector[RunType](dg['data']) + """Get detector data and apply user offsets to raw data. + + Parameters + ---------- + dg: + Data loaded with Mantid and converted to Scipp. + sample_offset: + Sample offset. + detector_bank_offset: + Detector bank offset. + """ + data = dg['data'] + sample_pos = data.coords['sample_position'] + sample_pos = sample_pos + sample_offset.to(unit=sample_pos.unit, copy=False) + pos = data.coords['position'] + pos = pos + detector_bank_offset.to(unit=pos.unit, copy=False) + return RawDetector[RunType]( + dg['data'].assign_coords(position=pos, sample_position=sample_pos) + ) + + +def assemble_detector_data( + detector: CalibratedDetector[RunType], +) -> DetectorData[RunType]: + """Dummy assembly of detector data, detector already contains neutron data.""" + return DetectorData[RunType](detector) def get_monitor_data( @@ -66,7 +94,7 @@ def get_monitor_data( def data_to_tof( - da: RawDetectorData[ScatteringRunType], + da: DetectorData[ScatteringRunType], ) -> TofData[ScatteringRunType]: """Dummy conversion of data to time-of-flight data. The data already has a time-of-flight coordinate.""" @@ -138,6 +166,7 @@ def get_detector_ids_from_sample_run(data: TofData[SampleRun]) -> DetectorIDs: providers = ( + assemble_detector_data, get_detector_data, get_detector_ids_from_sample_run, get_monitor_data, diff --git a/src/ess/isissans/sans2d.py b/src/ess/isissans/sans2d.py index e0598b03..42eea835 100644 --- a/src/ess/isissans/sans2d.py +++ b/src/ess/isissans/sans2d.py @@ -5,7 +5,7 @@ import sciline import scipp as sc from ess.sans import providers as sans_providers -from ess.sans.types import MaskedData, SampleRun, ScatteringRunType, TofData +from ess.sans.types import DetectorMasks, RawDetector, SampleRun from .data import load_tutorial_direct_beam, load_tutorial_run from .general import default_parameters @@ -22,15 +22,19 @@ """Sample holder mask""" -def detector_edge_mask(sample: TofData[SampleRun]) -> DetectorEdgeMask: +# It may make more sense to depend on CalibratedDetector here, but the current +# x and y limits are before setting the beam center, so we use RawDetector +def detector_edge_mask(sample: RawDetector[SampleRun]) -> DetectorEdgeMask: mask_edges = ( sc.abs(sample.coords['position'].fields.x) > sc.scalar(0.48, unit='m') ) | (sc.abs(sample.coords['position'].fields.y) > sc.scalar(0.45, unit='m')) return DetectorEdgeMask(mask_edges) +# It may make more sense to depend on CalibratedDetector here, but the current +# x and y limits are before setting the beam center, so we use RawDetector def sample_holder_mask( - sample: TofData[SampleRun], low_counts_threshold: LowCountThreshold + sample: RawDetector[SampleRun], low_counts_threshold: LowCountThreshold ) -> SampleHolderMask: summed = sample.hist() holder_mask = ( @@ -43,31 +47,31 @@ def sample_holder_mask( return SampleHolderMask(holder_mask) -def mask_detectors( - da: TofData[ScatteringRunType], - edge_mask: DetectorEdgeMask, - holder_mask: SampleHolderMask, -) -> MaskedData[ScatteringRunType]: - """Apply pixel-specific masks to raw data. +def to_detector_masks( + edge_mask: DetectorEdgeMask, holder_mask: SampleHolderMask +) -> DetectorMasks: + """Gather detector masks. + + Unlike :py:func:`ess.sans.masking.to_detector_mask`, this function does not rely on + mapping using a table of mask filenames. Instead it directly returns a dictionary + if multiple masks. Parameters ---------- - da: - Raw data. edge_mask: Mask for detector edges. holder_mask: Mask for sample holder. """ - da = da.copy(deep=False) + masks = {} if edge_mask is not None: - da.masks['edges'] = edge_mask + masks['edges'] = edge_mask if holder_mask is not None: - da.masks['holder_mask'] = holder_mask - return MaskedData[ScatteringRunType](da) + masks['holder_mask'] = holder_mask + return DetectorMasks(masks) -providers = (detector_edge_mask, sample_holder_mask, mask_detectors) +providers = (detector_edge_mask, sample_holder_mask, to_detector_masks) def Sans2dWorkflow() -> sciline.Pipeline: diff --git a/src/ess/loki/general.py b/src/ess/loki/general.py index 3ddb18d2..842f3b69 100644 --- a/src/ess/loki/general.py +++ b/src/ess/loki/general.py @@ -11,7 +11,10 @@ from ..sans.common import gravity_vector from ..sans.types import ( + BeamCenter, + CalibratedDetector, CorrectForGravity, + DetectorData, DetectorEventData, DetectorPixelShape, DimsToKeep, @@ -25,7 +28,6 @@ NonBackgroundWavelengthRange, PixelShapePath, RawDetector, - RawDetectorData, RawMonitor, RawMonitorData, RawSample, @@ -113,61 +115,60 @@ def get_detector_data( return RawDetector[ScatteringRunType](da) +def calibrate_detector( + detector: RawDetector[ScatteringRunType], + beam_center: BeamCenter, + source_position: SourcePosition[ScatteringRunType], + sample_position: SamplePosition[ScatteringRunType], +) -> CalibratedDetector[ScatteringRunType]: + return CalibratedDetector[ScatteringRunType]( + detector.assign_coords( + position=detector.coords['position'] - beam_center, + source_position=source_position, + sample_position=sample_position, + gravity=gravity_vector(), + ) + ) + + def get_monitor_data( monitor: NeXusMonitor[RunType, MonitorType], + source_position: SourcePosition[RunType], ) -> RawMonitor[RunType, MonitorType]: return RawMonitor[RunType, MonitorType]( - nexus.extract_monitor_data(monitor).assign_coords(position=monitor['position']) + nexus.extract_monitor_data(monitor).assign_coords( + position=monitor['position'], source_position=source_position + ) ) -def _add_variances_and_coordinates( - da: sc.DataArray, - source_position: sc.Variable, - sample_position: sc.Variable | None = None, -) -> sc.DataArray: +def _add_variances(da: sc.DataArray) -> sc.DataArray: out = da.copy(deep=False) if out.bins is not None: content = out.bins.constituents['data'] if content.variances is None: content.variances = content.values - # Sample position is not needed in the case of a monitor. - if sample_position is not None: - out.coords['sample_position'] = sample_position - out.coords['source_position'] = source_position - out.coords['gravity'] = gravity_vector() return out def assemble_detector_data( - detector_data: RawDetector[ScatteringRunType], + detector: CalibratedDetector[ScatteringRunType], event_data: DetectorEventData[ScatteringRunType], - source_position: SourcePosition[ScatteringRunType], - sample_position: SamplePosition[ScatteringRunType], -) -> RawDetectorData[ScatteringRunType]: +) -> DetectorData[ScatteringRunType]: grouped = nexus.group_event_data( - event_data=event_data, detector_number=detector_data.coords['detector_number'] - ) - detector_data.data = grouped.data - return RawDetectorData[ScatteringRunType]( - _add_variances_and_coordinates( - da=detector_data, - source_position=source_position, - sample_position=sample_position, - ) + event_data=event_data, detector_number=detector.coords['detector_number'] ) + detector.data = grouped.data + return DetectorData[ScatteringRunType](_add_variances(da=detector)) def assemble_monitor_data( monitor_data: RawMonitor[RunType, MonitorType], event_data: MonitorEventData[RunType, MonitorType], - source_position: SourcePosition[RunType], ) -> RawMonitorData[RunType, MonitorType]: meta = monitor_data.drop_coords('event_time_zero') da = event_data.assign_coords(meta.coords).assign_masks(meta.masks) - return RawMonitorData[RunType, MonitorType]( - _add_variances_and_coordinates(da=da, source_position=source_position) - ) + return RawMonitorData[RunType, MonitorType](_add_variances(da=da)) def _convert_to_tof(da: sc.DataArray) -> sc.DataArray: @@ -178,7 +179,7 @@ def _convert_to_tof(da: sc.DataArray) -> sc.DataArray: def data_to_tof( - da: RawDetectorData[ScatteringRunType], + da: DetectorData[ScatteringRunType], ) -> TofData[ScatteringRunType]: return TofData[ScatteringRunType](_convert_to_tof(da)) @@ -206,6 +207,7 @@ def detector_lab_frame_transform( providers = ( detector_pixel_shape, detector_lab_frame_transform, + calibrate_detector, get_detector_data, get_monitor_data, get_sample_position, diff --git a/src/ess/loki/io.py b/src/ess/loki/io.py index dd0e2699..39cb8c54 100644 --- a/src/ess/loki/io.py +++ b/src/ess/loki/io.py @@ -51,8 +51,7 @@ def load_nexus_detector( def load_nexus_monitor( - file_path: Filename[RunType], - monitor_name: NeXusMonitorName[MonitorType], + file_path: Filename[RunType], monitor_name: NeXusMonitorName[MonitorType] ) -> NeXusMonitor[RunType, MonitorType]: return NeXusMonitor[RunType, MonitorType]( nexus.load_monitor( @@ -71,8 +70,7 @@ def load_detector_event_data( def load_monitor_event_data( - file_path: Filename[RunType], - monitor_name: NeXusMonitorName[MonitorType], + file_path: Filename[RunType], monitor_name: NeXusMonitorName[MonitorType] ) -> MonitorEventData[RunType, MonitorType]: da = nexus.load_event_data(file_path=file_path, component_name=monitor_name) return MonitorEventData[RunType, MonitorType](da) diff --git a/src/ess/sans/beam_center_finder.py b/src/ess/sans/beam_center_finder.py index 14f086d3..367cb8b7 100644 --- a/src/ess/sans/beam_center_finder.py +++ b/src/ess/sans/beam_center_finder.py @@ -2,7 +2,6 @@ # Copyright (c) 2023 Scipp contributors (https://github.com/scipp) import uuid -from typing import NewType import numpy as np import sciline @@ -10,35 +9,19 @@ from ess.reduce.uncertainty import UncertaintyBroadcastMode from scipp.core import concepts -from .conversions import ( - ElasticCoordTransformGraph, - calibrate_positions, - compute_Q, - detector_to_wavelength, - mask_wavelength, -) -from .i_of_q import bin_in_q +from .conversions import ElasticCoordTransformGraph from .logging import get_logger -from .normalization import ( - iofq_denominator, - normalize_q, - process_wavelength_bands, - solid_angle, -) from .types import ( BeamCenter, - CalibratedMaskedData, - DetectorPixelShape, DimsToKeep, IofQ, - LabFrameTransform, MaskedData, NormWavelengthTerm, QBins, + RawDetector, ReturnEvents, SampleRun, WavelengthBands, - WavelengthBins, WavelengthMask, ) @@ -51,10 +34,7 @@ def _xy_extrema(pos: sc.Variable) -> sc.Variable: return sc.concat([x_min, x_max, y_min, y_max], dim='extremes') -def beam_center_from_center_of_mass( - data: MaskedData[SampleRun], - graph: ElasticCoordTransformGraph, -) -> BeamCenter: +def beam_center_from_center_of_mass(workflow: sciline.Pipeline) -> BeamCenter: """ Estimate the beam center via the center-of-mass of the data counts. @@ -67,16 +47,22 @@ def beam_center_from_center_of_mass( Parameters ---------- - data: - The data to find the beam center of. - graph: - Coordinate transformation graph for elastic SANS. + workflow: + The reduction workflow to compute MaskedData[SampleRun]. Returns ------- : The beam center position as a vector. """ + try: + beam_center = workflow.compute(BeamCenter) + except sciline.UnsatisfiedRequirement: + beam_center = sc.vector([0.0, 0.0, 0.0], unit='m') + workflow = workflow.copy() + workflow[BeamCenter] = beam_center + data = workflow.compute(MaskedData[SampleRun]) + graph = workflow.compute(ElasticCoordTransformGraph) dims_to_sum = set(data.dims) - set(data.coords['position'].dims) if dims_to_sum: @@ -114,7 +100,7 @@ def beam_center_from_center_of_mass( n_beam = incident_beam / sc.norm(incident_beam) com_shift = com - sc.dot(com, n_beam) * n_beam xy = [com_shift.fields.x.value, com_shift.fields.y.value] - return _offsets_to_vector(data=summed, xy=xy, graph=graph) + return beam_center + _offsets_to_vector(data=summed, xy=xy, graph=graph) def _offsets_to_vector(data: sc.DataArray, xy: list[float], graph: dict) -> sc.Variable: @@ -134,13 +120,9 @@ def _offsets_to_vector(data: sc.DataArray, xy: list[float], graph: dict) -> sc.V def _iofq_in_quadrants( xy: list[float], - data: sc.DataArray, + workflow: sciline.Pipeline, + detector: sc.DataArray, norm: sc.DataArray, - graph: dict, - q_bins: int | sc.Variable, - wavelength_bins: sc.Variable, - transform: sc.Variable, - pixel_shape: sc.DataGroup, ) -> dict[str, sc.DataArray]: """ Compute the intensity as a function of Q inside 4 quadrants in Phi. @@ -149,16 +131,10 @@ def _iofq_in_quadrants( ---------- xy: The x,y offsets in the plane normal to the beam. - data: - The sample data. + detector: + The raw detector. norm: The denominator data for normalization. - graph: - Coordinate transformation graph. - q_bins: - Bin edges for Q. - wavelength_bins: - The binning in wavelength to use for computing the intensity as a function of Q. Returns ------- @@ -171,33 +147,10 @@ def _iofq_in_quadrants( phi_bins = sc.linspace('phi', -pi, pi, 5, unit='rad') quadrants = ['south-west', 'south-east', 'north-east', 'north-west'] - providers = [ - compute_Q, - bin_in_q, - normalize_q, - iofq_denominator, - mask_wavelength, - detector_to_wavelength, - solid_angle, - calibrate_positions, - process_wavelength_bands, - ] - params = {} - params[UncertaintyBroadcastMode] = UncertaintyBroadcastMode.upper_bound - params[ReturnEvents] = False - params[WavelengthBins] = wavelength_bins - params[QBins] = q_bins - params[DetectorPixelShape[SampleRun]] = pixel_shape - params[LabFrameTransform[SampleRun]] = transform - params[ElasticCoordTransformGraph] = graph - params[BeamCenter] = _offsets_to_vector(data=data, xy=xy, graph=graph) - params[DimsToKeep] = () - params[WavelengthMask] = None - params[WavelengthBands] = None - - pipeline = sciline.Pipeline(providers, params=params) - pipeline[MaskedData[SampleRun]] = data - calibrated = pipeline.compute(CalibratedMaskedData[SampleRun]) + graph = workflow.compute(ElasticCoordTransformGraph) + workflow = workflow.copy() + workflow[BeamCenter] = _offsets_to_vector(data=detector, xy=xy, graph=graph) + calibrated = workflow.compute(MaskedData[SampleRun]) with_phi = calibrated.transform_coords( 'phi', graph=graph, keep_intermediate=False, keep_inputs=False ) @@ -215,11 +168,13 @@ def _iofq_in_quadrants( for i, quad in enumerate(quadrants): # Select pixels based on phi sel = (phi >= phi_bins[i]) & (phi < phi_bins[i + 1]) - pipeline[MaskedData[SampleRun]] = data[sel] - pipeline[NormWavelengthTerm[SampleRun]] = ( + workflow[RawDetector[SampleRun]] = detector[sel] + # MaskedData would be computed automatically, but we did it above already + workflow[MaskedData[SampleRun]] = calibrated[sel] + workflow[NormWavelengthTerm[SampleRun]] = ( norm if norm.dims == ('wavelength',) else norm[sel] ) - out[quad] = pipeline.compute(IofQ[SampleRun]) + out[quad] = workflow.compute(IofQ[SampleRun]) return out @@ -291,26 +246,12 @@ def _cost(xy: list[float], *args) -> float: return out -BeamCenterFinderQBins = NewType('BeamCenterFinderQBins', sc.Variable) -"""Q binning used for the beam center finder""" - -BeamCenterFinderTolerance = NewType('BeamCenterFinderTolerance', float | None) -"""Tolerance used for the beam center finder""" - -BeamCenterFinderMinimizer = NewType('BeamCenterFinderMinimizer', str | None) -"""Minimizer used for the beam center finder""" - - def beam_center_from_iofq( - data: MaskedData[SampleRun], - graph: ElasticCoordTransformGraph, - wavelength_bins: WavelengthBins, - norm: NormWavelengthTerm[SampleRun], - q_bins: BeamCenterFinderQBins, - transform: LabFrameTransform[SampleRun], - pixel_shape: DetectorPixelShape[SampleRun], - minimizer: BeamCenterFinderMinimizer, - tolerance: BeamCenterFinderTolerance, + *, + workflow: sciline.Pipeline, + q_bins: int | sc.Variable, + minimizer: str | None = None, + tolerance: float | None = None, ) -> BeamCenter: """ Find the beam center of a SANS scattering pattern using an I(Q) calculation. @@ -327,12 +268,8 @@ def beam_center_from_iofq( Parameters ---------- - data: - The DataArray containing the detector data. - graph: - Coordinate transformation graph for elastic SANS. - wavelength_bins: - The binning in the wavelength dimension to be used. + workflow: + The reduction workflow to compute I(Q). q_bins: The binning in the Q dimension to be used. minimizer: @@ -420,17 +357,39 @@ def beam_center_from_iofq( logger.info('Using minimizer: %s', minimizer) logger.info('Using tolerance: %s', tolerance) + keys = ( + RawDetector[SampleRun], + MaskedData[SampleRun], + NormWavelengthTerm[SampleRun], + ElasticCoordTransformGraph, + ) + results = workflow.compute(keys) + detector = results[RawDetector[SampleRun]] + data = results[MaskedData[SampleRun]] + norm = results[NormWavelengthTerm[SampleRun]] + graph = results[ElasticCoordTransformGraph] + # Flatten positions dim which is required during the iterations for slicing with a # boolean mask - pos_dims = data.coords['position'].dims + pos_dims = detector.coords['position'].dims new_dim = uuid.uuid4().hex - data = data.flatten(dims=pos_dims, to=new_dim) + detector = detector.flatten(dims=pos_dims, to=new_dim) dims_to_flatten = [dim for dim in norm.dims if dim in pos_dims] if dims_to_flatten: norm = norm.flatten(dims=dims_to_flatten, to=new_dim) + workflow = workflow.copy() + # Avoid reloading the detector + workflow[RawDetector[SampleRun]] = detector + workflow[UncertaintyBroadcastMode] = UncertaintyBroadcastMode.upper_bound + workflow[ReturnEvents] = False + workflow[DimsToKeep] = () + workflow[WavelengthMask] = None + workflow[WavelengthBands] = None + workflow[QBins] = q_bins + # Use center of mass to get initial guess for beam center - com_shift = beam_center_from_center_of_mass(data, graph) + com_shift = beam_center_from_center_of_mass(workflow) logger.info('Initial guess for beam center: %s', com_shift) coords = data.transform_coords( @@ -445,7 +404,7 @@ def beam_center_from_iofq( res = minimize( _cost, x0=[com_shift.fields.x.value, com_shift.fields.y.value], - args=(data, norm, graph, q_bins, wavelength_bins, transform, pixel_shape), + args=(workflow, detector, norm), bounds=bounds, method=minimizer, tol=tolerance, diff --git a/src/ess/sans/conversions.py b/src/ess/sans/conversions.py index 0f71185b..1e4bb823 100644 --- a/src/ess/sans/conversions.py +++ b/src/ess/sans/conversions.py @@ -11,8 +11,6 @@ from .common import mask_range from .types import ( - BeamCenter, - CalibratedMaskedData, CleanQ, CleanQxy, CleanWavelength, @@ -153,24 +151,11 @@ def monitor_to_wavelength( ) -def calibrate_positions( - detector: MaskedData[ScatteringRunType], beam_center: BeamCenter -) -> CalibratedMaskedData[ScatteringRunType]: - """ - Calibrate pixel positions. - - Currently the only applied calibration is the beam-center offset. - """ - detector = detector.copy(deep=False) - detector.coords['position'] = detector.coords['position'] - beam_center - return detector - - # TODO This demonstrates a problem: Transforming to wavelength should be possible # for RawData, MaskedData, ... no reason to restrict necessarily. # Would we be fine with just choosing on option, or will this get in the way for users? def detector_to_wavelength( - detector: CalibratedMaskedData[ScatteringRunType], + detector: MaskedData[ScatteringRunType], graph: ElasticCoordTransformGraph, ) -> CleanWavelength[ScatteringRunType, Numerator]: return CleanWavelength[ScatteringRunType, Numerator]( @@ -227,7 +212,6 @@ def compute_Qxy( providers = ( sans_elastic, sans_monitor, - calibrate_positions, monitor_to_wavelength, detector_to_wavelength, mask_wavelength, diff --git a/src/ess/sans/normalization.py b/src/ess/sans/normalization.py index ca164dc6..17a041d8 100644 --- a/src/ess/sans/normalization.py +++ b/src/ess/sans/normalization.py @@ -5,19 +5,21 @@ from scipp.core import concepts from .types import ( - CalibratedMaskedData, + CalibratedDetector, CleanDirectBeam, CleanMonitor, CleanSummedQ, CleanSummedQxy, CleanWavelength, Denominator, + DetectorMasks, DetectorPixelShape, EmptyBeamRun, Incident, IofQ, IofQxy, LabFrameTransform, + MaskedSolidAngle, NormWavelengthTerm, Numerator, ProcessedWavelengthBands, @@ -33,7 +35,7 @@ def solid_angle( - data: CalibratedMaskedData[ScatteringRunType], + data: CalibratedDetector[ScatteringRunType], pixel_shape: DetectorPixelShape[ScatteringRunType], transform: LabFrameTransform[ScatteringRunType], ) -> SolidAngle[ScatteringRunType]: @@ -81,6 +83,13 @@ def solid_angle( ) +def mask_solid_angle( + solid_angle: SolidAngle[ScatteringRunType], + masks: DetectorMasks, +) -> MaskedSolidAngle[ScatteringRunType]: + return MaskedSolidAngle[ScatteringRunType](solid_angle.assign_masks(masks)) + + def _approximate_solid_angle_for_cylinder_shaped_pixel_of_detector( pixel_position: sc.Variable, cylinder_axis: sc.Variable, @@ -219,7 +228,7 @@ def iofq_norm_wavelength_term( def iofq_denominator( wavelength_term: NormWavelengthTerm[ScatteringRunType], - solid_angle: SolidAngle[ScatteringRunType], + solid_angle: MaskedSolidAngle[ScatteringRunType], uncertainties: UncertaintyBroadcastMode, ) -> CleanWavelength[ScatteringRunType, Denominator]: """ @@ -452,4 +461,5 @@ def normalize_qxy( normalize_qxy, process_wavelength_bands, solid_angle, + mask_solid_angle, ) diff --git a/src/ess/sans/types.py b/src/ess/sans/types.py index 5be48ded..19f5528f 100644 --- a/src/ess/sans/types.py +++ b/src/ess/sans/types.py @@ -217,6 +217,10 @@ class SolidAngle(sciline.Scope[ScatteringRunType, sc.DataArray], sc.DataArray): """Solid angle of detector pixels seen from sample position""" +class MaskedSolidAngle(sciline.Scope[ScatteringRunType, sc.DataArray], sc.DataArray): + """Same as :py:class:`SolidAngle`, but with pixel masks applied""" + + class NeXusDetector(sciline.Scope[RunType, sc.DataGroup], sc.DataGroup): """Detector data, loaded from a NeXus file, containing not only neutron events but also pixel shape information, transformations, ...""" @@ -240,13 +244,15 @@ class MonitorEventData( class RawDetector(sciline.Scope[ScatteringRunType, sc.DataArray], sc.DataArray): - """Raw detector data component extracted from :py:class:`NeXusDetector`""" + """Raw detector component extracted from :py:class:`NeXusDetector`""" -class RawDetectorData(sciline.Scope[ScatteringRunType, sc.DataArray], sc.DataArray): - """Raw event data where variances and necessary coordinates - (e.g. sample and source position) have been added, and where optionally some - user configuration was applied to some of the coordinates.""" +class CalibratedDetector(sciline.Scope[ScatteringRunType, sc.DataArray], sc.DataArray): + """Calibrated version of raw detector""" + + +class DetectorData(sciline.Scope[ScatteringRunType, sc.DataArray], sc.DataArray): + """Calibrated detector with added raw event data""" class TofData(sciline.Scope[ScatteringRunType, sc.DataArray], sc.DataArray): @@ -266,12 +272,6 @@ class MaskedData(sciline.Scope[ScatteringRunType, sc.DataArray], sc.DataArray): """Raw data with pixel-specific masks applied""" -class CalibratedMaskedData( - sciline.Scope[ScatteringRunType, sc.DataArray], sc.DataArray -): - """Raw data with pixel-specific masks applied and calibrated pixel positions""" - - class NormWavelengthTerm(sciline.Scope[ScatteringRunType, sc.DataArray], sc.DataArray): """Normalization term (numerator) for IofQ before scaling with solid-angle.""" diff --git a/tests/isissans/sans2d_reduction_test.py b/tests/isissans/sans2d_reduction_test.py index 5f03b47a..1b81a04c 100644 --- a/tests/isissans/sans2d_reduction_test.py +++ b/tests/isissans/sans2d_reduction_test.py @@ -1,7 +1,5 @@ # SPDX-License-Identifier: BSD-3-Clause # Copyright (c) 2023 Scipp contributors (https://github.com/scipp) -from collections.abc import Callable - import pytest import sciline import scipp as sc @@ -71,6 +69,7 @@ def make_params() -> dict: params[UncertaintyBroadcastMode] = UncertaintyBroadcastMode.upper_bound params[ReturnEvents] = False params[DimsToKeep] = () + params[BeamCenter] = sc.vector([0, 0, 0], unit='m') return params @@ -86,7 +85,6 @@ def sans2d_providers(): isis.data.load_tutorial_run, isis.data.transmission_from_background_run, isis.data.transmission_from_sample_run, - sans.beam_center_finder.beam_center_from_center_of_mass, ) ) @@ -121,6 +119,7 @@ def test_pipeline_can_compute_background_subtracted_IofQ_in_wavelength_bands(): def test_pipeline_wavelength_bands_is_optional(): params = make_params() pipeline = sciline.Pipeline(sans2d_providers(), params=params) + pipeline[BeamCenter] = sans.beam_center_from_center_of_mass(pipeline) noband = pipeline.compute(BackgroundSubtractedIofQ) assert pipeline.compute(WavelengthBands) is None band = sc.linspace('wavelength', 2.0, 16.0, num=2, unit='angstrom') @@ -134,6 +133,7 @@ def test_workflow_is_deterministic(): params = make_params() params[UncertaintyBroadcastMode] = UncertaintyBroadcastMode.drop pipeline = sciline.Pipeline(sans2d_providers(), params=params) + pipeline[BeamCenter] = sans.beam_center_from_center_of_mass(pipeline) # This is Sciline's default scheduler, but we want to be explicit here scheduler = sciline.scheduler.DaskScheduler() graph = pipeline.get(IofQ[SampleRun], scheduler=scheduler) @@ -179,15 +179,6 @@ def test_pipeline_can_compute_intermediate_results(): assert result.dims == ('spectrum',) -# TODO See scipp/sciline#57 for plans on a builtin way to do this -def as_dict(funcs: list[Callable[..., type]]) -> dict: - from typing import get_type_hints - - return dict( - zip([get_type_hints(func)['return'] for func in funcs], funcs, strict=True) - ) - - def pixel_dependent_direct_beam( filename: DirectBeamFilename, shape: RawDetector[SampleRun] ) -> DirectBeam: @@ -203,9 +194,9 @@ def pixel_dependent_direct_beam( def test_pixel_dependent_direct_beam_is_supported(uncertainties): params = make_params() params[UncertaintyBroadcastMode] = uncertainties - providers = as_dict(sans2d_providers()) - providers[DirectBeam] = pixel_dependent_direct_beam - pipeline = sciline.Pipeline(providers.values(), params=params) + pipeline = sciline.Pipeline(sans2d_providers(), params=params) + pipeline.insert(pixel_dependent_direct_beam) + pipeline[BeamCenter] = sc.vector([0, 0, 0], unit='m') result = pipeline.compute(BackgroundSubtractedIofQ) assert result.dims == ('Q',) @@ -217,27 +208,31 @@ def test_beam_center_from_center_of_mass_is_close_to_verified_result(): params = make_params() providers = sans2d_providers() pipeline = sciline.Pipeline(providers, params=params) - center = pipeline.compute(BeamCenter) + center = sans.beam_center_from_center_of_mass(pipeline) # This is the result obtained from Mantid, using the full IofQ # calculation. The difference is about 3 mm in X or Y, probably due to a bias # introduced by the sample holder, which the center-of-mass approach cannot ignore. assert sc.allclose(center, MANTID_BEAM_CENTER, atol=sc.scalar(3e-3, unit='m')) +def test_beam_center_from_center_of_mass_independent_of_set_beam_center(): + params = make_params() + providers = sans2d_providers() + pipeline = sciline.Pipeline(providers, params=params) + pipeline[BeamCenter] = sc.vector([0.1, -0.1, 0], unit='m') + center = sans.beam_center_from_center_of_mass(pipeline) + assert sc.allclose(center, MANTID_BEAM_CENTER, atol=sc.scalar(3e-3, unit='m')) + + def test_beam_center_finder_without_direct_beam_reproduces_verified_result(): params = make_params() - params[sans.beam_center_finder.BeamCenterFinderQBins] = sc.linspace( - 'Q', 0.02, 0.3, 71, unit='1/angstrom' - ) del params[DirectBeamFilename] providers = sans2d_providers() - providers.remove(sans.beam_center_finder.beam_center_from_center_of_mass) - providers.append(sans.beam_center_finder.beam_center_from_iofq) pipeline = sciline.Pipeline(providers, params=params) - pipeline[sans.beam_center_finder.BeamCenterFinderMinimizer] = None - pipeline[sans.beam_center_finder.BeamCenterFinderTolerance] = None pipeline[DirectBeam] = None - center = pipeline.compute(BeamCenter) + center = sans.beam_center_finder.beam_center_from_iofq( + workflow=pipeline, q_bins=sc.linspace('Q', 0.02, 0.3, 71, unit='1/angstrom') + ) assert sc.allclose(center, MANTID_BEAM_CENTER, atol=sc.scalar(2e-3, unit='m')) @@ -250,52 +245,54 @@ def low_counts_mask( params = make_params() params[sans2d.LowCountThreshold] = sc.scalar(80.0, unit='counts') - params[sans.beam_center_finder.BeamCenterFinderQBins] = sc.linspace( - 'Q', 0.02, 0.3, 71, unit='1/angstrom' - ) del params[DirectBeamFilename] providers = sans2d_providers() providers.remove(sans2d.sample_holder_mask) - providers.remove(sans.beam_center_finder.beam_center_from_center_of_mass) - providers.append(sans.beam_center_finder.beam_center_from_iofq) providers.append(low_counts_mask) pipeline = sciline.Pipeline(providers, params=params) - pipeline[sans.beam_center_finder.BeamCenterFinderMinimizer] = None - pipeline[sans.beam_center_finder.BeamCenterFinderTolerance] = None pipeline[DirectBeam] = None - center = pipeline.compute(BeamCenter) + q_bins = sc.linspace('Q', 0.02, 0.3, 71, unit='1/angstrom') + center = sans.beam_center_finder.beam_center_from_iofq( + workflow=pipeline, q_bins=q_bins + ) assert sc.allclose(center, MANTID_BEAM_CENTER, atol=sc.scalar(5e-4, unit='m')) def test_beam_center_finder_works_with_direct_beam(): params = make_params() - params[sans.beam_center_finder.BeamCenterFinderQBins] = sc.linspace( - 'Q', 0.02, 0.3, 71, unit='1/angstrom' + providers = sans2d_providers() + pipeline = sciline.Pipeline(providers, params=params) + q_bins = sc.linspace('Q', 0.02, 0.3, 71, unit='1/angstrom') + center_with_direct_beam = sans.beam_center_finder.beam_center_from_iofq( + workflow=pipeline, q_bins=q_bins + ) + assert sc.allclose( + center_with_direct_beam, MANTID_BEAM_CENTER, atol=sc.scalar(2e-3, unit='m') ) + + +def test_beam_center_finder_independent_of_set_beam_center(): + params = make_params() providers = sans2d_providers() - providers.remove(sans.beam_center_finder.beam_center_from_center_of_mass) - providers.append(sans.beam_center_finder.beam_center_from_iofq) pipeline = sciline.Pipeline(providers, params=params) - pipeline[sans.beam_center_finder.BeamCenterFinderMinimizer] = None - pipeline[sans.beam_center_finder.BeamCenterFinderTolerance] = None - center_with_direct_beam = pipeline.compute(BeamCenter) + pipeline[BeamCenter] = sc.vector([0.1, -0.1, 0], unit='m') + q_bins = sc.linspace('Q', 0.02, 0.3, 71, unit='1/angstrom') + center_with_direct_beam = sans.beam_center_finder.beam_center_from_iofq( + workflow=pipeline, q_bins=q_bins + ) assert sc.allclose( center_with_direct_beam, MANTID_BEAM_CENTER, atol=sc.scalar(2e-3, unit='m') ) def test_beam_center_finder_works_with_pixel_dependent_direct_beam(): + q_bins = sc.linspace('Q', 0.02, 0.3, 71, unit='1/angstrom') params = make_params() - params[sans.beam_center_finder.BeamCenterFinderQBins] = sc.linspace( - 'Q', 0.02, 0.3, 71, unit='1/angstrom' - ) providers = sans2d_providers() - providers.remove(sans.beam_center_finder.beam_center_from_center_of_mass) - providers.append(sans.beam_center_finder.beam_center_from_iofq) pipeline = sciline.Pipeline(providers, params=params) - pipeline[sans.beam_center_finder.BeamCenterFinderMinimizer] = None - pipeline[sans.beam_center_finder.BeamCenterFinderTolerance] = None - center_pixel_independent_direct_beam = pipeline.compute(BeamCenter) + center_pixel_independent_direct_beam = ( + sans.beam_center_finder.beam_center_from_iofq(workflow=pipeline, q_bins=q_bins) + ) direct_beam = pipeline.compute(DirectBeam) pixel_dependent_direct_beam = direct_beam.broadcast( @@ -306,14 +303,12 @@ def test_beam_center_finder_works_with_pixel_dependent_direct_beam(): ).copy() providers = sans2d_providers() - providers.remove(sans.beam_center_finder.beam_center_from_center_of_mass) - providers.append(sans.beam_center_finder.beam_center_from_iofq) pipeline = sciline.Pipeline(providers, params=params) - pipeline[sans.beam_center_finder.BeamCenterFinderMinimizer] = None - pipeline[sans.beam_center_finder.BeamCenterFinderTolerance] = None pipeline[DirectBeam] = pixel_dependent_direct_beam - center = pipeline.compute(BeamCenter) + center = sans.beam_center_finder.beam_center_from_iofq( + workflow=pipeline, q_bins=q_bins + ) assert sc.identical(center, center_pixel_independent_direct_beam) diff --git a/tests/isissans/zoom_reduction_test.py b/tests/isissans/zoom_reduction_test.py index 8bc6b980..24b12a1f 100644 --- a/tests/isissans/zoom_reduction_test.py +++ b/tests/isissans/zoom_reduction_test.py @@ -5,6 +5,7 @@ from ess import isissans as isis from ess import sans from ess.sans.types import ( + BeamCenter, CorrectForGravity, Filename, Incident, @@ -63,13 +64,13 @@ def zoom_providers(): isis.data.load_tutorial_run, isis.data.transmission_from_background_run, isis.data.transmission_from_sample_run, - sans.beam_center_finder.beam_center_from_center_of_mass, ) ) def test_can_create_pipeline(): pipeline = sciline.Pipeline(zoom_providers(), params=make_params()) + pipeline[BeamCenter] = sc.vector([0, 0, 0], unit='m') pipeline = sans.with_pixel_mask_filenames( pipeline, isis.data.zoom_tutorial_mask_filenames() ) @@ -78,6 +79,7 @@ def test_can_create_pipeline(): def test_pipeline_can_compute_IofQ(): pipeline = sciline.Pipeline(zoom_providers(), params=make_params()) + pipeline[BeamCenter] = sc.vector([0, 0, 0], unit='m') pipeline = sans.with_pixel_mask_filenames( pipeline, isis.data.zoom_tutorial_mask_filenames() ) @@ -87,6 +89,7 @@ def test_pipeline_can_compute_IofQ(): def test_pipeline_can_compute_IofQxQy(): pipeline = sciline.Pipeline(zoom_providers(), params=make_params()) + pipeline[BeamCenter] = sc.vector([0, 0, 0], unit='m') pipeline = sans.with_pixel_mask_filenames( pipeline, isis.data.zoom_tutorial_mask_filenames() ) diff --git a/tests/loki/common.py b/tests/loki/common.py index a066eba5..092b490e 100644 --- a/tests/loki/common.py +++ b/tests/loki/common.py @@ -55,7 +55,7 @@ def make_params(no_masks: bool = True) -> dict: return params -def loki_providers_no_beam_center_finder() -> list[Callable]: +def loki_providers() -> list[Callable]: from ess.isissans.io import read_xml_detector_masking return list( @@ -66,10 +66,3 @@ def loki_providers_no_beam_center_finder() -> list[Callable]: loki.io.dummy_load_sample, ) ) - - -def loki_providers() -> list[Callable]: - return [ - *loki_providers_no_beam_center_finder(), - sans.beam_center_finder.beam_center_from_center_of_mass, - ] diff --git a/tests/loki/directbeam_test.py b/tests/loki/directbeam_test.py index e4b3295f..a9f6045d 100644 --- a/tests/loki/directbeam_test.py +++ b/tests/loki/directbeam_test.py @@ -6,7 +6,13 @@ import sciline import scipp as sc from ess import loki, sans -from ess.sans.types import DimsToKeep, QBins, WavelengthBands, WavelengthBins +from ess.sans.types import ( + BeamCenter, + DimsToKeep, + QBins, + WavelengthBands, + WavelengthBins, +) from scipp.scipy.interpolate import interp1d sys.path.insert(0, str(Path(__file__).resolve().parent)) @@ -30,6 +36,7 @@ def test_can_compute_direct_beam_for_all_pixels(): ) providers = loki_providers() pipeline = sciline.Pipeline(providers, params=params) + pipeline[BeamCenter] = sc.vector([0, 0, 0], unit='m') I0 = _get_I0(qbins=params[QBins]) results = sans.direct_beam(workflow=pipeline, I0=I0, niter=4) @@ -59,6 +66,7 @@ def test_can_compute_direct_beam_with_overlapping_wavelength_bands(): providers = loki_providers() pipeline = sciline.Pipeline(providers, params=params) + pipeline[BeamCenter] = sc.vector([0, 0, 0], unit='m') I0 = _get_I0(qbins=params[QBins]) results = sans.direct_beam(workflow=pipeline, I0=I0, niter=4) @@ -84,6 +92,7 @@ def test_can_compute_direct_beam_per_layer(): params[DimsToKeep] = ['layer'] providers = loki_providers() pipeline = sciline.Pipeline(providers, params=params) + pipeline[BeamCenter] = sc.vector([0, 0, 0], unit='m') I0 = _get_I0(qbins=params[QBins]) results = sans.direct_beam(workflow=pipeline, I0=I0, niter=4) @@ -111,6 +120,7 @@ def test_can_compute_direct_beam_per_layer_and_straw(): params[DimsToKeep] = ('layer', 'straw') providers = loki_providers() pipeline = sciline.Pipeline(providers, params=params) + pipeline[BeamCenter] = sc.vector([0, 0, 0], unit='m') I0 = _get_I0(qbins=params[QBins]) results = sans.direct_beam(workflow=pipeline, I0=I0, niter=4) diff --git a/tests/loki/iofq_test.py b/tests/loki/iofq_test.py index ec7d2d1c..f14f92f3 100644 --- a/tests/loki/iofq_test.py +++ b/tests/loki/iofq_test.py @@ -13,7 +13,6 @@ BackgroundSubtractedIofQ, BackgroundSubtractedIofQxy, BeamCenter, - CalibratedMaskedData, CleanSummedQ, CleanWavelengthMasked, CorrectForGravity, @@ -21,6 +20,7 @@ DimsToKeep, IofQ, IofQxy, + MaskedData, NeXusDetectorName, Numerator, QBins, @@ -37,13 +37,13 @@ sys.path.insert(0, str(Path(__file__).resolve().parent)) from common import ( loki_providers, - loki_providers_no_beam_center_finder, make_params, ) def test_can_create_pipeline(): pipeline = sciline.Pipeline(loki_providers(), params=make_params()) + pipeline[BeamCenter] = sc.vector([0, 0, 0], unit='m') pipeline.get(BackgroundSubtractedIofQ) @@ -52,6 +52,7 @@ def test_can_create_pipeline_with_pixel_masks(): pipeline = sans.with_pixel_mask_filenames( pipeline, loki.data.loki_tutorial_mask_filenames() ) + pipeline[BeamCenter] = sc.vector([0, 0, 0], unit='m') pipeline.get(BackgroundSubtractedIofQ) @@ -67,6 +68,7 @@ def test_pipeline_can_compute_IofQ(uncertainties, qxy: bool): pipeline = sans.with_pixel_mask_filenames( pipeline, loki.data.loki_tutorial_mask_filenames() ) + pipeline[BeamCenter] = sans.beam_center_from_center_of_mass(pipeline) if qxy: result = pipeline.compute(BackgroundSubtractedIofQxy) assert result.dims == ('Qy', 'Qx') @@ -103,6 +105,7 @@ def test_pipeline_can_compute_IofQ_in_event_mode(uncertainties, target): params = make_params() params[UncertaintyBroadcastMode] = uncertainties pipeline = sciline.Pipeline(loki_providers(), params=params) + pipeline[BeamCenter] = sans.beam_center_from_center_of_mass(pipeline) reference = pipeline.compute(target) pipeline[ReturnEvents] = True result = pipeline.compute(target) @@ -141,6 +144,7 @@ def test_pipeline_can_compute_IofQ_in_wavelength_bands(qxy: bool): 11, ) pipeline = sciline.Pipeline(loki_providers(), params=params) + pipeline[BeamCenter] = _compute_beam_center() result = pipeline.compute( BackgroundSubtractedIofQxy if qxy else BackgroundSubtractedIofQ ) @@ -159,6 +163,7 @@ def test_pipeline_can_compute_IofQ_in_overlapping_wavelength_bands(qxy: bool): [edges[:-2], edges[2::]], dim='wavelength' ).transpose() pipeline = sciline.Pipeline(loki_providers(), params=params) + pipeline[BeamCenter] = _compute_beam_center() result = pipeline.compute( BackgroundSubtractedIofQxy if qxy else BackgroundSubtractedIofQ ) @@ -171,6 +176,7 @@ def test_pipeline_can_compute_IofQ_in_layers(qxy: bool): params = make_params() params[DimsToKeep] = ['layer'] pipeline = sciline.Pipeline(loki_providers(), params=params) + pipeline[BeamCenter] = _compute_beam_center() result = pipeline.compute( BackgroundSubtractedIofQxy if qxy else BackgroundSubtractedIofQ ) @@ -180,8 +186,7 @@ def test_pipeline_can_compute_IofQ_in_layers(qxy: bool): def _compute_beam_center(): pipeline = sciline.Pipeline(loki_providers(), params=make_params()) - center = pipeline.compute(BeamCenter) - return center + return sans.beam_center_from_center_of_mass(pipeline) def test_pipeline_can_compute_IofQ_merging_events_from_multiple_runs(): @@ -195,7 +200,7 @@ def test_pipeline_can_compute_IofQ_merging_events_from_multiple_runs(): loki.data.loki_tutorial_background_run_60248(), loki.data.loki_tutorial_background_run_60393(), ] - pipeline = sciline.Pipeline(loki_providers_no_beam_center_finder(), params=params) + pipeline = sciline.Pipeline(loki_providers(), params=params) pipeline[BeamCenter] = _compute_beam_center() pipeline = sans.with_sample_runs(pipeline, runs=sample_runs) @@ -209,7 +214,7 @@ def test_pipeline_can_compute_IofQ_by_bank(): params = make_params() del params[NeXusDetectorName] - pipeline = sciline.Pipeline(loki_providers_no_beam_center_finder(), params=params) + pipeline = sciline.Pipeline(loki_providers(), params=params) pipeline[BeamCenter] = _compute_beam_center() pipeline = sans.with_banks(pipeline, banks=['larmor_detector']) @@ -227,7 +232,7 @@ def test_pipeline_can_compute_IofQ_merging_events_from_multiple_runs_by_bank(): loki.data.loki_tutorial_background_run_60248(), loki.data.loki_tutorial_background_run_60393(), ] - pipeline = sciline.Pipeline(loki_providers_no_beam_center_finder(), params=params) + pipeline = sciline.Pipeline(loki_providers(), params=params) pipeline[BeamCenter] = _compute_beam_center() pipeline = sans.with_sample_runs(pipeline, runs=sample_runs) @@ -248,9 +253,7 @@ def test_pipeline_IofQ_merging_events_yields_consistent_results(): N = 3 params = make_params() center = _compute_beam_center() - pipeline_single = sciline.Pipeline( - loki_providers_no_beam_center_finder(), params=params - ) + pipeline_single = sciline.Pipeline(loki_providers(), params=params) pipeline_single[BeamCenter] = center sample_runs = [loki.data.loki_tutorial_sample_run_60339()] * N @@ -285,7 +288,7 @@ def test_beam_center_from_center_of_mass_is_close_to_verified_result(): pipeline = sans.with_pixel_mask_filenames( pipeline, loki.data.loki_tutorial_mask_filenames() ) - center = pipeline.compute(BeamCenter) + center = sans.beam_center_from_center_of_mass(pipeline) reference = sc.vector([-0.0291487, -0.0181614, 0], unit='m') assert sc.allclose(center, reference) @@ -293,6 +296,7 @@ def test_beam_center_from_center_of_mass_is_close_to_verified_result(): def test_phi_with_gravity(): params = make_params() pipeline = sciline.Pipeline(loki_providers(), params=params) + pipeline[BeamCenter] = _compute_beam_center() pipeline[CorrectForGravity] = False data_no_grav = pipeline.compute( CleanWavelengthMasked[SampleRun, Numerator] @@ -316,7 +320,7 @@ def test_phi_with_gravity(): # Exclude pixels near y=0, since phi with gravity could drop below y=0 and give a # difference of almost 2*pi. y = sc.abs( - pipeline.compute(CalibratedMaskedData[SampleRun]) + pipeline.compute(MaskedData[SampleRun]) .coords['position'] .fields.y.flatten(to='pixel') )