OSA's Digital Library

Optics Express

Optics Express

  • Editor: C. Martijn de Sterke
  • Vol. 18, Iss. 5 — Mar. 1, 2010
  • pp: 4103–4117
« Show journal navigation

The relationship between wave and geometrical optics models of coded aperture type x-ray phase contrast imaging systems

Peter R.T. Munro, Konstantin Ignatyev, Robert D. Speller, and Alessandro Olivo  »View Author Affiliations


Optics Express, Vol. 18, Issue 5, pp. 4103-4117 (2010)
http://dx.doi.org/10.1364/OE.18.004103


View Full Text Article

Acrobat PDF (285 KB)





Browse Journals / Lookup Meetings

Browse by Journal and Year


   


Lookup Conference Papers

Close Browse Journals / Lookup Meetings

Article Tools

Share
Citations

Abstract

X-ray phase contrast imaging is a very promising technique which may lead to significant advancements in medical imaging. One of the impediments to the clinical implementation of the technique is the general requirement to have an x-ray source of high coherence. The radiation physics group at UCL is currently developing an x-ray phase contrast imaging technique which works with laboratory x-ray sources. Validation of the system requires extensive modelling of relatively large samples of tissue. To aid this, we have undertaken a study of when geometrical optics may be employed to model the system in order to avoid the need to perform a computationally expensive wave optics calculation. In this paper, we derive the relationship between the geometrical and wave optics model for our system imaging an infinite cylinder. From this model we are able to draw conclusions regarding the general applicability of the geometrical optics approximation.

© 2010 Optical Society of America

1. Introduction

It is hoped that x-ray Phase Contrast Imaging (XPCi) will provide a generational improvement in the effectiveness of mammography [1

1. R. Lewis, “Medical phase contrast x-ray imaging: current status and future prospects,” Phys. Med. Biol. 49(16), 3573–3583 (2004). [CrossRef] [PubMed]

]. To our knowledge, the only in vivo mammography program is in progress in Trieste, Italy, using the SYRMEP beam line [2

2. E. Castelli, F. Arfelli, D. Dreossi, R. Longo, T. Rokvic, M. Cova, E. Quaia, M. Tonutti, F. Zanconati, A. Abrami, V. Chenda, R. Menk, E. Quai, G. Tromba, P. Bregant, and F. de Guarrini, “Clinical mammography at the SYRMEP beam line,” Nucl. Instrum. Meth. A 572(1), 237–240 (2007). [CrossRef]

]. This program has provided mammograms of improved spatial resolution and detail visibility compared with conventional mammography. It cannot, however, be considered a viable tool for clinical screening due to its reliance on a synchrotron source.

An alternative XPCi technique employing laboratory sources was suggested by Olivo et. al [3

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

, 4

4. A. Olivo and R. Speller, “Phase contrast imaging.”, International Patent WO/2008/029107, (2008).

] in 2007. This technique is known as coded aperture XPCi and has since been under continuous development within the radiation physics group at UCL (see references [5

5. 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) (2009). [CrossRef]

, 6

6. 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]

, 3

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

] for example). This technique has been demonstrated experimentally and validated theoretically in the aforementioned references. We are now building a pre-prototype coded aperture XPCi system in order to demonstrate the efficacy of the technique using in vitro human breast tissue samples. In order to design the system and verify the experiments, it is necessary to model the entire imaging system, including the interaction of x-rays with tissue. The small refractive index contrast of tissue combined with the unpolarised x-ray source mean that a full electromagnetic calculation for the scattered x-rays can be avoided. Furthermore, the short wavelength of x-rays relative to typical cell structure dimensions means that a geometrical optics model is often sufficient. This is important as a rigorous scalar calculation of the scattered field would require prohibitively large computational resources. In this paper, we thus attempt to establish conditions under which a geometrical optics approximation can be employed to model a coded aperture XPCi system.

2. Wave optics model

Ui(P)=C𝓐exp(ik(ξ2+ψ2)zso+zod2zsozod)exp(ik(ξx′+ψyzod))dξdψ
(1)

where

x′=x+xszod/zso
C=iλzsozodexp(ik(zso+zod))exp(ik(xs22zso+x2+y22zod))
k=2π/λ
(2)

and 𝓐 represents the transmitting regions of the sample apertures. In addition, x, ξ and z are defined in Fig. 1 and (ξ, ψ, z) and (x,y, z) form right handed coordinate systems. The integration over ψ can be performed by noting that the apertures have no dependence upon ψ. We must thus evaluate:

Γ(y)=exp(ik(zso+zod2zsozodψ2yzodψ))dψ
(3)

evaluating this integral with limits at ±∞ contradicts the Fresnel approximation used to obtain Eq. (1) and should be solved using the theory of distributions [15

15. J. Arsac, Fourier transforms and the theory of distributions (Prentice-Hall, 1966).

]. This problem can however be avoided by noting that the kernel of the integral in Eq. (3) is a rapidly oscillating function which lends itself to asymptotic evaluation by the method of stationary phase. According to the method of stationary phase [14, Pgs. 29–34], an integral of the form:

I=f(x)exp(ikg(x))dx
(4)
Fig. 1. Schematic diagram of imaging system including reference frames used in the paper. Note that (x̄,ȳ,z), (ξ, ψ, z) and (x, y, z) all form right handed coordinate systems. The imaging system is assumed to have no y dependence. Note that dL is defined by the displacement between the detector apertures and the projection of the sample apertures onto the detector apertures.

where g(x) has a single first order stationary point, x 0, such that g′(x 0) = 0, g″(x 0) ≠ 0, can be approximated as:

I02πkg(x0)f(x0)exp(i[kg(x0)+π4sgn(g(x0))])
(5)

in the limit of large k. Applying this approximation to Eq. (3) we find that

Γ(y)izsozodλzso+zodexp(iky2zso2zod(zso+zod))
(6)

the role of this term is to ensure energy conservation and give the incident field the correct phase relationship with y. This result is also obtainable using Fourier theory applied to distributions [15

15. J. Arsac, Fourier transforms and the theory of distributions (Prentice-Hall, 1966).

] which reveals that Eq. (6) is in fact the solution to Eq. (3) [17

17. A. Erdélyi, ed., Tables of integral transforms : based, in part, on notes left by Harry Bateman and compiled by the staff of the Bateman Manuscript Project, vol. I (New York; London: McGraw-Hill, 1954).

]. This enables us to write Eq. (1) as:

Ui(P)=CΓ(y)T(ξ)exp(ikξ2zso+zod2zsozod)exp(ikξx′zod)dξ
(7)

where we now introduce the periodic function T(ξ) to represent the transmission function of the sample aperture. It is now easy to include the effect of a phase object with phase function ϕ(ξ) by following an approach similar to that of Arfelli et. al [18

18. F. Arfelli, M. Assante, V. Bonvicini, A. Bravin, G. Cantatore, E. Castelli, L. Dalla Palma, M. Di Michiel, R. Longo, A. Olivo, S. Pani, D. Pontoni, P. Poropat, M. Prest, A. Rashevsky, G. Tromba, A. Vacchi, E. Val-lazza, and F. Zanconati, “Low-dose phase contrast x-ray medical imaging,” Phys. Med. Biol. 43(10), 2845–2852 (1998). [CrossRef] [PubMed]

]. The total field at the detector apertures may be found according to:

U(P)=Ui(P)+T(ξ)[exp((ξ))1]exp(ikξ2zso+zod2zsozod)exp(ikξx′zod)dξ
(8)

where is the extent of the object.

3. Efficient evaluation of wave optics field

We now turn our attention to how the expression in Eq. (7) may be efficiently evaluated. As T(ξ) is a periodic function with period L, it can be represented as a complex Fourier series written in general as:

T(ξ)=n=Cnexp(i2πnξL)
(9)

which upon substitution into Eq. (7) yields

Ui(P)=CΓ(y)iλzsozodzso+zodexp(ikx′2zso2zod(zso+zod))n=exp(iπλ(nL)2zsozodzso+zod)Cnexp(i2πnLx′zsozso+zod)
(10)

As also suggested by Engelhardt et. al [10

10. M. Engelhardt, C. Kottler, o. Bunk, C. David, C. Schroer, J. Baumann, m. Schuster, and F. Pfeiffer, “The fractional Talbot effect in differential x-ray phase-contrast imaging for extended and polychromatic x-ray sources,” J. Microsc. 232, 145–157 (2008). [CrossRef] [PubMed]

], the Fast Fourier Transform (FFT) can be used to efficiently evaluate this expression. In particular, starting with the definition of the discrete Fourier transform [19

19. G. Arfken, Mathematical Methods for Physicists, 3rd ed. (Academic Press, Boston, 1985).

]:

G(fp)=12Nκ=02N1g(tκ)exp(i2πfptκ),tκ=κT2N,p=0,1,,2N1
g(tκ)=p=02N1G(fp)exp(i2πfptκ),fp=pT,κ=0,1,,2N1
(11)

by allowing xzso/(zso + zod) to take on values κL/(2N) the summation in Eq. (10) may be evaluated, for a finite number of terms, by constructing a vector of the form:

w=[w0w1w2wNwN+1wN+2wN+1]
(12)

where

wn=Cnexp(iπλ(nL)2zsozodzso+zod)
(13)

and finally taking the Fourier transform of the vector in Eq. (12). Noting also that the coefficients Cn may also be evaluated using the FFT, Eq. (10) may be evaluated very efficiently.

The second term in Eq. (8) must, in general, be evaluated numerically unless the object has a phase function permitting analytic evaluation. It was found that Gaussian Quadrature integration [19

19. G. Arfken, Mathematical Methods for Physicists, 3rd ed. (Academic Press, Boston, 1985).

] provided accurate results.

4. Geometrical optics model

Olivo and Speller [6

6. 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]

, 20

20. 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]

] have previously used geometrical optics to model the coded aperture XPCi system. Their approach used a “forward” technique where photons emitted by the source were traced through the system. Photons could be blocked by an aperture, refracted by a sample or both. The number of photons reaching a particular pixel represent the signal detected by that pixel. We now consider the ray optics approach in a more formal manner in order to relate it to the wave optics approach. For the remainder of this section we consider only non-trivial rays which are transmitted by the sample aperture. We consider here a first order geometrical optics. It has been shown by Keller [13

13. J. B. Keller, “Geometrical Theory of Diffraction,” J. Opt. Soc. Am. A 52(2), 116–130 (1962). [CrossRef]

] and later by James [14

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

] that geometrical optics may be extended to include higher order terms which represent what is usually termed diffraction. Here we consider only the first order terms of the geometrical optics approximation. The trajectory of a light ray is described by the expression [21

21. M. Born and E. Wolf, Principles of Optics, seventh ed. (Cambridge University Press, Cambridge, 1999).

]:

ndrds=𝓢
(14)

where r is the position vector of a point on the ray, s the length of the ray, n the refractive index of the medium and 𝓢 defines a wave front of constant phase, ie, 𝓢 (ξ) = constant. It is evident from this that we assume rays are deflected in the ξ direction only. Consider a phase object as depicted in Fig. 2. We define the phase function, ϕ(ξ), as

ϕ(ξ)=kzob/2zob/2n(ξ,z)dz
(15)

where n(ξ, z) is the refractive index at position (ξ, z) and we have assumed that rays make only small angles, θi, with the z-axis. The angle by which the ray is deflected in then given by:

α(ξ)=1kdξ
(16)

With reference to Fig. 1, we can say that a ray emitted at angle θi to the z-axis will intercept the ξ-axis at position ξ = zso tan(θi) and, if deflected by an object, will intercept the x-axis at position

x=zsotan(θi)+zodtan(θi+α(ξ))
(17)

The phase of the ray at the detector apertures is calculated by taking into account the phase introduced by the object and the distance travelled in free space according to the Fresnel approximation. The amplitude of the ray must be such that energy is conserved. In particular, the time average power propagating in a small pencil of rays emanating from the source must remain constant. The ratio between ray amplitudes at z = zod and z = 0 is thus given by:

U(z=zod)U(z=0)=dξdx=11+zodzso+zodsec2(α)dαdξ
(18)
Fig. 2. Diagram illustrating the phase function ϕ(ξ) of a phase object of extent zob in the z direction.

5. Modelling a finite size source

Secs. (2)–(4) show how to calculate the field incident upon the detector apertures. The detected signal is found by integrating the intensity of x-rays transmitted by the detector apertures and incident upon a particular pixel. In general, the pth transmitting region of the detector apertures is given by [p LMLM/4+dL,pLM+LM/4+dL] where dL is the displacement of the detector apertures relative to the projection of the sample apertures as shown in Fig. 1 and M = (zod + zso)/zso is the system magnification. We assume that the pixels are aligned as shown in Fig. 1 such that a single pixel entirely covers a single transmitting region of the detector apertures. Before calculating the signal detected by each pixel, we introduce a source of finite size in the x̄ direction. The brightness is described by P(x̄) which we will take to have a Gaussian profile. We can then take the signal of the pth pixel to be given by:

S(p)=(p1/4)LM+dL(p+1/4)LM+dL[P(x̄)I(x+x̄zodzso)dx̄]dx
(19)

where, for mathematical convenience we have assumed that the source brightness profile limits the effective source size rather than the limits of integration. By making the substitution P(x̄) = exp(-(x̄/σ)2), Eq. (19) may be expressed as:

S(p)=σI(x)K(x)dx
(20)

where

K(x)=π2{erf[zsozodσ(x(p1/4)LMdL)]erf[zsozodσ(x(p+1/4)LMdL)]}
(21)

and erf(z) is the error function

erf(z)=2π0zexp(t2)dt
(22)

Equation 20 shows that K(x) may effectively be considered as a pixel sensitivity function. Figure 3 shows plots of K(x) for a variety of source Full Widths at Half Maximum (FWHM). This shows how a broad source leads to a broad K(x) thus diminishing the sensitivity of the system to fine variations in the intensity caused by phase variations in the object.

Fig. 3. Plot of K(x) for some values of source FWHM, given in the legend in μm. Values of M = 1.25 and L = 40μm were used.

6. XPCi model of a dielectric fibre

6.1. Wave optics model

ϕ(ξ)=kcylinderδdz=2R2(ξξ0)2
(23)

We consider first the wave optics calculation. We must evaluate the second term of Eq. (8) after substituting Eq. (23) into it. The bounds of integration are found by taking the intersection of the transmitting part of the sample aperture and the cylinder. Without loss of generality we consider three cases depicted in Fig. 4 where the transmitting part of the sample aperture is assumed to be centered upon ξ = 0. Note that the cases depicted in Fig. 4 do not limit the cylinder radius, all that is important is where the cylinder boundaries lie relative to the transmitting regions of the apertures. A cylinder covering more than one sample aperture could be modelled using a combination of the cases depicted in Fig. 4. In practice our system employs a series of apertures to simultaneously image a wide field of view. For clarity, we consider here a single sample/detector aperture pair and scan the object to obtain its image. Images obtained in this way will be equivalent to those obtained in practice only when photons are not scattered between differing pre-sample/detector aperture pairs. Only a simple extension is required to model the practical system as is shown at the end of Sec. (6.2). An analysis of when this approximation is valid is given in Sec. (6.3).

Fig. 4. Three cases which must be considered in order to evaluate the integral in Eq. (8).

The integration may be evaluated numerically however it is instructive to analytically evaluate by approximation. We start by writing Eq. (8) as

U(P)=Ui+CΓ(y)[U1U2]
(24)

where

U1=Ωexp[ik(ξ2M2zodx′ξzod2δR2(ξξ0)2)]dξ
U2=Ωexp[ik(ξ2M2zodx′ξzod)]dξ
(25)

where Ω = [ξ 0R, ξ 0 +R] ∩ [−/2,)/2] and η is the fill factor of the sample apertures. We now attempt to find asymptotic solutions, for large k, to the integrals in Eq. (25) by again applying the stationary phase approximation. In this case we must also consider the end points of the integrals. It is shown by James [14, Pgs. 29–34] that when the integration in Eq. (4) is evaluated over the interval [a, b] and g′(x) is non-zero and finite at the end points, a term

I1[fexp(ikg)ikg′]ab
(26)

must be included. It is shown by Murray [22, Pg. 77] that if one of the end points is a stationary point, then Eq. (26) should be omitted and a factor of 1/2 included in the term of Eq. (5). We now turn our attention to evaluating leading terms in the asymptotic expansions of the integrals in Eq. (25). We start by defining the functions g 1 (ξ) and g 2 (ξ) and finding their derivatives as:

g1=Mξ22zodx′ξzod2δR2(ξξ0)2g2=Mξ22zodx′ξzod
g′1=x′zod+zod+2δ(ξξ0)R2(ξξ0)2g′2=x′zod+zod
g″1=Mzod+2δ(ξξ0)2(R2(ξξ0)2)3/2+2δR2(ξξ0)2g″2=Mzod
(27)

it is easy to verify that when M, δ and zod are limited to those values experienced in practice, g1 (ξ 1,0) = 0 has a unique solution for every value of x′. This solution must in general be calculated numerically. This may be done efficiently by evaluating γ = [γi] = [g1(ξi)+x′/zod] where ξi = [ξi] is a discretisation of the domain [ξ 0R + ε, ξ 0 + Rε] for some small ε. This corresponds to case (2) in Fig. 4 where the entire cylinder is illuminated and thus rays are refracted to all values of x′. In cases (1) and (3), the bounds of integration are affected by the sample apertures which in turn affects the values of x′ to which rays are refracted. The stationary point of g 1 for a particular x′ may be found by interpolation with γ as the abscissa. This enables the leading term in the expansion to be calculated as

I1,02πkg″1(ξ1,0)exp(i[kg1(ξ1,0)+π4])
(28)

examination of g2 shows that g 2 has a single stationary point at ξ 2,0 = x′ /M. In the case that x′/M is within the bounds of integration of U 2, the following term is contributed by the stationary point ξ 2,0:

I2,0zod/Mexp(ikx′22Mzod)
(29)

I2,1[exp(ikg2)ikg′2]ab
(30)

which, in the special case where b = −a, becomes

I2,1zodexp(ik(Ma)2/(2zod))ik(x′2(Ma)2)(i2x′sin(kxMa/zod)2Macos(kxMa/zod))
(31)

6.2. Relationship between wave and geometrical optics models

We show here how the wave and geometrical optics solutions are related. As θi and α in Eq. (17) are small we can write Eq. (17) as

x′Mξ+zod2(ξξ0)δR2(ξξ0)2
(32)

where Eq. (23) has been substituted into Eq. (16) to find α. This expression is identical to g1 = 0 in Eq. (27), the stationary phase condition for integral U 1. It is then easy to verify that assuming identical incident field conditions, substituting Eqs. (23) and (16) into Eq. (18) results in the same magnitude as CΓ(y)I 1,0. Furthermore, substitution of the phase contributions from the phase object and the Fresnel approximation for free space propagation result in the same phase as in CΓ(y)I 1,0. This shows that the leading term in the asymptotic expansion of U 1 gives the same field as the geometrical optics approximation to the refracted field. Examination of CΓ(y)l 2,0 = 1/(zso + Zod) exp (ik((y 2 + (xxs)2)/(2(zso + zod)) + zso + zod)) shows that this is the geometrical optics field of the light which reaches the detector apertures without being refracted by the cylinder or blocked by the sample apertures. Closer examination of CΓ(y)I 2,1 shows that this is the field due to diffraction at the edges of Ω. Note that this quantity becomes infinite at the edges of the geometrical projection of Ω onto the detector apertures. This non-physical result can be remedied by modifying the stationary phase solution [23

23. R. Buchal and J. Keller, “Boundary layer problems in diffraction theory,” Commun. Pur. Appl. Math. 13, 85–114 (1960). [CrossRef]

] however this is beyond the scope of this work.

Table 6.2 shows how the intensity and complex amplitude, for the geometrical and wave optics models respectively, are calculated in each region defined in Fig. 5. Note that the geometrical optics solution provides the intensity of the field whilst the wave optics solution provides the complex amplitude of the field.

Fig. 5. Diagram illustrating the three regions which must be considered when analysing the field incident upon the detector apertures.

Table 1. Summary of components contributing to the field and intensity in the different regions for the case of wave and geometrical optics solutions respectively. Note that objects illuminated by adjacent sample apertures may be modelled by replacing each term I 1,0 and I 1,0 2 with summations, Σi I i 1,0 and Σi(I i 1,0)2, over all objects i respectively.

table-icon
View This Table

6.3. Examples and analysis

The validity of the presented model and technique depend on the tendency for photons to be scattered between adjacent sample/detector aperture pairs. It is however possible to develop a minimum bound upon the separation of apertures required to maintain validity. If a cylinder of radius R is placed with its centre at ξ = 0 in the imaging system of Fig. 1, its edge will be projected onto the position x′ = MR in the space of the detector. We are interested in knowing how quickly the field scattered by the cylinder decays away from x′ = MR. Assuming that the edge of the cylinder is illuminated, photons are refracted to values of x′ approaching ∞ and are described by the term I 1,0 defined in Eq. (28). Photons reaching a position x′ ≫ MR must be incident upon the cylinder for a value of ξ very close to, but not exceeding R. By writing ξ = Rε, ε > 0, in Eqs. (27) it is easy to find a simple analytic expression giving I 1,0 for xMR as ε tends to 0. It is then simple to show that I 1,0 2 will reduce by two orders of magnitude at a position x′ = RM + Δx′ where:

Δx′3=Mzod(1+2δzodRM)(4δ2Rzod3)
(33)

Figure 6 shows contours of Δx′ for values of R and δ encountered in practice. Δx′ may be considered the minimum separation of adjacent sample/detector aperture pairs to ensure detector apertures principally detect photons originating from their associated sample aperture. The above analysis considers only a point source. A source of finite width may be considered by noting the definition of x′ in Eqs. (2) and thus adding (W/2)zod/zso to Δx′, where W is the detector FWHM.

Fig. 6. Contours of Δx′ for a variety of values of R and δ. This diagram effectively shows the minimum separation required between adjacent sample/detector aperture pairs to ensure detector apertures principally detect photons originating from their associated sample aperture.

Figure 8 shows the intensity incident upon the detector apertures as calculated by the wave optics and geometrical optics solutions for a point source illuminating a cylinder. The cylinder has a value of δ = 10−7, a radius of 5μm and was situated with its axis at ξ = − 5μm. As is expected, the wave optics intensity exhibits oscillations resulting from interference between different field components. The geometrical optics solution is physically impossible as the sharp edge occurring at x = 0 would require the field to contain infinite spatial frequencies. Consideration of the angular spectrum of a propagating aperiodic field shows that such a field would require evanescent waves which, in our case, would have negligible magnitude such a distance from the sample apertures.

Figure 9 compares the directly calculated wave optics intensity to that calculated using the stationary phase approximation. As explained in Sec. (6.2), singularity anomalies arise in this solution which have been neglected. This plot shows that apart from these anomalies, the approximate solution agrees well with the directly calculated intensity. One can use the components which comprise the approximate solution to determine when the geometrical optics and wave optics solutions converge. This is however made difficult by the singularity anomalies present in the approximate solution and so we have opted to use a more pragmatic approach as outlined below.

One can envisage that when a source of finite width is employed, the XPCi signals predicted by the geometrical and wave optics models should converge. There are two explanations for why this should be the case. The first explanation observes that the point source intensity incident upon the detector apertures is convolved with the magnified source profile which in our case is Gaussian. This is equivalent to applying a low pass filter to the intensity distribution causing the oscillations in the wave optics intensity and the sharp transition in the geometrical optics intensity to be smoothed. The second explanation considers that, as shown in Sec. (5), a source of finite width may be modelled by a system employing a point source and equivalent detector apertures which cause the pixels to have a spatially dependent sensitivity as described by Eq. (21). Figure 3 shows that as the source broadens, so does the width of the equivalent detector aperture sensitivity function. Because of energy conservation, one would expect the geometrical and wave optics XPCi signals to converge as the source broadens. In particular, consider the plot shown in Fig. 7. This shows the difference between the intensities, incident upon the detector apertures, predicted by wave and geometrical optics. This signal has a zero mean value as required by conservation of energy. The coded aperture XPCi signal thus depends upon the domain over which the field intensity is integrated by the detector pixel. As the sensitive part of each detector pixel increases, or equivalently, as the source broadens, the geometrical and wave optics signals thus tend to converge.

Fig. 7. Plot of the difference between intensities calculated using the wave optics (full expression evaluated numerically) and geometrical optics approximations. The sensitive region of the pixel is shaded. Simulation parameters were the same as in Fig. 8.

Fig. 8. Plot of the intensity of the field incident upon the detector apertures for the geometrical and wave optics (full expression evaluated numerically) solutions. Simulation parameters used were R = 5μm and ε = 10−7, all other parameters were as described in Sec. (6.3).
Fig. 9. Plot of the intensity of the field incident upon the detector apertures as calculated using the exact and approximate wave optics formulations. Simulation parameters were the same as in Fig. 8.

Before proceeding to calculate ε it is useful to note that some approximations can provide further insight into the problem. In the case of ξ 0 = −R, g 1 in Eq. (27) can be well approximated by g1x′ξ/zod2δ2ξR for x′ > −MR, but not too close to −MR. This approximate form leads to a solution of ξ 1,0 = −2δ 2 Rzod 2/x2 for the stationary point of g 1. Substitution of ξ 1,0 back into the approximate forms of g 1 and g1 show that both of these functions have a dependence upon δ 2 R rather than each of these independently. This suggests that it is reasonable to expect ε for a particular source FWHM to be constant for constant values of δ 2 R. This is indeed the case as was verified by a large number of simulations, a small selection of which are shown in Fig. 11. This significantly simplifies the task of determining the source size for which the geometrical and wave optics signals converge. Figure 12 is a contour plot of e as a function of source FWHM and δ 2 R. The important conclusion which we can draw from this is that for our particular choice of zod and zso, as we expect a source to have a FWHM of around 50μm, the geometrical optics model will provide results consistent with those of the wave optics model. This result will make it feasible to model much larger objects.

Fig. 10. Plot wave optics (WO, full expression evaluated numerically) and geometrical optics (GO) XPCi signal traces for a cylinder of radius 5μm and δ = 10−6, for three different values of source FWHM. The cylinder is scanned from ξ 0 = −L/4 − R to ξ 0 = 0. The signals have been normalised by the signal for the object free case.
Fig. 11. Plots of ε against δ 2 R for three different values of δ and a point source (left) and a source of FWHM 50μm.
Fig. 12. Contour plot of the error between the normalised XPCi signals as calculated by geometrical and wave optics models. Source FWHM is on the vertical axis and δ 2 R is on the horizontal axis.

7. Conclusions

In this paper we have outlined the two most widely used techniques for modelling XPCi systems: wave and geometrical optics. We have used the theory developed to model the image of an infinite cylinder in a coded aperture XPCi system. This problem has practical significance as it can be tested experimentally. For this particular problem, we show how the geometric and wave optics models are related. we then show how this theory can be used to develop a guide for when the two techniques can be trusted to give consistent results.

Acknowledgements

We would like to thank Dr Elizabeth Skelton, Imperial College London, for correcting our implementation of the stationary phase approximation.

This work was supported by the Wellcome Trust (085856/Z/08/Z). A. Olivo is supported by a Career Acceleration Fellowship awarded by the UK Engineering and Physical Sciences Research Council (EP/G004250/1).

References and links

1.

R. Lewis, “Medical phase contrast x-ray imaging: current status and future prospects,” Phys. Med. Biol. 49(16), 3573–3583 (2004). [CrossRef] [PubMed]

2.

E. Castelli, F. Arfelli, D. Dreossi, R. Longo, T. Rokvic, M. Cova, E. Quaia, M. Tonutti, F. Zanconati, A. Abrami, V. Chenda, R. Menk, E. Quai, G. Tromba, P. Bregant, and F. de Guarrini, “Clinical mammography at the SYRMEP beam line,” Nucl. Instrum. Meth. A 572(1), 237–240 (2007). [CrossRef]

3.

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

4.

A. Olivo and R. Speller, “Phase contrast imaging.”, International Patent WO/2008/029107, (2008).

5.

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) (2009). [CrossRef]

6.

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]

7.

T. Gureyev and S. Wilkins, “On x-ray phase imaging with a point source,” J. Opt. Soc. Am. A 15(3), 579–585 (1998). [CrossRef]

8.

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]

9.

M. Engelhardt, J. Baumann, M. Schuster, C. Kottler, F. Pfeiffer, O. Bunk, and C. David, “High-resolution differential phase contrast imaging using a magnifying projection geometry with a microfocus x-ray source,” Appl. Phys. Lett. 90(22), 224101 (pages 3) (2007). [CrossRef]

10.

M. Engelhardt, C. Kottler, o. Bunk, C. David, C. Schroer, J. Baumann, m. Schuster, and F. Pfeiffer, “The fractional Talbot effect in differential x-ray phase-contrast imaging for extended and polychromatic x-ray sources,” J. Microsc. 232, 145–157 (2008). [CrossRef] [PubMed]

11.

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

12.

A. Peterzol, A. Olivo, L. Rigon, S. Pani, and D. Dreossi, “The effects of the imaging system on the validity limits of the ray-optical approach to phase contrast imaging,” Med. Phys. 32(12), 3617–3627 (2005). [CrossRef]

13.

J. B. Keller, “Geometrical Theory of Diffraction,” J. Opt. Soc. Am. A 52(2), 116–130 (1962). [CrossRef]

14.

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

15.

J. Arsac, Fourier transforms and the theory of distributions (Prentice-Hall, 1966).

16.

A. Olivo and R. Speller, “Experimental validation of a simple model capable of predicting the phase contrast imaging capabilities of any x-ray imaging system,” Phys. Med. Biol. 51(12), 3015–3030 (2006). [CrossRef] [PubMed]

17.

A. Erdélyi, ed., Tables of integral transforms : based, in part, on notes left by Harry Bateman and compiled by the staff of the Bateman Manuscript Project, vol. I (New York; London: McGraw-Hill, 1954).

18.

F. Arfelli, M. Assante, V. Bonvicini, A. Bravin, G. Cantatore, E. Castelli, L. Dalla Palma, M. Di Michiel, R. Longo, A. Olivo, S. Pani, D. Pontoni, P. Poropat, M. Prest, A. Rashevsky, G. Tromba, A. Vacchi, E. Val-lazza, and F. Zanconati, “Low-dose phase contrast x-ray medical imaging,” Phys. Med. Biol. 43(10), 2845–2852 (1998). [CrossRef] [PubMed]

19.

G. Arfken, Mathematical Methods for Physicists, 3rd ed. (Academic Press, Boston, 1985).

20.

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]

21.

M. Born and E. Wolf, Principles of Optics, seventh ed. (Cambridge University Press, Cambridge, 1999).

22.

J. Murray, Asymptotic analysis (Springer Verlag, 1984). [CrossRef]

23.

R. Buchal and J. Keller, “Boundary layer problems in diffraction theory,” Commun. Pur. Appl. Math. 13, 85–114 (1960). [CrossRef]

OCIS Codes
(050.1960) Diffraction and gratings : Diffraction theory
(080.0080) Geometric optics : Geometric optics
(110.7440) Imaging systems : X-ray imaging
(340.7430) X-ray optics : X-ray coded apertures

ToC Category:
Imaging Systems

History
Original Manuscript: July 13, 2009
Revised Manuscript: November 23, 2009
Manuscript Accepted: December 16, 2009
Published: February 17, 2010

Virtual Issues
Vol. 5, Iss. 6 Virtual Journal for Biomedical Optics

Citation
Peter R. Munro, Konstantin Ignatyev, Robert D. Speller, and Alessandro Olivo, "The relationship between wave and geometrical optics models of coded aperture type x-ray phase contrast imaging systems," Opt. Express 18, 4103-4117 (2010)
http://www.opticsinfobase.org/oe/abstract.cfm?URI=oe-18-5-4103


Sort:  Author  |  Year  |  Journal  |  Reset  

References

  1. R. Lewis, "Medical phase contrast x-ray imaging: current status and future prospects," Phys. Med. Biol. 49(16), 3573-3583 (2004). [CrossRef] [PubMed]
  2. E. Castelli, F. Arfelli, D. Dreossi, R. Longo, T. Rokvic, M. Cova, E. Quaia, M. Tonutti, F. Zanconati, A. Abrami, V. Chenda, R. Menk, E. Quai, G. Tromba, P. Bregant, and F. de Guarrini, "Clinical mammography at the SYRMEP beam line," Nucl. Instrum. Meth. A 572(1), 237-240 (2007). [CrossRef]
  3. 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]
  4. A. Olivo and R. Speller, "Phase contrast imaging," International Patent WO/2008/029107, (2008).
  5. 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) (2009). [CrossRef]
  6. 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]
  7. T. Gureyev and S. Wilkins, "On x-ray phase imaging with a point source," J. Opt. Soc. Am. A 15(3), 579-585 (1998). [CrossRef]
  8. 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]
  9. M. Engelhardt, J. Baumann, M. Schuster, C. Kottler, F. Pfeiffer, O. Bunk, and C. David, "High-resolution differential phase contrast imaging using a magnifying projection geometry with a microfocus x-ray source," Appl. Phys. Lett. 90(22), 224101 (2007). [CrossRef]
  10. M. Engelhardt, C. Kottler, O. Bunk, C. David, C. Schroer, J. Baumann, M. Schuster, and F. Pfeiffer, "The fractional Talbot effect in differential x-ray phase-contrast imaging for extended and polychromatic x-ray sources," J. Microsc. 232, 145-157 (2008). [CrossRef] [PubMed]
  11. 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]
  12. A. Peterzol, A. Olivo, L. Rigon, S. Pani, and D. Dreossi, "The effects of the imaging system on the validity limits of the ray-optical approach to phase contrast imaging," Med. Phys. 32(12), 3617-3627 (2005). [CrossRef]
  13. J. B. Keller, "Geometrical Theory of Diffraction," J. Opt. Soc. Am. A 52(2), 116-130 (1962). [CrossRef]
  14. G. James, Geometrical theory of diffraction for electromagnetic waves (Peter Peregrinus Ltd., 1976).
  15. J. Arsac, Fourier transforms and the theory of distributions (Prentice-Hall, 1966).
  16. A. Olivo and R. Speller, "Experimental validation of a simple model capable of predicting the phase contrast imaging capabilities of any x-ray imaging system," Phys. Med. Biol. 51(12), 3015-3030 (2006). [CrossRef] [PubMed]
  17. A. Erdélyi, ed., Tables of integral transforms: based, in part, on notes left by Harry Bateman and compiled by the staff of the Bateman Manuscript Project, vol. I (New York; London: McGraw-Hill, 1954).
  18. F. Arfelli, M. Assante, V. Bonvicini, A. Bravin, G. Cantatore, E. Castelli, L. Dalla Palma, M. Di Michiel, R. Longo, A. Olivo, S. Pani, D. Pontoni, P. Poropat, M. Prest, A. Rashevsky, G. Tromba, A. Vacchi, E. Vallazza, and F. Zanconati, "Low-dose phase contrast x-ray medical imaging," Phys. Med. Biol. 43(10), 2845-2852 (1998). [CrossRef] [PubMed]
  19. G. Arfken, Mathematical Methods for Physicists, 3rd ed. (Academic Press, Boston, 1985).
  20. 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]
  21. M. Born and E. Wolf, Principles of Optics, seventh ed. (Cambridge University Press, Cambridge, 1999).
  22. J. Murray, Asymptotic analysis (Springer Verlag, 1984). [CrossRef]
  23. R. Buchal and J. Keller, "Boundary layer problems in diffraction theory," Commun. Pur. Appl. Math. 13, 85-114 (1960). [CrossRef]

Cited By

Alert me when this paper is cited

OSA is able to provide readers links to articles that cite this paper by participating in CrossRef's Cited-By Linking service. CrossRef includes content from more than 3000 publishers and societies. In addition to listing OSA journal articles that cite this paper, citing articles from other participating publishers will also be listed.


« Previous Article  |  Next Article »

OSA is a member of CrossRef.

CrossCheck Deposited