CN111090830B - 一种高轨非合作目标在轨光压辨识方法 - Google Patents
一种高轨非合作目标在轨光压辨识方法 Download PDFInfo
- Publication number
- CN111090830B CN111090830B CN201911247650.5A CN201911247650A CN111090830B CN 111090830 B CN111090830 B CN 111090830B CN 201911247650 A CN201911247650 A CN 201911247650A CN 111090830 B CN111090830 B CN 111090830B
- Authority
- CN
- China
- Prior art keywords
- orbit
- solar
- satellite
- relative
- light pressure
- 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
Links
- 238000000034 method Methods 0.000 title claims abstract description 16
- 230000033001 locomotion Effects 0.000 claims abstract description 57
- 230000015572 biosynthetic process Effects 0.000 claims abstract description 11
- 230000000694 effects Effects 0.000 claims abstract description 5
- 230000001133 acceleration Effects 0.000 claims description 20
- 230000014509 gene expression Effects 0.000 claims description 10
- 238000004364 calculation method Methods 0.000 claims description 8
- 238000005259 measurement Methods 0.000 claims description 8
- 230000009471 action Effects 0.000 claims description 5
- 230000008859 change Effects 0.000 claims description 4
- 239000011159 matrix material Substances 0.000 claims description 4
- 230000007704 transition Effects 0.000 claims description 2
- 230000003287 optical effect Effects 0.000 abstract description 2
- 238000009825 accumulation Methods 0.000 description 2
- 238000004422 calculation algorithm Methods 0.000 description 2
- 238000001914 filtration Methods 0.000 description 2
- 230000005855 radiation Effects 0.000 description 2
- 238000006243 chemical reaction Methods 0.000 description 1
- 230000000295 complement effect Effects 0.000 description 1
- 238000010276 construction Methods 0.000 description 1
- 230000010006 flight Effects 0.000 description 1
- 238000005286 illumination Methods 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 230000008569 process Effects 0.000 description 1
- 238000004088 simulation Methods 0.000 description 1
- 238000006467 substitution reaction Methods 0.000 description 1
- 238000012795 verification Methods 0.000 description 1
Images
Classifications
-
- 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
- G06F17/12—Simultaneous equations, e.g. systems of linear equations
-
- 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)
- Mathematical Physics (AREA)
- Pure & Applied Mathematics (AREA)
- Computational Mathematics (AREA)
- Mathematical Optimization (AREA)
- Mathematical Analysis (AREA)
- Theoretical Computer Science (AREA)
- Data Mining & Analysis (AREA)
- Algebra (AREA)
- Databases & Information Systems (AREA)
- Software Systems (AREA)
- General Engineering & Computer Science (AREA)
- Computing Systems (AREA)
- Operations Research (AREA)
- Control Of Position, Course, Altitude, Or Attitude Of Moving Bodies (AREA)
- Navigation (AREA)
Abstract
本发明公开了一种高轨非合作目标在轨光压辨识方法,首先,得到卫星之间的相对位置、相对速度;其次,求得太阳矢量在轨道坐标系下的面内角和面外角;其次,引入太阳光压摄动差对相对运动的影响,推导得到含太阳光压作用的相对运动方程;最后,建立两星相对位置与光压摄动系数之间的关系,解算相对运动方程中光压摄动系数。本发明提供的高轨非合作目标在轨光压辨识方法解决了编队卫星高精度相对运动建模问题,解决了在非合作目标的面质比不可知情况下的光压摄动系数估计问题,大大提高了相对运动模型精度以及在轨制导律的控制精度,为提高编队飞行的控制精度提供有力支撑。
Description
技术领域
本发明涉及卫星制导、导航与控制技术领域,具体涉及一种高轨非合作目标在轨光压辨识方法。
背景技术
卫星受到太阳光压照射时,太阳辐射能量一部分被吸收,另一部分被反射,这种能量转换使卫星受到力的作用,称为太阳辐射压力,简称光压。卫星受到的光压加速度主要取决于卫星光照面积与卫星质量的比值。对于近距离编队飞行,对相对运动模型产生影响的轨道摄动主要是太阳光压摄动差。
高轨卫星一个轨道周期为24小时,当高轨编队卫星之间的面值比差异较大时,光压力差的长时间累积影响显著,用CW方程建模预测的相对运动轨迹在光压作用下异化严重,传统的CW方程建立的相对运动模型已经不能满足高精度非合作目标相对运动状态的精度要求,需要构造新的方程建立高精度相对运动模型,准确表示在光压作用下非合作目标相对运动状态。
发明内容
本发明提供了一种高轨非合作目标在轨光压辨识方法,通过分析光压摄动影响规律,构造含光压摄动的相对运动方程,基于星上相对测距测角设备信息和相对导航算法,解决非合作目标光压摄动系数的估计问题,从而实现高精度的相对运动建模。
本发明提供的高轨非合作目标在轨光压辨识方法包含以下过程:
步骤1根据星上相对导航设备所提供的测距测角信息得到卫星之间相对位置、相对速度;
步骤2根据太阳矢量在地球静止卫星轨道坐标系中的变化关系,计算描述太阳位置关系的面内角和面外角;
步骤3利用编队卫星上相对测量设备,获得编队卫星相对位置速度信息,根据两星相对的光压摄动加速度模型,构造包含光压作用的相对运动方程;
步骤4根据两星伴飞状态下间隔半个轨道周期相对导航输出的编队卫星的高精度相对位置信息,解算相对运动方程中光压摄动系数。
可选地,所述的步骤2中,建立面内角θ和面外角β与轨道系下太阳矢量的关系计算公式为:
其中,Sox、Soy、Soz分别为太阳矢量在轨道系下X、Y、Z方向的位置坐标。
可选地,所述的步骤3中,构建的含光压作用的相对运动方程的计算公式为:
其中,太阳光压系数k=1,p=4.65×10-6N/m2,/>表示两颗卫星的面质比之差,A为卫星受太阳照射面积,m为卫星质量,θ和β表示太阳矢量在卫星轨道系下的面内角和面外角,n为卫星轨道角速度,/>为卫星的初始相对位置和速度。
可选地,所述的步骤4中,解算步骤3中构建的相对运动方程中光压摄动系数的计算公式为:
M=N cotβ (16)
与现有技术相比,本发明实现了高轨非合作目标在轨光压辨识,避免了非合作目标相对运动建模过程中由于光压差摄动带来的误差,同时解决了地面无法对非合作目标面质比进行估计的难题,数学仿真验证表明该辨识方法可以满足高精度相对运动建模要求。
附图说明
图1为太阳矢量在轨道坐标系中示意图。
具体实施方式
当编队卫星之间的面值比差异较大时,两星受到的太阳光压力存在108~107量级的压力差,高轨一个轨道周期为24小时,光压力差经长时间累计影响显著,由光压力差引起的相对运动轨迹异化超过200m/轨。为实现非合作目标的高精度相对运动建模,本发明提供了一种高轨非合作目标在轨光压辨识方法;首先,得到卫星之间的相对位置、相对速度;其次,求得太阳矢量在轨道坐标系下的面内角和面外角;其次,引入太阳光压摄动差对相对运动的影响,推导得到含太阳光压作用的相对运动方程;最后,建立两星相对位置与光压摄动系数之间的关系,解算相对运动方程中光压摄动系数。
以下结合附图以及具体实施例,对本发明做进一步说明阐述。
1.获得卫星之间的相对位置、相对速度
相对量测方程由星载测量设备的输出,例如微波雷达量测的相对距离,相对视线方位角和相对视线仰角,对于其它种类的测量设备,测量原理相似,只是精度不同,得到测量方程为:
其中,ρ为两星的相对距离,α为目标星在伴飞卫星轨道坐标下的俯仰角,γ为目标星在伴飞卫星轨道坐标下的方位角,x、y、z分别是目标星伴飞卫星轨道坐标下的三轴相对位置,V是观测噪声,t为相对运动时间。
设计导航滤波算法,根据姿态测量系统提供的相对参考轨道系的姿态和相对导航敏感器的输出完成滤波计算。最终输出得到两星的相对位置相对速度信息。
2.求太阳矢量在轨道坐标系下的面内角和面外角
如图1所示,在静止轨道卫星轨道坐标中,通过两个角度,即面内角θ和面外角β描述太阳矢量S在轨道系中任意位置关系,太阳矢量S方向定义为卫星指向太阳。
面内角θ定义为太阳矢量S在轨道面内的投影矢量与轨道坐标系Z轴正向的夹角,取值范围为(0,2π),顺时针增大。面内角θ物理意义为太阳相对于静止轨道卫星,一天内均匀绕卫星旋转一圈。
面外角β定义为太阳矢量S与追踪星轨道坐标系轨道面的夹角,即太阳矢量S与轨道面法向(Y轴负向)夹角的余角,取值范围-is~is,其中is为黄赤交角,常值23.45°。面外角β物理意义为在一年太阳直射点在地球的南北回归线间来回一周,太阳光相对地球赤道面(即静止轨道卫星的轨道面)的夹角(即面外角)变化周期为一年的长周期运动,变化范围为-23.45°~23.45°。
面内角θ和面外角β与轨道系下太阳矢量的关系计算公式为:
其中,Sox、Soy、Soz分别为太阳矢量S在轨道系下X、Y、Z方向的位置坐标。
3.引入太阳光压摄动差对相对运动的影响,推导得到含太阳光压作用的相对运动方程
根据太阳矢量S在轨道系的表达式,得到单星轨道系下太阳光压摄动加速度方程,建立两星相对光压摄动加速度模型,从而得到包含太阳光压的相对运动状态方程,然后求到含光压摄动的相对运动的解析解表达式。
1)单星轨道系下太阳光压摄动加速度方程
利用面内角θ和面外角β,轨道坐标系中太阳矢量S可表示为:
进而,在轨道坐标系中太阳光压加速度表示为:
其中,k=1,p=4.65×10-6N/m2,为卫星面质比,n为卫星轨道角速度。A为卫星受太阳照射的面积,m为卫星质量,t为相对运动时间。对于不同太阳高度角来说,对卫星的光压影响仅仅是幅值E,F不同,如春分点和秋分点(面外角最小),轨道面内加速度为幅值达到最大,Y方向加速度为零;夏至和冬至(面外角最大),Y方向加速度最大,所以本发明的光压模型具有普适性。
2)相对运动光压摄动加速度建模
由光压矢量模型得到地球同步轨道单颗卫星太阳光压摄动加速度模型为:
a(t)=[E sin(nt+θ) F E cos(nt+θ)]T (6)
地球同步轨道上两颗近距离卫星,光压摄动加速度表达式中,仅仅是卫星的面质比不同,因此相对运动受到的光压摄动加速度为:
u(t)=[M sin(nt+θ) N M cos(nt+θ)] (7)
3)构建含光压作用的相对运动方程
包含太阳光压的相对运动状态方程为:
可得:
得到含太阳光压的相对运动解析解(简称为S_CW方程):
4.通过两星相对位置与光压摄动系数之间的关系,求解相对运动方程中光压摄动系数
根据S_CW方程中Y方向相对位置表达式(11)分析可知,由所述函数表达式可知,通过观测不同时刻Y的值,就可以辨识得到N的值,M=N cotβ,由于面外角β已知,可通过N得到M,从而辨识得到在轨的光压摄动在各分量上的增量模型,实现在轨辨识的目的,简化了在轨求解光压模型的复杂计算。
相对导航记录初始时刻Y方向相对位置:
y(0)=y0 (13)
半个轨道周期后再一次记录Y方向相对位置:
式(13)、(14)联立求得N的表达式:
根据M=N cotβ,又已知面外角β,通过N得到M表达式:
M=N cotβ (16)
即解算相对运动方程中光压摄动系数的计算公式为:
M=N cotβ (16)
本发明提供的高轨非合作目标在轨光压辨识方法,先采用面内角θ和面外角β描述太阳矢量S在轨道系中任意位置关系,然后根据面内角θ和面外角β关系,建立光压摄动差表达式,代入相对运动状态方程,求解含光压摄动差的相对运动解析解,根据解析解表达式,通过间隔半个轨道周期的Y方向的相对导航数据直接求解到光压摄动系数。本发明提供的高轨非合作目标在轨光压辨识方法满足高精度相对运动建模要求。
尽管本发明的内容已经通过上述优选实施例作了详细介绍,但应当认识到上述的描述不应被认为是对本发明的限制。在本领域技术人员阅读了上述内容后,对于本发明的多种修改和替代都将是显而易见的。因此,本发明的保护范围应由所附的权利要求来限定。
Claims (1)
1.一种高轨非合作目标在轨光压辨识方法,其特征在于,包含以下过程:
步骤1根据星上相对导航设备所提供的测距测角信息得到卫星之间相对位置、相对速度;
步骤2根据太阳矢量在地球静止卫星轨道坐标系中的变化关系,计算描述太阳位置关系的面内角和面外角;
步骤3利用编队卫星上相对测量设备,获得编队卫星相对位置速度信息,根据两星相对的光压摄动加速度模型,构造包含光压作用的相对运动方程;
步骤4根据两星伴飞状态下间隔半个轨道周期相对导航输出的编队卫星的高精度相对位置信息,解算相对运动方程中光压摄动系数;
所述的步骤2中,建立面内角θ和面外角β与轨道系下太阳矢量的关系计算公式为:
其中,Sox、Soy、Soz分别为太阳矢量在轨道系下X、Y、Z方向的位置坐标;
所述的步骤3中,包含以下过程:
步骤S3.1根据太阳位置关系的面内角和面外角,得到单星轨道系下的太阳光压摄动加速度方程:
利用面内角θ和面外角β,轨道坐标系中太阳矢量S可表示为:
进而,在轨道坐标系中太阳光压加速度表示为:
步骤S3.2建立两星相对光压摄动加速度模型,从而得到包含太阳光压的相对运动状态方程:
由光压矢量模型得到地球同步轨道单颗卫星太阳光压摄动加速度模型为:
a(t)=[Esin(nt+θ)F Ecos(nt+θ)]T (6)
地球同步轨道上两颗近距离卫星,光压摄动加速度表达式中,仅仅是卫星的面质比不同,因此相对运动受到的光压摄动加速度为:
u(t)=[Msin(nt+θ)N Mcos(nt+θ)] (7)
步骤S3.3构建造含光压作用的相对运动方程,求解得到含太阳光压的相对运动解析解:
包含太阳光压作用的相对运动状态方程为:
可得:
求解得到含太阳光压的相对运动解析解:
其中,太阳光压系数k=1,p=4.65×10-6N/m2,/>表示两颗卫星的面质比之差,A为卫星受太阳照射面积,m为卫星质量,θ和β表示太阳矢量在卫星轨道系下的面内角和面外角,n为卫星轨道角速度,/>为卫星的初始相对位置和速度;
所述的步骤4中,解算步骤3中构建的相对运动方程中光压摄动系数的计算公式为:
M=Ncotβ (16)
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201911247650.5A CN111090830B (zh) | 2019-12-09 | 2019-12-09 | 一种高轨非合作目标在轨光压辨识方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201911247650.5A CN111090830B (zh) | 2019-12-09 | 2019-12-09 | 一种高轨非合作目标在轨光压辨识方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN111090830A CN111090830A (zh) | 2020-05-01 |
CN111090830B true CN111090830B (zh) | 2023-06-02 |
Family
ID=70396349
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201911247650.5A Active CN111090830B (zh) | 2019-12-09 | 2019-12-09 | 一种高轨非合作目标在轨光压辨识方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN111090830B (zh) |
Families Citing this family (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN112068130B (zh) * | 2020-09-01 | 2024-04-09 | 上海卫星工程研究所 | 基于整星二维扫描运动的静止轨道微波成像方法及系统 |
Citations (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN102305630A (zh) * | 2011-05-17 | 2012-01-04 | 哈尔滨工业大学 | 基于扩展卡尔曼滤波的sar卫星自主定轨方法 |
CN104021241A (zh) * | 2014-05-13 | 2014-09-03 | 北京空间飞行器总体设计部 | 一种导航卫星复杂模型的地球反照光压摄动建模方法 |
CN105865459A (zh) * | 2016-03-31 | 2016-08-17 | 北京理工大学 | 一种考虑视线角约束的小天体接近段制导方法 |
-
2019
- 2019-12-09 CN CN201911247650.5A patent/CN111090830B/zh active Active
Patent Citations (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN102305630A (zh) * | 2011-05-17 | 2012-01-04 | 哈尔滨工业大学 | 基于扩展卡尔曼滤波的sar卫星自主定轨方法 |
CN104021241A (zh) * | 2014-05-13 | 2014-09-03 | 北京空间飞行器总体设计部 | 一种导航卫星复杂模型的地球反照光压摄动建模方法 |
CN105865459A (zh) * | 2016-03-31 | 2016-08-17 | 北京理工大学 | 一种考虑视线角约束的小天体接近段制导方法 |
Non-Patent Citations (2)
Title |
---|
J_2摄动和太阳光压对绳系卫星系统的影响分析;徐方暖等;《机械科学与技术》(第01期);全文 * |
地球J_2摄动下卫星编队飞行自主相对定轨研究;宋琛等;《计算机仿真》(第03期);全文 * |
Also Published As
Publication number | Publication date |
---|---|
CN111090830A (zh) | 2020-05-01 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN107607947B (zh) | 基于卡尔曼滤波的星载雷达成像参数在线估计方法 | |
CN102175241B (zh) | 一种火星探测器巡航段自主天文导航方法 | |
CN103674032B (zh) | 融合脉冲星辐射矢量和计时观测的卫星自主导航系统及方法 | |
CN103323026B (zh) | 星敏感器和有效载荷的姿态基准偏差估计与修正方法 | |
CN109613583B (zh) | 基于单星与地面站测向及联合测时差的无源目标定位方法 | |
CN102168981B (zh) | 一种深空探测器火星捕获段自主天文导航方法 | |
CN101692001B (zh) | 一种借力飞行轨道上深空探测器的自主天文导航方法 | |
CN103900576B (zh) | 一种深空探测自主导航的信息融合方法 | |
CN107132542B (zh) | 一种基于光学和多普勒雷达的小天体软着陆自主导航方法 | |
CN104374388A (zh) | 一种基于偏振光传感器的航姿测定方法 | |
CN104567880A (zh) | 一种基于多源信息融合的火星最终接近段自主导航方法 | |
CN103389099B (zh) | 基于x射线脉冲星的航天器姿态和位置测量系统及方法 | |
CN104655152A (zh) | 一种基于联邦滤波的机载分布式pos实时传递对准方法 | |
CN107727101B (zh) | 基于双偏振光矢量的三维姿态信息快速解算方法 | |
CN103148856B (zh) | 一种基于自适应尺度变化的借力飞行探测器自主天文导航方法 | |
CN103852082A (zh) | 一种星间测量和陀螺的姿态轨道一体化滤波估计方法 | |
CN107300700B (zh) | 敏捷合成孔径雷达卫星聚束模式姿态机动需求计算方法 | |
CN103645489A (zh) | 一种航天器gnss单天线定姿方法 | |
CN104457705A (zh) | 基于天基自主光学观测的深空目标天体初定轨方法 | |
CN102168989B (zh) | 一种pos方位精度和姿态精度的地面测试方法 | |
CN108917764A (zh) | 一种双星编队仅测距相对导航方法 | |
CN102706363A (zh) | 一种高精度星敏感器的精度测量方法 | |
CN102514734B (zh) | 基于日地系统Halo轨道探测器构型与姿态指向的姿态递推方法 | |
CN104764449A (zh) | 一种基于星历修正的捕获段深空探测器自主天文导航方法 | |
CN104729510A (zh) | 一种空间目标相对伴飞轨道确定方法 |
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 |