CN103714563B - 一种曲线形农田作业区域边界建模方法 - Google Patents

Abstract

本发明公开了一种曲线形农田作业区域边界建模方法，采集曲线形农田作业区域边界的关键点，找出所有关键点中横坐标值最大的点和横坐标值最小的点，利用三次样条函数建立上边界模型和下边界模型。本发明可以将曲线形农田作业区域的边界以及农田作业区域内的障碍物顶点精准地记录并绘制出来，从而为农业植保机的作业人员提供准确的农田作业区域，为农业植保机的路径和航迹规划以及计算农田作业面积提供可靠的依据。

Description

一种曲线形农田作业区域边界建模方法

技术领域

本发明涉及一种曲线形农田作业区域边界建模方法。

背景技术

随着科技的不断发展，采用现代化机械来替换人工劳动已成为各行各业的流行趋势。我国作为一个传统的农业大国，耕地面积十分广阔，然而，目前国内农田作业方面还是采用传统的人工作业方式。其中农药喷洒采取的是人工喷洒的方式，这种方式不仅效率低，而且对作业人员的身心有极大的伤害，因此，急需先进的技术来改变这一现象。农业植保机应运而生，然而，现今所有的植保机都是靠人工进行操作，对农田的喷洒区域及边界完全靠操作员的眼睛进行判断，难免产生误判、漏喷及多喷的现象。因此，如果操作员能够获得准确的农田作业区域的边界，在无人直升机进行作业时，实时的将路径和航迹在显示器上显示出来，操作员便可根据航迹判断出飞机的位置，从而很大程度的解决误判、漏喷及多喷的问题，且大大提高了农业喷洒的安全性。而这一问题的关键便在于要准确的获取农田的边界信息。因此，急需一种好的农田作业区域边界的建模方法来解决这一实际问题。

北京农业信息技术研究中心的发明专利“采集农田关键顶点测绘成图的方法”提供了一种采集农田关键顶点测绘成图的方法，通过采集四种交界顶点、两种轮廓顶点一级一种辅助顶点后，根据关键顶点间的拓扑关系实时自动校验并查看测绘结果。包括以下步骤：S₁：获取GPS位置信息；S₂：勾勒待测区域轮廓；S₃：测绘农田关键分界顶点，并标注地块名称和注释；S₄：实时校验、提示与分割农田；S₅：上传测绘数据；S₆：获取矢量地图。

其中，在S₂中，勾勒待测区域轮廓包括：顺序标定和无序标定，所述顺序标定，以轮廓顶点的编号为顺序，依次标定待测区域轮廓；所述无序标定，是指自动标定包含全部轮廓顶点的面积最大的多边形为区域轮廓。

现有方法是以绘出农田的矢量地图为目的，而并非以画出农田作业区域边界为目的，因此，所绘出的农田轮廓包含两个部分：作业区域和非作业区域。在将此农田轮廓图用于农业植保机喷洒农药的路径和航迹规划等场合时，会出现非作业区域的多喷现象，浪费农药；在用于计算农田作业面积时，会出现作业面积的计算不准确现象；现有方法未对农田区域边界建立数学模型，从而无法为农田作业区域的路径和航迹规划提供可靠依据；没有将农田作业区域以内可能存在的障碍物顶点（如电线杆、树、信号发射塔等）标记出来，在将此农田轮廓图用于农业植保机喷洒农药的路径和航迹规划等场合时，可能导致农业植保机出现撞机或坠毁危险现象发生。

发明内容

本发明所要解决的技术问题是，针对现有技术不足，提供一种曲线形农田作业区域边界建模方法，将曲线形农田作业区域的边界以及农田作业区域内的障碍物顶点精准地记录并绘制出来，从而为农业植保机的作业人员提供准确的农田作业区域，为农业植保机的路径和航迹规划以及计算农田作业面积提供可靠的依据。

为解决上述技术问题，本发明所采用的技术方案是：一种曲线形农田作业区域边界建模方法，该方法为：

1）从曲线形农田作业区域边界的某一个关键点D₁开始，以顺时针方向采集曲线形农田作业区域边界的所有关键点，设共有m个关键点，分别为D₁(x₁,y₁),D₂(x₂,y₂),D₃(x₃,y₃),…,D_m(x_m,y_m)；建立直角坐标系，将所述m个关键点（每隔一个采样周期（GPS本身的采样频率的倒数）就自动采集一个点，当回到初始位置，即关键点D₁所在位置时结束采点）设置在直角坐标系的第一象限内；

2）比较上述m个关键点的横坐标，找出横坐标值最大的关键点和最小的关键点；若有多个关键点的横坐标值均为最小值，且该多个关键点是第一个关键点D₁与最后一个关键点D_m之间的关键点，则将次序最先的关键点作为横坐标值最小的关键点，若第一个关键点D₁和最后一个关键点D_m的横坐标值都是最小值，则将最后一个关键点作为横坐标值最小的关键点；若有多个关键点的横坐标值均为最大值，且该多个关键点是第一个关键点D₁与最后一个关键点D_m之间的关键点，则将次序最后的关键点作为横坐标值最大的关键点，若第一个关键点D₁和最后一个关键点D_m的横坐标值都是最大值，则将第一个关键点D₁作为横坐标值最大的关键点；将横坐标值最小的关键点D_i和最大的关键点D_j分别作为所有的关键点在直角坐标系中的左端和右端；

3）以D_i为起始点，D_j为终点；以顺时针方向，令D_i=B₀,D_i+1=B₁,…,D_j=B_n，将从D_i到D_j的各点依次即为B₀(x₀,y₀),K,B_n(x_n,y_n)；建立以下上边界函数模型S_t(x)：

$S_{t\left(x\right)}=M_t\frac{{(x_{t+1}-x)}^3}{{6h}_t}+M_{t+1}\frac{{(x-x_t)}^3}{{6h}_t}+(y_t-\frac{M_t{h}_t^2}{6})\frac{x_{t+1}-x}{h_t}+(y_{t+1}-\frac{M_{t+1}{h}_t^2}{6})\frac{{x-x}_t}{h_t};$

其中，x∈[x_t,x_t+1,],t=0,1,...,n-1；x_t、x_t+1分别为B₀、B_n之间的点B_t、B_t+1的横坐标，且B_t+1为B_t顺时针方向上的下一点；系数M_t、M_t+1通过追赶法求解；h_t=x_t+1-x_t；y_t、y_t+1分别为点B_t、B_t+1的纵坐标；

4）以D_i为起始点，D_j为终点；以逆时针方向，令D_i=C₀,D_i-1=C₂,…,D_j=C_r，建立以下上边界函数模型S_r(x)：

$S_{r\left(x\right)}=M_r\frac{{(x_{r+1}-x)}^3}{{6h}_r}+M_{r+1}\frac{{(x-x_r)}^3}{{6h}_r}+(y_r-\frac{M_r{h}_r^2}{6})\frac{x_{r+1}-x}{h_r}+(y_{r+1}-\frac{M_{r+1}{h}_r^2}{6})\frac{{x-x}_r}{h_r};$

其中，x∈[x_r,x_r+1,],r=0,1,...,m-n；x_r、x_r+1分别为C₀、C_m-n之间的点C_r、C_r+1的横坐标，且C_r+1为C_r逆时时针方向上的下一点；系数M_r、M_r+1通过追赶法求解；h_r=x_r+1-x_r；y_r、y_r+1分别为点C_r、C_r+1的纵坐标。

当所述曲线形农田作业区域内存在障碍物时，采集曲线形农田作业区域内障碍物顶点的坐标Z(x₀,y₀)，建立以下障碍物警戒线函数方程：

一级障碍物警戒线的函数方程为(x-x₀)²+(y-y₀)²=a²；

二级障碍物警戒线的函数方程为(x-x₀)²+(y-y₀)²=b²；

三级障碍物警戒线的函数方程为(x-x₀)²+(y-y₀)²=c²；

其中，0<a<b<c<50。

通过追赶法求解系数M_t、M_t+1的步骤为：

1）令β₁=b₁,y₁=d₁；其中， $b_1=2+{\mathrm{\&mu;}}_1+\frac{{\mathrm{\&mu;}}_1{h}_0}{h_1},b_2=2,...,b_{n-2}=2,$ $b_{n-1}=2+{\mathrm{\&lambda;}}_{n-1}+\frac{{\mathrm{\&lambda;}}_{n-1}{h}_{n-1}}{h_{n-2}};$ ${\mathrm{\&mu;}}_p=\frac{h_{p-1}}{h_{p-1}+h_p};$ $d_p=6\frac{f[x_p,x_{p+1}]-f[x_{p-1},x_p]}{h_{p-1}+h_p}=6\frac{\frac{y_{p+1}-y_p}{x_{p+1}-x_p}-\frac{y_p-y_{p-1}}{x_p-x_{p-1}}}{h_{p-1}+h_p};$ λ_p=1-μ_p；p=1,2,...,n-1；

2）计算β_q=b_q-l_qc_q-1,y_q=d_q-l_qy_q-1；其中，q=2,3,…,n-1； $a_2={\mathrm{\&mu;}}_2,$ $a_3={\mathrm{\&mu;}}_3,...,a_{n-1}={\mathrm{\&mu;}}_{n-1}-\frac{{\mathrm{\&lambda;}}_{n-1}{h}_{n-1}}{h_{n-2}};$ $c_1={\mathrm{\&lambda;}}_1-\frac{{\mathrm{\&mu;}}_1{h}_0}{h_1},c_2={\mathrm{\&lambda;}}_2,c_3={\mathrm{\&lambda;}}_3,...,c_{n-2}={\mathrm{\&lambda;}}_{n-2};$

3）利用下式求解M_n-1：

4）利用下式计算M_n-2,M_n-3,L,M₁：s=n-2,n-3,…,1；

5）通过下式计算M₀和M_n： $M_0=M_1-\frac{h_0}{h_1}(M_2-M_1),$ $M_n=M_{n-1}+\frac{h_{n-1}}{h_{n-2}}(M_{n-1}-M_{n-2}).$

与现有技术相比，本发明所具有的有益效果为：本发明可以将曲线形农田作业区域的边界以及农田作业区域内的障碍物顶点精准地记录并绘制出来，从而为农业植保机的作业人员提供准确的农田作业区域，为农业植保机的路径和航迹规划以及计算农田作业面积提供可靠的依据。

附图说明

图1为本发明一实施例曲线型边界模型示意图；

图2为本发明一实施例有障碍物的农田作业区域障碍物警戒线模型。

具体实施方式

以下结合附图详细说明本发明的具体实施方式。

1）假设采点工作者以第一个关键点开始，以顺时针方向依次采集m个关键点（以GPS的周期T_s作为采集农田作业区域边界点的采样周期），分别为D₁(x₁,y₁),D₂(x₂,y₂),D₃(x₃,_y3),…,D_m(x_m,y_m)，如图1所示；

2）比较这m个点的横坐标，找出其中横坐标最小和最大的点。假设横坐标最小的点为D_i(x_i,y_i)，横坐标最大的点为D_j(x_j,y_j)。(若遇到多个点的横坐标都为最小值，则把取点次序最先的点作为最小点，在此有个特例：如果D₁和D_m的横坐标都相等，但是在此取D_m作为最小点；若遇到多个点的横坐标都为最大值，则把取点次序最后的点作为最大点，在此有个特例：如果D_m和D₁的横坐标都相等，但是在此取D₁作为最大点)，那么D_i(x_i,y_i)即为左端，D_j(x_j,y_j)即为右端；

3）假设由左端D_i与右端D_j所确定的上边界有|j-i|+1个点，令n=|j-i|（包括左端D_i与右端D_j）,为便于接下来的计算，将上边界上各点（包括左端D_i与右端D_j）从左至右按顺时方向依次记为B₀(x₀,y₀),K,B_n(x_n,y_n)（即令D_i=B₀,D_i+1=B₁,…,D_j=B_n若当j<i，则在顺时针方向移点的过程中，如果D_n出现在D_j的左边（即x_n<x_j）,令D_n的下一个点为D₁,如果D_n=B_k(k<n)那么D₁=B_k+1,D₂=B_k+2,K,D_j=B_n。

则上边界函数S(x)满足：

（1）S(x_t)=y_t(j=0,1,…,n)；

（2）S(x)在每个小区间[x_t,x_t+1](t=0,1,…,n-1)上是3次多项式；

（3）S(x)在[x₀,x_n]上有连续2阶导数；

（4）端点条件：强迫第一个和第二的三次多项式的三阶导数相等，对最后一个和倒数第二个三次多项式也做同样处理。

按求解3次样条函数的三弯矩方程方法，设S″(x_t)=M_t，t=0,1,...,n，

记h_t=x_t+1-x_t,t=0,1,...,n-1。

（1）整理三次样条函数的推导过程得μ_tM_t-1+2M_t+λ_tM_t+1=d_t，t=1,2,...,n-1,其中， ${\mathrm{\&mu;}}_t=\frac{h_{t-1}}{h_{t-1}+h_t},$ λ_t=1-μ_t, $d_t=6f[x_{t-1},x_t,x_{t+1}]=6\frac{f[x_t,x_{t+1}]-f[x_{t-1},x_t]}{h_{t-1}+h_t}=6\frac{\frac{y_{t+1}-y_t}{x_{t+1}-x_t}-\frac{y_t-y_{t-1}}{x_t-x_{t-1}}}{h_{t-1}+h_t}.$

（2）由端点条件知， $\begin{array}{c}\frac{M_n-M_{n-1}}{h_{n-1}}=\frac{M_{n-1}-M_{n-2}}{h_{n-2}}\&DoubleRightArrow;{M_n=M}_{n-1}+\frac{h_{n-1}}{h_{n-2}}(M_{n-1}-M_{n-2})\\ \frac{M_{\text{1}}-M_{\text{0}}}{h_0}=\frac{M_2-M_1}{h_1}\&DoubleRightArrow;M_0=M_1-\frac{h_0}{h_1}(M_2-M_1)\end{array},$

由μ_tM_t-1+2M_t+λ_tM_t+1=d_t,t=1,2,...,n-1知μ₁M₀+2M₁+λ₁M₂=d₁

∴ $(2+{\mathrm{\&mu;}}_1+\frac{{\mathrm{\&mu;}}_1{h}_0}{h_1})M_1+({\mathrm{\&lambda;}}_1-\frac{{\mathrm{\&mu;}}_1{h}_0}{h_1})M_2=d_1,$

由μ_tM_t-1+2M_t+λ_tM_t+1=d_t,t=1,2,...,n-1知μ_n-1M_n-2+2M_n-1+λ_n-1M_n=d_n-1,

∴ $({\mathrm{\&mu;}}_{n-1}-\frac{{\mathrm{\&lambda;}}_{n-1}{h}_{n-1}}{h_{n-2}})M_{n-2}+(2+{\mathrm{\&lambda;}}_{n-1}+\frac{{\mathrm{\&lambda;}}_{n-1}{h}_{n-1}}{h_{n-2}})M_{n-1}=d_{n-1}.$

（3）结合（1）和（2）得矩阵

求出三弯矩阵的系数阵和右端项。

（4）另 $b_1=2+{\mathrm{\&mu;}}_1+\frac{{\mathrm{\&mu;}}_1{h}_0}{h_1},b_2=2,...,b_{n-2}=2,b_{n-1}=2+{\mathrm{\&lambda;}}_{n-1}+\frac{{\mathrm{\&lambda;}}_{n-1}{h}_{n-1}}{h_{n-2}},$

$c_1={\mathrm{\&lambda;}}_1-\frac{{\mathrm{\&mu;}}_1{h}_0}{h_1},c_2={\mathrm{\&lambda;}}_2,c_3={\mathrm{\&lambda;}}_3,...,c_{n-2}={\mathrm{\&lambda;}}_{n-2},$

$a_2={\mathrm{\&mu;}}_2,a_3={\mathrm{\&mu;}}_3,...,a_{n-1}={\mathrm{\&mu;}}_{n-1}-\frac{{\mathrm{\&lambda;}}_{n-1}{h}_{n-1}}{h_{n-2}},$

则（3）中的矩阵变为 $\left(\begin{array}{ccccccc}{b}_1& c_1& & & & & \\ a_2& b_2& c_2& & & & \\ & a_3& b_3& c_3& & & \\ & & \&CenterDot;& \&CenterDot;& \&CenterDot;& & \\ & & & \&CenterDot;& \&CenterDot;& \&CenterDot;& \\ & & & & a_{n-2}& b_{n-2}& c_{n-2}\\ & & & & & a_{n-1}& b_{n-1}\end{array}\right)\left(\begin{array}{c}{M}_1\\ M_2\\ M_3\\ \&CenterDot;\\ \&CenterDot;\\ M_{n-2}\\ M_{n-1}\end{array}\right)=\left(\begin{array}{c}{d}_1\\ d_2\\ d_3\\ \&CenterDot;\\ \&CenterDot;\\ d_{n-2}\\ d_{n-1}\end{array}\right)$

（5）用追赶法求出M_t,t=1,2,…,n-1。

具体过程如下：

①β₁=b₁,y₁=d₁；

②对于i=2,3,…,n-1，计算β_i=b_i-l_ic_i-1,y_i=d_i-l_iy_i-1；

③ $M_{n-1}=\frac{y_{n-1}}{{\mathrm{\&beta;}}_{n-1}};$

④对于i=n-2,n-3,…,1，计算

⑤计算 $M_0=M_1-\frac{h_0}{h_1}(M_2-M_1),$ $M_n=M_{n-1}+\frac{h_{n-1}}{h_{n-2}}(M_{n-1}-M_{n-2});$

（6）按下列公式求出S_(x)在各个区间[x_t,x_t+1]上表达式S_t(x)，其为 $S_{t\left(x\right)}=M_t\frac{{(x_{t+1}-x)}^3}{{6h}_t}+M_{t+1}\frac{{(x-x_t)}^3}{{6h}_t}+(y_t-\frac{M_t{h}_t^2}{6})\frac{x_{t+1}-x}{h_t}+(y_{t+1}-\frac{M_{t+1}{h}_t^2}{6})\frac{{x-x}_t}{h_t}$ $(x\&Element;[x_t,x_{t+1,}],t=\mathrm{0,1},...,n-1)$

4）下边界函数模型

左端D_i与右端D_j所确定的下边界上有n+1=m-|j-i|+1个点（包括左端D_i与右端D_j）,为便于后来计算，将下边界上各点（包括左端D_i与右端D_j）从左至右按逆时针方向依次记为C₀(x₁,y₁),K,C_n(x_n,y_n)（即令D_i=C₀,D_i-1=C₂,…,D_j=C_n，若当j>i，则在逆时针方向移点的过程中，如果D₁出现在D_j的左边（即X₁<X_j），那么令D₁的下一个点为D_n；如果D₁=C_k(k<n),那么D_n=C_k+1，D_n-1=C_k+2,...,D_j=C_n）。

下边界函数计算与上边界函数计算过程一样。

如果农田中有障碍物，则由采点工作者采集到障碍物点的坐标为Z(x₀,y₀)，假设一级警戒线为离障碍物点a米的圆，二级警戒线为离障碍物点b米的圆，三级警戒线为离障碍物点c米的圆（50>c>b>a>0），如图2所示，图2中三个圆从内至外分别为一级警戒线、二级警戒线、三级警戒线。

则一级警戒线的函数方程为(x-x₀)²+(y-y₀)²=a²；

二级警戒线的函数方程为(x-x₀)²+(y-y₀)²=b²；

三级警戒线的函数方程为(x-x₀)²+(y-y₀)²=c²。

Claims (3)

1.一种曲线形农田作业区域边界建模方法，其特征在于，该方法为：

1)从曲线形农田作业区域边界的某一个关键点D₁开始，以顺时针方向采集曲线形农田作业区域边界的所有关键点，设共有m个关键点，分别为D₁(x₁,y₁),D₂(x₂,y₂),D₃(x₃,y₃),…,D_m(x_m,y_m)；建立直角坐标系，将所述m个关键点设置在直角坐标系的第一象限内；

2)比较上述m个关键点的横坐标，找出横坐标值最大的关键点和最小的关键点；若有多个关键点的横坐标值均为最小值，且该多个关键点是第一个关键点D₁与最后一个关键点D_m之间的关键点，则将次序最先的关键点作为横坐标值最小的关键点，若第一个关键点D₁和最后一个关键点D_m的横坐标值都是最小值，则将最后一个关键点作为横坐标值最小的关键点；若有多个关键点的横坐标值均为最大值，且该多个关键点是第一个关键点D₁与最后一个关键点D_m之间的关键点，则将次序最后的关键点作为横坐标值最大的关键点，若第一个关键点D₁和最后一个关键点D_m的横坐标值都是最大值，则将第一个关键点D₁作为横坐标值最大的关键点；将横坐标值最小的关键点D_i和最大的关键点D_j分别作为所有的关键点在直角坐标系中的左端和右端；若有多个关键点的横坐标值均为最小值，其中一个关键点为第一个关键点或最后一个关键点，其余关键点为第一个关键点与最后一个关键点之间的关键点，则取从第一个关键点到最后一个关键点的闭区间内次序最先的点作为最小点；若有多个关键点的横坐标值均为最大值，其中一个关键点为第一个关键点或最后一个关键点，其余关键点为第一个关键点与最后一个关键点之间的关键点，则取从第一个关键点到最后一个关键点的闭区间内次序最后的点作为最大点；

3)以D_i为起始点，D_j为终点；以顺时针方向，令D_i＝B₀,D_i+1＝B₁,…,D_j＝B_n，将从D_i到D_j的各点依次即为B₀(x₀,y₀),...,B_n(x_n,y_n)；建立以下上边界函数模型S_t(x)：

$S_{t\left(x\right)}=M_t\frac{{(x_{t+1}-x)}^3}{6h_t}+M_{t+1}\frac{{(x-x_t)}^3}{6h_t}+(y_t-\frac{M_t{h_t}^2}{6})\frac{x_{t+1}-x}{h_t}+(y_{t+1}-\frac{M_{t+1}{h_t}^2}{6})\frac{x-x_t}{h_t};$

其中，x∈[x_t,x_t+1,],t＝0,1,...,n-1；x_t、x_t+1分别为B₀、B_n之间的点B_t、B_t+1的横坐标，且B_t+1为B_t顺时针方向上的下一点；系数M_t、M_t+1通过追赶法求解；

h_t＝x_t+1-x_t,t＝0,1,...,n-1；y_t、y_t+1分别为点B_t、B_t+1的纵坐标；

4)以D_i为起始点，D_j为终点；以逆时针方向，令D_i＝C₀,D_i-1＝C₁,…,D_j＝C_r，建立以下上边界函数模型S_r(x)：

$S_{r\left(x\right)}=M_r\frac{{(x_{r+1}-x)}^3}{6h_r}+M_{r+1}\frac{{(x-x_r)}^3}{6h_r}+(y_r-\frac{M_r{h_r}^2}{6})\frac{x_{r+1}-x}{h_r}+(y_{r+1}-\frac{M_{r+1}{h_r}^2}{6})\frac{x-x_r}{h_r};$

其中，x∈[x_r,x_r+1,],r＝0,1,...,m-n；x_r、x_r+1分别为C₀、C_m-n之间的点C_r、C_r+1的横坐标，且C_r+1为C_r逆时时针方向上的下一点；系数M_r、M_r+1通过追赶法求解；h_r＝x_r+1-x_r；y_r、y_r+1分别为点C_r、C_r+1的纵坐标。

2.根据权利要求1所述的曲线形农田作业区域边界建模方法，其特征在于，当所述曲线形农田作业区域内存在障碍物时，采集曲线形农田作业区域内障碍物顶点的坐标Z(x₀,y₀)，建立以下障碍物警戒线函数方程：

一级障碍物警戒线的函数方程为(x-x₀)²+(y-y₀)²＝a²；

二级障碍物警戒线的函数方程为(x-x₀)²+(y-y₀)²＝b²；

三级障碍物警戒线的函数方程为(x-x₀)²+(y-y₀)²＝c²；

其中，0＜a＜b＜c＜50。

3.根据权利要求1或2所述的曲线形农田作业区域边界建模方法，其特征在于，通过追赶法求解系数M_t、M_t+1的步骤为：

1)令β₁＝b₁,y₁＝d₁；其中，b₂＝2，...，b_n-2＝2， λ_p＝1-μ_p；p＝1，2，...，n-1；

2)计算β_q＝b_q-l_qc_q-1,y_q＝d_q-l_qy_q-1；其中，q＝2,3,…,n-1；a₂＝μ₂，c₂＝λ₂，c₃＝λ₃，...，c_n-2＝λ_n-2；

3)利用下式求解M_n-1：

4)利用下式计算M_n-2,M_n-3,…,M₁：s＝n-2，n-3，...，1；

5)通过下式计算M₀和M_n： $M_n=M_{n-1}+\frac{h_{n-1}}{h_{n-2}}(M_{n-1}-M_{n-2}).$