Quantitative evaluation of in-line x-ray phase contrast for a class of edge objects

Yakov I. Nesterets

doi:10.1364/OE.22.005937

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 I_abs, the SNR of the ideal observer, SNR_I, can be calculated as follows (we assume weak contrast due to diffraction/refraction),

S N R_{I}^{2} = \int d u \frac{| \hat{I} (u) - {\hat{I}}_{a b s} (u) |^{2} | {\hat{P}}_{s y s} (u) |^{2}}{W_{n} (u)},

where

\hat{I} (u)

and

{\hat{I}}_{a b s} (u)

are the Fourier transforms of the corresponding images,

{\hat{P}}_{s y s} (u) = {\hat{P}}_{s r c} (u) {\hat{P}}_{d e t} (u)

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 W_n(u) is the noise power spectrum (NPS) in the images.

In the case of the NPWMF, the corresponding SNR, SNR_npw, can be calculated as follows,

S N R_{n p w}^{2} = \frac{{[\int d u | \hat{I} (u) - {\hat{I}}_{a b s} (u) |^{2} | {\hat{P}}_{s y s} (u) |^{2}]}^{2}}{\int d u | \hat{I} (u) - {\hat{I}}_{a b s} (u) |^{2} | {\hat{P}}_{s y s} (u) |^{2} W_{n} (u)} .

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,

S N R_{w h t}^{2} = W_{n}^{- 1} (0) \int d u | \hat{I} (u) - {\hat{I}}_{a b s} (u) |^{2} | {\hat{P}}_{s y s} (u) |^{2},

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, W_n(u) = W_n(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]:

I_{z} ≅ {\bar{I}}_{0} {1 - 2 B * G_{a m} - 2 φ * G_{p h}},

where B and φ are the amplitude attenuation and the phase distributions induced by the features in the object and the asterisk denotes convolution. Here G_am and G_ph are the amplitude and phase in-line propagators that, in the one-dimensional (1D) case used in the subsequent analysis, are written as follows:

G_{a m} (x) = Re [G (x)], G_{p h} (x) = Im [G (x)], G (x) = {(λ z)}^{- 1 / 2} \exp [i π (\frac{x^{2}}{λ z} - \frac{1}{4})],

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,

{\bar{I}}_{0}

, in Eq. (4).

For comparison, by neglecting the diffraction effects (this effectively corresponds to setting $G (x) = δ_{D} (x)$ in Eq. (5), where δ_D(x) is the Dirac delta function), the intensity distribution I_z,abs in a hypothetical pure absorption image is written as follows,

I_{z, a b s} ≅ {\bar{I}}_{0} {1 - 2 B} .

The Fourier counterpart of Eq. (4) is written as follows [10]:

{\hat{I}}_{z} ≅ {\bar{I}}_{0} {δ_{D} - 2 \hat{B} {\hat{G}}_{a m} - 2 \hat{φ} {\hat{G}}_{p h}},

where the in-line amplitude and phase transfer functions, in the 1D case, are written below,

{\hat{G}}_{a m} (u) = Re [\hat{G} (u)], {\hat{G}}_{p h} (u) = Im [\hat{G} (u)], \hat{G} (u) = \exp [- i π λ z u^{2}] .

We compare two images collected at some finite distance from the object (z > 0): a phase-contrast image I_z (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:

I_{d i f f} \equiv I_{z} - I_{z, a b s} = {\bar{I}}_{0} {2 B - 2 B * G_{a m} - 2 φ * G_{p h}},

and the corresponding Fourier transform is

{\hat{I}}_{d i f f} = {\bar{I}}_{0} {2 \hat{B} [1 - {\hat{G}}_{a m}] - 2 \hat{φ} {\hat{G}}_{p h}} = {\bar{I}}_{0} {2 \hat{B} [1 - \cos (π λ z u^{2})] + 2 \hat{φ} \sin (π λ z u^{2})} .

In the following we restrict our consideration to a class of monomorphous objects, whose attenuation function is proportional to the phase distribution,

B = - γ^{- 1} φ,

where γ ≡ δ/β is defined by the X-ray complex refractive index of the object,

n = 1 - δ + i β

. Equation (11) allows one to rewrite Eq. (10),

{\hat{I}}_{d i f f} = {\bar{I}}_{0} 2 \hat{φ} {- 2 γ^{- 1} \sin^{2} (π λ z u^{2} / 2) + \sin (π λ z u^{2})} .

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, P_obj(x;σ_obj), whose width is characterised by a parameter σ_obj:

φ (x; σ_{o b j}) = (φ_{i d} * P_{o b j}) (x; σ_{o b j}) .

We assume that the kernel P_obj 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,

φ_{i d} (x) = {\begin{cases} A x^{α}, x > 0, \\ 0, x \leq 0, \end{cases}

where A and α are some constants. The family of “edge” objects described by Eq. (14) includes the following important cases (see Fig. 1(a)):

Fig. 1 (a) Schematic representation of several “edge” objects from the family described by Eqs. (13) and (14) and the corresponding phase contrast calculated using Eq. (12) with γ ⁻¹ = 0, for large Fresnel number (b) and small Fresnel number (c), $N_{F} = σ_{o b j}^{2} / (λ z)$ .

Download Full Size | PDF

1. a “hyperbolic” edge, for which α = −1/2 and $A_{- 1 / 2} = ϕ h / \sqrt{l}$ , 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 $A_{0} = ϕ t$ , where t is the height of the step edge;
3. an edge of a cylinder, for which α = ½ and $A_{1 / 2} = ϕ \sqrt{2 r}$ , where r is the radius of the cylinder;
4. a wedge, for which α = 1 and $A_{1} = ϕ \tan β$ , 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 $γ^{- 1} = 0$ (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, $\int_{0}^{+ \infty} d x x^{q - 1} e^{\pm i m x} = Γ (q) m^{- q} e^{\pm i q π / 2}$ , where Γ(x) is the Gamma function, the Fourier transform of the phase function defined in Eq. (14) can formally be written as follows,

{\hat{φ}}_{i d}^{} (u) = A Γ (1 + α) {(- 2 π i u)}^{- 1 - α} .

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):

S N R_{I}^{2} = L {\bar{I}}_{0}^{2} A^{2} C_{α} \int_{- \infty}^{+ \infty} d u \frac{| u^{- 1 - α} {- 2 γ^{- 1} \sin^{2} (π λ z u^{2} / 2) + \sin (π λ z u^{2})} |^{2} | {\hat{P}}_{t o t} (u) |^{2}}{W_{n} (u)},

S N R_{n p w}^{2} = L {\bar{I}}_{0}^{2} A^{2} C_{α} \frac{{[\int_{- \infty}^{+ \infty} d u | u^{- 1 - α} {- 2 γ^{- 1} \sin^{2} (π λ z u^{2} / 2) + \sin (π λ z u^{2})} |^{2} | {\hat{P}}_{t o t} (u) |^{2}]}^{2}}{\int_{- \infty}^{+ \infty} d u | u^{- 1 - α} {- 2 γ^{- 1} \sin^{2} (π λ z u^{2} / 2) + \sin (π λ z u^{2})} |^{2} | {\hat{P}}_{t o t} (u) |^{2} W_{n} (u)},

S N R_{w h t}^{2} = \frac{L {\bar{I}}_{0}^{2} A^{2} C_{α}}{W_{n} (0)} \int_{- \infty}^{+ \infty} d u | u^{- 1 - α} {- 2 γ^{- 1} \sin^{2} (π λ z u^{2} / 2) + \sin (π λ z u^{2})} |^{2} | {\hat{P}}_{t o t} (u) |^{2},

where

C_{α} \equiv 4 Γ^{2} (1 + α) / {(2 π)}^{2 (1 + α)}

,

| {\hat{P}}_{t o t} (u) |^{2} = | {\hat{P}}_{s y s} (u) |^{2} | {\hat{P}}_{o b j} (u) |^{2}

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 $I_{d i f f} \equiv I - I_{a b s}$ 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 ${\bar{I}}_{0}$ 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 W_n, the NPS in the recorded images, is affected by the modulation transfer function (MTF) of the detector, $| {\hat{P}}_{d e t} (u) |$ (see, for example, Section 8.2.6 in [8]):

W_{n} (u) = W_{n, w h t} | {\hat{P}}_{d e t} (u) |^{2},

where W_n,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,

W_{n, w h t} = {\bar{I}}_{0}

. Equation (19) allows one to rewrite Eqs. (16)-(18) as follows,

S N R_{I}^{2} = L {\bar{I}}_{0}^{} A^{2} C_{α} \int_{- \infty}^{+ \infty} d u | u^{- 1 - α} {- 2 γ^{- 1} \sin^{2} (π λ z u^{2} / 2) + \sin (π λ z u^{2})} |^{2} | {\hat{P}}_{s r c} (u) {\hat{P}}_{o b j} (u) |^{2},

S N R_{n p w}^{2} = L {\bar{I}}_{0}^{} A^{2} C_{α} \frac{{[\int_{- \infty}^{+ \infty} d u | u^{- 1 - α} {- 2 γ^{- 1} \sin^{2} (π λ z u^{2} / 2) + \sin (π λ z u^{2})} |^{2} | {\hat{P}}_{t o t} (u) |^{2}]}^{2}}{\int_{- \infty}^{+ \infty} d u | u^{- 1 - α} {- 2 γ^{- 1} \sin^{2} (π λ z u^{2} / 2) + \sin (π λ z u^{2})} |^{2} | {\hat{P}}_{t o t} (u) {\hat{P}}_{d e t} (u) |^{2}},

S N R_{w h t}^{2} = L {\bar{I}}_{0}^{} A^{2} C_{α} \int_{- \infty}^{+ \infty} d u | u^{- 1 - α} {- 2 γ^{- 1} \sin^{2} (π λ z u^{2} / 2) + \sin (π λ z u^{2})} |^{2} | {\hat{P}}_{t o t} (u) |^{2} .

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 ${\hat{P}}_{i} (i = s r c, o b j, d e t)$ (more precisely, all modulation transfer functions $| {\hat{P}}_{i} |$ ) 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),

P_{1} (x; σ) = {(\sqrt{2 π} σ)}^{- 1} \exp [- 1 / 2 {(x / σ)}^{2}], {\hat{P}}_{1} (u; σ) = {\hat{P}}_{1} (σ u) = \exp [- 2 π^{2} {(σ u)}^{2}],

or exponentials (that corresponds to a Lorentz point response function),

P_{2} (x; σ) = \frac{{(π σ)}^{- 1}}{1 + {(x / σ)}^{2}}, {\hat{P}}_{2} (u; σ) = {\hat{P}}_{2} (σ u) = \exp (- 2 π | σ u |) .

Under the above assumptions on the MTFs, it is sufficient to analyse separately the following integral,

F_{α} (λ z, γ, σ) \equiv \int_{- \infty}^{+ \infty} d u | u^{- 1 - α} {- 2 γ^{- 1} \sin^{2} (π λ z u^{2} / 2) + \sin (π λ z u^{2})} |^{2} | \hat{P} (σ u) |^{2} .

This integral can be rewritten in an equivalent form,

F_{α} (λ z, γ, σ) = {(λ z)}^{1 / 2 + α} {F^{'}}_{α} (γ, N_{F}),

where

N_{F} \equiv σ^{2} / (λ z)

is the Fresnel number associated with the imaging parameters and

{F^{'}}_{α} (γ, N_{F}) \equiv \int_{- \infty}^{+ \infty} d v | v^{- 1 - α} {- 2 γ^{- 1} \sin^{2} (π v^{2} / 2) + \sin (π v^{2})} |^{2} | \hat{P} (N_{F}^{1 / 2} v) |^{2}

is 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

{F^{'}}_{α} (γ, N_{F}) ≅ {F^{'}}_{α} (N_{F})

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):

{F^{'}}_{- 1 / 2} (N_{F}) ≅ 4^{- 1} \ln [1 + {(2 π N_{F})}^{- 2}],

{F^{'}}_{0} (N_{F}) ≅ π {(2 π N_{F}^{})}^{1 / 2} [\sqrt{1 + \sqrt{1 + {(2 π N_{F})}^{- 2}}} - \sqrt{2}],

{F^{'}}_{1 / 2} (N_{F}) ≅ π arc \tan [{(2 π N_{F})}^{- 1}] - (π / 2) (2 π N_{F}) \ln [1 + {(2 π N_{F})}^{- 2}],

{F^{'}}_{1} (N_{F}) ≅ 3^{- 1} {(2 π)}^{2} {(2 π N_{F}^{})}^{3 / 2} [\sqrt{1 + {(2 π N_{F})}^{- 2}} \sqrt{1 + \sqrt{1 + {(2 π N_{F})}^{- 2}}} + \sqrt{2}] .

Analysis of Eq. (28) shows that the function ${F^{'}}_{α}$ has two asymptotes, at large and small Fresnel numbers, respectively. In the former case, up to some constant factors, the function ${F^{'}}_{α}$ is proportional to ${(2 π N_{F})}^{α - 3 / 2}$ , while in the latter case the function ${F^{'}}_{α}$ asymptotically reaches its maximum value (finite for α ≥ 0 and infinite for α < 0). This is illustrated in Fig. 2 where the dependencies of the integral ${F^{'}}_{α}$ 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 ${F^{'}}_{α}$ , which is observed in Fig. 2 and in Eq. (28), is independent of the functional form of the MTFs.

Fig. 2 Dependencies of the integral ${F^{'}}_{α}$ on the Fresnel number, given by Eq. (28).

Download Full Size | PDF

4.1 The near-field imaging regime

This imaging regime occurs when the Fresnel number is much larger than one, N_F >> 1 (in fact, the formulas for the SNRs obtained in this section are applicable if N_F ≥ 1). In this case the general Eq. (25) can be simplified as follows,

\begin{array}{l} {F_{α} (λ z, γ, σ) |}_{N_{F} > > 1} ≅ \int_{- \infty}^{+ \infty} d u | u^{- 1 - α} {- 2 γ^{- 1} {(π λ z u^{2} / 2)}^{2} + π λ z u^{2}} |^{2} | \hat{P} (σ u) |^{2} \\ = b_{0} (α) \frac{{(λ z)}^{2}}{σ^{3 - 2 α}} - b_{1} (α) γ^{- 1} \frac{{(λ z)}^{3}}{σ^{5 - 2 α}} + b_{2} (α) γ^{- 2} \frac{{(λ z)}^{4}}{σ^{7 - 2 α}} \\ = σ^{1 + 2 α} {b_{0} (α) N_{F}^{- 2} - b_{1} (α) γ^{- 1} N_{F}^{- 3} + b_{2} (α) γ^{- 2} N_{F}^{- 4}}, \end{array}

where the positive parameters b_i(α) (i = 0,1,2) are defined as follows,

\begin{array}{l} b_{0} (α) \equiv π^{2} \int_{- \infty}^{+ \infty} d v | \hat{P} (v) |^{2} {(v^{2})}^{1 - α}, \\ b_{1} (α) \equiv π^{3} \int_{- \infty}^{+ \infty} d v | \hat{P} (v) |^{2} {(v^{2})}^{2 - α}, \\ b_{2} (α) \equiv \frac{π^{4}}{4} \int_{- \infty}^{+ \infty} d v | \hat{P} (v) |^{2} {(v^{2})}^{3 - α} . \end{array}

In particular, for a Gaussian OTF given by Eq. (23), one can utilise the following integral, $\int_{- \infty}^{+ \infty} d u | u |^{β} \exp (- A^{2} u^{2}) = \frac{Γ [(β + 1) / 2]}{A^{β + 1}}, β > - 1$ , where Γ(x) is the Gamma function, to obtain the following simple analytical expressions for the parameters b_i(α) (i = 0,1,2):

b_{0} (α) = π^{2} \frac{Γ (3 / 2 - α)}{{(2 π)}^{3 - 2 α}}, b_{1} (α) = π^{3} \frac{Γ (5 / 2 - α)}{{(2 π)}^{5 - 2 α}}, b_{2} (α) = \frac{π^{4}}{4} \frac{Γ (7 / 2 - α)}{{(2 π)}^{7 - 2 α}} .

For instance,

b_{0} (0) \approx 0.035, b_{1} (0) \approx 0.0042, b_{2} (0) \approx 2.1 \times 10^{- 4} .

In the same way, for an exponential OTF given by Eq. (24), one can utilise the following integral, $\int_{- \infty}^{+ \infty} d u | u |^{β} \exp (- A | u |) = \frac{2 Γ (β + 1)}{A^{β + 1}}, β > - 1$ , and the parameters b_i(α) (i = 0,1,2) are calculated as follows:

b_{0} (α) = 2 π^{2} \frac{Γ (3 - 2 α)}{{(4 π)}^{3 - 2 α}}, b_{1} (α) = 2 π^{3} \frac{Γ (5 - 2 α)}{{(4 π)}^{5 - 2 α}}, b_{2} (α) = \frac{π^{4}}{2} \frac{Γ (7 - 2 α)}{{(4 π)}^{7 - 2 α}} .

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 γ ⁻¹ and γ ⁻²) 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,

{S N R_{I}^{2} |}_{N_{F} > > 1} \approx L {\bar{I}}_{0}^{} A^{2} C_{α} b_{0} (α) \frac{{(λ z)}^{2}}{σ_{s r c, o b j}^{3 - 2 α}},

where

σ_{s r c, o b j}

incorporates the effect of the finite source and of the feature blurring. In the case of Gaussian distributions for P_src and P_obj one has:

σ_{s r c, o b j}^{2} = M^{- 2} {(M - 1)}^{2} σ_{s r c}^{2} + σ_{o b j}^{2}

. Here M = R / R₁ is the geometrical magnification of the imaging system, R = R₁ + R₂ is the total source-to-detector distance, R₁ and R₂ 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 P_src and P_obj, the smearing parameter

σ_{s r c, o b j}^{}

is defined as follows:

σ_{s r c, o b j}^{} = M^{- 1} (M - 1) σ_{s r c}^{} + σ_{o b j}^{}

.

Similarly, the SNR of the non-ideal observer is:

{S N R_{n p w}^{2} |}_{N_{F} > > 1} \approx L {\bar{I}}_{0}^{} A^{2} C_{α} b_{0} (α) {(λ z)}^{2} \frac{σ_{t o t, d e t}^{3 - 2 α}}{{[σ_{t o t}^{3 - 2 α}]}^{2}} .

Here, in the case of Gaussian distributions for P_src, P_obj and P_det, one has

σ_{t o t}^{2} = M^{- 2} {(M - 1)}^{2} σ_{s r c}^{2} + σ_{o b j}^{2} + M^{- 2} σ_{d e t}^{2}

, and

σ_{t o t, d e t}^{2} = σ_{t o t}^{2} + M^{- 2} σ_{d e t}^{2}

, and, in the case of Lorentz distributions for P_src, P_obj and P_det, one has:

σ_{t o t}^{} = M^{- 1} (M - 1) σ_{s r c}^{} + σ_{o b j}^{} + M^{- 1} σ_{d e t}^{}

, and

σ_{t o t, d e t}^{} = σ_{t o t}^{} + M^{- 1} σ_{d e t}^{}

.

Finally, the SNR of the naïve observer (that treats the noise as white) is:

{S N R_{w h t}^{2} |}_{N_{F} > > 1} \approx L {\bar{I}}_{0}^{} A^{2} C_{α} b_{0} (α) \frac{{(λ z)}^{2}}{σ_{t o t}^{3 - 2 α}} .

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 (α),

z \equiv R_{1} R_{2} / (R_{1} + R_{2}) = R M^{- 2} (M - 1) .

On the other hand, dependencies of the SNR in Eqs. (33)-(35) on the smearing parameters

σ_{s r c, o b j}^{}

,

σ_{t o t}^{}

and

σ_{t o t, d e t}^{}

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)),

{{\hat{I}}_{d i f f} (u; σ_{t o t}) |}_{N_{F} > > 1} ≅ 2 π {\bar{I}}_{0} λ z u^{2} {\hat{φ}}_{i d} (u) {\hat{P}}_{t o t} (σ u) = σ_{t o t}^{α - 1} {{\hat{I}}_{d i f f} (σ_{t o t} u; 1) |}_{N_{F} > > 1} .

Hence the difference image is

{I_{d i f f} (x; σ_{t o t}) |}_{N_{F} > > 1} = σ_{t o t}^{α - 2} {I_{d i f f} (x / σ_{t o t}; 1) |}_{N_{F} > > 1} .

Equation (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

σ_{t o t}^{α - 2}

. 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, N_F << 1 (strictly speaking, the formulas presented in this section correspond to N_F = 0). In this case Eq. (27) can be simplified as follows,

\begin{array}{l} {F^{'}}_{α} (γ, N_{F} < < 1) ≅ \int_{- \infty}^{+ \infty} d v | v^{- 1 - α} {- 2 γ^{- 1} \sin^{2} (π v^{2} / 2) + \sin (π v^{2})} |^{2} \\ = c_{0} (α) - c_{1} (α) γ^{- 1} + c_{2} (α) γ^{- 2}, \end{array}

where the positive parameters c_i(α) (i = 0,1,2) are calculated as follows,

\begin{array}{l} c_{0} (α) \equiv \int_{- \infty}^{+ \infty} d u \frac{\sin^{2} (π u^{2})}{{(u^{2})}^{1 + α}} = \frac{π^{1 + α} Γ (3 / 4 - α / 2)}{(1 + 2 α) Γ (3 / 4 + α / 2)}, - 1 / 2 \leq α < 3 / 2, \\ c_{0} (- 1 / 2) = \infty, c_{0} (0) = π, c_{0} (1 / 2) = π^{2} / 2, c_{0} (1) = 4 π^{2} / 3, \end{array}

\begin{array}{l} c_{1} (α) \equiv 4 \int_{- \infty}^{+ \infty} d u \frac{\sin (π u^{2}) \sin^{2} (π u^{2} / 2)}{{(u^{2})}^{1 + α}}, \\ c_{1} (- 1 / 2) = π / 2, c_{1} (0) = 2 π (2^{1 / 2} - 1), c_{1} (1 / 2) = 2 π \ln 2, c_{1} (1) = \frac{8 π^{2}}{3} (1 - 2^{- 1 / 2}), \end{array}

\begin{array}{l} c_{2} (α) \equiv 4 \int_{- \infty}^{+ \infty} d u \frac{\sin^{4} (π u^{2} / 2)}{{(u^{2})}^{1 + α}}, \\ c_{2} (- 1 / 2) = \infty, c_{2} (0) = π (2^{3 / 2} - 1), c_{2} (1 / 2) = π^{2} / 2, c_{2} (1) = \frac{4 π^{2}}{3} (2^{1 / 2} - 1) . \end{array}

Taking into account Eq. (38), the corresponding integral F_α, Eq. (25), is written as follows:

{F_{α} (λ z, γ, σ) |}_{N_{F} < < 1} = {(λ z)}^{1 / 2 + α} [c_{0} (α) - c_{1} (α) γ^{- 1} + c_{2} (α) γ^{- 2}] .

It should be emphasized that as long as the Fresnel number is small,

σ^{2} / (λ z) < < 1

, the function F_α is almost independent of σ, at least for non-negative α.

Importantly, in this imaging regime (although with different definitions of N_F for the ideal and non-ideal observers) the SNRs for all observers considered here are the same,

{S N R_{i}^{2} |}_{N_{F} < < 1} ≅ L {\bar{I}}_{0}^{} A^{2} C_{α} {(λ z)}^{1 / 2 + α} [c_{0} (α) - c_{1} (α) γ^{- 1} + c_{2} (α) γ^{- 2}], i = I, n p w, w h t .

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, z_max = R/4, is achieved at the optimum magnification M_opt = 2. On the other hand, if the source-to-object distance is fixed (this is typical for a synchrotron-based imaging system), R₁ = const, then the effective propagation distance z monotonically increases with the object-to-detector distance R₂, and its maximum value, z_max = R₁, is achieved at R₂ →∞.

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, ${{\hat{I}}_{d i f f} (u; \sqrt{λ z}) |}_{N_{F} < < 1} ≅ 2 {\bar{I}}_{0} {\hat{φ}}_{i d} (u) {- 2 γ^{- 1} \sin^{2} (π λ z u^{2} / 2) + \sin (π λ z u^{2})} = {(λ z)}^{\frac{1 + α}{2}} {{\hat{I}}_{d i f f} (\sqrt{λ z} u; 1) |}_{N_{F} < < 1} .$ Hence the difference image is

{I_{d i f f} (x; \sqrt{λ z}) |}_{N_{F} < < 1} = {(λ z)}^{α / 2} {I_{d i f f} (x / \sqrt{λ z}; 1) |}_{N_{F} < < 1} .

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.

Fig. 3 Dependencies of various SNRs (a-c) and of the corresponding Fresnel numbers (d) on the magnification of the imaging system, for different types of the edge (α = 0, ½ and 1). Gaussian distributions are assumed for the source, the edge smearing (σ_obj = 1μm) and the detector PSF (σ_det = 10μm). Solid and dashed lines correspond to two different source sizes, σ_src = 5μm and 50μm, respectively.

Download Full Size | PDF

Table 1. Edge length L = 1mm, image fluence ${\bar{I}}_{0}$ = 10 ph/μm², wavelength λ = 0.05nm, total distance R = 2m, linear phase shift ϕ = −5.5 rad/mm and γ = −1061 (Nylon in water), sheet thickness t = 100μm, cylinder radius r = 50μm, wedge slope tanβ = 1. Gaussian PSFs.

View Table

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 SNR_wht 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 SNR_wht in Fig. 3 and Table 1) which corresponds to our earlier analysis in [2,4] turns out to be much closer to the NPWMF (SNR_npw 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).

σ_src, μm	σ_obj, μm	σ_det, μm	SNR	α = 0		α = 1/2		α = 1
σ_src, μm	σ_obj, μm	σ_det, μm	SNR	M_opt	SNR _max	M_opt	SNR _max	M_opt	SNR _max
5	1	10	SNR_wht	2.747	6.916	2.562	1.455	2.295	0.497
			SNR_I	1.276	26.459	1.408	3.369	1.639	0.735
			SNR_npw	2.388	9.788	2.301	1.860	2.175	0.567
50	1	10	SNR_wht	1.245	1.143	1.341	0.397	1.584	0.245
			SNR_I	1.028	3.953	1.075	0.553	1.501	0.253
			SNR_npw	1.184	1.524	1.260	0.454	1.522	0.252

Quantitative evaluation of in-line x-ray phase contrast for a class of edge objects

Abstract

1. Introduction

2. Figure-of-merit

3. In-line phase contrast

4. Generalized edge

4.1 The near-field imaging regime

4.2 The far-field imaging regime

5. Discussion and conclusions

References and links

Cited By

Figures (3)

Tables (1)

Equations (48)

Optics Express