CN102564486A

CN102564486A - Correction method for slow deviation faults of sensor

Info

Publication number: CN102564486A
Application number: CN2011104335424A
Authority: CN
Inventors: 高羽; 周琴
Original assignee: Shanghai Dianji University
Current assignee: Shanghai Dianji University
Priority date: 2011-12-21
Filing date: 2011-12-21
Publication date: 2012-07-11

Abstract

The present invention provides a correction method for sensor slow bias faults, comprising step a: establishing a state equation and an observation equation for the measured value of the sensor; step b: extending the state equation and observation equation in time; step c: using The Kalman filter of the extended state estimates the bias, which is then subtracted from the measured value. The invention can estimate the slowly changing deviation of the sensor in real time, so that when the slowly changing deviation occurs in the sensor, the deviation can be estimated and subtracted from the measured value of the sensor.

Description

A Correction Method for Slow Offset Fault of Sensor

技术领域 technical field

本发明属于传感器故障诊断及校正技术领域，涉及一种可用于校正传感器慢偏故障的方法。The invention belongs to the technical field of sensor fault diagnosis and correction, and relates to a method for correcting slow-bias faults of sensors.

背景技术 Background technique

传感器是一种检测装置，它能检测到被测量的信息并将其按照一定的规律变换成为相应的电信号或其他所需形式的信息输出，以满足信息的传输、处理、存储、显示、记录和控制等要求，是实现自动检测和自动控制的首要环节。随着信息处理技术以及微处理机和计算机技术的高速发展，传感器已广泛应用于航空航天、军事、医疗卫生、工农业生产以及环境保护监测等各个领域。尤其是近几十年来，随着材料科学领域的突破性进展，出现了越来越多新型的传感器，如光纤、超声波和微机械传感器等等。A sensor is a detection device that can detect the measured information and convert it into a corresponding electrical signal or other required form of information output according to certain rules, so as to meet the requirements of information transmission, processing, storage, display, and recording. It is the primary link to realize automatic detection and control. With the rapid development of information processing technology and microprocessor and computer technology, sensors have been widely used in various fields such as aerospace, military, medical and health, industrial and agricultural production, and environmental protection monitoring. Especially in recent decades, with breakthroughs in the field of material science, more and more new sensors have emerged, such as optical fiber, ultrasonic and micromechanical sensors, etc.

作为信息采集系统的前端单元，传感器的准确性和稳定性是整个系统良好运行的先决条件。对于一些复杂的大型动态系统，如飞机，电站等，需要依赖准确可靠的传感器去获得测量信号来判断和监视系统状态，以进行有效控制。可以说，传感器的准确性和可靠性极大程度的决定了整个系统的性能、可靠性和安全性。As the front-end unit of the information collection system, the accuracy and stability of the sensor is a prerequisite for the good operation of the entire system. For some complex large-scale dynamic systems, such as airplanes and power stations, it is necessary to rely on accurate and reliable sensors to obtain measurement signals to judge and monitor the system status for effective control. It can be said that the accuracy and reliability of the sensor greatly determine the performance, reliability and safety of the entire system.

传感器的测量精度等信息通常由传感器生产厂家在产品出厂时给定或用户使用前通过其他方法标定。这两种标定精度通常是在一定条件下获得的，而实际应用中环境条件复杂多变，很难始终保持其与标定时的条件一致，当传感器的工作条件与它的标定条件不一致时，传感器精度往往不能真实反映传感器测量值的实际准确程度。很多传感器的工作环境非常恶劣，例如风力发电机状态监测系统中的传感器就工作在低温、大风或者高盐雾条件下，因此非常容易出现各种各样的故障。如果故障传感器的输出数据一直作为控制系统的输入时，可能会导致整个系统运行性能严重下降，甚至出现灾难性的事故。Information such as the measurement accuracy of the sensor is usually given by the sensor manufacturer when the product leaves the factory or calibrated by other methods before the user uses it. These two calibration accuracies are usually obtained under certain conditions, but the environmental conditions in practical applications are complex and changeable, and it is difficult to keep them consistent with the calibration conditions. When the working conditions of the sensor are inconsistent with its calibration conditions, the sensor Accuracy is often not a true reflection of the actual accuracy of sensor measurements. The working environment of many sensors is very harsh. For example, the sensors in the wind turbine condition monitoring system work under the conditions of low temperature, strong wind or high salt spray, so it is very prone to various failures. If the output data of the faulty sensor is always used as the input of the control system, it may lead to a serious decline in the performance of the entire system, or even a catastrophic accident.

为了保证测量的精度和可靠性，传统的方法主要是通过一个较高精度的传感器测量值作为基准，用来校正其余传感器测量值。这种方法的缺点在于，首先精度较高的传感器通常价格昂贵，使得成本大大增加；其次，即便在不需顾及成本的应用场合中，复杂的环境因素往往也会导致高精度的传感器精度下降，这样其基准作用也就无从谈起。另外一种实际中极为常用的方法是对系统中的传感器定期进行校准和维护，在定期校准带来大量的人力物力损耗的同时，人们还是无从得知两次校准之间传感器的输出是否准确，而校准周期的确定也是另外一个难题。因此怎样能够准确有效的在传感器出现测量偏差时，对这个偏差进行校正，具有十分重要的实用价值。In order to ensure the accuracy and reliability of the measurement, the traditional method mainly uses a higher-precision sensor measurement value as a reference to correct the other sensor measurement values. The disadvantage of this method is that firstly, sensors with higher precision are usually expensive, which greatly increases the cost; secondly, even in applications where cost is not required, complex environmental factors often lead to a decrease in the accuracy of high-precision sensors. In this way, its benchmark role is out of the question. Another method that is very commonly used in practice is to regularly calibrate and maintain the sensors in the system. While regular calibration brings a lot of manpower and material resources, people still have no way of knowing whether the output of the sensor between two calibrations is accurate. The determination of the calibration cycle is another problem. Therefore, how to accurately and effectively correct the deviation when the sensor has a measurement deviation has very important practical value.

发明内容 Contents of the invention

本发明的目的在于克服现有技术的不足，提供一种传感器慢偏故障的校正方法，该方法能够实时估计出传感器的慢变偏差，从而能在传感器出现缓慢变化的偏差时，将该偏差估计出来，从传感器测量值中减去。这样就实现了测量偏差和测量真值的分离。The purpose of the present invention is to overcome the deficiencies of the prior art, and provide a correction method for sensor slow-bias faults, which can estimate the slowly-varying deviation of the sensor in real time, so that the deviation can be estimated out and subtracted from the sensor measurement. In this way, the separation of measurement deviation and measurement true value is realized.

为达到上述目的，本发明采用如下技术方案：To achieve the above object, the present invention adopts the following technical solutions:

一种传感器慢偏故障的校正方法，包含如下步骤：A method for correcting a sensor slow bias fault, comprising the following steps:

步骤a：对传感器的测量值建立状态方程和观测方程；Step a: establish a state equation and an observation equation for the measured value of the sensor;

步骤b：对所述状态方程和观测方程进行时间上的扩展；Step b: expanding the state equation and observation equation in time;

步骤c：利用扩展状态的卡尔曼滤波估计偏差，再将偏差值从测量值中减去。Step c: use the Kalman filter of the extended state to estimate the deviation, and then subtract the deviation value from the measured value.

进一步地，所述步骤a可以按如下步骤实施：Further, the step a can be implemented as follows:

步骤a1：将含有传感器测量偏差的线性动态系统的状态空间模型描述为：Step a1: Describe the state-space model of a linear dynamical system with sensor measurement bias as:

x_t+1＝Fx_t+w_t (1)x _t+1 ＝Fx _t +w _t (1)

y_t＝Φx_t+b_t+v_t (2)y _t ＝Φx _t +b _t +v _t (2)

b_t+1＝Ab_t+η_t (3)b _t+1 =Ab _t +η _t (3)

其中，下标t表示时间，x_t+1为N×1维的状态向量，F为状态转移矩阵，y_t为p×1维的观测向量，Φ为测量矩阵，b_t为M×1维的传感器测量偏差向量，A为b_t的状态转移矩阵，w_t为状态噪声，v_t为测量噪声，η_t为偏差状态噪声，且状态噪声w_t，测量噪声v_t及偏差状态噪声η_t均为相互独立的零均值的高斯白噪声，且满足：Among them, the subscript t represents time, x _t+1 is the N×1-dimensional state vector, F is the state transition matrix, y _t is the p×1-dimensional observation vector, Φ is the measurement matrix, and b _t is the M×1-dimensional The sensor measurement deviation vector, A is the state transition matrix of b _t , w _t is the state noise, v _t is the measurement noise, η _t is the deviation state noise, and the state noise w _t , the measurement noise v _t and the deviation state noise η _t are independent Gaussian white noise with zero mean, and satisfy:

E(w_t)＝0E(w _t )＝0

$E E. (({w w}_{j j} {w w}_{t t}^{T T})) = = R R {δ δ}_{jt jt} - - - - - - ((44))$

E(v_t)＝0E(v _t )＝0

$E E. (({v v}_{j j} {v v}_{t t}^{T T})) = = Q Q {δ δ}_{jt jt} - - - - - - ((55))$

E(η_t)＝0E(η _t )=0

$E E. (({η η}_{j j} {η η}_{t t}^{T T})) = = G G {δ δ}_{jt jt} - - - - - - ((66))$

式中R、Q和G分别是状态、测量和偏差状态噪声的方差；where R, Q, and G are the variances of the state, measurement, and bias state noise, respectively;

步骤a2，将所述偏差向量扩展为状态向量，获得状态方程和观测方程：Step a2, expanding the deviation vector into a state vector to obtain a state equation and an observation equation:

所述状态方程为：The equation of state is:

${\overset{~ ~}{x x}}_{t t + + 11} = = \overset{~ ~}{F f} {\overset{~ ~}{x x}}_{t t} + + {\overset{~ ~}{w w}}_{t t}$

所述观测方程为：The observation equation is:

${y the y}_{t t} = = \overset{~ ~}{Φ Φ} {\overset{~ ~}{x x}}_{t t} + + {v v}_{t t}$

其中， ${\tilde{x}}_{t + 1} = [\begin{matrix} x_{t + 1} \\ b_{t + 1} \end{matrix}],$ $\tilde{F} = [\begin{matrix} F & 0 \\ 0 & A \end{matrix}],$ ${\tilde{w}}_{t} = [\begin{matrix} w_{t} \\ η_{t} \end{matrix}],$ $\tilde{Φ} = [\begin{matrix} Φ & I \end{matrix}],$ I为单位矩阵。in, ${\tilde{x}}_{t + 1} = [\begin{matrix} x_{t + 1} \\ b_{t + 1} \end{matrix}],$ $\tilde{f} = [\begin{matrix} f & 0 \\ 0 & A \end{matrix}],$ ${\tilde{w}}_{t} = [\begin{matrix} w_{t} \\ η_{t} \end{matrix}],$ $\tilde{Φ} = [\begin{matrix} Φ & I \end{matrix}],$ I is the identity matrix.

进一步地，所述步骤b可以包含：Further, the step b may include:

步骤b1，将状态变量进行时间上的扩展，即Step b1, expand the state variable in time, namely

${X x}_{t t + + 11} = = [[{x x}_{t t + + 11}^{11},, {x x}_{t t}^{11},, . . . . . .,, {x x}_{t t - - K K + + 22}^{11},, {x x}_{t t + + 11}^{22},, {x x}_{t t}^{22},, . . . . . .,, {x x}_{t t - - K K + + 22}^{22},, . . . . . . {x x}_{t t + + 11}^{N N},, {x x}_{t t}^{N N},, . . . . . .,, {x x}_{t t - - K K + + 22}^{N N}]]$

其中上角标为维数，where the upper corner is labeled dimension,

步骤b2，将传感器测量值进行时间上的扩展，即Step b2, expand the measured value of the sensor in time, namely

${Y Y}_{t t + + 11} = = [[{y the y}_{t t + + 11}^{11},, {y the y}_{t t}^{11},, . . . . . .,, {y the y}_{t t - - K K + + 22}^{11},, {y the y}_{t t + + 11}^{22},, {y the y}_{}^{t t},, . . . . . .,, {y the y}_{t t - - K K + + 22}^{22},, . . . . . . {y the y}_{t t + + 11}^{p p},, {y the y}_{t t}^{p p},, . . . . . .,, {y the y}_{t t - - K K + + 22}^{p p}]]$

其中上角标为维数，where the upper corner is labeled dimension,

步骤b3，使状态转移矩阵为：Step b3, make the state transition matrix as:

${\overset{&OverBar; &OverBar;}{F f}}_{ppf ppf} = = [\begin{matrix} {F f}_{ppf ppf}^{((11))} & 00 & . . . . . . & 00 \\ 00 & {F f}_{ppf ppf}^{((22))} & . . . . . . & 00 \\ . . . . . . & . . . . . . & . . . . . . & . . . . . . \\ 00 & . . . . . . & 00 & {F f}_{ppf ppf}^{((N N))} \end{matrix}]$

${F f}_{ppf ppf} = = [\begin{matrix} h h ((00)) & h h ((11)) & . . . . . . & h h ((K K - - 11)) \\ 11 & 00 & . . . . . . & 00 \\ . . . . . . & . . . . . . & . . . . . . & . . . . . . \\ 00 & . . . . . . & 11 & 00 \end{matrix}]$

步骤b4，使测量矩阵为：Step b4, make the measurement matrix as:

${Φ Φ}_{ppf ppf} = = [\begin{matrix} φ φ ((1,1 1,1)) & 00_{11 \times \times ((K K - - 11))} & φ φ ((1,2 1,2)) & 00_{11 \times \times ((K K - - 11))} & . . . . . . & φ φ ((11,, N N)) & 00_{11 \times \times ((K K - - 11))} \\ 00 & φ φ ((1,1 1,1)) & 00_{11 \times \times ((K K - - 11))} & φ φ ((1,2 1,2)) & . . . . . . & φ φ ((11,, N N)) & 00_{11 \times \times ((K K - - 22))} \\ . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . \\ 00 & 00 & . . . . . . & φ φ ((1,1 1,1)) & 00 & . . . . . . & φ φ ((11,, N N)) \\ . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . \\ φ φ ((p p,, 11)) & 00_{11 \times \times ((K K - - 11))} & φ φ ((p p,, 22)) & 00_{11 \times \times ((K K - - 11))} & . . . . . . & φ φ ((p p,, N N)) & 00_{11 \times \times ((K K - - 11))} \\ . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . \\ 00 & 00 & . . . . . . & φ φ ((p p,, 11)) & 00 & . . . . . . & φ φ ((p p,, N N)) \end{matrix}]$

步骤b5，扩展后的状态方程和观测方程分别为：In step b5, the expanded state equation and observation equation are respectively:

${X x}_{t t + + 11} = = {\overset{&OverBar; &OverBar;}{F f}}_{ppf ppf} {X x}_{t t}$

Y_t＝Φ_ppfX_t+V_t Y _t ＝Φ _ppf X _t +V _t

其中，V_t为扩展后的测量噪声，且有V_t＝[v_t，…，v_t-K+1]^T。Wherein, V _t is the expanded measurement noise, and V _t =[v _t , . . . , v _t-K+1 ] ^T .

进一步地，所述步骤c可以包含：Further, the step c may include:

步骤c1，假设传感器的测量偏差与可测物理量c_t之间存在线性关系，即：Step c1, assuming that there is a linear relationship between the measurement deviation of the sensor and the measurable physical quantity c _t , namely:

b_t≈α_tc_t b _t ≈α _t c _t

状态方程和观测方程分别为：The state equation and observation equation are respectively:

${\overset{~ ~}{X x}}_{t t + + 11} = = {\overset{~ ~}{F f}}_{c c} {\overset{~ ~}{X x}}_{t t}$

${Y Y}_{t t} = = {\overset{~ ~}{Φ Φ}}_{c c} {\overset{~ ~}{X x}}_{t t} + + {V V}_{t t}$

其中in

${\overset{~ ~}{X x}}_{t t + + 11} = = {[[{x x}_{t t + + 11},, {x x}_{t t},, . . . . . .,, {x x}_{t t - - K K + + 22},, {α α}_{t t + + 11},, {α α}_{t t},, . . . . . .,, {α α}_{t t - - K K + + 22}]]}^{T T}$

Y_t＝[y_t，…，y_t-K+1]^T Y _t =[y _t ,...,y _t-K+1 ] ^T

${\overset{~ ~}{F f}}_{c c} = = [\begin{matrix} {F f}_{ppf ppf} & 00 \\ 00 & {F f}_{ppf ppf} \end{matrix}]$

${\overset{~ ~}{Φ Φ}}_{c c} = = [\begin{matrix} Φ Φ & 00 & . . . . . . & 00 & {c c}_{t t} & 00 & . . . . . . & 00 \\ 00 & Φ Φ & . . . . . . & 00 & 00 & {c c}_{t t - - 11} & . . . . . . & 00 \\ 00 & . . . . . . & . . . . . . & . . . . . . & 00 & . . . . . . & . . . . . . & . . . . . . \\ 00 & . . . . . . & 00 & Φ Φ & 00 & . . . . . . & 00 & {c c}_{t t - - K K + + 11} \end{matrix}]$

c_t≠c_t-1≠…≠c_t-K+1.c _t ≠c _t-1 ≠…≠c _t-K+1 .

步骤c2，假设t＝0时的状态

和方差P₀已知，

{\hat{\tilde{X}}}_{t + 1} = {[{\hat{x}}_{t + 1}, {\hat{x}}_{t}, . . ., {\hat{x}}_{t - K + 2}, {\hat{α}}_{t + 1}, {\hat{α}}_{t}, . . ., {\hat{α}}_{t - K + 2}]}^{T},

对t＝1，2，…(t为自然数)，Step c2, assume the state at t=0

and variance _P0 known,

{\hat{\tilde{x}}}_{t + 1} = {[{\hat{x}}_{t + 1}, {\hat{x}}_{t}, . . ., {\hat{x}}_{t - K + 2}, {\hat{α}}_{t + 1}, {\hat{α}}_{t}, . . ., {\hat{α}}_{t - K + 2}]}^{T},

To t=1, 2, ... (t is a natural number),

对一阶信号，令For a first-order signal, let

$h h ((k k)) = = \frac{44 K K - - 66 k k - - 44}{K K ((K K - - 11))}$

对二阶信号，令For second-order signals, let

$h h ((k k)) = = \frac{{99 K K}^{22} + + ((- - 2727 - - 3636 k k)) K K + + {3030 k k}^{22} + + 4242 k k + + 1818}{{K K}^{33} - - 33 {K K}^{22} + + 22 K K}$

步骤c21，预测：Step c21, predict:

${\overset{^^}{\overset{~ ~}{X x}}}_{t t + + 11}^{- -} = = {\overset{~ ~}{F f}}_{c c} {\overset{^^}{\overset{~ ~}{X x}}}_{t t}$

${P P}_{t t + + 11}^{- -} = = {\overset{~ ~}{F f}}_{c c} {P P}_{t t} {\overset{~ ~}{F f}}_{c c}^{T T}$

为t+1时刻扩展状态预测值，

为扩展的状态转移矩阵，

为t时刻的扩展状态估计值，

is the extended state prediction value at time t+1,

is the extended state transition matrix,

is the extended state estimate at time t,

为t+1时刻状态估计方差预测值，

is the state estimation variance prediction value at time t+1,

P_t为t时刻状态估计方差，P _t is the state estimation variance at time t,

步骤c22，更新Step c22, update

${K K}_{t t + + 11} = = {P P}_{t t + + 11}^{- -} {\overset{~ ~}{H h}}_{c c}^{T T} (({\overset{~ ~}{H h}}_{c c} {P P}_{t t + + 11}^{- -} {\overset{~ ~}{H h}}_{c c}^{T T} + + {R R}_{t t + + 11}))$

K_t+1为t+1时刻新息加权值，

为扩展后的测量矩阵，R_t+1为测量方差阵，K _t+1 is the weighted value of new information at time t+1,

is the expanded measurement matrix, R _t+1 is the measurement variance matrix,

${\overset{^^}{\overset{~ ~}{X x}}}_{t t + + 11} = = {\overset{^^}{\overset{~ ~}{X x}}}_{t t + + 11}^{- -} + + {K K}_{t t + + 11} (({Y Y}_{t t + + 11} - - {\overset{~ ~}{H h}}_{c c} {\overset{^^}{\overset{~ ~}{X x}}}_{t t + + 11}^{- -}))$

${P P}_{t t + + 11} = = ((I I - - {K K}_{t t + + 11} {\overset{~ ~}{H h}}_{c c})) {P P}_{t t + + 11}^{- -}$

Y_t+1为扩展后的测量值，P_t+1为t+1时刻的状态估计方差阵；Y _t+1 is the expanded measured value, P _t+1 is the variance matrix of state estimation at time t+1;

步骤c23，计算出测量偏差Step c23, calculate the measurement deviation

${\overset{^^}{b b}}_{t t + + 11} = = {\overset{^^}{α α}}_{t t + + 11} {c c}_{t t + + 11}$

步骤c3，将偏差值从测量值中减去。Step c3, subtracting the deviation value from the measured value.

本发明的有益效果是：The beneficial effects of the present invention are:

通过对状态方程和测量方程进行扩展，并利用扩展状态的卡尔曼滤波估计偏差，能在进行系统状态估计的同时估计出测量偏差，进而可以从测量值中将其除去，因此估计出的系统状态与真实值非常接近，完全消除了传感器测量偏差对状态估计结果的影响。By extending the state equation and the measurement equation, and using the extended state Kalman filter to estimate the deviation, the measurement deviation can be estimated while the system state is estimated, and then it can be removed from the measured value, so the estimated system state It is very close to the real value, and completely eliminates the influence of the sensor measurement deviation on the state estimation result.

附图说明 Description of drawings

附图1所示为仿真中传感器的工作温度变化图。Accompanying drawing 1 shows the operating temperature change diagram of the sensor in the simulation.

附图2所示为系统的真实状态与传感器的测量值。Figure 2 shows the real state of the system and the measured values of the sensors.

附图3所示为根据现有技术和根据本发明方法进行状态估计时，所得到的状态估计结果的区别。Accompanying drawing 3 shows the difference between the state estimation results obtained when the state estimation is performed according to the prior art and according to the method of the present invention.

附图4所示为根据现有技术和根据本发明方法估计测量偏差时，所得到的测量偏差的区别。Figure 4 shows the difference in the measurement deviation obtained when the measurement deviation is estimated according to the prior art and according to the method of the present invention.

附图5所示为根据本发明方法进行测量偏差估计时，得到的测量偏差估计结果。Accompanying drawing 5 shows the measurement deviation estimation result obtained when the measurement deviation estimation is performed according to the method of the present invention.

附图6所示为根据本发明方法进行测量偏差估计时，得到的测量偏差估计标准差。Figure 6 shows the estimated standard deviation of the measurement deviation obtained when the measurement deviation is estimated according to the method of the present invention.

具体实施方式 Detailed ways

下面结合实例详细说明本发明的具体实施。The specific implementation of the present invention will be described in detail below in conjunction with examples.

在一种实施方式，一个含有传感器测量偏差的线性动态系统的状态空间模型可以描述为：In one embodiment, a state-space model of a linear dynamical system with sensor measurement bias can be described as:

x_t+1＝Fx_t+w_t (1)x _t+1 ＝Fx _t +w _t (1)

y_t＝Φx_t+b_t+v_t (2)y _t ＝Φx _t +b _t +v _t (2)

b_t+1＝Ab_t+η_t (3)b _t+1 =Ab _t +η _t (3)

E(w_t)＝0E(w _t )＝0

E(v_t)＝0E(v _t )＝0

E(η_t)＝0E(η _t )=0

式中R、Q和G分别是状态、测量和偏差状态噪声的方差。where R, Q, and G are the variances of state, measurement, and bias state noise, respectively.

将偏差向量扩展为状态向量时，式(1)-(3)可写为：When the deviation vector is extended to a state vector, equations (1)-(3) can be written as:

${\overset{~ ~}{x x}}_{t t + + 11} = = \overset{~ ~}{F f} {\overset{~ ~}{x x}}_{t t} + + {\overset{~ ~}{w w}}_{t t} - - - - - - ((77))$

${y the y}_{t t} = = \overset{~ ~}{Φ Φ} {\overset{~ ~}{x x}}_{t t} + + {v v}_{t t} - - - - - - ((88))$

为了通过扩展状态的卡尔曼滤波估计出测量偏差，需要利用一种全新的建模方法对状态方程建模，即建立状态方程的多项式预测模型。In order to estimate the measurement deviation through the Kalman filter of the extended state, it is necessary to use a new modeling method to model the state equation, that is, to establish a polynomial prediction model of the state equation.

下面简述状态方程的一种多项式预测建模方法。A polynomial predictive modeling method for the state equation is briefly described below.

1、对状态变量扩展，1. To expand the state variable,

将状态变量进行时间上的扩展，即Expand the state variables in time, namely

${X x}_{t t + + 11} = = [[{x x}_{t t + + 11}^{11},, {x x}_{t t}^{11},, . . . . . .,, {x x}_{t t - - K K + + 22}^{11},, {x x}_{t t + + 11}^{22},, {x x}_{}^{t t},, . . . . . .,, {x x}_{t t - - K K + + 22}^{22},, . . . . . . {x x}_{t t + + 11}^{N N},, {x x}_{t t}^{N N},, . . . . . .,, {x x}_{t t - - K K + + 22}^{N N}]] - - - - - - ((99))$

其中上角标为维数。where the upper corner is labeled dimension.

2、对传感器测量值扩展，2. To expand the measured value of the sensor,

将传感器测量值进行时间上的扩展，即Extend the sensor measurement value in time, that is,

${Y Y}_{t t + + 11} = = [[{y the y}_{t t + + 11}^{11},, {y the y}_{t t}^{11},, . . . . . .,, {y the y}_{t t - - K K + + 22}^{11},, {y the y}_{t t + + 11}^{22},, {y the y}_{}^{t t},, . . . . . .,, {y the y}_{t t - - K K + + 22}^{22},, . . . . . . {y the y}_{t t + + 11}^{p p},, {y the y}_{t t}^{p p},, . . . . . .,, {y the y}_{t t - - K K + + 22}^{p p}]] - - - - - - ((1010))$

其中上角标为维数。where the upper corner is labeled dimension.

3、使状态转移矩阵为，3. Let the state transition matrix be,

${\overset{&OverBar; &OverBar;}{F f}}_{ppf ppf} = = [\begin{matrix} {F f}_{ppf ppf}^{((11))} & 00 & . . . . . . & 00 \\ 00 & {F f}_{ppf ppf}^{((22))} & . . . . . . & 00 \\ . . . . . . & . . . . . . & . . . . . . & . . . . . . \\ 00 & . . . . . . & 00 & {F f}_{ppf ppf}^{((N N))} \end{matrix}] - - - - - - ((1111))$

${F f}_{ppf ppf} = = [\begin{matrix} h h ((00)) & h h ((11)) & . . . . . . & h h ((K K - - 11)) \\ 11 & 00 & . . . . . . & 00 \\ . . . . . . & . . . . . . & . . . . . . & . . . . . . \\ 00 & . . . . . . & 11 & 00 \end{matrix}] - - - - - - ((1212))$

4、使测量矩阵为，4. Let the measurement matrix be,

${Φ Φ}_{ppf ppf} = = [\begin{matrix} φ φ ((1,1 1,1)) & 00_{11 \times \times ((K K - - 11))} & φ φ ((1,2 1,2)) & 00_{11 \times \times ((K K - - 11))} & . . . . . . & φ φ ((11,, N N)) & 00_{11 \times \times ((K K - - 11))} \\ 00 & φ φ ((1,1 1,1)) & 00_{11 \times \times ((K K - - 11))} & φ φ ((1,2 1,2)) & . . . . . . & φ φ ((11,, N N)) & 00_{11 \times \times ((K K - - 22))} \\ . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . \\ 00 & 00 & . . . . . . & φ φ ((1,1 1,1)) & 00 & . . . . . . & φ φ ((11,, N N)) \\ . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . \\ φ φ ((p p,, 11)) & 00_{11 \times \times ((K K - - 11))} & φ φ ((p p,, 22)) & 00_{11 \times \times ((K K - - 11))} & . . . . . . & φ φ ((p p,, N N)) & 00_{11 \times \times ((K K - - 11))} \\ . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . \\ 00 & 00 & . . . . . . & φ φ ((p p,, 11)) & 00 & . . . . . . & φ φ ((p p,, N N)) \end{matrix}] - - - - - - ((1313))$

系统的状态方程和观测方程可写为：The state equation and observation equation of the system can be written as:

${X x}_{t t + + 11} = = {\overset{&OverBar; &OverBar;}{F f}}_{ppf ppf} {X x}_{t t} - - - - - - ((1414))$

Y_t＝Φ_ppfX_t+V_t (15)Y _t ＝Φ _ppf X _t +V _t (15)

在进行传感器测量偏差估计时，假设与传感器的测量偏差与一可测物理量c_t之间存在一近似线性关系，例如传感器的测量偏差与其工作温度存在这样的关系，即：When estimating the measurement deviation of the sensor, it is assumed that there is an approximately linear relationship between the measurement deviation of the sensor and a measurable physical quantity _ct , for example, there is such a relationship between the measurement deviation of the sensor and its operating temperature, namely:

b_t≈α_tc_t (16)b _t ≈α _t c _t (16)

则式(7)和(8)可写为：Then formulas (7) and (8) can be written as:

${\overset{~ ~}{X x}}_{t t + + 11} = = {\overset{~ ~}{F f}}_{c c} {\overset{~ ~}{X x}}_{t t} - - - - - - ((1717))$

${Y Y}_{t t} = = {\overset{~ ~}{Φ Φ}}_{c c} {\overset{~ ~}{X x}}_{t t} + + {V V}_{t t} - - - - - - ((1818))$

其中in

${\overset{~ ~}{X x}}_{t t + + 11} = = {[[{x x}_{t t + + 11},, {x x}_{t t},, . . . . . .,, {x x}_{t t - - K K + + 22},, {α α}_{t t + + 11},, {α α}_{t t},, . . . . . .,, {α α}_{t t - - K K + + 22}]]}^{T T} - - - - - - ((1010))$

Y_t＝[y_t，…，y_t-K+1]^T (20) _Yt = [ _yt , ..., yt _-K+1 ] ^T (20)

${\overset{~ ~}{F f}}_{c c} = = [\begin{matrix} {F f}_{ppf ppf} & 00 \\ 00 & {F f}_{ppf ppf} \end{matrix}] - - - - - - ((21 twenty one))$

${\overset{~ ~}{Φ Φ}}_{c c} = = [\begin{matrix} Φ Φ & 00 & . . . . . . & 00 & {c c}_{t t} & 00 & . . . . . . & 00 \\ 00 & Φ Φ & . . . . . . & 00 & 00 & {c c}_{t t - - 11} & . . . . . . & 00 \\ 00 & . . . . . . & . . . . . . & . . . . . . & 00 & . . . . . . & . . . . . . & . . . . . . \\ 00 & . . . . . . & 00 & Φ Φ & 00 & . . . . . . & 00 & {c c}_{t t - - K K + + 11} \end{matrix}] - - - - - - ((22 twenty two))$

当c_t≠c_t-1≠…≠c_t-K+1时，易知

和

满足控制系统的可观测性条件，即可以用扩展的卡尔曼滤波估计偏差。When c _t ≠c _t-1 ≠…≠c _t-K+1 , it is easy to know

and

The observability condition of the control system is met, that is, the bias can be estimated with the extended Kalman filter.

综上，基于多项式预测模型的传感器偏差估计方法具体步骤如下：In summary, the specific steps of the sensor bias estimation method based on the polynomial prediction model are as follows:

假设t＝0时的状态

和方差P₀已知，

{\hat{\tilde{X}}}_{t + 1} = {[{\hat{x}}_{t + 1}, {\hat{x}}_{t}, . . ., {\hat{x}}_{t - K + 2}, {\hat{α}}_{t + 1}, {\hat{α}}_{t}, . . ., {\hat{α}}_{t - K + 2}]}^{T},

对t＝1，2，…Suppose the state at t=0

and variance _P0 known,

{\hat{\tilde{x}}}_{t + 1} = {[{\hat{x}}_{t + 1}, {\hat{x}}_{t}, . . ., {\hat{x}}_{t - K + 2}, {\hat{α}}_{t + 1}, {\hat{α}}_{t}, . . ., {\hat{α}}_{t - K + 2}]}^{T},

For t = 1, 2, ...

通常对一阶信号，Usually for first-order signals,

$h h ((k k)) = = \frac{44 K K - - 66 k k - - 44}{K K ((K K - - 11))} - - - - - - ((23 twenty three))$

对二阶信号，For second-order signals,

$h h ((k k)) = = \frac{{99 K K}^{22} + + ((- - 2727 - - 3636 k k)) K K + + {3030 k k}^{22} + + 4242 k k + + 1818}{{K K}^{33} - - 33 {K K}^{22} + + 22 K K} - - - - - - ((24 twenty four))$

1、预测：1. Forecast:

${\overset{^^}{\overset{~ ~}{X x}}}_{t t + + 11}^{- -} = = {\overset{~ ~}{F f}}_{c c} {\overset{^^}{\overset{~ ~}{X x}}}_{t t} - - - - - - ((2525))$

${P P}_{t t + + 11}^{- -} = = {\overset{~ ~}{F f}}_{c c} {P P}_{t t} {\overset{~ ~}{F f}}_{c c}^{T T} - - - - - - ((2626))$

为t+1时刻扩展状态预测值，

为扩展的状态转移矩阵，

为t时刻的扩展状态估计值， is the extended state prediction value at time t+1,

is the extended state transition matrix,

is the extended state estimate at time t,

为t+1时刻状态估计方差预测值，

is the state estimation variance prediction value at time t+1,

2、更新：2. Update:

${K K}_{t t + + 11} = = {P P}_{t t + + 11}^{- -} {\overset{~ ~}{H h}}_{c c}^{T T} (({\overset{~ ~}{H h}}_{c c} {P P}_{t t + + 11}^{- -} {\overset{~ ~}{H h}}_{c c}^{T T} + + {R R}_{t t + + 11})) - - - - - - ((2727))$

K_t+1为t+1时刻新息加权值，

is the expanded measurement matrix, R _t+1 is the measurement variance matrix,

${\overset{^^}{\overset{~ ~}{X x}}}_{t t + + 11} = = {\overset{^^}{\overset{~ ~}{X x}}}_{t t + + 11}^{- -} + + {K K}_{t t + + 11} (({Y Y}_{t t + + 11} - - {\overset{~ ~}{H h}}_{c c} {\overset{^^}{\overset{~ ~}{X x}}}_{t t + + 11}^{- -})) - - - - - - ((2828))$

${P P}_{t t + + 11} = = ((I I - - {K K}_{t t + + 11} {\overset{~ ~}{H h}}_{c c})) {P P}_{t t + + 11}^{- -} - - - - - - ((2929))$

3、计算出测量偏差：3. Calculate the measurement deviation:

${\overset{^^}{b b}}_{t t + + 11} = = {\overset{^^}{α α}}_{t t + + 11} {c c}_{t t + + 11} - - - - - - ((3030))$

下面给出本发明的一个具体例子：Provide a concrete example of the present invention below:

设一个传感器测量系统模型如式(31)(32)所描述：Suppose a sensor measurement system model is described by formula (31) (32):

x_t＝0.1t (31) _xt = 0.1t (31)

y_t＝x_t+b_t+v_t (32)y _t =x _t +b _t +v _t (32)

其中，x_t表示系统状态，y_t表示传感器测量输出，b_t表示传感器测量偏差，假设传感器测量偏差与工作温度的关系由式(33)表示，TE_t代表t时刻可测量的传感器工作温度。Among them, x _t represents the system state, y _t represents the sensor measurement output, b _t represents the sensor measurement deviation, assuming that the relationship between sensor measurement deviation and operating temperature is expressed by Equation (33), TE _t represents the measurable sensor operating temperature at time t.

b_t＝0.001t·TE_t (33)b _t =0.001t·TE _t (33)

上述h(k)的计算公式是从一个复杂的方程推导出来的，该方程具有参数L和K，不同的L和K的取值，对应不同的表达形式，在公式23中，N＝1，L＝1；在公式24中，N＝1，L＝2。具体实施时，若要求出在不同L和K的取值下的h(k)的值，可参阅文献“P.Heinonen，Y.Neuvo.FIR-Median HybridFilters with Predictive FIR Substructures.IEEE Trans.on Acoustics，Speech andSignal Processing.1998，36(6)：892-899”The above calculation formula of h(k) is derived from a complex equation, which has parameters L and K, and different values of L and K correspond to different expressions. In formula 23, N=1, L=1; in Equation 24, N=1, L=2. During specific implementation, if you want to obtain the value of h(k) under different values of L and K, please refer to the literature "P.Heinonen, Y.Neuvo.FIR-Median HybridFilters with Predictive FIR Substructures.IEEE Trans.on Acoustics , Speech and Signal Processing.1998, 36(6): 892-899"

如果直接采用扩展状态的卡尔曼滤波算法，则测量偏差无法满足可观测性条件，因此无法获得准确的估计结果。而采用本发明提出的算法，通过控制传感器工作的温度并测量出其随时间变化的温度值，那么通过本发明的模型和方法就能够同时准确估计出系统状态和传感器的测量偏差，计算机仿真结果如图1-图6所示。If the extended state Kalman filter algorithm is used directly, the measurement deviation cannot satisfy the observability condition, so accurate estimation results cannot be obtained. And adopt the algorithm that the present invention proposes, by controlling the working temperature of the sensor and measure its temperature value that changes with time, then by the model and method of the present invention, the measurement deviation of the system state and the sensor can be accurately estimated at the same time, computer simulation results As shown in Figure 1-Figure 6.

仿真中假设人为的使传感器的工作温度按照图1所示的规律变化，图2所示为系统的真实状态与传感器测量值，由于仿真中的测量矩阵为单位阵，传感器的测量值与真实状态在均值上是相等的，但从图2中可以明显看出，由于与温度有关的测量偏差的存在，传感器的测量输出随时间逐渐偏离真实值。In the simulation, it is assumed that the operating temperature of the sensor is artificially changed according to the law shown in Figure 1. Figure 2 shows the real state of the system and the measured value of the sensor. Since the measurement matrix in the simulation is a unit matrix, the measured value of the sensor and the real state They are equal in mean value, but it is obvious from Figure 2 that the measured output of the sensor gradually deviates from the true value over time due to the existence of temperature-related measurement bias.

图3和图4所示为根据现有技术未校准测量偏差进行状态估计与根据本发明方法法进行状态估计的结果比较。从图3中可以看出，采用未经处理的测量数据直接进行状态估计，得到的估计结果由于受到测量偏差的影响，远远偏离了真实值；而采用本发明方法，如图4所示，由于可以同时估计出测量偏差，进而可以从测量值中将其除去，因此估计出的系统状态与真实值非常接近，完全消除了传感器测量偏差对状态估计结果的影响。Fig. 3 and Fig. 4 show the results of state estimation according to the uncalibrated measurement deviation of the prior art compared with the state estimation according to the method of the present invention. It can be seen from Fig. 3 that if the unprocessed measurement data is directly used for state estimation, the obtained estimation result is far away from the true value due to the influence of measurement deviation; while using the method of the present invention, as shown in Fig. 4, Since the measurement deviation can be estimated at the same time, and then can be removed from the measured value, the estimated system state is very close to the real value, and the influence of the sensor measurement deviation on the state estimation result is completely eliminated.

图5和图6所示为传感器测量偏差的估计结果。从图中可以看出，采用本发明算法，在进行系统状态估计的同时可以获得传感器测量偏差的准确估计，如果不采用本发明的多项式预测建模方法，在未知系统状态和测量偏差动态特性的情况下是无法同时估计出二者的。Figures 5 and 6 show the estimated results of the sensor measurement bias. As can be seen from the figure, the algorithm of the present invention can obtain an accurate estimate of the sensor measurement deviation while estimating the system state. If the polynomial predictive modeling method of the present invention is not used, the dynamic characteristics of the unknown system state and measurement deviation It is impossible to estimate both at the same time.

以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员，在不脱离本发明构思的前提下，还可以做出若干改变、改进或润饰也应视为本发明的保护范围。The above are only preferred implementations of the present invention, and it should be pointed out that for those of ordinary skill in the art, without departing from the concept of the present invention, some changes, improvements or modifications can also be made, which should also be regarded as the present invention. scope of protection.

Claims

1. A correction method for sensor slow bias fault, characterized in that, comprising the steps:

Step a: establish a state equation and an observation equation for the measured value of the sensor;

Step b: expanding the state equation and observation equation in time;

Step c: use the Kalman filter of the extended state to estimate the deviation, and then subtract the deviation value from the measured value.

2. The correcting method of sensor slow offset fault according to claim 1, is characterized in that, described step a comprises:

Step a1: Describe the state-space model of a linear dynamical system with sensor measurement bias as:

x _t+1 ＝Fx _t +w _t (1)

y _t ＝Φx _t +b _t +v _t (2)

b _t+1 =Ab _t +η _t (3)

Among them, the subscript t represents time, x _t+1 is the N×1-dimensional state vector, F is the state transition matrix, y _t is the p×1-dimensional observation vector, Φ is the measurement matrix, and b _t is the M×1-dimensional The sensor measurement deviation vector, A is the state transition matrix of b _t , w _t is the state noise, v _t is the measurement noise, η _t is the deviation state noise, and the state noise w _t , the measurement noise v _t and the deviation state noise η _t are independent Gaussian white noise with zero mean, and satisfy:

E(w _t )＝0

E E. (({w w}_{j j} {w w}_{t t}^{T T})) = = R R {δ δ}_{jt jt} - - - - - - ((44))

E(v _t )＝0

E E. (({v v}_{j j} {v v}_{t t}^{T T})) = = Q Q {δ δ}_{jt jt} - - - - - - ((55))

E(η _t )=0

E E. (({η η}_{j j} {η η}_{t t}^{T T})) = = G G {δ δ}_{jt jt} - - - - - - ((66))

where R, Q, and G are the variances of the state, measurement, and bias state noise, respectively;

Step a2, expanding the deviation vector into a state vector to obtain a state equation and an observation equation:

The state equation is:

{\overset{~ ~}{x x}}_{t t + + 11} = = \overset{~ ~}{F f} {\overset{~ ~}{x x}}_{t t} + + {\overset{~ ~}{w w}}_{t t}

The observation equation is:

{y the y}_{t t} = = \overset{~ ~}{Φ Φ} {\overset{~ ~}{x x}}_{t t} + + {v v}_{t t}

in,

{\tilde{x}}_{t + 1} = [\begin{matrix} x_{t + 1} \\ b_{t + 1} \end{matrix}],

\tilde{f} = [\begin{matrix} f & 0 \\ 0 & A \end{matrix}],

{\tilde{w}}_{t} = [\begin{matrix} w_{t} \\ η_{t} \end{matrix}],

\tilde{Φ} = [\begin{matrix} Φ & I \end{matrix}],

I is the identity matrix.

3. The correcting method of sensor slow bias fault according to claim 1, is characterized in that, described step b comprises:

Step b1, expand the state variable in time, namely

{X x}_{t t + + 11} = = [[{x x}_{t t + + 11}^{11},, {x x}_{t t}^{11},, . . . . . .,, {x x}_{t t - - K K + + 22}^{11},, {x x}_{t t + + 11}^{22},, {x x}_{t t}^{22},, . . . . . .,, {x x}_{t t - - K K + + 22}^{22},, . . . . . . {x x}_{t t + + 11}^{N N},, {x x}_{t t}^{N N},, . . . . . .,, {x x}_{t t - - K K + + 22}^{N N}]]

where the upper corner is labeled dimension,

Step b2, expand the measured value of the sensor in time, namely

{Y Y}_{t t + + 11} = = [[{y the y}_{t t + + 11}^{11},, {y the y}_{t t}^{11},, . . . . . .,, {y the y}_{t t - - K K + + 22}^{11},, {y the y}_{t t + + 11}^{22},, {y the y}_{t t}^{22},, . . . . . .,, {y the y}_{t t - - K K + + 22}^{22},, . . . . . . {y the y}_{t t + + 11}^{p p},, {y the y}_{t t}^{p p},, . . . . . .,, {y the y}_{t t - - K K + + 22}^{p p}]]

where the upper corner is labeled dimension,

Step b3, make the state transition matrix as:

{\overset{&OverBar; &OverBar;}{F f}}_{ppf ppf} = = [\begin{matrix} {F f}_{ppf ppf}^{((11))} & 00 & . . . . . . & 00 \\ 00 & {F f}_{ppf ppf}^{((22))} & . . . . . . & 00 \\ . . . . . . & . . . . . . & . . . . . . & . . . . . . \\ 00 & . . . . . . & 00 & {F f}_{ppf ppf}^{((N N))} \end{matrix}]

{F f}_{ppf ppf} = = [\begin{matrix} h h ((00)) & h h ((11)) & . . . . . . & h h ((K K - - 11)) \\ 11 & 00 & . . . . . . & 00 \\ . . . . . . & . . . . . . & . . . . . . & . . . . . . \\ 00 & . . . . . . & 11 & 00 \end{matrix}]

Step b4, make the measurement matrix as:

{Φ Φ}_{ppf ppf} = = [\begin{matrix} φ φ ((1,1 1,1)) & 00_{11 \times \times ((K K - - 11))} & φ φ ((1,2 1,2)) & 00_{11 \times \times ((K K - - 11))} & . . . . . . & φ φ ((11,, N N)) & 00_{11 \times \times ((K K - - 11))} \\ 00 & φ φ ((1,1 1,1)) & 00_{11 \times \times ((K K - - 11))} & φ φ ((1,2 1,2)) & . . . . . . & φ φ ((11,, N N)) & 00_{11 \times \times ((K K - - 22))} \\ . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . \\ 00 & 00 & . . . . . . & φ φ ((1,1 1,1)) & 00 & . . . . . . & φ φ ((11,, N N)) \\ . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . \\ φ φ ((p p,, 11)) & 00_{11 \times \times ((K K - - 11))} & φ φ ((p p,, 22)) & 00_{11 \times \times ((K K - - 11))} & . . . . . . & φ φ ((p p,, N N)) & 00_{11 \times \times ((K K - - 11))} \\ . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . & . . . . . . \\ 00 & 00 & . . . . . . & φ φ ((p p,, 11)) & 00 & . . . . . . & φ φ ((p p,, N N)) \end{matrix}]

In step b5, the expanded state equation and observation equation are respectively:

{X x}_{t t + + 11} = = {\overset{&OverBar; &OverBar;}{F f}}_{ppf ppf} {X x}_{t t}

Y _t ＝Φ _ppf X _t +V _t

Wherein, V _t is the expanded measurement noise, and V _t =[v _t , . . . , v _t-K+1 ] ^T .

4. The correcting method of sensor slow offset fault according to claim 1, is characterized in that, described step c comprises:

Step c1, assuming that there is a linear relationship between the measurement deviation of the sensor and the measurable physical quantity c _t , namely:

b _t ≈α _t c _t

The state equation and observation equation are respectively:

{\overset{~ ~}{X x}}_{t t + + 11} = = {\overset{~ ~}{F f}}_{c c} {\overset{~ ~}{X x}}_{t t}

{Y Y}_{t t} = = {\overset{~ ~}{Φ Φ}}_{c c} {\overset{~ ~}{X x}}_{t t} + + {V V}_{t t}

in

{\overset{~ ~}{X x}}_{t t + + 11} = = {[[{x x}_{t t + + 11},, {x x}_{t t},, . . . . . .,, {x x}_{t t - - K K + + 22},, {α α}_{t t + + 11},, {α α}_{t t},, . . . . . .,, {α α}_{t t - - K K + + 22}]]}^{T T}

Y _t =[y _t ,...,y _t-K+1 ] ^T

{\overset{~ ~}{F f}}_{c c} = = [\begin{matrix} {F f}_{ppf ppf} & 00 \\ 00 & {F f}_{ppf ppf} \end{matrix}]

{\overset{~ ~}{Φ Φ}}_{c c} = = [\begin{matrix} Φ Φ & 00 & . . . . . . & 00 & {c c}_{t t} & 00 & . . . . . . & 00 \\ 00 & Φ Φ & . . . . . . & 00 & 00 & {c c}_{t t - - 11} & . . . . . . & 00 \\ 00 & . . . . . . & . . . . . . & . . . . . . & 00 & . . . . . . & . . . . . . & . . . . . . \\ 00 & . . . . . . & 00 & Φ Φ & 00 & . . . . . . & 00 & {c c}_{t t - - K K + + 11} \end{matrix}]

c _t ≠c _t-1 ≠…≠c _t-K+1 .

Step c2, assume the state at t=0

and variance _P0 known,

{\hat{\tilde{x}}}_{t + 1} = {[{\hat{x}}_{t + 1}, {\hat{x}}_{t}, . . ., {\hat{x}}_{t - K + 2}, {\hat{α}}_{t + 1}, {\hat{α}}_{t}, . . ., {\hat{α}}_{t - K + 2}]}^{T},

To t=1, 2, ... (t is a natural number),

For a first-order signal, let

h h ((k k)) = = \frac{44 K K - - 66 k k - - 44}{K K ((K K - - 11))}

For second-order signals, let

h h ((k k)) = = \frac{{99 K K}^{22} + + ((- - 2727 - - 3636 k k)) K K + + {3030 k k}^{22} + + 4242 k k + + 1818}{{K K}^{33} - - 33 {K K}^{22} + + 22 K K}

Step c21, predict:

{\overset{^^}{\overset{~ ~}{X x}}}_{t t + + 11}^{- -} = = {\overset{~ ~}{F f}}_{c c} {\overset{^^}{\overset{~ ~}{X x}}}_{t t}

{P P}_{t t + + 11}^{- -} = = {\overset{~ ~}{F f}}_{c c} {P P}_{t t} {\overset{~ ~}{F f}}_{c c}^{T T}

is the extended state prediction value at time t+1,

is the extended state transition matrix,

is the extended state estimate at time t,

is the state estimation variance prediction value at time t+1,

P _t is the state estimation variance at time t,

Step c22, update

{K K}_{t t + + 11} = = {P P}_{t t + + 11}^{- -} {\overset{~ ~}{H h}}_{c c}^{T T} (({\overset{~ ~}{H h}}_{c c} {P P}_{t t + + 11}^{- -} {\overset{~ ~}{H h}}_{c c}^{T T} + + {R R}_{t t + + 11}))

K _t+1 is the weighted value of new information at time t+1,

is the expanded measurement matrix, R _t+1 is the measurement variance matrix,

{\overset{^^}{\overset{~ ~}{X x}}}_{t t + + 11} = = {\overset{^^}{\overset{~ ~}{X x}}}_{t t + + 11}^{- -} + + {K K}_{t t + + 11} (({Y Y}_{t t + + 11} - - {\overset{~ ~}{H h}}_{c c} {\overset{^^}{\overset{~ ~}{X x}}}_{t t + + 11}^{- -}))

{P P}_{t t + + 11} = = ((I I - - {K K}_{t t + + 11} {\overset{~ ~}{H h}}_{c c})) {P P}_{t t + + 11}^{- -}

Y _t+1 is the expanded measured value, P _t+1 is the variance matrix of state estimation at time t+1;

Step c23, calculate the measurement deviation

{\overset{^^}{b b}}_{t t + + 11} = = {\overset{^^}{α α}}_{t t + + 11} {c c}_{t t + + 11}

Step c3, subtracting the deviation value from the measured value.