Abstract

Wavelet analysis is a powerful tool for signal processing and mechanical equipment fault diagnosis due to the advantages of multiresolution analysis and excellent local characteristics in time-frequency domain. Wavelet total variation (WATV) was recently developed based on the traditional wavelet analysis method, which combines the advantages of wavelet-domain sparsity and total variation (TV) regularization. In order to guarantee the sparsity and the convexity of the total objective function, nonconvex penalty function is chosen as a new wavelet penalty function in WATV. The actual noise reduction effect of WATV method largely depends on the estimation of the noise signal variance. In this paper, an improved wavelet total variation (IWATV) denoising method was introduced. The local variance analysis on wavelet coefficients obtained from the wavelet decomposition of noisy signals is employed to estimate the noise variance so as to provide a scientific evaluation index. Through the analysis of the numerical simulation signal and real-word failure data, the results demonstrated that the IWATV method has obvious advantages over the traditional wavelet threshold denoising and total variation denoising method in the mechanical fault diagnose.

1. Introduction

The research of faulty characteristics extraction method is the basis of mechanical equipment fault diagnosis, which is directly related to the accuracy of fault diagnosis results [1, 2]. Traditional signal processing methods, such as Fast Fourier Transform (FFT) and time-frequency analysis, are based on the assumption that the vibration signal is stable and linear. However, when the mechanical equipment is in the state of failure, the dynamic behavior of the mechanical equipment often shows complex nonlinear and nonstationary characteristics [3, 4]. Therefore, the nonstationary signal processing methods such as Empirical Mode Decomposition (EMD) and wavelet transform method become the research hotspot. EMD is introduced by the characteristic of adaptive signal decomposition. However, it lacks theoretical basis and there are some problems in its own algorithm such as the phenomenon of model mixing and the end effect [5]. Currently, the wavelet analysis method is the most important tool to deal with the nonstationary signals [6, 7]. Specially, the parameters of the wavelet transform analysis are adjustable. Hence, the profile and the details of the signal can be detected. For this reason, wavelet analysis is also called “mathematical microscope.” Wavelet threshold filtering is employed based on the assumption that the larger amplitude of coefficient is produced by the useful signal. Because the wavelet transform has good time-domain and frequency-domain local characteristics, it has been widely used in various fields, such as nanoring driven in logic gates and memory mass devices [8], wind farm power prediction [9], and mechanical fault diagnosis [10].

Although the wavelet analysis has achieved good results in mechanical fault diagnosis, it also has some limitations. When the noise signal in the wavelet domain is able to achieve sparse representation, the wavelet denoising method based on the threshold processing has a very strong effectiveness, but it also will produce spurious noise spikes and the pseudo-Gibbs oscillations [11]. Noise spikes are due to noisy wavelet coefficients exceeding the threshold, whereas pseudo-Gibbs artifacts are due to nonzero coefficients being erroneously set to zero [12]. Therefore, the choice of thresholding is very important in the algorithm of wavelet denoising. The total variation minimization has been successfully applied to recover thresholded coefficients [13, 14]. In particular, the significant coefficients are maintained while the insignificant coefficients are then optimized to minimize a TV objective function. However, the single total variation denoising method often produces undesirable staircase artifacts [15]. Therefore, researchers have combined the traditional wavelet denoising method with the total variation method recently, which is called wavelet total variation (WATV) [16]. WATV method takes advantage of the multiresolution analysis of wavelet analysis and preserves the advantages of the total variation method in dealing with the spurious noise spikes and the pseudo-Gibbs oscillations. It is worth mentioning that a novel nonconvex regularizer was introduced in the WATV algorithm as it can better recover signal [17]. The literature indicated that the WATV method has some advantages over the traditional soft-thresholding denoising method, hard-thresholding denoising method, and total variation denoising method. However, the actual effect of WATV method is affected by the accuracy of noise variance estimation, which limits its performance.

This paper proposes an improved WATV method, namely, IWATV, based on the original theory of WATV. The main idea is to estimate the noise variance of the wavelet coefficients by using the local variance analysis [18], which is obtained by wavelet decomposition of noisy signal. This method avoids the interference of the human experience to the noise variance estimation. In order to verify the validity of the proposed method, numerical simulation signal, the faulty bearing data in public data set, and gearbox vibration data in the field of the metallurgical industry are analyzed. The results show that the proposed method can effectively identify the fault characteristics and has obvious advantages compared to other traditional methods.

2. Theoretical Descriptions

2.1. Basic Principles of WATV Denoising Algorithm

Assume that a superposition of the additive white Gaussian noise of finite length signal can be expressed aswhere represents the actual observation signal, represents the noise signal with the variance , and is a useful signal (wanted signal) to be separated from the observation signal. Wavelet transforming on the signal can obtain the corresponding wavelet coefficients , where represents the operator of wavelet transform, the subscript of , respectively, represents the scale and time. In this paper, we use undecimated wavelet transform, which satisfies Parseval frame condition; that is, .

The total variation of signal () can be expressed as . The first-order difference matrix () is defined as

Generally, norm () is expressed as . When , the formula is set as norm; namely, . When , the formula is set as norm; namely, . Then the total variation (TV) of signal is defined as .

WATV method obtains wavelet coefficients from the following convex optimization problem by the total variation regularization:

Through the inverse transform of wavelet coefficients obtained by optimization, the estimate of the original signal is obtained. In formula (3), the penalty term is the total variation of the estimated signal . Among them, , are the regularization parameters. is a parameter penalty function. The function is chosen considering the characteristic of both nonconvex and sparsity. Typical parametric penalty function is arctangent function and its expression is listed as follows:

In order to guarantee that the penalty function satisfies the characteristics of sparse and nonconvex regularization, the range of the penalty parameter should be . In formula (3), the choice of , , and has a great effect on the actual noise reduction. In general, the value is chosen as or . The main purpose is to ensure the maximum degree of formula (3) sparse, while ensuring that the optimization function is convex. When , formula (3) has a problem of a single wavelet analysis of the noise reduction, and when , formula (3) is converted to a total variation of noise reduction. Generally, the values of and are set as follows:

In (5), the weight of the wavelet and the total variation regularization are controlled by and it is usually limited to . In this paper, wavelet threshold is set to be , where represents the noise variance of the wavelet scale. In undecimated wavelet transform, we can get , where represents the variance of the noise. According to the total variation denoising theory, , so (5) can be converted to

2.2. Improved Wavelet Total Variation (IWATV) Denoising Algorithm

Considering that the multiresolution of wavelet analysis theory is able to describe nonstationary signal and the characteristic of total variation regularization, the improved total variation wavelet denoising approach is introduced. According to the distribution of signal and noise under different scales by convex optimization, the contradiction between noise reduction and feature loss can be effectively resolved, which provides a powerful tool for fault characteristics extraction and signal reconstruction. From (3), it can be seen that the selection of parameters () is very important for the original WATV noise reduction method. It is directly related to the effect of the actual noise reduction. Basically, the choice of parameters depends on the noise variance estimation by the actual measurement signal. In this paper, the local variance analysis is used to solve the problem of variance estimation about the noise signal so as to avoid the interference by human experience.

In this analysis, the signals are assumed being corrupted by additive gauss white noise (AWGN) with known variance. Let represents the orthonormal wavelet coefficients of the “clean” signal. The wavelet coefficients of the noisy signal are given bywhere is AWGN due to orthonormality of the chosen wavelet transform. Let represents the variance of the wavelet coefficients of wanted signal. Moreover, the variance of wavelet coefficients can be regarded as the average energy of the signal in the corresponding scale. The value of the variance varies with the size of the scale, and the variance of different spatial positions (such as edges and smooth regions) is also different in the same detail. Thus, it is unreasonable to estimate the noise variance in the whole details. To address this problem, the method of local variance analysis is proposed by Mihçak et al. [18]. The estimation of the variance field is the crux of the proposed denoising method. For each data point , an estimation of the noise variance is formed based on a local neighborhood . We use a square window centered at . Assuming that the correlation between variances of neighboring coefficients is high, the variance of the signal is estimated by Maximum Likelihood (ML) instead of estimating in the entire subbands. The calculation equation is as follows:where is the number of coefficients in , is the standard variance of wavelet noise in the scale of , and is th wavelet coefficients of noisy signal.

The flowchart of the proposed method is shown in Figure 1. The main steps of the improved WATV method in this paper are as follows.

Step 1. Select the appropriate wavelet basis function and the wavelet coefficients by wavelet transforming to the noisy signal.

Step 2. The local variance analysis was applied to the wavelet coefficients aiming to estimate the noise variance () of the noisy signal by (8).

Step 3. Calculate the key parameters () of the WATV method with (6).

Step 4. Calculate the new wavelet coefficients () by using (3), and obtain signal after denoising by using the inverse wavelet transform to the new wavelet coefficients.

3. The Verification of Simulated Signal

In order to verify the noise reduction performance of the improved WATV (IWATV) method proposed in this paper, the following simulation signals are designed as follows:where is composed of the linear superposition of the harmonic signal and impulse signal . is the additive gauss white noise (AGWN) with the variance of 0.6. The signal length is set as 4096 and sampling frequency is selected as 4096. Figure 2 describes the simulation signal and its corresponding component. The actual noise reduction effect of the wavelet soft-threshold, WATV, and the method proposed in this paper is shown in the Figure 3. As can be seen in Figure 3, the proposed method in this paper has better capacity in maintaining the pulse characteristic of the original simulation signal and reducing the reconstruction error to the maximum extend comparing with other denoising methods.

4. Application to Bearing Fault Diagnosis

The proposed improved total variation wavelet denoising method is applied to the signal containing frequency modulated signal in this section. The fault bearing data from Case Western Reserve University Data Center Website [19] is used to verify the proposed method. The website provides a large collection of experimental bearing data related to faulty bearings, and among them the data for faulty bearings are acquired by electric spark machining. To date, the bearing faulty data of Case Western Reserve University has been accepted as a public data set to verify the validity of fault diagnosis methods.

The schematic diagram and a picture of the device [19] are shown in Figure 4, respectively. The test bearings support the motor shaft. Single point faults were introduced to the test bearings using electrodischarge machining with fault diameters of 21 mils. Vibration data was collected by using accelerometers, which were attached to the housing with magnetic bases. Accelerometers were placed at the 12 o’clock position at both the drive end and fan end of the motor housing. The motor speed is 1730 r/min. Therefore, the rotation frequency of the theoretical calculation is  Hz. The sampling frequency of the signal is 12 kHz and the sampling point is 15000. The fault feature frequency of inner ring is  Hz. The detailed parameter and fault feature frequency are listed in Tables 1 and 2, respectively. Figure 5 shows a time-domain waveform of the fault bearing, where the presence of an impact component can be inspected clearly.

To verify the validity of the method, the classical wavelet denoising method and the proposed method based on the improved WATV (IWATV) were introduced into analysis of the faulty bearing data. The function of Db5 is chosen as the wavelet basis function and the decomposition level is set as 3. The frequency domain of filtered bearing fault signal by wavelet denoising method is shown in Figure 6(a). It is important to note that the fault frequency 155.5 Hz () and the second harmonic frequency component 314 Hz () are visible. However, the third harmonic frequency 468 Hz () is difficult to inspect because its amplitude is smaller than others. The result provided by envelop analysis is also employed in Figure 6(b). Envelop analysis generates some unknown frequency components, which can affect the determination of characteristic spectral lines. Then, the result obtained by the proposed method is plotted in Figure 7. It is obvious that the rotation frequency including its frequency multiplication such as , and frequency multiplication of inner ring fault () can be clearly detected. Comparing with the results in Figure 6, the proposed method has significant advantages over the wavelet denoising methods for the faulty feature extraction.

5. Experimental Verification of Gear Fault Diagnosis

Gears are the key parts of the power transmission box and the rotating movement. The gear transmission is the most common transmission structural in the mechanical equipment. In the practical vibration monitoring and diagnosis, sensors are installed on the gear box (including bearing seat), and then the vibration signal of the gear transmission can be acquired. Because of the influence of the characteristics of the gear vibration and the complicated and varying environment, the vibration signal of faulty gear not only is complex but also has the characteristic of strong nonstationary and nonlinear. The reason can be summarized by following aspects: the gear meshing effect, amplitude modulation, phase modulation, a variety of signal components deriving from the complex interference excitation, transmission path and dynamic coupling, time varying operating conditions, and so on. In order to verify the effectiveness of the proposed method in this paper, the measured data from the gear transmission system in the metallurgical industry was analyzed.

The data is from a gearbox for F2 mill, in which the meshing state of the gears has a direct influence on the normal operation of the whole unit. The transmission structure and the layout of the measuring points are shown in Figure 8. A CSI2130 dual channel vibration meter is used to test the transmission. Sampling frequency is set to 2.56 kHz and the sampling length is 5120. Moreover, the vibration signals of the vertical direction and horizontal direction are collected at each point. The main parameters of the transmission system are shown in Table 3. It can be calculated that the rotation frequency of the I axis is  Hz, and the gear meshing frequency of I and II axes is  Hz.

The time-domain waveform of the vibration signal at measuring point #1 in horizontal direction was represented in Figure 9. After careful observation of the time-domain waveform, preliminary conclusions can be obtained that the shock characteristic was obvious and the measured signal was interfered by the noisy component. Thus, the traditional signal processing methods such as Fast Fourier Transform (FFT) may not work in feature extraction. Subsequently, the wavelet denoising method and the proposed method in this paper are employed to analyze the signal.

The discrete wavelet transform can decompose the signal into different scales. Thus, the vibration signal is decomposed by using db5 orthogonal wavelet basis function and the decomposition level is set as 3. The low frequency signal is obtained by reconstructing the approximate coefficient of the third layer. Figure 10 is a view of a reconstructed signal in the frequency spectrum analysis. Although the second harmonic frequency of rotation frequency can be identification in Figure 10(b), the spectral line has been interfered by others and the gear meshing frequency cannot be detected properly. Increasing the number of decomposition levels cannot outstand the failure frequency as well.

The faulty signal is then filtered using the method proposed in this paper; the time-domain graph and frequency spectrum analysis after denoising are shown in Figure 11. As mentioned above, the rotation frequency of I axis is 10.8 Hz, while the meshing frequency of gear is 238.3 Hz in theory. It should be noted that, due to the influence of fluctuations in motor speed and frequency resolution, the actual frequency components of collected signal cannot correspond exactly with the calculated frequency. From Figure 11, it can be seen that the rotation frequency is 10.63 Hz (), the second harmonic frequency is 20.94 Hz (), and its third harmonic frequency is 31.56 Hz () which can be detected by the proposed method. Besides, the gear meshing frequency 230 Hz () and its second harmonic frequency 460 Hz () can be detected in the frequency spectrum. It is worth mentioning that the phenomenon of sideband about the gear meshing frequency (230 Hz) can also be inspected, which was proved to be the gear failure. Comparing with Figures 10(b) and 11(b), the analysis results demonstrated that the proposed method has significant advantages in gear fault diagnosis.

Then, the gear box was opened for inspection. As shown in Figure 12, the examination result indicated that there are eight serious surfaces flaking on gear #1 and gear #2, having a partial break. The results indicated that the proposed method is capable of effectively identifying the fault characteristics of gearbox.

6. Conclusions

In this paper, the method of local variance analysis was utilized to achieve the adaptive estimation of the noise variance, based on the improved wavelet denoising method (IWATV) using nonconvex penalty function. The results of the numerical simulation signal indicated that the proposed method can effectively reconstruct the original useful signal. The fault bearing data from Case Western Reserve University Data Center Website and gearbox fault signal in metallurgical industries were analyzed. The results demonstrated that the proposed method can reduce the noisy components effectively for extracting the fault features comparing with the typical wavelet denoising and WATV method.

Competing Interests

The authors declare that there is no conflict of interests regarding the publication of this article.

Acknowledgments

This work was supported by the National Natural Science Foundation of China (nos. 51475339 and 51405353), the Key Laboratory of Metallurgical Equipment and Control of Education Ministry, Wuhan University of Science and Technology (no. 2015B11), and the State Key Laboratory of Refractories and Metallurgy, Wuhan University of Science and Technology (no. ZR201603). The authors also appreciate the free download of the original bearing failure data and one photo picture provided by Case Western Reserve University Bearing Data Center Website.