Abstract
In order to quantify the effect of phase contrast on X-ray image formation, the theory of statistical decision making has been applied to a binary classification task between two signals known exactly, namely, a phase-contrast image (that combines both the absorption and phase contrast) and the corresponding hypothetical pure absorption image that would be obtained under the same imaging conditions but without diffraction/refraction effects. The signal-to-noise ratio (SNR) for two widely used observers, including the ideal observer (also known as the prewhitening matched filter) and a non-ideal observer (the non-prewhitening matched filter) has been estimated in the case of in-line phase-contrast imaging, thus providing a figure-of-merit for the optimisation of the imaging conditions. A broad class of edge objects has been investigated and simple analytical expressions for the corresponding SNRs have been obtained and discussed.
© 2014 Optical Society of America
1. Introduction
In the recent years, several papers have been published on quantitative assessment of X-ray phase-contrast image quality, in particular, in the context of the effect of phase contrast on the detectability of small weakly absorbing features inside objects [1–7]. Several figures of merit (FoMs) have been proposed for the quantitative characterisation of the image quality, including contrast, contrast-to-noise ratio (CNR), CNR per unit dose etc. Some authors [5,6] utilized a Bayesian approach to the imaging tasks [8], with the ideal (and/or non-ideal) observer’s SNR as the FoM. This latter approach is very generic and covers many important cases, including (but not restricted to) the detectability task (both in the raw [6] and phase-retrieved images [5]) and the reconstruction task (for example, extraction of the projected thicknesses of the different materials constituting the object, from the X-ray images).
As soon as an imaging task is formulated and a quantitative FoM is chosen, one is usually interested in optimising imaging conditions in order to maximise the FoM. These imaging conditions may include, for example, the X-ray energy (or energies), the imaging geometry, the source and the detector characteristics.
In the last two decades, several phase-contrast imaging methods have been rapidly developing and have attracted the attention of the bio-medical community. A recent review of these methods can be found in [9]. In the present paper, our analysis will be restricted to the propagation-based phase-contrast imaging. Nevertheless, the mathematical methods used here can find applications in quantitative analysis of other phase-contrast imaging modalities.
The paper is organized as follows. In Section 2 an imaging task is formulated and the corresponding FoMs are introduced. In Section 3 a linear approximation for propagation-based image formation is reviewed and relevant expressions for intensity distributions are introduced. In Section 4 simple analytical expressions for the SNRs corresponding to three different observers are derived for a broad class of edge objects. Discussion of the results and conclusions are given in Section 5.
2. Figure-of-merit
When assessing the quality of a phase-contrast image, usually one is trying to find out whether the phase contrast expected for the feature of interest is visible/detectable in the image, in the presence of noise. This is equivalent to deciding whether this phase-contrast image is distinguishable from the corresponding pure absorption image of the same feature, under the same imaging conditions. In the theory of statistical decision making [8] this problem can be formulated as a binary classification task between two signals known exactly, namely, a phase-contrast image (that combines both the absorption and the phase contrast) and the corresponding hypothetical pure absorption image that would be obtained under the same imaging conditions but without diffraction/refraction effects. By choosing a proper decision maker (observer) which returns a single decision variable for each input image, one then is able to calculate the corresponding signal-to-noise ratio (SNR), which is a measure of the overlap of the decision variable outcomes that would be obtained for the ensemble of noisy images in the above two classes. In this paper, we estimate the SNR for two widely used observers [8], namely the ideal observer (also known as the prewhitening matched filter, PWMF) and a non-ideal observer, the non-prewhitening matched filter (NPWMF), which is known to better match a real human observer. We restrict our consideration to the quantitative evaluation of the diffraction/refraction contrast in a phase-contrast image due to small weakly absorbing features inside the object.
Designating the phase-contrast images with I and the corresponding hypothetical pure absorption images with Iabs, the SNR of the ideal observer, SNRI, can be calculated as follows (we assume weak contrast due to diffraction/refraction),
where and are the Fourier transforms of the corresponding images, is the optical transfer function (OTF) of the imaging system (which combines the effect of the finite source size as well as the detector properties) and Wn(u) is the noise power spectrum (NPS) in the images.In the case of the NPWMF, the corresponding SNR, SNRnpw, can be calculated as follows,
Also, it is worth noting that the SNR introduced in our earlier papers [2,4] can, by analogy with Eqs. (1) and (2), be effectively presented as follows,which corresponds to an observer that ignores correlations of noise in the projection images, thus assuming that the noise is white (hence the abbreviation). Equation (3) formally follows from either Eq. (1) or Eq. (2) assuming constant NPS, Wn(u) = Wn(0).3. In-line phase contrast
From the experimental point of view, in-line (also known as propagation-based) phase-contrast imaging is the simplest of the currently used X-ray phase-contrast imaging methods, as it does not involve any optical elements between the source, object and detector. The phase modulation, induced by the object in the transmitted X-ray beam, results in modulation of the intensity in the images due to X-ray diffraction associated with the free-space propagation.
The intensity distribution in an ideal (i.e. not affected by the resolution properties of the imaging system) in-line image located at some distance z from the object, can be well approximated by the following formula, also known as the weak object approximation [10]:
where B and φ are the amplitude attenuation and the phase distributions induced by the features in the object and the asterisk denotes convolution. Here Gam and Gph are the amplitude and phase in-line propagators that, in the one-dimensional (1D) case used in the subsequent analysis, are written as follows:where λ is the X-ray wavelength. We assume the attenuation and the phase shift due to the bulk of the object to be quasi-uniform, so that it results in a uniform attenuation of the incident beam. This contribution to the image formation can be conveniently incorporated into a separate constant factor, , in Eq. (4).For comparison, by neglecting the diffraction effects (this effectively corresponds to setting in Eq. (5), where δD(x) is the Dirac delta function), the intensity distribution Iz,abs in a hypothetical pure absorption image is written as follows,
The Fourier counterpart of Eq. (4) is written as follows [10]:where the in-line amplitude and phase transfer functions, in the 1D case, are written below,We compare two images collected at some finite distance from the object (z > 0): a phase-contrast image Iz (which includes both the absorption and phase contrast) and a hypothetical pure absorption image (which excludes phase contrast). The difference between these two images is completely due to diffraction of the wave exiting the object. Using Eqs. (4) and (6), the difference (diffraction) image, between the phase-contrast and absorption images, is written as follows:
and the corresponding Fourier transform isIn the following we restrict our consideration to a class of monomorphous objects, whose attenuation function is proportional to the phase distribution,
where γ ≡ δ/β is defined by the X-ray complex refractive index of the object, . Equation (11) allows one to rewrite Eq. (10),4. Generalized edge
Consider the following important class of “blurred edge” features whose phase shift function, φ(x;σobj), can be presented as a convolution of the phase function, φid(x), of an “ideal” (sharp) edge and some kernel function, Pobj(x;σobj), whose width is characterised by a parameter σobj:
We assume that the kernel Pobj is normalised in such a way that its integral over x is unity.Consider the following general form for the phase function of an “ideal” edge,
where A and α are some constants. The family of “edge” objects described by Eq. (14) includes the following important cases (see Fig. 1(a)):- 1. a “hyperbolic” edge, for which α = −1/2 and , where ϕ is the linear phase shift of the material of the edge and h is the thickness of the edge at a point x = l;
- 2. a step edge, for which α = 0 and , where t is the height of the step edge;
- 3. an edge of a cylinder, for which α = ½ and , where r is the radius of the cylinder;
- 4. a wedge, for which α = 1 and , where β is the wedge angle.
Figures 1(b) and 1(c) show the examples of the spatial distribution of the phase contrast calculated using Eqs. (13)-(14) and Eq. (12) with (i.e. for the case of pure phase objects), for the selected values of α, in the case of large and small Fresnel numbers, respectively.
Using the following integral, , where Γ(x) is the Gamma function, the Fourier transform of the phase function defined in Eq. (14) can formally be written as follows,
Substituting Eqs. (12) and (15) into Eqs. (1)-(3) one obtains the following expressions for the SNR of the ideal and non-ideal observers (note that these are finite only for −1< α <3/2):
where , and L is the length of the edge feature in the direction orthogonal to both x and z.At this point, it should be emphasized that Eqs. (1)-(3) implicitly assume that the difference image represents a spatial distribution of the photon fluence, i.e. the number of photons per unit area. For this reason the meaning of the parameter in Eqs. (16)-(18) is the photon fluence that would be observed in the image plane in the absence of the feature of interest. Equations (16)-(18) are applicable to in-line images with noise of any origin (e.g. additive electronic noise of the detector etc.). However, in our subsequent analysis we assume that noise in the projection images is associated with the photon counting statistics (Poisson noise) only. In this case Wn, the NPS in the recorded images, is affected by the modulation transfer function (MTF) of the detector, (see, for example, Section 8.2.6 in [8]):
where Wn,wht is the NPS of uncorrelated (white) noise in the image plane. It should be noted that as long as the contrast due to the features is small, this latter NPS, in the case of Poisson noise, is equal to the photon fluence in the image plane, . Equation (19) allows one to rewrite Eqs. (16)-(18) as follows,Equations (20)-(22) indicate a major difference between the ideal observer and non-ideal observers. Indeed, the SNR of the ideal observer depends on the properties of the object and the source and is not affected by the detector properties (the ideal observer performs pre-whitening of the images thus removing any correlations in the noise). In contrast, a non-ideal observer does not perform pre-whitening of the images and therefore its SNR is affected not only by the properties of the object and the source but also by the properties of the detector.
For the sake of simplicity, we assume that all optical transfer functions (more precisely, all modulation transfer functions ) have the same functional form. Moreover we will consider only such MTFs the product of which preserves their original functional form, for example, Gaussians (that corresponds to a Gaussian point response function),
or exponentials (that corresponds to a Lorentz point response function),Under the above assumptions on the MTFs, it is sufficient to analyse separately the following integral,
This integral can be rewritten in an equivalent form,where is the Fresnel number associated with the imaging parameters andis a dimensionless quantity. In the case of large γ (γ ≥ 100), which is typical for materials composed of light elements and for hard X-rays, the function is almost independent of γ and, assuming a Gaussian MTF, Eq. (23), it has simple analytical representations for α = -½, 0, ½ and 1 (i.e. in the cases of “hyperbolic” edge, straight edge, cylinder edge and wedge, respectively):Analysis of Eq. (28) shows that the function has two asymptotes, at large and small Fresnel numbers, respectively. In the former case, up to some constant factors, the function is proportional to , while in the latter case the function asymptotically reaches its maximum value (finite for α ≥ 0 and infinite for α < 0). This is illustrated in Fig. 2 where the dependencies of the integral on the Fresnel number, calculated using Eq. (28), are plotted. Below we investigate these two limiting cases in more detail. In particular, it will be shown that the asymptotic behaviour of the function , which is observed in Fig. 2 and in Eq. (28), is independent of the functional form of the MTFs.
4.1 The near-field imaging regime
This imaging regime occurs when the Fresnel number is much larger than one, NF >> 1 (in fact, the formulas for the SNRs obtained in this section are applicable if NF ≥ 1). In this case the general Eq. (25) can be simplified as follows,
where the positive parameters bi(α) (i = 0,1,2) are defined as follows,In particular, for a Gaussian OTF given by Eq. (23), one can utilise the following integral, , where Γ(x) is the Gamma function, to obtain the following simple analytical expressions for the parameters bi(α) (i = 0,1,2):
For instance,In the same way, for an exponential OTF given by Eq. (24), one can utilise the following integral, , and the parameters bi(α) (i = 0,1,2) are calculated as follows:
Analysis of Eqs. (29), (31) and (32) shows that in the near-field regime the effect of absorption (the terms in Eq. (29) proportional to γ −1 and γ −2) is usually negligibly small and will be neglected in our subsequent analysis of the SNR. The SNR of the ideal observer can then be written as follows,
where incorporates the effect of the finite source and of the feature blurring. In the case of Gaussian distributions for Psrc and Pobj one has: . Here M = R / R1 is the geometrical magnification of the imaging system, R = R1 + R2 is the total source-to-detector distance, R1 and R2 is the source-to-object and object-to-detector distance, respectively. Note that hereafter all images as well as the point response functions are re-projected back to the object plane (this is found to be convenient in the case of cone-beam imaging geometry). This is also consistent with the definition of the effective free space propagation distance given below, in Eq. (36). In the case of Lorentz distributions for Psrc and Pobj, the smearing parameter is defined as follows: .Similarly, the SNR of the non-ideal observer is:
Here, in the case of Gaussian distributions for Psrc, Pobj and Pdet, one has , and , and, in the case of Lorentz distributions for Psrc, Pobj and Pdet, one has: , and .Finally, the SNR of the naïve observer (that treats the noise as white) is:
It is worth mentioning that Eq. (35), in the case of a step edge (α = 0), coincides (up to a constant factor) with the expression for the SNR obtained in [4].Equations (33)-(35) indicate that in the near-field imaging regime the SNR of the ideal observer as well as of the non-ideal observers is proportional to the effective propagation distance z, regardless of the edge type (α),
On the other hand, dependencies of the SNR in Eqs. (33)-(35) on the smearing parameters , and are different for different edge types (α). This behaviour results, in general, in different optimum geometrical configurations (e.g. different optimum magnifications if the total distance R is fixed) for different edges (see Section 5 for more detailed analysis).Besides, in the near-field regime, the Fourier transform of the difference image, Eq. (12) multiplied by the OTF of the imaging system, can be simplified as follows (see also Eq. (15)),
Hence the difference image isEquation (37) indicates that the difference image for any value of σtot can be obtained, by proper rescaling, from one difference image, calculated e.g. for σtot = 1. Also, the extreme values of the difference image (the phase contrast) are proportional to . These results indicate that, in the near-field imaging regime, the spatial resolution in the images of the edge objects is defined by the parameter σtot, in agreement with our earlier findings in [2,4].Importantly, the unnumbered equation preceding Eq. (37) also predicts that in the case of a step edge (α = 0) and a wedge (α = 1) the difference (phase-contrast) image is proportional to the derivative of the total point response function and to the total point response function, respectively, in the direction perpendicular to the edge. These results open a new way for direct measurements of the point response function of the imaging system, by using simple phase objects, a step edge or a wedge.
4.2 The far-field imaging regime
This imaging regime occurs when the Fresnel number is much smaller than one, NF << 1 (strictly speaking, the formulas presented in this section correspond to NF = 0). In this case Eq. (27) can be simplified as follows,
where the positive parameters ci(α) (i = 0,1,2) are calculated as follows, Taking into account Eq. (38), the corresponding integral Fα, Eq. (25), is written as follows:It should be emphasized that as long as the Fresnel number is small, , the function Fα is almost independent of σ, at least for non-negative α.Importantly, in this imaging regime (although with different definitions of NF for the ideal and non-ideal observers) the SNRs for all observers considered here are the same,
It is worth mentioning that Eq. (43), in the case of a step edge (α = 0), coincides (up to a constant factor) with the expression for the SNR obtained in [4].Equation (43) indicates that the only geometrical parameter of the imaging system that affects the SNR in the far-field regime is the effective free space propagation distance, z. Moreover, for the α-values of interest (−1/2 < α < 3/2), the SNR increases with the propagation distance. Thus the maximum SNR in the phase-contrast images of the generalised edge is achieved at the maximum possible propagation distance. For instance, if the total source-to-detector distance R is fixed (this is typical for a laboratory-based imaging system) then the maximum propagation distance, zmax = R/4, is achieved at the optimum magnification Mopt = 2. On the other hand, if the source-to-object distance is fixed (this is typical for a synchrotron-based imaging system), R1 = const, then the effective propagation distance z monotonically increases with the object-to-detector distance R2, and its maximum value, zmax = R1, is achieved at R2 →∞.
It should be mentioned that in the far-field regime, the Fourier transform of the difference image, Eq. (12) multiplied by the OTF of the imaging system, can be simplified as follows, Hence the difference image is
Equation (44) indicates that the difference image, for any value of the width of the first Fresnel zone, (λz)1/2, can be obtained, by proper rescaling of a reference difference image, calculated e.g. for (λz)1/2 = 1. Also, it shows that the dependence of the extreme values of the difference image on the X-ray wavelength and the effective propagation distance is different for different edge types (α). This is a new result which complements our results obtained earlier in [2,4] for a step edge corresponding to α = 0. Equation (44) also suggests that in the far-field imaging regime, the resolution in the images of the edge objects is defined by the width of the first Fresnel zone, (λz)1/2, which is in agreement with our earlier findings [2,4].5. Discussion and conclusions
The formalism proposed in Section 4, in particular the general Eqs. (20)-(22) and their special limiting equivalents for the cases of large and small Fresnel numbers, Eqs. (33)-(35) and (43), represent the main results of this paper that allow one to evaluate quantitatively the additional information content that is introduced in the X-ray projection images by in-line phase contrast, for a broad class of edge objects. The SNRs utilised throughout this paper present a measure of distinguishability between the phase-contrast images and their pure absorption equivalents. It is well known that in-line phase contrast results in enhancing the boundaries between different features in the object thus improving the image quality by means of a “natural” (as opposed to image post-processing) edge enhancement. For this reason one would usually be interested in optimising the imaging conditions so as to maximise the SNR. One example of such optimisation is the optimisation of the X-ray energy (or the energy band) for a fixed radiation dose delivered to the object. Another important aspect is the optimisation of the imaging geometry, for example, in terms of the geometrical magnification of the imaging system provided that the total source-to-detector distance is fixed. An example of this latter optimisation is shown in Fig. 3 and in Table 1. Here the SNR values for three different observers are plotted together versus the magnification. The results correspond to certain values of the parameters of the source, σsrc, the edge blurring, σobj, the detector resolution, σdet, and the total source-to-detector distance, R. Gaussian distributions are assumed for all point response functions.
Analysis of Fig. 3 indicates that the optimum magnification that maximises the SNR depends not only on the above mentioned imaging parameters but also on the edge type (described by α in Eq. (14)) and on the observer. For instance, the optimum magnification for the ideal observer is always less than 2. For comparison, the optimum magnification that maximises SNRwht is bigger than 2 if σdet > σsrc, less than 2 if σdet < σsrc and equal to 2 if σdet = σsrc. This is true for both Gaussian and Lorentz point response functions. The dependence of the optimum magnification on the imaging system parameters for the NPWMF is even more complicated (for instance, it is not in general equal to 2 in the case σdet = σsrc). Table 1 contains optimum magnifications and the corresponding maximum values of three observers and three different edges, for certain imaging conditions. One can see that, as expected, the ideal observer has the maximum SNR, as a result of noise pre-whitening. Note however that a real human observer is unable to perform this noise pre-whitening and their performance is better matched by the NPWMF [8]. It is worth noting that the naïve observer (see SNRwht in Fig. 3 and Table 1) which corresponds to our earlier analysis in [2,4] turns out to be much closer to the NPWMF (SNRnpw in Fig. 3 and Table 1) than to the ideal observer.
Figure 3 and Table 1 also show that the optimum magnification and the corresponding SNR differ significantly with the edge type. One clear tendency is that, under the same imaging conditions, by increasing the parameter α the optimum magnification tends to be closer to 2, regardless of the observer. Another important observation is that the maximum SNR decreases with α (at least for the geometrical parameters of different edges used in our calculations). This can be qualitatively explained by the fact that the in-line phase contrast is very sensitive to the edge ‘sharpness’ or curvature. This ‘sharpness’ generally decreases with α, at least in the range from 0 to 1 (see Fig. 1).
It should be emphasized that the FoMs presented by Eqs. (1)-(3) in Section 2 are generic and can be applied to quantification of phase contrast in other imaging modalities, at least numerically. This would allow a direct quantitative comparison of different imaging modalities. This comparison is however beyond the scope of this paper.
References and links
1. A. Pogany, D. Gao, and S. W. Wilkins, “Contrast and resolution in imaging with a microfocus x-ray source,” Rev. Sci. Instrum. 68(7), 2774–2782 (1997). [CrossRef]
2. Ya. I. Nesterets, S. W. Wilkins, T. E. Gureyev, A. Pogany, and A. W. Stevenson, “On the optimization of experimental parameters for x-ray in-line phase-contrast imaging,” Rev. Sci. Instrum. 76(9), 093706 (2005). [CrossRef]
3. E. Pagot, S. Fiedler, P. Cloetens, A. Bravin, P. Coan, K. Fezzaa, J. Baruchel, and J. Härtwig, “Quantitative comparison between two phase contrast techniques: diffraction enhanced imaging and phase propagation imaging,” Phys. Med. Biol. 50(4), 709–724 (2005). [CrossRef] [PubMed]
4. T. E. Gureyev, Y. I. Nesterets, A. W. Stevenson, P. R. Miller, A. Pogany, and S. W. Wilkins, “Some simple rules for contrast, signal-to-noise and resolution in in-line x-ray phase-contrast imaging,” Opt. Express 16(5), 3223–3241 (2008). [CrossRef] [PubMed]
5. C. Y. Chou and M. A. Anastasio, “Influence of imaging geometry on noise texture in quantitative in-line X-ray phase-contrast imaging,” Opt. Express 17(17), 14466–14480 (2009). [CrossRef] [PubMed]
6. M. A. Anastasio, C. Y. Chou, A. M. Zysk, and J. G. Brankov, “Analysis of ideal observer signal detectability in phase-contrast imaging employing linear shift-invariant optical systems,” J. Opt. Soc. Am. A 27(12), 2648–2659 (2010). [CrossRef] [PubMed]
7. P. C. Diemoz, A. Bravin, M. Langer, and P. Coan, “Analytical and experimental determination of signal-to-noise ratio and figure of merit in three phase-contrast imaging techniques,” Opt. Express 20(25), 27670–27690 (2012). [CrossRef] [PubMed]
8. H. H. Barrett and K. J. Myers, Foundations of Image Science (John Wiley & Sons, Inc., 2004).
9. A. Bravin, P. Coan, and P. Suortti, “X-ray phase-contrast imaging: from pre-clinical applications towards clinics,” Phys. Med. Biol. 58(1), R1–R35 (2013). [CrossRef] [PubMed]
10. J. P. Guigay, “Fourier transform analysis of Fresnel diffraction patterns and in-line holograms,” Optik (Stuttg.) 49, 121–125 (1977).