准南煤田水西沟火区地表温度及植被覆盖度时空变化特征

杨 洁，曾 强

(1.新疆大学资源与环境科学学院，新疆 乌鲁木齐 830046；2.新疆大学干旱生态环境研究所，新疆 乌鲁木齐 830046；3.绿洲生态教育部重点实验室，新疆 乌鲁木齐 830046)

地下煤火是指在自然环境或人为因素下，地下破碎煤体与空气接触后因氧化聚热引发燃烧并不断发展形成的大面积煤田火灾[1-5]。据不完全统计，我国新疆、内蒙古、山西、宁夏、陕西等产煤大省(区)的煤田火灾达200多处，造成的直接经济损失数以亿计，同时燃烧对地下及近地表造成的破坏也很大，主要体现为煤火燃烧所引起的地表大面积的断裂、塌陷，对土地、土壤资源及地表植被造成直接影响；煤火的热效应会使火区及其周围的地表温度上升，加速土壤涵养的水分蒸发，导致地表植被生存条件破坏甚至死亡[6]。遥感数据被广泛应用于煤火监测和分析：JIANG等[7]、HUO等[8]基于Landsat TM/ETM影像，反演内蒙古乌达火区地表温度，分析了煤火蔓延方向；MARTHA等[9]、CHATTERJEE等[10]利用Landsat遥感数据计算出Jharia火区地表温度，并通过反射红外数据和煤火横向蔓延识别出了火区位置；KUENZER等[11]利用高分辨率多光谱遥感影像和全色Quick bird数据反演乌达火区2005—2012年在灭火施工期间火区地表温度变化以及土地覆盖度变化；康高峰等[12]通过遥感技术对新疆的奇台进行了地下煤火的监测，研究和实验都表明利用遥感和GIS手段可以初步圈定煤田火区的范围；陈雅君等[13]、陈浩等[14]、夏安全等[15]基于Landsat 8影像分别反演得到了福州市、新沂市、济南市的地表温度并分别对研究区的土地类型进行了解译；杨菊等[16]选取5期遥感影像数据，以像元二分模型估算提取植被覆盖度，得到赤水丹霞植被覆盖度时空变化规律。

目前，针对不同研究区域地表温度或植被覆盖因素的研究做了许多分析，算法方法相对成熟，但表征火区方面的研究还很少，而新疆地区煤火自燃灾害严重，掌握火区地表温度和植被覆盖变化情况及趋势具有一定的现实意义。本文尝试融合地表温度和植被覆盖情况表征火区基本特征，掌握研究区历年监测数据，并应用GM(1,1)灰色预测模型进一步描述火区基本特征的变化趋势，以期为运用遥感方法监测火区，推算火区基本特征时空变化提供技术支撑。

1 研究区概况及数据来源

水西沟位于新疆吉木萨尔县西南20 km，东经88°55′～88°58′之间，北纬43°55′～43°57′之间，火区面积74 975 m2，燃烧年损失储量21.39万t，隶属于准南煤田(图1)。导致水西沟煤田火灾的主要原因是该地区地质活动较为剧烈，地层中的水平煤层经过多次地质运动，大多变为倾斜煤层，煤层露头受到雨水冲刷和风蚀作用，并与空气中氧气发生氧化作用，积热增温，使煤火发生概率大大增加。火区燃烧煤层厚度分别为10.49 m和3.23 m，燃烧深度为40～60 m[3]。

图1 研究区地理位置Fig.1 Geographical location of the study area

考虑季节、云和雾的影响，选取1998—2018年(2012年缺失)覆盖研究区夏季、云量少且图像清晰的TM影响和OLI/TIRS影像，作为反演温度和植被盖度的数据源。

2 方法基本原理

2.1 地表温度反演

考虑到精度及参数以及Landsat数据热红外通道数目的限制，采用覃志豪等[17]改进的单窗算法，地表真实温度可由式(1)计算求得。

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

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

(1)

式中：LST为最终反演得到的真实地表温度，℃；a、b为拟合系数；C=ε×τ；D=(1-τ)×(1+(1-τ)×τ)；Ta为大气平均作用温度；ε为地表比辐射率；τ为大气透过率；Tsensor为传感器在卫星高度上观测热辐射强度相对应的地表温度，包含有大气和地表对热辐射传导的影响，需将其演算去除才能得到真实地表温度，计算见式(2)。

Lλ=Gain×DN+Offset

(2)

式中：Lλ为辐射亮度值，W/(m2·μm·sr)；Gain为遥感图像的增益；Offset为遥感图像的偏移；DN为遥感图像的灰度值。经过辐射定标后得到辐射亮度值Lλ，由三部分组成：大气向上辐射亮度L↑，地面的真实辐射亮度经过大气层之后到达卫星传感器的能量，大气向下辐射到达地面后反射的能量L↓，表达式见式(3)。通过NASA官网查询所需要参数，可通过式(4)求得温度为T的黑体在热红外波段的辐射亮度B(TS)。

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

(3)

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

(4)

获取温度为T的黑体在热红外波段的辐射亮度后，根据Planck公式反函数求得地表真实温度，见式(5)。

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

(5)

式中，K1、K2为常数，可通过遥感影像头文件查询。

2.2 温度异常区的提取

按照统一的人工阈值法[18]确定煤火提取区域的温度阈值，低于分割点温度阈值的称为背景区，高于分割点温度阈值的称为煤火异常区。温度阈值具体由式(6)～式(8)确定。

(6)

(7)

T阈=Tm+2Tδ

(8)

式中：TN为地表温度图像中任一像元温度值；N为地表温度图像像元总数；Tm为地表温度平均值；Tδ为地表温度标准偏差；T阈为分割点阈值。

2.3 植被覆盖度反演

基于与归一化植被指数NDVI之间存在的较高相关性，采取像元二分模型[19]估算植被覆盖度，通过式(9)和式(10)确定。

(9)

(10)

式中：NDVIsoil为完全是裸土或无植被覆盖区域的NDVI值；NDVIveg为完全被植被所覆盖的像元的NDVI值；NIR为近红外波段的反射值；R为红光波段的反射值。

2.4 GM(1,1)灰色预测模型

灰色预测模型(grey models)是根据过去及现在已知或非知的信息，建立从过去引申到未来的GM模型，以发现和掌握系统的发展规律，对系统的未来做出科学的定量预测。其基本原理是对原始数据累加淡化数据序列随机性，提高数据序列的内在规律，建立动态微分方程[16]。

设X(0)为非负序列，X(0)={X(0)(i)，i=1，2，…，n}为某一预测对象的非负单调原始数据列，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为待识别参数。

(14)

(15)

(16)

(17)

3 结果与分析

3.1 地表温度时空变化

将经过预处理的遥感数据代入单窗算法反演模型，统计研究期间地表温度信息，得到各期地表温度的最小值、最大值、变化值、平均值、标准差和温度阈值，结果见表1。

由表1可知，温度随时间持续发生动态变化，温度最高值出现在2017年(53.17 ℃)，变化值为16.47 ℃；温度最小值出现在2018年(7.08 ℃)，最大值仅为2 186 ℃，变化值为14.78 ℃。1998—2018年期间温度最小值、最大值及其平均值分别减小18.54 ℃、22.28 ℃、3.74 ℃和19.73 ℃。研究期间温度最小值、最大值、平均值及温度阈值变化趋势基本一致，如图2所示。

由图2可知，1998—2006年四个温度信息均先增加后降低，峰值均出现在2003年区间附近；2007—2013年温度信息呈现为增加后减小再增加；2014—2018年温度信息呈现为减小后增大再减小。总体看来，火区地表温度呈现整体上升趋势，除2018年在原先上升趋势下表现为骤减。

表1 地表温度信息统计结果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 温度异常区变化

按照上述算法及流程，利用ArcGIS软件对地表数据进行处理，结果见图3。由图3可知，研究区温度异常范围主要不连续分布在西中东三个区域(经度等分为西中东三个区域)，随着年份不断增加，温度异常区面积因煤火燃烧状态呈波动变化。1998年火区三个异常区面积都不大，分布较为零散但均有存在，至2018年三个异常区面积不断扩展，西部向东南方向扩展，中部向西北方向，并且不断向东西方向延伸，东部面积却明显减少。这说明温度异常区存在新生火区，危害程度增加，同时也有部分火区消亡。

3.3 植被覆盖度

根据水西沟火区植被实地考察资料，研究区以稀疏草地、裸地等低植被覆盖指数土地类型为主，因此将研究区植被覆盖度等分为5个等级。其中，0%～20%为一级；20%～40%为二级；40%～60%为三级；60%～80%为四级；80%～100%为五级。一级为低植被覆盖区，即重点研究区域，主要包括裸地、水域和戈壁区域。反演得到一级植被覆盖面积和温度异常提取面积见表2。由表2数据绘制研究区一级植被覆盖面积及温度异常区面积关系折线图，见图4。

由表2和图4可知，1998—2018年研究区一级植被统计面积从6.89 km2增加至23.14 km2，温度异常区统计面积从0.71 km2增加至0.83 km2，二者变化较为一致。 一级植被覆盖面积最大值则出现在2017年，为24.15 km2，一级植被覆盖面积的最小值则出现在2001年，为0.65 km2。 2009—2018年温度异常区、一级植被面积较2000—2008年均有增加，截至2018年分别增加了0.12 km2和16.25 km2。

图3 温度异常区分布图Fig.3 Distribution of temperature anomalies

表2 一级植被覆盖面积及温度异常区面积统计表Table 2 Statistics of first-level vegetation coverage area and temperature anomaly area

图4 温度异常区和一级植被覆盖面积关系图Fig.4 Relationship between temperature anomaly areaand first-level vegetation coverage area

3.4GM(1,1)模型预测分析

(18)

(19)

由式(18)和式(19)分别计算可得2018—2038年水西沟火区温度异常区面积和一级植被覆盖度面积，见表3。

表3 水西沟火区2018—2038年温度异常区面积及一级植被覆盖度面积预测表Table 3 Forecast of temperature anomaly area and first-level vegetation coverage areain Shuixigou fire area from 2018 to 2038

图5 温度异常区面积及一级植被覆盖度面积预测趋势图Fig.5 Prediction trend of temperature anomaly areaand first-level vegetation coverage area

由表3可知，水西沟火区温度异常区面积和一级植被覆盖度面积预测值，如图5所示。由图5可知，水西沟火区2018—2038年期间，随着温度异常区面积的持续增加，一级植被面积覆盖度面积减小后持续增加。

4 结 论

1) 1998—2018年，准南煤田水西沟火区温度要素不断变化，其中，温度最小值、最大值、平均值及温度阈值变化趋势相似，相互之间影响较大。火区东西中三个区域有增有减，火区西部向东南方向扩展，中部向西北方向，并且不断向东西方向延伸，东部面积有明显减少。

2) 一级植被覆盖度面积受温度异常区面积的扩大减小影响较大，总体波动变化，呈现增长趋势。2009—2018年期间温度异常区和一级植被面积较前十年增长速度均有增加，截至2018年分别增加了0.12 km2和16.25 km2。灰色预测模型结果表明在预测期内，随着温度异常区面积的增大，一级植被覆盖度面积也增加。

3) 卫星遥感影像反演地表温度和植被覆盖度可融合表征火区基本特征，动态监测准南煤田火灾的变化情况，以期为相关部门治理火灾提供科学依据。