From da45a3270aac8335717d143ed2daaee8d78b0c1d Mon Sep 17 00:00:00 2001 From: Carter Francis Date: Wed, 8 May 2024 16:16:47 -0500 Subject: [PATCH] Refactor: Remove changes from RLV --- .../crystallography/_diffracting_vector.py | 44 ++++++++++++++++++- .../reciprocal_lattice_vector.py | 39 ++-------------- diffsims/generators/simulation_generator.py | 5 +-- .../test_diffracting_vector.py | 14 ++++++ .../test_reciprocal_lattice_vector.py | 14 ------ 5 files changed, 62 insertions(+), 54 deletions(-) diff --git a/diffsims/crystallography/_diffracting_vector.py b/diffsims/crystallography/_diffracting_vector.py index 4ede0284..fb2a1b13 100644 --- a/diffsims/crystallography/_diffracting_vector.py +++ b/diffsims/crystallography/_diffracting_vector.py @@ -19,6 +19,7 @@ from diffsims.crystallography import ReciprocalLatticeVector import numpy as np from orix.vector.miller import _transform_space +from orix.quaternion import Rotation class DiffractingVector(ReciprocalLatticeVector): @@ -90,7 +91,18 @@ def __init__(self, phase, xyz=None, hkl=None, hkil=None, intensity=None): self._intensity = np.array(intensity) def __getitem__(self, key): - dv_new = super().__getitem__(key) + new_data = self.data[key] + dv_new = self.__class__(self.phase, xyz=new_data) + + if np.isnan(self.structure_factor).all(): + dv_new._structure_factor = np.full(dv_new.shape, np.nan, dtype="complex128") + + else: + dv_new._structure_factor = self.structure_factor[key] + if np.isnan(self.theta).all(): + dv_new._theta = np.full(dv_new.shape, np.nan) + else: + dv_new._theta = self.theta[key] if np.isnan(self.intensity).all(): dv_new._intensity = np.full(dv_new.shape, np.nan) else: @@ -105,6 +117,36 @@ def __getitem__(self, key): return dv_new + @property + def basis_rotation(self): + """ + Returns the lattice basis rotation. + """ + return Rotation.from_matrix(self.phase.structure.lattice.baserot) + + def rotate_with_basis(self, rotation): + """Rotate both vectors and the basis with a given `Rotation`. + This differs from simply multiplying with a `Rotation`, + as that would NOT update the basis. + + Parameters + ---------- + rot : orix.quaternion.Rotation + A rotation to apply to vectors and the basis. + """ + + if rotation.size != 1: + raise ValueError("Rotation must be a single rotation") + # rotate basis + new_phase = self.phase.deepcopy() + br = new_phase.structure.lattice.baserot + # In case the base rotation is set already + new_br = br @ rotation.to_matrix().squeeze() + new_phase.structure.lattice.setLatPar(baserot=new_br) + # rotate vectors + vecs = ~rotation * self.to_miller() + return ReciprocalLatticeVector(new_phase, xyz=vecs.data) + @property def intensity(self): return self._intensity diff --git a/diffsims/crystallography/reciprocal_lattice_vector.py b/diffsims/crystallography/reciprocal_lattice_vector.py index e0492de2..12f190f2 100644 --- a/diffsims/crystallography/reciprocal_lattice_vector.py +++ b/diffsims/crystallography/reciprocal_lattice_vector.py @@ -24,7 +24,6 @@ import numba as nb import numpy as np from orix.vector import Miller, Vector3d -from orix.quaternion import Rotation from orix.vector.miller import ( _check_hkil, _get_highest_hkl, @@ -78,8 +77,6 @@ class ReciprocalLatticeVector(Vector3d): Indices of reciprocal lattice vector(s), often preferred over ``hkl`` in trigonal and hexagonal lattices. Default is ``None``. This, ``xyz``, or ``hkl`` is required. - rotation : orix.quaternion.Rotation, optional - Rotation to apply to the vectors. Default is ``None``. Examples -------- @@ -126,8 +123,8 @@ def __init__(self, phase, xyz=None, hkl=None, hkil=None): self._structure_factor = np.full(self.shape, np.nan, dtype="complex128") def __getitem__(self, key): - new_data = self.data[key] - rlv_new = self.__class__(self.phase, xyz=new_data) + miller_new = self.to_miller().__getitem__(key) + rlv_new = self.from_miller(miller_new) if np.isnan(self.structure_factor).all(): rlv_new._structure_factor = np.full( @@ -152,36 +149,6 @@ def __repr__(self): phase_name = self.phase.name return f"{name} {shape}, {phase_name} ({symmetry})\n" f"{data}" - @property - def basis_rotation(self): - """ - Returns the lattice basis rotation. - """ - return Rotation.from_matrix(self.phase.structure.lattice.baserot) - - def rotate_with_basis(self, rotation): - """Rotate both vectors and the basis with a given `Rotation`. - This differs from simply multiplying with a `Rotation`, - as that would NOT update the basis. - - Parameters - ---------- - rot : orix.quaternion.Rotation - A rotation to apply to vectors and the basis. - """ - - if rotation.size != 1: - raise ValueError("Rotation must be a single rotation") - # rotate basis - new_phase = self.phase.deepcopy() - br = new_phase.structure.lattice.baserot - # In case the base rotation is set already - new_br = br @ rotation.to_matrix().squeeze() - new_phase.structure.lattice.setLatPar(baserot=new_br) - # rotate vectors - vecs = ~rotation * self.to_miller() - return ReciprocalLatticeVector(new_phase, xyz=vecs.data) - @property def hkl(self): """Miller indices. @@ -1171,7 +1138,7 @@ def from_min_dspacing(cls, phase, min_dspacing=0.7, include_zero_vector=False): new = cls(phase, hkl=hkl).unique() if include_zero_vector: new_data = np.vstack((new.hkl, np.zeros(3, dtype=int))) - new = ReciprocalLatticeVector(phase, hkl=new_data) + new = cls(phase, hkl=new_data) return new @classmethod diff --git a/diffsims/generators/simulation_generator.py b/diffsims/generators/simulation_generator.py index 956cdf99..0a701e3c 100644 --- a/diffsims/generators/simulation_generator.py +++ b/diffsims/generators/simulation_generator.py @@ -24,7 +24,6 @@ from orix.quaternion import Rotation from orix.crystal_map import Phase -from diffsims.crystallography import ReciprocalLatticeVector from diffsims.crystallography._diffracting_vector import DiffractingVector from diffsims.utils.shape_factor_models import ( linear, @@ -191,7 +190,7 @@ def calculate_diffraction2d( # Rotate using all the rotations in the list vectors = [] for p, rotate in zip(phase, rotation): - recip = ReciprocalLatticeVector.from_min_dspacing( + recip = DiffractingVector.from_min_dspacing( p, min_dspacing=1 / reciprocal_radius, include_zero_vector=with_direct_beam, @@ -335,7 +334,7 @@ def calculate_diffraction1d( def get_intersecting_reflections( self, - recip: ReciprocalLatticeVector, + recip: DiffractingVector, rot: np.ndarray, wavelength: float, max_excitation_error: float, diff --git a/diffsims/tests/crystallography/test_diffracting_vector.py b/diffsims/tests/crystallography/test_diffracting_vector.py index ef3f4d30..bf4ccdf6 100644 --- a/diffsims/tests/crystallography/test_diffracting_vector.py +++ b/diffsims/tests/crystallography/test_diffracting_vector.py @@ -65,3 +65,17 @@ def test_flat_polar(self, ferrite_phase): r, t = dv.to_flat_polar() assert np.allclose(r, [np.sqrt(2), np.sqrt(2), 0.70710678, 0.70710678]) assert np.allclose(t, [np.pi / 4, np.pi / 4, -np.pi / 4, -np.pi / 4]) + + def test_get_lattice_basis_rotation(self, ferrite_phase): + """Rotation matrix to align the lattice basis with the Cartesian + basis is correct. + """ + rlv = DiffractingVector(ferrite_phase, hkl=[[1, 1, 1], [2, 0, 0]]) + rot = rlv.basis_rotation + assert np.allclose(rot.to_matrix(), np.eye(3)) + + def test_rotation_with_basis_raises(self, ferrite_phase): + rlv = DiffractingVector(ferrite_phase, hkl=[[1, 1, 1], [2, 0, 0]]) + rot = Rotation.from_euler([[90, 90, 0], [90, 90, 1]], degrees=True) + with pytest.raises(ValueError): + rlv.rotate_with_basis(rotation=rot) diff --git a/diffsims/tests/crystallography/test_reciprocal_lattice_vector.py b/diffsims/tests/crystallography/test_reciprocal_lattice_vector.py index 5e226466..fe02cd07 100644 --- a/diffsims/tests/crystallography/test_reciprocal_lattice_vector.py +++ b/diffsims/tests/crystallography/test_reciprocal_lattice_vector.py @@ -155,20 +155,6 @@ def test_get_hkil(self, silicon_carbide_phase): assert np.allclose(rlv.i, [0, 0, 0, 0, 0, 0]) assert np.allclose(rlv.l, [3, 2, 1, -1, -2, -3]) - def test_get_lattice_basis_rotation(self, ferrite_phase): - """Rotation matrix to align the lattice basis with the Cartesian - basis is correct. - """ - rlv = ReciprocalLatticeVector(ferrite_phase, hkl=[[1, 1, 1], [2, 0, 0]]) - rot = rlv.basis_rotation - assert np.allclose(rot.to_matrix(), np.eye(3)) - - def test_rotation_with_basis_raises(self, ferrite_phase): - rlv = ReciprocalLatticeVector(ferrite_phase, hkl=[[1, 1, 1], [2, 0, 0]]) - rot = Rotation.from_euler([[90, 90, 0], [90, 90, 1]], degrees=True) - with pytest.raises(ValueError): - rlv.rotate_with_basis(rotation=rot) - def test_gspacing_dspacing_scattering_parameter(self, ferrite_phase): """Length and scattering parameter properties give desired values.