## Nonlinear phase retrieval from single-distance radiograph |

Optics Express, Vol. 18, Issue 25, pp. 25771-25785 (2010)

http://dx.doi.org/10.1364/OE.18.025771

Acrobat PDF (2450 KB)

### Abstract

Phase contrast in the object plane of a phase object is retrieved from intensity contrast at a *single* object-detector distance. Expanding intensity contrast and phase shift in the detector plane in powers of object-detector distance, phase retrieval is performed beyond the solution to the linearized transport-of-intensity equation. The expansion coefficients are determined by the entire paraxial wave equation. The Laplacian of the phase shift in the object plane thus is written as a local expression linear in the intensity contrast and nonlinear in the phase shift in the object plane. A perturbative approach to this expression is proposed and tested with simulated phantom data.

© 2010 Optical Society of America

## 1. Introduction

*Z*and of low number density

*ρ̂*(polymers, soft biological tissue, etc.), the attenuation of highly energetic X-rays is weak yielding a poor absorption contrast. Sizable phase shifts are, however, induced to the coherent X-ray wave fronts. The latter induce intensity contrast that can be used for 2D and 3D imaging.

*z*. Retrieving phase information from intensity at a single distance is a precursor for 2D and 3D in-situ imaging. Specifically, we aim at a higher flexibility in setting

*z*to values beyond the edge-enhancement regime (from typically a few centimeters for X-ray energies above 20 keV to several ten centimeters). This flexibility is required in experiments where the sample is suspended by a large environment not accessible to the detector. As can be argued on the basis of the contrast-transfer function [1

1. A. Papoulis, “Ambiguity function in Fourier optics,” J. Opt. Soc. Am. **64**, 779 (1974). [CrossRef]

3. J.-P. Guigay, “The ambiguity function in diffraction and isoplanatic imaging by partially coherent beams,” Opt. Commun. **26**, 136 (1978). [CrossRef]

*z*small object scales, such as the fine-structure of edges, increasingly contribute to the fringe pattern discerned by the detector. The promise to resolve such structures at larger

*z*represents another motivation for our work.

*δ*= 1 – Re(

*n*) is the real refractive index decrement introduced by elastic scattering, and

*β*refers to the real absorption index. Assuming the projection approximation to be valid (negligible scattering away from the ray path), the effects induced by the object onto a monchromatic, coherent wave field propagating along, say, the

*x*

_{2}direction, see Fig. 1, is described by the transmission function

*ψ*

_{z}_{=0}(

*x, y*) as where

*λ*denotes the wavelength, the integration over

*x*

_{2}is terminated at the object boundaries, and

*x,y*are coordinates transverse to the ray path along

*x*

_{2}. The 2D phase map

*ϕ*

_{z}_{=0}(

*x, y*) is interpreted as the projection of the electron density. Functions

*B*and

*ϕ*are radiographs of the object which contain some information on the object that can be interpreted. Since –

*B*(

*x,y*) +

*iϕ*

_{z}_{=0}(

*x, y*) essentially is the Radon transform of

*n*(

*x*

_{1}

*, x*

_{2}

*, x*

_{3}) a stack of such radiographs, obtained by a stepwise rotation of the object about axis

*x*

_{3}, is input for computed tomographic 3D reconstruction of

*n*(

*x*

_{1}

*, x*

_{2}

*, x*

_{3}), see for example [4

4. A. V. Bronnikov, “Theory of quantitative phase-contrast computed tomography,” J. Opt. Soc. Am. A **19**, 472 (2002). [CrossRef]

*I*generated by the phase-distorted wavefront through interference with an undistorted reference beam [5–9

9. F. Zernike, “Phase-contrast, a new method for microscopic observation of transparent objects. Part I,” Physica **9**, 686 (1942). [CrossRef]

10. C. David, B. Nöhammer, H. H. Solak, and E. Ziegler, “Differential x-ray phase contrast imaging using a shearing interferometer,” Appl. Phys. Lett. **81**, 3287 (2002). [CrossRef]

11. 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**, 258 (2006). [CrossRef]

*α*= −1

*/k∂*, the Schlieren technique (or crystal diffraction enhanced imaging) analyses the angular deviation Δ

_{x}ϕ*α*of the transmitted X-rays. The quantitative analysis of the images requires input from dynamical X-ray diffraction theory (finite Darwin width) [12

12. E. Förster, K. Goetz, and P. Zaumseil, “Double crystal diffractometry for the characterization of targets for laser fusion experiments,” Krist. Tech. **1**, 937 (1980). [CrossRef]

14. T. J. Davis, D. Gao, T. E. Gureyev, A. W. Stevenson, and W. Wilkins, “Phase-contrast imaging of weakly absorbing materials using hard X-rays,” Nature **373**, 595 (1995). [CrossRef]

15. A. Snigirev, I. Snigireva, V. Kohn, S. Kuznetsov, and I. Schelokov, “On the possibilities of xray phase contrast microimaging by coherent highenergy synchrotron radiation,” Rev. Sci. Instrum. **66** (12), 5486 (1995). [CrossRef]

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

*I*, which develops due to interference of different, directionally displaced parts of the transmitted and propagated wavefront

_{z}*ψ*(

_{z}*x, y*), orginating by plane-wave illumination of the object, is observed by a detector placed at distance

*z*behind the object. Within any transverse plane

*z*≥ 0 only a non-homogeneous gradient of the phase (

*∂*≠ const and

_{x}ϕ_{z}*∂*≠ const) generates an intensity contrast

_{y}ϕ_{z}17. D. Paganin, *Coherent X-Ray Optics*, Oxford University Press (2006). [CrossRef]

19. D. Paganin and K. A. Nugent, “Noninterferometric Phase Imaging with Partially Coherent Light,” Phys. Rev. Lett. **80**, 2586 (1998). [CrossRef]

_{⊥}denotes the two-dimensional Nabla operator in the tranverse plane. In the limit of constant absorption (

*I*

_{z}_{=0}≡ const) Eq. (3) takes the following form Relaxing the assumption of constant absorption, the contrast-transfer function (CTF) considers a linearized version of the transmission function of Eq. (2) in

*B*and

*ϕ*

_{z}_{=0}. In [20

20. L. D. Turner, B. Dhal, J. Hayes, A. Mancuso, K. Nugent, D. Paterson, R. Scholten, C. Tran, and A. Peele, “X-ray phase imaging: Demonstration of extended conditions for homogeneous objects,” Opt. Express **12**, 2960 (2004). [CrossRef] [PubMed]

22. T. E. Gureyev, D. M. Paganin, G. R. Myers, Ya. I. Nesterets, and S. W. Wilkins, “Phase-and-amplitude computer tomography,” Appl. Phys. Lett. **89**, 034102 (2006). [CrossRef]

20. L. D. Turner, B. Dhal, J. Hayes, A. Mancuso, K. Nugent, D. Paterson, R. Scholten, C. Tran, and A. Peele, “X-ray phase imaging: Demonstration of extended conditions for homogeneous objects,” Opt. Express **12**, 2960 (2004). [CrossRef] [PubMed]

21. T. E. Gureyev, Ya.I. Nesterets, D.M. Paganin, A. Pogany, and S.W. Wilkins, “Linear algorithms for phase retrieval in the Fresnel region. 2. Partially coherent illumination,” Opt. Commun. **259**, 569 (2006). [CrossRef]

22. T. E. Gureyev, D. M. Paganin, G. R. Myers, Ya. I. Nesterets, and S. W. Wilkins, “Phase-and-amplitude computer tomography,” Appl. Phys. Lett. **89**, 034102 (2006). [CrossRef]

*linear*relation between phase in the object plane and intensity a distance

*z*away from the object) appealing to a multiplicative separation of a slowly varying factor from a small factor in the transmission function was shown. In Sec. 3.3 we point out the differences of TIE- and CTF-based phase retrieval to the nonlinear method proposed in the present paper.

23. M. Op de Beeck, D. Van Dyck, and W. Coene, “Wave function reconstruction in HRTEM: the parabola method,” Ultramicroscopy **64**, 167 (1996) and refs. therein. [CrossRef]

24. P. Cloetens, Contribution to Phase Contrast Imaging, Reconstruction and Tomography with Hard Synchrotron Radiation, PhD thesis at Vrije Universiteit Brussel (1999) and references therein. [PubMed]

*B*and

*ϕ*

_{z}_{=0}by intensity measurements performed at several distances

*z*yielding remarkably accurate results with nonlinear effects taken into account iteratively. Judging by the linear approximation only, large spatial object scales are poorly retrieved by the CTF based approach because the Fourier transform of

_{i}*ϕ*

_{z}_{=0}has a sine-function coefficient that vanishes at the origin (zeros of the sine function at nonvanishing frequencies are regularized by summation over various values of

*z*). A mixed TIE and CTF based approach is also applied, see [25

25. X. Wu and H. Liu, “A new theory of phase-contrast x-ray imaging based on Wigner distributions,” Med. Phys. **31**, 2378 (2004). [CrossRef] [PubMed]

27. M. Langer, F. Peyrin, P. Cloetens, and J.-P. Guigay, “Quantitative comparison of direct phase retrieval algorithms in in-line phase tomography,” Med. Phys. **35**, 4556 (2008). [CrossRef] [PubMed]

*z*function

*g*varies in an approximately linear way in

_{z}*z*, and Eq. (4) simplifies as: [Substituting

*g*(

_{z}*x, y*) =

*g*

_{1}(

*x, y*)

*z*and

*ϕ*(

_{z}*x, y*) =

*ϕ*

_{z}_{=0}(

*x, y*) +

*ϕ*

_{1}(

*x, y*)

*z*into Eq. (4), consistently ignoring nonvanishing powers in

*z*on the right-hand side, and multiplying by

*z*yields Eq. (5).] In [4

4. A. V. Bronnikov, “Theory of quantitative phase-contrast computed tomography,” J. Opt. Soc. Am. A **19**, 472 (2002). [CrossRef]

*δ*(

*x*

_{1}

*, x*

_{2}

*, x*

_{3}) directly in terms of the data

*g*where

_{z;θ}*θ*refers to an angle at which the object is rotated about the

*x*

_{3}direction, see Fig. 1. To retrieve

*ϕ*

_{z}_{=0}radiographs need to be collected at a single distance

*z*only which is a prerequisite for the characterization of processes [27

27. M. Langer, F. Peyrin, P. Cloetens, and J.-P. Guigay, “Quantitative comparison of direct phase retrieval algorithms in in-line phase tomography,” Med. Phys. **35**, 4556 (2008). [CrossRef] [PubMed]

29. M. Boone, Y. De Witte, M. Dierick, J. Van den Bulcke, J. Vlassenbroeck, and L. Van Hoorebeke, “Practical use of the modified Bronnikov algorithm in micro-CT,” Nucl. Instrum. Methods Phys. Res. B **267**, 1182 (2009). [CrossRef]

*z*only. The edge-enhancement regime with increasing distance

*z*emerges first in the intensity contrast due to a large density of interference points introduced by the strongly curved wavefront in the vicinity of boundaries between regions within the object that exhibit a large difference in

*δ*.

*z*, simple algorithm for phase retrieval) but works systematically beyond the linear approximation of Eq. (5). We demonstrate the potential usefulness of the extended algorithm by analysing simulated phantom data.

4. A. V. Bronnikov, “Theory of quantitative phase-contrast computed tomography,” J. Opt. Soc. Am. A **19**, 472 (2002). [CrossRef]

*g*as well as

_{z}*ϕ*are expanded in powers series in

_{z}*z*. The coefficients of the power series expansion are determined by the paraxial wave equation [Eq. (11)]. We give explicit expressions up to the next-to-leading order and discuss how they can be treated perturbatively. In Sec. 4 we perform first tests of our phase-retrieval algorithm applying it to simulated phantom data. Section 5 summarizes our work. An appendix contains the relation between

*g*and

_{z}*ϕ*

_{z}_{=0}at next-to-next-to leading order.

## 2. Tomographic reconstruction algorithm

*δ*(

*x*

_{1}

*, x*

_{2}

*, x*

_{3}) by inverse Radon transformation is

*g*. Although we do not perform tomographic reconstruction in the present paper we would like to sketch conventions for the tomographic set-up as used in [4

_{z;θ}**19**, 472 (2002). [CrossRef]

**19**, 472 (2002). [CrossRef]

*δ̂*(

*s*,

*θ*,

*ω*) of the object function

*δ*(

*x*

_{1}

*, x*

_{2}

*, x*

_{3}) to the 2D Radon transform

*x*

_{1}

*, x*

_{2}and

*x*

_{3}denote the Cartesian coordinates in the object frame. The coordinates in the plane transverse to the direction of propagation are

*x*and

*y*. In this plane the 2D Radon transformation is understood as an integration over lines parameterized as

*x*sin

*ω*+

*y*cos

*ω*=

*s*(see Fig. 1). In order to apply the Radon inversion formula for the reconstruction of object function

*δ*(

*x*

_{1}

*, x*

_{2}

*, x*

_{3}) from the data

*g*(

_{z;θ}*x,y*), Eqs. (6) and (7) thus requires a one-to-one relation between

*g*. To leading (linear) order in a Taylor expansion of

_{z;θ}*g*in

_{z;θ}*z*this relation is provided by Eq. (5). Notice that Eq. (4) continues to hold approximately for weakly absorbing objects, that is, if |∇

_{⊥}

*I*

_{z}_{=0;}

*| is finite but |∇*

_{θ}_{⊥}

*I*

_{z}_{=0;}

*| ≪ |∇*

_{θ}_{⊥}

*I*

_{z}_{>0;}

*|. In this case two measurements,*

_{θ}*I*

_{z}_{=0;}

*and*

_{θ}*I*

_{z}_{>0;}

*, are required to determine*

_{θ}*g*and hence, by virtue of Eqs. (6), (5), and (7), object function

_{z;θ}*δ*(

*x*

_{1}

*, x*

_{2}

*, x*

_{3}).

## 3. Intensity contrast beyond linearity in *z*

### 3.1. Powers in z, paraxial wave equation, and perturbation theory

*z*should be increased when applying the single-distance in-line propagation technique: As discussed in the Introduction, the intensity pattern

*I*then contains discernable information on the exiting phase

_{z;θ}*ϕ*

_{z}_{=0;}

*(*

_{θ}*x,y*) and thus on the object function

*δ*(

*x*

_{1}

*, x*

_{2}

*, x*

_{3}). The model of Eq. (5), assuming a linear

*z*dependence of

*g*(

_{z;θ}*x,y*), however, breaks down at larger

*z*, and thus a more general ansatz is needed: Also, function

*ϕ*(

_{z;θ}*x,y*) is expanded in powers of

*z*: Note that {

*g*} = {

_{l;θ}*g*

_{1;θ},

*g*

_{2;θ},...} and {

*ϕ*} = {

_{m;θ}*ϕ*

_{0;θ},

*ϕ*

_{1;θ},..} denote the

*z*-independent coefficients of the expansion of the

*z*-dependent functions

*g*and

_{z;θ}*ϕ*, respectively. Substituting ansatz (8) and ansatz (9) into Eq. (4) and comparing coefficients up to linear order in

_{z;θ}*z*yields Equation (4) follows from the imaginary part of the paraxial equation for the wavefield

*I*

_{z}_{=0;θ}≡ const, the real part of Eq. (11) is given as In the limit

*z*→ 0 we have from Eq. (12) for the coefficient

*ϕ*

_{1;θ}in Eq. (9) Combining Eqs. (8), (9), (10), and (13) (truncating at

*l*

_{max}= 2 and

*m*

_{max}= 1, next-to-leading-order) finally yields The result in Eq. (14) can also be obtained along the lines of [4

**19**, 472 (2002). [CrossRef]

*I*from the wavefield

_{z}*ψ*, Fresnel-propagated from

_{z}*z*= 0, and by neglecting variations in

*x*and

*y*of the attenuation factor. The counting in powers of

*z*then is due to a Taylor expansion of the Fourier transform of the Fresnel propagator. In such a derivation of Eq. (14) it is, however, impossible to keep track of how the nonlinear corrections arise. This becomes clear when directly appealing to Eq. (10). Namely, the first two terms in the square brackets arise from the unpropagated phase while the third term is due to phase-propagation.

*l*

_{max}– 2 derivatives of Eq. (12) with respect to

*z*are needed to express the coefficients in the expansions of Eqs. (8) and (9) in terms of polynomials of transverse derivatives acting on

*ϕ*

_{0;θ}. For later discussion we notice that the nonlinear correction in Eq. (14) represents a sum of total derivatives in

*x*and

*y*. The expansion of

*ϕ*

_{0;θ}. We thus are required to demand a weak variation of

*ϕ*

_{0;θ}in dependence of

*x,y*. (Otherwise the coefficients of powers in

*z*would be unacceptably large when increasing the order of the expansion in powers of

*z*.)

*ϕ*

_{z}_{=0;θ}linking the latter to the data

*g*. Let us now discuss how the solution to Eq. (14) can be approached. In principle, one could try to solve Eq. (14) numerically subjecting it to appropriate boundary conditions. Notice that due to the presence of additional powers of ∇

_{z;θ}_{⊥}

*ϕ*

_{0;θ}in Eq. (14) as compared to Eq. (5) an ambiguity ∇

_{⊥}

*ϕ*

_{0;θ}→ ∇

_{⊥}

*ϕ*

_{0;θ}+

*v⃗*,

*v⃗*a constant 2D vector no longer is present. Also, a reduction of the influence of optical vortices [18, 19

19. D. Paganin and K. A. Nugent, “Noninterferometric Phase Imaging with Partially Coherent Light,” Phys. Rev. Lett. **80**, 2586 (1998). [CrossRef]

*ϕ*

_{0;θ}: Phase retrieval is not only constrained by the transport-of-intensity equation (imaginary part of the paraxial wave equation representing Fresnel diffraction theory incompletely) but also by the real part of the paraxial wave equation. Thus, expanding beyond leading order increasingly accounts for the constraints of the full Fresnel theory.

*l*

_{max}= 1 and

*m*

_{max}= 0 in Eqs. (8) and (9), respectively) to find a useful approximation. The nonlinear terms in square brackets are then approximated in terms of the solution

*ϕ*

_{0;θ}of Eq. (5) easily obtained by Fourier transformation (perturbation theory).

*c*

_{0}(

*x, y*) vanish in Eq. (15). For example, setting

*l*

_{max}= 2, as assumed in deriving Eq. (14), the term

*k*(

*g*

_{1;θ}(

*x,y*) +

*g*

_{2;θ}(

*x,y*)

*z*) while the terms in square brackets is by definition linear in

*z*. Ideally, this linear term cancels the linear term arising in

*l*

_{max}= 2 →

*l*

_{max}= 3 the quadratic term in

*z*since by inverting the Laplacian in a truncation of the right-hand side of Eq. (14) at the leading order one, indeed, ends up with an estimate for

*ϕ*

_{0;θ}that is independent of

*z*. Due to its perturbative nature, this estimate, when substituted in the corrections to next-to-leading and next-to-next-to-leading order of Eq. (34), does, however, not completely cancel the respective powers in

*z*in the term

### 3.2. Numerical implementation

*x,y*dimensionless as

*g*is dimensionless its

_{z;θ}*z*dependence is determined by the dimensionless variables

*kz*and

*a*,

*b*, ··· are object scales.

*ϕ*

_{0;θ}in terms of the data we appeal to Fourier transformation. To leading order we have

*g*

_{1;θ}=

*g*/

_{z}*z*. Thus in perturbation theory NLO

_{1}reads as follows: where ℱ (ℱ

^{−1}) denotes (inverse) Fourier transformation,

*ξ*

_{1},

*ξ*

_{2}are the Fourier conjugates to

*x,y*, respectively, and a summation convention for repeated indices is understood (

*i*= 1,2). NLO

_{2}simply is and NLO

_{3}is given as Let us now discuss the expressions NLO

_{1}(

*x, y*), NLO

_{2}(

*x, y*), and NLO

_{3}(

*x, y*) in more detail. Up to the factor ∇

_{⊥}

*ϕ*

_{0;θ}in NLO

_{1}, which relates the Laplacian of

*ϕ*

_{0;θ}to the data

*g*

_{1;θ}in a

*nonlocal way*, the right-hand side of Eq. (17) is local in

*g*

_{1;θ}. Loosely speaking, ∇

_{⊥}

*ϕ*

_{0;θ}is obtained by integrating the leading-order equation

_{2}is completely local in

*g*

_{1;θ}while NLO

_{3}is nonlocal. Each of the terms NLO

_{1}(

*x, y*), NLO

_{2}(

*x, y*), and NLO

_{3}(

*x, y*), can be decomposed as where 〈NLO

*〉 ≡ ∫*

_{i}*dxdy*NLO

*(*

_{i}*x, y*). Since NLO

_{1}(

*x, y*) + NLO

_{2}(

*x, y*) and NLO

_{3}(

*x, y*) separately are sums of partial derivatives we have Since NLO

_{3}is

*nonlocal*in

*g*, this term

_{z;θ}*smeares g*

_{1;θ}and averages out large amplitudes.

*α*is a real constant. Requiring that the retrieved phase is insensitive to a variation of

*α*fixes

*α*’s value. It is straightforward to show that using the regularization of Eq. (21) in applying

*F*(

*x,y*) one has where

_{1}〉 and 〈NLO

_{2}〉 individually can be large but their sum is nil.). Deviations from this situation are introduced by numerical artifacts introduced through discrete Fourier transformations and thus should be subtracted. The exact phase, which is expressed by a line integral over the object function

*δ*(

*x*

_{1}

*, x*

_{2}

*, x*

_{3}), however has a nonvanishing mean (

*δ*(

*x*

_{1}

*, x*

_{2}

*, x*

_{3}) is sign definite). This ambiguity is due to the fact that an arbitrary constant can be added to without changing Eq. (14) for

*ϕ*

_{0;θ}. To compare exact with retrieved phase we thus subtract means in both cases.

### 3.3. Comparison with CTF-based single-distance retrieval

*ψ*

_{0},

*to the Fourier transform of intensity ℱ*

_{θ}*I*

_{z,θ}at distance

*z*: where

*r⃗*≡ (

*x,y*). Assuming a pure-phase object, one has

*I*

_{0}

*,*is homogeneous. In exploiting Eq. (25) one makes the assumption that

_{θ}*r⃗*and for all

*ξ⃗*at fixed

*z*and fixed

*k*. If this assumption is met then one may write Substituting Eq. (26) into Eq. (25), one has [30] Eq. (27) provides for an

*algebraic*solution of the phase-retrieval problem in Fourier space provided the zeros of the sine function are regularized in a natural way. Notice that going beyond Eq. (26) in taking higher powers of

*ϕ*

_{0;θ}into account algebraic inversion in the sense of Eq. (27) no longer is possible. This is because in Fourier space quadratic and higher orders enter in Eq. (25) in a highly nonlocal way.

^{−1}act on Eq. (27), that Therefore, by assuming the validity of Eq. (26) the inversion of Eq. (27) corresponds to an infinite summation of powers of

*ϕ*

_{0;θ}(

*r⃗*). We stress that our expansion (8) contains these powers of

*z*× Eq. (14) at order

*z*and the term

*z*× Eq. (34) at order

*z*

^{3}. Notice that in addition to these terms nonlinear corrections in

*ϕ*

_{0;θ}enter at these orders. This is the reason why the above assumption can increasingly be relaxed taking increasing powers of

*z*into account in Eq. (8). The series Eq. (8) expands in powers of transverse derivatives

*including*nonlinearities in

*ϕ*

_{0;θ}whereas expansion (29) takes into account all powers of derivatives at

*linear order*in

*ϕ*

_{0;θ}only. For weak phase objects expansion (29) is an excellent approximation in treating the case where all frequencies contribute to the overall phase shift such that the latter is smaller than unity but breaks down as soon as a part of the spectrum contributes relative phase shifts comparable or larger than unity. Since the nonlinear expansion of Eq. (8) for pragmatic reasons, needs to be truncated at a

*finite*order in

*z*large frequencies are retrieved in a worse way than they are by inverting Eq. (27) within that equation’s regime of validity.

## 4. Phase retrieval for simulated phantom data

*ϕ*

_{0;θ}.

*α*we define a function Φ(

*θ*,

*α*) as follows where

*ϕ*

_{0;θ}(

*x,y*) is the retrieved phase at regularization parameter value

*α*. The value of

*α*is chosen such that Φ(

*θ*,

*α*) is least sensitive to a variation in

*α*. For the phantoms considered below there occurs a wide region (several orders of magnitude with central value

*α*∼ 10

_{c}^{−6}) of insensitivity for Φ(

*θ*,

*α*).

_{θ,αc}as follows: where

*ϕ*

_{0;θ}either is the leading-order (LO) or includes the next-to-leading-order (LO + NLO) result. The quantity 〈Ψ

*〉 is defined as the mean of Ψ*

_{θ}*over all pixels.*

_{θ}*θ*= 0) at an X-ray energy of

*E*= 30keV, at

*δ*= 10

^{−7}on a

*δ*= 0 background, and for thickness

*d*= 0.256mm. A central disk with constant phase shift was introduced to avoid numerical problems arising from unresolved spokes segments. To avoid extreme phase jumps we have applied a normalized Gaussian blurring to the Siemens star to yield an exact phase map from which intensity is computed by free-space propagation. As Figs. 2(e) and 2(f) indicates, function Ψ

*is smaller when taking next-to-leading corrections into account. To investigate the dependence on distance and resolution of the retrieved phase we perform line cuts through the phase map as indicated in Fig. 2(a). The results are shown in Fig. 3. Keeping the standard deviation of the Gaussian blurring constant at 1.5 pixels, results do not change substantially for smaller values of*

_{θ}*δ*than 10

^{−7}. (This means that at a higher resolution the averaged-over length scale becomes shorter.) For

*δ*> 10

^{−6}, however, strong deviations of retrieved from exact phase maps occur both in the leading-order and the next-to-leading-order result suggesting that a violation of our assumption on weak phase variation takes place. Recall that strong phase variations yield large coefficients in the power series of Eq. (15) and thus worsen the convergence properties of the expansion.

## 5. Summary and outlook

*z*only.

*z*. In Fresnel theory the expansion coefficients are determined in terms of increasing numbers of transverse derivatives of the phase map in the object plane. The resulting partial differential equation starts to be nonlinear at next-to-leading order and can be approached in a perturbative way. Using simulated phantom data, we have demonstrated the promise of this method.

*δ*) the requirement of small values of derivatives (convergence of the expansion in Eq. (8)) is better met in this case (less nervousness on small transverse scales). Therefore, we expect the present approach to have applications in the tomography of samples introducing large phase shifts.

*z*(Eq. (8)) yield even better convergence properties of the expansion.

## Appendix: Next-to-next-to leading corrections to
∇ ⊥ 2 ϕ 0 ; θ

*l*

_{max}= 3 and

*m*

_{max}= 2, respectively. Inserting this truncation into Eq. (4) and comparing coefficients up to quadratic order in

*z*yields

*z*we obtain for

*z*→ 0: Including the next-to-next-to-leading order, we thus have

## Acknowledgments

24. P. Cloetens, Contribution to Phase Contrast Imaging, Reconstruction and Tomography with Hard Synchrotron Radiation, PhD thesis at Vrije Universiteit Brussel (1999) and references therein. [PubMed]

## References and links

1. | A. Papoulis, “Ambiguity function in Fourier optics,” J. Opt. Soc. Am. |

2. | J.-P. Guigay, “Fourier transform analysis of Fresnel diffraction patterns and in-line holograms,” Optik |

3. | J.-P. Guigay, “The ambiguity function in diffraction and isoplanatic imaging by partially coherent beams,” Opt. Commun. |

4. | A. V. Bronnikov, “Theory of quantitative phase-contrast computed tomography,” J. Opt. Soc. Am. A |

5. | M. Born and E. Wolf, |

6. | A. Momose, “Demonstration of phase-contrast X-ray computed tomography using an X-ray interferometer,” Nucl. Instrum. Methods Phys. Res. A |

7. | F. Beckmann, U. Bonse, F. Busch, O. Günnewig, and T. Biermann, “A novel system for X-ray phase-contrast microtomography,” HASYLAB Annual Report |

8. | U. Bonse and M. Hart, “An X-ray interferometer,” Appl. Phys. Lett. |

9. | F. Zernike, “Phase-contrast, a new method for microscopic observation of transparent objects. Part I,” Physica |

10. | C. David, B. Nöhammer, H. H. Solak, and E. Ziegler, “Differential x-ray phase contrast imaging using a shearing interferometer,” Appl. Phys. Lett. |

11. | F. Pfeiffer, T. Weitkamp, O. Bunk, and C. David, “Phase retrieval and differential phase-contrast imaging with low-brilliance X-ray sources,” Nat. Phys. |

12. | E. Förster, K. Goetz, and P. Zaumseil, “Double crystal diffractometry for the characterization of targets for laser fusion experiments,” Krist. Tech. |

13. | V. A. Bushuev, V. N. Ingal, and E. A. Belyaevskaya, “Dynamical Theory of Images Generated by Noncrystalline Objects for the Method of Phase-Dispersive Introscopy,” Kristallografiya |

14. | T. J. Davis, D. Gao, T. E. Gureyev, A. W. Stevenson, and W. Wilkins, “Phase-contrast imaging of weakly absorbing materials using hard X-rays,” Nature |

15. | A. Snigirev, I. Snigireva, V. Kohn, S. Kuznetsov, and I. Schelokov, “On the possibilities of xray phase contrast microimaging by coherent highenergy synchrotron radiation,” Rev. Sci. Instrum. |

16. | S. W. Wilkins, T. E. Gureyev, D. Gao, A. Pogany, and A. W. Stevenson, “Phase-contrast imaging using polychromatic hard X-rays,” Nature |

17. | D. Paganin, |

18. | T. E. Gureyev, A. Roberts, and K. A. Nugent, “Phase retrieval with the transport-of-intensity equation: matrix solution with use of Zernike polynomials,” J. Opt. Soc. Am. A |

19. | D. Paganin and K. A. Nugent, “Noninterferometric Phase Imaging with Partially Coherent Light,” Phys. Rev. Lett. |

20. | L. D. Turner, B. Dhal, J. Hayes, A. Mancuso, K. Nugent, D. Paterson, R. Scholten, C. Tran, and A. Peele, “X-ray phase imaging: Demonstration of extended conditions for homogeneous objects,” Opt. Express |

21. | T. E. Gureyev, Ya.I. Nesterets, D.M. Paganin, A. Pogany, and S.W. Wilkins, “Linear algorithms for phase retrieval in the Fresnel region. 2. Partially coherent illumination,” Opt. Commun. |

22. | T. E. Gureyev, D. M. Paganin, G. R. Myers, Ya. I. Nesterets, and S. W. Wilkins, “Phase-and-amplitude computer tomography,” Appl. Phys. Lett. |

23. | M. Op de Beeck, D. Van Dyck, and W. Coene, “Wave function reconstruction in HRTEM: the parabola method,” Ultramicroscopy |

24. | P. Cloetens, Contribution to Phase Contrast Imaging, Reconstruction and Tomography with Hard Synchrotron Radiation, PhD thesis at Vrije Universiteit Brussel (1999) and references therein. [PubMed] |

25. | X. Wu and H. Liu, “A new theory of phase-contrast x-ray imaging based on Wigner distributions,” Med. Phys. |

26. | T. E. Gureyev, A. Pogany, D. M. Paganin, and S. W. Wilkins, “Linear algorithms for phase retrieval in the Fresnel region,” Opt. Commun. |

27. | M. Langer, F. Peyrin, P. Cloetens, and J.-P. Guigay, “Quantitative comparison of direct phase retrieval algorithms in in-line phase tomography,” Med. Phys. |

28. | A. Groso, R. Abela, and M. Stampanoni, “Implementation of a fast method for high resolution phase contrast tomography,” Opt. Express |

29. | M. Boone, Y. De Witte, M. Dierick, J. Van den Bulcke, J. Vlassenbroeck, and L. Van Hoorebeke, “Practical use of the modified Bronnikov algorithm in micro-CT,” Nucl. Instrum. Methods Phys. Res. B |

30. | J.-P. Guigay, R. H. Wade, and C. Delpha, “Optical diffraction of Lorentz microscope images,” Proceedings of the 25th meeting of the Electron Microscopy and Analysis Group, W.C. Nixon ed. (The Institute of Physics, London, 1971), pp. 238–239. |

31. | M. Langer, |

**OCIS Codes**

(100.2960) Image processing : Image analysis

(100.3010) Image processing : Image reconstruction techniques

(100.5070) Image processing : Phase retrieval

(340.7440) X-ray optics : X-ray imaging

**ToC Category:**

Image Processing

**History**

Original Manuscript: September 16, 2010

Revised Manuscript: October 11, 2010

Manuscript Accepted: October 14, 2010

Published: November 23, 2010

**Citation**

Julian Moosmann, Ralf Hofmann, Andrei Bronnikov, and Tilo Baumbach, "Nonlinear phase retrieval from single-distance radiograph," Opt. Express **18**, 25771-25785 (2010)

http://www.opticsinfobase.org/oe/abstract.cfm?URI=oe-18-25-25771

Sort: Year | Journal | Reset

### References

- A. Papoulis, "Ambiguity function in Fourier optics," J. Opt. Soc. Am. 64, 779 (1974). [CrossRef]
- J.-P. Guigay, "Fourier transform analysis of Fresnel diffraction patterns and in-line holograms," Optik 49, 121 (1977).
- J.-P. Guigay, "The ambiguity function in diffraction and isoplanatic imaging by partially coherent beams," Opt. Commun. 26, 136 (1978). [CrossRef]
- A. V. Bronnikov, "Theory of quantitative phase-contrast computed tomography," J. Opt. Soc. Am. A 19, 472 (2002). [CrossRef]
- M. Born, and E. Wolf, Principles of Optics, 6th ed., Pergamon Press, Oxford, New York (1980).
- A. Momose, "Demonstration of phase-contrast X-ray computed tomography using an X-ray interferometer," Nucl. Instrum. Methods Phys. Res. A 352, 622 (1995). [CrossRef]
- F. Beckmann, U. Bonse, F. Busch, O. Günnewig, and T. Biermann, "A novel system for X-ray phase-contrast microtomography," HASYLAB Annual Report II, 691 (1995).
- U. Bonse, and M. Hart, "An X-ray interferometer," Appl. Phys. Lett. 6, 155 (1965). [CrossRef]
- F. Zernike, "Phase-contrast, a new method for microscopic observation of transparent objects. Part I," Physica 9, 686 (1942). [CrossRef]
- C. David, B. Nöhammer, H. H. Solak, and E. Ziegler, "Differential x-ray phase contrast imaging using a shearing interferometer," Appl. Phys. Lett. 81, 3287 (2002). [CrossRef]
- 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, 258 (2006). [CrossRef]
- E. Förster, K. Goetz, and P. Zaumseil, "Double crystal diffractometry for the characterization of targets for laser fusion experiments," Krist. Tech. 1, 937 (1980). [CrossRef]
- V. A. Bushuev, V. N. Ingal, and E. A. Belyaevskaya, "Dynamical Theory of Images Generated by Noncrystalline Objects for the Method of Phase-Dispersive Introscopy," Kristallografiya 41, 808 (1996).
- T. J. Davis, D. Gao, T. E. Gureyev, A. W. Stevenson, and W. Wilkins, "Phase-contrast imaging of weakly absorbing materials using hard X-rays," Nature 373, 595 (1995). [CrossRef]
- A. Snigirev, I. Snigireva, V. Kohn, S. Kuznetsov, and I. Schelokov, "On the possibilities of xray phase contrast microimaging by coherent high energy synchrotron radiation," Rev. Sci. Instrum. 66(12), 5486 (1995). [CrossRef]
- S. W. Wilkins, T. E. Gureyev, D. Gao, A. Pogany, and A. W. Stevenson, "Phase-contrast imaging using polychromatic hard X-rays," Nature 384, 335 (1996). [CrossRef]
- D. Paganin, Coherent X-Ray Optics, Oxford University Press (2006). [CrossRef]
- T. E. Gureyev, A. Roberts, and K. A. Nugent, "Phase retrieval with the transport-of-intensity equation: matrix solution with use of Zernike polynomials," J. Opt. Soc. Am. A 12, 1942 (1995).
- D. Paganin, and K. A. Nugent, "Noninterferometric Phase Imaging with Partially Coherent Light," Phys. Rev. Lett. 80, 2586 (1998). [CrossRef]
- L. D. Turner, B. Dhal, J. Hayes, A. Mancuso, K. Nugent, D. Paterson, R. Scholten, C. Tran, and A. Peele, "X-ray phase imaging: Demonstration of extended conditions for homogeneous objects," Opt. Express 12, 2960 (2004). [CrossRef] [PubMed]
- T. E. Gureyev, Ya. I. Nesterets, D. M. Paganin, A. Pogany, and S. W. Wilkins, "Linear algorithms for phase retrieval in the Fresnel region. 2. Partially coherent illumination," Opt. Commun. 259, 569 (2006). [CrossRef]
- T. E. Gureyev, D. M. Paganin, G. R. Myers, Ya. I. Nesterets, and S. W. Wilkins, "Phase-and-amplitude computer tomography," Appl. Phys. Lett. 89, 034102 (2006). [CrossRef]
- M. Op de Beeck, D. Van Dyck, and W. Coene, "Wave function reconstruction in HRTEM: the parabola method," Ultramicroscopy 64, 167 (1996) and refs. therein. [CrossRef]
- P. Cloetens, Contribution to Phase Contrast Imaging, Reconstruction and Tomography with Hard Synchrotron Radiation, PhD thesis at Vrije Universiteit Brussel (1999) and references therein. [PubMed]
- X. Wu, and H. Liu, "A new theory of phase-contrast x-ray imaging based on Wigner distributions," Med. Phys. 31, 2378 (2004). [CrossRef] [PubMed]
- T. E. Gureyev, A. Pogany, D. M. Paganin, and S. W. Wilkins, "Linear algorithms for phase retrieval in the Fresnel region," Opt. Commun. 231, 53 (2004). [CrossRef]
- M. Langer, F. Peyrin, P. Cloetens, and J.-P. Guigay, "Quantitative comparison of direct phase retrieval algorithms in in-line phase tomography," Med. Phys. 35, 4556 (2008). [CrossRef] [PubMed]
- A. Groso, R. Abela, and M. Stampanoni, "Implementation of a fast method for high resolution phase contrast tomography," Opt. Express 14, 8103 (2006). [CrossRef] [PubMed]
- M. Boone, Y. De Witte, M. Dierick, J. Van den Bulcke, J. Vlassenbroeck, and L. Van Hoorebeke, "Practical use of the modified Bronnikov algorithm in micro-CT," Nucl. Instrum. Methods Phys. Res. B 267, 1182 (2009). [CrossRef]
- J.-P. Guigay, R. H. Wade, and C. Delpha, "Optical diffraction of Lorentz microscope images," Proceedings of the 25th meeting of the Electron Microscopy and Analysis Group, W. C. Nixon ed. (The Institute of Physics, London, 1971), pp. 238-239.
- M. Langer, Phase Retrieval in the Fresnel Region for Hard X-ray Tomography, PhD thesis at L’Insitut National des Sciences Appliquees de Lyon (2008) and references therein.

## 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.