Skip to content

Commit

Permalink
implemented review, added to changelog
Browse files Browse the repository at this point in the history
  • Loading branch information
@din14970
din14970 committed Nov 23, 2020
1 parent d24eb78 commit ea7cb65
Show file tree
Hide file tree
Showing 3 changed files with 59 additions and 45 deletions.
7 changes: 5 additions & 2 deletions CHANGELOG.md
View file
Original file line number Diff line number Diff line change
Expand Up @@ -6,11 +6,14 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0

## [Unreleased]
### Changed
- get_grid_beam_directions, now works based off of meshes
- `get_grid_beam_directions`, now works based off of meshes
- the arguments in the `DiffractionGenerator` constructor and the `DiffractionLibraryGenerator.get_diffraction_library` function have been shuffled so that the former captures arguments related to "the instrument/physics" while the latter captures arguments relevant to "the sample/material".

### Added
- API reference documentation via Read The Docs: https://diffsims.rtfd.io
- New module: "sphere_mesh_generators"
- New module: `sphere_mesh_generators`
- beam precession is now supported in simulating electron diffraction patterns
- more shape factor functions have been added
- This project now keeps a Changelog

### Removed
Expand Down
61 changes: 35 additions & 26 deletions diffsims/generators/diffraction_generator.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -51,40 +51,48 @@
}


def _z_sphere_precession(phi, r_spot, wavelength, theta):
def _z_sphere_precession(theta, r_spot, wavelength, phi):
"""
Returns the z-coordinate in reciprocal space of the Ewald sphere at
distance r from the origin along the x-axis
distance r_spot from the origin
Parameters
----------
phi : float
The angle the beam is currently precessed to in degrees.
If the beam were projected on the x-y plane, the angle
is the angle between this projection with the x-axis.
theta : float
The azimuthal angle in degrees
r_spot : float
The distance to the point on the x-axis where we calculate z in A^-1
The projected length of the reciprocal lattice vector onto the plane
perpendicular to the optical axis in A^-1
wavelength : float
The electron wavelength in A^-1
theta : float
The electron wavelength in A
phi : float
The precession angle in degrees (angle between beam and optical axis)
Returns
-------
z : float
The height of the ewald sphere at the point r in A^-1
Notes
-----
* The azimuthal angle is the angle the beam is currently precessed to.
It is the angle between the projection of the beam and the projection of
the relevant diffraction spot both onto the x-y plane.
* In the derivation of this formula we assume that we will always integrate
over a full precession circle, because we do not explicitly take into
consideration x-y coordinates of reflections.
"""
phi = np.deg2rad(phi)
r = 1 / wavelength
theta = np.deg2rad(theta)
r = 1 / wavelength
phi = np.deg2rad(phi)
return -np.sqrt(
r ** 2 * (1 - np.sin(theta) ** 2 * np.sin(phi) ** 2)
- (r_spot - r * np.sin(theta) * np.cos(phi)) ** 2
) + r * np.cos(theta)
r ** 2 * (1 - np.sin(phi) ** 2 * np.sin(theta) ** 2)
- (r_spot - r * np.sin(phi) * np.cos(theta)) ** 2
) + r * np.cos(phi)


def _shape_factor_precession(
z_spot, r_spot, wavelength, precession_angle, function, max_excitation, **kwargs
z_spot, r_spot, wavelength, phi, shape_function, max_excitation, **kwargs
):
"""
The rel-rod shape factors for reflections taking into account
Expand All @@ -98,34 +106,35 @@ def _shape_factor_precession(
An array representing the distance of spots from the z-axis in A^-1
wavelength : float
The electron wavelength in A
precession_angle : float
phi : float
The precession angle in degrees
function : callable
shape_function : callable
A function that describes the influence from the rel-rods. Should be
in the form func(excitation_error: np.ndarray, max_excitation: float,
**kwargs)
max_excitation : float
Maximum excitation angle to consider if precession_angle = 0. With
precession, it is rather a parameter to describe the "width" of the
rel-rods.
Parameter to describe the "extent" of the rel-rods.
Other parameters
----------------
** kwargs: passed directly to function
** kwargs: passed directly to shape_function
Notes
-----
We calculate excitation_error as z_spot - z_sphere so that it is
* We calculate excitation_error as z_spot - z_sphere so that it is
negative when the spot is outside the ewald sphere and positive when inside
conform W&C chapter 12, section 12.6
* We assume that the sample is a thin infinitely wide slab perpendicular
to the optical axis, so that the shape factor function only depends on the
distance from each spot to the Ewald sphere parallel to the optical axis.
"""
shf = np.zeros(z_spot.shape)
# loop over all spots
for i, (z_spot_i, r_spot_i) in enumerate(zip(z_spot, r_spot)):

def integrand(phi):
z_sph = _z_sphere_precession(phi, r_spot_i, wavelength, precession_angle)
return function(z_spot_i - z_sph, max_excitation, **kwargs)
z_sph = _z_sphere_precession(phi, r_spot_i, wavelength, phi)
return shape_function(z_spot_i - z_sph, max_excitation, **kwargs)

# average factor integrated over the full revolution of the beam
shf[i] = (1 / (360)) * quad(integrand, 0, 360)[0]
Expand Down Expand Up @@ -175,7 +184,7 @@ class DiffractionGenerator(object):
Angle about which the beam is precessed. Default is no precession.
approximate_precession : boolean
When using precession, whether to precisely calculate average
excitation errors and intensities or use an approximation.
excitation errors and intensities or use an approximation. See notes.
shape_factor_model : function or string
A function that takes excitation_error and
`max_excitation_error` (and potentially kwargs) and returns
Expand All @@ -190,7 +199,7 @@ class DiffractionGenerator(object):
* A full calculation is much slower and is not recommended for calculating
a diffraction library for precession diffraction patterns.
* When using precession and approximate_precession=True, the shape factor
model defaults to lorentzian; shape_factor_model is ignored. Only with
model defaults to Lorentzian; shape_factor_model is ignored. Only with
approximate_precession=False the custom shape_factor_model is used.
"""

Expand Down
36 changes: 19 additions & 17 deletions diffsims/utils/shape_factor_models.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -81,9 +81,9 @@ def sinc(excitation_error, max_excitation_error, minima_number=5):
num = np.sin(fac * excitation_error)
denom = fac * excitation_error
return np.nan_to_num(
np.abs(np.divide(num, denom, out=np.zeros_like(num), where=denom != 0)),
nan=1,
)
np.abs(np.divide(num, denom, out=np.zeros_like(num), where=denom != 0)),
nan=1,
)


def sin2c(excitation_error, max_excitation_error, minima_number=5):
Expand All @@ -105,7 +105,7 @@ def sin2c(excitation_error, max_excitation_error, minima_number=5):
-------
intensity : array-like or float
"""
return sinc(excitation_error, max_excitation_error, minima_number)**2
return sinc(excitation_error, max_excitation_error, minima_number) ** 2


def atanc(excitation_error, max_excitation_error, minima_number=5):
Expand All @@ -131,15 +131,15 @@ def atanc(excitation_error, max_excitation_error, minima_number=5):
"""
fac = np.pi * minima_number / np.abs(max_excitation_error)
return np.nan_to_num(
np.arctan(fac*excitation_error)/(fac*excitation_error),
nan=1,
)
np.arctan(fac * excitation_error) / (fac * excitation_error),
nan=1,
)


def lorentzian(excitation_error, max_excitation_error):
"""
Lorentzian intensity profile that should approximate
the two-beam rocking curve. After [1].
the two-beam rocking curve. This is equation (6) in reference [1].
Parameters
----------
Expand All @@ -160,16 +160,18 @@ def lorentzian(excitation_error, max_excitation_error):
"""
# in the paper, sigma = pi*thickness.
# We assume thickness = 1/max_exitation_error
sigma = np.pi/max_excitation_error
fac = sigma/(np.pi*(sigma**2*excitation_error**2+1))
sigma = np.pi / max_excitation_error
fac = sigma / (np.pi * (sigma ** 2 * excitation_error ** 2 + 1))
return fac


def lorentzian_precession(excitation_error, max_excitation_error,
r_spot, precession_angle):
def lorentzian_precession(
excitation_error, max_excitation_error, r_spot, precession_angle
):
"""
Intensity profile factor for a precessed beam assuming a Lorentzian
intensity profile for the un-precessed beam. After [1].
intensity profile for the un-precessed beam. This is equation (10) in
reference [1].
Parameters
----------
Expand All @@ -195,8 +197,8 @@ def lorentzian_precession(excitation_error, max_excitation_error,
----------
[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.
"""
sigma = np.pi/max_excitation_error
u = sigma**2*(r_spot**2*precession_angle**2 - excitation_error**2)+1
z = np.sqrt(u**2 + 4*sigma**2 + excitation_error**2)
fac = (sigma/np.pi)*np.sqrt(2*(u+z)/z**2)
sigma = np.pi / max_excitation_error
u = sigma ** 2 * (r_spot ** 2 * precession_angle ** 2 - excitation_error ** 2) + 1
z = np.sqrt(u ** 2 + 4 * sigma ** 2 + excitation_error ** 2)
fac = (sigma / np.pi) * np.sqrt(2 * (u + z) / z ** 2)
return fac

0 comments on commit ea7cb65

Please sign in to comment.