==============================================
XAFS: Post-edge Background Subtraction
==============================================
.. :noindex: module:: larch.xafs
Background subtraction is one of the most important data processing steps
in EXAFS analysis, converting the measured :math:`\mu(E)` into the
:math:`\chi(k)` ready for quantitative analysis. In Larch, this step is
performed by the :func:`autobk` function, which has many options and
subtleties. This section is devoted to background subtraction with the
:func:`autobk` function and the underlying algorithm it uses.
The :func:`autobk` function
=============================
.. function:: autobk(energy, mu=None, group=None, rbkg=1.0, ...)
Determine the post-edge background function, :math:`\mu_0(E)`, and
corresponding :math:`\chi(k)`.
:param energy: 1-d array of x-ray energies, in eV
:param mu: 1-d array of :math:`\mu(E)`
:param group: output group (and input group for ``e0`` and ``edge_step``).
:param rbkg: distance (in :math:`\rm\AA`) for :math:`\chi(R)` above
which the signal is ignored. Default = 1.
:param e0: edge energy, in eV. If `None`, it will be determined.
:param edge_step: edge step. If `None`, it will be determined.
:param pre_edge_kws: dictionary containing keyword arguments to pass to :func:`pre_edge`.
:param nknots: number of knots in spline. If `None`, it will be determined.
:param kmin: minimum :math:`k` value [0]
:param kmax: maximum :math:`k` value [``None``, full data range].
:param kweight: :math:`k` weight for FFT. [1]
:param dk: FFT window window parameter. [0.1]
:param win: FFT window function name. ['hanning']
:param nfft: array size to use for FFT [2048]
:param kstep: :math:`k` step size to use for FFT [0.05]
:param k_std: optional :math:`k` array for standard :math:`\chi(k)` [``None``]
:param chi_std: optional :math:`\chi` array for standard :math:`\chi(k)` [``None``]
:param nclamp: number of energy end-points for clamp [2]
:param clamp_lo: weight of low-energy clamp [1]
:param clamp_hi: weight of high-energy clamp [1]
:param calc_uncertaintites: Flag to calculate uncertainties in :math:`\mu_0(E)` and :math:`\chi(k)` [``True``]
:param err_sigma: sigma level for uncertainties in mu_0(E) and chi(k) [1]
:returns: ``None``.
Follows the First Argument Group convention, using group members named ``energy`` and ``mu``.
The following data is put into the output group:
================= ===============================================================
attribute meaning
================= ===============================================================
bkg array of :math:`\mu_0(E)` (not normalized)
chie array of :math:`\chi(E)` values.
k array of :math:`k` values, on uniform grid.
chi array of :math:`\chi(k)` values, the EXAFS.
*delta_chi* array of uncertainty in :math:`\chi(k)`.
*delta_bkg* array of uncertainty in :math:`\mu_0(E)`.
autobk_details Group of arrays with autobk details.
================= ===============================================================
Here, the arrays ``group.k``, ``group.chi`` will be the same length,
giving :math:`\chi(k)` from 0 to a maximum k value set by ``kmax`` or
the range of available data. The arrays ``group.bkg`` and
``group.chie`` will correspond to the input ``energy`` array.
The parameters ``e0`` and ``edge_step`` are needed to determine
:math:`\chi(k)`, and can be specified here. If they are not, but the
input ``group`` has attributes ``e0`` and ``edge_step`` -- say, from
running :func:`pre_edge` -- those values will be used. If values
cannot be found in either place, they will be determined by calling
:func:`pre_edge` from within :func:`autobk`, and the output parameters
from :func:`pre_edge` will be added to the output ``group`` here.
The ``pre_edge_kws`` argument can be used here to hold a dictionary of
parameters to pass to :func:`pre_edge`.
If ``calc_uncertainties`` is set to ``True``, the outputs
``group.delta_chi`` and ``group.delta_bkg``, holding the uncertainties
in :math:`\chi(k)` and :math:`\mu_0(E)`, respectively. These
calculations are relatively costly, and can take 20 times longer to
determine :math:`\chi(k)` and :math:`\mu_0(E)` themselves.
The ``group.autobk_details`` group will contain the following attributes:
================= ===========================================================
attribute meaning
================= ===========================================================
spline_pars Parameters used for determining :math:`\mu_0(E)` spline.
init_bkg Initial value for :math:`\mu_0(E)`
init_chi Initial value for :math:`\chi(k)`
knots_e Spline knot energies
knots_y Spline knot values
init_knots_y Initial Spline knot values
================= ===========================================================
The AUTOBK Algorithm
======================
The background subtraction method used is the **AUTOBK** algorithm
(see :cite:`autobk`), in which a spline function is matched to the
low-*R* components of the resulting :math:`\chi(k)`.
For reference, :math:`k = \sqrt{2m_e (E-E_0)/\hbar^2}` is the wavenumber of
the ejected photo-electron, where :math:`E_0` is the absorption threshold
energy (the 0 of photo-electron energy). For :math:`k` in units of
:math:`\rm \AA^{-1}` and :math:`E` in units of eV, :math:`k \approx
\sqrt{(E-E_0)/3.81}`. With this conversion of energy to wavenumber,
:math:`\chi(k)` is defined from
:math:`\chi(E) = \frac{\mu(E)-\mu_0(E)}{\Delta\mu}`
where :math:`\mu_0(E)` represents the idealized x-ray absorption of a
bare atom, embedded in the molecular or solid environment, but without
scattering of the outgoing photo-electron that gives rise to the EXAFS, and
:math:`\Delta\mu` is the edge step in :math:`\mu(E)`.
The quantity :math:`\mu_0(E)` cannot be independently measured. Instead,
we empirically determine it here given the dxta for :math:`\mu(E)` by
fitting a spline (a piece-polynomial function that is easily adjusted even
while its smoothness is controlled) to :math:`\mu(E)`. An advantage to
what could easily be described as an *ad hoc* approach is that the data for
:math:`\mu(E)` can include many systematic drifts and dependencies that are
slowly varying with energy.
The definition above normalizes to the edge step :math:`\Delta\mu`
(``edge_step``) instead of the energy-dependent :math:`\mu_0(E)`, as is
often described in introductory texts on EXAFS. The reason for this is
essentially the same -- so that we do not need to carefully take care of
slow energy drifts in the measured :math:`\mu(E)` or worry about having an
absolute measure of :math:`\mu(E)`.
:math:`R_{\rm bkg}` and spline for :math:`\mu_0(E)`
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
A spline is a general mathematical function, and so using one to match
:math:`\mu(E)` that we will subtract from the same :math:`\mu(E)` could
easily match *too* well, and erase much of the data we're interested in.
Therefore, we have to carefully consider two aspects:
1. what portions of the spectrum we want it to match, and
2. how flexible the spline can be.
What we want is a :math:`\mu_0(E)` that removes the low frequency (low-*R*)
portions of the :math:`\chi` spectra, and recognizing this can tell us how
to determine both of these.
First, since we know the EXAFS oscillations do not extend far below the
near-neighbor distance and since atoms are essentially never closer than 1
:math:`\rm\AA`, and generally more like 2 :math:`\rm\AA`, we can say that
we want to remove the frequencies of :math:`\chi` below some distance,
:math:`R_{\rm bkg}` (``rbkg``). That is, we want the spline function to
be adjusted so that the components of :math:`\chi(R)` below :math:`R_{\rm
bkg}` are minimized. We can use a default of 1 :math:`\rm\AA`, and
recommend that the value be roughly half the expected near-neighbor
distance.
To do this, we create a spline function consisting of :math:`N_{\rm knots}`
adjustable points (*knots* in the spline, where the second derivatives are
allowed to change), evenly spaced in :math:`k` (that is, quadratic in
energy). We do a fit that adjusts the values of the function at each knot,
computes :math:`\mu_0(E)` by simple spline interpolation, calculates
:math:`\chi(k)`, does a Fourier transform to :math:`\chi(R)`, and then
seeks to minimize the components of :math:`\chi(R)` below :math:`R_{\rm
bkg}`.
This leads to the second question of how flexible to make the spline is
answered from basic signal processing through the `Nyquist-Shannon sampling
theorem `_,
which tells us how many adjustable parameters to use for our spline:
.. math::
N_{\rm knots} = 1 + \frac{2R_{\rm bkg} \Delta k}{\pi}
where :math:`\Delta k` is the :math:`k` range of the data. This value can
be over-ridden by specifing ``nknots``, but this is strongly recommended
against.
In this way we satisfy both considerations above: limiting the number of
knots prevents the spline from being able to follow the frequencies of
:math:`\chi(k)` we care about most, and limiting the portion of the
spectrum to be minimized to the very low-:math:`R` components, we do not
even try to match those frequencies. Thus, we are ensured of a
:math:`\mu_0(E)` that only removes the low-:math:`R` components of
:math:`\chi`.
Selecting Fourier transform parameters
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
In order to determine :math:`\mu_0(E)`, we must do a Fourier transform.
These are discussed in more detail in the next section, but here we give
the parameters affecting the Fourier transform used, and their default
values.
========== ====================== ============================================
argument meaning default, recommended values
========== ====================== ============================================
kmin :math:`k_{\rm min}` 0. should be below 1.0
kmax :math:`k_{\rm max}` highest :math:`k` value. Useful data range.
kweight :math:`kw` 1. should be 0 or 1 (but not > 2)
win window function name Hanning, Parzen works well too.
dk window sill size 0. should be 0 or 1
nfft FFT array size 2048. No reason to change this.
kstep :math:`k` step size 0.05. No reason to change this.
========== ====================== ============================================
The most important point is that the Fourier transform used to fit
and visualize :math:`\chi(k)` is not necessarily the one to use here. In
particular, the ``kweight`` parameter should be 0 or 1 here, as we want to
emphasize the low-:math:`k` portion of :math:`\chi(k)`. The default
parameter values for :func:`autobk` are generally sufficient, and need only
minor adjustments except in unusual cases.
Using a standard :math:`\chi(k)`
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
We said above that we want to minimize the low-:math:`R` components of
:math:`\chi`. In principle, there should be some leakage of the first
shell to fairly low-:math:`R`. This is especially noticeable for a short
ligand (say, below 1.75 :math:`\rm\AA` or so) of a low-Z atom (notably, C,
N, and O). Normally, this is not a serious problem, but it does point out
that :math:`\mu_0(E)` might want to leave some first-shell leakage at
low-:math:`R`.
To accout for this, you can proved a spectrum of :math:`\chi(k)` for a
standard that is meant to be close to the spectrum being analyzed. This is
done by providing 2 arrays: ``k_std`` and ``chi_std``, which need to be the
same lenth. By providing these, the best values for :math:`\mu_0(E)` will
be those that minimize the Fourier transorm of the difference of
:math:`\chi(k)` for the data and standard. This can be especially helpful
to give more consistent background for a set of related spectra.
End-point clamps for the spline
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
Because the spline is chosen to match the low-:math:`R` components of
:math:`\chi`, the process does not even look at how the resulting
:math:`\mu_0(E)` looks. This can lead to some unstable behavior at the
end-points of the spline, causing :math:`\chi(k)` to diverge at the ends.
This can be remedied by adding *clamps* to the end-points of :math:`chi(k)`
to the fit. That is, a few points at the low- and high-ends of the
:math:`\chi(k)` array are appended to the array to be minimized in the fit,
which has the tendency to push the end-points of :math:`\chi(k)` to zero.
There are three parameters influencing these end-points clamps. ``nclamp``
sets the number of data points of :math:`\chi(k)` at the beginning and end
of the :math:`k` range to add to the fit. ``clamp_lo`` sets the weighting
factor applied to the first ``nclamp`` points for the low-:math:`k` end fo
the data range, and ``clamp_hi`` sets the weighting factor applied to the
last ``nclamp`` points at the high-:math:`k` end.
There is not an easy way to determine ahead of time what the clamp
parameters will be. Typically, ``nclamp`` need only be 1 to 5 data points
(with 2 being the default). Setting ``nclamp`` to 0 will remove the clamps
completely. Values for ``clamp_lo`` and ``clamp_hi`` are floating point
numbers, and should range from about 1 (the default) to 10 or 20 for a very
strong clamp that will almost certainly force :math:`\chi(k)` to be 0 at
the end-points.
Recommendations
~~~~~~~~~~~~~~~~~
Here, we give a few recommendations on what parameters most affect the
background subtraction. These parameters are (in roughly increasing
order):
0. ``rbkg``: This is the main parameter that sets how flaexible the
spline function can be.
1. ``e0``: It can be hard to know what :math:`E_0` should be for any XAFS
spectra. As this sets the value of :math:`k=0`, it will affect the
location of the knots, and can have a profound affect on the results. If
you're unhappy with the :math:`\mu_0(E)` and :math:`\chi(k)`, playing with
the ``e0`` parameter.
2. ``kweight``: Increasing this emphasizes the fit at high *k*. Trying
both values of 0 and 1 is always a good idea.
3. ``clamp_lo`` and ``clamp_hi``: If :math:`\mu_0(E)` swoops up or down at
the endpoints too much, increasing these values can help greatly.
Examples
============
Here, we give a few examples using :func:`autobk`, and outline some of its
features.
Basic example
~~~~~~~~~~~~~~~~~~
The simplest example of reading data and doing background subtraction would
using all the default inputs would be:
.. literalinclude:: ../../examples/xafs/doc_autobk1.lar
with the resulting outputs looking like this:
.. subfigstart::
.. _fig_xafs2a:
.. figure:: ../_images/xafs_autobk1a.png
:target: ../_images/xafs_autobk1a.png
:width: 100%
:align: center
:math:`\mu(E)` and background :math:`\mu_0(E)`
.. _fig_xafs2b:
.. figure:: ../_images/xafs_autobk1b.png
:target: ../_images/xafs_autobk1b.png
:width: 100%
:align: center
resulting :math:`k\chi(k)`
.. subfigend::
:width: 0.45
:label: fig-xafs2
Example of simple usage of :func:`autobk` for Cu metal.
Thus demonstrating that we can process data on Cu metal.
For a slightly more challenging example, we try the As K-edge of scorodite
(:math:`\rm FeAsO4 \cdot (nH_2O)`). Here, we start by reading in data and
constructing :math:`\mu(E)`, then run :func:`autobk`. We take most of the
default values for pararmeters, but pass explicit parameters such as
``nnorm`` to :func:`pre_edge` using ``pre_edge_kws``. Setting these
parameters is not critically important for this data set, and is mostly
here to demonstrate the method for passing arguments to :func:`pre_edge`.
.. literalinclude:: ../../examples/xafs/doc_autobk2.lar
The resulting outputs looks OK:
.. subfigstart::
.. _fig-xafs3a:
.. figure:: ../_images/xafs_autobk2a.png
:target: ../_images/xafs_autobk2a.png
:width: 100%
.. _fig-xafs3b:
.. figure:: ../_images/xafs_autobk2b.png
:target: ../_images/xafs_autobk2b.png
:width: 100%
.. subfigend::
:width: 0.45
:label: fig-xafs3
Example usage of :func:`autobk` for :math:`\rm FeAsO4 \cdot (nH_2O)`.
A close examimation of :math:`k^2\chi(k)` suggests we might be able to do
better, which brings us to the next section.
Adjusting spline clamps
~~~~~~~~~~~~~~~~~~~~~~~~~
For the scorodite example above, the XAFS diverges from 0 at the
high-:math:`k` end of the spectra. Playing with the spline clamps
discussed above, we try a few options to study the effect of adding a
clamp:
.. literalinclude:: ../../examples/xafs/doc_autobk3.lar
results in the following comparison:
.. _xafs_fig4:
.. figure:: ../_images/xafs_autobk3.png
:target: ../_images/xafs_autobk3.png
:width: 65%
Influence of spline clamps of :math:`\chi(k)`
Using a Standard :math:`\chi(k)`
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
Since the EXAFS contribution from the first shell is not perfectly isolated
in :math:`R`, it can leak into portion of :math:`\chi(R)` below
:math:`R_{\rm bkg}`. It may be desirable to not simply reduce the portion
of :math:`\chi(R)` below :math:`R_{\rm bkg}`, but to account for the
expected leakage from the first shell. To do this, you can provide a
*standard* :math:`\chi(k)` spectrum that will simulate the first shell by
specify arrays for `k_std` and `chi_std` arrays to the :func:`autobk`
function.
As an example:
.. literalinclude:: ../../examples/xafs/doc_autobk4.lar
.. rubric:: References
.. bibliography:: ../larch.bib
:filter: cited and ({'xafs/autobk'} >= docnames)
:style: authorlist