杨 萍，郭志成
（兰州理工大学 机电工程学院，甘肃 兰州 730050）
花椒（Zanthoxylum bungeagumMaxim）是“八大调味品”之一，且具有一定的中药价值。随着我国花椒市场需求的增加，花椒树的培养面积日益扩大，人工采摘过程中所面临的诸多问题也日渐明显。首先，花椒树高2 ～5 m，目前花椒的采摘全靠人工，采摘工作量大。花椒树表皮多带有皮刺，对采摘人员容易造成伤害。其次，由于花椒种植环境复杂、不规范，特别是甘肃省花椒种植地多为山地，人工采摘花椒效率非常低，所以研发1 款能够自主采摘花椒的机器人迫在眉睫。
目前，国内基于视觉的花椒采摘机器人的研发处于初级阶段。其难点主要是花椒果实质量较轻且生长方向不确定；采摘时主要以串为单位而每串花椒的果实颗粒的生长较分散，这就给识别和定位花椒提出挑战。文献［1］在吸气式花椒采摘机械臂的基础上设计了基于红绿图像差异的方法进行图像识别，然后通过计算花椒果实串的形心和半径最终确定花椒的采摘点。文献［2］设计的花椒采摘视觉系统利用的是花椒图像在RGB 和HSV 颜色空间中各颜色分量的特点，采用阈值分割方法将花椒目标与背景分割，然后提取目标轮廓，结合形心偏差计算方法和惯性主轴最终确定出采摘点。
笔者针对识别、定位花椒过程中，所面临的离散型花椒串难以定位质心和花椒母枝生长方向不确定难以定位采摘点等问题。提出了摘花椒的单目视觉系统，其作用不但能够准确识别出花椒果实，同时也能够获取采摘点的深度信息和定位采摘点三维世界坐标。最后通过实验验证算法具有一定的可行性。
1 花椒果实的识别1.1 识别算法阈值分割法［3］是经典且常用的图像分割算法，由于成熟花椒呈现较明显的红色，并且花椒种植农田间的背景多为树木、树叶、土地和少量天空、杂草，这与成熟目标花椒呈红色调的特点具有较明显的颜色差异。所以通过设置RGB 颜色空间下的红色R值将目标提取出来。设分割后图像函数为h(x,y)、原图像为f(x,y)、R为阈值，则阈值分割算法满足下列表达式：
为了观察不同R值的分割结果，将R值从0.3以间隔为0.01 开始取值，最终确定R分量的阈值完成识别任务。经过多次实验确定R分量的阈值为0.6，图像识别结果如图1（b）所示。
模糊C 均值算法［4］是多元统计分析中的1 种无监督分割模式算法，其基本思想是用迭代的方式计算隶属度函数和簇中心同时达到最优的过程。让目标函数达到最小：
式(2)中，m为加权指数；若m值过大则聚类效果差，m值过小则算法变为硬C 均值聚类算法。uij为子集xj与子集xi之间的隶属关系，xi、xj为数据组X={x1,x2,…,xc}的子集。
为了使算法有很好的区分效果，保证模糊数组中每点到聚类中心的价值函数Jm达到最小，需要确定模糊组C的个数，若C值过大则目标被分到不同的类，若C值过小则聚类结果不够精确。因此为保证聚类结果的准确性，经过多次实验最终确定聚类数C=3,聚类结果如图1 所示，当C=2 时算法提取的结果中包括一些树枝和树叶。当C=3 和C=4 时算法能将花椒提取出来。当C=5 时由于将花椒聚成不同类，所以花椒的提取效果变差。为了减少计算量所以选择C=3，从图中可以看出算法能很好地将花椒目标分割出来。
图1 模糊C 均值算法Fig.1 fuzzy c-means algorithm
K-means 聚类算法［5］与FCM 算法同属于无监督学习算法，是把n个数据划分到K个簇，通过迭代计算来让函数J达到最小。
当数据点被划分到相似类的时候为rnk=1，否则rnk=0，μk为同一簇的均值。算法就是让同一类中各数据间的相似性较高；不同类间的相似性较低。K-means 聚类算法关键是聚类数K的选取，通过K=2,K=3,K=4 进行实验，如下图2 所示，由于花椒田背景比较复杂同时包括：农田、农作物、天空、房屋等多种背景，从图中可以看出当K=2 时目标无法被提取，当K=3，K=4 时提取效果基本一致，为提高算法速度并通过多次实验最终确定K=3。
图2 K 值选择Fig.2 K value selection
最大熵分割算法是以图像的熵为准则，用统计测量方法对目标图像进行分割，熵定义为：
用阈值q(0＜=q＜K-1)将图像分割为C0和C1区域，其概率密度函数为：
上式中p0(q),p1(q)代表的是近景和远景的累计概率，和为1。p0(q)和p1(q)对应的熵如下：
在该阈值下，图像总熵为:
图3（e）为最大熵阈值算法的分割结果。
通过对比上述算法，如图3 所示，可知K-means聚类算法可以更完整提取图像中的花椒。
图3 算法对比Fig.3 Algorithm comparison
1.2 花椒果实的提取由于K-means 聚类算法无法自动判别目标和背景。所以根据成熟花椒红色调较明显的特点，先把花椒图像从RGB 转换到Lab 颜色空间［6］，再选Lab 颜色空间下a分量最大的子类为目标,从而完成花椒的自动提取。因图像在拍摄、传输和处理过程中，存在噪声和失真等干扰，为提高图像的清晰度和便于分割，应先对花椒图片进行滤波、降噪等预处理。
步骤1：用摄像机对目标进行图像采集。
步骤2：使用模数转换器将采集到的图片信息模拟量转换成数字量。
步骤3：将得到的数字化图片进行滤波，去噪，尺寸修改等预处理工作。
步骤4：将花椒图像从RGB 颜色空间转换至Lab 颜色空间。
步骤5：运用K-means 聚类算法对花椒图像进行识别，设置K=3。
步骤6：将提取出的花椒结果转化成2 值图像。
步骤7：将2 值化的图像进行2 次降噪。
步骤8：运用Canny 算子提取2 值图像边缘。
步骤9：运用2 值填充方法对下个图进行2 次填充。
步骤10：运用1 阶质心矩算法定位花椒质心，图中蓝色星型点为质心位置。
步骤11：输出结果。
如图4 所示，图中的图片来源于甘肃省定西市临洮县龙门镇魏家坪农田中的花椒照片。从图4 中可以看出图片背景包括花椒、花椒树叶、农田、房 屋、其他农作物等，算法均能很好地将花椒目标提取出来。
图4 目标提取Fig.4 Target extraction
1.3 花椒质心的提取根据花椒串的生长特点具有离散性的特点，可知花椒果实的二值图像会有2 种情况；一是多个独立的连通区域，二是能够形成完整的连通区域。针对第1 种情况本文提出了根据求平面内质点系的质心算法，确定平面内离散花椒颗粒的质心。
针对第1 种情况需将花椒果实二值图像等效为平面内的质点系。首先将已提取出的花椒果实图像二值化求出连通区域，然后将此二值图像内的各独立花椒果实颗粒，如图6(b)所示，等效为平面内的质点系如图5 所示，求出二值图像内花椒颗粒的质心坐标xi、yi，和单个独立连通区域内的所有像素和作为质量mi，如图6(b)所示。然后根据平面内质点系的质心算法，将xi、yi、mi带入公式（12）和公式（13）求出花椒串的质心坐标。图6(c)中红色圆点为最终求得的花椒串质心。
图5 质点系的质心Fig.5 Center of mass of particle system
若花椒果实串可以形成1 个完整的连通区域，则可以直接计算花椒果实串的质心坐标，如图6（d—e）所示，图中红色圆点为质心。
图6 提取目标花椒之心Fig.6 Extraction of the heart of the target Zanthoxylum bangeanum
2 花椒结果母枝的提取和采摘点的定位2.1 花椒结果母枝的提取通过分析花椒图像在Lab 颜色空间下a分量直方图可知，花椒和花椒结果母枝相对于叶子和背景在直方图中灰度值的对比度较大［7］。如下图7 所示，从图中可以看出a分量灰度直方图具有明显的双峰特点。
图7 花椒图像的a 分量直方图Fig.7 A component histogram of Zanthoxylum bangeanum image
所以，可先采用Otsu 算法与最大熵阈值的线性组合［8］计算出花椒a分量的阈值将叶子和背景去除如图7，得到仅有花椒和花椒结果母枝如图8（b），采用形态学腐蚀和膨胀处理滤掉图像中的噪声和空洞图8（c），将K-means 聚类算法提取出的图像，如图8(d)所示，转换为二值图像，如图8(f)所示，然后经过1 次膨胀运算，再将其与已分割出的花椒串的二值图像，如图8(e)所示，进行图像相减，实现花椒结果母枝区域的提取［9］如图8(g)。
图8 提取结果母枝Fig.8 Extraction of fruiting mother branch
2.2 定位采摘点分析花椒的生长特点可知，花椒采摘点所在的结果母枝直线的延长线通过花椒的质心，或与花椒质心的距离最短。因此，笔者用点到直线的最短距离为限定条件确定最终采摘点所在的线段［10］，式（16）中函数是花椒质心到各直线距离，目标采摘点在min(Di)所对应的直线段上，取该线段的中点作为采摘点［11］。如图9 示意图所示，图中xC、yc为质心像素坐标。
先将花椒结果母枝进行Canny 边缘提取，然后将边缘图像进行霍夫直线检测，如式(15),并记录线段端点坐标［12］。
式中x、y—图像像素坐标；
ρ—原点到直线的距离；
θ—直线法线和X轴的夹角。
设上式(14)检测得出的线段端点为p1=(xi1，yi1)和p2=(xi2，yi2)，设两点直线方程为Li(x,y)=0，如图9中的红色线段。
根据求解点到直线距离，求出各个直线到质心的距离Di，当D最小时直线为所需要直线。
图9 示意图Fig.9 Sketch map
采摘点提取结果如图10，其中红色直线的中点为花椒的采摘点。
图10 定位采摘点Fig.10 Locating of the picking point
2.3 获取采摘点的深度信息根据相机成像的基本原理：
式中f为焦距，u为物距，v为像距。
相机在保证内参数和外参数不变的情况下，将相机沿光轴方向移动，对静止物体进行2 次拍照。根据透镜的2 次成像原理，如图11 所示，经过几何关系推导可得出式(18)和式(19)，式中k为大于零的实数。
图11 透镜2 次成像Fig.11 Twice imaging of lens
根据薄透镜成像原理可知物体在成像时的放缩率是相同的，首先假设相机对目标物体拍照形成的曲边形的面积为Ai，相机沿光轴移动后对物体拍照图像上形成的面积为Bi，从图12 中可以看出Ai相似于Bi。设Ai的极坐标方程为：
则Bi的极坐标方程为：
图12 极坐标下任意曲边面积Fig.12 Area of any curved edge in polar coordinates
根据曲线积分求得Ai和Bi区域面积为：
可得SBi=k2SBi，令SA=S像1，SB=S像2，则有：S像2=k2S像1。所以下式成立;
从式（25）中可以看出u为照相机到目标物体的深度距离［13］。
基于上述理论对目标花椒进行同光轴2次拍照，由于拍摄的每组照片需要沿着相机的光轴移动距离d，但实际拍摄的角度和精确的距离难以控制，所以本文将相机固定在高精度直线导轨上且用机械的方式限定相机的前后移动位置。提取出2 张结果母枝图像如图13（c，d），分别用Surf 算法对两图像进行特征点提取见图13（e，f），并对两图像进行特征点匹配见图13（g），提取20 个特征最明显的特征点形成凸包如图13（h，i），计算出结果母枝凸包区域的面积S像1和S像2带入式（25），相机沿光轴移动距离d=100 mm，求出花椒的深度距离，如图13 所示，镜头到花椒的真实距离为105 mm，算法计算的距离为107.39 mm，经过多次实验验证 如表1 所示，从表中可以看出误差范围满足采摘的 需求。
图13 花椒母枝的深度计算Fig.13 Calculation of the depth of the mother branch of Zanthoxylum bungeanum
表1 深度计算误差分析Table 1 Error analysis of depth calculation
2.4 获取采摘点的三维世界坐标根据相机针孔成像的基本原理可知1 幅图像中的物体坐标与三维世界坐标之间的转化需要经过4个坐标系的转换，这4 个坐标系分别是：世界坐标系、相机坐标系、图像物理坐标系、图像像素坐标系。如图14 所示，图中可以得出物体在图像中的坐标与在三维世界中的坐标的转换公式（26）。
图14 坐标转换Fig.14 Transformation of coordinate
上式（26）中等式右端第1 个矩阵与第2 个矩阵的乘积称为相机的内参数矩阵，矩阵［R|T］称为相机的旋转平移矩阵又叫相机外部参数矩阵，等式最右边为物体的三维空间坐标向量，等式左边Zc为采摘点在相机坐标下的坐标值，x，y为采摘点在图片内的坐标值。笔者利用Matlab2017A 的App 中Camera calibrator 对相机做标定，得出相机的内部参数矩阵为：
以机器人的基坐标为基准，根据相机的摆放位置和偏转角度得出相机的外部参数矩阵为：
设采摘点的三维世界坐标为（Xw,Yw,Zw）带入式（26），深度以表1 中第1 组数据为例，求出采摘点的真实坐标为（-145.953，-211.607，-218.111）mm， 实际测量坐标为（145，210，216）mm，满足采摘需求。
3 实验结果本文花椒图像拍摄于甘肃省定西市临洮县龙门镇魏家坪，背景包含农田、天空、房屋和其他农作物。拍摄相机为Microsoft 公司生产的LifeCam studio 型号的相机，分辨率为1 920×1 080 像素。本文实验图像的拍摄来自于摄像机对3 颗花椒树上不同位置花椒的拍摄。实验结果，如图15 和表2 所示，从图中可以看出算法能够定位出采摘点的位置，但是当图像背景中包含许多花椒树的枝干时，会导致采摘点的定位失败。
图15 采摘点定位实验Fig.15 Picking point positioning experiment
由表2 可以看出，在算法具有一定的可行性。
表2 实验结果Table 2 Experiment result
4 结论笔者对花椒采摘机器人的视觉分割及定位技术进行了研究，首先通过对比4 种分割算法的性能，确定利用提取效果最好的K-means 聚类算法识别出花椒果实，根据花椒串具有离散性的特点提出利用平面内质点系模型计算其质心。然后结合Otsu 算法与最大熵阈值的线性组合提取出花椒结果母枝，利用单目摄像机获取结果母枝的深度数据。最后通过限定霍夫直线到质心最短距离确定出采摘点在图像中的坐标，再根据三维世界坐标与图像坐标的转换定位采摘点的三维世界坐标。