Fitting SANS data

Fitting SANS data#

Previously, some small angle neutron scattering (SANS) data has been simulated and reduced, and can now be analysed with easyCore. Before the analysis can begin, it is necessary to load the experimental data and check that it looks reasonable. The data can be loaded with np.loadtxt as the data has been stored in a simple space-separated column file.

import numpy as np

q, i, di = np.loadtxt('../4-reduction/sans_iofq_3pulses.dat', unpack=True)

---------------------------------------------------------------------------
FileNotFoundError                         Traceback (most recent call last)
Cell In[1], line 3
      1 import numpy as np
----> 3 q, i, di = np.loadtxt('../4-reduction/sans_iofq_3pulses.dat', unpack=True)

File /opt/hostedtoolcache/Python/3.10.13/x64/lib/python3.10/site-packages/numpy/lib/npyio.py:1373, in loadtxt(fname, dtype, comments, delimiter, converters, skiprows, usecols, unpack, ndmin, encoding, max_rows, quotechar, like)
   1370 if isinstance(delimiter, bytes):
   1371     delimiter = delimiter.decode('latin1')
-> 1373 arr = _read(fname, dtype=dtype, comment=comment, delimiter=delimiter,
   1374             converters=converters, skiplines=skiprows, usecols=usecols,
   1375             unpack=unpack, ndmin=ndmin, encoding=encoding,
   1376             max_rows=max_rows, quote=quotechar)
   1378 return arr

File /opt/hostedtoolcache/Python/3.10.13/x64/lib/python3.10/site-packages/numpy/lib/npyio.py:992, in _read(fname, delimiter, comment, quote, imaginary_unit, usecols, skiplines, max_rows, converters, ndmin, unpack, dtype, encoding)
    990     fname = os.fspath(fname)
    991 if isinstance(fname, str):
--> 992     fh = np.lib._datasource.open(fname, 'rt', encoding=encoding)
    993     if encoding is None:
    994         encoding = getattr(fh, 'encoding', 'latin1')

File /opt/hostedtoolcache/Python/3.10.13/x64/lib/python3.10/site-packages/numpy/lib/_datasource.py:193, in open(path, mode, destpath, encoding, newline)
    156 """
    157 Open `path` with `mode` and return the file object.
    158 
   (...)
    189 
    190 """
    192 ds = DataSource(destpath)
--> 193 return ds.open(path, mode, encoding=encoding, newline=newline)

File /opt/hostedtoolcache/Python/3.10.13/x64/lib/python3.10/site-packages/numpy/lib/_datasource.py:533, in DataSource.open(self, path, mode, encoding, newline)
    530     return _file_openers[ext](found, mode=mode,
    531                               encoding=encoding, newline=newline)
    532 else:
--> 533     raise FileNotFoundError(f"{path} not found.")

FileNotFoundError: ../4-reduction/sans_iofq_3pulses.dat not found.

With the data read in, we can produce a quick plot.

import matplotlib.pyplot as plt

plt.errorbar(q, i, di)
plt.yscale('log')
plt.xlabel('$q$/Å')
plt.ylabel('I(q)')
plt.show()

We now want to consider the mathematical model to be used in the analysis. There are SANS models a range of different systems, see for instance the models in SasView. However, initially, we will assume that our data has arisen from a spherical system.

The mathematical model for a sphere is

(3)#\[ I(q) = \frac{\text{scale}}{V} \bigg(\frac{3 V \Delta \rho [\sin{(qr)} - qr \cos{(qr)}]}{(qr)^3}\bigg)^2 + \text{bkg}, \]

where \(\text{scale}\) is a scale factor, \(V\) is the volume of the sphere, \(\Delta \rho\) is the difference between the solvent and particle scattering length density, \(r\) is the radius of the sphere, \(\text{bkg}\) is a uniform background, and \(q\) is the q-vector that the intensity is being calculated for.

Task

The mathematical model described in Eqn. (3) has five parameters. What simple mathematical simplification can be performed to reduce this to four?

Click to show solution

The volume of a sphere is related to the radius of the sphere as

(4)#\[ V = \frac{4}{3} \pi r^3. \]

Therefore, the parameter \(V\) can be replaced with Eqn. (4).

Task

Four parameters is a suitable number for modelling. Therefore, we should write a function that implements your reduced dimensionality version of Eqn. (3).

Click below to show code solution

Task

Create four Parameter objects, for \(\text{scale}\), \(\Delta \rho\), \(r\), and \(\text{bkg}\). Each should have an initial value and a uniform prior distribution based on the values given in Table 1, except for the \(\text{scale}\) which should be fixed to a value of 1.4 × 10^-8.

Table 1 Parameter values for the spherical model.#
Parameter	Initial Value	Min	Max
\(\text{scale}\)	1.4 × 10^-8	N/A	N/A
\(\Delta \rho\)	3	0	10
\(r\)	80	10	1000
\(\text{bkg}\)	3.0 × 10^-3	1.0 × 10^-3	1.0 × 10^-2

Click below to show code solution

It is now possible to compare our model at the initial estimates to the simulated data.

plt.errorbar(q, i, di, marker='.', ls='', color='C0')
plt.plot(q, sphere(q), 'k', zorder=10)
plt.yscale('log')
plt.xlabel('$q$/Å')
plt.ylabel('I(q)')
plt.show()

Task

Using easyCore, obtain maximum likelihood estimates for the four parameters of the model from comparison with the data.

Click below to show code solution

We can then plot the model and the data together as before and print the values of the parameters along with their uncertainties.

plt.errorbar(q, i, di, marker='.', ls='', color='C0')
plt.plot(q, sphere(q), 'k-', zorder=10)
plt.yscale('log')
plt.xlabel('$q$/Å')
plt.ylabel('I(q)')
plt.show()

delta_rho.value

r.value

bkg.value