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

Two-term scattering phase function for photon transport to model subdiffuse reflectance in superficial tissues

Open Access Open Access

Abstract

For Monte Carlo simulations of light transport in a variety of diffuse scattering applications, a single-scattering two-term phase function with five adjustable parameters is sufficiently flexible to separately control the forward and backward components of scattering. The forward component dominates light penetration into a tissue and the resulting diffuse reflectance. The backward component controls early subdiffuse scatter from superficial tissues. The phase function consists of a linear combination of two phase functions [Reynolds and McCormick, J. Opt. Soc. Am. 70, 1206 (1980) [CrossRef]  ] that were derived from the generating function for Gegenbauer polynomials. The two-term phase function (TT) accommodates strongly-forward anisotropic scattering with enhanced backscattering and is a generalization of the two-term, three-parameter Henyey-Greenstein phase function. An analytical inverse of the cumulative distribution function for scattering is provided for implementation in Monte Carlo simulations. Explicit TT equations are given for the single-scattering metrics g1, g2, γ, and δ. Scattering data from previously published bio-optical data are shown to fit better with the TT than other phase function models. Example Monte Carlo simulations illustrate the use of the TT and its independent control of subdiffuse scatter.

© 2023 Optica Publishing Group under the terms of the Optica Open Access Publishing Agreement

1. Introduction

1.1 Research motivation

An analytical scattering phase function (pf) model can be helpful for analyzing experimental measurements of backscattered reflectance from a turbid material in order to specify the optical properties of the material (the inverse problem), and subsequently predicting light transport behavior within the turbid material (the forward problem). In biological tissues, backscattered light is an opportunity for noninvasive diagnostic monitoring, and determining the tissue optical properties allows specification of light transport in the tissue and quantification of a dose of treatment light.

The nano-scale and micro-scale structures of a biological tissue scatter photons in a distinct angular pattern at each scattering event. This angular pattern can be analyzed to yield a set of single-scattering metrics such as the anisotropy of scattering traditionally denoted by $g_1$ and the backscatter fraction $b$. As will be shown, additional metrics such as $\gamma$ and $\delta$ can be obtained from higher-order angular moments of a pf. These metrics then become input to the categorization of single-scattering behavior and its relationship to the nano- and micro-scale structures of the tissue. The metrics can help categorize a single-scattering pf. The pf can be used in a Monte Carlo simulation to predict an experimental measurement like reflectance and transmission. Hence, the relation between single-scattering metrics and experimental measurements can be studied.

Sometimes experimental goniometric measurements of the angular pattern of single scattering is missing at certain angles such as near 0$^{\circ }$, −180$^{\circ }$, and $\pm 90^{\circ }$. An analytical pf can fill in such missing data.

The nano- and micro-scale structures of a tissue can be simulated by a scattering model, as for example by Lorenz-Mie theory of scattering by spheres or by other shapes. Then the scattering behavior, metrics, and categorization of well-defined but simulated “tissues" can be specified so the relationship between the nano- and micro-scale structures of a tissue and its observable single-scattering behavior can be explored. Ultimately the categorization can aid a decision in a medical, agricultural, or industrial biomaterials application.

Figure 1 summarizes the motivations for this investigation. The first is to further develop the two-term (TT) single-scattering pf by increasing the number of metrics (e.g., $g_2,\ g_3, \ \gamma,\ \beta$) that can be used to characterize the scattering properties. The second motivation is to apply those equations to investigate sub-diffusive scattering in tissues. The figure shows that a math description of light transport theory can generate simulated measurements. Such simulations can contribute to the categorization of a tissue on the basis of its multiple scattering behavior. While laboratory specimens can be assessed by single-scattering transmission experiments, clinical reflectance measurements on patients done noninvasively involve multiple scattering. Such measurements can be done either in the subdiffuse regime where only a few scattering events have occurred, or in the diffuse regime where photon trajectories are fully randomized and photon propagation occurs by diffusion directed by photon concentration gradients.

 figure: Fig. 1.

Fig. 1. Overview of how a math description of light scattering by a tissue can lead to categorization of the tissue type, which then aids medical decisions. Additionally, the math description can be used in light transport theory to simulate experiments and generate simulated measurements, which then can contribute to the analysis. Hence knowledge of multiple scattering behavior in both the subdiffuse and diffuse regimes can participate in characterizing a tissue. The dashed line indicates how a simulation can augment experimental data by simulating missing data. In red, a tissue site on a patient or biomaterial of interest can be characterized by the analysis module to yield an assessment of the tissue’s nano- and micro-scale structural status to inform clinical or industrial decisions.

Download Full Size | PDF

1.2 History of scattering phase functions

Monte Carlo calculations simulate how photons are scattered in a medium. For azimuthally-symmetric scattering a pf $f(\mu )$ [sr$^{-1}$], $-1 \le \mu \le 1$, is defined to satisfy the normalization condition

$$2\pi\int_{{-}1}^1 f(\mu) d\mu = 1$$
for $\mu \equiv \cos \theta$ in terms of the scattering polar angle $\theta$, i.e., the angle of photon deflection. For any pf, two desired properties are that it be non-negative and that it can be sampled with random numbers so that a lookup table is not needed. If a pf model requires a series expansion in $\mu$ then convergence issues arise unless that expansion can be expressed in a closed-form, analytical generating function. As an example, consider a generating function for the classical Legendre polynomials $P_n(\mu )$ that satisfies the equation
$$(1+g^2-2g\mu)^{{-}1/2} = \sum_{n=0}^\infty P_n(\mu)g^n ,$$
for $-1 \le g \le 1$. This generating function, however, does not satisfy Eq. (1) but Henyey and Greenstein (HG) [1] published their one-parameter generating function for Legendre polynomials,
$$(4\pi)^{{-}1}(1-g^2)(1 + g^2 -2g\mu)^{{-}3/2} = (4\pi)^{{-}1} \sum_{n=0}^\infty (2n+1)P_n(\mu)g^n,$$
that does. Later Reynolds and McCormick (RM) [2] developed a generalization of the HG pf in terms of Gegenbauer polynomials $C_n^{(\alpha )}(\mu )$ that also satisfies Eq. (1). That pf provides an additional free parameter $\alpha$, besides $g$, to enable more flexibility in modeling strongly forward scattering.

A major advantage of such closed-form pfs, obtained from generating functions, is that Monte Carlo sampling of the new direction after scattering can be explicitly given [3,4] in terms of a random number $\eta$, $0 \le \eta \le 1$. Such is not the case for phase functions such as Lorenz-Mie scattering [5], for example, where pre-computed lookup tables are needed.

Another advantage of the HG and RM phase functions is that directional moments of the pf, as derived here, can yield analytical metrics to describe the scattering properties of an unknown medium. The analytical equation for the backscattering fraction of a phase function, for example, is a directional moment that has been shown to be very important [68] for light transport applications when interpreting reflectance data. By analogy, the forward-scattering fraction is important for interpreting transmittance data.

The one-term RM pf (OT) has been used in a wide variety of applications, including diffusers in lighting and display applications [9], spatially resolved reflectance in a subdiffusive regime [1012], and light scattering in human blood [2,1315] and tissue [16].

The OT $f(\mu )$ also has been used to analyze refractive index fluctuations in biological tissue [17] because it is the phase function associated with the Whittle-Matérn correlation function [3]. The OT also has been obtained from the small-angle-scattering approximation for the radiative transfer equation, which has led to its use for analyzing the propagation of light in media with large-scale inhomogeneities [1820], or when seeking the reflectance for grazing angles of incidence on a turbid medium [21,22].

A disadvantage of the OT $f(\mu )$ is that it is monotonic, so specifying a strongly-forward scattering pf also causes a very low backscatter, often lower than experimentally seen in scattering media, including tissues. To accommodate some backscattering in a pf with strongly-forward scattering, a linear combination of two phase functions, one directed primarily in the forward direction with parameters $\alpha _f$ and $g_f$ and one directed primarily in the backward direction with $\alpha _b$ and $g_b$, provides separate controls for forward and backward scattering. Such a combination can be implemented with either a three-parameter model of two one-term HG (OTHG) phase functions [23], a four-parameter model of one OT and one OTHG [24], or a five-parameter two-term pf. Although there are several two-term pf models, this paper considers the five-parameter, two-term phase function (TT) derived from the generating function for Gegenbauer polynomials. A special case of the TT is the two-term Henyey-Greenstein (TTHG) with $\alpha _f = \alpha _b = 1/2$. Recently a TT with $\alpha _f = \alpha _b = 1$ has been shown to be excellent for interstellar dust [25].

The general TT recently has been selected as the preferred pf to study light scattering by colloidal nanospheres [26,27] and also for modeling microalgae and mineral hydrosols [28] and polydisperse spherical particles [29]. Harmel et al. [28] have shown for optical oceanography applications the TT provides a better fit to Lorenz-Mie scattering than the one-term, two-parameter Fournier-Forand (FF) phase function model [8,30,31] that also is an excellent closed-form analytical model that does not require a series expansion to be evaluated. The two parameters of the FF model depend on the ratio of the circumference of a spherical particle to the wavelength of the light in the medium and the refractive index at the particle surface. Thus the two FF parameters are physical properties, whereas the two parameters $\alpha$ and $g$ of the OT and the five parameters of the TT are not.

Ma et al. [29] have investigated the influence of the TT scattering phase function on the radiative hemispherical reflectance and transmittance using Monte Carlo calculations for a layer of partially disperse particle sizes with a TT fit to Mie scattering.

1.3 Objective of paper

The general TT will be shown to be a flexible and accurate descriptor of experimentally observed single-scattering in tissues. The objective of this research is to take advantage of few-scattered photons for subdiffusive radiative transfer. Although there may be biomedical or other applications for which not all the generality of the five-parameter TT is needed, it will be shown that for subdiffuse photon transport the generality is needed.

In this paper, the terms “forward scatter" and “backward scatter" are used to indicate the probabilities of forward and backward scatter or to indicate forward or backward scattered photons.

Section 2 presents the OT $f(\mu )$ including the Monte Carlo sampling of the OT for use in simulations and the new metrics for the RM pf. Section 3 gives the corresponding results for the TT pf. Section 4 illustrates how the TT can match previously published goniometric data for $f(\mu )$ for cells and blood. Section 5 presents Monte Carlo simulations that include the flux escaping near a source, reflectance measurements from tissue, and distinguishing early subdiffuse versus diffuse measurements to predict how practical measurements will vary with changes in a TT. Section 6 concludes with a summary and overview of potential future applications of the TT for interpreting biological reflectance measurements from superficial tissues.

2. One-term phase function

2.1 OT phase function definition

The OT $f(\mu )$ can be expressed as [2],

$$f(\mu) = K(1+g^2 - 2g\mu)^{-(\alpha + 1)}, \quad \alpha >{-}1/2,\ -1 \le g \le 1,$$
with the normalization
$$K = \pi^{{-}1}\alpha g(1-g^2)^{2\alpha} [(1+g)^{2\alpha} - (1-g)^{2\alpha}]^{{-}1} \,.$$

The OTHG is a special case of the OT for $\alpha = 0.5$. The parameter $g$ for $g \in (0,1]$ creates a “forward bias” scattering parameter, and similarly if $g \in [-1,0)$ there is a backward scattering bias. The $\alpha$, $\alpha > -0.5$, acts to enhance the effects of $g$.

2.2 Influence of the OT parameters

The OT scattering function has two parameters $g$ and $\alpha$ to control its angular dependence. Figure 2 illustrates the influence of the parameter $\alpha$ on the scattering function $f(\mu )$ while holding $g$ constant at 0.90. The dynamic range between backward scatter and forward scatter increases as $\alpha$ increases. Figure 2(A) plots $f(\mu )$ as $\alpha$ is varied over 0.5, 1, 2, and 3. Higher $\alpha$ decreases the intensity of backscatter by 5 orders of magnitude (see Fig. 2(C)). Higher $\alpha$ linearly increases the intensity of forward scatter (see Fig. 2(D)).

 figure: Fig. 2.

Fig. 2. One-term phase function. As the parameter $\alpha$ increases, the backscatter falls exponentially while the forward increases linearly. (A) The scattering function $f(\mu )$ is shown as $\alpha$ is varied over 0.5, 1, 2, and 3, while $g$ is constant at 0.90. (B) Close-up view of (A) to show forward scatter in the 0-8$^{\circ }$ range. (C) Backscatter at 180$^{\circ }$ drops exponentially vs $\alpha$. (D) Forward scatter at 0$^{\circ }$ increases linearly vs $\alpha$.

Download Full Size | PDF

Such flexibility in the dynamic range between forward and backward scattering allows for a description of the scattering by biological materials. The details of forward scattering are important for the penetration of light into a tissue and strongly affect the diffuse reflectance. Penetration to a focus within a biological tissue during microscopy is an example. The details of backcattering are important for subdiffuse reflectance from the superficial layer of a tissue.

2.3 OT Monte Carlo sampling

To obtain the equation with which to determine the post-scattering direction $\mu$ with respect to its initial direction $\mu =0$, for a single random number $\eta$ uniformly distributed over $[\,0,1]$, first the integral

$$\eta = 2\pi\int_{{-}1}^{\,\mu}f(\mu') d\mu' \,$$
is evaluated. An algebraic rearrangement of the result then gives the new direction $\mu$ for a Monte Carlo sampling in terms of $\eta$ [3],
$$\mu = \frac{1+g^2}{2g} -\frac{1}{2g} \left[ \frac{\eta}{(1-g)^{2\alpha}} + \frac{1-\eta}{(1+g)^{2\alpha}}\right]^{{-}1/\alpha}, \ -1 \le \mu \le 1$$
that is used here instead of an expression given earlier [4],
$$\mu = \frac{1+g^2}{2g} - \frac{1}{2g\sqrt[\alpha]{\varrho}}\,$$
where
$$\varrho = \frac{2\alpha g\eta}{K} + \frac{1}{(1+g)^{2\alpha}}\,$$
and $K$ is in Eq. (5). Both Eqs. (7) and (8a) yield the same sampled value of $\mu$.

For the OTHG obtained by setting $\alpha = 0.5,$ from Eq. (7) it follows that a random number $\eta$ will give the sampled value of $\mu$ as [32]

$$\mu = \frac{1+g^2}{2g} - \frac{1}{2g}\left[\frac{1-g^2}{1-g+ 2g\eta}\right]^2, \ -1 \le \mu \le 1 \,.$$

2.4 OT metrics for $\mu \in [-1,1]$

Directional moments of a phase function $f(\mu )$ are metrics that yield information about the shape of the phase function as a function of $\mu$. These metrics serve to distinguish one type of scattering medium from another, but some metrics are more useful than others. The most important metric is the expectation value for the new scattering direction projected onto the original one, i.e., the mean cosine of the scattering angle.

The Legendre polynomial $P_n(\mu )$ moments of a phase function are important in some diffusion theory applications and are given by the ratio

$$g_{\,n} = \frac{2\pi \int_{{-}1}^1 P_n(\mu) f(\mu) d\mu}{2\pi \int_{{-}1}^1 f(\mu) d\mu}, \quad n= 0,\ 1,\ 2, \ldots\;,$$
with $g_0 \equiv 1$. For $n=1$ with $P_1(\mu ) = \mu$, the value of $g_1$ is the mean or expectation value of $\mu$ [2],
$$g_1 = \frac{2g\alpha L - (1+g^2)}{2g(\alpha - 1)}, \quad \alpha \ne 1$$
where $L$ is
$$L = \frac{(1+g)^{2\alpha} + (1-g)^{2\alpha}}{(1+g)^{2\alpha} - (1-g)^{2\alpha}}.$$
(The numerator factor $L$ compensates for the denominator factor ($\alpha$-1) as $\alpha$ transitions from <1 to >1 such that $g_1$ transitions smoothly and only becomes ill-defined precisely at $\alpha$ = 1. Similarly, this occurs for other metrics.) The OTHG with $\alpha = 1/2$ is just
$$g_1 = \frac{1 + g^2}{2g} + \left(\frac{1-g^2}{2g}\right)^2 \ln\left(\frac{1-g}{1+g}\right) .$$

Figure 3(A) plots the anisotropy $g_1$ versus $g$ for a range of $\alpha$ values. Increasing $g$ causes $g_1$ to increase toward 1.0. Higher $\alpha$ accelerates this rate of increase in $g_1(g)$. Figure 3(B) plots $g_1$ versus $\alpha$ for a range of $g$ values.

 figure: Fig. 3.

Fig. 3. (A) Anisotropy $g_1$ versus $g$ for a range of $\alpha$ values. (B) $g_1$ versus $\alpha$ for a range of $g$ values

Download Full Size | PDF

Equation (11) gives the metrics $g_n$ for phase functions other than the Henyey-Greenstein, which do not have a direct physical meaning yet are very useful. For a positive $g$ the value of $g_1$ is the “forward bias” of the scattering directions; similarly, a negative $g$ indicates a backscattering bias. For computations with the OT, the usual reduced scattering coefficient $\mu _s'$ of traditional diffusion theory with the OTHG, where $g_1=g$ and $\mu _s' = \mu _s(1-g),$ should be replaced by

$$\mu_s' = \mu_s(1-g_{\,1})\,$$
with the $g_1$ of Eq. (11) and not the parameter $g$. For $n=2$, Eq. (10) gives the metric
$$g_{\,2} = \frac{3(1+g^2)g_{\,1}}{2g(2-\alpha)} - \frac{1+\alpha}{2-\alpha}, \quad \alpha \ne 2$$
with $g_{\,2} = g^2$ for $\alpha = 0.5$. For $n=3$,
$$g_{\,3} = \frac{5}{2(3-\alpha)}\left[ \frac{1+g^2}{g} g_{\,2} - \alpha L + \frac{1+g^2}{2g} \right] - \frac{3}{2} g_{\,1} , \quad \alpha \ne 3$$
with $g_{\,3} = g^3$ for $\alpha = 0.5$. Equations (11), (15), and (16) enable us to generalize the helpful parameters $\gamma$ and $\delta$ that have been used in diffusion theory modeling [11,33,34] to help characterize the reflectance from a scattering medium,
$$\gamma = \frac{1 - g_{\,2}}{1 - g_{\,1}} $$
$$\delta = \frac{1 - g_{\,3}}{1 - g_{\,1}}. $$

2.5 OT metrics for $\mu \in [-1,0]$ and $\mu \in [0,1]$

The backscattering fraction $b$ is very important for interpreting reflectance measurements from a target material. It is defined by

$$b = \frac{2\pi\int_{{-}1}^0 f(\mu)d\mu}{2\pi \int_{{-}1}^1 f(\mu) d\mu}$$
and for the one-term $f(\mu )$ it is
$$b = \frac{(1-g)^{2\alpha} }{(1+g)^{2\alpha} - (1-g)^{2\alpha}} \left [\frac{(1+g)^{2\alpha}}{(1+ g^2)^\alpha} - 1 \right]$$
or
$$b = \frac{U - 1}{V-1}$$
with
$$U = \frac{(1+g)^{2\alpha}}{(1+g^2)^\alpha} \quad {\mathrm{and}} \quad V = \frac{(1+g)^{2\alpha}}{(1-g)^{2\alpha}}\,.$$

For $\alpha = 0.5$ the Henyey-Greenstein backscattering fraction is obtained [7,35],

$$b = \frac{1-g}{2g}\left [ \frac{1+g}{(1+g^2)^{1/2}} - 1 \right ]\,.$$

Figure 4(A) plots the backscatter fraction $b$ versus $g$ for a range of $\alpha$ values. Increasing $g$ causes $b$ to drop. Higher $\alpha$ accelerates this rate of drop in $b$. Figure 4(B) plots $b$ versus $\alpha$ for a range of $g$ values.

 figure: Fig. 4.

Fig. 4. (A) Backscatter fraction $b$ versus $g$ for a range of $\alpha$ values. (B) $b$ versus $\alpha$ for a range of $g$ values

Download Full Size | PDF

Higher-order backscattering angular moments can be defined by the metrics

$$b_n = ({-}1)^n \frac{2\pi \int_{{-}1}^0 P_n(\mu) f(\mu) d\mu}{2\pi\int_{{-}1}^1 f(\mu) d\mu}, \quad n= 0,\ 1,\ \ldots\;,$$
with $b_0 \equiv b$ of Eq. (19). The factor $(-1)^n$ is chosen to make the $b_n$ nonnegative. For example,
$${} b_1 = \frac{2\pi K(1+g^2)}{(2g)^2 \alpha (\alpha - 1)} \left[\frac{1}{(1+g^2)^\alpha} - \frac{1}{(1+g)^{2\alpha}} \right], \quad \alpha \ne 1\,.$$

Now consider just the fraction of the backscattering that projects onto the backward direction. A new reflectance metric $h_R$ is

$$h_R ={-} \frac{2\pi \int_{{-}1}^0 \mu f(\mu) d\mu}{2\pi\int_{{-}1}^0 f(\mu) d\mu} \equiv b_1/b_0 \,.$$

The $h_R$ is of practical use where the fraction of total backscattering that directly scatters into a narrow solid angle of collection may characterize the apparent size of the scatterers in a superficial tissue layer. Hence it is a good metric for interpreting subdiffuse scattering applications.

In a similar way, forward scattering angular moments can be defined by the metrics

$$f_n = \frac{2\pi \int_0^1 P_n(\mu) f(\mu) d\mu}{2\pi\int_{{-}1}^1 f(\mu) d\mu}, \quad n= 0,\ 1,\; \ldots\;,$$
and from Eqs. (11) and (23) it follows that $f_n = g_n - b_n$ with $f_0 = 1-b_0$. The values of $f_0$ and $f_1$ follow immediately from the equations for $b_0$ and $b_1$ by replacing $g$ by $-g$.

The mean or expectation value of $\mu$ for just the forward-scattered photons,

$$h_T = \frac{2\pi \int_0^1 \mu f(\mu) d\mu}{2\pi\int_0^1 f(\mu) d\mu} \equiv f_1/f_0,$$
is a new transmittance metric that is of practical use where strongly-directed forward scattering allows a focused beam of photons to reach the focus within a tissue.

Figure 5(A) plots $h_R$ versus $g$ for a range of $\alpha$ values. Increasing $g$ causes $h_R$ to decrease. Higher $\alpha$ accelerates this rate of drop in $h_R$. Figure 5(B) plots $h_R$ versus $\alpha$ for a range of $g$ values. Figure 5(C) plots $h_T$ versus $g$ for a range of $\alpha$ values. Figure 5(D) plots $h_T$ versus $\alpha$ for a range of $g$ values.

 figure: Fig. 5.

Fig. 5. (A) Fraction of backward scattering that projects onto the backward direction, metric $h_R$, plotted versus $g$ for a range of $\alpha$ values. (B) $h_R$ plotted versus $\alpha$ for a range of $g$ values. (C) Fraction of forward scatter that projects onto the forward direction, metric $h_T$, plotted versus $g$ for a range of $\alpha$ values. (D) $h_T$ plotted versus $\alpha$ for a range of $g$ values.

Download Full Size | PDF

In Section 3, when the TT metrics $h_R$ and $h_T$ are introduced, comprised of a forward-oriented OT and a backward-oriented OT, both the $h_R$ of the forward-oriented OT and the $h_T$ of the backward-oriented OT contribute to the narrow solid angle of backscattering from scatterers located within the superficial region of a biological tissue.

3. Two-term phase function

3.1 TT phase function definition

The TT requires summing two OT phase functions $f_f(\mu )$ and $f_b(\mu )$ to obtain a $f(\mu )$ that satisfies Eq. (1) and to obtain both an enhanced backward-directed component as well as the strongly-forward scattering component. The TT has four parameters $g_f$, $\alpha _f$, $g_b$, $\alpha _b$, and can be denoted with two coupling parameters $C_f$ and $C_b$ to control its angular dependence. The constraint on the coupling coefficient is such that, with $0 \le C_f \le 1$ and $C_b \equiv (1-C_f)$, Eq. (1) is preserved to give a total of five free parameters. Henceforth, the subscripts $b$ and $f$ are used in equations and numerical results are reported with $C \equiv C_f$.

The forward scattering $g$ is specified by $g_f \in [0,1]$ and the backscattering $g$ by $-g_{\,b} \in [0,1]$ for graphical purposes and the ease of comparing the magnitude of the forward bias $g_f$ to $|\! -g_b|$. Note that the factor $K$ in Eqs. (4) and (5) now becomes $K_f$ or $K_b$. The TT $f(\mu )$ then can be written as

$$\begin{aligned} f(\mu) &= C_b f_{\,b}(\mu) + C_f f_{f}(\mu)\\ &= \sum_{i = \{b\,\&f\}} C_i f_{\,i}(\mu) \end{aligned}$$
for $f_i(\mu )$, $i = b$ and $f$. (The second expression of Eq. (27) will be used in later equations.)

Figure 6 illustrates the forward and backward scatter components of the scattering functions using polar plots of the strength of scatter [$sr^{-1}$] versus the angle of scatter [$^{\circ }]$. In each figure, the black line is a Henyey-Greenstein function while the colored lines are Reynolds-McCormick functions.

 figure: Fig. 6.

Fig. 6. Polar plots emphasizing the backward scatter of the forward and backward components of the two-term scattering function. (A) Forward scatter component. $\alpha _f$ is varied (0.001, 0.1, 0.3, 0.5, 0.7, 0.9, 0.95) and $g_f$ = 0.90. (B) Backward scatter component. $\alpha _b$ is similarly varied and $g_f$ = 0.90. Black curve indicates Henyey-Greenstein scattering ($\alpha _f$ or $\alpha _b$ = 0.5). Other colors vary $\alpha _f$ and $\alpha _b$.

Download Full Size | PDF

3.2 Influence of the TT parameters

Figure 7 illustrates the influence of the parameter $g_b$ (varied as 0.1, 0.7, and 0.9) on the scattering function $f(\mu )$, while the other parameters are held constant ($\alpha _f$ = 1, $g_f = 0.8$, $\alpha _b$ = 1, $C = 0.98$). The high $\alpha _f$ and $g_f$ provide a strong forward scatter and a low backscatter. The $\alpha _b$ and $g_b$ provide a backscatter peak that exceeds the backscatter due to $\alpha _f$ and $g_f$. The figure also shows the backscatter pattern for the three values of $g_b$ becoming more narrowly directed as $g_b$ increases.

 figure: Fig. 7.

Fig. 7. Two-term $f(\mu )$ for $\alpha _f = 1,\ g_f = 0.8$, and $C = 0.98$. The forward scatter $f_f(\mu )$ (black line), backward scatter $f_b(\mu )$ (red line), and the total scatter $f(\mu )$ (green line) are shown. Note the backward scattering of $f_f(\mu )$ and the forward scattering of $f_b(\mu )$ both contribute to the total backscattering of $f(\mu )$. (A) $g_b$ = 0.1 yields a broad backscatter, (B) $g_b$ = 0.7 yields a more directional backscatter, and (C) $g_b$ = 0.9 yields a sharply peaked backscatter.

Download Full Size | PDF

3.3 TT Monte Carlo sampling

With Eqs. (1) and (27) the new direction of scattering for TT Monte Carlo calculations can be sampled by

$$\begin{aligned} \mu &= C_b \, \mu_{b} + C_f \,\mu_{f}, \quad -1 \le \mu \le 1\\ &= \sum_{i = \{b\,\&f\}} C_i \mu_{i} \end{aligned}$$
with this $\mu$ now the sum of the two sampled one-term $\mu$ values of Eq. (7). Note this requires the sampling of a separate random number $\eta _f$ for the forward $f_f(\mu )$ and $\eta _{\,b}$ for the backward $f_{\,b}(\mu )$. (See Appendix.)

3.4 TT metrics for $\mu \in [-1,1]$

The classic diffusion theory $P_n(\mu )$ moments of the TT $f(\mu )$ are given by the ratio

$$g_n = \frac{2\pi\int_{{-}1}^1 P_n(\mu) f(\mu) d\mu}{2\pi\int_{{-}1}^1 f(\mu) d\mu}, \quad n = 1,\ 2,\ \ldots\ .$$

From Eqs. (1), (10), and (27), in terms of the one-term $P_n(\mu )$ parameters $g_{n,i}$, the result is

$$g_n = \sum_{i = \{b\,\&f\}} C_i \,g_{n,i} .$$

Because $g_1$ represents the net flow of photons in the forward direction for the TT representation of the scattering, just as it does for the OT, the reduced scattering coefficient $\mu _s' = \mu _s(1-g_1)$ used in diffusion theory should be

$$\mu_s' = \mu_s[1 - \sum_{i = \{b\,\&f\}} C_i \,g_{1,i}]\,.$$

3.5 TT metrics for $\mu \in [-1,0]$ and $\mu \in [0,1]$

Backscattering TT metrics are defined similarly to Eq. (22) and are

$$\begin{aligned} b_n &= ({-}1)^n \frac{2\pi \int_{{-}1}^0 P_n(\mu)f(\mu) d\mu}{2\pi\int_{{-}1}^1 f(\mu) d\mu}, \quad n=0, 1,\ 2,\ \ldots \,\\ &= ({-}1)^n \sum_{i = \{b\,\&f\}} C_i 2\pi \int_{{-}1}^0 P_n(\mu) f_i(\mu) d\mu\\ &= \sum_{i = \{b\,\&f\}}C_i \,b_{n,i}\ , \end{aligned}$$
where the $b_{n,i} = b_{n,b}$ or $b_{n,f}$ of Eq. (22) for $i = b$ or $f$, respectively, and with $b_{0,i}=b_i$ the backscattering fraction of Eq. (19) for $i=b$ or $f$.

In a similar way to Eq. (24), the TT value corresponding to the sum of two OT components of $h_{R,i} = h_{R,b}$ or $h_{R,f}$ can be used to investigate the backward flow of backscattered light. The result is

$$h_R ={-} \frac{2\pi \int_{{-}1}^0 \mu f(\mu) d\mu}{2\pi \int_{{-}1}^0 f(\mu) d\mu} = \sum_{i = \{b\,\& f\}}C_i \frac{b_{1,i}}{b_i}\left[\sum_{i = \{b\,\& f\}}C_i b_i \right]^{{-}1} \,.$$

The corresponding forward scattering TT metrics can be defined by

$$\begin{aligned} f_n &= \frac{2\pi \int_0^1 P_n(\mu)f(\mu) d\mu}{2\pi\int_{{-}1}^1 f(\mu) d\mu}, \quad n=0,\ 1,\ 2,\ \ldots\\ &= \sum_{i = \{b\,\&f\}}C_i \,f_{n,i}\, , \end{aligned}$$
where $f_{n,i} = f_{n,b}$ or $f_{n,f}$ of Eq. (25) for $i = b$ or $f$, respectively. In a similar way, the corresponding TT obtained from the sum of two OT components can be used to investigate the forward flow of light. The result is
$$\begin{aligned} h_T &= \frac{2\pi \int_0^1 \mu f(\mu) d\mu}{2\pi \int_0^1 f(\mu) d\mu}\\ &= \sum_{i = \{b\,\&f\}}C_i \frac{f_{1,i}}{f_{i}}\left[\sum_{i = \{b\,\&f\}}C_if_i \right] ^{{-}1} \,. \end{aligned}$$

4. TT matches experimental data

Experimental measurements of the angular pattern of light scattering from cells or thin tissues indicate a very forward-directed $f(\mu )$, with a small elevation of backscatter centered near $\mu$ = −1 or $\theta = 180^{\circ }$. For example, Fig. 8 shows the scatter from cells measured by Mourant et al. [36].

 figure: Fig. 8.

Fig. 8. Scattering by cells (Mourant et al. [36]). Data is red line. (A) Two-term Henyey-Greenstein (TTHG), dashed line, does not fit the data. (B) Two-term (TT) fits the data.

Download Full Size | PDF

First consider a two-term Henyey-Greenstein scatter function (TTHG). Figure 8(A) shows data (red line) plus the fit (black dashed line) by the TTHG. The limited dynamic range of the TTHG function makes it difficult to fit both the strong forward scatter and the peak of backward scatter centered around $\theta = 180^{\circ }$. In the effort to match the strong forward scatter, the HG function generates a rather high backscatter (cyan line) that exceeds the data. Therefore, the backward-directed HG function (green line) cannot lower the backward fit to match the data.

Figure 8(B) uses the TT to fit Mourant’s cellular data. The fit closely approximates the measured data (red line). The curves show the forward contribution (cyan line), the backward contribution (green line), and the total scatter (black dashed line).

The optical scattering of blood was reported by Yaroslavsky et al. [4], who measured the angular scattering from a 100-$\mu$m-thick layer of blood. Figure 9 shows their data (red circles) and the TT fit to the data (red line). The black x’s are reported data that are ignored here because they are likely in error, perhaps since measurements as $\mu$ approaches zero are difficult. The forward (blue line) and backward (green line) components of the TT are shown, indicating that the backward component is needed to fit the backward peak centered around $\mu$ = −1. This example also illustrates how the scattering theory module of Fig. 1 can simulate missing experimental data.

 figure: Fig. 9.

Fig. 9. Scattering function ($f(\mu )$ [sr$^{-1}$] from a 100-$\mu$m-thick layer of blood Yaroslavsky et al. [4]. The red circles are the data fit by the TT, and the black x’s are data ignored by the fitting. The blue and green lines are the contributions of the forward and backward components of the TT, respectively. The red line is the total TT, and predicts the backscatter behavior in the angular range near $\mu$ = −1 where data are difficult to measure.

Download Full Size | PDF

These examples of scattering functions from biological materials illustrate the need for the TT. The backward term allows fitting the peak of backscatter, and the $\alpha _f$ of the forward OT provides the wider dynamic range needed to accommodate the very forward-directed scattering of tissues.

5. Monte Carlo simulations

The TT Monte Carlo program of this paper was adapted by a small change (see Appendix) in the Monte Carlo simulation that was the core of the MCML program [38].

5.1 Local escaping flux density near a source

Figure 10 shows a point spread function simulation of the fluence within a semi-infinite tissue and the flux density escaping the tissue as a function of the lateral distance $r$ in response to a very narrow pencil beam of collimated light at $r=0$. The simulations used a matched boundary condition at the tissue surface. The Monte Carlo simulation uses the sampling of TT described in Section 3.3. The parameter $g_b$ was varied as 0, .1, .3, .5, .7, .9, i.e., from an isotropic value (0) to a very backward-directed value (0.9). The other parameters ($\alpha _f,\ \alpha _b,\ C$) were kept constant. The absorption and scattering coefficients were kept constant ($\mu _a$ = 1 cm$^{-1}$ and $\mu _s$ = 100 cm$^{-1}$). The simulations kept the overall value of the anisotropy, $g_1$, constant at 0.90, which caused the diffusion of light distant from the point source of light entry to behave the same for all simulations, as shown for $r$ greater than about 0.06 cm. However, as $g_b$ increased, the $g_f$ that controlled the forward component of the TT was also increased to keep $g_1$ constant. The consequence was that the escaping flux density near the source decreased as $g_b$ increased, which at first seems counter-intuitive. Should not a more backward-directed backscatter cause escaping flux density near the source to increase? The explanation is that $g_f$ is scaled by the high value of $C$ (0.90) so the increase in $g_f$ outweighed the increase in $g_b$.

 figure: Fig. 10.

Fig. 10. (A) Monte Carlo simulation of the point spread function of light within an example tissue that has a scattering anisotropy of $g_1 = 0.9$. (B) The lateral profile of escaping flux density $R(r)$ [cm$^{-2}$]. (C) R $\leq$ 0.6 mm equals the integral of $R(r)$ over $r=0$ to 0.06 cm, which has a contribution from subdiffuse superficially scattered light and a contribution from diffuse multiply scattered light. As the backscatter component becomes increasing isotropic ($g_b$ approaching 0), the magnitude of backscatter reaches a maximum when $g_1$ is constant.

Download Full Size | PDF

When conducting experimental measurements of escaping flux density in response to a narrow pencil beam collimated source, one can capture the flux density as a function of $r$ from the source, which depends on the absorption ($\mu _a$) and the reduced scattering ($\mu _s(1-g_1)$) when distant from the source. Then the behavior near the source could be analyzed to specify the value of $g_b$. If a broader beam is used, convolution of the point spread function of Fig. 10 over the beam area will yield a result for flux density escaping both inside and just outside the beam, which can be analyzed to specify $g_b$.

Figure 10(C) shows the integral of escaping flux density over $r = 0$ to 0.0600 cm, $R_{\leq {\mathrm {0.6mm}}}$, which decreases as $g_b$ increases. The behavior suggests separate contributions to $R_{\leq {\mathrm {0.6mm}}}$ from the superficial subdiffuse scatter and the deeper multiply scattered light.

5.2 Simulated practical measurements in vivo

In practical experimental measurements in vivo, diffuse reflectance measurements are relatively easy to acquire. Both the total reflectance $R_d$ and the lateral point spread function of escaping flux density $R(r)$ [mm$^{-2}$], where $r$ is the distance from a point source, depend on the number density of scatterers $\rho$ [1/$\mu$m$^3$], the size of the scatterers $dia$ [$\mu$m], the photon wavelength $\lambda$ in the medium, and the refractive indices of the scatterer and surrouding medium. The $\rho$ and $dia/\lambda$ specify the scattering coefficient $\mu _s$ [mm$^{-1}$], and $dia/\lambda$ specifies the anisotropy of scattering $g_1$. The combination $\mu _s(1-g_1)$, called the reduced scattering coefficient, together with the absorption coefficient $\mu _a$ [cm$^{-1}$] govern the values of $R_d$ and $R(r)$ at large r.

However, as indicted by Fig. 10, the first several scattering events that occur within the superficial tissue layers, i.e., the epithelium or in skin the epidermis and the upper papillary dermis, can contribute to reflectance at the spot of light delivery. This superficial contribution is referred to as subdiffuse reflectance, and is sensitive to the backward component of the single-scattering function $f(\mu$) for $-1\le \mu \le 0$. The difference between the reflectance at the delivery site expected from multiply scattered light and the experimentally observed reflectance can be attributed to the first several scattering events in the subdiffuse domain. Mapping this difference signal over a tissue surface may provide a map of superficial tissue structure, indicating regions of optically perturbing structures.

The TT scattering function offers flexibility in characterizing the scattering function of a tissue. The total $g_1$ of the TT governs the diffuse $R_d$ and $R(r)$ at distant $r$. The backward-directed component of the TT can dominate over the backward portion of the forward-directed component of the TT, and hence dominate the subdiffuse reflectance. In particular, the backscatter $b$ and the direct backscatter $h_R$ can map to the number density of scatterers $\rho$ and the apparent size of scatterers $dia$ for a given $\lambda$ using the Mie sphere diameter as a model.

5.3 Distinguishing subdiffuse from diffuse scattering in superficial tissues

Much pathology arises in the superficial layer (epithelial or epidermal) of a tissue, which can best be seen with the few-scattering events in subdiffusive light scattering. Imaging with polarized light, in particular the difference image equal to a co-polarized image minus a cross-polarized image, yields an image based on only the still polarized photons, which is a superficial tissue image [37]. So it is possible to isolate the subdiffuse scatter from the superficial tissue layer. Such subdiffuse polarized reflectance usually involves the first 1 to 20 scattering events, depending on the tissue. Can the two-term pf offer a theoretical approach toward characterizing such subdiffuse scatter?

A TT Monte Carlo simulation was modified to record the photons escaping as reflectance $R(n)$ that have scattered once $(n=1)$, twice $(n=2)$, etc., after illumination of the surface. The simulation was run with the forward component as a Henyey-Greenstein (HG) function ($g_f=0.95$, $\alpha _f = 0.5$), which would dominate the penetration of light into the tissue and hence dominate the amount and point spread function of multiply scattered diffuse light reflectance. In a series of 8 simulations, the $g_b$ was set to values of 0.01 to 0.99. For each $g_b$, the value of $C$ was adjusted to keep the overall scattering anisotropy $g_1$ constant at 0.90. With $g_f=0.95$ and $g_1=0.90$, it was possible for a positive $g_b$ to add backscatter to the overall pf.

Figure 11 shows how varying $g_b$ influenced the early scattering events. But after about 10 scattering events, the curves all became similar as diffuse reflectance began to dominate. This example illustrates how the two-term pf allows the early subdiffuse scatter to be characterized by $g_b, \alpha _b$, while the $g_f, \alpha _f$ pair dominates the diffuse scatter, and using $C$ to keep $g_1$ constant.

 figure: Fig. 11.

Fig. 11. R(n) is the fraction of delivered photons that escape as reflectance after n scattering events from a semi-infinite tissue. Illustration of subdiffuse backscatter, controlled by $g_b$ with $\alpha _b = 1.0$ while $g_f$ = 0.95 with $\alpha _f$ = 0.50, which dominates the forward scatter and the diffuse reflectance. The values of $g_b$ were adjusted from 0.01 to 0.99 and for each $g_b$ the value of $C$ was adjusted to keep the overall scattering anisotropy $g_1$ constant at 0.90. Escaping reflected photons undergoing less than $\sim$10 scattering events are sensitive to $g_b$, and for $n>10$ $R(n)$ is insensitive to $g_b$. (a) Henyey-Greenstein, $g_1 = 0.90$, $C=1$. (Five repeated simulations, since $\alpha _b$ is ignored.) (b) $g_b$ = 0.10, $C$ adjusted. (c) $g_b$ = 0.01, $C$ adjusted. (Absorption coefficient $\mu _a$ = 1 cm$^{-1}$, scattering coefficient $\mu _s$ = 100 cm$^{-1}$, $g_1$ = 0.90.)

Download Full Size | PDF

6. Summary and discussion

6.1 Summary

The two-term phase function (TT) linearly combines two monotonic, one-term Gegenbauer-type phase functions, biased in opposite directions, to give five free parameters. The forward-directed component of the TT, $f_f(\mu )$, dominates the penetration of light into a tissue, causing multiple scattering to dominate the diffuse reflectance. The backward-directed component of the TT, $f_b(\mu )$, can dominate subdiffuse backscatter.

Compared to the two-term Henyey-Greenstein model, the two free parameters of the TT forward component provide an enhanced dynamic range for fitting the magnitude and shape of the forward scattering. The back-scattering component of the TT with its two free parameters, on the other hand, can be used to overcome the back-scattering portion of the forward-scattering component to obtain a match to early backscatter from superficial tissues.

In this paper the TT has been further developed by increasing the number of metrics that can be used to characterize the scattering properties. The new metric $h_R$ may be useful for characterizing superficial tissues, especially for experiments with subdiffuse scattering. Figure 1 provides an overview road map for the use of the TT as envisioned for forward and inverse light scattering applications with a light source incident on biological laboratory samples, and maybe eventually also for clinical applications involving tissue.

The first application of the equations illustrates how experimental data can be fit to the TT if there is difficulty in making measurements at certain angles. Then Monte Carlo results illustrate how the flux density escaping very near to a localized surface source will vary as $g_b$ is varied while keeping $g_1$ constant. Also shown is how the subdiffuse scatter involving $\le$10 scattering events is sensitive to the TT parameters $g_b$ and $\alpha _b$ while holding $g_1$ constant.

6.2 Discussion

Experimentally one way to observe the subdiffuse scattering events, and thereby possibly use the TT metrics, is to deliver light with an optical fiber while collecting the reflectance with the same fiber, but also collecting light from a nearby fiber. Both fibers collect diffuse light but the source fiber also collects the early backscatter. The difference in the signals of the two fibers isolates the subdiffuse scatter.

A second experimental method is to use a polarized-light camera that collects co-polarized and cross-polarized light. Both images contain equal amounts of diffuse reflectance, but the co-polarized image also contains early, still-polarized backscatter. The difference in images isolates the subdiffuse scatter.

A third experimental method is to use optical coherence tomography (OCT) to measure reflected light. The early signals are subdiffuse.

A fourth implementation of the TT is for goniometry, where experiments measure the angular dependence of scatter. The narrow-angle backscatter can be sensitive to superficial subdiffuse scatter. The TT was shown to augment the fitting of goniometric data for near-grazing and near-backward scattering angles, where measurements are difficult to make and tend to be less accurate as in the blood data of Fig. 9.

The TT may also be very helpful for fitting phase functions arising in non-biological applications where the illumination tends to be uniform over a surface, as in atmospheric and oceanic optics.

7. Appendix

7.1. Code for one-term (OT) and two-term (TT) phase functions

The one-term (OT) and two-term (TT) phase functions are calculated using MATLAB notation: where $\mu$ = mu = −1:1e-6:1. By reversing the order of mu to yield invmu, the fb component is pointed backwards.

7.2. Monte Carlo sampling of two-term phase function (TT)

The TT Monte Carlo program of this paper was adapted by a small change in the standard Monte Carlo simulation that is the core of the mcml.c program [38]. The TT code was validated by comparison of reflectance for simulations using Henyey-Greenstein (HG) scattering ($\alpha$ = 0.5, $C$ = 1.0) versus mcml.c simulations that also use HG.

First create the six dependent parameters that depend on the 5 independent TT parameters $g_f, \alpha _f, g_b, \alpha _b, C$ (using C-code notation):

Then in the routine Monte Carlo sampling, the following code is used, which uses the parameter C:

Then costheta and sintheta are used by the Monte Carlo simulation to specify the polar angle of photon deflection. This implementation is equivalent to using two random numbers, $\eta _f$ and $\eta _b$ to scale the $f_f$ and $f_b$ components of the TT scattering function, as stated in Eq. (28).

Acknowledgements

The reviewers’ careful critiques were much appreciated and improved the manuscript.

Disclosures

The authors declare no conflicts of interest.

Data availability

No data were generated or analyzed in the presented research.

References

1. L. G. Henyey and J. L. Greenstein, “Diffuse radiation in the galaxy,” Astrophys. J. 93, 70–83 (1941). [CrossRef]  

2. L. O. Reynolds and N. J. McCormick, “Approximate two-parameter phase function for light scattering,” J. Opt. Soc. Am. 70(10), 1206–1212 (1980). [CrossRef]  

3. V. Turzhitsky, A. Rodosevich, J. D. Rogers, A. Taflove, and V. Backman, “A predictive model of backscattering at subdiffusion length scales,” Biomed. Opt. Express 1(3), 1034–1046 (2010). [CrossRef]  

4. A. N. Yaroslavsky, I. V. Yaroslavsky, T. Goldbach, and H.-J. Schwarzmaier, “Influence of the scattering function approximation on the optical properties of blood determined from the integrating sphere measurements,” J. Biomed. Opt. 4(1), 47–53 (1999). [CrossRef]  

5. G. Mie, “Beiträge zur Optik der über Medien, speziell kolloidaler Metallösungen,” Ann. Phys. 330(3), 377–445 (1908). [CrossRef]  

6. C. D. Mobley, Light and Water: Radiative Transfer in Natural Waters, p. 438 (Academic, (1994).

7. C. D. Mobley, L. K. Sundman, and E. Boss, “Phase function effects on oceanic light fields,” Appl. Opt. 41(6), 1035–1050 (2002). [CrossRef]  

8. C. Mobley, E. Boss, and C. Roeler, “Ocean Optics Web Book,” (2020), (see Radiative Transfer Theory, level 2), [retrieved Nov., 2022] https://www.oceanopticsbook.info.

9. S. Leyre, Y. Meuret, G. Durinck, J. Hofkens, G. Deconinck, and P. Hanselaer, “Estimation of the effective phase function of diffusing materials with the inverse adding-doubling method,” Appl. Opt. 53(10), 2117–2125 (2014). [CrossRef]  

10. N. Bodenschatz, P. Krauter, A. Liemert, and A. Kienle, “Quantifying phase function influence in subdiffusively backscattered light,” J. Biomed. Opt. 21(3), 035002 (2016). [CrossRef]  

11. P. Naglič, F. Pernuš, B. Likar, and M. Bürmen, “Estimation of optical properties by spatially resolved reflectance spectroscopy in the subdiffusive regime,” J. Biomed. Opt. 21(9), 095003 (2016). [CrossRef]  

12. M. Ivančič, P. Naglič, F. Pernuš, B. Likar, and M. Bürmen, “Virtually increased acceptance angle for efficient estimation of spatially resolved reflectance in the subdiffusive regime: a Monte Carlo study,” Biomed. Opt. Express 8(11), 4872–4886 (2017). [CrossRef]  

13. A. N. Yaroslavsky, I. V. Yaroslavsky, T. Goldbach, and H-J. Schwarzmaier, “The optical properties of blood in the near infrared spectral range,” Proc. SPIE 2678, 314–324 (1996). [CrossRef]  

14. A. N. Yaroslavsky, I. V. Yaroslavsky, T. Goldbach, and H-J. Schwarzmaier, “Different phase function approximations to determine optical properties of blood: A comparison,” Proc. SPIE 2982, 324–330 (1997). [CrossRef]  

15. M. Friebel, A. Roggan, G. Müller, and M. Meinke, “Determination of optical properties of human blood in the spectral range 250 to 1100 nm using Monte Carlo simulations with hematocrit-dependent effective scattering phase functions,” J. Biomed. Opt. 11(3), 034021 (2006). [CrossRef]  

16. A. L. Post, S. L. Jacques, H. J. C. M. Sterenborg, et al., “Modeling subdiffusive light scattering by incorporating the tissue phase function and detector numerical aperture,” J. Biomed. Opt. 22(5), 050501 (2017) and Erratum in 24(6), 069801-1 (2019). [CrossRef]  

17. J. D. Rogers, A. J. Radosevich, J. Yi, and V. Backman, “Modeling light scattering in tissue as continuous random media using a versatile refractive index correlation function,” IEEE J. Select. Topics Quantum Electron. 20(2), 173–186 (2014). [CrossRef]  

18. V. V. Marinyuk, D. B. Rogozkin, and S. V. Sheberstov, “Propagation of a light beam in an absorbing medium with large-scale inhomogeneities,” Opt. Spectrosc. 117(1), 102–110 (2014). [CrossRef]  

19. V. V. Marinyuk and S. V. Sheberstov, “Impact of the scattering phase function on the reflectance of a turbid medium with large-scale inhomogeneities,” Appl. Opt. 56(32), 9105–9113 (2017). [CrossRef]  

20. V. V. Marinyuk and S. V. Sheberstov, “Effect of the single-scattering phase function on light transmission through disordered media with large inhomogeneities,” J. Phys. 788, 012045 (2017). [CrossRef]  

21. V. V. Marinyuk and S. V. Sheberstov, “Total reflectance of a random medium with the Reynolds-McCormick phase function at grazing incidence of light,” J. Phys.: Conf. Ser. 1238(1), 012043 (2019) IOP Publishing. [CrossRef]  

22. V. V. Marinyuk, V. S. Remizovich, and S. V. Sheberstov, “Angular reflectance of a highly forward scattering medium at grazing incidence of light,” J. Opt. Soc. Am. A 37(3), 501–510 (2020). [CrossRef]  

23. G. W. Kattawar, “A three-parameter analytic phase function for multiple scattering calculations,” J. Quant. Spectrosc.Radiat. Transfer 15(9), 839–849 (1975). [CrossRef]  

24. K. W. Calabro and W. Cassarly, “Modeling scattering in turbid media using the Gegenbauer phase function,” Proc. SPIE 9333, 93330F (2015). [CrossRef]  

25. M. Baes, P. Camps, and A. U. Kapoor, “A new analytical scattering phase function for intersteller dust,” Astron. & Astrophys. 659, A149 (2022). [CrossRef]  

26. J. Wang, C. Xu, A. M. Nilsson, D. L. A. Fernandes, and G. A. Niklasson, “A novel phase function describing light scattering of layers containing colloidal nanospheres,” Nanoscale 11(15), 7404–7413 (2019). [CrossRef]  

27. J. Wang, C. Xu, H. Xiong, Y. Han, A. M. Nilsson, M. Strömberg, T. Edvinsson, and G. A. Niklasson, “Extraction of backscattering and absorption coefficients of magnetite nanosphere composites from light-scattering measurements: Implications for optomagnetic sensing,” ACS Appl. Nano Mater. 3(11), 11172–11183 (2020). [CrossRef]  

28. T. Harmel, J. Agagliate, M. Hieronymi, and P. Gernez, “Two-term Reynolds-McCormick phase function parameterization better describes light scattering by microalgae and mineral hydrosols,” Opt. Lett. 46(8), 1860–1863 (2021). [CrossRef]  

29. L. Ma, L. Hu, C. Jia, C. Wang, and L. Liu, “Quantitative evaluation of the phase function effects on light scattering and radiative transfer in dispersed systems,” Photonics 9(8), 584 (2022). [CrossRef]  

30. G. Fournier and J. L. Forand, “Analytic phase function for ocean water,” Proc. SPIE 2258, 194–201 (1994). [CrossRef]  

31. G. Fournier and M. Jonas, “Computer-based underwater imaging analysis,” Proc. SPIE 3761, 62–70 (1999). [CrossRef]  

32. M. Keijzer, S. L. Jacques, S. A. Prahl, and A. J. Welch, “Light distribution in artery tissue: Monte Carlo simulations for finite diameter laser beams,” Lasers Surg. Med. 9(2), 148–154 (1989). [CrossRef]  

33. F. Bevilacqua and C. Depeursinge, “Monte Carlo study of diffuse reflectance at source-detector separations close to one transport mean free path,” J. Opt. Soc. Am. A 16(12), 2935–2945 (1999). [CrossRef]  

34. N. Bodenschatz, P. Krauter, A. Liemert, J. Wiest, and A. Kienle, “Model-based analysis on the influence of spatial frequency selection in spatial frequency domain imaging,” Appl. Opt. 54(22), 6725–6731 (2015). [CrossRef]  

35. H. C. van de Hulst, Multiple Light Scattering, Vol. 2, p. 488 (Academic1980).

36. J. R. Mourant, M. Canpolat, C. Brocker, O. Espona-Ramos, T. M. Johnson, A. Matanock, K. Stetter, and J. P. Freyer, “Light scattering from cells: the contribution of the nucleus and the effects of proliferative status,” J. Biomed. Opt. 5(2), 131 (2000). [CrossRef]  

37. S. L. Jacques, J. C. Ramella-Roman, and K. Lee, “Imaging skin pathology with polarized light,” J. Biomed. Opt. 7(3), 329–340 (2002). [CrossRef]  

38. L. Wang, S. L. Jacques, and L. Zheng, “MCML-Monte Carlo modeling of light transport in multi-layered tissues,” Methods Programs Biomed. 47(2), 131–146 (1995). [CrossRef]  

Data availability

No data were generated or analyzed in the presented research.

Cited By

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

Alert me when this article is cited.


Figures (11)

Fig. 1.
Fig. 1. Overview of how a math description of light scattering by a tissue can lead to categorization of the tissue type, which then aids medical decisions. Additionally, the math description can be used in light transport theory to simulate experiments and generate simulated measurements, which then can contribute to the analysis. Hence knowledge of multiple scattering behavior in both the subdiffuse and diffuse regimes can participate in characterizing a tissue. The dashed line indicates how a simulation can augment experimental data by simulating missing data. In red, a tissue site on a patient or biomaterial of interest can be characterized by the analysis module to yield an assessment of the tissue’s nano- and micro-scale structural status to inform clinical or industrial decisions.
Fig. 2.
Fig. 2. One-term phase function. As the parameter $\alpha$ increases, the backscatter falls exponentially while the forward increases linearly. (A) The scattering function $f(\mu )$ is shown as $\alpha$ is varied over 0.5, 1, 2, and 3, while $g$ is constant at 0.90. (B) Close-up view of (A) to show forward scatter in the 0-8$^{\circ }$ range. (C) Backscatter at 180$^{\circ }$ drops exponentially vs $\alpha$. (D) Forward scatter at 0$^{\circ }$ increases linearly vs $\alpha$.
Fig. 3.
Fig. 3. (A) Anisotropy $g_1$ versus $g$ for a range of $\alpha$ values. (B) $g_1$ versus $\alpha$ for a range of $g$ values
Fig. 4.
Fig. 4. (A) Backscatter fraction $b$ versus $g$ for a range of $\alpha$ values. (B) $b$ versus $\alpha$ for a range of $g$ values
Fig. 5.
Fig. 5. (A) Fraction of backward scattering that projects onto the backward direction, metric $h_R$, plotted versus $g$ for a range of $\alpha$ values. (B) $h_R$ plotted versus $\alpha$ for a range of $g$ values. (C) Fraction of forward scatter that projects onto the forward direction, metric $h_T$, plotted versus $g$ for a range of $\alpha$ values. (D) $h_T$ plotted versus $\alpha$ for a range of $g$ values.
Fig. 6.
Fig. 6. Polar plots emphasizing the backward scatter of the forward and backward components of the two-term scattering function. (A) Forward scatter component. $\alpha _f$ is varied (0.001, 0.1, 0.3, 0.5, 0.7, 0.9, 0.95) and $g_f$ = 0.90. (B) Backward scatter component. $\alpha _b$ is similarly varied and $g_f$ = 0.90. Black curve indicates Henyey-Greenstein scattering ($\alpha _f$ or $\alpha _b$ = 0.5). Other colors vary $\alpha _f$ and $\alpha _b$.
Fig. 7.
Fig. 7. Two-term $f(\mu )$ for $\alpha _f = 1,\ g_f = 0.8$, and $C = 0.98$. The forward scatter $f_f(\mu )$ (black line), backward scatter $f_b(\mu )$ (red line), and the total scatter $f(\mu )$ (green line) are shown. Note the backward scattering of $f_f(\mu )$ and the forward scattering of $f_b(\mu )$ both contribute to the total backscattering of $f(\mu )$. (A) $g_b$ = 0.1 yields a broad backscatter, (B) $g_b$ = 0.7 yields a more directional backscatter, and (C) $g_b$ = 0.9 yields a sharply peaked backscatter.
Fig. 8.
Fig. 8. Scattering by cells (Mourant et al. [36]). Data is red line. (A) Two-term Henyey-Greenstein (TTHG), dashed line, does not fit the data. (B) Two-term (TT) fits the data.
Fig. 9.
Fig. 9. Scattering function ($f(\mu )$ [sr$^{-1}$] from a 100-$\mu$m-thick layer of blood Yaroslavsky et al. [4]. The red circles are the data fit by the TT, and the black x’s are data ignored by the fitting. The blue and green lines are the contributions of the forward and backward components of the TT, respectively. The red line is the total TT, and predicts the backscatter behavior in the angular range near $\mu$ = −1 where data are difficult to measure.
Fig. 10.
Fig. 10. (A) Monte Carlo simulation of the point spread function of light within an example tissue that has a scattering anisotropy of $g_1 = 0.9$. (B) The lateral profile of escaping flux density $R(r)$ [cm$^{-2}$]. (C) R $\leq$ 0.6 mm equals the integral of $R(r)$ over $r=0$ to 0.06 cm, which has a contribution from subdiffuse superficially scattered light and a contribution from diffuse multiply scattered light. As the backscatter component becomes increasing isotropic ($g_b$ approaching 0), the magnitude of backscatter reaches a maximum when $g_1$ is constant.
Fig. 11.
Fig. 11. R(n) is the fraction of delivered photons that escape as reflectance after n scattering events from a semi-infinite tissue. Illustration of subdiffuse backscatter, controlled by $g_b$ with $\alpha _b = 1.0$ while $g_f$ = 0.95 with $\alpha _f$ = 0.50, which dominates the forward scatter and the diffuse reflectance. The values of $g_b$ were adjusted from 0.01 to 0.99 and for each $g_b$ the value of $C$ was adjusted to keep the overall scattering anisotropy $g_1$ constant at 0.90. Escaping reflected photons undergoing less than $\sim$10 scattering events are sensitive to $g_b$, and for $n>10$ $R(n)$ is insensitive to $g_b$. (a) Henyey-Greenstein, $g_1 = 0.90$, $C=1$. (Five repeated simulations, since $\alpha _b$ is ignored.) (b) $g_b$ = 0.10, $C$ adjusted. (c) $g_b$ = 0.01, $C$ adjusted. (Absorption coefficient $\mu _a$ = 1 cm$^{-1}$, scattering coefficient $\mu _s$ = 100 cm$^{-1}$, $g_1$ = 0.90.)

Equations (38)

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

2 π 1 1 f ( μ ) d μ = 1
( 1 + g 2 2 g μ ) 1 / 2 = n = 0 P n ( μ ) g n ,
( 4 π ) 1 ( 1 g 2 ) ( 1 + g 2 2 g μ ) 3 / 2 = ( 4 π ) 1 n = 0 ( 2 n + 1 ) P n ( μ ) g n ,
f ( μ ) = K ( 1 + g 2 2 g μ ) ( α + 1 ) , α > 1 / 2 ,   1 g 1 ,
K = π 1 α g ( 1 g 2 ) 2 α [ ( 1 + g ) 2 α ( 1 g ) 2 α ] 1 .
η = 2 π 1 μ f ( μ ) d μ
μ = 1 + g 2 2 g 1 2 g [ η ( 1 g ) 2 α + 1 η ( 1 + g ) 2 α ] 1 / α ,   1 μ 1
μ = 1 + g 2 2 g 1 2 g ϱ α
ϱ = 2 α g η K + 1 ( 1 + g ) 2 α
μ = 1 + g 2 2 g 1 2 g [ 1 g 2 1 g + 2 g η ] 2 ,   1 μ 1 .
g n = 2 π 1 1 P n ( μ ) f ( μ ) d μ 2 π 1 1 f ( μ ) d μ , n = 0 ,   1 ,   2 , ,
g 1 = 2 g α L ( 1 + g 2 ) 2 g ( α 1 ) , α 1
L = ( 1 + g ) 2 α + ( 1 g ) 2 α ( 1 + g ) 2 α ( 1 g ) 2 α .
g 1 = 1 + g 2 2 g + ( 1 g 2 2 g ) 2 ln ( 1 g 1 + g ) .
μ s = μ s ( 1 g 1 )
g 2 = 3 ( 1 + g 2 ) g 1 2 g ( 2 α ) 1 + α 2 α , α 2
g 3 = 5 2 ( 3 α ) [ 1 + g 2 g g 2 α L + 1 + g 2 2 g ] 3 2 g 1 , α 3
γ = 1 g 2 1 g 1
δ = 1 g 3 1 g 1 .
b = 2 π 1 0 f ( μ ) d μ 2 π 1 1 f ( μ ) d μ
b = ( 1 g ) 2 α ( 1 + g ) 2 α ( 1 g ) 2 α [ ( 1 + g ) 2 α ( 1 + g 2 ) α 1 ]
b = U 1 V 1
U = ( 1 + g ) 2 α ( 1 + g 2 ) α a n d V = ( 1 + g ) 2 α ( 1 g ) 2 α .
b = 1 g 2 g [ 1 + g ( 1 + g 2 ) 1 / 2 1 ] .
b n = ( 1 ) n 2 π 1 0 P n ( μ ) f ( μ ) d μ 2 π 1 1 f ( μ ) d μ , n = 0 ,   1 ,   ,
b 1 = 2 π K ( 1 + g 2 ) ( 2 g ) 2 α ( α 1 ) [ 1 ( 1 + g 2 ) α 1 ( 1 + g ) 2 α ] , α 1 .
h R = 2 π 1 0 μ f ( μ ) d μ 2 π 1 0 f ( μ ) d μ b 1 / b 0 .
f n = 2 π 0 1 P n ( μ ) f ( μ ) d μ 2 π 1 1 f ( μ ) d μ , n = 0 ,   1 , ,
h T = 2 π 0 1 μ f ( μ ) d μ 2 π 0 1 f ( μ ) d μ f 1 / f 0 ,
f ( μ ) = C b f b ( μ ) + C f f f ( μ ) = i = { b & f } C i f i ( μ )
μ = C b μ b + C f μ f , 1 μ 1 = i = { b & f } C i μ i
g n = 2 π 1 1 P n ( μ ) f ( μ ) d μ 2 π 1 1 f ( μ ) d μ , n = 1 ,   2 ,     .
g n = i = { b & f } C i g n , i .
μ s = μ s [ 1 i = { b & f } C i g 1 , i ] .
b n = ( 1 ) n 2 π 1 0 P n ( μ ) f ( μ ) d μ 2 π 1 1 f ( μ ) d μ , n = 0 , 1 ,   2 ,   = ( 1 ) n i = { b & f } C i 2 π 1 0 P n ( μ ) f i ( μ ) d μ = i = { b & f } C i b n , i   ,
h R = 2 π 1 0 μ f ( μ ) d μ 2 π 1 0 f ( μ ) d μ = i = { b & f } C i b 1 , i b i [ i = { b & f } C i b i ] 1 .
f n = 2 π 0 1 P n ( μ ) f ( μ ) d μ 2 π 1 1 f ( μ ) d μ , n = 0 ,   1 ,   2 ,   = i = { b & f } C i f n , i ,
h T = 2 π 0 1 μ f ( μ ) d μ 2 π 0 1 f ( μ ) d μ = i = { b & f } C i f 1 , i f i [ i = { b & f } C i f i ] 1 .
Select as filters


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