CN110555398B - 一种基于滤波最优平滑确定故障首达时刻的故障诊断方法 - Google Patents

Abstract

本发明公开了一种基于滤波最优平滑确定故障首达时刻的故障诊断方法。本发明中历史数据首先要通过Kalman滤波迭代获得，迭代步数的上限为总的采样次数，下限为滤波过程中发现状态产生故障或漂移的时刻。对历史数据中的观测值和状态估计值的误差求取平均范数，设定故障或漂移的阈值为3倍误差平均2‑范数。通过对历史数据的最优平滑估计使得估计精度得到提升，从而更加精确的判断故障或漂移产生的首达时刻，并可将该方法应用在船舶航行发生漂移后GPS导航定位系统中。

Description

一种基于滤波最优平滑确定故障首达时刻的故障诊断方法

技术领域

本发明属于故障诊断领域，具体涉及一种基于Kalman滤波最优平滑确定故障首达时刻的故障诊断方法。

背景技术

故障诊断领域中，由于系统某些状态故障和漂移是缓变的，在故障首达时刻附近变化量较少，又由于实际工况下系统的噪声和观测噪声影响，故障或漂移可能被其掩盖，难以被察觉，并且对于系统可观性来讲，状态的重构在理论上又具有延迟性，而从工程应用观点来看，确定系统故障的首达时刻尤其重要，寻求一种确定故障首达时刻的方法具有较大价值。传统的Kalman滤波是建立在模型精确和随机干扰信号统计特性已知的基础上的，然而随机干扰信号过大可能会使传统的Kalman滤波算法失去最优性，估计精度大大降低。而Kalman滤波最优平滑可以提高估计精度，从而更早的发现系统故障的首达时刻。

发明内容

本发明针对现有技术的不足，设计一种基于Kalman滤波最优平滑确定故障首达时刻的故障诊断方法。首先，本发明通过Kalman滤波算法对系统状态进行估计，在估计过程中状态某一分量发生故障或参数漂移，经历一定时刻的迭代后，发生故障的状态分量超过故障报警阈值，开始报警。然后，滤波过程中获得的历史数据开始对系统状态进行平滑估计，当状态分量超过平滑故障报警阈值后获得系统的故障首达时刻，并将此方法应用在船舶GPS导航系统对船舶航行的精确定位中。

本发明包括以下各步骤：

步骤(1)通过系统的状态方程和观测方程使用Kalman滤波算法对系统状态进行估计，并持续迭代：

步骤(1-1)Kalman滤波过程：

给定如下离散线性系统：

其中，k为离散时间，系统在时刻k的状态为X_k∈Rⁿ；Y_k∈Rⁿ为对应观测信号；W_k∈R^r为输入的白噪声；V_k∈R^m为观测噪声；称式(1)为状态方程和观测方程。称A为状态转移矩阵，C为观测矩阵。

递推Kalman滤波器如下：

式(2)为状态一步预测式。

式(3)为新息表达式；式(4)为状态更新式。

K_k+1＝P_k+1|kA^T[AP_k+1|kA^T+R]^-1 (5)

式(5)为滤波增益矩阵求取式。

P_k+1|k＝AP_k|kA^T+Q (6)

式(6)为一步预测协方差阵。

P_k+1|k+1＝[I_n-K_k+1A]P_k+1|k (7)

式(7)为协方差更新式。

其中状态估计值的初始值和初始协方差阵为：

步骤(1-2)某一状态分量发生故障：

系统状态某一分量f开始发生故障，由于故障首达时刻未知，将其设为r。此时系统的状态方程中状态转移部分加入缓变矩阵X_fault：

X_k+1＝AX_k-X_fault+W_k (9)

其中X¹到X^f-1和X^f+1到Xⁿ为正常的状态分量，X^f为发生故障的状态分量。

此时系统的观测方程不变。

步骤(2)滤波过程发现状态发生故障或漂移：

在迭代过程中，状态真实值无法直接获取，只能通过传感器的观测获得状态的观测值，对历史估计值和观测值之间的误差求取2-范数均值，设定报警阈值为3倍2-范数均值。

其中N为迭代当前时刻。若N时刻状态估计值和观测值的误差大于报警阈值，则认为系统状态已经发生了故障或漂移。

步骤(3)对系统状态进行平滑(Smooth)：

平滑过程如下：

J_k＝P_k|kA^T[P_k|k+1]^-1 (12)

式(12)为最优平滑增益矩阵求取式。式(12)的边界条件为：k＝N-1时，P_k+1|N＝P_N|N。

式(13)为最优平滑状态更新式。式(13)的边界条件为：k＝N-1时，

式(14)为最优平滑协方差阵更新式。

步骤(4)确定故障的首达时刻：

在后验分析中，等待获得更多的观测数据，这时对状态的估计可以通过最优平滑进一步提高精度。在最优平滑过程中，判断平滑后当前时刻状态是否超过(11)报警阈值，如果超出则可确定故障首达时间。

本发明的有益效果：在Kalman滤波可以有效地跟踪线性系统的基础上，在滤波过程中初步发现状态故障或漂移并获得历史数据后，用该历史数据反向估计，利用平滑可以进一步提高系统的估计精度，从而更早的发现超出阈值的故障状态。最终将其应用到船舶的匀速直线航行的跟踪中，相较于传统检测手段，本发明有效提高了故障首达时刻的检测精度。

附图说明

图1是本发明的算法实现流程图；

图2是Kalman滤波和最优平滑示意图；

图3是船舶的跟踪轨迹图；

图4是观测值和滤波估计值和滤波+平滑估计值与真实值之间的误差对比；

图5是滤波估计值和滤波+平滑估计值与观测值之间的误差对比。

具体实施方式

以下结合附图对本发明作进一步说明。

如图1所示，本发明提出基于Kalman滤波加平滑提高估计精度的故障诊断方法，包括以下各步骤：

1、通过系统的状态方程和观测方程使用Kalman滤波算法对系统状态进行估计，并持续迭代：

步骤(1-1)Kalman滤波过程：

给定如下离散线性系统：

其中，k为离散时间，系统在k时刻的状态为X_k∈Rⁿ；Y_k∈Rⁿ为对应观测信号；W_k∈R^r为输入的白噪声；V_k∈R^m为观测噪声；称式(1)为状态方程和观测方程。称A为状态转移矩阵，C为观测矩阵。

假定船舶驶出码头沿直线方向航行。以码头的出发处为坐标原点，设采样时间T₀，用l_k表示传播在采样时刻kT₀处的真实位置，用y_k表示在时刻kT₀处GPS定位的观测值，则有观测模型：

y_k＝l_k+v_k (2)

式中，v_k表示GPS定位误差(观测噪声)，可假设它是零均值、方差为

的白噪声，方差

是通过大量GPS观测试验数据用统计方法获取。记在时刻kT₀处船舶速度为

加速度为a_k，由匀加速运动公式有：

其中加速度a_k由机动加速度计u_k和随机加速度w_k两部分合成，即：

a_k＝u_k+w_k (5)

式中，u_k为船舶动力系统的控制信号，它是人为输出的已知机动控制信号；w_k是由海风和洋流等外界扰动因素造成的随机加速度，假设它的均值为零、方差为

的独立于v_k的白噪声。定义在采样时刻kT₀处系统的状态x_k为船舶的位置和速度，即：

可得到传播运动的状态方程：

观测方程为：

在不考虑机动目标自身的动力因素时，匀速直线运动的船舶系统为四维系统。状态包含水平方向的位置和速度和纵向的位置和速度，则该船舶系统方程可以用下式表示。

假定传播在二维水平面上运动，初始位置为(-100m,200m),水平运动速度为2m/s，垂直方向的运动速度为20m/s,GPS接收机的扫描周期为T＝1s，总的采样次数为100次，观测噪声的均值为0，方差为100。

递推Kalman滤波器如下：

式(11)为状态一步预测式。

式(12)为新息表达式；式(13)为状态更新式。

K_k+1＝P_k+1|kA^T[AP_k+1|kA^T+R]^-1 (14)

式(14)为滤波增益矩阵求取式。

P_k+1|k＝AP_k|kA^T+Q (15)

式(15)为一步预测协方差阵。

P_k+1|k+1＝[I_n-K_k+1A]P_k+1|k (16)

式(16)为协方差更新式。

其中状态估计值的初始值和初始协方差阵为：

步骤(1-2)某一状态分量发生故障：

系统状态某一分量f开始发生故障，由于故障首达时刻未知，将其设为r。此时系统的状态方程中状态转移部分加入缓变矩阵X_fault：

X_k+1＝AX_k-X_fault+W_k (18)

其中X¹到X^f-1和X^f+1到Xⁿ为正常的状态分量，X^f为发生故障的状态分量。此时观测方程不变。

对系统只进行Kalman滤波估计的跟踪轨迹在图2和图3中有体现。

2、滤波过程发现状态发生故障或漂移：

在迭代过程中，状态真实值是无法直接获取的，只能通过传感器的观测获得状态的观测值，对历史估计值和观测值之间的误差求取2-范数均值，设定报警阈值为3倍2-范数均值。

其中N为当前迭代时刻。若N时刻状态估计值和观测值的误差大于报警阈值，则认为系统状态已经发生了故障或漂移。

3、对系统状态进行平滑(Smooth)：

平滑过程如下：

J_k＝P_k|kA^T[P_k|k+1]^-1 (21)

式(21)为最优平滑增益矩阵求取式。式(21)的边界条件为：k＝N-1时，P_k+1|N＝P_N|N。

式(22)为最优平滑状态更新式。式(22)的边界条件为：k＝N-1时，

式(23)为最优平滑协方差阵更新式。

图3所表示的内容即为对系统进行Kalman滤波估计并再此基础进行最优平滑估计的跟踪轨迹，可以看出由于观测噪声较大Kalman滤波对状态估计的精度不高，但在此基础上的最优平滑过程可以进一步提高估计精度。

4、确定故障的首达时刻：

在后验分析中，等待获得更多的观测数据，这时对状态的估计可以通过最优平滑进一步提高精度。在最优平滑过程中，判断平滑后当前时刻状态是否超过(20)报警阈值，如果超出则可确定故障首达时间。

图4和图5分别为测值和滤波估计值和滤波+平滑估计值与真实值之间的误差对比。图5为滤波估计值和滤波+平滑估计值与观测值之间的误差对比。可以看出与仅进行Kalman滤波的误差对比，滤波加最优平滑后的误差更小。

表1 仅Kalman滤波和Kalman滤波+最优平滑的误差均值和均方差对比(1)

	误差均值	均方差
			仅Kalman滤波	14.459	29.177
Kalman滤波+最优平滑	11.571	23.111

在该仿真试验中，真实故障首达时刻r＝79,仅滤波发现故障首达时刻r_f＝96,而经过Kalman滤波+最优平滑后发现故障首达时刻r_f+s＝81。结合上表可以进一步证明Kalman滤波+最优平滑可以有效提高估计精度从而更精确的发现状态故障或漂移的首达时刻。

为进一步确定方法的正确性，再进行蒙特卡洛仿真，以下为实验结果：

表2 仅Kalman滤波和Kalman滤波+最优平滑的误差均值和均方差对比(2)

	误差均值	均方差
			仅Kalman滤波	17.616	35.469
Kalman滤波+最优平滑	14.029	29.313

在该仿真试验中，真实故障首达时刻r＝77,仅滤波发现故障首达时刻r_f＝84,而经过Kalman滤波+最优平滑后发现故障首达时刻r_f+s＝78。

表3 仅Kalman滤波和Kalman滤波+最优平滑的误差均值和均方差对比(3)

	误差均值	均方差
			仅Kalman滤波	7.929	10.868
Kalman滤波+最优平滑	8.567	16.115

在该仿真试验中，真实故障首达时刻r＝87,仅滤波发现故障首达时刻r_f＝92,而经过Kalman滤波+最优平滑后发现故障首达时刻r_f+s＝88。

表4 仅Kalman滤波和Kalman滤波+最优平滑的误差均值和均方差对比(4)

	误差均值	均方差
			仅Kalman滤波	6.678	11.383
Kalman滤波+最优平滑	5.912	10.092

在该仿真试验中，真实故障首达时刻r＝51,仅滤波发现故障首达时刻r_f＝68,而经过Kalman滤波+最优平滑后发现故障首达时刻r_f+s＝56。

Claims (1)

1.一种基于滤波最优平滑确定故障首达时刻的故障诊断方法，其特征在于该方法具体包括以下各步骤：

步骤(1)通过系统的状态方程和观测方程使用Kalman滤波算法对系统状态进行估计，并持续迭代：

步骤(1-1)Kalman滤波过程：

给定如下离散线性系统：

其中，k为离散时间，系统在时刻k的状态为X_k∈Rⁿ；Y_k∈Rⁿ为对应观测信号；W_k∈R^r为输入的白噪声；V_k∈R^m为观测噪声；称式(1)为状态方程和观测方程；A为状态转移矩阵，C为观测矩阵；

递推Kalman滤波器如下：

式(2)为状态一步预测式；

式(3)为新息表达式；式(4)为状态更新式；

K_k+1＝P_k+1|kA^T[AP_k+1|kA^T+R]^-1 (5)

式(5)为滤波增益矩阵求取式；

P_k+1|k＝AP_k|kA^T+Q (6)

式(6)为一步预测协方差阵；

P_k+1|k+1＝[I_n-K_k+1A]P_k+1|k (7)

式(7)为协方差更新式；

其中状态估计值的初始值和初始协方差阵为：

步骤(1-2)某一状态分量发生故障：

系统状态某一分量f开始发生故障，由于故障首达时刻未知，将其设为r；此时系统的状态方程中状态转移部分加入缓变矩阵X_fault：

X_k+1＝AX_k-X_fault+W_k (9)

其中X¹到X^f-1和X^f+1到Xⁿ为正常的状态分量，X^f为发生故障的状态分量；此时观测方程不变；

步骤(2)滤波过程发现状态发生故障或漂移：

在迭代过程中，状态真实值无法直接获取，只能通过传感器的观测获得状态的观测值，对历史估计值和观测值之间的误差求取2-范数均值，设定报警阈值为3倍2-范数均值；

其中N为迭代当前时刻；若N时刻状态估计值和观测值的误差大于报警阈值，则认为系统状态已经发生了故障或漂移；

步骤(3)对系统状态进行平滑：

平滑过程如下：

J_k＝P_k|kA^T[P_k|k+1]^-1 (12)

式(12)为最优平滑增益矩阵求取式；式(12)的边界条件为：k＝N-1时，P_k+1|N＝P_N|N；

式(13)为最优平滑状态更新式；式(13)的边界条件为：k＝N-1时，

式(14)为最优平滑协方差阵更新式；

步骤(4)确定故障的首达时刻：

在后验分析中，等待获得更多的观测数据，这时对状态的估计可以通过最优平滑进一步提高精度；在最优平滑过程中，判断平滑后当前时刻状态是否超过(11)报警阈值，如果超出则可确定故障首达时间。

Images