楊 潔，曾 強

(1.新疆大學(xué)資源與環(huán)境科學(xué)學(xué)院，新疆 烏魯木齊 830046；2.新疆大學(xué)干旱生態(tài)環(huán)境研究所，新疆 烏魯木齊 830046；3.綠洲生態(tài)教育部重點實驗室，新疆 烏魯木齊 830046)

地下煤火是指在自然環(huán)境或人為因素下，地下破碎煤體與空氣接觸后因氧化聚熱引發(fā)燃燒并不斷發(fā)展形成的大面積煤田火災(zāi)[1-5]。據(jù)不完全統(tǒng)計，我國新疆、內(nèi)蒙古、山西、寧夏、陜西等產(chǎn)煤大省(區(qū))的煤田火災(zāi)達200多處，造成的直接經(jīng)濟損失數(shù)以億計，同時燃燒對地下及近地表造成的破壞也很大，主要體現(xiàn)為煤火燃燒所引起的地表大面積的斷裂、塌陷，對土地、土壤資源及地表植被造成直接影響；煤火的熱效應(yīng)會使火區(qū)及其周圍的地表溫度上升，加速土壤涵養(yǎng)的水分蒸發(fā)，導(dǎo)致地表植被生存條件破壞甚至死亡[6]。遙感數(shù)據(jù)被廣泛應(yīng)用于煤火監(jiān)測和分析：JIANG等[7]、HUO等[8]基于Landsat TM/ETM影像，反演內(nèi)蒙古烏達火區(qū)地表溫度，分析了煤火蔓延方向；MARTHA等[9]、CHATTERJEE等[10]利用Landsat遙感數(shù)據(jù)計算出Jharia火區(qū)地表溫度，并通過反射紅外數(shù)據(jù)和煤火橫向蔓延識別出了火區(qū)位置；KUENZER等[11]利用高分辨率多光譜遙感影像和全色Quick bird數(shù)據(jù)反演烏達火區(qū)2005—2012年在滅火施工期間火區(qū)地表溫度變化以及土地覆蓋度變化；康高峰等[12]通過遙感技術(shù)對新疆的奇臺進行了地下煤火的監(jiān)測，研究和實驗都表明利用遙感和GIS手段可以初步圈定煤田火區(qū)的范圍；陳雅君等[13]、陳浩等[14]、夏安全等[15]基于Landsat 8影像分別反演得到了福州市、新沂市、濟南市的地表溫度并分別對研究區(qū)的土地類型進行了解譯；楊菊等[16]選取5期遙感影像數(shù)據(jù)，以像元二分模型估算提取植被覆蓋度，得到赤水丹霞植被覆蓋度時空變化規(guī)律。

目前，針對不同研究區(qū)域地表溫度或植被覆蓋因素的研究做了許多分析，算法方法相對成熟，但表征火區(qū)方面的研究還很少，而新疆地區(qū)煤火自燃災(zāi)害嚴重，掌握火區(qū)地表溫度和植被覆蓋變化情況及趨勢具有一定的現(xiàn)實意義。本文嘗試融合地表溫度和植被覆蓋情況表征火區(qū)基本特征，掌握研究區(qū)歷年監(jiān)測數(shù)據(jù)，并應(yīng)用GM(1,1)灰色預(yù)測模型進一步描述火區(qū)基本特征的變化趨勢，以期為運用遙感方法監(jiān)測火區(qū)，推算火區(qū)基本特征時空變化提供技術(shù)支撐。

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

水西溝位于新疆吉木薩爾縣西南20 km，東經(jīng)88°55′～88°58′之間，北緯43°55′～43°57′之間，火區(qū)面積74 975 m2，燃燒年損失儲量21.39萬t，隸屬于準南煤田(圖1)。導(dǎo)致水西溝煤田火災(zāi)的主要原因是該地區(qū)地質(zhì)活動較為劇烈，地層中的水平煤層經(jīng)過多次地質(zhì)運動，大多變?yōu)閮A斜煤層，煤層露頭受到雨水沖刷和風(fēng)蝕作用，并與空氣中氧氣發(fā)生氧化作用，積熱增溫，使煤火發(fā)生概率大大增加?；饏^(qū)燃燒煤層厚度分別為10.49 m和3.23 m，燃燒深度為40～60 m[3]。

圖1 研究區(qū)地理位置Fig.1 Geographical location of the study area

考慮季節(jié)、云和霧的影響，選取1998—2018年(2012年缺失)覆蓋研究區(qū)夏季、云量少且圖像清晰的TM影響和OLI/TIRS影像，作為反演溫度和植被蓋度的數(shù)據(jù)源。

2 方法基本原理

2.1 地表溫度反演

考慮到精度及參數(shù)以及Landsat數(shù)據(jù)熱紅外通道數(shù)目的限制，采用覃志豪等[17]改進的單窗算法，地表真實溫度可由式(1)計算求得。

LST=[a×(1-C-D)+(b×(1-C-D)+

C+D)×Tsensor-D×Ta]/C

(1)

式中：LST為最終反演得到的真實地表溫度，℃；a、b為擬合系數(shù)；C=ε×τ；D=(1-τ)×(1+(1-τ)×τ)；Ta為大氣平均作用溫度；ε為地表比輻射率；τ為大氣透過率；Tsensor為傳感器在衛(wèi)星高度上觀測熱輻射強度相對應(yīng)的地表溫度，包含有大氣和地表對熱輻射傳導(dǎo)的影響，需將其演算去除才能得到真實地表溫度，計算見式(2)。

Lλ=Gain×DN+Offset

(2)

式中：Lλ為輻射亮度值，W/(m2·μm·sr)；Gain為遙感圖像的增益；Offset為遙感圖像的偏移；DN為遙感圖像的灰度值。經(jīng)過輻射定標后得到輻射亮度值Lλ，由三部分組成：大氣向上輻射亮度L↑，地面的真實輻射亮度經(jīng)過大氣層之后到達衛(wèi)星傳感器的能量，大氣向下輻射到達地面后反射的能量L↓，表達式見式(3)。通過NASA官網(wǎng)查詢所需要參數(shù)，可通過式(4)求得溫度為T的黑體在熱紅外波段的輻射亮度B(TS)。

Lλ= [εB(TS)+(1-ε)L↓]τ+L↑

(3)

B(TS)=[Lλ-L↑-τ(1-ε)L↓]/τε

(4)

獲取溫度為T的黑體在熱紅外波段的輻射亮度后，根據(jù)Planck公式反函數(shù)求得地表真實溫度，見式(5)。

Ts=K2/ln(K1/B(Ts)+1)-273

(5)

式中，K1、K2為常數(shù)，可通過遙感影像頭文件查詢。

2.2 溫度異常區(qū)的提取

按照統(tǒng)一的人工閾值法[18]確定煤火提取區(qū)域的溫度閾值，低于分割點溫度閾值的稱為背景區(qū)，高于分割點溫度閾值的稱為煤火異常區(qū)。溫度閾值具體由式(6)～式(8)確定。

(6)

(7)

T閾=Tm+2Tδ

(8)

式中：TN為地表溫度圖像中任一像元溫度值；N為地表溫度圖像像元總數(shù)；Tm為地表溫度平均值；Tδ為地表溫度標準偏差；T閾為分割點閾值。

2.3 植被覆蓋度反演

基于與歸一化植被指數(shù)NDVI之間存在的較高相關(guān)性，采取像元二分模型[19]估算植被覆蓋度，通過式(9)和式(10)確定。

(9)

(10)

式中：NDVIsoil為完全是裸土或無植被覆蓋區(qū)域的NDVI值；NDVIveg為完全被植被所覆蓋的像元的NDVI值；NIR為近紅外波段的反射值；R為紅光波段的反射值。

2.4 GM(1,1)灰色預(yù)測模型

灰色預(yù)測模型(grey models)是根據(jù)過去及現(xiàn)在已知或非知的信息，建立從過去引申到未來的GM模型，以發(fā)現(xiàn)和掌握系統(tǒng)的發(fā)展規(guī)律，對系統(tǒng)的未來做出科學(xué)的定量預(yù)測。其基本原理是對原始數(shù)據(jù)累加淡化數(shù)據(jù)序列隨機性，提高數(shù)據(jù)序列的內(nèi)在規(guī)律，建立動態(tài)微分方程[16]。

設(shè)X(0)為非負序列，X(0)={X(0)(i)，i=1，2，…，n}為某一預(yù)測對象的非負單調(diào)原始數(shù)據(jù)列，n為序列長度，X(1)為X(0)的一次累加生成一次累加序列，見式(11)。

X(1)={X(1)(k),k=1,2,…,n}

(11)

式中，X(1)(k)由式(12)確定。

(12)

對X(1)可建立白化形式的微分方程，見式(13)。

(13)

式中，e、f為待識別參數(shù)。

(14)

(15)

(16)

(17)

3 結(jié)果與分析

3.1 地表溫度時空變化

將經(jīng)過預(yù)處理的遙感數(shù)據(jù)代入單窗算法反演模型，統(tǒng)計研究期間地表溫度信息，得到各期地表溫度的最小值、最大值、變化值、平均值、標準差和溫度閾值，結(jié)果見表1。

由表1可知，溫度隨時間持續(xù)發(fā)生動態(tài)變化，溫度最高值出現(xiàn)在2017年(53.17 ℃)，變化值為16.47 ℃；溫度最小值出現(xiàn)在2018年(7.08 ℃)，最大值僅為2 186 ℃，變化值為14.78 ℃。1998—2018年期間溫度最小值、最大值及其平均值分別減小18.54 ℃、22.28 ℃、3.74 ℃和19.73 ℃。研究期間溫度最小值、最大值、平均值及溫度閾值變化趨勢基本一致，如圖2所示。

由圖2可知，1998—2006年四個溫度信息均先增加后降低，峰值均出現(xiàn)在2003年區(qū)間附近；2007—2013年溫度信息呈現(xiàn)為增加后減小再增加；2014—2018年溫度信息呈現(xiàn)為減小后增大再減小?？傮w看來，火區(qū)地表溫度呈現(xiàn)整體上升趨勢，除2018年在原先上升趨勢下表現(xiàn)為驟減。

表1 地表溫度信息統(tǒng)計結(jié)果Table 1 Statistical results of surface temperature information 單位：℃

圖2 1998—2018年溫度最小值、最大值、平均值及閾值變化折線圖Fig.2 Line chart of minimum,maximum,average andthreshold changes from 1998 to 2018

3.2 溫度異常區(qū)變化

按照上述算法及流程，利用ArcGIS軟件對地表數(shù)據(jù)進行處理，結(jié)果見圖3。由圖3可知，研究區(qū)溫度異常范圍主要不連續(xù)分布在西中東三個區(qū)域(經(jīng)度等分為西中東三個區(qū)域)，隨著年份不斷增加，溫度異常區(qū)面積因煤火燃燒狀態(tài)呈波動變化。1998年火區(qū)三個異常區(qū)面積都不大，分布較為零散但均有存在，至2018年三個異常區(qū)面積不斷擴展，西部向東南方向擴展，中部向西北方向，并且不斷向東西方向延伸，東部面積卻明顯減少。這說明溫度異常區(qū)存在新生火區(qū)，危害程度增加，同時也有部分火區(qū)消亡。

3.3 植被覆蓋度

根據(jù)水西溝火區(qū)植被實地考察資料，研究區(qū)以稀疏草地、裸地等低植被覆蓋指數(shù)土地類型為主，因此將研究區(qū)植被覆蓋度等分為5個等級。其中，0%～20%為一級；20%～40%為二級；40%～60%為三級；60%～80%為四級；80%～100%為五級。一級為低植被覆蓋區(qū)，即重點研究區(qū)域，主要包括裸地、水域和戈壁區(qū)域。反演得到一級植被覆蓋面積和溫度異常提取面積見表2。由表2數(shù)據(jù)繪制研究區(qū)一級植被覆蓋面積及溫度異常區(qū)面積關(guān)系折線圖，見圖4。

由表2和圖4可知，1998—2018年研究區(qū)一級植被統(tǒng)計面積從6.89 km2增加至23.14 km2，溫度異常區(qū)統(tǒng)計面積從0.71 km2增加至0.83 km2，二者變化較為一致。 一級植被覆蓋面積最大值則出現(xiàn)在2017年，為24.15 km2，一級植被覆蓋面積的最小值則出現(xiàn)在2001年，為0.65 km2。 2009—2018年溫度異常區(qū)、一級植被面積較2000—2008年均有增加，截至2018年分別增加了0.12 km2和16.25 km2。

圖3 溫度異常區(qū)分布圖Fig.3 Distribution of temperature anomalies

表2 一級植被覆蓋面積及溫度異常區(qū)面積統(tǒng)計表Table 2 Statistics of first-level vegetation coverage area and temperature anomaly area

圖4 溫度異常區(qū)和一級植被覆蓋面積關(guān)系圖Fig.4 Relationship between temperature anomaly areaand first-level vegetation coverage area

3.4GM(1,1)模型預(yù)測分析

(18)

(19)

由式(18)和式(19)分別計算可得2018—2038年水西溝火區(qū)溫度異常區(qū)面積和一級植被覆蓋度面積，見表3。

表3 水西溝火區(qū)2018—2038年溫度異常區(qū)面積及一級植被覆蓋度面積預(yù)測表Table 3 Forecast of temperature anomaly area and first-level vegetation coverage areain Shuixigou fire area from 2018 to 2038

圖5 溫度異常區(qū)面積及一級植被覆蓋度面積預(yù)測趨勢圖Fig.5 Prediction trend of temperature anomaly areaand first-level vegetation coverage area

由表3可知，水西溝火區(qū)溫度異常區(qū)面積和一級植被覆蓋度面積預(yù)測值，如圖5所示。由圖5可知，水西溝火區(qū)2018—2038年期間，隨著溫度異常區(qū)面積的持續(xù)增加，一級植被面積覆蓋度面積減小后持續(xù)增加。

4 結(jié) 論

1) 1998—2018年，準南煤田水西溝火區(qū)溫度要素不斷變化，其中，溫度最小值、最大值、平均值及溫度閾值變化趨勢相似，相互之間影響較大。火區(qū)東西中三個區(qū)域有增有減，火區(qū)西部向東南方向擴展，中部向西北方向，并且不斷向東西方向延伸，東部面積有明顯減少。

2) 一級植被覆蓋度面積受溫度異常區(qū)面積的擴大減小影響較大，總體波動變化，呈現(xiàn)增長趨勢。2009—2018年期間溫度異常區(qū)和一級植被面積較前十年增長速度均有增加，截至2018年分別增加了0.12 km2和16.25 km2?；疑A(yù)測模型結(jié)果表明在預(yù)測期內(nèi)，隨著溫度異常區(qū)面積的增大，一級植被覆蓋度面積也增加。

3) 衛(wèi)星遙感影像反演地表溫度和植被覆蓋度可融合表征火區(qū)基本特征，動態(tài)監(jiān)測準南煤田火災(zāi)的變化情況，以期為相關(guān)部門治理火災(zāi)提供科學(xué)依據(jù)。