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

Two-stage equalization for ultra-fast RSOP and inter symbol interference compensation based on a simplified Kalman filter

Open Access Open Access

Abstract

We propose a two-stage equalization based on a simplified Kalman filter, which is used to solve the rapid rotation of the state of polarization (RSOP) that is caused by lightning strikes on optical cables and the extra inter symbol interference (ISI) introduced in the system. By analyzing the special expression of matrix coefficient in the Kalman filter under polarization demultiplexing, the simplified idea of a Kalman filter is provided, and its updating process is transformed into a kind of multiple-input-multiple-output (MIMO) structure algorithm. At the same time, the second stage finite impulse response filter is used to solve the ISI that is difficult to be solved by a Kalman filter. The performance of the proposed algorithm was tested in a coherent system of 28Gbaud PDM-QPSK/16QAM. The results confirm that on the basis of lower complexity than a Kalman filter, the proposed algorithm reduces its complexity by more than 30% compared to traditional MIMO equalization algorithm under the premise of linear operation, and which also can handle RSOP of 20 Mrad/s. When the system suffers from the extra ISI due to the limited device bandwidth, the optical signal to noise ratio of the proposed algorithm is about 4 dB lower than the Kalman filter at the same bit error rate.

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

1. Introduction

The polarization division multiplexing (PDM) system that expands capacity by transmitting two polarization signals has been widely applied. The existing 100 G systems generally use PDM system [1]. However, after the introduction of PDM, polarization related signal damage occurs which mainly including polarization mode dispersion (PMD), rotation of state of polarization (RSOP), and polarization dependent loss (PDL) [2].

The current mainstream equalization algorithm used in commercial solutions is the classic constant modulus algorithm (CMA) based on multiple-input-multiple-output (MIMO) structure [3]. However, when the Optical Fiber Composite Overload Ground Wire (OPGW) encounter extreme weather such as thunderstorm season, the fiber transmission link will experience signal interruption due to the rapid changes in the state of polarization (SOP) caused by the Faraday rotation effect [4,5], which further leads to bit error (BER) at the receiver. Previous experiments have shown that RSOP caused by lightning strikes can reach Mrad/s [6]. In this case, conventional algorithms such as CMA cannot handle fast RSOP and leading to performance degradation or even failure [7].

At present, there have been many studies on fast polarization demultiplexing, which can be roughly divided into three types: adaptive filtering algorithms based on MIMO structure [8,9]; Kalman filter structure, including extended Kalman filter [1014], radius guided Kalman filter [15]; Stokes space [1618]. The problem with adaptive filtering algorithms based on MIMO structure is that they cannot cope with fast RSOP; The Kalman filter structure currently has high polarization demultiplexing ability, but its structure is too complex [12]. The nonlinear operation process of its matrix is not conducive to commercial implementation, and Kalman cannot handle inter symbol interference (ISI) outside of PMD and dispersion well; The Stokes space algorithm also has the problem of complex algorithm structure.

By analyzing the special expression of matrix coefficients in Kalman filter under polarization demultiplexing, this paper provides a simplified idea of Kalman filter, transforms its updating process into a kind of MIMO structure algorithm, and in order to compensate the residual dispersion and ISI of the system, two finite impulse response (FIR) filters are cascaded subsequently. Through the coherent system of 28Gbaud PDM-QPSK and 28Gbaud PDM-16QAM, we have confirmed that the simplified Kalman filter can cope with the fast RSOP of more than 20 Mrad/s, When the system suffers from the extra ISI due to the limited device bandwidth, optical signal to noise ratio (OSNR) of the proposed algorithm is about 4 dB lower than Kalman filter at the same BER. At the same time, the simplified Kalman joint FIR filter structure is simple addition, subtraction, multiplication and division linear operation as the MIMO algorithm, and reduces the complexity by more than 30% compared with the complex Kalman filter, it is more conducive to commercial implementation.

2. Algorithm principle

2.1 Principle of simplified Kalman filter

There are different damages in optical transmission systems, such as amplified self-emission noise (ASE), chromatic dispersion (CD), frequency offset (FO), and phase noise (PN). For PDM, the signal distortion caused by polarization effects mainly composed of RSOP, PMD, and PDL will be considered a key part of the overall damage [19]. Compared to static PDL damage, this paper focuses more on fast RSOP caused by lightning strikes and ISI introduced by the system. In order to better study the polarization damage of signals, it is necessary to establish a model.

For Kalman filter, the main purpose is to obtain the damage inverse matrix of the signal and select the Kalman state vector. The existing literature generally selects the compensation matrix corresponding to three parameters according to the RSOP damage matrix form [13]:

$${J_{RSOP}} = \left( {\begin{array}{{cc}} {cos\kappa {e^{j\alpha }}}&{ - sin\kappa {e^{j\beta }}}\\ {sin\kappa {e^{ - j\beta }}}&{cos\kappa {e^{ - j\alpha }}} \end{array}} \right)$$
where $\kappa$, and $\alpha ,\beta $ denote the azimuth and phase angles .Since it contains trigonometric function and exponent, which is not conducive to the simplification of calculation, the compensation matrix and Kalman state vector are selected as:
$${J^{ - 1}} = \left( \begin{array}{c} a - jb\;\;\;\;\;\;\;\;\;\;c + jd\\ - c + jd\;\;\;\;\;\;\;a + jb \end{array} \right)$$
$${x_k} = {({a,b,c,d} )^T}$$
where $a,b,c,d$ are the coefficients of the polarization damage matrix, k represents the symbol number. The original Kalman filtering process can be described as follows: firstly, compensate the signal based on the initial value of the state vector:
$$\left( \begin{array}{l} {E_{xout}}\\ {E_{yout}} \end{array} \right) = \left( \begin{array}{c} a - jb\;\;\;\;\;\;\;\;\;\;c + jd\\ - c + jd\;\;\;\;\;\;\;a + jb \end{array} \right)\left( \begin{array}{l} {E_{xin}}\\ {E_{yin}} \end{array} \right)$$
where ${E_{xin}},{E_{yin}}$ is the input signal of Kalman filter and ${E_{xout}},{E_{yout}}$ is the output signal of filter. Secondly select the error function, where the QPSK signal can be expressed as:
$$\left\{ \begin{array}{l} ex = ({|{E_{xout}}(k )|{\;^2} - {R^2}} )\\ ey = ({|{E_{yout}}(k )|{\;^2} - {R^2}} )\end{array} \right.$$
where $R$ is 1, $ex,ey$ is the difference between the output signal modulus and the standard modulus $R$.${E_{xout}},{E_{yout}}$ is the output of the Kalman filter, Kalman's observation matrix ${h_k}$ can be represented as:
$${h_k} = \left( \begin{array}{l} 0\\ 0 \end{array} \right) - \left( \begin{array}{c} ex\\ ey \end{array} \right)$$
The form of Jacobian matrix $H$ is:
$$H = \left( {\begin{array}{{cc}} {\frac{{\partial ex}}{{\partial a}}}&{\frac{{\partial ex}}{{\partial b}}\;\;\;\;\frac{{\partial ex}}{{\partial c}}\;\;\;\;\frac{{\partial ex}}{{\partial d}}}\\ {\frac{{\partial ey}}{{\partial a}}}&{\frac{{\partial ey}}{{\partial b}}\;\;\;\;\frac{{\partial ey}}{{\partial c}}\;\;\;\;\frac{{\partial ey}}{{\partial d}}} \end{array}} \right) = \left( {\begin{array}{{cc}} {xa}&{xb\;\;\;\;xc\;\;\;\;xd}\\ {ya}&{yb\;\;\;\;yc\;\;\;\;yd} \end{array}} \right)$$
where $xa,xb,xc,xd$ represents the derivative of the observation matrix $ex$ to the state vector respectively, $ya,yb,yc,yd$ represent the derivative of the observation matrix $ey$ to the state vector respectively, and then update the covariance matrix, Kalman gain, state vector:
$${P_k} = {P_k} + Q$$
$${K_k} = {P_k}{H^T}{(H{P_k}{H^T} + R)^{ - 1}}$$
$${x_k} = {x_k} + {K_k}{h_k}$$
$${P_k} = (I - {K_k}H){P_k}$$
where $Q$ represents the covariance matrix of state process noise. The main difficulty in implementing Kalman in communication lies in its complex matrix operations. The algorithm calculation process of MIMO structure is completely linear, so it occupies the mainstream in digital signal processing (DSP). If the Kalman operation process can be linearized, it will greatly reduce its computational difficulty.

In general, the RSOP damage of the system is very small, usually at the Krad/s level. Literature shows that when lightning occurs around the optical cable, the Faraday rotation effect can cause rapid changes in the SOP of the optical signal with RSOP up to 8 Mrad/s [6]. Although this is a very large value, its variation is very small for the two symbols in a coherent communication system.

When RSOP is 10 Mrad/s, the damage matrix of two adjacent symbols can be expressed as follows:

$${J_{Symbo{l_1}}} = \left( {\begin{array}{{cc}} {\cos {k_1}{e^{j{\alpha_1}}}}&{ - \sin {k_1}{e^{j{\beta_1}}}}\\ {\sin {k_1}{e^{ - j{\beta_1}}}}&{\cos {k_1}{e^{ - j{\alpha_1}}}} \end{array}} \right)$$
$${J_{Symbo{l_2}}} = \left( {\begin{array}{{cc}} {\cos {k_2}{e^{j{\alpha_2}}}}&{ - \sin {k_2}{e^{j{\beta_2}}}}\\ {\sin {k_2}{e^{ - j{\beta_2}}}}&{\cos {k_2}{e^{ - j{\alpha_2}}}} \end{array}} \right)$$
where $Symbo{l_1},Symbo{l_2}$ respectively represent two adjacent symbols, ${k_1},{\alpha _1},{\beta _1}$ and ${k_2},{\alpha _2},{\beta _2}$ denote the azimuth and phase angles for polarization damage matrix coefficients of two symbols. When RSOP is 20 Mrad/s, due to the Faraday rotation effect, $k$ will change rapidly with $\alpha ,\beta $ relatively small changes. In order to obtain the maximum matrix coefficient change, the exponential function in the damage matrix can be set to 1. For ${k_2}$, the RSOP rotation speed can be used to obtain:
$${k_2} - {k_1} = {V_{RSOP}} \ast {T_{Symbol}}$$
where ${T_{Symbol}}$ and ${V_{RSOP}}$ represents the symbol interval time and RSOP speed respectively. When the system communication rate is 100 G, for the PDM-QPSK signal, the symbol interval is 40ps. At this time, the maximum value of the damage matrix coefficient change can be calculated:
$$\begin{aligned} Max|{({\cos {k_2} - \cos {k_1}} )} |&= Max\left|{ - 2\sin (\frac{{{k_2} + {k_1}}}{2})\ast \sin (\frac{{{k_2} - {k_1}}}{2})} \right|\\ &= 2\sin \left( {\frac{{{V_{RSOP}} \ast {T_{Symbol}}}}{2}} \right) \approx 8e - 4 \end{aligned}$$
It can be seen that even in large RSOP cases, the maximum difference in the damage matrix coefficients between two symbols is very small. Therefore, for the initial state selection of Kalman filter, the values of the process noise covariance matrix Q of the state vector and the mean squared error matrix P of the state variable estimate are generally set to be very small. For the QPSK signal, the error between signal normalization amplitude and unit circle will greatly exceed 8e-4 due to OSNR which is why the initial value of the noise covariance matrix R of the Kalman observation matrix is much larger than Q and P. Two conclusions can be drawn from this:

(1) The matrix elements of P must be a minimal quantity in the iterative process of Kalman;

(2) $R \gg P,Q$.

When the Kalman filter is approximately considered to complete the polarization demultiplexing, the two signals have no crosstalk and the variables of their observation matrix are relatively independent. For the Kalman state vector $a,b,c,d$, it is evolved from the three parameter RSOP model. In fact, when the power of two signals fluctuate, then the value of ${a^2} + {b^2} + {c^2} + {d^2}$ is not completely equal to 1, so it can be approximately considered that its update process is relatively independent.

Based on the above analysis, there is a theoretical basis for selecting the initial value of Kalman. We can take $P = diag({\mu a,\mu b,\mu c,\mu d} )$, $Q = ({\mu ,\mu ,\mu ,\mu } )$ and $R = diag(C,C)$ where C is much greater than $\mu $. Therefore, the initial iteration process of Kalman and the first update process of covariance matrix can be written as:

$${P_k} = {P_k} + Q = diag\left( \begin{array}{l} \mu a + \mu \\ \mu b + \mu \\ \mu c + \mu \\ \mu d + \mu \end{array} \right) = diag\left( \begin{array}{l} \mu a\\ \mu b\\ \mu c\\ \mu d \end{array} \right)$$
The calculation formula for traditional Kalman gain is:
$${K_k} = {P_k}{H^T}{(H{P_k}{H^T} + R)^{ - 1}}$$
This can be simplified in two ways:

  • (1) The meaning of Kalman gain:

For the derivation of the original Kalman:

$$\left\{ \begin{array}{l} \hat{x}k = \hat{x}{k^ - } + {K_k}({yk - H\hat{x}{k^ - }} )\\ yk = Hxk + {v_k} \end{array} \right.$$
where ${v_k}$ is the Gaussian noise vector of the observation process, $\hat{x}{k^ - },\hat{x}k$ are the prior and posterior estimates of the state vector respectively, $xk$ is true value of the state vector, and $yk$ is the observation matrix, ${K_k}$ is actually a trade-off between the observed value and the predicted value. When the observed value of the system is very reliable and the mean squared error ${P_k}$ of the state estimation is large, $H{P_k}{H^T} \gg R$, then the Kalman gain will tend to ${H^{ - 1}}$, and $\hat{x}k$ will tend to ${H^{ - 1}}yk$, the state estimation is more dependent on the observed value; For the polarization demultiplexing, according to the previous analysis $R \gg P,Q$, the coefficient change of adjacent symbols is very small, the difference between the current state vector and the previous state vector is very small, then obviously $\hat{x}k$ will prefer to trust $\hat{x}{k^ - }$, in other words, $H{P_k}{H^T} \ll R$, the Kalman gain tends to ${P_k}{H^T}{R^{ - 1}}$ and because $R$ is a constant coefficient diagonal matrix in the polarization demultiplexing, the expression of Kalman gain has a more concise form:
$${K_k} \approx {P_k}{H^T}{R^{ - 1}} = {P_k}{H^T}{\left( \begin{array}{l} C\textrm{ 0}\\ 0\textrm{ C} \end{array} \right)^{ - 1}} = N{P_k}{H^T}$$
where N can be considered as the gain factor introduced by $R$ for ${P_k}{H^T}$.
  • (2) Mathematical analysis
Substituting $H$ and ${P_k}$, $(H{P_k}{H^T} + R)$ can be expressed as:
$$(H{P_k}{H^T} + R) = \left( \begin{array}{l} \mu x + C\textrm{ }\mu xy\\ \mu xy\textrm{ }\mu y + C \end{array} \right)$$
$$\left\{ \begin{array}{l} \mu x = x{a^2} \cdot \mu a + x{b^2} \cdot \mu b + x{c^2} \cdot \mu c + x{d^2} \cdot \mu d\\ \mu xy = xa \cdot ya \cdot \mu a + xb \cdot yb \cdot \mu b + xc \cdot yc \cdot \mu c + xd \cdot yd \cdot \mu d\textrm{ }\\ \mu y = y{a^2} \cdot \mu a + y{b^2} \ast \mu b + y{c^2} \cdot \mu c + y{d^2} \cdot \mu d \end{array} \right.$$
In the previous analysis, we clarified that the elements in ${P_k}$ are much smaller than $R$.In other words there is $\mu x,\mu y,\mu xy \ll C$:
$${K_k} = {P_k}{H^T}{\left( \begin{array}{l} \mu x + C\textrm{ }\mu xy\\ \mu xy\textrm{ }\mu y + C \end{array} \right)^{ - 1}} \approx {P_k}{H^T}{\left( \begin{array}{l} C\textrm{ }0\\ 0\textrm{ }C \end{array} \right)^{ - 1}} = N{P_k}{H^T}$$
Based on the above approximation of Kalman gain, the specific form of ${K_k}$ can be written:
$${K_k} = N\left( \begin{array}{l} \mu a \cdot xa\textrm{ }\mu a \cdot ya\\ \mu b \cdot xb\textrm{ }\mu b \cdot yb\\ \mu c \cdot xc\textrm{ }\mu c \cdot yc\\ \mu d \cdot xd\textrm{ }\mu d \cdot yd \end{array} \right)$$
where $xa,xb,xc,xd,ya,yb,yc,yd$ respectively represents the derivative of the observation matrix error function to the state vector. For QPSK signal, it can be expressed as:
$$\left( \begin{array}{l} xa\\ xb\\ xc\\ xd \end{array} \right) = \left( \begin{array}{c} {E_{xin}}{E_{xout}}^ \ast{+} {E_{xout}}{E_{xin}}^ \ast \\ j({{E_{xin}}{E_{xout}}^ \ast{-} {E_{xout}}{E_{xin}}^ \ast } )\\ {E_{yin}}{E_{xout}}^ \ast{+} {E_{xout}}{E_{yin}}^ \ast \\ j({{E_{yin}}{E_{xout}}^ \ast{-} {E_{xout}}{E_{yin}}^ \ast } )\end{array} \right) = \left( \begin{array}{c} 2Re({E_{xin}}{E_{xout}}^ \ast )\\ - 2Im({E_{xin}}{E_{xout}}^ \ast )\\ 2Re({E_{yin}}{E_{xout}}^ \ast )\\ - 2Im({E_{yin}}{E_{xout}}^ \ast ) \end{array} \right)$$
$$\left( \begin{array}{l} ya\\ yb\\ yc\\ yd \end{array} \right) = \left( \begin{array}{c} {E_{yin}}{E_{yout}}^ \ast{+} {E_{yout}}{E_{yin}}^ \ast \\ - j({{E_{yin}}{E_{yout}}^ \ast{-} {E_{yout}}{E_{yin}}^ \ast } )\\ - ({{E_{xin}}{E_{yout}}^ \ast{+} {E_{yout}}{E_{xin}}^ \ast } )\\ j({{E_{xin}}{E_{yout}}^ \ast{-} {E_{yout}}{E_{xin}}^ \ast } )\end{array} \right) = \left( \begin{array}{c} 2Re({E_{yin}}{E_{yout}}^ \ast )\\ 2Im({E_{yin}}{E_{yout}}^ \ast )\\ - 2Re({E_{xin}}{E_{yout}}^ \ast )\\ - 2Im({E_{xin}}{E_{yout}}^ \ast ) \end{array} \right)$$
where $Re$ and $Im$ represent the real and imaginary parts respectively. After obtaining the Kalman gain, the updated expression can be obtained:
$${x_k} = {x_k}+ {K_k}{h_k}$$
$$x = \left( \begin{array}{l} a\\ b\\ c\\ d \end{array} \right) = \left( \begin{array}{l} a^{\prime}\\ b^{\prime}\\ c^{\prime}\\ d^{\prime} \end{array} \right) + N\left( \begin{array}{c} \mu a(xa \cdot ex + ya \cdot ey)\\ \mu b(xb \cdot ex + yb \cdot ey)\\ \mu c(xc \cdot ex + yc \cdot ey)\\ \mu d(xd \cdot ex + yd \cdot ey) \end{array} \right)$$
$${P_k} = (I - {K_k}H){P_k}$$
It can be seen that Eq. (28) is update of ${P_k}$. If it is linked to the MIMO structure algorithm, it is not difficult to find that the element of ${P_k}$ actually corresponds to the update step in the MIMO structure algorithm which is constantly changing in the Kalman filter. Based on the Kalman gain, Jacobian matrix, the update process of ${P_k}$ can be obtained as follows:
$$\left( \begin{array}{l} 1 - \mu a(x{a^2} + y{a^2})\textrm{ }\mu a(xaxb + yayb)\textrm{ }\mu a(xaxc + yayc)\textrm{ }\mu a(xaxd + yayd)\\ \mu b(xaxb + yayb)\textrm{ }1 - \mu b(x{b^2} + y{b^2})\textrm{ }\mu b(xcxb + ycyb)\textrm{ }\mu b(xdxb + ydyb)\\ \mu c(xaxc + yayc)\textrm{ }\mu c(xcxb + ycyb)\textrm{ }1 - \mu c(x{c^2} + y{c^2})\textrm{ }\mu c(xcxd + ycyd)\\ \mu d(xaxd + yayd)\mu d(xdxb + ydyb)\textrm{ }\mu d(xcxd + ycyd)\textrm{ }1 - \mu d(x{d^2} + y{d^2}) \end{array} \right){P_k}$$
The update of covariance matrix in Kalman filter is very complex, and the simplification of ${P_k}$ update can be analyzed from two aspects:
  • (1) Combining MIMO algorithm:
From the update process of $a,b,c,d$ in Eq. (27), it can be seen that it is very similar to the process of CMA update coefficient. At this point, $\mu a,\mu b,\mu c,\mu d$ is considered as the update step size in CMA, and ${P_k} = (I - {K_k}H){P_k}$ is considered as the update of step size. Therefore, it is obvious that the update of step size should be relatively independent:
$$\left\{ \begin{array}{l} \mu a = \mu a({1 - N\mu a(x{a^2} + y{a^2})} )\\ \mu b = \mu b({1 - N\mu b(x{b^2} + y{b^2})} )\\ \mu c = \mu c({1 - N\mu c(x{c^2} + y{c^2})} )\\ \mu d = \mu d({1 - N\mu d(x{d^2} + y{d^2})} )\end{array} \right.$$
At this time, Kalman filter is simplified into a MIMO like algorithm that can dynamically update the step size.
  • (2) Mathematical analysis:
Based on the previous analysis, the elements in the matrix of $P,Q$ must be very small during the iteration process. Due to the addition of constant 1, the elements on the ${P_k}$ diagonal are obviously larger than their non diagonal elements, so it can be approximated that the ${P_k}$ diagonal elements dominate. Therefore, the update process of ${P_k}$ can be represented as:
$${P_k} \approx diag\left( \begin{array}{l} 1 - N\mu a(x{a^2} + y{a^2})\\ 1 - N\mu b(x{b^2} + y{b^2})\\ 1 - N\mu c(x{c^2} + y{c^2})\\ 1 - N\mu d(x{d^2} + y{d^2}) \end{array} \right){P_k} = diag\left( \begin{array}{l} \mu a({1 - N\mu a(x{a^2} + y{a^2})} )\\ \mu b({1 - N\mu b(x{b^2} + y{b^2})} )\\ \mu c({1 - N\mu c(x{c^2} + y{c^2})} )\\ \mu d({1 - N\mu d(x{d^2} + y{d^2})} )\end{array} \right)$$
Therefore, the overall process of simplified Kalman filter (S-KF) can be expressed as shown in Fig. 1.

 figure: Fig. 1.

Fig. 1. S-KF compared to KF algorithm process.

Download Full Size | PDF

It can be seen that the simplified Kalman filter is more like a single tap MIMO structure algorithm, but compared with the MIMO algorithm, the simplified Kalman filter has the following advantages:

(1) The update step size corresponding to the four state vectors is not fixed which changes with the process of Kalman iteration, and the principle of step size change is based on Kalman filter, so that it always faces the direction of small error.

(2) Compared with the traditional Kalman filter, it greatly reduces the computational complexity and makes commercial application possible, not only because its linear computing structure similar to CMA, but also for the greatly reduced computational complexity compared with CMA due to single symbol computing.

2.2 Simplified Kalman filter joint FIR

From the above analysis, the simplified Kalman filter can be regarded as a single tap algorithm similar to MIMO structure. In fact, the signal not only suffers from polarization damage, but also suffers from ISI due to the bandwidth of system devices.

Kalman filter can be modeled according to different damages. Zhang et al. established damage models for RSOP, PMD and CD [13,14], but it is generally difficult to model the ISI introduced by device bandwidth. The general idea of handling ISI is use a multi tap FIR filter, so in order to compensate this part of ISI, we can cascade two multi tap FIR filters later to achieve the compensation of ISI. After the first stage simplified Kalman filter, the signal becomes two parallel signals, and then it is input to two N-order FIR filters for ISI and certain CD compensation. Then the first stage simplified Kalman filter continuously inputs data to the FIR filter. The FIR filter intercepts data and generates balanced data, the overall process is carried out synchronously as a whole (Fig. 2).

 figure: Fig. 2.

Fig. 2. S-KF joint FIR filter.

Download Full Size | PDF

It can be seen that the overall equalization algorithm is composed of two parts. The first part is based on the simplified Kalman filter to compensate for the polarization damage of the signal, and the second part is used to compensate for the signal damage except for the polarization damage, mainly ISI, which is equivalent to a two-stage CMA structure, and has extremely strong ability to respond to rapid RSOP because of learning from the Kalman filter idea.

After the first stage of polarization demultiplexing, the signal is input into two subsequent N-order FIR filters. If the number of taps selected is L, the output of the FIR filter can be expressed as:

$$\left\{ \begin{array}{l} {E_{xout}}(k )= {h_{xx}}^k{E_{xmid}}(l )\\ {E_{yout}}(k )= {h_{yy}}^k{E_{ymid}}(l )\end{array} \right.$$
where ${E_{xmid}},{E_{ymid}}$ is polarization demultiplexed data output by the simplified Kalman filter, ${E_{xout}},{E_{yout}}$ is the final output of the equalization.${h_{xx}},{h_{yy}}$ is the tap coefficient of the FIR filter, if $L$ is 11, then the range of l is:$K - 5:K + 5$, and the error function of the FIR filter can be expressed as:
$$\left\{ \begin{array}{l} ex = ({|{E_{xout}}(k )|{\;^2} - {R^2}} )\\ ey = ({|{E_{yout}}(k )|{\;^2} - {R^2}} )\end{array} \right.$$
For QPSK signal, R is 1. For 16QAM signal, MMA algorithm is used. If the update step size is $\mu$, the update formula of tap coefficient obtained by gradient descent is:
$$\left\{ \begin{array}{l} h_{xx}^{k + 1} = h_{xx}^k + \mu ey(k ){E_{xout}}(k ){E_{xmid}}^ \ast (l )\\ h_{yy}^{k + 1} = h_{yy}^k + \mu ey(k ){E_{yout}}(k ){E_{ymid}}^ \ast (l )\end{array} \right.$$
The signal compensates ISI at the second stage, thus completing all the equalization processes.

3. Simulation results

3.1 Simulation platform setup

To verify the performance of the proposed simplified Kalman joint FIR filter equalization algorithm for fast RSOP and ISI, we established a 28Gbaud PDM-QPSK and 28Gbaud PDM-16QAM coherent system as simulation platforms which can be seen in Fig. 3. The light emitted by the originating laser is divided into two polarized beams by a polarization beam splitter, then the two polarized beams are modulated into two orthogonal signals by IQ modulators, and pulse shaping is performed in the transmitting, the roll off coefficient is 0.1. It should be noted that the optical signal is modulated by four electrical signals, and then the signal is transmitted through the optical fiber of 80 km. The dispersion coefficient of the SMF is 16.9ps/nm/km, the polarization mode dispersion coefficient is 0.2ps/sqrt(km), the nonlinear coefficient is 1.1/W/km, and the attenuation coefficient is 0.2 dB/km. The introduction of RSOP is to multiply the sampled complex signal with the formula in the literature to simulate the polarization damage encountered in the optical fiber link [19]. The optical signal is coherently received by the local oscillator laser in the receiver, introducing the frequency offset between the transmitting laser and the local oscillator laser and the phase noise introduced by the laser linewidth, 100kHz carrier phase noise is introduced due to laser linewidth and 100 MHz frequency offset is introduced between the lasers at transmitter and receiver. After coherent reception, four electrical signals are generated by the photodetector (PD) and then merge into two electrical signals enter the digital signal processing (DSP) module. The 3 dB bandwidth of IQ modulator and photodetector will introduce additional ISI as simulation parameters to compare the algorithm performance. The proposed algorithm works after dispersion compensation and before frequency offset and phase noise estimation.

 figure: Fig. 3.

Fig. 3. Simulation platform.

Download Full Size | PDF

3.2 Different algorithm parameter settings

For the proposed simplified Kalman filter, the initial value of its state vector is set as ${x_k} = {(1,0,0,0)^T}$. According to the previous analysis, the values of the Covariance matrix P and Q of the Kalman filter should be very small during the depolarization process, so set $P = diag(2{e^{ - 5}},2{e^{ - 5}},2{e^{ - 5}},2{e^{ - 5}})$, $Q = diag(2{e^{ - 5}},2{e^{ - 5}},2{e^{ - 5}},2{e^{ - 5}})$,where the gain factor is chosen as $N = 10$. For the simplified Kalman filter, the number of taps of the N-order FIR filter is set to 15, and the update step is set to 3e-4.

For the CMA algorithm, the number of taps is set to 21, and the error update step is set to 3e-4. For the original Kalman filter, we set the state vector as:

$${x_k} = {({\tau _1},{\tau _2},{\tau _3},\alpha ,\beta ,\kappa )^T}$$
The selection of covariance matrix is $P = diag(2{e^{ - 5}},2{e^{ - 5}},2{e^{ - 5}},2{e^{ - 5}},2{e^{ - 5}},2{e^{ - 5}})$, $Q = diag(2{e^{ - 5}},$ $2{e^{ - 5}},2{e^{ - 5}},2{e^{ - 5}},2{e^{ - 5}},2{e^{ - 5}})$. ${\tau _1},{\tau _2},{\tau _3}$ is the polarization mode dispersion compensation factor.

Based on the above parameters, we simulated the coherent system of 28Gbaud PDM-QPSK and 28Gbaud PDM-16QAM with a symbol length of 2^16. The damaged signals were equalized using simplified Kalman filter, Kalman filter, and CMA/MMA, respectively. And then after frequency offset estimation and carrier phase recovery, the final bit error rate (BER) was used as the evaluation criterion for algorithm performance.

Figures 4(a) and 5(a) show the BER of three algorithms in different RSOP size with OSNR of 15 dB for PDM-QPSK and 22 dB for PDM-16QAM. Among them, RSOP is set to (0:19) Mrad/s. For QPSK systems, the BER of CMA has exceeded the FEC threshold when RSOP is 1.5 Mrad/s. Kalman and simplified Kalman can still polarization demultiplexing when RSOP is 20 Mrad/s. For the 16QAM system, the Kalman exceeds the 3.8e-3 FEC threshold, while the simplified Kalman filter can still remain below the FEC threshold. The main reason is that the Kalman filter unable to model the ISI introduced by the device bandwidth because of its structure, while the simplified Kalman filter has a good equalization effect on ISI because of its subsequent addition of FIR filter; Fig. 4(b) and Fig. 5(b) show the change of BER of Kalman filter and simplified Kalman structure with OSNR under fixed RSOP where RSOP is 2 Mrad/s, 10 Mrad/s and 20 Mrad/s respectively. When FEC threshold is 3.8e-3, OSNR tolerance of simplified Kalman for PDM-QPSK and PDM-16QAM systems is about 4 dB higher than Kalman, which shows that simplified Kalman can have higher RSOP margin.

 figure: Fig. 4.

Fig. 4. The performance of different algorithms under PDM-QPSK signals, (a) BER vs. RSOP.

Download Full Size | PDF

 figure: Fig. 5.

Fig. 5. The performance of different algorithms under PDM-16QAM signals, (a) BER vs. RSOP, (b) OSNR vs. RSOP.

Download Full Size | PDF

In order to verify the tolerance of simplified Kalman filter to ISI, we set the 3 dB bandwidth of IQ modulator and PD at the transmitter of the simulation system to 0.8 times and 0.7 times of the Baud of the signal respectively to introduce additional ISI, and test the restriction of the algorithm to the bandwidth.

Figure 6(a) and Fig. 7(a) show the change of BER of different algorithms with different RSOP when the 3 dB bandwidth of IQ modulator and PD is 2 × 0.8 × 28 Gbaud and 2 × 0.7 × 28 Gbaud respectively. It can be seen that the simplified Kalman filter has a higher tolerance for bandwidth. For PDM-QPSK signals, BER of Kalman filter has exceeded the FEC threshold when the bandwidth is 2 × 0.7 × 28 Gbaud, and the simplified Kalman filter also have good performance. For PDM-16QAM signals, The Kalman filter exceeds the FEC threshold when the bandwidth is 2 × 0.8 × 28 Gbaud, and the simplified Kalman filter is below the FEC threshold within the RSOP range. Figure 6(b) and Fig. 7(b) are the simplified Kalman filter and Kalman filter performance of BER with OSNR under different bandwidths. It can be seen that in the PDM-QPSK signal, under the same BER, the OSNR tolerance of the simplified Kalman filter is 3 dB higher than that of the Kalman filter, and the OSNR tolerance will be higher as the bandwidth decreases; In the PDM-16QAM system, when the bandwidth is 2 × 0.7 × 28 Gbaud, the OSNR tolerance of the simplified Kalman is 6 dB higher than that of the Kalman, which shows the excellent performance of the simplified Kalman filter in dealing with bandwidth constraints.

 figure: Fig. 6.

Fig. 6. The performance of different algorithms under PDM-QPSK signals, (a) BER vs. RSOP with different bandwidths, (b) OSNR vs. RSOP with different bandwidths.

Download Full Size | PDF

 figure: Fig. 7.

Fig. 7. The performance of different algorithms under PDM-16QAM signals, (a) BER vs. RSOP with different bandwidths, (b) OSNR vs. RSOP with different bandwidths.

Download Full Size | PDF

The simplified Kalman filter has a structure similar to the MIMO structure, but its update step size will change constantly, and its change degree will accelerate with the increase of RSOP. Figure 8 shows the changes in the four update step sizes of the simplified Kalman filter under PDM-16QAM signals. It can be seen that the step size changes will become more and more severe as the RSOP increases, but the range of changes is roughly the same. This is similar to the previous discussion. To quickly track, it is not necessary to pursue a large update step size, but more importantly, the step size provided at each moment is adapted to the effect of the RSOP on the signal at that time is more important. The simplified Kalman filter achieves this perfectly and has excellent depolarization ability.

 figure: Fig. 8.

Fig. 8. Change of update step size of simplified Kalman filter under different RSOP of PDM-16QAM signals: (a) RSOP = 2 Mrad/s, (b) RSOP = 10 Mrad/s, (c) RSOP = 20 Mrad/s.

Download Full Size | PDF

3.3 Algorithm complexity analysis

This section will analyze the complexity of the simplified Kalman filter algorithm and the complexity of CMA under the PDM-QPSK signal, because these two algorithms are linear operations and very suitable for commercial implementation. We will use the number of multipliers and adders required in the calculation process as a reference. For EKF denote window length as 16, slide step as 4. Looking up table (LUT) represents exponential operation (Table 1).

Tables Icon

Table 1. Algorithm complexity analysis (single symbol)

N is the number of taps of the CMA algorithm. S-KF(NO FIR) indicates that there is no subsequent FIR filter. It is evident that the complexity of S-KF(NO FIR) is much lower than EKF. The complexity of S-KF will be analyzed below. When the number of CMA taps is 11, the complexity of the simplified Kalman filter structure is 30% lower than CMA and 40% lower than EKF. When the number of CMA taps is 21, the complexity of the simplified Kalman filter structure is 40% lower than CMA and equivalent to EKF. The more the number of CMA taps, the lower the complexity, and the limit is 50%. It can be seen that the simplified Kalman filter algorithm has lower complexity than CMA, and has the fast depolarization ability of matching Kalman filter.

4. Conclusions

In this paper, we propose a two-stage equalizer based on the simplified Kalman filter, which aims to solve the problems of matrix operation of Kalman filter complexity and nonlinear calculation that is not conducive to commercial implementation, as well as the problem of insufficient ability of traditional MIMO algorithm to solve RSOP. Through our reasonable simplification, we successfully linearized the calculation process in the Kalman filter. While having the powerful depolarization ability comparable to the Kalman filter, the complexity is reduced by 30% -40% compared with the traditional CMA, and will be lower with the increase of the number of taps. The performance of this algorithm was tested in coherent systems of 28Gbaud PDM-QPSK and 28Gbaud PDM-16QAM. The results show that under the premise of linear operation, the algorithm can handle RSOP of 20 Mrad/s; When the system suffers from additional ISI due to limited equipment bandwidth, the OSNR of the proposed algorithm is about 4 dB lower than that of the Kalman filter at the same bit error rate. The proposed algorithm provides a feasible and less complex practical application approach for dealing with fast RSOP.

Funding

China Southern Power Grid.

Disclosures

The authors declare no conflicts of interest.

Data availability

Data underlying the results presented in this paper are not publicly available at this time but may be obtained from the authors upon reasonable request.

References

1. V. Houtsma, D. V. Veen, and E. Harstead, “Recent Progress on Standardization of Next-Generation 25, 50, and 100 G EPON,” J. Lightwave Technol. 35(6), 1228–1234 (2017). [CrossRef]  

2. C. Xie, “Polarization and nonlinear impairments in fiber communication systems,” in Enabling Technologies for High Spectral-Efficiency Coherent Optical Communication Networks, X. Zhou and C. Xie, eds. (Wiley, 2016).

3. K. Zhong, X. Zhou, J. Huo, C. Yu, C. Lu, and A. Lau, “Digital signal processing for short-reach optical communications: A review of current technologies and future trends,” J. Lightwave Technol. 36(2), 377–400 (2018). [CrossRef]  

4. F. Pittalà, C. Stone, D. Clark, M. Kuschnerov, C. Xie, and A. Haddad, “Laboratory Measurements of SOP Transients due to Lightning Strikes on OPGW Cables,” in Optical Fiber Communication Conference (Optica Publishing Group, 2018), paper M4B-5.

5. P. M. Krummrich, D. Ronnenberg, W. Schairer, D. Wienold, F. Jenau, and M. Herrmann, “Demanding response time requirements on coherent receivers due to fast polarization rotations caused by lightning events,” Opt. Express 24(11), 12442–12457 (2016). [CrossRef]  

6. D. Charlton, S. Clarke, D. Doucet, M. O’Sullivan, D. L Peterson, D. Wilson, G. Wellbrock, and M. Bélanger, “Field measurements of SOP transients in OPGW, with time and location correlation to lightning strikes,” Opt. Express 25(9), 9689–9696 (2017). [CrossRef]  

7. M. Kuschnerov, F. N. Hauske, K. Piyawanno, B. Spinnler, M. S. Alfiad, A. Napoli, and B. Lank, “DSP for coherent single-carrier receivers,” J. Lightwave Technol. 27(16), 3614–3622 (2009). [CrossRef]  

8. P. Yi, D. Li, H. Song, M. Cheng, D. Liu, and L. Deng, “Experimental Investigation on Low-Complexity Adaptive Equalizer Including RSOP Tracking and Phase Recovery for 112 Gb/s PDM-QPSK Transmission System,” IEEE Photonics J. 13(2), 1–15 (2021). [CrossRef]  

9. J. Li, T. Zeng, X. Li, L. Meng, M. Luo, Z. He, and S. Yu, “Real-time fast polarization tracking based on polarization phase locking least mean square algorithm,” Opt. Express 27(16), 22116–22126 (2019). [CrossRef]  

10. W. Yi, Z. Zheng, N. Cui, X. Zhang, L. Qiu, N. Zhang, L. Xi, W. Zhang, and X. Tang, “Joint equalization scheme of ultra-fast RSOP and large PMD compensation in presence of residual chromatic dispersion,” Opt. Express 27(15), 21896–21913 (2019). [CrossRef]  

11. N. Cui, Z. Zheng, X. Zhang, W. Yi, R. Guo, W. Zhang, X. Tang, H. Xu, and L. Xi, “Joint blind equalization of CD and RSOP using a time-frequency domain Kalman filter structure in Stokes vector direct detection system,” Opt. Express 27(8), 11557–11570 (2019). [CrossRef]  

12. W. Yi, Z. Zheng, N. Cui, L. Qiu, X. Zhang, N. Zhang, W. Zhang, and L. Xi, “Joint Equalization Scheme of Ultra-fast RSOP and Large PMD in Presence of Residual Chromatic Dispersion,” in Optical Fiber Communication Conference (Optica Publishing Group, 2019), paper Th2A-51.

13. Z. Zheng, N. Cui, H. Xu, X. Zhang, W. Zhang, L. Xi, Y. Fang, and L. Li, “Window-split structured frequency domain Kalman equalization scheme for large PMD and ultra-fast RSOP in an optical coherent PDM-QPSK system,” Opt. Express 26(6), 7211–7226 (2018). [CrossRef]  

14. Y. Feng, L. Li, J. Lin, H. Xu, W. Zhang, X. Tang, L. Xi, and X. Zhang, “Joint tracking and equalization scheme for multi-polarization effects in coherent optical communication systems,” Opt. Express 24(22), 25491–25501 (2016). [CrossRef]  

15. Y. Yang, G. Cao, K. Zhong, X. Zhou, Y. Yao, A. P. Lau, and C. Lu, “Fast polarization-state tracking scheme based on radius-directed linear Kalman filter,” Opt. Express 23(15), 19673–19680 (2015). [CrossRef]  

16. Z. Yu, X. Yi, J. Zhang, M. Deng, H. Zhang, and K. Qiu, “Modified Constant Modulus Algorithm With Polarization Demultiplexing in Stokes Space in Optical Coherent Receiver,” J. Lightwave Technol. 31(19), 3203–3209 (2013). [CrossRef]  

17. N. J. Muga and A. N. Pinto, “Extended Kalman Filter vs. Geometrical Approach for Stokes Space-Based Polarization Demultiplexing,” J. Lightwave Technol. 33(23), 4826–4833 (2015). [CrossRef]  

18. W. Hou, X. Zhai, Y. Wang, Q. Lu, and X. Zhang, “A Stokes-space-rotation Based Randomly Fast RSOP Tracking Algorithm Using Extended Kalman Filter,” in 26th Optoelectronics and Communications Conference (Optica Publishing Group, 2021), paper T4B-3.

19. N. Cui, X. Zhang, Z. Zheng, H. Xu, W. Zhang, X. Tang, L. Xi, Y. Fang, and L. Li, “Two-parameter-SOP and three-parameter-RSOP fiber channels: problem and solution for polarization demultiplexing using Stokes space,” Opt. Express 26(16), 21170–21183 (2018). [CrossRef]  

Data availability

Data underlying the results presented in this paper are not publicly available at this time but may be obtained from the authors upon reasonable request.

Cited By

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

Alert me when this article is cited.


Figures (8)

Fig. 1.
Fig. 1. S-KF compared to KF algorithm process.
Fig. 2.
Fig. 2. S-KF joint FIR filter.
Fig. 3.
Fig. 3. Simulation platform.
Fig. 4.
Fig. 4. The performance of different algorithms under PDM-QPSK signals, (a) BER vs. RSOP.
Fig. 5.
Fig. 5. The performance of different algorithms under PDM-16QAM signals, (a) BER vs. RSOP, (b) OSNR vs. RSOP.
Fig. 6.
Fig. 6. The performance of different algorithms under PDM-QPSK signals, (a) BER vs. RSOP with different bandwidths, (b) OSNR vs. RSOP with different bandwidths.
Fig. 7.
Fig. 7. The performance of different algorithms under PDM-16QAM signals, (a) BER vs. RSOP with different bandwidths, (b) OSNR vs. RSOP with different bandwidths.
Fig. 8.
Fig. 8. Change of update step size of simplified Kalman filter under different RSOP of PDM-16QAM signals: (a) RSOP = 2 Mrad/s, (b) RSOP = 10 Mrad/s, (c) RSOP = 20 Mrad/s.

Tables (1)

Tables Icon

Table 1. Algorithm complexity analysis (single symbol)

Equations (35)

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

J R S O P = ( c o s κ e j α s i n κ e j β s i n κ e j β c o s κ e j α )
J 1 = ( a j b c + j d c + j d a + j b )
x k = ( a , b , c , d ) T
( E x o u t E y o u t ) = ( a j b c + j d c + j d a + j b ) ( E x i n E y i n )
{ e x = ( | E x o u t ( k ) | 2 R 2 ) e y = ( | E y o u t ( k ) | 2 R 2 )
h k = ( 0 0 ) ( e x e y )
H = ( e x a e x b e x c e x d e y a e y b e y c e y d ) = ( x a x b x c x d y a y b y c y d )
P k = P k + Q
K k = P k H T ( H P k H T + R ) 1
x k = x k + K k h k
P k = ( I K k H ) P k
J S y m b o l 1 = ( cos k 1 e j α 1 sin k 1 e j β 1 sin k 1 e j β 1 cos k 1 e j α 1 )
J S y m b o l 2 = ( cos k 2 e j α 2 sin k 2 e j β 2 sin k 2 e j β 2 cos k 2 e j α 2 )
k 2 k 1 = V R S O P T S y m b o l
M a x | ( cos k 2 cos k 1 ) | = M a x | 2 sin ( k 2 + k 1 2 ) sin ( k 2 k 1 2 ) | = 2 sin ( V R S O P T S y m b o l 2 ) 8 e 4
P k = P k + Q = d i a g ( μ a + μ μ b + μ μ c + μ μ d + μ ) = d i a g ( μ a μ b μ c μ d )
K k = P k H T ( H P k H T + R ) 1
{ x ^ k = x ^ k + K k ( y k H x ^ k ) y k = H x k + v k
K k P k H T R 1 = P k H T ( C  0 0  C ) 1 = N P k H T
( H P k H T + R ) = ( μ x + C   μ x y μ x y   μ y + C )
{ μ x = x a 2 μ a + x b 2 μ b + x c 2 μ c + x d 2 μ d μ x y = x a y a μ a + x b y b μ b + x c y c μ c + x d y d μ d   μ y = y a 2 μ a + y b 2 μ b + y c 2 μ c + y d 2 μ d
K k = P k H T ( μ x + C   μ x y μ x y   μ y + C ) 1 P k H T ( C   0 0   C ) 1 = N P k H T
K k = N ( μ a x a   μ a y a μ b x b   μ b y b μ c x c   μ c y c μ d x d   μ d y d )
( x a x b x c x d ) = ( E x i n E x o u t + E x o u t E x i n j ( E x i n E x o u t E x o u t E x i n ) E y i n E x o u t + E x o u t E y i n j ( E y i n E x o u t E x o u t E y i n ) ) = ( 2 R e ( E x i n E x o u t ) 2 I m ( E x i n E x o u t ) 2 R e ( E y i n E x o u t ) 2 I m ( E y i n E x o u t ) )
( y a y b y c y d ) = ( E y i n E y o u t + E y o u t E y i n j ( E y i n E y o u t E y o u t E y i n ) ( E x i n E y o u t + E y o u t E x i n ) j ( E x i n E y o u t E y o u t E x i n ) ) = ( 2 R e ( E y i n E y o u t ) 2 I m ( E y i n E y o u t ) 2 R e ( E x i n E y o u t ) 2 I m ( E x i n E y o u t ) )
x k = x k + K k h k
x = ( a b c d ) = ( a b c d ) + N ( μ a ( x a e x + y a e y ) μ b ( x b e x + y b e y ) μ c ( x c e x + y c e y ) μ d ( x d e x + y d e y ) )
P k = ( I K k H ) P k
( 1 μ a ( x a 2 + y a 2 )   μ a ( x a x b + y a y b )   μ a ( x a x c + y a y c )   μ a ( x a x d + y a y d ) μ b ( x a x b + y a y b )   1 μ b ( x b 2 + y b 2 )   μ b ( x c x b + y c y b )   μ b ( x d x b + y d y b ) μ c ( x a x c + y a y c )   μ c ( x c x b + y c y b )   1 μ c ( x c 2 + y c 2 )   μ c ( x c x d + y c y d ) μ d ( x a x d + y a y d ) μ d ( x d x b + y d y b )   μ d ( x c x d + y c y d )   1 μ d ( x d 2 + y d 2 ) ) P k
{ μ a = μ a ( 1 N μ a ( x a 2 + y a 2 ) ) μ b = μ b ( 1 N μ b ( x b 2 + y b 2 ) ) μ c = μ c ( 1 N μ c ( x c 2 + y c 2 ) ) μ d = μ d ( 1 N μ d ( x d 2 + y d 2 ) )
P k d i a g ( 1 N μ a ( x a 2 + y a 2 ) 1 N μ b ( x b 2 + y b 2 ) 1 N μ c ( x c 2 + y c 2 ) 1 N μ d ( x d 2 + y d 2 ) ) P k = d i a g ( μ a ( 1 N μ a ( x a 2 + y a 2 ) ) μ b ( 1 N μ b ( x b 2 + y b 2 ) ) μ c ( 1 N μ c ( x c 2 + y c 2 ) ) μ d ( 1 N μ d ( x d 2 + y d 2 ) ) )
{ E x o u t ( k ) = h x x k E x m i d ( l ) E y o u t ( k ) = h y y k E y m i d ( l )
{ e x = ( | E x o u t ( k ) | 2 R 2 ) e y = ( | E y o u t ( k ) | 2 R 2 )
{ h x x k + 1 = h x x k + μ e y ( k ) E x o u t ( k ) E x m i d ( l ) h y y k + 1 = h y y k + μ e y ( k ) E y o u t ( k ) E y m i d ( l )
x k = ( τ 1 , τ 2 , τ 3 , α , β , κ ) T
Select as filters


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