朱超 吳凡 劉長斌 趙健翔 林麗麗 田雪瑩 苗騰

摘要： 針對當(dāng)前三維點(diǎn)云處理方法在玉米植株點(diǎn)云中識別雄穗相對困難的問題，提出一種基于超體素聚類和局部特征的玉米植株點(diǎn)云雄穗分割方法。首先通過邊連接操作建立玉米植株點(diǎn)云無向圖，利用法向量差異計算邊權(quán)值，并采用譜聚類方法將植株點(diǎn)云分解為多個超體素子區(qū)域；隨后結(jié)合主成分分析方法和點(diǎn)云直線特征提取植株頂部的子區(qū)域；最后利用玉米植株點(diǎn)云的平面局部特征在頂部子區(qū)域中識別雄穗點(diǎn)云。對3種點(diǎn)云密度的15株成熟期玉米植株點(diǎn)云進(jìn)行測試，采用F1分?jǐn)?shù)作為分割精度判別指標(biāo)，試驗(yàn)結(jié)果與手動分割真值相比，當(dāng)點(diǎn)云密度為0.8、1.3和1.9個點(diǎn)/cm時，雄穗點(diǎn)云分割的平均F1分?jǐn)?shù)分別為0.763、0.875和0.889，分割精度隨點(diǎn)云密度增加而增高。結(jié)果表明，本研究提出的基于超體素聚類和局部特征的玉米植株點(diǎn)云雄穗分割方法具備在玉米植株點(diǎn)云中提取雄穗的能力，可為玉米高通量表型檢測、玉米三維重建等研究和應(yīng)用提供技術(shù)支持。

關(guān)鍵詞： 玉米雄穗；三維點(diǎn)云分割；表型檢測；超體素聚類；局部特征；主成分分析

中圖分類號： S24???????????? 文獻(xiàn)標(biāo)志碼： A??????????? 文章編號： 202102-SA001

引用格式：朱超, 吳凡, 劉長斌, 趙健翔, 林麗麗, 田雪瑩, 苗騰. 基于超體素聚類和局部特征的玉米植株點(diǎn)云雄穗分割[J]. 智慧農(nóng)業(yè)(中英文), 2021, 3 (1): 75-85.

Citation:ZHU Chao, WU Fan, LIU Changbin, ZHAO Jianxiang, LIN Lili, TIAN Xueying, MIAO Teng. Tassel segmentation of maize point cloud based on super voxels clustering and local features[J]. Smart Agriculture, 2021, 3 (1): 75-85. (

1? 引? 言

玉米是重要糧食作物，其產(chǎn)量對保障全球糧食供應(yīng)至關(guān)重要［1］。玉米雄穗相關(guān)的表型參數(shù)是玉米育種的重要農(nóng)藝性狀［2-4］。傳統(tǒng)的雄穗表型參數(shù)提取方法主要是人工測量，方法效率低、主觀性強(qiáng)、破壞性強(qiáng)。近年來，計算機(jī)視覺技術(shù)的發(fā)展為作物表型參數(shù)的自動測量提供了條件［5,6］。

基于二維圖像的視覺技術(shù)率先被應(yīng)用于雄穗表型提取中。Gage等［7］開發(fā)了一種基于圖像的玉米雄穗表型分析系統(tǒng)，可對離體的雄穗進(jìn)行標(biāo)準(zhǔn)圖像獲取，進(jìn)而提取穗長、穗寬等表型參數(shù)。一些研究工作聚焦于在玉米植株或冠層圖像中識別雄穗。Ye等［8］結(jié)合方向梯度直方圖特征和支持向量機(jī)分類器對雄穗?yún)^(qū)域進(jìn)行識別，再采用顯著性檢測方法實(shí)現(xiàn)雄穗像素分割。Kurtulmus和Kavdr［9］結(jié)合顏色特征和支持向量機(jī)分類器進(jìn)行雄穗定位，并采用層次聚類方法實(shí)現(xiàn)雄穗分割。Lu等［10］采用多個角度的圖像提升雄穗分割的精度，并實(shí)現(xiàn)了對雄穗顏色、分支數(shù)等參數(shù)的粗略測量。Lu等［11］構(gòu)建了雄穗識別深度學(xué)習(xí)網(wǎng)絡(luò)TasselNet，實(shí)現(xiàn)了在玉米冠層圖像中對雄穗的準(zhǔn)確計數(shù)。二維圖像只是三維空間的一個投影，利用二維圖像可以較好地進(jìn)行雄穗計數(shù)工作，但對穗長、穗寬等參數(shù)的提取結(jié)果并不準(zhǔn)確。

隨著三維激光掃描儀［12］、飛行時間相機(jī)［13］，以及激光雷達(dá)［14］等三維傳感器的普及，基于三維點(diǎn)云的視覺技術(shù)被用于作物表型研究，如作物種類識別［15］、作物長勢動態(tài)監(jiān)控［16］、冠層三維重建［17,18］，以及各類表型參數(shù)的測量［19-22］等，這些研究證明了從三維空間提取作物表型具有更高的準(zhǔn)確性［23,24］。但是，當(dāng)前基于三維點(diǎn)云的玉米雄穗表型研究仍然較少。韓東等［25］利用立體視覺技術(shù)對離體的雄穗進(jìn)行三維點(diǎn)云重建，同時采用基于密度的聚類算法對雄穗進(jìn)行分割，進(jìn)而提取表型參數(shù)。

目前，在玉米植株或群體點(diǎn)云中對雄穗表型參數(shù)進(jìn)行提取仍具挑戰(zhàn)，其中，對雄穗點(diǎn)云的準(zhǔn)確分割是技術(shù)難點(diǎn)。針對該問題，本研究提出一種基于超體素聚類和局部特征的玉米植株點(diǎn)云雄穗分割方法，可直接在玉米植株點(diǎn)云上對雄穗點(diǎn)云進(jìn)行分割和識別，為玉米高通量表型檢測、三維幾何重建等研究和應(yīng)用提供技術(shù)基礎(chǔ)。

2? 材料與方法

2.1試驗(yàn)材料與點(diǎn)云數(shù)據(jù)獲取

田間試驗(yàn)于2019年5月至10月在沈陽農(nóng)業(yè)大學(xué)玉米試驗(yàn)田（116.16°E，39.56°N）進(jìn)行。選取種植行距60 cm、株距25 cm的成熟期玉米植株為研究對象，使用FreeScan X3手持激光掃描儀進(jìn)行玉米植株三維點(diǎn)云獲取。將田間成熟期玉米植株移植到盆中，并轉(zhuǎn)移到室內(nèi)進(jìn)行點(diǎn)云獲取。由于FreeScan X3激光掃描儀需要借助激光感光片進(jìn)行三維定位，因此，在采集數(shù)據(jù)時，先將感光片粘貼在輔助支架上，并將支架圍在被測植株周圍，保證激光線能同時掃描到植株和感光片，數(shù)據(jù)采集過程如圖1（a）所示。利用CloudCompare軟件對獲取的原始點(diǎn)云數(shù)據(jù)進(jìn)行預(yù)處理，手動刪除盆、移動支架的點(diǎn)云數(shù)據(jù)，僅保留玉米植株點(diǎn)云數(shù)據(jù)。

為驗(yàn)證本研究分割算法對不同密度點(diǎn)云植株的分割效果，選取15株不同形態(tài)的成熟期玉米植株點(diǎn)云作為測試樣本，每個玉米植株點(diǎn)云均被采樣到3個點(diǎn)云密度（如圖1（b）），分別為0.8、1.3和1.9個/cm，點(diǎn)云密度表示每厘米存在多少個點(diǎn)；3個密度的植株點(diǎn)云平均數(shù)量約為6000、15,000和30,000個點(diǎn)，雄穗點(diǎn)云平均數(shù)量約為140、360和890個點(diǎn)。利用CloudCompare軟件對被測玉米植株點(diǎn)云進(jìn)行手動分割獲得雄穗點(diǎn)云的分割真值。

2.2雄穗點(diǎn)云分割方法

本研究方法以玉米植株散亂點(diǎn)云為輸入，通過點(diǎn)云超體素聚類、植株頂部子區(qū)域提取和雄穗點(diǎn)云識別3個主要步驟，實(shí)現(xiàn)對雄穗點(diǎn)云的識別和分割。

2.2.1點(diǎn)云超體素聚類

本研究首先通過點(diǎn)云超體素聚類，將植株點(diǎn)云分解為多個點(diǎn)云子區(qū)域，這些子區(qū)域由具有幾何相似性的局部區(qū)域點(diǎn)云整合形成，并且這些子區(qū)域較好地保留了植株點(diǎn)云器官間的邊界信息，便于后續(xù)進(jìn)行雄穗點(diǎn)云的分割。本研究采用譜聚類算法實(shí)現(xiàn)點(diǎn)云超體素聚類，假設(shè)輸入的三維玉米點(diǎn)云為P =｛pi|i=1，2，…，n｝，點(diǎn)數(shù)量為n，對P進(jìn)行超體素聚類的具體步驟如下。

（1）邊連接操作。對P中的所有點(diǎn)均進(jìn)行邊連接操作形成無向圖。對于任意一點(diǎn)pi，首先通過K近鄰空間搜索方法［26］獲得P中與點(diǎn)pi歐式距離最小的k個近鄰點(diǎn)，形成點(diǎn)集Ai =｛aj| j=1，2，…，k，aj ∈ P｝。之后遍歷Ai中的所有點(diǎn)，判斷每個點(diǎn)是否能與pi連接形成邊。如果Ai中的點(diǎn)aj在P中的k個近鄰點(diǎn)不包含pi，則pi與aj不能連接成邊；反之則可連接形成1條邊，且該邊的權(quán)重值e通過公式（1）進(jìn)行計算：

θ_ij=acos ???????????????????????? (n_i?n_j )

e={(0.1@0.9)┤ ((@)&(θ_ij<0.9@θ_ij≥0.9))??? （1）其中，ni、nj分別表示點(diǎn)pi、aj的單位法向量；θij為ni和nj之間的夾角，rad。本研究利用點(diǎn)云處理庫（Point Cloud Library，PCL）［27］中的法向量估算模塊計算點(diǎn)云法向量。該模塊首先采用主成分分析法［28］估算法向量，初始的法向量存在定向不一致的問題，導(dǎo)致相鄰點(diǎn)云的法向量方向可能相反；該模塊通過視點(diǎn)位置對法向量方向進(jìn)行調(diào)整，給定一個視點(diǎn)三維坐標(biāo)，調(diào)整所有點(diǎn)云法向量的朝向，使其均指向視點(diǎn)方向。

（2）構(gòu)建歸一化拉普拉斯矩陣。構(gòu)建植株點(diǎn)云P的鄰接矩陣W，該矩陣為n×n大小的方陣，其第i行、k列的矩陣元素wik表示點(diǎn)pi與pk之間的連接關(guān)系，如果pi與pk之間無邊連接，則wik=0；如果pi與pk之間有邊連接，則wik等于該邊的權(quán)重值。構(gòu)建n×n大小的度矩陣D，該矩陣第i行、k列的矩陣元素值dik通過公式（2）計算。

d_ik={(∑_(k=1)^n?〖w_ik (&i=k)〗@0(& ????????????????????????????????????????????? ??????????? ??????????????????????????????????????????????????????????? ??????????????????????????????? ??????????????????????????? (&i≠k)))┤?? （2）

從公式（2）中可以看出，度矩陣D的對角線元素為非0值，其余元素均為0。得到鄰接矩陣W和度矩陣D之后，通過公式（3）計算歸一化拉普拉斯矩陣L。

L=1/√(D^ )×(D-W)×1/√(D^ )?? （3）

（3）點(diǎn)云聚類。采用奇異值分解方法計算歸一化拉普拉斯矩陣L的特征值和特征向量，并選出最小的K個特征值對應(yīng)的特征向量ui（i =1, 2, …, K），每個特征向量ui均為n維列向量，K為用戶輸入的類別數(shù)量參數(shù)，本研究設(shè)K=80。將上述K個特征向量組合成n×K大小的特征矩陣F =［u1, u2, …, uK］，并將F按行標(biāo)準(zhǔn)化。如果將F中的每個K維行向量看作1個數(shù)據(jù)樣本，則F矩陣中共有n個數(shù)據(jù)樣本。采用K-means聚類方法將這n個數(shù)據(jù)樣本分為K個類別，其中第i個數(shù)據(jù)的類別即為植株點(diǎn)云P中第i個點(diǎn)pi的類別。最后將類別相同的點(diǎn)劃分為同一子區(qū)域，則植株點(diǎn)云最終被劃分為K個子區(qū)域點(diǎn)集。用符號?表示子區(qū)域點(diǎn)集，第m個子區(qū)域用符號?m（m=1，2，…，K）表示。本研究以點(diǎn)云法向量作為特征，將法向量相似的鄰近點(diǎn)云聚類成同一子區(qū)域。

不同株型植株的點(diǎn)云超體素聚類結(jié)果如圖2所示，從圖中可以看出，植株的每個器官均被分解成多個邊界光滑的子區(qū)域，且每個子區(qū)域的整體形態(tài)相對簡單，便于用簡單的幾何特征進(jìn)行表示和識別。例如，葉片上的子區(qū)域近似平面，莖上的子區(qū)域近似圓柱。但在個別器官交界處的位置，會存在分屬不同器官的點(diǎn)云被誤劃分為同一子區(qū)域的現(xiàn)象，如部分葉領(lǐng)點(diǎn)云與莖形成同一子區(qū)域、部分雄穗點(diǎn)云與頂端葉片形成同一子區(qū)域。

2.2.2植株頂部子區(qū)域提取

雄穗位于玉米植株的頂部區(qū)域，因此，可通過點(diǎn)云在植株中的相對位置對雄穗點(diǎn)云進(jìn)行大體定位。一些研究利用深度相機(jī)獲取植株點(diǎn)云數(shù)據(jù)［29,30］，在數(shù)據(jù)獲取時可通過調(diào)整相機(jī)位置使植株點(diǎn)云所在坐標(biāo)系的某個坐標(biāo)軸與玉米植株的生長方向重合，從而通過點(diǎn)云坐標(biāo)值即可初步對植株頂部的子區(qū)域進(jìn)行識別。但利用手持激光三維掃描儀、多角度圖像重建方法獲取的點(diǎn)云數(shù)據(jù)屬于散亂點(diǎn)云，其與所在空間坐標(biāo)系坐標(biāo)軸之間的相對位置是隨機(jī)的，因此增大了對植株頂部區(qū)域的提取難度。本研究采用如下步驟對植株頂部區(qū)域進(jìn)行定位。

（1）擬合植株點(diǎn)云的主成分向量。采用主成分分析法提取植株點(diǎn)云P的主成分方向。首先求植株的中心點(diǎn)pc。

p_c=1/n ∑_(i=1)^n?p_i （4）

之后利用公式（5）構(gòu)建協(xié)方差矩陣C。

C=1/n ∑_(i=1)^n?〖(p_i-p_c )×(〗 p_i-p_c )^T （5）其中，矩陣C為3×3對稱實(shí)方陣。

采用奇異值分解算法對C進(jìn)行特征值分解，得到3個非負(fù)特征值λ1、λ2、λ3 （λ1≥λ2≥λ3）以及對應(yīng)的3個特征向量I1、I2、I3。特征向量I1、I2、I3分別為植株點(diǎn)云的第1、2、3主成分向量。

（2）提取候選的頂部子區(qū)域集合。以植株點(diǎn)云第1主成分向量I1為方向向量，過植株中心點(diǎn)pc的直線參數(shù)方程L（t）為：

L(t)=I_1 t+p_c??? （6）

利用公式（7）計算P中每個點(diǎn)pi在直線L（t）上的投影坐標(biāo)ti。

t_i=(p_i-p_c)?I_1 （7）

找到所有點(diǎn)投影坐標(biāo)的最大值t1和最小值t2，之后利用公式（8）和公式（9）計算2個閾值t3和t4用于提取植株兩端區(qū)域的點(diǎn)云。

t_3=0.9t_1+0.1t_2? （8）

t_4=0.9t_2+0.1t_1? （9）

定義子區(qū)域集合V1和V2，并遍歷所有點(diǎn)云子區(qū)域。如果子區(qū)域中存在投影坐標(biāo)>t3的點(diǎn)云，則將該子區(qū)域加入集合V1。如果子區(qū)域中存在投影坐標(biāo)

（3）識別頂部子區(qū)域。成熟期玉米植株底部區(qū)域的莖為圓柱形，具有較明顯的直線線性特征。與之相反，植株頂部區(qū)域包含了雄穗，呈現(xiàn)一定的掃帚狀，直線線性特征較弱?；谏鲜鎏攸c(diǎn)，本研究通過計算子區(qū)域中點(diǎn)云的直線線性特征從V1和V2中選出頂部子區(qū)域集合。

定義1個子區(qū)域點(diǎn)集的直線線性特征為集合內(nèi)每個點(diǎn)的直線線性特征的平均值。為計算子區(qū)域點(diǎn)集?中一點(diǎn)p的直線線性特征，首先在?中查找與p點(diǎn)距離小于r1的所有點(diǎn)，組成點(diǎn)集B1。其中r1為用戶可調(diào)節(jié)參數(shù)，本研究設(shè)r1=10 cm；之后采用主成分分析法，計算B1點(diǎn)集的3個特征值λ1、λ2、λ3 （λ1≥λ2≥λ3），則p點(diǎn)的直線線性特征f1通過公式（10）計算。

f_1=((λ_1-λ_2))/λ_1 ? （10）

利用上述方法分別計算V1、V2中每個子區(qū)域的直線線性特征，并求V1中所有子區(qū)域直線線性特征的平均值fv1以及V2中所有子區(qū)域直線線性特征的平均值fv2，比較fv1和fv2的大小，較小值對應(yīng)的子區(qū)域集合為植株頂部區(qū)域，較大值對應(yīng)的區(qū)域集合為植株底部區(qū)域。V1、V2集合的子區(qū)域直線線性特征熱力圖如圖3（b）所示，從圖中可以看出，頂部區(qū)域的雄穗子區(qū)域具有非常小的直線線性特征值，導(dǎo)致該候選區(qū)域的直線線性特征的平均值大幅度降低。而莖則具有較大的直線線性特征值，從而使底部區(qū)域的特征值保持在較高的數(shù)值。

2.2.3雄穗點(diǎn)云識別

上節(jié)提取的頂部子區(qū)域集合（圖4（a））包含了3類子區(qū)域，分別為只包含雄穗點(diǎn)的雄穗子區(qū)域、只包含葉片點(diǎn)的葉片子區(qū)域以及同時包含雄穗點(diǎn)和葉片點(diǎn)的混合子區(qū)域。為實(shí)現(xiàn)雄穗點(diǎn)云的精確分割，需要在頂部區(qū)域中識別出雄穗子區(qū)域，并在混合子區(qū)域中識別出雄穗點(diǎn)云。由于葉片比雄穗具有更明顯的平面特征，因此本研究通過子區(qū)域內(nèi)點(diǎn)云的平面特征來識別雄穗點(diǎn)云和葉片點(diǎn)云，具體步驟如下。

（1）提取子區(qū)域內(nèi)點(diǎn)云的平面特征。為計算子區(qū)域點(diǎn)集?中一點(diǎn)p的平面特征，在?中找到與p點(diǎn)距離小于r2的所有點(diǎn)，組成點(diǎn)集B2，其中r2為用戶可調(diào)節(jié)參數(shù)；之后計算B點(diǎn)集的3個特征值λ1、λ2、λ3 （λ1≥λ2≥λ3），則p點(diǎn)的平面特征f2（p）通過公式（11）計算。

f_2 (p)=((λ_2-λ_3))/λ_1 （11）

子區(qū)域點(diǎn)集?的平面特征f2（?）為集合內(nèi)每個點(diǎn)的平面特征的平均值。植株頂部點(diǎn)云平面特征的熱力圖如圖4（b）所示，特征值越大，表示平面特征越弱。

（2）點(diǎn)云分類。計算所有子區(qū)域平面特征的平均值?(f_2 )和方差S2。如果子區(qū)域平面特征> ?(f_2 )，則該子區(qū)域被分類為葉片子區(qū)域；如果子區(qū)域平面特征≤ ?(f_2 )，則該子區(qū)域被分類為雄穗子區(qū)域，其內(nèi)點(diǎn)云被分類成雄穗點(diǎn)；如果子區(qū)域的平面特征f2（?）滿足 ?(f_2 ) < f2(?)< ?(f_2 )+ 0.5S2，則該區(qū)域?yàn)榛旌献訁^(qū)域。在混合子區(qū)域中，將點(diǎn)云平面特征f2（p）小于f2的點(diǎn)分類成雄穗點(diǎn)云。點(diǎn)云分類過程如圖4（c）所示。

（3）去噪。在混合子區(qū)域中，葉片邊緣的點(diǎn)云的平面特征也較小，因此在分類時容易被分類成雄穗點(diǎn)。但通常情況下，這些誤分的點(diǎn)與整個雄穗點(diǎn)云的距離較遠(yuǎn)，因此本研究采用基于統(tǒng)計分析的點(diǎn)云去噪方法［30］從雄穗點(diǎn)云中排除這些誤分點(diǎn)。

最終雄穗分割的結(jié)果和真值的結(jié)果分別如圖4（d）和4（e）所示。從圖中可以看出，本研究方法分割出了大部分的雄穗點(diǎn)云。但是，本研究的點(diǎn)云分類方法存在對子區(qū)域類別的誤識別現(xiàn)象，從而降低雄穗點(diǎn)云的分割精度。例如，圖4（a）中的子區(qū)域1為混合子區(qū)域，但其內(nèi)部葉片點(diǎn)云要遠(yuǎn)多于雄穗點(diǎn)云，因此該子區(qū)域平面特征值更多地體現(xiàn)了葉片點(diǎn)云的特征，導(dǎo)致該子區(qū)域被誤分類為葉片子區(qū)域，使得最終的分割結(jié)果會缺失該子區(qū)域內(nèi)的雄穗點(diǎn)。圖4（a）中的混合子區(qū)域2，其內(nèi)部葉片點(diǎn)云和雄穗點(diǎn)云之間的數(shù)量相近，因此可通過平面特征值正確識別子區(qū)域類別，進(jìn)而識別出該子區(qū)域內(nèi)的多數(shù)雄穗點(diǎn)云。

2.3分割精度評估

本研究使用精確度（Precision）、召回率（Recall）、F1分?jǐn)?shù)（F1 score）對分割精度進(jìn)行評估。對于雄穗點(diǎn)云，用TP表示屬于雄穗的點(diǎn)云被正確分割為雄穗的點(diǎn)數(shù)量，用FP表示不屬于雄穗的點(diǎn)云被誤分割成雄穗的點(diǎn)數(shù)量，用FN表示屬于雄穗的點(diǎn)云被誤分割為其他器官點(diǎn)云集合的點(diǎn)數(shù)量，則精確度、召回率和F1分?jǐn)?shù)的計算如公式（12）~（14）所示。

p_r=T_p/(T_p+F_p )?? （12）

R_e=T_p/(T_p+F_N )?? （13）

F_1=2×(P_r+R_e)/(P_r+R_e ) （14）其中，Pr表示精確度；Re表示召回率；F1表示F1分?jǐn)?shù)。一般來說，分割精確度和召回率之間是矛盾的，因此常采用F1分?jǐn)?shù)來對分割算法效果進(jìn)行綜合考量，越高的F1分?jǐn)?shù)說明分割效果越好。

3? 結(jié)果與分析

3.1分割效率和精度

在配置為2.2 GHZ CPU、DDR32 G內(nèi)存的筆記本工作站上進(jìn)行了算法測試。本研究在進(jìn)行點(diǎn)云超體素聚類時，利用奇異值分解方法對n×n大小的歸一化拉普拉斯矩陣進(jìn)行特征向量分解，導(dǎo)致整個算法的時間復(fù)雜度達(dá)到O（n3）。因此隨著點(diǎn)云密度和點(diǎn)云數(shù)量的增長，算法的運(yùn)行時間會顯著增加。當(dāng)點(diǎn)云密度為0.8個/cm時，成熟期玉米植株點(diǎn)云個數(shù)平均在6000個左右，此時算法的平均總運(yùn)行時間為2.26 s；點(diǎn)云密度為1.3個/cm時，植株點(diǎn)云數(shù)量平均約為15,000，平均總運(yùn)行時間為6.85 s；當(dāng)點(diǎn)云密度增大到1.9個/cm時，植株平均點(diǎn)云數(shù)量增加到30,000個左右，此時算法的總運(yùn)行時間增長到27.37 s。本研究方法的具體運(yùn)行時間如表1所示。

本研究算法的分割精度如表2所示。點(diǎn)云密度為0.8個/cm時，雄穗的平均F1分?jǐn)?shù)、平均精確率和平均召回率分別為0.763、0.746和0.832；點(diǎn)云密度為1.3時，上述3個指標(biāo)分別為0.875、0.898和0.867；點(diǎn)云密度為1.9個/cm時，3個分割指標(biāo)分別為0.889、0.950和0.849。從F1分?jǐn)?shù)指標(biāo)可以看出，隨著點(diǎn)云密度增加，本研究的分割精度逐漸增加，但點(diǎn)云密度1.3個/cm和1.9個/cm之間的分割效果差距不大。

點(diǎn)云密度和點(diǎn)云個數(shù)的增加，主要影響雄穗分割步驟中的點(diǎn)云平面特征計算結(jié)果。當(dāng)點(diǎn)云個數(shù)較少時，頂部葉片點(diǎn)云也相對稀疏，葉片點(diǎn)的平面特征和雄穗點(diǎn)的平面特征差距不明顯，造成點(diǎn)云類別誤判的情況增大，分割出的雄穗點(diǎn)云中包含較多葉片點(diǎn)云，分割精確度不高；當(dāng)點(diǎn)云個數(shù)增加時，葉片點(diǎn)云的平面特征更加明顯，因此，雄穗點(diǎn)和葉片點(diǎn)的區(qū)分度加大，分割出的雄穗點(diǎn)云中包含的葉片點(diǎn)云減少，分割精確度顯著提高。從表1的精確度指標(biāo)可以看出，隨著點(diǎn)云密度的增大，精確度顯著提高，而召回率則變化不大，F(xiàn)1分?jǐn)?shù)的增長主要來自于精確度的增加。

點(diǎn)云密度為1.9個/cm的玉米植株雄穗分割可視化結(jié)果和F1分?jǐn)?shù)如圖5所示。從圖中可以看出，大多數(shù)測試樣本均出現(xiàn)了雄穗分割不完整的情況，這說明本研究方法存在欠分割問題。當(dāng)雄穗完全高于頂部葉片或與葉片距離較遠(yuǎn)時，本研究方法分割精度非常高（如圖5（m）~圖5（o））；與之相反，當(dāng)雄穗被葉片包裹住或距離較近時（如圖5（a）~圖5（c）），分割精度較低。造成這種現(xiàn)象的原因是，如果雄穗與葉片距離較近，兩者之間的邊界點(diǎn)云的法向量和平面局部特征非常相似，難以通過點(diǎn)云超體素聚類以及平面特征進(jìn)行準(zhǔn)確分割。

3.2算法參數(shù)

影響本研究算法的參數(shù)主要為點(diǎn)云超體素聚類步驟中的K參數(shù)、植株頂部子區(qū)域提取中的r1參數(shù)和雄穗點(diǎn)云識別中的r2參數(shù)。K參數(shù)會影響點(diǎn)云超體素聚類得到的子區(qū)域個數(shù)，K值越大，子區(qū)域個數(shù)越大。本研究在后續(xù)植株頂部子區(qū)域提取和雄穗點(diǎn)云識別操作中，均以子區(qū)域作為點(diǎn)云集合，計算點(diǎn)的直線線性特征和平面特征。K值越小，每個器官被分成的子區(qū)域個數(shù)越少，每個子區(qū)包含多個器官點(diǎn)云的概率會增大，造成直線線性特征和平面特征對莖、葉片、雄穗點(diǎn)的區(qū)別度下降，降低分割精度。但K值越大，超體素聚類的計算時間會顯著增加。本研究中，通過試驗(yàn)確定參數(shù)K=80，該值對于不同株型、不同密度的玉米成熟期點(diǎn)云均適用，可將每個器官分解成多個邊界光滑且形態(tài)簡單的子區(qū)域，便于利用直線線性特征和平面特征對點(diǎn)云類別進(jìn)行判斷。

r1參數(shù)會影響點(diǎn)的直線線性特征計算結(jié)果。如果r1值較小，底部子區(qū)域莖點(diǎn)云在計算一個點(diǎn)的直線特征時，搜索得到的鄰域點(diǎn)僅包含莖圓柱面上的一個側(cè)面，使計算得到的直線特征值減少。相反，r1值越大，莖點(diǎn)云的直線特征值越大，底部子區(qū)域和頂部子區(qū)域的直線特征對比越明顯，但過大的r1值會增多鄰域點(diǎn)個數(shù)，降低特征值的計算效率。本研究設(shè)r1=10.0 cm，由于玉米莖的直徑通常不會大于10.0 cm，因此該值可以使每個莖點(diǎn)的鄰域點(diǎn)集合均為一個較完整的圓柱狀點(diǎn)云集合，使莖點(diǎn)云的直線特征值較大，同時，也不會過多地降低計算效率。r2參數(shù)會影響點(diǎn)的平面特征計算結(jié)果，進(jìn)而改變雄穗點(diǎn)云的識別結(jié)果。玉米葉片的平面特征具有局部性，即整個葉片表面呈現(xiàn)較明顯的彎曲特征，但局部小區(qū)域是較平展的。隨著r2值的增大，每個點(diǎn)搜索得到的鄰域點(diǎn)范圍會增大，則葉片點(diǎn)的平面性特征值會增大，降低葉片點(diǎn)和雄穗點(diǎn)之間的平面特征對比度，因此r2值不宜過大。由于頂部子區(qū)域包含的葉片均為上部葉的葉尖區(qū)域，該區(qū)域葉片寬度通常不會大于6.0 cm，同時頂部葉片經(jīng)常沿著葉脈卷曲，因此r2值要小于3.0 cm，這樣才能盡可能地保證葉片點(diǎn)云的平面特征較明顯。但過小的r2值會使鄰域點(diǎn)個數(shù)過少，很難描述局部區(qū)域的幾何特性，降低平面特征計算的準(zhǔn)確性，同時，雄穗的每個分支在極小區(qū)域內(nèi)也可能會呈現(xiàn)一定的平面特征，因此要保證r2值不能小于雄穗分支的寬度值。玉米雄穗每個分支的寬度通常會小于1.0 cm，因此r2值要大于1.0 cm。

本研究在設(shè)定r2具體參數(shù)值時，針對每一個點(diǎn)云密度的植株，均從r2? = 1.0 cm開始進(jìn)行參數(shù)選擇測試，按照間隔0.1 cm逐漸增大r2，直到r2 = 3.0 cm。選擇分割效果較好的r2參數(shù)進(jìn)行計算，當(dāng)多個參數(shù)值分割效果差距不大時，選擇其中的最小值作為r2參數(shù)，從而盡可能地減少鄰域點(diǎn)數(shù)量，提高平面特征計算效率。本研究通過該方式確定r2值，當(dāng)點(diǎn)云密度= 0.8個/cm時，r2設(shè)為2.4 cm；點(diǎn)云密度= 1.3個/cm時，r2設(shè)為1.8 cm；點(diǎn)云密度=1.9個/cm時，r2設(shè)為1.0 cm。隨著點(diǎn)云密度的減少，r2參數(shù)值明顯增大，這是由于點(diǎn)云密度較少時，為了找到足夠多的鄰域點(diǎn)描述局部平面特征，需要增大鄰域搜索半徑。當(dāng)給定一個新密度的植株點(diǎn)云時，可以利用本研究給出的r2值作為參考提高參數(shù)選擇效率。例如，如果植株點(diǎn)云密度為1.6個/cm，則r2值可在點(diǎn)云密度為1.9和1.3個/cm的參數(shù)值之間進(jìn)行測試，即r2值從1.0 cm開始進(jìn)行參數(shù)選擇，按照間隔0.1 cm逐漸增大r2到1.8 cm。

3.3展望

本研究算法目前不能較好地處理雄穗被葉片包裹住或距離較近的情況，此時葉片和穗之間的邊界點(diǎn)云的法向量和平面局部特征非常相似，造成分割精度下降。Jin等［31］開發(fā)了基于中值歸一化向量的區(qū)域生長莖葉分割算法，利用歐式距離作為特征，采用啟發(fā)式區(qū)域增長策略從器官底部的種子點(diǎn)開始，逐步將器官點(diǎn)云進(jìn)行完整分割，該方法證明了距離特征在玉米器官分割時的有效性。在后續(xù)的工作中，以本研究分割出的穗點(diǎn)云作為種子點(diǎn)，結(jié)合基于距離條件限定的區(qū)域增長方法，有望優(yōu)化最終的穗點(diǎn)云分割結(jié)果。

本研究算法只能對單株玉米植株點(diǎn)云進(jìn)行雄穗分割。為獲得單株點(diǎn)云，需將玉米植株移植入盆并采用激光掃描儀手動進(jìn)行數(shù)據(jù)獲取，數(shù)據(jù)獲取效率低。MVS Pheno平臺可快速、自動地獲取不用角度的玉米植株圖像，并利用運(yùn)動恢復(fù)結(jié)構(gòu)方法獲得質(zhì)量較好地玉米單株點(diǎn)云［32］。下一步可將本研究方法與MVS Pheno平臺進(jìn)行結(jié)合，加快雄穗表型檢測速度。

4? 結(jié)? 論

本研究提出一種基于超體素聚類和點(diǎn)云局部特征的玉米植株三維點(diǎn)云雄穗分割方法，驗(yàn)證了在玉米成熟期植株點(diǎn)云上識別、分割雄穗的可能性。通過分析本研究試驗(yàn)結(jié)果，得出以下結(jié)論。

（1）本研究方法的分割精度隨植株點(diǎn)云密度的增大而增加。當(dāng)點(diǎn)云密度為0.8、1.3和1.9個/cm時，雄穗點(diǎn)云分割的平均F1分?jǐn)?shù)分別為0.763、0.875和0.889。

（2）本研究算法的時間復(fù)雜度為O（n3），隨著點(diǎn)云密度和點(diǎn)云數(shù)量的增長，算法的運(yùn)行時間會顯著增加。

（3）綜合考慮分割精度和運(yùn)行時間，本研究對點(diǎn)云密度為1.3個/cm、點(diǎn)云的平均數(shù)量約為15,000的成熟期玉米植株點(diǎn)云處理效果最好，分割F1分?jǐn)?shù)可達(dá)到0.875，運(yùn)行時間為6.85 s。

結(jié)果表明，本研究方法具有在玉米植株散亂點(diǎn)云中提取雄穗點(diǎn)云的能力，可為玉米高通量表型檢測、玉米三維重建等研究和應(yīng)用提供技術(shù)手段。

參考文獻(xiàn)：

[1]??? TESTER M, LANGRIDGE P. Breeding technologies to increase crop production in a changing world[J]. Science, 2010, 327: 818-822.

[2]??? DUNCAN W G, WILLIAMS W A, LOOMIS R S. Tassels and the productivity of maize[J]. Crop Science, 1967, 7: 7-9.

[3]??? MEGHJI M R, DUDLEY J W, LAMBERT R J, et al. Inbreeding depression, inbred and hybrid grain yields, and other traits of maize genotypes representing three eras[J]. Crop Science, 1984, 24(3): 545-549.

[4]??? BROWN P J, UPADYAYULA N, MAHONE G S, et al. Distinct genetic architectures for male and female inflorescence traits of maize[J]. PLoS Genet, 2011,7(11): ID e1002383.

[5]??? ZHANG Y, ZHANG N. Imaging technologies for plant high-throughput phenotyping: A review[J]. Frontiers of Agricultural Science and Engineering, 2018, 5(4): 406-419.

[6]??? FERNANDA D M, GEMMA M, CAROLINA R A, et al. Yielding to the image: How phenotyping reproductive growth can assist crop improvement and production[J]. Plant Science, 2018, 282: 73-82.

[7]??? GAGE J L, MILLER N D, SPALDING E P, et al. TIPS: A system for automated image-based phenotyping of maize tassels[J]. Plant Methods, 2017, 13(1): 27-32.

[8]??? YE M, CAO Z, YU Z. An image-based approach for automatic detecting tasseling stage of maize using spatio-temporal saliency[C]// SPIE Conference on Multispectral Image Processing and Pattern Recognition. Washington, USA: SPIE Digital library, 2013.

[9]??? KURTULMU F, KAVDR S. Detecting corn tassels using computer vision and support vector machines[J]. Expert Systems with Applications, 2014, 41(16): 7390-7397.

[10]? LU H, CAO Z, XIAO Y, et al. Fine-grained maize tassel trait characterization with multi-view representations[J]. Computers and Electronics in Agriculture, 2015, 118: 143-158.

[11]? LU H, CAO Z, XIAO Y, et al. TasselNet: Counting maize tassels in the wild via local counts regression network[J]. Plant Methods, 2017, 13(1):1-17.

[12]? MAKANZA R, ZAMAN-ALLAH M, CAIRNS J E, et al. High-throughput method for ear phenotyping and kernel weight estimation in maize using ear digital imaging[J]. Plant Methods, 2018, 15(1): ID 49.

[13]? 王傳宇, 郭新宇, 吳升, 等. 基于計算機(jī)視覺的玉米果穗三維重建方法[J]. 農(nóng)業(yè)機(jī)械學(xué)報, 2014, 45(9): 274-279.WANG C, GUO X, WU S, et al. Three dimensional reconstruction of maize ear based on computer vision[J]. Transactions of the CSAM, 2014, 45(9): 274-279.

[14]? 溫維亮, 郭新宇, 楊濤, 等. 玉米果穗點(diǎn)云分割方法研究[J]. 系統(tǒng)仿真學(xué)報, 2017, 29(12): 3030-3034, 3041.WEN W, GUO X, YANG T, et al. Point cloud segmentation method of maize ear[J]. Journal of System Simulation, 2017, 29(12): 3030-3034, 3041.

[15]? 溫維亮, 王勇健, 許童羽, 等. 基于三維點(diǎn)云的玉米果穗幾何建模[J]. 中國農(nóng)業(yè)科技導(dǎo)報, 2016, 18(5): 88-93.WEN W, WANG Y, XU T, et al. Geometric modeling of maize ear based on three-dimensional point cloud[J]. Journal of Agricultural Science and Technology, 2016, 18(5): 88-93.

[16]? YU Z, ZHOU H, LI C. An image based automatic recognition method for the flowering stage of maize[C]// International Symposium on Multispectral Image Processing and Pattern Recognition. Washington, USA: SPIE Digital library, 2019.

[17]? BRICHET N, FOURNIER C, TURC O, et al. A robot-assisted imaging pipeline for tracking the growths of maize ear and silks in a high-throughput phenotyping platform[J]. Plant Methods, 2017, 13(1): ID 96.

[18]? JIA H, QU M H, WANG G, et al. Dough-stage maize (Zea mays L.) ear recognition based on multiscale hierarchical features and multifeature fusion[J]. Mathematical Problems in Engineering, 2020, 2020(2): ID 5582598.

[19]? PAULUS S, SCHUMANN H, KUHLMANN H, et al. High-precision laser scanning system for capturing 3D plant architecture and analysing growth of cereal plants[J]. Biosystems Engineering, 2014, 121: 1-11.

[20]? PAULUS S, DUPUIS, J, RIEDEL, S, et al. Automated analysis of barley organs using 3D laser scanning: An approach for high throughput phenotyping[J]. Sensors, 2014, 14(7): 12670-12686.

[21]? CH?N? Y, ROUSSEAU D, LUCIDARME P, et al. On the use of depth camera for 3D phenotyping of entire plants[J]. Computers and Electronics in Agriculture, 2012, 82: 122-127.

[22]? BUSEMEYER L, MENTRUP D, M?LLER K, et al. BreedVision－A multi-sensor platform for non-destructive field-based phenotyping in plant breeding[J]. Sensors, 2013, 13(3): 2830-2847.

[23]? LIN Y. LiDAR: An important tool for next-generation phenotyping technology of high potential for plant phenomics?[J]. Computers and Electronics in Agriculture, 2015, 119: 61-73.

[24]? GARRIDO M, PARAFOROS D S, REISER D, et al. 3D maize plant reconstruction based on georeferenced overlapping lidar point clouds[J]. Remote Sensing, 2015, 7(12): 17077-17096.

[25]? 韓東, 楊貴軍, 楊浩, 等. 基于立體視覺的玉米雄穗三維信息提取[J]. 農(nóng)業(yè)工程學(xué)報, 2018, 34(11): 166-173.HAN D, YANG G, YANG H, et al. Three dimensional information extraction from maize tassel based on stereoscopic vision[J]. Transactions of the CSAE, 2018, 34(11): 166-173.

[26]? MILLER F P, VANDOME A F, MCBREWSTER J. KD-TREE[M]. San Francisco, California, USA: Alpha Press, 2009.

[27]? RUSU R B, COUSINS S. 3D is here: Point cloud library (PCL)[C]// IEEE International Conference on Robotics & Automation. New Jersey, USA: IEEE Press, 2011.

[28]? GROTH D, HARTMANN S, KLIE S, et al. Principal components analysis[J]. Methods in Molecular Biology, 2013, 930: 527-547.

[29]? XIANG L, BAO Y, TANG L, et al. Automated morphological traits extraction for sorghum plants via 3D point cloud data analysis[J]. Computers and Electronics in Agriculture, 2019, 162: 951-961.

[30]? 仇瑞承, 張漫, 魏爽, 等. 基于 RGB-D相機(jī)的玉米莖粗測量方法[J]. 農(nóng)業(yè)工程學(xué)報, 2017, 33(S1): 170-176.QIU R, ZHANG M, WEI S, et al. Method for measurement of maize stem diameters based on RGB-D camera[J]. Transactions of the CSAE, 2017, 33(S1): 170-176.

[31]? JIN S, SU Y, WU F, et al. Stem-leaf segmentation and phenotypic trait extraction of individual maize using terrestrial LiDAR data[J]. IEEE Transactions on Geoscience and Remote Sensing, 2019, 57(3): 1336-1346.

[32]? WU S, WEN W, WANG Y, et al. MVS-Pheno: A portable and low-cost phenotyping platform for maize shoots using multiview stereo 3D reconstruction[J]. Plant Phenomics. 2020, 2020(2): ID 1848437.

Tassel Segmentation of Maize Point Cloud Based on Super Voxels Clustering and Local Features

ZHU Chao1,2, WU Fan1,2, LIU Changbin3, ZHAO Jianxiang1,2, LIN Lili1,2,

TIAN Xueying1,2, MIAO Teng1,2*

（1.College of Information and Electrical Engineering, Shenyang Agricultural University, Shenyang 110866, China;

2.Liaoning Engineering Research Center for Information Technology in Agriculture, Shenyang 110866, China;

3.Beijing PAIDE Science and Technology Development Company Limited, Beijing 100097,China）

Abstract: Accurate and high-throughput maize plant phenotyping is vital for crop breeding and cultivation research. Tassel-related phenotypic parameters are important agronomic traits. However, fully automatic and fine tassel organ segmentation of maize shoots from three-dimensional (3D) point clouds is still challenging. To address this issue, a tassel point cloud segmentation method based on point cloud super voxels clustering and local geometric features was proposed in this study. Firstly, the undirected graph of the maize plant point cloud was established, the edge weights were calculated by using the difference of normal vectors, and the spectral clustering method was used to cluster the point cloud to form multiple super voxel sub-regions. Then, the principal component analysis method was used to find the two end regions of the plant and based on the observation of the straight direction of the bottom stem regions, the top and bottom regions were distinguished by the point cloud linear features. Finally, the tassel points were identified based on the plane local features of the point cloud. The sub-regions of the top region of the plant were classified into leaf regions, tassel regions, and mixed regions by plane local features of the point cloud, the tassel points in the tassel sub-region, and the mixed region were the finally segmented tassel point clouds. In this study, 15 mature maize plants with 3 point cloud densities were tested. Compared with the ground truth segmented manually, the average F1 scores of the tassel segmentation were 0.763, 0.875 and 0.889 when the point cloud density was 0.8/cm, 1.3/cm, and 1.9/cm, respectively. The segmentation accuracy of this method increased with the increase of plant point cloud density. The increase of point cloud density and the number of point clouds mainly affected the calculation results of point cloud plane features in tassel segmentation. When the number of point clouds was small, the top leaf point cloud was relatively sparse. Therefore, the difference between the plane feature of the leaf point and the plane feature of the tassel point was not obvious, which led to the increase of the misclassification of the point cloud. However, the time complexity of the algorithm was O(n3), so the increase in the density and number of point clouds would lead to a significant increase in the running time. Considering the segmentation accuracy and running time, the research obtained the best effect on the mature maize plants with a point cloud density of 1.3/cm and an average number of 15,000. The segmentation F1 score reached 0.875 and the running time was 6.85 s. The results showed that this method could extract tassels from maize plant point cloud, and provided technical support for the research and application of high-throughput phenotyping and three-dimensional reconstruction of maize.

Key words: maize tassel; 3D point cloud segmentation; phenotyping; super voxels clustering; local features; principal component analysis