高光谱遥感的特点是在获得目标2维空间信息分布的同时获取描述其属性的精确1维光谱特征，实现对目标的空间属性和物质属性的定量化分析(童庆禧等，2016；Pu，2017)。但受到高光谱成像硬件设备的限制以及现实物理条件与成本的制约，当前高光谱图像的空间分辨率还比较低，无法满足同时获取所观测目标的高分辨率空间信息和光谱信息的应用需求。而空谱融合可以充分融合不同传感器数据在空间维和光谱维的信息，得到同时具有高空间分辨率和高光谱分辨率的融合图像(肖亮等，2020)。此外，现有方法对待融合图像空间分辨率之比很大时的空谱融合问题还考虑不足，待融合图像的空间分辨率之比较小，大多为2或者4。而Sentinel-2可以提供10 m空间分辨率的多光谱图像，3倍的空间分辨率之比有利于与现有融合方法进行对比分析。因此，深入研究高分五号(GF-5)高光谱图像的空谱融合方法，实现国产卫星数据的有效应用，具有十分重要的实用意义。

当前针对高光谱图像的空谱融合方法仍处于发展阶段，主要分为3类：将传统多光谱图像空谱融合方法应用到高光谱图像的融合方法、基于模型的融合方法和基于深度学习的融合方法。

将传统多光谱图像空谱融合方法应用到高光谱图像主要是将成分替换法(cmponent substitution, CS)和多分辨率分析法(multiresolution analysis, MRA)这两类中的方法应用于高光谱图像空谱融合。CS中应用于高光谱图像融合最具代表性的方法是Gram-Schmidt(GS)融合方法。GS变换通过将模拟的低分辨率全色图像和高光谱图像进行GS正变换获得各个GS分量，再用直方图匹配后的全色图像替换GS变换的第1分量，与其余GS分量一起进行GS逆变换获得高空间分辨率高光谱图像(Dalla Mura等，2015)。该算法可以对任意波段数的遥感图像进行融合，简单而高效地提升图像的空间分辨率，但是GS方法融合结果的好坏主要取决于使用高光谱图像拟合的模拟图像。高空间分辨率图像与被替换的模拟图像之间的相关性越大，融合图像的光谱失真越小(Yokoya等，2017)。因此，GS变换应用于高光谱图像融合时仍有很大的提升空间。MRA中具有代表性的方法是基于亮度调节的平滑滤波方法(smoothing filter-based intensity modulation, SFIM) (王密等，2018)和拉普拉斯金字塔融合方法(generalized Laplacian pyramid, GLP) (Vivone等，2019)。SFIM基于简化的太阳辐射和地面反射模型，在逐像素的基础上将上采样的低分辨率图像乘以高分辨率图像与其低通滤波后图像的比值，从而提高了图像的清晰度。由于SFIM使用了盒式滤波器，所以融合图像的清晰度提升有限。在GLP中，将高分辨率图像与其低通滤波后图像之间的差异乘以增益因子来获得低分辨率波段的空间细节信息，但是应用于高光谱图像融合时会有较严重的光谱失真现象。

基于模型的融合方法主要分为基于统计估计的融合方法和基于混合光谱分解的融合方法。基于统计估计的融合方法通过观测得到的低分辨率高光谱图像和高空间分辨率多光谱图像求解图像降质的逆过程。在求解过程中引入先验知识降低问题复杂度，再通过最大后验概率(maximum a posteriori, MAP)求解。在MAP方法中加入随机混合模型(stochastic mixing model, SMM)即MAPSMM方法(Xie等，2020)。该方法实现多光谱图像对高光谱所有波段空间分辨率的提升，但其依赖于先验假设，在实际运用中计算量大、运算时间长是其主要缺点。基于混合光谱分解的融合方法假设高空间分辨率的高光谱图像可以表示为由低空间分辨率高光谱图像所构成的端元和高空间分辨率多光谱图像所确定的丰度联合实现。Zhang等人(2016)利用耦合非负矩阵分解(coupled nonnegative matrix factorization，CNMF)实现了高光谱与多光谱图像的融合，该算法结合了传感器观测模型和线性光谱混合模型，以无约束的欧几里得距离作为代价函数，并通过乘法更新规则更新光谱混合模型中的端元光谱矩阵和端元丰度矩阵，取得了良好的融合效果。但是，此方法需要预先知道有关端元的先验知识，如果缺少端元信息，在端元丰度矩阵提取时会产生很大的误差，影响融合效果。

基于深度学习的融合方法由于其非线性、深度特征学习的能力，在空谱融合领域展现出优越的性能。Palsson等人(2017)使用一个3D卷积神经网络(three-dimensional convolutional neural network, 3DCNN)来学习低空间分辨率高光谱图像和高空间分辨率多光谱图像空间域的映射关系，实现了高光谱图像和多光谱图像的空谱融合。Yang等人(2018)提出了2个分支结构的CNN，分别学习低空间分辨率高光谱图像(hyperspectral image, HS)和高空间分辨率多光谱图像(high-resolution multispectral image, HMS)的光谱和空间特征，实现高光谱图像空谱融合。目前应用于高光谱图像空谱融合的网络比较简单，融合效果一般。且基于深度学习的方法建立在大量数据集训练的基础上，而当前高光谱图像来源较少，与之匹配的高空间分辨率图像更少，这限制了深度学习在高光谱图像融合中的应用。

针对目前的空谱融合方法应用到高光谱图像融合时，会存在空间细节信息提升明显但光谱失真，或者光谱保真但空间细节信息提升不足的问题，本文通过GS变换与非下采样轮廓波变换(nonsubsampled contourlet transform, NSCT)提出一种波段自适应细节注入的高光谱图像空谱融合方法。首先，针对两个多波段图像不便于直接融合的问题，本文提出一种波段自适应的分组融合策略，即以相关系数为标准对待融合图像进行分组，将多光谱图像每个波段与对应的多个高光谱图像波段进行分组融合。其次，传统GS融合方法用平均权重系数模拟低分辨率图像造成融合图像光谱失真，本文通过最小均方误差(minimizing mean square error, MMSE)估计获取权重系数模拟低空间分辨率的多光谱图像(low-resolution multispectral image, LMS)作为第1分量进行GS正变换，改善融合图像光谱失真的问题。最后，为了能同时注入更多的空间细节信息，本文利用NSCT的平移不变性以及良好的空间、频域特性(Meng等，2017)，通过NSCT分解模拟图像、空间细节信息图像和多光谱图像，把空间、光谱信息融入到重构的高空间分辨率多光谱图像中，提升图像融合的效果。

提出的波段自适应细节注入融合方法主要从两个方面对融合图像的光谱和空间细节信息进行提升。一方面，针对传统GS方法应用到高光谱图像融合时光谱失真比较严重的问题，参考BDSD(band-dependent spatial-detail)模型(Vivone等，2015)求取模拟图像权重系数的思路，提出了一种波段自适应的低分辨率图像模拟方法，可以获得光谱信息更丰富的模拟图像。另一方面，为了保证融合图像的细节纹理信息，采用NSCT将模拟图像、空间细节信息图像和多光谱图像的空间以及部分光谱信息注入到重构的高空间分辨率图像中，再把该图像作为第1分量与其余GS分量进行逆变换，得到最终的高空间分辨率高光谱图像。波段自适应的细节注入融合方法具体流程如图 1所示，图中⊖为求差值运算。

GS变换是线性代数常用的多维正交变换。在可内积的空间中，GS变换可以为任意一组相互独立的向量找到一组正交基。设$\left\{ {{\mathit{\boldsymbol{x}}_1}, {\mathit{\boldsymbol{x}}_2}, \cdots, {\mathit{\boldsymbol{x}}_n}} \right\}$是一组相互独立的向量，设$\left\{ {{\mathit{\boldsymbol{y}}_1}, {\mathit{\boldsymbol{y}}_2}, \cdots, {\mathit{\boldsymbol{y}}_n}} \right\}$为GS变换构造的正交向量，变换过程为

$ \begin{gathered} \boldsymbol{y}_{1}=\boldsymbol{x}_{1} \\ \boldsymbol{y}_{2}=\boldsymbol{x}_{2}-\frac{\boldsymbol{y}_{1} \cdot \boldsymbol{x}_{1}}{\boldsymbol{y}_{1} \cdot \boldsymbol{y}_{1}} \boldsymbol{y}_{1} \\ \vdots \\ \boldsymbol{y}_{n}=\boldsymbol{x}_{n}-\frac{\boldsymbol{y}_{1} \cdot \boldsymbol{x}_{n}}{\boldsymbol{y}_{1} \cdot \boldsymbol{y}_{1}} \boldsymbol{y}_{1}-\frac{\boldsymbol{y}_{2} \cdot \boldsymbol{x}_{n}}{\boldsymbol{y}_{2} \cdot \boldsymbol{y}_{2}} \boldsymbol{y}_{2}-\cdots-\frac{\boldsymbol{y}_{n-1} \cdot \boldsymbol{x}_{n}}{\boldsymbol{y}_{n-1} \cdot \boldsymbol{y}_{n-1}} \boldsymbol{y}_{n-1} \end{gathered} $

(1)

一方面，使用GS变换对图像进行正交变换，消除各波段之间的相关性，且包含的信息量相近，改进了主成分分析方法中信息过分集中的问题。另一方面，GS变换应用于遥感图像融合时不受波段数的限制。因此，GS变换是一种比较适用于高光谱遥感影像空谱融合的方法。

NSCT是一种多尺度、多方向的图像分解方法(如图 2)，使用了双滤波器组结构，包括多尺度性质的非下采样金字塔(nonsubsampled pyramid, NSP)分解结构和多方向分解性质的非下采样方向滤波器组(nonsubsampled dirctional filter bank, NSDFB)。首先使用NSP对原图像进行多尺度分解，将输入的2维信号分解成一个低通子带图像和一个带通方向子带图像。每级的NSP分解继续在低通子带图像上迭代进行，得到不同尺度的子带系数。对于各尺度的带通方向子带图像，利用NSDFB进行多方向分解，最后得到不同尺度下不同方向的子带系数。由于NSCT整个分解过程没有上采样和下采样，因此所有的分解子带都和原图像的大小相同，从而使得NSCT获得了平移不变性。

对某尺度下的子带图像进行$k$级方向分解，可得到${2^k}$个方向子带图像；原图像经过$L$级NSCT分解，可以得到1个低频子带图像和$\sum\limits_{l = 1}^{L - 1} {{2^{{k_l}}}} $个带通方向子带图像，其中$k_l$为尺度$l$下的方向分解级数。

具有良好平移特性的多尺度NSCT变换可以使同一尺度变换下的高频分量表现多个方向的细节信息，更好地表现图像的细节特征(Li等，2018)。因此，本文选用GS变换与NSCT结合的融合方法，在保真高光谱遥感图像光谱信息的同时，尽可能地提升图像的空间细节信息。

多光谱图像与高光谱图像的波谱范围通常有很大的重叠。由于基于GS变换的图像融合方法对于待融合图像波谱重叠范围内的波段会有更好的融合效果(Meng等，2019)，而在多光谱图像和高光谱图像波谱重叠范围内，多光谱图像的一个波段对应着多个高光谱图像波段，这种多对多的关系不利于后续融合操作。因此，本文在进行融合前先对待融合的高光谱图像和多光谱图像进行分组操作。

对于多光谱波段范围内的高光谱波段，它们之间的分组关系是确定的，直接将高光谱图像的波段分配到对应多光谱图像的组中。不过，相邻多光谱波段的波谱范围并不是连续的(如图 3所示，4条彩色曲线是Sentinel-2的4个10 m空间分辨率波段的光谱响应函数，绿色直线为GF-5所有波段的中心波长)，对于它们之间的高光谱波段，通过计算相关系数求解其对应的分组会更精确。对于从近红外到短波红外范围的高光谱波段，一般没有多光谱波段的覆盖(高光谱卫星的波段范围一般是400~2 500 nm，多光谱卫星的波段范围一般是400~1 000 nm)，采用相关系数作为分组标准。首先对高空间分辨率多光谱图像HMS按照高光谱图像的尺寸大小进行下采样，得到下采样后的多光谱数据$\widetilde {\mathit{\boldsymbol{MS}}}$。然后计算高光谱图像$\mathit{\boldsymbol{HS}}$每个波段和$\widetilde {\mathit{\boldsymbol{MS}}}$各个波段之间的相关系数

$ \rho(\boldsymbol{H S}, \widetilde{\boldsymbol{M S}})=\frac{{Cov}(\boldsymbol{H} \boldsymbol{S}, \widetilde{\boldsymbol{M S}})}{\sigma_{1} \times \sigma_{2}} $

(2)

式中，$Cov(\mathit{\boldsymbol{HS}}, \widetilde {\mathit{\boldsymbol{MS}}}) = E\left\{ {\left({\mathit{\boldsymbol{HS}} - {\mu _1}} \right)\left({\widetilde {\mathit{\boldsymbol{MS}}} - {\mu _2}} \right)} \right\}$，$E\left\{ \cdot \right\}$为数学期望，${{\mu _1}}$、${{\mu _2}}$分别为$\mathit{\boldsymbol{HS}}$、$\widetilde {\mathit{\boldsymbol{MS}}}$的均值，${{\sigma _1}}$、${{\sigma _2}}$分别为$\mathit{\boldsymbol{HS}}$、$\widetilde {\mathit{\boldsymbol{MS}}}$的标准差。最后得到一个${\mathit{\boldsymbol{HS}}}$每个波段与${\widetilde {\mathit{\boldsymbol{MS}}}}$所有波段的相关系数表，根据表 1中最大值将HS与对应的HMS分为一组。以本文使用的数据集为例，HMS有4个波段，HS有330个波段，得到一张330×4的相关系数表，根据每1行的最大值可以将HS所有波段与对应的HMS波段分为一组，示例分组结果如表 1所示，分组过程如图 4所示。

基于GS变换的图像融合方法中非常重要的一步就是模拟生成低分辨率图像用于GS正变换，生成模拟图像的好坏直接影响融合结果。理论上，重叠范围内的多光谱波段的辐射能量等于对应高光谱波段辐射能量之和。但是传感器的光谱响应曲线并不是理想的。考虑到传感器的光谱响应特性，可以由对应波谱范围内的高光谱波段的辐射值近似估算多光谱波段的辐射值。在先前的一些多光谱图像与全色图像的空谱融合算法中，假设低分辨率的全色图像($\mathit{\boldsymbol{PA}}{\mathit{\boldsymbol{N}}_L}$)和多光谱图像各个波段存在多元一次线性回归关系，即

$ \boldsymbol{P} \boldsymbol{A} \boldsymbol{N}_{L}=\omega_{1} \boldsymbol{B}+\omega_{2} \boldsymbol{G}+\omega_{3} \boldsymbol{R}+\omega_{4}\boldsymbol{ N I R}+a $

(3)

式中，$\mathit{\boldsymbol{PAN}}$、$\mathit{\boldsymbol{B}}$、$\mathit{\boldsymbol{G}}$、$\mathit{\boldsymbol{R}}$、$\mathit{\boldsymbol{NIR}}$分别代表低分辨率全色、蓝色、绿色、红色、近红外波段的辐射值。由于全色图像波谱范围与多光谱图像不完全重合，因此存在常数$a$。类似地，高光谱图像与多光谱图像融合时，对于多光谱图像某一波段的低分辨率图像$(\boldsymbol{M S}_{\mathrm{L}}^{(k)})$，可以用对应的多个高光谱波段($\boldsymbol{H S}_{i}$)表示

$ \boldsymbol{M S}_{\mathrm{L}}^{(k)}=\sum\limits_{i=1}^{N} \omega_{k, i} \boldsymbol{H S}_{i}, \quad k=1, 2, \cdots, M $

(4)

式中，$\omega_{k, i}$是用多个高光谱波段拟合对应多光谱波段的系数。为了改善传统GS融合方法中通过平均权重系数拟合低分辨率图像的不足，本文借鉴BDSD模型求取权重系数时使用的MMSE估计，求取$\boldsymbol{M S}_{\mathit{L}}^{(k)}$和$\boldsymbol{H S}_{i}$对应之间的拟合系数${\omega _k}$，即

$ \omega_{k}=\left(\boldsymbol{H}_{d}^{\mathrm{T}} \boldsymbol{H}_{d}\right)^{-1} \boldsymbol{H}_{d}^{\mathrm{T}}\left(\boldsymbol{H} \boldsymbol{S}-\boldsymbol{H} \boldsymbol{S}^{\mathrm{LP}}\right) $

(5)

式中，${\mathit{\boldsymbol{H}}_d} = \left[ {\mathit{\boldsymbol{H}}{\mathit{\boldsymbol{S}}_1}, \mathit{\boldsymbol{H}}{\mathit{\boldsymbol{S}}_2}, \cdots, \mathit{\boldsymbol{H}}{\mathit{\boldsymbol{S}}_n}, \mathit{\boldsymbol{M}}{\mathit{\boldsymbol{S}}^{{\rm{LP}}}}} \right]$，${\mathit{\boldsymbol{M}}{\mathit{\boldsymbol{S}}^{{\rm{LP}}}}}$为经过低通滤波和下采样处理后的多光谱图像，${\mathit{\boldsymbol{H}}{\mathit{\boldsymbol{S}}^{{\rm{LP}}}}}$为低通滤波处理后的高光谱图像。

通过建立${\mathit{\boldsymbol{MS}}_{\rm{L}}^{\left(k \right)}}$与${\mathit{\boldsymbol{H}}{\mathit{\boldsymbol{S}}_i}}$之间的回归方程，采用最小均方误差估计求解拟合系数，进而拟合${\mathit{\boldsymbol{MS}}_{\rm{L}}^{\left(k \right)}}$。将${\mathit{\boldsymbol{MS}}_{\rm{L}}^{\left(k \right)}}$作为第1分量，与${\mathit{\boldsymbol{H}}{\mathit{\boldsymbol{S}}_i}}$按照式(1)进行GS正变换，可以得到正交分量$\left({\mathit{\boldsymbol{G}}{\mathit{\boldsymbol{S}}_1}, \mathit{\boldsymbol{G}}{\mathit{\boldsymbol{S}}_2}, \cdots, \mathit{\boldsymbol{G}}{\mathit{\boldsymbol{S}}_n}} \right)$。其中，第$n$分量计算为

$ \begin{aligned} &\boldsymbol{G S}_{n}(i, j)=\left(\boldsymbol{H} \boldsymbol{S}_{n}(i, j)-\mu_{n}\right)- \\ &\sum\limits_{m=1}^{n-1} {Cov}\left(\boldsymbol{H} \boldsymbol{S}_{n}, \boldsymbol{G} \boldsymbol{S}_{m}\right) \times \boldsymbol{G} \boldsymbol{S}_{m}(i, j) \end{aligned} $

(6)

式中，$\boldsymbol{G S}_{n}$表示GS变换的分量；$\boldsymbol{H S}_{n}$表示高光谱第$n$波段的DN(digital number)值；${{\mu _n}}$为高光谱波段平均值；${Cov}\left(\boldsymbol{H} \boldsymbol{S}_{n}, \boldsymbol{G} \boldsymbol{S}_{m}\right)$为${\mathit{\boldsymbol{H}}{\mathit{\boldsymbol{S}}_n}}$与${\mathit{\boldsymbol{G}}{\mathit{\boldsymbol{S}}_n}}$的协方差。

将${\mathit{\boldsymbol{MS}}_{}^{\left(k \right)}}$与对应的${\mathit{\boldsymbol{MS}}_{\rm{L}}^{\left(k \right)}}$作差，得到空间细节信息${\mathit{\boldsymbol{\delta }}^{\left(k \right)}}$，即

$ \boldsymbol{\delta}^{(k)}=\boldsymbol{M} \boldsymbol{S}^{(k)}-\boldsymbol{M} \boldsymbol{S}_{\mathrm{L}}^{(k)} $

(7)

为了提高最终融合图像的空间细节分辨能力，本文通过NSCT变换将更多的空间细节信息保留到融合图像中。

将含有空间信息的${\mathit{\boldsymbol{\delta }}^{\left(k \right)}}$、含有光谱信息的${\mathit{\boldsymbol{MS}}_{\rm{L}}^{\left(k \right)}}$以及高空间分辨率的${\mathit{\boldsymbol{MS}}^{\left(k \right)}}$分别进行NSCT分解，分离出高低频系数$\left\{ {{L_i}, {H_{i, k}}} \right\}$。其中，$L_i$是分解第$i$层的低频系数；${{H_{i, k}}}$是第$i$层$k$方向的高频系数。由于低频分量包含了主要能量，决定图像的轮廓信息，因此将${\mathit{\boldsymbol{MS}}_{\rm{L}}^{\left(k \right)}}$和由NSCT分解得到的低频系数采用加权融合的方法，即

$ L_{\mathrm{F}}=\gamma_{1} L_{\mathrm{MS}_{\mathrm{L}}}+\gamma_{2} L_{\mathrm{MS}} $

(8)

式中，$L_{\mathrm{F}}$为融合后的低频系数，$L_{\mathrm{MS}_{\mathrm{L}}}$、$L_{\mathrm{MS}}$分别为${\mathit{\boldsymbol{MS}}_{\rm{L}}^{\left(k \right)}}$和${\mathit{\boldsymbol{MS}}^{\left(k \right)}}$的低频系数，${\gamma _1}$、${\gamma _2}$分别为${\mathit{\boldsymbol{MS}}_{\rm{L}}^{\left(k \right)}}$和${\mathit{\boldsymbol{MS}}^{\left(k \right)}}$的权值系数，计算为

$ \left\{\begin{array}{l} \gamma_{1}=\mathrm{e}^{-\frac{\left(L_{\mathrm{M S}}-\mu\right)^{2}}{2 \sigma^{2}}} \\ \gamma_{2}=1-\gamma_{1} \end{array}\right. $

(9)

式中，$\mu $为系数${L_{{\rm{MS}}}}$的均值，$\sigma $为系数${L_{{\rm{MS}}}}$的方差。

对于${\mathit{\boldsymbol{\delta }}^{\left(k \right)}}$的高频系数和${\mathit{\boldsymbol{MS}}^{\left(k \right)}}$的高频系数采用模值取大的方法，即

$ \begin{array}{cc} \boldsymbol{H}_{i, k}(x, y)= \\ \begin{cases}H_{i, k}^{\mathrm{MS}}(x, y) & \left|H_{i, k}^{\mathrm{MS}}(x, y)\right| \geqslant\left|H_{i, k}^{\delta}(x, y)\right| \\ H_{i, k}^{\delta}(x, y) & \left|H_{i, k}^{\mathrm{MS}}(x, y)\right|<\left|H_{i, k}^{\delta}(x, y)\right|\end{cases} \end{array} $

(10)

式中，$\boldsymbol{H}_{i, k}(x, y)$为融合后的高频系数，$H_{i, k}^{\mathrm{MS}}(x, y)$、$H_{i, k}^{\rm{ \mathsf{ δ} }}(x, y)$分别为$\mathit{\boldsymbol{M}}{\mathit{\boldsymbol{S}}^{\left(k \right)}}$和${\mathit{\boldsymbol{\delta }}^{\left(k \right)}}$的高频系数。

将融合后的高低频系数进行NSCT重构，得到新的高空间分辨率同时包含一定光谱信息的多光谱图像$\mathit{\boldsymbol{MS}}_{{\rm{new}}}^{\left(k \right)}$。将$\mathit{\boldsymbol{MS}}_{{\rm{new}}}^{\left(k \right)}$和式(6)得到的GS变换其他分量$\left({\mathit{\boldsymbol{G}}{\mathit{\boldsymbol{S}}_2}, \mathit{\boldsymbol{G}}{\mathit{\boldsymbol{S}}_3}, \cdots, \mathit{\boldsymbol{G}}{\mathit{\boldsymbol{S}}_n}} \right)$一起进行GS逆变换，可得到最终的融合图像。GS逆变换计算为

$ \begin{aligned} &\boldsymbol{H} \boldsymbol{S}_{n}(i, j)=\left(\boldsymbol{G} \boldsymbol{S}_{n}(i, j)+\mu_{n}\right)+ \\ &\sum\limits_{m=1}^{N-1} {Cov}\left(\boldsymbol{H} \boldsymbol{S}_{n}, \boldsymbol{G} \boldsymbol{S}_{m}\right) \times \boldsymbol{G} \boldsymbol{S}_{m}(i, j) \end{aligned} $

(11)

选用2019年5月1日获取的GF-5卫星影像和2019年5月28日获取的Sentinel-2卫星影像作为实验数据。实验所用多光谱数据为Sentinel-2影像10 m空间分辨率的4个波段(Band 2, Band 3, Band 4和Band 8)；使用的高光谱数据为高分五号30 m空间分辨率的330个波段，前150个波段的光谱分辨率为5 nm，后180个波段的光谱分辨率为10 nm(刘银年等，2020；孙允珠等，2018)。在图像融合前，对两幅图像进行有理多项式系数(rational polynomial coefficients, RPC)校正、几何纠正和配准等预处理操作，并用数字高程模型(digital elevation model, DEM)数据辅助进行正射校正、几何精校正和图像精确配准。

为了评价本文提出融合方法的有效性与可靠性，将融合结果与已经应用到高光谱图像融合的经典算法——GS、SFIM、GLP、CNMF和MAPSMM进行比较。本文根据不同地物类型组成的场景，对融合图像的光谱保真度、空间细节信息进行对比分析。6种方法的融合结果如图 5—图 7所示，融合图像使用了《GF-5产品使用说明》中推荐的Band 62，Band 38和Band 15组成真彩色进行展示，图像大小为240×240像素。

区域1(见图 5)的主要地物是建筑、道路、裸地和植被。从整体来看，GS、本文方法和CNMF都很好地注入了MS图像的细节信息，空间分辨率有明显的提升；而GLP、MAPSMM和SFIM都出现了虚影，空间分辨率虽然有提升，但是细节注入明显没有前3种方法好。从局部来看，6种方法对于地物细节都有提升，例如原始HS下方红顶厂房轮廓并不明显，而融合后的图像可以看到两个红色棚顶的边缘以及中间的空隙。对于大面积的连续区域如裸地和植被，GS、本文方法和CNMF过渡平滑，没有多余的冗杂信息；而GLP、MAPSMM和SFIM均出现了马赛克状的斑块。对于图像中部和下部的小建筑群，GS、本文方法和CNMF很好地注入了MS中的细节信息，而GLP、MAPSMM和SFIM也注入了细节信息，但是部分建筑旁边出现了黑色的噪点(如图 5红圈所示)。

区域2(见图 6)的主要地物是建筑、农田、裸地和河流。整体上，6种方法对空间分辨率都有提升，其中GS、本文方法和CNMF的地物细节更加清晰，而GLP、MAPSMM和SFIM都出现了虚影。6种方法的光谱都接近原始HS图像。局部上，对于蓝色顶棚这类大建筑，各个方法都注入了细节信息。但是对于小建筑群，GLP、MAPSMM和SFIM出现了明显的虚影，个别地区出现了黑色的噪点(如图 6红圈所示)。对于大面积的连续区域如农田和裸地，GS、本文方法和CNMF过渡平滑，而GLP、MAPSMM和SFIM均出现了多余的条纹纹理。GS、本文方法和CNMF中河流的轮廓更加锐利，而GLP、MAPSMM和SFIM中河流的边缘较为模糊，有些河段甚至出现了虚影。

区域3(图 7)的主要地物是湖泊、山脉、植被和建筑。从整体上看，GS、本文方法和CNMF的空间分辨率提升明显，而GLP、MAPSMM和SFIM与之相比更模糊，个别地区出现了黑色的噪点(如图 7红圈所示)。此外，6种方法与原始HS光谱差异都比较大。由于GF-5数据还未商用，数据较为稀少，而与之匹配的云量低的Sentinel-2数据更加稀少，这使得参与融合的图像时间差异大(相差27天，5月正是植被返青期，光谱变化较大)，因而山脉上植被覆盖差异较大，导致各个方法生成的融合图像光谱失真比较严重。

从局部来看，对于大面积的湖泊，GS与本文方法的湖岸轮廓清晰、纹理信息更丰富，但是水体受到MS的影响有白色斑块；CNMF与GLP的湖岸轮廓也很清晰，但CNMF的水体有蓝紫色的畸变，GLP的水体有墨绿色的畸变；MAPSMM和SFIM较为模糊，水体中出现马赛克状的斑块，水体也有墨绿色的畸变。对于山脉以及山脉上的植被，6种方法都为山脊、山谷等注入了更清晰的纹理，但是都产生了较为严重的颜色失真。对于小建筑群，GS和本文方法展示了更清晰的轮廓和细节，其他方法则较为模糊。

从定性分析角度，本文算法能较好提升高光谱图像的空间分辨率，同时具有较高的光谱保真度，融合图像中地物的边缘和纹理都比较清晰，图像的色调与高光谱图像基本一致。

为了更准确、客观地评价不同融合方法对空间分辨率的提升和对光谱的保真效果，本文选取了标准差(standard deviation, STD)、信息熵(entropy)、通用图像质量指数(universal image quality index, UIQI)、相关系数(correlation coefficient, CC)、相对无量纲全局误差(the erreur relative globale adimensionnelle de Synthèse, ERGAS)及光谱角映射(spectral angle mapper, SAM) 6个指标进行评价(Vivone等，2015, 2018)。其中，标准差反映融合图像的空间细节信息，UIQI、相关系数、ERGAS、SAM反映融合图像的光谱信息，信息熵反映融合图像整体的信息量，3类研究区的具体评价参数如表 2—表 4所示。

标准差表示图像像素灰度值相对于均值的离散程度，标准差越大表明图像中灰度级分布越离散，图像的空间细节信息也越丰富。在3类区域中，GS的标准差均为最大，符合定性分析中GS融合图像轮廓清晰、纹理丰富的视觉表现。本文方法的标准差在除区域2以外的区域比GS稍差；而在区域2各个方法对空间细节的提升都比较明显，标准差相差不大，这与定性分析中的情况相符合。CNMF的标准差在6种方法里处于中等水平，其融合图像的空间细节信息不如GS和本文方法，但是比GLP、MAPSMM和SFIM好。而对于GLP、MAPSMM和SFIM，其标准差在6种方法中排名靠后。从前面的融合结果也可以印证：这些方法的融合图像没有GS、本文方法、CNMF清晰，地物细节模糊、轮廓有虚影，一些局部区域甚至出现了马赛克状的斑块。

信息熵是评估图像包含信息量的指标，信息熵通过计算灰度分布的聚集程度来体现信息的丰富程度，信息熵越大表明图像信息越丰富。GS的信息熵在4类区域中均是最大的，其融合图像中的各类地物也是6种方法中最清晰的，图像信息丰富。在区域1和区域2，各个方法的信息熵相差不大；在区域3，本文方法的信息熵仅次于GS，这表明本文方法在山脉、湖泊和农田等大面积地物融合时可以获得更好的细节信息。CNMF、GLP、MAPSMM和SFIM虽然在有些区域的信息熵接近本文方法甚至略微超过本文方法，但是从目视解译可以看出，这些方法中的地物轮廓模糊、纹理粗糙，甚至有噪声点的存在，由于信息熵并未考虑图像中像元的位置关系，轮廓的虚影和噪点会造成信息熵虚大的问题。

通用图像质量指数UIQI从相关性失真、亮度失真和对比度失真3个方面衡量融合图像质量，UIQI的理想值为1，UIQI越接近1表示融合效果越好。在3类区域中，本文方法的UIQI相比于GS有明显提升，表明本文方法相比于GS在光谱保真度上有明显提高。本文方法的UIQI在区域1和区域2均排名前列，说明本文方法对于受时相变化影响不大的城镇区域有着很好的光谱保真度。对于区域3这种以植被为主的场景，本文方法的UIQI略低于GLP、MAPSMM和SFIM，但高于GS和CNMF，说明本文方法对于受时相变化影响大的植被密集区域的光谱保真度一般。

相关系数CC评价融合图像和原始图像之间的相关性。CC的最佳值为1，相关系数越大说明融合图像从原始图像获取的信息越多，融合图像光谱失真也越小。本文方法的CC相比于GS有明显提升，表明本文方法相比于GS在光谱保真度上有明显提高。在3类区域中，本文方法的CC均排名前列，说明本文方法从原始高光谱图像中获取的信息多，且光谱失真较小，取得了很好的融合效果。

相对无量纲全局误差ERGAS评价融合图像所有波段的光谱变化情况。ERGAS的理想值为零，值越小，表示融合图像的光谱质量越好。本文方法的ERGAS明显小于GS，表明本文方法的光谱保真度相比GS有明显的改善。本文方法的ERGAS在区域1和区域2的数值较低，说明本文方法对于受时相变化影响不大的城镇区域有着很好的光谱保真度。对于区域3这类受时相变化影响大的区域，本文方法ERGAS略高于GLP、MAPSMM和SFIM，但在对于地物空间细节的提升上优于这3种方法；同时本文方法的ERGAS比GS、NSCT这些空间分辨率提升明显的方法要低，说明本文方法在提升空间分辨率的同时，具有良好的光谱保真度。

光谱角映射SAM是用于测量光谱失真的指标。该指标把融合图像和原始高光谱图像同一个像素点的不同光谱带的像素值看做是一对向量，计算这对向量之间的角度即为光谱角度。SAM的理想值是零，该指标越小表示融合图像质量越好，融合图像的光谱失真越小。在3类区域中，本文方法的SAM都比GS小，说明本文方法在GS基础上光谱保真度有了很大的提升。本文方法的SAM大于GLP、MAPSMM和SFIM，表明从SAM指标来看本文的光谱保真度没有这3种方法好，但是本文在空间分辨率的提升上更有优势。而对比同样空间分辨率提升较为优秀的CNMF，本文的SAM基本都小于CNMF，光谱保真度比其更优秀。

综上，从定量分析的角度，本文方法在GS的基础上光谱保真度有较明显的提升，尤其是对于受时相影响较小的城镇区域，本文方法在保证地物轮廓清晰、纹理丰富且没有噪点的同时，光谱十分接近原始的高光谱图像。而对于受时相变化影响大的植被密集区域，本文方法的光谱保真度虽然略逊色于GLP、MAPSMM和SFIM，但是本文方法可以提供更清晰、细节信息更丰富的地物，且没有噪点出现。因此，本文方法是一种兼具高空间分辨率和高光谱保真度的高光谱图像融合方法。

光谱曲线评价是评价高光谱图像融合质量的一个重要指标。根据不同研究区域的主要地物，选取各个研究区的代表性地物进行光谱曲线评价(图 8)，3个研究区内主要地物(建筑、道路、山体、植被、水体、农田、裸地)的光谱曲线如图 8所示。通过比较融合图像和原始高光谱图像的光谱曲线形状以及光谱曲线的数值差异来评价融合图像的光谱质量。从图 8可以看出，本文方法在各种地物中的光谱曲线与原始高光谱图像的光谱曲线十分接近，每个波段的数值与高光谱图像的相差很小，说明本文方法的融合结果保留了原始高光谱图像的光谱特性。在农田场景，CNMF的光谱曲线虽然与原始高光谱图像的光谱曲线走势一致，但是部分波段数值出现了较为明显的偏离。在水体场景，各方法在1 000~2 000 nm的光谱曲线相比原始光谱曲线要低。

此外，针对GLP、MAPSMM和SFIM融合图像中明显的噪点，本文以研究区2左下角噪点(如图 9所示，图 9右下角为放大显示)为例，对噪点像元进行光谱曲线评价。对比图 9可以看出，GLP(图 9(f))、MAPSMM(图 9(g))和SFIM(图 9(h))红框中心像元周围出现了明显的黑圈，融合图像上显示出明显的斑点。CNMF(图 9(e))红框中心像元虽然没有出现噪点，但是对比其他方法可以发现，中心像元左边的红色像元发生了失真，像元变成了黑色。GS(图 9(c))和本文方法(图 9(d))没有出现噪点。进一步对噪点进行光谱曲线评价如图 10所示。

GS和本文方法(图 10(b))在500~700 nm和1 500~1 700 nm的光谱曲线出现了异常大的值，其他波段和高光谱图像的光谱曲线几乎保持一致。对比GS和本文方法的光谱曲线可以看出，本文方法相比GS更接近于原始高光谱图像的光谱曲线。CNMF(图 10(c))在500~700 nm出现了异常大的值，在1 500~2 500 nm光谱曲线的值也比原始高光谱的值大。GS、本文方法和CNMF光谱曲线在局部出现了异常值，但是光谱曲线的整体形态和高光谱的光谱曲线仍然很接近。GLP(图 10(d))在500~800 nm出现了异常大的值，而在1 000~2 500 nm光谱值全部变成零，GLP的光谱曲线发生了很大的扭曲，其噪点附近出现了明显的光谱失真。MAPSMM(图 10(e))和SFIM(图 10(f))在1 400~2 500 nm出现了65 535的异常极大值，光谱曲线发生了非常大的扭曲，它们的噪点附近出现了明显的光谱失真。

通过这些融合图像的噪点位置，对比融合前的HS和MS可以发现：这些噪点出现在MS亮度很大而HS亮度很小的区域，如图 9 (a)(b)所示。GLP、MAPSMM和SFIM对于这些亮度差别很大的点，在融合图像中会出现噪点，且噪点的光谱会严重扭曲。而本文方法、GS和CNMF对于这些亮度差别很大的点，依然可以保证融合图像空间上不出现噪点，融合图像的光谱曲线与原始高光谱图像的光谱曲线基本保持一致。

从光谱评价的角度来看，GLP、MAPSMM和SFIM虽然整体光谱保真度较高，但是空间细节信息提升不明显，而且融合图像出现的噪点在光谱上出现了严重的扭曲，这对于高光谱图像的后续应用来说是无法接受的。而本文方法在提升空间分辨率的同时兼顾光谱保真度，同时光谱曲线更接近原始高光谱图像的光谱曲线。而且，即使待融合的HS和MS的亮度差别很大，融合图像空间上没有出现噪点，光谱曲线与原始光谱基本保持一致，取得了较好的融合效果。

提出了一种波段自适应细节注入的高分五号与Sentinel-2遥感影像空谱融合方法。在融合前，对多光谱图像波谱范围之外的高光谱波段，采用了波段自适应的融合策略，利用相关系数作为波段分组标准，将高光谱图像与多光谱图像的波段进行分组，解决了两个多波段图像不便于直接融合的问题，同时也间接提升了融合图像的光谱保真度。分组后，本文在GS融合的基础上，一方面，使用最小均方误差估计计算用于GS正变换的模拟图像，改善了使用平均权重系数拟合低分辨率多光谱图像造成的光谱失真问题。另一方面，通过NSCT重构高空间分辨率图像，进一步增强了融合图像的空间细节信息。最后，不同地物类型的实验结果表明，本文方法对于受时相影响较小的城镇区域，在提升空间分辨率的同时有较好的光谱保真度，且不会出现噪点；对于受时相变化影响大的植被密集区域，本文方法融合图像有较好的清晰度和地物细节信息，且没有噪点出现。本文目前选择了加权融合NSCT低频系数和模值取大融合高频系数的融合规则，在后续研究中，可以进一步选择更合适的融合规则，生成空间细节更好同时光谱保真度更高的融合图像。