《1　前言》

1　前言

遥感图像具有丰富的地貌信息，而在现有硬件条件下很难达到遥感图像高精度的要求，所以常常采用合成多帧图像的方法来重建高分辨率图像，即超分辨率重建[1~3]。但超分辨率方法往往计算量很大，需要融合多帧（至少 4 帧以上图像），不适于星上处理。笔者运用亚像元技术[4]，通过对采样点在水平和垂直方向分别错半个像元的两帧图像进行插值、去噪处理来提高图像的分辨率。其采样方式是模拟两排在空间位置错半个像元的 CCD 探测器。

图 1 表示了相差半个像元的两幅图像互相交错的空间结构示意图。其中 p 为任意像元，p/2 为采样间隔，这就提高了采样密度，从而使得图像的分辨率在理论上提高了 2 倍。如果两幅图像大小均为 m ×n ，那么对它们进行高分辨率重建的空间结构示意图如图 2 所示。

《图1》

图1　错半个像元的两帧图像的空间结构图

Fig.1　The spatial structure of two images offsetting half a pixel

《图2》

图2　重建后图像的结构示意图

Fig.2　The structure of reconstructed image

如图 2 所示， 到 之间的距离可看作单位像元，则 到 0 之间是半个像元的距离。0 点就是要估计的像元点，也就是待插值的像素点。文章利用二分树复小波插值算法，通过对两帧错半个像元的图像进行梅花点阵的交错采样方式，最终形成一幅更高分辨率的遥感图像。

《2　常用的插值方法》

2　常用的插值方法

1）最近邻法[5]。输出像素的值指定为当前点所属像素的值，而不考虑其他像素；

2）双线性法。输出像素的值是最近的 2 × 2 邻域内像素的加权平均值，利用周围4个点确定一个平面，并且在一个矩形栅格上进行一阶插值，它需要用到双线性函数；

3）双三次法。输出像素的值是最近邻的 4 × 4 邻域内像素值的加权平均值。

这 3 种插值方法的基本原理如图 3 所示，“0”是插值点，它由周围的 4 个节点，按照加权平均值的方式，进行估算。上述常用的插值方法可以快速、方便地对图像进行插值。但是，由于常用的插值方法只是对图像的像素灰度值直接进行处理，所以仅仅考虑到邻近点像素的关联，并未考虑到图像中高频细节的完整特征。而小波插值的方法利用多分辨率理论，可以完整地估计出插值点的信息。

《图3》

图3　单个插值点的邻域图

Fig.3　The adjacent region structure of interpolated point

《3　小波变换》

3　小波变换

利用小波的多分辨率的性质，把图像的高、低频分开。利用公式[6]：

式中，   为图像灰度； 为尺度函数； 为尺度系数； 为小波函数； 为小波系数。这样， 就可以表示某尺度 J  分解后的图像低频部 分，而  则可以分别表示水平、竖直和对角线方向上的高频细节之和。但是针对遥感图像的地貌特征复杂、信息量丰富的特点，必须采用一种能够充分反映地貌各个方向特征信息的一种插值方法。常用的离散小波变换（discrete wavelet translation，DWT）不仅缺少平移不变性，而且只能根据分解后有限的3 个方向（水平、垂直和对角线方向）特征信息，很难重构出真实地貌多方向性的遥感图像。针对这一局限，N. G. Kingsbury 提出了二分树复数小波[7]的概念，其一维结构如图 4 所示。

《图4》

图4　一维二分树复小波结构图

Fig.4　The complex wavelet structure of one dimension

图 4 中两个分树 A 和 B 分别对应复小波的实部和虚部，两棵树信息互补，这使得复小波具有平移不变性。拓展到二维后，4 个分树组合得到 6 个复数分量对应 6个方向的高频信息。因此，二分树复小波具有近似位移不变性、多种方向可选性、完美的重建等特点[8]。

《4　亚像元图像重建算法流程》

4　亚像元图像重建算法流程

如图 5 所示，对亚像元图像以 X  进行复小波分解，分别经过高频滤波器 H  和低频滤波器 L 分解后分为高频信号 D  和低频信号 A。用复小波进行插值估计，估计小波系数。由于原始图像上本身有噪声，插值又加入了另一部分噪声，因此应该在插值后对图像去噪。采用基于层间的小波去噪方法[9]，对高频部分进行去噪处理。对插值、去噪声后的高频信号 H ′，低频 L′ 进行重构，重构后的图像为 X ′。本文算法试图把遥感图像的插值和去噪过程整合成一个系统，利用复小波进行一次分解、重构，这样既可以减轻工作量，又可以减少插值引入的噪声。

《图5》

图5　重建算法过程

Fig.5　Process of reconstruction algorithm

《5　二分树复小波插值》

5　二分树复小波插值

二维双树复小波变换 DT -CWT 变换的复尺度函数和各个方向上的小波函数[10]为

和 ，d ∈ B  分别为尺度系数和小波系数； k 和 分别为尺度或小波函数在水平和竖直方向上的位移， 为二维 DT – CWT 的尺度函数矩阵。

分别为各方向上的二维 DT – CWT 小波矩阵。

为二维 DT – CWT 的尺度函数矩阵，， ，，， 和 分别为各方向上的二维 DT -CWT 小波矩阵。 和各  为尺度和小波系数矢量， d  表示各个方向。对低分辨率图像序列的二维 DT -CWT 分解后，再将各个方向的图像叠加起来，得到关于未知系数 和各 的方程组

其中

式中， ｜SJ ｜和 ｜Sj ｜分别为有效的尺度系数和小波系数个数。

对尺度函数系数进行估计：在小波变换中，绝对值大的小波系数对应着图像中突出特征，主要是纹理部分、区域的边界或噪声，而分解后的 6 个方向的细节子图分别反映了各方向的特征。对系数估计，就是用欠采样图像中已知的像素点来估测小波系数，最后由估计出的系数进行图像重建，即完成插值过程。首先利用最小二乘法估计 （代表尺度系数矢量 ）， 由近似式, = 

即

 

式中， （ > 0）为平衡因子。若 过大，会使得到的解不够理想；若 过小，则会引出过多噪声。从而有估计 ，再利用 估计小波系数 bJ 。在式（4）的 DT – CWT 分解中，拟合部分为低分辨率图像序列的近似，因此有

对小波函数系数估计：J 尺度上 15° 方向细节  近似等于公式（6）估计系数后剩下的残差

因此由式（7）可估计小波系数矢量 。

进一步，由 J 尺度上 45°方向的细节 与前面系数估计后剩下的残差的近似相等，即由式可估计 。

则其他尺度上的小波系数也可得到同样的估计。尺度系数矢量和各小波系数矢量的估计结果均为复值。

再用式（9）重建高分辨率图像高频在 6 个方向上的细节。

到此，完成了把图像按照 DT – CWT 的分解，这样就把图像分成一个低频部分和 6 个方向上的高频部分。

《6　小波的基于层间的自适应去噪》

6　小波的基于层间的自适应去噪

由于在小波插值过程中，噪声很可能成为比较大的小波系数，被误认为是有用信息估计到图像中，再加上估计小波系数的过程中可能引入无用信息，因而插值会带来噪声，所以有必要进行去噪处理。虽然各层小波系数间是独立的，但其模是相关的。把层间小波系数的相关性考虑到去噪中，可以避免“过扼杀”现象，能够保留图像中完整的细节信息。考虑到这一特性，我们利用一种自适应的双变量紧缩模型[11]对分解后的小波系数进行去噪。

去噪算法的流程如下：

1）计算噪声方差的估计 。已知图像经过小波分解后各方向上的小波系数 = ，，这里 和 σ 2 分别为噪声方差和局部方差。

2）对每一个小波系数（k = 1，…，m）：

a. 计算 。从观测模型得到 ，则  的经验估计为

b. 用式（12）计算 σ 。M  是矩形窗区域 N（k）的大小。由此可得 σ 的估计：

其中，代表

 

c. 用式（13）计算每一个原始图像的小波系数的估计 Bj 。而 代表 的父层的小波系数（即同意空间位置的邻近尺度层的小波系数）。则有 Bj 代表 j  尺度下估计的小波系数，见式（13）：

这种层间的对应关系如图 6 所示。

《图6》

图6　空间上对应同一点在不同尺度的层间关系

Fig.6　The relationship between scale and its father scale corresponding same point

由图 6 可以看出，这种层间相关的去噪方法，不仅考虑到了层内小波系数之间的影响，而且充分考虑了层间小波系数的相关性。

最后，对融合后的低频和 6 个方向上的高频部分进行重构，由式（1） 重建内插后的高分辨率图像。

《7　实验结果的评价标准》

7　实验结果的评价标准

笔者采用图像复原中广泛使用的峰值信噪比（PSNR）作为评价准则，主要用于对原始图像与重建图像进行比较，其公式如下：

式中，   为原始图像； 为内插后的图像；M，N 分别为图像的长度和宽度。由于 PSNR 能够客观、稳定地反映图像的质量特征，因而常被采用。

《8　实验结果与对比分析》

8　实验结果与对比分析

本文分别对两帧 128 ×128 错半个像元的靶标图像和两帧 256 ×256 的遥感图像进行了实验，硬件配置为 Intel PIV 2.94 GHz，内存为 512 Mb 的 PC。此外，还实验了几种常用的插值算法以及 DWT 插值和 AdapShrink[12,13]去噪方法，并分别计算其PSNR 用于对比分析。

如图 7 所示，对靶标图像进行二分树复小波分解的 6 个方向的高频细节。

《图7》

图7　二分树复小波分解的高频子图

Fig.7　The high frequency subimages of DT -CWT decomposition

实验中发现，采用二分树复小波对靶标进行 4 层分解后，效果较好。对遥感图像进行 6 层分解，重构结构较好。超过 6 层分解后，重构发现图像的细节过分突出，边缘保持不好，甚至发生扭曲。而这种对更多方向的可选择性（优于 DWT 仅仅对 3 个方向的分解），满足了对复杂内容图像的重构。各种插值去噪声重构方法的实验图比对结果，如图 8（靶标图）、图 9（ 城区图 -真实的遥感图像）所示。

《图8》

图8　不同方法插值及去噪结果（靶标图）

Fig.8　The result image of different interpolation methods（target image）

《图9》

图9　不同方法插值及去噪结果（城区图）

Fig.9　The result image of different interpolation methods（city image）

表1 是使用几种插值方法得到的峰值信噪比（PSNR）（噪声方差 σn =30 ）的比较结果。

《表1》

表1　各种插值方法的 PSNR 比较。

Table 1　PSNR of various methods

由表 1 可以看出，用小波插值的算法对两帧错半个像元的图像进行重建、去噪，峰值信噪比相对于常用的插值方法都要高，而且比用 DWT 插值结合AdapShrink 去噪的方法高。这主要是因为本文算法考虑到了遥感图像复杂的纹理信息，采用了具有多方向选择性的二分树复小波。结合基于层间的自适应阈值法，有效地去除了原有图像中的噪声和插值引入的噪声。另外，本文提出算法复杂度为 O（nlogn），其中图像大小为 n ×n 像素。

图 10 是退化的遥感图像与经过文章的算法处理后的遥感图像的频谱图。

《图10》

图10　处理前后的频谱对比图

Fig.10　Spectrum figure comparison of degraded image and restored image

通过频谱图的比较可见，经过本文算法处理的遥感图像的主要的高频信息被有效地集中到了几个区域，这使得噪声得以抑制和减少，因而图像变得更清晰。

《9　结语》

9　结语

提出了一种在亚像元的基础上运用二分树复小波插值和基于层间的小波去噪方法。该方法在插值的过程中能够充分地反映遥感图像中复杂的地貌特征，结合基于层间的自适应阈值法，有效地去除了原有图像中的噪声和插值引入的噪声。实验结果能达到对图像进行平滑、去模糊、保留阶梯函数的不连续性的目的。因此，这种小波插值方法与传统方法相比，能更好地提高遥感图像质量。