CN114092807A - 一种基于时段累计蒸散发的实际灌溉面积识别方法 - Google Patents

Abstract

本发明提供一种基于时段累计蒸散发的实际灌溉面积识别方法，步骤1、计算农田范围内像元逐日灌溉蒸散发量；步骤2、分析不同时长的累计灌溉蒸散发量概率分布直方图，确定灌溉识别的最短时长窗口T*；步骤3、利用T*对应的时段累计灌溉蒸散发量概率分布直方图，分析确定灌溉与非灌溉像元的阈值ET*；步骤4、用选定的阈值ET*对按T*时长累加的灌溉蒸散发数据进行划分，累计灌溉蒸散发量大于ET*的像元为灌溉像元，否则为非灌溉像元；步骤5、当时长窗口T大于T*时，重复步骤3、4，得到T时长内的实际灌溉面积。本发明提供高植被覆盖和降水影响实际灌溉面积影响识别精度，给出确定基于蒸散发量识别灌溉面积的最短时长窗口的方法。

Description

一种基于时段累计蒸散发的实际灌溉面积识别方法

技术领域

本发明涉及农业实际灌溉面积识别技术领域，特别涉及一种基于时段累计蒸散发的实际灌溉面积识别方法。

背景技术

在全球气候变化和社会经济发展的双重影响下,水资源短缺问题日益凸显，成为社会经济可持续发展的主要瓶颈之一。农业用水占我国用水总量的62％，灌溉用水占农业用水的90％以上。实际灌溉面积是农业用水统计和现代农业管理的重要支撑信息，对监测灌区灌溉系统、优化灌区水资源利用、预防干旱、优化作物种植结构和水资源调控配置都具有重要意义。传统的灌溉面积监测方法存在着信息获取困难，监测站点少获取信息不全面，费时耗力，更新时间长，不能按灌溉轮次及时提供灌溉面积数据的问题，不能满足现代化农业灌区管理所需的数据支撑。卫星遥感技术作为地球信息科学的前沿技术，被认为是现阶段最有效的对地观测技术以及最全面的获取地面信息的手段之一。具有监测覆盖范围广、时效性高、数据客观、经济效益高等优势。应用遥感数据进行灌溉面积识别可以改善识别的时效性和可获取性。

现行的基于遥感技术的灌溉面积识别方法大体可以分为两类，一类是基于遥感种植结构的灌溉面积识别方法，一类是基于遥感土壤含水量的灌溉面积识别方法。前者是基于种植结构与灌溉的相关性进行分析识别。适用于特定的气候条件(干旱条件)，主要田间水分补给由灌溉提供的区域，通过结合种植结构与作物长势进行分析可以得到精度较高的的实际灌溉面积识别结果，但是时效性较差，不能分析出灌溉的灌次。

基于遥感土壤含水量的灌溉面积识别方法，是基于土壤含水量变化对灌溉的响应。虽然土壤含水量变化对灌溉的响应过程迅速、精度较高，但是响应过程容易受到其它水分补充(如降水)的影响，且土壤含水量反演精度受到植被生长的影响。这些因素都影响了基于土壤含水量变化识别灌溉面积的方法的精度。

蒸散发是全球水循环中的重要环节，区域的大部分水资源都会通过蒸散发消耗进入大气。对于农田水量平衡，大部分农田的补充水量都会以蒸散发的形式消耗。故在农业节水管理中蒸散发通常被认为是重要的管理支撑数据。遥感技术的发展使得利用遥感数据反演区域蒸散发量得以实现，估算的区域蒸散发量通常被认为是作物的耗水量(需水量)，被用于推算灌溉需水量、评价灌溉效率和作物水分监测与估产。蒸散发在农田节水管理中的应用都基于蒸散发对灌溉有较为显著的响应过程。蒸散发对灌溉的响应有潜力应用于实际灌溉面积识别。

发明内容

本发明针对现有技术的缺陷，提供了一种实际灌溉面积识别方法。对灌溉引起的田间实际蒸散发进行时段累加，根据区域灌溉蒸散发累加结果的概率分布直方图体现的聚类特征，获得灌溉识别的最短时长窗口。通过对最短时长窗口对应的累计灌溉蒸散发概率分布直方图的聚类特征分析，得到灌溉与未灌溉的判别阈值，使用该阈值识别可判别最短时长窗口期内农田区域的实际灌溉面积与空间分布。类似的，可用本方法判别任何大于最短时长窗口期的农田实际灌溉面积与空间分布。

为了实现以上发明目的，本发明采取的技术方案如下：

一种基于时段累计蒸散发的实际灌溉面积识别方法，包括以下步骤：

步骤1、利用遥感影像反演得到农田范围内像元逐日总蒸散发量。将站点实测日降水量插值匹配到蒸散发像元网格，得到像元网格日降水量，再计算由降水产生的像元逐日蒸散发量即有效降水。从农田逐像元总蒸散发中扣除由降水产生的蒸散发量，得到农田灌溉产生的蒸散发量。

步骤2、自小向大设定不同时长T，分析不同时长窗口下，像元灌溉蒸散发累加结果的概率分布直方图。随着时长T的增大，当灌溉蒸散发量的概率直方图表现为两个聚类时，选定该时长T*为灌溉识别的最短时长窗口。

步骤3、在T*对应的时段累计灌溉蒸散发量概率分布直方图中，选取两个聚类分界点对应的累计灌溉蒸散发量ET*，作为判别灌溉与非灌溉像元的阈值。

步骤4、用选定的阈值ET*对按T*时长累加的灌溉蒸散发数据进行划分，累计灌溉蒸散发量大于ET*的像元为灌溉像元，否则为非灌溉像元。灌溉像元的面积总和为T*时长内的实际灌溉面积。

步骤5、当时长窗口T大于T*时，可重复步骤3，得到T对应的灌溉判别阈值ET*。重复步骤4，得到T时长内的实际灌溉面积。

进一步地，所述步骤1包括以下子步骤：

步骤11、利用卫星遥感影像，反演得到农田范围内的像元逐日总蒸散发量；

步骤12、将农田区域的雨量站实测日降水量进行空间插值，得到与步骤11蒸散发像元网格匹配的网格日降水量，再选取适合当地农田的有效降水计算公式，计算有降水产生的像元逐日蒸散发量；

步骤13、从农田像元逐日总蒸散发量中扣除降水产生的蒸散发量，得到农田像元逐日灌溉产生的蒸散发量。

进一步地，所述步骤2包括以下子步骤：

步骤21、设定时长T为1、2、3至n日，对步骤13得到逐日灌溉蒸散发量，分别进行T时长的累加，对每个T时长对应的累计蒸散发量E_T，统计得到像元累计蒸散发量的概率直方图。

步骤22、对比T为1、2、3至n日的累计蒸散发量概率直方图，当随着T的增大，概率直方图变为两个聚类时，T为灌溉识别的最短时间窗口T*。

进一步地，所述步骤3包括以下子步骤：

步骤31、基于步骤22得到的最短时间窗口T*，选取对应的累计蒸散发概率直方图，用来分析判别像元是否灌溉的阈值；

步骤32、在步骤31选定的概率直方图中，选取两个聚类之间的波谷，波谷处对应的蒸散发量为判别像元是否灌溉的阈值ET*。

进一步地，所述步骤4包括以下子步骤：

步骤41、利用步骤32中选定的阈值ET*，对农田区域按T时长累加的灌溉蒸散发量进行判别，当像元的累计灌溉蒸散发量大于ET*时，判定像元在T时长中有灌溉，否则，无灌溉；

步骤42、对步骤41识别的灌溉像元进行统计分析，可得到该农田区域在T时段内的实际灌溉面积及空间分布。

进一步地，所述步骤5包括以下子步骤：

步骤51、当时长窗口T大于T*时，重复步骤31、32，得到T时长对应的灌溉判别阈值ET*；

步骤52、重复步骤41、42，可得到该农田区域在T时段内的实际灌溉面积及空间分布。

与现有技术相比，本发明的优点在于：

1、利用遥感反演的农田蒸散发来识别像元灌溉情况，可以考虑作物生长的影响，结果受植被覆盖的影响较小，在植被茂密的时期也能保证较好的识别精度；

2、利用蒸散发量的时空可累加性，可结合区域特点选取有效降水的计算方法，通过有效降水的扣除，解决降水影响识别精度的问题。

3、本发明给出了利用累计蒸散发量可识别像元是否灌溉的最短时长窗口的分析方法。

附图说明

图1是本发明实施例区域农田实际灌溉面积识别的流程图；

图2a是本发明实施例区域灌溉蒸散发空间分布计算结果示意图之一；

图2b是本发明实施例区域灌溉蒸散发空间分布计算结果示意图之二；

图3a是本发明实施例区域时段累计5日灌溉蒸散发量概率分布直方图；

图3b是本发明实施例区域时段累计10日灌溉蒸散发量概率分布直方图；

图3c是本发明实施例区域时段累计15日灌溉蒸散发量概率分布直方图；

图3d是本发明实施例区域时段累计20日灌溉蒸散发量概率分布直方图；

图3e是本发明实施例区域时段累计25日灌溉蒸散发量概率分布直方图；

图3f是本发明实施例区域时段累计30日灌溉蒸散发量概率分布直方图；

图4是本发明实施例区域灌溉识别结果空间分布图。

具体实施方式

为使本发明的目的、技术方案及优点更加清楚明白，以下根据附图并列举实施例，对本发明做进一步详细说明。

如图1所示，一种基于时段累计蒸散发的实际灌溉面积识别方法，包括以下步骤：

步骤1、利用遥感影像反演得到农田范围内像元逐日总蒸散发量。将站点实测日降水量插值匹配到蒸散发像元网格，得到像元网格日降水量，再计算由降水产生的像元逐日蒸散发量。从农田逐像元总蒸散发中扣除由降水产生的蒸散发量，得到农田灌溉产生的蒸散发量。

步骤2、自小向大设定不同时长T，分析不同时长窗口下，像元灌溉蒸散发累加结果的概率分布直方图。随着时长T的增大，当灌溉蒸散发量的概率直方图表现为两个聚类时，选定该时长T*为灌溉识别的最短时长窗口。

步骤3、在T*对应的时段累计灌溉蒸散发量概率分布直方图中，选取两个聚类分界点对应的累计灌溉蒸散发量ET*，作为判别灌溉与非灌溉像元的阈值。

步骤4、用选定的阈值ET*对按T*时长累加的灌溉蒸散发数据进行划分，累计灌溉蒸散发量大于ET*的像元为灌溉像元，否则为非灌溉像元。灌溉像元的面积总和为T*时长内的实际灌溉面积。

步骤5、当时长窗口T大于T*时，可以T为累计时长，重复步骤3，得到T对应的灌溉判别阈值ET*。重复步骤4，得到T时长内的实际灌溉面积。

以下是对上述流程的具体分析：

步骤1、利用遥感影像反演得到农田范围内像元逐日总蒸散发量。将站点实测日降水量插值匹配到蒸散发像元网格，得到像元网格日降水量，再计算由降水产生的像元逐日蒸散发量。从农田逐像元总蒸散发中扣除由降水产生的蒸散发量，得到农田灌溉产生的蒸散发量。

步骤11、利用MODIS和Sentinel-2卫星遥感影像联合，以实施例研究区农田范围为掩膜，选用改进的PML蒸散发反演模型计算得到农田范围内的像元逐日总蒸散发量ET_总；

步骤12、将农田区域的雨量站实测日降水量进行空间插值，得到与步骤11蒸散发像元网格匹配的网格日降水量，对实施例研究区选用公式(1)计算计算有降水产生的像元逐日蒸散发量：

其中，ET_C为根据气象数计算得到区域作物需水量；

步骤13、从农田像元逐日总蒸散发量中扣除降水产生的蒸散发量，得到农田像元逐日灌溉产生的蒸散发量ET_灌。示例见图2a和图2b。

ET_灌＝ET_总-P_有效 (2)

步骤2、自小向大设定不同时长T，分析不同时长窗口下，像元灌溉蒸散发累加结果的概率分布直方图。随着时长T的增大，当灌溉蒸散发量的概率直方图表现为两个聚类时，选定该时长T*为灌溉识别的最短时长窗口。具体包括：

步骤21、设定时长T为1、2、3至n日，对步骤13得到逐日灌溉蒸散发量，分别进行T时长的累加，对每个T时长对应的累计灌溉蒸散发量ET_灌，统计得到像元累计蒸散发量的概率直方图。T为5、10、15、20、25、30日的示例见图3a到图3f；

步骤22、对比T为1、2、3至n日的累计ET_灌概率直方图，以图3a到图3f为例，当T增大到20日时，概率直方图变为两个聚类，则20日为实施例农田区域灌溉识别的最短时间窗口。

步骤3、在实施例农田区域20日累计灌溉蒸散发量概率分布直方图中，选取两个聚类分界点对应的累计灌溉蒸散发量ET*，作为判别灌溉与非灌溉像元的阈值。

步骤31、基于步骤22得到的最短时间窗口T*，选取20日累计蒸散发概率直方图，用来分析判别像元是否灌溉的阈值；

步骤32、在步骤31选定的概率直方图中，选取两个聚类之间的波谷，波谷处对应的累计蒸散发量(61mm)为判别像元是否灌溉的阈值ET*。

步骤4、用选定的阈值ET*对按T*时长累加的灌溉蒸散发数据进行划分，累计灌溉蒸散发量大于ET*的像元为灌溉像元，否则为非灌溉像元。灌溉像元的面积总和为T*时长内的实际灌溉面积。

步骤41、以61mm为阈值，对实施例农田区域20日累计灌溉蒸散发量进行判别，当像元的累计灌溉蒸散发量大于61mm时，判定像元在此20日内有灌溉，否则，无灌溉；

步骤42、对步骤41识别的灌溉像元进行统计分析，可得到该实施例农田区域在此20日内的实际灌溉面积及空间分布。

步骤5、当时长窗口T大于T*时，可以T为累计时长，重复步骤3，得到T对应的灌溉判别阈值ET*。重复步骤4，得到T时长内的实际灌溉面积。

步骤51、当时长窗口T大于T*时，重复步骤31、32，得到T时长对应的灌溉判别阈值ET*。以图3为例，T为25日时，ET*为79mm；T为30日时，ET*为105mm；

步骤52、重复步骤41、42，可得到该农田区域在T时段内的实际灌溉面积及空间分布。图4为实施例3月1日-30日内的实际灌溉面积空间分布图。

本领域的普通技术人员将会意识到，这里所述的实施例是为了帮助读者理解本发明的实施方法，应被理解为本发明的保护范围并不局限于这样的特别陈述和实施例。本领域的普通技术人员可以根据本发明公开的这些技术启示做出各种不脱离本发明实质的其它各种具体变形和组合，这些变形和组合仍然在本发明的保护范围内。

Claims (6)

1.一种基于时段累计蒸散发的实际灌溉面积识别方法，其特征在于，包括以下步骤：

步骤1、利用遥感影像反演得到农田范围内像元逐日总蒸散发量；

将站点实测日降水量插值匹配到蒸散发像元网格，得到像元网格日降水量，再计算由降水产生的像元逐日蒸散发量即有效降水；从农田逐像元总蒸散发中扣除由降水产生的蒸散发量，得到农田灌溉产生的蒸散发量；

步骤2、自小向大设定不同时长T，分析不同时长窗口下，像元灌溉蒸散发累加结果的概率分布直方图；随着时长T的增大，当灌溉蒸散发量的概率直方图表现为两个聚类时，选定该时长T*为灌溉识别的最短时长窗口；

步骤3、在T*对应的时段累计灌溉蒸散发量概率分布直方图中，选取两个聚类分界点对应的累计灌溉蒸散发量ET*，作为判别灌溉与非灌溉像元的阈值；

步骤4、用选定的阈值ET*对按T*时长累加的灌溉蒸散发数据进行划分，累计灌溉蒸散发量大于ET*的像元为灌溉像元，否则为非灌溉像元；灌溉像元的面积总和为T*时长内的实际灌溉面积；

步骤5、当时长窗口T大于T*时，重复步骤3，得到T对应的灌溉判别阈值ET*；重复步骤4，得到T时长内的实际灌溉面积。

2.根据权利要求1所述的一种基于时段累计蒸散发的实际灌溉面积识别方法，其特征在于，所述步骤1包括以下子步骤：

步骤11、利用卫星遥感影像，反演得到农田范围内的像元逐日总蒸散发量；

步骤12、将农田区域的雨量站实测日降水量进行空间插值，得到与步骤11蒸散发像元网格匹配的网格日降水量，再选取适合当地农田的有效降水计算公式，计算有降水产生的像元逐日蒸散发量；

步骤13、从农田像元逐日总蒸散发量中扣除降水产生的蒸散发量，得到农田像元逐日灌溉产生的蒸散发量。

3.根据权利要求2所述的一种基于时段累计蒸散发的实际灌溉面积识别方法，其特征在于，所述步骤2包括以下子步骤：

步骤21、设定时长T为1、2、3至n日，对步骤13得到逐日灌溉蒸散发量，分别进行T时长的累加，对每个T时长对应的累计蒸散发量E_T，统计得到像元累计蒸散发量的概率直方图；

步骤22、对比T为1、2、3至n日的累计蒸散发量概率直方图，当随着T的增大，概率直方图变为两个聚类时，T为灌溉识别的最短时间窗口T*。

4.根据权利要求3所述的一种基于时段累计蒸散发的实际灌溉面积识别方法，其特征在于，所述步骤3包括以下子步骤：

步骤31、基于步骤22得到的最短时间窗口T*，选取对应的累计蒸散发概率直方图，用来分析判别像元是否灌溉的阈值；

步骤32、在步骤31选定的概率直方图中，选取两个聚类之间的波谷，波谷处对应的蒸散发量为判别像元是否灌溉的阈值ET*。

5.根据权利要求4所述的一种基于时段累计蒸散发的实际灌溉面积识别方法，其特征在于，所述步骤4包括以下子步骤：

步骤41、利用步骤32中选定的阈值ET*，对农田区域按T时长累加的灌溉蒸散发量进行判别，当像元的累计灌溉蒸散发量大于ET*时，判定像元在T时长中有灌溉，否则，无灌溉；

步骤42、对步骤41识别的灌溉像元进行统计分析，得到该农田区域在T时段内的实际灌溉面积及空间分布。

6.根据权利要求5所述的一种基于时段累计蒸散发的实际灌溉面积识别方法，其特征在于，所述步骤5包括以下子步骤：

步骤51、当时长窗口T大于T*时，重复步骤31、32，得到T时长对应的灌溉判别阈值ET*；

步骤52、重复步骤41、42，得到该农田区域在T时段内的实际灌溉面积及空间分布。

Images