砂岩铀矿光谱特征匹配技术研究

1砂岩铀矿光谱特征匹配技术研究1张杰林1,2岑长华1陆书宁3韩晓青1（1．核工业北京地质研究院，北京100029；２．中国矿业大学，北京100083；3．国家卫星应用中心，北京100055）[摘要]本文根据成像光谱数据的特点，采用高精度光谱数据预处理技术和基于相关系数测度的光谱匹配技术，实现了岩石矿物光谱特征波形对比分析及诊断光谱信息提取。在此基础上，以云南腾冲铀矿区和内蒙古海拉尔盆地为例，进行试验研究，取得实效。[关键词]成像光谱；光谱匹配；砂岩铀矿中图分类号：TP790引言高分辨率定量化成像光谱技术是当代遥感技术的前沿之一。它将遥感技术由定性分析阶段推向定量化分析识别这一新的高度[1，2，3]。成像光谱技术具有图谱合一的特征，通过成像光谱测量所获得的地球表面图像包含了大量的空间、辐射及光谱三重信息[2]，形成了高光谱图像数据立方体。因此，如何充分地有效利用这一海量数据进行地面目标的定量化、自动化识别，从而为铀矿勘查新技术、新方法的研究服务，便成为当前铀资源勘查领域必须解决的基础关键技术问题之一。本文就是在这种背景下，重点开展砂岩铀矿光谱特征匹配技术的研究。1光谱数据预处理技术针对地面光谱测量提供的大量的不同时间、不同区域的光谱数据，采用有效的预处理技术，可以实现高质量的数据管理与维护，为进一步的数据处理、分析与应用奠定基础。该项技术主要包括野外实测光谱曲线的数据库管理，包络线消除技术及光谱微分技术等。1.1实测光谱曲线的数据库管理本次研究工作中，基于ENVI软件系统，将云南腾冲、内蒙海拉尔等地区的有关铀成矿岩石及背景地物的实测光谱数据，实现数据库管理，建立了砂岩铀矿及背景光谱数据库（表1、图1）。通过数据库管理，可以方便查询及实时显示任何一条光谱曲线。入库的光谱数据可以是野外实测数据，也可以从成像光谱图像上采集。光谱数据的输出格式有两种：文本格式和图形格式。另外，系统的主要功能还包括光谱数据重采样、特征数据提取与波形分析、混合光谱分解等。1核工业遥感信息与图像分析技术重点实验室基金项目编号：遥U66.６2表1研究区测量目标一览表Table1Thelistofmeasuringobjectinthestudyareas测量地区典型目标背景地物云南腾冲盆地含砾砂岩（矿化）、细粒砂岩（矿化）、中粗粒片麻质砂砾岩（矿化）、中粗粒砂砾岩（矿化）黑云母斑状花岗岩、白马草（50号矿区）、白马草（团田）、粗粒砂岩、含砾砂岩、泥岩、泥质砂岩、细砂岩、中粒砂岩内蒙古海拉尔盆地含碳泥岩(矿化)、细砂岩（矿化）、砂质泥岩（矿化）、细砂岩（矿化）中粒砂岩、含碳泥岩、砂砾岩、中细粒砂岩、尖叶草、阔叶草、短叶草.5001000150020002500020406080100123456图1云南腾冲地区铀成矿岩石及背景地物光谱数据库Fig.1Spectrumdatabaseofsandstone-typeuraniumdepositsandbackgroundvegetationandgroundobjectsinTengchongarea,Yunnanprovince1——含砾砂岩（矿化）；2——细粒砂岩（矿化）；3——片麻质砂砾岩（矿化）；4——中粗粒砂砾岩（矿化）；5——白马草（50号矿区）；6——水汽吸收带。1.2包络线消除技术为了客观分析波谱曲线的吸收特征，利用包络线消除技术，将波谱曲线中与目标物质成分密切相关的典型的吸收峰提取出来，用统一的基线来对比每一个吸收峰，进行光谱波形的分析研究。波长(nm)反射率(%)3包络线定义为逐点直线连接那些突出的峰值点，并使折线在峰值点上的外角180°[1]。包络线消除归一化后(图2)，那些峰值点上的相对值为1，非峰值点则1，这样可以很容易地测定吸收峰的各项参数，显示地物成分的分布，并对光谱特征匹配技术研究提供基础。1.3光谱微分技术采用光谱微分技术，主要是为了有效地消除在不同测量时间获得的光谱数据之间的系统误差、消弱大气散射和吸收对目标光谱的影响，以便提取不同的光谱参数（波长位置、深度和宽度等），进行光谱匹配分析。对反射光谱计算不同阶次的导数，可以迅速确定光谱曲线的拐点及最大最小反射率和波长位置(图3)。其计算公式如下:00.20.40.60.813007001100150019002300波长(nm)包洛线消除前反射率曲线包络线包洛线消除后光谱吸收特征曲线图2包络线消除前后光谱对比图Fig.2ComparasionThespectrumcomparemapofrawdataandcontinuumremovingcurve400600800100012000204060400.00800.001200.00-0.200.000.200.40(a)(b)图3内蒙海拉尔盆地砂岩的反射光谱曲线(a)和微分光谱曲线(b)对比分析图Fig.3Thesandstonespectrumcomparemapofreflectancecurve（A）andderivativecurve（B）A、一阶差商计算公式[1，4，5]：波长(nm)波长(nm)反射率(%)反射率导数值反射率(%)4λλλλ∆×−=−+2][)()()('11iiiRRR（1）其中：)('iRλ——波长λi处的反射率一阶差商值；)(1+iRλ、)(1−iRλ——波长λi+1和λi-1处的反射率值；λ∆——波长间隔。B、二阶差商的计算公式[1，4，5]为：2)()()()('']2[11λλλλλ∆+−=−+iiiiRRRR（2）其中：)(''iRλ——波长λi处的反射率二阶差商值；)(1+iRλ、)(iRλ、)(1−iRλ——波长λi+1、λi和λi-1处的反射率值；λ∆——波长间隔。2光谱匹配模式分析光谱匹配技术主要是指实测光谱与标准光谱（参考光谱）的对比，以计算二者之间的相似性或差异性，定量化分析地物吸收光谱特征，从而进一步研究岩石矿物的物质组成、结构构造及表面状况等。匹配模式主要包括：二值编码匹配和基于相关系数测度的光谱匹配模式。2.1光谱二值编码匹配光谱二值编码技术是一种较为常用的光谱匹配方法，其基本原理就是对预处理的光谱数据首先指定一个阀值（T），大于此阀值的数据赋1，否则赋0，这样实测光谱可用简单的0-1序列来表述，完成编码后，可利用距离度量的算法来进行匹配识别[2，6-8]。(a)（b）图4云南腾冲381矿岩芯(中粗粒砂砾岩)原始反射率曲线及二值编码曲线图Fig.4The0-1codingcurvemapofmedium-grainedgranulitein381miningareaofYunnanTengchonga——原始反射率曲线b——二值编码曲线波长(nm)波长(nm)反射率%编码值012000210022002300240002040608010020002100220023002400反射率(%)编码值5图4就是利用二值编码技术获得的云南腾冲381矿岩芯(中粗粒砂砾岩)的编码曲线图。实验分析表明：二值编码技术使光谱匹配过程快速、有效，有助于提高光谱数据分析处理的效率，但由于这种技术在处理编码过程中会丢掉许多光谱信息的细节，因此该项技术只适合于较粗略地识别岩石矿物光谱，而不能应用于精细的高光谱数据的匹配识别研究中。正基于此，本次工作采用相关系数测度进行光谱匹配的分析处理研究。2.2基于相关系数测度的光谱匹配技术光谱匹配技术主要是依据距离来度量样本的相似性。对于成像光谱数据匹配而言，相关系数测度同时兼顾平均值和协方差两个统计参量，是计算成像光谱超维空间中两类统计距离的理想测度。为了消除各种干扰因素的影响，保证各变量在分布上的同等地位，光谱数据匹配时，要首先对数据进行标准化变换。2.2.1数据变换及匹配算法(1)利用标准差进行数据标准化变换[4]iiisxxx−=′（3）其中：∑==niixnx11（4）∑=−−=niiixxns12)(11（5）Xi——波长i处的反射率值用标准差处理后的变量其均值为0，标准差（均方差）为1。(2)、利用相关系数进行匹配∑∑∑===−−−−=nknkjkjikinkjkjikiijxxxxxxxxc11221)(.)())((（6）其中：Ｋ——样本数；Xki——第K个样本实测光谱数据在波长i处的反射率值；Xkj——第K个样本标准光谱数据在波长i处的反射率值；6ix、jx的计算同x。2．2．2算法实现本次工作中利用VC设计完成了光谱匹配软件的设计，其主要功能包括数据的显示，匹配区选择，光谱数据的多参数匹配及匹配结果的动态显示与存储，为典型目标波形对比分析提供了强有力的工具。(1)程序的主要类和类结构根据本程序的功能，编程中共设计了8个主要类，分别是：cCMath1Doc文档类：负责数据的存储、管理及各种处理结果计算等；dCMainFrame框架窗口类，负责数据的读入，视图窗口的建立；eCMathView视图类，负责显示、管理反射率曲线图形；fCOneView视图类，负责显示、管理一阶导数曲线图形；gCTwoView视图类，负责显示、管理二阶导数曲线图形；hCDiffView视图类，负责显示、管理标准差曲线图形；iCselectDialog对话框类，负责数值选择匹配区域；jCresultDialog对话框类，负责显示保存的匹配结果。相应的类结构图如图5所示。图5光谱匹配程序的类结构图Fig.5FlowchartshowingthestructureofthespectrummatchingprogramCMathCMath1DocCMathFrameCMathDialogCMathViewCMainFrameCChildFrameCSelectDialogCResultDialogCOneViewCTwoViewCDiffView7Cmath-项目名称；CMathFrame-框架窗口类；CChildFrame-子框架窗口类；CMathDialog-对话框类。(2)程序功能本程序可以计算两条光谱曲线的一阶导数、二阶导数及标准差，并在指定区间或最大范围内（无指定区间时）计算原始光谱反射率值、一阶导数、二阶导数及标准差的匹配系数。在计算过程中，可以分别显示原始光谱反射率曲线图、一阶导数、二阶导数和标准差曲线，并将匹配结果以文本文件的形式存盘。(3)对输入的光谱数据文件的要求光谱数据文件格式为：波长值空格反射率值回车换行符数据文件内可有空行或注释行，程序将自动忽略，继续读下一行数据。另外，波长数据不必按顺序排列，程序将自动排序，波长间隔也可以是任意的，程序将自动用插值法将波长间隔采样为1nm。再者，两个文件的波长范围应有重叠，否则将显示无重叠区，无法进行匹配计算。每个数据文件的最大波长值与最小波长值的差值应小于等于5000，但可以通过修改Cmath1Doc.cpp中的#defineMAX-POINT5000参数，改变这个值。每行的光谱数据可以是整型，也可以是浮点型，但文件中不得有波长值相同的行，否则将影响匹配的精度。2.2.3计算结果本次工作中，采用云南腾冲铀矿区和内蒙古海拉尔盆地的岩芯光谱测量数据进行实验研究。光谱匹配数据为工作区同种岩性（砂岩）不同地点的测量数据。由匹配结果(图6、图7)可以看出匹配结果较好，但同时表明：匹配区间不同，匹配结果略有差异（表2）。由此看出，匹配区的选择在整个匹配计算中非常重要。因此在光谱匹配计算中，应尽量选择地物诊断光谱特征区数据作为匹配区间，这样才能客观反映匹配结果。这正体现了高光谱遥感最佳测量波段选择的重要性。表2不同匹配区间匹配系数计算结果表Table2Differentmatchingcoefficientsofdifferentranges匹配系数匹配区间（nm）云南腾冲砂岩岩芯内蒙古海拉尔砂岩岩芯500～7480.9995720.999957350～10000.9958220.998592379～12910.9841530.992076835～13200.8717380.947737682～12380.7235820.9599858350500