## Abstract

Here a transport model is used to simulate amplitude-only imaging and intensity-based quantitative phase imaging in a turbid medium. We derive an optical transfer function for propagation through a scattering medium. We also show that, as expected, scattering leads to a degradation in the spatial resolution in both forms of imaging, while the magnitude of the phase retrieved using a solution of the transport-of-intensity equation decreases as the optical density of the scattering medium increases.

© 2008 Optical Society of America

## 1. Introduction

Understanding the propagation of light through scattering media is important in a range of contexts, most notably biomedical imaging. It is well-known that images of objects embedded in turbid surroundings are degraded and a variety of strategies and imaging modalities have been developed to improve the signal-to-noise ratio and enhance contrast. Central to the development of new imaging techniques is understanding the interaction between light and the scattering medium and how this impacts on the propagation of information about the specimen under investigation. The most widely used approach to studying propagation of light through random media involves solutions of the radiative transfer equation [1–3]. Methods used to solve this equation involve either stochastic, the Monte-Carlo method [4] for example, or deterministic approaches involving approximations to the radiative transfer equation [5].

We approach modeling the image formation process through the study of the coherence properties of the field transmitted through the medium. The coherence properties of a wavefield are embodied in the correlation functions. It is well known that the spatial coherence properties of a light field evolve as it propagates through a scattering medium such as biological tissue [6, 7]. The decrease in the degree of coherence has been shown to influence the quality of biological imaging methods [8, 9].

Both the Generalized Radiance (GR) and the Spatial Coherence Function (SCF) are sensitive to both phase and amplitude changes in a multiply scattered field [10]. Knowledge of these changes in an optical field is critical to a complete understanding of techniques such as microscopy [11], optical lithography [12] and the development of new imaging methods in the optical [13], X-ray [14] and other regions of the spectrum.

Most studies into imaging through scattering media have concentrated on measuring intensity distributions [15]. Many samples of interest, however, exhibit poor absorption contrast but possess significant refractive index or thickness variations. A number of phase imaging modalities [16] have consequently been developed for imaging such objects. One of these, Quantitative Phase Microscopy (QPM) [17] has been shown to be able to reproduce phase images of transparent specimens in the absence of scattering by random inhomogeneities in the surrounding medium and the technique can be implemented using conventional bright field optics. Here we investigate the influence of scattering on both conventional amplitude imaging and quantitative phase imaging utilizing the transport-of-intensity equation. The approach we adopt involves a transport-based method for finding the generalized radiance of initially coherent light after it has propagated through a scattering medium. Using this we derive an optical transfer function that describes the impact of scattering on the spatial frequency content of the transmitted field. Implications for imaging amplitude objects as well as utilizing QPM in a turbid medium are discussed.

## 2. The correlation functions

Quasi-monochromatic, scalar, partially coherent, paraxial wavefields are completely characterized by various correlation functions. One of these, the mutual optical intensity, *G*(*r*⃑_{1},*r*⃑_{2}), is defined [18] as:

where the brackets 〈〉 denotes an average over all realizations of the field and *r*⃑_{1} and *r*⃑_{2} are the coordinates of two points in the field in the plane z. Another correlation function called the spatial coherence function, *J*(*x*⃑, *s*⃑, *z*⃑), is defined in terms of the mean, *x*=⃑(*r*⃑_{1}+*r*⃑_{2})/2 and the difference, *s*⃑=*r*⃑_{1}-*r*⃑_{2}, coordinates [18] as:

Equivalently, the field can be described using a phase space representation called the generalized radiance, Wigner distribution function or phase-space density [19] which is written as:

where *u*⃗ is a variable conjugate to the difference coordinate *s*⃑. Finally, the ambiguity function [20] is written in the form of a four-dimensional Fourier transform of the generalized radiance [21]:

It is straightforward to show [20] that in the paraxial approximation in a homogeneous medium the ambiguity function after propagation through a distance *z* can be written:

where *k* is the wavenumber in the medium. Thus the ambiguity function is sheared in a four-dimensional space as the field propagates.

## 3. Model

We consider the situation shown in Fig. 1. A field (coherent or partially coherent) described by a spatial coherence function, *J*
_{0}(*x*⃗, *s*⃗) (and ambiguity function X_{0}(*q*⃗,*s*⃗) is incident upon a thin object with a complex transmission coefficient *O*(*x*⃗)=*A*(*x*⃗)exp(*iφ*(*x*⃗)). The light then passes through a turbid medium of thickness *ℓ*. In order to study the influence of scattering on image formation we wish to determine the intensity of the field leaving the scattering medium (i.e. at *z*=0) and, in order to investigate its impact on phase imaging utilizing the transport-of-intensity equation, we also wish to find the derivative of the intensity at this plane with respect to *z*.

The formulation used here to model the propagation of light through turbid media is a Boltzmann-like transport equation proposed by Cheng and Raymer [6]. Under the assumption that most of the scattered light propagates in the forward direction and the correlation between near-forward and wide-angle scattered light is weak, it is possible to write down a differential equation for the near-forward scattered component of the spatial coherence function *J _{S}*(

*x*⃗,

*s*⃗,0) and it has been shown [6, 22] that this equation has a solution of the form:

where

and

$$\times \left\{\mathrm{erf}\left(\frac{{k}_{\mathrm{med}}{\theta}_{0}}{2\mid \overrightarrow{q}\mid}\overrightarrow{q}\u2022\overrightarrow{s}\right)-\mathrm{erf}\left[\frac{{k}_{\mathrm{med}}{\theta}_{0}}{2}\left(\frac{\overrightarrow{q}\u2022\overrightarrow{s}}{\mid \overrightarrow{q}\mid}-\frac{\ell \mid \overrightarrow{q}\mid}{{k}_{\mathrm{med}}}\right)\right]\right\}$$

where *θ _{0}* is the width of a Gaussian fitted to the Mie differential scattering cross-section for the given scattering particles [23],

*σ*is the scattering cross section integrated over the angle defined by the collection angle of the optical system used,

_{N}*N*is the number density of scattering particles,

*k*is the wavenumber in the medium and

_{med}*µ*is the total extinction coefficient.

_{T}Taking the Fourier transform of Eq. (6) over the spatial coordinate, *x*⃗, leads to the ambiguity function [20] (Eq. (4)):

where *X _{0}* is the ambiguity function of the input field. Hence the term

*S*(

*q*⃗,

*s*⃗) acts as an optical transfer function describing the influence of only scattering on the propagation. This observation forms the key conclusion of the present paper. The shearing of the ambiguity function due to propagation in a homogeneous medium is described by the second factor. In the remainder of the paper, we investigate implications for imaging both amplitude and phase objects.

The spatial coherence function of light transmitted through a thin object with complex transmission function, *O*(*x*)=*A*(*x*)exp(-*iφ*(*x*)), is written

where *J*
_{0}(*x*⃗,*s*⃗), is the spatial coherence function of the field incident on the object.

Here we focus on the case where the length of the scattering region is short compared with typical distances over which diffraction effects become apparent or that the object is perfectly re-imaged using an optical system of unity magnification. In this case, the second, shearing, term in Eq. (9) can be ignored and assuming uniform, coherent illumination, i.e. *J*
_{0}(*x*⃗, *s*⃗)=constant, and weak phase gradients, the ambiguity function describing the transmitted light takes the form:

since *J*
_{obj}(*x*⃗,*s*⃗)≈*A*
^{2}(*x*⃗)exp(*is*⃗•*∇*φ(*x*⃗)). Using Eq. (5) the ambiguity function a distance *z* from the output face of the scattering medium is given by:

$$\times \int \mathrm{dx}{A}^{2}(\overrightarrow{x})\mathrm{exp}\left(i\left(\overrightarrow{s}-\frac{z}{k}\overrightarrow{q}\right)\u2022\nabla \phi (\overrightarrow{x})\right)\mathrm{exp}(-i\overrightarrow{q}\u2022\overrightarrow{x}).$$

The intensity a distance *z* from the output face of the scattering medium is given by:

$$\times \mathrm{exp}\left(-i\frac{z}{k}\overrightarrow{q}\u2022\nabla \phi (\overrightarrow{x}\text{'})\right)\mathrm{exp}(i\overrightarrow{q}\u2022(\overrightarrow{x}-\overrightarrow{x}\text{'}))$$

and at the *z*=0 plane (i.e. at the output face) by:

$$=\frac{1}{{(2\pi )}^{2}}\int d\overrightarrow{q}S(\overrightarrow{q},0){\tilde{I}}_{0}(\overrightarrow{q})\mathrm{exp}(i\overrightarrow{q}\u2022\overrightarrow{x}),$$

where

is the Fourier transform of the intensity transmission coefficient of the object. Hence, in the case of an amplitude-only object, it is apparent that the presence of scattering leads to the well-known reduction in the spatial resolution of the system and the term *S*(*q*⃗,0) acts as the amplitude transfer function.

Similarly, the intensity derivative:

$$+\frac{1}{{(2\pi )}^{2}k}\int \overrightarrow{d}q(i\overrightarrow{q}\u2022\overrightarrow{x}\text{'})\overrightarrow{q}\u2022{\nabla}_{s}S(\overrightarrow{q},0){\tilde{I}}_{0}(\overrightarrow{q})$$

where ∇* _{s}* denotes the gradient with respect to the second coordinate,

*I*̃

_{0}(

*q*⃗) is defined in Eq. (16) and

*G*̃(

*q*⃗)=

*iq*⃗•∫

*dx*⃗

*x*′

*A*

^{2}(

*x*⃗′)∇φ(

*x*⃗′)exp(-

*iq*⃗•

*x*⃗′) is the Fourier transform of ∇•(

*A*

^{2}(

*x*⃗′)∇

*φ*(

*x*⃗′)). In the case of a phase-only object and uniform, coherent illumination,

*I*̃

_{0}(

*q*⃗)=

*I*δ(

_{0}*q*⃗) and the second term in Eq. (17) is zero so

In the absence of scattering, *I _{1}*(

*q*⃗,0)=0,

*S*(

*q*⃗,0)=1 and Eq. (20) reduces to:

which is the well-known transport-of-intensity equation [24,25].

Note that the expression for the intensity (Eq. (15)) and the axial intensity gradient (Eq. (18)) are written in terms of *S*(*q*⃗,0) which is given by:

Hence, the presence of scattering will firstly reduce the spatial resolution obtainable with phase retrieval methods based on the transport-of-intensity equation and secondly, reduce the magnitude of the apparent phase of the object when such methods are used to quantitatively determine the phase excursion. Note that this will occur even in the absence of noise or other imaging imperfections.

## 4. Transfer functions

Given the interpretation of *S*(*q*⃗,*s*⃗) as being an optical transfer-like function describing scattering, it is instructive to plot this function. In all simulations, the refractive index of the spherical particles relative to the medium has been chosen to be 1.2 and the free-space wavelength 632.8 nm. The near-forward scattering cross-section, *σN* is also assumed to be equal to the total scattering cross-section. A typical transfer function is shown in Fig. 2. In this case the concentration of spheres is 1.6×10^{11}/m^{3} which corresponds to a mean free path of 10 mm, i.e. equal to the thickness of the medium. Hence, this medium has an optical density, OD=*NσTℓ*, of 1. Figure 3 shows *S*(*q*⃗,0) for media containing spherical particles of diameter (a) 5 µm and (b) 20µm. As expected, the transfer functions have a maximum at *q*=0 and can be regarded as consisting of two components associated with the ballistic and scattered components of the field. The scattered component contributes a component that is peaked around *q*=0. The width of this peak depends on the size of the scattering spheres: larger spheres produce a broader peak consistent with the narrower differential scattering cross-section while smaller spheres produce a narrower scattered contribution to *S*(*q*⃗,0). It can be seen that the ballistic contribution decreases as the optical density increases which would correspond to a decrease in spatial resolution for images containing spatial frequencies greater than approximately 2 mm^{-1} for 5 µm spheres and 8 mm^{-1} for 20 µm spheres. These values correspond to the 1/*e*-half-widths of Gaussians fitted to OD=5 data.

## 5. Simulated images

In this section imaging of a one-dimensional transmission object sandwiched between two slices of scattering medium is discussed. Calculations were performed using the tools described in the previous sections in Eqns (8–10).

A spatially coherent Gaussian input field with an intensity of half-width a was assumed. Such a field has *a* spatial coherence function given by

An amplitude-only and a phase-only object were considered. The objects in both cases were considered to be sandwiched between two layers of scattering media of identical thickness and scatterer concentration.

As a test amplitude object, a one-dimensional double-slit with slit width of 200 µm and center separation 300 µm was used (Fig. 4). The light incident on the medium had a uniform intensity profile and was spatially coherent. The turbid medium consisted of spheres of diameter 5 µm (Fig. 4(a)) and 20 µm (Fig. 4(b)) of relative refractive index 1.2 and the object was sandwiched between two layers of thickness 5 mm. The wavelength was chosen to be 632.8 nm. The intensity of the output field is shown in Fig. 4. In both cases, effects due to diffraction have been suppressed producing an in-focus image of the structure. The influence of scattering on spatial resolution is apparent, but can be seen to be greater in the case of scattering by smaller spheres. Again, the broader angular spread of the differential scattering cross-section in this case is responsible.

To consider the influence of scattering on non-interferometric phase retrieval, a non-absorbing phase object with a parabolic phase profile of width 200 µm and phase excursion *π* radians was placed in the centre of a scattering medium of total thickness 10 mm. The light incident on the medium had a field width of 2 mm and was spatially coherent. The intensity arising from free-space propagation of the resulting coherence function through distances of Δ*z*=±1 mm on either side of an in-focus image of the sample was then calculated. These images were then used with an algorithm [26] based on the transport-of-intensity equation to determine the phase. The resulting calculated phase is shown in Fig. 5 for scattering spheres of diameter 5 µm (Fig. 5(a)) and 20 µm (Fig. 5(b)). It can be seen that accompanying a decrease in the value of the recovered phase from that obtainable in the absence of scattering, there is a decrease in spatial resolution manifested as a reduction in sharpness of the edges of the object. It can also be seen that phase imaging through scattering media containing smaller spheres is much more vulnerable to both loss of spatial resolution and the accuracy of the measured phase excursion. Hence, although phase images can be obtained through scattering media, the recovered phase will not be quantitative. Note that in the presence of noise, these recovered phase images would be expected to deteriorate even further [27], particularly with the decrease in intensity accompanying scattering. Also note that the transport model used here to model changes in coherence implicitly assumes that sufficient random movement of scatterers occurs during typical camera exposure times so that adequate ensemble averaging occurs. The validity of this assumption will depend on the precise experimental situation.

## 6. Conclusion

A basic theory to explain the imaging properties in terms of the ambiguity function and an optical transfer function-like scattering function that describes the impact of scattering on image formation has been presented. This theory, based on a wave-transport model, was used to simulate simple amplitude and non-interferometric phase imaging through scattering media containing different concentrations of spheres with two different diameters. The theory and the discussion of the variation of the coherence properties of the field due to the presence of transmission objects suggests that, as expected, the degree of spatial coherence affects the quality of images obtained. In particular, both spatial resolution and the quantitative nature of phase imaging, were more adversely affected by scattering by smaller particles. This model provides a foundation on which the quantitative investigation of imaging through scattering media can be studied in more detail. Here we have focussed on conventional amplitude imaging as well as non-interferometric phase imaging, but there are clearly a large number of other specific imaging configurations to which this approach could also be applied.

## Acknowledgments

This research was supported under the Australian Research Council’s Discovery funding scheme (project number DP0557505).

## References and links

**1. **S. R. Arridge and J. C. Hebden, “Optical imaging in medicine: II. Modelling and reconstruction,” Phys. Med. Biol. **42**, 841–853 (1997). [CrossRef] [PubMed]

**2. **R. Chandrasekhar, *Radiative Transfer*, (Oxford, 1950).

**3. **A. Ishimaru, *Wave propagation and scattering in random media Volume 1: Single scattering and transport theory* (New York: Academic1978).

**4. **J. J. Duderstadt and L. J. Hamilton, *Nuclear Reactor Analysis* (Wiley, 1976).

**5. **H. W. Lewis, “Multiple scattering in an infinite medium,” Phys. Rev. **78**, 526–529 (1950). [CrossRef]

**6. **C.-C. Cheng and M. G. Raymer, “Propagation of transverse optical coherence in random multiple scattering media,” Phys. Rev. A **62**, 023811 (2000). [CrossRef]

**7. **A. Wax and J. E. Thomas, “Measurement of smoothed Wigner phase-space distribution for small-angle scattering in a turbid medium,” J. Opt. Soc. Am. A **15**, 1896–1908 (1998). [CrossRef]

**8. **F. Dubois, L. Joannes, and J.-C. Legros, “Improved three-dimensional imaging with a digital holography microscope with a source of partial spatial coherence,” Appl. Opt. **38**, 7085–7094 (1999). [CrossRef]

**9. **F. Dubois, M.-L. N. Requena, C. Minetti, O. Monnom, and E. Istasse, “Partial spatial coherence effects in digital holographic microscopy with a laser source,” Appl. Opt. **43**, 1131–1139 (2004). [CrossRef] [PubMed]

**10. **C.-C. Cheng and M. G. Raymer, “Long range saturation of spatial decoherence in wave-field transport in random multiple scattering media,” Phys. Rev. Lett. **82**, 4807–4810 (1999). [CrossRef]

**11. **F. Dubois, M.-L. N. Requena, and C. Minetti, “Partial spatial coherence effects in digital holographic microscopy with a laser source,” Appl. Opt. **43**, 1131–1139 (2004). [CrossRef] [PubMed]

**12. **N. A. Beaudry and T. D. Milster, “Effects of object roughness on partially coherent image formation,” Opt. Lett. **25**, 454–456 (2000). [CrossRef]

**13. **D. M. Marks, R. A. Stack, and D. J. Brady, “Astigmatic coherence sensor for digital imaging,” Opt. Lett. **25**, 1726–1728 (2000). [CrossRef]

**14. **A. Momose, “Phase-sensitive imaging and phase tomography using X-ray interferometers,” Opt. Express **11**, 2303–2314 (2003). [CrossRef] [PubMed]

**15. **M. Alrubaiee, M. Xu, S. K. Gayen, M. Brito, and R. R. Alfano, “Three-dimensional optical tomographic imaging of scattering objects in tissue-simulating turbid media using independent component analysis,” Appl. Phys. Lett. **87**, 19112 (2005). [CrossRef]

**16. **E. D. Barone-Nugent, A. Barty, and K. A. Nugent, “Quantitative phase-amplitude microscopy I: optical microscopy,” J. Microsc. **206**, 194–203 (2002). [CrossRef] [PubMed]

**17. **A. Barty, K.A. Nugent, D. Paganin, and A. Roberts, “Quantitative optical phase microscopy,” Opt. Lett. , **23**, 817–819 (1998). [CrossRef]

**18. **L. Mandel and E. Wolf, *Optical Coherence and Quantum Optics* (Cambridge, 1995).

**19. **M. J. Bastiaans, “Application of the Wigner distribution function to partially coherent light,” J. Opt. Soc. Am. A **3**, 1227–1238 (1986). [CrossRef]

**20. **K.-H. Brenner and J. Ojeda-Castaneda, “Ambiguity function and Wigner distribution function applied to partially coherent imagery,” Opt. Acta **31**, 213–223 (1984). [CrossRef]

**21. **M. J. Bastiaans and T. Alieva, “Wigner distribution moments in fractional fourier transform systems,” J. Opt. Soc. Am. A **19**, 1763–1773 (2002). [CrossRef]

**22. **C. K. Aruldoss, N. M. Dragomir, and A. Roberts, “Non-interferometric characterization of partially coherent scalar wavefields and application to scattered light,” J. Opt. Soc. Am. A **24**, 3189–3197 (2007). [CrossRef]

**23. **C. F. Bohren and D. R. Huffman, *Absorption and Scattering of Light by Small Particles* (Wiley1998). [CrossRef]

**24. **M. R. Teague, “Deterministic phase retrieval: a Green’s function solution,” J. Opt. Soc. Am. **73**1434–1441 (1983). [CrossRef]

**25. **T. E. Gureyev, A. Roberts, and K. A. Nugent, “Partially coherent fields, the transport-of-intensity equation, and phase uniqueness,” J. Opt. Soc. Am. A **12**, 1942–1946 (1995). [CrossRef]

**26. **D. Paganin and K. A. Nugent, “Non-interferometric phase imaging using partially coherent light,” Phys. Rev. Lett. **80**, 2586–2589.(1998). [CrossRef]

**27. **D. Paganin, A. Barty, P. J. McMahon, and K. A. Nugent, “Quantitative phase-amplitude microscopy. III. The effects of noise,” J. Microsc. **214**, 51–61 (2004). [CrossRef] [PubMed]