王啟元 肖 武,2 李素萃 張文凱 張 偉

(1.中國礦業(yè)大學(xué)(北京)土地復(fù)墾與生態(tài)重建研究所，北京市海淀區(qū)，100083；2.浙江大學(xué)公共管理學(xué)院，浙江省杭州市，310058)

摘 要利用時序TM遙感影像數(shù)據(jù)，在對比單波段閾值法、譜間關(guān)系法、歸一化水體指數(shù)法和改進歸一化水體指數(shù)法4種水體信息提取精度的基礎(chǔ)上，選擇改進歸一化水體指數(shù)(MNDWI)法并結(jié)合GIS空間疊加技術(shù)，對潘謝礦區(qū)2005-2017年采沉陷水域進行提取并分析其時空演變特征。結(jié)果顯示，2005-2017年，淮南潘謝礦區(qū)采煤沉陷水域面積增至4948.65 hm2，總增量達到4003.52 hm2，年均增長333.63 hm2，年均增長率為20.1%；分析研究區(qū)內(nèi)煤炭開采作為主要貢獻因素對沉陷水域面積變化影響及預(yù)測。通過線性擬合分析，預(yù)計到2020年潘謝礦區(qū)沉陷水域面積5968.12 hm2；2021-2030年沉陷水域面積年均增約為350 hm2，至2030其總量將達到9468.12 hm2。利用時序TM遙感影像對高潛水位礦區(qū)沉陷水體的時序動態(tài)變化監(jiān)測，分析了其變化原因，為后期采煤沉陷區(qū)的土地復(fù)墾與治理提供了科學(xué)依據(jù)。

東部高潛水位礦區(qū)煤炭開采引起的地表沉陷積水，嚴重影響了當(dāng)?shù)氐霓r(nóng)業(yè)發(fā)展和居民生活。國內(nèi)礦區(qū)生態(tài)環(huán)境方面的研究大多集中在生態(tài)環(huán)境修復(fù)、沉降監(jiān)測以及環(huán)境影響評估等方面，對于礦區(qū)采煤沉陷積水的研究主要集中在沉陷積水動態(tài)監(jiān)測以及范圍確定等方面。快速且精準獲取沉陷水域時空變化信息是高潛水位煤糧復(fù)合區(qū)進行環(huán)境整治與沉陷地復(fù)墾工作的前提。與傳統(tǒng)手段相比較，遙感監(jiān)測方法具有快速、精確、宏觀性強和低費用等特點，在沉陷水域的長時序動態(tài)監(jiān)測上具有無可比擬的優(yōu)勢。

國內(nèi)外常見的水體信息提取方法有單波段閾值法、水體指數(shù)法、多波段譜間關(guān)系法等。Frazier等針對水體在近紅外波段表現(xiàn)為強吸收性而植被和干土壤在該波段具有強反射性特點，提出單波段閾值法對水體進行提取。Kloiber等用TM影像的2波段加上3波段大于4波段加上5波段，提取水體信息。McFeeters提出歸一化水體指數(shù)(NDWI)，該指數(shù)提取水體時植被和土壤信息能夠被抑制。徐涵秋提出了改進歸一化水體指數(shù)(MNDWI)來提取水體信息。彭蘇萍等基于多時相TM影像，提取了淮南礦區(qū)積水沉陷面積的擴展變化信息；張瑞婭等將采區(qū)的原始地形歸入地表沉陷分析，確定了采煤沉陷的積水區(qū)域范圍；楊光華等利用高分遙感影像，對濟寧市采煤導(dǎo)致的塌陷積水耕地信息用人機交互解譯的方式進行了提取測算；孟磊等用TM數(shù)據(jù)分析淮南潘謝礦區(qū)1976-2006年的水體變化情況并分析了其演變的景觀生態(tài)效應(yīng)。利用遙感技術(shù)手段對水體信息進行提取，必須考慮到在光譜特征方面水體與其他地物的差異性，針對不同的研究區(qū)域的不同自然地理條件，在選擇水體提取方法時亦有所差異。本文以淮南潘謝礦區(qū)為研究區(qū)，比較4種遙感水體信息提取方法，結(jié)合礦區(qū)塌陷范圍實測資料及GIS空間疊加分析技術(shù)，實現(xiàn)研究區(qū)2005-2017年沉陷水域提取，揭示了潘謝礦區(qū)近十多年沉陷水域的時空演變特征極，并分析煤炭開采作為主要原因?qū)ζ浔O(jiān)測時段內(nèi)的影響及趨勢預(yù)測，對未來采煤沉陷水域變化趨勢預(yù)測以及生態(tài)環(huán)境治理具有現(xiàn)實意義。

1 研究區(qū)及數(shù)據(jù)源

1.1 研究區(qū)概況

潘謝礦區(qū)位于安徽省淮南市西北部，地形平坦開闊，地處淮河沖積平原，整體地勢西北高、東南低，地表高程一般為21.5～25.4 m，平均標高在23.1 m左右。所在區(qū)域?qū)偌撅L(fēng)溫暖帶半濕潤氣候，季節(jié)性表現(xiàn)為夏季炎熱，冬季寒冷。6～8月降雨量約占全年的40%，年平均降雨量865.5 mm，年平均蒸發(fā)量1610.44 mm，蒸發(fā)量大于降雨量，潮濕系數(shù)近似0.5。區(qū)內(nèi)溝渠縱橫，水利設(shè)施完善，當(dāng)?shù)貪撍宦裆罴s為1.5 m。因此地下煤炭開采沉陷后，地表極易形成積水，如圖1所示。截至2014年，淮南礦區(qū)累計原煤產(chǎn)量約6.7億t，根據(jù)生產(chǎn)規(guī)劃，到2020年，礦區(qū)累計原煤產(chǎn)量約11億t；2021-2030年，原煤產(chǎn)量穩(wěn)定在7000萬t/a。

圖1 淮南潘榭礦區(qū)遙感影像

1.2 數(shù)據(jù)選擇及預(yù)處理

為了研究不同水體提取方法在本區(qū)域內(nèi)的精度與適用性，選取研究區(qū)內(nèi)丁集礦、顧橋礦和顧北礦作為實驗區(qū)域?qū)Ρ确治?種水體提取方法中的最優(yōu)法，選用2016年9月2日Landsat8 OLI影像作為數(shù)據(jù)源，然后整體提取潘謝礦區(qū)6期沉陷水域信息，見表1。所有影像來源于地理空間數(shù)據(jù)云，空間分辨率為30 m，寬幅185 km×185 km，影像包含紅、綠、藍、近紅外和短波紅外等9個波段。由于研究所采用的遙感影像獲取時均已經(jīng)過系統(tǒng)的輻射校正與幾何校正等相關(guān)預(yù)處理工作，本文對數(shù)據(jù)的預(yù)處理工作主要包括輻射定標、大氣校正以及影像裁剪等。2005-2017年6期遙感影像日期均在5月份左右，未到研究區(qū)多降水季節(jié)，降水量對于水體提取結(jié)果的影響不予考慮。

表1 研究區(qū)遙感數(shù)據(jù)基本情況

2 研究方法

基于TM影像的水體信息自動提取方法，考慮影像光譜特性的提取方法目前應(yīng)用最多，其中單波段閾值法、多波段譜間關(guān)系法、水體指數(shù)法這3種方法最為廣泛。本文使用4種水體信息提取方法，選擇研究區(qū)內(nèi)丁集、顧橋、顧北3個礦作為實驗區(qū)，對實驗區(qū)內(nèi)所有水體進行提取，經(jīng)過精度評價找到適合研究區(qū)的最優(yōu)水體信息提取方法，再利用優(yōu)選后的方法對研究區(qū)內(nèi)2005年、2008年、2009年、2013年、2015年、2017年共6期水體進行分階段提取，分析煤炭開采對沉陷水域變化的影響以及趨勢預(yù)測。

2.1 TM影像水體信息提取方法

目前利用遙感影像進行水體信息提取的方法主要有單波段(近紅外)閾值法、多波段譜間關(guān)系法和水體指數(shù)法等。單一的水體信息提取方法無法通用，需要根據(jù)研究數(shù)據(jù)情況以及區(qū)域地理條件擇優(yōu)選擇。

2.1.1 單波段閾值法

單波段閾值法利用近紅外波段水體的強吸收性以及植被和干土壤的強反射性特點，選取適當(dāng)?shù)拈撝祵⑺w信息提取出來。此方法原理簡單，但水體與非水體之間存在的過渡區(qū)域容易被忽略，同時在區(qū)分細小水體和混淆在水體中的陰影上很困難，從而產(chǎn)生誤提、漏提的現(xiàn)象。利用水體信息和其他地物信息在遙感圖像不同波段上表現(xiàn)出來的不同亮度值來設(shè)定閾值，提取水體信息。依據(jù)實驗，水體在TM5(近紅外)波段具備很強的吸收作用，亮度值表現(xiàn)較低，而在該波段其他地物亮度值和水體有明顯區(qū)別。使用ENVI5.2軟件，對2016年9月2日Landsat8 OLI影像近紅外波段灰度圖像使用密度分割法確定閾值，設(shè)定水體提取閾值，利用波段運算工具，提取水體信息。

2.1.2 多波段譜間關(guān)系法

波譜間關(guān)系法可根據(jù)水體與地物波譜曲線的特征，找出其中的變化規(guī)律并設(shè)定邏輯判別式提取水體，這樣能夠?qū)⑺w與陰影較好地區(qū)分出來，而且提取精度較高。水體具備TM2(藍)加TM3(綠)大于TM4(紅)加TM5(近紅外)的特征，經(jīng)過試驗發(fā)現(xiàn)，設(shè)定一定的閾值來得到提取水體的效果更好，設(shè)定模型如下：

(1)

通過多次反復(fù)實驗，選取K1及K2值，使用波段運算工具，得到水體信息提取圖。

2.1.3 歸一化水體指數(shù)法

水體的反射率從可見光波段到短紅外波段逐級減弱，尤其在近紅外和短波熱紅外范圍內(nèi)顯示強烈的吸收性，因此可采用水體在可見光波段和近紅外波段(NIR)的光譜差異構(gòu)建水體指數(shù)，而植被在近紅外波段的反射率最強，采用綠波段(GREEN)和近紅外波段構(gòu)建指數(shù)可以有效地抑制植被和土壤信息，NDWI指數(shù)法如下：

(2)

同樣對于NDWI指數(shù)設(shè)定閾值上下限1，K2〗，可以更好地獲得水體信息提取效果。

2.1.4 改進歸一化水體指數(shù)法

土壤、建筑物因素的干擾在NDWI指數(shù)中被忽略，而水體綠光和近紅外的波段波譜特征與土壤、建筑物幾乎一致。徐秋涵等利用中紅外波段(MIR)，提出MNDWI指數(shù)用于水體信息的提取，計算公式如下：

(3)

2.2 遙感影像與采煤沉陷邊界空間疊加

采用水體信息提取方法中精度最高的改進歸一化水體指數(shù)法分別對研究區(qū)內(nèi)2005-2017年(2005年5月8日、2008年5月15日、2009年6月3日、2013年5月21日、2015年5月19日、2017年4月30日)共6期遙感影像分別進行體提取，在Arcgis軟件中，使用空間疊加分析模塊，將水體提取結(jié)果與淮南潘榭礦區(qū)2015年地表沉降范圍疊加分析，實現(xiàn)對潘榭礦區(qū)采煤沉陷水域的提取。

3 結(jié)果與分析

3.1 水體信息提取方法比較及精度分析

經(jīng)過前期預(yù)處理，對實驗區(qū)OLI遙感影像應(yīng)用各水體信息自動提取方法，得到提取結(jié)果，見表2。

表2 不同方法水體信息提取結(jié)果對比(黑色為水體)

表2從左到右依次為單波段閾值法、多波段譜間關(guān)系法、NDWI指數(shù)法和MNDWI水體指數(shù)的水體信息提取結(jié)果，黑色部分為水體信息。結(jié)合原始影像目視解譯分析上述提取結(jié)果，MNDWI指數(shù)法對水體信息的提取效果最好，對水體斑塊完整性及連續(xù)性提取上表現(xiàn)出色，并且誤提現(xiàn)象最少；單波段閾值法對水體的提取效果次之，雖然也存在誤提現(xiàn)象，但在水體信息連續(xù)性提取上表現(xiàn)出色，相比與其他兩種自動提取方法效果好；水體指數(shù)法效果和譜間關(guān)系法相對較差，水體信息邊界提取不連貫，存在較多誤提現(xiàn)象。

使用ENVI5.3軟件對實驗區(qū)遙感原始影像目視解譯，均勻選擇100組樣本，包含像元5562個。經(jīng)統(tǒng)計得到水體提取混淆矩陣和Kappa統(tǒng)計，見表3。由表3可知，4種水體提取方法精度均達到90%以上，但MNDWI法提取研究區(qū)水體信息的用戶精度、總體精度以及Kappa相關(guān)系數(shù)為各方法中最高。

表3 不同水體方法提取結(jié)果精度評價

3.2 歷年水體信息提取

基于MNDWI指數(shù)法，用2005年、2008年、2009年、2013年、2015年和2017年TM/OLI影像潘謝礦區(qū)整體水域進行提取，結(jié)果如圖2所示。

由圖2可以看到，2005-2017年，潘謝礦區(qū)地表全域水體擴張明顯，尤其在2010年之后增速加快，區(qū)內(nèi)水體斑塊總量不斷增加，個體積水斑塊多數(shù)由小變大，且斑塊形狀變化巨大。

圖2 監(jiān)測年份潘謝礦區(qū)地表水域分布

3.3 采煤沉陷水域監(jiān)測

由于本文主要側(cè)重于分析開采沉陷導(dǎo)致的水體變化，而礦區(qū)內(nèi)仍存在大面積的自然水體，為了區(qū)分未受擾動的天然水體與開采沉陷水體，實現(xiàn)對研究區(qū)內(nèi)采煤沉陷水域的提取監(jiān)測，將淮南潘謝礦區(qū)2015年采煤沉陷邊界與水域提取結(jié)果結(jié)合，利用Arcgis10.3軟件進行空間疊加運算，獲得2005-2017年共6期潘謝礦區(qū)沉陷區(qū)范圍內(nèi)的水域結(jié)果，如圖3所示，并在SPSS軟件中統(tǒng)計分析潘謝礦區(qū)監(jiān)測時段的沉陷水域面積變化情況，見表4。

分析結(jié)果表明，淮南潘謝礦區(qū)沉陷水域面積由2005年的945.13 hm2增加至2017年的4948.65 hm2，累積增長量為4003.52 hm2，年均增長333.63 hm2，年均增長率為20.1%。沉陷水體總量明顯呈上升趨勢。其中，2005-2009年沉陷水域面積增加1459.3 hm2，增幅達154%，變化幅度相對最大；2009-2013年沉陷水域面積增加603.32 hm2，增幅為25.1%；2013-2017年沉陷水域面積增量為1940.9 hm2，增幅為64.5%，雖然增幅未至最大，但其沉陷水域面積的絕對增量在監(jiān)測時段內(nèi)為最大值。

表4 2005-2017年沉陷水域面積及變化值

綜合圖3和表4可知，2005年淮南潘謝礦區(qū)由于已采礦井?dāng)?shù)量較少，還未出現(xiàn)大范圍的采煤沉陷積水，只有潘一礦、潘二礦、潘三礦、顧北礦、張集礦出現(xiàn)小范圍的沉陷積水；2005-2009年，潘一礦、潘二礦、潘三礦、顧北礦、張集礦沉陷水域面積增加，其中謝集礦、張集礦擴張最為明顯；2009-2013年期間，除朱集礦外，研究區(qū)內(nèi)各礦沉陷水域面積均有所變化，顧橋礦水域增加明顯，2013年其內(nèi)部產(chǎn)生一大面積沉陷積水區(qū)；2013-2017年，研究區(qū)境內(nèi)所有沉陷水域面積均擴有增加，且分布均勻。

圖3 2005-2017年潘謝礦區(qū)沉陷水域分布

圖4 累積原煤產(chǎn)量與沉陷水域面積相關(guān)性分析

淮南潘謝礦區(qū)2005-2017年累積原煤產(chǎn)量與沉陷水域面積關(guān)系如圖4所示，二者相關(guān)系數(shù)r=0.94。隨著累積煤炭產(chǎn)量的增加以及礦區(qū)沉陷水域面積不斷擴大，根據(jù)淮南礦業(yè)集團生產(chǎn)規(guī)劃，到2020年潘謝礦區(qū)累積原煤產(chǎn)量預(yù)計達10.85億t；2020-2030年，原煤產(chǎn)量穩(wěn)定在7000萬t/a，累計原煤產(chǎn)量達17.85億t。根據(jù)上述擬合結(jié)果預(yù)計到2020年潘謝礦區(qū)沉陷水域面積總量預(yù)計增至5968.12 hm2；2021-2030年期間，區(qū)內(nèi)沉陷水域面積年均增量預(yù)計為350 hm2，到2030年沉陷水域面積將達到9468.12 hm2。

4 結(jié)論

本研究在對比分析各種水體信息方法的基礎(chǔ)上，以淮南潘榭礦區(qū)為例，對研究區(qū)內(nèi)所有水體進行了提取，結(jié)合潘謝礦區(qū)2015年采煤沉陷邊界，利用GIS空間疊加，實現(xiàn)對研究區(qū)2005-2017年地表沉陷水體的動態(tài)監(jiān)測，分析了區(qū)內(nèi)煤炭開采對沉陷水域面積變化的影響并做出預(yù)測，得出以下結(jié)論：

(1)依據(jù)研究所用數(shù)據(jù)特點以及研究區(qū)地理條件，對比了單波段閾值法、譜間關(guān)系法、歸一化水體指數(shù)法及改進歸一化水體指數(shù)法4種水體信息提取方法的優(yōu)劣，選擇提取效果最好、精度最高的改進歸一化水體指數(shù)法作為潘謝礦區(qū)水體信息提取方法，其總體精度達96.03%，Kappa系數(shù)為0.88。

(2)2005-2017年潘謝礦區(qū)沉陷水域的監(jiān)測結(jié)果顯示，2005-2017年，區(qū)內(nèi)沉陷水域面積增至4948.65 hm2，總增量為4003.52 hm2，年均增長333.63 hm2，年均增長率為20.1%。其中，2005-2009年沉陷水域面積增加1459.3 hm2，增幅達154%，變化幅度相對最大，2013-2017年沉陷水域面積增量為1940.9 hm2，沉陷水域面積的絕對增量在監(jiān)測時段內(nèi)為最大值。

(3)通過對監(jiān)測時段內(nèi)的累積原煤產(chǎn)量與沉陷水域面積的相關(guān)性定量分析可知，隨著煤炭開采量的累積增加，潘謝礦區(qū)沉陷水域面積隨之不斷上升，二者具有很強的相關(guān)性，相關(guān)系數(shù)為0.94。依據(jù)生產(chǎn)規(guī)劃和擬合結(jié)果，預(yù)計到2020年，潘謝礦區(qū)沉陷水域面積總量將增至5968.12 hm2，2021-2030年，其沉陷水域面積年均增量預(yù)計為350 hm2，到2030年沉陷水域面積將達到9468.12 hm2。