Improved Goldstein Interferogram Filter Based on Local Fringe Frequency Estimation

The quality of an interferogram, which is limited by various phase noise, will greatly affect the further processes of InSAR, such as phase unwrapping. Interferometric SAR (InSAR) geophysical measurements’, such as height or displacement, phase filtering is therefore an essential step. In this work, an improved Goldstein interferogram filter is proposed to suppress the phase noise while preserving the fringe edges. First, the proposed adaptive filter step, performed before frequency estimation, is employed to improve the estimation accuracy. Subsequently, to preserve the fringe characteristics, the estimated fringe frequency in each fixed filtering patch is removed from the original noisy phase. Then, the residual phase is smoothed based on the modified Goldstein filter with its parameter alpha dependent on both the coherence map and the residual phase frequency. Finally, the filtered residual phase and the removed fringe frequency are combined to generate the filtered interferogram, with the loss of signal minimized while reducing the noise level. The effectiveness of the proposed method is verified by experimental results based on both simulated and real data.


Introduction
As an all-weather all-time remote sensing technique, the synthetic aperture radar (SAR) is of great importance in many fields, such as natural hazards' monitoring, ocean investigation, geographic mapping, and so on [1][2][3]. Synthetic aperture radar interferometry (InSAR), a further development of the traditional SAR technology, employs two or more SAR antennas to retrieve the height profile or the deformation of the ground surface [4][5][6]. Due to its high measurement accuracy, the InSAR technique has been applied in a wide range of areas, such as forest height estimation [7], electronic warfare [8], ocean measurements [9] and geolocation determination of ground targets [10]. The main InSAR processing procedures include image registration, interferogram generation, phase unwrapping, etc. [4,11]. In practice, the quality of the generated interferogram is limited by phase noise due to co-registration errors, thermal noise, temporal decorrelation, baseline decorrelation, electromagnetic interference, and so on [12,13]. To reduce the phase unwrapping difficulty and improve the precision of the unwrapped phase, phase filtering has become an essential step for InSAR data processing [14][15][16]. An ideal phase filter should be able to reduce phase residues significantly while preserving the fringe details well [17]. Since failure in edge preservation will cause serious errors in subsequent procedures, fringe preservation is also a very important topic.
In recent years, various filtering methods were proposed to improve the quality of the interferometric phase. The simplest ones are the mean filter in the spatial domain and the low-pass filter in the frequency domain. As an extension, in the spatial domain, the fringe adaptive smoothing approach, proposed by Lee et al. (commonly referred to as the Lee filter) [18], has been successfully applied for interferometric phase processing with direction-dependent windows. Another extension is the slope-compensated mean filter based on local fringe frequency estimation (topography adaptive filter) [19], and its effectiveness is greatly affected by the size and shape of the filtering window [19].
In the frequency domain, the coherent signal in the filtering patch accumulates and forms a dominant peak, while the power of uncorrelated noises disperses stochastically in different directions [20,21]. Based on this characteristic, many filtering algorithms have been further developed [14,[22][23][24][25][26]. The traditional Goldstein filter, as a low-pass filtering method, smoothes the intensity of Fourier-transformed samples in overlapped interferogram patches [22]. It is widely used for InSAR because of its notable noise suppression capability and fast operation [14]. However, one disadvantage of this method is that it destroys phase continuity in those dense fringe regions. To keep more texture details in the interferogram, Baran et al. proposed an adaptive Goldstein filtering method with the parameter alpha varying according to the coherence of the filtering window [23]. To improve the estimation accuracy of alpha, Song et al. proposed modified Goldstein filters based on empirical mode decomposition and the adaptive-neighbourhood technique [24,25]. All of these modified versions have been shown to be more effective than the traditional Goldstein filter. However, although these methods can reduce most types of phase noise, their ability to preserve fringes and edges is limited. Since the filter response in each overlapping window can be essentially considered as a low-pass filter, the high frequency components of fringes are suppressed [26]. Therefore, these methods may result in loss of fine details in an interferogram, especially in areas with dense fringes and complex textures.
To address the aforementioned problem, in this work, an enhanced Goldstein filtering method is proposed to preserve the fringes by removing the local fringe frequency before phase filtering. Firstly, a mean-filter moving patch is utilized to prefilter the pixels before fringe frequency estimation. The size of the prefilter window is adjusted with the mean coherence value and the PSD (phase standard deviation) in each filtering patch. In this way, the areas with low coherence or high noise level are filtered with larger windows. Since the phase precision decreases sharply when the size of mean filtering approaches to or exceeds the critical averaging look numbers [27], the size of the prefiltering window is limited to prevent the distortion of fringes. Therefore, the subsequent frequency estimation, achieved by Fourier transforms, will be more efficient and accurate. Then, the detected fringe frequency in the estimation window is removed from the original noisy phase, and the remaining part is filtered by the modified Goldstein filter. The filtering parameter α, changing according to the mean coherence and the re-estimated residual phase frequency, not only prevents the areas of high coherence (less noise level) from being over-filtered, but also allows stronger filtering in low coherence (high noise level) regions. Finally, the local fringe frequency and the filtered residual phase are combined to generate the filtered interferogram. Compared to the original Goldstein interferogram filter, there are three changes proposed: prefilter with an adaptive window size, the application of fringe frequency estimation and the optimization of the parameter α. As shown in the simulation results, the new method can suppress noise more effectively while preserving the fringe well, even for fringes with strong curvatures.
The remainder of this paper is organized as follows. The improved Goldstein filter based on local frequency estimation is introduced in Section 2, and evaluation results based on the simulated dataset and collected real data are provided in Section 3. Conclusions are drawn in Section 4.

Analysis of the Goldstein Filter
The Goldstein filter smooths the intensity of Fourier-transformed samples of small, overlapping interferogram patches. The spectrum of the filtered interferogram can be expressed as [22]: where f x and f y respectively represent the spatial frequencies in range and azimuth, α is the filter parameter, Z f x , f y represents the Fourier spectrum of each filtering window and S {·} is a smoothing operator, which is normally achieved by a low-pass filter. The patches are defined as small windows of the interferogram, with overlaps to maintain continuity at the boundaries. The parameter α, taking a value in the range of [0, 1], indicates the desired effectiveness level of the filtering operation. For which means no filtering applied. The filtering effect will become more significant with the increase of α. A large value of α will result in a loss of resolution in the filtered phases, while a small value will reduce the ability of noise suppression. In general, it is difficult to choose an appropriate value for α, and for the original Goldstein filter, α = 0.5 is normally used to ensure a balance between noise suppression and phase preservation [24,25]. To find the best value for α, in [23], the relationship between α and the mean value of the absolute coherence is derived as: where γ is the mean coherence value of the effective patch (patch minus overlap). Clearly, this choice prevents areas of high coherence from being over-filtered and meanwhile allows strong filtering in areas of low coherence, which can effectively reduce the loss of resolution in areas of high coherence.

Combination of Goldstein Filter and Local Frequency Estimation
The Goldstein filter discussed above has the following drawbacks. In each filtering window, the smoothing operator removes high frequency components of the noisy interferogram. As a result, the dense fringe, composed of high frequencies in the filtering window, is also removed. Therefore, the filtering process may destroy the fringe frequency and result in the loss of fine details, especially in areas of complex textures.
To overcome this problem, we can incorporate the fringe frequency estimation technique into the Goldstein filter, which suppresses the noise of slope-compensated phase after estimating and removing the local fringe frequency in the filtering window. The overall flowchart for the improved method is shown in Figure 1, with the following three major steps: 1.
The proposed adaptive mean filter is applied to ensure the accuracy of fringe frequency estimation. The prefilter window size, limited by the critical averaging look number, is varying according to the mean coherence value and PSD.

2.
Fringe frequency estimation using Fourier transform is performed after adaptive mean prefiltering. Note that the estimated principal phase component is removed from the original noisy phase rather than the prefiltered phase. Hence, the prefiltering operation improves the accuracy of fringe frequency estimation and does not reduce the resolution of the interferogram.

3.
The Goldstein filter is utilized to smooth the residual noisy phase with modified parameter α dependent on both the coherence map and residual phase frequency. The filtered residual phase and the removed fringe frequency are ultimately combined to derive the filtered interferogram.
There are three improvements to the original Goldstein filter: prefiltering with size-varied windows, estimating fringe frequency before the Goldstein filter and optimizing the Goldstein filter parameter α. In the following, we provide more details for each part.

Size-Varied Windows Prefilter
Initially, to improve the accuracy of fringe frequency estimation, the prefiltering operation is implemented. To deal with phase noise, an adaptive mean filter is applied to prefilter the estimation patches. Thus, the complex interferogram of the centre pixel in the filtering window becomes: where S (x, y) denotes the original complex interferogram. The range and azimuth radius of the mean filtering patch are respectively limited by: where min {·} takes the minimum value of its two parameters, γ is the mean coherence value in the estimation window, σ denotes the PSD in the estimation window and · rounds down its parameter to the nearest integer. n r and n a represent the critical averaging look numbers in the range and azimuth directions, respectively. The PSD is calculated by [25]: where ϕ (x, y) represents the phase value in the estimation window and ϕ (x, y) denotes the linear phase ramp in a moving patch. A larger PSD means a rough phase with more noise. In the range and azimuth directions, the critical averaging look numbers calculated by a priori information of the SAR imaging mode and control point altitude are respectively expressed as [27]: where λ is the radar operating wavelength, r is the slant range between the master antenna and the target on the ground, B is the length of baseline, R r is the slant range resolution, θ is the look angle, β represents the angle of ground inclination, r g is the projection of r on the ground and H denotes the satellite flight altitude. As can be seen, the areas of lower coherence or more noise will be prefiltered with larger windows, while the low noise region will be prefiltered with smaller windows.

Principal Phase Component Estimation
Next, fringe frequency estimation is performed after the above prefiltering operation. According to the model of the noisy phase described in [21], the power spectrum is characterized by a principal narrow-band component and a broadband component. The principal phase component corresponds to the fringe frequency of the real phase, while the broadband component includes the residual phase component and the phase noise. To derive the local fringe frequency of the filtering window, the maximum likelihood (ML) method [15] is employed using Fourier transforms [28][29][30]. The local frequency of a (2P + 1) × (2Q + 1) window in the interferogram can be expressed as: where S (x, y) represents the prefiltered complex interferogram and (x 0 , y 0 ) denotes the centre pixel in the local fringe frequency estimation window. To speed up the optimization of Equation (9), FFT (fast Fourier transform) is usually employed. Then, the slope-compensated pixel in a (2P + 1) × (2Q + 1) patch can be obtained as: where x ∈ [x 0 − P, x 0 + P], y ∈ [y 0 − Q, y 0 + Q], S (x, y) represents the original noise phase before prefiltering and S (x, y) denotes the slope-compensated complex phase. Since the estimated principal phase component is subtracted from the original noisy phase, the prefiltering operation improves the accuracy of fringe frequency estimation without reducing the resolution of the interferogram.

Residual Noisy Phase Filter
In this part, the residual noisy phase, containing the noise and the residual phase component, is further discussed. In an ideal interferogram without noise, the residue phase without fringe frequency approaches zero in the frequency domain. Nevertheless, the frequencies of the residual phase are not close to zero in practice due to different types of phase noise. Therefore, a change to the Goldstein parameter α is required depending on both the residual phase frequency and coherence map. Firstly, the dominant frequency spectrum amplitude of the residual noisy phase in a (2P + 1) × (2Q + 1) window can be obtained as: where f x_res and f y_res respectively represent the range and azimuth frequency component of the residual noisy phase in the (2P + 1) × (2Q + 1) window. In Equation (11), a larger frequency means greater noise intensity. Therefore, we can combine the dominant frequency spectrum amplitude and coherence value to modify the Goldstein filter parameter into: This modification ensures that the incoherent or noisy areas are filtered more than coherent or low-noise ones. The filtering operation is implemented by combining Equations (1) and (12). To express it more clearly, we define the filtered residual complex phase asŜ (i, j). Accordingly, the derived interferogram filtered by our method can be written as: where (i, j) denotes the position in the interferogram. f i ,f j represents the estimated frequencies derived by Equation (9) in a local window centred on (i, j). In Equation (13), the processed interferometric phase is composed of two parts: the estimated fringe frequency and the filtered residual phase. This enhanced Goldstein filter will effectively prevent the local fringe frequency from being suppressed by the low-pass filter. Therefore, slope compensation can effectively reduce phase noise and at the same time maintain the fringe details. Figure 2 shows the detailed procedure of our method with a patch of 31 × 31 pixels. The original noisy patch and the corresponding real phase are respectively presented in Figure 2a,b. Then, the phase after prefiltering is obtained, as shown in Figure 2c. Subsequently, the dominant frequency spectrum amplitude is derived in Figure 2d, and the local fringe component is produced in Figure 2e. By subtracting the local fringe component from the original noisy phase, the residual phase is obtained as Figure 2f. The residual noisy phases, composed of residual phase components and phase noise, are further filtered by using Equations (1) and (12), as shown in Figure 2g. Finally, the local fringe component and filtered residual phase constitute the filtered interferogram (Figure 2h), which is closer to the real phase ( Figure 2b). As can be seen, after these modifications, texture details of the effective phases are preserved while the regions with strong noise are filtered more, which provides a more accurate filter for the interferogram.
Sensors 2016, 16,1976 6 of 17 greater noise intensity. Therefore, we can combine the dominant frequency spectrum amplitude and coherence value to modify the Goldstein filter parameter into: This modification ensures that the incoherent or noisy areas are filtered more than coherent or low-noise ones. The filtering operation is implemented by combining Equations (1) and (12). To express it more clearly, we define the filtered residual complex phase as Accordingly, the derived interferogram filtered by our method can be written as: where   , i j denotes the position in the interferogram.  ,i j f f represents the estimated frequencies derived by Equation (9) in a local window centred on   , i j . In Equation (13), the processed interferometric phase is composed of two parts: the estimated fringe frequency and the filtered residual phase. This enhanced Goldstein filter will effectively prevent the local fringe frequency from being suppressed by the low-pass filter. Therefore, slope compensation can effectively reduce phase noise and at the same time maintain the fringe details. Figure 2 shows the detailed procedure of our method with a patch of 31 × 31 pixels. The original noisy patch and the corresponding real phase are respectively presented in Figure 2a,b. Then, the phase after prefiltering is obtained, as shown in Figure 2c. Subsequently, the dominant frequency spectrum amplitude is derived in Figure 2d, and the local fringe component is produced in Figure 2e. By subtracting the local fringe component from the original noisy phase, the residual phase is obtained as Figure 2f. The residual noisy phases, composed of residual phase components and phase noise, are further filtered by using Equations (1) and (12), as shown in Figure 2g. Finally, the local fringe component and filtered residual phase constitute the filtered interferogram (Figure 2h), which is closer to the real phase ( Figure 2b). As can be seen, after these modifications, texture details of the effective phases are preserved while the regions with strong noise are filtered more, which provides a more accurate filter for the interferogram.

Results and Analysis
In this section, the performance of the proposed method is studied using simulated and real data. The effect of the three proposed modifications to the standard Goldstein filter is demonstrated by simulated data; both the simulated and the real dataset are used to show the effectiveness of the proposed method in phase noise reduction and fringe preservation in comparison with several existing methods.

Comparison with Our Modifications
The interferogram is derived from the simulated complex SAR image pairs through interferometry. The SAR images with 150 × 150 samples are initially simulated according to the SAR geometry and the imaging area [31]. Subsequently, image co-registration is performed to derive the noisy interferogram, as shown in Figure 3a. The corresponding coherence map calculated by the coherence estimator [32] with 3 × 3 windows is presented in Figure 3b. For comparison, the true phase of the mountainous terrain is given in Figure 3c using the method in [33]. The corresponding phase error of the original noisy interferogram is shown in Figure 3d.

Results and Analysis
In this section, the performance of the proposed method is studied using simulated and real data. The effect of the three proposed modifications to the standard Goldstein filter is demonstrated by simulated data; both the simulated and the real dataset are used to show the effectiveness of the proposed method in phase noise reduction and fringe preservation in comparison with several existing methods.

Comparison with Our Modifications
The interferogram is derived from the simulated complex SAR image pairs through interferometry. The SAR images with 150 × 150 samples are initially simulated according to the SAR geometry and the imaging area [31]. Subsequently, image co-registration is performed to derive the noisy interferogram, as shown in Figure 3a. The corresponding coherence map calculated by the coherence estimator [32] with 3 × 3 windows is presented in Figure 3b. For comparison, the true phase of the mountainous terrain is given in Figure 3c using the method in [33]. The corresponding phase error of the original noisy interferogram is shown in Figure 3d. We filter the simulated data using five algorithms: the reference Goldstein interferogram filter [23], the improved filter with Modification 1 only (adding the adaptive prefilter to the reference Goldstein filter), the improved filter with Modification 2 only (removing fringe frequency before applying the reference Goldstein filter), the improved filter with Modification 3 only (with the We filter the simulated data using five algorithms: the reference Goldstein interferogram filter [23], the improved filter with Modification 1 only (adding the adaptive prefilter to the reference Goldstein filter), the improved filter with Modification 2 only (removing fringe frequency before applying the reference Goldstein filter), the improved filter with Modification 3 only (with the Goldstein filter parameter α dependent on coherence and noise frequency) and the proposed method. The filtered results and the corresponding error images are respectively shown in Figure 4a-e. All of these methods use an 11 × 11 window for Goldstein phase noise filtering. Comparing the proposed Modification 1 with the reference Goldstein filter, as shown in Figure 4a,b, the residues are dramatically reduced due to the noise suppression capability of the adaptive mean prefiltering operation, leading to a more accurate fringe frequency estimation result based on Figure 4b. Figure 4c indicates that Modification 2 has a great advantage of preserving fringe continuity. The proposed Modification 3 makes Goldstein filter parameter α dependent on the absolute value of noise frequency, as well as coherence, as shown in Figure 4d, where although several dense fringes are not preserved well, such a modification has effectively increased the smoothness of the filtering result. The proposed method, combining the advantages of Modification 1, Modification 2 and Modification 3 (see Figure 4e), provides promising performance in both noise reduction and fringe preservation.
Sensors 2016, 16,1976 8 of 17 Goldstein filter parameter  dependent on coherence and noise frequency) and the proposed method. The filtered results and the corresponding error images are respectively shown in Figure 4a-e. All of these methods use an 11 × 11 window for Goldstein phase noise filtering.
Comparing the proposed Modification 1 with the reference Goldstein filter, as shown in Figure 4a,b, the residues are dramatically reduced due to the noise suppression capability of the adaptive mean prefiltering operation, leading to a more accurate fringe frequency estimation result based on Figure 4b. Figure 4c indicates that Modification 2 has a great advantage of preserving fringe continuity. The proposed Modification 3 makes Goldstein filter parameter  dependent on the absolute value of noise frequency, as well as coherence, as shown in Figure 4d, where although several dense fringes are not preserved well, such a modification has effectively increased the smoothness of the filtering result. The proposed method, combining the advantages of Modification 1, Modification 2 and Modification 3 (see Figure 4e), provides promising performance in both noise reduction and fringe preservation.  In order to evaluate the performance of our modifications in a quantitative way, the phase residue number, the edge preservation index (EPI) and the mean-square errors (MSE) between the filtered and the true phase are calculated. The number of residues determines the effectiveness of eliminating phase noise, and the EPI, indicating the capability of fringe and edge preservation, is given by [24]: In order to evaluate the performance of our modifications in a quantitative way, the phase residue number, the edge preservation index (EPI) and the mean-square errors (MSE) between the filtered and the true phase are calculated. The number of residues determines the effectiveness of eliminating phase noise, and the EPI, indicating the capability of fringe and edge preservation, is given by [24]: where ϕ o (i, j) and ϕ s (i, j) represent the edge phase before and after filtering, respectively. In this experiment, ϕ o (i, j) is the true phase value. A better fringe means an EPI closer to one. Clearly, EPI > 1 indicates more noise or false fringe details on the filtered interferogram, while EPI < 1 means the fringe curvature of the filtered phase is decreased in comparison with the "true" phase. MSE, which is the mean-square error of the filtered phase, is defined as [15]: where E [·] denotes the statistics expectation,φ (i, j) represents the filtered InSAR phase and ϕ real (i, j) is the ideal phase value. The filtered interferograms are evaluated using the aforementioned criteria. The results are shown in Table 1. It can be seen that Modification 1 is mainly for reducing the residues, while Modification 2 and Modification 3 have a good performance in fringe and edge preservation. Furthermore, all of those modifications have smaller MSE values than the reference Goldstein filter, which indicates an increased accuracy by the proposed method. These results clearly demonstrate that the proposed filtering method has a much better performance in reducing noise and preserving fringes than the reference Goldstein interferogram filter.

Comparison with Other Filters
In this part, the performance of our method is evaluated with the simulated dataset. As a comparison, the reference Goldstein filter [23], the topography adaptive filter [19] and the Lee filter [18] are also used, and the corresponding filtered results of different methods are shown in Figure 5a-c, where 11 × 11 phase denoising windows are used for all cases.
As shown in Figure 5a, the reference Goldstein filter is less capable of edge preservation in low-coherence regions compared with the true interferogram (Figure 3c). In Figure 5b, the reference topography adaptive filter induces phase errors along the fringe because of the inaccurate estimation of fringe frequency. In Figure 5c, several residues are introduced by the Lee filter, since a high noise level will cause serious problems for the orientation determination of the Lee filter, and a wrong window direction will result in an inaccurate filtered result. Figure 5d again demonstrates the superiority of the proposed method.
The calculated quantitative results are shown in Table 2. We can see that the filtered interferogram of the proposed method has the least residue number among these methods, indicating a very high noise reduction result. Moreover, the EPI and MSE obtained by our method have almost reached those of the real phase. This quantitative evaluation has again verified that the proposed method can not only eliminate phase noise more effectively, but also preserve the fringe well.  Then, in order to validate the robustness of our method in filtering the transitional region and phase jumping region, two cross-sections (A in the azimuth direction and B in the range direction) where phase distortion often occurs are extracted from Figure 3a. By comparing the cross-sections, the advantage of our method can be revealed.
As show in Figure 6a, the filtered phase obtained by our method is much closer to the real phase than others. In Figure 6b, the Goldstein filter and the Lee filter cannot achieve a very good result in the phase jumping area with a high level of noise and produce a shift in the fringe peak. In both figures, the topography adaptive filter leads to some variations in the transition zone. As can be seen, the proposed method has the best edge preservation performance among all four methods considered.
(a)   Then, in order to validate the robustness of our method in filtering the transitional region and phase jumping region, two cross-sections (A in the azimuth direction and B in the range direction) where phase distortion often occurs are extracted from Figure 3a. By comparing the cross-sections, the advantage of our method can be revealed.
As show in Figure 6a, the filtered phase obtained by our method is much closer to the real phase than others. In Figure 6b, the Goldstein filter and the Lee filter cannot achieve a very good result in the phase jumping area with a high level of noise and produce a shift in the fringe peak. In both figures, the topography adaptive filter leads to some variations in the transition zone. As can be seen, the proposed method has the best edge preservation performance among all four methods considered.
Sensors 2016, 16,1976 12 of 17 (d)  Then, in order to validate the robustness of our method in filtering the transitional region and phase jumping region, two cross-sections (A in the azimuth direction and B in the range direction) where phase distortion often occurs are extracted from Figure 3a. By comparing the cross-sections, the advantage of our method can be revealed.
As show in Figure 6a, the filtered phase obtained by our method is much closer to the real phase than others. In Figure 6b, the Goldstein filter and the Lee filter cannot achieve a very good result in the phase jumping area with a high level of noise and produce a shift in the fringe peak. In both figures, the topography adaptive filter leads to some variations in the transition zone. As can be seen, the proposed method has the best edge preservation performance among all four methods considered.
(a) Finally, the density function of phase errors (filtered phase minus the real phase) within the range [ , )    is shown in Figure 7, and we can draw the same conclusion about our proposed method.

Real Data Experiment
In this part, the interferogram of the SIR C-SAR data (C-band, Etna Volcano of Italy) is used to investigate the performance of the proposed method. A typical area of the experimental results with a 400 × 400 sample dimension is selected.
In Figure 8, the entire phase image is presented in the left column, and the enlarged area in the white rectangle is on the right. The original noisy image is provided in Figure 8a-c. The fringes in Figure 8b, buried in significant noise in the white rectangle, represent the steep terrain of Mount Etna. As shown in Figure 8c, the enlarged area contains more texture details. The filtering results produced by the reference Goldstein filter, the topography adaptive filter, the Lee filter and our  Finally, the density function of phase errors (filtered phase minus the real phase) within the range [−π, π) is shown in Figure 7, and we can draw the same conclusion about our proposed method. Finally, the density function of phase errors (filtered phase minus the real phase) within the range [ , )    is shown in Figure 7, and we can draw the same conclusion about our proposed method.

Real Data Experiment
In this part, the interferogram of the SIR C-SAR data (C-band, Etna Volcano of Italy) is used to investigate the performance of the proposed method. A typical area of the experimental results with a 400 × 400 sample dimension is selected.
In Figure 8, the entire phase image is presented in the left column, and the enlarged area in the white rectangle is on the right. The original noisy image is provided in Figure 8a-c. The fringes in Figure 8b, buried in significant noise in the white rectangle, represent the steep terrain of Mount Etna. As shown in Figure 8c, the enlarged area contains more texture details. The filtering results produced by the reference Goldstein filter, the topography adaptive filter, the Lee filter and our

Real Data Experiment
In this part, the interferogram of the SIR C-SAR data (C-band, Etna Volcano of Italy) is used to investigate the performance of the proposed method. A typical area of the experimental results with a 400 × 400 sample dimension is selected.
In Figure 8, the entire phase image is presented in the left column, and the enlarged area in the white rectangle is on the right. The original noisy image is provided in Figure 8a-c. The fringes in Figure 8b, buried in significant noise in the white rectangle, represent the steep terrain of Mount Etna. As shown in Figure 8c, the enlarged area contains more texture details. The filtering results produced by the reference Goldstein filter, the topography adaptive filter, the Lee filter and our method with a constant denoising window size of 17 × 17 are shown in Figure 8d-o. As can be seen in Figure 8e,f, the reference Goldstein filter has reduced the phase noise effectively, but the fringes are ambiguous and not continuous in dense fringe areas. Comparing Figure 8g-i with Figure 8d-f, we see that the topography adaptive filter is more capable of fringe preservation than the reference Goldstein filter in dense fringe areas. However, the estimated fringe frequency by the topography adaptive method is not accurate, and the fringe details of the resultant interferogram are lost. In Figure 8j-l, the Lee filter shows a better performance in detail preservation than the reference Goldstein filter and topography adaptive filter. Nevertheless, its windows' direction is difficult to determine, especially in areas with a high noise level, leading to broken fringes in some regions. Figure 8m-o shows the filtering result using the proposed method, which has effectively improved the accuracy of fringe frequency estimation and simultaneously smoothed noise. In the enlarged areas of Figure 8, the proposed method gives the best result in preserving texture details and has the most continuous fringe in the steep terrain.
Sensors 2016, 16,1976 14 of 17 method with a constant denoising window size of 17 × 17 are shown in Figure 8d-o. As can be seen in Figure 8e,f, the reference Goldstein filter has reduced the phase noise effectively, but the fringes are ambiguous and not continuous in dense fringe areas. Comparing Figure 8g-i with Figure 8d-f, we see that the topography adaptive filter is more capable of fringe preservation than the reference Goldstein filter in dense fringe areas. However, the estimated fringe frequency by the topography adaptive method is not accurate, and the fringe details of the resultant interferogram are lost. In Figure 8j-l, the Lee filter shows a better performance in detail preservation than the reference Goldstein filter and topography adaptive filter. Nevertheless, its windows' direction is difficult to determine, especially in areas with a high noise level, leading to broken fringes in some regions. Figure 8m-o shows the filtering result using the proposed method, which has effectively improved the accuracy of fringe frequency estimation and simultaneously smoothed noise. In the enlarged areas of Figure 8, the proposed method gives the best result in preserving texture details and has the most continuous fringe in the steep terrain.  Table 3 lists the results for phase residue number and phase standard deviation. In terms of residues in the interferogram, the improvements by the reference Goldstein filter, the topography adaptive filter and the Lee filter are 97.41%, 96.17% and 93.98%, respectively, while our proposed one is 99.05%, which is the best result. Besides, the PSD for the proposed method is the smallest, giving the best smoothing effect. Overall, we can conclude that the proposed interferogram filter based on local frequency estimation has consistently outperformed the existing ones.

Conclusions
In this paper, a modified Goldstein filtering method to reduce the phase noise has been proposed. As demonstrated by experimental results based on both simulated data and real data, it can suppress noise effectively while preserving the fringe details well. In detail, we have the following four major findings:  Table 3 lists the results for phase residue number and phase standard deviation. In terms of residues in the interferogram, the improvements by the reference Goldstein filter, the topography adaptive filter and the Lee filter are 97.41%, 96.17% and 93.98%, respectively, while our proposed one is 99.05%, which is the best result. Besides, the PSD for the proposed method is the smallest, giving the best smoothing effect. Overall, we can conclude that the proposed interferogram filter based on local frequency estimation has consistently outperformed the existing ones.

Conclusions
In this paper, a modified Goldstein filtering method to reduce the phase noise has been proposed. As demonstrated by experimental results based on both simulated data and real data, it can suppress noise effectively while preserving the fringe details well. In detail, we have the following four major findings:

•
The adaptive prefiltering operation based on phase standard deviation and coherence can effectively improve the accuracy of local fringe frequency estimation for areas incoherent or with a high level of noise without reducing the resolution of the interferogram.

•
The fringe frequency estimation and slope compensation before applying the Goldstein filter can significantly enhance its performance in edge preservation.

•
The modified Goldstein parameter α, varying with coherence and the dominant frequency component in the residual noise phase, provides a promising result in noise reduction.

•
Fringe frequency compensation and residual phase filtering are combined to reduce the number of phase residues significantly while preserving the fringe details well, even for fringes with strong curvatures.
As a result, the filtered interferogram can keep its fringe frequency components and also benefit from the reduced noise level provided by the Goldstein filter, leading to an improved performance by the proposed method.