基于Landsat 8數據的弓長嶺礦區遙感蝕變異常信息提取

宋 坤 王恩德 付建飛 姚玉增 孫 健

(1.東北大學資源與土木工程學院,遼寧 沈陽 110819;2.中國海洋大學海底科學與探測技術教育部重點實驗室,山東 青島 266100)

遼寧鞍本地區弓長嶺鐵礦是我國唯一由鞍山式貧鐵礦經后期熱液改造形成的大型磁鐵富礦,也是唯一具有工業開采價值的富鐵礦[1]。國內外諸多學者在弓長嶺鐵礦成因方面進行了深入研究,其找礦方法以地質、地球物理及地球化學方法為主,對于遙感方法的應用較少。自20世紀70年代以來,遙感技術以其快速、無損的優勢被廣泛應用在地質領域[2]。相對于較傳統的地質找礦方法在礦產勘探工作中的高成本、低效率等不足來說,遙感技術可在較大范圍內提取研究區蝕變礦化信息,為礦產資源潛力評估提供幫助,為找礦勘探工作提供方向[3]。本研究以弓長嶺鐵礦為例,選取Landsat8OLI遙感數據,以羥基礦物和鐵染礦物在特定波段的吸收或反射特征為理論依據,應用主成分分析法和波段比值法對遙感異常蝕變信息進行提取。通過總結并分析蝕變信息與鐵礦體之間的關系,為鞍本地區及其他同類型礦床的找礦工作提供依據。

1 研究區概況

1.1 研究區地質背景

研究區位于弓長嶺及其周邊區域,大致以弓長嶺為中心區域,向西至遼陽市東側,向東至南芬鐵礦區,呈70 km×50 km的矩形,面積約3 500 km2(圖 1)。弓長嶺鐵礦資源十分豐富,具有儲量大、分布廣等特點。研究區位于華北克拉通北緣東段的遼東地區,具有下部太古宙基底與上部古生界沉積蓋層的雙層結構。基底巖性主要由鞍山群茨溝組、大峪溝組和櫻桃園組變質巖組成[4]。我國賦存于太古宙變質巖系的鐵礦最早發現于這一帶,稱之為鞍山式鐵礦,以貧礦多、少有富礦、儲量大為主要特點。據統計,遼寧省鐵礦資源儲量占全國總量的1/4左右,且鞍山式鐵礦占遼寧省鐵礦總量的96%。

圖1 研究區區域地質特征Fig.1 Regional geological characteristics of the study area

1.2 礦床地質特征

鞍本地區以弓長嶺礦區規模最大,也最為典型。弓長嶺礦區主要為沉積變質型鐵礦,磁鐵富礦體一部分分布在各層條帶狀鐵礦床中,還有部分含鐵礦體產在蝕變圍巖或混合巖中[5]。這表明富鐵礦體可能是條帶狀磁鐵礦體經熱液改造形成[6]。礦集區呈NW—SE向展布,弓長嶺一、二、三礦區以及獨木礦區為鐵礦床的主要分布區域(圖2)。鞍山群茨溝組是該地區的主要賦礦地層,二礦區茨溝組地層出露最為完整。弓長嶺背斜為控礦構造[7]。礦區內礦石類型主要為磁鐵石英巖,其中二礦區內含有富鐵礦石(圖2)。出露的巖漿巖主要為晚太古代混合花崗巖,弓長嶺礦區呈孤島狀分布在混合花崗巖中[8]。

圖2 弓長嶺礦區地質特征Fig.2 Geological characteristics of Gongchangling mining area

2 遙感數據處理

2.1 數據源

Landsat 8衛星于2013年成功發射,衛星上搭載了OLI陸地成像儀和TIRS熱紅外傳感器,其空間分辨率和光譜分辨率等方面與之前的Landsat系列衛星基本保持一致。Landsat 8與Landsat 7衛星波段對應關系參見表1。Landsat 8遙感影像數據分辨率較高,且影像涉及范圍廣,在遙感異常信息提取中被廣泛使用;同時,OLI和ETM+數據波段分配大體一致,且OLI波段范圍比ETM+更窄,可更好地表現出蝕變礦物波的譜特征[9]。

表1 Landsat 8與Landsat 7衛星波段參數Table 1 Satellite bands parameters of Landsat 8 and Landsat 7

通過對TM和ETM+數據波段特征及其應用的分析研究,總結出OLI傳感器中不同波段的光譜特征及其適用領域,即Band 5、Band 6和Band 7 3個波段與蝕變異常提取相關,其中,對碳酸鹽以及黏土礦物等有明顯吸收特征的光譜波段為Band 7。

研究區Landsat 8 OLI圖像的成像時間為2015-04-20,衛星軌道號為119-31,經反復比對分析,該時段圖像地物輪廓清晰、層次感強,圖上無云層,地面無積雪,植被覆蓋率低,地表巖石、土壤裸露較好,成像質量最佳,有效滿足了本次研究的要求(圖3)。研究區圖像范圍為東經 123°04′40″～123°58′01″,北緯40°53′14″～41°25′54″。

圖3 研究區原始遙感影像Fig.3 Original remote sensing image of the study area

2.2 遙感影像預處理

在遙感成像過程中,遙感影像因受其搭載平臺、傳感器誤差以及大氣環境、太陽高度角、地形等不同因素影響,得到的測量值與地物本身的光譜反射率存在誤差。為了消除或修正輻射畸變,得到真實的地物反射率,須對原始影像數據進行輻射校正,通常包括傳感器校正和大氣校正兩類[10]。

對傳感器誤差的校正是通過輻射定標進行,將傳感器中的電壓或數字量化值(Digital Number,DN)換算為大氣外層的表面反射率。大氣校正是通過將大氣頂層的輻射亮度值轉換為地表反射率消除或修正誤差,從而獲取地物的真實反射率(圖4)。上述處理均使用ENVI4.7軟件的輻射定標和FLAASH大氣校正模塊,參照遙感數據的參數設置完成校正。

圖4 大氣校正前后遙感圖像及其波譜曲線Fig.4 Remote sensing images and their pop curves before and after atmospheric correction

3 研究方法

不同的礦物、巖石及礦化蝕變信息都具有其特有的光譜特征[11]。依據蝕變礦物的不同波譜特征,可以有效提取遙感蝕變異常信息[12]。通常情況下,蝕變異常信息提取的主要相關礦物離子包括鐵離子、羥基或碳酸根等離子(表2)。

表2 對巖石反射光譜特征起主導作用的離子和基團的光譜特征[11]Table 2 Spectral characteristics of ions and groups that dominate the spectral characteristics of rock reflections

分析表2可知:

(1)含鐵(Fe2+、Fe3+)離子。以含Fe3+的褐鐵礦、赤鐵礦等蝕變礦物為主,其中含Fe2+的蝕變礦物相對較少。如圖5所示,Fe3+礦物在0.40～0.55μm和0.85～0.95μm波段表現為明顯的吸收特征,對應OLI 2(0.45～0.51μm)波段、OLI 3(0.53～0.59μm)波段和OLI 5(0.85～0.88μm)波段的強吸收特征,以及在OLI 4(0.64～0.67μm)波段的強反射特征。

圖5 含鐵類蝕變礦物反射光譜曲線Fig.5 Reflectance spectral curves of iron dyeing altered minerals

(2)含羥基(OH-)或碳酸根(CO32-)離子。該類離子蝕變礦物主要為黏土礦物,如圖6所示,在1.4 μm及2.2～2.4μm波段存在吸收特征,對應OLI 7(2.11～2.29μm)波段產生低反射值,在 OLI 6(1.57～1.65μm)波段產生高反射值。含碳酸根礦物主要存在5個特征吸收譜帶,分別在2.35μm、2.55μm波段處表現出較強吸收,在1.9、2.0、2.16 μm波段處表現出相對較弱吸收,對應OLI 7(2.11～2.29μm)波段形成吸收特征。

圖6 含羥基類蝕變礦物反射光譜曲線Fig.6 Reflectance spectral curves of hydroxyl altered minerals

3.1 波段運算法

波段運算法包括波段的加減運算與比值運算。波段加法通過將多幅圖像亮度值相加求平均,達到圖像視覺增強效果、邊緣增強效果。圖像減法運算則可以去背景,突出研究對象。蝕變圍巖與周圍巖石之間的亮度值存在差異,可通過波段加減法運算增強礦化蝕變信息[13]。

為了增強蝕變信息,還可以運用波段比值法。根據每個象元在兩個不同波段上的亮度比值形成新的比值圖像[14],可消除地物反射或避免光照引起的陰影干擾,突出重要信息。

比值法通常包括簡單比值法、差和比值法以及截取比值法等。結合弓長嶺礦區實際情況,突出異常信息,本研究運用簡單比值法對圖像進行處理[15]。簡單比值法運算公式為

式中,DNm(x,y),DNn(x,y)分別為像元 (x,y)在m和n波段上的亮度值;Rmn(x,y)為輸出的亮度比值。

OLI遙感影像通常可對鐵氧化物、氫氧化物類,羥基礦物類,碳酸鹽、泥化類3類蝕變礦物進行識別[16]。在遙感異常信息提取過程中增強蝕變信息的主要OLI波段比值為:①Band 4/Band 2,主要增強鐵氧化物、氫氧化物類蝕變;②Band 6/Band 4,主要增強硅化類蝕變;③Band 6/Band 7,主要增強碳酸鹽化、泥化類蝕變。

本研究利用ENVI 4.7軟件中的“band math”工具進行上述相關波段進行運算,結果如圖7所示。

圖7 波段比值法圖像處理結果Fig.7 Image processing results by band ratio method

3.2 主成分分析法

主成分分析法(Principal Components Analysis,PCA)是較常用的蝕變信息提取方法。通過對反應蝕變信息的多個波段進行特征統計的多維正交線性變換(K-L變換),生成1組新的組分圖像[17]。主要目的是集中多波段中有用的信息到盡可能少的新組分圖像中,確保地物信息之間互不干擾再去相關性[18]。同時,選用計算結果中的新波段作為多波段原始圖像,可實現增強重要信息的目的[19]。

主成分分析法在遙感蝕變異常信息提取中應用廣泛,被稱為CROSTA法。根據前人研究成果,對于OL I數據,常用OLI 2、OLI 5、OLI 6、OLI 7這4個波段提取羥基蝕變異常,用OLI 2、OLI 4、OLI 5、OLI 6這4個波段提取鐵染遙感蝕變異常。

鐵染異常主分量的判斷準則為含鐵染異常礦物OLI 4與OLI 2及OLI 5的系數符號相反、OLI 4與OLI 6的系數符號相同[20]。同時,還可以根據不同波段在不同主分量上的載荷情況來輔助確定異常主成分(該主分量在OLI 4波段載荷較大,在OLI 2波段載荷較小)。本研究將與上述準則對應的主分量稱為鐵染異常主分量。

羥基異常主分量的判斷準則為含羥基異常礦物的OLI 6系數應與OLI 7及OLI 5的系數符號相反,OLI 2一般與OLI 6系數符號相同[21-22],且該主分量在OLI 6、OLI 7波段載荷較大。本研究將與上述準則對應的主分量稱為羥基異常主分量。

3.2.1 鐵染蝕變異常信息提取

根據上述PCA方法提取原理,對研究區圖像的OLI 2、OLI4、OLI 5、OLI6波段進行主成分分析,得到4個主分量,并對其進行統計分析,結果見表3。

表3 OLI 2、4、5、6波段提取鐵染異常特征向量矩陣及各主分量Table 3 Eigenvector matrix and eigenvalues of each principal component of iron dyeing anomaly with OLI 2,4,5 and 6

表3主要體現了4個波段的亮度信息,其特征向量值在4個波段中均為負值,不符合判斷準則。同時,根據上述鐵染異常主分量判斷準則,發現并沒有完全符合條件的主成分。但是第4主成分(PC4)最接近判斷準則,其在 OLI 4上的載荷系數(0.609 854)與在OLI 2上的載荷系數(-0.782 576)符號相反,而在OLI 5上的載荷系數(0.002 123)雖然是正,但數值很小接近于0。考慮到獲取圖像和預處理時可能會帶入一些誤差,因此可以將PC4視為滿足判斷準則的主分量。從含鐵礦物光譜特征來看,含鐵蝕變巖石在OLI 2強吸收、OLI 4強反射,因此含鐵染信息的主分量在OLI 2和OLI4顯示為強吸收和強反射且在這兩個波段上具有較大的載荷,可以看出PC4主分量滿足判斷準則,代表了鐵染蝕變信息。

通過分析PC4灰度圖像發現,與實際礦點對比發現異常區域為低值(顯示為暗),因此對PC4數據取相反數,使得鐵染異常區顯示為高值(亮),取反后的-PC4灰度圖像如圖8所示。

圖8 鐵染異常主分量DN值取反圖像Fig.8 Principal component DN value inverse image of iron dyeing anomaly

3.2.2 羥基蝕變異常信息提取

根據上述PCA方法提取原理,對研究區圖像的OLI 2、OLI5、OLI6、OLI7波段進行主成分分析,得到4個主分量,并對其進行了統計分析,結果見表4。

表4 OLI 2、5、6、7 提取羥基異常特征向量矩陣及各主分量Table 4 Eigenvector matrix and eigenvalues of each principal component of hydroxyl anomaly with OLI 2,5,6 and 7

由表4可知:根據上述羥基異常主分量判斷準則,該主分量在OLI 6、OLI 7波段載荷較大。同時,可以發現 PC4主分量符合判定準則,其 OLI 6(0.617 392)系數與OLI 7(-0.543 045)及OLI 5(-0.331 771)系數符號相反,OLI2(0.462447)與OLI 6(0.617 392)系數符號相同。PC4特征向量載荷系數絕對值較大的是OLI 6和OLI 7,分別為0.617 392和0.543 045,二者符號相反,而其他波段載荷較小,表明PC4主成分的信息主要來自OLI 6和OLI 7,與黏土礦物和碳酸鹽礦物存在有關,且與OLI 6高反射和OLI7強吸收相符,突出了泥化蝕變信息。因此,本研究選擇PC4作為羥基信息。結合研究區實際情況發現灰度圖像中亮的區域為羥基異常區(圖 9)。

圖9 羥基異常主分量圖像Fig.9 Principal component image of hydroxyl anomaly

4 弓長嶺礦區蝕變異常信息提取

4.1 鐵染異常

首先對OLI 2、4、5、6 主成分分析得到的鐵染異常主分量-PC4數據進行數學統計(表5),利用門限法分級,多次試驗確定最優分割(表6),得到如圖10所示的鐵異常分級結果。

表5 鐵染異常主分量-PC4數據數學統計信息Table 5 Mathematical statistics information of the principal component-PC4 data of iron dyeing anomaly

表6 鐵染異常主分量異常分級Table 6 Classification of principal component of iron dyeing anomaly

圖10 弓長嶺地區鐵異常圖Fig.10 Iron anomaly map of Gongchangling Area

4.2 羥基異常

對OLI 2、5、6、7主成分分析得到的羥基異常主分量PC4數據進行數學統計(表7),利用門限法分級,多次試驗確定最優分割(表8),得到如圖11所示的羥基異常分級結果。

表7 羥基異常主分量PC4數據數學統計信息Table 7 Mathematical statistics information of the principal component PC4 data of hydroxyl anomaly

表8 羥基異常主分量異常分級Table 8 Classification of principal component of hydroxyl anomaly

圖11 弓長嶺地區羥基異常圖Fig.11 Hydroxyl anomaly map of Gongchangling Area

針對弓長嶺礦區范圍進行詳細的蝕變異常分級,觀察其在小范圍內(礦區)是否能表達出更多、更有規律的信息。將主成分分析得到的鐵染和羥基異常主分量進行裁剪,只保留弓長嶺礦區圖像。對裁剪得到的圖像數據重新進行數學統計,然后利用門限法重新分級[23]。蝕變信息特征中大部分礦區以外的背景值、假異常等干擾因素可通過處理去除,結果如圖12和圖13所示。

圖12 弓長嶺礦區鐵異常圖Fig.12 Iron anomaly map of Gongchangling Mning Area

圖13 弓長嶺礦區羥基異常圖Fig.13 Hydroxyl anomaly map of Gongchangling Mining Area

4.3 分析驗證

主成分分析法對蝕變信息的提取效果比較理想,研究區內主要礦區都有異常顯示,能明顯的與非異常區分。鐵和羥基的蝕變異常區域重疊較好,也吻合了研究區內礦化區域的鐵和羥基含量均較高的特征。為了進一步驗證遙感蝕變異常信息與區域礦化信息的吻合程度,本研究將鐵和羥基蝕變異常進行疊加處理,并將研究區內已知礦點投影到異常疊加圖上,分析蝕變異常與礦點的吻合程度,結果如圖14所示。

圖14 蝕變異常疊加圖和部分已知礦點投影Fig.14 Alteration anomaly superposition map and partial known ore occurrences projection map

由圖14可知:投影的7個超大型礦床和6個大型礦床中有5個超大型礦床和4個大型礦床位于強蝕變異常區范圍內,吻合程度較好,2個大型礦床處于輕微蝕變異常區范圍內,可能是礦床圍巖蝕變作用較弱或為隱伏礦床所致。2個超大型礦床不在明顯的蝕變異常區內,吻合程度較差,分別為大臺溝礦床和思山嶺礦床。大臺溝礦床2011年詳查顯示埋深達到1 100～1 400m,含鐵巖系分布在太古宙鞍山群,上覆遼河群及近水平的青白口系、震旦系等沉積蓋層[24]。思山嶺鐵礦埋深達400 m以上,賦礦巖系為鞍山群茨溝組,上覆沉積蓋層為青白口系和震旦系巖層[25]。兩個礦床的變質蝕變信息均被沉積蓋層覆蓋,影響蝕變信息提取結果。除了以上兩個超大埋深的隱伏礦床外,其余礦床與蝕變異常區域基本吻合,證明主成分分析法能有效地提取礦化蝕變異常信息,在同類型礦床中具有較好的應用前景。

5 結 論

(1)以Landsat 8遙感數據為基礎,通過提取弓長嶺礦區蝕變礦物信息,在同類型BIF鐵礦床中建立了一種快速準確的遙感找礦勘查方法。

(2)采用主成分分析法和波段比值法分別對研究區以及弓長嶺礦區進行蝕變異常信息提取,在研究區內主成分分析法可以有效地抑制干擾信息,提取鐵和羥基蝕變信息效果較好。弓長嶺礦區鐵異常主要分布在礦區中心且向南延展。羥基異常范圍較鐵染異常面積大,且礦區南部有較多弱異常分布,其中有部分信息為地表覆蓋下的假異常。

(3)提取的蝕變異常信息與研究區內已知礦點對比發現,異常區與礦點吻合較好。蝕變信息集中區域為對找礦有指示意義的區域,反映出本研究方法對于鞍本地區及其他同類型礦床研究有一定的參考價值。但該方法在提取地表覆蓋層較厚的隱伏礦床蝕變信息過程中,對假異常的識別和蝕變弱信息的增強仍存在不足。下一步應結合鞍本地區的其他礦床信息進行補充和完善,進一步提升該方法的普適性。