Higher order statistics of the Mueller matrix in a fiber with an arbitrary length impacted by PMD

Junhe Zhou; Junhe Zhou; Qinsong Hu

doi:10.1364/OE.404223

1. Introduction

Polarization mode dispersion (PMD) in optical fibers has been thoroughly studied from the 1990s to the 2010s [1–11]. The mature theory provides a general picture of the statistical nature of the signal propagation inside the fiber with PMD [1–11]. The signal polarization state keeps changing due to the random birefringence and this gives rise to the averaging effect over the Poincare sphere in a sufficiently long fiber [2,4]. Many theoretical results are based on the above conclusion, and they are essential for the system optimization [12–13] and the sensing data analysis [14]. One profound example of the theoretical derivation from the uniform polarization state distribution is the well-known Manakov equation, which governs the propagation of two polarized signals in a long transmission fiber [15].

While the statistics in fibers with PMD has been paid special attention to, there are still some unreached territories. First of all, the existing statistical dynamic study usually focuses on the statistical dynamic behavior of the PMD vector rather than the statistical dynamic information of transfer matrix (such as the Jones matrix and the Mueller matrix). The matrix elements are usually assumed to be purely random with the related Stokes vector distributed uniformly on the Poincare sphere [2,4,16]. While the above assertion is true in a fiber which is sufficiently longer than the PMD coherent length, it is not the case in short fibers. In the fiber systems designed for sensing applications, the implementation of the short fibers is not unusual [17]. However, no special treatment has been proposed for the short fiber transfer matrix statistics analysis, particularly the statistical behavior of the matrix elements before reaching their final “stable” uniformly distributed random state [2,4]. Secondly, the current discussions on the matrix elements statistics are usually limited to its 2nd order [2]. Although a particular type of the 4th order Mueller matrix elements moments is studied in [16], it does not involve the general 4th order moments including the cross-correlation terms between the Mueller matrix elements, which are essential for the back-scattered signal transfer matrix characterization. While the gain/power of the back-scattered signal can be characterized by the 2nd order moments of the forward propagation Mueller matrix elements [2], their variations must be analyzed based on the 4th order moments. Such gain/power variation can be important for practical applications like the backwardly pumped Raman amplifiers [18] and the Brillouin gain fluctuations in the distributed sensors [19].

Henceforth, despite the outstanding pioneering researches [1–19] in this field, there is still a great demand for a mathematical tool for the treatment of the transfer matrix statistics, particularly the dynamic statistical behavior in short fibers and the higher order statistics of the transfer matrix elements.

In this paper, we proposed a novel coupled stochastic differential equation (SDE) approach to analyze the Mueller matrix elements statistics. Based on the fundamental SDEs, one is not only able to derive their related ordinary differential equations (ODEs), but also to derive the SDEs for the higher order moments. By this systematic approach, it is possible to obtain the ODEs of the moments with arbitrary orders along with their analytical solutions. Since the derivation is not based on the long fiber assumption, the dynamic behavior of the Mueller matrix elements is included in the results, which can be applied to a fiber with an arbitrary length. The long length limit of the analytical solutions to the 2nd order moments agree well with the published results. The dynamic behavior of the fiber statistics and the 4th order moments (including the cross-correlation terms) are calculated by the analytical formulas and the Monte Carlo simulations, and an excellent agreement has been demonstrated.

2. Theory

The Jones matrix [4] and the Mueller matrix [2,19–20] are two approaches to represent the transfer matrix in a random birefringent fiber. Here, we use the Mueller matrix as the main mathematical tool for the statistical analysis.

2.1 Fundamental propagation model in the presence of PMD

Since the fiber is assumed to be lossless and the overall power does not vary during the propagation, the reduced Stokes vector (with three components) is used. The Stokes vector evolution can be characterized by the following equation [1–2]

(1)$$\begin{array}{l} \frac{{d{\textbf s}}}{{dz}} = {{\mathrm{\boldsymbol {\beta}}} } \times {\textbf s}\\ {{\mathrm{\boldsymbol {\beta}}} } = {({{\beta_1},{\beta_2},{\beta_3}} )^T}, \end{array}$$

where T stands for the matrix transpose, s the Stokes vector, β the local birefringence vector.

Here we may denote the operation β× by a matrix [1]

(2)$${{\mathrm{\boldsymbol {\beta}}} } \times{=} {{\textbf T}_\beta } = \left( {\begin{array}{ccc} 0&{ - {\beta_3}}&{{\beta_2}}\\ {{\beta_3}}&0&{ - {\beta_1}}\\ { - {\beta_2}}&{{\beta_1}}&0 \end{array}} \right),$$

The matrix T_β is sometimes referred to as the differential Mueller matrix [21–22]. The input and the output Stokes vectors are related to each other by [1–2]

(3)$${\textbf s}(z )= {\textbf R}(z ){\textbf s}(0 ),$$

It is straight forward to have

(4)$$\frac{{d{\textbf R}}}{{dz}} = {{\textbf T}_\beta }{\textbf R},$$

Matrix R is also referred to as the forward propagation Mueller matrix [2], and its relationship with the differential Mueller matrix is also mentioned in [22]. The round-tip propagation matrix for the back-scattered signal is T = MR^TMR with M as the diagonal matrix diag (1,1,-1) [2]. Clearly, the average value of the T elements corresponds the 2nd order moments of the matrix R elements, while the variance (the 2nd order moments) of the T elements corresponds to the 4th order moments of the R elements.

2.2 Derivation of the fundamental SDEs

The local birefringence vector β can be modeled as the three dimensional Brownian motion [5,12,13]

(5)$$\begin{array}{l} \left\langle {{\mathrm{\boldsymbol {\beta}}} } \right\rangle = 0\\ \left\langle {{{\mathrm{\boldsymbol {\beta}}} }({{z_1}} ){{{\mathrm{\boldsymbol {\beta}}} }^T}({{z_2}} )} \right\rangle = {\sigma ^2}{\textbf I}\delta ({{z_1} - {z_2}} )\\ \left\langle {{{\mathrm{\boldsymbol {\beta}}} }(z )dz{{{\mathrm{\boldsymbol {\beta}}} }^T}(z )dz} \right\rangle = {\sigma ^2}{\textbf I}dz, \end{array}$$

where σ² stands for the variance of the birefringence vector element, I the 3×3 identity matrix. It is worth mentioning that σ² should be frequency dependent [5,18] and equals ω²D_p²/3, where ω is the angular frequency of the signal, D_p the differential group delay [18]. Since we are interested in the statistical property of the Mueller matrix rather than the PMD vector, the frequency dependent terms are not separated in Eq. (5) and the unit for σ² should be (1/m). Hence, we have

(6)$$\left\langle {{{\textbf T}_\beta }dz{{\textbf T}_\beta }dz} \right\rangle ={-} 2{\sigma ^2}{\textbf I}dz$$

Therefore, one may derive the fundamental SDE from Eq. (4) by converting it from the stratonovich sense to the Ito sense (detailed in the appendix)

(7)$$d{\textbf R} = {{\textbf T}_\beta }dz{\textbf R} - {\sigma ^2}dz{\textbf R},$$

2.3 Derivation of the ODEs for the 2nd order and the 4th order statistics

The ODE for the 2nd order moments can be derived from Eq. (7) based on the important matrix property, i.e.$({{\textbf A} \otimes {\textbf B}} )({{\textbf C} \otimes {\textbf D}} )= ({{\textbf {AC}}} )\otimes ({{\textbf {BD}}} )$.

(8)$$\begin{array}{l} d({{\textbf R} \otimes {\textbf R}} )= d{\textbf R} \otimes {\textbf R} + {\textbf R} \otimes d{\textbf R} + d{\textbf R} \otimes d{\textbf R}\\ = ({{{\textbf T}_\beta }dz{\textbf R} - {\sigma^2}{\textbf R}dz} )\otimes {\textbf R} + {\textbf R} \otimes ({{{\textbf T}_\beta }dz{\textbf R} - {\sigma^2}{\textbf R}dz} )+ {{\textbf T}_\beta }dz{\textbf R} \otimes {{\textbf T}_\beta }dz{\textbf R}\\ = ({{{\textbf T}_\beta }dz \otimes {{\textbf I}_3}} )({{\textbf R} \otimes {\textbf R}} )+ ({{{\textbf I}_3} \otimes {{\textbf T}_\beta }dz} )({{\textbf R} \otimes {\textbf R}} )- 2{\sigma ^2}({{\textbf R} \otimes {\textbf R}} )dz + {\textbf K}({{\textbf R} \otimes {\textbf R}} )dz, \end{array}$$

where I_n stands for the identity matrix with the size of n×n, and ${\otimes}$ denotes Kronecker matrix product, with

\begin{array}{l} {\textbf K} = \frac{{\left\langle {{{\textbf T}_\beta }dz \otimes {{\textbf T}_\beta }dz} \right\rangle }}{{dz}}\\ = \left( {\begin{array}{ccccccccc} 0&0&0&0&{{\sigma^2}}&0&0&0&{{\sigma^2}}\\ 0&0&0&{ - {\sigma^2}}&0&0&0&0&0\\ 0&0&0&0&0&0&{ - {\sigma^2}}&0&0\\ 0&{ - {\sigma^2}}&0&0&0&0&0&0&0\\ {{\sigma^2}}&0&0&0&0&0&0&0&{{\sigma^2}}\\ 0&0&0&0&0&0&0&{ - {\sigma^2}}&0\\ 0&0&{ - {\sigma^2}}&0&0&0&0&0&0\\ 0&0&0&0&0&{ - {\sigma^2}}&0&0&0\\ {{\sigma^2}}&0&0&0&{{\sigma^2}}&0&0&0&0 \end{array}} \right), \end{array}

Taking the average on Eq. (8), one has

(9)$$\begin{array}{l} \frac{{d\left\langle {{\textbf R} \otimes {\textbf R}} \right\rangle }}{{dz}} = {{\textbf K}_0}\left\langle {{\textbf R} \otimes {\textbf R}} \right\rangle \\ {{\textbf K}_0} = {\textbf K} - 2{\sigma ^2}{{\textbf I}_9}, \end{array}$$

The equation for the 4th order moments can be derived based on the equation for the 2nd order moments

(10)$$\begin{array}{l} d({({{\textbf R} \otimes {\textbf R}} )\otimes ({{\textbf R} \otimes {\textbf R}} )} )\\ = d({{\textbf R} \otimes {\textbf R}} )\otimes ({{\textbf R} \otimes {\textbf R}} )+ ({{\textbf R} \otimes {\textbf R}} )\otimes d({{\textbf R} \otimes {\textbf R}} )+ d({{\textbf R} \otimes {\textbf R}} )\otimes d({{\textbf R} \otimes {\textbf R}} )\\ = ({({({{{\textbf T}_\beta }dz \otimes {{\textbf I}_3}} )+ ({{{\textbf I}_3} \otimes {{\textbf T}_\beta }dz} )} )({{\textbf R} \otimes {\textbf R}} )} )\otimes ({{\textbf R} \otimes {\textbf R}} )\\ + ({{\textbf R} \otimes {\textbf R}} )\otimes ({({({{{\textbf T}_\beta }dz \otimes {{\textbf I}_3}} )+ ({{{\textbf I}_3} \otimes {{\textbf T}_\beta }dz} )} )({{\textbf R} \otimes {\textbf R}} )} )\\ + ({({{{\textbf K}_0} \otimes {{\textbf I}_9}} )+ ({{{\textbf I}_9} \otimes {{\textbf K}_0}} )} )({{\textbf R} \otimes {\textbf R}} )\otimes ({{\textbf R} \otimes {\textbf R}} )dz\\ + \left\langle {({({{{\textbf T}_\beta }dz \otimes {{\textbf I}_3}} )+ ({{{\textbf I}_3} \otimes {{\textbf T}_\beta }dz} )} )\otimes ({({{{\textbf T}_\beta }dz \otimes {{\textbf I}_3}} )+ ({{{\textbf I}_3} \otimes {{\textbf T}_\beta }dz} )} )} \right\rangle ({{\textbf R} \otimes {\textbf R}} )\otimes ({{\textbf R} \otimes {\textbf R}} ), \end{array}$$

Taking average on both sides of Eq. (10), we have

(11)$$\begin{array}{l} d\left\langle {({{\textbf R} \otimes {\textbf R}} )\otimes ({{\textbf R} \otimes {\textbf R}} )} \right\rangle \\ = {{\textbf L}_0}\left\langle {({{\textbf R} \otimes {\textbf R}} )\otimes ({{\textbf R} \otimes {\textbf R}} )} \right\rangle dz\\ {{\textbf L}_0} = ({({{{\textbf K}_0} \otimes {{\textbf I}_9}} )+ ({{{\textbf I}_9} \otimes {{\textbf K}_0}} )} )\\ + \frac{{\left\langle {({({{{\textbf T}_\beta }dz \otimes {{\textbf I}_3}} )+ ({{{\textbf I}_3} \otimes {{\textbf T}_\beta }dz} )} )\otimes ({({{{\textbf T}_\beta }dz \otimes {{\textbf I}_3}} )+ ({{{\textbf I}_3} \otimes {{\textbf T}_\beta }dz} )} )} \right\rangle }}{{dz}}, \end{array}$$

Similar to the procedures to find K₀, one may use Eq. (11) in combination with Eq. (5) to get the expression of L₀ (added in Dataset 1 [23]).

One may further derive the ODE for the 8th order moments as

(12)$$\begin{array}{l} d\left\langle {({({{\textbf R} \otimes {\textbf R}} )\otimes ({{\textbf R} \otimes {\textbf R}} )} )\otimes ({({{\textbf R} \otimes {\textbf R}} )\otimes ({{\textbf R} \otimes {\textbf R}} )} )} \right\rangle \\ = {{\textbf N}_0}\left\langle {({({{\textbf R} \otimes {\textbf R}} )\otimes ({{\textbf R} \otimes {\textbf R}} )} )\otimes ({({{\textbf R} \otimes {\textbf R}} )\otimes ({{\textbf R} \otimes {\textbf R}} )} )} \right\rangle \\ {{\textbf N}_0} = ({{{\textbf L}_0} \otimes {{\textbf I}_{81}}} )+ ({{{\textbf I}_{81}} \otimes {{\textbf L}_0}} )+ \frac{{d({{\textbf {QQ}}} )}}{{dz}}\\ {\textbf Q} = ({({{{\textbf T}_\beta }dz \otimes {{\textbf I}_3}} )+ ({{{\textbf I}_3} \otimes {{\textbf T}_\beta }dz} )} )\otimes {{\textbf I}_9} + {{\textbf I}_9} \otimes ({({{{\textbf T}_\beta }dz \otimes {{\textbf I}_3}} )+ ({{{\textbf I}_3} \otimes {{\textbf T}_\beta }dz} )} ), \end{array}$$

During the evaluation of the 4th order moments for the round-tip transfer matrix, one will need the information for the 8th order moments of the forward propagation matrix elements. (The expression of N₀ is added in Dataset 2 [24]).

Since Eqs. (9), (11), and (12) are the linear ODEs with constant coefficients, they can be solved analytically as

(13)$$\begin{array}{l} \left\langle {{\textbf R}(z )\otimes {\textbf R}(z )} \right\rangle \\ = \exp ({{{\textbf K}_0}z} )\left\langle {{\textbf R}(0 )\otimes {\textbf R}(0 )} \right\rangle \\ \left\langle {({{\textbf R}(z )\otimes {\textbf R}(z )} )\otimes ({{\textbf R}(z )\otimes {\textbf R}(z )} )} \right\rangle \\ = \exp ({{{\textbf L}_0}z} )\left\langle {({{\textbf R}(0 )\otimes {\textbf R}(0 )} )\otimes ({{\textbf R}(0 )\otimes {\textbf R}(0 )} )} \right\rangle \\ \left\langle {({({{\textbf R}(z )\otimes {\textbf R}(z )} )\otimes ({{\textbf R}(z )\otimes {\textbf R}(z )} )} )\otimes ({({{\textbf R}(z )\otimes {\textbf R}(z )} )\otimes ({{\textbf R}(z )\otimes {\textbf R}(z )} )} )} \right\rangle \\ = \exp ({{{\textbf N}_0}z} )\left\langle {({({{\textbf R}(0 )\otimes {\textbf R}(0 )} )\otimes ({{\textbf R}(0 )\otimes {\textbf R}(0 )} )} )\otimes ({({{\textbf R}(0 )\otimes {\textbf R}(0 )} )\otimes ({{\textbf R}(0 )\otimes {\textbf R}(0 )} )} )} \right\rangle , \end{array}$$

where exp denotes the matrix exponential function.

2.4 Asymptotic behavior of the moments: validity check

To verify the correctness of the derived analytical formulas, i.e. Eq. (13), we check the asymptotic behavior of the moments when the length of the fiber approaches infinity. Some of these results (particularly the 2nd order moments) have been derived previously and widely adopted in the analysis.

It is easy to find the eigen values of matrix K₀ as σ²(-3,-3,-3,-3,-3,-1,-1,-1,0)^T. The eigen vector u₉ related to the 9th eigen value 0 can be obtained as

(14)$${{\textbf u}_9} = {\left( {\begin{array}{ccccccccc} {\frac{1}{{\sqrt 3 }}}&0&0&0&{\frac{1}{{\sqrt 3 }}}&0&0&0&{\frac{1}{{\sqrt 3 }}} \end{array}} \right)^T}$$

Obviously, one has

(15)$$\begin{array}{l} \mathop {\lim }\limits_{z \to \infty } \exp ({{{\textbf K}_0}z} )= {\textbf U}\left( {\begin{array}{ccc} {\exp ({{\lambda_1}z} )}&{}&{}\\ {}& \ddots &{}\\ {}&{}&{\exp ({{\lambda_9}z} )} \end{array}} \right){{\textbf U}^T}\\ = {\textbf U}\left( {\begin{array}{ccc} 0&{}&{}\\ {}& \ddots &{}\\ {}&{}&1 \end{array}} \right){{\textbf U}^T} = {{\textbf u}_9}{{\textbf u}_9}^T, \end{array}$$

Therefore, we have the 2nd order moments of matrix R element when z approaches infinity as follows:

(16)$$\begin{array}{l} \mathop {\lim }\limits_{z \to \infty } \left\langle {{R_{ij}}^2(z )} \right\rangle = \frac{1}{3}\\ \mathop {\lim }\limits_{z \to \infty } \left\langle {{R_{ij}}(z ){R_{i^{\prime}j^{\prime}}}(z )} \right\rangle = 0({i \ne i^{\prime}\;or\;j \ne j^{\prime}} ), \end{array}$$

where R_ij stands for the element of R on the ith row and jth column. The first equation agrees with the results in [16]. Based on Eq. (16), one may have the asymptotic behavior of round-tip transfer matrix as

(17)$$\begin{array}{l} \mathop {\lim }\limits_{z \to \infty } \left\langle {{\textbf {MR}}{{(z )}^T}{\textbf {MR}}(z )} \right\rangle \\ = \mathop {\lim }\limits_{z \to \infty } \left\langle {\left( {\begin{array}{ccc} 1&{}&{}\\ {}&1&{}\\ {}&{}&{ - 1} \end{array}} \right)\left( {\begin{array}{ccc} {{R_{11}}}&{{R_{21}}}&{{R_{31}}}\\ {{R_{12}}}&{{R_{22}}}&{{R_{32}}}\\ {{R_{13}}}&{{R_{23}}}&{{R_{33}}} \end{array}} \right)\left( {\begin{array}{ccc} 1&{}&{}\\ {}&1&{}\\ {}&{}&{ - 1} \end{array}} \right)\left( {\begin{array}{ccc} {{R_{11}}}&{{R_{12}}}&{{R_{13}}}\\ {{R_{21}}}&{{R_{22}}}&{{R_{23}}}\\ {{R_{31}}}&{{R_{32}}}&{{R_{33}}} \end{array}} \right)} \right\rangle \\ = \mathop {\lim }\limits_{z \to \infty } \left\langle {\left( {\begin{array}{ccc} {{R_{11}}}&{{R_{21}}}&{{R_{31}}}\\ {{R_{12}}}&{{R_{22}}}&{{R_{32}}}\\ { - {R_{13}}}&{ - {R_{23}}}&{ - {R_{33}}} \end{array}} \right)\left( {\begin{array}{ccc} {{R_{11}}}&{{R_{12}}}&{{R_{13}}}\\ {{R_{21}}}&{{R_{22}}}&{{R_{23}}}\\ { - {R_{31}}}&{ - {R_{32}}}&{ - {R_{33}}} \end{array}} \right)} \right\rangle \\ = \left( {\begin{array}{ccc} {1/3}&{}&{}\\ {}&{1/3}&{}\\ {}&{}&{ - 1/3} \end{array}} \right), \end{array}$$

which is in exact agreement with the results in [2].

Similarly, we may get the asymptotic behavior of the 4th order moments of matrix R when z approaches infinity. The eigen values of L₀ are -10σ² (9 eigen values with this quantity), -6σ² (21 eigen values with this quantity), -3σ² (30 eigen values with this quantity), -σ² (18 eigen values with this quantity), and 0 (3 eigen values with this quantity) respectively. Clearly, all of the eigen vectors vanish when z approaches infinity except for the ones associated with the zero eigen values. Therefore, it is possible to evaluate $\mathop {\lim }\limits_{z \to \infty } \exp ({{{\textbf L}_0}z} )$ to find the related asymptotic behavior of the 4th order moments, which is summarized in Eq. (18)

(18)$$\begin{array}{l} \mathop {\lim }\limits_{z \to \infty } \left\langle {{R_{ij}}^4(z )} \right\rangle = \frac{1}{5}\\ \mathop {\lim }\limits_{z \to \infty } \left\langle {{R_{ij}}^2(z ){R_{i^{\prime}j^{\prime}}}^2(z )} \right\rangle = \frac{2}{{15}}({i^{\prime} \ne j,j^{\prime} \ne j} )\\ \mathop {\lim }\limits_{z \to \infty } \left\langle {{R_{ij}}^2(z ){R_{i^{\prime}j^{\prime}}}^2(z )} \right\rangle = \frac{1}{{15}}({i = i^{\prime},j \ne j^{\prime}ori \ne i^{\prime},j = j^{\prime}} )\\ \mathop {\lim }\limits_{z \to \infty } \left\langle {{R_{ij}}(z ){R_{ij}}(z ){R_{ji}}(z ){R_{jj}}(z )} \right\rangle ={-} \frac{1}{{30}}({i \ne j} )\\ \mathop {\lim }\limits_{z \to \infty } \left\langle {{R_{ij}}(z ){R_{i^{\prime}j^{\prime}}}(z ){R_{{i_1}{j_1}}}(z ){R_{{i_1}{j_1}}}(z )} \right\rangle = 0,\\ ({i,j,i^{\prime},j^{\prime},{i_1},{j_1},{i_1}^{\prime},{j_1}^{\prime}\;\textrm{not}\;\textrm{the}\;\textrm{above}\;\textrm{case}} )\end{array}$$

The first equation in Eq. (19) matches the results in [16], while other forms of the 4th order moments have not been reported before.

Henceforth, the validity of the proposed formulas is verified by checking the results in [2] and [16]. The quantitative asymptotic results which had not appeared in the published literatures will be checked via the Monte Carlo simulations and will be presented in the next section.

3. Results and discussions

The proposed analytical model has been verified by the numerical simulations. The local birefringence vector β has been assumed to be the three dimensional Brownian motion. Monte Carlo simulations are conducted by dividing the fiber into many small sections with each section possessing a constant local birefringence vector, which has its three elements to be normally distributed. 5000 realizations of Monte Carlo simulations are calculated to measure the average values including the 2nd order and the 4th order moments. The numerical simulations are based on the original polarization state propagation equation, i.e. Eq. (1). The analytical solutions are obtained via the derived formulas, i.e. Eq. (13). The agreement between the two methods is prominent.

Figure 1 demonstrates the results by the numerical and the analytical methods for the 2nd order moments of the forward propagation Mueller matrix elements. The x-axis uses the normalized propagation distance (σ²L) as the quantity measure.

Fig. 1. The 2nd order moments of the forward propagation matrix elements: results obtained by the Monte Carlo simulations and the analytical formulas.

Download Full Size | PDF

The exact agreement between the analytical prediction and the numerical validation confirms the accuracy of the derived formulas. The maximum discrepancy between the two curves for < R₁₁²> is 5%, which can be reduced if the number of realizations further increases. It can be seen from the figure that the dynamics of the Mueller matrix elements moments is limited within the range of [0 2] for σ²L. As a realistic example, we may assume that the fiber has the differential group delay D_p=0.1 ps²/$\sqrt {km} $, and the signal wavelength is 1550nm. It can be calculated that σ²=4.93/m, while σ²L = 2 suggests L is about 0.2 m. If the differential group delay D_p=0.01 ps²/$\sqrt {km} $, L is about 20 m. The results are in agreement with the previous reports [2,4] that the PMD induced dynamics of the transfer matrix elements occurs only within a short fiber length at the beginning of the propagation. After this, the moments reach a relatively “stable” random state, whose probability density function (pdf) is not related to the fiber position z [4]. While < R_ij²> reaches 1/3, the cross correlation terms approach zero. The 2nd order moments can be used for the forward signal power variation characterization. Furthermore, as suggested by Eq. (17), the 2nd order moments can be used to evaluate the average value of the elements for the round-tip propagation matrix, which will be discussed shortly after.

The 4th order moments of the forward propagation Mueller matrix elements can be useful for the evaluation of the gain fluctuation in the forward pumped amplifiers and the characterization of the 2nd order moments of back-scattered matrix elements, and hence, it is quite important to have this quantity precisely determined. Figures 2–3 demonstrate the 4th order moments of the forward propagation matrix elements which have not been reported before. Again, the analytical predictions match the numerical simulations, and the validity of the analytical formulas is confirmed again. The maximum discrepancy between the two curves for < R₁₁⁴> is 7%, which is comparable with the case of the 2nd order moments. The dynamics of the 4th order moments resembles that of the 2nd order moments, which converge to a stable value after σ²L > 2. With the 2nd order and the 4th order moments of the Mueller matrix elements, it is possible to evaluate the variance and the standard deviation of the round-tip propagation matrix elements. Two interesting features can be noticed from the figures and Eq. (18). Firstly, the 4th order moments ${R_{ij}}^2(z ){R_{i^{\prime}j^{\prime}}}^2(z )$ have different expectations, while i, j, i’ and j’ take different values. Secondly, the moments of ${R_{ij}}(z ){R_{ij}}(z ){R_{ji}}(z ){R_{jj}}(z )$are non-zero and have a negative expectation as z approaches infinity. These two features can be useful during the round-tip propagation matrix analysis and other applications.

Fig. 2. The 4th order moments of the forward propagation matrix elements: results obtained by the Monte Carlo simulations and the analytical formulas.

Download Full Size | PDF

Fig. 3. The 4th order moments of the forward propagation matrix elements: results obtained by the Monte Carlo simulations and the analytical formulas.

Download Full Size | PDF

Figures 4–5 demonstrate the previously mentioned average value and the 2nd order moments of the round-tip propagation matrix T elements, which refer to the statistical dynamics of the back-scattered signals. Again, the analytical formulas and the numerical results agree well, which demonstrates the capability of the proposed method to deal with the dynamics and the power fluctuations for the back-scattered signals. In Fig. 4, the expected value of T₁₁ is 1/3 while z approaches zero, which is in exact agreement with the results in [2] derived from a completely different approach. The value of T₁₁ can be used to characterize the back-scattered signal such as the Rayleigh scattering and the Brillouin gain efficiency in the presence of PMD [2,19]. As is known, the back scattered signal fluctuates and such fluctuation can be characterized by the 2nd order moments. The calculated < T₁₁²> in Fig. 5 and Eq. (18) suggests its value asymptotically approaches 7/15 and the variance of T₁₁ reaches a constant value as

(19)$$\begin{array}{l} {\mathop{\rm var}} ({{T_{11}}} )\\ = \left\langle {{T_{11}}^2} \right\rangle - {\left\langle {{T_{11}}} \right\rangle ^2}\\ = \frac{{16}}{{45}}, \end{array}$$

The interesting outcome suggests the back-scattered signal from one specific point of fiber has a constant fluctuation due to the PMD randomization effect. However, the Raman gain and the Brillouin gain are impacted by the integral of the round-tip matrix along the fiber, and the averaging effect along the fiber will mitigate the fluctuation [19].

Fig. 4. The comparison of average value of the round-tip propagation matrix elements (related to the 2nd order moments of the forward propagation matrix elements) obtained by the Monte Carlo simulations and the analytical formulas.

Download Full Size | PDF

Fig. 5. The comparison of 2nd order moments of the round-tip propagation matrix elements (related to the 4th order moments of the forward propagation matrix elements) obtained by the Monte Carlo simulations and the analytical formulas.

Download Full Size | PDF

A further study on the 4th order moments of the round-tip propagation matrix T elements is shown in Fig. 6. The perfect match of the curves indicates the accuracy of the proposed theory. These quantities can help to determine the pdf of the elements albeit with the numerous computational efforts to obtain the other higher order moments.

Fig. 6. The comparison of 4th order moments of the round-tip propagation matrix elements (related to the 8th order moments of the forward propagation matrix elements) obtained by the Monte Carlo simulations and the analytical formulas.

Download Full Size | PDF

4. Summary

In summary, we have proposed a comprehensive analytical theory to study the higher order statistics of the Mueller matrix in optical fibers with arbitrary length affected by the random PMD. The ODEs with constant coefficients are derived based on the coupled SDEs. Analytical formulas are found afterwards, which agree perfectly with the Monte Carlo simulation results. Asymptotic behaviors of the moments are discussed in detail with some of them exactly agree with the published literatures. The proposed method can be extremely useful in the fiber communication system analysis and the sensing applications.

Appendix

In this appendix, the detailed procedures for the conversion of the fundamental SDE from the Stratonovich sense to the Ito sense is presented. The Stratonovich sense equation reads as

(20)$$d{\textbf R} = {{\textbf T}_\beta }dz{\textbf R},$$

On the right hand side of the equation, the Stratonovich sense SDE takes the value in the middle of the integration interval, which means

(21)$$d{\textbf R} = {{\textbf T}_\beta }dz{\textbf R}\left( {z + \frac{{dz}}{2}} \right),$$

The Ito sense SDE takes the value on the left side of the integration interval. Hence, we may convert Eq. (21) by omitting the higher order terms O(dz) and using Eq. (6)

(22)$$\begin{array}{l} d{\textbf R} = {{\textbf T}_\beta }dz\left( {{\textbf R}(z )+ \frac{{d{\textbf R}}}{2}} \right)\\ = {{\textbf T}_\beta }dz{\textbf R}(z )+ \frac{1}{2}\left\langle {{{\textbf T}_\beta }dz{{\textbf T}_\beta }dz} \right\rangle {\textbf R}(z )\\ = {{\textbf T}_\beta }dz{\textbf R}(z )- {\sigma ^2}{\textbf I}dz, \end{array}$$

which is the Ito sense SDE, i.e. Eq. (7).

Funding

National Natural Science Foundation of China (61775168).

Disclosures

The authors declare no conflicts of interest.

References

1. G. J. Foschini and C. D. Poole, “Statistical theory of polarization dispersion in single mode fibers,” J. Lightwave Technol. 9(11), 1439–1456 (1991). [CrossRef]

2. F. Corsi, A. Galtarossa, and L. Palmieri, “Analytical treatment of polarization-mode dispersion in single-mode fibers by means of the backscattered signal,” J. Opt. Soc. Am. A 16(3), 574–583 (1999). [CrossRef]

3. J. P. Gordon and H. Kogelnik, “PMD fundamentals: polarization mode dispersion in optical fibers,” Proc. Natl. Acad. Sci. U. S. A. 97(9), 4541–4550 (2000). [CrossRef]

4. A. Vannucci and A. Bononi, “Statistical Characterization of the Jones Matrix of Long Fibers Affected by Polarization Mode Dispersion (PMD),” J. Lightwave Technol. 20(5), 811–821 (2002). [CrossRef]

5. A. Galtarossa and C. R. Menyuk, Polarization Mode Dispersion (Optical and Fiber Communications Reports Book 1), (Springer, 2005).

6. L. Palmieri, “Polarization properties of spun single-mode fibers,” J. Lightwave Technol. 24(11), 4075–4088 (2006). [CrossRef]

7. A. Mecozzi, “Theory of polarization mode dispersion with linear birefringence,” Opt. Lett. 33(12), 1315 (2008). [CrossRef]

8. A. Mecozzi, “A Theory of Polarization-Mode Dispersion of Spun Fibers,” J. Lightwave Technol. 27(7), 938–943 (2009). [CrossRef]

9. L. Chen, Z. Zhang, and X. Bao, “Moment-Generating Function Method Used to Evaluate the Performance of a Linear Optical Communication System,” J. Lightwave Technol. 27(16), 3399–3409 (2009). [CrossRef]

10. W. T. Ireeta, V. Musara, L. Wu, and A. W. R. Leitch, “Statistical characterization of the output state of polarization in deployed single-mode fibers,” Opt. Lett. 35(12), 2049–2051 (2010). [CrossRef]

11. M. Gazeau, “Numerical simulation of nonlinear pulse propagation in optical fibers with randomly varying birefringence,” J. Opt. Soc. Am. B 30(9), 2443–2451 (2013). [CrossRef]

12. L. Zhang, A. A. E. Hajomer, X. Yang, and W. Hu, “Error-free secure key generation and distribution using dynamic Stokes parameters,” Opt. Express 27(20), 29207–29216 (2019). [CrossRef]

13. T. Wang, P. Huang, S. Wang, and G. Zeng, “Polarization-state tracking based on Kalman filter in continuous-variable quantum key distribution,” Opt. Express 27(19), 26689–26700 (2019). [CrossRef]

14. X. Wang, J. Hu, F. Wang, Y. Yong, Y. Zhang, M. Xue, X. Zhang, and S. Pan, “Multi-vibration detection by probe pulses with ergodic SOPs in a POTDR system,” Opt. Express 26(22), 28349–28362 (2018). [CrossRef]

15. D. Marcuse, C. R. Menyuk, and P. K. A. Wai, “Application of the Manakov-PMD equations to studies of signal propagation in fibers with randomly-varying birefringence,” J. Lightwave Technol. 15(9), 1735–1746 (1997). [CrossRef]

16. S. Said and N. L. Bihan, “Higher-order statistics of Stokes parameters in a random birefringent medium,” Waves in Random and Complex Media 18(2), 275–292 (2008). [CrossRef]

17. K. Hotate, K. Abe, and K. Y. Song, “Suppression of Signal Fluctuation in BrillouinOptical Correlation Domain Analysis System Using Polarization Diversity Scheme,” IEEE Photonics Technol. Lett. 18(24), 2653–2655 (2006). [CrossRef]

18. Q. Lin and G. P. Agrawal, “Vector theory of stimulated Raman scattering and its application to fiber-based Raman amplifiers,” J. Opt. Soc. Am. B 20(8), 1616–1631 (2003). [CrossRef]

19. A. Zadok, E. Zilka, A. Eyal, L. Thévenaz, and M. Tur, “Vector analysis of stimulated Brillouin scattering amplification in standard single-mode fibers,” Opt. Express 16(26), 21692–21707 (2008). [CrossRef]

20. H. Dong, P. Shum, M. Yan, J. Q. Zhou, G. X. Ning, Y. D. Gong, and C. Q. Wu, “Generalized Mueller matrix method for polarization mode dispersion measurement in a system with polarization-dependent loss or gain,” Opt. Express 14(12), 5067–5072 (2006). [CrossRef]

21. V. Devlaminck, “Physical model of differential Mueller matrix for depolarizing uniform media,” J. Opt. Soc. Am. A 30(11), 2196–2204 (2013). [CrossRef]

22. R. Ossikovski and O. Arteaga, “Statistical meaning of the differential Mueller matrix of depolarizing homogeneous media,” Opt. Lett. 39(15), 4470–4473 (2014). [CrossRef]

23. J. Zhou and Q. Hu, “Supplement L0,” figshare (2020) [retrieved 22 September 2020], https://doi.org/10.6084/m9.figshare.12988577 .

24. J. Zhou and Q. Hu, “Supplement N0,” figshare (2020) [retrieved 22 September 2020], https://doi.org/10.6084/m9.figshare.12988919 .

Higher order statistics of the Mueller matrix in a fiber with an arbitrary length impacted by PMD

Abstract

1. Introduction

2. Theory

2.1 Fundamental propagation model in the presence of PMD

2.2 Derivation of the fundamental SDEs

2.3 Derivation of the ODEs for the 2nd order and the 4th order statistics

2.4 Asymptotic behavior of the moments: validity check

3. Results and discussions

4. Summary

Appendix

Funding

Disclosures

References

Supplementary Material (2)

Cited By

Figures (6)

Equations (23)

Optics Express

Name	Description
Dataset 1	Supplement L0
Dataset 2	The expression of N0.