白 洋，宋唐雷，賈玉娜，李孟倩，康會濤

(1.華北理工大學礦業(yè)工程學院，河北 唐山 063210；2.河北省地質(zhì)礦產(chǎn)勘查開發(fā)局第七地質(zhì)大隊，河北 三河 065201)

尾礦庫是礦山開采過程中形成的最大危險源之一，嚴重威脅著尾礦庫周圍的生態(tài)環(huán)境及周邊居民的生命財產(chǎn)安全[1-3]。因此，實時、快速、高效地掌握尾礦庫的分布情況，對提高礦產(chǎn)資源循環(huán)利用能力具有十分重要的現(xiàn)實意義。傳統(tǒng)的尾礦庫監(jiān)測方法多以實地調(diào)查、勘測為主，耗費大量的人力、物力，且監(jiān)測結(jié)果存在一定的主觀性和局限性，無法滿足大面積尾礦監(jiān)測。隨著遙感技術(shù)在空間分辨率、光譜分辨率和時間分辨率方面的顯著提升[4-6]，利用遙感技術(shù)探測范圍廣、數(shù)據(jù)獲取快捷、受自然條件約束小等優(yōu)勢，可大大提高尾礦庫監(jiān)測的成效[7-8]。閻永忠等[9]利用遙感技術(shù)和GIS技術(shù)系統(tǒng)對尾礦庫進行監(jiān)測分析；申彥科等[10]依據(jù)TM影像總結(jié)出尾礦庫增容會對周圍的生態(tài)環(huán)境產(chǎn)生的影響；方雪娟等[11]對陳貴鎮(zhèn)小型尾礦庫周圍的4類地物進行了危害分級分析；南竣祥等[12]對宣威市煤礦開采情況展開遙感解譯，提升了尾礦庫遙感監(jiān)測技術(shù)的作用；張云英等[13]利用高分一號影像對尾礦庫進行目視解譯最佳波段組合，為運用遙感影像進行尾礦庫研究提供了依據(jù)；凌子燕[14]利用GF-1遙感衛(wèi)星，采取支持向量機(support vector machine，SVM)建立遙感反演模型，提高了尾礦庫遙感監(jiān)測精度；曹蘭杰等[15]對典型地區(qū)鐵尾礦樣本進行了光譜測量，并利用SAM(spectral angle mapper)方法和SFF(spectral feature fitting)方法進行匹配，優(yōu)化了河北省兩種礦床類型鐵尾礦的識別窗口；郝利娜等[16]從尾砂的光譜、紋理特征等因素考慮，根據(jù)WorldView-2遙感影像對鄂東南尾礦庫規(guī)模進行了尾礦識別因子研究；馬國超等[17]集合高分遙感衛(wèi)星數(shù)據(jù)和無人機遙感數(shù)據(jù)，利用其快速排查，精準監(jiān)控的優(yōu)勢，極大提高了礦山安全監(jiān)測效率；汪金花等[18-19]對比分析了差異化粒徑、干濕狀態(tài)下的鐵尾砂樣本光譜特征，擬合確定了鐵尾礦高光譜識別的有效窗口，經(jīng)光譜特征參量及光譜匹配分析獲得了鐵尾礦特征識別波段，得到了鐵尾礦多源信息提取的可行性。多源遙感數(shù)據(jù)具有一定的優(yōu)勢，但同樣存在局限性，如多光譜數(shù)據(jù)空間信息豐富但無法區(qū)分“同物異譜”或“同譜異物”，而高光譜數(shù)據(jù)光譜信息豐富但數(shù)據(jù)量大，信息冗余，存在Hughes現(xiàn)象。因此，如何利用高光譜數(shù)據(jù)和多光譜數(shù)據(jù)的優(yōu)勢進行尾礦庫信息提取尤為重要。

鑒于此，本文以唐山市司家營尾礦庫為例，利用實測尾砂波譜信息分別與多光譜Landsat8 OLI影像數(shù)據(jù)和高光譜珠海一號影像數(shù)據(jù)的端元波譜進行匹配，并采用光譜角分類的方法進行尾礦信息填圖，獲取尾礦庫的位置信息。

1 研究區(qū)概況與數(shù)據(jù)源

本文研究區(qū)域是河北省唐山市灤州市司家營尾礦庫，位于灤州市城南10 km處，屬于高硅鞍山型鐵尾礦。礦體南北長10 km，東西2 km，礦石資源總儲量為23.48億t，屬于我國整裝特大型鐵礦床之一。司家營尾礦庫距選礦廠約12 km，尾礦庫匯水面積為2.5 km2，平均坡度4.1°。尾礦庫由東、西兩個溝組成，分別建有初期壩和副壩。尾礦堆積壩最終堆積標高為128 m，總壩高93 m，總庫容約8 714萬m3[20]。

2 數(shù)據(jù)源介紹與預處理

2.1 數(shù)據(jù)源介紹

Landsat8衛(wèi)星于2013年2月發(fā)射成功，是目前世界范圍內(nèi)應用領域最廣的民用對地觀測衛(wèi)星，本次數(shù)據(jù)選用采集于2019年10月29日的一幅OLI影像，共9個波段，包括8個30 m的多光譜波段和1個15 m的全色波段。珠海一號高光譜衛(wèi)星于2018年4月發(fā)射成功，其采集的數(shù)據(jù)空間分辨率達到10 m，光譜分辨率達到2.5 nm，共256個波段，均在400～1 000 nm之間。本文使用的60個尾砂樣本均取自河北省灤州市司家營尾礦庫，由ASD Field Spec Pro光譜儀對各樣本進行測量，每個樣本均測量5次，每5 min測一次，白板進行校正，最終得到60個樣本的平均波譜曲線。

2.2 數(shù)據(jù)預處理

根據(jù)數(shù)據(jù)源的具體情況，對Landsat8 OLI影像數(shù)據(jù)進行了輻射定標、大氣校正等預處理，提高了影像的定位精度、空間分辨率及光譜分辨率，有利于進行尾礦庫的精準識別。以研究區(qū)的Landsat8 OLI影像文件和珠海一號影像文件為基準，對去除包絡線后的波譜曲線進行重新采樣，確保實測波譜與影像波譜范圍相同，方便后續(xù)開展實測波譜與端元波譜的匹配。

3 光譜匹配實驗

影像數(shù)據(jù)經(jīng)最小噪聲分離(minimum noise fractions，MNF)達到降維和去噪目的，通過MNF特征值曲線選擇信息量較大的幾個分量，進行純凈像元指數(shù)(PPI)計算，通過PPI求得的端元信息，在N維可視化中找到緊密抱團的端元，并繪制其特征波譜曲線，經(jīng)與實測的尾礦自定義波譜庫進行匹配，并進行尾礦庫填圖。

3.1 端元波譜提取

端元波譜是地物識別分類中十分重要的參考光譜，端元波譜的質(zhì)量直接影響地物識別結(jié)果的精度。預處理后的影像進行MNF變換后可以得到各波段的特征信息，如圖1和圖2所示。從圖1和圖2中可以看出沿X軸方向各特征數(shù)包含的特征信息逐漸遞減。圖1為Landsat8 OLI影像波段經(jīng)過MNF變換后的特征圖，其中前3個分量的特征值比較大，第4分量以后特征曲線的特征值幾乎都在5以下，因此，選取前3個特征波段進行PPI計算。圖2為珠海一號影像的特征曲線，其中前4個波段涵蓋了大部分的圖像信息，從第5分量開始曲線基本保持平穩(wěn)，后面波段的影像信息幾乎全部為噪聲，所以后續(xù)PPI計算時只選取了前4個波段，以大幅度減少處理工作量和處理時間。

圖1 Landsat8 OLI MNF變換特征曲線Fig.1 MNF transformation characteristic curve of Landsat8 OLI

圖2 珠海一號MNF變換特征曲線Fig.2 MNF transformation characteristic curve of Zhuhai-1

純凈像元指數(shù)(pixel purity index，PPI)是在數(shù)據(jù)降維處理后進行純凈像元提取的過程，PPI越大，說明與之對應的像元越純凈，對應的PPI結(jié)果圖就越亮。在影像進行MNF變換之后，經(jīng)過多次實驗對Landsat8 OLI影像進行PPI計算采用了10 000次迭代，閾值系數(shù)設置為2.5，珠海一號影像迭代次數(shù)為10 000次，閾值系數(shù)設置為2.5，如圖3和圖4所示，其中，高亮區(qū)域為PPI提取后的極值像元的分布位置，這些純凈像元多聚集在尾礦庫中央的位置，還有一些零散的分布在尾礦庫的外圍。

圖3 Landsat8 OLI影像PPI結(jié)果圖Fig.3 PPI results of Landsat8 OLI images

圖4 珠海一號影像PPI結(jié)果圖Fig.4 PPI results of Zhuhai-1 images

將MNF和PPI計算后的結(jié)果輸入到N維可視化中，選取不同的波段進行旋轉(zhuǎn)，找到多個獨立的點集群，如圖5(a)和圖6(a)所示，圖中有明顯聚集的4個團簇，這4個團簇代表了4種“純凈像元”，以此定義為端元，圖5(b)和圖6(b)為得到的4種端元的波譜曲線，通過現(xiàn)場采集的尾砂樣本制作的自定義波譜庫可以識別出這4種像元代表的物質(zhì)和分布的區(qū)域。

圖5 Landsat8 OLI影像N維可視化圖Fig.5 N-dimensional visualization of Landsat8 OLI images

圖6 珠海一號影像N維可視化圖Fig.6 N-dimensional visualization of Zhuhai-1 images

3.2 光譜匹配算法

在端元提取出來之后，還不能明確各端元代表哪種地物類型，需要通過波譜分析對端元進行識別，同時要修改波譜曲線分辨率的縮放系數(shù)，使端元波譜與波譜庫具有相同的分辨率，此次研究分別使用光譜角度匹配(spectral angle match，SAM)和波譜特征擬合(spectral feature fitting，SFF)2種方法進行識別。

3.2.1 光譜角度匹配法

光譜角度匹配(spectral angle match，SAM)是將待驗證光譜和參考光譜均作為高維向量，通過計算二者的夾角大小來確定匹配程度[21]，夾角越小說明相似程度越高，是同類地物的可能性越大。

計算夾角通常通過兩個向量的余弦來求得，公式見式(1)。

(1)

式中：X*為像元光譜的向量；Xi為參考光譜的向量。

3.2.2 波譜特征擬合法

波譜特征擬合(spectral feature fitting，SFF)是一種利用波譜曲線的吸收特征對像元波譜和參考波譜進行匹配的方法，使用最小二乘法逐次迭代以評價兩者相似程度[22]。SFF會對輸入的像元波譜進行包絡線去除、刪去背景信息、增強特征波段的吸收特性，比較吸收谷的位置和深度，將參考波譜與像元波譜相匹配。

4 匹配結(jié)果與分析

4.1 端元匹配結(jié)果

N維可視化后的端元波譜需要通過波譜分析與已知物質(zhì)的波譜進行匹配，方可得知哪一個端元最能代表需填圖的地物類型，并用最終得分的高低來衡量匹配的可信度，分值越高匹配度越好。本研究使用兩種方法參與權(quán)重分析，分別是光譜角度匹配(SAM)和波譜特征擬合(SFF)，表1為Landsat8 OLI數(shù)據(jù)得到的4種端元與實測光譜的識別結(jié)果。從表1中可以看出n-D Class #1與尾礦實測波譜的匹配程度最高達到1.618，因此，確定端元n-D Class #1為尾礦波譜。表2為珠海一號影像數(shù)據(jù)得到的4種端元與實測光譜的識別結(jié)果。從表2可以看出n-D Class #1與尾礦實測波譜的匹配程度最高達到2.069，因此，確定端元n-D Class #1為尾礦波譜。從綜合波譜匹配的結(jié)果來看珠海一號影像數(shù)據(jù)的匹配精度要優(yōu)于Landsat8 OLI影像。

表1 Landsat8 OLI影像端元波譜與實測波譜的匹配得分表Table 1 Matching score of end element spectrum and measured spectrum of Landsat8 OLI images

表2 珠海一號影像端元波譜與實測波譜的匹配得分表Table 2 Matching scores of end element spectrum and measured spectrum of Zhuhai-1 images

4.2 尾礦識別填圖

針對Landsat8 OLI影像數(shù)據(jù)進行尾礦地物識別，確定影像分類閾值分別為0.10、0.13、0.15，圖7中深色區(qū)域為識別出的尾礦庫，其中閾值為0.10時，在非尾礦庫區(qū)域識別出的深色區(qū)域最少，即錯分最少；閾值為0.15時，尾礦位置的識別面積最大，此時識別尾礦位置的像元數(shù)量最多。針對珠海一號影像數(shù)據(jù)識別的結(jié)果如圖8所示，最優(yōu)閾值分別為0.28、0.29、0.30，都取得了較為精準的識別結(jié)果，基本覆蓋整個尾礦庫。對比Landsat8 OLI影像和珠海一號影像數(shù)據(jù)尾礦識別結(jié)果，提取的尾礦空間分布基本一致，珠海一號影像數(shù)據(jù)填圖結(jié)果優(yōu)于Landsat8 OLI的填圖結(jié)果。

圖7 Landsat8 OLI數(shù)據(jù)尾礦識別填圖Fig.7 Tailings identification mapping of Landsat8 OLI data

圖8 珠海一號數(shù)據(jù)尾礦識別填圖Fig.8 Tailings identification mapping of Zhuhai-1 data

5 結(jié) 論

本文對實測高光譜數(shù)據(jù)重采樣后，利用多源遙感影像進行尾礦庫識別研究，得到結(jié)論如下所述。

1) 依據(jù)影像波段范圍對實測波譜進行重采樣，在降低數(shù)據(jù)維數(shù)的同時，還能確保在相同波譜范圍開展波譜匹配實驗。

2) 對比高光譜數(shù)據(jù)尾礦識別結(jié)果，Landsat8 OLI影像能夠有效識別尾礦信息且空間位置基本一致，從精度來說，高光譜影像識別尾礦庫的解譯效果更接近真實地物。