基于多光譜遙感的礦區(qū)積水區(qū)水深反演

吳 雙 武瑞杰 吳 侃 王 瑞，4 周大偉

（1.江蘇省資源環(huán)境信息工程重點實驗室，江蘇 徐州 221116；2.中國礦業(yè)大學(xué)環(huán)境與測繪學(xué)院，江蘇 徐州 221116；3.武漢大學(xué)資源與環(huán)境科學(xué)學(xué)院，湖北 武漢 430079；4.江西理工大學(xué)應(yīng)用科學(xué)學(xué)院，江西 贛州 341000）

淮南礦區(qū)經(jīng)過長期開采，地面發(fā)生大面積坍塌、沉陷，生態(tài)環(huán)境和地表景觀受到破壞，形成了眾多沉陷積水區(qū)，導(dǎo)致土地資源被浪費，嚴(yán)重影響了礦區(qū)周圍居民的正常生活［1-2］。了解礦區(qū)沉陷積水面積和積水深度可以為礦區(qū)生態(tài)環(huán)境評價和治理以及后期土地復(fù)墾等工作提供基礎(chǔ)數(shù)據(jù)支撐。遙感水深反演是水深測量的重要方法之一，相比于傳統(tǒng)的聲納測深和機載激光測深，具有探測范圍廣、快速靈活、重復(fù)觀測、不受地理位置和人為條件限制、費用低廉等優(yōu)點［3-6］。

自20世紀(jì)70年代起，國內(nèi)外眾多學(xué)者基于遙感手段在水深測量方面進行了大量研究，近年來，隨著遙感技術(shù)不斷發(fā)展，遙感水深測量的精度和實用性都得到了很大提高。陸天啟等［7］基于SPOT-6高分辨率遙感影像，對南海東鑼島海域利用綠-紅波段比值模型得出了具有借鑒意義的水深反演結(jié)果；張磊等［8］基于Landsat 8遙感影像，利用12個反演因子構(gòu)建了多波段反演模型，在內(nèi)蒙古渤海灣水庫得到了絕對誤差為0.68 m的水深反演結(jié)果；盛琳等［9］利用分辨率為2 m的Worldview-2遙感影像，基于逐步回歸的非線性模型對南海島環(huán)礁中東島水深進行了反演，精度得到大幅提高，但會受到波浪影響；李麗［10］利用WorldView-2衛(wèi)星數(shù)據(jù)對西沙群島水深進行了反演，均方根誤差在1.25 m以內(nèi)?？傮w上，水深反演研究主要集中在近海、湖泊，水庫等區(qū)域，礦區(qū)開采沉陷積水區(qū)水深反演研究成果較少，代表性成果有：鄧軍［11］對徐州九里塌陷湖建立了多元回歸模型，得到了平均誤差為0.32 m的反演的相關(guān)結(jié)果；武彥斌等［12］利用近紅外波段以2 m為閾值分區(qū)構(gòu)建了反演模型，在淮南積水沉陷區(qū)2 m深度以內(nèi)得到了平均相對誤差為13.367%、2～6 m深度平均相對誤差為15.02%的結(jié)果。

本研究以淮南潘一礦為例，基于2005年Landsat 5 TM多光譜遙感影像，選取擬合優(yōu)度高的反演因子，構(gòu)建精度較高的水深反演模型，并應(yīng)用在相同地區(qū)其他年份的遙感影像，通過與礦區(qū)沉陷預(yù)測預(yù)報系統(tǒng)（MSPS）軟件計算的水深結(jié)果進行比較分析，探討開采沉陷積水區(qū)水深反演的可行性，為礦區(qū)生態(tài)環(huán)境動態(tài)監(jiān)測及治理提供依據(jù)。

1 研究區(qū)概況

淮南潘謝礦區(qū)地處洪淮平原，內(nèi)有泥河流域，該地區(qū)雨量充沛，地下水補給量大，潛水位高［13］。潘一礦自20世紀(jì)90年代建成投產(chǎn)后，隨著開采工作不斷推進，造成了數(shù)量眾多、分布面積大的地表沉陷積水區(qū)，積水區(qū)水源主要來自于淺層水以及大氣降水的匯流。研究區(qū)域位于淮南潘一礦采后面積最大、水深最大的沉陷積水區(qū)，如圖1所示。

該沉陷積水區(qū)與海洋、湖泊以及水庫等不同：潮汐影響較?。恍纬赡晗尴鄬^短，深度相對較淺，由于主要的水源是大氣降水，所以受到季節(jié)和天氣的影響較大；塌陷積水區(qū)水生植被較少，因而遙感影像的光譜信息與水深值相關(guān)性較高。

2 數(shù)據(jù)獲取及處理

2.1 數(shù)據(jù)獲取

2005年水深實測值來自于2005年4月中旬淮南礦業(yè)集團對潘謝礦區(qū)塌陷積水區(qū)進行的現(xiàn)場調(diào)查數(shù)據(jù)，2010年水深預(yù)計值來自于礦區(qū)沉陷預(yù)測預(yù)報系統(tǒng)（MSPS）軟件結(jié)合淮南礦區(qū)地質(zhì)條件、采煤計劃等實際情況基于概率積分法得到的豐水期（每年8—9月）穩(wěn)定態(tài)預(yù)計數(shù)據(jù)。為了保證遙感影像與實際水深之間的同步性，遙感影像分別選取2005年4月13日與2010年8月1日采集的Landsat 5TM多光譜遙感影像，具體參數(shù)見表1。

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

（1）輻射定標(biāo)。本研究采用的兩期遙感影像為Leve 1T級別的遙感影像，已經(jīng)經(jīng)過系統(tǒng)的輻射校正和基于地面控制點的幾何校正，并且通過DEM進行了地形校正。在此基礎(chǔ)上，利用兩期Landsat 5遙感影像的拍攝日期和太陽高度角等參數(shù)，對每個波段都進行了輻射定標(biāo)，使兩期影像各波段的DN值轉(zhuǎn)換為具有可比性的絕對輻射亮度值［14］。

（2）水體提取。淮南礦區(qū)沉陷積水區(qū)與海洋、水庫等具有明顯水跡線的水域不同，為了保證水深反演不受其他類型地物的影響，需要在研究區(qū)內(nèi)進行水體提取。本研究采用歸一化差分水指數(shù)（Normal?ized Difference Water Index，NDWI）［15］，并基于 ENVI 5.3軟件進行波段運算，設(shè)置水體閾值，構(gòu)建水體掩膜進行水體提取。

（3）實測信息獲取。本研究采用的2005年水深實測值為區(qū)域水深等值線，在1，3，5 m 3條水深等值線上均勻、隨機地選取了45個測試點，并記錄每一點的經(jīng)緯度坐標(biāo)以及水深值。

（4）預(yù)計值獲取。利用礦區(qū)沉陷預(yù)測預(yù)報系統(tǒng)（MSPS）軟件，根據(jù)淮南礦區(qū)開采沉陷現(xiàn)狀并參考《淮南礦區(qū)開采沉陷參數(shù)及預(yù)計軟件研究》項目成果，選取預(yù)計參數(shù)，按照概率積分法的基本原理得到穩(wěn)定態(tài)預(yù)計下沉數(shù)據(jù)，結(jié)合2005年水深實測值以及2010年當(dāng)?shù)厮痪€等地質(zhì)資料，生成積水區(qū)預(yù)計水深等值線，利用ArcGIS軟件，構(gòu)建TIN數(shù)據(jù)，再轉(zhuǎn)成柵格形式的水下DEM數(shù)據(jù)。

3 水深反演模型建立

本研究主要依據(jù)水體在不同深度的光譜特性構(gòu)建水深反演模型。如圖2所示，電磁波在水體深水區(qū)由于被吸收了更多的能量，輻射出的能量相對較弱，從而表現(xiàn)在遙感影像上輻射亮度較小，淺水區(qū)輻射出的能量相對較強，所以遙感影像的輻射亮度值會較大，故可據(jù)此利用遙感影像反演水體深度［14-16］。此外，由于光譜傳感器接收的水體遙感光譜信號不僅包含了水體水深信息，還包含了水體懸浮物，水體葉綠素含量以及水底等多個層次的光譜信息，并且不同波段對水體深度信息的敏感度也不同，所以隨著水體深度變大，各波段輻射亮度值并非簡單呈現(xiàn)線性下降趨勢，而是有著不同幅度的波動，故而僅使用單個波段或者雙波段比值模型進行水深反演，可能會影響精度［17］。多波段模型綜合利用了傳感器多個波段的光譜信息，盡可能獲取更多的水深信息，可以在一定程度上提高反演模型的擬合優(yōu)度［18］。

3.1 確定反演因子

水深反演模型構(gòu)建需要確定良好的反演因子［19］。由圖2可知，在波長大于1 μm時，水體反射的相對能量極其微弱，于是選擇波長在1 μm以內(nèi)的B1藍色波段、B2綠色波段、B3紅色波段和B4近紅外波段進行反演因子研究以及模型構(gòu)建。

利用ENVI 5.3軟件及其二次開發(fā)IDL平臺編寫程序，批量提取45個水深實測值對應(yīng)地理坐標(biāo)下的遙感影像像元灰度值，并輸出為CSV文件，將水深實測值與各波段的像元輻射亮度值一一對應(yīng)。以水深實測值為因變量y，各波段的像元輻射亮度值為自變量x，分別建立了單波段線性、指數(shù)、對數(shù)模型以及雙波段指數(shù)模型，并進行了擬合優(yōu)度計算，如表2所示［6］。

由表2可知：單波段模型中藍色波段與近紅外波段擬優(yōu)度最佳，各類型模型擬優(yōu)度均達到0.8以上，綠色波段模型擬合優(yōu)度低于紅色波段，可能是由于沉陷積水區(qū)水質(zhì)特性（懸浮泥沙較多），水深光譜特性發(fā)生了“紅移”現(xiàn)象所致。而雙波段比值模型擬合優(yōu)度最佳的組合為藍色波段與近紅外波段的比值，擬合優(yōu)度為0.720 7，剩余比值模型擬合優(yōu)度欠佳，整體精度不如單波段比值模型。本研究選取反演效果較好的單波段因子進行建模，由于單波段各類型擬合優(yōu)度相差不大，為了方便計算，統(tǒng)一選取b1、b2、b3和b4線性因子構(gòu)建的多波段回歸模型公式為

式中，y為水深反演值，m；b為反演因子；p為反演因子對應(yīng)的擬合系數(shù)；q為常數(shù)項擬合系數(shù)；n為反演因子個數(shù)。

基于Visual Studio 2015平臺，編寫多元線性回歸算法，代入對應(yīng)的水深實測值和各波段輻射亮度值，經(jīng)過運算得到以下模型：

經(jīng)過取樣檢驗，部分水深反演值出現(xiàn)負(fù)值情況，并且通過與水深實測值進行對比分析，擬合優(yōu)度為0.931 9，平均絕對誤差為0.33 m，平均相對誤差為18.32%，造成誤差較大的原因主要是由于水深較淺區(qū)域出現(xiàn)了部分負(fù)值，故而對上述模型進行了改進，將初步水深反演值繼續(xù)迭代回歸，得到多元非線性正值模型（圖3），模型公式為

式中，y為最終水深反演值，m；x為初步水深反演值，m。

由圖3可知：經(jīng)過改進的多元非線性回歸模型使得水深均為正值，在一定程度上保證了數(shù)值的準(zhǔn)確性，擬合優(yōu)度為0.970 7，比未校正前的初步反演模型有所提高。

3.3 誤差分析

利用ENVI 5.3軟件波段運算模塊將水深反演模型應(yīng)用在2005年研究區(qū)遙感影像上，并在ArcGIS中生成水深等值面，結(jié)果如圖4所示。為了進一步驗證反演精度，利用ENVI 5.3軟件及其二次開發(fā)IDL平臺編寫程序，批量提取了水深柵格影像中對應(yīng)水深的實測值和反演值，并計算了絕對誤差與相對誤差。經(jīng)過統(tǒng)計計算，平均絕對誤差為0.21 m，平均相對誤差為10.92%，精度較高，表明反演結(jié)果可靠性強。

由圖4可知：反演結(jié)果良好，兩片表現(xiàn)為黑色的深水區(qū)域與實測水深圖基本吻合，水體中間也存在部分顯示為白色斑塊的異常區(qū)域，可能是受到植被或者漂浮物影響，以及反演模型本身存在一定程度的不穩(wěn)定性造成的。水深反演部分區(qū)域表現(xiàn)為不連續(xù)的現(xiàn)象，主要是因為Landsat 5 TM影像單個像素面積為30 m×30 m，單像素水深值為900 m2區(qū)域內(nèi)的平均水深，從而導(dǎo)致在水深變化區(qū)域會出現(xiàn)斑駁現(xiàn)象。

4 水深反演模型應(yīng)用分析

根據(jù)2005年研究區(qū)實測值建立的水深反演模型擬合優(yōu)度高，誤差為分米級，由于遙感數(shù)據(jù)的重復(fù)觀測特性以及可比性，本研究將該模型應(yīng)用在2010年同一地區(qū)的遙感影像中，并將反演值與同期水深預(yù)計值進行比較分析。

4.1 模型應(yīng)用

由于已經(jīng)對兩期遙感影像做過輻射定標(biāo)處理，像元DN值已轉(zhuǎn)換為絕對輻射亮度，使得不同時間獲取的同一區(qū)域的遙感影像具有一定的可比性。因此，通過2005年研究區(qū)遙感數(shù)據(jù)建立的水深反演模型在一定程度上可以應(yīng)用在2010年相同研究區(qū)域的遙感數(shù)據(jù)中。基于ENVI 5.3軟件波段運算模塊，將反演模型應(yīng)用在2010年Landsat 5遙感影像，可得到研究區(qū)2010年水深反演圖（圖5（a）），并與2010年水深預(yù)計圖進行對比分析，結(jié)果見圖5（b）。

由圖5可知：反演結(jié)果較好，2010年水深反演圖出現(xiàn)了與2005年近似的異常白色斑塊區(qū)域以及部分區(qū)域不連續(xù)的斑駁現(xiàn)象，但是異常區(qū)域明顯增多，除去植被或者漂浮物的影響以外，還有由于時間間隔過長，以及所處季節(jié)不同等原因引起的水質(zhì)變化等現(xiàn)象。

4.2 誤差分析

對2010年水深反演值和水深預(yù)計值進行變化監(jiān)測分析，發(fā)現(xiàn)2010年積水區(qū)域共計20 784個像元。表3為水深差值統(tǒng)計結(jié)果，經(jīng)過統(tǒng)計可得，水深差值小于1 m的像元占比最多，為47.30%，也存在小部分差值較大，相差大于6 m的區(qū)域，但僅占到0.97%，水深差值總體集中在3 m以內(nèi)，共計83.73%。為了更直觀地分析研究區(qū)水深差值分布，繪制了如圖6所示的水深差值分布圖。

由圖6可知：水深差值近似呈正態(tài)分布，無明顯的系統(tǒng)誤差，水深反演值與水深預(yù)計值變化規(guī)律一致，在一定程度上說明二者可以起到互相驗證的作用，兩種方法都可以得到在一定精度范圍內(nèi)的水深值。為了進一步驗證精度，選取測試點進行統(tǒng)計分析，2010年水深反演值與水深預(yù)計值的平均絕對誤差為0.63 m，平均相對誤差為15.24%。相比于2005年，誤差有所增大，這是由于模型本身的誤差與預(yù)計值的誤差疊加所致。研究表明：本研究構(gòu)建的模型，即利用多光譜遙感對其包含的水體深度信息進行多元非線性回歸得到的水深反演模型，可以為礦區(qū)沉陷積水區(qū)動態(tài)監(jiān)測以及生態(tài)環(huán)境評價提供參考。

5 結(jié) 論

結(jié)合2005年淮南潘一礦水深實測值，基于多光譜遙感影像建立了非線性多元回歸模型，將其應(yīng)用到2010年同一地區(qū)，并與基于概率積分法的水深預(yù)計值進行了比較分析，得出如下結(jié)論；

（1）通過多波段線性回歸模型建立了遙感數(shù)據(jù)的像元絕對輻射亮度值與水深值之間的定量關(guān)系，擬合優(yōu)度為97.07%，平均絕對誤差為0.21 m，平均相對誤差為10.92%，水深反演值與水深實測值的水深等值面符合度較高，模型可靠性強。

（2）將水深反演模型應(yīng)用于2010年同一地區(qū)的遙感數(shù)據(jù)，并與基于概率積分法的水深預(yù)計值進行了比較分析，平均絕對誤差為0.63 m，平均相對誤差為15.24%。由于該年度沒有實測值，所以無法判斷哪種方法的精度更高，但是二者水深等值面變化規(guī)律相似，水深差值呈正態(tài)分布，可以互相驗證和補充，從某種意義來說，兩種方法均可以為積水區(qū)水深測量提供數(shù)據(jù)支撐。在沒有充分掌握礦區(qū)資料的情況下可以側(cè)重于使用遙感反演方法；在礦區(qū)資料充分可靠的情況下，可以選擇礦區(qū)沉陷預(yù)測預(yù)報系統(tǒng)（MSPS）軟件進行水深預(yù)計，或者將二者結(jié)合起來互相補充得到精度較高的成果。

（3）本研究還存在一定的不足，后續(xù)工作中需要進一步研究的內(nèi)容有：①水深反演模型的精度在很大程度上取決于水深實測值的精度，試驗中的實測數(shù)據(jù)為水深等值線，只有間隔2 m的離散程度較大的實測值，后期可以通過補充與豐富實測值來優(yōu)化模型；②將水深反演模型應(yīng)用在其他年份或者其他地區(qū)時，可以添加調(diào)整因子，并利用少量實測值對模型進行改進；③通過獲取2010年水深實測值，對遙感方法和基于概率積分法的預(yù)計方法進行精度評判，進一步建立遙感結(jié)合概率積分法的高精度水深反演模型。