pyxem / orix Goto Github PK

MTEX' Vector3d has these properties (https://github.com/mtex-toolbox/mtex/blob/develop/geometry/%40vector3d/vector3d.m#L47), we should too.

I just found myself wanting the spherical coordinates of Vector3d, i.e. the radius r = sqrt(x2 + y2 + z**2), polar coordinate theta = atan2(y, x) and azimuthal coordinate phi/rho = arccos(z / r). The computations are simple, so I think these should be attributes of Vector3d computed on the fly.

Will make a PR shortly if there is agreement.

Adding a orix.crystal_map.Phase to a orix.crystal_map.PhaseList is currently done like so

>>> from orix.crystal_map import Phase, PhaseList
>>> pl = PhaseList(names=["a", "b"], space_group=[10, 20])
>>> pl
Id  Name  Space group  Point group  Proper point group       Color
 0     a         P2/m          2/m                 112    tab:blue
 1     b        C2221          222                 222  tab:orange
>>> pl["c"] = "m-3m"
>>> pl
Id  Name  Space group  Point group  Proper point group       Color
 0     a         P2/m          2/m                 112    tab:blue
 1     b        C2221          222                 222  tab:orange
 2     c         None         m-3m                 432   tab:green

I didn't think this that much through upon implementing the PhaseList.__setitem__() method, as now I would have to set the space group in another line...

I think removing this __setitem__() method and adding a PhaseList.append(Phase()) method is better.

Our grid generations are not uniform on SO3, this can be corrected by using the Harr measure, specifically applying a cos() term within the grid_rzxz functionality. Although this code is going to orix it seems correct to fix this particular bit now.

I'm going to release a patch version of the v0.2 series to force matplotlib to be <3.3. I've chosen this solution (rather than a backport) as:

a) I think the way we were using the (now deprecated) arguments did serve a specific purpose within the code
b) Supporting the older release cycle (a tiny bit) isn't the worst idea in the world but ...
c) as we want people moving forward we have a bit more flexibility in approach (eg - requiring constraints on other packages)
d) The demos will fail with v0.2.3 and matplotlib > 3.2, and we as developers can easily make the avoidable

I used function/method signature type hints in the initial commits in #66, but decided to leave them for another PR since we haven't decided whether to use them or not (#66 (review)).

Why I want to use them:

• The editor I use, PyCharm, picks up these signatures type hints. It then tells me when I'm trying to pass something of the "wrong" type or use a return argument in the wrong way. I find this quite useful when coding.
• I understand what a function takes and returns quicker with type hints, instead of looking through the docstring.

So, we don't gain that much by introducing them, but we gain more than we lose, in my opinion.

Just because it's not being used any more I think we should delete the gh-pages branch and I don't think this should cause any issues - does anyone see a potential issue I'm missing?

I want a simple HDF5 reader/writer for the CrystalMap class into v0.3 (#58), so we can start to put stuff from pyxem in it, and read the data back into pyxem and kikuchipy.

I like the h5ebsd format specified here: Jackson, M.A., Groeber, M.A., Uchic, M.D. et al. h5ebsd: an archival data format for electron back-scatter diffraction data sets. Integrating Materials 3, 44–55 (2014). https://doi.org/10.1186/2193-9772-3-4. I plan to follow this as closely as it suits our CrystalMap class. Will add manufacturer and version datasets in the top group. This way we can easily create readers into e.g. MTEX's EBSD.load of this format.

Are there any interest in having readers returning anything else than a CrystalMap? If not, we should consider having only orix.load() figuring out which format the file is in, and returning a CrystalMap using the proper reader (like HyperSpy's clever hs.load()).

I would like to have the possibility to use the Phase and PhaseList classes in kikuchipy, not always via a CrystalMap. I therefore plan to make the reader/writer modular, so that these two classes are written to more or less standard data groups and sets.

Does anyone have any other requirements for the reader/writer? I'm not familiar with everything pyxem would want to write, so if there are anything that does not fit into a CrystalMap, please let me know.

At present our algorithm is:

1. generate a full covering of SO3 (in euler angles, rzxz - ie. Bunge)
2. convert to axis-angle representation
3. cut out rotations bigger than the omega_max of the fundamental zone in question
4. cut out the rotations using the full RF routine

A known perk of using rzxz is that it allows you to restrict the domain of considered rotations as the final rotation is along a high symmetry axis. This reduces the number of rotations that need considering in (2) by 6 for hexagonal, 4 for cubic etc.

If profiling suggests this would be a big time saving it could be worth implementing.

tests should be packaged for release

Let's set up a bare API reference with Sphinx via Read the Docs (RtD) at https://orix.readthedocs.io, as is done by e.g. HyperSpy (https://github.com/hyperspy/hyperspy/tree/RELEASE_next_minor/doc) and kikuchipy (https://github.com/kikuchipy/kikuchipy/tree/master/doc).

I'm fine with setting this up after v0.4.0 is out the door.

This requires

1. readthedocs.yml config file in the top directory
2. a orix/doc/ directory with make files (Makefile and make.bat), .rst files (index.rst and reference.rst/api.rst at a minimum) and a _static/ directory with any images, favicons etc.
3. someone with admin access to this repo sets up a GitHub webhook to RtD (https://docs.readthedocs.io/en/latest/webhooks.html#github). (It seems I don't have that yet.)
4. now merges into specified branches/tags (like master and some older v0.3.x for example) will prompt RtD to rebuild the docs!

We will need to update docstring parameters, links to other packages (so-called intersphinx), examples etc. for everything to look nice and linked correctly. We can do this before pushing the RtD doc stuff by building the docs locally (https://docs.readthedocs.io/en/latest/intro/getting-started-with-sphinx.html#quick-start).

• I also think we should add a doc extra feature requirements in setup.py with sphinx, sphinx-rtd-theme (if we use that theme) and other sphinx related packages, that can be used as pip install --editable .[doc] to install doc related dependencies. This way people can build the docs locally themselves if they want to.

• I suggest the standard RtD Sphinx theme (https://docs.readthedocs.io), however there are many others shown here (https://sphinx-themes.org/) and elsewhere. I guess the pyxem website uses the sphinx-bootstrap theme.

• We can tell RtD to build a PDF of the API as well (even an epub file!).

• Note that since RtD will be hosting the docs, https://orix.readthedocs.io will have a small ad at the bottom of the side bar.

orix-0.4.0 will contain significant functionality that we will depend on in https://github.com/pyxem/diffsims and https://github.com/pyxem/pyxem . This will make backward compatibility a more pressing concern for us and we are therefore considering an orix-1.0.0 release. Before that, we would like to review how crystal geometry and crystallographic axes are handled and make all conventions (and options for conventions) clear in the documentation.

There will probably be some more minor releases before we get to 1.0.0 - in this issue, I would like to list all features/developments that we would like to incorporate in orix-1.0.0 so that we can then split their addition up into minor release milestones.

The list so far includes:

• Standardization & documentation of crystal axis conventions
• Inverse pole figure plotting options
• cubochoric orientation sampling
• equal angle orientation sampling

Two questions:

1. Would we prefer CircleCI + CodeCov to our current Travis + Coveralls - or do we not see much difference?

2. Do we think it's worth automating releases to avoid manual errors, make life easier, and/or enable nightly builds?

I think this is a good time to (1) formalize some guides for the import structure and (2) simplify the API.

I see there was an attempt to simplify the API via orix.objects and orix.symmetry, is this still the goal?

I think there are generally 2,5 (2 and 3 go together) options for an import structure:

1. Import everything

>>> import orix
>>> orix.vector.Vector3d.xvector()
Vector3d (1,)
[[1 0 0]]

2. Import modules/submodules/sub-packages

>>> from orix import vector
>>> vector.Vector3d.xvector()
Vector3d (1,)
[[1 0 0]]

3. Import functions and classes

>>> from orix.vector import Vector3d
>>> Vector3d.xvector()
Vector3d (1,)
[[1 0 0]]

I think we shouldn't make anything available from orix (except __version__ etc.), and force users to either import modules, or functions and classes from modules.

We should also make the API as intuitive as possible. This is an example of how it is now for the orix.quaternion module, which I think is unintuitive

>>> from orix.quaternion import Quaternion
>>> from orix.quaternion.orientation import Orientation
>>> from orix.quaternion.orientation_region import OrientationRegion
>>> from orix.quaternion.symmetry import Symmetry

It should be like this

>>> from orix.quaternion import Quaternion, Orientation, OrientationRegion, Symmetry

In general I'm a fan of these guiding principles:

• To use a module it must be imported explicitly
• Functions and classes can either be called from an imported module, or directly imported from the module
• A function or class in the public API should only be importable directly from a module

scikit-image does it like this: https://github.com/scikit-image/scikit-image/blob/master/skimage/__init__.py. I like their import structure a lot, and would like us to use that.

Practically speaking I think what is required is to

• move class definitions Vector3d, Quaternion, Scalar (and Object3d if that's to be part of the public API) away from init.py to their own files
• import what's supposed to be available from module/init.py into that file, and perhaps specify a all list to really show what is part of the public API and not, like scikit-image does

These changes would create this (public) API

>>> from orix.crystal_map import CrystalMap, CrystalMapProperties, Phase, PhaseList
>>> from orix.io import load, save, loadang, loadctf
>>> from orix.plot import CrystalMapPlot, AxAnglePlot, RodriguesPlot, RotationPlot
>>> from orix.quaternion import Quaternion, Rotation, Misorientation, Orientation, OrientationRegion, Symmetry
>>> from orix.sampling import get_sample_fundamental, get_sample_local, uniform_SO3_sample
>>> from orix.scalar import Scalar
>>> from orix.vector import Vector3d, AxAngle, Homochoric, NeoEuler, Rodrigues, SphericalRegion

There a few: from x import *

We should replace these with explicit imports for namespace clarity.

I'm afraid I forgot this use case in PR #67.

This fails:

>>> from orix.io import load, save
>>> cm = load("austenite_ferrite.ang")
>>> cm
Phase   Orientations       Name  Symmetry       Color
    1   5657 (48.4%)  austenite       432    tab:blue
    2   6043 (51.6%)    ferrite       432  tab:orange
Properties: iq, dp
Scan unit: um
>>> save("austenite.h5", cm["austenite"])
>>> cm_austenite = load("austenite.h5")
[...]
  File "/home/hakon/kode/orix/orix/crystal_map/crystal_map.py", line 295, in y
    return self._y[self.is_in_data]
IndexError: boolean index did not match indexed array along dimension 0; dimension is 5657 but corresponding boolean dimension is 11700

What happens is that the datasets crystal_map/header/nx (and ny, nz) and crystal_map/data/x (and y, z) are has the value 5657 and a size of 5657, respectively, while the dataset crystal_map/data/is_in_data has the full map size of 11700.

Changing this line in orix.io.plugins.orix_hdf5.crystalmap2dict():

    z, y, x = [0 if i is None else i for i in crystal_map._coordinates.values()]

to

    z, y, x = [
        0 if i is None else i for i in [crystal_map._z, crystal_map._y, crystal_map._x]
    ]

fixes the issue, which I will do, while adding a test or two, in a PR immediately.

On a side note, for the readers of this issue's information, masking/slicing/subindexing a crystal map does not remove the information from the object, it just "hides" it by setting the masked out pixels to False in CrystalMap.is_in_data. Thus, when writing to file, both masked and unmasked values are included. This is in contrast to MTEX, where selecting a subset of the pixels removes the unwanted pixels. I chose this design because it made working with the data simpler.

When computing with symmetry operations I sometimes want the rotation matrix instead of the quaternion. Should we add some Quaternion.to_om etc. as in the Rowenhorst et al. (2015) paper? If yes, how should the API be? Quaternion.qu2om() looks a bit strange, and Rotation.to_euler() is already in use... I prefer the above Quaternion.to_om(), Quaternion.to_ax() etc.

And, is there a reason for to_euler() to be a Rotation method and not a Quaternion method?

Specifically, currently the code works with "an orientation as a passive rotation from the crystal reference frame to the specimen reference" (see preprint) while many other adopt the opposite convention (ie -
"orientation matrix transforms the specimen reference frame into the crystal reference frame") - this should be a fairly simple change to implement.

Eventually it should point to the orix page of the pyxem-website, but we're not quite there yet.

After #50 (c5daae0), but perhaps also before, I do not get the results from the docstrings examples for Misorientation.set_symmetry() and Orientation.set_symmetry():

Docstring:

>>> from orix.quaternion.orientation import Misorientation, Orientation
>>> from orix.quaternion.symmetry import C4, C2
>>> data = np.array([[0.5, 0.5, 0.5, 0.5], [0, 1, 0, 0]])
>>> Misorientation(data).set_symmetry(C4, C2)
Misorientation (2,) 4, 2
[[-0.7071  0.     -0.7071  0.    ]
 [ 0.      0.7071 -0.7071  0.    ]]
>>> Orientation(data).set_symmetry(C4)
Orientation (2,) 4
[[-0.7071  0.     -0.7071  0.    ]
 [ 0.     -0.7071 -0.7071  0.    ]]

Actual:

>>> from orix.quaternion.orientation import Misorientation, Orientation
>>> from orix.quaternion.symmetry import C4, C2
>>> data = np.array([[0.5, 0.5, 0.5, 0.5], [0, 1, 0, 0]])
>>> Misorientation(data).set_symmetry(C4, C2)
Misorientation (2,) 4, 2
[[-0.7071  0.7071  0.      0.    ]
 [ 0.      1.      0.      0.    ]]
>>> Orientation(data).set_symmetry(C4)
Orientation (2,) 1, 4
[[-0.7071  0.     -0.7071  0.    ]
 [ 0.      1.      0.      0.    ]]

What is wrong, the docstring or the code?

We should also make a Orientation.__repr__(), since it is with reference to C1 and only has one associated symmetry.

I'm motivated to create a CrystallographicMap class with Rotation objects and quality metric/property maps to couple these to EBSD pattern analysis in KikuchiPy.

See start of discussion on this class' implementation here: #47 (comment).

It should have existing functionality available in pyXem (https://github.com/pyxem/pyxem/blob/master/pyxem/signals/crystallographic_map.py).

Disclaimer: I haven't used orix before, and don't know all the CrystallographicMap class' potential use cases in orix, or things that need to change (if any) in the package overall (or how it affects other packages depending on it).

Briefly, I would like something similar to the following

>>> import orix
>>> import matplotlib.pyplot as plt
>>> data = orix.io.load_h5ebsd('/my/file.h5')
>>> data.__class__
orix.CrystallographicMap
>>> data.orientations
# The Rotation object
>>> data.prop
# Listing all available property maps
>>> plt.imshow(data.prop['image_quality'])  # Gives me a nice 2D image quality map

A new property can easily be added with data.prop['my_new_property'].

A natural source of inspiration for this class is MTEX' EBSD class (https://github.com/mtex-toolbox/mtex/tree/develop/EBSDAnalysis/%40EBSD).

Since the respository has moved to pyxem/texpy documentation and install in travis are outdated

This can follow the same style as was used for pyxem and diffsims.

@pc494 any thoughts on how long it might take us to get coverage to 100% - it would be nice to be able to say that we hit that in 0.2.0

Hi all!

Are there any plans to create a class also containing quality metrics (CI, IQ etc.) in addition to orientations? I'm asking for mainly two reasons:

• I'm starting to create quality metrics maps (Krieger-Lassen's image quality and soon average dot product map akin to EMsoft's) in KikuchiPy, and want it to be easy to create navigation masks and to navigate patterns in these instead of a virtual image or just the mean intensity of each pattern. I had in mind a sort of Map class, where I can easily do s.plot(navigator=m.iq) e.g. However, it would be natural to couple these maps to orientation data, which would fit in orix...
• Also, I plan to write a reader for EMsoft's dictionary indexing result file (HDF5), both orientations and quality maps. Should I add a load_emsoft function along with the other load functions under io/__init__.py, or does anyone have any other plans or preferanses?

I see the use of orix as the natural step forward to couple EBSD patterns directly with orientation data (from indexing in other softwares, for now), quality metric maps etc. in KikuchiPy. I'm therefore also interested in just generally hearing your future plans for orix.

We should make a proper pypi and conda release of this package, which needs:

• Continuous integration updated
• Coveralls setup
• Tests packaged within orix
• make license GPL-3
• Github and pypi release
• conda release

Not sure if this is the place to ask, but I would like to get a quick visual on the rotations generated by the functions in rotation_list_generators, as shown in pyxem/diffsims#103 . Is there already some kind of IPF plot function in diffsims/orix, or is this being considered? Right now, I feel like I've created a very crappy workaround by creating a fictitious ebsd map from the orientations, which I then import into MTEX:

rot_list_cubic = diffsims.generators.rotation_list_generators.get_fundamental_zone_grid(225, 10)
euler = np.array(rot_list_cubic)
ln = euler.shape[0]
phase_id = np.ones((ln, 1))
score_xy = np.vstack([np.ones(ln), range(ln), np.ones(ln)]).T
results_array = np.hstack([phase_id, euler, score_xy])
np.savetxt("rotlist.txt", results_array, delimiter="\t", newline="\r\n")

Also, do you have (or would you consider opening) a Gitter channel for pyxem/orix/diffsims? I feel like I have questions that are more discussion type, not really feature requests or bug reports.

At present the main orix workflow is designed with certain (navigation space) assumed.

The orientation workflow (demo-2) makes no internal use of the place (in navigation space) from which rotations come.

The misorientation workflow (demo-3) use a one pixel shift to say between which two pixels a misorientation is formed. However the main calculation (used for clustering) is naive to navigation space.

A new, workflow we may want to consider is creating two aligned datasets, < x , y | zxz > (this notation is from hyperspy) from which we then measure the misorientation (error if one dataset is considered the ground truth). I believe that this differs from the above cases as the final array that is needed is much smaller, and so the code could be much quicker.

This issue is a place for people with thoughts to offer them (mainly @EirikOpheim)

This ties into my question in #90 of where the projection from Cubochoric coordinates to a Homochoric vector could be placed.

EMsoft (Fortran) places this projection in their Lambert module (https://github.com/EMsoft-org/EMsoft/blob/develop/Source/EMsoftLib/Lambert.f90#L952), which also includes a stereographic projection (https://github.com/EMsoft-org/EMsoft/blob/develop/Source/EMsoftLib/Lambert.f90#L1630). These projections are used to e.g. orientation sampling, and projecting part of an EBSD/ECP master pattern onto a detector. I mentioned my desire to implement such a projection in #62. I believe it could be implemented in orix in such a general way so that it can be useful for other use cases apart from this projection.

MTEX (Matlab) has several spherical projections (https://mtex-toolbox.github.io/SphericalProjections.html, https://github.com/mtex-toolbox/mtex/tree/develop/plotting/sphericalProjections) in their geometry/geometry_tools part. These are used to great effect by McAuliffe and Britton to perform spherical integration of intensities along specific Kikuchi bands in EBSD patterns (code: https://github.com/tmcaul/SphericalAngularDF/blob/master/SphericalAnalysis/detector.m, arXiv: https://arxiv.org/abs/2005.10581).

The stereographic projection in #62 is "tied to" Matplotlib. Would it be possible to add it to a projection module, and use that for plotting? Then the cubochoric to homochoric projection could be placed in that module as well.

@dnjohnstone suggested in a review comment #66 (comment) that we could start looking at the space group/symmetry stuff in diffpy.structure, e.g. the SpaceGroup class: https://github.com/diffpy/diffpy.structure/blob/master/src/diffpy/structure/spacegroupmod.py#L194. After briefly looking at the SpaceGroup class now when working on #67 to get symmetrically equivalent atom positions in a unit cell (via diffpy.structure.symmetryutilities.expandPosition), this is what an object itself contains:

>>> from diffpy.structure.spacegroups import GetSpaceGroup
>>> sg225 = GetSpaceGroup(225)
>>> print(sg225.__class__)
diffpy.structure.spacegroupmod.SpaceGroup
>>> print(sg225.__dict__)
{'number': 225,
 'num_sym_equiv': 192,
 'num_primitive_sym_equiv': 48,
 'short_name': 'Fm-3m',
 'point_group_name': 'PGm3barm',
 'crystal_system': 'CUBIC',
 'pdb_name': 'F 4/m -3 2/m',
 'symop_list': [<diffpy.structure.spacegroupmod.SymOp at 0x7f01817c7dd0>,
  <diffpy.structure.spacegroupmod.SymOp at 0x7f01817c7e10>,
  <diffpy.structure.spacegroupmod.SymOp at 0x7f01817c7e50>,
  <diffpy.structure.spacegroupmod.SymOp at 0x7f01817c7e90>,
  [...]
  <diffpy.structure.spacegroupmod.SymOp at 0x7f01817d1e50>]}

Definition of rotation matrices and rotation vectors for the space groups: https://github.com/diffpy/diffpy.structure/blob/master/src/diffpy/structure/spacegroupmod.py

Definition of the space groups: https://github.com/diffpy/diffpy.structure/blob/master/src/diffpy/structure/mmlibspacegroups.py

How does the orix.quaternion.symmetry.Symmetry class compare to the diffpy.structure.spacegroupmod.SymOp class? It seems like Symmetry has more functionality, @pc494? I haven't used the Symmetry class yet, so therefore I cannot say what we gain by starting using the SymOp and SpaceGroup classes.

Thought it would be best to discuss this here than in the issue #46.

EMsoft's dot product files contain all their Euler angles and quality metrics in an HDF5 file in the h5ebsd format introduced by Jackson et al., 2015. Both EDAX TSL, Bruker and Oxford can write angles/metrics to these formats as well. Briefly, the format requires a structure in which the top HDF group contains at least two string data sets, Manufacturer and Version, and one group, e.g. Scan 1 (EMsoft/TSL) or Scan 0 (Bruker). Within these scan groups, the formats differ somewhat between vendors, but we could first determine vendor from Manufacturer, and then read the scans appropriately.

Examples of these files:

• EDAX TSL h5ebsd files can be found within the data released alongside Jackson et al. 2019.
• A Bruker Nano h5ebsd file was released by Wilkinson and Collins on Zenodo
• The Oxford format: https://github.com/oinanoanalysis/h5oina

There is also the question of how to handle quality metric maps, or any maps in general. So far the load_ang and load_ctf returns Rotation objects. In the short term, I guess we can also return a dictionary with all the quality maps entered by name?

Should we create a load_h5ebsd function (located in io/__init__.py)? Or should we create separate load_emsoft, load_h5ebsd_edax, load_h5ebsd_bruker etc.?

Since I have a setup for reading patterns from the EDAX and Bruker's h5ebsd files in KikuchiPy (https://github.com/kikuchipy/kikuchipy/blob/master/kikuchipy/io/plugins/h5ebsd.py), I can take responsibility for at least starting on this (or these) functions. Will try to not add any additional dependencies apart from h5py.

As soon as we make a decision on one/multiple readers, I will make a work-in-progress PR here.

texpy and pytex are both already packages on pypi so we'll need to rename for packaging

The property equivalent (of Misorientation, but inhereted by Orientation) is broken for the non-grain exchange case. Inspecting the copy+paste snippet below it's clear that for Gl._tuples != Gr._tuples we will always get an UnboundLocal error (on "orientations")

    def equivalent(self):
        """Equivalent misorientations
        Returns
        -------
        Misorientation
        """
        Gl, Gr = self._symmetry
        if Gl._tuples == Gr._tuples:  # Grain exchange case
            orientations = Orientation.stack([self, ~self]).flatten()
        equivalent = Gr.outer(orientations.outer(Gl))
        return self.__class__(equivalent).flatten()

See EMsoft-org/EMsoft#120 (comment). Will fix this when I've gotten an answer there.

At present the two (poorly named) distance functions (distance_1 and distance_2) are somewhat limited by RAM/speed. Enhancing the speed (of the slow function) or refactoring (for memory) would be a very productive use of developers time.

I can get all symmetry operations for space group 229 with diffpy.structure:

>>> from diffpy.structure.spacegroups import GetSpaceGroup
>>> sg229 = GetSpaceGroup(229)
>>> list_sym_ops = sg229.symop_list[:sg229.num_primitive_sym_equiv]
>>> sym_ops = np.stack([op.R for op in list_sym_ops])
>>> print(sym_ops.shape)
(48, 3, 3)
>>> print(sym_ops[0])
[[1. 0. 0.]
 [0. 1. 0.]
 [0. 0. 1.]]

Is it possible to get the same list from orix.quaternion.symmetry.Symmetry?

Since orix-0.3.0 is now released it seems a good time to write a list of functionality on the road map for 0.4.0.

The following will be part of orix-0.4.0:

• Porting and improving the creation of orientation space grids from diffsims to orix #71
  Adding stereographic projection #62
• Implementation of to_matrix() etc. methods in Rotation #79
• Developer improvements (#88)
• Add space group to Phase class (#104)

Does anyone want anything else on the list? @hakonanes @pc494

This is a continuation from pyxem/pyxem#304.

As I understand it, we would rather implement rotation lists in texpy than in pyxem directly. Should this live under texpy/grid?

I have a "pyxem" implementation at https://github.com/shogas/pyxem/blob/rotation_list/pyxem/utils/sim_utils.py. I think the base structure is fine, but I will rewrite it to fit the texpy style if it should be part of texpy.

Do we want to add a dependency to https://github.com/diffpy/diffpy.structure? My current implementation actually only needs the transformation matrix from the lattice (diffpy.structure.lattice.Lattice.stdbase), and it might be good to keep texpy independent on a specific structure representation.

Because error checking skills always improve by an order of magnitude after pressing submit....

We should:

• update setup.py so license says GPL-3
• update setup.py more generally to include a "long_description" link to readme and classifiers
• decide on and implement branch structuring of main repo
• check for any other silly mistakes
• sort zenodo doi etc
• conda release

As mentioned pyxem/pyxem#596 we have 'enough' for a new minor version of orix, this issue serves to collate developer request for things that it would be good to have in by the release point - @hakonanes

From my perspective

• Clean up the deprecations from the 0.2 series
• Improve the documentation for OrientationRegion
• Address #57
• Address the review comment on #60
• Add reader/writer for CrystalMap #59
• check the "grain-exchange" (GE) case of the .equivalent method, add path for S1==S2 no GE

I'll add more things as they occur to me.

I believe there is a small bug in Misorientation.set_symmetry() for octahedral symmetries. The bug is that for Euler angles (-45,0,0) and (45,0,0) the distance is reported to be 90 degrees, while it should be 0 degrees. For (0,0,0) and (90,0,0) it works fine. Even for (-45.0001,0,0) and (44.9999,0,0), so it is just for a special case that this goes wrong.

There is also another issue that I am not quite sure about for this method. The Euler angles (10,10,0) and (-80,10,0) have a non-zero distance. I would expect them to have zero distance but their distance is the same as for (0,0,0) and (10,10,0). I am not sure if this is correct, as at least the diffraction patterns, as simulated in DiffractionLibraryGenerator.get_diffraction_library() are equal for the two for an FCC structure.

To reproduce, simply run:

ori = Orientation.from_euler([np.radians(np.asarray([45,0,0]))]).set_symmetry(O)
ori2 = Orientation.from_euler(np.radians(np.asarray([-45,0,0]))).set_symmetry(O)
misori_base = Misorientation(~ori * ori2)
misori = misori_base.set_symmetry(O, O)
np.rad2deg(misori_base.angle.data)

I am not quite sure how to resolve the issue, but I am guessing that both (45,0,0) and (-45,0,0) are defined as inside the symmetry region of the O symmetry.

I thought it might be helpful to create a list of other resources that we might learn something from, consider using, or just find interesting in relation to orix. Partly because I'm clearing out stuff and had a list that is better here...

• MTEX - https://mtex-toolbox.github.io/

• numpy-quaternion - https://github.com/moble/quaternion

• NIST texture - https://github.com/usnistgov/texture

See initial discussion and motivation in issue #87. I'll submit a WIP PR tomorrow.

Questions:

1. Should we allow initializing a Phase only by a point group (and not providing the space group), and setting Phase.space_group to None by default? Simplest solution is to require the space group, and let the point group always be a derived property.
2. After initialization, should we allow setting the point group, or should this always be derived from the space group? (Difficult if it is None.) If allowed, should we set the Phase.space_group to None if the current space group is not derived from the new point group the user wants to set?
3. What symmetry info should the __repr__() method of Phase, PhaseList and CrystalMap show? As it is now they show the point group name. Should we change it to show "symmetry" with both the space group name, and point group name in parantheses? Or just the space group?
4. Should the properties be point/space_group or point/spacegroup?

I'll only add the SpaceGroup object to the Phase class, i.e. not extract any of its properties for easier access than via phase.spacegroup.this_property. If we want to use the SpaceGroup.symop_list (with SymOp objects) internally, we can extract the symmetry elements (rotation matrices and translation vectors) like so

>>> from orix.quaternion.rotation import Rotation
>>> from orix.vector import Vector3d
>>> from diffpy.structure.spacegroups import sg35
>>> sg35.short_name
'Cmm2'
>>> Rotation.from_matrix([i.R for i in sg35.symop_list])
Rotation (8,)
[[1.     0.     0.     0.    ]
 [0.     0.     0.     1.    ]
 [0.5774 0.5774 0.     0.5774]
 [0.5774 0.     0.5774 0.5774]
 [1.     0.     0.     0.    ]
 [0.     0.     0.     1.    ]
 [0.5774 0.5774 0.     0.5774]
 [0.5774 0.     0.5774 0.5774]]
>>> Vector3d([i.t for i in sg100.symop_list])
Vector3d (8,)
[[0.  0.  0. ]
 [0.  0.  0. ]
 [0.  0.  0. ]
 [0.  0.  0. ]
 [0.5 0.5 0. ]
 [0.5 0.5 0. ]
 [0.5 0.5 0. ]
 [0.5 0.5 0. ]]

Not sure if we still need this branch @hakonanes and would like to delete if not, just having a tidy up...

Quote from #17:

The Euler angles (10,10,0) and (-80,10,0) have a non-zero distance. I would expect them to have zero distance but their distance is the same as for (0,0,0) and (10,10,0). I am not sure if this is correct, as at least the diffraction patterns, as simulated in DiffractionLibraryGenerator.get_diffraction_library() are equal for the two for an FCC structure.

End quote

Reproduction of the first result is again possible with:

from orix.quaternion.symmetry import O
from orix.quaternion.orientation import Misorientation, Orientation
ori = Orientation.from_euler([np.radians(np.asarray([10,10,0]))]).set_symmetry(O)
ori2 = Orientation.from_euler(np.radians(np.asarray([-80,10,0]))).set_symmetry(O)
misori_base = Misorientation(~ori * ori2)
misori = misori_base.set_symmetry(O, O)
np.rad2deg(misori_base.angle.data)

As per #58.

We should be able to pass lattice parameters (a, b, c, alpha, beta, gamma) and atoms (element, site coordinates, site occupation, Debye Waller factors) when defining a Phase object used either by itself or as part of the PhaseList object in a CrystalMap object. For this we should use the diffpy.structure.Structure class (https://www.diffpy.org/diffpy.structure/package.html#module-diffpy.structure).

From before a Phase object only has a name, color and crystal symmetry. I propose we put the name as the Structure title, keep the space group as the orix symmetry, and the color as it is. We'll make the structure available simply via a Phase.structure attribute? This way all interactions with the Structure object happens via the diffpy.structure.Structure API.

One question: does diffpy allow storing the symmetry along with the Structure object? I couldn't seem to find it. And does that make sense?

This will give the following syntax:

>>> cm
Phase   Orientations       Name  Symmetry       Color
    1   5657 (48.4%)  austenite       432    tab:blue
    2   6043 (51.6%)    ferrite       432  tab:orange
>>> cm.phases["austenite"].structure.lattice.a  # Lattice parameters
>>> cm.phases["austenite"].structure.atoms[0]  # Atoms

Are we happy with this?

This will force users to get to know the diffpy.structure API. Their API is clearly explained in the above link, but they have almost no examples. We should keep this in mind.

EMsoft provides the orientation similarity metric (OSM). From their wiki:

[...] a map for which the intensity in each pixel measures how many of the top N matches that pixel has in common with its nearest neighbors. For a well annealed material, neighbouring pixels have very similar top match lists, so the OSM value will be high; near grain boundaries, one or more of the neighbors will have a completely different top match list, so the OSM value will be substantially smaller. As a result, Orientation Similarity maps generally show a clear outline of the microstructural components (grains, phases, voids, etc.).

This should be implemented as a method in the CrystalMap class, I think.

EMsoft implementation: EBSDgetOrientationSimilarityMap (https://github.com/EMsoft-org/EMsoft/blob/f5e6c8a5cb186478f5d34597708d4962f245a794/Source/EMsoftHDFLib/commonmod.f90#L335).

I think it should be implemented as two tools, the second using the first:

1. Compare the best n matches for one map position
2. Compare the best n in one map position to the best n matches in the four nearest neighbours, or whatever number of neighbours you want

Pulled out of #17

Quote: EirikOpheim

I believe there is a small bug in Misorientation.set_symmetry() for octahedral symmetries. The bug is that for Euler angles (-45,0,0) and (45,0,0) the distance is reported to be 90 degrees, while it should be 0 degrees. For (0,0,0) and (90,0,0) it works fine. Even for (-45.0001,0,0) and (44.9999,0,0), so it is just for a special case that this goes wrong.

To reproduce, simply run:

ori = Orientation.from_euler([np.radians(np.asarray([45,0,0]))]).set_symmetry(O)
ori2 = Orientation.from_euler(np.radians(np.asarray([-45,0,0]))).set_symmetry(O)
misori_base = Misorientation(~ori * ori2)
misori = misori_base.set_symmetry(O, O)
np.rad2deg(misori_base.angle.data)

I am not quite sure how to resolve the issue, but I am guessing that both (45,0,0) and (-45,0,0) are defined as inside the symmetry region of the O symmetry.

End quote

This belies a general issue. Often (sometimes with some explanation) the literature will leave degenerate rotations in. For example when working with euler angles, inclusive 0,m * pi ranges are used. In theory, either 0 or m * pi would suffice, but in practice people tend not to worry about it. Similarly, axis angles tend include pi in the range, which leaves both [1,0,0] pi and [-1,0,0] pi, a degenerate pair in.

In orix workflows this problem should stop being a problem once we run the clustering step.

A cursory glance suggests this is due to matplotlib going to v3.3 (from v3.2)

This is really two related issues (some previous discussion in #90)

1. There are lots of ways to grid in SO3, it seems reasonable to add new ones in as kwarg specified routines with uniform_SO3_sample - specifically cubochoric seems like a decent idea.

2. Deciding exactly what resolution should mean. In general the problem arises as uniform sampling is equal volume, but most users will want an angle based report. My belief is that the density of points is always given by dividing N (number of points) by the "volume" dimensionless of the area covered (8*pi for all of SO3).

Currently we use "resolution" to mean:

The distance from a rotation to 4 of its neighbours. This is because our implementation creates 6 direct neighbours, four of which are a fixed distance apart (ie - independent of where within SO3 we are). The other two neighbours are equidistant but there distance will vary. This is so that:
dV = sin beta dalpha dbeta dgamma

remains constant (using alpha, beta and gamma, left to right as in https://iopscience.iop.org/article/10.1088/0965-0393/23/8/083501/pdf)

So assuming a small change in each angle, near the origin we have

dV = sin(beta) dalpha dbeta dgamma

and at some other beta = Beta

dV = sin(Beta) dBeta dalpha dgamma

Setting the two equal and dividing through we find:

dbeta = (sin(Beta)/sin(beta) dBeta

so this size (in radians) changes. Our currently implementation is such that the average value is the same as dalpha, but half the time the two "beta neighbours" will be closer than the characteristic distance and half the time they will be further away.

The alternative is to talk about the "characteristic disorientation", as discussed in https://iopscience.iop.org/article/10.1088/0965-0393/24/8/085013

This issue arises because we make use of multiple Python versions, pyxem should have solution soon enough.