Investigating the missing-wedge problem in small-angle X-ray scattering tensor tomography across real and reciprocal space

Nielsen, L.C.; Tänzer, T.; Rodriguez-Fernandez, I.; Erhart, P.; Liebi, M.

doi:10.1107/S1600577524006702

research papers

JOURNAL OF
SYNCHROTRON
RADIATION

ISSN: 1600-5775

Volume 31| Part 5| September 2024| Pages 1327-1339

https://doi.org/10.1107/S1600577524006702

Open

access

Investigating the missing-wedge problem in small-angle X-ray scattering tensor tomography across real and reciprocal space

Leonard C. Nielsen,^a Torne Tänzer,^b,^c Irene Rodriguez-Fernandez,^b Paul Erhart ^a and Marianne Liebi ^a,^b,^c ^*

^aDepartment of Physics, Chalmers University of Technology, Gothenburg, Sweden, ^bPhoton Science Division, Paul Scherrer Institute (PSI), Villigen, Switzerland, and ^cInstitute of Materials, École Polytechnique Fédérale de Lausanne (EPFL), Lausanne, Switzerland
^*Correspondence e-mail: marianne.liebi@psi.ch

Edited by S. D. Kelly, Advanced Photon Source, USA (Received 22 April 2024; accepted 8 July 2024; online 28 August 2024)

Small-angle-scattering tensor tomography is a technique for studying anisotropic nanostructures of millimetre-sized samples in a volume-resolved manner. It requires the acquisition of data through repeated tomographic rotations about an axis which is subjected to a series of tilts. The tilt that can be achieved with a typical setup is geometrically constrained, which leads to limits in the set of directions from which the different parts of the reciprocal space map can be probed. Here, we characterize the impact of this limitation on reconstructions in terms of the missing wedge problem of tomography, by treating the problem of tensor tomography as the reconstruction of a three-dimensional field of functions on the unit sphere, represented by a grid of Gaussian radial basis functions. We then devise an acquisition scheme to obtain complete data by remounting the sample, which we apply to a sample of human trabecular bone. Performing tensor tomographic reconstructions of limited data sets as well as the complete data set, we further investigate and validate the missing wedge problem by investigating reconstruction errors due to data incompleteness across both real and reciprocal space. Finally, we carry out an analysis of orientations and derived scalar quantities, to quantify the impact of this missing wedge problem on a typical tensor tomographic analysis. We conclude that the effects of data incompleteness are consistent with the predicted impact of the missing wedge problem, and that the impact on tensor tomographic analysis is appreciable but limited, especially if precautions are taken. In particular, there is only limited impact on the means and relative anisotropies of the reconstructed reciprocal space maps.

Keywords: tensor tomography; small-angle scattering; nanostructure; missing wedge problem; optimization.

1. Introduction

Small-angle X-ray scattering tensor tomography (SAXSTT) is a promising method for probing anisotropic nanostructures of macroscopic samples in a volume-resolved manner (Liebi et al., 2015 , 2018 ; Schaff et al., 2015 ; Guizar-Sicairos et al., 2020 ). It has been applied to the study of a variety of biological materials, including bone, tendon and myelin (Georgiadis et al., 2021 ; Casanova et al., 2023 ; Grünewald et al., 2023 ; Silva Barreto et al., 2024 ). In the absence of very strong real-space uniformity and reciprocal space symmetry constraints (Stribeck et al., 2006 ; Skjønsfjell et al., 2016 ), SAXSTT requires a more general acquisition scheme than traditional scalar tomography, such as rotating the sample while subjecting the axis of rotation to a series of tilts (Schaff et al., 2015; Liebi et al., 2015, 2018), carrying out measurements over part of a sphere of rotation. Such acquisition schemes are generally geometrically constrained to a tilt of up to 45°, since the rotation stage will obstruct the beam at greater tilt angles. Nielsen et al. (2023 ), using this measurement scheme with simulated data, observed that the degree of correlation with the original reciprocal space maps (RSMs) approached lower values than what should be theoretically attainable in terms of the RSM representations used in the simulations and reconstructions, even at very low noise levels. This can likely be attributed in part to the so-called missing wedge problem, a common data incompleteness problem in tomography (e.g. Liu et al., 2018 ). While a limited investigation into the effect of reduced data was carried out by Liebi et al. (2018), a thorough examination of the effects of data incompleteness is still outstanding. A deeper understanding of data incompleteness in SAXSTT is a crucial component of the development of approaches to counteract this incompleteness, similar to those used in other tomography methods, as in, for example, Trampert et al. (2018 ), Ding et al. (2019 ) and Moebel & Kervrann (2020 ).

Here, to investigate these effects under real experimental conditions, we present a scheme utilizing sample remounting to yield two incomplete data sets, each measured using the 0°–45° tilt scheme, which when combined form a complete data set. The scheme was applied in measurements on a sample of human trabecular bone. For the reconstruction, we employed an RSM representation which uses local Gaussian radial basis functions on a spherical grid to interpolate measured data. This reconstruction is closely related to the spherical integral geometric tensor tomograph (SIGTT) approach (Nielsen et al., 2023) but replaces the model for the RSM with local functions, which avoids artifacts due to the spherical harmonic Gibbs phenomenon (Gelb, 1997 ). Both models have in common that the only symmetry enforced is Friedel symmetry, and thus allow for reconstruction of complex textures. In addition, the use of local radial basis functions permits the problem of SAXSTT to be analyzed as a set of scalar tomography problems, which allows the application of the framework of standard tomographic analysis. By carrying out reconstructions from the two separate data sets, as well as the combined data set, and comparing the reconstructions, this work seeks to investigate whether imperfect reconstructions in limited-angle small-angle X-ray scattering tensor tomography can indeed be attributed to missing wedges. In addition, we aim to provide insight into the impact of this effect on SAXSTT analysis. We conclude that the differences between complete and partial data sets are consistent with the predicted effects of the missing wedge problem. These effects impact typical SAXSTT analysis in a non-trivial but manageable way, and suggest strategies for mitigation. In particular, two important scalar quantities, the means and relative anisotropies of the reconstructed RSMs, show relatively little impact from the missing wedges.

2. Theory

The RSM measured by small-angle X-ray scattering (SAXS) in a small volume may be written as

$[{\rm{RSM}}({\bf{q}},{\bf{r}}) = \int\!\!\int\!\!\int {\rm{d}}V\ \left\{ \tilde{\rho}\,\left({\bf{r}}^{\prime}-{\bf{r}}\right) \exp\left[-i{\bf{q}}\cdot\left({\bf{r}}^{\prime}-{\bf{r}}\right)\right] \right\}, \eqno(1)]$

where $[\tilde{\rho}\,({\bf{r}}^{\prime}-{\bf{r}})]$ is the auto-correlation function of the electron density over the small volume, r is the position of the center of the volume, r′ is the point of integration within the volume and q is the reciprocal space vector. For a small scattering angle, such that the scattered intensity travels approximately the same path as the transmitted intensity, and assuming that the total amount of scattering is small enough to not significantly influence the transmission, we can probe the RSM by measuring the small-angle scattering intensity with a small beam, and correcting by the transmitted intensity as

$[\int\limits_{p_{0}}^{\,p_{\rm{end}}} {\rm{d}}p\ {\rm{RSM}}({\bf{q}},j,k,p) \,\propto\, {{I_{\rm{S}}({\bf{q}},j,k)}\over{I_{\rm{T}}(\,j,k)}}, \eqno(2)]$

where I_S(q, j, k) is the measured scattering intensity at reciprocal space coordinates (q, θ, ϕ) and real-space coordinates (j, k). Here, (j, k) are two Cartesian coordinates which give the position of the beam relative to the sample in the plane orthogonal to the incident beam direction, and I_T is the transmitted intensity (Liebi et al., 2015). The location of the RSM is given in the experimental system coordinates (j, k, p), where p is the coordinate of the direction in which the X-ray beam travels, see Fig. 1(a).

Figure 1
SAXSTT measurement setup and probed directions. (a) Full setup with detector, goniometer with sample mounted on a pin, X-ray source and measurement coordinate system. (b) Initial sample mounting during first measurement and initial orientation of sample coordinate system. (c) Initial sample mounting during second measurement. (d) Points on the sphere of projection probed during first measurement. (e) Points probed during second measurement. The points on the sphere of projection give the coordinates of the projection direction p in the sample coordinate system spanned by x, y and z.

The subset of RSM(q), which is possible to measure at a given sample orientation under the small-angle approximation, lies on a great circle given by

$[C(\varphi,\alpha,\beta) = \cos(\varphi)\,{\bf{q}}_{0}(\alpha,\beta)+ \sin(\varphi)\,{\bf{q}}_{90}(\alpha,\beta), \eqno(3)]$

where φ is an angle on the detector, q₀(α, β) and q₉₀(α, β) are two unit vectors in the sample coordinate system aligned with the 0° and 90° direction of the detector, respectively, (α, β) are two angles which give the sample orientation as a sequence of rotations about two axes $[\hat{\boldalpha}]$ (the main tomographic rotation axis) and $[\hat{\boldbeta}]$ (the tilt axis used for tensor tomography) orthogonal to the direction of the impinging beam, see Fig. 1(a). In Fig. 2, the general relationship between directions on a sphere (the direction of the impinging beam) and their unique orthogonal great circle (probed part of reciprocal space) is illustrated.

Figure 2
Examples of vectors defining directions on a sphere, and their respective orthogonal great circles. The three vectors labeled v₁ (blue), v₂ (red) and v₃ (yellow) each have a unique orthogonal great circle C₁, C₂ and C₃.

The sample is mounted on a rotation stage such that the first rotation R_β also rotates $[\hat{\boldalpha}]$ , and a sequential rotation of the sample may therefore be described by the composite rotation R_βR_α. At each tilt and rotation, a raster scan over a square grid spanned by the coordinates (j, k) is carried out, yielding a two-dimensional image mapping in which each pixel corresponds to a measured diffraction pattern. We choose the sample coordinate system such that it coincides with the experimental coordinate system when α, β = 0, see the initial sample mounting in Fig. 1(b). To simplify this analysis without loss of generality, we will parameterize the measured reciprocal space vector as $[\hat{\boldalpha}]$ = $[{\bf{q}}_{90}(0,0)]$ and $[\hat{\boldbeta}]$ = $[{\bf{q}}_{0}(0,0)]$ . Then, the coordinate system of the sample will be subject to the composite rotation $[R_{{\bf{q}}_{0}}(\beta)\,R_{{\bf{q}}_{90}}(\alpha)]$ , and the direction of the impinging beam in the sample coordinate system thus changes according to

$[{\bf{p}}(\alpha,\beta) = R^{T}_{{\bf{q}}_{90}}(\alpha) \, R^{T}_{{\bf{q}}_{0}}(\beta) \, {\bf{p}}(0,0), \eqno(4)]$

where $[{\bf{p}}(0,0)]$ is the direction of the beam in the sample coordinate system prior to any rotation or tilt of the sample. We can understand the values taken by $[{\bf{p}}(\alpha,\beta)]$ as points on a sphere of projection, which is a unit sphere consisting of all unique measurement directions up to Friedel symmetry, see Fig. 1(d). Equation (4) also applies to $[{\bf{q}}_{90}]$ and $[{\bf{q}}_{0}]$ .

A composite rotation of the projection vector about the two axes $[R_{{\bf{q}}_{0}}(\beta)\,R_{{\bf{q}}_{90}}(\alpha)]$ may be described as a single rotation around a third axis $[\hat{\imath}]$ , which lies on C(φ, α, β). Consequently, this direction is a rotational invariant with respect to said composite rotation. This means that $[{\rm{RSM}}(\|{\bf{q}}\|\,\hat{\imath},{\bf{r}})]$ can be regarded as a scalar quantity for the purpose of tomography, and standard tomographic analysis can therefore, in principle, be applied to the reconstruction of this component. Although carrying out the experiment in practice requires a third rotation axis when remounting, as it is not possible to tilt the sample by more than 45°, it is possible to specify all points on the sphere of projection using rotations about only two orthogonal axes. Any component of the rotation that occurs about the axis of projection can be discounted in this analysis, since it does not change the information contained in the projection. The following line of reasoning is therefore also valid when combining data from the the two measurements. According to the projection-slice theorem, the Fourier transform of a projection along p constitutes a slice orthogonal to p in Fourier space (e.g. Garces et al., 2011 ). Tomographic reconstruction can therefore be understood as the problem of interpolating between slices in Fourier space. This implies that a set of sufficiently densely placed projections along a great semicircle on the sphere of projection must be measured for a reconstruction of good quality of any given point on the RSM. This leads to the so-called missing wedge problem, where the absence of projections along any section of this great semicircle leads to a missing wedge in the Fourier transform of the reconstruction, and thus a blurring in that direction.

In Fig. 3, projections of reconstructions from data set 1 [Fig. 1(d)] along the y-direction of the absorbance as well as the RSM amplitude in three different directions are shown, along with the discrete Fourier transform of each projection. The absorbance of each projection can be defined as

$[a(\,j,k) = -\log\left[{{I_{\rm{T}}(\,j,k)}\over{I_{0}(\,j,k)}}\right],]$

where I_T is the transmitted intensity at each point in the raster scan, and I₀ is the incident intensity, where we use the projection-wise approximation I₀( j,k) ≃ $[\max{I_{\rm{T}}(\,j,k)}]$ , since each measurement includes some air, which has a very small absorbance. The absorbance is a scalar quantity and therefore very accurately reconstructed without missing wedges. Hence, it is included for comparison. This projection direction, along the y-axis (or equivalently for data set 1, the main tomographic axis $[\hat{\boldalpha}]$ ), is not part of the measurement of data set 1, as seen in Fig. 1(b), which is why it is useful in illustrating missing wedges. The x- and z-components of the RSM shown in Figs. 3(b) and 3(d) can be regarded as measurements which are missing from the data; the y-component in Fig. 3(c) would not actually be measured along this projection direction but does not suffer from missing wedges, because it is measured from every orthogonal direction [the measurements which lie on the line where y = 0 in Fig. 1(d)]. We can observe how the amplitudes in Figs. 3(b) and 3(d), which are those of RSM components orthogonal to $[\hat{\boldalpha}]$ , are smeared out in the directions orthogonal to the RSM component compared with Figs. 3(a) and 3(c). Equivalently, the discrete Fourier transforms are attenuated in the directions of smearing, relative to the more symmetric Fourier transforms in Figs. 3(e) and 3(g). The attenuated segments of the discrete Fourier transforms correspond to missing orthogonal projections, per the projection-slice theorem.

Figure 3
Illustration of the missing wedge problem in SAXSTT. Panels (a)–(d) show a projection along the y-direction of the component of reconstructions of (a) the absorbance and (b), (c) and (d) the RSM aligned with the the x-, y- and z-directions, respectively. Panels (e)–(h) show the amplitudes, normalized by the zero frequency component, of the discrete Fourier transforms of the projections as functions of the discrete frequency ω. Smearing of (b) and (d) occur in the z- and x-directions, respectively, but no smearing can be seen for (a) or (c).

A set of reconstruction constraints similar to those given by the projection-slice theorem exist for the general case of three-dimensional projections in the form of John's equation (e.g. Ma et al., 2017 ), which has been generalized to the case of arbitrary-rank symmetric tensor fields, resulting in additional smoothness constraints on the components of the tensor field (Sharafutdinov, 2012 ; Nadirashvili et al., 2016 ). Therefore, to treat this problem more generally, and not just for discrete, precisely measured RSM components, we need to assume that the RSM does not change too quickly across real and reciprocal space dimensions. Then, given the existence of the invariant axis $[\hat{\imath}]$ , we can define a sampling quality factor on the reciprocal space sphere based on the density of sampling on the sphere of projection. To accomplish this, we define a Friedel symmetric sampling density $[\rho[{\bf{p}}(\alpha,\beta)]]$ as

$[\rho\big[{\bf{p}}(\alpha,\beta)\big] = \left\{ \matrix{ \eqalign{&1,\cr&\cr&} \hfill & \eqalign{ & {\rm{if\ the\ direction\ }} {\bf{p}}(\alpha^{\prime},\beta^{\,\prime})\ {\rm{of\ the }} \cr & {\rm{nearest\ measurement\ satisfies\ }} \cr& \Delta\big[{\bf{p}}(\alpha,\beta),{\bf{p}}(\alpha^{\prime},\beta^{\,\prime})\big]\ \lt\ \delta, } \hfill \cr&\cr& \cr 0,\hfill & {\rm{otherwise,}}\hfill } \right.]$

with $[{\bf{p}}(\alpha,\beta)]$ being defined by equation (4), and where $[\Delta({\bf{v}},{\bf{u}})]$ is the Friedel symmetric great-circle distance, defined as

$[\Delta({\bf{v}},{\bf{u}}) = \arccos\left({{{|{\bf{v}}\cdot{\bf{u}}|} \over {\|{\bf{v}}\|\|{\bf{u}}\|}}}\right), \eqno(5)]$

where v and u are two vectors. The precise choice of the threshold parameter δ depends on assumptions about both the real and reciprocal space continuity of the sample, the size of the sample, as well as the chosen reconstruction method. We chose δ under the assumption that our sampling density in well sampled regions was sufficient to obtain a good tomographic reconstruction. This was chosen over more quantitative thresholds such as the great-circle distance implied by the Nyquist–Shannon sampling theorem applied to tomographic reconstruction (Natterer & Wübbeling, 2001 ) for several reasons. First, we carry out the reconstruction under smoothness and sparsity constraints, which reduce the required number of samples. Second, John's equation for tensor tomography complicates the assumptions based on which standard sampling factors are calculated, since they will also depend on smoothness in reciprocal space (Nadirashvili et al., 2016), and consequently using standard measures would give a misleading impression of certainty about the required number of samples. Therefore, to avoid over-complicating the analysis, we prefer to use $[\rho\,[{\bf{p}}(\alpha,\beta)]]$ as a relative measure of sampling density.

Because we sample points on the reciprocal space sphere orthogonal to the respective points on the sphere of projection, we can compute a quality factor by evaluating the Funk–Radon transform, i.e. the normalized integral over the orthogonal great circle (Funk, 1913 ), of the sampling density ρ. In other words, the quality factor can be computed as

$[F[\rho](\theta,\phi) = {{1}\over{2\pi}} \int_{0}^{\,2\pi} {\rm{d}}\tau\ \rho\big[C(\tau,\theta,\phi)\big], \eqno(6)]$

where (θ, ϕ) specifies a direction in reciprocal space and C(τ, θ, ϕ) is defined as in equation (3). Although we defined C(τ, θ, ϕ) as giving the directions in reciprocal space measured, given a real-space direction, the symmetry of this relationship means that we can invert it to give real-space directions to measure, given a reciprocal space direction that we wish to probe. The relationship between directions on a unit sphere and their unique orthogonal great circle is shown in Fig. 2. Applying equation (6), the value of the quality factor at, for example, v₁ would be defined by the integral of ρ over all directions in C₁. Using this quality factor, which lies in the range [0, 1], where 0 means the lowest possible quality and 1 means the highest possible quality, we may predict the quality of the tomographic reconstruction at any point in reciprocal space.

In order to obtain a reconstruction suitable for the analysis of the missing wedge problem, we want to represent the RSM on the sphere using smooth local functions on the sphere which can be projected into our measurement basis. The locality of the representation is important, since non-local artifacts (such as the spherical harmonic Gibbs phenomenon) could otherwise affect the evaluation of the reconstructions in an unpredictable fashion. We also want to reduce the measured data into azimuthal bins in order to keep the data size manageable, and this reduction can be represented by integrating $[I_{\rm{S}}[\|{{\bf{q}}}\|C(\varphi,\alpha,\beta),j,k]]$ [equation (2)] over segments of φ [as in, for example, Bunk et al. (2009 )]. Schaff et al. (2015) utilized the existence of an invariant point in the RSM for any given rotation, referring to this point as a `virtual axis'. However, Schaff et al. (2015) employed only limited reduction of the measured RSM, necessitating extensive sorting of measurements according to their nearest virtual axis, and processing of a large number of separate tomographic problems, followed by subsequent composition and analysis of the separate reconstructions. De Falco et al. (2021 ) also utilized rotational axis invariance, and carried out a reconstruction using the component of the SAXS measurement orthogonal to the main axis of rotation to study a subpopulation of mineral particles within a sample. These approaches utilize the separability of measurements in order to simplify the tomographic reconstruction problem. However, azimuthal binning reduces this separability, unless the bins are made very small, which would work against the purpose of reducing data size. Thus, rather than aiming to carry out reconstructions at specific points in reciprocal space and subsequently fitting a function to these points, we define a grid of Gaussian basis functions on the unit sphere. This basis set forms a local representation of spherical functions which is used to interpolate measured data into a smooth function on the sphere (Fornberg & Piret, 2008 ). Gaussian radial basis functions have an advantage over the spherical harmonic representation used by Nielsen et al. (2023) because they do not suffer from the Gibbs phenomenon or other non-local artifacts (Gelb, 1997). We define a set of projection matrices from the spherical RSM to the detector by left-multiplication as

$[G_{i,nm} = {{1}\over{{\cal{N\!}}_{n}\,\varphi^{\,\prime}_{m}}} \int\limits^{\varphi_{m+1}}_{\varphi_{m}} {\rm{d}}\tau\ \exp\left\{ {{ \Delta\left[\left(1,\theta_{n},\phi_{n}\right),C\left(\tau,\alpha_{i},\beta_{i}\right)\right]^{2} }\over{ 2\sigma^{2} }} \right\}, \eqno(7)]$

where $[{\cal{N}}_{\!n}]$ is a normalization factor, [φ_m, φ_m+1) parameterizes the mth detector segment on the unit circle, $[\varphi^{\,\prime}_{m}]$ = |φ_m − φ_m+1|, (α_i, β_i) gives the sample orientation, (1, θ_n, ϕ_n) is a unit vector expressed in spherical coordinates giving the location of the mode of basis function n, σ parameterizes the width of each basis function, Δ(u, v) for any two vectors (u, v) is the great-circle distance defined by equation (5), and finally C(τ, α_i, β_i) is defined by equation (3), with τ being an integration variable that parameterizes the integration over each segment. The normalization factor $[{\cal N\!}_{n}]$ , which evens out irregularities in the distribution of grid points, is given by the sum of all rows in an auto-projection matrix $[G_{nn^{\,\prime}}]$ , which can be expressed in a similar form as equation (7) but evaluated only at one point rather than integrated over,

$[{\cal N}_{n} = \sum\limits_{n^{\prime}}G_{nn^{\prime}} = \sum\limits_{n^{\prime}} \exp\left\{ {{ \Delta\big[\left(1,\theta_{n},\phi_{n}\right),\left(1,\theta_{n^{\prime}},\phi_{n^{\prime}}\right)\big]^{2} }\over{ 2\sigma^{2} }} \right\},]$

where n and n′ both run the indices of all basis functions. In this work, the basis functions have been distributed on a modified Kurihara mesh (Kurihara, 1965 ), with an approximately equal distribution over the unit sphere. See Fig. S2 of the supporting information for an illustration of the Gaussian kernel representation. The modified Kurihara mesh depends on an integer scale parameter s which determines the number of basis functions on the hemisphere, according to N = 2s². The width parameter was chosen based on a simple heuristic for smooth and non-oscillatory interpolation, $[\sigma]$ = $[{{\pi}/{2s}}]$ . We chose s = 9 as the scale parameter, thus yielding N = 162 basis functions and width parameter $[\sigma]$ = $[{{\pi}/{18}}]$ . These choices yield smoothly interpolated functions without oscillations, and a density of basis functions greater than the density of detector segments, but smaller than the density of projection directions. For more details on the modified Kurihara mesh and the basis functions, see Supplementary Note S2. Since Gaussians do not have compact support, i.e. they do not fall off to zero, the kernels are local only in a non-strict sense – the vast majority of the amplitude of each basis function is located within a small area around its mode, assuming the standard deviation σ is at least a few times smaller than $[{{\pi}/{2}}]$ . Because of this locality property, we expect the reliability of the coefficient of the basis function located at (θ_n, ϕ_n), considered across real space, to be related to the value of the quality factor F[ρ](θ_n, ϕ_n) given by equation (6).

Completing the description of the forward model requires the definition of a John transform matrix for tensor tomography, which is treated in greater detail by Nielsen et al. (2023). The John transform is also known as the X-ray transform. The resulting set of matrices P_i together define a transform of a tensor field in three-dimensional space into a tensor field in projection space, with each i indicating a projection direction, similarly to how G_i [equation (7)] describes a transform between detector space and spherical space. This allows us to describe the system of equations to be solved for each projection i as

$[{\bf P}_{i}{\bf{X}}{\bf G}_{i} = {\bf{D}}_{i}, \eqno(8)]$

where D_i is matrix of data measured from a single projection.

3. Methods

3.1. Formalism

In order to improve the rate of convergence of the system in equation (8), we compute a series of weight and preconditioning matrices. Each weight matrix is computed as

$[{\bf W}_{i} = ({\bf P}_{i}{\bf{U}}{\bf G}_{i})^{\circ(-1)},]$

where U is a matrix filled with the value 1 everywhere and A°⁽⁻¹⁾ denotes a relaxed element-wise multiplicative inverse of A,

$[\left[{\bf A}^{\circ(-1)}\right]_{ij} \ = \left\{ \matrix{ A_{ij}^{-1},\hfill & {\rm{if}}\ A_{ij}\,\geq\,\epsilon,\hfill \cr \epsilon^{-1},\hfill & {\rm{otherwise,}}\hfill } \right.]$

for some predefined ε > 0. Similarly, each preconditioning matrix is computed as

$[{\bf C}_{i} = ({\bf P}^{T}_{i}{\bf{V}} {\bf G}^{T}_{i})^{\circ(-1)},]$

where V is a matrix filled with the value 1 everywhere. We may now write the system to be solved for each projection i as

$[{\bf C}_{i}\odot\left\{{\bf P}^{T}_{i}\left[{\bf W}_{i}\odot\left({\bf P}_{i}{\bf{X}}{\bf G}_{i}\right)\right]\,{\bf G}_{i}^{T}\right\} = {\bf C}_{i}\odot\left[{\bf P}^{T}_{i}\left({\bf W}_{i}\odot{\bf{D}}_{i}\right)\,{\bf G}_{i}^{T}\right],]$

where ⊙ denotes the Hadamard or elementwise product. This is analogous to the weights and preconditioner used in the SIRT algorithm for scalar tomography [see, for example, Gregor & Fessler (2015 )], which was utilized by Schaff et al. (2015) in the separate reconstructions about each virtual axis. This weight and preconditioner pair serves to normalize the gradient by accounting for the number of voxels that contribute to each pixel, and the number of projections that contribute to each voxel. This normalization is done for each detector segment and each RSM basis function, accounting also for the detector-to-sphere mapping of equation (7). This system is then solved through least-squares Nestorov-accelerated gradient descent, subject to coefficient-wise total variation and L₁ norm regularization, optimized in the Huber approximation of each (Huber, 1964 ), using existing implementations in the mumott package (Nielsen et al., 2023, 2024 ). The reconstruction of the full data set took 630 s, and the reconstruction of each partial data set took 380 s, on a workstation using an Nvidia RTX 3060 GPU, an 8-core AMD Ryzen 7 3700X CPU and 64 GB DDR4 2666 MHz RAM.

3.2. Computations

The integral in equation (7) was computed by quadrature utilizing the adaptive Simpson's rule (e.g. Lyness, 1969 ), terminating when the largest change in a matrix element, relative to the largest element in the matrix, fell below 10⁻⁵. Analysis of the orientation (Fig. 8), the Funk–Radon transform [equation (6)] and computation of the scalar quantities in Fig. 7 requires transforming the spherical function representation from a local Gaussian kernel representation to a spherical harmonic representation. This is done by Driscoll–Healy quadrature (Driscoll & Healy, 1994 ), sampling the function by evaluating the representation on a dense curvilinear grid. The mean amplitude, and the relative anisotropy, are defined as in Nielsen et al. (2023), i.e. as the spherical mean and the spherical standard deviation normalized by the mean. These figures of merit are similar to the `symmetric intensity' and `degree of orientation' used by, for example, Bunk et al. (2009). The fiber symmetry factor is given by

$[S({\bf{a}}) = {{ \left({\sum_{\ell\,=\,1} \sum_{m\,=\,-\ell}^{\ell} F[{\bf{a}}]_{\ell}^{m}\,\hat{Y\,}\!_{\ell}^{\,m}(\theta,\phi)} \right)^{1/2} }\over{ \left({\sum_{\ell\,=\,1} \sum_{m\,=\,-\ell}^{\ell} F[{\bf{a}}]_{\ell}^{m}} \right)^{1/2} \left( {\sum_{\ell\,=\,1}\sum_{m\,=\,-\ell}^{\ell}\hat{Y\,}\!_{\ell}^{\,m}(\theta,\phi)} \right)^{1/2} }}, \eqno(9)]$

where $[{\bf{a}}]$ is the spherical harmonic representation of a RSM, F[·] is the Funk–Radon transform, $[\hat{Y\,}\!_{\ell}^{\,m}]$ is a spherical harmonic basis function and (θ, ϕ) is the orientation of the RSM, as given by the minimal eigenvector of its rank-2 tensor representation. This figure of merit evaluates how similar the RSM is to an ideal ring function (which has all of its amplitude at a great circle consisting of the points orthogonal to its orientation). Consequently, it quantifies the extent to which an RSM exhibits the equatorial symmetry expected from diffuse mineral scattering in bone.

The orientation error is computed as Δ(v, u) [equation (5)], where v and u are two orientation vectors. The orientation error is thus the angle subtended by the two orientation vectors, accounting for the Friedel symmetry of orientation vectors.

3.3. Implementation

The version of mumott used in this work can be found at https://doi.org/10.5281/zenodo.10708583 (Nielsen et al., 2024). New versions of mumott are continuously made available at https://zenodo.org/doi/10.5281/zenodo.7919448. The John transform in mumott is implemented using a bilinear interpolation algorithm which supports multiple channels per voxel and pixel, written using the CUDA API of the Python package Numba (Lam et al., 2015 ). The algorithm employed is based on the work of Xu et al. (2010 ) and Palenstijn et al. (2011 ). Other computations were carried out using the Python packages NumPy and SciPy (Harris et al., 2020 ; Virtanen et al., 2020 ). Two-dimensional plots were created using the package Matplotlib (Hunter, 2007 ). The color maps used in this work are from the package ColorCET (Kovesi, 2015 , 2020 ). The experimental setup render in Fig. 1 was created using Blender (Blender Online Community, 2018 ). All other 3D renders in this work were created using ParaView (Ahrens et al., 2005 ). The two data sets were aligned using the cross-correlation algorithm of Guizar-Sicairos et al. (2008 ), and rotated by modifying the set of vectors used to define the reconstruction geometry in mumott. The rotations were first determined by eye and then refined by comparing absorptivity reconstructions in ParaView.

4. Experiment

The sample chosen for this study was trabecular bone fixed and embedded in polymethyl methylacryate (PMMA). A cube was extracted from the bulk and subsequently milled into a cylinder of diameter 1.2 mm and height 1.2 mm using a custom-made lathe system (Holler et al., 2020 ). The sample was measured at the cSAXS beamline of the Swiss Light Source (SLS) at the Paul Scherrer Institut (PSI), Switzerland. The X-ray energy was set to 12.4 keV using a Si(111) double-crystal monochromator, and the scattering patterns were recorded on a Pilatus 2M detector placed at a sample-to-detector distance of 2.17 m. A flight tube, approximately 2 m in length, was placed in between the sample and detector to reduce the air scattering. A 1.5 mm steel beamstop inside the flight tube blocked the directly transmitted beam. The fluorescence signal from the beamstop, proportional to the intensity of the impinging X-rays [I_T in equation (2), was measured by a Cyberstar (Oxford Danfysik). This allowed the relative X-ray transmission through the sample to be measured. The sample was measured with a beam that had a full width at half-maximum of 12 µm × 24 µm as measured by a knife-edge scan. The raster scan used a step size of 25 µm in both the vertical and horizontal directions, with continuous fly scanning in the vertical direction. The experimental setup is illustrated in Fig. 1(a). Two sets of small-angle X-ray scattering tensor tomography measurements were carried out, each consisting of 224 scanning SAXS images. The two sets of measurements are shown on the sphere of projection in Figs. 1(d) and 1(e), where each marker indicates the direction of the X-ray beam [given by p(α, β) in equation (4)] in the sample coordinate system.

During the first set of SAXSTT measurements, the base of the cylinder sample was glued to the end of a PMMA needle, using a hot water-soluble glue (Norland Blocking Adhesive 107). Before the second SAXSTT experiment, the sample was glued with UV-glue (Norland Optical Adhesive 81) to a second pin, before detaching the first pin by placing the sample in hot water. The second pin was placed at approximately 90° to the first pin, measured around the axis of the initial direction of the beam, see Fig. 1. In total, 1716960 scattering images were measured. The measurement of the first data set took 1218 min, and the measurement of the second data set took 1364 min.

Figs. 4(a)–4(c) show the directions of measurement on the unit hemisphere of projection while Figs. 4(d)–4(f) show the quality factor F defined by equation (6) on the reciprocal space hemisphere. Note that Friedel symmetry is accounted for in the hemispheric representation. The reciprocal space quality factors follow the expected symmetry, where measurements along the entirety of a great semicircle result in a quality factor of 1 at the point orthogonal to this semicircle. The lowest obtained quality factor is 0.5, since the lowest possible coverage (for data set 1) of a great semicircle occurs when the semicircle lies at a single longitude and varies only in latitude. Such a semicircle is still covered by measurements at a fixed longitude, the latitude (tilt, for data set 1) of which span the range [−45°, 45°] — thus, in the worst-case scenario, half of the semicircle's total arc length of 180° is covered.

Figure 4
Points on the hemisphere of projection and theoretical quality factors. (a) Probed points on the sphere of projection in first measurement. (b) Probed points in second measurement. (c) Combined points from both measurements. (d) Quality factor in reciprocal space from first data set. (e) Quality factor from second data set. (f) Quality factor from combined data sets. The dotted lines show great circles at longitudes 0°, ±30° and ±60° with the y-axis as the meridian. The dashed lines show small circles with elevations of 0°, ±30° and ±60° with the x-axis as the equator.

For the reconstruction and analysis, a q-range of 0.597–0.607 nm⁻¹ was used, corresponding to a d-spacing range of 10.36–10.53 nm. This range was used due to artifacts in the second measurement at lower q-ranges, possibly due to the water-soluble glue penetrating the outer layer during the remounting of the sample; see Supplementary Note S3, as well as Supplementary Figs. S4 and S5, for details.

5. Results

The results of comparing a reconstruction of the full data set, which combines data sets 1 and 2, with reconstructions that include, respectively, only data set 1 or 2, are shown in Fig. 5. Valid voxels for comparison, i.e. voxels containing trabecular bone sample, were identified based on the mean amplitude of each RSM in the full dataset reconstruction. The location of each marker in Figs. 5(a) and 5(b) corresponds to the mode of a RSM basis function, see equation (7), while the color of the marker corresponds to the error computed from comparing the coefficients of that basis function to the corresponding coefficients of the full dataset reconstruction. The markers are overlaid on the reciprocal space quality factor. In Figs. 5(c) and 5(d), the error for the amplitude at each point on the reciprocal space sphere is shown. The distribution of errors in the reciprocal space amplitude follows the quality factor closely, with errors above approximately 0.1 occurring exclusively in the region where the quality factor is smaller than 1. The errors for the basis set coefficients in Figs. 5(a) and 5(b) are larger than the errors of the amplitude in Figs. 5(c) and 5(d), which is especially apparent when comparing the upper region of Fig. 5(a) with the same region in Fig. 5(c). This is explained by the fact that the basis set functions are not orthogonal but overlap. This means that some variations in the basis set coefficients cancel out when the amplitude of each RSM function is evaluated for the calculation of the error.

Figure 5
Error distribution for each RSM basis function. (a) Coefficient errors for data set 1. (b) Coefficient errors for data set 2. (c) RSM errors for data set 1. (d) RSM errors for data set 2. The markers have been placed over the corresponding theoretical quality factor from Fig. 4

. The errors in (a) and (b) were calculated by computing the overall Pearson correlation coefficient for each basis function coefficient when comparing the partial and full data sets. In (c) and (d) the correlation coefficients were computed for the RSM amplitude at each coordinate. The correlation factors for the coefficients and the RSM amplitudes are not the same, since the Gaussian radial basis functions of the basis set overlap. The dotted lines show great circles at longitudes 0°, ±30° and ±60° with the y-axis as the meridian. The dashed lines show small circles with elevations of 0°, ±30° and ±60° with the x-axis as the equator.

In Fig. 6 the reconstructed RSMs can be seen in a spherical function glyph render for data set 1 only [Fig. 6(a)], data set 2 only [Fig. 6(b)] and the full data set [Fig. 6(c)]. Each rendered glyph shows the RSM reconstructed in that voxel, colored by its amplitude, and scaled by the Funk–Radon transform of the amplitude. Because the scattering at this q-range is dominated by diffuse equatorial mineral scattering, deforming each glyph by the Funk–Radon transform allows its shape to visually indicate the orientation of the underlying nanostructure.

Figure 6
Spherical function glyph render of reconstructions. (a) Partial data set 1 and a reciprocal space sphere showing the error distribution compared with the full data set. (b) Partial data set 2 and a reciprocal space sphere showing the error distribution compared with the full data set. (c) Full data set. (d) Error distribution in reciprocal space for data set 1. (e) Error distribution in reciprocal space for data set 2. The color of each spherical function indicates the RSM amplitude, whereas the shape indicates the orientation. The shape is computed from the Funk–Radon transform of the RSM amplitude. The insets show two sets of RSMs that the partial reconstructions each have difficulty reconstructing, compared with the full data.

As illustrated in the insets [Figs. 6(d) and 6(e)], which show the error distribution on the RSM, data set 1 has the best sampling and therefore the most reliable reconstruction along the y-axis, i.e. along the main tomographic axis (Fig. 1). Data set 2 has the smallest error along the x-axis. The effect of this on the reconstructed 3D RSM is illustrated in the enlarged views below each render. The respective upper enlarged views shows RSMs that are better reconstructed by data set 1, as data set 2 has difficulty reconstructing the amplitude near the y-axis, leading to increased asymmetry in the equatorial scattering due to missing wedges. The lower enlarged views show RSMs which are better reconstructed by data set 2, as data set 1 has difficulty reconstructing amplitudes near the x-axis, introducing additional texture in the equatorial scattering. Both data sets have some difficulty reconstructing amplitudes that lie along the z-axis, but the difficulty is overall greater for data set 1, as indicated by the distribution on the spherical inset [Fig. 6(d)], compared with the spherical inset [Fig. 6(e)], which shows the amplitude error from Figs. 5(c) and 5(d), respectively, rendered on a spherical surface.

It is likely that the primary explanation for this larger error is that the measurements close to the x-axis on the sphere of projection which result in the large errors near the z-axis must pass through the thickest part of the sample. This means that the transmission is small, around 4%, compared with values of 20–50% for thinner parts of the sample. Consequently, noise in the transmission will have a relatively large impact on these measurements, see equation (2).

Three scalar quantities for each of the three reconstructions can be seen in Fig. 7: the mean RSM amplitude, the relative anisotropy (similar to quantities often referred to as degree of orientation) and a fiber symmetry factor. The fiber symmetry factor quantifies the degree to which the scattering is equatorial, see Section 3 and equation (9) for details. These quantities are of interest in evaluating the RSMs, and therefore their similarity between partial and full data set reconstructions are of importance in evaluating the impact of the missing wedge problem. The mean amplitude in the top row shows no large variations, except for slightly higher values at the edges of protrusions in data set 2, which may be due to the leeching of water-soluble glue into the sample during remounting, see Supplementary Note S3. The mean amplitude is an important scalar value which is used for q-resolved reconstruction and further analysis of nanostructure information contained in the SAXS curve (Liebi et al., 2021 ; Casanova et al., 2023; Silva Barreto et al., 2024). The relative anisotropy is also very similar for all three reconstructions, with almost no discernible differences. Somewhat greater differences can be seen in the fiber symmetry factor, especially in the right-hand-side interface region where the insets in Fig. 6 are located. The full data set has a high fiber symmetry factor in this area except at the very center of this interface, whereas the partial data sets appear to have a lower factor around the edges. Thus, the fiber symmetry factor is more sensitive to missing wedges than the ordinary relative anisotropy. This can also be seen in Fig. 6, where additional texture within the ring of the equatorial scattering appears as an artifact of the missing wedge. For quantitative plots of the distribution of the quantities, see Supplementary Note S4.

Figure 7
Volume renders of scalar quantities. Panels (a)–(c) show the mean amplitude of the RSMs for partial data sets 1 and 2, as well as for the full data set. Panels (d)–(f) show the relative anisotropy of the RSMs. Panels (g)–(i) show the fiber symmetry factor of the RSMs.

One of the most important properties that can be retrieved from a SAXSTT measurement is the local orientation, and it is therefore of interest to see how much uncertainty the missing wedge problem introduces in determining this. Fig. 8 shows glyph renders of the orientation error for each partial reconstruction. The orientation error is defined in Section 3, using equation (5). Most orientations are determined to within an error of no more than 10°, as seen qualitatively for the blue and green colors in Figs. 8(a) and 8(b) and quantitatively in the density plots Figs. 8(c) and 8(d). The enlarged areas shown in green and orange rectangles show a region in the trabecular bone where differently oriented domains are intersecting. This is the region where the orientation error is the largest in both partial data sets. Comparing with Figs. 7(d)–7(i), it can be seen that the larger orientation errors lie in regions where both the relative anisotropy and the fiber symmetry are small. This means that the orientation is less well defined, and may include multiple orientations within a voxel. This is consistent with the region containing an interface of domains of different orientation. The same tendency is seen in the overall RSM error which is largely similar to the distribution of the orientation errors (see Supplementary Note S1).

Figure 8
Orientation errors. (a) Glyph render of orientations and errors of partial data set 1 (b) Glyph render of orientations and errors of partial data set 2. (c) Probability density plot of orientation errors of partial data set 1. (d) Probability density plot of orientation errors of partial data set 2. The color of the glyph indicates the orientation error in that voxel compared with the full data set, with each glyph being scaled by the relative anisotropy in the partial reconstruction. The insets highlight an interface area where the different tendencies of the orientation errors for (a) and (b) can be seen, with (a) showing larger errors for orientations closer to the y-axis and (b) showing larger errors orientations closer to the x-axis. The density plots show the orientation errors for high (greater than 0.6) and low (less than 0.6) relative anisotropy, showing that the error is greater in low relative anisotropy regions.

Fig. 9 illustrates a single slice from the enlarged region in Fig. 8 with larger errors, as well as low values of relative anisotropy and fiber symmetry. The Funk–Radon transform of the RSM shows that in this interface region multiple orientations are present inside single voxels. The reconstruction with a model which does not impose strong symmetries, such as the grid of Gaussian radial basis functions used here, opens up the possibility to extract multiple orientations in each voxel. However, comparing datasets 1 and 2, as well as the full data set illustrates that the missing wedge problem influences the accuracy of the reconstructed RSM in the partial data set reconstructions. The partial data set reconstructions in Figs. 8(a) and 8(b) have more voxels with apparent multi-orientation, and with a greater relative amplitude in the secondary orientation when compared with the full data set reconstruction in Fig. 8(c). This is likely due to missing-wedge smearing of certain parts of the RSM amplitude across real space. Thus, while multi-orientation analysis can be used to precisely localize this interface in a full-data reconstruction, the missing wedge problem makes this localization much more difficult in partial-data reconstructions.

Figure 9
Multiple orientations in interface region. (a) Data set 1, renders of Funk–Radon transform of anisotropic part of the RSM. (b) Data set 2. (c) Full data set. (d) Location of interface region in the sample. Maxima in the Funk–Radon transform indicate the orientation of each RSM, and voxels with multiple local maxima appear to have multiple orientations. There are more apparent multi-orientation voxels in the partial reconstructions in (a) and (b), which is likely due to missing-wedge smearing of certain parts of the RSM amplitude across real space.

6. Conclusions

In this work we have devised a scheme for complete acquisition of SAXSTT data, and applied it to the analysis of a sample of trabecular bone. Reconstructing incomplete as well as complete data sets and comparing them across both real and reciprocal space, we conclude that the understanding of data incompleteness in terms of the missing-wedge problem, as indicated by the computed quality factor, is consistent with the observed errors in the reconstruction. Analyzing the orientations as well as scalar quantities, we find that the impact of the missing-wedge problem in a typical SAXSTT analysis is limited, but appreciable in edge and interface areas. In particular, the impact on mean RSM amplitude and relative anisotropy is very limited, except for apparent artifacts in the mean. Moreover, we observe that the impact of errors can be reduced by choosing the sample orientation during acquisition in a way that takes into account the missing wedge problem, i.e. by orienting the sample such that as much scattering as possible is close to the main axis of rotation. Prior understanding of the nanostructure and expected RSM of a sample, such as acquired by scanning SAXS, is crucial in this process.

This understanding could also be employed in various measures to reduce the impact of the missing wedge problem, e.g. by enforcing a particular RSM symmetry. Such symmetries can be encoded in the SAXSTT basis set [as in Liebi et al. (2018), which used a spherical harmonic model that enforced rotational symmetry about an axis]. One disadvantage of encoding symmetries in the basis set is that more complex textures (such as the multi-orientation investigated in this work) cannot be captured. However, symmetries can also be selectively enforced (based on a robust quantity, such as the relative anisotropy) in a post-processing step, or encouraged through regularization.

The further exploration of these possibilities and their impact on reconstruction quality is an interesting avenue for future research. Finally, we remark that the complete acquisition scheme devised in this work is likely to be useful for specialized applications, such as the analysis of interface regions with overlapping domains of multiple orientations, or the reconstruction of especially complicated RSMs.

Supporting information

Link https://doi.org/10.5281/zenodo.10995088
Data, as well as reconstructions and Paraview data used in this work. Also included is a Jupyter notebook that shows how the reconstruction and analysis in this work is carried out.

Supplementary figures and notes referred to in the main text. DOI: https://doi.org/10.1107/S1600577524006702/vy5028sup1.pdf

Acknowledgements

We acknowledge the Paul Scherrer Institut, Villigen, Switzerland for provision of synchrotron radiation beamtime at the beamline cSAXS of the SLS. We extend our thanks to Manuel Guizar-Sicairos at the Paul Scherrer Institut and EPFL, for discussions on the experimental strategy to solve the missing wedge problem, as well as for co-authoring the experimental proposal. We finally acknowledge Andreas Menzel at the cSAXS beamline at the Paul Scherrer Institut for assisting with the planning and administration of the experiment.

Conflict of interest

The authors declare no competing financial or non-financial interest.

Data availability

The data used in this work, along with code demonstrating the analysis and reconstructions which can be viewed in ParaView, is available at https://doi.org/10.5281/zenodo.10995088.

Funding information

The following funding is acknowledged: H2020 European Research Council (grant No. ERC-2020-StG 949301 to Marianne Liebi); Vetenskapsrådet (grant No. VR 2018-041449 to Marianne Liebi); funding from the Chalmers Initiative for Advancement of Neutron and Synchrotron Techniques.

References

Ahrens, J. P., Geveci, B. & Law, C. C. (2005). The Visualization Handbook, edited by C. D. Hansen & C. R. Johnson, pp. 717–731. Oxford: Elsevier Butterworth–Heinemann. Google Scholar
Blender Online Community (2018). Blender – a 3D modelling and rendering package. Blender Foundation, Stichting Blender Foundation, Amsterdam. https://www.blender.org. Google Scholar
Bunk, O., Bech, M., Jensen, T., Feidenhans'l, R., Binderup, T., Menzel, A. & Pfeiffer, F. (2009). New J. Phys. 11, 123016. Web of Science CrossRef Google Scholar
Casanova, E. A., Rodriguez-Palomo, A., Stähli, L., Arnke, K., Gröninger, O., Generali, M., Neldner, Y., Tiziani, S., Dominguez, A. P., Guizar-Sicairos, M., Gao, Z., Appel, C., Nielsen, L. C., Georgiadis, M., Weber, F. E., Stark, W., Pape, H.-C., Cinelli, P. & Liebi, M. (2023). Biomaterials, 294, 121989. CrossRef PubMed Google Scholar
De Falco, P., Weinkamer, R., Wagermaier, W., Li, C., Snow, T., Terrill, N. J., Gupta, H. S., Goyal, P., Stoll, M., Benner, P. & Fratzl, P. (2021). J. Appl. Cryst. 54, 486–497. Web of Science CrossRef CAS IUCr Journals Google Scholar
Ding, G., Liu, Y., Zhang, R. & Xin, H. (2019). Sci. Rep. 9, 12803. CrossRef PubMed Google Scholar
Driscoll, J. & Healy, D. (1994). Adv. Appl. Math. 15, 202–250. CrossRef Web of Science Google Scholar
Fornberg, B. & Piret, C. (2008). J. Comput. Phys. 227, 2758–2780. CrossRef Google Scholar
Funk, P. (1913). Math. Ann. 74, 278–300. CrossRef Google Scholar
Garces, D. H., Rhodes, W. T. & Peña, N. M. (2011). J. Opt. Soc. Am. A, 28, 766–769. CrossRef Google Scholar
Gelb, A. (1997). Math. C. 66, 699–717. CrossRef Google Scholar
Georgiadis, M., Schroeter, A., Gao, Z., Guizar-Sicairos, M., Liebi, M., Leuze, C., McNab, J. A., Balolia, A., Veraart, J., Ades-Aron, B., Kim, S., Shepherd, T., Lee, C. H., Walczak, P., Chodankar, S., DiGiacomo, P., David, G., Augath, M., Zerbi, V., Sommer, S., Rajkovic, I., Weiss, T., Bunk, O., Yang, L., Zhang, J., Novikov, D. S., Zeineh, M., Fieremans, E. & Rudin, M. (2021). Nat. Commun. 12, 2941. Web of Science CrossRef PubMed Google Scholar
Gregor, J. & Fessler, J. A. (2015). IEEE Trans. Comput. Imaging, 1, 44–55. CrossRef PubMed Google Scholar
Grünewald, T. A., Johannes, A., Wittig, N. K., Palle, J., Rack, A., Burghammer, M. & Birkedal, H. (2023). IUCrJ, 10, 189–198. Web of Science CrossRef PubMed IUCr Journals Google Scholar
Guizar-Sicairos, M., Georgiadis, M. & Liebi, M. (2020). J. Synchrotron Rad. 27, 779–787. Web of Science CrossRef IUCr Journals Google Scholar
Guizar-Sicairos, M., Thurman, S. T. & Fienup, J. R. (2008). Opt. Lett. 33, 156–158. Web of Science PubMed Google Scholar
Harris, C. R., Millman, K. J., van der Walt, S. J., Gommers, R., Virtanen, P., Cournapeau, D., Wieser, E., Taylor, J., Berg, S., Smith, N. J., Kern, R., Picus, M., Hoyer, S., van Kerkwijk, M. H., Brett, M., Haldane, A., del Río, J. F., Wiebe, M., Peterson, P., Gérard-Marchant, P., Sheppard, K., Reddy, T., Weckesser, W., Abbasi, H., Gohlke, C. & Oliphant, T. E. (2020). Nature, 585, 357–362. Web of Science CrossRef CAS PubMed Google Scholar
Holler, M., Ihli, J., Tsai, E. H. R., Nudelman, F., Verezhak, M., van de Berg, W. D. J. & Shahmoradian, S. H. (2020). J. Synchrotron Rad. 27, 472–476. CrossRef IUCr Journals Google Scholar
Huber, P. J. (1964). Ann. Math. Stat. 35, 73–101. CrossRef Web of Science Google Scholar
Hunter, J. D. (2007). Comput. Sci. Eng. 9, 90–95. Web of Science CrossRef Google Scholar
Kovesi, P. (2015). arXiv:1509.03700. Google Scholar
Kovesi, P. (2020). Colorcet, https://colorcet.com/index.html. Google Scholar
Kurihara, Y. (1965). Mon. Weather Rev. 93, 399–415. CrossRef Google Scholar
Lam, S. K., Pitrou, A. & Seibert, S. (2015). Proceedings of the Second Workshop on the LLVM Compiler Infrastructure in HPC (LLVM2015), 15 November 2015, Austin, TX, USA, Article No. 7. New York: ACM. Google Scholar
Liebi, M., Georgiadis, M., Kohlbrecher, J., Holler, M., Raabe, J., Usov, I., Menzel, A., Schneider, P., Bunk, O. & Guizar-Sicairos, M. (2018). Acta Cryst. A74, 12–24. Web of Science CrossRef IUCr Journals Google Scholar
Liebi, M., Georgiadis, M., Menzel, A., Schneider, P., Kohlbrecher, J., Bunk, O. & Guizar-Sicairos, M. (2015). Nature, 527, 349–352. Web of Science CrossRef CAS PubMed Google Scholar
Liebi, M., Lutz-Bueno, V., Guizar-Sicairos, M., Schönbauer, B. M., Eichler, J., Martinelli, E., Löffler, J. F., Weinberg, A., Lichtenegger, H. & Grünewald, T. A. (2021). Acta Biomaterialia, 134, 804–817. Web of Science CrossRef CAS PubMed Google Scholar
Liu, J., Guan, Y., Chen, L., Bai, H., Wei, W., Tian, Y. & Liu, G. (2018). Microsc. Microanal. 24, 138–139. Google Scholar
Lyness, J. N. (1969). J. ACM, 16, 483–495. CrossRef Google Scholar
Ma, J., Wu, S., Qi, H., Li, B., Yan, H., Zhou, L. & Xu, Y. (2017). Sci. Rep. 7, 4920. CrossRef PubMed Google Scholar
Moebel, E. & Kervrann, C. (2020). J. Struct. Biol. X, 4, 100013. Web of Science PubMed Google Scholar
Nadirashvili, N. S., Sharafutdinov, V. A. & Vlăduţ, S. G. (2016). Inverse Probl. 32, 105013. CrossRef Google Scholar
Natterer, F. & Wübbeling, F. (2001). Mathematical Methods in Image Reconstruction. Society for Industrial and Applied Mathematics. Google Scholar
Nielsen, L., Carlsen, M., Wang, S., Baroni, A., Tänzer, T., Liebi, M. & Erhart, P. (2024). mumott – a Python library for the analysis of multi-modal tensor tomography data. https://doi.org/10.5281/zenodo.10708583. Google Scholar
Nielsen, L. C., Erhart, P., Guizar-Sicairos, M. & Liebi, M. (2023). Acta Cryst. A79, 515–526. Web of Science CrossRef IUCr Journals Google Scholar
Palenstijn, W., Batenburg, K. & Sijbers, J. (2011). J. Struct. Biol. 176, 250–253. Web of Science CrossRef CAS PubMed Google Scholar
Schaff, F., Bech, M., Zaslansky, P., Jud, C., Liebi, M., Guizar-Sicairos, M. & Pfeiffer, F. (2015). Nature, 527, 353–356. Web of Science CrossRef CAS PubMed Google Scholar
Sharafutdinov, V. A. (2012). Integral Geometry of Tensor Fields. Berlin: Walter de Gruyter. Google Scholar
Silva Barreto, I., Pierantoni, M., Nielsen, L. C., Hammerman, M., Diaz, A., Novak, V., Eliasson, P., Liebi, M. & Isaksson, H. (2024). Acta Biomaterialia, 174, 245–257. Web of Science CrossRef CAS PubMed Google Scholar
Skjønsfjell, E. T., Kringeland, T., Granlund, H., Høydalsvik, K., Diaz, A. & Breiby, D. W. (2016). J. Appl. Cryst. 49, 902–908. Web of Science CrossRef IUCr Journals Google Scholar
Stribeck, N., Camarillo, A. A., Nöchel, U., Schroer, C., Kuhlmann, M., Roth, S. V., Gehrke, R. & Bayer, R. K. (2006). Macro Chem. Phys. 207, 1139–1149. CrossRef CAS Google Scholar
Trampert, P., Wang, W., Chen, D., Ravelli, R. B., Dahmen, T., Peters, P. J., Kübel, C. & Slusallek, P. (2018). Ultramicroscopy, 191, 1–10. CrossRef CAS PubMed Google Scholar
Virtanen, P., Gommers, R., Oliphant, T. E., Haberland, M., Reddy, T., Cournapeau, D., Burovski, E., Peterson, P., Weckesser, W., Bright, J., van der Walt, S. J., Brett, M., Wilson, J., Millman, K. J., Mayorov, N., Nelson, A. R. J., Jones, E., Kern, R., Larson, E., Carey, C. J., Polat, İ., Feng, Y., Moore, E. W., VanderPlas, J., Laxalde, D., Perktold, J., Cimrman, R., Henriksen, I., Quintero, E. A., Harris, C. R., Archibald, A. M., Ribeiro, A. H., Pedregosa, F., van Mulbregt, P., Vijaykumar, A., Bardelli, A. P., Rothberg, A., Hilboll, A., Kloeckner, A., Scopatz, A., Lee, A., Rokem, A., Woods, C. N., Fulton, C., Masson, C., Häggström, C., Fitzgerald, C., Nicholson, D. A., Hagen, D. R., Pasechnik, D. V., Olivetti, E., Martin, E., Wieser, E., Silva, F., Lenders, F., Wilhelm, F., Young, G., Price, G. A., Ingold, G., Allen, G. E., Lee, G. R., Audren, H., Probst, I., Dietrich, J. P., Silterra, J., Webber, J. T., Slavič, J., Nothman, J., Buchner, J., Kulick, J., Schönberger, J. L., de Miranda Cardoso, J. V., Reimer, J., Harrington, J., Rodríguez, J. L. C., Nunez-Iglesias, J., Kuczynski, J., Tritz, K., Thoma, M., Newville, M., Kümmerer, M., Bolingbroke, M., Tartre, M., Pak, M., Smith, N. J., Nowaczyk, N., Shebanov, N., Pavlyk, O., Brodtkorb, P. A., Lee, P., McGibbon, R. T., Feldbauer, R., Lewis, S., Tygier, S., Sievert, S., Vigna, S., Peterson, S., More, S., Pudlik, T., Oshima, T., Pingel, T. J., Robitaille, T. P., Spura, T., Jones, T. R., Cera, T., Leslie, T., Zito, T., Krauss, T., Upadhyay, U., Halchenko, Y. O. & Vázquez-Baeza, Y. (2020). Nat. Methods, 17, 261–272. Web of Science CrossRef CAS PubMed Google Scholar
Xu, W., Xu, F., Jones, M., Keszthelyi, B., Sedat, J., Agard, D. & Mueller, K. (2010). J. Struct. Biol. 171, 142–153. CrossRef PubMed Google Scholar