Laue three dimensional neutron diffraction

This article presents a measurement technique and data analysis tool to perform 3D grain distribution mapping and indexing of oligocrystalline samples using neutrons: Laue three-dimensional neutron diffraction (Laue3DND). The approach builds on forward modelling used for correlation and multiple fitting of the measured diffraction spots relative to individual grains. This enables not only to identify individual grains, but also their position and orientation in the sample. The feasibility and performance of the Laue3DND approach are tested using multi-grain synthetic datasets from cubic (α-Fe) and tetragonal (YBaCuFeO5) symmetries. Next, experimental results from two data sets measured at the FALCON instrument of Helmholtz-Zentrum Berlin are presented: A cylindrical alpha iron (α-Fe) reference sample with 5 mm diameter and 5 mm height, as well as a 2 mm3 layered perovskite (YBaCuFeO5). Using Laue3DND, we were able to retrieve the position and orientation of 97 out of 100 grains from a synthetic α-Fe data set, as well as 24 and 9 grains respectively from the α-Fe and YBaCuFeO5 sample measured at FALCON. Results from the synthetic tests also indicate that Laue3DND is capable of indexing 10 out of 10 grains for both symmetries in two extreme scenarios: using only 6 Laue projections and using 360 projections with extremely noisy data. The precision achieved in terms of spatial and orientation resolution for the current version of the method is 430 μm and 1° respectively. Based on these results obtained, we are confident to present a tool that expands the capabilities of standard Laue diffraction, providing the number, position, orientation and relative size of grains in oligocrystalline samples.


Introduction
Understanding the link between a material structure at different length scales and its emerging macroscopic properties is a general theme of material science. In the case of crystalline materials, retrieving 3D spatial information about the individual crystallites (grains) in the bulk non-destructively has been the motivation behind several methods developed with X-rays and widely applied for corresponding valuable studies since the millennium [1][2][3][4][5] . These methods make possible the retrieval of grain maps from millimetric samples with sub-micron resolution. Sub-micron resolutions are outside the reach of neutron grain mapping techniques, but neutrons bear the promise to retrieve grain maps with sub-millimeter resolution of larger sample volumes due to their better penetration characteristics for many structural engineering materials and dense crystals. Grain mapping of bulky engineering samples is particularly interesting for large directionally grown pieces, such as nickel-based turbine blades, and large samples undergoing phase transformations during use, such as iron-based shape memory alloys. In the field of solid-state physics, a 3D grain mapping characterization tool with neutrons enables the utilization of imperfect crystals for diffraction studies, (e.g.) by providing 3D information about the position of the largest grain within a sample so it can be cut out. In both metallurgy and solid-state physics, neutrons are better suited than X-rays for in-situ testing with bulkier sample environments. A first approach and proof-of-principle experiment with neutrons (nDCT) was reported recently 6 , using full sample illumination and a neutron beam with a narrow energy spectrum. The result was a grain boundary topological 3D reconstruction of 13 grains from an aluminum sample, measured with cold neutrons at the ICON beamline 7 of Paul Scherrer Institut (PSI). More recently, a substantially more efficient Time-of-Flight (ToF) approach utilizing a pulsed neutron source has been introduced 8 utilizing the SENJU beamline 9 at J-PARC (ToF3DND) in conjunction with a timepix transmission imaging detector 10 . ToF3DND enabled indexing and reconstruction of more than a hundred grains under full sample illumination. Both realizations of grain mapping with neutrons were based on partial or full wavelength resolved measurements enabling diffraction analyses and retrieval of 3D grain maps. Therefore, introducing a white beam Laue diffraction technique facilitates grain mapping at most neutron sources, in particular at continuous sources, where energy selection implies selecting only a part of available flux. Thus we present Laue3DND, which draws partially on intense continuous white beam flux but on the other hand on complex computational efforts in a forward modelling approach.
Firstly, the necessary crystallographic and geometrical concepts to build the forward model are explained, structured in the sample, laboratory and detector reference systems. Once the Laue patterns can be simulated, the solver is introduced and structured in seeding, single grain fitting and global fitting. These are the set of iterative algorithms implemented to find the best possible match between the simulated and measured spots. The indexing procedure is then finished, so the analysis of the diffracted intensities can be carried out. In order to test the performance of Laue3DND a series of tests are conducted using synthetic data sets, which provide benchmarks for robustness, precision and limitations of the method. Next, the experimental setup of the FALCON beamline of Helmholtz Zentrum Berlin (HZB) is detailed, followed by the experimental results from the α-Fe and YBaCuFeO 5 oligocrystalline samples. Finally, the results obtained, the current performance of the method and the future challenges and improvements are laid out in the discussion and conclusions.

Forward model
The forward model is the tool that allows simulating diffraction patterns from crystal and beamline parameters as they would be during an experiment. All the crystallographic and geometrical calculations are performed in the forward model as shown in Fig. 1, so later the solver can compare the position of the predicted spots (P) with the measured spots. The neutron beam traverses the center of the back scattering detector through hole in the back scattering detector and meets the beam stop in front of the forward scattering detector. Ideally, the sample is placed in the center of the rotation axis ω, as shown in the figure. b) Flowchart of the forward model. From a given grain orientationr, the valid scattering wavelengths λ are calculated, followed by the reciprocal scattering vectorḠ r and the direction of the diffracted vectork /|k | which isL in the sample reference system. The position of the grainX and the detector positionD and its orientationN are then used together withL to calculate the position of diffraction spotP. This process is repeated for every hkl plane and ω rotation step. The black labeled part of the flowchart (Seeding and fitting) is run every time the forward model is used, while the green part is only used during the first part of the algorithm (Seeding).
Ultimately, the forward model simulates the Laue diffraction pattern measured at one or more given detectors (D andN), for a given grain orientation (r) and a given position of the grain in space (X). In order to do this, three different reference systems are used throughout the forward modelling: the sample, laboratory and detector reference systems.

2/15
whereL samp hkl,ω is the unit vector with the direction of the diffracted beam in the SRS. A scheme of the sample reference system is shown in Fig. 2a.   Figure 2. a) Scheme of the sample reference system, wherek hkl,ω represents the direction of the incoming beam and k hkl,ω =Ḡ hkl,ω +k hkl,ω is the direction of the diffracted beam andḠ hkl,ω is a reciprocal lattice vector for the plane hkl at a given rotation angle ω. b) Scheme of the laboratory reference system. Given the detector position and orientation (D,N), the diffraction spotP hkl,ω can be calculated for every ω, as the intersection of the line formed by the scattering vectorL hkl,ω and the grain positionX ω with the detector plane.
First, the forward model requires a grain orientation defined in the sample reference system. We define the orientation of the grain asr, a three component vector in Rodrigues space 11 . Defining the orientation in Rodrigues formulation has two main benefits with respect to Euler angles: any orientation can be defined with only three components (which simplifies the optimization process) and the orientations are uniquely defined, given that ther vector lies in the fundamental zone of the given crystal symmetry. The Rodrigues vectorr is used to calculate the rotation matrix of the grain orientation by where U is the rotation matrix of the grain orientation and r 1 , r 2 and r 3 are the three components ofr. Once the orientation of the grain is defined, the next step is to calculate all the neutron wavelengths which can potentially fulfill the Bragg condition by where B is the transformation matrix between the Cartesian lattice and the reciprocal space,Ḡ hkl are the Miller indices of the reflecting hkl planes, Γ ω is the right-hand rotation matrix around the z-axis for the angle ω, U is the rotation matrix for the grain orientation, and λ hkl,ω is the resulting diffracted wavelength. With the formulation given in equation (3) only the first component of the product Γ ω UBḠ hkl is used for the calculation, as explained in Ref. 12. The next step in the forward model is to obtain all the reciprocal lattice vectors for the given structure, orientation and diffraction planes, which satisfy the Bragg condition. From the neutron wavelengths found in equation (3), we remove those which are not present in the incident neutron spectrum in order to find the valid reciprocal scattering vectors, bȳ whereḠ r hkl,ω are the reciprocal scattering vectors. Finally, one can calculate the direction of the diffracted beam for each reciprocal scattering vector bȳ

Laboratory reference system
The laboratory reference system (LRS) is used to calculate the intersection point (P hkl,ω ) with the detector of the line formed by the diffraction vector (L hkl,ω ) and the center-of-mass (CMS) of the grain (X ω ). That intersection point is the position of the diffraction spot on the detector in the LRS. We first transformL samp hkl,ω into de LRS bȳ whereL lab hkl,ω is the direction of the diffracted beam in the LRS, as shown in Fig. 2b. The position of the CMS of the grain (X ω ) is initially assumed to be at the origin of the coordinate system (O) as a first approximation. The position of the predicted diffraction spotP lab hkl,ω on the detector, given bȳ whereX ω is the CMS of the grain in the LRS and where t is the modulus required for the unit vectorL lab hkl,ω to reach the detector from X ω ,D is the detector position andN is the orientation.
The orientation of the detector plane is given bȳ where A f and A b are 3D right hand rotation matrices for the forward (9) and backward (10) detectors respectively, which account for the three tilt directions of the detectors. Following this procedure for each of the detectors the positionP lab hkl,ω of every diffraction spot is obtained for every hkl, ω and a given orientation.

Detector reference system
To calculate the difference between experiment and calculation the modelled diffraction spots need to be projected from the grain in the LRS onto the detector in the detector reference system (DRS) bȳ whereP det hkl,ω is the diffraction spot position in the DRS, after subtracting the position of the detectorD fromP lab hkl,ω and applying the inverse rotation matrix A −1 from (9) for the forward detector or (10) for the backward detector.
Once the diffraction spot position is calculated in the DRS (P det hkl,ω in Fig. 3a), a comparison with the positions of the experimental diffraction spots (peaks) can be performed. The blue lines describe the seeding process, while the red lines describe the fitting process. Both algorithms are iterated during the solver execution to find the optimal parameters for every grain fit. The green label indicates the starting values for the solver.

Solver
The solver (Fig. 3b) comprises of the set of algorithms dedicated to find the best possible fit between the experimental data and the grain-by-grain modelled solutions. The core of the solver is an algorithm which attempts to assign each calculated spot to the CMS of every segmented peak. We argue that a k-nearest-neighbors (KNN) algorithm is the most efficient assignment solution for this task, despite providing non-unique assignments for many spots 13 . Classic assignment algorithms, like the Hungarian or Munkres 14 algorithm, provide a unique assignment but the complexity of these is O(n 3 ), while the KNN has a complexity of O(nd) where n is the number of spots and d is the dimensionality of the position (in this case d=2). In the case of an experimental spot being assigned to several predicted spots (or vice versa) we discard all assignments except the one with minimum Euclidean distance. Once the KNN is finished, a cost function is computed by using the Euclidean distance of every assignment. Since the experimental data include diffraction spots from many different grains, the algorithm has to be robust enough to identify and optimize a correct orientation with a large number of outliers. Different strategies are followed in order to deal with incorrectly fitted spots, but the underlying principle is to give a higher weight to assignments with the lower Euclidean distance over those with larger Euclidean distances.

Seeding
Seeding is the overall search of the orientation space in order to find orientations which have higher probability of being close to a grain's orientation in the sample. As explained in Fig. 3b, the solver first generates a number m of diffraction patterns according to m corresponding divisions of the Rodrigues 11 orientation space. The value of m is a compromise between computation time and the size of the fundamental zone. For our experiments, the value of m ranges between 15000 and 45000, which corresponds to a maximum step-width in orientation space between of 3.8 • and 2.7 • (for a cubic system). Initially, the forward model assumes a grain position at the origin of the laboratory reference system, and the detectors' position assumptions are based on approximate distance measurements when configuring the experiment. The resulting diffraction patterns calculated on such detector planes are then compared to the experimental . The m calculated orientations are then sorted by the median Euclidean distance between the peaks and the predicted spots. The median of the distances proves to be useful in this case because the goal is to find best suited orientation candidates that have a large number of low-distance assignments without being affected by long-distance outliers, which would be the case if (e.g.) the mean distance was used. Once the m orientations have been sorted by their median distance, the best match is selected for further fitting.

5/15
Single grain fitting A constrained optimization algorithm is used initially to fit each one of the individual grains, in contrast to the seeding process where only direct comparison was performed to sort the orientations by the median. The variables are the detectors' positions and orientations (D,N) and the grain's position and orientation (X ω ,r) constituting 6 input parameters for the grain plus 6 parameters per detector used. Since we are using two detectors in our experiments every single grain fit involves a total of 18 variables. The constrains forD andN are defined by the user and relate to measurement errors. In the cases presented here the tolerances are 10 mm detector misalignment in every direction and 2 • for every tilt. With respect toX ω , a 3D space larger than the volume of the sample is set as boundary condition. The best orientation obtained during the seeding process is chosen (r i ) for the first iteration of the optimization. The volume of the voxel in Rodrigues space which serves as boundary condition for the optimization is defined by the user. In the case of a cubic system we use a volume of 0.066x0.066x0.066 in Rodrigues' space, which corresponds approximately to 7.55 • for every Euler angle. Once the boundary conditions for all the parameters are set and the KNN assigns the predicted spots to every peak, we need to estimate the goodness of the fit. Instead of using the median of the euclidean distances, the cost function used in this case is where j is the index of N total assignments and dist j is the euclidean distance of the j th assignment in millimeters. The value of 0.25 added to the distance corresponds to the thickness of the scintillator layer of the detector (0.25 mm), which is a good estimate of the resolution of the detector system. By adding 0.25 to the cost function we give a similar weight to all the assignments which have a distance smaller than the resolution of the detector system. Once the cost function is minimized through the constrained optimization algorithm, a criterion to segment the correctly fitted spots from the outliers is applied. Plotting a histogram of the distances resulting from all assignments produces a rightskewed histogram with a long tail, as shown in Fig. 4. This is the result of optimizing with the chosen cost function (equation (12)), which gives a higher weight to assignments with small distances and quickly lowers the weight to the assignment when distance increases. This gives the algorithm a higher incentive on reducing the distance of the best assignments even more, while not being affected significantly by an increase in the distance of a bad assignment. Bad assignments are predicted diffraction spots that do not correlate to corresponding measured peaks. This happens because planes with low structure factors or diffracting neutron wavelengths with low intensity in the incident beam might not provide enough diffracted intensity for the peak to be segmented successfully by the watershed algorithm. Two overlapping distributions which can be described using a Gaussian fit arise when applying the logarithm to the assignments' euclidean distances, as shown in Fig. 4. This approach is not only used for defining the threshold and hence removal of bad assignments but also for validation of the orientation fit. Once the logarithms of the assignments' distances is calculated, the algorithm checks if a combination of two Gaussian curves can be fitted to it. The distribution including the lower 6/15 distances is expected to be the correctly assigned, while the distribution including the higher distances is expected to be the outliers. If such Gaussian mixture model converges, the median of the good assignment distances is checked to be lower than a threshold value chosen by the user. We consider that a grain was found if there are more than three diffraction spots assigned per angle, and the median of the assignments' distance is lower than 2 mm. If only one Gaussian distribution can be fitted, the grain is accepted as correct if it fulfills the same conditions regarding number of spots per angle and median of the assignments' distance. Once a grain has been identified as valid, the output parameters,D,N,X ω andr are stored, and the peaks correctly assigned are removed from the list of peaks to be fitted. The fitting process is then repeated for the next grain starting with the next best median distance obtained from the seeding process, until no seed orientation fulfills the described fitting criterion. Then the seeding process starts again creating a new list of m median distances with the peaks yet to be assigned to a grain. If there is no grain found within the first fitting iteration afterwards, the search finishes as seen in Fig. 3b.

Global fitting
Global fitting refines all the grain and detector parameters found during the single grain fitting procedure. Diffraction spots from all the predicted grains are generated using the orientations and positions found during the single grain fitting, and are optimized simultaneously by comparison with the peaks through a constrained optimization algorithm. The global fitting algorithm uses again 6 variables per detector plus 6 variables per grain, which in the case of the α-Fe sample led to a constrained optimization function with 156 variables. The global fitting creates a competition between the predicted grains to have the best fit for all diffraction spots simultaneously, therefore some of the assigned measured spots are moved to a different grain than the initial one. The output of this process is the final result of the developed indexing procedure and provides a list of grain orientations and positions with the corresponding assigned diffraction spots classified by hkl and ω.

Analysis
The output of the indexing process provides: • The number of grains found This information can be used to evaluate the validity of the fitted grains in various ways, by computing how many peaks have been assigned to every grain or what is the mean and median assignment distance of every found grain. Moreover, given that the peaks were cropped and stored before the indexing process, the information from the indexing can now be combined with the shape and intensity of every peak for further analysis, such as the estimation of relative volume of every grain in the sample.
In order to estimate the relative size of the grains from the diffracted intensities, we use only the diffraction data from a single family of planes. In that way, the neutron wavelengths scattered for a given Bragg angle are equivalent and therefore their intensities can be compared directly. Since the intensity of the diffraction spot is in first order proportional to the volume, we compute the average summed intensity among all spots for every different grain. The value is a relative measure of the volume of every grain, that can be used for comparison. The relative volume is calculated by

7/15
and where RV i is the relative volume of the i th grain, N is the number of grains, V i is the average intensity (I i, j ) of the diffraction spots ( j) from grain i restricted to a single family of planes (hkl) and a narrow θ interval, θ ∈ [θ 1 , θ 2 ].
The position and orientation of the grains (relative to each other) fitted during the indexing process can be combined with the estimated volume of the grains for the reconstruction of 3D grain maps.

Performance
In order to test the performance of the tool under well-known conditions, synthetic (simulated) data sets with well known grain positions and orientations are used for fitting. These tests give an idea of the number of grains that can be successfully indexed by the code under realistic experimental conditions for the cubic and tetragonal symmetries. The parameters defined for the benchmarks are: These parameters are set to the corresponding values unless explicitly said otherwise. We aimed to simulate that measurements often have bright spots and other artifacts, while sometimes experimental peaks are not intense enough to be segmented. For example: with a 5% noise we remove 5% of the simulated peaks and we add the same amount of peaks randomly placed in the synthetic data set. This is done for every grain and every ω angular step individually.
Since the position and orientation of the synthetic grains is defined a-priori, we can evaluate the spatial and angular precision of the solver and fitting criteria chosen. Based on a detector resolution of 0.25 mm and sample-to-detector distances of 160 mm, the maximum spatial and angular deviations accepted are 0.43 mm and 1 • respectively. As can be seen in Fig. 5a the solver can handle up to 50% of missing spots substituted by noise without missing any grain. Even with a rate of 75% of peaks substituted with noise, all ten grains can be found for the tetragonal symmetry, while nine grains could be found for the cubic symmetry. When increasing the spot substitution to 90%, no grain could be found successfully with neither the cubic nor the tetragonal symmetry. The robustness of Laue3DND with respect to noisy data sets relies on a large angular range and number of angular steps, which means better statistics. Figure 5b shows the performance of the code with synthetic data sets using the parameters described previously, but trying to find the minimum amount of angular steps required to find 100% of the simulated grains. For a sample with 10 grains and a data set with 360 projections around 360 • , only six projections are required to index the 10 grains successfully. Nevertheless, a higher number of grains in the sample would likely increase this requirement.
Finally, a test was performed following the parameters described previously, but simulating 100 grains of Fe instead of 10 within the same sample volume. Although the computation time was increased significantly up to two days, the solver was capable of successfully finding 97 grains out of 100 grains generated. Based on these results, we consider that complications arising from the indexing of a 100-grain sample might not be due to the indexing algorithm, but most likely the limitations of the watershed algorithm to deal with peak overlap.

Experimental setup
The experiments presented in this work were performed at the E11 beam port of Helmholtz-Zentrum Berlin (HZB) with the FALCON 15 Laue diffractometer installed. FALCON is composed of a back-diffraction and a forward-diffraction detector, with sizes of 400x400 mm, 4000x4000 px CCDs and a 6 LiF-ZnS based scintillator with 250 µm thickness. The thickness of the scintillator screen establishes a compromise between light output and resolution 16 . This scintillator provides a good light output at the expense of a resolution limited to 250 µm. The thermal neutron spectrum has half of the highest intensity for the wavelengths of 0.8 Å and 3.2 Å. The center of the thermal neutron beam traverses the rotation axis of the sample holder, placed in between two far field detectors at 160 mm from the sample, as depicted in Fig. 1a.
During the experiment, the whole sample is illuminated by the white neutron beam, diffracting simultaneously from all crystallites in the sample. Images are acquired in rotation steps of ∆ω within the largest possible angular range (241 • in our experiments) in order to provide data from as many diffraction spots (peaks) per crystallite as possible. The identification and segmentation of the individual peaks is obtained through a watershed algorithm 17 . From this we obtain a list of all the peaks CMS coordinates, and a corresponding set of cropped peaks from the experimental data, which are an ensemble of pixels showing the 2D peak profile. The list of CMS coordinates is used in the solver for the indexing procedure by comparison with the simulated spots, while the cropped peaks can be later used for analysis of the intensity.
Once the peaks CMS and profiles have been extracted, the forward model then predicts different diffraction patterns, so that the solver can find fits with the experimental data.

Fe oligocrystal
The first sample was a Fe oligocrystalline cylinder of 5 mm diameter and 5 mm height. The sample was measured using both forward and backward detectors with ∆ω=1 • over 241 • and 10 seconds of exposure time per angular step. Figure 6 shows two examples of backward (a) and forward (b) scattered Laue data from the oligocrystalline Fe sample. The peaks of the 24 grains found are already color coded and indexed in the image. It can be seen that some spots could not been indexed. The reason is foremost that the process for thresholding from the logarithm of the assignments will inevitably disregard some of the valid but not well fitted spots and secondly, that larger grains will generate more visible diffraction spots than smaller ones. This might only generate few visible spots from the largest structure factors for small grains, being not sufficient for the algorithm to successfully identify a grain. Hence, there might be grains in the sample yet to be found which do not pass the criteria discussed in the solver section (at least three good assignments per angle and a median assignment distance smaller than 2 mm). Figure 6 shows the path followed by every identified diffraction spot in back (c) and forward (d) scattering over the 241 steps of ω for a specific grain (grain 12). Table 1 provides the percentage of extracted diffraction peaks that can be assigned to the individual grains identified. The algorithm has been able to assign 40% of the diffraction peaks (segmented with the watershed algorithm) with statistical significance. This value can be increased by lowering the criteria for a valid grain found, for instance requiring two or only one spot per scattering angle to consider a grain found, or accepting a large distance between predicted spots and experimental peaks as valid, at the cost of higher chances of making wrong assignments. Although the acceptance criteria could be less conservative, we consider that these values give enough statistical certainty to validate the indexing method. Figure 7 is a representation of the solution found by the solver, including the position of the CMS of every grain found in the sample and its relative size calculated from the integrated intensities of the diffraction spots using equation (13). Note how the positions of the grains found inside the α-Fe sample represent a cylindrical geometry with 5 mm diameter and 5 mm height, which are the sample's dimensions. Taking into account that the boundaries for the grain CMS positions was set to 1 cm 3 , the 9/15 Figure 6. Representation of the fitted solution against the experimental data. Top: Diffraction planes indexed in back scattering (a) and forward scattering (b) for all 24 grains. Bottom: Scatter plot with the predicted and peaks positions for every omega in back scattering (c) and forward scattering (d) for grain 12. The assignments will appear blue or orange depending on which of the two spots appears to be in front from the current viewpoint.
code was able to find all the grains inside the sample's volume without strong constraints. This underlines the accuracy of the applied method.

YBaCuFeO 5 oligocrystal
The second sample is an oligocrystal of the high-temperature multiferroic candidate YBaCuFeO 5 with layered perovskite structure 18,19 , in which several grains with a common c-axis and slightly different orientation in the ab plane were formed during the process of crystal growth 20 . The main objective concerning this sample was to identify the number of grains and the respective misalignment of the 0 0 l direction, as well as the relative contribution of each domain to the diffracted signal (i.e. size distribution). The sample was measured using the forward detector with ∆ω=1 • over 241 • and 230 seconds of exposure time per angular step. Table 2 summarizes some statistics evaluated from the fit and highlights stronger variations as compared to the previous Fe reference sample. A consideration to be taken into account in this case is the fact that the misorientation between the 0 0 l directions of the different crystallites is quite small. Hence, it has to be expected that spot overlap from the reflections of this plane, more than for the others, might lead to segmentation of two actual spots as one and hence biased spot assignments.
In Fig. 8a peaks overlap in the horizontal line around the center of the detector underline that the different grains within the sample have very similar orientation with respect to the 0 0 l direction. The orientations of the two other main axis of the crystal lattice have a wider spread, generating less overlap and thus making them easier to identify and to distinguish. Figure 8b shows a 3D representation of the grain's positions within the sample and their relative sizes, calculated proportionally to the average integrated intensities of all diffraction spots from a single family. The box around the plotted prisms represents the search space defined as boundary conditions for the positions of the grains, which is slightly larger than the sample's size.

Discussion and Conclusions
Laue3DND has been presented using experimental data for α-Fe and YBaCuFeO 5 oligocrystals and corresponding synthetic data sets. This method has shown to be capable of indexing 97/100 grains on a cubic synthetic sample under realistic conditions with spatial and angular resolution of 0.43 mm and 1 • respectively. On experimental data sets, Laue3DND has been able to retrieve 24 grains from a α-Fe sample and 9 grains from a YBaCuFeO 5 . No significant differences in the code performance have been found between the cubic and tetragonal symmetries. One of the strengths of this method is the robustness towards noise and incomplete data sets, as shown in Fig. 5a, in which all the grains are found successfully even when 50% of the peaks in the data set are substituted for random noise. A very exciting case for data sets with few projections as shown in Fig. 5b is the possibility to generate 4D grain maps (obtaining a 3D grain map every few minutes), in order to study the evolution of the grain macrostructure of a sample under thermal or mechanical loading.
It is important to underline the fact that this 3D grain mapping method is several times more efficient than previous neutron methods 6,8 . Since this method does not require to select the energy of the incoming neutrons, the necessary neutron instrumentation is required and the effective neutron flux on the sample is increased between one and two orders of magnitude in comparison with other methods. In practice, this means that our new method represents a remarkable reduction in the measurements time required for neutron 3D grain mapping. The presented samples required exposure times of 10 seconds (Fe) and 215 seconds (YBaCuFeO 5 ) per angular step respectively, compared to 250 total seconds (Al 6 ) and 1 hour (Fe 8 ) for previous neutron methods.
Given the percentage of spots which have been assigned correctly for both samples, it is expected that some grains can still be found. Note that the criteria used for the code to confirm the existence of a grain are based on the amount of peak-spot 11/15 Figure 7. Iron sample of 5 mm diameter and 5 mm height, with low opacity overlaid with coloured cubes representing the different grains and their and their orientations. The volume of the cubes is based on the integrated intensity of the diffraction spots of every grain and is proportional to the relative volume of the grain following equation (13 assignments per angle and the distance of these assignments. Since the thresholding method ( Fig. 4) will inevitably discard some right (and include some wrong) assignments, small grains with a low number of visible diffraction spots might not be accepted.
The current spatial resolution of Laue3DND (0.43 mm) is limited by the resolution of the detector itself, as well as the sample to detector distance. The experiments and simulations provided show an orientation resolution comparable to that of other X-ray and neutron methods (1 • ), arguably thanks to the far field setup, which is also a limiting factor for the spatial resolution. A first step to bring the spatial resolution closer to other neutron methods (0.1 mm) would be to set at least one of the two detectors in near field position.
The relative volume of the crystallites has been estimated by comparison of integrated intensities of equivalent reflections (13). Further quantitative analysis of the diffracted signal will be explored in future work, in order to recover the shape of the grain boundaries 21,22 or obtain information about the wavelength distribution of the incoming beam.
The maximum sample volume measure with Laue3DND so far has been the α-Fe sample with 5 mm diameter and 5 mm height. However, the FALCON beamline can accommodate samples up to 2 cm 3 , limited by the beam diameter. Samples larger than 2 cm in one of their dimensions could be oriented vertically and translated after every rotation. Larger volumes could be analyzed at imaging beamlines using large detectors for forward scattering, but the increase in beam size will reduce the collimation of the beam and the orientation resolution of the measurement.
The minimum grain size accessible with neutrons is physically limited by the coherent scattering cross section of the materials, the signal-to-noise ratio and the spatial resolution of the instrument used. We consider the practical limit for grain indexing with neutrons to be in the range of tens of micrometers for the strongest coherent nuclei, and in the order of hundreds 12/15 Figure 8. a) Representation of the fitted solution against the experimental data for the Perovskite: Diffraction planes indexed using only the forward scattering detector. b) Plot with the 9 grains found by the algorithm inside the layered perovskite sample, represented by coloured prisms representing the different grains and their orientations. The volume of the cubes is based on the integrated intensity of the diffraction spots of every grain and is proportional to the relative volume of the grain following equation (13).
of micrometers for others.
The number of grains that can potentially be indexed is harder to estimate since there are more factors in play. Peak overlap is the biggest concern when trying to index a large amount of grains, and it can be tackled by improving the angular resolution of the detection systems. One can do that by reducing the thickness of the scintillator screens at the expense of light output, in the case of scintillator based camera setups like FALCON. This should be taken into consideration for further Laue instrumentation development, since typical neutron imaging detectors already use a variety of scintillator screens suited to the resolution requirements.
The forward model, solver and analysis tools have been so far used only to analyze neutron data, but adapting them to X-ray Laue data would not require fundamental changes. From the instrumentation and sample point of view, some changes are already easily implemented in the current form of the code, such as changes on the wavelength spectrum, detector geometry or sample sizes.
The study of lower symmetry crystals would not require fundamental changes on the Laue3DND indexing algorithms, except for monoclinic and triclinic systems. In order to index crystals from these two systems, crystal misorientation defined with quaternions would have to be coded as an alternative to the Rodrigues formulation. This is because the fundamental zone for triclinic and monoclinic symmetries corresponds to the entire Rodrigues space 23 , and therefore the seeding time is infinite. Any other symmetry can be searched with the Rodrigues formulation in the present form, although the seeding time will increase with the asymmetry of the crystal.
The code presented in this work was written by Marc Raventós, Søren Schmidt and Stavros Samothrakitis for MATLAB, and can be found in the GitHub repository Laue3DND doi: 10.5281/zenodo.1553164.