Neutron Powder Diffraction#

In this tutorial demonstrates how neutron-scattering data can be loaded, visualized, and manipulated with generic functionality from scipp as well as neutron-specific functionality from scippneutron. It focuses on reducing data from the ORNL POWGEN neutron diffractometer.

[1]:
import numpy as np
import scipp as sc
import scippneutron as scn

Loading Nexus files#

Loading Nexus files requires Mantid. See, e.g., Installation on how to install scipp and Mantid with conda.

We start by loading two files: the sample and the vanadium runs.

[2]:
raw_sample = scn.load_with_mantid(
    scn.data.get_path('PG3_4844_event.nxs'),
    load_pulse_times=False,
    mantid_args={'LoadMonitors': True},
)
raw_vanadium = scn.load_with_mantid(scn.data.get_path('PG3_4866_event.nxs'), load_pulse_times=False)
Downloading file 'PG3_4844_event.nxs' from 'https://public.esss.dk/groups/scipp/scippneutron/5/PG3_4844_event.nxs' to '/home/runner/.cache/scippneutron/5'.
FrameworkManager-[Notice] Welcome to Mantid 6.8.0
FrameworkManager-[Notice] Please cite: http://dx.doi.org/10.1016/j.nima.2014.07.029 and this release: http://dx.doi.org/10.5286/Software/Mantid6.8
CheckMantidVersion-[Notice] A new version of Mantid(6.9.0) is available for download from https://download.mantidproject.org
DownloadInstrument-[Notice] Downloading 5 files from the instrument repository
Load-[Notice] Load started
InstrumentDefinitionParser-[Notice] Geometry cache is not available
InstrumentDefinitionParser-[Notice] Creating cache in /home/runner/.mantid/instrument/geometryCache/POWGEN243edfd42ebf5d0140e5f79fd93cf4f644e0815f.vtp
Load-[Notice] Load successful, Duration 1.70 seconds
Downloading file 'PG3_4866_event.nxs' from 'https://public.esss.dk/groups/scipp/scippneutron/5/PG3_4866_event.nxs' to '/home/runner/.cache/scippneutron/5'.
ExtractSpectra-[Notice] ExtractSpectra started
ExtractSpectra-[Notice] ExtractSpectra successful, Duration 0.00 seconds
DeleteWorkspace-[Notice] DeleteWorkspace started
DeleteWorkspace-[Notice] DeleteWorkspace successful, Duration 0.00 seconds
DeleteWorkspace-[Notice] DeleteWorkspace started
DeleteWorkspace-[Notice] DeleteWorkspace successful, Duration 0.00 seconds
DeleteWorkspace-[Notice] DeleteWorkspace started
DeleteWorkspace-[Notice] DeleteWorkspace successful, Duration 0.00 seconds
Load-[Notice] Load started
Load-[Notice] Load successful, Duration 4.06 seconds
DeleteWorkspace-[Notice] DeleteWorkspace started
DeleteWorkspace-[Notice] DeleteWorkspace successful, Duration 0.00 seconds

The optional mantid_args dict is forwarded to the Mantid algorithm used for loading the files – in this case LoadEventNexus – and can be used to control, e.g., which part of a file to load. Here we request loading monitors, which Mantid does not load by default. The resulting dataset looks as follows:

[3]:
raw_sample
[3]:
  • data
    scipp
    DataArray
    (spectrum: 24794, tof: 1)
    DataArrayView
    binned data [len=0, len=0, ..., len=0, len=0]
  • sample
    scipp
    Variable
    ()
    PyObject
    <mantid.api._api.Sample object at 0x7ff05434e960>
  • instrument_name
    str
    ()
    POWGEN
  • start_time
    scipp
    Variable
    ()
    string
    𝟙
    2011-08-12T15:50:17
  • end_time
    scipp
    Variable
    ()
    string
    𝟙
    2011-08-12T17:22:05
  • ChopperStatus1
    scipp
    DataArray
    (time: 2)
    float64
    𝟙
    4.0, 4.0
  • ChopperStatus2
    scipp
    DataArray
    (time: 2)
    float64
    𝟙
    4.0, 4.0
  • ChopperStatus3
    scipp
    DataArray
    (time: 2)
    float64
    𝟙
    4.0, 4.0
  • CurrentSP
    scipp
    DataArray
    (time: 2)
    float64
    𝟙
    300.0, 300.0
  • EnergyRequest
    scipp
    DataArray
    (time: 2)
    float64
    meV
    287.955, 287.955
  • LKSRampRate
    scipp
    DataArray
    (time: 2)
    float64
    𝟙
    0.0, 0.0
  • LambdaRequest
    scipp
    DataArray
    (time: 2)
    float64
    Å
    0.533, 0.533
  • Phase1
    scipp
    DataArray
    (time: 1794)
    float64
    µs
    8166.720, 8165.158, ..., 8163.853, 8163.853
  • Phase2
    scipp
    DataArray
    (time: 1793)
    float64
    µs
    8335.626, 8334.088, ..., 8332.859, 8332.859
  • Phase3
    scipp
    DataArray
    (time: 1777)
    float64
    µs
    1.400e+04, 1.400e+04, ..., 1.400e+04, 1.400e+04
  • PhaseRequest1
    scipp
    DataArray
    (time: 2)
    float64
    µs
    8164.075, 8164.075
  • PhaseRequest2
    scipp
    DataArray
    (time: 2)
    float64
    µs
    8332.893, 8332.893
  • PhaseRequest3
    scipp
    DataArray
    (time: 2)
    float64
    µs
    1.400e+04, 1.400e+04
  • SampleTemp
    scipp
    DataArray
    (time: 467)
    float64
    𝟙
    299.352, 299.446, ..., 300.0, 300.0
  • Speed1
    scipp
    DataArray
    (time: 2)
    float64
    Hz
    60.0, 60.0
  • Speed2
    scipp
    DataArray
    (time: 2)
    float64
    Hz
    60.0, 60.0
  • Speed3
    scipp
    DataArray
    (time: 2)
    float64
    Hz
    60.0, 60.0
  • SpeedRequest1
    scipp
    DataArray
    (time: 2)
    float64
    Hz
    60.0, 60.0
  • SpeedRequest2
    scipp
    DataArray
    (time: 2)
    float64
    Hz
    60.0, 60.0
  • SpeedRequest3
    scipp
    DataArray
    (time: 2)
    float64
    Hz
    60.0, 60.0
  • TolRequest
    scipp
    DataArray
    (time: 2)
    float64
    𝟙
    20.0, 20.0
  • currentsample
    scipp
    DataArray
    (time: 2)
    float64
    𝟙
    4.0, 4.0
  • fernsstatus
    scipp
    DataArray
    (time: 2)
    float64
    𝟙
    3.0, 3.0
  • frequency
    scipp
    DataArray
    (time: 330473)
    float64
    Hz
    0.0, 60.024, ..., 60.002, 59.999
  • proton_charge
    scipp
    DataArray
    (time: 330473)
    float64
    pC
    1.484e+07, 1.484e+07, ..., 1.487e+07, 1.484e+07
  • samplerequest
    scipp
    DataArray
    (time: 2)
    float64
    𝟙
    4.0, 4.0
  • S1HCenter
    scipp
    DataArray
    (time: 2)
    float64
    mm
    0.0, 0.0
  • S1HCenterOffset
    scipp
    DataArray
    (time: 2)
    float64
    𝟙
    0.0, 0.0
  • S1HWidth
    scipp
    DataArray
    (time: 2)
    float64
    mm
    10.0, 10.0
  • S1VCenter
    scipp
    DataArray
    (time: 2)
    float64
    mm
    5.0, 5.0
  • S1VCenterOffset
    scipp
    DataArray
    (time: 2)
    float64
    𝟙
    0.0, 0.0
  • S1VHeight
    scipp
    DataArray
    (time: 2)
    float64
    mm
    30.0, 30.0
  • commErrs
    scipp
    DataArray
    (time: 2)
    float64
    𝟙
    0.0, 0.0
  • guide
    scipp
    DataArray
    (time: 2)
    float64
    mm
    -55.463, -55.463
  • s1b
    scipp
    DataArray
    (time: 2)
    float64
    mm
    20.0, 20.0
  • s1l
    scipp
    DataArray
    (time: 2)
    float64
    mm
    5.0, 5.0
  • s1r
    scipp
    DataArray
    (time: 2)
    float64
    mm
    -5.0, -5.0
  • s1t
    scipp
    DataArray
    (time: 2)
    float64
    mm
    -10.0, -10.0
  • vGuide
    scipp
    DataArray
    (time: 2)
    float64
    𝟙
    2.0, 2.0
  • veto_pulse_time
    scipp
    DataArray
    (time: 1)
    float64
    𝟙
    0.0
  • gd_prtn_chrg
    scipp
    Variable
    ()
    float64
    µAh
    1171.953902925
  • run_start
    scipp
    Variable
    ()
    string
    𝟙
    2011-08-12T15:50:17
  • run_title
    scipp
    Variable
    ()
    string
    𝟙
    diamond cw0.533 4.22e12 60Hz [10x30]
  • file_notes
    scipp
    Variable
    ()
    string
    𝟙
    NONE
  • run_number
    scipp
    Variable
    ()
    string
    𝟙
    4844
  • experiment_identifier
    scipp
    Variable
    ()
    string
    𝟙
    IPTS-2767
  • duration
    scipp
    Variable
    ()
    float64
    s
    5508.0
  • running
    scipp
    DataArray
    (time: 1)
    bool
    𝟙
    True
  • Filename
    scipp
    Variable
    ()
    string
    𝟙
    /home/runner/.cache/scippneutron/5/PG3_4844_event.nxs
  • scipp
    DataGroup
    (tof: 200001)
      • scipp
        DataGroup
        (tof: 200001)
          • data
            scipp
            DataArray
            (tof: 200001)
            float64
            counts
            25.0, 10.0, ..., 0.0, 0.0
          • instrument_name
            str
            ()
            POWGEN

Extract the actual events.

[4]:
sample = raw_sample['data']
vanadium = raw_vanadium['data']
[5]:
sc.plot(sample.hist(spectrum=500, tof=400))
[5]:
../_images/tutorials_powder-diffraction_8_0.svg

Instrument view#

Scipp provides a simple 3D instrument view inspired by Mantid’s own instrument view, which can be used to take a quick look at the neutron counts on the detector panels in 3D space or using various cylindrical and spherical projections

Note

The support and handling of event data is currently being rewritten. The interface for manipulating event data (binned) data is not representative for the final API.

[6]:
scn.instrument_view(sample.hist())
[6]:

Plot against scattering angle \(\theta\) using bin#

This is not an essential step and can be skipped.

Plotting raw data directly yields a hard-to-interpret figure. We can obtain something more useful by binning the spectrum axis based on its \(\theta\) value:

[7]:
sample.coords['two_theta'] = scn.two_theta(sample)
vanadium.coords['two_theta'] = scn.two_theta(vanadium)
two_theta = sc.linspace(dim='two_theta', unit='rad', start=0.0, stop=np.pi, num=2001)

We concatenate events lists from different spectra that fall into a given \(2\theta\) range into longer combined lists:

[8]:
theta_sample = sample.bin(two_theta=two_theta)
[9]:
theta_sample
[9]:
Show/Hide data repr Show/Hide attributes
scipp.DataArray (273.59 MB)
    • tof: 1
    • two_theta: 2000
    • sample_position
      ()
      vector3
      m
      [0. 0. 0.]
      Values:
      array([0., 0., 0.])
    • source_position
      ()
      vector3
      m
      [ 0. 0. -60.]
      Values:
      array([ 0., 0., -60.])
    • tof
      (tof [bin-edge])
      float64
      µs
      19.0, 1.669e+04
      Values:
      array([ 19. , 16694.30078125])
    • two_theta
      (two_theta [bin-edge])
      float64
      rad
      0.0, 0.002, ..., 3.140, 3.142
      Values:
      array([0.00000000e+00, 1.57079633e-03, 3.14159265e-03, ..., 3.13845106e+00, 3.14002186e+00, 3.14159265e+00])
    • (tof, two_theta)
      DataArrayView
      binned data [len=0, len=0, ..., len=0, len=0]
      dim='event',
      content=DataArray(
                dims=(event: 17926980),
                data=float32[counts],
                coords={'tof':float64[µs]})
[10]:
sc.plot(theta_sample.hist(tof=400))
[10]:
../_images/tutorials_powder-diffraction_17_0.svg

Coordinate transformation#

Note: We are back to working with ``sample``, not ``theta_sample``.

scippneutron provides building blocks for scipp.transform_coords to convert between coordinates related to time-of-flight. The loaded raw data has dimension tof, and we convert to interplanar lattice spacing (dspacing):

[11]:
dspacing_graph = {
    **scn.conversion.graph.beamline.beamline(scatter=True),
    **scn.conversion.graph.tof.elastic_dspacing('tof'),
}
[12]:
dspacing_vanadium = vanadium.transform_coords('dspacing', graph=dspacing_graph)
dspacing_sample = sample.transform_coords('dspacing', graph=dspacing_graph)
dspacing_sample
[12]:
Show/Hide data repr Show/Hide attributes
scipp.DataArray (412.88 MB)
    • spectrum: 24794
    • dspacing: 1
    • L1
      ()
      float64
      m
      60.0
      Values:
      array(60.)
    • L2
      (spectrum)
      float64
      m
      2.560, 2.556, ..., 4.251, 4.252
      Values:
      array([2.56027902, 2.55590928, 2.55268677, ..., 4.25113991, 4.25116338, 4.25188 ])
    • Ltotal
      (spectrum)
      float64
      m
      62.560, 62.556, ..., 64.251, 64.252
      Values:
      array([62.56027902, 62.55590928, 62.55268677, ..., 64.25113991, 64.25116338, 64.25188 ])
    • dspacing
      (spectrum, dspacing [bin-edge])
      float64
      Å
      0.001, 0.557, ..., 0.003, 2.343
      Values:
      array([[6.33939029e-04, 5.57008885e-01], [6.32685444e-04, 5.55907426e-01], [6.31487177e-04, 5.54854572e-01], ..., [2.67322756e-03, 2.34882447e+00], [2.67090407e-03, 2.34678294e+00], [2.66651403e-03, 2.34292565e+00]])
    • incident_beam
      ()
      vector3
      m
      [ 0. 0. 60.]
      Values:
      array([ 0., 0., 60.])
    • position
      (spectrum)
      vector3
      m
      [ 1.17451004 -1.01106149 -2.03796699], [ 1.18147634 -0.95946649 -2.05334117], ..., [1.81428985 0.09565841 3.84338287], [1.81375055 0.1499371 3.84269584]
      Values:
      array([[ 1.17451004, -1.01106149, -2.03796699], [ 1.18147634, -0.95946649, -2.05334117], [ 1.18844265, -0.90787149, -2.06871534], ..., [ 1.81482915, 0.04137972, 3.8440699 ], [ 1.81428985, 0.09565841, 3.84338287], [ 1.81375055, 0.1499371 , 3.84269584]])
    • sample_position
      ()
      vector3
      m
      [0. 0. 0.]
      Values:
      array([0., 0., 0.])
    • scattered_beam
      (spectrum)
      vector3
      m
      [ 1.17451004 -1.01106149 -2.03796699], [ 1.18147634 -0.95946649 -2.05334117], ..., [1.81428985 0.09565841 3.84338287], [1.81375055 0.1499371 3.84269584]
      Values:
      array([[ 1.17451004, -1.01106149, -2.03796699], [ 1.18147634, -0.95946649, -2.05334117], [ 1.18844265, -0.90787149, -2.06871534], ..., [ 1.81482915, 0.04137972, 3.8440699 ], [ 1.81428985, 0.09565841, 3.84338287], [ 1.81375055, 0.1499371 , 3.84269584]])
    • source_position
      ()
      vector3
      m
      [ 0. 0. -60.]
      Values:
      array([ 0., 0., -60.])
    • spectrum
      (spectrum)
      int32
      1, 2, ..., 24793, 24794
      Values:
      array([ 1, 2, 3, ..., 24792, 24793, 24794], dtype=int32)
    • tof
      (dspacing [bin-edge])
      float64
      µs
      19.0, 1.669e+04
      Values:
      array([ 19. , 16694.30078125])
    • two_theta
      (spectrum)
      float64
      rad
      2.491, 2.504, ..., 0.442, 0.442
      Values:
      array([2.49144445, 2.50372965, 2.51564288, ..., 0.44118918, 0.44157918, 0.44231324])
    • (spectrum, dspacing)
      DataArrayView
      binned data [len=0, len=0, ..., len=0, len=0]
      dim='event',
      content=DataArray(
                dims=(event: 17926980),
                data=float32[counts],
                coords={'tof':float64[µs], 'dspacing':float64[Å]})

Neutron monitors#

This is an optional section. The next section does not use the monitor-normalized data produced here. This section could thus be skipped.

If available, neutron monitors are stored as attributes of a data array:

[13]:
mon = raw_sample['monitors']['monitor1']['data']
mon
[13]:
Show/Hide data repr Show/Hide attributes
scipp.DataArray (4.58 MB)
    • tof: 200001
    • position
      ()
      vector3
      m
      [ 0. 0. -1.]
      Values:
      array([ 0., 0., -1.])
    • sample_position
      ()
      vector3
      m
      [0. 0. 0.]
      Values:
      array([0., 0., 0.])
    • source_position
      ()
      vector3
      m
      [ 0. 0. -60.]
      Values:
      array([ 0., 0., -60.])
    • spectrum
      ()
      int32
      1
      Values:
      array(1, dtype=int32)
    • tof
      (tof [bin-edge])
      float64
      µs
      0.0, 1.0, ..., 2.000e+05, 2.000e+05
      Values:
      array([0.00000e+00, 1.00000e+00, 2.00000e+00, ..., 1.99999e+05, 2.00000e+05, 2.00001e+05])
    • (tof)
      float64
      counts
      25.0, 10.0, ..., 0.0, 0.0
      σ = 5.0, 3.162, ..., 0.0, 0.0
      Values:
      array([25., 10., 12., ..., 0., 0., 0.])

      Variances (σ²):
      array([25., 10., 12., ..., 0., 0., 0.])

The monitor could, e.g., be used to normalize the data. To do so, both data and monitor need to be converted to a unit that accounts for differing flight paths, e.g., wavelength or energy:

[14]:
wavelength_graph = {
    **scn.conversion.graph.beamline.beamline(scatter=True),
    **scn.conversion.graph.tof.elastic_wavelength('tof'),
}
wavelength_graph_monitor = {
    **scn.conversion.graph.beamline.beamline(scatter=False),
    **scn.conversion.graph.tof.elastic_wavelength('tof'),
}
[15]:
sample_lambda = sample.transform_coords('wavelength', graph=wavelength_graph)
mon = mon.transform_coords('wavelength', graph=wavelength_graph_monitor)

The sample data is in event-mode, i.e., is not histogrammed. Event data can be divided by a histogram (such as mon in this case), using a specialized function for scaling (see Binned data). First we rebin the monitor since the original binning is very fine:

[16]:
edges = sc.linspace(dim='wavelength', unit='Angstrom', start=0.0, stop=1.0, num=1001)
mon = sc.rebin(mon, wavelength=edges)
mon
[16]:
Show/Hide data repr Show/Hide attributes
scipp.DataArray (25.97 KB)
    • wavelength: 1000
    • Ltotal
      ()
      float64
      m
      59.0
      Values:
      array(59.)
    • position
      ()
      vector3
      m
      [ 0. 0. -1.]
      Values:
      array([ 0., 0., -1.])
    • sample_position
      ()
      vector3
      m
      [0. 0. 0.]
      Values:
      array([0., 0., 0.])
    • source_position
      ()
      vector3
      m
      [ 0. 0. -60.]
      Values:
      array([ 0., 0., -60.])
    • spectrum
      ()
      int32
      1
      Values:
      array(1, dtype=int32)
    • wavelength
      (wavelength [bin-edge])
      float64
      Å
      0.0, 0.001, ..., 0.999, 1.0
      Values:
      array([0. , 0.001, 0.002, ..., 0.998, 0.999, 1. ])
    • (wavelength)
      float64
      counts
      194.397, 214.709, ..., 206.055, 193.213
      σ = 13.943, 14.653, ..., 14.355, 13.900
      Values:
      array([194.39748464, 214.70889565, 215.31141101, 221.66343536, 195.61509218, 184.13926377, 210.34429449, 204.1798771 , 202.92165623, 221.50386493, 224.95938667, 219.86073623, 204.61024521, 226.93908 , 215.93423304, 224.13926377, 225.8965026 , 214.2833745 , 191.96711652, 228.3597542 , 261.99711711, 196.54447826, 218.12956985, 205.67312928, 227.49128811, 199.62282203, 200.59478551, 198.8810429 , 193.91196231, 217.1034974 , 214.63251595, 218.10441622, 201.85497045, 194.31141101, 194.24852695, 238.97969333, 229.6740481 , 198.44779189, 212.49325218, 195.1527594 , 225.12576813, 204.21852636, 205.06000118, 191.28822146, 209.86361912, 180.46901738, 199.70889565, 214.97196348, 190.09865044, 223.87816001, 201.05607304, 210.87043015, 233.6923907 , 207.69055306, 191.67889506, 231.44294493, 239.33368175, 222.01454088, 192.32687072, 212.6402458 , 195.73012114, 225.54447826, 197.52705448, 196.09773162, 217.38097968, 217.17411132, 226.09576755, 214.34717738, 232.00288289, 203.68374202, 226.62282203, 193.41202552, 208.68374202, 226.40325042, 213.7254006 , 214.05607304, 212.33944753, 198.02607245, 217.73601335, 233.8810429 , 198.79693334, 201.08122666, 180.75619359, 196.50110846, 208.58601041, 213.16245332, 214.35294317, 193.81147423, 224.55993796, 202.69540002, 198.89349328, 222.45656699, 207.79012231, 197.07257799, 221.69920173, 187.82877157, 200.5367484 , 197.27760871, 205.9661977 , 208.93698951, 216.43337744, 220.74662609, 223.72527418, 199.4294493 , 223.81619476, 214.3752139 , 200.87043015, 228.7349681 , 251.89270089, 229.2407971 , 219.9478551 , 218.57828056, 203.70024698, 197.32687072, 211.2562568 , 193.94288172, 225.39552057, 236.89453853, 209.77675311, 221.08018141, 204.5367484 , 219.03576637, 200.38595307, 210.24668931, 210.0338023 , 235.7437432 , 187.63723648, 246.6740481 , 219.43325101, 215.77767193, 219.50190086, 218.51932463, 212.79496927, 212.32883479, 200.12092117, 209.2291391 , 221.0851548 , 223.21079651, 215.00956749, 228.4294493 , 227.5367484 , 227.7592029 , 212.62282203, 214.41202552, 214.31429391, 217.14214666, 215.23411249, 212.79484285, 212.01454088, 221.18092235, 232.90122314, 181.77767193, 217.4720267 , 205.37324983, 209.17398489, 223.88576343, 224.66473346, 205.60042486, 216.65177737, 198.45067478, 215.73994149, 205.049262 , 224.77282497, 214.28533857, 208.01245039, 209.28533857, 191.8810429 , 199.19809327, 220.83178089, 231.61325453, 187.59963247, 191.35779013, 232.12275881, 213.3280424 , 209.99699068, 189.32791597, 224.0838567 , 215.55246097, 201.96711652, 205.95362089, 180.87724119, 183.07048749, 226.97300873, 203.44858428, 208.40337685, 201.02410838, 210.60343418, 205.16834553, 199.07729852, 213.62282203, 198.18275999, 187.79810502, 194.88760107, 181.41595366, 195.53582958, 203.08281145, 188.23921231, 200.05319014, 213.31049219, 220.01061274, 227.93607068, 201.79496927, 257.74086031, 208.17202082, 186.76221222, 204.96711652, 187.88091647, 219.24668931, 230.79496927, 195.56858663, 223.18288642, 194.63343477, 186.46233277, 200.44858428, 206.97981976, 203.43245861, 187.00760342, 201.61876746, 202.31546558, 190.83951075, 204.93882715, 202.01466731, 196.90895299, 184.13362441, 189.50085561, 214.62282203, 193.70024698, 206.17018318, 199.28717621, 212.92754844, 188.14594837, 208.44399018, 211.46404399, 237.98466672, 163.84540296, 171.76601394, 207.71361618, 228.56846021, 208.50896474, 191.87790715, 197.91313398, 193.98821558, 216.36656523, 194.54853283, 202.96515245, 216.92244862, 208.5321543 , 176.71545383, 191.40141278, 201.70889565, 208.31337508, 206.41386316, 204.39944871, 194.25809444, 225.00601864, 211.00196407, 188.53687483, 195.36251066, 202.03643233, 185.89780071, 214.57146953, 169.81147423, 204.96895416, 184.79313163, 184.8810429 , 200.29870777, 170.37166505, 204.76956279, 196.23804063, 185.56833378, 216.38804357, 200.35699773, 170.555091 , 185.96711652, 199.31834847, 187.61588457, 175.9127547 , 193.16454382, 217.85313281, 190.23373321, 205.41817059, 187.94760225, 194.01035989, 187.53136191, 195.1014069 , 204.8810429 , 164.28730264, 188.10717268, 199.24865338, 199.0742892 , 197.87749404, 197.42434948, 181.70889565, 175.32018612, 176.09760519, 171.23110317, 196.31442033, 189.93122372, 192.30288877, 184.76484226, 202.15785922, 184.24263474, 189.20804004, 187.05895593, 149.40153921, 203.50728735, 199.11214607, 187.96135074, 185.5530005 , 191.88614272, 202.09117344, 187.89152921, 160.84920467, 172.00367528, 192.38804357, 164.02607245, 188.7099409 , 181.71152569, 160.08135309, 198.33643821, 150.33852871, 191.51067596, 184.28638382, 174.36355591, 177.96711652, 181.08711887, 179.26281497, 177.78356414, 181.7349681 , 190.18000353, 171.67718385, 163.46771927, 168.55600982, 174.7032563 , 209.5571815 , 197.54919879, 175.03956808, 162.50950428, 182.78945635, 173.93435947, 157.69527359, 169.03944165, 177.41123313, 160.22533739, 179.92073741, 180.18551645, 178.83466378, 184.86767369, 172.15931758, 154.75527476, 178.21865279, 169.45669342, 155.24551763, 169.13926377, 177.47085502, 182.8810429 , 172.60932639, 163.41018788, 157.92232219, 178.56883949, 173.24133663, 189.16885125, 170.49116168, 158.3055188 , 169.04558672, 167.56202846, 176.54435183, 178.37809679, 163.20594955, 174.3393211 , 152.35607891, 143.06405574, 170.1865617 , 168.0338023 , 164.30643762, 186.97208991, 177.23975185, 154.57460527, 168.46914381, 168.85760049, 167.79091471, 155.90356649, 183.00405457, 160.22220164, 170.05319014, 159.72644585, 159.00091882, 175.43220576, 162.08791127, 170.8810429 , 148.8810429 , 167.36276352, 169.18785981, 156.39656582, 163.57238835, 162.885637 , 174.8810429 , 159.62374085, 181.07809091, 194.37074623, 172.4515936 , 174.8810429 , 160.63172356, 167.90410603, 161.10729911, 174.61049807, 179.2474817 , 184.95888096, 167.29908706, 176.28926671, 146.26737526, 175.95704331, 157.61233571, 141.59085737, 159.97760284, 153.75394284, 157.79354474, 163.90461174, 148.33146482, 167.92349387, 160.29018553, 160.20084973, 181.2136794 , 157.50804592, 136.37834965, 187.73012114, 145.81828526, 144.24894006, 167.16467025, 166.0020905 , 154.93214254, 161.70233748, 162.15238012, 182.53829937, 147.12760578, 156.25876041, 157.79496927, 174.19546323, 144.53975772, 178.75644644, 160.52208109, 145.96109788, 149.05619946, 181.88367293, 152.1698965 , 152.40312399, 162.36159184, 157.13926377, 159.49220693, 135.58692923, 152.7850225 , 163.59817411, 163.2858781 , 158.47449648, 169.07913616, 142.56269442, 177.72448178, 145.22533739, 182.70195819, 181.34258328, 160.6314707 , 147.45932345, 143.547614 , 142.10194644, 149.69552645, 150.54289347, 152.61195643, 148.18301285, 139.92336744, 155.1439843 , 145.8810429 , 141.38804357, 157.5125136 , 150.7370586 , 142.31141101, 161.28245568, 145.82630179, 133.39589985, 154.36539355, 153.31768251, 161.47649437, 161.07821734, 146.45460292, 117.98821558, 150.92884654, 152.98428744, 164.63853459, 160.30355473, 134.22926553, 142.57748816, 156.08673959, 148.01178442, 153.05385611, 164.43755843, 138.97171063, 149.4382244 , 159.22533739, 124.81933051, 136.82839229, 154.91785452, 149.45865749, 126.89349328, 149.31939372, 150.43470936, 147.19650848, 152.6064435 , 151.97522566, 135.85250067, 131.02448766, 136.90974538, 157.950738 , 161.47136073, 160.63501956, 143.02055952, 167.39748464, 152.4995575 , 140.79522213, 142.96699009, 137.32816882, 133.63920055, 145.93398019, 129.34099849, 158.64274941, 149.12275881, 163.13926377, 134.86453794, 141.77478904, 153.11908353, 140.19980448, 126.97614448, 149.90122314, 153.30393402, 130.72907589, 135.42146659, 155.4482976 , 116.65808269, 168.4962615 , 139.53294669, 156.52404516, 128.65950722, 136.31651083, 164.44177325, 149.40128635, 154.52534327, 125.95190967, 153.70892948, 166.69749052, 145.22400546, 141.95704331, 148.32908764, 152.31141101, 159.33669107, 145.48355826, 159.11398371, 142.85313281, 144.34049278, 147.02556674, 147.9682882 , 154.25441916, 141.22533739, 135.22533739, 143.85196113, 142.9770633 , 130.30146424, 144.87606951, 151.24944577, 144.11474228, 131.75752551, 151.67442739, 141.28730264, 131.86729441, 130.4578651 , 152.92925965, 162.36618594, 153.45944988, 130.63159713, 152.80374438, 146.22312046, 141.06418217, 128.38870954, 151.22533739, 138.22533739, 126.76810443, 143.48355826, 122.48355826, 139.5312693 , 136.28454617, 139.52195467, 157.80700655, 170.78895064, 143.40299756, 140.18142806, 142.85041017, 149.98346122, 119.41452913, 131.95282849, 151.53950487, 136.13926377, 149.11636089, 137.45051453, 128.77190615, 153.85091588, 139.09225251, 149.55363264, 121.28834789, 120.79325806, 135.08949605, 150.72213843, 140.31609772, 136.14097498, 124.21551704, 145.54827997, 144.61300168, 142.36393519, 122.24668931, 135.31207698, 154.73091353, 148.56896592, 128.30054541, 135.07492134, 156.13888448, 129.26946575, 133.91250184, 133.84920467, 139.02135191, 148.8397636 , 141.87857312, 137.05319014, 159.63226309, 115.51156095, 149.0783776 , 138.18614859, 137.24982505, 120.80231984, 148.5388389 , 154.68336274, 140.94050454, 134.90870013, 126.05632589, 162.96293552, 149.833079 , 141.10115404, 123.34847549, 143.29151746, 145.33130457, 131.10428979, 141.38072683, 131.77634 , 127.98428744, 151.4157008 , 143.20712122, 152.43141337, 151.98533269, 147.22533739, 134.09066773, 140.8243039 , 147.66029961, 142.50741378, 150.93407279, 117.29230985, 144.84695392, 133.60815472, 147.92979919, 129.55939842, 137.46429684, 137.03617948, 164.29778895, 135.10872364, 145.63644409, 126.13128106, 141.39184528, 133.05319014, 142.12564171, 134.9398724 , 126.59253476, 135.55471172, 160.86871894, 115.66932757, 129.72945517, 143.29493988, 152.55192143, 146.50256682, 137.69682455, 128.86048338, 126.22533739, 134.31834847, 137.73458881, 133.84172767, 147.56887331, 143.88154861, 158.09944283, 134.00693746, 151.75514834, 153.35123195, 136.5166946 , 127.58968569, 137.02645173, 153.76237248, 124.53620887, 138.13926377, 152.47687366, 155.6490209 , 160.27238247, 162.58379348, 156.35043956, 158.19913852, 164.29724942, 140.1055879 , 152.24029138, 146.1055879 , 136.64822851, 143.74257152, 160.13926377, 179.3129958 , 150.22137543, 150.73984888, 158.42719854, 140.59662315, 132.42447591, 155.44126754, 169.65491311, 146.4185837 , 153.15535561, 137.04308311, 163.62193703, 147.71057304, 169.05908235, 144.60004558, 144.81017612, 150.62871424, 191.72068007, 151.05319014, 138.05319014, 149.05319014, 128.9019817 , 149.6051454 , 150.34366235, 155.83963717, 168.00915439, 160.78644703, 157.98352887, 175.42247801, 138.3501867 , 153.9794743 , 175.48824497, 161.3501867 , 142.71032018, 130.97421423, 147.39890917, 174.70747112, 154.6944812 , 150.32135779, 128.51685485, 163.57957866, 166.16809268, 158.68047985, 167.20014994, 176.77691336, 150.20056305, 149.39343007, 172.53899916, 151.50788567, 140.78686014, 152.18025638, 164.02252359, 135.21501133, 139.0838567 , 121.50829878, 155.0162521 , 157.9486475 , 132.8647908 , 154.41373673, 152.13926377, 175.38123254, 142.31141101, 143.93201612, 146.57001117, 137.95302257, 163.58518419, 139.51682102, 166.52341302, 159.17068889, 163.24013113, 162.34796978, 156.57146953, 175.39932228, 170.66682396, 145.74545441, 169.01653878, 164.43781129, 153.30589809, 145.31141101, 136.02949488, 167.22349975, 176.25270794, 140.166921 , 154.8810429 , 156.94129693, 154.3664388 , 161.77098733, 156.11528182, 162.33539296, 154.34245685, 171.12418335, 154.34157185, 164.4657552 , 157.8527197 , 188.77282497, 152.41072742, 172.49245979, 143.97396138, 190.8013672 , 142.87853929, 156.43926965, 158.58354063, 176.14749933, 147.43103408, 175.09814472, 167.21459822, 171.26579047, 161.31141101, 158.30133781, 152.4828923 , 154.07466849, 153.31157127, 176.09614683, 167.50453089, 171.33288936, 166.67718385, 154.69815648, 174.9444665 , 164.70889565, 171.10755197, 167.8810429 , 155.72432153, 154.8255094 , 159.13692041, 166.35620534, 166.00956749, 171.53557673, 192.35503366, 172.10520861, 156.80149363, 175.31445416, 152.52718091, 168.15418394, 154.22533739, 142.26360737, 156.18706741, 164.06275764, 154.89599689, 172.37834965, 170.98504601, 157.77703979, 195.69514716, 161.74596013, 153.7925921 , 157.99280968, 166.24146306, 187.81828526, 148.09982211, 190.09855783, 185.85243302, 161.20857958, 172.05319014, 185.51343242, 173.84461056, 175.51343242, 188.73522096, 168.00709771, 190.58585016, 200.11648732, 173.54029726, 170.16558907, 178.50988356, 191.25467201, 174.87857312, 157.4536841 , 185.62583135, 180.79797859, 179.97012584, 191.97291613, 177.47474934, 184.60476611, 181.32344829, 185.22834671, 164.48431683, 156.62206346, 179.5726412 , 194.13926377, 181.52246038, 175.31236366, 208.85547616, 214.117659 , 175.85246684, 193.44335804, 184.45498221, 180.67300285, 201.8810429 , 170.94852107, 171.02290288, 170.75648027, 181.52074916, 167.58822734, 204.82915086, 168.50126872, 178.08737173, 176.44896357, 172.31141101, 207.65741672, 210.55233454, 202.22533739, 185.64011937, 192.46797213, 204.78109436, 206.42121373, 198.23921231, 184.61537886, 186.31071123, 173.76379701, 171.33615153, 171.95796214, 175.49914439, 194.07520802, 181.58521802, 200.94066479, 182.58521802, 199.2520758 , 181.547614 , 186.18272617, 191.13812591, 166.72590632, 194.82642822, 215.49442386, 211.44452972, 174.18773338, 200.35988063, 174.21732086, 187.37261769, 180.70417512, 189.8668813 , 186.56979214, 217.5367484 , 222.28324807, 212.46939666, 183.11407632, 193.79199378, 187.01400135, 195.76069509, 175.17353796, 169.83102232, 205.21903207, 173.90762106, 199.99689808, 192.48038869, 198.64943401, 210.94050454, 172.62909352, 197.96084503, 238.74954281, 201.84039574, 185.99214372, 184.04219812, 182.37559319, 231.95308135, 186.34743024, 198.45853106, 193.11873807, 189.90942488, 191.50051015, 167.1394828 , 209.57109024, 182.1343242 , 205.44921642, 196.21602275, 198.74817705, 183.38782454, 208.99068536, 201.17674135, 214.70610536, 192.87318662, 173.70103937, 171.52889212, 191.96385435, 196.57748816, 214.66948782, 210.91126252, 210.85867956, 222.07023463, 192.89554996, 203.77792478, 198.06443503, 195.85855313, 189.03735116, 210.31020551, 209.50845903, 184.24117638, 196.03843023, 206.05549967, 193.21288701])

      Variances (σ²):
      array([194.39748464, 214.70889565, 215.31141101, 221.66343536, 195.61509218, 184.13926377, 210.34429449, 204.1798771 , 202.92165623, 221.50386493, 224.95938667, 219.86073623, 204.61024521, 226.93908 , 215.93423304, 224.13926377, 225.8965026 , 214.2833745 , 191.96711652, 228.3597542 , 261.99711711, 196.54447826, 218.12956985, 205.67312928, 227.49128811, 199.62282203, 200.59478551, 198.8810429 , 193.91196231, 217.1034974 , 214.63251595, 218.10441622, 201.85497045, 194.31141101, 194.24852695, 238.97969333, 229.6740481 , 198.44779189, 212.49325218, 195.1527594 , 225.12576813, 204.21852636, 205.06000118, 191.28822146, 209.86361912, 180.46901738, 199.70889565, 214.97196348, 190.09865044, 223.87816001, 201.05607304, 210.87043015, 233.6923907 , 207.69055306, 191.67889506, 231.44294493, 239.33368175, 222.01454088, 192.32687072, 212.6402458 , 195.73012114, 225.54447826, 197.52705448, 196.09773162, 217.38097968, 217.17411132, 226.09576755, 214.34717738, 232.00288289, 203.68374202, 226.62282203, 193.41202552, 208.68374202, 226.40325042, 213.7254006 , 214.05607304, 212.33944753, 198.02607245, 217.73601335, 233.8810429 , 198.79693334, 201.08122666, 180.75619359, 196.50110846, 208.58601041, 213.16245332, 214.35294317, 193.81147423, 224.55993796, 202.69540002, 198.89349328, 222.45656699, 207.79012231, 197.07257799, 221.69920173, 187.82877157, 200.5367484 , 197.27760871, 205.9661977 , 208.93698951, 216.43337744, 220.74662609, 223.72527418, 199.4294493 , 223.81619476, 214.3752139 , 200.87043015, 228.7349681 , 251.89270089, 229.2407971 , 219.9478551 , 218.57828056, 203.70024698, 197.32687072, 211.2562568 , 193.94288172, 225.39552057, 236.89453853, 209.77675311, 221.08018141, 204.5367484 , 219.03576637, 200.38595307, 210.24668931, 210.0338023 , 235.7437432 , 187.63723648, 246.6740481 , 219.43325101, 215.77767193, 219.50190086, 218.51932463, 212.79496927, 212.32883479, 200.12092117, 209.2291391 , 221.0851548 , 223.21079651, 215.00956749, 228.4294493 , 227.5367484 , 227.7592029 , 212.62282203, 214.41202552, 214.31429391, 217.14214666, 215.23411249, 212.79484285, 212.01454088, 221.18092235, 232.90122314, 181.77767193, 217.4720267 , 205.37324983, 209.17398489, 223.88576343, 224.66473346, 205.60042486, 216.65177737, 198.45067478, 215.73994149, 205.049262 , 224.77282497, 214.28533857, 208.01245039, 209.28533857, 191.8810429 , 199.19809327, 220.83178089, 231.61325453, 187.59963247, 191.35779013, 232.12275881, 213.3280424 , 209.99699068, 189.32791597, 224.0838567 , 215.55246097, 201.96711652, 205.95362089, 180.87724119, 183.07048749, 226.97300873, 203.44858428, 208.40337685, 201.02410838, 210.60343418, 205.16834553, 199.07729852, 213.62282203, 198.18275999, 187.79810502, 194.88760107, 181.41595366, 195.53582958, 203.08281145, 188.23921231, 200.05319014, 213.31049219, 220.01061274, 227.93607068, 201.79496927, 257.74086031, 208.17202082, 186.76221222, 204.96711652, 187.88091647, 219.24668931, 230.79496927, 195.56858663, 223.18288642, 194.63343477, 186.46233277, 200.44858428, 206.97981976, 203.43245861, 187.00760342, 201.61876746, 202.31546558, 190.83951075, 204.93882715, 202.01466731, 196.90895299, 184.13362441, 189.50085561, 214.62282203, 193.70024698, 206.17018318, 199.28717621, 212.92754844, 188.14594837, 208.44399018, 211.46404399, 237.98466672, 163.84540296, 171.76601394, 207.71361618, 228.56846021, 208.50896474, 191.87790715, 197.91313398, 193.98821558, 216.36656523, 194.54853283, 202.96515245, 216.92244862, 208.5321543 , 176.71545383, 191.40141278, 201.70889565, 208.31337508, 206.41386316, 204.39944871, 194.25809444, 225.00601864, 211.00196407, 188.53687483, 195.36251066, 202.03643233, 185.89780071, 214.57146953, 169.81147423, 204.96895416, 184.79313163, 184.8810429 , 200.29870777, 170.37166505, 204.76956279, 196.23804063, 185.56833378, 216.38804357, 200.35699773, 170.555091 , 185.96711652, 199.31834847, 187.61588457, 175.9127547 , 193.16454382, 217.85313281, 190.23373321, 205.41817059, 187.94760225, 194.01035989, 187.53136191, 195.1014069 , 204.8810429 , 164.28730264, 188.10717268, 199.24865338, 199.0742892 , 197.87749404, 197.42434948, 181.70889565, 175.32018612, 176.09760519, 171.23110317, 196.31442033, 189.93122372, 192.30288877, 184.76484226, 202.15785922, 184.24263474, 189.20804004, 187.05895593, 149.40153921, 203.50728735, 199.11214607, 187.96135074, 185.5530005 , 191.88614272, 202.09117344, 187.89152921, 160.84920467, 172.00367528, 192.38804357, 164.02607245, 188.7099409 , 181.71152569, 160.08135309, 198.33643821, 150.33852871, 191.51067596, 184.28638382, 174.36355591, 177.96711652, 181.08711887, 179.26281497, 177.78356414, 181.7349681 , 190.18000353, 171.67718385, 163.46771927, 168.55600982, 174.7032563 , 209.5571815 , 197.54919879, 175.03956808, 162.50950428, 182.78945635, 173.93435947, 157.69527359, 169.03944165, 177.41123313, 160.22533739, 179.92073741, 180.18551645, 178.83466378, 184.86767369, 172.15931758, 154.75527476, 178.21865279, 169.45669342, 155.24551763, 169.13926377, 177.47085502, 182.8810429 , 172.60932639, 163.41018788, 157.92232219, 178.56883949, 173.24133663, 189.16885125, 170.49116168, 158.3055188 , 169.04558672, 167.56202846, 176.54435183, 178.37809679, 163.20594955, 174.3393211 , 152.35607891, 143.06405574, 170.1865617 , 168.0338023 , 164.30643762, 186.97208991, 177.23975185, 154.57460527, 168.46914381, 168.85760049, 167.79091471, 155.90356649, 183.00405457, 160.22220164, 170.05319014, 159.72644585, 159.00091882, 175.43220576, 162.08791127, 170.8810429 , 148.8810429 , 167.36276352, 169.18785981, 156.39656582, 163.57238835, 162.885637 , 174.8810429 , 159.62374085, 181.07809091, 194.37074623, 172.4515936 , 174.8810429 , 160.63172356, 167.90410603, 161.10729911, 174.61049807, 179.2474817 , 184.95888096, 167.29908706, 176.28926671, 146.26737526, 175.95704331, 157.61233571, 141.59085737, 159.97760284, 153.75394284, 157.79354474, 163.90461174, 148.33146482, 167.92349387, 160.29018553, 160.20084973, 181.2136794 , 157.50804592, 136.37834965, 187.73012114, 145.81828526, 144.24894006, 167.16467025, 166.0020905 , 154.93214254, 161.70233748, 162.15238012, 182.53829937, 147.12760578, 156.25876041, 157.79496927, 174.19546323, 144.53975772, 178.75644644, 160.52208109, 145.96109788, 149.05619946, 181.88367293, 152.1698965 , 152.40312399, 162.36159184, 157.13926377, 159.49220693, 135.58692923, 152.7850225 , 163.59817411, 163.2858781 , 158.47449648, 169.07913616, 142.56269442, 177.72448178, 145.22533739, 182.70195819, 181.34258328, 160.6314707 , 147.45932345, 143.547614 , 142.10194644, 149.69552645, 150.54289347, 152.61195643, 148.18301285, 139.92336744, 155.1439843 , 145.8810429 , 141.38804357, 157.5125136 , 150.7370586 , 142.31141101, 161.28245568, 145.82630179, 133.39589985, 154.36539355, 153.31768251, 161.47649437, 161.07821734, 146.45460292, 117.98821558, 150.92884654, 152.98428744, 164.63853459, 160.30355473, 134.22926553, 142.57748816, 156.08673959, 148.01178442, 153.05385611, 164.43755843, 138.97171063, 149.4382244 , 159.22533739, 124.81933051, 136.82839229, 154.91785452, 149.45865749, 126.89349328, 149.31939372, 150.43470936, 147.19650848, 152.6064435 , 151.97522566, 135.85250067, 131.02448766, 136.90974538, 157.950738 , 161.47136073, 160.63501956, 143.02055952, 167.39748464, 152.4995575 , 140.79522213, 142.96699009, 137.32816882, 133.63920055, 145.93398019, 129.34099849, 158.64274941, 149.12275881, 163.13926377, 134.86453794, 141.77478904, 153.11908353, 140.19980448, 126.97614448, 149.90122314, 153.30393402, 130.72907589, 135.42146659, 155.4482976 , 116.65808269, 168.4962615 , 139.53294669, 156.52404516, 128.65950722, 136.31651083, 164.44177325, 149.40128635, 154.52534327, 125.95190967, 153.70892948, 166.69749052, 145.22400546, 141.95704331, 148.32908764, 152.31141101, 159.33669107, 145.48355826, 159.11398371, 142.85313281, 144.34049278, 147.02556674, 147.9682882 , 154.25441916, 141.22533739, 135.22533739, 143.85196113, 142.9770633 , 130.30146424, 144.87606951, 151.24944577, 144.11474228, 131.75752551, 151.67442739, 141.28730264, 131.86729441, 130.4578651 , 152.92925965, 162.36618594, 153.45944988, 130.63159713, 152.80374438, 146.22312046, 141.06418217, 128.38870954, 151.22533739, 138.22533739, 126.76810443, 143.48355826, 122.48355826, 139.5312693 , 136.28454617, 139.52195467, 157.80700655, 170.78895064, 143.40299756, 140.18142806, 142.85041017, 149.98346122, 119.41452913, 131.95282849, 151.53950487, 136.13926377, 149.11636089, 137.45051453, 128.77190615, 153.85091588, 139.09225251, 149.55363264, 121.28834789, 120.79325806, 135.08949605, 150.72213843, 140.31609772, 136.14097498, 124.21551704, 145.54827997, 144.61300168, 142.36393519, 122.24668931, 135.31207698, 154.73091353, 148.56896592, 128.30054541, 135.07492134, 156.13888448, 129.26946575, 133.91250184, 133.84920467, 139.02135191, 148.8397636 , 141.87857312, 137.05319014, 159.63226309, 115.51156095, 149.0783776 , 138.18614859, 137.24982505, 120.80231984, 148.5388389 , 154.68336274, 140.94050454, 134.90870013, 126.05632589, 162.96293552, 149.833079 , 141.10115404, 123.34847549, 143.29151746, 145.33130457, 131.10428979, 141.38072683, 131.77634 , 127.98428744, 151.4157008 , 143.20712122, 152.43141337, 151.98533269, 147.22533739, 134.09066773, 140.8243039 , 147.66029961, 142.50741378, 150.93407279, 117.29230985, 144.84695392, 133.60815472, 147.92979919, 129.55939842, 137.46429684, 137.03617948, 164.29778895, 135.10872364, 145.63644409, 126.13128106, 141.39184528, 133.05319014, 142.12564171, 134.9398724 , 126.59253476, 135.55471172, 160.86871894, 115.66932757, 129.72945517, 143.29493988, 152.55192143, 146.50256682, 137.69682455, 128.86048338, 126.22533739, 134.31834847, 137.73458881, 133.84172767, 147.56887331, 143.88154861, 158.09944283, 134.00693746, 151.75514834, 153.35123195, 136.5166946 , 127.58968569, 137.02645173, 153.76237248, 124.53620887, 138.13926377, 152.47687366, 155.6490209 , 160.27238247, 162.58379348, 156.35043956, 158.19913852, 164.29724942, 140.1055879 , 152.24029138, 146.1055879 , 136.64822851, 143.74257152, 160.13926377, 179.3129958 , 150.22137543, 150.73984888, 158.42719854, 140.59662315, 132.42447591, 155.44126754, 169.65491311, 146.4185837 , 153.15535561, 137.04308311, 163.62193703, 147.71057304, 169.05908235, 144.60004558, 144.81017612, 150.62871424, 191.72068007, 151.05319014, 138.05319014, 149.05319014, 128.9019817 , 149.6051454 , 150.34366235, 155.83963717, 168.00915439, 160.78644703, 157.98352887, 175.42247801, 138.3501867 , 153.9794743 , 175.48824497, 161.3501867 , 142.71032018, 130.97421423, 147.39890917, 174.70747112, 154.6944812 , 150.32135779, 128.51685485, 163.57957866, 166.16809268, 158.68047985, 167.20014994, 176.77691336, 150.20056305, 149.39343007, 172.53899916, 151.50788567, 140.78686014, 152.18025638, 164.02252359, 135.21501133, 139.0838567 , 121.50829878, 155.0162521 , 157.9486475 , 132.8647908 , 154.41373673, 152.13926377, 175.38123254, 142.31141101, 143.93201612, 146.57001117, 137.95302257, 163.58518419, 139.51682102, 166.52341302, 159.17068889, 163.24013113, 162.34796978, 156.57146953, 175.39932228, 170.66682396, 145.74545441, 169.01653878, 164.43781129, 153.30589809, 145.31141101, 136.02949488, 167.22349975, 176.25270794, 140.166921 , 154.8810429 , 156.94129693, 154.3664388 , 161.77098733, 156.11528182, 162.33539296, 154.34245685, 171.12418335, 154.34157185, 164.4657552 , 157.8527197 , 188.77282497, 152.41072742, 172.49245979, 143.97396138, 190.8013672 , 142.87853929, 156.43926965, 158.58354063, 176.14749933, 147.43103408, 175.09814472, 167.21459822, 171.26579047, 161.31141101, 158.30133781, 152.4828923 , 154.07466849, 153.31157127, 176.09614683, 167.50453089, 171.33288936, 166.67718385, 154.69815648, 174.9444665 , 164.70889565, 171.10755197, 167.8810429 , 155.72432153, 154.8255094 , 159.13692041, 166.35620534, 166.00956749, 171.53557673, 192.35503366, 172.10520861, 156.80149363, 175.31445416, 152.52718091, 168.15418394, 154.22533739, 142.26360737, 156.18706741, 164.06275764, 154.89599689, 172.37834965, 170.98504601, 157.77703979, 195.69514716, 161.74596013, 153.7925921 , 157.99280968, 166.24146306, 187.81828526, 148.09982211, 190.09855783, 185.85243302, 161.20857958, 172.05319014, 185.51343242, 173.84461056, 175.51343242, 188.73522096, 168.00709771, 190.58585016, 200.11648732, 173.54029726, 170.16558907, 178.50988356, 191.25467201, 174.87857312, 157.4536841 , 185.62583135, 180.79797859, 179.97012584, 191.97291613, 177.47474934, 184.60476611, 181.32344829, 185.22834671, 164.48431683, 156.62206346, 179.5726412 , 194.13926377, 181.52246038, 175.31236366, 208.85547616, 214.117659 , 175.85246684, 193.44335804, 184.45498221, 180.67300285, 201.8810429 , 170.94852107, 171.02290288, 170.75648027, 181.52074916, 167.58822734, 204.82915086, 168.50126872, 178.08737173, 176.44896357, 172.31141101, 207.65741672, 210.55233454, 202.22533739, 185.64011937, 192.46797213, 204.78109436, 206.42121373, 198.23921231, 184.61537886, 186.31071123, 173.76379701, 171.33615153, 171.95796214, 175.49914439, 194.07520802, 181.58521802, 200.94066479, 182.58521802, 199.2520758 , 181.547614 , 186.18272617, 191.13812591, 166.72590632, 194.82642822, 215.49442386, 211.44452972, 174.18773338, 200.35988063, 174.21732086, 187.37261769, 180.70417512, 189.8668813 , 186.56979214, 217.5367484 , 222.28324807, 212.46939666, 183.11407632, 193.79199378, 187.01400135, 195.76069509, 175.17353796, 169.83102232, 205.21903207, 173.90762106, 199.99689808, 192.48038869, 198.64943401, 210.94050454, 172.62909352, 197.96084503, 238.74954281, 201.84039574, 185.99214372, 184.04219812, 182.37559319, 231.95308135, 186.34743024, 198.45853106, 193.11873807, 189.90942488, 191.50051015, 167.1394828 , 209.57109024, 182.1343242 , 205.44921642, 196.21602275, 198.74817705, 183.38782454, 208.99068536, 201.17674135, 214.70610536, 192.87318662, 173.70103937, 171.52889212, 191.96385435, 196.57748816, 214.66948782, 210.91126252, 210.85867956, 222.07023463, 192.89554996, 203.77792478, 198.06443503, 195.85855313, 189.03735116, 210.31020551, 209.50845903, 184.24117638, 196.03843023, 206.05549967, 193.21288701])

Normalizing by (that is dividing by) the monitor would introduce correlations between different detector pixels if the monitor has variances. This is because the same monitor bin is applied to many detector bins. Scipp would reject such an operation. To work around this, we drop the variances of the monitor. In practice, we have to carefully examine the uncertainties to find out if they really can be neglected!

[17]:
mon = sc.values(mon)

We intend to normalize each event to the relative monitor counts (compared to the total monitor counts). We use sum to compute the total monitor counts and obtain the relative counts using division:

[18]:
mon /= mon.sum('wavelength')

The sample data is event data in bins and the monitor is a histogram. Multiplication and division operations for such cases are supported by modifying the weights (values) for each event using the operators of the bins property, in combination with the sc.lookup helper, a wrapper for a discrete “function”, given by the monitor:

[19]:
sample_over_mon = sample_lambda.bins / sc.lookup(func=mon, dim='wavelength')
sample_over_mon
[19]:
Show/Hide data repr Show/Hide attributes
scipp.DataArray (412.88 MB)
    • spectrum: 24794
    • wavelength: 1
    • L1
      ()
      float64
      m
      60.0
      Values:
      array(60.)
    • L2
      (spectrum)
      float64
      m
      2.560, 2.556, ..., 4.251, 4.252
      Values:
      array([2.56027902, 2.55590928, 2.55268677, ..., 4.25113991, 4.25116338, 4.25188 ])
    • Ltotal
      (spectrum)
      float64
      m
      62.560, 62.556, ..., 64.251, 64.252
      Values:
      array([62.56027902, 62.55590928, 62.55268677, ..., 64.25113991, 64.25116338, 64.25188 ])
    • incident_beam
      ()
      vector3
      m
      [ 0. 0. 60.]
      Values:
      array([ 0., 0., 60.])
    • position
      (spectrum)
      vector3
      m
      [ 1.17451004 -1.01106149 -2.03796699], [ 1.18147634 -0.95946649 -2.05334117], ..., [1.81428985 0.09565841 3.84338287], [1.81375055 0.1499371 3.84269584]
      Values:
      array([[ 1.17451004, -1.01106149, -2.03796699], [ 1.18147634, -0.95946649, -2.05334117], [ 1.18844265, -0.90787149, -2.06871534], ..., [ 1.81482915, 0.04137972, 3.8440699 ], [ 1.81428985, 0.09565841, 3.84338287], [ 1.81375055, 0.1499371 , 3.84269584]])
    • sample_position
      ()
      vector3
      m
      [0. 0. 0.]
      Values:
      array([0., 0., 0.])
    • scattered_beam
      (spectrum)
      vector3
      m
      [ 1.17451004 -1.01106149 -2.03796699], [ 1.18147634 -0.95946649 -2.05334117], ..., [1.81428985 0.09565841 3.84338287], [1.81375055 0.1499371 3.84269584]
      Values:
      array([[ 1.17451004, -1.01106149, -2.03796699], [ 1.18147634, -0.95946649, -2.05334117], [ 1.18844265, -0.90787149, -2.06871534], ..., [ 1.81482915, 0.04137972, 3.8440699 ], [ 1.81428985, 0.09565841, 3.84338287], [ 1.81375055, 0.1499371 , 3.84269584]])
    • source_position
      ()
      vector3
      m
      [ 0. 0. -60.]
      Values:
      array([ 0., 0., -60.])
    • spectrum
      (spectrum)
      int32
      1, 2, ..., 24793, 24794
      Values:
      array([ 1, 2, 3, ..., 24792, 24793, 24794], dtype=int32)
    • tof
      (wavelength [bin-edge])
      float64
      µs
      19.0, 1.669e+04
      Values:
      array([ 19. , 16694.30078125])
    • two_theta
      (spectrum)
      float64
      rad
      2.491, 2.504, ..., 0.442, 0.442
      Values:
      array([2.49144445, 2.50372965, 2.51564288, ..., 0.44118918, 0.44157918, 0.44231324])
    • wavelength
      (spectrum, wavelength [bin-edge])
      float64
      Å
      0.001, 1.056, ..., 0.001, 1.028
      Values:
      array([[0.00120148, 1.05567339], [0.00120156, 1.05574713], [0.00120162, 1.05580152], ..., [0.00116986, 1.02789183], [0.00116986, 1.02789145], [0.00116984, 1.02787999]])
    • (spectrum, wavelength)
      DataArrayView
      binned data [len=0, len=0, ..., len=0, len=0]
      dim='event',
      content=DataArray(
                dims=(event: 17926980),
                data=float32[counts],
                coords={'tof':float64[µs], 'wavelength':float64[Å]})

Finally, we can plot the data, which needs to be histogrammed before being plotted. By default, .hist() uses the coordinates of the binned data to define histogram edges, which, in this case, would give a single bin along the 'wavelength' dimension. For a better representation of the data, we supply a finer binning, yielding a more meaningful figure:

[20]:
sc.plot(sample_over_mon.hist(wavelength=400).hist(spectrum=500).transpose())
[20]:
../_images/tutorials_powder-diffraction_35_0.svg
[21]:
del sample_lambda
del sample_over_mon
del sample
del vanadium

From events to histogram#

Note: We are continuing here with data that has not been normalized to the monitors.

We histogram the event data:

[22]:
dspacing = sc.arange(dim='dspacing', start=0.3, stop=2.0, step=0.001, unit='Angstrom')
hist = sc.Dataset(
    data={
        'sample': dspacing_sample.hist(dspacing=dspacing),
        'vanadium': dspacing_vanadium.hist(dspacing=dspacing),
    }
)
sc.show(hist['spectrum', 0:3]['dspacing', 0:7])
samplesample(dims=('spectrum', 'dspacing'), shape=(3, 7), unit=counts, variances=True)variances spectrumdspacingvalues spectrumdspacing vanadiumvanadium(dims=('spectrum', 'dspacing'), shape=(3, 7), unit=counts, variances=True)variances spectrumdspacingvalues spectrumdspacing two_t..two_theta(dims=('spectrum',), shape=(3,), unit=rad, variances=False)values spectrum spect..spectrum(dims=('spectrum',), shape=(3,), unit=None, variances=False)values spectrum scatt..scattered_beam(dims=('spectrum',), shape=(3,), unit=m, variances=False)values spectrum posit..position(dims=('spectrum',), shape=(3,), unit=m, variances=False)values spectrum Ltota..Ltotal(dims=('spectrum',), shape=(3,), unit=m, variances=False)values spectrum L2L2(dims=('spectrum',), shape=(3,), unit=m, variances=False)values spectrum sourc..source_position(dims=(), shape=(), unit=m, variances=False)values sampl..sample_position(dims=(), shape=(), unit=m, variances=False)values incid..incident_beam(dims=(), shape=(), unit=m, variances=False)values L1L1(dims=(), shape=(), unit=m, variances=False)values dspacingdspacing(dims=('dspacing',), shape=(8,), unit=Å, variances=False)values dspacing
[23]:
sc.plot(hist['sample'].hist(spectrum=500).transpose())
[23]:
../_images/tutorials_powder-diffraction_39_0.svg

Summing (focussing) and normalizing#

After conversion to 'dspacing', generic sum and / operations can be used to “focus” and normalize the diffraction data to the vanadium run:

[24]:
summed = sc.sum(hist, 'spectrum')
sc.plot(summed)
[24]:
../_images/tutorials_powder-diffraction_41_0.svg
[25]:
normalized = summed['sample'] / summed['vanadium']
sc.plot(normalized)
[25]:
../_images/tutorials_powder-diffraction_42_0.svg

Focussing with \(\theta\) dependence in event-mode#

The approach used above combines reflections from all crystallographic planes and is therefore of limited use. We can use bin to focus each of multiple groups of spectra into a distinct output spectrum. Here we define groups based on a range of scattering angles – a simple \(\theta\)-dependent binning. This also demonstrates how we can postpone histogramming until after the focussing step.

[26]:
two_theta = sc.linspace(dim='two_theta', unit='rad', start=0.0, stop=np.pi, num=16)

focussed_sample = dspacing_sample.bin(two_theta=two_theta)
focussed_vanadium = dspacing_vanadium.bin(two_theta=two_theta)
norm = focussed_vanadium.hist(dspacing=dspacing)

Similarly as when normalizing by monitor, we have to drop the variances of the normalization term. Otherwise, the normalization would be broadcast to focussed_sample which would introduce correlations. In practice, we have to make sure that the uncertainties of the vanadium measurement can be neglected!

[27]:
norm = sc.values(norm)
[28]:
focussed_sample.bins /= sc.lookup(func=norm, dim='dspacing')
normalized = focussed_sample.hist(dspacing=dspacing)

The normalized output looks as follows:

[29]:
normalized.plot(vmin=sc.scalar(0), vmax=sc.scalar(0.4))
[29]:
../_images/tutorials_powder-diffraction_49_0.svg

As a bonus, we can use slicing and a dict-comprehension to quickly create of plot comparing the spectra for different scattering angle bins:

[30]:
# compute centers of theta bins
angle = normalized.coords['two_theta'].values
angle = 0.5 * (angle[1:] + angle[:-1])
results = {
    f'{round(angle[group], 3)} rad': normalized['dspacing', 300:500]['two_theta', group]
    for group in range(2, 6)
}
sc.plot(results)
[30]:
../_images/tutorials_powder-diffraction_51_0.svg