API reference#
This reference manual details the public modules, classes, and functions in diffsims, as generated from their docstrings. Many of the docstrings contain examples. See the user guide for other ways to use diffsims.
Caution
diffsims is in an alpha stage, so there may be breaking changes with each release.
The list of top modules:
Generation of reciprocal lattice vectors (crystal plane, reflector, g, hkl) for a crystal structure. |
|
Generation of diffraction simulations and libraries, and lists of rotations. |
|
Diffraction, structure and vector libraries. |
|
Diffraction simulations. |
|
Calculation of scattering factors and structure factors. |
|
Diffraction utilities used by the other modules. |
crystallography#
Generation of reciprocal lattice vectors (crystal plane, reflector, g, hkl) for a crystal structure.
|
Reciprocal lattice vectors \((hkl)\) for use in electron diffraction analysis and simulation. |
|
[Deprecated] Reciprocal lattice point (or crystal plane, reflector, g, etc.) with Miller indices, length of the reciprocal lattice vectors and other relevant structure_factor parameters. |
|
Return symmetrically equivalent Miller indices. |
|
Return the highest Miller indices hkl of the plane with a direct space interplanar spacing greater than but closest to a lower threshold. |
|
Return a list of planes from a set of highest Miller indices. |
ReciprocalLatticeVector#
Methods
|
Calculate angles between reciprocal lattice vectors, possibly using symmetrically equivalent vectors to find the smallest angle under symmetry. |
|
Populate |
|
Populate |
|
Cross product between reciprocal lattice vectors producing zone axes \([uvw]\) or \([UVTW]\) in the direct lattice. |
|
Get a deepcopy of the vectors. |
|
Draw great or small circles with a given |
|
Dot product of the vectors with other reciprocal lattice vectors. |
|
Outer dot product of the vectors with other reciprocal lattice vectors. |
|
A new instance with these reciprocal lattice vectors in a single column. |
|
Create a set of unique reciprocal lattice vectors from three highest indices. |
|
Create a set of unique reciprocal lattice vectors with a a direct space interplanar spacing greater than a lower threshold. |
|
Create a new instance from a |
|
Get vectors delineating great or small circle(s) with a given |
Get unique sets of \({hkl}\) for the vectors and the indices of vectors in each set. |
|
Table with indices, structure factor values and multiplicity of each set of \({hkl}\). |
|
|
A new instance with these reciprocal lattice vectors reshaped. |
Sanitise the |
|
|
Plot vectors in the stereographic projection. |
|
A new instance with these reciprocal lattice vectors where singleton dimensions are removed. |
|
A new instance from a sequence of reciprocal lattice vectors. |
|
Unique vectors symmetrically equivalent to the vectors. |
Return the vectors as a |
|
|
Return the azimuth \(\phi\), polar \(\theta\), and radial \(r\) spherical coordinates defined as in the ISO 31-11 standard :cite:`weisstein2005spherical`. |
|
A new instance with the navigation shape of these reciprocal lattice vectors transposed. |
|
The unique vectors. |
- class diffsims.crystallography.ReciprocalLatticeVector(phase, xyz=None, hkl=None, hkil=None)[source]#
Bases:
Vector3d
Reciprocal lattice vectors \((hkl)\) for use in electron diffraction analysis and simulation.
All lengths are assumed to be given in Å or inverse Å.
This class extends
orix.vector.Vector3d
to reciprocal lattice vectors \((hkl)\) specifically for diffraction experiments and simulations. It is thus different fromorix.vector.Miller
, which is a general class for Miller indices both in reciprocal and direct space. It supports relevant methods also supported in Miller, like obtaining a set of vectors from a minimal interplanar spacing.Create a set of reciprocal lattice vectors from \((hkl)\) or \((hkil)\).
The vectors are stored internally as cartesian coordinates in
data
.- Parameters:
phase (orix.crystal_map.Phase) – A phase with a crystal lattice and symmetry.
xyz (numpy.ndarray, list, or tuple, optional) – Cartesian coordinates of indices of reciprocal lattice vector(s)
hkl
. Default isNone
. This,hkl
, orhkil
is required.hkl (numpy.ndarray, list, or tuple, optional) – Indices of reciprocal lattice vector(s). Default is
None
. This,xyz
, orhkil
is required.hkil (numpy.ndarray, list, or tuple, optional) – Indices of reciprocal lattice vector(s), often preferred over
hkl
in trigonal and hexagonal lattices. Default isNone
. This,xyz
, orhkl
is required.
Examples
>>> from diffpy.structure import Atom, Lattice, Structure >>> from orix.crystal_map import Phase >>> from diffsims.crystallography import ReciprocalLatticeVector >>> phase = Phase( ... "al", ... space_group=225, ... structure=Structure( ... lattice=Lattice(4.04, 4.04, 4.04, 90, 90, 90), ... atoms=[Atom("Al", [0, 0, 1])], ... ), ... ) >>> rlv = ReciprocalLatticeVector(phase, hkl=[[1, 1, 1], [2, 0, 0]]) >>> rlv ReciprocalLatticeVector (2,), al (m-3m) [[1. 1. 1.] [2. 0. 0.]]
- property allowed#
Return whether vectors diffract according to diffraction selection rules assuming kinematic scattering theory.
Integer vectors are assumed.
- Returns:
allowed – Boolean array.
- Return type:
Examples
>>> from diffpy.structure import Atom, Lattice, Structure >>> from orix.crystal_map import Phase >>> from diffsims.crystallography import ReciprocalLatticeVector >>> phase = Phase( ... "al", ... space_group=225, ... structure=Structure( ... lattice=Lattice(4.04, 4.04, 4.04, 90, 90, 90), ... atoms=[Atom("Al", [0, 0, 1])], ... ), ... ) >>> rlv = ReciprocalLatticeVector( ... phase, hkl=[[1, 0, 0], [2, 0, 0]] ... ) >>> rlv.allowed array([False, True])
- angle_with(other, use_symmetry=False)[source]#
Calculate angles between reciprocal lattice vectors, possibly using symmetrically equivalent vectors to find the smallest angle under symmetry.
- Parameters:
other (ReciprocalLatticeVector) – Other vectors of compatible shape to the vectors.
use_symmetry (bool, optional) – Whether to consider equivalent vectors to find the smallest angle under symmetry. Default is
False
.
- Returns:
The angle between the vectors, in radians. If
use_symmetry=True
, the angles are the smallest under symmetry.- Return type:
- property azimuth: ndarray#
Azimuth spherical coordinate, i.e. the angle \(\phi \in [0, 2\pi]\) from the positive z-axis to a point on the sphere, according to the ISO 31-11 standard :cite:`weisstein2005spherical`.
- Return type:
azimuth
- calculate_structure_factor(scattering_params='xtables')[source]#
Populate
structure_factor
with the complex kinematical structure factor \(F_{hkl}\) for each vector.- Parameters:
scattering_params (str) – Which atomic scattering factors to use, either
"xtables"
(default) or"lobato"
.
Examples
See
ReciprocalLatticeVector
for the creation ofrlv
>>> rlv ReciprocalLatticeVector (2,), al (m-3m) [[1. 1. 1.] [2. 0. 0.]]
A unit cell with all asymmetric atom positions is required to calculate structure factors
>>> rlv.phase.structure [Al 0.000000 0.000000 1.000000 1.0000] >>> rlv.sanitise_phase() >>> rlv.phase.structure [Al 0.000000 0.000000 0.000000 1.0000, Al 0.000000 0.500000 0.500000 1.0000, Al 0.500000 0.000000 0.500000 1.0000, Al 0.500000 0.500000 0.000000 1.0000]
>>> rlv.calculate_structure_factor() >>> rlv.structure_factor array([8.46881663-1.55569638e-15j, 7.04777513-8.63103525e-16j])
Default atomic scattering factors are from the International Tables of Crystallography Vol. C Table 4.3.2.3. Alternative scattering factors are available from Lobato and Van Dyck Acta Cryst. (2014). A70, 636-649 https://doi.org/10.1107/S205327331401643X
>>> rlv.calculate_structure_factor("lobato") >>> rlv.structure_factor array([8.44934816-1.55212008e-15j, 7.0387957 -8.62003862e-16j])
- calculate_theta(voltage)[source]#
Populate
theta
with the Bragg angle \(theta_B\) in radians.Assumes
phase.structure
lattice parameters and Debye-Waller factors are expressed in Ångströms.- Parameters:
voltage (float) – Beam energy in V.
Examples
See
ReciprocalLatticeVector
for the creation ofrlv
>>> rlv ReciprocalLatticeVector (2,), al (m-3m) [[1. 1. 1.] [2. 0. 0.]]
>>> rlv.calculate_theta(20e3) >>> rlv.theta array([0.0184036 , 0.02125105]) >>> rlv.calculate_theta(200e3) >>> rlv.theta array([0.00537583, 0.00620749])
- property coordinate_format#
Vector coordinate format, either
"hkl"
or"hkil"
.- Return type:
Examples
See
ReciprocalLatticeVector
for the creation ofrlv
>>> rlv ReciprocalLatticeVector (2,), al (m-3m) [[1. 1. 1.] [2. 0. 0.]] >>> rlv.coordinate_format 'hkl' >>> rlv.coordinate_format = "hkil" >>> rlv ReciprocalLatticeVector (2,), al (m-3m) [[ 1. 1. -2. 1.] [ 2. 0. -2. 0.]]
- property coordinates#
Miller or Miller-Bravais indices.
- Returns:
coordinates – Miller indices if
coordiante_format
is"hkl"
or Miller-Bravais indices if it is"hkil"
.- Return type:
Examples
See
ReciprocalLatticeVector
for the creation ofrlv
>>> rlv ReciprocalLatticeVector (2,), al (m-3m) [[1. 1. 1.] [2. 0. 0.]] >>> rlv.coordinates array([[1., 1., 1.], [2., 0., 0.]]) >>> rlv.coordinate_format = "hkil" >>> rlv.coordinates array([[ 1., 1., -2., 1.], [ 2., 0., -2., 0.]])
- cross(other)[source]#
Cross product between reciprocal lattice vectors producing zone axes \([uvw]\) or \([UVTW]\) in the direct lattice.
- Parameters:
other (ReciprocalLatticeVector) – Other vectors of compatible shape to the vectors.
- Returns:
Direct lattice vector(s) \([uvw]\) or \(UVTW\), depending on whether the vector’s
coordinate_format
ishkl
orhkil
, respectively.- Return type:
- dim = 3#
Return the number of dimensions for this object.
- dot(other)[source]#
Dot product of the vectors with other reciprocal lattice vectors.
- Parameters:
other (ReciprocalLatticeVector) – Other vectors of compatible shape to the vectors.
- Return type:
- dot_outer(other)[source]#
Outer dot product of the vectors with other reciprocal lattice vectors.
The dot product for every combination of the vectors are computed.
- Parameters:
other (ReciprocalLatticeVector) – Other vectors of compatible shape to the vectors.
- Return type:
- draw_circle(projection: str = 'stereographic', figure: Figure | None = None, opening_angle: float | ndarray = 1.5707963267948966, steps: int = 100, reproject: bool = False, axes_labels: List[str] | None = None, hemisphere: str | None = None, show_hemisphere_label: bool | None = None, grid: bool | None = None, grid_resolution: Tuple[float, float] | None = None, figure_kwargs: Dict | None = None, reproject_plot_kwargs: Dict | None = None, return_figure: bool = False, **kwargs: Any) Figure | None #
Draw great or small circles with a given
opening_angle
to to the vectors in the stereographic projection.A vector must be present in the current hemisphere for its circle to be drawn.
- Parameters:
projection – Which projection to use. The default is
"stereographic"
, the only current option.figure – Which figure to plot onto. Default is
None
, which creates a new figure.opening_angle – Opening angle(s) around the vector(s). Default is \(\pi/2\), giving a great circle. If an array is passed, its size must be equal to the number of vectors.
steps – Number of vectors to describe each circle, default is 100.
reproject – Whether to reproject parts of the circle(s) visible on the other hemisphere. Reprojection is achieved by reflection of the circle(s) parts located on the other hemisphere in the projection plane. Ignored if hemisphere is
"both"
. Default isFalse
.axes_labels – Reference frame axes labels, defaults to
[None, None, None]
.hemisphere – Which hemisphere(s) to plot the vectors in, defaults to
None
, which means"upper"
if a new figure is created, otherwise adds to the current figure’s hemispheres. Options are"upper"
,"lower"
, and"both"
, which plots two projections side by side.show_hemisphere_label – Whether to show hemisphere labels
"upper"
or"lower"
. Default isTrue
ifhemisphere
is"both"
, otherwiseFalse
.grid – Whether to show the azimuth and polar grid. Default is whatever
axes.grid
is set to inmatplotlib.rcParams
.grid_resolution – Azimuth and polar grid resolution in degrees, as a tuple. Default is whatever is default in
stereographic_grid
.figure_kwargs – Dictionary of keyword arguments passed to
matplotlib.pyplot.subplots()
.reproject_plot_kwargs – Keyword arguments passed to
matplotlib.axes.Axes.plot()
to alter the appearance of parts of the circle(s) visible on the other hemisphere ifreproject=True
. These lines are dashed by default. Values used for circle(s) on the current hemisphere are used unless values are passed here.return_figure – Whether to return the figure (default is
False
).**kwargs – Keyword arguments passed to
matplotlib.axes.Axes.plot()
to alter the circles’ appearance.
- Returns:
The created figure, returned if
return_figure=True
.- Return type:
fig
Notes
This is a somewhat customizable convenience method which creates a figure with axes using
StereographicPlot
, however, it is meant for quick plotting and prototyping. This figure and the axes can also be created using Matplotlib directly, which is more customizable.
- property dspacing#
Direct lattice interplanar spacing \(d = 1 / g\).
- Return type:
Examples
See
ReciprocalLatticeVector
for the creation ofrlv
>>> rlv ReciprocalLatticeVector (2,), al (m-3m) [[1. 1. 1.] [2. 0. 0.]]
Lattice parameters are given in \(Å\)
>>> rlv.phase.structure.lattice Lattice(a=4.04, b=4.04, c=4.04, alpha=90, beta=90, gamma=90)
so \(d\) is given in \(Å\)
>>> rlv.dspacing array([2.33249509, 2.02 ])
- flatten()[source]#
A new instance with these reciprocal lattice vectors in a single column.
- Return type:
- classmethod from_highest_hkl(phase, hkl)[source]#
Create a set of unique reciprocal lattice vectors from three highest indices.
- Parameters:
phase (orix.crystal_map.Phase) – A phase with a crystal lattice and symmetry.
hkl (numpy.ndarray, list, or tuple) – Three highest reciprocal lattice vector indices.
Examples
>>> from diffpy.structure import Atom, Lattice, Structure >>> from orix.crystal_map import Phase >>> from diffsims.crystallography import ReciprocalLatticeVector >>> phase = Phase( ... "al", ... space_group=225, ... structure=Structure( ... lattice=Lattice(4.04, 4.04, 4.04, 90, 90, 90), ... atoms=[Atom("Al", [0, 0, 1])], ... ), ... ) >>> rlv = ReciprocalLatticeVector.from_highest_hkl(phase, [3, 3, 3]) >>> rlv ReciprocalLatticeVector (342,), al (m-3m) [[ 3. 3. 3.] [ 3. 3. 2.] [ 3. 3. 1.] ... [-3. -3. -1.] [-3. -3. -2.] [-3. -3. -3.]]
Vectors are included regardless of whether they are kinematically allowed or not
>>> rlv.allowed.all() False
- classmethod from_miller(miller)[source]#
Create a new instance from a
Miller
instance.- Parameters:
miller (orix.vector.Miller) – Reciprocal lattice vectors \((hk(i)l)\).
- Return type:
Examples
>>> from diffpy.structure import Atom, Lattice, Structure >>> from orix.crystal_map import Phase >>> from orix.vector import Miller >>> from diffsims.crystallography import ReciprocalLatticeVector >>> phase = Phase( ... "al", ... space_group=225, ... structure=Structure( ... lattice=Lattice(4.04, 4.04, 4.04, 90, 90, 90), ... atoms=[Atom("Al", [0, 0, 1])], ... ), ... ) >>> miller = Miller(hkl=[[1, 1, 1], [2, 0, 0]], phase=phase) >>> miller Miller (2,), point group m-3m, hkl [[1. 1. 1.] [2. 0. 0.]] >>> rlv = ReciprocalLatticeVector.from_miller(miller) >>> rlv ReciprocalLatticeVector (2,), al (m-3m) [[1. 1. 1.] [2. 0. 0.]] >>> rlv.to_miller() Miller (2,), point group m-3m, hkl [[1. 1. 1.] [2. 0. 0.]]
- classmethod from_min_dspacing(phase, min_dspacing=0.7)[source]#
Create a set of unique reciprocal lattice vectors with a a direct space interplanar spacing greater than a lower threshold.
- Parameters:
phase (orix.crystal_map.Phase) – A phase with a crystal lattice and symmetry.
min_dspacing (float, optional) – Smallest interplanar spacing to consider. Default is 0.7, in the unit used to define the lattice parameters in
phase
, which is assumed to be Ångström.
Examples
>>> from diffpy.structure import Atom, Lattice, Structure >>> from orix.crystal_map import Phase >>> from diffsims.crystallography import ReciprocalLatticeVector >>> phase = Phase( ... "al", ... space_group=225, ... structure=Structure( ... lattice=Lattice(4.04, 4.04, 4.04, 90, 90, 90), ... atoms=[Atom("Al", [0, 0, 1])], ... ), ... ) >>> rlv = ReciprocalLatticeVector.from_min_dspacing(phase) >>> rlv ReciprocalLatticeVector (798,), al (m-3m) [[ 5. 2. 2.] [ 5. 2. 1.] [ 5. 2. 0.] ... [-5. -2. 0.] [-5. -2. -1.] [-5. -2. -2.]] >>> rlv.dspacing.min() 0.7032737300610338
Vectors are included regardless of whether they are kinematically allowed or not
>>> rlv = ReciprocalLatticeVector.from_min_dspacing(phase, 1) >>> rlv.size 256 >>> rlv.allowed.all() False
- get_circle(opening_angle: float | ndarray = 1.5707963267948966, steps: int = 100) Vector3d #
Get vectors delineating great or small circle(s) with a given
opening_angle
about each vector.Used for plotting plane traces in stereographic projections.
- Parameters:
opening_angle – Opening angle(s) around the vector(s). Default is \(\pi/2\), giving a great circle. If an array is passed, its size must be equal to the number of vectors.
steps – Number of vectors to describe each circle, default is 100.
- Returns:
Vectors delineating circles with the
opening_angle
about the vectors.- Return type:
circles
Notes
A set of
steps
number of vectors equal to each vector is rotated twice to obtain a circle: (1) About a perpendicular vector to the current vector atopening_angle
and (2) about the current vector in a full circle.
- get_hkl_sets()[source]#
Get unique sets of \({hkl}\) for the vectors and the indices of vectors in each set.
- Returns:
hkl_sets – Dictionary with (h, k, l) as keys and a tuple with
numpy.ndarray
with integers of the vectors (possibly multi-dimensional) in each set. The keys (h, k, l) are rounded to six decimals so that applying integer values (h, k, l) as dictionary keys work.- Return type:
defaultdict
Examples
See
ReciprocalLatticeVector
for the creation ofrlv
>>> rlv ReciprocalLatticeVector (2,), al (m-3m) [[1. 1. 1.] [2. 0. 0.]] >>> hkl_sets = rlv.get_hkl_sets() >>> hkl_sets defaultdict(<class 'tuple'>, {(2.0, 0.0, 0.0): (array([1]),), (1.0, 1.0, 1.0): (array([0]),)}) >>> hkl_sets[2, 0, 0] (array([1]),) >>> rlv[hkl_sets[2, 0, 0]] ReciprocalLatticeVector (1,), al (m-3m) [[2. 0. 0.]]
- property gspacing#
Reciprocal lattice vector spacing \(g\).
- Return type:
Examples
See
ReciprocalLatticeVector
for the creation ofrlv
>>> rlv ReciprocalLatticeVector (2,), al (m-3m) [[1. 1. 1.] [2. 0. 0.]]
Lattice parameters are given in \(Å\)
>>> rlv.phase.structure.lattice Lattice(a=4.04, b=4.04, c=4.04, alpha=90, beta=90, gamma=90)
so \(g\) is given in \(Å^-1\)
>>> rlv.gspacing array([0.42872545, 0.4950495 ])
- property h#
First reciprocal lattice vector index.
- Return type:
Examples
See
ReciprocalLatticeVector
for the creation ofrlv
>>> rlv ReciprocalLatticeVector (2,), al (m-3m) [[1. 1. 1.] [2. 0. 0.]] >>> rlv.h array([1., 2.])
- property has_hexagonal_lattice#
Whether the crystal lattice is hexagonal/trigonal.
- Return type:
Examples
See
ReciprocalLatticeVector
for the creation ofrlv
>>> rlv ReciprocalLatticeVector (2,), al (m-3m) [[1. 1. 1.] [2. 0. 0.]] >>> rlv.has_hexagonal_lattice False
- property hkil#
Miller-Bravais indices.
- Return type:
Examples
See
ReciprocalLatticeVector
for the creation ofrlv
>>> rlv ReciprocalLatticeVector (2,), al (m-3m) [[1. 1. 1.] [2. 0. 0.]] >>> rlv.hkil array([[ 1., 1., -2., 1.], [ 2., 0., -2., 0.]])
- property hkl#
Miller indices.
- Return type:
Examples
See
ReciprocalLatticeVector
for the creation ofrlv
>>> rlv ReciprocalLatticeVector (2,), al (m-3m) [[1. 1. 1.] [2. 0. 0.]] >>> rlv.hkl array([[1., 1., 1.], [2., 0., 0.]])
- property i#
Third reciprocal lattice vector index in 4-index Miller-Bravais indices, equal to \(-(h + k)\).
- Return type:
Examples
See
ReciprocalLatticeVector
for the creation ofrlv
>>> rlv ReciprocalLatticeVector (2,), al (m-3m) [[1. 1. 1.] [2. 0. 0.]] >>> rlv.i array([-2., -2.])
- property k#
Second reciprocal lattice vector index.
- Return type:
Examples
See
ReciprocalLatticeVector
for the creation ofrlv
>>> rlv ReciprocalLatticeVector (2,), al (m-3m) [[1. 1. 1.] [2. 0. 0.]] >>> rlv.k array([1., 0.])
- property l#
Third reciprocal lattice vector index, or fourth index in 4-index Miller Bravais indices.
- Return type:
Examples
See
ReciprocalLatticeVector
for the creation ofrlv
>>> rlv ReciprocalLatticeVector (2,), al (m-3m) [[1. 1. 1.] [2. 0. 0.]] >>> rlv.l array([1., 0.])
- property multiplicity#
Number of symmetrically equivalent directions per vector.
- Returns:
mult
- Return type:
Examples
See
ReciprocalLatticeVector
for the creation ofrlv
>>> rlv ReciprocalLatticeVector (2,), al (m-3m) [[1. 1. 1.] [2. 0. 0.]] >>> rlv.multiplicity array([8, 6])
- property polar: ndarray#
Polar spherical coordinate, i.e. the angle \(\theta \in [0, \pi]\) from the positive z-axis to a point on the sphere, according to the ISO 31-11 standard :cite:`weisstein2005spherical`.
- Return type:
polar
- print_table()[source]#
Table with indices, structure factor values and multiplicity of each set of \({hkl}\).
Examples
See
ReciprocalLatticeVector
for the creation ofrlv
>>> rlv ReciprocalLatticeVector (2,), al (m-3m) [[1. 1. 1.] [2. 0. 0.]] >>> rlv.print_table() h k l d |F|_hkl |F|^2 |F|^2_rel Mult 1 1 1 2.332 nan nan nan 8 2 0 0 2.020 nan nan nan 6 >>> rlv.sanitise_phase() >>> rlv.calculate_structure_factor() >>> rlv.print_table() h k l d |F|_hkl |F|^2 |F|^2_rel Mult 1 1 1 2.332 8.5 71.7 100.0 8 2 0 0 2.020 7.0 49.7 69.3 6
- property radial: ndarray#
Return the radial spherical coordinate, i.e. the distance from a point on the sphere to the origin, according to the ISO 31-11 standard :cite:`weisstein2005spherical`.
- Returns:
Radial spherical coordinate.
- Return type:
radial
- reshape(*shape)[source]#
A new instance with these reciprocal lattice vectors reshaped.
- Parameters:
*shape (int) – Multiple integers designating the new shape.
- Return type:
- sanitise_phase()[source]#
Sanitise the
phase
inplace for calculation of structure factors.The phase is sanitised when it’s
structure
has an expanded unit cell with all symmetrically atom positions filled, and the atoms have theirelement
set to a string, e.g. “Al”.Examples
See
ReciprocalLatticeVector
for the creation ofrlv
>>> rlv ReciprocalLatticeVector (2,), al (m-3m) [[1. 1. 1.] [2. 0. 0.]] >>> rlv.phase.structure [Al 0.000000 0.000000 1.000000 1.0000] >>> rlv.sanitise_phase() >>> rlv.phase.structure [Al 0.000000 0.000000 0.000000 1.0000, Al 0.000000 0.500000 0.500000 1.0000, Al 0.500000 0.000000 0.500000 1.0000, Al 0.500000 0.500000 0.000000 1.0000]
- scatter(projection: str = 'stereographic', figure: Figure | None = None, axes_labels: List[str] | None = None, vector_labels: List[str] | None = None, hemisphere: str | None = None, reproject: bool = False, show_hemisphere_label: bool | None = None, grid: bool | None = None, grid_resolution: Tuple[float, float] | None = None, figure_kwargs: Dict | None = None, reproject_scatter_kwargs: Dict | None = None, text_kwargs: Dict | None = None, return_figure: bool = False, **kwargs: Any) Figure | None #
Plot vectors in the stereographic projection.
- Parameters:
projection – Which projection to use. The default is
"stereographic"
, the only current option.figure – Which figure to plot onto. Default is
None
, which creates a new figure.axes_labels – Reference frame axes labels, defaults to
[None, None, None]
.vector_labels – Vector text labels, which by default are not added.
hemisphere – Which hemisphere(s) to plot the vectors in, defaults to
None
, which means"upper"
if a new figure is created, otherwise adds to the current figure’s hemispheres. Options are"upper"
,"lower"
, and"both"
, which plots two projections side by side.reproject – Whether to reproject vectors onto the chosen hemisphere. Reprojection is achieved by reflection of the vectors located on the opposite hemisphere in the projection plane. Ignored if
hemisphere
is"both"
. Default isFalse
.show_hemisphere_label – Whether to show hemisphere labels
"upper"
or"lower"
. Default isTrue
ifhemisphere
is"both"
, otherwiseFalse
.grid – Whether to show the azimuth and polar grid. Default is whatever
axes.grid
is set to inmatplotlib.rcParams
.grid_resolution – Azimuth and polar grid resolution in degrees, as a tuple. Default is whatever is default in
stereographic_grid
.figure_kwargs – Dictionary of keyword arguments passed to
matplotlib.pyplot.subplots()
.reproject_scatter_kwargs – Dictionary of keyword arguments for the reprojected scatter points which is passed to
scatter()
, which passes these on tomatplotlib.axes.Axes.scatter()
. The default marker style for reprojected vectors is"+"
. Values used for vector(s) on the visible hemisphere are used unless another value is passed here.text_kwargs – Dictionary of keyword arguments passed to
text()
, which passes these on tomatplotlib.axes.Axes.text()
.return_figure – Whether to return the figure (default is
False
).**kwargs – Keyword arguments passed to
scatter()
, which passes these on tomatplotlib.axes.Axes.scatter()
.
- Returns:
The created figure, returned if
return_figure=True
.- Return type:
fig
Notes
This is a somewhat customizable convenience method which creates a figure with axes using
StereographicPlot
, however, it is meant for quick plotting and prototyping. This figure and the axes can also be created using Matplotlib directly, which is more customizable.See also
- property scattering_parameter#
Scattering parameter \(0.5 \cdot g\).
- Return type:
Examples
See
ReciprocalLatticeVector
for the creation ofrlv
>>> rlv ReciprocalLatticeVector (2,), al (m-3m) [[1. 1. 1.] [2. 0. 0.]]
Lattice parameters are given in \(Å\)
>>> rlv.phase.structure.lattice Lattice(a=4.04, b=4.04, c=4.04, alpha=90, beta=90, gamma=90)
so the scattering parameters are given in \(Å^-1\)
>>> rlv.scattering_parameter array([0.21436272, 0.24752475])
- squeeze()[source]#
A new instance with these reciprocal lattice vectors where singleton dimensions are removed.
- Return type:
- classmethod stack(sequence)[source]#
A new instance from a sequence of reciprocal lattice vectors.
- Parameters:
sequence (iterable of ReciprocalLatticeVector) – One or more sets of compatible reciprocal lattice vectors.
- Return type:
- property structure_factor#
Kinematical structure factors \(F\).
- Returns:
structure_factor – Complex array. Filled with
None
ifcalculate_structure_factor()
hasn’t been called yet.- Return type:
Examples
See
ReciprocalLatticeVector
for the creation ofrlv
>>> rlv ReciprocalLatticeVector (2,), al (m-3m) [[1. 1. 1.] [2. 0. 0.]]
Kinematical structure factors are by default not calculated
>>> rlv.structure_factor array([nan+0.j, nan+0.j])
A unit cell with all asymmetric atom positions is required to calculate structure factors
>>> rlv.phase.structure [Al 0.000000 0.000000 1.000000 1.0000] >>> rlv.sanitise_phase() >>> rlv.phase.structure [Al 0.000000 0.000000 0.000000 1.0000, Al 0.000000 0.500000 0.500000 1.0000, Al 0.500000 0.000000 0.500000 1.0000, Al 0.500000 0.500000 0.000000 1.0000]
>>> rlv.calculate_structure_factor() >>> rlv.structure_factor array([8.46881663-1.55569638e-15j, 7.04777513-8.63103525e-16j])
- symmetrise(return_multiplicity=False, return_index=False)[source]#
Unique vectors symmetrically equivalent to the vectors.
- Parameters:
- Returns:
ReciprocalLatticeVector – Flattened symmetrically equivalent vectors.
multiplicity (numpy.ndarray) – Multiplicity of each vector. Returned if
return_multiplicity=True
.idx (numpy.ndarray) – Index into the vectors for the symmetrically equivalent vectors. Returned if
return_index=True
.
Examples
See
ReciprocalLatticeVector
for the creation ofrlv
>>> rlv ReciprocalLatticeVector (2,), al (m-3m) [[1. 1. 1.] [2. 0. 0.]] >>> rlv.symmetrise() ReciprocalLatticeVector (14,), al (m-3m) [[ 1. 1. 1.] [-1. 1. 1.] [-1. -1. 1.] [ 1. -1. 1.] [ 1. -1. -1.] [ 1. 1. -1.] [-1. 1. -1.] [-1. -1. -1.] [ 2. 0. 0.] [ 0. 2. 0.] [-2. 0. 0.] [ 0. -2. 0.] [ 0. 0. 2.] [ 0. 0. -2.]] >>> _, mult, idx = rlv.symmetrise( ... return_multiplicity=True, return_index=True ... ) >>> mult array([8, 6]) >>> idx array([0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1])
- property theta#
Twice the Bragg angle.
- Returns:
theta – Filled with
None
ifcalculate_theta()
hasn’t been called yet.- Return type:
Examples
See
ReciprocalLatticeVector
for the creation ofrlv
>>> rlv ReciprocalLatticeVector (2,), al (m-3m) [[1. 1. 1.] [2. 0. 0.]]
Bragg angles are by default not calculated
>>> rlv.theta array([nan, nan])
>>> rlv.calculate_theta(20e3) >>> rlv.theta array([0.0184036 , 0.02125105])
- to_miller()[source]#
Return the vectors as a
Miller
instance.- Return type:
Examples
See
ReciprocalLatticeVector
for the creation ofrlv
>>> rlv ReciprocalLatticeVector (2,), al (m-3m) [[1. 1. 1.] [2. 0. 0.]] >>> rlv.to_miller() Miller (2,), point group m-3m, hkl [[1. 1. 1.] [2. 0. 0.]]
- to_polar(degrees: bool = False) Tuple[ndarray, ndarray, ndarray] #
Return the azimuth \(\phi\), polar \(\theta\), and radial \(r\) spherical coordinates defined as in the ISO 31-11 standard :cite:`weisstein2005spherical`.
- Parameters:
degrees – If
True
, the given angles are returned in degrees. Default isFalse
.- Returns:
azimuth – Azimuth angles in radians (
degrees=False
) or degrees (degrees=True
).polar – Polar angles in radians (
degrees=False
) or degrees (degrees=True
).radial – Radial values.
- transpose(*axes)[source]#
A new instance with the navigation shape of these reciprocal lattice vectors transposed.
If
ndim
is originally 2, then order may be undefined. In this case the first two dimensions will be transposed.- Parameters:
*axes (int, optional) – Transposed axes order. Only navigation axes need to be defined. May be undefined if the vectors only contain two navigation dimensions.
- Return type:
- unique(use_symmetry=False, return_index=False)[source]#
The unique vectors.
- Parameters:
- Returns:
ReciprocalLatticeVector – Flattened unique vectors.
idx (numpy.ndarray) – Indices of the unique data in the (flattened) array.
- property unit#
Unit reciprocal lattice vectors.
- Return type:
- property x: ndarray#
Return or set the x coordinates.
- Parameters:
value (np.ndarray) – The new x coordinates.
- property xyz: Tuple[ndarray, ndarray, ndarray]#
Return the coordinates as three arrays, useful for plotting.
ReciprocalLatticePoint#
- class diffsims.crystallography.ReciprocalLatticePoint(phase, hkl)[source]#
Bases:
object
[Deprecated] Reciprocal lattice point (or crystal plane, reflector, g, etc.) with Miller indices, length of the reciprocal lattice vectors and other relevant structure_factor parameters.
Notes
Deprecated since version 0.5: Class
ReciprocalLatticePoint
is deprecated and will be removed in version 0.6. UseReciprocalLatticeVector
instead.- property allowed#
Return whether planes diffract according to structure_factor selection rules assuming kinematical scattering theory.
- calculate_structure_factor(method=None, voltage=None)[source]#
Populate self.structure_factor with the structure factor F for each plane.
- Parameters:
method (str, optional) – Either “kinematical” for kinematical X-ray structure factors or “doyleturner” for structure factors using Doyle-Turner atomic scattering factors. If None (default), kinematical structure factors are calculated.
voltage (float, optional) – Beam energy in V used when method=doyleturner.
- calculate_theta(voltage)[source]#
Populate self.theta with the Bragg angle \(theta_B\) for each plane.
- Parameters:
voltage (float) – Beam energy in V.
- property dspacing#
Return
np.ndarray
of direct lattice interplanar spacings.
- classmethod from_highest_hkl(phase, highest_hkl)[source]#
Create a CrystalPlane object populated by unique Miller indices below, but including, a set of higher indices.
- Parameters:
phase (orix.crystal_map.phase_list.Phase) – A phase container with a crystal structure and a space and point group describing the allowed symmetry operations.
highest_hkl (np.ndarray, list, or tuple of int) – Highest Miller indices to consider (including).
- classmethod from_min_dspacing(phase, min_dspacing=0.5)[source]#
Create a CrystalPlane object populated by unique Miller indices with a direct space interplanar spacing greater than a lower threshold.
- Parameters:
phase (orix.crystal_map.phase_list.Phase) – A phase container with a crystal structure and a space and point group describing the allowed symmetry operations.
min_dspacing (float, optional) – Smallest interplanar spacing to consider. Default is 0.5 Å.
- property gspacing#
Return
np.ndarray
of reciprocal lattice point spacings.
- property h#
Return
np.ndarray
of Miller index h.
- property k#
Return
np.ndarray
of Miller index k.
- property l#
Return
np.ndarray
of Miller index l.
- property multiplicity#
Return either int or
np.ndarray
of int.
- property scattering_parameter#
Return
np.ndarray
of scattering parameters s.
- property shape#
Return tuple.
- property size#
Return int.
- property structure_factor#
Return
np.ndarray
of structure factors F or None.
- symmetrise(antipodal=True, unique=True, return_multiplicity=False)[source]#
Return planes with symmetrically equivalent Miller indices.
- Parameters:
antipodal (bool, optional) – Whether to include antipodal symmetry operations. Default is True.
unique (bool, optional) – Whether to return only distinct indices. Default is True. If True, zero-entries, which are assumed to be degenerate, are removed.
return_multiplicity (bool, optional) – Whether to return the multiplicity of indices. This option is only available if unique is True. Default is False.
- Returns:
ReciprocalLatticePoint – Planes with Miller indices symmetrically equivalent to the original planes.
multiplicity (np.ndarray) – Multiplicity of the original Miller indices. Only returned if return_multiplicity is True.
Notes
Should be the same as EMsoft’s CalcFamily in their symmetry.f90 module, although not entirely sure. Use with care.
- property theta#
Return
np.ndarray
of twice the Bragg angle.
Functions#
- diffsims.crystallography.get_equivalent_hkl(hkl, operations, unique=False, return_multiplicity=False)[source]#
Return symmetrically equivalent Miller indices.
- Parameters:
hkl (orix.vector.Vector3d, np.ndarray, list or tuple of int) – Miller indices.
operations (orix.quaternion.symmetry.Symmetry) – Point group describing allowed symmetry operations.
unique (bool, optional) – Whether to return only unique Miller indices. Default is False.
return_multiplicity (bool, optional) – Whether to return the multiplicity of the input indices. Default is False.
- Returns:
new_hkl (orix.vector.Vector3d) – The symmetrically equivalent Miller indices.
multiplicity (np.ndarray) – Number of symmetrically equivalent indices. Only returned if return_multiplicity is True.
- diffsims.crystallography.get_highest_hkl(lattice, min_dspacing=0.5)[source]#
Return the highest Miller indices hkl of the plane with a direct space interplanar spacing greater than but closest to a lower threshold.
- Parameters:
lattice (diffpy.structure.Lattice) – Crystal lattice.
min_dspacing (float, optional) – Smallest interplanar spacing to consider. Default is 0.5 Å.
- Returns:
highest_hkl – Highest Miller indices.
- Return type:
np.ndarray
- diffsims.crystallography.get_hkl(highest_hkl)[source]#
Return a list of planes from a set of highest Miller indices.
- Parameters:
highest_hkl (orix.vector.Vector3d, np.ndarray, list, or tuple of int) – Highest Miller indices to consider.
- Returns:
hkl – An array of Miller indices.
- Return type:
np.ndarray
generators#
Generation of diffraction simulations and libraries, and lists of rotations.
diffraction_generator#
Electron diffraction pattern simulation.
- class diffsims.generators.diffraction_generator.AtomicDiffractionGenerator(accelerating_voltage, detector, reciprocal_mesh=False)[source]#
Bases:
object
Computes electron diffraction patterns for an atomic lattice.
- Parameters:
accelerating_voltage (float, 'inf') – The accelerating voltage of the microscope in kV
detector (list of 1D float-type arrays) – List of mesh vectors defining the (flat) detector size and sensor positions
reciprocal_mesh (bool, optional) – If True then detector is assumed to be a reciprocal grid, else (default) it is assumed to be a real grid.
- calculate_ed_data(structure, probe, slice_thickness, probe_centre=None, z_range=200, precessed=False, dtype='float64', ZERO=1e-14, mode='kinematic', **kwargs)[source]#
Calculates single electron diffraction image for particular atomic structure and probe.
- Parameters:
structure (Structure) – The structure for upon which to perform the calculation
probe (instance of probeFunction) – Function representing 3D shape of beam
slice_thickness (float) – Discretisation thickness in the z-axis
probe_centre (ndarray (or iterable), shape [3] or [2]) – Translation vector for the probe. Either of the same dimension of the space or the dimension of the detector. default=None focusses the probe at [0,0,0]
zrange (float) – z-thickness to discretise. Only required if sample is not thick enough to fully resolve the Ewald-sphere. Default value is 200.
precessed (bool, float, or (float, int)) – Dictates whether beam precession is simulated. If False or the float is 0 then no precession is computed. If <precessed> = (alpha, n) then the precession arc of tilt alpha (in degrees) is discretised into n projections. If n is not provided then default of 30 is used.
dtype (str or numpy.dtype) – Defines the precision to use whilst computing diffraction image.
ZERO (float > 0) – Rounding error permitted in computation of atomic density. This value is the smallest value rounded to 0. Default is 1e-14.
mode (str) – Only <mode>=’kinematic’ is currently supported.
kwargs (dictionary) – Extra key-word arguments to pass to child simulator. For kinematic: GPU (bool): Flag to use GPU if available, default is True. pointwise (bool): Flag to evaluate charge pointwise on voxels rather than average, default is False.
- Returns:
Diffraction data to be interpreted as a discretisation on the original detector mesh.
- Return type:
ndarray
- class diffsims.generators.diffraction_generator.DiffractionGenerator(accelerating_voltage, scattering_params='lobato', precession_angle=0, shape_factor_model='lorentzian', approximate_precession=True, minimum_intensity=1e-20, **kwargs)[source]#
Bases:
object
Computes electron diffraction patterns for a crystal structure.
Calculate reciprocal lattice of structure. Find all reciprocal points within the limiting sphere given by \(\\frac{2}{\\lambda}\).
For each reciprocal point \(\\mathbf{g_{hkl}}\) corresponding to lattice plane \((hkl)\), compute the Bragg condition \(\\sin(\\theta) = \\frac{\\lambda}{2d_{hkl}}\)
The intensity of each reflection is then given in the kinematic approximation as the modulus square of the structure factor. \(I_{hkl} = F_{hkl}F_{hkl}^*\)
- Parameters:
accelerating_voltage (float) – The accelerating voltage of the microscope in kV.
scattering_params (str) – “lobato”, “xtables” or None. If None is provided then atomic scattering is not taken into consideration.
precession_angle (float) – Angle about which the beam is precessed in degrees. Default is no precession.
shape_factor_model (func or str) – A function that takes excitation_error and max_excitation_error (and potentially kwargs) and returns an intensity scaling factor. If None defaults to shape_factor_models.linear. A number of pre-programmed functions are available via strings.
approximate_precession (bool) – When using precession, whether to precisely calculate average excitation errors and intensities or use an approximation.
minimum_intensity (float) – Minimum intensity for a peak to be considered visible in the pattern (fractional from the maximum).
kwargs – Keyword arguments passed to shape_factor_model.
Notes
When using precession and approximate_precession=True, the shape factor model defaults to Lorentzian; shape_factor_model is ignored. Only with approximate_precession=False the custom shape_factor_model is used.
- calculate_ed_data(structure, reciprocal_radius, rotation=(0, 0, 0), with_direct_beam=True, max_excitation_error=0.01, shape_factor_width=None, debye_waller_factors={})[source]#
Calculates the Electron Diffraction data for a structure.
- Parameters:
structure (diffpy.structure.structure.Structure) – The structure for which to derive the diffraction pattern. Note that the structure must be rotated to the appropriate orientation and that testing is conducted on unit cells (rather than supercells).
reciprocal_radius (float) – The maximum radius of the sphere of reciprocal space to sample, in reciprocal Angstroms.
rotation (tuple) – Euler angles, in degrees, in the rzxz convention. Default is (0, 0, 0) which aligns ‘z’ with the electron beam.
with_direct_beam (bool) – If True, the direct beam is included in the simulated diffraction pattern. If False, it is not.
max_excitation_error (float) – The cut-off for geometric excitation error in the z-direction in units of reciprocal Angstroms. Spots with a larger distance from the Ewald sphere are removed from the pattern. Related to the extinction distance and roungly equal to 1/thickness.
shape_factor_width (float) – Determines the width of the reciprocal rel-rod, for fine-grained control. If not set will be set equal to max_excitation_error.
debye_waller_factors (dict of str:value pairs) – Maps element names to their temperature-dependent Debye-Waller factors.
- Returns:
The data associated with this structure and diffraction setup.
- Return type:
- calculate_profile_data(structure, reciprocal_radius=1.0, minimum_intensity=0.001, debye_waller_factors={})[source]#
Calculates a one dimensional diffraction profile for a structure.
- Parameters:
structure (diffpy.structure.structure.Structure) – The structure for which to calculate the diffraction profile.
reciprocal_radius (float) – The maximum radius of the sphere of reciprocal space to sample, in reciprocal angstroms.
minimum_intensity (float) – The minimum intensity required for a diffraction peak to be considered real. Deals with numerical precision issues.
debye_waller_factors (dict of str:value pairs) – Maps element names to their temperature-dependent Debye-Waller factors.
- Returns:
The diffraction profile corresponding to this structure and experimental conditions.
- Return type:
library_generator#
Diffraction pattern library generator and associated tools.
- class diffsims.generators.library_generator.DiffractionLibraryGenerator(electron_diffraction_calculator)[source]#
Bases:
object
Computes a library of electron diffraction patterns for specified atomic structures and orientations.
- get_diffraction_library(structure_library, calibration, reciprocal_radius, half_shape, with_direct_beam=True, max_excitation_error=0.01, shape_factor_width=None, debye_waller_factors={})[source]#
Calculates a dictionary of diffraction data for a library of crystal structures and orientations.
Each structure in the structure library is rotated to each associated orientation and the diffraction pattern is calculated each time.
Angles must be in the Euler representation (Z,X,Z) and in degrees
- Parameters:
structure_library (difffsims:StructureLibrary Object) – Dictionary of structures and associated orientations for which electron diffraction is to be simulated.
calibration (float) – The calibration of experimental data to be correlated with the library, in reciprocal Angstroms per pixel.
reciprocal_radius (float) – The maximum g-vector magnitude to be included in the simulations.
half_shape (tuple) – The half shape of the target patterns, for 144x144 use (72,72) etc
with_direct_beam (bool) – Include the direct beam in the library.
max_excitation_error (float) – The extinction distance for reflections, in reciprocal Angstroms.
shape_factor_width (float) – Determines the width of the shape functions of the reflections in Angstroms. If not set is equal to max_excitation_error.
debye_waller_factors (dict of str:value pairs) – Maps element names to their temperature-dependent Debye-Waller factors.
- Returns:
diffraction_library – Mapping of crystal structure and orientation to diffraction data objects.
- Return type:
DiffractionLibrary
- class diffsims.generators.library_generator.VectorLibraryGenerator(structure_library)[source]#
Bases:
object
Computes a library of diffraction vectors and pairwise inter-vector angles for a specified StructureLibrary.
- get_vector_library(reciprocal_radius)[source]#
Calculates a library of diffraction vectors and pairwise inter-vector angles for a library of crystal structures.
- Parameters:
reciprocal_radius (float) – The maximum g-vector magnitude to be included in the library.
- Returns:
vector_library – Mapping of phase identifier to phase information in dictionary format.
- Return type:
DiffractionVectorLibrary
rotation_list_generators#
Provides users with a range of gridding functions
- diffsims.generators.rotation_list_generators.get_beam_directions_grid(crystal_system, resolution, mesh='spherified_cube_edge')[source]#
Produces an array of beam directions, within the stereographic triangle of the relevant crystal system. The way the array is constructed is based on different methods of meshing the sphere [Cajaravelli2015] and can be specified through the mesh argument.
- Parameters:
crystal_system (str) – Allowed are: ‘cubic’,’hexagonal’,’trigonal’,’tetragonal’, ‘orthorhombic’,’monoclinic’,’triclinic’
resolution (float) – An angle in degrees representing the worst-case angular distance to a first nearest neighbor grid point.
mesh (str) – Type of meshing of the sphere that defines how the grid is created. Options are: uv_sphere, normalized_cube, spherified_cube_corner (default), spherified_cube_edge, icosahedral, random.
- Returns:
rotation_list
- Return type:
list of tuples
- diffsims.generators.rotation_list_generators.get_fundamental_zone_grid(resolution=2, point_group=None, space_group=None)[source]#
Generates an equispaced grid of rotations within a fundamental zone.
- Parameters:
resolution (float, optional) – The characteristic distance between a rotation and its neighbour (degrees)
point_group (orix.quaternion.symmetry.Symmetry, optional) – One of the 11 proper point groups, defaults to None
space_group (int, optional) – Between 1 and 231, defaults to None
- Returns:
rotation_list – Grid of rotations lying within the specified fundamental zone
- Return type:
list of tuples
- diffsims.generators.rotation_list_generators.get_grid_around_beam_direction(beam_rotation, resolution, angular_range=(0, 360))[source]#
Creates a rotation list of rotations for which the rotation is about given beam direction.
- Parameters:
beam_rotation (tuple) – A desired beam direction as a rotation (rzxz eulers), usually found via get_rotation_from_z_to_direction.
resolution (float) – The resolution of the grid (degrees).
angular_range (tuple) – The minimum (included) and maximum (excluded) rotation around the beam direction to be included.
- Returns:
rotation_list
- Return type:
list of tuples
Examples
>>> from diffsims.generators.zap_map_generator import get_rotation_from_z_to_direction >>> beam_rotation = get_rotation_from_z_to_direction(structure, [1, 1, 1]) >>> grid = get_grid_around_beam_direction(beam_rotation, 1)
- diffsims.generators.rotation_list_generators.get_list_from_orix(grid, rounding=2)[source]#
Converts an orix sample to a rotation list
- Parameters:
grid (orix.quaternion.rotation.Rotation) – A grid of rotations
rounding (int, optional) – The number of decimal places to retain, defaults to 2
- Returns:
rotation_list – A rotation list
- Return type:
list of tuples
- diffsims.generators.rotation_list_generators.get_local_grid(resolution=2, center=None, grid_width=10)[source]#
Generates a grid of rotations about a given rotation
- Parameters:
resolution (float, optional) – The characteristic distance between a rotation and its neighbour (degrees)
center (euler angle tuple or orix.quaternion.rotation.Rotation, optional) – The rotation at which the grid is centered. If None (default) uses the identity
grid_width (float, optional) – The largest angle of rotation away from center that is acceptable (degrees)
- Returns:
rotation_list
- Return type:
list of tuples
sphere_mesh_generators#
- diffsims.generators.sphere_mesh_generators.beam_directions_grid_to_euler(vectors)[source]#
Convert list of vectors representing zones to a list of Euler angles in the bunge convention with the constraint that phi1=0.
- Parameters:
vectors (numpy.ndarray (N, 3)) – N 3-dimensional vectors to convert to Euler angles
- Returns:
grid – Euler angles in bunge convention corresponding to each vector in degrees.
- Return type:
numpy.ndarray (N, 3)
Notes
The Euler angles represent the orientation of the crystal if that particular vector were parallel to the beam direction [001]. The additional constraint of phi1=0 means that this orientation is uniquely defined for most vectors. phi1 represents the rotation of the crystal around the beam direction and can be interpreted as the rotation of a particular diffraction pattern.
- diffsims.generators.sphere_mesh_generators.get_cube_mesh_vertices(resolution, grid_type='spherified_corner')[source]#
Return the (x, y, z) coordinates of the vertices of a cube mesh on a sphere. To generate the mesh, a cube is made to surround the sphere. The surfaces of the cube are subdivided into a grid. The vectors from the origin to these grid points are normalized to unit length. The grid on the cube can be generated in three ways, see grid_type and reference [Cajaravelli2015].
- Parameters:
- Returns:
points_in_cartesian – Rows are x, y, z where z is the 001 pole direction
- Return type:
numpy.ndarray (N,3)
Notes
The resolution determines the maximum angle between first nearest neighbor grid points, but to get an integer number of points between the cube face center and the edges, the number of grid points is rounded up. In practice this means that resolution is always an upper limit. Additionally, where on the grid this maximum angle will be will depend on the type of grid chosen. Resolution says something about the maximum angle but nothing about the distribution of nearest neighbor angles or the minimum angle - also this is fixed by the chosen grid.
In the normalized grid, the grid on the surface of the cube is linear. The maximum angle between nearest neighbors is found between the <001> directions and the first grid point towards the <011> directions. Points approaching the edges and corners of the cube will have a smaller angular deviation, so orientation space will be oversampled there compared to the cube faces <001>.
In the spherified_edge grid, the grid is constructed so that there are still two sets of perpendicular grid lines parallel to the {100} directions on each cube face, but the spacing of the grid lines is chosen so that the angles between the grid points on the line connecting the face centers (<001>) to the edges (<011>) are equal. The maximum angle is also between the <001> directions and the first grid point towards the <011> edges. This grid slightly oversamples the directions between <011> and <111>
The spherified_corner case is similar to the spherified_edge case, but the spacing of the grid lines is chosen so that the angles between the grid points on the line connecting the face centers to the cube corners (<111>) is equal. The maximum angle in this grid is from the corners to the first grid point towards the cube face centers.
References
[Cajaravelli2015] (1,2,3,4)O. S. Cajaravelli, “Four Ways to Create a Mesh for a Sphere,” https://medium.com/@oscarsc/four-ways-to-create-a-mesh-for-a-sphere-d7956b825db4.
- diffsims.generators.sphere_mesh_generators.get_icosahedral_mesh_vertices(resolution)[source]#
Return the (x, y, z) coordinates of the vertices of an icosahedral mesh of a cube, see [Cajaravelli2015]. Method was adapted from meshzoo [Meshzoo].
- Parameters:
resolution (float) – The maximum angle in degrees between neighboring grid points. Since the mesh is generated iteratively, the actual maximum angle in the mesh can be slightly smaller.
- Returns:
points_in_cartesian – Rows are x, y, z where z is the 001 pole direction
- Return type:
numpy.ndarray (N,3)
References
[Meshzoo]The meshzoo.sphere module.
- diffsims.generators.sphere_mesh_generators.get_random_sphere_vertices(resolution, seed=None)[source]#
Create a mesh that randomly samples the surface of a sphere
- Parameters:
- Returns:
points_in_cartesian – Rows are x, y, z where z is the 001 pole direction
- Return type:
numpy.ndarray (N,3)
References
- diffsims.generators.sphere_mesh_generators.get_uv_sphere_mesh_vertices(resolution)[source]#
Return the vertices of a UV (spherical coordinate) mesh on a unit sphere [Cajaravelli2015]. The mesh vertices are defined by the parametrization:
\[ \begin{align}\begin{aligned}x = sin(u)cos(v)\\y = sin(u)sin(v)\\z = cos(u)\end{aligned}\end{align} \]- Parameters:
resolution (float) – An angle in degrees. The maximum angle between nearest neighbor grid points. In this mesh this occurs on the equator of the sphere. All elevation grid lines are separated by at most resolution. The step size of u and v are rounded up to get an integer number of elevation and azimuthal grid lines with equal spacing.
- Returns:
points_in_cartesian – Rows are x, y, z where z is the 001 pole direction
- Return type:
numpy.ndarray (N,3)
zap_map_generator#
- diffsims.generators.zap_map_generator.corners_to_centroid_and_edge_centers(corners)[source]#
Produces the midpoints and center of a trio of corners
- diffsims.generators.zap_map_generator.generate_directional_simulations(structure, simulator, direction_list, reciprocal_radius=1, **kwargs)[source]#
Produces simulation of a structure aligned with certain axes
- Parameters:
structure (diffpy.structure.structure.Structure) – The structure from which simulations need to be produced.
simulator (DiffractionGenerator) – The diffraction generator object used to produce the simulations
direction_list (list of lists) – A list of [UVW] indices, eg. [[1,0,0],[1,1,0]]
reciprocal_radius (float) – Default to 1
- Returns:
direction_dictionary – Keys are zone axes, values are simulations
- Return type:
- diffsims.generators.zap_map_generator.generate_zap_map(structure, simulator, system='cubic', reciprocal_radius=1, density='7', **kwargs)[source]#
Produces a number of zone axis patterns for a structure
- Parameters:
structure (diffpy.structure.structure.Structure) – The structure to be simulated.
simulator (DiffractionGenerator) – The simulator used to generate the simulations
system (str) – ‘cubic’, ‘hexagonal’, ‘trigonal’, ‘tetragonal’, ‘orthorhombic’, ‘monoclinic’. Defaults to ‘cubic’.
reciprocal_radius (float) – The range of reciprocal lattice spots to be included. Default to 1.
density (str) – ‘3’ for the corners or ‘7’ (corners + midpoints + centroids). Defaults to 7.
kwargs – Keyword arguments to be passed to simulator.calculate_ed_data().
- Returns:
zap_dictionary – Keys are zone axes, values are simulations
- Return type:
Examples
Plot all of the patterns that you have generated
>>> zap_map = generate_zap_map(structure,simulator,'hexagonal',density='3') >>> for k in zap_map.keys(): >>> pattern = zap_map[k] >>> pattern.calibration = 4e-3 >>> plt.figure() >>> plt.imshow(pattern.get_diffraction_pattern(),vmax=0.02)
- diffsims.generators.zap_map_generator.get_rotation_from_z_to_direction(structure, direction)[source]#
Finds the rotation that takes [001] to a given zone axis.
- Parameters:
structure (diffpy.structure.structure.Structure) – The structure for which a rotation needs to be found.
direction (array like) – [UVW] direction that the ‘z’ axis should end up point down.
- Returns:
euler_angles – ‘rzxz’ in degrees.
- Return type:
Notes
This implementation works with an axis arrangement that has +x as left to right, +y as bottom to top and +z as out of the plane of a page. Rotations are counter clockwise as you look from the tip of the axis towards the origin
libraries#
Diffraction, structure and vector libraries.
diffraction_library#
- class diffsims.libraries.diffraction_library.DiffractionLibrary(*args, **kwargs)[source]#
Bases:
dict
Maps crystal structure (phase) and orientation to simulated diffraction data.
- identifiers#
A list of phase identifiers referring to different atomic structures.
- Type:
list of strings/ints
- structures#
A list of diffpy.structure.Structure objects describing the atomic structure associated with each phase in the library.
- Type:
list of diffpy.structure.Structure objects.
- diffraction_generator#
Diffraction generator used to generate this library.
- Type:
- get_library_entry(phase=None, angle=None)[source]#
Extracts a single DiffractionLibrary entry.
- Parameters:
- Returns:
library_entries – Dictionary containing the simulation associated with the specified phase and orientation with associated properties.
- Return type:
structure_library#
- class diffsims.libraries.structure_library.StructureLibrary(identifiers, structures, orientations)[source]#
Bases:
object
A dictionary containing all the structures and their associated rotations in the .struct_lib attribute.
- identifiers#
A list of phase identifiers referring to different atomic structures.
- Type:
list of strings/ints
- structures#
A list of diffpy.structure.Structure objects describing the atomic structure associated with each phase in the library.
- Type:
list of diffpy.structure.Structure objects.
- orientations#
A list over identifiers of lists of euler angles (as tuples) in the rzxz convention and in degrees.
- Type:
- classmethod from_crystal_systems(identifiers, structures, systems, resolution, equal='angle')[source]#
Creates a structure library from crystal system derived orientation lists
- Parameters:
identifiers (list of strings/ints) – A list of phase identifiers referring to different atomic structures.
structures (list of diffpy.structure.Structure objects.) – A list of diffpy.structure.Structure objects describing the atomic structure associated with each phase in the library.
systems (list) – A list over indentifiers of crystal systems
resolution (float) – resolution in degrees
equal (str) – Default is ‘angle’
- Raises:
NotImplementedError: – “This function has been removed in version 0.3.0, in favour of creation from orientation lists”
- classmethod from_orientation_lists(identifiers, structures, orientations)[source]#
Creates a structure library from “manual” orientation lists
- Parameters:
identifiers (list of strings/ints) – A list of phase identifiers referring to different atomic structures.
structures (list of diffpy.structure.Structure objects.) – A list of diffpy.structure.Structure objects describing the atomic structure associated with each phase in the library.
orientations (list of lists of tuples) – A list over identifiers of lists of euler angles (as tuples) in the rzxz convention and in degrees.
- Return type:
vector_library#
- class diffsims.libraries.vector_library.DiffractionVectorLibrary(*args, **kwargs)[source]#
Bases:
dict
Maps crystal structure (phase) to diffraction vectors.
The library is a dictionary mapping from a phase name to phase information. The phase information is stored as a dictionary with the following entries:
- ‘indices’np.array
List of peak indices [hkl1, hkl2] as a 2D array.
- ‘measurements’np.array
List of vector measurements [len1, len2, angle] in the same order as the indices. Lengths in reciprocal Angstrom and angles in radians.
- identifiers#
A list of phase identifiers referring to different atomic structures.
- Type:
list of strings/ints
sims#
Diffraction simulations.
diffraction_simulation#
- class diffsims.sims.diffraction_simulation.DiffractionSimulation(coordinates, indices=None, intensities=None, calibration=None, offset=(0.0, 0.0), with_direct_beam=False)[source]#
Bases:
object
Holds the result of a kinematic diffraction pattern simulation.
- Parameters:
coordinates (array-like, shape [n_points, 2]) – The x-y coordinates of points in reciprocal space.
indices (array-like, shape [n_points, 3]) – The indices of the reciprocal lattice points that intersect the Ewald sphere.
intensities (array-like, shape [n_points, ]) – The intensity of the reciprocal lattice points.
calibration (float or tuple of float, optional) – The x- and y-scales of the pattern, with respect to the original reciprocal angstrom coordinates.
offset (tuple of float, optional) – The x-y offset of the pattern in reciprocal angstroms. Defaults to zero in each direction.
- property calibrated_coordinates#
Coordinates converted into pixel space.
- Type:
ndarray
- property calibration#
The x- and y-scales of the pattern, with respect to the original reciprocal angstrom coordinates.
- property coordinates#
The coordinates of all unmasked points.
- Type:
ndarray
- property direct_beam_mask#
If with_direct_beam is True, returns a True array for all points. If with_direct_beam is False, returns a True array with False in the position of the direct beam.
- Type:
ndarray
- get_as_mask(shape, radius=6.0, negative=True, radius_function=None, direct_beam_position=None, in_plane_angle=0, mirrored=False, *args, **kwargs)[source]#
Return the diffraction pattern as a binary mask of type bool
- Parameters:
shape (2-tuple of ints) – Shape of the output mask (width, height)
radius (float or array, optional) – Radii of the spots in pixels. An array may be supplied of the same length as the number of spots.
negative (bool, optional) – Whether the spots are masked (True) or everything else is masked (False)
radius_function (Callable, optional) – Calculate the radius as a function of the spot intensity, for example np.sqrt. args and kwargs supplied to this method are passed to this function. Will override radius.
direct_beam_position (2-tuple of ints, optional) – The (x,y) coordinate in pixels of the direct beam. Defaults to the center of the image.
in_plane_angle (float, optional) – In plane rotation of the pattern in degrees
mirrored (bool, optional) – Whether the pattern should be flipped over the x-axis, corresponding to the inverted orientation
- Returns:
mask – Boolean mask of the diffraction pattern
- Return type:
- get_diffraction_pattern(shape=(512, 512), sigma=10, direct_beam_position=None, in_plane_angle=0, mirrored=False)[source]#
Returns the diffraction data as a numpy array with two-dimensional Gaussians representing each diffracted peak. Should only be used for qualitative work.
- Parameters:
shape (tuple of ints) – The size of a side length (in pixels)
sigma (float) – Standard deviation of the Gaussian function to be plotted (in pixels).
direct_beam_position (2-tuple of ints, optional) – The (x,y) coordinate in pixels of the direct beam. Defaults to the center of the image.
in_plane_angle (float, optional) – In plane rotation of the pattern in degrees
mirrored (bool, optional) – Whether the pattern should be flipped over the x-axis, corresponding to the inverted orientation
- Returns:
diffraction-pattern – The simulated electron diffraction pattern, normalised.
- Return type:
numpy.array
Notes
If don’t know the exact calibration of your diffraction signal using 1e-2 produces reasonably good patterns when the lattice parameters are on the order of 0.5nm and a the default size and sigma are used.
- property indices#
- property intensities#
The intensities of all unmasked points.
- Type:
ndarray
- plot(size_factor=1, direct_beam_position=None, in_plane_angle=0, mirrored=False, units='real', show_labels=False, label_offset=(0, 0), label_formatting={}, ax=None, **kwargs)[source]#
A quick-plot function for a simulation of spots
- Parameters:
size_factor (float, optional) – linear spot size scaling, default to 1
direct_beam_position (2-tuple of ints, optional) – The (x,y) coordinate in pixels of the direct beam. Defaults to the center of the image.
in_plane_angle (float, optional) – In plane rotation of the pattern in degrees
mirrored (bool, optional) – Whether the pattern should be flipped over the x-axis, corresponding to the inverted orientation
units (str, optional) – ‘real’ or ‘pixel’, only changes scalebars, falls back on ‘real’, the default
show_labels (bool, optional) – draw the miller indices near the spots
label_offset (2-tuple, optional) – the relative location of the spot labels. Does nothing if show_labels is False.
label_formatting (dict, optional) – keyword arguments passed to ax.text for drawing the labels. Does nothing if show_labels is False.
ax (matplotlib Axes, optional) – axes on which to draw the pattern. If None, a new axis is created
**kwargs – passed to ax.scatter() method
- Return type:
ax,sp
Notes
spot size scales with the square root of the intensity.
- rotate_shift_coordinates(angle, center=(0, 0), mirrored=False)[source]#
Rotate, flip or shift patterns in-plane
- property size#
- class diffsims.sims.diffraction_simulation.ProfileSimulation(magnitudes, intensities, hkls)[source]#
Bases:
object
Holds the result of a given kinematic simulation of a diffraction profile.
- Parameters:
magnitudes (array-like, shape [n_peaks, 1]) – Magnitudes of scattering vectors.
intensities (array-like, shape [n_peaks, 1]) – The kinematic intensity of the diffraction peaks.
hkls ([{(h, k, l): mult}] {(h, k, l): mult} is a dict of Miller) – indices for all diffracted lattice facets contributing to each intensity.
- get_plot(annotate_peaks=True, with_labels=True, fontsize=12)[source]#
- Plots the diffraction profile simulation for the
calculate_profile_data method in DiffractionGenerator.
- Parameters:
annotate_peaks (boolean) – If True, peaks are annotaed with hkl information.
with_labels (boolean) – If True, xlabels and ylabels are added to the plot.
fontsize (integer) – Fontsize for peak labels.
structure_factor#
Calculation of scattering factors and structure factors.
- diffsims.structure_factor.find_asymmetric_positions(positions, space_group)[source]#
Return the asymmetric atom positions among a set of positions when considering symmetry operations defined by a space group.
- Parameters:
positions (list) – A list of cartesian atom positions.
space_group (diffpy.structure.spacegroupmod.SpaceGroup) – Space group describing the symmetry operations.
- Returns:
Asymmetric atom positions.
- Return type:
- diffsims.structure_factor.get_atomic_scattering_parameters(element, unit=None)[source]#
Return the eight atomic scattering parameters a_1-4, b_1-4 for elements with atomic numbers Z = 1-98 from Table 12.1 in [DeGraef2007], which are themselves from [Doyle1968] and [Smith1962].
- Parameters:
- Returns:
a (numpy.ndarray) – The four atomic scattering parameters a_1-4.
b (numpy.ndarray) – The four atomic scattering parameters b_1-4.
References
[DeGraef2007]De Graef, M. E. McHenry, “Structure of Materials,” Cambridge University Press (2007).
[Doyle1968] (1,2,3)P. A. Doyle, P. S. Turner, “Relativistic Hartree-Fock X-ray and electron scattering factors,” Acta Cryst. 24 (1968), doi: https://doi.org/10.1107/S0567739468000756.
[Smith1962]G. Smith, R. Burge, “The analytical representation of atomic scattering amplitudes for electrons,” Acta Cryst. A15 (1962), doi: https://doi.org/10.1107/S0365110X62000481.
- diffsims.structure_factor.get_doyleturner_atomic_scattering_factor(atom, scattering_parameter, unit_cell_volume)[source]#
Return the atomic scattering factor f for a certain atom and scattering parameter using Doyle-Turner atomic scattering parameters [Doyle1968].
Assumes structure’s Debye-Waller factors are expressed in Ångströms.
This function is adapted from EMsoft.
- Parameters:
atom (diffpy.structure.atom.Atom) – Atom with element type, Debye-Waller factor and occupancy number.
scattering_parameter (float) – The scattering parameter s for these Miller indices describing the crystal plane in which the atom lies.
unit_cell_volume (float) – Volume of the unit cell.
- Returns:
f – Scattering factor for this atom on this plane.
- Return type:
- diffsims.structure_factor.get_doyleturner_structure_factor(phase, hkl, scattering_parameter, voltage, return_parameters=False)[source]#
Return the structure factor for a given family of Miller indices using Doyle-Turner atomic scattering parameters [Doyle1968].
Assumes structure’s lattice parameters and Debye-Waller factors are expressed in Ångströms.
This function is adapted from EMsoft.
- Parameters:
phase (orix.crystal_map.phase_list.Phase) – A phase container with a crystal structure and a space and point group describing the allowed symmetry operations.
hkl (numpy.ndarray or list) – Miller indices.
scattering_parameter (float) – Scattering parameter for these Miller indices.
voltage (float) – Beam energy in V.
return_parameters (bool, optional) – Whether to return a set of parameters derived from the calculation as a dictionary. Default is False.
- Returns:
structure_factor (float) – Structure factor F.
params (dict) – A dictionary with (key, item) (str, float) of parameters derived from the calculation. Only returned if return_parameters=True.
- diffsims.structure_factor.get_element_id_from_string(element_str)[source]#
Get periodic element ID for elements \(Z\) = 1-98 from one-two letter string.
- diffsims.structure_factor.get_kinematical_atomic_scattering_factor(atom, scattering_parameter)[source]#
Return the kinematical (X-ray) atomic scattering factor f for a certain atom and scattering parameter.
Assumes structure’s Debye-Waller factors are expressed in Ångströms.
This function is adapted from EMsoft.
- Parameters:
atom (diffpy.structure.atom.Atom) – Atom with element type, Debye-Waller factor and occupancy number.
scattering_parameter (float) – The scattering parameter s for these Miller indices describing the crystal plane in which the atom lies.
- Returns:
f – Scattering factor for this atom on this plane.
- Return type:
- diffsims.structure_factor.get_kinematical_structure_factor(phase, hkl, scattering_parameter)[source]#
Return the kinematical (X-ray) structure factor for a given family of Miller indices.
Assumes structure’s lattice parameters and Debye-Waller factors are expressed in Ångströms.
This function is adapted from EMsoft.
- Parameters:
phase (orix.crystal_map.phase_list.Phase) – A phase container with a crystal structure and a space and point group describing the allowed symmetry operations.
hkl (numpy.ndarray or list) – Miller indices.
scattering_parameter (float) – Scattering parameter for these Miller indices.
- Returns:
structure_factor – Structure factor F.
- Return type:
- diffsims.structure_factor.get_refraction_corrected_wavelength(phase, voltage)[source]#
Return the refraction corrected relativistic electron wavelength in Ångströms for a given crystal structure and beam energy in V.
This function is adapted from EMsoft.
- Parameters:
phase (orix.crystal_map.Phase) – A phase container with a crystal structure and a space and point group describing the allowed symmetry operations.
voltage (float) – Beam energy in V.
- Returns:
wavelength – Refraction corrected relativistic electron wavelength in Ångströms.
- Return type:
pattern#
detector_functions#
- diffsims.pattern.detector_functions.add_dead_pixels(pattern, n=None, fraction=None, seed=None)[source]#
Adds randomly placed dead pixels onto a pattern
- Parameters:
pattern (numpy.ndarray) – The diffraction pattern at the detector
n (int) – The number of dead pixels, defaults to None
fraction (float) – The fraction of dead pixels, defaults to None
seed (int or None) – seed value for the random number generator
- Returns:
corrupted_pattern – The pattern, with dead pixels included
- Return type:
- diffsims.pattern.detector_functions.add_detector_offset(pattern, offset)[source]#
Adds/subtracts a fixed offset value from a pattern
- Parameters:
pattern (numpy.ndarray) – The diffraction pattern at the detector
offset (float or numpy.ndarray) – Added through the pattern, broadcasting applies
- Returns:
corrupted_pattern – The pattern, with offset applied, pixels that would have been negative are instead 0.
- Return type:
np.ndarray
- diffsims.pattern.detector_functions.add_gaussian_noise(pattern, sigma, seed=None)[source]#
Applies gaussian noise at each pixel within the pattern
- Parameters:
pattern (numpy.ndarray) – The diffraction pattern at the detector
sigma (float) – The (absolute) deviation of the gaussian errors
seed (int or None) – seed value for the random number generator
- Return type:
corrupted_pattern
- diffsims.pattern.detector_functions.add_gaussian_point_spread(pattern, sigma)[source]#
Blurs intensities across space with a gaussian function
- Parameters:
pattern (numpy.ndarray) – The diffraction pattern at the detector
sigma (float) – The standard deviation of the gaussian blur, in pixels
- Returns:
blurred_pattern – The blurred pattern (deterministic)
- Return type:
- diffsims.pattern.detector_functions.add_linear_detector_gain(pattern, gain)[source]#
Multiplies the pattern by a gain (which is not a function of the pattern)
- Parameters:
pattern (numpy.ndarray) – The diffraction pattern at the detector
gain (float or numpy.ndarray) – Multiplied through the pattern, broadcasting applies
- Returns:
corrupted_pattern – The pattern, with gain applied
- Return type:
- diffsims.pattern.detector_functions.add_shot_and_point_spread(pattern, sigma, shot_noise=True, seed=None)[source]#
Adds shot noise (optional) and gaussian point spread (via a convolution) to a pattern
- Parameters:
pattern (numpy.ndarray) – The diffraction pattern at the detector
sigma (float) – The standard deviation of the gaussian blur, in pixels
shot_noise (bool) – Whether to include shot noise in the original signal, default True
seed (int or None) – seed value for the random number generator (effects the shot noise only)
- Returns:
detector_pattern – A single sample of the pattern after accounting for detector properties
- Return type:
See also
add_shot_noise
adds only shot noise
add_gaussian_point_spread
adds only point spread
- diffsims.pattern.detector_functions.add_shot_noise(pattern, seed=None)[source]#
Applies shot noise to a pattern
- Parameters:
pattern (numpy.ndarray) – The diffraction pattern at the detector
seed (int or None) – seed value for the random number generator
- Returns:
shotted_pattern – A single sample of the pattern after accounting for shot noise
- Return type:
Notes
This function will (as it should) behave differently depending on the pattern intensity, so be mindful to put your intensities in physical units
- diffsims.pattern.detector_functions.constrain_to_dynamic_range(pattern, detector_max=None)[source]#
Force the values within pattern to lie between [0,detector_max]
- Parameters:
pattern (numpy.ndarray) – The diffraction pattern at the detector after corruption
detector_max (float) – The maximum allowed value at the detector
- Returns:
within_range_pattern – The pattern, with values >=0 and =< detector_max
- Return type:
utils#
Diffraction utilities used by the other modules.
atomic_diffraction_generator_utils#
Back end for computing diffraction patterns with a kinematic model.
- diffsims.utils.atomic_diffraction_generator_utils.get_diffraction_image(coordinates, species, probe, x, wavelength, precession, GPU=True, pointwise=False, **kwargs)[source]#
Return kinematically simulated diffraction pattern
- Parameters:
coordinates (numpy.ndarray [float], (n_atoms, 3)) – List of atomic coordinates
species (numpy.ndarray [int], (n_atoms,)) – List of atomic numbers
probe (diffsims.ProbeFunction) – Function representing 3D shape of beam
x (list [numpy.ndarray [float] ], of shapes [(nx,), (ny,), (nz,)]) – Mesh on which to compute the volume density
wavelength (float) – Wavelength of electron beam
precession (a pair (float, int)) – The float dictates the angle of precession and the int how many points are used to discretise the integration.
dtype ((str, str)) – tuple of floating/complex datatypes to cast outputs to
ZERO (float > 0, optional) – Rounding error permitted in computation of atomic density. This value is the smallest value rounded to 0.
GPU (bool, optional) – Flag whether to use GPU or CPU discretisation. Default (if available) is True
pointwise (bool, optional) – Optional parameter whether atomic intensities are computed point-wise at the centre of a voxel or an integral over the voxel. default=False
- Returns:
DP – The two-dimensional diffraction pattern evaluated on the reciprocal grid corresponding to the first two vectors of x.
- Return type:
numpy.ndarray [dtype[0]], (nx, ny, nz)
- diffsims.utils.atomic_diffraction_generator_utils.grid2sphere(arr, x, dx, C)[source]#
Projects 3d array onto a sphere
- Parameters:
arr (np.ndarray [float], (nx, ny, nz)) – Input function to be projected
x (list [np.ndarray [float]], of shapes [(nx,), (ny,), (nz,)]) – Vectors defining mesh of <arr>
dx (list [np.ndarray [float]], of shapes [(3,), (3,), (3,)]) – Basis in which to orient sphere. Centre of sphere will be at C*dx[2] and mesh of output array will be defined by the first two vectors
C (float) – Radius of sphere
- Returns:
out – If y is the point on the line between i*dx[0]+j*dx[1] and C*dx[2] which also lies on the sphere of radius C from C*dx[2] then: out[i,j] = arr(y). Interpolation on arr is linear.
- Return type:
np.ndarray [float], (nx, ny)
- diffsims.utils.atomic_diffraction_generator_utils.precess_mat(alpha, theta)[source]#
Generates rotation matrices for precession curves.
- Parameters:
alpha (float) – Angle (in degrees) of precession tilt
theta (float) – Angle (in degrees) along precession curve
- Returns:
R – Rotation matrix associated to the tilt of alpha away from the vertical axis and a rotation of theta about the vertical axis.
- Return type:
numpy.ndarray [float], (3, 3)
atomic_scattering_params#
discretise_utils#
Utils for converting lists of atoms to a discretised volume
- diffsims.utils.discretise_utils.do_binning(x, loc, Rmax, d, GPU)[source]#
Utility function which takes in a mesh, atom locations, atom radius and minimal grid-spacing and returns a binned array of atom indices.
- Parameters:
x (list [np.ndarray [float]], of shape [(nx,), (ny,), ...]) – Dictates the range of the box over which to bin atoms.
loc (np.ndarray, (n, 3)) – Atoms to bin.
Rmax (float > 3) – Maximum radius of an atom (rounded up to 3).
d (list of float > 0) – The finest permitted binning.
GPU (bool) – If True then constrains to memory of GPU rather than RAM.
- Returns:
subList (np.ndarray [int]) – subList[i0,i1,i2] is a list of indices [j0, j1, …, jn, -1,…] such that the only atoms which are contained in the box: [x[0].min()+i0*r,x[0].min(),+(i0+1)*r] x [x[1].min()+i1*r,x[1].min(),+(i1+1)*r]…
r (np.ndarray [float]) – Size of each bin.
mem (int) – Upper limit of memory in bytes.
- diffsims.utils.discretise_utils.get_atoms(Z, returnFunc=True, dtype='f8')[source]#
This function returns an approximation of the atom with atomic number Z using a list of Gaussians.
- Parameters:
- Returns:
obj1, obj2 – Continuous atom is represented by: .. math:: ymapsto sum_i a[i]*exp(-b[i]*|y|^2)
- Return type:
numpy.ndarray or function
This is data table 3 from ‘Robust Parameterization of Elastic and Absorptive Electron Atomic Scattering Factors’ by L.-M. Peng, G. Ren, S. L. Dudarev and M. J. Whelan, 1996
- diffsims.utils.discretise_utils.get_discretisation(loc, Z, x, GPU=False, ZERO=None, dtype=('f8', 'c16'), pointwise=False, FT=False, **kwargs)[source]#
- Parameters:
loc (numpy.ndarray, (n, 3)) – Atoms to bin
Z (str, int, or numpy.ndarray [str or int], (n,)) – atom labels, either string or atomic masses.
x (list [numpy.ndarray [float]]) – Dictates mesh over which to discretise. Volume will be discretised at points [x[0][i],x[1][j],…]
GPU (bool, optional) – If True (default) then attempts to use the GPU.
ZERO (float > 0) – Approximation threshold
dtype ((str, str), optional) – Real and complex data precisions to use, default=(‘float64’, ‘complex128’)
pointwise (bool, optional) – If True (default) then computes pointwise atomic potentials on mesh points, else averages the potential over cube of same size as the discretisation.
FT (bool, optional) – If True then computes the Fourier transform directly on the reciprocal mesh, otherwise (default) computes the volume potential
- Returns:
out – Discretisation of atoms defined by loc/Z on mesh defined by x.
- Return type:
numpy.ndarray, (x[0].size, x[1].size, x[2].size)
- diffsims.utils.discretise_utils.rebin(x, loc, r, k, mem)[source]#
Bins each location into a grid subject to memory constraints.
- Parameters:
x (list [np.ndarray [float]], of shape [(nx,), (ny,), ...]) – Dictates the range of the box over which to bin atoms.
loc (np.ndarray, (n, 3)) – Atoms to bin.
r (float or [float, float, float]) – Mesh size (in each direction).
k (int) – Integer such that the radius of the atom is <= k*r. Consequently, each atom will appear in approximately 8k^3 bins.
mem (int) – Upper limit of number of bytes permitted for mesh. If not possible then raises a MemoryError.
- Returns:
subList – subList[i0,i1,i2] is a list of indices [j0, j1, …, jn, -1,…] such that the only atoms which are contained in the box: [x[0].min()+i0*r,x[0].min(),+(i0+1)*r] x [x[1].min()+i1*r,x[1].min(),+(i1+1)*r]… are the atoms with locations loc[j0], …, loc[jn].
- Return type:
np.ndarray [int]
fourier_transform#
Created on 31 Oct 2019
Module provides optimised fft and Fourier transform approximation.
@author: Rob Tovey
- diffsims.utils.fourier_transform.convolve(arr1, arr2, dx=None, axes=None)[source]#
Performs a centred convolution of input arrays
- Parameters:
arr1 (numpy.ndarray) – Arrays to be convolved. If dimensions are not equal then 1s are appended to the lower dimensional array. Otherwise, arrays must be broadcastable.
arr2 (numpy.ndarray) – Arrays to be convolved. If dimensions are not equal then 1s are appended to the lower dimensional array. Otherwise, arrays must be broadcastable.
dx (float > 0, list of float, or None , optional) – Grid spacing of input arrays. Output is scaled by dx**max(arr1.ndim, arr2.ndim). default=`None` applies no scaling
axes (tuple of ints or None, optional) – Choice of axes to convolve. default=`None` convolves all axes
- diffsims.utils.fourier_transform.fast_abs(x, y=None)[source]#
Fast computation of abs of an array
- Parameters:
x (numpy.ndarray) – Input
y (numpy.ndarray or None, optional) – If y is not None, used as preallocated output
- Returns:
y – Array equal to abs(x)
- Return type:
numpy.ndarray
- diffsims.utils.fourier_transform.fast_fft_len(n)[source]#
Returns the smallest integer greater than input such that the fft can be computed efficiently at this size
- Parameters:
n (int) – minimum size
- Returns:
N – smallest integer greater than n which permits efficient ffts.
- Return type:
int
- diffsims.utils.fourier_transform.fftshift_phase(x)[source]#
Fast implementation of fft_shift: fft(fftshift_phase(x)) = fft_shift(fft(x))
Note two things: - this is an in-place manipulation of the (3D) input array - the input array must have even side lengths. This is softly guaranteed by fast_fft_len but will raise error if not true.
- diffsims.utils.fourier_transform.from_recip(y)[source]#
Converts Fourier frequencies to spatial coordinates.
- Parameters:
y (list [numpy.ndarray [float]], of shape [(nx,), (ny,), …]) – List (or equivalent) of vectors which define a mesh in the dimension equal to the length of x
- Returns:
x – List of vectors defining a mesh such that for a function, f, defined on the mesh given by y, ifft(f) is defined on the mesh given by x. 0 will be in the middle of x.
- Return type:
list [numpy.ndarray [float]], of shape [(nx,), (ny,), …]
- diffsims.utils.fourier_transform.get_DFT(X=None, Y=None)[source]#
Returns discrete analogues for the Fourier/inverse Fourier transform pair defined from grid X to grid Y and back again.
- Parameters:
X (list [numpy.ndarray [float]], of shape [(nx,), (ny,), …], optional) – Mesh on real space
Y (list [numpy.ndarray [float]], of shape [(nx,), (ny,), …], optional) – Corresponding mesh on Fourier space
other (If either X or Y is None then it is inferred from the) –
- Returns:
DFT (function(f, axes=None)) – If f is a function on X then DFT(f) is the Fourier transform of f on Y. axes parameter can be used to specify which axes to transform.
iDFT (function(f, axes=None)) – If f is a function on Y then iDFT(f) is the inverse Fourier transform of f on X. axes parameter can be used to specify which axes to transform.
- diffsims.utils.fourier_transform.get_recip_points(ndim, n=None, dX=inf, rX=0, dY=inf, rY=1e-16)[source]#
Returns a minimal pair of real and Fourier grids which satisfy each given requirement.
- Parameters:
ndim (int) – Dimension of domain
n (int, list of length ndim, or None , optional) – Sugested number of pixels (per dimension). default=`None` infers this from other parameters. If enough other constraints are given to define a discretisation then this will be shrunk if possible.
dX (float > 0 or list of float of length ndim, optional) – Maximum grid spacing (per dimension). default=`numpy.inf` infers this from other parameters
rX (float > 0 or list of float of length ndim, optional) – Minimum grid range (per dimension). default=`None` infers this from other parameters. In this case, range is maximal span, i.e. diameter.
dY (float > 0 or list of float of length ndim) – Maximum grid spacing (per dimension) in Fourier domain. default=`None` infers this from other parameters
rY (float > 0 or list of float of length ndim) – Minimum grid range (per dimension) in Fourier domain. default=`None` infers this from other parameters. In this case, range is maximal span, i.e. diameter.
- Returns:
x (list [numpy.ndarray [float]], of shape [(nx,), (ny,), …]) – Real mesh of points, centred at 0 with at least n pixels, resolution higher than dX, and range greater than rX.
y (list [numpy.ndarray [float]], of shape [(nx,), (ny,), …]) – Fourier mesh of points, centred at 0 with at least n pixels, resolution higher than dY, and range greater than rY.
- diffsims.utils.fourier_transform.plan_fft(A, n=None, axis=None, norm=None, **_)[source]#
Plans an fft for repeated use. Parameters are the same as for pyfftw’s fftn which are, where possible, the same as the numpy equivalents. Note that some functionality is only possible when using the pyfftw backend.
- Parameters:
A (numpy.ndarray, of dimension d) – Array of same shape to be input for the fft
n (iterable or None, len(n) == d, optional) – The output shape of fft (default=`None` is same as A.shape)
axis (int, iterable length d, or None, optional) – The axis (or axes) to transform (default=`None` is all axes)
overwrite (bool, optional) – Whether the input array can be overwritten during computation (default=False)
planner ({0, 1, 2, 3}, optional) – Amount of effort put into optimising Fourier transform where 0 is low and 3 is high (default=`1`).
threads (int, None) – Number of threads to use (default=`None` is all threads)
auto_align_input (bool, optional) – If True then may re-align input (default=`True`)
auto_contiguous (bool, optional) – If True then may re-order input (default=`True`)
avoid_copy (bool, optional) – If True then may over-write initial input (default=`False`)
norm ({None, 'ortho'}, optional) – Indicate whether fft is normalised (default=`None`)
- Returns:
plan (function) – Returns the Fourier transform of B, plan() == fftn(B)
B (numpy.ndarray, A.shape) – Array which should be modified inplace for fft to be computed. If possible, B is A.
Example
A = numpy.zeros((8,16)) plan, B = plan_fft(A)
B[:,:] = numpy.random.rand(8,16) numpy.fft.fftn(B) == plan()
B = numpy.random.rand(8,16) numpy.fft.fftn(B) != plan()
- diffsims.utils.fourier_transform.plan_ifft(A, n=None, axis=None, norm=None, **_)[source]#
Plans an ifft for repeated use. Parameters are the same as for pyfftw’s ifftn which are, where possible, the same as the numpy equivalents. Note that some functionality is only possible when using the pyfftw backend.
- Parameters:
A (numpy.ndarray, of dimension d) – Array of same shape to be input for the ifft
n (iterable or None, len(n) == d, optional) – The output shape of ifft (default=`None` is same as A.shape)
axis (int, iterable length d, or None, optional) – The axis (or axes) to transform (default=`None` is all axes)
overwrite (bool, optional) – Whether the input array can be overwritten during computation (default=False)
planner ({0, 1, 2, 3}, optional) – Amount of effort put into optimising Fourier transform where 0 is low and 3 is high (default=`1`).
threads (int, None) – Number of threads to use (default=`None` is all threads)
auto_align_input (bool, optional) – If True then may re-align input (default=`True`)
auto_contiguous (bool, optional) – If True then may re-order input (default=`True`)
avoid_copy (bool, optional) – If True then may over-write initial input (default=`False`)
norm ({None, 'ortho'}, optional) – Indicate whether ifft is normalised (default=`None`)
- Returns:
plan (function) – Returns the inverse Fourier transform of B, plan() == ifftn(B)
B (numpy.ndarray, A.shape) – Array which should be modified inplace for ifft to be computed. If possible, B is A.
- diffsims.utils.fourier_transform.to_recip(x)[source]#
Converts spatial coordinates to Fourier frequencies.
- Parameters:
x (list [numpy.ndarray [float]], of shape [(nx,), (ny,), …]) – List (or equivalent) of vectors which define a mesh in the dimension equal to the length of x
- Returns:
y – List of vectors defining a mesh such that for a function, f, defined on the mesh given by x, fft(f) is defined on the mesh given by y
- Return type:
list [numpy.ndarray [float]], of shape [(nx,), (ny,), …]
generic_utils#
Created on 31 Oct 2019
Generic tools for all areas of code.
@author: Rob Tovey
- class diffsims.utils.generic_utils.GLOBAL_BOOL(val)[source]#
Bases:
object
An object which behaves like a bool but can be changed in-place by set or by calling as a function.
- diffsims.utils.generic_utils.to_mesh(x, dx=None, dtype=None)[source]#
- Generates dense meshes from grid vectors, broadly:
to_mesh(x)[i,j,…] = (x[0][i], x[1][j], …)
- Parameters:
x (list [numpy.ndarray], of shape [(nx,), (ny,), …]) – List of grid vectors
dx (list [numpy.ndarray] or None, optional) – Basis in which to expand mesh, default is the canonical basis
dtype (str or None, optional) – String representing the numpy type of output, default inherits from x
- Returns:
X – X[i,j,…, k] = x[0][i]*dx[0][k] + x[1][j]*dx[1][k] + …
- Return type:
numpy.ndarray [dtype], (x[0].size, x[1].size, …, len(x))
kinematic_simulation_utils#
Created on 1 Nov 2019
Back end for computing diffraction patterns with a kinematic model.
@author: Rob Tovey
- diffsims.utils.kinematic_simulation_utils.get_diffraction_image(coordinates, species, probe, x, wavelength, precession, GPU=True, pointwise=False, **kwargs)[source]#
Return kinematically simulated diffraction pattern
- Parameters:
coordinates (numpy.ndarray [float], (n_atoms, 3)) – List of atomic coordinates
species (numpy.ndarray [int], (n_atoms,)) – List of atomic numbers
probe (diffsims.ProbeFunction) – Function representing 3D shape of beam
x (list [numpy.ndarray [float] ], of shapes [(nx,), (ny,), (nz,)]) – Mesh on which to compute the volume density
wavelength (float) – Wavelength of electron beam
precession (a pair (float, int)) – The float dictates the angle of precession and the int how many points are used to discretise the integration.
dtype ((str, str)) – tuple of floating/complex datatypes to cast outputs to
ZERO (float > 0, optional) – Rounding error permitted in computation of atomic density. This value is the smallest value rounded to 0.
GPU (bool, optional) – Flag whether to use GPU or CPU discretisation. Default (if available) is True
pointwise (bool, optional) – Optional parameter whether atomic intensities are computed point-wise at the centre of a voxel or an integral over the voxel. default=False
- Returns:
DP – The two-dimensional diffraction pattern evaluated on the reciprocal grid corresponding to the first two vectors of x.
- Return type:
numpy.ndarray [dtype[0]], (nx, ny, nz)
- diffsims.utils.kinematic_simulation_utils.grid2sphere(arr, x, dx, C)[source]#
Projects 3d array onto a sphere.
- Parameters:
arr (np.ndarray [float], (nx, ny, nz)) – Input function to be projected
x (list [np.ndarray [float]], of shapes [(nx,), (ny,), (nz,)]) – Vectors defining mesh of <arr>
dx (list [np.ndarray [float]], of shapes [(3,), (3,), (3,)]) – Basis in which to orient sphere. Centre of sphere will be at C*dx[2] and mesh of output array will be defined by the first two vectors.
C (float) – Radius of sphere.
- Returns:
out – If y is the point on the line between i*dx[0]+j*dx[1] and C*dx[2] which also lies on the sphere of radius C from C*dx[2] then: out[i,j] = arr(y). Interpolation on arr is linear.
- Return type:
np.ndarray [float], (nx, ny)
- diffsims.utils.kinematic_simulation_utils.precess_mat(alpha, theta)[source]#
Generates rotation matrices for precession curves.
- Parameters:
alpha (float) – Angle (in degrees) of precession tilt
theta (float) – Angle (in degrees) along precession curve
- Returns:
R – Rotation matrix associated to the tilt of alpha away from the vertical axis and a rotation of theta about the vertical axis.
- Return type:
numpy.ndarray [float], (3, 3)
lobato_scattering_params#
probe_utils#
Created on 5 Nov 2019
@author: Rob Tovey
- class diffsims.utils.probe_utils.ProbeFunction(func=None)[source]#
Bases:
object
Object representing a probe function.
- Parameters:
func (function) – Function which takes in an array, r, of shape [nx, ny, nz, 3] and returns an array of shape [nx, ny, nz]. r[…,0] corresponds to the x coordinate, r[…, 1] to y etc. If not provided (or None) then the __call__ and FT methods must be overrided.
- __call__(x, out=None, scale=None)[source]#
Returns func(x)*scale. If out is provided then it is used as preallocated storage. If scale is not provided then it is assumed to be 1. If x is a list of arrays it is converted into a mesh first.
- Parameters:
x (numpy.ndarray, (nx, ny, nz, 3) or list of arrays of shape) – [(nx,), (ny,), (nz,)] Mesh points at which to evaluate the probe density.
out (numpy.ndarray, (nx, ny, nz), optional) – If provided then computation is performed inplace.
scale (numpy.ndarray, (nx, ny, nz), optional) – If provided then the probe density is scaled by this before being returned.
- Returns:
out – An array equal to probe(x)*scale.
- Return type:
numpy.ndarray, (nx, ny, nz)
- FT(y, out=None)[source]#
Returns the Fourier transform of func on the mesh y. Again, if out is provided then computation is inplace. If y is a list of arrays then it is converted into a mesh first. If this function is not overridden then an approximation is made using func and the fft.
- Parameters:
y (numpy.ndarray, (nx, ny, nz, 3) or list of arrays of shape) – [(nx,), (ny,), (nz,)] Mesh of Fourier coordinates at which to evaluate the probe density.
out (numpy.ndarray, (nx, ny, nz), optional) – If provided then computation is performed inplace.
- Returns:
out – An array equal to FourierTransform(probe)(y).
- Return type:
numpy.ndarray, (nx, ny, nz)
- class diffsims.utils.probe_utils.BesselProbe(r)[source]#
Bases:
ProbeFunction
Probe function given by a radially scaled Bessel function of the first kind.
- Parameters:
r (float) – Width of probe at the surface of the sample. More specifically, the smallest 0 of the probe.
- __call__(x, out=None, scale=None)[source]#
If X = sqrt(x[…,0]**2+x[…,1]**2)/r then returns J_1(X)/X*scale. If out is provided then this is computed inplace. If scale is not provided then it is assumed to be 1. If x is a list of arrays it is converted into a mesh first.
- Parameters:
x (numpy.ndarray, (nx, ny, nz, 3) or list of arrays of shape) – [(nx,), (ny,), (nz,)] Mesh points at which to evaluate the probe density. As a plotting utility, if a lower dimensional mesh is provided then the remaining coordinates are assumed to be 0 and so only the respective 1D/2D slice of the probe is returned.
out (numpy.ndarray, (nx, ny, nz), optional) – If provided then computation is performed inplace.
scale (numpy.ndarray, (nx, ny, nz), optional) – If provided then the probe density is scaled by this before being returned.
- Returns:
out – An array equal to probe(x)*scale. If ny=0 or nz=0 then array is of smaller dimension.
- Return type:
numpy.ndarray, (nx, ny, nz)
- FT(y, out=None)[source]#
If Y = sqrt(y[…,0]**2 + y[…,1]**2)*r then returns an indicator function on the disc Y < 1, y[2]==0. Again, if out is provided then computation is inplace. If y is a list of arrays then it is converted into a mesh first.
- Parameters:
y (numpy.ndarray, (nx, ny, nz, 3) or list of arrays of shape) – [(nx,), (ny,), (nz,)] Mesh of Fourier coordinates at which to evaluate the probe density. As a plotting utility, if a lower dimensional mesh is provided then the remaining coordinates are assumed to be 0 and so only the respective 1D/2D slice of the probe is returned.
out (numpy.ndarray, (nx, ny, nz), optional) – If provided then computation is performed inplace.
- Returns:
out – An array equal to FourierTransform(probe)(y). If ny=0 or nz=0 then array is of smaller dimension.
- Return type:
numpy.ndarray, (nx, ny, nz)
scattering_params#
Scattering Paramaters as Tabulated in “Advanced Computing in Electron Microscopy - Second Edition (2010) - Earl.J.Kirkland” ISBN 978-1-4419-6532-5 Pages 253-260 Appendix C
This transcription comes from scikit-ued (MIT license) - https://pypi.org/project/scikit-ued/
shape_factor_models#
- diffsims.utils.shape_factor_models.atanc(excitation_error, max_excitation_error, minima_number=5)[source]#
Intensity with arctan(s)/s profile that closely follows sin(s)/s but is smooth for s!=0.
- Parameters:
excitation_error (array-like or float) – The distance (reciprocal) from a reflection to the Ewald sphere
max_excitation_error (float) – The distance at which a reflection becomes extinct
minima_number (int) – The minima_number’th minima in the corresponding sinx/x lies at max_excitation_error from 0
- Returns:
intensity
- Return type:
array-like or float
- diffsims.utils.shape_factor_models.binary(excitation_error, max_excitation_error)[source]#
Returns a unit intensity for all reflections
- diffsims.utils.shape_factor_models.linear(excitation_error, max_excitation_error)[source]#
Returns an intensity linearly scaled with by the excitation error
- diffsims.utils.shape_factor_models.lorentzian(excitation_error, max_excitation_error)[source]#
Lorentzian intensity profile that should approximate the two-beam rocking curve. This is equation (6) in reference [1].
- Parameters:
- Returns:
intensity_factor – Vector representing the rel-rod factor for each reflection
- Return type:
array-like or float
References
[1] L. Palatinus, P. Brázda, M. Jelínek, J. Hrdá, G. Steciuk, M. Klementová, Specifics of the data processing of precession electron diffraction tomography data and their implementation in the program PETS2.0, Acta Crystallogr. Sect. B Struct. Sci. Cryst. Eng. Mater. 75 (2019) 512–522. doi:10.1107/S2052520619007534.
- diffsims.utils.shape_factor_models.lorentzian_precession(excitation_error, max_excitation_error, r_spot, precession_angle)[source]#
Intensity profile factor for a precessed beam assuming a Lorentzian intensity profile for the un-precessed beam. This is equation (10) in reference [1].
- Parameters:
excitation_error (array-like or float) – The distance (reciprocal) from a reflection to the Ewald sphere
max_excitation_error (float) – The distance at which a reflection becomes extinct
r_spot (array-like or float) – The distance (reciprocal) from each reflection to the origin
precession_angle (float) – The beam precession angle in radians; the angle the beam makes with the optical axis.
- Returns:
intensity_factor – Vector representing the rel-rod factor for each reflection
- Return type:
array-like or float
References
[1] L. Palatinus, P. Brázda, M. Jelínek, J. Hrdá, G. Steciuk, M. Klementová, Specifics of the data processing of precession electron diffraction tomography data and their implementation in the program PETS2.0, Acta Crystallogr. Sect. B Struct. Sci. Cryst. Eng. Mater. 75 (2019) 512–522. doi:10.1107/S2052520619007534.
- diffsims.utils.shape_factor_models.sin2c(excitation_error, max_excitation_error, minima_number=5)[source]#
Intensity with sin^2(s)/s^2 profile, after Howie-Whelan rel-rod
- Parameters:
- Returns:
intensity
- Return type:
array-like or float
sim_utils#
- diffsims.utils.sim_utils.acceleration_voltage_to_relativistic_mass(acceleration_voltage)[source]#
Get relativistic mass of electron as function of acceleration voltage.
- Parameters:
acceleration_voltage (float) – In Volt
- Returns:
mr – Relativistic electron mass
- Return type:
Example
>>> import diffsims.utils.sim_utils as sim_utils >>> mr = sim_utils.acceleration_voltage_to_relativistic_mass(200000) # 200 kV
- diffsims.utils.sim_utils.acceleration_voltage_to_velocity(acceleration_voltage)[source]#
Get relativistic velocity of electron from acceleration voltage.
Example
>>> import diffsims.utils.sim_utils as sim_utils >>> v = sim_utils.acceleration_voltage_to_velocity(200000) # 200 kV >>> round(v) 208450035
- diffsims.utils.sim_utils.acceleration_voltage_to_wavelength(acceleration_voltage)[source]#
Get electron wavelength from the acceleration voltage.
- diffsims.utils.sim_utils.beta_to_bst(beam_deflection, acceleration_voltage)[source]#
Calculate Bs * t values from beam deflection (beta).
- Parameters:
beam_deflection (NumPy array) – In radians
acceleration_voltage (float) – In Volts
- Returns:
bst – In Tesla * meter
- Return type:
NumPy array
Examples
>>> import numpy as np >>> import diffsims.utils.sim_utils as sim_utils >>> data = np.random.random((100, 100)) # In radians >>> acceleration_voltage = 200000 # 200 kV (in Volt) >>> bst = sim_utils.beta_to_bst(data, 200000)
- diffsims.utils.sim_utils.bst_to_beta(bst, acceleration_voltage)[source]#
Calculate beam deflection (beta) values from Bs * t.
- Parameters:
bst (NumPy array) – Saturation induction Bs times thickness t of the sample. In Tesla*meter
acceleration_voltage (float) – In Volts
- Returns:
beta – Beam deflection in radians
- Return type:
NumPy array
Examples
>>> import numpy as np >>> import diffsims.utils.sim_utils as sim_utils >>> data = np.random.random((100, 100)) # In Tesla*meter >>> acceleration_voltage = 200000 # 200 kV (in Volt) >>> beta = sim_utils.bst_to_beta(data, acceleration_voltage)
- diffsims.utils.sim_utils.diffraction_scattering_angle(acceleration_voltage, lattice_size, miller_index)[source]#
Get electron scattering angle from a crystal lattice.
Returns the total scattering angle, as measured from the middle of the direct beam (0, 0, 0) to the given Miller index.
Miller index: h, k, l = miller_index Interplanar distance: d = a / (h**2 + k**2 + l**2)**0.5 Bragg’s law: theta = arcsin(electron_wavelength / (2 * d)) Total scattering angle (phi): phi = 2 * theta
- diffsims.utils.sim_utils.et_to_beta(et, acceleration_voltage)[source]#
Calculate beam deflection (beta) values from E * t.
- Parameters:
et (NumPy array) – Electric field times thickness t of the sample.
acceleration_voltage (float) – In Volts
- Returns:
beta – Beam deflection in radians
- Return type:
NumPy array
Examples
>>> import numpy as np >>> import diffsims.utils.sim_utils as sim_utils >>> data = np.random.random((100, 100)) >>> acceleration_voltage = 200000 # 200 kV (in Volt) >>> beta = sim_utils.et_to_beta(data, acceleration_voltage)
- diffsims.utils.sim_utils.get_atomic_scattering_factors(g_hkl_sq, coeffs, scattering_params)[source]#
Calculate atomic scattering factors for n atoms.
- Parameters:
g_hkl_sq (numpy.ndarray) – One-dimensional array of g-vector lengths squared.
coeffs (numpy.ndarray) – Three-dimensional array [n, 5, 2] of coefficients corresponding to the n atoms.
scattering_params (str) – Type of scattering factor calculation to use. One of ‘lobato’, ‘xtables’.
- Returns:
scattering_factors – The calculated atomic scattering parameters.
- Return type:
- diffsims.utils.sim_utils.get_electron_wavelength(accelerating_voltage)[source]#
Calculates the (relativistic) electron wavelength in Angstroms for a given accelerating voltage in kV.
- diffsims.utils.sim_utils.get_holz_angle(electron_wavelength, lattice_parameter)[source]#
Converts electron wavelength and lattice paramater to holz angle :param electron_wavelength: In nanometers :type electron_wavelength: scalar :param lattice_parameter: In nanometers :type lattice_parameter: scalar
- Returns:
scattering_angle – Scattering angle in radians
- Return type:
scalar
Examples
>>> import diffsims.utils.sim_utils as sim_utils >>> lattice_size = 0.3905 # STO-(001) in nm >>> wavelength = 2.51/1000 # Electron wavelength for 200 kV >>> angle = sim_utils.get_holz_angle(wavelength, lattice_size)
- diffsims.utils.sim_utils.get_intensities_params(reciprocal_lattice, reciprocal_radius)[source]#
Calculates the variables needed for get_kinematical_intensities
- Parameters:
reciprocal_lattice (diffpy.Structure.Lattice) – The reciprocal crystal lattice for the structure of interest.
reciprocal_radius (float) – The radius of the sphere in reciprocal space (units of reciprocal Angstroms) within which reciprocal lattice points are returned.
- Returns:
unique_hkls (array-like) – The unique plane families which lie in the given reciprocal sphere.
multiplicites (array-like) – The multiplicites of the given unqiue planes in the sphere.
g_hkls (list) – The g vector length of the given hkl in the sphere.
- diffsims.utils.sim_utils.get_interaction_constant(accelerating_voltage)[source]#
Calculates the interaction constant, sigma, for a given accelerating voltage.
- diffsims.utils.sim_utils.get_kinematical_intensities(structure, g_indices, g_hkls_array, debye_waller_factors=None, scattering_params='lobato', prefactor=1)[source]#
Calculates peak intensities.
The peak intensity is a combination of the structure factor for a given peak and the position the Ewald sphere intersects the relrod. In this implementation, the intensity scales linearly with proximity.
- Parameters:
structure (diffpy.structure.Structure) – The structure for which to derive the structure factors.
g_indices (numpy.ndarray) – Indicies of spots to be considered.
g_hkls_array (numpy.ndarray) – Coordinates of spots to be considered.
debye_waller_factors (dict) – Maps element names to their temperature-dependent Debye-Waller factors.
scattering_params (str) – “lobato”, “xtables” or None
prefactor (float or numpy.ndarray) – Multiplciation factor for structure factor.
- Returns:
peak_intensities – The intensities of the peaks.
- Return type:
- diffsims.utils.sim_utils.get_points_in_sphere(reciprocal_lattice, reciprocal_radius)[source]#
Finds all reciprocal lattice points inside a given reciprocal sphere. Utilised within the DiffractionGenerator.
- Parameters:
reciprocal_lattice (diffpy.Structure.Lattice) – The reciprocal crystal lattice for the structure of interest.
reciprocal_radius (float) – The radius of the sphere in reciprocal space (units of reciprocal Angstroms) within which reciprocal lattice points are returned.
- Returns:
spot_indices (numpy.array) – Miller indices of reciprocal lattice points in sphere.
cartesian_coordinates (numpy.array) – Cartesian coordinates of reciprocal lattice points in sphere.
spot_distances (numpy.array) – Distance of reciprocal lattice points in sphere from the origin.
- diffsims.utils.sim_utils.get_scattering_params_dict(scattering_params)[source]#
Get scattering parameter dictionary from name.
- Parameters:
scattering_params (string) – Name of scattering factors. One of ‘lobato’, ‘xtables’.
- Returns:
scattering_params_dict – Dictionary of scattering parameters mapping from element name.
- Return type:
- diffsims.utils.sim_utils.get_unique_families(hkls)[source]#
Returns unique families of Miller indices, which must be permutations of each other.
- diffsims.utils.sim_utils.get_vectorized_list_for_atomic_scattering_factors(structure, debye_waller_factors, scattering_params)[source]#
Create a flattened array of coeffs, fcoords and occus for vectorized computation of atomic scattering factors.
Note: The dimensions of the returned objects are not necessarily the same size as the number of atoms in the structure as each partially occupied specie occupies its own position in the flattened array.
- Parameters:
structure (diffpy.structure.Structure) – The atomic structure for which scattering factors are required.
debye_waller_factors (dist) – Debye-Waller factors for atoms in the structure.
scattering_params (string) – The type of scattering params to use. “lobato”, “xtables”, and None are supported.
- Returns:
coeffs (numpy.ndarray) – Coefficients of atomic scattering factor parameterization for each atom.
fcoords (numpy.ndarray) – Fractional coordinates of each atom in structure.
occus (numpy.ndarray) – Occupancy of each atomic site.
dwfactors (numpy.ndarray) – Debye-Waller factors for each atom in the structure.
- diffsims.utils.sim_utils.is_lattice_hexagonal(latt)[source]#
Determines if a diffpy lattice is hexagonal or trigonal. :param latt: The diffpy lattice object to be determined as hexagonal or not. :type latt: diffpy.Structure.lattice
- Returns:
is_true – True if hexagonal or trigonal.
- Return type:
- diffsims.utils.sim_utils.scattering_angle_to_lattice_parameter(electron_wavelength, angle)[source]#
Convert scattering angle data to lattice parameter sizes.
- Parameters:
electron_wavelength (float) – Wavelength of the electrons in the electron beam. In nm. For 200 kV electrons: 0.00251 (nm)
angle (NumPy array) – Scattering angle, in radians.
- Returns:
lattice_parameter – Lattice parameter, in nanometers
- Return type:
NumPy array
Examples
>>> import diffsims.utils.sim_utils as sim_utils >>> angle_list = [0.1, 0.1, 0.1, 0.1] # in radians >>> wavelength = 2.51/1000 # Electron wavelength for 200 kV >>> lattice_size = sim_utils.scattering_angle_to_lattice_parameter( ... wavelength, angle_list)
- diffsims.utils.sim_utils.simulate_kinematic_scattering(atomic_coordinates, element, accelerating_voltage, simulation_size=256, max_k=1.5, illumination='plane_wave', sigma=20, scattering_params='lobato')[source]#
Simulate electron scattering from arrangement of atoms comprising one elemental species.
- Parameters:
atomic_coordinates (array) – Array specifying atomic coordinates in structure.
element (string) – Element symbol, e.g. “C”.
accelerating_voltage (float) – Accelerating voltage in keV.
simulation_size (int) – Simulation size, n, specifies the n x n array size for the simulation calculation.
max_k (float) – Maximum scattering vector magnitude in reciprocal angstroms.
illumination (string) – Either ‘plane_wave’ or ‘gaussian_probe’ illumination
sigma (float) – Gaussian probe standard deviation, used when illumination == ‘gaussian_probe’
scattering_params (string) – Type of scattering factor calculation to use. One of ‘lobato’, ‘xtables’.
- Returns:
simulation – ElectronDiffraction simulation.
- Return type:
ElectronDiffraction
- diffsims.utils.sim_utils.tesla_to_am(data)[source]#
Convert data from Tesla to A/m
- Parameters:
data (NumPy array) – Data in Tesla
- Returns:
output_data – In A/m
- Return type:
NumPy array
Examples
>>> import numpy as np >>> import diffsims.utils.sim_utils as sim_utils >>> data_T = np.random.random((100, 100)) # In tesla >>> data_am = sim_utils.tesla_to_am(data_T)
vector_utils#
- diffsims.utils.vector_utils.get_angle_cartesian(a, b)[source]#
Compute the angle between two vectors in a cartesian coordinate system.
- Parameters:
a (array-like with 3 floats) – The two directions to compute the angle between.
b (array-like with 3 floats) – The two directions to compute the angle between.
- Returns:
angle – Angle between a and b in radians.
- Return type:
- diffsims.utils.vector_utils.get_angle_cartesian_vec(a, b)[source]#
Compute the angles between two lists of vectors in a cartesian coordinate system.
- Parameters:
a (np.array()) – The two lists of directions to compute the angle between in Nx3 float arrays.
b (np.array()) – The two lists of directions to compute the angle between in Nx3 float arrays.
- Returns:
angles – List of angles between a and b in radians.
- Return type:
np.array()
- diffsims.utils.vector_utils.vectorised_spherical_polars_to_cartesians(z)[source]#
Converts an array of spherical polars into an array of (x,y,z) = r(cos(psi)sin(theta),sin(psi)sin(theta),cos(theta))
- Parameters:
z (np.array) – With rows of r : the radius value, r = sqrt(x**2+y**2+z**2) psi : The azimuthal angle generally (0,2pi]) theta : The elevation angle generally (0,pi)
- Returns:
xyz – With rows of x,y,z
- Return type:
np.array