一种基于模型参数估计的星光矢量蒙气差修正方法与流程

1.本发明属于天文导航、大气光学领域，具体涉及一种基于模型参数估计的星光矢量蒙气差修正方法。


背景技术：

2.星敏感器是一种以恒星作为参考，经星图捕获、星点质心提取、星图识别/跟踪、姿态解算等处理步骤，获得空间姿态的嵌入式传感设备。作为目前绝对姿态测量精度最高的传感器，在大气层外及临近空间的航行体上已得到广泛的应用。然而对于航行于大气层内，尤其在近地面空间的载体平台来说，使用星敏感器作为姿态测量器件，首先需要克服星光大气折射(蒙气差)导致的测量误差问题。
3.对于蒙气差修正技术，常用方法是获得观测星矢量/参考星矢量的天顶角以及观测环境的大气参数，根据蒙气差经验模型计算每个星光矢量的蒙气差经验值，再对星光矢量进行修正。所以这种方法就依赖于大气参数传感器(如温度、湿度、气压传感器等)，并受传感器误差及经验公式误差的影响(在修正观测星矢量时，还需依赖外部姿态传感器，并受该传感器准确度的影响。如具有漂移特点的陀螺仪等)。


技术实现要素：

4.本发明的目的在于克服以往星光矢量蒙气差修正方法的不足，提出一种不依赖其他外部传感器并自适应获得蒙气差经验模型参数的方法。该方法利用星光矢量间角距不变的特点，建立蒙气差参数估计量测方程，再采用扩展卡尔曼滤波(extended kalman filter,ekf)算法，建立蒙气差参数滤波方程，实现不依赖外部传感器而获得随大气变化的蒙气差经验模型参数，可以简单有效地进行参考星矢量蒙气差修正，进而提高大气内星敏感器姿态解算准确度。
5.本发明采用的技术方案为：一种基于模型参数估计的星光矢量蒙气差修正方法，包括以下步骤：
6.步骤一：星敏感器拍摄观测星图，获得成像面上星点坐标，根据星敏感器内参、畸变模型及小孔成像模型获得星敏感器坐标系下的若干观测星矢量，再进行星图识别，得到观测星与导航星表上参考星的匹配关系；
7.步骤二：根据星敏感器拍摄时间和地点以及地球自转等信息，基于国际天文联合会(international astronomical union,iau)基础天文学标准(standards of fundamental astronomy,sofa)提供的高精度坐标系转换链，可将步骤一获得的参考星从导航星表采用的天球坐标系转换到星敏感器观测位置的地理坐标系，即获得地理坐标系下无蒙气差的参考星矢量，进而可以获得参考矢量的真天顶距(无蒙气差的参考星矢量和天顶矢量的夹角)；
8.步骤三：选择满足工程需要的蒙气差经验模型，该模型将星光矢量的蒙气差大小表示为关于天顶距的三阶函数，各阶系数都是关于大气环境特征的函数，在本方法中设为
待估计的参数。在步骤二基础上，计算所有参考星矢量的蒙气差大小并对其沿仰角增加的方向做蒙气差修正，所得修正的参考星矢量都以待估计参数表示；
9.步骤四：根据星光矢量夹角不变原理，星敏感器坐标系下任意两观测矢量的夹角等于观测点地理坐标系下两对应的修正后参考矢量的夹角。所以可构造蒙气差参数估计量测方程，即以所有观测矢量对的夹角与对应参考矢量对的夹角之差表示观测变量列，该观测变量列可以构造为其关于待估计参数的雅可比矩阵和参数代差的乘积形式。
10.步骤五：在步骤四的基础上采用扩展卡尔曼滤波(extended kalman filter,ekf)算法解算蒙气差模型参数：1、初始化蒙气差模型参数及误差协方差矩阵；2、计算步骤四所述观测变量列及观测变量列雅可比矩阵；3、计算卡尔曼滤波增益；4、计算参数代差；5、更新蒙气差模型参数及误差协方差矩阵；
11.步骤六：重复上述步骤一到步骤五，直至获得稳定的参数估计结果即可认为蒙气差模型参数估计值可用于修正参考星矢量。代入步骤三得到修正的参考星矢量。
12.其中，所述步骤一中，星敏感器拍摄观测星图，所得第i颗观测星在星敏感器坐标系下的观测星矢量记为v
i
。所述星图识别方法可采用子图同构类算法(如三角形算法、金字塔算法等)，也可采用模式识别类算法(如栅格算法、采用径向和环向特征的识别算法等)，还可以采用其他类型的识别算法(如基于神经网络的星图识别算法等)。所述星图识别算法不唯一，需结合实际星敏条件选择合适的算法，但需保证识别结果正确。
13.其中，所述步骤二中，参考星坐标系转换时，星敏感器拍摄时间为曝光中心时刻，表示方法可选utc时间，精确到毫秒；地点为曝光中心时刻星敏感器的经度、纬度、高度，表示方法可选wgs84椭球参考系，精确到秒和米。地球自转信息包括地球当日极移和dut1时间，可由国际地球自转服务(international earth rotation service,iers)提供。基于国际天文联合会(international astronomical union,iau)基础天文学标准(standards of fundamental astronomy，sofa)提供的高精度坐标系转换链进行参考星坐标系转换，还需依赖导航星表上参考星的赤经、赤纬、自行、视差、径向速度等参数。观测点地理坐标系即为东北天坐标系，所得地理坐标系下第i颗参考星无蒙气差的矢量记为w
i
，天顶矢量u＝[0，0，1]
t
，真天顶距为：
[0014]
z
i
＝arccos(u
·
w
i
)
ꢀꢀ
(1)
[0015]
其中，所述步骤三中，蒙气差经验模型如下：
[0016]
ρ＝a1tanz a2tan3z a3tan5z
ꢀꢀ
(2)
[0017]
其中，a1，a2，a3即为待估计参数，ρ为星光矢量的蒙气差大小，z为星光矢量的真天顶距；计算所有参考星矢量的蒙气差大小ρ
i
并对其沿仰角增加的方向做蒙气差修正，所得修正的参考星矢量都是关于待估计参数a1，a2，a3的函数，修正后的参考星矢量表示为
[0018][0019]
其中，i表示单位矩阵，n
i
表示参考星矢量w
i
增加仰角ρ
i
时的转轴矢量：
[0020]
n
i
＝w
i
×
u
ꢀꢀ
(4)
[0021]
n
i
∧表示由矢量元素构成的反对称矩阵，n
i
(1)、n
i
(2)、n
i
(3)分别表示矢量n
i
的第1、2、3个元素：
[0022][0023]
其中，所述步骤四中，根据星光矢量夹角不变原理，星敏感器坐标系下任意两观测矢量的夹角等于观测点地理坐标系下两对应的修正后参考矢量的夹角，可得：
[0024][0025]
若有n对匹配的星光矢量，则可构造个观测变量r
ij
：
[0026][0027]
待估计参数可写为：
[0028][0029]
可构造蒙气差参数估计量测方程：
[0030]
r
k
＝h
k
d
k
ꢀꢀ
(9)
[0031]
其中，表示所有观测变量r
ij
构成的观测变量列，
[0032]
表示r
k
关于待估计参数a1，a2，a3的雅可比矩阵，表示待估计参数代差。
[0033]
其中，所述步骤五中，采用扩展卡尔曼滤波(ekf)算法解算蒙气差模型参数的滤波过程具体如下：
[0034]
(1)初始化参数a1，a2，a3，d0，p0；其中，d0为初始的d
k
(k＝0)；p0为初始误差协方差矩阵：
[0035]
(2)读取一幅星图中星点坐标按上述方法及公式完成对量测方程的观测变量r
k
以及系数矩阵h
k
的计算。
[0036]
(3)计算卡尔曼滤波增益k
k
为：
[0037][0038]
式中，e为星敏感器的测量噪声方差强度阵。
[0039]
(4)由观测变量r
k
更新d
k
为：
[0040]
d
k
＝d
k
‑1 k
k
[r
k
‑
h
k
d
k
‑1]
ꢀꢀ
(11)
[0041]
式中，d
k
‑1代表第k
‑
1张星图迭代出的蒙气差模型参数a1，a2，a3与第k
‑
2张星图迭代出的参数a1，a2，a3之间的差值。
[0042]
(5)更新参数a1，a2，a3的值：
[0043]
(a1，a2，a3)
k 1
＝(a1，a2，a3)
k
d
k
ꢀꢀ
(12)
[0044]
(6)计算下一张星图的协方差矩阵p
k 1
并重置d
k 1
，返回(1)：
[0045]
p
k 1
＝(i
‑
k
k
h
k
)p
k
ꢀꢀ
(13)
[0046]
其中，所述步骤六中，重复上述滤波过程直至获得稳定的参数估计结果即可认为完成了蒙气差模型参数估计过程，此后每帧星图所迭代得到的蒙气差模型参数a1，a2，a3可带入公式(3)完成当前参考星矢量蒙气差修正，进而使得星敏感器姿态解算更加准确。
[0047]
本发明与现有技术相比的优点在于：
[0048]
(1)本发明无需依赖其他外部传感器，如温度、气压、湿度、姿态传感器等，降低了星光矢量蒙气差补偿实施难度及系统成本；
[0049]
(2)本发明可根据环境变化自主估计蒙气差模型参数，避免了外部传感器及经验模型公式引入误差的问题；
附图说明
[0050]
图1为本发明一种基于模型参数估计的星光矢量蒙气差修正方法实现流程图；
[0051]
图2为参考星矢量蒙气差修正示意图。
具体实施方式
[0052]
下结合附图以及具体实施方式进一步说明本发明。
[0053]
一种基于模型参数估计的星光矢量蒙气差修正方法，其流程如图1所示，包括以下步骤：
[0054]
步骤一：星敏感器拍摄观测星图，获得成像面上星点坐标，根据星敏感器内参、畸变模型及小孔成像模型获得星敏感器坐标系下的若干观测星矢量，再进行星图识别，得到观测星与导航星表上参考星的匹配关系；
[0055]
其中，所述步骤一中，星敏感器拍摄观测星图，所得第i颗观测星成像坐标为x
i
，y
i
，星敏感器内参包括主点x0，y0、焦距f、畸变参数p1，p2，p3，q1，q2，q3，畸变校正后的星点成像坐标为：
[0056][0057]
其中：
[0058]
r2＝(x
i
‑
x0)2 (y
i
‑
y0)2ꢀꢀ
(2)
[0059]
所得第i颗观测星在星敏感器坐标系下的观测星矢量为：
[0060][0061]
所述星图识别方法可采用子图同构类算法(如三角形算法、金字塔算法等)，也可采用模式识别类算法(如栅格算法、采用径向和环向特征的识别算法等)，还可以采用其他类型的识别算法(如基于神经网络的星图识别算法等)。所述星图识别算法不唯一，需结合实际星敏条件选择合适的算法，但需保证识别结果正确。
[0062]
步骤二：根据星敏感器拍摄时间和地点以及地球自转等信息，基于国际天文联合会(international astronomical union，iau)基础天文学标准(standards of fundamental astronomy,sofa)提供的高精度坐标系转换链，可将步骤一获得的参考星从导航星表采用的天球坐标系转换到星敏感器观测位置的地理坐标系，即获得地理坐标系下无蒙气差的参考星矢量，进而可以获得参考矢量的真天顶距(无蒙气差的参考星矢量和天顶矢量的夹角)；
[0063]
其中，所述步骤二中，参考星坐标系转换时，星敏感器拍摄时间为曝光中心时刻，表示方法可选utc时间，精确到毫秒；地点为曝光中心时刻星敏感器的经度、纬度、高度，表示方法可选wgs84椭球参考系，精确到秒和米。地球自转信息包括地球当日极移和dut1时间，可由国际地球自转服务(international earth rotation service，iers)提供。基于国际天文联合会(iau)基础天文学标准(sofa)提供的高精度坐标系转换链进行参考星坐标系转换，还需依赖导航星表上参考星的赤经、赤纬、自行、视差、径向速度等参数。观测点地理坐标系即为东北天坐标系，所得地理坐标系下第i颗参考星无蒙气差的矢量记为w
i
，天顶矢量u＝[0，0，1]
t
，真天顶距为：
[0064]
z
i
＝arccos(u
·
w
i
)
ꢀꢀ
(4)
[0065]
步骤三：选择满足工程需要的蒙气差经验模型如式(5)，该模型将星光矢量的蒙气差大小ρ表示为关于天顶距z的三阶函数，系数a1，a2，a3都是关于大气环境特征的函数，在本方法中设为待估计的参数。
[0066]
ρ＝a1tanz a2tan3z a3tan5z
ꢀꢀ
(5)
[0067]
在步骤二基础上，将所有参考星真天顶距z
i
带入式(5)计算所有参考星矢量的蒙气差大小ρ
i
并如图2所示，对其沿仰角增加的方向做蒙气差修正，所得修正的参考星矢量都是关于待估计参数a1，a2，a3的函数，修正后的参考星矢量表示为
[0068][0069]
其中，i表示单位矩阵，n
i
表示参考星矢量w
i
增加仰角ρ
i
时的转轴矢量：
[0070]
n
i
＝w
i
×
u
ꢀꢀ
(7)
[0071]
n
i
∧表示由矢量元素构成的反对称矩阵，n
i
(1)、n
i
(2)、n
i
(3)分别表示矢量n
i
的第1、2、3个元素：
[0072][0073]
步骤四：根据星光矢量夹角不变原理，星敏感器坐标系下任意两观测矢量的夹角等于观测点地理坐标系下两对应的修正后参考矢量的夹角。所以可构造蒙气差参数估计量
测方程，即以所有观测矢量对的夹角与对应参考矢量对的夹角之差表示观测变量列，该观测变量列可以构造为其关于待估计参数的雅可比矩阵和参数代差的乘积形式。
[0074]
其中，所述步骤四中，根据星光矢量夹角不变原理，星敏感器坐标系下任意两观测矢量的夹角等于观测点地理坐标系下两对应的修正后参考矢量的夹角，可得：
[0075][0076]
若有n对匹配的星光矢量，则可构造个观测变量r
ij
：
[0077][0078]
待估计参数可写为：
[0079][0080]
可构造蒙气差参数估计量测方程：
[0081]
r
k
＝h
k
d
k
ꢀꢀ
(12)
[0082]
其中，表示所有观测变量r
ij
构成的观测变量列，
[0083]
表示r
k
关于待估计参数a1，a2，a3的雅可比矩阵，
[0084]
表示待估计参数代差。
[0085]
步骤五：在步骤四的基础上采用扩展卡尔曼滤波(extended kalman filter，ekf)算法解算蒙气差模型参数：1、初始化蒙气差模型参数及误差协方差矩阵；2、计算步骤四所述观测变量列及观测变量列雅可比矩阵；3、计算卡尔曼滤波增益；4、计算参数代差；5、更新蒙气差模型参数及误差协方差矩阵；
[0086]
其中，所述步骤五中，采用扩展卡尔曼滤波(ekf)算法解算蒙气差模型参数的滤波过程具体如下：
[0087]
(1)初始化参数a1，a2，a3，d0，p0；其中，d0为初始的d
k
(k＝0)；p0为初始误差协方差矩阵：
[0088]
(2)读取一幅星图中星点坐标按上述方法及公式完成对量测方程的观测变量r
k
以及系数矩阵h
k
的计算。
[0089]
(3)计算卡尔曼滤波增益k
k
为：
[0090][0091]
式中，e为星敏感器的测量噪声方差强度阵。
[0092]
(4)由观测变量r
k
更新d
k
为：
[0093]
d
k
＝d
k
‑1 k
k
[r
k
‑
h
k
d
k
‑1]
ꢀꢀ
(14)
[0094]
式中，d
k
‑1代表第k
‑
1张星图迭代出的蒙气差模型参数a1，a2，a3与第k
‑
2张星图迭代出的参数a1，a2，a3之间的差值。
[0095]
(5)更新参数a1，a2，a3的值：
[0096]
(a1，a2，a3)
k 1
＝(a1，a2，a3)
k
d
k
ꢀꢀ
(15)
[0097]
(6)计算下一张星图的协方差矩阵p
k 1
并重置d
k 1
，返回(1)：
[0098]
p
k 1
＝(i
‑
k
k
h
k
)p
k
ꢀꢀ
(16)
[0099]
步骤六：重复上述步骤一到步骤五，直至获得稳定的参数估计结果即可认为蒙气差模型参数估计值可用于修正参考星矢量。代入步骤三得到修正的参考星矢量。
[0100]
其中，所述步骤六中，重复上述滤波过程直至获得稳定的参数估计结果即可认为完成了蒙气差模型参数估计过程，此后每帧星图所迭代得到的蒙气差模型参数a1，a2，a3可带入公式(3)完成当前参考星矢量蒙气差修正，进而使得星敏感器姿态解算更加准确。
[0101]
实施效果
[0102]
为验证本发明提出的星光矢量蒙气差修正方法的有效性，在国家天文台兴隆观测站进行了外场观星实验，实验具体内容如下：
[0103]
在观星前，将星敏感器安装在转台上，星敏感器初始视轴指向约为天顶方向。在观星时，设置转台分别转动0
°
，10
°
，20
°
，30
°
，40
°
，50
°
，得到对应的六组星图数据各2000帧。对星图数据采用本发明提出的方法进行参考星矢量蒙气差修正，采用星内角统计偏差δr作为评价指标：
[0104][0105]
得到修正前后对比结果如表1：
[0106]
表1修正前后对比结果表
[0107]

技术特征：
1.一种基于模型参数估计的星光矢量蒙气差修正方法，其特征在于：所述方法包括以下步骤：步骤一，由星敏感器拍摄星图获得星敏感器坐标系下的观测星矢量，再经过星图识别获得观测星与参考星匹配关系；步骤二，将参考星从天球坐标系转换到星敏感器观测时的地理坐标系，并计算参考星矢量的真天顶距；步骤三，根据蒙气差经验模型和参考星矢量真天顶距计算星光矢量蒙气差值，再对参考星矢量沿仰角增加的方向做蒙气差修正，保留蒙气差模型参数为待估计值，得到修正的参考星矢量；步骤四，根据星光矢量夹角不变原理，构造蒙气差参数估计量测方程，并计算方程中观测变量列及其关于参数的雅可比矩阵；步骤五，采用扩展卡尔曼滤波算法解算并更新蒙气差模型参数和误差协方差矩阵；步骤六，重复上述步骤一到步骤五，直至获得稳定的参数估计结果并用于修正参考星矢量。2.如权利要求1所述的一种基于模型参数估计的星光矢量蒙气差修正方法，其特征在于：所述步骤一中，星敏感器拍摄观测星图，所得第i颗观测星在星敏感器坐标系下的观测星矢量记为v
i
。3.如权利要求1所述的一种基于模型参数估计的星光矢量蒙气差修正方法，其特征在于：所述步骤二中，参考星坐标系转换时，星敏感器拍摄时间为曝光中心时刻，表示方法可选utc时间，精确到毫秒；地点为曝光中心时刻星敏感器的经度、纬度、高度，表示方法可选wgs84椭球参考系，精确到秒和米；地球自转信息包括地球当日极移和dut1时间，可由国际地球自转服务(international earth rotation service，iers)提供；基于国际天文联合会(international astronomical union，iau)基础天文学标准(standards of fundamental astronomy，sofa)提供的高精度坐标系转换链进行参考星坐标系转换，还需依赖导航星表上参考星的赤经、赤纬、自行、视差、径向速度等参数。观测点地理坐标系即为东北天坐标系，所得地理坐标系下第i颗参考星无蒙气差的矢量记为w
i
，天顶矢量u＝[0，0，1]
t
，第i颗参考星的真天顶距为：z
i
＝arccos(u
·
w
i
)
ꢀꢀꢀꢀ
(1)。4.如权利要求1所述的一种基于模型参数估计的星光矢量蒙气差修正方法，其特征在于，所述步骤三中，蒙气差经验模型如下：ρ＝a1tanz a2tan3z a3tan5z
ꢀꢀꢀꢀ
(2)其中，a1，a2，a3即为待估计参数，ρ为星光矢量的蒙气差大小，z为星光矢量的真天顶距；计算所有参考星矢量的蒙气差大小ρ
i
并对其沿仰角增加的方向做蒙气差修正，所得修正的参考星矢量都是关于待估计参数a1，a2，a3的函数，修正后的参考星矢量表示为的函数，修正后的参考星矢量表示为其中i表示单位矩阵，n
i
表示参考星矢量w
i
增加仰角ρ
i
时的转轴矢量：n
i
＝w
i
×
u
ꢀꢀꢀꢀ
(4)n
i∧
表示由矢量元素构成的反对称矩阵，n
i
(1)、n
i
(2)、n
i
(3)分别表示矢量n
i
的第1、2、3个元素：
5.如权利要求1所述的一种基于模型参数估计的星光矢量蒙气差修正方法，其特征在于：所述步骤四中，根据星光矢量夹角不变原理，星敏感器坐标系下任意两观测矢量的夹角等于观测点地理坐标系下两对应的修正后参考矢量的夹角，可得：若有n对匹配的星光矢量，则可构造个观测变量r
ij
：待估计参数可写为：可构造蒙气差参数估计量测方程：r
k
＝h
k
d
k
(9)其中，表示所有观测变量r
ij
构成的观测变量列，表示r
k
关于待估计参数a1，a2，a3的雅可比矩阵，表示待估计参数代差。6.如权利要求1所述的一种基于模型参数估计的星光矢量蒙气差修正方法，其特征在于：所述步骤五中，采用扩展卡尔曼滤波(ekf)算法解算蒙气差模型参数的滤波过程具体如下：(1)初始化参数a1，a2，a3，d0，p0；其中，d0为初始的d
k
(k＝0)；p0为初始误差协方差矩阵；(2)读取一幅星图中星点坐标按上述方法及公式完成对量测方程的观测变量r
k
以及系数矩阵h
k
的计算；
(3)计算卡尔曼滤波增益k
k
为：式中，e为星敏感器的测量噪声方差强度阵；(4)由观测变量r
k
更新d
k
为：d
k
＝d
k
‑1 k
k
[r
k
‑
h
k
d
k
‑1]
ꢀꢀꢀꢀ
(11)式中，d
k
‑1代表第k
‑
1张星图迭代出的蒙气差模型参数a1，a2，a3与第k
‑
2张星图迭代出的参数a1，a2，a3之间的差值；(5)更新参数a1，a2，a3的值：(a1，a2，a3)
k 1
＝(a1，a2，a3)
k
d
k
ꢀꢀꢀ
(12)(6)计算下一张星图的协方差矩阵p
k 1
并重置d
k 1
，返回(1)：p
k 1
＝(i
‑
k
k
h
k
)p
k
ꢀꢀꢀꢀ
(13)。7.如权利要求1所述的一种基于模型参数估计的星光矢量蒙气差修正方法，其特征在于：所述步骤六中，重复上述滤波过程直至获得稳定的参数估计结果即可认为完成了蒙气差模型参数估计过程，此后每帧星图所迭代得到的蒙气差模型参数a1，a2，a3可带入公式(3)完成当前参考星矢量蒙气差修正，进而使得星敏感器姿态解算更加准确。
技术总结
本发明公开了一种基于模型参数估计的星光矢量蒙气差修正方法，该方法可用于大气内星敏感器参考星矢量蒙气差修正。该方法将关于外部大气特征的蒙气差模型参数设为待估计值，根据星光矢量夹角不变原理，构造了蒙气差参数估计量测方程，通过扩展卡尔曼滤波算法循环迭代实时星图数据，估计出稳定且跟随环境变化的蒙气差模型参数，能够实时对参考星矢量进行蒙气差修正，进而提高大气内星敏感器姿态测量准确度。本发明实现了蒙气差模型参数的自主估计而不依赖其他外部传感器，采用扩展卡尔曼滤波算法实时估计出的蒙气差模型参数，不受传感器及经验模型误差影响。经验模型误差影响。经验模型误差影响。


技术研发人员：江洁 张晗 杨季三 张广军
受保护的技术使用者：北京航空航天大学
技术研发日：2021.03.25
技术公布日：2021/6/29