普朗斑岩型铜矿区矿化蚀变特征矿物填图及找矿潜力分析：来自资源一号02D 遥感卫星的证据

2022-04-01赵志芳张新乐张瑞丝李文昌曹晓民

沉积与特提斯地质 2022年1期

赵志芳，张新乐，陈琪，张瑞丝，李文昌，曹晓民

（1.云南大学地球科学学院，云南昆明 650050；2. 自然资源部三江成矿作用及资源勘查利用重点实验室，云南昆明 650050；3.云南省国产高分卫星遥感地质工程研究中心，云南昆明 650050；4.中国地质调查局成都地质调查中心，四川成都 610081；5.云南省地质调查院，云南昆明 650051）

0 前言

普朗铜矿是近年来西南三江地区取得的重要找矿勘查突破，目前已规模化开发利用。前人在普朗铜矿矿床地质特征、成岩成矿时代、含矿斑岩及构造背景、成矿物质来源、矿床成因模式、矿化蚀变特征、地球物理地球化学遥感等找矿标志等方面均取得了一系列研究认识。已有研究表明：普朗斑岩型铜矿矿化与花岗闪长斑岩、石英二长玢岩关系密切，上述两类岩体为普朗铜矿主要含矿岩体（谭康华等，2005；曾普胜等，2003，2004a，2004b，2006；范玉华等，2006；李文昌等，2011；Xia et al.，2021）。普朗斑岩型铜矿具有显著的蚀变分带特征，即从中心向外围主要发育钾化硅化带—绢英岩化带—青磐岩化带，且高品位铜矿矿化富集主要赋存于钾化硅化带（钾长石、石英）、绢英岩化带（绢云母、石英）中。普朗铜矿于2019 年投产使用，开采设计规模达1250 万吨。目前勘查程度较高的主矿体已作为首采区进行规模开采，但外围开发潜力仍不明确，开展找矿勘查研究紧迫性强。由于首采区外围海拔多处于4000 米以上，地形陡峻，气候严寒，高山草甸、针叶林等植被覆盖较厚，地质露头较少，传统地质调查手段探测存在困难。加之岩矿测试分析有限，使得首采区外围储量探明及深部找矿勘查难以深化推进。

高光谱遥感技术尤其是国产高光谱遥感技术的飞速发展，无疑为条件艰苦地区全覆盖地质调查提供了有力手段。以往研究表明，高光谱遥感技术在矿物识别方面具有显著优势，已被成功应用到找矿勘查领域（Bierwirth et al.，2002；甘甫平，2009；王润生等，2010；Bedini et al.，2011；董新丰等，2018；连琛芹等，2020）。但以往受昂贵的国外Hyperion 遥感数据及小区域范围的机载成像CASI、SASI、TASI高光谱数据源制约，高光谱矿物填图多停留在科学研究、应用示范层面（Cudahy et al.，2002；甘甫平等，2002；Kruse et al.，2003；Pignatti et al.，2009；侯珂，2012；宋晚郊，2013；Amin et al.，2014；Kumari et al.，2014；张川等，2015；Oskouei and Babakan，2016；叶发旺等，2018；Jain et al.，2019）。国内外学者主要探索建立的最小噪声分离变换（Minimum Noise Fraction，MNF）高光谱噪声去除、纯净像元法（Pixel Purity Index，PPI）、矿物分层识别、基于光谱吸收的矿物端元探测及光谱角（Spectral Angle Mapper，SAM）等方法，具有较好地分离高光谱冗余噪声、识别绢云母、绿泥石等蚀变矿物的能力；构建的独立成分分析（Independent Component Correlation Algorithm，ICA）及基于端元矿物标准波谱的匹配噪声分离等方法则在一定程度上能较好消除植被影响（李亮，2014；Tayebi et al.，2015；董新丰等，2020；连琛芹等，2021）。随着资源一号02D高光谱卫星于2019 年9月12 日成功发射，因其高光谱载荷具有幅宽大，短波红外光谱分辨率达20nm、波谱数量达166 的技术优势及高性价比等特点，在我国大范围精细化矿物填图及其找矿勘查等业务化应用中具有良好应用前景。目前资源一号02D 星仅在哈密遥感实验场地层、岩性、蚀变矿物识别方法研究（李娜等，2020），尚未见相关植被覆盖区高光谱矿物识别应用研究报道。

基于此，本文以普朗矿山为研究区，基于资源一号02D高光谱遥感技术，采用比值植被指数阈值大于2 区分植被覆盖区及非植被覆盖区，基于实测波谱曲线，建立符合普朗铜矿区的分层次蚀变矿物波谱库，采用匹配滤波（Matched Filtering，MF）的方法，优选绢云母、绿泥石及绿泥石蚀变信息提取方法并予以示范应用，野外验证。本文探索建立了植被覆盖背景下的混合像元蚀变矿物高效精细化提取方法，并初步摸清普朗矿山首采区外围蚀变矿物分布状况，为矿山后续找矿勘查提供技术支撑。

1 研究区概况

普朗斑岩型铜矿区位于云南省西北部，隶属于云南省迪庆藏族自治州香格里拉市，位于县城北东方向36km 处（图1）。普朗斑岩铜矿属印支期成矿，是我国近年来新发现的超大型矿床（吴练荣等，2020）。

研究区位于云南省普朗铜多金属矿集区，义敦岛弧南端，西界为格咱断裂带，东部和南部为甘孜-理塘结合带，是西南“三江”古特提斯洋演化阶段形成的重要地质构造单元之一，即印支期甘孜-理塘洋向西俯冲于格咱-中甸微陆块所诱发的大规模构造-岩浆作用的产物，也是三江特提斯成矿域内重要的斑岩-矽卡岩型铜钼多金属矿集区（图2）。

图2 普朗铜矿地质简图（据吴练荣，2021修编）Fig.2 Geological sketch map of Pulang copper mining area（modified from Wu et al.，2021）

区内地层主要由上三叠统碎屑岩、碳酸盐岩和火山岩及次火山岩组成。区内晚三叠世普朗铜矿成矿作用发生于普朗复式斑岩体内，在岩体中心形成了以细脉浸染状矿石为主的脉状矿体。该斑岩体由3 个期次岩浆侵入形成，其中第1 阶段（（221.0 ±1.0）Ma）主要形成石英闪长玢岩，第2 阶段（（211.8 ±0.5）Ma）主要形成石英二长斑岩，第3阶段（（206.3 ±0.7）Ma）形成花岗闪长斑岩，均属于印支期构造-岩浆活动的产物（刘学龙等，2018）。

普朗矿区斑岩体蚀变强烈，具典型的斑岩型蚀变分带，从内到外可划分为钾化硅化带、绢英岩化带、青磐岩化带（张少颖等，2020；陈琪等，2021；吴练荣，2021）。这些蚀变特征矿物也是普朗铜矿的重要示矿标志。

钾化硅化带（KSi）：钾化硅化带大体上与工业矿体相对应，位于岩体中心部位，具有强钾交代和硅化的特征。蚀变岩中发育有大量的钾长石和石英，并伴随有少量的黑云母脉，蚀变矿物总量约占20%～25%。

绢英岩化带（SiSe）：绢英岩化带通常分布在钾化-硅化带的边缘地带。两个蚀变分带之间无明显的界线，但在其交界处绢英岩化带呈穿插、叠加关系。绢英岩化带与工业矿体基本上对应，紧靠岩体中心部位。绢英岩化带内具有强烈的绢云母化、硅化，钾交代作用少量发育，多被绢云母化、硅化逐渐取代。绢英岩化带亦是矿区的重要含矿部位之一，发育有大量黄铁矿化和黄铜矿化。黄铜矿化品位比较低，品位约为0.2% ～0.5%，在绢英岩化带与钾化硅化带交界处常形成富集的工业矿体。

青磐岩化带（ChEp）：青磐岩化带通常远离岩体中心部位，主要发育在石英闪长玢岩中。蚀变范围发育广，但强度低。青磐岩化也会叠加发育在钾化硅化带、绢英岩化带上，此时矿化亦随之增强。

2 数据源及预处理

2.1 数据源

资源一号02D卫星于2019 年9 月12 日成功发射，卫星搭载的两台相机，可有效获取115km 幅宽的9 谱段多光谱数据以及60km 幅宽的166 谱段高光谱数据，其中全色谱段分辨率可达2.5m、多光谱为10m、高光谱优于30米，高光谱载荷可见近红外和短波红外光谱分辨率分别达到10nm和20nm（表1）。

表1 资源一号02D卫星参数Table 1 Parameters of Ziyuan-1 02D satellite

2.2 数据预处理

原始数据影像信息储存的格式为传感器记录的数字量化值，需要进行辐射定标将其转化为具有符合实际物理意义的辐射亮度值：

式（1）中R 表示辐射亮度值，Gain 和Offset 为传感器自身的增量和偏移，DN 表示传感器记录的灰度值。

对辐射定标后的影像进行大气校正，将辐射亮度值恢复为地表实际反射率。利用RPC 模型进行正射校正，纠正各种因素引起的几何变形，对影像进行地理坐标定位，使影像恢复真实大地坐标信息（图3、图4）。

图4 大气校正前（a）和后（b）波谱曲线Fig.4 Spectral curve before（a）and after（b）atmospheric correction

2.3 植被覆盖程度划分

研究区有大面积覆盖植被，需要去除或抑制植被信息才能较好获得地质异常信息。传统处理中，往往将植被覆盖区进行掩膜去除，但该方法会将隐含于植被中的矿化弱蚀变矿物等信息也一并去除，可能造成矿化蚀变矿物信息的缺失。为较好弥补以往工作不足，本文构建了分层次蚀变信息提取方法。即首先区分植被覆盖区及非植被覆盖区，并实际测定不同区域蚀变矿物波谱，构建植被覆盖区及非植被覆盖区分层次蚀变矿物波谱特征，进而采用匹配滤波等提取蚀变矿物信息。

本文采用植被指数划分植被覆盖程度。常用的植被指数法包括：归一化植被指数（NDVI）、比值植被指数（RVI）、土壤调整植被指数（SAVI）和增强型植被指数（EVI）（王泽等，2021）。其中比值植被指数在中高植被覆盖区域划分中效果显著。结合研究区中高植被覆盖程度实际，本文采用比值植被指数进行植被覆盖程度划分，具体公式为：

式（2）中，R 为红光波段，NIR 为近红外波段。经反复试验，研究区RVI ＝2 ～8 为植被覆盖区，RVI值接近1 则为建筑、裸地、水体等。本次依据实验结果，将＜2 的区域划为非植被覆盖区，≥2 的区域划为植被覆盖区（图5）。

图5 非植被覆盖区（a）与植被覆盖区（b）划分Fig.5 Non-vegetation-covered area（a）and vegetation-covered area（b）

3 蚀变矿物信息提取

3.1 实测波谱分析

针对普朗斑岩型铜矿矿化蚀变矿物为绢云母、绿泥石、绿帘石，本次实验采用中国地质调查局南京中心研制的地物波谱仪CSD350A（扫描范围为350 ～2500nm，光谱分辨率高达2nm）对普朗矿区86件样品进行绢云母、绿泥石、绿帘石等波谱实测（图6）。扫描时将每件样品分为3 ～4 面，每一面扫描3～4 条谱线。去除实际误差较大的谱线后共1039条谱线。

图6 普朗地区蚀变矿物实测波谱采样位置Fig. 6 Sampling locations of ground measured alteration minerals in Pulang area

根据CSD350A 携带的矿物光谱分析专家系统MSA测定数据运算，分植被覆盖区及非植被覆盖区，分别筛选出样品中含量大于60%的绢云母、绿泥石和绿帘石矿物波谱曲线，计算绢云母、绿泥石、绿帘石等矿物的平均波谱数值，构建植被覆盖区及非植被覆盖区绢云母、绿泥石、绿帘石波谱曲线（图7、图8、图9）。

图7 普朗矿区非植被覆盖区（a）及植被覆盖区（b）绢云母实测波谱特征Fig.7 Ground measured spectral curve of sericite at non-vegetation-covered area（a）and at vegetation-covered area（b）in Pulang area

图8 普朗矿区非植被覆盖区（a）及植被覆盖区（b）绿泥石实测波谱特征Fig.8 Ground measured spectral curve of chlorite at non-vegetation-covered area（a）and at vegetation-covered area（b）in Pulang area

图9 普朗地区非植被覆盖区（a）及植被覆盖区（b）绿帘石实测波谱特征Fig.9 Ground measured spectral curve of epidote at non-vegetation-covered area（a）and at vegetation-covered area（b）in Pulang area

对比植被覆盖区实测波谱曲线、非植被覆盖区实测矿物波谱曲线、美国地质调查局USGS 标准矿物光谱库波谱曲线，蚀变矿物反射率、吸收深度与吸收宽度等存在一定差异，但特征蚀变矿物在某些吸收位置保持了一致性（图10）。例如，绢云母在2205nm、2355nm左右吸收峰（张川等，2017），绿泥石在2252nm、2340nm左右吸收峰（梁树能等，2014）以及绿帘石在2252nm、2334nm 左右吸收峰（熊燕云等，2019）。基于此，结合羟基矿物的吸收位置具有特殊性，本实验初步优选2000 ～2500nm波段，并将其重采样到资源一号02D高光谱数据上进行后续运算。

图10 蚀变特征矿物非植被覆盖区实测波谱曲线（a）、植被覆盖区实测波谱曲线（b）与USGS波谱曲线（c）Fig.10 Ground measured spectrum curve of non-vegetation-covered area（a），vegetation-covered area（b）and USGS spectrum library（c）for alteration minerals

3.2 最小噪声分离变换

高光谱遥感数据往往包含有噪声信息，这对蚀变信息的提取会产生较大干扰，因此采用最小噪声分离变换（MNF）对植被覆盖区及非植被覆盖区分别进行降噪处理，为下一步的混合像元分解做好准备。

MNF 本质上是两次主成分变换，第一次变换（基于估计的噪声协方差矩阵）用于分离和重新调节数据中的噪声，该运算使变换后的影像噪声数据之间具有最小的方差且波段性弱；第二步是对噪声白化数据（Noise-whitened）的标准主成分变换（李海涛等，2007）。

首先，利用高通滤波器，对高光谱影像数据进行滤波处理，得到噪声协方差矩阵CN，将其对角化为矩阵DN：

式（3）中，DN为CN的特征值按照降序排列的对角矩阵；U 为由特征向量组成的正交矩阵，I 为单位矩阵，P为变换矩阵。

其次，对噪声数据进行标准主成分变换：

式（4）中，CD为影像X的协方差矩阵，CD-adj为经过P变换后的矩阵，DD-adj为CD-adj的特征值按照降序排列的对角矩阵，V为由特征向量组成的正交矩阵。

变换后结果如图11、图12 所示。

图11 MNF处理结果（a. 非植被覆盖区；b. 植被覆盖区）Fig.11 MNF transformation results（a. non-vegetation-covered area；b. vegetation-covered area）

图12 MNF处理波谱曲线（a. 非植被覆盖区；b. 植被覆盖区）Fig.12 MNF transformation result curves（a. non-vegetation-covered area；b. vegetation-covered area）

3.3 混合像元分解

受传感器空间分辨率的限制，在同一像元内包含可各种各样的地物类型，形成了混合像元。为便于计算，往往在数据处理中将混合像元视为线性混合，即确定某一空间内具有各自丰度的端元混合到一起，因此也往往按线性方式进行混合像元分解。

（1）线性混合分解模型（Linear Mixture Model，简称LMM）

线性混合模型（LMM）忽略地物间的多次散射（Adams et al.，1986），认为混合光谱是端元光谱及其丰度的线性组合：

式（5）中，n 为混合像元内端元的总数，fk为在k 点位的端元光谱，Rk为在k 点位的丰度，ε 为残差项。该模型通常用最小二乘法进行求解得到端元丰度。为保证求解结果不失去物理意义，常附加约束条件：端元丰度总和为1 约束及非负性约束。

（2）匹配滤波器模型

匹配滤波（MF）最初的设计是在信号处理领域，被认为是最优的线性检测放大。由于它不需要已知所有的端元波谱，在高光谱图像处理领域得到了广泛的应用（林娜等，2011）。该方法的优点在于，它将已选定感兴趣区作为端元波谱，未知波谱作为背景值，最大化地突出已知端元波谱的信息，同时抑制背景的干扰。这种方法提供了一种快速探测需求波谱，而并不需要了解影像内所有波谱的信息的技术（陶秋香，2007）。该方法的取值结果-1 ～1，＜0 的部分为背景值，＞0 且趋近于1 则接近于标准测定波谱。

为较好探索植被覆盖的影响，本次实验将区分植被覆盖区和非植被覆盖区、不区分植被覆盖区和非植被覆盖区两种情况进行匹配滤波信息提取处理。

本次实验以“平均值＋标准差”作为阈值进行蚀变矿物异常信息划分，实现蚀变矿物填图。

蚀变矿物异常信息提取结果详见表2、表3、表4 和图13、图14。

表2 未划分植被与非植被覆盖区MF处理结果Table 2 MF results of undivided vegetation

表3 非植被覆盖区MF处理结果Table 3 Result of MF in non vegetation areas

表4 植被覆盖区MF结果Table 4 Result of MF in vegetation areas

图13 未划分植被与非植被覆盖区MF蚀变信息提取结果（a. 绢云母；b. 绿泥石；c. 绿帘石；d. 蚀变矿物异常综合图）Fig.13 MF extraction results of alteration minerals in undivided vegetation-covered and non-vegetation-covered area（a. sericite；b.chlorite；c. epidote；d. comprehensive alteration minerals）

图14 划分植被与非植被覆盖区蚀变矿物异常信息提取结果（a. 绢云母；b. 绿泥石；c. 绿帘石；d. 综合蚀变矿物）Fig.14 MF extraction results of alteration minerals of divided vegetation-covered and non-vegetation-covered area（a. sericite；b.chlorite；c. epidote；d. comprehensive alteration minerals）

4 结果与讨论

4.1 实验结果分析

分析上述试验结果表明：划分了植被覆盖区和非植被覆盖区，进行匹配滤波数据处理可较好揭示植被覆盖区蚀变特征矿物异常信息（图15）。

图15 研究区蚀变矿物异常信息分布及成矿预测靶区分布图（a.未区分植被与非植被覆盖区；b.区分植被与非植被覆盖区）Fig.15 Distribution of alteration minerals anomaly information and metallogenic prediction targets in Pulang mining area（a.undivided vegetation and non-vegetation cover；b. divided vegetation and non-vegetation cover）

进一步剖析典型遥感异常信息表明，无论是否划分植被覆盖程度，在6、23 号野外观测点分别有绢云母、绿泥石蚀变特征矿物异常显示；但对于5 号、24 号和25 号野外观测点来说，当未划分植被覆盖区和非植被覆盖区进行匹配滤波数据处理时，绢云母蚀变矿物信息不能得到有效识别；而在划分植被覆盖区和非植被覆盖区后，开展匹配滤波数据处理时，绢云母蚀变矿物异常信息均得到较好反映。可见，区分植被覆盖程度的蚀变矿物填图具有显著优势（表5）。

表5 分层次蚀变矿物信息提取结果野外验证情况对比Table 5 Contrast of verification for hierarchical extraction results of alteration mineral information

收集前人研究成果及勘探资料进一步验证表明，采用资源一号02D 岩矿敏感的2000 ～2500nm高光谱数据、划分植被覆盖区及非植被覆盖区的分层次蚀变矿物提取方法，可较好地探测研究区尤其是植被覆盖区的蚀变矿物分布状况。即，本文采用分层次蚀变矿物提取方法具有较强可靠性。特别地，在普朗首采区东矿段植被覆盖区，提取的绢云母遥感异常面状浓集发育，这与野外地质勘查认识一致（图16）。普朗首采区东矿段成、控矿主要受岩体边部的节理、裂隙构造控制，岩浆上涌沿主矿体边部的次级构造充填富集，形成绢云母及铜、铅锌多金属矿化富集，其中矿化以脉状矿产出为主（吴练荣等，2020；王国强等，2018）。

图16 研究区首采区东段成矿预测靶区岩心样镜下鉴定照片Fig.16 Microscopic identification photos of rocks and minerals in the east section of the first mining area

结果分析亦表明：普朗矿区首采区外围，仍具有与首采区类似、绢云母和绿泥石、绿帘石蚀变遥感异常信息叠加发育组合的特点。据此，在首采区外围具有上述蚀变矿物异常浓集的区域，圈定找矿靶区3 处，即Ba1、Ba2、Ba3。进一步找矿勘查亦证实，首采区外围尤其是东部地区Ba1，植被覆盖区内绢云母浓集，部署的1 个钻孔见矿其岩芯样中铜品位达0.1% ～0.4%，预示了区内资源潜力巨大（图17）。

图17 钻孔D1岩心Fig.17 Photograph of No.1 drilling core

4.2 讨论

（1）区内绢云母蚀变特征矿物可能是普朗矿区斑岩型铜矿化重要稳定的示矿近矿标志。

高光谱遥感提取的绢云母富集特征，在野外尤其是东矿段均得到了较好验证，高光谱遥感信息提取能更细化的揭示首采区外围地带蚀变分带（绢英岩化、绿泥石化、绿帘石化）特征。结果同时显示：无论在普朗首采区还是外围尤其是东矿段，区内绢云母遥感异常均具有面状浓集发育特点。结合钻孔及镜下鉴定等勘探成果分析，普朗矿区斑岩型铜矿化过程中，绢云母化蚀变特征矿物均较发育，推测其参与了斑岩型铜矿化的重要过程，绢云母化蚀变特征矿物可能是区内斑岩型铜矿化最重要稳定的示矿近矿标志。

（2）无论是植被覆盖区还是非植被覆盖区，蚀变特征矿物提取结果可能指示了普朗矿区斑岩型铜矿具有多期次成矿作用叠加复合成矿的特点。

本次高光谱遥感蚀变矿物异常信息提取结果表明，在首采区主矿体部位或外围植被覆盖区或非植被覆盖区，均具有面状绢云母和点状绿泥石、绿帘石遥感异常信息叠加发育特点，可能反映的是多期次成矿作用进而矿化富集的指示。因此，圈定的具有面状绢云母和点状绿泥石、绿帘石遥感异常富集特征的成矿预测靶区Ba1、Ba2、Ba3，应具有较好的找矿条件。其中，矿区东矿段绢云母化和绿泥石绿帘石叠加发育的地段，应指示了斑岩型铜矿矿化富集中心。