CN103791899A - 一种基于成像传感器视差的卫星姿态颤振探测方法 - Google Patents

Abstract

本发明涉及一种基于成像传感器视差的卫星姿态颤振探测方法，包括以下步骤：1)获取视差影像，并对其进行预处理，所述的视差影像包括相互间具有观测视差的参考影像和待匹配影像；2)对视差影像进行精确密集匹配处理，获得两幅影像的视差偏差值；3)根据视差偏差值计算卫星姿态颤振值。与现有技术相比，本发明可在不依赖地面控制点和附加硬件设备的基础上，对星体的姿态颤振进行探测和分析，得到真实的星体在垂轨和沿轨方向的姿态颤振值。

Description

一种基于成像传感器视差的卫星姿态颤振探测方法

技术领域

本发明涉及卫星遥感成像技术，尤其是涉及一种基于成像传感器视差的卫星姿态颤振探测方法。

背景技术

高分辨率遥感成像卫星自身的姿态稳定性和姿态控制性能，是决定高分辨率卫星成像质量和几何定位精度的重要因素之一。卫星星体的姿态数据，作为传感器严格成像模型中的外方位参数，在几何模型的构建、定位精度的提高、有理函数模型的生成等方面起着决定性的作用。卫星星体的姿态颤振作为卫星运行过程中的普遍现象，对卫星姿态数据的精度有着直接的影响。在高性能定轨和定姿设备的支持下，国外的星载姿态传感器能够对星体的姿态颤振进行精确的探测，从而获得高精度的姿态数据，其高分辨率卫星遥感影像的直接定位精度可达到10米级。但由于硬件技术等原因，国产姿态传感器存在着采样频率低、定姿精度较差和可靠性不足等问题，无法对卫星星体的颤振进行探测，从而对几何定位精度和成像质量产生了较大影响。同时，星载精密定姿传感器也是西方发达国家对中国禁运的主要设备之一。如何在现有硬件条件下，对国产高分辨率遥感卫星姿态颤振进行探测，提高姿态数据精度，成为我国当前遥感测绘领域一个迫切的任务。

发明内容

本发明的目的就是为了克服上述现有技术存在的缺陷而提供一种基于成像传感器视差的卫星姿态颤振探测方法，可以在无地面控制点和附加硬件设备的情况下，对遥感卫星的姿态颤振进行探测。

本发明的目的可以通过以下技术方案来实现：

一种基于成像传感器视差的卫星姿态颤振探测方法，包括以下步骤：

1)获取视差影像，并对其进行预处理，所述的视差影像包括相互间具有观测视差的参考影像和待匹配影像；

2)对视差影像进行精确密集匹配处理，获得两幅影像的视差偏差值；

3)根据视差偏差值计算卫星姿态颤振值。

步骤1)中，所述的预处理包括影像增强和视差影像预配准。

步骤2)中，所述的精确密集匹配具体为：对两幅影像中相同地物对应的各影像行，以等间隔均匀取点的方式进行逐点匹配，首先采用基于灰度互相关的方法进行粗匹配，获取整像素匹配点，再以粗匹配结果作为初值，使用基于二维曲面峰值拟合的相位相关的方法进行精确匹配，获取亚像素匹配点，对整行影像所有匹配点的结果取平均值，作为整行影像的匹配结果，得到各行影像的视差偏差值。

步骤3)中，所述的计算卫星姿态颤振值具体为：

301)使用频谱分析的方法由步骤2)获得的视差偏差值得到卫星姿态颤振的频率和幅值，当颤振只含有一种频率分量时，执行步骤302)，当颤振含有多种频率分量时，执行步骤303)；

302)通过频谱分析获得视差偏差值的幅值

和频率ω，并由相位谱获取其初相位

将t时刻的视差偏差值与影像几何偏移值表示为：

f(t)=A_fcos(ωt+f₀)

f(t+Δt)＝A_fcos(ωt+f₀+Δf)

式中，A_f、f₀分别为影像几何偏移值的幅值和初相位，Δt为视差影像的成像时间间隔，Δf=ωΔt；

根据

，计算得：

从而获得影像几何偏移值f(t)；

303)通过频谱分析获得视差偏差值各频率分量的幅值

频率和初相位

影像几何偏移值的各频率等于视差偏差值的频率；

对式

进行离散傅里叶变换，获得：

式中，为虚数单位，

k为离散傅里叶变换后频域下对应t时刻的虚数，N为参与离散傅里叶变换的数据总数；

则

对f(t)进行频谱分析，获取影像几何偏移值各个频率分量的幅值A_fi和频率ω_i，由相位谱获取各个频率分量的初相位f_0i；

从而获得影像几何偏移值f(t)：

M为频率分量个数；

304)根据影像几何位移值f(t)计算星体在垂轨方向和沿轨方向的姿态颤振值

(t)和

(t)：

${\∂g}_{\mathrm{cross}}\left(t\right)=\frac{f_x\left(t\right)}{f_c}$

${\∂g}_{\mathrm{along}}\left(t\right)=\frac{f_y\left(t\right)}{f_c}$

式中，f_x(t)和f_y(t)分别为f(t)在像平面行方向和列方向的分量，f_c为相机主距。

与现有技术相比，本发明在不依赖地面控制点和附加硬件设备的基础上，根据高分卫星多光谱传感器各波段线阵成像视差的特点，对星体的姿态颤振进行了探测。实验表明该方法在国产高分辨率测绘卫星上取得了良好的效果。使用本发明方法在国产高分辨率测绘卫星上得到了应用，使用其多光谱影像，探测到了频率为0.63Hz，垂轨方向最大为2.6角秒，沿轨方向最大为1.5角秒的姿态颤振，其对多光谱影像的影响为，垂轨方向最大1.0像素，沿轨方向最大0.6像素。

附图说明

图1为颤振探测原理示意图；

图2为本发明的流程示意图；

图3为视差影像的匹配流程示意图；

图4为各影像行的等间隔密集匹配示意图；

图5为某国产测绘卫星多光谱传感器结构示意图；

图6为国产测绘卫星传感器每片CCD阵列的视差构成；

图7为国产测绘卫星河南某地区域影像；

图8为各片CCD影像的沿轨方向视差图；

其中，(8a)为CCD1的沿轨方向视差图，(8b)为CCD2的沿轨方向视差图，(8c)为CCD3的沿轨方向视差图；

图9为各片CCD影像的垂轨方向视差图；

其中，(9a)为CCD1，(9b)为CCD2，(9c)为CCD3；

图10为各片CCD影像的视差偏差值；

其中，(10a)、(10b)分别为CCD1的沿轨、垂轨方向视差偏差值，(10c)、(10d)分别为CCD2的沿轨、垂轨方向视差偏差值，(10e)、(10f)分别为CCD3的沿轨、垂轨方向视差偏差值；

图11为影像的几何偏移值；

其中，(11a)为沿轨方向，(11b)为垂轨方向；

图12为影像的几何偏移值的拟合结果；

其中，(12a)为沿轨方向，(12b)为垂轨方向；

图13为卫星的姿态颤振值；

其中，(13a)为垂轨方向颤振值，(13b)为沿轨方向颤振值。

具体实施方式

下面结合附图和具体实施例对本发明进行详细说明。本实施例以本发明技术方案为前提进行实施，给出了详细的实施方式和具体的操作过程，但本发明的保护范围不限于下述的实施例。

1、基于传感器视差的颤振探测原理

具有视差观测的线阵CCD推扫式传感器系统是指传感器的多条CCD线阵在不同的观测时刻连续对同一地物进行成像，即多条CCD线阵在对同一地物成像时，在沿轨方向的成像时刻存在着时间差，如图1所示。对于多条CCD线阵在推扫过程中形成的视差影像，在影像的重叠区域，理论上在像平面相同地物应投影在相同的位置上，即相同的地物在配准后的影像中应具有相同的像素坐标。但是在推扫过程中，在多条CCD线阵依次成像的时间间隔里，星体的姿态颤振会引起传感器外方位元素的变化，从而造成相同地物在多条CCD线阵形成的视差影像中，其几何位置会有微小的偏移，如图1所示。

因此，针对具有视差观测的线阵CCD推扫式传感器的这一特点，可以对视差影像对进行逐行的密集匹配，根据匹配的结果，计算各成像时刻对应行的在沿轨和垂轨方向的视差值。理论上，当卫星的姿态稳定不存在颤振时，各行影像的视差值应为一常量。当卫星的姿态受到颤振影响时，其各影像行存在着颤振引起的几何偏移，各行影像的视差值不再为常量，存在视差偏差。对此成像时刻t对应的各行影像的视差偏差值与姿态颤振在像平面引起的几何偏移f(t)的关系可以用下式表示：

式中，Δt为视差影像对应的传感器线阵的成像时间间隔，对于固定的传感器线阵，其值为常数。

通过对

的量测，可以对具有视差传感器系统的线阵CCD推扫式遥感成像卫星的姿态颤振进行探测，借助频谱分析和傅里叶变换等方法，可以由

得到颤振的真实频率和幅值。

2、颤振探测流程

如图2所示，基于成像传感器视差的姿态颤振的探测方法的主要包括三个步骤：

(1)视差影像预处理

视差影像是指由具有视差的观测系统获得的相互间具有观测视差的影像对。视差影像预处理包括影像的增强和视差影像预配准。影像增强目的是为了增强影像的特征，并消除影响对之间的线性灰度差异，提高后续匹配的可靠性，这里使用的影像增强方法为基于Wallis滤波的影像增强方法。具有视差的两幅影像按照成像时间的先后，将成像时间早的影像设为参考影像，将成像时间晚的影像设为待匹配影像，预配准目的是为了消除了两幅影像间大的坐标差异，减少影像密集匹配时误匹配点的可能性。

(2)影像精确密集匹配

视差影像对间的精确密集匹配，是获取影像视差偏差值的关键。匹配结果的精度将直接影响颤振探测的精度。精确密集匹配具体为：对两幅影像中相同地物对应的各影像行，以等间隔均匀取点的方式进行逐点匹配，首先采用基于灰度互相关的方法进行粗匹配，获取整像素匹配点，再以粗匹配结果作为初值，使用基于二维曲面峰值拟合的相位相关的方法进行精确匹配，获取亚像素匹配点，对整行影像所有匹配点的结果取平均值，作为整行影像的匹配结果，得到各行影像的视差偏差值。

基于灰度互相关的匹配方法，其基本原理为比较两幅模板影像的相似性。其使用归一化互相关系数作为相似性的度量值，表达式如下：

$\mathrm{Cov}(u,v)=\frac{{\mathrm{\Σ}}_{x,y}[M(x,y)-{\stackrel{\‾}{M}}_{u,v}][S(x,y)-{\stackrel{\‾}{S}}_{u,v}]}{{\left\{{\mathrm{\Σ}}_{x,y}{[M(x,y)-{\stackrel{\‾}{M}}_{x,y}]}^2{\mathrm{\Σ}}_{x,y}{[S(x,y)-{\stackrel{\‾}{S}}_{x,y}]}^2\right\}}^{0.5}}---\left(2\right)$

式中M(x，y)为参考影像模板矩阵，S(x，y)为待匹配影像模板矩阵，Cov(u，v)为两幅模板影像在坐标(u，v)处的归一化互相关系数。当两模板影像对应的归一化互相关系数达到最大时，其两模板对应的中心点即为整像素匹配点。

基于二维曲面峰值拟合的相位相关的匹配方法，其基本原理为：对两幅模板影像在频域下进行卷积，获得其互功率谱。对互功率谱进行傅里叶逆变换，再使用二维的Sinc函数对互功率谱逆变换后的峰值进行拟合，由拟合的结果获取其亚像素的匹配值。式(3)给出了峰值拟合的表达式。

$r(n_1,n_2)\≅\mathrm{\α}\·\frac{\sin \left\{\mathrm{\π}(n_1+{\mathrm{\δ}}_1)\right\}}{\mathrm{\π}(n_1+{\mathrm{\δ}}_1)}\·\frac{\sin \left\{\mathrm{\π}(n_2+{\mathrm{\δ}}_2)\right\}}{\mathrm{\π}(n_2+{\mathrm{\δ}}_2)}---\left(3\right)$

式中r(n₁，n₂)为互功率谱逆变换的峰值，等式右边为二维Sinc函数的表达式。当r(n₁，n₂)达到最大值时，其在拟合曲面中的对应点，即为亚像素匹配点。

对相同地物对应的各影像行，以等间隔均匀取点的方式进行逐点匹配，计算待匹配影像中的各点在参考影像的像平面坐标下的位置差值。对整行影像所有匹配点的结果取平均值，作为整行影像的匹配结果，得到各行影像的视差偏差值。视差影像匹配的具体流程与各影像行的等间隔密集匹配方式如图3和图4所示。

(3)由视差偏差值计算姿态颤振值

使用频谱分析的方法由匹配获得的视差偏差值得到星体姿态颤振的频率和幅值。再结合传感器视差间隔Δt，由视差偏差值获取真实的颤振在各行影像引起的几何偏移值。最后，计算姿态颤振值在横滚、俯仰和航向方向的分量。

由影像匹配得到的视差偏差值，其本质上是间隔为Δt的像平面各行影像几何偏移值f(t)的差分。因此，在求取真实姿态颤振值之前，需要先由视差偏差值计算像平面各行影像的几何偏移值。

根据傅里叶变换的性质，我们可以将颤振引起的像平面各行影像的几何偏移值分解为多个余弦表达式相叠合的形式：

$f\left(t\right)=\underset{i=1}{\overset{M}{\mathrm{\Σ}}}{A}_{\mathrm{fi}}\cos ({\mathrm{\ω}}_{i}t+f_{0i})---\left(4\right)$

式(4)中，A_fi为第i分量的幅值，ω_i为第i分量的频率，f_0i为第i分量的初相位，M为频率分量个数。

同理由式(1)和同频合振动的性质，可知视差偏差值也可以表达为类似的形式：

式(5)中，A_i为第i分量的幅值，ω_i为第i分量的频率，为第i分量的初相位。

在此基础上，本发明采用两种不同思路的像平面颤振偏移值的计算方法。当颤振本身只含有一种频率分量时，可以在直接空域反算颤振偏移值的频率、幅值和相位。当颤振本身含有多种频率分量时，则需要在频域下反算几何偏移值的参数。

单一频率分量时影像几何偏移值的计算如下：

单频率情况的处理方法的基本思路是利用同频余弦振动的叠合与分解的特性。通过频谱分析获得视差偏差值的幅值

和频率ω。再由相位谱获取其初相位

根据同频合振动的性质，影像几何偏移值的频率等于视差偏差值的频率，即ω。t时刻的视差偏差值与影像几何偏移值可以由下式表达：

f(t)=A_fcos(ωt+f₀)                         (7)

f(t+Δt)=A_fcos(ωt+f₀+Δf)                 (8)

式中Δf=ωΔt；

由式(7)和式(8)可得：

即可得影像几何偏移值的振幅A_f：

由同频合振动的初相位特性可得：

由式(11)可得影像几何偏移值的初相位f₀：

多个频率分量时影像几何偏移值的计算如下：

多频率分量情况的处理方法的基本思路是利用傅里叶变换的平移特性，在频域下计算差分前的影像几何偏移值。通过频谱分析获得视差偏差值各频率分量幅值

频率ω_i、初相位根据同频合振动的性质，影像几何偏移值的各频率等于视差偏差值的频率，即ω_i。

对式(1)两边同时做离散傅里叶变换(DFT)，得到下式：

式中F表示离散傅里叶变换。根据离散傅里叶变换的性质，对式(12)变形后可得下式：

式中

为虚数单位，

k为DFT后频域下对应t时刻的虚数，N为参与DFT的数据总数。

由式(14)可得：

多频率分量情况下获得的影像几何偏移值f(t)，仍为多个频率分量的叠合。此时针对其各个频率分量，再对f(t)进行频谱分析，由频谱获取其各个频率分量的幅值A_fi和频率ω_i，由相位谱获取其各个频率分量的初相位f_0i。

这种方法亦适用于单频率分量情况的计算。

星体在垂轨和沿轨方向的姿态颤振值

(t)和

(t)与颤振引起的影像几何位移值f(t)有以下关系：

$\begin{array}{c}{\∂g}_{\mathrm{cross}}\left(t\right)=\frac{f_x\left(t\right)}{f_c}\\ {\∂g}_{\mathrm{along}}\left(t\right)=\frac{f_y\left(t\right)}{f_c}\end{array}---\left(16\right)$

式中，f_x(t)和f_y(t)分别为f(t)在像平面行方向和列方向的分量，f_c为相机主距。

以下分别使用国产测绘卫星的数据对其姿态颤振进行了探测，并对其结果进行了分析。

3.1实验数据

国产测绘卫星的传感器视差由其多光谱传感器构成。其多光谱传感器包含三片CCD，使用透反结合的方式如图5所示安置于焦平面上。每片CCD的焦平面上平行安置了4个波段的CCD线阵，如图6所示。每片CCD上平行安置的CCD线阵的成像时刻具有时间间隔，其构成了颤振探测所需的微小视差。

本实验所使用的实验数据为河南某地区的多光谱零级未拼接的影像数据。影像地面分辨率为5.8m，包含G、B、R和NR四个波段影像。单片CCD影像幅宽为3072像素。整景影像快视图如图7所示。

3.2实验结果

对三片CCD影像的颤振探测结果如表1以及图8～图13所示。实验结果表明本发明方法在国产高分辨率测绘卫星上取得了良好的效果。

Claims (4)

1.一种基于成像传感器视差的卫星姿态颤振探测方法，其特征在于，包括以下步骤：

1)获取视差影像，并对其进行预处理，所述的视差影像包括相互间具有观测视差的参考影像和待匹配影像；

2)对视差影像进行精确密集匹配处理，获得两幅影像的视差偏差值；

3)根据视差偏差值计算卫星姿态颤振值。

2.根据权利要求1所述的一种基于成像传感器视差的卫星姿态颤振探测方法，其特征在于，步骤1)中，所述的预处理包括影像增强和视差影像预配准。

3.根据权利要求1所述的一种基于成像传感器视差的卫星姿态颤振探测方法，其特征在于，步骤2)中，所述的精确密集匹配具体为：对两幅影像中相同地物对应的各影像行，以等间隔均匀取点的方式进行逐点匹配，首先采用基于灰度互相关的方法进行粗匹配，获取整像素匹配点，再以粗匹配结果作为初值，使用基于二维曲面峰值拟合的相位相关的方法进行精确匹配，获取亚像素匹配点，对整行影像所有匹配点的结果取平均值，作为整行影像的匹配结果，得到各行影像的视差偏差值。

4.根据权利要求1所述的一种基于成像传感器视差的卫星姿态颤振探测方法，其特征在于，步骤3)中，所述的计算卫星姿态颤振值具体为：

301)使用频谱分析的方法由步骤2)获得的视差偏差值得到卫星姿态颤振的频率和幅值，当颤振只含有一种频率分量时，执行步骤302)，当颤振含有多种频率分量时，执行步骤303)；

302)通过频谱分析获得视差偏差值的幅值

和频率ω，并由相位谱获取其初相位将t时刻的视差偏差值与影像几何偏移值表示为：

f(t)=A_fcos(ωt+f₀)

f(t+Δt)=A_fcos(ωt+f₀+Δf)

式中，A_f、f₀分别为影像几何偏移值的幅值和初相位，Δt为视差影像的成像时间间隔，Δf=ωΔt；

根据

计算得：

从而获得影像几何偏移值f(t)；

303)通过频谱分析获得视差偏差值各频率分量的幅值

频率ω_i和初相位

影像几何偏移值的各频率等于视差偏差值的频率；

对式

进行离散傅里叶变换，获得：

式中，

为虚数单位，

k为离散傅里叶变换后频域下对应t时刻的虚数，N为参与离散傅里叶变换的数据总数；

则

对f(t)进行频谱分析，获取影像几何偏移值各个频率分量的幅值A_fi和频率ω_i，由相位谱获取各个频率分量的初相位f_0i；

从而获得影像几何偏移值f(t)：

M为频率分量个数；

304)根据影像几何位移值f(t)计算星体在垂轨方向和沿轨方向的姿态颤振值

(t)和

(t)：

${\∂g}_{\mathrm{cross}}\left(t\right)=\frac{f_x\left(t\right)}{f_c}$

${\∂g}_{\mathrm{along}}\left(t\right)=\frac{f_y\left(t\right)}{f_c}$

式中，f_x(t)和f_y(t)分别为f(t)在像平面行方向和列方向的分量，f_c为相机主距。

Images