王浩云 閆茹琪 周小莉 馬仕航 胡皓翔 徐煥良

(1.南京農業(yè)大學信息科技學院， 南京 210095； 2.南京農業(yè)大學農業(yè)工程博士后流動站， 南京 210031)

0 引言

作物表型受基因和環(huán)境因素影響，反映植物結構及組成、植物生長發(fā)育過程等全部物理、生理、生化特征。植物果實外形參數很大程度上反映了作物生長環(huán)境、作物產量等，是作物表型的一個重要參數。因傳統(tǒng)的接觸式測量[1-2]存在效率低、不準確等諸多問題，非接觸式測量方法得到飛速發(fā)展。早期，非接觸式測量多基于二維圖像處理，文獻[3]采集葡萄不同時期的果實圖像，利用機器視覺及圖像處理，計算了葡萄果實不同生長時期的表型參數。文獻[4]利用圖像處理的方法，通過閾值分割得到了核桃仁的分割圖像，采用像素統(tǒng)計方法得到核桃仁的面積，實現了核桃仁大小分級。文獻[5]通過分析提取蘋果RGB圖像中單色、波長差等信息，利用形態(tài)學方法處理圖像，提取了目標區(qū)域幾何、灰度和色調頻度等特征，采用多元線性和偏最小二乘回歸預測了蘋果的質量和糖度等指標。文獻[6]運用圖像處理方法提取了黃瓜的瓜長、把長、橫徑差、弓形高度等特征，實現了黃瓜精確分級等。這些方法雖然可以應用于果實分級、檢測及監(jiān)測果實生長態(tài)勢，但單幅二維圖像中所包含三維信息缺失，即使采用機器視覺方法進行處理，其獲取全面的三維信息過程復雜，并伴隨著相應的誤差。

近年來，立體相機、深度相機和激光雷達等三維數據采集設備[7-8]的興起，點云數據的出現為作物三維重建、三維測量等提供了新方法。文獻[9]使用Kinect傳感器獲取角果期油菜分枝圖像，使用基于KD-tree搜索的最近點迭代算法實現點云配準，實現了油菜枝的三維重建。文獻[10]提出一種基于深度信息分割聚類的草莓冠層結構形態(tài)三維重建算法，為草莓采摘機器人果實識別過程中枝葉空間結構關系的構建提供了技術支持。文獻[11]使用三維掃描儀獲取植物葉片的點云數據，并通過點云數據進行預處理，對處理后的點云數據經網格化后得到了高精度的植物葉片模型。文獻[12]使用Kinect相機與旋轉臺獲取蘋果各個角度的點云數據，通過點云配準及貼皮建立了蘋果果實三維模型。文獻[13]利用體著色三維重建的方法獲取變態(tài)根點云數據，并提取變態(tài)根的構型參數進行曲面建模，實現了胡蘿卜和馬鈴薯根系的三維可視化。

基于點云數據的植物三維重建和測量雖然有了長足發(fā)展，但方法多為從不同角度獲取點云信息后進行配準拼接得到完整的作物。這種方法在室內單株植物摘取果實較適用，但在獲取戶外果實的多角度點云信息相對困難，相反獲取單角度的局部點云數據較簡單。為了實現從局部點云信息中直接獲取果實外形指標，本文以蘋果作為研究對象，對基于局部點云的蘋果外形指標估測方法進行研究。采用橢球曲面方程的方法構建蘋果幾何模型，并計算蘋果幾何模型的高度、直徑、體積，使用Kinect V2相機從任意角度獲取點云數據，采用直通濾波法去除點云數據的背景，并用包圍盒算法精簡點云得到蘋果局部點云數據后，采用粒子群算法將蘋果局部點云數據與蘋果模型進行空間匹配，并以遺傳算法求解蘋果最優(yōu)匹配模型的參數，利用蘋果最優(yōu)匹配模型參數估測與其匹配的真實蘋果的外形指標。

1 實驗與方法

1.1 蘋果外形指標測定

實驗選取紅富士蘋果作為研究對象，測量了250個紅富士蘋果的直徑、高度及體積3個參數作為其外形指標。采用游標卡尺測量蘋果直徑、高度。游標卡尺的測量精度為0.01 mm。蘋果體積測量使用排水法，1 000 mL燒杯2個，500、100 mL標準量筒各1個，量筒精度為1 mL。

1.2 點云數據采集與處理

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

開發(fā)軟件使用Matlab 2016b，其工具箱中發(fā)布的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所示。蘋果原始點云包含實驗臺等冗余的背景信息和不穩(wěn)定的噪聲點，不利于后續(xù)數據處理。采用目前最簡單有效的直通濾波法[14]去除背景環(huán)境數據得到蘋果局部點云。由于拍攝設備較為精密，經處理后的蘋果局部點云中數據量仍較大，但是一部分數據點所包含的信息量并不能為后期蘋果局部點云與蘋果模型匹配作貢獻，即處理后的蘋果局部點云中仍包含大量冗余數據點，采用包圍盒算法[15-16]對蘋果局部點云數據進行精簡壓縮。包圍盒算法首先對點云數據建立一個最小的立方體包圍盒，然后將這個包圍盒分割成大小相等的立方體，將點云數據按照三維坐標歸入到一個小立方體中，對每個小立方體找出其中最接近立方體中心的點，并移除其余點從而精簡點云數據。原始點云經直通濾波算法去除背景、又經包圍盒算法精簡后得到只包含蘋果部分的局部點云，如圖3b所示。

1.3 蘋果幾何模型構建及尺寸計算

1.3.1基于橢球曲面方程的蘋果模型

蘋果外形的最小外接曲面近似于一個橢球面，且上下兩端有不同程度的凹陷。采用基于橢球變形的蘋果果實造型方法[17]，通過添加合適的凹凸干擾函數將橢球曲面變形實現蘋果外形的造型。蘋果外層形態(tài)的參數方程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個參數的取值決定了蘋果模型的尺寸，為了使每個模型區(qū)分度更大，確定了蘋果模型的外形參數，分別為蘋果幾何模型的高度(記為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)繼續(xù)改變式(1)中u、w的值，找到蘋果幾何模型在Z=0，X=0時，與Y軸的交點C、D，點C與點D在Y軸方向的差值作為蘋果幾何模型的一個直徑D1；同理找到蘋果幾何模型在Z=0，Y=0時，與X軸的交點，將其在X軸方向的差值作為蘋果幾何模型的另一個直徑D2。

(3)利用積分法對曲面方程式(1)進行積分得到蘋果模型體積作為蘋果模型體積指標V。

1.4 蘋果外形指標估測方法

基于Kinect相機實采的蘋果點云數據是從任意角度獲取的蘋果部分點云，由于信息不完整，從中直接獲取蘋果的全部尺寸參數難度較大。本實驗尋找與真實蘋果尺寸最為接近的蘋果幾何模型并用其尺寸估測真實蘋果外形指標參數。為了尋找與真實蘋果尺寸最為接近的蘋果幾何模型，需要將蘋果局部點云數據與離散的蘋果幾何模型進行空間匹配并尋找匹配效果最優(yōu)的蘋果模型。

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)更新最優(yōu)值。如果位置更新后粒子的距離和E小于原來粒子的距離和，則將更新后粒子位置作為該粒子的局部最優(yōu)解pr。如果該位置更新后粒子的距離和小于原來的全局最小距離和，則將該粒子位置作為本次迭代的全局最優(yōu)解pg，且令兩點集間的最小距離和Emin=E。

(6)根據公式更新粒子速度和位置

(8)

式中c1——自我學習因子，取值為1.3

c2——群體學習因子，取值為1.7

rand()——[0,1]的隨機數

(7)判斷是否到達最大迭代次數，是則執(zhí)行步驟(8)，否則循環(huán)回到步驟(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基于遺傳算法的蘋果幾何模型參數尋優(yōu)

利用粒子群算法實現了蘋果局部點云與離散的蘋果幾何模型的匹配，并可得到蘋果局部點云與該模型的最小距離和。蘋果幾何模型的大小由a、b、c、p1、p25個參數決定，不同的參數組合都會得到不同尺寸的蘋果模型，同一個蘋果局部點云與不同的蘋果幾何模型匹配，效果不同，如圖8所示?？芍?，當匹配到的蘋果幾何模型尺寸過大或過小均不能與蘋果局部點云完全重合，只有匹配到的蘋果幾何模型尺寸與真實蘋果尺寸相近時才可以與蘋果局部點云完全重合。尋找與真實蘋果尺寸最為接近的蘋果幾何模型，即尋找與其局部點云匹配效果最優(yōu)的蘋果模型。本實驗采用遺傳算法[22]求解最優(yōu)蘋果幾何模型參數。

基于遺傳算法的蘋果幾何模型參數尋優(yōu)步驟：

(1)輸入蘋果局部點云M。

(2)隨機初始化種群，種群規(guī)模記為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，是則執(zhí)行步驟(9)，否則返回到步驟(5)。

(9)將最佳個體的染色體解碼得到蘋果幾何模型的5個參數a、b、c、p1、p2，并輸出。

圖11 蘋果外形參數實測結果統(tǒng)計Fig.11 Apple volume measurement results statistics

根據上述步驟得到與真實蘋果尺寸最為接近的蘋果幾何模型的5個參數a、b、c、p1、p2，根據參數建立相應的蘋果幾何模型，然后計算蘋果幾何模型的外形尺寸，即該蘋果局部點云對應的真實蘋果的外形參數估測值?；诰植奎c云的蘋果外形指標估測方法流程圖如圖10所示。

圖10 基于局部點云的蘋果外形指標估測方法流程圖Fig.10 Flow chart of apple shape index estimation method based on local point cloud

2 結果與分析

2.1 蘋果尺寸實際測量結果

實驗選取了250個紅富士蘋果作為實驗樣本，通過排水法測量蘋果的體積，利用游標卡尺測量蘋果的高度和2個相互垂直的直徑。圖11為蘋果外形參數實測結果統(tǒng)計圖。利用數學統(tǒng)計學方法對測量結果進行分析，結果如表1所示。由表1和圖11可知，蘋果的體積分布范圍為200～700 mL，蘋果的高度分布范圍為65～100 mm，直徑分布范圍為70～110 mm，實驗所用蘋果外形參數分布范圍廣，說明實驗樣本大小不一，類型較豐富，避免了樣本單一性，使實驗更加具有說服性、普遍性。

2.2 基于局部點云的蘋果外形指標估測結果

分別獲取了250個蘋果在3個角度下(圖2)的點云數據，共計750幅點云圖像。將處理后的250個蘋果3個角度下的點云數據作為基于局部點云的蘋果外形指標估測方法的輸入，分別估測其對應的蘋果果實外形指標(蘋果的體積V、高度h、直徑D1和D2)，并將估測值與蘋果果實實測值進行線性回歸分析，如圖12～14所示。圖中橫軸真實值為實驗人工測量的蘋果外形參數；縱軸估測值為使用基于局部點云的蘋果外形指標估測方法估測的蘋果外形指標；R2表示回歸直線對觀測值的擬合程度，最大值為1，其值越接近1，說明回歸直線對觀測值的擬合程度越好；RMSE表示均方根誤差(亦稱標準誤差)，用來衡量觀察值同真實值之間的偏差。表2為蘋果外形指標估測結果。同時，利用數學統(tǒng)計法計算了各個蘋果外形指標的平均估測誤差，如表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下蘋果外形指標估測結果最好，而在田間作業(yè)時，機器人獲取角度2的果實數據最為方便，說明該方法可以為田間作業(yè)中獲取生長期的果實外形指標提供理論基礎。

(6)由表3可知，蘋果體積估測結果的平均誤差不大于16.16 mL，高度估測結果的平均誤差不大于2.92 mm，直徑估測結果的平均誤差不大于2.35 mm，估測結果的平均誤差較小，在允許誤差范圍內，按照國內現行的蘋果分級標準GB/T 10651—1989《鮮蘋果》標準可知蘋果直徑是蘋果分級的重要指標，該方法可以用于蘋果分級中，有實用價值。實驗結果表明，本文方法預測結果相對穩(wěn)定，結果符合實際，具有實用性。

2.3 基于彩色圖像的蘋果外形指標估測結果

Kinect相機為RGB-D相機，可以同時采集彩色圖像、深度圖像和點云數據。本實驗在采集250個蘋果局部點云數據的同時采集了其在3個角度下的彩色圖像。在角度1、角度3下的彩色圖像只包含蘋果果實直徑，無法得到高度信息，且由2.2節(jié)分析可知在角度2下基于局部點云的蘋果外形指標估測效果最好。因此，基于彩色圖像的蘋果外形指標估測使用角度2下的彩色圖像進行對比實驗。實驗通過數字圖像處理對彩色圖像進行處理，得到蘋果圖像分割圖，并利用最小外接矩形檢測蘋果。通過外接矩形的4個頂點坐標計算出外接矩陣的長和寬，該長和寬是蘋果果實的高和直徑的像素點數，根據張正友標定方法得到每個像素點代表的實際距離為0.56 mm。將每個像素數乘以像素點大小得到蘋果的高度和最大直徑。利用SVR支持向量機回歸模型擬合體積參數，其中200個數據作為訓練集，50個數據作為測試集。根據數學統(tǒng)計學方法計算結果的平均誤差。將角度2下2種蘋果外形指標估測結果進行對比，結果如表4所示。

表4 估測平均誤差對比Tab.4 Estimated average error comparison

由表4可知，基于彩色圖像的蘋果外形參數估測結果平均誤差均大于基于局部點云的蘋果外形指標估測結果。說明，基于局部點云的蘋果外形指標估測更加準確，符合實際。

3 結論

(1)采用直通濾波法去除原始點云的背景信息，得到只包含蘋果部分的點云數據，然后使用包圍盒算法精簡點云，得到了準確的蘋果局部點云。

(2)采用遺傳算法尋找蘋果最優(yōu)匹配模型并求解模型參數，其中匹配方法采用基于粒子群算法的蘋果點云空間匹配，實現了蘋果局部點云的位置變換，使空間中處于不同位置的蘋果局部點云與蘋果幾何模型重合。

(3)本文算法估測的蘋果外形指標結果與真實值具有較高的相關性，其線性回歸分析的R2均高于0.7，且蘋果體積估測的平均誤差不大于16.16 mL，高度估測的平均誤差不大于2.92 mm，直徑估測的平均誤差不大于2.35 mm，估測結果的平均誤差較小，在允許誤差范圍內，說明該方法具有可行性、實用性。

(4)本文方法為機器人田間獲取生長期果實的外形指標提供了方法，且在角度2下的蘋果外形指標估測效果最好。