李浩然,刘子菁,邱巧勇,何永宁,赵文化,魏永顺
(1.博罗县气象局,广东博罗 516000;2.广州气象卫星地面站,广东广州 510640)
森林是植物分布密集的区域,对地球生态系统的平衡有着举足轻重的作用。森林火灾使得森林蓄积下降,引发水土流失并危害人类健康,给人民的生命财产安全带来巨大的威胁[1]。全世界每年平均发生森林火灾近30万起,森林防火的重要性毋庸置疑[2]。森林过火面积是对森林火灾进行描述的重要参考指标之一,国内外开展了许多这方面的研究。Pereira[3]使用AVHRR的单通道阈值法来进行火点的判识;赵文化等[4]也使用EOS/MODIS数据来进行森林火灾的监测与面积估算等方面的研究。然而极轨气象卫星最大空间分辨率为250 m,这极大地影响了使用该类卫星进行过火区域面积计算的精度。贾永红[5]指出,当两种数据的空间分辨率差距太大时,很难将两种不同分辨率的数据插值在同一网格中进行计算。
高分一号(GF-1)卫星是中国高分辨率对地观测系统中的首发星,高分数据广泛应用于灾害预警、林火监测等方面。杨斌[6]基于 GF-1与Landsat-8不同分辨率的数据分别提取归一化植被指数(normalized vegetation index,NDVI),并用来计算四川茂县山林的植被覆盖度。随着人工智能和图像识别的兴起与发展,图像识别也应用在了更广泛的领域。邹成明等[7]基于图像色彩直方分布提出一种图片相似度比较和特征计算算法;郭仪权[8]使用图像哈希算法比较连续时刻的图片来实现目标追踪;黄嘉恒等[9]对比3种不同的哈希算法在图像相似度中的识别效果。因此,将图像识别的方法引入过火区域的判识与计算,有着广阔的应用前景。
选取广东省佛山市高明区作为研究区域,该区域地处东经 112°45′—112°92′,北纬 27°73′—30°02′。佛山地区全年降水量接近2 000 mm,但降水分布不均,其中12月降水量往往不足50 mm。而冬季植物蒸发蒸腾量的作用,使得季节性、区域性气象干旱也时有发生[10]。高明区森林火灾发生于2019年12月5日14:00(北京时,下同)左右,林火区域位于高明区新城村南侧的凌云山,火场区域为马尾松、灌木林交错覆盖的山地地形。12月8日11:30,凌云山外线已基本看不到明火点,火灾形势得到有效控制。
GF-1号高分辨率极轨卫星所携带的多光谱相机星下点分辨率最高可达2 m。此外还搭载了1台分辨率为8 m的多光谱相机、4台分辨率16 m的多光谱相机,其光谱范围覆盖可见光与近红外波段,可以满足火情监测与过火面积估算的需要。GF-6号的轨道参数和GF-1号基本相似。
本研究分别采用归一化植被指数方法和图像识别方法对过火区域面积进行计算。在使用归一化植被指数方法计算时,先在非过火区域选定NDVI变化阈值,再在过火区域内用阈值判断过火像元,最终计算出过火面积;在使用图像识别方法进行计算时,提取过火前后两张图片的特征生成图片指纹,最终比较两张图片中不同的比例来计算过火面积。最后将两种方法得到的过火面积与真值进行对比。
当树木燃烧后,树叶中的叶绿素遭到破坏,使得过火区域在近红外波段的反射率下降,光谱特征曲线发生改变。利用过火区域和未过火区域在近红外通道的反射率差异来进行过火区的判识是最常用的原理[11]。NDVI参数能够反映观测区域红光通道与近红外通道的反射率,可用来评判植被长势的重要因子之一。NDVI的表达式为
其中,RRed代表红光通道的反射率;RNir代表近红外通道的反射率。通过分析过火区域前后的NDVI差异,确定NDVI阈值,即可以判断出过火区域的像元数量以及位置。在该次研究中,将GF-1号卫星WFV相机通道4的数据作为RRed,通道5的数据作为RNir,来计算选取区域内各点的NDVI值。该次研究中通过对2019年12月8日GF-1B的1景图像进行分析,大致确定火场位置,并选定 112°71′E—112°79′E,22°90′N—23°00′N区域作为火场区域,将其定义为区域A。通过判别式(2)来确定过火区域:
其中,NDVI1208和NDVI1124分别为12月8日和11月24日区域A各点NDVI值。TNDVI为判别阈值,可利用区域A附近的下垫面类型较为相似的非起火区域来确定,以消除由于大气和植被生长造成的NDVI值变化的影响。本研究选取A区域西北方向约20 km的羚羊峡的一片山林植被作为非过火区域 B,范围为 112°58′E—112°61′E,23°09′N—23°12′N(图 1)。将区域 B过火时间前后的NDVI变化最大值作为阈值,以此来确定自然状态下NDVI变化的极值。当区域A内像元的NDVI变化大于阈值时,则认为该像元为过火像元。
图1 选取的过火区与非过火区
首先选取2019年11月24日的GF-1B和12月8日GF-6号的2轨卫星数据,对PMS数据进行辐射定标校正,FLAASH大气校正与正射校正,以保证两颗卫星的误差在1个像元以内。区域B内过火时间前后NDVI最大减少值为0.202,将0.202作为阈值,计算火灾前后区域A内各点NDVI的值(图2)并判断是否为过火像元。
图2 火灾前(a)和后(b)NDVI值
在本研究中共将184 255个像素点判识为过火像素点。由于受到分辨率的限制,因此单一使用像元个数乘以像元面积的算法会使得估算的过火面积大于实际过火面积。因此通过计算植被覆盖度来减小裸露的地表在过火像素点计算中造成的影响。通过计算公式得到过火面积公式:
其中,S为过火区域面积(m2);n为像素数;Si为每个像元的面积(m2);NDVImix为区域内某一像元土壤与植被NDVI的混合值,NDVIs为区域内NDVI土壤端的值,一般为区域内NDVI极小值,NDVIv为区域内NDVI植被端的值,一般为区域内NDVI极大值。在本研究中,Si为 64 m2,NDVIs和 NDVIv分别为 0.105与 0.866,通过计算得出的过火面积为11.737 km2。
随着人工智能技术的发展,使用计算机自动识别图像并判断图像相似度的方法得到了越来越广泛的应用。使用人工智能进行图像识别的优点在于直接从图像本身来寻找过火后的图像变化,可以不受卫星仪器波长通道设置和空间分辨率的影响,对不同种类的卫星在同一地区的图像进行图像识别,找到连续两张图像中的不同,进而计算出过火面积。
本研究使用感知哈希算法,通过将图像转换为灰度图像,获取图片的hash值,再比较两张图片hash值(指纹)的不同位的个数(汉明距离)来判断两张图片的相似程度。两张图片的指纹汉明距离越小,则证明两张图片越相似。
感知哈希算法对图片相似度的处理过程如下:(1)遍历像素点。通过遍历图中32×32的像素点,可以简化离散余弦变换的计算。(2)将图片转换为灰度图像并计算图片中所有像素的灰度平均值。(3)计算hash值。根据离散余弦变换得到的N×N矩阵,将其设置成二进制的hash值。将二进制数字组合在一起,即为这张图片的指纹。
在哈希算法中,可通过以下的方法来计算二维离散余弦变换(DCT):
当 u≥0,v<N时,α(u)=α(v)。
其中,g为N×N图像像素点,一般为32×32;G为矩阵中阈矩阵,α为余弦系数矩阵。通过对两张图片进行处理,找到其指纹不同的地方,便可以计算出两张图片中不同的面积比例。
为了尽量剔除准备自然生长的影响,尽可能选取相隔时间短的火灾前后的晴空图像进行对比。选取2019年11月24日的GF-1和12月8日GF-6号的2景图像,研究区域与第2章相同,同样进行地理投影校正,以保证其误差在一个像元以内(图3)。
图3 GF-1号(a)和6号(b)通道合成图像
将两幅图片使用上述感知哈希算法进行相似度判识,保留高分辨率图像所带有的低频信息,得到了过火区的二值图像(图4)。感知哈希算法计算过火前后两张图片相似比例为14.76%,所选择区域面积为 73.618 km2,因此感知哈希算法计算的过火面积为10.866 km2。
图4 过火区域二值图
该次研究选择佛山市高明区的山火作为研究对象,以火灾前后的高分卫星数据为基础,通过未过火区域NDVI的变化情况确定阈值进而判断并确定火灾区域,并计算出过火区域面积。通过3种方法的哈希相似度比较,将火灾前后的遥感图像进行对比并找到图像中的不同部分。通过计算两张图片的汉明距离来计算两张图片中不同部分的比例,从所使用的遥感图像总面积来计算过火区域的面积。该次火灾后,佛山市林业局公布的火灾过火面积为10.667 km2,通过NDVI的变化计算得出的过火面积为11.737 km2,与实际的误差为10.03%;通过感知哈希算法计算的过火面积为10.866 km2,与实际的误差为2.05%。说明感知哈希算法计算出的过火面积可靠性较好,可以应用于火灾后过火面积的估算。
本研究所提出的结合高分辨率卫星和图像识别的方法可以充分结合分辨率较高和人工智能计算速度快、精度高的特点,在提高遥感手段进行林火监测能力的同时,保证了过火区域面积计算的精度。同时,该方法可将 GF-4、Himawari-8、Landsat-8等不同轨道、不同传感器的卫星数据与图像结合起来。既发挥气象、环境卫星的高频次观测,又发挥对地观测高分卫星的高空间分辨率的特点,做到卫星遥感的大数据融合观测。可以为森林火灾的防灾减灾提供更加及时准确的服务信息。