中国科技核心期刊      中国指挥与控制学会会刊     军事装备类重点期刊
Command Control and Simulation

Command Control and Simulation >

2023 , Vol. 45 >Issue 4: 154 - 160

DOI: https://doi.org/10.3969/j.issn.1673-3819.2023.04.023

Radar & Communication & Navigation

A navigation method based on the inertial navigation information, satellite positioning and scene matching points

  • WANG Jianshu ,
  • YIN Kuiying ,
  • YIN Jinrong ,
  • HUANG Zhaoyou ,
  • WANG Zhongbao
Expand
  • Nanjing Research Institute of Electronics Technology, Nanjing 210000, China

Received date: 2022-11-01

  Revised date: 2022-12-10

  Online published: 2023-08-17

Fold

Abstract

Focus on the poor navigation precision when the aircraft suffers from the loss of satellite signals, the navigation method of the Synthetic Aperture Radar (SAR) platform based on Inertial Navigation System (INS), satellite positioning and scene matching is proposed. First, the corresponding equations of measurements are formed by the models of the INS, satellite positioning and scene matching points, respectively. And then applying the Kalman Filter (KF), the corresponding filter processes are established according to the different equations of measurements. The Extended Kalman Filter (EKF) is used to linearize the equations of satellite positioning and scene matching measurements. The proposed method can be executed when the satellites are not enough to provide the positioning information of the aircraft. In numeral simulations, when the satellite signals are partially missing, the proposed method can obtain lower absolute positioning errors than the INS and satellite-based method, and the proposed method obtains better precisian than the INS and satellite-based method when the satellite signals always exist. The simulation results prove the effectiveness and the superior performance of the proposed method.

Key words： inertial navigation system; satellite positioning; scene matching; Kalman filter; synthetic aperture radar

Cite this article

WANG Jianshu , YIN Kuiying , YIN Jinrong , HUANG Zhaoyou , WANG Zhongbao . A navigation method based on the inertial navigation information, satellite positioning and scene matching points[J]. Command Control and Simulation, 2023 , 45(4) : 154 -160 . DOI: 10.3969/j.issn.1673-3819.2023.04.023

飞行器在航行过程中的定位导航是实际工程中需要解决的重要问题。在目前大多数场景中,惯性导航系统(Inertial Navigation System, INS)与卫星导航是应用最广泛的。INS具有较强的抗干扰能力,能全天候工作,但精度不高,且误差会随时间积累[1]。卫星导航主要有全球定位系统(Global Positioning System, GPS)与北斗卫星导航系统等,精度较高,但会受制于时政与干扰等因素。目前将INS与卫星导航组合应用,可以充分利用各导航系统的优势,提高定位精度。然而,当卫星信号失效时,由于仅INS作用,该导航方法获得的导航信息误差会不断增大,致使飞行器面临风险[2-3]。
对于合成孔径雷达(Synthetic Aperture Radar, SAR)[4]载机平台,充分利用实时测得的地面SAR图像与内置参考图像的景象匹配[5-6]结果,可以增强传统导航方法的可靠性。陈圣义等[7]研究了景象匹配点与载机的几何关系,构造了景象匹配与INS的联合定位方法。Huang等[8]研究了6维INS误差模型下利用景象匹配的SAR载机平台位置实时解算方法。李家松等[9]研究了GPS拒止条件下的景象匹配导航方法。Liao等[10]研究了SINS、SAR与GPS的联合定位方法。上述定位导航方法均需要在各量测同时获得的条件下进行,在实际场景下,景象匹配点往往具有较低的采样率,这限制了这些方法的应用。另外,当平台连接的卫星数量较少以致无法给出定位结果时,部分利用GPS信息的导航方法将失效。
本文提出一种基于INS、卫星定位与景象匹配的SAR平台导航方法。首先根据INS、卫星定位与景象匹配的数学模型,构造了对应的量测方程,然后应用卡尔曼滤波(Kalman Filter, KF)[11],对不同的量测推导了对应的滤波算法。其中对于卫星数据与景象匹配数据,通过扩展卡尔曼滤波(Extended Kalman Filter, EKF)[12]进行了线性化处理。该方法不要求卫星数据直接提供定位信息,故在卫星数较少时也可以应用。仿真实验中,当卫星信号缺失时,本文方法获得比INS与卫星组合定位方法更小的位置误差,在包含卫星信号时,也能提供比INS与卫星组合定位方法更优的精度。该结果证明了本文方法的有效性与良好性能。

1 数学模型

图1显示了SAR平台与卫星以及景象匹配点之间的几何关系。图1a)中,q(t)=[x,y,z]T与pG(t)=[xG,yG,zG]T分别表示SAR平台与卫星在世界坐标系中t时刻的位置坐标,rG(t)为该卫星测得SAR平台的伪距。为简单起见,忽略卫星钟差等因素的影响,该测量的伪距有如下关系
rG(t)=‖q(t)-pG(t)‖2+nG(t)
图1 SAR平台与卫星和景象匹配点几何关系示意图

Fig.1 The geometric relation between the SAR platform and the scene-matching point

式中,‖·‖2为向量的二范数,nG(t)为伪距测量误差,并假设服从高斯分布nG(t)∈N(0, δ G 2)。图1b)显示了世界坐标系下SAR平台与景象匹配点的关系,xoy平面为成像区域的当地水平面,oz为当地的高度方向。在t1时刻,SAR平台获得一个地面SAR成像与参考图像的景象匹配点pM(t1),以及该点的斜距rM(t1),在不计测量误差的情况下其几何关系可以表示为
rM(t1)=‖q(t1)-pM(t1)‖2
实际情况下,景象匹配点pM(t1)的xy坐标由SAR图像与参考图像局部匹配测得,z坐标由对应xy坐标的数字高程提供,假设匹配点坐标的测量误差服从0均值高斯分布,且协方差为diag{ δ p x 2, δ p y 2, δ p z 2},其中diag{·}表示构造对角矩阵的运算。斜距rM(t1)由SAR测得,假设测量误差服从0均值方差为 δ r M 2 高斯分布。同理可见示意图上t2与t3时刻SAR平台与景象匹配的几何关系。

2 基于INS卫星定位与景象匹配组合导航方法

本文提出一种INS、卫星定位与景象匹配组合导航方法。首先根据数学模型,构建状态方程与量测方程,然后通过KF与EKF,对不同量测构造不同的滤波流程,实时估计SAR平台的导航信息。

2.1 状态方程

令状态矢量为X=[x, x ˙, x ¨,y, y ˙, y ¨,z, z ˙, z ¨]T,式中x,y,z为SAR平台的位置坐标, x ˙, y ˙, z ˙分别为SAR平台在x,y,z方向的速度分量, x ¨, y ¨, z ¨分别为SAR平台在x,y,z方向的加速度分量。则状态方程可以写为
Xk+1|k=FkXk|k+Vk
式中,Xk+1|k与Xk|k分别为第k+1和k时刻滤波后的状态,Fk为第k时刻状态转移矩阵,且Fk= F ~ 0 3 × 3 0 3 × 3 0 3 × 3 F ~ 0 3 × 3 0 3 × 3 0 3 × 3 F ~,其中 F ~= 1 T k T k 2 / 2 0 1 T k 0 0 1,Tk为第k时刻和第k+1时刻的时间差,Vk包含了噪声与其他输入信号,假设其满足均值为0,协方差为Qk的高斯分布,即Vk~N(09×1,Qk)。

2.2 量测方程

INS、卫星定位与景象匹配系统独立获得相关信息的测量值,一般具有不同的采样率与采样时间,故需要分别构造量测方程和滤波流程。INS获得的数据转化到世界坐标系下,可以获得SAR平台的加速度ax、ay、az,以之作为量测,可以构造线性的量测方程。卫星提供的目标伪距作为量测,然后根据卫星坐标构造非线性量测方程。景象匹配系统中,根据其景象匹配点pM与该点的斜距rM作为量测,构造非线性的量测方程。以下对这三种情形分别讨论。
1)INS量测情形
该情形下,量测为SAR平台的加速度ax、ay、az。令量测向量为Z=[ax,ay,az]T,则量测方程可以写为
Zk+1|k=Hk+1Xk+1|k+Wk+1
式中,Zk+1|k为第k+1时刻的量测向量,Hk+1为量测矩阵,满足Hk+1= H ~ 0 1 × 3 0 1 × 3 0 1 × 3 H ~ 0 1 × 3 0 1 × 3 0 1 × 3 H ~,其中 H ~=[0,0,1],Wk+1为零均值协方差为R的高斯白噪声,即满足Wk+1~N(03×1,R)。
2)卫星量测情形
该情形下,量测为卫星测量的伪距rG。先考虑一颗卫星的情形,令量测向量为Z=rG,则量测方程可以写为
Zk+1|k=f(Xk+1|k, p G k + 1)+Wk+1
式中,f(Xk+1|k, p G k + 1)=‖ H -Xk+1|k- p G k + 1‖2为量测的非线性转化函数,其中 H -= H ~ 0 1 × 3 0 1 × 3 0 1 × 3 H ~ 0 1 × 3 0 1 × 3 0 1 × 3 H ~, H ~=[1,0,0], p G k + 1为该时刻卫星的真实位置,对应的斜距为 r G k + 1。该情形下,Wk+1~N(0,R),其中Wk+1的方差为R= δ G 2。在多颗卫星的情形下,量测方程(5)可以按行合并成一个表达式。假设卫星数为N,令Z=[ r G ( 1 ), r G ( 2 ),…, r G ( N )]表示N颗卫星测得的伪距,则量测方程(5)可以写为
Zk+1|k=[f(Xk+1|k, p G k + 1 ( 1 )),f(Xk+1|k, p G k + 1 ( 2 )),…,f(Xk+1|k, p G k + 1 ( N ))]T+Wk+1
其中, p G k + 1 ( 1 ), p G k + 1 ( 2 ),…, p G k + 1 ( N ) 表示N颗卫星的位置,为表达简便,使用Wk+1~N(0N×1,R)表示维数扩展后的噪声向量。
3)景象匹配量测情形
该情形下,景象匹配点与斜距的几何关系互相耦合,且为非线性关系,本文提出一种近似的量测方程构造方法。令量测向量为Z=rM,则量测方程可以写为
Zk+1=f(Xk+1|k, p M k + 1)+Wk+1
式中,f(Xk+1|k, p M k + 1)=‖ H -Xk+1|k- p M k + 1‖2为量测的非线性转化函数,其中 H -= H ~ 0 1 × 3 0 1 × 3 0 1 × 3 H ~ 0 1 × 3 0 1 × 3 0 1 × 3 H ~, H ~=[1,0,0], p M k + 1为该时刻景象匹配点真实位置,对应的斜距为 r M k + 1。令 p M k + 1= p ~ M+Δ p M k + 1,其中 p ~ M为该时刻景象匹配点测得位置,Δ p M k + 1为真实值与测量值的误差。对f(Xk+1|k, p M k + 1)在 p M k + 1= p ~ M处进行一阶泰勒展开,可得
f(Xk+1|k, p M k + 1)≈‖ H -Xk+1|k- p ~ M‖2+ J ~Δ p M k + 1
式中, J ~= ∂ f ( X k + 1 | k , p M k + 1 ) ∂ p M k + 1| p M k + 1 = p ~ M=-( H -Xk+1|k- p ~ M)T/‖ H -Xk+1|k- p ~ M‖2为3维行向量。将式(8)代入(7)并忽略近似项可得新的量测方程
Zk+1=h(Xk+1|k)+Uk+1
式中,h(Xk+1|k)=‖ H -Xk+1|k- p ~ M‖2,Uk+1=Wk+1+ J ~Δ p M k + 1。由于景象匹配点位置精度为δp=[ δ p x, δ p y, δ p z]T,假设Δ p M k + 1~N(03×1,diag{ δ p x 2, δ p y 2, δ p z 2}),则 J ~Δ p M k + 1~N(03×1, J ~diag{ δ p x 2, δ p y 2, δ p z 2} J ~ T)。由于Δ p M k + 1与Wk+1不相关,进而可得Uk+1=Wk+1+ J ~Δ p M k + 1~N(03×1,R+ J ~diag{ δ p x 2, δ p y 2, δ p z 2} J ~ T)。

2.3 滤波算法流程

对上述三种不同量测方程使用不同的卡尔曼滤波流程。
1)INS量测情形
对量测方程(4),使用标准卡尔曼滤波进行处理。具体地,计算状态的一步预测:
Xk+1|k=FkXk|k
计算新息
β=Zk+1-Hk+1Xk+1|k
计算状态协方差的一步预测
PK+1|k=FkPk|k F T k+Qk
估计新息协方差矩阵
S=Hk+1PK+1|k H T k + 1+R
计算增益
K=PK+1|k H T k + 1 S-1
更新状态变量
Xk+1|k+1=Xk+1|k+Kβ
更新状态协方差矩阵
PK+1|k+1=PK+1|k-KSKT
上述式(10)至式(16)为量测为加速度情形下的滤波流程。
2)卫星量测情形
对该情形的非线性量测方程(6)使用EKF流程。对量测的预测向量一阶泰勒展开,代替其非线性表示。计算函数f(X, p G k + 1)的雅克比矩阵可得
J k + 1 , p G k + 1= ∂ f ( X , p G k + 1 ) ∂ X| X = X k + 1 | k= H - T( H -Xk+1|k- p G k + 1)/‖ H -Xk+1|k- p G k + 1‖2
则量测矩阵为
H ⌢=[ J k + 1 , p G k + 1 ( 1 ), J k + 1 , p G k + 1 ( 2 ),…, J k + 1 , p G k + 1 ( N )]T
新息为
β=Zk+1-f(Xk+1|k, p G k + 1)
估计新息协方差矩阵
S= H ⌢PK+1|k H ⌢ T+R
计算增益
K=PK+1|k H ⌢ TS-1
其他滤波流程与标准卡尔曼滤波一致。则量测方程(6)的滤波流程分别为式(10)、(18)、(19)、(20)、(21)、(15)与(16)。
3)景象匹配量测情形
对该情形的非线性量测方程(9)使用EKF流程。计算函数h(Xk+1|k)的雅克比矩阵可得
Jk+1= ∂ h ( X ) ∂ X| X = X k + 1 | k= H - T( H -Xk+1|k- p ~ k + 1)/‖ H -Xk+1|k- p ~ k + 1‖2
新息为
β=Zk+1-h(Xk+1|k)
估计新息协方差矩阵
S=Jk+1PK+1|k J T k + 1+R+ J ~diag{ δ p x 2, δ p y 2, δ p z 2} J ~ T
计算增益
K=PK+1|k J T k + 1 S-1
其他滤波流程与标准卡尔曼滤波一致。则量测方程(9)的滤波流程分别为式(10)、(22)、(23)、(24)、(25)、(15)与(16)。
4)滤波流程
分别记上述INS量测情形、卫星量测和景象匹配量测情形的滤波器为滤波器1、滤波器2和滤波器3,则本文方法的滤波流程图如图2所示。本文的滤波方法根据此时刻获得量测的不同动态地进行选择,由于各量测的采样率差异,3组滤波器一般是独立地运行,但由于均作用于该飞行器运动的状态方程,故能不断地提供高精度的导航信息。
图2 滤波流程图

Fig.2 Flow chart of filtering

3 仿真实验

假设SAR平台在t∈[0,1 000]的时间内,由0~200 s匀速运动,201 s~500 s匀加速运动,501 s~710 s匀加速运动与711 s~1 000 s的匀速运动构成一段飞行轨迹,如图3所示。假设INS测量加速度与卫星测量伪距的采样率均为fS=100 Hz,INS测得加速度精度为δINS=[ δ a x, δ a y, δ a z]T=[0.05,0.05,0.05]T m/s2。为便于仿真,假设同时有3颗静止轨道卫星连接SAR平台,卫星坐标分别为 p G ( 1 )=[0,0,35 790]T, p G ( 2 )=[20 000,0,35 790]T与 p G ( 3 )=[0,20 000,35 790]T,测量伪距精度均为δG=20 m。假设机载SAR每1 s给出一个景象匹配点位置p(t),精度为δp=[3,3,10]T m,以及该位置的斜距r(t),精度为δr=10 m。
图3 SAR平台飞行轨迹

Fig.3 The flight track of the SAR platform

当量测为INS测得加速度时,量测噪声协方差设为R=diag{ δ I N S 2}。当量测为卫星测得伪距时,量测噪声协方差设为R=diag{ δ G 2, δ G 2, δ G 2}。当量测为景象匹配点与斜距时,量测噪声协方差在每一次迭代构造,且为 δ r 2+ J ~diag{ δ p x 2, δ p y 2, δ p z 2} J ~ T。使用Singer模型描述飞机运动,取状态噪声协方差为Qk= Q ~ k 0 3 × 3 0 3 × 3 0 3 × 3 Q ~ k 0 3 × 3 0 3 × 3 0 3 × 3 Q ~ k,式中 Q ~ k=q T k 5 / 20 T k 4 / 8 T k 3 / 6 T k 4 / 8 T k 3 / 3 T k 2 / 2 T k 3 / 6 T k 2 / 2 T k,q为适合飞机机动大小的超参数,本文取为q=0.000 1。在上述设置下,同时考虑INS联合卫星定位的方法[13](记为INS/卫星)与本文方法的性能。
首先假设卫星信号在500 s~800 s的时间内丢失,考虑该情形下有无景象匹配点量测的滤波算法效果。图4显示了一次仿真两种方法的位置误差,可以看出,两种方法在0~500 s时间内误差的绝对值均小于5 m;在500 s~800 s时间内,由于卫星信号丢失,INS/卫星方法只能依靠INS信息进行导航,误差不断积累变大,而本文方法由于存在景象匹配点,很好地避免了该问题;在800 s以后,由于卫星信号恢复,INS/卫星方法迅速修正位置信息,误差减小。该结果表明,在卫星信号失效时,本文提出的方法仍能有效地获取导航信息,优于INS/卫星方法。
图4 两种方法的位置误差

Fig.4 The positioning errors of the two methods

接着考虑两种方法在存在卫星信号下的导航精度性能。使用均方根误差(Root-Mean-Square Error, RMSE)来表征测量精度,其计算方法为
RMSE= 1 M ∑ m = 1 M ( θ ^ m - θ ) 2
式中θ代表待估计参数的真实值, θ ^ m代表第m次独立实验对θ的估计值,M为独立实验次数,本文取M=500。图5显示了该情形下SAR平台位置估计的RMSE,可以看出,在x、y、z三个维度上,本文方法在大部分时间段具有更低的RMSE,优于INS/卫星方法,表明了景象匹配信息增加后,能够有效改善SAR平台的导航精度。
图5 两种方法的位置RMSE

Fig.5 The positioning RMSEs of the two methods

4 结束语

本文提出一种基于INS、卫星与景象匹配的SAR平台组合导航方法。该方法对不同的量测构造了对应的滤波流程,能实时提供导航信息。仿真实验表明,在卫星信号部分缺失的情况下,本文方法能获得比INS与卫星组合定位方法更小的位置误差;在包含卫星信号时,也能提供比INS与卫星组合定位方法更优的位置精度。仿真结果证明了本文方法的有效性与良好性能。

References

Publishing order | Descend order by publishing year | Descend order by cited within
[1]
SAVAGE P G. Strapdown inertial navigation integration algorithm design part 1: attitude algorithms[J]. Journal of Guidance, Control, and Dynamics, 1998, 21(1): 19-28.

DOI

[2]
TANIL C, KHANAFSEH S, JOERGER M, et al. An INS monitor to detect GNSS spoofers capable of tracking vehicle position[J]. IEEE Transactions on Aerospace and Electronic Systems, 2018, 54(1): 131-143.

DOI

[3]
GUO Y, WU M P, TANG K H, et al. Covert spoofing algorithm of UAV based on GPS/INS-integrated navigation[J]. IEEE Transactions on Vehicular Technology, 2019, 68(7): 6 557-6 564.

DOI

[4]
SOUMEKH M. Synthetic aperture radar signal processing with MATLAB algorithms[M]. New York: Wiley, 1999.

[5]
LI Y J, YU Y, PAN Q, et al. Scene matching based on spatial relation constraint in suitable-matching area[C]// 2009 IEEE International Conference on Intelligent Computing and Intelligent Systems. Shanghai, 2009: 598-603.

[6]
郭勤. 景象匹配技术发展概述[J]. 红外与激光工程, 2007, 36(S2): 57-61.

GUO Q. Scene matching technology development outline[J]. Infrared and Laser Engineering, 2007, 36(S2): 57-61.

[7]
陈圣义, 刘晓春, 雷志辉, 等. 结合景象匹配与惯导信息的SAR平台定位方法[J]. 国防科技大学学报, 2017, 39(4): 105-109.

CHEN S Y, LIU X C, LEI Z H, et al. SAR platform location method based on scene matching and inertial navigation system information[J]. Journal of National University of Defense Technology, 2017, 39(4): 105-109.

[8]
HUANG Z Y, YIN K Y. Generalized pseud-range equation for SAR platform positioning under the 6-D INS error model[J]. IETE Journal of Research, 2021: 1-8.

[9]
李家松, 李明磊, 魏大洲, 等. GPS拒止条件下的景象匹配导航方法研究[J]. 电子技术应用, 2022, 48(3): 88-93.

LI J S, LI M L, WEI D Z, et al. Research on scene matching navigation method under GPS-denied environments[J]. Application of Electronic Technique, 2022, 48(3): 88-93.

[10]
LIAO M Y, LIU J C, MENG Z Y, et al. A SINS/SAR/GPS fusion positioning system based on sensor credibility evaluations[J]. Remote Sensing, 2021, 13(21): 4 463.

DOI

[11]
GREWAL M S, ANDREWS A P. Kalman Filtering: Theory and Practice Using Matlab[M]. Hoboken, NJ, USA: John Wiley & Sons, Inc., 2014.

[12]
BRONSHTEIN I N, SEMENDYAYEV K A, MUSIOL G, et al. Geometry[M]//Handbook of Mathematics. Berlin, Heidelberg: Springer Berlin Heidelberg, 2015: 129-268.

[13]
张鑫, 冯双记. 机载INS/GNSS深组合导航系统发展现状[J]. 光学与光电技术, 2021, 19(6): 88-96.

ZHANG X, FENG S J. Development status of airborne INS/GNSS deep integrated navigation system[J]. Optics & Optoelectronic Technology, 2021, 19(6): 88-96.

Options
Outlines

/

〈 〉