the Creative Commons Attribution 4.0 License.
the Creative Commons Attribution 4.0 License.
The baseline wander correction based on the improved ensemble empirical mode decomposition (EEMD) algorithm for grounded electrical source airborne transient electromagnetic signals
Yuan Li
Song Gao
Saimin Zhang
Hu He
Pengfei Xian
Chunmei Yuan
The grounded electrical source airborne transient electromagnetic (GREATEM) system is an important method for obtaining subsurface conductivity distribution as well as outstanding detection efficiency and easy flight control. However, there are the superposition of desired signals and various noises for the GREATEM signal. The baseline wander caused by the receiving coil motion always exists in the process of data acquisition and affects measurement results. The baseline wander is one of the main noise sources, which has its own characteristics such as being low frequency, large amplitude, nonperiodic, and nonstationary and so on. Consequently, it is important to correct the GREATEM signal for an inversion explanation. In this paper, we propose improving the method of ensemble empirical mode decomposition (EEMD) by adaptive filtering (EEMDAF) based on EEMD to suppress baseline wander. Firstly, the EEMDAF method will decompose the electromagnetic signal into multistage intrinsic mode function (IMF) components. Subsequently, the adaptive filter will process higherindex IMF components containing the baseline wander. Lastly, the denoised signal will be reconstructed. To examine the performance of our introduced method, we processed the simulated and field signal containing the baseline wander by different methods. Through the evaluation of the signaltonoise ratio (SNR) and meansquare error (MSE), the result indicates that the signal using the EEMDAF method can get a higher SNR and lower MSE. Comparing correctional data using the EEMDAF and the waveletbased method in the anomaly curve profile images of the response signal, it is proved that the EEMDAF method is practical and effective for the suppression of the baseline wander in the GREATEM signal.
The grounded electrical source airborne transient electromagnetic (GREATEM) system consists of two parts: the ground transmitter and air receiver system. This method takes advantage of the airborne electromagnetic method (AEM) and the magnetotelluric method (MT), which has large detection depth, a higher signal resolution, and outstanding detection efficiency (Mogi, 1998; Smith, 2001).
There are the superposition of desired signals and various noises for the GREATEM signal. The noise may be divided into stationary white noise and nonstationary noise. According to the various noise sources, the noise is usually classified as sferics noise, human electromagnetic noise, and motioninduced noise (Bouchedda et al., 2010; Buselli et al., 1998; Macnae et al., 1984). The sferics noise is mainly caused by the charge–discharge in the atmosphere, and the frequency is within 1 kHz. Human electromagnetic noise is caused by 50 or 60 Hz industrial frequencies and its odd harmonics. Motioninduced noise comes from the receiving coil motion and has its own characteristics, such as being low frequency, large amplitude, nonperiodic and nonstationary. The signal baseline wander caused by motion noise is one of the major interferences with the GREATEM signal. This phenomenon always exists in the process of data acquisition and affects the measurement results. Severe baseline wander leads to inferior resistivity image formation and a lower reliability of inversion explanation in the measured response signal. After removing the above noises, the processed data will be stacked and averaged in the next stage.
Because the receiver system is mounted on aircraft such as the rotorwing unmanned aerial vehicle (UAV) and airship, the GREATEM system is different from AEM. First, during the flight, the amplitude of the mounted coil swing is smaller because the vibration and speed of the aircraft are weaker than the airborne electromagnetic system. Hence, the fake anomalous amplitude of GREATEM caused by baseline wander is smaller than that of the AEM. Second, there is narrower frequency distribution of baseline wander for the GREATEM signal. The frequency distribution of the motioninduced noise is within 1 kHz for the AEM, while the frequency distribution of the baseline is mostly within a few Hertz for GREATEM in the actual measurement. Third, due to the use of miniaturized aircraft for GREATEM, the maximum flight loads are much less than for the AEM. It is impossible to install the complex mechanical structure to suppress baseline wander in the receiver system. Therefore, it is necessary to develop algorithm processing for the suppression of baseline wander.
In the method to suppress baseline wander, on the one hand the mechanical correctional structure and the hardware filter can be installed; on the other hand the digital filter and fitting can be used for data processing. Some studies have focussed on the correctional method to suppress motioninduced noise on the transient electromagnetic system. The Fugro company developed the time domain airborne electromagnetic system, where hardware compensation devices and a notch filter with a centre frequency of 0.5 Hz are installed to correct coil motion in the data acquisition system. Buselli et al. (1998) proposed that a highpass filter with a cutoff frequency of 10 Hz be used to reduce this noise. Lemire et al. (2001) put forward the spline interpolation and Lagrange optimization method to reject lowfrequency noise. Wang et al. (2013) introduced a waveletbased baseline drift correction method using a sym8 wavelet and 10 decomposition layers. The waveletbased method is based on multiresolution decomposition analysis. But because it is difficult for wavelet decomposition to choose an optimal wavelet basis function and layer levels, such methods have poor adaptability. Flandrin et al. (2004) put forward the detrend method based on ensemble empirical mode decomposition (EEMD). This method consists of two steps. First, the trend should be regarded as a baseline estimate which is expected to be captured by intrinsic mode functions (IMFs) of a large index. Second, the reconstructed signal should be amounted to the partial IMFs from lower index to middle index without the higherindex components directly. Liu et al. (2017) focussed on the EEMD method to distinguish and suppress motioninduced noise in a GREATEM system. Because the components containing motioninduced noise are excluded from the reconstructed signal directly, the reconstructed signal was distorted by this EEMD method.
Huang et al. (1998) proposed the empirical mode decomposition (EMD), and then Wu et al. (2009) put forward EEMD. The EMD and EEMD methods are scaleadaptive timedomain methods, which are applied to nonlinear and nonstationary signal decomposition. For nonstationary signal processing, it is necessary to propose the shorttime Fourier transform (STFT) and wavelet transform generally. The main method of STFT is to divide the signal into short time intervals where the signal is approximately stationary and then to perform the Fourier transform of the signal in each time interval to get the frequency distribution. And the main method of the wavelet transform is to utilize a variablescale sliding window where the specific data are approximately stationary on signal. The width of the window is variable for the time and frequency domain. However, because it is difficult to choose an optimal wavelet basis function and the layer levels of wavelet decomposition by the signal itself, this method has poor adaptability. Therefore, the requirement of signal characteristics of the above method is stationary in a specific window in the same way as the Fourier transform.
Different from previous methods, the major advantage of the EEMD is that the decomposition is derived from the signal itself. Therefore, the EEMD analysis is an adaptive decomposition in contrast to the traditional methods, where the decomposition functions are fixed in a specific window throughout processing. In addition, the characteristics of the signal itself are not affected in the sifting process.
According to the characteristics of baseline wander for the GREATEM signal, the EEMD adaptive filtering (EEMDAF) method consists of three steps.

The signal is decomposed into the Nlevel IMF components and the residual component by the EEMD method.

It is prudent to use an adaptive lowpass filter for higherindex IMFs to get a baseline wander estimate.

The denoised signal can be obtained by subtracting baseline wander from the noisy signal.
In the later sections, the correctional result shows that, compared with the waveletbased method and EEMD without the higherindex components, the EEMDAF method is practical and effective for the suppression of the baseline wander in the GREATEM signal.
2.1 EMD method
The signal S(t) is decomposed into Nlevel IMF components and a residual component by the EMD method. The EMD involves the adaptive decomposition of signal S(t) by means of the sifting process. The term IMF is adopted because it represents the oscillation mode embedded in the data. The sifting process is defined by the following steps

Identify levels of decomposition N, and ${r}_{j\mathrm{1}}\left(t\right)=S\left(t\right)$ as residual parameter.

Extract IMF_{j}:
 a.
all extrema of r_{j−1}(t);
 b.
interpolate local maxima and minima as the upper and lower envelopes separately by a cubic spline line, and compute “envelope” E_{min}(t) and E_{max}(t);
 c.
compute the average component $m\left(t\right)=({E}_{\mathrm{min}}\left(t\right)+{E}_{\mathrm{max}}(t\left)\right)/\mathrm{2}$;
 d.
extract the detail component ${D}_{i}\left(t\right)=x\left(t\right)m\left(t\right)$;
 e.
iterate steps (a) to step (d) on the detail component D(t) until the stopping criterion satisfies the threshold, dd<ε; once criterion is achieved, the D(t) is regarded as the effective IMF_{j}, which can also be regarded as the zero mean generally; calculate the stopping criterion:
$$\begin{array}{}\text{(1)}& \mathrm{dd}=\sum {\displaystyle \frac{\mathrm{}{D}_{i\mathrm{1}}\left(t\right){D}_{i}\left(t\right){\mathrm{}}^{\mathrm{2}}}{{D}_{i\mathrm{1}}(t{)}^{\mathrm{2}}}}.\end{array}$$
 a.

Update the residual: ${r}_{j}\left(t\right)={r}_{j\mathrm{1}}\left(t\right){\mathrm{IMF}}_{j}\left(t\right)$; the residual is deemed as the input for a new round of iterations.

Repeat steps (2) and (3) until the value of j is equal to N; the stopping criterion threshold ε is set to between 0.2 and 0.3; the result of the sifting procedure is that S(t) will be decomposed into IMF_{j}(t), $j=\mathrm{1},\mathrm{\dots},N$ and residual r_{N}(t):
$$\begin{array}{}\text{(2)}& S\left(t\right)=\sum _{j=\mathrm{1}}^{N}{\mathrm{IMF}}_{j}\left(t\right)+{r}_{N}\left(t\right).\end{array}$$
2.2 EEMD method
The EEMD method is an improved method based on the EMD method. Compared with the EMD method, the EEMD method resolves the mode mixing problem and achieves better performance by adding white noise to the original signal (Wu et al., 2009). For the EEMD method, the first step is to produce an ensemble of datasets by adding the finite amplitude σ of Gaussian distribution white noise to the original data. The σ stands for standard deviation of white noise. Then EMD method is applied to each realization of datasets to get IMF_{i}(t) with NE times repeatedly. The next step is that the expected ${\widehat{\mathrm{IMF}}}_{j}$ is obtained by averaging the respective components in each realization to compensate for the effect of the addition of Gaussian white noise.
where NE is the ensemble numbers. Finally, the result of the sifting procedure is that S(t) will be decomposed into ${\widehat{\mathrm{IMF}}}_{j}\left(t\right)$, $j=\mathrm{1},\mathrm{\dots},N$ and residual r_{N}(t).
where σ is set to between 0.05 and 0.2 and NE is set to between 100 to 400. In this paper, we set σ and NE to 0.1 and 200 respectively.
2.3 EEMDAF method
The EEMD method is equivalent to a sifting filter which sifts out each IMF component from signal S(t) according to oscillations from fast to slow. The lowerindex IMF component mainly contains fast oscillations; meanwhile the higherindex IMF component mainly contains slow oscillations. The baseline wander is expected to be captured by higherindex IMFs. The simple removal of several higherindex IMFs may introduce significant distortions in the reconstructed signal.
Thus, the baseline wander is distributed over the desired components in several higherindex IMFs. To suppress the baseline wander, this method introduces a group of adaptive lowpass filters to process the several higherindex IMFs successively. The sum of the output of these filters is regarded as the reconstructed baseline estimate. Finally, the denoised signal can be obtained by subtracting an estimated baseline from the noisy signal.
First of all, we suppose the signal S(t) contained severe baseline wander. After processing by EEMD, S(t) will be decomposed into IMFs which be referred to as a_{k}(t).
where N is the number of IMFs. Then, it is important to find out how much IMFs contribute to the baseline wander. Denote this number value as M. The a_{k}(t) is processed from the higher to lower index by lowpass filter h_{k}(t). The output of filter is b_{k}(t).
where ∗ denotes the convolution. The h_{k}(t) is the Butterworth lowpass filter whose cutoff frequency is ω_{k}. As the IMF index decreases, fewer slowoscillations component, but more signal components, are contained in each IMF. So, we design a group of adaptive lowpass filters whose cutoff frequency is decreased as the IMF index decreases. In other words, the first processing of filters is that the last IMF, a_{N}(t), convolves with the first lowpass filter whose cutoff frequency is ω_{N}(t). And the cutoff frequency decreases along with the k index filter decreases.
where α is set to between 0.1 and 0.99 and $k=N,\mathrm{\dots},\mathrm{2},\mathrm{1}$. By this means the filter output b_{k}(t) containing lowfrequency components are extracted from each IMF. Because the algorithm has to be adaptive, the output can be used to determine the value of M as the condition of the reconstructed signal. According to the analysis of the procedure above, the amplitude of the baseline should gradually decrease as the IMF index decreases on filter output b_{k}(t). As a result, to determine the value of M, we consider using the evaluation coefficient function P_{k} as a stopping criterion where the SD(b_{k}) stands for standard deviation b_{k}.
where $k=\mathrm{1},\mathrm{2},\mathrm{\dots},N$. The operator flip is the flipped function that the data rearrange in the opposite direction. The evaluation coefficient threshold δ is set to between 0.05 and 0.1. If P_{k}<δ, the value $M=N+\mathrm{1}k$. In this process, we set ω_{N}, α, and δ to 10, 0.9, and 0.01 respectively. The sum of the output of filtered IMFs whose index is from M+1 to N is regarded as the reconstructed baseline estimate.
Finally, to obtain a reconstructed denoised signal, the baseline estimate is subtracted from the original signal.
3.1 Simulation data
In the GREATEM system, the transmitter injects a bipolar square wave current into the ground; meanwhile the receiver and frontmounted coil were installed on an aircraft in response to the vertical component of the induced electromotive force in a horizontal layered earth model (Nabighian et al., 1988). Response signals are related to the size and depth of the underground conductor, the length and current of transmitter, the equivalent area of receiving coil, the horizontal offset, the flight altitude, and so on. These parameters can be used to calculate the time domain response as a clean signal in the horizontal layer earth model for simulation. In Fig. 1, the model parameters as follows: the length of the transmitter line TX is 1000 m on the ground, the transmitter current is 10 A with 50 % duty cycles at 25 Hz, the equivalent area of receiving coil RX is 1000 m^{2}, the horizontal offset is 50 m, the flight altitude is 35 m, and the sample rate of the receiver is 32 kHz. In this paper, we consider a threelayer earth model where parameters are shown in Table 1. In the end, we calculated the corresponding time domain signal and the vertical response–decay curve on a threelayer earth model.
Because of nonperiodic and nonstationary characteristics of baseline wander, it is difficult to synthesize these noises from a simulation on the computer. The simulated signal is obtained by superimposing the clean signal on the field baseline wander measured by the inertial navigation system. Figure 2a is a simulated noisy signal which is obtained by adding baseline wander to a clean signal with the duration is 10 s. Figure 2b is the field baseline wander noise which is measured by the inertial navigation system with the duration is 10 s.
3.2 Performance of the correction and analysis
In this paper, in order to quantitatively assess the denoised quality between our method and other methods, we propose the signaltonoise ratio (SNR) and meansquare error (MSE) to evaluate the correctional methods in Eqs. (11) and (12), in which S(n) is a synthetic clean signal, $\widehat{S\left(n\right)}$ is a denoised signal, and L is the number of samples. The higher the SNR, the better the correctional effect; the lower the MSE, the better the fitting result.
There are comparisons of SNR and MSE obtained by three correctional methods for the noisy signal whose duration is 60 s. The correctional results are shown in Table 2. In the table, the term noisy signal means the simulated noisy signal before correction. The SNR value shows that the three methods have a remarkable improvement in signal quality. It is proved that EEMDAF and the waveletbased method had better correctional performance than EEMD. Quantitatively, the SNR of the EEMDAF method is significantly close to the SNR achieved by the waveletbased method. It is obvious that the EEMDAF achieves a correction performance similar to the waveletbased method.
For further analysis, the response–decay curve is related to the conductivity of underground geological bodies in the data processing of GREATEM. Besides the SNR comparison above, the original data are used for stacking, averaging, and extracting a secondary field to build a onebyone test point along with the survey path. Then the number of time gates is 24 per test point, where the width of time gates increases approximately logarithmically.
To generate the anomaly curve profile image, we process the simulated noisy data and correctional data by stacking, averaging, extracting, and gating. The anomaly curve profile generated from the clean signal responses is shown in Fig. 3a, where the 24 parallel lines of time gates are represented along with the test point. And Fig. 3b is an anomaly curve profile generated from the noisy signal, where the fake anomaly from Gate 14 to 24 is identified clearly due to baseline wander affecting the horizontal layer model. From Gate 20 to 24, they are mixed with each other. After the processing using the waveletbased method, the anomaly curve profile is shown in Fig. 3c, where the fake anomaly is not accurately represented and the curves are almost parallel to each other. Figure 3d is the anomaly curve profile using the EEMDAF method; it is obvious that the paralleled curves between the gates are better than in the above method.
The comparison of SNR and MSE profiles produced by the datasets by different methods is illustrated in Fig. 4 along with the test point, where SNR and MSE of the noisy signal are marked as reference (black solid curve). In Fig. 4a, the black solid curve shows that the stacking and averaging may produce the improved SNR for a noisy signal. Quantitatively, the EEMDAF and waveletbased method yield an SNR which is significantly higher than the value achieved by the EEMD method. It is observed that there are fluctuations in the SNR using the waveletbased method (red solid curve); meanwhile there are stabilities in the SNR using the EEMDAF method (blue solid curve). And in Fig. 4b, the MSE curves indicate the same conclusion. Results from the comparison of the figures also show that the EEMDAF method significantly outperforms the waveletbased one for the suppression of nonstationary baseline wander.
In October 2018, a field experimental GREATEM survey was performed to detect infiltration water in the refuse landfill of Longquanyi District, Chengdu, in China. The GREATEM system was developed by Chengdu University of Technology. The electrical source transmitter was fixed on the ground; meanwhile the receiver system was mounted on the sixrotor UAV. The survey area and flight paths of the receiver were shown in Fig. 5b. The length of the transmitter line was 1100 m on the ground, the transmitter waveform was a bipolar square wave, and the current was 20 A with 50 % duty cycles at 5 Hz. The receiver system made use of 24 bit analoguetodigital converter whose sample rate was 32 kHz, and the equivalent area of the receiving coil was 1000 m^{2}. The transmitter line was set in the middle of flight paths and almost perpendicular to them. Each length of the flight path was 800 m, and the intervals were 80 m to each other. The flight speed of the UAV was 2.5 m s^{−1}, and the height was 50 m from the ground.
The amplitude of the response will decrease with the transmitter–receiver offset increase. We choose the measured data of part of the flight path L4 for our processing, and the duration of 60 s of data is shown in Fig. 6a. The baseline wander is observed to be significant on the measured data. And Fig. 6b is the correctional result from the EEMDAF. By comparison of the signals before and after processing, the baseline wander is effectively suppressed by the EEMDAF method.
Besides the comparison of the time domain data above, we produced an anomaly curve profile image from the original measured data and the correctional data, which are processed by the waveletbased method and the EEMDAF method respectively. The number of time gates is 18, and the widths increase approximately logarithmically per test point. Figure 7a shows the anomaly curve profile generated from the measured raw data. The correctional data using the waveletbased and the EEMDAF method are shown in Fig. 7b and c respectively. Based on the survey area for refuse landfill, the geological structure can be regarded as layered earth; there may be partial regions where the infiltration water was leaked.
In Fig. 7a, the higher amplitude of response anomaly curves is reflected at 220, 270, and 300 m in the flight survey path. Therefore, the baseline wander exists in the original signal and affects exploration elevation and the anomaly result on inversion. The fake anomalies from Gate 10 to 15 and their interference with one another from Gate 16 to 18 are obvious. In Fig. 7b and c, after using two correction methods, the fake anomalies are suppressed on Gate 10 to 15, and their interference with one another is improved from Gate 16 to 18. In addition, Fig. 7c shows that there is no interference between the gates. However, there is partial interference on Gate 16 to 18 in Fig. 7b. By contrast with the waveletbased method, there is no interference between the last three channels for datasets using the EEMDAF method on an anomaly curve profile. Regarding the decay time of curves, the EEMDAF method holds decay time 4.5 ms more than the waveletbased method to improve the exploration elevation on the survey path. In a comparison of Fig. 7b and c, the results reveal that the performance of the EEMDAF method is significantly superior to the waveletbased method to suppress baseline wander. In a word, the results confirm that the EEMDAF method is an effective and practical correctional method.
Motioninduced noise is usually referred to as baseline wander, which is an inevitable noise and always exists for the GREATEM system. The noise caused by the receiving coil motion has its own characteristics such as being low frequency, large amplitude, nonperiodic, and nonstationary. This phenomenon affects the accuracy of measurements severely, leading to an inferior exploration elevation and a fake anomaly result on inversion. Therefore, we proposed the improved EEMDAF method for baseline wander correction. The noisy signal is decomposed into Nlevel IMF components and a residual component by the EEMD method, and the baseline wander is generally distributed over several higherindex IMFs; then a group of adaptive lowpass filters process these IMFs successively. The baseline estimate is reconstructed by the sum of these filter output. Finally, the denoised signal can be obtained by subtracting an estimated baseline wander from the noisy signal.
First of all, through a comparison of different denoised methods in this paper, the SNR and MSE results show that the denoised performance of the EEMDAF method is significantly superior to the other methods. And the same conclusion can be reached from the anomaly curve profile image. Furthermore, in field data processing, the baseline wander is effectively suppressed by EEMDAF and the waveletbased method. Because there is no interference between the last few gates, the comparison of the anomaly curve profile image reveals that the EEMDAF method is significantly better than the waveletbased method. These results also indicate that the EEMDAF method is a practical as well as an effective method for the suppression of the baseline wander in the GREATEM signal.
In this paper, the data are not publicly accessible because funder terms require the original geological data to be kept confidential without cooperative licensing agreements.
First, YL and SG designed the method model and developed the code. Second, the author SZ designed the field experiments and carried them out with HH and PX on the survey area. Third, YL and CY performed the simulations and processed data. Finally, YL prepared the paper with contributions from all coauthors.
The authors declare that they have no conflict of interest.
The authors thank the members of the project committee for their help. We also acknowledge the referees and the editor for their constructive remarks and suggestions.
This research has been supported by the “Research on fixed wing time domain airborne electromagnetic measurement technology system (grant no. 2017YFC0601900)”.
This paper was edited by Luis Vazquez and reviewed by two anonymous referees.
Bouchedda, A., Chouteau, M., Keating, P., and Smith, R.: Sferics noise reduction in timedomain electromagnetic systems: application to MegaTEM II signal enhancement, Explor. Geophys., 41, 225–239, https://doi.org/10.1071/EG09007, 2010.
Buselli, G., Hwang, H. S., and Pik, P. J.: AEM noise reduction with remote referencing, Explor. Geophys., 29, 71–76, https://doi.org/10.1071/EG998071, 1998.
Flandrin, P., Goncalves, P., and Rilling, G.: Detrending and denoising with empirical mode decomposition, 12th European Signal Processing Conference, Vienna, Austria, September 2004.
Liu, F., Li, J., Liu, L., Huang, L., and Fang, G.: Application of the EEMD method for distinction and suppression of motioninduced noise in grounded electrical source airborne TEM system, J. Appl. Geophys., 139, 109–116, https://doi.org/10.1016/j.jappgeo.2017.02.013, 2017.
Huang, N. E., Shen, Z., Long, S. R., Wu, M. C., Shin, H. H., Zheng, Q., Yen, N. C., Tung, C. C., and Liu, H. H.: The Empirical Mode Decomposition and the Hilbert spectrum for nonlinear and nonstationary time series analysis, P. R. Soc. Lond. A, 454, 903–995, https://doi.org/10.1098/rspa.1998.0193, 1998.
Lemire, D.: Baseline asymmetry, Tau projection, Bfield estimation and automatic halfcycle rejection, Technical Report, THEM Geophysics Inc., PhD, THEM Geophysics Inc, Canada, 13 pp., 2001.
Macnae, J. C., Lamontagne, Y., and West, G. F.: Noise processing techniques for timedomain EM system, Geophysics, 49, 934–948, https://doi.org/10.1190/1.1826986, 1984.
Mogi, T., Tanaka, Y., Kusunoki, K., Morikawa, T., and Jomori, N.: Development of grounded electrical source airborne transient EM(GREATEM), Explor. Geophys., 29, 61–64, https://doi.org/10.1071/EG998061, 1998.
Nabighian, M. N.: Electromagnetic methods in applied geophysics, Volume 1: Theory, Society of Exploration Geophysicists Publications, United Kingdom, 1988.
Smith, R. S., Annan, A. P., and McGowan, P. D.: A comparison of data from airborne, semiairborne and ground electromagnetic systems, Geophysics, 66, 13791385, https://doi.org/10.1190/1.1487084, 2001.
Wang, Y., Ji, Y., Li, S., Lin, J., Zhou, F., and Yang, G.: A waveletbased baseline drift correction method for grounded electrical source airborne transient electromagnetic signals, Explor. Geophys., 44, 229–237, https://doi.org/10.1071/EG12078, 2013.
Wu, Z. and Huang, N. E.: Ensemble empirical mode decomposition: a noiseassisted data analysis method, Adv. Adapt. Data Anal., 1, 1–41, https://doi.org/10.1142/S1793536909000047, 2009.