[发明专利]暗弱小天体光学导航特征成像模拟方法

权利要求书

1.暗弱小天体光学导航特征成像模拟方法，其特征在于：包括如下步骤，

步骤一、给定暗弱小天体探测器任务远距离接近段长曝光下产生的拖曳星图的光学系统参数；

给定暗弱小天体探测器任务远距离接近段长曝光下产生的拖曳星图的光学系统参数，所述光学系统参数包括CCD相机的光圈半径A、视场角FOV、像素尺寸px/py、拍摄曝光时间T、光学效率η、量子效率Q和ADC增益K；

步骤二、给定暗弱小天体探测器拖曳星图的背景噪声参数，并根据给定的背景噪声参数生成噪声背景；

给定暗弱小天体探测器拖曳星图的背景噪声参数，所述背景噪声参数包括读数噪声偏置值λ_bias、暗电流噪声泊松参数λ_dark、成像背景亮度Sky_level，并根据给定的背景噪声参数生成噪声背景；

步骤三、根据目标小天体所在位置确定光轴指向，检索星表确定在该指向下能够进入视场范围的导航恒星列表，并计算所述导航星在CCD平面上的位置；

步骤四、给定PSF能量扩散半径，根据背景导航星及目标小天体的视星等调用高斯扩散公式得到各星点散焦后的能量；

步骤五、随机生成一系列光轴指向扰动，计算由该扰动产生的星图拖曳轨迹的能量，生成拖曳星图，即实现模拟暗弱小天体探测器任务远距离接近段长曝光下产生的拖曳星图；

步骤六，将步骤五得到暗弱小天体探测器拖曳星图用于验证轨迹图像星点质心提取的图像处理算法的有效性，进行图像处理算法的测试；

步骤七，改变光学系统参数，通过步骤一至步骤五模拟不同光学系统参数下的暗弱小天体探测器拖曳星图，进而选定满足实际深空探测任务的光学系统参数，实现实际深空探测任务；

步骤一的具体实现方法为，

CCD相机的感应器件将入射光子转化为电子并将生成的电子贮存起来，随后，读入器件测量每个像素中的电荷量在电容器两端产生模拟电压，模数转换器将每个像素中的电压数字化，产生一个计数写入电脑磁盘中；最后将电脑磁盘中的计数转化为灰度，生成灰度图像；

每个像素点接收到的入射光子数决定该像素点处信号的大小；当光学系统入瞳处的光通量为F时，CCD处接受到的光子数P为：

其中，为平均波长，h为普朗克常数，c为光速，A为光圈半径，η为光学效率；

入射光子数P与生成光电子数C的比率为光学系统的量子效率Q：

C＝Q*P (2)

光电子数C与计数值N之间的关系由ADC的增益K所确定：

N＝K*C (3)

星载光学系统生成的图像为8位灰度图像，有256级灰度，计数N和灰度DN间的转换关系为线性关系：

DN＝a*N+b (4)

步骤二的具体实现方法为，

总噪声包括读数噪声、暗电流噪声、成像背景亮度噪声：

S_噪＝I_background+I_read+I_dark (5)

由于深空中星际尘埃、太空碎片及其他物质的影响，拍摄的拖曳星图中存在大量杂散光，背景并不是完全的黑色；模拟星图的背景由用户输入参数背景亮度Sky_level和曝光时间T所决定，服从泊松分布，背景信号强度写为：

I_background～P(Sky_level*T) (6)

读数噪声是读数器件自身的电压所带来的噪声，所述读数噪声服从关于偏置值的高斯分布，其信号强度写为：

I_read～N(λ_bias) (7)

暗电流噪声是由CCD的结构内热造成的，其在单位时间内的噪声水平与CCD温度相关，属于散粒噪声的一种，服从泊松分布，其信号强度写为；

I_dark～P(λ_dark) (8)

根据给定的背景噪声参数生成噪声背景；

步骤三的具体实现方法为，

定义光学相机的视轴方向为(α₀,δ₀)，则出现在光学相机视场内的导航目标所需满足的赤经赤纬条件如下：

其中，(α_i,δ_i)为导航恒星的赤经赤纬，FOV为视场角大小；

考虑天球赤道坐标系和相机坐标系之间的坐标变换；由于导航星距离遥远，将相机坐标系的原点O和天球赤道坐标系的原点O'视为重合，忽略坐标系平移所带来的误差,仅需要考虑旋转变换；导航恒星在天球赤道坐标系下的坐标为(α_i,δ_i)，天球赤道坐标系O-UVW和相机坐标系O'-XYZ间的转换关系为：

[X Y Z]^T＝M[U V W]^T (10)

[U V W]^T＝[cosα_icosδ_i sinα_icosδ_i sinδ_i]^T (11)

所述相机坐标系O'-XYZ的Z轴为视轴正向，X轴和Y轴方向与CCD相机平面方向一致；

其中，j角与CCD相机安装方向有关，取j＝90°，θ为视轴与OW轴的夹角，O为O'X轴与OU轴的夹角，θ＝90°-δ₀，O＝90°+α₀；因此，M矩阵简化为：

最后，考虑相机坐标系和CCD成像平面坐标系之间的坐标变换；根据光学相机的成像模型，得到相机坐标系和CCD平面坐标系之间的几何关系：

上式中，(x,y)为导航目标在CCD成像平面投影的坐标，f为相机焦距；

同时，根据CCD相机视场角大小、传感器尺寸和焦距间的关系知：

f＝2pxNx^-1tan(FOV/2)＝2pyNy^-1tan(FOV/2) (14)

其中，Nx,Ny是CCD平面尺寸，FOV为视场角大小；px,py为像素的高和宽；

最终得导航目标在CCD成像平面投影的坐标位置的表达式：

所述导航目标在CCD成像平面投影的坐标位置即为导航星在CCD平面上的位置；

步骤四的具体实现方法为，

在天文观测中，视星等与单色光流间的关系为：

m_v＝2.5log₁₀F_v+m_c (16)

其中，m_v为信源视星等，F_v为信源单色光流，m_c为常数；

由公式(17)得出，信源的成像强度主要由信源的星等决定，信源视亮度与信号强度之间的关系为：

其中，I为某像素点上的光能，I_m为0等星所对应的光能，mag为某像素点对应的星等；

由于深空导航目标为点光源，需要考虑在CCD焦平面上造成的点扩散问题，在输出的数字灰度图上实际表现为包含若干像素点的光斑；采用Gauss点扩散函数能够分别计算出目标及导航星在CCD面的投影位置周围像素收集的能量：

上式中，(m,n)为所计算的像素点坐标，I(m,n)为所计算的像素点收集的光能，I₀为目标点光源的光能，(x₀,y₀)为目标点光源在CCD面的投影坐标；

根据公式(18)能够分别计算各星点散焦后的能量；

步骤五的具体实现方法为：

由于光轴扰动，在单次曝光时间T内相机光轴发生k次变化，每次扰动后的天球赤道坐标系和相机坐标系之间转换矩阵记为M_k，扰动后导航目标投影的位置坐标：

上式中，Nx,Ny是CCD平面尺寸，FOV为视场角大小，M_k1,M_k2,M_k3为M_k各行，X_i为导航恒星在天球赤道坐标系下的位置；

由于两次扰动之间的时间间隔较短，故将两次光轴方向变化间的运动视为是匀速的，在任意时刻t有：

M_t＝M_k·(I+ω·Δt) (20)

其中，M_t为任意时刻t的转换矩阵，ω为光轴的方向变换角速度，Δt为时间步长；

在任意时刻t，将M_t作为天球赤道坐标系和相机坐标系间的状态转移矩阵M，按照步骤三和步骤四计算任意时刻t星点的位置和能量；将不同时刻下的星点能量叠加即得到暗弱小天体探测器拖曳星图；

步骤六具体实现方法为，

步骤6.1：对长曝光下产生的拖曳星图进行去噪预处理；

考虑到图像噪声的非高斯性及低信噪比性，采用阈值小波变化的方法进行去噪；小波去噪方法的思路是：进行小波变换后，信号主要集中在低频，而噪声集中在高频上，故给定适当的阈值除去噪声；

采用软阈值小波方法进行去噪，软阈值是指若小波变换系数比设置的阈值小则将该系数置零，若小波变换系数比该阈值大则将该系数向阈值方向收敛，小波变换系数表达形式为：

其中，ω_i,j为小波变换系数，λ为设定的阈值；

步骤6.2：通过连通域判定轨迹所在的兴趣区域ROI，通过兴趣轨迹选取准则选取兴趣轨迹，采用逐像素细化算法One Pass Thinning Algorithm对兴趣轨迹进行细化并提取出轨迹骨架；

长曝光拖曳图像模型表示式为：

I_Blur＝K_blur·I_static+E (22)

其中，I_Blur是长曝光的拖曳星图，I_static是需要被复原的静态星图，E为噪声项，K_blur是每一像素的模糊转移矩阵；当星图中各像素发生的模糊是均匀时，则式(22)变换为卷积的形式，如式(23)所示：

I_Blur＝k*I_static+e (23)

其中，k被称为模糊核，与静态星图I_static的每一像素发生卷积运算使图像产生模糊，其物理含义是反映每一像素点的运动特性；

由于长曝光拖曳星图的拖曳是点状目标的运动模糊，故星点拖曳轨迹本身能够反映地曝光时间内的运动，而模糊核的本质就是像素点的运动特性，故采用星点拖曳轨迹来估计拖曳图像的模糊核；

通过连通域判定轨迹所在的兴趣区域ROI，通过兴趣轨迹选取准则选取兴趣轨迹；

所述兴趣轨迹选取准则选为亮度兴趣轨迹选取准则，即选取亮度最亮的轨迹为兴趣轨迹，兴趣轨迹的选取准则为：

max average(G_{ROI_i}) (24)

其中，G_{ROI_i}是第i个兴趣区的灰度矩阵；

OPTA算法是一种基于形态学的单循环并行模板匹配算法，其优点是运算简单，计算速度很快，且能够较好地保留图像的连通性和形状特征；采用OPTA算法对兴趣轨迹进行细化并提取出轨迹骨架，具体实现方法如下：按从左到右，从上到下的方式，从最左上角的像素点开始抽取出每个灰度值不为0的像素点的邻域，将该邻域与消除模板中的任意一个消除模板匹配，若不能匹配则保留像素点的值，否则要再将该邻域与保留模板匹配，若能匹配则仍然保留像素点的值，若不能匹配则删去该像素点，即将像素点的灰度值置零；对所有像素点进行消除模板和保留模板的匹配，保留下的灰度值不为0的像素点即为提取出的轨迹骨架；

步骤6.3：以步骤6.2中提取出的轨迹骨架为模糊核的估计初值，通过稀疏先验信息构建图像复原算法的代价函数，对长曝光拖曳星图进行盲复原，得到长曝光拖曳星图复原图像；

复原图像问题需估计出合适的模糊核及复原图像，使得图像复原的代价函数尽量最小，代价函数形式为：

然而，复原图像问题本质上是NP，即non-deterministic polynomial问题，即由于信息的匮乏，会导致结果不唯一甚至病态；为解决星图模糊复原的问题，需要根据导航星图的实际特点引入相关先验信息建立正则惩罚项；

深空导航星图的重要特征是星点能量集中在极少数像素点中，通过步骤一去噪声预处理后，绝大部分像素点为黑色背景，故利用图像的灰度稀疏作为先验信息；同理，由于有值像素点较少，故图像的边缘像素点也就更少，所以图像的灰度梯度也是稀疏的；稀疏先验信息通过L₀范数构建，即

最终采用的星图复原代价函数的形式为：

其中，项是为保证估计出的模糊核的简单性；

将步骤6.2中提取出的轨迹骨架为模糊核的估计初值，以代价函数最小为指标估计复原图像；由估计出的复原图像更新模糊核；如此循环迭代直至达到给定的循环的终止条件，此时估计出的复原图像即为最终得到的长曝光拖曳星图的复原图像；

步骤6.4：对步骤6.3得到的复原图像进行连通域判定，并提取复原图像中各星点的灰度光心，完成长曝光拖曳星图特征提取；

对步骤6.3得到的长曝光拖曳星图复原图进行连通域判定，去除连通域小于预设像素数量的目标；

步骤6.4中去除连通域小于预设像素数量的目标选去除连通域小于2个像素的目标；

然后，根据导航星覆盖像素的像素值确定复原星图导航星光心所在的位置，即实现提取复原图像中各星点的灰度光心，完成长曝光拖曳星图特征提取，并验证轨迹图像星点质心提取的图像处理算法的有效性，即进行图像处理算法的测试；

步骤6.4中计算灰度光心的方法采用灰度加权算法，灰度加权算法的公式为：