Rietveld refinement

Rietveld refinement is a technique described by Hugo Rietveld for use in the characterisation of crystalline materials. The neutron and X-ray diffraction of powder samples results in a pattern characterised by reflections (peaks in intensity) at certain positions. The height, width and position of these reflections can be used to determine many aspects of the material's structure.

The Rietveld method uses a least squares approach to refine a theoretical line profile until it matches the measured profile. The introduction of this technique was a significant step forward in the diffraction analysis of powder samples as, unlike other techniques at that time, it was able to deal reliably with strongly overlapping reflections.

The method was first implemented in 1967,^[1] and reported in 1969^[2] for the diffraction of monochromatic neutrons where the reflection-position is reported in terms of the Bragg angle, 2θ. This terminology will be used here although the technique is equally applicable to alternative scales such as x-ray energy or neutron time-of-flight. The only wavelength and technique independent scale is in reciprocal space units or momentum transfer Q, which is historically rarely used in powder diffraction but very common in all other diffraction and optics techniques. The relation is

${\displaystyle Q={\frac {4\pi \sin \left(\theta \right)}{\lambda }}.}$

The most common powder X-ray diffraction (XRD) refinement technique used today is based on the method proposed in the 1960s by Hugo Rietveld.^[2] The Rietveld method fits a calculated profile (including all structural and instrumental parameters) to experimental data. It employs the non-linear least squares method, and requires the reasonable initial approximation of many free parameters, including peak shape, unit cell dimensions and coordinates of all atoms in the crystal structure. Other parameters can be guessed while still being reasonably refined. In this way one can refine the crystal structure of a powder material from PXRD data. The successful outcome of the refinement is directly related to the quality of the data, the quality of the model (including initial approximations), and the experience of the user.

The Rietveld method is an incredibly powerful technique which began a remarkable era for powder XRD and materials science in general. Powder XRD is at heart a very basic experimental technique with diverse applications and experimental options. Despite being slightly limited by the one-dimensionality of PXRD data and limited resolution, powder XRD's power is astonishing. It is possible to determine the accuracy of a crystal structure model by fitting a profile to a 1D plot of observed intensity vs angle. Rietveld refinement requires a crystal structure model, and offers no way to come up with such a model on its own. However, it can be used to find structural details missing from a partial or complete ab initio structure solution, such as unit cell dimensions, phase quantities, crystallite sizes/shapes, atomic coordinates/bond lengths, micro strain in crystal lattice, texture, and vacancies.^[3][4]

Before exploring Rietveld refinement, it is necessary to establish a greater understanding of powder diffraction data and what information is encoded therein in order to establish a notion of how to create a model of a diffraction pattern, which is of course necessary in Rietveld refinement. A typical diffraction pattern can be described by the positions, shapes, and intensities of multiple Bragg reflections. Each of the three mentioned properties encodes some information relating to the crystal structure, the properties of the sample, and the properties of the instrumentation. Some of these contributions are shown in Table 1, below.

The structure of a powder pattern is essentially defined by instrumental parameters and two crystallographic parameters: unit cell dimensions, and atomic content and coordination. So, a powder pattern model can be constructed as follows:

1. Establish peak positions: Bragg peak positions are established from Bragg's law using the wavelength and d-spacing for a given unit cell.
2. Determine peak intensity: Intensity depends on the structure factor, and can be calculated from the structural model for individual peaks. This requires knowledge of the specific atomic coordination in the unit cell and geometrical parameters.
3. Peak shape for individual Bragg peaks: Represented by functions of the FWHM (which vary with Bragg angle) called the peak shape functions. Realistically ab initio modelling is difficult, and so empirically selected peak shape functions and parameters are used for modelling.
4. Sum: The individual peak shape functions are summed and added to a background function, leaving behind the resultant powder pattern.

It is easy to model a powder pattern given the crystal structure of a material. The opposite, determining the crystal structure from a powder pattern, is much more complicated. A brief explanation of the process follows, though it is not the focus of this article.

To determine structure from a powder diffraction pattern the following steps should be taken. First, Bragg peak positions and intensities should be found by fitting to a peak shape function including background. Next, peak positions should be indexed and used to determine unit cell parameters, symmetry, and content. Third, peak intensities determine space group symmetry and atomic coordination. Finally, the model is used to refine all crystallographic and peak shape function parameters. To do this successfully, there is a requirement for excellent data which means good resolution, low background, and a large angular range.

For general application of the Rietveld method, irrespective of the software used, the observed Bragg peaks in a powder diffraction pattern are best described by the so-called peak shape function (PSF). The PSF is a convolution of three functions: the instrumental broadening ${\displaystyle \Omega (\theta )}$, wavelength dispersion ${\displaystyle \Lambda (\theta )}$, and the specimen function ${\displaystyle \Psi (\theta )}$, with the addition of a background function, ${\displaystyle b(\theta )}$. It is represented as follows:

${\displaystyle PSF(\theta )=\Omega (\theta )\otimes \Lambda (\theta )\otimes \Psi (\theta )+b(\theta )}$,

where ${\displaystyle \otimes }$ denotes convolution, which is defined for two functions ${\displaystyle f}$ and ${\displaystyle g}$ as an integral:

${\displaystyle f(t)\otimes g(t)=\int _{-\infty }^{\infty }f(\tau )g(t-\tau )d\tau =\int _{-\infty }^{\infty }g(\tau )f(t-\tau )d\tau }$

The instrumental function depends on the location and geometry of the source, monochromator, and sample. Wavelength function accounts for the distribution of the wavelengths in the source, and varies with the nature of the source and monochromatizing technique. The specimen function depends on several things. First is dynamic scattering, and secondly the physical properties of the sample such as crystallite size, and microstrain.

A short aside: unlike the other contributions, those from the specimen function can be interesting in materials characterization. As such, average crystallite size, ${\displaystyle \tau }$, and microstrain, ${\displaystyle \varepsilon }$, effects on Bragg peak broadening, ${\displaystyle \beta }$ (in radians), can be described as follows, where ${\displaystyle k}$ is a constant:

${\displaystyle \beta ={\frac {\lambda }{\tau \cdot \cos \theta }}}$ and ${\displaystyle \beta =\kappa \cdot \varepsilon \cdot \tan \theta }$.

Returning to the peak shape function, the goal is to correctly model the Bragg peaks which exist in the observed powder diffraction data. In the most general form, the intensity, ${\displaystyle Y(i)}$, of the ${\displaystyle i^{\text{th}}}$ point (${\displaystyle 1\leq i\leq n}$, where ${\displaystyle n}$ is the number of measured points) is the sum of the contributions ${\displaystyle y_{k}}$ from the ${\displaystyle m}$ overlapped Bragg peaks (${\displaystyle 1\leq k\leq m}$), and the background, ${\displaystyle b(i)}$, and can be described as follows:

${\displaystyle Y(i)=b(i)+\sum _{k=1}^{m}{I_{k}[y_{k}(x_{k})]}}$

where ${\displaystyle I_{k}}$ is the intensity of the ${\displaystyle k^{\text{th}}}$ Bragg peak, and ${\displaystyle x_{i}=2\theta _{i}-2\theta _{k}}$. Since ${\displaystyle I_{k}}$ is a multiplier, it is possible to analyze the behaviour of different normalized peak functions ${\displaystyle y(x)}$ independently of peak intensity, under the condition that the integral over infinity of the PSF is unity. There are various functions that can be chosen to do this with varying degrees of complexity. The most basic functions used in this way to represent Bragg reflections are the Gauss, and Lorentzian functions. Most commonly though, is the pseudo-Voigt function, a weighted sum of the former two (the full Voigt profile is a convolution of the two, but is computationally more demanding). The pseudo-Voigt profile is the most common and is the basis for most other PSF's. The pseudo-Voigt function can be represented as:

${\displaystyle y(x)=V_{p}(x)=n*G(x)+(1-n)*L(x)}$,

where

${\displaystyle G(x)={\frac {C_{G}^{\frac {1}{2}}}{{\sqrt {\pi }}H}}e^{-C_{G}x^{2}}}$

and

${\displaystyle L(x)={\frac {C_{L}^{\frac {1}{2}}}{{\sqrt {\pi }}H'}}\left(1+C_{L}x^{2}\right)^{-1}}$

are the Gaussian and Lorentzian contributions, respectively.

Thus,

${\displaystyle V_{p}(x)=\eta {\frac {C_{G}^{\frac {1}{2}}}{{\sqrt {\pi }}H}}e^{-C_{G}x^{2}}+(1-\eta ){\frac {C_{L}^{\frac {1}{2}}}{{\sqrt {\pi }}H'}}(1+C_{L}x^{2})^{-1}.}$

where:

• ${\displaystyle H}$ and ${\displaystyle H'}$ are the full widths at half maximum (FWHM)
• ${\displaystyle x={\frac {2\theta _{i}-2\theta _{k}}{H_{k}}}}$ is essentially the Bragg angle of the ${\displaystyle i^{\text{th}}}$ point in the powder pattern with its origin in the position of the ${\displaystyle k^{\text{th}}}$ peak divided by the peak's FWHM.
• ${\displaystyle C_{G}=4\ln 2}$, ${\displaystyle C_{L}=4}$ and ${\textstyle {\frac {C_{G}^{\frac {1}{2}}}{{\sqrt {\pi }}H}}}$ and ${\textstyle {\frac {C_{L}^{\frac {1}{2}}}{{\sqrt {\pi }}H'}}}$ are normalization factors such that ${\textstyle \int _{-\infty }^{\infty }G(x)dx=1}$ and ${\textstyle \int _{-\infty }^{\infty }L(x)dx=1}$ respectively.
• ${\displaystyle H^{2}=U\tan ^{2}\theta +V\tan \theta +W}$, known as the Caglioti formula,^[6] is the FWHM as a function of ${\displaystyle \theta }$ for Gauss, and pseudo-Voigt profiles. ${\displaystyle U}$, ${\displaystyle V}$, and ${\displaystyle W}$ are free parameters.
• ${\textstyle H'={\frac {X}{\cos \theta }}+Y\tan \theta }$ is the FWHM vs. ${\displaystyle 2\theta }$ for the Lorentz function. ${\displaystyle X}$ and ${\displaystyle Y}$ are free variables
• ${\displaystyle \eta =\eta _{0}+\eta _{1}2\theta +\eta _{2}\theta ^{2}}$, where ${\displaystyle 0\leq \eta \leq 1}$ is the pseudo-Voigt mixing parameter, and ${\displaystyle \eta _{0,1,2}}$ are free variables.

The pseudo-Voigt function, like the Gaussian and Lorentz functions, is a centrosymmetric function, and as such does not model asymmetry. This can be problematic for non-ideal powder XRD data, such as those collected at synchrotron radiation sources, which generally exhibit asymmetry due to the use of multiple focusing optics.

The Finger–Cox–Jephcoat function is similar to the pseudo-Voigt, but has better handling of asymmetry, which is treated in terms of axial divergence. The function is a convolution of pseudo-Voigt with the intersection of the diffraction cone and a finite receiving slit length using two geometrical parameters, ${\displaystyle S/L}$, and ${\displaystyle H/L}$, where ${\displaystyle S}$ and ${\displaystyle H}$ are the sample and the detector slit dimensions in the direction parallel to the goniometer axis, and ${\displaystyle L}$ is the goniometer radius.^[7]

The shape of a powder diffraction reflection is influenced by the characteristics of the beam, the experimental arrangement, and the sample size and shape. In the case of monochromatic neutron sources the convolution of the various effects has been found to result in a reflex almost exactly Gaussian in shape. If this distribution is assumed then the contribution of a given reflection to the profile ${\displaystyle y_{i}}$ at position ${\displaystyle 2\theta _{i}}$ is:

${\displaystyle y_{i}=I_{k}\exp \left[{\frac {-4\ln \left(2\right)}{H_{k}^{2}}}\left(2\theta _{i}-2\theta _{k}\right)^{2}\right]}$

where ${\displaystyle H_{k}}$ is the full width at half peak height (full-width half-maximum), ${\displaystyle 2\theta _{k}}$ is the center of the reflex, and ${\displaystyle I_{k}}$ is the calculated intensity of the reflex (determined from the structure factor, the Lorentz factor, and multiplicity of the reflection).

At very low diffraction angles the reflections may acquire an asymmetry due to the vertical divergence of the beam. Rietveld used a semi-empirical correction factor, ${\displaystyle A_{s}}$ to account for this asymmetry:

${\displaystyle A_{s}=1-\left[{\frac {P\left(2\theta _{i}-2\theta _{k}\right)^{2}}{\tan \theta _{k}}}\right]}$

where ${\displaystyle P}$ is the asymmetry factor and ${\displaystyle s}$ is ${\displaystyle +1}$, ${\displaystyle 0}$, or ${\displaystyle -1}$ depending on the difference ${\displaystyle 2\theta _{i}-2\theta _{k}}$ being positive, zero, or negative respectively.

At a given position more than one diffraction peak may contribute to the profile. The intensity is simply the sum of all reflections contributing at the point ${\displaystyle 2\theta _{i}}$.

For a Bragg peak ${\displaystyle (hkl)}$, the observed integrated intensity, ${\displaystyle I_{hkl}}$, as determined from numerical integration is

${\displaystyle I_{hkl}=\sum _{i=1}^{j}(Y_{i}^{obs}-b_{i})}$,

where ${\displaystyle j}$ is the total number of data points in the range of the Bragg peak. The integrated intensity depends on multiple factors, and can be expressed as the following product:

${\displaystyle I_{hkl}=K\times p_{hkl}\times L_{\theta }\times P_{\theta }\times A_{\theta }\times T_{hkl}\times E_{hkl}\times |F_{hkl}|^{2}}$

where:

• ${\displaystyle K}$: scale factor
• ${\displaystyle p_{hkl}}$: multiplicity factor, which accounts for symmetrically equivalent points in the reciprocal lattice
• ${\displaystyle L_{\theta }}$: Lorentz multiplier, defined by diffraction geometry
• ${\displaystyle P_{\theta }}$: polarization factor
• ${\displaystyle A_{\theta }}$: absorption multiplier
• ${\displaystyle T_{hkl}}$: preferred orientation factor
• ${\displaystyle E_{hkl}}$: extinction factor (often neglected as it is usually insignificant in powders)
• ${\displaystyle F_{hkl}}$: structure factor as determined by the crystal structure of the material

The width of the diffraction peaks are found to broaden at higher Bragg angles. This angular dependency was originally represented by

${\displaystyle H_{k}^{2}=U\tan ^{2}\theta _{k}+V\tan \theta _{k}+W}$

where ${\displaystyle U}$, ${\displaystyle V}$, and ${\displaystyle W}$ are the half-width parameters and may be refined during the fit.

In powder samples there is a tendency for plate- or rod-like crystallites to align themselves along the axis of a cylindrical sample holder. In solid polycrystalline samples the production of the material may result in greater volume fraction of certain crystal orientations (commonly referred to as texture). In such cases the reflex intensities will vary from that predicted for a completely random distribution. Rietveld allowed for moderate cases of the former by introducing a correction factor:

${\displaystyle I_{\text{corr}}=I_{\text{obs}}\exp \left(-G\alpha ^{2}\right)}$

where ${\displaystyle I_{\text{obs}}}$ is the intensity expected for a random sample, ${\displaystyle G}$ is the preferred orientation parameter and ${\displaystyle \alpha }$ is the acute angle between the scattering vector and the normal of the crystallites.

The principle of the Rietveld method is to minimize a function ${\displaystyle M}$ which analyzes the difference between a calculated profile ${\displaystyle y^{\text{calc}}}$ and the observed data ${\displaystyle y^{\text{obs}}}$. Rietveld defined such an equation as:

${\displaystyle M=\sum _{i}W_{i}\left\{y_{i}^{\text{obs}}-{\frac {1}{c}}y_{i}^{\text{calc}}\right\}^{2}}$

where ${\displaystyle W_{i}}$ is the statistical weight and ${\displaystyle c}$ is an overall scale factor such that ${\displaystyle y^{\text{calc}}=cy^{\text{obs}}}$.

The fitting method used in Rietveld refinement is the non-linear least squares approach. A detailed derivation of non-linear least squares fitting will not be given here. Further detail can be found in Chapter 6 of Pecharsky and Zavalij's text 12 . There are a few things to note however. First, non-linear least squares fitting has an iterative nature for which convergence may be difficult to achieve if the initial approximation is too far from correct, or when the minimized function is poorly defined. The latter occurs when correlated parameters are being refined at the same time, which may result in divergence and instability of the minimization. This iterative nature also means that convergence to a solution does not occur immediately for the method is not exact. Each iteration depends on the results of the last which dictate the new set of parameters used for refinement. Thus, multiple refinement iterations are required to eventually converge to a possible solution.

Using non-linear least squares minimization, the following system is solved:

${\displaystyle {\begin{pmatrix}Y_{i}^{\text{calc}}=kY_{i}^{\text{obs}}\\\vdots \\Y_{n}^{\text{calc}}=kY_{n}^{\text{obs}}\end{pmatrix}}}$

where ${\displaystyle Y_{i}^{\text{calc}}}$ is the calculated intensity and ${\displaystyle Y_{i}^{\text{obs}}}$ is the observed intensity of a point ${\displaystyle i}$ in the powder pattern, ${\displaystyle k}$, is a scale factor, and ${\displaystyle n}$ is the number of measured data points. The minimized function is given by:

${\displaystyle \Phi =\sum _{i=1}^{n}w_{i}(Y_{i}^{\text{obs}}-Y_{i}^{\text{calc}})^{2}}$

where ${\displaystyle w_{i}}$ is the weight, and ${\displaystyle k}$ from the previous equation is unity (since ${\displaystyle k}$ is usually absorbed in the phase scale factor). The summation extends to all ${\displaystyle n}$ data points. Considering the peak shape functions and accounting for the overlapping of Bragg peaks because of the one-dimensionality of XRD data, the expanded form of the above equation for the case of a single phase measured with a single wavelength becomes:

${\displaystyle \Phi =\sum _{i=1}^{n}w_{i}{\biggl (}Y_{i}^{\text{obs}}-{\Bigl (}b_{i}+K\sum _{j=1}^{m}I_{j}y_{j}(x_{j}){\Bigr )}{\biggr )}^{2}}$

where:

• ${\displaystyle b_{i}}$ is the background at the ${\displaystyle i^{\text{th}}}$ data point.
• ${\displaystyle K}$ is the phase scale factor.
• ${\displaystyle m}$ is the number of Bragg reflections contributing to the intensity of the ${\displaystyle i^{\text{th}}}$ reflection.
• ${\displaystyle I_{j}}$ is the integrated intensity of the ${\displaystyle j^{\text{th}}}$ Bragg peak.
• ${\displaystyle y_{i}(x_{i})}$ is the peak shape function.

For a material that contains several phases (${\displaystyle p}$), the contribution from each is accounted for by modifying the above equation as follows:

${\displaystyle \Phi =\sum _{i=1}^{n}w_{i}{\biggl (}Y_{i}^{\text{obs}}-{\Bigl (}b_{i}+\sum _{l=1}^{p}K_{l}\sum _{j=1}^{m}I_{l,j}y_{l,j}(x_{l,j}){\Bigr )}{\biggr )}^{2}}$

It can easily be seen from the above equations that experimentally minimizing the background, which holds no useful structural information, is paramount for a successful profile fitting. For a low background, the functions are defined by contributions from the integrated intensities and peak shape parameters. But with a high background, the function being minimized depends on the adequacy of the background and not integrated intensities or peak shapes. Thus, a structure refinement cannot adequately yield structural information in the presence of a large background.

It is also worth noting the increased complexity brought forth by the presence of multiple phases. Each additional phase adds to the fitting, more Bragg peaks, and another scale factor tied to corresponding structural parameters, and peak shape. Mathematically they are easily accounted for, but practically, due to the finite accuracy and limited resolution of experimental data, each new phase can lower the quality and stability of the refinement. It is advantageous to use single phase materials when interested in finding precise structural parameters of a material. However, since the scale factors of each phase are determined independently, Rietveld refinement of multi phase materials can quantitatively examine the mixing ratio of each phase in the material.

Generally, the background is calculated as a Chebyshev polynomial. In GSAS and GSAS-II they appear as follows. Again, background is treated as a Chebyshev polynomial of the first kind ("Handbook of Mathematical Functions", M. Abramowitz and IA. Stegun, Ch. 22), with intensity given by:

${\displaystyle I_{i}=\sum _{j=1}^{N}P_{j}T'_{j-1}}$

where ${\displaystyle T'_{j-1}}$ are the coefficients of the Chebyshev polynomial taken from Table 22.3, pg. 795 of the Handbook. The coefficients have the form:

${\displaystyle T'_{n}=\sum _{m=0}^{i-1}C_{m}X^{m}}$

and the values for ${\displaystyle C_{m}}$ are found in the Handbook. The angular range (${\displaystyle 2\theta }$) is converted to ${\displaystyle X}$ to make the Chebyshev polynomial orthogonal by

${\displaystyle X={\frac {2}{\tau }}-1}$

And, the orthogonal range for this function is –1 to +1.

Now, given the considerations of background, peak shape functions, integrated intensity, and non-linear least squares minimization, the parameters used in the Rietveld refinement which put these things together can be introduced. Below are the groups of independent least squares parameters generally refined in a Rietveld refinement.

• Background parameters: usually 1 to 12 parameters.
• Sample displacement: sample transparency, and zero shift corrections. (move peak position)
• Multiple peak shape parameters.
  • FWHM parameters: i.e. Caglioti parameters (see section 3.1.2)
  • Asymmetry parameters (FCJ parameters)
• Unit cell dimensions
  • one to six parameters (a, b, c, α, β, γ), depending on the crystal family/system, for each present phase.
• Preferred orientation, and sometimes absorption, porosity, and extinction coefficients, which can be independent for each phase.
• Scale factors (for each phase)
• Positional parameters of all independent atoms in the crystal model (generally 0 to 3 per atom).
• Population parameters
  • Occupation of site positions by atoms.
• Atomic displacement parameters
  • Isotropic and anisotropic (temperature) parameters.

Each Rietveld refinement is unique and there is no prescribed sequence of parameters to include in a refinement. It is up to the user to determine and find the best sequence of parameters for refinement. It is worth noting that it is rarely possible to refine all relevant variables simultaneously from the beginning of a refinement, nor near the end since the least squares fitting will be destabilized or lead to a false minimum. It is important for the user to determine a stopping point for a given refinement. Given the complexity of Rietveld refinement it is important to have a clear grasp of the system being studied (sample, and instrumentation) to ensure that results are accurate, realistic, and meaningful. High data quality, a large enough ${\displaystyle 2\theta }$ range, and a good model – to serve as the initial approximation in the least squares fitting – are necessary for a successful, reliable, and meaningful Rietveld refinement.

Since refinement depends on finding the best fit between a calculated and experimental pattern, it is important to have a numerical figure of merit quantifying the quality of the fit. Below are the figures of merit generally used to characterize the quality of a refinement. They provide insight to how well the model fits the observed data.

Profile residual (reliability factor):

${\displaystyle R_{p}=\sum _{i}^{n}{\frac {|Y_{i}^{\text{obs}}-Y_{i}^{\text{calc}}|}{\sum _{i}^{n}Y_{i}^{\text{obs}}}}\times 100\%}$

Weighted profile residual:

${\displaystyle R_{wp}=\left(\sum _{i}^{n}{\frac {w_{i}(Y_{i}^{\text{obs}}-Y_{i}^{\text{calc}})^{2}}{\sum _{i}^{n}w_{i}(Y_{i}^{\text{obs}})^{2}}}\right)^{\frac {1}{2}}\times 100\%}$

Bragg residual:

${\displaystyle R_{B}=\sum _{j}^{m}{\frac {|I_{j}^{\text{obs}}-I_{j}^{\text{calc}}|}{\sum _{i}^{n}I_{j}^{\text{obs}}}}\times 100\%}$

Expected profile residual:

${\displaystyle R_{\text{exp}}=\left({\frac {n-p}{\sum _{i}^{n}w_{i}(Y_{i}^{\text{obs}})^{2}}}\right)^{\frac {1}{2}}\times 100\%}$

Goodness of fit:

${\displaystyle \mathrm {X} ^{2}=\sum _{i}^{n}{\frac {(Y_{i}^{\text{obs}}-Y_{i}^{\text{calc}})^{2}}{n-p}}=\left({\frac {R_{wp}}{R_{\text{exp}}}}\right)}$

It is worth mentioning that all but one (${\displaystyle R_{B}}$) figure of merit include a contribution from the background. There are some concerns about the reliability of these figures, as well there is no threshold or accepted value which dictates what represents a good fit.^[8] The most popular and conventional figure of merit used is the goodness of fit which should approach unity given a perfect fit, though this is rarely the case. In practice, the best way to assess quality is a visual analysis of the fit by plotting the difference between the observed and calculated data plotted on the same scale.

• Pecharsky, Vitalij K.; Zavalij, Peter Y. (2009). Fundamentals of powder diffraction and structural characterization of materials (2nd ed.). New York: Springer. ISBN 978-0-387-09579-0. OCLC 314182615.
• V. Emond (2018). Optimizing and Analyzing X-Ray Powder Diffraction of Orthosilicate Cathodes using a Combined Synchrotron X-Ray Diffraction and Absorption Spectroscopy Setup (MS thesis). University of Guelph. hdl:10214/13005.