CN110288647A - 一种基于高分辨率卫星数据监测灌区灌溉面积方法 - Google Patents
一种基于高分辨率卫星数据监测灌区灌溉面积方法 Download PDFInfo
- Publication number
- CN110288647A CN110288647A CN201910553283.5A CN201910553283A CN110288647A CN 110288647 A CN110288647 A CN 110288647A CN 201910553283 A CN201910553283 A CN 201910553283A CN 110288647 A CN110288647 A CN 110288647A
- Authority
- CN
- China
- Prior art keywords
- area
- spectrum
- ndvi
- remote sensing
- value
- 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.)
- Granted
Links
- 238000000034 method Methods 0.000 title claims abstract description 46
- 230000002262 irrigation Effects 0.000 title claims abstract description 45
- 238000003973 irrigation Methods 0.000 title claims abstract description 45
- 238000012544 monitoring process Methods 0.000 title claims abstract description 19
- 238000013461 design Methods 0.000 title claims abstract description 14
- 238000001228 spectrum Methods 0.000 claims abstract description 109
- 230000003595 spectral effect Effects 0.000 claims abstract description 29
- 230000003287 optical effect Effects 0.000 claims abstract description 28
- 238000005070 sampling Methods 0.000 claims abstract description 25
- 238000011160 research Methods 0.000 claims abstract description 20
- 238000009826 distribution Methods 0.000 claims abstract description 16
- 238000005259 measurement Methods 0.000 claims abstract description 12
- 230000003044 adaptive effect Effects 0.000 claims abstract description 7
- 239000000284 extract Substances 0.000 claims abstract description 6
- 238000004445 quantitative analysis Methods 0.000 claims abstract description 4
- 230000011218 segmentation Effects 0.000 claims abstract description 4
- 238000004364 calculation method Methods 0.000 claims description 9
- 238000012937 correction Methods 0.000 claims description 6
- 238000000605 extraction Methods 0.000 claims description 4
- 230000005855 radiation Effects 0.000 claims description 3
- 241000209140 Triticum Species 0.000 description 9
- 235000021307 Triticum Nutrition 0.000 description 9
- 240000008042 Zea mays Species 0.000 description 8
- 235000005824 Zea mays ssp. parviglumis Nutrition 0.000 description 8
- 235000002017 Zea mays subsp mays Nutrition 0.000 description 8
- 235000005822 corn Nutrition 0.000 description 8
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 description 8
- 239000003621 irrigation water Substances 0.000 description 3
- 238000005457 optimization Methods 0.000 description 3
- 238000003066 decision tree Methods 0.000 description 2
- 238000010586 diagram Methods 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 238000005516 engineering process Methods 0.000 description 2
- 238000001556 precipitation Methods 0.000 description 2
- 230000002265 prevention Effects 0.000 description 2
- 240000007594 Oryza sativa Species 0.000 description 1
- 235000007164 Oryza sativa Nutrition 0.000 description 1
- 244000131316 Panax pseudoginseng Species 0.000 description 1
- 235000005035 Panax pseudoginseng ssp. pseudoginseng Nutrition 0.000 description 1
- 235000003140 Panax quinquefolius Nutrition 0.000 description 1
- 238000004458 analytical method Methods 0.000 description 1
- 239000012141 concentrate Substances 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 230000018109 developmental process Effects 0.000 description 1
- 238000001704 evaporation Methods 0.000 description 1
- 230000008020 evaporation Effects 0.000 description 1
- 239000004744 fabric Substances 0.000 description 1
- 235000008434 ginseng Nutrition 0.000 description 1
- 229940023495 irrigation product Drugs 0.000 description 1
- 239000010977 jade Substances 0.000 description 1
- 230000007774 longterm Effects 0.000 description 1
- VMXUWOKSQNHOCA-UKTHLTGXSA-N ranitidine Chemical compound [O-][N+](=O)\C=C(/NC)NCCSCC1=CC=C(CN(C)C)O1 VMXUWOKSQNHOCA-UKTHLTGXSA-N 0.000 description 1
- 238000002310 reflectometry Methods 0.000 description 1
- 235000009566 rice Nutrition 0.000 description 1
- 238000012706 support-vector machine Methods 0.000 description 1
- 230000002123 temporal effect Effects 0.000 description 1
- 238000012360 testing method Methods 0.000 description 1
- 238000012876 topography Methods 0.000 description 1
- 238000013519 translation Methods 0.000 description 1
- 230000000007 visual effect Effects 0.000 description 1
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/22—Matching criteria, e.g. proximity measures
-
- G06T5/80—
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/10—Segmentation; Edge detection
- G06T7/11—Region-based segmentation
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/10—Segmentation; Edge detection
- G06T7/136—Segmentation; Edge detection involving thresholding
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/60—Analysis of geometric attributes
- G06T7/62—Analysis of geometric attributes of area, perimeter, diameter or volume
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/10—Image acquisition modality
- G06T2207/10032—Satellite or aerial image; Remote sensing
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/30—Subject of image; Context of image processing
- G06T2207/30181—Earth observation
- G06T2207/30188—Vegetation; Agriculture
Landscapes
- Engineering & Computer Science (AREA)
- Theoretical Computer Science (AREA)
- Physics & Mathematics (AREA)
- Computer Vision & Pattern Recognition (AREA)
- General Physics & Mathematics (AREA)
- Data Mining & Analysis (AREA)
- Geometry (AREA)
- Life Sciences & Earth Sciences (AREA)
- Artificial Intelligence (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Bioinformatics & Computational Biology (AREA)
- Evolutionary Biology (AREA)
- Evolutionary Computation (AREA)
- General Engineering & Computer Science (AREA)
- Image Processing (AREA)
Abstract
本发明公开一种基于高分辨率卫星数据监测灌区灌溉面积方法：步骤1、卫星遥感数据的获取和预处理；步骤2、实测样点数据；步骤3、根据步骤2获得的真实样点数据提取端元光谱；步骤4、采用光谱匹配法中的统计算法和光谱波形特征算法计算光谱相似度，通过三个指标对研究区主要作物端元光谱与目标光谱匹配程度进行定量分析；步骤5、采用OTSU自适应阈值算法计算SSV分割阈值来判断是否为灌溉区域，小于该阈值即为灌溉区域，从而识别研究区的灌溉面积空间分布情况，最终获得灌溉区域范围。本发明方法适用于高分辨率卫星遥感数据，满足提取小地块灌溉信息的需求，能够提高灌溉面积监测结果的精度。
Description
技术领域
本发明涉及一种基于高分辨率卫星数据监测灌区灌溉面积方法，尤其是一种基于高分辨率卫星数据，通过光谱匹配方法像元尺度应用，并引入OTSU自适应阈值算法，构建的高分辨率灌溉面积遥感监测方法，属于水利遥感技术领域。
背景技术
准确监测灌溉面积及其空间分布对于预防干旱、优化灌区水资源利用和水资源调控配置具有重要意义。我国作为一个农业大国，灌溉面积占总耕地面积的近50％，在全球气候变暖、干旱等极端事件增加、经济发展用水需求猛增情况下，水资源短缺形式日趋严峻。准确监测灌溉面积时空分布及其变化，是掌握灌溉用水量、灌溉用水去向的关键，对监测灌溉系统运转、预防干旱、优化作物种植结构、提高灌溉用水效率、优化灌区水资源配置均具有重要意义。
国内外农业灌溉面积监测研究在决策树算法、时空螺旋曲线与变化矢量分析、物候特征、支持向量机、自动农田算法、光谱匹配、k-means和ISOCLASS算法、自动机器学习算法、植被指数与干旱指数反演等方法方面取得诸多进展。其中，光谱匹配技术(SMTs)能够提取长时间序列多期遥感影像的光谱反射率或特征指数变化特征，在周期性农业灌溉面积监测方面具有独特优势。Thenkabai最先应用这种方法开展土地类型识别与灌溉面积识别研究。2006年，世界水资源管理研究所(IWMI)通过使用多种数据源，先对全球范围进行分类，而后使用光谱匹配技术对全球范围内的土地利用类型、灌溉面积分布进行定类与归类，发布了世界第一份1km空间分辨率的全球灌溉区域图(GIAM)。2012年，IWMI又针对非洲和亚洲进行了重点研究，绘制了250m空间分辨率的灌溉与雨养区分布图。2016年，Ambika将这种方法应用到了印度，使用了250m空间分辨率的MODIS NDVI数据(MOD13Q1)，先通过匹配计算分类，再使用决策树通过NDVI阈值划分灌溉区域，获得了印度2000-2015年的灌溉区域分布图，相比IWMI 250m空间分辨率的灌溉产品提升了精度。Teluguntla结合定量光谱匹配技术(QSMTs)与农田自动分类算法,使用连续14年的250m空间分辨率MODIS NDVI产品数据对澳大利亚的灌溉信息与种植强度进行了提取。
我国具有气候多样、农业种植结构复杂、田块破碎、灌区信息化水平不高等国情，农业统计缺乏空间信息，目前可公开获取的灌溉遥感产品缺乏且空间分辨率与精度难以满足我国实际需求，我国灌溉面积遥感监测对数据空间分辨率要求更高，灌溉面积遥感监测技术需要取得进一步突破。本发明针对我国国情，应用较高分辨率的卫星遥感数据，通过光谱匹配方法像元尺度应用，并引入自适应阈值算法，构建了高分辨率灌溉面积遥感监测新方法。
发明内容
本发明的目的在于提供一种基于高分辨率卫星数据监测灌区灌溉面积方法，以解决现有技术中农业统计缺乏空间信息，可公开获取的灌溉遥感产品缺乏且空间分辨率与精度难以满足我国国情的问题，构建一种高分辨率灌溉面积遥感监测新方法。
本发明采用如下技术手段以实现上述目的，一种基于高分辨率卫星数据监测灌区灌溉面积方法，具体步骤如下：
步骤1、卫星遥感数据的获取和预处理
其中，所述的预处理包括：对遥感数据进行辐射定标、几何校正、大气校正，并计算NDVI得到研究区各像元的NDVI时间序列曲线，用于提取灌溉信息。
步骤2、实测样点数据
结合历史灌溉情况与步骤1获取的卫星遥感数据，选择具有典型代表意义的灌溉区域与非灌溉区域作为真实样点数据；特别的，在选取灌溉区域时应避免附近有建筑物与道路的干扰，避开农田的边界处。
其中，所述的真实样点数据包括：样点经纬度位置、作物类型及灌溉情况等信息。
其中，所述的选择具有典型代表意义的灌溉区域与非灌溉区域作为真实样点数据，具体为：每个作物或地类上取10-20个点。
步骤3、根据步骤2获得的真实样点数据提取端元光谱；所述的端元光谱C＝(s1，s2，…,sn)为由实地采样点位置所在卫星像元提取得到的已知作物类型与灌溉情况的样本光谱；其中，C为端元光谱；s1…sn为各实地采样点位置所在卫星像元对应的多期NDVI值曲线；
步骤4、采用光谱匹配法中的统计算法和光谱波形特征算法计算光谱相似度，通过三个指标对研究区主要作物端元光谱(样本光谱)与目标光谱(待测像元光谱)匹配程度进行定量分析；所述的三个指标为：形状定量、距离定量及形状距离综合定量；具体过程如下：
所述的研究区主要作物端元光谱，为步骤3根据真实样点数据提取的端元光谱C＝(s1，s2，…,sn)，所述的目标光谱为从卫星影像各像元中获取的待确定灌溉与否的光谱Xi,j＝(b1,b2,…,bn)，b1….bn为卫星影像各像元对应的多期NDVI曲线；
S41、形状定量：光谱关联相似度SCS；SCS通过如下公式计算：
这里Xi,j为目标光谱的NDVI值，C为端元光谱的NDVI值，i,j分别代表第i行，第j列，n为数据集层数，μi,j为目标光谱NDVI均值，μs为端元光谱NDVI均值；σX为目标光谱的NDVI标准偏差，σC为端元光谱的NDVI标准偏差。
SCS是衡量端元光谱与目标光谱NDVI时间序列曲线形状相似度的指标。SCS值越高，端元光谱与目标光谱的NDVI时间序列曲线就越相似。当SCS值为0时，相似度最小；当SCS为1时，相似度最大，即两种光谱完全一致。所有SCS值在0-1范围之外的像元即可认为与端元光谱相似度差异过大，直接剔除。
S42、距离定量：欧氏距离EDS
EDS是用于衡量端元光谱与目标光谱在光谱空间中距离的指标，在光谱空间中距离越近则越相似。EDS计算公式如下:
通常欧式距离越大则代表着两种光谱差异性越大，反之则差异性越小。为了方便计算将上述公式获得结果通过归一化处理，将其标准范围处理至0到1。具体公式计算如下所示：
EDSnormal＝(EDS-m)/(M-m) (3)
式中，EDSnormal值的范围为0-1，它从光谱特征空间距离上测量端元光谱与目标光谱NDVI时间序列曲线间相关性的大小。M，m分别代表最大、最小的欧式距离。通过对EDS的计算可以把SCS运算结果中满足条件像元点进行计算，理想状态下，0表示二者完全一致，1表示二者完全不相关。此外，EDS受像元数量的影响，一旦有新的像元参与计算，该式中的M和m即可能会发生变化，对不同研究对象有较好的变动适应性。
S43、形状距离综合定量：光谱相似值SSV
SSV综合了形状定量SCS与距离定量EDS的特点，从NDVI时间序列曲线形状和光谱特征空间距离两方面衡量了端元光谱与目标光谱间的形似度，其计算公式如下:
式中，SSV值越小，光谱之间越相似，其范围通常介于0-1.414之间。不同作物之间，NDVI时间序列曲线的差异较大，SSV值高；对于同一种作物，灌溉区域的NDVI时间序列曲线比非灌溉区域具有更高的一致性，SSV值更低。
步骤5、采用OTSU自适应阈值算法计算SSV分割阈值来判断是否为灌溉区域，小于该阈值即为灌溉区域，从而识别研究区的灌溉面积空间分布情况，最终获得灌溉区域范围。
具体算法如下：
当遥感影像数据量级为L(G＝0,1,……,L-1)时，初始阈值为h，将图像分为目标区域T与背景区域B两部分，pi代表SSV为i的像元出现的概率，ni代表SSV为i的像元数量，N代表像元总数，T、B两部分概率为：
对应的T、B两部分的均值为：
则存在阈值h0即为最佳阈值：
h0＝Argmax(PTPB(μT-μB)2) (7)。
本发明的一种基于高分辨率卫星数据监测灌区灌溉面积方法，有别于IWMI采用的聚类光谱匹配，该新方法针对高分辨率遥感影像，将光谱匹配方法应用于像元尺度，保证所有像元的时序光谱曲线参与匹配计算，减少聚类过程造成的误差，并引入OTSU自适应阈值算法自动确定灌溉面积提取阈值。OTSU算法的引入可保证每个像元计算相似度信息的同时，自动确定相似度合理变化范围，对于不同年份与不同数据情况有良好变动适应性，无需人为干预即可提取灌溉区域。该新方法适用于高分辨率卫星遥感数据，满足提取小地块灌溉信息的需求，能够提高灌溉面积监测结果的精度。可广泛适用于高空间分辨率的灌溉区域遥感识别提取，可以更有效区分城镇、村庄、道路、渠道、田块等边界信息，更有效识别小田块灌溉分布，在作物种植强度及其灌溉区域分布空间信息提取方面更符合我国国情。
附图说明
图1所示为本发明方法流程框图；
图2所示为本发明实施例研究区内选取的225个真实样点数据分布图；
图3a为本发明实施例玉米灌溉区(6月至10月)9个实测样本点的NDVI时间序列曲线；
图3b本发明实施例小麦灌溉区(10月至6月)10个实测样本点的NDVI时间序列曲线；
图4a-图4f本发明实施例像元尺度光谱匹配算法得到的主要作物玉米和小麦灌溉区域识别相似性参数SCS、EDS和SSV分布图；
图5所示为OTSU法计算直方图，其中图5a为小麦OTSU法获得阈值，图5b为玉米OTSU法获得阈值；
图6所示为本发明实施例灌溉面积分布图；
图7所示为本发明方法与IWMI方法识别效果对比图，图中(a)为GF-1原图，(b)为本发明结果，(c)为IWMI-Irrigated Area Map Asia(2000-2010)结果；
具体实施方式
下面结合附图和实施例对本发明的技术方案做进一步的说明。
本发明实施例选择我国西北半干旱区典型渠灌灌区即东雷二期抽黄灌区为研究区，对其2018年的小麦与玉米灌溉面积开展了遥感识别提取，通过野外采样精度验证，结果总体精度为88.27％(Kappa＝0.8308)。
东雷二期抽黄灌区，位于东经109°10′-110°10′，北纬34°41′-35°00′，横跨陕西省富平县、大荔县和蒲城县，属黄土波状台原区，海拔介于385m-600m，地势东南低，西北高。灌区属温带大陆性季风气候，春季温暖干旱，夏季多暴雨，并有短暂大风，冬季严寒干燥，雨雪稀少，年平均降水量519mm-552mm，干旱年份仅为360mm-390mm，且降水时空分配不均衡，50％以上雨量集中在7月-9月，年蒸发量约为1700mm-2000mm。灌区灌溉方式为抽取黄河水进行渠道引水漫灌。设计灌溉面积126.5万亩，分春灌、夏灌和冬灌。主要作物为小麦和玉米，占灌区种植总面积的90％以上。灌溉作物生长周期及灌溉时间详见表1。
表1
所述的一种基于高分辨率卫星数据监测灌区灌溉面积方法，如图1所示，具体包括如下步骤：
步骤1、卫星遥感数据的获取和预处理
S11、所述的数据源为2018年GF-1号16m多光谱数据，来源为中国资源卫星公共网站平台(http://www.cresda.com/CN/)，详细参数见表2GF-1多光谱传感器信息。使用2018年研究区云覆盖量低于15％的影像，并编译成由92个波段组成的数据集(共有23景数据，每景4个波段)。此数据集用于准确获取不同类别的时序光谱特征信息。
表2
S12、对遥感数据进行辐射定标、几何校正、大气校正等预处理，并计算NDVI。所述的NDVI能够较好地反映植被生长状态，进行归一化处理后可以更好地限制数据范围，有利于后续的光谱匹配计算。
步骤2、实测样点数据
结合历史灌溉情况与步骤1获取的卫星遥感影像，选择具有典型代表意义的灌溉区域与非灌溉区域，在选取灌溉区域时应注意避免附近有建筑物与道路的干扰，避开农田的边界处，以防在取点时样本点中同时包含农田、草地或其他干扰地物。最终在研究区内共选取225个真实样点数据(见图2)，主要记录样点位置、作物类型及灌溉情况等信息。
步骤3、根据步骤2获得的真实样点信息，从步骤1所获得的卫星遥感影像像元中提取得到端元光谱(见图3a、b)；
步骤4、采用光谱匹配法中的统计算法和光谱波形特征算法计算光谱相似度，通过三个指标对研究区主要作物端元光谱与目标光谱匹配程度进行定量分析。所述的三个指标包括：形状定量、距离定量及形状距离综合定量。
具体算法如下：
S41、形状定量：光谱关联相似度SCS通过如下公式计算：
这里Xi,j为目标光谱的NDVI值，C为端元光谱的NDVI值，i,j分别代表第i行，第j列，n为数据集层数，μi,j为目标光谱NDVI均值，μs为端元光谱NDVI均值；σX为目标光谱的NDVI标准偏差，σC为端元光谱的NDVI标准偏差。
SCS是衡量端元光谱与目标光谱NDVI时间序列曲线形状相似度的指标。SCS值越高，端元光谱与目标光谱的NDVI时间序列曲线就越相似。当SCS值为0时，相似度最小；当SCS为1时，相似度最大，即两种光谱完全一致。所有SCS值在0-1范围之外的像元即可认为与端元光谱相似度差异过大，直接剔除。
S42、距离定量：欧氏距离EDS
EDS是用于衡量端元光谱与目标光谱在光谱空间中距离的指标，在光谱空间中距离越近则越相似。EDS计算公式如下:
通常欧式距离越大则代表着两种光谱差异性越大，反之则差异性越小。为了方便计算将上述公式获得结果通过归一化处理，将其标准范围处理至0到1。具体公式计算如下所示：
EDSnormal＝(EDS-m)/(M-m) (3)
式中，EDSnormal值的范围为0-1，它从光谱特征空间距离上测量端元光谱与目标光谱NDVI时间序列曲线间相关性的大小。M，m分别代表最大、最小的欧式距离。通过对EDS的计算可以把SCS运算结果中满足条件像元点进行计算，理想状态下，0表示二者完全一致，1表示二者完全不相关。此外，EDS受像元数量的影响，一旦有新的像元参与计算，该式中的M和m即可能会发生变化，对不同研究对象有较好的变动适应性。
尽管使用欧氏距离计算光谱匹配度是一种简单直观的方法，但其存在一定的局限性。由于不同年份的气候情况不同可能导致NDVI光谱曲线会发生上下平移，因此，本该被判断同一种类的可能被判定成不同种类。这说明，仅使用光谱距离的特征描述往往还不足以较好的完成光谱匹配的工作，需要进一步对端元光谱与目标光谱的相似程度进行量化。
S43、形状距离综合定量：光谱相似值SSV
SSV综合了形状定量SCS与距离定量EDS的特点，从NDVI时间序列曲线形状和光谱特征空间距离两方面衡量了端元光谱与目标光谱间的形似度，其计算公式如下:
式中，SSV值越小，光谱之间越相似，其范围通常介于0-1.414之间。不同作物之间，NDVI时间序列曲线的差异较大，SSV值高；对于同一种作物，灌溉区域的NDVI时间序列曲线比非灌溉区域具有更高的一致性，SSV值更低。
根据像元尺度光谱匹配算法得到的主要作物玉米和小麦灌溉区域识别相似性参数SCS、EDS和SSV分布如图4a-图4f所示。
步骤5、采用OTSU自适应阈值算法计算SSV分割阈值来判断是否为灌溉区域，小于该阈值即为灌溉区域，从而识别研究区的灌溉面积空间分布情况，最终获得灌溉区域范围。
以本发明实施例实地采样获取的灌溉区域玉米和小麦的NDVI时间序列曲线分别作为端元光谱，利用上述衡量光谱相似度的三个指标计算研究区每个像元的SSV值。同类型的光谱相似度高，SSV值低；不同类的光谱相似度低，SSV值更高。需要确定合理阈值进行判定，SSV值小于该阈值时与端元光谱为同一类别。因此引入OTSU算法，其是无参量的一种自适应阈值选取方法，该方法可以给定一个临时阈值，然后计算该临时阈值两侧范围的像元SSV的方差值，当方差达到最大，所对应的像元SSV为最佳阈值，即目标区域与背景区域差异最大化，而后将相似度大于或等于该阈值的归并成一类，小于该阈值的归并成另一类，得到相应的二值化图像。具体算法如下：
当影像数据量级为L(G＝0,1,……,L-1)时，初始阈值为h，将图像分为目标区域T与背景区域B两部分，pi代表SSV为i的像元出现的概率，ni代表SSV为i的像元数量，N代表像元总数，T、B两部分概率为：
对应的T、B两部分的均值为：
则存在阈值h0即为最佳阈值：
h0＝Argmax(PTPB(μT-μB)2) (7)。
本发明实施例使用OTSU算法，获得了两种作物对应的相似度阈值(见图5)，结果表明小麦灌溉区最佳分割阈值为0.3985(图5a)，玉米灌溉区最佳分割阈值为0.4639(图5b)。
本发明实施例最终获得如图6所示的灌溉面积分布图，研究区内主要作物根据种植强度可分为双季作物灌溉区、单季小麦作物灌溉区、单季玉米作物灌溉区。
进一步的，将本发明方法与现有国际研究技术与成果(水管理研究所灌溉数据产品)进行比较，如图7所示，图中(a)为GF-1原图，(b)为本发明结果，(c)为IWMI-IrrigatedArea Map Asia(2000-2010)结果；如图所示，本发明方法使用GF-1像元尺度匹配方法比IWMI-Irrigated Area Map Asia(2000-2010)的识别效果更好，可以有效的区分出建筑物、道路与田块信息，而IWMI结果无法获取灌溉分布细节，较多的建筑物与道路存在混淆错分的情况，通过目视可以看出，IWMI结果并不能很好的反映出灌区内的实际情况，本发明的基于像元尺度的光谱匹配方法更适合应用于高精度的研究区域研究，更适合农业种植结构复杂、田块破碎、灌区信息化水平不高的地区。
Claims (6)
1.一种基于高分辨率卫星数据监测灌区灌溉面积方法，其特征在于：该方法具体步骤如下：
步骤1、卫星遥感数据的获取和预处理
步骤2、实测样点数据
结合历史灌溉情况与步骤1获取的卫星遥感数据，选择具有典型代表意义的灌溉区域与非灌溉区域作为真实样点数据；
步骤3、根据步骤2获得的真实样点数据提取端元光谱；所述的端元光谱C＝(s1，s2，…,sn)为由实地采样点位置所在卫星像元提取得到的已知作物类型与灌溉情况的样本光谱；其中，C为端元光谱；s1…sn为各实地采样点位置所在卫星像元对应的多期NDVI值曲线；
步骤4、采用光谱匹配法中的统计算法和光谱波形特征算法计算光谱相似度，通过三个指标对研究区主要作物端元光谱与目标光谱匹配程度进行定量分析；所述的三个指标为：形状定量、距离定量及形状距离综合定量；所述的研究区主要作物端元光谱为步骤3根据真实样点数据提取的端元光谱C＝(s1，s2，…,sn)，所述的目标光谱为从卫星影像各像元中获取的待确定灌溉与否的光谱Xi,j＝(b1,b2,…,bn)，b1….bn为卫星影像各像元对应的多期NDVI曲线；
步骤5、采用OTSU自适应阈值算法计算SSV分割阈值来判断是否为灌溉区域，小于该阈值即为灌溉区域，从而识别研究区的灌溉面积空间分布情况，最终获得灌溉区域范围。
2.根据权利要求1所述的一种基于高分辨率卫星数据监测灌区灌溉面积方法，其特征在于：步骤1所述的预处理包括：对遥感数据进行辐射定标、几何校正、大气校正，并计算NDVI得到研究区各像元的NDVI时间序列曲线，用于提取灌溉信息。
3.根据权利要求1所述的一一种基于高分辨率卫星数据监测灌区灌溉面积方法，其特征在于：步骤2所述的真实样点数据包括：样点经纬度位置、作物类型及灌溉情况。
4.根据权利要求1所述的一种基于高分辨率卫星数据监测灌区灌溉面积方法，其特征在于：步骤2所述的选择具有典型代表意义的灌溉区域与非灌溉区域作为真实样点数据，具体为：每个作物或地类上取10-20个点。
5.根据权利要求1所述的一种基于高分辨率卫星数据监测灌区灌溉面积方法，其特征在于：所述步骤4的具体过程如下：
S41、形状定量：光谱关联相似度SCS；SCS通过如下公式计算：
Xi,j为目标光谱的NDVI值，C为端元光谱的NDVI值，i,j分别代表第i行，第j列，n为数据集层数，μi,j为目标光谱NDVI均值，μs为端元光谱NDVI均值；σX为目标光谱的NDVI标准偏差，σC为端元光谱的NDVI标准偏差；
SCS是衡量端元光谱与目标光谱NDVI时间序列曲线形状相似度的指标，SCS值越高，端元光谱与目标光谱的NDVI时间序列曲线就越相似；当SCS值为0时，相似度最小；当SCS为1时，相似度最大，即两种光谱完全一致；所有SCS值在0-1范围之外的像元即可认为与端元光谱相似度差异过大，直接剔除；
S42、距离定量：欧氏距离EDS
EDS计算公式如下:
通常欧式距离越大则代表着两种光谱差异性越大，反之则差异性越小，为了方便计算将上述公式获得结果通过归一化处理，将其标准范围处理至0到1；具体公式计算如下所示：
EDSnormal＝(EDS-m)/(M-m) (3)
式中，EDSnormal值的范围为0-1，它从光谱特征空间距离上测量端元光谱与目标光谱NDVI时间序列曲线间相关性的大小；M，m分别代表最大、最小的欧式距离；
S43、形状距离综合定量：光谱相似值SSV
SSV综合了形状定量SCS与距离定量EDS的特点，从NDVI时间序列曲线形状和光谱特征空间距离两方面衡量了端元光谱与目标光谱间的形似度，其计算公式如下:
式中，SSV值越小，光谱之间越相似，其范围通常介于0-1.414之间；不同作物之间，NDVI时间序列曲线的差异较大，SSV值高；对于同一种作物，灌溉区域的NDVI时间序列曲线比非灌溉区域具有更高的一致性，SSV值更低。
6.根据权利要求1所述的一种基于高分辨率卫星数据监测灌区灌溉面积方法，其特征在于：所述步骤5的具体过程如下：
当遥感影像数据量级为L(G＝0,1,……,L-1)时，初始阈值为h，将图像分为目标区域T与背景区域B两部分，pi代表SSV为i的像元出现的概率，ni代表SSV为i的像元数量，N代表像元总数，T、B两部分概率为：
对应的T、B两部分的均值为：
则存在阈值h0即为最佳阈值：
h0＝Argmax(PTPB(μT-μB)2) (7)。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910553283.5A CN110288647B (zh) | 2019-06-25 | 2019-06-25 | 一种基于高分辨率卫星数据监测灌区灌溉面积方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910553283.5A CN110288647B (zh) | 2019-06-25 | 2019-06-25 | 一种基于高分辨率卫星数据监测灌区灌溉面积方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN110288647A true CN110288647A (zh) | 2019-09-27 |
| CN110288647B CN110288647B (zh) | 2021-05-07 |
Family
ID=68005599
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201910553283.5A Active CN110288647B (zh) | 2019-06-25 | 2019-06-25 | 一种基于高分辨率卫星数据监测灌区灌溉面积方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN110288647B (zh) |
Cited By (10)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110738196A (zh) * | 2019-11-13 | 2020-01-31 | 云南农业大学 | 一种基于区域墒情监测和遥感数据的实时灌溉预报系统 |
| CN112418314A (zh) * | 2020-11-23 | 2021-02-26 | 广东弓叶科技有限公司 | 一种光谱相似度匹配分类中阈值设定方法及装置 |
| CN113221789A (zh) * | 2021-05-21 | 2021-08-06 | 广东省科学院广州地理研究所 | 一种田块灌溉设施识别方法及装置 |
| CN113378848A (zh) * | 2021-07-08 | 2021-09-10 | 湖南元想科技有限公司 | 南荻收割区或弃割区分布提取方法、系统及存储介质 |
| CN113470175A (zh) * | 2021-06-11 | 2021-10-01 | 北京大学 | 基于光学梯形模型的灌溉区域制图方法 |
| CN113793376A (zh) * | 2021-09-14 | 2021-12-14 | 安徽理工大学 | 一种灌溉水体提取方法 |
| CN114254266A (zh) * | 2021-12-22 | 2022-03-29 | 中国水利水电科学研究院 | 智慧水利数字化场景构建中的一种灌溉信息遥感监测方法 |
| CN114780904A (zh) * | 2022-06-17 | 2022-07-22 | 中国科学院、水利部成都山地灾害与环境研究所 | 一种端元自适应的山地植被覆盖度遥感反演方法 |
| CN115126686A (zh) * | 2022-08-31 | 2022-09-30 | 山东中聚电器有限公司 | 用于植保无人机载隔膜泵控制系统 |
| CN117011363A (zh) * | 2023-08-17 | 2023-11-07 | 中国水利水电科学研究院 | 一种基于modis地表温度日较差的灌区灌溉面积提取阈值确定方法 |
Citations (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN1837787A (zh) * | 2006-04-19 | 2006-09-27 | 南京大学 | 一种无损伤精准测定棉花生物物理参数的方法 |
| CN102012528A (zh) * | 2010-11-23 | 2011-04-13 | 北京理工大学 | 一种稀疏植被地区的高光谱遥感油气勘探方法 |
| CN104424390A (zh) * | 2013-09-06 | 2015-03-18 | 中国水利水电科学研究院 | 灌溉面积监测方法及装置 |
| EP2980669A2 (en) * | 2014-08-01 | 2016-02-03 | AGCO Corporation | Determining field characterisitics using optical recognition |
| CN106960190A (zh) * | 2017-03-22 | 2017-07-18 | 中国水利水电科学研究院 | 基于高分辨率图像监督中分辨率图像的灌溉面积统计方法 |
| AU2017204615A1 (en) * | 2016-07-27 | 2018-02-15 | Indicina Pty Ltd | Irrigation system and method |
| CN107864827A (zh) * | 2017-12-18 | 2018-04-03 | 无锡同春新能源科技有限公司 | 基于北斗导航定位的无人机调配海水稻灌溉水的灌溉系统 |
| CN207263153U (zh) * | 2017-10-13 | 2018-04-20 | 胡勇平 | 一种农业机械作业面积测量仪 |
-
2019
- 2019-06-25 CN CN201910553283.5A patent/CN110288647B/zh active Active
Patent Citations (8)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN1837787A (zh) * | 2006-04-19 | 2006-09-27 | 南京大学 | 一种无损伤精准测定棉花生物物理参数的方法 |
| CN102012528A (zh) * | 2010-11-23 | 2011-04-13 | 北京理工大学 | 一种稀疏植被地区的高光谱遥感油气勘探方法 |
| CN104424390A (zh) * | 2013-09-06 | 2015-03-18 | 中国水利水电科学研究院 | 灌溉面积监测方法及装置 |
| EP2980669A2 (en) * | 2014-08-01 | 2016-02-03 | AGCO Corporation | Determining field characterisitics using optical recognition |
| AU2017204615A1 (en) * | 2016-07-27 | 2018-02-15 | Indicina Pty Ltd | Irrigation system and method |
| CN106960190A (zh) * | 2017-03-22 | 2017-07-18 | 中国水利水电科学研究院 | 基于高分辨率图像监督中分辨率图像的灌溉面积统计方法 |
| CN207263153U (zh) * | 2017-10-13 | 2018-04-20 | 胡勇平 | 一种农业机械作业面积测量仪 |
| CN107864827A (zh) * | 2017-12-18 | 2018-04-03 | 无锡同春新能源科技有限公司 | 基于北斗导航定位的无人机调配海水稻灌溉水的灌溉系统 |
Non-Patent Citations (2)
| Title |
|---|
| 沈静: "遥感技术在灌溉面积监测上的应用研究", 《中国优秀硕士学位论文全文数据库》 * |
| 黄涛等: "遥感技术在灌区灌溉面积及作物种植结构中的应用", 《科技信息》 * |
Cited By (13)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN110738196A (zh) * | 2019-11-13 | 2020-01-31 | 云南农业大学 | 一种基于区域墒情监测和遥感数据的实时灌溉预报系统 |
| CN112418314A (zh) * | 2020-11-23 | 2021-02-26 | 广东弓叶科技有限公司 | 一种光谱相似度匹配分类中阈值设定方法及装置 |
| CN113221789A (zh) * | 2021-05-21 | 2021-08-06 | 广东省科学院广州地理研究所 | 一种田块灌溉设施识别方法及装置 |
| CN113470175A (zh) * | 2021-06-11 | 2021-10-01 | 北京大学 | 基于光学梯形模型的灌溉区域制图方法 |
| CN113378848A (zh) * | 2021-07-08 | 2021-09-10 | 湖南元想科技有限公司 | 南荻收割区或弃割区分布提取方法、系统及存储介质 |
| CN113378848B (zh) * | 2021-07-08 | 2024-01-30 | 湖南元想科技有限公司 | 南荻收割区或弃割区分布提取方法、系统及存储介质 |
| CN113793376B (zh) * | 2021-09-14 | 2023-12-05 | 安徽理工大学 | 一种灌溉水体提取方法 |
| CN113793376A (zh) * | 2021-09-14 | 2021-12-14 | 安徽理工大学 | 一种灌溉水体提取方法 |
| CN114254266A (zh) * | 2021-12-22 | 2022-03-29 | 中国水利水电科学研究院 | 智慧水利数字化场景构建中的一种灌溉信息遥感监测方法 |
| CN114780904B (zh) * | 2022-06-17 | 2022-09-27 | 中国科学院、水利部成都山地灾害与环境研究所 | 一种端元自适应的山地植被覆盖度遥感反演方法 |
| CN114780904A (zh) * | 2022-06-17 | 2022-07-22 | 中国科学院、水利部成都山地灾害与环境研究所 | 一种端元自适应的山地植被覆盖度遥感反演方法 |
| CN115126686A (zh) * | 2022-08-31 | 2022-09-30 | 山东中聚电器有限公司 | 用于植保无人机载隔膜泵控制系统 |
| CN117011363A (zh) * | 2023-08-17 | 2023-11-07 | 中国水利水电科学研究院 | 一种基于modis地表温度日较差的灌区灌溉面积提取阈值确定方法 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN110288647B (zh) | 2021-05-07 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN110288647A (zh) | 一种基于高分辨率卫星数据监测灌区灌溉面积方法 | |
| Pan et al. | Mapping cropping intensity in Huaihe basin using phenology algorithm, all Sentinel-2 and Landsat images in Google Earth Engine | |
| CN108846832B (zh) | 一种基于多时相遥感影像与gis数据的变化检测方法及系统 | |
| Li et al. | Estimating bamboo forest aboveground biomass using EnKF-assimilated MODIS LAI spatiotemporal data and machine learning algorithms | |
| Ruelland et al. | Comparison of methods for LUCC monitoring over 50 years from aerial photographs and satellite images in a Sahelian catchment | |
| Scudiero et al. | Remote sensing is a viable tool for mapping soil salinity in agricultural lands | |
| CN110472184A (zh) | 一种基于Landsat遥感数据的多云雨雾地区水稻识别方法 | |
| CN113128401B (zh) | 基于光学和雷达遥感数据的区域实际灌溉面积监测方法 | |
| Dubovyk et al. | Spatial targeting of land rehabilitation: A relational analysis of cropland productivity decline in arid Uzbekistan | |
| Bretreger et al. | Monitoring irrigation using landsat observations and climate data over regional scales in the Murray-Darling Basin | |
| CN108592888A (zh) | 一种居民地提取方法 | |
| Jiang et al. | Mapping interannual variability of maize cover in a large irrigation district using a vegetation index–phenological index classifier | |
| EP4081959A1 (en) | Agricultural decision support system | |
| Ji et al. | Object-based mapping of plastic greenhouses with scattered distribution in complex land cover using Landsat 8 OLI images: a case study in Xuzhou, China | |
| Deng et al. | Winter wheat planting area extraction technique using multi-temporal remote sensing images based on field parcel | |
| CN102663396A (zh) | 一种水稻乳熟期自动检测的方法 | |
| Guermazi et al. | Water irrigation management using remote sensing techniques: a case study in Central Tunisia | |
| Warner et al. | Remote sensing of land cover change | |
| Gorte | Land-use and catchment characteristics | |
| Chen et al. | Water requirement for irrigation of complicated agricultural land by using classified airborne digital sensor images | |
| Gomez et al. | Sentinel-2 images to assess soil surface characteristics over a rainfed Mediterranean cropping system | |
| Tiwari et al. | In-season crop-area mapping for wheat and rice in Afghanistan and Bangladesh | |
| CN109543729A (zh) | 基于特征参数聚类的时间序列数据土地覆被分类方法 | |
| Yang et al. | Extraction of multiple cropping information at the Sub-pixel scale based on phenology and MODIS NDVI time-series: a case study in Henan Province, China | |
| Cheema et al. | Quantification of land use changes in complex cropping of irrigated Indus basin, Pakistan using MODIS vegetation time series data. |
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 |