一种适用于弹道导弹ins/cns/gnss组合导航系统的最优数据融合方法

文档序号:1672196 发布日期:2019-12-31 浏览:30次 >En<

阅读说明:本技术 一种适用于弹道导弹ins/cns/gnss组合导航系统的最优数据融合方法 (Optimal data fusion method suitable for ballistic missile INS/CNS/GNSS combined navigation system ) 是由 陈熙源 柳笛 刘晓 石春凤 马振 于 2019-09-24 设计创作,主要内容包括:本发明公开了一种适用于弹道导弹INS/CNS/GNSS组合导航系统的最优数据融合方法,该方法包括以下步骤:构造INS/CNS/GNSS组合导航系统模型;在广义高阶CKF的时间更新阶段和量测更新阶段分别引入自适应渐消因子和最大相关熵准则进行INS/CNS子系统和INS/GNSS子系统的局部状态估计;根据最小方差原理和容积准则融合INS/CNS子系统和INS/GNSS子系统的局部估计得到全局最优状态估计。本发明可以同时抑制过程建模误差和非高斯量测噪声对状态估计的影响,提高弹道导弹INS/CNS/GNSS组合导航的自适应性和鲁棒性,获得全局最优的状态估计值。(The invention discloses an optimal data fusion method suitable for a ballistic missile INS/CNS/GNSS combined navigation system, which comprises the following steps: constructing an INS/CNS/GNSS integrated navigation system model; respectively introducing a self-adaptive fading factor and a maximum correlation entropy criterion to carry out local state estimation of the INS/CNS subsystem and the INS/GNSS subsystem in a time updating stage and a measurement updating stage of the generalized high-order CKF; and according to the minimum variance principle and the volume criterion, local estimation of the INS/CNS subsystem and the INS/GNSS subsystem is fused to obtain global optimal state estimation. The method can simultaneously inhibit the influence of process modeling errors and non-Gaussian measurement noise on state estimation, improve the adaptivity and robustness of ballistic missile INS/CNS/GNSS combined navigation, and obtain the globally optimal state estimation value.)

一种适用于弹道导弹INS/CNS/GNSS组合导航系统的最优数据 融合方法

技术领域

本发明涉及一种数据融合方法,具体涉及一种适用于弹道导弹INS/CNS/GNSS组合导航系统的最优数据融合方法。

背景技术

传统的基于EKF、UKF的联邦多传感器数据融合方法使用过程噪声协方差的上界而不是过程噪声协方差本身来消除局部状态估计之间的相关性,因此得到的全局状态估是次优的。此外,该类方法要求非线性系统模型必须达到足够的精度。然而,在弹道导弹飞行的过程中,一方面,其弹载的INS/CNS/GNSS组合导航系统很容易受到复杂环境的影响,进而导致INS/CNS/GNSS组合导航系统的量测信息受到非高斯噪声的影响,在这种情况下如果继续使用传统的基于高斯假设的非线性滤波方法将会使组合导航系统的导航精度严重下降;另一方面,由于真实的动态系统模型非常复杂,我们建立的系统过程模型只能是真实模型的理论近似,因而存在建模误差,这种建模误差也会导致组合导航系统的导航精度的下降。

所以,需要一个新的技术方案来解决上述问题。

发明内容

发明目的:为了克服现有技术中存在的组合导航系统的导航精度差的不足,提供一种适用于弹道导弹INS/CNS/GNSS组合导航系统的最优数据融合方法,其实现在过程模型含有不确定性以及量测噪声为非高斯噪声的情况下依然可以获取全局最优的状态估计。

技术方案:本发明提供一种适用于弹道导弹INS/CNS/GNSS组合导航系统的最优数据融合方法,包括以下步骤:

S1:构造INS/CNS/GNSS组合导航系统滤波模型;

S2:在广义高阶CKF的时间更新阶段和量测更新阶段分别引入自适应渐消因子和最大相关熵准则进行INS/CNS子系统和INS/GNSS子系统的局部状态估计;

S3:根据最小方差原理和容积准则融合INS/CNS子系统和INS/GNSS子系统的局部估计得到全局最优状态估计。

进一步的,所述步骤S1中构造INS/CNS/GNSS组合导航系统滤波模型包括以下步骤:

S1-1:设置INS/CNS/GNSS组合导航系统的状态向量为x=[φx,φy,φz,δvx,δvy,δvz,x,y,z,εx,εy,εz,Δx,Δy,Δz]T,其中φx,φy,φz表示发射点惯性系下姿态失准角,δvx,δvy,δvz表示发射点惯性系下的速度误差,x,y,z表示发射点惯性系下的位置误差,εx,εy,εz表示弹体坐标系下的陀螺常值漂移,Δx,Δy,Δz表示弹体坐标系下的加速度计常值偏置,T为转置符号;

S1-2:根据INS/CNS/GNSS组合导航系统的15维状态向量x建立系统的状态方程:

xk=f(xk-1)+vk-1

其中,f(·)为非线性系统函数,xk-1和xk分别表示k-1时刻和k时刻的状态向量,vk-1表示过程噪声,vk-1的协方差为

Figure BDA0002212558140000021

S1-3:分别建立INS/GNSS子系统和INS/CNS子系统的量测方程:

在INS/GNSS子系统中,将INS与GNSS输出的位置的差值和速度的差值作为量测信息,建立该子系统的量测方程:

z1,k=H1,kxk1,k

其中,z1,k表示INS/GNSS子系统k时刻的量测向量,H1,k表示INS/GNSS子系统k时刻的测量矩阵,ω1,k表示INS/GNSS子系统k时刻的量测噪声,其方差是

在INS/CNS子系统中,将INS与CNS输出的姿态角的差值作为量测信息,建立该子系统的量测方程:

z2,k=H2,kxk2,k

其中,z2,k表示INS/CNS子系统k时刻的量测向量,H2,k表示INS/CNS子系统k时刻的测量矩阵,ω2,k表示INS/CNS子系统k时刻的量测噪声,其方差为

Figure BDA0002212558140000023

进一步的,所述步骤S2具体包括以下步骤:

S2-1:因为INS/CNS子系统和INS/GNSS子系统采用同样的滤波过程进行局部状态估计,为了避免重复,在此仅具体说明INS/GNSS子系统的滤波过程。因此,将步骤S1-3中的z1,k、z2,k统一写成zk,H1,k和H2,k统一写成Hk,ω1,k和ω2,k统一写成ωk,R1,k和R2,k统一写成Rk,根据INS/GNSS子系统的局部状态估计要求设置核宽度γ的值,设置初始的状态向量、状态误差协方差、渐消因子分别为

Figure BDA0002212558140000024

τk(0)=1,S0|0=chol(P0|0),chol(·)表示乔列斯分解运算;

S2-2:根据公式

Figure BDA0002212558140000031

计算出传播容积点

Figure BDA0002212558140000032

其中i=1,...,2n2+1,n表示状态向量的维数,

Figure BDA0002212558140000033

[·]i表示集合[·]的第i列;

S2-3:预测k时刻状态向量

Figure BDA0002212558140000034

误差协方差阵

Figure BDA0002212558140000035

和计算Sk|k-1=chol(Pk|k-1),其中

S2-4:根据

Figure BDA0002212558140000037

和Sk|k-1产生新的容积点并预测k时刻的量测信息

Figure BDA0002212558140000039

S2-5:根据INS/GNSS子系统的量测方程和步骤S2-3中的

Figure BDA00022125581400000310

和Pk|k-1构建如下方程:

Figure BDA00022125581400000311

其中,I表示单位矩阵,

Figure BDA00022125581400000312

Mp,k|k-1=chol(Pk|k-1)、Mr,k=chol(Rk)和

在方程

Figure BDA00022125581400000314

的两边左乘

Figure BDA00022125581400000315

得到:

Figure BDA00022125581400000316

Figure BDA00022125581400000317

则上述方程可重写为:Dk=BkXk+ek

S2-6:对量测噪声方差进行更新:

Figure BDA0002212558140000041

其中,diag(·)表示矩阵的对角化运算,m为量测信息的维度,

Figure BDA0002212558140000042

Figure BDA0002212558140000043

di,k表示Dk的第i个元素,bi,k表示Bk的第i行元素;

S2-7:计算状态信息和量测信息之间的互协方差矩阵Pxz,k|k-1并令i=0;

S2-8:令

Figure BDA0002212558140000045

Figure BDA0002212558140000046

并设定卡方分布

Figure BDA0002212558140000047

的值,其中θ表示卡方分布的自由度,α表示卡方分布的分位数,如果

Figure BDA0002212558140000048

则计算

Figure BDA0002212558140000049

Figure BDA00022125581400000410

返回步骤S2-2继续执行下一滤波周期;

S2-9:

如果

Figure BDA00022125581400000411

计算

Figure BDA00022125581400000412

其中,i=(0,1,...),i=i+1,根据得到自适应渐消因子τk更新ak(i)的值,判断ak(i)是否满足

Figure BDA00022125581400000413

若不满足,则返回步骤S2-9继续执行;若满足,则得到最优的自适应渐消因子τk的值,然后继续执行后续步骤;

S2-10:使用自适应渐消因子τk更新预测的状态预测协方差:

Figure BDA00022125581400000414

计算卡尔曼增益:

Figure BDA00022125581400000415

估计后验状态和后验协方差矩阵:

Figure BDA00022125581400000416

返回步骤S2-2继续执行下一滤波周期。

进一步的,所述步骤S3包括以下步骤:

S3-1:在分别对INS/GNSS子系统和INS/CNS子系统执行完滤波过程后便可获得相应的局部后验状态估计

Figure BDA00022125581400000417

为了便于区分,令INS/GNSS子系统获得局部后验状态估计为

Figure BDA0002212558140000051

INS/CNS子系统获得局部后验状态估计为

Figure BDA0002212558140000052

53-2:根据最小方差原理和容积准则得到全局最优状态估计:

Figure BDA0002212558140000053

其中令β=[β1,β2]T

Figure BDA0002212558140000054

P11为INS/GNSS子系统执行局部状态估计得到的Pk|k,P22为INS/CNS子系统执行局部状态估计得到的Pk|k,P12表示

Figure BDA0002212558140000055

Figure BDA0002212558140000056

估计误差的互协方差矩阵,P21表示

Figure BDA0002212558140000057

Figure BDA0002212558140000058

估计误差的互协方差矩阵。

进一步的,P12和P21通过容积准则近似得到:

Figure BDA0002212558140000059

Figure BDA00022125581400000510

其中,

Figure BDA00022125581400000511

分别表示INS/GNSS子系统与INS/CNS子系统执行步骤S2-2得到的传播容积点

Figure BDA00022125581400000512

Figure BDA00022125581400000513

分别表示INS/GNSS子系统与INS/CNS子系统执行步骤S2-3得到的状态向量预测值,

Figure BDA00022125581400000515

分别表示INS/GNSS子系统与INS/CNS子系统执行步骤S2-4得到的新的容积点Xi,k|k-1

Figure BDA00022125581400000517

分别表示INS/GNSS子系统与INS/CNS子系统执行步骤S2-4得到的量测信息预测值

Figure BDA00022125581400000518

Figure BDA00022125581400000520

分别表示INS/GNSS子系统与INS/CNS子系统执行局部状态估计过程得到的滤波增益

Figure BDA00022125581400000522

I表示n维的单位矩阵。

有益效果:本发明与现有技术相比,在广义高阶CKF的时间更新阶段和量测更新阶段分别引入自适应渐消因子和最大相关熵准则进行INS/GNSS子系统和INS/CNS子系统的局部状态估计,最后根据最小方差原理和容积准则融合INS/CNS子系统和INS/GNSS子系统的局部估计,使得弹道导弹INS/CNS/GNSS组合导航系统在系统模型存在过程建模误差和量测噪声为非高斯分布的情况下,该组合导航系统依然可以获得全局最优的状态估计,本发明可以同时抑制过程建模误差和非高斯量测噪声对状态估计的影响,提高弹道导弹INS/CNS/GNSS组合导航的自适应性和鲁棒性,从而保证了INS/CNS/GNSS组合导航系统的导航精度。

附图说明

图1为本发明方法的流程图。

具体实施方式

下面结合附图和具体实施例,进一步阐明本发明。

如图1所示,本发明提供一种适用于弹道导弹INS/CNS/GNSS组合导航系统的最优数据融合方法,包括以下步骤:

S1:构造INS/CNS/GNSS组合导航系统滤波模型:

S1-1:设置INS/CNS/GNSS组合导航系统的状态向量为:

x=[φx,φy,φz,δvx,δvy,δvz,x,y,z,εx,εy,εz,Δx,Δy,Δz]T (1)

其中φx,φy,φz表示发射点惯性系下姿态失准角,δvx,δvy,δvz表示发射点惯性系下的速度误差,x,y,z表示发射点惯性系下的位置误差,εx,εy,εz表示弹体坐标系下的陀螺常值漂移,Δx,Δy,Δz表示弹体坐标系下的加速度计常值偏置,T为转置符号。

S1-2:根据INS/CNS/GNSS组合导航系统的15维状态向量x建立系统的状态方程:

xk=f(xk-1)+vk-1 (2)

其中,f(·)为非线性系统函数,xk-1和xk分别表示k-1时刻和k时刻的状态向量,vk-1表示过程噪声,vk-1的协方差为

Figure BDA0002212558140000061

S1-3:分别建立INS/GNSS子系统和INS/CNS子系统的量测方程:

在INS/GNSS子系统中,将INS与GNSS输出的位置和速度分别做差,并将其作为量测信息,建立该子系统的量测方程:

z1,k=H1,kxk1,k (3)

其中,z1,k表示INS/GNSS子系统k时刻的量测向量,H1,k表示INS/GNSS子系统k时刻的测量矩阵,ω1,k表示INS/GNSS子系统k时刻的量测噪声,其方差是

Figure BDA0002212558140000062

在INS/CNS子系统中,该子系统的量测方程将INS与CNS输出的姿态角的差值作为量测信息,建立该子系统的量测方程:

z2,k=H2,kxk2,k (4)

其中,z2,k表示INS/CNS子系统k时刻的量测向量,H2,k表示INS/CNS子系统k时刻的测量矩阵,ω2,k表示INS/CNS子系统k时刻的量测噪声,其方差为

Figure BDA0002212558140000063

S2:在广义高阶CKF的时间更新阶段和量测更新阶段分别引入自适应渐消因子和最大相关熵准则进行INS/CNS子系统和INS/GNSS子系统的局部状态估计,其具体过程如下:

S2-1:因为INS/GNSS子系统和INS/CNS子系统分别通过局部滤波器1和局部滤波器2采用同样的滤波过程进行局部状态估计,为了避免重复,在此仅具体说明INS/GNSS子系统的滤波过程:,因此,本实施例将步骤S1中的z1,k、z2,k统一写成zk,H1,k和H2,k统一写成Hk,ω1,k和ω2,k统一写成ωk,R1,k和R2,k统一写成Rk。根据INS/GNSS子系统的局部状态估计要求设置核宽度γ的值,设置初始的状态向量、状态误差协方差、渐消因子分别为

Figure BDA0002212558140000071

τk(0)=1,S0|0=chol(P0|0),chol(·)表示乔列斯分解运算。

S2-2:根据公式(5)计算出传播容积点

Figure BDA0002212558140000072

Figure BDA0002212558140000073

其中i=1,...,2n2+1,n表示状态向量的维数,[·]i表示集合[·]的第i列。

S2-3:预测k时刻状态向量和误差协方差阵:

Figure BDA0002212558140000075

其中

Figure BDA0002212558140000077

S2-4:根据

Figure BDA0002212558140000078

阳Sk|k-1产生新的容积点:

其中,Sk|k-1=chol(Pk|k-1)。

S2-5:预测k时刻的量测信息:

S2-6:根据INS/GNSS子系统的量测方程和和Pk|k-1构建如下方程:

Figure BDA00022125581400000712

其中,I表示单位矩阵,Mp,k|k-1=chol(Pk|k-1)、Mr,k=chol(Rk)和

Figure BDA0002212558140000082

在方程

Figure BDA0002212558140000083

的两边左乘

Figure BDA0002212558140000084

得到:

Figure BDA0002212558140000085

Figure BDA0002212558140000086

则公式(11)可重写为:

Dk=BkXk+ek (12)

S2-7:对量测噪声方差进行更新:

Figure BDA0002212558140000087

其中,diag(·)表示矩阵的对角化运算,m为量测信息的维度,

Figure BDA0002212558140000089

di,k表示Dk的第i个元素,bi,k表示Bk的第i行元素。

S2-8:计算状态信息和量测信息之间的互协方差矩阵Pxz,k|k-1

Figure BDA00022125581400000810

令i=0,

Figure BDA00022125581400000811

Figure BDA00022125581400000812

并设定卡方分布

Figure BDA00022125581400000813

的值,其中θ表示卡方分布的自由度,α表示卡方分布的分位数;

如果

Figure BDA00022125581400000814

则计算:

Figure BDA00022125581400000816

Figure BDA00022125581400000817

Figure BDA0002212558140000091

Sk|k=chol(Pk|k) (19)

继续执行k=k+1,并返回步骤S2-2中的公式(5)所在位置开始执行下一滤波周期;

如果

Figure BDA0002212558140000092

计算:

Figure BDA0002212558140000093

根据得到自适应渐消因子τk更新ak(i)的值,判断ak(i)是否满足

Figure BDA0002212558140000094

若不满足,i=i+1,则返回公式(20)更新渐消因子τk的值;若满足,则得到最优的自适应渐消因子τk的值,然后继续执行下列步骤。

S2-9:使用得到的自适应渐消因子τk更新预测的状态预测协方差:

Figure BDA0002212558140000095

计算卡尔曼增益:

S2-10:估计后验状态和后验协方差矩阵:

Figure BDA0002212558140000097

Figure BDA0002212558140000098

Sk|k=chol(Pk|k) (25)

执行k=k+1,并返回步骤S2-2中的公式(5)所在位置开始执行下一滤波周期。

S3:根据最小方差原理和容积准则融合INS/GNSS子系统和INS/CNS子系统的局部估计得到全局最优状态估计:

S3-1:在分别对INS/GNSS子系统和INS/CNS子系统执行完以上滤波过程后便可获得相应的局部后验状态估计

Figure BDA0002212558140000099

为了便于区分,令INS/GNSS子系统获得局部后验状态估计为

Figure BDA00022125581400000910

INS/CNS子系统获得局部后验状态估计为

Figure BDA00022125581400000911

S3-2:根据最小方差原理和容积准则得到全局最优状态估计:

Figure BDA00022125581400000912

其中令β=[β1,β2]T

Figure BDA0002212558140000101

P11为INS/GNSS子系统执行局部状态估计得到的Pk|k,P22为INS/CNS子系统执行局部状态估计得到的Pk|k,P12表示

Figure BDA0002212558140000102

估计误差的互协方差矩阵,P21表示

Figure BDA0002212558140000104

估计误差的互协方差矩阵,P12和P21通过容积准则近似得到:

Figure BDA0002212558140000106

Figure BDA0002212558140000107

其中,分别表示INS/GNSS子系统与INS/CNS子系统执行步骤S2的过程中得到的传播容积点

Figure BDA0002212558140000109

Figure BDA00022125581400001010

Figure BDA00022125581400001011

分别表示INS/GNSS子系统与INS/CNS子系统执行步骤S2的过程中得到的状态向量预测值,

Figure BDA00022125581400001012

Figure BDA00022125581400001013

分别表示INS/GNSS子系统与INS/CNS子系统执行步骤S2的过程中得到的新的容积点Xi,k|k-1

Figure BDA00022125581400001014

Figure BDA00022125581400001015

分别表示INS/GNSS子系统与INS/CNS子系统执行步骤S2的过程中得到的量测信息预测值

Figure BDA00022125581400001017

分别表示INS/GNSS子系统与INS/CNS子系统执行步骤S2的过程中得到的滤波增益

Figure BDA00022125581400001019

I表示n维的单位矩阵。

15页详细技术资料下载
上一篇:一种医用注射器针头装配设备
下一篇:自动驾驶车辆的定位方法、装置、电子设备及可读介质

网友询问留言

已有0条留言

还没有人留言评论。精彩留言会获得点赞!

精彩留言,会给你点赞!

技术分类