Direct x-ray scattering signal measurements in edge-illumination/beam-tracking imaging and their interplay with the variance of the refraction signals Open Access

X ray dark field (DF) imaging has attracted significant interest over recent years, thanks to its ability to detect the presence of features below the spatial resolution of an imaging system.^1–5 It belongs to the class of “phase-sensitive” x-ray imaging methods,⁶ in which image contrast arises from the unit decrement of the real part (δ), rather than from the imaginary part (β) of the refractive index. “Ultra-small angle x-ray scatter” (USAXS) is often used as a synonym, with the caveat that scattering angles involved in DF are in the region of microradians or sub-microradian, while USAXS usually refers to angles of a fraction of a degree. Following its introduction at synchrotron facilities in crystal-based experiments,^1–3 DF became increasingly popular after Pfeiffer et al.'s demonstration that the same signal could be accessed with conventional laboratory sources by using a three-grating setup.⁴ Thereafter, DF with laboratory sources was demonstrated with Edge-Illumination⁵ (EI) and its simplified, single mask “beam tracking” (BT) implementation,⁷ as well as other methods such as speckle-based imaging.^8,9 Following initial applications to, e.g., welding¹⁰ and composite materials,¹¹ DF has significantly evolved and is currently being used for in vivo lung imaging of human patients.¹² A recent study also showed significant potential in the detection of concealed explosives,¹³ and efforts to implement DF on a clinical CT gantry are under way.¹⁴ Among the many phase-based techniques, we focus here on EI¹⁵ (and on its BT implementation, which earlier work demonstrated to provide equivalent signals¹⁶) as it is achromatic¹⁷ and can be implemented with low spatial coherence.¹⁸ 

A straightforward way to understand the nature of the DF signal is by referring to x-ray refraction. This is greatest at sample boundaries (as the refraction angle is proportional to the first derivative of δ in space⁶) and after some propagation, directional changes translate into a change in the (transverse) position at which the x-rays are detected. These changes are usually referred to as “phase contrast” when they are resolved by the imaging system, and as DF when they are not. Refraction, modeled through ray-tracing, provides an approximation of phase-based mechanisms; wave-optics tools (e.g., Fresnel–Kirchhoff diffraction integrals) are required for a more accurate description.⁶ Differences are significant only when the object is coherently illuminated; under limited coherence conditions (e.g., those of a laboratory x-ray source), the refraction (i.e., ray-tracing) and wave-optics descriptions provide compatible results.¹⁹ The former is simpler and lends itself to straightforward implementation in Monte Carlo approaches.

We implemented both ray-tracing Monte Carlo and wave-optics models without resorting to the thin object approximation. Concurrently, we demonstrate that the same results are obtained using classical scattering theory when multiple scattering is negligible. Details on the two models are provided in the methods section.

Scattering measurements can refer to wide angles on the order of several degrees (molecular; crystalline), small angles on the order of arc seconds (complex biological materials; soft matter), and to ultra-small angles on the scale of microradians (microscopic or sub-microscopic structures).²¹ In wide and small angle scattering experiments, the primary beam is often stopped, as the characteristic peaks are of such low intensity that, if the primary beam were integrated alongside them, detectors would saturate quickly. In USAXS as measured in EI/BT and other phase methods, the higher order peaks of $I q$ are not measured, and the broadening of the primary beam becomes the object of analysis.

Nonetheless, by using the wave-optics model developed for EI/BT (see methods), we can demonstrate both adherence to and departure from Eq. (1) under conditions of single and multiple scattering, respectively. Figures 1(a) and 1(b) display the intensity field resulting from a 15 keV plane wave interacting with sphere ensembles of constant packing density in 200 and 400 μm thick slabs. Figure 1(c) shows that the (radially integrated) intensity resulting from interactions with the 200-μm-thick sample approximate Eq. (1) well (aside from the difficult-to-sample deep troughs), while the intensity oscillations tend to disappear in the case of multiple scattering caused by the thicker sphere ensemble. From the model's experimental validation (see below), these results demonstrate the link between scattering theory and EI/BT-based USAXS measurements. Figure 1(d) shows the radial intensity from the same sample as panel (b) modeled with Monte Carlo ray tracing, which only models photon intensities rather than amplitudes, and thus does not reproduce the oscillations caused by coherence effects but, rather, traces the envelope of the function. However, as shown by the dashed yellow line in Fig. 1(c), these same oscillations are washed out by multiple scattering in thick samples like those analyzed in this work. This demonstrates the link, and underlying consistence under the conditions of limited coherence and/or multiple scattering between classic scattering theory, wave-optics, and Monte Carlo modeling.

There are four parameters which may affect the DF signal in our (multiple-scattering) experiment:

1. the beamlet width;

2. the characteristic feature size of the sample;

3. the sample thickness; and

4. the X-ray energy.

Clearly, the sample material would also affect the DF signal but in our experiment is kept constant. The interplay between the first three parameters (at 60 keV) is explored in Fig. 2, while energy is discussed in the supplementary material.

Figures 2(a)–2(d) make it evident that the DF signal increases with sample thickness (the number of traversed layers containing microspheres)^27,28 and with decreasing sphere size, in part because, as explained above, effort was made to keep constant the fill factor by volume which led to a total number of spheres inversely proportional to the cube of their diameter. While a similar dependence on sample thickness and sphere size has been observed before,²⁹ a first novelty observed here is the independence from the size of the beamlets, despite the fact that resolution in EI/BT is defined by the aperture width.³⁰ In addition from the heatmaps in Figs. 2(a)–2(d), this can be clearly observed in the graphs of Figs. 2(e)–2(h). This is an important result because, together with EI/BT's independence from the spacing between adjacent beamlets (shown in previous work^22,30 and further discussed below), it demonstrates the independence of the retrieved scattering signal from system-specific parameters such as the autocorrelation length (defined as ξ = λd/p, where λ is the x-ray wavelength, d the sample-to-detector distance, and p the period of the mask or grating employed).

Another important result evident from Figs. 2(e)–2(h) is the agreement between model and experiment across the widely explored experimental conditions. Taken together with the agreement between model and scattering theory shown above, this demonstrates the ability of EI/BT to perform direct scattering measurements. For increased generality, a subset of samples were also scanned with a conventional lab system (Fig. S1) i.e., with an extended, polychromatic source. The quantitative agreement observed between model and experiment persists and demonstrates that the same scattering measurements can be performed with conventional lab sources. However, in this case the polychromaticity of the measurements requires the adoption of appropriate spectral weighting approaches when modeling.¹⁷ 

Measurements were also taken at 90 and 120 keV to demonstrate EI/BT's robustness against increasing x-ray energy, which opens the way to applications to higher $Z$ or denser materials while remaining quantitative, as proven by the good agreement between model and experiment observed in all cases (Fig. S2). This also provided an opportunity to study the behavior of the DF signal vs energy, with the proven reliability of the simulation allowing us to fill in more energies than the three (60; 90; 120 keV) that were experimentally tested. Results are shown in Fig. S3. The best fit indicates proportionality to E⁻³ (Ref. 7), although the previously reported^22,27 E⁻⁴ dependence also leads to a reasonably good fit with a reduction in the coefficient of determination $R 2$ from 0.998 to 0.988, meaning an E⁻⁴ dependence cannot be excluded. It should also be noted that $E − 3.89$ dependence was previously reported;³¹ hence, this could be an area that requires further investigation. The independence of the EI/BT DF signal from the beamlet size becomes apparent when one considers the way in which it is retrieved in Eq. (2b). By subtracting the beamlet width without the sample, its influence is eliminated. Conversely, in grating-based methods, the DF signal is measured as a reduction in fringe visibility^29,32–34 V^S/V^R, extracted either via fringe-scanning methods³⁵ or Fourier analysis,³⁶ with V = (I_max − I_min)/(I_max + I_min) with I_max, I_min the maximum and minimum intensities observed in a sinusoidal pattern and indices S and R having the same meaning as above. This makes the measured signal dependent on the system's autocorrelation length, ξ, regardless of whether absorption rather than phase gratings is used.³⁷ This significant difference between the two methods is exemplified in Fig. 3.

The key aspect to note in Fig. 3 is that, as predicted by Strobl³³ and Lynch,²⁹ in grating-based methods the DF signal “kicks in” and remains stable when ξ becomes larger than the average size of the scatterers, as can be seen by the curve “plateauing” around ξ = 10 μm in Fig. 3(b). It should be noted that while most GI literature plots the DF signal as $V S / V F$⁠, it is plotted here as 1 − Vr/V s to highlight the fact that the DF signal goes to zero for small values of ξ, and stabilizes to its final (asymptotic) value only when ξ exceeds the size of the scatterers. Conversely, with EI/BT, the scattering value is effectively retrieved at extremely small values of ξ, and is unaffected by it. What happens around ξ = 1.4 μm is that the specific geometry of the system causes excessive beamlet overlap, which prevents a “moment-type” retrieval as it becomes impossible to integrate the beamlets over their full extent. Other methods (e.g., Gaussian perturbation³⁴) can be used to perform an effective DF retrieval at (moderately) larger ξ values; however, the whole point of EI/BT is to keep the beamlet sufficiently separated from each other so as to allow a retrieval approach that focuses on the individual beamlet (considered as a separate, independent entity), rather than on the analysis of the pattern as a whole. Indeed, as beamlet overlap increases, the modulation ultimately becomes sinusoidal, and EI/BT becomes indistinguishable from grating-based approaches.³⁸ With this in mind, it is important not to misinterpret the graph in Fig. 3(a): the decrease in the retrieved DF signal for ξ > 1.4 does not indicate a dependence of the retrieved DF on ξ, but rather that EI/BT systems with such high ξ values should not be built, as they would effectively not be EI/BT systems anymore since they do not allow the individual analysis of each beamlet. In other words, in the described settings an EI/BT system must have a ξ < 1.4 and, in this ξ range, scatter retrieval is not affected by the specific value of ξ. This also means that the concept of ξ does not really apply to EI/BT, and indeed e.g., increasing d and decreasing p do not have the same effect on an EI/BT system.

Models for grating-based DF were also used to determine a sample's “scattering power” as a function of the scatterer's diameter, for example through a “DF extinction coefficient” (see, e.g., Fig. 4 in Lynch 2011²⁹). The models described in this paper allow doing the same for EI/BT, revealing that in this case a monotonic function is obtained (see Fig. S4), because of the discussed independence from specific imaging system parameters. For completeness, the ability of the same (wave-optics) model to reproduce the results obtained using a standard phase gratings setup is also shown (Fig. S5).

It has been argued that the dependence of “visibility-based” DF retrieval on ξ can be used to determine the size of the scatterers, although this requires sweeping ξ, which can be experimentally challenging (the most practical method probably being based on changing the distance between two phase gratings³⁹) This is not possible in EI/BT experiments in which, so long as the beamlets are kept physically separated, a constant value directly linked to the sample's scattering power, which is unaffected by system parameters, is returned. Although models exist that allow retrieving the scatterers' size so long as more information is extracted from the scattering function on top of its variance,²⁵ the fact that the refraction signal is sensitive to objects larger that the system's spatial resolution can be exploited for this purpose.

In analogy with Figs. 2(a)–2(d), Figs. 4(a)–4(d) show the (experimental) behavior of the variance of the refraction signal (VR) as a function of the scatter size and thickness for various beamlet sizes. The signal is obtained by applying a standard deviation filter with a nine-by-nine pixel window, on a flat-field corrected refraction image, squaring it, and averaging the same ROI as for the DF measurements. As can be seen, the behavior is comparable to that of DF vs sample thickness (i.e., number of scattering layers), but opposite vs scatterer size, with the maximum values observed for the larger spheres. This can be understood intuitively as features below the spatial resolution of the system do not produce distinct refraction signals. Also, the VR signal is not independent from the beamlet size, as made clear also by the graphs of Figs. 4(e)–4(h), which also show that the VR signal can be modeled with the same tools used for the DF signal, obtaining a comparable agreement with the experimental results. This allows us to use the model to explore the VR signal of monodisperse spheres as a function of a varying beamlet size, with a denser sampling than would be practically feasible in a real experiment. The results [Figs. 4(i)–4(l)] show a clear plateauing of the signal when the size of the beamlet reaches that of the scatterers, similarly to (visibility-based) DF vs ξ in grating-based methods. As a consequence, it could be used in the same way—e.g., by using a pre-sample mask with small apertures and translating the mask in aperture increments, the separate images then being summed will create virtual beamlets of larger widths.

We are effectively proposing the VR as a new contrast channel, which provides an “equivalent” to DF for randomly distributed features larger than the spatial resolution of the system, with the same characteristics of, e.g., increasing linearly with sample thickness for a given distribution of scatterers. This had already been used heuristically to determine the degree of porosity in composite materials⁴⁰ and is formalized here as an additional, quantitative contrast channel. The agreement with the theoretical models allows using the latter to explore the behavior of VR vs e.g., sample thickness and x-ray energy, with exemplar results reported in Fig. S6. In particular, the linear behavior of VR vs sample thickness demonstrates its suitability to tomographic implementations.

The simultaneous availability of DF and VR signals through a single scan of the same sample thus offers a method for the quantitative characterization of the sample's microstructure, above and below the system's spatial resolution, with direct access to the sample's scattering function unaffected by specific parameters of the imaging system.

Synchrotron experiments were carried out at beamline ID17 of the European Synchrotron Radiation Facility (ESRF). The radiation source is a wiggler placed at ∼145 m from the experimental hutch. The beam is monochromatized by a bent Laue crystal, and energies of 60, 90, and 120 keV were selected. An EI/BT setup was installed in the experimental hutch, a side view schematic of which is shown in Fig. 5. The beam (in yellow) is vertically thin but extends to several cm into the plane of the drawing so as to cover the entire imaged samples.

A tunable Huber slit (A₁) was used as a pre-sample mask, so that samples could be inspected with variable (vertical) beam sizes—sizes of 10, 20, 30, and 40 μm were used. A much larger (300 μm) slit (A₂) was placed in front of an indirect imaging system with an effective pixel size of 22.4 μm (PCO Edge 5.5 sCMOS detector⁴¹) mounted on a vertical translation stage. This allows rapid switching between the EI and the BT configuration. By letting the beam entirely through A₂ (as in the schematic of Fig. 5), the beamlet is directly resolved by the detector pixels (BT mode). By moving A₂ downward or upward so that the beam is partially intercepted by one of its edges, the EI configuration is obtained. In this latter case, multiple sample scans, with slightly different degrees of beam or edge overlap, are required to enable the retrieval of the beam's vertical width; in principle, a much larger (e.g., >300 μm) detector pixel could be used, but in this case the 22.4 μm pixels were binned together for simplicity. In both cases, the sample is scanned vertically through the beam to cover the desired region-of-interest.

These were performed with a standard laboratory EI system, previously described in Astolfo 2022.⁴² Briefly, a COMET MXR-160/11 tungsten anode x-ray source had its focal spot limited to 70 μm in the horizontal direction by a slit, with the vertical dimension left at the original ∼400 μm. For the experiments presented in this paper, it was operated at 60 kVp and 7.7 mA. Two gratings (fabricated by Microworks, Karlsruhe, Germany) with 21 and 28 μm parallel, vertical apertures arranged in 75 and 98 μm periods were used as pre-sample and detector masks, respectively. Both had ∼300-μm-thick gold septa on a 1-mm-thick graphite substrate, and they were placed 1.60 and 2.06 m downstream of the source, respectively. The samples were scanned horizontally (i.e., orthogonally to the gratings' slits) immediately downstream of the pre-sample mask with a Newport (Irvine, CA) translator at a speed of 39.1 μm s⁻¹, and a 100-μm pixel photon-counting detector (XC-FLITE FX2, Direct Conversion, Danderyd, Sweden—now part of Varex), acquiring at 1 frame per second, was placed 4 cm downstream of the detector mask.

The rigorous multi-slice wave optical simulation code developed by Munro⁴³ was adapted to simulate the EI/BT synchrotron experimental setup, as well as the simulations highlighting differences and similarities between EI/BT and grating interferometry. Planar wavefronts are propagated via the Fresnel propagator, with sample effects being applied through multiplication of the wavefront by the samples' complex transmission functions. Due to the large distances between source, sample, and detector, each “layer” of spheres could be simulated as a single projection. To keep computation times manageable, only small sample areas (0.6 × 2.1 mm²) were simulated; this is a particular advantage when considering that the field of view available in such simulations becomes more costly with increasing photon energy in order to preserve robustness against aliasing. Finally, upon reaching the detector, beam divergence is accounted for by application of the Fresnel scaling theorem.

The McXtrace neutron/x-ray scattering engine developed by Knudsen et al.⁴⁴ was used to simulate the synchrotron and lab-based experimental setups. Refraction is implemented in three dimensions via Snell's law, whereby photons refract or reflect as they encounter sample boundaries.²⁴ Absorption is calculated by adjustment of photon “weights” according to the Beer–Lambert law, and in this way, all photon paths that are computed contribute to the final image. In the polychromatic lab setting, photons are assigned energies according to the CdTe efficiency-weighted spectrum produced by TASMICS.⁴⁵ 

Virtual samples of randomly arranged and tightly packed spheres were created for use in both the simulation frameworks described above via the Advancing Front algorithm described by Valera et al.⁴⁶ The average volumetric fill of simulated spheres was 42 ± 8%.

Melted gel wax was poured into a 2 × 2 × 1 cm³ mold, with a 1 × 1 × 0.5 cm³ cuboid void extending inward from one of the wide faces, and allowed to cool. Mass-controlled amounts of calibrated microspheres were added to the voids and manually vibrated until their distribution was approximately uniform. Melted gel wax was then poured into the void to seal the microspheres in place, with the excess gel wax being scraped away with a glass microscope slide moving parallel to the existing gel wax surface to minimize any variation in gel wax (i.e., background material) thickness. As the prepared samples cooled, they were monitored for any signs of air having been trapped in the gel wax, as such pockets of air would produce their own phase signals—any such exhibiting samples were discarded to ensure that only the calibrated microspheres contributed to the final image contrasts.

See the supplementary material for figures describing simulated data, which (Fig. S1) matches lab-based, polychromatic DF data; Figs. S2 and S3 match multiple monochromatic, synchrotron DF measurements and describe the dependence of DF with x-ray energy; Fig. S4 describes the scattering power per unit thickness due to different feature sizes in EI and BT; Fig. S5 demonstrates the wave-optical model compatibility with DF measurements in a grating interferometry system; and Fig. S6 describes the dependence of VR with sample thickness and x-ray energy.

This work was supported by the EPSRC (Grant No. EP/T005408/1). A.O. and C.K.H. are supported by the Royal Academy of Engineering under the “Chairs in Emerging Technologies” and “Research Fellowships” schemes, respectively. P.R.T.M. was supported by the Royal Society under the “University Research Fellowships” scheme. The authors thank Dr. Heinz Amenitsch (Elettra Synchrotron Light Source, Trieste, Italy) for useful discussions.

D. Bate is a Nikon employee; M. Endrizzi, P.R.T. Munro, and A. Olivo are named inventors on patents owned by UCL protecting the imaging technology used to obtain the presented results.

Ian Buchanan: Conceptualization (equal); Data curation (equal); Formal analysis (lead); Investigation (lead); Project administration (equal); Software (lead); Visualization (lead); Writing – original draft (lead); Writing – review & editing (lead). Silvia Cipiccia: Data curation (equal); Methodology (equal); Writing – review & editing (equal). Carlo Peiffer: Data curation (equal); Writing – review & editing (equal). Carlos Navarrete-León: Data curation (equal); Writing – review & editing (equal). Alberto Astolfo: Conceptualization (equal); Writing – review & editing (equal). Tom Partridge: Writing – review & editing (equal). Michela Esposito: Investigation (equal); Visualization (equal); Writing – review & editing (equal). Luca Fardin: Data curation (equal); Writing – review & editing (equal). Alberto Bravin: Data curation (equal); Writing – review & editing (equal). Charlotte Klara Hagen: Writing – review & editing (equal). Marco Endrizzi: Writing – review & editing (equal). Peter Robert Thomas Munro: Writing – review & editing (equal). David Bate: Supervision (equal); Writing – review & editing (equal). Alessandro Olivo: Conceptualization (lead); Data curation (equal); Investigation (equal); Methodology (lead); Project administration (lead); Resources (lead); Supervision (lead); Writing – original draft (equal); Writing – review & editing (equal).

The data that support the findings of this study are available from the corresponding author upon reasonable request.