diff --git a/doc/user_guide/eels.rst b/doc/user_guide/eels.rst index b1fd92540..1aff5e059 100644 --- a/doc/user_guide/eels.rst +++ b/doc/user_guide/eels.rst @@ -409,6 +409,49 @@ edge functionalities: * :py:meth:`~.models.EELSModel.fix_fine_structure` * :py:meth:`~.models.EELSModel.free_fine_structure` +.. _eels.absolute_quantification: + +Absolute quantification +^^^^^^^^^^^^^^^^^^^^^^^ + +The fitted :py:attr:`~.components.EELSCLEdge.intensity` parameter is the +scaling factor applied to the theoretical cross-section (in +barns/eV/atom). When the model is convolved with a **raw** low-loss +spectrum (i.e. in counts, not normalised), the fitted intensity is +directly proportional to the areal density :math:`N` of the element: + +.. math:: + + \text{intensity} = N \times 10^{-10} + +where :math:`N` is in atoms/nm². Therefore, to obtain the areal density: + +.. code-block:: python + + >>> N_B = m.components.B_K.intensity.value * 1e10 + >>> N_N = m.components.N_K.intensity.value * 1e10 + +.. warning:: + + If the core-loss and low-loss spectra were acquired with different + dwell times (or total acquisition times), the fitted intensity must + be corrected by the ratio of the dwell times before converting to + areal density. The dwell time can often be found in the signal + metadata (e.g. + :py:attr:`~hyperspy.api.signals.BaseSignal.metadata`\ ``.Acquisition_instrument.TEM.Detector.EELS.dwell_time``). + The corrected conversion is: + + .. math:: + + N = \text{intensity} \times \frac{t_{\text{core}}}{t_{\text{low}}} \times 10^{10} + + where :math:`t_{\text{core}}` and :math:`t_{\text{low}}` are the + dwell times of the core-loss and low-loss spectra, respectively. + + If the low-loss spectrum has been normalised (e.g. to sum to unity) + before convolution, the zero-loss intensity factor is removed and + the simple relationship above no longer holds. + .. _eels.fine_structure: Fine structure analysis using a spline function diff --git a/examples/model_fitting/EELS_curve_fitting.py b/examples/model_fitting/EELS_curve_fitting.py index ac3ca162e..a6f64e86d 100644 --- a/examples/model_fitting/EELS_curve_fitting.py +++ b/examples/model_fitting/EELS_curve_fitting.py @@ -30,6 +30,16 @@ m.enable_fine_structure() m.multifit(kind="smart") +# %% +# Print the fitted intensities and convert them to areal density +# (atoms/nm²). This assumes the low-loss spectrum was passed *raw* +# (counts) and that both spectra share the same dwell time. + +m.quantify() +for edge in m.edges: + N = edge.intensity.value * 1e10 + print(f"{edge.name}: {N:.3e} atoms/nm²") + # %% # Plot the model fit result m.plot() diff --git a/exspy/_docstrings/model.py b/exspy/_docstrings/model.py index 0b2354180..81e990ee9 100644 --- a/exspy/_docstrings/model.py +++ b/exspy/_docstrings/model.py @@ -34,7 +34,13 @@ If an EELSSpectrum is provided, it will be assumed that it is a low-loss EELS spectrum, and it will be used to simulate the effect of multiple scattering by convolving it with the EELS - spectrum. + spectrum. For absolute quantification (areal density in + atoms/nm²), the low-loss spectrum should be the *raw* (i.e. + unnormalised) spectrum in counts. If the core-loss and low-loss + spectra were acquired with different dwell times, the fitted + ``intensity`` must be corrected by the ratio of dwell times + before conversion to areal density. See + :attr:`~exspy.components.EELSCLEdge.intensity` for details. auto_background : bool If True, and if spectrum is an EELS instance adds automatically a powerlaw to the model and estimate the parameters by the diff --git a/exspy/_misc/eels/base_gos.py b/exspy/_misc/eels/base_gos.py index a6bbfd874..c7ddbfa04 100644 --- a/exspy/_misc/eels/base_gos.py +++ b/exspy/_misc/eels/base_gos.py @@ -131,6 +131,23 @@ def as_dictionary(self, fullcopy=True): return dic def integrateq(self, onset_energy, angle, E0): + """Integrate the GOS over q to obtain the energy-differential cross-section. + + Parameters + ---------- + onset_energy : float + The edge onset energy in eV. + angle : float + The effective collection semi-angle in rad. + E0 : float + The electron beam energy in keV. + + Returns + ------- + scipy.interpolate.BSpline + A spline representing the energy-differential cross-section + in barns/eV/atom as a function of energy loss in eV. + """ a0 = scipy.constants.value("Bohr radius") R = scipy.constants.value("Rydberg constant times hc in eV") diff --git a/exspy/_misc/eels/hydrogenic_gos.py b/exspy/_misc/eels/hydrogenic_gos.py index 701a391b7..a1a17551e 100644 --- a/exspy/_misc/eels/hydrogenic_gos.py +++ b/exspy/_misc/eels/hydrogenic_gos.py @@ -154,6 +154,23 @@ def __init__(self, element_subshell): _logger.info(info_str) def integrateq(self, onset_energy, angle, E0): + """Integrate the GOS over q to obtain the energy-differential cross-section. + + Parameters + ---------- + onset_energy : float + The edge onset energy in eV. + angle : float + The effective collection semi-angle in rad. + E0 : float + The electron beam energy in keV. + + Returns + ------- + scipy.interpolate.BSpline + A spline representing the energy-differential cross-section + in barns/eV/atom as a function of energy loss in eV. + """ energy_shift = onset_energy - self.onset_energy self.energy_shift = energy_shift gamma = 1 + E0 / 511.06 diff --git a/exspy/components/_eels_cl_edge.py b/exspy/components/_eels_cl_edge.py index e7631b214..a6e606959 100644 --- a/exspy/components/_eels_cl_edge.py +++ b/exspy/components/_eels_cl_edge.py @@ -101,10 +101,34 @@ class EELSCLEdge(Component): onset_energy : Parameter The edge onset position intensity : Parameter - The factor by which the cross section is multiplied, what in - favourable cases is proportional to the number of atoms of - the element. It is a component.Parameter instance. - It is fixed by default. + The factor by which the cross section is multiplied. The + internal cross-section is given in barns per eV per atom, + so ``intensity`` is proportional to the areal density of + atoms. When the model is convolved with a *raw* low-loss + spectrum (i.e. counts), the fitted ``intensity`` equals + :math:`N \times 10^{-10}`, where :math:`N` is the areal + density in atoms/nm². Consequently: + + .. code-block:: python + + N_atoms_per_nm2 = edge.intensity.value * 1e10 + + If the core-loss and low-loss spectra were acquired with + different dwell times (or total acquisition times), the + ratio of the dwell times must be accounted for: + + .. code-block:: python + + N_atoms_per_nm2 = edge.intensity.value * (t_core / t_low) * 1e10 + + where ``t_core`` and ``t_low`` are the dwell times of the + core-loss and low-loss spectra, respectively. Note that + if the low-loss spectrum is normalised (e.g. summed to + unity) before convolution, the zero-loss intensity factor + is removed and the simple relationship above no longer + holds. + + It is a component.Parameter instance. It is fixed by default. effective_angle : Parameter The effective collection semi-angle. It is automatically calculated by ``set_microscope_parameters``. It is a diff --git a/exspy/models/eelsmodel.py b/exspy/models/eelsmodel.py index e5797e8ad..8633754af 100644 --- a/exspy/models/eelsmodel.py +++ b/exspy/models/eelsmodel.py @@ -673,7 +673,21 @@ def _fit_edge(self, edgenumber, start_energy=None, **kwargs): def quantify(self): """Prints the value of the intensity of all the independent - active EELS core loss edges defined in the model + active EELS core loss edges defined in the model. + + The printed ``intensity`` is the fitted scaling factor of the + cross-section (in barns/eV/atom). When the model has been + convolved with a *raw* low-loss spectrum, this value is + proportional to the areal density :math:`N` (atoms/nm²) via + :math:`\\text{intensity} = N \\times 10^{-10}`. See the + :attr:`~.components.EELSCLEdge.intensity` attribute for the + full conversion formula, including the dwell-time correction + required when the core-loss and low-loss acquisition times + differ. + + See Also + -------- + :attr:`exspy.components.EELSCLEdge.intensity` """ elements = {}