Coverage for local_installation/dynasor/post_processing/electron_scattering_factors.py: 97%
58 statements
« prev ^ index » next coverage.py v7.3.2, created at 2025-04-16 06:13 +0000
1import json
2from re import search
3from importlib.resources import files
4from typing import List, Dict, Optional
5from warnings import warn
7from numpy import exp, abs, pi
8from pandas import DataFrame, concat
9from .weights import Weights
12class ElectronScatteringFactors(Weights):
13 r"""This class generates sample weights corresponding to electron scattering factors.
14 The scattering factors are parametrized as sums of exponentials
15 of the following form
17 .. math::
19 f(q) = \sum_{i=1}^k a_i \exp(-b_i s^2),
21 with :math:`s=\sin(\theta)/\lambda`, where and :math:`a_i` and :math:`b_i`
22 are fitting parameters. Note that there are two parametrizations for the elemental
23 scattering factors; one valid up to :math:`s = 2.0\,\mathrm{Å}^{-1}`
24 and the other valid up to :math:`s = 6.0\,\mathrm{Å}^{-1}`.
25 By default, the parametrisation valid up to :math:`s = 6.0\,\mathrm{Å}^{-1}`
26 is used. Ionic scattering factors only have a single parametrization.
28 The parametrizations are based on the following publications:
30 * Elemental scattering factors are based on Table 1 and Table 3
31 from L.-M. Peng, G. Ren, S. L. Dudarev and M. J. Whelan,
32 Acta Crystallographica Section A **52**, 257-276 (1996);
33 `doi: 10.1107/S0108767395014371 <https://doi.org/10.1107/S0108767395014371>`_.
34 * Ionic scattering factors are based on Table 1 from
35 Lian-Mao. Peng, Acta Crystallographica Section A **54**, 481-485 (1998);
36 `doi: 10.1107/S0108767398001901 <https://doi.org/10.1107/S0108767398001901>`_.
38 Parameters
39 ----------
40 atom_types
41 List of atomic species for which to retrieve scattering lengths.
42 """
44 def __init__(
45 self,
46 atom_types: List[str],
47 parametrisation: str = 'smax6'
48 ):
49 self._parametrisation = parametrisation
50 scattering_factors = self._read_scattering_factors(self._parametrisation)
51 # Select the relevant species
52 scattering_factors = scattering_factors[
53 scattering_factors.index.isin(atom_types)] # Atom species is index
54 self._scattering_factors = scattering_factors
56 # Check if any of the fetched form factors are missing,
57 # indicating that it is missing in the experimental database.
58 for s in atom_types:
59 row = scattering_factors[scattering_factors.index == s]
60 if row.empty:
61 raise ValueError('Missing tabulated values '
62 f'for requested species {s}.')
64 weights_coh = scattering_factors.to_dict(orient='index')
65 supports_currents = False
66 super().__init__(weights_coh, None, supports_currents=supports_currents)
68 def get_weight_coh(self, atom_type, q_norm):
69 """Get the coherent weight for a given atom type and q-vector norm."""
70 return self._compute_f(self._weights_coh[atom_type],
71 q_norm,
72 self._parametrisation,
73 self._get_charge(atom_type))
75 def _get_charge(self, atom_type: str):
76 """
77 Extracts the ionic charge from the `atom_type`.
78 """
79 match = search(r'(\d+)?([+-])', atom_type)
80 if match:
81 magnitude = int(match.group(1)) if match.group(1) else 1
82 charge = magnitude * (1 if match.group(2) == '+' else -1)
83 return charge
84 return None # If no charge is found
86 def _compute_f(
87 self,
88 coefficients: Dict,
89 q_norm: float,
90 parametrisation: str,
91 charge: Optional[int] = None):
92 r"""Compute electronic scattering factors f(q).
93 The scattering factors are parametrized as a sum
94 of exponentials of the form
96 .. math::
98 f(q) = \sum_{i=1}^k a_i * exp(-b_i * s**2)
100 Ions are also offset by
102 .. math::
104 (m_0 e^2) / (8 \pi^2 \hbar^2) C / s^2 \approx 0.023934 C / s^2
106 where :math:`C` is the ionic charge.
108 Parameters
109 ----------
110 coefficients
111 Parametrization parameters, read from the corresponding source file.
112 q_norm
113 The |q|-value at which to evaluate the form factor.
114 parametrisation
115 If the Peng parametrisation valid up to :math:`s = 6.0\,\mathrm{Å}^{-1}`
116 or :math:`s = 2.0\,\mathrm{Å}^{-1}` should be used.
117 charge
118 Integer representing the ionic charge. None if no charge to avoid division by zero.
119 """
120 s_max = 6.0 if parametrisation == 'smax6' else 2.0
122 s = q_norm / (4 * pi) # q in dynasor is q = 4 pi sin(theta) / lambda.
123 s_squared = s*s # s = sin(theta) / lambda in the Waasmaier paper.
124 if abs(s) > s_max:
125 warn(f'Peng parametrization is not reliable for q'
126 f' above {(s_max * 4 * pi):.2f} rad/Å'
127 f' (corresponding to s={s_max} 1/Å).'
128 ' Parametrisations for ions are accurate up to s=6.0 1/Å')
129 nmax = 5
131 if charge is not None:
132 f = 0.023934 * charge / s_squared
133 else:
134 f = 0.0
135 for i in range(1, nmax+1):
136 f += coefficients[f'a{i}'] * exp(-coefficients[f'b{i}'] * s_squared)
137 return f
139 def _read_scattering_factors(self, parametrisation: str) -> DataFrame:
140 r"""
141 Extracts the parametrization for the form factors :math:`f(q)`,
142 for both elemental systems and ionic species.
144 Parameters
145 ----------
146 parametrisation
147 If the Peng parametrisation valid up to :math:`s = 6.0\,\mathrm{Å}^{-1}`
148 or :math:`s = 2.0\,\mathrm{Å}^{-1}` should be used.
149 """
150 if parametrisation == 'smax2':
151 data_file = files(__package__) / \
152 'form-factors/electron-parameters-kmax2-peng-1996.json'
153 elif parametrisation == 'smax6': 153 ↛ 157line 153 didn't jump to line 157, because the condition on line 153 was never false
154 data_file = files(__package__) / \
155 'form-factors/electron-parameters-kmax6-peng-1996.json'
156 else:
157 raise ValueError(f'Unknown parametrisation {parametrisation}')
159 data_file_ions = files(__package__) / \
160 'form-factors/electron-parameters-ions-peng-1998.json'
162 frames = []
163 for df in [data_file, data_file_ions]:
164 with open(df) as fp:
165 coefficients = json.load(fp)
166 frames.append(DataFrame.from_dict(coefficients))
168 scattering_factors = concat(frames)
169 scattering_factors.index.names = ['species']
170 return scattering_factors