柳炳利, 張 立, 黃昕蕾, 郭 科, 李 程, 趙 麗, 王玉蘭, 陳 聆
(1.中國地質(zhì)科學(xué)院 地球物理地球化學(xué)勘查研究所,河北 廊坊 065000;2.成都理工大學(xué) 數(shù)學(xué)地質(zhì)四川省重點實驗室,成都 610059)
?
基于小波變換的巖心高光譜特征提取及波譜匹配方法
柳炳利1,2, 張 立2, 黃昕蕾2, 郭 科2, 李 程2, 趙 麗2, 王玉蘭2, 陳 聆2
(1.中國地質(zhì)科學(xué)院 地球物理地球化學(xué)勘查研究所,河北 廊坊 065000;2.成都理工大學(xué) 數(shù)學(xué)地質(zhì)四川省重點實驗室,成都 610059)
探討巖心礦物高光譜特征提取和快速波譜匹配的方法。將小波變換應(yīng)用于巖心高光譜的特征提取,以多尺度小波系數(shù)能量值為特征表示方法,以光譜角度匹配為礦物波譜匹配方法。針對湖北省大冶市雞冠咀銅金礦的10條實測光譜,采用光譜角度匹配方法對實測波譜與標(biāo)準(zhǔn)波譜的多尺度小波系數(shù)能量值進(jìn)行匹配,得到的礦物匹配結(jié)果與鉆孔柱狀圖中的實際礦物完全吻合,驗證了該方法的可靠性。多尺度小波系數(shù)能量值可以作為巖心高光譜的波譜特征,光譜角度匹配法在巖心高光譜的波譜匹配中同樣具備適宜性。
巖心高光譜;小波變換;特征提取;波譜匹配
隨著便攜式波譜儀的推廣,巖礦光譜測試工作正逐漸為大多數(shù)遙感地質(zhì)工作者所熟悉,并形成除航天遙感、航空遙感之外的“近感”遙感[1],由此將巖礦光譜的應(yīng)用引向更為廣闊的領(lǐng)域,其中典型的代表就是巖心高光譜測量在礦床勘探中所起到的重要作用[2]。
特征提取是利用巖心高光譜進(jìn)行礦物識別的關(guān)鍵技術(shù)之一,特征提取的實質(zhì)是從海量的原始數(shù)據(jù)中抽取出最能描述該模式的少部分特征,是一種能有效解決“維數(shù)災(zāi)難”問題的預(yù)處理方法[3]。對輸入特征空間進(jìn)行某種變換或者組合,然后產(chǎn)生新的特征空間的過程就是特征提取。國內(nèi)外眾多學(xué)者在特征提取的研究中取得了成果,如Benediktsson[4]將邊界特征提取技術(shù)用于遙感圖像的分類。Nakariyakul[5]在農(nóng)業(yè)生產(chǎn)分析中應(yīng)用了高光譜數(shù)據(jù)的特征提取方法,Du[6]對光譜數(shù)據(jù)的編碼技術(shù)進(jìn)行了細(xì)致的研究。Li Jiang[7]采用小波方法進(jìn)行特征提取,在混合像元分解與端元提取方面取得了較好的效果。Kuo[8]對高光譜的圖像降維時,應(yīng)用了特征提取方法。
湖北省大冶市雞冠咀銅金礦床屬長江中下游地區(qū)Fe、Cu、Au成礦帶西段的矽卡巖型礦床,發(fā)育有17種主要礦物和6種次要礦物。為了利用巖心高光譜進(jìn)行準(zhǔn)確的礦物識別并進(jìn)行深部礦產(chǎn)預(yù)測,使用便攜式地物光譜儀ASD FieldSpec?4對巖心進(jìn)行了高光譜測量,取得了大量的波譜數(shù)據(jù)。同時,在礦區(qū)選取單體較大、相對純凈的礦物進(jìn)行測量,取得的波譜曲線作為標(biāo)準(zhǔn)波譜建庫,為實現(xiàn)有效的波譜匹配做準(zhǔn)備。
在研究同一個問題時,如果使用不同的小波基函數(shù)必然會產(chǎn)生不同的結(jié)果。小波基函數(shù)可以有很多種,在光譜數(shù)據(jù)處理中,選擇不同的小波基處理同一數(shù)據(jù),得到的結(jié)果也許就不相同,因此,尋找適用于巖心高光譜處理的最優(yōu)小波基顯得尤為重要。
在實際數(shù)據(jù)處理中可供選擇的小波基一般有Daubechies(dbN)系列、Symlets(symN)系列、Mexican小波、Morlet小波等。這些小波都是具有一定緊支撐性、平滑性以及對稱性的正交小波[9-10]。
最優(yōu)小波基的選取,是根據(jù)Mallat算法對原始信號f1僅作單尺度小波分解與重構(gòu)處理,然后計算重構(gòu)信號f2與原始信號f1之間的誤差。誤差計算公式為
其中N為信號的采樣點數(shù)。
針對湖北雞冠咀銅金礦26號勘探線鉆孔KZK10的光譜曲線在不同小波基下分解與重構(gòu)的誤差進(jìn)行分析。
據(jù)表1與表2發(fā)現(xiàn),使用db4小波基進(jìn)行小波分解和重構(gòu)的誤差最小,因此在處理巖心高光譜數(shù)據(jù)時選擇db4作為最優(yōu)小波基,并將分解層數(shù)確定為3層。
2.1 小波特征提取流程
特征提取可以通過小波分解與重構(gòu)來實現(xiàn),圖1為小波分解的示意圖。
小波變換通過特征空間轉(zhuǎn)換來研究光譜信號,光譜信號在新“坐標(biāo)系”內(nèi)的投影就構(gòu)成了小波系數(shù),可以通過選取小波系數(shù)的某些統(tǒng)計信息作為特征(如對高頻、低頻系數(shù)求平均值、絕對平均值等),這種方法得到的特征維數(shù)很低并且易于實現(xiàn)。圖2是小波特征提取的流程圖。
表1 sym系列小波分解重構(gòu)之后的誤差
表2 db系列小波分解重構(gòu)之后的誤差
圖1 小波分解示意圖Fig.1 Diagram showing wavelet decomposition
圖2 小波特征提取流程圖Fig.2 Flow chart showing extraction of wavelet characteristics
2.2 巖心高光譜特征提取
地物光譜的吸收特征,必須設(shè)定一些衡量其特征的參數(shù)。通常情況下,吸收位置(P)、吸收深度(H)、吸收總面積(A)、光譜絕對反射值(R)、吸收寬度(W)、斜率(K)以及對稱度(S)等就能描述該地物光譜的特征(圖3)。大冶銅金礦巖心高光譜常見礦物曲線如圖4所示。
圖3 地物光譜特征示意圖Fig.3 Sketch showing spectral characteristics of ground objects
圖4 常見蝕變礦物光譜曲線Fig.4 Spectral curves of common altered minerals
在對大冶銅金礦巖心高光譜的特征提取過程中,選取db4為小波基函數(shù),分解層數(shù)確定為3層,將原始信號X分解成低頻部分CA1與高頻部分CD1,將低頻部分CA1繼續(xù)分解為低頻部分CA2與高頻部分CD2,而高頻部分不再分解;最后將低頻部分CA2分解為低頻部分CA3與高頻部分CD3。得到的原始信號就是第3層的低頻部分和每層的高頻部分的疊加,其數(shù)學(xué)表達(dá)式為
X=CA3+CD3+CD2+CD1
以各尺度能量值作為小波特征向量(能量值為各尺度小波系數(shù)模的平方),分別求取標(biāo)準(zhǔn)波譜與實測波譜的特征如表3和表4。
3.1 波譜匹配方法
光譜匹配是根據(jù)2條光譜曲線的相似程度來判別礦物類型[11]。常用的光譜匹配方法主要有相關(guān)系數(shù)法、光譜角度匹配法、距離函數(shù)法(歐幾里得距離和馬哈拉諾比斯距離等)、最小距離匹配法等[12-13]。
本次采用的是光譜角度匹配法(spectral angel mapper, 簡稱SAM)。1992年,Kruse等人提出了光譜角度匹配的概念[14]。通過計算待判光譜向量和標(biāo)準(zhǔn)光譜向量之間的廣義夾角,光譜角越小說明這2條光譜的相似程度越高。
假設(shè)有一個待判光譜數(shù)據(jù)為X=(x1,x2, …,xn)T,其對應(yīng)的標(biāo)準(zhǔn)光譜數(shù)據(jù)Y=(y1,y2, …,yn)T,該廣義夾角的數(shù)學(xué)表達(dá)式為
表3 標(biāo)準(zhǔn)波譜特征提取結(jié)果
表4 鉆孔KZK10部分礦物光譜曲線的特征提取結(jié)果
3.2 巖心波譜匹配
在分別計算標(biāo)準(zhǔn)波譜與實測波譜的能量特征向量的基礎(chǔ)上,使用SAM方法計算特征向量的夾角,依據(jù)夾角最小原則進(jìn)行匹配。并將匹配結(jié)果與KZK10鉆孔柱狀圖對比,以實測位置定位,參考相應(yīng)巖性進(jìn)行分析,結(jié)果表明礦物匹配結(jié)果是可靠的(表5)。
本文以湖北省大冶市雞冠咀銅金礦的實測光譜為基礎(chǔ),采用小波變換的方法提取各尺度能量特征,采用SAM方法進(jìn)行波譜匹配。從匹配結(jié)果來看,特征向量夾角最大者0.001432°,說明匹配效果較好。通過與巖性對比發(fā)現(xiàn),匹配結(jié)果是準(zhǔn)確的。
表5 光譜特征匹配結(jié)果
小波變換的特征表示指標(biāo)較多,可以用均值、絕對均值、模極大值等多種指標(biāo),匹配方法也可以選取距離、相關(guān)系數(shù)等,本文僅以能量作為特征向量、SAM作為匹配方法進(jìn)行了嘗試,后續(xù)可繼續(xù)進(jìn)行系統(tǒng)研究。
[1] 李興國,李躍華.毫米波近感技術(shù)基礎(chǔ)[M].北京:北京理工大學(xué)出版社,2009.
Li X G, Li Y H. Millimeter Wave Proximity Sensing Technology [M]. Beijing: Beijing Institute of Technology Press, 2009. (In Chinese)
[2] 王潤生,熊盛青,聶洪峰,等.遙感地質(zhì)勘查技術(shù)與應(yīng)用研究[J].地質(zhì)學(xué)報,2011,85(11):1699-1743.
Wang Y S, Xiong S Q, Nie H F,etal. Remote sensing technology and its application in geological exploration[J]. Acta Geologica Sinica, 2011, 85(11): 1699-1743. ( In Chinese)
[3] 蘇紅軍,杜培軍.高光譜數(shù)據(jù)特征選擇與特征提取研究[J].遙感技術(shù)與應(yīng)用,2006,21(4):288-293.
Su H J, Du P J. Study on feature selection and extraction of hyperspectral data[J]. Remote Sensing Technology and Application, 2006, 21(4): 288-293. ( In Chinese)
[4] Benedktsson J A, Sveinsson J R, Arnask K. Classification and feature extraction of AVIRIS data[J]. IEEE Trans on Geoscience and Remote Sensing, 1995, 33(5): 1194-1205.
[5] Nakariyakul S, Casasent D. Hyperspectral ratio feature selection: Agricultural product inspection example[J]. Proceedings of SPIE - The International Society for Optical Engineering, 2004, 5587: 133-143.
[6] Du P J, Fang T, Tang H,etal. Encoding methods of spectral vector in hyperspectral remote sensing image[J]. Journal of Shanghai University: English Edition, 2005, 9(1): 52-57.
[7] Li J. Liner Unmixing of Hyperspectral Signals Via Wavelet Feature Extraction[D]. Mississippi: Mississippi State University, 2002.
[8] Kuo B C, Landgrebe D. Improved statistics estimation and feature extraction for hyperspectral data classification[D]. West Lafayette, USA: Purdue University, 2001.
[9] 姚建紅.基于小波變換的地震勘探信號處理技術(shù)研究[D].大慶:大慶石油學(xué)院檔案館,2009.
Yao J H. Research on the Seismic Signal Process Based on Wavelet Transform[D]. Daqing: The Archive of Daqing Petroleum University, 2009. ( In Chinese)
[10] 李正周. Matlab數(shù)字信號處理與應(yīng)用[M].北京:清華大學(xué)出版社, 2008.
Li Z J. Matlab Digital Signal Processing and Application[M]. Beijing: Tsinghua University Press, 2008. ( In Chinese)
[11] 黃文鈺,尚海興.基于小波變換的高光譜遙感影像光譜匹配技術(shù)研究[J].西北水電,2013(1):7-9.
Huang W Y, Shang H X. Study on technology for spectrum matching of high spectrum remote sensing image based on wavelet change [J]. Northwest Hydropower, 2013(1): 7-9. ( In Chinese)
[12] 唐攀科.成像光譜相似礦物識別及其礦物填圖的不確定性研究[D]. 北京:中國地質(zhì)大學(xué)檔案館,2006.
Tang P K. Study on the Similar Mineral Identification and Uncertainty of Mineral Mapping of Imaging Spectrometry[D]. Beijing: The Archive of China University of Geosciences, 2006. ( In Chinese)
[13] 杜培軍,陳云浩,方濤,等.基于光譜特征的高光譜遙感影像檢索[J].光譜學(xué)與光譜分析,2005,25(8):1171-1175.
Du P J, Chen Y H, Fang T,etal. Spectral feature based on hyperspectral RS image retrieval[J]. Spectroscopy and Spectral Analysis, 2005, 25(8): 1171-1175. (In Chinese)
[14] Kruse F A, Lefkoff A B, Boardman J W,etal. The spectral image processing system (SIPS) software for integrated analysis of AVIRIS data [C]//Summaries of the 4th Annual JPL Airborne Geoscience Workshop. Pasadena: JPL Pub, 2014: 23-25.
Extraction of hyperspectral features from drilling core and its spectrum matching based on wavelet transform
LIU Bing-li1, 2, ZHANG Li2, HUANG Xin-lei2, GUO Ke2, LI Cheng2, ZHAO Li2, WANG Yu-lan2, CHEN Ling2
1.InstituteofGeophysicalandGeochemicalExploration,ChineseAcademyofGeoscience,Langfang065000,China; 2.KeyLaboratoryofGeomathematicsofSichuanProvince,ChengduUniversityofTechnology,Chengdu610059,China
Method for extraction of hyperspectral features from drilling cores and its spectral matching based on wavelet transform is investigated. Multi scale wavelet coefficients are used to build the features and spectral angle mapper (SAM) method is used to construct the matching method. 10 measured spectral in Jiguanzui copper gold deposit is analyzed. Experiment results indicate that the matching result is identical with the mineral in the borehole diagram. The proposed method shows a better performance on mineral identification by the combination of wavelet transforms and spectral angel mapper.
hyperspectral core; wavelet transform; feature extraction; spectrum matching
10.3969/j.issn.1671-9727.2016.05.14
1671-9727(2016)05-0630-05
2015-07-01。
中國地質(zhì)調(diào)查局項目(12120114002001);國家自然科學(xué)基金應(yīng)急管理項目(41541023)。
柳炳利(1981-),男,博士,講師,主要從事數(shù)學(xué)地質(zhì)研究, E-mail:liubingli-82@163.com。
P628.1
A