CN116125503A - 一种高精度卫星轨道确定及预报算法 - Google Patents

一种高精度卫星轨道确定及预报算法 Download PDF

Info

Publication number
CN116125503A
CN116125503A CN202211231390.4A CN202211231390A CN116125503A CN 116125503 A CN116125503 A CN 116125503A CN 202211231390 A CN202211231390 A CN 202211231390A CN 116125503 A CN116125503 A CN 116125503A
Authority
CN
China
Prior art keywords
orbit
satellite
forecast
track
data
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.)
Pending
Application number
CN202211231390.4A
Other languages
English (en)
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.)
Shaanxi Aerospace Technology Application Research Institute Co ltd
Original Assignee
Shaanxi Aerospace Technology Application Research Institute Co ltd
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 Shaanxi Aerospace Technology Application Research Institute Co ltd filed Critical Shaanxi Aerospace Technology Application Research Institute Co ltd
Priority to CN202211231390.4A priority Critical patent/CN116125503A/zh
Publication of CN116125503A publication Critical patent/CN116125503A/zh
Pending legal-status Critical Current

Links

Images

Classifications

    • 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/01Satellite radio beacon positioning systems transmitting time-stamped messages, e.g. GPS [Global Positioning System], GLONASS [Global Orbiting Navigation Satellite System] or GALILEO
    • G01S19/13Receivers
    • G01S19/24Acquisition or tracking or demodulation of signals transmitted by the system
    • G01S19/27Acquisition or tracking or demodulation of signals transmitted by the system creating, predicting or correcting ephemeris or almanac data within the receiver
    • 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/24Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00 specially adapted for cosmonautical navigation
    • 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/01Satellite radio beacon positioning systems transmitting time-stamped messages, e.g. GPS [Global Positioning System], GLONASS [Global Orbiting Navigation Satellite System] or GALILEO
    • G01S19/13Receivers
    • G01S19/35Constructional details or hardware or software details of the signal processing chain
    • G01S19/37Hardware or software details of the signal processing chain
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00Computer-aided design [CAD]
    • G06F30/20Design optimisation, verification or simulation
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/14Force analysis or force optimisation, e.g. static or dynamic forces
    • YGENERAL 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
    • Y02TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02ATECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
    • Y02A90/00Technologies having an indirect contribution to adaptation to climate change
    • Y02A90/10Information and communication technologies [ICT] supporting adaptation to climate change, e.g. for weather forecasting or climate simulation

Landscapes

  • Engineering & Computer Science (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Signal Processing (AREA)
  • Automation & Control Theory (AREA)
  • Astronomy & Astrophysics (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • General Engineering & Computer Science (AREA)
  • Power Engineering (AREA)
  • Position Fixing By Use Of Radio Waves (AREA)

Abstract

本发明公开一种高精度卫星轨道确定及预报算法,包括以下步骤:轨道确定,处理外测与GNSS内测数据,计算卫星的精密轨道,为注入参数、轨道预报以及轨道控制提供轨道根数;轨道预报,使用精密轨道进行外推计算星历预报,在此基础上进行测控设备跟踪条件预报、设备引导预报、星下点预报、地月影预报等相关预报计算。本发明能够更快捷、更精确的确定大量卫星的瞬时轨道参数,并对未来某一段时间进行较为精确的轨道预报,从而实现对卫星精确有效控制、规避碰撞风险,同时也为卫星的遥测数据接收、卫星平台状态监控、载荷数据传输等工程任务提供可靠、有效、快捷的技术支撑。

Description

一种高精度卫星轨道确定及预报算法
技术领域
本发明属于卫星应用领域,特别涉及一种高精度卫星轨道确定及预报算法。
背景技术
随着国家商业航天产业的不断发展,更多数量、更多功能的卫星星座相继发射并投入使用,在轨卫星的数量得到了井喷式的爆发,卫星的运行轨道参数也不尽相同,尤其在中低轨道,卫星的轨道资源极其优先,卫星布设也较为拥挤,这极大的增加了卫星之间相互碰撞、干扰的风险。因此急需探索一种高精度卫星轨道确定及预报算法,能够更快捷、更精确的确定大量卫星的瞬时轨道参数,并对未来某一段时间进行较为精确的轨道预报,从而实现对卫星的精确有效控制、规避碰撞风险,同时也为卫星的遥测数据接收、卫星平台状态监控、载荷数据传输等工程任务提供可靠、有效、快捷的技术支撑。
发明内容
为解决上述技术问题,本发明公开了一种高精度卫星轨道确定及预报算法,包括以下步骤:
S1,数据预处理,对卫星轨道测量数据进行预处理,产生轨道计算所需格式的数据文件;
S2,初轨计算,利用轨道测量数据计算产生卫星初始轨道信息,供轨道改进使用;
S3,轨道改进,基于轨道动力学方法,利用测量数据对初轨进行迭代改进,生成精密轨道;
S4,星历预报,根据卫星轨道根数,外推卫星在未来一段时间的状态;
S5,测站跟踪预报,基于星历预报结果,结合测控站设备地理坐标,计算测站跟踪预报;
S6,设备引导预报,基于星历预报结果,结合测控站设备地理坐标,计算测控站设备天线跟踪卫星的引导角度;
S7,星下点预报,基于星历预报结果,计算卫星星下点位置;
S8,地月影预报,基于星历预报结果,计算卫星进入地球或者月球阴影区的时间段。
优选的,所述S1中,所述数据预处理的具体步骤包括:
⑴对流层延迟修正;
⑵电离层延迟修正;
⑶数据时标修正;
⑷设备延迟修正;
⑸野值剔除;
⑹数据格式转换。
优选的,所述S2中,对于不同的数据源,根据定轨方法类型分类处理计算初轨。
优选的,根据定轨方法类型分类处理计算初轨的方法:
⑴外测数据定轨,基于拉普拉斯方法进行初轨计算;
⑵GNSS内测数据定轨,利用多项式拟合方法进行初轨计算;
⑶预报星历定轨,提取指定历元轨道数据进行初轨计算。
优选的,所述S3中,所述轨道改进的具体步骤包括:
⑴观测数据建模;
⑵数据选用与权重的设置;
⑶动力学模型的选用设置;
⑷解算参数选取设置。
优选的,所述S3中,采用批处理或序贯处理估值方法可以有效改进定轨计算的稳定性。
优选的,两种估值方法的具体步骤:
⑴批处理
①获取全部观测数据;
②构建观测方程;
③计算协方差及残差;
④法方程求解;
⑤历元最优估值;
⑵序贯处理
①当前历元状态方程;
②状态估值;
③修正参考轨道;
④下一历元状态方程;
⑤最优估值。
优选的,所述S4中,星历预报计算中设置动力学模型参数,需要考虑米级以上的摄动影响。
优选的,所述S5中,所述测站跟踪预报的具体步骤包括:
⑴进站时间预报;
⑵过顶时间计算;
⑶出站时间预报。
优选的,所述S8中,所述地月影预报的具体步骤包括:
⑴太阳位置计算;
⑵阴影预报计算。
本发明的有益效果在于:实现一种高精度卫星轨道确定及预报方法,给定力学模型,实现支持外测和GNSS内测数据单独或者联合使用,进行轨道改进,突破低轨航天器的轨道定轨与预报精度的技术能力提升,形成团队核心竞争力和自主知识产权,并为大量卫星在轨管理提供业务支撑和配套成果转化。
附图说明
图1为高精度卫星轨道确定及预报算法的流程示意图;
图2为本发明轨道确定功能接口示意图;
图3为本发明轨道确定功能模块图;
图4为本发明数据预处理流程图;
图5为本发明初轨计算处理流程图;
图6为本发明轨道改进处理流程图;
图7为本发明轨道确定的两种估值算法流程图;
图8为本发明轨道预报功能接口示意图;
图9为本发明轨道预报功能模块图;
图10为本发明星历预报流程图;
图11为本发明测站跟踪预报流程图;
图12为本发明设备引导预报流程图;
图13为本发明星下点示意图;
图14为本发明星下点预报流程图;
图15为本发明地月影预报流程图;
图16为本发明航天器进出地影的平面图。
具体实施方式
为使本发明的目的、技术方案及优点更加清楚,下面将结合本发明的附图,对本发明的技术方案进行清楚完整的描述。
本发明公开一种高精度卫星轨道确定及预报算法,如图1所示,包括以下步骤:
S1、S2、S3为本发明轨道确定算法的步骤内容,处理外测与GNSS内测数据,计算卫星的精密轨道,为注入参数、轨道预报以及轨道控制提供轨道根数。
精密定轨是通过含有观测误差的各类观测数据,使用精确的动力学拟合卫星轨道的过程。显然,定轨计算的精度除观测几何外,很大程度受限于测量数据的精度,观测模型的精度以及动力学模型的精度。精密定轨涉及到各类观测数据的处理,包括精确的观测建模与各类测量误差的修正,不同轨道类型的卫星所受到的动力学影响也不一样,因此在定轨计算中所考虑动力学模型也不同。另外,不管是动力学模型还是观测模型,都存在一定的误差,其中尤其以动力学模型部分更为显著,比如低轨卫星所受到的大气阻尼力通过大气模型进行建模不可避免的会存在一定的误差。综上所述,在轨道计算过程中,需要选择合适的数学模型,并且对测量数据进行优化处理,从而提升计算结果精度。精确动力学建模技术,卫星定轨计算的精度取决于动力学模型和观测数据,对于当前精度的轨道计算,需要考虑的动力学模型包括,地球质点引力,地球非球形引力,太阳光学,大气阻尼力,后牛顿效应力等。
图2为轨道确定功能接口示意图,通过接口传递的各类信息内容如表1所示:
表1轨道确定功能接口信息定义
序号 接口方向 信息类型 收发方向 传递方式
1. 用户 轨道确定命令 接收 操作界面
2. 轨道数据库 轨道测量数据 接收 数据库访问
3. 轨道数据库 定轨计算结果 发送 数据库访问
轨道确定功能模块如图3所示。
S1,数据预处理,对卫星轨道测量数据进行预处理,产生轨道计算所需格式的数据文件;
数据预处理流程如图4所示,***启动时,初始化数据库接口;按照用户发出的轨道确定命令,从数据库中获取对应的轨道测量数据;执行数据修正、野值剔除、数据格式转换等一系列数据处理过程;完成数据预处理后,输出数据预处理结果。
数据预处理的具体步骤包括:
⑴对流层延迟修正;
基于测站处测量的温度、湿度及压强数据,通过大气折射修正模型对卫星测距与测速完成对流层折射修正。
⑵电离层延迟修正;
基于电离层格网数据对卫星测距与测速完成电离层延迟修正。
⑶数据时标修正;
对各类测量数据的时标进行修正。
⑷设备延迟修正;
对地基测距数据的应答机延迟进行修正。
⑸野值剔除;
对外测数据、GNSS内测数据按照特定的标准进行野值剔除。
⑹数据格式转换;
在轨道测量数据修正处理后进行格式转换,生成定轨计算所需格式的数据内容。
S2,初轨计算,利用轨道测量数据计算产生卫星初始轨道信息,供轨道改进使用;
利用外测数据、内测数据或是预报星历文件产生卫星初始轨道信息,供轨道改进使用。
初轨计算处理流程如图5所示,***启动时,按照指定的定轨方法类型进行分类处理;分别进行外测数据定轨、GNSS内测数据定轨、预报星历定轨等不同类型的初轨计算过程;生成初轨计算结果并输出。
对于不同的数据源,根据定轨方法类型分类处理计算初轨,方法包括:
⑴外测数据定轨,基于拉普拉斯方法进行初轨计算;
基于拉普拉斯方法,利用测距、方位角与高度角完成卫星初轨计算。初轨计算中采用的动力学模型包括地球质点引力,J2摄动,以及日月摄动,适用于地球低轨与中高轨卫星。当前测控***的地基测距精度大约为2m,测速精度为0.01m/s,方位角、高度角的精度为0.01°,对于中低轨航天器可以实现优于500m的初轨精度,满足后续轨道改进的初值需求。
⑵GNSS内测数据定轨,利用多项式拟合方法进行初轨计算;
利用星载GNSS自定位数据,通过多项式拟合方法计算卫星的初轨信息。
对于地轨卫星(轨道高度1000km以下),当前星载GNSS自定位数据的精度优于12m,速度精度优于0.05m/s,初轨精度取决于星历预报的精度。
⑶预报星历定轨,提取指定历元轨道数据进行初轨计算;
从预报星历文件中提取指定历元的卫星位置、速度信息,或是插值计算指定历元的位置速度信息,生成卫星初轨信息。
S3,轨道改进,基于轨道动力学方法,利用测量数据对初轨进行迭代改进,生成精密轨道;
轨道改进基于轨道动力学模型,利用测量数据对初轨进行迭代改进,以生成精密轨道。不同轨道类型的差异主要体现在摄动力模型选择与不同。针对特定类型的轨道,需要在配置文件中设置相应的动力模型使用和数据使用设置。
轨道改进依赖初轨计算提供的先验位置、速度信息,通常采用精确的动力学模型,考虑的力学因素包括地球质点、地球非球形、日月及大行星引力、大气阻尼力和太阳光压等,给定力模型,本发明支持外测和GNSS内测数据单独或者联合使用,进行轨道改进。对于GNSS内测数据,由于自定位精度优于15m,轨道改进通过动力学拟合可以实现优于12m的轨道精度;对于地基测量,受限于测量数据的精度和跟踪弧段,对于低轨卫星,采用一天4圈弧段的外测数据位置精度优于200米,速度精度优于0.05m/s;对于中高轨卫星,采用三站24小时分时跟踪,位置精度优于200米,速度精度优于0.05m/s。
轨道改进处理流程如图6所示,***启动时,执行一系列进程初始化过程;基于初轨计算结果及测量数据,执行迭代改进;迭代改进收敛完成后,生成精密定轨结果并输出。
轨道改进的具体步骤包括:
⑴观测数据建模;
对外测、内测以及星间测量等观测数据进行精确观测建模。
⑵数据选用与权重的设置;
通过配置文件,设置参数定轨计算的某种数据类型,或是几类数据类型的组合,观测数据的弧段范围,各类观测数据使用的权重。
⑶动力学模型的选用设置;
根据卫星轨道类型,设置使用指定的动力学模型,可选的动力学模型包括:非球形引力的阶数,固体潮、海潮摄动,大气阻尼力,太阳光压,日月引力及其他大天体摄动。对于不同轨道类型的卫星,可以根据配置文件设置使用不同的动力学模型,低轨与高轨卫星的在动力学模型方面的差异主要体现在大气阻尼力、地球非球形的阶数,以及第三体摄动。在不考虑效率的前提下,可以将高阶非球形、大气阻尼力以及第三体摄动都设置为使用,即可同时满足低轨与中高轨的轨道计算需求。
⑷解算参数选取设置;
根据需要设置解算参数,可选解算参数包括:①位置、速度;②大气阻尼力;③太阳光压系数;④RTN方向常值经验力;⑤RTN方向周期性经验力。
定轨的基本原理是,利用含有测量误差的测量数据,结合卫星受到的动力学约束,获取卫星状态的最佳估值。对于精密定轨,传统的叫法为轨道改进,但由于可以在卫星定轨的同时确定一些与轨道相关的动力学或是测量上的物理参数,对传统意义下单纯的轨道改进进行了扩充,现在称为精密定轨,采用合适的估值方法可以有效改进定轨计算的稳定性。
采用批处理或序贯处理两种估值方法有效改进定轨计算的稳定性,两种估值方法的具体步骤:
⑴批处理
批处理利用全部观测数据求某一时刻的轨道参数,适用于事后处理。特点是数据较多、求解的精度高,但是存储量大,需要将所有的观测数据存储起来,待数据采样结束后,才可以进行解算,其步骤如下:
①获取全部观测数据;
②构建观测方程;
③计算协方差及残差;
④法方程求解;
⑤历元最优估值;
⑵序贯处理
序贯处理是一有观测数据就解算出此时卫星轨道,没有存储量大的问题,因而在快速、实时定轨中应用广泛。主要问题是递推容易发散,随着时间推移,计算的轨道逐渐偏离真实轨道,需要引入受力模型以补偿算法误差。其步骤如下:
①当前历元状态方程;
②状态估值;
③修正参考轨道;
④下一历元状态方程;
⑤最优估值;
两种估值方法的具体处理流程如图7所示。
通常,卫星运动的状态微分方程可写成如下形式:
Figure BDA0003880330510000101
其中X表示卫星的状态量
Figure BDA0003880330510000102
和其他待估参数
Figure BDA0003880330510000103
Figure BDA0003880330510000104
将上式中X在给定参考状态
Figure BDA0003880330510000105
处进行展开,记
Figure BDA0003880330510000106
略去高阶项后可以表示为一个线性方程,
Figure BDA0003880330510000107
其中,
Figure BDA0003880330510000108
其解可以表示为
x(t)=Φ(t0,t)x(t0)
Figure BDA0003880330510000109
上式中Φ(t0,t0)称为状态转移矩阵,I表示单位矩阵。
观测量与卫星状态量存在一定的联系,可以表示为
Y=G(t,X)+ε
其中Y表示观测量,ε表示测量噪声。
同样,对上式在参考状态处展开,略去高阶项有
Figure BDA0003880330510000111
式中,
Figure BDA0003880330510000112
表示观测量对观测历元处状态量的观测偏导数,H表示对改进历元状态量的观测偏导数。对于时间序列的观测量,公式可以表示为
yi=Hix0i
状态方程的估计通常有两种方式实现:批处理和序贯处理。
所谓批处理就是将需要处理的所有观测数据一次解算,一般取观测的权矩阵W满足:
Figure BDA0003880330510000113
略去推导过程,根据线性无偏最小方差估计可得批处理的估值
Figure BDA0003880330510000114
Figure BDA0003880330510000115
其中Ri为观测数据噪声水平的平方,
Figure BDA0003880330510000116
为待估状态量的先验协方差。
由此计算待估历元的最优估值
Figure BDA0003880330510000117
为:
Figure BDA0003880330510000118
相应协方差矩阵为:
Figure BDA0003880330510000119
序贯处理是一种递推算法。假设已知tk,(k=0,1,2,…,n)时刻的估值
Figure BDA00038803305100001110
和协方差矩阵
Figure BDA00038803305100001111
则tk+1时刻的状态量为:
Figure BDA00038803305100001112
其中,
Figure BDA00038803305100001113
公式构成了序贯处理的递推公式。tk+1时刻修正后的状态估计可以表示为,
Figure BDA0003880330510000121
S4到S8为本发明轨道预报算法的步骤内容,轨道预报,使用精密轨道进行外推计算星历预报,再此基础上进行测控设备跟踪条件预报(三点预报:入站、过顶和出站)、设备引导预报、星下点预报、地月影预报等相关预报计算。
低轨航天器预报精度受限于大气阻尼力的计算,现有的大气模型存在10-15%的模型误差,磁暴期间的大气密度误差可达数倍,对于低轨卫星(轨道高度1000km以下),12小时预报精度优于100m,中高轨卫星7天轨道预报精度优于200米。
如前所述,轨道预报包括星历预报、测站跟踪预报、设备引导预报,星下点预报、地月影预报等。本发明轨道预报部分的代码全部自主研发,采用C++语言进行编码,设计合理的数据结构,采用面向对象的编程思想,重新梳理和优化了算法流程,尤其是对变量的存储和访问进行优化,矩阵和向量运算增加涵盖了线性代数的绝大部分基础功能,处理微分方程和数值计算提供了良好的运算基础,数据存储采用双精度浮点运算,效率更高精度更优。轨道预报24小时,计算时间优于30秒;单星多站8弧段数据定轨计算优于1分种。
图8为轨道预报功能接口示意图,通过接口传递的各类信息内容如表2所示:
表2轨道预报功能接口信息定义
序号 接口方向 信息类型 发方向 传递方式
1. 用户 轨道预报命令 接收 操作界面
2. 轨道数据库 卫星轨道根数 接收 数据库访问
3. 轨道数据库 预报计算结果 发送 数据库访问
轨道预报功能模块如图9所示。
S4,星历预报,根据卫星轨道根数,外推卫星在未来一段时间的状态;
轨道根数支持开普勒六根数和两行轨道根数,其中开普勒六根数的计算模型为HPOP高精度轨道预报模型。两行轨道根数的预报模型采用NOARD的SGP4/SDP4模型,二者都是目前业内普遍采用的计算方法,其精度和效率都经过大量的实际工程测试验证,能够很好的支持测运控业务的底层算法业务。
星历预报处理流程如图10所示,***启动时,初始化数据库连接;按照用户发出的轨道预报命令,从数据库中获取对应的卫星轨道根数;设置动力学模型参数;计算卫星星历预报;完成预报计算后,保存星历预报结果;
对于卫星的运动规律,在地心天球坐标系下考虑,采用位置矢量r和速度矢量
Figure BDA0003880330510000131
来描述。
运动方程即为初值问题P
Figure BDA0003880330510000132
右函数a可分为a0和aε两部分,a0为二体问题下航天器的加速度,aε为其它各种摄动力加速度之和,可表示如下
Figure BDA0003880330510000133
考虑在米级以上的摄动影响,卫星受影响包括:地球非球型引力;日月引力;太阳光压摄动;大气阻力摄动;固体潮摄动等。
⑴地球非球形引力摄动
在地固坐标系中,地球非球形引力摄动的位函数为:
Figure BDA0003880330510000134
式中:
Figure BDA0003880330510000135
Figure BDA0003880330510000136
分别为归一化的勒让德和缔合缔合勒让德函数
Figure BDA0003880330510000141
为归一化的地球引力场带谐项系数
Figure BDA0003880330510000142
Figure BDA0003880330510000143
为归一化的地球引力场田谐项系数
GME为地球引力常数
RE为平均地球赤道半径
为经纬度
r为地心距
有关归一化过程如下:
Figure BDA0003880330510000144
Figure BDA0003880330510000145
Figure BDA0003880330510000146
Figure BDA0003880330510000147
Figure BDA0003880330510000148
Figure BDA0003880330510000149
的递推关系如下:
Figure BDA00038803305100001410
Figure BDA00038803305100001411
Figure BDA00038803305100001412
Figure BDA00038803305100001413
的递推关系如下:
Figure BDA00038803305100001414
Figure BDA0003880330510000151
Figure BDA0003880330510000152
Figure BDA0003880330510000153
由此可得地球非球形摄动加速度如下:
Figure BDA0003880330510000154
式中:
Figure BDA0003880330510000155
Figure BDA0003880330510000156
Figure BDA0003880330510000157
Figure BDA0003880330510000158
Figure BDA0003880330510000159
Figure BDA00038803305100001510
其中:
Figure BDA00038803305100001511
Figure BDA00038803305100001512
勒让德多项式的偏导数使用下面的推导公式计算:
Figure BDA0003880330510000161
Figure BDA0003880330510000162
Figure BDA0003880330510000163
扇谐项和田谐项中缔合勒让德多项式偏导数递推公式为:
Figure BDA0003880330510000164
Figure BDA0003880330510000165
Figure BDA0003880330510000166
Figure BDA0003880330510000167
以上是航天器在地球固定坐标系中的加速度,还需转换到相应的历元地心天球坐标系。
⑵N体引力摄动
太阳、月球和太阳系的大行星对飞行器都有引力作用,其效应可以近似描述为点质量摄动P。在历元地心天球坐标系里,可以使用下面的公式计算N体摄动:
Figure BDA0003880330510000168
式中:GMBjB为第j个摄动体的引力常数;Rj为第j个摄动体在地心惯性系中的位矢;Δj为航天器相对第j个摄动体的位矢。美国喷气推进实验室(Jet PropulsionLaboratory,JPL)生成的历表提供了日月和大行星位置的精确值。
⑶太阳辐射压摄动
太阳辐射压力摄动是由于太阳辐射到航天器上对航天器产生压力造成的P。也简称光压摄动。因此,只有当航天器位置处于地球和月球的阴影之外时才有光压摄动。
本体部分的太阳辐射摄动加速度计算如下:
Figure BDA0003880330510000171
其中:
F:阴影因子,地球和月球的阴影分为本影、半影和伪半影,阴影因子可以根据日、月、地球和航天器的位置以及柱形或锥形地影模型来计算。
ρSR:在航天器附近太阳光压强常数。
η:航天器受照表面的反射系数。η的取值范围为:0<η<1,一般可取η=0.5。η可作为太阳辐射压力摄动的校正因子,在轨道估值过程中进行估值。
m:航天器质量。
ΔS:航天器至太阳的矢量。
A:垂直于
Figure BDA0003880330510000172
的航天器横截面积。
不同的其太阳帆板的安装角度和对太阳定向的控制方式可能是不同的。太阳辐射对太阳帆板压力的方向也就不同。因此,要根据各航天器太阳帆板的具体情况计算太阳辐射摄动加速度。太阳帆板摄动加速度的计算:
Figure BDA0003880330510000173
其中:
AP:太阳帆板的面积。
β:太阳帆板对太阳辐射的反射系数,可作为被估值量。
⑷大气阻力摄动
围绕地球飞行的低轨航天器,其所处的环境并非真空,大气对航天器会产生阻力P。大气阻力摄动可以描述为:
Figure BDA0003880330510000181
式中:CBDB为大气阻力系数;A为飞行器的截面积在垂直于轨道的平面上的投影;m为飞行器的质量;ρBaB为飞行器处的大气密度;Vr为飞行器相对于大气的速度。
对于有太阳帆板的航天器,同时需要考虑太阳帆板产生的大气阻力,其公式与上面的公式类似:
Figure BDA0003880330510000182
式中:ABpB为太阳帆板的面积;CBDPB为适用于太阳能帆板的大气阻力系数;γ为太阳能帆板法向和飞行器速度的夹角;|APcosγ|为太阳能帆板在垂直于轨道的平面上的有效面积。
上述的大气密度ρa可以使用大气密度模型计算,这些模型包括:指数模型,Harris-Priester,Jacchia 77,DTM(Drag Temperature Model)和MSIS系列等,可见有关资料。
⑸固体潮摄动
由于地球不是一个刚体,在月球和太阳引力位的作用下,地球的陆地部分会发生弹性形变P。这种形变称固体潮。固体潮使地球外壳的起伏振幅达到20~30厘米。固体潮使得地球内部质量分布随时间而变化,从而使得地球的引力场也随时间而变化,这一影响称动力潮汐。由固体潮引起的地球引力位的这种变化称固体潮附加位。通常,引潮位是用谐波函数表示的。不同阶的引潮位对应的LOVE数和负荷形变系数的数值不同。目前对于固体潮只考虑二阶引潮位。由固体潮引起的地球引力势的变化可以表达为地球引力势球谐函数系数的变化:
Figure BDA0003880330510000191
式中:
Figure BDA0003880330510000192
为日月引力引起的地球引力位系数的变化;HBkB,Θk为k分潮波的振幅和相位;
Figure BDA0003880330510000193
为k阶Love数;δ0m为Kronecker算子;ae为平均地球赤道半径。
固体潮附加位引起的摄动加速度的计算与地球非球形引力摄动加速度的计算类似,这里不再描述。
⑹相对论影响
由于广义相对论效应,飞行器在地球质心为原点的局部惯性坐标系的运动方程将不同于仅考虑牛顿引力场时的运动方程P。这种差异可看作飞行器受到了一附加摄动。近地飞行器的广义相对论摄动可以描述为
Figure BDA0003880330510000194
式中:c为光速;rE,
Figure BDA0003880330510000195
为飞行器在地心惯性系中的位置和速度矢量;GMe为地球的引力常数;β,γ为相对论参数,对Einstein广义相对论其值为1。
S5,测站跟踪预报,基于星历预报结果,结合测控站设备地理坐标,计算测站跟踪预报;
基于星历预报结果,结合测控站设备地理坐标,计算测站跟踪预报即三点预报,预报内容包括进站、出站和过顶等三个时间点的测站时间、方位角、高度角、距离和距离变化率。
测站跟踪预报处理流程如图11所示,***启动时,初始化数据库连接;按照用户发出的轨道预报命令,从数据库中获取对应的卫星星历预报;依次计算各角度过顶时间;计算过顶时间对应的距离等其它参数;完成预报计算后,保存测站跟踪预报结果。
测站跟踪预报的具体步骤包括:
⑴进站时间预报;
按照跟踪高度截止角计算卫星进出站时间和方位角,并计算距离和距离变率。
⑵过顶时间计算;
计算卫星过顶时间,并计算高度角和方位角,并计算距离和距离变率。
⑶出站时间预报;
按照跟踪高度截止角计算卫星进出站时间和方位角,并计算距离和距离变率。
S6,设备引导预报,基于星历预报结果,结合测控站设备地理坐标,计算测控站设备天线跟踪卫星的引导角度;
基于星历预报结果,结合测控站设备地理坐标,计算测控站设备天线跟踪卫星的引导角度。通常需要在天线所在的当地大地水平面建立坐标系(ENU,东北天坐标系),用来描述卫星在天线坐标系中的运动状态,状态的描述可由三个基本物理量表征方位角(Az,取值范围0~360°)、俯仰角(El,取值范围0~90°)、距离(R)。
设备引导预报的处理流程如图12所示,***启动时,初始化数据库连接;按照用户发出的轨道预报命令,从数据库中获取对应的卫星星历预报;在设备跟踪时间段内,间隔固定时间计算对应的设备天线方位角及俯仰角;完成预报计算后,保存设备引导预报结果。
S7,星下点预报,基于星历预报结果,计算卫星星下点位置;
基于星历预报结果,计算卫星星下点位置。星下点位置可以反映卫星相对于地球表面的运行轨迹,为卫星开展应用业务提供决策依据,如对某个区域的覆盖情况和重访周期等。卫星星下点可由三个基本物理量表征经度(Lon,取值范围-180~180°,东经为正),纬度(Lat,取值范围-90~90,北纬为正),高度(Alt)。将某段时间内卫星的星下点连成一片,会在地面形成一条轨迹,即星下点轨迹。假设地球为旋转椭球体,过卫星质心做地球的法线,定义该法线与椭球面的交点为星下点,即卫星在椭球面上的垂直投影定义为星下点,如图13所示。
星下点预报的处理流程如图14所示,***启动时,初始化数据库连接;按照用户发出的轨道预报命令,从数据库中获取对应的卫星星历预报;在指定的预报时间段内,间隔固定时间计算对应的星下点位置;完成预报计算后,保存星下点预报结果。
S8,地月影预报,基于星历预报结果,计算卫星进入地球或者月球阴影区的时间段;
基于星历预报结果,计算卫星进入地球或者月球阴影区的时间段。有些卫星载荷正常工作需要特定的光照条件,因此必须要考虑卫星当前相对于太阳、月球的相对位置关系,给出阴照面和阳照面的预报结果,为任务计划的编排提供决策依据。
地月影预报的处理流程如图15所示,***启动时,初始化数据库连接;按照用户发出的轨道预报命令,从数据库中获取对应的卫星星历预报;在指定的预报时间段内,计算太阳所在位置;在指定的预报时间段内,分别计算地球、月球所在位置;在指定的预报时间段内,计算卫星进出地球或月球阴影区的时刻;完成预报计算后,保存地月影预报结果。
地月影预报具体步骤包括:
⑴太阳位置计算;
读取行星历表计算太阳位置。
⑵阴影预报计;
根据太阳位置,计算卫星进出地影的时间。
地月影计算原理:由于太阳很远,可认为太阳在无穷远处(即忽略视差),地球上的太阳光线为平行光.从图16上可看出,测站位置矢量R与太阳方向L^(单位矢量)的夹角大于90°+|β0|,即在第二象限,此条件如下:
Figure BDA0003880330510000211
在地面由90°+|β0|确定的边界称为日界线,日界线把地球表面分为两个区域,即“白天”和“黑夜”,当测站为黑夜时才能进行光学观测。
光学观测还需要航天器表面被太阳照亮,此时设备才能观测到反射的光线,这一条件是由航天器与太阳的相对位置关系决定的,也就是航天器必须在地影外,对于这一条件,只需考虑简单的柱形地影即可,见图16。
由图16可知,地影边界上满足下述条件:
Figure BDA0003880330510000221
式中cosψ的右端取负号,是由于ψ>90°,当航天器和太阳方向的夹角小于ψ时,航天器不在地影内,可以被太阳照亮,太阳与地心的连线与椭球面的交点被称为日下点,地面上以日下点为中心半径为ψ的区域就是航天器可以被照亮的区域,该区域边界就称为地影线,同样可用轨道半长径代替r来估计无地影区域的大小,即
Figure BDA0003880330510000222
对于航天器的光学观测而言,前面提到的三个条件都必须满足,对于一个航天器,可观测范围相对测站是不变的,而日界线和地影线随太阳位置的变化而由东向西移动,对于其他观测,如激光观测、无线电观测等,上述后两个条件理论上不是必需的。

Claims (10)

1.一种高精度卫星轨道确定及预报算法,其特征在于,包括以下步骤:
S1,数据预处理,对卫星轨道测量数据进行预处理,产生轨道计算所需格式的数据文件;
S2,初轨计算,利用轨道测量数据计算产生卫星初始轨道信息,供轨道改进使用;
S3,轨道改进,基于轨道动力学方法,利用测量数据对初轨进行迭代改进,生成精密轨道;
S4,星历预报,根据卫星轨道根数,外推卫星在未来一段时间的状态;
S5,测站跟踪预报,基于星历预报结果,结合测控站设备地理坐标,计算测站跟踪预报;
S6,设备引导预报,基于星历预报结果,结合测控站设备地理坐标,计算测控站设备天线跟踪卫星的引导角度;
S7,星下点预报,基于星历预报结果,计算卫星星下点位置;
S8,地月影预报,基于星历预报结果,计算卫星进入地球或者月球阴影区的时间段。
2.根据权利要求1所述的高精度卫星轨道确定及预报算法,其特征在于,所述S1中,所述数据预处理的具体步骤包括:
⑴对流层延迟修正;
⑵电离层延迟修正;
⑶数据时标修正;
⑷设备延迟修正;
⑸野值剔除;
⑹数据格式转换。
3.根据权利要求1所述的高精度卫星轨道确定及预报算法,其特征在于,所述S2中,对于不同的数据源,根据定轨方法类型分类处理计算初轨。
4.根据权利要求3所述的高精度卫星轨道确定及预报算法,其特征在于,根据定轨方法类型分类处理计算初轨的方法:
⑴外测数据定轨,基于拉普拉斯方法进行初轨计算;
⑵GNSS内测数据定轨,利用多项式拟合方法进行初轨计算;
⑶预报星历定轨,提取指定历元轨道数据进行初轨计算。
5.根据权利要求1所述的高精度卫星轨道确定及预报算法,其特征在于,所述S3中,所述轨道改进的具体步骤包括:
⑴观测数据建模;
⑵数据选用与权重的设置;
⑶动力学模型的选用设置;
⑷解算参数选取设置。
6.根据权利要求1所述的高精度卫星轨道确定及预报算法,其特征在于,所述S3中,采用批处理或序贯处理估值方法可以有效改进定轨计算的稳定性。
7.根据权利要求6所述的高精度卫星轨道确定及预报算法,其特征在于,两种估值方法的具体步骤:
⑴批处理
①获取全部观测数据;
②构建观测方程;
③计算协方差及残差;
④法方程求解;
⑤历元最优估值;
⑵序贯处理
①当前历元状态方程;
②状态估值;
③修正参考轨道;
④下一历元状态方程;
⑤最优估值。
8.根据权利要求1所述的高精度卫星轨道确定及预报算法,其特征在于,所述S4中,星历预报计算中设置动力学模型参数,需要考虑米级以上的摄动影响。
9.根据权利要求1所述的高精度卫星轨道确定及预报算法,其特征在于,所述S5中,所述测站跟踪预报的具体步骤包括:
⑴进站时间预报;
⑵过顶时间计算;
⑶出站时间预报。
10.根据权利要求1所述的高精度卫星轨道确定及预报算法,其特征在于,所述S8中,所述地月影预报的具体步骤包括:
⑴太阳位置计算;
⑵阴影预报计算。
CN202211231390.4A 2022-10-09 2022-10-09 一种高精度卫星轨道确定及预报算法 Pending CN116125503A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202211231390.4A CN116125503A (zh) 2022-10-09 2022-10-09 一种高精度卫星轨道确定及预报算法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202211231390.4A CN116125503A (zh) 2022-10-09 2022-10-09 一种高精度卫星轨道确定及预报算法

Publications (1)

Publication Number Publication Date
CN116125503A true CN116125503A (zh) 2023-05-16

Family

ID=86305149

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202211231390.4A Pending CN116125503A (zh) 2022-10-09 2022-10-09 一种高精度卫星轨道确定及预报算法

Country Status (1)

Country Link
CN (1) CN116125503A (zh)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116299584A (zh) * 2023-05-18 2023-06-23 中科星图测控技术股份有限公司 一种卫星的实时测量跟踪并发处理方法和***
CN118276123A (zh) * 2024-05-29 2024-07-02 天津云遥宇航科技有限公司 基于sgp4算法的掩星卫星高度与相位监控***及方法

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116299584A (zh) * 2023-05-18 2023-06-23 中科星图测控技术股份有限公司 一种卫星的实时测量跟踪并发处理方法和***
CN118276123A (zh) * 2024-05-29 2024-07-02 天津云遥宇航科技有限公司 基于sgp4算法的掩星卫星高度与相位监控***及方法

Similar Documents

Publication Publication Date Title
CN116125503A (zh) 一种高精度卫星轨道确定及预报算法
Kai et al. Autonomous navigation for a group of satellites with star sensors and inter-satellite links
Hwang et al. GPS‐Based Orbit Determination for KOMPSAT‐5 Satellite
CN111427002A (zh) 地面测控天线指向卫星的方位角计算方法
Iwata Precision attitude and position determination for the Advanced Land Observing Satellite (ALOS)
CN112713922A (zh) 一种多波束通讯卫星的可见性快速预报算法
Hesar et al. Small body gravity field estimation using LIAISON supplemented optical navigation
Shtark et al. Regional positioning using a low Earth orbit satellite constellation
CN112649006A (zh) 一种太阳同步圆轨道的轨道规划方法
Zeitlhöfler Nominal and observation-based attitude realization for precise orbit determination of the Jason satellites
CN111125874B (zh) 一种动平台高精度测轨预报方法
Iiyama et al. Autonomous Distributed Angles-Only Navigation and Timekeeping in Lunar Orbit
Rutkowska et al. SLR technique used for description of the Earth elasticity
Vallado et al. Accurate orbit determination from short-arc dense observational data
CN115392540A (zh) 一种用于月球轨道交会制导的快速预报方法
Watkins Tracking station coordinates and their temporal evolution as determined from laser ranging to the LAGEOS satellite
Vailado et al. Accurate orbit determination from short-arc dense observational data
Mikrin et al. Circumlunar spacecraft navigation using the measurements from global navigation satellite systems glonass, gps, galileo and beidou
Arbinger et al. Impact of orbit prediction accuracy on low earth remote sensing flight dynamics operations
Zhou Onboard orbit determination using GPS measurements for low Earth orbit satellites
Keil Kalman filter implementation to determine orbit and attitude of a satellite in a molniya orbit
Martin-Mur et al. Deep-space navigation using optical communications systems
Duong et al. Orbit determination by range-only data
Straka et al. Navigation of satellite measurements without ground control points
Jacobs et al. The extragalactic and solar system celestial frames: Accuracy, stability, and interconnection

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