基于同步提取变换与Vold-Kalman滤波的燃机阶次跟踪方法研究*

张文海，胡明辉*，江志农，冯 坤

(1.北京化工大学 发动机健康监控及网络化教育部重点实验室，北京 100029；2.北京化工大学 高端机械装备健康监控与自愈化北京市重点实验室，北京 100029)

0 引 言

目前，燃气轮机已在我国航空、船舶及石油化工等行业得到了广泛应用，但其在运行过程中故障频发。当燃机发生故障时，振动是最敏感的参数之一。例如转子不平衡、支承不同心及叶片断裂等故障发生后，由于对转子状态或流场的影响，测得的机匣振动信号必然发生变化。

因此，振动分析已经被作为燃机状态监测的重要手段。然而，目前我国主要通过机匣的低频振动速度总值来监测燃机振动状态，该参数可从振动能量角度较好地用于联锁保护，却难以实现振动故障的早期预警，亦难以向故障诊断提供直接的数据支撑。

由于燃气轮机常见的不平衡、不对中、不同心、叶片断裂等故障，其主要振动特征均为跟转速相关的谐频振动即阶次成分，因此，燃气轮机振动信号中的阶次分量是其振动故障预警和故障分析的主要依据。在燃机的运行过程中，包括稳定工况和启停车等非稳定工况，有必要对阶次振动进行跟踪提取。基于卡尔曼滤波的阶次跟踪(vold-kalman filter based order tracking，VKFOT)是旋转机械振动阶次分量提取方面应用较广的方法之一，其不仅可以用于计算阶次分量幅值，还可以得到各阶次分量的时域波形，为后续相位分析等阶次分量深入挖掘提供足够的数据支撑。但该方法需要振动信号及与振动同步的转速数据同时作为输入[1]。

然而，由于燃机转速波动较大，同时很多现场的转速数据和振动数据往往通过不同设备进行采集，即很难获得与振动完全同步的精确转速数据。在机匣处测得的振动信号中，工频成分与直接测量得到的转速频率存在差异。若将该转速与振动信号作为VKFOT方法的输入，则必然会因为转频偏差而导致提取到的阶次幅值误差较大，对于高阶成分误差更大。

1995年，Vold和Leuridan首次将经典最优滤波理论中Kalman滤波器应用于阶次跟踪上[2]，但该方法需要振动信号及与振动同步且相同采样率的转速数据同时作为滤波器输入，否则便会引发提取到的频率分量与目标不一致问题。针对转速与振动同步性较差或无转速输入，学者们首先提出通过设置合理的滤波器带宽来补偿的方法，进而确保阶次提取结果的准确性。TUMA[3]通过理论分析，给出了滤波带宽与共振频率的关系；WANG K S等人[4]通过建模仿真得到了共振频率计算滤波的优化带宽；冯珂等人[5]通过VKF跟踪的多次循环计算来进行滤波带宽的优化选取。以上方法均可用于补偿目标阶次偏差带来的VKF跟踪误差，但仅适用于结构较为简单、共振频率易于获得的场合，循环优选也不能从根本上解决目标阶次偏移的问题。另一方面，为获得同步的转速数据，将VKFOT方法应用于无转速信号的场合，文献应用了能量重心法频谱校正、短时傅里叶变换(STFT)及短时chirp-Fourier变换(STCFT)等获得瞬时转速[6-10]，但由于以上时频变换方法建立在局部平稳性假设上，不适用于燃机变负荷过程转速变化剧烈时的转速提取。于刚[11]在2017年提出了同步提取变换(SET)方法，极大地提高了非平稳非线性信号的瞬时频率分析精度，为快变转速场合转速提取提供了可能。

针对上述问题，为获取与振动同步的转速数据，笔者提出一种基于同步提取变换与Vold-Kalman滤波的燃机振动阶次跟踪方法。

1 方法基础理论

1.1 同步提取变换(SET)

同步提取变换是基于短时傅里叶变化(STFT)的时频分析方法，用于提取瞬时频率(IF)。

主要实现过程分如下3步。

(1)STFT时频谱计算

计算信号x(t)=A·eiω0t的标准STFT时频谱：

(1)

式中：g(u-t)—可移动的紧支窗函数。

根据帕塞瓦尔定理，对式(1)进行修正，进一步展开可得：

(2)

由式(2)可知，STFT是由一系列与谐波信号频率ω0一致的时频系数组成。在频率ω=ω0处，时频系数幅值最大，距离ω0越远时频系数幅值越小。

(2)STFT中的瞬时频率估计

根据时频谱的相位信息，可得到瞬时频率：

(3)

式中：∂tGe(t,ω)—Ge(t,ω)对时间t的一阶偏导数。

式(3)表示在二维时频平面中，只要存在Ge(t,ω)≠0，ω0(t,ω)即为瞬时频率ω0。

(3)瞬时频率精确提取

利用同步提取算子将瞬时频率处的发散能量去除，得到SET谱：

Te(t,ω)=Ge(t,ω)·δ(ω-ω0(t,ω))

(4)

式中：δ(ω-ω0(t,ω))—同步提取算子。

式(4)通过δ函数(单位脉冲函数)，仅保留时频脊线上能量最大(ω=ω0)的时频系数，从而得到时频聚焦良好的时频分析效果。

原始信号与SET谱的重构关系如下：

(5)

式中：Re{·}—取实部。

因此，r(t)为SET谱中信号的瞬时频率轨迹，即实现了在高分辨率SET谱提取瞬时频率。

1.2 Vold-Kalman阶次跟踪方法原理

VKFOT的基本思想是使用转速数据构建跟踪阶次的结构方程作为估计值，以传感器测量信号形成的数据方程作为观测值，通过令估计误差与观测误差协方差最小的方法，估计出该状态下的跟踪阶次信号。

(1)数值方程

设一个振动信号波形x(n)，该信号与其复包络信号x′(n)、幅值A(n)三者的关系[12]如下：

x(n)=A(n)cos[Θ(n)+φ]=
Re{x′(n)exp(jΘ(n))}

(6)

式中：φ—初相位；Θ(n)—参考轴的角位移。

对于VKFOT数值方程而言，阶次跟踪就是求解复包络x′(n)的过程，通过式(6)可以将幅值、复包络与波形进行相互转换。

VKFOT数值方程向量形式如下：

y-Cx=η

(7)

式中：y—原振动信号；x—目标阶次波形，二者均为复包络形式；C(n)—由转速计算求得的角位移；η—误差项，代表其他谐波分量与随机噪声成分的总和。

由于在数值方程中引入角位移，该方法被称为基于角位移的VKF_OT[13]。

(2)结构方程

用包含异构项ε(n)的低阶多项式表示的结构方程如下(异构项ε(n)形容目标阶次受到噪声干扰而引起的幅值频率等变化)：

(8)

结构方程向量形式为：

Ax=ε

(9)

(3)求解过程

VKFOT的求解过程是通过结构方程与数据方程联立，将其转化为最小均方根问题。

通过最小二乘法建立异构项与误差项的平方和，称为损失函数J：

J=r2εTε+ηTη

(10)

式中：r—权重因子，用于形容平衡误差项η与异构项ε占比。

权重因子r与跟踪带宽Bw的关系为：

(11)

式中：Δf=Bw/fs，fs—信号采样频率。

由于J的展开式中二次项系数为正值，求J对于xH的一阶导数为0时J取到最小值。求解可得：

x=(r2ATA+E)-1CHy

(12)

求解该方程即可得到跟踪谐波分量x′(n)。

2 SET与VKFOT结合的阶次跟踪方法

笔者提出的SET与VKFOT结合的阶次跟踪方法实现流程，如图1所示。

图1 SET与VKFOT结合的阶次跟踪方法

具体步骤如下：

(1)由于加速度信号中主要成分为处于高频段的叶片通过频率等成分，成分多、波动大，难以直接提取，而工频成分不明显，对加速度信号进行积分得到工频成分占主导的速度信号；

(2)对(1)中得到的速度信号进行同步提取变换，得到包含多个IF分量的SET时频图；

(3)结合采集所得的同步性较差转速信号，从(2)中SET时频图中确定工频脊线明显、无噪声成分的频率范围作为工频提取带宽，从SET时频图中提取工频信号，即可得到与振动同步的转速信号；

(4)根据跟踪目标阶次的大小，选择振动加速度信号或速度信号，以及SET提取得到的转速数据作为VKF滤波器输入，分别构造数值方程和结构方程；根据跟踪目标设定目标阶次O，如跟踪一阶分量则取O=1；根据时频图中目标阶次处无其他噪声成分最小范围确定跟踪带宽Bw；根据式(10)计算权重因子r；

(5)联立数值方程、结构方程、权重因子，构造损失函数并求其最小值，求解得到目标阶次时域波形。

3 实际案例分析

为验证该方法在燃机振动阶次分量提取中的有效性及优势，笔者选取了某型燃气轮机实际试车数据进行分析；使用该方法和原始VKFOT方法分别提取目标阶次时域波形，采用计算阶次跟踪(COT)的方法，分别将二者提取所得的目标阶次幅值与原信号阶次谱中读出的幅值进行对比分析。

(1)数据基本情况

某型双转子燃气轮机启动及平稳运行过程级联图如图2所示。

图2 启动过程级联图A—振动加速度；v—转速；f—频率

振动信号采样频率为fs=51 200 Hz，信号长度T=12 s。图2中标出的明显阶次成分频率8 326 Hz，为低压转子的57.3阶(低压转子直接测量转速为8 718 r/min)，对应低压压气机某级动叶通过频率(57阶)，但按照测量转速计算得到的阶次为57.3，与实际阶次数57存在较大偏差。下面对这一明显阶次成分进行阶次提取。

(2)基于振动信号的转速提取

笔者对振动加速度信号进行积分得到速度信号，进行同步提取变换后得到包含低压工频的SET时频，如图3所示。

图3 SET时频图T—时间

从图3可以看出：低压转子工频脊线明显，周围约6 Hz处无噪声成分。以直接测量得到的转速为中心频率，带宽6 Hz，从SET时频图中提取得到低压转子工频。

提取得到的转速与直接测量转速对比图如图4所示。

图4 实测转速与SET提取转速对比

从图4中可以看出：通过SET谱提取的转速与直接测量所得转速存在一定偏差。当跟踪高阶成分时，该误差会呈倍数放大。

(3)VKFOT目标阶次时域波形提取

将振动信号，SET谱提取到的工频信号作为VKF滤波器输入，设定目标阶次O为57阶，由于通过SET谱提取的转频成分避免了阶次偏移的问题，根据时频图中目标阶次处无其他噪声成分最小范围，确定滤波带宽Bw分别为0.85个阶次，最大程度上保留有效阶次成分，滤除噪声干扰。

得到目标阶次成分时域波形如图5所示。

图5 57阶目标阶次时域波形

(4)结果对比分析

笔者单独使用VKFOT进行分析，以原始转速及振动信号作为滤波器输入，设定相同的目标阶次与滤波器参数，得到目标阶次成分时域波形。

分别对原始振动信号、SET-VKFOT得到的57阶时域波形、单独VKFOT得到的时域波形以原始转速为参考做阶次谱，分别读取目标阶次成分(即57阶)幅值，得到阶次幅值趋势对比图，如图6所示。

图6 阶次幅值趋势对比

由图6可知：SET-VKFOT方法得到的阶次幅值趋势与原阶次谱幅值趋势基本吻合，而单独使用VKFOT则误差较大。计算图中两种方法与原阶次谱幅值趋势的绝对误差，SET-VKFOT误差在15 m/s2以内，而VKFOT误差最大达到了160 m/s2以上。

由上述对比分析可知：在变负荷及稳定工况下，SET-VKFOT方法提取阶次分量幅值误差较小，相对单独使用VKFOT方法有明显优势。

4 结束语

针对燃气轮机转速波动较大、测量转速与振动难以同步，进而导致难以提取振动信号中目标阶次幅值变化趋势的问题，笔者提出了一种同步提取变换与Vold-Kalman结合的燃气轮机振动阶次跟踪方法；以机匣振动速度信号为输入，通过同步提取变换方法获得与振动同步的转速数据；以振动及与其同步的转速数据为输入，使用基于Vold-Kalman滤波的阶次跟踪方法提取振动信号中目标阶次时域波形。

笔者使用该方法对实际案例数据进行了分析，结果表明，该方法能够避免转速不精确、不同步带来的转频偏差问题，能够准确提取目标阶次时域波形。由此可见，该方法适用于在燃气轮机运行过程中跟踪振动阶次分量，可以为振动状态监测及预警提供数据支撑。