一种抗重叠的植物点云单叶片分割方法与流程
本发明涉及农业工程、植物表型学、计算机图形学领域，尤其是涉及一种抗叶片重叠的植物点云单叶片分割方法。
背景技术：
植物表型是受基因和环境因素决定或影响的,反映植物结构及组成、植物生长发育过程及结果的全部物理、生理、生化特征和性状。随着植物功能基因组学和作物分子育种研究的深入,传统的表型观测已经成为制约其发展的主要瓶颈,而高通量的植物表型研究是解决这一困境的有效途径。植物表型研究是对植物复杂性状(植物的生长、发育、耐受、抗性、结构、生理、生态、产量等)的综合评价，是对构成复杂性状的基本参数(根系形态、生物量、叶片特性、果实特性、产量相关性状、光合效率、生物和非生物胁迫响应等)的直接测量。对绝大多数植物而言，叶片构成了植物主要的表面形态和结构，能直接反映植物的生长状态。叶片特性(形态、纹理、颜色等)通常蕴含着影响生长的生物压力因素(植物疾病和害虫)或者非生物压力因素(干旱)。因此，自动叶片分割、以及在其基础上的识别、分类方法能够及时提供直接的生长状态监控以保障植物产量。自20世纪90年代初以来，研究人员开始使用多种成像工具监测和分析植物的生长情况，基于图像数据的叶片自动分割、识别以及分类方法是其中的研究重点与难点。植物通常具有复杂的空间结构，对单个叶片分割方法带来了极高的要求；而且，植物的品种千差万别，在某类植物上有效的算法通常无法用于其他类别的植物。所以设计准确及高效的通用单叶片分割算法依然是亟待解决的难题和挑战。
以数据类型分，植物叶片的自动分割方法一般基于2d图像或者3d点云。基于2d图像的叶片分割算法一般采用传统的图像处理、模式识别以及机器学习方法，该类方法普遍能对结构较为简单的植物(如拟南芥、烟草以及小麦等)实现单个叶片分割。然而，基于2d图像的叶片分割算法通常受限于复杂的环境因素以及理想的成像角度，当图像中存在大面积相互遮挡的叶片时，2d的叶片分割算法几乎束手无策，因此只能在少数莲座型以及标准植物上达到满意的效果。
3d点云弥补了2d图像深度信息不足的缺陷，很大程度上避免了因叶片相互遮挡而难以准确分割的窘境。然而，3d植物点云虽然蕴含着充足的表征信息，但是植物的叶片是成簇状分布的，相互重叠遮挡的现象极为严重，对拥挤的植物冠层点云进行单叶片自动分割是植物表型领域的大难题。
技术实现要素：
本发明的目的是：彻底解决复杂冠层中难以准确分割单叶片的难题。
为了达到上述目的，本发明的技术方案是提供了一种抗重叠的植物点云单叶片分割方法，其特征在于，包括以下步骤：
步骤s1：对植物点云数据进行预处理，滤除植物点云数据中存在的检测噪声以及非植物叶片的信息，从而获得植物冠层点云；
步骤s2：对步骤s1中获得的植物冠层点云进行3d联合滤波算子处理，3d联合滤波算子是由基于空间半径的离群点滤波器和基于曲面边界的滤波器共同组成的，3d联合滤波算子在三维空间中对点云进行2d图像中的腐蚀操作，通过对点云中叶片相互重叠遮挡的部分进行边缘腐蚀，将重叠的叶片分开；
经过3d联合滤波算子的处理，将植株冠层点云分为两部分，第一部分是所有叶片的叶片中心点云部分，是经过3d联合滤波算子处理后留下的部分；第二部分是所有叶片的叶片外缘点云部分，是经过3d联合滤波算子处理后滤除的部分；
步骤s3：将3d联合滤波算子去除的叶片外缘点云部分进行面片过分割，然后添加回步骤s2已经分割后的叶片中心点云，以标记后的叶片中心向外进行基于面片的区域生长，生长后的结果即为最终的植物冠层单叶片分割结果。
优选地，在所述步骤s2进行多次3d联合滤波算子的处理，以应对冠层中叶片复杂的重叠情况；
则在所述步骤3中，按照由内向外的顺序对每次滤除的叶片外缘点云部分进行面片过分割和面片区域生长，以达到最佳的效果。
优选地，所述3d联合滤波算子从下述5种滤波器中选择至少两种滤波器：
滤波器i：基于空间区域的滤波器，可滤除点云在空间坐标系x/y/z轴方向上存在的非叶片信息；
滤波器ii：基于半径的离群点滤波器，该滤波器以查询点xq为球心，r为球体的半径，即当球体内的点数超过阈值nthreshold时，xq在新点云中保留；
滤波器iii：基于统计的k-邻近点滤波器；
滤波器iv：颜色滤波器，根据rgb通道中g值的所占的比重区分叶片和非叶片信息；
滤波器v：降采样方法，即创建一个三维体素栅格，然后每个体素内用体素中所有点的重心来近似，从而降低点云复杂度。
优选地，所述步骤s2中对获得的植物冠层点云进行3d联合滤波过程，3d联合滤波算子由基于空间半径的离群点滤波器和曲面边界滤波器级联而成，包括以下步骤：
步骤s2-1：采用半径滤波器对点云进行滤波，其原理是以点云中任一点xq为球心，r为球体的半径，即当球体内的点数超过nthreshold时，xq在新点云中保留，半径滤波器适合于剔除稀疏离群点和边界点，用于分离交叉重叠的叶片区域；
步骤s2-2：采用边界滤波器对点云进行滤波，其原理是通过pca提取曲面边缘点，实现对曲面边缘点进行遍历，达到类似边缘腐蚀算子的效果；
步骤s2-3：若由双目立体视觉系统和kinect传感器获取的植物点云，仅需要一次3d联合滤波，即执行一次步骤s2-1及步骤s2-2；若由多视角成像方式获取的点云，则需要3～4次联合滤波，即重复步骤s2-1及步骤s2-2三至四次。
优选地，所述步骤s2-2中采用边界滤波器对点云进行滤波过程：
步骤s2-2-1：找出点云中查询点xi周围最近的k个点，建立点集
步骤s2-2-2：通过pca估计点集
步骤s2-2-3：查询点xi与点集
θj＝arccos<(xj-xi)uv,u>其中，-π≤θj≤π，j∈{1,...,k}(1)
步骤s2-2-4：对投影向量(xj-xi)uv与u的夹角集合θ＝{θj|j∈{1,…,k}}进行递增排序，若最大角度差max(θj+1-θj)>θthreshold，则xi是一个边缘点，在曲面边界滤波中滤除；
步骤s2-2-5：当遍历完点云中所有的点后，滤除所有边缘点就完成了一次边缘滤波。
优选地，为了达到满意的边缘滤波效果，在经过步骤s2-2-5后剩下的点云中重复多次步骤s2-2-1至步骤s2-2-4。
优选地，所述步骤s3中面片的区域生长包括以下步骤：
步骤s3-1：叶片中心区域点云χc预分割，经过一次或者多次联合滤波后，保留的冠层点云大部分为叶片中心区域，通过3d区域生长的方式对叶片中心区域点云χc进行预分割与标记；
步骤s3-2：对滤除的叶片外缘点云部分χb进行过分割；
步骤s3-3：基于面片标记生长的单叶片分割，采用从内而外的方式，将标记的叶片中心区域向面片化的叶片外缘部分进行基于面片的区域生长。
优选地，所述步骤s3-1中，叶片中心区域点云χc预分割包括以下步骤：
步骤s3-1-1：从滤波后的叶片中心区域点云χc中随机选择一个未标记的点，赋予一个新的叶片标记；
步骤s3-1-2：搜索步骤s3-1-1所选择的点一定半径内是否存在其他未标记的点，如果有则将这些点都赋予这个点相同标记，并将这些新搜索到的点压入一个队列中；
步骤s3-1-3：每次从步骤s3-1-2得到的队列中取出一个点，继续重复步骤s3-1-1及步骤s3-1-2直到队列为空，则认为对一片叶片的标记已经结束；
步骤s3-1-4：在叶片中心区域点云χc中继续随机寻找下一个未标记的点，如果有就赋予其一个新的叶片标记，再重复步骤s3-1-1、步骤s3-1-2和步骤s3-1-3，直到叶片中心区域点云χc中全部点都被标记为止。
优选地，所述步骤s3-2中，对滤除的叶片外缘点云部分χb进行过分割包括以下步骤：
步骤s3-2-1：计算待过分割部分所有点的空间特征；
步骤s3-2-2：基于空间特征生成一定数量的种子点，进行粗略的面片生成；
步骤s3-2-3：使用局部k均值聚类对粗面片区域进行改进，生成精细面片。
优选地，所述步骤s3-3中，基于面片标记生长的单叶片分割包括以下步骤：
步骤s3-3-1：遍历每个完成标记的叶片中心区域点云，若其中的点与某个小面片邻接，则将该叶片标号赋值给该小面片包含的所有点并渲染相应的标号色彩，直至所有叶片标记区域均不再扩大为止；
步骤s3-3-2：如果叶片外缘点云部分χb中还有面片未被标记，则该类型面片未与任何已标记叶片区域相邻接，所以该面片代表一个新叶片区域，对其赋予一个独立的叶片标记并渲染新色彩；
步骤s3-3-3：对每个新叶片区域进行基于广度优先搜索的面片标记生长，将中与其邻接但未标记的其他片面赋予相同的标记，直至所有叶片外缘点云部分χb面片完成分配。
与现有技术相比，本发明具有以下优点：
1)本发明设计了一种新型的3d联合滤波算子对叶片先分离后分割，该滤波算子可以有效地分离以各种角度相互遮挡的叶片。
2)本发明通过面片过分割算法和面片标记生长，抑制了分割噪声，使得分离的叶片中心得以向外生长，得以对充满遮挡的稠密冠层点云进行有效的单叶片分割。
3)本发明适用于双目立体视觉成像系统产生的植物点云，结构光传感器扫描的植物点云以及多视角立体成像系统产生的植物点云，具有对数据来源的适普性。
4)本发明所提出的算法可以为叶面积估计带来方便，对于植物学研究和农业工程是一个非常有效的工具，在孔雀竹芋冠的实验表明所提出方法可以得到平均误差仅为0.47％的叶面积估计结果。
5)本发明所提出的算法可以准确高效地对多类植物和多类点云进行单叶片分割。在4种不同的植物点云上进行了单叶片分割实验，分割的平均点覆盖率达到了96％，平均准确率达到了99.66％，平均f-measure达到了99.66％，分割平均速度达到了12.92秒/株。
附图说明
图1为三种植物点云成像平台，其中，(a)为双目立体视觉系统；(b)为kinect传感器；(c)为多视角立体视觉成像系统。
图2为本发明方法的技术路线。
图3为对四种植物点云采用的预处理方法。第一行为绿萝点云的预处理过程，它使用了滤波器i、ii、iii；第二行为龟背竹点云的预处理过程，它通过滤波器i和iii滤除了大量噪声和背景区域；第三行为孔雀竹芋点云依次通过滤波器i、iv、ii以及v进行预处理；第四行为常春藤点云通过滤波器i、ii以及v进行预处理。
图4为实际植物点云中交叉重叠与共面重叠情况的示例。左半部分是交叉重叠的示例，选取的是包含四片叶片(nos.1～4)的常春藤局部点云，其中，(a)和(c)分别为该点云的俯视图和斜视图，(b)与(d)分别对应(a)与(c)的四片叶片经过网格化渲染后的立体形态展示，以反映交叉重叠叶片相互之间的空间关系，箭头表示为叶片表面的正法向量。右半部分是共面重叠的示例，选取的是包含两片共面的龟背竹叶片(nos.5～6)的局部点云，其中(e)和(g)分别是该点云的俯视图和侧视图，(f)和(h)分别对应(e)与(g)的网格化渲染形态展示。
图5为基于空间半径的离群点滤波器的原理与展示。该图展示了该滤波器作用于两片相互接触的叶片点云的情况。其中xi为当前处理点，在空间半径为r的范围内的邻域点用其他灰度颜色点标记，其数量为ni。邻域点数量ni与nthreshold进行比较后可决定是否在滤波后去除当前点xi。最右侧是该点云经过半径滤波器后的效果，滤除点用更深的灰度颜色进行了标记。
图6为曲面边界滤波器的原理与展示。在(a)中的两片龟背竹叶片点云中我们选取了a点进行边界滤波器的计算，经过(b)中的邻域点pca计算、投影和最大角度差计算后，得到点a最大角度差为为122.04°，大于规定的阈值θthreshold＝90°因此被判定为边界点。在(c)中我们选取了b点进行边界滤波器的计算，经过邻域点pca计算、投影和最大角度差计算后，得到点b最大角度差为为49.32°，小于规定的阈值θthreshold＝90°因此被判定为非边界点。(d)是(a)中点云经过了一次边界滤波器后的结果，其中的某一灰度的点为检测出的曲面边缘点，注意到在右边叶片的内部有个别曲率较大的点也被滤波器认为是边缘点，这个现象是允许的，因为这些曲率突然变大的点也通常是由成像噪声而引起的。
图7为半径滤波器与曲面边界滤波器对不同种类叶片重叠的滤波作用展示。图中左半部分从两个视角展示了半径滤波器对常春藤点云中两片交叉重叠的叶片滤波后的叶片分离效果，右半部分从两个视角展示了曲面边界滤波器对龟背竹点云中共面重叠的叶片滤波后的叶片分离效果。(a)与(b)为交叉重叠叶片的原始点云的两个视角展示；(c)与(d)中将半径滤波器判定的离群点标记为深色点；(e)与(f)分别对应(a)与(b)点云去除深色点的剩余部分。(g)与(h)为共面重叠叶片的原始点云的两个视角展示；(i)与(j)中将曲面边缘滤波器判定的边缘点标记为深色点，集中于叶片的外圈轮廓部分；(k)与(l)分别对应(g)与(h)点云中去除了外圈深色点的剩余部分。
图8为完整单叶片分割算法在常春藤点云上的实际分割效果。由于重叠情况严重，使用了三次3d联合滤波过程，在之后的面片生长中进行了连续三次生长。
图9为对四种植物冠层点云的算法分割定性效果展示。第一行为绿萝点云，第二行为龟背竹点云，第三行为孔雀竹芋点云，第四行为常春藤点云。第一列为植物冠层的原始点云，第二列为人工分割的真实值，第三列为所提出算法对点云各自的分割结果的角度i展示，第四列是所提出算法对点云各自的分割结果的角度ii展示，第五列是所提出算法的分割结果角度iii展示。
图10为点级别分割结果的定量描述。(a)为由实际灰度色彩渲染的常春藤单叶片真实值；(b)为算法对(a)所在的整体点云使用所提出算法得到的该叶片分割结果，其中每一个叶片实例都会使用了一种不同的灰度色进行标记，灰度色的深浅不代表算法标记能力的强弱；(c)为用三种灰度色表示的算法分割定性定量双评价结果，中等深灰度色标记了(a)、(b)完全重合的部分，即代表了算法正确分割的点云部分，最深灰度色区域标记了算法遗漏的叶片部分，最浅灰度色区域标记了算法错误分割的部分。
图11为绿萝点云中所有叶片的点级别分割结果的评价展示。每个叶片对应三个点云，其中a列为由原色彩渲染的单个叶片点云真实值；b列为我们算法得到的分割区域，其中每一个叶片实例都会使用了一种不同的灰度色进行标记，灰度色的深浅不代表算法标记能力的强弱；c列中点云的三种不同灰度色分别对应图10(c)中定义的三种分割评价情况，能够反映本专利方法的性能。
图12为龟背竹点云中所有叶片的点级别分割结果的定性与定量展示。每个叶片对应三个点云，其中a列为由原色彩渲染的单个叶片点云真实值；b列为我们算法得到的分割区域，其中每一个叶片实例都会使用了一种不同的灰度色进行标记，灰度色的深浅不代表算法标记能力的强弱；c列中点云的三种不同灰度色分别对应图10(c)中定义的三种分割评价情况，能够反映本专利方法的性能。
图13为孔雀竹芋点云中所有叶片的点级别分割结果的评价展示。每个叶片对应三个点云，其中a列为由原色彩渲染的单个叶片点云真实值；b列为我们算法得到的分割区域，其中每一个叶片实例都会使用了一种不同的灰度色进行标记，灰度色的深浅不代表算法标记能力的强弱；c列中点云的三种不同灰度色分别对应图10(c)中定义的三种分割评价情况，能够反映本专利方法的性能。
图14为常春藤点云中所有叶片的点级别分割结果的评价展示。每个叶片对应三个点云，其中a列为由原色彩渲染的单个叶片点云真实值；b列为我们算法得到的分割区域，其中每一个叶片实例都会使用了一种不同的灰度色进行标记，灰度色的深浅不代表算法标记能力的强弱；c列中点云的三种不同灰度色分别对应图10(c)中定义的三种分割评价情况，能够反映本专利方法的性能。
图15为本发明的流程图。
具体实施方式
下面结合具体实施例，进一步阐述本发明。应理解，这些实施例仅用于说明本发明而不用于限制本发明的范围。此外应理解，在阅读了本发明讲授的内容之后，本领域技术人员可以对本发明作各种改动或修改，这些等价形式同样落于本申请所附权利要求书所限定的范围。
图1为三维点云成像工具，其中(a)为双目立体视觉成像系统，(b)为微软的kinectv2变结构光传感器，(c)为多视角立体视觉成像系统。双目立体视觉成像系统由两个网络摄像机、一个三脚架以及一台笔记本电脑组成，双目立体视觉成像技术通过利用摄像机从不同角度、位置获得二维图片，然后利用计算机进行立体匹配算法处理，最终获得景物的三维立体信息。传感器采集系统由一个kinect变结构光传感器、一个三脚架以及一台笔记本电脑组成，通过投射的红外线脉冲并测量反射回来的红外模式来获得深度信息，从而实现对植株的三维重建。基于运动重建(structurefrommotion，简称sfm)的多视角立体视觉成像系统由一台智能手机、一个三脚架、一个电动转盘以及一台笔记本组成，该系统通过智能手机的后置摄像机对电动转盘上的作物拍摄不同角度下的多张图像，将图像导入visualsfm软件获得物体的三维点云。
图2为本发明方法的技术路线，所述的植株点云单叶片分割方法适用于双目立体视觉成像系统采集的植株点云，kinect传感器获取的植株点云以及多视角立体视觉成像系统获取的植株点云。现以四个植株点云实例阐述本发明的实施，
步骤s1、对四种植株点云进行预处理。
绿萝点云进行了三段滤波，使用了滤波器i、ii、iii。龟背竹进行了两段段滤波，通过滤波器i和iii滤除了大量噪声和背景区域。孔雀竹芋点云进行了四段滤波，依次通过滤波器i、iv、ii以及v进行预处理。常春藤点云进行三段滤波，通过滤波器i、ii以及v进行预处理。四种植物的滤波结果如图3所示。
滤波器i：基于空间区域的滤波器，可滤除点云在空间坐标系x/y/z轴方向上存在的非叶片信息，比如如果以xoy平面为水平面，z方向的截断可以去除传感器采集到的地面点集。
滤波器ii：基于半径的离群点滤波器，该滤波器的原理是：以查询点xq为球心，r为球体的半径，即当球体内的点数超过阈值nthreshold时，xq在新点云中保留。该滤波器能够滤除点云中的稀疏部分与离群噪声，特别适合于移除传感器插值点以及直径较小的茎干部分。
滤波器iii：基于统计的k-邻近点滤波器，其原理是将输入数据中每个点到k个临近点的距离分布情况进行计算，如果平均距离在标准范围之外，就被定义为离群点而不出现在处理后的点云中。
滤波器iv：颜色滤波器，植物点云中存在部分区域的非叶片信息难以通过空间几何关系进行滤除，但我们可以借助于叶片大多为绿色的先验知识，根据rgb通道中g值的所占的比重区分叶片和非叶片信息。
滤波器v：降采样方法，即创建一个三维体素栅格，然后每个体素内用体素中所有点的重心来近似，从而降低点云复杂度。通过降采样可以使得不同成像方式生成的点云拥有相似的稠密度，从而为算法的参数调谐带来便利。
步骤s2、对步骤s1中获得的植物冠层点云进行3d联合滤波过程，包括以下步骤：
步骤s2-1、采用半径滤波器对点云进行滤波，如图5所示。对于绿萝点云，以点云中任一点xq为球心，r＝0.01m为球体的半径，即当球体内的点数超过阈值nthreshold＝40时，xq在新点云中保留。对于龟背竹点云，以点云中任一点xq为球心，r＝0.01m为球体的半径，即当球体内的点数超过阈值nthreshold＝15时，xq在新点云中保留。对于孔雀竹芋点云，以点云中任一点xq为球心，r＝0.005m为球体的半径，即当球体内的点数超过阈值nthreshold＝13时，xq在新点云中保留。对于常春藤点云，以点云中任一点xq为球心，r＝0.004m为球体的半径，即当球体内的点数超过阈值nthreshold＝8时，xq在新点云中保留。半径滤波器的滤波效果如图7左列所示。
步骤s2-2、采用边界滤波器对点云进行滤波，如图6所示，包括以下步骤：
步骤s2-2-1、找出点云中查询点云中点xi周围最近的k＝20个点，建立点集
步骤s2-2-2、通过pca估计点集
步骤s2-2-3、计算(xj-xi)uv与成分向量u的夹角θj，(xj-xi)uv表示向量xj-xi在u、v组成的法平面上的投影向量；
步骤s2-2-4、对投影向量(xj-xi)uv与成分向量u的夹角集合θ＝{θj|j∈{1,...,k}}进行递增排序，若最大角度差max(θj+1-θj)>θthreshold＝90°，则点xi是一个边缘点，应该在曲面边界滤波中滤除；
步骤s2-2-5、当遍历完点云中所有的点后，滤除所有边缘点就完成了一次边缘滤波。为了达到满意的边缘滤波效果，可以在剩下的点云中重复多次步骤s2-2-1至步骤s2-2-4。边界滤波器的滤波效果如图7右列所示。
步骤s3、按照由内向外的顺序对每次滤除的叶片外缘点云进行面片过分割和面片区域生长，该过程如图8所示，包括以下步骤：
步骤s3-1、对叶片中心区域进行预分割，包括以下步骤：
步骤s3-1-1、从滤波后的点云χc中随机选择一个未标记的点，赋予一个新的叶片标记；
步骤s3-1-2、搜索步骤s3-1-1所选择的点一定半径内是否存在其他未标记的点。如果有则将这些点都赋予这个点相同标记，并将这些新搜索到的点压入一个队列中。对于绿萝、孔雀竹芋点云，搜索半径为0.004m；对于龟背竹、常春藤点云，搜索半径为0.006m。
步骤s3-1-3、每次从步骤s3-1-2得到的队列中取出一个点，继续重复步骤s3-1-1及步骤s3-1-2直到队列为空，则认为对一片叶片的标记已经结束；
步骤s3-1-4、在点云χc中继续随机寻找下一个未标记的点，如果有就赋予其一个新的叶片标记，再重复步骤s3-1-1、步骤s3-1-2和步骤s3-1-3，直到点云中全部点都被标记为止。
步骤s3-2、对滤除的点云部分χb进行过分割，包括以下步骤：
步骤s3-2-1、计算待过分割部分所有点的空间特征，包括平滑度与法向量等；
步骤s3-2-2、基于空间特征生成一定数量的种子点，进行粗略的面片生成；
步骤s3-2-3、使用局部k均值聚类对粗面片区域进行改进，生成精细面片。
步骤s3-3、基于面片标记生长的单叶片分割，包括以下步骤：
步骤s3-3-1、遍历每个完成标记的叶片中心区域，若其中的点与某个小面片邻接，则将该叶片标号赋值给该小面片包含的所有点并渲染相应的标号色彩，直至所有叶片标记区域均不再扩大为止；
步骤s3-3-2、如果滤除的点云部分χb中还有面片未被标记，则该类型面片未与任何已标记叶片区域相邻接，所以该面片代表一个新叶片区域，对其赋予一个独立的叶片标记并渲染新色彩；
步骤s3-3-3、对每个新叶片区域进行基于广度优先搜索的面片标记生长，将中与其邻接但未标记的其他片面赋予相同的标记，直至所有χb面片完成分配。
四种作物点云的单个叶片分割结果如图9所示。
- 还没有人留言评论。精彩留言会获得点赞！