耿泽栋, 杨万能, 李峰, 戴雄泽, 欧立军, 邹学校, 冯慧
(1.华中农业大学植物科学技术学院, 武汉 430070; 2.华中农业大学园林学院, 武汉 430070; 3.湖南农业科学院, 长沙 410125; 4.湖南农业大学园艺园林学院, 长沙 410128; 5.华中农业大学信息学院, 武汉 430070)
我国是辣椒生产、消费和出口大国[1]。辣椒是我国人民生活中不可或缺的蔬菜之一,其经济价值越来越受到人们广泛的重视。“十二五”时期,我国十分重视辣椒产业,国家“863”计划项目、科技支撑计划项目等都对辣椒的研究和发展提供了强有力的支持[2]。当前,辣椒产业发展迅速,在遗传育种方面都取得了突破性的成绩,获得了高产优质的优良品种[3-4]。对辣椒农艺性状的分析有利于加深对辣椒性状特征的了解,为辣椒的选育和功能基因组分析提供基础。目前,我国对辣椒农艺性状的采集主要是通过人工方式,分析的速度较慢,严重浪费人力物力。
数字图像处理技术是对图像进行深度加工,进行有目的的改变,强调原始图像与加工图像之间的变换。图像处理技术具有成像清晰、适用性广、处理精度高等特点[5]。随着计算机技术的提高和成本的下降,图像处理技术在农业方面的应用越来越广泛,越来越向专业化、系统化、智能化的方向发展,其可以为作物农艺性状的采集和分析提供强有力的工具[6]。应用数字图像处理技术对辣椒农艺性状进行采集分析可以提高采集效率,节省人力资源,促进辣椒的研究。
大量研究表明,辣椒的表型特征可以反映辣椒性状的遗传规律,可以通过研究辣椒的农艺性状制定出相应的筛选机制,从而获得高品质的辣椒品种。曲晓斌等[7]对线辣椒农艺性状的研究表明,果横径的增加会导致单株结果数呈下降趋势,果纵径与单株结果数之间的相关性不大。果横径的增加还会导致单果果肉重,平均单果重呈明显上升趋势[8-9]。文锦芬等[10]和张建云[11]对辣椒径尖的研究表明,辣椒顶部表面附生着多种病原菌,严重时将导致辣椒病害。王得元[12]研究表明,辣椒果形指数(果实长宽比)受基因控制,稳定遗传。因此,可以通过挑选厚度较大、果形指数高的辣椒来选择高品质的辣椒的品种。
与人工提取相比较,基于图像处理可以更加快速地提取辣椒农艺性状。杨一璐等[13]基于叶绿素荧光图像对辣椒叶片氮含量进行预测,建立了辣椒叶片荧光特征参数对氮含量的预测模型,方便对辣椒叶片含氮量进行科学合理的预测。袁开放等[14]通过图像处理技术设计了牛角椒的大小分选系统,取代传统的人工作业,实现了机器对牛角椒进行大小筛选,机器选择正确率在95%以上。王建玺等[15]利用图像处理技术,改进了主动轮廓模型(c-v模型),提高了对辣椒病斑图像分割的效果。Chupawa等[16]利用图像处理技术对辣椒的种子进行快速质量检测。综上所述,图像处理技术可以为辣椒性状的研究提供强有力的工具,降低人工成本,提高采集效率。
目前辣椒直径、果型指数、心室个数等性状提取的方式主要是人工提取,而一些利用图像处理技术对辣椒大小、病斑等性状的研究缺乏对辣椒内部性状,如圆形度、果肉厚度、横向褶皱等综合的检测和提取。本研究采用机器视觉技术,详细分析多种辣椒内部性状特征,包含辣椒横截面的果长、果宽、果肉厚度、果形指数(果长/果宽)、心室数、横向褶皱及辣椒纵截面的中线、肩宽、胸宽、腹宽、尖宽、中线曲率等性状,对辣椒农艺性状的提取更加全面、快速。
本文通过图像分析技术,基于机器视觉技术编写了一套完整的对辣椒横切和纵切截面的农艺性状的采集和分析程序。对辣椒横切和纵切截面图进行处理分析,可以快速、准确得出辣椒的果肩径、直径、心室个数、果肉厚度等性状,实现了辣椒农艺性状高通量、自动化提取。
本文使用的编程软件:Labview2015(64-bit)[17]、Microsoft Visual Studio2015(VS2015)[18]、OpenCV3.1.0[19]。
根据辣椒的横切、纵切截面图特点和鉴别分类要求,建立了辣椒性状的自动化检测流程:背景板上放好比色板和直尺,摆放辣椒切片进行图像采集,通过Labview编写的程序对辣椒横纵截面图像进行分割、识别,进而提取辣椒的内部特征;在提取过程通过对各像素点的阈值分析,将各辣椒截面的图像数据转变为简单、易操作的像素集合;对截面进行目标识别和种类划分,提取出二值图、单个心室图、轮廓图等图像,分析图像得到辣椒截面直径、圆形度、辣椒壁厚度、心室数、横向褶皱等性状数据。计算并分析所有截面的性状数据后,将其自动保存至TXT文档,供后续分析研究。对比于人工检测,基于机器视觉的图像处理方法更加快速、简便,可实现对辣椒截面的自动化、高通量性状检测。
辣椒横截面的提取性状如图1所示,首先做出截面的最小外接矩形长Dx和宽Dy,再根据矩形的对角线确定辣椒横截面的质心(O点),黑色区域H1、H2为辣椒心室,红线部分为提取辣椒厚度的部分(心室连接处不计算横截面厚度),绿色虚线部分为辣椒的横向褶皱。
注:O—截面质心;Dx—最小外接矩形长;Dy—最小外接矩形宽;H1和H2—心室。红线区域为提取辣椒厚度,绿色区域为横向褶皱。Note: O—Facet center; Dx—Minimum circumscribed rectangle length, Dy:Minimum circumscribed rectangle width; H1 and H2:Ventricle. Red domain is pepper thickness, and green domain is lateral pleate.图1 横截面提取指标Fig.1 Crosscutting section map extraction index
1.3.1横截面图像分割 本研究采用阈值分割法[20]将灰度图像转换成二值图像。它的作用是减少图像中数据量,凸显出目标的轮廓[21]。由于背景较为固定,辣椒则有红色、绿色、黄色等颜色,故本文中采用OTSU算法进行分割,如图2所示。二值化处理的优劣程度对性状提取的准确性和合理性具有十分重要的作用[22]。
1.3.2直径比Dx/Dy分析 直径是反映圆特征的最直观数据。对于不规则的圆,分别取其最小外接矩形的长和宽。对于辣椒来说,直径的大小可以直接反映其粗细程度。如图2D所示,首先对横截面进行填充,然后取填充图最小外接矩形的长和宽。对横截面直径比Dx/Dy,大于1表示细长的截面,等于1表示圆形的截面[23]。
1.3.3圆形度DC分析 圆形度DC是不规则圆形相对于规则圆的偏离情况,表现了辣椒截面的不规则程度。本文通过Labview算法模块IMAQ Particle Analysis VI中的Area提取填充图的面积记为A;通过Perimeter提取填充图的周长记为L。圆形度计算公式(1)如下。
DC=4πA/L2
(1)
其中,DC为圆形度,A为外轮廓所包含面积,L为外轮廓周长。
1.3.4果肉厚度T分析 对于果肉横截面,内轮廓如果存在室间隔,则无法计算果肉厚度,需要去掉此部分,具体计算分为以下两个过程。首先提取出横截面的内轮廓(图2E、F),计算内轮廓上每一点与内轮廓凸包上每一点之间的距离,以最小值作为内轮廓上该点到凸包的距离。然后对此距离进行判定,给定一个阈值(此处为10个像素),只有距离小于阈值的部分才可计算厚度。
具体过程如图3所示,质心以及整幅图像的四个角点将图像区域分为四个区域,并分别计算四个区域的面积,其中图像高记为height,宽记为width。之所以进行区分是因为质心点和内轮廓上的点与图像边缘的交点不在同一个区域,四个区域的计算方法不同。首先根据面积判断内轮廓上的点在哪一块区域内,以内轮廓上的点为中心,分别与相邻的角点以及质心相连,应用海伦公式计算S1、S2、S3三个区域的面积之和,该面积与四个区域中某一区域的面积相等,则内轮廓上的点在该区域上。如在区域1内,计算该直线与y=0的交点,如在区域2内,计算该直线与x=width-1的交点,如在区域3内,计算该直线与y=height-1的交点,如果在区域4内计算该直线与x=0的交点。计算出交点后,统计从质心点到该交点的路径间,内轮廓与外轮廓的交点之间的距离,则为该处的厚度。
A:原始图;B:二值化优化;C:单个辣椒二值图;D:填充图;E:内轮廓及凸包;F:厚度轮廓;G:心室图;H:褶皱图。A: Original image; B: Binary optimization; C: Single pepper binary image; D: Fill image; E: Inner contour and conver hull; F: Thickness image; G: Venticular image; H: Pleated graphics.图2 辣椒横截面图像处理Fig.2 Image processing of pepper slice
1.3.5心室个数H取面积最大两块区域均值的50%,记为心室个数的阈值A,公式(2)如下。
A=(a+b)/4
(2)
式中,A为心室个数阈值,a为最大区域面积,b为第二大区域面积。
若剩下的区域面积超过这个阈值,则认为是心室,小于这个阈值则去除。本文通过Labview软件中的IMAQ Particle Analysis算法模块提取心室图中的粒子数,记为心室个数。对二值化心室图像(图2G)进行分析,计算出联通区域的个数及其面积、图像在X轴方向的中心坐标、在Y轴方向的中心坐标等[24]。
1.3.6横向褶皱VWN和VWD横向褶皱是指辣椒横截面外轮廓产生的弯曲,如图2H所示。计算横向褶皱,用外轮廓凸包减去外轮廓填充图,得到褶皱图形,分别计算褶皱外接矩形的长和宽,如果其中一项超过厚度的1/2,则认为有2个交点,如果其中一项与厚度的1/2相等,则认为有1个交点。用外轮廓凸包(图2H外圈黄线)和内轮廓凸包(图2H内圈黄线)将截面轮廓(图2H中的红线)夹起来,做出内外侧轮廓凸包的中线(图2H蓝色实线),蓝线长度(L)为内轮廓凸包周长和外轮廓凸包周长和的1/2,测量黄线间区域面积S为外轮廓填充面积减去内轮廓凸包的填充面积。横向褶皱的计算公式如下。
图3 果肉厚度示意图Fig.3 Skech of pulp thickness
VWN=N/L
(3)
VWD=S/L
(4)
式中,VWN为平均褶皱数量,VWD为平均褶皱面积,N为交点个数,L为内外轮廓凸包周长和的一半(蓝线长度),S为两黄线间的面积。
如图4所示,首先提取辣椒区域图像,然后分割出单个辣椒图像,提取辣椒骨架,并确定辣椒的两个端点(红色D1、D2),即绿色萼片果肩界线(绿色虚线)的中点D1及果尖D2。不同种类辣椒端点的位置不同,如在线椒型辣椒中D1与D2大多是突出的,而在甜椒型辣椒中大多是凹陷的。对辣椒截面图进行分割处理,提取出辣椒纵截面的性状特征。
通过以下几个数据可以得到辣椒的高度、宽度和弯曲程度等基本性状信息,基本还原出辣椒纵截面的形状。采集的性状特征如图4所示:果肩界线N(绿色虚线)为辣椒果柄处绿色萼片轮廓;果柄中线M(白色虚线)为果肩界线端点D1处的切线的垂线;中线长Y从D1到D2外围轮廓的中线(蓝色实线)长,随辣椒的弯曲而弯曲;肩宽X1为端点D1切线沿果柄中线向下平移5 mm,其与轮廓线的交点间的距离;胸宽X2为辣椒果实偏上1/3处的宽度;腹宽X3为辣椒果实偏下1/3处的宽度;尖宽X4为端点D2沿中线上移5 mm处中线切线的垂线与轮廓线交点间的距离。
本文程序对辣椒纵截面进行特征提取时,首先对整体图片进行二值化,得到辣椒区域,然后从图像左上方到右下方进行扫描得到单个辣椒纵截面的外接矩形,并对该外接矩形四周各加40像素值进行扩大剪切,提取出单个辣椒的二值化图像。由于两个辣椒截面摆放的距离较近而导致在剪切单个辣椒二值图时会包含相邻的辣椒部分,需要去除与边界相连的区域,得到单个辣椒的二值图。通过Labview模块扫描得到单个辣椒纵截面的外接矩形,消除接触图像边界的粒子;调用VS编写的dll(采用不断细化的方式,最后只剩下单个像素即为骨架)对辣椒纵截面二值图提取骨架。果柄部分是绿色的,所以采用超G分割提取辣椒柄部,计算方法如公式(5)所示,阈值设为20,N≥20则认为此处是果柄。然后利用二值化图像减去柄部图像,得到去除果柄的辣椒图像。
注:N—果肩界线;D1—顶部端点;D2—尖部端点;M—果柄中线;Y—中线;X1—肩宽;X2—胸宽;X3—腹宽;X4—尖宽;R—中线曲率。Note: N—Fruit shoulder line; D1—Top end point; D2—Tip end point; M—Fruit handle midline line; Y—Centerline; X1—Shoulder width; X2—Chest width; X3—Belly width; X4—Tip width; R—Midline curvature.图4 纵截面提取指标Fig.4 Longitudinal section map extraction index
N=2G-R-B
(5)
式中,G为绿色分量,R为红色分量,B为蓝色分量。
本文通过中线形状反映中线曲率。中线曲率即中线上每点(1/20Y间隔)切线与果柄中线M的夹角度数,如图4中紫色虚线与白色虚线的夹角R是针对中线上点的切线方向对弧长的转动率,可以用来衡量辣椒中线的不平坦程度,是对辣椒纵截面弯曲程度更详细的表示,曲率越大,则表示中线的弯曲程度越大,辣椒的弯曲程度越大。
本文提取辣椒中线长上每点的坐标并计算中线上所有像素点的个数,用像素个数除以20的值作为检索值对中线坐标进行检索,将中线进行20等分,可以找到中线上每一个20等分点的坐标。用图像的宽度减去每一个20等分点的横坐标可以得到共20 个中线到图像右边界的距离值,记为宽度W。对这些值做波形分析,可以还原辣椒的中线形状,它的弯曲程度可以反映中线曲率的变化趋势。
如图5所示,当DC值越接近1(0.994、0.988)时,辣椒横截面越接近圆形;DC值越小(0.934),表明图形的周长越大,面积越小,与圆形的差距越大,辣椒横截面越偏离圆形[25]。
如图6所示,统计从质心点到内外壁上交点的路径距离可以得到一系列内轮廓与外轮廓的交点之间的距离,即为果肉某一点的厚度。将所有得到的数据进行排列,得到果肉厚度变化曲线,可以形象的反映果肉厚度的变化过程。
心室内容较少可直接准确计算心室,心室内容较多无法准确计算。如图7所示,本文对单个二值图中心进行画圆处理,如果辣椒横截面中间果肉或籽粒较多,需去除中间的辣椒籽的部分,只提取出心室的模块图;如果辣椒横截面的心室较多,需对提取的心室图进行筛选。经计算分析,两个样品1心室个数为2,两个样品2心室个数为3。
本程序对辣椒纵截面进行分析可以得到一系列辣椒二值化图像(图8),用于性状提取。
对纵截面得到的中线大小做归一化处理,将20个宽度值放在同一坐标系下,将长度大小不同的辣椒中线变换为统一标准的折线图。消除图片大小和辣椒大小对还原中线形状的影响,方便不同规格辣椒之间的中线形状的比较。将不同辣椒的中线归一化(图9),可以明显观察到不同辣椒的弯曲程度,即中线形状。
图5 圆形度Fig.5 Circularity
图6 果肉厚度Fig.6 Pulp thickness
图7 辣椒心室提取Fig.7 Processing of pepper ventricular
2.6.1横截面特征性状参数分析 本研究通过上述方法检测了40个横截图像(800个辣椒横截面),58个纵截图像(422个辣椒纵截面),得到大量辣椒的截面表型数据,如表1所示。其中心室个数和VWN的CV值较高,其他性状的CV值较低,表明该品种的辣椒之间差异度较小。
图9 中线形状分析Fig.9 Analysis of centerline shape
2.6.2纵截面特征性状参数分析 提取到辣椒纵截面中线长Y、肩宽X1、胸宽X2、腹宽X3、尖宽X4等,如表2所示,其中肩宽CV值最大,尖宽CV值较大,表明该品种辣椒的肩部区别较大,中线长度、肩宽、腹宽基本一致。
图片分辨率是度量图像内像素个数多少的指标[26]。不同分辨率图像中同一区域所包含的像素个数不同,图片的分辨率更高,则表示该图像的内存占据更大,包含的数据更大,相同区域中包含的像素个数更多,能展示更多的细节。图片的实际大小可以通过分辨率和像素来计算,如公式(6)所示。
L=2.54A/DPI
(6)
表1 辣椒横截面图像的特征性状参数Table 1 Characteristics of pepper transect images
表2 辣椒纵截面图像的特征性状参数Table 2 Characteristic parameters of vertical cutting image of pepper
式中,L为图片长度(cm),A为像素个数,DPI为每英寸上的像素个数。
如图2中通过5 mm对应的像素个数的多少来判断肩宽、尖宽的位置。在不同分辨率的图像中,5 mm所代表的像素个数不同,对于96 DPI,5 mm像素个数为19,对于600 DPI,5 mm像素个数为118。所以在不同分辨率的图像中通过取不同的像素个数来表示同一长度。
图像处理技术在农业中的应用十分广泛,其本质就是使用不同的图像处理平台对农作物进行处理。孙宏佳[27]基于机器视觉利用Labview设计了花生种子自动识别系统,主要使用了IMAQ Vision开发环境,可以高效、准确地提取花生种子的性状特征。Labview程序具有强大的外部程序接口能力,本文通过Labview调用多个DLL,丰富了程序的运算功能。李莲等[28]利用基于卷积神经网络的深度学习识别方法,识别辣椒的颜色形状,正确率有86.67%。但是该方法仅识别辣椒的特色和形状,不能获得更详细的辣椒性状特征。与之相比,本文对辣椒的多种性状进行提取,可以获取更完善的辣椒性状特征。
对辣椒形态学性状的研究表明,辣椒性状可以反映其种植遗传规律[7-12]。本文提供了一种可以高通量、自动化的辣椒表型性状提取方法,提取了辣椒横截面的果长、果宽、果肉厚度、果形指数(果长/果宽),还增加了心室数、横向褶皱及辣椒纵截面的中线、肩宽、胸宽、腹宽、尖宽、中线曲率等性状。相比人工采集数据,图像处理技术更加快捷方便,效率大大提高。与传统的人工测量辣椒的农艺性状相比,本文采用数字图像处理技术采集辣椒横纵截面的性状特征。人工采集性状因为采集者的主观判断和测量方式不同,无法保证每次采集的标准一致,可能会存在误差,导致采集的性状出现偏差。而本文采用计算机图像处理技术,采集数据有统一的标准,消除了主观误差,确保性状采集的准确性和一致性。且该程序具有统一的性状采集标准,不会受到人为因素的干扰。简单高效,便于数据的储存和读取,可广泛地应用于辣椒内部性状的提取。