# Fault Diagnosis for Gearbox Based on Improved Empirical Mode Decomposition

Authors: Ling Zhao; Darong Huang; Yi Qin
Journal: Shock and Vibration
 (2015)

---

## Abstract

The application of the improved empirical mode decomposition (EMD) theory in gearbox fault diagnosis has been studied in this paper, and the transient features of gearbox vibration signals are shown. Based on using EMD, an improved algorithm of orthogonal empirical mode decomposition (OEMD) is put forward and is applied to extract the fault feature. Finally, fault diagnosis application in a gearbox is used as an example to prove the feasibility of the proposed method.


---

## Body

## 1. Introduction

Gearbox is the most important part of rotating machinery, which covers a broad range of mechanical equipment and plays a significant role in industrial applications. It generally operates under tough working environment and is therefore subject to faults, which may cause machinery to break down and decrease machinery service performance such as manufacturing quality and operation safety. One of the commonly used strategies in fault diagnosis of rotating machinery is adopting effective signal processing techniques to analyze the response signals and to reveal fault characteristics. However, traditional signal processing techniques, including time-domain and frequency-domain analysis, are based on the assumption that the process generating signals is stationary and linear. When transmission component faults occur in machinery, the dynamic behavior of the component is shown as nonlinear. The vibration signals have some characteristics, such as being nonstationary and nonperiodic, and the transient characteristics are especially obvious. Therefore, the need of signal processing techniques fit to the nonstationary and nonlinear signal is increasing in recent years.Equipment vibration signals not only include vibration signals, but also contain many background signals and noise signals, which contain large amounts of energy and are related with their own properties of each transmission component [1]. And when the bands of these signals overlap, it is difficult to extract the active fault information using the traditional analysis methods of time and frequency domain [2]. Therefore, in order to get accurate fault features, the noise signals and background signals must be removed from the vibration signals. For a gearbox, the nonlinear characteristics of its signals are very prominent and the filtering effects are not sufficient by means of the traditional linear filtering. Furthermore, the traditional method can easily cause signal distortion.Empirical mode decomposition (EMD) is one of the most powerful time-frequency analysis techniques and is based on the local characteristic of a signal. With EMD, the signal can be decomposed into a set of complete and almost orthogonal components called intrinsic mode function (IMF). EMD has a unique instantaneous ability of identification and the local characteristics scale which is based on signals that can decompose the complex vibration signals into a finite number of IMFs, which can reflect the local characteristics of the nonstationary signals [3]. The IMFs indicate the natural oscillatory mode imbedded in the signal and serve as the basis functions, which are determined by the signal itself, rather than predetermined kernels. Therefore, it is a self-adaptive signal processing technique and suitable for nonlinear and nonstationary processes. Since EMD was introduced in 1998, it has been extensively studied and widely applied to numerous fields, such as signal detection, bridge fault detection, and medical signal detection [4–6]. These applications can simultaneously demonstrate its unique advantages.However, currently the Hilbert Huang transform (HHT) signal analysis scheme is only an immature method, and it is inevitable that its theoretical basis and theoretical framework need to be further refined [7]. To a large degree, the analysis quality of this method depends on the EMD decomposition quality, but there are various factors that influence the EMD decomposition quality. One of the common causes is the poor performance of the EMD decomposition, which is influenced by the component features and noise [8]. Now many new methods, such as the mirroring closed continuation method, the continuation method based on neural networks, the continuation method based on support vector machines, the waveform feature matching the continuation method, and the mask signal method [9–13], are put forward by many scholars to solve the problem of the mode mixing or end effect; however, limitations still exist for each method that exists.For the multicomponent signals and the signals containing noise, mode mixing and false components may occur, if these signals are decomposed directly by means of EMD. Therefore, in this paper a method of extracting transient features of a nonstationary signal without mode mixing is proposed based on the combination of orthogonal empirical mode decomposition (OEMD) and Hilbert transform.


## 2. Empirical Mode Decomposition

EMD is used to decompose the signals into IMF with different scale characteristics, and these functions need to satisfy two conditions [14] as follows. (1) The number of extreme points is equal to the number of the zero-crossing points or their difference is one. (2) The mean value of the envelope, which is constituted by the maximum and the minimum values of dates, equals zero for any point. The signal filtering process is represented as follows: the mean value of the upper and lower envelope is determined in accordance with all the minimum points and maximum points of the analyzed data sequence x(t). Namely, (1)ht=xt-mt,where x(t) stands for the analyzed data sequence, m(t) represents the mean value of the upper and lower envelope, and h(t) is regarded as the new analyzed data sequence x(t). Then the first-order IMF is obtained by repeating the above operations until h(t) satisfies the two conditions mentioned above. c1(t) denotes the component with the highest frequency of the signal data sequence and it can be represented as follows:(2)c1t=ht.Quasiresidual component function r(t) is achieved by isolating the IMF as follows: (3)rt=xt-c1t,where r(t) is not always a stationary data sequence, so a repeated operation like above is needed so that the second IMF is obtained. The process as expressed above should be repeated until the final data sequence rn(t) cannot be decomposed, while rn(t) stands for the tendency or mean value of the analyzed data sequence x(t). Figure 1 shows the steps of the original EMD method.Figure 1
Flow chart of EMD.The detailed algorithm of the EMD is designed as below.Step 1.
Setr0(t)=x(t), i=1.Step 2.
Theith IMF is obtained by the following.(a)
Leth0(t)=ri(t), j=1.(b)
Find out the local extreme points ofhj-1(t).(c)
The maximum and minimum value points ofhj-1(t) are interpolated so that the upper and lower envelope are formed.(d)
Compute the mean valuemj-1(t) of the upper and lower envelope.(e)
Considerhj(t)=hj-1(t)-mj-1(t).(f)
If the judgment criteria of IMFs are contented, then setci(t)=hj(t); otherwise set j=j+1 and then return to (b).Step 3.
Considerrit=ri-1t-cit.Step 4.
If the number of the extreme points ofri(t) is less than two, let i=i+1 and then turn to Step 2; otherwise ri(t) is the residual component and the decomposition is complete.From the process above we can conclude that each mode has the identical number of extreme points within the entire length of the entire signals. The local mean value is defined by the upper and lower envelope of the signals. On this basis, the different mode components can be distinguished according to the characteristic scales of the signals (the time span between the adjacent extreme points of the signals). IMFs are no longer limited into the narrow-band signals, and simultaneously they can be presented as amplitude modulations and frequency modulations. Generally, associated with the addition of the order of IMFs, the characteristics scale of IMFs will increase. And consequently the multiorder IMFs sequence, whose time scales are varying from small to large, is achieved. From a spectrum angle, the signals are filtered from high frequency to low frequency.We compare EMD with classical time-frequency analysis methods, such as short time Fourier transform (STFT) and wavelets as follows.STFT adopts the same window for the whole signal, which makes it produce constant resolution for all frequencies. Although STFT can overcome the disadvantages of FFT-based methods in processing nonstationary signals, it still cannot obtain a good frequency resolution using wide windows and good time resolution (narrow window) at the same time. Therefore, STFT is suitable for the analysis of stationary signals instead of real nonstationary signals.Comparing with STFT, wavelets can be utilized to analyze multiscale signals through dilation and translation and extract time-frequency characteristics of the signals effectively. That means wavelets are more suitable than STFT for analyzing nonstationary signals. However, wavelets are nonadaptive and have their own disadvantage that their analysis results depend on the choice of the wavelet base function. This may lead to a subjective and a priori assumption on the characteristics of the signal.Different from wavelets and short time Fourier transform (STFT), EMD is a self-adaptive signal processing method. It is based on the local characteristic time scales of a signal and could decompose the signal into a set of IMFs. The IMFs represent the natural oscillatory mode embedded in the signal and work as the basis functions, which are determined by the signal itself, rather than predetermined kernels. Although the EMD method shows outstanding performance in processing nonlinear and nonstationary signals, the algorithm itself has some weaknesses. For example, mode mixing sometimes occurs between IMFs; the IMFs are not strictly orthogonal to each other. In conclusion, each time-frequency analysis method suffers from various problems. It is hard to say that one can always exceed the others for any case.


## 3. Mode Mixing Analysis

As discussed, after EMD, a complex signals(t) can be decomposed into several IMFs cjt and a residual component rit, which is expressed by (4)st=∑j=1icjt+ri(t).An IMF is defined as any function satisfying the following conditions. (1) In the whole dataset, the number of extrema and the number of zero crossings must either be equal or differ by at most one. (2) At any point, the mean value of the envelope defined by the local maxima and the envelope defined by the local minima is zero.Theoretically, EMD is perfect tool to separate monocomponents. But when high frequency components in signal create intermittence, the decomposition by EMD produces mode mixing [15]. And if the frequencies of two individual components lie close to each other, often within an octave, EMD cannot separate them too. The situation in which the IMF contains different time scales appearing after the signals have been through screening for mechanical vibration signals is called mode mixing, and the mechanical vibration signals stand for the signals containing multicomponents and complex noise. In other words, the various mode components cannot be effectively isolated according to the characteristics scales, making a certain IMF contain multiple modes, and consequently the intrinsic nature of the signals cannot be clearly reflected.There are two causes of the problem:(1) the IMF contains signals of widely disparate scales or (2) signals of a similar scale reside in different IMF components. To decompose (5)st=sht+sltby the EMD, where sht and slt denote the high frequency component and the low frequency component, respectively, the problem occurs when parts of the components of  sht are completely immersed in the corresponding parts of slt. Thus, the sifting process detects the extrema of slt in the parts of sht that are immersed. However, the sifting process also detects the extrema of sht in other parts of the signal st. This phenomenon can create intermittency in the extrema of the signal; that is, the extrema detected in the sifting process belong to different signals. As a consequence, the IMF is a signal comprised of different scales; some parts are from signal of high frequency component and the other parts are from signal of low frequency component.For the signals formed byxt=sin⁡100πt+sin⁡180πt+sin⁡(440πt) in 0~1 sec with the sample frequency 512 Hz, that means the signal is composed of three components with frequencies 50 Hz, 90 Hz, and 220 Hz, respectively. The original signal and frequency components can be achieved, as shown in Figure 2. As can be seen from Figure 1, x(t) can be decomposed into three simple harmonic motions, several small amplitude IMFs caused by some algorithm errors, and its residues. The IMF components of the original signal and decomposed signal can be achieved; the first four IMFs can be seen from Figure 3.Figure 2
Time-domain waveform of simulation signalx(t) and the frequency components.Figure 3
EMD effect and Hilbert spectrum of IMF1–IMF4.When the frequencies ratio of the two signal components varies from 0.5 to 2, it can easily be interpreted that the signal is made up of some IMFs with certain degrees of modulation according to the standard algorithms of EMD, which deviate from the nature of the original signals. Since the frequency ratio 90/50 = 1.8 < 2 lies within the same octave, which means normal EMD cannot separate the two modes, the mode mixing clearly appears in the IMF component IMF2. IMF2 contains not only 90 Hz component but also part of 50 Hz component during the intermittent period. Although IMF3 contains the intermittent 50 Hz component, the result of IMF2 is severely distorted.


## 4. Orthogonal Empirical Mode Decomposition

Essentially, empirical mode decomposition is the nonstationary signals that are treated axisymmetrically, separated the intrinsic mode functions according to the frequency in descending order. The intrinsic mode functions need to satisfy the two following conditions:(1) the number of extreme points equals the number of the zero-crossing points or their difference is one; (2) the mean value of the envelope which constituted the maximum and the minimum values of the dates equals zero for any point.It is difficult to meet the second condition mentioned above, so various criteria are proposed by different researchers to achieve ideal results. Aimed at the second condition, the criterion is proposed as follows: the ratio of the signals’ local mean curve and the energy of the signal is smaller than a threshold, where the threshold is determined by the characteristics of the components contained in the signals and the decomposition accuracy of the signals.Generally, the signals need to filter the IMF by EMD. Actually, the empirical mode decomposition is a process of adaptive filtering, obtaining the intrinsic mode functions by the adaptive band-pass filter [16], and the specific decomposition process is shown below.First, let the signal be detectedx(t) with a frequency range of [0,fc]. Then the frequency f1 can be found by searching within the frequency range above, making x(t) pass through the ideal band-pass filter. The filter is(6)H1f=1,f1≤f<fc,0,else.Then x(t) can meet the definition of IMF, therefore, obtaining the first IMF and letting it be c1(t). Second, the frequency f2 is searched within [0,f1], obtaining the second IMF or c2(t). The frequency f3 can be detected within [0,f2], achieving c3(t). Then the IMFs c1(t),c2(t),…,cn(t) and the residual function r(t) can be obtained by repeating the operation above, and the decomposition process ends when fn=0, r(t) is a monotonic function, or the energy of r(t) is small enough:(7)Hd(f)=1,0≤f<fc,0,else.Let fn=0; then r(t)=0. In order to get the frequency-modulated (FM) signals and amplitude-modulated (AM) signals, Hi(f) should be the filter with the widest bandpass.It can be seen that when obtaining thekth intrinsic mode function ck(t), the band-pass filter used can be presented as follows:(8)Hk(f)=1,fk≤f<fk-1,0,else.Whenk=1, f0=fc. Let Ck(f)=Hk(f)X(f), where X(f)=FTx(t), FT denotes the Fourier transform, and Hk(f) represents thekth filter. Then thekth IMF ck(t) can be expressed as IFTCk(f), so that ck(t)=IFTCk(f), where IFT indicates the inverse Fourier transform.For any two intrinsic modes function componentscj(t) and ck(t), according to the Parseval theorem, the following formulation can be obtained:(9)∫-∞∞cj(t)ck*(t)dt=∫-∞∞Ci(f)Ck*(f)df=∫-∞∞HifHkfCk*fdf,where j<k. Consider fj≥fk-1, so Hj(f)Hk(f)=0, and then the following formula is successfully obtained:(10)∫-∞∞cj(t)ck*(t)dt=0.Obviously, any two componentscj(t) and ck(t) are mutually orthogonal and its one component has nothing to do with the other. Moreover, its conversion does affect the represented function itself. Similarly, we can achieve that r(t) is orthogonal to any ck(t):(11)∑i=1nci(f)+Rf=∑i=1nHifXf+HdfXf=Xf,so (12)xt=∑i=1ncit+rt.All components, obtained by the decomposition mentioned above, are mutually orthogonal, so the decomposition process above is called OEMD. Meanwhile, formula (12) indicates the completeness of the method of OMED, where the sequential search method is adopted to search the frequency fi and the detailed algorithm flowchart can be seen in [7].It can be seen from the decomposition process above that every IMF can be obtained by applying OEMD to the signal; meanwhile, IMF is a result which can be obtained by filtering the original signal. Therefore, the resulting problem is not very serious and the previous components’ boundary effects will not have a later effect on the obtained component. The phenomenon of divergence and inward pollution, which appear in the data at both ends when the general EMD is adopted, will not happen, overcoming the end effect problem that exists in EMD.


## 5. Hilbert Transform

For an arbitrary time series,X(t), we can always have its Hilbert transform Y(t), as (13)Yt=1πP∫X(t′)t-t′dt′,where P indicates the Cauchy principal value. This transform exists for all functions of class LP. With this definition X(t) and Y(t) form a complex conjugate pair, so we can have an analytic signal, Z(t), as(14)Zt=Xt+Yt=ateiθt,in which(15)at=X2t+Y2t,(16)θt=arctanYtXt.In principle, there are infinitely many ways of defining the imaginary part, but the Hilbert transform provides a unique way of defining the imaginary part so that the result is an analytic function. A brief tutorial on the Hilbert transform with the emphasis on its physical interpretation was proposed by Bendat and Piersol. Essentially (13) defines the Hilbert transform as the convolution of X(t); therefore, it emphasizes the local properties of X(t). It is the best local fit of an amplitude and phase varying trigonometric function to X(t); in (14), the polar coordinate expression further clarifies the local nature. To conduct the Hilbert transform, there is considerable controversy in defining the instantaneous frequency as(17)ω=dθtdt.Theoretically, some limitations on the data are necessary, for the instantaneous frequency given in (16) is a single value function of time. At any given time, there is only one frequency value, so it can only represent one component, “monocomponent.” Unfortunately, no clear definition of the “monocomponent” signal was given to judge whether a function is or is not “monocomponent.” But we can use bandwidth to judge it.After obtaining the intrinsic mode function components, we will have no difficulties in applying the Hilbert transform to each component and computing the instantaneous frequency according to (16). After performing the Hilbert transform on each IMF component, we can express the signal. There are two steps: firstly, the signal is decomposed into several IMFs and a residual function by means of OEMD; secondly, the instantaneous frequency and the instantaneous amplitude are calculated for every IMF by Hilbert transform.Functionv(t) can be achieved by applying Hilbert transform to the signal u(t), and v(t) can be shown as follows:(18)vt=1π∫-∞∞uττ-tdτ.Then the analytical signal is constructed as follows:(19)zit=ut+jvt=aitejϕit.Thus, the amplitude function and phase function are all obtained and can be expressed by the following formulae (20) and (21), respectively:(20)ait=ui2t+vi2t,(21)ϕit=arctanvituit.At the same time, the instantaneous frequency can be obtained as follows:(22)fit=12πωit=12π·dϕitdt.Therefore, the instantaneous frequency and amplitude are used to depict the frequency of the signal instead of the power spectrum. So the Hilbert spectrum can be denoted as follows:(23)Hω,t=Re∑i=1naitej∫ωitdt,where Re means that the real part is adopted and the residual component ri(t) is ignored.However, Hilbert marginal spectrumh(ω) can be expressed as follows:(24)hω=∫0THω,tdt,where T indicates the total length of the signal. With time and frequency conversion, H(ω,t) can describe the signals’ conversion law over the entire frequency range. The instantaneous amplitude and instantaneous frequency are the variables of time, constituting the three-dimensional spectrum of time, frequency, and amplitude, that is, the Hilbert spectrum. However, with the frequency conversion, h(ω) reflects the signals’ transform situation over the entire frequency range.


## 6. Experimental Verification

### 6.1. Simulation Experiment

The method mentioned in this paper is applied to analyze the aforementioned simulated signal. The IMFs components and Hilbert spectrum are shown in Figure4.Figure 4
Main IMF and Hilbert spectrum of IMFs obtained by OEMD.As can be seen from Figure4, the sinusoidal signals with a frequency of 90 Hz and 50 Hz have been decomposed successfully, and they correspond to IMF2 and IMF3. Comparing Figure 3 with Figure 4, we can conclude that the method applied in this paper can not only eliminate the influence cause by noise but also avoid the shortcomings of the nonsingle component of IMF.


### 6.2. Vibration Feature Extraction Experiments of Gear Fault

The experiments are conducted by means of a gearbox failure detection experimental device in the Mechanical Transmission National Key Laboratory of Chongqing University. The gearbox is constituted by two pairs of gear pairings. The meshing frequencies of gear 1 and gear 2 are both 372 Hz, and the meshing frequencies of gear 3 and gear 4 are 228.5 Hz. Gear 1 with rationalf frequency of f1=14 Hz is installed on axis I; gear 2 and gear 3 with their rotation frequency of f2=20.4 Hz are both installed on axis II; gear 4 with rotation frequency f3=7.2 Hz is installed on axis III. The test point of the gear vibration signals is chosen on the bearing seat of axis I, and the sample frequency is 5 kHz. In order to get rid of the impacts of the higher harmonic wave, the maximum upper frequency of the analog filter is 1000 Hz.When the gearbox is in the fault state, its time-domain waveform is shown as Figure5. Not only are the impact compositions caused by the fault of the gearbox difficult to perceive, but also the failed gear position and failure mode cannot be determined because of the influence of noise. As can be seen from Figure 6 the main components of the original signals are shown as below. The frequency of the meshing vibration of the gear pair 1-2 is about 230 Hz, their frequency multiplication is about 466 Hz, and the mesh frequency of gear pair 3-4 is about 378 Hz.Figure 5
Time-domain signals of gearbox and description of gearbox signal spectrums.Figure 6
Main components of IMF achieved by means of EMD.The results obtained after processing the gearbox signals by means of EMD directly are shown as Figure6. IMF1 is the high frequency component of the sample signal, and some of it is shown in IMF2. IMF3 and IMF4 are the components of periodic vibration of gear. From the figure we cannot find the characteristic frequency to diagnose fault. This indicates that, because of the interference of the noise signals, the signals are difficult to meet the requirements of the analysis signal characteristics.The signals are dealt with in accordance with the method mentioned in this paper; the results of OEMD are shown in Figure7. It is clear that the signal mode mixing has been eliminated and the influence of the white noise has been reduced. It can be seen from Figure 7 that the edge frequency of the intervals with its rotational frequency of 20 Hz emerges surrounding the 458.5 Hz frequency multiplication and 369.3 Hz mesh frequency, which indicates that faults exit on axis II. Ultimately this proves that a broken tooth fault occurred on gear 2.Figure 7
IMF1–4 obtained after OEMD processing.


## 6.1. Simulation Experiment

The method mentioned in this paper is applied to analyze the aforementioned simulated signal. The IMFs components and Hilbert spectrum are shown in Figure4.Figure 4
Main IMF and Hilbert spectrum of IMFs obtained by OEMD.As can be seen from Figure4, the sinusoidal signals with a frequency of 90 Hz and 50 Hz have been decomposed successfully, and they correspond to IMF2 and IMF3. Comparing Figure 3 with Figure 4, we can conclude that the method applied in this paper can not only eliminate the influence cause by noise but also avoid the shortcomings of the nonsingle component of IMF.


## 6.2. Vibration Feature Extraction Experiments of Gear Fault

The experiments are conducted by means of a gearbox failure detection experimental device in the Mechanical Transmission National Key Laboratory of Chongqing University. The gearbox is constituted by two pairs of gear pairings. The meshing frequencies of gear 1 and gear 2 are both 372 Hz, and the meshing frequencies of gear 3 and gear 4 are 228.5 Hz. Gear 1 with rationalf frequency of f1=14 Hz is installed on axis I; gear 2 and gear 3 with their rotation frequency of f2=20.4 Hz are both installed on axis II; gear 4 with rotation frequency f3=7.2 Hz is installed on axis III. The test point of the gear vibration signals is chosen on the bearing seat of axis I, and the sample frequency is 5 kHz. In order to get rid of the impacts of the higher harmonic wave, the maximum upper frequency of the analog filter is 1000 Hz.When the gearbox is in the fault state, its time-domain waveform is shown as Figure5. Not only are the impact compositions caused by the fault of the gearbox difficult to perceive, but also the failed gear position and failure mode cannot be determined because of the influence of noise. As can be seen from Figure 6 the main components of the original signals are shown as below. The frequency of the meshing vibration of the gear pair 1-2 is about 230 Hz, their frequency multiplication is about 466 Hz, and the mesh frequency of gear pair 3-4 is about 378 Hz.Figure 5
Time-domain signals of gearbox and description of gearbox signal spectrums.Figure 6
Main components of IMF achieved by means of EMD.The results obtained after processing the gearbox signals by means of EMD directly are shown as Figure6. IMF1 is the high frequency component of the sample signal, and some of it is shown in IMF2. IMF3 and IMF4 are the components of periodic vibration of gear. From the figure we cannot find the characteristic frequency to diagnose fault. This indicates that, because of the interference of the noise signals, the signals are difficult to meet the requirements of the analysis signal characteristics.The signals are dealt with in accordance with the method mentioned in this paper; the results of OEMD are shown in Figure7. It is clear that the signal mode mixing has been eliminated and the influence of the white noise has been reduced. It can be seen from Figure 7 that the edge frequency of the intervals with its rotational frequency of 20 Hz emerges surrounding the 458.5 Hz frequency multiplication and 369.3 Hz mesh frequency, which indicates that faults exit on axis II. Ultimately this proves that a broken tooth fault occurred on gear 2.Figure 7
IMF1–4 obtained after OEMD processing.


## 7. Conclusions

The mode mixing appearing in the process of EMD is brought about by the IMF containing various time scales. The interference signals are eliminated by means of OEMD, avoiding the mode mixing and reducing the influence of the false frequency and noise. The results achieved by processing the gearbox vibration signals show that the method is practical and effective.



---

*Source: 215806-2015-09-21.xml*