1. Introduction

Recently, there has been a particularly strong interest in understanding the properties of diffuse backscattering Mueller matrix from scattering media. The Mueller matrix is well known to describe many optical elements and materials [1–3]. While for a highly scattering medium, any information about particle size, refractive index, particle shape etc. can be extracted from the various matrix elements because the description of its optical characteristics with the Mueller matrix is complete. The Mueller matrix is considered to be an “optical fingerprint” of the scattering medium under investigation.

Experimentally, the determination of the Mueller matrix is based on the measurement of the polarization states of the scattered light for the specified polarization state of the incident light, and subsequent calculation of the corresponding Mueller matrix components. By varying the polarization of the source and 49 independent intensity measurements, Heilscher et al. used the steady state backscattering Mueller matrix patterns to differentiate cancerous from noncancerous cell suspensions [4]. Jiao et al. measured the Mueller matrix of hard or soft biological tissues by a double-beam polarization-sensitive optical coherence tomography imaging technique, and extracted the polarization parameters, such as magnitude and orientation of birefringence and diattenuation [5–6]. Manhas et al. used the measured backscattering Mueller matrix to determine the optical rotation in chirl turbid media by a polar decomposition approach [7]. Berezhnyy et al. presented a complete spatial-temporal polarization pattern description of a Mueller matrix of scattering media and the polarization patterns were used to separate the polarimetric contributions of different scattering paths [8]. In addition, Justin et al. reported a novel optical polarization-imaging system which can be used to generate the full 16-element Mueller matrix in less than 70-sec [9].

Theoretically, the Mueller matrix of a scattering medium can be obtained by solving the Maxwell’s equations with the corresponding boundary conditions. However, this requires stipulating the appropriate assumptions concerning the scattering process as well as the properties of the medium, and involves intense numerical calculations. Therefore, directly solving Maxwell’s equations for the general case of multiple scattering media appears to be a computationally unrealistic task. To characterize the diffuse backscattering Mueller matrix from scattering media, a mathematical model is required to interpret data from a multiple scattering Mueller matrix data in a highly scattering medium. Hence, several research groups have developed numerical and approximately analytical methods that characterize the diffuse backscattering Mueller matrix from scattering media. Bartel et al., Yao et al., and Kaplan et al. developed a Monte Carlo model based on the Stokes-Mueller formulation for a multiple scattering Mueller matrix [11–13]. However, Monte Carlo simulation has limitations of computational efficiency. Wang et al. presented a single-scattering model and compared it with a Monte Carlo model. Agreement was satisfactorily achieved between the two methods in the non-diffusing regime, but the single-scattering model becomes inaccurate in the diffuse regime [14]. Raković et al. was the first to present a double scattering method to explain the diffuse backscattering Mueller matrix in turbid media. They were able to reproduce the azimuthal symmetry in all 16 matrix elements, but obtained only poor agreement in the radial dependence [15–17]. Moreover, use of the double scattering approximation method has weakly accurate the revealing of the polarized patterns dependent on scattering coefficients [18].

In this paper we completely characterize the spatially-resolved diffuse backscattering Mueller matrix from a highly scattering medium with an effective Mueller matrix (EMM), which is based on triple scattering approximation method reported by author in reference [18]. We numerically calculate the matrix expressions for the polystyrene sphere suspensions and obtain the two-dimensional intensity pattern maps of each matrix element. The calculated results are compared with the Monte Carlo simulation results performed by the authors using a code developed by Ramella-Roman [13].

In Section 2, we provide a description of the effective Mueller matrix based on triple scattering approximation. In Section 3, the calculated results for sphere suspensions are compared with Monte Carlo simulated results.

2. Theory

Each light beam can be represented by its four-component Stokes vector, S=(S₀, S₁, S₂, S₃)^T, which is a mathematical representation of the polarization state, where S ₀=〈|E_x|²+|E_y|²〉, S ₁=〈|E_x|²-|E_y|²〉, S ₂=〈2E _x E_ycosδ〉, S ₃=〈2E_xE_ysinδ〉, E_x and E_y are two orthogonal complex electric field components perpendicular to the propagation direction, δ denotes the phase difference with respect to Ex and Ey. And the brackets <> refers to ensemble average. Some of the Stokes components have an obvious interpretation. For example, S ₀ is equal to the total intensity.

Let S_in be the Stokes vector of the incident laser beam that is injected into a highly scattering medium surface, and S_out be the Stokes vector of diffuse backscattering light that escapes from the upper surface of the medium. As light is described by a four-component vector, this interaction with the scattering medium can be described as a multiplication of the Stokes vector with a 4x4 matrix: S_out=M S_in. This sixteen-element matrix is called the Mueller matrix of the scattering medium, as follows:

This 4×4 matrix operator gives a comprehensive description of scattering properties of a medium in the spatial domain. Rather than being just one number, each of the 16 components of the Mueller matrix is, in fact, a two-dimensional array of numbers, corresponding to different spatial locations across the surface of the medium. On the other hand, the matrix is only determined by the properties of the scattering medium and optical paths. It allows determining various characteristics of the medium such as an average size of the scattering particles, particle shape, scatterer spatial distribution, refractive index of the medium, scattering coefficient and anisotropy factor, etc.

In general, the Mueller matrix intensity patterns are formed as a result of mixing different orders of scattering. While for highly scattering media, the diffuse backscattering intensity patterns are essentially formed as a result of multiple scattering of light, which can be described by triple scattering approximate method. The formulations for triple scattering approximation are based on the assumption that the scattering medium is a slab structure and the scattering of polarized light is incoherent and the successive scattering events occur at the lower half-space of the medium. Considering that the major contribution to the patterns comes from the triple scattering of light, thus, an effective backscattering Mueller matrix from a scattering medium can takes the following form [18]:

with

$C_1(\theta ,\rho )=-{\mu }_t\rho \left[\cos \theta +\sin \theta -2\sin \theta \sin (\theta ⁄2)\right],C_2\left(\theta \right)=\cos \theta {\sin }^2\theta \sin (\theta ⁄2)$

$G\left(\theta \right)=M(\pi ⁄4-\theta ⁄2)M(3\pi ⁄4-\theta ⁄2)M\left(\theta \right)$

where i≤4 and j≤4 denote the matrix row and column, respectively; µ_t and µ_s are the extinction coefficient and the scattering coefficient, respectively; M (θ) is a single scattering matrix [2,19]. In this study we consider only isotropically distributed spherical particles, which reduces the Mueller matrix to four independent elements:

Mie theory is used to determine the scattering amplitudes from which the m_ij are derived [19]. R (Φ) is a 4×4 rotation matrix, which transforms between the scattering plane and the reference plane; Φ is defined as positive for clockwise rotation [2]:

Using Eq. (2), an effective Mueller matrix based on the triple scattering approximation, we numerically calculate the intensity pattern maps for the 16 matrix elements in polar coordinates and compare the calculated results with Monte Carlo simulated results. In the meantime, the results calculated from the double scattering approximation are also provided.

3. Results

The numerical calculation results based on scattering approximate model and Monte Carlo simulation results are compared for the case of the backscattering Mueller matrix from aqueous suspensions of polystyrene spheres. All 16 matrix elements displayed are normalized to the maximum intensity of the M11element that the amplitudes range from -1 to +1. And each image displayed is 2cm×2cm.

Within the approximate model calculations, The particle size, the refractive index of polystyrene sphere and water, and the particle concentration and wavelength of the laser are chosen identical to Monte Carlo simulating parameters, as follows: indices of refraction of the medium n₀=1.33 and of the polystyrene spheres scatterers n=1.59 at a wavelength of λ=0.6328µm. And Mie theory can be used to calculate the scattering cross section C_s of the individual particles. Given the number density N_s, that is, the particle concentration, the scattering coefficients of turbid media can be determined as µ_s=N_sC_s. For statistical results, large quantities of photon packets (10⁸) are traced in the Monte Carlo simulations. The thickness of the aqueous suspensions of polystyrene spheres is 2.5cm. Firstly, for the aqueous suspension of polystyrene spheres of particle diameter 0.486 µm with scattering coefficients 14.8 cm^-1, the numerical calculations of triple scattering model are given as two-dimensional images of the surface in Fig. 1(a), and the calculations of double scattering model are shown in Fig. 1(b). As a comparison, Monte Carlo simulated results for the same aqueous suspension of polystyrene spheres are also shown in Fig. 1(c).

 

Fig. 1. Comparison of the calculations of the effective Mueller matrix and Monte Carlo simulations for the suspension of 0.486 µm particles with scattering coefficients 14.8 cm^-1. (a) The calculations of triple scattering approximation. (b) The calculations of double scattering approximation. (c) Monte Carlo simulations.

Download Full Size | PDF

As can be seen in Fig. 1(a) and Fig. 1(b), the elements M14 and M41 are always zero. The elements M11 and M44 are independent of the azimuth, whereas the other elements have different azimuthal intensity variations. And several pairs of elements exist that differ only in the orientation of the intensity pattern, for example M13 can be obtained by rotating M12 by 45 degrees, the same relation holds for the pairs M31 and M21, and M22 and M33. Obviously, seven elements are only independent in the effective Mueller matrix [16]. Comparing with Monte Carlo simulated results in Fig. 1(c), the symmetric structure properties of patterns calculated by the triple and double scattering methods for the individual matrix element and the matrix as a whole are similar to the corresponding simulated diffuse backscattering Mueller matrix. And the azimuthal intensity variations in each row and each column also show good overall agreement with Monte Carlo simulated results. It can be concluded that the effective Mueller matrix for the double and triple scattering models can recover the basic symmetries and azimuthal structures of the diffuse backscattering Mueller matrix from a highly medium.

However, it is interesting to note that the matrix elements calculated by the two approximation methods shown Fig. 1(a) and Fig. 1(b) have obviously dissimilar radial extent. Comparing with Monte Carlo simulations shown Fig. 1(c), the calculations of the effective Mueller matrix based on triple scattering shown in Fig. 1(a) appear better agreement with the Monte Carlo simulations than that of the double scattering shown in Fig. 1(b). The matrix elements calculated by the triple scattering model are further compared quantitatively in radial dependence. Figure 2 shows the M44 element radial dependence along the horizontal direction in M44 image shown in Fig. 1(a). From Fig. 2, it is shown that, the values of the M44 element along the horizontal direction are always negative. M44 increases with increasing the radial distance. The numerical calculated results (line) reflect the same trend in radial dependence with the Monte Carlo simulations (dot). Similarly, Figure 3 shows the M21 element radial dependence along the horizontal direction in M21 image shown in Fig. 1a. It also can be seen that the calculations of the effective Mueller matrix based on triple scattering approximation (line) also have the same variations with Monte Carlo simulations (dot). It should be pointed out that, similar to the M44 and M21 elements, the other effective Mueller matrix elements also show the same symmetry relationships and radial dependence as Monte Carlo simulations. These comparisons demonstrate the effective Mueller matrix can not only reproduce the azimuthal symmetry in all 16 matrix elements, but also have a good agreement with Monte Carlo simulated results in the radial dependence.

 

Fig. 2. Comparison of the calculations of the effective Mueller matrix element M44 (line) and Monte Carlo simulations (dot) dependence on radial instance for the suspension of 0.486 µm particles with scattering coefficients 14.8 cm^-1.

Download Full Size | PDF

 

Fig. 3. Comparison of the calculations of the effective Mueller matrix element M21 (line) and Monte Carlo simulations (dot) dependence on radial instance for the suspension of 0.486 µm particles with scattering coefficients 14.8 cm^-1.

Download Full Size | PDF

Further, by varying the scatterer size, we calculate the effective Mueller matrix and simulate corresponding diffuse backscattering Mueller matrix with Monte Carlo method and compare the characteristics of patterns variations. For the elements dependent on azimuth, as an example, we observe the pattern variations of the element M12 for three sphere suspensions with different scatterer sizes.

 

Fig. 4. Comparison of the calculations of the effective Mueller matrix element M12 and Monte Carlo simulations for the suspensions of different particle diameter with scattering coefficients 14.8 cm^-1.

Download Full Size | PDF

 

Fig. 5. Comparison of the calculations of the effective Mueller matrix element M44 (line) and Monte Carlo simulations (dot) as a function of the particle diameter for the suspensions of different particle diameter with scattering coefficients 14.8 cm^-1

Download Full Size | PDF

Figure 4 shows the results. It can be seen that the patterns of the element M12 expand or shrink around the point of incidence of the laser beam as the scatterer sizes vary; this feature is reproduced by Monte Carlo simulations. For the elements independent on azimuth, such as M44 element, we extract the mean value of M44 elements as a function of scaterer sizes. Figure 5 shows the effect of scaterer sizes on the element M44 mean value for five suspensions of scattering coefficients scattering coefficient 14.8cm^-1. We found that, as the scaterer sizes increase, the mean values of M44 increase. The calculations have a good agreement with Monte Carlo simulated results. The results demonstrate that the pattern variations of the effective Mueller matrix can accurately reveal the variations of the scaterer sizes. Furthermore, notwithstanding the changes in scattering coefficients or anisotropy factor g of the suspensions, the pattern features of the effective Mueller matrix coincide with Monte Carlo simulations well.

Overall the numerical calculations of the effective Mueller matrix agree closely with Monte Carlo simulated results, despite some small differences. The discrepancies may be also due to the effective Mueller matrix in which the contribution of lower-order and high-order influences on scattered light are neglected. We will further apply this framework to study the contribution of different orders to the diffuse backscattering Mueller matrix and validity in future works.

4. Conclusions

By using an effective Mueller matrix based on triple scattering approximation, we characterize the patterns of diffuse backscattering Mueller matrix from a highly scattering medium. Both the calculations of the effective Mueller matrix and Monte Carlo simulations have excellent agreements in azimuthal symmetry relations and radial dependence for different aqueous suspensions of polystyrene spheres. And as the scatterer size varies, the effective Mueller matrix pattern variations have same trend with Monte Carlo simulations.

It has been demonstrated that the effective Mueller matrix can be used to quantitatively predict and interpret an experimentally-determined backscattering Mueller matrix from a highly scattering medium. The method provides a tool to investigate the effects of different particle sizes and optical properties on observable diffuse backscattering Mueller matrixes from highly scattering media, and it is potentially possible to be used in biomedical diagnostics such as blood cell diagnosis.