基于三維點云的植物葉片重建及其面積估算*

陳 柱，楊 君*

(1.武漢科技大學冶金自動化檢測技術(shù)教育部工程研究中心，湖北 武漢 430081；2.武漢科技大學信息科學與工程學院，湖北 武漢 430081)

植物葉片大小是監(jiān)測植物生長及預(yù)測植物生長的重要參數(shù)之一[1]，也是植物生長管理的重要依據(jù)；葉片面積是植物生長發(fā)育過程中的一個重要評價指標。 葉片面積測量由最初的回歸系數(shù)法[2]、方格紙法等發(fā)展到基于立體視覺三維點云的建模方法[3-8]。 基于立體視覺的三維點云重建技術(shù)為監(jiān)測植物生長和測算葉片面積提供了新方法。 文獻[3]利用Kinect v2 傳感器進行三維重建。 文獻[4]在機械臂上安裝攝像頭采集植株圖像，采用多視角立體成像(multi-view stero，MVS)算法生成點云。 MVS算法可應(yīng)用于任何種類或形態(tài)的植物，為植物自動化三維重建提供了一種方法，但文獻[4]存在圖像采集方法復(fù)雜的問題。 文獻[5]利用3DSOM 軟件進行植株三維重建，利用相機、三腳架、旋轉(zhuǎn)平臺等設(shè)備搭建圖像獲取系統(tǒng)，最終可獲得多角度拍攝的植株圖像，但當圖像數(shù)量較多時，存在重建模型效果較差且設(shè)備復(fù)雜的問題。 文獻[6]概述了運動恢復(fù)結(jié)構(gòu)算法和多視角立體算法的原理，同時指出算法重建結(jié)果可能會出現(xiàn)孔洞等問題。 文獻[7]利用3D 激光掃描儀從多個角度捕獲植物的表面點云數(shù)據(jù)，同時提出用于粗略初始配準的特征點匹配算法并建立了植物的三維模型，使用三角剖分算法計算得到重建網(wǎng)格的表面積和體積。 文獻[8]通過CNN網(wǎng)絡(luò)，從單幅RGB 圖像中得到待重建場景的深度圖、分割圖及平面參數(shù)，然后將平面深度圖與相機標定獲取的相機參數(shù)相結(jié)合，得到像素點的三維點云深度信息。

與上述方法相比，多視角拍攝二維圖像的運動恢復(fù)結(jié)構(gòu)(SFM)算法[9]更有優(yōu)勢，因為算法只需用到一臺相機且能自動處理相機標志定位問題。 文獻[10]采用運動恢復(fù)結(jié)構(gòu)算法獲取稀疏點云，并把點云分割、表面重建等算法結(jié)合，最終得到較精確的葉片面積。 SFM 算法在處理植株這類復(fù)雜結(jié)構(gòu)對象時，稀疏點云的三維結(jié)構(gòu)信息不足以完成植物的三維重建。 為解決上述問題，學者們把SFM 算法和MVS 算法結(jié)合，從而得到三維信息足夠多的稠密點云。

為處理葉片重疊導致葉面積測量不準確的問題，本文提出了一種植物葉片三維重建補償方法。該方法首先利用相機多角度拍攝植物圖像；其次結(jié)合SFM 算法、CMVS 算法與PMVS 算法處理圖像生成點云數(shù)據(jù)；然后對點云數(shù)據(jù)進行去噪、分割、填補、網(wǎng)格化處理；最后計算葉片面積。

1 方法與過程

本研究首先采用運動恢復(fù)結(jié)構(gòu)(structure from motion，SFM)算法對多角度拍攝的二維圖像進行處理，獲得葉片的三維稀疏點云數(shù)據(jù)，再采用聚類多視角立體(clustering multi-view stereo，CMVS)算法和基于面片的多視角立體(Patch-Based Multi-View Stereo，PMVS)算法處理三維稀疏點云得到三維稠密點云；其次通過顏色特征空間閾值方法去除三維稠密點云的噪點；然后采用K-means++聚類算法對葉片的三維稠密點云進行分割；分割后用Harris 角點檢測算法尋找分割重疊葉片的角點，通過角點對缺口葉片進行填補；最后用滾球算法對葉片進行網(wǎng)格化處理并計算葉片面積。 過程如圖1 所示，算法分為四個步驟，依次為輸入圖像、點云重建(SFM 算法，CMVS/PMVS 算法)、稠密點云優(yōu)化(點云去噪，點云分割，點云填補)、網(wǎng)格化及面積計算。

圖1 總體算法流程圖

1.1 實驗材料

實驗材料選取普通的樹葉葉片，實驗葉片分為四類:①葉片之間無重疊；②葉片小部分與另一葉片重疊；③葉片大部分與另一葉片重疊；④葉片交叉重疊。 本文僅考慮兩片葉片重疊的情況。 其中類1 有五片葉片，編號為a1、b1、c1、d1、e1；類2、3、4 均有兩片葉片，編號依次為a2、b2、a3、b3、a4、b4。 如圖2所示。

圖2 葉片重疊情況分類

圖3 為多角度拍攝的50 幅類1 葉片圖像。 每幅圖像尺寸為3 024×3 024，圖像為JPG 格式。 四種類別葉片各采集50 幅圖像。

圖3 50 幅類1 葉片圖像

1.2 點云重建

點云重建由稀疏點云重建與稠密點云重建兩個步驟組成。

1.2.1 稀疏點云重建

SFM 算法處理二維圖像形成三維稀疏點云。圖4 為SFM 算法流程。

圖4 SFM 算法流程圖

首先輸入圖像，采用尺度不變特征轉(zhuǎn)換(Scale invariant feature transform，SIFT)算法[11-12]對二維圖像進行特征提取與特征匹配；其次通過幾何驗證篩選圖像匹配對完成圖像配準；然后從圖像匹配對特征點的像素坐標開始，利用對極幾何原理，計算得到基礎(chǔ)矩陣；在基礎(chǔ)矩陣奇異分解運算后再進行三角化運算，得到特征點在世界坐標系下的坐標；最后使用光束法平差(Bundle adjustment，BA)最小化重投影誤差。

重投影誤差為真實三維空間點在圖像平面上的投影與重投影的差值。 重投影誤差圖如圖5 所示，軌跡j上的3D 點Xj投影至相機i的2D 圖像平面上形成的點I(Ci，Xj)，即為三維空間點在圖像平面上的投影點；虛擬三維空間點qij為重投影點。

圖5 投影誤差圖

當存在m個軌跡與n個拍攝視角時，重投影誤差累計和g(C，X)為式(1):

式中:當3D 點Xj在相機i平面上存在投影時，ωij為1，反之為0。 通過迭代的方式使累積和逐漸收斂到最小值，最終得到最優(yōu)的相機姿態(tài)與稀疏點云。

1.2.2 稠密點云重建

稠密點云由稀疏點云獲得法向量并向外擴散形成。 稠密點云重建分為兩步:

①在緊密性、圖像大小、覆蓋范圍的約束下，CMVS 算法將稀疏點云聚類到不同的圖像簇。

②在局部光度一致性、全局可視化一致性的約束下，PMVS 算法將稀疏點云匹配、擴展、過濾[13]得到稠密點云。

本文通過VisualSFM 軟件和MeshLab 軟件實現(xiàn)SFM 算法和CMVS/PMVS 算法。

圖6 為三類重疊葉片的稀疏點云重建結(jié)果、稠密點云重建結(jié)果和拍攝葉片時相機的空間位置。 圖6 左半部分為稀疏點云重建結(jié)果，對三類葉片圖像進行稀疏點云重建后，分別得到15 567、14 832、17 963 個 特 征 點， 生 成1 738 455、1 584 126、1 984 136個稀疏點。 圖6 右半部分為稠密點云重建結(jié)果，生成的稠密點云較稀疏點云能更好的還原葉片真實形態(tài)。

圖6 稀疏點云結(jié)果(左)與稠密點云結(jié)果(右)

1.3 點云去噪

在重建過程中，葉片存在白色背景或葉片表面反光導致葉片邊緣存在噪點，如圖6 葉片點云邊緣所示。 本研究以點云顏色與噪點顏色存在差異為依據(jù)濾除噪聲。

本文采用MATLAB Color Threshold 工具箱提取點云的HSV 顏色特征閾值[10]，歸一化后三個通道的閾值范圍依次為[0，0.225]、[0.395，0.969]、[0.067，0.585]，然后利用閾值范圍生成掩膜矩陣。

我們先將點云數(shù)據(jù)的顏色矩陣作為初始對象生成新的顏色矩陣，以掩膜矩陣為索引生成所有點新的坐標矩陣；后將新的顏色矩陣和新的坐標矩陣寫為點云格式，生成去噪后的葉片稠密點云。 去噪結(jié)果如圖7 所示，圖7 與圖6 對比，去噪后葉片稠密點云邊緣噪聲點顯著減少。 三類葉片點云數(shù)剩余1 553 845、1 321 255、1 626 512 個。

圖7 葉片點云去噪結(jié)果

1.4 點云分割及填補

本文采用無監(jiān)督聚類K-means++算法[14]來實現(xiàn)點云分割，對三類重疊葉片點云聚類分割的結(jié)果如圖8 所示，完整葉片點云為深色，缺口葉片點云為淺色。

圖8 三類重疊葉片分割圖

類1 葉片能正常分割，不用填補，其他重疊葉片均需填補。

填補缺口葉片需找到缺口葉片點云的角點。 角點是連接物體輪廓線的點，它可以展現(xiàn)圖像的圖形特征并幫助減少運算量。 本文采用Harris 角點檢測算法[15]尋找待填補葉片點云的角點。 Harris 角點檢測算法選取某像素的領(lǐng)域窗口，當窗口在各個方向上滑動且窗口內(nèi)像素平均灰度值產(chǎn)生較大變化時，則認為選取的像素點為角點。

以類2 葉片為例，填補算法步驟為:

①采用Harris 角點檢測算法找到類2b2葉片點云處于缺口邊緣的角點e21、e22，連接e21、e22形成線段l2，如圖9(a)所示；

②移動線段l2至同一枝葉的另一完整葉片點云a2上，線段l2構(gòu)成的平面將完整葉片點云截取為兩部分s21、s22，規(guī)定s21為大面積，s22為小面積，如圖9(b)、圖9(c)所示；

③根據(jù)缺口葉片待填補大小選擇s21或s22填補至缺口葉片點云b2，如圖9(d)所示，選擇s22；

④若填補后仍有空洞，則采用漫水填充算法[16]填充孔洞，填補算法完成，如圖9(e)所示。

由圖9(e)可看出填補算法對類2b2葉片進行了有效填補，使b2葉片的空缺部分得到填補。 類3、類4 葉片同理，如圖10、圖11 所示。

圖9 類2 b2 葉片點云填補圖

圖10 類3 b3 葉片點云填補圖

圖11 類4 b4 葉片點云填補圖

1.5 網(wǎng)格化及面積計算

1.5.1 網(wǎng)格化

本文采用滾球算法對稠密點云進行三角網(wǎng)格化，使用自適應(yīng)半徑的方法使球在搜索過程中更容易遇到新的點，可避免僅使用單一半徑進行數(shù)據(jù)點搜索而留下孔洞。

滾球算法步驟為:首先，在點云中隨機選取一點，尋找與其距離最近的兩個點并將三點連接形成三角形；其次，將三角形的三條邊作為網(wǎng)格的邊界并列為備選邊，選擇三角形的一條邊作為軸，使用自適應(yīng)半徑的方法使軸旋轉(zhuǎn)生成球；然后，球在點云表面滾動且接觸到新的點時，將該點和上述備選邊連接形成新三角形并列入重建網(wǎng)格，將新形成的邊列為備選邊，同時將上一步生成的邊取消備選；最后，再從備選邊中選擇下一邊重復(fù)上述操作，直至所有點都參與滾動并生成三維三角網(wǎng)格。

類1a1葉片通過滾球算法得到的三維三角網(wǎng)格如圖12 所示，葉片三角網(wǎng)格與真實葉片接近。

1.5.2 面積計算

本文在網(wǎng)格化后對所有三角網(wǎng)格面積求和，葉片面積是上述所有三角網(wǎng)格面積的累加。

以網(wǎng)格中△BCD 為例，如圖13 所示，根據(jù)△BCD 三點的坐標可確定各邊邊長為:

圖13 BCD 三角網(wǎng)格圖

式中:lBC、lBD、lCD為△BCD 的三條邊長，(xB，yB，zB)、(xC，yC，zC)、(xD，yD，zD)分別為點B、C、D 的坐標。

根據(jù)海倫凱勒公式可得△BCD 的面積為:

式中:S△BCD為△BCD的面積，p為半周長，值為(lBC+lCD+lBD)/2。

按上述方法可求得網(wǎng)格中所有三角形的面積，最后累加所有三角形的面積可得近似的葉片面積，累加公式如式(6):

式中:S葉為葉片面積，SΔ為單個三角網(wǎng)格的面積，t為三角網(wǎng)格序號，n為網(wǎng)格總數(shù)。

2 結(jié)果與分析

以掃描法為標準，將本文方法與方格紙法進行對比。 掃描法測量葉片面積時，將單個葉片按壓平放至掃描儀1 ∶1 掃描。 對掃描的圖像進行灰度處理，統(tǒng)計葉片像素數(shù)占紙張像素數(shù)的比例，已知紙張面積，通過比例計算得到葉面積。 方格紙法選取方格面積為0.25 cm2方格紙。 統(tǒng)計單個葉片占方格紙的數(shù)量，將滿格和大于半格計為0.25 cm2，小于半格計為0 cm2。 圖14、圖15 依次為方格紙法和掃描法處理類1a1葉片得到的結(jié)果。

圖14 類1 a1 葉片方格紙法結(jié)果

圖15 類1 a1 葉片掃描法結(jié)果

掃描法的精度接近真實值，故本文把掃描法的結(jié)果作為對比標準。

方格紙、掃描法和本文方法計算所有編號葉片的面積，結(jié)果如表1。

表1 葉面積計算結(jié)果 單位:cm2

由表1 數(shù)據(jù)可看出，本文方法相對于方格紙法更接近精度較高的掃描法。

表2 第三列是對重疊葉片的面積，第四列是對重疊葉片填補的面積。

表2 本文方法填補前后葉面積 單位:cm2

由表2 數(shù)據(jù)可知，三類重疊的b 葉片在填補處理后面積增加了3.810 cm2、27.337 cm2、18.116 cm2。填補后的結(jié)果更接近掃描法，驗證了填補處理的有效性。

3 結(jié)論

本研究采用三維點云重建的方法測算葉片面積，對重疊葉片進行填補來處理面積測量不準確的問題。 填補處理能有效解決由于重疊導致的葉面積測量不準確的問題，與方格紙法相比，本文方法更接近掃描法，并有效地還原真實葉片面積。