Home

Light Transport in Tissue Derivation of the diffusion equation The radiative transport equation is [8] $$(\hat\mathbf{s}\cdot\nabla)L(\mathbf{r},\hat\mathbf{s}) = -\... ...{s},\hat\mathbf{s}') L(\mathbf{r},\hat\mathbf{s}') \, d\omega'$$ (4.1) where $L(\mathbf{r},\hat\mathbf{s})$ is the radiance (W cm-2sr-1) at position r in the direction of the unit vector $\hat\mathbf{s}$. The differential solid angle $d\omega'$ has the unit vector $\hat\mathbf{s}'$ as an outward normal. The phase function $p(\hat\mathbf{s},\hat\mathbf{s}')$ represents the fraction of light scattered from the direction $\hat\mathbf{s}'$ into the direction $\hat\mathbf{s}$. The phase function is normalized such that it is unity when integrated over all directions. The scattering medium is an infinite slab with finite thickness and infinite breadth. The inward normal to the top of the slab is chosen as the positive z-direction. The top surface is illuminated with a monodirectional flux $\pi F_0(r)$ $$L_{\mathrm{inc}}(r,\hat\mathbf{s}) = \pi F_0(r) \delta((\hat... ...at\mathbf{s}_0\cdot\hat\mathbf{z}))= \pi F_0 \delta(\mu-\mu_0)$$ (4.2) where $\hat\mathbf{z}$ is a unit vector in the direction of the z-axis and $\hat\mathbf{s}_0$ is a unit vector in the direction of the incident flux and $\mu$ is the cosine of the angle $\hat\mathbf{s}$ makes with the z-axis. The delta function $\delta(\mu-\mu_0)$ is discussed in Appendix C. Historically, the factor of $\pi$ is included so that an isotropic diffuse source F0(r) will result in a flux of $\pi F_0(r)$ at a surface, $$\int_{2\pi\, \mu\ge0}F_0(r)(\hat\mathbf{s}\cdot\hat\mathbf{z... ... = F_0(r) \int_0^1\mu\,d\mu\int_0^{21\pi}\,d\varphi=\pi F_0(r)$$ (4.3) If the irradiance E0(r) represents collimated light, then $E_0(r)=\pi F_0(r)$. Notice that $L_{\mathrm{inc}}(\mathbf{r},\hat\mathbf{s})$ includes contributions from a cone containing all vectors having $\hat\mathbf{s}\cdot\hat\mathbf{z} = \mu_0$, and not just one particular azimuthal angle. This simplifies the mathematics by eliminating any azimuthal dependence of radiance in the slab. The radiance $L(\mathbf{r},\hat\mathbf{s})$ is a function of only the position r and the angle $\cos^{-1}(\hat\mathbf{s}\cdot\hat\mathbf{z})$. The phase function is modelled as a delta-Eddington phase function $$P_\mathrm{delta-E}(\cos\theta)= {1\over4\pi}\left\lbrace 2f\delta(1-\cos\theta)+(1-f)(1+3g'\cos\theta)\right\rbrace$$ (4.4) where $\cos\theta = (\hat\mathbf{s}\cdot\hat\mathbf{s}')$ is the cosine of the angle between the incident and the scattered light. The first term accounts for strong scattering in the forward direction and the second term approximates a more diffuse type of scattering. The parameter f determines the fraction of light scattered into the forward peak (the delta function) and g' denotes the degree of asymmetry in the diffuse portion of the scattering. The factor of $1/4\pi$ is included for normalization. Substituting Equation (4.4) into (4.1) yields $$(\hat\mathbf{s}\cdot\nabla)L(\mathbf{r},\hat\mathbf{s}) = -\... ...{s}') [1+3g'(\hat\mathbf{s}\cdot\hat\mathbf{s}')] \, d\omega'$$ (4.5) where the reduced scattering coefficient is $\mu_s'=\mu_s(1-f)$ and the transport coefficient is $\mu_t'=\mu_s'+\mu_a$. The transport $\mu_t'$ and reduced scattering $\mu_s'$ coefficients are less than the corresponding total attenuation $\mu_t$ and scattering coefficients $\mu_s$ respectively and represent effective total attenuation and scattering. The scattering coefficient is reduced because light scattered into the forward peak of Equation (4.4) is indistinguishable from unscattered light. The fraction of light not scattered into the forward peak is (1-f), and consequently, the scattering coefficient is reduced by a factor of (1-f). Thus the reduced scattering coefficient $\mu_s'$ represents the amount of light scattered out of the collimated portion of the radiance and into the diffuse portion of the radiance. The radiance is divided into collimated and diffuse components $$L(\mathbf{r},\hat\mathbf{s})=L_{\mathrm{coll}}(\mathbf{r},\hat\mathbf{s})+L_d(\mathbf{r},\hat\mathbf{s})$$ (4.6) The collimated radiance includes both the light scattered into a direction parallel to the incident beam and any unscattered light. Because the collimated radiance includes light scattered forward, the beam is attenuated not by the usual total attenuation coefficient $\mu_t$ but by the transport coefficient $\mu_t'$ The collimated light is attenuated at a rate proportional to the transport coefficient. $$(\hat\mathbf{s}_0\cdot\nabla)L_{\mathrm{coll}}(\mathbf{r},\h... ...athbf{s}) =-\mu_t'L_{\mathrm{coll}}(\mathbf{r},\hat\mathbf{s})$$ (4.7) The amount of collimated light entering the slab is given by Equation (4.2), less the amount of light lost to specular reflection from the surface $$L_{\mathrm{coll}}(\mathbf{r},\hat\mathbf{s})=(1-r_s)\pi F_0(r)\delta(\mu-\mu_0)$$ (4.8) where rs is the specular reflection coefficient given by the usual Fresnel reflection formula for specular reflection for an angle of incidence $\cos\theta_0=\hat\mathbf{s}_0\cdot\hat\mathbf{z}$ $$r_s = {1\over 2}\left[ {\sin^2(\theta_0-\theta_t)\over\sin^2... ...n^2(\theta_0-\theta_t)\over\tan^2(\theta_0+\theta_t)} \right]$$ (4.9) The incident and transmitted angles $\theta_0$ and $\theta_t$ are related by Snell's law $$n_{\mathrm{tissue}}\sin\theta_t = n_{\mathrm{air}}\sin\theta_0$$ where nair and ntissue are the indices of refraction of air and tissue. The solution of Equation (4.7) subject to the initial condition of Equation (4.8) is $$L_{\mathrm{coll}}(\mathbf{r},\hat\mathbf{s})=(1-r_s)\pi F_0(r)\exp(-\mu_t' z/\mu_0) \delta(\mu-\mu_0)$$ (4.10) where $z/\mu_0$ is the distance incident light travels in tissue to reach a depth z in the slab. Substituting Equation (4.6) into (4.5) and simplifying using Equations (4.7) and (4.10) yields $$(\hat\mathbf{s}\cdot\nabla)L_d(\mathbf{r},\hat\mathbf{s}) -\mu_t' L_d(\mathbf{r},\hat\mathbf{s}) + {\mu_s'\over4\pi}\int_{4... ...athbf{r},\hat\mathbf{s}') [1+3g'(\hat\mathbf{s}\cdot\hat\mathbf{s}')]\,d\omega'$$ = $${\mu_s'\over4\pi}(1-r_s)\pi F_0(r)\exp(-\mu_t' z/\mu_0)[1+3g'\mu_0\mu]$$ + (4.11) where $\hat\mathbf{s}\cdot\hat\mathbf{s}'$ has been rewritten in terms of the spherical angles determining $\hat\mathbf{s}$ and $\hat\mathbf{s}'$ If $\mu=\hat\mathbf{s}\cdot\hat\mathbf{z}$ and $\mu'=\hat\mathbf{s}'\cdot\hat\mathbf{z}$ then $$\hat\mathbf{s}\cdot\hat\mathbf{s}' = \mu\mu'-\sqrt{1-\mu^2}\sqrt{1-\mu'^2}\cos(\varphi-\varphi')$$ (4.12) The Eddington or diffusion approximation characterizes the diffuse radiance as a sum of a diffuse radiant fluence $\varphi _d(\mathbf {r})$ and a diffuse radiant flux per unit area Fd(r) These are defined as the first two moments of the radiance $L_d(\mathbf{r},\hat\mathbf{s})$ $$\varphi_d(\mathbf{r})=\int_{4\pi} L_d(\mathbf{r},\hat\mathbf... ...4\pi} L_d(\mathbf{r},\hat\mathbf{s}') \hat\mathbf{s}\,d\omega'$$ (4.13) The diffuse radiance can then be expressed as $$L_d(\mathbf{r},\hat\mathbf{s}) = {1\over4\pi}\varphi_d(\mathbf{r})+{3\over4\pi}\mathbf{F}_d(\mathbf{r})\cdot\hat\mathbf{s}$$ (4.14) The factors of $1/4\pi$ and $3/4\pi$ in Equation (4.14) result from normalization. Equation (4.14) represents the first two terms of the Taylor expansion for the diffuse radiance $L_d(\mathbf{r},\hat\mathbf{s})$, where $\varphi _d(\mathbf {r})$ represents the isotropic and Fd(r) the anisotropic contribution to the diffuse radiance. Recalling the solid angle integration identities (see Appendix C) $$\int_{4\pi}\hat\mathbf{s}\cdot\mathbf{A} \,d\mathbf{\omega} ... ...) \,d\mathbf{\omega} = {4\pi\over3}(\mathbf{A}\cdot\mathbf{B})$$ (4.15) and substituting Equation (4.14) into (4.11) and using (4.15) to simplify yields the following equation for the diffuse radiance $$(\hat\mathbf{s}\cdot\nabla)\varphi_d(\mathbf{r})+3(\hat\mathbf{s}\cdot\nabla)(\mathbf{F}_d(\mathbf{r})\cdot\hat\mathbf{s}) -\mu_a\varphi_d(\mathbf{r})-3(\mu_t'-g\mu_s')(\mathbf{F}_d(\mathbf{r})\cdot\hat\mathbf{s})$$ = (4.16) $$\mu_s'(1-r_s) \pi F_0(r)\exp(-\mu_t'z/\mu_0)(1+3g'\mu\mu_0)$$ + Integration of Equation (4.17) over all directions and using relations (4.15) results in the following equation for the diffuse flux $$\nabla\cdot\mathbf{F}_d(\mathbf{r}) =-\mu_a\varphi_d(\mathbf{r})+\mu_s'(1-r_s)\pi F_0(r)\exp(-\mu_t'z/\mu_0)$$ (4.17) Here the L. H. S. represents the net change in the diffuse radiant flux. This equals the intensity lost through absorption of the diffuse radiant fluence plus that gained through scattering of collimated light into the diffuse portion of the radiance. Multiplying Equation (4.17) by $\hat\mathbf{s}$ and integrating over all angles yields an energy flux equation for the diffuse radiance which states that the change in the diffuse fluence equals the diffuse flux lost to absorption plus that gained from collimated light. $$\nabla\varphi_d(\mathbf{r})=-3\mu_{tr}'\mathbf{F}_d(\mathbf{... ...\mu_s'(1-r_s)\pi F_0(r)\exp(-\mu_t'z/\mu_0)\mu_o\hat\mathbf{z}$$ (4.18) where $\mu_{\mit tr}'=\mu_a+\mu_s'(1-g')$ is the reduced transport coefficient. The reduced transport coefficient incorporates forward scattered light from the second term of the delta-Eddington phase function into the collimated beam in much the same way that such light for the first term was incorporated into expressions for the reduced scattering coefficient. Taking the divergence of Equation (4.18) and solving for $\nabla\cdot\mathbf{F}_d(\mathbf{r})$ yields $$\nabla\cdot\mathbf{F}_d(\mathbf{r}) =-{1\over3\mu_{tr}'}\nab... ...u_s'\mu_t'\over\mu_{tr}'}(1-r_s)\pi F_0(r)\exp(-\mu_t'z/\mu_0)$$ (4.19) Equating Equation (4.17) and (4.19) leads to a Helmholtz equation $$\nabla^2\varphi_d(\mathbf{r})-3\mu_{tr}'\mu_a\varphi_d(\math... ...mu_s'(\mu_{tr}'+\mu_t'g')(1-r_s)\pi F_0(r)\exp(-\mu_t'z/\mu_0)$$ (4.20) which is the well-known diffusion equation. The R. H. S. of this equation represents the collimated irradiation source.

S. A. Prahl."Light Transport in Tissue," PhD thesis, University of Texas at Austin, 1988.