本文内容来源于《测绘学报》2024年第11期（审图号GS京(2024)2421号）

勾万祥, 李崇辉, 佟帅, 郑勇, 杨原, 陈智兴, 詹银虎

基金项目

国家自然科学基金(42374011)(42074013); 河南省优秀青年基金(242300421150)(242300421151)

摘要

白昼星图具有高背景噪声、低信噪比等特征，传统方法很难高效提取星点，且由于星点为弱小目标，对噪声极为敏感，成像随机不规则，单帧提取星点质心稳健性较差。为此，本文提出一种附加恒星视运动修正的星点提取方法，该方法首先根据测站、星历、观测时间预测恒星概略位置，通过相机参数转换模型得出星点在像平面坐标系中的概略位置；然后，利用各帧星图拍摄时刻精确计算出星点相对运动量，将各帧星图星点预测区域进行配准叠加；最后，在星点预测位置区域对星点进行精确提取。试验表明：该方法星点提取成功率可达100%，较传统算法提升163%以上，且算法执行效率高，耗时缩短至传统方法的23%；除此之外，多帧叠加增强了星图信噪比，克服了星点成像随机误差，星点质心提取精度较单帧星图提取方法提升72%以上。

关键词

白昼星图; 视运动; 星点提取; 质心精度; 天文大地测量; 多帧叠加

作者简介

第一作者：勾万祥（1990—），男，硕士，研究方向为天文导航、图像处理等。　 E-mail：memgwa201611@163.com

通讯作者:李崇辉　E-mail：lichonghui6501@126.com

本文引用格式

勾万祥, 李崇辉, 佟帅, 郑勇, 杨原, 陈智兴, 詹银虎. 附加恒星视运动修正的白昼星图星点提取方法[J]. 测绘学报, 2024, 53(11): 2166-2177 doi:10.11947/j.AGCS.2024.20230521

GOU Wanxiang, LI Chonghui, TONG Shuai, ZHENG Yong, YANG Yuan, CHEN Zhixing, ZHAN Yinhu. A method for star extraction from daytime star maps with additional correction of eastern elongation[J]. Acta Geodaetica et Cartographica Sinica, 2024, 53(11): 2166-2177 doi:10.11947/j.AGCS.2024.20230521

阅读全文

天文大地测量是通过测量恒星方位来确定站点天文坐标及方位的方法，在空间基准建立、惯导设备标定、垂线偏差测量等领域发挥着重要作用[1]。近些年，随着以铟镓砷（InGaAs）、碲镉汞（HgCdTe）等新型红外传感器的投入使用，天文测量逐步从传统的夜间测量走向全天时测量[2]。然而，受白昼复杂天空背景影响，白昼星图往往呈现噪声强、信噪比低、星点目标小等突出特点[3]，给白昼星图星点的高精度提取带来了一定的困难挑战[4]。

白昼星图星点的高精度提取是决定天文定位、定向精度的关键因素[5]。围绕白昼条件下准确探测提取星点目标，各国学者分别提出了基于单帧星图的星点提取算法，如背景预测法[6]、形态学法[7]、局部对比度法[8]、小波变换法[9]等，以及基于多帧星图的星点提取算法，如关联校验法[10]、一维处理法[11]、帧间差分法[12]等。这些算法虽然一定程度上可以提升星图信噪比，剔除伪星点，但大多数算法来源或继承于红外弱小目标探测，对于星点提取的精度关注较少，直接应用于高精度的全天时天文大地测量普遍存在以下3个问题：①在提升星图信噪比的同时，容易改变星点能量分布，使星点边界产生弱化、偏移等现象，导致星点提取精度下降[13-14]；②对于漂浮物、灰尘、杂散光等产生的随机伪星点，无法有效辨别区分，导致星图提取过程含有伪星点，需额外设计判定条件予以剔除[15-17]；③在近红外波段，星点能量较弱，受大气湍流、随机噪声等影响，各帧星图星点成像随机不规则，依靠单帧星图提取星点的精度和稳健性较差[18-19]。基于此，文献[20]、文献[21]、文献[22]分别提出了不同星点识别策略的多帧星图叠加法，此类方法虽然可以对静态星图实现信噪比增强，提高星点提取精度，然而实际应用中近地空间传感器会受地球自转影响产生恒星视运动现象，多帧星图直接叠加易产生星点拖尾现象，星点提取难以达到理想效果，且对于随机伪星点也未提出有效的去除方法。文献[23-24]提出基于星点位置预测的星点提取方法，该方法可以有效剔除伪星点，提高星点提取的成功率，但无法克服单帧星图星点成像随机性误差，恒星提取精度和稳健性较差，且随着星图信噪比降低，提取成功率也会受到影响。

为此，本文提出一种附加恒星视运动修正的白昼星图星点提取方法，该方法利用测站位置、恒星星历、星图拍摄时刻，求解待提取星点在星图中的概略位置，并利用不同帧星图的时间序列，精确求取不同帧星图中星点的相对运动量，基于此进行多帧星图星点预测区域配准叠加，最后在星点预测区域进行星点质心高精度提取。该方法逻辑简单、执行效率高，在大幅提高星点提取效率和准确率的同时，可以有效提升星点提取的精度和稳健性。

1　经典的多帧叠加法原理

白昼条件下，红外星图成像噪声种类复杂、繁多，但大体可分为固定背景噪声如云层、漂浮物、图形噪声等[25]和随机噪声如热噪声、散粒噪声、暗电流噪声等[26]。据此可将红外星图图像分解为红外星点、红外背景和成像随机噪声，具体模型如下[11]

(1)

式中，T（i，j）为红外星点；B（i，j）为红外背景；N（i，j）为成像随机噪声。设随机噪声为均值0、标准差δ的高斯噪声，则星图信噪比可定义为

(2)

根据误差原理可知，n帧星图叠加后噪声标准差为

(3)

由式（2）可知，叠加后的星图信噪比R′SNR为

(4)

由式（4）可知，当n帧星图叠加后，星图信噪比将明显增强，达到原星图的倍。如图1所示。

图1

图1   多帧叠加效果

Fig.1   Effect of multi-frame superposition

如图1（a）所示，白昼星图星点成像能量较弱易受大气湍流、随机噪声等影响，各帧星图星点成像形状随机不规则，不论采用何种质心提取算法，依托单帧星图提取的质心精度和稳健性均不理想。而将图1（a）前7帧星图叠加后，如图1（a）第8帧星图所示，星点形状明显平滑规则，星点提取精度和稳健性明显增强。对比图1（b）、（c）可以看出，单帧星图各种随机噪声大且能量分布不均匀，星图信噪比较低，而叠加后星点能量明显增强，随机噪声明显减弱，星图信噪比得到明显提升。

2　附加恒星视运动修正的多帧叠加法

虽然传统的多帧叠加法可以有效提升白昼星图信噪比，但由于近地空间传感器始终伴随地球处于运动状态，不同时间拍摄的恒星处于视场中不同位置，若直接将图像进行叠加必然引入星点视运动误差，导致叠加失败。对于天文大地测量而言，虽然不同拍摄时刻的恒星位置无法精确求得，但只要能精确求出恒星在不同时刻位置的相对变化量，即可准确实现星图叠加，增强信噪比。

2.1　恒星位置预测原理

恒星运动具有一定规律性，天文学家通过长期观测制作出包含恒星赤经赤纬（a，δ）信息的星历，通过查询星历可以准确计算出恒星在任意观测位置、任意时刻的地平位置[27]。

如图2所示，以测站O（λ，φ）为中心建立地平坐标系OXYZ，则测站天顶Z，北天极P和恒星σ可构成定位三角形ZPσ。通过定位三角形ZPσ可将测站经纬度（λ，φ）、恒星赤经赤纬（a，δ）、恒星方位角高度角（As，hs）建立联系。依据球面三角形知识可知，恒星方位角高度角（As，hs）可由式（5）计算[24]

(5)

式中，恒星时角t=s+λ-a，s为钟面时。由于恒星位置（a，δ）可通过恒星星历精确计算，钟面时s误差为ms量级，产生的计算误差可忽略不计，因此利用式（5）预测恒星地平位置的误差公式可表示为

(6)

式中，（Δλ,Δφ）为测站位置误差。在定位三角形ZPσ代入五元素公式可推得式（6）中误差项（ΔAs,Δhs）为

(7)

图2

图2   坐标系转换

Fig.2   Coordinate system transformation

如图3所示，恒星高度角hs是以测站为中心的大圆弧来表征，恒星方位角As是以恒星所在等高圈小圆中的弧来表征，而实际相机视场中水平方位是以测站为中心大圆弧来表征。因此式（7）中（ΔAs，Δhs）需等价转换为以测站中心的大圆弧角度（ΔFA，ΔFV），根据球面几何知识有

(8)

图3

图3   恒星方位转换

Fig.3   Star orientation transformation

顾及、、、，则对于测站位置误差产生的星点位置预测弧度角误差（ΔFA，ΔFV）有

(9)

实际应用中，使用式（5）预测星点成像位置精度还取决于大气折射模型误差hρ、相机初始定向误差（Δh，Δv）。因此采用式（5）对恒星位置预测产生的水平和竖直方向的视场角误差（ΔH，ΔV）可表示为

(10)

式中，经过模型改正后的大气折射误差影响优于10″；相机初始定向误差，一般可控制在5″量级；测站概略经纬度通常由GNSS或惯导计算得出，误差优于10″[23]。若相机视场为1.2°×1.2°，对应像元为1024×1024像素，则每个像元对应视场角约为4.2″×4.2″，因此使用式（5）模型对星点成像位置的预测误差可以控制在10个像素以内，约占为原图像的3.8×10-4。

由此可知，采用恒星位置预测方法在指定范围进行恒星提取，既可以有效剔除伪星点大幅提高恒星提取准确率、成功率，也可以明显降低图像运算量进而提升星图处理效率。

2.2　改进的多帧叠加法

图4为恒星相对位置。如图4所示，设测站O（λ，φ）在t1、t2时刻拍摄恒星位置分别为

(11)

图4

图4   恒星相对位置

Fig.4   The relative positions of star

则对于测站O（λ，φ）而言，任意时刻、任何方向的恒星预测位置（As，hs）均可通过式（5）计算得出；由于t1、t2时刻间隔较短，恒星高度角并不会发生明显变化，因此大气折射改正误差hρ可视为不变。故对于t1至t2时刻恒星实际位置相对运动量（ΔA，Δh）可表示为

(12)

式中，ΔA′=At2-At1、Δh′=ht2-ht1表示恒星预测位置的相对运动量，可以通过式（5）精确计算。因此恒星在t1至t2时刻实际位置运动量主要取决于ΔAt2-ΔAt1和Δht2-Δht1。对于ΔAt2-ΔAt1和Δht2-Δht1有

(13)

由于测站近天顶方位恒星无法参与测站方位解算，且在1 min时间内恒星最多运动15′，所以对于t1、t2时刻拍摄的恒星方位可以有以下限定条件

(14)

结合式（13）、式（14）可以进一步推导-和-在视场内对应的弧度角（ΔFA,ΔFV）有以下关系

(15)

对于测站而言，（Δλ，Δφ）误差不会超过10″，因此由式（15）可知和引起的角度误差（ΔFA，ΔFV）分别小于0.8″、0.08″，该量级对星图配准影响可以忽略不计。因此t1至t2时刻恒星实际位置运动量与预测位置运动量方位角、高度角可表示为

(16)

因此在进行星图配准时，只需利用式（16）计算恒星预测位置的运动量，即可实现星图的精确配准，其中横向配准精度优于0.8″，纵向配准精度优于0.08″。

2.3　算法流程

根据上述分析可知，基于恒星星历、测站位置、观测时间等先验信息计算的恒星概略位置一般优于35″（不超过10个像素），考虑白昼星点成像不会超过10×10像素，因此若以预测位置为中心取50像素为边长的正方形区域为预测区域，则恒星必然在预测区域内；根据各帧星图拍摄时间，精确计算星图中恒星运动量，并依此进行预测区域星图平移叠加，即可实现星图的精确配准叠加，具体操作流程如图5所示。

图5

图5   本文算法流程

Fig.5   Flowchart of the proposed algorithm

具体流程归纳如下。

（1）建立与测站地平坐标系相关联的图像坐标系，并精确标定相机坐标转换模型参数[28]，然后对相机进行初始定向。

（2）依据测站位置、观测时刻，计算目标恒星在测站地平坐标系下的概略视位置（高度角HS，方位角AS）；驱动相机转至目标恒星方位（HS，AS），连续拍摄2N+1帧星图，精确记录每帧星图拍摄时刻为t1、t2、…、t2N+1，利用式（5）计算拍摄时刻恒星视位置为（H1，L1）、（H2，L2）、…、（H2N+1，L2N+1）。

（3）根据相机坐标转换参数计算各时刻恒星视位置（H1，L1）、（H2，L2）、…、（H2N+1，L2N+1）对应的图像坐标系下坐标（x1，y1）、（x2，y2）、…、（x2N+1，y2N+1），在图像坐标系中以（x1，y1）、（x2，y2）、…、（x2N+1，y2N+1）为中心，取边长为50像素的正方形区域作为目标星点预测区域。

（4）以第N帧时刻恒星预测坐标（xN，yN）为基准，分别计算各帧星图星点预测位置相对运动量（x1-xN，y1-yN）、（x2-xN，y2-yN）、…、（x2N+1-xN，y2N+1-yN），并将所有星图星点预测区域以第N帧星图时刻为基准进行位移配准。

（5）将2N+1帧预测区域星图进行叠加，采用自适应阈值分割法提取星点，并进行星点质心细分计算[29]，将最终质心坐标作为第N帧时刻星图的星点质心坐标。

（6）利用第N帧星图星点质心坐标和相机坐标系转换参数，计算恒星在第N帧星图拍摄时刻的精确地平方位。

3　试验与分析

3.1　试验设备

为验证本文算法的有效性，利用全天时红外测星系统进行白昼恒星观测。该系统由全天时红外全站仪、遮光罩、原子钟、GNSS接收机、数据处理系统等设备构成。其中GNSS接收机搭配原子钟进行授时守时，时间精度优于0.01 ms；全天时红外全站仪采用H波段为观测波段，在正午时段可以实现H波段0等星有效观测，具体光学参数见表1。

表1   红外全站仪光学系统参数

Tab.1  Optical system parameters of infrared total station

系统参数	参数值
视场角大小	1°20′
通光口径/mm	40
光学透过率/（%）	56
等效焦距/mm	240
工作波段/nm	1500~1700
传感器型号	QHY990
量子响应效率/（%）	77
像元尺寸/μm	5.0
像元总数	1392×1052
暗电流/（像素/s）	≤227e-

新窗口打开| 下载CSV

测量开始前，首先使用GNSS接收机接收测站概略位置，并同步系统时间；然后使用已知方位边对全天时红外全站仪进行精确定向；最后通过数据处理系统根据测站位置、恒星编号、预计观测时间计算出待测恒星概略方位，驱动红外全站仪寻星定位并进行星图拍摄，星图拍摄时同步记录拍摄时间。具体全天时红外测星系统如图6所示。

图6

图6   全天时红外测星系统

Fig.6   All-time infrared star measurement system

3.2　试验数据

2023年8月28日11时25分许，郑州某楼顶利用红外全站仪寻星至HIP 84345（H波段-3.22等，方位96°46′，65°44′），HIP 80704（H波段-1.85等，方位47°48′，13°13′），HIP 54061（H波段-0.61等，方位19°23′，59°50′），设置曝光时间为10 ms，每隔160 ms拍摄一帧星图，连续拍摄100帧星图，其中第1、50、100帧星图如图7（a）、（b）、（c）所示。

图7

图7   不同时刻拍摄星图

Fig.7   Star map taken at different times

对星图采用传统自适应阈值分割方法进行星点提取，并对比分析各帧星图可知：

（1）图7中，红色框为提取出的目标星点，蓝色框为伪星点，其中图7（a）、（b）采用传统算法虽能提取出目标星点，但星图中夹杂大量伪星点，常规方法很难区分；图7（c）采用传统方法只能提取出1个伪星点，无法提取出目标星点。说明白昼星图受传感器内部噪声、盲像元及外部背景灰度不均、漂浮物等影响，普遍存在高噪声、低信噪比等特点，传统阈值分割方法很难准确提取出目标星点。

（2）本次拍摄星图第1帧与第100帧拍摄时间相差约17 s，从各帧星图星点坐标可知恒星不同时刻处于图像坐标系中不同位置，其中图7（a）拍摄的HIP 84345恒星，星点位移达29个像素，平均移动速度为1.7像素/s；图7（b）拍摄的HIP 80704恒星，星点位移达56个像素，平均移动速度为3.3像素/s；图7（c）拍摄的HIP 54061恒星，星点位移也有8个像素，平均移动速度为0.5像素/s。说明受恒星视运动影响，恒星在图像坐标系具有明显的移动现象，若将不同时刻的星图直接叠加必然造成配准误差，导致星点提取精度下降。

（3）分别对HIP 84345第1、50、100帧星图星点区域进行放大，观察红外全站仪实际拍摄的星点形状及能量分布，如图8所示。由图8可知，在近红外波段恒星能量相对较弱，受各种随机噪声影响，各帧星图星点成像大小不一、形状不定，如第1帧星点为单峰型、第50帧星点为扁平型、第100帧星点为双峰型，即使星点处于同一位置，采用同一种星点质心提取算法也会产生0.5个像素左右的误差，这在高精度天文测量中是不容忽视的。

图8

图8   各帧星图星点形状

Fig.8   Star shape of different frame star map

3.3　试验结果

以第50帧时刻恒星预测位置为基准，利用式（16）将100帧星图全部归算至第50帧星图像平面坐标系中，归算精度优于0.01个像素，这对于分析本文算法精度影响可以忽略不计。在完成100帧星图配准的基础上，分别采用以下3种方法对100帧星图星点进行提取。

（1）传统方法：通过调整阈值和连通域参数，采用窗口大小为30×30像素的局部自适应算法对各帧星图进行阈值分割，采用灰度加权法进行质心计算，记录星点提取成功率和各帧星点坐标。

（2）位置预测法：根据每帧星图拍摄时间，通过式（5）计算恒星概略地平位置，并通过图像参数模型，归算至红外全站仪视场中像素位置。考虑到恒星大小小于10×10像素，预测误差不超过10个像素，取恒星预测位置为中心，50个像素为边界的正方形区域作为星点预测区域，并进行星点提取，记录星点成功率和各帧星点坐标。

（3）本文方法：将配准后的100帧星图按照每7帧为一组，随机组合出100组星图。取恒星预测位置为中心，50个像素为边界的正方形区域作为星点预测区域，将星点预测区域图像进行叠加并进行星点提取，记录星点成功率和各帧星点坐标。

根据各种方法提取星点具体坐标分布绘制散点图如图9所示。

图9

图9   3种方法提取星点位置分布

Fig.9   Extracting star position distribution by three methods

由图9可知，3种方法提取出的星点位置分布存在明显不同，其中传统算法和位置预测法单帧星点提取互差最大达2个像素；而本文方法星点提取最大互差仅为0.8个像素左右。同时本文所提方法与位置预测法提取坐标均值相差较小，HIP 84345、HIP 80704、HIP 54061 3颗恒星分别相差（0.02，-0.18）、（-0.13，-0.09）、（0.13，-0.06）个像素，说明利用多帧星点质心位置取平均方法可以有效克服白昼星点成像随机性误差，所提取星点位置具有更好的稳健性；而传统方法由于星点提取成功帧数较少，参与求平均的帧数也较少，星点质心随机性误差明显提升，其中HIP 84345、HIP 80704恒星分别相差达（-0.47，-0.26）、（-0.87，-0.25）个像素，充分说明依托单帧提取星点质心产生的误差在高精度天文测量中是不容忽视的。

分别统计3种方法提取星点的精度、成功率及时效见表2。

表2   3种方法提取星点结果

Tab.2  The results of extracting star points by three methods

采用方法	恒星编号	X坐标均值/像素	Y坐标均值/像素	均方差/像素	所用时间/ms	提取成功帧数	成功率/（%）
	HIP84345	83.04	141.40	0.57	2859	38	38
传统方法	HIP 80704	308.45	254.51	0.63	3043	7	14
	HIP 54061	—		—	2810	0	0
	HIP84345	82.55	141.32	0.43	438	100	100
位置预测法	HIP 80704	307.71	254.35	0.46	442	100	100
	HIP 54061	172.28	253.28	0.54	441	87	87
	HIP84345	82.57	141.14	0.12	656	100	100
本文方法	HIP 80704	307.58	254.26	0.11	681	100	100
	HIP 54061	172.41	253.22	0.15	652	100	100

新窗口打开| 下载CSV

分析表2数据可知，在提取成功率方面，采用传统方法，虽然经过不断调整阈值和连通域参数，但对于HIP 84345、HIP 80704恒星最终也仅获得38%、7%的提取成功率，而位置预测法和本文所提方法均可达100%，这主要因为部分伪星点无法准确辨别剔除。随着拍摄场景、方位不同，星图构成可能产生明显变化，不同方法的星点提取成功率也会产生明显变化，尤其对于能量较弱的恒星。如图10所示，对于星等较高的HIP 54061恒星，依托传统方法根本无法提取出有效星点，而采用位置预测法在指定区域进行提取，可以明显提高成功率，但由于星图信噪比过低，且受大气湍流影响，星点呈现闪烁现象，提取成功率也有所下降，仅为87%，而本文方法克服了单帧星点信噪比过低影响，有效提升了星点信噪比，提取成功率仍保持100%，较传统方法提升163%以上。

图10

图10   HIP 54061恒星能量分布

Fig.10   Energy distribution of HIP 54061 star

从提取时间来看，由于传统方法需要对全图像进行分割、计算、提取，涉及区域大，计算目标多，3颗恒星分别为2859、3043、2810 ms；而采用位置预测法，每帧星图仅计算预测区域约为原图的1.7×10-3，星点提取运算时间明显缩短，时间为438、442、441 ms，约为传统方法的15%，而采用本文方法，虽然进行7帧预测区域叠加，但计算区域较小且为加法运算，实际时间较位置预测法有所增加但远小于传统方法，时间为656、681、652 ms，约为传统方法的23%。

从质心提取精度来看，由于传统方法提取成功率较低，因此计算精度较单帧位置预测法有所下降，但总体相差不大，对于HIP 84345、HIP 80704恒星分别为0.57、0.63，0.43、0.46 pixel，说明单帧星点成像不规则导致星点提取精度受偶然性误差影响明显。而对于HIP 84345、HIP 80704、HIP 54061不同星等的3颗恒星，本文方法始终保持较高质心提取精度，分别为0.12、0.11、0.15 pixel，较传统单帧提取算法0.43、0.46、0.54分别提升了72%、76%、72%，说明本文提出的附加恒星视运动修正的白昼星图星点提取方法，可以有效克服了单帧星点成像的随机性误差，对各种信噪比下的星图均能维持较高的质心提取精度。

4　结论

白昼星图星点高精度提取是实现全天时天文测量的关键技术之一。本文提出了一种附加恒星视运动修正的白昼星图星点提取方法，试验结果表明：

（1）基于星点概略位置预测的星点提取策略可以精准识别待测星点，有效剔除伪星点，大幅提高星图星点提取成功率，所测3颗恒星提取成功率均为100%，较传统方法提升163%以上。

（2）基于星点预测区域的叠加策略大幅减少星图运算区域，有效提升了算法的执行效率，叠加7帧星图的运算时效约为传统方法的23%，且算法逻辑简单易于编程实现。

（3）附加恒星视运动修正的多帧叠加方法可以实现序列星图的精准叠加，在提高星点信噪比的同时，大幅提升了星点提取的精度和稳健性，整体星点提取精度较单帧提取方法提升72%以上。

初审：陈品馨

复审：宋启凡

终审：金   君