NOISE REDUCTION AND INTERPRETATION OF ICE-PENETRATING RADAR DATA IN ANTARCTIC ICE SHEET BASED ON VARIATIONAL MODE DECOMPOSITION

: Ice-penetrating radar is an effective method for studying the subglacial bedrock and ice information within the Antarctic ice sheet. Because of the low conductivity of ice and the relatively uniform composition of ice sheets in the polar region, ice-penetrating radar can penetrate deeper part of the ice sheet and collect the following data. However, it is still necessary to suppress the noise from radar data to obtain more accurate and effective data. In this paper, the entirely non-recursive Variational Mode Decomposition (VMD) is applied to the data noise reduction of ice-penetrating radar data. VMD is a decomposition method of adaptive and quasi-orthogonal signals, which decomposes airborne radar data into multiple frequency-limited quasi-orthogonal Intrinsic Mode Functions (IMFs). The IMFs containing noise are then removed according to the information distribution in the IMFs component and the remaining IMFs are reconstructed. We implements the method to process the real ice-penetrating radar data, which effectively eliminates the interference noise in the data, improves the signal-to-noise ratio and characterizes the internal layer structure of ice. It is verified that the method can be applied to the noise reduction processing of polar ice-penetrating radar data successful, which provides a better basis for data interpretation. Finally, we present the internal structure within the ice sheet based on VMD denoised radar profile.


INTRODUCTION
Ice-penetrating radar is a geophysical detection method based on electromagnetic wave theory, which uses radar echo to sound the characteristics of ice and snow medium (Bingham al., 2007). However, the original airborne ice radar echo signals usually include: coupling signal between the transmitting and receiving antennas, reflection signal on the ice surface and internal layer within the ice, echo signal from bedrock interface, and system thermal noise and other clutter interference signals. Noise and clutter can interfere with the extraction of effective signals and even mask the echoes of weakly reflective targets. For the acquired radar data, filtering methods in clutter suppression and noise reduction processing are generally implemented in the frequency domain. The Fourier transform assumes that the signal is steady-state, and that the Fourier transform-based filter denoising method is generally effective when the signal to be processed is stationary and its spectral characteristics are significantly different from noise. However, the signals to be processed are often non-stationary nonlinear, and they usually appear as a combination of harmonic components in the frequency domain. In this case, it is difficult to find the true spectral distribution of such signals. Therefore, the Empirical Mode Decomposition (EMD) method is introduced, which is a powerful tool for analysing non-stationary signals and nonlinear signals. Moreover, EMD has higher time-frequency resolution and better noise reduction capability than traditional timefrequency analysis methods such as Fast Fourier Transform (STFT) (Partyka al., 1999), Wavelet Transform (WT) (Chakraborty al., 2010), S Transform (ST) (Odebeatu et al., 2012), so it has also been introduced into seismic data processing. Unfortunately, the EMD method itself has some shortcomings, such as lack of strict mathematical foundation, low algorithm efficiency, modal aliasing problem, and poor noise resistance (Li et al., 2008). Although there are many improvements in the subsequent methods, such as EEMD (Wang et al., 2014), BEMD (Zhang et al., 2018a), CEEMD (Zhang et al., 2018b), there are still problems in the lack of mathematical theory and sensitivity to noise. Dragomiretskiy, K. et al. (2014) proposed (Variational mode decomposition) VMD which can better overcome the shortcomings of the EMD method (Dragomiretskiy et al., 2014;2015), and it immediately received popular attention, and the method was quickly introduced into many fields.
In this paper, the variational mode decomposition (VMD) is applied to the data noise reduction of ice-penetrating radar. VMD is a decomposition method of adaptive and quasiorthogonal signals, which decomposes airborne radar data into multiple frequency-limited quasi-orthogonal eigenmode functions (IMFs). The IMFs containing the noise is then removed based on the information distribution in the IMFs component and the remaining IMFs are reconstructed. It implements the method to deal with the real ice-penetrating radar data, which effectively removes the interference noise in the data, improves the signal-to-noise ratio, and fine characterizes the internal layer structure of ice. It is verified that the method can be applied to the noise reduction processing of polar ice-penetrating radar data very well, which provides a better basis for data interpretation.

VMD principle
VMD is an estimation of signal decomposition based on classical Wiener filtering, Hilbert transform and mixing variational problem (Dragomiretskiy et al., 2015;Xue et al., 2017). The purpose of VMD is to decompose a complex analytical signal into eigenmode functions IMFs with specific sparsity property. IMFs are defined as amplitude-modulated-frequency-modulated (AM-FM) signals ( ) k u t with the limited bandwidth property: where the phase  ( ) k t a nondecreasing function, and the envelope ( ) k A t and the instantaneous frequency both are nonnegative and vary much slower than the  ( ) k t .
Assuming that each modal function IMF is a finite bandwidth around the respective centre frequency  ( ) k t , by estimating the bandwidth of each modal function, a variational constraint model that minimizes the sum of the modal function bandwidths is established: where k u is the k th mode.  k is a centre frequency around which k u is mostly compact. The bandwidth of each mode is determined by the squared H1 norm of its baseband-shifted Hilbert complemented analytic signal with only positive frequencies. A quadratic penalty and Lagrangian multipliers are introduced to address the formula (2). The augmented Lagrangian L is set using the following equation: where  is the balancing parameter of the data-fidelity constraint and  is the Lagrangian multiplier.
Alternate direction method of multipliers (ADMM) approach is used to solve the variational problem in (3) and produce different decomposed modes and the centre frequency during each shifting operation. Each mode obtained from solutions in spectral domain can be represented as: VMD procedure mainly includes the following steps.
2. Centre frequencies' update. The centre frequencies  1 n k are updated as the centre of gravity of the corresponding mode's power spectrum as shown in (6).
3. Dual ascent update. For all   0 , the Lagrangian multiplier  1n is updated by (7) as dual ascent to enforce exact signal reconstruction until Where  represents the parameter of noise margin. When the signal contains strong noise, in order to achieve good denoising effect, we let   0 .
The detailed complete algorithm of VMD can be found in Dragomiretskiy et al. (2014). The fact that wiener filtering is embedded for modes update makes the VMD algorithm more robust to noise.

Comparison of VMD and EMD decomposition tests for synthetic data
In order to verify the validity of the VMD algorithm, we first use this method to denoise the synthetic data. The synthetic data consists of four basic sub-data, all of which are non-stationary signals. Sub-data 1, within 0~0.4 s, is a sine wave with amplitude value A = 1, frequency f = 50 Hz, A = 0 within 0.4~1 s; sub-data 2, within 0.3~0.7 s, is a sine wave with amplitude value A = 1, frequency f = 100 Hz, and the rest part A = 0; subdata 3, within 0.6~1 s, is a sine wave with amplitude value A = 0.5, frequency f = 200 Hz, A = 0 within 0~0.6 s; sub-data 4 is the random noise. We select the EMD method for comparison.
The synthetic data is processed by the VMD method, and k is taken as 4, and the obtained result is shown in Figure 1. By comparing the test results from the VMD and EMD methods on the synthetic data, it is found that the VMD method not only have a strict mathematical foundation, but also have the following advantages: (1) algorithm is more efficient and faster than traditional EMD methods; (2) the decomposition effect of the non-stationary model is better, and the result is more accurate, which can better overcome the modal aliasing phenomenon of the EMD method; (3) the noise immunity of the VMD is better than the EMD methods. This section uses the EMD method to decompose the synthetic data in Figure 1, and obtains the results shown in Figure 2.

ICE-PENETRATING RADAR FORWARD SIMULATION AND VMD TEST
We use the VMD method to test the simulated ice-penetrating radar data and to verify the noise reduction effect of the VMD method on the ice-penetrating radar signal,. First, a simple ice sheet model is established. Then, the forward modelling is performed to obtain the simulated ice-penetrating data. Finally, the VMD method is used to decompose and reconstruct the analog signal.
The model consists of air, ice cover and subglacial bedrock. It consists of three internal layers, one ice-bedrock interface, one air ice interface, and two internal layers. With regard to the computer performance at hand, a small model is realized. The size of the model was set to 45 m × 40 m, but the ice surface was 0 m as a reference surface, the air layer was 5 m thick. The electrical parameter setting is mainly described by MacGregor et al. (2007). In order to increase the noise ratio in the simulation data, the random medium is considered. The dielectric constant model and the conductivity model as shown in Figure 3  The forward modeling is carried out for the above model using the time domain finite difference method (FDTD). The antenna center frequency is 300MHz, the time window 400ns, the sampling interval 0.101ns, the number of collecting traces 121, and the distance interval 0.375m. The height of the acquisition position is 5m from the ground, and the offset of the antenna is 0.1m.
Setting k equal to 5, VMD is performed to obtain the main noise components. The combined profiles of the components show that the VMD method eliminates a large amount of noise, and through the elimination and recombination, the signal-to-noise ratio of the whole data is significantly improved. The data below the ice surface is taken as test data ( Figure 5a) and processed by VMD. VMD is performed to obtain the five IMF components in Figure 5b-f, where IMF1, IMF3, and IMF4 are the main noise components, and the main useful information is stored in IMF2 and IMF5. Therefore, the composite profile of the component IMF2 and the IMF5 is obtained as shown in Figure 6. From the comparison of Figure 5a with Figure 6 we can see, the latter is apparently clearer. The VMD method eliminates a large amount of noise. Through the elimination and recombination, the signal-to-noise ratio of the whole data is significantly improved.

INTERPRETATION OF ANTARCTIC ICE-PENETRATING RADAR DATA
The airborne radar data is from the 32nd China Antarctic scientific research survey in 2015. There are 1000 totally traces, each of which has 2000 sampling point. The sampling frequency is 50MHz. Before the VMD processing, DC removal, band-pass filtering, and averaging subtraction are carried out. Although the decomposition data used has been processed by conventional processing, its signal-to-noise ratio is still low, and there is a large amount of interference noise, which causes the isochronous layer in the ice being blurred. First, the 405th trace is taken as the test trace to be processed by the VMD method. By comparison, it can be seen that IMF2, IMF3 and IMF5 are components mainly containing noise, while IMF1 and IMF4 components mainly contain ice-bedrock interface and ice internal layer reflection information. The components IMF1 and IMF4 with high signal to noise ratio are selected to obtain the composition profile. For comparison purposes, the raw profile is placed beside the composite profile. Obviously, it can be seen that the composite section reconstructed by VMD is clearer, the interface between the isochronous layer and the bedrock is more prominent, and the noise is greatly reduced, which effectively improves the signal-to-noise ratio. Furthermore, it is explained that the VMD method can be very well used for noise reduction of ice-penetrating radar data.
Based on the information provided by the composition profile, we can see five distinct and continuous electromagnetic wave reflection interfaces resulting from differences in conductivity or permittivity (Figure 7). Set the depth of the first reflective interface to 0 m. According to the time-depth conversion, their average depth from top to bottom is 110m, 245 m, 425 m and 850 m. In general, the first reflective interface corresponds to the snow surface. Below the snow surface, there will be a strong reflection signal at the ice-ice interface due to the electrical difference between snow (Antarctic snow density of about 300 kg·m -3 ) and ice. According to the propagation speed of electromagnetic waves in the Antarctic snow, the depth of the second interface is 110 m, which is consistent with the average snow depth data of the area, so the second interface is determined to be the interface of ice and snow. The third and fourth reflective interfaces are time-depth converted, and their depths are all less than 700 m. According to the foregoing, it is currently assumed that the two interfaces are caused by the difference in the two critical densities of ice, but the specific causes of this need to be supported by the study of ice cores in the area. The fifth reflective interface is an ice-bedrock interface. Figure 7. Internal layer from the airborne ice-penetrating radar data profile The VMD algorithm uses a completely non-recursive approach to achieve frequency domain decomposition of signals and noise. It has the advantages of strong decomposition capability, anti-noise performance and fast processing speed. The real icepenetrating radar data processing example shows that the method can effectively improve the signal-to-noise ratio of the ice-penetrating radar data in the polar regions. It provides an effective means for obtaining more abundant information within the ice sheet and subglcial bedrock, which is beneficial for later interpretation and can be widely used in the noise reduction of ice-penetrating radar data processing. From the processed data, four layers and five interfaces are clearly interpreted by combing the ice sheet knowledge.