CN116734864A - 一种常值观测偏差条件下航天器自主相对导航方法 - Google Patents
一种常值观测偏差条件下航天器自主相对导航方法 Download PDFInfo
- Publication number
- CN116734864A CN116734864A CN202311017608.0A CN202311017608A CN116734864A CN 116734864 A CN116734864 A CN 116734864A CN 202311017608 A CN202311017608 A CN 202311017608A CN 116734864 A CN116734864 A CN 116734864A
- Authority
- CN
- China
- Prior art keywords
- relative
- time
- measurement
- constant
- deviation
- 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.)
- Granted
Links
- 238000000034 method Methods 0.000 title claims abstract description 78
- 238000005259 measurement Methods 0.000 claims abstract description 111
- 230000033001 locomotion Effects 0.000 claims abstract description 55
- 238000004364 calculation method Methods 0.000 claims abstract description 48
- 238000001914 filtration Methods 0.000 claims abstract description 21
- 239000011159 matrix material Substances 0.000 claims description 55
- 239000013598 vector Substances 0.000 claims description 35
- 230000007704 transition Effects 0.000 claims description 20
- 230000001186 cumulative effect Effects 0.000 claims description 3
- 230000000694 effects Effects 0.000 abstract description 5
- 230000002401 inhibitory effect Effects 0.000 abstract 1
- 238000010586 diagram Methods 0.000 description 12
- 230000005764 inhibitory process Effects 0.000 description 4
- 238000004088 simulation Methods 0.000 description 4
- 238000001514 detection method Methods 0.000 description 3
- 230000010354 integration Effects 0.000 description 3
- 230000003287 optical effect Effects 0.000 description 3
- 230000006978 adaptation Effects 0.000 description 2
- 230000006870 function Effects 0.000 description 2
- 239000013307 optical fiber Substances 0.000 description 2
- 230000000644 propagated effect Effects 0.000 description 2
- 238000004590 computer program Methods 0.000 description 1
- 230000002708 enhancing effect Effects 0.000 description 1
- 230000005484 gravity Effects 0.000 description 1
- 230000003252 repetitive effect Effects 0.000 description 1
- 239000004065 semiconductor Substances 0.000 description 1
- 230000001629 suppression Effects 0.000 description 1
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01C—MEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
- G01C21/00—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00
- G01C21/20—Instruments for performing navigational calculations
-
- B—PERFORMING OPERATIONS; TRANSPORTING
- B64—AIRCRAFT; AVIATION; COSMONAUTICS
- B64G—COSMONAUTICS; VEHICLES OR EQUIPMENT THEREFOR
- B64G1/00—Cosmonautic vehicles
- B64G1/22—Parts of, or equipment specially adapted for fitting in or to, cosmonautic vehicles
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01C—MEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
- G01C21/00—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00
- G01C21/24—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00 specially adapted for cosmonautical navigation
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
- Y02T90/00—Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Remote Sensing (AREA)
- Radar, Positioning & Navigation (AREA)
- Mathematical Physics (AREA)
- Mathematical Optimization (AREA)
- Pure & Applied Mathematics (AREA)
- Computational Mathematics (AREA)
- Mathematical Analysis (AREA)
- Theoretical Computer Science (AREA)
- Data Mining & Analysis (AREA)
- Databases & Information Systems (AREA)
- Algebra (AREA)
- Automation & Control Theory (AREA)
- Software Systems (AREA)
- General Engineering & Computer Science (AREA)
- Operations Research (AREA)
- Aviation & Aerospace Engineering (AREA)
- Astronomy & Astrophysics (AREA)
- Computing Systems (AREA)
- Navigation (AREA)
Abstract
本公开关于一种常值观测偏差条件下航天器自主相对导航方法。该方法包括:获取从航天器的初始的相对运动状态;根据非线性轨道动力学模型计算从航天器的相对运动状态的数值积分;基于卡尔曼滤波计算模型以相对测量的常值观测偏差作为估计对象,并以相对测量的随机误差作为校正条件;通过非线性轨道动力学模型和卡尔曼滤波计算模型求得修正的相对运动状态;并进行迭代计算从而得到最终的相对运动状态。本公开通过采用非线性轨道动力学模型,获取高精度数值解。并通过闭环标校相对测量常值观测偏差并实时补偿相对测量的扩展卡尔曼滤波计算策略,实现相对运动状态的快速精确估计,使相对导航方法对相对测量的常值观测偏差具有快速检测和抑制效果。
Description
技术领域
本公开实施例涉及航天器自主导航技术领域,尤其涉及一种常值观测偏差条件下航天器自主相对导航方法。
背景技术
目前航天器自主相对导航技术主要是基于主航天器观测从航天器的相对距离和相对角度,经常采用相对运动的希尔(Hill-ClohessyWiltshire,HCW)轨道动力学模型,难以对相对测量的随机误差产生抑制效果。这主要受限于两个缺陷:其一,HCW轨道动力学模型本质为从航天器地心引力的一阶近似,在高精度导航滤波计算中存在可观性差、计算精度低的问题;其二,相对测量通常包含的常值观测偏差未能在相对导航方程中估计补偿,从而引起相对导航状态误差偏大。
关于上述技术方案,申请人发现存在以下问题,相对导航方法中存在的相对测量常值观测偏差问题,无法对相对测量的常值观测偏差进行快速检测,从而无法起到抑制效果。
因此,有必要改善上述相关技术方案中存在的一个或者多个问题。
需要说明的是,在上述背景技术部分公开的信息仅用于加强对本公开的背景的理解,因此可以包括不构成对本领域普通技术人员已知的现有技术的信息。
发明内容
本公开实施例的目的在于提供一种常值观测偏差条件下航天器自主相对导航方法,进而至少解决上述相关技术方案中存在的一个或者多个问题。
本发明的目的采用以下技术方案实现:
本发明提供了一种常值观测偏差条件下航天器自主相对导航方法,包括:
以主航天器中心做参考坐标系,获取从航天器的初始的相对运动状态;其中,所述相对运动状态包括相对位置矢量和相对速度矢量;
根据非线性轨道动力学模型计算所述从航天器的所述相对运动状态的数值积分;
基于卡尔曼滤波计算模型以相对测量的常值观测偏差作为估计对象,并以相对测量的随机误差作为校正条件;
通过所述非线性轨道动力学模型和卡尔曼滤波计算模型求得修正的相对运动状态;
通过对所述非线性轨道动力学模型和卡尔曼滤波计算模型进行迭代计算从而得到最终的相对运动状态。
可选地,所述基于卡尔曼滤波计算模型以相对测量的常值观测偏差作为估计对象,并以相对测量的随机误差作为校正条件的步骤,包括:
所述相对测量的常值观测偏差包括:相对位置误差、相对速度误差、相对距离的常值观测偏差、相对俯仰角的常值观测偏差和相对方位角的常值观测偏差。
可选地,所述基于卡尔曼滤波计算模型以相对测量的常值观测偏差作为估计对象,并以相对测量的随机误差作为校正条件的步骤,包括:
所述相对测量的随机误差包括:相对距离的随机误差、相对方位角的随机误差和相对俯仰角的随机误差。
可选地,所述根据非线性轨道动力学模型计算所述从航天器的所述相对运动状态的数值积分的步骤,包括:
所述非线性轨道动力学模型的方程表示为:
其中,k+1表示为当前时刻,k表示当前时刻的前一时刻;为相对运动状态,且/>为相对位置矢量,为相对速度矢量;/>为主航天器位置矢量,且/>为主航天器的轨道半径;/>为从航天器位置矢量;/>为地球引力常数;0*×*为元素都为0的矩阵;I表示除主对角线元素都为1外其余元素都为0的矩阵;
通过公式/>求得,其中,/>为主航天器的轨道运动角速率。
可选地,所述基于卡尔曼滤波计算模型以相对测量的常值观测偏差作为估计对象,并以相对测量的随机误差作为校正条件的步骤,包括:
根据所述相对运动状态的雅克比矩阵计算所述常值观测偏差的从前一时刻到当前时刻的状态转移矩阵;
基于所述状态转移矩阵确定从前一时刻到当前时刻的所述常值观测偏差的状态时间预测;
计算从前一时刻到当前时刻的所述状态转移矩阵的状态协方差时间预测,并基于所述状态协方差时间预测和所述相对测量的随机误差生成对应的增益矩阵;
计算当前时刻的相对测量的随机误差的估计值;
基于当前时刻的所述相对测量的随机误差的估计值、所述增益矩阵和所述状态时间预测更新当前时刻的相对测量的常值观测偏差;
更新当前时刻的状态协方差,累积当前时刻的相对测量的常值观测偏差以用于计算下一时刻的相对测量的随机误差的估计值;
将当前时刻的相对测量的常值观测偏差清零,并重新计算下一时刻的相对测量的常值观测偏差。
可选地,所述根据所述相对运动状态的雅克比矩阵计算所述常值观测偏差的从前一时刻到当前时刻的状态转移矩阵的步骤,包括:
所述相对运动状态的雅克比矩阵为:
其中,为相对运动状态偏导;
通过公式/>求得,其中,T通过公式/>求得;
从时刻到/>时刻的状态转移矩阵/>为:
其中,表示时间更新的周期;
从时刻到/>时刻的所述常值观测偏差的状态时间预测/>为:
其中,为/>时刻的相对测量的常值观测偏差。
可选地,所述计算从前一时刻到当前时刻的所述状态转移矩阵的状态协方差时间预测,并基于所述状态协方差时间预测生成对应的增益矩阵的步骤,包括:
从时刻到/>时刻的状态协方差时间预测/>为:
其中,为/>时刻的状态协方差;/>为/>时刻的过程噪声矩阵;
时刻的过程噪声矩阵/>为:
其中,通过公式/>求得,其中,/>表示原始过程噪声矩阵,/>表示过程噪声分布矩阵/>;
时刻的增益矩阵/>为:
其中,通过公式/>求得,其中,E[*]为求数学期望计算公式;表示 />时刻的相对测量的随机误差;
其中,通过公式/>求得,
其中,;为相对位置矢量/>在x分量;/>为相对位置矢量/>在y分量;/>为相对位置矢量/>在z分量;表示相对距离的测量值;A表示相对方位角的测量值;E表示相对俯仰角的测量值。
可选地,所述计算当前时刻的相对测量的随机误差的估计值的步骤,包括:
时刻的相对测量的随机误差的估计值/>为:
其中,相对距离的随机误差的估计值;/>相对方位角的随机误差的估计值;/>相对俯仰角的随机误差的估计值;/>为/>时刻的相对位置矢量,/>为的第1个元素,/>为/>的第2个元素,/>为/>的第3个元素;
其中,为/>时刻的相对距离的估计值,并通过公式/>求得,其中,/>为/>时刻的相对距离的测量值;/>为/>时刻的相对距离的常值观测偏差;
其中,为/>时刻的相对方位角的估计值,并通过公式/>求得,其中,/>为/>时刻的相对方位角的测量值;/>为/>时刻的相对方位角的常值观测偏差;
其中,为/>时刻的相对俯仰角的估计值,并通过公式/>求得,其中,/>为/>时刻的相对俯仰角的测量值;/>为/>时刻的相对俯仰角的常值观测偏差;
时刻的相对测量的常值观测偏差/>的更新方程为:
。
可选地,所述更新当前时刻的状态协方差的步骤,包括:
时刻的状态协方差/>的更新方程为:
。
可选地,所述累积当前时刻的相对测量的常值观测偏差以用于计算下一时刻的相对测量的随机误差的估计值的步骤,包括:
时刻的相对测量的常值观测偏差的累积方程为:
其中,为/>时刻的相对距离的常值观测偏差;/>为更新后/>时刻的相对测量的常值观测偏差/>的第7个元素;/>为/>时刻的相对方位角的常值观测偏差;/>为更新后/>时刻的相对测量的常值观测偏差/>的第8个元素;为/>时刻的相对俯仰角的常值观测偏差;/>为更新后/>时刻的相对测量的常值观测偏差/>的第9个元素。
本公开的实施方式中,通过采用非线性轨道动力学模型,获取高精度数值解,并通过闭环标校相对测量常值观测偏差并实时补偿相对测量的扩展卡尔曼滤波计算策略,实现相对运动状态的快速精确估计,使相对导航方法对相对测量的常值观测偏差具有快速检测和抑制效果。
应当理解的是,以上的一般描述和后文的细节描述仅是示例性和解释性的,并不能限制本公开。
附图说明
此处的附图被并入说明书中并构成本说明书的一部分,示出了符合本公开的实施例,并与说明书一起用于解释本公开的原理。显而易见地,下面描述中的附图仅仅是本公开的一些实施例,对于本领域普通技术人员来讲,在不付出创造性劳动的前提下,还可以根据这些附图获得其他的附图。
图1示出本公开示例性实施例中常值观测偏差条件下航天器自主相对导航方法的流程示意图;
图2示出本公开示例性实施例中基于卡尔曼滤波计算模型的流程示意图;
图3示出本公开示例性实施例中常值观测偏差条件下航天器自主相对导航方法的逻辑示意图;
图4示出本公开示例性实施例中从航天器相对主航天器相对运动轨迹的示意图;
图5示出本公开示例性实施例中1组相对导航中相对距离的仿真计算结果的示意图;
图6示出本公开示例性实施例中1组相对导航中相对速度的仿真计算结果的示意图;
图7示出本公开示例性实施例中1组相对导航中常值观测偏差的仿真计算结果的示意图;
图8示出本公开示例性实施例中100组相对导航中相对距离的仿真计算结果的示意图;
图9示出本公开示例性实施例中100组相对导航中相对速度的仿真计算结果的示意图;
图10示出本公开示例性实施例中100组相对导航中常值观测偏差的仿真计算结果的示意图;
图11示出本公开示例性实施例中一种存储介质的示意图。
具体实施方式
现在将参考附图更全面地描述示例实施方式。然而,示例实施方式能够以多种形式实施,且不应被理解为限于在此阐述的范例;相反,提供这些实施方式使得本公开将更加全面和完整,并将示例实施方式的构思全面地传达给本领域的技术人员。所描述的特征、结构或特性可以以任何合适的方式结合在一个或更多实施方式中。
此外,附图仅为本公开的示意性图解,并非一定是按比例绘制。图中相同的附图标记表示相同或类似的部分,因而将省略对它们的重复描述。附图中所示的一些方框图是功能实体,不一定必须与物理或逻辑上独立的实体相对应。可以采用软件形式来实现这些功能实体,或在一个或多个硬件模块或集成电路中实现这些功能实体,或在不同网络和/或处理器装置和/或微控制器装置中实现这些功能实体。
本示例实施方式中首先提供了一种常值观测偏差条件下航天器自主相对导航方法,参考图1中所示,该方法可以包括下述步骤:
步骤S101:以主航天器中心做参考坐标系,获取从航天器的初始的相对运动状态。
步骤S102:根据非线性轨道动力学模型计算从航天器的相对运动状态的数值积分。
步骤S103:基于卡尔曼滤波计算模型以相对测量的常值观测偏差作为估计对象,并以相对测量的随机误差作为校正条件。
步骤S104:通过非线性轨道动力学模型和卡尔曼滤波计算模型求得修正的相对运动状态。
步骤S105:通过对非线性轨道动力学模型和卡尔曼滤波计算模型进行迭代计算从而得到最终的相对运动状态。
其中,相对运动状态包括相对位置矢量和相对速度矢量。
需要理解的是,相对测量的常值观测偏差包括:相对位置误差、相对速度误差、相对距离的常值观测偏差、相对俯仰角的常值观测偏差和相对方位角的常值观测偏差。
还需要理解的是,相对测量的随机误差包括:相对距离的随机误差、相对方位角的随机误差和相对俯仰角的随机误差。
还需要理解的是,考虑相对测量常值观测偏差对相对运动状态的影响,建立了9个状态相对导航方程,不同于传统相对导航方程通常选择相对位置和相对速度为状态变量,本发明专利选择相对位置误差、相对速度误差和相对测量常值观测偏差为***状态,考虑了常值观测偏差对相对运动状态的影响;同时不同于传统相对导航方程通常选择相对测量为测量值,本发明专利选择相对测量的随机误差为测量值,从而使得相对导航滤波***状态方程和测量方程均表现为线性形式,具有计算简洁的特点。
还需要理解的是,在扩展卡尔曼滤波计算***状态转移矩阵和测量矩阵的框架下,设计了一种相对运动状态更新、相对测量常值观测偏差估计、相对测量的随机误差补偿同步进行的相对导航计算策略。在主从航天器相对运动完全非线性轨道动力学模型高精度数值积分约束条件下,该相对导航计算策略满足***状态无偏估计的可观性约束,确保相对导航状态和协方差均可收敛至其理论值。
还需要理解的是,提出的相对测量常值观测偏差条件下航天器自主相对导航方法,适用于近地航天器近程空间自主相对导航,不依赖于地面轨道计算支持,相对导航状态收敛精度高,计算速度快,相对导航方法对相对测量常值观测偏差具有一定的快速检测和抑制效果。
还需要理解的是,在步骤S101中,假定主航天器为圆轨道,预先可获取其轨道半径为rc。以地心惯性坐标系为基准建立当地垂直当地水平坐标系即LVLH(local-vertical/local-horizontal)坐标系。在主航天器LVLH坐标系中,近程范围内从航天器相对主航天器的相对位置和相对速度矢量为。x,y,z表示相对位置矢量,初始值是航天器直接测量得到的,通过本专利提出的方法获取更精确的相对位置矢量。相对速度矢量同理。主航天器对从航天器进行相对测量,包含的相对测距、方位角和俯仰角通常具有随机误差和常值偏差的误差特性。为此,建立9个状态的相对导航模型,采用相对测量常值偏差标校补偿和相对导航同步并行的计算策略,提高主航天器在相对测量常值观测偏差条件下自主相对导航的性能。
还需要理解的是,在步骤S104中,相对运动状态的修正可表示为:
(1)
也就相当于:
(2)
为通过非线性轨道动力学模型得到的相对运动状态;/>为通过卡尔曼滤波计算模型得到的求得相对位置误差和相对速度误差;/>为修正后的相对运动状态。
通过上述方法,通过采用非线性轨道动力学模型,获取高精度数值解,并通过闭环标校相对测量常值观测偏差并实时补偿相对测量的扩展卡尔曼滤波计算策略,实现相对运动状态的快速精确估计,使相对导航方法对相对测量的常值观测偏差具有快速检测和抑制效果。
下面,将参考图1至图6对本示例实施方式中的上述方法的各个步骤进行更详细的说明。
在一个实施例中,参考图3所示,步骤S102中还可以包括:
非线性轨道动力学模型的方程表示为:
(3)
(4)
其中,k+1表示为当前时刻,k表示当前时刻的前一时刻;为相对运动状态,且/>为相对位置矢量,/>为相对速度矢量;/>为主航天器位置矢量,且/>为主航天器的轨道半径;/>为从航天器位置矢量;/>为地球引力常数;0*×*为元素都为0的矩阵;I表示除主对角线元素都为1外其余元素都为0的矩阵;
通过公式/>求得,其中,/>为主航天器的轨道运动角速率。
需要理解的是,除了、/>为常数,公式中的其它变量均是时间t的函数。每个不同的时间t,对应不同的r。后文中的时间t省略了。
相对导航滤波计算中,***状态变量取为相对位置误差、相对速度误差以及相对距离、相对方位角和相对俯仰角的常值观测偏差,即:
这些误差均可以同时通过仿真设定得到。
在一个实施例中,参考图2所示,步骤S103中还可以包括:
步骤S201:根据相对运动状态的雅克比矩阵计算常值观测偏差的从前一时刻到当前时刻的状态转移矩阵。
步骤S202:基于状态转移矩阵确定从前一时刻到当前时刻的常值观测偏差的状态时间预测。
步骤S203:计算从前一时刻到当前时刻的状态转移矩阵的状态协方差时间预测,并基于状态协方差时间预测和所述相对测量的随机误差生成对应的增益矩阵。
步骤S204:计算当前时刻的相对测量的随机误差的估计值。
步骤S205:基于当前时刻的相对测量的随机误差的估计值、增益矩阵和状态时间预测更新当前时刻的相对测量的常值观测偏差。
步骤S206:更新当前时刻的状态协方差,累积当前时刻的相对测量的常值观测偏差以用于计算下一时刻的相对测量的随机误差的估计值。
步骤S207:将当前时刻的相对测量的常值观测偏差清零,并重新计算下一时刻的相对测量的常值观测偏差。
需要理解的是,在步骤S207中将当前时刻的相对测量的常值观测偏差清零的过程,可以通过公式进行表达,并代入到下一时刻中。
在一个实施例中,在步骤S201中:
相对运动状态的雅克比矩阵为:
(6)
其中,为相对运动状态偏导;
通过公式/>求得,其中,T通过公式/>求得。
需要理解的是,假设相对测量噪声为高斯白噪声,建立连续***状态方程表示:
(7)
式中,连续***状态矩阵表示为:
(8)
过程噪声分布矩阵表示为:
(9)
过程噪声矢量的均值和方差为:
(10)/>
在步骤S201中的从时刻到/>时刻的状态转移矩阵/>为:
(11)
其中,表示时间更新的周期。
需要理解的是,离散化状态方程表示为:
(12)
从中可推导出从时刻到/>时刻的状态转移矩阵/>,状态转移矩阵/>为取6阶截断近似表示的矩阵。
在一个实施例中,参考图3所示,在步骤S202中:
从时刻到/>时刻的常值观测偏差的状态时间预测/>为:
(13)
其中,为/>时刻的相对测量的常值观测偏差。
在一个实施例中,参考图3所示,在步骤S203中:
从时刻到/>时刻的状态协方差时间预测/>为:
(14)
其中,为/>时刻的状态协方差;/>为/>时刻的过程噪声矩阵;
时刻的过程噪声矩阵/>为:
(15)
其中,通过公式/>求得,其中,/>表示原始过程噪声矩阵,/>表示过程噪声分布矩阵/>。
取值如下:
(16)
时刻的增益矩阵/>为:
(17)
其中,通过公式/>求得,其中,E[*]为求数学期望计算公式;表示 />时刻的相对测量的随机误差;/>
其中,通过公式/>求得,
其中,;为相对位置矢量/>在x分量;/>为相对位置矢量/>在y分量;/>为相对位置矢量/>在z分量;表示相对距离的测量值;A表示相对方位角的测量值;E表示相对俯仰角的测量值。
需要理解的是,测量模型表示:
(18)
式中,,/>分别为相对距离的随机误差、相对方位角的随机误差和相对俯仰角的随机误差。
在一个实施例中,参考图3所示,在步骤S204中:
时刻的相对测量的随机误差的估计值/>为:/>
(19)
其中,相对距离的随机误差的估计值;/>相对方位角的随机误差的估计值;/>相对俯仰角的随机误差的估计值;/>为/>时刻的相对位置矢量,/>为的第1个元素,/>为/>的第2个元素,/>为/>的第3个元素;
其中,为/>时刻的相对距离的估计值,并通过公式/>求得,其中,/>为/>时刻的相对距离的测量值;/>为/>时刻的相对距离的常值观测偏差;
其中,为/>时刻的相对方位角的估计值,并通过公式/>求得,其中,/>为/>时刻的相对方位角的测量值;/>为/>时刻的相对方位角的常值观测偏差;
其中,为/>时刻的相对俯仰角的估计值,并通过公式/>求得,其中,/>为/>时刻的相对俯仰角的测量值;/>为/>时刻的相对俯仰角的常值观测偏差。
在一个实施例中,参考图3所示,在步骤S205中:
时刻的相对测量的常值观测偏差/>的更新方程为:
(20)
在一个实施例中,参考图3所示,步骤S206中:
时刻的状态协方差/>的更新方程为:
(21)
需要理解的是,表示初始协方差矩阵,在公式(16)中需要利用。当k=0时刻时,计算/>,即/>,需要利用/>,即/>。后续当k>0时刻,就不再利用/>了,并依次迭代计算。
在一个实施例中,参考图3所示,在步骤S206中:
时刻的相对测量的常值观测偏差的累积方程为:
(22)
其中,为/>时刻的相对距离的常值观测偏差;/>为更新后/>时刻的相对测量的常值观测偏差/>的第7个元素;/>为/>时刻的相对方位角的常值观测偏差;/>为更新后/>时刻的相对测量的常值观测偏差/>的第8个元素;为/>时刻的相对俯仰角的常值观测偏差;/>为更新后/>时刻的相对测量的常值观测偏差/>的第9个元素。/>表示相对距离测量常值观测偏差的剩余估计值,/>表示相对方位角测量的常值观测偏差的剩余估计值,/>表示相对俯仰角测量的常值观测偏差的剩余估计值。
本示例实施方式中根据上述实施方式进行综合,公开一个更加具体的实施方式,参考图3所示,具体步骤如下:
步骤1:相对导航滤波计算初始化。确定主航天器轨道半长轴,设定初始相对运动状态、初始协方差矩阵、过程噪声矩阵和测量噪声矩阵。
步骤2:从时刻数值积分相对运动非线性轨道动力学方程至/>时刻,即:
步骤3:计算时刻雅可比矩阵/>:
步骤4:计算从时刻到/>时刻状态转移矩阵/>:
步骤5:计算时刻过程噪声矩阵/>:/>
步骤6:状态时间预测:
步骤7:状态协方差时间预测:
步骤8:增益矩阵计算:
步骤9:相对测量的修正:
步骤10:测量信息计算:
步骤11:状态测量更新:
步骤12:状态协方差测量更新:
/>
步骤13:相对测量的常值观测偏差累积:
步骤14:相对运动状态校正:
步骤15:***状态清零:
步骤16:循环步骤2至步骤15,直至相对导航计算结束,完成相对测量常值观测条件下相对导航计算。
通过上述具体的实施方式,还提供了仿真模拟算例:
主从航天器初始轨道参数分别设定如下:
表1 主从航天器初始轨道参数(历元:2023-01-01 0:0:0.00)
根据表1中初始轨道参数,表示卫星轨道六根数,主航天器和从航天器空间运动的轨道动力学模型均采用高精度轨道预报(The High—precisionorbit propagator,HPOP)模型,则相对运动轨迹如图4所示。
其中,相对测量的常值观测偏差和随机误差/>分别设定如下:
相对导航滤波计算中,初始位置和速度误差分别设定如下:
相对导航滤波计算中,初始状态协方差矩阵设定如下:
参考图5、图6和图7所示,共同组成了1组相对导航仿真计算结果。参考图8、图9和图10所示,共同组成了100组蒙特卡洛相对导航仿真计算结果。根据图5至图10可知,对于主从航天器,基于HPOP动力学模型和9状态导航模型,考虑相对测距100m常值观测偏差和100m随机误差、相对方位角0.3度常值观测偏差和0.3度随机误差,相对俯仰角0.3度常值观测偏差和0.3度随机误差,采用同步估计相对测量常值观测偏差和相对运动状态的相对导航计算策略,则0.5个轨道周期后,相对导航位置误差小于30m,速度误差小于0.005m/s,相对测距和相对测角常值观测偏差的估计精度优于90%。
本示例实施方式中还提供了一种终端设备,包括:
一个或多个处理器;
存储器;
一个或多个应用程序,其中所述一个或多个应用程序被存储在所述存储器中并被配置为由所述一个或多个处理器执行,所述一个或多个应用程序配置用于:执行根据上述实施方式中任意一项所述的航天器自主相对导航方法。
其中终端设备执行操作的具体方式已经在有关航天器自主相对导航方法的实施方式中进行了详细描述,此处将不做详细阐述说明。
进一步的,本示例实施方式中,还提供了一种实现上述方法的虚拟装置。该虚拟装置可以包括分别实现上述方法中各个步骤的虚拟模块。
其中各个模块执行操作的具体方式已经在有关该方法的实施例中进行了详细描述,此处将不做详细阐述说明。
应当注意,尽管在上文详细描述中提及了用于动作执行的设备的若干模块或者单元,但是这种划分并非强制性的。实际上,根据本公开的实施方式,上文描述的两个或更多模块或者单元的特征和功能可以在一个模块或者单元中具体化。反之,上文描述的一个模块或者单元的特征和功能可以进一步划分为由多个模块或者单元来具体化。作为模块或单元显示的部件可以是或者也可以不是物理单元,即可以位于一个地方,或者也可以分布到多个网络单元上。可以根据实际的需要选择其中的部分或者全部模块来实现本公开方案的目的。本领域普通技术人员在不付出创造性劳动的情况下,即可以理解并实施。
在本公开的示例性实施例中,还提供了一种计算机可读存储介质,其上存储有计算机程序,该程序被处理器执行时可以实现上述任意一个实施方式中航天器自主相对导航方法的步骤。在一些可能的实施方式中,本发明的各个方面还可以实现为一种程序产品的形式,其包括程序代码,当所述程序产品在终端设备上运行时,所述程序代码用于使所述终端设备执行本说明书上述控制方法部分中描述的根据本发明各种示例性实施方式的步骤。
参考图11所示,描述了根据本发明的实施方式的用于实现上述方法的程序产品110,其可以采用便携式紧凑盘只读存储器(CD-ROM)并包括程序代码,并可以在终端设备,例如个人电脑上运行。然而,本发明的程序产品不限于此,在本文件中,可读存储介质可以是任何包含或存储程序的有形介质,该程序可以被指令执行***、装置或者器件使用或者与其结合使用。
所述程序产品可以采用一个或多个可读介质的任意组合。可读介质可以是可读信号介质或者可读存储介质。可读存储介质例如可以为但不限于电、磁、光、电磁、红外线或半导体的***、装置或器件,或者任意以上的组合。可读存储介质的更具体的例子(非穷举的列表)包括:具有一个或多个导线的电连接、便携式盘、硬盘、随机存取存储器(RAM)、只读存储器(ROM)、可擦式可编程只读存储器(EPROM或闪存)、光纤、便携式紧凑盘只读存储器(CD-ROM)、光存储器件、磁存储器件或者上述的任意合适的组合。
所述计算机可读存储介质可以包括在基带中或者作为载波一部分传播的数据信号,其中承载了可读程序代码。这种传播的数据信号可以采用多种形式,包括但不限于电磁信号、光信号或上述的任意合适的组合。可读存储介质还可以是可读存储介质以外的任何可读介质,该可读介质可以发送、传播或者传输用于由指令执行***、装置或者器件使用或者与其结合使用的程序。可读存储介质上包含的程序代码可以用任何适当的介质传输,包括但不限于无线、有线、光缆、RF等等,或者上述的任意合适的组合。
可以以一种或多种程序设计语言的任意组合来编写用于执行本发明操作的程序代码,所述程序设计语言包括面向对象的程序设计语言—诸如Java、C++等,还包括常规的过程式程序设计语言—诸如“C”语言或类似的程序设计语言。程序代码可以完全地在用户计算设备上执行、部分地在用户设备上执行、作为一个独立的软件包执行、部分在用户计算设备上部分在远程计算设备上执行或者完全在远程计算设备或服务器上执行。在涉及远程计算设备的情形中,远程计算设备可以通过任意种类的网络,包括局域网(LAN)或广域网(WAN),连接到用户计算设备,或者,可以连接到外部计算设备(例如利用因特网服务提供商来通过因特网连接)。
本领域技术人员在考虑说明书及实践这里公开的发明后,将容易想到本公开的其它实施方案。本申请旨在涵盖本公开的任何变型、用途或者适应性变化,这些变型、用途或者适应性变化遵循本公开的一般性原理并包括本公开未公开的本技术领域中的公知常识或惯用技术手段。说明书和实施例仅被视为示例性的,本公开的真正范围和精神由所附的权利要求指出。
Claims (10)
1.一种常值观测偏差条件下航天器自主相对导航方法,其特征在于,包括:
以主航天器中心做参考坐标系,获取从航天器的初始的相对运动状态;其中,所述相对运动状态包括相对位置矢量和相对速度矢量;
根据非线性轨道动力学模型计算所述从航天器的所述相对运动状态的数值积分;
基于卡尔曼滤波计算模型以相对测量的常值观测偏差作为估计对象,并以相对测量的随机误差作为校正条件;
通过所述非线性轨道动力学模型和卡尔曼滤波计算模型求得修正的相对运动状态;
通过对所述非线性轨道动力学模型和卡尔曼滤波计算模型进行迭代计算从而得到最终的相对运动状态。
2.根据权利要求1所述的航天器自主相对导航方法,其特征在于,所述基于卡尔曼滤波计算模型以相对测量的常值观测偏差作为估计对象,并以相对测量的随机误差作为校正条件的步骤,包括:
所述相对测量的常值观测偏差包括:相对位置误差、相对速度误差、相对距离的常值观测偏差、相对俯仰角的常值观测偏差和相对方位角的常值观测偏差。
3.根据权利要求2所述的航天器自主相对导航方法,其特征在于,所述基于卡尔曼滤波计算模型以相对测量的常值观测偏差作为估计对象,并以相对测量的随机误差作为校正条件的步骤,包括:
所述相对测量的随机误差包括:相对距离的随机误差、相对方位角的随机误差和相对俯仰角的随机误差。
4.根据权利要求1-3任一项的所述的航天器自主相对导航方法,其特征在于,所述根据非线性轨道动力学模型计算所述从航天器的所述相对运动状态的数值积分的步骤,包括:
所述非线性轨道动力学模型的方程表示为:
;
;
其中,k+1表示为当前时刻,k表示当前时刻的前一时刻;为相对运动状态,且/>为相对位置矢量,/>为相对速度矢量;为主航天器位置矢量,且/>为主航天器的轨道半径;/>为从航天器位置矢量;/>为地球引力常数;0*×*为元素都为0的矩阵;I表示除主对角线元素都为1外其余元素都为0的矩阵;
通过公式/>求得,其中,n为主航天器的轨道运动角速率。
5.根据权利要求4所述的航天器自主相对导航方法,其特征在于,所述基于卡尔曼滤波计算模型以相对测量的常值观测偏差作为估计对象,并以相对测量的随机误差作为校正条件的步骤,包括:
根据所述相对运动状态的雅克比矩阵计算所述常值观测偏差的从前一时刻到当前时刻的状态转移矩阵;
基于所述状态转移矩阵确定从前一时刻到当前时刻的所述常值观测偏差的状态时间预测;
计算从前一时刻到当前时刻的所述状态转移矩阵的状态协方差时间预测,并基于所述状态协方差时间预测和所述相对测量的随机误差生成对应的增益矩阵;
计算当前时刻的相对测量的随机误差的估计值;
基于当前时刻的所述相对测量的随机误差的估计值、所述增益矩阵和所述状态时间预测更新当前时刻的相对测量的常值观测偏差;
更新当前时刻的状态协方差,累积当前时刻的相对测量的常值观测偏差以用于计算下一时刻的相对测量的随机误差的估计值;
将当前时刻的相对测量的常值观测偏差清零,并重新计算下一时刻的相对测量的常值观测偏差。
6.根据权利要求5所述的航天器自主相对导航方法,其特征在于,所述根据所述相对运动状态的雅克比矩阵计算所述常值观测偏差的从前一时刻到当前时刻的状态转移矩阵的步骤,包括:
所述相对运动状态的雅克比矩阵为:
;
其中,为相对运动状态偏导;
通过公式/>;求得,其中,T通过公式/>求得;
从时刻到/>时刻的状态转移矩阵/>为:
;
其中,表示时间更新的周期;
从时刻到/>时刻的所述常值观测偏差的状态时间预测/>为:
;
其中,为/>时刻的相对测量的常值观测偏差。
7.根据权利要求6所述的航天器自主相对导航方法,其特征在于,所述计算从前一时刻到当前时刻的所述状态转移矩阵的状态协方差时间预测,并基于所述状态协方差时间预测生成对应的增益矩阵的步骤,包括:
从时刻到/>时刻的状态协方差时间预测/>为:
;
其中,为/>时刻的状态协方差;/>为/>时刻的过程噪声矩阵;
时刻的过程噪声矩阵/>为:
;
其中,通过公式/>求得,其中,/>表示原始过程噪声矩阵,/>表示过程噪声分布矩阵/>;
时刻的增益矩阵/>为:
;
其中,通过公式/>求得,其中,E[*]为求数学期望计算公式;/>表示/>时刻的相对测量的随机误差;
其中,通过公式/>求得,
其中,;/>为相对位置矢量/>在x分量;/>为相对位置矢量/>在y分量;/>为相对位置矢量/>在z分量;/>表示相对距离的测量值;A表示相对方位角的测量值;E表示相对俯仰角的测量值。
8.根据权利要求7所述的航天器自主相对导航方法,其特征在于,所述计算当前时刻的相对测量的随机误差的估计值的步骤,包括:
时刻的相对测量的随机误差的估计值/>为:
;
其中,相对距离的随机误差的估计值;/>相对方位角的随机误差的估计值;相对俯仰角的随机误差的估计值;/>为/>时刻的相对位置矢量,/>为/>的第1个元素,/>为/>的第2个元素,/>为/>的第3个元素;
其中,为/>时刻的相对距离的估计值,并通过公式/>求得,其中,/>为/>时刻的相对距离的测量值;/>为/>时刻的相对距离的常值观测偏差;
其中,为/>时刻的相对方位角的估计值,并通过公式/>求得,其中,/>为/>时刻的相对方位角的测量值;/>为/>时刻的相对方位角的常值观测偏差;
其中,为/>时刻的相对俯仰角的估计值,并通过公式/>求得,其中,/>为/>时刻的相对俯仰角的测量值;/>为/>时刻的相对俯仰角的常值观测偏差;
时刻的相对测量的常值观测偏差/>的更新方程为:
。
9.根据权利要求8所述的航天器自主相对导航方法,其特征在于,所述更新当前时刻的状态协方差的步骤,包括:
时刻的状态协方差/>的更新方程为:
。
10.根据权利要求8所述的航天器自主相对导航方法,其特征在于,所述累积当前时刻的相对测量的常值观测偏差以用于计算下一时刻的相对测量的随机误差的估计值的步骤,包括:
时刻的相对测量的常值观测偏差的累积方程为:
;
其中,为/>时刻的相对距离的常值观测偏差;/>为更新后/>时刻的相对测量的常值观测偏差/>的第7个元素;/>为/>时刻的相对方位角的常值观测偏差;/>为更新后/>时刻的相对测量的常值观测偏差/>的第8个元素;/>为/>时刻的相对俯仰角的常值观测偏差;/>为更新后/>时刻的相对测量的常值观测偏差/>的第9个元素。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202311017608.0A CN116734864B (zh) | 2023-08-14 | 2023-08-14 | 一种常值观测偏差条件下航天器自主相对导航方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202311017608.0A CN116734864B (zh) | 2023-08-14 | 2023-08-14 | 一种常值观测偏差条件下航天器自主相对导航方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN116734864A true CN116734864A (zh) | 2023-09-12 |
CN116734864B CN116734864B (zh) | 2023-11-28 |
Family
ID=87904758
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202311017608.0A Active CN116734864B (zh) | 2023-08-14 | 2023-08-14 | 一种常值观测偏差条件下航天器自主相对导航方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN116734864B (zh) |
Cited By (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN117367436A (zh) * | 2023-12-08 | 2024-01-09 | 中国西安卫星测控中心 | 一种星间相对测量线性时变误差的实时估计方法 |
Citations (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN108692729A (zh) * | 2018-05-04 | 2018-10-23 | 北京空间飞行器总体设计部 | 一种空间非合作目标相对导航协方差自适应修正滤波方法 |
CN108717198A (zh) * | 2018-05-04 | 2018-10-30 | 北京空间飞行器总体设计部 | 一种空间非合作目标相对导航***误差补偿修正滤波方法 |
CN109631913A (zh) * | 2019-01-30 | 2019-04-16 | 西安电子科技大学 | 基于非线性预测强跟踪无迹卡尔曼滤波的x射线脉冲星导航定位方法及*** |
WO2020087846A1 (zh) * | 2018-10-31 | 2020-05-07 | 东南大学 | 基于迭代扩展卡尔曼滤波融合惯性与单目视觉的导航方法 |
CN112325886A (zh) * | 2020-11-02 | 2021-02-05 | 北京航空航天大学 | 一种基于重力梯度仪和陀螺仪组合的航天器自主定姿*** |
US20210048542A1 (en) * | 2019-08-16 | 2021-02-18 | California Institute Of Technology | Systems and Methods for Robust and Accurate Relative Navigation |
CN115856977A (zh) * | 2022-12-26 | 2023-03-28 | 上海航天控制技术研究所 | 一种基于差分gnss的相对导航方法 |
-
2023
- 2023-08-14 CN CN202311017608.0A patent/CN116734864B/zh active Active
Patent Citations (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN108692729A (zh) * | 2018-05-04 | 2018-10-23 | 北京空间飞行器总体设计部 | 一种空间非合作目标相对导航协方差自适应修正滤波方法 |
CN108717198A (zh) * | 2018-05-04 | 2018-10-30 | 北京空间飞行器总体设计部 | 一种空间非合作目标相对导航***误差补偿修正滤波方法 |
WO2020087846A1 (zh) * | 2018-10-31 | 2020-05-07 | 东南大学 | 基于迭代扩展卡尔曼滤波融合惯性与单目视觉的导航方法 |
CN109631913A (zh) * | 2019-01-30 | 2019-04-16 | 西安电子科技大学 | 基于非线性预测强跟踪无迹卡尔曼滤波的x射线脉冲星导航定位方法及*** |
US20210048542A1 (en) * | 2019-08-16 | 2021-02-18 | California Institute Of Technology | Systems and Methods for Robust and Accurate Relative Navigation |
CN112325886A (zh) * | 2020-11-02 | 2021-02-05 | 北京航空航天大学 | 一种基于重力梯度仪和陀螺仪组合的航天器自主定姿*** |
CN115856977A (zh) * | 2022-12-26 | 2023-03-28 | 上海航天控制技术研究所 | 一种基于差分gnss的相对导航方法 |
Non-Patent Citations (2)
Title |
---|
刘勇;徐世杰;徐鹏;马骏;: "基于改进型两步卡尔曼滤波的相对导航方法", 中国空间科学技术, no. 03 * |
李轶;张善从;: "基于视线测量的航天器相对导航滤波方法研究", 仪器仪表学报, no. 06 * |
Cited By (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN117367436A (zh) * | 2023-12-08 | 2024-01-09 | 中国西安卫星测控中心 | 一种星间相对测量线性时变误差的实时估计方法 |
CN117367436B (zh) * | 2023-12-08 | 2024-02-23 | 中国西安卫星测控中心 | 一种星间相对测量线性时变误差的实时估计方法 |
Also Published As
Publication number | Publication date |
---|---|
CN116734864B (zh) | 2023-11-28 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
RU2701194C2 (ru) | Способ оценки навигационного состояния в условиях ограниченной возможности наблюдения | |
CN108225337B (zh) | 基于sr-ukf滤波的星敏感器和陀螺组合定姿方法 | |
CN109163721A (zh) | 姿态测量方法及终端设备 | |
CN106772524B (zh) | 一种基于秩滤波的农业机器人组合导航信息融合方法 | |
Gross et al. | A comparison of extended Kalman filter, sigma-point Kalman filter, and particle filter in GPS/INS sensor fusion | |
CN116734864B (zh) | 一种常值观测偏差条件下航天器自主相对导航方法 | |
WO2008120145A1 (en) | Method and system for orientation sensing | |
Zhang et al. | Robust H‐infinity CKF/KF hybrid filtering method for SINS alignment | |
CN112146655A (zh) | 一种BeiDou/SINS紧组合导航***弹性模型设计方法 | |
Xiong et al. | Adaptive iterated extended Kalman filter for relative spacecraft attitude and position estimation | |
CN108508463B (zh) | 基于Fourier-Hermite正交多项式扩展椭球集员滤波方法 | |
Zhang et al. | Attitude determination using gyros and vector measurements aided with adaptive kinematics modeling | |
CN114323007A (zh) | 一种载体运动状态估计方法及装置 | |
Wang et al. | A robust backtracking CKF based on Krein space theory for in-motion alignment process | |
CN114139109A (zh) | 一种目标跟踪方法、***、设备、介质及数据处理终端 | |
CN112782732A (zh) | 一种基于粒子群算法的导航信号解析方法及计算机存可读介质 | |
CN116878510A (zh) | 基于因子图优化的多源融合室内定位方法、装置及*** | |
Avzayesh et al. | Improved-Performance Vehicle’s State Estimator Under Uncertain Model Dynam | |
Wang et al. | A line-of-sight rate estimation method for roll-pitch gimballed infrared seeker | |
Fiedler et al. | A probabilistic moving horizon estimation framework applied to the visual-inertial sensor fusion problem | |
CN114705223A (zh) | 多移动智能体在目标跟踪中的惯导误差补偿方法及*** | |
Girrbach et al. | Towards in-field and online calibration of inertial navigation systems using moving horizon estimation | |
Hao et al. | Research on data fusion for SINS/GPS/magnetometer integrated navigation based on modified CDKF | |
KR101860810B1 (ko) | 비선형 측정모델을 가지는 항법시스템 및 방법 | |
Girrbach et al. | Adaptive compensation of measurement delays in multi-sensor fusion for inertial motion tracking using moving horizon estimation |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |