CN104252703A - 基于小波预处理和稀疏表示的卫星遥感图像超分辨率重构方法 - Google Patents

Abstract

本发明的基于小波预处理和稀疏表示的卫星遥感图像超分辨率重构方法属于卫星遥感图像处理的技术领域。该方法的应用条件是已知相同观测地区、不同时间分辨率的高分辨率遥感图像和低分辨率遥感图像，对其它观测时间的低分辨率遥感图像进行超分辨率重构，提高低分辨率遥感图像的空间分辨率。具体分为字典训练和对低分辨率遥感图像进行重构两个步骤。本发明考虑到遥感图像的物候变化，构建了含不同特征信息的小波域的字典，基于三对子波段字典的训练，结合稀疏表示，有效地实现低分辨率遥感图像的超分辨率重构，更好地获取了图像的细节特征；本发明的方法可有效提高低分辨率遥感图像的重构质量，为低分辨率遥感图像的后期应用提供基础。

Description

基于小波预处理和稀疏表示的卫星遥感图像超分辨率重构方法

技术领域

本发明属于卫星遥感图像处理的技术领域，利用小波预处理技术，提出一种基于稀疏表示的卫星遥感图像超分辨率重构方法。该方法能够提高原始遥感图像的空间分辨率，提供比原始图像更多的细节信息，有助于遥感图像的后续处理。

背景技术

遥感影像的空间分辨率是影响目标识别和精确判读的一个重要因素。为了提高低分辨率遥感图像的空间分辨率，可以采用超分辨率重构方法进行图像增强。超分辨率重构方法是一种指利用同一场景的一帧或多帧低分辨率图像来重构一帧或多帧高分辨率图像的算法。

近年来，基于稀疏表示的图像超分辨率重构引起了很多学者的关注。目前，国内外的图像超分辨率重构算法主要是针对常规图像进行重构，具有代表性的方法有：Yang等首先提出基于稀疏表示的超分辨率算法，先求解低分辨率图像块在对应字典下的稀疏表示，再通过线性组合的方法利用稀疏表示系数获得相应的高分辨率图像块(参见Yang J,Wright J,Huang T S.et al.,“Image super-resolutionvia sparse representation,”IEEE Transactions on Image Processing,19(11),2861-2873(2010)；Yang J,Wright J,Huang T,et al.,“Image super-resolution assparse representation of raw image patches,”IEEE Conference on Computer Visionand Pattern Recognition,1-8(2008))。Zeyde和Elad等对Yang的算法提出改进，在图像重构过程中利用正交匹配追踪(OMP)算法求解稀疏表示系数，减少了算法的计算量(参见Zeyde R,Elad M,Protter M.,On single image scale-up usingsparse-representations[M],Curves and Surfaces.Springer Berlin Heidelberg,711-730(2012)；Elad M.,Sparse and Redundant Representations:Chapter 15.4.2TheSuper-Resolution Algorithm[M],Springer Press,346-357(2010))。2014年，在Zeyde的算法的基础上，Nazzal采用小波域的字典训练算法，提高了重构质量和计算效率(参见Nazzal M,Ozkaramanli H.,“Wavelet domain dictionary learning-basedsingle image superresolution,”Signal,Image and Video Processing,1-11(2014))。对于遥感图像的超分辨率重构算法的研究有：2013年，赵永光等采用基于稀疏表示的超分辨率重构算法来提高MODIS图像的空间分辨率，但该算法只能在一定程度上提高原MODIS图像的分辨率(参见赵永光,黄波,汪超亮.基于超分辨率重建的多时相MODIS与Landsat反射率融合方法[J].遥感学报,2013,17(3):590-608)。遥感图像序列中长时间的间隔导致了序列图像之间存在变化。宋慧慧分析了遥感图像之间的物候变化，认为该变化主要表现在低频信息的变化(参见宋慧慧.基于稀疏表示的图像超分辨率重建算法研究[D].合肥:中国科学技术大学,2011)。

综上所述，国内外学者利用基于稀疏表示的超分辨率重构方法主要是针对常规图像进行处理，取得了相关的研究成果。利用稀疏表示对于遥感影像超分辨率重构的研究属于初步的探讨，对于遥感图像空间分辨率提高的能力有限。

发明内容

本发明要解决的技术问题是，提供一种有效的基于小波预处理技术和稀疏表示的卫星遥感影像超分辨率重构的算法。该方法的应用条件是已知相同观测地区、不同时间分辨率的高分辨率遥感图像和低分辨率遥感图像，对其它观测时间的低分辨率遥感图像进行超分辨率重构，提高低分辨率遥感图像的空间分辨率。

为解决本发明要解决的技术问题，给出的具体技术方案如下。

一种基于小波预处理和稀疏表示的卫星遥感图像超分辨率重构方法，是已知相同观测地区、不同时间分辨率的高分辨率遥感图像和低分辨率遥感图像，对其它观测时间的低分辨率遥感图像进行超分辨率重构；分为字典训练和对低分辨率遥感图像进行重构两个步骤；

所述的字典训练的过程是，

1)对已知低分辨率遥感图像上采样，使其和已知高分辨率遥感图像具有相同的图像尺寸；

2)对已知高分辨率遥感图像、低分辨率遥感图像分别进行小波分解提取特征，获得低频、水平、和垂直高频三对子波段；

3)在已知高分辨率遥感图像、低分辨率遥感图像的三对子波段图像上提取图像块作为训练样本；

4)利用冗余字典训练算法对已知低分辨率图像块进行训练，得到已知低分辨率图像块对应的字典和稀疏表示系数；

5)计算已知高辨率卫星遥感图像块对应的字典；

所述的对低分辨率遥感图像进行重构的过程是，

1)对待重构的低分辨率遥感图像进行小波分解，提取特征，获得低频、水平、垂直和对角线高频四个子波段；

2)在待重构低分辨率遥感图像的低频、水平、和垂直高频三个子波段图像上提取图像块作为训练样本；

3)求解图像块在已知低分辨率遥感图像块对应的字典中的稀疏表示系数；

4)根据求得的已知高分辨率遥感图像块对应的字典，计算待重构遥感图像在低频、水平、和垂直高频三个子波段上对应的高分辨率遥感图像，结合待重构遥感图像的对角线高频子波段，通过小波逆变换合成高辨率遥感图像重构结果。

在字典训练的过程的第4)步中，所述的冗余字典训练算法，是K-SVD冗余字典训练算法；具体的是通过(1)式对已知低分辨率遥感图像块进行训练，得到三个尺寸大小为n²×m_a的低分辨率遥感图像块对应的字典和稀疏表示系数其中m_a是冗余字典的列数或者原子数；

$D_M^y,\left\{q_y^i\right\}=\arg \min \underset{i}{\mathrm{\Σ}}{\left|\right|M_y^i-D_M^y{q}_y^i\left|\right|}_2,s.t.\underset{i}{\mathrm{\Σ}}{\left|\right|q_y^i\left|\right|}_0\≤k_0,\∀i,y=f,h,v---\left(1\right)$

式(1)中，k₀为稀疏表示系数中非零元素的个数，f，h和v分别表示低频、水平、和垂直高频，i＝1,...,P，表示第i个训练样本，P表示训练样本的个数。

在字典训练的过程的第5)步中，所述的计算已知高辨率卫星遥感图像块对应的字典通过伪逆表达式(2)求得；

$D_L^y=L_y{Q_y}^T{\left(Q_y{Q_y}^T\right)}^{-1},y=f,h,v---\left(2\right)$

式(2)中，L_y是由作为其列组成，Q_y是字典训练的过程的第4)步求得的作为其列组成，其中T代表矩阵的转置。

在对低分辨率遥感图像进行重构过程的第3)步中，所述的求解图像块在已知低分辨率遥感图像块对应的字典中的稀疏表示系数，是利用OMP算法求解图像块在已知低分辨率遥感图像块对应的字典中的稀疏表示系数

在对低分辨率遥感图像进行重构过程的第4)步中，所述的通过小波逆变换合成高辨率遥感图像重构结果，是利用式(3)求得对应的分辨率增强后的图像块

${\hat{T}}_y^i=D_L^y\×q_y^{\′i},y=f,h,v---\left(3\right)$

再对所有图像块进行以上过程然后合并得到最后结合T_d通过小波逆变换合成重构结果

更完整具体的技术方案可以叙述如下：

已知相同观测地区、不同时间分辨率的高分辨率遥感图像L和低分辨率遥感图像M，获得其他观测时间的低分辨率遥感图像T的超分辨率重构结果

本发明的流程图如图1所示，具体包括如下过程：

步骤一：字典训练

1)对已知低分辨率遥感图像M上采样，使其和已知高分辨率遥感图像L具有相同的图像尺寸；

2)对已知高、低分辨率遥感图像分别进行小波分解来提取特征，获得低频、水平、和垂直高频三对子波段(L_f,M_f,L_h,M_h,L_v,M_v)，其中f，h和v分别表示低频，水平和垂直；

3)在已知高、低分辨率遥感图像的子波段图像上提取大小为n×n，同时块之间重叠大小为(n-4)×(n-4)的图像块作为训练样本；其中i＝1,...,P，表示第i个训练样本，P表示训练样本的个数；

4)利用K-SVD冗余字典训练算法，通过式(1)对已知低分辨率遥感图像块进行训练，得到三个尺寸大小为n²×m_a的低分辨率遥感图像块对应的字典和稀疏表示系数其中m_a是冗余字典的列数或者原子数；

$D_M^y,\left\{q_y^i\right\}=\arg \min \underset{i}{\mathrm{\Σ}}{\left|\right|M_y^i-D_M^y{q}_y^i\left|\right|}_2,s.t.\underset{i}{\mathrm{\Σ}}{\left|\right|q_y^i\left|\right|}_0\≤k_0,\∀i,y=f,h,v---\left(1\right)$

式(1)中，k₀为稀疏表示系数中非零元素的个数；

5)已知高分辨率遥感图像块对应的字典可以通过伪逆表达式(式(2))求得；

$D_L^y=L_y{Q_y}^T{\left(Q_y{Q_y}^T\right)}^{-1},y=f,h,v---\left(2\right)$

式(2)中，L_y是由作为其列组成，Q_y是4)求得的作为其列组成，T代表矩阵的转置；

步骤二：对低分辨率遥感图像进行重构

1)对待重构的低分辨率遥感图像T进行小波分解来提取特征，获得低频、水平、垂直和对角线高频四个子波段(T_f,T_h,T_v,T_d)，其中f，h，v,d分别表示低频，水平,垂直和对角线；

2)在待重构低分辨率遥感图像的低频、水平、和垂直高频三个子波段图像上提取大小为n×n，同时块之间重叠大小为(n-4)×(n-4)的图像块作为训练样本；

3)利用OMP算法求解图像块在已知低分辨率遥感图像块对应的字典中的稀疏表示系数

4)利用式(3)求得对应的分辨率增强后的图像块对所有图像块进行以上过程然后合并得到最后结合T_d通过小波逆变换合成重构结果

${\hat{T}}_y^i=D_L^y\×q_y^{\′i},y=f,h,v---\left(3\right).$

本发明的有益效果：

考虑到遥感图像的物候变化，本发明构建了三对包含不同特征信息(低频、水平和垂直高频子波段)的小波域的字典，基于三对子波段字典的训练，结合稀疏表示，有效地实现低分辨率遥感图像的超分辨率重构，更好地获取了图像的细节特征。该算法的研究结果可以有效提高低分辨率遥感图像的重构质量，为低分辨率遥感图像的后期应用提供基础。

附图说明

图1是本发明基于小波预处理和稀疏表示的卫星遥感图像超分辨率重构方法的流程图。

图2是本发明实施例1的2001年5月24日MODIS遥感图像。

图3是本发明实施例1的2001年5月24日Landsat ETM+遥感图像。

图4是本发明实施例1的2001年8月12日MODIS遥感图像。

图5是本发明实施例1的2001年8月12日Landsat ETM+遥感图像。

图6是本发明实施例1的2001年7月11日MODIS遥感图像。

图7是本发明实施例1的2001年7月11日Landsat ETM+遥感图像。

图8是本发明实施例1的2001年7月11日MODIS超分辨率重构遥感图像。

具体实施方式

实施例1：采用30m高空间分辨率的Landsat ETM+和250m低空间分辨率的MODIS遥感数据作为测试数据。MODIS遥感数据的时间分辨率1天，LandsatETM+遥感数据的时间分辨率是16天。两类数据来源于NASA’s Reverb ECHO的网站(http://reverb.echo.nasa.gov/reverb/)数据的采集时间分别为：2001年5月24日，2001年7月11日和2001年8月12日。观测地区位于(54°N,104°W)附近。本发明从2001年5月24日和2001年8月12日观测的MODIS和LandsatETM+遥感图像(图2-图5)中选择高、低分辨率训练样本对，采用本发明提出的方法对2001年7月11日的MODIS遥感图像(图6)进行超分辨率重构，并用2001年7月11日的Landsat ETM+遥感图像(图7)为作为目标图像进行验证。具体包括以下步骤：

步骤一：字典训练

1)对已知低分辨率MODIS图像M通过ENVI软件进行Bicubic双线性插值上采样，使其和已知高分辨率Landsat ETM+图像L具有相同的图像尺寸，处理后图像样本的尺寸均为400×400像素；

2)对已知高分辨率Landsat ETM+和低分辨率MODIS图像分别进行小波分解来提取特征，小波基选择db1小波，获得低频、水平、和垂直高频三对子波段(L_f,M_f,L_h,M_h,L_v,M_v)，其中f，h和v分别表示低频，水平和垂直；

3)在已知高分辨率Landsat ETM+和低分辨率MODIS图像的子波段图像上提取大小为22×22像素，同时块之间重叠大小为18×18像素的图像块作为训练样本；其中i＝1,...,2025，表示第i个训练样本；

4)利用K-SVD冗余字典训练算法，K-SVD算法进行1次迭代；通过式(1)对已知低分辨率MODIS图像块进行训练，得到三个尺寸大小为484×1800的低分辨率MODIS图像块对应的字典和稀疏表示系数

$D_M^y,\left\{q_y^i\right\}=\arg \min \underset{i}{\mathrm{\Σ}}{\left|\right|M_y^i-D_M^y{q}_y^i\left|\right|}_2,s.t.\underset{i}{\mathrm{\Σ}}{\left|\right|q_y^i\left|\right|}_0\≤k_0,\∀i,y=f,h,v---\left(1\right)$

式(1)中，k₀为稀疏表示系数中非零元素的个数，k₀＝6；

5)已知高分辨率Landsat ETM+图像块对应的字典可以通过伪逆表达式(式(2))求得；

$D_L^y=L_y{Q_y}^T{\left(Q_y{Q_y}^T\right)}^{-1},y=f,h,v---\left(2\right)$

式(2)中，L_y是由作为其列组成，Q_y是4)求得的作为其列组成；

步骤二：对低分辨率MODIS图像进行重构

1)对待重构的低分辨率MODIS图像T进行小波分解来提取特征，获得低频、水平、垂直和对角线高频四个子波段(T_f,T_h,T_v,T_d)，其中f，h，v，d分别表示低频，水平，垂直和对角线；

2)在待重构低分辨率的MODIS图像的低频、水平、和垂直高频三个子波段图像上提取大小为22×22像素，同时块之间重叠大小为18×18像素的图像块作为训练样本；

3)利用OMP算法求解图像块在已知低分辨率MODIS图像块对应的字典中的稀疏表示系数

4)利用式(3)求得对应的分辨率增强后的图像块对所有图像块进行以上过程然后合并得到最后结合T_d通过小波逆变换合成重构结果

${\hat{T}}_y^i=D_L^y\×q_y^{\′i},y=f,h,v---\left(3\right)$

实验结果：将2001年7月11日的MODIS超分辨率重构遥感图像(图8)和作为目标图像的2001年7月11日的Landsat ETM+遥感图像(图7)进行定量计算，通过峰值信噪比(PSNR)、结构相似性指数(SSIM)和视觉图像来评估本发明重构遥感图像的质量。

对于两个尺寸均为m×n的图像X和Y，图像的最大值为MAX，μ_x,μ_y,σ_x,σ_y,σ_xy分别为图像X和Y的均值，标准差和共变异数，C₁,C₂都为常数。峰值信噪比(PSNR)和结构相似性指数(SSIM)计算公式如下：

$\mathrm{PSNR}=10\log \left(\frac{{\mathrm{mn}\left(\mathrm{MAX}\right)}^2}{\underset{i=0}{\overset{m-1}{\mathrm{\Σ}}}\underset{j=0}{\overset{n-1}{\mathrm{\Σ}}}{(X(i,j)-Y(i,j))}^2}\right)---\left(4\right)$

$\mathrm{SSIM}=\frac{(2{\mathrm{\μ}}_x{\mathrm{\μ}}_y+C_1)(2{\mathrm{\σ}}_{\mathrm{xy}}+C_2)}{({{\mathrm{\μ}}_x}^2+{{\mathrm{\μ}}_y}^2+C_1)({{\mathrm{\σ}}_x}^2+{{\mathrm{\σ}}_y}^2+C_2)}---\left(5\right)$

当图像块尺寸为22×22像素，块之间重叠大小为18×18像素时，本发明提出方法重构图像的PSNR为29.7051db,SSIM为0.2694。研究结果表明通过采用基于小波预处理和稀疏表示的卫星遥感图像超分辨率重构方法，可以有效提高MODIS低分辨率遥感图像的空间分辨率。由于MODIS和Landsat ETM+遥感数据分辨率时间分辨率不同，当某一观测地区当日没有Landsat ETM+高分辨率遥感图像时，可以采用本发明对当日MODIS低分辨率遥感图像进行超高分辨率重构，能够提供观测地区较多细节信息。

Claims (5)

1.一种基于小波预处理和稀疏表示的卫星遥感图像超分辨率重构方法，是已知相同观测地区、不同时间分辨率的高分辨率遥感图像和低分辨率遥感图像，对其它观测时间的低分辨率遥感图像进行超分辨率重构；分为字典训练和对低分辨率遥感图像进行重构两个步骤；

所述的字典训练的过程是，

1)对已知低分辨率遥感图像上采样，使其和已知高分辨率遥感图像具有相同的图像尺寸；

2)对已知高分辨率遥感图像、低分辨率遥感图像分别进行小波分解提取特征，获得低频、水平、和垂直高频三对子波段；

3)在已知高分辨率遥感图像、低分辨率遥感图像的三对子波段图像上提取图像块作为训练样本；

4)利用冗余字典训练算法对已知低分辨率图像块进行训练，得到已知低分辨率图像块对应的字典和稀疏表示系数；

5)计算已知高辨率卫星遥感图像块对应的字典；

所述的对低分辨率遥感图像进行重构的过程是，

1)对待重构的低分辨率遥感图像进行小波分解，提取特征，获得低频、水平、垂直和对角线高频四个子波段；

2)在待重构低分辨率遥感图像的低频、水平、和垂直高频三个子波段图像上提取图像块作为训练样本；

3)求解图像块在已知低分辨率遥感图像块对应的字典中的稀疏表示系数；

4)根据求得的已知高分辨率遥感图像块对应的字典，计算待重构遥感图像在低频、水平、和垂直高频三个子波段上对应的高分辨率遥感图像，结合待重构遥感图像的对角线高频子波段，通过小波逆变换合成高辨率遥感图像重构结果。

2.根据权利要求1所述的基于小波预处理和稀疏表示的卫星遥感图像超分辨率重构方法，其特征在于，

在字典训练的过程的第4)步中，所述的冗余字典训练算法，是K-SVD冗余字典训练算法；具体的是通过(1)式对已知低分辨率遥感图像块进行训练，得到三个尺寸大小为n²×m_a的低分辨率遥感图像块对应的字典和稀疏表示系数其中m_a是冗余字典的列数或者原子数；

$D_M^y,\left\{q_y^i\right\}=\arg \min \underset{i}{\mathrm{\Σ}}{\left|\right|M_y^i-D_M^y{q}_y^i\left|\right|}_2,s.t.\underset{i}{\mathrm{\Σ}}{\left|\right|q_y^i\left|\right|}_0\≤k_0,\∀i,y=f,h,v---\left(1\right)$

式(1)中，k₀为稀疏表示系数中非零元素的个数，f，h和v分别表示低频、水平、和垂直高频，i＝1,...,P，表示第i个训练样本，P表示训练样本的个数。

3.根据权利要求1所述的基于小波预处理和稀疏表示的卫星遥感图像超分辨率重构方法，其特征在于，

在字典训练的过程的第5)步中，所述的计算已知高辨率卫星遥感图像块对应的字典通过伪逆表达式(2)求得；

$D_L^y=L_y{Q_y}^T{\left(Q_y{Q_y}^T\right)}^{-1},y=f,h,v---\left(2\right)$

式(2)中，L_y是由作为其列组成，Q_y是字典训练的过程的第4)步求得的作为其列组成，其中T代表矩阵的转置。

4.根据权利要求1所述的基于小波预处理和稀疏表示的卫星遥感图像超分辨率重构方法，其特征在于，

在对低分辨率遥感图像进行重构过程的第3)步中，所述的求解图像块在已知低分辨率遥感图像块对应的字典中的稀疏表示系数，是利用OMP算法求解图像块在已知低分辨率遥感图像块对应的字典中的稀疏表示系数

5.根据权利要求1所述的基于小波预处理和稀疏表示的卫星遥感图像超分辨率重构方法，其特征在于，

在对低分辨率遥感图像进行重构过程的第4)步中，所述的通过小波逆变换合成高辨率遥感图像重构结果，是利用式(3)求得对应的分辨率增强后的图像块 ${\hat{T}}_f^i,{\hat{T}}_h^i,{\hat{T}}_v^i;$

${\hat{T}}_y^i=D_L^y\×q_y^{\′i},y=f,h,v---\left(3\right)$

再对所有图像块进行以上过程然后合并得到最后结合T_d通过小波逆变换合成重构结果