CN106441291B - 一种基于强跟踪sdre滤波的组合导航***及导航方法 - Google Patents

一种基于强跟踪sdre滤波的组合导航***及导航方法 Download PDF

Info

Publication number
CN106441291B
CN106441291B CN201610855565.7A CN201610855565A CN106441291B CN 106441291 B CN106441291 B CN 106441291B CN 201610855565 A CN201610855565 A CN 201610855565A CN 106441291 B CN106441291 B CN 106441291B
Authority
CN
China
Prior art keywords
state
matrix
space
time
module
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201610855565.7A
Other languages
English (en)
Other versions
CN106441291A (zh
Inventor
赵良玉
任珊珊
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Beijing Institute of Technology BIT
Original Assignee
Beijing Institute of Technology BIT
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Beijing Institute of Technology BIT filed Critical Beijing Institute of Technology BIT
Priority to CN201610855565.7A priority Critical patent/CN106441291B/zh
Publication of CN106441291A publication Critical patent/CN106441291A/zh
Application granted granted Critical
Publication of CN106441291B publication Critical patent/CN106441291B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01CMEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
    • G01C21/00Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00
    • G01C21/10Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00 by using measurements of speed or acceleration
    • G01C21/12Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00 by using measurements of speed or acceleration executed aboard the object being navigated; Dead reckoning
    • G01C21/16Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00 by using measurements of speed or acceleration executed aboard the object being navigated; Dead reckoning by integrating acceleration or speed, i.e. inertial navigation
    • G01C21/165Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00 by using measurements of speed or acceleration executed aboard the object being navigated; Dead reckoning by integrating acceleration or speed, i.e. inertial navigation combined with non-inertial navigation instruments
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S19/00Satellite radio beacon positioning systems; Determining position, velocity or attitude using signals transmitted by such systems
    • G01S19/38Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system
    • G01S19/39Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system the satellite radio beacon positioning system transmitting time-stamped messages, e.g. GPS [Global Positioning System], GLONASS [Global Orbiting Navigation Satellite System] or GALILEO
    • G01S19/42Determining position
    • G01S19/48Determining position by combining or switching between position solutions derived from the satellite radio beacon positioning system and position solutions derived from a further system
    • G01S19/49Determining position by combining or switching between position solutions derived from the satellite radio beacon positioning system and position solutions derived from a further system whereby the further system is an inertial position system, e.g. loosely-coupled

Landscapes

  • Engineering & Computer Science (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Automation & Control Theory (AREA)
  • Position Fixing By Use Of Radio Waves (AREA)
  • Navigation (AREA)

Abstract

本发明公开了一种基于强跟踪SDRE滤波的组合导航***及导航方法,其中,所述***包括GPS接收机、惯性导航子***、和滤波器,其中,所述滤波器包括空间模块、SDRE分解模块、离散化处理模块和强跟踪滤波模块,并且,所述强跟踪滤波模块又包括初始化模块、一步预测子模块、增益矩阵获得子模块、实际状态估计子模块和误差方差阵更新子模块。所述方法包括以下步骤:先进行SDRE分解,然后进行离散化处理,最后进行强跟踪滤波,得到组合导航数据,并进行数据输出。本发明所述的***和方法引入了自适应渐消因子,能够进行准确滤波,得到精确的导航数据。

Description

一种基于强跟踪SDRE滤波的组合导航***及导航方法
技术领域
本发明涉及导航领域,尤其涉及组合导航***及导航方法,特别地,涉及一种基于强跟踪SDRE滤波的组合导航***及导航方法。
背景技术
基于GPS技术和惯性导航的组合导航已被广泛应用于军民诸多领域,该技术对GPS定位数据和惯性导航定位数据进行估计、修正、融合,从而达到优势互补以提高载体导航定位精度的目的。在GPS和惯性导航组合导航中以导航参数作为估计对象,可直观描述导航参数的动态变化过程,准确反映真实状态的演变情况。
在组合导航***中,实际***总是存在不同程度的非线性。同时,实际***还存在高斯或非高斯的随机噪声干扰。在这种情况下,对非线性***进行有效甚至最优状态估计十分重要,需要采用非线性滤波技术。在面对非线性***时,一般采用扩展卡尔曼滤波(EKF)进行处理。
但这种方法具有高阶截断误差,并且必须求解***的Jacobian矩阵,同时还存在着强非线性***下失去滤波稳定性的可能。滤波技术的基础是建立准确的模型,如果***模型不精确,会影响滤波精度。
发明内容
为了克服上述问题,本发明人进行了锐意研究,使用GPS接收机、惯性导航(INS)子***(简称惯导***)组成组合导航***,其次使用SDRE技术处理滤波模型,用残差序列处处正交的原理设计滤波器,设计出一种基于强跟踪SDRE滤波的组合导航***及导航方法,这样,可以得到高精度的位置和速度信息,实现精确的组合导航,从而完成本发明。
本发明一方面提供了一种基于强跟踪SDRE滤波的组合导航***,其中,所述***包括GPS接收机1、惯性导航子***2、和滤波器3,其中,所述滤波器3包括空间模块31、SDRE分解模块32、离散化处理模块33、和强跟踪滤波模块34,并且,所述强跟踪滤波模块34包括初始化模块340、一步预测子模块341、增益矩阵获得子模块342、实际状态估计子模块343和误差方差阵更新子模块344。
本发明另一方面提供了一种基于强跟踪SDRE滤波的导航方法,优选利用本发明第一方面所述***进行,其中,所述方法包括以下步骤:
步骤1、利用SDRE分解模块32对空间模块31内的状态空间和量测空间进行状态相关分解,分别得到状态分解空间和量测分解空间;
步骤2、利用离散化处理模块33对SDRE分解模块32内的状态分解空间和量测分解空间进行离散化处理,分别得到状态离散空间和量测离散空间;
步骤3、利用强跟踪滤波模块34对离散化处理模块33内的状态离散空间和量测离散空间进行强跟踪滤波处理,得到组合导航数据,并进行数据输出。
附图说明
图1示出本发明所述组合导航***的结构示意图;
图2示出本发明所述组合导航方法的流程图;
图3示出在实验例中飞行器所做的运动轨迹;
图4a示出利用实施例和对比例所述方法得到的在东向的位置误差;
图4b示出利用实施例和对比例所述方法得到的在北向的位置误差;
图4c示出利用实施例和对比例所述方法得到的在天向的位置误差;
图5a示出利用实施例和对比例所述方法得到的在东向的速度误差;
图5b示出利用实施例和对比例所述方法得到的在北向的速度误差;
图5c示出利用实施例和对比例所述方法得到的在天向的速度误差。
附图标记
1-GPS接收机;2-惯性导航子***;21-加速度计;22-陀螺仪;3-滤波器;31-空间模块;32-SDRE分解模块;33-离散化处理模块;34-强跟踪滤波模块;340-初始化模块;341-一步预测子模块;342-增益矩阵获得子模块;343-实际状态估计子模块;344-误差方差阵更新子模块。
具体实施方式
下面通过附图对本发明进一步详细说明。通过这些说明,本发明的特点和优点将变得更为清楚明确。
本发明一方面提供了一种基于强跟踪SDRE滤波的组合导航***,所述导航***包括GPS接收机1、惯性导航子***2和滤波器3。其中,所述GPS接收机1用于输出通过GPS技术得到的当前地理位置和当前速度信息;所述惯性导航子***2包括加速度计21和陀螺仪22,用于检测并输出通过惯性导航技术得到的当前地理位置、当前速度和当前姿态角信息;所述滤波器3用于对GPS接收机1和惯性导航子***2输出的数据进行数据处理,得到组合导航信息并进行输出。
根据本发明一种优选的实施方式,所述滤波器3包括空间模块31、SDRE分解模块32、离散化处理模块33和强跟踪滤波模块34。
在进一步优选的实施方式中,在所述空间模块31内集成有状态空间和量测空间。
在更进一步优选的实施方式中,所述状态空间和量测空间分别如式(31-1)和式(31-2)所示:
z=h(x)+V 式(31-2)
其中,在式(31-1)中:
其中,φe、φn和φu分别为数学平台误差角,实时进行更新;ve、vn和vu分别为载体的东向、北向、天向的速度;λ和h分别为载体所在的纬度,经度和高度;εe、εn和εu分别为陀螺仪常值漂移,Δe、Δn和Δu分别为加速度计常值漂移;
其中,分别表示陀螺仪在东向、北向和天向的随机误差,分别表示加速度计在东向、北向和天向的随机误差;其中,所述随机误差即为***噪声;
G为***噪声转移矩阵,其中,
其中,I表示单位矩阵。
根据本发明一种优选的实施方式,如式(31-2)所示的状态空间如下建立:以地理系(东、北、天)为导航解算的基本坐标系,考虑飞行高度h,并假设地球为椭球体,惯性导航(INS)状态方程如下:
(1)指北惯性导航***姿态误差角方程:
其中,在式(31-1-1)~式(31-1-3)中,RM是当地子午面内主曲率半径;RN是卯酉面内主曲率半径;ωie是地球自转角速率;λ和h分别为惯导***(INS)测得的载体所在的纬度、经度和高度;φe、φn和φe为数学平台误差角;ve、vn和vu分别为惯导***(INS)测得的载体的东向、北向和天向的速度;εe、εn和εu分别为陀螺仪在东向、北向和天向的常值漂移;δve为惯导***(INS)所测的东向速度与GPS接收机所测的东向速度之差,δvn为惯导***(INS)所测的东向速度与GPS接收机所测的北向速度之差;为惯导***(INS)所测的纬度与GPS接收机所测的纬度之差,δh为惯导***(INS)所测的高度与GPS接收机所测的高度之差。其中,地球长半径Re=6378245,地球扁率e=1/298.257。
(2)速度方程为:
其中,在式(31-1-4)~式(31-1-6)中,fe、fn和fu分别为加速度计的测量值在东向、北向和天向三个方向的加速度分量;Δe、Δn和Δu为加速度计在东向、北向和天向的常值漂移。
(3)位置方程:
(4)惯性元件误差方程:假设惯性元件的误差是常值噪声,则其误差方程为:
在进一步优选的实施方式中,将上述(1)所述的指北惯性导航***姿态误差角方程、上述(2)所述的速度方程、上述(3)所述的位置方程以及上述(4)所述的惯性元件误差方程联立起来,取15阶状态变量:即可得到如式(31-1)所示的状态空间,其中,f(x)为对应的式(31-1-1)~式(31-1-9)所示的9个基本导航参数方程和式(31-1-10)与式(31-1-11)所示的惯性元件误差方程。
根据本发明一种优选的实施方式,在式(31-2)中:
量测信息z=[vge,vgn,vgu,pge,pgn,pgu]T,其中,vge表示GPS接收机测得的载体在东向的速度,vgn表示GPS接收机测得的载体在北向的速度,vgu表示GPS接收机测得的载体在天向的速度,pge表示GPS接收机测得的载体在东向的位置,pgn表示GPS接收机测得的载体在北向的位置,pgu表示GPS接收机测得的载体在天向的位置;
V为量测噪声,V=[Vvge,Vvgn,Vvgu,Vpge,Vpgn,Vpgu]T,其中,Vvge表示GPS接收机测得的载体在东向速度的随机误差,Vvgn表示GPS接收机测得的载体在北向速度的随机误差,Vvgu表示GPS接收机测得的载体在天向速度的随机误差,Ppge表示GPS接收机测得的载体在东向位置的随机误差,Ppgn表示GPS接收机测得的载体在北向位置的随机误差,Ppgu表示GPS接收机测得的载体在天向位置的随机误差。
在进一步优选的实施方式中,h(x)是线性方程,可以写成h(x)=Hx,其中,H为量测矩阵,其中:
根据本发明一种优选的实施方式,所述SDRE分解模块32用于对空间模块31内如式(31-1)和式(31-2)所示的状态空间和量测空间进行状态相关分解,分别得到状态分解空间和量测分解空间。
在进一步优选的实施方式中,所述状态分解空间和量测分解空间分别如式(32-1)和式(32-2)所示:
z=Hx+V 式(32-2)。
其中,SDRE是用状态相关分解方法将非线性状态空间(式(31-1))分解成状态与转移矩阵相乘的形式,即f(x)=F(x)x。其中,式(31-2)所示的量测空间中,h(x)是线性方程,因此不需要SDRE分解,因此,SDRE分解是对式(31-1)所述状态空间进行SDRE分解,得到式(32-1)。在本发明中,式(32-2)虽然不是经式(31-2)SDRE分解得到,但是,为了统一表述,将式(32-2)表述为量测分解空间。
在本发明中,状态转移矩阵F(x)为:
其中,FN是9维矩阵,其中的非零元素为:
F(4,2)=-fu F(4,3)=fn
F(5,1)=fu F(5,3)=-fe
F(6,1)=-fn F(6,2)=fe
F(9,6)=1;
Fs为:
根据本发明一种优选的实施方式,所述离散化处理模块33用于对SDRE分解模块32内的状态分解空间和量测分解空间进行离散化处理,得到状态离散空间和量测离散空间。
在进一步优选的实施方式中,所述状态离散空间和量测离散空间分别如式(33-1)和式(33-2)所示:
xk+1=Φk+1/kxkkwk 式(33-1);
zk+1=Θk+1xk+1+Vk+1 式(33-2)。
其中,Φk+1/k是离散后的状态转移矩阵,Гk是离散后的噪声转移矩阵,Θk+1是离散后的量测矩阵,其中,离散化后:
根据泰勒展开,Φk+1/k=I+f(xk)Δt+f(xk)Δt2/2!+f(xk)Δt2/3!+…,
Γk=GΔt+GΔt2/2!+GΔt2/3!+…。
在更进一步优选的实施方式中,***噪声向量wk和量测噪声向量Vk的方差阵分别为Qk和Rk,并且wk和Vk满足:
E[wk]=0,表示wk的均值为0;
其中,Cov[wk,wj]表示wk和wj的协方差,用于衡量wk和wj的总体误差;
E[Vk]=0,表示Vk的均值为0;
Cov[Vk,Vj]表示Vk和Vj的协方差,用于衡量Vk和Vj的总体误差;
Cov[wk,Vj]表示wk和Vj的协方差,用于衡量wk和Vj的总体误差。
其中,wk表示在k时刻的***噪声向量,wj表示j时刻的***噪声向量,Vk表示k时刻的量测噪声向量,Vj表示j时刻的量测噪声向量,δkj为Kronecker-δ函数,表示不同时刻噪声相互独立。
根据本发明一种优选的实施方式,所述强跟踪滤波模块34用于对离散化处理模块33内的如式(33-1)所示的状态离散空间和如式(33-2)所示的量测离散空间进行强跟踪滤波处理,得到组合导航信息。
在进一步优选的实施方式中,强跟踪滤波模块34包括一步预测子模块341、增益矩阵获得子模块342、实际状态估计子模块343和误差方差阵更新子模块344。
根据本发明一种优选的实施方式,所述一步预测子模块用于进行状态的一步预测和误差方差阵的一步预测。
其中,所述状态的一步预测是指根据当前时刻(k时刻)的状态信息对下一时刻(k+1时刻)的状态信息进行预测,所述误差方差阵的一步预测是指根据当前时刻的误差方阵对下一时刻的误差方差阵进行预测。
在进一步优选的实施方式中,所述状态的一步预测如式(341-1)所示:
其中,在式(341-1)中,表示根据k时刻的状态预测k+1时刻的状态预测值,xk表示k时刻的状态,Δt表示滤波周期。
在更进一步优选的实施方式中,所述误差方差阵的一步预测如式(341-2)所示:
其中,在式(341-2)中,Pk+1/k表示根据k时刻的误差方差阵预测得到的k+1时刻的误差方差阵预测值,βk为自适应渐消因子,Pk表示k时刻的误差方差阵,Φk+1/k表示离散后的状态转移矩阵,表示Φk+1/k的转置矩阵,Qk表示wk的方差阵,即***噪声方差阵,k表示当前时刻,k+1表示下一时刻,其中,k=0,1,2,3,…,L,其中,L由滤波时间和滤波周期决定,L=滤波时间/滤波周期。
在本发明中,为了增强滤波器对模型不确定的鲁棒性,在进行误差方差阵的一步预测时,引入自适应渐消因子βk,提高了一步预测状态的不确定度。其主要是基于新息残差正交化原理,即其中,新息是实际量测值与一步预测的量测值之差如果一步预测值准确,那么,各时刻的新息是正交的,即新息是均值为0的高斯白噪声,一旦不准确,则不再正交。
根据本发明一种优选的实施方式,所述自适应渐消因子通过下式(a)~式(e)进行调整,(可参见:周东华,席裕庚张钟俊.非线性***带次优渐消因子的扩展卡尔曼滤波[J].控制与决策,1990(05):1-6.):
其中,在式(a)中,表示新息序列,表示的转置矩阵,ρ表示遗忘因子、为人为设定、且ρ=0~1,Vk-1表示k-1时刻的量测噪声;
其中,在式(b)中,Nk没有实际意义,其用来表示表示的,其中,Θk+1表示离散后的量测矩阵,表示Θk+1的转置矩阵,Qk表示wk的方差阵,Rk+1表示量测噪声向量Vk+1的方差阵;
其中,在式(c)中,Mk没有实际意义,其用来表示Θk+1表示离散后的量测矩阵,表示Θk+1的转置矩阵,Φk+1/k表示离散后的状态转移矩阵,表示Φk+1/k的转置矩阵,Pk表示k时刻的误差方差阵;
其中,在式(d)中,tr(Nk)表示Nk的迹,tr(Mk)表示Mk的迹,所述迹即为矩阵的迹,即矩阵对角线上各元素的和;
其中,通过上式进行调整得到的自适应渐消因子能够提高滤波的准确性,并且得到准确的状态值。具体地,获得自适应渐消因子βk的最终目的是在一步预测状态值与实际状态值值相差较大的情况下,增大滤波增益矩阵,从而加大对量测量的信任,得到更精确的状态估计值。
根据本发明一种优选的实施方式,所述增益矩阵获得子模块342用于获得k+1时刻的增益矩阵,如式(342-1)所示:
其中,在式(342-1)中,Pk+1/k表示根据k时刻的误差方差阵预测得到的k+1时刻的误差方差阵,Θk+1是离散后的量测矩阵,为Θk+1的转置矩阵,Rk+1为量测噪声向量Vk+1的方差阵。
根据本发明一种优选的实施方式,所述实际状态估计子模块343用于估计k+1时刻的实际状态值,即用于获得k+1时刻的实际状态估计值,具体如式(343-1)所示:
其中,表示根据k时刻的状态预测的k+1时刻的预测状态值,Kk+1表示k+1时刻的增益矩阵,zk+1为k+1时刻的量测数据(根据GPS给出的k+1时刻的信息获得),表示一步预测状态值的量测值。
根据本发明一种优选的实施方式,所述误差方差阵更新子模块344用于对k+1时刻的误差方差阵进行更新,其中,利用k+1时刻的增益矩阵Kk+1和根据k时刻的误差方差阵预测得到的k+1时刻的误差方差阵Pk+1/k获得k+1时刻的误差方差阵Pk+1
在进一步优选的实施方式中,k+1时刻的误差方差阵的更新如式(344-1)进行:
Pk+1=(I-Kk+1Θk+1)Pk+1/k 式(344-1)
其中,为了令滤波器不断地运行下去,直至***过程结束,必须进行下一时刻(k+1时刻)的误差方差阵的更新。
根据本发明一种优选的实施方式,本发明所述滤波器还包括初始化模块340,用于对状态和误差方差阵进行初始化,同时进行时间更新。
在进一步优选的实施方式中,所述初始化模块340赋予状态初始值x0,即惯导***提供的初始值。
在更进一步优选的实施方式中,所述初始化模块340赋予初始误差方差阵P0,其表示惯导***的不确定性。其中,根据误差统计值,预估初始误差方差阵,并且,为了使滤波更快到达准确位置,可以将P0设置为误差统计值的10倍,其中所述误差统计值为统计出来的传感器的误差值。
本发明另一方面还提供了一种基于强跟踪SDRE滤波的组合导航方法。
根据本发明一种优选的实施方式,如图2所示,所述导航方法包括以下步骤:
步骤1、利用SDRE分解模块32对空间模块31内的状态空间和量测空间进行状态相关分解,分别得到状态分解空间和量测分解空间。
根据本发明一种优选的实施方式,在步骤1中,所述状态空间和量测空间分别如式(31-1)和式(31-2)所示:
z=h(x)+V 式(31-2)。
根据本发明一种优选的实施方式,在步骤1中,SDRE分解模块32进行状态相关分解后得到的状态分解空间和量测分解空间分别如式(32-1)和式(32-2)所示:
z=Hx+V 式(32-2)。
其中,SDRE分解是用状态相关分解方法将非线性状态空间(式(31-1))分解成状态与转移矩阵相乘的线性***形式,即f(x)=F(x)x。其中,在式(31-2)所示的量测空间中,h(x)是线性方程,因此不需要SDRE分解,因此,SDRE分解是对式(31-1)所述状态空间进行SDRE分解,得到式(32-1)。在本发明中,式(32-2)虽然不是经式(31-2)SDRE分解得到,但是,为了统一表述,将式(32-2)表述为量测分解空间。
根据本发明一种优选的实施方式,在步骤1之前进行步骤1’:
步骤1’、利用初始化模块340对状态和误差方差阵进行初始化并进行时间更新。
在进一步优选的实施方式中,初始化模块340对状态进行初始化,得到状态初始值x0和初始误差方差阵P0,即分别为惯导***提供的初始值与惯导***的不确定性。
步骤2、利用离散化处理模块33对SDRE分解模块32内的状态分解空间和量测分解空间进行离散化处理,分别得到状态离散空间和量测离散空间。
根据本发明一种优选的实施方式,经离散化处理模块33处理得到的状态离散空间和量测离散空间分别如式(33-1)和式(33-2)所示:
xk+1=Φk+1/kxkkwk 式(33-1);
zk+1=Θk+1xk+1+Vk+1 式(33-2)。
其中,Φk+1/k是离散后的状态转移矩阵,Гk是离散后的噪声转移矩阵,Θk+1是离散后的量测矩阵。
在进一步优选的实施方式中,***噪声向量wk和量测噪声向量Vk的方差阵分别为Qk和Rk,并且wk和Vk满足:
E[wk]=0,E[Vk]=0,
其中,wk表示在k时刻的***噪声向量,wj表示j时刻的***噪声向量,Vk表示k时刻的量测噪声向量,Vj表示j时刻的量测噪声向量,δkj为Kronecker-δ函数,表示不同时刻噪声相互独立。
步骤3、利用强跟踪滤波模块34对离散化处理模块33内的状态离散空间和量测离散空间进行强跟踪滤波处理,得到组合导航数据,并进行数据输出。
根据本发明一种优选的实施方式,在步骤3中,所述组合导航数据包括:φe(k+1)、φn(k+1)、φe(k+1)、Δe(k+1)、Δn(k+1)、Δu(k+1)、ve(k+1)、vn(k+1)、vu(k+1)、λ(k+1)和h(k+1)。
在进一步优选的实施方式中,如图2所示,所述步骤3包括以下子步骤:
步骤3-1、利用一步预测子模块341对下一时刻(k+1时刻)的状态和误差方阵进行一步预测,获得下一时刻的状态预测值和误差方差阵预测值,分别以表示。
步骤3-2、利用步骤3-1获得的Pk+1/k,经增益矩阵获得子模块342处理,得到下一时刻的增益矩阵,以Kk+1表示;
步骤3-3、利用步骤3-1获得的和步骤3-2获得的Kk+1,经实际状态估计子模块343获得下一时刻的实际状态值,以xk+1表示,即组合导航信息,并进行信息输出。
根据本发明一种优选的实施方式,在步骤3-1中,所述状态的一步预测是指根据当前时刻(k时刻)的状态信息对下一时刻(k+1时刻)的状态信息进行预测,其如式(341-1)进行:
其中,在式(341-1)中,表示根据k时刻的状态预测的k+1时刻的状态预测值,xk表示k时刻的状态,Δt表示滤波周期。
在进一步优选的实施方式中,所述误差方差阵的一步预测是指根据当前时刻的误差方阵对下一时刻的误差方差阵进行预测,其如式(341-2)进行:
其中,在式(341-2)中,Pk+1/k表示根据k时刻的误差方差阵预测得到的k+1时刻的误差方差阵预测值,βk为自适应渐消因子,Pk表示k时刻的误差方差阵,Φk+1/k表示离散后的状态转移矩阵,表示Φk+1/k的转置矩阵,Qk表示wk的方差阵,k表示当前时刻,k+1表示下一时刻,其中,k=0,1,2,3,…,L,其中,L由滤波时间和滤波周期决定,L=滤波时间/滤波周期。
根据本发明一种优选的实施方式,在步骤3-2中,利用步骤3-1获得的Pk+1/k经增益矩阵获得子模块342如式(342-1)所示的处理,得到k+1时刻的增益矩阵:
其中,在式(342-1)中,Pk+1/k表示根据k时刻的误差方差阵预测得到的k+1时刻的误差方差阵,Θk+1是离散后的量测矩阵,为Θk+1的转置矩阵,Rk+1为量测噪声向量Vk+1的方差阵。
根据本发明一种优选的实施方式,在步骤3-3中,利用步骤3-1获得的以及步骤3-2获得的Pk+1/k在实际状态估计子模块343内进行如式(343-1)所示处理,得到k+1时刻的实际状态估计值xk+1
其中,表示根据k时刻的状态预测的k+1时刻的预测状态,Kk+1表示k+1时刻的增益矩阵,zk+1为k+1时刻的量测数据(根据GPS给出的k+1时刻的信息获得),表示一步预测状态值的量测值。
根据本发明一种优选的实施方式,为了令滤波器不断地运行下去,直至***过程结束,在进行步骤3-3的同时必须进行下一时刻(k+1时刻)的误差方差阵的更新,即与步骤3-3同时进行步骤3-4:
步骤3-4、利用步骤3-1获得的Pk+1/k以及步骤3-2获得的Kk+1经误差方差阵更新子模块344进行如式(344-1)所示进行处理:
Pk+1=(I-Kk+1Θk+1)Pk+1/k 式(344-1)
其中,所述误差方差阵更新子模块344用于对k+1时刻的误差方差阵进行更新,其中,利用k+1时刻的增益矩阵Kk+1和根据k时刻的误差方差阵预测得到的k+1时刻的误差方差阵Pk+1/k获得k+1时刻的误差方差阵Pk+1
本发明的有益效果:
(1)本发明所述***结构简单,运行速度快;
(2)在本发明所述***中引入自适应渐消因子,经过***处理可以得到高精度的状态值;
(3)本发明所述方法采用SDRE分解,无需进行线性化处理,也不用求解Jacobian矩阵,减少了处理量,对最终工程应用具有实际意义;
(4)本发明所述方法提高了对***不确定情况下的鲁棒性以及对飞行状态的突变情况下的跟踪能力,改善了导航精度,提高了导航***的可靠性;
(5)本发明所述方法简单易行,并且结果准确。
实施例
在本实施例中,设定飞行器做如图3所示的机动运动,飞行轨迹中含有爬升、下降、转弯和变速运动。惯导***(INS)的采样周期为0.01s,GPS接收机的采样周期及滤波周期为1s,仿真模拟时间为1932s。
根据上述所述对空间模块31内的状态空间和量测空间依次进行SDRE分解、离散化处理,然后进行强滤波跟踪,具体过程如下:
(1)在零时刻,通过惯导***测量得到x0,同时经过处理得到初始误差方差阵P0
其中,根据误差统计值,预估初始误差方差阵,并且,为了使滤波更快到达准确位置,可以将P0设置为误差统计值的10倍,其中所述误差统计值为统计出来的传感器的误差值。
(2)利用一步预测子模块341对第一时刻的状态和误差方阵进行一步预测,分别获得第一时刻的状态预测值和误差方差阵预测值P1/0
其中:
其中:x0由步骤(1)得到;Δt为设定值,在本实施例中设为0.01s;f(x0)=F(x0)x0,其中,x0由步骤(1)得到,F(x0)是将零时刻的惯导***输出的参数值输入状态转移矩阵F(x)中得到;
其中,β0为自适应渐消因子,Φ1/0=I+F(x0)T,其中,在F(x0)内输入零时刻的惯导信息,Q0为人为设置。
(3)利用步骤(2)获得的P1/0,经增益矩阵获得子模块342处理,得到第一时刻的增益矩阵K1
其中,其中,P1/0经上述(2)获得,Θ1=H,R1为人为设置。
(4)利用步骤(2)获得的和步骤(3)获得的K1,经实际状态估计子模块343获得下一时刻的实际状态值x1
其中,其中,由步骤(2)获得,K1由步骤(3)获得,z1表示在第一时刻的GPS检测数据,
(5)利用步骤(3)获得的K1、步骤(2)获得的P1/0经误差方差阵更新子模块344处理,得到第一时刻的误差方差阵P1
其中,P1=(I-K1Θ1)P1/0,其中,步骤(5)的目的是为了令滤波器不断地运行下去,并且,步骤(5)与步骤(4)同时进行。
将获得的x1与P1,以及惯性导航通过惯导解算后的第一时刻的信息输入滤波器,进行数据更新,重复步骤(2)~(5)的过程。
对比例
与实施例不同于:采用SDRE滤波进行处理,而不是强跟踪SDRE滤波,其在步骤3中的自适应渐消因子β为1。
实验例
分别采用实施例所述方法以及对比例所述方法对飞行器进行仿真,同时,假设现有空间模型是精确空间模型,在500s~800s内保持滤波模型不变,飞行数据中引入状态突变量Δx=[01×3,1,1,1,10/Re,10/Re,10,01×6]T,使这段时间的真实空间为:
xk+1=f(xk)+Δx
对比例和实施例与真实空间在东、北、天的位置与真实位置之间的位置误差结果分别如图4a~图4c所示,在东、北、天的速度与真实速度之间的速度误差结果如图5a~图5c所示。
其中:在图4a~图4c中可以看出,在500s~800s的时间段内,实施例在这一时间段内滤波后的位置误差明显小于对比例,保持在5m之内,实施例的滤波结果更接近于真实状态;在图5a~图5c中可以看出,在出现误差的时间段内,实施例的误差也均小于对比例,小于1m/s,再次证明实施例的滤波结果更接近于真实状态;
因此,在平稳飞行时间内,采用本发明所述方法,速度误差小于1m/s,位置误差小于5m。
以上结合了优选的实施方式对本发明进行了说明,不过这些实施方式仅是范例性的,仅起到说明性的作用。在此基础上,可以对本发明进行多种替换和改进,这些均落入本发明的保护范围内。

Claims (7)

1.一种基于强跟踪SDRE滤波的组合导航***,其特征在于,所述***包括:
GPS接收机(1),用于输出通过GPS技术得到的当前地理位置和当前速度信息;
惯性导航子***(2),其包括加速度计(21)和陀螺仪(22),用于检测并输出通过惯性导航技术得到的当前地理位置、当前速度和当前姿态角信息;和
滤波器(3),用于对GPS接收机(1)和惯性导航子***(2)输出的数据进行数据处理,得到组合导航信息并进行输出;
所述滤波器(3)包括
空间模块(31),在其内集成有状态空间和量测空间;
SDRE分解模块(32),用于对空间模块(31)内的状态空间和量测空间进行状态相关分解,分别得到状态分解空间和量测分解空间;
离散化处理模块(33),用于对SDRE分解模块(32)内的状态分解空间和量测分解空间进行离散化处理,分别得到状态离散空间和量测离散空间;和
强跟踪滤波模块(34),用于对离散化处理模块(33)内的状态离散空间和量测离散空间进行强跟踪滤波处理,得到组合导航信息;
所述状态空间和量测空间分别如式(31-1)和式(31-2)所示:
z=h(x)+V 式(31-2);
其中,在式(31-1)中:
其中,φe、φn和φu分别为数学平台误差角,实时进行更新;ve、vn和vu分别为载体的东向、北向、天向的速度;λ和h分别为载体所在的纬度,经度和高度;εe、εn和εu分别为陀螺仪常值漂移,Δe、Δn和Δu分别为加速度计常值漂移;
其中,分别表示陀螺仪在东向、北向和天向的随机误差,分别表示加速度计在东向、北向和天向的随机误差;其中,所述随机误差即为***噪声;
G为***噪声转移矩阵,其中,I表示单位矩阵;
在式(31-2)中:
量测信息z=[vge,vgn,vgu,pge,pgn,pgu]T,其中,vge表示GPS接收机测得的载体在东向的速度,vgn表示GPS接收机测得的载体在北向的速度,vgu表示GPS接收机测得的载体在天向的速度,pge表示GPS接收机测得的载体在东向的位置,pgn表示GPS接收机测得的载体在北向的位置,pgu表示GPS接收机测得的载体在天向的位置;
V为量测噪声,V=[Vvge,Vvgn,Vvgu,Vpge,Vpgn,Vpgu]T,其中,Vvge表示GPS接收机测得的载体在东向速度的随机误差,Vvgn表示GPS接收机测得的载体在北向速度的随机误差,Vvgu表示GPS接收机测得的载体在天向速度的随机误差,Ppge表示GPS接收机测得的载体在东向位置的随机误差,Ppgn表示GPS接收机测得的载体在北向位置的随机误差,Ppgu表示GPS接收机测得的载体在天向位置的随机误差;
h(x)是线性方程,h(x)=Hx,其中,H为量测矩阵,其中:
所述状态分解空间和量测分解空间分别如式(32-1)和式(32-2)所示:
z=Hx+V 式(32-2);
其中,在式(32-1)中,F(x)为状态转移矩阵,其中,FN是9维矩阵,其中的非零元素为:
F(4,2)=-fu F(4,3)=fn
F(5,1)=fu F(5,3)=-fe
F(6,1)=-fn F(6,2)=fe
F(9,6)=1;
其中,ωie是地球自转角速率;和h分别为惯导***测得的载体所在的纬度和高度;ve、vn分别为惯导***测得的载体的东向、北向的速度;δve为惯导***所测的东向速度与GPS接收机所测的东向速度之差,δvn为惯导***所测的东向速度与GPS接收机所测的北向速度之差;为惯导***所测的纬度与GPS接收机所测的纬度之差,δh为惯导***所测的高度与GPS接收机所测的高度之差;fe、fn和fu分别为加速度计的测量值在东向、北向和天向三个方向的加速度分量;RM是当地子午面内主曲率半径;RN是卯酉面内主曲率半径;其中,地球长半径Re=6378245,地球扁率e=1/298.257;
Fs为:
H为量测矩阵,其中:
所述状态离散空间和量测离散空间分别如式(33-1)和式(33-2)所示:
xk+1=Φk+1/kxkkwk 式(33-1),
zk+1=Θk+1xk+1+Vk+1 式(33-2);
其中,Φk+1/k是离散后的状态转移矩阵,Гk是离散后的噪声转移矩阵,Θk+1是离散后的量测矩阵;wk表示在k时刻的***噪声向量;xk表示k时刻的状态;xk+1表示k+1时刻的实际状态估计值;zk+1为k+1时刻的量测数据;
所述强跟踪滤波模块(34)包括一步预测子模块(341)、增益矩阵获得子模块(342)、实际状态估计子模块(343)、误差方差阵更新子模块(344);
所述状态的一步预测如式(341-1)所示进行处理:
其中,在式(341-1)中,表示根据k时刻的状态预测k+1时刻的状态预测值,xk表示k时刻的状态,Δt表示滤波周期;
所述误差方差阵的一步预测如式(341-2)所示进行处理:
其中,在式(341-2)中,Pk+1/k表示根据k时刻的误差方差阵预测得到的k+1时刻的误差方差阵的预测值,βk为自适应渐消因子,Pk表示k时刻的误差方差阵,Φk+1/k表示离散后的状态转移矩阵,表示Φk+1/k的转置矩阵,Qk表示wk的***噪声方差阵,k表示当前时刻,k+1表示下一时刻,其中,k=0,1,2,3,…,L,其中,L由滤波时间和滤波周期决定,L=滤波时间/滤波周期;
所述实际状态估计值xk+1的获得如式(343-1)所示:
其中,表示根据k时刻的状态预测的k+1时刻的预测状态值,Kk+1表示k+1时刻的增益矩阵,zk+1为k+1时刻的量测数据、根据GPS给出的k+1时刻的信息获得,表示根据一步预测状态值的量测值。
2.根据权利要求1所述的***,其特征在于,所述强跟踪滤波模块(34)包括
一步预测子模块(341),用于进行状态的一步预测和误差方差阵的一步预测,分别得到和Pk+1/k
增益矩阵获得子模块(342),用于利用Pk+1/k获得k+1时刻的增益矩阵,以Kk+1表示;
实际状态估计子模块(343),用于利用和Kk+1获得k+1时刻的实际状态估计值,以xk+1表示;
误差方差阵更新子模块(344),利用Kk+1和Pk+1/k对k+1时刻的误差方差阵进行更新,即获得Pk+1
3.根据权利要求1所述的***,其特征在于,
所述增益矩阵Kk+1的获得如式(342-1)所示:
其中,在式(342-1)中,Pk+1/k表示根据k时刻的误差方差阵的预测得到的k+1时刻的误差方差阵,Θk+1是离散后的量测矩阵,为Θk+1的转置矩阵,Rk+1为量测噪声向量Vk+1的量测噪声方差阵;
和/或
k+1时刻的误差方差阵的更新如式(344-1)进行:
Pk+1=(I-Kk+1Θk+1)Pk+1/k 式(344-1)。
4.根据权利要求1所述的***,其特征在于,所述滤波器还包括初始化模块(340),用于对状态和误差方差阵进行初始化,分别得到状态初始值和初始误差方差阵,分别以x0和P0表示。
5.一种基于强跟踪SDRE滤波的组合导航方法,所述方法利用权利要求1至4所述的组合导航***进行,其中,所述方法包括以下步骤:
步骤1、利用SDRE分解模块(32)对空间模块(31)内的状态空间和量测空间进行状态相关分解,分别得到状态分解空间和量测分解空间;
步骤2、利用离散化处理模块(33)对SDRE分解模块(32)内的状态分解空间和量测分解空间进行离散化处理,分别得到状态离散空间和量测离散空间;
步骤3、利用强跟踪滤波模块(34)对离散化处理模块(33)内的状态离散空间和量测离散空间进行强跟踪滤波处理,得到组合导航数据,并进行数据输出。
6.根据权利要求5所述的方法,其中,
所述步骤3包括以下子步骤:
步骤3-1、利用一步预测子模块(341)对下一时刻、k+1时刻的状态和误差方阵进行一步预测,获得下一时刻的状态预测值和误差方差阵预测值,分别以和Pk+1/k表示;
步骤3-2、利用步骤3-1获得的Pk+1/k,经增益矩阵获得子模块(342)处理,得到下一时刻的增益矩阵,以Kk+1表示;
步骤3-3、利用步骤3-1获得的和步骤3-2获得的Kk+1,经实际状态估计子模块(343)获得下一时刻的实际状态值,以xk+1表示,即组合导航信息,并进行信息输出。
7.根据权利要求6所述的方法,其中,在进行步骤3-3的同时必须进行下一时刻、k+1时刻的误差方差阵的更新,即与步骤3-3同时进行步骤3-4:
步骤3-4、利用步骤3-1获得的Pk+1/k以及步骤3-2获得的Kk+1经误差方差阵更新子模块(344)进行如式(344-1)所示进行处理:
Pk+1=(I-Kk+1Θk+1)Pk+1/k 式(344-1),
其中,所述误差方差阵更新子模块(344)用于对k+1时刻的误差方差阵进行更新,其中,利用k+1时刻的增益矩阵Kk+1和根据k时刻的误差方差阵预测得到的k+1时刻的误差方差阵Pk+1/k获得k+1时刻的误差方差阵Pk+1
CN201610855565.7A 2016-09-27 2016-09-27 一种基于强跟踪sdre滤波的组合导航***及导航方法 Active CN106441291B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201610855565.7A CN106441291B (zh) 2016-09-27 2016-09-27 一种基于强跟踪sdre滤波的组合导航***及导航方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201610855565.7A CN106441291B (zh) 2016-09-27 2016-09-27 一种基于强跟踪sdre滤波的组合导航***及导航方法

Publications (2)

Publication Number Publication Date
CN106441291A CN106441291A (zh) 2017-02-22
CN106441291B true CN106441291B (zh) 2019-06-21

Family

ID=58169569

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201610855565.7A Active CN106441291B (zh) 2016-09-27 2016-09-27 一种基于强跟踪sdre滤波的组合导航***及导航方法

Country Status (1)

Country Link
CN (1) CN106441291B (zh)

Families Citing this family (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106931966B (zh) * 2017-02-24 2019-07-26 西北工业大学 一种基于泰勒高阶余项拟合的组合导航方法
CN108168574B (zh) * 2017-11-23 2022-02-11 东南大学 一种基于速度观测的8位置捷联惯导***级标定方法
CN109000642A (zh) * 2018-05-25 2018-12-14 哈尔滨工程大学 一种改进的强跟踪容积卡尔曼滤波组合导航方法
CN111366151A (zh) * 2018-12-26 2020-07-03 北京信息科技大学 一种极地范围船舶导航的信息融合方法
CN112082548B (zh) * 2020-09-10 2022-04-26 中国人民解放军海军航空大学 一种无人机惯导与gps混合高度测量方法
CN113630106A (zh) * 2021-08-02 2021-11-09 杭州电子科技大学 基于强跟踪滤波的高阶扩展卡尔曼滤波器设计方法

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104062672A (zh) * 2013-11-28 2014-09-24 哈尔滨工程大学 基于强跟踪自适应Kalman滤波的SINSGPS组合导航方法
CN105737828A (zh) * 2016-05-09 2016-07-06 郑州航空工业管理学院 一种基于强跟踪的相关熵扩展卡尔曼滤波的组合导航方法

Family Cites Families (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102608631B (zh) * 2011-10-28 2014-09-03 北京航空航天大学 基于模糊逻辑的自适应强跟踪ukf定位滤波算法
PT3009860T (pt) * 2014-10-16 2020-03-23 Gmv Aerospace And Defence S A Método para calcular um limite de erro de um filtro de kalman baseado em solução de posição gnss

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104062672A (zh) * 2013-11-28 2014-09-24 哈尔滨工程大学 基于强跟踪自适应Kalman滤波的SINSGPS组合导航方法
CN105737828A (zh) * 2016-05-09 2016-07-06 郑州航空工业管理学院 一种基于强跟踪的相关熵扩展卡尔曼滤波的组合导航方法

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
Robust INS/GPS Sensor Fusion for UAV Localization Using SDRE Nonlinear Filtering;Abdelkrim Nemra,et al.;《IEEE SENSORS JOURNAL》;20100430;第10卷(第4期);789-798
基于强跟踪SDRE滤波的GPS_INS组合导航;任珊珊等;《弹箭与制导学报》;20171031;第37卷(第5期);43-46,55
非线性***带次优渐消因子的扩展卡尔曼滤波;周东华等;《控制与决策》;19901231(第5期);1-6

Also Published As

Publication number Publication date
CN106441291A (zh) 2017-02-22

Similar Documents

Publication Publication Date Title
CN106441291B (zh) 一种基于强跟踪sdre滤波的组合导航***及导航方法
Jung et al. Inertial attitude and position reference system development for a small UAV
Marantos et al. UAV state estimation using adaptive complementary filters
Kingston et al. Real-time attitude and position estimation for small UAVs using low-cost sensors
Madyastha et al. Extended Kalman filter vs. error state Kalman filter for aircraft attitude estimation
CN104655131B (zh) 基于istssrckf的惯性导航初始对准方法
Youn et al. Combined quaternion-based error state Kalman filtering and smooth variable structure filtering for robust attitude estimation
CN109931955B (zh) 基于状态相关李群滤波的捷联惯性导航***初始对准方法
CN103941273B (zh) 机载惯性/卫星组合导航***的自适应滤波方法与滤波器
No et al. Attitude estimation method for small UAV under accelerative environment
Oh Multisensor fusion for autonomous UAV navigation based on the Unscented Kalman Filter with Sequential Measurement Updates
CN105737823A (zh) 基于五阶ckf的gps/sins/cns组合导航方法
CN105806363A (zh) 基于srqkf的sins/dvl水下大失准角对准方法
Hu et al. Robust unscented Kalman filter-based decentralized multisensor information fusion for INS/GNSS/CNS integration in hypersonic vehicle navigation
de La Parra et al. Low cost navigation system for UAV's
JP5164645B2 (ja) カルマンフィルタ処理における繰り返し演算制御方法及び装置
CN104713559A (zh) 一种高精度sins模拟器的设计方法
CN111189442A (zh) 基于cepf的无人机多源导航信息状态预测方法
Liu et al. Interacting multiple model UAV navigation algorithm based on a robust cubature Kalman filter
Liu et al. An improved GNSS/INS navigation method based on cubature Kalman filter for occluded environment
Davari et al. Multirate adaptive Kalman filter for marine integrated navigation system
CN110873577B (zh) 一种水下快速动基座对准方法及装置
Eldesoky et al. Improved position estimation of real time integrated low-cost navigation system using unscented kalman filter
Xia et al. Low-cost MEMS-INS/GNSS integration using quaternion-based nonlinear filtering methods for USV
Khoder et al. A quaternion scaled unscented kalman estimator for inertial navigation states determination using ins/gps/magnetometer fusion

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant
EE01 Entry into force of recordation of patent licensing contract

Application publication date: 20170222

Assignee: BEIJING BOYING TONGHANG TECHNOLOGY Co.,Ltd.

Assignor: BEIJING INSTITUTE OF TECHNOLOGY

Contract record no.: X2023110000104

Denomination of invention: A Combined Navigation System and Navigation Method Based on Strong Tracking SDRE Filter

Granted publication date: 20190621

License type: Common License

Record date: 20230904

EE01 Entry into force of recordation of patent licensing contract