CN106645947A - Time-frequency analysis method based on nonlinear mode decomposition and adaptive optimal kernel - Google Patents

Abstract

The invention discloses a time-frequency analysis method based on nonlinear mode decomposition and self-adaptive optimal kernel. This method combines the advantages of the nonlinear mode decomposition analysis method and the adaptive optimal kernel analysis method. Firstly, the nonlinear mode decomposition algorithm is used to decompose the multi-component non-stationary signal into a group of physically meaningful nonlinear mode components. The algorithm has strong noise robustness, which reduces the requirements for the anti-noise performance of subsequent analysis methods; after that, the time-frequency analysis method of adaptive optimal kernel is used to make the kernel function adaptively change with the change of the signal, effectively The cross-terms are effectively suppressed and the time-frequency focusing ability is improved. This new analysis method inherits the advantages of the nonlinear mode decomposition analysis method and the adaptive optimal kernel analysis method, and its performance is excellent.

Description

A Time-Frequency Analysis Method Based on Nonlinear Mode Decomposition and Adaptive Optimal Kernel

technical field

The invention relates to the field of time-frequency analysis of non-stationary signals, in particular to a time-frequency analysis method based on nonlinear mode decomposition and adaptive optimal kernel.

Background technique

It is self-evident that a good time-frequency analysis method plays a role in the estimation and analysis of non-stationary signals. The signal itself has many attributes, and for signal estimation, the frequency domain characteristics are very important attributes.

At present, there are many methods of signal processing and estimation. The most classic one is of course the Fourier transform, which can reflect the frequency characteristics of the signal. The Fourier transform has a good frequency domain resolution ability for stationary signals. The time information is lost, that is, the time at which each spectrum appears is not known, so the Fourier transform is not suitable for dealing with non-stationary signals. For non-stationary signals in the actual environment, commonly used linear time-frequency representation methods include short-time Fourier transform, wavelet transform, S-transform, etc., but their time-frequency accuracy needs to be improved; while Wigner-Ville distribution, Cohen-like time-frequency Although bilinear time-frequency distributions such as bilinear distributions have high time-frequency accuracy, they inevitably introduce cross terms in the calculation process.

Empirical Mode Decomposition (EMD) can decompose complex signals into a series of Intrinsic Mode Functions (IMF). Although EMD has a good performance in suppressing the cross term of multi-component signals, its disadvantages are also obvious, such as poor anti-noise performance and modal aliasing due to the presence of discontinuous signals or noise. In order to overcome the shortcomings of EMD, a noise-assisted analysis method called ensemble empirical mode decomposition (Ensemble Empirical Mode Decomposition, EEMD) has emerged. Although the anti-noise performance of EEMD is better than that of EMD, it is still not up to expectations.

Time-Frequency Representation (TFR) is a very effective method for instantaneous frequency estimation, especially for multi-component signals. General TFR methods only have fixed window functions or kernel functions, so they only have a better processing effect on a certain type of signal. Adaptive Optimal Kernel Time-Frequency Representation (AOK TFR) is a signal-based instantaneous frequency estimation method, which uses a radial Gaussian kernel function that can change as the signal changes. Therefore, AOK has better time-frequency focusing performance and suppression Cross-entry capabilities.

EMD+AOK can reduce the influence of cross items when analyzing multi-component signals, but it is too sensitive to noise. EEMD+AOK can improve the anti-noise performance to a certain extent, but it is far from enough.

Nonlinear Mode Decomposition (NMD) can decompose a signal into a series of physically meaningful nonlinear mode components, and has strong noise robustness. It is a new type of algorithm that integrates multiple methods such as time-frequency analysis, data substitution test and harmonic identification. However, the nonlinear mode component obtained by NMD decomposition is not a single component like IMF, and needs to be optimized by using the kernel function of the AOK method that can adapt to the change of the signal.

The new time-frequency representation method based on NMD and AOK not only utilizes the effective decomposition performance of NMD for multi-component signals, but also inherits the excellent time-frequency focusing performance of AOK and the ability to effectively suppress cross terms, and solves the problems of EMD+AOK and EEMD+AOK has poor anti-noise performance and still contains a certain amount of cross-terms.

Contents of the invention

Purpose of the invention: in order to solve the above-mentioned technical problems, the present invention provides a kind of time-frequency analysis method based on nonlinear mode decomposition and self-adaptive optimal kernel, adopts nonlinear mode decomposition (Nonlinear Mode Decomposition, NMD) to decompose signal to be processed into a A series of nonlinear mode components to suppress noise interference; then the obtained signal components are processed by the adaptive optimal kernel (AOK) analysis method to suppress cross-terms, because the adaptive optimal kernel analysis method (AOK) can automatically track and analyze As the signal changes, the kernel function can change adaptively according to the signal change, and the kernel function generated for different signals is always optimal, so its time-frequency analysis performance is better.

Technical solution: In order to achieve the above-mentioned technical effects, the technical solution provided by the present invention is:

A time-frequency analysis method based on nonlinear mode decomposition and adaptive optimal kernel includes steps:

(1) Define the signal to be processed as s(t), the sampling frequency of s(t) is f _s , and the data length is N; the signal to be processed s(t) is decomposed into a group of nonlinear modes by nonlinear mode decomposition method Components, namely:

Among them, c _i (t) is the i-th nonlinear mode component of s(t), and n(t) represents noise;

(2) each nonlinear mode component that step (1) decomposes adopts the self-adaptive optimal kernel analysis method to analyze, including steps:

(2-1) Construct radial Gaussian kernel function

In the formula, Indicates the radial Gaussian function used to control the radial angle Direction extension function;

(2-2) To obtain the optimal kernel that changes adaptively with the signal as the target problem, the optimization problem model is constructed as follows:

The constraints for setting the optimization problem model are:

in, is the ambiguity function of the i-th nonlinear mode component c _i (t) in polar coordinates, and β is the volume of the optimal kernel; The expression in Cartesian coordinate system is:

Among them, A _i (t, θ, τ) is The mapping quantity in Cartesian coordinate system, s ^* (t), w ^* (u) is the conjugate complex form of s(t), c _i (t), w(u); the window function w(u) is a symmetrical rhombus window function centered on t and the width is 2T, and | When u|>T, w(u)=0, variables τ and θ are parameters of the general fuzzy domain {τ, u}, |τ|<2T;

(2-3) Solve the optimization problem model according to A _i (t, θ, τ), and obtain the optimal kernel function Φ _(i)opt (t, θ, τ) of the nonlinear mode component c _i (t);

(2-4) The time-frequency expression of the adaptive optimal kernel to calculate the nonlinear mode component c _i (t) is:

In the formula, is the energy value of nonlinear mode component c _i (t) at time t;

(3) According to the time-frequency representation of the adaptive optimal kernel of all nonlinear mode components obtained in step (2), the result of calculating the time-frequency analysis is:

Further, the method for decomposing the signal s(t) to be processed into a group of nonlinear mode components by using the nonlinear mode decomposition method in the step (1) includes steps:

(1-1) Calculate the wavelet transform expression W _s (ω, t) of the signal s(t):

in, is the Fourier transform of s(t), namely s ⁺ (t) is the positive frequency part of the s(t) signal, is the wavelet function, and Each other is a Fourier transform pair, and satisfies the condition ψ ^* (t), are ψ(t), respectively, The conjugate complex number of ; ω _ψ represents the wavelet peak frequency, ω _ψ = 1, f ₀ is the resolution parameter, which is used to weigh the resolution of time and frequency in the transformation process;

(1-2) Determine whether the wavelet transform calculated in step (1-1) is the best time-frequency representation of s(t); if the judgment result is no, then calculate the windowed Fourier transform of signal s(t) Expression G _s (ω, t):

Among them, g(t) is the window function of the windowed Fourier transform, is the Fourier transform of g(t), which satisfies the following conditions:

and

(1-3) Find out all the ridge curves represented by the best time-frequency of the signal s(t), and use the ridge method to reconstruct the harmonic components, where the hth harmonic component is:

x ^(h) (t)=A ^(h) (t)cosφ ^(h) (t), h∈[1,2,…,N]

In the formula, A ^(h) (t), φ ^(h) (t) are the amplitude and phase of the hth harmonic component respectively; v ^(h) (t) is the frequency of the hth harmonic component, y ^{( h)} (t)≡φ′ ^(h) (t); φ′ ^(h) (t) is the derivative of φ ^(h) (t) to time t; N represents the highest order of the harmonic component, that is, the data length;

(1-4) Utilize anti-noise substitution test method to discriminate the authenticity of the extracted harmonic components, screen out all real harmonic components,

(1-5) adding all the real harmonic components to obtain a nonlinear mode component c ₁ (t);

(1-6) Subtract c ₁ (t) obtained by nonlinear mode decomposition from the original signal s(t), and repeat steps 1-1 to 1-6 for the residual components to obtain all the nonlinear mode components c _i (t);

Finally, the original target signal s(t) can be expressed as:

In the formula, n(t) represents noise.

Further, in the step (1-3), the method of finding all ridge curves represented by the best time-frequency of signal s(t) is:

At time t _i , find h maximum points and connect the found h maximum points to form the ridge curve at time t _i as follows:

In the above formula, i=1, 2, ..., N, N is the data length; H _s (ω, t) is the best time-frequency representation of signal s(t), which is W _s (ω, t) or G _s (ω, t); H _s (ω, t) find out all the ridge point connections at all times, which constitute N ridge curves.

Further, reconstructing the h harmonic component x ^(h) (t)=A ^(h) (t)cosφ ^(h) (t) with the ridge method in the step (1-3) includes the following steps:

If the best time-frequency representation of the signal s(t) to be processed is wavelet transform, then

If the best time-frequency representation of the signal s(t) to be processed is the windowed Fourier transform, then

In the formula, with is an improved discretized impact factor obtained by parabolic interpolation.

Further, the method for screening real harmonic components in the step (1-4) includes steps:

(1-4-1) Construct the identification statistic D(α _A , α _ν ), which is used to measure the amplitude and frequency order of each harmonic component; the expression of D(α _A , α _v ) is:

In the formula, with represent the spectral entropy of A ^(h) (t) and v ^(h) (t) respectively; α _A , α _v are weight coefficients;

(1-4-2) Create N _s Fourier transform surrogate data for the original signal s(t); wherein, the expression of the jth Fourier transform surrogate data is:

(1-4-3) Calculate the optimal time-frequency expression formula for each Fourier transform surrogate data, and extract each harmonic component from it, and calculate the identification statistics corresponding to each Fourier transform surrogate data as follows:

Define significance level indicators In the formula In order to satisfy D _s > D ₀ the number of substitute data; D ₀ is the order degree of the hth harmonic of the original signal;

Set the significant level index as p, when there are x alternative data satisfying D _s > D ₀ and x≥p×N _s , it is determined that the corresponding harmonic component is not noise; otherwise, it is determined that the corresponding harmonic component is noise;

(1-4-4) Calculate the comprehensive measure of the correlation between harmonics

in,

In the formula, w _A , w _φ , w _v represent ρ ^(h) assigns equal weights for magnitude and phase consistency ρ ^(h) ≡ρ ^(h) (1, 1, 0);

(1-4-5) Define thresholds for comprehensive metrics

(1-4-6) When a harmonic component satisfies ρ ^(h) ≥ ρ _min and the significance level index p ≥ 95%, it is determined that the harmonic component has passed the test and is a real harmonic component.

Further, the method for judging whether a substitute data satisfies D _s > D ₀ in the step (1-4-3) is:

Select any three groups of parameters (α _A , α _v ) with different values, and calculate the identification statistics corresponding to the three groups of parameters respectively. If there is any identification statistics greater than D ₀ , it is determined that the substitute data satisfies D _s >D ₀ .

Beneficial effect: compared with the prior art, the present invention has the following advantages:

The NMD method not only has a good suppression effect on Gaussian white noise, but also can effectively suppress other types of noise signals, and has strong noise robustness and wide adaptability; however, the nonlinear mode components obtained by NMD decomposition are not It is not a single component like IMF. It is necessary to use the kernel function of the AOK method to adapt to the change of the signal for optimization. The AOK method can effectively and adaptively track the change of the non-stationary signal. The new time-frequency representation method based on NMD and AOK in the present invention not only utilizes the effective decomposition performance of NMD for multi-component signals, but also inherits the excellent time-frequency focusing performance and the ability to effectively suppress cross-terms of AOK, and solves the problem of EMD+ AOK and EEMD+AOK have poor anti-noise performance and still contain a certain amount of cross-terms.

Description of drawings

Fig. 1 is a flowchart of the present invention.

detailed description

Taking a multi-component non-stationary simulation signal as an example, the implementation and advantages of the present invention will be described in detail below.

Suppose s(t) is a multicomponent signal with Gaussian white noise:

s( _t )=sd(t)+n(t)

s _d (t)＝cos(20πt)+sin(200πt)+sin(400πt)+sin(100π(t-0.5) ² )

In the formula, n(t) represents Gaussian white noise; s _d (t) is an ideal multi-component signal. Set the sampling frequency to 1kHz, the sampling time to 1s, and the data length to 1000. Set the window length 2T=128, and the kernel function volume is limited to β=5.

Flow process of the present invention is shown in Figure 1, comprises the following steps:

Step A: Prepare the signal s(t) to be processed, whose sampling frequency is f _s and data length is N;

Step B: performing NMD analysis on the signal s(t);

Step B-1: Calculate the wavelet transform (Wavelet Transform, WT) W _s (ω, t) of the signal s(t),

in, is the Fourier transform of s(t), namely s ⁺ (t) is the positive frequency part of the s(t) signal, ψ(t) is a wavelet function, and Each other is a Fourier transform pair, and satisfies the condition ψ ^* (t), are ψ(t), respectively, conjugate complex numbers of . Indicates the wavelet peak frequency.

ω _ψ = 1

Among them, f ₀ is a resolution parameter, which is used to weigh the resolution of time and frequency in the transformation process (usually the default is 1). Higher time resolution results in lower frequency resolution, and vice versa.

Step B-2: Check whether the wavelet transform W _s (ω, t) obtained in Step B-1 is the best time-frequency representation of signal s(t), if not, use Windowed Fourier Transform (Windowed Fourier Transform , WFT) G _s (ω, t), WFT is defined as:

Among them, g(t) is the window function of WFT, is the Fourier transform of g(t), satisfying the conditions: The Gaussian window is selected as the window function of WFT, and its expression is:

Step B-3 Find all the ridge curves of the best time-frequency representation of the signal s(t) here is the ridge curve of the hth harmonic. The so-called ridge curve is a curve formed by connecting some local maximum points on the time-frequency diagram. Each maximum value point is called a ridge point.

At time t _i , using the following algorithm, h maximum points can be found:

In the above formula, i=1, 2, ..., N, where N is the data length. H _s (ω, t) is the best time-frequency representation of signal s(t), that is, W _s (ω, t) or G _s (ω, t).

Connect the ridge points at all times found in H _s (ω, t) to form H ridge curves.

Step B-4: Use the ridge method to reconstruct the hth harmonic component as x ^(h) (t) = A ^(h) (t)cosφ ^(h) (t), where A ^(h) (t), φ ^(h) (t) is the magnitude and phase of the hth harmonic component respectively, ν ^(h) (t) is the frequency of the hth harmonic component, ν ^(h) (t)≡φ′ ^(h) ( t), φ′ ^(h) (t) is the derivative of φ ^(h) (t) to time t.

If the optimal time-frequency of the signal s(t) to be processed is expressed as W _s (ω, t), then the hth harmonic component x ^(h) (t) of the signal s(t) can be calculated by the following formula .

If the optimal time-frequency expression of the signal s(t) to be processed is G _s (ω, t), then the harmonic correlation component x ^(h) (t) of order h of s(t) is calculated as:

In the formula, with is an improved discretized impact factor obtained by parabolic interpolation.

Step B-5: Use the anti-noise substitution test method to identify the authenticity of the extracted harmonic components, screen out all the real harmonic components, and stop the decomposition when three consecutive harmonic components are judged to be false process. Include the following steps:

(1) Extract the harmonic components from TFR and calculate the corresponding identification statistics D(α _A , α _v );

The order degree of the extracted amplitude A ^(h) (t) and frequency ν ^(h) (t) of each harmonic component can be calculated by its spectral entropy with To measure quantitatively, the identification statistics D and spectral entropy Q are defined as follows:

Wherein, α _A , α _v are weight coefficients for calculating D(α _A , α _v ).

(2) Create N _s Fourier transform substitute data for the original signal s(t);

make The amplitude of remains unchanged, and the phase becomes N _s random phases uniformly distributed on [0, 2π) The inverse Fourier transform corresponding to each phase of this random distribution is a Fourier transform surrogate data s′ _j (t) of s(t).

Here j=1, 2, . . . , N _s .

(3) Calculate the TFR associated with each surrogate data, and extract each harmonic component from it, so as to calculate the identification statistics corresponding to each surrogate data

Define significance level indicators In the formula is the number of substitute data where D _s >D ₀ ; D ₀ is the order degree of the hth harmonic of the original signal.

Assuming that N _s surrogate data are created and the significance level index is set to p, that is, there are at least p×N _s surrogate data satisfying D _s > D ₀ before the component is considered not to be noise, and the decomposition process continues.

The method for judging whether a substitute data satisfies D _s > D ₀ is as follows: select any three groups of parameters (α _A , α _v ) with different values, and calculate the corresponding identification statistics of the three groups of parameters. If there is any identification statistics If the value is greater than D ₀ , it is determined that the substitute data satisfies D _s >D ₀ .

(4) Calculate the comprehensive measure of the correlation between harmonics

in,

In the formula, w _A , w _φ , w _v represent , where ρ ^(h) ≡ ρ (h ^{) (} 1, 1, 0) is used by default to assign equal weights to amplitude and phase consistency, and not to assign weights to frequency consistency.

(5) In order to reduce the misjudgment of the real harmonic component, define the threshold of the comprehensive measurement value

(6) When a harmonic component satisfies ρ ^(h) ≥ ρ _min and the significance level index p ≥ 95%, it is considered that this harmonic component has passed the test and is a real harmonic component.

Step B-6: Add all real harmonic components to obtain a nonlinear mode component c ₁ (t).

Step B-7: Subtract c ₁ (t) obtained by nonlinear mode decomposition from the original signal s(t), and repeat steps B-1 to B-6 for the residual components to obtain all the individual nonlinear mode components c _i (t).

Finally, the original target signal s(t) can be expressed as:

In the formula, n(t) represents noise.

Step C: Perform AOK analysis on each nonlinear mode component c _i (t) decomposed by the NMD method, the process is as follows:

Step C-1: First select the radial Gaussian kernel function, which is defined as follows:

In the formula, Indicates the radial Gaussian function used to control the radial angle Direction extension function.

Step C-2: In order to obtain the optimal kernel that changes adaptively with the signal, the following constrained optimization problem should be solved.

The constraints are:

in, is the ambiguity function in polar coordinates, For each nonlinear mode component corresponding to β is the volume of the optimal kernel.

The definition A _i (t; θ, τ) in the Cartesian coordinate system is:

where, s ^* (t), w ^* (u) is the conjugate complex form of s(t), c _i (t), w(u); the window function w(u) is a symmetrical rhombus window function centered on t and the width is 2T, and | When u|>T, w(u)=0; variables τ and θ are parameters of the general fuzzy domain {τ, u}, and |τ|<2T.

Step C-3: put Substituting (4) and (5) in step C-2, an optimal kernel function can be obtained by solving this constrained optimization problem It varies with time like the short-term fuzzy function.

Step C-4: Calculate the time-frequency representation AOK TFR _i of the adaptive optimal kernel of a certain nonlinear pattern component c _i (t) at the current time point (time t):

In the formula, is the energy value of a certain component c _i (t) of the signal at a certain moment t.

Step C-5: Repeat steps C-1 to C-4 to obtain adaptive optimal kernel time-frequency representations of all decomposed nonlinear mode components.

Step D: The obtained AOK TFR _i of all signal components (ie ) to get the final time-frequency analysis result P _NMD-AOK (t, ω):

From the above formula, the result of NMD+AOK time-frequency analysis can be finally obtained.

The technical scheme provided by the present invention is compared with prior art, can obtain following analysis result:

Only use the AOK method to produce a certain amount of cross-terms, and the anti-noise performance is very poor; EMD+AOK and EEMD+AOK methods can suppress the cross-terms to a certain extent, but are still too sensitive to noise; and the NMD adopted in the present invention The +AOK method is obviously superior to several other methods in terms of removing noise and suppressing cross-terms.

The above is only a preferred embodiment of the present invention, it should be pointed out that for those of ordinary skill in the art, without departing from the principle of the present invention, some improvements and modifications can also be made. It should be regarded as the protection scope of the present invention.

Claims (6)

1. A time-frequency analysis method based on nonlinear mode decomposition and self-adaptive optimal kernel is characterized by comprising the following steps:

(1) defining the signal to be processed as s (t), the sampling frequency of s (t) is f_sThe data length is N; decomposing the signal s (t) to be processed into a group of nonlinear mode components by adopting a nonlinear mode decomposition method, namely:

$s\left(t\right)=\underset{i}{\&Sigma;}{c}_i\left(t\right)+n\left(t\right)$

wherein, c_i(t) is the ith nonlinear mode component of s (t), and n (t) represents noise;

(2) analyzing each nonlinear mode component decomposed in the step (1) by adopting a self-adaptive optimal kernel analysis method, wherein the method comprises the following steps:

(2-1) construction of radial Gaussian Kernel function

In the formula,representing radial angles for controlling radial Gaussian functionsA spread function of direction;

(2-2) taking the optimal kernel which changes along with the signal self-adaption as a target problem, and constructing an optimization problem model as follows:

the constraint conditions for setting the optimization problem model are as follows:

wherein,is the ith nonlinear mode component c_i(t) ambiguity function in polar coordinates, β being the volume of the optimal kernel;the expression in the rectangular coordinate system is:

$A_i(t,\mathrm{\&theta;},\mathrm{\&tau;})={\&Integral;}_{-\mathrm{\&infin;}}^{+\mathrm{\&infin;}}{c}_i^{*}(u-\frac{\mathrm{\&tau;}}{2})w^{*}(u-t-\frac{\mathrm{\&tau;}}{2})\&CenterDot;c_i(u-\frac{\mathrm{\&tau;}}{2})w(u-t-\frac{\mathrm{\&tau;}}{2})e^{j\mathrm{\&theta;}u}du$

wherein A is_i(t, θ, τ) isAmount of mapping in rectangular coordinate system, s^*(t)，w^*(u) is s (t), c_i(t), complex conjugate forms of w (u); the window function w (u) is a symmetrical diamond window function with the center of T and the width of 2T, when | u | is greater than T, w (u) is 0, variables tau and theta are parameters of a general fuzzy domain { tau, u }, and | tau | is less than 2T;

(2-3) according to A_i(t, theta, tau) solving the optimization problem model to obtain a nonlinear mode component c_i(t) optimum kernel function Φ_(i)opt(t，θ，τ)；

(2-4) calculating the nonlinear mode component c_i(t) the time-frequency representation of the adaptive optimal kernel is:

$P_{\left(i\right)AOK}(t,\mathrm{\&omega;})=\frac{1}{4\mathrm{\&pi;}}{\&Integral;}_{-\mathrm{\&infin;}}^{+\mathrm{\&infin;}}{\&Integral;}_{-\mathrm{\&infin;}}^{+\mathrm{\&infin;}}{A}_i(t,\mathrm{\&theta;},\mathrm{\&tau;}){\mathrm{\&Phi;}}_{\left(i\right)opt}(t,\mathrm{\&theta;},\mathrm{\&tau;})e^{-j\mathrm{\&theta;}t-j\mathrm{\&tau;}\mathrm{\&omega;}}d\mathrm{\&theta;}d\mathrm{\&tau;}$

in the formula, P_(i)AOK(t, ω) is a nonlinear mode component c_i(t) an energy value at time t;

(3) according to the time-frequency representation of the self-adaptive optimal kernels of all the nonlinear mode components obtained in the step (2), calculating the result of time-frequency analysis as follows:

$P_{NMD-AOK}(t,\mathrm{\&omega;})=\underset{i}{\&Sigma;}{P}_{{\left(i\right)}_{AOK}}(t,\mathrm{\&omega;}).$

2. the time-frequency analysis method based on nonlinear mode decomposition and adaptive optimal kernel as claimed in claim 1, wherein the method for decomposing the signal to be processed s (t) into a set of nonlinear mode components by using nonlinear mode decomposition in step (1) comprises the steps of:

(1-1) wavelet transform expression W of calculation signal s (t)_s(ω，t)：

$W_s(\mathrm{\&omega;},t)\&equiv;{\&Integral;}_{-\mathrm{\&infin;}}^{+\mathrm{\&infin;}}{s}^{+}\left(u\right){\mathrm{\&psi;}}^{*}\&lsqb;\frac{\mathrm{\&omega;}(u-t)}{{\mathrm{\&omega;}}_{\mathrm{\&psi;}}}\&rsqb;\frac{\mathrm{\&omega;}du}{{\mathrm{\&omega;}}_{\mathrm{\&psi;}}}=\frac{1}{2\mathrm{\&pi;}}{\&Integral;}_0^{\mathrm{\&infin;}}{e}^{i\mathrm{\&xi;}t}\hat{s}\left(\mathrm{\&xi;}\right){\widehat{\mathrm{\&psi;}}}^{*}\left(\frac{{\mathrm{\&omega;}}_{\mathrm{\&psi;}}\mathrm{\&xi;}}{\mathrm{\&omega;}}\right)d\mathrm{\&xi;}$

Wherein,is the Fourier transform of s (t), i.e.s⁺(t) is the positive frequency portion of the s (t) signal,ψ (t) is a wavelet function, anAre mutually Fourier transform pairs and satisfy the condition Respectively by psi (t),the conjugate complex number of (a); omega_ψWhich represents the peak frequency of the wavelet, ω_ψ＝1，f₀the resolution parameter is used for balancing the time and frequency resolution in the transformation process;

(1-2) judging whether the wavelet transform obtained by the calculation in the step (1-1) is the optimal time-frequency representation of s (t); if the judgment result is no, calculating a windowing Fourier transform expression G of the signal s (t)_s(ω，t)：

$G_s(\mathrm{\&omega;},t)\&equiv;\underset{-\mathrm{\&infin;}}{\overset{\mathrm{\&infin;}}{\&Integral;}}{s}^{+}\left(u\right)g(u-t)e^{-i\mathrm{\&omega;}(u-t)}dt=\frac{1}{2\mathrm{\&pi;}}\underset{0}{\overset{\mathrm{\&infin;}}{\&Integral;}}{e}^{i\mathrm{\&xi;}t}\hat{s}\left(\mathrm{\&xi;}\right)\hat{g}(\mathrm{\&omega;}-\mathrm{\&xi;})d\mathrm{\&xi;}$

Wherein g (t) is a windowed function of a windowed Fourier transform,is the Fourier transform of g (t), and satisfies the following conditions:

and is

(1-3) finding all ridge curves represented by the optimal time frequency of a signal s (t), and reconstructing harmonic components by using a ridge method, wherein the h-th harmonic component is as follows:

x^(h)(t)＝A^(h)(t)cosφ^(h)(t)，h∈[1，2，…，N]

in the formula, A^(h)(t)、φ^(h)(t) is the amplitude and phase of the h-th harmonic component, respectively; v is^(h)(t) is the frequency of the h-th harmonic component, v^(h)(t)≡φ^′(h)(t)；φ^′(h)(t) is phi^(h)(t) a derivative over time t; n represents the highest order of the harmonic component, i.e., the data length;

(1-4) identifying the truth of the extracted harmonic component by using an anti-noise alternative inspection method, screening out all real harmonic components,

(1-5) adding all the real harmonic components to obtain a nonlinear mode component c₁(t)；

(1-6) subtracting c obtained by nonlinear mode decomposition from the original signal s (t)₁(t) and repeating steps 1-1 to 1-6 on the residual component to obtain all the respective nonlinear mode components c_i(t)；

Finally, the original target signal s (t) can be expressed as:

$s\left(t\right)=\underset{i}{\&Sigma;}{c}_i\left(t\right)+n\left(t\right)$

in the formula, n (t) represents noise.

3. The time-frequency analysis method based on nonlinear mode decomposition and adaptive optimal kernel as claimed in claim 2, wherein all ridge curves of optimal time-frequency representation of signal s (t) are found in the steps (1-3)The method comprises the following steps:

at time t_iFinding h maximum value points, and connecting the h maximum value points to form a moment t_iThe ridge curves of (a) are:

${\mathrm{\&omega;}}_p\left(t\right)=\underset{\mathrm{\&omega;}\&Element;\&lsqb;{\mathrm{\&omega;}}_{-}\left(t_i\right),{\mathrm{\&omega;}}_{+}\left(t_i\right)\&rsqb;}{\arg max}\left|H_s(\mathrm{\&omega;},t)\right|$

in the above formula, i is 1, 2, …, N is data length; h_s(ω, t) is the best time-frequency representation of signal s (t), i.e. W_s(ω, t) or G_s(ω，t)；H_sAnd (omega, t) finding out ridge point connecting lines at all the time points, namely forming N ridge curves.

4. The time-frequency analysis method based on nonlinear mode decomposition and adaptive optimal kernel as claimed in claim 3, wherein the h-th harmonic component x is reconstructed by ridge method in the steps (1-3)^(h)(t)＝A^(h)(t)cosφ^(h)(t) comprises the steps of:

if the best time-frequency representation of the signal s (t) to be processed is a wavelet transform, then

If the best time-frequency representation of the signal to be processed s (t) is a windowed Fourier transform, then

$\left\{\begin{array}{c}{v}^{\left(h\right)}\left(t\right)={\mathrm{\&omega;}}_p^{\left(h\right)}\left(t\right)+\mathrm{\&delta;}{v}_d^{\left(h\right)}\left(t\right)\\ A^{\left(h\right)}{e}^{{\mathrm{i\&phi;}}^{\left(h\right)}\left(t\right)}=\frac{2G_s\left({\mathrm{\&omega;}}_p^{\left(h\right)}\right(t),t)}{\hat{g}\&lsqb;{\mathrm{\&omega;}}_p^{\left(h\right)}\left(t\right)-v^{\left(h\right)}\left(t\right)\&rsqb;}\end{array}\right.$

In the formula,seed of a plantIs an improved discretization impact factor obtained by parabolic interpolation.

5. The time-frequency analysis method based on nonlinear mode decomposition and adaptive optimal kernel as claimed in claim 4, wherein the method for screening true harmonic components in the steps (1-4) comprises the steps of:

(1-4-1) construction of the recognition statistic D (α)_A，α_ν) For measuring the order of magnitude and frequency of the respective harmonic components, D (α)_A，α_v) The expression of (a) is:

$D({\mathrm{\&alpha;}}_A,{\mathrm{\&alpha;}}_v)\&equiv;{\mathrm{\&alpha;}}_{A}Q\&lsqb;{\hat{A}}^{\left(h\right)}\left(\mathrm{\&xi;}\right)\&rsqb;+{\mathrm{\&alpha;}}_{v}Q\&lsqb;{\hat{v}}^{\left(h\right)}\left(\mathrm{\&xi;}\right)\&rsqb;,$

$Q\&lsqb;f\left(x\right)\&rsqb;\&equiv;-\&Integral;\frac{|f\left(x\right){|}^2}{\&Integral;|f\left(x\right){|}^{2}dx}\ln \frac{|f\left(x\right){|}^2}{\&Integral;|f\left(x\right){|}^{2}dx}dx$

in the formula,andrespectively represent A^(h)(t) and v^(h)(t) spectral entropy α_A，α_vIs a weight coefficient;

(1-4-2) creation of N for the original signal s (t)_sFourier transform the substitution data; wherein, the expression of the jth Fourier transform substitution data is as follows:

$s_j^{\&prime;}\left(t\right)=\frac{1}{2\mathrm{\&pi;}}\&Integral;e^{i\mathrm{\&xi;}t}\left|\hat{s}\left(\mathrm{\&xi;}\right)\right|e^{{\mathrm{i\&phi;}}_{{\mathrm{\&xi;}}_j}}d\mathrm{\&xi;},j=1,2,...,N_s;$

(1-4-3) calculating an optimal time-frequency expression formula of each Fourier transform alternative data, extracting each subharmonic component from the optimal time-frequency expression formula, and calculating identification statistics corresponding to each Fourier transform alternative data as follows:

$D_j({\mathrm{\&alpha;}}_A,{\mathrm{\&alpha;}}_v)\&equiv;{\mathrm{\&alpha;}}_{A}Q\&lsqb;{\hat{A}}^{\left(j\right)}\left(\mathrm{\&xi;}\right)\&rsqb;+{\mathrm{\&alpha;}}_{v}Q\&lsqb;{\hat{v}}^{\left(j\right)}\left(\mathrm{\&xi;}\right)\&rsqb;,j=1,2,...,N_s$

defining a significance level indicatorIn the formulaTo satisfy D_s＞D₀The number of substitute data of (a); d₀The degree of order of h harmonic of the original signal;

setting the significance level index as p, when presentx number of alternative data satisfy D_s＞D₀And x is more than or equal to p × N_sThen, the corresponding harmonic component is judged not to be noise; otherwise, judging that the corresponding harmonic component is noise;

(1-4-4) calculating a comprehensive metric of correlation between harmonics

${\mathrm{\&rho;}}^{\left(h\right)}(w_A,w_{\mathrm{\&phi;}},w_v)={\left(q_A^{\left(h\right)}\right)}^{w_A}{\left(q_{\mathrm{\&phi;}}^{\left(h\right)}\right)}^{w_{\mathrm{\&phi;}}}{\left(q_v^{\left(h\right)}\right)}^{w_v}$

Wherein,

$q_A^{\left(h\right)}\&equiv;\exp \{-\frac{\sqrt{<{\&lsqb;A^{\left(h\right)}\left(t\right)<A^{\left(1\right)}\left(t\right)>-A^{\left(1\right)}\left(t\right)<A^{\left(h\right)}\left(t\right)>\&rsqb;}^2>}}{<A^{\left(1\right)}\left(t\right)A^{\left(h\right)}\left(t\right)>}\},$

$q_{\mathrm{\&phi;}}^{\left(h\right)}\&equiv;a|<\exp \{i\&lsqb;{\mathrm{\&phi;}}^{\left(h\right)}\left(t\right)-{\mathrm{h\&phi;}}^{\left(1\right)}\left(t\right)\&rsqb;\}>|,$

$q_v^{\left(h\right)}\&equiv;\exp \{-\frac{\sqrt{<{\&lsqb;v^{\left(h\right)}\left(t\right)-{\mathrm{hv}}^{\left(1\right)}\left(t\right)\&rsqb;}^2>}}{<v^{\left(h\right)}\left(t\right)>}\}.$

in the formula, w_A，w_φ，w_vRepresentsWeight of rho^(h)Equal weight rho is distributed for amplitude and phase consistency^(h)≡ρ^(h)(1，1，0)；

(1-4-5) defining a threshold value for the composite metric value

(1-4-6) when a harmonic component satisfies ρ^(h)≥ρ_minAnd when the significance level index p is 95%, the harmonic component is determined to be the true harmonic component through the inspection.

6. The time-frequency analysis method based on nonlinear mode decomposition and adaptive optimal kernel as claimed in claim 5, characterized in thatJudging whether a substitute data satisfies D in the above step (1-4-3)_s＞D₀The method comprises the following steps:

any three sets of parameters with different values are selected (α)_A，α_v) Respectively calculating the identification statistic values corresponding to the three groups of parameters, and if any one identification statistic value is larger than D₀Then the substitute data is judged to satisfy D_s＞D₀。