基于遥感指数的岩性信息提取研究

汪清川，边 宇，孙 昂，周 冲，吕大庆，李向英，奚砚涛※

（1.中国矿业大学资源与地球科学学院，江苏 徐州 221116；2.中国地质调查局自然资源航空物探遥感中心，北京 100083；3.珠海欧比特宇航科技股份有限公司，广东 珠海 519080）

近年来利用多光谱遥感影像提取地表的一些岩性信息逐渐成为了遥感和地质研究的热点，国内外的学者利用遥感数据进行大量的有关岩性识别的研究[1－7]。由于不同的岩矿物质在遥感影像上呈现出不同的波谱曲线，因此可以基于不同岩矿石的理化性质及其成分差异，提取它们的遥感波谱信息，有针对性地构建相应的遥感指数，从而对相关岩矿石的岩性信息进行识别[8－10]。早期的遥感数据由于其波段数量较少，且光谱分辨率较低，在岩性信息提取研究中有一定的局限性；但近年来遥感影像的分辨率逐步提高，在岩矿物质识别精度上有了很大的提高。本文中用到的ASTER 传感器，有3 个独立的子系统，分别处于近红外、短波红外、热红外波段，共计14 个波段，其波段设置比其它多光谱数据具有更好的连续性，光谱分辨率也有很大的提高，且其传感器是专门为地质应用和火山监测而设计，在国内外岩性信息提取中有着广泛的应用，国内许多学者已通过ASTER 数据进行相关的地质岩性填图工作[11－14]。国外学者Rowan 等人基于ASTER 遥感数据对一些重要的岩石单元进行识别与提取，成功将方解石质岩石与白云石质岩石区分开，将侵入岩中常见的铁质白云母与花岗质片麻岩及中生代花岗岩中所含的铝质白云母分离出来[15]。Ninomiya等人针对ASTER 的TIR 谱域定义了CI、QI 和MI矿物比值指数，分别用于识别碳酸盐类岩石、石英质沉积岩和硅酸盐类岩石[16]。国内学者甘甫平等提出了基于特征谱带的高光谱遥感矿物谱系识别方法，初步建立了“矿物大类－族－种－亚种”的矿物识别分层谱系[17]。王润生等对岩石矿物的光谱变异特征和岩石矿物的混合光谱特征进行了较系统的研究[18]。这些研究为进一步的地质矿产信息提取奠定了理论基础[19－20]。

昆仑造山带为青藏高原北部的重要地质单元，经历了漫长的地质演化过程和多次造山运动，形成了典型的复合造山带。研究区域布喀达坂峰是昆仑山脉中段的最高峰，该区域地质条件复杂，气候恶劣，山势险峻，对该区域的地质调查与研究受制于恶劣的自然条件很难进行。本文以ASTER 数据对该地区进行地表岩性信息的提取，取得了良好的效果。

1 研究区域和数据源

1.1 研究区域地质概况

图1 研究区域遥感影像图（R（b8）、G（b3）、B（b1））Fig.1 Remote sensing image of the study area（R（b8）、G（b3）、B（b1））

研究区位于昆仑山中段阿尔格山东端与博卡雷克塔山西头交接处，在青海省玉树藏族自治州治多县境内，为新疆、青海的界山，地处东经90.9°，北纬36.0°，是昆仑山脉中段的最高峰，海拔达6860 m，耸立于群峰之上，雪原绵延，气势磅礴（图1）。布喀达坂峰气候干燥，峰区有巨大的冰帽冰川，全年最低温度达－30℃，山势险峻，春秋冬季风大。区内出露的主要岩石地层有：早古生代纳赤台群哈拉巴依沟组，以浅变质岩屑长石（杂）砂岩、长石岩屑（杂）砂岩为主夹粉砂质板岩、千枚岩，局部夹玄武岩；泥盆纪布拉克巴什组主要为一套海陆交互相陆源碎屑岩及一套海相复理石夹海相中基性火山岩的岩石组合；二叠纪浩特洛哇组以碎屑岩为主夹少量安山岩及生物碎屑灰岩，树维门科组下段岩性以浅灰白色亮晶、微晶生物碎屑灰岩及亮晶核形石灰岩为主，上段岩性以生物微屑岩为主；二叠纪马尔争组下段为一套碎屑岩、碳酸盐岩夹火山岩的地层，中段由一套砂板岩复理石组成，岩石种类较单一，以灰色中细粒岩屑长石砂岩、岩屑长石粉砂岩为主；中新世地层沉积类型繁多，分布广泛，包括冲积、冲洪积、冰碛、风积、湖积、沼泽堆积等。

1.2 遥感数据

ASTER 传感器是一个高分辨率多光谱成像仪，共有三个子系统，包括了从可见光到热红外3个谱段，其中可见光－近红外共3 个波段，其空间分辨率为15 m；短波红外有6 个波段，空间分辨率为30 m；近红外有6 个波段，空间分辨率为90 m。其具体参数如表1 所示。本文中利用的ASTERL1T 遥感影像数据，获取时间为2015 年4 月份，来源于美国地质勘探局官网，该时间图像云量遮盖较小，适合做岩性信息提取分析。

表1 ASTER数据相关参数Table 1 Related parameters of ASTER data

2 岩性信息提取

2.1 ASTER 数据预处理

2.1.1 串扰校正

由于ASTERR 传感器的自身缺陷，导致其频段4 的信号泄漏到相邻的频段5 和频段9，使得B5 和B9 波段反射率比实际反射率高，从而影响岩性信息提取的精度。因此，要对ASTER 遥感影像数据进行串扰校正，提高解译精度。

2.1.2 辐射定标与大气校正

在进行大气校正之前，需要对原始数据进行辐射定标。由于ASTER 遥感影像数据是以无量纲的数字量化值（DN 值）记录信息的，通过辐射定标可以把DN 值与辐射亮度值，反射率值和温度值等物理量进行转化。ENVI 可以从元数据文件中获取每个波段的Gain 和Offset 参数，自动进行定标。辐射定标完成后，为了充分利用短波红外波段的高分辨率特性，将热红外与短波红外波段进行波段融合，将分辨率统一至30 m。最后选择大气校正方式，本文选择FALASH 模块对数据进行大气校正，来消除大气中水蒸气，氧气，二氧化碳等气体对地物反射的影响。

2.1.3 干扰信息掩膜

干扰信息掩膜主要是为了消除影像上雪、植被、云层等对岩性信息提取时的干扰，由于它们的灰度值较高，分布区域较大，直接参与运算时，容易出现异常值，创建目标掩膜可以有效扣除这些干扰信息。在创建目标掩膜时，先基于归一化NDVI 指数提取植被信息，在NDVI 灰度图中亮度值高的白色区域就是植被覆盖的区域，然后建立掩膜，确定阈值，判定DN 值大于0.25 为植被。如图2 为植被掩膜结果，再根据比值band3／band4 增强雪体信息，使用与植被同样的方法对积雪区进行掩膜处理，通过逐层掩膜的方法来剔除植被和积雪等信息，最终掩膜结果如图3 所示。

使用掩膜处理一方面可以减少计算量，使有用信息得到增强，节省数据处理时间；另一方面可以有选择地去除原始影像中的干扰信息，提高岩性信息提取精度。

2.2 岩性信息识别方法及过程

图2 植被掩膜结果Fig.2 Results of vegetation mask

图3 植被和雪体掩膜结果Fig.3 Vegetation and snow mask results

岩石是由一种或多种矿物组成的固体聚集体，因此岩石的光谱是多种矿物的混合光谱。ASTER的VNIR（波长范围0.52 μm 至0.8 μm）3 个波段具有有关过渡金属吸收的重要信息来源，且其第二波段是叶绿素的主要吸收波段。SWIR（1.60 μm 至2.43 m）的六个谱带能够显示许多硅酸盐、碳酸盐、水合物和氢氧化物矿物的分子吸收特征。此外，ASTER 的五个TIR 光谱带还可用于计算表面温度和发射率，以识别提取重要的成岩矿物，例如石英和长石，并在TIR 波长区域显示基本的分子吸收特征。

2.2.1 碳酸盐，碳酸盐矿物在VNIR－SWIR 谱域矿物光谱特征

碳酸盐岩主要由碳酸盐矿物方解石、白云石、菱镁矿、菱铁矿等组成。碳酸盐矿物是由碳酸根离子与金属阳离子构成的化合物，其光谱特征的成因主要是金属阳离子的电子跃迁和CO32－的振动。碳酸盐矿物在可见光－近红外波段和短波红外波段都表现出较为明显的CO32－吸收特征。如图4（a）所示，方解石的吸收集中在2.33 μm 和2.34 μm 之间，白云石的吸收在2.31 μm 和2.3 μm 之间。但受到金属阳离子的影响，不同碳酸盐矿物的吸收中心波长位置存在差异。方解石的吸收中心为2.3465 μm，白云石的吸收中心为2.3039 μm，表明方解石和白云石可以通过在2.33 μm 和2.45 m 之间的吸收变化来区分和识别。

硅酸盐根据硅氧四面体的联结方式可分为5种类型，分别是岛状硅酸盐、链状硅酸盐、层状硅酸盐、架状硅酸盐和环状硅酸盐。图4（b）中列举出了4 种硅酸盐代表矿物在VNIR－SWIR 范围的光谱曲线，从图中光谱特征曲线上可以看出，角闪石2.3 μm和2.4 μm 附近出现两个吸收峰，原因是Mg－OH 的伸缩振动。粘土类矿物高岭石和白云母，高岭石有两个吸收峰，分别在2.16－2.17 μm 和2.2 μm 处，前者吸收强度较弱，后者较强，原因是Al－OH 基团伸缩振动；白云母在2.20 μm、2.35 μm、2.45 μm 处吸收特征明显，原因是Al－OH 基团振动。橄榄岩在VNIR 谱域的吸收特征是铁离子的电子跃迁引起的，橄榄石1.0 μm 附近出现强吸收带，主要是Fe2＋跃迁作用产生的，Fe3＋作用特征相对较弱。

2.2.2 碳酸盐，碳酸盐及长石矿物在TIR 谱域矿物光谱特征

岩石在热红外光谱中的发射光谱特征主要与硅酸盐矿物和碳酸盐矿物的选择性发射有关，因此使用热红外波段数据对碳酸盐岩矿物和硅酸盐矿物的提取效果显著。如图5（a），石英在8－10 μm 波段表现典型的双低发射率谱带，分别对应ASTER第10 波段和第12 波段。如图5（b）所示，长石矿物的低发射率特征谱带位置集中在8μm 至11μm，例如正长石、条纹长石、微斜长石等碱性长石矿物，更长石和钠长石等斜长石矿物，该谱域低发射率是Si－O 基团伸缩振动产生的，谱带中心具有随着SiO2含量的降低与铁镁质矿物含量的增加而向长波方向移动的特点。如图5（c），碳酸盐矿物方解石和白云石在TIR 波长区域具有典型发射率特征，分别为11.40 μm 和11.35 μm。主要原因是碳酸根离子基振动。

图4 碳酸盐矿物（a）和硅酸盐矿物（b）VNIR－SWIR谱域光谱特征曲线Fig.4 VNIR－SWIR spectral characteristic curve of carbonate minerals（a）and silicate minerals（b）

2.2.3 基于遥感指数的岩性提取

岩性信息识别的方法很多，本文采用波段比值法进行岩性的识别提取。ASTER 波段比值又称为矿物指数，是一种经常用来提取岩性信息的方法，通常选择目标矿物光谱的最小反射波段或最大反射波段作为比值波段进行代数运算，构建不同的遥感指数；当波段间差值相近但斜率不同时，增强不同岩性之间的微小差异，来更好的识别某种岩性信息。

（1）碳酸盐指数（SWIR 范围），计算公式如下：

分析碳酸盐矿物在SWIR 波段范围内的反射光谱特征曲线，可以看出碳酸盐矿物（如方解石）在SWIR 谱域因碳酸根离子振动，在2.3 μm 附近有显著的吸收特征，出现吸收谷，对应ASTER 波段8。波段6 和波段9 的反射率都高于波段8，所以选择它们建立方解石指数（b6／b8）×（b9／b8）来增大方解石与其它岩性的差异，来突出该类岩石的信息。同时，指数（b6＋b9）／（b7＋b8）可用于识别含绿泥石、绿帘石或者角闪石类矿物，也能指示含碳酸盐矿物的岩体的分布。比值指数（b7＋b9）／b8 也能用于碳酸盐岩的提取。三个指数的提取效果接近，得出指数图像十分相像。图6 为CLI 灰度图像，可以看到二叠纪树维门科组上段像元亮度值较高，该区域方解石含量高。

（2）碳酸盐指数（CI），计算公式如下：

图5 硅酸盐矿物（a）、长石矿物（b）、碳酸盐矿物（c）TIR谱域光谱特征曲线Fig.5 TIR spectral characteristic curves of Silicate minerals（a），feldspar minerals（b）and carbonate minerals（c）

图6 CLI灰度图像Fig.6 Grayscale image of CLI

通过TIR 发射光谱特征曲线确定了碳酸盐类矿物方解石和白云石在TIR 波长区域的吸收位置分别为11.40 μm 和11.35 μm，对应ASTER 第14波段；方解石和白云石在TIR 波长区域ASTER13波段有最高发射率，比值b13／b14 能有效增强这些碳酸盐矿物的信息，有效地突出含方解石、白云石等碳酸盐矿物的碳酸盐岩。图7 为CI 灰度图像。

（3）改良镁铁质指数（MI），计算公式如下：

图7 CI灰度图像Fig.7 Grayscale image of CI

镁铁质矿物即镁铁硅酸盐矿物，是氧化镁与氧化亚铁和二氧化硅结合形成的，比如橄榄石、辉石等。由于氧化镁与氧化亚铁与硅酸岩中的二氧化硅含量呈负关系，因此只有在二氧化硅含量低的情况下，才会出现镁铁质矿物，所以在MI 指数灰度图8上，高亮地区表示二氧化硅含量低，从而反映研究区硅酸盐信息。

（4）石英指数（QI），计算公式如下：

分析石英的发射光谱特征曲线，石英在8－10 μm 波段表现典型的双低发射率谱带，分别对应ASTER 第10 波段和第12 波段，分别与高发射率作比值处理建立石英指数能有效突出石英含量高的岩石信息；且因为与长石矿物的低发射率相反，所以石英指数能反映富含石英而缺乏长石的岩石。在QI 指数灰度图9 中，像元亮度值高说明该区域富含石英。

（5）碱性长石指数（AFI），计算公式如下：

根据长石矿物光谱曲线特征分析得出，正长石、条纹长石等碱性长石band12 出现低发射率谱带，band13 的发射率要高于band12，将band13 与band12 做比值运算，能增强碱性长石信息。因此，构造碱性长石指数AFI 的，可以突出碱性长石信息。图10 反映了研究区域内的AFI 灰度图像。

（6）斜长石指数（PLI），计算公式如下：

图8 MI灰度图像Fig.8 MI Grayscale image

图9 QI灰度图像Fig.9 Grayscale image of QI

图10 AFI灰度图像Fig.10 Grayscale image of AFI

根据长石矿物光谱曲线特征分析得出，斜长石矿物的低发射率波长位置与碱性长石相反，band11和band13 处形成吸收谷，而band12 的发射率值高于band11 与band13，所以比值band12／ band11 和band12／band13 都能一定程度增强斜长石的信息；因而，选择比值b12×b12 ／（b11×b13）构造PLI，很大程度地增强这些长石矿物信息。图11 为研究区域内的PLI 比值灰度图像。

3 岩性信息提取结果分析

本文通过将三种矿物指数分别置于红（R）绿（G）蓝（B）分量上，进行RGB 假彩色合成图像，然后通过合成图像色彩上的差别来识别造岩矿物，进而目视解译出岩性。图12 是由R（方解石矿物指数CLI）、G（指数（b7＋b9）／b8）、B（指数（b6＋b9）／（b7＋b8））假彩色合成的。从图中可以发现，指数图像颜色单调，整体呈暗色调，只有一个区域出现了明显的白色；按颜色组合进行解译，单调的白色区域碳酸盐特征明显，从矿物指数图像上可以看出比较清晰的界线。对比研究区地质图，该区发育有生物碎屑灰岩，界线与指数图像基本吻合。黑色区域是植被覆盖区域和雪地掩膜处理后表现，其余区域颜色单一，难以目视解译识别划分岩性。

图11 PLI灰度图像Fig.11 Grayscale image of PLI

图12 光谱指数假彩色合成影像Fig.12 False color composite image of spectral index

图13 是由R（石英指数QI）、G（碳酸盐指数CI）和B（铁镁质指数MI）彩色合成的矿物指数图像。为了便于解译分析描述，将提取结果划分多个区域，标注上字母代表对应区域。图中可以发现，A区东昆仑南活动陆缘带有条带状的蓝紫色区域，含石英矿物和镁铁矿物，应为硅酸盐矿物；该区域有雪和植被覆盖，岩性难以识别提取，根据地质关系推测覆盖区域与蓝紫色区域岩石类型相同。CI 碳酸盐指数对应绿色波段，所以N、K、L、R 绿色区域指示的为含碳酸盐矿物，该处岩石的白云石和方解石矿物含量较丰富。矿物指数图像上能够看得出清晰的岩性界线，对比研究区地质图，该区域发育有生物碎屑灰岩，界线提取结果比较吻合。QI 石英指数对应红色波段，以红色为主色调的F、J、M 区域应为高石英含量的岩石，应为石英质沉积岩；对比地质图，该处地层岩性主要为碎屑岩。H 区域下方的紫色区域与其同属一个地层，而指数图像上显示紫色，原因是该区岩石含有铁镁质矿物。C 区颜色斑杂，为杂色碎屑岩。B、D 区界线分明，参照提取出的岩性，B 区为岩屑长石砂岩，D 区与A 区同为花岗岩。

图14 是由R（碱性长石指数AFI）、G（斜长石指数PLI）和B（石英指数QI）彩色合成的矿物指数图像，从图中可以看出，图像主体主要呈现紫色和绿色，绿色区域岩性主要为沉积岩与第四纪的沉积物，紫色区域岩体主要矿物为碱性长石和石英，对比已有研究区的地质图，B 区为花岗岩，A 区为泥盆纪碎屑岩，两种岩石主要组成矿物都为长石和石英。

不同的矿物指数提取的侧重点不同，能突出不同矿物丰度信息，综合上述几个矿物指数图像，对比已有地质图绘制了研究区的岩性分布解译图（图15）。

4 结论

图13 光谱指数假彩色合成影像（R（QI）、G（CI）、B（MI））Fig.13 False color composite image with spectral index（R（QI），G（CI），B（MI））

图14 光谱指数假彩色合成影像（（AFI）、G（PLI）、B（QI））Fig.14 False color composite image of spectral index（（AFI），G（PLI），B（QI））

图15 研究区岩性解译图Fig.15 Lithologic interpretation map of the study area

本文利用ASTER 数据，研究分析不同矿物岩石的光谱特征曲线，通过构建不同矿物的遥感指数对布喀达坂峰地区的矿物岩石信息进行提取，再将不同矿物指数的灰度图像进行RGB 假彩色合成，对研究区的岩性信息进行解译，划分岩石类型。该方法很好地提取了相关岩石的岩性信息，特别是对灰岩的识别提取效果较好。通过目视解译可以在遥感图像上勾勒出不同岩性地质单元，从已有的地质图和野外验证来看，可以发现该方法提取出的岩性界线与地质图上的地质界线较为吻合，所解译出来的岩石地层分布与组级岩石地层单位基本匹配。虽然采用矿物指数方法提取裸岩层区域的岩性信息会取得比较好的效果，但是受ASTER 数据的空间分辨率和岩石露头宽度的影响，目前尚无法识别单位或者厚度小于ASTER 空间分辨率的单一岩性信息。另外在信息提取的过程中对植被和雪做了掩膜处理，当使用遥感指数法提取岩性信息时，某些矿物的信息会发生混淆，对遥感岩性解译造成影响。所以该方法也存在一定的局限性，在应用的时候应根据具体情况进行相应的调整和选择。