一种深空光通信扩展信标捕获跟踪方法

摘要

一种深空光通信扩展信标捕获跟踪方法，它涉及扩展信标图像的处理方法。它解决传统方法对扩展信标进行捕获和跟踪时精度低、误差大的问题。本发明把实际信标图像的噪声合理的近似为附加高斯白噪声，采用最小二乘法捕获信标图像，确定信标图像的中心位置，使系统进入跟踪模式，然后基于离散傅立叶变换和极大似然算法对信标图像的旋转角度和平移量进行计算，根据平移量和旋转角度的计算值，可以计算出信标图像的中心，根据该中心位置更新光通信终端天线的指向。本发明通过缩小探测器视场、提高测量精度、采用天线扫描结合像素扫描的方法进行扩展信标捕获，节省捕获时间，信标的成功捕获概率为98％以上；同时，本发明的控制光学天线光学天线的跟踪误差控制在5％以内。

说明书

技术领域

本发明涉及扩展信标图像的处理方法。

背景技术

用激光进行深空探测科学数据回传是深空通信的发展方向，用激光进行深空远距离信息传输面临的严峻挑战是精确的对准要求，现有基于卫星光通信系统设计的深空光通信系统以地面发射的上行激光束作为信标，但是以上行激光束作为信标，不但需要采用大功率的激光器，而且易受到大气的影响和天气的制约，不能满足深空光通信链路全天候运行的要求，因此，我们考虑采用自然天体图像(可视图像或红外图像)作为信标。由于链路在运行过程中，飞行器上光学终端总是对准地球方向，因此最常选用的自然天体为地球及其附近的天体。由于这些天体在信标探测器上成像总是扩展到多个像素，因此，我们称之为扩展信标。采用自然天体图像作为信标，无须地面发射上行光束作为对准基准，因此可使飞行器上光学终端保持相对上行链路的独立性，使飞行器上光学终端免受大气的影响和天气的制约，同时采用扩展天体图像作为信标，其光强相对稳定，且可以提供较高的跟踪速率。由于采用扩展天体图像作为信标，传统的基于上行激光束作为信标的捕获跟踪方法将不再适应本方案。

采用地球图像作为信标，通常根据地球星历表和飞行器星历表确定初始对准方向，由于飞行器星历表和地球星历表误差，飞行器上微动力学环境及热力学环境的的影响，初始对准方向与实际对准方向总是存在一定的误差，因此需要在初始对准方向确定以后，对地球图像进行捕获和跟踪，而捕获和跟踪地球图像的关键是精确的确定地球图像的中心位置，根据地球图像中心位置和地面站的关系，确定飞行器上光通信终端接收天线的瞄准方向。在这个过程中，对地球图像的捕获和跟踪算法的研究是其中的核心问题。传统方法总是直接根据测得实际图像进行计算获得其中心位置，这样做误差很大、精度低。

发明内容

为了解决传统方法对扩展信标进行捕获和跟踪时精度低、误差大的问题，本发明提供了一种深空光通信扩展信标捕获跟踪方法。深空光通信扩展信标的捕获用于深空光通信链路的建立，对扩展信标的捕获，采用天线扫描结合像素扫描的方式进行信标的搜索，扩展信标的跟踪用于通信过程中保持深空光通信链路不中断。

本发明的深空光通信扩展信标捕获跟踪方法按以下步骤进行：

一、按以下步骤捕获扩展信标位置：①在深空激光通信链路建立和运行的过程中，根据自然天体和飞行器星历表确定惯性坐标系下的光通信终端初始对准方向，该方向指向所要捕获的扩展信标；②将上述惯性坐标系中的瞄准方向矢量通过坐标变换转换到星上俯仰坐标系中，使所述通信光束的控制在飞行器的星上俯仰坐标系中进行；③根据航天器的相关技术参数，结合统计学原理确定不确定区域平面视场角；④在上述不确定区域内按照预定角间隔进行天线扫描，在扫描过程中，使天线在扫描轨迹上各位置点停留一段时间，信标探测器对所对准的区域成像，然后将获得的实际扩展信标图像和预存在系统内的标准扩展信标图像的二维矩阵分别先进行傅立叶-梅林变换并取二者的幅度矩阵，再利用二次插值将二者的幅度矩阵转换到对数极坐标系，并按下述公式(1)计算上述二者幅度矩阵的误差系数，

$>>L>=>>\Sigma >>i>=>1>>M>>>\Sigma >>j>=>1>>N>{}^{}$

上式中，L为上述二者幅度矩阵的误差系数，M为插值到对数极坐标系中矩阵的行数，N为插值到对数极坐标系中矩阵的列数，M_S(i，j)为上述标准扩展信标图像在对数极坐标系中的傅立叶-梅林变换的幅度，M_F(i，j)为实际信标图像在对数极坐标系中的傅立叶-梅林变换的幅度；⑤控制探测器的光通信终端天线对准上述误差系数取最小值时的位置，于是完成捕获过程；

二、在上述捕获过程结束后进入扩展信标的跟踪过程，跟踪的具体步骤如下：I、第一次跟踪时，将上述预存在系统内的标准扩展信标图像或上述捕获完成后跟踪过程初始时刻的信标图像作为这次跟踪的参考图像，然后利用探测器获得此时发生了相对运动的扩展信标的实测图像；II、按以下步骤计算在直角坐标系中实测图像相对于参考图像的旋转角：A、将实测图像看成是参考图像和附加高斯白噪声之和，然后对跟踪探测到的实际图像和参考图像分别进行离散傅立叶变换并取模，获得上述两图像的幅度谱；B、将实际图像和参考图像的幅度谱利用二次插值转换到极坐标系下获得上述两图像在极坐标系下的幅度矩阵，并按下述公式(2)计算附加高斯白噪声在极坐标系下的幅度矩阵为，

G(m，n)＝M_R(m，n)-M₁(m，n-k)           (2)

上式中，G(m，n)为附加高斯白噪声在极坐标系中的幅度二维矩阵，M_R(m，n)为极坐标系中实测图像的幅度矩阵，M₁(m，n-k)为极坐标系中参考图像的幅度矩阵，k为极坐标系中实测图像相对于参考图像角度平移量(极坐标系中角度的平移表现在直角坐标中为角度旋转)；C、定义上述附加高斯白噪声在极坐标系下的幅度矩阵的各像素点服从均值为0的高斯分布且各点噪声值相互独立，然后按下述公式(3)计算所述角度平移量k的极大似然估计，

$>ver>>k>^>>=>->>N>>2>\pi >>>·>>>>\Sigma >>m>=>0>>M>>>\Sigma >>n>=>0>>N>>n>|>\omega >>(>m>,>n>)>>|>\xi >>>>\Sigma >>m>=>0>>M>>>\Sigma >>n>=>0>>N>{}^{}$

上式中，为所述角度平移量k的极大似然估计，|ω(m，n)|为上述M_R(m，n)的离散傅立叶变换与上述M₁(m，n-k)的共轭离散傅立叶变换相乘的模，ξ为ω(m，n)的傅立叶相位角；D、按下述公式(4)计算直角坐标中实际图像相对于参考图像旋转角，

$>>>\theta >0>>=>>>2>\pi >>N>>·ver>>k>^>>->->->>(>4>)>>>s>$

上式中，θ₀为直角坐标系中实测图像相对于参考图像的旋转角，N为插值后的图像矩阵列数，它的值与公式(1)中的N一致；III、按以下步骤计算直角坐标系中实测图像相对于参考图像的平移量：a、将上述旋转角θ₀代入下述公式(5)计算旋转后的参考图像的二维矩阵，

s₂(m，n)＝s(m·cosθ₀+n·sinθ₀-x₀，-m·sinθ₀+n·cosθ₀-y₀)   (5)   上式中，s₂(m，n)为转动后的参考图像二维矩阵，s(m，n)为初始参考图像二维矩阵，x₀为直角坐标系中实测图像相对于参考图像的横坐标平移量，y₀为直角坐标系中实测图像相对于参考图像的纵坐标平移量；b、根据下述公式(6)计算直角坐标系中实测图像相对于参考图像的横坐标平移量和纵坐标平移量，

$>\left(\begin{array}{}>>->>1>>2>\pi >>>>\Sigma >>m>=>0>>M>>>\Sigma >>n>=>0>>N>>m>|>\upsilon >>(>m>,>n>)>>|>\beta >>(>m>,>n>)>>=>>x>0>>>\Sigma >>m>=>0>>M>>>\Sigma >>n>=>0>>N>>m>>m>M>>|>\upsilon >>(>m>,>n>)>>|>+>>y>0>>>\Sigma >>m>=>0>>M>>>\Sigma >>n>=>0>>N>>m>>n>N>>|>\upsilon >>(>m>,>n>)>>|>\; >>->>1>>2>\pi >>>>\Sigma >>m>=>0>>M>>>\Sigma >>n>=>0>>N>>n>|>\upsilon >>(>m>,>n>)>>|>\beta >>(>m>,>n>)>>=>>x>0>>>\Sigma >>m>=>0>>M>>>\Sigma >>n>=>0>>N>>n>>m>M>>|>\upsilon >>(>m>,>n>)>>|>+>>y>0>>>\Sigma >>m>=>0>>M>>>\Sigma >>n>=>0>>N>>n>>n>N>>|>\upsilon >>(>m>,>n>)>>|>\; >\end{array}\right)>->->->>(>6>)>>>s>$

上式中，υ(m，n)为旋转后的参考图像和实测图像的傅立叶变换的共轭矩阵的积，|υ(m，n)|为υ(m，n)的模，β(m，n)为υ(m，n)的傅立叶相位角矩阵；IV、根据旋转角度和平移量，按下述公式(7)确定实测图像的中心位置，

                 x_t＝x′cosθ₀+y′sinθ₀-x₀          (7)

                 y_t＝-x′sinθ₀+y′cosθ₀-y₀

上式中，(x′，y′)为这次跟踪所用参考图像的中心位置，(x_t，y_t)为实测图像的中心位置；V、将探测器的光通信终端天线对准此时实测图像的中心位置，于是就完成了第一次跟踪过程；VI、当进行下一次跟踪时，将第一次跟踪所测的扩展信标实测图像作为这次跟踪的参考图像，然后利用探测器获得此时发生了相对运动的扩展信标的实测图像，再重复执行第II至IV步，并将探测器的光通信终端天线对准这次跟踪实测图像的中心位置；VII、在以后的每次跟踪过程中都以上一次跟踪获得的实测图像作为本次跟踪的参考图像，并根据获得的本次跟踪实测图像的中心位置调整探测器光通信终端天线的对准方向。

工作原理：本发明根据飞行器技术参数，结合统计学原理确定了扫描的不确定区域平面视场角，并利用最小二乘法计算参考图像和实际测量图像的误差系数L，在天线扫描完成以后，控制天线对准L最小的位置点，这时，信标图像进入探测器视场，根据系统内预存的参考图像来确定实际信标图象的中心，即可认为捕获到该信标图像。在对数极坐标系中，图像的伸缩和旋转都可以看成是在对数极坐标系中两个轴上的平移，所以由傅立叶变换的性质，坐标平移只影响其傅立叶变换相位的变化，而不影响其幅度的变化，因此，其幅度应具有不变性，所以对于参考图像和实测图像而言，在对数极坐标系中傅立叶变换的幅度应具有最大相关性。在信标的捕获完成以后，系统进入跟踪模式，在跟踪阶段，系统的主要工作是根据信标的运动来调整天线的指向，因此跟踪过程的主要任务是对信标的运动进行补偿。本发明把信标探测器探测到的图像看成参考图像和高斯白噪声之和的形式，而探测的实际图像相对于参考图像有一定的平移量和旋转角，于是在跟踪时为了确定实际图像的中心位置必须分别计算它们之间的平移量和旋转角，这时本发明引入极坐标系中的角度平移量k，并利用统计学的方法和傅立叶变换计算上述平移量和旋转角。

发明效果：1)本发明的方法通过缩小探测器视场、提高测量精度、采用天线扫描结合像素扫描的方法进行扩展信标捕获，节省捕获时间，信标的成功捕获概率为98％以上；2)本发明能控制光学天线的跟踪误差在5％以内。

附图说明

图1是具体实施方式捕获过程的参考信标图像，图2是具体实施方式捕获过程的实测信标图像，图3是具体实施方式跟踪过程的参考图像，图4是具体实施方式跟踪过程的实测信标图像。

具体实施方式

本具体实施方式以自然天体地球为所要探测的扩展信标，并按以下步骤对其进行捕获跟踪：一、按以下步骤捕获地球的位置：①在深空激光通信链路建立和运行的过程中，根据自然天体和飞行器星历表确定惯性坐标系下的光通信终端初始对准方向，该方向指向所要捕获的扩展信标；②将上述惯性坐标系中的瞄准方向矢量通过坐标变换转换到星上俯仰坐标系中，使所述通信光束的控制在飞行器的星上俯仰坐标系中进行；③根据航天器的相关技术参数，结合统计学原理确定不确定区域平面视场角；④在上述不确定区域内按照预定角间隔进行天线扫描，在扫描过程中，使天线在扫描轨迹上各位置点停留一段时间，信标探测器对所对准的区域成像，然后将获得的实际地球图像和预存在系统内的标准地球图像的二维矩阵分别先进行傅立叶-梅林变换并取二者的幅度矩阵，再利用二次插值将二者的幅度矩阵转换到对数极坐标系，并按下述公式(1)计算上述二者幅度矩阵的误差系数，

$>>L>=>>\Sigma >>i>=>1>>M>>>\Sigma >>j>=>1>>N>{}^{}$

上式中，L为上述二者幅度矩阵的误差系数，M为插值到对数极坐标系中矩阵的行数，N为插值到对数极坐标系中矩阵的列数，M_S(i，j)为上述标准地球图像在对数极坐标系中的傅立叶-梅林变换的幅度，M_F(i，j)为实际地球图像在对数极坐标系中的傅立叶-梅林变换的幅度；⑤控制探测器的光通信终端天线对准上述误差系数取最小值时的位置，于是完成捕获过程，然后缩小信标探测器的视场、提高测量精度、获取此时地球的实际图像，并按下述公式(8)和(9)计算所探测地球图像的中心位置，

$>>x>=>>>>\Sigma >>m>=>0>>>M>->1>>>>\Sigma >>n>=>0>>>N>->1>>>m>·>S>>(>m>,>n>)>>>>>\Sigma >>m>=>0>>>M>->1>>>>\Sigma >>n>=>0>>>N>->1>>>S>>(>m>,>n>)>>>>->->->>(>8>)>>>s>$

$>>y>=>>>>\Sigma >>m>=>0>>>M>->1>>>>\Sigma >>n>=>0>>>N>->1>>>n>·>S>>(>m>,>n>)>>>>>\Sigma >>m>=>0>>>M>->1>>>>\Sigma >>n>=>0>>>N>->1>>>S>>(>m>,>n>)>>>>->->->>(>9>)>>>s>$

上式中，S(m，n)为上述标准地球图像，x为探测到的实际地球图像中心点的横坐标，y为探测到的实际地球图像中心点的纵坐标；

二、在上述捕获过程结束后进入扩展信标的跟踪过程，跟踪的具体步骤如下：I、第一次跟踪时，将上述捕获完成后跟踪过程初始时刻的地球图像作为这次跟踪的参考图像，然后利用探测器获得此时发生了相对运动的地球实测图像；II、按以下步骤计算直角坐标系中实测图像相对于参考图像的旋转角：A、将实测图像看成是参考图像和附加高斯白噪声之和，然后对跟踪探测到的实际图像和参考图像分别进行离散傅立叶变换并取模，获得上述两图像的幅度谱；B、将实际图像和参考图像的幅度谱利用二次插值转换到极坐标系下获得上述两图像在极坐标系下的幅度矩阵，并按下述公式(2)计算附加高斯白噪声在极坐标系下的幅度矩阵为，

G(m，n)＝M_R(m，n)-M₁(m，n-k)           (2)

上式中，G(m，n)为附加高斯白噪声在极坐标系中的幅度二维矩阵，M_R(m，n)为极坐标系中实测图像的幅度矩阵，M₁(m，n-k)为极坐标系中参考图像的幅度矩阵，k为极坐标系中实测图像相对于参考图像角度平移量(极坐标系中角度的平移表现在直角坐标中为角度旋转)；C、定义上述附加高斯白噪声在极坐标系下的幅度矩阵的各像素点服从均值为0的高斯分布且各点噪声值相互独立，然后按下述公式(3)计算所述角度平移量k的极大似然估计，

$>ver>>k>^>>=>->>N>>2>\pi >>>·>>>>\Sigma >>m>=>0>>M>>>\Sigma >>n>=>0>>N>>n>|>\omega >>(>m>,>n>)>>|>\xi >>>>\Sigma >>m>=>0>>M>>>\Sigma >>n>=>0>>N>{}^{}$

上式中，为所述角度平移量k的极大似然估计，|ω(m，n)|为上述M_R(m，n)的离散傅立叶变换与上述M₁(m，n-k)的共轭离散傅立叶变换相乘的模，ξ为ω(m，n)的傅立叶相位角；D、按下述公式(4)计算直角坐标中实际图像相对于参考图像旋转角，

$>>>\theta >0>>=>>>2>\pi >>N>>·ver>>k>^>>->->->>(>4>)>>>s>$

上式中，θ₀为直角坐标系中实测图像相对于参考图像的旋转角，N为插值后的图像矩阵列数，它的值与公式(1)中的N一致；III、按以下步骤计算直角坐标系中实测图像相对于参考图像的平移量：a、将上述旋转角θ₀代入下述公式(5)计算旋转后的参考图像的二维矩阵，

s₂(m，n)＝s(m·cosθ₀+n·sinθ₀-x₀，-m·sinθ₀+n·cosθ₀-y₀)    (5)

上式中，s₂(m，n)为转动后的参考图像二维矩阵，s(m，n)为初始参考图像二维矩阵，x₀为直角坐标系中实测图像相对于参考图像的横坐标平移量，y₀为直角坐标系中实测图像相对于参考图像的纵坐标平移量；b、根据下述公式(6)计算直角坐标系中实测图像相对于参考图像的横坐标平移量和纵坐标平移量，

$>\left(\begin{array}{}>>->>1>>2>\pi >>>>\Sigma >>m>=>0>>M>>>\Sigma >>n>=>0>>N>>m>|>\upsilon >>(>m>,>n>)>>|>\beta >>(>m>,>n>)>>=>>x>0>>>\Sigma >>m>=>0>>M>>>\Sigma >>n>=>0>>N>>m>>m>M>>|>\upsilon >>(>m>,>n>)>>|>+>>y>0>>>\Sigma >>m>=>0>>M>>>\Sigma >>n>=>0>>N>>m>>n>N>>|>\upsilon >>(>m>,>n>)>>|>\; >>->>1>>2>\pi >>>>\Sigma >>m>=>0>>M>>>\Sigma >>n>=>0>>N>>n>|>\upsilon >>(>m>,>n>)>>|>\beta >>(>m>,>n>)>>=>>x>0>>>\Sigma >>m>=>0>>M>>>\Sigma >>n>=>0>>N>>n>>m>M>>|>\upsilon >>(>m>,>n>)>>|>+>>y>0>>>\Sigma >>m>=>0>>M>>>\Sigma >>n>=>0>>N>>n>>n>N>>|>\upsilon >>(>m>,>n>)>>|>\; >\end{array}\right)>->->->>(>6>)>>>s>$

上式中，υ(m，n)为旋转后的参考图像和实测图像的傅立叶变换的共轭矩阵的积，|υ(m，n)|为υ(m，n)的模，β(m，n)为υ(m，n)的傅立叶相位角矩阵；IV、根据旋转角度和平移量，按下述公式(7)确定实测图像的中心位置，

       x_t＝x′cosθ₀+y′sinθ₀-x₀

       y_t＝x′sinθ₀+y′cosθ₀-y₀                          (7)

上式中，(x′，y′)为这次跟踪所用参考图像的中心位置，即按上述公式(8)和(9)所计算的实际地球图像的中心点位置，(x_t，y_t)为实测图像的中心位置；V、将探测器的光通信终端天线对准此时实测图像的中心位置，于是就完成了第一次跟踪过程；VI、当进行下一次跟踪时，将第一次跟踪所测的扩展信标实测图像作为这次跟踪的参考图像，然后利用探测器获得此时发生了相对运动的扩展信标的实测图像，再重复执行第II至IV步，并将探测器的光通信终端天线对准这次跟踪实测图像的中心位置；VII、在以后的每次跟踪过程中都以上一次跟踪获得的实测图像作为本次跟踪的参考图像，并根据获得的本次跟踪实测图像的中心位置调整探测器光通信终端天线的对准方向。

上述第③步根据飞行器技术参数，结合统计学原理确定不确定区域平面视场角的方法是：所述通信光束在星上俯仰坐标系中的方向可用方位角和俯仰角来表示，先定义上述方位角和俯仰角的误差都服从均值为0、方差为σ的高斯分布，方差σ可根据飞行器的技术参数确定；再按下述公式(10)计算不确定区域内成功捕获的概率，

$>>P>=sup>>\int >>->>>\theta >u>>2>>>>>\theta >u>>2>sup>sup>>\int >>->>>\theta >u>>2>>>>>\theta >u>>2>sup>>f>>(>>ϵ>h>>,>>ϵ>v>>)>>d>>ϵ>h>>>dϵ>v>>->->->>(>10>)>>>s>$

上式中，P为不确定区域内链路成功建立的概率，ε_h为方位角的误差，ε_v为俯仰角的误差，f(ε_h，ε_v)为上述方位角和俯仰角误差的联合概率密度，θ_u为视场角；然后计算上述公式(10)中不确定区域内链路成功建立的概率P≥98％时的视场角θ_u，即为不确定区域平面视场角。

在图像信标捕获过程中，视场扫描要求覆盖整个不确定区域，因此，相邻两个位置的角间隔d必须满足

$>>d>\le >>\theta >fov>>/>>2>>->->->>(>11>)>>>s>$

上式中，θ_fov是信标探测器的视场。

在捕获地球的过程中，如图1所示，预存的参考信标图像为一幅257×250的地球图像轮廓；如图2所示，探测到的一幅512×512的实际地球图像为上述参考图像与一附加高斯白噪声之和，图像各像素上的噪声独立同分布，其均值为0，方差为σ，信噪比为1。下表1为根据不同的方法计算的地球图像中心。

                                            表1

  实际信标图像中心  (pixel)  捕获算法计算的图像中心  (pixel)  直接根据信标图像计算的图像中心  (pixel)  x  y  x  y  x  y  251.05  254.27  251.05  254.27  250.78  254.75

根据上表可以看出，本发明采用最小二乘算法对信标图像进行捕获，可以精确的确定信标图像的中心，其误差非常微小，几乎为0，而直接采用实际信标图像来计算其中心，因附加高斯白噪声的影响，其误差较大，在0.5个像素左右。当精确的确定信标图像的中心后，可根据地面站和地球图像中心位置的相对关系来精确确定光通信终端的天线指向。

在跟踪阶段，需要对地球图像的旋转和平移运动进行补偿。参考图像如图3所示，运动后的实测地球图像如图4所示。根据本发明结合了离散傅立叶变换和极大似然算法的跟踪方法求出实际信标图像的平移量和旋转角如表2所示。

                                            表2

  实际平移量(pixel)  平移量计算值(pixel)  实际旋转角度(°)  旋转角计算值(°)  x  y  x₀  y₀  θ  θ₀  -3.00  -4.00  -2.92  -3.87  5.00  5.04

根据上表可以看出，计算得到的信标图像的平移值和旋转角度与实际值相差非常小，在X方向，误差仅为2.7％，在Y方向，其误差为3.3％，旋转角度误差仅为0.8％，完全可以满足深空激光通信光束精确对准的要求。根据计算得到的平移量和旋转角，由公式(7)可得运动后信标图像的中心如表3所示，

                                      表3

  参考信标图像中心(pixel)  计算的所测运动后的实际信标图像的中心(pixel)  x  y  x_t  y_t  124.05  123.27  137.32  115.77

由上表可知，根据平移量和旋转角度的计算值，可精确的确定信标图像的中心，从而动态的更新光通信终端天线的指向。