Welcome to Smart Agriculture 中文
Special Issue--Key Technologies and Equipment for Smart Orchard

Three-Dimensional Virtual Orchard Construction Method Based on Laser Point Cloud

  • FENG Han , 1 ,
  • ZHANG Hao 1 ,
  • WANG Zi 1 ,
  • JIANG Shijie 1 ,
  • LIU Weihong 1 ,
  • ZHOU Linghui 2 ,
  • WANG Yaxiong 2 ,
  • KANG Feng 2 ,
  • LIU Xingxing 1 ,
  • ZHENG Yongjun , 1, 3
Expand
  • 1. College of Engineering, China Agricultural University, Beijing 100083, China
  • 2. College of Engineering, Beijing Forestry University, Beijing 100083, China
  • 3. Modern Agricultural Equipment and Facilities Engineering Research Center, Beijing 100083, China

Received date: 2022-07-05

  Online published: 2022-11-04

Highlights

To solve the problems of low level of digitalization of orchard management and relatively single construction method, a three-dimensional virtual orchard construction method based on laser point cloud was proposed in this research. First, the hand-held 3D point cloud acquistion equipment (3D-BOX) combined with the lidar odometry and mapping (SLAM-LOAM) algorithm was used to complete the acquisition of the point cloud data set of orchard; then the outliers and noise points of the point cloud data were removed by using the statistical filtering algorithm, which was based on the K-neighbor distance statistical method. To achieve this, a distance threshold model for removing noise points was established. When a discrete point exceeded, it would be marked as an outlier, and the point was separated from the point cloud dataset to achieve the effect of discrete point filtering. The VoxelGrid filter was used for down sampling, the cloth simulation filtering (CSF) cloth simulation algorithm was used to calculate the distance between the cloth grid points and the corresponding laser point cloud, and the distinction between ground points and non-ground points was achieved by dividing the distance threshold, and when combined with the density-based spatial clustering of applications with noise (DBSCAN) clustering algorithm, ground removal and cluster segmentation of orchard were realized; finally, the Unity3D engine was used to build a virtual orchard roaming scene, and convert the real-time GPS data of the operating equipment from the WGS-84 coordinate system to the Gauss projection plane coordinate system through Gaussian projection forward calculation. The real-time trajectory of the equipment was displayed through the LineRenderer, which realized the visual display of the motion trajectory control and operation trajectory of the working machine. In order to verify the effectiveness of the virtual orchard construction method, the test of orchard construction method was carried out in the Begonia fruit and the mango orchard. The results showed that the proposed point cloud data processing method could achieve the accuracy of cluster segmentation of Begonia fruit trees and mango trees 95.3% and 98.2%, respectively. Compared with the row spacing and plant spacing of fruit trees in the actual mango orchard, the average inter-row error of the virtual mango orchard was about 3.5%, and the average inter-plant error was about 6.6%. And compared the virtual orchard constructed by Unity3D with the actual orchard, the proposed method can effectively reproduce the actual three-dimensional situation of the orchard, and obtain a better visualization effect, which provides a technical solution for the digital modeling and management of the orchard.

Cite this article

FENG Han , ZHANG Hao , WANG Zi , JIANG Shijie , LIU Weihong , ZHOU Linghui , WANG Yaxiong , KANG Feng , LIU Xingxing , ZHENG Yongjun . Three-Dimensional Virtual Orchard Construction Method Based on Laser Point Cloud[J]. Smart Agriculture, 2022 , 4(3) : 12 -23 . DOI: 10.12133/j.smartag.SA202207002

1 引 言

在果园面积、水果产量逐年增长的趋势下,果园的规范化、标准化管理变得尤为重要,利用三维重建技术构建出一个包含大量农业生产信息的三维虚拟果园,不仅能更加立体、直观地显示真实果园和果树三维形态数据,而且可以结合Unity3D引擎实现虚拟果园漫游场景以及果园作业机械工作轨迹的可视化,极大地方便工作人员对果园的管理,有效地促进智慧果园、无人果园的发展1
由于果园地形的复杂性和果树生长的随机性,实现真实环境下果园的三维模型重建成为难点。目前主要采用视觉和激光雷达(Light Detection and Ranging,LiDAR)进行三维重建技术研究。但是视觉传感器多适用于单个物体或者小型模型的构建,对于建筑物集群、森林、果园等大规模的三维模型重建存在易受外部环境影响、构建误差大等缺点2-6。激光雷达利用其自身穿透性强、直线传播距离远和扫描角度大等优点,被广泛应用于文物保护、建筑物外廓提取、森林建模等领域7-9。对于果园这种大规模的三维模型构建,Wang等10采用移动式激光雷达获取果园点云图,结合轻量级激光雷达惯性测量单元(LiDAR-IMU)状态估计器和旋转约束优化算法来重建果园。Hu等11提出了一种基于多站激光点云数据的树干参数估计和树干模型重建方法,通过建立的模型实现树干几何参数的准确提取,可用于森林等大尺度模型的可视化和生物量估算。为了提高三维模型的构建精度,Indirbai等12基于激光雷达扫描获取的点云库,结合分层最小分割算法和超体素聚类算法,实现了树冠和树干的检测以及单株植株的分割,并利用激光点云库重建出高精度的树木模型。翟晓晓等13基于激光雷达点云数据,通过图优化的即时定位与地图构建(Simultaneous Localization and Mapping,SLAM)方法实现了树木的三维重建。Zheng和Fu14利用激光雷达批量优化框架,研究出了实现批量优化的理论导数,并将批量优化框架融合到激光测距与测绘(LiDAR Odometry and Mapping,LOAM)框架中,降低了激光建图过程中的累计误差。Zhang和Singh15针对不同时间接收到的点云存在运动畸变导致点云配准结果出错等问题,提出了一种基于6自由度的2轴单线激光雷达实时构建激光里程计与建图的方法,但是不具备后端优化的功能。
上述研究表明,激光雷达可以有效获取果园等大型规模特征的三维点云数据,但是对点云数据进行处理的准确性仍有待提高,获取大型果园中果树的数量和位置特征描述也不详细。本研究以海棠果园和芒果园作为研究对象,提出了一种基于激光点云数据的三维虚拟果园构建方法。采用统计滤波算法、布料滤波算法(Cloth Simulation Filtering,CSF)和聚类算法(Density-Based Spatial Clustering of Applications with Noise,DBSCAN)对手持式三维点云采集设备(3D-BOX)获取果园点云数据集进行离散点去除、地面去除以及聚类分析,基于降采样处理和Unity3D引擎完成虚拟果园的构建,并将构建的果园与作业机械二次开发,实现作业机械运动轨迹控制与作业轨迹可视化。

2 设备与方法

2.1 三维数据采集

本研究选取两处地点作为虚拟果园构建对象,一处是中国农业大学工学院西侧海棠果园,试验时间是2021年3月1日,另一处是在广西省百色市田阳县芒果庄园,试验时间是2021年4月8日。使用激光雷达对这两处果园进行扫描完成三维点云的原始数据采集。
使用大连航佳机器人科技有限公司的手持式三维点云采集设备(3D-BOX)进行三维点云数据采集,如图1所示。该设备分为数据采集、客户端显示两部分。数据采集设备是由Velodyne VLP-16激光雷达、惯性测量单元(Inertial Measurement Unit,IMU)、微型计算机组成;微型计算机连接激光雷达和IMU,驱动传感器运行,采集三维点云数据并执行建图,同时通过无线网络与客户端进行数据交互。
图1 手持式三维点云采集设备

Fig. 1 Hand-held 3D point cloud acquistion equipment

选用北京雷格瑞科技有限公司提供的全球导航卫星系统(Global Navigation Satellite System,GNSS)定位装置LGR-BD982接收器模块作为GPS(Global Positioning System)数据采集设备。

2.2 点云预处理

利用3D-BOX中激光雷达完成对果园的三维立体扫描后,结合LOAM算法实现对果园的三维果园点云的预处理。该预处理过程主要包括特征点提取、查找对应边缘线与局部面片、激光雷达运动估计和建图四个方面,构建出三维果园原始点云模型。原始三维点云数据模型中,由于激光雷达的精度、数据采集时的速度和稳定性以及被测物体易受环境因素的影响,导致点云数据有大量噪声点,同时非被测物体的遮挡还会产生一些离群点,因此需要对此类数据点进行去除。此外,庞大的点云数据对后续的可视化处理会产生一定的影响,所以在保证果树点云结构表达清晰的情况下,对整体点云数据集进行降采样,减少点云数量。在点云数据后处理环节,主要采用去除离群点、去除地面、聚类以及点云降采样算法。本研究的模型构建方法流程图如图2所示。
图2 三维果园构建流程图

Fig. 2 Flowchart of three-dimensional orchard construction

2.2.1 特征点提取

特征点提取过程如下:通过曲率值c的大小对特征点进行分类,当c为连续数据点集中的极大值时,将该点设为边缘点;当c为连续数据点集中的极小值时,将该点设为平面点。边缘特征点的选取从最大曲率c开始进行,平面特征点则与之相反。假设激光雷达扫描点云数据集{L}中第k帧点云数据为 p k p k中的某一数据点为i i p k,则在相同扫描帧中数据点i邻域点的集合为S,并定义点i处局部曲面的c函数如公式(1)所示。
c = 1 S X k , i L j S , j i X k , i L - X k , j L

2.2.2 查找对应边缘线与局部面片

边缘线的选取方法是:假设边缘点 i ε ˙ k,在不同扫描线上查找i点的最近邻点,然后根据曲率值对这两个点进行判断,若是边缘点,则依据这两点确定一条直线,获得边缘线。对应局部面片的选取方法是:假设平面点 i H ˙ k,那么i对应的局部面片可通过不共线的三个点确定,在同一扫描线上选取i的两个最近邻点,然后在不同扫描线上选取一个对应i点的最近邻点,从而确定出平面点i对应的局部面片,该处最近邻点的平面点选取也要对其曲率值进行计算,判断是否为平面点。
上述 ε ˙ k H ˙ k定义如下:假设第k-1帧扫描起始时间为 t k - 1,扫描结束时间为 t k,扫描结束获得的当前帧点云数据为 p k - 1,利用空间变换矩阵将当前帧点云数据映射到 t k时刻点云数据坐标系下,转换后的点云数据为 p ˙ k - 1,利用 t k时刻获取的第k帧点云数据 p k与转换后的点云数据 p ˙ k - 1进行配准,完成粗略位姿估计,即完成 t k - 1t时刻的点云配准,从而完成激光雷达的位姿估计。在第k帧扫描获取的点云数据 p k中,提边缘点、平面点分别组成边缘点集 ε k、平面点集 H k,并在 p ˙ k - 1中查找边缘点集 ε k对应 ε ˙ k的边缘线以及平面点集 H k对应 H ˙ k的局部面片。
所以本研究需要判断相邻两帧点云中的最邻近点,分别计算边缘点集 ε k点到 ε ˙ k对应边缘线和平面点集 H k点到 H ˙ k对应局部面片的最小距离计算特征点与查找到的对应边缘线、局部面片的最小距离。计算方程如公式(2)公式(3)所示。
d ε = ε k , i L - p ˙ k - 1 , j L × ε k , i L - p ˙ k - 1 , l L p ˙ k - 1 , j L - p ˙ k - 1 , l L
d H = H ˙ k , i L - p ˙ k - 1 , j L T T                                       T = p ˙ k - 1 , j L - p ˙ k - 1 , l L × p ˙ k - 1 , j L - p ˙ k - 1 , m L
其中, t k时刻边缘点 i ε ˙ k t k - 1时刻最近邻点 j , l p ˙ k - 1 L j , l确定边缘线; d ε为边缘点i到对应边缘线的最小距离; t k时刻平面点 i H ˙ k t k - 1时刻最近邻点 j , l , m p ˙ k - 1 L d H为平面点i到其对应局部面片的距离。

2.2.3 激光雷达运动估计

假设激光雷达匀速运动, t k , t时间段激光雷达的位姿变换矩阵为 T k L = t x , t y , t z , θ x , θ y , θ zi为第k次帧点云数据 p k中的任意一点,则 t k , t k , i时段的激光雷达位姿 T k , i L可通过插值得到;点云数据 p k中的边缘点和平面点进行坐标转换得到,如公式(4)所示。
T ( ( k , i ) ) L = ( t ( ( k , i ) ) - t k ) / ( t - t k ) T k L X ˙ k , i L = R k , i L X k , i L + τ k , i L                                        R k L ( t ) = e θ ¯ k ( t ) )                                                               e θ ¯ k ( t ) ) = I + n ̂ s i n θ + n ̂ 2 ( 1 - c o s θ )                 
其中, X K , J L为边缘点或平面点; X ¯ K , J L为坐标转换后相应的边缘点或平面的坐标; R k , i L为旋转矩阵; τ k , i L为平移矩阵; θ = T k L 4 6 n ̂ n的斜对称矩阵; n为旋转方向的单位向量,计算如公式(5)所示。
n   =   T k , i L 4 6 T k , i L 4 6
边缘点与边缘线、平面点与局部面片的误差函数如公式(6)所示。
f ε ( X k , i L , T ( k + 1 ) L ) = d ε , i ε k    f H ( X k , i L , T ( k + 1 ) L ) = d H , i H k f ε ( X k , i L ) = d                               
误差函数 f与特征点是相对应的,d表示该点到对应点的距离,可计算 f的雅可比矩阵 J,更新激光雷达的位姿 T k , i L。通过非线性的方式使d趋近于0求解 f ( T k + 1 L )λ可由最小二乘法确定,函数如公式(7)所示。
T k L T k L - J T J + λ d i a g J T J - 1 J T d
其中, f为误差函数; J 为雅可比矩阵。

2.2.4 激光雷达构图

激光雷达LOAM算法三维建图过程是在世界坐标系W下进行,且建图的频率比运动估计的频率低一个数量级。假设在激光雷达坐标系L下, t k , t k + 1时间段转换后的点云数据为 p ˙ k,激光雷达位姿估计为 T k L t k + 1 t k时刻扫描k帧点云的地图为 Q k - 1,且激光雷达在该时刻的位姿为 T k W t k,建图算法可根据 T k L t k + 1 T k W t k转换到世界坐标系下得到 T k W t k + 1,同时将点云数据 p ˙ k转换到世界坐标系下为 Q ˙ k,通过 T k W t k + 1 Q ˙ k Q k - 1对齐,并对激光雷达位姿 T k W t k + 1不断进行优化。

2.3 点云数据后处理

经过对激光点云的预处理过后,可以基本展现出果园树的模型。在激光雷达扫描后的点云地图中,包含了果树、地面点云和大量的噪声点。经过点云预处理后的点云数据,虽然果树的轮廓可以凸显出来,但是大量的地面点云不仅会增加计算量,而且非常不利于后续对点云地图的聚类分析处理。因此,本研究将对果园点云模型进行后处理。

2.3.1 去除离群点

由于3D-BOX获取的点云为散乱点云,采用基于点云库(The Point Cloud Library,PCL)的统计滤波器16, 17去除离群点。该方法对点云集合 P中的任意点 p i,查找离该点最近的k个点 p j , p j + 1 , p j + k - 1,分别计算出点 p i与点 p j , p j + 1 , p j + k - 1的欧式距离并求平均值,如公式(8)所示。
a i = 1 k j = 1 j + k - 1 P i - P j 2
得到所有K邻近点的平均值 a 1 a N统计结果近似服从高斯分布,所有点与k个邻近点距离的均值 μ和标准差 σ决定了其形状,记为 A N μ , σ 2,如公式(9)公式(10)
μ = 1 N i = 1 N a i
σ = 1 N i = 1 N a i - μ 2
使用K邻近点距离统计方法建立去除噪声点的距离阈值模型,如公式(11)
t h r e s h o l d = μ + a × σ
其中,a表示为标准差的倍数。
遍历点云集合中所有点 p i的K邻近点距离平均值,将其分别与距离阈值进行比较,若不在 - t h r e s h o l d , + t h r e s h o l d范围内,即点 p ik个邻近点的距离超出平均距离,且大于a倍的标准差时,则标记为离群点,将该点从点云数据集中分离出。该处理流程如图3
图3 去除离群点流程图

Fig. 3 Flowchart of removing outlier

2.3.2 地面点云去除

由激光雷达扫描获取的原始数据集,点云数量庞大且杂乱无章,采用CSF算法18, 19进行地面点云去除操作,保留非地面点云。其工作原理是将点云倒置,通过布料模拟滤波算法20,确定出布料的最终形状,使其作为区分地面点和非地面点的基础。布料模拟滤波算法是针对激光点云基于布料模拟的一种改进。改进后的布料滤波算法网格点位置和外部作用力之间的关系表示为公式(12)
X ( t + Δ t ) = 2 X ( t ) - X ( t - Δ t ) + G m Δ t 2
其中, m是网格点的质量,通常设为1;Δt为时间步长;G为常数。
网格点的位移向量 d 计算方法如公式(13)所示。
d = 1 2 b p i - p 0 n
其中, p 0表示待移动网格点当前位置; p i表示为 p 0邻近网格点的位置; n 表示为垂直方向的单位向量 n = 0,0 , 1 Tb被用作判别网格点是否移动,若网格点为可移动点,则b为1,否则为0。该过程流程如图4
图4 布料滤波算法流程图

Fig. 4 Flowchart of cloth simulation filtering algorithm(CSF)

2.3.3 DBSCAN聚类分析

DBSCAN21是基于密度聚类的经典算法之一,可将点密度足够大的区域分割成簇,并能在含有噪声点的数据中找到任意形状的簇。本研究中以Eps聚类半径和MinPts核心点邻域内最少样本数两个参数作为区域划分的重要前提,将点密度高的区域划分成簇,形成满足密度相连点的最大集合。采用DBSCAN算法对果园树木进行聚类分析,通过聚类得出果园树木的总数。

2.3.4 点云降采样

将处理完成的点云地图导入Unity3D引擎进行二次开发,点云数量过多将会增大项目开发的内存占用,因此还需对点云数据进行降采样,在保留果园点云特征信息最大化的条件下降低点云密度,减少文件占用内存。本研究基于体素滤波方法对PCL点云库聚类后的点云数据进行降采样处理。

2.4 虚拟果园设计

果园漫游场景将通过虚拟的方式显示实际果园情况,以果园点云地图为基础,在Unity3D引擎22中进行二次开发,对果园点云模型进行处理,开发三维虚拟果园。主要包含场景漫游、轨迹显示、坐标交互以及文件管理等功能。构建虚拟果园所使用到的软件包括点云数据处理库PCL、Unity3D开发引擎、点云文件格式转换工具MeshLab、作业机械喷雾机模型处理工具3DMax以及脚本语言C#等。

2.4.1 GPS数据处理

(1)经纬度提取。该功能采用C#中的StreamReader类在文件流中实现原始数据逐行读取。本研究在读取一行数据时采用两种数据提取的方式,一种是按特殊字符进行查询,再对字符串进行截取;另一种是根据目标对象在字符串中所处的位置进行查询,对指定位置的字符串进行截取。截取后的经纬度字符串转换为双浮点型。为方便查询,将转换为双浮点型的经纬度数据先存入二维数组,再添加到词典临时存储。
(2)坐标转换。GPS数据采集设备获取的数据属于WGS-84大地坐标,是包含了经纬度和海拔高度的一种球面坐标。本研究将在二维平面进行轨迹绘制,故需对球面坐标进行处理,将大地坐标进行高斯投影坐标换算,转换为高斯投影平面坐标。坐标求解如公式(14)公式(15)所示。
x = X + t 2 N c o s 2 θ B l 2 + Δ x y = N c o s B l + N 6 c o s 3 B ( 1 - t 2 + η 2 ) l 3 + Δ y
t = t a n   B N = a 1 - e 2 s i n 2   B l = L - L 0 η 2 = e 2 1 - e 2 c o s 2   B
其中,所有经度纬度都是换算为弧度的经度纬度,(°); X对应纬度的子午线弧长; B为坐标点的纬度, N为卯酉圈半径; a为长半轴; e为椭球第一偏心率; L为点的经度,(°); L 0为中央子午线的经度,(°)。
本研究在WGS-84坐标系向高斯投影坐标系转换处理中采用3°分带,经坐标转换后获得的y值前两位数字为代号,同时防止y值为负,在y值的结果上加500 km处理,即整体对y值进行了相同距离的平移,因此对轨迹显示并无影响。

2.4.2 漫游场景设计

以重建的三维虚拟果园为基础,将作业机械加载到虚拟果园环境中,进而模拟果园作业机械工作情况。通过作业机械实现真实世界与虚拟世界的映射,建立作业机械工作轨迹绘制与作业过程模拟等功能。在虚拟果园环境中,添加平移、旋转和缩放等必要的查看功能,并在显示界面右上方地点显示栏接入了百度地图普通定位应用程序接口API(Application Programming Interface)和系统时间,实现果园漫游场景的实时呈现。在虚拟果园的其他功能中添加了测量工具、时间地点显示以及历史数据保存等,建立了一个比较完善果园作业机械显示场景。基于Unity3D 开发平台使用System.IO.Ports类设计了GPS数据采集交互界面,具备实时读取GPS串口数据、存储数据、经纬度提取、坐标转换处理。
Unity3D能够进行层级式开发环境、可视化编辑环境等操作,并且可使用Direct3D、OpenGL、APIs等图形引擎和PhysX物理引擎,完成开发调试功能,开发人员在开发过程中可随时查看三维虚拟果园漫游状况。
将点云进行纹理化处理主要在Unity3D中实现,主要过程如下:将目标点云信息与图像纹理信息进行匹配,得到具有真实颜色信息的点云模型;着色后的点云可以直接生成密集的三角网模型,对三角网模型进行渲染时进行RGB颜色匹配,得到具有真实颜色的三角网模型;对于较稀疏的三角网模型,进行纠正与纹理映射,最后得到具有纹理真实感的三维果园模型。

2.4.3 轨迹可视化

GPS数据采集设备采集的GPS数据处理是轨迹绘制的关键环节,处理的结果将直接影响后续的轨迹可视化。采集的GPS数据转换为高斯投影平面坐标数据后,使用该数据基于Unity3D平台实现轨迹的可视化。将果园点云模型导入Unity3D坐标原点位置为质心位置,绘制轨迹过程将以质心位置为起始点,经度数据沿Z轴方向,纬度数据沿X轴方向进行绘制;在实际采集GPS数据时,固定了采集的起始位置,因此在绘制轨迹时要先计算出质心坐标与固定起点两者之间的偏移量,然后读取平面坐标数据进行绘制。轨迹绘制部分依赖于Unity3D中的LineRenderer线渲染器,可对线的材质、颜色、宽度等属性进行设置。本研究采用默认材质和颜色,宽度设置为0.05 m,并设置绘制轨迹坐标系为世界坐标系。

3 结果分析与功能验证

利用3D-BOX扫描获取的海棠果园与芒果园全局点云数据,点云数据采集完成后,安装GPS数据采集设备于作业机械上,在该区域运行采集轨迹数据。图5为芒果园实景图。
图5 芒果园实景图

Fig. 5 View of mango orchard

3.1 算法验证结果分析

3.1.1 点云预处理

将激光雷达三维扫描的海棠果园与芒果园原始点云数据进行特征点提取、查找对应边缘线与局部面片、激光雷达运动估计和建图等初步预处理操作后,得到预处理后果园点云数据集。原始果园点云与预处理后点云如图6所示。
图6 点云预处理前后效果对比图

Fig. 6 Comparison of the effects before and after the point cloud pretreatment

3.1.2 点云后处理

(1)统计滤波去除离群点。在Ubuntu系统下基于PCL点云库进行统计滤波对预处理后的果园点云数据集进行去除离群点。设置K邻近搜索点数阈值r为40,标准差 σ的倍数a为0.6,此处两个数值r和a的设置为经验值。在海棠果园数据处理过程中,去除离群点滤波前该数据集共有4,746,821个点,滤波之后数据集为4,113,008个点。在芒果园数据处理过程中,滤波前该数据集共有4,303,845个点,滤波之后数据集共有3,548,440个点,芒果园果树的轮廓已经凸显,但仍存在少量杂点。由于实际环境杂草较多且高度不一,导致激光雷达扫描的地面点云数据较为稀疏,使用统计滤波去除了部分地面点云,效果如图7所示,同时验证了该阈值参数下的统计滤波可用于其他场景。
图7 统计滤波后点云效果图

Fig. 7 Point cloud effect after statistical filtering

(2)CSF地面去除。采用CSF算法对PCL点云库进行地面点云数据去除操作,海棠果园去除效果如图8(a)所示,绿色为果树点云,红色为地面点云。芒果园试验环境坡度相对较大,滤波后地面点云稀疏,去除效果如图8(b)所示。在对地面点云数据进行去除后,过滤后的树冠点云特征变化不明显。
图8 CSF去除地面点云效果图

Fig. 8 Results of ground point cloud removal using CSF

(3)聚类分析。通过分别统计实际该海棠果园与芒果园树林树木,进行聚类分析如下:设置聚类半径Eps为1,密度阈值MinPts为100;海棠果园树木为43棵,其中两棵非同一品种且树高不同。经聚类同样可得到的43个簇。统计实际环境芒果树共有57棵树;经聚类得到56个簇,使用正检率与准确率来对果园激光点云进行聚类分割进度分析,正检率(A)和准确率(R)可由公式(16)得出:
A = n N × 100 % R = n M × 100 %
其中,n为正确分割的树木数量,棵;N为分割得到的树木数量,棵;M为人工计数的实际树木数量。果园DBSCAN聚类准确率结果如表1所示。
表1 果园DBSCAN聚类准确率

Table1 Orchard DBSCAN clustering accuracy

果园种类 实际树木数量/棵 分割得到树木数量/棵 正确分割树木数量/棵 正检率/% 准确率/%
芒果果园 43 43 41 100 95.3
海棠果园 57 56 56 98.2 98.2
在海棠果园中,由于树林株间、行间树木均存在郁闭状况,因此DBSCAN聚类并没有完全将一整棵树聚类,准确率只有95.3%,但聚类结果与实际树林统计棵树一致。对于芒果园,树林株间、行间树木也均存在郁闭状况,且株间芒果树郁闭较为严重,因此DBSCAN聚类并没有完全将每棵树聚类,正检率和准确率为98.2%,聚类分割效果还有待提高。聚类效果如图9所示。
图9 果园 DBSCAN聚类效果图

Fig. 9 Results of orchard DBSCAN clustering

(4)降采样处理。在保证树木点云特征基本一致的情况下,采取了降采样处理,体素尺寸设置0.05×0.05×0.05,降采样前海棠果园点云数目为2,316,788个,降采样后为610,076个。芒果园降采样前点云数目为2,216,632个,降采样后为855,130个,进而建立出降采样处理效果图。将降采样处理前和处理后的点云数据合并展示,红色点为降采样处理前的数据点,绿色点为降采样处理后的数据点。降采样后效果图如图10所示。
图10 降采样处理前后数据合并显示效果图

Fig. 10 Merge display of before and after data downsampling

(5)精确度检测。因果园地面有一定坡度且地面不平整,地面杂草较多,重建后的果园在去除地面点云后果树高度会存在误差,所以本研究主要通过测量果树的株距和行距来检测重建精确度,如图11所示。由于海棠果园中夹杂着两颗非同一品种的树,所以在测量位置时会存在较大的误差。本研究在芒果园的实际环境中,选取4处位置进行测量,这4处位置分别是扫描第一行果树的起始位置、中间位置、第一行末尾位置以及扫描结束位置,测量数据如图12所示,两者进行比较。结果表明,平均行间误差约为3.5%,平均株间误差约为6.6%,因株间果树郁闭严重,误差较大,最大误差处于10%左右。通过比较实际芒果树与建模后株距和行距精确度,如图13所示,株距精确度大于88.8%,行距的精确度大于93.6%。
图11 芒果树行距、株距实际测量

Fig. 11 Mango tree row spacing and plant spacing actual measurement

图12 芒果树点云行距、株距扫描测量结果

Fig. 12 Scan results of mango tree point cloud spacing and plant spacing

图13 芒果园构建精确度

Fig. 13 Construction accuracy of mango orchard

3.1.3 三维果园构建效率

本研究果园三维模型重建工作环境CPU为i7-9750H,8 G运行内存,Ubuntu18.04系统下进行编译运行。三维果园建园过程中会产生大量的点云,算法的运行效率会随着点云的数量增加而逐渐降低。点云预处理过程中主要采用激光雷达结合LOAM算法,所消耗的时间较短。因此本研究重点比较海棠果园与芒果园在后处理过程中点云数量以及其运行时间,如表2所示。通过累计上述三个时间段的时长,海棠果园的点云后处理与芒果园的点云后处理时间都在30 s内,具有较高的处理效率。
表2 不同算法运行时间表

Table 2 Run schedule of the algorithms

点云数据 点云处理过程 预处理后点云 统计滤波 CSF布料滤波 降采样
海棠果园 点云数量/个 4,746,821 4,113,008 2,316,788 670,076
运行时间/s 约5 22.4 1.2 0.5
芒果园 点云数量/个 4,303,845 3,548,440 2,216,632 855,130
运行时间/s 约5 20.9 1.1 0.5

3.2 虚拟果园系统功能测试验证

为验证虚拟果园系统所具备的功能有效运行,本研究在广西省百色市田阳县芒果庄园和中国农业大学海棠果园进行了测试。主要对果园场景漫游以及果园作业机械运行轨迹可视场景进行验证。通过登录界面进入三维虚拟果园系统主界面,在系统主界面根据用户需求选择功能场景,点击对应按钮可进入相应场景。点击选择轨迹文件按钮后即可绘制轨迹,该处轨迹文件数据是经过转换的高斯投影平面坐标数据。通过对虚拟果园系统所涉及的功能进行测试以及系统所具备的功能进行展示,如图14所示,海棠果园和芒果园作业机械与作业轨迹可视化效果较好。测试结果表明,虚拟果园系统的各功能均能良好运行。
图14 虚拟果园功能测试结果

Fig. 14 Functional test results of virtual orchards

4 结 论

本研究提出了一种基于激光雷达的三维虚拟果园构建方法,并通过海棠果园和芒果园对方法进行验证。主要结论如下:
(1)构建的果园模型具有较高的准确性。本文选取实际环境下的芒果园果树行距、株距进行测量,与重建后的结果进行对比,平均行间误差约为3.5%,平均株间误差约为6.6%。
(2)所提出的方法具有普适性。本研究采用DBSCAN算法分别对海棠果园和芒果园进行聚类分析,对海棠果园与芒果园聚类准确率分别达到了95.3%和98.2%;并且两种果园的点云后处理过程时间均在30 s内。
(3)基于Unity3D引擎开发三维虚拟果园,所添加的各项功能够流畅运行。通过对虚拟果园系统所具备的功能进行测试,本研究构建的三维虚拟果园能够满足作业机械轨迹绘制,并可进行历史作业过程的模拟;具备果园场景的全方位立体显示、距离测量、位置显示等功能。
(4)本研究的DBSCAN聚类算法,在株间郁闭严重的情况下聚类分割效果有待改进;后续可以优化聚类分割算法,提升分割精度,并将聚类分割作为虚拟果园的一项功能,实现在虚拟果园中的聚类分析。
1
郑永军, 陈炳太, 吕昊暾, 等. 中国果园植保机械化技术与装备研究进展[J]. 农业工程学报, 2020, 36(20): 110-124.

ZHENG Y, CHEN B, LYU H, et al. Research progress of orchard plant protection mechanization technology and equipment in China[J]. Transactions of the CSAE, 2020, 36(20): 110-124.

2
NI Z, BURKS T F. Three-dimensional dense reconstruction of plant or tree canopy based on stereo vision[J]. Agricultural Engineering International: CIGR Journal, 2018, 20(2): 248-262.

3
BERVEGLIERI A, TOMMASELLI A, LIANG X, et al. Vertical optical scanning with panoramic vision for tree trunk reconstruction[J]. Sensors, 2017, 17(12): ID 2791.

4
郑立华, 麦春艳, 廖崴, 等. 基于Kinect相机的苹果树三维点云配准[J]. 农业机械学报, 2016, 47(5): 9-14.

ZHENG L, MAI C, LIAO W, et al. 3D point cloud registration of apple trees based on Kinect camera[J]. Transactions of the CSAM, 2016, 47(5): 9-14.

5
BERVEGLIERI A, TOMMASELLI A, LIANG X, et al. Vertical optical scanning with panoramic vision for tree trunk reconstruction[J]. Sensors, 2017, 17(12): ID 2791.

6
NEWCOMBE R A, FITZGIBBON A, IZADI S, et al. KinectFusion: Real-time dense surface mapping and tracking[C]// 10th IEEE International Symposium on Mixed and Augmented Reality. Piscataway, New York, USA: IEEE, 2011:127-136.

7
罗秋, 李国柱, 杨勇, 等. 地面三维激光扫描技术在文物三维重建中的应用——以金殿为例[J]. 城市勘测, 2020(3): 112-116.

LUO Q, LI G, YANG Y, et al. Application of ground 3D laser scanning technology in the 3D reconstruction of cultural relics: A case study of the Golden Temple[J]. Urban Survey 2020(3): 112-116.

8
杨俊, 王志鹏, 祝凯, 等. 激光扫描与无人机一体化控制下的复杂体育场馆高精度三维重建[J]. 测绘通报, 2020(12): 21-26.

YANG J, WANG Z, ZHU K, et al. High-precision three-dimensional reconstruction of complex sports venues under the integrated control of laser scanning and unmanned aerial vehicles[J]. Bulletin of Surveying and Mapping, 2020(12): 21-26.

9
王国利, 李群, 杨学博, 等. 基于机载激光雷达点云的森林场景建模[J]. 北京建筑大学学报, 2021, 37(2): 39-46.

WANG G, LI Q, YANG X, et al. Modeling of forest scene based on airborne lidar point cloud[J]. Journal of Beijing University of Civil Engineering and Architecture, 2021, 37(2): 39-46.

10
WANG K, ZHOU J, ZHANG W H, et al. Mobile LiDAR scanning system combined with canopy morphology extracting methods for tree crown parameters evaluation in orchards[J]. Sensors, 2021,21(2): ID 33419182.

11
HU C, PAN S, ZHANG H, et al. Trunk model establishment and parameter estimation for a single tree using multistation terrestrial laser scanning[J]. IEEE Access, 2020, 8: 102263-102277.

12
INDIRABAI I, NAIR M V H, JAISHANKER R N, et al. Terrestrial laser scanner based 3D reconstruction of trees and retrieval of leaf area index in a forest environment[J]. Ecological Informatics, 2019, 53: ID 100986.

13
翟晓晓, 邵杰, 张吴明, 等. 基于移动LiDAR点云的树木三维重建[J]. 中国农业信息, 2019, 31(5): 84-89.

ZHAI X, SHAO J, ZHANG W, et al. 3D reconstruction of trees based on mobile LiDAR point clouds[J]. China Agriculture Information 2019, 31(5): 84-89.

14
ZHENG L, FU Z. BALM: Bundle adjustment for LiDAR mapping[J/OL]. arXiv 2020: 8.

15
ZHANG J, SINGH S. LOAM: Lidar odometry and mapping in real-time. in robotics[C]// Robotics: Science and Systems Conference(RSS). New York City, USA: RSS, 2014.

16
RUSU R B, MARTON Z C, BLODOW N, et al. Towards 3D point cloud based object maps for household environments[J]. Robotics And Autonomous Systems, 2008, 56(11): 927-941.

17
罗方燕. PCL库点云统计去噪算法的应用研究[J]. 现代计算机(专业版), 2016(26): 63-66.

LUO F. Application of statistical denoising algorithm of PCL library point cloud[J]. Modern Computer (Professional Edition), 2016(26): 63-66.

18
ZHANG W, QI J, WAN P, et al. An easy-to-use airborne LiDAR data filtering method based on cloth simulation[J]. Remote Sensing, 2016, 8(6): ID 501.

19
PROVOT X.Deformation constraints in a mass-spring model to describe rigid cloth behaviour[C]//Graphics interface. Mississauga, Canada: Canadian Information Processing Society, 1995: 147.

20
肖冰. 激光和多波束点云水陆地形数据整合研究[D]. 泰安: 山东农业大学, 2018.

XIAO B. Research on the integration of land and water terrain data of laser and multi-beam point cloud[D]. Taian: Shandong Agricultural University, 2018.

21
ESTER M, KRIEGEL HP, SANDER J, et al. A density-based algorithm for discovering clusters in large spatial databases with noise[C]// 2nd International Conference on Knowledge Discovery and Data Mining. Broadway, New York, USA: ACM Digital Library, 1996: 226-231.

22
胡凡成. 基于Unity 3D的实时数据驱动数字化车间研究[D]. 长沙: 湖南大学, 2018.

HU F. Real-time data-driven digital shop floor research based on Unity 3D[D]. Changsha: Hunan University, 2018.

Outlines

/