Abstract
The vibration signal of rotating machinery compound faults acquired in actual fields has the characteristics of complex noise sources, the strong background noise, and the nonlinearity, causing the traditional blind source separation algorithm not be suitable for the blind separation of rotating machinery coupling fault. According to these problems, an extraction method of multisource fault signals based on wavelet packet analysis (WPA) and fast independent component analysis (FastICA) was proposed. Firstly, according to the characteristic of the vibration signal of rotating machinery, an effective denoising method of wavelet packet based on average threshold is presented and described to reduce the vibration signal noise. In the method, the thresholds of every node of the best wavelet packet basis are acquired and averaged, and then the average value is used as a global threshold to quantize the decomposition coefficient of every node. Secondly, the mixed signals were separated by using the improved FastICA algorithm. Finally, the results of simulations and real rotating machinery vibration signals analysis show that the method can extract the rotating machinery fault characteristics, verifying the effectiveness of the proposed algorithm.
1. Introduction
Rotating machinery plays an increasingly important role in modern industry and intelligent manufacturing. The real-time monitoring of the working state of rotating machinery can not only avoid the occurrence of disasters but also increase obvious economic benefits [1, 2]. The rotor system is one of the important joints of rotating machinery equipment in industrial production. It is widely used in different fields. The normal operation of the rotor system is directly related to the working performance of the whole mechanical equipment. Therefore, the rotor system is a good research object of the exploration and mining of new fault diagnosis methods. In practical engineering, the rotor vibration signal caused by local defects usually has the characteristics of nonlinear, nonstationary, low signal-to-noise ratio, and unclear fault characteristics [3, 4]. It is difficult to make effective diagnosis by directly using spectrum analysis. In addition, some conventional diagnostic methods such as spectral kurtosis and time-frequency analysis have their own limitations. Therefore, it is a practical problem to explore an effective fault diagnosis method for the rotor system in engineering practice [5, 6].
In order to effectively identify the feature information contained in the fault signal of rotating machinery and reveal its inherent characteristics, many fault feature extraction methods of rotating machinery have been proposed, such as empirical mode decomposition (EMD) [7, 8], mathematical morphology filtering [9, 10], wavelet decomposition [11, 12], adaptive filtering [13, 14], matching pursuit [15, 16], cyclostationary signal analysis [17, 18], Wiener filter [19], Kalman filter [20, 21], and stochastic resonance [22, 23] that are widely used in early fault diagnosis of rotating machinery. The EMD proposed by Huang et al. [7] is a nonstationary signal analysis method, which can find the hidden characteristic information in the signal, and has been widely used in the extraction and noise reduction of the impact signal of rotating machinery. Xin et al. [24] proposed a fault feature extraction and diagnosis method for vibration signal based on iterative empirical wavelet transform (EWT) and sparse filter (SF). Pang et al. [25] proposed an adaptive filtering algorithm based on mathematical morphology for rolling bearing fault diagnosis. Kang et al. [20] proposed an improvement method associated with the Kalman filter to estimate the bearing dynamic coefficients. Cheng et al. [26] proposed a novel deconvolution algorithm called adaptive multipoint optimal minimum entropy deconvolution adjusted (AMOMEDA) for extracting fault-related features from noisy vibration signals. The condition monitoring and diagnosis technology of mechanical equipment has made many encouraging achievements and has been gradually applied to practice. However, there are still many problems in the research of fault diagnosis mechanism, signal detection, feature extraction, and fault reasoning, which need to be further studied. For example, the detection method of fault characteristic signal of rotating machinery is not universal. One of the important problems is that the working condition of the equipment is often complex, and the fault characteristic signal is often superimposed by various interference signals. How to effectively remove the interference signal and extract the effective signal has become a key problem. The emergence of blind signal separation (BSS) theory can provide a feasible way to solve the above problems. Blind signal separation refers to the process of recovering the source signal only from the observed signal according to the statistical characteristics of the source signal when the parameters of the source signal and transmission channel are unknown.
BSS plays an increasingly important role in the field of digital signal processing and has been widely used in communication [27], speech processing [28], fault diagnosis [29, 30], seismic exploration [31], biomedicine [32, 33], image processing [34], radar [35], and economic data analysis [36]. In blind signal separation, the typical algorithms commonly used include the fast fixed-point algorithm [37], natural gradient algorithm [38], Equivariant Adaptive Separation via Independence (EASI) algorithm [39, 40], and Joint Approximation Diagonalization of Eigen-matrices (JADE) algorithm [41, 42], etc. Grotas et al. [43] developed the constrained maximum likelihood (ML) estimator of the Laplacian matrix for this graph BSS problem with Gaussian-distributed states. Eitner et al. [44] used two blind source separation algorithms to estimate the modal parameters of a reduced-scale rocket nozzle using only measurements of deformation. Ge and Jiang [9] proposed a framework based on joint blind source separation (JBSS) in order to solve the jamming suppression problem in the noise environment for the distributed radar with single transmitter and multiple receivers, where the multiple jamming enter into all the receivers through the main beam of the antennas. Belaid et al. [45] proposed a new multiscale decomposition algorithm which enables the blind separation of convolutely mixed images. All the above algorithms show good separation performance when separating noiseless mixed signals. However, in the separation of noisy signals, there will be a large error, especially in the case of low signal-to-noise ratio, and even draw a completely wrong conclusion, because these algorithms are derived without considering the noise model. In the process of rotating machinery running, the signal measured by vibration sensor installed on mechanical equipment inevitably contains signal noise. If the noise is ignored and the blind source separation algorithm is directly used to separate the mixed vibration signals, it may produce large errors or draw wrong conclusions. Therefore, before the blind separation of mechanical vibration signals, it is very important to reduce the noise and improve the signal-to-noise ratio.
Therefore, aiming at the problem of weak fault feature signal extraction of rotating machinery under the influence of noise, a fault separation method based on wavelet packet analysis (WPA) filter and improved fast independent component analysis (FastICA) algorithm is proposed. Firstly, wavelet packet filter is used to denoise the mixed signal under noise interference, and then the improved FastICA algorithm is used to separate the denoised signal, and then the weak fault signal is effectively extracted.
The contents of the following sections are as follows: Section 2 introduces the basic theory and model of blind source separation. Section 3 presents the principle of the wavelet packet denoising algorithm. Section 4 presents implementation algorithm steps of the improved WPA-FastICA Method. Simulated and experimental verifications are conducted in Sections 5 and 6. Finally, the conclusions and outlook are both given in Section 7.
2. The Basic Theory and Model of Blind Source Separation
Blind source separation (BSS) refers to the process of extracting and recovering (separating) the original signals which cannot be directly observed from multiple observed mixed signals. The “blind” means that the source signal is unobservable, and the hybrid system characteristics are unknown or only a small amount of prior knowledge is known. Independent component analysis (ICA) is the main method to solve the problem of blind source separation, which is based on the independence of the source signal. The solution steps of ICA are as follows: establish the objective function of the observed signal according to the principle of statistical independence and then use the optimization algorithm to separate the observed signal, so that the output component is as independent or completely independent as possible in statistics. In this way, the signal can be separated.
The mathematical model of fast ICA problem is described as follows: where is the observation signal with -dimension, is the mixing matrix of source signal aliasing and unknown, is the source signal, and is the noise.
It can be expressed as a vector
The purpose of fast ICA is to find a separation matrix to separate the source signal from , and the mathematical model of the separation process is expressed as follows: which is Where is the estimated value of .To separate the source signal is to find the separation matrix to make the following true: If where if is a unit matrix, then is an effective separation of .
3. The Principle of the Wavelet Packet Denoising Algorithm
Wavelet packet analysis is a more detailed analysis and reconstruction method of signal from wavelet analysis. In wavelet analysis, the signal is actually decomposed into low-frequency rough parts and high-frequency details, and then only the low-frequency details are decomposed for the second time, instead of the high-frequency parts [46, 47]. But in practice, we are only interested in some undetermined time (point) or frequency domain (point) signal, only need to extract the information of these specific time and frequency points. In this case, the fixed distribution of the time-frequency window of orthogonal wavelet transform is not an optimal choice. Wavelet packet analysis not only decomposes the low-frequency part but also makes a secondary decomposition of the high-frequency part, which can make a more detailed description of the high-frequency part of the signal, and has a stronger ability to analyze the signal.
In order to simplify the problem without losing generality, considering that the collected vibration signal contains noise signal and has statistical characteristics, and considering that the measured vibration signal data is discrete data in time, the equation for vibration signal is established as follows: where is the sampling signal, is the mixing matrix of source signal aliasing and unknown, is the source signal, and is the noise signal.
Suppose is the orthogonal low-pass real filter corresponding to the orthogonal scaling function , and is the high pass filter corresponding to the orthogonal wavelet function , where ; that is, the two coefficients satisfy the orthogonal relationship, and then the recursively defined function is called the wavelet packet determined by the orthogonal scaling function . So, there is the following function: Therefore, for any nonnegative integer Nez, there exists the following equation: If the orthogonal wavelet decomposition algorithm of the Mallat multiresolution analysis algorithm is extended to the wavelet packet decomposition algorithm, the fast decomposition and reconstruction formula of wavelet packet can be obtained. It is expressed as follows: where the coefficients and are the projections of an approximation function at scale on subspaces and , respectively. Therefore, wavelet packet decomposition can not only decompose the low-frequency components of the signal but also can decompose the high-frequency components. Figure 1 is the schematic diagram of wavelet packet decomposition, in which represents low-frequency coefficients, represents high-frequency coefficients, and numbers indicate the number of levels of wavelet packet decomposition.
(a) The time domain signals
(b) The frequency domain signals
A group of normal orthogonal bases extracted from wavelet library which can form are called wavelet packet bases. Different wavelet packet bases have different abilities of time-frequency localization, and wavelet packet decomposition of a signal can be realized by many wavelet packet bases. Therefore, in order to better express the characteristics of signals, it is necessary to find an optimal wavelet packet basis. People usually define an information cost function that satisfies both additive and concentration properties to select the optimal wavelet basis. The common information cost functions include Shannon entropy, norm , logarithmic energy entropy, and threshold entropy. The commonly used definition of Shannon entropy is as follows: where .According to these information cost functions, the wavelet packet sequence which minimizes the information cost function is obtained. The best wavelet packet basis can be obtained by using the fast search method from bottom to top.
In wavelet analysis, threshold denoising is a common method. Similarly, in wavelet packet analysis, the threshold method can also be used to process the coefficients of wavelet packet decomposition. In practice, how to select the threshold is the key to the algorithm design. The usual threshold selection is based on Donoho’s hard threshold and soft threshold. The mathematical expression is as follows: Where is the wavelet decomposition coefficient before processing, is its symbol, is the coefficient after processing, , is the maximum number of layers of wavelet decomposition, is the threshold value, when , and it is the hard threshold method; when , it is the soft threshold method.
Generally, the threshold selection criteria include unbiased likelihood estimation, fixed threshold criterion, hybrid criterion, and minimax criterion. Here, the threshold AA is selected as where is the number of high-frequency coefficients in the corresponding decomposition level, and is the noise standard deviation. Generally, its value is the median value of the absolute value of wavelet packet decomposition coefficient in each layer divided by 0.6745. The formula can be expressed as follows: where is the median value. In addition, soft threshold denoising can be divided into global threshold denoising and layered independent threshold denoising. The global threshold means that the same threshold is used to process each layer coefficient of wavelet packet decomposition, while the layered independent threshold means that different thresholds are used when processing each layer coefficient. The results show that the global denoising threshold value is better than the independent threshold value.
In this paper, an average threshold method is proposed to extract the vibration signal of rotating machinery. Firstly, the noise standard deviation of each node after wavelet packet decomposition is obtained by using formula (14), and the threshold value of each node is calculated by formula (13). Then, the threshold value of all nodes is averaged as the global threshold of wavelet packet denoising. The formula can be expressed as follows: where , , is the number of nodes after wavelet packet decomposition, .
Therefore, according to the wavelet packet analysis method, the following steps can be adopted to filter the noisy signal collected by the sensor. The specific algorithm steps are as follows.
Step 1. The level of wavelet packet decomposition is determined, and the vibration signal of rotating machinery under strong background noise is decomposed by wavelet packet.
Step 2. Shannon entropy is chosen as the information cost function to determine the optimal orthogonal wavelet packet basis.
Step 3. For the coefficients on each node of the optimal wavelet packet basis, the average threshold is obtained according to formula (15), and the threshold is processed by using the soft threshold method.
Step 4. The coefficients of the optimal wavelet packet basis after thresholding are reconstructed to realize the filtering of weak vibration signal under strong noise background.
4. The Algorithm of Blind Source Separation
4.1. Problem Description of Blind Source Separation
The ICA method is a kind of signal processing method to extract independent signal sources when the observed data are known and the mixing mode of signal sources is unknown. This method takes the non-Gaussian source signals as the research object, and under the assumption of their statistical independence, blind separation is carried out for the multichannel observed mixed signals, so as to separate the independent source signals hidden in the mixed signals perfectly The classical ICA model is as follows: where is the -dimensional observation signal; is the -dimensional source signal; is the mixed matrix, which is generally assumed to be .
By inverse equation (16), we can get the following results: where is the estimation of source signal and is the unmixing matrix.
The basic idea of ICA is to adjust the weight vector to make the components of the output signal as independent as possible. In this paper, the FastICA algorithm with simple calculation and fast convergence speed is used to extract the source signal.
The iterative formula of the separation matrix is as follows: In the formula, is the whitened signal of the observation signal ; the nonlinear function is , and .
In practice, the observed signal is often mixed with various noises. The model is as follows: where is the additive white Gaussian noise; so, the estimated is the noisy signal, where is the unbiased estimation of ; is the noise signal. When it is sparse, the standard ICA model can carry out sparse coding of uncorrelated components and then select the contraction function to remove the noise in it. The probability density function model of the signal needs to be selected. In this paper, the double exponential function is selected to approximate the fault signal. The function is as follows: where is the sparse distribution and is the scale value.
can be obtained by maximum likelihood estimation as follows: where is the contraction function, which is defined as where is the variance of the signal containing noise.
4.2. Algorithm Steps of the Improved WPA-FastICA Method
In this paper, the specific steps of the multisource weak signal extraction method are as follows:
Step 1. Using sensor to collect -channel linear aliasing noise observation signal . The wavelet packet decomposition of each observation signal is to use as the wavelet basis function to decompose the noisy signal in five layers and use Shannon entropy to get the reconstructed denoised signal.
Step 2. The denoised signal is deaveraged to make and then whitened to obtain the signal . The -path signal in is obtained after deaveraging and whitening.
Step 3. Let be equal to the number of independent source signals to be estimated.
Step 4. The initial weight vector is generated randomly, so that .
Step 5. .
Step 6. Orthogonal normalization: .
Step 7. If does not converge, return to step (5).
Step 8. If converges, then . If , return to step (4).
Step 9. The solution matrix is , and according to , the unmixing signal is .
Step 10. is calculated from equation (21), and pseudoinverse transformation is performed: .
Finally, the separation signal is .
5. Simulation
5.1. Evaluation Index of Separation Performance of the Method
In the fast ICA algorithm, the main indexes to evaluate the separation performance of the algorithm usually include correlation coefficient, normalized mean square error (NMSE), and scatter plot of source signal and separated signal. At the same time, when the separation results with the same accuracy are obtained, the number of iterations is compared to evaluate the performance of the three algorithms.
5.1.1. Correlation Coefficient
It describes the parameters of correlation between two signals. Suppose that the -th signal in the source signal is , and the estimated signal obtained after the algorithm is run is , then the calculation formula of the correlation coefficient between and is as follows [48]:
where is the variance.
The closer the absolute value of the correlation coefficient between the two signals is 1, the closer the signal is, the better the separation effect of the algorithm is.
5.1.2. Normalized Mean Square Error (NMSE)
Normalized mean square error measure is used to measure the accuracy of hybrid matrix estimation. Its expression is as follows [49]: In equation (18), is the number of rows of , is the number of columns of , and is the estimated mixed matrix element. The smaller the value of NMSE, the higher the accuracy of mixed matrix estimation.
5.1.3. Number of Iterations
It refers to the number of times that the formula in the algorithm is repeatedly calculated in order to approximate the desired target or result. When the ideal results can be obtained, the less iterations, the shorter the running time and the faster the convergence speed.
5.2. Simulation of Mixed Signal Separation under Different Carrier Frequencies and Noises
In order to effectively simulate the fault characteristics of rotating machinery, the fault signal features under the actual working conditions of rotating machinery are selected to simulate the unbalance and rub impact faults. The simulation functions of these four types of source signals are as follows: where simulates the unbalance fault signal and simulates the characteristic signal of rub impact fault.
A mixing matrix is generated randomly, and the two source signals are mixed according to equation (6). where is the observable mixed signal and is the mixing matrix, where The time-frequency waveform of the mixed signal generated by matrix a mixing is shown in Figure 2; the time-frequency waveform of two source signals is shown in Figure 1.
(a) The time domain signals
(b) The frequency domain signal
It can be seen from Figure 1 that the simulated unbalance and rub impact fault signals are concentrated at 50 kHz and 100 kHz, respectively. Therefore, when wavelet packet is used for noise reduction, wavelet is selected for 5-layer wavelet packet decomposition, and the second to sixth characteristic wavelet packet is selected for reconstruction, and the signal with frequency of 42.5-127.5 Hz is obtained as wavelet packet denoising signal. The interference noise is generated automatically by MATLAB.
5.3. Result Analysis
In order to compare with other separation methods, we add classic FastICA and SOBI separation methods to separate mixed signals and compare with WPA-FastICA method proposed in this paper.
The multisource unbalance and rub impact fault signals mixed with different intensities of white Gaussian noise are separated by using the extraction method in this paper. The source signals are mixed by random matrix, and then 0 dB, 10 dB, and 20 dB random white Gaussian noises are added, respectively, to obtain the observation signal shown in Figure 3. In order to intuitively compare the characteristics of the vibration signals before and after separation, it is necessary to analyze and compare the spectrum characteristics of the separated vibration signals under different noises.
(a) The time domain signals
(b)
(c)
(d)
Figure 4 is the time-frequency waveform after separating the simulation aliasing fault signal by using the classic FastICA algorithm. It can be seen from Figure 4 that the classical FastICA algorithm has a very good separation effect when the signal-to-noise ratio is greater than 20, but when the signal-to-noise ratio is less than 20, the separation result is wrong. Figure 5 is the time-frequency waveform after separating the simulation aliasing fault signal by using the classical SOBI algorithm. It can be seen from Figure 5 that the separation result of the classical SOBI algorithm is wrong no matter under any Gaussian noise interference.
(a) The separated time domain signals with different noises by FastICA
(b) The separated frequency domain signals with different noises by FastICA
(a) The separated time domain signals with different noises by SOBI
(b) The separated frequency domain signals with different noises by SOBI
It can be seen from the time-domain waveform in Figure 6 that the algorithm proposed in this paper is used to separate the sampled signal under different noise interference, and the noise can be well filtered from the time-domain waveform. Through the analysis of the frequency domain waveform characteristics, rub impact and unbalance fault can also be clearly identified when the signal-to-noise ratio is 0. The simulation fault is consistent with the separated fault signal, which shows that the WPA-FastICA algorithm can be used to separate the coupling faults of the rotor system.
(a) The separated time domain signals with different noises by WPA-FastICA
(b) The separated frequency domain signals with different noises by WPA-FastICA
Therefore, through Table 1 and Figures 4–6, compared with the three methods, the extraction effect of the proposed method is obviously better than that of the other two methods, and the direct separation by the SOBI method has the worst effect.
6. Applications
In order to effectively verify the effectiveness of the proposed method, a fault simulation experimental platform is built. The experimental table is shown in Figure 7(a). The experimental platform is 134 mm long, 51 mm wide, 120 mm high, and weighs about 50 kg. The rotor is composed of two shafts about 49 mm long and connected by coupling. There are three large mass disks and four small mass disks installed on the rotor, which are used to simulate the blades of rotating machinery in the industrial field. 13 sensors are installed on the bearing. The front 12 sensors are installed in the direction and direction of the same section in a mutually perpendicular way, as shown in Figure 7(b). They are used to measure the vibration displacement signal of the rotor. After the vibration signals output by the sensors in the and directions are combined, the rotor axis trajectory diagram of this section can be obtained. In the experimental process, the inverter is used to control the operation of the variable frequency motor, and then the rotor is driven by the variable frequency motor, and the motor power is about 1.1 KW. The minimum speed of the test bed can be set arbitrarily, and the maximum speed can reach 12000 rpm, which is far more than the critical speed of the system, which can meet the experimental requirements. Moreover, the speed up and down of the rotor can be adjusted through the keyboard of the frequency converter. The control mode of the motor can be manual control or computer control.
(a) The rotor test bed
(b) Installation position of sensors
The rotor test bed is an experimental device used to simulate the vibration of rotating machinery. Various rotor fault simulation experiments can be done on the rotor test bench, such as rotor misalignment, imbalance, pedestal looseness, dynamic and static rubbing, oil film whirl, and oil film oscillation. In addition, it can also be used for lubrication theory, nonelectric measurement, data acquisition, signal analysis, rotor bearing system dynamics, condition monitoring, and fault diagnosis. In this experiment, by adding bolts on the mass disk of the rotor test bed and adding rubbing device on the bearing, the simulation experiment of rotor unbalance and dynamic static rubbing fault is realized. The experimental speed is 3000rd/min, and the sampling frequency is 5000 Hz.
The time-frequency waveform of rotor vibration signal collected by two sensors is shown in Figure 8. The classical second-order blind identification (SOBI) and separation algorithm are used to separate the sampled signal directly, and the time-frequency waveform of the separated signal is shown in Figure 9. The classical fast independent component analysis (FastICA) algorithm is used to separate the sampled signal directly, and the time-frequency waveform of the separated signal is shown in Figure 10. The WPA-FastICA algorithm proposed in this paper is used to separate the sampled signal directly, and the time-frequency waveform of the separated signal is shown in Figure 11.
It can be seen from the time domain waveform in Figure 8 that the sampled signal is interfered by strong noise, and the fault characteristics of the rotor system cannot be identified from the time domain. However, from the characteristics of the frequency domain signal, we can clearly see the first harmonic of 50 Hz and the weak second harmonic of 100 Hz. From the analysis of frequency domain characteristics, only slight unbalance fault can be identified.
As can be seen from the time-domain waveform in Figure 9, the time-domain waveform separated by the classical SOBI algorithm contains strong noise; so, the fault characteristics of the rotor system cannot be identified from the time-domain waveform. However, from the characteristics of frequency-domain signals, it is obvious that one of the signals contains the first harmonic of 50 Hz and the second harmonic of 100 Hz, and the other signal is random and disordered. From the analysis of frequency domain characteristics, only slight unbalance fault can be identified.
It can be seen from the time-domain waveform in Figure 10 that under strong noise interference, the classical FastICA is used to directly separate the sampling signal, and the fault type of the rotor system cannot be identified from the time-domain waveform. However, through the frequency domain waveform characteristics, rub impact and unbalance fault can be identified, but the noise interference is particularly obvious.
As can be seen from the time-domain waveform in Figure 11 that under the strong noise interference, the algorithm proposed in this paper is used to separate the sampled signal. From the time-domain waveform, it can be seen that the noise is well filtered. According to the frequency domain waveform characteristics, rub impact and unbalance fault can be clearly identified. The simulated fault of the test bed is consistent with the separated fault signals, which shows that the WPA FastICA algorithm can be used to separate the coupling faults of the rotor system.
7. Conclusions
A blind source separation algorithm based on WPA-FastICA is proposed to solve the problem of fast feature extraction for weak fault of rotating machinery. Our conclusions are as follows: (1)In the background of noise interference, the vibration signal of rotating machinery has nonlinear and nonstationary characteristics, which will affect the effective extraction of fault features. In order to solve this problem, an algorithm combining wavelet packet filtering and blind source separation are introduced to effectively extract the components related to fault features, eliminate irrelevant components, and enhance the effective extraction of weak fault features(2)If the traditional blind source separation (BSS) algorithm is used to separate the observed signals directly under noise interference, the separation result error is very different, especially when the signal-to-noise ratio is lower than 15 dB, even the wrong result may be obtained(3)The improved WPA-FastICA blind source separation algorithm proposed in this paper can effectively separate the simulated rotor vibration signals, and its performance has been greatly improved compared with the traditional method under different noise interferences(4)The improved WPA-FastICA blind source separation algorithm can separate the real rotor fault signal through the time and frequency domain analysis of the separated signal, which is consistent with the experimental set-up fault
Data Availability
The data used to support the findings of this study are included within the article.
Conflicts of Interest
The authors declare no conflict of interest.
Acknowledgments
The authors thank all those who helped in the course of this research. This project was also supported by the National Natural Science Foundation of China (Grant No. 51675253) and Key Technologies R&D Program of Henan Province (Grant No. 172102210097).