[发明专利]一种自校准秩滤波方法

摘要

本发明提供一种自校准秩滤波方法，步骤如下一建立包含未知输入的非高斯、非线性离散系统；二滤波初始化；三对非高斯、非线性离散系统进行时间更新；四结合量测信息，对状态一步预测和一步预测误差方差矩阵进行量测更新；五进行迭代计算；本发明针对非高斯、非线性系统，建立了含有未知输入项的系统模型，使未知输入的影响可以通过数学语言进行描述；基于秩采样理论，本发明不仅可以处理高斯分布非线性系统，同时得以对不服从高斯分布的系统进行状态估计，填补了传统方法在该方面研究的空白；本发明通过引入自校准技术，抵消了状态方程中未知输入的不利影响，减少了滤波发散的现象，提升了滤波精度，增强了系统的鲁棒性。

主权项

一种自校准秩滤波方法，其特征在于：它包含以下五个步骤：步骤一：建立包含未知输入的非高斯、非线性离散系统针对工程实际中所遇到的状态方程受未知输入影响的非高斯、非线性离散系统，特别是具有加性系统噪声和量测噪声的系统，其状态方程和量测方程表示为Xk＝f(Xk‑1)+bk‑1+Wk‑1···········(1)Zk＝h(Xk)+Vk·············(2)式中，Xk表示系统的状态向量，Zk表示系统量测向量，f(·)和h(·)为非线性向量函数,bk表示未知输入，Wk与Vk分别为系统噪声向量和量测噪声向量，其方差矩阵分别为Qk和Rk，并且满足E[Wk]=0Cov[Wk,Wj]=E[WkWjT]=QkδkjE[Vk]=0Cov[Vk,Vj]=E[VkVjT]=RkδkjCov[Wk,Vj]=E[WkVjT]=0...(3)]]>式中，Cov[·]为协方差，E[·]为数学期望，δkj为δ函数，当k＝j时，δkj＝1，当k≠j时，δkj＝0；步骤二：滤波初始化对由式(1)和式(2)所组成系统的状态值和误差方差矩阵进行初始化对系统状态和误差方差矩阵分别赋予初值X^0=E[X0]...(4)]]>P0=E[(X0-X^0)(X0-X^0)T]...(5)]]>步骤三：对非高斯、非线性离散系统进行时间更新设k‑1时刻的状态估计值和误差方差矩阵分别为和Pk‑1，基于它们对系统进行时间更新，即计算k时刻的状态一步预测和一步预测误差方差矩阵Pk/k‑1；基于秩滤波的一般性递推公式，在时间更新过程中首先需要计算秩采样点集{χi}当k＝1,2时χk-1,i=X^k-1+up1(Pk-1)i,i=1,...,nX^k-1-up1(Pk-1)i-n,i=n+1,...,2nX^k-1+up2(Pk-1)i-2n,i=2n+1,...,3nX^k-1-up2(Pk-1)i-3n,i=3n+1,...,4n...(6)]]>Xk/k-1i=f(χk-1,i)...(7)]]>当k≥3时χk-2,i=X^k-2+up1(Pk-2)i,i=1,...,nX^k-2-up1(Pk-2)i-n,i=n+1,...,2nX^k-2+up2(Pk-2)i-2n,i=2n+1,...,3nX^k-2-up2(Pk-2)i-3n,i=3n+1,...,4n...(8)]]>χk-1,i=X^k-1+up1(Pk-1)i,i=1,...,nX^k-1-up1(Pk-1)i-n,i=n+1,...,2nX^k-1+up2(Pk-1)i-2n,i=2n+1,...,3nX^k-1-up2(Pk-1)i-3n,i=3n+1,...,4n...(9)]]>Xk/k‑1,i＝f(χk‑1,i)+χk‑1,i‑f(χk‑2,i)········(10)式中，n为状态向量Xk的维度，标准正态偏量表示Pk‑1平方根的第i列向量；在式(10)中，用到了自校准技术的基本假设：由于工程中相邻两步滤波之间的时间间隔很短，所以它们的未知输入也往往相等及近似相等，即有bk≈bk‑1；由式(7)和式(10)计算系统状态一步预测值X^k/k-1=14nΣi=14nXk/k-1,i...(11)]]>由式(7)、式(10)和式(11)计算一步预测误差方差矩阵Pk/k-1=1ωΣi=14n{(Xk/k-1,i-X^k/k-1)(Xk/k-1,i-X^k/k-1)T}+Qk-1...(12)]]>式中，ω为协方差权重系数，其计算公式为ω=2(up12+up22)=3.0105...(13);]]>步骤四：结合量测信息，对状态一步预测和一步预测误差方差矩阵进行量测更新首先，基于状态一步预测和一步预测误差方差矩阵Pk/k‑1，根据秩滤波的一般性递推公式对系统状态进行重采样，得到新的秩采样点集{χk/k‑1,i}χk/k-1,i=X^k/k-1+up1(Pk/k-1)i,i=1,...,nX^k/k-1-up1(Pk/k-1)i-n,i=n+1,...,2nX^k/k-1+up2(Pk/k-1)i-2n,i=2n+1,...,3nX^k/k-1-up2(Pk/k-1)i-3n,i=3n+1,...,4n...(14)]]>然后，根据{χk/k‑1,i}计算量测估计值Zk/k‑1,i＝h(χk/k‑1,i)···········(15)Z^k/k-1=14nΣi=14nZk/k-1,i...(16)]]>由式(16)和式(17)计算量测误差方差矩阵PZZPZZ=1ωΣi=14n{(Zk/k-1,i-Z^k/k-1)(Zk/k-1,i-Z^k/k-1)T}+Rk...(17)]]>由式(11)、式(15)、式(16)和式(17)计算误差协方差矩阵PXZPXZ=1ωΣi=14n{(χk/k-1,i-X^k/k-1)(Zk/k-1,i-Z^k/k-1)T}.....(18)]]>则，由量测更新得到最终的状态估计值和状态估计误差方差矩阵X^k=X^k/k-1+Kk(Zk-Z^k/k-1)...(19)]]>Pk=Pk/k-1-KkPZZKkT...(20)]]>式中，根据秩滤波的一般性递推公式，步骤五：迭代计算根据k时刻的状态估计值和误差方差矩阵Pk，重复步骤(二)和(三)，进而得到k+1时刻的状态估计值和误差方差矩阵，往复迭代，直至滤波过程结束。