原文链接：http://tecdat.cn/?p=22467

本文展示了如何使用扩展卡尔曼滤波器进行故障检测。本文使用扩展的卡尔曼滤波器对一个简单的直流电机的摩擦力进行在线估计。估计的摩擦力的重大变化被检测出来，并表明存在故障。

电机模型

电机被模拟成具有阻尼系数c，_转动惯量_J，由一个扭矩u驱动。电机的角速度w和加速度˙w是测量输出。

为了使用扩展卡尔曼滤波器估计阻尼系数c，为阻尼系数引入一个辅助状态，并将其导数设为零。

因此，模型状态，x = [w;c]，和测量，y，方程式为：

使用近似值˙x=xn+1-xnTs将连续时间方程转换为离散时间，其中Ts为离散采样周期。这就得到了离散时间模型方程。

指定电机参数

J = 10; % 惯性 Ts = 0.01; % 采样时间

指定初始状态。

x = \[... 0; ... % 角速度 1\]; % 摩擦力% 以摩擦为状态的电机的状态更新方程 % 状态更新方程式 x1 = \[... x0(1)+Ts/J*(u-x0(1)*x0(2)); ... x0(2)\];% 以摩擦为状态的电机的测量方程


% 输出。

% y - 电机测量元素\[角速度；角加速度\]。

y = \[... x(1); ... (u-x(1)*x(2))/J\];

电机经历状态（过程）噪声干扰，q，和测量噪声干扰，r。噪声项是相加的。

过程和测量噪声的平均值为零，E[q]=E[r]=0，协方差Q=E[qq']和R=E[rr']。摩擦状态有很高的过程噪声干扰。这反映了我们希望摩擦力在电机正常运行期间会有变化，并且滤波器能跟踪这种变化。加速和速度状态的噪声很低，但速度和加速度测量的噪声相对较大。

指定过程噪声协方差。

\[... 1e-6 0; ... % 角速度

0 1e-2\]; % 摩擦力

指定测量噪声协方差。

\[... 1e-4 0; ... % 速度测量 0 1e-4\]; % 加速度测量

创建一个扩展的卡尔曼过滤器

创建一个扩展的卡尔曼滤波器来估计模型的状态。我们特别关注阻尼状态，因为这个状态值的急剧变化表明存在故障事件。
创建一个扩展卡尔曼滤波器对象，并指定状态转换和测量函数的雅各布系数。

扩展卡尔曼滤波器的输入参数是之前定义的状态转换和测量函数。初始状态值x0、初始状态协方差、过程和测量噪声协方差也是扩展卡尔曼滤波器的输入。在这个例子中，精确的雅各布函数可以从状态转换函数f和测量函数h中得到。

% 输出 Jac - 在x处计算出的状态雅各布系数


% 雅各布系数

Jac = \[

    1-Ts/J\*x(2) -Ts/J\*x(1); ...

    0 1\];
% 电机模型测量方程的雅各布系数
% 输出

Jac - 在 x 处计算的测量雅各布系数

% 雅各布系数 J = \[ ... 1 0; -x(2)/J -x(1)/J\];

Simulation仿真

为了模拟工厂，创建一个环路，在电机中引入一个故障（虚构的电机剧烈变化）。在模拟回路中，使用扩展的卡尔曼滤波器来估计电机状态，并特别跟踪摩擦状态，检测摩擦力何时发生统计意义上的变化。

电机被模拟成一个脉冲序列，反复加速和减速。这种类型的电机操作对于生产线上的采摘机器人来说是典型的。

在模拟电机时，加入与构建扩展卡尔曼滤波器时使用的Q和R噪声协方差值相似的过程和测量噪声。对于摩擦，使用一个小得多的噪声值，因为除了故障发生时，摩擦大多是恒定的。在模拟过程中人为地诱发故障。

Qv = chol(Q); % 过程噪声的标准偏差 Qv(end) = 1e-2; % 较小的摩擦噪声 Rv = chol(R); % 测量噪声的标准偏差

使用状态更新方程对模型进行仿真，并在模型状态中加入过程噪声。仿真十秒钟后，强制改变电机的摩擦力。使用模型测量功能来模拟电机传感器，并在模型输出中加入测量噪声。

for ct = 1:numel(t) % 模型输出更新 y = y+Rv*randn(2,1); % 添加测量噪声 % 模型状态更新 xSig(:,ct) = x0; % 诱发摩擦力的变化 if t(ct) == 10 x1(2) = 10; % 步骤变化 x1n = x1+Qv*randn(nx % 加入过程噪声Significant friction change at 10.450000

为了从电机测量值中估计电机状态，使用扩展卡尔曼滤波器的预测和纠正命令。

% 使用扩展卡尔曼滤波器进行状态估计 x_corr = correct(ekf,y,u(ct),J,Ts); % 根据当前测量结果修正状态估计。 predict(ekf,u(ct),J,Ts); % 根据当前状态和输入预测下一个状态。

为了检测摩擦力的变化，使用4秒的移动窗口来计算估计的摩擦力平均值和标准偏差。在最初的7秒之后，锁定计算的平均值和标准差。这个最初计算出的平均值是摩擦力的预期无故障平均值。7秒后，如果估计的摩擦力与预期的无故障平均值相差超过3个标准差，这就意味着摩擦力有了明显的变化。为了减少噪音和估计摩擦力的影响，在与3个标准差的界限比较时，使用估计摩擦力的平均值。

% 计算估计平均值和标准偏差。 else % 存储计算出的平均数和标准差，不需要 %重新计算。 fMean(ct) = fMean(ct-1) % 使用预期的摩擦力平均值和标准偏差来检测 %摩擦力变化。 estFriction = mean(xSigEst(2, if fChanged(ct) && ~fChanged(ct-1) % 检测摩擦变化信号的上升沿|fChanged|

使用估计的状态来计算估计的输出。计算测量输出和估计输出之间的误差，并计算出误差统计。误差统计可用于检测摩擦力的变化。这一点将在后面详细讨论。

kurtosis(ySigEst(1,idx)-ySig(1,idx)); kurtosis(ySigEst(2,idx)-ySig(2,idx))\];