Reintroduce option for 3D datasets

[smason747]

I would like to open a discussion on the reintroduction of CrystalMap objects to handle 3D EBSD data. It looks like this base level of functionality was removed several years ago #392 due to incompatibility with plotting functions.

In talking with @argerlt, we considered one option for allowing 3D data without breaking current functionality (avoiding Matplotlib warnings/errors as mentioned in #392) could be to create an additional CrystalMap3D class that extends the standard CrystalMap class while adding in a z variable as a new property. From there, a user could choose to return 2D slices of the larger 3D dataset, returned as traditional CrystalMap objects, where normal functionality could be used.

For example, taking slices on each of the orthogonal axes:

xmap3d[100].plot('ci') xmap3d[:,100,:].plot('ci') xmap3d[:,:,100].plot('ci')

This approach would be greatly helpful to my workflow (often involving serial sectioned 3D-EBSD datasets) by allowing me to process orientation data, apply rotations, etc., all of which already works in ORIX, agnostic of a prop.z property. Currently, any orientation processing I do has to be done either in individual slices, or by passing 1D quat arrays in and out of ORIX to other tools.

I have written up a rough, proof-of-concept CrystalMap3D object extending CrystalMap adding only the z property on my fork, but wanted to see if there were other ideas on how to approach this problem before going any further.

Please let me know what you think. Thanks!

[CSSFrancis]

A couple of leading questions. What does a 3D (x, y, time) dataset look like and how is that different from a 3D (x, y, z) crystal map. What about a crystal map that is 2D but you happen to have mutiple phases at some point? What about a more flexible crystal map that has dimensions where you can have mutiple orientations for each position in x/y as well as some score associated.

How should this be plotted?

[smason747]

Hmm, I may not have described my thoughts on the issue re plotting well. As far as the distinction between temporal data and static data (either 3D in (x, y, time) or 4D in (x, y, z, time)), I think similar principles could hold, where when you want to be able to plot a (2D) figure, it can be "sliced" out from the larger bulk object, indexing along whichever parameter is of interest, and returned back to a standard CrystalMap object.

This brings the main concern to adding flexibility in the data object structure itself, allowing for a z (or a time, multi-phase, or multi-orientation) parameter, with the option to add increased plotting capabilities, like those mentioned in discussion #489, down the line after basic functionality is implemented.

[hakonanes]

What about a more flexible crystal map that has dimensions where you can have mutiple orientations for each position in x/y as well as some score associated.

This is already supported.

[hakonanes]

What does a 3D (x, y, time) dataset look like and how is that different from a 3D (x, y, z) crystal map.

I see them as different data structures. A 3D crystal map is useful for morphological analysis in 3D, while a temporally correlated dataset is not. A new class storing registered crystal maps with the same origin and step sizes but possibly varying shapes would be better suited.

What about a crystal map that is 2D but you happen to have mutiple phases at some point?

We may look to support this. As mentioned, we currently support multiple rotations and properties per point. We can introduce multiple phase IDs per point too.

[argerlt]

@CSSFrancis ,

A couple of leading questions. What does a 3D (x, y, time) dataset look like and how is that different from a 3D (x, y, z) crystal map.

I think, in this instance, it's the same, right? ideally some day you would want to be able to rotate in spatial coordinates, but that requires a lot of plotting support that doesn't exist yet, and should not be tackled in a single PR.

What about a crystal map that is 2D but you happen to have multiple phases at some point?

Since it's a single class with a single PhaseList, this should still be fine. There could be some complications when extracting a regular CrystalMap from CrystalMap3d, but that's a solvable problem.

Unless you mean multiple phases at a single voxel? in which case, that seems like a problem to tackle in CrystalMap.

What about a more flexible crystal map that has dimensions where you can have multiple orientations for each position in x/y as well as some score associated.

Seems like a solid use case for #576 in the Rotatoins class:)
But also, I don't think that's a 3d-specific feature, that's just a potentially helpful addition for adding to the base CrystalMap class. I wouldn't try to support it in 3d-only.

How should this be plotted

This is the crux of why i suggested to @smason747 that CrystalMap3d should be it's own seperate class. Many of the current plotting functions break if you go 3d, and a few don't really have 3d-equivalents. So, just use all the old functions for only 2D maps, and make a new, smaller set of 3d-specific plotting functions. If you want to plot a slice with a 2D plotter, you just export the slice first to a 2D CrystalMap first.

Saying this out loud, I'm now realizing the correct setup might be _CrystalMapBase, then CrystalMap and CrtystalMap3D` that inherit form it....

Also, I think trying too hard to get the best visualization tools in the first attempt might tank this effort. I would suggest we just do the easy first-step of a 2D fly-through plot with a z-vector slider, and leave anything more complicated involivng pyQT or pyVista or plotly for a future effort

@smason747

Hmm, I may not have described my thoughts on the issue re plotting well. As far as the distinction between temporal data and static data (either 3D in (x, y, time) or 4D in (x, y, z, time)), I think similar principles could hold, where when you want to be able to plot a (2D) figure, it can be "sliced" out from the larger bulk object, indexing along whichever parameter is of interest, and returned back to a standard CrystalMap object.

Would you actually want this? i think adding a time axis might be a step too much.
I guess if we ever add spatial rotation, it would be nice to have time separate from space, but in that case I think i'd just create a new CrystalMap3d for every discrete time step.

[smason747]

Would you actually want this? i think adding a time axis might be a step too much.

I think if time were to be added, it would likely need to be added to the base CrystalMap object (which would extend to xmap3d as well), since @CSSFrancis brought up the concept of an (x, y, time) object.

I would suggest we just do the easy first-step of a 2D fly-through plot with a z-vector slider
I like this suggestion. This is something that could be done with matplotlib.widgets and not need to add any further dependencies like PyVista or plotly. I can work on adding this to my fork this week.

@CSSFrancis and @hakonanes, did either of you have any thoughts on @argerlt's suggestion of a _CrystalMapBase object? I am less familiar with structuring data objects and would appreciate your feedback.

[hakonanes]

But also, I don't think that's a 3d-specific feature, that's just a potentially helpful addition for adding to the base CrystalMap class. I wouldn't try to support it in 3d-only.

As mentioned above, this is already supported.

[hakonanes]

Hi @smason747,

Thank you for bringing this up! I think supporting 3D crystal maps is a great idea. This is well within the purpose of the existing CrystalMap class, so the additions in your fork should be added to that class.

Now, we have to make sure that all existing functionality works as-is and as expected. So there will be many edge-cases we have to handle. E.g., we should start by adding a "layer" or "depth" parameter to the current 2D crystal map Matplotlib plotter, with a default of the first layer.

I suggest to start by adding your additions to the current crystal map class. You can then open a draft PR. We can iron out the edge-cases and new methods and attributes there. What do you say?

[smason747]

Sounds good! I'll get to work on moving functionality over to the current CrystalMap class and go from there. Thank you all for your feedback.

[argerlt]

This is well within the purpose of the existing CrystalMap class, so the additions in your fork should be added to that class.

Okay, I've done a bit of a 180 on this.

Short version, I agree with @hakonanes now. I'd suggest ignoring my idea of sub-classing. Instead, I'd start a PR for supporting 3D in CrystalMap , and add a widget-based plotting extension to CrystalMapPlot that triggers if z is not None.

Longer version:
After reading @hakonanes comment, I went through every method in CrystalMap, which I listed below so we can discuss if we want (though it might be easier if you just start a PR and we work though these as you get to them):

    __init__
    id(self) -> np.ndarray:
    size(self) -> int:
    shape(self) -> tuple:
    ndim(self) -> int:
    x(self) -> Union[None, np.ndarray]:
    y(self) -> Union[None, np.ndarray]:
    dx(self) -> float:
    dy(self) -> float:
    row(self) -> Union[None, np.ndarray]:
    col(self) -> Union[None, np.ndarray]:
    phase_id(self) -> np.ndarray:
    phases(self) -> PhaseList:
    phases_in_data(self) -> PhaseList:
    rotations(self) -> Rotation:
    rotations_per_point(self) -> int:
    rotations_shape(self) -> tuple:
    orientations(self) -> Orientation:
    is_indexed(self) -> np.ndarray:
    all_indexed(self) -> np.ndarray:
    prop(self) -> CrystalMapProperties:
*    _coordinates(self) -> dict:
*    _all_coordinates(self) -> dict:
*    _step_sizes(self) -> dict:
*    _coordinate_axes(self) -> dict:
*    __getattr__(self, item):
*    __setattr__(self, name, value):
*    __getitem__(self, key: Union[str, slice, tuple, int, np.ndarray]):
*    __repr__(self) -> str:
    deepcopy(self) -> "CrystalMap":
    empty(cls,shape: Union[None, int, tuple] = None, step_sizes: Union[None, int, tuple] = None,) -> "CrystalMap":
*    get_map_data(self, item: Union[str, np.ndarray], decimals: Optional[int] = None, fill_value: Union[int, float, None] = np.nan,) -> np.ndarray:
    plot(...) -> plt.Figure:
*    def _data_slices_from_coordinates(self, only_is_in_data: bool = True) -> tuple:
*    def _data_shape_from_coordinates(self, only_is_in_data: bool = True) -> tuple:

I think all of these make sense in 3D. The starred ones are going to require some if/then logic to handle 2D vs 3D, but that's very doable.

How should this be plotted?

This is the big issue with 3D crystal maps IMO, but there is a separate class just for plotting, CrystalMapPlot. As of right now, today, I think adding an extension to CrystalMapPlot for fly-through maps with widgets that triggers when z is not None as well as an option for the viewing face (xy, yz, or zx), would make for great viewing tool and would be a concise and complete first feature addition.

Some day, possibly years from now, I want orix to support a FeatureMap class as well as being able to plot non-gridded data (this has been discussed repeatedly, #48, #531, #214, #15, etc). I suggested @smason747 sub-class CrystalMap with those features in mind, but:

1. @smason747 's changes wouldn't make of those additions more difficult, and in fact might make them easier
2. The challenges related to those new features are primarily associated with plotting, so if/when we get around to adding those features, we should probably only subclass CrystalMapPlot anyway.
3. they don't exist yet, and we can deal with the implications of them if/when they do.

In the same vein, it would be nice to have 3D KAM maps, or orientation isopleths, or 3D selected grain maps, but those are all things that would fit nicely into a CrystalMapPlot3d class if/when it ever came up.

With all that in mind, I would also drop the slice2d function It was a shortsighted suggestion on my part, and requires a 3D Voronoi tessellation to do correctly. We can tackle that if/when it comes up in the future, as it is surprisingly non-trivial.

Apologies for the wall of text.

[hakonanes]

I agree with the separation of 2D and 3D crystal map plotting functionality in two different classes.

An alternative to PyVista is vedo. Something to consider.

[CSSFrancis]

A bit of a round about way of getting here but what about just using hyperspy for 3D plotting?

I realize that orix specifically got to this point by trying to remove pyxem, diffsims and hyperspy as a dependancy, but if we are already adding (probably an optional) dependancy for 3D plotting then hyperspy isn't really the same plotting tool it was 5 or even 2 years ago. The benefits are that there is a large knowledge base to lean on, with the major drawback that it isn't a true 3D viewer.

I do like vedo although I'm not sure about it's performance. I've also played around with pygfx which is a bit newer but pretty awesome when considering speed. Only draw back is that it requires a gpu.

It also will render in a jupyter notebook which is a +.

[argerlt]

@CSSFrancis, I'll preface this by saying, I'm taking your workshop on Sunday specifically because I don't know enough about hyperspy, and it keeps showing up as a great drop-in solution to a lot of my problems.

That said, doesn't hyperspy work under the assumption all it's data is gridded? Ie, it's slicing through hypercubes? This would be fine for some data, but locks you out of spatial affine transforms, hex grid data, Voronoi pixel interpolation, some types of grain plotting, 3d contour maps, etc.

Which, again, is a bunch of features we don't currently support, but would like to.

Being able to export a CryatalMap to hyperspy more easily would sure be nice though, especially if we had an example showing off some of the advantageous reasons why you might do so for 3D data....

[CSSFrancis]

@argerlt I was hoping you were going to M&M, we can always dicuss things a bit more there. Hyperspy defines axes which can be:

• Uniform: Bascially what you would expect
• Non-Uniform: These have arbitary spacing much like how orix currently operates.
• Functional: Some axis defined by a function.

You're correct that all of the underlying data is gridded. I've thought a bit about what it looks like if that isn't the case.

Well a lot:

hyperspy/hyperspy#3440
hyperspy/hyperspy#3055

I don't know if hyperspy is the right answer here, it's probably not. The main benefit of hyperspy is that if you need some additional plotting tool it's pretty easy to get it added to hyperspy and harder to add it to a different plotting package.

[argerlt]

well, a lot.

This really got me good.

Yes, we should talk in person next week if there's time, i'd like that.