Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Fix simulations of rotated structure #60

Merged
merged 8 commits into from
Jan 30, 2020
9 changes: 9 additions & 0 deletions diffsims/generators/diffraction_generator.py
Original file line number Diff line number Diff line change
Expand Up @@ -22,6 +22,7 @@

import numpy as np
from math import pi
from transforms3d.euler import euler2mat

from diffsims.sims.diffraction_simulation import DiffractionSimulation
from diffsims.sims.diffraction_simulation import ProfileSimulation
Expand Down Expand Up @@ -83,6 +84,7 @@ def __init__(self,
def calculate_ed_data(self,
structure,
reciprocal_radius,
rotation=(0,0,0),
with_direct_beam=True):
"""Calculates the Electron Diffraction data for a structure.

Expand All @@ -95,6 +97,9 @@ def calculate_ed_data(self,
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.
Expand All @@ -118,6 +123,10 @@ def calculate_ed_data(self,
spot_indicies, cartesian_coordinates, spot_distances = get_points_in_sphere(
recip_latt, reciprocal_radius)

ai,aj,ak = np.deg2rad(rotation[0]),np.deg2rad(rotation[1]),np.deg2rad(rotation[2])
R = euler2mat(ai,aj,ak,axes='rzxz')
cartesian_coordinates = np.matmul(R,cartesian_coordinates.T).T

# Identify points intersecting the Ewald sphere within maximum
# excitation error and store the magnitude of their excitation error.
r_sphere = 1 / wavelength
Expand Down
4 changes: 1 addition & 3 deletions diffsims/generators/library_generator.py
Original file line number Diff line number Diff line change
Expand Up @@ -30,7 +30,6 @@
from diffsims.libraries.vector_library import DiffractionVectorLibrary

from diffsims.utils.sim_utils import get_points_in_sphere
from diffsims.utils.sim_utils import simulate_rotated_structure
from diffsims.utils.vector_utils import get_angle_cartesian_vec


Expand Down Expand Up @@ -101,8 +100,7 @@ def get_diffraction_library(self,
intensities = np.empty(num_orientations, dtype='object')
# Iterate through orientations of each phase.
for i, orientation in enumerate(tqdm(orientations, leave=False)):
matrix = euler2mat(*np.deg2rad(orientation), 'rzxz')
simulation = simulate_rotated_structure(diffractor, structure, matrix, reciprocal_radius, with_direct_beam)
simulation = diffractor.calculate_ed_data(structure,reciprocal_radius,rotation=orientation,with_direct_beam=with_direct_beam)

# Calibrate simulation
simulation.calibration = calibration
Expand Down
40 changes: 0 additions & 40 deletions diffsims/utils/sim_utils.py
Original file line number Diff line number Diff line change
Expand Up @@ -346,46 +346,6 @@ def simulate_kinematic_scattering(atomic_coordinates,

return intensity


def simulate_rotated_structure(diffraction_generator, structure, rotation_matrix, reciprocal_radius, with_direct_beam):
"""Calculate electron diffraction data for a structure after rotating it.

Parameters
----------
diffraction_generator : DiffractionGenerator
Diffraction generator used to simulate diffraction patterns
structure : diffpy.structure.Structure
Structure object to simulate
rotation_matrix : ndarray
3x3 matrix describing the base rotation to apply to the structure, applied on the left of the vector
reciprocal_radius : float
The maximum g-vector magnitude to be included in the simulations.
with_direct_beam : bool
Include the direct beam peak

Returns
-------
simulation : DiffractionSimulation
The simulation data generated from the given structure and rotation.
"""
# Convert left multiply (input) to right multiply (diffpy)
stdbase = structure.lattice.stdbase
stdbase_inverse = np.linalg.inv(stdbase)
rotation_matrix_diffpy = stdbase_inverse @ rotation_matrix @ stdbase

lattice_rotated = diffpy.structure.lattice.Lattice(
*structure.lattice.abcABG(),
baserot=rotation_matrix_diffpy)
# Don't change the original
structure_rotated = diffpy.structure.Structure(structure)
structure_rotated.placeInLattice(lattice_rotated)

return diffraction_generator.calculate_ed_data(
structure_rotated,
reciprocal_radius,
with_direct_beam)


def get_points_in_sphere(reciprocal_lattice, reciprocal_radius):
"""Finds all reciprocal lattice points inside a given reciprocal sphere.
Utilised within the DiffractionGenerator.
Expand Down