基于改进傅里叶模态分解和频带熵的滚动轴承故障诊断方法

刘俊锋，俞翔，万海波

1 海军工程大学 动力工程学院，湖北 武汉 430033

2 海军工程大学 舰船与海洋学院，湖北 武汉 430033

0 引 言

受强背景噪声、设备转速、摩擦和负载的影响，滚动轴承振动信号往往呈现信噪比低、非线性、非平稳等特点，故导致其故障特征难以识别[1-2]。为解决此类问题，相关学者将短时傅里叶变换(short-time Fourier transform，STFT )[3]、Winger-Ville分布[4]、小波变换[5]、希尔伯特-黄变换(Hilbert-Huang transform, HHT)[6]、变分模态分 解(variation mode decomposition，VMD)[7-8]等方法应用于滚动轴承故障诊断领域。然而，由于Heisenberg 测不准原理的限制，短时傅里叶变换的分辨率较低；Winger-Ville 分布则会产生无法消除的二次交叉项干扰；HHT 是一种自适应的高效时频分析方法，包括经验模态分解（empirical mode decomposition，EMD）和希尔伯特谱分析两个部分，但EMD 存在模态混叠和端点效应等问题[9]。

在此基础上，相关学者提出了大量改进方案。Wu 等[10]提出了集合经验模态分解（ensemble empirical mode decomposition，EEMD)方法，Smith[11]提出了局部均值分解(local mean decomposition，LMD)方法，程军圣等[12]提出了局部特征尺度分解（local characteristic-scale decomposition，LCD）方法。但是，这类算法存在因极值点拟合所导致的模态混叠和端点效应以及缺乏严格的数学基础等问题。虽然VMD 方法可以通过搜寻约束变分模型的最优解来实现信号的自适应分解，但其受预设模态数量和惩罚参数的影响较大，且计算复杂、耗时，故效率较低[13]。

Singh 等[14]基于傅里叶理论提出了一种基于零相位滤波器组的自适应时频分析方法，即傅里叶分解方法（Fourier decomposition method，FDM），其通过高频到低频或低频到高频的边界频率搜索，可以将一个有限长度的非线性和非平稳数据分解为若干个傅里叶本征模态函数（Fourier intrinsic mode function，FIMF）之和。FDM 具有自适应性、局部性、完备性和正交性等优点，但在实际强背景噪声下，其存在边界频率偏移、信号过分解、有效FIMF 选取困难等问题。郑近德等[15]通过改进FDM 的边界频率搜索方法，提出了自适应经验傅里叶分解（adaptive empirical Fourier decomposition，AEFD），但仍然存在耗时长、敏感分量选取困难等问题。

当滚动轴承发生表面损伤型故障时，其冲击作用将诱发系统的高频固有振动成分，而低频部分在强背景噪声下的信噪比则较低[16]。为此，本文拟提出一种基于改进傅里叶模态分解和频带熵的滚动轴承故障诊断方法：首先，利用频带熵（frequency band entropy，FBE）分析确定敏感频带的中心频率并确定敏感区间边界；然后，在敏感频带区间内对信号进行带限傅里叶模态分解；接着，根据FIMFs 与原信号的FBE 区域从属关系，选取可以反映故障特征的敏感FIMFs；最后，对选取的FIMFs 进行包络谱分析和故障特征提取，并开展滚动轴承仿真和实验，用以验证该方法的有效性和精确性。

1 理论分析

1.1 改进傅里叶模态分解

改进傅里叶模态分解（modified fourier mode decomposition，MFMD）方法是在FDM 的基础上，通过引入FBE 分析等方法来确定敏感频带，进而设置频率搜索的初始边界并在各区间内进行优化的带限傅里叶模态分解。MFMD 可以将任意能量有限的非线性和非平稳信号x(t)自适应分解为一系列FIMF之和，即x(t)=yi(t)+rI(t)，其中：t为时间；yi(t)为 FIMFs，i=1, 2, ···，I，其中I为模态个数；rI(t)为残留信号。

对任意满足Dirichlet 条件的有限长的非线性和非平稳的零均值实信号x(t)（t∈[t0,t0+T]，其中t0为起始时间，T为周期）：

1）构造x(t)的周期延拓并对其进行快速傅里叶变换。令xT(t)=x(t-kT)，其中k∈[-∞,∞]，k为周期个数。使得x(t)=xT(t)w(t)，其中：当t0≤t≤t0+T时，w(t)=1， 其 他 则w(t)=0。因 此，xT(t)的傅里叶级数展开的复数形式为

式中：n∈[-∞,∞]，n为波数；j为虚数单位；角频率ω0=2π/T；cn=xT(t)exp(-jnω0t)dt。通过对xT(t)进行快速傅里叶变换，即可得到其复系数F(f)=∫x(t)e-iftdt，其中f为频率。

2） 频率边界的准确性将直接影响傅里叶模态的分解结果，而FDM 边界频率搜索受背景噪声的影响较大，为了取得理想的分解效果，本文提出一种基于FBE 的边界频率搜索方法，其确定频率边界集的步骤如下：

（1） 首先，对原始信号进行FBE 分析和短时傅里叶变换，选取区域熵值最小点作为中心频率，而敏感频带边界则由离中心频率最近的频带熵包络极大值点所确定。

（2） 然后，依据敏感频带边界，划分整个待搜索频带内的初始边界频率集{Bs}=[fs-1,fs)=[0,Fs/2)，其中：s=1, 2,···，S，S为初始区间个数；最小值f0=0；最大值fS=Fs/2 ， 其中Fs为频率上限。

（3） 根据初始边界频率集，在各区间内进行边界频率二次搜索，判断标准为在满足瞬时幅值ai(t)≥0 和 瞬时频率fi(t)≥0的情况下获取最小数量的解析FIMFs。设定最终优化的频率边界集{Bi}=[fi-1,fi)=[0,Fs/2)， 其中最小值f0=0，最大值fI=Fs/2，然后依据边界频率对其复系数的实部Re{F(f)}进行自适应分割。

3） 对信号在区间Bi=[fi-1,fi)内进行逆快速傅里叶变换，得到每个区间Bi内的解析FIMF分量Ii(t)=ai(t)exp(jφi(t))， 其中：ai(t)为瞬时幅值；φi(t)为瞬时相位。

因此，原始信号可以表示为

其离散形式为

式中，x[n]，ai[n]， φi[n]分 别为xT(t)，ai(t) ， φi(t)的离散形式，其中n=1,2,···,N，N为离散信号长度。

4） 每个FIMF 的瞬时幅值ai(t)和瞬时频率fi(t)均为有关时间的函数，故定义三维时频能量分布 {t,fi(t),ai(t)}为傅里叶希尔伯特谱，记为H(f,t)， 其边际希尔伯特谱h(f)为

1.2 频带熵

频带熵是一种时频分析和信息熵相结合的信号分析方法[17]，基于幅值谱熵的频带熵计算方法如下：

1） 首先，对信号y(z)（z=1, 2, ···,Z，其中Z为信号长度）做短时傅里叶变换，得到其时频分布TER为

式中：M为频率点的个数；C=Z/L，为傅里叶变换次数，其中L为步长；rM.C为信号在第C个窗口对应局部时间内频率M成分的估计值。

2）将第q个频率分量的幅值沿时间的变化定义为Xfq=(rq,1,rq,2,···,rq,C)，则单个频率分量的频带熵为

式中：m=1, 2, ···,C；q=1, 2, ···,M；pm,q为 第q个频率分量在整个频谱中的占比；Hsq为第q个频率分量的频带熵值；Fm为 频率分量Xfq沿时间轴的谱分布。

3） 计算各频率分量的频带熵值，即可得到全频带的频带熵分布Hsf为

如果频率分量Xfq随着时间平缓变化或规律变化，该频率分量的频带熵值将较小，反之则较大，故可在故障诊断中用于寻找设备的共振频率[18]，并为自适应滤波参数设置提供参考。

2 仿真信号分析

本节将通过仿真信号分析验证改进傅里叶分解的可行性，并通过轴承故障仿真信号分析阐明基于MFMD 和FBE 的故障诊断方法的优越性。

2.1 改进傅里叶模态分解

设定模拟信号为

式中：模拟信号x(t)由3 个调幅调频的时变模态信号x1(t)，x2(t)，x3(t)叠加而成；采样频率Fs为4 096 Hz；n(t)为信噪比为5 dB 的高斯白噪声。x(t)的时域波形和频谱如图1 所示。

图1 模拟信号的时域波形与频谱Fig. 1 Time domain waveform and spectrum of analog signal

分别对x(t)进行FDM 和MFMD，其中FDM 采用低频到高频的频率边界搜索。由于MFMD 本质上是一种带限傅里叶分解方法，故首先依据信号频谱确定仿真信号的初始分解频率边界为[40，90，200]，并在初始边界下进行优化傅里叶模态分解。选取算法处理所得的前5 个分量，分别命名为FIMF1，FIMF2，FIMF3，FIMF4，FIMF5，结果如图2 所示。受噪声信号n(t)影响，FDM 搜索得到的边界频率有所偏移，其在低频部分无法获取准确的信号模态信息，故分解所得的FIMFs 与预期结果的误差较大，而MFMD 则获取了3 个与其模态成分相对应的FIMFs，其分解效果较为理想。

图2 FDM 和MFMD 方法的信号分解结果Fig. 2 Signal decomposition results of FDM and MFMD methods

进一步通过傅里叶希尔伯特谱分析其时频分布特性，结果如图3 所示：FDM 分解的FIMFs 整体波动较大，在强背景噪声下存在严重的边界频率识别误差；MFMD 方法兼具傅里叶分解的完备性和正交性，其获取的模态分量符合预期，精确度较高。

图3 模拟信号的时频分布Fig. 3 Time-frequency distribution of analog signals

2.2 基于MFMD 和FBE 的故障特征提取

为了验证MFMD 和FBE 故障特征提取方法的可行性和优越性，设定滚动轴承内圈的故障仿真信号为

式中：x(t)为模拟信号；Ae为第e次冲击的调幅信号，其中e=1, 2, ···,E，E为最大冲击次数；τe为第e次冲击的微小波动；z(t)为高斯白噪声，信噪比为-10 dB；A0为冲击幅值；fr=28 Hz，为转频；h(t)为调频信号；B=500，为系统衰减系数；fn=4 000 Hz，为结构共振频率。系统采样频率Fs=12 000 Hz，分析点数为12 000，内圈故障频率f1=1/T=80 Hz。

仿真信号的时域波形及其包络谱如图4 所示，由于强背景噪声的影响，频谱包络图中无法直接提取故障特征频率、倍频或转频。FDM 分析结果如图5 所示，由于过分解、边界频率偏移和有效FIMF 选取困难等问题，故难以提取滚动轴承的故障特征。

图4 轴承故障仿真信号的时域波形与频谱Fig. 4 Time domain waveform and spectrum of bearing fault simulation signal

图5 轴承故障仿真信号的时频分布Fig. 5 Time-frequency distribution of bearing fault simulation signal

因此，本文提出基于MFMD 和FBE 的故障特征提取方法，其流程如图6 所示。

图6 故障诊断流程Fig. 6 Fault diagnosis process

通过对原始信号进行FBE 分析，即可得到窗长分别为16，32，64，128，256 时的频带熵分布，如图7 所示，可知其熵值最小点在4 000 Hz 处，这表明轴承固有频率也在该值附近。当旋转机械发生表面损伤故障时，其冲击将诱发系统的高频固有振动成分，从而放大故障特征，因此，本文将选取以熵值最小点作为中心频率的区域为敏感频带区间，其频率区间边界则由离中心频率最近的频带熵极大值点所确定。本文选取的敏感区间为3 500～4 500 Hz，并在敏感区间内进行信号优化傅里叶模态分解，结果如图8 所示。从时域波形上观察，FIMF4 的周期性冲击故障特征非常明显。

图7 原始信号不同窗长下的频带熵分析Fig. 7 Frequency band entropy analysis of original signal with different window lengths

图8 敏感区间内的信号分解结果Fig. 8 Signal decomposition results in the sensitive interval

为了选取敏感FIMF，进一步对各分量做FBE 分 析，如 图9（a）所 示，结 果 表 明：FIMF3，FIMF4 的频带熵分布在固有频率4 000 Hz 处与原始信号存在较高的区域从属关系，其中FIMF4 的从属特征更明显。为了验证敏感FIMFs 的选取准确性并提取故障特征，对FIMFs 进行包络谱分析，如 图9（b）所 示，结 果 表 明：FIMF1，FIMF2，FIMF3，FIMF4 均表征出基频为f1=1/T=80 Hz 的故障特征频率，其中FIMF4 的故障特征最为清晰，这也符合基于FBE 分析的敏感FIMFs 选取判断结果。

图9 FIMFs 的FBE 和包络谱分析Fig. 9 FBE and envelope spectrum analysis of FIMFs

为了验证本方法在强背景噪声下提取故障特征的优越性，本文将对比应用离散小波变换和经验模态分解处理该故障仿真信号，结果如图10 和图11 所示。在强背景噪声下，离散小波变换和EMD 均存在有效分量选取困难和故障特征受噪声信号影响明显等问题，因此，基于MFMD 和FBE的故障特征提取算法具有较高的可行性和优越性。

图10 小波包络谱的分析结果Fig. 10 Results of wavelet envelope spectrum analysis

图11 EMD 分解和包络谱的分析结果Fig. 11 Results of EMD decomposition and envelope spectrumanalysis

3 实验信号分析

本节将针对轴承故障模拟平台的实验轴承进行故障诊断，如图12 所示。轴承型号为NSK7010C，外径为80 mm，内径为50 mm，接触角为15°，滚动体直径为8.7 mm，滚动体共计19 个。轴承故障为人工模拟故障，利用激光在外圈加工了一个平行于轴承轴线，宽0.5 mm、深0.5 mm 的槽。电机转频为50 Hz，测试过程的采样频率为65 536 Hz，轴承外圈故障的理论特征频率为412 Hz。

图12 轴承故障模拟平台与模型Fig. 12 Bearing fault simulation platform and model

基于该故障模拟平台，采集实验轴承内侧的径向振动加速度数据，分析点数为65 536。首先，对原始振动信号进行短时傅里叶变换（STFT）和FBE 分析，如图13 和图14 所示，结果表明：轴承在低频部分的振动分量较复杂，其能量分散在较宽的频带范围内，故难以选取轴承固有频率所在的频带区域；振动信号在20 000 Hz 左右的频带范围内具有较高的能量，且频带熵在此范围内具有区域下降趋势，这符合固有频率的判断标准。由此，可判断实验轴承在20 000 Hz 左右具有某一阶固有频率。

图13 轴承振动信号的短时傅里叶变换Fig. 13 Short time Fourier transform of bearing vibration signal

图14 轴承振动信号的频带熵分析Fig. 14 Frequency band entropy analysis of bearing vibration signal

根据FBE 分析，即可进一步确定敏感频带区间为19 000～21 000 Hz，在此区间内对信号进行优化傅里叶模态分解，可以获得4 个FIMFs，如图15 所示。为了选取敏感FIMFs，对各分量做FBE 分 析，如 图16(a) 所 示，结 果 表 明：FIMF1，FIMF2，FIMF3，FIMF4 的频带熵分布在20 000 Hz附近均与原始信号存在区域从属关系，其中FIMF2 的从属特征最明显。为了验证敏感FIMFs的选取准确性，对FIMFs 做包络谱分析，如图16(b) 所示，结果表明：FIMF1，FIMF2，FIMF3，FIMF4 均表征出基频为416 Hz 的故障特征频率，与理论结果的412 Hz 基本吻合；FIMF2 的故障特征最为清晰，可以观察到一阶和二阶故障特征频率，符合上文基于FBE 分析的敏感FIMFs 选取判断结果。由此可见，基于MFMD 和FBE 的机械故障特征提取算法具有较高的可行性和精确性。

图15 敏感区间内的信号分解结果Fig. 15 Signal decomposition results in the sensitive interval

图16 FIMFs 的FBE 和包络谱分析Fig. 16 FBE and envelope spectrum analysis of FIMFs

4 结 论

本文提出了一种基于MFMD 与FBE 的旋转机械故障诊断方法，适用于多分量和强背景噪声下的滚动轴承早期故障诊断，主要结论如下：

1） 针对傅里叶分解在整个频率范围内搜索频率边界耗时长、精确性差、抗噪性能差和过分解等问题，提出了基于FBE 的故障特征敏感频带选取方法和基于初始敏感频率区间的MFMD 算法。仿真和实验信号分析结果证明：MFMD 的效果优于FDM、小波变换和EMD，具有较高的有效性和精确性，同时兼具傅里叶分解的自适应性、局部性、正交性和完备性等优势。

2） 针对敏感分量选取困难的问题，提出了FBE 和包络谱分析相结合的敏感FIMFs 选取方法。首先，依据FIMFs 与原信号的频带熵分布在固有频率附近的相似关系，选取敏感FIMFs；然后，通过包络谱分析即可提取故障特征，进而对敏感分量选取进行验证。