DRAFT: Diffraction model with intragranular misorientation

[MACarlsen]

Hi All,

I was reading a very nice paper-draft this wednesday, which was using this package to simulate diffraction from deformed microstructures. But I couldn't help to thing that there must be a better way to make such a simulation. So I tried to implement a concept where the sample consists of a collection of gaussian-shaped "grains" with a narrow anisotropic gaussian ODF each.

If you make a the right approximations, the diffraction peaks produced by such a grain is also a "gaussian" in the detector coordinates, so you can use some concept from "gaussian splatting" to do fast evaluation. Since this package already has a few different model for rendering diffraction peaks, I figured why not add another one.

There is still a lot of work to be done testing and optimizing, but if you are interested and can give me some pointers about how to implement it, I would try to do it.

The model is also fully differentiable, so one can dream about doing proper gaussian-splatting fitting...

In the jupyternotebook from the front-page of the documentation you can run this:

from xrd_simulator.gaussian_crystal_model import GaussianGrainish, GaussianPolycrystal
import torch

def make_random_tensor(axis_1, axis_2):
    random_direction = np.random.normal(size=3)
    random_direction = random_direction/np.linalg.norm(random_direction)
    tensor = axis_1**2 * np.eye(3) + (axis_2**2-axis_1**2) * np.outer(random_direction, random_direction)
    return tensor

N = mesh.number_of_elements
grain_list = []
for ii in range(N):

    shape_tensor = torch.eye(3) * mesh.eradius[ii]**2

    misorientation_tensor = make_random_tensor(
        np.random.uniform(0.002, 0.0005),
        np.random.uniform(0.002, 0.0005),
    )

    grain = GaussianGrainish(
        phase=quartz, #  For now it assumes all gaussians are the same phase, but it just needs a wrapper for multiphase
        position=mesh.espherecentroids[ii], # 3 vector centroid real-space position
        shape_tensor=shape_tensor, # 3-by-3 symmetric shape tensor where the eigenvalues are the radii-squared.
        orientation=orientation[ii], # 3-by-3 rotation matrix.
        misorientation_tensor=misorientation_tensor, # 3-by-3 misorientation tensor where the eigenvalues are the misorientaion spread in radians squared.
                                                     # misorientation vectors live in laboratory coordinates.
        strain_tensor=polycrystal.strain_lab[ii] # 3-by-3 symmetric strain tensor.
    )

    grain_list.append(grain)

gauss_polycrystal = GaussianPolycrystal(grain_list)


angle = np.radians(1.0)
axis = np.array([0, 1 / np.sqrt(2), -1 / np.sqrt(2)])

f = gauss_polycrystal.render_detector_frame(
    detector=detector,    
    xray_propagation_direction=np.array([1.0, 0.0, 0.0]),
    wavelength=0.28523,
    sample_orientation=R.from_rotvec(0.5*axis * angle).as_matrix(),
    sample_rotation_during_exposure = 0.5 * axis * angle,
)

fig, ax = plt.subplots(1, 1, figsize=(12, 12))
ax.imshow(f, cmap="gray", vmax=100000)
plt.show()

max_misorientation=np.radians(1.0),
sample_rotation_during_exposure = np.array([0, 1 / np.sqrt(2), -1 / np.sqrt(2)])*np.radians(1.0), # This is of course with a gaussian time-window

[AxelHenningsson]

This looks amazing! ❤️

is the idea to allow each tet in a mesh to be associated to such a Guassian ODF expansion? I feel like this could be overlaid the classical rendering. I.e each splat is simply scaled with origin tet - beam intersection volume.

Exciting! Let us know if we can help!

Cheers
Axel

[MACarlsen]

Hi Axel,

I have two suggestions. I'd like to implement number one first and then see if there is time for number two, but I understand if you don't want to pollute the codebase with two separate workflows.

1. Full gaussian: As it is right now, I model also the real-space shape of the grains as gaussians. In that setting the peaks are always gaussians on the detector, even for big grains, but of course the shapes are very wrong in the near-field low misorientation setting (DCT) unless you fit the grains by many gaussians. This could then be a parallel implementation of the renderer which doesn't interfere with the current version, but can be used when smearing effects are important. This is almost already implemented.

2. Replacement for _find_solutions_to_tangens_half_angle_equation Instead of looking for exact solutions of the Laue equation during a finite rotation interval, we can calculate based on a spread-out Bragg condition. Even for grains with very little misorientation this has some advantages: a) No divergence of the Lorenz-factor near the rotation axis, b) works without continuous rotation, c) easier to add polychromaticity/incident beam angular divergence effects.

The second solution of course requires a lot more work, but I think it's manageable, since you are anyways already convolving the projected tet with the (gaussian) detector point spread, the propagated beam divergence could just be added to this point-spread kernel.

The second solution requires a lot more work and modifying a lot of the existing classes (and if I'm not home with a flu, I only have a couple of hours a week to play with this most likely, so the timescale would be a year or so :-) )

One issue I see is that I can only do continuous rotations by approximating it with a Gaussian time window, or by numerical integration. (a top hat is pretty well fit by 3 gaussians for most applications I guess)

If you want to understand what I'm doing, you really only need to look at
lauel._get_diffraction_arcsegment()

[AxelHenningsson]

I am happy for both of these ideas to go into xrd_simulator!

I am not opposed to multi-workflows, as long as it is clear in the docs how things work. It also makes sense to me to put this bigger goal on a two stage development track.

Thanks for looking into this Mads, we are excited to see what comes out of this.

Best,
Axel