基于mnf变换的蚀变信息提取方法

基于mnf变换的蚀变信息提取方法

围岩雕刻是寻找热液矿床的重要标志之一。星载热发射和反射辐射仪ASTER(AdvancedSpaceborneThermalEmissionandReflectionRadiometer)数据在0.52～0.86μm可见光长范围有三个波段(15m分辨率);在1.6～2.43μm的短波红外波段范围有6个波段(30m分辨率);在8.125～11.65μm热红外波段有5个波段(90m分辨率)(Abramsetal.,2000)。蚀变矿物在近红外、短波红外和热红外波段有吸收特征(Abramsetal.,1983;Spatzetal.,1995),ASTER在短波红外的6个波段可用于岩性和矿物的识别与探测。在ASTER多光谱中,为了突出不同矿物间的差异,一般是通过比值、主成分变换等方法来完成。如Mars和Rowan用比值法进行泥化蚀变和绢英岩化蚀变填图(Rowanetal.,2003);Crósta、Tangestani等分别用主成分分析法进行了蚀变矿物填图(Tangestanietal.,2001;Cróstaetal.,2003,2009;Marsetal.,2006;Perryetal.,2009)。目前,国内多数遥感蚀变信息提取也是主要基于这两种方法(毛晓长等,2005;丛丽娟等,2007;陈江等,2007;李培军等,2007;黄照强等,2009;杨长保等,2009;代晶晶等,2010)。对于植被发育、地形复杂的研究区来说,植被和阴影对蚀变信息提取的影响较大。比值法可以消除电磁波传输过程中乘性因子,可以很好的消除谱带吸收或反射特征较弱的背景因子(如土壤),然而它受植被、阴影、积雪和云团等影响非常严重,因此主要适用于地形起伏较小的干旱地区的信息增强。主成分分析是针对降低波段数据间的相关性提出的。它以波段间方差最大为原则,实现降维、数据压缩和信息分离的目的,具有变换前后总信息量不变的特点。对于植被发育、地形复杂地区,当噪声方差和信号方差没有明显差别时,主成分分析中的噪声分量并不随着阶次的增加而降低。而Green等(Greenetal.,1988)将改进后的主成分变换称为最大噪声分量变换(maximumnoisefractiontransformation,MNF)可以使变换后的结果按图像质量从高到低的顺序排列。在前人工作的基础上,本文提出了一种基于MNF变换和矿物标识的植被覆盖区的斑岩铜矿蚀变信息提取方法。方法的基本原理是先将已知光谱转换为图像像元,作为蚀变信息提取的矿物标识,然后对处理后的图像进行MNF变换,并根据变换后的各主成分进行RGB合成,形成突入目标的假彩色合成图像。1含矿岩体蚀变研究普朗-红山铜多金属矿带位于云南中甸地区,属于热带高原地带,植被发育,既有典型的植物覆盖区松树林,又有平原牧区的杂草和灌木,如图1所示。研究区位于义顿岛弧南段的格咱火山-岩浆弧,为典型的热液蚀变斑岩型铜矿床。夹持于娘央-夏隆瓦和弥里躺断裂之间,沿中甸岛弧的复背斜及次级背斜、亚拉夏向斜、红山背斜、冈错向斜、普朗向斜、热绒背斜轴部发育。斑(玢)岩的产出受区域构造控制明显,岩体数量大小多达46个,呈岩株、岩枝和岩墙状产出,总体呈NNW向延伸。其中较大者从南到北有懒中、普郎、松诺、萨决和促纳卡等岩体(曾普胜等,2006;尹光侯等,2009;李文昌等,2009),可以作为典型的植被覆盖区蚀变信息提取方法验证点来讨论。岩体地表露头较大,斑岩型矿化发育。含矿岩体蚀变分带明显,矿床的斑(玢)岩体的蚀变与国内外常见的斑岩型矿床基本相似,故按通常的分带方式划分为钾化硅化带、绢英岩化带、泥化蚀变带和青盘岩化带如图2,中心是钾化硅化带如钾长石、黑云母等,其周围有多种蚀变带(曹殿华等,2006;余海军等,2009;王安建等,2009)。这些含水蚀变矿物具有诊断吸收光谱。如宽的绢英岩化带,通常是由于黄铁矿的氧化产生的褐铁矿化,可以由伊利石或者白绢云母的诊断光谱判断。窄的泥化带可以由高岭土和明矾石的波谱判断。外围的青磐岩化由于围岩组份的变化,其矿物组合变化较大,但是常见组合是绿帘石、绿泥石和碳酸岩等矿物。围岩和侵入岩都会产生铜矿化并且有热液蚀变。由于青磐岩化的矿物组合变化较大,而钾化硅化带很少暴露在外,因此对于ASTER数据的分析主要集中在绢英岩化和泥化蚀变带。绢英岩化蚀变光谱特征包括伊利石和白云母的光谱,其Al-OH键在2.20μm(ASTER波段6)附近有强吸收特征,在2.38μm(ASTER波段8)附近有次吸收特征。ASTER数据可以探测到由于Fe2+置换Al从而使得低值向长波波段移动(Rowanetal.,2003)。Fe3+的强吸收使得矿物(如褐铁矿)在可见光下的反射率较低。泥化岩石的Al-OH键在2.20μm附近也有吸收特征,但是高岭石和明矾石的光谱特征和白云母及伊利石的光谱特征明显不同(Lowelletal.,1997;Rowanetal.,2003)。高岭石在2.17μm附近有第二特征,明矾石在2.17μm有吸收(ASTER波段5),在2.20μm有小吸收特征,这两种矿物在2.38μm附近均有亚吸收特征。ASTER数据在短波红外的波段设置可以识别这些波段差异(Rowanetal.,2003),如图3。2方法2.1大气上界反射率函数aster本文使用ASTER数据为L1B产品,成像时间为2005年11月30日,过境时间为4点1分。影像覆盖区域为云南中甸地区,覆盖范围为北纬27.632255°～28.275776°,东经99.222843°～99.994819°ASTER图像数据的大气上界视反射率,是图像辐射校正后得到的基本参数。假设天空辐照度各向同性和地面朗伯面反射,并忽略大气的折射、湍流和偏振,由遥感方程,地物表面反射辐亮度是太阳辐照度和表面反射率的函数。本次获取的ASTER数据是1A级数据,通过头文件发现所得的为辐射定标数据,像元值为辐亮度值,可直接转换成大气上界反射率。公式为:RTOA=(π×Lrad×d2)/(ESUNi×cos(z))(1)其中π=3.14159,RTOA是大气上界反射率,Lrad是辐亮度值,ESUNi是大气顶部的平均太阳辐射照度(对应各波段,可由通道响应函数计算得出)。z为太阳天顶角(z=90°-solarelevationangle),可以在ASTER的头文件中查到。d是日-地距离(天文单位),可以用公式2计算(Achardetal.,1994;Evaetal.,1998):d=(1-0.01672×cos(RADIANS(0.9856)×(JulianDay-4)))(2)JulianDay为卫星过境日期的儒略日。根据1∶20万地质图对ASTER数据进行几何纠正,误差控制在1个像元。2.2aster特征分析在提取信息结果的验证方面,以前都是利用地面露头或者实验室内观测波谱进行验证,本文通过将波谱库的标识矿物的波谱采样为像元波谱作为标识,可以快速验证不同信息提取方法的有效性。已知目标识别矿物为绢云母、高岭石和明矾石。根据ASTER波段响应函数对标准波谱数据重采样,获得ASTER传感器成像波段的地表反射率,如公式(3)。r=∫λ2λ1fλrλdλ∫λ2λ1fλdλ(3)r=∫λ1λ2fλrλdλ∫λ1λ2fλdλ(3)其中r为波段反射率,rλ为对应波长地物反射率,fλ为传感器的波段响应函数在波长λ处的响应,波长范围为λ1到λ2。在ASTER遥感影像上构建新像元,将标准矿物波段反射率转换为图像像元波段反射率,如公式(4)。xni=si(4)xni=si(4)其中xni为ASTER影像第i个像元在第n个波段上的反射率,si矿物在ASTER标准波段的反射率值。新像元为纯净像元,为蚀变识别的参考端元,如图4。2.3mnf变换结果MNF变换是由Green等于1988年提出的(Greenetal.,1988),目的是为了解决主成分分析中存在的图像质量不随主成分数值的增加而递减的问题。假设当前图像由p个波段组成,像素值Z(x)由信号S(x)与噪声N(x)组成,公式为Z(x)=S(x)+N(x)(5)Ζ(x)=S(x)+Ν(x)(5)其中,x为坐标;ZT(x)={Z1(x),Z(x),...,Zp(x)},S(x)与N(x)分别为互不相关的信号和噪声分量。因此,Cov{Z(x)}=∑=∑δ+∑N。其中∑δ和∑N分别表示S(x)与N(x)的协方矩阵;定义第i波段的噪声分量为Var{Ni(x)}/Var{Zi(x)},也就是该波段的噪声方差与波段方差的比值。MNF变换可以用线性公式表示为Yi(x)=aTiZ(x),i=1,...,p(6)Yi(x)=aiΤΖ(x),i=1,...,p(6)上述公式在求解过程中满足如下条件:在于Yi(x)(j=1,2,...,i)正交的线性变换过程中,噪声分量保持最大。运用主成分的推导过程,可以得出ai是正规化后的矩阵∑N∑-1的左侧特征向量,满足ai∑ai=I;与ai相对于的特征值μi是Yi(x)中的噪声分量,满足μ1≥μ2≥...μp。因此,MNF变换后的成分所代表的图像质量将逐渐提升。遥感影像的2.0～2.4μm波段属于中红外波段,可用于区分主要岩石类型。水体反射率总体较低(4%～5%),并随波长的增大逐渐降低,0.6μm处约2%～3%,过了0.75μm,水体几乎成为全吸收体。由于植被叶细胞结构的影响,其在中红外波段受到绿色植物含水量的影响,吸收率大增,反射率下降。因此在近红外的遥感影像下植被和水体都是暗目标,其反射率都非常低。由以上分析,将ASTER图像数据短波红外波谱范围(1.65～2.395μm)的6个波段进行MNF变换,其结果将干扰信息降低到最小,岩性信息提高到最大。从图5可以明显看出,前三个波段的特征值远远大于后面三个波段的特征值,前三个波段显示的地物影像清晰,集中地物的绝大部分光谱信息。后面三个波段,包含地物光谱信息减少,噪声增大。由标识矿物可以看出,目标矿物是亮色的。根据图6可以看出MNF2主要反映绢云母蚀变,MNF3主要反映高岭石蚀变,MNF4主要反映明矾石蚀变。3主成分分析及数据处理为了更进一步的验证分析本方法对植被发育、地形复杂地区的有效性,下面将分别采取主成分分析法和比值法提取泥化蚀变,将MNF变换的结果与之进行对比。主成分分析法是利用一个特征向量算子作用于输入图像,图像中的方差信息进行重新分配,利用其不同主成分浓缩信息的侧重不同,通过对特征向量载荷的判别分析,从某一个或某几个主成分中得到特征信息,它的基础是输人图象的相关性大小。Crósta和Moore在1989年用主成分分析法提取了TM数据的氧化铁(Cróstaetal.,1989)。随后又将其成功应用在ASTER数据中,用ASTER1,4,6,7波段填高岭石,ASTER1,3,5,6填白/绢云母(Cróstaetal.,2003;Marsetal.,2006)。进行主成分分析后,分析特征向量矩阵,确定包含目标矿物信息的主成分。其判断标准由Loughlin提出(Loughlinetal.,1991):ASTER波段的特征向量对主成分有高负载荷且具有相反的贡献值。对实验区的ASTER数据1,4,6,7波段做主成分分析后,其特征向量矩阵如表1所示。从图3分析可知,高岭石在ASTER4和ASTER7有高反射,在ASTER1和ASTER6有强吸收,从表1可以看出,波段1对PC3有高负载荷,为负贡献(-0.672),而波段4有高负载荷,为正贡献(0.730),因此PC3符合条件,为高岭石矿物信息增强图像。PC3图像中的亮目标包含高岭石,其高岭石丰度图如图7(b)所示。对实验区的ASTER数据1,3,5,6波段做主成分分析后,其特征向量矩阵如表2所示。由表可以看出,波段6对PC4有高负载荷,为负贡献(-0.872),而波段5有高负载荷,为正贡献(0.409),因此PC4符合条件,为白/绢云母矿物信息增强图像,其矿物丰度图如图8(b)所示。波段比值法作为一种提取波谱信息的有效手段,广泛应用于地质遥感中(Rowanetal.,1974,1977)。比值法以地物波谱的斜率为依据,增大波谱对比度,突出特征信息。光谱指数可以用于比值对提取铁氧化物、碳酸盐和铁镁质矿物等。根据研究表明,(b4+b6)/b5可以提取高岭石蚀变,(b5+b7)/b6可以提取白/绢云母蚀变。由此得到的高岭石矿物丰度图如图7(c)所示,图8(c)为白/绢云母蚀变矿物丰度图。图7结果显示,PCA(ASTER1,4,6,7)的PC3和MNF变换的MN3获得的高岭石丰度图结果很类似,比值(b4+b6)/b5获得的丰度图较差;从图7我们可以看出,MNF变换的MNF2的丰度图结果最好,而PCA(ASTER1,3,5,7)的PC4有噪声,比值(b5+b7)/b6噪声最大。由于研究区位于热带高寒区,受植被、阴影、积雪和云团等影响严重,这些干扰因素对蚀变信息提取的影响较大。波段比值法是根据代数运算的原理,当波段间差值相近但斜率不同时,利用反射波段与吸收波段的比值处理增强各种岩性之间的波谱差异,但是干枯树叶纤维素的吸收特征和含Al-OH、Fe和Mg-OH的矿物类似,其在2.10(ASTERband5)和2.30μm(ASTERband8)附近有吸收特征,因此用比值法容易混淆干枯树叶和矿物。主成分分析能够把原来多波段影像中的有用信息集中到数目尽可能少的特征影像组中,并使不同波段的影像互不相关,从而达到减少数据量的效果。主成分分析对噪声比较敏感,主成分变换后的影像并不能消除一定的噪声影响,如图7、图8所示,而比值法明显受噪声影响较大。而MNF变换不但可以将目标等有效信息集中到尽可能少的低维数据中,而且可以有效地分离噪声。由于植被和水体在短波红外波段反射率极低,可以视为黑体,本文只将ASTER短波红外波段进行MNF变换,因此效果良好。通过以上分析,将MNF2、MNF3和MNF4的丰度图叠加与ASTER3上,得到普朗-红山成矿带绢英岩化和泥化蚀变信息蚀变图,如图9所示。可以发现普朗-红山区整体表现为岩体蚀变信息较强,泥化蚀变和绢英岩化蚀变整体呈NNW向延伸,受构造控制明显。蚀变在区域上具有反映成矿背景的明显特征。实验区蚀变遥感异常展布方向、组合类型等均显示与测区区域地质背景及成矿系列有着密切关系。蚀变遥感异常展布总体呈NW-NNW向,与实验区总体构造方向一致。研究区位于高山侵蚀区,绢英岩化蚀变岩石抗侵蚀能力较弱,是一种最常见的重要蚀变,粘土化强度不大。异常带与矿化产出的曲嘎寺组、图姆沟组等地层岩性具一致展布性,且多表现为与岩体相伴出露,反映了较强的构造岩体蚀变矿化背景。蚀变遥感异常强浓度区往往分布在地表有岩体出露或部分出露地带,可能是浅埋或出露岩体风化较强、地表成矿元素富集导致岩石蚀变所致;而蚀变遥感异常强度弱的区域落在隐伏岩体区域。故推测蚀变遥感异常强度可反映岩体埋深。通过蚀变遥感异常点群与研究区已有矿床(点)产出进行吻合度分析发现,矿化蚀变遥感异常与岩石蚀变带吻合率高,可以将蚀变遥感异常视为矿化蚀变类型的直接指示标志之一。其中蚀变遥感异常与已知的矿(化)点信息对比,与斑岩体铜矿有关的绢英岩化带及泥化带无一例外都落在羟基异常一级或多级异常重叠区,且羟基异常往往浓度强、面积大且集中分布。蚀变遥感异常可较好地缩小地球化学异常面源靶区,可起到优化化探异常的作用。经地表地质调查证实,蚀变遥感异常发育的地带,往往都落在了1:5万地球化学异常图区内,且由于蚀变遥感异常具有点源特点,可较好地缩小地球化学异常的面源靶区范围,故在测区可用蚀变遥感异常这一点源信息进行近矿指示。研究区蚀变遥感异常部分地区出现蚀变遥感异常不连续、局部地区与矿化出现不吻合等的现象,这从干扰因素分析及地形地貌条件分析中可以找到答案,故蚀变遥感异常的认识需借助于其它因素进行全面研究。研究区蚀变遥感异常多出现在向阳坡,在阴影区则蚀变遥感异常缺失,造成蚀变遥感异常分布陡然断开、连续性差。究其原因,应是由于测区处于地形深切割区,植被发育,为蚀变遥感异常提取过程中最重要的干扰因子,从主成分分析法和比值法可以看出。MNF变换