基于资源一号02D卫星高光谱数据的甘肃头吊泉-南大滩地区蚀变矿物填图及金矿化带识别

孙 雨，刘家军，翟德高，柳振江，张方方，赵英俊，刘鹏飞，王子涛

(1.中国地质大学(北京)地球科学与资源学院,北京 100083；2.核工业北京地质研究院，北京 100029；3.遥感信息与图像分析技术国家级重点实验室，北京 100029；4.中国地质大学(北京)地质过程与矿产资源国家重点实验室,北京 100083)

0 引言

高光谱技术能够根据指纹性光谱特征进行地物鉴别，在岩性识别和蚀变矿物提取等地质领域取得了较好的应用效果，是当前国内外遥感应用研究的热点之一。卫星高光谱遥感由于视域宽广，一次成像即可获取大面积、价格低廉的高光谱数据，经过数据处理后能够为地质矿产勘查和基础地质调查提供重要遥感信息。搭载有高光谱传感器的卫星主要有美国EO-1、我国的GF-5和资源一号02D卫星，以往学者针对EO-1的Hyperion卫星高光谱数据处理和应用开展了大量研究工作(Hubbard and Crowley，2005；Zadeh et al.，2014；Ayoobi and Tangestani，2018；Hu et al.，2019；Souza et al.，2021)，针对我国GF-5卫星高光谱数据也开展了一些地质应用研究(董新丰等，2020；Ye et al.，2020；李娜等，2021)。当前，全世界民用在轨能够提供短波红外谱段卫星高光谱数据的只有我国的资源一号02D卫星，有关这一新型高光谱数据的地质应用研究还较少(李娜等，2020；魏红艳等，2020；李根军等，2021)，迫切需要在数据处理方法和应用效果方面加强研究交流，推动“天上卫星好用”和“地上数据用好”两方面更好的协调发展。

本文应用甘肃省瓜州县头吊泉-南大滩地区的资源一号02D高光谱数据，探索了该高光谱数据的地质应用前景。采用改进的SAM填图算法提取该地区的褐铁矿、绿泥石、方解石、短波云母等7种蚀变矿物，分析蚀变矿物的地质分布特征，并对新发现的金矿化蚀变带进行野外查证，证实存在地表金矿化异常，为该区找矿勘查提供了重要线索。

1 研究区地质背景

图1 甘肃头吊泉-南大滩地区地质简图

2 卫星高光谱数据源

2019年9月12日11时26分，资源一号02D卫星在我国太原卫星发射中心成功发射升空进入预定轨道。作为我国自主建造并成功运行的首颗民用高光谱业务卫星，该卫星是国家民用空间基础设施中新型对地观测卫星发展的又一重要成果。该星可有效获取宽幅高光谱等遥感数据，将与后续系列卫星组网形成全球领先的业务化对地光谱探测能力，为自然资源领域的卫星高光谱数据业务化应用提供了重要的数据保障(杨诗瑞，2019)。

资源一号02D卫星轨道高度为778 km，搭载的高光谱传感器载荷可获取幅宽60 km、空间分辨率优于30 m的高光谱(AHSI)数据。资源一号02D高光谱数据谱段共166个，其中可见近红外谱段光谱分辨率优于10 nm，短波红外谱段光谱分辨率优于20 nm。

本次使用的是1景2020年10月12日获取的资源一号02D AHSI高光谱遥感数据，产品级别为L1A级，产品号为L1A0000179333，包含GeoTiff数据文件、轨道参数RPB文件及定标及光谱响应RAW文件等。头吊泉-南大滩地区位于整景高光谱数据的西北侧(图2)，面积约450 km2。

图2 头吊泉-南大滩地区高光谱影像位置图

3 高光谱数据处理方法

3.1 高光谱数据预处理

资源一号02D卫星的L1A级高光谱原始数据为带RPC校正参数的DN值数据，需要进行数据预处理才能用于蚀变矿物填图，数据预处理由基础预处理和大气校正两部分构成。在基础预处理方面，应用当前热门的python编程语言，开发了具有自主知识产权的资源一号02D高光谱数据批量基础预处理软件模块，具备解压缩、正射校正、辐射定标、波段组合、BSQ转BIL、波段去除功能，实现一键化对指定文件夹下所有高光谱原始数据进行批量的基础预处理。软件模块的主要功能函数如图3所示。

图3 高光谱数据基础预处理模块主要功能函数

批量基础预处理软件模块的处理思路是首先应用ASTER GDEM V3高程数据进行正射校正，得到具有高斯克吕格-CGCS2000投影坐标系统的高光谱数据。然后，应用原始数据自带的定标参数进行辐射定标，得到具有物理意义的辐射亮度数据，并将VNIR和SWIR数据文件合并成具有166波段的一个数据文件，再将BSQ格式转换为BIL格式。随后，去除高光谱数据的水汽吸收扰动明显的波段和目视效果差的传感器起止波段(共24个)，最终得到142个波段的高光谱数据。

资源一号02D高光谱数据的单景tar压缩包大小为1.31 G，应用本次开发高光谱预处理模块的处理时长为121 s，使用ENVI 5.3商业软件的处理时长为579 s，本次开发的模块在耗时上约为商业软件的1/5，可有效助力于区域大面积的资源一号02D高光谱数据快速预处理。

在采用开发的算法模块进行基础预处理后，再设置好各项参数后应用ENVI 5.3中FLAASH算法模块进行大气校正，得到可用于矿物填图的高光谱反射率数据。

首先，是对学生进行分层。根据学生之前的学习成绩、班主任建议等因素，综合评估，把学生分为A班，B班，原则上比例为1∶1。A班侧重提高，B班侧重基础。学生分层工作中最重要的是心理辅导，要向学生宣传正确的教育理念，疏导学生心理问题，避免学生出现自卑等负面情绪。

3.2 高光谱数据降维增值

高光谱数据拥有上百个波段，存在两个处理关键问题。一是，如果将142个波段全部进行计算将极大降低处理效率，而且蚀变矿物并不是在所有波段都具有明显吸收特征；二是，在蚀变矿物反映敏感的短波红外谱段(燕守勋等，2003)，其光谱采样间隔在16 nm左右，与GF-5高光谱数据的8 nm相比光谱分辨率有所降低，需要考虑在光谱维度进行升采样以提高光谱分辨率从而可以更加精细的识别矿物。

针对第一个问题，本文采取波段筛选组合的高光谱数据降维处理方法，有针对性地对蚀变矿物有明显吸收反射特征的波段进行挑选后组合，使矿物特征更加集中(郑向涛，2017)。对于VNIR高光谱数据，选择499～997 nm谱段；对于SWIR高光谱数据，选择2014～2451 nm谱段。

三次样条函数作为数学分析中常用的插值方法，广泛应用于离散点插值计算(张成良和刘小泉，2006；刘翔舸等，2010；殷瑞娟等，2013)。杨萍等(2007)应用三次样条插值法根据8波段数值求解出106个波段(间隔5nm)的光谱响应值，认为该方法可以非常近似地模拟出光滑的地物光谱曲线。梁晨(2016)对同型号不同分辨率仪器实测数据进行了三次样条插值，解决了光谱测量点的匹配问题。因此，本文针对第二个问题，引入了三次样条函数对短波红外谱段在相邻波段间进行插值运算，实现不破坏光谱固有信息的光谱维增值，将光谱采样间隔由16 nm升采样到了2 nm。

3.3 蚀变填图端元选取

本文应用降维增值后的高光谱数据，使用ENVI 5.3软件的光谱沙漏(Spectral Hourglass)流程算法(李光辉等，2013)，得到蚀变矿物在VNIR谱段和SWIR谱段的图像候选端元光谱。然后，基于专家光谱知识对候选端元光谱进行优选，得到VNIR谱段褐铁矿的图像端元光谱(图4a)和SWIR谱段绿泥石、方解石、白云石、短波云母、中短波云母、中长波云母的图像端元光谱(图4b)。为了进一步凸显不同蚀变矿物的吸收和反射特征，对SWIR谱段图像候选端元光谱进行包络线去除，得到该谱段用于蚀变矿物填图的图像端元光谱(图5)。其中，短波云母端元第一特征吸收峰在2199 nm附近，中短波云母端元光谱第一特征吸收峰在2206 nm附近，中长波云母端元光谱第一特征吸收峰在2216 nm附近。

图4 优选后的蚀变矿物图像端元光谱Fig.4 Selected image endmember spectra of alteration minerals

图5 SWIR谱段数据包络线去除后的蚀变矿物图像端元光谱Fig.5 Image endmember spectra of alteration minerals in SWIR data after continuum removal

3.4 蚀变矿物填图与制图

在VNIR谱段，褐铁矿等蚀变矿物的光谱吸收峰较为宽缓。针对这一特点，将优选后的图像端元光谱作为参考光谱，以各个像元的图像光谱为待匹配光谱。在SWIR谱段，云母等蚀变矿物的光谱吸收峰较为尖锐，而且往往存在不止一个光谱吸收峰，因此将包络线去除后的图像端元光谱作为参考光谱，以各个像元包络线去除后的图像光谱为待匹配光谱。在选择好参考光谱和待匹配光谱后，采用本文建立的改进的光谱角填图方法进行蚀变矿物填图。

改进的光谱角填图方法是针对传统光谱角填图(Spectral Angle Mapping)方法仅度量全局谱形相似度的缺陷，在对SAM方法增加一系列特征吸收-反射峰位判别条件后，逐像元计算待匹配光谱与参考光谱的相似度。本文应用当前流行的python程序开发语言，开发了改进的光谱角填图算法，可以批量进行高光谱蚀变矿物填图。

高光谱蚀变矿物填图计算完成后，得到矿物匹配得分灰度图。像元灰度值为弧度值，灰度值越小则与参考光谱匹配程度越高，进行低值端阈值切割并填充对应的颜色，得到各种单矿物分布图，并以空间坐标信息为链接，将所有蚀变矿物全部叠加到岩性-构造地质图上，在添加比例尺等图饰信息后，编制头吊泉-南大滩高光谱遥感地质综合图(图6)。然后，应用成像时间为2021年1月29日的亚米级(0.8 m)空间分辨率的GF-2国产遥感数据制作了真彩色高分影像图，在叠加了褐铁矿和中短波云母蚀变矿物后，编制了南大滩-西红滩地段蚀变高清分布图(图7)。

图6 头吊泉-南大滩地区高光谱遥感地质综合图Fig.6 Comprehensive map showing hyperspectral geological information of the Toudiaoquan-Nandatan area

图7 南大滩-西红滩地段金矿化蚀变带分布图Fig.7 Distribution of gold mineralization alteration belt in the Nandatan-Xihongtan section

4 结果与讨论

4.1 蚀变矿物分布特征

4.2 金矿化蚀变带展布特征

许多类型的矿产都与热液蚀变关系密切，是开展矿产资源勘查的有利信号(朱永峰和安芳，2010，杜斌等，2021)。蚀变矿物是成矿热液蚀变活动的直观证据，以往多光谱遥感技术提取的羟基等矿物族级别的蚀变矿物为金属矿产勘查提供了重要的遥感判据(王钦军等，2017；Yousefi et. al，2018；刘小雨等，2020)，高光谱遥感更是提供了传统地质方法难以获取的无缝像元尺度的蚀变单矿物信息(孙雨等，2015；连琛芹等，2020；吴德文等，2021)，为金矿、铜矿、铀矿、钨钼矿等金属矿种的找矿勘查圈定了范围更小的找矿远景区(汪重午等，2014；贺金鑫等，2017；刘德长等，2017；叶发旺等，2021)。

众多研究表明，金矿化与硅化、黄铁矿化、褐铁矿化、绢云母化、绿泥石化等蚀变关系密切。范潇等(2015)对焦家金矿田进行了构造蚀变测量，蚀变岩带的地球化学测试表明绢英岩化蚀变阶段是矿化元素富集的主要阶段。矿物填图中使用的“云母”系指具有不同结构和粒度的绢云母、白云母族矿物，通常统称为“绢云母”。“绢云母”的二八面体晶体结构主要受Tschemark替代过程影响，从而造成矿物结构中的Al-OH吸收峰位置发生漂移，一般介于2185～2225 nm之间。邵雪维等(2021)对胶东新城金矿深部蚀变岩进行了光谱测量，研究结果表明金矿化岩心段的绢云母Al-OH吸收峰位置在2205～2208 nm。其他矿种围岩蚀变中的绢云母矿物Al-OH吸收峰偏移也有类似的规律。如郭娜等(2012)对西藏甲玛斑岩铜多金属矿进行了蚀变矿物光谱测量，结果表明矿体主要与Al-OH吸收峰位置处于2205～2208 nm的绢云母共存。周岩等(2019)对福建紫金山铜金矿区进行了岩心扫描，发现与铜钼矿体关系密切的绢云母Al-OH吸收峰位置大于2205 nm，以此找矿标志发现了肉眼难以识别的新矿体。徐清俊等(2016)对新疆白杨河铀矿进行了全孔岩心光谱测量与分析，证实铀矿化与赤铁矿化、高铝绢云母(Al-OH吸收峰位置在2192～2203 nm)和中铝绢云母(Al-OH吸收峰位置在2204～2212 nm)有关，三者共同存在的地段是铀成矿有利地段。

金矿化普遍与黄铁绢英岩化关系密切，即与黄铁矿(以及进一步蚀变形成的褐铁矿)、绢云母、石英矿物密不可分。在谱段范围400～2500 nm的资源一号02D高光谱数据上，代表硅化的石英和代表黄铁矿化的黄铁矿均不具有明显光谱特征，无法直接探测识别出来；褐铁矿、绢云母等蚀变矿物具有光谱特征，可以直接探测识别出来。因此，金矿化在高光谱蚀变图上可表现为褐铁矿+绢云母。在详细分析高光谱蚀变矿物分布图结合地质岩性背景，在研究区南部识别出了1条可能的金矿化蚀变带。该矿化蚀变带分布于南大滩-西红滩地段，发育褐铁矿和中短波云母蚀变矿物(图7)。蚀变带总体近东西向展布，断续延伸约15 km，所处地质体为海西中期及晚期花岗岩。

4.3 金矿化蚀变带验证情况

针对高光谱技术在南大滩-西红滩地段识别出的金矿化蚀变带，综合采用野外地质踏勘、ASD光谱测量和地球化学样品采集及分析测试等技术手段，验证了该带地表存在金矿化现象。

野外地质踏勘发现该蚀变带位于第四系大滩的北侧，蚀变原岩为褐红色中粗粒块状花岗岩，发育石英脉，可见明显的硅化及褐铁矿化现象(图8a、8b)，绢云母化难以肉眼辨识。通过对蚀变岩的ASD光谱测量结果分析，显示同时存在褐铁矿和中短波云母矿物(图8c)，褐铁矿体现在实测光谱曲线中880 nm附近的宽缓吸收峰和765 nm附近的反射峰，短波云母体现在实测光谱曲线中2304 nm附近的第一吸收峰、2350 nm附近的第二吸收峰和2445 nm附近的第三吸收峰，不仅与USGS标准光谱库中褐铁矿、绢云母矿物的标准光谱曲线特征吻合，也与用于高光谱蚀变矿物提取的图像端元光谱特征吻合较好，证实高光谱提取的蚀变矿物确实存在，高光谱填图结果可信。

图8 南大滩-西红滩地段金矿化蚀变带综合验证结果Fig.8 Comprehensive verification results for alteration belt of gold mineralization at Nandatan-Xihongtan section

在南大滩-西红滩地段东部的西红滩地区(图7)开展了地表捡块地球化学样品采集和室内分析测试，测试结果(图8d)显示6个样品的金元素含量高于0.3×10-6，最高值可达1.839×10-6，这还仅是地表岩石的金元素含量，深部岩石的金含量还会增高。多数样品的二氧化硅含量高于80%甚至超过90%，远超花岗岩平均值(66%)，表明蚀变岩中存在较强的硅化现象。总FeO含量也在2%上下浮动，为含铁矿物黄铁矿、褐铁矿的显示。此外，银元素也存在高值异常。值得关注的是，南大滩-西红滩地段的西侧还存在相似的岩性-蚀变地段(图7)，后续将开展相关矿化查证工作。

5 结论

(1)针对资源一号02D遥感数据特点开展了高光谱数据预处理、降维增值和SAM光谱匹配方法改进，开发了基于python语言的资源一号02D高光谱数据批量基础预处理、改进的SAM匹配的软件模块。在甘肃头吊泉-南大滩地区进行了高光谱蚀变矿物填图，成功提取了褐铁矿、绿泥石、方解石、白云石、短波云母、中短波云母、中长波云母共7种蚀变矿物。

(2)资源一号02D高光谱数据在该地区提取的蚀变矿物分布具有明显规律性。一是蚀变矿物与地质体岩性关系密切，多数与断裂构造并无明显直接联系；二是蚀变矿物并不是均匀分布在地质体中，而是集中在特定地质体的局部地段。方解石和白云石碳酸盐类蚀变矿物的分布区可以指示碳酸盐类沉积/变质地层的出露范围，蚀变的定向展布揭示了地层的走向。

(3)根据资源一号02D高光谱数据提取的蚀变矿物，结合地质背景信息分析，新发现1条长约15 km的金矿化蚀变带。东部的西红滩地段经野外地质调查、光谱测量、地化采样测试等综合验证存在地表金矿化异常，西部地段有待后续开展相关查证工作。资源一号02D高光谱遥感数据为地质找矿提供了成功案例，具有较好的地质深入应用潜力。

致谢：感谢中国资源卫星应用中心提供的资源一号02D高光谱数据，感谢论文评审专家对本文提出的宝贵意见。