CN113128401B - 基于光学和雷达遥感数据的区域实际灌溉面积监测方法 - Google Patents

基于光学和雷达遥感数据的区域实际灌溉面积监测方法 Download PDF

Info

Publication number
CN113128401B
CN113128401B CN202110418925.8A CN202110418925A CN113128401B CN 113128401 B CN113128401 B CN 113128401B CN 202110418925 A CN202110418925 A CN 202110418925A CN 113128401 B CN113128401 B CN 113128401B
Authority
CN
China
Prior art keywords
irrigation
remote sensing
image
radar
area
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Expired - Fee Related
Application number
CN202110418925.8A
Other languages
English (en)
Other versions
CN113128401A (zh
Inventor
杨永民
李璐
顾涛
吴迪
徐磊
杨鹏
黄诗峰
杨昆
江威
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
China Institute of Water Resources and Hydropower Research
Original Assignee
China Institute of Water Resources and Hydropower Research
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by China Institute of Water Resources and Hydropower Research filed Critical China Institute of Water Resources and Hydropower Research
Priority to CN202110418925.8A priority Critical patent/CN113128401B/zh
Publication of CN113128401A publication Critical patent/CN113128401A/zh
Application granted granted Critical
Publication of CN113128401B publication Critical patent/CN113128401B/zh
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06VIMAGE OR VIDEO RECOGNITION OR UNDERSTANDING
    • G06V20/00Scenes; Scene-specific elements
    • G06V20/10Terrestrial scenes
    • G06V20/188Vegetation
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00Pattern recognition
    • G06F18/20Analysing
    • G06F18/21Design or setup of recognition systems or techniques; Extraction of features in feature space; Blind source separation
    • G06F18/213Feature extraction, e.g. by transforming the feature space; Summarisation; Mappings, e.g. subspace methods
    • G06F18/2135Feature extraction, e.g. by transforming the feature space; Summarisation; Mappings, e.g. subspace methods based on approximation criteria, e.g. principal component analysis
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00Pattern recognition
    • G06F18/20Analysing
    • G06F18/21Design or setup of recognition systems or techniques; Extraction of features in feature space; Blind source separation
    • G06F18/214Generating training patterns; Bootstrap methods, e.g. bagging or boosting
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00Pattern recognition
    • G06F18/20Analysing
    • G06F18/24Classification techniques
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06NCOMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
    • G06N20/00Machine learning
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06QINFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
    • G06Q50/00Information and communication technology [ICT] specially adapted for implementation of business processes of specific business sectors, e.g. utilities or tourism
    • G06Q50/02Agriculture; Fishing; Forestry; Mining
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/10Segmentation; Edge detection
    • G06T7/136Segmentation; Edge detection involving thresholding
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00Image analysis
    • G06T7/60Analysis of geometric attributes
    • G06T7/62Analysis of geometric attributes of area, perimeter, diameter or volume
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06VIMAGE OR VIDEO RECOGNITION OR UNDERSTANDING
    • G06V10/00Arrangements for image or video recognition or understanding
    • G06V10/20Image preprocessing
    • G06V10/26Segmentation of patterns in the image field; Cutting or merging of image elements to establish the pattern region, e.g. clustering-based techniques; Detection of occlusion
    • G06V10/267Segmentation of patterns in the image field; Cutting or merging of image elements to establish the pattern region, e.g. clustering-based techniques; Detection of occlusion by performing operations on regions, e.g. growing, shrinking or watersheds
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10032Satellite or aerial image; Remote sensing
    • G06T2207/10036Multispectral image; Hyperspectral image
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/10Image acquisition modality
    • G06T2207/10032Satellite or aerial image; Remote sensing
    • G06T2207/10044Radar image
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06TIMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00Indexing scheme for image analysis or image enhancement
    • G06T2207/30Subject of image; Context of image processing
    • G06T2207/30181Earth observation
    • G06T2207/30188Vegetation; Agriculture

Landscapes

  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Data Mining & Analysis (AREA)
  • Evolutionary Computation (AREA)
  • General Engineering & Computer Science (AREA)
  • Artificial Intelligence (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Evolutionary Biology (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Business, Economics & Management (AREA)
  • Health & Medical Sciences (AREA)
  • General Health & Medical Sciences (AREA)
  • Multimedia (AREA)
  • Software Systems (AREA)
  • Medical Informatics (AREA)
  • Mathematical Physics (AREA)
  • Geometry (AREA)
  • Agronomy & Crop Science (AREA)
  • Animal Husbandry (AREA)
  • Marine Sciences & Fisheries (AREA)
  • Mining & Mineral Resources (AREA)
  • Computing Systems (AREA)
  • Economics (AREA)
  • Human Resources & Organizations (AREA)
  • Marketing (AREA)
  • Primary Health Care (AREA)
  • Strategic Management (AREA)
  • Tourism & Hospitality (AREA)
  • General Business, Economics & Management (AREA)
  • Image Processing (AREA)
  • Investigating Or Analysing Materials By Optical Means (AREA)

Abstract

本发明提出一种基于光学和雷达遥感数据的区域实际灌溉面积监测方法,将对土壤和植被含水量变化敏感的短波红外指数和雷达后向散射系数变化结合,基于卫星遥感数据,可对区域灌溉活动进行动态探测,实现灌溉面积的准确提取。

Description

基于光学和雷达遥感数据的区域实际灌溉面积监测方法
技术领域
本发明涉及农业水资源监测管理技术领域,为一种基于光学和雷达时序遥感数据的区域实际灌溉面积动态监测方法,可服务于灌区水资源管理、农业用水效率评估、农业节水、区域人类取用水活动影响评估、灌区现代化管理等方面的应用领域。
背景技术
区域实际灌溉面积的动态监测是农业水资源管理、农业用水效率评估、农业节水、区域人类活动影响评估、灌区现代化管理、流域水资源管理等应用领域的基础。农业是我国第一用水大户,农业用水占用水总量的60%以上,在北方干旱/半干旱地区占比达70%以上。目前,我国灌溉面积达11.1亿亩,居世界第一,其中耕地灌溉面积10.2亿亩,占全国耕地总面积的50.3%,科学准确的灌溉面积调查数据是实施水资源管理的重要基础。人类活动的水文效应是国内外关注的焦点之一,人类通过取用水活动直接或间接地影响着陆地水循环,人类活动特别是灌溉用水是影响陆面水文循环最难刻画的部分,大部分的水文模型对灌溉事件的刻画仍然采用预设的土壤水分消退阈值法或灌溉制度进行参数化,存在较大的不确定性。此外,目前的公布的区域灌溉面积资料存在较大的不确定性,且难以反映灌溉时间和灌溉区域的分布,因此,准确的提取区域实际灌溉面积是进一步深入评估人类取用水活动对区域陆地水循环影响研究的前提。区域实际灌溉面积的动态监测结果能够为水利、农业、林业等专业技术人员提供重要下垫面状况信息,为辅助决策提供关键支持信息。
遥感作为地面观测的重要手段,为区域实际灌溉面积的动态监测提供了重要的技术手段。现有的基于卫星遥感的灌溉面积调查方法主要包括以下几类:(1)植被指数阈值法,目前大部分的研究工作主要集中在可见光和近红外谱段范围,如垂直植被指数方法。此类方法不确定性大,阈值确定较难。(2)基于遥感地表温度方法,这类方法基于地表温度和下垫面土壤水分之间的关系来探测灌溉活动。此类方法一方面受制于目前热红外卫星遥感的分辨率,目前民用热红外最高的分辨率在60米左右。此外,优于100米热红外遥感的遥感数据的重访周期较长。(3)基于雷达的方法,这类方法一种是结合雷达后向散射系数的灌溉事件探测,另外一种是结合雷达反演土壤水分进行探测。此类方法目前收到雷达后向散射系数变化的多重因素影响,仍然有待进一步研究。(4)结合地表蒸散和土壤水分遥感反演的土壤水分探测方法,此类方法需要结合SMAP等土壤水分遥感反演结果,空间分辨率较粗。
区域实际灌溉面积的动态监测对区域水资源管理等方面的应用至关重要。目前,传统的基于光学植被指数阈值法和基于遥感地表温度的灌溉面积提取方法对灌溉活动导致的土壤和植被含水量的变化响应较差,往往提取的灌溉范围及面积存在较大不确定性。
发明内容
针对目前基于遥感的灌区大范围灌溉面积提取方法存在的严重不足,本发明提出一种基于光学和雷达遥感数据的区域实际灌溉面积监测方法。
本发明的目的是这样实现的:
一种基于光学和雷达遥感数据的区域实际灌溉面积监测方法,包括以下步骤:
步骤1研究区基础资料收集:收集研究区域内农作物(主要种植农作物及主要农作物)的物候期及农事历信息;所述农事历信息包括农作物的常年灌溉时间;
步骤2研究区卫星遥感数据收集:收集研究区域光学卫星遥感数据和雷达卫星遥感数据,光学卫星遥感数据包含短波红外波波段,所述短波红外波谱范围为1.3-3.0μm;
步骤3卫星遥感数据预处理:光学遥感数据的预处理包括辐射定标和大气校正,光学遥感影像数据处理主要包括几何校正、正射校正、图像增强、影像融合、影像镶嵌等处理;雷达遥感数据的预处理包括辐射定标、多视处理、滤波、正射校正;
步骤4基于面向对象分类法的研究区农作物种植结构遥感提取:利用预处理后的光学卫星遥感数据计算光谱指数,采集农作物的训练样本数据,基于面向对象分类方法实现研究区作物种植结构的遥感提取;
步骤5短波红外指数差异图生成:基于预处理后的光学遥感影像数据,计算短波红外指数,根据监测时间范围计算短波红外指数差异图,采用均值比算子运算;
步骤6基于光学的灌溉探测方法:基于短波红外指数差异图,使用自动阈值分割方法实现区域灌溉识别;
步骤7雷达后向散射系数差异图生成:基于预处理后的雷达遥感影像数据,计算后向散射系数,根据检测时间需求,选择两个时相的遥感数据采用均值比算子计算差异图;
步骤8基于雷达的灌溉探测方法:基于雷达后向散射系数差异图,使用自动阈值分割方法实现区域灌溉识别;
步骤9光学和雷达灌溉探测结果判别与融合:对步骤6和步骤8的灌溉初步探测结果进行判别与融合计算;
步骤10实际灌溉面积提取与制图:对步骤9计算的灌溉探测结果进行面积统计计算,并进行空间制图表达。
进一步的方案,步骤2中光学卫星选择Landsat8或Sentinel 2,雷达卫星选择Sentinel1;Landsat 8和Sentinel1组合满足分辨率为30米的监测所需,Sentinel2和Sentinel1组合满足20分辨率监测所需。
进一步的,步骤4的具体操作为:
4-1结合光学卫星遥感数据计算遥感光谱指数:光谱指数包括:归一化植被指数(NDVI)、归一化水体指数(NDWI)、改进型归一化差值水体指数(MNDWI)、裸土指数(BSI);
4-2基于主成分分析的图像变换:选择三个主成分组成多波段图像,实现图像数据的降维;
4-3波段组合及训练样本采集:使用步骤4-2中生成的多波段图像结合步骤4-1中的光谱指数计算结果组成待分类的多波段图像;结合步骤1收集的研究区基础资料,选择研究区主要农作物的训练样本(为保证遥感图像分类的准确性)每类作物采集样本数目大于12个(且尽可能保证分布均匀),基于采样点位置数据,提取对应点位的待分类的多波段图像信息,形成训练样本集合;
4-4基于超像素分割算法实现图像的分割计算:基于超像素的分割将图像分割成一系列子区域,每个子区域内部具有波段特征的一致性;采用SNIC超像素算法进行图像分割计算;
4-5机器学习图像分类算法训练:在4-3训练样本采集的基础上,进行机器学习分类模型方法的训练,采用随机森林机器学习方法,对采集的训练样本数据进行训练,生成图像分类器;
4-6面向超像素的图像分类及作物种植结构提取:在图像的分割计算和图像分类算法训练的基础上,使用训练的随机森林图像分类器对图像分割的超像素图像数据进行分类计算;在图像分类的基础上,提取出研究区的主要作物类型的分布,并统计各类作物的面积。
进一步的,步骤5的具体操作为:
5-1、对步骤3预处理后的光学遥感影像数据计算短波红外指数,公式如下:
Figure BDA0003027060420000031
Figure BDA0003027060420000032
gvmimax=max(gvmimax,seasonal,gvmiwater)
其中,gvmi为全球植被水分指数,gvmimax为全球植被水分指数的最大值;NIR为近红外波段的反射率;SWIR为短波红外的反射率;gvmimax,seasonal为像元点gvmi的最大值,gvmiwater为距离像元点最近的水体像元的gvmi值;
5-2、在短波红外指数计算的基础上,进一步计算监测时段的短波红外指数时序差异图:结合超像素分割算法的区域图像分割计算结果,计算每个分割子区域监测时段起始两期图像的WSCI均值,对每个子区域计算监测起始WSCI的差值,生成短波红外指数时序差异图:
ΔWSCI=WSCIt2-WSCIt1
其中,ΔWSCI为监测时段WSCI的变量值,WSCIt2和WSCIt1分别为监测时段末和监测时段初的WSCI值。
进一步的,步骤6的具体操作为:
6-1.将按照研究区的范围生成的一系列规则格网,基于生成的规则格网将短波红外指数时序差异图分割为一系列的子图像;
6-2.对于每个分割的子图像计算其像元的平均值,以此像元平均值作为分割阈值实现子图像的潜在灌溉探测识别;
6-3.将所有子图像的潜在灌溉探测识别结果进行合并,生成最终的潜在灌溉光学遥感探测结果。
进一步的,步骤7的具体操作为:
7-1对步骤3预处理后的雷达卫星遥感数据,计算后向散射系数;
7-2基于超像素分割算法的区域图像分割计算结果,计算每个分割子区域监测时段起始两期图像的后向散射系数均值,对每个子区域计算监测起始后向散射系数的差值,生成后向散射系数时序差异图:
ΔVV=VVt2-VVt1
ΔVH=VHt2-VHt1
其中,ΔVV和ΔVH为监测时段VV极化和VH极化散射系数的变量值,VVt2和VVt1分别为监测时段末和监测时段初的VV极化的后向散射系数值,VHt2和VHt1分别为监测时段末和监测时段初的VH极化的后向散射系数值。
进一步的,步骤8的具体操作为:
8-1、对VV和VH的后向散射系数时序差异图进行合并处理:
ΔVHVV=ΔVH+ΔVV
8-2、将按照研究区的范围生成一系列规则格网,基于生成的规则格网将VV和VH雷达后向散射系数时序差异图分割为一系列的子图像;
8-3、对于每个分割的子图像计算其像元的平均值,以此像元平均值作为分割阈值实现子图像的潜在灌溉探测识别;
8-4、最后将所有子图像的潜在灌溉探测识别结果进行合并,生成最终的潜在灌溉雷达遥感探测结果。
进一步的,步骤9的具体操作为:
基于超像素分割算法的区域图像分割计算结果,对每个分割的区域基于光学和雷达的灌溉探测结果进行判别:
如果两者的探测结果都为潜在灌溉,则此区域标记为灌溉;对于两者探测结果存在分歧的超像素区域,结合步骤4研究区作物种植结构提取结果对存在分歧的超像素区域进行判别:
如果存在分歧超像素区域的不是农作物类型,则剔除;如果存在分歧的超像素区域属于农作物类型,结合步骤1研究区农作物灌溉时间进行判别:
如果在监测时段内存在分歧的超像素区域属于农作物类型为灌溉期,则判定为灌溉,相反,则判定为非灌溉;
基于以上判定,对所有判定为灌溉的区域进行合并,生成最终的灌溉光学和雷达遥感探测结果。
本发明的优点和有益效果是:
本发明提出一种基于光学和雷达遥感数据的区域实际灌溉面积监测方法,特别涉及一种对土壤和植被含水量变化敏感的短波红外指数和雷达后向散射系数变化的灌溉探测方法,结合卫星星座观测数据(重访周期优于15天的),可对灌溉活动调研影响进行动态探测,实现灌溉面积的准确提取。
本发明提出的模型方法通过有效结合高频次光学和雷达卫星星座观测,可实现对短期灌溉活动导致的土壤水分和植被含水量变化特征的准确探测,有效的杜绝了传统灌溉探测方法的不确定性。本发明提出的方法可促进基于卫星遥感的灌溉探测,实现区域实际灌溉面积的提取,具有重要的应用价值。本发明所提出的模型方法可服务于农业灌区水资源管理、农业用水效率评估、农业节水、区域人类取用水活动影响评估、灌区现代化管理等方面的应用领域。
附图说明
下面结合附图和实施例对本发明作进一步说明。
图1为本发明基于光学和雷达遥感数据的区域实际灌溉面积监测方法的流程图。
图2为本发明应用于河北磁右灌区探测案例。
具体实施方式
实施例1:
(一)研究区基础资料收集
根据研究所需,确定研究区的地理位置及边界。收集区域内主要种植农作物及主要农作物的物候期,并收集区域主要农作物农事历信息。调查研究区域主要农作物农事历信息包括主要农作物的常年的灌溉时间。本实施例将基于此调查信息制作研究区主要农作物常年灌溉时间的查找表。基于此查找表可实现对应作物主要灌溉时间的快速查询。
(二)研究区光学和雷达卫星遥感数据收集
根据研究区域需要收集研究区光学和雷达卫星遥感数据收集。本发明提及的方法使用对土壤和植被含水量非常敏感的短波红外信息做支撑,因此光学卫星数据需包含短波红外波段。短波红外波谱范围为1.3-3.0μm。光学卫星可选择Landsat8或Sentinel 2,雷达卫星选择Sentinel 1。Landsat 8和Sentinel1组合满足分辨率为30米的灌溉探测所需,Sentinel2和Sentinel1组合满足20分辨率灌溉探测所需。本实施例将以Sentinel2和Sentinel1组合为例。
(三)光学和雷达卫星遥感数据预处理
本发明的技术方法所需的遥感数据预处理包括光学遥感数据预处理和雷达遥感数据预处理。光学遥感数据预处理包括辐射定标和大气校正。辐射定标是将传感器接收到的遥感数据,转换为辐射亮度或天顶反射率,根据所获取数据的传感器类型和相应的定标系数,对光学遥感数据进行辐射定标计算。大气校正处理旨在于消除大气的吸收、散射等效应对地表反射率的影响,消除由大气影响所造成的辐射误差,计算地物的反射率。采用ENVI/FLAASH对ASTER传感器采集的数据进行大气校正处理。依据ENVI/FLAASH输入要求和卫星遥感数据获取的头文件信息对遥感数据进行大气校正处理,获得地表反射率数据。本实施例采用Sentinel2遥感数据,可选择ESA提供的共享网站获得相应的地表反射率数据。
雷达遥感数据的预处理主要包括斑点滤波、辐射定标、多视处理、地理编码、后向散射系数提、影像裁剪。经过预处理计算的区域的雷达后向散射系数。Sentinel1的预处理可使用SNAP软件进行操作。对于Sentinel1的预处理操作流程可参考[http://eoscience.esa.int/landtraining2017/files/materials/D2P1__I.pdf]
(四)基于面向对象分类方法的研究区作物种植结构提取
在光学遥感数据预处理的基础上,本步骤将结合光学遥感数据实现研究区作物种植结构的提取。
①结合光学卫星遥感数据计算遥感光谱指数。
光谱指数包括:归一化植被指数(NDVI)、归一化水体指数(NDWI)、改进型归一化差值水体指数(MNDWI)、裸土指数(BSI)。
基于Sentinel2卫星遥感数据的归一化植被指数(NDVI)计算如下:
NDVI=(Band8-Band4)/(Band8+Band4)
其中,Band4和Band8为Sentinel2第4波段和第8波段的反射率值。
归一化水体指数(NDWI)计算如下:
NDWI=(Band8-Band12)/(Band8+Band12)
其中,Band8和Band12为Sentinel2第8波段和第12波段的反射率值。
改进型归一化差值水体指数计算如下:
MNDWI=(Band3-Band12)/(Band3+Band12)
其中,Band3和Band12为Sentinel2第3波段和第12波段的反射率值。
裸土指数(BSI)计算如下:
BSI=((Band12+Band4)-(Band8+Band2))/((Band12+Band4)-(Band8+Band2))
其中,Band2,Band4,Band8,Band12为Sentinel2第2、4、8和第12波段的反射率值。
②基于主成分分析的图像变换。
由于Sentinel2的图像波段数目较多,各波段间具有较强的相关性,因此使用主成分分析(PCA)方法对Sentinel2光谱数据进行图像变化操作,选择第1、2、3主成分组成三波段图像,实现图像数据的降维。
③波段组合及训练样本采集。
使用②主成分分析生成的多波段图像结合①中遥感光谱指数计算结果组成待分类的多波段图像。结合步骤(一)收集的研究区基础资料,选择研究区主要农作物的训练样本,为保证遥感图像分类的准确性,每类作物采集样本数目大于12个,且尽可能保证分布均匀。基于采样点位置数据,提取对应点位的待分类图像的多波段信息,形成训练样本集合。
④基于超像素分割算法实现图像的分割计算。
基于超像素的分割可将图像分割成一系列子区域,每个子区域内部具有波段特征的一致性。本发明采用简单非迭代聚类(Simple Non-Iterative Clustering,SNIC)超像素算法开展图像的分割计算。SNIC分割算法能够实现多波段图像的分割,并大幅缩短运算时间,具有较强的可操作性。关于SNIC分割算法的介绍可参考[Achanta R,SusstrunkS.Superpixels and Polygons Using Simple Non-iterative Clustering[C]//IEEEConference on Computer Vision&Pattern Recognition.IEEE,2017.]。基于超像素分割算法实现图像的分割计算,保存图像分割计算数据,用于后续的面向对象的图像分类及区域灌溉的探测。
⑤机器学习图像分类算法训练。
在③训练样本采集的基础上,进行机器学习分类模型方法的训练。本实施例使用python的sklearn机器学习模型库下使用其随机森林机器学习方法,对采集的训练样本数据进行训练,生成图像分类器,用于下一步的图像分类。
⑥面向超像素的图像分类及作物种植结构提取。
在图像的分割计算和图像分类算法训练的基础上,使用训练的随机森林图像分类器对图像分割的超像素图像数据进行分类计算。在图像分类的基础上,提取出研究区的主要作物类型的分布,并统计各类作物的面积。
(五)短波红外指数时序差异图生成
短波红外波段(1.3-3.0μm)对土壤水分和植被含水量较为敏感,基于短期内(10天内)的短波红外指数的变化对灌溉信息进行探测识别。
首先,对步骤(三)处理后的光学遥感影像数据计算短波红外指数,本实施例采用Ceccato等(2002)提出全球植被水分指数(global vegetation moisture index,gvmi)构建刻画下垫面土壤水分供给状况的短波红外指数(WSCI,water supply conditionindex):
Figure BDA0003027060420000081
Figure BDA0003027060420000082
gvmimax=max(gvmimax,seasonal,gvmiwater)
其中,gvmi为全球植被水分指数,gvmimax为全球植被水分指数的最大值;NIR为近红外波段的反射率,对应于Sentinel2的第8波段;SWIR为短波红外的反射率,对应于Sentinel2的第12波段;gvmimax,seasonal为像元点gvmi的最大值,gvmiwater为距离像元点最近的水体像元的gvmi值。
然后,在短波红外指数计算的基础上,进一步计算监测时段的短波红外指数时序差异图。由于基于像元的图像差值计算往往导致计算结果的随机噪声较大,无法满足灌溉信息探测的需要。对此,本实施例结合超像素分割计算结果,实现超像素的时序差值计算。结合步骤(四)④基于超像素分割算法的区域图像分割计算结果,计算每个分割子区域监测时段起始和结束两期图像的WSCI均值,对每个子区域计算监测起始WSCI的差值,生成短波红外指数时序差异图:
ΔWSCI=WSCIt2-WSCIt1
其中,ΔWSCI为监测时段WSCI的变量值,WSCIt2和WSCIt1分别为监测时段末和监测时段初的WSCI值。通过结合超像素分割的短波红外指数差值计算一方面可实现像元点随机噪声的抑制,另一方面实现潜在灌溉地块信息的增强。
(六)基于光学遥感数据的潜在灌溉探测
基于步骤(五)生成的短波红外指数时序差异图,本实施例使用局部阈值分割方法开展基于光学遥感数据的潜在灌溉探测。由于受区域内灌溉先后顺序及作物类型的影响,基于固定的全局阈值则不能有效的兼顾区域各处的情况,将会对灌溉探测的结果产生影响。对此,本实施采用局部阈值分割法进行灌溉地块的探测。首先将按照研究区的范围生成6km×6km(此数值可以根据区域情况具体调整)的一系列规则格网,基于生成的规则格网将短波红外指数时序差异图分割为一系列的子图像。然后,对于每个分割的子图像计算其像元的平均值,以此像元平均值作为分割阈值实现子图像的潜在灌溉探测识别。最后将所有子图像的潜在灌溉探测识别结果进行合并,生成最终的潜在灌溉光学遥感探测结果。
(七)雷达后向散射系数时序差异图生成
雷达后向散射系数对土壤和植被的含水量敏感,基于Sentinel1的雷达后向散射系数常被应用于区域土壤水分的估计。由于灌溉活动最终会导致地块的土壤水分和植被含水量的增加,因此可以基于雷达后向散射系数时序差异图对灌溉信息进行探测。首先,对步骤(三)处理后的雷达遥感后向散射系数,使用和步骤(五)相似的处理方法,基于超像素算法计算后向散射系数的时序差异图。结合步骤(四)④基于超像素分割算法的区域图像分割计算结果,计算每个分割子区域监测时段起始两期图像的后向散射系数均值,对每个子区域计算监测起始后向散射系数的差值,生成后向散射系数时序差异图:
ΔVV=VVt2-VVt1
ΔVH=VHt2-VHt1
其中,ΔVV和ΔVH为监测时段VV极化和VH极化散射系数的变量值,VVt2和VVt1分别为监测时段末和监测时段初的VV极化的后向散射系数值,VHt2和VHt1分别为监测时段末和监测时段初的VH极化的后向散射系数值。结合超像素分割的雷达后向散射系数差值计算可极大的实现雷达随机噪声的抑制,并可实现潜在灌溉信息的增强。
(八)基于雷达遥感数据的灌溉探测
基于步骤(七)生成的VV和VH雷达后向散射系数时序差异图,使用局部阈值分割法开展基于雷达数据的潜在灌溉探测。由于雷达遥感数据包括VV和VH的灌溉探测结果,本发明将对两者的后向散射系数时序差异图进行合并处理:
ΔVHVV=ΔVH+ΔVV
与步骤(六)类似,本发明使用局部阈值分割方法开展基于雷达遥感数据的灌溉探测。对此,本发明采用局部阈值分割法进行灌溉地块的探测。首先将按照研究区的范围生成6km×6km(此数值可以根据区域情况具体调整)的一系列规则格网,基于生成的规则格网将VV和VH雷达后向散射系数时序差异图分割为一系列的子图像。然后,对于每个分割的子图像计算其像元的平均值,以此像元平均值作为分割阈值实现子图像的潜在灌溉探测识别。最后将所有子图像的潜在灌溉探测识别结果进行合并,生成最终的潜在灌溉雷达遥感探测结果。
(九)基于光学和雷达遥感数据的灌溉探测结果判别与合并
本步骤将在步骤(六)和(八)开展基于光学和雷达遥感数据的潜在灌溉探测结果判别与融合。首先,基于步骤(四)④基于超像素分割算法的区域图像分割计算结果,对每个分割的区域基于光学和雷达的灌溉探测结果进行判别,如果两者的探测结果都为潜在灌溉,则此区域标记为灌溉。对于两者探测结果存在分歧的超像素区域,本发明将结合土地利用和植被指数进行进一步的判别。结合步骤(四)研究区作物种植结构提取结果对存在的分歧的超像素区域进行判别,如果存在分歧超像素区域的不是农作物类型,则剔除。如果存在分歧的超像素区域属于农作物类型,结合步骤(一)研究区主要农作物常年灌溉时间的查找表进行判别,如果在监测时段内存在分歧的超像素区域属于农作物类型为灌溉期,则判定为灌溉,相反,则判定为非灌溉。基于以上判定,对所有判定为灌溉的区域进行合并,生成最终的灌溉光学和雷达遥感探测结果。
(十)区域实际灌溉面积的统计与制图
对步骤(九)计算的灌溉探测结果进行面积统计计算,结合步骤(四)基于面向对象分类方法的研究区作物种植结构提取结果,进一步对研究区各主要作物种植面积进行统计。最后,对基于遥感光学和雷达的区域实际灌溉探测结果,结合区域基础地理信息数据实现区域的实际灌溉面积的制图输出。本发明提出的一种基于光学和雷达时序遥感数据的区域实际灌溉面积动态监测方法具体步骤见图1.
为验证本发明提出方法的可靠性,将本发明应用于河北省磁右灌区灌溉面积的探测之中,图2为河北磁右灌区2021年3月底至四月初灌溉范围,灌溉面积为13.7万亩。经实地调研检验,本发明探测的结果准确,与当地的实际灌溉状况基本一致。
最后应说明的是,以上仅用以说明本发明的技术方案而非限制,尽管参照较佳布置方案对本发明进行了详细说明,本领域的普通技术人员应当理解,可以对本发明的技术方案进行修改或者等同替换,而不脱离本发明技术方案的精神和范围。

Claims (7)

1.一种基于光学和雷达遥感数据的区域实际灌溉面积监测方法,其特征在于:包括以下步骤:
步骤1研究区基础资料收集:收集研究区域内农作物的物候期及农事历信息;所述农事历信息包括农作物的常年灌溉时间;
步骤2研究区卫星遥感数据收集:收集研究区域光学卫星遥感数据和雷达卫星遥感数据,光学卫星遥感数据包含短波红外波波段,所述短波红外波谱范围为1.3-3.0μm;
步骤3卫星遥感数据预处理:光学遥感数据的预处理包括辐射定标和大气校正;雷达遥感数据的预处理包括辐射定标、多视处理、滤波、正射校正;
步骤4基于面向对象分类法的研究区作物种植结构遥感提取:利用预处理后的光学卫星遥感数据计算光谱指数,采集农作物的训练样本数据,基于面向对象分类方法实现研究区作物种植结构的遥感提取;
步骤5短波红外指数时序差异图生成:基于预处理后的光学遥感影像数据,计算短波红外指数,根据监测时间范围计算短波红外指数时序差异图,采用均值比算子运算;步骤5的具体操作为:
5-1、对步骤3预处理后的光学遥感影像数据计算短波红外指数,公式如下:
Figure FDA0003295335270000011
Figure FDA0003295335270000012
gvmimax=max(gvmimax,seasonal,gvmiwater)
其中,gvmi为全球植被水分指数,gvmimax为全球植被水分指数的最大值;NIR为近红外波段的反射率;SWIR为短波红外的反射率;gvmimax,seasonal为像元点gvmi的最大值,gvmiwater为距离像元点最近的水体像元的gvmi值;
5-2、在短波红外指数计算的基础上,进一步计算监测时段的短波红外指数时序差异图:结合超像素分割算法的区域图像分割计算结果,计算每个分割子区域监测时段起始两期图像的WSCI均值,对每个子区域计算监测起始WSCI的差值,生成短波红外指数时序差异图:
ΔWSCI=WSCIt2-WSCIt1
其中,ΔWSCI为监测时段WSCI的变量值,WSCIt2和WSCIt1分别为监测时段末和监测时段初的WSCI值;
步骤6基于光学的灌溉探测方法:基于短波红外指数时序差异图,使用自动阈值分割方法实现区域灌溉识别;
步骤7雷达后向散射系数时序差异图生成:基于预处理后的雷达遥感影像数据,计算后向散射系数,根据监测时间需求,选择两个时相的遥感数据采用均值比算子计算时序差异图;
步骤8基于雷达的灌溉探测方法:基于雷达后向散射系数时序差异图,使用自动阈值分割方法实现区域灌溉识别;
步骤9光学和雷达灌溉探测结果判别与融合:对步骤6和步骤8的灌溉初步探测结果进行判别与融合计算;
步骤10实际灌溉面积提取与制图:对步骤9计算的灌溉探测结果进行面积统计计算,并进行空间制图表达。
2.根据权利要求1所述的基于光学和雷达遥感数据的区域实际灌溉面积监测方法,其特征在于:步骤2中光学卫星选择Landsat8或Sentinel 2,雷达卫星选择Sentinel 1;Landsat 8和Sentinel1组合满足分辨率为30米的监测所需,Sentinel2和Sentinell组合满足20米分辨率监测所需。
3.根据权利要求1所述的基于光学和雷达遥感数据的区域实际灌溉面积监测方法,其特征在于:步骤4的具体操作为:
4-1结合光学卫星遥感数据计算遥感光谱指数:光谱指数包括:归一化植被指数、归一化水体指数、改进型归一化差值水体指数、裸土指数;
4-2基于主成分分析的图像变换:选择三个主成分组成多波段图像,实现图像数据的降维;
4-3波段组合及训练样本采集:使用步骤4-2中生成的三波段图像结合步骤4-1中的光谱指数计算结果组成待分类的多波段图像;结合步骤1收集的研究区基础资料,选择研究区主要农作物的训练样本每类作物采集样本数目大于12个,基于采样点位置数据,提取对应点位的待分类的多波段图像信息,形成训练样本集合;
4-4基于超像素分割算法实现图像的分割计算:基于超像素的分割将图像分割成一系列子区域,每个子区域内部具有波段特征的一致性;采用SNIC超像素算法进行图像分割计算;
4-5机器学习图像分类算法训练:在4-3训练样本采集的基础上,进行机器学习分类模型方法的训练,采用随机森林机器学习方法,对采集的训练样本数据进行训练,生成图像分类器;
4-6面向超像素的图像分类及作物种植结构提取:在图像的分割计算和图像分类算法训练的基础上,使用训练的随机森林图像分类器对图像分割的超像素图像数据进行分类计算;在图像分类的基础上,提取出研究区的主要作物类型的分布,并统计各类作物的面积。
4.根据权利要求1所述的基于光学和雷达遥感数据的区域实际灌溉面积监测方法,其特征在于:步骤6的具体操作为:
6-1.将按照研究区的范围生成的一系列规则格网,基于生成的规则格网将短波红外指数时序差异图分割为一系列的子图像;
6-2.对于每个分割的子图像计算其像元的平均值,以此像元平均值作为分割阈值实现子图像的潜在灌溉探测识别;
6-3.将所有子图像的潜在灌溉探测识别结果进行合并,生成最终的潜在灌溉光学遥感探测结果。
5.根据权利要求1所述的基于光学和雷达遥感数据的区域实际灌溉面积监测方法,其特征在于:步骤7的具体操作为:
7-1对步骤3预处理后的雷达卫星遥感数据,计算后向散射系数;
7-2基于超像素分割算法的区域图像分割计算结果,计算每个分割子区域监测时段起始两期图像的后向散射系数均值,对每个子区域计算监测起始后向散射系数的差值,生成后向散射系数时序差异图:
ΔVV=VVt2-VVt1
ΔVH=VHt2-VHt1
其中,ΔVV和ΔVH为监测时段VV极化和VH极化散射系数的变量值,VVt2和VVt1分别为监测时段末和监测时段初的VV极化的后向散射系数值,VHt2和VHt1分别为监测时段末和监测时段初的VH极化的后向散射系数值。
6.根据权利要求5所述的基于光学和雷达遥感数据的区域实际灌溉面积监测方法,其特征在于:步骤8的具体操作为:
8-1、对VV和VH的后向散射系数时序差异图进行合并处理,合并后的值ΔVHVV为:
ΔVHVV=ΔVH+ΔVV
8-2、将按照研究区的范围生成一系列规则格网,基于生成的规则格网将VV和VH雷达后向散射系数时序差异图分割为一系列的子图像;
8-3、对于每个分割的子图像计算其像元的平均值,以此像元平均值作为分割阈值实现子图像的潜在灌溉探测识别;
8-4、最后将所有子图像的潜在灌溉探测识别结果进行合并,生成最终的潜在灌溉雷达遥感探测结果。
7.根据权利要求1所述的基于光学和雷达遥感数据的区域实际灌溉面积监测方法,其特征在于:步骤9的具体操作为:
基于超像素分割算法的区域图像分割计算结果,对每个分割的区域基于光学和雷达的灌溉探测结果进行判别:
如果两者的探测结果都为潜在灌溉,则此区域标记为灌溉;对于两者探测结果存在分歧的超像素区域,结合步骤4研究区作物种植结构提取结果对存在分歧的超像素区域进行判别:
如果存在分歧超像素区域的不是农作物类型,则剔除;如果存在分歧的超像素区域属于农作物类型,结合步骤1研究区农作物灌溉时间进行判别:
如果在监测时段内存在分歧的超像素区域属于农作物类型为灌溉期,则判定为灌溉,相反,则判定为非灌溉;
基于以上判定,对所有判定为灌溉的区域进行合并,生成最终的灌溉光学和雷达遥感探测结果。
CN202110418925.8A 2021-04-19 2021-04-19 基于光学和雷达遥感数据的区域实际灌溉面积监测方法 Expired - Fee Related CN113128401B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110418925.8A CN113128401B (zh) 2021-04-19 2021-04-19 基于光学和雷达遥感数据的区域实际灌溉面积监测方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110418925.8A CN113128401B (zh) 2021-04-19 2021-04-19 基于光学和雷达遥感数据的区域实际灌溉面积监测方法

Publications (2)

Publication Number Publication Date
CN113128401A CN113128401A (zh) 2021-07-16
CN113128401B true CN113128401B (zh) 2021-11-26

Family

ID=76777671

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110418925.8A Expired - Fee Related CN113128401B (zh) 2021-04-19 2021-04-19 基于光学和雷达遥感数据的区域实际灌溉面积监测方法

Country Status (1)

Country Link
CN (1) CN113128401B (zh)

Families Citing this family (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113642191B (zh) * 2021-08-25 2022-03-22 中国水利水电科学研究院 一种基于短波红外的遥感蒸散模型构建方法
CN114114255B (zh) * 2021-11-29 2022-11-29 商丘师范学院 记录黄河生态资源变化状态的遥感检测方法
CN115376016B (zh) * 2022-08-16 2023-04-11 水利部交通运输部国家能源局南京水利科学研究院 一种基于植被水分指数与蒸散发结合的稻田实际灌溉面积识别方法
CN115100534B (zh) * 2022-08-22 2022-12-06 西南林业大学 森林优势树种识别方法、装置、设备及存储介质
CN115965866B (zh) * 2022-12-28 2023-10-27 广州市嘉卉园林绿化建筑工程有限公司 一种立体绿化水肥遥感监测方法及系统

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109345555A (zh) * 2018-10-15 2019-02-15 中科卫星应用德清研究院 基于多时相多源遥感数据进行水稻识别的方法
CN110703244A (zh) * 2019-09-05 2020-01-17 中国科学院遥感与数字地球研究所 基于遥感数据识别城区内水体的方法和装置
CN110852262A (zh) * 2019-11-11 2020-02-28 南京大学 基于时间序列高分一号遥感影像的农业用地提取方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110334880B (zh) * 2019-07-16 2023-05-26 安徽理工大学 一种有效灌溉面积监测方法

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109345555A (zh) * 2018-10-15 2019-02-15 中科卫星应用德清研究院 基于多时相多源遥感数据进行水稻识别的方法
CN110703244A (zh) * 2019-09-05 2020-01-17 中国科学院遥感与数字地球研究所 基于遥感数据识别城区内水体的方法和装置
CN110852262A (zh) * 2019-11-11 2020-02-28 南京大学 基于时间序列高分一号遥感影像的农业用地提取方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
Combined use of optical and radar satellite data for the detection of tillage and irrigation operations: Case study in Central Morocco;R. Hadria,et al.;《Agricultural Water Management》;20090326;全文 *
农田灌溉遥感监测技术的发展与前景;张威等;《节水灌溉》;20190405(第04期);全文 *

Also Published As

Publication number Publication date
CN113128401A (zh) 2021-07-16

Similar Documents

Publication Publication Date Title
CN113128401B (zh) 基于光学和雷达遥感数据的区域实际灌溉面积监测方法
CN110472184B (zh) 一种基于Landsat遥感数据的多云雨雾地区水稻识别方法
Prabhakara et al. Evaluating the relationship between biomass, percent groundcover and remote sensing indices across six winter cover crop fields in Maryland, United States
Sun et al. Classification mapping and species identification of salt marshes based on a short-time interval NDVI time-series from HJ-1 optical imagery
Leinenkugel et al. Characterisation of land surface phenology and land cover based on moderate resolution satellite data in cloud prone areas—A novel product for the Mekong Basin
Conrad et al. Satellite based calculation of spatially distributed crop water requirements for cotton and wheat cultivation in Fergana Valley, Uzbekistan
CN109685081B (zh) 一种遥感提取撂荒地的联合变化检测方法
CN114821362B (zh) 一种基于多源数据的水稻种植面积提取方法
CN110288647A (zh) 一种基于高分辨率卫星数据监测灌区灌溉面积方法
Cao et al. Mapping paddy rice using Landsat time series data in the Ganfu Plain irrigation system, Southern China, from 1988− 2017
Mtibaa et al. Land cover mapping in cropland dominated area using information on vegetation phenology and multi-seasonal Landsat 8 images
Akbari et al. Crop and land cover classification in Iran using Landsat 7 imagery
Pang et al. Pixel-level rice planting information monitoring in Fujin City based on time-series SAR imagery
Dong et al. Using RapidEye imagery to identify within-field variability of crop growth and yield in Ontario, Canada
Mokarram et al. RELATIONSHIP BETWEEN LAND COVER AND VEGETATION INDICES. CASE STUDY: EGHLID PLAIN, FARS PROVINCE, IRAN.
Jiang et al. Mapping interannual variability of maize cover in a large irrigation district using a vegetation index–phenological index classifier
CN115641504A (zh) 一种基于作物物候特征与决策树模型的田块边界自动化遥感提取方法
Choudhary et al. Rice growth vegetation index 2 for improving estimation of rice plant phenology in costal ecosystems
TAO et al. Fusing multi-source data to map spatio-temporal dynamics of winter rape on the Jianghan Plain and Dongting Lake Plain, China
Sprintsin et al. Relationships between stand density and canopy structure in a dryland forest as estimated by ground-based measurements and multi-spectral spaceborne images
Gorte Land-use and catchment characteristics
Huang et al. Comparing the effects of temporal features derived from synthetic time-series NDVI on fine land cover classification
Fleischmann et al. Multi‐temporal AVHRR digital data: An approach for landcover mapping of heterogeneous landscapes
Diem et al. Assessing the applicability of Fusion Landsat-MODIS data for mapping agricultural land use-A case study in An Giang Province
CN109543729A (zh) 基于特征参数聚类的时间序列数据土地覆被分类方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20211126

CF01 Termination of patent right due to non-payment of annual fee