Framework for Reconstruction of Pseudo Quad Polarimetric Imagery from General Compact Polarimetry

: Pseudo quad polarimetric (quad-pol) image reconstruction from the hybrid dual-pol (or compact polarimetric (CP)) synthetic aperture radar (SAR) imagery is a category of important techniques for radar polarimetric applications. There are three key aspects concerned in the literature for the reconstruction methods, i.e., the scattering symmetric assumption, the reconstruction model, and the solving approach of the unknowns. Since CP measurements depend on the CP mode conﬁgurations, different reconstruction procedures were designed when the transmit wave varies, which means the reconstruction procedures were not uniﬁed. In this study, we propose a uniﬁed reconstruction framework for the general CP mode, which is applicable to the mode with an arbitrary transmitted ellipse wave. The uniﬁed reconstruction procedure is based on the formalized CP descriptors. The general CP symmetric scattering model-based three-component decomposition method is also employed to ﬁt the reconstruction model parameter. Finally, a least squares (LS) estimation method, which was proposed for the linear π /4 CP data, is extended for the arbitrary CP mode to estimate the solution of the system of non-linear equations. Validation is carried out based on polarimetric data sets from both RADARSAT-2 (C-band) and ALOS-2/PALSAR (L-band), to compare the performances of reconstruction models, methods, and CP modes. model constant- decomposi- outperforms the other methods the circular mode for the ocean surface. pared the iterative methods, the LS estimator for reconstruction of the urban area data combined with the decomposition-based model. the LS-based method, the linear CP mode gives lower estimates for the cross-pol term as compared with the and modes for land areas. that the constant 4 iterative method is not suitable for ocean surface reconstruction, in the circularly polarized mode. Comparatively, the constant- N method with N estimated from the decomposition has the best reconstruction circular mode to reconstruct the scattering coefﬁcient of ocean surface. The combination use of the linear CP mode and the LSND method is better than any other combinations of the CP modes and the reconstruction algorithms for the urban areas which do not possess obvious rotation angles. For the forested and tilted urban areas, reconstruction results obtained from the circular mode are better, especially when LSND is used.


Introduction
Synthetic Aperture Radar (SAR) uses electromagnetic waves to characterize target geometrical structures and dielectric properties. When the incident wave interacts with the object, the polarization state of the reradiated wave is changed. Target information is recorded in the scattered wave, and the scattering process is a function of the incident and scattered fields. Radar measurement is polarization dependent [1]. For the single-polarization transmit and single-polarization receive case, which corresponds to a single-pol system, target backscatter is characterized by a single scattering coefficient. For the single-polarization transmit and orthogonal-polarization simultaneous receive case, which corresponds to the dual-pol system, backscatter is characterized by a complex scattering vector. If we alternatively transmit orthogonal polarizations and use orthogonal polarizations to simultaneously receive backscatter, we get the scattering matrix, and this corresponds to the full-pol (FP) system. Compared to the full-pol system, the dual-pol system has advantages on imaging swath, power consumption, system complexity, and data volume [2,3], but meanwhile it cannot complete characterize the backscattering natures of scatterers. In recent years, the dual-pol system is discussed in two cases, i.e., the conventional dual-pol imaging mode (HH/HV or VH/VV) and the hybrid dual-pol imaging mode (i.e., the transmitted wave is not H or V polarized, generally known as compact polarimetry). Some studies have shown that the compact polarimetric (CP) mode performs An overview of the current related techniques is presented. The current iterationbased techniques, which were developed for either the circular or the linear CP modes, are also extended for the general CP mode. The contents of this study are summarized in Figure 1. In Section 2, we introduce the test data used for illustration. In Section 3, the formalism of the general CP descriptors in terms of the covariance matrix is introduced. In Section 4, we first summarize the reconstruction models, and then based on the formalized CP covariance matrix the iterative solving approach is extended for the general CP SAR images. In Section 5, the LS model function is proposed, and the approximation to the model parameter N is also given. In Section 6, multiple polarimetric data sets are used to show the reconstruction performances of both the iterative and the LS methods, as well as the abilities of different CP modes for revealing the quad-pol information. Finally, the paper is summarized in Section 7.
Remote Sens. 2021, 13, x FOR PEER REVIEW 3 of 20 An overview of the current related techniques is presented. The current iterationbased techniques, which were developed for either the circular or the linear CP modes, are also extended for the general CP mode. The contents of this study are summarized in Figure 1. In Section 2, we introduce the test data used for illustration. In Section 3, the formalism of the general CP descriptors in terms of the covariance matrix is introduced. In Section 4, we first summarize the reconstruction models, and then based on the formalized CP covariance matrix the iterative solving approach is extended for the general CP SAR images. In Section 5, the LS model function is proposed, and the approximation to the model parameter is also given. In Section 6, multiple polarimetric data sets are used to show the reconstruction performances of both the iterative and the LS methods, as well as the abilities of different CP modes for revealing the quad-pol information. Finally, the paper is summarized in Section 7.

Test Data
RADARSAT-2 data acquired over Fuzhou, China, and Wallerfing, Germany, as well as ALOS-2/PALSAR data acquired over San Francisco are used for verification. Corresponding Pauli-basis images are shown in Figure 2. The Fuzhou area has complex terrain types, including dispersed residential areas, mountain and hills, agriculture fields, wetlands, and sea surface. The Fuzhou data set was collected on 20 October 2013, with incidence angle ranging from 34.4º to 36.0º. This image has 6140 × 3332 pixels, and the pixel spacing is 4.73 × 4.81 m 2 . The Wallerfing area consists of large agriculture fields. Data sets were collected regularly over this area to monitor the crop growth in 2014. The data sets acquired on 28 May, 21 June, 15 July, and 8 August 2014 are included for demonstration. All the Wallerfing data sets were collected with incidence angle ranging from 40.2º to 41.6º on the descending pass. The image has 2001 × 2001 pixels, and the pixel spacing is 4.73 × 5.12 m 2 . The farmland mainly consists of 5 crop types, i.e., barley, corn, potato, sugarbeat, and wheat. The ground truth of crop types is also given in Figure 2. The San Francisco data was collected on 21 August 2018 with the center incidence angle of 33.871º. The image has 2389 × 2640 pixels, and the pixel spacing is 2.8 × 3.2 m 2 . The transversal lines and the outlined areas in Figure 2c,d will be used for analysis in experiments.

Test Data
RADARSAT-2 data acquired over Fuzhou, China, and Wallerfing, Germany, as well as ALOS-2/PALSAR data acquired over San Francisco are used for verification. Corresponding Pauli-basis images are shown in Figure 2. The Fuzhou area has complex terrain types, including dispersed residential areas, mountain and hills, agriculture fields, wetlands, and sea surface. The Fuzhou data set was collected on 20 October 2013, with incidence angle ranging from 34.4º to 36.0º. This image has 6140 × 3332 pixels, and the pixel spacing is 4.73 × 4.81 m 2 . The Wallerfing area consists of large agriculture fields. Data sets were collected regularly over this area to monitor the crop growth in 2014. The data sets acquired on 28 May, 21 June, 15 July, and 8 August 2014 are included for demonstration. All the Wallerfing data sets were collected with incidence angle ranging from 40.2º to 41.6º on the descending pass. The image has 2001 × 2001 pixels, and the pixel spacing is 4.73 × 5.12 m 2 . The farmland mainly consists of 5 crop types, i.e., barley, corn, potato, sugarbeat, and wheat. The ground truth of crop types is also given in Figure 2. The San Francisco data was collected on 21 August 2018 with the center incidence angle of 33.871º. The image has 2389 × 2640 pixels, and the pixel spacing is 2.8 × 3.2 m 2 . The transversal lines and the outlined areas in Figure 2c

Formalism of the CP Descriptors
In this section, we will briefly introduce the formalized scattering vector and covariance matrix [18].
A general form of the electromagnetic field is represented by a transverse ellipse which is described by two parameters, i.e., the ellipse orientation angle and the ellipticity angle [1], as follows.

Formalism of the CP Descriptors
In this section, we will briefly introduce the formalized scattering vector and covariance matrix [18].
A general form of the electromagnetic field is represented by a transverse ellipse which is described by two parameters, i.e., the ellipse orientation angle θ and the ellipticity angle χ [1], as follows. where E i θ, χ = 1 . The CP system measures a projection of the complex scattering matrix S onto a transmitted wave. Then for an arbitrary transmit wave E i θ, χ , the received backscattered wave E r θ, χ is For an arbitrary S, the received compact (or traditional dual) polarimetric signal is totally dependent on θ and χ. For the commonly considered linear π/4, left circular, right circular, horizontal and vertical transmitted waves, the values of θ, χ are π/4, 0 , −π/2 π/2 , π/4 , −π/2 π/2 , −π/4 , 0, 0 and π/2, 0 , respectively. Note that the circular polarization is not affected by wave orientation angles, and thus θ ∈ −π/2 π/2 . E r θ, χ in (2) is received in the linear-polarization orthogonal basis and can be easily transformed to the circular polarization basis through a unitary transformation [1,6,21]. The dual-polarization measurement is independent of the elliptical basis of radar receivers, and thus we only discuss the linear H/V-polarization received CP data. We notice that E HC contains S HH and E VC contains S VV when none of E ix θ, χ and E iy θ, χ is 0. Then, the scattering vector (2) can be formulated as where diag(·) is the diagonal matrix of a vector, and are the wave component ratios. E 1 and E 2 are the formalized elements for E r θ, χ to characterize the scattered wave. This formalism is only for the general compact polarimetric mode. It is not applicable to the conventional dual polarizations, i.e., the S HV /S VV or S HH /S VH measurements.
From the defined scattering vector k 1 , the second order product, named as the CP covariance matrix following that in full polarimetry, is constructed for the partially polarized backscattered waves, as follows.
where H denotes matrix conjugate transpose, and · denotes ensemble average. C 2 is the basis in this study to analyze the general CP matrix for the multi-polarization reconstruction.
The formalized scattering vector in (3) and the corresponding CP covariance matrix defined in (4) provide a unified method representing CP data. For a monostatic polarimetric SAR and reciprocal scatterers, the scattering matrix is symmetric, i.e., S HV = S VH . In this case, expansion of C 2 in terms of the scattering coefficients of the medium is given in (5). If scattering reflection symmetry is further assumed for the CP covariance matrix, which means the components involving products of co-polarized and cross-polarized terms are much smaller than the others and thus negligible, C 2 only contains the co-polarized and cross-polarized terms, as shown in (6). When the transmit wave is linearly polarized, e.g., Remote Sens. 2021, 13, 530 6 of 21 in the linear π/4 mode, C 2 has the same form as the wave covariance matrix [7,8,12] except that the matrix power is doubled. When the transmit wave is circularly polarized, C 2 has a different form from the wave covariance matrix as that shown in (4) in [8]. This formalism provides a unified description for the backscattered wave, which facilitates the analysis of CP imagery. Under this formalism, the system of non-linear equations for the reconstruction is easily formulated as a function of the transmit wave, which will be discussed in Sections 4 and 5.
x S HV S * VV + y * S HH S * HV x * S VV S * HV + y S HV S * HH 2Re y * S VV S * HV (5)

A General Form of Reconstruction Models
In the literature, there are mainly five reconstruction models. The reconstruction models in [7][8][9][10][11] have the same form as follows.
where ρ is the co-polarized correlation defined by ρ = S HH S * VV |S HH | 2 |S VV | 2 (8) and N is the model parameter derived from either theoretical analysis or empirical curve fitting. The Souyris et al. model [7] is derived from a pseudo-deterministic relationship by linking ρ and the cross-and co-polarization ratio. By assuming the backscattered signal is either fully polarized or fully depolarized, N is obtained as 4. The Nord et al. model [8] is approximated from a mathematical inequality, and the equality is asserted by a relatively small difference of the two sides of the inequality. Both Collins's [9,10] and Li's [11] models are empirical models but with different regression equations, where Collins et al. used the exponential model to fit N by considering the incidence angle, and Li et al. fitted N as a function of the polarization ratio. A backscattered process is a mixture of several elemental scattering processes. By assuming that the backscattered energy is a sum of the surface scattering power P s , double-bounce scattering power P d , and volume scattering power P v , we established a decomposition-based power-weighted reconstruction model [12], as shown in (9), which use Re( S HH S * VV ) as an indicator to discriminate between the surface and double-bounce reflections (refer to [12] for more details).
where sg = sign(Re( S HH S * VV )) is a signum function, and N = 4(2P d + P v )/P v ). This power-weighted reconstruction model has an extra parameter sg which is either 1 or −1. In the following sections, we use (9) as a general formulation of the reconstruction models, and the model parameters are summarized in Table 1. The model parameter N is updated during the iterative reconstruction procedure in the studies [8,11,12], which are referred to as the variable N methods. While parameter N is a constant in the methods [7,9], which are referred to as the constant N methods.

Unified Experssion of the Model Parameters
In previous studies, the reason that different CP modes need different reconstruction procedures arises from the fact that different equations are necessary to calculate the copolarized correlation coefficient ρ. See Equations (7) and (11) in [8] for instance. However, with the unified CP covariance matrix defined in (6), the co-polarized coherence has a general expression, given by x and y are parameters depending on the transmit wave. With this expression, the reconstruction model in (9) is applicable to all CP modes. Next, we consider the parameter sg, which only exists in the decomposition-based reconstruction model [12]. Our previous studies [12,14] tested that for the linear π/4 mode, the total agreement between sign(Re( C 12 )) and sign(Re( S HH S * VV )) is quite good. For the general case, from the off-diagonal term of (6) it is observed that C 12 = S HH S * VV + |S HV | 2 xy * , which indicates the sign of S HH S * VV is only affected by |S HV | 2 e j2δ , where δ = angle E iy E * ix is the relative phase of the transmit wave, which is a known number. For the area dominated by a single scattering mechanism, |Re( S HH S * VV )| is usually larger than |S HV | 2 . Therefore, for the general reconstruction algorithm, we still use Re( C 12 ) to approximate the decision whether backscatter is dominated by surface or double-bounce scattering [12][13][14]. Let sg c = sign(Re( C 12 )). By using the test data shown in Section 2, variation of the total agreement between sign(Re( C 12 )) and sign(Re( S HH S * VV )) with the transmit wave parameter δ is shown in Figure 3. It shows that for the C-band data when the relative phase δ is ±π/4 or ±3π/4, the total agreement is as high as 99.92% (99.96% for the L-band data), and when the relative phase δ is ±π/2 (e.g., one possible mode is the circular modes), the total agreement is lowest but still with a good percentage of 88%. Areas dominated by double-bounce scattering or double-bounce scattering taking a relatively high amount, such as the urban areas, would have a lower overall agreement. When the transmit wave is linearly polarized (i.e., δ = 0 or ±π), the overall agreement is 93.65% for C-band data and 93.83% for L-band data. The general form of the reconstruction models with model parameters derived from (6) constructs the unified reconstruction model. of 88%. Areas dominated by double-bounce scattering or double-bounce scattering taking a relatively high amount, such as the urban areas, would have a lower overall agreement. When the transmit wave is linearly polarized (i.e., = 0 or ±π), the overall agreement is 93.65% for C-band data and 93.83% for L-band data. The general form of the reconstruction models with model parameters derived from (6) constructs the unified reconstruction model.

The Decomposition-Based Variable N
The three-component model-based decomposition [16] is written as

The Decomposition-Based Variable N
The three-component model-based decomposition [16] is written as where S VV T ; f s , f d and f v are the decomposed parameters to be determined corresponding to the surface, double-bounce, and volume scattering models, which are given by where β and α are model parameters with arg(β) ≈ 0 for surface scattering and arg(α) ≈ ±π for double-bounce scattering. In Freeman and Durden's 3-component decomposition, b is set to 1/3. However, for the purpose of the pseudo quad-pol image reconstruction, we set b = 0 [12,14]. According to the general CP expression in (6), scattering models in (12) are synthesized as follows. [ Then the general CP covariance matrix can be expanded as where f sc , f dc , and f vc are expansion coefficients. This decomposition should be performed under the constraint that the ratio of the total backscattered powers of the FP and the CP data is known. Otherwise, the CP decomposed powers cannot be related to the FP decomposed powers. Only when |x| 2 = |y| 2 , which corresponds to the CP modes with either θ = ±π/4 or χ = ±π/4, the decomposed powers can be expressed as a function of |S HV | 2 . In this case, the total powers of the FP covariance matrix and the formalized CP covariance matrix are equal by assuming reflection symmetry. We assume that the CP decomposed volume scattering power is equal to that of FP, i.e., f vc = f v = 2 |S HV | 2 (refer to (3) in [14]). Then, (14) can be solved following the Freeman Remote Sens. 2021, 13, 530 9 of 21 and Durden procedure. When Re( C 12 ) is positive, surface scattering is dominant and we let α = −1. When Re( C 12 ) is negative, double-bounce scattering is dominant and we let β = 1. It is easily derived that the decomposed powers, represented as a function of |S HV | 2 , are . P sc , P dc , and P vc are the surface, double-bounce, and volume scattering powers of CP data, respectively. Note that this solution is only valid for the case when |x| 2 = |y| 2 , i.e., |x| = |y| = 1. Then N can be updated by (16), which is embedded in the iteration procedure.

The Iterative Solving Approach
By assuming reflection symmetry, we can obtain 3 equations from (6) for each pixel. However, we have four unknows, i.e., |S HH | 2 , |S VV | 2 , |S HV | 2 , and S HH S * VV . When the reconstruction model is further considered as the fourth equation, |S HV | 2 can be solved by an iterative approach, as the procedures introduced in [7][8][9][10][11][12][13]. These studies were only discussed for the typical linear-or circular-transmit CP mode.
Based on the unified CP descriptor and the reconstruction model in (9), the iterative approach for solving |S HV | 2 for the mode with |x| = |y| is summarized as below.
Step 0: Initialization Step 1: Iteration where i = 0, 1, 2 . . . is the iteration number. For the methods in [7][8][9][10][11], sg c = sg = 1, as listed in Table 1. For the decomposition-based reconstruction model [12], sg c = sign(Re( C 12 )). N varies and updated according to different models. For the Nord et al. [8] and the Li et al. [11] models, N can be updated easily through the constructed FP covariance matrix at each iteration. For the decomposition-based method in [12], the model parameter N is updated by the three-component decomposed-powers, which are a function of the cross-pol term |S HV | 2 (i) . The updating method has been introduced in Section 4.3. Given an estimated value for |S HV | 2 c , the pseudo quad-pol data is constructed as Remote Sens. 2021, 13, 530 10 of 21 The iteration usually converges after several steps. However, the estimated crosspolarized term is often overestimated [7][8][9]14]. In the next section, we propose a least squares error estimation framework to optimize |S HV | 2 for all CP modes.

Least Squares Objective Function
The LS estimation is performed in a local window, and can be applied to both multilook and single-look data. However, if the multi-look data is used for the LS estimation, the reconstructed pseudo quad-pol image may be over-filtered. Therefore, we use the single-look complex data to fit the reconstruction model. Assuming that an optimal crosspolarized term |S HV | 2 exists, then for each pixel in the window, the system of non-linear equations is where ρ is the local co-polarized correlation coefficient, as defined in (10). Substituting the first 3 equations in (20) to the last one, a least squares error objective function is obtained for |S HV | 2 , as follows.
where i is the i-th pixel in the local window, n is the window pixel number, and |S HV | 2 c is the LS estimate for the problem. The local window can be either square or non-square, if an edge-aligned non-square window as that in [22] is used, image texture information will be better preserved. The LS solution |S HV | 2 c is for the center pixel of the window. After optimization, the reconstructed FP covariance matrix is obtained by (19).
In the LS model function, since pixels in a local window are used to find the best fitting parameter, N should also be an average estimate from the local window. In this study, we do not consider the empirical models, because they are related to specific observation scenarios. The variable N method is not applicable to the LS estimation because N needs to be a fixed number. The Souyris et al. model [7] can be integrated into the LS objective function since the model parameters are constant, with parameter settings as sg c = 1 and N = 4. For the decomposition-based model, sg c = sign(Re( C 12 )) and N = 4(2P dc + P vc )/P vc , where P dc and P vc are the decomposed powers from the CP data. In [14], we proposed a new method to give a priori estimate for N based on the linear π/4 mode. In [19], the decomposition was extended to a general case. Therefore, in this paper, the decomposed powers are used to estimate N.

The Decomposition-Based Approximation to the Constant N
N is a parameter indicating the ratio between double-bounce and volume scattering components. The iterative approach updates N based on the intermediate result of |S HV | 2 . For the LS estimation, N should be a constant. In [20], a CP decomposition was proposed for the circular CP mode, where the single scattering mechanism is modeled by the symmetric scattering type angle α s . In [19], we generalize this decomposition to the general CP case, but leave one degree of freedom, i.e., F p in (14) in [19], for the volume scattering model, and take it as an argument for the decomposed powers. In this study, when approximating N, let F p = 1. The reasons are as follows. (1) The decomposition-based reconstruction model was proposed with F p = 1 [12]. (2) The volume scattering model has the same degree of polarization of 1/3 for all the CP modes which are with balanced transmit wave channel amplitudes (i.e., |E ix | = E iy ). Suppose the decomposed powers from the α s -based decomposition are P dc and P vc , respectively, for the double-bounce and volume scattering components. Then an approximation to N is obtained by (16).

Solution Constraint
The cross-polarized term |S HV | 2 should be optimized under constraints. In [14], we apply 2 conditions for |S HV | 2 , which can be easily extended to the general CP case.
(1) The polarization intensity should be larger than 0, and in general the crosspolarized intensity is smaller than the co-polarized intensities. Based on (6), it follows that |S HV | 2 should be estimated in the interval 0 min We tested this on several RADARSAT-2 and ALOS-2/PALSAR data sets for the linear π/4 mode, and on average 98.8% pixels are within this interval.
(2) In the model-based target decomposition, the decomposed volume scattering power P v = 6 |S HV | 2 should be smaller than the total backscattered energy span = |S HH | 2 + 2 |S HV | 2 + |S VV | 2 , from which it follows that 4 |S HV | 2 < |S HH | 2 + |S VV | 2 . Substituting |S HH | 2 and |S VV | 2 by the CP covariance elements C 11 and C 22 in (6), we obtain |S HV | 2 < C 11 + C 22 4 + |x| 2 + |y| 2 By using the same data sets as those used in the last test for the first condition, 98.2% pixels satisfy this condition.
In summary, the cross-polarized term should be estimated in the interval 0 min We use the Fuzhou data for illustration, percentage of pixels that satisfy this constraint is shown in Figure 4 for various CP modes. It shows that when the transmit wave is balanced in amplitude, i.e., θ = ±π/4 or χ = ±π/4, the linear CP mode has the lowest agreement and the circular mode the highest, but the difference is very small. When the transmit wave is imbalanced in amplitude, i.e., θ = ±π/4 and χ = ±π/4, a decreased satisfactory is observed, and the percentage decreases with channel imbalance increasing. When the channel imbalance increases, the neglected terms in (5) by assuming reflection symmetry gradually becomes comparable to the cross-pol term, thus resulting a decreased percentage. However, in any case, the pixels falling in the interval take part in more than 96.5%. transmit wave is imbalanced in amplitude, i.e., ≠ ±π/4 and ≠ ±π/4, a decreased satisfactory is observed, and the percentage decreases with channel imbalance increasing. When the channel imbalance increases, the neglected terms in (5) by assuming reflection symmetry gradually becomes comparable to the cross-pol term, thus resulting a decreased percentage. However, in any case, the pixels falling in the interval take part in more than 96.5%.

Experiments
The LS-based methods are applied to the single-look complex data, and the estimation is performed within a 7 × 7 square window. The iterative methods are applied to the multi-look data where a 7 × 7 sliding window is used for speckle reduction.
is set to 1 for the symmetric scattering type-based decomposition.
The pseudo quad-pol imagery reconstruction can be implemented via the combination of the reconstruction models and the LS-based or the iteration-based solving approach. The reconstruction model can be with either constant , which is a pre-defined parameter, or variable , which is updated during the iteration procedure. In Section 4.1, we formulated the reconstruction model and presented both the empirical and theoretical model parameters in Table 1. The empirical model is closely related to the observation scenarios and can only be applied to areas with a single terrain type, so we do not include the empirical model in the experiments. The Nord et al. variable model needs initial values for the iterative approach and is sensitive to the initialization [9], which cannot be embedded in the LS estimator. Thus, we only consider two reconstruction models, i.e., the = 4 model [7] and the decomposition-based model [12]. These two models can be combined with both the LS estimator and the iterative approach. For the iterative approach,

Experiments
The LS-based methods are applied to the single-look complex data, and the estimation is performed within a 7 × 7 square window. The iterative methods are applied to the multi-look data where a 7 × 7 sliding window is used for speckle reduction. F p is set to 1 for the symmetric scattering type-based decomposition.
The pseudo quad-pol imagery reconstruction can be implemented via the combination of the reconstruction models and the LS-based or the iteration-based solving approach. The reconstruction model can be with either constant N, which is a pre-defined parameter, or variable N, which is updated during the iteration procedure. In Section 4.1, we formulated the reconstruction model and presented both the empirical and theoretical model parameters in Table 1. The empirical model is closely related to the observation scenarios and can only be applied to areas with a single terrain type, so we do not include the empirical model in the experiments. The Nord et al. variable N model needs initial values for the iterative approach and is sensitive to the initialization [9], which cannot be embedded in the LS estimator. Thus, we only consider two reconstruction models, i.e., the N = 4 model [7] and the decomposition-based model [12]. These two models can be combined with both the LS estimator and the iterative approach. For the iterative approach, different updating strategies can be applied to the model parameter N, either constant or variable with different initial values. Hence, from the explanation above, we consider 6 reconstruction algorithms, that is, the LS estimators with Souyris's and Yin' model parameters, the iterative approaches with constant Souyris's and Yin' model parameters, and the iterative approaches with Yin's model in which N is variable and initialized with N = 4 as well as initialized by the CP α s -based decomposition method. The above algorithms are denoted in turn as LSN4, LSND, IterNC4, IterNCD, IterNV4, and IterNVD, respectively.
Experiments are conducted on the following aspects: (1) performances of the iterative and LS-based methods, (2) reconstruction accuracies of the reconstruction models, (3) reconstruction accuracies under different CP modes, and (4) comparison of the performances of CP data, pseudo quad-pol data, and FP data for multi-temporal agriculture field classification. Three CP modes, i.e., the linear π/4 mode (θ = π/4, χ = 0), an elliptical mode (θ = π/4, χ = π/8), and the left circular mode (θ = π/4, χ = π/4), are used for demonstration. For these CP modes, the transmit waves are balanced in amplitude, i.e., |x| = |y|, so the variable-N algorithm is applicable. The linear π/4 and circular modes are commonly considered imaging configurations. They are included so as to get comparable results with those found in the open literature.
Reconstructed results for the line segments in Figure 2 are given in Figure 5. It shows that in general the results from the decomposition-based model is better than those from the N = 4 model, no matter which solving approaches is used. For the variable N methods, we again verify that the initialization affects the reconstructed results, especially for the circular CP mode, as can be observed from the last two rows in Figure 5. Initial values given by the α s -based decomposition generates better results than the initial guess with N = 4. For the iterative methods with constant N, it shows that the N = 4 model under the circular mode does not performs well for the ocean surface, as can be observed in the IterNC4 plots that the blue profiles for the ocean surface has larger variations. Results show that the iterative method with constant-N estimated by the α s -based decomposition outperforms the other methods under the circular mode for the ocean surface. Compared to the iterative methods, the LS estimator is superior for reconstruction of the urban area data when combined with the decomposition-based model. In the LS-based method, the linear CP mode gives lower estimates for the cross-pol term as compared with the elliptical and circular CP modes for land areas. From Figure 5, it shows that the constant N = 4 iterative method is not suitable for ocean surface reconstruction, especially in the circularly polarized mode. Comparatively, the constant-N method with N estimated from the α s -based decomposition (IterNCD) has the best reconstruction accuracy for ocean surface under the circular mode. While for land areas, the LS-based method LSND has the best result. However, for the urban area without obvious rotation, all methods tend to overestimate the cross-pol term. Further, results in Figure 5 also show that the decomposition-based model fit the urban areas better in both iterative and LS-based solving approaches, which is because the typical reflection asymmetry model, i.e., the helix scattering models, also satisfies the decomposition-based reconstruction model. However, the reflection asymmetric models do not agree with the N = 4 model. Evaluation in terms of root mean square errors (RMSE) and percentages of the pixels that deviate from the real values by 5% of the total range are given in Table 2. |S HV | 2 directly affects the reconstruction performances and CPD is associated with the identification of scattering mechanisms. RMSE is calculated for both |S HV | 2 and CPD, and the pixel percentage of deviation is calculated for CPD. In Table 2, for each assessment index, the first two methods that are with best estimation results are highlighted. Reconstruction are performed under the aforementioned 3 CP modes for the test data sets. In total, each method is quantified 27 times. Results show that for the LSND, LSN4, IterNCD, IterNC4, IterNVD, and IterNV4 methods, the frequencies that those methods perform best are 21/27, 1/27, 10/27, 8/27, 13/27, and 2/27, respectively. The LS estimator with the decompositionbased model can provide the best overall results, and the variable-N method initialized by parameters from the α s -based decomposition is in the second position.
In the LSND method, for the linear π/4 mode, the elliptical model with (θ, χ) = (π/4, π/8), and the circular mode, on average 0.98%, 0.6% and 0% pixels, respectively, cannot find minimums in the solution constraint. With the LSN4 method, for the 3 CP modes, 1.17%, 1.14%, and 1.83% pixels cannot find minimums within the constraint. This indicates that the decomposition-based model is more accurate to fit the CP data under different modes. face under the circular mode. While for land areas, the LS-based method LSND has the best result. However, for the urban area without obvious rotation, all methods tend to overestimate the cross-pol term. Further, results in Figure 5 also show that the decomposition-based model fit the urban areas better in both iterative and LS-based solving approaches, which is because the typical reflection asymmetry model, i.e., the helix scattering models, also satisfies the decomposition-based reconstruction model. However, the reflection asymmetric models do not agree with the = 4 model.  Evaluation in terms of root mean square errors (RMSE) and percentages of the pixels that deviate from the real values by 5% of the total range are given in Table 2. 〈| | 〉 directly affects the reconstruction performances and CPD is associated with the identification of scattering mechanisms. RMSE is calculated for both 〈| | 〉 and CPD, and the pixel percentage of deviation is calculated for CPD. In Table 2 Table 2. Quantitative evaluation for log |S HV | 2 and CPD by using the root mean square error (RMSE) and the percentage of pixels deviating from the real values by 5% of the total range. Three CP modes illustrated are (θ = π/4, χ = 0), (θ = π/4, χ = π/8), and (θ = π/4, χ = π/4). The reconstructed accuracy of CPD is another important factor to evaluate the reconstruction performance. On average, the LSND method performs best in phase reconstruction. By using the L-band ALOS-2/PALSAR data for illustration, CPDs estimated by LSND for the three CP modes are shown in Figure 6a. It shows that under the linear π/4 mode, most pixels distribute along the diagonal line, indicating a superior reconstruction result. The reconstructed |S HV | 2 is evaluated in Figure 6b in terms of the relative error. Similar analysis results can be found in [14] (see Figures 5 and 6) for the C-band RADARSAT-2 data over the same test site for the Linear CP mode. Figure 6b shows that for the urban area without obvious rotation, the linear π/4 mode has smaller relative errors compared to the other two modes, but still the relative error of city areas is larger than those of ocean surface and forested areas. For the rotated urban area, the linear π/4 mode has the larger relative errors and in contrast the circular mode performs best. Since for this test set, the rotated urban area and the forested areas are characterized by a similar cross-pol and co-pol ratio, i.e., ratio between the cross-pol intensity and the sum of co-pol intensities, it is expected that the circular mode would also perform better for the forested areas. Relative errors of the reconstructed log |S HV | 2 for the areas outlined in Figure 2c are listed in Table 3. The quantitative result is consistent with the explanation for the results shown in Figure 6b.

RS-2 WF data
= π/8  In the LSND method, for the linear π/4 mode, the elliptical model with ( , ) = (π/4, π/8), and the circular mode, on average 0.98%, 0.6% and 0% pixels, respectively, cannot find minimums in the solution constraint. With the LSN4 method, for the 3 CP modes, 1.17%, 1.14%, and 1.83% pixels cannot find minimums within the constraint. This indicates that the decomposition-based model is more accurate to fit the CP data under different modes.
The reconstructed accuracy of CPD is another important factor to evaluate the reconstruction performance. On average, the LSND method performs best in phase reconstruction. By using the L-band ALOS-2/PALSAR data for illustration, CPDs estimated by LSND for the three CP modes are shown in Figure 6a. It shows that under the linear π/4 mode, most pixels distribute along the diagonal line, indicating a superior reconstruction result. The reconstructed 〈| | 〉 is evaluated in Figure 6b in terms of the relative error. Similar analysis results can be found in [14] (see Figures 5 and 6) for the C-band RADARSAT-2 data over the same test site for the Linear CP mode. Figure 6b shows that for the urban area without obvious rotation, the linear π/4 mode has smaller relative errors compared to the other two modes, but still the relative error of city areas is larger than those of ocean surface and forested areas. For the rotated urban area, the linear π/4 mode has the larger relative errors and in contrast the circular mode performs best. Since for this test set, the rotated urban area and the forested areas are characterized by a similar cross-pol and copol ratio, i.e., ratio between the cross-pol intensity and the sum of co-pol intensities, it is expected that the circular mode would also perform better for the forested areas. Relative errors of the reconstructed log〈| | 〉 for the areas outlined in Figure 2c are listed in Table 3. The quantitative result is consistent with the explanation for the results shown in Figure 6b. Results in Figures 5 and 6b, and Table 3 show that a same method performs differently for different terrain types when the CP mode varies. Since from Table 2, it is shown that apart from the LSND method, the IterNCD and IterNVD methods have relatively better overall reconstruction accuracies, in Table 4 we further show the relative errors of 〈| | 〉 for those methods. In combination with the profiles in Figure 5, it confirms that the IterNCD method is favorable for the circular mode to reconstruct the scattering coefficient of ocean surface. The combination use of the linear CP mode and the LSND method is better than any other combinations of the CP modes and the reconstruction algorithms . ALOS-2/PALSAR San Francisco data. Table 3. Relative errors of the reconstructed log |S HV | 2 for the outlined areas in Figure 2c by using the LSND method. The CP modes are with θ = π/4. Results in Figures 5 and 6b, and Table 3 show that a same method performs differently for different terrain types when the CP mode varies. Since from Table 2, it is shown that apart from the LSND method, the IterNCD and IterNVD methods have relatively better overall reconstruction accuracies, in Table 4 we further show the relative errors of |S HV | 2 for those methods. In combination with the profiles in Figure 5, it confirms that the IterNCD method is favorable for the circular mode to reconstruct the scattering coefficient of ocean surface. The combination use of the linear CP mode and the LSND method is better than any other combinations of the CP modes and the reconstruction algorithms for the urban areas which do not possess obvious rotation angles. For the forested and tilted urban areas, reconstruction results obtained from the circular mode are better, especially when LSND is used. Table 4. Relative errors of the reconstructed log |S HV | 2 for the outlined areas in Figure 2c by using (a) the IterNCD method, and (b) the IterNVD method. The CP modes are with θ = π/4.  Table 2 shows that for all methods, the deviation percentage which evaluates the agreement of reconstructed and real CPDs is smallest in the linear π/4 mode. In combination with the scatter plots in Figure 6a, we found that the π/4 mode outperforms the other two CP modes in the reconstruction of the phase term. RMSE of CPD is larger for the π/4 mode is due to the period of 2π radians phase. When the real CPD is close to ±π, a small perturbation on the estimated |S HV | 2 will lead to dramatic changes of ±2π in the estimated CPD (e.g., from −π to π or from π to −π), as can be observed in Figure 6a that in the π/4 mode a certain amount of pixels lies in the corners of (−π, π) and (π, −π). This results in larger RMSEs. For the other two test sets in Figure 2, distributions of the CPD estimated by LSND under the 3 CP modes are quite similar to those in Figure 6a. The omitted terms in (5) by assuming reflection symmetry has different impacts on CPD when the CP mode varies. In Figure 7, assuming that |S HV | 2 is perfectly reconstructed, we show variation of the deviation percentages of the reconstructed CPD with the varying CP modes. It shows that for all test data sets, when the transmit wave channel amplitude is balanced, the estimated CPD is least affected by the assumption of reflection symmetry in the linear mode. The case that the channel amplitudes of the transmit wave are imbalanced is also tested. We found that the more the CP mode is linearly polarized, the closer the reconstructed CPD is to the real value. This implies that if applications are based on the pseudo quad-pol images and the algorithms applied subsequently are CPD-based, results obtained from the linear π/4 mode may be expected to be closer to that of the FP mode. For example, if CPD is used to discriminate surface and double-bounce scattering such as in the application of Freeman-Durden's decomposition, the overall agreements in |CPD|>π/2 between the FP data (ALOS-2/PALSAR SF data) and the pseudo-FP data reconstructed by LSND are 88.75% and 84.72%, respectively, for the linear π/4 and circular modes, and by IterNVD the agreements are 89.35% and 81.34%, respectively.

CP
Another application is carried out for example, which is the crop type classification by using the RS-2 Wallerfing data. The iterative Wishart classifier [23] is used to perform the classification, with 5% pixels randomly selected as the training samples based on ground data. For each data set, same training samples are used for all the subsequent experimental implementation. Although several studies [3,4,24] have been conducted on the comparison of the performances of CP, pseudo quad-pol, and FP data for terrain type classification, there is rare study giving a clear illustration on the variation of the overall classification accuracies under different circumstances. In Figure 8, by using the Wallerfing data acquired by RADARSAT-2 on May 28, 2014, we tested the classification performances of the reconstructed data, the original 2 × 2 CP covariance matrix data, and the original 3 × 3 FP covariance matrix data. The reconstructed |S HV | 2 is assumed to be rX min with r ∈ 0.1 1 , where X min is the upper bound for the cross-pol term in estimation. By this means, the pseudo FP data is accordingly reconstructed. It shows that compared to the original CP data, the reconstructed data with the accurately estimated |S HV | 2 greatly improve the classification accuracy, especially for the linear CP mode, in which case the classification accuracy is improved by 5.3%. For the circular mode, when the estimated |S HV | 2 is too large, the classification accuracy will deteriorate dramatically.
In Figure 8, we did not give the classification rates in the cases of |S HV | 2 c = 0.9X min and |S HV | 2 c = X min for the circular mode, because in both cases the term |S HH − S VV | 2 /2 of the reconstructed coherency matrix is close to 0 or even negative for a certain amount of pixels, for instance, 15% pixels in the r = 0.9 case and 23% in the r = 1 case. However, if the estimation is accurate, Figure 8 shows that in general the circular mode is better for agriculture field classification. Another application is carried out for example, which is the crop type classification by using the RS-2 Wallerfing data. The iterative Wishart classifier [23] is used to perform the classification, with 5% pixels randomly selected as the training samples based on ground data. For each data set, same training samples are used for all the subsequent experimental implementation. Although several studies [3,4,24] have been conducted on the comparison of the performances of CP, pseudo quad-pol, and FP data for terrain type classification, there is rare study giving a clear illustration on the variation of the overall classification accuracies under different circumstances. In Figure 8, by using the Wallerfing data acquired by RADARSAT-2 on May 28, 2014, we tested the classification performances of the reconstructed data, the original 2 × 2 CP covariance matrix data, and the original 3 × 3 FP covariance matrix data. The reconstructed 〈| | 〉 is assumed to be X with ∈ 0.1 1 , where X is the upper bound for the cross-pol term in estimation. By this means, the pseudo FP data is accordingly reconstructed. It shows that compared to the original CP data, the reconstructed data with the accurately estimated 〈| | 〉 greatly improve the classification accuracy, especially for the linear CP mode, in which case the classification accuracy is improved by 5.3%. For the circular mode, when the estimated 〈| | 〉 is too large, the classification accuracy will deteriorate dramatically. In Figure 8, we did not give the classification rates in the cases of 〈| | 〉 = 0.9X and 〈| | 〉 = X for the circular mode, because in both cases the term 〈| − | 〉/ 2 of the reconstructed coherency matrix is close to 0 or even negative for a certain amount of pixels, for instance, 15% pixels in the = 0.9 case and 23% in the = 1 case. However, if the estimation is accurate, Figure 8 shows that in general the circular mode is better for agriculture field classification.
Next we consider the performances of the pseudo quad-pol data sets reconstructed via the test methods for multi-phase crop classification. The overall classification accuracies of the 6 methods, and those of the original CP as well as the FP data are shown in Figure 9. We observe that the circular mode performs best, which is in accordance with the analysis in Figure 8. However, for the multi-phase data sets, due to different phenology periods and disturbance of environmental changes, the best method evaluated in  Next we consider the performances of the pseudo quad-pol data sets reconstructed via the test methods for multi-phase crop classification. The overall classification accuracies of the 6 methods, and those of the original CP as well as the FP data are shown in Figure 9.
We observe that the circular mode performs best, which is in accordance with the analysis in Figure 8. However, for the multi-phase data sets, due to different phenology periods and disturbance of environmental changes, the best method evaluated in terms of overall classification accuracy varies. However, in general, data reconstructed by the iterative approach with the decomposition-based model provides comparatively better and stable classification accuracies if the Wishart classifier is applied. It also shows that classification carried out on the reconstructed data is better than that on the original CP data.

Summary and Discussion
In this study, based on the formalized CP descriptors which was one of the recent developed works, we unified the current pseudo quad-pol imagery reconstruction models into two frameworks, which are the LS error model and the iterative approach, to estimate the FP scattering coefficients. By inducing this formalism, the system of nonlinear equations is parameterized by the transmit wave coefficients, and the most important is the copolarized coherence has a unified form for the general CP mode. By this means, the two unknown solving approaches can be combined with any reconstruction models (either empirical or theoretical models) for an arbitrary CP mode. The decomposition-based

Summary and Discussion
In this study, based on the formalized CP descriptors which was one of the recent developed works, we unified the current pseudo quad-pol imagery reconstruction models into two frameworks, which are the LS error model and the iterative approach, to estimate the FP scattering coefficients. By inducing this formalism, the system of nonlinear equations is parameterized by the transmit wave coefficients, and the most important is the copolarized coherence has a unified form for the general CP mode. By this means, the two unknown solving approaches can be combined with any reconstruction models (either empirical or theoretical N models) for an arbitrary CP mode. The decomposition-based variable N method is also extended for the CP mode which is with balanced transmitted wave channel amplitudes.
Experiments were carried out on 6 data sets to demonstrate the performances of reconstruction models and solving approaches. On average, the decomposition-based model is more adaptive and provides better results than the N = 4 model. Results showed that by using a same reconstruction method the linear mode has better reconstruction accuracy in the cross-pol term than the circular mode for urban areas. While the circular mode always performs better for areas not dominated by double-bounce scattering, especially for ocean surface. If evaluated in terms of CPD, the linear π/4 mode always outperforms the circular Remote Sens. 2021, 13, 530 20 of 21 mode. We also carried out same experiments on several other data sets in addition to the data sets presented in this paper, and similar conclusions were obtained.
Most reconstruction algorithms, including the methods summarized and extended in this study, were proposed based on reflection symmetry. If backscatter of scatterers does not follow the reflection symmetric assumption, then the reconstructed results will have a larger estimation error. This is caused by the omitted terms, which are assumed to be 0 by reflection symmetry but actually are comparable to the cross-pol term, thus leading to larger relative errors. The reconstruction model parameter N is the key factor for good reconstruction. If N is accurate, both the LS-based and the constant-N iterative methods can achieve a perfect reconstruction. At this point the variable-N method is no longer preferable. However, except the empirical N, there are only 2 theoretical methods that can give a prior estimation for N, i.e., the N = 4 and the N from the α s -based decomposition. The N = 4 model applies well for areas dominated by volume scattering, but this application situation is very restricted. N estimated from the α s -based decomposition is more adaptive to the observation scenarios. However, since the CP mode only measures partial information of the backscattering space, the difference of the decomposed powers from CP and FP data is inevitable. Improvement on the accuracy of N still needs further studies, e.g., by developing semi-empirical models.  Data Availability Statement: Data is available upon request for scientific purpose.