CN106247926B - 基于单元格扫描和gps轨迹插值的农机耕作面积测算方法 - Google Patents

Abstract

基于单元格扫描和GPS轨迹插值的农机耕作面积测算方法，涉及基于GPS轨迹信息的农机耕作面积测算技术。本发明是为了解决现有的农机耕作面积测算方法测算精度低的问题。本发明通过对上传经纬度点进行插值处理标记已耕地的面积，利用插值使面积测算的精度达到一个很高的水平。然后根据坐标值判断出XY方向上的最值，确定出该段耕地路径的一个外切矩形；以农机设备的耕作幅宽为边长作为单元格，以此单元格从外切矩形的一端进行累加扫描，如果单元格内存在插值点则进行累加记录，最后通过累加单元格的数量计算得到最终耕地面积。本发明适用于农机耕作面积测算场合。

Description

基于单元格扫描和GPS轨迹插值的农机耕作面积测算方法

技术领域

本发明涉及基于GPS轨迹信息的农机耕作面积测算技术。

背景技术

土地面积测算是地籍测量中一项十分重要的工作内容。它为调整土地利用结构，合理分配土地，制定国民经济计划、农业区划以及土地利用规划等提供数据基础。土地面积测算包括行政管辖区、宗地、农业耕作等面积的测算。目前土地面积测算的方法有两种，即解析法面积测算和图解法面积测算。解析法由于运算复杂而导致计算速度较慢，而图解法由于设备和人为的因素往往精度不满足要求。

发明内容

本发明是为了解决现有的农机耕作面积测算方法测算精度低的问题，从而提供一种基于单元格扫描和GPS轨迹插值的农机耕作面积测算方法。

基于单元格扫描和GPS轨迹插值的农机耕作面积测算方法，其特征是：

步骤一、利用安装在农机设备上的GPS定位装置，每隔时间t采集一次农机设备的经纬度坐标，并将采集的经纬度坐标按照高斯投影变换法转换为平面直角坐标系下的XY坐标；

步骤二、将步骤一的时间序列作为原始时间序列，若相邻的坐标间的距离大于预设值M米，则进行分段处理，M为正数；

经分段处理后得到k段坐标点集C_i，i＝1,...,k，再将原始时间序列的坐标点按照等分的方式插入至k段坐标点，使相邻插值点间的距离小于或等于给定值L；

步骤三、判断每一段坐标点的方位点，找到XY方向最靠近外侧的四个坐标，即：X_min、X_max、Y_min、Y_max，为了使边界的坐标点也包含在单元格中，则拓展这四个坐标值为X_min-L、X_max+L、Y_min-L、Y_max+L；将这四个坐标通过组合形成四个边界点，以这四个边界点为基准作XY轴的平行线，构建外切矩形R_i；

步骤四、以农机设备的耕作幅宽为边长，在每段点集的外切矩形R_i中构建正方形单元格，用此正方形单元格铺满所述外切矩形R_i，遍历每个单元格，判断此单元格中是否存在插值后的坐标点；

若存在，则进行计数累加，最后得到每个点集C_i含有坐标点单元格的总个数A_i(i＝1,...,k)，对每一分段的A_i进行累加得到时段内整体作业单元格数再利用S＝N*L²即得到该农机时段内的耕地作业面积S。

步骤二中的给定值L取值为：农机设备耕作幅宽长度。

步骤四中，遍历每个单元格，判断此单元格中是否存在插值后的坐标点；若存在，则进行计数累加，最后得到每个点集C_i含有坐标点单元格的总个数A_i，i＝1,...,k，具体为：

从矩形边界的左上方开始扫描，从起始坐标开始以幅宽为长度进行累加，以确定所在单元格的边界；对于每个单元格遍历矩形框内的坐标点，若判断该单元格中存在坐标点则进行计数累加；当判断到达矩形右边界时坐标进行换行，即将单元格顶点横坐标设置为起始值而纵坐标在原来的基础上减去幅宽，然后继续重复坐标累加、判断；当坐标值达到右下边界时，跳出循环。

本发明有效利用现代信息技术中的地理信息系统，通过单元格扫描的方法对农机设备的耕地面积进行测算，由于采用了微积分的计算方式，所以是一种高效、便利、精确度高的解决方案。本发明通过对上传经纬度点的分段处理，将大片随机不连续的耕地区域划分成较为简单的图形面积的计算，提高了运算效率。

附图说明

图1是本发明的基于单元格扫描和GPS轨迹插值农机耕作面积测算方法流程图；

图2是分带投影示意图；

图3是高斯-克吕格投影示意图；

图4是确定单元格扫描示意图。

具体实施方式

具体实施方式一、基于单元格扫描和GPS轨迹插值的农机耕作面积测算方法，

步骤一、利用安装在农机设备上的GPS定位装置，每隔时间t采集一次农机设备的经纬度坐标，并将采集的经纬度坐标按照高斯投影变换法转换为平面直角坐标系下的XY坐标；

步骤二、将步骤一的时间序列作为原始时间序列，若相邻的坐标间的距离大于预设值M米，则进行分段处理，M为正数；

经分段处理后得到k段坐标点集C_i，i＝1,...,k，再将原始时间序列的坐标点按照等分的方式插入至k段坐标点，使相邻插值点间的距离小于或等于给定值L；

步骤三、判断每一段坐标点的方位点，找到XY方向最靠近外侧的四个坐标，即：X_min、X_max、Y_min、Y_max，为了使边界的坐标点也包含在单元格中，则拓展这四个坐标值为X_min-L、X_max+L、Y_min-L、Y_max+L；将这四个坐标通过组合形成四个边界点，以这四个边界点为基准作XY轴的平行线，构建外切矩形R_i；

步骤四、以农机设备的耕作幅宽为边长，在每段点集的外切矩形R_i中构建正方形单元格，用此正方形单元格铺满所述外切矩形R_i，遍历每个单元格，判断此单元格中是否存在插值后的坐标点；

若存在，则进行计数累加，最后得到每个点集C_i含有坐标点单元格的总个数A_i(i＝1,...,k)，对每一分段的A_i进行累加得到时段内整体作业单元格数再利用S＝N*L²即得到该农机时段内的耕地作业面积S。

步骤二中的给定值L取值为：农机设备耕作幅宽长度。

步骤四中，遍历每个单元格，判断此单元格中是否存在插值后的坐标点；若存在，则进行计数累加，最后得到每个点集C_i含有坐标点单元格的总个数A_i(i＝1,...,k)，具体为：

从矩形边界的左上方开始扫描，从起始坐标开始以幅宽为长度进行累加，以确定所在单元格的边界；对于每个单元格遍历矩形框内的坐标点，若判断该单元格中存在坐标点则进行计数累加；当判断到达矩形右边界时坐标进行换行，即将单元格顶点横坐标设置为起始值而纵坐标在原来的基础上减去幅宽，然后继续重复坐标累加、判断；当坐标值达到右下边界时，跳出循环。

其中，在步骤一中，通过安装在农机设备上的GPS定位装置每隔固定时间上传农机的经纬度坐标，并将上传的经纬度按照高斯投影变换转换为平面直角坐标系下的XY坐标，高斯投影全称为高斯-克吕格投影变换，是一种等角横切椭圆柱投影变换，经过变换后的坐标才能进行距离上的运算处理，进而计算出面积。

具体实施例：

设GPS设备某时刻上传的经纬度分别为longitude和latitude(单位为度)。

计算1°所对应的弧度为iPI＝π/180.0。

如图2所示，采用分带投影的方法，使投影边缘的变形不致过大：ZoneWide＝6即6°带。

坐标系参数选用54年北京坐标系参数：

椭球体长半轴长度：a＝6378137.0

椭球体短半轴长度：b＝6356752.3142

将经纬度分别转换为弧度longitude1和latitude1：

longitude1＝longitude*iPI

latitude1＝latitude*iPI

查找带号：

ProjNo＝longitude/ZoneWide

longitude0＝ProjNo*ZoneWide+ZoneWide/2

得到中央经线：longitude2＝longitude0*iPI

椭球第一偏心率：

椭球第二偏心率：

计算卯酉圈曲率半径：

T＝tan²(latitude1)

C＝e₁*cos²(latitude1)

A＝(longitude1-longitude0)*cos(latitude1)

计算子午线弧长：

S＝a*((1-e₂/4-3*e₂ ²/64-5*e₂ ³/256)*latitude1-(3*e₂/8+3*e₂ ²/32+45*e₂ ³/1024)*sin(2*latitude1)+(15*e₂ ²/256+45*e₂ ³/1024)*sin(4*latitude1)-(35*e₂ ³/3072)*sin(6*latitude1))

带内大地坐标：

xval＝R*(A+(1-T+C)*A³/6+(5-18*T+T²+72*C-58*e)*A⁵/120)

yval＝S+R*tan(latitude1)*(A²/2+(5-T+9*C+4*C²)*A⁴/24+(61-58*T+T²+600*C-330*e)*A⁶/720)

计算横坐标方向上的偏移量：X₀＝1000000*(ProjNo+1)+500000。

计算得到直角坐标：X＝xval+X₀，Y＝yval(单位为米)。

进一步的，在步骤二中，插值点的间隔可以设置为农机设备幅宽L的k倍(k<1),这样插值点的密集程度可以根据计算要求的精度和单元格选取的边长来决定；理论上在一定范围内，插值越密集，计算的精度越高。

进一步的，在步骤四中，对单元格进行遍历，判断是否有插值后的坐标点落在此单元格中，若存在则记录下来，最后得到含有坐标点单元格的个数n，具体过程为：

如图4所示，在不规则图形abmn中存在插值后的坐标点，以这些坐标点确定所在扫描矩形，从矩形边界的左上方开始扫描，从起始坐标开始以幅宽为长度进行累加，以确定所在单元格的边界；对于每个单元格遍历矩形框内的坐标点，若判断该单元格中存在坐标点则进行计数累加；当判断到达矩形右边界时坐标进行换行，即把单元格顶点横坐标设置为起始值而纵坐标在原来的基础上减去幅宽L，然后继续重复坐标累加、判断；当坐标值达到右下边界时，跳出循环。对于第i段计算得到的单元格数A_i进行累加得到上传点所在的单元格总数N,即面积S＝N*L²(单位：平方米)。

本发明对比现有技术，有如下有益效果：传统的解析算法由于需要经过较为复杂的解析运算，所以运算速度较慢，对于大量数据的处理效率较低；而图解算法由于对于设备和人为操作等诸多因素的要求，使得精确度往往较低。针对计算精度问题，本发明利用微积分的思想对上传的坐标进行插值处理，这样可以以任意小的单元格进行积分累加，可以将精确度控制在理想水平。同时由于分段操作避免了大量无上传坐标点区域的判断，同时由于没有经过复杂的运算极大地提高了计算效率。

以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。

Claims (3)

1.基于单元格扫描和GPS轨迹插值的农机耕作面积测算方法，其特征是：

步骤一、利用安装在农机设备上的GPS定位装置，每隔时间t采集一次农机设备的经纬度坐标，并将采集的经纬度坐标按照高斯投影变换法转换为平面直角坐标系下的XY坐标；

步骤二、将步骤一的时间序列作为原始时间序列，若相邻的坐标间的距离大于预设值M米，则进行分段处理，M为正数；

经分段处理后得到k段坐标点集C_i，i＝1,...,k，再将原始时间序列的坐标点按照等分的方式插入至k段坐标点，使相邻插值点间的距离小于或等于给定值L；

步骤三、判断每一段坐标点的方位点，找到XY方向最靠近外侧的四个坐标，即：X_min、X_max、Y_min、Y_max，为了使边界的坐标点也包含在单元格中，则拓展这四个坐标值为X_min-L、X_max+L、Y_min-L、Y_max+L；将这四个坐标通过组合形成四个边界点，以这四个边界点为基准作XY轴的平行线，构建外切矩形R_i；

步骤四、以农机设备的耕作幅宽为边长，在每段点集的外切矩形R_i中构建正方形单元格，用此正方形单元格铺满所述外切矩形R_i，遍历每个单元格，判断此单元格中是否存在插值后的坐标点；

若存在，则进行计数累加，最后得到每个点集C_i含有坐标点单元格的总个数A_i，i＝1,...,k，对每一分段的A_i进行累加得到时段内整体作业单元格数再利用S＝N*L²即得到该农机时段内的耕地作业面积S。

2.根据权利要求1所述的基于单元格扫描和GPS轨迹插值的农机耕作面积测算方法，其特征在于步骤二中的给定值L取值为：农机设备耕作幅宽长度。

3.根据权利要求2所述的基于单元格扫描和GPS轨迹插值的农机耕作面积测算方法，其特征在于步骤四中，遍历每个单元格，判断此单元格中是否存在插值后的坐标点；若存在，则进行计数累加，最后得到每个点集C_i含有坐标点单元格的总个数A_i，i＝1,...,k，具体为：

从矩形边界的左上方开始扫描，从起始坐标开始以幅宽为长度进行累加，以确定所在单元格的边界；对于每个单元格遍历矩形框内的坐标点，若判断该单元格中存在坐标点则进行计数累加；当判断到达矩形右边界时坐标进行换行，即将单元格顶点横坐标设置为起始值而纵坐标在原来的基础上减去幅宽，然后继续重复坐标累加、判断；当坐标值达到右下边界时，跳出循环。