Vold-Kalman滤波和高阶能量分离在时变工况行星齿轮箱故障诊断中的应用研究

秦嗣峰，冯志鹏，LIANG Ming

（1.北京科技大学机械工程学院，北京100083；2.渥太华大学机械工程系，安大略渥太华K1N 6N5，加拿大）

Vold-Kalman滤波和高阶能量分离在时变工况行星齿轮箱故障诊断中的应用研究

秦嗣峰1，冯志鹏1，LIANG Ming2

（1.北京科技大学机械工程学院，北京100083；2.渥太华大学机械工程系，安大略渥太华K1N 6N5，加拿大）

时变工况下的行星齿轮箱振动信号具有明显的时变调制特点，常规的频谱分析方法难以识别齿轮故障特征频率。提出了基于Vold-Kalman滤波和能量分离的时频分析方法，识别行星齿轮箱的时变特征频率，诊断齿轮故障。与传统的时频分析方法相比，基于Vold-Kalman滤波和能量分离的时频表示具有良好的时频分辨率，而且没有交叉项干扰，能够有效提取非平稳信号中的时变频率成分。分析了行星齿轮箱时变工况下的实验信号，准确地诊断了齿轮故障，验证了该方法的有效性。

故障诊断；行星齿轮箱；高阶能量分离；Vold-Kalman滤波；时频分析

引 言

行星齿轮箱在风力发电机、直升机和车辆等装备中应用广泛，实际运行中，其工况经常发生变化，产生非平稳的振动响应。目前的研究主要集中在恒定工况方面［1-2］，而针对时变工况方面的研究并不多见，其中，Williams和Zalubas［3］考虑直升机传动系统振动信号的非平稳特点，应用Wigner-Ville分布检测行星齿轮箱故障。Bartelmus和Zimroz［4-5］发现行星齿轮箱故障对载荷变化敏感，应用齿轮啮合振动幅值与工况之间的线性依赖关系监测时变工况下的行星齿轮箱运行状态。

行星齿轮箱故障诊断的关键是根据故障特征频率及其幅值的变化识别故障原因。行星齿轮箱啮合频率和齿轮故障特征频率与转速直接相关。工况（转速和负荷）的变化将引起转速波动，啮合频率和齿轮故障特征频率也会随之发生变化。因此，如何从非平稳信号中有效地识别齿轮故障特征频率及其时变特征，是时变工况下行星齿轮箱故障诊断的关键问题。

对于时变工况下的行星齿轮箱的非平稳振动信号，时频分析是一种有效的分析手段。但是，常见的时频分析方法各有优缺点［6-7］。例如：线性时频表示（如短时Fourier变换和小波变换）虽然不存在交叉项干扰问题，但是其时频分辨率受Heisenberg不确定性原理的限制，不能同时达到最佳；Wigner-Ville分布虽然时频分辨率高，却存在固有的交叉项干扰问题，不适合分析复杂多分量信号；以Wigner-Ville分布为基础的双线性时频分布（包括Cohen类分布和仿射类分布）通过各种核函数进行平滑处理，虽然抑制了交叉项干扰，但是会造成信号自项畸变，降低时频分辨率；Hilbert-Huang变换虽然对信号的时变具有自适应性，而且具有良好的局部时频聚集能力，但是对信号中的奇异点敏感，容易产生模式混淆，得到虚假的本质模式函数，影响分析结果。

行星齿轮箱振动信号具有明显的调制特点，在时变工况下，将表现出时变调制特征。非线性信号处理中的能量算子方法为分析这种时变调制信号提供了有效的途径。能量算子通过信号的瞬时幅值及其微分的非线性组合估计信号的“能量”，在其基础上提出的能量分离算法可以计算任意调制信号的幅值包络和瞬时频率［8-12］。能 量算 子 和 能 量分离在信号及其导数的基础上进行计算，算法简单、时间分辨率高，对信号的瞬态变化具有良好的适应性。但是，该方法仅适用于单分量信号。对于成分复杂的行星齿轮箱振动信号，需要将其分解为单分量再进行分析。

近年来提出的Vold-Kalman滤波方法［13-15］在最小化结构和数据方程误差的基础上，可以有效分解复杂多分量信号的频率成分，为解决频率交错重叠的行星齿轮箱振动信号单分量分解问题提供了一种有效手段。

针对时变工况下行星齿轮箱故障诊断面临的复杂非平稳振动信号分析问题，考虑Vold-Kalman滤波方法在单分量分解以及能量分离方法在瞬时频率计算方面的独特优势，提出了基于二者的时频分析方法。

1 Vold-Kalman滤波

Vold-Kalman滤波可以将复杂多分量信号分解为单分量。与传统的滤波方法相比，该方法直接从时域中提取感兴趣的信号分量，避免了由时域至频域变换带来的相位偏差；与经验模式分解方法相比，Vold-Kalman滤波的中心频率可以根据瞬时频率自适应调节，能够有效分离在时频域内邻近甚至交叉的信号分量［14-15］，避免 了 由经验 模 式分 解 方法带 来的模式混淆现象。

任意调制信号可以表示为

式中k为阶次，Ak（t）为第k个分量的幅值包络，Θk（t）为载波信号。

式中ω（τ）为瞬时角频率（τ）dτ为瞬时相位。

与载波信号相比，幅值包络以较低频率缓慢变化，可以用低阶多项式表示。对于离散信号

式中 ▽为差分算子，s为差分阶次，εk（n）为非齐次项。设多项式阶次为2，根据式（3）得

当n=1时，Ak（0）-2Ak（1）+Ak（2）=εk（1）。由于实际信号中，Ak（0）不存在，因此将Ak（0）视为0，可得-2Ak（1）+Ak（2）=εk（1）。同样，当n=N时，Ak（N+1）不存在，将Ak（N+1）视为0，可得Ak（N-1）-2Ak（N）=εk（N），式中N为采样点数。

将式（4）展开成矩阵形式

不同阶次的多项式具有相同的矩阵形式

式中M为N×N矩阵。

实际测试信号y（n）由各阶分量与噪声或误差δ（n）组成，则

式中y（n）是实测振动信号，δ（n）为噪声或误差项。

为了获得某一感兴趣的信号分量xk（n），由式（7）得

其矩阵形式为

根据常规的时频分析（如短时Fourier变换）的时频脊线，可以估计感兴趣信号分量的瞬时频率ωk（n），从而获得其载波矩阵C。

根据式（6）和（9），应用最小二乘法使非齐次项和噪声或误差项的平方和最小［16］，即

式中C*为C共轭转置；A*为A共轭转置；r为加权因子，决定了滤波器的频率分辨力。r越大，滤波器的分辨力越高，但幅值包络的收敛速度越慢；相反，r越小，滤波器的分辨力越低，但幅值包络的收敛速度越快。为了保证幅值包络的收敛性及计算的可行性，r的取值不能过大，其最大值如表1所示［17］。

表1 加权因子最大值Tab.1 The maximum value of the weighting factor

由式（10）求得各信号分量的幅值包络矩阵A为

由载波矩阵C和幅值包络矩阵A可以重构感兴趣的各信号分量

2 高阶能量分离算法

能量算子和能量分离方法完全由信号数据驱动，无需构造任何基函数，能够适应信号的瞬时变化，时间分辨率高，适合检测信号幅值和频率的非平稳变化特征，结合时频分析方法，能够分析处理任意时变调制信号。高阶能量算子不仅具有普通二阶能量算子的各种优点，而且对噪声干扰的鲁棒性更好［11，12］，因此 本 文应用高 阶 能量 算 子进行 能 量分离，计算信号的幅值包络和瞬时频率。

对于任意连续信号x（t），高阶能量算子定义 为［10］

式中k为高阶能量算子的阶次，为任意整数；x（k）为x（t）的第k阶导数

当k为负数时，x（k）（t）由微分变为积分，将削弱高阶能量算子检测瞬态成分的能力，因此本文只关注k为正数时的高阶能量算子。通常，高阶能量算子可由低阶能量算子通过递归推导得到

由式（13），（14）可知，当k=2时，γ2（x）=（˙x）2-x¨x，即为二阶Teager能量算子［8］；零阶能量算子为γ0（x）=x-x2；一阶微分能量算子对于任何信号都为零。在高阶能量算子中，三阶和四阶能量算子比较实用，分别定义为

其中

对于简谐信号x（t）=Acos（ωt+θ），可以得到高阶能量算子输出的递归方程

初始条件为E0=-A2，E1=0。由式（18）递推，可得

对于简谐信号x（t）=Acos（ωt+θ），由四阶和二阶能量算子可以估计信号的频率和幅值

对于调幅-调频信号

若幅值包络a（t）和瞬时频率ω（t）变化的速率和幅度相对于载波频率来说不大，则依然可以由四阶和二阶能量算子估计信号的时变频率和包络

3 分析过程

（1）应用短时Fourier变换等常规时频分析方法对信号进行初步时频分析，判断感兴趣的信号成分。根据时频脊线估计感兴趣信号成分的瞬时频率。

（2）根据已获得的瞬时频率曲线，构造其载波信号，从而获得感兴趣信号成分的载波矩阵。

（3）根据感兴趣信号成分的载波矩阵，应用Vold-Kalman滤波将感兴趣的信号成分分解为单分量，满足能量分离方法的要求。

（4）应用高阶能量分离方法计算各单分量信号成分xi（t）的四阶“能量”Ei（t）和瞬时频率fi（t）。

（5）根据各单分量信号成分的四阶“能量”和瞬时频率，构造信号的高阶时频“能量”分布

式中Ei（t）为第i个单分量在t时刻的四阶“能量”，fi（t）为第i个单分量在t时刻的瞬时频率，δ（·）为Dirac函数。

4 实验信号分析

4.1 实验说明

本实验在Ottawa大学机器监控实验室完成。图1为该实验室的风电机组行星齿轮传动实验台，其中行星齿轮箱和定轴齿轮箱均为两级结构，齿轮参数分别见表2和3所示。该实验台的具体传动路线为：电机→定轴齿轮箱→第1级行星齿轮箱→第2级行星齿轮箱。为了模拟行星齿轮箱的齿轮故障，在第1级太阳轮的某个轮齿上加工了剥落损伤，如图2所示。在行星齿轮箱箱体顶部安装了加速度传感器，采集行星齿轮箱的振动信号。在自然停机过程中采集了振动信号，采样频率为20 k Hz，行星齿轮箱第2级输出轴承受的负荷为16.284 N·m。为了计算的可行性，仅取60至40 Hz降速过程内的振动信号进行分析。

图1 风电机组行星齿轮传动实验系统Fig.1 Wind turbine planetary gearbox experimental system

表2 行星齿轮箱参数Tab.2 Planetarygearboxconfigurationparameter

表3 定轴齿轮箱参数Tab.3 Fixed-shaftgearboxconfigurationparameter

图2 太阳轮损伤Fig.2 Sun gear damage

根据行星齿轮箱和定轴齿轮箱的参数及电机转速，可计算出第1级行星齿轮箱齿轮特征频率与电机转速之间的关系：

啮合频率

太阳轮旋转频率

行星架旋转频率

太阳轮局部故障特征频率

行星轮局部故障特征频率

齿圈局部故障特征频率

4.2 实验信号分析

4.2.1 正常信号

图3为正常行星齿轮箱降速过程中驱动电机的转速变化曲线。在降速过程中，驱动电机的最高转速为60 Hz，相应地，第1级行星齿轮箱的最高啮合频率为222.222 Hz。由于高阶倍频处特征频率的变化规律与基频处相同，不失一般性，以下分析的频率范围为0～400 Hz，覆盖第1级行星齿轮箱最高啮合频率的3/2倍，振动信号中包含该齿轮箱的健康状态信息。

图3 正常状态电机转速Fig.3 Motor speed of normal gearbox

图4（a）为正常信号的短时Fourier变换谱图。从图中可以看出信号成分主要包括啮合频率fm、啮合频率与太阳轮故障特征频率之差fm-fs、以及啮合频率、太阳轮故障特征频率的2倍频和太阳轮旋转频率的组合fm±2fs-fsr，但是时频分辨率较低。

图4 正常信号Fig.4 Normal gearbox signal

应用Vold-Kalman滤波分离上述4个主要信号分量，然后应用高阶能量分离算法计算各分量的4阶能量和瞬时频率，最后应用各分量的瞬时频率和4阶能量构造高阶时频能量分布，结果如图4（b）所示。该方法得到的时频分布的分辨率高，而且没有交叉项干扰。可以清晰地识别出啮合频率fm、啮合频率与太阳轮故障特征频率之差fm-fs、以及啮合频率、太阳轮故障特征频率的2倍频和太阳轮旋转频率的组合fm±2fs-fsr，但它们的幅值均较小。齿轮制造和安装误差以及微小的损伤在所难免，导致振动信号中出现上述频率成分，这些现象符合正常信号的理论预期结果。

4.2.2 太阳轮故障信号

图5为太阳轮故障状态下驱动电机的转速变化曲线。在降速过程中，驱动电机的最高转速仍为60 Hz，因此，以下分析频率范围仍为0～400 Hz。

图5 太阳轮故障状态电机转速Fig.5 Motor speed of faulty sun gearbox

图6（a）为太阳轮故障信号的短时Fourier变换谱图。从图中可以看出信号成分主要包括啮合频率fm、啮合频率与太阳轮故障特征频率及其倍频的组合fm±nfs、以及啮合频率、太阳轮故障特征频率及其倍频和太阳轮旋转频率及其倍频的组合fm±nfs-nfsr等7个频率成分，但是时频分辨率较低。

应用Vold-Kalman滤波和高阶能量分离方法构造的高阶时频能量分布如图6（b）所示。啮合频率fm、啮合频率与太阳轮故障特征频率之差fmfs、以及啮合频率、太阳轮故障特征频率的2倍频和太阳轮旋转频率的组合fm±2fs-fsr依然存在，但在1～2 s阶段，啮合频率与太阳轮故障特征频率之差fm-fs的能量明显增强。此外，还出现了啮合频率与太阳轮故障特征频率及其3倍频之和fm+fs，fm+3fs以及啮合频率、太阳轮故障特征频率和太阳轮旋转频率2倍频的组合fm+fs-2fsr等3个频率成分。与正常信号相比，该信号的高阶时频能量分布中出现了更多的频率成分，且这些频率均和太阳轮故障直接相关，说明太阳轮出现了故障，符合实验中的实际情况。

图6 太阳轮故障信号Fig.6 Faulty sun gearbox signal

5 结 论

风力发电机组、直升机和车辆等装备中的行星齿轮箱的运行工况经常变化，振动信号具有非平稳性。针对时变工况下行星齿轮箱振动信号的时变调制特点，发挥Vold-Kalman滤波单分量分解和高阶能量分离时间分辨率高的优势，提出了一种具有良好的时频聚集能力、无交叉项干扰的时频分析方法。应用该方法分析了行星齿轮箱时变工况下的实验信号，准确识别信号中的时变频率成分，正确诊断了齿轮故障，验证了该方法的有效性。

［1］ 雷亚国，何正嘉，林京，等.行星齿轮箱故障诊断技术的研究进展［J］.机械工程学报，2011，47（19）：59—67. Lei Yaguo，He Zhengjia，Lin Jing，et al.Research advances of fault diagnosis technique for planetary gearboxes［J］.Journal of Mechanical Engineering，2011，47（19）：59—67.

［2］ Lei Yaguo，Lin Jing，Zuo M J，et al.Condition monitoring and fault diagnosis of planetary gearboxes：a review［J］.Measurement，2014，48：292—305.

［3］ Williams W J，Zalubas E J.Helicopter transmission fault detection via time-frequency，scale and spectral methods［J］.Mechanical Systems and Signal Processing，2000，14（4）：545—559.

［4］ Bartelmus W，Zimroz R.Vibration condition monitoring of planetary gearbox under varying external load［J］.Mechanical Systems and Signal Processing，2009，23：246—257.

［5］ Bartelmus W，Zimroz R.A new feature for monitoring the condition of gearboxes in non-stationary operation conditions［J］.Mechanical Systems and Signal Processing，2009，23：1 528—1 534.

［6］ 冯志鹏，范寅夕，Liang Ming，等.行星齿轮箱故障诊断的非平稳振动信号分析方法［J］.中国电机工程学报，2013，33（17）：105—110.Feng Zhipeng，Fan Yinxi，Liang Ming，et al.A nonstationary vibration signal analysis method for fault diagnosis of planetary gearboxes［J］.Proceedings of the CSEE，2013，33（17）：105—110.

［7］ Feng Zhipeng，Liang Ming，Chu Fulei.Recent advances in time-frequency analysis methods for machinery fault diagnosis：A review with application examples［J］.Mechanical Systems and Signal Processing，2013，38：165—205.

［8］ Kaiser J F.On a simple algorithm to calculate the ‘energy’of a signal［A］.Proceedings of IEEE International Conference on Acoustics，Speech and Signal Processing［C］.Albuquerque，USA.1990，1：381—384.

［9］ Maragos P，Kaiser J F，Quatieri T F.Energy separation in signal modulations with application to speech analysis［J］.IEEE Transactions on Signal Processing，1993，41（10）：3 024—3 051.

［10］Maragos P，Potamianos A.Higher order differential energy operators［J］.IEEE Signal Processing，1995，2 （8）：152—154.

［11］Diop E H S，Boudraa A O，Salzenstein F.A joint 2D AM-FM estimation based on higher order Teager-Kaiser energy operators［J］.Signal Image and Video Processing，2011，5：61—68.

［12］Salzenstein F，Boudraa A O，Chonavel T.A new class of multi-dimensional Teager-Kaiser and higher order operators based on directional derivatives［J］.Multidimensional Systems and Signal Processing，2013，24 （3）：543—572.

［13］Vold H，Herlufsen H.Order tracking at extreme slew rates，using Kalman tracking filters［J］.SAE Paper Number 931288，1993.

［14］Vold H，Mains M，Blough J.Theoretical foundations for high performance order tracking with the Vold-Kalman tracking filter［J］.SAE Paper Number 972007，1997.

［15］Vold H，Herlufsen H，Mains M.Multi axle order tracking with the Vold-Kalman tracking filter［J］. Journal of Sound and Vibration，1997，13（5）：30—34.

［16］Feldbauer Ch，Holdrich R.Realization of a Vold-Kalman tracking filter—A least squares problem［A］. Proceedings of the COST G-6 Conference on Digital Audio Effects（DAFX-000）［C］.Verona Italy，December 7-9，2000：1—4.

［17］Tuma J.Setting the pass bandwidth in the Vold-Kalman order tracking filter［A］.Twelfth International Congress on Sound and Vibration［C］.Lisbon.2005，Paper 719.

Application of Vold-Kalman filter and higher order energy separation to fault diagnosis of planetary gearbox under time-varying conditions

QIN Si-feng1，FENG Zhi-peng1，LIANG Ming2
（1.School of Mechanical Engineering，University of Science and Technology Beijing，Beijing 100083，China；2.Department of Mechanical Engineering，University of Ottawa，Ottawa，ON K1N 6N5，Canada）

The vibration signals of planetary gearboxes under nonstationary running conditions have significant time-varying modulation feature.Conventional spectral analysis methods are unable to identify the gear fault characteristic frequencies from such nonstationary signals.A method based on Vold-Kalman filter and higher order energy separation is proposed to analyze the vibration signals of planetary gearboxes under nonstationary conditions in time-frequency domain，thus to identify the timevarying characteristic frequencies and diagnose the gear faults.The time-frequency representation derived from Vold-Kalman filter and higher order energy separation provides nice time-frequency resolution and is free from cross term interference，and thus it can effectively pinpoint the time-varying constituent frequency components of nonstationary signals.The proposed method is validated with analysis of lab experimental signals of a planetary gearbox under time-varying running conditions.

fault diagnosis；planetary gearbox；higher order energy separation；Vold-Kalman filter；time-frequency analysis

TH165+.3；TH132.425

A

1004-4523（2015）05-0839-07

10.16385/j.cnki.issn.1004-4523.2015.05.20

秦嗣峰（1989—），男，硕士研究生。电话：18866226821；E-mail：qsf0312@163.com

冯志鹏（1973—），男，博士，教授，博士生导师。电话：13520978607；E-mail：fengzp@ustb.edu.cn

2014-05-13<； class="emphasis_bold">；修订日期：2；

2014-11-05

国家自然科学基金资助项目（11272047，51475038）和教育部新世纪优秀人才支持计划项目（NCET-12-0775）