Artificial source electromagnetic exploration denoising method and system based on STL-FMD decomposition algorithm, terminal and medium

Abstract

The invention discloses an artificial source electromagnetic exploration denoising method and system based on an STL-FMD decomposition algorithm, a terminal and a medium, and the method comprises the steps: obtaining observation signals of an original reference station and a noisy observation station, and carrying out the down-sampling, and obtaining a first observation signal of the reference station and the observation station; performing STL decomposition on the first observation signal of the observation station and eliminating a trend component to obtain a second observation signal of the observation station; noise extraction is realized by using the signal variance ratio of the first observation signal of the reference station and the second observation signal of the observation station for synchronous observation, and a noisy observation signal of the observation station is obtained; performing low-pass FMD decomposition on the noisy observation signal to obtain a noise signal; noise signals are removed, and effective signals are obtained; and replacing the noisy part in the second observation signal of the noisy observation station with an effective signal to obtain a denoised time domain signal. Strong electromagnetic noise and baseline drift are effectively eliminated by means of combining STL and FMD decomposition, and the denoising effect and the signal-to-noise ratio of observation signals are improved.

Description

Technical Field

[0001] The present invention relates to the technical field of artificial source electromagnetic exploration data processing, and in particular to an artificial source electromagnetic exploration denoising method, system, terminal and medium based on a STL-FMD decomposition algorithm. Background Art

[0002] The artificial source electromagnetic method (CSEM) is a geophysical exploration method with controllable source, known signal and strong anti-interference ability. It has been widely used in energy and mineral resource exploration, geological survey and disaster prevention and control. In actual field observations, strong electromagnetic interference is always one of the main factors affecting the effect of electromagnetic exploration; at the same time, due to the complexity and uncontrollability of noise, a single denoising method often cannot achieve satisfactory results. Therefore, how to effectively suppress strong electromagnetic interference has always been a research hotspot and difficulty in this field.

[0003] In the field of electromagnetic exploration data processing, empirical mode decomposition (EMD) can separate different components in the original signal into a series of intrinsic mode functions (IMFs) according to the local time scale characteristics of the signal, and perform feature analysis in the time-frequency domain to achieve signal-to-noise identification and separation. This type of method is also considered to be one of the most effective tools for multi-component analysis of many signals, and is widely used to suppress strong interference, baseline drift and other noise. However, EMD and its improved methods (EEMD, LMD, etc.) have problems such as mode aliasing, endpoint effect, over-envelope and under-envelope, and susceptibility to noise. The variational mode decomposition algorithm (VMD) transforms mode decomposition into a variational solution problem, which effectively makes up for the shortcomings of methods such as EMD; however, the number of modal components K and the penalty factor α in the algorithm parameters determine the decomposition performance. Therefore, when VMD is used without any prior knowledge, there is still a large decomposition error, which reduces the denoising effect. Summary of the invention

[0004] In view of the shortcomings of the existing artificial source electromagnetic exploration denoising technology, the present invention provides an artificial source electromagnetic exploration denoising method, system, terminal and medium based on the STL-FMD decomposition algorithm, which utilizes STL (Seasonal-Trend decomposition procedure using Loess, a time series decomposition algorithm based on local weighted regression) decomposition to effectively eliminate the long-period nonlinear baseline drift (trend component) and large-scale step noise of the observation data, and then utilizes low-pass FMD (Feature Mode Decomposition, characteristic mode decomposition) decomposition to obtain the low-frequency noise of the processing signal, ultimately achieving effective suppression of strong electromagnetic noise and improving the signal-to-noise ratio of the observation data.

[0005] In a first aspect, a method for denoising artificial source electromagnetic exploration based on an STL-FMD decomposition algorithm is provided, comprising the following steps:

[0006] S1: Obtain the observation signals of the original reference station and the noisy measurement station and perform downsampling processing to obtain the first observation signals of the reference station and the noisy measurement station;

[0007] S2: Perform STL decomposition on the first observation signal of the noisy station and eliminate the trend component therein to obtain the second observation signal of the noisy station;

[0008] S3: noise identification and positioning are realized by using the signal variance ratio of the first observation signal of the reference station with high signal-to-noise ratio synchronous observation and the second observation signal of the noisy measurement station, and the noisy observation signal of the noisy measurement station is extracted;

[0009] S4: Perform low-pass FMD decomposition on the noisy observation signal of the noisy station to obtain a noise signal;

[0010] S5: subtract the noise signal from the noisy observation signal of the noisy measuring station to obtain the effective signal after denoising;

[0011] S6: Replace the noisy part of the second observation signal of the noisy measuring station with the valid signal to obtain a complete time domain signal after denoising.

[0012] Furthermore, the downsampling process in step S1 includes:

[0013] The Nyquist sampling theorem is used to downsample the original observed signal, that is, to ensure that f s ≥2f max , and the error between the downsampled main frequency amplitude curve and the original main frequency amplitude curve is less than 1%; where f s is the sampling frequency, f max is the highest frequency of the original observed signal.

[0014] Furthermore, step S2 specifically includes:

[0015] S2.1, STL decomposes the inner loop, uses the local weighted regression algorithm to fit the periodic component, and then solves the trend component;

[0016] S2.2, STL decomposes the outer loop, obtains the residual, and calculates the robust weight in each loop, updates the weight of the local weighted regression algorithm, and uses it for the next inner loop calculation;

[0017] S2.3. When the iteration result meets the convergence condition, the decomposition is completed and the final trend component T is obtained. t ; The convergence condition is defined as:

[0018]

[0019] in, and Represents the trend component or cycle component of continuous iteration;

[0020] S2.4, the first observation signal X of the noisy station t Eliminate the trend component T t , the second observation signal of the noisy station is obtained and recorded as XN t , namely XN t =X t -T t .

[0021] Furthermore, STL decomposes the inner loop: Assume that the trend component and the period component in the kth iteration result are represented by T t k , The initial value of the trend component is set to 0, that is, T t 0 =0; the decomposition process is as follows:

[0022] S2.1.1. Detrending: Use the first observation signal X of the noisy station at the current moment t Subtract the trend component after the last inner cycle to get the detrended sequence, expressed as X t -T t k ;

[0023] S2.1.2, periodic subsequence smoothing: Divide the sequence obtained in S2.1.1 into multiple subsequences, perform local weighted regression algorithm smoothing on each subsequence, and combine and restore the smoothed subsequences. The obtained periodic subsequences are recorded as

[0024] S2.1.3. Periodic subsequences Low-pass filtering: For the periodic subsequence obtained in S2.1.2 Perform sliding average processing and then perform local weighted regression algorithm smoothing to obtain sequence;

[0025] S2.1.4. Detrending periodic subseries: using As a periodic component

[0026] S2.1.5. Deperiodicization of the original sequence: As the sequence after the deperiodic component;

[0027] S2.1.6, Trend smoothing: The sequence obtained in S2.1.5 is smoothed by a local weighted regression algorithm to obtain the trend component T t k+1 ;

[0028] STL decomposition outer loop: According to S2.1.6, the trend component T t k+1 , calculate the residual Furthermore, the robust weight ρ is calculated using the residual t , the calculation formula is as follows:

[0029] in

[0030] In the formula, R t represents the iterative residual of this round of cycle, that is, B(·) represents the bisquare weight function; median represents the median, and || represents the absolute value of the data;

[0031] Get the robust weight ρ t After that, the weights of the local weighted regression algorithm are updated in the inner loop, that is, the weights of the local weighted regression algorithm are multiplied by ρ t , completing the next iteration of the inner loop.

[0032] Furthermore, step S3 specifically includes:

[0033] S3.1. The signal variance ratio of the first observation signal of the reference station with high signal-to-noise ratio synchronous observation and the second observation signal of the noisy station is calculated by the following formula:

[0034]

[0035] Where, ROV i is the signal variance ratio corresponding to the i-th window, var represents the signal variance, e r represents the first observation signal of the reference station, el represents the second observation signal of the noisy station, (t i ,t i+1 ) represents the data segment corresponding to the i-th window, and Win_Nums is the number of windows;

[0036] S3.2. Use the density-based clustering algorithm (DBSCAN) to perform cluster analysis on the signal variance ratio obtained in S3.1 to complete the identification and positioning of the noise and extract the noisy observation signal x t .

[0037] Furthermore, step S4 specifically includes:

[0038] S4.1, parameter initialization: set the initialization FMD parameters, decomposition modulus Mode_Nums, filter length L, maximum number of iterations Max_IterNums, low-pass filter cutoff frequency f c ;

[0039] S4.2, FIR filter design: Design M Hanning window low-pass filter groups, M = 5-10, initialize the number of iterations j = 1;

[0040] S4.3, noisy observation signal x t Filter decomposition: The low-pass filter designed by S4.2 is used to decompose the noisy observation signal x t Filter and obtain the M modal components of the jth iteration Its expression is as follows:

[0041]

[0042] In the formula, represents the mth modal component in the jth iteration, represents the filter coefficient of the mth low-pass filter in the jth iteration (m=1,2,…M), and “*” represents the convolution operation;

[0043] S4.4. Filter coefficient update: Using the noisy observation signal x t , modal components and the estimated signal period Implement filter coefficient update and set j=j+1;

[0044] S4.5, iteration judgment: judge whether the current iteration has reached the maximum number of iterations; if so, execute the next step S4.6; otherwise, return to S4.3;

[0045] S4.6, Modal component selection: Calculate the correlation coefficients of all modal components and construct the correlation matrix CC K×K , select the two modal components with the largest correlation coefficient and calculate the failure period by the estimated Calculate their related kurtosis; then select the modal component with larger related kurtosis as the decomposed modal component, and set M = M-1;

[0046] S4.7, FMD decomposition result judgment: judge whether the current mode number reaches the set mode number Mode_Nums. If not, return to step S4.3; otherwise, stop iteration and obtain the final FMD decomposition result; the final low-frequency noise signal can be obtained, recorded as u t .

[0047] Furthermore, the correlation coefficient between any two modal components in step S4.6 is calculated as follows:

[0048]

[0049] In the formula, u p 、u q represents the pth and qth decomposition moduli, represents the mean of the p-th and q-th decomposition moduli, n is the sampling point, and N is the signal length.

[0050] In the second aspect, a denoising system for artificial source electromagnetic exploration based on the STL-FMD decomposition algorithm is provided, comprising the following steps:

[0051] A downsampling module is used to obtain the observation signals of the original reference station and the noisy measurement station and perform downsampling processing to obtain the first observation signals of the reference station and the noisy measurement station;

[0052] An STL decomposition module is used to perform STL decomposition on the first observation signal of the noisy measuring station and eliminate the trend component therein to obtain the second observation signal of the noisy measuring station;

[0053] A noise identification and positioning module is used to realize noise identification and positioning by using the signal variance ratio of the first observation signal of the reference station with high signal-to-noise ratio synchronous observation and the second observation signal of the noisy measurement station, and extract the noisy observation signal of the noisy measurement station;

[0054] The FMD decomposition module is used to perform low-pass FMD decomposition on the noisy observation signal of the noisy measuring station to obtain the noise signal;

[0055] A denoising module is used to subtract the noise signal from the noisy observation signal of the noisy measuring station to obtain a valid signal after denoising;

[0056] The time domain signal splicing module is used to replace the noisy part of the second observation signal of the noisy measuring station with the effective signal to obtain the complete time domain signal after denoising.

[0057] In a third aspect, an electronic terminal is provided, comprising:

[0058] a memory having a computer program stored thereon;

[0059] The processor is used to load and execute the computer program to implement the artificial source electromagnetic exploration denoising method based on the STL-FMD decomposition algorithm as described above.

[0060] In a fourth aspect, a computer-readable storage medium is provided, on which a computer program is stored. When the computer program is executed by a processor, the artificial source electromagnetic exploration denoising method based on the STL-FMD decomposition algorithm as described above is implemented.

[0061] The invention proposes a denoising method, system, terminal and medium for artificial source electromagnetic exploration based on STL-FMD decomposition algorithm. For large-scale electromagnetic interference and long-period nonlinear baseline drift in artificial source electromagnetic method (CSEM) observation data, the original data is downsampled under the premise of following the Nyquist sampling theorem and not affecting the main frequency amplitude curve, so as to eliminate part of the random interference and improve the speed of subsequent processing; the baseline drift and large-scale step noise of the observation signal are eliminated by STL decomposition; the noise is quickly and accurately identified and located by using the variance ratio of the high signal-to-noise ratio reference signal and the station signal observed synchronously, and the noisy signal is extracted and post-processed by DBSCAN as a technical means; for the small-scale noise retained after STL processing, the low-pass FMD decomposition is used to obtain the low-frequency noise signal, so as to realize the extraction of the effective signal; finally, the noise-free signal and the processed signal are spliced ​​in the time domain, and the subsequent processing can be carried out. The advantage of the present invention is that it uses the time domain signal variance ratio to achieve rapid and accurate identification and positioning of complex and uncontrollable CSEM strong electromagnetic interference. It is the first time in the field of electromagnetic exploration data processing to use the combination of STL and low-pass FMD decomposition to effectively eliminate strong electromagnetic noise and baseline drift, thereby improving the denoising effect and signal-to-noise ratio of the observation signal, and opening up a new research path for CSEM signal-noise separation. BRIEF DESCRIPTION OF THE DRAWINGS

[0062] In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the drawings required for use in the embodiments or the description of the prior art will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For ordinary technicians in this field, other drawings can be obtained based on these drawings without paying creative work.

[0063] Figure 1 It is a flow chart of a method for denoising artificial source electromagnetic exploration based on STL-FMD decomposition algorithm provided by an embodiment of the present invention;

[0064] Figure 2 It is the time domain waveform and main frequency curve of the measured original reference station and the strong interference station provided by the embodiment of the present invention; wherein (a) is the time domain waveform of the measured reference station, (b) is the time domain waveform of the strong interference station, and (c) is the main frequency curve of the measured reference station and the strong interference station;

[0065] Figure 3 1 is a comparison of the time domain waveform and main frequency curve before and after the downsampling process of the strong interference station provided by the embodiment of the present invention; wherein (a) is the time domain waveform before the downsampling process of the strong interference station, (b) is the time domain waveform after the downsampling process of the strong interference station, and (c) is the main frequency curve before and after the downsampling process of the strong interference station;

[0066] Figure 41 is a comparison of the time domain waveform and the main frequency curve before and after the STL decomposition of the strong interference station provided by the embodiment of the present invention; wherein (a) is the time domain waveform before the STL decomposition of the strong interference station, (b) is the time domain waveform of the trend component after the STL decomposition of the strong interference station, (c) is the time domain waveform after the STL decomposition of the strong interference station with the trend component removed, and (d) is the main frequency curve before and after the STL decomposition of the strong interference station;

[0067] Figure 5 is the noise identification result of the strong interference station provided by the embodiment of the present invention; wherein (a) is the DBSCAN clustering result, and (b) is the noisy observation data of the strong interference station;

[0068] Figure 6 It is a comparison of the time domain waveform and main frequency curve before and after the low-pass FMD decomposition processing of the strong interference measurement station provided by the embodiment of the present invention; wherein (a) is the time domain waveform before the low-pass FMD decomposition processing of the strong interference measurement station, (b) is the time domain waveform of the low-frequency noise signal after the low-pass FMD decomposition processing of the strong interference measurement station, (c) is the time domain waveform after the time domain signal of the strong interference measurement station is spliced, and (d) is the main frequency curve before and after the low-pass FMD decomposition processing of the strong interference measurement station. DETAILED DESCRIPTION

[0069] To make the purpose, technical solution and advantages of the present invention clearer, the technical solution of the present invention will be described in detail below. Obviously, the described embodiments are only part of the embodiments of the present invention, rather than all the embodiments. Based on the embodiments of the present invention, all other implementation methods obtained by ordinary technicians in this field without creative work belong to the scope of protection of the present invention.

[0070] As a new vibration signal processing and detection method, Feature Mode Decomposition (FMD) takes advantage of the superiority of signal correlation kurtosis and decomposes different modes through the designed adaptive finite impulse response (FIR) filter. It also takes into account the impulsiveness and periodicity of the fault signal, and the decomposition result has better robustness. The time series decomposition algorithm based on local weighted regression (Seasonal-Trend decomposition procedure using Loess, STL) is a general and robust time series decomposition method that can decompose the time series into three terms: quantity trend component (baseline drift), period component (seasonal component) and residual, which can effectively eliminate the long-period nonlinear baseline drift of the time series.

[0071] Aiming at the problem that the existing artificial source electromagnetic exploration denoising technology has poor denoising effect, the present invention applies the STL and FMD decomposition algorithms to the field of artificial source electromagnetic exploration signal processing for the first time, realizes the effective separation of signal and noise in a strong interference environment, and improves the denoising effect and signal-to-noise ratio of the observation data. The present invention is further described in detail below in conjunction with the accompanying drawings and specific embodiments, but it is not intended to limit the present invention.

[0072] The embodiment of the present invention provides a method for denoising artificial source electromagnetic exploration based on the STL-FMD decomposition algorithm, using the wide field electromagnetic method (WFEM) as the exploration method, such as Figure 1 As shown, the following steps are included:

[0073] S1. Downsampling processing of original observation signals: Obtain the observation signals of the original reference station and the noisy measurement station and perform downsampling processing to obtain the first observation signals of the reference station and the noisy measurement station.

[0074] Specifically, the Nyquist sampling theorem is used to downsample the original observed signal, that is, to ensure that f s ≥2f max , and the error between the downsampled main frequency amplitude curve and the original main frequency amplitude curve is less than 1%; where f s is the sampling frequency, f max is the highest frequency of the original observed signal.

[0075] Figure 2 The time domain waveform and main frequency curve of the measured original reference station and the strong interference station (i.e., the noisy station) provided in the embodiment of the present invention, the time domain signal of the reference station (see Figure 2 In (a), except for some pulse interference, there is no large-scale step noise and it is relatively stable. Its main frequency amplitude curve (see Figure 2 (c)) is continuous, smooth, and has no distortion points, which meets the basic requirements of the reference signal; while the station signal (see Figure 2 (b) is obviously affected by large-scale step noise such as square wave and triangle wave. The main frequency amplitude curve (see Figure 2 (c) The low frequency band (<32Hz) is severely distorted. Figure 3 The time domain waveform and main frequency curve comparison of the strong interference station before and after downsampling processing provided in the embodiment of the present application show that some random pulse interference is eliminated after downsampling processing (see Figure 3 (b)), and has no significant effect on the main frequency amplitude curve (see Figure 3 (c)).

[0076] S2, STL decomposition: Perform STL decomposition on the first observation signal of the noisy measurement station, and convert the first observation signal X t Represented as X t =T t +St +R t , where T t represents the trend component (baseline drift), S t represents the periodic component (seasonal component), R t Represents the residual, t represents the time; the trend component in the first observation signal of the noisy measuring station is eliminated to obtain the second observation signal of the noisy measuring station.

[0077] The STL decomposition includes two parts: inner loop and outer loop, which are implemented using MATLAB functions. The inner loop uses a local weighted regression algorithm to fit the periodic component S t , and then solve the trend component T t ; Obtain the residual R through the outer loop t , and calculate the robust weight in each loop for the next inner loop calculation. The specific steps include:

[0078] S2.1, STL decomposition inner loop: Assume that the trend component and cycle component in the k-th iteration result are represented by T t k , The initial value of the trend component is set to 0, that is, T t 0 =0; the decomposition process is as follows:

[0079] S2.1.1. Detrending: Use the first observation signal X of the noisy station at the current moment t Subtract the trend component after the last inner cycle to get the detrended sequence, expressed as X t -T t k ;

[0080] S2.1.2, periodic subsequence smoothing: Divide the sequence obtained in S2.1.1 into multiple subsequences, perform local weighted regression algorithm smoothing on each subsequence, and combine and restore the smoothed subsequences. The obtained periodic subsequences are recorded as

[0081] S2.1.3. Periodic subsequences Low-pass filtering: For the periodic subsequence obtained in S2.1.2 Perform sliding average processing and then perform local weighted regression algorithm smoothing to obtain sequence;

[0082] S2.1.4. Detrending periodic subseries: using As a periodic component

[0083] S2.1.5. Deperiodicization of the original sequence: As the sequence after the deperiodic component;

[0084] S2.1.6, Trend smoothing: The sequence obtained in S2.1.5 is smoothed by a local weighted regression algorithm to obtain the trend component T t k+1 ;

[0085] S2.2, STL decomposition outer loop: According to S2.1.6, the trend component T t k+1 , calculate the residual Furthermore, the robust weight ρ is calculated using the residual t , the calculation formula is as follows:

[0086]

[0087] In the formula, R t represents the iterative residual of this round of cycle, that is B(·) represents the bisquare weight function; median represents the median, and || represents the absolute value of the data;

[0088] Get the robust weight ρ t After that, the weights of the local weighted regression algorithm are updated in the inner loop, that is, the weights of the local weighted regression algorithm are multiplied by ρ t , complete the next iteration of the inner loop;

[0089] S2.3. When the iteration result meets the convergence condition, the decomposition is completed and the final trend component T is obtained. t ; The convergence condition is defined as:

[0090]

[0091] in, and Represents the trend component or cycle component of continuous iteration;

[0092] S2.4, the first observation signal X of the noisy station t Eliminate the trend component T t , and the second observation signal is recorded as XN t , namely XN t =X t -T t .

[0093] Figure 4 The time domain waveform and main frequency curve comparison of the strong interference station STL decomposition before and after the embodiment of the present application. Figure 4 (a)) is decomposed by STL to obtain its trend components (see Figure 4(b)), further, by subtracting the trend component from the downsampled data, we can obtain a processed signal that eliminates the trend and large-scale step noise (see Figure 4 In (c), by comparing the main frequency amplitude curves before and after processing, it can be seen that the low-frequency data is not significantly improved after eliminating the baseline drift and large-scale step noise (see Figure 4 In (d), further processing is required to remove low-frequency interference.

[0094] S3, noise identification and positioning: The signal variance ratio of the first observation signal of the reference station with high signal-to-noise ratio synchronous observation and the second observation signal of the noisy measurement station is used to realize rapid and accurate noise identification and positioning, and the noisy observation signal of the noisy measurement station is extracted. Specifically, the following steps are included:

[0095] S3.1. The signal variance ratio of the first observation signal of the reference station with high signal-to-noise ratio synchronous observation and the second observation signal of the noisy station is calculated by the following formula:

[0096]

[0097] Where, ROV i is the signal variance ratio corresponding to the i-th window, var represents the signal variance, e r represents the first observation signal of the reference station, el represents the second observation signal of the noisy station, (t i ,t i+1 ) represents the data segment corresponding to the i-th window, and Win_Nums is the number of windows;

[0098] S3.2. Use the density-based clustering algorithm (DBSCAN) to perform cluster analysis on the signal variance ratio obtained in S3.1. The discrete points are noise points. The noise is identified and located, and the noisy observation signal x is extracted. t .

[0099] Figure 5 This is the strong interference station noise identification result provided by the embodiment of the present application. DBSCAN clustering is performed using the variance ratio of the reference signal and the station signal (see Figure 5 In the subsequent processing, only the noisy data (see Figure 5 (b)) to preserve the original information to the greatest extent.

[0100] S4: Low-pass FDM decomposition of noisy observation signals of noisy stations: Perform low-pass FMD decomposition on the noisy observation signals to obtain the noise signal. This is implemented using MATLAB functions, the source program is improved, and the filter is set to low-pass filtering. The FMD decomposition process is as follows:

[0101] S4.1, parameter initialization: set the initialization FMD parameters, decomposition modulus Mode_Nums, filter length L, maximum number of iterations Max_IterNums, low-pass filter cutoff frequency f c .

[0102] In this embodiment, Mode_Nums=7, filter length L=25, maximum number of iterations Max_IterNums=30, low-pass filter cutoff frequency f c =100Hz.

[0103] S4.2, FIR filter design: Design M Hanning window low-pass filter groups, M = 5-10, and initialize the number of iterations j = 1. In this embodiment, M = 10.

[0104] S4.3, noisy observation signal x t Filter decomposition: A low-pass filter designed with S4.2 (cutoff frequency is f c ) for the noisy observation signal x t Filter and obtain the M modal components of the jth iteration Its expression is as follows:

[0105]

[0106] In the formula, represents the mth modal component in the jth iteration, represents the filter coefficient of the mth low-pass filter in the jth iteration (m=1, 2, ...M), and “*” represents the convolution operation.

[0107] S4.4. Filter coefficient update: Using the noisy observation signal x t , modal components and the estimated signal period Implement filter coefficient update and set j=j+1. Specifically include:

[0108] S4.4.1. Signal period estimation: For each modal component Perform Hilbert transform to obtain the analytical signal Computing analytical signals The time when the autocorrelation spectrum reaches the local maximum after passing through the zero point is the signal iteration period.

[0109] S4.4.2. Filter coefficient update: for modal components The relevant kurtosis is maximized and constrained decomposition is performed, where the relevant kurtosis is calculated as follows:

[0110]

[0111] Among them, f m =[f m (1)f m (2)…f m (L)],R XWX =x t T W H x t , R XX =x t T x t , f m represents the filter coefficient vector of the mth filter, f m (L) represents the filter coefficient of the Lth point of the mth filter.

[0112]

[0113] Among them, R XWX represents the weighted correlation matrix; R XX represents the correlation matrix; W H represents the weighting matrix; H represents the shift order of the correlation kurtosis, and the shift order is set to 1; N is the signal length; L is the filter length; T represents the conjugate transposed sign.

[0114] Furthermore, solving the generalized eigenvalue problem R XWX f m =λR XX f m , λ is its generalized eigenvalue, and the maximum value λ among the generalized eigenvalues ​​is selected max The corresponding eigenvector updates the filter coefficients.

[0115] S4.5, iteration judgment: judge whether the current iteration has reached the maximum number of iterations; if so, execute the next step S4.6; otherwise return to S4.3.

[0116] S4.6, Modal component selection: Calculate the correlation coefficients of all modal components and construct the correlation matrix CC K×K , select the two modal components with the largest correlation coefficient and calculate the failure period by the estimated Calculate their related kurtosis; then select the modal component with larger related kurtosis as the decomposed modal component, and set M = M-1.

[0117] The correlation coefficient between any two modal components is calculated as follows:

[0118]

[0119] In the formula, u p 、u q represents the pth and qth decomposition moduli, represents the mean of the p-th and q-th decomposition moduli, n is the sampling point, and N is the signal length.

[0120] S4.7, FMD decomposition result judgment: judge whether the current mode number reaches the set mode number Mode_Nums. If not, return to step S4.3; otherwise, stop iteration and obtain the final FMD decomposition result; the final low-frequency noise signal can be obtained, recorded as u t .

[0121] S5, noisy signal removal: The noisy observation signal x t Subtract the noise signal u t , get the effective signal after denoising, expressed as Y t =x t -u t .

[0122] S6. Evaluation of time domain signal splicing and processing effect: The second observation signal X t The noisy part in is replaced by the effective signal Y t , obtain the complete time domain signal after denoising. Perform FFT transformation on the processed result to evaluate the continuity of its main frequency amplitude in the frequency domain.

[0123] Figure 6 The comparison results of the time domain waveform and main frequency curve before and after the low-pass FMD decomposition processing of the strong interference station provided in the embodiment of the present application are shown in FIG. t (See Figure 6 In (a), low-pass FMD decomposition is performed to obtain its low-frequency noise u t (See Figure 6 (b)); after removing the low-frequency noise, it is further spliced ​​with the noise-free data to obtain the complete time domain signal after processing (see Figure 6 (c)); and then the main frequency amplitude curve before and after processing (see Figure 6 The signal processed by the present invention eliminates the large-scale step noise of the original observation signal of the observation station, the main frequency amplitude curve is continuous and smooth, and the quality is significantly improved, indicating that the artificial source electromagnetic exploration denoising technology method provided by the present invention is feasible and has a good denoising effect.

[0124] The above embodiment proposes a method and system for denoising artificial source electromagnetic exploration based on the STL-FMD decomposition algorithm. For large-scale electromagnetic interference and long-period nonlinear baseline drift in artificial source electromagnetic method (CSEM) observation data, the original data is downsampled under the premise of following the Nyquist sampling theorem and not affecting the main frequency amplitude curve to eliminate part of the random interference and improve the speed of subsequent processing; the baseline drift and large-scale step noise of the observation signal are eliminated by STL decomposition; the variance ratio of the high signal-to-noise ratio reference signal and the station signal observed synchronously is used to realize fast and accurate identification and positioning of noise, and the noisy signal is extracted by DBSCAN as a technical means for post-processing; for the small-scale noise retained after STL processing, the low-pass FMD decomposition is used to obtain the low-frequency noise signal, thereby realizing the extraction of the effective signal; finally, the noise-free signal and the processed signal are spliced ​​in the time domain for subsequent processing.

[0125] The embodiment of the present invention further provides an artificial source electromagnetic exploration denoising system based on the STL-FMD decomposition algorithm, comprising the following steps:

[0126] A downsampling module is used to obtain the observation signals of the original reference station and the noisy measurement station and perform downsampling processing to obtain the first observation signals of the reference station and the noisy measurement station;

[0127] An STL decomposition module is used to perform STL decomposition on the first observation signal of the noisy measuring station and eliminate the trend component therein to obtain the second observation signal of the noisy measuring station;

[0128] A noise identification and positioning module is used to realize noise identification and positioning by using the signal variance ratio of the first observation signal of the reference station with high signal-to-noise ratio synchronous observation and the second observation signal of the noisy measurement station, and extract the noisy observation signal of the noisy measurement station;

[0129] The FMD decomposition module is used to perform low-pass FMD decomposition on the noisy observation signal of the noisy measuring station to obtain the noise signal;

[0130] A denoising module is used to subtract the noise signal from the noisy observation signal of the noisy measuring station to obtain a valid signal after denoising;

[0131] The time domain signal splicing module is used to replace the noisy part of the second observation signal of the noisy measuring station with the effective signal to obtain the complete time domain signal after denoising.

[0132] It should be understood that the functional unit modules in the embodiments of the present invention may be concentrated in one processing unit, or each unit module may exist physically separately, or two or more unit modules may be integrated in one unit module, and may be implemented in the form of hardware or software.

[0133] An embodiment of the present invention further provides an electronic terminal, including:

[0134] a memory having a computer program stored thereon;

[0135] The processor is used to load and execute the computer program to implement the artificial source electromagnetic exploration denoising method based on the STL-FMD decomposition algorithm described in the above embodiment.

[0136] An embodiment of the present invention further provides a computer-readable storage medium having a computer program stored thereon. When the computer program is executed by a processor, the method for denoising artificial source electromagnetic exploration based on the STL-FMD decomposition algorithm described in the above embodiment is implemented.

[0137] It can be understood that the same or similar parts of the above embodiments can be referenced to each other, and the contents not described in detail in some embodiments can refer to the same or similar contents in other embodiments.

[0138] Those skilled in the art will appreciate that the embodiments of the present application may be provided as methods, systems, or computer program products. Therefore, the present application may adopt the form of a complete hardware embodiment, a complete software embodiment, or an embodiment in combination with software and hardware. Moreover, the present application may adopt the form of a computer program product implemented in one or more computer-usable storage media (including but not limited to disk storage, CD-ROM, optical storage, etc.) that include computer-usable program code.

[0139] The present application is described with reference to the flowcharts and / or block diagrams of the methods, devices (systems), and computer program products according to the embodiments of the present application. It should be understood that each process and / or box in the flowchart and / or block diagram, as well as the combination of the processes and / or boxes in the flowchart and / or block diagram, can be implemented by computer program instructions. These computer program instructions can be provided to a processor of a general-purpose computer, a special-purpose computer, an embedded processor, or other programmable data processing device to generate a machine, so that the instructions executed by the processor of the computer or other programmable data processing device generate instructions for implementing the processes in the flowchart and / or block diagram. Figure 1 A process or multiple processes and / or boxes Figure 1 A device that provides the functions specified in a block or multiple blocks.

[0140] These computer program instructions may also be stored in a computer-readable memory capable of directing a computer or other programmable data processing device to operate in a specific manner, so that the instructions stored in the computer-readable memory produce an article of manufacture comprising an instruction device, which implements the process Figure 1 A process or multiple processes and / or boxes Figure 1 A function specified in one or more boxes.

[0141] These computer program instructions can also be loaded onto a computer or other programmable data processing device so that a series of operating steps are executed on the computer or other programmable device to produce a computer-implemented process, thereby providing instructions for implementing the process. Figure 1 A process or multiple processes and / or boxes Figure 1 The steps for the functions specified in one or more boxes.

[0142] Although the embodiments of the present invention have been shown and described above, it is to be understood that the above embodiments are exemplary and are not to be construed as limitations of the present invention. A person skilled in the art may change, modify, replace and vary the above embodiments within the scope of the present invention.

Claims

1. A method for denoising artificial source electromagnetic exploration based on STL-FMD decomposition algorithm, characterized in that: The steps include: S1: Obtain the observation signals of the original reference station and the noisy measurement station and perform downsampling processing to obtain the first observation signals of the reference station and the noisy measurement station; S2: Perform STL decomposition on the first observation signal of the noisy station and eliminate the trend component therein to obtain the second observation signal of the noisy station; S3: noise identification and positioning are realized by using the signal variance ratio of the first observation signal of the reference station with high signal-to-noise ratio synchronous observation and the second observation signal of the noisy measurement station, and the noisy observation signal of the noisy measurement station is extracted; S4: Perform low-pass FMD decomposition on the noisy observation signal of the noisy station to obtain a noise signal; S5: subtract the noise signal from the noisy observation signal of the noisy measuring station to obtain the effective signal after denoising; S6: Replace the noisy part of the second observation signal of the noisy measuring station with the valid signal to obtain a complete time domain signal after denoising.

2. The method for denoising artificial source electromagnetic exploration based on the STL-FMD decomposition algorithm according to claim 1 is characterized in that: The downsampling process in step S1 includes: The Nyquist sampling theorem is used to downsample the original observed signal, that is, to ensure that f s ≥2f max , and the error between the downsampled main frequency amplitude curve and the original main frequency amplitude curve is less than 1%; where f s is the sampling frequency, f max is the highest frequency of the original observed signal.

3. The method for denoising artificial source electromagnetic exploration based on the STL-FMD decomposition algorithm according to claim 1 is characterized in that: Step S2 specifically includes: S2.1, STL decomposes the inner loop, uses the local weighted regression algorithm to fit the periodic component, and then solves the trend component; S2.2, STL decomposes the outer loop, obtains the residual, and calculates the robust weight in each loop, updates the weight of the local weighted regression algorithm, and uses it for the next inner loop calculation; S2.

3. When the iteration result meets the convergence condition, the decomposition is completed and the final trend component T is obtained. t ; The convergence condition is defined as: in, and Represents the trend component or cycle component of continuous iteration; S2.4, the first observation signal X of the noisy station t Eliminate the trend component T t , the second observation signal of the noisy station is obtained and recorded as XN t , namely XN t =X t -T t .

4. The method for denoising artificial source electromagnetic exploration based on the STL-FMD decomposition algorithm according to claim 3 is characterized in that: STL decomposition inner loop: Assume that the trend component and cycle component in the kth iteration result are represented by T t k , The initial value of the trend component is set to 0, that is, T t 0 =0; The decomposition process is as follows: S2.1.

1. Detrending: Use the first observation signal X of the noisy station at the current moment t Subtract the trend component after the last inner cycle to get the detrended sequence, expressed as X t -T t k ; S2.1.2, periodic subsequence smoothing: Divide the sequence obtained in S2.1.1 into multiple subsequences, perform local weighted regression algorithm smoothing on each subsequence, and combine and restore the smoothed subsequences. The obtained periodic subsequences are recorded as S2.1.

3. Periodic subsequences Low-pass filtering: For the periodic subsequence obtained in S2.1.2 Perform sliding average processing and then perform local weighted regression algorithm smoothing to obtain sequence; S2.1.

4. Detrending periodic subseries: using As a periodic component S2.1.

5. Deperiodicization of the original sequence: As the sequence after the deperiodic component; S2.1.6, Trend smoothing: The sequence obtained in S2.1.5 is smoothed by local weighted regression algorithm to obtain the trend component STL decomposition outer loop: According to S2.1.6, the trend component T t k+1 , calculate the residual Use the residual to calculate its robust weight ρ t , the calculation formula is as follows: in In the formula, R t represents the iterative residual of this round of cycle, that is B(·) represents the bisquare weight function; median represents the median, and || represents the absolute value of the data; Get the robust weight ρ t After that, the weights of the local weighted regression algorithm are updated in the inner loop, that is, the weights of the local weighted regression algorithm are multiplied by ρ t , completing the next iteration of the inner loop.

5. The method for denoising artificial source electromagnetic exploration based on STL-FMD decomposition algorithm according to claim 1, characterized in that: Step S3 specifically includes: S3.

1. The signal variance ratio of the first observation signal of the reference station with high signal-to-noise ratio synchronous observation and the second observation signal of the noisy station is calculated by the following formula: Where, ROV i is the signal variance ratio corresponding to the i-th window, var represents the signal variance, e r represents the first observation signal of the reference station, el represents the second observation signal of the noisy station, (t i ,t i+1 ) represents the data segment corresponding to the i-th window, and Win_Nums is the number of windows; S3.

2. Use density-based clustering algorithm to cluster the signal variance ratio obtained in S3.1 to complete the identification and location of noise and extract the noisy observation signal x t .

6. The method for denoising artificial source electromagnetic exploration based on STL-FMD decomposition algorithm according to claim 1, characterized in that: Step S4 specifically includes: S4.1, parameter initialization: set the initialization FMD parameters, decomposition modulus Mode_Nums, filter length L, maximum number of iterations Max_IterNums, low-pass filter cutoff frequency f c ; S4.2, FIR filter design: Design M Hanning window low-pass filter groups, M = 5-10, initialize the number of iterations j = 1; S4.3, noisy observation signal x t Filter decomposition: The low-pass filter designed by S4.2 is used to decompose the noisy observation signal x t Filter and obtain the M modal components of the jth iteration Its expression is as follows: In the formula, represents the mth modal component in the jth iteration, represents the filter coefficient of the mth low-pass filter in the jth iteration, and "*" represents the convolution operation; S4.

4. Filter coefficient update: Using the noisy observation signal x t , modal components and the estimated signal period Implement filter coefficient update and set j=j+1; S4.5, iteration judgment: judge whether the current iteration has reached the maximum number of iterations; if so, execute the next step S4.6; otherwise, return to S4.3; S4.6, Modal component selection: Calculate the correlation coefficients of all modal components and construct the correlation matrix CC K×K , select the two modal components with the largest correlation coefficient and calculate the failure period by the estimated Calculate their related kurtosis; then select the modal component with larger related kurtosis as the decomposed modal component, and set M = M-1; S4.7, FMD decomposition result judgment: judge whether the current mode number reaches the set mode number Mode_Nums. If not, return to step S4.3; otherwise, stop iteration and obtain the final FMD decomposition result; the final low-frequency noise signal can be obtained, recorded as u t .

7. The method for denoising artificial source electromagnetic exploration based on the STL-FMD decomposition algorithm according to claim 6 is characterized in that: The correlation coefficient between any two modal components in step S4.6 is calculated as follows: In the formula, u p 、u q represents the pth and qth decomposition moduli, represents the mean of the p-th and q-th decomposition moduli, n is the sampling point, and N is the signal length.

8. An artificial source electromagnetic exploration denoising system based on STL-FMD decomposition algorithm, characterized in that: The steps include: A downsampling module is used to obtain the observation signals of the original reference station and the noisy measurement station and perform downsampling processing to obtain the first observation signals of the reference station and the noisy measurement station; An STL decomposition module is used to perform STL decomposition on the first observation signal of the noisy measuring station and eliminate the trend component therein to obtain the second observation signal of the noisy measuring station; A noise identification and positioning module is used to realize noise identification and positioning by using the signal variance ratio of the first observation signal of the reference station with high signal-to-noise ratio synchronous observation and the second observation signal of the noisy measurement station, and extract the noisy observation signal of the noisy measurement station; The FMD decomposition module is used to perform low-pass FMD decomposition on the noisy observation signal of the noisy measuring station to obtain the noise signal; A denoising module is used to subtract the noise signal from the noisy observation signal of the noisy measuring station to obtain a valid signal after denoising; The time domain signal splicing module is used to replace the noisy part of the second observation signal of the noisy measuring station with the effective signal to obtain the complete time domain signal after denoising.

9. An electronic terminal, characterized in that: include: a memory having a computer program stored thereon; A processor is used to load and execute the computer program to implement the artificial source electromagnetic exploration denoising method based on the STL-FMD decomposition algorithm as described in any one of claims 1 to 7.

10. A computer-readable storage medium having a computer program stored thereon, characterized in that: When the computer program is executed by a processor, the artificial source electromagnetic exploration denoising method based on the STL-FMD decomposition algorithm as described in any one of claims 1 to 7 is implemented.

Images