Polarimetric SAR Speckle Reduction by Hybrid Iterative Filtering

Speckle filtering in synthetic aperture radar polarimetry (PolSAR) is essential for the extraction of significant information. In this study, the authors introduced a hybrid iterative filtering scheme. The proposed iterative filter is initialized by a polarimetric filter ensuring a high speckle reduction level. Then, to enhance the spatial details, the iterative filter is applied for few iterations. The key parameter <inline-formula> <tex-math notation="LaTeX">$b_{k}$ </tex-math></inline-formula> which traduced the variability of the pixels has been studied. An expression that takes into account the following important numerical considerations is proposed. 1) The variability of the pixel is measured by coefficient variation (CV) rather than the variance. 2) The statistics are computed using non-local neighborhood instead of local square one. 3) A new normalizing function (i. e. hyperbolic tangent) is implemented. 4) The dynamic of range of <inline-formula> <tex-math notation="LaTeX">$b_{k}$ </tex-math></inline-formula> is increased by multiplying the CVs of the filtered and the originals images. Comparison with various state of the art PolSAR filters demonstrated the effectiveness of the proposed iterative method. Simulated, one-look and multilook real PolSAR data were used for validation.


I. INTRODUCTION
The synthetic aperture radar polarimetry (PolSAR) plays an important role in the analysis and characterization of the Earth's surface. However, due to the coherent record of echoes, PolSAR data are affected by speckle noise which is one of the major problems of the SAR imagery that should be filtered correctly.
The multi-looking process is a widely applied speckle filtering technique [1]. However, its main weakness lies in the spatial resolution degradation. Based on the multiplicative noise model and the minimum mean square error (MMSE), the filtered Mueller matrix was estimated [2]. The model-based PolSAR (MBPolSAR) filter [3] processes the covariance matrix according to the multiplicative-additive speckle noise model [4]. To preserve spatial details and polarimetric information, research in PolSAR image filtering has mainly focused on the selection of pixels having similar The associate editor coordinating the review of this manuscript and approving it for publication was Jenny Mahoney. scattering properties. The refined Lee filter [5] selects the most homogeneous window from eight edge aligned windows. The intensity-driven adaptive-neighborhood (IDAN) filter uses a region growing technique to produce an adaptive neighborhood [6]. The scattering model-based filter (SMBF) selects homogeneous pixels of the same scattering mechanisms [7]. In the polarimetric Lee sigma filter, homogeneous pixels are selected within refined sigma ranges [8].
The dimension of the homogenous window is increased by selecting non-local (NL) pixels. The NL filtering includes the pretest filter [9], the stochastic based filter [10], the NL-SAR filter [11], the NL-discriminative similarity measure NL-DSM filter [12], the NL-distributed Lee filter [13], the patch-based NLM adaptive speckle filter based on constant false alarm rate (CFAR) edge detector [14]. The capability of polarimetric non-local filters for speckle reduction in PolSAR measurements is investigated in [15].
Variational methods based on total variation (TV) regularization, introduced in [16], constitutes a new trend for PolSAR speckle reduction. A TV-based variational model, VOLUME 8, 2020 This work is licensed under a Creative Commons Attribution 4.0 License. For more information, see https://creativecommons.org/licenses/by/4.0/ called WisTV, was introduced in [17]. In [18], WisNLTV filter hybridized the NL and TV techniques. Deep Learning techniques have been successful in computer vision applications and showed promising performances when applied to PolSAR filtering [19]- [21].
Recently, a PolSAR image speckle reduction algorithm based on a new definition of similarity coefficient between coherency matrices is proposed [22]. A polarimetric despeckling filter based on a generalization of the structure tensor has introduced [23]. F. Shang et al. introduced a despeckling filter based for fully polarimetric synthetic aperture radar (PolSAR) based on the degree of polarization (DoP) information [24]. The mean shift technique, which is commonly used in filtering, segmentation, clustering, and object tracking in the field of pattern recognition, has been generalized in PolSAR filtering [25]. A comparison of the basic PolSAR filters can be found in [26] and [27].
An iterative MMSE (IMMSE) filter has been introduced to filter SAR images [28] and [29]. Hamrouni et al. extended the use of the IMMSE to PolSAR images [30]. The IMMSE filter showed promising results in both SAR and PolSAR filtering.
The objective of this study is the enhancement of the PolSAR speckle filtering performance by hybrid filtering. The first step is to apply one PolSAR filter ensuring high speckle reduction level serving as reference image. Generally, this operation is accompanied by blurring spatial details to a certain degree. Then by applying the proposed iterative filter, the initially highly filtered homogeneous areas are maintained filtered and the initially blurred details are enhanced.
This paper is organized as follows: Section II introduces SAR polarimetry. In section III, the hybrid iterative filter is introduced. Section IV presents results obtained using simulated and real PolSAR images. Finally, section V gives the conclusions of this paper.

II. SAR POLARIMETRY
Generally, PolSAR data are represented by 2 × 2 S scattering matrix which is expressed as [ The subindices h and v are the horizontal and vertical orthogonal polarization, respectively. For reciprocal backscattering case, S hv = S vh . The scattering matrix can be expressed in a vector form in linear basis or in Pauli basis where the superscript '' t '' indicates matrix transpose. The span which represents the total power of a pixel is defined as [1] | · | is the module operation. The complex correlation coefficient is where the indices i and j represent the hh, hv and vv polarizations. The covariance matrix C and coherency matrix T are defined as [1] where E is the mathematical expectation and '' * '' denotes the complex conjugation. The eigen-decompsition of the coherency matrix leads to the definition the three eigenvalues λ, the entropy H , the anisotropy A and the alpha angle α.

A. RELATED WORK
To achieve robust speckle filtering, a PolSAR filter should ensure three constraints [8]: (i) Preserve the polarimetric information.
(ii) Reduce speckle in extended homogeneous areas. The optimal scenario is to average all pixels (i. e.X (i) = C(i)). WhereX (i) is the filtered covariance matrix. (iii) Retain the spatial details in structured areas. The filtered pixelX (i) should maintain the value of the original covariance matrix (i. e.X (i) = C(i)). In general, a linear MMSE filtered pixel verifiesX (i) ∈ C(i),C(i) . The essence of the iterative method is to scanX (i) in C(i),C(i) by the following iterative filtering procedure [30]: Equation (9) can be written as a geometric sequence as, The parameter b k (i) which is the step size controls the convergence of the geometric sequence. The principle of the iterative method is to apply few iterations (i. e. N iterations), then the homogeneous areas should be maintained filtered as in the initially applied filter (i. e.X k (i) ≈C(i)) and the heterogeneous areas should be enhanced (i. eX k (i) ≈ C(i)). To guarantee these results, b k (i) should demonstrate three important numerical properties which are [30]: ( * ) 0 < b k (i) < 1, (to ensure the convergence of the iterative procedure i.e.X k (i) ∈ C(i),C(i) ) ( * * ) b k (i) ≈ 0 in homogeneous areasX k (i) ≈C(i), and ( * * * ) b k (i) ≈ 1 in heterogeneous areasX k (i) ≈ C(i).
The property ( * * * ) represents the ideal case. In heterogeneous areas, any value that satisfies 0 < b k (i) < 1 will The pixel (i) of one diagonal element of the covariance matrix is affected by a multiplicative noise as y is one diagonal element of the covariance matrix. x is the noise-free image to be estimated, and ν is the speckle noise with unit mean and standard deviation σ ν [1].ȳ(i) and var (y (i)) were the mean and the variance of y (i), respectively. In [30], the parameter b k (i) is expressed as wherex k (i) is the filtered pixel at the k th iteration.

B. THE PROPOSED HYBRID ITERATIVE FILTER
From numerical point of view, it can be seen from (8), (9) and (10) that the iterative filter is controlled by two parameters: i) The initial filtered imageX 0 (i) and ii) the parameter b k (i). This study contributes to the polarimetric filtering improvement regarding these two points.

1) THE INITIAL FILTERED IMAGEX 0 (i)
It has been demonstrated in [30], that the speckle reduction level of the filtered imageX k (i) decreases as the number of iteration increases. For this reason, to ensure high speckle reduction level, the equivalent number of looks inX 0 (i) must be sufficiently high. In [30], the boxcar filter is implemented as initial filter (i. e.X 0 (i) =C(i)). This choice permits high speckle reduction level but spatial details are considerably blurred. So that the iterative procedure requires high number of iterations to recover the spatial information which reduced the speckle filtering level in homogeneous areas. To improve the filtering performance, the authors proposed an original filtering scheme by hybridizing two polarimetric filters. The iterative filter is initialized by more sophisticated filter ensuring high speckle reduction level. Then, to enhance the spatial details the iterative filter is applied for few iterations. In this study, the Mean Shift Filter has been implemented as initially applied filter [25].

2) THE PARAMETER b k i
For one look image, b k (i) (12) can be expressed as, By substitutingx k (i) by its meanx k (i), the new expression becomes where where std () is the standard deviation. Hence, the coefficient variation CV(i) is employed to estimate the variability of the pixels rather than the variance in (12). In heterogeneous areas, the CV(i) is high and b k (i) ≈ 1 while in homogeneous area, the CV(i) is low and b k (i) ≈ 0. The parameter b k (i) in (14) can be expressed as where From numerical point of view, it can be seen that • The parameter b k (i) which regulates the speed of convergence is controlled by the CV(i) (i.e. faster convergence in spatial details).
• The objective of f is to normalize the CV(i) between 0 and 1, to ensure the convergence of the iterative method (numerical property ( * )) However, it can be observed that (12) presents the following numerical limitations: 1) f (x) ≤ 0.5 1 (see Fig. 1) which violates the third property needed for b k (i) (i. e. b k (i) ≈ 1 in heterogeneous areas). 2) In Fig. 1, it can be observed that the increase of f is relatively slow which reduces the discrimination of homogeneous areas from heterogeneous ones.
3) The parameter b k (i) is estimated from the filtered image. As a consequence, the CVs have small values and the dynamic range of b k (i) is small. This reduces the discrimination of homogeneous areas from heterogeneous ones and delays the convergence of the process in heterogeneous areas. 4) The parameter b k (i) is estimated using local statistics in a square window. This choice could reduce the accuracy of the estimated statistics. For example, non-similar VOLUME 8, 2020 pixels could be included for the computing of the mean value. To address the limitations 1) and 2), the function f is substituted by a normalizing function g having better numerical properties. The function g is given by where tanh is the hyperperbolic tangent function. It can be observed from Fig. 1 that the function g alleviated the limitations 1) and 2).
To increase the dynamic range of b k (i) and address the numerical limitation 3), the authors implemented the following expression where n is a tuning integer. In (19), In addition, in the initially multilooked data, CVx k (i) and CV y (i) have small values and so b k (i). To fix this problem, the normalized in a homogenous area. It can be easily verified that the expression (19) satisfied the three numerical constraints ( * ), ( * * ) and ( * * * ).
To address the limitation 4), the authors adopted a strategy analog to the non-local filtering [9]. Each pixel in the search area gets a weight computed from comparing its patch with the centered pixel patch. The ENL versus EPD-ROA using different percents of retained pixels (i. e. 25%, 50% and 75 %) has been studied in [29]. To obtain a good compromise between low computing time and good filtering performances, the authors opted to preserve half pixels (i.e 50%) with strongest patches to estimate the local statistics of the pixels as in [29].

3) COMPUTATION OF b k (i): CHOICE OF THE REFERENCE IMAGE
To ensure the preservation of the polarimetric information, it has been suggested that all the elements of the covariance matrices should be filtered equally similarly to multi-look process [8]. Accordingly, the iterative procedure is applied to all elements of the covariance matrix using the same parameter b k (i) computed from a reference image. In the literature, the common used reference image is the span [5] and [8]. However, this choice demonstrates several limitations. First, the equivalent number of looks of the span ENL(span) depends on how the hh, hv and vv channels contribute to the span {i. e. 1<ENL(span)<3 (see appendix)}. As a result, the parameter b k (i) could have different values in different homogeneous areas. Therefore, by using of the span image as a reference image, the homogeneous areas are filtered differently. Secondly, when one channel has low backscattering intensities, its contribution to the span is trivial. Then, the filtering process smoothed its spatial details. To fix all these limitations, the authors proposed to estimate the parameter b k (i) from hh, hv and vv channels separately and not combined in the span.
The maximum is opted to emphasize spatial detail preservation. Hence, when one channel presented high spatial variability, b k (i) = 1. The entire filtered covariance matrix is equal to the original one (i. eX k (i) ≈ C(i)).
In this study, the hh, vv and hv channels are used to estimate b k (i) in (20). This choice exploits the complimentarily between these intensity channels. However, since the diagonal elements of the coherency matrix (i, e. the Pauli basis elements) are more sensitive to the physical scattering mechanism, they could be used to compute the weight b k (i) rather than diagonal elements of the covariance matrix [8].

C. SPECKLE REDUCTION AND POINT TARGETS PRESERVATION 1) SPECKLE REDUCTION
In extended homogeneous areas, CVx k (i) ≈ 0 and CV y (i) = 1 so that the parameter b k (i) ≈ 0 andX N (i) = X 0 (i). As a result, after few iterations (i. e. N iterations), the extended homogeneous areas are maintained highly filtered.

2) PRESERVATION OF POINT TARGETS
Point targets are characterized by strong backscattering reflection. Their preservation by the filtering process is crucial [8]. Three cases can happen • First case: The original filter preserved the point target: In (9)). So that, the point target is preserved by the iterative filter.
• Second case: The original filter degraded the point target: In this case, CVx k (i) and CV y (i) are high. So that b k (i) ≈ 1 and the point target is retrieved in the first iteration.
• Third case: The point is completely smoothed. In the point target, CVx k (i) is low (smoothed point target) and CV y (i) is high (original image) whereas in homogeneous areas, CVx k (i) ≈ 0 and CV y (i) = 1. We let b k1 (i), b k2 (i) and b k3 (i) be the parameters b k (i) of the homogeneous areas, the deleted point target and the degraded point target (i. e. the second case), respectively.
As a result, the point target is recovered using higher number of iterations than of the second case. This scenario has been studied using simulated data.

D. PRACTICAL IMPLEMENTATION OF THE HYBRID ITERATIVE POLARIMETRC FILTER
i. ComputeX 0 image by applying a filter ensuring a high speckle reduction level. For a given pixel of the hh image ii. Compute the parameter CV 2 y0 (i) in (19). iii. Define the search window . Compute the similarity weights (search window , patch ) [9]. Then, discard half pixels having the lowest similarity coefficient values. Compute CVx k (i) using the retained pixels. iv. From the selected pixel of the original image y, discard half pixels and compute CV y (i) using the same process in iii. v. Compute b k,hh (i) using (19). Repeat steps ii-v for hv and vv images. vi. Compute b k (i) using (11). vii. Update the filtered pixel using (20). viii. Apply the process for all pixels of the image.
ix. Repeat iii to viii N iterations.

IV. EXPERIMENTAL RESULTS
In this section, the performance of the proposed method is evaluated on simulated and real PolSAR data, in comparison with the improved Lee sigma filter [8], the IDAN filter [6], and the MSF filter [9]. These filters are implemented with the PolSARpro v6.0 software [32].
To illustrate the effectiveness of the proposed hybrid iterative method, both simulated and real PolSAR images were tested. The first data is Les-Landes image. It is a single look PolSAR data having rectangular zones containing trees with different ages (see Fig. 2. a). The second are the San Francisco Bay PolSAR L-band four-look data. The image contains various typical scattering mechanisms: ocean (surface scattering) city area (high returns from double bounce and point scatterings) and park areas (volume scattering) (see Fig. 2. b). Both real data were collected by the National Aeronautics and Space Administration Jet Propulsion Laboratory (NASA/JPL) Airborne SAR (AIRSAR).
In the simulated data set, two kinds of speckled extended homogeneous areas corresponding to surface scattering and volume scattering have been introduced. In addition, two deterministic targets (i. e. no speckled points and lines) corresponding to diplane scattering and characterized by high backscattering power has been introduced (see Figs. 2. c and d). The procedure of simulating 1-look PolSAR data is given in [8]. For a given covariance matrix C (i. e. ground truth), compute C 1/2 (i. e. C 1/2 (C 1/2 ) * t = C). Then, simulate a complex random vector u having 5000 samples that are complex and have normal distribution, with zero mean and identity covariance matrix. This can be achieved by separately generating the real and imaginary parts of each component of u from a normal distribution with zero mean and 0.5 in variance. The complex single-look vector (2) is given by

B. EVALUATION CRITERIA
To evaluate the performances of the implemented PolSAR speckle filtering algorithms quantitatively in terms of speckle reduction level, the equivalent number of looks ENL is chosen. .
VOLUME 8, 2020 To illustrate the algorithm validity in structural feature retaining, various spatial resolution indexes are used. For simulated data, the mean square error is used.
where M is the number of processed pixels. For real PolSAR images, the Edge Preservation Degree based on the Ratio of Averages (EPD-ROA) is used [33]. The EPD-ROA in horizontal direction is: where m and n are the xy coordinates of the pixel in the selected zone, respectively. EPD-ROA V is obtained by substituting in (24) the indexes (m,n+1) by (m+1,n). For original image, EPD-ROA=1. When the EPD-ROA is closer to one, it means better ability of spatial detail preservation To evaluate PolSAR data preservation, the elements of the covariance matrix, the eigenvalues and Cloude-Pottier parameters [31], i. e. entropy H , the anisotropy A and alpha angle α, were assessed. The span and Pauli decomposition images are considered for visual criterion of intensity PolSAR data filtering.

C. STUDY OF THE PROPOSED ITERATIVE METHOD 1) CONVERGENCE STUDY
Eq. (9) can be expressed aŝ where the function G is It is possible from the fixed point method formula x k+1 = G(x k ) = starting withX 0 (i) =C(i) and for k ≥0, to get closer and closer approximations of the root of the function f ( [34]. To ensure the convergence of the iterative process, |G'(x)|< 1, x>0, then 0 < b k (i) < 2 [34]. Since the sequence {X k (i)} could be negative for 1 < b k (i) < 2, the iterative process converges to the original covariance matrix C(i) for 0 < b k (i) < 1. Fig. 3 displays the ENL of I zone as a function of EPD-ROA H of J zone using the studied filters. As expected, for high number of iterations, the iterative process converged to the original data where ENL=1 and EPD-ROA=1.

2) THE WEIGHT PARAMETER b k
The parameter b k (i) (19) is the key parameter of the proposed iterative filter. It permits an automatic identification of the nature of the processed pixel (i.e. homogeneous or heterogeneous pixel).
These results are traduced by a stable ENL in homogenous areas and an improvement of the EPD-ROA of the spatial detail (see curve for n = 3 in Fig. 3).
For infinite iterations, since 0 < b k (i) < 1, the iterative process converges to the original pixel y(i). In this case, in homogeneous areas, CVx k (i) = CV y (i) = 1. Then, b k (i) = tanh (1) n > 0. As a consequence, the convergence of the iterative in homogeneous areas is relatively high. This result is noticeable in Fig. 3 where the ENL-EPDROA curve decreases rapidly for high number of iterations (see curve for n = 2).
In (19), the parameter n is a tuning integer which controls the scan resolution of the iterative method in [y(i),ȳ(i)]. High values of n decreased the weight b k (i) which delays the convergence of the filtering process. Fig. 3 displays the ENL of I zone as a function of EPD-ROA H of J zone for various values of the tuning integer n. For n>1, the proposed expression (19) outperformed the expression (12). It can be observed also that the increase of parameter n improved the filtering performances but reduced the speed of the convergence of the iterative method. It this study, the selected value is n = 2.

3) THE NUMBER OF ITERATIONS N
In Fig. 3, the studied filters are tuned using the window size parameter. It can be observed that the curve ENL-EPDROA decreases with negative exponential low. By analogy to the previous study, the ENL and EPD-ROA of the iterative method for various iterations were plotted in Fig.3. Compared to the previous study, the ENL decrease is considerably improved (higher ENL and EPD-ROA than boxcar and Lee sigma filter, IDAN and MSF). The main attractiveness of the proposed filter is the enhancement of the spatial details while practically maintaining the initial high speckle 89608 VOLUME 8, 2020   The MSF ensured high speckle reduction level. However, the introduced points and lines were suppressed completely. For the IDAN filter, these details were blurred and the speckle persisted. The Lee sigma filter ensured higher spatial detail preservation at the cost of speckle reduction. The hybrid iterative method maintained the high speckle reduction level of the reference image (i. e. MSF) and retrieved the originally suppressed features to a great extent. Quantitative results in Table 1 agree with visual inspections. The MSF produced the best speckle reduction level with ENL=251 at the detriment of spatial detail preservation (MSE=69058). The IDAN degraded also the spatial details (MSE=68773) with lower speckle reduction level (ENL=26). The Lee sigma filter ensured higher spatial detail preservation (ENL=22361) at the cost of speckle reduction (ENL=13). The hybrid method ensured the best compromise between speckle reduction and detail preservation (ENL=190, MSE=10332).
Figs. 5, 6 and 7 display the entropy, the alpha angle and the anisotropy of the filtered images, respectively. The MSF filter reduced the biases largely (see the mean value of entropy and alpha angle of the volume scattering), but suppressed the polarimetric information of the points and the lines (H = 0 and α = 90 • ). The IDAN filter smoothed the polarimetric information and introduced biases. The Lee sigma filter ensured better polarimetric information preservation. However, the smoothing action around the points and the lines is visible and biases in the extended homogeneous areas are also noticeable due to low speckle reduction level [33]. The hybrid method ensured the best compromise between polarimetric information preservation and bias compensation.     blurred spatial details considerably. The MSF ensured better spatial details preservations but the line is smoothed slightly (see rectangles). The Lee sigma filter preserved better these lines. The hybrid iterative filter corrected the blurring effects introduced by the reference image (i. e. MSF Filter) and produced the highest spatial resolution. Table 2 displays the ENL of I zone and the EPD-ROA of J zone. It can be seen that the hybrid iterative method maintained the high speckle reduction level of the initially applied MSF filter and improved the spatial detail preservation (ENL=37.49≈38.45 and EPD-ROA H = 0.82 0.75). The IDAN and Lee sigma filter produced the same speckle reduction level (ENL=17.67 and 17.80 respectively) with EPD-ROA H = 0.77 and 0.78 respectively. As a result, the hybrid iterative filter ensured the best balance between speckle reduction and spatial detail preservation Figs. 9, 10 and 11 display the entropy, the alpha angle and the anisotropy images of the implemented filters of H zone. The line which is characterized by simple reflection   (H = 0 and α = 0 • ) is smoothed by the MSF and IDAN whereas the hybrid iterative filter enhanced it. The Lee sigma visually showed also high spatial detail preservation particularly in the entropy and anisotropy images whereas the alpha image showed a certain blurring effect. On the other hand, the ability of the MSF filter and the proposed method in terms of speckle reduction and bias compensation in the extended homogeneous areas are comparable. The IDAN and Lee sigma showed lower speckle reduction and bias compensation. Fig. 12 represents the probability density functions (pdf) of the covariance matrix-related parameters of I zone using Lee sigma filter and the proposed hybrid iterative method. It can be observed that the proposed method better reduced the speckle than Lee sigma filter since the pdfs are more compacted. It can be seen also that proposed method better compensated the biases that are related to the parameters of the covariance matrix such as the coherences, phasedifferences, the eigenvalues and the Cloude-Pottier parameters. For example, for underestimated parameters (e. g. the entropy) the pdfs are translated to the right (i. e. to the mean values) ensuring better bias compensation. The opposite case can be observed for overestimated parameters such as the anisotropy [35]. Fig. 13 displays the filtered Pauli color coded images of the real PolSAR data. The IDAN filter blurred spatial details considerably. The Lee sigma filter ensured better spatial details preservations but dark lines in the park are smoothed slightly (see circles). The MSF preserved better these features but it smoothed features in the rectangles. The hybrid iterative filter corrected the blurring effects introduced by the reference image and produced the highest spatial resolution. Quantitative results in table 3 confirmed the visual inspections where the hybrid iterative filter maintained the ENL of the initially applied filter (ENL=31) and enhanced the spatial detail preservation (EPD-ROA H = 0.60 > 0.56). It can be seen also that the hybrid iterative filter outperformed IDAN and Lee sigma in terms of speckle reduction (ENL= 31 > 30 and 27) and spatial detail preservation (EPD-ROA H = 0.60 > 0.58 and 0.56).

3) FOUR LOOKS DATA
Figs. 14 and 15 display the entropy and alpha angle images of the implemented filters of B and C zones. The ability of the MSF filter and the proposed method to reduce speckle and compensate the biases in the extended homogeneous areas are comparable (see rectangles). However, a notable difference between both methods in fine details is observed. The point (see circle) and line (see triangular) are smoothed VOLUME 8, 2020  by the MSF filter. The hybrid iterative filter retrieved these features. The IDAN and Lee sigma filters show this spatial information with a certain degree of smoothness. In C zone, which represents the city zone characterized by double bounce reflection (H = 0 and α = 90 • ), the MSF filtered image exhibited high entropy and low alpha angle values (see triangles). Consequently, the MSF filter produced erroneous information in these areas. The inaccuracy of MSF filter is due to the effect of mixing different scattering media [36]. The Lee sigma filter preserved the polarimetric information   to some extent. In the IDAN filter, mixing media effect is also noticeable. The hybrid iterative filter gives the same entropy and alpha angle information as the Lee sigma filter, but it demonstrates higher spatial resolution. Fig. 16 displays the anisotropy images of the implemented filters of D and E zones. From D zone, it can be observed that the IDAN produced erroneous anisotropy information in the ocean. The hybrid iterative method produced the same anisotropy information as the MSF with a notable correction of the initially blurred details see (squares). Compared to the Lee sigma filter, the hybrid iterative filter ensured better smoothing in the ocean area and better spatial detail preservation (see squares). Concerning the E zone, it can be seen that the hybrid iterative filter conserved the information produced by the reference filter (i. e. MSF) whereas the IDAN filter produced erroneous information and the Lee sigma filter did not suppressed speckle in this homogeneous extended area.

4) COMPUTATIONAL COMPLEXITY
For many of the new SAR systems, such as Gaofen-3, TerraSAR-X and Cosmo-SkyMed, their image dimensions can be higher than 10 000 × 10 000 pixels. Consequently, in addition to the filtering performance, the computational VOLUME 8, 2020 complexity of the PolSAR filter is a great concern for practical implementations.
The conventional local and non iterative state-of-the-art filters such as refined Lee, Lee sigma, as well as the boxcar filters are the most computationally efficient method among other better speckle filtering algorithms for voluminous datasets.
Let the size of the image be N x ×N y , the search window size be w 1 ×w 1 and the patch size be w 2 ×w 2 . The computational complexity of the Lee sigma filter is O(N x ×N y ×w 2 1 ). The computational complexity of the NLM filter is O(N x ×N y ×w 2 1 ×w 2 2 ). For the NLM pretest filter [9], w 1 = 11 and w 2 = 3. The computational complexity of the hybrid iterative method is O(N x ×N y ×w 2 1 ×w 2 2 ×N) where, w 1 = 11, w 2 = 3 and N = 3. The computational complexity of other recent state-of the-art PolSAR filters can be found in [14] and [27].

E. DISCUSSION
Since the speckle reduction level of the iterative filter did not exceed that of the original filter, the initial filter should ensure a very high speckle reduction level. In this study, the iterative method is initialized by the MSF filter. However, other state of the art filters could be used. Then, by running few iterations, the highly de-speckled extended areas are maintained filtered and the spatial details are enhanced. The initial filter can be tuned to obtain better balance between speckle reduction and spatial detail preservation. However, due to the ability to retrieve the initially blurred spatial details, initializing the iterative filter by a highly filtered image provides better performance. It is noticed that if the hybrid iterative filter is initialized by the clean image, the proposed filter preserves it since the second part of the right side of (9) is zero (i. e.X (i) = C(i) in heterogeneous zones and b k (i) = 0 in the noise free heterogeneous areas).

V. CONCLUSION
In this study, the enhancement of PolSAR speckle filtering is addressed. An original hybrid filtering scheme is proposed. The proposed iterative filter is initialized by the MSF filter. Then, to enhance the spatial details the proposed iterative filter is applied for few iterations. The key parameter b k of the iterative method has been improved taking into account the following numerical considerations. 1) The coefficient variation (CV) is considered instead of the variance. 2) The statistics are computed using non-local neighborhood instead of local square one. 3) The hyperbolic tangent is applied as a normalizing function. 4) The dynamic of range of b k is increased by multiplying the CVs of the filtered and the originals images. Comparisons with various state of the art PolSAR filters (i.e. MSF, IDAN and Lee sigma) demonstrated the effectiveness of the proposed iterative filtering method.
Dr. Abdelfattah was an Elected Member with the University Council of Carthage, from 2011 to 2017. He has served as a member with the Executive Committee for the IEEE Tunisia Section, from 2013 to 2015. He was also an Elected Member with the Scientific Council of the Agence Universitaire de la Francophonie (AUF), from 2016 to 2018, and a member with the Commission Régionale des Experts, AUF. He serves as the Co-Chairing for the Mediteranean and MENA Geoscience and Remote Sensing Symposium (M2GARSS).
SAMY EL MAHDY received the Bachelor of Science degree in geology from Ain Shams University, in 1987, and the Master of Science degree in remote sensing and the Ph.D. degree in application of GIS and geomatics engineering (geotechnical engineering) from the University of Putra, in 2006 and 2012, respectively. He has more than 25 years of experience in research with particular interest in quantitative remote sensing, GIS, and geosciences. He is currently a Geospatial Scientist with the American University of Sharjah. He has successfully published more than 35 research articles in journals and conference proceedings. He received the Gold Medal from the Malaysian National Space Agency, in 2011, and the First Place from the Rashid bin Humaid Award for Scientific Research, in 2017.