卡尔曼滤波与组合导航原理(十二)扩展卡尔曼滤波:EKF、二阶EKF、迭代EKF

news2025/1/10 10:53:10

文章目录

    • 一、多元向量的泰勒级数展开
    • 二、扩展Kalman滤波
    • 三、二阶滤波
    • 四、迭代EKF滤波

一、多元向量的泰勒级数展开

{ y 1 = f 1 ( X ) = f 1 ( x 1 , x 2 , ⋯ x n ) y 2 = f 2 ( X ) = f 2 ( x 1 , x 2 , ⋯ x n ) ⋮ y m = f m ( X ) = f m ( x 1 , x 2 , ⋯ x n ) \left\{\begin{array}{c} y_{1}=f_{1}(\boldsymbol{X})=f_{1}\left(x_{1}, x_{2}, \cdots x_{n}\right) \\ y_{2}=f_{2}(\boldsymbol{X})=f_{2}\left(x_{1}, x_{2}, \cdots x_{n}\right) \\ \quad \vdots \\ y_{m}=f_{m}(\boldsymbol{X})=f_{m}\left(x_{1}, x_{2}, \cdots x_{n}\right) \end{array}\right. y1=f1(X)=f1(x1,x2,xn)y2=f2(X)=f2(x1,x2,xn)ym=fm(X)=fm(x1,x2,xn)

上面线性方差组,写成向量形式就是:
Y = f ( X ) \boldsymbol{Y} = \boldsymbol{f} \left( \boldsymbol{X}\right) Y=f(X)
泰勒展开式为:
Y = f ( X ( 0 ) ) + ∑ i = 1 ∞ 1 i ! ( ∇ T ⋅ δ X ) i f ( X ( 0 ) ) \boldsymbol{Y}=\boldsymbol{f}\left(\boldsymbol{X}^{(0)}\right)+\sum_{i=1}^{\infty} \frac{1}{i !}\left(\nabla^{\mathrm{T}} \cdot \delta \boldsymbol{X}\right)^{i} \boldsymbol{f}\left(\boldsymbol{X}^{(0)}\right) Y=f(X(0))+i=1i!1(TδX)if(X(0))
其中 δ X = X − X ( 0 ) \delta \boldsymbol{X}=\boldsymbol{X}-\boldsymbol{X}^{(0)} δX=XX(0) ∇ = [ ∂ ∂ x 1 ∂ ∂ x 2 ⋯ ∂ ∂ x n ] T \nabla=\left[\begin{array}{llll}\frac{\partial}{\partial x{1}} & \frac{\partial}{\partial x{2}} & \cdots & \frac{\partial}{\partial x_{n}}\end{array}\right]^{\mathrm{T}} =[x1x2xn]T

一阶导代表梯度、二阶导代表曲率半径

一阶展开项的详细展开式

标量的求导运算对向量求导是对向量里的每一个元素求导, m m m 列多元函数对 n n n 行的自变量求导:
1 1 ( ∇ T ⋅ δ X ) ′ f ( X ( 0 ) ) = ( ∂ ∂ x 1 δ x 1 + ∂ ∂ x 2 δ x 2 + ⋯ + ∂ ∂ x n δ x n ) f ( X ) ∣ X = X ( 0 ) = [ ∂ f 1 ( X ) ∂ x 1 δ x 1 + ∂ f 1 ( X ) ∂ x 2 δ x 2 + ⋯ + ∂ f 1 ( X ) ∂ x n δ x n ∂ f 2 ( X ) ∂ x 1 δ x 1 + ∂ f 2 ( X ) ∂ x 2 δ x 2 + ⋯ + ∂ f 2 ( X ) ∂ x n δ x n ∂ f m ( X ) ∂ x 1 δ x 1 + ∂ f m ( X ) ∂ x 2 δ x 2 + ⋯ + ∂ f m ( X ) ∂ x n δ x n ] X = X ( 0 ) = [ ∂ f 1 ( X ) ∂ x 1 ∂ f 1 ( X ) ∂ x 2 ⋯ ∂ f 1 ( X ) ∂ x n ∂ f 2 ( X ) ∂ x 1 ∂ f 2 ( X ) ∂ x 2 ⋯ ∂ f 2 ( X ) ∂ x n ⋮ ⋮ ⋱ ⋮ ∂ f m ( X ) ∂ x 1 ∂ f m ( X ) ∂ x 2 ⋯ ∂ f m ( X ) ∂ x n ] δ [ δ x 1 δ x 2 ⋮ δ x n ] X = X ( 0 ) = J ( f ( X ) ) δ X ∣ X = X ( 0 ) \begin{array}{l} \frac{1}{1}\left(\nabla^{\mathrm{T}} \cdot \delta \boldsymbol{X}\right)^{\prime} \boldsymbol{f}\left(\boldsymbol{X}^{(0)}\right)=\left.\left(\frac{\partial}{\partial x_{1}} \delta x_{1}+\frac{\partial}{\partial x_{2}} \delta x_{2}+\cdots+\frac{\partial}{\partial x_{n}} \delta x_{n}\right) \boldsymbol{f}(\boldsymbol{X})\right|_{\boldsymbol{X}=\boldsymbol{X}^{(0)}} \\ =\left[\begin{array}{c} \frac{\partial f_{1}(\boldsymbol{X})}{\partial x_{1}} \delta x_{1}+\frac{\partial f_{1}(\boldsymbol{X})}{\partial x_{2}} \delta x_{2}+\cdots+\frac{\partial f_{1}(\boldsymbol{X})}{\partial x_{n}} \delta x_{n} \\ \frac{\partial f_{2}(\boldsymbol{X})}{\partial x_{1}} \delta x_{1}+\frac{\partial f_{2}(\boldsymbol{X})}{\partial x_{2}} \delta x_{2}+\cdots+\frac{\partial f_{2}(\boldsymbol{X})}{\partial x_{n}} \delta x_{n} \\ \frac{\partial f_{m}(\boldsymbol{X})}{\partial x_{1}} \delta x_{1}+\frac{\partial f_{m}(\boldsymbol{X})}{\partial x_{2}} \delta x_{2}+\cdots+\frac{\partial f_{m}(\boldsymbol{X})}{\partial x_{n}} \delta x_{n} \end{array}\right]_{\boldsymbol{X}=\boldsymbol{X}^{(0)}}=\left[\begin{array}{cccc} \frac{\partial f_{1}(\boldsymbol{X})}{\partial x_{1}} & \frac{\partial f_{1}(\boldsymbol{X})}{\partial x_{2}} & \cdots & \frac{\partial f_{1}(\boldsymbol{X})}{\partial x_{n}} \\ \frac{\partial f_{2}(\boldsymbol{X})}{\partial x_{1}} & \frac{\partial f_{2}(\boldsymbol{X})}{\partial x_{2}} & \cdots & \frac{\partial f_{2}(\boldsymbol{X})}{\partial x_{n}} \\ \vdots & \vdots & \ddots & \vdots \\ \frac{\partial f_{m}(\boldsymbol{X})}{\partial x_{1}} & \frac{\partial f_{m}(\boldsymbol{X})}{\partial x_{2}} & \cdots & \frac{\partial f_{m}(\boldsymbol{X})}{\partial x_{n}} \end{array}\right]^{\delta}\left[\begin{array}{c} \delta x_{1} \\ \delta x_{2} \\ \vdots \\ \delta x_{n} \end{array}\right]_{X=X^{(0)}} \\ =\left.J(\boldsymbol{f}(\boldsymbol{X})) \delta \boldsymbol{X}\right|_{X=X^{(0)}} \quad \\ \end{array} 11(TδX)f(X(0))=(x1δx1+x2δx2++xnδxn)f(X) X=X(0)= x1f1(X)δx1+x2f1(X)δx2++xnf1(X)δxnx1f2(X)δx1+x2f2(X)δx2++xnf2(X)δxnx1fm(X)δx1+x2fm(X)δx2++xnfm(X)δxn X=X(0)= x1f1(X)x1f2(X)x1fm(X)x2f1(X)x2f2(X)x2fm(X)xnf1(X)xnf2(X)xnfm(X) δ δx1δx2δxn X=X(0)=J(f(X))δXX=X(0)
其中 ∫ J ( f ( X ) ) = ∂ f ( X ) ∂ X T \int J(\boldsymbol{f}(\boldsymbol{X}))=\frac{\partial \boldsymbol{f}(\boldsymbol{X})}{\partial \boldsymbol{X}^{\mathrm{T}}} J(f(X))=XTf(X) 可称之为雅各比矩阵, m m m 列多元函数对 n n n 行的自变量求一阶导得到的矩阵。

二阶展开项详细表达式
1 2 ( ∇ T ⋅ δ X ) 2 f ( X ( 0 ) ) = 1 2 ( ∂ ∂ x 1 δ x 1 + ∂ ∂ x 2 δ x 2 + ⋯ + ∂ ∂ x n δ x n ) 2 f ( X ) ∣ X = X ( 0 ) = 1 2 tr ⁡ ( [ ∂ 2 ∂ x 1 2 ∂ 2 ∂ x 1 ∂ x 2 ⋯ ∂ 2 ∂ x 1 ∂ x n ∂ 2 ∂ x 2 ∂ x 1 ∂ 2 ∂ x 2 2 ⋯ ∂ 2 ∂ x 2 ∂ x n ⋮ ⋮ ⋱ ⋮ ∂ 2 ∂ x n ∂ x 1 ∂ 2 ∂ x n ∂ x 2 ⋯ ∂ 2 ∂ x n 2 ] [ δ x 1 2 δ x 1 δ x 2 ⋯ δ x 1 δ x n δ x 2 δ x 1 δ x 2 2 ⋯ δ x 2 δ x n ⋮ ⋮ ⋱ ⋮ δ x n δ x 1 δ x n δ x 2 ⋯ δ x n 2 ] ) f ( X ) ∣ X = X ( 0 ) = 1 2 tr ⁡ ( ∇ ∇ T ⋅ δ X δ X T ) f ( X ) ∣ X = X ( 0 ) = 1 2 tr ⁡ ( ∇ ∇ T ⋅ δ X δ X T ) ∑ j = 1 m e j f j ( X ) ∣ X = X ( 0 ) = 1 2 ∑ j = 1 m e j tr ⁡ ( ∇ ∇ T ⋅ δ X δ X T ) f j ( X ) ∣ X = X ( 0 ) = 1 2 ∑ j = 1 m e j tr ⁡ ( ∇ ∇ T f j ( X ) ∣ X = X ( 0 ) ⋅ δ X δ X T ) = 1 2 ∑ j = 1 m e j tr ⁡ ( H ( f j ( X ) ∣ X = X ( 0 ) ⋅ δ X δ X T ) \begin{array}{l} \frac{1}{2}\left(\nabla^{\mathrm{T}} \cdot \delta \boldsymbol{X}\right)^{2} \boldsymbol{f}\left(\boldsymbol{X}^{(0)}\right)=\left.\frac{1}{2}\left(\frac{\partial}{\partial x_{1}} \delta x_{1}+\frac{\partial}{\partial x_{2}} \delta x_{2}+\cdots+\frac{\partial}{\partial x_{n}} \delta x_{n}\right)^{2} \boldsymbol{f}(\boldsymbol{X})\right|_{X=X^{(0)}} \\ =\left.\frac{1}{2} \operatorname{tr}\left(\left[\begin{array}{cccc} \frac{\partial^{2}}{\partial x_{1}^{2}} & \frac{\partial^{2}}{\partial x_{1} \partial x_{2}} & \cdots & \frac{\partial^{2}}{\partial x_{1} \partial x_{n}} \\ \frac{\partial^{2}}{\partial x_{2} \partial x_{1}} & \frac{\partial^{2}}{\partial x_{2}^{2}} & \cdots & \frac{\partial^{2}}{\partial x_{2} \partial x_{n}} \\ \vdots & \vdots & \ddots & \vdots \\ \frac{\partial^{2}}{\partial x_{n} \partial x_{1}} & \frac{\partial^{2}}{\partial x_{n} \partial x_{2}} & \cdots & \frac{\partial^{2}}{\partial x_{n}^{2}} \end{array}\right]\left[\begin{array}{cccc} \delta x_{1}^{2} & \delta x_{1} \delta x_{2} & \cdots & \delta x_{1} \delta x_{n} \\ \delta x_{2} \delta x_{1} & \delta x_{2}^{2} & \cdots & \delta x_{2} \delta x_{n} \\ \vdots & \vdots & \ddots & \vdots \\ \delta x_{n} \delta x_{1} & \delta x_{n} \delta x_{2} & \cdots & \delta x_{n}^{2} \end{array}\right]\right) \boldsymbol{f}(\boldsymbol{X})\right|_{X=\mathbf{X}^{(0)}} \\ \begin{array}{l} =\left.\frac{1}{2} \operatorname{tr}\left(\nabla \nabla^{\mathrm{T}} \cdot \delta \boldsymbol{X} \delta \boldsymbol{X}^{\mathrm{T}}\right) \boldsymbol{f}(\boldsymbol{X})\right|_{X=X^{(0)}}=\left.\frac{1}{2} \operatorname{tr}\left(\nabla \nabla^{\mathrm{T}} \cdot \delta \boldsymbol{X} \delta \boldsymbol{X}^{\mathrm{T}}\right) \sum_{j=1}^{m} \boldsymbol{e}_{j} f_{j}(\boldsymbol{X})\right|_{X=X^{(0)}} \\ =\left.\frac{1}{2} \sum_{j=1}^{m} \boldsymbol{e}_{j} \operatorname{tr}\left(\nabla \nabla^{\mathrm{T}} \cdot \delta \boldsymbol{X} \delta \boldsymbol{X}^{\mathrm{T}}\right) f_{j}(\boldsymbol{X})\right|_{X=X^{(0)}}=\frac{1}{2} \sum_{j=1}^{m} \boldsymbol{e}_{j} \operatorname{tr}\left(\left.\nabla \nabla^{\mathrm{T}} f_{j}(\boldsymbol{X})\right|_{X=X^{(0)}} \cdot \delta \boldsymbol{X} \delta \boldsymbol{X}^{\mathrm{T}}\right) \end{array} \\ =\frac{1}{2} \sum_{j=1}^{m} \boldsymbol{e}_{j} \operatorname{tr}\left(\mathscr{H}\left(\left.f_{j}(\boldsymbol{X})\right|_{X=X^{(0)}} \cdot \delta \boldsymbol{X} \delta \boldsymbol{X}^{\mathrm{T}}\right) \quad\right. \\ \end{array} 21(TδX)2f(X(0))=21(x1δx1+x2δx2++xnδxn)2f(X) X=X(0)=21tr x122x2x12xnx12x1x22x222xnx22x1xn2x2xn2xn22 δx12δx2δx1δxnδx1δx1δx2δx22δxnδx2δx1δxnδx2δxnδxn2 f(X) X=X(0)=21tr(TδXδXT)f(X) X=X(0)=21tr(TδXδXT)j=1mejfj(X) X=X(0)=21j=1mejtr(TδXδXT)fj(X) X=X(0)=21j=1mejtr(Tfj(X) X=X(0)δXδXT)=21j=1mejtr(H(fj(X)X=X(0)δXδXT)

Y = f ( X ( 0 ) ) + J ( f ( X ( 0 ) ) ) δ X + 1 2 ∑ i = 1 m e i tr ⁡ ( H ( f i ( X ( 0 ) ) ) ⋅ δ X δ X T ) + O 3 \boldsymbol{Y}=\boldsymbol{f}\left(\boldsymbol{X}^{(0)}\right)+J\left(f\left(X^{(0)}\right)\right) \delta \boldsymbol{X}+\frac{1}{2} \sum_{i=1}^{m} e_{i} \operatorname{tr}\left(\mathscr{H}\left(f_{i}\left(X^{(0)}\right)\right) \cdot \delta \boldsymbol{X} \delta \boldsymbol{X}^{\mathrm{T}}\right)+O^{3} Y=f(X(0))+J(f(X(0)))δX+21i=1meitr(H(fi(X(0)))δXδXT)+O3

其中 H \mathscr{H} H 为海森矩阵, m m m 列多元函数对 n n n 行的自变量求二阶导得到的矩阵,m 个向量函数有 m 个海森矩阵。

二、扩展Kalman滤波

遇到非线性的状态空间模型时,对状态方程和量测方程做线性化处理,取泰勒展开的一阶项。

状态空间模型(加性噪声,噪声和状态直接是相加的)

这里简单的用加性噪声模型表示,一般形式为 X k = f ( X k − 1 , W k − 1 , k − 1 ) \boldsymbol{X}_{k}=\boldsymbol{f}\left(\boldsymbol{X}_{k-1}, \boldsymbol{W}_{k-1}, k-1\right) Xk=f(Xk1,Wk1,k1) ,非线性更强,更一般,可以表示噪声和状态之间复杂的关系

{ X k = f ( X k − 1 ) + Γ k − 1 W k − 1 Z k = h ( X k ) + V k \left\{\begin{array}{l}\boldsymbol{X}_{k}=\boldsymbol{f}\left(\boldsymbol{X}_{k-1}\right)+\boldsymbol{\Gamma}_{k-1} \boldsymbol{W}_{k-1} \\ \boldsymbol{Z}_{k}=\boldsymbol{h}\left(\boldsymbol{X}_{k}\right)+\boldsymbol{V}_{k}\end{array}\right. {Xk=f(Xk1)+Γk1Wk1Zk=h(Xk)+Vk

选择 k − 1 k-1 k1 时刻参考值 X k − 1 n \boldsymbol{X}_{k-1}^{n} Xk1n ,参考点和真实值偏差 δ X k − 1 = X k − 1 − X k − 1 n \delta \boldsymbol{X}_{k-1}=\boldsymbol{X}_{k-1}-\boldsymbol{X}_{k-1}^{n} δXk1=Xk1Xk1n

状态一步预测:上一时刻参考点带入状态方程 X k / k − 1 n = f ( X k − 1 n ) ​ \boldsymbol{X}_{k / k-1}^{n}=\boldsymbol{f}\left(\boldsymbol{X}_{k-1}^{n}\right)​ Xk/k1n=f(Xk1n) ,预测值的偏差 δ X k = X k − X k / k − 1 n ​ \delta \boldsymbol{X}_{k}=\boldsymbol{X}_{k}-\boldsymbol{X}_{k / k-1}^{n}​ δXk=XkXk/k1n

量测一步预测 Z k / k − 1 n = h ( X k / k − 1 n ) \boldsymbol{Z}_{k / k-1}^{n}=\boldsymbol{h}\left(\boldsymbol{X}_{k / k-1}^{n}\right) Zk/k1n=h(Xk/k1n) ,偏差 δ Z k = Z k − Z k / k − 1 n \delta \boldsymbol{Z}_{k}=\boldsymbol{Z}_{k}-\boldsymbol{Z}_{k / k-1}^{n} δZk=ZkZk/k1n

展开得状态偏差方程
X k ≈ f ( X k − 1 n ) + J ( f ( X k − 1 n ) ) ( X k − 1 − X k − 1 n ) + Γ k − 1 W k − 1 X k − f ( X k − 1 n ) ≈ Φ k / k − 1 n ( X k − 1 − X k − 1 n ) + Γ k − 1 W k − 1 δ X k = Φ k / k − 1 n δ X k − 1 + Γ k − 1 W k − 1 \begin{array}{ll} & \boldsymbol{X}_{k} \approx \boldsymbol{f}\left(\boldsymbol{X}_{k-1}^{n}\right)+\boldsymbol{J}\left(\boldsymbol{f}\left(\boldsymbol{X}_{k-1}^{n}\right)\right)\left(\boldsymbol{X}_{k-1}-\boldsymbol{X}_{k-1}^{n}\right)+\boldsymbol{\Gamma}_{k-1} \boldsymbol{W}_{k-1} \\ & \boldsymbol{X}_{k}-\boldsymbol{f}\left(\boldsymbol{X}_{k-1}^{n}\right) \approx \boldsymbol{\Phi}_{k / k-1}^{n}\left(\boldsymbol{X}_{k-1}-\boldsymbol{X}_{k-1}^{n}\right)+\boldsymbol{\Gamma}_{k-1} \boldsymbol{W}_{k-1} \\ & \delta \boldsymbol{X}_{k}=\boldsymbol{\Phi}_{k / k-1}^{n} \delta \boldsymbol{X}_{k-1}+\boldsymbol{\Gamma}_{k-1} \boldsymbol{W}_{k-1} \end{array} Xkf(Xk1n)+J(f(Xk1n))(Xk1Xk1n)+Γk1Wk1Xkf(Xk1n)Φk/k1n(Xk1Xk1n)+Γk1Wk1δXk=Φk/k1nδXk1+Γk1Wk1
展开得量测偏差方程
Z k ≈ h ( X k / k − 1 n ) + J ( h ( X k / k − 1 n ) ) ( X k − X k / k − 1 n ) + V k δ Z k = H k n δ X k + V k \begin{array}{ll} & \boldsymbol{Z}_{k} \approx \boldsymbol{h}\left(\boldsymbol{X}_{k / k-1}^{n}\right)+\boldsymbol{J}\left(\boldsymbol{h}\left(\boldsymbol{X}_{k / k-1}^{n}\right)\right)\left(\boldsymbol{X}_{k}-\boldsymbol{X}_{k / k-1}^{n}\right)+\boldsymbol{V}_{k} \\ & \delta \boldsymbol{Z}_{k}=\boldsymbol{H}_{k}^{n} \delta \boldsymbol{X}_{k}+\boldsymbol{V}_{k} \end{array} Zkh(Xk/k1n)+J(h(Xk/k1n))(XkXk/k1n)+VkδZk=HknδXk+Vk
得偏差状态空间估计模型
{ δ X k = Φ k ∣ k − 1 n δ X k − 1 + Γ k − 1 W k − 1 δ Z k = H k n δ X k + V k \left\{\begin{array}{l}\delta \boldsymbol{X}_{k}=\boldsymbol{\Phi}_{k \mid k-1}^{n} \delta \boldsymbol{X}_{k-1}+\boldsymbol{\Gamma}_{k-1} \boldsymbol{W}_{k-1} \\ \delta \boldsymbol{Z}_{k}=\boldsymbol{H}_{k}^{n} \delta \boldsymbol{X}_{k}+\boldsymbol{V}_{k}\end{array}\right. {δXk=Φkk1nδXk1+Γk1Wk1δZk=HknδXk+Vk
偏差量是线性的了,可以用Kalman滤波进行处理:

滤波更新方程 δ X ^ k = Φ k / k − 1 n δ X ^ k − 1 + K k n ( δ Z k − H k n Φ k / k − 1 n δ X ^ k − 1 ) {\color{red}\delta \hat{\boldsymbol{X}}_{k}}=\boldsymbol{\Phi}_{k / k-1}^{n} {\color{green}\delta \hat{\boldsymbol{X}}_{k-1}}+\boldsymbol{K}_{k}^{n}\left(\delta \boldsymbol{Z}_{k}-\boldsymbol{H}_{k}^{n} \boldsymbol{\Phi}_{k / k-1}^{n} {\color{green}\delta \hat{\boldsymbol{X}}_{k-1}}\right) δX^k=Φk/k1nδX^k1+Kkn(δZkHknΦk/k1nδX^k1)
⟶ δ X ^ k = X ^ k − X k / k − 1 n δ X ^ k − 1 = X ^ k − 1 − X k − 1 n X ^ k = X k / k − 1 n + K k n ( Z k − Z k / k − 1 n ) + ( I − K k n H k n ) Φ k / k − 1 n ( X ^ k − 1 − X k − 1 n ) ⟶ X k − 1 n = X ^ k − 1 = X k / k − 1 n + K k n ( Z k − Z k / k − 1 n ) \begin{array}{l} \underset{\delta \hat{\boldsymbol{X}}_{k-1}=\hat{\boldsymbol{X}}_{k-1}-\boldsymbol{X}_{k-1}^{n}}{\underset{\delta \hat{\boldsymbol{X}}_{k}=\hat{\boldsymbol{X}}_{k}-\boldsymbol{X}_{k / k-1}^{n}}{\longrightarrow}} \hat{\boldsymbol{X}}_{k}=\boldsymbol{X}_{k / k-1}^{n}+\boldsymbol{K}_{k}^{n}\left(\boldsymbol{Z}_{k}-\boldsymbol{Z}_{k / k-1}^{n}\right)+\left(\boldsymbol{I}-\boldsymbol{K}_{k}^{n} \boldsymbol{H}_{k}^{n}\right) \boldsymbol{\Phi}_{k / k-1}^{n}\left(\hat{\boldsymbol{X}}_{k-1}-\boldsymbol{X}_{k-1}^{n}\right) \\ \stackrel{\boldsymbol{X}_{k-1}^{n}=\hat{\boldsymbol{X}}_{k-1}}{\longrightarrow}=\boldsymbol{X}_{k / k-1}^{n}+\boldsymbol{K}_{k}^{n}\left(\boldsymbol{Z}_{k}-\boldsymbol{Z}_{k / k-1}^{n}\right) \end{array} δX^k1=X^k1Xk1nδX^k=X^kXk/k1nX^k=Xk/k1n+Kkn(ZkZk/k1n)+(IKknHkn)Φk/k1n(X^k1Xk1n)Xk1n=X^k1=Xk/k1n+Kkn(ZkZk/k1n)
δ X ^ k − 1 = X ^ k − 1 − X k − 1 n {\color{green}\delta \hat{\boldsymbol{X}}_{k-1}=\hat{\boldsymbol{X}}_{k-1}-\boldsymbol{X}_{k-1}^{n}} δX^k1=X^k1Xk1n δ X ^ k = X ^ k − X k / k − 1 n {\color{red}\delta \hat{\boldsymbol{X}}_{k}=\hat{\boldsymbol{X}}_{k}-\boldsymbol{X}_{k / k-1}^{n}} δX^k=X^kXk/k1n 带入得:
X ^ k = X k / k − 1 n + K k n ( Z k − Z k / k − 1 n ) + ( I − K k n H k n ) Φ k / k − 1 n ( X ^ k − 1 − X k − 1 n ) \hat{\boldsymbol{X}}_{k}=\boldsymbol{X}_{k / k-1}^{n}+\boldsymbol{K}_{k}^{n}\left(\boldsymbol{Z}_{k}-\boldsymbol{Z}_{k / k-1}^{n}\right)+{\color{pink} \left(\boldsymbol{I}-\boldsymbol{K}_{k}^{n} \boldsymbol{H}_{k}^{n}\right) \boldsymbol{\Phi}_{k / k-1}^{n}\left(\hat{\boldsymbol{X}}_{k-1}-\boldsymbol{X}_{k-1}^{n}\right)} X^k=Xk/k1n+Kkn(ZkZk/k1n)+(IKknHkn)Φk/k1n(X^k1Xk1n)
此公式就相当于对状态进行操作。前面部分就是Kalman滤波方程。后面多出的粉色部分包含了参考值,令 X k − 1 n = X ^ k − 1 \boldsymbol{X}_{k-1}^{n}=\hat{\boldsymbol{X}}_{k-1} Xk1n=X^k1 ,参考点选成上一时刻的估计值,后面一部分就为 0 0 0 ,得
X k / k − 1 n + K k n ( Z k − Z k / k − 1 n ) \boldsymbol{X}_{k / k-1}^{n}+\boldsymbol{K}_{k}^{n}\left(\boldsymbol{Z}_{k}-\boldsymbol{Z}_{k / k-1}^{n}\right) Xk/k1n+Kkn(ZkZk/k1n)
EKF滤波公式汇总
{ X ^ k / k − 1 = f ( X ^ k − 1 ) P k / k − 1 = Φ k / k − 1 P k − 1 Φ k / k − 1 T + Γ k − 1 Q k − 1 Γ k − 1 T Φ k / k − 1 = J ( f ( X ^ k − 1 ) ) K k = P k / k − 1 H k T ( H k P k / k − 1 H k T + R k ) − 1 H k = J ( h ( X ^ k / k − 1 ) ) X ^ k = X ^ k / k − 1 + K k [ Z k − h ( X ^ k / k − 1 ) ] P k = ( I − K k H k ) P k / k − 1 \left\{\begin{array}{ll}\hat{\boldsymbol{X}}_{k / k-1}=\boldsymbol{f}\left(\hat{\boldsymbol{X}}_{k-1}\right) & \\ \boldsymbol{P}_{k / k-1}={\color{red}\boldsymbol{\Phi}_{k / k-1}} \boldsymbol{P}_{k-1} {\color{red}\boldsymbol{\Phi}_{k / k-1}^{\mathrm{T}}}+\boldsymbol{\Gamma}_{k-1} \boldsymbol{Q}_{k-1} \boldsymbol{\Gamma}_{k-1}^{\mathrm{T}} & {\color{red}\boldsymbol{\Phi}_{k / k-1}=\boldsymbol{J}\left(\boldsymbol{f}\left(\hat{\boldsymbol{X}}_{k-1}\right)\right)} \\ \boldsymbol{K}_{k}=\boldsymbol{P}_{k / k-1} {\color{green}\boldsymbol{H}_{k}^{\mathrm{T}}}\left({\color{green}\boldsymbol{H}_{k}} \boldsymbol{P}_{k / k-1} {\color{green}\boldsymbol{H}_{k}^{\mathrm{T}}}+\boldsymbol{R}_{k}\right)^{-1} & {\color{green}\boldsymbol{H}_{k}=\boldsymbol{J}\left(\boldsymbol{h}\left(\hat{\boldsymbol{X}}_{k / k-1}\right)\right)} \\ \hat{\boldsymbol{X}}_{k}=\hat{\boldsymbol{X}}_{k / k-1}+\boldsymbol{K}_{k}\left[\boldsymbol{Z}_{k}-\boldsymbol{h}\left(\hat{\boldsymbol{X}}_{k / k-1}\right)\right] & \\ \boldsymbol{P}_{k}=\left(\boldsymbol{I}-\boldsymbol{K}_{k} {\color{green}\boldsymbol{H}_{k}}\right) \boldsymbol{P}_{k / k-1} & \end{array}\right. X^k/k1=f(X^k1)Pk/k1=Φk/k1Pk1Φk/k1T+Γk1Qk1Γk1TKk=Pk/k1HkT(HkPk/k1HkT+Rk)1X^k=X^k/k1+Kk[Zkh(X^k/k1)]Pk=(IKkHk)Pk/k1Φk/k1=J(f(X^k1))Hk=J(h(X^k/k1))

  • 状态转移矩阵是状态方程在 k − 1 k-1 k1 处的一阶偏导组成的雅可比矩阵。
  • 设计矩阵是量测方程在 k − 1 k-1 k1 处的一阶偏导组成的雅可比矩阵。

三、二阶滤波

状态方程的二阶泰勒展开

X k = f ( X k − 1 n ) + Φ k / k − 1 n δ X + 1 2 ∑ i = 1 n e i tr ⁡ ( D f i − ⋅ δ X δ X T ) + Γ k − 1 W k − 1 \boldsymbol{X}_{k}=\boldsymbol{f}\left(\boldsymbol{X}_{k-1}^{n}\right)+\boldsymbol{\Phi}_{k / k-1}^{n} \delta \boldsymbol{X}+\frac{1}{2} \sum_{i=1}^{n} \boldsymbol{e}_{i} \operatorname{tr}\left(\boldsymbol{D}_{f i}^{-} \cdot \boldsymbol{\delta} \boldsymbol{X} \delta \boldsymbol{X}^{\mathrm{T}}\right)+\boldsymbol{\Gamma}_{k-1} \boldsymbol{W}_{k-1} Xk=f(Xk1n)+Φk/k1nδX+21i=1neitr(DfiδXδXT)+Γk1Wk1
状态预测:上面的方程两边求均值得:
X ^ k / k − 1 = E [ f ( X k − 1 n ) + Φ k / k − 1 n δ X + 1 2 ∑ i = 1 n e i tr ⁡ ( D f i ⋅ δ X δ X T ) + Γ k − 1 W k − 1 ] = f ( X k − 1 n ) + Φ k / k − 1 n E [ δ X ] + 1 2 ∑ i = 1 n e i tr ⁡ ( D f i ⋅ E [ δ X δ X T ] ) + Γ k − 1 E [ W k − 1 ] = f ( X k − 1 n ) + 1 2 ∑ i = 1 n e i tr ⁡ ( D f i P k − 1 ) \begin{aligned} \hat{\boldsymbol{X}}_{k / k-1} & =\mathrm{E}\left[\boldsymbol{f}\left(\boldsymbol{X}_{k-1}^{n}\right)+\boldsymbol{\Phi}_{k / k-1}^{n} \delta \boldsymbol{X}+\frac{1}{2} \sum_{i=1}^{n} \boldsymbol{e}_{i} \operatorname{tr}\left(\boldsymbol{D}_{f i} \cdot \boldsymbol{\delta} \boldsymbol{X} \delta \boldsymbol{X}^{\mathrm{T}}\right)+\boldsymbol{\Gamma}_{k-1} \boldsymbol{W}_{k-1}\right] \\ & =\boldsymbol{f}\left(\boldsymbol{X}_{k-1}^{n}\right)+\boldsymbol{\Phi}_{k / k-1}^{n} \mathrm{E}[\delta \boldsymbol{X}]+\frac{1}{2} \sum_{i=1}^{n} \boldsymbol{e}_{i} \operatorname{tr}\left(\boldsymbol{D}_{f i} \cdot \mathrm{E}\left[\delta \boldsymbol{X} \delta \boldsymbol{X}^{\mathrm{T}}\right]\right)+\boldsymbol{\Gamma}_{k-1} \mathrm{E}\left[\boldsymbol{W}_{k-1}\right] \\ & =\boldsymbol{f}\left(\boldsymbol{X}_{k-1}^{n}\right)+\frac{1}{2} \sum_{i=1}^{n} \boldsymbol{e}_{i} \operatorname{tr}\left(\boldsymbol{D}_{f i} \boldsymbol{P}_{k-1}\right) \end{aligned} X^k/k1=E[f(Xk1n)+Φk/k1nδX+21i=1neitr(DfiδXδXT)+Γk1Wk1]=f(Xk1n)+Φk/k1nE[δX]+21i=1neitr(DfiE[δXδXT])+Γk1E[Wk1]=f(Xk1n)+21i=1neitr(DfiPk1)
状态预测既和状态的传递有关系,也和方差 P P P 阵有关系,比较麻烦。

状态预测的误差
X ~ k / k − 1 = X k − X ^ k / k − 1 = Φ k / k − 1 n δ X + 1 2 ∑ i = 1 n e i tr ⁡ ( D f i ⋅ δ X δ X T ) + Γ k − 1 W k − 1 − 1 2 ∑ i = 1 n e i tr ⁡ ( D f i P k − 1 ) = Φ k / k − 1 n δ X + 1 2 ∑ i = 1 n e i tr ⁡ ( D f i ( δ X δ X T − P k − 1 ) ) + Γ k − 1 W k − 1 \begin{aligned} \tilde{\boldsymbol{X}}_{k / k-1} & =\boldsymbol{X}_{k}-\hat{\boldsymbol{X}}_{k / k-1} \\ & =\boldsymbol{\Phi}_{k / k-1}^{n} \delta \boldsymbol{X}+\frac{1}{2} \sum_{i=1}^{n} \boldsymbol{e}_{i} \operatorname{tr}\left(\boldsymbol{D}_{f i} \cdot \delta \boldsymbol{X} \delta \boldsymbol{X}^{\mathrm{T}}\right)+\boldsymbol{\Gamma}_{k-1} \boldsymbol{W}_{k-1}-\frac{1}{2} \sum_{i=1}^{n} \boldsymbol{e}_{i} \operatorname{tr}\left(\boldsymbol{D}_{f i} \boldsymbol{P}_{k-1}\right) \\ & =\boldsymbol{\Phi}_{k / k-1}^{n} \delta \boldsymbol{X}+\frac{1}{2} \sum_{i=1}^{n} \boldsymbol{e}_{i} \operatorname{tr}\left(\boldsymbol{D}_{f i}\left(\delta \boldsymbol{X} \delta \boldsymbol{X}^{\mathrm{T}}-\boldsymbol{P}_{k-1}\right)\right)+\boldsymbol{\Gamma}_{k-1} \boldsymbol{W}_{k-1} \end{aligned} X~k/k1=XkX^k/k1=Φk/k1nδX+21i=1neitr(DfiδXδXT)+Γk1Wk121i=1neitr(DfiPk1)=Φk/k1nδX+21i=1neitr(Dfi(δXδXTPk1))+Γk1Wk1
状态预测的方差:比较复杂,最后算下来和四阶矩有关系
P k / k − 1 = E [ X ~ k / k − 1 X ~ k / k − 1 T ] = E [ [ Φ k ∣ k − 1 n δ X + 1 2 ∑ i = 1 n e i tr ⁡ ( D f i ( δ X δ X T − P k − 1 ) ) + Γ k − 1 W k − 1 ] × [ Φ k l k − 1 n δ X + 1 2 ∑ i = 1 n e i tr ⁡ ( D f i ( δ X δ X T − P k − 1 ) ) + Γ k − 1 W k − 1 ] T ] = Φ k k / k − 1 n P k − 1 ( Φ k / k − 1 n ) T + 1 4 E [ ∑ i = 1 n e i tr ⁡ ( D f i ( δ X δ X T − P k − 1 ) ) ( ∑ i = 1 n e i tr ⁡ ( D f i ( δ X δ X T − P k − 1 ) ) ) T ] + Γ k − 1 Q k − 1 Γ k − 1 T = Φ k / k − 1 n P k − 1 ( Φ k / k − 1 n ) T + 1 4 ∑ i = 1 n ∑ j = 1 n e i e j T E [ tr ⁡ ( D f i ( δ X δ X T − P k − 1 ) ) ⋅ tr ⁡ ( D j j ( δ X δ X T − P k − 1 ) ) ] + Γ k − 1 Q k − 1 Γ k − 1 T = Φ k / k − 1 n P k − 1 ( Φ k / k − 1 n ) T + 1 2 ∑ i = 1 n ∑ j = 1 n e i e j T tr ⁡ ( D f i P k − 1 D f j P k − 1 ) + Γ k − 1 Q k − 1 Γ k − 1 T \begin{array}{l} \boldsymbol{P}_{k / k-1}=\mathrm{E}\left[\tilde{\boldsymbol{X}}_{k / k-1} \tilde{\boldsymbol{X}}_{k / k-1}^{\mathrm{T}}\right] \\ =\mathrm{E}\left[\left[\boldsymbol{\Phi}_{k \mid k-1}^{n} \delta \boldsymbol{X}+\frac{1}{2} \sum_{i=1}^{n} \boldsymbol{e}_{i} \operatorname{tr}\left(\boldsymbol{D}_{f i}\left(\delta \boldsymbol{X} \delta \boldsymbol{X}^{\mathrm{T}}-\boldsymbol{P}_{k-1}\right)\right)+\boldsymbol{\Gamma}_{k-1} \boldsymbol{W}_{k-1}\right]\right. \\ \left.\times\left[\boldsymbol{\Phi}_{k l k-1}^{n} \delta \boldsymbol{X}+\frac{1}{2} \sum_{i=1}^{n} \boldsymbol{e}_{i} \operatorname{tr}\left(\boldsymbol{D}_{f i}\left(\delta \boldsymbol{X} \delta \boldsymbol{X}^{\mathrm{T}}-\boldsymbol{P}_{k-1}\right)\right)+\boldsymbol{\Gamma}_{k-1} \boldsymbol{W}_{k-1}\right]^{\mathrm{T}}\right] \\ =\boldsymbol{\Phi}_{k k / k-1}^{n} \boldsymbol{P}_{k-1}\left(\boldsymbol{\Phi}_{k / k-1}^{n}\right)^{\mathrm{T}} \\ +\frac{1}{4} \mathrm{E}\left[\sum_{i=1}^{n} \boldsymbol{e}_{i} \operatorname{tr}\left(\boldsymbol{D}_{f i}\left(\delta \boldsymbol{X} \delta \boldsymbol{X}^{\mathrm{T}}-\boldsymbol{P}_{k-1}\right)\right)\left(\sum_{i=1}^{n} \boldsymbol{e}_{i} \operatorname{tr}\left(\boldsymbol{D}_{f_{i}}\left(\delta \boldsymbol{X} \delta \boldsymbol{X}^{\mathrm{T}}-\boldsymbol{P}_{k-1}\right)\right)\right)^{\mathrm{T}}\right]+\boldsymbol{\Gamma}_{k-1} \boldsymbol{Q}_{k-1} \boldsymbol{\Gamma}_{k-1}^{\mathrm{T}} \\ =\boldsymbol{\Phi}_{k / k-1}^{n} \boldsymbol{P}_{k-1}\left(\boldsymbol{\Phi}_{k / k-1}^{n}\right)^{\mathrm{T}} \\ +\frac{1}{4} \sum_{i=1}^{n} \sum_{j=1}^{n} \boldsymbol{e}_{i} \boldsymbol{e}_{j}^{\mathrm{T}} \mathrm{E}\left[\operatorname{tr}\left(\boldsymbol{D}_{f i}\left(\delta \boldsymbol{X} \delta \boldsymbol{X}^{\mathrm{T}}-\boldsymbol{P}_{k-1}\right)\right) \cdot \operatorname{tr}\left(\boldsymbol{D}_{j j}\left(\delta \boldsymbol{X} \delta \boldsymbol{X}^{\mathrm{T}}-\boldsymbol{P}_{k-1}\right)\right)\right]+\boldsymbol{\Gamma}_{k-1} \boldsymbol{Q}_{k-1} \boldsymbol{\Gamma}_{k-1}^{\mathrm{T}} \\ =\boldsymbol{\Phi}_{k / k-1}^{n} \boldsymbol{P}_{k-1}\left(\boldsymbol{\Phi}_{k / k-1}^{n}\right)^{\mathrm{T}}+\frac{1}{2} \sum_{i=1}^{n} \sum_{j=1}^{n} \boldsymbol{e}_{\boldsymbol{i}} \boldsymbol{e}_{j}^{\mathrm{T}} \operatorname{tr}\left(\boldsymbol{D}_{f i} \boldsymbol{P}_{k-1} \boldsymbol{D}_{f j} \boldsymbol{P}_{k-1}\right)+\boldsymbol{\Gamma}_{k-1} \boldsymbol{Q}_{k-1} \boldsymbol{\Gamma}_{k-1}^{\mathrm{T}} \\ \end{array} Pk/k1=E[X~k/k1X~k/k1T]=E[[Φkk1nδX+21i=1neitr(Dfi(δXδXTPk1))+Γk1Wk1]×[Φklk1nδX+21i=1neitr(Dfi(δXδXTPk1))+Γk1Wk1]T]=Φkk/k1nPk1(Φk/k1n)T+41E[i=1neitr(Dfi(δXδXTPk1))(i=1neitr(Dfi(δXδXTPk1)))T]+Γk1Qk1Γk1T=Φk/k1nPk1(Φk/k1n)T+41i=1nj=1neiejTE[tr(Dfi(δXδXTPk1))tr(Djj(δXδXTPk1))]+Γk1Qk1Γk1T=Φk/k1nPk1(Φk/k1n)T+21i=1nj=1neiejTtr(DfiPk1DfjPk1)+Γk1Qk1Γk1T
量测方程的二阶展开
Z k / k − 1 = h ( X k / k − 1 n ) + H k n δ X + 1 2 ∑ i = 1 m e i i r ( D h i ⋅ δ X δ X T ) + V k Z_{k / k-1}=h\left(X_{k / k-1}^{n}\right)+\boldsymbol{H}_{k}^{n} \delta \boldsymbol{X}+\frac{1}{2} \sum_{i=1}^{m} e_{i} i r\left(D_{h i} \cdot \delta \boldsymbol{X} \delta \boldsymbol{X}^{\mathrm{T}}\right)+V_{k} Zk/k1=h(Xk/k1n)+HknδX+21i=1meiir(DhiδXδXT)+Vk
量测预测
Z ^ k / k − 1 = h ( X k / k − 1 n ) + 1 2 ∑ i = 1 m e i tr ⁡ ( D h i P k / k − 1 ) \hat{\boldsymbol{Z}}_{k / k-1}=\boldsymbol{h}\left(\boldsymbol{X}_{k / k-1}^{n}\right)+\frac{1}{2} \sum_{i=1}^{m} e_{i} \operatorname{tr}\left(\boldsymbol{D}_{h i} \boldsymbol{P}_{k / k-1}\right) Z^k/k1=h(Xk/k1n)+21i=1meitr(DhiPk/k1)
量测预测方差阵
P z z , k k − 1 = H k n P k / k − 1 ( H k n ) T + 1 2 ∑ i = 1 m ∑ j = 1 m e i e e j t ( D h i P k / k − 1 D h j P k / k − 1 ) + R k \boldsymbol{P}_{z z, k k-1}=\boldsymbol{H}_{k}^{n} \boldsymbol{P}_{k / k-1}\left(\boldsymbol{H}_{k}^{n}\right)^{\mathrm{T}}+\frac{1}{2} \sum_{i=1}^{m} \sum_{j=1}^{m} \boldsymbol{e}_{i}^{\mathrm{e}} \mathrm{e}_{j}^{\mathrm{t}}\left(\boldsymbol{D}_{h i} \boldsymbol{P}_{k / k-1} \boldsymbol{D}_{hj} \boldsymbol{P}_{k / k-1}\right)+\boldsymbol{R}_{k} Pzz,kk1=HknPk/k1(Hkn)T+21i=1mj=1meieejt(DhiPk/k1DhjPk/k1)+Rk
状态/量测互协方差阵
P x z , k k − 1 = P k / k − 1 H k n \boldsymbol{P}_{x z, k k-1}=\boldsymbol{P}_{k / k-1} \boldsymbol{H}_{k}^{n} Pxz,kk1=Pk/k1Hkn
最后,得二阶滤波公式:红色部分是二阶多出来的,运算量巨大收益微小

在这里插入图片描述

四、迭代EKF滤波

{ X ^ k / k − 1 = f ( X ^ k − 1 ) P k / k − 1 = Φ k / k − 1 P k − 1 Φ k / k − 1 T + Γ k − 1 Q k − 1 Γ k − 1 T Φ k / k − 1 = J ( f ( X ^ k − 1 ) ) K k = P k / k − 1 H k T ( H k P k / k − 1 H k T + R k ) − 1 H k = J ( h ( X ^ k / k − 1 ) ) X ^ k = X ^ k / k − 1 + K k [ Z k − h ( X ^ k / k − 1 ) ] P k = ( I − K k H k ) P k / k − 1 \left\{\begin{array}{ll}\hat{\boldsymbol{X}}_{k / k-1}=\boldsymbol{f}\left(\hat{\boldsymbol{X}}_{k-1}\right) & \\ \boldsymbol{P}_{k / k-1}={\color{red}\boldsymbol{\Phi}_{k / k-1}} \boldsymbol{P}_{k-1} {\color{red}\boldsymbol{\Phi}_{k / k-1}^{\mathrm{T}}}+\boldsymbol{\Gamma}_{k-1} \boldsymbol{Q}_{k-1} \boldsymbol{\Gamma}_{k-1}^{\mathrm{T}} & {\color{red}\boldsymbol{\Phi}_{k / k-1}=\boldsymbol{J}\left(\boldsymbol{f}\left(\hat{\boldsymbol{X}}_{k-1}\right)\right)} \\ \boldsymbol{K}_{k}=\boldsymbol{P}_{k / k-1} {\color{green}\boldsymbol{H}_{k}^{\mathrm{T}}}\left({\color{green}\boldsymbol{H}_{k}} \boldsymbol{P}_{k / k-1} {\color{green}\boldsymbol{H}_{k}^{\mathrm{T}}}+\boldsymbol{R}_{k}\right)^{-1} & {\color{green}\boldsymbol{H}_{k}=\boldsymbol{J}\left(\boldsymbol{h}\left(\hat{\boldsymbol{X}}_{k / k-1}\right)\right)} \\ \hat{\boldsymbol{X}}_{k}=\hat{\boldsymbol{X}}_{k / k-1}+\boldsymbol{K}_{k}\left[\boldsymbol{Z}_{k}-\boldsymbol{h}\left(\hat{\boldsymbol{X}}_{k / k-1}\right)\right] & \\ \boldsymbol{P}_{k}=\left(\boldsymbol{I}-\boldsymbol{K}_{k} {\color{green}\boldsymbol{H}_{k}}\right) \boldsymbol{P}_{k / k-1} & \end{array}\right. X^k/k1=f(X^k1)Pk/k1=Φk/k1Pk1Φk/k1T+Γk1Qk1Γk1TKk=Pk/k1HkT(HkPk/k1HkT+Rk)1X^k=X^k/k1+Kk[Zkh(X^k/k1)]Pk=(IKkHk)Pk/k1Φk/k1=J(f(X^k1))Hk=J(h(X^k/k1))

非线性系统展开点越接近你想研究的那个点,展开的精度越高。普通EKF的展开点是上一时刻 k − 1 k-1 k1,考虑迭代修正展开点。做完一次滤波之后,有 k k k 时刻的新量测值,可以用新量测值对 k − 1 k-1 k1 时刻做反向平滑,用反向平滑的值作为新的上一时刻的状态估计值,进行Kalman滤波。

在这里插入图片描述

预滤波
{ X ^ k / k − 1 ∗ = f ( X ^ k − 1 ) P k / k − 1 ∗ = Φ k / k − 1 ∗ P k − 1 ( Φ k / k − 1 ∗ ) T + Γ k − 1 Q k − 1 Γ k − 1 T Φ k / k − 1 ∗ = J ( f ( X ^ k − 1 ) ) K k ∗ = P k / k − 1 ∗ ( H k ∗ ) T [ H k ∗ P k / k − 1 ∗ ( H k ∗ ) T + R k ] − 1 H k ∗ = J ( h ( X ^ k / k − 1 ∗ ) ) X ^ k ∗ = X ^ k / k − 1 ∗ + K k ∗ [ Z k − h ( X ^ k / k − 1 ∗ ) ] \left\{\begin{array}{ll} \hat{\boldsymbol{X}}_{k / k-1}^{*}=\boldsymbol{f}\left(\hat{\boldsymbol{X}}_{k-1}\right) \\ \boldsymbol{P}_{k / k-1}^{*}=\boldsymbol{\Phi}_{k / k-1}^{*} \boldsymbol{P}_{k-1}\left(\boldsymbol{\Phi}_{k / k-1}^{*}\right)^{\mathrm{T}}+\boldsymbol{\Gamma}_{k-1} \boldsymbol{Q}_{k-1} \boldsymbol{\Gamma}_{k-1}^{\mathrm{T}} & \boldsymbol{\Phi}_{k / k-1}^{*}=\boldsymbol{J}\left(\boldsymbol{f}\left(\hat{\boldsymbol{X}}_{k-1}\right)\right) \\ \boldsymbol{K}_{k}^{*}=\boldsymbol{P}_{k / k-1}^{*}\left(\boldsymbol{H}_{k}^{*}\right)^{\mathrm{T}}\left[\boldsymbol{H}_{k}^{*} \boldsymbol{P}_{k / k-1}^{*}\left(\boldsymbol{H}_{k}^{*}\right)^{\mathrm{T}}+\boldsymbol{R}_{k}\right]^{-1} & \boldsymbol{H}_{k}^{*}=\boldsymbol{J}\left(\boldsymbol{h}\left(\hat{\boldsymbol{X}}_{k / k-1}^{*}\right)\right) \\ \hat{\boldsymbol{X}}_{k}^{*}=\hat{\boldsymbol{X}}_{k / k-1}^{*}+\boldsymbol{K}_{k}^{*}\left[\boldsymbol{Z}_{k}-\boldsymbol{h}\left(\hat{\boldsymbol{X}}_{k / k-1}^{*}\right)\right] & \end{array}\right. X^k/k1=f(X^k1)Pk/k1=Φk/k1Pk1(Φk/k1)T+Γk1Qk1Γk1TKk=Pk/k1(Hk)T[HkPk/k1(Hk)T+Rk]1X^k=X^k/k1+Kk[Zkh(X^k/k1)]Φk/k1=J(f(X^k1))Hk=J(h(X^k/k1))
反向一步平滑:可以用RTS算法
X ^ k − 1 / k = X ^ k − 1 + P k − 1 ( Φ k / k − 1 ∗ ) T ( P k / k − 1 ∗ ) − 1 ( X ^ k ∗ − X ^ k / k − 1 ∗ ) \hat{\boldsymbol{X}}_{k-1 / k}=\hat{\boldsymbol{X}}_{k-1}+\boldsymbol{P}_{k-1}\left(\boldsymbol{\Phi}_{k / k-1}^{*}\right)^{\mathrm{T}}\left(\boldsymbol{P}_{k / k-1}^{*}\right)^{-1}\left(\hat{\boldsymbol{X}}_{k}^{*}-\hat{\boldsymbol{X}}_{k / k-1}^{*}\right) X^k1/k=X^k1+Pk1(Φk/k1)T(Pk/k1)1(X^kX^k/k1)
重做状态一步预测 X ^ k − 1 / k \hat{\boldsymbol{X}}_{k-1 / k} X^k1/k 再做泰勒级数展开
X ^ k / k − 1 = f ( X ^ k − 1 ) = f ( X ^ k − 1 / k ) + f ( X ^ k − 1 ) − f ( X ^ k − 1 / k ) ≈ f ( X ^ k − 1 / k ) + J ( f ( X ^ k − 1 / k ) ) ( X ^ k − 1 − X ^ k − 1 / k ) \begin{aligned} \hat{\boldsymbol{X}}_{k / k-1} & =\boldsymbol{f}\left(\hat{\boldsymbol{X}}_{k-1}\right)=\boldsymbol{f}\left(\hat{\boldsymbol{X}}_{k-1 / k}\right)+\boldsymbol{f}\left(\hat{\boldsymbol{X}}_{k-1}\right)-\boldsymbol{f}\left(\hat{\boldsymbol{X}}_{k-1 / k}\right) \\ & \approx \boldsymbol{f}\left(\hat{\boldsymbol{X}}_{k-1 / k}\right)+\boldsymbol{J}\left(\boldsymbol{f}\left(\hat{\boldsymbol{X}}_{k-1 / k}\right)\right)\left(\hat{\boldsymbol{X}}_{k-1}-\hat{\boldsymbol{X}}_{k-1 / k}\right) \end{aligned} X^k/k1=f(X^k1)=f(X^k1/k)+f(X^k1)f(X^k1/k)f(X^k1/k)+J(f(X^k1/k))(X^k1X^k1/k)
重做量测一步预测
Z ^ k / k − 1 = h ( X ^ k / k − 1 ) = h ( X ^ k ∗ ) + h ( X ^ k / k − 1 ) − h ( X ^ k ∗ ) ≈ h ( X ^ k ∗ ) + J ( h ( X ^ k ∗ ) ) ( X ^ k / k − 1 − X ^ k ∗ ) \begin{aligned} \hat{\boldsymbol{Z}}_{k / k-1} & =\boldsymbol{h}\left(\hat{\boldsymbol{X}}_{k / k-1}\right)=\boldsymbol{h}\left(\hat{\boldsymbol{X}}_{k}^{*}\right)+\boldsymbol{h}\left(\hat{\boldsymbol{X}}_{k / k-1}\right)-\boldsymbol{h}\left(\hat{\boldsymbol{X}}_{k}^{*}\right) \\ & \approx \boldsymbol{h}\left(\hat{\boldsymbol{X}}_{k}^{*}\right)+\boldsymbol{J}\left(\boldsymbol{h}\left(\hat{\boldsymbol{X}}_{k}^{*}\right)\right)\left(\hat{\boldsymbol{X}}_{k / k-1}-\hat{\boldsymbol{X}}_{k}^{*}\right) \end{aligned} Z^k/k1=h(X^k/k1)=h(X^k)+h(X^k/k1)h(X^k)h(X^k)+J(h(X^k))(X^k/k1X^k)
迭代滤波:在此基础上做Kalman滤波
{ X ^ k / k − 1 = f ( X ^ k − 1 / k ) + Φ k / k − 1 ( X ^ k − 1 − X ^ k − 1 / k ) P k / k − 1 = Φ k / k − 1 P k − 1 Φ k / k − 1 T + Γ k − 1 Q k − 1 Γ k − 1 T Φ k / k − 1 = J ( f ( X ^ k − 1 / k ) ) K k = P k / k − 1 H k T ( H k P k / k − 1 H k T + R k ) − 1 H k = J ( h ( X ^ k ∗ ) ) X ^ k = X ^ k / k − 1 + K k [ Z k − h ( X ^ k ∗ ) − H k ( X ^ k / k − 1 − X ^ k ∗ ) ] P k = ( I − K k H k ) P k / k − 1 \left\{\begin{array}{ll} \hat{\boldsymbol{X}}_{k / k-1}=\boldsymbol{f}\left(\hat{\boldsymbol{X}}_{k-1 / k}\right)+\boldsymbol{\Phi}_{k / k-1}\left(\hat{\boldsymbol{X}}_{k-1}-\hat{\boldsymbol{X}}_{k-1 / k}\right) & \\ \boldsymbol{P}_{k / k-1}=\boldsymbol{\Phi}_{k / k-1} \boldsymbol{P}_{k-1} \boldsymbol{\Phi}_{k / k-1}^{\mathrm{T}}+\boldsymbol{\Gamma}_{k-1} \boldsymbol{Q}_{k-1} \boldsymbol{\Gamma}_{k-1}^{\mathrm{T}} & \boldsymbol{\Phi}_{k / k-1}=\boldsymbol{J}\left(\boldsymbol{f}\left(\hat{\boldsymbol{X}}_{k-1 / k}\right)\right) \\ \boldsymbol{K}_{k}=\boldsymbol{P}_{k / k-1} \boldsymbol{H}_{k}^{\mathrm{T}}\left(\boldsymbol{H}_{k} \boldsymbol{P}_{k / k-1} \boldsymbol{H}_{k}^{\mathrm{T}}+\boldsymbol{R}_{k}\right)^{-1} & \boldsymbol{H}_{k}=\boldsymbol{J}\left(\boldsymbol{h}\left(\hat{\boldsymbol{X}}_{k}^{*}\right)\right) \\ \hat{\boldsymbol{X}}_{k}=\hat{\boldsymbol{X}}_{k / k-1}+\boldsymbol{K}_{k}\left[\boldsymbol{Z}_{k}-\boldsymbol{h}\left(\hat{\boldsymbol{X}}_{k}^{*}\right)-\boldsymbol{H}_{k}\left(\hat{\boldsymbol{X}}_{k / k-1}-\hat{\boldsymbol{X}}_{k}^{*}\right)\right] \\ \boldsymbol{P}_{k}=\left(\boldsymbol{I}-\boldsymbol{K}_{k} \boldsymbol{H}_{k}\right) \boldsymbol{P}_{k / k-1} & \end{array}\right. X^k/k1=f(X^k1/k)+Φk/k1(X^k1X^k1/k)Pk/k1=Φk/k1Pk1Φk/k1T+Γk1Qk1Γk1TKk=Pk/k1HkT(HkPk/k1HkT+Rk)1X^k=X^k/k1+Kk[Zkh(X^k)Hk(X^k/k1X^k)]Pk=(IKkHk)Pk/k1Φk/k1=J(f(X^k1/k))Hk=J(h(X^k))

还可以多次迭代,但收益小

本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如若转载,请注明出处:http://www.coloradmin.cn/o/630176.html

如若内容造成侵权/违法违规/事实不符,请联系多彩编程网进行投诉反馈,一经查实,立即删除!

相关文章

大家都说Java有三种创建线程的方式,并发编程中的惊天骗局

在Java中,创建线程是一项非常重要的任务。线程是一种轻量级的子进程,可以并行执行,使得程序的执行效率得到提高。Java提供了多种方式来创建线程,但许多人都认为Java有三种创建线程的方式,它们分别是继承Thread类、实现…

论文浅尝 | Dually Distilling KGE for Faster and Cheaper Reasoning

笔记整理:张津瑞,天津大学硕士,研究方向为知识图谱 链接:https://dl.acm.org/doi/10.1145/3488560.3498437 动机 知识图谱已被证明可用于各种 AI 任务,如语义搜索,信息提取和问答等。然而众所周知&#xff…

【C++】C++11常用新特性

✍作者:阿润菜菜 📖专栏:C 目录 一、统一的列表初始化二、 简化声明2.1 auto2.2 decltype2.3 nullptr 三、右值引用和移动语义 -- 重要3.1 区分左值引用和右值引用3.2 对比左值引用看看右值引用使用价值3.3 万能引用和完美转发(st…

基于word文档,使用Python输出关键词和词频,并将关键词的词性也标注出来

点击上方“Python爬虫与数据挖掘”,进行关注 回复“书籍”即可获赠Python从入门到进阶共10本电子书 今 日 鸡 汤 移船相近邀相见,添酒回灯重开宴。 大家好,我是Python进阶者。 一、前言 前几天在有个粉丝问了个问题,大概意思是这样…

一道北大强基题背后的故事(三)——什么样的题是好题?

早点关注我,精彩不错过! 上回我们针对这道北大强基题[((1 sqrt(5)) / 2) ^ 12]在答案的基础上给出了出题的可能思路,想一探究竟,相关内容请戳: 一道北大强基题背后的故事(二)——出题者怎么想的…

【Kubernetes入门】Service四层代理入门实战详解

文章目录 一、Service四层代理概念、原理1、Service四层代理概念2、Service工作原理3、Service原理解读4、Service四种类型 二、Service四层代理三种类型案例1、创建ClusterIP类型Service2、创建NodePort类型Service3、创建ExternalName类型Service 三、拓展1、Service域名解析…

Latex中图片排版(多个子图、横排、竖排、添加小标题)

1、两个子图横排 \begin{figure}[!t] \centering %\includegraphics[width3in]{fig5} \subfloat[subfig figure title]{\includegraphics[scale0.5]{superd2}} \subfloat[subfig figure title]{\includegraphics[scale0.5]{superd2}} \caption{title} \label{fig_6} \end{figu…

阿里发布的百亿级高并发系统(全彩版小册),涵盖了所有的高并发操作

高并发 提到“高并发”相信你们应该都不会感到陌生!此时你脑中应该会浮现好多有关高并发的:业务急剧增长、电商购物、电商秒杀、12306抢票、淘宝天猫各种活动等;都是需要用到高并发的,那么如何去设计一个高并发系统抵挡这些冲击呢…

Django的app里面的视图函数

我之前说过需要重点去了解view和model,下面是我的总结。 视图函数是存在view.py里面的,视图函数的主要功能是接收请求、返回响应。在建立应用程序后,先在URL配置文件中加一条配置项指明URL与视图函数的对应关系。然后按照实际需求在视图函数…

三次握手四次挥手过程剖析

【一】预备知识: 1.三次握手并不一定非得成功,最担心得其实就是最后一个akc(应答)丢失,但是还是有配套得解决方案,比如超时重传机制。 2.连接是需要被保存下来得,是需要被os管理起来得&#xf…

被一个gpio口搞死的一天

今天是新项目调试的第一天。 我起的很早,起早的原因很简单,我家楠哥要我送他上学,他说爸爸没有起到一个当爸爸的责任,他也想让爸爸送他上学,然后我就送了。 7点30起来,8点出发,然后回来看了一下…

IDEA把css/js压缩成一行min文件,idea实现右键压缩css和js文件

前言 发布时有些css和js文件较长多行,导致加载的时候略慢,所以想把指定的css或js压缩 实现 整合 yuicompressor-2.4.8.jar 下载地址1:https://github.com/yui/yuicompressor/releases 下载地址2:https://github.com/yui/yuicom…

小学生开“卷”AIGC,绝不能输在起跑线上

图片来源:由无界AI生成 OpenAI的研究报告称,未来,大量工作岗位将受到AI冲击,首当其冲的岗位是作家、数学家、网页设计师、记者、律师…… 自从ChatGPT问世以来,人类会被AI替代的讨论甚嚣尘上,焦虑情绪无处不…

chatgpt赋能python:Python中的倒序遍历:如何使用Python倒序遍历?

Python中的倒序遍历:如何使用Python倒序遍历? Python是一种高级编程语言,它非常适合数据科学、机器学习和人工智能等领域。Python的强大之处在于它有很多内置功能,其中包括倒序遍历。 在本篇文章中,我们将介绍如何使…

EmbodiedGPT|具身智能或将成为实现AGI的最后一公里

卷友们好,我是穆尧。 最近由Chatgpt所引爆的新一代人工智能的革命正在如火如荼的进行,几乎重塑了所有的互联网产品,如办公软件、浏览器插件、搜索引擎、推荐系统等。这样巨大的改变,让大家对通用人工智能又燃起了新的希望&#xf…

CTPN文本检测详解 面试版本

二.关键idea 1.采用垂直anchor回归机制,检测小尺度的文本候选框 2.文本检测的难点在于文本的长度是不固定,可以是很长的文本,也可以是很短的文本.如果采用通用目标检测的方法,将会面临一个问题:**如何生成好…

Autosar诊断实战系列01-手把手教你增加一路31Routine服务

本文框架 1.系列概述2. UDS Routine服务添加3. DcmDspRoutine配置3.1 DcmDspRoutineInfos配置3.2 DcmDspRoutines配置1.系列概述 在本系列笔者将结合工作中对诊断实战部分的应用经验进一步介绍常用UDS服务的进一步探讨及开发中注意事项, Dem/Dcm/CanTp/Fim模块配置开发及注意…

编译tolua——2、基础编译tolua

目录 1、编译工具和环境说明 2、基础编译tolua 大家好,我是阿赵。 继续来讲tolua的各个常用平台的编译。 这里使用官方的tolua_runtime-master项目来做编译 具体需要的编译软件和源码地址,在上一篇文章已经介绍过了,先把环境准备好&#xff…

飞桨AI4S污染物扩散快速预测模型,亮相全国数据驱动计算力学研讨会

5月19-21日,第一届全国数据驱动计算力学研讨会在大连召开。本次研讨会由中国力学学会主办,大连理工大学运载工程与力学学部承办,北京理工大学先进结构技术研究院协办。 会议共吸引了400多位来自全国各地高校与企业的老师与学生参会&#xff0…

DNSPod十问林洪祥:顶级带货主播,其实是数字人?

本期嘉宾 林洪祥 风平智能CEO 林洪祥,风平智能CEO。风平智能拥有全球领先的数字人AIGC预训练大模型技术,利用数字人AI知识大模型打造视频版ChatGPT,实现数字人名师、数字人医生、数字人保险客服、数字广告模特、数字人AI直播等,…