Analyzing the Impact of Different Filtering Methods on Satellite Altimetry Full Waveform Decomposition

Filtering is an essential step in the denoising of satellite altimetry full waveform data, since any deformation and distortion in the shape of the waveform can cause errors in range estimation and further waveform decomposition will also be adversely affected. This paper evaluated comprehensive performance of the popular filtering approaches like Gaussian filter, Taubin filter, Wavelet filter, and EMD based filter by simulated waveform data and ICESat/GLAS waveform. Firstly, according to the principle of each filter, the optimal parameters of filtering algorithm by ergodic tests were selected, then the Gaussian function using Levenberg-Marquardt method was used to fit full waveform to exact waveform parameters (i.e. peak amplitude, position, and half-width). Thirdly, through comparing SNR, RMSE of the pre and post filtering simulation waveform, and the consistency ratio, the average error of peak amplitude, position, and half-width in each Gaussian components of the fitted simulation waveform, verified the effectiveness of these filters and analyzed their influence on decomposition accuracy. Both the simulation experiments and ICESat/GLAS experimental results suggested that the Taubin filter had superior performance with the lowest peak position error, which turns out it has advantage in full waveform denoising and contributes to better full waveform decomposition. However, it introduces more parameters needed to be selected. The self-adaptive EMD based approach has the highest consistency, which shows EMD-based method is more suitable for the denoising of satellite altimetry full waveform decomposition.


INTRODUCTION
Space-borne laser altimetry is a powerful remote sensing technology, which can provide high resolution information on surface elevations and ground cover to the Earth science community by digitally recording the return laser pulse. NASA's Ice, Cloud, and land Elevation Satellite (ICESat/GLAS) was the first Earth-orbiting laser altimeter mission that operated from 2003-2009 [1] , its product data have been widely used in monitoring glaciers [2] , ice shelves [3] , vegetation, and land elevation [4] , etc. Building on this heritage, the next generation of space-borne Advanced Topographic Laser Altimeter System (ICESat-2/ATLAS) is seeking to monitor a diverse range of geophysical phenomena at unprecedented levels of precision and accuracy by photon-counting technolgy [5] . Comparing with America, though China started later in this aspect, currently ZY3-02 as the first experimental altimeter acquired 44 tracks effective data, and the combined adjustment of its laser and image data can actually improve the accuracy of uncontrolled elevation [6]. . GaoFen7 (GF-7) is the first civil laser three-dimensional mapping satellite in China, and will launch this year whose main objective is to aid high precision mapping.
Critical to the success of full waveform satellites is the ability to interpret the shape of the digitallyrecorded return laser pulse and extract accurate time for specific reflecting surfaces within the footprint (e.g., first and last), and thus to geolocate the desired reflecting surface directly. Brenner proposed Gaussian function to decompose the waveforms [7]. Specifically, from the real waveform signal, initial parameters need to be identified first for the fitting step, more precisely the number of peaks or modes together with width, amplitude and location for each mode. Due to the noisy nature of many waveforms, estimation of initial values from the raw signal results in a large number of modes with a low amplitude and a narrow width. Therefore, it is necessary to smooth the waveforms in order to avoid the noise-induced Gaussian components and help increase the fitting accuracy.
Many methods have been proposed for both noise removal and improvement of the signal-to-noise ratio (SNR) in full waveform decomposition. But some methods are not appropriate in all cases and need specific assumptions. Typically, in Low-Pass filters, the high frequency segments are assumed to contain noise. Wiener filter requires to pre-estimate a cross-correlation vector between input signal and expected signal, which is rather difficult to carry out because the desired signal is unknown before filtering. To Gaussian filter, the selection of an appropriate kernel width is a key but difficult work, and the denoising effect comparing with other filter like Taubin filter, Wavelet filter and EMD filter lack explicit comparison. In this paper, we report our experience of using those different kinds of filters. In order to keep the results objective, the principle of the filters are first described and their smoothing parameters are investigated and evaluated. And following the results of these methods applied to simulation and ICESat/GLAS data are presented and compared. At last, experimental results are discussed and concluding remarks are offered, with expectation to provide reference for the upcoming GF-7 data process.

Gaussian Filter
Gaussian filtering is realized by convolution between the Gaussian kernel function with the original signal [8] . The smoothing performance of it could be adjusted by Gaussian kernel width σ. The width of flat plane pulse is derived as the best Gaussian kernel width for a fixed width Gaussian filtering [9] . In the paper, we chose width of GLAS transmit pulse 6ns as the kernel width to determine the weight of each bin. When the bin length is 9, the smoothing effect is better through many tests.

Taubin Filter
On the basis of Gaussian filter, Taubin filter adds two scale factors and to adjust the smoothing degree [10] , which effectively avoid the defects of edge data shrink led by Gaussian filter. It uses formula eq.(1) for iterative processing of waveform data.
In the equation, denotes the raw waveform data before filtering; ′ means the waveform after the first smoothing; * is the adjacent point set of i，and j is the point of the set; ′′ stands for the waveform after the second smoothing; λ and μ mean the scale，and impact on ( ′ ) ，which could adopt the reciprocal of point number in the set [11] . In the paper, we chose equals 0.9057, equals -0.9072 by many experiments.

Wavelet filter
While the wavelet filter is quite different from the above two methods, the basic idea of it is to start with a basis wavelet function with regularity, locality and oscillation, and to get a function family through scaling and translation transformations [12].
In this paper, we adopted"Daubechies" basis (db4) with the level 3, and the soft threshold method for estimating wavelet coefficients.

EMD Filter
Empirical mode decomposition (EMD), introduced by Huang et al [13] is a method for decomposing complex, multi-component signal into several elementary Intrinsic Mode Functions (IMFs) as eq. (2) shows. Differing from Fourier and wavelet analysis which have predefined basis, EMD uses only scale and frequency characters of the original signal. EMD is a local, fully data-driven and selfadaptive analysis approach. In order to let its instantaneous frequency have a meaningful interpretation, the IMF has to satisfy two conditions: (1) in the whole data set, the number of extrema and the number of zero-crossings must either equal or differ at most by one; (2) at any point, the local average of upper and lower envelop is zero.
Where 0 ( ) means the smoothing signal after deleting the first j IMFs that are considered as noise. In this paper, if the Hurst exponent of a certain IMF is less than 0.5, then it is regarded as noise.

RUSULTS AND DISCUSSIONS
In this section, we used simulated waveform data and ICESat/GLAS waveform to test each filter performance. Due to the true Gaussian components are known in simulation waveform data, therefore, we compared average absolute error of peak value, peak position, half-width error between the fitted waveform and real waveform. While, the ICESat/GLAS data we just compared these waveform with same Gaussian number and their fitting accuracy.

Simulation tests and results
In order to analyze fitting correct rate and fitting accuracy, a group of the full-waveform echoes 5000 in total were simulated using addition of Gaussian functions and Gaussian white noise. The number of Gaussian function from fitting result was compared with the number of Gaussian function given in simulation signal to verify correctness of fitting result. If these two number are equal and the location error is less than 1ns, the fitting result is considered as consistency. The specific parameters are as follows:  [6,7,8,9,10,15,20,25,30,35,40] Figure 1 presents one of the typical waveforms denoising effect of different filters, we can see that Wavelet filter and EMD filter and keep the shape of waveform better, while the Gaussian filter and Taubin filter slightly weaken peak value. Figure 2 reveals SNR gain [14] (the SNR after waveform smoothing/ raw waveform SNR) all increased with smaller RMSE. When the maximum peak number is three and random SNR, fitting difference of consistency ratio (CR), and absolute error difference of peak value (A), peak location (T), peak half width (σ), and fit deviation (FitDev) of the total waveforms are shown in figure 3. It shows EMD filter has superior performance with biggest CR, and Taubin filter has the highest T precision. (G W T E are abbreviation of Gaussian filter,Wavelet filter,Taubin filter and EMD filter respectively).

ICESat/GLAS data experiments
In this article we selected GLAS01 waveform data on March 9, 2008, and the corresponding GLAS05 peak number, peak value, peak location, and sigma of each Gaussian component as validation parameters. First we chose the available data through ancillary science packet APID availability flag (apid_ADSci_flg), then the land waveform data by waveform type flag (i_waveformType), and sample padding with zero filling by sample padding flag (i_samp_pad) to acquire the 1000 pieces of input waveform data. During decomposition, some constraint conditions were set to fitted Gaussian parameters, like minimum peak value (minPeak>noiseAverage+3*Devnoise), and the peak width (2.5<σ<300). Figure 4 shows the number waveform ratio (number of waveforms could be decomposed successfully/1000, NumberWf Ratio), consistency ratio, and absolute ratio (consistency ratio/number ratio, Absolute Ratio). It depicts some waveforms cannot be decomposed after using Wavelet filter and EMD filter. The reason is that after waveform smoothing, the sigma is too narrow which is less than the width of transmit pulse and thus the peak is deleted during Levenberg-Marquardt optimization. Figure  5 illustrates the root mean squared error of fitting reduced near one time comparing than GLAS05 results, and Taubin filter has the smallest fit deviation.

CONCLUSIONS
The paper evaluated the influence of popular filtering methods on satellite full waveform decomposition. Besides Gaussian filter, the other three filters are all effective denoising methods for full waveform decomposition. And each method has its advantage over others, for example, smoothing parameters in Gaussian are relatively easy to set, while others introduce scale and thresholds coefficients, etc. The simulation experiments turn out that Taubin filter can acquire higher time accuracy, and EMD filter can get higher consistency ratio in waveform decomposition. In ICESat/GLAS waveform data tests, the GLAS05 validation results also show that Taubin filter has better performance. And the fit deviations all excel the results of GLAS05.