基于区域滤波的模糊星图复原方法

王 军，何 昕，魏仲慧，穆治亚，吕 游，何家维

(1.中国科学院 长春光学精密机械与物理研究所，吉林 长春 130033；2.中国科学院大学，北京 100049)

1 引 言

星敏感器是一种以恒星为测量对象的高精度、高可靠性的姿态敏感测量仪器，目前已经在深空探测以及天文导航等空间任务中得到了广泛应用。其工作原理为：由图像传感器对星空成像，经过星点提取和质心定位得到星点在传感器靶面上的位置和亮度信息，之后利用星图识别找到星点在星表中相对应的导航星，最后根据识别结果解算出星敏感器相对于惯性坐标系的三轴姿态[1-3]。星敏感器通常工作在静态条件下，多用于载体稳态飞行，即假设在曝光时间内，导航星和星敏感器相对静止，星点只在传感器靶面上的固定位置成像。当星敏感器工作于载体初始入轨、机动、调姿等动态场合时，在曝光时间内星点在传感器感光区内移动，最终在传感器靶面上形成一段轨迹图像[4-5]，导致图像信噪比和星点质心定位的精度降低，严重影响星敏感器的姿态测量精度[6]。因此动态性能作为星敏感器的关键指标之一，已经成为星敏感器研究领域的重点内容。在动态条件下，由于图像的信噪比降低以及星点运动模糊现象的存在，高精度的质心定位已经成为动态星敏感器研究的主要内容之一[7]。

目前，许多研究机构都已经提出了提高动态条件下星点质心定位精度的方法，这些方法主要可以分为两类：一类是硬件增强法，另一类是软件复原法。硬件增强法通过对星敏感器的硬件电路进行改进，减轻或抑制星点运动模糊的程度，从而提高星点质心定位的精度。Bezooijen等人[8]提出了一种延时积分法(Time Delayed Integration,TDI)，该方法通过设计特殊的硬件时序电路对星点的运动模糊进行抑制，一定程度上提高了图像的信噪比。但是该方法增加了硬件时序电路的复杂度，减小了传感器的有效面积。此外，该方法只消除了y方向上的运动模糊，在此基础上，仍需通过设计特殊的图像处理算法，进一步提高信噪比。Pasetti等人[9]在TDI技术的基础上，利用过采样和多像素合成的方法补偿运动的影响，但是该方法的运动参数仅限于水平、垂直和对角线方向，运动参数误差很大，只是一种近似估计的方法。软件复原法利用图像复原算法对模糊星图进行复原，消除星点的运动模糊现象，从而提高星点质心定位的精度。Sun等人[10]首先根据从陀螺仪获取的角速度信息建立星点质心的运动模型，进而得到模糊星图的退化函数，最终利用RL(Richardson Lucy)算法[11-12]对模糊星图进行复原。但是在实际应用的过程中，陀螺仪提供的角速度数据会随时间发生漂移，严重影响星点质心运动模型的精度，因此，Sun等人[13]采用扩展卡尔曼滤波(Extended Kalman Filter,EKF)的方法对陀螺仪数据进行修正，从而提高了星点质心运动模型的精度。此外，虽然Sun等人提出的RL算法对模糊星图具有很好的复原效果，但是一方面RL算法是一种无终止条件的迭代方法，往往只能根据经验选择迭代次数，如果选择不当会对最终的复原结果造成很大的影响。另一方面，RL算法对噪声相当敏感，在对模糊星图进行复原的过程中容易放大噪声，严重影响图像的复原效果，因此，Ma等人[14]在对模糊星图进行复原之前先采用多子区域生长法对图像进行预处理，但是该算法选用的模板尺寸具有局限性，容易造成星点能量的丢失。此外，Zhang[15]和Sun[16]分别采用自适应小波阈值法以及相关滤波和形态学滤波相结合的方法对模糊星图进行预处理，但是由于星图中星点所占的像元数相对较少，对模糊星图进行全局滤波往往会造成系统资源的浪费，导致星敏感器的数据更新率降低，最终影响星敏感器的实时性。

文中针对目前大多主流方法在提高动态条件下星点质心定位精度时存在的不足，分析了星敏感器在不同工作条件下星点质心的运动模型和模糊星图的退化函数，并介绍了一种基于运动模型的区域滤波算法，采用该方法对模糊星图进行预处理，并在此基础上利用复原算法对图像进行复原。

2 模糊星图的参数估计

2.1 星点质心的运动模型

星敏感器是一种具有角秒级测量精度的姿态测量仪器，其姿态测量模型如图1所示[17]。

图1 星敏感器的姿态测量模型Fig.1 Attitude measurement model of star sensor

在惯性坐标系oc-xcyczc中，星点的参考矢量uc可以用星点对应的赤经α和赤纬δ来表示，即：

(1)

在像空间坐标系o-xyz中，星点的观测矢量ws可以用星点在传感器靶面上的坐标p(x,y)和光学系统的焦距f来表示，即：

(2)

理论上，参考矢量uc和观测矢量ws满足下列的等式关系：

ws=Auc,

(3)

式中，A为像空间坐标系和惯性坐标系之间的旋转矩阵，也被称作星敏感器的姿态矩阵，利用姿态矩阵可以解算出星敏感器相对于惯性坐标系的三轴姿态。

假设在t0时刻，星点在传感器靶面上的质心坐标为(x(t0),y(t0))，对应的观测矢量为ws(t0)，在t0+Δt时刻(Δt<

(4)

由于在t0和t0+Δt时刻，星点在惯性坐标系中具有相同的参考矢量uc，则有：

(5)

(6)

式中，w=[wxwywz]T为星敏感器的角速度矢量，通常由与其捷联的惯性导航系统提供，w×为与w正交的三阶方阵。

假设在极短的时间间隔Δt内，星敏感器的角速度几乎不发生改变，则通过解算公式(5)，可以得到t0至t0+Δt时刻星点质心的运动模型：

(7)

由于时间间隔Δt在数值上远小于星敏感器的焦距f，则[x(t0)wyΔt-y(t0)wxΔt]/f=1成立，因此星点质心的运动模型可以近似表示为：

(8)

2.2 模糊星图的退化函数

在公式(8)的基础上，分析星敏感器在两种不同工作条件下星点质心的运动模型和模糊星图的退化函数，其中一种工作条件是星敏感器只绕z轴旋转，另一种工作条件是星敏感器同时绕x轴和y轴旋转。

(1)星敏感器只绕z轴旋转

当星敏感器只绕z轴旋转时，x轴和y轴的角速度满足wx=wy=0，将其代入到公式(8)中，可得：

(9)

解算公式(9)，则有：

x2(t)+y2(t)=c2,

(10)

(11)

(2)星敏感器同时绕x轴和y轴旋转

当星敏感器同时绕x轴和y轴旋转时，z轴的角速度满足wz=0，将其代入到公式(8)中，可得：

(12)

解算公式(12)，则有：

(13)

(14)

3 区域滤波

当星敏感器在动态条件下工作，曝光时间内星点在传感器感光区内移动时，星点因能量的分散而出现运动模糊现象，导致图像的信噪比降低并且星点的模糊区域很难被提取，因此文中提出了一种区域滤波算法对模糊星图进行预处理。首先，根据星点质心的运动模型得到星点质心的运动轨迹，然后确定星点的模糊区域，如图2所示。从图2可以看出，模糊区域是由星点质心运动轨迹窗口向四周各外扩一行/一列像元得到的。这是因为在静态条件下，由于质心计算的需要，星敏感器在像面调整时需要采用离焦技术使无穷远的星光信号弥散到大约3×3个像元上[18]，形成传感器靶面上的星点，星点的灰度呈高斯点扩散函数分布，中心像元的灰度值最高，并且星点的质心通常位于中心像元内[19-20]。因此为了尽可能地保留星点的能量，文中选用图2中模糊区域的确定方法。

图2 星点的模糊区域Fig.2 Blurred region of star

由于星图中星点所占的像元数相对较少， 并且模糊区域外的图像中不包含星点的其它信息，因此为了在有效提高图像信噪比的同时提高算法的处理速度，文中将模糊区域外图像的灰度值置0，即：

(15)

式中，I(x,y)为图像中(x,y)坐标处的灰度值，Rk为图像中第k个星点的模糊区域。

在模糊区域内，文中采用形态学滤波中的开操作对图像进行滤波[21]。开操作常用于消除与结构元素相比尺寸较小的亮细节， 而保持图像整体灰度值和大的亮区域基本不变，有利于避免对星点的能量造成损失，其原理为先利用结构元素b对图像f进行腐蚀，即：

(16)

再对腐蚀后的图像进行膨胀，最终得到开操作滤波后的图像f2：

f2=f1⊕b=

(17)

文中b选用元素值都为1，大小为3×3的方形模板。

4 模糊星图复原

图3 图像的退化/复原模型Fig.3 Degradation and restoration model of image

文中采用RL算法对模糊星图进行复原，RL算法是图像复原领域的经典算法之一，其假设图像服从Poission分布，采用极大似然法进行估计，是一种基于贝叶斯分析的迭代算法，需要的先验知识少，复原效果好。其迭代方程为：

(18)

式中，*代表卷积运算，⊕代表相关运算，g为模糊星图，h为公式(11)和公式(14)卷积运算后的结果，即h(x,y)=hz(x,y)*hxy(x,y)，f(k)和f(k+1)分别为迭代k次和k+1次后得到的复原图像。

从公式(18)中的迭代方程可以看出，RL算法的迭代过程不存在终止条件，往往只能根据经验选择迭代次数k。如果选择的k值过小，模糊星图复原的不充分，复原后的图像无法满足星点质心定位精度的要求，反之，如果选择的k值过大，不仅会增加复原过程的处理时间，降低星敏感器的数据更新率，影响星敏感器的实时性，还会导致噪声在迭代过程中被放大，严重影响星点质心定位的精度。因此，文中对RL算法进行了改进，利用连续三帧复原图像中星点的质心坐标和模糊星图的参数构建迭代方程的终止条件，即：

(19)

式中，(xk-1,yk-1)，(xk,yk)和(xk+1,yk+1)分别为第k-1，k和k+1次迭代后复原图像中星点的质心坐标，Δx和Δy分别为t0至t0+T时刻星点质心坐标的差值，ε(Δx)和ε(Δy)分别为关于Δx和Δy的方程，需要在权衡星点质心的定位精度和复原算法的处理时间之后选择合适的表达式。利用该方法不仅可以在复原图像满足终止条件时及时停止迭代过程，同时，由于使用了连续多帧复原图像，该方法还具有很强的鲁棒性。

5 实验结果与分析

在处理器为Intel Core i5-7300HQ，主频为 2.5 GHz的计算机上利用Matlab进行仿真实验，仿真过程中选用的星敏感器参数如表1所示，实验选用J2000星表，星敏感器的角速度的变化范围为0～2°/s。

表1 星敏感器参数

5.1 模糊星图去噪实验

利用计算机仿真得到的静态条件下的模拟星图如图4(a)所示，图像的分辨率为512×512，图中包含7个星点，在wx=wy=2°/s的条件下，根据星敏感器的参数对图4(a)进行模糊处理，得到动态条件下的模拟星图如图4(b)所示。在图4(a)和图4(b)中分别截取7号星点模糊前后的局部图像，图像的分辨率为128×128，模糊前后星点的灰度分布如图4(c)和图4(d)所示(为了便于显示将星点进行了放大)。对比图中星点的灰度分布可以看出，在动态条件下，星点的能量被分散，星点的灰度值明显降低。

图4 静态和动态条件下的模拟星图及其灰度分布Fig.4 Simulated star images and gray distribution under static and dynamic conditions

在图4(d)中依次加入方差为30、50、70的随机分布的高斯白噪声，分别采用自适应小波阈值法，开操作和区域滤波算法对图像进行去噪，并计算去噪后图像的峰值信噪比，结果如表2所示。实验结果表明，在不同噪声方差的条件下，3种算法都可以有效提高图像的峰值信噪比，但区域滤波算法对应的峰值信噪比明显高于其它两种算法。据此可知，与自适应小波阈值法和开操作相比，区域滤波算法对模糊星图具有更好的去噪效果。

表2 去噪算法在不同噪声方差下的峰值信噪比

为了验证区域滤波算法的处理速度，文中在实验的过程中分别记录了3种算法在不同噪声方差条件下的处理时间，如图5所示。实验结果表明，在不同噪声方差的条件下，区域滤波算法的处理速度明显快于自适应小波阈值法和开操作。由于自适应小波阈值法的复杂度较高，与其它两种算法相比其处理时间相对较长，而区域滤波算法只对星点模糊区域内的图像进行滤波，因此与全局滤波的开操作相比，其处理时间明显更短。

图5 去噪算法的处理时间Fig.5 Processing time of each method

5.2 模糊星图复原实验

在图4(d)中加入方差为30的随机分布的高斯白噪声，采用区域滤波算法对图像进行去噪，然后利用RL算法对去噪后的图像进行复原，得到的复原图像及其灰度分布如图6所示。从图6可以看出，复原图像中星点的能量分布更加集中，但是星点的灰度值与图4(c)相比有所降低，可能是在预处理过程中，去噪算法对星点的能量造成了损失所引起的。

图6 复原后的星图及其灰度分布Fig.6 Restored star image and gray distribution

为了进一步验证区域滤波算法对模糊星图复原效果的影响，文中采用去噪效果相对较好的开操作与其进行对比，分析两种去噪算法复原后图像中星点质心定位精度。在图4(b)中加入方差为30的随机分布的高斯白噪声，分别采用开操作和区域滤波算法对图像进行去噪，然后利用RL算法对去噪后的图像进行复原，最后利用质心法计算复原图像中各个星点的质心坐标，如表3所示。从表3可以看出，与开操作相比，采用区域滤波算法所对应的星点质心坐标更加接近真值。为了进一步对比两种去噪算法对复原图像中星点质心定位精度的影响，在图4(b)中依次加入方差为30，50，70的随机分布的高斯白噪声，分别采用开操作和区域滤波算法对图像进行去噪，然后利用RL算法对去噪后的图像进行复原，最后利用质心法计算复原图像中各个星点的质心坐标，并分别计算不同噪声方差条件下两种去噪算法所对应的星点的质心偏差|δx|+|δy|，如图7所示。

表3 复原星图中星点的质心坐标

图7 去噪算法在不同噪声方差条件下对应的星点的质心偏差Fig.7 Star centroid error for different methods under various noise variances

从图7可以看出，在同一噪声方差条件下，与开操作相比，区域滤波算法对应的星点的质心偏差明显更小，当噪声的方差发生变化时，区域滤波算法对应的星点质心偏差变化比较稳定，并且基本保持在0.1 pixel以内，而开操作对应的星点质心偏差变化幅度较大。因此，与开操作相比，采用区域滤波算法预处理后对应的复原图像中星点质心坐标的精度更高并且稳定性更好。

通常，当星敏感器在动态条件下工作时，星敏感器的角速度并非是恒定不变的，而是按照一定的规律不断变化的[25]。因此，为了验证本文方法在星敏感器角速度变化时的有效性，文中在0～2°/s的范围内随机选取了A～D 4组大小依次递增的角速度矢量，分别代表星敏感器在连续4帧中4种不同的运动状态，如表4所示。

表4 星敏感器的角速度

将图4(a)作为第一帧的原始星图，首先，在角速度矢量A的条件下，根据星敏感器的参数对原始星图进行模糊处理，并在模糊星图中加入噪声方差为30的随机分布的高斯白噪声，然后采用区域滤波算法对图像进行去噪，利用RL算法对去噪后的图像进行复原，最后利用质心法计算复原图像中每个星点的质心坐标，并与原始星图中对应星点的质心坐标进行对比，得到第一帧中7个星点的质心偏差|δx|+|δy|，将得到的复原图像作为下一帧的原始星图，改变角速度矢量，重复上述实验过程直至第四帧结束。将噪声方差依次改为50和70并重复上述实验，实验结果如图8所示。

图8 角速度变化条件下星点的质心偏差Fig.8 Star centroid errors under conditions of variable angular velocity

从图8可以看出，当星敏感器的角速度恒定时，噪声的方差越大，噪声对模糊星图的干扰越严重，图像复原后星点的质心偏差越大；在同一噪声方差条件下，星敏感器的角速度越大，星点的运动模糊越严重，图像复原后星点的质心偏差越大，尽管如此，当星敏感器的角速度在0～2°/s的范围内变化时，4帧复原图像中星点的质心偏差不超过0.1 pixel，可以满足星敏感器对高质心定位精度的要求。

6 结 论

针对动态条件下星图的信噪比降低导致星点质心定位精度降低的问题，提出了一种基于区域滤波的模糊星图复原方法。首先根据星敏感器的工作特性，建立不同工作条件下星点质心的运动模型，然后根据运动模型确定星点质心的运动轨迹，进而提取星点的模糊区域，再对模糊区域内外的图像分别进行预处理，最后利用复原算法对模糊星图进行复原。本文进行了模糊星图去噪实验。结果表明在2°/s的动态条件下，区域滤波算法的去噪效果和处理速度都明显优于两种经典的去噪算法，接着进行了模糊星图复原实验，结果表明采用区域滤波算法预处理后对应的复原图像中星点的质心偏差基本保持在0.1 pixel以内，最后在星敏感器角速度变化的条件下进行了多帧星图复原实验，结果表明当星敏感器的角速度在0～2°/s的范围内变化时，采用本文方法得到的复原图像中星点的质心偏差不超过0.1 pixel，可以满足星敏感器对高质心定位精度的要求。