单天婵，郑伟，陈洁

国家卫星气象中心,北京 100081

1 引言

过火区是指经火灾烧毁后尚未长成新林的土地（武晋雯 等，2020），是火灾监测的一项重要指标。快速、精准地获取过火区的位置，范围等信息有助于评估灾后损失，对灾后恢复措施制定和生态环境变化监测有重要作用（李静 等，2018）。遥感卫星具有覆盖范围广、成像周期短的特点，适合过火区的监测和提取（郑伟 等，2011）。中国对森林火灾的监管严格，通常救援工作能够在火灾发生后及时有效地开展（鹿德林，2019），近几年单一过火区范围普遍小于100 km2。与低空间分辨率遥感数据相比，高空间分辨率数据提取范围较小的过火区具有明显优势。目前已有一些研究用到我国高分数据提取过火区，其中高分一号（GF-1）宽幅相机WFV（Wide Field of View）数据具有空间分辨率高，幅宽大等特点，在该研究领域中使用广泛。祖笑锋等（2015）利用GF-1 WFV各波段反射率信息，结合5种光谱指数，构建森林过火区识别决策树模型提取过火区。孙桂芬等（2019）基于GF-1 WFV 数据采用多种光谱指数对内蒙古鄂伦春自治旗过火区进行了识别，取得了良好的效果。

上述方法主要是针对研究区采集过火像元样本，分析过火像元在影像中的光谱特性，确定分类阈值，对过火区进行提取。由于高空间分辨率数据普遍缺少对过火区敏感的短波红外等波段，因此若仅基于该数据构成的光谱指数提取过火区，可能会出现阴影、裸地或其他非过火像元与过火像元混淆的情况。同时由于过火区所在地理位置、遥感影像采集的时间和卫星传感器性能差异等原因，上述方法中提出的分类阈值往往只适用于某一特定区域，普适性较差。

为了解决以上问题，可结合火点信息提取过火区。火点像元在一定时间、空间上的累积可用于判断过火区所在位置和范围。目前越来越多对过火区提取的研究用到了火点产品。主要包括2种方法：一种是先以火点产品结合光谱指数，选择明确过火的像元作为种子点，然后按一定规则迭代生长种子点形成过火区（谭明艳 等，2007；Giglio 等，2009；Yang 等，2013；Alonso-Canas 和Chuvieco，2015；Giglio 等，2018）。另一种是首先根据光谱指数设定阈值，确定可能的过火像元，在此基础上设置缓冲区，去除缓冲区外的像元形成过火区（Merino-de-Miguel 等，2010；Boschetti等，2015）。

上述方法用到的火点信息通常来自中分辨率成像光谱仪MODIS（Moderate-resolution Imaging Spectroradiometer）火点产品，这类产品主要基于中红外通道的亮温识别火点（张婕 等，2016）。自我国极轨气象卫星发射以来，可见光红外扫描辐射计VIRR（Visible and Infra-Red Radiometer）、中分辨率光谱成像仪MERSI（Medium Resolution Spectral Imager）等传感器均具有可判识火点的中红外等通道，对火点监测的研究已取得一定成果（赵立成，2004；王钊，2011；董晓锐和魏迎，2014；杨蕾等，2019）。目前，风云三号D星（FY-3D）MERSI火点监测产品已具有较高的精度，与MODIS 火点产品在时间分辨率、空间分辨率上相近，特别是在中国区域，风云三号D星（FY-3D）MERSI 火点监测产品具有更高的精度（郑伟 等，2020）。因此，可利用FY-3D MERSI 火点监测产品提取过火区。

综上所述，本文综合利用火点像元和过火像元在空间、时间、光谱特征上的关系，提出一种基于GF-1 WFV 数据和FY-3D MERSI 火点产品自动提取过火区的方法，以人机交互方式获得的过火区参考真值作验证，并与神经网络分类法提取过火区的结果作对比，快速准确地提取小范围过火区，进一步提升中国卫星遥感数据火情监测的能力。

2 研究区与数据

2.1 研究区

本文选择四川省凉山彝族自治州木里藏族自治县和西昌市2 处过火区作为研究区（图1）。每年3、4 月是上述地区的火灾高发期，2019 年3 月30 日木里县雅砻江镇立尔村发生的森林火灾，造成严重的人员伤亡（饶月明 等，2020）。2020 年3月底，木里、西昌又相继发生森林火灾。木里藏族自治县位于凉山彝族自治州西北部，处于横断山脉中段，地形复杂，海拔差异大，属于典型的高山峡谷地区（饶月明 等，2020）。西昌过火区所在的西昌市区东临邛海盆地、西临安宁河谷平原，地势北高南低，以山地为主，山区沟谷深切，地形崎岖（胡卸文 等，2020）。本文以2020年3月底的2 处火灾形成的过火区作为研究区，其中木里研究区范围为27°45′N—28°9′N，101°12′E—101°36′E；西昌研究区范围为27°39′N—27°57′N，102°24′E—101°36′E。

图1 研究区行政分布图Fig.1 Administrative distribution map of the study area

2.2 数据

2.2.1 GF-1 WFV影像

本文选择GF-1 WFV 数据结合FY-3D MERSI每日火点产品提取过火区。GF-1 WFV数据空间分辨率为16 m，在卫星不侧摆的情况下重访周期可达4 d，幅宽为800 km，包括蓝光、绿光、红光和近红外4个波段，各波段参数如表1所示。

表1 GF-1 WFV波段参数Table 1 The band parameters of GF-1 WFV

从中国资源卫星应用中心-陆地观测卫星数据服务平台（http：//36.112.130.153：7777/DSSPlatform/index.html/［2020-08-13］）查找并获取GF-1 WFV数据。为了降低过火区提取的误差，选择数据时除了要考虑火灾发生的时间和地点，还应尽量选用云量较少的影像。根据以上要求，本文选用的过火前后的2 景GF-1 WFV 影像，获取日期分别是2020年2月21日和2020年5月8日。

2.2.2 FY-3D MERSI每日火点产品

本文选择FY-3D MERSI 每日火点产品获取火点信息。该产品主要包括每日全球范围内火点的经纬度、火点强度和可信度等信息。其中可信度有4种，分别为火点、云区火点、可能是火点和可能是噪声。其火点判识算法主要根据FY-3D MERSI 空间分辨率为1000 m 的中红外通道对高温热源的敏感特性，建立合适的阈值探测火点像元（郑伟 等，2020）。根据过火前后的GF-1 WFV 影像获取时间，获取2020 年2 月21 日—5 月8 日研究区内的FY-3D MERSI每日火点产品。

2.2.3 地表类型信息

选择中国科学院地理科学与资源研究所自然环境科学与数据中心提供的土地覆被遥感监测数据集LUCC（Land Use/Land Cover Change）作为地表类型信息，空间分辨率为30 m。该数据集以Landsat 遥感数据作为主信息源，通过人工目视解译将土地利用分为了耕地、林地、草地等6个一级地类和25 个二级地类（徐新良 等，2018）。2 个研究区的一级地类分布情况如图2所示。

3 研究方法

本方法提取过火区的过程包括2个部分。第一部分：过火区粗提取。数据预处理后，基于高分数据读取并叠加火点，结合火点的时间、空间、光谱特征筛选火点像元并扩充，确定过火区粗略范围。第二部分：过火区精提取。根据每种地表类型所含火点像元的数量，确定过火区含有的地表类型；采用阈值分割算法，在过火区粗略范围内对每种地表类型确定分类阈值，分类过火像元和非过火像元，最后剔除小斑块，获得过火区提取结果。该方法的技术流程如图3所示。

图3 技术流程图Fig.3 Technique flow chart

3.1 计算光谱指数

根据GF-1 WFV 数据包含的波段，对过火前后2 景影像分别计算归一化差值植被指数NDVI（Normalized Difference Vegetation Index）和 由McFeeters（1996）提出的归一化差异水体指数NDWI（Normalized Difference Water Index），计算公式如下：

式中，ρNIR、ρR和ρG分别表示近红外、红光和绿光波段的反射率。

火灾发生后，植被结构遭到破坏，NDVI 值会明显下降，因此选用NDVI 差值dNDVI来突出植被的变化，dNDVI计算公式如下：

式中，NDVIpre表示过火前影像的NDVI值，NDVIpost表示过火后影像的NDVI值。

3.2 过火区粗略提取

3.2.1 去除云、水体和其他非植被像元

本文仅考虑提取在森林、草地等植被地区的过火区，因此可根据地表类型信息和光谱信息去除云、水体、城市等非植被区域。

在GF-1 WFV 影像选择和预处理阶段，尽可能选择无云或少云的影像。在有少量的云时，需要先去除过火前后2幅影像的云像元。基于已有研究（李俊杰和傅俏燕，2020），本文利用云像元在红、绿、蓝3个波段的反射率特性，粗略提取并剔除影像上的云像元。

水体反射率从可见光到近红外波段依次降低，在近红外波段反射率几乎为零，而植被在近红外波段反射率较高（童李霞 等，2017）。McFeeters（1996）利用上述绿光波段与近红外波段的关系，提出NDWI。该指数可以在一定程度上抑制植被信息，突出水体信息，在水体提取中得到广泛应用。本文方法利用过火前后2 景影像的NDWI 值提取影像中的水体信息，同时满足式（4）和式（5）的像元视为水体像元，需要去除：

式中，NDWIpre表示过火前影像的NDWI 值，NDWIpost表示过火后影像的NDWI值。

在此基础上，利用地表类型信息提取植被区域。为了便于计算，将LUCC 数据重采样为16 m。根据LUCC 的分类级别和编号，仅保留林地和草地部分的像元，具体编号和名称如表2所示。

3.2.2 火点处理

3.2.2.1 火点块叠加

根据过火前后2幅影像的获取时间，读取并叠加研究区范围内2020 年2 月21 日—5 月8 日的FY-3D MERSI 每日火点产品。这里只选择可信度为火点和云区火点的火点信息。以空间分辨率为16 m的WFV 影像图为基础，根据产品所提供的火点经纬度，可以确定火点所处的像元位置，作为中心火点像元，并记录其探测时间t。t用一年中的第几天（Day of Year，DOY）表示。若一个像元多时次被监测为中心火点像元，则取最晚的探测时间作为t。

由于FY-3D MERSI 每日火点产品是以空间分别率为1 km 的数据制作生成，为了降低空间分辨率从1 km 到16 m 可能带来的影响，尽可能保持火点在空间上的连续性，减轻火点识别遗漏等问题，将每个中心火点像元进行膨胀，记为一个火点块。根据1 km 和16 m 的倍数关系，将膨胀核设为63×63。火点块内的所有像元均为火点像元，探测时间与其中心火点像元一致，均为t。若一个火点像元在空间上同时属于多个相邻火点块，则把它归于探测时间最晚的火点块。以西昌研究区为例，叠加后的火点块如图4（a）所示。

图4 西昌研究区火点在各阶段的处理效果Fig.4 The results of fire point in each step of Xichang study area

3.2.2.2 火点块筛选

根据火点块的时间特征进行筛选。一个过火区的形成通常是某一时间段内集中发生的火点造成的，远离该时间段内监测到的火点很大可能与形成该过火区无关，可以剔除。根据火点块的探测时间t，统计火点块数量最多的探测时间tmax，tmax即为火势最强的时间。由于火灾强度的变化一般是由弱到强，然后再减弱直至熄灭，且持续时间一般不超过30 d，因此只取时间为tmax前后15 d内的火点块，即保留满足式（6）的火点块：

以西昌研究区为例，筛选后的火点块分布如图4（b）所示。

在此基础上，根据火点块的光谱特征进行筛选。每个火点块内所有像元的最大dNDVI值作为该火点块的dNDVI，记为dNDVI_f。火灾后植被遭到破坏，过火后像元的NDVI 值会明显降低，因此过火像元的dNDVI值会较高。一部分火点像元在燃烧结束后即为过火像元，因此可根据火点块的dNDVI_f进行筛选。首先保留满足下式的火点块：

在此基础上，采用迭代阈值法，确定分割阈值dNDVI_f_th，保留满足下式的火点块：

以西昌研究区为例，筛选后的火点块分布如图4（c）所示。

在此基础上，根据火点像元的空间特征进行筛选。一般情况下，在空间上集中的火点更容易形成过火区。独立的火点有可能是由于异常高温等非过火原因造成的，并不能形成过火区。因此可以删除空间上孤立的火点块，进一步确定过火区在空间中的位置。以西昌研究区为例，剔除孤立火点块后累计火点块的分布如图4（d）所示。累计火点块中的像元为累计火点像元。

3.2.2.3 火点块扩充

在此基础上，扩充累计火点像元，确定过火区分类阈值确定范围和过火区粗略范围。用d表示像元到累计火点像元的距离。若像元满足式（9），则作为过火区阈值确定范围内的像元：

若像元满足式（10），则作为过火区粗略范围的像元：

上述两式中，一般设置dr＞dth，即过火区粗略范围比阈值确定范围要大。针对不同的地表类型，累计火点像元在dth范围内区域增长以后，即可包含一定数量的非过火像元，适用于后续确定过火像元与非过火像元的分割阈值。根据多次实验比较发现在本文研究区中dth取2 km 较为合适。而由于受云覆盖和卫星过境时间周期的影响，过火期间的部分火点可能无法探测到，导致过火区的实际范围可能远大于火点块叠加范围，因此dr的值要大于dth。根据多次实验比较发现在本文研究区中dr取4 km较为合适。以西昌研究区为例，过火区阈值确定范围和过火区粗略范围如图4（d）所示。

3.3 过火区精细提取

3.3.1 确定过火区地表类型

地表类型不同，过火前后的光谱变化可能会有差异。为了提高过火区的提取精度，本方法按地表类型确定过火像元和非过火像元的分割阈值。对于每种二级土地类型l，统计阈值确定范围内的像元个数Nl和统计累计火点像元个数nl，若满足式（11）则l属于过火区土地类型。

式中，pl=()nl/Nl× 100%。pl可以视为每种地表类型的火点像元个数在阈值确定范围内的百分比，当pl过小时，说明地表类型l的火点个数相对较少，可忽略不计。

3.3.2 初步精提取

与3.2.2 中确定dNDVI_f_th方法相同，采用图像分割阈值法中的迭代阈值法，基于dNDVI对每种地表类型确定分割阈值，分类过火像元与非过火像元。迭代法首先根据图像的全局直方图，将取阈值后得到的区域作为子图像，再对各子图像选灰度均值点作为新的阈值再分割图像，重复上述过程，直到新的阈值不再变化或变化很小为止。随着循环次数的增加，这种算法越来越细地考虑了图像的局部特征，可以获得精细的分割（李琳琳，2012），具体的实现步骤如下：

（1）首先求出图像的最大dNDVI和最小dNDVI，分别记为ZMAX和ZMIN，计算初始阈值T为

（2）根据阈值T将图像分割为前景和背景，分别求出两者的平均值Z0和ZB；

（3）求出新阈值T'为

若T=T′，则T所得即为阈值；否则将T'设为T，并从第（2）步开始重复，进行迭代计算。

在过火区阈值确定范围内，采用上述迭代阈值法确定分类阈值dNDVI_l_th。在过火区粗略范围内，地表类型为l的像元满足下式则为过火像元：

所有地表类型的过火像元合集则为过火区初步提取结果。

3.3.3 小斑块去除

过火区初步提取结果周边会有一些细小的斑块。为了使过火区更为完整，采用滤波算法，去除细小的斑块。统计过火像元3×3范围内的火点像元数量Nb，满足式（15）的像元仍为过火像元，否则为非过火像元，重复上述步骤直至不再有新的非过火像元产生。

4 结果处理与分析

根据上述方法，2 个研究区的过火区提取结果如图5 和图6 所示。同时基于过火前后影像，分别选择水体、云、裸土、植被、城市、过火区等类别的样本，用神经网络分类方法提取2个研究区的过火区作为对比。基于2 个研究区的过火后GF-1 WFV 影像，通过人机交互的方式判识过火区作为参考真值。将本方法和神经网络分类法的结果分别于参考真值比较，通过混淆矩阵验证过火区提取精度，精度评价见表3。结合相关概念（许文宁等，2011）和本实验情况，表3中，错分误差是指非过火像元错分为过火像元的误差，漏分误差是指过火像元被分为非过火像元的误差。用Kappa系数作为精度评价指标。

表3 研究区过火区分类精度Table 3 Accuracy of burned area mapping in study areas

图5 木里过火区分类结果图Fig.5 The classification result of burned area in Muli

图6 西昌过火区参考真值、本方法结果以及神经网络的分类结果Fig.6 The classification result of reference true value，results of the proposed method and results of neural network classification of burned area in Xichang

从图5 和图6 可以看出，本文方法结果与人工提取的过火区参考真值在空间上的分布和形状差异较小。从表3 可以看出，本文方法结果的Kappa系数整体较高，达到0.82，说明本文方法的分类精度整体较高。本方法木里研究区的Kappa系数为0.82，要低于西昌研究区的Kappa 系数0.87，说明本文方法对木里过火区的提取精度要低于西昌过火区，主要有2 个原因：一是与西昌过火区相比，木里过火区分布较为疏散，斑块较多。从图5 和图6可以看出，错分误差和漏分误差多出现在过火区斑块的边缘处，因此分布较为疏散的木里过火区的误差较大；二是木里的错分误差较大，达到19.38%，这是因为在木里过火区的分类结果周边有较多的非过火像元构成的细小斑块造成的错分。

与神经网络分类法结果的Kappa 系数相比，可以看出本方法的分类精度要高于神经网络分类法。神经网络分类法木里研究区结果的漏分误差较高，主要漏分了过火区斑块边缘部分；西昌研究区结果的错分误差较大，主要是错分了图像上非过火区的像元，如裸土、山体阴影等。这可能是因为用于神经网络分类的GF-1 WFV 数据仅有4 个波段，可能会混淆相似光谱特征的地物像元；同时神经网络法等监督分类法的分类结果往往与样本的选择有较大关系，分类过程有较多不确定性；此外，在选择样本时，需要根据图像中地物的复杂程度来考虑样本种类，每种样本的数量要相近，空间分布要均匀，因此耗费的人力、时间成本较高。

5 讨论

本文方法充分结合了GF-1 WFV影像和FY-3D MERSI 火点产品的优势和特点，得到较高的过火区提取精度。

FY-3D MERSI火点产品空间分辨率为1000 m，与GF-1 WFV 影像差异较大，在过火区提取中如何能够充分结合利用上述2种数据，取长补短，本方法从以下2个方面进行了考虑：（1）根据火点与过火区形成关系的时间特征、空间特征和光谱特征提取过火区，可以弥补GF-1 WFV 影像在时间和光谱分辨率上的不足。（2）为了更好的结合这2种数据，本方法根据数据特点，由低空间分辨率逐步到高空间分辨率转换：从空间分辨率为1000 m的火点产品切入，然后将空间分辨率为16 m 和1000 m的数据相结合，最后基于空间分辨率为16 m的影像获得结果。这种从低空间分辨率到高空间分辨率，从粗略到精细的分类方法，在充分利用2种数据的同时能最大程度地降低空间分辨率差异性所带来的影响。

6 结论

综上，本文基于GF-1 WFV 影像和FY-3D MERSI 火点，对小范围过火区分类提取，主要结论如下：（1）该方法可以充分将GF-1 WFV 影像和FY-3D MERSI 火点产品相结合，弥补了GF-1 WFV 影像在时间分辨率、光谱分辨率上的不足；（2）该方法降低了分类样本选择的不确定性和成本，能快速准确地对小范围过火区分类提取，进一步提升我国卫星遥感数据火情监测的能力。未来可利用静止卫星的火点产品，结合实地考察，以提升火点产品和参考真值的精度，对本方法改进完善。