CN104236553A - 一种全天时的星光折射卫星自主定位方法 - Google Patents

Abstract

本发明一种全天时的星光折射卫星自主定位方法，首先采用具有饱和抑制功能的CCD进行观测并获取星图，并对获得的星图进行图像处理，提取未折射星和折射星；用三角匹配算法实现该星图与导航星表的星图匹配，获得星图上每颗星的赤经、赤纬；从星图中为每颗折射星寻找离自身最近的恒星，作为该颗折射星的对应星，计算获得每颗折射星的折射角γ；根据大气模型和折射角，计算获得每颗折射星的切向高度h；选取观测星图上三颗折射星，根据他们的赤经、赤纬以及切向高度h，计算获得地心惯性坐标系下卫星的位置；本发明实现了卫星的全天时高精度天文自主导航，且数据更新率高、涉及设备简单、成本低、功耗低。

Description

一种全天时的星光折射卫星自主定位方法

技术领域

本发明属于天文导航领域，涉及一种全天时的星光折射卫星自主定位方法，尤其适用于具有全天时和高精度自主导航需求的卫星。

背景技术

星光折射天文导航是利用光学敏感器测量恒星星光在通过地球边缘大气层时所发生的折射，间接得到地平信息，获得卫星在地心坐标系中的位置，从而确定卫星轨道的方法。星光折射导航的概念自1982年提出至今，受到许多学者的关注，并取得了一定的研究成果，但Robert，Robert和Eliezer(Guidanceand Control Conference,Gatlinburg,TN,August 15-17,1983:359-367)，Whilte，Thurman和Barnes(Proceedings of the Forty-First Annual Meeting,the Institute of Navigation,Annapolis,MD,June 25-27,1985:83-89)，Lair和Duchon(Acta Astronautica,1988,17(10):1069-1079)，Ning和Fang(Aerospace Science and Technology 11(2007)：222-228)等这些主要研究者更多讨论的是不考虑空间环境影响下的星光折射导航位置解算原理和方法，由于地球的光照情况将严重影响星光折射的导航精度，甚至直接导致星光折射导航系统不能获得定位结果，所以之前的方法未能满足卫星自主定位对全天时和高精度的需求。

中国专利公开号CN 103616028A，公开日是2014年3月5日，名称为“一种基于单星敏感器的星光折射卫星自主导航方法”中公开了一种基于单星敏感器利用星光折射原理进行卫星自主导航的方法，该方法介绍了利用星光折射进行卫星自主导航的基本思路和初步仿真结果，并未解决该方法涉及到的地球强背景下恒星提取等关键技术，所以目前该方法只适用于地球光照条件极好的一小段时间内，不能满足卫星对全天时和高精度自主导航的需求。中国专利公开号CN 103630109A，公开日2014年3月12日，名称为“一种基于星光折射确定地心矢量的方法”中公开了一种基于多个星敏感器利用星光折射原理进行地心矢量确定的方法，该方法所需仪器设备多，技术难度大且精度低，不适用于具有全天时和高精度需求的卫星自主导航。

发明内容

本发明的技术解决问题是：克服现有技术的不足之处，提供一种全天时的星光折射卫星自主定位方法，实现了卫星的全天时高精度天文自主导航，且数据更新率高、涉及设备简单、成本低、功耗低。

本发明的技术解决方案是：一种全天时的星光折射卫星自主定位方法，其特征在于步骤如下：

1)采用具有饱和抑制功能的CCD对地球、及其周边大气和大气层外的恒星进行观测并获取星图；观测时，用于观测未经大气折射的未折射星的视场大于50平方度；

2)对步骤1)获得的星图进行图像处理；

21)对星图进行定标；

22)遍历整个星图，设定亮目标信噪比阈值，扣除亮目标，包括地球和亮恒星，再根据星图上剩余的能量分布拟合出背景图，并将该背景图从原始星图中扣除；

23)将扣除背景图后的星图进行滤波，设定恒星提取信噪比阈值，提取信噪比大于该阈值的目标，完成对恒星的提取；所述的恒星包括经大气折射的折射星和未经大气折射的未折射星；

24)通过恒星形心定位法对步骤23)提取的恒星进行定位，得到每颗恒星在星图上的位置信息；

3)根据步骤24)获得的未折射星的位置信息，利用三角匹配算法实现经过步骤2)处理后的星图与导航星表的星图匹配，获得星图上每颗星的赤经、赤纬；

4)在导航星表中，为每颗折射星寻找离其自身最近的恒星，作为该颗折射星的对应星，根据折射星的赤经、赤纬以及该颗折射星对应星的赤经、赤纬，计算由折射星及其对应星的赤经、赤纬形成的两个向量的夹角，获得每颗折射星的折射角γ；

5)根据大气模型γ＝2350.1074e^-0.10326788h和步骤4)获得的折射角γ，计算获得每颗折射星的切向高度h；

6)选取观测星图上三颗折射星，根据他们的赤经、赤纬以及切向高度h，计算获得地心惯性坐标系下卫星的位置；

61)根据三颗折射星的赤经、赤纬解算中心轴O_bE在天球上的位置(α_c,δ_c)；所述的中心轴O_bE为以三个折射星方向为母线的圆锥的中心轴；并计算出该圆锥的半锥角η；

62)旋转卫星惯性坐标系F_b，使卫星惯性坐标系F_b的Z_b轴与O_bE方向一致，旋转后的坐标系记为F_d＝(X_d,Y_d,Z_d)^T，其中上标T表示转置，则其转换关系为：

其中，表示绕Z_b轴逆时针旋转-α_c，表示绕Y_b轴逆时针旋转δ_c-90°；

根据转换关系转换三颗折射星的坐标，得出在F_d坐标系下的位置为(α_i，δ_i)(i＝1，2，3)；

63)根据以下三个方程求解三个变量ξ、α

ζ²+ξ²-2ζξcos(α-α_i)＝a_i ²，(i＝1，2，3)；

其中的h_i为第i颗折射星的切向高度；

64)根据以下方程获得F_d坐标系中的地心矢量e_d

e_d＝ζ(cosαX_d+sinαY_d)+ξctgηZ_d；

65)根据以下坐标转换关系获得地心在卫星惯性坐标系F_b中的矢量e_b为

66)获得卫星在地心惯性坐标系F_i下的坐标矢量为-e_b。

本发明与现有技术相比的有益效果是：

本发明采用了地球强背景下恒星提取技术，保证了星光折射导航可进行全天时的工作，否则当太阳位于星光折射敏感器的后方时，将大部分时间无法开展工作，其它太阳光照情况下，星光折射导航也只能获得有限范围下的折射星，极大地影响了导航精度，因此，在同等条件下，相比于其它方法，基于这项关键技术的星光折射敏感器可获得最大范围下的折射星，即提高了导航精度，也保证了较高的数据更新率。

附图说明

图1是本发明星光折射导航定位流程图；

图2是本发明恒星相对位置示意图；

图3是本发明卫星位置解算几何关系图；

图4是本发明卫星位置解算投影示意图。

具体实施方式

下面结合附图和实例对本发明作进一步的详细说明。

本发明提出一种全天时的星光折射卫星自主定位方法，工作流程见图1，首先基于饱和抑制功能的CCD获取观测图像，进行图像处理，包括地球强背景下恒星提取，接下来基于提取的未折射星和高精度导航星表完成星图识别，证认出折射星对应的恒星，算出折射角，再基于大气模型解算出选取的折射星的切向高度，最后根据折射星的折射角和位置矢量解算每个折射星组合的定位结果。

具体实现步骤如下：

步骤一：星图成像

由于视场中受光照的地球亮度最亮时超过-22星等，而实现星光折射导航视场中至少需有3颗折射星，由于当探测极限星等达到7.5星等时，才可保证运行在GEO轨道时有效视场中(大气层20km～80km)有3颗折射星的概率为96.5％以上，因此，只有使用具有饱和抑制功能的CCD才可实现同一视场中同时对高亮度地球和7.5等恒星的同时探测。

观测时采用具有饱和抑制功能的CCD获取星图，实现对高亮度地球和恒星的同时成像。饱和抑制功能可行性分析如下：

在空间观测受太阳光照的地球，计算亮地球背景产生的光子数为：

X＝ψ_bAq_btα²

上式中ψ_b为单位像元上的背景光子流，单位：m^-2s^-1(″)^-2；A为望远镜接收面积，单位：m^-2；q_b为探测器对天空背景的平均量子效率；t为曝光时间，单位：s；α²为单位像元角面积(″)²。

目前CCD的满阱电荷一般都能达到120000，并且具有抑制过载电荷量可超过10000倍的饱和抑制功能，即具有抑制超过1.2×10⁹电荷量的功能，大于亮地球背景引起的1.0×10⁸的电荷量，因此，CCD通过饱和抑制功能，在保证正常观测的情况下，可满足地球强背景下超大动态范围成像的需求。

由于实施步骤中的星图识别需使用未折射恒星来完成与高精度导航星表的匹配，因此，星图成像时需保证留有50平方度以上的视场来观测未折射恒星。

步骤二：图像处理

地球强背景下的恒星提取为星光折射导航的关键技术之一。由于适合用于星光折射导航的大气层范围为20km～80km，可见折射星的视场范围很小，GEO轨道时最大有效视场角仅为0.11°，并且由于星光折射导航至少需探测到3颗折射星，且折射星分布越均匀定位精度越高，因此，需使用大视场星敏感器，在高轨时可观测整个地球，然而，观测视场中的地球绝大部分时候都有光照的影响，受太阳照射的地球大气亮度最高可达9mag/arcsec²，所以需实现在高亮度的大气背景噪声下对7.5星等折射星的探测，即通过解决地球强背景下的恒星提取技术来实现最大范围的折射星获取，保证星光折射导航精度。

地球强背景下恒星提取的步骤如下：

1)对星图进行定标

去除星图中由于观测系统的不均匀性带来的误差；具体操作包括修正CCD本底，校正仪器像素-像素不均匀效应，图像大尺度改正，标记坏像素等。

2)扣除背景

CCD势阱电荷数字化后进入计算机的灰度值，并非仅仅由来自观测目标的光线产生，也包括天光背景等环境光源的贡献。为了探测暗弱目标，需要精确估算不同天区的天光背景噪声；大尺度的背景噪声可能会淹没局部的暗弱信号。计算天光背景时，首先遍历整帧图像，设定亮目标信噪比阈值，扣除亮目标，包括地球和亮恒星，再根据星图上剩余的能量分布拟合出背景图，并将该背景图从原始星图中扣除；

3)目标探测与提取

根据图形特征参数对扣除背景图后的星图进行识别，完成滤波处理，并设定恒星提取信噪比阈值，提取信噪比大于该阈值的目标，完成对恒星的提取；所述的恒星包括经大气折射的折射星和未经大气折射的未折射星；

4)恒星形心定位

由于象元匹配和衍射效应等因素，观测目标往往在CCD靶面上占据多个像素。首先需要对候选像素中的能量进行提取，得出目标的能量分布图。将临近的像素点聚类形成各个目标的能量分布图。再通过恒星形心定位法完成对恒星的高精度定位，生成目标在CCD靶面XY坐标系中的坐标。

步骤三：星图识别

由于规划的姿态保证了观测时留有50平方度以上的视场来观测未折射恒星，且星光折射敏感器的探测极限星等优于7.5星等，因此，获得的未折射恒星数目几乎都在50颗以上，且这些星的分布范围较大，利用三角匹配算法来实现与高精度导航星表的星图匹配，可保证95％以上的匹配成功率。完成星图匹配后，可获得星图上每颗星的赤经、赤纬。

步骤四：折射星折射角解算

星图识别成功后，可证认出每颗折射星对应的恒星，即获得每颗折射星折射前的赤经和赤纬，基于星图上提取的折射星赤经赤纬，计算两个向量的夹角，获得每颗折射星的折射角γ。折射角的测量精度直接取决于恒星位置测量精度。

步骤五：折射星切向高度解算

目前构建了范围适用、与星光大气折射规律相吻合的基于大气折射率变化的星光大气折射模型来进行仿真分析。该大气模型的拟合公式如下所示：

γ＝2350.1074e^-0.10326788h

其中h为大气高度，单位km；γ为大气折射角，单位为"。因此，把折射角输入到以上的大气模型中来即可获得该颗折射星的切向高度。

折射星切向高度的误差主要由折射角测量精度、大气模型及其误差决定。

步骤六：位置解算

位置解算算法的输入为三颗折射星的赤经、赤纬和切向高度，输出为地心惯性坐标系下卫星的位置。具体算法如下：

假设单位天球的球心在卫星质心O_b，地球为球形体，其半径为R_e；星敏感器观测到星光穿过地球边缘附近的三颗恒星S₁，S₂，S₃在天球上的视在位置用S_i(α_i,δ_i)(i＝1，2，3)表示；天球坐标(α,δ)与卫星惯性坐标系F_b的关系为

X_b＝cosδcosα，Y_b＝cosδsinα，Z_b＝sinδ

天球上的三点为S_i(i＝1,2,3)。下面由几何关系得到地心的方位，其中恒星相对位置示意图见图2。

假设三条星光不是全与地球表面相切，可以确定出O_bS_i(i＝1，2，3)为母线的圆锥A的中心轴O_bE的方向，假设O_bE在天球上的位置为(α_c,δ_c)，圆锥A的半锥角η为O_bS_i与O_bE的夹角。

旋转坐标系F_b，使Z_b轴与O_bE的方向一致，旋转后的坐标系记为F_d＝(X_d,Y_d,Z_d)^T，则其转换关系为：

以上公式中的表示绕Z_b轴逆时针旋转-α_c，表示绕Y_b轴逆时针旋转δ_c-90°。

图3为卫星位置解算几何关系图，设在F_d坐标系下的S_i的位置分别为(α_i，δ_i)(i＝1，2，3)。过地心O_e和O_e在O_bE上的垂足O_c分别作星光O_bS_i(i＝1，2，3)的垂线，垂足分别为E_i和C_i，其中O_eE_i即为折射星的切向高度h_i。设星光O_bS_i与过O_e且垂直于O_bE的平面交点为B_i(i＝1，2，3)，把O_e、O_c、B_i及B_i所在的圆在O_bX_bY_b平面上投影，O_e投影点记为O_c投影点即为O_b，B_i投影点记为D_i。则D_i在F_d坐标系下的天球上的位置为(α_i，0)(i＝1，2，3)，见图4的卫星位置解算投影示意图。

由几何关系，显然

$O_{e^{\′}}{D}_1=O_e{B}_1=\frac{O_e{E}_1}{\cos \mathrm{\η}},O_{e^{\′}}{D}_2=O_e{B}_2=\frac{O_e{E}_2}{\cos \mathrm{\η}}$

$O_{e^{\′}}{D}_3=O_e{B}_3=\frac{O_e{E}_3}{\cos \mathrm{\η}},O_b{D}_i=O_c{B}_i=\frac{O_c{C}_i}{\cos \mathrm{\η}}\stackrel{\mathrm{\Δ}}{=}\mathrm{\ξ},(i=\mathrm{1,2,3})$

记

$O_b{O}_{e^{\′}}=O_c{O}_e=\mathrm{\ζ},O_{e^{\′}}{D}_i=a_i,(i=\mathrm{1,2,3})$

在F_d坐标系中的位置矢量O_bO_e记为e_d，则

当a₁＝a₂＝a₃时，此时与O_b重合,O_e与O_c重合，则

e_d＝a₁ctgηZ_d

当a₁＝a₂＝a₃不成立时，则与O_b不重合,O_e与O_c不重合，记在F_d坐标系下的天球上的位置为(α，0)，则

ζ²+ξ²-2ζξcos(α-α_i)＝a_i ²

解上述三个方程(三个未知量)可以求得ξ、α。显然

$O_b{O}_c=\frac{O_c{C}_1}{\sin \mathrm{\η}}=\mathrm{\ξctg\η}$

在F_d坐标系中，O_cO_e方向的单位矢量记为e_c，则

e_c＝cosαX_d+sinαY_d

e_d＝ξctgηZ_d+ζe_c

从而得地心在F_b坐标系下的矢量e_b为

由给出的F_i与F_b的转换关系，以及上面得到的地心在F_b坐标系下的矢量e_b，可得卫星在地心惯性坐标系F_i中的坐标矢量r＝-C_bi ^Te_b。

实施例1：

本发明基于以上介绍的强背景下恒星提取图像处理技术，通过地基观测月亮及其周围恒星，验证地球强背景下恒星提取的可行性，试验结果表明：使用兴隆观测基地测光辅助望远镜，10秒曝光，在月亮边缘附近可识别12等星。再根据天基观测和地基观测情况的不同和相同之处，针对不同的观测背景亮度和观测设备，基于以下探测系统信噪比公式即可推算出星光折射敏感器可实现对地球强背景下7.5星等恒星的探测。

$s/N={\mathrm{\ψ}}_s\sqrt{\left({\mathrm{Aq}}_s^{2}t\right)/\left({\mathrm{\ψ}}_b{q}_b{\mathrm{\α}}^2\right)}$

上式中ψ_s为单位像元上的信号光子流，单位：m^-2s^-1；ψ_b为单位像元上的背景光子流，单位：m^-2s^-1(″)^-2；A为望远镜接收面积，单位：m^-2；q_s为探测器对信号的平均量子效率；q_b为探测器对天空背景的平均量子效率；t为曝光时间，单位：s，其中设背景光和信号光积分时间相同；α²为单位像元角面积(″)²。

实施例2：

卫星星光折射导航仿真系统在恒星测量精度3"(3σ)、大气模型误差为1％的情况下，通过仿真卫星GTO(20000km～36000km)过程中的星光折射导航功能，得出星光折射导航系统的定位精度优于1.2km。

本发明说明书中未作详细描述的内容属本领域技术人员的公知技术。

Claims (1)

1.一种全天时的星光折射卫星自主定位方法，其特征在于步骤如下：

1)采用具有饱和抑制功能的CCD对地球、及其周边大气和大气层外的恒星进行观测并获取星图；观测时，用于观测未经大气折射的未折射星的视场大于50平方度；

2)对步骤1)获得的星图进行图像处理；

21)对星图进行定标；

22)遍历整个星图，设定亮目标信噪比阈值，扣除亮目标，包括地球和亮恒星，再根据星图上剩余的能量分布拟合出背景图，并将该背景图从原始星图中扣除；

23)将扣除背景图后的星图进行滤波，设定恒星提取信噪比阈值，提取信噪比大于该阈值的目标，完成对恒星的提取；所述的恒星包括经大气折射的折射星和未经大气折射的未折射星；

24)通过恒星形心定位法对步骤23)提取的恒星进行定位，得到每颗恒星在星图上的位置信息；

3)根据步骤24)获得的未折射星的位置信息，利用三角匹配算法实现经过步骤2)处理后的星图与导航星表的星图匹配，获得星图上每颗星的赤经、赤纬；

4)在导航星表中，为每颗折射星寻找离其自身最近的恒星，作为该颗折射星的对应星，根据折射星的赤经、赤纬以及该颗折射星对应星的赤经、赤纬，计算由折射星及其对应星的赤经、赤纬形成的两个向量的夹角，获得每颗折射星的折射角γ；

5)根据大气模型γ＝2350.1074e^-0.10326788h和步骤4)获得的折射角γ，计算获得每颗折射星的切向高度h；

6)选取观测星图上三颗折射星，根据他们的赤经、赤纬以及切向高度h，计算获得地心惯性坐标系下卫星的位置；

61)根据三颗折射星的赤经、赤纬解算中心轴O_bE在天球上的位置(α_c,δ_c)；所述的中心轴O_bE为以三个折射星方向为母线的圆锥的中心轴；并计算出该圆锥的半锥角η；

62)旋转卫星惯性坐标系F_b，使卫星惯性坐标系F_b的Z_b轴与O_bE方向一致，旋转后的坐标系记为F_d＝(X_d,Y_d,Z_d)^T，其中上标T表示转置，则其转换关系为：

其中，表示绕Z_b轴逆时针旋转-α_c，表示绕Y_b轴逆时针旋转δ_c-90°；

根据转换关系转换三颗折射星的坐标，得出在F_d坐标系下的位置为(α_i，δ_i)(i＝1，2，3)；

63)根据以下三个方程求解三个变量ξ、α

ζ²+ξ²-2ζξcos(α-α_i)＝a_i ²，(i＝1，2，3)；

其中的h_i为第i颗折射星的切向高度；

64)根据以下方程获得F_d坐标系中的地心矢量e_d

e_d＝ζ(cosαX_d+sinαY_d)+ξctgηZ_d；

65)根据以下坐标转换关系获得地心在卫星惯性坐标系F_b中的矢量e_b为

66)获得卫星在地心惯性坐标系F_i下的坐标矢量为-e_b。