改进的图像重建迭代算法

张兆田; 张朋

摘要: 针对不完全数据图像重建问题,该文从积分方程角度提出一种改进的图像重建迭代算法(Im-proved Iterative Reprojection Reconstruction，IIRR)，并证明其在范数下收敛。该算法的图像重建收敛速度和重建图像误差依赖于引入的参数因子、已知的投影数据与图像的先验信息。重建图像可以表示为由已知投影数据所构成图像的级数形式。数值模拟结果表明，IIRR算法在不完全数据的情况下具有较好的图像重建能力。

关键词:

图像重建; 迭代算法; 不完全数据

Abstract: An Improved Iterative Reprojection Reconstruction (IIRR) from a view of in-tegral equation for incomplete image reconstruction is presented in this paper. IIRR is convergence under L~2 norm. The convergent velocity and the error of the reconstructed image depend on the relaxed parameter, the detected data and some prior information on image. The reconstructed image can be expressed in term of a series of image reconstructed by the detected data. Numerical simulation also shows that IIRR is effective to image reconstruction.

1. 引言

随着空间信息技术的发展进步，通信、导航、遥感等各种空间应用系统得以在全球范围内普及，正广泛、深入、潜移默化地改变着人类的生产生活方式，空间信息网络作为信息时代的国家公共基础设施，将成为保障“海洋远边疆、太空高边疆、网络新边疆”的重要支撑，是各国竞相争夺的战略制高点。高分辨遥感卫星具有高清分辨率、大覆盖范围、成像速度快、不受地域限制等特点，可为空间信息精细化获取提供重要服务。

从最早发射数百米空间分辨率的Landsat陆地观测系列卫星，到如今遥感卫星进入亚米级分辨率时代。美国的GeoEye, WorldView、法国Pleiades、中国的资源系列、高分系列、遥感系列等国内外高分辨率光学卫星相继成功发射，卫星影像最高分辨率达到0.3 m甚至更高，为全球对地观测提供高质量的数据产品。目前，多数高分率光学遥感卫星成像模式主要为相机线阵推扫成像或面阵相机推帧成像中的一种，获取大区域的长条带影像。但单一成像模越来越难以满足用户的多样化需求。为了实现对目标区域的精细化观测，获取现势性强、精度高的测绘遥感数据，各国竞相发展高分率多模式成像卫星。如2012年和2014年发射的SPOT-6/7双子卫星，采用同轨前、后视立体或前、下、后视三视立体成像，具有60 km的大幅宽，获取1.5 m空间分辨率全色影像和6 m空间分辨率多光谱影像^[1]；Pleiades-1A/B卫星具备快速机动与稳定控制能力，可短时间调整观测角度实现对不同目标观测^[2]，主要用于大面积区域测绘，以及矿业、工业、军事区域及自然灾害的监测等；我国2019年发射的高分七号卫星是中国首颗民用亚米级高分辨率立体测绘卫星，双线阵立体相机的空间分辨率达0.65 m，能够用于中国1∶10000立体测图及更大比例尺基础地理信息产品的更新^[3,4]；在商业卫星中，以北京三号卫星为代表，同时具备多目标成像、条带拼幅成像、立体成像和沿迹成像等模式，提高了卫星观测范围、获取信息效率和立体覆盖能力^[5]。

为了提升全球范围快速响应和空间信息时空连续支撑能力，2013年国家自然科学基金委立项“空间信息网络基础理论与关键技术”重大研究计划，以空间信息网络为核心开展交叉学科研究，支持了国家自然科学基金重大研究计划集成项目“天基信息网络在轨处理与实时传输综合集成演示验证”，推动了首颗互联网智能遥感卫星双清一号(珞珈三号01星)及其实时服务系统的研制，形成了开放式重大基础设施平台，为全国科研单位在遥感领域科研成果进行在轨验证提供服务。

2023年1月15日，双清一号在山西太原发射场成功发射，该星搭载了1台高分辨彩色面阵相机，可依据用户不同需求实现视频凝视、多角度立体和面阵推帧/面阵推扫成像。双清一号将作为国家自然科学基金委科学实验卫星平台为全国科研单位提供开放式卫星在轨实验验证环境；该星还将作为烟台东方航天港立项建设的东方慧眼星座工程首颗实验验证卫星，为后续东方慧眼星座工程的建设提供技术支撑。本文介绍了双清一号卫星的特点，公开了双清一号多模式成像样例数据集，研究了该数据集的构建流程与关键算法，该数据集是通过信号数据解码、Bayer插值、相对辐射校正、传感器校正、视频稳像和3维重建等处理步骤，遵循高精度遥感数据处理准则生成，空间分辨率为0.7 m，最后对数据集进行了可视化展示和定量化精度评价。

2. 双清一号性能与特点

双清一号互联网智能遥感试验卫星是一颗以0.7 m视频相机为主载荷，具备高性能实时处理、开放软件平台、在轨智能规划、在轨智能处理、灵活APP(Mobile Application)上注更新特征的智能遥感试验卫星。该星搭载了具有面阵推扫、面阵推帧和视频凝视多种成像模式的轻小型遥感相机，实现了光学遥感卫星单一相机的多模式体制成像，可依据用户不同需求提供单张图像数据、凝视视频数据、3维立体数据等多种类型的观测数据^[6,7]。表1列举了双清一号数据的技术指标。

表 1 双清一号数据技术指标

参数	指标数值
成像模式	视频凝视/面阵推帧/面阵推扫
空间分辨率	0.7 m@500 km
地面幅宽	双CMOS优于10 km
视频帧频	2～12 Hz
图像格式	Bayer彩色图像
量化位数	8 bit/12 bit
几何定位精度	优于30 m
相对辐射精度	优于3%

下载: 导出CSV

| 显示表格

(1) 视频凝视成像模式

该卫星平台机动能力强、指向精度和稳定度高，在此基础上可以进行敏捷的高清视频凝视成像。相机的光轴始终指向固定的地面目标点，在成像时间段内通过实时调整卫星姿态对目标点实现连续快速跟踪，满足对热点区域进行最多90 s的长时间连续观测。图1为对地凝视成像的工作示意图。

图 1 凝视视频成像模式

下载: 全尺寸图片幻灯片

(2) 面阵推帧成像模式

该卫星面阵推帧成像时，通过姿控控制光轴指向沿卫星飞行方向在地面匀速扫描，相对地速小于卫星飞行速度，通过控制相机可以进行推帧成像，每帧图像之间有一定重叠，多帧成像构成一个条带。图2给出了面阵推帧成像模式示意图。

图 2 面阵推帧成像模式

下载: 全尺寸图片幻灯片

(3) 面阵推扫成像模式

该卫星在面阵推扫成像时，保持在一个固定的姿态，采用亚米级彩色面阵数字TDI (Time Delay Integration)方法，进行双CMOS(Complementary Metal-Oxide Semiconductor)组合成像，可以等效为“线阵推扫”成像，实现高分辨率彩色面阵推扫功能获取连续影像，满足大区域覆盖要求。图3为面阵推扫成像的工作示意图。

图 3 面阵推扫成像模式

下载: 全尺寸图片幻灯片

3. 数据集构建

双清一号多模式成像数据集通过天基信息实时智能服务系统软件，对目标区域进行任务规划，采集遥感成像数据。一方面，可以将原始数据下传至地面接收站，将数据通过5G光纤传输至大数据中心，发起产品生产处理；另一方面，也可以直接通过在轨计算单元对原始信号数据进行解码、高精度定位和智能处理，将在轨处理产品传输至地面用户，整个处理过程均遵循高精度遥感数据处理准则。如图4所示，样例数据集生成流程包括成像数据获取、高精度在轨定标和产品生产3个关键步骤。

图 4 多模式成像样例数据集生成流程

下载: 全尺寸图片幻灯片

3.1 成像数据获取

高效获取成像数据是卫星发挥服务的关键，本文对卫星任务规划方案进行了实验分析。首先，通过网上开源的两行根数计算其卫星轨道数据，根据侧摆能力计算其星下点可覆盖区域。然后，在卫星星下点可覆盖区域，依次查看该区域内优先级高(包含动态添加的待拍摄区域、全球主要城市、机场、港口等)的目标，并在卫星过境时间查看目标区域对应时间的云层情况，若该区域内云量较少，则将该目标纳入待观测集合；否则查看下一个目标。当待观测集合目标数目已达当日可下传任务数量上限时，停止本次任务规划。最后，根据待观测集合内各目标的经纬度坐标等信息生成对应的成像指令来进行成像。

3.2 高精度在轨定标

卫星发射前对星上载荷参数通常会进行严格的实验室检校，然而由于卫星发射和运行过程中空间环境改变等因素的影响，星上载荷结构和状态发生改变，导致实验室检校参数在轨后与真实参数存在较大的误差，直接影响产品处理精度^[8]。因此，需要开展在轨相对辐射定标与几何定标，来满足高精度的影像处理和应用需求。

3.2.1 相对辐射定标

相对辐射定标首先分析影响探元光电响应变化的各种因素，并以最大后验概率理论为依托，通过引入约束条件和加权思想构建在轨自适应的系统辐射校正模型，利用定标数据或者实时获取的遥感数据，构建统计样本的智能筛选模型，并采用增量统计的策略，优化辐射定标参数的解算精度。

在地面辐射定标的基础上，将定标查找表上注更新至星上，星上根据查找表进行相对辐射校正，得到校正后影像，如图5所示。

图 5 在轨相对辐射定标技术路线

下载: 全尺寸图片幻灯片

(1) 线阵成像相对辐射定标

对于等效线阵推扫成像体制的相对辐射定标，本文采用在轨统计定标方法。利用卫星在轨后获取的常规影像数据，统计每个探元出现的灰度频次。基于大数定理，样本趋于无穷多时，随机因素(地形、光照、大气等)的影响会趋向于0，器件的响应非均匀性会得到有效改善。因此，在轨统计相对辐射定标方法以传感器的整体直方图作为标准校正每个探元的响应函数，即可获得相对辐射定标系数。

在轨统计相对辐射定标处理步骤主要包括：

(a) 计算各探元累积概率密度函数

选取多组具有代表性的样本条带数据，统计各像素灰度值出现的次数，根据式(1)所示，计算出每个探元传感器的累积概率密度函数 ${S_k}(i)$ 。

$\left. \begin{gathered} {P_{_k}}(i) = {m_k}(i)/M(i) \\ {S_k}(i) = \sum\limits_{j = 0}^k {{P_k}(i)} \\ \end{gathered} \right\}$

(1)

其中， ${m_k}(i)$ 为第i个探元灰度级等于 $k$ 的像素数目， $M(i)$ 为对应探元直方图的总像素数， ${P_k}(i)$ 为直方图中像素的灰度级为k对应的概率密度函数。

(b) 计算综合探元直方图

将每个探元的综合直方图进行合并就可以获得所有探元的综合直方图。对于直方图匹配算法来说，所有探元的综合直方图就是期望直方图，建立查找表的原理是使匹配处理后每个探元的综合直方图的概率密度函数和期望直方图的概率密度函数相同。

(c) 生成查找表。设期望直方图灰度级为的累积概率密度函数为 ${V_l}$ ，某探元灰度级为 $k$ 的累积概率密度函数为 ${S_k}$ ，那么如果满足以下关系

$\left. \begin{gathered} {V_l} \le {S_k} \le {V_{l + 1}} \\ \left| {{V_l} - {S_k}} \right| \le \left| {{V_{l + 1}} - {S_k}} \right| \\ \end{gathered} \right\}$

(2)

用灰度级 $l$ 代替灰度级 $k$ ，否则，就用灰度级 $l + 1$ 代替灰度级 $k + 1$ 。用这种方法处理所有探元就可以得到相应的灰度范围查找表，最后对统计过程中超过响应区间的查找表值予以剔除。

(2) 面阵成像相对辐射定标

不同于线阵推扫成像体制，面阵成像相对辐射定标是对像素排列成维矩阵的图像传感器进行相对辐射定标，其成像探元数是线阵探元数百倍以上，基于线阵统计的方法将难以适用。本文针对面阵视频和面阵推帧的成像特点，采用卫星对夜间海洋成像模式来标定面阵传感器暗电流误差，并对沙漠、云、雪等均匀场地物进行序列帧成像，实现面阵传感器的高精度增益定标，以完成卫星面阵传感器相机的在轨相对辐射定标。

(a) 暗电流标定。卫星在轨后对深空成像数据或夜间海洋成像多帧序列图像数据，标定传感器各个探元暗电流，如式(3)所示

${\rm{offset}}(i) = \sum\limits_{n = 1}^N {{\rm{DN}}(n,i)}$

(3)

(b) 增益标定。通过卫星对匀场景区域(均匀海洋、雪、沙漠等特殊均匀地物场景)成像获取视频多帧序列图像数据，对视频多帧序列图像数据进行暗电流校正后，以传感器单波段有效探元(除异常探元外)均值为相对辐射定标基准，实现传感器各波段探元间响应不一致性标定。

首先，对多帧序列进行暗电流去除并进行探元序列帧平均，如式(4)所示

${{\overline{\rm{DN}}_i}} = 1/N \times \sum\limits_{n = 1}^N {({\rm{DN}}(n,i) - {\rm{offset}}(i))}$

(4)

然后，计算各面阵探元的相对增益定标基准，如式(5)所示

$\overline {{\rm{DN}}} = 1/M \times \sum\limits_{i = 1}^M {{{\overline {{\rm{DN}}} }_i}}$

(5)

最后，生成各面阵探元传感器的增益系数，如式(6)所示

${\rm{gain}}(i) = \frac{{\overline {{\rm{DN}}} }}{{{{\overline {{\rm{DN}}} }_i}}}$

(6)

式(3)–式(6)中，N为总帧数，M为总探元数， ${\rm{offset}}(i)$ 为定后探元i ( $0 \le i < M$ )的暗电流标偏置， ${\rm{DN}}(n,i)$ 为第n帧第i个探元的灰度值， ${{\overline{\rm{DN}}_i}}$ 为探元i在暗电流校正后的序列帧灰度均值， $\overline {{\rm{DN}}}$ 为各面阵探元的相对增益定标基准， ${\rm{gain}}(i)$ 为标定后的增益系数。

3.2.2 在轨几何定标

在轨几何定标旨在消除卫星平台外部系统误差(如相机安装角在卫星发射过程中受空间力学环境影响相对于实验室检校值的形变、空间复杂热环境下引起的长周期性变形误差)和相机内部系统误差(主点主距检校误差及物镜光学畸变、CCD畸变)，从而保证影像产品的几何质量^[9,10]。

对于光学卫星影像，通常基于共线方程建立严格几何成像模型，以对各项系统误差构建严格的数学模型。在严格几何成像模型中各项系统误差参数均具有其严格的物理意义，但由于内外系统误差参数强相关、过渡参数化，导致各项系统误差参数难以精确标定。因此本文根据双清一号相机成像的特点，采用了一种基于指向角模型的几何定标模型，结合姿轨观测数据，利用地面控制信息解求内外定标模型参数，对各探元在相机坐标系下的指向值进行拟合，来实现高精度的在轨定标。

基于2维探元指向角模型的在轨几何定标的流程如图6所示，主要包括：

图 6 基于探元指向角模型的在轨几何定标流程

下载: 全尺寸图片幻灯片

(1) 几何定标模型构建

为了克服严密物理量测模型中参数高度相关、部分参数低显著的问题，基于多项式指向角的成像模型(线阵成像采用一元多项式模型，面阵成像采用二元多项式模型)可以有效地作为严密物理量测模型的替代。通过在轨几何内定标可以确定阵传感器上每个探元在参考相机坐标系下的精确光线指向，以克服相机内部畸变的影响。

(2) 控制点匹配。传统的几何定标场可以提供高精度数字正射影像(DOM)与数字高程模型(DEM)等参考数据，通常采用基于影像模拟的高精度匹配方法，从参考DOM与DEM影像中提取待定标影像上像方坐标所对应的地面控制点的物方坐标，从而实现控制点的自动量测。

(3) 定标参数的解算。目前光学遥感卫星在轨几何定标一般采用空间后方交会的方法，结合姿轨观测数据，利用地面控制信息解求内外定标模型参数。外定标参数与内定标参数虽然代表着不同的物理含义，但是数学上他们仍然是高度相关的，若同时将其作为未知数进行求解，虽不会影响方程最后的收敛解算，但求得的外定标参数与内定标参数将失去其原本该有的物理意义。因此，本文设计了定标参数的分步解算方法，在保证正确求解定标参数的前提下，最大限度地保留内外定标参数本身蕴含的物理意义，便于后期的分析应用。

3.3 产品生产

3.3.1 标准产品生产

在大温差、高辐照的空间环境下，星上处理平台的体积、功耗、性能受到严格限制，本文采用了面向星上兴趣区提取算法^[11]，实现对任务指定兴趣区标准影像产品的实时生产，为基于星上受限平台实现更丰富的实时应用和服务提供了一种可行的解决方案。

面向任务区域的星上提取算法以像平面像点的高精度几何定位算法为基础，根据目标区待定位的像平面点坐标 $(s,l)$ ，结合来自相机原始数据中的成像时间数据，得到像点所在原始成像行的精确时刻；根据成像时刻，对临近时间采样的姿态数据、轨道数据进行内插，得到成像时刻的姿态、轨道信息；根据姿态、轨道信息及上注的相机几何参数可构建该成像行的严格几何模型；根据严格几何模型和成像探元序号 $s$ ，即可求得像平面点对应的地理坐标 $(B,L,H)$ 。

光学卫星遥感影像的标准产品生产需要进行相对辐射校正和传感器校正。传感器校正是指针对传感器设计与成像特点，通过建立原始影像与校正后影像几何映射模型，实现多片CMOS拼接和多光谱影像波段配准，得到无缝完整无畸变影像及其高精度RFM，便于后续控制点搜索及快速正射校正提供初始的模型基础^[10,11]。因此，高精度传感器校正是进行高精度测绘应用的关键。

区别于传统以景为单位的影像传感器校正，兴趣区影像传感器校正在成像时间和成像空间范围上更加灵活，因此，对兴趣区影像进行传感器校正时，需要依据成像范围动态构建传感器模型。为了保证兴趣区传感器校正精度，在兴趣区成像范围内还需要要考虑平台震颤引起的姿态抖动、影像变形问题；在兴趣区覆盖多片CCD影像时需要考虑影像拼接问题；对于多光谱影像，还需要考虑多光谱影像配准问题；以及镜头畸变和CCD变形引起的影像畸变问题等。为了校正这些因素带来了的影响，本文拟基于虚拟CMOS稳态重成像原理在兴趣区成像范围内动态建立高精度的传感器校正模型，算法过程如算法1。利用原始严密成像模型和传感器校正模型的几何定位一致性，建立原始影像与虚拟稳态重成像的坐标映射关系，校正原始影像变形，并解决影像拼接、多光谱影像配准的问题，同时可获取兴趣区传感器校正影像的高精度有理函数模型，便于后续处理与应用。

表 1 传感器校正算法

输入：分片影像、成像时间、姿轨数据、地球自转参数、高程　dem等
输出：有理函数模型RFM和虚拟大影像
1 分别构建单片原始影像和待校正虚拟影像的严格几何成像模　型：　 $\left[ \begin{gathered} X \\ Y \\ Z \\ \end{gathered} \right] = \dfrac{1}{\lambda } \cdot R_{ {\rm{J2000} } }^{ {\rm{WGS84} } }(t)R_{ {\rm{body} } }^{ {\rm{J2000} } }(t)R_{ {\rm{senor} } }^{ {\rm{body} } }\left[ \begin{gathered} x + {x_0} + \Delta x \\ y + {y_0} + \Delta y \\ f \\ \end{gathered} \right] + \left[ \begin{gathered} {X_{ {\rm{GPS} } } }(t) \\ {Y_{ {\rm{GPS} } } }(t) \\ {Z_{ {\rm{GPS} } } }(t) \\ \end{gathered} \right]$
2 建立分片原始影像与待校正虚拟影像的坐标映射关系；
${(s,l)}_{ {\rm{ori} } }\iff {(x,y)}_{{\rm{vir}}},其中{(s,l)}_{ {\rm{ori} } }\Rightarrow (B,L,H)\Rightarrow {(x,y)}_{ {\rm{dst} } }$
3 影像重采样(CUDA并行)
3.1 对单片原始影像进行格网划分，每个格网(grid)包含　　 $(w \times h)$ 个像素
3.2 设定CUDA线程参数，包括：
(a) 每个CUDA格网(grid)中的块(block)数目： $({B_x},{B_y})$ ;
(b) 每个CUDA块(block)中的线程(thread)数目： $({T_x},{Y_y})$ ;
(c) 每个CUDA线程(thread)中处理的像素数目: $({P_x},{P_y})$
确保 $({B_x} \times {T_x} \times {P_x},{B_y} \times {T_y} \times {P_y}) \ge (w \times h)$
3.3 每个格网(grid)中的 $({B_x} \times {T_x},{B_y} \times {T_y})$ 个线程(thread) 　并行执行：
	每个线程处理 $({P_x},{P_y})$ 个像素
	每个像素
		(a) 计算虚拟影像像点 ${(x,y)_{{\rm{vir}}}}$ 所对应单片原始影像像　　点坐标 ${(s,l)_{{\rm{ori}}}}$ ;
		(b) 分别计算像点 ${(s,l)_{{\rm{ori}}} }$ 邻域像素的灰度值；
		(c) 利用多项式重采样方法，计算出 ${(s,l)_{{\rm{ori}}} }$ 灰度值.
	像素循环结束
影像格网循环结束

下载: 导出CSV

| 显示表格

3.3.2 视频物方稳像

针对卫星相机光轴指向变化导致卫星视频场景不断变化的问题，采用一种带有地理编码的卫星视频物方稳像方法。该方法采用逐帧运动估计和补偿的策略，在卫星视频帧序列中依次以相邻两帧的前一帧为主帧、后一帧为辅帧，估计辅帧相对于主帧的帧间运动参数，并对辅帧定向参数进行运动补偿，最后在物方空间对视频序列进行地理编码，实现带有地理编码的光学视频卫星物方稳像处理^[12]，流程如图7所示。

图 7 带有地理编码的卫星视频物方稳像流程

下载: 全尺寸图片幻灯片

(1) 帧间运动模型构建

对于光学视频卫星稳像处理而言，若以视频序列中的固定帧为主帧、其余帧为辅帧，则需要将所有的辅帧分别与主帧进行同名像点匹配。然而，受卫星平台抖动和姿态控制误差等因素的影响，部分辅帧与主帧之间的重叠度可能较小，以至于难以获得分布均匀的同名像点。此外，视频卫星的成像时间越长，成像视角的差异也随之增大，这必将导致视频帧间几何变形很大，给主帧和辅帧之间的同名像点匹配带来难度，以至于难以实现子像素级帧间运动估计。因此，本文采用逐帧运动估计策略，充分利用卫星视频帧的定向参数信息，在辅帧定向模型的基础上引入帧间运动参数，以消除帧间相对误差的影响。

(2) 帧间运动参数求解

卫星视频帧间运动参数的精确求解是帧间运动估计的另一个重要环节。在卫星视频帧定向模型和帧间运动模型的基础上，以卫星视频帧序列中相邻两帧的前一帧为主帧、后一帧为辅帧，可以明显减小主帧和辅帧之间的几何变形，为帧间同名像点精确匹配提供便利，但在求解帧间运动参数时，主帧和辅帧之间的交会角会非常小。若采用传统区域网平差方法求解帧间运动参数，交会角过小会使得高程误差被异常放大，导致平差不收敛。因此，本文通过引入DEM作为高程约束条件，求解辅帧相对于主帧的运动参数。

(3) 辅帧定向参数运动补偿

求解出当前辅帧相对于主帧的运动参数后，当前辅帧将作为视频序列中下一个相邻两帧运动估计的主帧。为了便于构建下一个相邻两帧的帧间运动模型，同时便于后续卫星视频帧的地理编码，需要对当前辅帧的定向参数进行运动补偿。首先，根据当前辅帧原始定向参数构建定向模型，并生成虚拟控制格网点。其次，利用当前辅帧相对于主帧的运动参数，根据帧间运动模型，对虚拟格网点的像点坐标进行帧间运动补偿。然后，利用运动补偿后的虚拟格网点重新求解辅帧定向参数，得到运动补偿后的定向参数。

(4) 视频序列物方地理编码

通过卫星视频帧间运动估计与帧间运动补偿处理，可以消除卫星平台抖动和卫星定轨测姿等误差的影响，恢复卫星视频帧间的精确几何关系。在此基础上，本文对视频序列中的每一帧进行物方地理编码处理，从而得到带有地理编码的流畅视频。

4. 数据集展示与分析

4.1 标准产品影像

双清一号采用彩色Bayer型探测系统，极大减少了数据量，降低了数传的压力，同时还保证了各颜色谱段的高分辨率，融合彩色影像后色彩纹理清晰，层次丰富。如图8所示，双清一号卫星原始影像是单波段的Bayer格式，在经过插值处理，可以计算出每个像素的RGB值，恢复出颜色鲜明的彩色图像。

图 8 Bayer插值前后影像对比

下载: 全尺寸图片幻灯片

双清一号卫星搭载了两片CMOS成像器件，采用光学拼接的设计，在影像边缘处由于光路的遮挡，会有明显的渐晕区，如图9(a)和图9(c)所示。此外由于探元响应不一致，相对辐射校正前影像还存在有明显的条带噪声，在相对辐射校正后，如图9(b)、图9(d)所示，包括边缘部分在内的影像的条带噪声全部得以消除，目视效果显著提升。

图 9 相对辐射校正前后影像对比(意大利比纳斯科)

下载: 全尺寸图片幻灯片

双清一号卫星搭载了两台高分辨率相机以有效增加成像视场，达到10 km幅宽，如果不能将两台相机进行有效拼接，且生成高精度的RFM模型，将直接影响双清一号卫星数据的有效应用。如图10所示，通过基于虚拟相机的双相机传感器校正方法，能有效地实现原始双相机影像的无缝拼接，且校正后影像的精度与原始单相机影像精度一致。

图 10 传感器校正结果(意大利比纳斯科)

下载: 全尺寸图片幻灯片

4.2 凝视视频数据

珞珈三号01星在实现凝视的过程中，需要不停地调整成像姿态，在成像过程中不可避免地会受到各方面因素的影响，导致视频不能稳定地成像，在结果上呈现晃动的现象，因此需要对其进行稳像处理。通过稳像处理可以发现将不同的影像背景进行对比，变动的只有前景部分。除此之外，由于成像角度变化过大，观测比较高的地物时能够发现成像角度的变化。图11展示了卫星在视频成像过程中，短时间序列下不同角度观测的地面目标。

图 11 视频稳像产品示意图(上海市黄浦江)

注：该图为视频稳像产品每隔8 s抽取的单帧影像

下载: 全尺寸图片幻灯片

帧间几何精度评定步骤主要包括^[12]：(1)根据主帧定向参数构建定向模型，并根据定向模型将主帧上的每一个同名像点投影至DEM 上，得到物方投影点；(2)根据辅帧定向参数构建定向模型，并根据定向模型将每一个物方投影点投影至像方，得到像方投影点；(3)根据辅帧上的像方投影点坐标及其匹配坐标，统计两者坐标残差中误差，作为帧间几何精度。

为了验证本文卫星视频稳像方法的正确性，利用2023年2月28日拍摄的上海市黄浦江卫星视频开展了实验验证。其中，实验数据的视频拍摄时长约为34 s、共计202 帧。实验数据的示意图如图12所示。

图 12 稳像前后的卫星视频帧间几何精度

下载: 全尺寸图片幻灯片

受视频成像期间卫星平台抖动和卫星定轨测姿等误差的影响，视频序列中相邻两帧之间存在着明显抖动，从图12可以看出，帧间几何误差最大近6像素。由此可见，原始卫星视频数据的帧间几何精度并未达到子像素级，这将严重制约卫星视频的高精度应用。利用视频几何稳像对辅帧进行运动估计和运动补偿后，可以有效消除卫星平台抖动和卫星定轨测姿等误差的影响，显著提高帧间几何精度，最大误差提升至低于0.3像素。

4.3 实景3维产品

实景3维作为真实、立体、时序化反映人类生产、生活和生态空间的时空信息，是国家重要的新型基础设施，可以为数字中国提供统一的空间定位框架和分析基础，是数字政府、数字经济重要的战略性数据资源和生产要素数据，在城市规划与管理、自然灾害评估等多个领域发挥重要作用。3维立体建模的关键技术包括区域网平差、密集匹配、构网及纹理映射等步骤^[5]。

(1) 区域网平差

该过程首先利用SIFT算子提取影像特征点，并完成影像之间的匹配。在匹配的连接点基础上，区域网平差过程可采用在RFM模型变换上叠加像方平移、仿射变换来实现。

(2) 密集匹配

通过两视图立体匹配、多视角点云融合等步骤，实现被摄区域稠密点云的生成。两视图立体匹配采用SGM算法，通过寻找具有图像一致性的区域，生成影像对稠密点云，然后使用一定的内插和滤波算法得到数字表面模型DSM，并确定像素点相对于相机平面的深度(3维点坐标)，最后对多视角的DSM以及点云进行融合处理，获取被摄区域的DSM模型及稠密点云。

(3) 构网及纹理映射

在DSM产品基础上构建Delaunay三角网，通过三角面片剖分、纹理映射生成高质量的实景3维模型。对于可能出现基于卫星影像构建的建筑物3维模型侧面纹理不够清晰，还需要利用地面影像纹理贴图或者纹理编辑和色彩处理等手段^[13]，实现对地物纹理特征的恢复。

双清一号具备高敏捷、高稳定、高精度获取影像能力，在视频凝视的成像模式下，可实现同轨对一个区域进行多次多角度观测。通过实景3维重建处理后，可构建满足不同层次的3维地理场景，这将在城市规划与管理、自然灾害评估等多个领域发挥重要作用。图13是土耳其哈塔伊地区多角度3维立体重建结果示意图。

图 13 多角度3维立体重建(土耳其哈塔伊地区)

下载: 全尺寸图片幻灯片

5. 结论

目前，多数高分率遥感卫星成像模式单一，在变化检测、动目标跟踪和3维建模等处理方面受到诸多限制，难以满足多样化的应用服务需求。因此，该文公开了双清一号多模式成像样例数据集，空间分辨率为0.7 m，该数据集是由天基信息实时智能服务系统经过信号数据解码、Bayer插值、相对辐射校正、几何定位、视频稳像和3维重建等处理步骤，遵循高精度遥感数据处理准则生成的数据集，覆盖了意大利、土耳其、阿联酋和中国等地区，包括单张高分辨率图像、凝视视频、3维重建等多模式观测数据及产品。同时，进一步对构建的数据集进行了可视化和质量分析。该数据集可广泛用于车流量管理、变化监测、3维模型构建等应用场景。

致谢感谢国家自然科学基金委员会重大研究计划集成项目支持，航天东方红卫星公司卫星研制支持，东方航天港卫星发射支持，以及珞珈三号团队其他老师与同学所做的大量工作！

Papoulis A. A new algorithm in spectra analysis and band-limited signal extrapolation. IEEE Trans. on Circuits Syst., 1975, CAS-22(9): 735-742.[2]Papoulis A. Signal Analysis. New York: McGraw-Hill, 1977: 205-215.[3]Gerchberg R W. Super-resolution through error energy reduction. Opt. Acta, 1974, 21(9): 709-720.[4]Gerchberg R W, Saxton W O. A practical algorithm for the determination of phase from image and diffraction plane pictures. Optik, 1972, 35(2): 237-246.[5]Chamzas C, Xu W Y. An improved version of Papoulis-Gerchberg algorithm on band-limited extrapolation. IEEE Trans. on Acoustics, Speech. Signal Processing, 1984, ASSP-32(2): 437-440.[6]Xu W Y, Chamzas C. On the extrapolation of band-limited functions with energy constraints[J].IEEE Trans. on Acoustic, Speech Signal Processing.1983, 31(5):1222-1234[7]Sanz J L, Huang T S. Discrete and continuous band-limited signal extrapolation. IEEE Trans.on Acoustic, Speech Signal Processing, 1983, ASSP-31(5): 1276-1285.[8]Sato T, Norton S J, et al.. Tomographic image reconstruction from limited projections using iterative revisions in image and transform space[J].Applied Optics.1981, 20(3):395-399[9]Peng H, Stark H. One-step image reconstruction from incomplete data in computer tomography.IEEE Trans. on Medical Imaging, 1989, MI-8(1): 16-31.[10]Virador P R G, et al.. 3D reconstruction in PET cameras with irregular sampling and depth of interaction[J].IEEE Trans. on Nuclear Science.2001, 48(4):1524-1529[11]张兆田,邱佩璋.一类广义Radon变换的反演:迭代重建再投影的外插算法.CT理论与应用研究,1993,2(3):20-23.

姓名
邮箱
手机号码
标题
留言内容
验证码

留言板