技术领域
[0001] 本
发明属于
图像处理技术领域，更进一步涉及遥感图像处理技术领域中的基于压缩感知的多光谱与全色图像超分辨融合方法。本发明可以应用在土地利用、城区识别和森林资源调查等遥感领域。本发明是将低空间
分辨率的多光谱图像与高空间分辨率的全色图像进行融合，同时利用压缩感知技术进行超分辨重建，以得到一幅同时具有高空间分辨率和光谱信息的融合图像。
背景技术
[0002] 多光谱图像和全色图像融合是指将从多光谱探测器获得的同一场景的多谱图像的信息综合利用，利用它们在
时空上的相关性及信息上的互补性，得到对景物更全面、清晰的描述，更有利于人眼的识别和机器的自动探测。多光谱图像具有丰富的光谱信息，但空间分辨率较低；全色图像具有丰富的空间细节信息，空间分辨率较高，但光谱信息较少。要获得同时具有高的光谱和空间分辨率的图像，可通过两个途径实现：一是采用更新型的
传感器，这样能获取更高
质量的图像，但对新型传感器的研制将耗费大量时间，且成本高；二是进行多光谱图像和全色图像的融合，这样能在现有传感器技术的
基础上获取更高质量的图像，且技术成本低。目前的多光谱图像融合技术主要有两类：一类是IHS色彩变换法、主成分分析法(PCA)等；另一类是
小波变换法、Laplacian金字塔分解等多尺度变换方法，但这些方法在多光谱图像和全色图像融合中存在着一个共同的问题：均采用“注入”的方法或运用分辨率变换提取全色图像的细节信息进行融合，利用全色图像简单的替代多光谱图像的
亮度分量或当
抽取全色图像的高频分量与多光谱图像的亮度分量的高频分量不一致时，会产生严重的
颜色失真和光谱扭曲；另一方面，融合后图像的空间分辨率受限于融合前全色图像的分辨率。
[0003] 湖南大学
申请的
专利“一种基于压缩传感理论的卫星遥感图像融合方法”(申请号201010283310.0，申请公布号CN 101996396A)中，该专利申请的方法是把压缩感知技术运用到多光谱图像的融合中，先向量化低空间分辨率多光谱图像和高空间分辨率全色图像；
然后构造高空间分辨率多光谱图像
块的稀疏表示字典；再用基追踪
算法求解高空间分辨率多光谱图像在该字典下的稀疏表示；再将稀疏表示与预先设定的过完备字典相乘得到高空间分辨率多光谱图像块向量；最后将向量表示转化为图像块得到融合图像。该方法能够较好的克服颜色失真和光谱扭曲，但是仍然存在的不足是，融合后的多光谱图像的分辨率取决于融合前全色图像的分辨率，只能得到与融合前全色图像相同分辨率的图像。
发明内容
[0004] 本发明的目的在于克服上述已有技术的缺点，提出一种基于压缩感知的多光谱图像与全色图像的超分辨融合方法，以提高空间分辨率，同时降低颜色失真和光谱扭曲。
[0005] 为实现上述目的，本发明在压缩感知的
框架下实现全色图像与多光谱图像的融合，其技术方案是首先将全色图像与多光谱图像作为压缩
采样后的观测值，将高分辨率多光谱图像与低分辨率多光谱图像构造的
观测矩阵，与高分辨率多光谱图像和高分辨率全色图像构造的观测矩阵进行联合；然后利用KSVD(K-Singular Value Decomposition)方法对高分辨率多光谱图像随机抽块进行字典学习，得到训练字典；最后利用压缩感知基追踪算法得到高分辨率多光谱图像。
[0006] 本发明的具体步骤如下：
[0007] (1)分别输入低分辨率多光谱图像和高分辨率全色图像；
[0008] (2)获得高分辨率联合采样值矩阵；
[0009] 2a)对低分辨率多光谱图像和高分辨率全色图像进行分块列向量化；
[0010] 2b)将分块列向量化的高分辨率全色图像矩阵连接到分块列向量化的低分辨率多光谱图像矩阵后，得到高分辨率联合采样值矩阵；
[0011] (3)对高分辨率的多光谱图像和低分辨率的多光谱图像构造观测矩阵，获得高-低观测矩阵；
[0012] (4)对高分辨率的多光谱图像和高分辨率的全色图像构造观测矩阵，获得高-高观测矩阵；
[0013] (5)将高-低观测矩阵连接到高-高观测矩阵后，得到高分辨率联合观测矩阵；
[0014] (6)利用KSVD方法对高分辨率多光谱图像进行字典学习，得到高分辨率过完备字典；
[0015] (7)利用基追踪算法求解高分辨率过完备字典下的高分辨率稀疏系数；
[0016] (8)重构高分辨率多光谱图像；
[0017] 8a)将高分辨率过完备字典与高分辨率稀疏系数相乘，得到列向量化的多光谱图像；
[0018] 8b)将列向量化的多光谱图像转化为图像块，恢复到原图像的
位置，得到超分辨率的多光谱图像。
[0019] 本发明与现有的技术相比具有以下优点：
[0020] 第一，本发明将多光谱图像融合引入到压缩感知领域来，克服了
现有技术中颜色失真和光谱的扭曲的
缺陷，使得本发明可以保持空间分辨率和光谱信息，有利于后期对图像的处理和识别。
[0021] 第二，本发明引入了矩阵联合和字典学习的方法，克服了现有技术中融合图像分辨率受限于融合前全色图像分辨率的缺陷，使得本发明得到的融合图像分辨率比现有技术的融合前全色图像分辨率明显提高。
[0022] 第三，本发明引入了压缩感知基追踪算法，克服了现有技术中传感器和融合过程中噪声的影响，使得本发明比现有技术的峰值
信噪比大大提高。
[0023] 本发明的技术过程和效果可结合以下
附图详细说明。
附图说明
[0025] 图2为本发明的仿真图。
[0026] 具体实施方法
[0027] 参照附图1，本发明具体实施方式如下：
[0028] 步骤1，分别输入低分辨率多光谱图像和高分辨率全色图像；
[0029] 分别读入待处理的低分辨多光谱图像和高分辨率的全色图像。
[0030] 本发明
实施例中输入的低分辨多光谱图像大小为125×125×4，分辨率为4m；高分辨率的全色图像大小为500×500，分辨率为1m。
[0031] 步骤2，获得高分辨率联合采样值矩阵；
[0032] 将低分辨率多光谱图像的红、绿、蓝、
近红外四个颜色通道图像分为5×5的小块然后列向量化，将第二块构成的列向量连接到第一块构成的列向量后，将第三块构成的列向量连接到第二块构成的列向量后，以此类推，低分辨率的多光谱图像变换为一个行数为100的矩阵，本发明实施例中低分辨率的多光谱图像变换为一个大小为100×625的矩阵。
[0033] 将高分辨率全色图像分为20×20的小块然后列向量化，将第二块构成的列向量连接到第一块构成的列向量后，将第三块构成的列向量连接到第二块构成的列向量后，以此类推，高分辨率的全色图像变成一个行数为400的矩阵，本发明实施例中高分辨率的多光谱图像变换为一个大小为400×625的矩阵。
[0034] 将高分辨率全色图像矩阵连接到低分辨率多光谱图像矩阵后，得到一个行数为500的高分辨率联合采样值矩阵，本发明实施例中高分辨率联合采样值矩阵大小为500×625。
[0035] 步骤3，对高分辨率的多光谱图像和低分辨率的多光谱图像构造观测矩阵，获得高-低观测矩阵；
[0036] 将低分辨率多光谱图像的红、绿、蓝、近红外四个颜色通道图像分为5×5的小块然后列向量化，将高分辨率多光谱图像的红、绿、蓝、近红外四个颜色通道图像分为40×40的小块然后列向量化；
[0037] 将高分辨率多光谱图像每8个
像素取平均得到低分辨率多光谱图像的1个像素，将高分辨率多光谱图像的8×8的一个小块取平均得到低分辨率多光谱图像的1×1的一个小块，高分辨率多光谱图像与低分辨率多光谱图像的
下采样比为8∶1；
[0038] 本发明实施例中按下式求解高-低观测矩阵：
[0039]
[0040] 为低分辨率的多光谱图像， 为高分辨率的多光谱图像，M为待求解的高低观测矩阵，vp为噪声向量，本发明实施例中高-低观测矩阵大小为100×6400。
[0041] 步骤4，对高分辨率的多光谱图像和高分辨率的全色图像构造观测矩阵，获得高-高观测矩阵；
[0042] 将高分辨率多光谱图像的红、绿、蓝、近红外四个颜色通道图像分40×40的小块然后列向量化，将高分辨率全色图像分为40×40的小块然后列向量化，将低分辨率全色图像分为20×20的小块然后列向量化；
[0043] 将高分辨率多光谱图像的红、绿、蓝、近红外四个颜色通道的每一个像素加权求和，得到高分辨率全色图像对应位置的一个像素，本发明实施例中按下式求解权值矩阵：
[0044] yPAN＝M2x+v2
[0045] YPAN为高分辨率的全色图像， 为高分辨率的多光谱图像，M2为待求解的权值矩阵，vp为噪声向量。本发明实施例中红、绿、蓝、近红外每一谱段所取的权值系数大小分别为0.1139，0.2315，0.2308，0.4239，权值矩阵大小为1600×6400。
[0046] 将高分辨率全色图像每2个像素取平均得到低分辨率全色图像的1个像素，将高分辨率全色图像的2×2的一个小块取平均得到低分辨率全色图像的1×1的一个小块，高分辨率多光谱图像与低分辨率多光谱图像的下采样比为2∶1，本发明实施例中按下式求解高分辨率全色图像下采样矩阵：
[0047]
[0048] 其中， 为高分辨率的全色图像， 为低分辨率的全色图像，M3为高分辨率全色图像下采样矩阵，v3为空间噪声，本发明实施例中高分辨率全色图像下采样矩阵大小为400×1600。
[0049] 将高分辨率全色图像下采样矩阵与权值矩阵相乘得到高-高观测矩阵，本发明实施例中高-高观测矩阵大小为400×6400。
[0050] 步骤5，将高-低观测矩阵连接到高-高观测矩阵后，得到高分辨率联合观测矩阵；
[0051] 将步骤3得到的高-低观测矩阵连接到步骤4得到的高-高观测矩阵后，得到高分辨率联合观测矩阵，本发明实施例中高分辨率联合观测矩阵大小为500×6400。
[0052] 步骤6，利用KSVD方法对高分辨率的多光谱图像进行字典学习，得到高分辨率过完备字典；
[0053] 对输入的与待融合图像统计特性相近且比全色图像高一倍分辨率的20幅高分辨率多光谱图像进行分块列向量化，随机抽取高分辨率多光谱图像的25000个图像块进行字典训练；
[0054] 本发明实施例中按照下式求解训练字典D：
[0055] Subject to
[0056] 其中，min{}为最小化， 为求解Y-DX的二范数，Y为输入的初始离散余弦字典，D为待求解的高分辨率训练字典，X为随机抽取高分辨率多光谱图像块，Subject to表示‖Xi‖0为求解 的限制条件， 为任意第i列，‖Xi‖0为Xi的0范数，T0为稀疏度控制系数，本发明实施例中训练字典D大小为6400×6400。
[0057] 步骤7，利用基追踪算法求解高分辨率过完备字典下的高分辨率稀疏系数；
[0058] 本发明实施例中按照下式求解高分辨率的稀疏系数：
[0059] min‖α‖1 Subject to
[0060] 其中，min‖α‖1为最小化α的一范数，α为待求解的高分辨率稀疏系数，Subject to表示 为求解min‖α‖1的限制条件， 为求解y-Φα的二范数，y为高分辨率联合采样值矩阵，Φ为高分辨率联合观测矩阵与高分辨率训练字典的乘积，ε为重构误差。
[0061] 本发明实施例中ε设置为1，高分辨率联合采样值矩阵大小为500×625，高分辨率联合观测矩阵大小为500×6400。
[0062] 步骤8，重构高分辨率多光谱图像；
[0063] 将高分辨率过完备字典与高分辨率稀疏系数相乘，得到列向量化的多光谱图像；
[0064] 将列向量化的多光谱图像转化为40×40的图像块，恢复到原图像的位置，得到超分辨率的多光谱图像。
[0065] 本发明实施例中高分辨率稀疏系数矩阵大小为6400×625，列向量化的多光谱图像大小为6400×625，超分辨率的多光谱图像大小为1000×1000×4，分辨率为0.5m。
[0066] 本发明的效果可以通过以下仿真实验作进一步说明。
[0067] 实 验 仿 真 环 境 为 MATLAB 7.10.0(R2010a)Service Pack 3，Intel(R)Pentium(R)4CPU 3.00GHz，WindowXP Professional，图2(a)和图2(b)是一组待融合图像，图2(c)是一幅参考图像，这三幅图像都是由Quickbrid公司拍摄的，图2(d)为采用本发明对图2(a)和图2(b)待融合图像进行融合后获得的超分辨多光谱图像。
[0068] 将图2(d)与图2(b)进行比较，图2(b)的边缘不清晰而且细节也保持的不好，如河岸、桥、树木、房顶等，而图2(d)边缘清晰，细节也保持的很好，在空间分辨率、光谱信息、视觉效果上也明显优于图2(b)，在分辨率上得到比融合前全色图像高一倍分辨率的多光谱图像。由此可以看出，本发明融合的超分辨率多光谱图像重构效果好。
[0069] 为了证明本发明的效果，分别用本发明的方法和现有技术小波变换法、HIS色彩变换法、Contourlet变换法、PCA变换法对图2(a)和图2(b)的待融合图像进行融合。通过常用评价指标公式分别计算融合结果图像和图2(c)参考图像的红R、绿G、蓝B、近红外NIR及全部波段的像素矩阵，得到相关系数CC、均方根误差RMSE、光谱弧度SAM、相对整体维数综合误差ERGAS和整体图像质量指数Q4如下表。
[0070]
[0071] 由上表可见，本发明的相关系数CC和整体图像质量指数Q4均大于现有技术的评价值，相关系数越接近1，表示融合图像与参考图像越接近，融合结果越好，整体图像质量指数用来客观评估光谱信息的保留程度，结果范围在[0，1]，整体图像质量指数越大越好，越与参考图像相似，相关系数和Q4的值越大表示融合效果越好；均方根误差RMSE、光谱弧度SAM和相对整体维数综合误差ERGAS均小于现有技术的评价值，均方根误差表示融合图像与参考图像误差的大小，均方根误差越小，图像融合的效果越好，光谱弧度表示光谱的曲折程度，越接近0，融合结果越好，相对整体维数综合误差越小融合图像和参考图像越接近，融合的效果就越好，理想情况为0，均方根误差RMSE、光谱弧度SAM、相对整体维数综合误差ERGAS越小表示融合效果越好。由此可以看出本发明的评价结果均优于现有技术的评价结果，本发明比已有技术有较好的客观评价效果。
