CN103489208B - 一种高光谱图像非线性光谱混合模型及地质填图应用研究 - Google Patents

Abstract

本发明属于高光谱遥感影像处理技术领域，具体为一种高光谱遥感影像光谱解混方法。针对线性光谱混合模型解混精度不高的问题，提出一种基于二次散射的高光谱遥感影像非线性光谱混合模型，并定义了一种模型中二次散射系数的确定方法。研究了基于该模型解混结果的遥感地质填图，提出了3种填图规则，通过与线性模型填图结果对比，利用二次散射模型进行解混能得到更精确的填图结果。

Description

一种高光谱图像非线性光谱混合模型及地质填图应用研究

技术领域

本发明属于高光谱遥感影像处理技术领域，具体为一种高光谱遥感影像光谱解混方法。

背景技术

20世纪80年代兴起了新型的对地观测技术-高光谱遥感技术。它把人们研究地表信息的能力由陆地推向太空，拓宽了人们的研究视野，极大地提高了人们宏观、准确、及时、综合地对地观测与监测能力。成像光谱技术把遥感波段从几个、几十个推向数百个、上千个，其覆盖的光谱范围也从可见近红外、短波红外，逐步向紫外和中红外拓展。高光谱遥感数据每个像元可以提供几乎连续的地物光谱曲线，使人们利用高光谱反演陆地细节成为可能。其应用领域已涵盖地球科学的各个方面，在地质找矿和制图、大气和环境监测、农业和森林调查、海洋生物和物理研究等领域发挥着越来越重要的作用。高光谱遥感技术已经成为当前遥感领域的前沿技术。

然而，在遥感成像系统、高光谱遥感成像光谱仪的设计中，光谱分辨率和空间分辨率是不可兼得的(参考对比文件1)。高光谱遥感光谱分辨率高，一般小于10nm，这使得其空间分辨率较低。较低的空间分辨率导致混合像元大量存在，纯净像元相对较少，给数据处理如精准匹配、光谱解混和目标检测与识别等技术带来了巨大的困难。

针对高光谱图像中混合像元大量存在的问题，混合像元分解技术应运而生。其目的是求解各种不同物质(即端元)在混合像元中所占的比例(丰度)。高光谱图像混合像元分解一般称为光谱解混(Spectral Unmixing)，因为高光谱图像混合像元的光谱为各个纯像元光谱以某种方式的混合。在进行混合像元分解之前需要进行端元提取，以确定各种纯像元的光谱信息，再利用各端元的光谱信息确定其在混合像元内的比例。

从理论上讲，光谱解混首先需要确定光谱混合模型，当前研究最广泛最深入的光谱混合模型是线性光谱混合模型(参考对比文件2)，然而它不能正确反映光谱的混合方式，使得解混精度不高，非线性光谱混合模型能更好的体现光谱混合方式，能得到更精确的解混结果(参考对比文件3)。

对比文件1：马艳华，高空间分辨率和高光谱分辨率遥感图像的融合[J]，红外，2003，10：11-16

对比文件2：Gautam Ghosh，Suresh Kumar and S.K.Saha，“HyperspectralSatellite Data in Mapping Salt-Affected Soils Using Linear Spectral UnmixingAnalysis”，JOURNAL OF THE INDIAN SOCIETY OF REMOTE SENSING，Volume 40，Number 1，129-136，2012

对比文件3：A.Halimi，Y.Altmann，N.Dobigeon，and J.-Y.Tourneret，“Nonlinearunmixing of hyperspectral images using a generalized bilinear model，”IEEETrans.Geosci.and Remote Sensing，49(11)：4153-4162，2011

发明内容

根据地物反射机理，发明一种基于二次散射的高光谱遥感影像非线性光谱混合模型，该模型考虑了光子在地物间发生二次散射作用，模型如下：

y＝c(Ma)+d(Ma)²+n

c+d＝1，c≥0，d≥0

$0\≤a_i\≤1,\underset{i=1}{\overset{R}{\mathrm{\Σ}}}{a}_i=1$

其中，y＝[y₁，…y_L]^T为L波段的高光谱图像，每个波段用矩阵的一行表示；M为L×R的混合矩阵，每一列为一端元(即纯地物)的光谱；R为端元的个数；a＝[a₁，…a_R]^T为丰度向量；n为误差项，一般认为0均值，方差为σ²的高斯噪声；d叫做二次散射系数。

考虑到不同混合像元的二次散射系数应该是不同的，提出二次散射模型的改进模型，模型如下：

y_ij＝c_ij(Ma)+d_ij(Ma)²+n，

0≤c_ij∈C ≤1，0≤d_ij∈D≤1，c_ij+d_ij＝1

其中，y_ij表示遥感图像第i行第j列混合像元的混合光谱向量，d_ij为该混合像元的二次散射系数。

提出一种基于丰度方差的确定二次散射系数的方法，定义的二次散射系数如下：

$d(i,j)=d_{\max }-\frac{d_{\max }\×\mathrm{var}(i,j)}{{\mathrm{var}}_{\max }},d(i,j)\∈[0,d_{\max }]$

其中，d_max是最大的二次散射系数(设为0.5)，var(i，j)为地物丰度方差，var_max为各混合像元中地物丰度方差的最大值。由此定义式可知，每个混合像元的二次散射系数为0到d_max之间。

根据该模型进行光谱解混，得到不同地物的丰度信息。提出了一种基于光谱解混结果的地质填图方法(见具体实施方式)。

附图说明

图1为两种地物二次散射情况

图2为地物分布复杂度：(a)两种地物均等分布，(b)不均等分布

图3为地质填图规则1流程图

图4为地质填图规则2流程图

图5为地质填图规则3流程图

图6为cuprite地区第100波段遥感影像

图7为USGS1995年绘制的cuprite地区地质填图

图8为二次散射模型解混结果

图9为利用规则1分别基于线性解混和非线性解混得到的地质填图。

图10为利用规则2分别基于线性解混和非线性解混得到的地质填图。

图11为利用规则3分别基于线性解混和非线性解混得到的地质填图。

具体实施方式

●模型推导：

假设某一混合像元内只存在A、B两种地物(图1)，两种地物的空间位置分布及比例都未知，设两种地物的光谱向量分别为m₁和m₂，其丰度向量分别为a₁和a₂，设M＝[m₁，m₂]，a＝[a₁，a₂]^T。

混合方式分线性和二次散射两大类：

线性情况：(a)光子入射到地物A后直接反射到传感器中；(b)光子在入射到地物B后直接反射到传感器中。

二次散射情况：光子在到达传感器前的二次撞击包含4种情况，(c)入射光子射到地物A，散射到地物A后反射到传感器；(d)入射光子射到地物A，散射到地物B后反射到传感器；(e)入射光子射到地物B，散射到地物A后反射到传感器；(f)入射光子射到地物B，散射到地物B后反射到传感器(图2)。

光子在入射到某种地物上后，将以概率c直接反射到传感器。以概率d(d＝1-c)散射到其他地物，散射的这部分光子又分别以概率a₁和概率a₂碰撞到地物A和地物B，之后反射到传感器中。该混合像元得到的光谱为：

y＝ca₁m₁+da₁m₁a₁m₁+da₁m₁a₂m₂+ca₂m₂+da₂m₂a₁m₁+da₂m₂a₂m₂+n

上式中的前6项分别对应图1中的(a)、(c)、(d)、(b)、(e)、(f)六种情况，第7项为误差项。上式可进一步化简为：

y＝ca₁m₁+ca₂m₂+2da₁m₁a₂m₂+d(a₁m₁)²+d(a₂m₂)²+n

＝c(a₁m₁+a₂m₂)+d(a₁m₁+a₂m₂)²+n

＝c(Ma)+d(Ma)²+n

●二次散射系数确定：

二次散射模型是半确定性模型，需要首先确定二次散射系数d(c＝1-d)。二次散射系数直观反映的是地物分布复杂程度，地物分布越简单二次散射系数越低，极限情况d＝0(线性混合)，相反地物分布越复杂，非线性程度越强，则二次散射系数越高。

对于遥感影像，完全靠人力实地调查确定地物复杂程度是不现实的。而只根据一幅遥感图像，缺乏必要的先验知识，很难确定混合像元地物复杂程度。二次散射情况主要发生在不同地物相交界的地方，交界地带有着最复杂的空间分布状况。因此混合像元内不同地物交界地带越多则地物越复杂，相应的二次散射系数越高。基于这种思想，可以设计确定混合模型二次散射系数的方法。

假设某混合区由两种纯地物组成，如果某个像元只存在其中一种地物，则该像元内没有不同地物间的边界，其地物分布最简单，可认为二次散射系数为0(线性情况)；假设像元内两种地物等量分布，则其交界地带最多，可认为此种情况有着最大的二次散射系数(图2(a))；如果地物不等量分布，比如地物A占10％，地物B占90％，此种情况地物间存在交界，但交界处相对较少，二次散射系数较低(图2(b))。

地物丰度的方差能表征地物复杂度情况，当像元中只有一种纯地物存在，方差最大，对应的地物复杂度最低，即二次散射系数最小(设为0)；当像元中地物分布均等时，各地物丰度方差为0，但地物分布最复杂，二次散射系数最大。可见，二次散射系数与地物丰度方差成反比。

我们定义即二次散射系数如下：

$d(i,j)=d_{\max }-\frac{d_{\max }\×\mathrm{var}(i,j)}{{\mathrm{var}}_{\max }},d(i,j)\∈[0,d_{\max }]$

其中，d_max是最大的二次散射系数(在本文实验中设为0.5)，var(i，j)为地物丰度方差，var_max为各混合像元中地物丰度方差的最大值。

●基于光谱解混的地质填图方法：

高光谱图像空间分辨率一般较低，每个像元内往往包含不止一种地物，形成混合像元，这使得传统的硬分类方法具有一定的局限性，因为无论将该像元分为哪一类都是错误的。高光谱图像光谱解混能给出图像的软分类结果，对于每个混合像元都能得到所有地物的比重。为此研究基于光谱解混的高光谱图像地质填图相较于直接分类填图更有实际意义。

基本思路是：首先对已知的高光谱图像进行波段去除等预操作，去除掉信噪比过低会严重影响解混效果的波段；然后对高光谱图像进行光谱解混，得到不同地物在混合像元内的比例分布；基于该比例分布对每个像元进行填图操作。这里由于非线性模型解混结果优于线性模型解混结果(见第四章)，因此选用基于二次散射非线性混合模型进行光谱解混。为了便于简化，对于每个混合像元单独进行填图操作，不考虑周围像元的类型。

提出了三种地质填图规则：

√填图规则1：将混合像元内比例最大的地物作为该像元的地物类型进行填图(图3)。

√填图规则2：首先对混合像元内所有地物的丰度从大到小排序，基于如下流程进行填图操作：首先判断最大比例是否大于80％，若是，则该像元标为该种地物类型，若否，则将最大比例的两种地物作为该像元类型进行填图(图4)。

√填图规则3：首先对混合像元内所有地物的丰度从大到小排序，基于如下流程进行填图操作：首先判断最大比例是否大于80％，若是，则该像元标为该种类型，若否，则将第一和第二比例相加，若大于80％，则将该像元标为此两种地物混合类型，若否，则继续加入第三种地物比例，将该像元类型标为此3种地物混合类型。由于将某一混合像元标为地物类型太多失去填图意义，这里选择最多标出3种地物类型(图5)。

仿真结果

选取美国内华达州南部Cuprite地区的AVIRIS数据，该地区富含各种矿物质，该图像大小为250*191，包含224个波段，波长为0.389～2.467微米，除去低信噪比和水吸收波段(波段1～2，104～113，148～167，221～224)，剩余188个波段数据。

图6为cuprite地区高光谱遥感影像第100波段。图7为1995年USGS绘制的地质填图。图8为基于二次散射模型解混得到的六种地物丰度图。图9为利用规则1分别基于线性解混和非线性解混得到的地质填图。图10为利用规则2分别基于线性解混和非线性解混得到的地质填图。图11为利用规则3分别基于线性解混和非线性解混得到的地质填图。

Claims (2)

1.一种基于光谱解混的地质填图方法，其特征在于，所述光谱解混是采用基于二次散射的非线性光谱混合模型；

填图规则1：将混合像元内比例最大的地物作为该像元的地物类型进行填图；

填图规则2：首先对混合像元内所有地物的丰度从大到小排序，判断最大比例是否大于80％，若是，则将该像元标为对应该比例的地物类型，若否，则将最大比例的两种地物作为该像元类型进行填图；

填图规则3：首先对混合像元内所有地物的丰度从大到小排序，首先判断最大比例是否大于80％，若是，则将该像元标为对应该比例的地物类型，若否，则将第一和第二比例相加，若大于80％，则将该像元标为对应该两种比例的地物混合类型，若否，则继续加入第三比例，将该像元标为对应该三种比例的地物混合类型；

所述非线性光谱混合模型考虑了光子在地物间发生二次散射作用，所述非线性光谱混合模型如下：

y＝c(Ma)+d(Ma)²+n

c+d＝1，c≥0，d≥0

$0\≤a_k\≤1,\underset{k=1}{\overset{R}{\Σ}}{a}_k=1$

其中，y＝[y₁，...y_L]^T为L波段的高光谱图像，每个波段用矩阵的一行表示；M为L×R的混合矩阵，每一列为一端元的光谱；R为端元的个数；a＝[a₁，...a_R]^T为丰度向量；n为误差项，是均值为0，方差为σ²的高斯噪声；d叫做二次散射系数；

当考虑到不同混合像元的二次散射系数是不同的情况下，所述非线性光谱混合模型如下：

y_ij＝c_ij(Ma)+d_ij(Ma)²+n，

0≤c_ij∈c≤1，0≤d_ij∈d≤1，c_ij+d_ij＝1

其中，y_ij表示遥感图像第i行第j列混合像元的混合光谱向量，d_ij为该混合像元的二次散射系数；

根据以下公式确定二次散射系数：

$d_{ij}=d_{max}-\frac{d_{\max }\×\mathrm{var}\left(i,j\right)}{{\mathrm{var}}_{\max }},d_{ij}\∈\[0,d_{\max }\]$

其中，d_max是最大的二次散射系数，Var(i，j)为地物丰度方差，Var_max为各混合像元中地物丰度方差的最大值；由此定义式可知，每个混合像元的二次散射系数为0到d_max之间。

2.如权利要求1所述的基于光谱解混的地质填图方法，其特征在于，d_max取值的为0.5。