CN104457626B - 一种基于激光雷达点云技术的植物叶面积指数测定方法 - Google Patents

Abstract

本发明涉及一种基于激光雷达点云技术的植物叶面积指数测定方法，与现有技术相比解决了叶面指数测定效率低、环境要求高的缺陷。本发明包括以下步骤：获取植物冠层的三维激光点云数据；点云数据的分割；植物边缘轮廓之内的点云图分类；对所有植物分别计算其半侧植物的体积；统计击中点的点数；统计击穿点的比例；利用反演计算出叶面积指数。本发明利用地面激光雷达测量系统生成的三维点云数据，结合计算机图形学技术，从点云图像中获取叶面积指数计算所需的相关参数并进行多元线性回归处理，从而计算叶面积指数。

Description

一种基于激光雷达点云技术的植物叶面积指数测定方法

技术领域

本发明涉及植物叶面积指数测定技术领域，具体来说是一种基于激光雷达点云技术的植物叶面积指数测定方法。

背景技术

叶面积指数是指一块土地上作物叶片的总面积与土地面积的比值，是生态学研究的关键参数之一，具有重要的理论价值与实际价值。目前树木的页面指数的测量方法大多为间接法，即利用仪器设备测量与叶面积指数相关的其他参数，再通过计算求得叶面积指数，而并非直接依照叶面积指数的定义进行测量。例如专利号为CN200810166822的文件提出利用统计落叶的相关参数来求取大面积林地的叶面积指数，但是这种方法无法在树木的生长期进行统计，且对测量环境要求苛刻(风速，树叶含水率)也很费时费力。

另外一种叶面积指数的测量方法主要是依靠间接光学模型测量法，即，主要研究空隙率，即冠层内太阳辐射未被截取的概率，从而出现了一系列基于空隙率分析的冠层叶面积指数分析仪器，例如专利号为CN201010136947的文件提出采用冠层分析仪。然而这种测量方法在每次测量前均需要做大量准备工作包括校对仪器、测量时间的选择(仅限于正午)、光照强度的限制、方位角校正等。每次测量准备时间长，冠层分析仪拍照之后需要在计算机中对图像进行后处理，无法实时获取测量结果，因此有必要设计一种新型的测量方法以满足智能农业中边测量边施药的新需求。

虽有少量采用直接法测量的技术，如专利号为CN201110169182和CN201210004979的文件，利用激光传感器测量树叶的层数，从而计算叶面积指数。但是这种方法存在激光的发射端发射激光穿过叶片到达接收端这一过程，由于机械结构的限制很难做得够大以方便测量高大的树木，在每次测量过程中都得按照每个树冠的形状调整测量时所走的路径，测量速度慢，效率不够高。

如何开发出一种测量效率高、适用环境广的叶面积测定方法已经成为急需解决的技术问题。

发明内容

本发明的目的是为了解决现有技术中叶面指数测定效率低、环境要求高的缺陷，提供一种基于激光雷达点云技术的植物叶面积指数测定方法来解决上述问题。

为了实现上述目的，本发明的技术方案如下：

一种基于激光雷达点云技术的植物叶面积指数测定方法，包括以下步骤：

获取植物冠层的三维激光点云数据，利用二维激光雷达和移动机器人沿着植物冠层的侧面进行扫描，获取每一帧图像拍摄时的三维坐标以及扫描目标点的三维空间坐标信息，将每一帧图像拼接起来，构成三维点云图像；

点云数据的分割，沿用采样时的坐标系，利用二维激光雷达和植物的位置关系，对植物边缘轮廓之外的点云图进行分类，将植物边缘轮廓之外的点云图从地面、背景杂物、天空无穷远处的点云中分离出来；植物边缘轮廓之内的属于地面、背景杂物、天空无穷远处的点云依然保留在点云图中；

植物边缘轮廓之内的点云图分类，对于植物边缘轮廓之内的所有点，以雷达行进直线到植物树干所在平面的倾角对应的距离为划分，小于此距离的标记为击中，大于此距离的标记为击穿；

对所有植物分别计算其半侧植物的体积V_half，通过切片的每个微元的面积S_i计算出每个切片的总面积S_frame，从而计算出半侧植物的体积V_half；

统计击中点的点数，对点云图中所有标记为击中的点进行计数，累加得到击中点的数目为I；

统计击穿点的比例，对点云图中所有标记为击穿的点进行统计，计算出击穿点的比例Z；

利用反演计算出叶面积指数，利用线性回归的方法求取同一地区同一植物在某段较短时间内的四个参数，计算同一地区同一植物的叶面积指数LAI。

所述的对所有植物分别计算其半侧植物的体积V_half包括以下步骤：

对于某一帧上所扫描的切片，计算切片的每个微元的面积，其计算公式如下：

其中S_OAB为二维激光雷达至植物的中心线所形成三角形的面积，S_OCD为二维激光雷达至植物的树冠最近处所形成三角形的面积，OA和OB为二维激光雷达与植物中心线之间的两根相邻的激光扫描线距离，OC和OD为二维激光雷达与植物树冠最近处之间的两根相邻的激光扫描线距离；θ为两根相邻的激光扫描线之间的夹角；

计算每个切片的总面积，其计算公式如下：

通过累加计算整个半侧植物的体积，其计算公式如下：

V_half＝ΣS_frame×Δs，其中Δs为激光雷达两次扫描帧之间前进的距离。

所述的统计击穿点的比例Z包括以下步骤：

对点云图中所有标记为击穿的点进行计数，累加得到击穿点的总数为z₀；

计算击穿点的比例Z，其计算公式如下：

Z＝z₀/(z₀+I)×100％。

所述的利用反演计算出叶面积指数包括以下步骤：

选择一定数量的样本植物，对样本植物进行人工取样测量或使用冠层分析仪测量，获取样本植物内单棵植物准确的叶面积指数LAI₀；

对样本植物建立多元线性回归模型，利用样本植物的叶面积指数LAI₀计算出体积因子A、击中因子B、击穿因子C、修正值D，其公式如下：

LAI₀＝A×V_half+B×I+C×Z+D，

其中，LAI₀为个体准确的叶面积指数、I为击中点的数目、Z为击穿点的比例、V_half为半侧植物的体积；

利用体积因子A、击中因子B、击穿因子C、修正值D，对于同一片区域内的其他植物分别计算单棵植物的叶面积指数LAI，其计算公式如下：

LAI＝A×V_half+B×I+C×Z+D，

其中V_half为需计算的单棵植物的半侧植物体积、I为需计算的单棵植物的击中点的数目、Z为需计算的单棵植物的击穿点的比例。

有益效果

本发明的一种基于激光雷达点云技术的植物叶面积指数测定方法，与现有技术相比利用地面激光雷达测量系统生成的三维点云数据，结合计算机图形学技术，从点云图像中获取叶面积指数计算所需的相关参数并进行多元线性回归处理，从而计算叶面积指数。

本发明提供了一种实时、准确、间接地估算树木叶面积指数的方法。通过提取单株树木的三维的激光点云图像中的与树木冠层相关的参数，来实现利用激光雷达点云数据进行叶面积指数的估算与反演，简洁高效(相比光学测量方法)且不会对树木结构和辐射特性造成不良影响(相比人工测量)，而且还可以永久性地记录树木的三维结构特征，有利于将来的深入研究。可应用于大范围内同一时间同一地区的同种植物上，测量速度快，无需针对每棵树单独调整，且不依赖外界光照，日夜均可工作。在测量过程结束后可以实时获取叶面积指数，在精准农业领域有着巨大的潜在应用价值。

附图说明

图1为本发明的方法流程图；

图2为本发明中二维激光雷达的测定原理图；

图3为叶面积指数与半侧植物体积的相关性展示图；

图4为叶面积指数与击中点数的相关性展示图；

图5为叶面积指数与击穿点比例的相关性展示图；

其中，1-二维激光雷达、2-植物。

具体实施方式

为使对本发明的结构特征及所达成的功效有更进一步的了解与认识，用以较佳的实施例及附图配合详细的说明，说明如下：

如图1所示，本发明所述的一种基于激光雷达点云技术的植物叶面积指数测定方法，包括以下步骤：

第一步，获取植物冠层的三维激光点云数据。利用二维激光雷达和移动机器人沿着植物冠层的侧面进行扫描，获取每一帧图像拍摄时的三维坐标以及扫描目标点的三维空间坐标信息。二维激光雷达和移动机器人对植物进行扫描，如图2所示，其拍摄的每一帧图像为透视扫描的切片图像，非日常视频中的帧图像。将每一帧图像(切片图像)拼接起来，构成三维点云图像。三维激光点云数据是利用地面二维激光扫描仪在平动中进行扫描所获取的树木冠层的点云，其中包含了扫描目标点的空间几何信息以及每个点的空间位置坐标信息，对每帧扫描图像进行拼接，从而获取对一棵树木扫描的完整点云图，作为后续步骤进行计算分析的数据源。构造三维点云图像在本发明中是通过现有技术中的二维激光雷达进行切片，再进行组合，现有技术中还有三维激光雷达，可以直接对植物进行三维点云图像的构成，其后续所涉及的步骤也与本发明不同。

第二步，点云数据的分割。沿用采样时的坐标系，利用二维激光雷达和植物的位置关系，对植物边缘轮廓之外的点云图进行分类。将植物边缘轮廓之外的点云图从地面、背景杂物、天空无穷远处的点云中分离出来，这部分点云图无关植物叶面积指数，因此无需进行处理，直接分离剔除。植物边缘轮廓之内的属于地面、背景杂物、天空无穷远处的点云依然保留在点云图中，待后续进行处理。

第三步，植物边缘轮廓之内的点云图分类。对于植物边缘轮廓之内的所有点，以雷达行进直线到植物树干所在平面的倾角对应的距离为划分，区分出击中和击穿。小于此距离的标记定义为击中，击中表示该点在植物冠层上；大于此距离的标记定义为击穿，击穿表示该点不属于植物冠层。在对于树木的叶面积指数的研究中，我们注意到了叶面积指数作为一个真实反映树木冠层属性的物理量，与树木冠层的别的物理量之间存在一定的相关性。比如，通常情况下，树木的冠层体积越大，含有树叶就越多，而且体积和面积是三次方对二次方的关系，会导致叶面积指数的增加。同时，冠层内树叶的疏密程度也显然与叶面积指数的大小存在相关性。而击中点数目与击穿点比例就可以很好地分别从总数和比例上体现树冠内树叶的疏密程度。因此，我们可以基于击中点数目、击穿点比例和冠层体积来进行植物叶面积的准确计算。

而击中点数目、击穿点比例和冠层体积的获取，则通过二维激光雷达切片出的结果进行计算，在实际应用中，当二维激光雷达在移动机器人沿着植物冠层的侧面进行扫描的过程中，就已得到每个切片的击中点数量、击穿点数量和微元面积(S_OAB和S_OCD的相关数据信息)，再利用下面的步骤进行计算。

第四步，对所有植物分别计算其半侧植物的体积V_half。通过切片的每个微元的面积S_i计算出每个切片的总面积S_frame，从而计算出半侧植物的体积V_half，以用于后期的反演计算。其具体步骤如下：

(1)对于某一帧上所扫描的切片，计算切片的每个微元的面积，其计算公式如下：

如图2所示，其中S_OAB为二维激光雷达至植物的中心线所形成三角形的面积，S_OCD为二维激光雷达至植物的树冠最近处所形成三角形的面积，OA和OB为二维激光雷达与植物中心线之间的两根相邻的激光扫描线距离，OC和OD为二维激光雷达与植物树冠最近处之间的两根相邻的激光扫描线距离；θ为两根相邻的激光扫描线之间的夹角。

(2)计算每个切片的总面积，其计算公式如下：

(3)通过累加计算整个半侧植物的体积，由于二维激光雷达所扫描的是植物的半侧，因此在此针对所有切片的总面积进行累加，从而得到整体半侧植物的体积。其计算公式如下：V_half＝ΣS_frame×Δs，其中Δs为激光雷达两次扫描帧之间前进的距离。

第五步，统计击中点的点数，对点云图中所有标记为击中的点进行计数，累加得到击中点的数目为I。由于点云图为多个切片构成的三维点云图，在此所统计的点云图中所有击中标记，为统计整个半侧植物的击中点总和，即单个植物的半侧击中点总和。

第六步，统计击穿点的比例，对点云图中所有标记为击穿的点进行统计，计算出击穿点的比例Z。与击中点点数统计同理，在此所统计的点云图中的击穿点比例，也是整个半侧植物的击穿点比例，即单个植物的半侧击穿点比例。其具体步骤如下：

(1)对点云图中所有标记为击穿的点进行计数，累加得到击穿点的总数为z₀，此时得到的为整个半侧植物的击穿点总数。

(2)计算击穿点的比例Z，此时得到的为整个半侧植物的击穿点比例，其计算公式如下：

Z＝z₀/(z₀+I)×100％。

第七步，利用反演计算出叶面积指数。利用线性回归的方法求取同一地区同一植物在某段较短时间内的四个参数(因子)，再通过这些参数计算同一地区同一植物的叶面积指数LAI。

(1)选择一定数量的样本植物，对样本植物进行人工取样测量或使用冠层分析仪测量，获取样本植物内单棵植物准确的叶面积指数LAI₀。人工取样测量方法和冠层分析仪测量方法获取叶面积指数LAI₀，为现有技术的方法，可以准确获取单棵植物的叶面积指数LAI₀。在此选择样本植物的数量依据解多元线性回归模型方程式来做，例如：求解A、B、C、D这四个值，从理论上来讲，4棵植物的叶面积指数LAI₀就可以完成求解过程，但实际应用中一般选择获取10棵左右树木的叶面积指数LAI₀，样本选择多一些，则相当于采用最小二乘法等方法,求个最接近的解即可。

(2)对样本植物建立多元线性回归模型，利用样本植物的叶面积指数LAI₀计算出体积因子A、击中因子B、击穿因子C、修正值D，其公式如下：

LAI₀＝A×V_half+B×I+C×Z+D，

其中，LAI₀为个体准确的叶面积指数、I为击中点的数目、Z为击穿点的比例、V_half为半侧植物的体积。在上一步中选择一定数量的样本植物，分别将其叶面积指数、击中点数、击穿比例和半侧植物的体积代入公式，组成方程组，从而求解出体积因子A、击中因子B、击穿因子C和修正值D。

(3)上述四个参数在实际应用中，不同植物均不相同。但是对于同一地区同一植物，则可以直接利用这四个参数，加上利用二维激光雷达对单棵植物所检测出的半侧植物的体积V_half、击中点数I、击穿比例Z，可以直接计算出所对应的单棵植物的叶面积指数。

利用体积因子A、击中因子B、击穿因子C、修正值D，对于同一片区域内的其他植物分别计算单棵植物的叶面积指数LAI，其计算公式如下：

LAI＝A×V_half+B×I+C×Z+D，

其中V_half为需计算的单棵植物的半侧植物体积、I为需计算的单棵植物的击中点的数目、Z为需计算的单棵植物的击穿点的比例。

叶面积指数与V_half、I、Z这3个变量分别存在一定的相关性，但是不明显，有些变量的相关性不强，有些变量的数据差距比较大。如图3所示，叶面积指数与V_half的相关性中，LAI＝0.0014434*V_half+0.20836。其相对的RMSE＝0.18241，R²＝0.47398。

其中，RMSE为均方根误差亦称标准误差，其定义为i＝1，2，3，…n。在有限测量次数中，均方根误差常用下式表示：√[∑di^2/n]＝Re，式中：n为测量次数；di为一组测量值与真值的偏差。R²拟合优度(Goodness of Fit)是指回归直线对观测值的拟合程度。

如图4所示，叶面积指数与I的相关性中，LAI＝5.0422*I+0.21424，其相对的RMSE＝0.18305，R²＝0.46819。

如图5所示，叶面积指数与Z的相关性中，LAI＝-0.0064961*I+0.61955，其相对的RMSE＝0.16579，R²＝0.59958。

但是在将三个参数充分考虑之后，LAI计算如下：

LAI＝0.30177*V_half+-0.005446*I+0.00098752*Z+0.44049，

其相对的RMSE＝0.15572，R²＝0.68426。

由此可以明显地看出，RMSE和R²均有明显提高，其中R²越大，说明相关性越好；RMSE越小，说明偏差越小，总体上计算的精确度有了相当大程度的提高。

实验表明，本发明提供了间接观测树木冠层叶面积指数的有效方法，无论树种和光照条件，都可以利用本方法进行快速测量。在实际应用中，先针对一定数量的植物采集准确的叶面积指数，从而计算出体积因子A、击中因子B、击穿因子C和修正值D，将这四个参数应用于二维激光雷达，二维激光雷达在移动机器人(拖拉机)带动下对大面积的植物进行扫描，二维激光雷达在扫描的同时采集并计算出半侧植物的体积V_half、击中点数I和击穿比例Z，直接联同体积因子A、击中因子B、击穿因子C和修正值D应用于线性回归模型的公式中，直接快速地计算出植物的叶面积指数。二维激光雷达扫描完一棵，则直接计算出一棵。克服了传统的手工测量方法对于树木结构的破坏，以及光学测量对天空光照条件的要求与图像后处理的滞后性，提供了一种实时测量的新方法，特别适用于大规模的果园环境，在智能农业中有着巨大的应用价值。

以上显示和描述了本发明的基本原理、主要特征和本发明的优点。本行业的技术人员应该了解，本发明不受上述实施例的限制，上述实施例和说明书中描述的只是本发明的原理，在不脱离本发明精神和范围的前提下本发明还会有各种变化和改进，这些变化和改进都落入要求保护的本发明的范围内。本发明要求的保护范围由所附的权利要求书及其等同物界定。

Claims (3)

1.一种基于激光雷达点云技术的植物叶面积指数测定方法，其特征在于，包括以下步骤：

11)获取植物冠层的三维激光点云数据，利用二维激光雷达和移动机器人沿着植物冠层的侧面进行扫描，获取每一帧图像拍摄时的三维坐标以及扫描目标点的三维空间坐标信息，将每一帧图像拼接起来，构成三维点云图像；

12)点云数据的分割，沿用采样时的坐标系，利用二维激光雷达和植物的位置关系，对植物边缘轮廓之外的点云图进行分类，将植物边缘轮廓之外的点云图从地面、背景杂物、天空无穷远处的点云中分离出来；植物边缘轮廓之内的属于地面、背景杂物、天空无穷远处的点云依然保留在点云图中；

13)植物边缘轮廓之内的点云图分类，对于植物边缘轮廓之内的所有点，以雷达行进直线到植物树干所在平面的倾角对应的距离为划分，小于此距离的标记为击中，大于此距离的标记为击穿；

14)对所有植物分别计算其半侧植物的体积V_half，通过切片的每个微元的面积S_i计算出每个切片的总面积S_frame，从而计算出半侧植物的体积V_half；

15)统计击中点的点数，对点云图中所有标记为击中的点进行计数，累加得到击中点的数目为I；

16)统计击穿点的比例，对点云图中所有标记为击穿的点进行统计，计算出击穿点的比例Z；

17)利用反演计算出叶面积指数，利用线性回归的方法求取同一地区同一植物在某段较短时间内的四个参数，计算同一地区同一植物的叶面积指数LAI；

所述的利用反演计算出叶面积指数包括以下步骤：

171)选择一定数量的样本植物，对样本植物进行人工取样测量或使用冠层分析仪测量，获取样本植物内单棵植物准确的叶面积指数LAI₀；

172)对样本植物建立多元线性回归模型，利用样本植物的叶面积指数LAI₀计算出体积因子A、击中因子B、击穿因子C、修正值D，其公式如下：

LAI₀＝A×V_half+B×I+C×Z+D，

其中，LAI₀为个体准确的叶面积指数、I为击中点的数目、Z为击穿点的比例、V_half为半侧植物的体积；

173)利用体积因子A、击中因子B、击穿因子C、修正值D，对于同一片区域内的其他植物分别计算单棵植物的叶面积指数LAI，其计算公式如下：

LAI＝A×V_haif+B×I+C×Z+D，

其中V_half为需计算的单棵植物的半侧植物体积、I为需计算的单棵植物的击中点的数目、Z为需计算的单棵植物的击穿点的比例。

2.根据权利要求1所述的一种基于激光雷达点云技术的植物叶面积指数测定方法，其特征在于：所述的对所有植物分别计算其半侧植物的体积V_half包括以下步骤：

21)对于某一帧上所扫描的切片，计算切片的每个微元的面积，其计算公式如下：

$S_{ABCD}=S_{OAB}-S_{OCD}=\frac{1}{2}\×\sin \mathrm{\θ}\×(OA\×OB-OC\×OD)=S_i,$

其中S_OAB为二维激光雷达至植物的中心线所形成三角形的面积，S_OCD为二维激光雷达至植物的树冠最近处所形成三角形的面积，OA和OB为二维激光雷达与植物中心线之间的两根相邻的激光扫描线距离，OC和OD为二维激光雷达与植物树冠最近处之间的两根相邻的激光扫描线距离；θ为两根相邻的激光扫描线之间的夹角；

22)计算每个切片的总面积，其计算公式如下：

$S_{frame}=\underset{i=1}{\overset{i_{max}}{\Σ}}{S}_i;$

23)通过累加计算整个半侧植物的体积，其计算公式如下：

V_half＝ΣS_frame×Δs，其中Δs为激光雷达两次扫描帧之间前进的距离。

3.根据权利要求1所述的一种基于激光雷达点云技术的植物叶面积指数测定方法，其特征在于，所述的统计击穿点的比例Z包括以下步骤：

31)对点云图中所有标记为击穿的点进行计数，累加得到击穿点的总数为z₀；

32)计算击穿点的比例Z，其计算公式如下：

Z＝z₀/(z₀+I)×100％。