Research on a Signal Separation Method Based on Vold-Kalman Filter of Improved Adaptive Instantaneous Frequency Estimation

The fault vibration signal of rotating machinery system under strong background noise has the characteristics of non-stationary, non-Gaussian and complex components. In view of these characteristics, an improved method of signal separation based on Vold-Kalman filter (VKF) of adaptive instantaneous frequency estimation is proposed. First, a method for adaptive multiridge extraction of peaks detection based on synchro-squeezing wavelet transform (SWT) is proposed as the high-precision adaptive instantaneous frequency (IF) estimation method. The high precision IF estimation is used as the instantaneous frequency parameter of VKF, so that the complex multi-component non-stationary signal can be separated directly in the time domain and transformed into a signal combination composed of multiple stationary single-component signals and signal residues. Secondly, an improved method is proposed combining the adaptive IF estimation method with order tracking analysis and diagonal slice of bispectrum. In the improved method, the corresponding IF estimation of each component signal is taken as the reference frequency of its order tracking and the order spectrum analysis of each component signal is carried out respectively. Meanwhile, the signal residual is analyzed by diagonal slice of bispectrum, so as to suppress Gaussian noise and effectively separate and extract fault features in the vibration signal. Finally, the method is verified on simulation data and experimental data under different conditions. The results show that the improved method has higher extraction accuracy than other traditional methods. It has the superiority and the great potential for practical applications.


I. INTRODUCTION
For rotating machinery, periodic shock vibration will occur at the fault location due to its own rotation characteristics, so the fault vibration signal is non-stationary [1]. At the same time, due to the impact of periodic shock, the vibration system often shows non-linearity, so that the output signal of the system has non-Gaussian characteristics [2]. For a rotating system with rotation speed variation, the vibration signal itself is non-stationary due to changes in speed, and its dynamic modulation sideband will vary with the speed [3]. For complex multi-component signals The associate editor coordinating the review of this manuscript and approving it for publication was Chuan Li. generated by complex rotating systems, such as a multi-stage transmission gearbox or a complex gearbox with multiple inputs and outputs, there are multiple fault frequencies, and the accurate extraction of the relevant features of each fault frequency lays a good foundation for its fault diagnosis.
Traditional fault feature extraction methods often treat the vibration signal as a locally approximately stationary Gaussian signal processing in a short period of time, but the traditional signal processing methods are weak when faced with the multi-component complex components, nonstationary and non-Gaussian characteristics of fault signals in rotating machinery systems. Therefore, according to the fault characteristics of mechanical equipment, how to extract effective fault features from non-stationary and non-Gaussian complex multi-component fault signals has become the focus of current research.
In recent years, many scholars have studied the field of the non-stationary, non-Gaussian and complex components of fault signals in rotating machinery systems. For example, Cheng et al. [4] successfully applied local mean decomposition (LMD) and order tracking to the fault diagnosis of gears; Yang et al. [5] further optimized the parameterized time-frequency analysis method and successfully applied it to the signal separation of multi-component FM signals; Hong et al. [6] proposed a novel vibration-based fault diagnostic algorithm for gearboxes under speed fluctuations without rotational speed measurement; Wang et al. [7] proposed an adaptive noise reduction method combined with sparse wavelet coefficients, and successfully recovered and extracted weak bearing fault features from multi-component signals; Wang et al. [8] successfully applied to the fault diagnosis of gearboxes by optimizing the parameters of variational modal decomposition (VMD); Hou et al. [9] proposed a rotation speed order tracking method based on generalized demodulation (GD) and applied it to wavelet transform fault detection; Feng et al. [10] successfully applied the time-frequency demodulation analysis method to the fault diagnosis of wind turbine planetary gearboxes at nonstationary speeds; Huang et al. [11] proposed a method of multi-frequency curve extraction based on the fast path and applied it to the bearing fault diagnosis under variable rotation speed; Li et al. [12] proposed a method combining Vold-Kalman filtering with compound multi-scale fuzzy entropy and Laplace operator. It was used for fault diagnosis of rolling bearings.
Although many scholars have studied and made positive contributions to this field, the scope of research characteristics targeted is not comprehensive. Because of limitations of the methods, the effect of vibration signal processing is not satisfactory. Therefore, it is necessary to take into account the characteristics of non-stationary, non-Gaussian and complex components of fault signals of rotating machinery system, and study the solutions from various aspects.
High-order cumulates and corresponding high-order spectra [13] are the main methods for processing non-Gaussian signals. For example, Zhou et al. [14] has deeply studied characteristics of bispectrum in gearbox vibration signals, which provides the reference direction about the application of bispectrum to gearbox fault diagnosis. Chen and Yang [15] deeply studied various algorithms of highorder spectra and made great contributions to the promotion and application of high-order spectra.
In view of the non-stationary characteristic, scholars have improved a series of time-frequency analysis methods (such as short-time Fourier transform [16], Wigner-Ville distribution [17] and Hilbert-Huang transform [18]) or time-scale analysis (such as wavelet transform and wavelet packet transform [19]). Although these methods have achieved good results in fault diagnosis of non-stationary signals, there are some deficiencies. For example, the time-frequency focusing of the short-time Fourier transform is poor, cross-interference terms appear in the Wigner-Ville distribution, and the accuracy of the wavelet transform or wavelet packet transform depends on the choice of the wavelet basis, etc.
In addition, order tracking analysis [20] is also an effective method for processing non-stationary signals. The order tracking signal is an equal-angle interval sampling signal, and the number of sampling points is guaranteed to be the same in each vibration period. When the reference shaft is selected, the order spectrum distribution is not affected by the change of the reference shaft speed. Moreover, this method can ensure that the entire period of the periodic signal is truncated, thereby highlighting the periodic component of the signal. For a single input gearbox, speed can be used as the reference frequency, but for a multi-input gearbox or machine, since the input speed is not unique, speed cannot be used as the reference frequency.
Signal separation [21] is an indispensable link in the process of multi-component signal processing, which lays a good foundation for the subsequent signal parameter estimation. How to accurately separate each component from the multi-component signal is a difficult problem in the multi-component signal analysis. According to different application fields and analysis methods, many algorithms suitable for signal separation are proposed, such as polynomial fitting method of signal components [22], single frequency assumption method of signal components [23], independent components Analysis method (ICA) [24], time-varying filter method [25], time-frequency analysis method [26], etc. The Vold-Kalman filtering method [27] can effectively decompose the frequency components of complex multi-component signals on the basis of minimizing the structure and data equation errors. It provides an effective means for solving the problem about the single-component decomposition of complex multi-component vibration signals.
Compared with traditional filtering methods, Vold-Kalman filtering directly extracts the signal components of interest from the time domain. It can avoid the phase deviation caused by the time-to-frequency transformation. Compared with the empirical mode decomposition method (EMD), the center frequency of the Vold-Kalman filter can be adaptively adjusted according to the instantaneous frequency. It can effectively separate signal components that are adjacent or even crossed in the time-frequency domain [28]. It can avoid the pattern aliasing caused by EMD.
However, VKF is based on instantaneous frequency estimation, and peaks detection is the most commonly used traditional instantaneous frequency estimation method for multi-component signals. At present, the commonly used peaks detection method is time-frequency transform peak searching method [29], that is, peak searching is carried out after time-frequency transformation of the signal. Since the size of the time-frequency coefficient represents the similarity between the current position of the signal and the time-frequency transformation function, a local maximum value will be formed at the spectral peak of the time-frequency coefficient, and a convex ridge will be formed vertically. The basic principle of time-frequency transform peak searching method [30] is to determine the position of peak by searching the ridge line in the time-frequency coefficient matrix, so as to complete the peaks detection.
By summarizing the research methods of complex multicomponent vibration signals, it is found that the problems to be solved are concentrated in three aspects. The first is how to obtain a high-precision instantaneous frequency estimate. The second is how to separate signals efficiently. The third is how to extract the characteristics of the separated signals effectively. At the same time, due to the non-stationary, non-Gaussian and complex components of the fault signal, the difficulty of achieving the effective extraction of fault features is increased.
For the above Problem, this paper has proposed a signal separation method based on Vold-Kalman filtering of improved adaptive instantaneous frequency estimation. The main contributions are as follows: (1) The method for adaptive multiridge extraction of peaks detection based on synchro-squeezing wavelet transform (SWT) is proposed as a method of high-precision adaptive instantaneous frequency (IF) estimation. It achieves high precision and fast adaptive IF estimation.
(2)It is proposed to use high-precision IF estimates as the instantaneous frequency parameters of Vold-Kalman filtering (VKF), so that the complex multi-component non-stationary signal can be separated directly in the time domain and transformed into a signal combination composed of multiple stationary single-component signals and signal residues.
(3)An improved method is proposed combining the adaptive IF estimation methods with order tracking and diagonal slice of bispectrum. In the improved method, the corresponding IF estimation of each component signal is taken as the reference frequency of its order tracking and the order spectrum analysis of each component signal is carried out respectively. Meanwhile, the signal residual is analyzed by diagonal slice of bispectrum, so as to suppress Gaussian noise and effectively separate and extract fault features in the vibration signal.
(4)It is proved that the improved method proposed in this paper is superior to the traditional speed order tracking analysis and traditional signal separation methods in suppressing noise and extracting signal fault features.
(5) It proves the superiority and practicability of the improved method proposed in this paper in solving the problem of fault feature extraction of complex multi-component non-stationary signals under the background of strong noise.
In the second section of this paper, the principle and algorithm flow involved in the improved method are introduced. In the third section, compared with the traditional spectral peak detection method, the advantages of the improved method in high-precision multi-ridge extraction were verified by simulation signals. Compared with EMD and traditional order tracking analysis, the feasibility and superiority of the improved method in signal separation and feature extraction are verified by simulation signals. In the fourth section, different experimental data are used to verify the practicability of the improved method. Finally, the work of this paper was summarized and prospected in the fifth section, which laid the foundation for the subsequent research work.

A. PRINCIPLE OF SIGNAL SEPARATION BASED ON VOLD-KALMAN FILTERING OF ADAPTIVE INSTANTANEOUS FREQUENCY ESTIMATION
In this paper, the method for adaptive multiridge extraction of peaks detection based on SWT is proposed as the method of adaptive instantaneous frequency estimation. The number and starting points of ridge extraction, conditions of the extraction and conditions of screening curves, etc., can be adaptively determined. This method can efficiently and accurately extract the time-frequency ridge of multi-component complex signals with high precision.
By reassigning the wavelet transform coefficients of the signal, Synchro-squeezing Wavelet Transform (SWT) [31] can concentrate the energy in the time-frequency spectrum of the signal to the vicinity of the real instantaneous frequency of the signal, so as to eliminate interference terms and improve the time-frequency resolution [32]. Compared with traditional wavelet transform, SWT has higher timefrequency resolution and time-frequency aggregation.
SWT is a time-frequency recombination analysis algorithm based on CWT. First, the continuous wavelet transform is applied to the signalx(t): where: a is the scale factor; b is the translation factor; t is time; ψ * (t) is the conjugate of the wavelet basis function ψ(t) . When the signal is a pure harmonic function, namely x(t) = A cos(t), and the frequency is ξ < 0, then the Fourier transform isψ(ξ ) = 0. According to Plancherel's theorem, its continuous wavelet coefficient is: If the frequency of the mother waveletψ(ξ ) is ξ = ω 0 , then the wavelet coefficient W x (a, b) will aggregate in the time scale plane a = ω 0 /ω. However, the actual wavelet coefficients diverge in the scale direction, which makes the time spectrum of wavelet transform relatively fuzzy. It is pointed out that although W x (a, b) is distributed near a, its oscillation behavior at point b is determined by the original frequency ω, which is independent of the value of a [33]. Therefore, for any (a, b) and W x (a, b) = 0, calculate its instantaneous frequency ω x (a, b) is: Finally, for the signal x(t), setting up the mapping relationship (a, b) → (ω x (a, b) , b), the wavelet coefficients W x (a, b) transforms from time-scale plane to time-frequency plane. The wavelet coefficients of the interval ω l − 1 2 ω, ω l + 1 2 ω around the arbitrary center frequency ω l , which belong to W x (ω x (a, b) , b), are compressed by SWT on time-frequency plane. In this way, each frequency component of the signal is compressed in the direction of the frequency domain. The disadvantage of the traditional wavelet transform with low resolution in the time-frequency domain is changed, so that each frequency component of the signal can be clearly displayed on the time-frequency graph. The discrete calculation formula of the SWT coefficient is: The principle of Vold-Kalman filtering (VKF) is quoted from Vold's improved method of Kalman filter [19]. This method can directly extract the signal components in the time domain. In this paper, the IF estimation is used as the reference instantaneous frequency parameter of each component, and the second-order Vold-Kalman filter is selected extract multiple components simultaneously.
The modulation signal can be expressed as: The amplitude envelope can be expressed as a low-order polynomial. For discrete signals, the amplitude envelope can be expressed as: where: s is the difference order; ∇ is a difference operator; ε k (n) is the non-homogeneous item. Set s = 2, according to equation (5): . For the actual signal,A k (0) = 0. When n = N , for the actual signal, which is where, δ(n) is the noise or error term. The matrix form is: The time-frequency ridgeline extracted by the instantaneous frequency estimation method can estimate ω k (n), the instantaneous frequency of the signal component of interest, to obtain its carrier matrix B. The reconstruction form of each component signal is: By resampling the discrete signals in time domain at equal angular intervals, the order tracking algorithm skillfully transforms the non-stationary signals into stationary signals. The traditional order tracking takes the frequency of the rotation shaft as the instantaneous frequency parameter to analyze the order tracking of the signal. But in this paper, the IF estimation of each component is used as the reference to analyze the order tracking of each component, so as to obtain the amplitude characteristics of the order tracking of the corresponding component signal.
Diagonal slice of bispectrum [34] for a zero-mean Gaussian process, the third-order cumulate and bispectrum are zero. Diagonal slice of bispectrum is a special case of bispectrum and its value is also zero. That is, the residual signal δ (n), its third-order cumulate is: Its bispectrum is: When ω 1 = ω 2 = ω, the diagonal slice of bispectrum is estimated as:Ŝ In mechanical fault diagnosis, the fault signal is often non-Gaussian, and the non-fault signal is often Gaussian. So the effect of Gaussian noise can be suppressed by diagonal slice of bispectrum analysis, thereby extracting the characteristics of the fault signal.

B. ALGORITHM FLOW CHART AND MAIN STEPS
Based on the above principles, the design algorithm flow is shown in Figure 1: Among the figure 1, the algorithm flow chart of the method for adaptive multiridge extraction of peaks detection based on SWT is shown in figure 2.
The main steps of the signal separation method based on the Vold-Kalman filtering of adaptive instantaneous frequency estimation are summarized as follows: (1) The signal x(t) is divided into N sections. (2) The time-frequency matrix TFR(a) is obtained by per- (3) The peak value information vmt a (m, n) and position information amt a (m, n) of wavelet ridge are extracted from the matrix TFR(a).
(4) The peak positions and number of amt a (:, 1) were calculated by the first-order derivative peaks detection method, which was used as the starting positions Xapmx and the number M.
(5) Set the adaptive extraction threshold (6) When the coordinate distance of adjacent starting position D 0 ≤ C 0 , the starting position with the corresponding small peak value in vmt a (:, 1) is set to zero and removed until the threshold condition is not satisfied.
is not an empty set, the ridge line is continuous. (11) If F(k) = ϕ and F(k − 1) = ϕ, then the ridge line is discontinuous. Search forward from starting points.
If F(k +1) = ϕ, continue to judge according to step (12) until F(k + d) = ϕ, the ridge line is discontinuous, where, d is the maximum number of consecutive empty sets allowed.
If f mn (a, b) = 0, then f (b) zeroing and culling. (15) Set the adaptive filtering threshold of adjacent ridges as The signal residual δ(n) is analyzed by diagonal slice of bispectrum to extract the amplitude and frequency characteristics of the residual signal.
(20) Diagnose vibration signals using instantaneous frequency estimation, order amplitude characteristics and amplitude-frequency characteristics of diagonal slice of bispectrum.
This improved method utilizes the method for adaptive multiridge extraction of peaks detection based on SWT to adaptively extract high-precision IF estimates of the signal, that is, to obtain the number of signal decompositions and the instantaneous frequency estimates of each single-component signal. It can reduce the SWT calculation time. It improves the accuracy of peaks detection and solves the problem that multi-ridge extraction is difficult to extract separately in the traditional peaks detection.
This improved method uses the Vold-Kalman filter algorithm to adaptively decompose the time series signal into multiple single-component signals according to IF estimation. At the same time, the sum of the single-component signals is recorded as the reconstructed signal. The difference is subtracted between the original signal and the VOLUME 8, 2020 reconstructed signal. The difference contains noise and other useful information components, which are recorded as the signal residual. The signal is decomposed into multiple single-component signals with actual physical significance. The disadvantages of frequency aliasing, complex components and noise interference caused by direct analysis of the original signal are avoided, and the problem of efficient adaptive separation of complex multi-component signals is solved.
This improved method uses the instantaneous frequency curves estimated obtained by IF estimation to track the order of the corresponding component signals. It solves the problem that for a multi-input gearbox or machine, the input speeds are not unique and cannot be used as the reference frequency. Since each component signal is a single harmonic signal, there is almost no noise interference. Therefore, the effect of order tracking analysis is obvious. There are no approximate orders and cross orders. It is sufficient to extract the amplitude of the corresponding order as the fault feature. The order information missing in the order spectrum can be determined by the relation between the instantaneous frequency curve and the rotation frequency, so as to determine the fault component. It solves the problem of fault diagnosis caused by the non-stationary characteristics of the vibration signal.
Because the signal residual contains not only noise, but also some other useful information, such as the instantaneous frequency of the vibration harmonic components cannot be extracted because of the phase cancellation. Since there is no corresponding IF estimation for signal residual, order tracking analysis cannot be carried out in the case of multiinput speeds, and some weak fault features may be missed. At the same time, based on the characteristic that the noise is Gaussian and the fault signal is non-Gaussian, this paper uses diagonal slice of bispectrum analysis to suppress noise and extract fault features. It solves the problem of fault diagnosis caused by the non-Gaussian characteristic of fault signal.

A. THE SUPERIORITY OF THE IMPROVED METHOD IN MULTI-RIDGE EXTRACTION OF HIGH PRECISION IF CURVES
In order to verify the superiority of the improved method in multi-ridge extraction of high precision instantaneous frequency curves, set the simulation signal 1 to: x(t) = sin(2πt 2 ) + sin(4πt 2 ) + sin(8πt 2 ) + noise (14) where, noise is randomly distributed noise. The sampling frequency is 5120Hz, and the sampling time is 4 seconds. Its instantaneous frequency of speed is ϕ (t) = 2t, and the instantaneous frequency of the three components for x(t) are respectively: ϕ 1 (t) = 2t, ϕ 2 (t) = 4t and ϕ 3 (t) = 8t. Therefore, the simulation signal 1 contains one, double and four times of the speed frequency information.
As shown in Figure 3, (a) (b) (c) is the time domain diagram of simulation signal 1 and its corresponding SWT time-frequency image when the noise is 0dB, 2dB and 4dB. All three time-frequency images show the one, double and four times ridges of the instantaneous frequency of speed, and they are clustered. However, as the intensity of noise increases, the time-frequency ridges become more and more blurred. Especially the time-frequency ridges at four times, the frequency are most seriously affected by noise, and the degree of dispersion for peak points increases. For traditional method of peaks detection, if the selected frequency extraction width is too large, it is easy to cause over-recognition when the initial ridgeline distance is close. But if the selected frequency extraction width is too small, under recognition occurs at the blur of ridgeline. These conditions increase the difficulty of accurate ridge extraction. Figure 3 (d) (e) (f) are the instantaneous frequency curves extracted by the traditional peaks detection method from (a) (b) (c) respectively. Figure 3 (g) (h) (i) are the instantaneous frequency curves extracted by the improved method from (a) (b) (c) respectively. The comparison between (d) (e) (f) and (g) (h) (i) found that the instantaneous frequency curve IF1 extracted by the traditional method has obvious over-recognition in the green circles of (d) (e) (f). It leads to abrupt and stepped frequency changes that become more pronounced as the noise increases. In the red circles of (d) (e) (f), the frequency of instantaneous frequency curves IF2 and IF4, extracted by the traditional method, have obvious underrecognition. It causes the frequency change to be so slow that appearing as 'burrs'. The phenomenon becomes more obvious as the noise increases. The improved method solves such problems of traditional peaks detection under different noise. The extracted curves are smoother.
As shown in Figure 4, where (a), (b) and (c) are the relative error (RE) curves of the instantaneous frequency curves in Figure 3 and theoretical values, respectively. In the figure, the word 'old' represents the traditional peak detection method and the word 'new' represents the improved method. The changes and comparison of the curves verify the analysis about the extraction effect of two methods; (d), (e) and (f) are the comparison about the average relative error (MRE) curves and the data statistical graphs of the instantaneous frequency between the traditional method and the improved method under different noise. These intuitive show the superiority of the improved method. The MRE curves and statistical graphs of IF1 show that the improved method is significantly better than the traditional method in both the value of the error and the amount of error distribution. Due to the under identification of the traditional method, the MRE curves of IF2 and IF4 show that the two methods are basically consistent. But the statistical graphs of IF2 show that the error distribution range of the improved method is smaller than that of the traditional method. The statistical graphs of IF4 data show that the improved method is superior to the traditional method in the distribution of small errors.
Therefore, in terms of suppressing noise, the improved method of adaptive multi-ridge extraction peaks detection is superior to the traditional peaks detection method.

B. THE FEASIBILITY AND SUPERIORITY OF THE IMPROVED METHOD IN SIGNAL SEPARATION AND FEATURE EXTRACTION
In order to verify the feasibility of the improved method in signal separation and feature extraction, the simulation signal 1 when noise=0dB was analyzed.
As shown in figure 5, (a) is the simulation signal 1 and the IF estimation curves when noise=0dB, the signal is a non-stationary signal with acceleration, the instantaneous frequency of each component corresponds to different colors, the same below. , respectively. Diagonal slice of bispectrum is sensitive to non-Gaussian components. Therefore, it can extract 10Hz (one-octave frequency), 30Hz (the sum of one-octave frequency and two-octave frequency) and 50Hz (the sum of one-octave frequency and four-octave frequency) in the signal residual. However, almost no information is extracted from the power spectrum. The frequency in the diagonal slice of bisectrum is not the stationary component, but the average result after bisectrum transform which treats the signal residual as the multi-segment short-time stationary signal. It can be used as auxiliary information to judge the vibration state comprehensively. The above results are in line with the expected results and verify the feasibility of the improved method.
To verify the superiority of the improved method in signal separation and feature extraction, the simulation signal 1 with different SNR is analyzed by traditional signal analysis methods such as EMD and speed order tracking, and the results are shown in figure 6.
In Figure 6, (a) (b) (c) are the simulated signal 1 and the IF estimation curves with different SNR, the harmonic signals are submerged by strong noise. The instantaneous frequency of each component corresponds to different colors, the same below. The speed order tracking and the order diagonal slice of bispectrum are carried out to the signals in (a)(b)(c), and the results are shown in (d)(e)(f). By comparison, it is found that with the increase of noise intensity, the useful information extracted from the rotational speed order spectrum becomes less and less. When noise=6dB, it is impossible to extract useful information from the results of speed order tracking analysis. The effect of bispectrum in suppressing Gaussian noise is obviously. Due to the large number of components of the simulated signal 1, the adjacent order and cross order interference are serious. It leads to unsatisfactory analysis results. At the same time, for the signal of multi-input speeds, there is a problem that the reference frequency cannot be determined. Therefore, for signals with strong noise and complex components, signal separation must be performed.
Perform EMD analysis on the simulated signal 1 in Figure 6 (c) to obtain instantaneous frequency estimates as shown in Figure 6 (g). Compared with the IF estimates in (c), the IF estimates in (g) have frequency aliasing and their physical significance is not clear. As shown in Figure 6 (h), the first eight intrinsic mode functions (IMF1 to IMF8) obtained after EMD analysis are taken, and it can be found that each component signal is a non-stationary signal. Due to the IF estimation of EMD has frequency aliasing, so diagonal slice of bispectrum analysis is performed on the IMF component signals, and the results are shown in Figure 6 (i). Although bispectrum suppresses noise, the frequency obtained by bispectrum is the average frequency of the transformation period. So the analysis results are difficult to meet the application needs.
From the above comparative analysis, it can be found that speed order tracking analysis and traditional signal separation methods have limitations in processing signals with strong noise and complex components. Therefore, it is necessary to conduct in-depth research on the characteristics of the signal and propose more effective signal separation methods.
Using the improved method proposed in this paper, the simulation signal in Figure 6 (c) is analyzed, and the Because the power spectrum cannot suppress Gaussian noise, its analysis results are greatly disturbed, which causes great difficulty in identifying the extracted features. Bispectrum suppresses the noise, and the extraction effect is good. It can be found that the diagonal slice of bispectrum extracts the information of 10Hz (one-octave frequency), 40Hz (fouroctave frequency), 80Hz (eight-octave frequency) and 120Hz (twelve-octave frequency) in the signal residual.
Comparing the analysis results of Fig. 7 (b) and (d) with Fig. 6 (f) and (i), the improved method can separate the components with definite physical meaning, so it is better than EMD in signal separation; the improved method can avoid the influence of neighboring order and Gaussian noise, so it is better than traditional rotational speed ratio analysis in feature extraction.
The relative error between the reconstructed harmonic signals obtained by different methods and the theoretical values is obtained. Since the median is not affected by the maximum or minimum value of the distribution sequence, it is highly representative of the distribution sequence. Therefore, the median of the relative error is obtained as shown in table 1. Through the comparison, found that medians of EMD are greater than 1. This is due to EMD is unable to realize the VOLUME 8, 2020 separation of noise. The combination of traditional peaks detection based on SWT and VKF can realize the separation of noise. But compared with the improved method, the difference of median between the traditional method and the improved method is more obvious with the increase of noise.
To sum up, it is proved that the improved method proposed in this paper has more advantages than the speed order tracking analysis and traditional signal separation methods.

C. COMPARISON OF THE IMPROVED METHOD WITH EMD FOR COMPLEX MULTI-COMPONENT SIGNALS
The simulated signal 1 simulates a rotating system with varying speeds, which verifies the feasibility and superiority of the improved method proposed in this paper. There are complex multi-component signals produced by complex rotating systems in production, such as the vibration signal of multi-stage transmission gearboxes or a complex gearbox with multiple inputs and multiple outputs. In order to verify that the improved method proposed in this paper is also feasible and advantageous for processing the vibration signal of a complex rotating system, set the simulation signal 2 to: where,  x 2 (t) = sin 3.2π (t − 7.5) + 0.11π (t − 7.5) 2 + sin 6πt + 0.08πt 2 + 4 sin (0.2πt) The sampling frequency is 40Hz. The sampling time is 25.5 seconds. The signal from 0 to7.5 seconds is x 1 (t), and the signal from 7.5 seconds to 25.5 seconds is x 2 (t). f (t) contains three components and the instantaneous frequencies are: where, ϕ 1 (t) indicates the state of uniform speed, ϕ 2 (t) indicates the state of speed increase, ϕ 3 (t) represents the state of speed fluctuations. ϕ 1 (t) and ϕ 2 (t) together show the transition process of speed from steady to rising speed. The instantaneous frequency curve has an intersection point at 7.5 seconds, which meets the design requirements of multi-component complex vibration signals.
As shown in Figure 8, (a) is the time-domain diagram and IF estimation of the simulated signal 2, indicating that there is no correlation between the signal components, which can be regarded as a non-stationary multi-component signal in the case of multiple inputs and multiple outputs. Because of the components are not related to each other but affect each other, so the speed order tracking analysis of the signal cannot obtain the expected effect, so effective signal separation must be performed for this type of signal.
Using EMD to decompose simulation signal 2, the IF estimation is shown in Figure 8 (b). Compared with the IF estimation in (a), the IF estimation in (b) is aliasing and has no clear physical meaning. As shown in Figure 8 (c), taking the first four intrinsic modal functions (IMF1 to IMF4) obtained after EMD analysis, it can be found that each component signal is a non-stationary signal. Diagonal slice of bispectrum analysis is performed on the IMFs, the result is shown in Figure 8 (d). Bispectrum analysis has a good effect on noise suppression. But due to frequency aliasing, IMF1 contains the frequency information of IMF2. It indicates that the signal is not sufficiently decomposed in the corresponding frequency band. The frequency information of IMF3 and IMF4 are consistent. It indicates that the signal is excessively decomposed in the corresponding frequency band.
Therefore, in view of the characteristics of complex multicomponent signals generated by complex rotating systems, it is necessary to study a more effective signal separation method. VOLUME 8, 2020 Using the improved method proposed in this paper, the simulation signal 2 in Figure 8 (a) is analyzed, and the result is shown in Figure 9. Among them, (a) is the single non-stationary component signals extracted by VKF, and the color of each component is corresponding to the color of the IF estimation, the same below. (b) is the order tracking analysis of the IF estimation of each component signal. The blue order analysis is a 0.5-order frequency multiplication, and the red order analysis is the coupling between the 1-order frequency multiplication and the 1.5-order frequency multiplication. (c) is the comparison between the reconstructed signal and the original signal and the signal residual. The comparison shows that VKF can be better separated and reconstructed the signal harmonic components. (d) is the diagonal slice of bispectrum analysis and power spectrum analysis of the signal residual in (c), and the extracted information can be used as a reference for fault identification.
Comparing the analysis results of Fig. 9 (a) and (b) with Fig. 8 (c) and (d), it proves that the improved method proposed in this paper is more advantageous than the traditional signal separation method for complex multi-component signals.

IV. EXPERIMENT A. THE DESIGN SCHEME OF EXPERIMENT
In order to verify the practicability of the improved method proposed in this paper, a closed power flow gearbox test bench was used to perform related experiments. The vibration signals of the gearbox were collected under normal and pitting of the tooth surface conditions. The bench is loaded by the internal force generated by the torsion bar. In the experiment, the pitting corrosion of tooth surface was observed and recorded every 60 minutes when the machine stopped. Data were collected once at the start-up and shutdown stage, each time for 60 seconds. This method can be used to simulate the non-stationary operation of the gearbox in the production environment. The gear can reach the limit of fatigue life quickly, and the natural pitting of tooth surface will appear.
The test gear has a transmission ratio of 1: 1 and a number of 18 teeth. The test bench is shown in Figure 10 (a). The model of piezoelectric sensor is CA-YD-186 (sensitivity is 10.41mV/m•s −2 ), the sampling frequency is 12000Hz, and layout positions are shown in Figure 10

B. COMPARISON OF THE IMPROVED METHOD AND EMD ON DIFFERENT SIGNALS
The signals within the 19th to 20th seconds of the speed-up phase are selected, as shown in Figure 11(a) (b) (c), which are the normal signal, the slight gear-pitting signal and the gear-pitting signal. As can be seen from Figure 10 (c), the speed increases from 500r/min to 600r/min, and the calculated instantaneous frequency of gear engagement ranges from 150Hz to 180Hz.
In Figure 11 (a), the one-octave frequency, two-octave frequency and five-octave frequency of the instantaneous frequency of the gear meshing signal are obtained. In Figure 11 (b)(c), the one-octave frequency, two-octave frequency, threeoctave frequency and five-octave frequency of the instantaneous frequency of the gear meshing signal are obtained. Each time-frequency curve in Figure 11 (a) (b) (c) has its specific physical meaning. The results of EMD analysis for the signals are shown in Figure 11 (d) (e) (f). The instantaneous frequency estimation of EMD has no clear physical meaning. At the same time, it cannot suppress the noise, so that the frequency curves of the high-frequency range fluctuate too much.
As shown in Figure 11 (g) (h) (i), taking the first four intrinsic mode functions (IMF1 to IMF4) obtained after EMD analysis, it can be found that each component signal is a non-stationary signal. Diagonal slice of bispectrum analysis was performed, and the results are shown in Figure 11 (j) (k) (l). Due to frequency aliasing, the frequency information of IMF1 to IMF4 partially overlapped, indicating that the signal was not sufficiently decomposed in the corresponding frequency band. Since the IMFs are non-stationary signals, the frequency value is the average value during the conversion period of the bispectral analysis. Therefore, the physical meaning represented by each frequency value cannot be judged, that is, the failure component corresponding to the frequency value cannot be judged.
Using the improved method proposed in this paper, the signals in Figure 11  twice that in (e). It proves that the gear of the vibration signal represented by (f) has a fault and the gear of the vibration signal represented by (e) has a slight fault, which is consistent with the experimental design. (g) (h) (i) are the comparison of the reconstructed signals with the original signals and the signal residuals. The comparison shows that VKF can better separate and reconstruct the harmonic components of the signals. (j) (k) (l) are performed by diagonal slice of bispectrum analysis and power spectrum analysis on the signal residual in (g) (h) (i), respectively. Bispectrum analysis has a good effect on noise suppression. The amplitude range of the bispectrum analysis in (l) is 100 times that in (j) and is 50 times that in (k). Form an obvious frequency cluster in (l). It proves that the vibration signal represented by (l) has a fault.
Comparing the analysis results of Figure 12 and Figure 11, it proves that the improved method proposed in this paper is superior to the traditional signal separation method.

C. COMPARISON OF THE IMPROVED METHOD AND SPEED ORDER TRACKING AND EMD ON THE COMPOSITE FAULT SIGNAL
The above experiment is designed to analyze and compare the fault-free case and the single fault case of gear pitting. The results show that the traditional signal separation method has limitations and deficiencies in the practical application of multi-component non-stationary signals with strong noise.  The improvement method proposed in this article is more practical.
Aiming at the complex multi-component non-stationary vibration signal of multi-input and multi-output, a composite fault experiment was designed for coexistence of gear pitting and bearing outer ring pitting, in which the simulated pitting pits of the bearing outer ring was made by EDM, as shown in Figure 13. Its model number is 32212. This composite fault experiment meets the requirement of multiple outputs for multiple fault vibration sources, which increases the difficulty for effectively extracting different fault characteristics.
Select the signal within the 16th to 17th seconds of the speed-up phase, as shown in Figure 14 (a). From Figure 10 (c), it can be found that the speed increases from 450r/ min to 500r/ min, and the instantaneous frequency of gear engagement is 135Hz to 150Hz. The outer ring fault frequency of the test bearing is 35.5Hz to 39.4Hz.
The IF estimation of the composite fault signal is shown in Figure14(a). Different colors represent different components, the same below. Four IF curves from low frequency to high frequency are: the one-octave frequency of the test bearing outer ring pitting fault frequency, and the one-octave frequency, two-octave frequency and three-octave frequency of the gear meshing frequency.
The speed order tracking is performed to the vibration signal in Figure 14  Comparing with the original signal, it can be found that after isometric resampling, the signal is transformed from a non-stationary signal in the time domain to a stationary signal in the angular domain. It has a certain impact period, but due to the different vibration periods of the two fault sources and the influence of strong noise, it is not possible to directly find specific laws from the angular domain diagram. The order spectrum is shown in Figure 14 (c), which are the traditional order spectrum and the order diagonal slice bispectrum. The comparison shows that the order bispectrum is superior to the traditional order spectrum in suppressing strong noise. However, the two order spectrums are based on the rotation speed. From the order values in (c), it can be found that the order spectrums form order clusters at the 18th, 36th, 54th, and 72th orders, which respectively correspond to one, two, three, and four times of the frequency of gear pitting fault. Because the outer ring pitting fault of the bearing is a weak component, it is easily overwhelmed by strong noise and other faults. Although the speed order tracking can be extracted the pitting fault characteristics of gear, but there are obvious deficiencies in extracting weak fault components. VOLUME 8, 2020 It is necessary to extract different components. It is also the demand characteristics of the composite fault.
Using EMD to separate the vibration signal in Figure 14 (a), the IF estimation of EMD is shown in Figure 14 (d). As shown in Figure 14 (e), the first four intrinsic mode functions (IMF1 to IMF4) obtained after EMD analysis are taken, and each component signal is a nonstationary signal. The diagonal slice of bispectrum analysis is performed on the IMFs, and the results are shown in iFigure 14 (f). Due to frequency aliasing, the frequency information of IMF1, IMF2, and IMF3 partially overlap. It indicates that the signal is insufficient decomposition in the corresponding frequency band. Since the IMFs are non-stationary signals, so the physical meaning of each frequency value cannot be judged.
Using the improved method proposed in this paper, the composite fault signal in Figure 14 (a) is analyzed, and the result is shown in Figure 15. Among them, (a) is the single non-stationary component signals extracted by VKF, and the color of each component is corresponding to the corresponding IF estimation, the same below. (b) is the order tracking analysis of each component signal by the IF estimates, and the results show that there is only first-order information, indicating that VKF has a good effect on the separation of noise and the signal of each component. The order spectrum of the blue component corresponds to the fault frequency of the bearing outer ring, and the rest of the order spectrum corresponds to the one, two and three times of the frequency of gear pitting fault. The analysis results in (b) are good separation of different components of the multi-component signal. (c) is the comparison between the reconstructed signal and the original signal and the signal residual. The comparison shows that VKF can better separate and reconstruct the harmonic components of the signal. As shown in (d), the signal residuals in (c) are performed by diagonal slice bispectrum analysis and power spectrum analysis respectively.
Comparing the analysis results of Figure 14 and Figure 15, in terms of signal separation and feature extraction of complex multi-component signals with multiple inputs and multiple outputs, it proves that the improved method proposed in this paper is more advantageous than the traditional order tracking analysis and the traditional signal separation method.

V. CONCLUSION AND RECOMMENDATIONS
This paper proposes a signal separation method based on Vold-Kalman filter of adaptive instantaneous frequency estimation, and analyzes the simulation data and experimental data under different conditions. Comparison of different signal separation methods, the results show that the improved method proposed in this paper has advantages and practicability in adaptively separating and extracting the fault characteristics of complex multi-component non-stationary signals and suppressing noise.
In future work, it is recommended to extend this improved method to other composite faults and vibration environments, such as composite faults of planetary gearboxes. At the same time, the research method in this paper laid the foundation that complex multi-component signals need to be separated and extracted non-stationary fault characteristics, but the actual application and effect still need to be verified in industrial machines.