《测绘学 》
构建与学术的桥梁 拉近与权威的距离
许丽娜, 何鲁晓
中国地质大学(武汉)地球物理与空间信息学院, 湖北 武汉 430074
收稿日期：2017-02-16；修回日期：2017-05-04
摘要：基于凸集投影（projection on convex sets，POCS）对GF-4数据进行超分辨率重建处理。利用能量连续降减法构造参考帧以减少迭代次数，提高算法效率；利用Vandewalle配准法对图像块进行运动估计从而获取精确的运动参数有利于图像的分块处理。通过对GF-4真实数据进行试验与分析，表明超分图像结果的清晰度与序列影像帧数呈正相关；信噪比与帧数呈负相关。应用5帧数据进行超分重建，清晰度（平均梯度）从7.803提高到14.386；信噪比会略有下降，从3.411下降到3.336。试验表明该方法可以有效提高图像清晰度，增加图像细节信息，提升图像判别能力，扩展数据的应用范围。
GF-4 Images Super Resolution Reconstruction Based on POCS
XU Lina, HE Luxiao
Abstract: The super resolution reconstruction of GF-4 is made by projection on convex sets (POCS). Papoulis-Gerchberg is used to construct reference frame which can reduce iteration and improve algorithm efficiency.Vandewalle is used to estimate motion parameter which is benefit to block process. Tested and analyzed by real GF-4 series images, it shows that sharpness of super resolution result is positive correlatie to frame amount, and signal to noise ratio (SNR) is negative correlate to frame amount. After processing by 5 frames, information entropy (IE) changes little; sharpness (average gradient) increases from 7.803 to 14.386; SNR reduces a little, from 3.411 to 3.336. The experiment shows that after super resolution reconstruction, sharpness and detail information of results can be greatly improved.
Key words: GF-4 super resolution POCS Papoulis-Gerchberg Vandewalle
地球静止轨道卫星具有对地观测位置相对固定、时间分辨率高、观测范围广等特点，是航天对地观测成像技术的重要组成部分，是对低轨观测的有益补充。GF-4是我国第1颗民用面阵静止轨道光学成像卫星，具有凝视成像、机动巡查、区域成像等多个模式，能够满足减灾、气象、地震、林业等多个领域的不同需求[1]。
GF-4采用面阵凝视成像模式，搭载有一台10 240×10 240的面阵CMOS探测器和一台1024×1024的面阵碲镉汞探测器。在可见光近红外通道，GF-4可以提供全色影像与多光谱影像，空间分辨率为50 m，重复成像周期为15 s，单帧影像区域大小为500 km×500 km；在中波红外通道，提供光谱范围为3.5 μm~4.1μm的中波红外影像，空间分辨率为400 m，重复成像周期为1 s，单帧区域大小为400 km×400 km。GF-4可以对地面进行连续观测，提供高时间分辨率的多帧序列影像。相较于低轨推扫式卫星，GF-4的空间分辨率依然较低，有必要利用其高时间分辨率的数据特点对其进行超分辨率重建，提高其空间分辨率，增强影像的可判读能力，扩大数据的使用范围。
超分辨率重建(super resolution，SR)的概念是由文献[2-3]于20世纪60年代提出的，指用若干张低分辨率图像重建生成一幅高分辨率图像，可以分为频率域与空间域两种。频率域方法主要是基于傅里叶变换的位移特性，通过在频率域消除频谱混叠现象来改善图像的空间分辨率的[4]。频率域方法原理简单，计算方便，但是难以建立复杂的运动模型与降质模型，造成其复原效果不佳[5]。空间域方法可以建立较为复杂的运动模型，使用范围广，具有较强的先验知识包含能力，使其逐渐取代了频率域方法。空间域方法主要包括非均匀样本内插法(non-uniform interpolation)[6]、迭代反向投影(iterative back projection, IBP)[7]、最大后验概率估计(maximum a posteriori, MAP)[8]、凸集投影[9-10]、混合MAP/POCS[11]多尺度细节增强[12]等。近年来，基于学习的超分辨率重建方法[13-14]在试验与实际应用中取得了较好的效果。
1 原理与方法1.1 凸集投影
凸集投影是一种空间域超分辨率算法[15]，原理直观，算法简单。超分辨率重建是一个病态求逆的过程，一般是通过正则化方法，利用图像的约束条件与先验信息来约束可行解的范围。低分序列影像的每一帧都可以认为是由高分图像经过一系列降质因素得到的，即每一帧低分图像都构成一个约束条件。低分图像的获取模型可以由下式表示
(1)
式中，gl(m1,m2)是第l帧低分图像；f(n1,n2)是高分图像；hl(m1,m2;n1,n2)是相应的降质函数，一般描述为点扩散函数；η1(m1,m2)是加性噪声。
凸集投影的基本原理就是求解闭合凸集的交集。假设有m个先验知识或约束条件，则对应有m个闭合凸集Ci,i=1, 2, …,m。对于每个闭合凸集Ci，记对应的投影算子为Pi，假设交集C0是非零解空间，那么理想高分图像f∈C0。因为影像超分重建结果包含于解空间，所以POCS算法得到的结果是不唯一的。其迭代过程可以表示为
式中，Ti是第i个先验知识所对应的松弛投影算子，定义为
式中，0＜λi＜2是松弛投影算子。
在考虑噪声影响时，需要在模型中加入噪声的先验知识。假设图像中的噪声是方差为σv2的高斯噪声，先验边界为cσv(c≥0)，决定了其统计置信度。对于图像中的各个像素应当满足以下公式
(4)
式中，r(y)(m1,m2)表示x(i1,i2)和g(m1,m2)的残差，与噪声统计特性相同。对于任意一个x(i1,i2)，它在Cm1,m2上的投影[16]可以表示如下
(5)
GF-4采用10 bit存储，将振幅也作为一个限制条件，如下式所示
图 1 总体技术路线图Fig. 1 Overall technology
图选项
1.2 能量连续降减法
能量连续降减法是一种基于频率域的外推方法[19-20]，也被称为Papoulis-Gerchberg。该方法能够消除下采样后信号的频率混叠现象。文献[21]提出了一维信号混叠的外推恢复算法。一维原始信号f(t)的傅里叶频谱为
(7)
对原始信号抽样得到的实际信号为
(8)
式中，rect为矩形函数；a用于调整采样频率。
(9)
用截止频率为ω0低通滤波，得到
(10)
于是
(11)
再将信号h0(t)加载在h1(t)上，得到
(12)
在对图像这种二维信号进行操作时，将低分序列影像的各个像素放置于对应的高分 格中作为已知点，其余点赋值为零。对图像进行低通滤波后(空间域操作时可以使用均值滤波)，把已知点重新加载覆盖到图像上，这样的操作作为一次迭代。这一方法的实质是将已知点的信息向周围点进行扩散，以填充未知点。以一定准则退出迭代后，即可产生一幅具有较高分辨率的图像。
1.3 Vandewalle配准法
Vandewalle配准法是一种频率域配准方法[22]。它使用图像的低频信息进行配准，因为噪声主要包含在高频信息中，而低频信息不包含走样信息[23]，具有更强的鲁棒性。记参考图像为f1(x)，待匹配图像为f2(x)，两者关系可表示为
式中，
。将f1(x)和f2(x)的傅里叶变换形式记作F1(u)和
(14)
因为|F1(u)|和|F2(u)|之间的关系与旋转角度有关，而与平移量无关，所以先进行旋转角度估计，再进行平移参数估计。在进行旋转角度估计时，定义一个关于旋转角度α的函数
(15)
由于|F(r,θ)|是一个离散信号，可以把h(α)每隔0.1°进行一次计算。此时只考虑ερ＜r＜ρ的圆形区域，其中ρ为图像半径，ε设为0.1。因为低频信息的值相对很大，所以可以忽略r＜ερ的值。
在进行平移参数估计时，平移参数是通过频率域中的一个线性相位移动来表示的，即
(16)
相位差∠(F2(u)/F1(u))是关于u的线性函数，斜率为2πΔx。平移参数Δx可以通过计算相位差平面斜率来获得。
1.4 超分辨率具体实现
对于参考帧构造，首先选取低分序列图像中的一帧P1，利用双线性内插法进行上采样，得到底图P′，使得其分辨率提高两倍。双线性内插法是一种重采样方法，使用投影位置周围4个像素作为采样数据，在x、y方向上分别进行一次插值，即
(17)
式中, F是输出图像；f是输入图像；(i+u,j+v)是(m,n)在输入影像上的映射坐标，i、j是整数坐标，u、v是0到1之间的小数坐标。同时建立一个大小与P′相同的高分 格Ptemp。分别计算其余各帧与P′的偏移量，并根据偏移量将P2、P3、…、Pn的像素填入Ptemp。将P′用以下规则更新
(18)
对上述步骤进行多次迭代操作后，即可得到参考帧Pref。
对于运动估计，首先对图像进行傅立叶变换，并提取其低频信息，利用低频信息计算图像之间的旋转参数与平移参数。根据旋转参数与平移参数计算P2、P3、…、Pn在Pref上的准确位置。
对于参考帧修正，首先应用高斯模型生成一个大小为5×5的模拟点扩散函数。高斯模型为
式中，u′、v′是PSF中心点的坐标值。归一化后，PSF如表 1所示。
表 1 点扩散函数Tab. 1 Point spread function
| 0.018 3 | 0.082 1 | 0.135 3 | 0.082 1 | 0.018 3 |
| 0.082 1 | 0.367 9 | 0.606 5 | 0.367 9 | 0.082 1 |
| 0.135 3 | 0.606 5 | 1.000 0 | 0.606 5 | 0.135 3 |
| 0.082 1 | 0.367 9 | 0.606 5 | 0.367 9 | 0.082 1 |
| 0.018 3 | 0.082 1 | 0.135 3 | 0.082 1 | 0.018 3 |
表选项
根据运动估计得到的具体坐标，提取5×5的图像块，以式(21)、(22) 计算像素的理想灰度值
(21)
(22)
误差E定义为
(23)
(24)
2 试验与分析2.1 试验数据
试验数据为5帧GF-4可见光近红外影像，拍摄于2016年3月6日，左上角坐标为100.440 9°E，37.382 8°N，右下角坐标为106.214 7°E，31.743 8°N。在图幅中，分布有山川、城镇、冰山、湖泊等多种地形地貌。从中截取2000×2000的图像作为输入数据，如图 2所示。
图 2 试验数据Fig. 2 Experiment data
图选项
2.2 评价标准
(25)
式中，M、N是图像的长度与宽度；f是测试图像。
信息熵代表图像信息量，数值越大越好。其定义为
(26)
式中，pi表示灰度值为i的像素数量占整幅图像的比例。
信噪比通过局部平均值与标准差法测量[25-26]，其值越大越好。信噪比的测量步骤为：① 求图像平均值M；② 对图像进行分块，求取局部均值与标准差；③ 得到局部标准差的最大值STDmax；④ 求解信噪比
(27)
2.3 结果分析
图 3分别是输入的序列图像、双线性内插图像和超分图像的局部区域。
图 3 超分辨率重建结果图Fig. 3 Super resolution results
图选项
表 2是图 3中各子图的评价参数。
表 2 评价参数Tab. 2 Evaluate parameters of Fig. 3
| 评价参数 | 图像 | band1 | band2 | band3 | band4 | band5 | mean |
| 平均梯度 | original | 8.374 | 4.155 | 7.319 | 11.309 | 7.858 | 7.803 |
| bilinear | 6.266 | 4.178 | 5.806 | 7.367 | 6.011 | 5.926 | |
| SR | 14.303 | 9.008 | 11.260 | 23.004 | 14.355 | 14.386 | |
| 信息熵 | original | 8.827 | 7.920 | 8.821 | 8.737 | 8.679 | 8.597 |
| bilinear | 8.827 | 7.920 | 8.820 | 8.734 | 8.679 | 8.596 | |
| SR | 8.831 | 7.921 | 8.825 | 8.742 | 8.683 | 8.600 | |
| 信噪比 | original | 3.258 | 4.397 | 3.622 | 2.846 | 2.932 | 3.411 |
| bilinear | 3.703 | 4.948 | 4.080 | 3.212 | 3.379 | 3.865 | |
| SR | 3.167 | 4.346 | 3.547 | 2.695 | 2.927 | 3.336 |
表选项
如表 2所示，原始图像，双线性内插图像，超分图像的信息熵没有明显差别。这是因为信息熵是对全图信息量的体现，而原始图像，双线性内插图像，超分图像三者的差异主要体现在图像细节上，而非图像的整体特性上。图 4是三者的灰度直方图分布，可以看到三者的灰度分布几乎是一致的，仅在灰度极值附近有所差异，而这一差异无法在信息熵这一指标上有明显反映。
图 4 灰度直方图比较Fig. 4 Comparison of three gray level histograms
图选项
平均梯度体现了图像的清晰度，就平均梯度而言，双线性内插图像最小，平均为5.926；原始图像居中，平均为7.803；超分图像最大，平均为14.386。这说明经过超分重建后，图像的清晰度得到明显提高，图像细节更为丰富。从图 3(b)、(c)的对比中可以明显看出超分图像具有更为清晰的地物特征。
超分辨率重建的效果与序列影像帧数有直接关系，图 5是不同帧数的超分辨率重建结果。
图 5 不同帧数的影像超分重建结果Fig. 5 Super resolution results processed by different input frame amount
图选项
表 3是图 5中各幅子图的评价参数，图 6是序列影像帧数与影像超分重建结果平均梯度的关系图，图 7是序列影像帧数与影像超分重建结果信噪比的关系图。
表 3 不同帧数影像超分重建结果的评价参数Tab. 3 Evaluate parameters of Fig. 4
| 评价参数 | 帧数 | band1 | band2 | band3 | band4 | band5 | mean |
| 平均梯度 | 2 | 7.156 | 4.789 | 6.466 | 10.300 | 10.296 | 7.801 |
| 3 | 9.428 | 6.504 | 8.412 | 14.084 | 14.092 | 10.504 | |
| 4 | 11.511 | 7.717 | 9.543 | 18.292 | 18.527 | 13.118 | |
| 5 | 14.303 | 9.008 | 11.260 | 23.004 | 14.354 | 14.386 | |
| 信息熵 | 2 | 8.828 | 7.922 | 8.822 | 8.738 | 8.733 | 8.609 |
| 3 | 8.829 | 7.921 | 8.823 | 8.739 | 8.739 | 8.610 | |
| 4 | 8.830 | 7.921 | 8.824 | 8.740 | 8.740 | 8.611 | |
| 5 | 8.831 | 7.921 | 8.825 | 8.742 | 8.683 | 8.600 | |
| 信噪比 | 2 | 3.487 | 4.458 | 3.849 | 3.017 | 3.017 | 3.566 |
| 3 | 3.369 | 4.445 | 3.734 | 2.901 | 2.901 | 3.470 | |
| 4 | 3.269 | 4.462 | 3.641 | 2.798 | 2.798 | 3.394 | |
| 5 | 3.167 | 4.346 | 3.547 | 2.695 | 2.927 | 3.336 |
表选项
图 6 帧数-平均梯度关系图Fig. 6 Relationship between frame amount and average gradient
图选项
图 7 帧数-信噪比关系图Fig. 7 Relationship between frame amount and SNR
图选项
3 结论
【引文格式】许丽娜，何鲁晓。基于凸集投影的高分四号卫星影像超分辨率重建[J]. 测绘学 ，2017，46(8)：1026-1033. DOI:
10.11947/j.AGCS.2017.20170070
声明：本站部分文章内容及图片转载于互联 、内容不代表本站观点，如有内容涉及侵权，请您立即联系本站处理，非常感谢！