Sans2d Data Reduction#

Introduction#

This notebook gives an overview of how to use the esssans package with Sciline, on the example of the data reduction of a Sans2d experiment. We begin with relevant imports:

[1]:

import scipp as sc
import plopp as pp
from ess import sans
from ess import isissans as isis
import ess.isissans.data  # noqa: F401
from ess.sans.types import *

Create and configure the workflow#

We begin by creating the Sans2d workflow object (this is a sciline.Pipeline which can be consulted for advanced usage). The Sans2d workflow uses Mantid to load files. This tutorial comes with files that do not require Mantid, so we use a slightly modified workflow that does not require Mantid. The workflow is otherwise identical to the full Mantid-based workflow:

[2]:

workflow = isis.sans2d.Sans2dTutorialWorkflow()
# For real data use:
# workflow = isis.sans2d.Sans2dWorkflow()

We can insert steps for configuring the workflow. 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. We also want to compute the beam center using a simple center-of-mass estimation:

[3]:

workflow.insert(isis.io.transmission_from_background_run)
workflow.insert(isis.io.transmission_from_sample_run)

The workflow can be visualized as a graph. For readability we show only sub-workflow for computing IofQ[Sample]. The workflow can actually compute the full BackgroundSubtractedIofQ, which applies and equivalent workflow to the background run, before a subtraction step:

[4]:

# left-right layout works better for this graph
workflow.visualize(IofQ[SampleRun], graph_attr={"rankdir": "LR"})

[4]:

../../_images/user-guide_isis_sans2d_7_0.svg

Note the red boxes which indicate missing input parameters. We can set these missing parameters, as well as parameters where we do not want to use the defaults:

[5]:

workflow[OutFilename] = 'reduced.nxs'

workflow[NeXusMonitorName[Incident]] = 'monitor2'
workflow[NeXusMonitorName[Transmission]] = 'monitor4'

workflow[isis.SampleOffset] = sc.vector([0.0, 0.0, 0.053], unit='m')
workflow[isis.MonitorOffset[Transmission]] = sc.vector([0.0, 0.0, -6.719], unit='m')

workflow[WavelengthBins] = sc.linspace(
    'wavelength', start=2.0, stop=16.0, num=141, unit='angstrom'
)

workflow[isis.sans2d.LowCountThreshold] = sc.scalar(100, unit='counts')

mask_interval = sc.array(dims=['wavelength'], values=[2.21, 2.59], unit='angstrom')
workflow[WavelengthMask] = sc.DataArray(
    sc.array(dims=['wavelength'], values=[True]),
    coords={'wavelength': mask_interval},
)

workflow[QBins] = sc.linspace(dim='Q', start=0.01, stop=0.6, num=141, unit='1/angstrom')
workflow[NonBackgroundWavelengthRange] = sc.array(
    dims=['wavelength'], values=[0.7, 17.1], unit='angstrom'
)
workflow[CorrectForGravity] = True
workflow[UncertaintyBroadcastMode] = UncertaintyBroadcastMode.upper_bound
workflow[ReturnEvents] = True

Configuring data to load#

We have not configured which files we want to load. In this tutorial, we use helpers to fetch the tutorial data which return the filenames of the cached files. In a real use case, you would set these parameters manually:

[6]:

workflow[DirectBeamFilename] = isis.data.sans2d_tutorial_direct_beam()
workflow[Filename[SampleRun]] = isis.data.sans2d_tutorial_sample_run()
workflow[Filename[BackgroundRun]] = isis.data.sans2d_tutorial_background_run()
workflow[Filename[EmptyBeamRun]] = isis.data.sans2d_tutorial_empty_beam_run()

Downloading file 'DIRECT_SANS2D_REAR_34327_4m_8mm_16Feb16.dat.h5' from 'https://public.esss.dk/groups/scipp/ess/sans2d/1/DIRECT_SANS2D_REAR_34327_4m_8mm_16Feb16.dat.h5' to '/home/runner/.cache/ess/sans2d/1'.
Downloading file 'SANS2D00063159.nxs.h5' from 'https://public.esss.dk/groups/scipp/ess/sans2d/1/SANS2D00063159.nxs.h5' to '/home/runner/.cache/ess/sans2d/1'.

The beam center is not set yet. Unless we have a previously computed value, we can do so now:

[7]:

center = sans.beam_center_from_center_of_mass(workflow)
center

[7]:

scipp.Variable (336 Bytes)

()

vector3

m

[ 0.09007484 -0.08263978 0. ]

Values:
array([ 0.09007484, -0.08263978,  0.        ])

Remember to update the workflow with the new center:

[8]:

workflow[BeamCenter] = center

Use the workflow#

Compute final result#

We can now compute the background-subtracted \(I(Q)\):

[9]:

result = workflow.compute(BackgroundSubtractedIofQ)
result.hist().plot(scale={'Q': 'log'}, norm='log')

[9]:

../../_images/user-guide_isis_sans2d_17_0.svg

As the result was computed in event-mode, we can also use a different \(Q\)-binning, without re-reducing the data:

[10]:

result.hist(Q=60).plot(scale={'Q': 'log'}, norm='log')

[10]:

../../_images/user-guide_isis_sans2d_19_0.svg

In the above we used an upper bound for the uncertainties of the normalization factors. We can also compute the result with dropped normalization-factor uncertainties. This is incorrect, but is useful for understanding whether the normalization factors significantly contribute to the uncertainty of the result:

[11]:

workflow[UncertaintyBroadcastMode] = UncertaintyBroadcastMode.drop
result_drop = workflow.compute(BackgroundSubtractedIofQ)
# Reset the UnsertaintyBroadcastMode to the old value
workflow[UncertaintyBroadcastMode] = UncertaintyBroadcastMode.upper_bound
sc.DataGroup(upper_bound=result, dropped=result_drop).hist().plot(norm='log')

[11]:

../../_images/user-guide_isis_sans2d_21_0.svg

Save reduced data to file#

ess.sans provides a function for saving the reduced data as an NXcanSAS file. It could be used directly with the result computed above, but we would have to provide the required metadata ourselves. Instead, we use Sciline to get all required information directly from the pipeline: (See also the File output docs.)

[12]:

from ess.sans.io import save_background_subtracted_iofq

workflow.bind_and_call(save_background_subtracted_iofq)

Compute intermediate results#

For inspection and debugging purposes we can also compute intermediate results. To avoid repeated computation (including costly loading of files) we can request multiple results at once, including the final result, if desired. For example:

[13]:

monitors = (
    WavelengthMonitor[SampleRun, Incident],
    WavelengthMonitor[SampleRun, Transmission],
    WavelengthMonitor[BackgroundRun, Incident],
    WavelengthMonitor[BackgroundRun, Transmission],
)
parts = (CleanSummedQ[SampleRun, Numerator], CleanSummedQ[SampleRun, Denominator])
iofqs = (IofQ[SampleRun], IofQ[BackgroundRun], BackgroundSubtractedIofQ)
keys = (*monitors, MaskedData[SampleRun], *parts, *iofqs)

results = workflow.compute(keys)

display(sc.plot({str(key): results[key] for key in monitors}, norm='log'))

display(
    isis.plot_flat_detector_xy(
        results[MaskedData[SampleRun]]['spectrum', :61440].hist(), norm='log'
    )
)

wavelength = workflow.compute(WavelengthBins)
display(
    results[CleanSummedQ[SampleRun, Numerator]]
    .hist(wavelength=wavelength)
    .transpose()
    .plot(norm='log')
)
display(results[CleanSummedQ[SampleRun, Denominator]].plot(norm='log'))
parts = {str(key): results[key].sum('wavelength') for key in parts}
display(sc.plot(parts, norm='log'))

iofqs = {str(key): results[key] for key in iofqs}
iofqs = {key: val if val.bins is None else val.hist() for key, val in iofqs.items()}
display(sc.plot(iofqs, norm='log'))