王浩云 闫茹琪 周小莉 马仕航 胡皓翔 徐焕良
(1.南京农业大学信息科技学院, 南京 210095; 2.南京农业大学农业工程博士后流动站, 南京 210031)
作物表型受基因和环境因素影响,反映植物结构及组成、植物生长发育过程等全部物理、生理、生化特征。植物果实外形参数很大程度上反映了作物生长环境、作物产量等,是作物表型的一个重要参数。因传统的接触式测量[1-2]存在效率低、不准确等诸多问题,非接触式测量方法得到飞速发展。早期,非接触式测量多基于二维图像处理,文献[3]采集葡萄不同时期的果实图像,利用机器视觉及图像处理,计算了葡萄果实不同生长时期的表型参数。文献[4]利用图像处理的方法,通过阈值分割得到了核桃仁的分割图像,采用像素统计方法得到核桃仁的面积,实现了核桃仁大小分级。文献[5]通过分析提取苹果RGB图像中单色、波长差等信息,利用形态学方法处理图像,提取了目标区域几何、灰度和色调频度等特征,采用多元线性和偏最小二乘回归预测了苹果的质量和糖度等指标。文献[6]运用图像处理方法提取了黄瓜的瓜长、把长、横径差、弓形高度等特征,实现了黄瓜精确分级等。这些方法虽然可以应用于果实分级、检测及监测果实生长态势,但单幅二维图像中所包含三维信息缺失,即使采用机器视觉方法进行处理,其获取全面的三维信息过程复杂,并伴随着相应的误差。
近年来,立体相机、深度相机和激光雷达等三维数据采集设备[7-8]的兴起,点云数据的出现为作物三维重建、三维测量等提供了新方法。文献[9]使用Kinect传感器获取角果期油菜分枝图像,使用基于KD-tree搜索的最近点迭代算法实现点云配准,实现了油菜枝的三维重建。文献[10]提出一种基于深度信息分割聚类的草莓冠层结构形态三维重建算法,为草莓采摘机器人果实识别过程中枝叶空间结构关系的构建提供了技术支持。文献[11]使用三维扫描仪获取植物叶片的点云数据,并通过点云数据进行预处理,对处理后的点云数据经网格化后得到了高精度的植物叶片模型。文献[12]使用Kinect相机与旋转台获取苹果各个角度的点云数据,通过点云配准及贴皮建立了苹果果实三维模型。文献[13]利用体着色三维重建的方法获取变态根点云数据,并提取变态根的构型参数进行曲面建模,实现了胡萝卜和马铃薯根系的三维可视化。
基于点云数据的植物三维重建和测量虽然有了长足发展,但方法多为从不同角度获取点云信息后进行配准拼接得到完整的作物。这种方法在室内单株植物摘取果实较适用,但在获取户外果实的多角度点云信息相对困难,相反获取单角度的局部点云数据较简单。为了实现从局部点云信息中直接获取果实外形指标,本文以苹果作为研究对象,对基于局部点云的苹果外形指标估测方法进行研究。采用椭球曲面方程的方法构建苹果几何模型,并计算苹果几何模型的高度、直径、体积,使用Kinect V2相机从任意角度获取点云数据,采用直通滤波法去除点云数据的背景,并用包围盒算法精简点云得到苹果局部点云数据后,采用粒子群算法将苹果局部点云数据与苹果模型进行空间匹配,并以遗传算法求解苹果最优匹配模型的参数,利用苹果最优匹配模型参数估测与其匹配的真实苹果的外形指标。
实验选取红富士苹果作为研究对象,测量了250个红富士苹果的直径、高度及体积3个参数作为其外形指标。采用游标卡尺测量苹果直径、高度。游标卡尺的测量精度为0.01 mm。苹果体积测量使用排水法,1 000 mL烧杯2个,500、100 mL标准量筒各1个,量筒精度为1 mL。
1.2.1点云数据采集
实验采用Kinect V2相机采集苹果局部点云数据。将Kinect相机倒挂于装有横臂、水平仪的三脚架上,将三脚架置于实验台上,在相机正下方的位置铺上黑色硬纸板作为背景。图1为Kinect V2相机安装实物图。点云数据采集过程中相机垂直于实验台,高度固定为65 cm,每次以单个苹果为研究对象,将苹果以任意角度置于实验台,苹果拍摄角度分为正立、侧躺、倒立,为方便表述,分别记为角度1、角度2、角度3,如图2所示。
图1 点云采集装置Fig.1 Point cloud acquisition device
图2 苹果拍摄角度Fig.2 Apple shooting angle
开发软件使用Matlab 2016b,其工具箱中发布的Image Acquisition Toolbox工具包可以将Kinect相机连接到计算机上并与Matlab进行搭接,可以使用Matlab直接连接Kinect相机采集彩色点云数据,保存为PLY格式文件。
1.2.2点云数据处理
图3 苹果点云数据Fig.3 Apple point cloud data
通过Kinect V2相机获得的苹果点云数据包含苹果的RGB颜色信息和X、Y、Z位置信息。点云的坐标系以Kinect深度相机为原点,单位m,精确度为0.001 m,如图3a所示。苹果原始点云包含实验台等冗余的背景信息和不稳定的噪声点,不利于后续数据处理。采用目前最简单有效的直通滤波法[14]去除背景环境数据得到苹果局部点云。由于拍摄设备较为精密,经处理后的苹果局部点云中数据量仍较大,但是一部分数据点所包含的信息量并不能为后期苹果局部点云与苹果模型匹配作贡献,即处理后的苹果局部点云中仍包含大量冗余数据点,采用包围盒算法[15-16]对苹果局部点云数据进行精简压缩。包围盒算法首先对点云数据建立一个最小的立方体包围盒,然后将这个包围盒分割成大小相等的立方体,将点云数据按照三维坐标归入到一个小立方体中,对每个小立方体找出其中最接近立方体中心的点,并移除其余点从而精简点云数据。原始点云经直通滤波算法去除背景、又经包围盒算法精简后得到只包含苹果部分的局部点云,如图3b所示。
1.3.1基于椭球曲面方程的苹果模型
苹果外形的最小外接曲面近似于一个椭球面,且上下两端有不同程度的凹陷。采用基于椭球变形的苹果果实造型方法[17],通过添加合适的凹凸干扰函数将椭球曲面变形实现苹果外形的造型。苹果外层形态的参数方程Q(u,w)的各分量为
(1)
其中
g(u,w)=g1(u,w)+g2(u,w)
(2)
g1(u,w)=-p1e-2u
g2(u,w)=-p2e2u
式中 (xQ,yQ,zQ)——曲面Q(u,w)上点的坐标值
a——椭球曲面与X轴交点到中心的距离
b——椭球曲面与Y轴交点到中心的距离
c——椭球曲面与Z轴交点到中心的距离
u、w——椭球曲面方程参数
g(u,w)——干扰函数
g1(u,w)——上端凹陷干扰函数
g2(u,w)——下端凹陷干扰函数
p1、p2——上、下凹陷程度控制参数
根据该果实造型方法所构建的苹果模型如图4所示。
图4 苹果几何模型Fig.4 Apple geometric model
1.3.2苹果几何模型尺寸计算
苹果曲面方程中有a、b、c、p1、p25个参数,这5个参数的取值决定了苹果模型的尺寸,为了使每个模型区分度更大,确定了苹果模型的外形参数,分别为苹果几何模型的高度(记为h),苹果模型的2个相互垂直的直径(分别记为D1、D2)以及苹果几何模型的体积(记为V)。图5为苹果几何模型在X=0时,沿Z轴方向的剖面图。图中A为苹果剖面图最高点,B为其最低点,C、D为图像与Y轴的交点。
图5 苹果几何模型剖面图Fig.5 Apple geometric model profile
苹果几何模型尺寸计算步骤如下:
(1)改变式(1)中u、w的值,找到最高点A和最低点B,点A与点B在Z轴方向的差值作为苹果几何模型的高度h。
(2)继续改变式(1)中u、w的值,找到苹果几何模型在Z=0,X=0时,与Y轴的交点C、D,点C与点D在Y轴方向的差值作为苹果几何模型的一个直径D1;同理找到苹果几何模型在Z=0,Y=0时,与X轴的交点,将其在X轴方向的差值作为苹果几何模型的另一个直径D2。
(3)利用积分法对曲面方程式(1)进行积分得到苹果模型体积作为苹果模型体积指标V。
基于Kinect相机实采的苹果点云数据是从任意角度获取的苹果部分点云,由于信息不完整,从中直接获取苹果的全部尺寸参数难度较大。本实验寻找与真实苹果尺寸最为接近的苹果几何模型并用其尺寸估测真实苹果外形指标参数。为了寻找与真实苹果尺寸最为接近的苹果几何模型,需要将苹果局部点云数据与离散的苹果几何模型进行空间匹配并寻找匹配效果最优的苹果模型。
1.4.1曲面点集的空间变换
设苹果局部点云为点集M={mi}(i=1,2,…,n),离散苹果几何模型点集为O={oj}(j=1,2,…,k)。苹果局部点云与离散苹果几何模型匹配即点集M与点集O两曲面点集之间的欧氏空间变换,因而需要寻找一个欧氏变换矩阵作用于其中一个点集,从而实现点集M到点集O的距离最小[18]。对点集M进行欧氏空间变换,公式为
(3)
(4)
(5)
T=[tx;ty;tz]
(6)
式中H——欧氏变换矩阵
R——3×3正交旋转矩阵
α、β、γ——点绕坐标轴X、Y、Z的旋转角
T——3×1平移矩阵
tx、ty、tz——沿坐标轴X、Y、Z的平移量
1.4.2基于粒子群算法的苹果点云空间匹配
为了使变换后的苹果局部点云点集M与离散的苹果模型点集O达到最大限度的重合,需寻找点集最佳变换位置。实验固定离散后的苹果几何模型点集O,移动苹果局部点云点集M,并选用粒子群算法[21]对点集M的移动位置进行评估,找到两点集之间距离最小时欧氏变换矩阵H。基于粒子群算法的苹果点云空间匹配方法实现了苹果局部点云的空间位置变换,使空间中处于不同位置的苹果局部点云与苹果几何模型重合。
基于粒子群算法的苹果点云空间匹配步骤:
(1)输入点集M、点集O。
(2)随机初始化粒子的位置xr=[txtytzαβγ](r=1,2,…,N)和速度vr=[1 1 1 1 1 1],其中N表示粒子数,N=50。
(4)计算点集M′与点集O之间的距离和E表示两个点集的重合度,其值越小重合度越好,距离和E计算公式为
(7)
(5)更新最优值。如果位置更新后粒子的距离和E小于原来粒子的距离和,则将更新后粒子位置作为该粒子的局部最优解pr。如果该位置更新后粒子的距离和小于原来的全局最小距离和,则将该粒子位置作为本次迭代的全局最优解pg,且令两点集间的最小距离和Emin=E。
(6)根据公式更新粒子速度和位置
(8)
式中c1——自我学习因子,取值为1.3
c2——群体学习因子,取值为1.7
rand()——[0,1]的随机数
(7)判断是否到达最大迭代次数,是则执行步骤(8),否则循环回到步骤(3),其中,最大迭代次数为200。
(8)输出两点集间的距离和Emin。
图6为匹配前苹果局部点云与离散的苹果几何模型点集在空间中的位置,图中红色点集为离散的苹果几何模型点集,蓝色为移动前苹果局部点云;根据上述步骤将苹果局部点云M与离散的苹果几何模型点集O进行匹配,匹配结果如图7所示。图7a为经粒子群算法后两点集匹配的结果,其中绿色点集为经粒子群算法移动后的苹果局部点云,即苹果局部点云经过粒子群算法后从原来蓝色点集位置变化到现在绿色点集位置;图7b为匹配结果放大图。
图6 匹配前位置关系Fig.6 Spatial positional relationship before matching
图8 苹果点云与不同几何模型匹配结果Fig.8 Apple point cloud with different geometric model matching results
图7 粒子群算法匹配结果Fig.7 Particle swarm optimization match results
1.4.3基于遗传算法的苹果几何模型参数寻优
利用粒子群算法实现了苹果局部点云与离散的苹果几何模型的匹配,并可得到苹果局部点云与该模型的最小距离和。苹果几何模型的大小由a、b、c、p1、p25个参数决定,不同的参数组合都会得到不同尺寸的苹果模型,同一个苹果局部点云与不同的苹果几何模型匹配,效果不同,如图8所示。可知,当匹配到的苹果几何模型尺寸过大或过小均不能与苹果局部点云完全重合,只有匹配到的苹果几何模型尺寸与真实苹果尺寸相近时才可以与苹果局部点云完全重合。寻找与真实苹果尺寸最为接近的苹果几何模型,即寻找与其局部点云匹配效果最优的苹果模型。本实验采用遗传算法[22]求解最优苹果几何模型参数。
基于遗传算法的苹果几何模型参数寻优步骤:
(1)输入苹果局部点云M。
(2)随机初始化种群,种群规模记为S=15;采用二进制编码法对染色体进行编码。其中,染色体代表了苹果几何模型的5个参数a、b、c、p1、p2,长度为28,如图9所示。其中前4位0001代表参数a,第5~8位0010代表参数b,第9~16位01010111代表参数c,第17~20位1010代表p1,后8位10111100代表p2。
(3)计算个体适应度,先将染色体解码得到5个参数,根据参数得到离散的苹果几何模型点集O,将点集M与点集O作为粒子群算法的输入,调用粒子群算法,将粒子群算法的输出结果Emin作为个体适应度fitness。其中染色体解码方式,首先根据图9将染色体分成5个二进制编码,然后将5个二进制编码转换为十进制,分别记为aret、bret、cret、p1ret、p2ret,然后根据下式解码。
(9)
图9 染色体二进制编码示意图
Fig.9 Chromosome binary coding diagram
(4)将最小fitness记为最佳适应度bestfitness,其对应个体为最佳个体bestchrom。
(5)对种群进行遗传操作,包括选择、交叉、变异。其中,选择操作采用轮盘赌法,交叉概率为0.75,变异概率为0.05。
(6) 调用粒子群算法计算新一代的个体适应度,并计算新一代最佳适应度newbestfitness和新一代最佳个体newbestchrom。
(7)如果newbestfitness小于bestfitness,则新一代的值newbestfitness、newbestchrom代替上一代的值bestfitness、bestchrom。
(8)判断是否到达最大迭代次数,若达到最大迭代次数5,是则执行步骤(9),否则返回到步骤(5)。
(9)将最佳个体的染色体解码得到苹果几何模型的5个参数a、b、c、p1、p2,并输出。
图11 苹果外形参数实测结果统计Fig.11 Apple volume measurement results statistics
根据上述步骤得到与真实苹果尺寸最为接近的苹果几何模型的5个参数a、b、c、p1、p2,根据参数建立相应的苹果几何模型,然后计算苹果几何模型的外形尺寸,即该苹果局部点云对应的真实苹果的外形参数估测值。基于局部点云的苹果外形指标估测方法流程图如图10所示。
图10 基于局部点云的苹果外形指标估测方法流程图Fig.10 Flow chart of apple shape index estimation method based on local point cloud
实验选取了250个红富士苹果作为实验样本,通过排水法测量苹果的体积,利用游标卡尺测量苹果的高度和2个相互垂直的直径。图11为苹果外形参数实测结果统计图。利用数学统计学方法对测量结果进行分析,结果如表1所示。由表1和图11可知,苹果的体积分布范围为200~700 mL,苹果的高度分布范围为65~100 mm,直径分布范围为70~110 mm,实验所用苹果外形参数分布范围广,说明实验样本大小不一,类型较丰富,避免了样本单一性,使实验更加具有说服性、普遍性。
分别获取了250个苹果在3个角度下(图2)的点云数据,共计750幅点云图像。将处理后的250个苹果3个角度下的点云数据作为基于局部点云的苹果外形指标估测方法的输入,分别估测其对应的苹果果实外形指标(苹果的体积V、高度h、直径D1和D2),并将估测值与苹果果实实测值进行线性回归分析,如图12~14所示。图中横轴真实值为实验人工测量的苹果外形参数;纵轴估测值为使用基于局部点云的苹果外形指标估测方法估测的苹果外形指标;R2表示回归直线对观测值的拟合程度,最大值为1,其值越接近1,说明回归直线对观测值的拟合程度越好;RMSE表示均方根误差(亦称标准误差),用来衡量观察值同真实值之间的偏差。表2为苹果外形指标估测结果。同时,利用数学统计法计算了各个苹果外形指标的平均估测误差,如表3所示。
表1 红富士苹果果实实测参数Tab.1 Measured parameter analysis of red Fuji apple fruit
(1)苹果尺寸大小不一、分布范围广,分析图12~14可知,在各个角度下估测的苹果各个外形指标与相应的真实值有较高的相关性,其线性回归拟合的R2均高于0.7,RMSE均在允许的偏差范围之内,证明该方法适用于各种尺寸的苹果在任意角度下拍摄的点云数据,说明方法具有普适性。
(2) 纵观图12~14中苹果体积参数估测结果,角度2下R2最大,且RMSE最小,角度3下R2最小且RMSE最大,综合R2与RMSE分析可知,在角度2下估测的体积参数结果最好,在角度3下体积估测值相对较差。
图12 角度1下苹果外形指标估测结果Fig.12 Estimated results of apple shape indicators under angle 1
图13 角度2下苹果外形指标估测结果Fig.13 Estimated results of apple shape indicators under angle 2
图14 角度3下苹果外形指标估测结果Fig.14 Estimated results of apple shape indicators under angle 3
表2 苹果外形指标估测结果Tab.2 Apple shape indicator estimation result
表3 苹果外形指标估测平均误差Tab.3 Average error of apple shape indicator estimation
(3)纵观图12~14中苹果高度估测结果,角度2下其R2最大且RMSE最小,相反,角度3下其R2最小且RMSE最大,综合R2与RMSE分析可知,在角度2下估测的高度参数结果最好,在角度3下高度估测值相对较差。
(4)纵观图12~14中苹果直径估测结果,角度2下两个直径的R2均较大且RMSE均较小,相反,角度3下两个直径的R2均较小且RMSE均较大,综合R2与RMSE分析可知,在角度2下估测的直径参数结果最好,在角度3下直径估测值相对较差。
(5)横观图12~14中估测的各个苹果外形指标结果,3个角度下体积参数的R2均高于0.9,两个直径参数R2均高于0.88,而高度的R2均偏小,说明体积估测结果最好,高度估测结果与其他指标估测结果相比较差;由以上分析可知,在角度2下苹果外形指标估测结果最好,而在田间作业时,机器人获取角度2的果实数据最为方便,说明该方法可以为田间作业中获取生长期的果实外形指标提供理论基础。
(6)由表3可知,苹果体积估测结果的平均误差不大于16.16 mL,高度估测结果的平均误差不大于2.92 mm,直径估测结果的平均误差不大于2.35 mm,估测结果的平均误差较小,在允许误差范围内,按照国内现行的苹果分级标准GB/T 10651—1989《鲜苹果》标准可知苹果直径是苹果分级的重要指标,该方法可以用于苹果分级中,有实用价值。实验结果表明,本文方法预测结果相对稳定,结果符合实际,具有实用性。
Kinect相机为RGB-D相机,可以同时采集彩色图像、深度图像和点云数据。本实验在采集250个苹果局部点云数据的同时采集了其在3个角度下的彩色图像。在角度1、角度3下的彩色图像只包含苹果果实直径,无法得到高度信息,且由2.2节分析可知在角度2下基于局部点云的苹果外形指标估测效果最好。因此,基于彩色图像的苹果外形指标估测使用角度2下的彩色图像进行对比实验。实验通过数字图像处理对彩色图像进行处理,得到苹果图像分割图,并利用最小外接矩形检测苹果。通过外接矩形的4个顶点坐标计算出外接矩阵的长和宽,该长和宽是苹果果实的高和直径的像素点数,根据张正友标定方法得到每个像素点代表的实际距离为0.56 mm。将每个像素数乘以像素点大小得到苹果的高度和最大直径。利用SVR支持向量机回归模型拟合体积参数,其中200个数据作为训练集,50个数据作为测试集。根据数学统计学方法计算结果的平均误差。将角度2下2种苹果外形指标估测结果进行对比,结果如表4所示。
表4 估测平均误差对比Tab.4 Estimated average error comparison
由表4可知,基于彩色图像的苹果外形参数估测结果平均误差均大于基于局部点云的苹果外形指标估测结果。说明,基于局部点云的苹果外形指标估测更加准确,符合实际。
(1)采用直通滤波法去除原始点云的背景信息,得到只包含苹果部分的点云数据,然后使用包围盒算法精简点云,得到了准确的苹果局部点云。
(2)采用遗传算法寻找苹果最优匹配模型并求解模型参数,其中匹配方法采用基于粒子群算法的苹果点云空间匹配,实现了苹果局部点云的位置变换,使空间中处于不同位置的苹果局部点云与苹果几何模型重合。
(3)本文算法估测的苹果外形指标结果与真实值具有较高的相关性,其线性回归分析的R2均高于0.7,且苹果体积估测的平均误差不大于16.16 mL,高度估测的平均误差不大于2.92 mm,直径估测的平均误差不大于2.35 mm,估测结果的平均误差较小,在允许误差范围内,说明该方法具有可行性、实用性。
(4)本文方法为机器人田间获取生长期果实的外形指标提供了方法,且在角度2下的苹果外形指标估测效果最好。