徐 越,张 震,朱 骏,谢 毅,王璇璇
(1.安徽理工大学 测绘学院,安徽 淮南 232001;2.浙江水利水电学院 测绘学院,浙江 杭州 310000)
遥感矿化蚀变信息的提取是遥感找矿的重要步骤,对于矿产的勘测与评价具有重要指示意义[1],如何从遥感影像中提取出更为准确的蚀变信息一直是遥感找矿中的关键问题[2]。近些年来,Landsat8 OLI、WorldView-2等影像在蚀变信息提取中得到了广泛应用,OLI数据在继承TM/ETM数据特点的基础上,增大了波段数量和波段光谱范围,提高了影像的辐射分辨率,增强了对岩性的识别[3],Landsat8在VNIR有5个波段,更适用于铁染蚀变信息的提取[4]。WorldView-2影像数据除了拥有较高的分辨率和丰富的空间信息外,其较为丰富的光谱波段与光谱信息将更加有利于信息的提取和遥感制图能力[5],使地物的识别能力比多光谱数据强,但存在波段宽度窄的问题,使地物识别能力受到限制。将两种不同的影像数据融合能使光谱信息更为精确,每种数据能在新形成的影像数据中发挥其优势,来提高最终蚀变信息提取的准确性。王守志等人[6]将Landsat8 OLI的多光谱波段和GF-1 PMS的全色波段进行数据融合,采用主成分分析法对研究区的铁染蚀变信息进行提取,并通过对比前后的提取结果来验证融合后的数据在蚀变信息提取中的优势。将两种影像进行协同作为另一种有效的方法,也能够使光谱信息更为准确。协同与融合不同,不是将光谱波段进行融合,而是将所有的光谱波段拼接起来,使光谱数大大增加,提高光谱的详细程度,对于研究提取蚀变信息有着重要的意义。唐超等人[7]采用这种方法将Landsat8 OLI数据和ASTER数据进行集成后,为提高集成数据的空间分辨率又与Landsat8数据的全色波段进行融合,对研究区进行铁染和羟基异常的提取,发现经集成融合后的数据能够提取出更多蚀变细节信息,得出该方法可用于提升矿化蚀变信息精度的结论。本研究将Landsat8数据和WorldView-2数据进行协同,用以浙江漓渚地区铁染蚀变信息的提取。
近年来的蚀变信息提取方法主要有光谱角匹配法(SAM)、主成分分析法(PCA)、波段比值法(BRM)。光谱角匹配法是通过参考波谱的相对重采样与所要处理的遥感影像波谱信息进行匹配,根据两者的相似程度进行地物的识别,从而达到识别地物的目的,在高光谱遥感地物识别中,光谱角匹配法又被称为光谱角填图法,当目标光谱与参考光谱矢量之间的广义夹角越小,两者相似性越大,当该夹角小于某个阈值时,就认为目标光谱与参考光谱匹配成功。也即目标地物被识别为已知参考地物[8]。贺金鑫等人[9]采用光谱角填图法,利用Hyperion 卫星高光谱遥感数据对内蒙古大营铀矿区进行矿化蚀变信息的提取,验证方法的有效性,并且确定含铁离子和粘土矿物类蚀变信息可以作为当地重要的找矿标志。主成分分析法是一种降维思想的数学统计方法,是大部分蚀变信息提取中使用较多的一种较为精确和实用的方法。该方法利用 K-L 变换,将给定的多维相关变量转换成新的彼此相互独立且不相关的变量,由于各变量之间互不相关,新生成的各变量之间信息没有重复与冗余,从而实现将同1幅影像内的各个波段中较为明显的蚀变信息集中到尽可能少的重新组分图像中[10]。刘桂萍等人[11]采用主成分分析法对ETM+遥感影像进行蚀变异常信息提取并进行强弱分级,通过野外查证表明遥感异常区域有较好的矿化蚀变显示,最后进行成矿预测,圈定成矿远景区。
本研究将光谱波段范围较宽的Landsat8影像与空间数据分辨率较高的WorldView-2影像数据进行协同处理,生成具有两种影像数据特征的协同影像,在研究区域使用光谱角匹配法和主成分分析法分别对单个数据及协同后的影像进行处理,得到研究区的铁染矿物蚀变信息的分布,结合当地铁矿分布和实地考察来验证协同影像提取的蚀变信息相对于单一影像更为准确的结论,为当地新矿源的寻找提供一种新的参考方法。
如图1所示,漓渚位于鉴湖水系源头,为绍兴县西南部半山区镇,东邻福全镇,南连兰亭,西接诸暨市店口镇,北与湖塘街道、福全镇仅一山之隔。全区以山地为主,其次是丘陵。研究区的地理坐标:120°26'6''~120°29'12''E,29°54'34''~29°57'16''N。浙江漓渚地区处于钦杭成矿带北东段[12-13],区内有着较为连续的出露地层,发育矽卡岩型铁或铁多金属矿床或矿(化)点,成矿作用与沿北东走向褶皱-断裂带发育的晚侏罗世末栅溪—广山侵入岩密切相关。漓渚地区主要以栅溪岩体和广山岩体为主,成岩构造环境具有从俯冲碰撞作用向后碰撞或后造山作用转换的特征,转换时间与华南地区构造环境转换时间吻合,为华南构造-岩浆演化提供新的约束,且这两类侵入岩在浙江分布广,均具有较好的成矿潜力及找矿前景[14]。
OLI传感器搭载的Landsat8卫星发射于2013年2月,在基本继承TM/ETM+数据特点的基础上,增加了光谱范围和信噪比,进一步对波段设置进行了优化调整。尤其是对NIR波段(0.845~0.885 μm)的调整,避开了大气吸收特征,较为显著地增强了OLI对地物的识别能力[15]。该数据具有9个波段,分辨率为30 m,其中,包括1个15 m的全色波段(0.500~0.680 μm),成像宽幅为185×185 km。Landsat8 OLI数据NIR波段(0.845~0.885 μm)和SWIR2波段(2.100~2.300 μm)波谱分辨率的提高,已分别成为含Fe3+类、Al-OH类和Mg-OH类蚀变矿物的诊断性谱段。本研究选用的Landsat8遥感影像数据获取时间为2018年4月,下载的影像只包含有效蚀变信息的1~7波段。原始影像云量较少,地物层次较为清晰,色调对比度好,能够较好地满足本次遥感矿化信息提取的工作需求。该数据可以从美国地质勘探局(https://earthexplorer.usgs.gov/)下载。影像波段参数如表1所示。
图1 研究区地理位置
WorldView-2卫星于2009年10月发射升空,运行于770 km高的太阳同步轨道,能够提供0.46 m的全色图像和1.85 m分辨率的多光谱图像[16]。其星载多光谱遥感器不仅具有蓝色波段(450~510 nm)、绿色波段(510~580 nm)、红色波段(630~690 nm)、近红外1波段(770~895 nm)4个标准谱段,还包括海岸波段(400~450 nm)、黄色波段(585~625 nm)、红色边缘波段(705~745 nm)和近红外2波段(860~1 040 nm)4个额外波段。影像波段参数如表2所示。
WorldView-2影像数据拥有较高的分辨率和丰富的空间信息,有利于信息的提取和遥感制图能力。本次购买使用的WorldView-2数据获取时间为2018年4月,多光谱图像空间分辨率为2 m。原始影像质量良好,无云及阴影遮盖,能够满足本次遥感矿化信息提取工作的要求。
表1 Landsat8影像波段参数(部分)
表2 WorldView-2影像波段参数
在对两种数据进行预处理后,进行影像协同前需要统一两种数据的分辨率,所以要对影像数据进行尺度转换。两者影像数据的分辨率不同,Landsat8遥感影像的分辨率为30 m,WorldView-2遥感影像的分辨率为2 m。不同分辨率的栅格数据进行空间尺度的统一,转换的方法主要有两种:一是以其中一类数据的分辨率为基准,将另一类数据转换为相同的分辨率;二是以两类数据的分辨率为上下限,选择其中某个适宜的尺度作为基准,将两类数据进行转换[17-21]。本研究采用第一种方法进行实验。以WorldView-2的2 m分辨率为基准,将30 m分辨率的Landsat8遥感影像重采样提高至2 m,再将Landsat8影像数据的7个波段与WorldView-2影像数据的8个波段集成,生成具有15个波段且分辨率为2 m的协同影像,该影像包含比Landsat8 OLI和WorldView-2都丰富的波段信息。波段参数如表3所示。
表3 协同影像波段参数
在生成协同影像后,采用光谱角匹配法和主成分分析法,对3种数据进行蚀变信息的提取,并将提取结果进行对比验证。蚀变信息提取过程的流程如图2所示。
铁染类是对以含铁离子矿物为特征矿物的蚀变类型的统称。铁离子的价态及矿物质的透明度和含水性等很大程度上决定其波谱特征,含铁矿物主要有褐铁矿、磁铁矿、赤铁矿、黄钾铁矾等。它们的光谱曲线均具有如下特点[22]:在0.4~0.5 μm波段处出现吸收的特征;在附近有一个强的吸收特征,意味着特别典型的矿物;从强吸收开始,随着波长增加反射率急速上升至附近达到最高值,随后呈缓慢下降趋势。
图3 采样后的蚀变光谱曲线
根据研究区的成矿特征,选取褐铁矿作为铁染蚀变的特征矿物的代表。褐铁矿(Limonite)是在温湿条件下黄铁矿的一种蚀变矿物。由数量不一的含水氧化硅、枯土、纤铁矿等混合而成。USGS光谱库和分别用3种数据重采样后的褐铁矿的蚀变光谱曲线如图3所示。
由图3 (a)可知,褐铁矿在0.5 μm和0.9 μm附近具有较为明显的吸收特征,是Fe3+类矿物的诊断性谱段。虽然从图3 (b)中的OLI波谱也能观察到这一特征,但同USGS光谱库的褐铁矿蚀变光谱中的吸收特征相比不太显著,而通过集成Landast8和WorldView-2数据生成图3 (d)的协同影像数据重采样波谱,在0.5 μm和0.9 μm附近相比OLI数据有着更为显著的光谱吸收特征。协同影像波谱具备这种更显著的光谱吸收特征是因为在0.5 μm和0.9 μm的波段范围内,图3(c)的WorldView-2数据对Landsat8数据缺失的部分信息进行了补充,这正是影像协同的关键和原理所在。
光谱角匹配法(SAM)通过参考波谱的相对重采样与所要处理的遥感影像波谱信息进行匹配,把光谱在N个波段的响应作为N维空间的矢量,计算两个光谱之间的广义夹角,来表示地物之间相似程度,进行地物提取。其表达式为
两者之间的广义夹角越小,则光谱相似程度越高,表明识别地物是参考地物的可能性也越高[23]。本研究的参考波谱为USGS岩石波谱,对影像进行光谱适配,选择合理的阈值进行密度分割,可以获得较为准确的蚀变信息二值图像,如图4(c)、图4(d)、图4(e)所示。
主成分分析法(PCA)可以尽量减少新生成变量之间信息的重复与冗余,去除波段间的多余信息,一般来说,前3个主成分包含所有波段95%以上的信息量。主成分分析法是目前用于蚀变信息提取使用次数较多、也较为精确一种方法。借鉴 Crosta 主成分分析法[24]的思想,对Landsat8数据的铁染蚀变信息提取,采用 OLI2、OLI4、OLI5、OLI6等 4 个波段组合进行主成分分析。由于WorldView-2数据最大的波长范围在0.860~1.040 μm,没有办法找到对应OLI6(1.560~1.660 μm)的波段,所以无法对WorldView-2数据进行铁染异常的主成分分析。选择协同影像数据的Band3、Band9、Band12、Band14进行蚀变信息提取。然后观察特征向量矩阵(见表4、表5),在其中选择出能代表铁染矿物的主成分向量。
选择OLI2与 OLI5系数符号相同,OLI5系数与OLI4 及OLI6的系数符号相反的本征向量PC3作为代表铁染矿物的主分量[25],利用门限化技术,根据PC3的标准差σ及其倍数进行等级划分、密度分割渲染,用不同层次的颜色代表不同的铁染蚀变强度,其等级划分标准见表6。以表7中的铁染异常等级进行划分,得到铁染类蚀变信息分级结果见图4(a) 。
表4 PCA[2 4 5 6]主成分特征向量矩阵
表5 PCA[3 9 12 14]主成分特征向量矩阵
表6 铁染强度等级划分标准
图4 铁染类蚀变信息分级结果
表7 铁染蚀变遥感异常分级
根据Crosta法异常所在组分判断原则,选定协同数据的PC4为主成分向量后,根据标准差及其倍数进行铁染蚀变强度等级划分。分级结果如图4 (b)所示。
从图5可以看出,铁染异常主要是沿着研究区的东北部分延展,以中部位置蚀变信息最为明显,在西北部和东南部有零散分布。WorldView-2数据因为缺少1.04 μm以上的波段,所以无法进行铁染蚀变信息的主成分分析,即在主成分分析的步骤中,在高于1.04 μm的部分,协同影像提取蚀变信息的精度没有提升,但对比图4(a)和图4(b)发现,协同影像和Landsat影像提取铁染蚀变信息的范围和位置较为接近,协同影像在细节提取方面比Landsat影像要更好,说明WorldView-2数据中低于1.04 μm的部分在协同后不仅没有给数据带来相对较大的误差影响,还提供了丰富的空间信息。通过图4也可以看出,协同影像较单一影像的SAM法提取结果要更为接近主成分分析法的提取结果。
对比单一影像和协同影像提取结果图4(c)、图4(d)、图4(e)存在差异的地方,随机挑选区域1和区域2进行研究。在区域1中,WorldView-2数据中没有表现出来明显铁染异常,而Landsat8数据和协同数据均表现出铁染异常;在区域2中,Landsat8数据和WorldView-2数据中均没有表现出明显铁染异常,而协同影像结果表现出铁染异常。进行野外考察,采集区域1、2内有铁染异常区域的岩石样本,发现铁染蚀变信息的存在,证明协同影像提取的蚀变信息相对于单一影像更为准确。
图5 野外实地考察照片
本研究利用协同影像集成Landsat8数据和WorldView-2数据优势,采用主成分分析法和光谱角匹配法对漓渚地区的铁染类矿物的蚀变信息进行提取,通过两种方法的提取结果对比可以发现,协同影像在集成了WorldView-2数据优势后,在细节提取方面比Landsat影像要更好。通过对比SAM法蚀变信息提取结果不同的分布区域,进行野外实地考察,依据铁矿矿点的分布和采集样本的验证,得到结论:协同影像提取的蚀变信息相对于单一影像更为准确。