2.1.

SCBM for Anisotropic Media

As shown in Fig. 1, SCBM approximates the anisotropic turbid medium to a mixture of solid spherical and infinitely long cylindrical scatterers embedded in a linearly birefringent medium. In this work, we assume that the surrounding medium is a uniaxial material with its extraordinary axis being along the direction of the cylinders.⁹ Parameters for the scatterers are the same as those for the SCSM.¹¹ They include number densities and sizes of both the spheres and cylinders, and the mean value and standard deviation of the angular distribution for the cylinders. Parameters for the surrounding medium include the value and axis direction of birefringence. Both the cylindrical orientation and the axis direction of birefringence can be adjusted in the three-dimensional (3-D) space. SCBM reduces to SCSM if the birefringence is set to zero, or to SBM if the scattering coefficient of the cylinders is set to zero.

Fig. 1

Schematics of the sphere–cylinder birefringence model (SCBM).

2.2.

Monte Carlo Algorithm

The Monte Carlo program for SCBM was developed based on the polarization-sensitive SCSM program,¹² which tracks the trajectory and polarization state of photons. In the SCSM program, after each scattering event, the photon moves along a new direction determined by the Mie theory^13,14 and may lose part of its energy because of the absorption. A statistical method was designed to determine whether the photon is scattered by the cylinders or the spherical particles. The simulation process continues until the photon is completely absorbed or moves out of detection range. For SCBM, as the polarized photons transmit in the anisotropic medium, they alternately experience the transmission process in the birefringent medium and the scattering process by the spherical or infinitely long cylindrical scatterers.

The Monte Carlo program for SCBM was also developed based on the polarization-sensitive SCSM program¹² and the SBM program,⁹ which describes the detailed algorithm of the birefringence module. In SCBM, we assumed that the birefringence effect is the property of the surrounding medium and does not affect the scattering phase function, but the birefringence does alter the polarization states of the photons as they propagate between two successive scattering events.^9,15 We added a new module to the SCSM Monte Carlo program to calculate the effects of birefringence.

In SCSM, the polarization states of the photons are altered by scattering only. If the medium around the scatterers is birefringent, a scattered photon undergoes retardation between successive scatterings. The new module calculates the changes of Stokes vector due to the birefringent medium:

2.3.

Validity Test: Comparison with SBM

The SCBM becomes SBM if the scattering coefficient of cylinders is set to zero. To test the validity of the new module, we set the scattering coefficient of the cylinders to zero in the SCBM program and compared the simulated Mueller matrix patterns with previously published results of SBM.¹⁶ Figure 2 shows the simulation results of both models using the same set of parameters for the spheres and birefringent medium.¹⁶ For the spheres, the diameter is 0.46 μm and the scattering coefficient is $10{\mathrm{cm}}^{-1}$. The average refractive index of the birefringent surrounding medium is 1.34, and its extraordinary axis is along the horizontal $x$ axis. The maximum retardation at the transport length distance is calculated by ${\mathrm{\Delta }}_{\mathrm{TR}}=2\pi \mathrm{\Delta }n{l}_{\mathrm{TR}}/{\lambda }^{\text{'}}=2\pi \mathrm{\Delta }n/[{\lambda }^{\text{'}}{\mu }_s(1-g)]$.¹⁶ In Fig. 2(a)–2(d), the birefringence value takes 0, $1.5\times {10}^{-6}$, $7.5\times {10}^{-6}$, and $4.5\times {10}^{-5}$, corresponding to ${\mathrm{\Delta }}_{\mathrm{TR}}=0$, 0.1, 0.5, and 3, respectively. Figure 2 shows a good agreement between the simulation results of SCBM and SBM.

Fig. 2

Monte Carlo–simulated backscattering Mueller matrices of the SBM: ${\mathrm{\Delta }}_{\mathrm{TR}}=0$ (a); ${\mathrm{\Delta }}_{\mathrm{TR}}=0.1$ (b); ${\mathrm{\Delta }}_{\mathrm{TR}}=0.5$ (c); and ${\mathrm{\Delta }}_{\mathrm{TR}}=3$ (d). Each image size is $1\times 1\mathrm{cm}$. Note that the $m11$ images used the same color maps but with the range from 0 to 1.

2.4.

Validity Test: Comparison with Mueller Matrix Experiments

2.4.1.

Experiment setup

Figure 3 shows a typical experimental setup for backscattering Mueller matrix measurement.¹⁰ A polarized He-Ne laser ($\lambda =633\mathrm{nm}$) is used as the light source. A quarter-wave plate (QW1) and a linear polarizer P1 control the incident polarization. The incident light reaches the sample through a small hole at the center of the mirror ($M$), which is tilted at 45 deg to the incident beam. The mirror directs the backscattered light from the sample toward the camera through a quarter-wave plate (QW2) and a linear polarizer (P2). The backscattered light is imaged by a 14-bit CCD camera (Canon APS-C). A laboratory reference frame is defined such that both the muscle fiber orientation and the vertical linear polarization direction are along the $y$ axis.

Fig. 3

Scheme of experimental setup. Light source (LS): 12-mW He-Ne laser; $P$, polarizer; QW, quarter-wave plate; CCD, imaging camera.

In the experiments, six different polarization states were used for both incidence and detection: horizontal linear ($H$), vertical linear ($V$), 45-deg linear ($P$), 135-deg linear ($M$), right circular ($R$), and left circular ($L$). A total of 36 reflectance images were captured with different combinations of the six incidence and detection polarization states. The Mueller matrix was calculated from these images:

Eq. (7)

$$\begin{gathered} M=\left[\begin{array}{cccc}m11& m12& m13& m14\\ m21& m22& m23& m24\\ m31& m32& m33& m34\\ m41& m42& m43& m44\end{array}\right] \\ =\frac{1}{2}\left[\begin{array}{ccc}HH+HV+VH+VV& HH+HV-VH-VV& PH+PV-MP-MM\\ HH-HV+VH-VV& HH-HV-VH+VV& PH-PV-MH+MV\\ HP-HM+VP-VM& HP-HM-VP+VM& PP-PM-MP+MM\\ HR-LL+VR-RL& HR-VR+VL-HL& PR-MR+ML-PL\end{array}\begin{array}{c}RH+RV-LH-LV\\ RH-RV-LH+LV\\ RP-RM-LP+LM\\ RR-RL-LR+LL\end{array}\right]. \end{gathered}$$

For convenience, each reflectance image is indexed by two capital letters: the first letter represents the incidence polarization state and the second letter represents the detection polarization state. For example, HV represents a reflectance image of horizontally polarized ($H$) incidence and vertically polarized ($V$) detection.

2.5.

Results from the Experiment and SCBM Simulations

Using the new program based on SCBM, we simulated the experimental results of the sample. All the parameters used in the simulation matched the experiments. As shown in Fig. 4, the simulated result using SCBM agrees well with the experimental result of the strained microsphere-silk polyacrylamide. Although Fig. 4(a) and 4(b) are still different quantitatively, the characteristic features in the two sets of Mueller matrix elements are very similar. The detailed analysis about the effect of birefringence on the backscattering Mueller matrix is given in the following section.

Fig. 4

The backscattering Mueller matrices of the strained microsphere-silk polyacrylamide with silk fibers aligned along the $y$ axis and the direction of strain (the extraordinary axis of birefringence) along the 45-deg angle on the $x\text{-}y$ plane (a) and the SCBM simulation (b). The size of each image is $1\times 1\mathrm{cm}$. Note that the $m11$ images used the same color maps but with the range from 0 to 1.

3.1.

Comparison with the Simulations of SCSM

The characteristic features of the Mueller matrix of sphere scattering media have been reported previously by other groups.^17,18 For the isotropic scattering media, the Mueller matrix has diagonal symmetry properties. The $m14$ and $m41$ elements are always blanks. The patterns of $m11$ and $m44$ are both circular. The $m12$ element is quatrefoil, and its intensity distributions along the $x$ and $y$ axes are nearly the same. The $m22$ element has a symmetric cross-like pattern with almost identical intensity distribution along the $x$ and $y$ axes. There are several pairs of elements of similar intensity patterns but different rotations. For example, $m13$ can be obtained by rotating $m12$ by 45 deg, and the same relation applies to the $m31/m21$, $m22/m33$, and $m24/m34$ pairs.

In our previous paper,¹¹ the Mueller matrix patterns of the sphere–cylinder scattering medium were described in detail. According to the Mie scattering theory for infinitely long cylinders,^13,14 polarization-maintaining photons tend to be scattered in the direction perpendicular to the cylindrical scatterers (i.e., the $x$ axis in Fig. 4). Therefore the intensity distributions of all the Mueller matrix elements around the $x$ axis are higher than those around the $y$ axis. For example, $m11$ has the typical rhombus profile with an elongation along the $x$ axis. $m12$ and $m21$ have similar shapes of quatrefoils, but the intensity along the $x$ axis is higher than that along the $y$ axis. The total intensity distributions of the $m13$ and $m31$ elements are weaker than those of the $m12$ and $m21$ elements. The $m22$ element still shows a cross-like pattern, but the intensity along the $x$ axis is higher than that along the $y$ axis. The total intensity distribution and size of the $m33$ element are weaker than those of the $m22$ element.

Compared to the sphere–cylinder scattering medium, the characteristic features of the Mueller matrix elements for the sphere–cylinder birefringence medium are very different. Monte Carlo–simulated backscattering Mueller matrix patterns of SCSM and SCBM are shown in Fig. 5. The Mueller matrix elements have been normalized to the $m11$ element to compensate for the radial decay of intensity. The extraordinary axis of the birefringent medium is along the $y$ axis in the laboratory frame. For easier comparison, parameters of the scatterers are the same as those in Ref. 11: The diameters of spheres and cylinders are 0.2 and 1.5 μm, and their refractive indices are 1.59 and 1.56, respectively. The refractive index of the surrounding medium is 1.33. The cylinders are aligned in the $x\text{-}y$ plane, and their orientations fluctuate around the $y$ axis following a Gaussian distribution with a standard deviation of 10 deg half width. The scattering coefficients of spheres and cylinders are 10 and $65{\mathrm{cm}}^{-1}$, respectively. In SCBM, since the birefringence value of tissues ranges from $1\times {10}^{-4}$ to $1\times {10}^{-2}$ (Ref. 19), we set $\mathrm{\Delta }n$ to 0, $5\times {10}^{-4}$, $1\times {10}^{-3}$, and $1\times {10}^{-2}$.

Fig. 5

Monte Carlo simulations of backscattering Mueller matrices. Sphere–cylinder scattering medium and a sphere–cylinder birefringence medium (a) with $\mathrm{\Delta }n=5\times {10}^{-4}$ (b); $\mathrm{\Delta }n=1\times {10}^{-3}$ (c); and $\mathrm{\Delta }n=1\times {10}^{-2}$ (d). The size of each image is $1\times 1\mathrm{cm}$. Note that the $m11$ images used the same color maps but with the range from 0 to 1.

Figure 5 shows the influence of linear birefringence on the spatial intensity distributions of the backscattering Mueller matrix elements. The detailed features are summarized as follows: (1) The $m11$ has the same typical rhombus profile as that in SCSM. Linear birefringence does not affect the $m11$ element, which relates only to unpolarized light intensity. (2) The $m12$ and $m21$ elements have the same shapes of quatrefoils, and the $m22$ also has the same cross-like pattern, as those in the SCSM simulations. For the present configuration, where both the extraordinary axis and the cylindrical scatterers are along the vertical ($y$) direction, $m12$, $m21$, and $m22$ are hardly affected by linear birefringence. This is because the Mueller matrix patterns represent mainly contributions from the polarization-maintaining photons. These elements relate only to the horizontal (0 deg) and vertical (90 deg) polarization components, which are either parallel or perpendicular to the extraordinary axis and subject to minimal effects due to the birefringent medium. (3) The intensity distributions of the $m13$, $m23$, $m33$, $m31$, $m32$, $m33$, and $m44$ elements are sensitive to birefringence $\mathrm{\Delta }n$, since circular or 45-deg and 135-deg linear polarization states are affected by retardation due to birefringence. As $\mathrm{\Delta }n$ increases, the $m13$, $m31$, $m33$, and $m44$ patterns stay in the same shape but shrink in size. The shapes of the $m23$ and $m32$ patterns change to quatrefoils stretched along the $x$ axis. (4) The $m14(m41)$), $m24(m42)$, and $m34(m43)$ patterns have similar shapes to those of $m13$, $m23$, and $m33$, respectively, but the contrasts are different. The transformation between circular and linear polarization states in the turbid medium is enhanced by linear birefringence.⁹

3.2.

Effects of the Value and Orientation of Birefringence on the Mueller Matrix Patterns

In Fig. 6, we compare the different responses of $m22$ and $m33$ to changes in the birefringence value $\mathrm{\Delta }n$ while keeping both the extraordinary axis and the cylindrical scatterers along the $y$ axis in the laboratory frame. As the value of birefringence $\mathrm{\Delta }n$ increases, the total intensity of $m22$ remains almost the same, but the $m33$ intensity decreases rapidly and approaches zero when $\mathrm{\Delta }n$ exceeds $1\times {10}^{-3}$. One may use $m33$ to detect the presence of birefringence.

Fig. 6

The total intensity of $m22$ and $m33$ with different birefringence value $\mathrm{\Delta }n$. The extraordinary axis of the birefringence and the cylindrical direction are both along the $y$ axis in the laboratory coordinate. The data are normalized by the maximum of $m22$.

As $\mathrm{\Delta }n$ increases, the intensity of the $m13$, $m23$, $m31$, $m32$, and $m44$ Mueller matrix elements quickly decrease to zero, and the terms related to 0-deg or 90-deg linear polarization ($m12$, $m21$) remain almost the same. These phenomena are consistent with the simulation results of SBM (see Fig. 9 in Ref. 9), in which the degree of polarization (DOP) of diffusely reflected light does not change significantly with $\mathrm{\Delta }n$ when the incident light is in horizontal (0 deg) or vertical (90 deg) polarization, but does decrease fast and approaches zero with increasing $\mathrm{\Delta }n$ if the incident light is 45-deg linear polarization or circle polarization. The decrease of DOP means more recorded photons are depolarized and the Mueller matrix loses their features. The intensity of the $m14$, $m24$, $m34$, $m41$, $m42$, and $m43$ elements will first increase then quickly decrease to zero as $\mathrm{\Delta }n$ increases from zero. This is also consistent with the simulation results of SBM (see Fig. 4 in Ref. 16).

In the above simulations, the extraordinary axis of birefringence has been set parallel to the orientation of the cylindrical scatterers, and to the $y$ axis of the laboratory frame. We may also assume that the two directions are different and simulate the Mueller matrix patterns for different orientation of the extraordinary axis. As shown in a previous work,¹¹ the long axis of the $m11$ pattern is always perpendicular to the direction of the fibers. For simplification, we set the direction of the cylindrical scatterers along the $y$ axis in the laboratory frame and the extraordinary axis of birefringence in the $x\text{-}y$ plane. All the patterns except for the $m11$ display asymmetric behavior, as shown in Fig. 4(b). A further study indicates that $m22$ is the most sensitive to the angle mismatch. As shown in Fig. 7, we changed the extraordinary axis from the $x$ axis to $y$ axis at intervals of 15 deg in the simulations while keeping the other parameters the same as in Fig. 5(c). The $m22$ pattern has a cross-like shape when the extraordinary axis is parallel or perpendicular to the cylinders, but is distorted when neither the fast axis nor the slow axis of birefringence matches the fiber orientation. Therefore, the $m22$ pattern may serve as an indicator for the mismatch between the directions of the birefringence axis and the fiber orientation. In real biological tissues, the birefringence axis may vary in the 3-D space. The corresponding features can also be simulated using the current SCBM program.

Fig. 7

$m22$ patterns from sphere–cylinder birefringence media with the cylinders aligned along the $y$ axis and the angle between the extraordinary axis and $x$ axis of 0, 15, 30, 45, 60, 75, and 90 deg.

3.3.

Comparison with Experimental Results on Skeletal Muscles

In a previous paper,¹⁰ we demonstrated that SCSM can be used to characterize the structural and optical properties of skeletal muscle. The Mueller matrix elements of both skeletal muscles and a model sample containing polystyrene microspheres and silk fibers show very similar features. However, more careful examinations reveal inconsistency between the experimentally obtained and the simulated results. For example, the third row and the third column of the simulated Mueller matrix using SCSM display clear patterns that are not apparent in the experimental results of skeletal muscle. Simulations using SCBM show that this discrepancy is due to the birefringence effect.

In the present work, the skeletal muscle samples are fresh bovine sternomandibularis muscle excised from the animals 2 h after slaughtering. The Mueller matrix elements of the skeletal muscle were compared with the simulations using SCBM, SCSM, and SBM as shown in Fig. 8.

Fig. 8

The backscattering Mueller matrices of fresh skeletal muscle (a), the SCBM simulation (b), the SCSM simulation (c), and the SBM simulation (d). The muscle fiber orientation, the extraordinary axis of the birefringent medium, and the cylinder orientation are all along the $y$ axis. Note that the $m11$ images used the same color maps but with the range from 0 to 1.

The following parameters are used in the SCBM. The extraordinary axis of the birefringent surrounding medium and the direction of the muscle fibers are both along the $y$ axis. The birefringence value $\mathrm{\Delta }n$ is $2\times {10}^{-3}$ (Ref. 20). The diameters of the spheres and infinitely long cylinders are 0.2 and 1.5 μm, respectively, and the refractive indices of scatterers are both 1.4. The average refractive index of the surrounding medium is taken as 1.33. The standard deviation of angular distribution of the cylinders is 10 deg. Using the previous published results¹⁰ as a reference, the ratio of the scattering coefficients of spheres and cylinders is set at $1/4$ to give the best match between the shapes of the experimental and simulated images. The total scattering coefficients along the $y$ axis is $50\pm 10{\mathrm{cm}}^{-1}$, which is close to the scattering coefficient of the skeletal muscle.¹⁰ In the SCSM, the birefringence value is set to zero, and in the SBM, the scattering coefficient of cylinders is set to zero. The other parameters are the same as those in the SCBM.

Compared with the experiments, simulations of the three anisotropic models (SCBM, SCSM, and SBM) had respective similarities and differences. The characteristic features in the Mueller matrix patterns of anisotropic turbid media showed contributions by both optical anisotropy due to birefringence and scattering anisotropy due to cylindrical scatterers.

As shown in Fig. 8, for the Mueller matrix elements of a muscle sample, $m11$ has the typical rhombus shape, and $m22$ has a cross-like pattern with the dominant distribution along the $x$ axis. The simulations using SCBM and SCSM regenerated similar features in these two elements. However, simulations using SBM resulted in totally different shapes. The $m11$ pattern is centrosymmetric and $m22$ has a symmetric cross-like pattern with identical intensity distributions along the $x$ and $y$ axes. Thus one has to include cylindrical scatterers in the model to regenerate these characteristic features in $m11$ and $m22$ of anisotropic tissues. These features of $m11$, $m12$, $m21$, and $m22$ are due to the anisotropy of the cylinders, since the polarization-maintaining photons tend to be scattered to the direction perpendicular to the cylindrical scatterers (the $x$ axis).¹¹

On the other hand, the third and fourth rows and columns of the experimentally obtained Mueller matrix for the muscle samples show only weak patterns. Simulations using both SCBM and SBM regenerated such features, but simulations using SCSM resulted in much clearer patterns, for example clear quatrefoils for the $m13$ and $m31$ elements. More detailed simulations using SCBM showed that the birefringence effect is responsible for smoothing out the distinctive features in the third and fourth rows and columns of the Mueller matrix elements.

The above analysis proves that simulations using SCBM result in better agreements with the experiments than SCSM and SBM. It should be pointed out that the Mueller matrix elements related to the third and fourth rows and columns do not show good agreements between the skeletal muscle experiments [Fig. 8(a) and also Ref. 7] and the SCBM simulations [Fig. 8(b)]. Simulation shows that these elements are all weak and noisy when birefringence exists, as shown in Fig. 8(b) and 8(d). The corresponding experimental patterns also display noisy and irregular patterns,⁷ which may also contain measurement errors and inhomogeneity effects of the tissue sample.

One has to consider contributions from both the cylindrical scatterers and the birefringence effect to explain the characteristic features of polarized photons scattering in complicated anisotropic turbid media such as skeletal muscles. These two terms on tissue anisotropy are perhaps corresponding to the anisotropy of form defined by the tissue microstructure and the anisotropy of material controlled by intrinsic anisotropic character of metabolic molecules.²⁰

In addition, all the experimentally obtained $m22$ elements showed cross-like patterns without distortion. From the above analysis in Sec. 3.2, we conclude that the extraordinary axis of the birefringent media is parallel to the direction of the fibrous structures in all the samples we have measured, including the skeletal muscles. This may be a hint that the two come from the same origin.