一种夜光遥感影像的农田火识别与去除方法

1.本发明属于农田火识别与去除技术领域，具体涉及一种夜光遥感影像的农田火识别与去除方法。


背景技术：

2.随着全球科技和经济的快速发展，人类活动日益活跃，对城市化、环境、政治、文化等方面均施加了复杂影响，因此实时监测与分析相应变化尤为重要。在遥感和对地观测领域，随着对地观测技术的发展，人类对地球的综合观测能力达到空前水平，经由卫星、飞机、无人机等多种非接触式手段获得的、特点各异的多光谱遥感数据具有客观性、实时性与易获得性，使其被广泛应用于监测和分析社会经济环境与自然环境变化的各个领域。在夜间无云条件下，遥感卫星捕捉到城镇灯光、渔船灯火、火点等可见光辐射源的影像即为夜光遥感影像，相比于日间遥感，夜光遥感具有独特的、反映人类社会活动的能力，在估算社会经济参数、评估城市化水平、评估受灾影响等方面具有独特优势。尤其是在欠发达地区，当地政府缺乏开展精确统计与分析的资源和能力，导致统计数据源匮乏、缺乏可信数据，使得相应的社会经济评估工作难以开展。夜光遥感影像作为可信度高、易获得、经济成本低的数据，能够以较低的成本为欠发达地区开展社会经济参数评估等工作，为这类地区今后的发展方向提供有利的数据支撑和重要的决策依据。
3.在这类欠发达地区，农业是最重要的产业之一，农业用地所占面积较高。以非洲为例，其70％的人口是农户，但农业现代化水平落后于世界其他地区，仍然保留着“刀耕火种”的耕作方法。农户们在播种前，放火烧光田地的植物以增加草木灰、提高土地肥力。这导致了在夜间灯光影像上能够观察到这些地区在播种时间里出现频次高、范围广的农田火点。这些农田火噪声和人工灯光混杂在一起，将会影响夜光影像在估算贫困率、电力普及率等社会经济参数时的精度。目前尚未出现一套已去除农田火噪声的夜光数据，虽然已有多代卫星传感器的数据用于火产品的生产，如活跃火(active fire,af)和夜间火(viirs night fires,vnf)产品，但是其算法一脉相承，基于单景影像的近红外、中红外的亮度温度确定火灾，而没有考虑某些特定类型火的季节性特征，如农田火，且现有技术通常对火灾的亮度温度识别阈值为300k以上。而研究发现，在耕地位置上发生辐亮度季节性增长的农田火像素，其亮度温度绝大部分在290k以下，且nasa发布的活跃火数据基本覆盖不了这些像素点，但这些农田火像素确实降低了夜光影像的质量。


技术实现要素：

4.本发明提出了一种夜光遥感影像的农田火识别与去除方法，，根据农田火的季节性特征、其发生地与耕地的高度耦合关系，以及农田火发生时的夜光辐亮度与亮度温度的变化，并采用时间上的最邻近重采样去除农田火的异常高值，重新生成多时相的夜光影像。
5.为了达到以上目的，本发明所采用的技术方案如下。
6.一种夜光遥感影像的农田火识别与去除方法，包括以下步骤：
7.步骤1、获取研究区域的每日夜光数据；
8.步骤2、对研究区域的每日夜光数据进行预处理获取批量数据集；
9.步骤3、对批量数据集进行筛选，获取高质有效的时间序列夜光辐亮度和时间序列亮度温度；
10.步骤4、以中心窗口和一定大小的邻域窗口的均值作为最终结果，提取研究区域每个像素点的时间序列夜光辐亮度和时间序列亮度温度；
11.步骤5、处理土地覆盖产品，使土地覆盖产品与研究区域的每日夜光数据范围对齐；
12.步骤6、基于耕地比例及时间序列夜光辐亮度和时间序列亮度温度，提取研究区域每个像素点的多元特征；
13.步骤7、采集训练样本；
14.步骤8、基于训练样本得到分类结果，得到农田火点的分布；
15.步骤9、识别出农田火点后，找到并去除农田火点在时间序列上的异常高值，并将滤除火点后的时序夜光辐亮度重新生成多时相的夜光影像。
16.进一步地，所述步骤1中，每日夜光数据为nasa发布的viirs夜光数据black marble中vnp46a1与vnp46a2。
17.进一步地，在步骤2中，对每日夜光数据进行格式转换获取批量数据集，且提取brdf纠正的夜光辐亮度、m13亮度温度、太阳天顶角、强制质量标志和云层掩码的数据图层。
18.进一步地，在步骤3中，筛选条件为：夜光辐亮度和亮度温度的dn值为1到65534的天数，太阳天顶角大于101
°
的天数，强制质量标志为“非常清晰”的天数，云层掩码为“无云”的天数。
19.进一步地，在步骤4中，以每个像素点为中心窗口，辐射周围多个领域窗口，以天为单位判断中心窗口与邻域窗口是否满足筛选条件，若中心窗口与多个邻域窗口均满足筛选条件，则对中心窗口与多个邻域窗口进行均值计算，获取每个像素点当天的夜光辐亮度和亮度温度，若中心窗口与多个邻域窗口不完全满足筛选条件，则将当前像素点当天的夜光辐亮度和亮度温度视为无效值，提取研究区域的所有像素点的高质量时间序列夜光辐亮度和时间序列亮度温度。
20.进一步地，在步骤5中，将对应范围内的globeland30重投影到与vnp46a1与vnp46a2夜光数据相同的坐标系，进而将土地覆盖产品根据vnp46a1与vnp46a2夜光数据的影像块范围进行镶嵌和裁剪。
21.进一步地，在步骤6中，将一年内所有天数的夜光影像堆叠在一起，基于每个像素点读取其时间序列夜光辐亮度、时间序列亮度温度和每个像素点对应的耕地比例，将像素点分为农田火点、稳定发光点和全黑点，其中，农田火点在全年大部分天数夜光辐亮度基本为0或略大于0，但在播种季节会产生较高夜光辐亮度；稳定发光点在全年均有较高夜光辐亮度且夜光辐亮度稳定；全黑点在全年夜光辐亮度近似为0且变化波动很小。
22.进一步地，多元特征为夜光辐亮度的最大值f1、夜光辐亮度大于一定阈值的比例f2、夜光辐亮度突变率f3、亮度温度最大值f4、夜光辐亮度最大值对应的亮度温度f5和夜光数据的每个像素中的耕地比例f6，识别农田火点、稳定发光点和全黑点，夜光辐亮度的最大值f1的计算公式为：
23.f1＝max(radiance.ts)
ꢀꢀ
(1)
24.其中，radiance.ts是指像素点的有效时间序列夜光辐亮度，全黑点的夜光辐亮度最大值在0.5nw
·
cm-2
·
sr-1
以下，发光点和农田火点均大于1nw
·
cm-2
.sr-1
；
25.夜光辐亮度大于某一阈值的比例f2的计算公式为：
[0026][0027]
其中，gtd为夜光辐亮度大于阈值的天数，vdc指所有的有效天数，对于稳定发光点，夜光辐亮度大于一定阈值的比例f2在0.9以上，对于农田火点和全黑点，夜光辐亮度大于一定阈值的比例f2会低于0.5；
[0028]
夜光辐亮度突变率f3的计算公式为：
[0029]
f3＝(max
in-max
out
)/max
out
ꢀꢀ
(3)
[0030]
其中，max
in
指火点时间窗口内的夜光辐亮度最大值，max
out
指火点时间窗口以外的夜光辐亮度最大值；特殊情况下：若max
out
为0，若max
in
也为0，说明没有突变，将夜光辐亮度突变率f3赋值为0；若max
in
不为0，说明存在变化，将夜光辐亮度突变率f3赋值为max
in
，全黑点和稳定发光点f3的绝对值在1以下，农田火点f3为1以上的较大值；
[0031]
亮度温度最大值f4的计算公式为：
[0032]
f4＝max(bt.ts)
ꢀꢀ
(4)
[0033]
其中，bt.ts指的是有效亮度温度的时间序列；
[0034]
农田火点与稳定发光点亮度温度最大值f4与夜光辐亮度最大值对应的亮度温度f5较大，而全黑点亮度温度最大值f4与夜光辐亮度最大值对应的亮度温度f5较小；
[0035]
夜光数据的每个像素中的耕地比例f6的计算公式为：
[0036][0037]
其中，plou
num
指的是一个夜光像素范围内，土地覆盖产品上对应的耕地像元的总数，而num指的是一个夜光像素范围内，土地覆盖产品上对应的像元总数，夜光数据的每个像素中的耕地比例f6与该夜光像素是否是耕地有关，如果是耕地，则夜光数据的每个像素中的耕地比例f6对应的值就越大，反之则越小，而一块像素是耕地与该块像素会发生农田火正相关。
[0038]
进一步地，在步骤8中，采用随机森林分类方法实现了分类算法，得到分类结果，得到农田火点的分布。
[0039]
进一步地，在步骤9中，滤除农田火点异常高值时，执行以下子步骤：
[0040]
步骤9.1，找到异常高值，对农田火像素一年内的时间序列夜光辐亮度进行搜索，找到农田火发生窗口内夜光辐亮度高于农田火发生窗口之外的夜光辐亮度最大值的天数；
[0041]
步骤9.2，最邻近采样，将步骤9.1搜寻到的天数对应的夜光辐亮度都赋值为其最邻近的非异常高值的天数对应的夜光辐亮度；
[0042]
步骤9.3，重新生成夜光影像，将经过步骤9.1和步骤9.2处理的夜光辐亮度重新生成多时相的夜光影像。
[0043]
本发明具有以下优点和有益效果：
[0044]
1.本发明提出的一种夜光遥感影像的农田火识别与去除方法，通过时间序列夜光
辐亮度和亮度温度，在日级时间分辨率上提取了农田火的季节特征，可以更好地识别出夜光影像上的农田火。
[0045]
2.本发明提出的一种夜光遥感影像的农田火识别与去除方法，顾及了农田火的发生位置与耕地的高耦合度，引入夜光像素的耕地比例这一特征，使得识别的农田火区别于森林、草地等其他野火。
[0046]
3.本发明提出的一种夜光遥感影像的农田火识别与去除方法，计算速度快，除了标记训练样本外，其他环节均可以进行全自动处理，减少人为干预，可以实现对一整个影像块一年的每日影像进行批量读取、识别与去除农田火点。
附图说明
[0047]
此处所说明的附图用来提供对本发明的进一步理解，构成本技术的一部分。本发明的示意性实施例及其说明用于解释本发明，并不构成对本发明的不当限定。在附图中：
[0048]
图1为本发明的一种夜光遥感影像的农田火识别与去除方法流程示意图；
[0049]
图2为本发明实施例中的夜光影像；
[0050]
图3为本发明实施例中的农田火点的特点示意图；
[0051]
图4为本发明实施例中的全黑点的特点示意图；
[0052]
图5为本发明实施例中的稳定发光点的特点示意图；
[0053]
图6为本发明实施例中农田火点异常高值去除前后的时间序列夜光辐亮度曲线。
具体实施方式
[0054]
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述。
[0055]
一种夜光遥感影像的农田火识别与去除方法，如图1所示，包括以下步骤：
[0056]
步骤1、获取研究区域的每日夜光数据；
[0057]
步骤2、对研究区域的每日夜光数据进行预处理获取批量数据集；
[0058]
步骤3、对批量数据集进行筛选，获取高质有效的时间序列夜光辐亮度和时间序列亮度温度；
[0059]
步骤4、以中心窗口和一定大小的邻域窗口的均值作为最终结果，提取研究区域每个像素点的时间序列夜光辐亮度和时间序列亮度温度；
[0060]
步骤5、处理土地覆盖产品，使土地覆盖产品与研究区域的每日夜光数据范围对齐；
[0061]
步骤6、基于耕地比例及时间序列夜光辐亮度和时间序列亮度温度，提取研究区域每个像素点的多元特征；
[0062]
步骤7、采集训练样本；
[0063]
步骤8、基于训练样本得到分类结果，得到农田火点的分布；
[0064]
步骤9、识别出农田火点后，找到并去除农田火点在时间序列上的异常高值，并将滤除火点后的时序夜光辐亮度重新生成多时相的夜光影像。
[0065]
本发明提供一种夜光遥感影像的农田火识别与去除方法，通过时光序列夜光辐亮度和时光序列亮度温度，在日级时间分辨率上提取了农田火的季节特征，使得可以更好地
识别出夜光影像上的农田火。
[0066]
本发明提供一种夜光遥感影像的农田火识别与去除方法，顾及了农田火的发生位置与耕地的高耦合度，引入夜光像素的耕地比例这一特征，使得识别的农田火区别于森林、草地等其他野火。
[0067]
通过步骤1能够获取所需的数据，在步骤1中，每日夜光数据为nasa发布的viirs夜光数据black marble中vnp46a1与vnp46a2，下载得到研究区域的vnp46a1与vnp46a2每日夜光数据以及涵盖了研究区域的globeland30土地覆盖产品，其中，研究区域的vnp46a1与vnp46a2每日夜光数据数据来源：https://ladsweb.modaps.eosdis.nasa.gov/search/，涵盖了研究区域的globeland30土地覆盖产品数据来源：http://www.globallandcover.com/。
[0068]
本发明以nasa发布的viirs夜光数据black marble中vnp46a1与vnp46a2夜光数据作为研究数据，vnp46a1与vnp46a2夜光数据具有日级时间分辨率、校正了大气辐射、地形、月光、杂散光的影响。
[0069]
南部非洲的大陆部分包括赞比亚、莫桑比克、南非、安哥拉等10个国家，是非洲重要的农业耕种地区，出现农田火的概率比较大。涵盖南部非洲10个大陆国家的vnp46影像块的区号分别是：h19v09、h20v09、h21v09、h19v10、h20v10、h21v10、h22v10、h19v11、h20v11、h21v11、h19v12、h20v12、h21v12。
[0070]
在本发明的一个实施例中，选择的研究年份为2018年，每个影像块各365天(数据缺失时可能少于365天)的vnp46a1和vnp46a2夜光数据。
[0071]
考虑到程序批处理的便利性，将数据进行格式转换，在步骤2中，对每日夜光数据进行格式转换获取批量数据集，将vnp46a1与vnp46a2数据从hdf5格式转换为geotiff格式，且提取brdf纠正的夜光辐亮度、m13的亮度温度、太阳天顶角、强制质量标志和云层掩码的数据图层，其中，夜光辐亮度和亮度温度作为输入数据源，采用太阳天顶角、强制质量标志和云层掩码作为质量筛选，所用图层的具体信息如表1所示：
[0072]
表1本发明中涉及的vnp46a1与vnp46a2夜光数据的图层信息
[0073][0074]
步骤3能够提高分类的可靠性，筛选出高质量数据，在步骤3中，将步骤2中转换的太阳天顶角、强制性质量标志和云层掩码作为质量筛选图层，使得每个像素筛选出高质有效的夜光辐亮度和亮度温度。基于数据的特性，本发明对批量数据集的筛选过程如下所示：
[0075]
用gdal开源库对vnp46数据转换出来的5个波段进行批读取，需要进行如下处理：vnp46a1与vnp46a2夜光数据各个波段都做了一些转换，将值转换为整型，压缩了影像存储空间，因此在读取时，根据用户手册提供的系数，将值转换回来；筛选条件为：夜光辐亮度和亮度温度的dn值为1到65534的天数，太阳天顶角大于101
°
的天数，强制质量标志为“非常清晰”的天数，云层掩码为“无云”的天数。
[0076]
在步骤4中，以每个像素点为中心窗口，辐射周围多个领域窗口，以天为单位判断中心窗口与邻域窗口是否满足筛选条件，若中心窗口与多个邻域窗口均满足筛选条件，则对中心窗口与多个邻域窗口进行均值计算，获取每个像素点当天的夜光辐亮度和亮度温度，若中心窗口与多个邻域窗口不完全满足筛选条件，则将当前像素点当天的夜光辐亮度和亮度温度视为无效值，提取研究区域的所有像素点的高质量时间序列夜光辐亮度和时间序列亮度温度。
[0077]
在步骤4中，为了后续提取多元特征的便利，基于步骤3的筛选条件，提取每个像素点的有效的时间序列夜光辐亮度和时间序列亮度温度。提取的每个像素点的时间序列夜光辐亮度和时间序列亮度温度需要以一定范围的均值作为最终结果，因为vnp46a2数据的夜光辐亮度存在像素偏移的问题。
[0078]
在本发明的一个实施例中，所采用的辐射窗口大小为3
×
3，本发明以每个像素点为中心窗口，辐射周围8个领域窗口，对中心窗口和8个领域窗口进行均值计算，获取每个像素点的时间序列夜光辐亮度和时间序列亮度温度温度均值。
[0079]
夜光辐亮度的均值处理公式为：
[0080][0081]
其中，ra
′
r,c
是指第r行第c列像素在均值处理后的夜光辐亮度，ra(i,j)为均值处理前的辐亮度；
[0082]
亮度温度的均值处理公式为：
[0083][0084]
其中，bt
′
r,c
是指均值处理后的亮度温度，bt(i,j)为均值处理前的亮度温度。
[0085]
在步骤5中，处理土地覆盖产品，使土地覆盖产品与研究区域的每日夜光数据范围对齐，将对应范围内的globeland30重投影到与vnp46a1与vnp46a2夜光数据相同的坐标系，即从utm分度带投影变为wgs1984投影，然后将土地覆盖产品根据vnp46a1与vnp46a2夜光数据的影像块范围进行镶嵌、裁剪。土地覆盖产品为每一类地物赋了不同的dn值，耕地类型的dn值为10，该产品的空间分辨率为30米，而夜光数据的空间分辨率为500米，计算每个夜光像素中的耕地比例。
[0086]
每个夜光像素中的耕地比例的计算公式：
[0087][0088]
其中，plo指的是某夜光像素中的耕地比例，plonum指的是该夜光像素范围对应的
土地覆盖产品上dn值为10的像素数目。
[0089]
在步骤6中，提取分类所需的多元特征，提取的多元特征如下：
[0090]
将一年内所有天数的影像堆叠在一起，基于每个像素点读取其时间序列夜光辐亮度、时间序列亮度温度和每个像素点对应的耕地比例。
[0091]
首先，在经过对月平均夜光影像的目视评估和对夜光辐亮度一年变化曲线的判别之后，本发明将夜光影像上的点分为三类：
[0092]
(1)农田火点：在全年大部分天数夜光辐亮度基本为0或略大于0，但在播种季节会产生较高夜光辐亮度的像素点区域。
[0093]
(2)稳定发光点：全年均有较高夜光辐亮度且夜光辐亮度稳定的像素点，这类像素点可被认为是人类活动较为频繁密集的城市化区域。
[0094]
(3)全黑点：全年夜光辐亮度近似为0且变化波动很小的像素点，这类像素点可被认为是自然植被、裸地等没有人类活动和农田火焚烧的区域。
[0095]
其次，根据研究区的农作物生长时间确定农田火出现的时间窗口。本发明以南部非洲地区为例，该地区主要种植的农作物有谷物(玉米、小麦、稻谷、小米等)、油料(花生、橄榄、油棕果)、棉花、蔬菜(辣椒、秋葵、洋葱、番茄等)、水果(芭蕉、山竹、芒果、西瓜等)、坚果(开心果、腰果、榛果、核桃等)。本发明基于该地区的主要农作物播种时间，将6月到11月确定为农田火点出现的时间窗口。
[0096]
本发明基于像素的一年时间序列信息，提取了以下6个特征(表2)：
[0097]
表2用于农田火识别的多元特征变量
[0098]
[0099]
(1)夜光辐亮度的最大值f1的计算公式为：
[0100]
f1＝max(radiance.ts)
ꢀꢀ(i)[0101]
其中，radiance.ts是指像素点的有效时间序列夜光辐亮度，全黑点的夜光辐亮度最大值在0.5nw
·
cm-2
·
sr-1
以下，明显小于稳定发光点和农田火点，发光点和农田火点均大于1nw
·
cm-2
·
sr-1
。夜光辐亮度的最大值f1能够将全黑点与发光点和农田火点区分开来。
[0102]
(2)夜光辐亮度大于一定阈值的比例f2，即辐亮度大于阈值的天数占所有有效天数的比例，夜光辐亮度大于某一阈值的比例f2的计算公式为：
[0103][0104]
其中，gtd为夜光辐亮度大于阈值的天数，vdc指所有的有效天数。对于稳定发光点，夜光辐亮度大于一定阈值的比例f2在0.9以上，对于农田火点和全黑点，夜光辐亮度大于一定阈值的比例f2会低于0.5。夜光辐亮度大于一定阈值的比例f2能够区分稳定发光点和农田火点和全黑点。
[0105]
本发明经过对2000多个稳定发光点、3000多个黑点和2000多个农田火点的时间序列曲线进行查看，最终将这一阈值确定为1nw
·
cm-2
·
sr-1
。
[0106]
(3)夜光辐亮度突变率f3，即夜光辐亮度在时间窗口内的最大值与在时间窗口外的最大值的归一化差异，夜光辐亮度突变率f3的计算公式为：
[0107]
f3＝(max
in-max
out
)/max
out
ꢀꢀ
(3)
[0108]
其中，max
in
指火点时间窗口内的夜光辐亮度最大值，max
out
指火点时间窗口以外的夜光辐亮度最大值；特殊情况下：若max
out
为0，若max
in
也为0，说明没有突变，将突变率赋值为0；若max
in
不为0，说明存在变化，将突变率赋值为max
in
。如果夜光辐亮度没有发生突然的跃迁增长，一直在平稳波动，那么f3的绝对值在1以下，但对于农田火点这样会发生突然性增长的辐亮度曲线，夜光辐亮度突变率f3则会产生1以上的较大值，因此可以将农田火点与全黑点和稳定发光点分开。
[0109]
(4)亮度温度最大值f4的计算公式为：
[0110]
f4＝max(bt.ts)
ꢀꢀ
(4)
[0111]
其中，bt.ts指的是有效亮度温度的时间序列。
[0112]
(5)夜光辐亮度最大值对应的亮度温度f5。
[0113]
亮度温度与地表温度和辐射能量都有关系，发生火点或辐亮度很高都可能会导致较高的亮度温度。因此农田火点与稳定发光点都会对应较高的亮度温度，而黑点的亮度温度会显著低于前两者，亮度温度最大值f4与夜光辐亮度最大值对应的亮度温度f5可以将黑点与农田火点与稳定发光点分开。
[0114]
(6)夜光数据的每个像素中的耕地比例f6：根据土地覆盖数据，计算夜光数据的每个像素中的耕地比例f6：
[0115][0116]
其中，plou
num
指的是一个夜光像素范围内，土地覆盖产品上对应的耕地像元的总数，而num指的是一个夜光像素范围内，土地覆盖产品上对应的像元总数，夜光数据的每个像素中的耕地比例f6与该夜光像素是否是耕地有关，如果是耕地，则夜光数据的每个像素
中的耕地比例f6对应的值就越大，反之则越小，而一块像素是耕地与该块像素会发生农田火正相关。
[0117]
在步骤7中，如图2-图5所示，采集训练样本，采集的样本基于npp/viirs的月度夜光影像、土地覆盖数据以及google earth高分辨率日间遥感影像。首先人工采集了影像上每个月都不发光的全黑点(在采集时，本发明认为全年夜光辐亮度都低于1nw
·
cm-2
·
sr-1
的点可以作为全黑点的备选)；在火点时间窗口内出现过较高夜光辐亮度，而在其余时间为黑点的农田火点(本发明认为火点窗口内出现5nw
·
cm-2
·
sr-1
以上的夜光辐亮度，而火点窗口外夜光辐亮度都低于1nw
·
cm-2
·
sr-1
，可以作为农田火的备选)；每个月都发光的稳定发光点(本发明认为，全年夜光辐亮度值基本都高于5nw
·
cm-2
·
sr-1
，可以作为稳定发光点的备选)各多个。
[0118]
本发明的一个实施例中，取全黑点、农田火点和稳定发光点各2000个，其次，会将这些点与土地覆盖数据以及google earth对照，以弥补人工采集可能存在的误差。最后提取通过对照的样本点的多元特征，作为接下来进行机器学习分类的训练样本。
[0119]
在步骤8中，如图2-图5所示，农田火识别方法流程如下：将夜光影像上的样本点预先分为3类(即农田火点、全黑点、稳定发光点)。训练集为采集的样本点提取的步骤6中的多元特征，获得多元特征矩阵。测试集为每个影像块上每一个有值输出像素的多元特征向量。本发明采用随机森林方法实现了分类算法，并采用网格搜索进行参数优化，本发明使用实例中随机森林的参数：决策树棵数为10，自举重采样，其他参数不做约束，采用程序包的默认参数。
[0120]
在步骤9中，如图6所示，为了去除农田火点的异常夜光辐亮度，滤除农田火点异常高值时，执行以下子步骤：
[0121]
步骤9.1，找到异常高值，对农田火像素一年内的时间序列夜光辐亮度进行搜索，找到农田火发生窗口(即每年的6月到11月)内夜光辐亮度高于农田火发生窗口之外(即每年的1月-5月以及12月)的夜光辐亮度最大值的天数；
[0122]
步骤9.2，最邻近采样，将第一步搜寻到的天数对应的夜光辐亮度值都赋值为其最邻近的非异常高值的天数对应的夜光辐亮度值。
[0123]
步骤9.3，重新生成夜光影像，将经过前两步处理的夜光辐亮度重新生成多时相的夜光影像。
[0124]
具体实施时，本发明提供的流程可采用软件技术实现自动运行。
[0125]
显然，本领域的技术人员可以对本发明进行各种改动和变型而不脱离本发明的精神和范围。这样，倘若本发明的这些修改和变型属于本发明权利要求及其等同技术的范围之内，则本发明也意图包含这些改动和变型在内。