Expand this Topic clickable element to expand a topic
Skip to content
Optica Publishing Group

A quantitative, non-interferometric X-ray phase contrast imaging technique

Open Access Open Access

Abstract

We present a quantitative, non-interferometric, X-ray differential phase contrast imaging technique based on the edge illumination principle. We derive a novel phase retrieval algorithm which requires only two images to be acquired and verify the technique experimentally using synchrotron radiation. The technique is useful for planar imaging but is expected to be important for quantitative phase tomography also. The properties and limitations of the technique are studied in detail.

© 2012 Optical Society of America

1. Introduction

Conventional radiographs record the projection, along the direction of X-ray propagation, of an object’s absorption properties. Image contrast thus arises due to spatial variation of an object’s thickness and constituent material absorption coefficients. Since the pioneering work of Bonse and Hart [1], a variety of techniques have been developed to acquire X-ray images possessing contrast due to sample induced phase shifts [27]. Such images are referred to as phase contrast images. Since then, techniques have been developed which are capable of quantifying the phase shift. Such techniques, whilst being useful for planar imaging, are also essential for performing quantitative phase tomography.

Existing quantitative X-ray phase contrast imaging (XPCI) techniques may be divided into three categories. The first category may be denoted analyser based imaging (ABI) and develops contrast using the rocking curve of an analyser crystal [2]. Using this technique, contributions to image contrast due to the phase and absorption are able to separated. Furthermore, researchers have developed algorithms for quantifying the phase shift using various simplifying approximations [811]. This method requires a beam of high spectral purity and small angular divergence thus limiting it to synchrotron sources or monochromated laboratory sources [2, 12].

The second method for performing quantitative XPCI is based upon a technique known as in-line holography [4, 13]. This technique requires a source of high spatial coherence and so must be performed using synchrotron radiation or a microfocal source. Algorithms for phase extraction have largely employed transport of intensity formalism and generally require one or more assumptions to be made. For example, Nugent et al. [14] introduced a phase retrieval algorithm for the special case of uniform irradiance. Algorithms were introduced for performing quantitative phase imaging of phase objects [15] and weakly absorbing objects [16]. Paganin et al. [17] introduced a formalism for phase extraction based upon the assumption of a homogeneous sample. Gureyev later combined the transport of intensity formalism with the first Born formalism to obtain new phase retrieval algorithms [18, 19]. These techniques have all been extended and we refer the reader to a recent review [20] for a thorough account.

The final method, known as grating interferometry, employs two or three gratings and the principle of Talbot self-imaging. A sample’s phase and absorption information is thus encoded in the relative shift and amplitude, respectively, of fringes formed through self-imaging. Detectors with pixel sizes small enough to resolve the fringes are, however, unsatisfactory for imaging large fields of view with clinically compatible exposure times. Phase stepping [21] and the moiré [22] configuration, developed for use in optical interferometry, have both been employed to overcome this problem. The generation of fringes requires a spatially coherent source or an array or such sources in order to achieve Talbot self-imaging. This method has thus been employed using synchrotron radiation [23, 24] and apertured laboratory sources [2527].

The principle employed in this paper is called edge illumination XPCI (EIXPCI) [28] and is depicted in Fig. 1(a). EIXPCI works by projecting an X-ray beam onto the edge of a sensitive region of a detector. When the arrangement in Fig. 1(a) is tiled as shown in Fig. 1(b), the coded aperture XPCI (CAXPCI) system results. Image formation is equivalent in both systems with the exception that the sample must be scanned through the beam in the EIXPCI case. The analysis presented in this paper thus applies equally to both cases. For the remainder of this paper we refer principally to EIXPCI assuming that the results also apply to CAXPCI.

 figure: Fig. 1

Fig. 1 (a) Demonstration of the edge illumination principle whereby pre-sample apertures, located a distance zso from the source, shape a beam which is incident upon detector apertures and detector, placed a distance zod from the pre-sample apertures. The centre of the detector apertures is offset from the centre of the beam by a distance ΔP. (b) A CAXPCI system made by tiling the setup used in (a). A1 and A2 are apertures with the same projected pitch which also matches the detector pixel width P. All parameters have the same meaning as in (a). A1 creates X-ray beams which are incident upon A2 which is used to achieve the partial pixel illumination condition.The imaging system continues uniformly out of the page.

Download Full Size | PDF

Superficially, EIXPCI shares similarities with systems from two categories: ABI and grating interferometry. On first inspection, EIXPCI seems to be in the family of grating based interferometers due to the presence of two sets of apertures. The technique is, however, non-interferometric since it does not employ Talbot’s self-imaging phenomenon [29]. In this regard the edge illumination technique is related to that reported recently by Huang et. al [30] who, contrary to EIXPCI, employ phase stepping and do not make use of edge illumination. On the other hand, EIXPCI was first developed in analogy with ABI [28] due to the manner in which the pixel edge creates sensitivity to the angle of photon refraction. It is no surprise then, that the quantitative adaptation of EIXPCI shares more in common with ABI than grating interferometry since an image of the object is taken for two orientations of the detector apertures relative to the object. Furthermore, phase stepping is not employed.

In this paper we begin by deriving an algorithm for extracting the gradient of the phase shift due to a sample imaged using the EIXPCI system. We then derive the limitations and important properties of the technique before showing experimental results obtained using synchrotron radiation. Although the theoretical and experimental results of this paper imply the use of highly coherent synchrotron sources, the technique is in fact amenable to use with laboratory sources using the coded aperture technique [3133].

2. Derivation of quantitative method

2.1. Outline of the problem

In the X-ray regime, the refractive index is normally expressed as n = 1 − δ + where δ and β are the refractive index decrements. δ can be as large as 10−6 whilst β can be as much as three orders of magnitude less than δ. A wavefront which propagates through a thin object will be perturbed by the object’s complex transmission function T(x,y) = exp(−(x,y) − μ(x,y)) where ϕ and μ, both real valued, are defined as:

ϕ(x,y)=k𝒪δ(x,y,z)dzμ(x,y)=k𝒪β(x,y,z)dz
where 𝒪 is the extent of the object and k is the wave number. μ results in attenuation of the wavefront amplitude, which is directly detectable. Please note that our definition of μ differs from the attenuation coefficient as traditionally used in X-ray imaging which is represented by the same symbol. For the sake of reference, the attenuation coefficient is given by 2(x,y,z). Generally only the gradient or Laplacian of ϕ can be measured using one of the previously mentioned quantitative XPCI techniques.

As mentioned previously, the proposed quantitative imaging technique is based upon the EIXPCI method [28]. The system is designed to perform XPCI with laboratory sources, however in the first instance we develop the theory of quantitative imaging assuming an idealised monochromatic point source. In practice, a synchrotron X-ray source can be assumed to closely approximate this condition. EIXPCI is generally performed using a system illustrated schematically in Fig. 1(b), which employs two sets of apertures. The set of apertures labeled A1 in Fig. 1(b) is called the pre-sample aperture as it is placed immediately upstream of the sample. The set of apertures labeled A2 is called the detector aperture as it is positioned as close to the flat panel detector as possible. The projected pitches of A1 and A2 match that of the flat panel pixel width.

The key to EIXPCI is the principle of pixel edge illumination [28, 34]. This condition is achieved by offsetting apertures A1 and A2 from each other as demonstrated in Fig. 1(b). The partial illumination condition is usually quantified by the illuminated pixel fraction (IPF) [29] which is defined as the total integrated X-ray intensity normalised by its maximum. The maximum pixel intensity is usually observed when the two apertures are in alignment. Although an EIXPCI system sensitive to phase gradients in two dimensions has been developed [35, 36], the system in Fig. 1(b) is sensitive to phase gradients in one direction only. As a result it is henceforth assumed, for the mathematical derivations only, that the imaging system and objects extend uniformly out of the page in Fig. 1(b). In practice, however, the stationary phase approximation [37, Pgs. 29–34] shows that to within a good approximation, a region of a sample generates contrast only in the row of pixels that it projects onto. This allows the assumption of a sample uniform in one direction to be made whilst at the same time forming an image of the sample.

2.2. Derivation for object with linear phase and absorption profile

For simplicity and without loss of generality, while deriving the equations for quantitative EIX-PCI we consider a single pair of apertures as demonstrated in Fig. 2. Quantitative EIXPCI is performed by imaging the same object using two complementary positions of A2. Both configurations are set to an IPF of 0.5, yet will result in inverted contrasts. For a non-absorbing, prism-like, object depicted in Fig. 2, the configuration denoted by I will result in a detected signal which is lower than the flat field whilst the I+ configuration will result in a detected signal which is greater than the flat field. If an absorbing object were employed, I and I+ would also be affected by absorption. We will derive the equations of the quantitative method by assuming a monochromatic point source. Within the paraxial (i.e., small angle) approximation, a point source results in an X-ray intensity which is uniform in the y direction of Fig. 2. Then, by assuming that the imaging system and object are also uniform in the y direction we need only consider variations of intensity in the x direction. If the complex amplitude of X-rays incident upon A2 is given by U(x) then, assuming a pixel of height P, I and I+ are given by:

I=0PMW0|U(x)|2dxdyI+=0P0MW|U(x)|2dxdy
where M is the system magnification given by (zso + zod)/zso and zso and zod are defined in Fig. 2.

 figure: Fig. 2

Fig. 2 The two complementary aperture configurations used to perform quantitative EIX-PCI. W is the width of the pre-sample apertures, A1, M = (zso + zod)/zso is the magnification, ξ is the spatial coordinate within the open region of the pre-sample aperture and Iand I+ represent the signals measured by the detector pixels in the left and right configurations respectively. A2 represents the detector aperture.

Download Full Size | PDF

We now turn to the task of calculating the complex amplitude U(x). We begin by dropping the y dependence of the complex transmission function thus re-writing it as T(x). Furthermore, in the object space we use the symbol ξ instead of x to distinguish between the object and detector spaces. We also assume that the phase and absorption functions may be described locally by a linear function giving the complex transmission function as:

T(ξ+ξs;ξi)=exp(iϕ(ξs)iϕξ|ξs(ξξi)μ(ξs)μξ|ξs(ξξi))
where the parameter ξs is used to represent a shift of the sample, the Taylor series expansion is taken about ξ = 0 and the parameter ξi is used to shift the Taylor series representation of the transmission function. The series expansion is centered on ξ = 0 as this coincides with the centre of the pre-sample aperture. Note that Eq. (3) means that we assume that the sample has constant phase and absorption gradients within the transmitting region of the pre-sample aperture. Then, by applying the paraxial approximation to the Fresnel-Kirchhoff diffraction integral the complex amplitude incident upon A2 may be found according to [38]:
U(x)=CW/2W/2T(ξξs;ξi)exp(ikξ2zso+zod2zsozod)exp(ikξxzod)dξ=Cexp(i(ϕξ|ξsξiϕ(ξs)))exp(μξ|ξsξiμ(ξs))W/2W/2exp(μξ|ξsξ)exp(ikξ2zso+zod2zsozodikξ(xzod+1kϕξ|ξs))dξ
where, since we are ultimately interested in |U(x)|, complex phase terms dependent upon y have been omitted and
C=U0iλzsozod(zso+zod)exp(ik(zso+zod))exp(ik(x22zod)).
Where U0 represents the amplitude of the incident wave emitted by the source. As in Ref. [38] we employ the stationary phase approximation (SPA) to yield the first order asymptotic solution to Eq. (4) as
U0~Cexp(i(ϕξ|ξsξiϕ(ξs)))exp(μξ|ξsξiμ(ξs))iλzsozodzso+zodexp[μξ|ξszsozodzso+zod(xzod+1kϕξ|ξs)]exp[ik2zsozodzso+zod(xzod+1kϕξ|ξs)2]
which is valid for x ∈ [−MW/2 − 1/k(∂ϕ/∂ξ)|ξszod, MW/2 − 1/k(∂ϕ/∂ξ)|ξszod]. Outside of this range the first order asymptotic solution is zero. Returning now to the ultimate objective of evaluating I and I+, substituting Eq. (6) into Eqs. (2) yields:
I=𝒜[exp(μξ|ξsW)]I+=𝒜[exp(μξ|ξsW)+]
where 𝒜 = |U0|2Pexp(2(∂μ/∂ξ)|ξs ξi − 2μ(ξs))/((zso + zod)2zso (μ/∂ξ)|ξs) and ℬ = exp(−2(∂μ/∂ξ)|ξs zsozod(1/k)(∂ϕ/∂ξ)|ξs/(zso + zod)). Any X-ray transmitted by A1 and not absorbed by the sample should be detected in either the I or I+ configurations, so long as it is not scattered beyond the sensitive region of the pixel of both configurations. Further, any redistribution of X-rays between the I and I+ configurations not due to absorption may be assumed to be due to a sample phase gradient. These assumptions, and the mathematical form of Eqs. (7), lead us to form the sum and difference of these two quantities as:
I++I=2𝒜sinh(μξ|ξsW)
I+I_=2𝒜[cosh(μξ|ξsW)]
Note that Eq. (8) is equivalent to a conventional absorption image and an identical result could be obtained by integrating the X-ray intensity over the equivalent region of the sample aperture. This also follows from conservation of energy within the paraxial limit. In the case that there is no sample or the sample is, not absorbing, Eq. (8) simplifies to:
I0=|U0|2MWP(zso+zod)2
which we denote by I0 as it is equivalent to the flat field in a traditional absorption image. In order to perform quantitative XPCI we must extract the quantity 1/k (∂ϕ/∂ξ)|ξs. We proceed by noting that an image, as opposed to a profile, may be obtained using one of two methods, either independently or in combination. One method makes use of periodic apertures and a flat panel detector as demonstrated in Fig. 1(b). The other method employs a technique known as dithering whereby the object is scanned relative to the imaging system. These approaches are equivalent and, in keeping with the current derivation, we choose to scan the object through a single aperture pair. We thus assume that we have access to I and I+ for object positions ξi = −Δξ, 0, Δξ. This allows us to determine (∂μ/∂ξ)|ξs as:
μξ|ξs=14Δξlog((I++I)|ξi=Δξ(I++I)|ξi=Δξ)
If total acquisition time is limited we may have access to only two dithering steps ξi = −Δξ, 0 in which case (∂μ/∂ξ)|ξs may be found slightly less accurately as:
μξ|ξs=12Δξlog((I++I)|ξi=0(I++I)|ξi=Δξ)
It is then easy to obtain the total absorption function at ξ = 0 as:
μ(ξs)=12log[I0(I++I)|ξi=0sinh(μξ|ξsW)μξ|ξsW]
12log[I0(I++I)|ξi=0]
Having evaluated (∂μ/∂ξ)|ξs and since W is known, (1/k)(∂ϕ/∂ξ)|ξs may be found according to:
(1/k)ϕξ|ξs=log[I+II++I|ξi=0sinh(μξ|ξsW)+cosh(μξ|ξsW)](zso+zod)2μξ|ξszsozod
By making the approximations sinh(x) ≈ x, cosh(x) ≈ 1 + x2/2 and log(1 + x) ≈ x for x < 0.5, Eq. (15) can be found in a simpler form, valid so long as μξ|ξsW<1/2, as:
(1/k)ϕξ|ξs=[I+II++I|ξi=0+12μξ|ξsW]W(zso+zod)2zsozod
Finally, we note that we have hitherto tacitly assumed that β is non-zero. When β is zero, the entire derivation is considerably simpler resulting in:
(1/k)ϕξ|ξs=I+II++I|ξi=0W(zso+zod)2zsozod
meaning that the phase gradient can be extracted without performing the central difference scheme in Eq. (11). Thus the main results of this paper have been derived, in particular, Eqs. (11), (13) and (15) give the gradient of the absorption function, the total absorption function and the gradient of the phase function respectively.

3. Examples and analysis

3.1. Limitations of the technique

The first restriction of the technique is that the object’s phase and absorption functions (i.e., ϕ and μ) should be well approximated locally by a linear function. Furthermore, this linear approximation should be reasonable over the range of the sample aperture opening (W). This requirement is, however, common to all grating techniques since a projected fringe samples the average phase and absorption properties over a region of a sample. Secondly, Eqs. (11) and (15) were derived under the assumption that −MW/2 − 1/k(∂ϕ/∂ξ)|ξs zod ≥ − MW and MW/2 − 1/k(∂ϕ/ξ)|ξs zodMW. This assumption is employed in the evaluation of Eqs. (7) but it also follows from energy conservation. In particular, if photons deviate outside of the union of the two exposed pixel regions in Fig. 2 then the quantity I+ + I will not be equivalent to the conventional absorption image of the object. These restrictions require that:

(1/k)|ϕξ|ξsMW2zod.
For the experimental system described in Sec. (3.2) this would require an object to have a gradient not exceeding 1.9×10−5/δ m−1. Given that δ is seldom greater than 10−6 this restriction will be satisfied in most practical applications. Furthermore, if a typical laboratory based system were used with zso = 1.6m and zod = 0.4m [36], the object slope should not exceed the less restrictive value of 4 × 10−5/δ, assuming the same value of W.

The next restriction concerns the gradient of the absorption function μ(ξ). For an object with properties described precisely by Eq. (3), the measured values of I+ and I will differ slightly from those predicted by Eqs. (7) due to the SPA employed in the evaluation of Eq. (6). I+ and I are, however, very well approximated by the expressions in Eqs. (7) and we demonstrate this by an example, employing a wedge-like object with thickness described by h(ξ) = 5ξ + 1.5 × 10−4 where all quantities are in m. δ was fixed at 10−6 whilst β was varied. The system parameters were zso = 1.6m, zod = 0.4m, W = 40μm, Δξ = 10μm and the source had a photon energy of 20keV. In this analysis, the phase extraction algorithm was applied to data found by direct numerical evaluation of Eq. (4), rather than applying the SPA. We then seek to quantify the error in (∂ϕ/∂ξ)|ξs which results from employing the SPA in deriving the phase and absorption extraction formulae.

For simplicity, we introduce the quantity ρ to denote the ratio (I+I)/(I+ + I), evaluated using the SPA. We then introduce the quantity ε to represent the perturbation to ρ, which becomes ρ(1+ε), when I+ and I are calculated by evaluating Eq. (4) directly. Figure 3 contains a plot of ε for a range of values of β and shows that the SPA indeed leads to a low value of ε across the range of values of β considered. ε decreases as β increases to 10−7. This is because the oscillations present in the actual X-ray intensity incident upon the detector aperture, but not represented by the SPA, become smaller due to the increased absorption. ε begins to increase again when the absorption becomes so high that the oscillation in the X-ray intensity achieves a magnitude on the same order as the average intensity predicted using the SPA.

 figure: Fig. 3

Fig. 3 Plot of the error, ε, introduced into the ratio (I+I)/(I+ + I) by the SPA as β varies for a wedge-like object.

Download Full Size | PDF

In the next example, we show that the value of (∂ϕ/∂ξ)|ξs, extracted from field data calculated by direct evaluation of Eq. (4), departs steadily from the true value as β increases, as shown in Fig. 4, even though ε decreases. In the absence of noise, the technique calculates (μ/∂ξ)|ξs very accurately. Thus, the error in the predicted value of (∂ϕ/∂ξ)|ξs is due to the absolute perturbation to ρ which occurs when field data calculated using direct evaluation of Eq. (4) is used. The absolute error in (∂ϕ/∂ξ)|ξs is thus the difference between the values found when ρ and ρ(1 +ε) are used as input to Eq. (15). By assuming that, for large x, cosh(x) ≈ sinh(x) and for small x, log(1+x) ≈ x, it can be shown that, to a good approximation, the absolute error in (∂ϕ/∂ξ)|ξs is given by:

ϕξ|ξsϕξ|ξs˜ε(zso+zod)4zsozod(1II+)(1/k)μξ|ξs
where the tilde is used to represent the phase gradient extracted from values of I and I+ found by evaluating Eq. (4) directly. Examining this equation it is evident that the error is introduced due to the large disparity between I+ and I which arises for large values of (∂μ/∂ξ)|ξs. This is demonstrated for this example by the plot of Fig. 5 which shows how (1 − I/I+)/((1/k)(∂μ/∂ξ)|ξs) varies with β for the wedge-like example considered. The acceptable value for (1 − I/I+)/((1/k)(∂μ/∂ξ)|ξs) will vary with the accuracy required for the extracted phase. For design purposes, if an upper bound is placed upon (∂μ/∂ξ)|ξs, the maximum allowable value of W can be found easily by assuming δ = 0 and substituting Eqs. (7) into the expression (1 − I/I+)/((1/k)(∂μ/∂ξ)|ξs).

 figure: Fig. 4

Fig. 4 Plot of the value of (1/k) (∂ϕ/∂ξ)|ξs calculated for the example wedge-like object for a range of values of β. The value is calculated by simulating the values of I+ and I by evaluating the diffraction integrals of Eq. (4). The actual value of (1/k) (∂ϕ/∂ξ|ξs is 5 × 10−6.

Download Full Size | PDF

 figure: Fig. 5

Fig. 5 Plot of (1 − I/I+)/((1/k) (∂μ/∂ξ)|ξs) for the wedge-like object for a range of values of β.

Download Full Size | PDF

As a final restriction we note that even in the limit of zero absorption, the estimated value for (∂ϕ/∂ξ)|ξs in Fig. 4 is slightly below the real value of 5×10−6. This comparatively small error is again due to the SPA which does not model diffraction by the pre-sample aperture. As a result, the estimated value of (∂ϕ/∂ξ)|ξs is seen to have a small dependence on wavelength. Figure 6 contains a plot of (∂ϕ/∂ξ)|ξs against photon energy for the wedge-like object considered previously with β = 10−8 and all other parameters as specified previously. In practice, this error will be less significant than in the example as the effects of diffraction are less significant when sources with finite spectral width and spatial extent are employed.

 figure: Fig. 6

Fig. 6 Plot of how (∂ϕ/∂ξ)|ξs varies with the source photon energy for the wedge-like object example considered in this section with β = 10−8.

Download Full Size | PDF

3.2. Experimental setup

Experimental verification of the phase extraction technique was performed on the SYRMEP bending magnet of the Elettra synchrotron radiation facility in operation in Trieste, Italy. The beamline is described in more detail elsewhere [39, 40], we give here only the details relevant to the experiment. The sample stage was located, in the experimental room, approximately 22m from the apparent X-ray source which has full width at half maximum dimensions of 0.280(horizontal)×0.080(vertical) mm2. This results in a usable beam of dimension approximately 120×4 mm2. A channel-cut Si (1,1,1) crystal monochromates the beam to nominal photon energy of 20 keV with a fractional bandwidth of 0.2%. A photon counting, linear array silicon microstrip detector known as PICASSO was employed [41]. The detector works in the so-called “edge-on” configuration and provides an array of 2368 pixels 50μm wide and 300μm high. An important property of the PICASSO detector is that global and channel specific thresholds can be adjusted in order to have negligible pixel cross talk [42, 43].

The experimental system is shown schematically in Fig. 7. The diagram shows the narrow dimension of the beam and thus the beam extends 120mm into the page. The system in Fig. 7 is sensitive to the phase gradients in the vertical direction and this is also the direction in which the sample was scanned. A 20μm slit was used as the pre-sample aperture, thus setting the value of W to 20μm. A single edge was sufficient to act as the detector aperture and this was simply translated vertically to switch between the I+ and I configurations. A single edge is able to be used since a single X-ray beam is used, in contrast to the plurality of beams in the system shown in Fig. 1(b). If the system in Fig. 1(b) were employed, a larger value of W would need to be employed to ensure that nearly all photons which propagate through the object are counted in either I+ or I. With the edge in place against the detector, each pixel recorded approximately 4.5×104 counts per second.

 figure: Fig. 7

Fig. 7 Schematic, side-on view of the experimental system employed to verify the phase extraction algorithm. The diagram is not to scale.

Download Full Size | PDF

Figure 8 shows a photograph of the sample which was used in the experiment. It was composed of six filaments of varying diameter and material properties as outlined in Table 1 in which each filament is labeled with a letter from A to F which is used throughout the remainder of this paper. The diameters were estimated from the experimental results shown in Fig. 9. The sample holder was positioned in the experimental setup of Fig. 7 such that the wires were approximately normal to the page and were thus roughly parallel with the beam and aperture edges. 770 vertical scan positions were taken for both the I+ and I configurations. All 770 scan positions were acquired for the I+ configuration before changing into the I configuration and returning the sample to its initial position. A scan step of 10μm was employed however several scan steps were passed, without acquiring a detector signal, when the wires were not within the beam. An exposure time of 1s was employed and each of I+ and I took approximately 35 minutes in total to acquire. A motorised sample stage (Newport IMS300V) was employed which has a unidirectional reputability of 0.5μm and a guaranteed on-axis accuracy of 10±5μm. The reputability of the sample position is important to ensure that the I+ and I images are correctly registered.

 figure: Fig. 8

Fig. 8 Photograph of the sample used in the experiment.

Download Full Size | PDF

 figure: Fig. 9

Fig. 9 Images of the sample taken in the I+ and I configurations. Note that some vertical scan positions, in between filaments, were omitted to reduce the total number of acquisitions. The image shows 16mm of the sample in the horizontal direction. Each row of pixels was normalised by the mean of the rows within the boxed region in each image. The filaments are arranged in the same order as in Fig. 8.

Download Full Size | PDF

Tables Icon

Table 1. Properties of wires making up the sample. PEEK stands for polyetheretherketone. The boron is supplied on a 5μm diameter tungsten core. The filaments were sourced from Goodfellow (Goodfellow Cambridge Ltd., Huntingdon, UK) and the nominal diameters were obtained from data sheets provided by Goodfellow.

3.3. Experimental results

Figure 9 shows images representing the raw data acquired during the experiment with the filaments arranged in the same order as in Fig. 8. Each row of pixels represents an acquisition from the detector for a particular vertical position of the sample. The I+ and I images were normalised independently by a row of pixel values found by averaging over a sample-free region, denoted by boxes in each image. The images exhibit a transition from predominantly phase contrast, in the case of boron (lower filament), through to strong absorption contrast in the case of titanium (upper filament). Note that the 5μm diameter tungsten core is clearly visible in the image of the boron filament. Figure 10 shows line plots of the intensity along a single column of both images in Fig. 9 and a plot of their sum. Note that the left side of the line plots correspond to the upper part of the images in Fig. 9.

 figure: Fig. 10

Fig. 10 Profile plots for a single pixel column for the absorption image ((I + I+)/2), I+ and I.

Download Full Size | PDF

Figure 11 contains images of both μ and ∂ϕ/∂ξ which have been found by using the data from Fig. 9 as input to Eqs. (14) and (15) respectively. Note that all filaments are easily resolved in the phase image, whilst only the three most absorbing filaments are resolvable in the absorption image. The absorption profiles for filaments D, E and F (i.e., the right most three) in Fig. 10 exhibit anomalies in the form of erroneous peaks. It is likely that these are due to small misalignments arising due to the I+ and I images being taken during separate scan sequences. It should, however, be noted that this anomaly impacts minimally upon extracted phase in Fig. 12. This demonstrates the robustness of the phase retrieval algorithm.

 figure: Fig. 11

Fig. 11 Images of μ (left) and (1/k)∂ϕ/∂ξ (right) derived from Fig. 9 according to Eqs. (14) and (15) respectively.

Download Full Size | PDF

 figure: Fig. 12

Fig. 12 Plots comparing the extracted value (broken line) for (1/k)∂ϕ/∂ξ for each wire and the theoretical value (solid line). The dots superimposed on the broken line represent the experimentally determined values of (1/k)∂ϕ/∂ξ. The wires are labeled with a letter which corresponds to a row in Table 1.

Download Full Size | PDF

Fig. 12 compares the extracted phase gradient with the analytic value for (∂ϕ/∂ξ)|ξs found by differentiating ϕ in Eq. (1) directly, assuming a cylindrical sample. The analytic value for the boron filament was calculated assuming a 5μm diameter tungsten core. Each of the analytically calculated phase gradients reach infinity at the filament edges. The extracted phase gradient does not, however, reach infinity. There are a few reasons for this, firstly, a degree of smoothing is inevitable since the X-ray beam used in the imaging technique has a finite width and thus senses a segment of the filament. Secondly, the extracted value of (∂ϕ/∂ξ)|ξs relies upon the extracted value of (∂μ/∂ξ)|ξs which is calculated by a central difference calculation. The central difference calculation is accurate only for smooth functions and breaks down when a sharp transition occurs, such as at the the filament edge. Despite these limitations, the extracted phase gradient is very close to the analytic phase gradient. Plot F of Fig. 12 contains two sharp peaks around ξ = 0 due to the 5μm tungsten core on which the boron is constructed. The core is below the spatial resolution of the imaging system and thus why it appears as a spike.

4. Conclusions

We have presented the derivation, analysis and experimental verification of a quantitative X-ray phase contrast imaging technique. The experiment was performed using synchrotron radiation, yet the technique has been adapted to be used with laboratory sources using the coded aperture technique [33]. The synchrotron experiments were important in establishing the absolute accuracy of the method by eliminating spectral averaging. The experimentally determined phase gradients matched the theoretically expected results very closely. The technique also works well in close proximity to the infinite phase gradients near the filament edges. This work provides the theoretical foundation for a number of experiments and theoretical developments, in both planar and three dimensional imaging, currently in progress within our laboratory.

Acknowledgments

This work was funded by the UK Engineering and Physical Sciences Research Council ( EP/G004250/1 and EP/I021884/1). K.I. was supported by the Wellcome Trust ( 085856/Z/08/Z) and P.M. is supported by a Discovery Early Career Research Award from the Australian Research Council ( DE120101331). F.L. is supported by the Prof. Giulio Brautti PhD Memorial Fellowship. We would like to thank the personnel from the Elettra synchrotron and the University of Trieste working on the SYRMEP beamline for assistance with obtaining the experimental results.

References and links

1. U. Bonse and M. Hart, “An X-ray interferometer,” Appl. Phys. Lett. 6(8), 155–156 (1965). [CrossRef]  

2. T. Davis, D. Gao, T. Gureyev, A. Stevenson, and S. Wilkins, “Phase-contrast imaging of weakly absorbing materials using hard X-rays,” Nature 373(6515), 595–598 (1995). [CrossRef]  

3. A. Momose, “Demonstration of phase-contrast X-ray computed tomography using an X-ray interferometer,” Nucl. Instrum. Meth. A 352(3), 622 – 628 (1995). [CrossRef]  

4. A. Snigirev, I. Snigireva, V. Kohn, S. Kuznetsov, and I. Schelokov, “On the possibilities of X-ray phase contrast microimaging by coherent high-energy synchrotron radiation,” Rev. Sci. Instrum. 66(12), 5486–5492 (1995). [CrossRef]  

5. D. Chapman, W. Thomlinson, F. Arfelli, N. Gmür, Z. Zhong, R. Menk, R. E. Johnson, D. Washburn, E. Pisano, and D. Sayers, “Mammography imaging studies using a Laue crystal analyzer,” The 9th National Conference on Synchrotron Radiation Instrumentation 67(9), 3360–3360 (1996).

6. P. Cloetens, R. Barrett, J. Baruchel, J.-P. Guigay, and M. Schlenker, “Phase objects in synchrotron radiation hard X-ray imaging,” J. Phys. D Appl. Phys. 29, 133–146 (1996). [CrossRef]  

7. D. Chapman, W. Thomlinson, R. Johnston, D. Washburn, E. Pisano, N. Gmür, Z. Zhong, R. Menk, F. Arfelli, and D. Sayers, “Diffraction enhanced X-ray imaging,” Phys. Med. Biol. 42(11), 2015–2025 (1997). [CrossRef]   [PubMed]  

8. D. Paganin, T. E. Gureyev, K. M. Pavlov, R. A. Lewis, and M. Kitchen, “Phase retrieval using coherent imaging systems with linear transfer functions,” Opt. Commun. 234, 87 – 105 (2004). [CrossRef]  

9. Y. I. Nesterets, T. Gureyev, D. Paganin, K. Pavlov, and S. W. Wilkins, “Quantitative diffraction-enhanced X-ray imaging of weak objects,” J. Phys. D Appl. Phys. 37(8), 1262–1274 (2004). [CrossRef]  

10. M. J. Kitchen, D. M. Paganin, K. Uesugi, B. J. Allison, R. A. Lewis, S. B. Hooper, and K. M. Pavlov, “X-ray phase, absorption and scatter retrieval using two or more phase contrast images,” Opt. Express 18(19), 19,994–20,012 (2010). [CrossRef]  

11. P. C. Diemoz, P. Coan, C. Glaser, and A. Bravin, “Absorption, refraction and scattering in analyzer-based imaging: comparison of different algorithms,” Opt. Express 18, 3494–3509 (2010). [CrossRef]   [PubMed]  

12. D. J. Vine, D. M. Paganin, K. M. Pavlov, J. Kraeusslich, O. Wehrhan, I. Uschmann, and E. Foerster, “Analyzer-based phase contrast imaging and phase retrieval using a rotating anode X-ray source,” Appl. Phys. Lett. 91(25), 254110 (2007). [CrossRef]  

13. S. Wilkins, T. Gureyev, D. Gao, A. Pogany, and A. Stevenson, “Phase-contrast imaging using polychromatic hard X-rays,” Nature 384, 335–338 (1996). [CrossRef]  

14. K. Nugent, T. Gureyev, D. Cookson, D. Paganin, and Z. Barnea, “Quantitative phase imaging using hard X-rays,” Phys. Rev. Lett. 77, 2961–2964 (1996). [CrossRef]   [PubMed]  

15. P. Cloetens, W. Ludwig, J. Baruchel, D. Van Dyck, J. Van Landuyt, J. Guigay, and M. Schlenker, “Holotomography: quantitative phase tomography with micrometer resolution using hard synchrotron radiation X-rays,” Appl. Phys. Lett. 75(19), 2912–2914 (1999). [CrossRef]  

16. T. Gureyev, C. Raven, A. Snigirev, I. Snigireva, and S. Wilkins, “Hard X-ray quantitative non-interferometric phase-contrast microscopy,” J. Phys. D Appl. Phys. 32(5), 563–567 (1999). [CrossRef]  

17. D. Paganin, S. Mayo, T. Gureyev, P. Miller, and S. Wilkins, “Simultaneous phase and amplitude extraction from a single defocused image of a homogeneous object,” J. Microsc.-Oxford 206, 33–40 (2002). [CrossRef]  

18. T. Gureyev, A. Pogany, D. Paganin, and S. Wilkins, “Linear algorithms for phase retrieval in the fresnel region,” Opt. Commun. 231, 53–70 (2004). [CrossRef]  

19. T. Gureyev, Y. Ne, D. Paganin, A. Pogany, and S. Wilkins, “Linear algorithms for phase retrieval in the fresnel region. 2. partially coherent illumination,” Opt. Commun. 259, 569–580 (2005). [CrossRef]  

20. K. Nugent, “The measurement of phase through the propagation of intensity: an introduction,” Contemp. Phys. 52, 55–69 (2011). [CrossRef]  

21. K. Creath, “Phase-measurement interferometry techniques,” Prog. Optics XXVI, 349–393 (1988). [CrossRef]  

22. A. Lohmann and D. Silva, “An interferometer based on the Talbot effect,” Opt. Commun. 2(9), 413–415 (1971). [CrossRef]  

23. T. Weitkamp, A. Diaz, C. David, F. Pfeiffer, M. Stampanoni, P. Cloetens, and E. Ziegler, “X-ray phase imaging with a grating interferometer,” Opt. Express 13(16), 6296–6304 (2005). [CrossRef]   [PubMed]  

24. A. Momose, S. Kawamoto, I. Koyama, Y. Hamaishi, K. Takai, and Y. Suzuki, “Demonstration of X-ray Talbot interferometry,” Jpn. J. Appl. Phys. 42(Part 2, No. 7B), L866–L868 (2003). [CrossRef]  

25. F. Pfeiffer, T. Weitkamp, O. Bunk, and C. David, “Phase retrieval and differential phase-contrast imaging with low-brilliance X-ray sources,” Nat. Phys. 2(4), 258–261 (2006). [CrossRef]  

26. P. Zhu, K. Zhang, Z. Wang, Y. Liu, X. Liu, Z. Wu, S. A. McDonald, F. Marone, and M. Stampanoni, “Low-dose, simple, and fast grating-based X-ray phase-contrast imaging,” Proc. Natl. Acad. Sci. U. S. A. 107(31), 576–581 (2010). [CrossRef]  

27. T. Gureyev, S. Mayo, D. Myers, Y. Nesterets, D. Paganin, A. Pogany, A. Stevenson, and S. Wilkins, “Refracting Rontgen’s rays: propagation-based X-ray phase contrast for biomedical imaging,” J. Appl. Phys. 105(10), 102005 (2009). [CrossRef]  

28. A. Olivo, F. Arfelli, G. Cantatore, R. Longo, R. Menk, S. Pani, M. Prest, P. Poropat, L. Rigon, G. Tromba, E. Vallazza, and E. Castelli, “An innovative digital imaging set-up allowing a low-dose approach to phase contrast applications in the medical field,” Med. Phys. 28(8), 1610–1619 (2001). [CrossRef]   [PubMed]  

29. P. Munro, K. Ignatyev, R. Speller, and A. Olivo, “Source size and temporal coherence requirements of coded aperture type X-ray phase contrast imaging systems,” Opt. Express 18(19), 19681–19692 (2010). [CrossRef]  

30. Z.-F. Huang, K.-J. Kang, L. Zhang, Z.-Q. Chen, F. Ding, Z.-T. Wang, and Q.-G. Fang, “Alternative method for differential phase-contrast imaging with weakly coherent hard X-rays,” Phys. Rev. A 79(1), 013815 (2009). [CrossRef]  

31. A. Olivo and R. Speller, “A coded-aperture technique allowing X-ray phase contrast imaging with conventional sources,” Appl. Phys. Lett. 91(7), 074106 (2007). [CrossRef]  

32. A. Olivo and R. Speller, “Modelling of a novel X-ray phase contrast imaging technique based on coded apertures,” Phys. Med. Biol. 52(22), 6555–6573 (2007). [CrossRef]   [PubMed]  

33. P. Munro, K. Ignatyev, R. Speller, and A. Olivo, “Phase and absorption retrieval using incoherent X-ray sources,” Proc. Natl. Acad. Sci. U. S. A. 109, 13,922–13,927 (2012). [CrossRef]  

34. A. Olivo and R. Speller, “Image formation principles in coded-aperture based X-ray phase contrast imaging,” Phys. Med. Biol. 53(22), 6461–6474 (2008). [CrossRef]   [PubMed]  

35. A. Olivo, S. E. Bohndiek, J. A. Griffiths, A. Konstantinidis, and R. D. Speller, “A non-free-space propagation X-ray phase contrast imaging method sensitive to phase effects in two directions simultaneously,” Appl. Phys. Lett. 94(4), 044108 (2009). [CrossRef]  

36. P. Munro, K. Ignatyev, R. Speller, and A. Olivo, “Design of a novel phase contrast X-ray imaging system for mammography,” Phys. Med. Biol. 55(14), 4169–4185 (2010). [CrossRef]   [PubMed]  

37. G. James, Geometrical theory of diffraction for electromagnetic waves (Peter Peregrinus Ltd., 1976).

38. P. Munro, K. Ignatyev, R. Speller, and A. Olivo, “The relationship between wave and geometrical optics models of coded aperture type X-ray phase contrast imaging systems,” Opt. Express 18(5), 4103–4117 (2010). [CrossRef]   [PubMed]  

39. F. Arfelli, A. Bravin, G. Barbiellini, G. Cantatore, E. Castelli, M. D. Michiel, P. Poropat, R. Rosei, M. Sessa, A. Vacchi, L. D. Palma, R. Longo, S. Bernstorff, A. Savoia, and G. Tromba, “Digital mammography with synchrotron radiation,” Rev. Sci. Instrum. 66(2), 1325–1328 (1995). [CrossRef]  

40. A. Abrami, F. Arfelli, R. Barroso, A. Bergamaschi, F. Bille, P. Bregant, F. Brizzi, K. Casarin, E. Castelli, V. Chenda, L. Palma, D. Dreossi, A. Fava, R. Longo, L. Mancini, R. Menk, F. Montanari, A. Olivo, S. Pani, A. Pillon, E. Quai, S. Kaiser, L. Rigon, T. Rokvic, M. Tonutti, G. Tromba, A. Vaseotto, C. Venanzi, F. Zanconati, A. Zanetti, and F. Zanini, “Medical applications of synchrotron radiation at the SYRMEP beamline of Elettra,” Nucl. Instrum. Meth. A 548, 221–227 (2005). [CrossRef]  

41. L. Rigon, F. Arfelli, A. Astolfo, A. Bergamaschi, D. Dreossi, R. Longo, R.-H. Menk, B. Schmitt, E. Vallazza, and E. Castelli, “A single-photon counting edge-on silicon detector for synchrotron radiation mammography,” Nucl. Instrum. Meth. A 608(1, Supplement 1), S62 – S65 (2009). [CrossRef]  

42. L. Rigon, F. Arfelli, A. Bergamaschi, R. C. Chen, D. Dreossi, R. Longo, R. H. Menk, B. Schmitt, E. Vallazza, and E. Castelli, “Evaluation of charge -sharing effects on the spatial resolution of the PICASSO detector,” Nucl. Instrum. Meth. A 617(1–3), 244–245 (2010). [CrossRef]  

43. F. C. Lopez, L. Rigon, R. Longo, F. Arfelli, A. Bergamaschi, R. C. Chen, D. Dreossi, B. Schmitt, E. Vallazza, and E. Castelli, “Development of a fast read-out system of a single photon counting detector for mammography with synchrotron radiation,” J. Instrum. 6, C12031 (2011). [CrossRef]  

44. B. Henke, E. Gullikson, and J. Davis, “X-Ray interactions: photoabsorption, scattering, transmission, and reflection at E = 50–30,000 eV, Z = 1–92,” Atom. Data Nucl. Data 54(2), 181 – 342 (1993). [CrossRef]  

Cited By

Optica participates in Crossref's Cited-By Linking service. Citing articles from Optica Publishing Group journals and other participating publishers are listed here.

Alert me when this article is cited.


Figures (12)

Fig. 1
Fig. 1 (a) Demonstration of the edge illumination principle whereby pre-sample apertures, located a distance zso from the source, shape a beam which is incident upon detector apertures and detector, placed a distance zod from the pre-sample apertures. The centre of the detector apertures is offset from the centre of the beam by a distance ΔP. (b) A CAXPCI system made by tiling the setup used in (a). A1 and A2 are apertures with the same projected pitch which also matches the detector pixel width P. All parameters have the same meaning as in (a). A1 creates X-ray beams which are incident upon A2 which is used to achieve the partial pixel illumination condition.The imaging system continues uniformly out of the page.
Fig. 2
Fig. 2 The two complementary aperture configurations used to perform quantitative EIX-PCI. W is the width of the pre-sample apertures, A1, M = (zso + zod)/zso is the magnification, ξ is the spatial coordinate within the open region of the pre-sample aperture and Iand I+ represent the signals measured by the detector pixels in the left and right configurations respectively. A2 represents the detector aperture.
Fig. 3
Fig. 3 Plot of the error, ε, introduced into the ratio (I+I)/(I+ + I) by the SPA as β varies for a wedge-like object.
Fig. 4
Fig. 4 Plot of the value of (1/k) (∂ϕ/∂ξ)|ξs calculated for the example wedge-like object for a range of values of β. The value is calculated by simulating the values of I+ and I by evaluating the diffraction integrals of Eq. (4). The actual value of (1/k) (∂ϕ/∂ξ|ξs is 5 × 10−6.
Fig. 5
Fig. 5 Plot of (1 − I/I+)/((1/k) (∂μ/∂ξ)|ξs) for the wedge-like object for a range of values of β.
Fig. 6
Fig. 6 Plot of how (∂ϕ/∂ξ)|ξs varies with the source photon energy for the wedge-like object example considered in this section with β = 10−8.
Fig. 7
Fig. 7 Schematic, side-on view of the experimental system employed to verify the phase extraction algorithm. The diagram is not to scale.
Fig. 8
Fig. 8 Photograph of the sample used in the experiment.
Fig. 9
Fig. 9 Images of the sample taken in the I+ and I configurations. Note that some vertical scan positions, in between filaments, were omitted to reduce the total number of acquisitions. The image shows 16mm of the sample in the horizontal direction. Each row of pixels was normalised by the mean of the rows within the boxed region in each image. The filaments are arranged in the same order as in Fig. 8.
Fig. 10
Fig. 10 Profile plots for a single pixel column for the absorption image ((I + I+)/2), I+ and I.
Fig. 11
Fig. 11 Images of μ (left) and (1/k)∂ϕ/∂ξ (right) derived from Fig. 9 according to Eqs. (14) and (15) respectively.
Fig. 12
Fig. 12 Plots comparing the extracted value (broken line) for (1/k)∂ϕ/∂ξ for each wire and the theoretical value (solid line). The dots superimposed on the broken line represent the experimentally determined values of (1/k)∂ϕ/∂ξ. The wires are labeled with a letter which corresponds to a row in Table 1.

Tables (1)

Tables Icon

Table 1 Properties of wires making up the sample. PEEK stands for polyetheretherketone. The boron is supplied on a 5μm diameter tungsten core. The filaments were sourced from Goodfellow (Goodfellow Cambridge Ltd., Huntingdon, UK) and the nominal diameters were obtained from data sheets provided by Goodfellow.

Equations (19)

Equations on this page are rendered with MathJax. Learn more.

ϕ ( x , y ) = k 𝒪 δ ( x , y , z ) d z μ ( x , y ) = k 𝒪 β ( x , y , z ) d z
I = 0 P M W 0 | U ( x ) | 2 d x d y I + = 0 P 0 M W | U ( x ) | 2 d x d y
T ( ξ + ξ s ; ξ i ) = exp ( i ϕ ( ξ s ) i ϕ ξ | ξ s ( ξ ξ i ) μ ( ξ s ) μ ξ | ξ s ( ξ ξ i ) )
U ( x ) = C W / 2 W / 2 T ( ξ ξ s ; ξ i ) exp ( i k ξ 2 z s o + z o d 2 z s o z o d ) exp ( i k ξ x z o d ) d ξ = C exp ( i ( ϕ ξ | ξ s ξ i ϕ ( ξ s ) ) ) exp ( μ ξ | ξ s ξ i μ ( ξ s ) ) W / 2 W / 2 exp ( μ ξ | ξ s ξ ) exp ( i k ξ 2 z s o + z o d 2 z s o z o d i k ξ ( x z o d + 1 k ϕ ξ | ξ s ) ) d ξ
C = U 0 i λ z s o z o d ( z s o + z o d ) exp ( i k ( z s o + z o d ) ) exp ( i k ( x 2 2 z o d ) ) .
U 0 ~ C exp ( i ( ϕ ξ | ξ s ξ i ϕ ( ξ s ) ) ) exp ( μ ξ | ξ s ξ i μ ( ξ s ) ) i λ z s o z o d z s o + z o d exp [ μ ξ | ξ s z s o z o d z s o + z o d ( x z o d + 1 k ϕ ξ | ξ s ) ] exp [ i k 2 z s o z o d z s o + z o d ( x z o d + 1 k ϕ ξ | ξ s ) 2 ]
I = 𝒜 [ exp ( μ ξ | ξ s W ) ] I + = 𝒜 [ exp ( μ ξ | ξ s W ) + ]
I + + I = 2 𝒜 sinh ( μ ξ | ξ s W )
I + I _ = 2 𝒜 [ cosh ( μ ξ | ξ s W ) ]
I 0 = | U 0 | 2 M W P ( z s o + z o d ) 2
μ ξ | ξ s = 1 4 Δ ξ log ( ( I + + I ) | ξ i = Δ ξ ( I + + I ) | ξ i = Δ ξ )
μ ξ | ξ s = 1 2 Δ ξ log ( ( I + + I ) | ξ i = 0 ( I + + I ) | ξ i = Δ ξ )
μ ( ξ s ) = 1 2 log [ I 0 ( I + + I ) | ξ i = 0 sinh ( μ ξ | ξ s W ) μ ξ | ξ s W ]
1 2 log [ I 0 ( I + + I ) | ξ i = 0 ]
( 1 / k ) ϕ ξ | ξ s = log [ I + I I + + I | ξ i = 0 sinh ( μ ξ | ξ s W ) + cosh ( μ ξ | ξ s W ) ] ( z s o + z o d ) 2 μ ξ | ξ s z s o z o d
( 1 / k ) ϕ ξ | ξ s = [ I + I I + + I | ξ i = 0 + 1 2 μ ξ | ξ s W ] W ( z s o + z o d ) 2 z s o z o d
( 1 / k ) ϕ ξ | ξ s = I + I I + + I | ξ i = 0 W ( z s o + z o d ) 2 z s o z o d
( 1 / k ) | ϕ ξ | ξ s M W 2 z o d .
ϕ ξ | ξ s ϕ ξ | ξ s ˜ ε ( z s o + z o d ) 4 z s o z o d ( 1 I I + ) ( 1 / k ) μ ξ | ξ s
Select as filters


Select Topics Cancel
© Copyright 2024 | Optica Publishing Group. All rights reserved, including rights for text and data mining and training of artificial technologies or similar technologies.