CN113361398A - 草场围栏识别方法、装置及存储介质 - Google Patents
草场围栏识别方法、装置及存储介质 Download PDFInfo
- Publication number
- CN113361398A CN113361398A CN202110625305.1A CN202110625305A CN113361398A CN 113361398 A CN113361398 A CN 113361398A CN 202110625305 A CN202110625305 A CN 202110625305A CN 113361398 A CN113361398 A CN 113361398A
- Authority
- CN
- China
- Prior art keywords
- image
- ndvi
- resolution
- time
- image data
- 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 52
- 238000000605 extraction Methods 0.000 claims abstract description 24
- 230000004927 fusion Effects 0.000 claims abstract description 17
- 238000007781 pre-processing Methods 0.000 claims abstract description 13
- 238000000513 principal component analysis Methods 0.000 claims description 27
- 230000002123 temporal effect Effects 0.000 claims description 17
- 230000008439 repair process Effects 0.000 claims description 14
- 238000002310 reflectometry Methods 0.000 claims description 8
- 238000003709 image segmentation Methods 0.000 claims description 5
- 238000012952 Resampling Methods 0.000 claims description 4
- 230000015572 biosynthetic process Effects 0.000 claims description 4
- 238000012847 principal component analysis method Methods 0.000 claims description 4
- 238000003786 synthesis reaction Methods 0.000 claims description 4
- 238000012937 correction Methods 0.000 claims description 3
- 230000000873 masking effect Effects 0.000 claims description 3
- 230000011218 segmentation Effects 0.000 description 10
- 238000012545 processing Methods 0.000 description 8
- 238000004364 calculation method Methods 0.000 description 6
- 238000009304 pastoral farming Methods 0.000 description 6
- 244000025254 Cannabis sativa Species 0.000 description 5
- 238000011156 evaluation Methods 0.000 description 4
- 239000011159 matrix material Substances 0.000 description 4
- 238000012986 modification Methods 0.000 description 4
- 230000004048 modification Effects 0.000 description 4
- 230000000007 visual effect Effects 0.000 description 4
- 238000010586 diagram Methods 0.000 description 3
- 238000005516 engineering process Methods 0.000 description 3
- 230000003203 everyday effect Effects 0.000 description 3
- 230000000670 limiting effect Effects 0.000 description 3
- 230000008569 process Effects 0.000 description 3
- 238000012795 verification Methods 0.000 description 3
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 description 3
- 238000004458 analytical method Methods 0.000 description 2
- 238000010276 construction Methods 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 230000003287 optical effect Effects 0.000 description 2
- 230000009467 reduction Effects 0.000 description 2
- 230000002829 reductive effect Effects 0.000 description 2
- 239000004576 sand Substances 0.000 description 2
- 230000003068 static effect Effects 0.000 description 2
- 101100425560 Rattus norvegicus Tle4 gene Proteins 0.000 description 1
- 101100327317 Saccharomyces cerevisiae (strain ATCC 204508 / S288c) CDC1 gene Proteins 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 238000006243 chemical reaction Methods 0.000 description 1
- 230000002354 daily effect Effects 0.000 description 1
- 230000006866 deterioration Effects 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 230000018109 developmental process Effects 0.000 description 1
- 230000002708 enhancing effect Effects 0.000 description 1
- 238000000556 factor analysis Methods 0.000 description 1
- 230000006870 function Effects 0.000 description 1
- 230000003993 interaction Effects 0.000 description 1
- 238000011835 investigation Methods 0.000 description 1
- 238000003064 k means clustering Methods 0.000 description 1
- 230000007774 longterm Effects 0.000 description 1
- 230000007246 mechanism Effects 0.000 description 1
- 230000036961 partial effect Effects 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/21—Design or setup of recognition systems or techniques; Extraction of features in feature space; Blind source separation
- G06F18/213—Feature extraction, e.g. by transforming the feature space; Summarisation; Mappings, e.g. subspace methods
- G06F18/2135—Feature extraction, e.g. by transforming the feature space; Summarisation; Mappings, e.g. subspace methods based on approximation criteria, e.g. principal component analysis
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/25—Fusion techniques
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
- Y02A40/00—Adaptation technologies in agriculture, forestry, livestock or agroalimentary production
- Y02A40/70—Adaptation technologies in agriculture, forestry, livestock or agroalimentary production in livestock or poultry
Landscapes
- Engineering & Computer Science (AREA)
- Data Mining & Analysis (AREA)
- Theoretical Computer Science (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Artificial Intelligence (AREA)
- Bioinformatics & Computational Biology (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Evolutionary Biology (AREA)
- Evolutionary Computation (AREA)
- Physics & Mathematics (AREA)
- General Engineering & Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Image Processing (AREA)
Abstract
本申请公开了一种草场围栏识别方法、装置及存储介质,该方法包括:获取目标区域内预定时间的高时间分辨率和高空间分辨率的卫星影像;对所述高时间分辨率和高空间分辨率的卫星影像进行预处理;计算预处理后影像的归一化植被指数NDVI,对得到的所述高时间分辨率和高空间分辨率的卫星影像的NDVI影像数据进行融合,得到高时空分辨率NDVI时间序列影像数据;对所述高时空分辨率NDVI时间序列影像数据进行特征提取;基于提取的特征,识别所述目标区域的草场围栏。本申请通过影像融合得到更加全面的高时空分辨率的遥感影像数据,能够实现大范围内草场围栏的准确识别、提取,兼顾高精度和大范围识别,且提取方法自动化程度高。
Description
技术领域
本申请涉及遥感影像处理技术领域,尤其涉及一种草场围栏识别方法、装置及存储介质。
背景技术
围栏作为积极有效的管理工具,其生态影响是广泛而复杂的,具有两面性,其短期功效可能优于长期后果,因此“围栏生态学”受到越来越多的关注。在我国,草原广泛分布于西北地区和青藏高原地区,普遍采用划区轮牧制度进行放牧。草场围栏是划区轮牧制度下运用围栏轮牧技术的一种产物,旨在保护草原,减轻一系列生态恶化问题。因此,草场围栏在我国分布典型,成为当地人类活动影响草原生态的主要形式,准确可靠的草场围栏空间识别和定量化分析,对于了解草场围栏分布范围和状况,分析围栏生态效应,实现草原地区的可持续发展具有重要的现实意义。
然而,大规模的围栏分布量化数据不足,成为当前研究围栏与生态***相互作用的局限。目前,直接从光学遥感影像中识别围栏仍然是一项艰巨的工作。围栏是一种线性基础设施,一般由编结/刺丝网和固定立柱构成,编结/刺丝网的直径不到10cm,固定立柱直径不到20cm。目前,主要的围栏识别方法可分为目视解译法和密度建模法。目视解译方法,使用非常高分辨率(VHR)遥感影像,通过目视解译手动绘制围栏,是目前识别围栏的一种主要方法。但这种方法既繁琐又耗时,同时VHR影像往往覆盖范围小,时间分辨率低,价格昂贵,并不能满足大范围围栏提取的实际需求。密度建模法,通过结合多种空间数据,构建围栏密度模型,从而确定某区域的围栏分布密度。这种密度建模方法提高了识别效率,但其局限性在于不能在空间上明确围栏的具***置,只能估算相对围栏密度,无法评估围栏对生态***的影响机制。
因此,需要新的技术以至少部分解决现有技术中存在的不能实现大范围围栏精准识别的局限。
发明内容
鉴于现有技术中存在的上述技术问题,本申请实施例提供一种草场围栏识别方法、装置及存储介质,能够解决现有可用的大范围围栏空间数据不足,使得现有的围栏识别方法不能兼顾高精度和大范围识别的技术问题。
为了解决上述技术问题,本申请实施例采用了如下技术方案:
一种草场围栏识别方法,包括:
获取目标区域内预定时间的高时间分辨率和高空间分辨率的卫星影像;
对所述高时间分辨率和高空间分辨率的卫星影像进行预处理;
计算预处理后影像的归一化植被指数NDVI,对得到的所述高时间分辨率和高空间分辨率的卫星影像的NDVI影像数据进行融合,得到高时空分辨率NDVI时间序列影像数据;
对所述高时空分辨率NDVI时间序列影像数据进行特征提取;
基于提取的特征,识别所述目标区域的草场围栏。
在一些实施例中,所述高时间分辨率的卫星影像为MODIS影像,所述高空间分辨率的卫星影像为Sentinel-2影像。
在一些实施例中,对所述高时间分辨率和高空间分辨率的卫星影像进行预处理,包括:
对所述Sentinel-2影像进行大气校正;
对所述MODIS影像进行重投影、重采样和波段提取;
使用三次修复算法对所述MODIS影像进行云像素修复。
在一些实施例中,计算预处理后影像的归一化植被指数NDVI,对得到的所述高时间分辨率和高空间分辨率的卫星影像的NDVI影像数据进行融合,得到高时空分辨率NDVI时间序列影像数据,包括:
计算预处理后影像的归一化植被指数NDVI;
选取前后两个日期的基于MODIS计算的NDVI影像数据和基于Sentinel-2计算的NDVI影像数据作为参考图像,将所述两个日期之间的基于MODIS计算的NDVI影像数据作为待融合图像,通过时空数据融合算法进行图像融合,得到所述两个日期之间的高时空分辨率NDVI影像数据。
在一些实施例中,所述归一化植被指数NDVI的计算公式为:
其中,NIR为近红外波段反射率,R为红光波段反射率;对于MODIS影像,NIR为所使用的MOD09GQ数据的b02波段,R为b01波段;对于Sentinel-2影像,NIR为B08波段,R为B04波段。
在一些实施例中,对所述高时空分辨率NDVI时间序列影像数据进行特征提取,包括:
基于主成分分析方法,对所述NDVI时间序列影像数据进行特征提取,得到主成分分析结果影像。
在一些实施例中,基于主成分分析方法,对所述NDVI时间序列影像数据进行特征提取,得到主成分分析结果影像,包括:
对所述NDVI时间序列影像数据进行波段合成,得到多波段影像;
对所述多波段影像进行主成分分析,得到主成分多波段影像。
在一些实施例中,基于提取的特征,识别所述目标区域的草场围栏,包括:
对所述主成分分析结果影像进行影像分割,得到初始的草场围栏识别结果;
对所述初始的草场围栏识别结果,依次进行线简化、土地覆盖掩膜以及利用形状规则性指数进行草场围栏提取,得到最终的草场围栏识别结果。
本申请实施例还提供一种草场围栏识别装置,包括:
获取模块,配置为获取目标区域内预定时间的高时间分辨率和高空间分辨率的卫星影像;
预处理模块,配置为对所述高时间分辨率和高空间分辨率的卫星影像进行预处理;
影像融合模块,配置为计算预处理后影像的归一化植被指数NDVI,对得到的所述高时间分辨率和高空间分辨率的卫星影像的NDVI影像数据进行融合,得到高时空分辨率NDVI时间序列影像数据;
特征提取模块,配置为对所述高时空分辨率NDVI时间序列影像数据进行特征提取;
识别模块,配置为基于提取的特征,识别所述目标区域的草场围栏。
本申请实施例还提供一种计算机可读存储介质,其上存储有计算机可执行指令,所述计算机可执行指令由处理器执行时,实现上述的草场围栏识别方法的步骤。
与现有技术相比,本申请实施例提供的草场围栏识别草场围栏识别方法、装置及存储介质,通过对现有的高时间分辨率和高空间分辨率的遥感卫星影像进行融合,得到更加全面的高时空分辨率的遥感影像数据,能够实现大范围内草场围栏的准确识别、提取,草场围栏识别兼顾高精度和大范围识别,且提取方法自动化程度高。
附图说明
图1为本申请实施例的草场围栏识别方法的流程图;
图2为本申请实施例的目标区域的地理位置和Sentinel-2影像;
图3为本申请实施例中仅使用Sentinel-2影像构建的NDVI时间序列曲线与使用本申请的方法构建的NDVI时间序列曲线的对比图。
图4为本申请实施例中目标区域的最终围栏识别结果示意图;
图5为本申请实施例的草场围栏识别装置的结构示意图。
具体实施方式
为使本领域技术人员更好的理解本申请的技术方案,下面结合附图和具体实施方式对本申请作详细说明。
应理解的是,可以对此处公开的实施例做出各种修改。因此,上述说明书不应该视为限制,而仅是作为实施例的范例。本领域的技术人员将想到在本申请的范围和精神内的其他修改。
包含在说明书中并构成说明书的一部分的附图示出了本申请的实施例,并且与上面给出的对本申请的大致描述以及下面给出的对实施例的详细描述一起用于解释本申请的原理。
通过下面参照附图对给定为非限制性实例的实施例的优选形式的描述,本申请的这些和其它特性将会变得显而易见。
还应当理解,尽管已经参照一些具体实例对本申请进行了描述,但本领域技术人员能够确定地实现本申请的很多其它等效形式,它们具有如权利要求所述的特征并因此都位于借此所限定的保护范围内。
当结合附图时,鉴于以下详细说明,本申请的上述和其他方面、特征和优势将变得更为显而易见。
此后参照附图描述本申请的具体实施例;然而,应当理解,所公开的实施例仅仅是本申请的实例,其可采用多种方式实施。熟知和/或重复的功能和结构并未详细描述以避免不必要或多余的细节使得本申请模糊不清。因此,本文所公开的具体的结构性和功能性细节并非意在限定,而是仅仅作为权利要求的基础和代表性基础用于教导本领域技术人员以实质上任意合适的详细结构多样地使用本申请。
由于采用围栏轮牧技术进行放牧的草原地区,植被的生长状况和覆盖率在不同的草场以及同一草场的不同时间具有差异性。具体而言,空间上,不同草场在同一日期存在植被生长状况的差异,时间上,不同时间同一草场的植被生长状况存在差异。因此,可以通过增强这种差异性实现草场边界的识别,进而提取出目标区域的草场围栏。
图1为本申请实施例提供的草场围栏识别方法的流程图。如图1所示,本申请实施例提供的草场围栏识别方法包括如下步骤。
S1:获取目标区域内预定时间的高时间分辨率和高空间分辨率的卫星影像。
具体地,以内蒙古自治区锡林郭勒草原为目标区域,对该草原的围栏进行识别。如图2所示,锡林郭勒草原位于内蒙古自治区锡林郭勒盟境内,介于43.25°N-44.25°N,115.76°E-117.12°E之间,具有优质的天然草场,草原总面积19.2万km2,占全盟总面积的97.8%。20世纪80年代以来,一系列政策极大地推动了锡林郭勒草原的围栏建设,据2019年不完全统计,锡林郭勒盟草场围栏面积1.9亿亩,占草原总面积的65%。
本实施例中,预定时间为草地生长期,通常为每年的5-9月,因此,可以通过遥感卫星采集待识别的目标区域处于草地生长期的卫星影像,其中,高时间分辨率卫星影像为MODIS影像,高空间分辨率影像为Sentinel-2影像。
通过遥感卫星获取的MODIS影像为包含每个日期的连续影像;而由于无云Sentinel-2影像并不是每个日期均出现,因此,获取的无云Sentinel-2影像的数量有限,并不是处于草地生长期的每个日期都有。
本实施例中,如表1所示,通过遥感卫星获取2018年草地生长期内(5-9月)覆盖目标区域的5幅无云Sentinel-2影像和所有日期的MODIS影像。Sentinel-2影像为L1C级大气顶部发射率数据;MODIS影像包含MOD09GA和MOD09GQ数据,为L2G地表反射率数据。MOD09GA数据中包含云波段(即state波段),可以用来识别MOD09GQ数据中的云。
表1获取的高时间分辨率和高空间分辨率的卫星影像及其参数
如表1所示,MODIS影像数据(高时间分辨率卫星影像)的时间分辨率为1天,Sentinel-2影像数据(高空间分辨率卫星影像)的空间分辨率为10-60m。
S2:对所述高时间分辨率(MODIS影像)和高空间分辨率(Sentinel-2影像)的卫星影像进行预处理。
具体地,通过预处理修复地表反射率影像中被云覆盖的像元,得到真实反射率影像。
不同影像的预处理方式不同,具体地:
S21:对于选取的Sentinel-2影像,使用Sen2Cor工具,对L1C级数据进行大气校正,获得地表反射率数据。
S22:对于获取的MODIS影像,借助MRT工具,分别对MOD09GQ和MOD09GA数据进行重投影、重采样和波段提取。
其中,重投影为批量重投影,批量重投影具体包括:将MODIS影像(即MOD09GA和MOD09GQ数据)原始的正弦投影转换为UTM zone 49N投影,以与Sentinel-2影像使用的投影保持一致。
重采样具体包括:将MOD09GA数据原始的空间分辨率1000m重采样至250m,以与MOD09GQ数据的空间分辨率保持一致。如表1所示,MOD09GA数据的原始空间分辨率是1000m,MOD09GQ数据的原始空间分辨率是250m,为了方便后续的云识别、云修复等操作,需要将空间分辨率统一,因此,将MOD09GA数据的原始空间分辨率1000重采样至250m。
波段提取具体包括:MOD09GQ数据提取的波段为红光b01波段与近红外b02波段;MOD09GA数据提取的波段为state波段,以便进行云识别,从而进行云修复。
在一些实施例中,预处理还包括格式转换,即将Sentinel-2影像和MODIS影像的格式统一,以便后续进行影像数据处理。
进一步地,对MODIS影像的预处理还包括:对MODIS影像进行云像素修复。
具体地,可以使用从MOD09GA数据中提取的state波段进行MODIS影像云识别,并基于GNSPI算法对MOD09GQ影像中相应区域的云缺失像元进行修复(即使用三次修复算法对MODIS影像进行云像素修复)。
通过云识别可以找到影像中受到云干扰的像素,但是,这些像素被云覆盖,不能反映地面真实的情况,是无效像素。本实施例中,先利用MOD09GA数据中的state波段对MOD09GQ数据中受到云干扰的无效像素进行识别标记;然后,通过三次修复法对MOD09GQ中被标记的受到云干扰的无效像素进行修复,还原地面真实的情况(即没被云覆盖时的真实影像)。
三次修复法可以对影像中被云污染(即被云盖住的)的像素进行修复,还原地面真实的情况。该方法使用待修复影像给定时间窗口内的有效像素实现修复,可以提高有效观测像素的利用率。三次修复法利用被污染像素周围(一定空间窗口)邻近日期(即前后一定时间窗口)的有效参考像素对它进行修复。具体包括如下三次修复:
(1)第一次修复:针对待修复像素所在日期前后都可以找到参考像素的情况,对云污染像素进行修复;(2)第二次修复:针对第一次修复中没被修复到的云污染像素,而且待修复日期前可以找到参考像素的情况,对云污染像素进行修复;(3)第三次修复:针对第二次修复中还没被修复到的云污染像素,而且待修复日期后可以找到参考像素的情况,对云污染像素进行修复。其中,空间窗口固定为待修复位置周围25像素(待修复位置为中心的5X5像素矩形范围);前后一定时间窗口固定为以待修复日期为中心的前后各48天。
通过上述云识别与云修复可以对上述MOD09GQ数据中的受到云干扰的无效像素进行修复,提高MOD09GQ数据中有效像素的利用率。
S3:计算预处理后影像的归一化植被指数NDVI,对得到的所述高时间分辨率和高空间分辨率的卫星影像的NDVI影像数据进行融合,得到高时空分辨率NDVI时间序列影像数据。
步骤S3具体包括:
S31:计算预处理后影像的归一化植被指数NDVI。
归一化植被指数NDVI的计算公式为:
式(1)中,NIR为近红外波段反射率,R为红光波段反射率。对于MODIS影像,NIR为所使用的MOD09GQ数据的b02波段(近红外波段),R为b01波段(红光波段);对于Sentinel-2影像,NIR为B08波段,R为B04波段。本实施例中,可以使用ENVI 5.3软件中的批量波段运算扩展工具完成上述NDVI计算。
归一化植被指数NDVI是反映植被变化的有效指标,可以用它来区分不同状态的植被。因此,基于获取的遥感影像计算NDVI,得到NDVI影像,影像中不同的NDVI值可以反映植被不同的生长情况,采用各天的NDVI构建的时间序列可以反映植被从返青,生长到枯黄的过程。
NDVI作为植被生长状态和覆盖度的最佳指示因子可以增强上述草场在时间和空间上的差异。因此,通过获取高时空分辨率的NDVI时间序列影像数据可以增强草场内部的同质性和草场之间的异质性,从而便于草场围栏的识别。
S32:选取前后两个日期的基于MODIS(具体为MOD09GQ)计算的NDVI影像数据(空间分辨率250m)和基于Sentinel-2计算的NDVI影像数据(空间分辨率10m)作为参考图像,将所述两个日期之间的基于MODIS计算的NDVI影像数据作为待融合图像,通过时空数据融合算法(ESTARFM)进行图像融合,得到两个日期之间的高时空分辨率(每天,且空间分辨率为10m)NDVI影像数据。
使用时空数据融合算法进行影像融合时,将基于Sentinel-2影像计算得到的NDVI影像数据与基于MODIS影像计算得到的NDVI影像数据两种影像作为输入,因此,输入的影像需要具有一致的空间分辨率和空间范围。
本实施例中,对MODIS影像计算得到的NDVI进行批量重采样到10m分辨率,以与Sentinel-2影像具有相同的空间分辨率。
上述的空间范围指待识别的目标区域的范围,可以根据实际需要自行设定,本实施例中,目标区域的范围即为Sentinel-2影像所覆盖的范围,因此,使用Sentinel-2影像的空间范围对重采样结果进行裁剪,以使Sentinel-2影像和MOD09GQ影像具有相同的空间范围。
由于本申请实施例是基于草场植被生长期的NDVI时间序列找到不同草场之间的差异,从而识别草场围栏,因此,时间序列越完整,草场之间的差异会越明显。时间序列完整的前提是需要每天都有可用的高空间分辨率的影像。Sentinel-2影像数据是具有高空间分辨率的数据,但是受到重访周期以及天气的限制,导致研究区草场植被生长期之间的可用影像(干净无云)有且仅有数个(例如本实施例中为5幅影像)。为了获取这5幅影像之间的每一天的可用影像,采用了ESTARFM算法,利用MODIS影像的高时间分辨率特征(每天),对相应日期的Sentinel-2影像进行预测,从而得到一定连续时间范围内的高空间分辨率(空间分辨率10m)的Sentinel-2影像,即具有高时空分辨率的Sentinel-2影像。
本实施例中,将上述5幅基于Sentinel-2的NDVI影像和对应日期基于MODIS的NDVI影像为参考数据,参考日期之间基于MODIS的NDVI影像为待融合数据,进行影像融合。即使用表1中的2018DOY 118、2018DOY 128、2018DOY 248、2018DOY 253和2018DOY 258这5个日期的Sentinel-2影像和MODIS影像得到的NDVI影像预测2018DOY118-2018DOY 258这一完整时间序列的Sentinel-2影像。
具体地,可以以上述5幅Sentinel-2影像日期为界限两两划分时间范围,第一个融合的时间范围是2018DOY 118-2018DOY 128,第二个融合的时间范围是2018DOY 128-2018DOY 248,以此类推,一共确定出4个需要融合的时间范围。将上述每个时间范围内的前后两个日期的Sentinel-2数据与对应的该天的MODIS数据作为参考影像输入,另外还需输入前后两个日期之间的所有的可用的MODIS数据作为待融合影像,最终得到该完整的时间范围内的Sentinel-2预测影像,增强了时间序列的完整性。
可以理解的是,本实施例中,可以对每个时间范围的影像数据进行融合后,再将个时间范围的影像进行合并,得到2018DOY 118-2018DOY 258这一完整时间序列的Sentinel-2影像。另一些实施例中,也可以直接选取2018DOY 118、2018DOY 258这两个日期的MODIS数据和Sentinel-2数据作为参考数据,以及该时间范围内的MODIS数据作为待融合数据,通过一次融合得到2018DOY 118-2018DOY 258这一完整时间序列的Sentinel-2影像。
在将数据进行融合后,对融合结果进行筛选,构建高时空分辨率NDVI时间序列影像数据,如图3所示,与仅使用研究期(上述表1中5个日期)Sentinel-2影像构建的NDVI时间序列影像数据相比,采用本实施例的融合结果构建的NDVI时间序列影像数据具有更好的完整性。
可以理解的是,由于上述用于预测的参考数据为基于Sentinel-2影像计算得到的NDVI影像数据和基于MODIS影像计算得到的NDVI影像数据,因此,该预测出的Sentinel-2影像为基于Sentinel-2影像计算得到的NDVI影像数据。
S4:对所述高时空分辨率NDVI时间序列影像数据进行特征提取。
具体地,基于主成分分析(PCA)方法,对步骤S3中得到的NDVI时间序列影像进行特征提取,得到主成分分析结果影像。
具体包括如下步骤:
S41:对所述NDVI时间序列影像数据进行波段合成,得到多波段影像。
所述NDVI时间序列影像的每个日期包含一个影像,即为单波段影像,因此,可以将获取的NDVI时间序列影像进行波段合成,得到多波段影像。
S42:对所述多波段影像进行主成分分析,得到主成分多波段影像。
使用ENVI中的工具进行主成分分析,进行主成分分析包含两种计算方式,一种是基于波段之间相关系数矩阵,一种是基于波段之间的协方差矩阵。本实施例中,基于相关系数矩阵进行主成分分析。主成分分析会对多波段像中的有效信息进行整合,生成主成分多波段影像,波段越靠前包含的信息量越大。即通过主成分分析可以提取影像中的主要特征。得到的主成分多波段影像即为PCA分析结果影像。
具体实施中,可以借助ENVI 5.3软件的Forward PCA Rotation New Statisticsand Rotate工具,将步骤S3得到的NDVI时间序列影像作为输入,计算方式选择相关系数矩阵,得到输出结果,本实施例中,输出结果保留信息量较大的前五个主成分波段,即PC1至PC5。具体实施中,可以根据实际需要确定需要保留的波段,本申请不具体限定。
S5:基于提取的特征,识别所述目标区域的草场围栏。
步骤S5具体包括如下步骤:
S51:对所述主成分分析结果影像进行影像分割,得到初始的草场围栏识别结果。
利用步骤S4中的主成分分析结果影像,通过面向对象影像分割技术,基于多尺度分割MRS方法对所述主成分分析结果影像进行影像分割,得到初始的草场围栏识别结果。
具体地,可以利用eCognition Developer 9中的Multiresolution Segmentation工具,将步骤S4中得到的主成分多波段影像作为输入,通过设置尺度、形状和紧致度等参数实现多尺度分割。多尺度分割MRS方法涉及的上述三个关键参数(尺度、形状和紧致度),可以使用eCognition Developer 9中的ESP2工具自动确定。本实施例中,将形状参数设置为0.01,紧致度参数设置为0.5,尺度参数设置为12。
S52:对所述初始的草场围栏识别结果,依次进行线简化、土地覆盖掩膜以及利用形状规则性指数进行草场围栏提取,得到最终的草场围栏识别结果。
步骤S51通过多尺度分割得到的分割结果是基于像素边界绘制的带有锯齿状的边界线,而现实中的围栏应该是平滑的,因此,通过线简化算法(DP算法),有助于将原始的草场边界变得更加平滑,更接近真实的草场边界。
本实施例中,使用QGIS 3.10中的简化工具,将距离阈值确定为30m,对步骤S51中得到的分割结果进行线简化。
土地覆盖掩膜处理包括:
通过使用步骤S3中得到的高时空分辨率NDVI时间序列影像数据,基于K-均值聚类方法进行土地覆盖分类,实现对水体和城镇建筑的掩膜处理。
具体地,通过叠加分割对象和NDVI时间序列影像数据,获取每个分割对象的NDVI时间序列曲线。根据分割对象的NDVI时间序列曲线上表现的不同土地覆盖类型物候特征,利用SPSS 26.0中的K-均值聚类工具,实现土地覆盖分类,土地覆盖分类的具体类别主要有水体、不透水面、沙地、草地等。在得到土地覆盖分类后,通过掩膜剔除分割对象中的水体、城镇建筑和部分沙地等其他地类,保留草地。
进一步地,基于掩膜后的分割对象,计算形状规则性指数S,计算公式为:
S=W1×Nmax+W2×Lmax+W3×Dprei+W4×Darea (2)
式(2)中,Nmax为最大凹口比,其中,AOi为分割对象的凹口数量,AOmax为目标区域内所有分割对象的凹口数量最大值;Lmax为分割对象的最大边长;Dprei为周长比,其中,preii为分割对象的周长,preih为分割对象外包多边形的周长;Darea为面积比,其中,areai为分割对象的面积,areac为分割对象外包圆的面积;W1至W4分别为上述指标的权重,可以通过对上述指标的主成分分析进行确定。
具体实施中,最大凹口比、最大边长、周长比和面积比可以结合ArcGIS 10.5软件和Python 3.7计算得到,主成分分析利用SPSS 26.0中的因子分析工具实现。权重可以通过运用SPSS中的主成分分析(PCA)工具,输入上面计算得到的Nmax、Lmax、Dprei、Darea四个指标,进行PCA分析,取各指标对应的PC1特征值确定的。
形状规则性指数的指标按地块进行统计,利用不同地块计算出的形状规则性指数构建形状规则性指数的统计表,其中,每个地块相当于表中的一行,各个指标相当于表的列,将该统计表输入SPSS中进行主成分分析,可以得到四个主成分,即PC1到PC4。基于分析结果构建分析结果表,每个主成分的特征值相当于表中的一列,各个指标相当于表中的行,因此,每个主成分与各指标通过特征值进行对应。PC1所包含的信息量最大,所以把各指标在PC1的特征值作为各指标的权重。
根据掩膜后分割对象的形状规则性指数S计算结果,通过判断形状规则性指数S是否满足预设阈值以剔除非草场对象,从而获得如图4所示的最终的草场围栏识别结果。
因为即使在草地内部,得到的边界中也存在错误的非草场边界,因此,本实施例中,在经过DP线简化(30m距离阈值)和土地覆盖掩膜(保留草地)后,通过经计算形状规则性指数S可以有效剔除非草场边界。
形状规则性指数S的取值范围通常为小于0,具体到每个目标区域,需要根据实际情况结合S的均值和标准差进行判定。本实施例中,当确定S满足S≤S+0.25×σ时,即可确定为草场边界。其中,和分别为σ为整个目标区域内所有形状规则性指数的均值和标准差。
为验证本申请草场围栏识别结果的准确性,本申请实施例通过比较采用本申请的方法与单一时相的Sentinel-2NDVI影像得到的围栏识别结果,并使用总体精度、漏分误差和错分误差三个精度评价指标对识别结果的准确性进行精度评价。其中,采用单一时相的Sentinel-2NDVI影像进行围栏识别是指仅仅用某一景(一天)的Sentinel-2影像进行识别。总体精度为验证区域内被正确提取的围栏长度占样本围栏长度的百分比;漏分误差为未被提取的围栏长度占样本围栏长度的百分比;错分误差为过度提取的围栏长度占样本围栏长度的百分比。
根据实地调查和ArcGIS 10.5在线地图中的非常高分辨率遥感影像,可以选取10个10×10km的验证区域,通过目视解译围栏手动绘制验证样本,精度评价数据如表2所示。
表2围栏识别结果的精度评价数据
输入影像数据 | 总体精度(OA) | 漏分误差(OE) | 错分(CE) |
NDVI时间序列影像 | 80.06% | 20.88% | 19.01% |
单一时相NDVI影像 | 68.62% | 39.13% | 23.61% |
从表2中可以看出,采用本申请的方法得到的围栏识别精度显著优于采用单一时相的Sentinel-2NDVI影像。
本申请实施例提供的草场围栏识别方法能够通过对现有的高时间分辨率和高空间分辨率的遥感卫星影像进行融合,得到更加全面的高时空分辨率的遥感影像数据,能够实现大范围内草场围栏的准确识别、提取,草场围栏识别兼顾高精度和大范围识别,且提取方法自动化程度高。
图5为本申请实施例提供的草场围栏识别装置的流程图。如图5所示,本申请实施例还提供了一种草场围栏识别装置,包括:
获取模块100,配置为获取目标区域内预定时间的高时间分辨率和高空间分辨率的卫星影像;
预处理模块200,配置为对所述高时间分辨率和高空间分辨率的卫星影像进行预处理;
影像融合模块300,配置为计算预处理后影像的归一化植被指数NDVI,对得到的所述高时间分辨率和高空间分辨率的卫星影像的NDVI影像数据进行融合,得到高时空分辨率NDVI时间序列影像数据;
特征提取模块400,配置为对所述高时空分辨率NDVI时间序列影像数据进行特征提取;
识别模块500,配置为基于提取的特征,识别所述目标区域的草场围栏。
本申请实施例提供的草场围栏识别装置对应于上述实施例的草场围栏识别方法,草场围栏识别方法实施例中的任何可选项也适用于本实施例,这里不再详述。
本申请实施例还提供了一种计算机可读存储介质,其上存储有计算机可执行指令,所述计算机可执行指令由处理器执行时,实现上述根据本申请的实施例中的草场围栏识别方法。
在一些实施例中,执行算机可执行指令处理器可以是包括一个以上通用处理设备的处理设备,诸如微处理器、中央处理单元(CPU)、图形处理单元(GPU)等。更具体地,该处理器可以是复杂指令集计算(CISC)微处理器、精简指令集计算(RISC)微处理器、超长指令字(VLIW)微处理器、运行其他指令集的处理器或运行指令集的组合的处理器。该处理器还可以是一个以上专用处理设备,诸如专用集成电路(ASIC)、现场可编程门阵列(FPGA)、数字信号处理器(DSP)、片上***(SoC)等。
在一些实施例中,计算机可读存储介质可以为存储器,诸如只读存储器(ROM)、随机存取存储器(RAM)、相变随机存取存储器(PRAM)、静态随机存取存储器(SRAM)、动态随机存取存储器(DRAM)、电可擦除可编程只读存储器(EEPROM)、其他类型的随机存取存储器(RAM)、闪存盘或其他形式的闪存、缓存、寄存器、静态存储器、光盘只读存储器(CD-ROM)、数字通用光盘(DVD)或其他光学存储器、盒式磁带或其他磁存储设备,或被用于储存能够被计算机设备访问的信息或指令的任何其他可能的非暂时性的介质等。
以上实施例仅为本申请的示例性实施例,不用于限制本申请,本申请的保护范围由权利要求书限定。本领域技术人员可以在本申请的实质和保护范围内,对本申请做出各种修改或等同替换,这种修改或等同替换也应视为落在本申请的保护范围内。
Claims (10)
1.一种草场围栏识别方法,其特征在于,包括:
获取目标区域内预定时间的高时间分辨率和高空间分辨率的卫星影像;
对所述高时间分辨率和高空间分辨率的卫星影像进行预处理;
计算预处理后影像的归一化植被指数NDVI,对得到的所述高时间分辨率和高空间分辨率的卫星影像的NDVI影像数据进行融合,得到高时空分辨率NDVI时间序列影像数据;
对所述高时空分辨率NDVI时间序列影像数据进行特征提取;
基于提取的特征,识别所述目标区域的草场围栏。
2.根据权利要求1所述的草场围栏识别方法,其特征在于,所述高时间分辨率的卫星影像为MODIS影像,所述高空间分辨率的卫星影像为Sentinel-2影像。
3.根据权利要求2所述的草场围栏识别方法,其特征在于,对所述高时间分辨率和高空间分辨率的卫星影像进行预处理,包括:
对Sentinel-2影像进行大气校正;
对MODIS影像进行重投影、重采样和波段提取;
使用三次修复算法对所述MODIS影像进行云像素修复。
4.根据权利要求2所述的草场围栏识别方法,其特征在于,计算预处理后影像的归一化植被指数NDVI,对得到的所述高时间分辨率和高空间分辨率的卫星影像的NDVI影像数据进行融合,得到高时空分辨率NDVI时间序列影像数据,包括:
计算预处理后影像的归一化植被指数NDVI;
选取前后两个日期的基于MODIS计算的NDVI影像数据和基于Sentinel-2计算的NDVI影像数据作为参考图像,将所述两个日期之间的基于MODIS计算的NDVI影像数据作为待融合图像,通过时空数据融合算法进行图像融合,得到所述两个日期之间的高时空分辨率NDVI影像数据。
6.根据权利要求1所述的草场围栏识别方法,其特征在于,对所述高时空分辨率NDVI时间序列影像数据进行特征提取,包括:
基于主成分分析方法,对所述NDVI时间序列影像数据进行特征提取,得到主成分分析结果影像。
7.根据权利要求6所述的草场围栏识别方法,其特征在于,基于主成分分析方法,对所述NDVI时间序列影像数据进行特征提取,得到主成分分析结果影像,包括:
对所述NDVI时间序列影像数据进行波段合成,得到多波段影像;
对所述多波段影像进行主成分分析,得到主成分多波段影像。
8.根据权利要求6所述的草场围栏识别方法,其特征在于,基于提取的特征,识别所述目标区域的草场围栏,包括:
对所述主成分分析结果影像进行影像分割,得到初始的草场围栏识别结果;
对所述初始的草场围栏识别结果,依次进行线简化、土地覆盖掩膜以及利用形状规则性指数进行草场围栏提取,得到最终的草场围栏识别结果。
9.一种草场围栏识别装置,其特征在于,包括:
获取模块,配置为获取目标区域内预定时间的高时间分辨率和高空间分辨率的卫星影像;
预处理模块,配置为对所述高时间分辨率和高空间分辨率的卫星影像进行预处理;
影像融合模块,配置为计算预处理后影像的归一化植被指数NDVI,对得到的所述高时间分辨率和高空间分辨率的卫星影像的NDVI影像数据进行融合,得到高时空分辨率NDVI时间序列影像数据;
特征提取模块,配置为对所述高时空分辨率NDVI时间序列影像数据进行特征提取;
识别模块,配置为基于提取的特征,识别所述目标区域的草场围栏。
10.一种计算机可读存储介质,其特征在于,其上存储有计算机可执行指令,所述计算机可执行指令由处理器执行时,实现根据权利要求1-8中任一项所述的草场围栏识别方法。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202110625305.1A CN113361398B (zh) | 2021-06-04 | 2021-06-04 | 草场围栏识别方法、装置及存储介质 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202110625305.1A CN113361398B (zh) | 2021-06-04 | 2021-06-04 | 草场围栏识别方法、装置及存储介质 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN113361398A true CN113361398A (zh) | 2021-09-07 |
CN113361398B CN113361398B (zh) | 2022-10-18 |
Family
ID=77532227
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202110625305.1A Active CN113361398B (zh) | 2021-06-04 | 2021-06-04 | 草场围栏识别方法、装置及存储介质 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN113361398B (zh) |
Citations (11)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN101930548A (zh) * | 2010-06-24 | 2010-12-29 | 北京师范大学 | 一种基于独立成分分析算法的高空间分辨率影像的识别方法 |
CN102609726A (zh) * | 2012-02-24 | 2012-07-25 | 中国科学院东北地理与农业生态研究所 | 利用面向对象技术融合高空间和高时间分辨率数据的遥感图像分类方法 |
CN103886334A (zh) * | 2014-04-08 | 2014-06-25 | 河海大学 | 一种多指标融合的高光谱遥感影像降维方法 |
CN103914678A (zh) * | 2013-01-05 | 2014-07-09 | 中国科学院遥感与数字地球研究所 | 基于纹理与植被指数的撂荒地遥感识别方法 |
CN104915674A (zh) * | 2014-10-24 | 2015-09-16 | 北京师范大学 | Landsat8和MODIS融合构建高时空分辨率数据识别秋粮作物的方法 |
CN106295498A (zh) * | 2016-07-20 | 2017-01-04 | 湖南大学 | 光学遥感图像目标区域检测装置与方法 |
CN108170926A (zh) * | 2017-12-12 | 2018-06-15 | 伊犁师范学院 | 一种河谷草地退化情况的信息数据采集及分析方法 |
US20190034725A1 (en) * | 2016-01-29 | 2019-01-31 | Global Surface Intelligence Limited | System and method for earth observation and analysis |
CN110008908A (zh) * | 2019-04-10 | 2019-07-12 | 中国科学院、水利部成都山地灾害与环境研究所 | 一种基于高分遥感影像的草原围栏提取方法 |
CN110879992A (zh) * | 2019-11-27 | 2020-03-13 | 内蒙古工业大学 | 基于迁移学习的草原地表覆盖物分类方法和*** |
US20210142447A1 (en) * | 2019-07-01 | 2021-05-13 | David P. Groeneveld | Method to Correct Satellite Data to Surface Reflectance Using Scene Statistics |
-
2021
- 2021-06-04 CN CN202110625305.1A patent/CN113361398B/zh active Active
Patent Citations (11)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN101930548A (zh) * | 2010-06-24 | 2010-12-29 | 北京师范大学 | 一种基于独立成分分析算法的高空间分辨率影像的识别方法 |
CN102609726A (zh) * | 2012-02-24 | 2012-07-25 | 中国科学院东北地理与农业生态研究所 | 利用面向对象技术融合高空间和高时间分辨率数据的遥感图像分类方法 |
CN103914678A (zh) * | 2013-01-05 | 2014-07-09 | 中国科学院遥感与数字地球研究所 | 基于纹理与植被指数的撂荒地遥感识别方法 |
CN103886334A (zh) * | 2014-04-08 | 2014-06-25 | 河海大学 | 一种多指标融合的高光谱遥感影像降维方法 |
CN104915674A (zh) * | 2014-10-24 | 2015-09-16 | 北京师范大学 | Landsat8和MODIS融合构建高时空分辨率数据识别秋粮作物的方法 |
US20190034725A1 (en) * | 2016-01-29 | 2019-01-31 | Global Surface Intelligence Limited | System and method for earth observation and analysis |
CN106295498A (zh) * | 2016-07-20 | 2017-01-04 | 湖南大学 | 光学遥感图像目标区域检测装置与方法 |
CN108170926A (zh) * | 2017-12-12 | 2018-06-15 | 伊犁师范学院 | 一种河谷草地退化情况的信息数据采集及分析方法 |
CN110008908A (zh) * | 2019-04-10 | 2019-07-12 | 中国科学院、水利部成都山地灾害与环境研究所 | 一种基于高分遥感影像的草原围栏提取方法 |
US20210142447A1 (en) * | 2019-07-01 | 2021-05-13 | David P. Groeneveld | Method to Correct Satellite Data to Surface Reflectance Using Scene Statistics |
CN110879992A (zh) * | 2019-11-27 | 2020-03-13 | 内蒙古工业大学 | 基于迁移学习的草原地表覆盖物分类方法和*** |
Non-Patent Citations (2)
Title |
---|
郭利彪等: ""基于数据机理的植被叶面积指数遥感反演研究"", 《遥感技术与应用》 * |
鲍甜甜等: ""GOOGLE EARTH 在野外地质工作中的使用方法"", 《中国水运》 * |
Also Published As
Publication number | Publication date |
---|---|
CN113361398B (zh) | 2022-10-18 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Ni et al. | An enhanced pixel-based phenological feature for accurate paddy rice mapping with Sentinel-2 imagery in Google Earth Engine | |
Li et al. | SinoLC-1: the first 1-meter resolution national-scale land-cover map of China created with the deep learning framework and open-access data | |
Wang et al. | Graph-based block-level urban change detection using Sentinel-2 time series | |
González-Yebra et al. | Methodological proposal to assess plastic greenhouses land cover change from the combination of archival aerial orthoimages and Landsat data | |
CN110991430B (zh) | 基于遥感图像的地物识别及覆盖率计算方法及*** | |
Szantoi et al. | A tool for rapid post-hurricane urban tree debris estimates using high resolution aerial imagery | |
Mtibaa et al. | Land cover mapping in cropland dominated area using information on vegetation phenology and multi-seasonal Landsat 8 images | |
Wang et al. | Deep segmentation and classification of complex crops using multi-feature satellite imagery | |
Zhou et al. | Deep learning-based local climate zone classification using Sentinel-1 SAR and Sentinel-2 multispectral imagery | |
CN113780307A (zh) | 一种区域年度最大蓝绿空间信息提取方法 | |
Vigneshwaran et al. | Comparison of classification methods for urban green space extraction using very high resolution worldview-3 imagery | |
CN114596489A (zh) | 一种高精度人类居住指数的多源遥感城市建成区提取方法 | |
CN110334623A (zh) | 一种基于Sentinel-2A卫星遥感影像提取崩岗信息的方法 | |
Nazmfar et al. | Classification of satellite images in assessing urban land use change using scale optimization in object-oriented processes (a case study: Ardabil city, Iran) | |
Bagan et al. | Improved subspace classification method for multispectral remote sensing image classification | |
Zhang et al. | Classification of desert grassland species based on a local-global feature enhancement network and UAV hyperspectral remote sensing | |
Shi et al. | An improved framework for assessing the impact of different urban development strategies on land cover and ecological quality changes-A case study from Nanjing Jiangbei New Area, China | |
Deng et al. | Comparison of 2D and 3D vegetation species mapping in three natural scenarios using UAV-LiDAR point clouds and improved deep learning methods | |
Feng et al. | S2EFT: Spectral-spatial-elevation fusion transformer for hyperspectral image and LiDAR classification | |
CN113361398B (zh) | 草场围栏识别方法、装置及存储介质 | |
Xie et al. | Improvement and application of UNet network for avoiding the effect of urban dense high-rise buildings and other feature shadows on water body extraction | |
Meng et al. | Large-scale and high-resolution paddy rice intensity mapping using downscaling and phenology-based algorithms on Google Earth Engine | |
CN115830441A (zh) | 一种作物识别方法、装置、***和介质 | |
Zhang et al. | A Mapping Approach for Eucalyptus Plantations Canopy and Single-Tree Using High-Resolution Satellite Images in Liuzhou, China | |
Álvarez | Monitoring Land Cover Changes using Landsat Data and Maximum Likelihood Classification: A Case Study of Hanoi, Vietnam |
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 |