1. Introduction

Existing models of dark-field imaging primarily include a phenomenological model [1–3], referred to here as convolutional blurring, as well as more rigorous wave-optics formulations [4–6]. While numerical wave-optics models provide accurate results, they are generally computationally intensive, requiring large high-resolution 2D meshes and significant execution time to propagate the radiation field from one plane of interest to another. They are also restricted to simulation of one realization of the random scattering medium at a time, so it may be the case that multiple runs are needed to yield statistically meaningful results [6]. While there is no doubt that wave-optics models can generate useful field profiles at various locations in an optical path, they are not necessarily well suited for providing insight into the basic physics of a problem, particularly when a scattering medium is involved. In contrast, the statistical optics model formulated in Part 1 [7] is based on fundamental principles associated with grating diffraction and propagation of the ensemble-average mutual intensity through the system to the detection plane where the Talbot fringe forms. For the case in which the scattering medium comprises a random ensemble of microspheres having an arbitrary concentration, we show that the resulting detected fringe visibility can be calculated in a numerically efficient manner while also providing physical insight into the underlying optical phenomena. The primary focus here centers on a monochromatic point-source configuration that can serve as the basis for modeling a more general polychromatic extended-source setup.

Part 2 of the paper is organized as follows. In Section 2 we describe a 28-keV interferometer setup that is used for our simulation examples. Section 3 provides a summary of the statistical optics model. Section 4, which forms the bulk of the paper, contains several simulation examples in which the detected fringe visibility is calculated as various parameters related to the object are varied. These results are validated in Section 5 with a numerical Fourier optics model that is computationally much more demanding but does provide support for the statistical optics formulation. Additionally, the Fourier optics model provides images of the Talbot fringe showing how the visibility is degraded by varying amounts of scatter-induced speckle. In Appendix A we extend the Fourier optics model to investigate a thick lossless medium that exhibits strong scattering and significant diffraction within the medium, and show that even in this case the statistical optics model based on the projection approximation provides accurate results.

2. System layout

Figure 1 shows a schematic representation of a typical three-grating Talbot-Lau shearing interferometer using a standard X-ray tube source. For purposes of simulation, we consider a specific interferometer configuration based on a system built in our lab with the following properties:

 

Fig. 1. Diagram of a Talbot-Lau X-ray interferometer. The object is shown in front of G$_1$, but it can also be placed after G$_1$.

Download Full Size | PDF

3. Summary of the dark-field fringe visibility model

In this section we review the key theoretical aspects of the statistical optics model as formulated in Part 1 [7]. We begin by noting that the Talbot fringe pattern originates from interference among the G$_1$ diffraction orders. The various spatial frequency terms of the fringe are indexed by an integer $q = 1, 2, \ldots, 2M$, where $M$ is the highest diffraction order of significance. Each fringe component is weighted by a corresponding coherence factor, $\gamma _q$, that depends on the statistical properties of the scattering object and its location within the interferometer. The net Talbot fringe is given by

The coherence coefficients in Eq. (1) are given by the product of two terms,

The second term in Eq. (5), $\boldsymbol {\mu }_{d} (\Delta x_q)$, pertains to the ensemble-average spatial coherence of the zeroth-order wavefront as observed in the detection plane. Upon traversing a lossless scattering medium, the zeorth-order beam acquires a random phase perturbation which is described statistically by the following field autocorrelation function evaluated in the output plane of the medium,

In the limit of low particle number density ($f \ll 1$) the normalized phase autocovariance function for microspheres of radius $R$ is well approximated by

Finally, the fringe signal associated with a detector pixel is found by convolving (along the $x$-direction) the Talbot intensity pattern of Eq. (1) with an ideal square-wave G$_2$ grating transmission function having full zero-to-one modulation. The visibility, $V_d$, of this detected signal is found by dividing the peak-to-peak amplitude by two times the average signal. The dark-field effect corresponds to a change in $V_d$ when a scattering object is introduced into the beam path.

4. Simulation examples

For simulation purposes, we focus on one type of scattering object with statistical properties that are particularly amenable to analytical analysis, namely a suspension of monodisperse microspheres (i.e., a random spatial distribution of small non-overlapping hard spheres having a fixed radius and a uniform refractive index that is different from the background). The microspheres are assumed to be nonabsorbing, so the medium is lossless. The relative concentration of microspheres within a given sample volume is an important parameter, quantified by the volume fraction $f$. By using the projection approximation, a 3D volume of scattering microspheres can be represented as a 2D random phase screen. The interferometer layout is our 28 keV configuration described in Section 2. The emitter is taken to be a monochromatic point source radiating at the design energy of 28 keV. The Talbot fringe is found by using $M = 15$ (i.e., by using G$_1$ diffraction orders between $\pm 15$).

In the following subsections, we provide examples of the dark-field signal ($V_d$) computed as a function of various parameters related to this phase screen representation of the object, along with its placement within the interferometer, as well as properties of the monodisperse collection of microspheres, including the particle radius, the standard deviation of the projected phase profile, the volume fraction of microspheres, and the object thickness. Again, for all of these examples we assume the medium is sufficiently thin that $\rho _q \sim 1$ (see Eq. (5)).

4.1 Visibility vs. size and phase of scattering microspheres: $V_d(R)$ and $V_d(\sigma _{\phi })$

As a first set of examples, we consider a dilute suspension of microspheres and examine the change in Talbot fringe visibility as a function of particle radius and the standard deviation $\sigma _{\phi }$ of the projected phase front, $\phi (\alpha,\beta )$, at the output plane of the object. A dilute sample has a small volume fraction of particles ($f << 1$), so the mean value of $\phi (\alpha,\beta )$ is close to zero, which implies the projected phase autocorrelation function and the projected phase autocovariance function are essentially identical and, as discussed in Appendix A of Part 1 [7], can be well approximated by a Gaussian function as given by Eq. (9). As an aside, we note that for dilute suspensions of randomly oriented particles having a shape other than spherical (e.g., cylindrical), the autocovariance function can be modified to $\exp [-(\Delta \alpha / R)^{2H}]$ where $H$ is the Hurst exponent $(0 < H < 1)$ that varies with particle shape and/or aspect ratio, and $R$ now represents an average particle dimension as discussed by Yashiro, et al. [4]. Here, for simplicity, we focus solely on spherical particles, so from Eqs. (6) and (9) the normalized autocorrelation function of the scattered field at the output of the diffusing screen becomes,

 

Fig. 2. Plots of the detected monochromatic visibility for a dilute suspension of monodisperse microspheres: (a) versus microsphere radius for three different values of phase standard deviation, and (b) versus phase standard deviation for three different values of microsphere radius. The phase screen is positioned adjacent to G$_1$.

Download Full Size | PDF

In Fig. 2(a) we see that the detected visibility signal goes up with increasing particle radius as would be expected because larger particles yield a broader coherence function, so the scattered diffraction orders become more correlated with increasing radius. For the curves corresponding to $\sigma _{\phi }$ = 0.5 and 1.0, the visibility displays a minimum non-zero plateau for small particles. This behavior can be understood by first recalling that $\exp (-\sigma _{\phi }^2)$ equals the fraction of specular transmitted power. For $\sigma _{\phi } < 2.0$ there is a specular component that yields a Talbot fringe having a minimum visibility regardless of the degree of correlation between the scattered portion of the diffraction orders. For sufficiently small particles, the scattered diffraction orders are completely uncorrelated, so the detected visibility remains pinned at a constant value that depends on $\sigma _{\phi }$. As the particles become larger, the coherence function broadens, so that a portion of the scattered radiation now contributes to the Talbot fringe, thereby improving its visibility. The particle radius at which this transition occurs depends on the location of the scattering object, but for an object adjacent to G$_1$ this radius is approximately equal to $p_2/4$ (i.e., one quarter of the Talbot fringe period, or 1.5 $\mu$m in this example). This result is consistent with the findings reported in [9]. For $\sigma _{\phi } > 2.0$ there is little or no specular radiation, so the visibility is zero for small particles, but a fringe can eventually form from the scattered radiation for sufficiently large particles, with the transition point depending on both $\sigma _{\phi }$ and the object location.

The plots in Fig. 2(b) show the visibility falls rapidly with increasing projected phase standard deviation, particularly for microspheres having a radius less than about 3-5 $\mu m$. For larger diameter microspheres, reduction of the visibility requires a somewhat stronger phase modulation characterized by larger values of $\sigma _{\phi }$. In general, $\sigma _{\phi } \propto \Delta n$ (index contrast between the particle and its background), while the phase variance $\sigma ^2_{\phi } \propto TR$ (sample thickness times particle radius).

4.2 Visibility vs. scattering object location: $V_d(z_s)$

Next we calculate the monochromatic detected visibility as a function of the scattering object position $z_s$, with results shown in Fig. 3 for three different values of the projected phase standard deviation ($\sigma _{\phi } =$ 0.5, 1 and 3 radians) and three particle radii ($R =$ 1, 3 and 10 $\mu$m). As in the previous subsection, these simulations are based on our 28 keV system layout operating at the design energy, under the assumption of a dilute suspension of monodisperse microspheres. From these results it is clear that for a given scattering object, the visibility can vary significantly based solely on the object’s position within the interferometer. This finding is consistent with previous work [10]. A similar effect is known to occur with the position sensitivity of a differential phase contrast measurement [11], although the visibility response displays a richer range of variation depending on the size and statistical properties of the scattering particles.

 

Fig. 3. Plots of the detected monochromatic visibility as a function of the object location for $\sigma _{\phi }$ equal to: (a) 0.5, (b) 1.0, and (c) 3.0 radians. Curves for microsphere radii of 1, 3 and 10 $\mu$m are shown in each plot. The microsphere suspension is assumed to be dilute ($f << 1$) so the phase autocovariance function is well-approximated as Gaussian with a correlation length equal to the particle radius.

Download Full Size | PDF

4.3 Visibility vs. volume fraction of scattering microspheres: $V_d(f)$

All of the simulation examples shown so far have been based on scattering from a dilute microsphere suspension in which the mean value of the projected phase is close to zero and the phase autocovariance function is well-approximated by a simple Guassian function. We now remove this constraint and utilize a more general model for a monodisperse ensemble of microspheres with arbitrary volume fraction, up to a maximum of $f_\textrm {max} \sim 0.64$ for randomly packed hard spheres. Again, we utilize Eq. (6) for the scattered field autocorrelation function, but now we must work with generalized expressions for the phase variance and the normalized phase autocovariance function as given by Eqs. (11) and (12), respectively.

To generate sample plots of $V_d(f)$, we begin with numerical evaluation of $\gamma _e$ for an initial value of volume fraction, $f$, and a specific particle radius, $R$. Then, for given values of sample thickness $T$, particle refractive index contrast $\Delta n$, and the wavenumber $k = 2\pi /\lambda$, we determine $\sigma _{\phi }$ and $c_{\phi }(\Delta \alpha )$, from which follows $\boldsymbol {\mu }_{s}(\Delta \alpha )$. For a given position of the scattering object within the interferometer, we can scale $\boldsymbol {\mu }_{s}(\Delta \alpha )$ to obtain the field autocorrelation function in the detection plane $\boldsymbol {\mu }_{d}(\Delta x)$, and subsequently calculate the Talbot fringe and its detected visibility. We then incrementally step the value of $f$ and repeat the calculations. The overall computational load is minimal, with a curve of $V_d(f)$ evaluated at 60 points taking about 7 seconds to execute in Matlab.

For the example here we consider the case in which the microsphere beads are made from silica (amorphous SiO$_2$, $\delta = 1-n$ = 5.8e-07 @ 28 keV) and are immersed in glycerin (C$_3$H$_8$O$_3$, $\delta _0 = 1 - n_0$ = 3.6e-07 @ 28 keV). This mimics the sample used by Gkoumas et al. for experimental measurements of the linear dark-field extinction coefficient versus the microsphere volume fraction [12]. Figure 4 shows the phase standard deviation and the detected visibility as functions of $f$ for sphere radii of 1.0, 3.0 and 10.0 $\mu$m; the sample thickness is 10 mm. In Fig. 4(a) the phase standard deviation exhibits a maximum at about $f = 0.33$ and then decreases for higher concentrations as the spheres become more closely packed. In Fig. 4(b) the detected visibility drops rapidly with increasing $f$, but the behavior for $f > 0.2$ is more complicated, with the visibility for $R = 1.0$ $\mu$m dipping slightly below $V_d = 0.2$ when $f \sim 0.3$, and then actually going back up somewhat as the concentration continues to increase. The corresponding plots of the dark-field extinction coefficient ($\epsilon _d = -\ln (V_d)/T$) versus $f$ are shown in Fig. 5. We see $\epsilon _d$ reaches a peak at a value of $f$ that increases with increasing microsphere radius. These findings are similar to those in Ref. [12] (Fig. 2) for $R = 0.93$ and 3.88 $\mu$m, although our simulation predicts somewhat higher values of $\epsilon _d$ compared to those reported in [12]; however, our system configuration is different, which may account for this discrepancy.

 

Fig. 4. Plots of: (a) the phase standard deviation $\sigma _{\phi }(f)$, and (b) the monochromatic detected visibility $V_d(f)$, as a function of particle volume fraction for three different values of microsphere radius. The scattering medium, with thickness $T$ = 10 mm, is positioned immediately in front of G$_1$. The microspheres are taken to be silica, immersed in a background of glycerin (E = 28 keV).

Download Full Size | PDF

 

Fig. 5. The visibility curves of Fig. 4(b) converted into plots of the dark-field extinction coefficient versus microsphere volume fraction.

Download Full Size | PDF

4.4 Visibility vs. particle radius: modified version of $V_d(R)$

In Section 4.1 we investigated the visibility versus particle radius for specified values of $\sigma _{\phi }$. However, in general $\sigma _{\phi }$ depends on the particle radius, so it is of practical interest to use our analysis for arbitrary $f$ and incorporate the $\sigma ^{2}_{\phi }(R)$ dependence into the visibility calculation while keeping other sample parameters fixed (such as the particle index contrast, particle volume fraction, and sample thickness) as would likely be the case in an experiment. Figure 6 shows the phase variance and the detected visibility versus particle radius for the case of silica beads immersed in glycerin with fixed values of volume fraction. The sample thickness is 10 mm.

 

Fig. 6. Plots of: (a) $\sigma _{\phi }^2(R)$ and (b) $V_d(R)$ for three different values of microsphere volume fraction. The sample is 10 mm thick and is located immediately in front of G$_1$. The microspheres are taken to be silica immersed in glycerin (E = 28 keV).

Download Full Size | PDF

In this case, we see that the phase variance is linearly proportional to the particle radius, and the detected visibility displays a minimum when the particle radius is in the vicinity of the interferometer’s primary correlation length ($\xi _c = 2.3$ $\mu$m). These results are consistent with previous wave-optics simulations [6,13]. For high particle concentrations, the visibility can vanish over a range of particle size.

4.5 Visibility vs. scattering object thickness: $V_d(T)$

Of practical importance is the dependence of the detected fringe visibility on thickness of the scattering medium. Figure 7 shows both the phase variance and the detected visibility versus thickness for the case of silica microspheres having a radius of 3.0 $\mu$m immersed in glycerin. The phase variance increases linearly with thickness. The detected visibility decreases in an approximate exponential fashion with thickness; however, closer inspection shows it is not a pure exponential because the shape of the monochromatic Talbot fringe evolves from a square wave for thin media to a sine wave as the thickness increases, and this affects the resulting visibility. However, a pure exponential decay does arise if the visibility is calculated by using the amplitude of the fundamental Fourier component of the detected signal instead of the amplitude of the detected signal itself. In practice, with a polychromatic extended source, the distinction between these two methods for evaluating the detected visibility may be minimal. Exponential visibility decay with thickness for a uniform medium (similar to Beer’s law for attenuation) means the dark-field extinction coefficient ($\epsilon _d = -\ln (V_d)/T$) is a constant. In the more general case of a spatially varying scattering medium, the 3D distribution of the extinction coefficient can be determined from computed tomography in a way entirely analogous to the standard case of the linear attenuation coefficient [2].

 

Fig. 7. Plots of: (a) $\sigma _{\phi }^2(T)$ and (b) $V_d(T)$ for three different values of microsphere volume fraction. The sample is located immediately in front of G$_1$. The microspheres (with $R$ = 3.0 $\mu$m) are taken to be silica immersed in glycerin (E = 28 keV).

Download Full Size | PDF

5. Fourier optics simulations

To validate the results of our statistical optics model, we have undertaken numerical Fourier optics simulations that seek to replicate previous findings for: (1) visibility versus the phase standard deviation of the diffusing screen (Fig. 2(b)), and (2) visibility versus the diffusing screen position within the interferometer (Fig. 3(b)). As was the case with the statistical optics models for these two examples, the screen is assumed to mimic a dilute suspension of scattering microspheres, which upon projection displays Gaussian phase statistics as well as a Gaussian phase autocovariance function $c_{\phi }(r_{\perp }) = \exp [-(r_{\perp }/R)^2]$ where $R$ is the correlation length of the wavefront phase profile, which here equals the microsphere radius. The object is represented by a thin random phase screen, thereby replicating the approach used for the statistical optics simulation examples.

5.1 Visibility vs. phase of diffusing screen: $V_d(\sigma _{\phi })$

We begin with a monochromatic, on-axis point source having an energy of 28 keV. The field incident on G$_1$ is therefore

The field leaving the G$_1$ and phase-screen combination is

To calculate the detected fringe visibility when phase stepping G$_2$, twenty periods of the Talbot fringe pattern are extracted from the center portion of the $I_T$ array and convolved along the $x$-direction with an ideal square-wave G$_2$ transmittance function. This detection convolution is a spatial averaging process performed on one realization of the Talbot pattern. However, since the speckled Talbot pattern is assumed to be spatially ergodic, the result should be equivalent to the ensemble average used in the statistical optics formulation. Of course, this equivalence requires that the area of the Talbot pattern subject to convolution be large enough for the spatial average to properly reflect an ensemble average. Our investigation indicates that twenty periods is adequate, although a larger area (requiring more computation time) would help further reduce small variations observed when the convolution is taken over different random realizations of the diffusing screen. This analysis is repeated for a sequence of $\sigma _{\phi }$ values, yielding a set of $V_d(\sigma _{\phi })$ data points which are shown as solid circles in Fig. 8, overlaid on the curves previously generated from the statistical optics model. The apparent good agreement supports the validity of our statistical optics model. The small difference between the circle data points and the curves is largely attributed to the slight departure from ergodicity mentioned above for the Fourier optics model.

 

Fig. 8. Plots of $V_d(\sigma _{\phi })$ for three separate dilute suspensions of monodispere microspheres ($R$ = 1, 3, and 10 $\mu$m). The corresponding diffuse phase screen is positioned adjacent to G$_1$. The graph lines are from the previous statistical optics analysis (Fig. 2(b)), while the solid circles are for the Fourier optics simulation.

Download Full Size | PDF

It is also instructive to see how the Talbot fringe images themselves change as $\sigma _{\phi }$ increases. In Fig. 9 we display example images for the case of $R = 3.0$ $\mu$m using four different values of $\sigma _{\phi }$. The nominal Talbot pattern without scattering ($\sigma _{\phi } = 0$) is shown in Fig. 9(a); it has a detected visibility of 0.97 (note: the X-ray Talbot fringe is an approximate square wave with a visibility of one, but the sharp edges of this fringe pattern have slight tails that, upon convolution with G$_2$, produce a detected triangular-wave signal that does not quite reach zero, hence $V_d = 0.97$). As the strength of the random phase screen increases, the degree of speckle goes up and the detected fringe visibility drops.

 

Fig. 9. Fourier optics simulation of monochromatic Talbot fringe images corrupted by speckle, which reduces the detected fringe visibility. They correspond to the following points on the $R = 3.0$ $\mu$m curve of Fig. 8: (a) $\sigma _{\phi } = 0$ $(V_d = 0.97)$, (b) $\sigma _{\phi } = 0.6$ $(V_d = 0.81)$, (c) $\sigma _{\phi } = 1.0$ $(V_d = 0.58)$, and (d) $\sigma _{\phi } = 2.0$ $(V_d = 0.14)$. All four images are plotted on the same linear intensity scale having a maximum value of six times the mean intensity.

Download Full Size | PDF

5.2 Visibility vs. scattering object location: $V_d(z_s)$

For additional validation, we apply a Fourier optics model to find the Talbot fringe visibility as a function of the diffusing screen location within the interferometer. The approach is similar to that of the previous section, but here we first evaluate the point-source spherical wavefront at either the diffusing screen or the G$_1$ grating, whichever comes first, and then multiply the field by the corresponding amplitude transmittance function. We then Fresnel propagate the field to the second object, multiply by its amplitude transmittance, and then finally propagate once again to the detection plane. Only one realization of the random diffusing screen is utilized. The results are shown in Fig. 10. Again, good correspondence with the statistical optics model is observed.

 

Fig. 10. Plots of $V_d(z_s)$ for a dilute suspension of microspheres, represented by a Gaussian random-phase screen having $\sigma _{\phi } = 1.0$. Three different microsphere radii are modeled ($R$ = 1, 3, and 10 $\mu$m). The graph lines are from the previous statistical optics analysis (Fig. 3(b)); the solid circles are for the Fourier optics simulation.

Download Full Size | PDF

We conclude this section with a few comments regarding the relative computational complexity of the two simulation approaches. On a high-end dual-CPU workstation (Intel Xeon Gold 6154 3.0 GHz processors), we find that the Fourier optics method requires about 45 seconds per visibility data point, while the statistical optics approach only takes 0.12 second per point. Also, the Fourier optics method requires significantly more memory because it uses multiple large 2D grids to represent the field at various locations along the propagation path; however, it does provide images of the speckled Talbot fringe, which are helpful for visual interpretation of the results.

6. Conclusion

Dark-field X-ray imaging with a Talbot-Lau interferometer has the ability to detect and image small-angle scattering from a sample with high sensitivity. As its importance for a variety of applications continues to grow, improved numerical simulation techniques will prove helpful in advancing the technology. In this paper we have demonstrated the versatility of our statistical optics model to calculate the detected fringe visibility for a scattering object consisting of randomly positioned dielectric microspheres. The model is based on fundamental properties of grating diffraction and propagation of the ensemble-average mutual intensity from the scattering object to the detection plane. It therefore provides physical insight into the underlying optical phenomena associated with the dark-field signal.

The results of several example monochromatic simulations based on a 28-keV system configuration are included. These simulations explore the dependence of the visibility signal on various object properties (microsphere radius, volume fraction, and object thickness) as well as for changes in the object location within the interferometer. We replicate a subset of these results via numerical Fourier optics simulation and find close agreement, which serves to validate the statistical optics model. Numerical implementation of the statistical optics model is found to execute significantly faster than a Fourier optics model. While the emphasis here is on dark-field simulation of known scattering objects, the numerical efficiency of our model suggests it could also be useful for fitting of experimental results to estimate properties of the scattering particles such as their average size, the refractive index contrast, or volume fraction. Lastly, here in Appendix A of Part 2 we extend the Fourier optics model to show that the statistical optics model (in which the projection approximation is used to find the phase perturbation introduced by a lossless object) remains valid even for a strongly scattering thick object for which diffraction within the medium becomes significant.

Appendix A. The projection approximation: range of validity

The basic idea of the projection approximation is illustrated in Fig. 11. Parallel rays propagate through an inhomogeneous scattering medium having a random spatial distribution of refractive index, for example a suspension of small randomly located microspheres. The scale sizes associated with the inhomogeneities are assumed to be much larger than the wavelength (similar to the assumption typically made for turbulent media, in contrast to turbid media in which the scale sizes are comparable to or smaller than the wavelength). In this approximation the rays incur no bending, but each ray does acquire an accumulated phase delay along its path. Therefore, this approach comprises a purely geometrical model with diffraction being neglected during propagation within the scattering medium itself.

 

Fig. 11. Illustration of the plane-wave projection approximation as applied to a random distribution of microspheres. Rays propagating parallel to the $z$-axis yield a 2D projected phase distribution at the output plane. We generally neglect the uniform background phase ($\phi _0 = kn_0T$, where $k$ is the free-space wavenumber) and retain only the spatial phase modulation $\phi (x,y)$ at $z = T$.

Download Full Size | PDF

From simple diffraction considerations, one might assume this model is accurate as long as the smallest Fresnel number for the particles in the distribution, calculated relative to the output plane, remains much greater than unity. In our case, the minimum value of $N_F$ arises for a sphere near the input plane, so that $N_F = R^2/\lambda T$. For example, if the sphere radius is $R = 1$ $\mu$m and $\lambda =$ 4.4 x 10$^{-5}$ $\mu$m (E = 28 keV), and we desire $N_F \ge 10$, then we may conclude that the maximum thickness of the scattering region for which the projection approximation remains valid is $T = R^2/(\lambda N_{F}) =$ 2.3 mm. However, the detected fringe visibility does not depend directly on the projected phase profile, but rather on the autocorrelation function of the scattered field leaving the medium (see Eqs. (1) and (5)). It has also been theoretically established for propagation through lossless turbulent media that the field autocorrelation function (or mutual intensity) calculated using the near-field geometrical projection model actually remains valid for distances well beyond what a simple Fresnel number estimate suggests should be the case ([16], Sec. 8.5.6). This finding stems from a wave-optics based propagation analysis for a statistically isotropic, lossless turbulent medium (including, but not necessarily specific to, the atmosphere) subject to the Rytov approximation, which like the Born approximation, assumes weak scattering. Even then, theories have been developed for the strong-fluctuation regime, and as noted by Goodman ([16], p. 405), “The results of these various (strong-fluctuation) theories yield a remarkable and important result. The predictions of all methods for finding the mutual coherence function of a propagating wave yield the same result, namely, the result we obtained with the Rytov approximation.” It therefore seems quite reasonable that these same findings are applicable to small-angle X-ray scattering. Thus the dark-field signal calculated using the projection approximation is valid for much greater object thicknesses than might have been originally assumed.

To test this assumption, we utilize the Fourier optics model of Sec. 5 but now simulate a thick random medium by using the beam propagation method (BPM) [17] and then compare the visibility results to those of our statistical optics model based on the plane-wave projection approximation. The two models are illustrated in Fig. 12. For the projection approximation, the entire object is replaced by a single phase screen located at the object output as shown in Fig. 12(a). On the other hand, the BPM approach, illustrated in Fig. 12(b), uses a sequence of $N$ uncorrelated phase screens with equal spacing $\Delta z$ along the propagation direction, and the field is Fresnel propagated from one screen to the next so that diffraction within the medium is taken into account. Each screen represents the geometric phase projection of the object over a slice thickness $\Delta z$. The total phase variance of the object, $\sigma _{\phi }^2$, is split equally among the individual screens so each screen has a phase variance of $\sigma _{\phi }^2/N$, consistent with the screens being uncorrelated (allowed because their separation $\Delta z$ is taken to be greater than the medium’s correlation length which is on the order of the average particle size). For these simulations the object is placed adjacent to, but in front of G$_1$, and is assumed to range in thickness from 0 to 200 mm. The medium consists of silica particles with $R$ = 1 $\mu$m and $f$ = 0.01, immersed in a background of glycerin, so at 28 keV the refractive index contrast is $\Delta n$ = 2.2 x 10$^{-7}$. For the BPM model $\Delta z$ = 1 mm, and the transverse object dimensions are approximately 150 $\mu$m x 150 $\mu$m (see Sec. 5 for mesh details). The detected visibility versus object thickness is shown in Fig. 12(c). Both models yield virtually identical results, even in the strong-scattering regime where the visibility becomes small. The BPM model took about two hours to complete, while the statistical optics model executed in roughly one minute.

The extent to which diffraction plays a role can be quantified from the intensity variations in the BPM model. Figure 13(a) shows probability density functions (PDF) of the intensity for various object thicknesses, computed in Matlab as properly scaled intensity histograms. The incident field has a delta-function PDF, but after just a few millimeters inside the medium a partially developed speckle pattern forms in accordance with the spread of intensities for the PDF at $z=10$ mm. As propagation in the medium continues, the PDF broadens further and the degree of speckle increases. For z < 50 mm the PDF is well represented by a log-normal distribution ([16], Sec. 8.4.4). For z > 50 mm the beam eventually becomes highly speckled, and by $z=200$ mm the intensity PDF appears to resemble a so-called $K$-distribution [18,19]. The degree of speckle is quantified by its contrast, which equals the standard deviation of the intensity fluctuations divided by the mean intensity ($\sigma _I / \bar {I}$), and is plotted in Fig. 13(b). The contrast rises approximately linearly until it reaches a value of one at $z =85$ mm, which is generally considered to be fully developed speckle ([20], Sec. 3.2.1). However, the contrast continues to rise more slowly up to a peak value of about 1.15 at $z = 150$ mm. With continued propagation in an even thicker medium, we expect the contrast would drop slightly and asymptotically approach one, corresponding to negative-exponential intensity statistics. This type of behavior is well-established in the theory of propagation in random media [21–23]. The primary point here is that even though we are deep into the diffraction-dominated regime with strong scattering, the projection approximation for a thick phase-only medium provides the correct solution for the dark-field signal. Determining specific conditions under which the projection approximation may begin to break down for purposes of dark-field analysis is beyond the scope of this paper, but clearly its range of applicability covers object thicknesses of practical interest.

 

Fig. 12. (a) Projection approximation for a thick phase-only object; a single planar phase screen is located at the object output. (b) Beam propagation method (BPM) in which the object is represented as a sequence of phase screens, and Fresnel propagation takes the field from the output of one screen to the input of the next. (c) Example comparison of the detected dark-field fringe visibility versus object thickness for both models. The medium consists of silica particles immersed in a background of glycerin.

Download Full Size | PDF

 

Fig. 13. BPM numerical simulation results: (a) intensity probability density functions for various object thicknesses, and (b) speckle contrast as a function of object thickness.

Download Full Size | PDF

Funding

Transportation Security Administration (HSTS04-17-C-CT7224).

Acknowledgments

We would like to gratefully acknowledge review of a preliminary manuscript by Prof. Joseph W. Goodman, as well as many helpful discussions with members of the research team at Stanford including Max Yuen, Bill Aitkenhead, Yao-Te Cheng, Paul Hansen, Ludwig Galambos, Niharika Gupta, Ching-Wei Chang, and George Herring. In addition, we thank the anonymous reviewers for their thoughtful insight and suggestions that helped improve the paper.

Disclosures

The authors declare no conflicts of interest.

Data availability

Data underlying the results presented in this paper are not publicly available at this time but may be obtained from the authors upon reasonable request.

Supplemental document

See Supplement 1 for supporting content.

References

1. Z.-T. Wang, K.-J. Kang, Z.-F. Huang, and Z.-Q. Chen, “Quantitative grating-based x-ray dark-field computed tomography,” Appl. Phys. Lett. 95(9), 094105 (2009). [CrossRef]  

2. M. Bech, O. Bunk, T. Donath, R. Feidenhans, C. David, and F. Pfeiffer, “Quantitative x-ray dark-field computed tomography,” Phys. Med. Biol. 55(18), 5529–5539 (2010). [CrossRef]  

3. M. Strobl, “General solution for quantitative dark-field contrast imaging with grating interferometers,” Sci. Rep. 4(1), 7243 (2015). [CrossRef]  

4. W. Yashiro, Y. Terui, K. Kawabata, and A. Momose, “On the origin of visibility contrast in x-ray Talbot interferometry,” Opt. Express 18(16), 16890–16901 (2010). [CrossRef]  

5. S. K. Lynch, V. Pai, J. Auxier, A. F. Stein, E. E. Bennett, C. K. Kemble, X. Xiao, W.-K. Lee, N. Y. Morgan, and H. H. Wen, “Interpretation of dark-field contrast and particle-size selectivity in grating interferometers,” Appl. Opt. 50(22), 4310–4319 (2011). [CrossRef]  

6. A. Malecki, G. Potdevin, and F. Pfeiffer, “Quantitative wave-optical numerical analysis of the dark-field signal in grating-based x-ray interferometry,” Europhys. Lett. 99(4), 48001 (2012). [CrossRef]  

7. J. P. Wilde and L. Hesselink, “Statistical optics modeling of dark-field scattering in X-ray grating interferometers: Part 1. Theory,” Opt. Express 29(25), 40891–40916 (2021). [CrossRef]  

8. S. Torquato, T. M. Truskett, and P. G. Debenedetti, “Is random close packing of spheres well defined?” Phys. Rev. Lett. 84(10), 2064–2067 (2000). [CrossRef]  

9. V. Revol, I. Jerjen, C. Kottler, P. Schütz, R. Kaufmann, T. Lüthi, U. Sennhauser, U. Straumann, and C. Urban, “Sub-pixel porosity revealed by x-ray scatter dark field imaging,” J. Appl. Phys. 110(4), 044912 (2011). [CrossRef]  

10. F. Prade, A. Yaroshenko, J. Herzen, and F. Pfeiffer, “Short-range order in mesoscale systems probed by x-ray grating interferometry,” Europhys. Lett. 112(6), 68002 (2015). [CrossRef]  

11. T. Donath, M. Chabior, F. Pfeiffer, O. Bunk, E. Reznikova, J. Mohr, E. Hempel, S. Popescu, M. Hoheisel, M. Schuster, J. Baumann, and C. David, “Inverse geometry for grating-based x-ray phase-contrast imaging,” J. Appl. Phys. 106(5), 054703 (2009). [CrossRef]  

12. S. Gkoumas, P. Villanueva-Perez, Z. Wang, L. Romano, M. Abis, and M. Stampanoni, “A generalized quantitative interpretation of dark-field contrast for highly concentrated microsphere suspensions,” Sci. Rep. 6(1), 35259 (2016). [CrossRef]  

13. A. Velroyen, M. Bech, A. Malecki, A. Tapfer, A. Yaroshenko, M. Ingrisch, C. C. Cyran, S. D. Auweter, K. Nikolaou, M. Reiser, and F. Pfeiffer, “Microbubbles as a scattering contrast agent for grating-based x-ray dark-field imaging,” Phys. Med. Biol. 58(4), N37–N46 (2013). [CrossRef]  

14. R. Gonzalez, R. Woods, and S. Eddins, Digital Image Processing Using MATLAB (Pearson Prentice Hall, 2004).

15. J. Schmidt, Numerical Simulation of Optical Wave Propagation with Examples in MATLAB, Press monograph (SPIE, 2010).

16. J. W. Goodman, Statistical Optics (John Wiley & Sons, 2015), 2nd ed.

17. M. Mansuripur, Classical Optics and Its Applications (Cambridge, 2009), chap. 32, 2nd ed.

18. E. Jakeman and P. N. Pusey, “Significance of K distributions in scattering experiments,” Phys. Rev. Lett. 40(9), 546–550 (1978). [CrossRef]  

19. G. Parry and P. N. Pusey, “K distributions in atmospheric propagation of laser light,” J. Opt. Soc. Am. 69(5), 796–798 (1979). [CrossRef]  

20. J. W. Goodman, Speckle Phenomena in Optics: Theory and Applications (SPIE, 2020), 2nd ed.

21. C. Macaskill and T. E. Ewart, “Computer simulation of two-dimensional random wave propagation,” IMA J. Appl. Math. 33(1), 1–15 (1984). [CrossRef]  

22. B. J. Uscinski, “Analytical solution of the fourth-moment equation and interpretation as a set of phase screens,” J. Opt. Soc. Am. A 2(12), 2077–2091 (1985). [CrossRef]  

23. E. Jakeman and K. D. Ridley, Modeling Fluctuations in Scattered Waves (Taylor & Francis, 2006), chap. 9.