get_rotation_from_z_to_direction + calculate_ed_data don't work together for monoclinic (and likely most non-orthonormal) unit cells

[JoonatanL]

I'm working on a monoclinic cocrystal. Running get_rotation_from_z_to_direction in zap_map_generators does not return the correct structure when u (in [uvw]) is nonzero. For example, for 100:

when plotted with calculate_ed_data as follows:

gives:

whereas

from CrystalMaker is correct. This is not present for rotations where u=0, such as 010:

in diffsims versus:

in CrystalMaker.

I wrote my own check for the Euler angles using scipy.Rotation.align_vectors, which gives a different orientation that seems equivalent but is not for the plots on calculate_ed_data.

For 100:
diffsims: 90, 99.526, -90 versus mine (which gives the correct diffraction): 90, 99.526,90.
For 010:
diffsims: 180, 90, -180 versus mine 0, 90, 0. (results are identical as we are rotating perpendicular to the monoclinic axis).

Working on this as I discover more of the problem.

[JoonatanL]

Just had a look, and I don't think align_vectors solves the problem, since 110 is still not right after, so wouldn't recommend the fix. Will delve into it a bit deeper.

edit: I wonder if it might be an issue with how the rotations are applied to the monoclinic cell?

[JoonatanL]

I've had a deeper look, and actually think the error comes in calculate_ed_data rather than in get_rotation, as the rotated cartesian_coordinates seem correct, so problem might be elsewhere. Will update the name.

[pc494]

Just be careful with the plotting, I would be inclined to use the properties of the simulations (specifically the indexes of the vectors to check) - for example I think the top one from diffsims might be okay?! Will leave it with you for now...

[hakonanes]

Thank you for pointing this out @JoonatanL. I haven't touched that part of diffsims, so I cannot say if this is a bug or not. But, we have a new Miller class in orix, initialized from a Phase with a diffpy structure (atoms, lattice) and a Symmetry (point group operations); I think this class could be used to handle these crystallographic operations in diffsims. I get your result with:

>>> from diffpy.structure import Lattice, Structure
>>> import numpy as np
>>> from orix.crystal_map import Phase
>>> from orix.quaternion import Rotation
>>> from orix.vector import AxAngle, Miller, Vector3d
>>> monoclinic = Phase(
...     point_group="2/m",
...     structure=Structure(lattice=Lattice(1, 2, 3, 90, 85, 90))
... )
>>> monoclinic
<name: . space group: None. point group: 2/m. proper point group: 112. color: tab:blue>
>>> uvw = Miller(uvw=(0, 1, 0), phase=monoclinic)
>>> uvw
Miller (1,), point group 2/m, uvw
[[0. 1. 0.]]
>>> uvw.data  # Vector with respect to cartesian coordinate system stored and used internally
array([[0., 2., 0.]])
>>> vz = Vector3d.zvector()
>>> axis = vz.cross(uvw).data
>>> angle = vz.angle_with(uvw).data
>>> r_axangle = Rotation.from_neo_euler(AxAngle.from_axes_angles(axis, angle))
>>> r_axangle
Rotation (1,)
[[ 0.7071 -0.7071  0.      0.    ]]
>>> np.degrees(r_axangle.to_euler())
array([[-0., 90.,  0.]])

Note that using Vector3d.yvector() instead of Miller(uvw=(0, 1, 0)) should work here, since we normalize the vector anyway.

[hakonanes]

@JoonatanL, what are your lattice parameters?

[JoonatanL]

@JoonatanL, what are your lattice parameters?
Lattice(a=8.7337, b=15.3838, c=11.5271, alpha=90, beta=99.256, gamma=90)

[hakonanes]

Hm, OK, I can't seem to get the correct answer with the above approach either...

[JoonatanL]

I've just been looking at the rotation matrices within the calculate_ed_data class, and they're wrong as well:

Deceptively similar, but when rotating the reciprocal 100 vector, they lead to different results. I think it's that stage that is going wrong, I'll try your method as well, but I'm wondering if it might be something like the a* and b* pointing the wrong way (as that would matter for monoclinic etc but not cubic.

[JoonatanL]

I think I may have figured it out. get_rotation_from_z is not entirely wrong I think. If you multiply the first two columns of the rotation matrix in calculate_ed_data by -1 (i.e. I guess flip the axes), everything becomes right again (not only 100, but 110, 210, etc.).

Might be why this has never shown up before, since I'm guessing it's basically -a pointing in the a direction and -b pointing in b. That would only cause problems in non-orthonormal systems. Will still double check.

[hakonanes]

OK, so there might be a misalignment of crystal axes between transforms3d.euler.euler2mat and the reciprocal lattice retreived from diffpy.structure.Structure.lattice?

[JoonatanL]

Yeah I imagine so, I'm not sure which one should be changed though, since I imagine this might not be the only function it comes up in. Any thoughts?

[hakonanes]

Is the issue here perhaps related to pyxem/orix#236? The issue there was not getting the expected symmetrically equivalent Miller indices from (311), which was "solved" by switching the b and c axes when setting the lattice parameters.

[pc494]

Both, been thinking about sweeping this up for the next release. It looks like it's two things:

1. get_rotation_to_z_direction
   should be rewritten into orix syntax.

2. calculate_ed_data
   needs to be checked against a sensible case study, @JoonatanL is there any chance you could send me a CrystalMaker pattern at a given rzxz angles (ideally for 100) and the diffsims results, then I can did into it.

Also, could you both check that what I've said makes sense?

Ta.

[hakonanes]

Rewriting using orix sounds very good, don't have any other comments besides that.

[hakonanes]

I don't know the chosen convention in CrystalMaker for aligning the crystal's lattice vectors (a, b, c) with the cartesian reference frame (e1, e2, e3), but diffsims uses diffpy.structure by default, which has the alignment e1 || a*, e2 || e1 x e3, e3 || c. This might be the reason for the misunderstanding here.

See orix' documentation for the definition chosen in orix, which uses diffpy.structure but makes sure to switch the alignment internally to e1 || a, e2 || e1 x e3, e3 || c* so that this is used every time.