Statistical measures of angular correlations in amorphous materials from electron nano-diffraction in the scanning/transmission electron microscope

Crystallography employing conventional large-volume diffraction has enabled the firm connections between structure and properties and structure and function that have solved many of the most difficult problems in materials science and biology. Disordered materials possess a large variety of local structural arrangements and pose a special challenge for crystallography. Often the local structures and symmetries that are responsible for observed phenomena in structurally complex disordered materials cannot be distinguished from conventional diffraction alone. In this article, we review analytical approaches for understanding local structure and symmetry from angular correlations in limited-volume diffraction patterns of amorphous materials, with a special focus on electron nano-diffraction. We discuss how these angular correlations can be interpreted in the context of dense, disordered, three-dimensional materials probed in a projection geometry and highlight the experimental challenges and considerations. New developments in this field are described whereby these angular correlations are statistically analyzed to probe the symmetry and variety of local structures, transformed to a real-space function that contains the 2-, 3- and 4-body particle correlations, and employed to develop reverse Monte Carlo models with more realistic higher-order correlations.


Introduction
Routine crystallographic structure determination has been one of the great underpinning methods that has enabled large advances in materials science, chemistry, physics, engineering, and biology in the past 100 years. The translational symmetry of crystals places such strict constraints on structure solutions from total scattering, that in structural biology it has become routine to solve structures with tens of thousands of non-hydrogen atoms from conventional diffraction measurements. In contrast, many materials are not crystals and have complex structures with elements of randomness as well as order and multiple levels of hierarchy that pose a challenge to crystallography. Examples of technological importance include molecular and metallic glasses, ion-irradiated amorphous solids, polymers, chars, activated and non-graphitizing carbons, glassy and jammed colloidal systems, chemically disordered alloys, and gels [1]. In these cases, the disorder imparts special properties-materials can be more homogenous, isotropic, tough, wear resistant, corrosion resistant, formable, hierarchical, desirably metastable, magnetically soft, or simply requiring less embodied energy to manufacture. In the extreme case of glasses, each atom occupies a unique position relative to others, and the unit cell is infinite. The large volume diffraction pattern from amorphous materials contains diffuse rings that reflect the isotropy of the system and average inter-particle separations. Structural models can be developed that are consistent with these diffracted intensities, but they are under-specified by the data, and thus non-unique. Often the local structure, or structures that are responsible for observed properties or behaviours, are not able to be identified from diffraction measurements alone. There is a large scope to develop a 'generalized crystallography' to better quantify the complexity and variation in structure that nature presents us with [2].
From the radial average of the large volume, isotropic diffraction pattern of an amorphous material, I(q) = ⟨I(q, ϕ)⟩ ϕ , where ⟨X⟩ ϕ = 1 2π´2 π 0 Xdϕ denotes averaging over the azimuthal angle ϕ, one can construct a reduced intensity function (or interference function), Φ(q), which is related by Fourier transform to the pair-correlation function (g(r)) of the material [3]. Here q = |⃗ q| = √ q 2 x + q 2 y = 4πsinθ/λ (for scattering angle 2θ) and r = √ x 2 + y 2 are the scattering vector magnitude and real-space co-ordinate, respectively.
In these equations for monatomic specimens, S(q) = I(q)/(N at f 2 (q)) is the structure factor, N at is the number of atoms in the diffracting volume, f (q) is the atomic scattering factor and ρ 0 is the average density [3]. The one-dimensional pair correlation function and related radial distribution function gives reliable, quantitative information about the average structure of the material and can be readily accessed using freely available software tools [4][5][6][7]. While these functions show the average interparticle distances, the symmetry of the local structure cannot be directly accessed.
Kam observed angular fluctuations in the scattered intensity from identical, sparse, randomly oriented, macro-molecules in suspension and proposed that these fluctuations could be employed to gain more insight into the symmetry of molecules than from the uniform intensity at a given q produced by conventional averaged diffraction experiments [8]. These angular fluctuations have been the focus of much recent, intensive study in the free-electron laser community [9][10][11][12][13][14] where, similarly, the structure of identical, isolated, arbitrarily oriented, bio-molecules is a major research focus. The electron microscopy community has also realized the potential of characterizing the angular fluctuations in limited-volume diffraction patterns to study short-to medium-range order in amorphous materials [15][16][17][18][19][20][21][22]. Here, the assembly of particles is dense-the central atom of one co-ordination polyhedron is the vertex atom of the neighboring polyhedron-and thus the interpretation of these angular fluctuations is more challenging. In both communities, increases in brightness, coherence, beam size tunability, aberration-reduction, and event-counting pixelated detectors have considerably facilitated experimental studies.
Several parameters have been proposed to quantify order from these fluctuations. These include the angular auto-and cross-correlation [8,15,22], the intensity cross-correlation [18], the intensity variance [20,23], and transforms to real-space of the average angular auto-correlation [9,11,12,14,24]. The ergodicity of many of the systems studied presents a unique opportunity as temporal and spatial sampling can give insight into the energetic landscape. Correlation lengths beyond the size of the probe have been directly measured from maps of order parameters from experiments where specimen stage or probe-beam mapping has been realized [18,[25][26][27][28][29]. There is still an ongoing challenge in interpreting these intensity fluctuations, but despite this there has been an explosion of experimental activity. This paper is organized as follows: section 2 summarizes the experimental conditions and quantitative tools for extracting structural parameters from these angular correlations employed by the authors and other researchers. Further discussions in section 3 showcase promising new directions being actively investigated by the authors for using these angular correlations (i) as input to advanced statistical analysis including matrix factorization algorithms, (ii) to measure higher-order real-space correlation functions and (iii) to refine more realistic structural models.

Quantifying angular symmetries
The angular cross-correlation function can be expressed in the form: [10,22,24]: This can be simplified to the angular auto-correlation function in cases where q = q ′ [10,30,31] to examine angular symmetry at different length scales: As shown, this angular auto-correlation function can be decomposed into a Fourier series to extract the Fourier coefficients, C n q , at each q and thus quantify the magnitude of each angular symmetry, indexed here by the integer n [10,25,31]. When averaged over a q-range corresponding to the nearest-neighbour shell, these Fourier coefficients display a characteristic fingerprint of different projected local symmetries in the same way as the three-dimensional bond-orientational order (BOO) parameters [32] fingerprint the symmetry of particle arrangements around a central atom [31]. Contributions from different candidate local orientational orders can be extracted by fitting [26]. Figure 1(a) shows the geometry of scanning electron nano-diffraction [18] (SEND -also often called 'nanobeam diffraction' or 'nanobeam electron diffraction' [33]) in a modern scanning/transmission electron microscope (S/TEM). A quasi-parallel electron beam is focussed onto the specimen and scanned in the spatial coordinates, x and y. A detector is placed in the far-field plane and a full diffraction pattern is collected at each point in the scanned array. This detector could be a scintillator coupled to a charge-coupled device, but is more latterly an electron counting, pixellated detector. Increased detector speeds have made the implementation of such '4D-STEM' [34] with data dimensions x, y, q x , q y , much more tractable.

Experimental geometry
The flexible condenser optics in a modern S/TEM allows SEND to be implemented with a large variety of convergence angles (deep sub-mrad < α < 1-5 mrad ) and corresponding beam diameters (10-50 nm > d p > sub-nm), depending on the specifics of the instrument optics and lens excitations [33]. The final four elements of the pre-specimen optics are composed of a condenser aperture and second condenser lens (CL2), a third condenser lens (CL3), a condenser mini-lens (CM) (or 'mini-condenser' or 'twin' lens) and the objective pre-field (OL) as shown in figure 1(b). Large convergence angles, typically used for high spatial resolution scanning transmission electron microscope (STEM) imaging are obtained using a large condenser aperture and the mini-condenser lens very weakly excited or switched off. In this case, the effect of CL3 and OL is to create a highly focused probe at the specimen plane (red ray path in figure 1(b)). A small condenser aperture can also be employed in this mode to reduce the convergence angle while maintaining the same lens excitations as for STEM (orange ray path in figure 1(b)). A particularly advantageous condition is achieved when the condenser aperture is small enough to be filled with a flat electron phase profile. In this case, the effect of aberrations is minimized and the diameter of the final electron probe is diffraction limited [33]. The diffraction pattern then consists of discs that are the image of the condenser aperture. The CM lens can be used to change the convergence angle of the final electron probe without changing the condenser aperture. The CL3 and CM lenses can be excited to form a focused probe with a smaller convergence angle and larger diameter than for normal STEM imaging (dashed dark blue ray path in figure 1(b)). In this situation, an additional cross-over is formed between CL3 and CM. In the extreme case, CM can be strongly excited to form a cross-over in the back focal plane of OL and a parallel probe beam [33,35] is formed at the specimen (dashed light blue ray path in figure 1(b)) to achieve the best angular resolution in the diffraction pattern. The presence of a probe spherical aberration corrector (CsC as shown with dashed blue box in figure 1 gives an extra degree of flexibility. In normal operation CsC is optically neutral, but when switched off, greater de-magnification of the probe beam is achieved, allowing parallel beam operation with a small diameter, high-acuity probe [35]. With such freedom of choice over probe beam convergence angle and diameter, the experimenter must tune these parameters to access desired correlation lengths and minimise the ambiguity of data interpretation. Obviously, the beam size must be large enough to encompass the structural correlations that are to be probed, but beyond this, the optimum choice depends on the analysis. For example, in the case of intensity variance, changing the probe size has been shown to be an effective way to estimate extended correlation lengths [19,36] as the variance decreases if two or more structurally correlated regions are included in the probed volume. When angular correlations are quantified using C n q , the ideal probe size (where C n q from structurally correlated regions dominated), equalled the correlation length of interest [13,37]. Even when employing an optimal probe size, in dense, disordered systems, where multiple short-range structures may be included within the probed volume, particles that are not structurally correlated may give rise to serendipitous, strong angular correlations [13,22,37]. This situation is shown in figure 1(c) for the case of a 4 nm thick Cu-Zr metallic glass specimen, where the correlation length decays rapidly after three co-ordination shells. Both structurally correlated (in red-from one coordination polyhedron) and uncorrelated atoms (in green-distributed in the specimen thickness) can have a transverse arrangement to give a strong six-fold symmetry. This 'structural noise' will add a spurious signal to the C n q as we demonstrate below in section 2.1. 3. In experiments where the magnitude of angular symmetries is quantified with C n q , non-zero odd symmetries have been noted [25,30,38] demonstrating an absence of Friedel symmetry. Friedel symmetry is the exact inversion symmetry in the diffraction pattern (I(q, ϕ) = I(q, ϕ + π)) that arises in kinematical diffraction conditions where the diffraction pattern is the magnitude squared of the Fourier transform of the projected structure. Friedel's law can be broken due to Ewald sphere curvature [30] in the case of specimens with an inverse thickness that is comparable to the curvature of the Ewald sphere at the q-range being considered [39]. Lens aberrations can also cause this lack of π-azimuthal symmetry [16,17]. Both of these contributions were shown to be negligible in the context of examining thin (≲ 4 nm at 300 keV) specimens in an aberration-corrected S/TEM [37,38]. Rather, this breakdown in π-rotational symmetry was found to be caused by dynamical diffraction from a non-centrosymmetric object, resulting in significant ratios of odd-n to even-n C n q magnitudes [37,38]. Thus even in this weak scattering regime, at the lower limit of what can be experimentally realized, the diffraction pattern of a disordered material is not related by Fourier transform to the projected structure of the object. Figure 2 demonstrates this phenomenon for the 4 nm thick metallic glass shown in figure 1(c). Here the atomic glass model was simulated using molecular dynamics [40] and the electron diffraction was simulated using well-known algorithms [41] implemented in Interactive Data Language [38]. The figure displays the simulated END pattern, the angular auto-correlation function, C(q, ∆), the C n q and the average C n q from the first diffraction peak, C n , for the kinematical case (figure 2(a-d)) and for the dynamical case (figure 2(e-h)), modelled using multislice simulations (slice thickness 1 Å).

Noise
The angular symmetries in electron nano-diffraction patterns from amorphous materials are often subtle and interpretation requires knowledge of the contribution of sources of noise and uncertainty and how these propagate though the analysis. The following analysis examines this for the angular auto-correlation function (equation (4)) and other parameters derived from this [12,19,24]. The angular auto-correlation function is calculated from a set of i intensities at different azimuthal angles I(q, ϕ) = I i . These intensities can be approximated as a mean intensity I 0 , divergences from the mean due to diffracted intensities from structurally correlated and uncorrelated atoms within the probed volume, (i) S i , and (ii) R i , respectively, and (iii) white experimental noise, n i , which may be related to the detector or instrument stability. Here, The angular auto-correlation C(q, ∆ k ) for a single angular separation ∆ k , is calculated from the intensity at (q, ϕ) and at (q, ϕ + ∆ k ) which we will denote I i and I i + k , respectively. In the discrete case is the number of pixels in the diffraction pattern at a radius of q. Thus we can calculate the components of the angular auto-correlation (equation (4)).
Here the separate angular auto-correlations in ⟨⟩ ϕ can be positive or negative and give rise to systematic correlations or anti-correlations at different ∆ k . When the angular auto-correlation is transformed into Fourier coefficients these contributions are additive, and will not average to zero if we average the C n q over a large number of diffraction patterns. Assuming the experimental additive noise for these stochastic events can be modelled as white noise, we have ⟨n i n i+k ⟩ ϕ = δ k ⟨n 2 i ⟩ ϕ , such that the angular auto-correlation of the noise at ∆ k = 0 matches its variance (δ k is a delta function). This noise spike will contribute to C 1 q . For weak scattering, the structurally correlated and randomly correlated signals from atomic scattering will be proportional to the input intensity. The unchanging additive noise variance then makes comparatively diminishing contributions as the flux is increased.
In measurements where angular correlations are compared to favoured atomic motifs [42] or C n q averaged over a given q-range are mapped [25,27,28] extreme caution must be exercised as 'structural noise' can give rise to spurious symmetries that do not correspond to any real correlated structural motif in the material [22,37]. In measurements where the C n are averaged from an ensemble of patterns and plotted as a function of n to obtain a fingerprint of local order [26,31], the 'structural noise' and the experimental noise add a constant offset across n with the correlated structural fingerprint superimposed. In all the measurements employing C n q , care must be exercised to use the same experimental geometry and exposure time to ensure that data sets are comparable. In measurements where the angular auto-correlation function itself is averaged from a number of diffraction patterns [12,22,24], the noise terms will eventually tend to zero (except the spike at ∆ k = 0) and significant hidden signatures of order will emerge in the final function, even if individual diffraction patterns are too noisy to clearly see angular correlations by eye.

Experimental details
The data presented in this article is drawn from many sources to illustrate the power and universality of employing angular correlations in limited-volume diffraction patterns to probe local structure and symmetry in disordered materials, but all these analytical methods can readily be applied to electron nanodiffraction. SEND measurements from Cu-Zr metallic glasses were obtained in a FEI Titan 3 80-300 FEGTEM operated at 300 keV. A 10 µm condenser aperture produced a diffraction-limited probe with 2.7 mrad convergence angle and 0.51 nm full-width-half-maximum. Diffraction patterns were collected using a Gatan UltraScan CCD camera (2048 × 2048) and binned to 512 × 512 [25]. For each specimen 2700 patterns were obtained from three separate 3 nm × 3 nm regions scanned with 0.1 nm steps. The thickness of the probed areas was kept constant and measured using electron energy loss spectroscopy to be ≈ 5 nm [25].
Scanning micro-small-angle x-ray diffraction (µSAXS) was performed on colloidal glasses composed of 300 nm diameter silica spheres using a 5 keV x-ray beam limited in size by a 20 µm near-field aperture and a Further details in [25]. b Further details in [26].
Dectris Pilatus 1 M diffraction camera [26]. Monodisperse colloidal dispersions have a simple crystalline phase diagram, featuring hexagonal close-packed (HCP), face-centred cubic (FCC), and body-centred cubic (BCC) crystals, and this model system has proved an excellent testbed for novel characterization methods [26]. In the case of both the metallic and colloidal glasses, the specimen thickness was limited to contain 10-20 nearest-neighbour structures in the probed depth to simplify the interpretation of the angular correlations [25,26]. Table 1 summarizes the size and resolution of the diffraction patterns employed in the analyses presented in section 3.

Transforming angular correlations in diffraction patterns to real-space functions
An alternative approach to obtaining structural information from angular correlation functions is to calculate the real-space pair-angle distribution function (PADF) [24] from the averaged angular crosscorrelation function. The PADF contains statistical information about three-and four-atom arrangements in a material and is a natural generalization of the pair-distribution analysis (PDF) that is commonly performed on disordered materials. To provide an explicit form of the PADF function, we first define the n-body correlation function as: where N is the total number of atoms in the sample, r is an atomic coordinate and the average ⟨⟩ is an ensemble average over many nanoscale regions in a bulk disordered sample (or structural states in a thermodynamic ensemble). The PADF is defined as: The PADF depends on the separation distance of two pairs of atoms (r and r ′ ) and the relative angle between those atomic pairs. The three-body terms arise when the two atom pairs share a common atom and the two-body term arises when the atom pairs are identical. The tilde symbol˜indicates that each n-body correlation function is modified by integrating out degrees of freedom to which the measurement is insensitive, namely the absolute position, the absolute orientation of the sample, and the distance between the pairs of atoms. The presence of two three-body terms generates an angular symmetry around 90 degrees. The four-body term obeys the same symmetry. When r and r ′ are equal to a bond distance, then Θ(r, r ′ , θ) provides information about the bond angle distribution. At larger distances, the PADF contains angular peaks whose locations are determined by common local structures in the material. An analysis of angular peak positions can be sufficient to detect the presence of local structures in the material. Simulations studies on amorphous and crystalline nickel [24] and an experimental study of gold nanocrystals [43] has demonstrated that three-body combinations can account for most peaks at short distances. In a study of activated carbon [43], five-fold ring defects and adatom defects that were not possible to distinguish from PDF analysis were identified by their angular peaks in the PADF.
For this article, the potential of the PADF was investigated to distinguish local order in colloidal glasses with simple archetypal local structures. One hundred µSAXS patterns of different local structures of 300 nm diameter silica spheres (icosahedral, HCP, FCC, BCC) were simulated in random orientations employing the kinematical approximation and without modeling absorption. The angular cross-correlation function was calculated for each pattern according to equation (3) and averaged from the ensemble of patterns. The angular correlations in C ′ (q, q ′ , ∆) can be transformed into the three-dimensional (3D) pair-angle distribution function (PADF-Θ(r, r ′ , θ)) employing a spherical harmonic expansion. This is analogous to the Fourier decomposition given in equation (4) except that the spherical harmonics invoke a 3D geometry. The curvature of the Ewald sphere is accounted for in the sampling of the harmonic expansion. The radial coordinates q and q ′ are transformed into real space variables r and r ′ with spherical Bessel transformations. For the full mathematical details see reference [24].

Modelling methods
The atomic structure of disordered solids is challenging to model with ab-initio or quantum-chemistry methods, since significant numbers of atoms (and electrons) are required to capture the structural heterogeneity. Density functional theory can efficiently simulate glassy structure for small homogeneous solids [44][45][46][47] but has failed to describe some important archetypal structures, including large continuous random networks [48]. In contrast, classical force-field approaches based upon Newtonian molecular dynamics or statistical-mechanical Monte Carlo algorithms can accurately approximate structure and properties (forgoing electronic structure) for many thousands of atoms. Mimetic modelling methods, which incorporate material synthesis pathways, can successfully predict structure, for example, in the case of stress induced transitions from graphitic to diamond-like carbon in thin film synthesis [49]. However, such approaches require intimate knowledge of material history, without which predictions can fail.
When the material synthesis history cannot be readily modelled, an alternative is to consider atomic structures that are statistically consistent with experimental observations like quantitative electron diffraction [3,6,50,51]. One such approach, the 'reverse Monte Carlo' (RMC) algorithm [52], seeks to minimize differences between model and experimental data [53], including pair correlations governed by the radial distribution function g(r), the structure factor S(q) and geometric constraints for the models, such as coordination numbers and average bond angles. Despite the utility of intuitive geometric constraints and topological controls, for example ring statistics, it is well known that experimental g(r) and S(q) data cannot reproduce important many-body correlations, and so an overly large ensemble of different atomic structures can be non-uniquely consistent with the constraints. It is for this reason that the hybrid reverse Monte Carlo algorithm (HRMC) was developed to handle disordered covalently bonded solids [54][55][56][57]. HRMC is able to model non-equilibrium structures with minimised energy configurations, using algorithm protocols such as forms of simulated annealing.
There have been a few examples of orientational order being incorporated into RMC models by employing experimental constraints. In the technique called 'fluctuation electron microscopy' the second moment in the intensity distribution, the variance, is calculated, which is a function of the 2-, 3-, and 4-body correlations [20]. Energy-optimised RMC simulations employing the variance have shown that the explicit constraint by higher-order correlations via the intensity variance reduces the region of solution-space for amorphous silicon to those that contain paracrystals, but that a suitable potential is also invaluable to guide the refinement [58,59]. Our HRMC simulations have shown that electron phase contrast images from aberration corrected electron microscopy could also provide additional constraints for paracrystals buried in an amorphous network [60]. Constraints based on high-resolution TEM images have been implemented for some time [61] and are useful for modelling extensive orientational order, such as graphitic stacking, basal plane curvature and other complexities in chars and pyrolysed carbons [62] [63].
Recent advances in related END and µSAXS experiments and theory has simplified the interpretation of these richly complicated data sets, leading to new measures of many-body correlations [24,25,31,38], which ought to be of great use to RMC and HRMC simulations. The direct experimental fitting of angular correlations to infer statistical proportions of locally ordered clusters by Liu et al [26] can be used in RMC to improve agreement with known measurements and create disordered structures littered with archetypal structural motifs, which compete to describe the overall BOO [32]. Although computationally costly, one can in principle forward-calculate four dimensional scanned diffraction data to fit angular symmetries, projected BOO patterns, fluctuation microscopy etc However, the BOO basis measurements of Liu et al [26] offer simpler avenues of thought. Here we show a first such application to RMC, whereby proportions of clusters ascribed to the BOO parameters are inserted into an amorphous network which is originally consistent with volume-averaged diffraction experiments, yet lacks known angular symmetries and distinct packing geometries. RMC relaxation is then performed, whilst maintaining many body correlations arising from the inserted clusters. The final model possesses the average orientational order derived from measurements, but post-relaxation, shows how this order can be incorporated in 3D into the amorphous structure.
Starting from a random configuration of 2000 spheres of radius 140 nm, RMC calculations simultaneously fitted the µSAXS g(r) and S(q) pair correlation data for 4 million steps, over which constraint weightings were lowered in four stages and step sizes were altered to maintain an efficient acceptance ratio. Random voids were then carved out of the amorphous structure, to make room for the inclusion of densely-packed archetypal FCC, BCC, and ICO clusters comprising 13, 15, and 13 spheres, respectively, in the proportions measured from experiment. The µSAXS BOO parameter fitting of Liu et al [26] indicated just 30% of essentially random bond orientation order and so 70% of spheres were removed from the amorphous structure. The network of spherical voids was created using a template of disordered spheres. To this end, Monte Carlo moves were used to melt a close-packed FCC lattice of 380 nm radii spheres, to circumvent random-sphere-insertion stagnation for such high packing fractions. Allowance for a slight amount of sphere void overlap was found necessary to achieve this, reducing the hard core minimum separation by 7%.
The BOO clusters with equal 280 nm bond lengths were then randomly inserted into the voids without overlap of the constituent ordered spheres, with the remaining experimentally fitted structure outside the voids. The bond angle coordination number distribution functions were calculated from this structure of jammed BOO archetypes as target constraints for a final RMC simulation, to maintain the characteristic many body correlations. Over 8 million steps, these and the experimental pair correlation constraints were applied to construct a statistical ensemble of models that agrees with the diffraction-based experimental data while retaining the many body angular and coordination information arising from the original BOO clusters.
Reverse Monte Carlo (RMC) [52] modelling was performed using the hybrid RMC algorithm code (HRMC version 2.0), as downloaded from the Computer Physics Communications database [57]. The calculations described in section 3.3 were designed to test the merits of including high order correlation information obtained from scanning µSAXS experiments [26] of silica spheres. As such, energy minimisation (pertaining to the 'hybrid aspect' of the code) was switched off, as were all other extraneous constraints, such as three-membered ring penalties. The sum square differences between diffraction data were minimised over several million steps, including that of the inferred bond-angle and coordination distributions (see section 3.3). This was achieved using the simulated annealing procedure available in HRMC 2.0, whereby the constraint magnitudes and Monte Carlo step sizes can be gradually decreased in four stages, to maintain an efficient probabilistic acceptance ratio for rapid minimisation of the global error.

Data-intensive analysis
Statistical analysis of many spatial samples from employing a scanned beam can give unique insights into local order and structural variability in amorphous materials. Here, we employ the Fourier coefficients C n q averaged over the first sharp diffraction peak C n as a set of order parameters that probe the arrangements in the nearest-neighbour polyhedra in metallic glasses in the sampled column of material and present statistical analyses of these to examine differences in local order and variability. The advantage of such measures is that they allow the data to 'speak' for itself, and do not require comparison with models or assumed structural motifs. The statistical measures employed were the mean of m spatial samples,C n , the covariance in the Fourier coefficients Cov(n, n ′ ) and Pearson correlation coefficient between the C n from the m scanned-points, Corr(m, m ′ ). We focus on 1 ≤ n ≤ 12 as although n in principle goes up to half the size of the image, in practise physically meaningful n from projected short-range structures finish at 12. Spatially averaged order parameters are denoted by: where elements C n 1 , C n 2 , · · · specify the order parameters at different spatial locations. Using this notation, the estimated covariance and variance among symmetry magnitudes are described by standard forms.

Cov(n, n
It is also instructive to measure the correlation (normalized covariance) between scanning probe location indices,  represents an average over a maximum of 12 symmetry indices for scanning index m. Figure 3 displays the average Fourier coefficients,C n and the covariance matrices, from Cu-rich and Zr-rich, Cu-Zr metallic glasses for 1 ≤ n ≤ 12. Previous work on the model system of monodisperse colloidal glasses used the idealC n calculated from the archetypal short-range structures from the liquid and the two-component crystalline phase diagram to decompose theC n into contributions from different kinds of order [26,31]. For the Cu-Zr system, the phase diagram in the glass-forming region contains many more crystalline phases and this decomposition is challenging. However, from inspection, the Cu-rich glass has much more persistent magnitudes in the higher symmetries, which is a strong indicator for local structures with higher point symmetry such as the icosahedral cluster [26,31]. The covariance matrices also support this conclusion. They show how one C n varies relative to another averaged over all spatial samples, and thus a high covariance between certain n could evidence a certain symmetric local structure. The covariance from the Cu-rich glass shows greater variance (diagonal elements) out to higher n and also larger covariance between even-n C n , which shows that the Zr-rich glass has less well-defined local order, with perhaps greater elements of disorder.
In figure 4 we display the matrices of the Pearson correlation coefficients for the Cu-rich and Zr-rich metallic glasses, noting that the normalization in equation (13) ensures values span -1 to 1 for perfect anti-correlation and perfect correlation, respectively. The diagonal elements in the m = m ′ line in the correlation coefficient matrices have a value of 1 while the bright and dark stripe features in the matrices arise from regions with particularly strong or weak diffraction. From the correlation coefficient matrices ( figure 4(a) and (b)) and the plots of their distributions (figure 4(c)) we see that the Cu-rich glass shows greater variation from place-to-place. This is an indicator that the Cu-rich glass has either more distinct or more varied local structures. Such an observation may give valuable insight into glass-forming ability, a key quantity in glass science. Theorists have demonstrated that larger numbers of more distinct local structures suppress the fluctuations in order in the under-cooled melt that may give rise to crystallization [64]. Here we see that the Cu-rich glass (the better glass-former) has larger variability in local structure, which may reflect the population of local structures frozen in from the melt at the glass transition. It is worth noting that the Pearson distributions for each of the three individual regions were consistent with the combined histograms in figure 4(c) and distinct between the two glasses. For example, the Cu-rich region means and medians were (0.42, 0.49, 0.48), and (0.44, 0.51, 0.40), respectively (cf. a mean of 0.42 and median of 0.44 in the combined data). Similarly, the Zr-rich means and medians were (0.63, 0.57, 0.59), and (0.66, 0.59, 0.62), respectively, consistent with the mean of 0.59 and median 0.62 in the combined data.
A final statistical measurement we propose is principle components analysis (PCA) that is widely used to reduce the dimensionality of a large dataset with many correlated variables, while retaining the largest possible proportion of variation present in the data. PCA has been widely employed in the electron microscopy community for both scanning spectral and diffraction-based datasets to identify principle components that can be mapped [65][66][67]. The PCA decomposition performed in this work is done via singular value decomposition using the implementation in HyperSpy [68,69]. Decompositions were performed on combined (from three 30 × 30 scans) datasets, where the central beam was masked to focus on the variation in scattered intensity, and the mean was subtracted. Prior to decomposition, all the data was cropped to have exactly the same dimensions with the direct beam at the center of the patterns. This step was needed to counteract slight difference of beam center across the different scans and samples. The final datasets were of the dimensions 3 × 30 × 30 × 452 × 452.
The first twelve components for the two metallic glasses are shown in figure 5 (a) and (b), respectively. While it is extremely challenging to interpret individual features of the components for amorphous materials as has been noted in section 2.1.2, interesting features can be observed. For example, there is a tendency towards more structured components in the Cu-rich sample, thus supporting the notion of a larger number of more symmetrical clusters in this specimen. Future work implementing non-negative matrix factorization [67,[70][71][72], which enforces a non-negativity constraint on the components, or other decomposition algorithms, which previously have successfully been applied to crystalline solids, holds the potential to provide more interpretable decomposition factors and loadings.
A major part of applying PCA in practice is to estimate the number of components needed to describe the data satisfactorily. Common practice is to inspect either the explained variance ratio or the cumulative explained variance ratio as a function of components in the search of a 'knee-point' [73,74]. The cumulative explained variance ratio plot quantifies to what extent the total variance is contained within the corresponding number of components. This plot indirectly may help one assess the level of redundancy present in the dataset, which by extension may be correlated to the structural variability within a given sample. This curve, and the closely related Scree plot, have successfully been used to determine the number of spectral components needed to describe S/TEM spectroscopy datasets and precession diffraction from simple crystalline materials [65,67,75]. In contrast, the 'knee-point' is not particularly prominent in the cumulative explained variance ratio from PCA of the full diffraction patterns of the metallic glasses, evidencing the complexity of interpreting the components that are found. However, in common with the other statistical analyses presented here, a clear difference between the glasses is detected with the Cu-rich glass requiring more components to describe a given variance than the Zr-rich glass. This further supports the hypothesis that the better glass-forming material will possess a richer and more varied set of local configurations frozen into the structure at the glass transition.

Pair-angle distribution function (PADF)
Simulated examples of PADF distributions for different nearest-neighbor colloidal clusters are shown in figure 6 and clearly demonstrate the power of the PADF to 'fingerprint' key local structures. Icosahedral, hexagonal close packed (HCP), face-centred cubic (FCC), and body-centred cubic (BCC) local structures taken from reference [26] were used to simulate a hundred x-ray diffraction patterns each in random orientations. The ensembles of diffraction patterns were used to generate angular intensity correlation functions according to equation (3) for each local structure. The intensity correlation functions were converted into the PADF distributions. Both the intensity correlation functions and the PADF distributions show distinctive angular peaks for each local structure. This portrays the power of angular correlations to provide a fingerprint of different local structures both in q-space and real space. However, the PADF provides extra advantage in that we can directly identify the three-and four-atom arrangements that contribute to each angular peak, as shown by the models in figure 6. The simulated structures are all unit polyhedra scaled to the size of the sphere, and have the same nearest-neighbour inter-particle spacing, giving rise to similar short-range order in the PDF. When you include the angular information from the PADF the structures are much more distinct, showing how angular correlations provide a clear and direct method of identifying dominant local structures in amorphous materials. PADF analysis is based on an assumption of kinematic scattering. In PDF analysis, multiple scattering changes peak amplitudes but not the peak positions [76]. The first experimental PADF results with electrons [43] indicated that is also true for the PADF. Further understanding of how multiple scattering impacts PADF analysis will be required for quantitative studies of PADF peak amplitudes measured with electrons.

Convergence and noise in PADF analysis
A key parameter for a successful PADF measurement is the number of diffraction patterns required for convergence to the ensemble average. Convergence is achieved when there is sufficient averaging of: (i) uniform, random orientations of the local environments; (ii) variations of the local structures; and (iii) measurement noise.
Orientational averaging necessitates the measurement of multiple diffraction patterns even if the sample is comprised of well-ordered, single-crystal nanoparticles. This is because the PADF function defined in equation (8) is sensitive to local angular structure, but not the global orientation of the sample. The number of diffraction patterns required for randomly sampled reproducible structures is known from the theory of single-particle tomography using serially collected data [77]. The formula for the number of particles required depends on the size of the structure of interest and the target resolution. For highly disordered samples that contain only short-range order, the number of patterns will be of the order of 100 to 1000, which has been confirmed in PADF simulation [24] and experiments [43]. For materials with mesoscale ordering such as activated carbons small, but visible, effects of incomplete convergence were observed across different subsamples of 900 diffraction patterns. Employing datasets of 2700 diffraction patterns in total, the PADF with distances from 1.5 Å to 3 Å [43] and an angular resolution of 12 degrees could be obtained. To access larger length scales and higher angular resolution, more data may be required. In summary, the early evidence is that for many samples of interest, the number of electron diffraction patterns required should be of the order of 10 3 -10 4 to observe the PADF statistics of the short-range order at atomic resolution.
Measurement noise is primarily an issue when the scattering is very weak, the beam size is much larger than the local structures of interest, or the specimen is extremely electron dose intolerant. Electron measurements are not greatly impacted by either of the first two problems. The elastic atomic cross-section for fast electrons is sufficiently large for angular intensity variation due to structure to dominate measurement noise. In other words, tens or hundreds of electrons can be measured per pixel in a single exposure achieving a positive signal-to-noise ratio up to high resolution. Electrons are readily focused to sub-nanometer length scales, so the beam focal size can be optimally matched to length scale of interest. Figure 7 show the results of incorporating measured bond-orientational order from angular correlations in diffraction patterns into RMC models. Part (a) of the figure shows the inserted BOO clusters in purple amidst the amorphous configuration in green, prior to the final RMC modelling. Note that all clusters were unrealistically aligned in the same direction, to visually highlight the subtle inclusion of the BOO clusters, however all statistical measures of the short range order reported here are isotropic. Figure 7 (b) shows one configuration after those 8 million final steps. One cannot discern from this view (or other visual sections and projections), if the BOO cluster information has been preserved-many-body statistical measures are required for that assessment, as follows.

Incorporating measured bond-orientational order in reverse Monte Carlo modelling
Satisfactory agreement between experiment and simulation was maintained between RMC runs, as shown in figures 7 (c) and (d), however it is not clear whether the final RMC calculation has removed all memory of the BOO cluster insertion. The retainment of this many-body experimental information is Figure 6. Simulated PADF calculations. Four local structures were studied: icosahedral, hexagonal close packed (HCP), body-centred cubic (BCC) and face-centred cubic (FCC). For each local structure one hundred diffraction patterns were simulated with the structure in random orientations. The top row shows a single diffraction pattern for each structure. The second row show a cross section of the intensity correlation function calculated from the ensemble of diffraction patterns. The third row shows a cross section r = r ′ of the PADF for each structure as computed from the intensity correlation function. Structure diagrams show a selection of 3-and 4-atom combinations that contribute to angular peaks in the PADF plots at the locations annotated with white dots. Grey dots indicate the symmetry equivalent positions. however demonstrated in figures 7 (e) and 7 (f). Prior to BOO insertion, the broad bond angle distribution shown by the black dots in figure 7 (e) portrays the highly disordered nature of the structure, with a peak at 60 • primarily associated with a high fraction of three membered rings that are created to efficiently fit the 1st nearest neighbour peak in g(r), whether physically realistic or not [54]. These findings are typical of RMC simulations without geometric or many-body constraints, as the error minimisation tends to maximise entropy due to the large degeneracy of structures consistent with the diffraction data. As figure 7(e) and (to some extent) figure 7(f) show, the inclusion of the BOO clusters, with subsequent bond-angle and coordination fitting, leads to structures with considerably more short range ordering.
As another measure of the BOO cluster information preservation, we can calculate the quadratic and third-order rotationally invariant spherical harmonic orientational order parameters, Q l and W l that are sensitive to the local symmetry in particle arrangements around a central particle [32]. These classic BOO parameters were calculated from one of the RMC structures in the statistical ensemble (after the final fitting), as shown in figure 8. In particular, the averaged spherical harmonic bond orientation order parameter of figure 8(b) exhibits increased angular correlations, indicating that the BOO cluster information has been retained. The retention of this order shows that a relaxed amorphous structure can be formed with 70% of   [32] in (a) shows little difference between the local order before and after cluster insertion, whereas the associated W l BOO parameter in (b) shows a marked increase in angular order. The error bars in each graph represent the standard error of the mean. particles in a highly-defined local order validating the analysis method and interpretation of the original measurement [26].
Similar approaches to match many body correlations between experiment and simulation can be conceived beyond this preliminary analysis. Prospects for modelling disordered structure at the atomic scale are particularly enticing, given the recent nano-beam diffraction advances described in section 2.1.1. Likewise, the PADF presented in section 3.2 disentangles many body correlations inherent within fluctuation microscopy and scanned diffraction data. The PADF contains angular short-range order information that can be calculated more efficiently than diffraction intensity variance or other medium-range order measures, which ought to be useful for refining structural models. Indeed, a recent RMC uniqueness study has shown the efficacy of additionally constraining known bond and coordination distributions for reproducing amorphous network statistics with desirably-low density functional energies [78].

Conclusions
Advances in experiment, analysis, interpretation and modelling using angular symmetries are opening new avenues for research into the structure of disordered solids, some of which have been outlined and demonstrated here. Experimental geometry can be optimized to facilitate the interpretation of the angular symmetries in limited-volume diffraction patterns from amorphous materials. Nevertheless, the interpretation of these symmetries is complicated for dense, disordered materials and the effects of specimen thickness, beam aberrations, dynamical diffraction and 'structural' and experimental noise need to be considered for each analysis. Conclusions about structure from single, rare, observations of a particular angular symmetry should be avoided, and statistical measures employed, as is appropriate for these materials that possess a spectrum of local particle arrangements.
Subtle angular symmetries can be uncovered through calculating the angular crosscorrelation and autocorrelation from the diffracted intensity. The magnitude of the symmetries can be further quantified by calculating the Fourier coefficients of the angular autocorrelation. If these are averaged over a q-range corresponding to the nearest-neighbour shell and many spatial samples, they provide a strong fingerprint of local bond-orientational order for metallic and colloidal glasses. For simple systems, these can be broken into components from different proportions of candidate local structures. Further statistical measures of these Fourier coefficients such as covariance, correlation and principal components analysis showed systematic trends for a Zr-Cu metallic glass, promising tantalizing new insights into the number and variability of local structures in amorphous materials. We note that these approaches offer an attractive data-centric approach that does not assume a favoured structural motif or motifs.
Analogous to the procedure for inverting the average diffracted intensity to the pair-correlation function, the angular cross-correlation function that additionally encodes the three and four-body correlations can be inverted to the real-space pair-angle correlation function. This PADF is a rich fingerprint of the local structure and often contains features that can be directly interpreted to arise from a given arrangement of particles in the first and second nearest-neighbor shells, as demonstrated for hard-sphere glasses.
These angular symmetries and the other parameters derived from them (Fourier coefficients, PADF, proportions of local structures) can be employed to efficiently constrain models created using reverse Monte Carlo and address some of the longstanding challenges in this field, such as non-uniqueness and tendency towards heightened disorder. Accurate modelling can supplement the structural interpretation of these new statistical measures derived from diffraction and harness or test the validity of structural interpretations, as we show here.