CN114782324A - 一种基于地块形态特征的农田作物行向遥感识别方法 - Google Patents
一种基于地块形态特征的农田作物行向遥感识别方法 Download PDFInfo
- Publication number
- CN114782324A CN114782324A CN202210329863.8A CN202210329863A CN114782324A CN 114782324 A CN114782324 A CN 114782324A CN 202210329863 A CN202210329863 A CN 202210329863A CN 114782324 A CN114782324 A CN 114782324A
- Authority
- CN
- China
- Prior art keywords
- land
- plot
- row direction
- crop row
- image
- 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 35
- 230000000877 morphologic effect Effects 0.000 title claims abstract description 17
- 238000011160 research Methods 0.000 claims abstract description 13
- 238000009826 distribution Methods 0.000 claims description 9
- 101100356682 Caenorhabditis elegans rho-1 gene Proteins 0.000 claims 1
- 244000037666 field crops Species 0.000 claims 1
- 238000005516 engineering process Methods 0.000 abstract description 6
- 238000012545 processing Methods 0.000 abstract description 4
- 238000001514 detection method Methods 0.000 description 14
- 240000008042 Zea mays Species 0.000 description 8
- 235000002017 Zea mays subsp mays Nutrition 0.000 description 8
- 239000002689 soil Substances 0.000 description 8
- 230000009466 transformation Effects 0.000 description 8
- 238000002474 experimental method Methods 0.000 description 7
- 235000005824 Zea mays ssp. parviglumis Nutrition 0.000 description 5
- 235000005822 corn Nutrition 0.000 description 5
- 238000004364 calculation method Methods 0.000 description 4
- 230000001788 irregular Effects 0.000 description 4
- 230000003287 optical effect Effects 0.000 description 4
- 241000196324 Embryophyta Species 0.000 description 3
- 235000016383 Zea mays subsp huehuetenangensis Nutrition 0.000 description 3
- 230000009286 beneficial effect Effects 0.000 description 3
- 230000008901 benefit Effects 0.000 description 3
- 238000003708 edge detection Methods 0.000 description 3
- 238000000605 extraction Methods 0.000 description 3
- 235000009973 maize Nutrition 0.000 description 3
- 230000000737 periodic effect Effects 0.000 description 3
- 230000000007 visual effect Effects 0.000 description 3
- 235000016068 Berberis vulgaris Nutrition 0.000 description 2
- 241000335053 Beta vulgaris Species 0.000 description 2
- 244000068988 Glycine max Species 0.000 description 2
- 235000010469 Glycine max Nutrition 0.000 description 2
- 238000004458 analytical method Methods 0.000 description 2
- 230000002457 bidirectional effect Effects 0.000 description 2
- 230000008859 change Effects 0.000 description 2
- 238000012937 correction Methods 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 238000011156 evaluation Methods 0.000 description 2
- 238000011835 investigation Methods 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000004048 modification Effects 0.000 description 2
- 238000012544 monitoring process Methods 0.000 description 2
- 238000007781 pre-processing Methods 0.000 description 2
- 238000001556 precipitation Methods 0.000 description 2
- 238000002310 reflectometry Methods 0.000 description 2
- 230000011218 segmentation Effects 0.000 description 2
- 238000012795 verification Methods 0.000 description 2
- 239000002028 Biomass Substances 0.000 description 1
- 208000003643 Callosities Diseases 0.000 description 1
- 241000208818 Helianthus Species 0.000 description 1
- 235000003222 Helianthus annuus Nutrition 0.000 description 1
- 240000005979 Hordeum vulgare Species 0.000 description 1
- 235000007340 Hordeum vulgare Nutrition 0.000 description 1
- 206010020649 Hyperkeratosis Diseases 0.000 description 1
- 241000209140 Triticum Species 0.000 description 1
- 235000021307 Triticum Nutrition 0.000 description 1
- 238000012271 agricultural production Methods 0.000 description 1
- 239000003153 chemical reaction reagent Substances 0.000 description 1
- 238000007796 conventional method Methods 0.000 description 1
- 238000013527 convolutional neural network Methods 0.000 description 1
- 230000001186 cumulative effect Effects 0.000 description 1
- 238000005520 cutting process Methods 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000010586 diagram Methods 0.000 description 1
- 238000005315 distribution function Methods 0.000 description 1
- 239000003814 drug Substances 0.000 description 1
- 229940079593 drug Drugs 0.000 description 1
- 230000002349 favourable effect Effects 0.000 description 1
- 238000003384 imaging method Methods 0.000 description 1
- 239000000463 material Substances 0.000 description 1
- 230000008569 process Effects 0.000 description 1
- 230000005855 radiation Effects 0.000 description 1
- 238000005070 sampling Methods 0.000 description 1
- 238000001228 spectrum Methods 0.000 description 1
- 239000007921 spray Substances 0.000 description 1
- 238000012360 testing method Methods 0.000 description 1
- 238000000844 transformation Methods 0.000 description 1
- 230000005068 transpiration Effects 0.000 description 1
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/0002—Inspection of images, e.g. flaw detection
-
- 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/13—Edge detection
-
- 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/20—Special algorithmic details
- G06T2207/20112—Image segmentation details
- G06T2207/20132—Image cropping
-
- 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
Abstract
本发明涉及一种基于地块形态特征的农田作物行向遥感识别方法,属于作物行向识别技术领域。解决了现有的作物行向识别方法通常是针对厘米级无人机影像或照片进行的,难以应用于大范围的作物行向制图,而常用的商业卫星和民用卫星的分辨率很难利用地块内部纹理结构特征来识别作物行向的问题。本发明的基于地块形态特征的农田作物行向遥感识别方法,放弃使用地块内部复杂结构特征,以地块为最小研究对象,利用图像处理技术对地块形状和长宽比等简单特征进行作物行方向的识别。本发明的方法对耕地地块形状大小差异大、随机性高的区域有更强的适用性,在大空间尺度的影像上具有较高的应用潜力。
Description
技术领域
本发明属于作物行向识别技术领域,具体涉及一种基于地块形态特征的农田作物行向遥感识别方法,尤其适用于地块平均面积大于0.5km2以上的耕地作物行向的识别。
背景技术
耕地作物种植一般采用大田种植的方式,种植时采用成行种植以充分增加作物受光照面积,利于作物呼吸蒸腾和气体交换,并且能够有效促进雨水充沛季节的田内排水。卫星遥感作为一门能够快速地获取大范围地面遥感信息的实用技术,在农业领域被广泛的应用于作物生长监测、作物分类、作物估产、农业灾害监测以及土壤水分反演等。遥感卫星主要依靠接收遥感信号生成遥感图像,但是遥感信号的反射会受到地面周期性结构的影响。在雷达图像中,垂直于雷达观测视角的行结构和平行于雷达观测视角的行结构相比,同极化后向散射系数相差可达10dB,使用雷达数据进行土壤湿度反演时,由土壤垄行结构引起的误差会超过土壤湿度反演的容限值,说明不同的作物行向会对模型反演结果的精度产生影响。在光学数据中由作物行向引起的误差比雷达数据中误差更直观、更明显。农田作物行的二向反射分布特性使反射率随太阳入射角及观测角变化而有规律的变化,在影像处理时,需要进行双向反射分布函数BRDF校正以减少由成像几何角度引起的反射率差异,有利于提高多源光学影像之间的辐射一致性。所以作物行向作为卫星遥感反演地表参数(生物量、作物产量、叶面积指数、含水量等)的一种重要影响因素,对其定量提取非常重要。
为了定量的提取作物行向,现有的方法一般采用:(1)人工获取;(2)基于厘米级的无人机影像识别技术;(3)基于米级的卫星遥感影像识别技术。人工获取包括实地测量和对高分辨率影像进行目视解译,这两种获取作物行向的方法费时费力、效率低下,特别是对实现大范围、逐地块、快速及时的作物行向识别的可行性不高。无人机在快速获取区域范围内农作物的表型信息方面具有很大优势,在复杂的农田环境条件下,机动灵活的无人机逐渐成为获取农作物行向的重要手段。在农业种植园内,使用基于无人机采集的尺寸为320像素*240像素的生长初期玉米作物图像,利用边缘检测和霍夫变换实现了基于作物行检测的无人机导航作业,在作物生长初期、杂草较少的高分辨率无人机影像中的作物行检测精度可达100%,减少了无人机喷撒药物时偏离规划航线的误差,提高了无人机的作业效率和作业精度。通过对多时相无人机影像的分析,发现识别时期会影响作物行的识别精度,并且确定拔节期是使用高分辨率影像识别玉米行并提取的最佳时间,在土壤和玉米行之间具有高对比度的高分辨率无人机影像中提取玉米行的精度可达100%,为农业管理活动提供了可靠数据。近年来,卷积神经网络被广泛应用于杂草检测、植被分类分割和作物产量估测等方面。SegNet分割网络和霍夫变换相结合得到用于作物行检测的CRowNet网络,对无人机影像中的甜菜行进行检测,在测试集上获得了93.58%检测识别率,此外在甜菜数据集上训练的模型也能够检测出玉米的行结构,证明所提出的算法模型具有较好的泛化性。在空间分辨率为0.075m的航空影像中,作物行结构作为地面上的一种周期性结构可以很容易的被观测到,考虑到频域识别周期性信息的便利性,使用傅里叶幅度谱和方位角扫描可以对作物行向和周期进行提取,并取得相关系数R=0.99,RMSE=2.17°的识别精度,证明在更大范围内使用超高分辨率的航空影像识别作物行向的可行性,提供了一种在大范围尺度下,替代使用无人机影像识别作物行向的有效方案。米级的光学遥感卫星影像在保留一定的地面细节观测能力的基础上,大大提高了地面的观测尺度。使用方向空间滤波器和数学形态学算法对空间分辨率为2米Formosat-2卫星影像中的作物纹理特征进行识别,针对向日葵、玉米等夏季作物和小麦、大麦和油菜等冬季作物分别实现精度为0.27<R2<0.99,7.15°<RMSE<43.02°,0<R2<0.96, 10.22°<RMSE<80°的作物行向检测,为大范围机械化农业生产管理提供了可靠信息。
无人机影像所拥有的丰富细节和纹理信息能有效支撑作物行向的提取,但厘米级的无人机影像空间覆盖范围受续航能力限制,只适合于田块尺度的观测,难以将此类技术实际应用于大区域。虽然米级的遥感卫星影像解决了观测尺度的问题,但是米级的商业遥感卫星影像通常是需要付费获取的,不利于推广。并且在进行大范围的米级遥感影像处理时,需要耗费巨量的计算资源。当使用分辨率较低的Sentinel-2遥感卫星影像以及其他10米以上分辨率的卫星影像时,影像的一个像元范围包括了数个作物行,无法体现作物行和作物行之间,作物行和裸土之间的差异性,因此耕地中的作物行在10米分辨率下的纹理特征较弱,针对地块内部的作物行结构识别具有较高的实现难度。
发明内容
本发明解决了现有的作物行向识别方法通常是针对厘米级无人机影像或照片进行的,难以应用于大范围的作物行向制图,而常用的商业卫星和民用卫星的分辨率在几十厘米到几十米之间,一个像元包含了一个或者多个作物行,导致很难利用地块内部纹理结构特征来识别作物行向的技术问题,提供一种基于地块形态特征的农田作物行向遥感识别方法。
为实现上述目的,采用以下技术方案:
本发明的一种基于地块形态特征的农田作物行向遥感识别方法,步骤如下:
步骤一、获取研究区范围的矢量文件;
步骤二、根据研究区范围的矢量文件,从卫星数据分发网站提取影像,对影像勾画地块矢量边界,并将每个地块的文件格式转换为二值图,且填充二值图的边缘;
步骤三、计算长宽比:
3.1)对步骤二得到的边缘填充后的二值图进行地块边界的轮廓逼近,得到轮廓点形成的曲线围成的图形周长,取周长的1%作为阈值,依次将轮廓点之间的距离与阈值作比较,大于阈值的点舍去,小于阈值的点保留,最后统计保留点的数量,得到顶点数d和所有顶点的坐标(xi,yi);
3.2)使用顶点数d以及顶点坐标(xi,yi)判断地块形状,当d=3时,认为此地块为三角形;当d=4,图形四个内角的值均处于80°~100°之间,认为此地块为矩形;当d=4,两组对角之间的差小于10°且相临的角之间大于10°,认为此地块为平行四边形;当d=4,有一组平行边,认为此地块为梯形;当d=4,有零组平行边,认为此地块为不规则四边形;当5≤d<10时,认为此地块为多边形;当10≤d时,认为此地块为复杂图形;
3.3)若地块为三角形、矩形和平行四边形,直接执行步骤四;
若地块为梯形、不规则四边形和多边形,先对地块进行最小面积矩形拟合,再执行步骤四;
步骤四、使用欧氏距离计算地块三个点的两两之间的距离,得到长边和短边的长度,进而计算地块的长宽比;若地块的长宽比大于长宽比阈值,则执行步骤五;若地块的长宽比不满足长宽比阈值,则放弃识别;
所述长宽比阈值为1~4.5;
步骤五、使用概率霍夫变换或者LSD算法标记地块长边并得到长边端点的坐标(xi1,xy1)和(xi2,xy2);
步骤六、先使用步骤五长边端点的坐标(xi1,xy1)和(xi2,xy2),求得反正切函数值,然后把得到的弧度转换为角度,得到地块长边的角度即为该地块的作物行向角度。
进一步的,所述步骤一中,影像的分辨率为10米或更低。
进一步的,所述步骤二中,使用ArcGIS软件对影像勾画地块矢量边界,并将每个地块的文件格式转换为二值图,使用开源计算机视觉库OpenCV中的填充边界函数填充二值图的边缘。
进一步的,所述步骤二中,二值图的每个边缘各填充200行/列的亮度值为 255的像素点。
进一步的,所述步骤3.3)中,先使用OpenCV中的最小面积矩形拟合函数对地块进行最小面积矩形拟合。
进一步的,所述步骤四中,长宽比阈值为1.4。
与现有技术相比,本发明的有益效果为:
本发明的基于地块形态特征的农田作物行向遥感识别方法,放弃使用地块内部复杂结构特征,以地块为最小研究对象,利用图像处理技术对地块形状和长宽比等简单特征进行作物行方向的识别。解决了无人机影像可表达范围小、难以在大尺度空间范围推广和10米以上分辨率的卫星影像中作物行的纹理特征较弱,难以使用地块内部的作物行结构特征识别作物行向的缺点。实验结果证明使用地块形态特征进行识别作物行方向的方法是有效的,并且证明在进行作物行向识别时,LSD算法(0.81<R2<0.99,21.46°<RMSE<1.78°)在精度评价方面均优于概率霍夫变换算法(0.65<R2<0.92,30.76°<RMSE<13.82°),尤其是对耕地地块形状大小差异大、随机性高的区域有更强的适用性。
本发明提出的作物行向识别方法在大空间尺度的影像上具有较高的应用潜力。
附图说明
为了更清楚地说明本发明实施例中的技术方案,下面将对实施例中所需要使用的附图作简单地介绍,显而易见地,下面描述中的附图仅仅是本发明的一些实施例,对于本领域普通技术人员来讲,在不付出创造性劳动的前提下,还可以根据这些附图获得其它的附图。
图1为本发明实施例1中的研究区概况,其中,(a)为研究区的地理位置; (b)为研究区范围图;(c)为研究区内地块分布的局部放大图;
图2为本发明实施例1中的研究区内的作物行向分布数量图;
图3为本发明实施例1的地块数据预处理结果图,其中,(a)为耕地地块; (b)为地块二值图;(c)为填充后的地块二值图;
图4为本发明实施例1中的作物行向识别方法流程图;
图5为本发明实施例1中的耕地地块形状举例,其中,(a)为三角形;(b) 为梯形;(c)为矩形;(d)为多边形。
图6为本发明实施例1中的不规则地块的最小面积矩形拟合及长边识别结果,其中,(a)不规则地块;(b)最小矩形拟合结果;(c)拟合后边缘检测结果; (d)拟合后直线检测结果;
图7为本发明实施例1的正确率拟合结果;
图8为本发明实施例1的识别率和拟合的程序识别率;
图9为本发明实施例1中的地块二值图中的边缘“阶梯”现象示意图。
具体实施方式
为了进一步理解本发明,下面对本发明优选实施方案进行描述,但是应当理解,这些描述只是为进一步说明本发明的特征和优点,而不是对本发明权利要求的限制。
本发明的基于地块形态特征的农田作物行向遥感识别方法,步骤如下:
步骤一、获取研究区范围的矢量文件;
步骤二、根据研究区范围的矢量文件,从卫星数据分发网站提取影像,对影像勾画地块矢量边界,并将每个地块的文件格式转换为二值图,且填充二值图的边缘;
步骤三、计算长宽比:
3.1)对步骤二得到的边缘填充后的二值图进行地块边界的轮廓逼近,得到轮廓点形成的曲线围成的图形周长,取周长的1%作为阈值,依次将轮廓点之间的距离与阈值作比较,大于阈值的点舍去,小于阈值的点保留,最后统计保留点的数量,得到顶点数d和所有顶点的坐标(xi,yi);
3.2)使用顶点数d以及顶点坐标(xi,yi)判断地块形状,当d=3时,认为此地块为三角形;当d=4,图形四个内角的值均处于80°~100°之间,认为此地块为矩形;当d=4,两组对角之间的差小于10°且相临的角之间大于10°,认为此地块为平行四边形;当d=4,有一组平行边,认为此地块为梯形;当d=4,有零组平行边,认为此地块为不规则四边形;当5≤d<10时,认为此地块为多边形;当10≤d时,认为此地块为复杂图形;
3.3)若地块为三角形、矩形和平行四边形,直接执行步骤四;
若地块为梯形、不规则四边形和多边形,先对地块进行最小面积矩形拟合,再执行步骤四;
步骤四、使用欧氏距离计算地块三个点(如a1、a2、a3)的两两之间的距离,得到长边(a1、a2)和短边(a2、a3)的长度,进而计算地块的长宽比;若地块的长宽比大于长宽比阈值,则执行步骤五;若地块的长宽比不满足长宽比阈值,则放弃识别;
步骤五、使用概率霍夫变换或者LSD算法标记地块长边并得到长边端点的坐标(xi1,xy1)和(xi2,xy2);
步骤六、先使用步骤五长边端点的坐标(xi1,xy1)和(xi2,xy2),求得反正切函数值,然后把得到的弧度转换为角度,得到地块长边的角度即为该地块的作物行向角度。
上述技术方案,步骤一中,影像的分辨率为10米或更低。
上述技术方案,步骤二中,使用ArcGIS软件对影像勾画地块矢量边界,并将每个地块的文件格式转换为二值图,使用开源计算机视觉库OpenCV中的填充边界函数填充二值图的边缘。
上述技术方案,步骤二中,二值图的每个边缘各填充200行/列的亮度值为 255的像素点。
上述技术方案,步骤3.3)中,先使用OpenCV中的最小面积矩形拟合函数对地块进行最小面积矩形拟合。
上述技术方案,步骤四中,长宽比阈值为1~4.5,最优为1.4。
在本发明中所使用的术语,一般具有本领域普通技术人员通常理解的含义,除非另有说明。为了使本领域的技术人员更好地理解本发明的技术方案,下面将结合实施例对本发明作进一步的详细介绍。
在以下实施例中,未详细描述的各种过程和方法是本领域中公知的常规方法。下述实施例中所用的材料、试剂、装置、仪器、设备等,如无特殊说明,均可从商业途径获得。
以下结合实施例进一步说明本发明。
实施例1
1.1研究区概况
研究区位于黑龙江省双鸭山市友谊县(131°27′~132°15′E,46°28~46°59′N),如图1中(a)所示,该地处于黑龙江、松花江、乌苏里江三江平原腹地,属中温带大陆气候,冬冷夏热四季分明,年平均气温2.5℃。年降水量500-600mm,降水集中在6-8月。研究区内地势平坦,土地肥沃,自然条件优越,现有耕地 9.2万公顷,耕地作物中,玉米和大豆种植面积占农作物总种植面积的绝大部分,玉米和大豆在5月初播种,10月收获。种植以农场种植为主,农户种植为辅,机械化程度较高。农场最大地块面积为3.0479km2,最小面积为0.0011km2,平均面积0.5398km2。
1.2数据源
1.2.1Sentinel-2卫星影像和地块边界
本发明使用的数据来源于Google Earth Engine。使用云平台下载的Sentinel-2MSI影像10米分辨率影像,影像包括蓝波段,绿波段,红波段和近红外波段四个波段数据。该影像拍摄于2020年7月23日,由经过辐射定标、大气校正和去云后的四幅影像裁剪拼接得到,如图1中(b)所示。发明用到的地块矢量边界是使用ArcGIS软件对Sentinel-2影像根据地块形状、颜色以及纹理特征进行目视解译勾画得到,如图1中(c)所示。
1.2.2GF-2卫星影像和作物行向的地面真值
数据来源于国产GF-2卫星,国产GF-2卫星是我国自主研发的亚米级高分辨率光学遥感卫星,搭载全色、多光谱两台相机,全色影像星下点分辨率为0.8m,多光谱影像星下点分辨率为3.2m。本发明采用2020年6月1日拍摄的GF-2影像数据,使用ArcGIS内置的COGO工具参照影像中地块的纹理特征,以正东方向为0°,沿逆时针方向旋转,正西方向为180°进行目视解译勾画获得作物行向角度真值,用于验证本发明算法识别作物行向的精度,研究区内500个地块的作物行向在各个角度区间的分布如图2。
1.3地块行向识别算法
如图4所示,步骤如下:
(1)地块数据预处理,包括地块矢量边界二值化和边界填充;
使用ArcGIS Python2.7提供的要素转栅格函数,批量实现将多个耕地地块 (图3(a))矢量文件转为8bit位深的TIFF文件格式的单地块二值图。从图3 (b)可以发现生成的TIFF文件的边缘紧邻地块,会导致最小面积矩形拟合结果越界,因此使用开源计算机视觉库OpenCV中的填充边界函数分别对TIFF图像的上下左右各填充200行/列的亮度值为255的像素点以消除对最小面积矩形拟合的影响,如图3(c)所示。
(2)长宽比计算,包括基于轮廓逼近的地块顶点计算方法、地块形状复杂度判别、最小面积矩形拟合以及地块长宽比计算;
①轮廓逼近,获取地块顶点数:对边缘填充后得到的地块二值图进行地块边界的轮廓逼近,得到轮廓点形成的曲线围成的图形周长,取其1%作为阈值,依次将轮廓点之间的距离与1%的周长阈值作比较,大于这个值的点舍去,小于这个值的点保留,最后统计保留点的数量得到顶点数d和所有顶点的坐标(xi,yi)。
②地块形状复杂度判别:使用轮廓逼近得到的顶点数d以及顶点坐标进行判断地块形状,当d=3时,认为此地块为三角形,如图5(a)。当d=4时所包含的图形种类较多,需要进一步判断:当d=4且图形四个内角的值均处于80°~100°之间,则认为是矩形,如图5(c);当d=4,两组对角之间的差小于10°并且相临的角之间大于10°,认为是平行四边形;当d=4,且有一组平行边,认为是梯形,如图5(b);其他为不规则四边形。当顶点数5≤d<10时,为多边形,如图 5(d)。当10≤d时,判断地块形状为复杂图形。
③长宽比计算:使用欧氏距离计算地块三个点的两两之间的距离,得到长边和短边的长度,进而计算地块的长宽比;若地块的长宽比大于长宽比阈值,则执行步骤五;若地块的长宽比不满足长宽比阈值,则放弃识别。
若地块为三角形、矩形和平行四边形的简单地块,使用顶点坐标中的连续三个顶点的坐标,计算两组相邻的两个顶点的距离。
对于梯形、不规则四边形以及多边形等具有复杂形状,为了避免识别误差和简化计算,算法将复杂形状拟合为矩形后,再进行后续操作。拟合使用OpenCV 中的最小面积矩形拟合函数对不规则地块(如图6中(a))进行最小面积矩形拟合,得到最小矩形拟合结果(如图6中(b)),然后再进行地块长宽比计算(使用拟合后的最小矩形坐标中的连续三个顶点的坐标,计算两组相邻的两个顶点的距离。)
(3)直线检测和角度计算。
概率霍夫变换和LSD算法作为检测直线时的常用算法被引入到本发明中。概率霍夫变换算法相较于经典霍夫变换算法无需对所有像素点逐一进行检测,而是随机选择少量边缘点进行参数空间中的小范围投票。同时可提前停止直线检测的过程,减少了算法运算开销,提高直线检测效率,表1列举了本发明使用概率霍夫变换的具体参数设置。
表1本发明中使用的概率霍夫变换的参数设置
表1中Rho代表以像素为单位的搜索距离精度,Theta代表以弧度为单位的角度精度,Threshold代表识别某部分为直线的累加平面阈值,MinLL代表能组成一条直线的最少的点数量,点数量不足的直线将被抛弃,即线段的最小长度。 MaxLG代表线段上两点之间的最近距离阈值。
LSD算法相比于概率霍夫变换在本发明中最大的优点是通过高斯采样将输入图像缩小为原来的80%,目的在于减弱甚至消除很多图像中出现的阶梯效应。并且LSD算法在检测长度差距大的直线时无需设置过多的参数,这也是在本发明中LSD算法优于概率霍夫变换的重要原因。
对于形状简单的矩形和平行四边形地块,如果地块长宽比满足设置的长宽比阈值,则直接进行概率霍夫变换或者LSD直线检测。对于形状复杂的地块,先进行最小面积矩形拟合(如图6(c)),然后再进行直线检测(如图6(d))。
角度计算首先使用检测得到长边两点坐标(xi1,xy1)和(xi2,xy2),通过两点求得反正切函数值,然后把得到的弧度转换为角度,最终得到地块长边的角度作为该地块的作物行向角度。
1.4长宽比与地块识别精度、角度精度的关系
通过观察遥感卫星影像发现,在我国北方的耕地一般是具有高长宽比的矩形地块,同时地块中的作物行向一般与地块长边的方向相同,但是这两者之间的关系并没有文献对其进行具体的定量描述。因此本发明设置多组不同地块长宽比阈值的实验进行探究作物行向和地块长边的关系。
(1)为了选取使概率霍夫变换检测直线正确率最高的阈值参数MinLL,在长宽比阈值为1.5的实验中,设置概率霍夫变换中的参数MinLL为165、175、 185、195、205、215的六组实验。并针对使用一次概率霍夫变换检测地块长边精度低的问题,设计使用两次概率霍夫变换,提高检测地块长边的精度。本发明设置两次概率霍夫变换的参数MinLL为(120,210)、(130,220)、(140,230)、 (150,240)、(160,250)、(170,260)六组对比实验。
(2)为了对比不同长宽比区间内的作物行向识别正确率和识别精度,设置长宽比的阈值区间为[1.5,2)、[2,2.5)、[2.5,3)、[3,3.5)、[3.5,4)、[4, 4.5)、[4.5,∞)的七组实验。使用在GF-2卫星影像上进行目视解译提取的地面作物行向角度真值作为参考,设置程序结果和地面真值之间的误差小于10°时,认为是程序识别结果是正确的。
(3)为了探究算法在不同长宽比阈值下的地块识别率和作物行向识别精度, 本发明设置长宽比阈值为1.5、2.0、2.5、3.0、3.5、4.0、4.5的对比实验,分析概率霍夫变换算法和LSD算法在不同阈值下检测地块长边的精度。
1.5精度评价指标
本发明使用根据纹理特征目视解译GF-2影像获得作物行向的地面角度真值作为参考,验证本发明算法识别作物行向的精度。为了定量化描述验证结果及算法精度,本发明使用识别率、正确率、决定系数R2以及均方根误差RMSE进行精度验证。
(1)识别率:表示满足长宽比阈值的地块占所有地块的比例。计算公式:
(2)正确率:算法给出正确角度(程序估计值-角度实测值<±10°)的地块数所占满足长宽比阈值的总地块数的比例。计算公式:
(3)决定系数(R2):
(4)均方根误差RMSE:
其中,yi为行向实测值,y′i为行向估计值。
1.6数据结果
图7为正确率拟合结果;
图8识别率和拟合的程序识别率;
表2各长宽比阈值的程序正确率
表3各长宽比阈值区间内的地块数量分布及算法正确率
1.7分析
当地块长宽比较小时,地块的长宽差距不明显,种植作物的规律性不强,所以在低阈值区间范围内的地块长边不一定是作物行的方向,造成结果正确率较低。当地块长宽比阈值大于2.5时,程序识别作物行向的准确率达到97.76%,证明地块中实际的作物种植基本以长边为方向进行种植。但是从表3可以发现,实验地块的长宽比分布是极其不均匀的,地块长宽比t<1.5的地块仅占总地块数的12.2%,其他各阈值区间分别占比为10.49%、10.28%、9.42%、8.35%、6.22%、 5.78%、37.26%,其中长宽比t>2.5的地块占67.89%,同时正确率为93.18%以上。这就造成在进行后面不同阈值的实验时,比如t>1.5、t>2.0的时候包含了t>2.5、 t>3.0、t>3.5等的结果,因此造成了在低阈值正确率过高,因此无法利用研究区内真实的地块长宽比分布数据确定平衡识别率和正确率的最佳阈值。为了从理论上确定一个最佳地块长宽比阈值来平衡地块识别率和识别作物行向的正确率。首先设置地块长宽比阈值以1开始,以0.1为步长进行长宽比阈值细化,属于每个细化后阈值的地块数为10。然后使用LSD算法识别作物行向的正确率数据确定正确率拟合公式,然后利用拟合公式对各个细化后的长宽比阈值进行拟合,得到各个长宽比阈值的识别作物行向的正确率。最后综合地块87.2%的识别率和 79.66%的正确率,确定识别作物行向的最佳长宽比阈值为1.4时,使作物行向的识别正确率最高达到69.47%。
综上,本发明以黑龙江省友谊县为研究区,使用空间分辨率为10m的 Sentinel-2卫星多光谱影像为数据源,将地块作为最小发明对象,尝试分析基于地块形态特征的作物行向识别方法的可行性。根据影像提取地块的矢量边界,融合轮廓检索算法、最小矩形拟合算法、边缘检测算法和直线检测算法等多种图像处理算法计算地块长边与短边的长度比值(长宽比),分析了地块长宽比与行向识别率、识别精度之间关系。结果表明:(1)地块长宽比与作物行方向识别率呈负相关关系,而与方向识别精度呈正相关关系,随着地块长宽比阈值增加,行向的识别率从82.0%降低到34.8%,行向识别均方根误差(RMSE)从21.46°降低到1.78°;(2)在不同长宽比阈值条件下,直线检测器LSD(Line Segment Detector)算法识别作物行向的平均精度(决定系数R2=0.93,RMSE=9.53°)高于概率霍夫变换(R2=0.81,RMSE=20.80°);(3)本发明构建的基于10米分辨率遥感数据的农田作物行向识别技术体系,可以对农田地块作物行向进行大范围识别。本发明方法可为使用遥感卫星影像识别作物行方向提供新的思路参考。
显然,上述实施方式仅仅是为清楚地说明所作的举例,而并非对实施例的限定。对于所属领域的普通技术人员来说,在上述说明的基础上还可以做出其它不同形式的变化或变动。这里无需也无法对所有实施例予以穷举。而由此所引申出的显而易见的变化或变动仍处于本发明创造的保护范围之中。
Claims (7)
1.一种基于地块形态特征的农田作物行向遥感识别方法,其特征在于,步骤如下:
步骤一、获取研究区范围的矢量文件;
步骤二、根据研究区范围的矢量文件,从卫星数据分发网站提取影像,对影像勾画地块矢量边界,并将每个地块的文件格式转换为二值图,且填充二值图的边缘;
步骤三、计算长宽比:
3.1)对步骤二得到的边缘填充后的二值图进行地块边界的轮廓逼近,得到轮廓点形成的曲线围成的图形周长,取周长的1%作为阈值,依次将轮廓点之间的距离与阈值作比较,大于阈值的点舍去,小于阈值的点保留,最后统计保留点的数量,得到顶点数d和所有顶点的坐标(xi,yi);
3.2)使用顶点数d以及顶点坐标(xi,yi)判断地块形状,当d=3时,认为此地块为三角形;当d=4,图形四个内角的值均处于80°~100°之间,认为此地块为矩形;当d=4,两组对角之间的差小于10°且相临的角之间大于10°,认为此地块为平行四边形;当d=4,有一组平行边,认为此地块为梯形;当d=4,有零组平行边,认为此地块为不规则四边形;当5≤d<10时,认为此地块为多边形;当10≤d时,认为此地块为复杂图形;
3.3)若地块为三角形、矩形和平行四边形,直接执行步骤四;
若地块为梯形、不规则四边形和多边形,先对地块进行最小面积矩形拟合,再执行步骤四;
步骤四、使用欧氏距离计算地块三个点的两两之间的距离,得到长边和短边的长度,进而计算地块的长宽比;若地块的长宽比大于长宽比阈值,则执行步骤五;若地块的长宽比不满足长宽比阈值,则放弃识别;
所述长宽比阈值为1~4.5;
步骤五、使用概率霍夫变换或者LSD算法标记地块长边并得到长边端点的坐标(xi1,xy1)和(xi2,xy2);
步骤六、先使用步骤五长边端点的坐标(xi1,xy1)和(xi2,xy2),求得反正切函数值,然后把得到的弧度转换为角度,得到地块长边的角度即为该地块的作物行向角度。
2.根据权利要求1所述的一种基于地块形态特征的农田作物行向遥感识别方法,其特征在于,所述步骤一中,影像的分辨率为10米或更低。
3.根据权利要求1所述的一种基于地块形态特征的农田作物行向遥感识别方法,其特征在于,所述步骤二中,使用ArcGIS软件对影像勾画地块矢量边界,并将每个地块的文件格式转换为二值图,使用开源计算机视觉库OpenCV中的填充边界函数填充二值图的边缘。
4.根据权利要求1所述的一种基于地块形态特征的农田作物行向遥感识别方法,其特征在于,所述步骤二中,二值图的每个边缘各填充200行/列的亮度值为255的像素点。
5.根据权利要求1所述的一种基于地块形态特征的农田作物行向遥感识别方法,其特征在于,所述步骤3.3)中,先使用OpenCV中的最小面积矩形拟合函数对地块进行最小面积矩形拟合。
6.根据权利要求1所述的一种基于地块形态特征的农田作物行向遥感识别方法,其特征在于,所述步骤四中,长宽比阈值为1.4。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202210329863.8A CN114782324B (zh) | 2022-03-31 | 一种基于地块形态特征的农田作物行向遥感识别方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202210329863.8A CN114782324B (zh) | 2022-03-31 | 一种基于地块形态特征的农田作物行向遥感识别方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN114782324A true CN114782324A (zh) | 2022-07-22 |
CN114782324B CN114782324B (zh) | 2024-05-17 |
Family
ID=
Cited By (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN115223059A (zh) * | 2022-08-31 | 2022-10-21 | 自然资源部第三航测遥感院 | 基于多元遥感影像的多云雾区域农作物种植模式提取方法 |
Cited By (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN115223059A (zh) * | 2022-08-31 | 2022-10-21 | 自然资源部第三航测遥感院 | 基于多元遥感影像的多云雾区域农作物种植模式提取方法 |
CN115223059B (zh) * | 2022-08-31 | 2022-12-09 | 自然资源部第三航测遥感院 | 基于多元遥感影像的多云雾区域农作物种植模式提取方法 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN109613513B (zh) | 一种顾及InSAR形变因子的光学遥感潜在滑坡自动识别方法 | |
Che et al. | Estimation of maize plant height and leaf area index dynamics using an unmanned aerial vehicle with oblique and nadir photography | |
CN106023133B (zh) | 一种基于多特征联合处理的高分辨率遥感影像水体提取方法 | |
Jin et al. | High-throughput measurements of stem characteristics to estimate ear density and above-ground biomass | |
Shi et al. | Mapping plastic greenhouses using spectral metrics derived from GaoFen-2 satellite data | |
CN112285710B (zh) | 一种多源遥感水库蓄水量估算方法与装置 | |
CN108982369B (zh) | 融合gf-1wfv和modis数据的地块尺度作物长势监测方法 | |
CN112131946B (zh) | 光学遥感影像植被和水体信息自动提取方法 | |
CN112183209A (zh) | 一种基于多维特征融合的区域农作物分类方法及系统 | |
Cao et al. | Use of unmanned aerial vehicle imagery and a hybrid algorithm combining a watershed algorithm and adaptive threshold segmentation to extract wheat lodging | |
CN114359728B (zh) | 一种河道类型遥感提取方法 | |
CN113505635A (zh) | 基于光学和雷达的冬小麦与大蒜混种区识别方法及装置 | |
CN111091079A (zh) | 基于tls的高寒脆弱区植被优势单株结构参数测定方法 | |
CN112418133A (zh) | 一种基于多源遥感影像的秸秆焚烧监测方法 | |
CN111353402B (zh) | 一种油棕林遥感提取方法 | |
Bai et al. | A fast and robust method for plant count in sunflower and maize at different seedling stages using high-resolution UAV RGB imagery | |
CN113063375B (zh) | 一种直线型耕作田埂的无人机遥感提取方法 | |
Che et al. | The effect of growth stage and plant counting accuracy of maize inbred lines on LAI and biomass prediction | |
Yang et al. | Fraction vegetation cover extraction of winter wheat based on RGB image obtained by UAV | |
CN113534083B (zh) | 基于sar的玉米留茬方式识别方法、装置和介质 | |
Xu et al. | Vegetation information extraction in karst area based on UAV remote sensing in visible light band | |
CN112036313A (zh) | 烟草种植面积检测方法、装置、设备及可读存储介质 | |
CN114782324B (zh) | 一种基于地块形态特征的农田作物行向遥感识别方法 | |
CN114782324A (zh) | 一种基于地块形态特征的农田作物行向遥感识别方法 | |
CN116310864A (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 |