0 引言

分辨率是遥感图像的最重要特征之一。由于许多实际条件的限制，获取的遥感图像分辨率不能达到实际需求。图像超分辨重建技术可以对获取的一帧或多帧低分辨率图像进行处理以得到更高分辨率的图像。

图像超分辨重建大致可分为三类：基于重建的方法、基于插值的方法和基于学习的方法。基于重建的方法^[1-3]一般用于多帧超分辨，多帧低分辨率图像间的运动估计是该方法的难点，运动估计的误差会对超分辨效果产生比较大的影响。基于插值的方法^[4-5]可以快速地对单帧图像进行超分辨，但合适的先验信息不易加入，放大倍数较大时性能下降。基于学习的方法近年来获得了巨大的进展，大致又可分为基于神经网络^[6-8]和基于稀疏表示与压缩感知的超分辨^[9-17]两大类:前者利用神经网络学习高和低分辨率图像之间的映射关系，该类方法在学习阶段计算量较大; 后者利用稀疏性作为先验信息，为求解欠定问题提供了新的视角和工具。

文献[9-13]将压缩感知和稀疏表示理论应用于图像的超分辨重建，获得了良好的超分辨效果。文献[9-12]利用K-奇异值分解 (K-Singular Value Decomposition, K-SVD) 方法进行字典学习，基于压缩感知理论进行超分辨重建。由于训练图像与超分辨图像纹理有可能存在较大差别，从而影响超分辨重建的效果。在此基础上，文献[13]将训练图像限定在超分辨图像本身，一定程度上解决了上述问题，但需要对每帧超分辨图像都进行字典学习，即字典是一次性的。

上述都是基于单字典学习。一般认为，利用单一字典表达复杂多样的纹理有一定的局限性。基于多字典学习获得了发展^[14-16]。文献[14]按照图像块梯度和的大小将图像块纹理分成7类，建立7个不同字典，重建阶段将图像块归类并在相应的类别中进行稀疏求解，但该分类方法没有考虑图像纹理的方向性。文献[15]对每组图像块提取K个纹理，形成K个子字典；重建阶段对每个低分辨率图像块进行K次稀疏求解，K值较大时计算量很大。考虑到小波变换可以将图像分解为低频频带和3个方向的高频频带，并且具有多分辨特性。利用小波变换对图像纹理进行分类是一种有效的方法。文献[16]在小波域低频部分进行稀疏求解，另3个频率采用双立方插值，显然高频成分的恢复不足。

在上述基础上，本文利用小波变换将复杂多样的纹理分解为4个不同的频带，对应建立4个字典来描述。基于压缩感知理论在4个不同频段分别进行稀疏求解，最后增加全局性限制以消除块状效应。

1 理论基础 1.1 压缩感知

假定信号x∈R^N，y=Φx为对该信号的一次测量，其中Φ∈R^M×N为投影算子。令x=Ψα，其中Ψ∈R^N×N为信号基。若‖α‖₀=S＜N，则称信号x相对于基Ψ是S稀疏的。定义μ(Φ, Ψ)=$\sqrt N \cdot \mathop {\max }\limits_{1 \le k,j \le N} $〈Φ_k, ψ_j〉为Φ与Ψ的相干性度量。

压缩感知理论^[17-19]认为，在满足一定条件下，可以通过较少的测量信号y完美重建原始信号x。完美重建所需的测量数目与Φ和Ψ的相干性、x的稀疏性密切相关。

1.2 字典学习

字典学习指冗余字典的构建过程。假定I_b=Dα_b，其中I_b为图像块矢量，D∈R^N×K为冗余字典，α_b∈R^K为稀疏表示系数，且有‖α_b‖₀=r$ \ll $K。K-SVD字典学习方法^[20]借鉴K均值聚类的思想，通过求解最优化问题

$\mathop {\min }\limits_{\mathit{\boldsymbol{D}},{\mathit{\boldsymbol{\alpha }}_b}} \left\{ {\left\| {{\mathit{\boldsymbol{I}}_b} - \mathit{\boldsymbol{D}}{\mathit{\boldsymbol{\alpha }}_b}} \right\|_2^2} \right\},满足{\left\| {{\mathit{\boldsymbol{\alpha }}_b}} \right\|_0} \le {T_0}$

(1)

同时获得冗余字典D和对应的稀疏表示α_b，其中T₀为稀疏限制的阈值。

1.3 稀疏编码

稀疏编码指稀疏表示系数的求解过程，即求解式 (2) :

$\mathop {\min }\limits_{{{\bf{ \pmb{\mathsf{ α}} }}_b}} {\left\| {{\mathit{\boldsymbol{\alpha }}_b}} \right\|_0},满足\left\| {{\mathit{\boldsymbol{I}}_b} - \mathit{\boldsymbol{D}}{\mathit{\boldsymbol{\alpha }}_b}} \right\|_2^2 \le \varepsilon $

(2)

近年来涌现出多种稀疏编码算法，如基追踪 (Basis Pursuit, BP) 算法、匹配追踪 (Matching Pursuits, MP) 算法、正交匹配追踪 (Orthogonal Matching Pursuit, OMP) 算法等^[21]，但BP算法的计算量较大，而OMP算法是对MP算法的改进，在分解的每一步对所选择的全部原子进行正交化处理，收敛速度更快^[22]。

2 小波域多字典超分辨重建

为克服单字典表达复杂纹理图像能力的不足，利用多个字典进行表示将会获得更好的效果。小波变换可以将图像分解成4个不同的频带，分别对应低频、垂直、水平和倾斜频带，并且可以进行多级的分解。本文将利用多级小波变换，在小波域利用压缩感知理论求解面向遥感图像的超分辨重建问题。

假设X_b、Y_b分别表示对应的第b个高分辨率和低分辨率图像块，H表示退化算子 (图像退化过程包含模糊和下采样)，依据压缩感知理论，X_b又可以利用字典D进行稀疏表示，则两者关系可表示为式 (3)：

${\mathit{\boldsymbol{Y}}_b} = {\mathit{\boldsymbol{H}}}{{\mathit{\boldsymbol{X}}}_b} = {\mathit{\boldsymbol{HD}}}{{\mathit{\boldsymbol{\alpha }}}_b}$

(3)

式 (3) 可通过式 (4) 所示的最优化问题求解：

$\mathop {\min }\limits_{{{\mathit{\boldsymbol{\alpha }}}_b}} {\left\| {{{\mathit{\boldsymbol{\alpha }}}_b}} \right\|_0},满足\left\| {{{\mathit{\boldsymbol{Y}}}_b} - {\mathit{\boldsymbol{HD}}}{{\mathit{\boldsymbol{\alpha }}}_b}} \right\|_2^2 \le \varepsilon $

(4)

进一步在小波域通过多字典求解，设W_S和W_A分别表示小波综合与分析算子，对应小波域4个频带分别训练获得4个字典D_k(k=H(表示水平方向), V (表示垂直方向), D (表示对角线方向), A (表示低频部分分解系数))，α_{b, k}表示第b图像块第k频带的稀疏表示，则最优化问题变为：

$\mathop {\min }\limits_{{{\mathit{\boldsymbol{\alpha }}}_{wb}}} {\left\| {{{\mathit{\boldsymbol{\alpha }}}_{wb}}} \right\|_0},满足\left\| {{{\mathit{\boldsymbol{Y}}}_b} - {\mathit{\boldsymbol{H}}}{{\mathit{\boldsymbol{W}}}_S}{{\mathit{\boldsymbol{D}}}_w}{{\mathit{\boldsymbol{\alpha }}}_{wb}}} \right\|_2^2 \le \varepsilon $

(5)

其中：D_w=[D_A, D_H, D_V, D_D]，α_wb=α_{b, A}，α_{b, H}，α_{b, V}，α_{b, D}]^T。

上述最优化问题可通过式 (6) 求解：

${{\mathit{\boldsymbol{\hat \alpha }}}_{wb}} = \mathop {\arg \min }\limits_{{\mathit{\boldsymbol{\alpha }}_{wb}}} \left\| {{\mathit{\boldsymbol{Y}}_b} - \mathit{\boldsymbol{H}}{\mathit{\boldsymbol{W}}_S}{\mathit{\boldsymbol{D}}_w}{\mathit{\boldsymbol{\alpha }}_{wb}}} \right\|_2^2 + {\mu _b}{\left\| {{\mathit{\boldsymbol{\alpha }}_{wb}}} \right\|_0}$

(6)

由于H和W_S两个算子不具有可交换性，为简化求解，将其分成如下两步 (式 (7) 和 (8))：

$\mathit{\boldsymbol{\hat X}} = \mathop {\arg \min }\limits_\mathit{\boldsymbol{X}} \left\{ {\left\| {\mathit{\boldsymbol{Y}} - \mathit{\boldsymbol{HX}}} \right\|_2^2} \right\}$

(7)

${{\mathit{\boldsymbol{\hat \alpha }}}_{wb}} = \mathop {\arg \min }\limits_{{{\bf{ \pmb{\mathsf{ α}} }}_{wb}}} \left\| {{{\mathit{\boldsymbol{\hat X}}}_b} - {\mathit{\boldsymbol{W}}_S}{\mathit{\boldsymbol{D}}_w}{\mathit{\boldsymbol{\alpha }}_{wb}}} \right\|_2^2 + {\mu _b}{\left\| {{\mathit{\boldsymbol{\alpha }}_{wb}}} \right\|_0}$

(8)

式 (7) 的求解需要加入适当的正则化项，这里采用Laplace正则化项。式 (8) 等价于：

${{\mathit{\boldsymbol{\hat \alpha }}}_{wb}} = \arg \mathop {\min }\limits_{{\mathit{\boldsymbol{\alpha }}_{wb}}} \left\| {{\mathit{\boldsymbol{W}}_A}{{\mathit{\boldsymbol{\hat X}}}_b} - {\mathit{\boldsymbol{D}}_w}{\mathit{\boldsymbol{\alpha }}_{wb}}} \right\|_2^2 + {\mu _b}{\left\| {{\mathit{\boldsymbol{\alpha }}_{wb}}} \right\|_0}$

(9)

根据以上分析，本文建立基于小波域多字典学习和压缩感知的超分辨重建算法，步骤如下：

步骤1 将训练图像变换到小波域，对4个分量都分成大小$\sqrt N \times \sqrt N $的块，利用K-SVD算法分别进行训练，获得4个字典D_A、D_H、D_V、D_D。

步骤2 求解${\mathit{\boldsymbol{\hat X}}} = \arg \mathop {\min }\limits_{\mathit{\boldsymbol{X}}} \left\{ {\left\| {{\mathit{\boldsymbol{Y}}} - {\mathit{\boldsymbol{HX}}}} \right\|_2^2 + \gamma {{\left\| {\Gamma {\mathit{\boldsymbol{X}}}} \right\|}_2}} \right\}$，式中，Γ为Laplace高通算子，获得全局限制解。

步骤3 利用OMP算法，求解$\mathop {\min }\limits_{{{\mathit{\boldsymbol{\alpha }}}_{k,b}}} {\left\| {{{\mathit{\boldsymbol{\alpha }}}_{k,b}}} \right\|_0}$，满足$\left\| {{{\left[ {{{\mathit{\boldsymbol{W}}}_a}{\mathit{\boldsymbol{\hat X}}}} \right]}_{k,b}} - {{\mathit{\boldsymbol{D}}}_k}{{\mathit{\boldsymbol{\alpha }}}_{k,b}}} \right\|_2^2 \le \varepsilon $，其中：k=H、V、D、A；b表示第b图像块。

步骤4 求解${{\mathit{\boldsymbol{\hat X}}}_b} = {{\mathit{\boldsymbol{W}}}_s}{{\mathit{\boldsymbol{D}}}_w}{{\mathit{\boldsymbol{\alpha }}}_{wb}}$，式中D_w和α_wb的含义同公式 (5)，生成高分辨率图像块${{\mathit{\boldsymbol{\hat X}}}_b}$。

步骤5 将各图像块${{\mathit{\boldsymbol{\hat X}}}_b}$拼接后进行逆小波变换形成一帧高分辨率图像。

3 实验结果 3.1 对比实验

为证明算法的有效性，本文进行了下面实验，实验所用图像来自USC-SIP遥感图像库。在字典训练中，对如图 1所示的不同纹理的6帧尺寸为1024×1024的遥感图像进行二级小波分解，对分解后的低频部分和3个高频部分进行分块，分块大小为2×2，然后对这4组图像块利用K-SVD算法进行字典学习，得到4个小波域字典。每个字典的原子数目本文参照文献[12]的结果，设为16，则每个字典是大小为4×16的矩阵。

利用训练好的字典，对其他遥感图像进行超分辨重建。文献[13]提出了基于压缩感知的单字典算法，与本文算法原理接近。文献[14]给出了一种基于稀疏表示和多字典对的算法。本文算法主要与这两种算法进行对比分析。对图 2所示的4帧图像进行实验，得到4组实验结果。首先对一帧尺寸为256×256的高分辨率图像 (High Resolution Image, HRI) 执行支持域为3×3、标准差为2的高斯模糊，接着1/2×1/2倍的下采样，得到一帧退化的低分辨率图像 (Low Resolution Image, LRI)，对LRI分别执行双线性插值、文献[13]算法、文献[14]算法和本文算法进行比较分析。

图 3给出了第一组实验结果, 图 4给出了该实验结果图像中截取的两块小区域的放大图像。

从图 4中可以看出，相比LRI的双线性插值，文献[13]和[14]算法的超分辨图像的分辨率有了很大提高，图像更加清晰，原来不能分辨的许多纹理已经能够分辨出来。文献[14]算法利用多字典对在一定程度上比单字典对表达图像纹理更有效，性能上比文献[13]算法有一定提高, 而本文算法利用小波多级分解进行纹理分类，能更有效地表达复杂的纹理，结果图像的分辨率有了较大提高，一些在图 4(c)和图 4(d)中模糊或者不能分辨的纹理在图 4(e)中能够明显分辨出来。

为定量比较算法的性能，表 1和表 2给出了上述4组实验结果的峰值信噪比 (Peak Signal-to-Noise Ratio, PSNR) 和结构相似度 (Structural SIMilarity, SSIM) 指标。相比双线性插值，文献[13]算法结果的PSNR指标和SSIM指标分别提高3 dB和0.04以上；文献[14]算法结果又相对于文献[13]算法结果分别提高0.4 dB和0.003以上；本文算法比双线性插值提高6 dB和0.057以上。

3.2 计算量分析

从算法原理上，文献[13]算法与本文算法最为接近，本文主要与该算法进行计算量的比较。根据算法原理，与文献[13]算法相比，本文算法的计算主要增加了全局求解、小波正反变换和多字典的建立三部分，如果将字典建立包含在内，计算量有适当的增加; 但本文算法所建字典可重复使用，减少了每次都重建字典的计算, 当字典一旦建立重复使用时，计算量又会大大减小。

为便于计算量的比较，对本文算法包含字典建立和不包含字典建立两种计算时间进行实验比较，而文献[13]算法的计算时间则必须包含字典建立部分。分别对128×128、256×256和512×512三种尺寸的10组遥感图像执行上述三种算法的计算时间取平均值进行比较，结果如图 5所示。可以看出，如果考虑字典建立时间，本文算法不到文献[13]算法的4倍，主要为建立字典的时间。如果不考虑字典建立时间，本文算法大约为文献[13]算法的计算时间的一半。

4 结语

针对遥感图像的分辨率增强以及单字典表达复杂图像纹理的局限性问题，本文提出了一种基于压缩感知理论的小波域多字典超分辨重建算法。基本思想是基于压缩感知理论，将单字典改为小波域多字典进行超分辨重建。利用小波变换的频带分解和多分辨特性，在小波域进行多字典的训练和稀疏求解。结果表明，本文算法重建图像的分辨率相比单字典算法和同类多字典算法有很大提高。主观评价上超分辨重建图像纹理更加丰富清晰, 客观评价指标PSNR和SSIM指标都有较大提高。考虑字典建立的时间，与单字典算法相比，本文算法计算量稍大。由于字典可重复使用，字典建立以后进行大批量遥感图像超分辨处理时，计算时间又大幅减少，为单字典算法的一半左右。对于小波域字典与感知矩阵的相干性分析、利用更有效的图像表达方式如Contourlet等进行多字典多尺度的超分辨重建、同时考虑相邻图像块间的相关性等是下一步的研究方向。