网站导航

  • 首 页
  • 期刊介绍
    • 期刊简介
    • 收录情况
    • 荣       誉
  • 编委会
    • 成员介绍
    • 主编简介
  • 作者中心
    • 投稿须知
    • 写作模板
    • 保密协议
    • 版面费交纳须知
    • EI摘要写作要求
    • 中国分类号
  • 文件法规
    • 伦       理
    • 同行评审
    • 撤       稿
    • 数字出版
  • 审者中心
    • 审稿须知
    • 审稿单
  • 联系我们
  • English
高阶容积卡尔曼滤波及其在目标跟踪中的应用
Download PDF  
文章快速检索  
  高级检索
引用本文
张龙, 崔乃刚, 杨峰, 路菲, 卢宝刚. 高阶容积卡尔曼滤波及其在目标跟踪中的应用[J]. 哈尔滨工程大学学报, 2016, 37(04): 573-578
ZHANG Long, CUI Naigang, YANG Feng, LU Fei, LU Baogang. High-degree cubature Kalman filter and its application in target tracking[J]. Journal of Harbin Engineering University, 2016, 37(04): 573-578.

DOI:10.11990/jheu.201412079
高阶容积卡尔曼滤波及其在目标跟踪中的应用
张龙1, 崔乃刚1, 杨峰1, 路菲1, 卢宝刚2    
1. 哈尔滨工业大学 航天工程系, 黑龙江 哈尔滨 150001;
2. 北京航天长征飞行器研究所, 北京 100076    
基金项目: 国家自然科学基金项目(61304236).
作者简介: 张龙(1987-), 男,博士研究生;崔乃刚(1965-), 男,教授,博士生导师
收稿日期: 2014-12-31
通信作者: 张龙, E-mail: zhanglong_hit@163.com
摘要:针对传统的容积卡尔曼滤波(CKF)估计精度有限的问题,提出了一种基于任意阶容积规则的高阶容积卡尔曼滤波(HCKF)方法并应用于机动目标跟踪问题。传统的CKF采用三阶球面-相径容积规则,可获得优于其他非线性滤波如不敏卡尔曼滤波(UKF)的估计精度和数值稳定性。为了进一步提高CKF的估计精度,在基于点的高斯近似滤波框架下,分别使用Genz积分方法和矩匹配法推导出任意阶的球面规则和相径规则,以此构造高阶球面-相径容积规则来计算高斯型积分,并建立高阶容积卡尔曼滤波算法。将提出的HCKF算法应用于机动目标跟踪问题中并进行数值仿真。仿真结果表明,相对于传统容积卡尔曼滤波,高阶容积卡尔曼滤波对目标位置和速度估计的精度分别提高了11%和24%,可获得更高的估计精度。
关键词: 高阶容积卡尔曼滤波    目标跟踪    非线性系统    贝叶斯估计    球面相径规则    容积规则    
High-degree cubature Kalman filter and its application in target tracking
ZHANG Long1, CUI Naigang1, YANG Feng1, LU Fei1, LU Baogang2    
1. Department of Astronautics Engineering, Harbin Institute of Technology, Harbin 150001, China;
2. Beijing Institute of Space Long March Vehicle, Beijing 100191, China
Abstract: A high-order cubature Kalman filter (HCKF) based on the arbitrary-order cubature rule was proposed and applied to the maneuvering target tracking problem resulting from the limited precision of the conventional cubature Kalman filter (CKF). The conventional CKF, which employs the third-order spherical-radial cubature rule, can achieve better estimation precision and numerical stability than the other nonlinear filters, such as Unscented Kalman Filter(UKF). To further improve the estimation precision of CKF, in the framework of point-based Gaussian approximation filters, Genz's method and the moment-matching method were employed to deduce the arbitrary-order spherical rule and the radial rule, respectively. By combining the two rules, the high-order spherical-radial cubature rule was developed for computing Gaussian weighted integrals, and HCKF was proposed based on the high-order spherical-radial cubature rule. The proposed HCKF algorithm was applied to the maneuvering target tracking to test its performance. The simulation result shows that HCKF can achieve better estimation precision than conventional CKF, with the improvement of position and velocity estimates by 11% and 24% respectively.
Key words: high-degree cubature Kalman filter    target tracking    nonlinear system    Bayesian filtering    spherical-radial rule    cubature rule    

近年,基于贝叶斯框架的非线性滤波和估计得到了深入研究,特别是基于高斯概率密度函数假设下的高斯近似滤波[1]。其中,扩展卡尔曼滤波(EKF)是目前应用最广泛的滤波方法,其通过对当前时刻状态估计的局部线性化来实现一阶近似,但是当系统非线性程度较高以及噪声较大时,EKF的估计性能随之下降[2]。为了提高EKF的估计性能,一类基于点的高斯近似滤波得到越来越多的关注,例如不敏卡尔曼滤波(UKF)、高斯-厄米特积分滤波(GHQF)以及容积卡尔曼滤波(CKF)等,相对于EKF,它们无需计算微分且可以获得更高阶的估计精度。在贝叶斯框架下,这些滤波方法仅在计算高斯权值积分的方式上不同,即采用的数值积分规则不同。UKF通过使用一组Sigma点经非线性函数传递后的加权和来近似状态的概率密度函数,其估计精度可以达到三阶[3, 4, 5]。但由于UKF的中心Sigma点权值为负值,容易导致滤波数值特性不稳定甚至发散。GHQF采用高斯-厄米特积分规则来近似高斯积分,可以获得任意阶的估计精度[6, 7]。但是GHQF由于积分点数和计算复杂度随着维数指数增长,容易发生维数爆炸问题,很难被应用到实际问题中[7]。针对高维系统的状态估计问题,Arasaratnam等基于三阶球面-相径容积规则提出了CKF[8, 9, 10],并得到广泛应用[11, 12, 13, 14]。因为CKF使用的所有点权值相同且为正数,因此其值稳定性优于UKF。但是CKF估计精度仍然有限且算法存在一些缺陷。例如,CKF无法准确计算一些简单多项式函数的高斯权值积分,如x12x22,其中x1和x2是高斯随机向量函数的两个元素。为了克服CKF的缺点,需要更高阶的CKF来取得更高精度的估计结果。文献[15, 16, 17]基于Genz积分方法和矩匹配法分别推导了任意阶球面规则和相径规则,建立了可获取任意阶精度的高阶容积规则并应用于二维空间目标跟踪问题中,得到了优于传统CKF的估计精度。

本文针对CKF估计精度有限的问题提出了高阶容积卡尔曼滤波,推导了五阶球面-相径容积规则来近似高斯权值积分,在五阶容积规则的基础上建立高阶容积卡尔曼滤波算法模型和计算流程,并应用于三维空间中的机动目标跟踪系统。

1 高阶球面-相径容积规则 1.1 高斯权值积分近似方法

考虑如下非线性离散时间动力学系统:

${x_k} = f({x_{k - 1}}) + {v_{k - 1}}$ (1)
${y_k} = h({x_k}) + {w_k}$ (2)
式中:状态向量xk∈Rn,测量向量yk∈Rm;vk-1和wk是相互独立的零均值系统高斯白噪声和测量高斯白噪声,方差分别为Qk-1和Rk。

对于高斯假设下的非线性动力学系统,将贝叶斯估计基本理论与任意阶容积规则相结合,可推导出任意阶容积卡尔曼滤波方法。任意阶CKF与其他一般的高斯近似滤波方法结构相同,但使用任意阶球面-相径容积规则来计算高斯权值积分式I=∫Rng(x)N(x;0,I)dx的近似值。 其中,g(x)为一般非线性函数,N(x;0,I)表示零均值、协方差为I的高斯分布。

在容积规则中,考虑如下积分[9]:

$I(g) = \int_{{{\text{R}}^n}} {g(x)\exp ( - x{\text{T}}x){\text{dx}}} $ (3)

令x=rs,则式(3)可转换到球面-相径的坐标系统中,有

$I(g) = \int_0^\infty {\int_{{U_n}} {g(rs){r^{n - 1}}\exp ( - {r^2}} } ){\text{d}}\sigma (s)dr$ (4)
式中:方向向量s=[s1 s2 … sn]T,n维球面Un满足Un={s∈Rn:s12+s22+…sn2=1},σ(·)为该球面上的积分微元。

容易看出,式(4)包含两种积分:球面积分$\int_{{U_n}} {{g_s}} d\sigma (s)$和相径积分$\int_0^\infty {{g_r}{{(r)}^{n - 1}}\exp ( - {r^2})dr} $,权重函数分别是${w_g}(s) = 1$和${w_g}(r) = {(r)^{n - 1}}\exp ( - {r^2})$。容积卡尔曼滤波的球面-相径容积规则就是基于这两种积分的结合[8, 9, 10],则式(4)可近似表达为

$\sum\limits_{i = 1}^{{N_r}} {\sum\limits_{j = 1}^{{N_s}} {{w_r},{w_{s,j}}g({r_i}{s_j})} } $ (5)
式中:ri和wr,i为计算相径积分的点和对应的权值,sj和ws,j为计算球面积分的点和对应的权值,Nr和Ns分别是相径积分和球面积分的点数。则高斯权值积分式可转换为

$\int_{{{\text{R}}^n}} {g(x)N(x;0,I} )dx \approx \frac{1}{{{\pi ^{n/2}}}}\sum\limits_{i = 1}^{{N_r}} {\sum\limits_{j = 1}^{{N_s}} {{w_{r,i}}} } {w_{s,j}}g(\sqrt 2 {r_i}{s_j})$ (6)

采用不同阶数的容积规则,并根据该规则选取不同容积点和权值,式(6)可得到不同的积分结果。目前,三阶球面-相径容积规则已经得到广泛的关注,下面给出五阶球面-相径容积规则。

1.2 五阶球面-相径容积规则

文献[15, 16]对球面规则进行了详细推导。由Genz积分方法[18, 19],五阶球面规则满足:

$\begin{gathered} {I_{{U_n}}},5({g_s}) = {\varpi _{s1}}\sum\limits_{j = 1}^{n(n - 1)/2} {(({g_s}(s_j^ + } ) + gs( - s_j^ + ) + {g_s}(s_j^ - ) + \hfill \\ {g_s}(s_j^ - )) + {\varpi _{s2}}\sum\limits_{j = 1}^n {({g_s}(} {e_j}) + {g_s}( - {e_j})) \hfill \\ \end{gathered} $ (7)
式中:ej为空间Rn的单位矢量,其第j个元素为1。点集sj+和sj-分别给出:
$\{ s_j^ - \} = \{ \sqrt {1/2} ({e_k} + {e_l}):k < l,k,l = 1,2, \cdot \cdot \cdot ,n\} $ (8)
$\{ s_j^ + \} = \{ \sqrt {1/2} ({e_k} + {e_l}):k < l,k,l = 1,2, \cdot \cdot \cdot ,n\} $ (9)
式(7)中的权值$\varpi $s1和$\varpi $s2分别为
${\varpi _{s1}} = {A_n}/n(n + 2)$ (10)
${\varpi _{s2}} = (4 - n){A_n}/n(n + 2)$ (11)
式中:${A_n} = 2\sqrt {{\pi ^n}} /\Gamma (n/2)$为单位球面的表面积,且函数Γz定义为$\Gamma (z) = \int_0^\infty {\exp ( - \lambda ){\lambda ^{z - 1}}d\lambda } $。 由矩匹配法(moment matching method),对于Nr=2的五阶相径规则,积分点和权值需要满足:

$\left\{ {\begin{array}{*{20}{c}} {{w_{r,1}}r_1^0 + {w_{r,2}}r_2^0 = \Gamma (n/2)/2} \\ {{w_{r,1}}r_1^2 + {w_{r,2}}r_2^2 = n\Gamma (n/2)/4} \\ {{w_{r,1}}r_1^4 + {w_{r,2}}r_2^4 = (n/2)(n/2) + 1)\Gamma (n/2)/2} \end{array}} \right.$ (12)

选择r1为自由变量并令r1=0,解式(12),可得五阶相径规则的点和权值为

$\left\{ {\begin{array}{*{20}{c}} {{r_1} = 0} \\ {{r_2} = \sqrt {n/2 + 1} } \end{array}} \right.$ (13)
$\left\{ {\begin{array}{*{20}{c}} {{w_{r,1}} = n\Gamma (n/2)(n + 2)} \\ {{w_{r,2}} = n\Gamma (n/2)(2(n + 2))} \end{array}} \right.$ (14)

结合式(6)、(7)、(13)和(14),可以得到满足分布为x,~N(x;0,I)的五阶球面-相径容积规则Nr=2,Ns=2n2:

$\begin{gathered} \int_{{R^n}} {g(x)N(x;0,I){\text{d}}x \approx \frac{1}{{{\pi ^{n/2}}}}} \sum\limits_{i = 1}^{{N_r}} {\sum\limits_{j = 1}^{{N_s}} {{w_{r,i}}{w_{s,j}}g(2{r_i}{s_j}) = } } \hfill \\ \frac{2}{{n + 2}}g(0) + \frac{1}{{{{(n + 2)}^2}}}\sum\limits_{j = 1}^{n(n - 1)/2} {[g(\beta \cdot s_j^ + } ) + g( - \beta \cdot s_j^ + )] + \hfill \\ \frac{1}{{{{(n + 2)}^2}}}\sum\limits_{j = 1}^{n(n - 1)/2} {[g(\beta \cdot s_j^ - } ) + g( - \beta \cdot s_j^ - )] + \hfill \\ \frac{{4 - n}}{{2{{(n + 2)}^2}}}\sum\limits_{j = 1}^n {[g(\beta \cdot {e_j}} ) + g( - \beta \cdot {e_j})] \hfill \\ \end{gathered} $ (15)
式中$\beta = \sqrt {n + 2} $。对于j=1,2,…,2n(n-1),ws,j=$\varpi $s1;对于j=2n(n-1)+1,2n(n-1)+2,…,2n2,ws,j=$\varpi $s2。由式(15)可以看出,五阶球面-相径容积规则使用的容积点个数为2n2+1。

对于更一般的高斯分布N(x;${\hat x}$,P),五阶球面-相径容积规则如下:

$\begin{gathered} \int_{{R^{_n}}} {g(x)N(x;\hat x,P){\text{d}}x} = \int_{{R^{_n}}} {g(sx + \hat x)N(x;0,I){\text{d}}x} \approx \hfill \\ \frac{1}{{{{(n + 2)}^2}}}\sum\limits_{j = 1}^{n(n - 1)/2} {[g(\beta S \cdot s_j^ + } + \hat x) + \hfill \\ [g( - \beta S \cdot s_j^ + + \hat x)] + \hfill \\ \frac{1}{{{{(n + 2)}^2}}}\sum\limits_{j = 1}^{n(n - 1)/2} {[g(\beta S \cdot s_j^ - } + \hat x) + \hfill \\ [g( - \beta S \cdot s_j^ - + \hat x)] + \hfill \\ \frac{{4 - n}}{{2{{(n + 2)}^2}}}\sum\limits_{j = 1}^n {[g(\beta S \cdot {e_j}} + \hat x) + \hfill \\ g( - \beta S \cdot {e_j} + \hat x)] + \frac{{2g(\hat x)}}{{n + 2}} \hfill \\ \end{gathered} $ (16)

2 高阶容积卡尔曼滤波

根据式(8)、(9)、(13)、(14)、(16)可得到高阶容积卡尔曼滤波HCKF(五阶)的具体实现步骤。HCKF算法步骤同三阶CKF算法基本相同,只是在积分点和权值选取上有区别,HCKF采用2n2+1个容积点。针对非线性系统(1)~(2),假设k时刻的状态向量x满足xk~N(xk;${\hat x}$k,Pk),则HCKF具体算法如下:

1)计算容积点xk,i(i=0,1,…,2n2)。

xk,i=${\hat x}$k+Skξi (17)
式中Sk是Pk的Cholesky分解,即Pk=SkSkT 。向量ξi为
${\xi _i} = \left\{ {\begin{array}{*{20}{l}} {{{[00 \cdot \cdot \cdot 0]}^{\text{T}}},i = 0} \\ {\beta S_j^ + ,i = 1,2, \cdot \cdot \cdot n(n - 1)/2} \\ { - \beta S_{i - n(n - 1)/2}^ + ,} \\ {i = n(n - 1)/2 + 1,n(n - 1)/2 + 2, \cdot \cdot \cdot ,n(n - 1)} \\ {\beta S_{i - n(n - 1)}^ - ,} \\ {i = n(n - 1) + 1,n(n - 1) + 2, \cdot \cdot \cdot ,3n(n - 1)/2} \\ { - \beta S_{i - 3n(n - 1)/2}^ - ,} \\ {i = 3n(n - 1)/2 + 1,3n(n - 1)/2 + 2, \cdot \cdot \cdot ,2n(n - 1)} \\ {\beta e_{i - 2n(n - 1)}^{},} \\ {i = 2n(n - 1) + 1,2n(n - 1) + 2, \cdot \cdot \cdot ,n(2n - 1)} \\ { - \beta e_{i - n(2n - 1)}^{},} \\ {i = n(2n - 1) + 1,n(2n - 1) + 2, \cdot \cdot \cdot ,2{n^2}} \end{array}} \right.$ (18)
式中:$\beta = \sqrt {n + 2} $,ei表示n维单位向量且其第i个元素为1。sj+和sj-分别为
$\left\{ {\begin{array}{*{20}{c}} {s_j^ + = \sqrt {1/2} ({e_p} + {e_q}),p < q,pq = 1,2, \cdot \cdot \cdot ,n} \\ {s_j^ - = \sqrt {1/2} ({e_p} + {e_q}),p < q,pq = 1,2, \cdot \cdot \cdot ,n} \end{array}} \right.$ (19)

2)计算经状态方程传递后的容积点χk+1/k,i:

${x_{k + 1/k,i}} = f({x_{k,i}})$ (20)

3)计算k+1时刻的状态预测值${{\hat x}_{k + 1/k}}$:

${{\hat x}_{k + 1/k}} = \sum\limits_{i = 0}^{2n2} {{w_i}{x_{k + 1/k,i}}} $ (21)
其中权值wi分别为

${w_i} = \left\{ {\begin{array}{*{20}{c}} {2/n + 2,i = 0} \\ {1/{{(n + 2)}^2},i = 1,2 \cdot \cdot \cdot ,2n(n - 1)} \\ {(4 - n)/{{(n + 2)}^2},} \\ {i = 2n(n - 1) + 1,2n(n - 1) + 2, \cdot \cdot \cdot ,2{n^2}} \end{array}} \right.$ (22)

4)估计k+1时刻的状态误差协方差阵Pk+1/k。

${\sum\limits_{i = 0}^{2n2} {{w_i}({x_{k + 1/k,i}} - {{\hat x}_{k + 1/k}}){{({x_{k + 1/k,i}} - {{\hat x}_{k + 1/k}})}^T} + {Q_k}} }$ (23)

5)计算更新后的状态容积点xk+1/k,i: 其中,${P_{k + 1/k}} = {S_{k + 1/k}}S_{k + 1/k}^{\text{T}}$。

${x_{k + 1/k,i}} = {S_{k + 1/k}}{\xi _i} + {{\hat x}_{k + 1/k}},i = 0,1, \cdot \cdot \cdot ,2{n^2}$ (24)
其中,${{P_{k + 1/k}} = {S_{k + 1/k}}S_{k + 1/k}^{\text{T}}}$。

6)计算经过测量方程传递的容积点yk+1,i:

yk +1,i = h(xk +1/ k,i ) (25)

7)计算k+1时刻的测量预测值${\hat y}$k+1:

${{\hat y}_{k + 1}} = \sum\limits_{i = 0}^{2n2} {{w_i}{y_{k + 1,i}}} $ (26)

8)估计k+1时刻的测量误差协方差阵Pk+1yy和一步预测互相关协方差阵Pk+1/kxy:

$P_{k + 1}^{yy} = \sum\limits_{i = 0}^{2n2} {{w_i}({y_{k + 1,i}}} - {{\hat y}_{k + 1}}){({y_{k + 1,i}} - {{\hat y}_{k + 1}})^{\text{T}}} + {R_k}$ (27)
$P_{k + 1}^{yy} = \sum\limits_{i = 0}^{2n2} {{w_i}({x_{k + 1/k,i}}} - {{\hat x}_{k + 1/k}}){({y_{k + 1,i}} - {{\hat y}_{k + 1}})^{\text{T}}}$ (28)

9)计算k+1时刻的滤波增益矩阵Kk+1:

${K_{k + 1}} = P_{k + 1/k}^{xy}{(P_{k + 1}^{yy})^{ - 1}}$ (29)

10)计算k+1时刻的状态估计值${\hat x}$k+1:

${{\hat x}_{k + 1}} = {{\hat x}_{k + 1/k}} + {K_{k + 1}}({y_{k + 1}} - {{\hat y}_{k + 1}})$ (30)

11)估计k+1时刻的状态误差协方差阵Pk+1。

${P_{k + 1}} = {P_{k + 1/k}} - {K_{k + 1}}P_{k + 1}^{yy}K_{k + 1}^{\text{T}}$ (31)

给定初始条件${\hat x}$0|0,P0|0,按照以上计算流程即可进行高阶容积滤波。

3 数值仿真

为了表明高阶容积卡尔曼滤波算法的优越性能,将高阶容积滤波(五阶CKF)应用到三维空间中的目标跟踪系统中,并与不敏卡尔曼滤波(UKF)和传统容积滤波(三阶CKF)进行对比。

定义地面坐标系o-xyz的原点位于地面,ox指向东向,oy指向北向,oz与ox和oy成右手系,指向天向。假设目标在坐标系o-xyz中等高度机动飞行,目标跟踪系统动力学模型如下[9, 16, 17]: ${X_k} = f({X_{K - 1}}) + {v_{k - 1}}\left[ {\begin{array}{*{20}{l}} 1&{\frac{{\sin (\Omega T)}}{\Omega }}&0&{\frac{{ - \cos (\Omega T)}}{\Omega }}&0&0&0 \\ 0&{\cos (\Omega T)}&0&{ - \sin (\Omega T)}&0&0&0 \\ 0&{\frac{{1 - \cos (\Omega T)}}{\Omega }}&1&{\frac{{\sin (\Omega T)}}{\Omega }}&0&0&0 \\ 0&{\sin (\Omega T)}&0&{\cos (\Omega T)}&0&0&0 \\ 0&0&0&0&1&T&0 \\ 0&0&0&0&0&0&0 \\ 0&0&0&0&0&0&1 \end{array}} \right]{X_{k - 1}} + {v_{k - 1}}$ 式中:${X_k} = {\left[ {x\dot xy\dot yz\dot z\Omega } \right]^{\text{T}}}$为k时刻的目标状态,${\left[ {xyz} \right]^{\text{T}}}$和${\left[ {\dot x\dot y\dot z} \right]^{\text{T}}}$分别为目标的位置和速度矢量;Ω为目标转弯角速度且未知;T为连续两次测量的时间间隔;vk-1为零均值的系统高斯白噪声,其协方差为Qk-1。Q=diag(M1,M1,M1,M2) 式中:$M1 = \left[ {\begin{array}{*{20}{c}} {{q_1}{T^3}/3}&{{q_1}{T^2}/2} \\ {{q_1}{T^2}/2}&{{q_1}T} \end{array}} \right]$,M2=q2T,q1和q2是与系统噪声强度相关的标量。

假设测量雷达固定在地面坐标系原点,获取对目标的相对距离r、高低角η和方位角θ信息。假设目标的测量向量为y=[r θ η]T,满足:$\begin{gathered} {y_k} = \left[ {\begin{array}{*{20}{c}} {{r_k}} \\ {{\theta _k}} \\ {{\eta _k}} \end{array}} \right] = h({X_k}) + {w_k} = \hfill \\ \left[ {\begin{array}{*{20}{c}} {\sqrt {x_k^2 + y_k^2 + z_k^2} } \\ {a\tan ({y_k}/{x_k})} \\ {a\tan ({z_k}/\sqrt {x_k^2 + y_k^2} } \end{array}} \right] + {w_k} \hfill \\ \end{gathered} $ 式中:w为测量噪声,其协方差阵R为$R = {\text{diag}}(\sigma _r^2,\sigma _\theta ^2,\sigma _\beta ^2)$ 式中σr、σθ、σβ是雷达测量噪声的标准差。

仿真进行100次蒙特卡洛打靶实验,并令每次打靶的初始状态估计${\hat X}$0|0由满足均值为X0、协方差为P0|0的高斯正态分布N(${\hat X}$0;${\hat X}$0|0,P0|0)随机生成。同时,采用UKF、三阶CKF与五阶CKF三种滤波方法对目标跟踪并对位置均方根误差、速度均方根误差和转速均方根误差进行对比。定义位置均方根误差为$\sqrt {\frac{{\sum\limits_{i = 1}^N {({{(x_k^i - \hat x_k^i)}^2} + {{(y_k^i - \hat y_k^i)}^2} + {{(z_k^i - \hat z_k^i)}^2})\frac{1}{2}} }}{2}} $式中:(xki,yki,zki)和(${\hat x_k^i}$,${\hat y_k^i}$,${\hat z_k^i}$)分别是第i次打靶中第k时刻的目标真实位置分量和估计值,N为打靶次数。同理可得速度均方根误差RMSEvel和转速均方根误差RMSEomg。在每一次打靶中所有滤波方法的初始状态相同,打靶时间为100 s。

图1为目标飞行的平面轨迹和雷达位置,“I”表示轨迹起点,“F”表示轨迹终点,“★”为雷达位置。图2,图3,图4分别给出了三种滤波方法对目标位置估计均方差、速度估计均方差和转速估计均方差的对比结果。由于0~40 s内三种滤波方法得到估计结果差别不大,故没有显示在图2,图3,图4中。由仿真结果可知,三阶CKF与UKF的估计精度较为接近并且三阶CKF精度稍高于UKF。五阶CKF的估计精度均高于UKF和三阶CKF,且数值稳定性更好。说明高阶球面-相径容积规则可以较大提高高斯权值积分的计算精度,进而提高滤波精度,体现了高阶球面-相径容积规则在估计精度方面的优越性。

图1 目标飞行的平面轨迹和雷达位置 Fig.1 The plane trajectory of the target and the radar's position
图选项

图2 目标位置的估计均方根误差 Fig.2 The RMSE of the target position estimation
图选项

图3 目标速度的估计均方根误差 Fig.3 The RMSE of the target velocity estimation
图选项

图4 目标转速的估计均方根误差 Fig.4 The RMSE of the target turn rate estimation
图选项

表1 仿真参数设置 Table 1 The parameters of simulation
参数项 取值
基本仿真参数T=1,q1=0.1 m2/s3,Ω=-3 (°)/s,q2=1.75×10-4/s3
测量噪声标准差 σr=10 m,σθ=$\sqrt {10} $mrad,ση=$\sqrt {10} $mrad
目标初始状态X0 X0=[1 000 m 300 m/s 1 000 m 0 m/s 2 000 m 0 m/s -3(°)/s]T
初始状态协方差矩阵P0|0 P0|0=diag(100 m2,10 m2/s2,100 m2,10 m2/s2,100 m2,10 m2/s2,100 (mrad/s)2)
表选项

图5对比了10次随机打靶中各滤波方法计算所消耗的CPU时间。其中三阶CKF采用的容积点数最少,为14个容积点,因此算法计算量最小、运行时间最短。而基于高阶球面-相径容积规则的五阶CKF由于在时间更新和测量更新计算中采用99个容积点,提高了计算复杂度,因此算法运算时间最长。

图5 3种滤波方法的执行时间 Fig.5 Execution time of the three methods
图选项

表2给出了三种滤波方法对目标位置、速度、转速的均方根误差均值以及算法的容积点数和平均执行时间。

表2 三种滤波方法的估计精度对比 Table 2 The estimation accuracy of three methods
参数项 UKF 三阶CKF 五阶CKF
位置精度/m 17.38 16.66 14.82
速度精度/(m·s-1) 8.55 7.95 6.04
转速精度/((°)·s-1) 0.99 0.92 0.86
积分点数 15 14 99
平均执行时间/s 0.042 0.036 0.242
表选项

可以清楚看到,五阶CKF的估计精度明显优于其他两种滤波方法。但是相对于UKF和三阶CKF,五阶CKF采用的容积点数较多,运算时间较长。由此可知,五阶CKF在取得较高估计精度的同时,计算量的代价相对较大,特别是在系统状态向量维数n较高的情况下,其使用的容积点数几乎是三阶CKF的n倍。因此,五阶CKF适用于滤波时间间隔较长或对滤波实时性要求不高的高精度滤波系统中。此次仿真中五阶CKF算法的运算时间与滤波时间间隔相差较大,可以满足实时性需求。

4 结论

本文针对传统容积卡尔曼滤波估计精度有限的问题,在基于点的高斯近似滤波框架下,提出高阶容积卡尔曼滤波算法HCKF。该算法采用Genz积分方法和矩匹配法分别推导出任意阶的球面规则和相径规则,以此构造高阶球面-相径容积规则。通过数学仿真可以得出以下结论:

1)对于三维空间中的机动目标跟踪问题,高阶CKF可获得优于UKF和传统CKF的估计精度和数值稳定性,表明高阶球面-相径容积规则可以较大提高高斯权值积分的计算精度,进而提高滤波精度,体现了高阶球面-相径容积规则在估计精度方面的优越性。

2)由于HCKF采用的容积点数量为2n2+1(n为状态向量维数),可知当状态向量维数较大时,采用的容积点较多,计算量相对CKF和UKF较大。因此,作者的下一个研究思路是通过降低状态向量维数或观测向量维数的方法在保持HCKF估计精度的同时明显减小算法计算量。

3)本文研究的HCKF是基于五阶球面-相径容积规则。实际上,采用不同阶数的球面规则和相径规则,就可以得到不同精度的容积规则。因此,一个可行的研究思路是设计更高阶的球面-相径容积滤波算法并对其分析。

参考文献
[1] ITO K, XIONG Kaiqi. Gaussian filters for nonlinear filtering problems[J]. IEEE transactions on automatic control, 2000, 45(5): 910-927.
[2] FORBES J R. Extended Kalman filter and sigma point filter approaches to adaptive filtering [C]//AIAA Guidance, Navigation, and Control Conference. Toronto, Canada, 2010: 7748-7763.
[3] JULIER S J, UHLMANN J K. Unscented filtering and nonlinear estimation [J]. Proceedings of the IEEE, 2004, 92(3): 401-422.
[4] LIU Xing, JIANG Shoushan. Research on target tracking based on unscented Kalman filter[J]. Sensors & transducers journal, 2013, 153(6): 13-21.
[5] 王宏建, 徐金龙, 么洪飞, 等. 基于二阶差分滤波器的水下目标纯方位角跟踪[J]. 哈尔滨工程大学学报, 2014, 35(1): 87-92. WANG Hongjian, XU Jinlong, YAO Hongfei, et al. Research on second order divided difference filter algorithm for underwater target bearing-only tracking [J]. Journal of Harbin Engineering University, 2014, 35(1): 87-92.
[6] ARASARATNAM I, HAYKIN S, ELLIOTT R J, Discrete-time nonlinear filtering algorithms using gauss-hermite quadrature[J]. Proceedings of the IEEE, 2007, 95(5): 953-977.
[7] JIA Bin, XIN Ming, CHENG Yang. Sparse Gauss-Hermite quadrature filter with application to spacecraft attitude estimation[J]. Journal of guidance, control, and dynamics, 2011, 34(2): 367-379.
[8] JIA Bin, XIN Ming, CHENG Yang. Sparse-grid quadrature nonlinear filtering[J]. Automatica, 2012, 48(2): 327-341.
[9] ARASARATNAM I, HAYKIN S. Cubature Kalman filters [J]. IEEE transactions on automatic control, 2009, 54(6): 1254-1269.
[10] ARASARATNAM I, HAYKIN S, HURD, T R. Cubature Kalman filtering for continuous-discrete systems: theory and simulations [J]. IEEE transactions on signal processing, 2010, 58(10): 4977-4993.
[11] FERNÁNDEZ-PRADES C, VILÀ-VALLS J. Bayesian nonlinear filtering using quadrature and cubature rules applied to sensor data fusion for positioning[C]//IEEE International Conference on Communications (ICC). Cape Town, South Africa, 2010: 1-5.
[12] PESONEN H, PICHÉ R, Cubature-based Kalman filters for positioning[C]//7th Workshop on Positioning Navigation and Communication (WPNC). Dresden, 2010: 45-49.
[13] LI Qiurong, SUN Feng. Strong tracking cubature Kalman filter algorithm for GPS/INS integrated navigation system[C]//IEEE International Conference on Mechatronics and Automation (ICMA). Takamatsu, 2013: 1113-1117.
[14] LI Chao, GE Quanbo. SCKF for MAV attitude estimation[C]//2011 International Conference on Machine Learning and Cybernetics. Guilin, China, 2011: 1313-1318.
[15] TAO Sun, MING Xin. Hypersonic entry vehicle state estimation using high-degree cubature Kalman filter[C]//AIAA Atmospheric Flight Mechanics Conference. Atlanta, 2014: 2383-2394.
[16] JIA Bin, XIN Ming, CHENG Yang. The high-degree cubature Kalman filter[C]//51st IEEE Conference on Decision and Control. Maui, USA, 2012: 4095-4100.
[17] JIA Bin, XIN Ming, CHENG Yang. High-degree cubature Kalman filter [J]. Automatica, 2013, 49(2): 510-518.
[18] GENZ, A. Fully symmetric interpolatory rules for multiple integrals over hyper-spherical surfaces[J]. Journal of computational and applied mathematics, 2003, 157(1): 187-195.
[19] GENZ A, MONAHAN J. A Stochastic algorithm for high-dimensional integrals over unbounded regions with Gaussian weight[J]. Journal of computational and applied mathematics, 1999, 112(1/2): 71-81.
DOI: 10.11990/jheu.201412079
0

文章信息

张龙, 崔乃刚, 杨峰, 路菲, 卢宝刚
ZHANG Long, CUI Naigang, YANG Feng, LU Fei, LU Baogang
高阶容积卡尔曼滤波及其在目标跟踪中的应用
High-degree cubature Kalman filter and its application in target tracking
哈尔滨工程大学学报, 2016, 37(04): 573-578
Journal of Harbin Engineering University, 2016, 37(04): 573-578
DOI: 10.11990/jheu.201412079

文章历史

收稿日期:2014-12-31

相关文章

工作空间