Polarized Adding Method of Discrete Ordinate Approximation for Ultraviolet-Visible and Near-Infrared Radiative Transfer

Kun Wu Feng Zhang fengzhang@fudan.edu.cn Wenwen Li Fengzi Bao Yi-ning Shi Key Laboratory of Meteorological Disaster, Ministry of Education (KLME)/ Joint International Research Laboratory of Climate and Environment Change (ILCEC)/ Collaborative Innovation Center on Forecast and Evaluation of Meteorological Disasters (CIC-FEMD) / Institute for Climate and Application Research (ICAR), Nanjing University of Information Science and Technology, Nanjing 210044, China Department of Atmospheric and Oceanic Sciences and Institutes of Atmospheric Sciences, Fudan University, Shanghai 200433, China CMA Earth System Modeling and Prediction Centre, China Meteorological Administration, Beijing 100081, China

Abstract

The polarization characteristics of atmospheric scattering are important and should not be ignored in radiative transfer simulations. In this study, a new vector radiative transfer model called the polarized adding method of discrete ordinate approximation (POLDDA) is proposed for use in remote sensing applications for ultraviolet-visible and near-infrared spectra. The single-layer radiative transfer process and inhomogeneous multi-layer connection are solved using the discrete ordinate method (DOM) and adding methods, respectively. By combining the advantages of DOM and the adding method, the Stokes vector (including the $I$ -, $Q$ -, $U$ -, and $V$ -components) calculated using the new method conforms to the results of PolRadtran/RT3, whether in a Rayleigh scattering atmosphere or the water cloud case. Moreover, the relative root-mean-square error (RMSE) values of the Stokes vector for the test cases between MYSTIC and the new method or RT3 prove the accuracy of the proposed method. Meanwhile, the new method has a higher computational efficiency than RT3, particularly for an atmosphere with a large scattering optical depth. Unlike RT3, the computation time of the proposed method does not increase with the optical depth of each layer.

keywords:

vector radiative transfer , ultraviolet-visible and near-infrared , adding method , discrete ordinate approximation

1 Introduction

In recent years, interest in the development of vector radiative transfer methods and instruments with polarization capabilities for aircraft, balloons, satellites, and ground-based platforms has increased rapidly. Radiation scattered by air molecules (Rayleigh scattering), aerosols, and cloud particles is polarized and shows different polarization characteristics based on the scattering event [1, 2]. The degree of linear and circular polarization produced by clouds and aerosol particles is more sensitive to the shape, size, and refractive index of the poly-dispersed small scattering particles. Polarization can provide additional information on atmospheric optical phenomena and their constituents [3], including the properties of atmospheric components and vertical distributions. Numerous studies [4, 5, 6, 7, 8, 9] have shown that approximately 10% of the errors introduced in atmospheric radiative transfer simulations and remote sensing are caused by ignoring the atmospheric polarization, particularly at short wavelengths [10, 5, 8, 11] owing to more scattering compared to long wavelengths [12]. Therefore, a scalar radiative transfer model cannot sufficiently describe the nature of radiation processes [13].

Various vector radiative transfer schemes have been developed, which commonly use the adding-doubling method [14, 15, 16, 17, 18, 19], discrete ordinates method [20, 21, 22, 23, 24], Monte-Carlo solutions method [25, 26, 27, 28], successive orders of scattering method [29, 30, 1], Invariant Imbedding Method [31], and spherical harmonic method [32]. Each method has its own advantages and disadvantages. For instance, the adding-doubling method deals with multiple scattering using the doubling method [17]. Though the adding method can simplify the layer-to-layer connection, it is time-consuming to compute an atmospheric layer with a thick optical depth, such as a sky with thick stratus. This is caused by the increase in the number of thin sub-layers (a number 2^N of identical thin layers), each of which is characterized by single scattering. The problem can be overcome using the discrete ordinate method (DOM). DOM solves the vector transfer equation using eigenvectors and eigenvalues [33, 22]. It can calculate the interior variation of the reflection and transmission function within a layer; therefore, it is accurate and efficient for thick scattering media [21], such as aerosols and clouds.

In modeling scalar solar radiation, there is the method [34] combines the advantages of the DOM and adding method, which via using DOM for a single-layer radiative transfer solution and the adding method employed for the inhomogeneous multi-layer connection in a plane-parallel scattering atmosphere. In this study, it will be expanded to vector radiative transfer to optimize the accuracy and computational efficiency. The remainder of this study is organized as follows. Section 2 describes the formulation and algorithm of the proposed method. In Section 3, the calculation accuracy and efficiency of the proposed method are evaluated for different scattering atmospheres. Finally, the conclusions and discussion are presented in Section 4.

2 Solution by discrete ordinates method for a single-layer

2.1 Solar radiation

The polarized radiative transfer equation with only a solar source in a plane-parallel medium is given by [10, 35, 36]

\begin{split}\mu\frac{d\mathbf{L}(\tau,\mu,\varphi)}{d\tau}=&\mathbf{L}(\tau,\mu,\varphi)-\frac{\omega}{4\pi}\int_{0}^{2\pi}\int_{-1}^{1}\mathbf{Z}(\mu,\varphi,\mu^{\prime},\varphi^{\prime})\mathbf{L}(\tau,\mu^{\prime},\varphi^{\prime})d\mu^{\prime}d\varphi^{\prime}\\ &-\frac{\omega}{4\pi}\mathbf{Z}(\mu,\varphi,\mu_{0},\varphi_{0})\mathbf{F_{0}}(\mu_{0},\varphi_{0})e^{-\tau/\mu_{0}},\end{split}

(1)

where $\mathbf{L}=[I,Q,U,V]^{\mathbb{T}}$ represents the Stokes vector, and the superscript $\mathbb{T}$ indicates a matrix transpose. $\tau$ indicates the optical depth, $\omega$ signifies the single scattering albedo, $\varphi$ indicates the azimuthal angle, and $\mu$ denotes the zenith angle cosine. $\mathbf{Z}$ is a 4 $\times$ 4 scattering matrix, defined as $\mathbf{Z}=\mathbf{C}(\pi-i_{2})\mathbf{P}(\mu,\varphi,\mu^{\prime},\varphi^{\prime})\mathbf{C}(-i_{1})$ , where $\mathbf{C}$ represents the rotational matrices with $i_{1}$ and $i_{2}$ defining the rotation angles, and $\mathbf{P}$ represents the scattering phase matrix, and each element of $\mathbf{P}$ can be expanded using the Legendre series. $\mathbf{F_{0}}$ represents the incident solar Stokes vector at the top of the atmosphere (TOA).

The radiation transfer equation can be split into 2M equations by expanding the phase matrix and $\mathbf{L}$ into the Fourier cosine and sine series. The discrete ordinate method is used to solve the equation and neglecting the superscript $m$ , Eq.(1) can be written as in matrix form as

\frac{d}{d\tau}{\left[{\begin{array}[]{*{20}{c}}{{{\bf{L}}^{+}}(\tau)}\\ {{{\bf{L}}^{-}}(\tau)}\end{array}}\right]}={\left[{\begin{array}[]{cc}\bf{X}_{11}&\bf{X}_{12}\\ \bf{X}_{21}&\bf{X}_{22}\end{array}}\right]}{\left[{\begin{array}[]{c}{{{\bf{L}}^{+}}}\\ {{{\bf{L}}^{-}}}\end{array}}\right]}+{\left[{\begin{array}[]{c}\bf{E}^{+}\\ \bf{E}^{-}\end{array}}\right]}e^{-\tau/\mu_{0}},

(2)

where $\mathbf{L}^{\pm}(\tau)=[I(\tau,\pm\mu_{1}),\cdot\cdot\cdot,I(\tau,\pm\mu_{N}),Q(\tau,\pm\mu_{1}),\cdot\cdot\cdot,Q(\tau,\pm\mu_{N}),\\ U(\tau,\pm\mu_{1}),\cdot\cdot\cdot,U(\tau,\pm\mu_{N}),V(\tau,\pm\mu_{1}),\cdot\cdot\cdot,V(\tau,\pm\mu_{N})]^{\mathbb{T}}$ . The details are given in Appendix A.

This can be solved using the eigenvalue method [23] as follows:

{\left[{\begin{array}[]{*{20}{c}}{{{\bf{L}}^{+}}(\tau)}\\ {{{\bf{L}}^{-}}(\tau)}\end{array}}\right]}=\bf{G}\circ\bf{K}{\left[{\begin{array}[]{c}\bf{C}_{1}\\ \bf{C}_{2}\end{array}}\right]}-\mu_{0}{\left[{\begin{array}[]{c}{{{\bf{E}}^{+}}}\\ {{{\bf{E}}^{-}}}\end{array}}\right]}e^{-\tau/\mu_{0}}

(3)

where $\mathbf{G}$ is composed of the eigenvectors of $\mathbf{X}$ , and $\mathbf{K}$ is composed of the eigenvalues of $\mathbf{X}$ . $[\bf{C}_{1},\bf{C}_{2}]^{\mathbb{T}}$ can be obtained from the boundary condition, that is, $\bf{L}^{-}(0)=0$ and $\bf{L}^{+}(\tau_{1})=0$ . By solving Eq.(6), $\bf{L}^{+}(0)$ and $\bf{L}^{-}(\tau_{1})$ can be expressed as follows: $\bf{L}^{+}(0)=\left[{\begin{array}[]{c}\mathbf{R}^{I}\\ \mathbf{R}^{Q}\\ \mathbf{R}^{U}\\ \mathbf{R}^{V}\\ \end{array}}\right]$ , $\bf{L}^{-}(\tau_{1})=\left[{\begin{array}[]{*{20}{c}}\mathbf{T}^{I}\\ \mathbf{T}^{Q}\\ \mathbf{T}^{U}\\ \mathbf{T}^{V}\\ \end{array}}\right]$ . $\mathbf{R}$ and $\mathbf{X}$ represent the direct reflection and transmission matrices, respectively, of a single layer with optical depth $\tau_{1}$ , and can be expressed as

	$\displaystyle\mathbf{R}^{I}=\left[R^{I\leftarrow I}(0,\mu_{1},\mu_{0}),\cdots,R^{I\leftarrow I}(0,\mu_{N},\mu_{0})\right]^{\mathbb{T}},$		(4)
	$\displaystyle\mathbf{R}^{Q}=\left[R^{Q\leftarrow I}(0,\mu_{1},\mu_{0}),\cdots,R^{Q\leftarrow I}(0,\mu_{N},\mu_{0})\right]^{\mathbb{T}},$
	$\displaystyle\mathbf{R}^{U}=\left[R^{U\leftarrow I}(0,\mu_{1},\mu_{0}),\cdots,R^{U\leftarrow I}(0,\mu_{N},\mu_{0})\right]^{\mathbb{T}},$
	$\displaystyle\mathbf{R}^{V}=\left[R^{V\leftarrow I}(0,\mu_{1},\mu_{0}),\cdots,R^{V\leftarrow I}(0,\mu_{N},\mu_{0})\right]^{\mathbb{T}},$
	$\displaystyle\mathbf{T}^{I}=\left[T^{I\leftarrow I}(\tau_{1},\mu_{1},\mu_{0}),\cdots,T^{I\leftarrow I}(\tau_{1},\mu_{N},\mu_{0})\right]^{\mathbb{T}},$
	$\displaystyle\mathbf{T}^{Q}=\left[T^{Q\leftarrow I}(\tau_{1},\mu_{1},\mu_{0}),\cdots,T^{Q\leftarrow I}(\tau_{1},\mu_{N},\mu_{0})\right]^{\mathbb{T}},$
	$\displaystyle\mathbf{T}^{U}=\left[T^{U\leftarrow I}(\tau_{1},\mu_{1},\mu_{0}),\cdots,T^{U\leftarrow I}(\tau_{1},\mu_{N},\mu_{0})\right]^{\mathbb{T}},$
	$\displaystyle\mathbf{T}^{V}=\left[T^{V\leftarrow I}(\tau_{1},\mu_{1},\mu_{0}),\cdots,T^{V\leftarrow I}(\tau_{1},\mu_{N},\mu_{0})\right]^{\mathbb{T}}.$

A connection exists between $I$ -, $Q$ -, $U$ - and $V$ -component transfer in a medium considering scattering. For the I-component, its change is not only from the I-component self, but also from the other three components. The superscript ’ $Q\leftarrow I$ ’ denotes the reflection or transmission of $I$ -component to $Q$ -component dimension. The other parameters are defined in a similar manner.

2.2 Diffuse radiation

The diffuse equation of solar-polarized radiative transfer can be given as

\mu\frac{d\mathbf{L}(\tau,\mu,\varphi)}{d\tau}=\mathbf{L}(\tau,\mu,\varphi)-\frac{\omega}{4\pi}\int_{0}^{2\pi}\int_{-1}^{1}\mathbf{Z}(\mu,\varphi,\mu^{\prime},\varphi^{\prime})\mathbf{L}(\tau,\mu^{\prime},\varphi^{\prime})d\mu^{\prime}d\varphi^{\prime}.

(5)

This solution is also required for the addition method. The details of the solutions of the diffuse radiation for the $I$ and $Q$ components are almost the same as those for the polarized thermal infrared radiative transfer described in [37]. The $U$ and $V$ components are obtained in a similar manner. Thus, the reflection and transmission matrices of diffuse radiation are given as


$\displaystyle\overline{\mathbf{R}}$	$\displaystyle=\left[{\begin{array}[]{cccc}\overline{\cal R}^{I\leftarrow I}&{\overline{\cal R}}^{I\leftarrow Q}&{\overline{\cal R}}^{I\leftarrow U}&{\overline{\cal R}}^{I\leftarrow V}\\ {\overline{\cal R}}^{Q\leftarrow I}&{\overline{\cal R}}^{Q\leftarrow Q}&{\overline{\cal R}}^{Q\leftarrow U}&{\overline{\cal R}}^{Q\leftarrow V}\\ {\overline{\cal R}}^{U\leftarrow I}&{\overline{\cal R}}^{U\leftarrow Q}&{\overline{\cal R}}^{U\leftarrow U}&{\overline{\cal R}}^{U\leftarrow V}\\ {\overline{\cal R}}^{V\leftarrow I}&{\overline{\cal R}}^{V\leftarrow Q}&{\overline{\cal R}}^{V\leftarrow U}&{\overline{\cal R}}^{V\leftarrow V}\\ \end{array}}\right]_{4N\times 4N}$	(6e)
$\displaystyle\overline{\mathbf{T}}$	$\displaystyle=\left[{\begin{array}[]{cccc}{\overline{\cal T}}^{I\leftarrow I}\qquad{\overline{\cal T}}^{I\leftarrow Q}\qquad{\overline{\cal T}}^{I\leftarrow U}\qquad{\overline{\cal T}}^{I\leftarrow V}\\ {\overline{\cal T}}^{Q\leftarrow I}\qquad{\overline{\cal T}}^{Q\leftarrow Q}\qquad{\overline{\cal T}}^{Q\leftarrow U}\qquad{\overline{\cal T}}^{Q\leftarrow V}\\ {\overline{\cal T}}^{U\leftarrow I}\qquad{\overline{\cal T}}^{U\leftarrow Q}\qquad{\overline{\cal T}}^{U\leftarrow U}\qquad{\overline{\cal T}}^{U\leftarrow V}\\ {\overline{\cal T}}^{V\leftarrow I}\qquad{\overline{\cal T}}^{V\leftarrow Q}\qquad{\overline{\cal T}}^{V\leftarrow U}\qquad{\overline{\cal T}}^{V\leftarrow V}\\ \end{array}}\right]_{4N\times 4N}.$	(6j)

where ${{\overline{\cal R}}^{I\leftarrow I}}={\left[{\begin{array}[]{*{20}{c}}{(1+\delta_{0,m}){a_{1}}{R^{I\leftarrow I}}({\mu_{1}},{\mu_{1}}){\mu_{1}}}&\cdots&{(1+\delta_{0,m}){a_{N}}{R^{I\leftarrow I}}({\mu_{1}},{\mu_{N}}){\mu_{N}}}\\ \cdots&\cdots&\cdots\\ {(1+\delta_{0,m}){a_{1}}{R^{I\leftarrow I}}({\mu_{N}},{\mu_{1}}){\mu_{1}}}&\cdots&{(1+\delta_{0,m}){a_{N}}{R^{I\leftarrow I}}({\mu_{N}},{\mu_{N}}){\mu_{N}}}\end{array}}\right]_{N\times N}}$ ,
${{\overline{\cal T}}^{I\leftarrow I}}={\left[{\begin{array}[]{*{20}{c}}{(1+\delta_{0,m}){a_{1}}{{\tilde{T}}^{I\leftarrow I}}({\mu_{1}},{\mu_{1}}){\mu_{1}}}&\cdots&{(1+\delta_{0,m}){a_{N}}{{\tilde{T}}^{I\leftarrow I}}({\mu_{1}},{\mu_{N}}){\mu_{N}}}\\ \cdots&\cdots&\cdots\\ {(1+\delta_{0,m}){a_{1}}{{\tilde{T}}^{I\leftarrow I}}({\mu_{N}},{\mu_{1}}){\mu_{1}}}&\cdots&{(1+\delta_{0,m}){a_{N}}{{\tilde{T}}^{I\leftarrow I}}({\mu_{N}},{\mu_{N}}){\mu_{N}}}\end{array}}\right]_{N\times N}}$ , and other matrices on the right hand side of Eq.(6) are defined similarly.

3 Four invariance principles and adding method for inhomogeneous multi-layer connection

In this section, the four invariance principles [10] are extended to the polarized solar radiative transfer process and applied to the new model. A schematic diagram is shown in Figure 1 (a) $\sim$ (d): More specifically, we consider a combination of two inhomogeneous layers with $\tau_{1}$ (first layer) and $\tau_{2}$ (second layer) as the optical depths. $R_{i}^{a\leftarrow b}(\mu,\mu^{\prime})$ and $T_{i}^{a\leftarrow b}(\mu,\mu^{\prime})$ ( $a=I,Q,U,V$ and $b=I,Q,U,V$ ) indicate the reflection and transmission functions of the $i$ -th layer. The superscript ^∗ denotes that the radiation originates from below. From Figure 1, the four invariance principles can be expressed as follows:

(a) The reflected matrices $\mathbf{A}$ ( $A^{I}(\mu)$ / $A^{Q}(\mu)$ / $A^{U}(\mu)$ / $A^{V}(\mu)$ ) of the Stokes vector at level 2 originate from two parts: the reflection of the direct solar beam by the second layer $\mathbf{R}_{2}(0,\mu)$ ( $R^{I}_{2}(0,\mu)$ / $R^{Q}_{2}(0,\mu)$ / $R^{U}_{2}(0,\mu)$ / $R^{V}_{2}(0,\mu)$ ) and the coupled reflection of the downward components $\mathbf{D}$ ( $D^{I}(\mu)$ / $D^{Q}(\mu)$ / $D^{U}(\mu)$ / $D^{V}(\mu)$ ) through the second layer;

(b) The transmitted matrices $\mathbf{D}$ of the Stokes vector at level 2 originate from two parts: direct transmission by the first layer $\mathbf{T}_{1}(\tau_{1},-\mu)$ ( $T^{I}_{1}(\tau_{1},-\mu)$ / $T^{Q}_{1}(\tau_{1},-\mu)$ / $T^{U}_{1}(\tau_{1},-\mu)$ / $T^{V}_{1}(\tau_{1},-\mu)$ ) and the coupled reflection of the upward components $\mathbf{A}$ through the first layer.

(c) The reflected matrices $\mathbf{R}_{1,2}(0,\mu)$ ( $R^{I}_{1,2}(0,\mu)$ / $R^{Q}_{1,2}(0,\mu)$ / $R^{U}_{1,2}(0,\mu)$ / $R^{V}_{1,2}(0,\mu)$ ) of the Stokes vector at the top of the two-layer are composed of two parts: the reflection by the first layer $\mathbf{R}_{1}(0,\mu)$ ( $R^{I}_{1}(0,\mu)$ / $R^{Q}_{1}(0,\mu)$ / $R^{U}_{1}(0,\mu)$ / $R^{V}_{1}(0,\mu)$ ), and the total transmission including the direct beam transmissions of $\exp(-\frac{\tau_{1}}{\mu})$ and $\mathbf{T}_{1}^{*}(\mu,\mu^{\prime})$ of the upward components $\mathbf{A}$ through the first layer.

(d) The transmitted matrices $\mathbf{T}_{1,2}(\tau_{1}+\tau_{2},-\mu)$ ( $T^{I}_{1,2}(\tau_{1}+\tau_{2},-\mu)$ / $T^{Q}_{1,2}(\tau_{1}+\tau_{2},-\mu)$ / $T^{U}_{1,2}(\tau_{1}+\tau_{2},-\mu)$ / $T^{V}_{1,2}(\tau_{1}+\tau_{2},-\mu)$ ) of the Stokes vector at the bottom of the two-layer are composed of two parts: transmission of the direct solar beam $\exp(-\frac{\tau_{1}}{\mu_{0}})$ through the second layer and the total coupled transmission, including $\exp(-\frac{\tau_{2}}{\mu})$ and $\mathbf{T}_{2}(\mu,\mu^{\prime})$ of the downward components $\mathbf{D}$ through the second layer.

Based on the above (a $\sim$ d) statements, we can write the four invariance principles in vector solar radiative transfer as


	$\displaystyle\mathbf{A}=\mathbf{R}_{2}e^{-\frac{\tau_{1}}{\mu_{0}}}+2\int_{0}^{1}\overline{\mathbf{R}}_{2}\mathbf{D}\mu^{\prime}d\mu^{\prime},$		(7a)
	$\displaystyle\mathbf{D}=\mathbf{T}_{1}+2\int_{0}^{1}\overline{\mathbf{R}}^{*}_{1}\mathbf{A}\mu^{\prime}d\mu^{\prime},$		(7b)
	$\displaystyle\mathbf{R}_{1,2}=\mathbf{R}_{1}+\mathbf{A}e^{-\frac{\tau_{1}}{\mu_{0}}}+2\int_{0}^{1}\overline{\mathbf{T}}^{*}_{1}\mathbf{A}\mu^{\prime}d\mu^{\prime},$		(7c)
	$\displaystyle\mathbf{T}_{1,2}=\mathbf{T}_{2}e^{-\frac{\tau_{1}}{\mu_{0}}}+\mathbf{D}e^{-\frac{\tau_{2}}{\mu_{0}}}+2\int_{0}^{1}\overline{\mathbf{T}}_{2}\mathbf{D}\mu^{\prime}d\mu^{\prime}.$		(7d)

In Eq.(7), $\mathbf{A}=[\mathbf{A}^{I},\mathbf{A}^{Q},\mathbf{A}^{U},\mathbf{A}^{V}]^{\mathbb{T}}$ , $\mathbf{D}=[\mathbf{D}^{I},\mathbf{D}^{Q},\mathbf{D}^{U},\mathbf{D}^{V}]^{\mathbb{T}}$ ,
$\mathbf{R}_{1,2}=[\mathbf{R}_{1,2}^{I},\mathbf{R}_{1,2}^{Q},\mathbf{R}_{1,2}^{U},\mathbf{R}_{1,2}^{V}]^{\mathbb{T}}$ , $\mathbf{T}_{1,2}=[\mathbf{T}_{1,2}^{I},\mathbf{T}_{1,2}^{Q},\mathbf{T}_{1,2}^{U},\mathbf{T}_{1,2}^{V}]^{\mathbb{T}}$ .
Others matrix in Eq.(7) are detailed in the Appendix B.

$N$ -node Gaussian integration was used to handle the integration; thus, Eq.(7) can be written as follows:


	$\displaystyle\left[\begin{array}[]{cccc}\mathbf{A}^{I}\\ \mathbf{A}^{Q}\\ \mathbf{A}^{U}\\ \mathbf{A}^{V}\\ \end{array}\right]=\left[\begin{array}[]{cccc}\mathbf{R}_{2}^{I}\\ \mathbf{R}_{2}^{Q}\\ \mathbf{R}_{2}^{U}\\ \mathbf{R}_{2}^{V}\\ \end{array}\right]e^{-\frac{\tau_{1}}{\mu_{0}}}+\left[\begin{array}[]{cccc}\overline{\cal R}_{2}^{I\leftarrow I}&\overline{\cal R}_{2}^{I\leftarrow Q}&\overline{\cal R}_{2}^{I\leftarrow U}&\overline{\cal R}_{2}^{I\leftarrow V}\\ \overline{\cal R}_{2}^{Q\leftarrow I}&\overline{\cal R}_{2}^{Q\leftarrow Q}&\overline{\cal R}_{2}^{Q\leftarrow U}&\overline{\cal R}_{2}^{Q\leftarrow V}\\ \overline{\cal R}_{2}^{U\leftarrow I}&\overline{\cal R}_{2}^{U\leftarrow Q}&\overline{\cal R}_{2}^{U\leftarrow U}&\overline{\cal R}_{2}^{U\leftarrow V}\\ \overline{\cal R}_{2}^{V\leftarrow I}&\overline{\cal R}_{2}^{V\leftarrow Q}&\overline{\cal R}_{2}^{V\leftarrow U}&\overline{\cal R}_{2}^{V\leftarrow V}\\ \end{array}\right]\left[\begin{array}[]{cccc}\mathbf{D}^{I}\\ \mathbf{D}^{Q}\\ \mathbf{D}^{U}\\ \mathbf{D}^{V}\\ \end{array}\right]$		(8q)
	$\displaystyle\left[\begin{array}[]{cccc}\mathbf{D}^{I}\\ \mathbf{D}^{Q}\\ \mathbf{D}^{U}\\ \mathbf{D}^{V}\\ \end{array}\right]=\left[\begin{array}[]{cccc}\mathbf{T}_{1}^{I}\\ \mathbf{T}_{1}^{Q}\\ \mathbf{T}_{1}^{U}\\ \mathbf{T}_{1}^{V}\\ \end{array}\right]+\left[\begin{array}[]{cccc}\overline{\cal R}_{1}^{,I\leftarrow I}&\overline{\cal R}_{1}^{,I\leftarrow Q}&\overline{\cal R}_{1}^{,I\leftarrow U}&\overline{\cal R}_{1}^{,I\leftarrow V}\\ \overline{\cal R}_{1}^{,Q\leftarrow I}&\overline{\cal R}_{1}^{,Q\leftarrow Q}&\overline{\cal R}_{1}^{,Q\leftarrow U}&\overline{\cal R}_{1}^{,Q\leftarrow V}\\ \overline{\cal R}_{1}^{,U\leftarrow I}&\overline{\cal R}_{1}^{,U\leftarrow Q}&\overline{\cal R}_{1}^{,U\leftarrow U}&\overline{\cal R}_{1}^{,U\leftarrow V}\\ \overline{\cal R}_{1}^{,V\leftarrow I}&\overline{\cal R}_{1}^{,V\leftarrow Q}&\overline{\cal R}_{1}^{,V\leftarrow U}&\overline{\cal R}_{1}^{,V\leftarrow V}\\ \end{array}\right]\left[\begin{array}[]{cccc}\mathbf{A}^{I}\\ \mathbf{A}^{Q}\\ \mathbf{A}^{U}\\ \mathbf{A}^{V}\\ \end{array}\right]$		(8ah)
	$\displaystyle\left[\begin{array}[]{cccc}\mathbf{R}^{I}_{1,2}\\ \mathbf{R}^{Q}_{1,2}\\ \mathbf{R}^{U}_{1,2}\\ \mathbf{R}^{V}_{1,2}\\ \end{array}\right]=\left[\begin{array}[]{cccc}\mathbf{R}_{1}^{I}\\ \mathbf{R}_{1}^{Q}\\ \mathbf{R}_{1}^{U}\\ \mathbf{R}_{1}^{V}\\ \end{array}\right]+\left[\begin{array}[]{cccc}\overline{\cal T}_{1}^{,I\leftarrow I}&\overline{\cal T}_{1}^{,I\leftarrow Q}&\overline{\cal T}_{1}^{,I\leftarrow U}&\overline{\cal T}_{1}^{,I\leftarrow V}\\ \overline{\cal T}_{1}^{,Q\leftarrow I}&\overline{\cal T}_{1}^{,Q\leftarrow Q}&\overline{\cal T}_{1}^{,Q\leftarrow U}&\overline{\cal T}_{1}^{,Q\leftarrow V}\\ \overline{\cal T}_{1}^{,U\leftarrow I}&\overline{\cal T}_{1}^{,U\leftarrow Q}&\overline{\cal T}_{1}^{,U\leftarrow U}&\overline{\cal T}_{1}^{,U\leftarrow V}\\ \overline{\cal T}_{1}^{,V\leftarrow I}&\overline{\cal T}_{1}^{,V\leftarrow Q}&\overline{\cal T}_{1}^{,V\leftarrow U}&\overline{\cal T}_{1}^{,V\leftarrow V}\\ \end{array}\right]\left[\begin{array}[]{cccc}\mathbf{A}^{I}\\ \mathbf{A}^{Q}\\ \mathbf{A}^{U}\\ \mathbf{A}^{V}\\ \end{array}\right]$		(8ay)
	$\displaystyle\left[\begin{array}[]{cccc}\mathbf{T}^{I}_{1,2}\\ \mathbf{T}^{Q}_{1,2}\\ \mathbf{T}^{U}_{1,2}\\ \mathbf{T}^{V}_{1,2}\\ \end{array}\right]=\left[\begin{array}[]{cccc}\mathbf{T}_{2}^{I}\\ \mathbf{T}_{2}^{Q}\\ \mathbf{T}_{2}^{U}\\ \mathbf{T}_{2}^{V}\\ \end{array}\right]e^{-\frac{\tau_{1}}{\mu_{0}}}+\left[\begin{array}[]{cccc}\overline{\cal T}_{2}^{I\leftarrow I}&\overline{\cal T}_{2}^{I\leftarrow Q}&\overline{\cal T}_{2}^{I\leftarrow U}&\overline{\cal T}_{2}^{I\leftarrow V}\\ \overline{\cal T}_{2}^{Q\leftarrow I}&\overline{\cal T}_{2}^{Q\leftarrow Q}&\overline{\cal T}_{2}^{Q\leftarrow U}&\overline{\cal T}_{2}^{Q\leftarrow V}\\ \overline{\cal T}_{2}^{U\leftarrow I}&\overline{\cal T}_{2}^{U\leftarrow Q}&\overline{\cal T}_{2}^{U\leftarrow U}&\overline{\cal T}_{2}^{U\leftarrow V}\\ \overline{\cal T}_{2}^{V\leftarrow I}&\overline{\cal T}_{2}^{V\leftarrow Q}&\overline{\cal T}_{2}^{V\leftarrow U}&\overline{\cal T}_{2}^{V\leftarrow V}\\ \end{array}\right]\left[\begin{array}[]{cccc}\mathbf{D}^{I}\\ \mathbf{D}^{Q}\\ \mathbf{D}^{U}\\ \mathbf{D}^{V}\\ \end{array}\right]$		(8bp)

From (8a-b), we obtain


	$\displaystyle\left[\begin{array}[]{cccc}\mathbf{A}^{I}\\ \mathbf{A}^{Q}\\ \mathbf{A}^{U}\\ \mathbf{A}^{V}\\ \end{array}\right]=\overline{\mathbf{X}}_{2}^{1}\mathbf{Y}_{2}^{1},$		(9e)
	$\displaystyle\left[\begin{array}[]{cccc}\mathbf{D}^{I}\\ \mathbf{D}^{Q}\\ \mathbf{D}^{U}\\ \mathbf{D}^{V}\\ \end{array}\right]=\left[\begin{array}[]{cccc}\mathbf{T}_{1}^{I}\\ \mathbf{T}_{1}^{Q}\\ \mathbf{T}_{1}^{U}\\ \mathbf{T}_{1}^{V}\\ \end{array}\right]+\left[\begin{array}[]{cccc}\overline{\cal R}_{1}^{,I\leftarrow I}&\overline{\cal R}_{1}^{,I\leftarrow Q}&\overline{\cal R}_{1}^{,I\leftarrow U}&\overline{\cal R}_{1}^{,I\leftarrow V}\\ \overline{\cal R}_{1}^{,Q\leftarrow I}&\overline{\cal R}_{1}^{,Q\leftarrow Q}&\overline{\cal R}_{1}^{,Q\leftarrow U}&\overline{\cal R}_{1}^{,Q\leftarrow V}\\ \overline{\cal R}_{1}^{,U\leftarrow I}&\overline{\cal R}_{1}^{,U\leftarrow Q}&\overline{\cal R}_{1}^{,U\leftarrow U}&\overline{\cal R}_{1}^{,U\leftarrow V}\\ \overline{\cal R}_{1}^{,V\leftarrow I}&\overline{\cal R}_{1}^{,V\leftarrow Q}&\overline{\cal R}_{1}^{,V\leftarrow U}&\overline{\cal R}_{1}^{,V\leftarrow V}\\ \end{array}\right]\overline{\mathbf{X}}_{2}^{1}\mathbf{Y}_{2}^{1}$		(9r)

where $\overline{\mathbf{X}}_{i}^{j}$ and $\mathbf{Y}_{i}^{j}$ are given in Appendix B.

Substituting (9) into (8c-d), the direct reflection and transmission can be expressed as follows:


	$\displaystyle\left[\begin{array}[]{cccc}\mathbf{R}^{I}_{1,2}\\ \mathbf{R}^{Q}_{1,2}\\ \mathbf{R}^{U}_{1,2}\\ \mathbf{R}^{V}_{1,2}\\ \end{array}\right]=\left[\begin{array}[]{cccc}\mathbf{R}_{1}^{I}\\ \mathbf{R}_{1}^{Q}\\ \mathbf{R}_{1}^{U}\\ \mathbf{R}_{1}^{V}\\ \end{array}\right]+\left[\begin{array}[]{cccc}\overline{\cal T}_{1}^{,I\leftarrow I}&\overline{\cal T}_{1}^{,I\leftarrow Q}&\overline{\cal T}_{1}^{,I\leftarrow U}&\overline{\cal T}_{1}^{,I\leftarrow V}\\ \overline{\cal T}_{1}^{,Q\leftarrow I}&\overline{\cal T}_{1}^{,Q\leftarrow Q}&\overline{\cal T}_{1}^{,Q\leftarrow U}&\overline{\cal T}_{1}^{,Q\leftarrow V}\\ \overline{\cal T}_{1}^{,U\leftarrow I}&\overline{\cal T}_{1}^{,U\leftarrow Q}&\overline{\cal T}_{1}^{,U\leftarrow U}&\overline{\cal T}_{1}^{,U\leftarrow V}\\ \overline{\cal T}_{1}^{,V\leftarrow I}&\overline{\cal T}_{1}^{,V\leftarrow Q}&\overline{\cal T}_{1}^{,V\leftarrow U}&\overline{\cal T}_{1}^{,V\leftarrow V}\\ \end{array}\right]\overline{\mathbf{X}}_{2}^{1}\mathbf{Y}_{2}^{1},$		(10m)
	$\displaystyle\left[\begin{array}[]{cccc}\mathbf{T}^{I}_{1,2}\\ \mathbf{T}^{Q}_{1,2}\\ \mathbf{T}^{U}_{1,2}\\ \mathbf{T}^{V}_{1,2}\\ \end{array}\right]=\left[\begin{array}[]{cccc}\mathbf{T}_{2}^{I}\\ \mathbf{T}_{2}^{Q}\\ \mathbf{T}_{2}^{U}\\ \mathbf{T}_{2}^{V}\\ \end{array}\right]e^{-\frac{\tau_{1}}{\mu_{0}}}+\left[\begin{array}[]{cccc}\overline{\cal T}_{2}^{I\leftarrow I}&\overline{\cal T}_{2}^{I\leftarrow Q}&\overline{\cal T}_{2}^{I\leftarrow U}&\overline{\cal T}_{2}^{I\leftarrow V}\\ \overline{\cal T}_{2}^{Q\leftarrow I}&\overline{\cal T}_{2}^{Q\leftarrow Q}&\overline{\cal T}_{2}^{Q\leftarrow U}&\overline{\cal T}_{2}^{Q\leftarrow V}\\ \overline{\cal T}_{2}^{U\leftarrow I}&\overline{\cal T}_{2}^{U\leftarrow Q}&\overline{\cal T}_{2}^{U\leftarrow U}&\overline{\cal T}_{2}^{U\leftarrow V}\\ \overline{\cal T}_{2}^{V\leftarrow I}&\overline{\cal T}_{2}^{V\leftarrow Q}&\overline{\cal T}_{2}^{V\leftarrow U}&\overline{\cal T}_{2}^{V\leftarrow V}\\ \end{array}\right]\left[\begin{array}[]{cccc}\mathbf{T}_{1}^{I}\\ \mathbf{T}_{1}^{Q}\\ \mathbf{T}_{1}^{U}\\ \mathbf{T}_{1}^{V}\\ \end{array}\right]+$		(10ad)
	$\displaystyle\quad\left[\begin{array}[]{cccc}\overline{\cal T}_{2}^{I\leftarrow I}&\overline{\cal T}_{2}^{I\leftarrow Q}&\overline{\cal T}_{2}^{I\leftarrow U}&\overline{\cal T}_{2}^{I\leftarrow V}\\ \overline{\cal T}_{2}^{Q\leftarrow I}&\overline{\cal T}_{2}^{Q\leftarrow Q}&\overline{\cal T}_{2}^{Q\leftarrow U}&\overline{\cal T}_{2}^{Q\leftarrow V}\\ \overline{\cal T}_{2}^{U\leftarrow I}&\overline{\cal T}_{2}^{U\leftarrow Q}&\overline{\cal T}_{2}^{U\leftarrow U}&\overline{\cal T}_{2}^{U\leftarrow V}\\ \overline{\cal T}_{2}^{V\leftarrow I}&\overline{\cal T}_{2}^{V\leftarrow Q}&\overline{\cal T}_{2}^{V\leftarrow U}&\overline{\cal T}_{2}^{V\leftarrow V}\\ \end{array}\right]\left[\begin{array}[]{cccc}\overline{\cal R}_{1}^{,I\leftarrow I}&\overline{\cal R}_{1}^{,I\leftarrow Q}&\overline{\cal R}_{1}^{,I\leftarrow U}&\overline{\cal R}_{1}^{,I\leftarrow V}\\ \overline{\cal R}_{1}^{,Q\leftarrow I}&\overline{\cal R}_{1}^{,Q\leftarrow Q}&\overline{\cal R}_{1}^{,Q\leftarrow U}&\overline{\cal R}_{1}^{,Q\leftarrow V}\\ \overline{\cal R}_{1}^{,U\leftarrow I}&\overline{\cal R}_{1}^{,U\leftarrow Q}&\overline{\cal R}_{1}^{,U\leftarrow U}&\overline{\cal R}_{1}^{,U\leftarrow V}\\ \overline{\cal R}_{1}^{,V\leftarrow I}&\overline{\cal R}_{1}^{,V\leftarrow Q}&\overline{\cal R}_{1}^{,V\leftarrow U}&\overline{\cal R}_{1}^{,V\leftarrow V}\\ \end{array}\right]\overline{\mathbf{X}}_{2}^{1}\mathbf{Y}_{2}^{1}$		(10am)

By applying the invariance principles to diffuse radiation, we obtain


	$\displaystyle\left[\begin{array}[]{cccc}\overline{\cal R}_{1,2}^{I\leftarrow I}&\overline{\cal R}_{1,2}^{I\leftarrow Q}&\overline{\cal R}_{1,2}^{I\leftarrow U}&\overline{\cal R}_{1,2}^{I\leftarrow V}\\ \overline{\cal R}_{1,2}^{Q\leftarrow I}&\overline{\cal R}_{1,2}^{Q\leftarrow Q}&\overline{\cal R}_{1,2}^{Q\leftarrow U}&\overline{\cal R}_{1,2}^{Q\leftarrow V}\\ \overline{\cal R}_{1,2}^{U\leftarrow I}&\overline{\cal R}_{1,2}^{U\leftarrow Q}&\overline{\cal R}_{1,2}^{U\leftarrow U}&\overline{\cal R}_{1,2}^{U\leftarrow V}\\ \overline{\cal R}_{1,2}^{V\leftarrow I}&\overline{\cal R}_{1,2}^{V\leftarrow Q}&\overline{\cal R}_{1,2}^{V\leftarrow U}&\overline{\cal R}_{1,2}^{V\leftarrow V}\\ \end{array}\right]=\left[\begin{array}[]{cccc}\overline{\cal R}_{1}^{I\leftarrow I}&\overline{\cal R}_{1}^{I\leftarrow Q}&\overline{\cal R}_{1}^{I\leftarrow U}&\overline{\cal R}_{1}^{I\leftarrow V}\\ \overline{\cal R}_{1}^{Q\leftarrow I}&\overline{\cal R}_{1}^{Q\leftarrow Q}&\overline{\cal R}_{1}^{Q\leftarrow U}&\overline{\cal R}_{1}^{Q\leftarrow V}\\ \overline{\cal R}_{1}^{U\leftarrow I}&\overline{\cal R}_{1}^{U\leftarrow Q}&\overline{\cal R}_{1}^{U\leftarrow U}&\overline{\cal R}_{1}^{U\leftarrow V}\\ \overline{\cal R}_{1}^{V\leftarrow I}&\overline{\cal R}_{1}^{V\leftarrow Q}&\overline{\cal R}_{1}^{V\leftarrow U}&\overline{\cal R}_{1}^{V\leftarrow V}\\ \end{array}\right]$		(11i)
	$\displaystyle\quad+\left[\begin{array}[]{cccc}\overline{\cal T}_{1}^{,I\leftarrow I}&\overline{\cal T}_{1}^{,I\leftarrow Q}&\overline{\cal T}_{1}^{,I\leftarrow U}&\overline{\cal T}_{1}^{,I\leftarrow V}\\ \overline{\cal T}_{1}^{,Q\leftarrow I}&\overline{\cal T}_{1}^{,Q\leftarrow Q}&\overline{\cal T}_{1}^{,Q\leftarrow U}&\overline{\cal T}_{1}^{,Q\leftarrow V}\\ \overline{\cal T}_{1}^{,U\leftarrow I}&\overline{\cal T}_{1}^{,U\leftarrow Q}&\overline{\cal T}_{1}^{,U\leftarrow U}&\overline{\cal T}_{1}^{,U\leftarrow V}\\ \overline{\cal T}_{1}^{,V\leftarrow I}&\overline{\cal T}_{1}^{,V\leftarrow Q}&\overline{\cal T}_{1}^{,V\leftarrow U}&\overline{\cal T}_{1}^{,V\leftarrow V}\\ \end{array}\right]\overline{\mathbf{X}}_{2}^{1}$		(11n)
	$\displaystyle\quad\times\left[\begin{array}[]{cccc}\overline{\cal R}_{2}^{I\leftarrow I}&\overline{\cal R}_{2}^{I\leftarrow Q}&\overline{\cal R}_{2}^{I\leftarrow U}&\overline{\cal R}_{2}^{I\leftarrow V}\\ \overline{\cal R}_{2}^{Q\leftarrow I}&\overline{\cal R}_{2}^{Q\leftarrow Q}&\overline{\cal R}_{2}^{Q\leftarrow U}&\overline{\cal R}_{2}^{Q\leftarrow V}\\ \overline{\cal R}_{2}^{U\leftarrow I}&\overline{\cal R}_{2}^{U\leftarrow Q}&\overline{\cal R}_{2}^{U\leftarrow U}&\overline{\cal R}_{2}^{U\leftarrow V}\\ \overline{\cal R}_{2}^{V\leftarrow I}&\overline{\cal R}_{2}^{V\leftarrow Q}&\overline{\cal R}_{2}^{V\leftarrow U}&\overline{\cal R}_{2}^{V\leftarrow V}\\ \end{array}\right]\left[\begin{array}[]{cccc}\overline{\cal T}_{1}^{I\leftarrow I}&\overline{\cal T}_{1}^{I\leftarrow Q}&\overline{\cal T}_{1}^{I\leftarrow U}&\overline{\cal T}_{1}^{I\leftarrow V}\\ \overline{\cal T}_{1}^{Q\leftarrow I}&\overline{\cal T}_{1}^{Q\leftarrow Q}&\overline{\cal T}_{1}^{Q\leftarrow U}&\overline{\cal T}_{1}^{Q\leftarrow V}\\ \overline{\cal T}_{1}^{U\leftarrow I}&\overline{\cal T}_{1}^{U\leftarrow Q}&\overline{\cal T}_{1}^{U\leftarrow U}&\overline{\cal T}_{1}^{U\leftarrow V}\\ \overline{\cal T}_{1}^{V\leftarrow I}&\overline{\cal T}_{1}^{V\leftarrow Q}&\overline{\cal T}_{1}^{V\leftarrow U}&\overline{\cal T}_{1}^{V\leftarrow V}\\ \end{array}\right],$		(11w)
	$\displaystyle\left[\begin{array}[]{cccc}\overline{\cal T}_{1,2}^{I\leftarrow I}&\overline{\cal T}_{1,2}^{I\leftarrow Q}&\overline{\cal T}_{1,2}^{I\leftarrow U}&\overline{\cal T}_{1,2}^{I\leftarrow V}\\ \overline{\cal T}_{1,2}^{Q\leftarrow I}&\overline{\cal T}_{1,2}^{Q\leftarrow Q}&\overline{\cal T}_{1,2}^{Q\leftarrow U}&\overline{\cal T}_{1,2}^{Q\leftarrow V}\\ \overline{\cal T}_{1,2}^{U\leftarrow I}&\overline{\cal T}_{1,2}^{U\leftarrow Q}&\overline{\cal T}_{1,2}^{U\leftarrow U}&\overline{\cal T}_{1,2}^{U\leftarrow V}\\ \overline{\cal T}_{1,2}^{V\leftarrow I}&\overline{\cal T}_{1,2}^{V\leftarrow Q}&\overline{\cal T}_{1,2}^{V\leftarrow U}&\overline{\cal T}_{1,2}^{V\leftarrow V}\\ \end{array}\right]=\left[\begin{array}[]{cccc}\overline{\cal T}_{2}^{I\leftarrow I}&\overline{\cal T}_{2}^{I\leftarrow Q}&\overline{\cal T}_{2}^{I\leftarrow U}&\overline{\cal T}_{2}^{I\leftarrow V}\\ \overline{\cal T}_{2}^{Q\leftarrow I}&\overline{\cal T}_{2}^{Q\leftarrow Q}&\overline{\cal T}_{2}^{Q\leftarrow U}&\overline{\cal T}_{2}^{Q\leftarrow V}\\ \overline{\cal T}_{2}^{U\leftarrow I}&\overline{\cal T}_{2}^{U\leftarrow Q}&\overline{\cal T}_{2}^{U\leftarrow U}&\overline{\cal T}_{2}^{U\leftarrow V}\\ \overline{\cal T}_{2}^{V\leftarrow I}&\overline{\cal T}_{2}^{V\leftarrow Q}&\overline{\cal T}_{2}^{V\leftarrow U}&\overline{\cal T}_{2}^{V\leftarrow V}\\ \end{array}\right]$		(11af)
	$\displaystyle\quad\times\left[\begin{array}[]{cccc}\overline{\cal T}_{1}^{I\leftarrow I}&\overline{\cal T}_{1}^{I\leftarrow Q}&\overline{\cal T}_{1}^{I\leftarrow U}&\overline{\cal T}_{1}^{I\leftarrow V}\\ \overline{\cal T}_{1}^{Q\leftarrow I}&\overline{\cal T}_{1}^{Q\leftarrow Q}&\overline{\cal T}_{1}^{Q\leftarrow U}&\overline{\cal T}_{1}^{Q\leftarrow V}\\ \overline{\cal T}_{1}^{U\leftarrow I}&\overline{\cal T}_{1}^{U\leftarrow Q}&\overline{\cal T}_{1}^{U\leftarrow U}&\overline{\cal T}_{1}^{U\leftarrow V}\\ \overline{\cal T}_{1}^{V\leftarrow I}&\overline{\cal T}_{1}^{V\leftarrow Q}&\overline{\cal T}_{1}^{V\leftarrow U}&\overline{\cal T}_{1}^{V\leftarrow V}\\ \end{array}\right]+\left[\begin{array}[]{cccc}\overline{\cal T}_{2}^{I\leftarrow I}&\overline{\cal T}_{2}^{I\leftarrow Q}&\overline{\cal T}_{2}^{I\leftarrow U}&\overline{\cal T}_{2}^{I\leftarrow V}\\ \overline{\cal T}_{2}^{Q\leftarrow I}&\overline{\cal T}_{2}^{Q\leftarrow Q}&\overline{\cal T}_{2}^{Q\leftarrow U}&\overline{\cal T}_{2}^{Q\leftarrow V}\\ \overline{\cal T}_{2}^{U\leftarrow I}&\overline{\cal T}_{2}^{U\leftarrow Q}&\overline{\cal T}_{2}^{U\leftarrow U}&\overline{\cal T}_{2}^{U\leftarrow V}\\ \overline{\cal T}_{2}^{V\leftarrow I}&\overline{\cal T}_{2}^{V\leftarrow Q}&\overline{\cal T}_{2}^{V\leftarrow U}&\overline{\cal T}_{2}^{V\leftarrow V}\\ \end{array}\right]$		(11ao)
	$\displaystyle\quad\times\left[\begin{array}[]{cccc}\overline{\cal T}_{1}^{,I\leftarrow I}&\overline{\cal T}_{1}^{,I\leftarrow Q}&\overline{\cal T}_{1}^{,I\leftarrow U}&\overline{\cal T}_{1}^{,I\leftarrow V}\\ \overline{\cal T}_{1}^{,Q\leftarrow I}&\overline{\cal T}_{1}^{,Q\leftarrow Q}&\overline{\cal T}_{1}^{,Q\leftarrow U}&\overline{\cal T}_{1}^{,Q\leftarrow V}\\ \overline{\cal T}_{1}^{,U\leftarrow I}&\overline{\cal T}_{1}^{,U\leftarrow Q}&\overline{\cal T}_{1}^{,U\leftarrow U}&\overline{\cal T}_{1}^{,U\leftarrow V}\\ \overline{\cal T}_{1}^{,V\leftarrow I}&\overline{\cal T}_{1}^{,V\leftarrow Q}&\overline{\cal T}_{1}^{,V\leftarrow U}&\overline{\cal T}_{1}^{,V\leftarrow V}\\ \end{array}\right]\overline{\mathbf{X}}_{2}^{1}\left[\begin{array}[]{cccc}\overline{\cal R}_{2}^{I\leftarrow I}&\overline{\cal R}_{2}^{I\leftarrow Q}&\overline{\cal R}_{2}^{I\leftarrow U}&\overline{\cal R}_{2}^{I\leftarrow V}\\ \overline{\cal R}_{2}^{Q\leftarrow I}&\overline{\cal R}_{2}^{Q\leftarrow Q}&\overline{\cal R}_{2}^{Q\leftarrow U}&\overline{\cal R}_{2}^{Q\leftarrow V}\\ \overline{\cal R}_{2}^{U\leftarrow I}&\overline{\cal R}_{2}^{U\leftarrow Q}&\overline{\cal R}_{2}^{U\leftarrow U}&\overline{\cal R}_{2}^{U\leftarrow V}\\ \overline{\cal R}_{2}^{V\leftarrow I}&\overline{\cal R}_{2}^{V\leftarrow Q}&\overline{\cal R}_{2}^{V\leftarrow U}&\overline{\cal R}_{2}^{V\leftarrow V}\\ \end{array}\right]$		(11ax)
	$\displaystyle\quad\times\left[\begin{array}[]{cccc}\overline{\cal T}_{1}^{I\leftarrow I}&\overline{\cal T}_{1}^{I\leftarrow Q}&\overline{\cal T}_{1}^{I\leftarrow U}&\overline{\cal T}_{1}^{I\leftarrow V}\\ \overline{\cal T}_{1}^{Q\leftarrow I}&\overline{\cal T}_{1}^{Q\leftarrow Q}&\overline{\cal T}_{1}^{Q\leftarrow U}&\overline{\cal T}_{1}^{Q\leftarrow V}\\ \overline{\cal T}_{1}^{U\leftarrow I}&\overline{\cal T}_{1}^{U\leftarrow Q}&\overline{\cal T}_{1}^{U\leftarrow U}&\overline{\cal T}_{1}^{U\leftarrow V}\\ \overline{\cal T}_{1}^{V\leftarrow I}&\overline{\cal T}_{1}^{V\leftarrow Q}&\overline{\cal T}_{1}^{V\leftarrow U}&\overline{\cal T}_{1}^{V\leftarrow V}\\ \end{array}\right]$		(11bc)

For a light beam incident from below, $\mathbf{R}^{*}_{1,2}(\mu_{0})$ , $\mathbf{T}^{*}_{1,2}(\mu_{0})$ , $\overline{\cal R}^{*}_{1,2}$ , and $\overline{\cal T}^{*}_{1,2}$ can be obtained in a similar manner.

In the downward path calculation, we apply it to multiple layers (from layer 1 to layer $k$ ).


	$\displaystyle\left[\begin{array}[]{cccc}\mathbf{T}^{I}_{1,k}\\ \mathbf{T}^{Q}_{1,k}\\ \mathbf{T}^{U}_{1,k}\\ \mathbf{T}^{V}_{1,k}\\ \end{array}\right]=\left[\begin{array}[]{cccc}\mathbf{T}_{k}^{I}\\ \mathbf{T}_{k}^{Q}\\ \mathbf{T}_{k}^{U}\\ \mathbf{T}_{k}^{V}\\ \end{array}\right]e^{-\frac{\tau_{1,k-1}}{\mu_{0}}}+\left[\begin{array}[]{cccc}\overline{\cal T}_{k}^{I\leftarrow I}&\overline{\cal T}_{k}^{I\leftarrow Q}&\overline{\cal T}_{k}^{I\leftarrow U}&\overline{\cal T}_{k}^{I\leftarrow V}\\ \overline{\cal T}_{k}^{Q\leftarrow I}&\overline{\cal T}_{k}^{Q\leftarrow Q}&\overline{\cal T}_{k}^{Q\leftarrow U}&\overline{\cal T}_{k}^{Q\leftarrow V}\\ \overline{\cal T}_{k}^{U\leftarrow I}&\overline{\cal T}_{k}^{U\leftarrow Q}&\overline{\cal T}_{k}^{U\leftarrow U}&\overline{\cal T}_{k}^{U\leftarrow V}\\ \overline{\cal T}_{k}^{V\leftarrow I}&\overline{\cal T}_{k}^{V\leftarrow Q}&\overline{\cal T}_{k}^{V\leftarrow U}&\overline{\cal T}_{k}^{V\leftarrow V}\\ \end{array}\right]\left[\begin{array}[]{cccc}\mathbf{T}_{1,k-1}^{I}\\ \mathbf{T}_{1,k-1}^{Q}\\ \mathbf{T}_{1,k-1}^{U}\\ \mathbf{T}_{1,k-1}^{V}\\ \end{array}\right]+$		(12q)
	$\displaystyle\quad\left[\begin{array}[]{cccc}\overline{\cal T}_{k}^{I\leftarrow I}&\overline{\cal T}_{k}^{I\leftarrow Q}&\overline{\cal T}_{k}^{I\leftarrow U}&\overline{\cal T}_{k}^{I\leftarrow V}\\ \overline{\cal T}_{k}^{Q\leftarrow I}&\overline{\cal T}_{k}^{Q\leftarrow Q}&\overline{\cal T}_{k}^{Q\leftarrow U}&\overline{\cal T}_{k}^{Q\leftarrow V}\\ \overline{\cal T}_{k}^{U\leftarrow I}&\overline{\cal T}_{k}^{U\leftarrow Q}&\overline{\cal T}_{k}^{U\leftarrow U}&\overline{\cal T}_{k}^{U\leftarrow V}\\ \overline{\cal T}_{k}^{V\leftarrow I}&\overline{\cal T}_{k}^{V\leftarrow Q}&\overline{\cal T}_{k}^{V\leftarrow U}&\overline{\cal T}_{k}^{V\leftarrow V}\\ \end{array}\right]\left[\begin{array}[]{cccc}\overline{\cal R}_{1,k-1}^{,I\leftarrow I}&\overline{\cal R}_{1,k-1}^{,I\leftarrow Q}&\overline{\cal R}_{1,k-1}^{,I\leftarrow U}&\overline{\cal R}_{1,k-1}^{,I\leftarrow V}\\ \overline{\cal R}_{1,k-1}^{,Q\leftarrow I}&\overline{\cal R}_{1,k-1}^{,Q\leftarrow Q}&\overline{\cal R}_{1,k-1}^{,Q\leftarrow U}&\overline{\cal R}_{1,k-1}^{,Q\leftarrow V}\\ \overline{\cal R}_{1,k-1}^{,U\leftarrow I}&\overline{\cal R}_{1,k-1}^{,U\leftarrow Q}&\overline{\cal R}_{1,k-1}^{,U\leftarrow U}&\overline{\cal R}_{1,k-1}^{,U\leftarrow V}\\ \overline{\cal R}_{1,k-1}^{,V\leftarrow I}&\overline{\cal R}_{1,k-1}^{,V\leftarrow Q}&\overline{\cal R}_{1,k-1}^{,V\leftarrow U}&\overline{\cal R}_{1,k-1}^{,V\leftarrow V}\\ \end{array}\right]\overline{\mathbf{X}}_{k}^{1,k-1}\mathbf{Y}_{k}^{1.k-1}$		(12z)
	$\displaystyle\left[\begin{array}[]{cccc}\overline{\cal R}_{1,k}^{,I\leftarrow I}&\overline{\cal R}_{1,k}^{,I\leftarrow Q}&\overline{\cal R}_{1,k}^{,I\leftarrow U}&\overline{\cal R}_{1,k}^{,I\leftarrow V}\\ \overline{\cal R}_{1,k}^{,Q\leftarrow I}&\overline{\cal R}_{1,k}^{,Q\leftarrow Q}&\overline{\cal R}_{1,k}^{,Q\leftarrow U}&\overline{\cal R}_{1,k}^{,Q\leftarrow V}\\ \overline{\cal R}_{1,k}^{,U\leftarrow I}&\overline{\cal R}_{1,k}^{,U\leftarrow Q}&\overline{\cal R}_{1,k}^{,U\leftarrow U}&\overline{\cal R}_{1,k}^{,U\leftarrow V}\\ \overline{\cal R}_{1,k}^{,V\leftarrow I}&\overline{\cal R}_{1,k}^{,V\leftarrow Q}&\overline{\cal R}_{1,k}^{,V\leftarrow U}&\overline{\cal R}_{1,k}^{,V\leftarrow V}\\ \end{array}\right]=\left[\begin{array}[]{cccc}\overline{\cal R}_{k}^{,I\leftarrow I}&\overline{\cal R}_{k}^{,I\leftarrow Q}&\overline{\cal R}_{k}^{,I\leftarrow U}&\overline{\cal R}_{k}^{,I\leftarrow V}\\ \overline{\cal R}_{k}^{,Q\leftarrow I}&\overline{\cal R}_{k}^{,Q\leftarrow Q}&\overline{\cal R}_{k}^{,Q\leftarrow U}&\overline{\cal R}_{k}^{,Q\leftarrow V}\\ \overline{\cal R}_{k}^{,U\leftarrow I}&\overline{\cal R}_{k}^{,U\leftarrow Q}&\overline{\cal R}_{k}^{,U\leftarrow U}&\overline{\cal R}_{k}^{,U\leftarrow V}\\ \overline{\cal R}_{k}^{,V\leftarrow I}&\overline{\cal R}_{k}^{,V\leftarrow Q}&\overline{\cal R}_{k}^{,V\leftarrow U}&\overline{\cal R}_{k}^{,V\leftarrow V}\\ \end{array}\right]$		(12ai)
	$\displaystyle\quad+\left[\begin{array}[]{cccc}\overline{\cal T}_{k}^{I\leftarrow I}&\overline{\cal T}_{k}^{I\leftarrow Q}&\overline{\cal T}_{k}^{I\leftarrow U}&\overline{\cal T}_{k}^{I\leftarrow V}\\ \overline{\cal T}_{k}^{Q\leftarrow I}&\overline{\cal T}_{k}^{Q\leftarrow Q}&\overline{\cal T}_{k}^{Q\leftarrow U}&\overline{\cal T}_{k}^{Q\leftarrow V}\\ \overline{\cal T}_{k}^{U\leftarrow I}&\overline{\cal T}_{k}^{U\leftarrow Q}&\overline{\cal T}_{k}^{U\leftarrow U}&\overline{\cal T}_{k}^{U\leftarrow V}\\ \overline{\cal T}_{k}^{V\leftarrow I}&\overline{\cal T}_{k}^{V\leftarrow Q}&\overline{\cal T}_{k}^{V\leftarrow U}&\overline{\cal T}_{k}^{V\leftarrow V}\\ \end{array}\right]\overline{\mathbf{X}}_{k}^{1,k-1}\left[\begin{array}[]{cccc}\overline{\cal R}_{1,k-1}^{,I\leftarrow I}&\overline{\cal R}_{1,k-1}^{,I\leftarrow Q}&\overline{\cal R}_{1,k-1}^{,I\leftarrow U}&\overline{\cal R}_{1,k-1}^{,I\leftarrow V}\\ \overline{\cal R}_{1,k-1}^{,Q\leftarrow I}&\overline{\cal R}_{1,k-1}^{,Q\leftarrow Q}&\overline{\cal R}_{1,k-1}^{,Q\leftarrow U}&\overline{\cal R}_{1,k-1}^{,Q\leftarrow V}\\ \overline{\cal R}_{1,k-1}^{,U\leftarrow I}&\overline{\cal R}_{1,k-1}^{,U\leftarrow Q}&\overline{\cal R}_{1,k-1}^{,U\leftarrow U}&\overline{\cal R}_{1,k-1}^{,U\leftarrow V}\\ \overline{\cal R}_{1,k-1}^{,V\leftarrow I}&\overline{\cal R}_{1,k-1}^{,V\leftarrow Q}&\overline{\cal R}_{1,k-1}^{,V\leftarrow U}&\overline{\cal R}_{1,k-1}^{*,V\leftarrow V}\\ \end{array}\right]$		(12ar)
	$\displaystyle\quad\times\left[\begin{array}[]{cccc}\overline{\cal T}_{k}^{,I\leftarrow I}&\overline{\cal T}_{k}^{,I\leftarrow Q}&\overline{\cal T}_{k}^{,I\leftarrow U}&\overline{\cal T}_{k}^{,I\leftarrow V}\\ \overline{\cal T}_{k}^{,Q\leftarrow I}&\overline{\cal T}_{k}^{,Q\leftarrow Q}&\overline{\cal T}_{k}^{,Q\leftarrow U}&\overline{\cal T}_{k}^{,Q\leftarrow V}\\ \overline{\cal T}_{k}^{,U\leftarrow I}&\overline{\cal T}_{k}^{,U\leftarrow Q}&\overline{\cal T}_{k}^{,U\leftarrow U}&\overline{\cal T}_{k}^{,U\leftarrow V}\\ \overline{\cal T}_{k}^{,V\leftarrow I}&\overline{\cal T}_{k}^{,V\leftarrow Q}&\overline{\cal T}_{k}^{,V\leftarrow U}&\overline{\cal T}_{k}^{,V\leftarrow V}\\ \end{array}\right]$		(12aw)

where $\overline{\mathbf{X}}_{i}^{*j}$ is provided in Appendix B.

Similarly, in an upward path calculation, we apply the multi-layer from the surface to layer $k+1$ :


	$\displaystyle\left[\begin{array}[]{cccc}\mathbf{R}^{I}_{k,N}\\ \mathbf{R}^{Q}_{k,N}\\ \mathbf{R}^{U}_{k,N}\\ \mathbf{R}^{V}_{k,N}\\ \end{array}\right]=\left[\begin{array}[]{cccc}\mathbf{R}_{k}^{I}\\ \mathbf{R}_{k}^{Q}\\ \mathbf{R}_{k}^{U}\\ \mathbf{R}_{k}^{V}\\ \end{array}\right]+\left[\begin{array}[]{cccc}\overline{\cal T}_{k}^{,I\leftarrow I}&\overline{\cal T}_{k}^{,I\leftarrow Q}&\overline{\cal T}_{k}^{,I\leftarrow U}&\overline{\cal T}_{k}^{,I\leftarrow V}\\ \overline{\cal T}_{k}^{,Q\leftarrow I}&\overline{\cal T}_{k}^{,Q\leftarrow Q}&\overline{\cal T}_{k}^{,Q\leftarrow U}&\overline{\cal T}_{k}^{,Q\leftarrow V}\\ \overline{\cal T}_{k}^{,U\leftarrow I}&\overline{\cal T}_{k}^{,U\leftarrow Q}&\overline{\cal T}_{k}^{,U\leftarrow U}&\overline{\cal T}_{k}^{,U\leftarrow V}\\ \overline{\cal T}_{k}^{,V\leftarrow I}&\overline{\cal T}_{k}^{,V\leftarrow Q}&\overline{\cal T}_{k}^{,V\leftarrow U}&\overline{\cal T}_{k}^{,V\leftarrow V}\\ \end{array}\right]\overline{\mathbf{X}}_{k}^{k+1,N}\mathbf{Y}_{k}^{k+1,N}$		(13m)
	$\displaystyle\left[\begin{array}[]{cccc}\overline{\cal R}_{k,N}^{I\leftarrow I}&\overline{\cal R}_{k,N}^{I\leftarrow Q}&\overline{\cal R}_{k,N}^{I\leftarrow U}&\overline{\cal R}_{k,N}^{I\leftarrow V}\\ \overline{\cal R}_{k,N}^{Q\leftarrow I}&\overline{\cal R}_{k,N}^{Q\leftarrow Q}&\overline{\cal R}_{k,N}^{Q\leftarrow U}&\overline{\cal R}_{k,N}^{Q\leftarrow V}\\ \overline{\cal R}_{k,N}^{U\leftarrow I}&\overline{\cal R}_{k,N}^{U\leftarrow Q}&\overline{\cal R}_{k,N}^{U\leftarrow U}&\overline{\cal R}_{k,N}^{U\leftarrow V}\\ \overline{\cal R}_{k,N}^{V\leftarrow I}&\overline{\cal R}_{k,N}^{V\leftarrow Q}&\overline{\cal R}_{k,N}^{V\leftarrow U}&\overline{\cal R}_{k,N}^{V\leftarrow V}\\ \end{array}\right]=\left[\begin{array}[]{cccc}\overline{\cal R}_{k}^{I\leftarrow I}&\overline{\cal R}_{k}^{I\leftarrow Q}&\overline{\cal R}_{k}^{I\leftarrow U}&\overline{\cal R}_{k}^{I\leftarrow V}\\ \overline{\cal R}_{k}^{Q\leftarrow I}&\overline{\cal R}_{k}^{Q\leftarrow Q}&\overline{\cal R}_{k}^{Q\leftarrow U}&\overline{\cal R}_{k}^{Q\leftarrow V}\\ \overline{\cal R}_{k}^{U\leftarrow I}&\overline{\cal R}_{k}^{U\leftarrow Q}&\overline{\cal R}_{k}^{U\leftarrow U}&\overline{\cal R}_{k}^{U\leftarrow V}\\ \overline{\cal R}_{k}^{V\leftarrow I}&\overline{\cal R}_{k}^{V\leftarrow Q}&\overline{\cal R}_{k}^{V\leftarrow U}&\overline{\cal R}_{k}^{V\leftarrow V}\\ \end{array}\right]$		(13v)
	$\displaystyle\quad+\left[\begin{array}[]{cccc}\overline{\cal T}_{k}^{,I\leftarrow I}&\overline{\cal T}_{k}^{,I\leftarrow Q}&\overline{\cal T}_{k}^{,I\leftarrow U}&\overline{\cal T}_{k}^{,I\leftarrow V}\\ \overline{\cal T}_{k}^{,Q\leftarrow I}&\overline{\cal T}_{k}^{,Q\leftarrow Q}&\overline{\cal T}_{k}^{,Q\leftarrow U}&\overline{\cal T}_{k}^{,Q\leftarrow V}\\ \overline{\cal T}_{k}^{,U\leftarrow I}&\overline{\cal T}_{k}^{,U\leftarrow Q}&\overline{\cal T}_{k}^{,U\leftarrow U}&\overline{\cal T}_{k}^{,U\leftarrow V}\\ \overline{\cal T}_{k}^{,V\leftarrow I}&\overline{\cal T}_{k}^{,V\leftarrow Q}&\overline{\cal T}_{k}^{,V\leftarrow U}&\overline{\cal T}_{k}^{,V\leftarrow V}\\ \end{array}\right]\overline{\mathbf{X}}_{k}^{k+1,N}\left[\begin{array}[]{cccc}\overline{\cal R}_{k+1,N}^{I\leftarrow I}&\overline{\cal R}_{k+1,N}^{I\leftarrow Q}&\overline{\cal R}_{k+1,N}^{I\leftarrow U}&\overline{\cal R}_{k+1,N}^{I\leftarrow V}\\ \overline{\cal R}_{k+1,N}^{Q\leftarrow I}&\overline{\cal R}_{k+1,N}^{Q\leftarrow Q}&\overline{\cal R}_{k+1,N}^{Q\leftarrow U}&\overline{\cal R}_{k+1,N}^{Q\leftarrow V}\\ \overline{\cal R}_{k+1,N}^{U\leftarrow I}&\overline{\cal R}_{k+1,N}^{U\leftarrow Q}&\overline{\cal R}_{k+1,N}^{U\leftarrow U}&\overline{\cal R}_{k+1,N}^{U\leftarrow V}\\ \overline{\cal R}_{k+1,N}^{V\leftarrow I}&\overline{\cal R}_{k+1,N}^{V\leftarrow Q}&\overline{\cal R}_{k+1,N}^{V\leftarrow U}&\overline{\cal R}_{k+1,N}^{V\leftarrow V}\\ \end{array}\right]$		(13ae)
	$\displaystyle\quad\times\left[\begin{array}[]{cccc}\overline{\cal T}_{k}^{I\leftarrow I}&\overline{\cal T}_{k}^{I\leftarrow Q}&\overline{\cal T}_{k}^{I\leftarrow U}&\overline{\cal T}_{k}^{I\leftarrow V}\\ \overline{\cal T}_{k}^{Q\leftarrow I}&\overline{\cal T}_{k}^{Q\leftarrow Q}&\overline{\cal T}_{k}^{Q\leftarrow U}&\overline{\cal T}_{k}^{Q\leftarrow V}\\ \overline{\cal T}_{k}^{U\leftarrow I}&\overline{\cal T}_{k}^{U\leftarrow Q}&\overline{\cal T}_{k}^{U\leftarrow U}&\overline{\cal T}_{k}^{U\leftarrow V}\\ \overline{\cal T}_{k}^{V\leftarrow I}&\overline{\cal T}_{k}^{V\leftarrow Q}&\overline{\cal T}_{k}^{V\leftarrow U}&\overline{\cal T}_{k}^{V\leftarrow V}\\ \end{array}\right]$		(13aj)

where $\mathbf{R}_{N}=\left[\begin{array}[]{cccc}\textbf{R}_{N}^{I}\\ \textbf{R}_{N}^{Q}\\ \textbf{R}_{N}^{U}\\ \textbf{R}_{N}^{V}\\ \end{array}\right]$ is the surface albedo and $\overline{\cal R}_{N}=\left[\begin{array}[]{cccc}\overline{\boldsymbol{\mathcal{R}}}_{N}^{I\leftarrow I}&\overline{\boldsymbol{\mathcal{R}}}_{N}^{I\leftarrow Q}&\overline{\boldsymbol{\mathcal{R}}}_{N}^{I\leftarrow U}&\overline{\boldsymbol{\mathcal{R}}}_{N}^{I\leftarrow V}\\ \overline{\boldsymbol{\mathcal{R}}}_{N}^{Q\leftarrow I}&\overline{\boldsymbol{\mathcal{R}}}_{N}^{Q\leftarrow Q}&\overline{\boldsymbol{\mathcal{R}}}_{N}^{Q\leftarrow U}&\overline{\boldsymbol{\mathcal{R}}}_{N}^{Q\leftarrow V}\\ \overline{\boldsymbol{\mathcal{R}}}_{N}^{U\leftarrow I}&\overline{\boldsymbol{\mathcal{R}}}_{N}^{U\leftarrow Q}&\overline{\boldsymbol{\mathcal{R}}}_{N}^{U\leftarrow U}&\overline{\boldsymbol{\mathcal{R}}}_{N}^{U\leftarrow V}\\ \overline{\boldsymbol{\mathcal{R}}}_{N}^{V\leftarrow I}&\overline{\boldsymbol{\mathcal{R}}}_{N}^{V\leftarrow Q}&\overline{\boldsymbol{\mathcal{R}}}_{N}^{V\leftarrow U}&\overline{\boldsymbol{\mathcal{R}}}_{N}^{V\leftarrow V}\\ \end{array}\right]$ is the reflection property of the surface.

Thus, the internal Stokes parameters at level $k+1$ can be determined as follows:


	$\displaystyle\left[\begin{array}[]{cccc}\mathbf{A}^{I}_{k+1}\\ \mathbf{A}^{Q}_{k+1}\\ \mathbf{A}^{U}_{k+1}\\ \mathbf{A}^{V}_{k+1}\\ \end{array}\right]=\overline{\mathbf{X}}_{1,k}^{k+1,N}\mathbf{Y}_{1,k}^{k+1,N},$		(14e)
	$\displaystyle\left[\begin{array}[]{cccc}\mathbf{D}^{I}_{k+1}\\ \mathbf{D}^{Q}_{k+1}\\ \mathbf{D}^{U}_{k+1}\\ \mathbf{D}^{V}_{k+1}\\ \end{array}\right]=\left[\begin{array}[]{cccc}\mathbf{T}_{1,k}^{I}\\ \mathbf{T}_{1,k}^{Q}\\ \mathbf{T}_{1,k}^{U}\\ \mathbf{T}_{1,k}^{V}\\ \end{array}\right]+\left[\begin{array}[]{cccc}\overline{\cal R}_{1,k}^{,I\leftarrow I}&\overline{\cal R}_{1,k}^{,I\leftarrow Q}&\overline{\cal R}_{1,k}^{,I\leftarrow U}&\overline{\cal R}_{1,k}^{,I\leftarrow V}\\ \overline{\cal R}_{1,k}^{,Q\leftarrow I}&\overline{\cal R}_{1,k}^{,Q\leftarrow Q}&\overline{\cal R}_{1,k}^{,Q\leftarrow U}&\overline{\cal R}_{1,k}^{,Q\leftarrow V}\\ \overline{\cal R}_{1,k}^{,U\leftarrow I}&\overline{\cal R}_{1,k}^{,U\leftarrow Q}&\overline{\cal R}_{1,k}^{,U\leftarrow U}&\overline{\cal R}_{1,k}^{,U\leftarrow V}\\ \overline{\cal R}_{1,k}^{,V\leftarrow I}&\overline{\cal R}_{1,k}^{,V\leftarrow Q}&\overline{\cal R}_{1,k}^{,V\leftarrow U}&\overline{\cal R}_{1,k}^{,V\leftarrow V}\\ \end{array}\right]\overline{\mathbf{X}}_{1,k}^{k+1,N}\mathbf{Y}_{1,k}^{k+1,N}.$		(14r)

An adjusted absorption and scattering atmosphere may be considered to incorporate the forward peak contribution into multiple scattering. The fraction of the scattered energy residing in the forward peak $f$ is separated from the phase function. This is the $\delta$ -M adjustment [38] used in the proposed model.

The intensity at an arbitrary zenith angle for satellite applications can be obtained by replacing $\mu_{i}$ in Eq.(4) and the satellite zenith angle $\mu_{sat}$ . The entire calculation process is the same for $I(\tau,\pm\mu_{sat})$ , $Q(\tau,\pm\mu_{sat})$ , $U(\tau,\pm\mu_{sat})$ , and $V(\tau,\pm\mu_{sat})$ . We call this new method the polarized discrete ordinate adding approximation (POLDDA) for vector solar radiative transfer.

In the PolRadtran/RT3 model, which is a code based on the adding-doubling method, the adding method is always used together with the doubling method because the multiple scattering process is not considered in the single-layer solution and the calculation time increases with an increase in the optical thickness. However, in POLDDA, no doubling process is required and the adding method can be used directly after a single-layer solution (considering multiple scattering processes). Therefore, in theory, the calculation time does not change with the variation in the optical depth. This point is demonstrated in Section 4.3.

4 Comparison Results

In this section, the calculation accuracy and efficiency of POLDDA are evaluated by comparing with the Monte Carlo model (MYSTIC) [39] and PolRadtran/RT3, where the MYSTIC model is considered the benchmark. Test cases for both single-layer and multi-layer atmospheres are included. POLDDA and RT3 were executed with 16 streams (half-sphere) for all Rayleigh cases, and 32 streams/64 streams (half-sphere) for the water cloud case.

4.1 Single-layer test cases

In the first set of cases, one layer including molecules of atmospheric constituents with different depolarization factors and no surface reflection was tested. Another case exists, which includes Lambertian surface reflection.

4.1.1 Case 1 - Rayleigh scattering with different depolarization factor

Case 1 (including cases 1-1, 1-2, 1-3) was used to verify the correct treatment of the anisotropy of the molecules through the Rayleigh depolarization factor. It contains one layer of non-absorbing molecules (single scattering albedo $\omega=1$ ) and an optical thickness of 0.5. Surface reflection was not considered (the surface albedo was 0). Stokes vectors were calculated for various depolarization factors and sun positions. The definition of the viewing zenith angle is with respect to the downward normal instead of the upward normal, so viewing zenith angles are 0^∘ -80^∘ (down-looking) at the bottom and 100^∘-180^∘ (up-looking) at the top with 5^∘ increments. The results at viewing zenith angle near 90^∘ (direction near horizontal) are not shown. The viewing azimuthal angle definition is clockwise and ranged from 0^∘ to 360^∘ (5^∘ increments).

Test case 1-1 had the simplest setup with a zero Rayleigh depolarization factor. The Stokes vectors at the bottom and TOA were calculated for solar zenith angles $\theta_{0}=$ 0^∘ and $\phi_{0}=$ 65^∘, respectively. The left plots shown in Figure 2 present the results calculated using POLDDA. The absolute errors and relative differences between POLDDA and RT3 compared to the MYSTIC model are shown. The $U$ -component is zero because the sun was located at the zenith. From Fig. 2, the relative differences by POLDDA are smaller than 0.1% for the $I$ -component, which is almost the same as that of RT3 in both the upward and downward directions. The related bias are almost less than 2% for the $Q$ -component at the top and bottom of the atmosphere using both POLDDA and RT3. In one point, the abnormally large related bias is shown at viewing azimuth angles of 180^∘ (top) and 0^∘ (bottom). It is because the values of $Q$ -component are unsuitable as a denominator for calculating relative error if they are near zero. From the absolute error figures, a large error does not exist, which is a good illustration. When we removed abnormally large relative differences from Fig.2, the relative difference plot was consistent with the absolute difference plot (as shown in Fig. 2). The same phenomenon also appeared in subsequent cases, but we did not remove the deviation again. The relative root mean square errors (RMSE) for all the test cases are listed in Table 1. We calculated the RMSE for the $I$ -component between MYSTIC and POLDDA or RT3 using RMSE= $\frac{\sqrt{\sum_{i=1}^{N}(I_{MYSTIC}^{i}-I_{o}^{i})^{2}}}{\sqrt{\sum_{i=1}^{N}(I_{MYSTIC}^{i})^{2}}}$ (where $I_{o}$ refers to the other models) for the radiation field including all down-and up-looking directions. The RMSEs of POLDDA and RT3 are 0.0177% for the $I$ -component and 0.0226% for the $Q$ -component, which also demonstrates the accuracy of POLDDA in this case.

Test cases 1-2 and 1-3 consider two cases in which the depolarization factor is not zero; it is 0.03 in case 1-2 and 0.1 in case 1-3. The solar positions ( $\theta_{0}=$ 30^∘, $\phi_{0}$ =65^∘) are the same in cases 1-2 and 1-3. The $U$ -component is nonzero because the solar zenith angle is not zero. Figures 3 and 4 present the Stokes vector results at the TOA and bottom, respectively. We show the relative bias, and as shown in Fig. 3, the relative bias for the $I$ -component of POLDDA is between -0.05% and 0.05%. The largest bias of the $Q$ -component reached a high of 2%, ranging from 120^∘-150^∘ and near 170^∘ and 285^∘ viewing azimuth angles, for both POLDDA and RT3. The $U$ -component calculated by POLDDA and RT3 did not show a significant difference, except near viewing azimuth angles of 0^∘ and 180^∘ owing to the $U$ -component near zero. The large bias of the $I$ , $Q$ , and $U$ components yielded by POLDDA at the bottom (Fig. 4) is almost the same as that of RT3. For Case 1-2, the RMSEs of POLDDA were close to those of RT3 (Table 1).

In Case 1-3, although the value of the depolarization factor increased, the bias of POLDDA (Figures 5 and 6) also differed slightly from that of RT3. The RMSE of POLDDA are 0.0166%, 0.025%, and 0.0229% (Table 1) for the $I$ -, $Q$ -, and $U$ -components, respectively, which are similar to those of RT3.

4.1.2 Case 2 - Rayleigh atmosphere with Lambertian surface

Surface reflection was considered in Case 2 with a Lambertian surface albedo of 0.3, as well as for a one-layer medium of 0.1 without absorption. The optical depth of the layer was set as 0.1. The solar position was $\theta_{0}=$ 50^∘, with $\phi_{0}$ =0^∘, and the Rayleigh depolarization factor is the same as that in Case 1-3. The Stokes vectors ( $I$ -, $Q$ -, $U$ -, and $V$ -components) were calculated for viewing zenith angles of 0-80^∘ (100-180^∘) at the bottom (top) with 5^∘ increments, and viewing azimuth angles of 0-180^∘ with 5^∘ increments.

Figure 7 shows the results of Case 2. The relative bias of POLDDA is less than 0.05% for the $I$ -component. A relative bias greater than 1% occurs when viewing azimuth angles are approximately 170 degrees at the top of the layer for the $Q$ -component and approximately 30 degrees at the bottom. For the $U$ -component, over 5% relative bias is observed in the viewing azimuth angles ranging from 0 to 10^∘ at the bottom and 100^∘-110^∘ at the top of the layer. Significant differences were not found in the $I$ -, $Q$ -, and $U$ -components calculated using POLDDA and RT3. The RMSEs of POLDDA and RT3 were also similar to those of nonzero Lambertian surface reflection.

4.2 Test cases for realistic atmosphere with multi-layer

In this set of cases, the coupling between POLDDA model layers was examined using a multi-layer in a plane-parallel atmosphere. We used the addition method to handle multi-layer connections. Unlike RT3, the doubling process for POLDDA is not required; only the addition process is used with single-layer reflection and transmission functions.

The U.S. standard atmosphere [40] was used for the multi-layer cases. The model atmosphere comprised 30 layers from 0 km to 30 km in an altitude with a thickness of 1 km. The surface albedo was set to 0 for Cases 3 and 4 and 0.2 for Case 5. The Rayleigh scattering depolarization factor is 0.03. The sun?¡¥s position is $\theta_{0}=$ 60^∘ and $\phi_{0}=$ 0^∘. For Cases 3 and 4, the radiance field was calculated at the surface (0 km) and TOA (30 km) for viewing zenith angles from 0^∘ to 85^∘ (up-looking), and 95^∘ to 180^∘ (down-looking) in 5^∘ increments. The viewing azimuth angle ranged from 0^∘ to 180^∘ in 5^∘ increments. For Case 5, the radiance field was calculated at the surface (0 km) and TOA (30 km) for viewing zenith angles from 0^∘ to 90^∘ (up-looking) and 90^∘ to 180^∘ (down-looking) at a 90^∘ viewing azimuth angle.

4.2.1 Case 3 - Only Rayleigh scattering with standard atmosphere

The radiance field was calculated as 450 nm. Only Rayleigh scattering was considered. The profiles of the scattering optical thicknesses are depicted in Figure 8(a) and can be downloaded from the International Polarized Radiative Transfer (IPRT) website.

Figure 9 shows the $I$ -, $Q$ -, and $U$ -components at the output altitudes for a multi-layer atmosphere. For the one-layer cases, good agreement between the three models was observed for pure Rayleigh scattering with zero depolarization factor. In the inhomogeneous multi-layer, the three models are in good agreement, except near 0^∘-10^∘ at the bottom and 100^∘-110^∘ viewing azimuth angles at TOA for the $U$ -component. The bias of POLDDA is less than 0.2% for the $I$ -component and 0.5% for the $Q$ -component. The RMSE of POLDDA for the $I$ -, $Q$ -, and $U$ -components are only 0.017%, 0.0266%, and 0.0198%, respectively, against MYSTIC (Table 1).

4.2.2 Case 4 - Rayleigh scattering and absorption with standard atmosphere

In this case, it was determined whether the absorption effect was appropriately considered. The radiance field was calculated at a wavelength of 325 nm, where the ozone absorption is strong. The profiles of scattering optical thicknesses and absorption optical thicknesses are illustrated in Figures 8 (b-c), which can also be downloaded from the IPRT website.

Figure 10 shows the $I$ -, $Q$ -, and $U$ -component results at the top and on the surface. Similar to Case 3, little difference exists as calculated by POLDDA and RT3 for $I$ -, $Q$ -, and $U$ -components. The difference between RMSE of POLDDA and RT3 decreases compared to the values in Case 4 owing to the addition of gas absorption, particularly for the $I$ -component.

4.2.3 Case 5 - Standard atmosphere with cloud layer

In this case, a cloud layer was added between 2 and 3 km. The radiance was measured at 800 nm. The cloud optical thickness is 5, which is larger than the atmosphere that only includes air molecules. The profiles of the Rayleigh scattering optical thicknesses with no molecular absorption are shown in Figure 8(d). The single-scattering albedo for the cloud layer is 0.999979 and the cloud phase function was calculated using Mie scattering. POLDDA and RT3 were executed with 32 and 64 streams (half-sphere), respectively.

Figures 11 and 12 show the $I$ -, $Q$ -, $U$ -, and $V$ -component results for RT3 (64 streams), and the relative differences between POLDDA and RT3 at the top and surface, respectively. The relative differences between POLDDA and RT3 are less than 0.01%, 0.025%, 0.05%, and 0.05% for the $I$ -, $Q$ -, $U$ -, and $V$ -components at TOA (Fig. 11), respectively, and are notably close to zero as $|cos\theta|$ decreases. At the surface (Fig. 12), the relative difference between POLDDA and RT3 is close to zero when $|cos\theta|>0.1$ . Thus, in the case of water clouds, the results of POLDDA and RT3 agree.

4.3 Computational efficiency

Although the single-layer solution of POLDDA is complicated by the consideration of multiple scattering compared to RT3, POLDDA also saves time in the operation of the doubling process. To illustrate the efficiency of POLDDA, we compare the efficiency of POLDDA with that of RT3. Figure 13(a) shows the computational times calculated using the POLDDA and RT3 models with increasing optical depth. A two-layer atmosphere with only molecular scattering is considered, and the optical depth of each layer is in the range from 0.001 to 100. The half-stream was 16, as calculated using POLDDA and RT3. Admittedly, the RT3 model has better computational efficiency when the optical depth is less than 0.0001 because it ignores multiple scattering effects when solving the equation of a single-layer atmosphere. However, the time run by RT3 increases with optical depth. The computational time increases by more than four times when the optical depth reaches 100. Unlike RT3, the computing time of POLDDA does not increase with increasing optical depth of each layer. The operation-time advantage of POLDDA is gradually apparent when the optical depth exceeds 0.0005. Following this trend, POLDDA saves more computational time than RT3 in aerosols or clouded atmospheres (optically thicker than 1). The computational efficiencies of POLDDA and RT3 with different half streams are also considered (Fig. 13b). The layer optical depth was set as 1. As shown in Fig. 13(b), the computational time increases with an increasing number of streams for both POLDDA and RT3, whereas the computational time for RT3 is longer than that for POLDDA with an increasing number of streams.

5 Conclusions and discussion

In this study, the polarized adding method of discrete ordinate approximation (POLDDA) was developed for ultraviolet-visible and near-infrared spectra. The single-layer polarized radiative transfer equation and inhomogeneous multilayer connection are solved using the discrete ordinate method and adding method, respectively. From the accuracy evaluation results in a multiple-layer standard atmosphere with Rayleigh scattering and a cloud layer, POLDDA proved to conform to the results of PolRadtran/RT3. The RMSE values of the Stokes vectors between POLDDA and RT3 against MYSTIC were found to be similar, which further confirms the good accuracy of POLDDA. POLDDA also has a high computational efficiency, particularly for an atmosphere with an optical depth of over 0.0005 compared to RT3. Moreover, the computation time of POLDDA does not increase with the optical depth of each layer.

Appendix A

From Eq.(1), we expand the phase matrix and $\mathbf{L}$ into the Fourier cosine and sine series,


$\displaystyle\left[\begin{array}[]{c}I(\tau,\mu,\varphi)\\ Q(\tau,\mu,\varphi)\\ U(\tau,\mu,\varphi)\\ V(\tau,\mu,\varphi)\\ \end{array}\right]$	$\displaystyle=$	$\displaystyle\sum^{2M-1}_{m=0}\left[\begin{array}[]{c}I^{m}(\tau,\mu)\cdot cos(\varphi-\varphi_{0})\\ Q^{m}(\tau,\mu)\cdot cos(\varphi-\varphi_{0})\\ U^{m}(\tau,\mu)\cdot sin(\varphi-\varphi_{0})\\ V^{m}(\tau,\mu)\cdot sin(\varphi-\varphi_{0})\\ \end{array}\right]$	(A.1i)
$\displaystyle\mathbf{Z}$	$\displaystyle=$	$\displaystyle\left[\begin{array}[]{cccc}z_{11}&z_{12}&z_{13}&z_{14}\\ z_{21}&z_{22}&z_{23}&z_{24}\\ z_{31}&z_{32}&z_{33}&z_{34}\\ z_{41}&z_{42}&z_{43}&z_{44}\end{array}\right]$	(A.1n)

where $z_{xy}(\mu,\varphi,\mu^{\prime},\varphi^{\prime})=\sum^{2M-1}_{m=0}\left\{\begin{array}[]{lr}z^{m}_{xy}(\mu,\varphi,\mu^{\prime},\varphi^{\prime})\cdot cosm(\varphi-\varphi^{\prime})\\ z^{m}_{xy}(\mu,\varphi,\mu^{\prime},\varphi^{\prime})\cdot sinm(\varphi-\varphi^{\prime})\\ \end{array}\right.$ Thus, the radiation transfer equation can be split into 2M equations.

The discrete ordinate method is used to solve the equation. To achieve this, a Gaussian quadrature was used to handle the integration in Eq.(1) as

\int_{-1}^{1}\mathbf{Z}^{m}(\mu,\mu^{\prime})\mathbf{L}^{m}(\tau,\mu^{\prime})d\mu^{\prime}=\sum^{N}_{j=-N,j\neq 0}a_{j}\mathbf{Z}^{m}(\mu,\mu_{j})\mathbf{L}^{m}(\tau,\mu_{j})

(A.2)

where $N$ denotes the number of half-streams, and $\mu_{j}=-\mu_{-j}$ and $a_{j}=a_{-j}$ $(j=1,2,...,N)$ denote the quadrature angles and weights, respectively.

Upon substituting Eqs.(A1)-(A2) into Eq.(1) and neglecting the superscript $m$ , Eq.(1) can be written as follows:

\begin{split}\mu\frac{d\mathbf{L}(\tau,\mu_{i})}{d\tau}=&\mathbf{L}(\tau,\mu_{i})-\frac{\omega}{4}(1+\delta_{0,m})\sum^{N}_{j=-N,j\neq 0}a_{j}\mathbf{Z}(\mu_{i},\mu_{j})\mathbf{L}(\tau,\mu_{j})\\ &-\frac{\omega}{4\pi}\mathbf{Z}(\mu_{i},-\mu_{0})\mathbf{F_{0}}e^{-\tau/\mu_{0}}.\end{split}

(A.3)

where $N$ indicates the number of half-streams, and $\mu_{j}=-\mu_{-j}$ and $a_{j}=a_{-j}$ $(j=1,2,...,N)$ denote the quadrature angles and weights, respectively. They can be written in matrix form as Eq.(2).

Appendix B

In Eq.(7),

\begin{array}[]{lcrc}\textbf{A}^{I}=\left[\begin{array}[]{rc}R^{I}(\tau_{1},\mu_{1}),\cdots,R^{I}(\tau_{1},\mu_{N})\end{array}\right]_{N\times 1}^{\mathbb{T}},\textbf{A}^{Q}=\left[\begin{array}[]{rc}R^{Q}(\tau_{1},\mu_{1}),\cdots,R^{Q}(\tau_{1},\mu_{N})\end{array}\right]_{N\times 1}^{\mathbb{T}},\\ \textbf{A}^{U}=\left[\begin{array}[]{rc}R^{U}(\tau_{1},\mu_{1}),\cdots,R^{U}(\tau_{1},\mu_{N})\end{array}\right]_{N\times 1}^{\mathbb{T}},\textbf{A}^{V}=\left[\begin{array}[]{rc}R^{V}(\tau_{1},\mu_{1}),\cdots,R^{V}(\tau_{1},\mu_{N})\end{array}\right]_{N\times 1}^{\mathbb{T}},\\ \textbf{D}^{I}=\left[\begin{array}[]{rc}T^{I}(\tau_{1},-\mu_{1}),\cdots,T^{I}(\tau_{1},-\mu_{N})\end{array}\right]_{N\times 1}^{\mathbb{T}},\textbf{D}^{Q}=\left[\begin{array}[]{rc}T^{Q}(\tau_{1},-\mu_{1}),\cdots,T^{Q}(\tau_{1},-\mu_{N})\end{array}\right]_{N\times 1}^{\mathbb{T}},\\ \textbf{D}^{U}=\left[\begin{array}[]{rc}T^{U}(\tau_{1},-\mu_{1}),\cdots,T^{U}(\tau_{1},-\mu_{N})\end{array}\right]_{N\times 1}^{\mathbb{T}},\textbf{D}^{V}=\left[\begin{array}[]{rc}T^{V}(\tau_{1},-\mu_{1}),\cdots,T^{V}(\tau_{1},-\mu_{N})\end{array}\right]_{N\times 1}^{\mathbb{T}}\end{array}

and others matrix are $\mathbf{R}_{1}=[\mathbf{R}_{1}^{I},\mathbf{R}_{1}^{Q},\mathbf{R}_{1}^{U},\mathbf{R}_{1}^{V}]^{\mathbb{T}}$ , $\mathbf{T}_{1}=[\mathbf{T}_{1}^{I},\mathbf{T}_{1}^{Q},\mathbf{T}_{1}^{U},\mathbf{T}_{1}^{V}]^{\mathbb{T}}$ ,
$\mathbf{R}_{2}=[\mathbf{R}_{2}^{I},\mathbf{R}_{2}^{Q},\mathbf{R}_{2}^{U},\mathbf{R}_{2}^{V}]^{\mathbb{T}}$ , $\mathbf{T}_{2}=[\mathbf{T}_{2}^{I},\mathbf{T}_{2}^{Q},\mathbf{T}_{2}^{U},\mathbf{T}_{2}^{V}]^{\mathbb{T}}$ ,

\overline{\mathbf{R}}_{2}=\left[\begin{array}[]{cccc}\overline{\cal R}_{2}^{I\leftarrow I}&\overline{\cal R}_{2}^{I\leftarrow Q}&\overline{\cal R}_{2}^{I\leftarrow U}&\overline{\cal R}_{2}^{I\leftarrow V}\\ \overline{\cal R}_{2}^{Q\leftarrow I}&\overline{\cal R}_{2}^{Q\leftarrow Q}&\overline{\cal R}_{2}^{Q\leftarrow U}&\overline{\cal R}_{2}^{Q\leftarrow V}\\ \overline{\cal R}_{2}^{U\leftarrow I}&\overline{\cal R}_{2}^{U\leftarrow Q}&\overline{\cal R}_{2}^{U\leftarrow U}&\overline{\cal R}_{2}^{U\leftarrow V}\\ \overline{\cal R}_{2}^{V\leftarrow I}&\overline{\cal R}_{2}^{V\leftarrow Q}&\overline{\cal R}_{2}^{V\leftarrow U}&\overline{\cal R}_{2}^{V\leftarrow V}\\ \end{array}\right]

\overline{\mathbf{R}}^{*}_{1}=\left[\begin{array}[]{cccc}\overline{\cal R}_{1}^{*,I\leftarrow I}&\overline{\cal R}_{1}^{*,I\leftarrow Q}&\overline{\cal R}_{1}^{*,I\leftarrow U}&\overline{\cal R}_{1}^{*,I\leftarrow V}\\ \overline{\cal R}_{1}^{*,Q\leftarrow I}&\overline{\cal R}_{1}^{*,Q\leftarrow Q}&\overline{\cal R}_{1}^{*,Q\leftarrow U}&\overline{\cal R}_{1}^{*,Q\leftarrow V}\\ \overline{\cal R}_{1}^{*,U\leftarrow I}&\overline{\cal R}_{1}^{*,U\leftarrow Q}&\overline{\cal R}_{1}^{*,U\leftarrow U}&\overline{\cal R}_{1}^{*,U\leftarrow V}\\ \overline{\cal R}_{1}^{*,V\leftarrow I}&\overline{\cal R}_{1}^{*,V\leftarrow Q}&\overline{\cal R}_{1}^{*,V\leftarrow U}&\overline{\cal R}_{1}^{*,V\leftarrow V}\\ \end{array}\right]

\overline{\mathbf{T}}^{*}_{1}=\left[\begin{array}[]{cccc}\overline{\cal T}_{1}^{*,I\leftarrow I}&\overline{\cal T}_{1}^{*,I\leftarrow Q}&\overline{\cal T}_{1}^{*,I\leftarrow U}&\overline{\cal T}_{1}^{*,I\leftarrow V}\\ \overline{\cal T}_{1}^{*,Q\leftarrow I}&\overline{\cal T}_{1}^{*,Q\leftarrow Q}&\overline{\cal T}_{1}^{*,Q\leftarrow U}&\overline{\cal T}_{1}^{*,Q\leftarrow V}\\ \overline{\cal T}_{1}^{*,U\leftarrow I}&\overline{\cal T}_{1}^{*,U\leftarrow Q}&\overline{\cal T}_{1}^{*,U\leftarrow U}&\overline{\cal T}_{1}^{*,U\leftarrow V}\\ \overline{\cal T}_{1}^{*,V\leftarrow I}&\overline{\cal T}_{1}^{*,V\leftarrow Q}&\overline{\cal T}_{1}^{*,V\leftarrow U}&\overline{\cal T}_{1}^{*,V\leftarrow V}\\ \end{array}\right]

\overline{\mathbf{T}}_{2}=\left[\begin{array}[]{cccc}\overline{\cal T}_{2}^{I\leftarrow I}&\overline{\cal T}_{2}^{I\leftarrow Q}&\overline{\cal T}_{2}^{I\leftarrow U}&\overline{\cal T}_{2}^{I\leftarrow V}\\ \overline{\cal T}_{2}^{Q\leftarrow I}&\overline{\cal T}_{2}^{Q\leftarrow Q}&\overline{\cal T}_{2}^{Q\leftarrow U}&\overline{\cal T}_{2}^{Q\leftarrow V}\\ \overline{\cal T}_{2}^{U\leftarrow I}&\overline{\cal T}_{2}^{U\leftarrow Q}&\overline{\cal T}_{2}^{U\leftarrow U}&\overline{\cal T}_{2}^{U\leftarrow V}\\ \overline{\cal T}_{2}^{V\leftarrow I}&\overline{\cal T}_{2}^{V\leftarrow Q}&\overline{\cal T}_{2}^{V\leftarrow U}&\overline{\cal T}_{2}^{V\leftarrow V}\\ \end{array}\right]

In Eq.(9),

\begin{array}[]{lcrc}\overline{\mathbf{X}}_{i}^{j}=\{\mathbf{E}-\left[\begin{array}[]{cccc}\overline{\cal R}_{i}^{I\leftarrow I}&\overline{\cal R}_{i}^{I\leftarrow Q}&\overline{\cal R}_{i}^{I\leftarrow U}&\overline{\cal R}_{i}^{I\leftarrow V}\\ \overline{\cal R}_{i}^{Q\leftarrow I}&\overline{\cal R}_{i}^{Q\leftarrow Q}&\overline{\cal R}_{i}^{Q\leftarrow U}&\overline{\cal R}_{i}^{Q\leftarrow V}\\ \overline{\cal R}_{i}^{U\leftarrow I}&\overline{\cal R}_{i}^{U\leftarrow Q}&\overline{\cal R}_{i}^{U\leftarrow U}&\overline{\cal R}_{i}^{U\leftarrow V}\\ \overline{\cal R}_{i}^{V\leftarrow I}&\overline{\cal R}_{i}^{V\leftarrow Q}&\overline{\cal R}_{i}^{V\leftarrow U}&\overline{\cal R}_{i}^{V\leftarrow V}\\ \end{array}\right]\left[\begin{array}[]{cccc}\overline{\cal R}_{j}^{*,I\leftarrow I}&\overline{\cal R}_{j}^{*,I\leftarrow Q}&\overline{\cal R}_{j}^{*,I\leftarrow U}&\overline{\cal R}_{j}^{*,I\leftarrow V}\\ \overline{\cal R}_{j}^{*,Q\leftarrow I}&\overline{\cal R}_{j}^{*,Q\leftarrow Q}&\overline{\cal R}_{j}^{*,Q\leftarrow U}&\overline{\cal R}_{j}^{*,Q\leftarrow V}\\ \overline{\cal R}_{j}^{*,U\leftarrow I}&\overline{\cal R}_{j}^{*,U\leftarrow Q}&\overline{\cal R}_{j}^{*,U\leftarrow U}&\overline{\cal R}_{j}^{*,U\leftarrow V}\\ \overline{\cal R}_{j}^{*,V\leftarrow I}&\overline{\cal R}_{j}^{*,V\leftarrow Q}&\overline{\cal R}_{j}^{*,V\leftarrow U}&\overline{\cal R}_{j}^{*,V\leftarrow V}\\ \end{array}\right]\}\\ \mathbf{Y}_{i}^{j}=\{\left[\begin{array}[]{cccc}\mathbf{R}_{i}^{I}\\ \mathbf{R}_{i}^{Q}\\ \mathbf{R}_{i}^{U}\\ \mathbf{R}_{i}^{V}\\ \end{array}\right]e^{-\frac{\tau_{j}}{\mu_{0}}}+\left[\begin{array}[]{cccc}\overline{\cal R}_{i}^{I\leftarrow I}&\overline{\cal R}_{i}^{I\leftarrow Q}&\overline{\cal R}_{i}^{I\leftarrow U}&\overline{\cal R}_{i}^{I\leftarrow V}\\ \overline{\cal R}_{i}^{Q\leftarrow I}&\overline{\cal R}_{i}^{Q\leftarrow Q}&\overline{\cal R}_{i}^{Q\leftarrow U}&\overline{\cal R}_{i}^{Q\leftarrow V}\\ \overline{\cal R}_{i}^{U\leftarrow I}&\overline{\cal R}_{i}^{U\leftarrow Q}&\overline{\cal R}_{i}^{U\leftarrow U}&\overline{\cal R}_{i}^{U\leftarrow V}\\ \overline{\cal R}_{i}^{V\leftarrow I}&\overline{\cal R}_{i}^{V\leftarrow Q}&\overline{\cal R}_{i}^{V\leftarrow U}&\overline{\cal R}_{i}^{V\leftarrow V}\\ \end{array}\right]\left[\begin{array}[]{cccc}\mathbf{T}_{j}^{I}\\ \mathbf{T}_{j}^{Q}\\ \mathbf{T}_{j}^{U}\\ \mathbf{T}_{j}^{V}\\ \end{array}\right]\}\end{array}

and $\mathbf{E}$ is a 4N $\times$ 4N identity matrix.

In Eq.(12),

\begin{array}[]{lcrc}\overline{\mathbf{X}}_{i}^{*j}=\{\mathbf{E}-\left[\begin{array}[]{cccc}\overline{\cal R}_{i}^{*,I\leftarrow I}&\overline{\cal R}_{i}^{*,I\leftarrow Q}&\overline{\cal R}_{i}^{*,I\leftarrow U}&\overline{\cal R}_{i}^{*,I\leftarrow V}\\ \overline{\cal R}_{i}^{*,Q\leftarrow I}&\overline{\cal R}_{i}^{*,Q\leftarrow Q}&\overline{\cal R}_{i}^{*,Q\leftarrow U}&\overline{\cal R}_{i}^{*,Q\leftarrow V}\\ \overline{\cal R}_{i}^{*,U\leftarrow I}&\overline{\cal R}_{i}^{*,U\leftarrow Q}&\overline{\cal R}_{i}^{*,U\leftarrow U}&\overline{\cal R}_{i}^{*,U\leftarrow V}\\ \overline{\cal R}_{i}^{*,V\leftarrow I}&\overline{\cal R}_{i}^{*,V\leftarrow Q}&\overline{\cal R}_{i}^{*,V\leftarrow U}&\overline{\cal R}_{i}^{*,V\leftarrow V}\\ \end{array}\right]\left[\begin{array}[]{cccc}\overline{\cal R}_{j}^{I\leftarrow I}&\overline{\cal R}_{j}^{I\leftarrow Q}&\overline{\cal R}_{j}^{I\leftarrow U}&\overline{\cal R}_{j}^{I\leftarrow V}\\ \overline{\cal R}_{j}^{Q\leftarrow I}&\overline{\cal R}_{j}^{Q\leftarrow Q}&\overline{\cal R}_{j}^{Q\leftarrow U}&\overline{\cal R}_{j}^{Q\leftarrow V}\\ \overline{\cal R}_{j}^{U\leftarrow I}&\overline{\cal R}_{j}^{U\leftarrow Q}&\overline{\cal R}_{j}^{U\leftarrow U}&\overline{\cal R}_{j}^{U\leftarrow V}\\ \overline{\cal R}_{j}^{V\leftarrow I}&\overline{\cal R}_{j}^{V\leftarrow Q}&\overline{\cal R}_{j}^{V\leftarrow U}&\overline{\cal R}_{j}^{V\leftarrow V}\\ \end{array}\right]\}\end{array}

Acknowledgments

The code of POLDDA is available from the corresponding author upon reasonable request. This study was supported by the National Natural Science Foundation of China (42222506 and 42105081) and China Postdoctoral Science Foundation (2023M730618).

References

[1] M. Duan, Q. Min, D. Lu, A polarized radiative transfer model based on successive order of scattering, Advances in Atmospheric Sciences 27 (2010) 891–900.
[2] Q. Yin, C. Song, Fundamental definition of two-stream approximation for radiative transfer in scattering atmosphere, IEEE Transactions on Geoscience and Remote Sensing 60 (2022) 2003614.
[3] J. E. Hansen, L. D. Travis, Light scattering in planetary atmosphere, Space Science Reviews 16 (1974) 527–610.
[4] M. I. Mishchenko, A. A. Lacis, L. D. Travis, Errors induced by the neglect of polarization in radiance calculations for rayleigh-scattering atmospheres, Journal of Quantitative Spectroscopy and Radiative Transfer 51 (1994) 491–510.
[5] A. A. Lacis, J. Chowdhary, M. I. Mishchenko, B. Cairns, Modeling errors in diffuse-sky radiation: vector versus scalar treatment, Journal of the Atmospheric Sciences 25 (1998) 135–138.
[6] L. Oikarinen, Polarization of light in uv-visible limb radiance measurements, Journal of Geophysical Research 106 (2001) 1533–1544.
[7] R. C. Levy, L. A. Remer, Y. J. Kaufman, Effects of neglecting polarization on the modis aerosol retrieval over land, IEEE Transactions on Geoscience and Remote Sensing 42 (2004) 2576–2583.
[8] D. M. Stam, J. W. Hovenier, Errors in calculated planetary phase functions and albedos due to neglecting polarization, Astronomy and Astrophysics 444 (2005) 275–286.
[9] V. Natraj, R. Spurr, H. Boesch, Y. Jiang, Y. Yung, Evaluation of errors in neglecting polarization in the forward modeling of o2 a band measurements from space, with relevance to co2 column retrieval from polarization sensitive instruments, Journal of Quantitative Spectroscopy and Radiative Transfer 103 (2007) 245–259.
[10] S. Chandrasekhar, Radiative Transfer, Oxford University Press, 1950.
[11] Y. Gao, M. Z. Duan, X. Y. Huang, Preliminary comparisons of the typical polarized radiative transfer models: precision and efficiency, Journal of Remote Sensing 14 (2010) 839–851.
[12] Y. Cai, F. Zhang, H. Lin, J. Li, H. Zhang, W. Li, S. Hu, Optimized alternate mapping correlated k-distribution method for atmospheric longwave radiative transfer, Journal of Advances in Modeling Earth Systems 15 (2023) e2022MS003419.
[13] M. I. Mishchenko, Coauthors, Accurate monitoring of terrestrial aerosols and total solar irradiance: introducing the glory mission, Bulletin of the American Meteorological Society 88 (2007) 677–691.
[14] J. E. Hansen, Multiple scattering of polarized light in planetary atmospheres. part ii. sunlight reflected by terrestrial water clouds, Journal of the Atmospheric Sciences 28 (1971) 1400–1426.
[15] G. Plass, G. Kattawar, F. Catchings, Matrix operator theory of radiative transfer. i: Rayleigh scattering, Applied Optics 12 (1973) 314–29.
[16] J. F. D. Haan, P. B. Bosma, J. W. Hovenier, The adding method for multiple scattering computations of polarized light, Astronomy and Astrophysics 183 (1987) 371–391.
[17] K. F. Evans, G. L. Stephens, A new polarized atmospheric radiative transfer model, Journal of Quantitative Spectroscopy and Radiative Transfer 46 (1991) 413–423.
[18] W. G. Bai, P. Zhang, W. J. Zhang, G. Ma, C.Qi, A model for accurately calculating hyper-spectral, middle-shortwave infrared radiative transfer for remote sensing, Science China Earth Sciences 47 (2017) 1483–1492.
[19] W. G. Bai, P. Zhang, W. J. Zhang, J. Li, G. Ma, C.Qi, H. Liu, Jacobian matrix for near-infrared remote sensing based on vector radiative transfer model, Science China Earth Sciences 63 (2020) 1353–1365.
[20] K. Stamnes, P. Conklin, A new multi-layer discrete ordinate approach to radiative transfer in vertically inhomogeneous atmospheres, Journal of Quantitative Spectroscopy and Radiative Transfer 31 (1984) 273–282.
[21] F. Weng, A multi-layer discrete-ordinate method for vector radiative transfer in a vertically-inhomogeneous, emitting and scattering atmosphere-i. theory, Journal of Quantitative Spectroscopy and Radiative Transfer 47 (1992) 19–33.
[22] F. M. Schulz, K. Stamnes, F. Weng, Vdisort: An improved and generalized discrete ordinate radiative transfer model for polarized (vector) radiative transfer, Journal of Quantitative Spectroscopy and Radiative Transfer 61 (1999) 105–122.
[23] C. E. Siewert, A discrete-ordinates solution for radiative-transfer models that include polarization effects, Journal of Quantitative Spectroscopy and Radiative Transfer 64 (2000) 227–254.
[24] Y. Ota, A. Higurashi, T. Nakajima, T. Yokota, Matrix formulations of radiative transfer including the polarization effect in a coupled atmosphere-ocean system, Journal of Quantitative Spectroscopy and Radiative Transfer 111 (2010) 878–94.
[25] G. W. Kattawar, G. N. Plass, Radiance and polarization of multiple scattered light from haze and clouds, Appl. Opt. 7 (1968) 1519–1527.
[26] D. Collins, W. Blattner, M. Wells, H. Horak, Backward monte-carlo calculations of polarization characteristics of the radiation emerging from a spherical shell atmosphere, Applied Optics 11 (1972) 2684–2705.
[27] B. Y. Wu, D. R. Lu, Simulation of twilights after el chic eruption with the monte-carlo method, Journal of the Atmospheric Sciences 13 (1989) 204–213.
[28] R. L. Roberti, C. Kummerow, Monte carlo calculations of polarized microwave radiation emerging from cloud structures, J. Geophy. Res. 104 (1999) 2093–2104.
[29] W. Irvine, Multiple scattering in planetary atmospheres, Icarus 25 (1975) 175–204.
[30] Q. L. Min, M. Z. Duan, A successive order of scattering model for solving vector radiative transfer in the atmosphere, Journal of Quantitative Spectroscopy and Radiative Transfer 87 (2004) 243–259.
[31] M. I. Mishchenko, The fast invariant imbedding method for polarized light: computational aspects and numerical results for rayleigh scattering, Journal of Quantitative Spectroscopy and Radiative Transfer 43 (1990) 163–171.
[32] K. F. Evans, The spherical harmonic discrete ordinate method for three-dimensional atmospheric radiative transfer, Journal of the Atmospheric Sciences 55 (1998) 429–446.
[33] T. Nakajima, M. Tanaka, Matrix formulations for the transfer of solar radiation in a plane-parallel scattering atmosphere, Journal of Quantitative Spectroscopy and Radiative Transfer 35 (1986) 13–21.
[34] F. Zhang, Z. Shen, J. Li, X. Zhou, L. Ma, Analytical delta-four-stream doubling-adding method for radiative transfer parameterizations, Journal of the Atmospheric Sciences 70 (2013) 794–808.
[35] Q. Liu, F. Weng, A microwave polarimetric two-stream radiative transfer model, Journal of the Atmospheric Sciences 59 (2002) 2396–2402.
[36] K. N. Liou, S. C. Ou, Y. Takano, Q. Liu, A polarized delta-four-stream approximation for infrared and microwave radiative, Journal of the Atmospheric Sciences 62 (2005) 2542–2554.
[37] W. Li, F. Zhang, F. Bao, K. Wu, J. Li, P. Zhang, W. Han, Polarized discrete ordinate adding approximation for infrared and microwave radiative transfer, Journal of Quantitative Spectroscopy and Radiative Transfer 293 (2022) 108368.
[38] R. J. Spurr, Vlidort: A linearized pseudo-spherical vector discrete ordinate radiative transfer code for forward model and retrieval studies in multilayer multiple scattering media, Journal of Quantitative Spectroscopy and Radiative Transfer 102 (2006) 316–342.
[39] C. Emde, V. Barlakas, C. Cornet, E. F., S. Korkin, Y. Ota, L. C. Labonnote, A. Lyapustin, A. Macke, B. Mayer, M. Wendisch, Iprt polarized radiative transfer model intercomparison project - phase a, Journal of Quantitative Spectroscopy and Radiative Transfer 164 (2015) 8–36.
[40] R. A. McClatchey, R. W. Fenn, J. E. A. Selby, F. E. Volz, J. S.Garing, Optical properties of the atmosphere, Air Force Rep. AFCRL-71-0279.

Table 1: Relative root mean square errors in percent between MYSTIC and POLDDA, RT3 for all cases in this study.

Case		POLDDA	RT3
Case1-1	I	0.01773589%	0.01773660%
	Q	0.02265512%	0.02264818%
Case1-2	I	0.01649993%	0.01650125%
	Q	0.02561625%	0.02561332%
	U	0.02343666%	0.02343717%
Case1-3	I	0.01660634%	0.01660778%
	Q	0.02499878%	0.02499902%
	U	0.02290941%	0.02290969%
Case2	I	0.008561691%	0.008559853%
	Q	0.03571449%	0.03571167%
	U	0.02874306%	0.02874456%
Case3	I	0.01723639%	0.01723684%
	Q	0.02661624%	0.02660797%
	U	0.01979122%	0.01979095%
Case4	I	0.01260397%	0.01260261%
	Q	0.02472782%	0.02472788%
	U	0.02089377%	0.02090551%

Refer to caption — Figure 1: Schematic diagram of the principles of invariance in vector solar radiative transfer.