基于無人機(jī)航拍圖像的大田玉米冠層結(jié)構(gòu)建模

朱冰琳 李 敏 劉扶桑 賈奧博 毛 秀 郭 焱

(中國農(nóng)業(yè)大學(xué)土地科學(xué)與技術(shù)學(xué)院， 北京 100193)

0 引言

高效選育高產(chǎn)的作物新品種是解決世界糧食安全問題的重要途徑[1]。影響現(xiàn)代育種的關(guān)鍵因素之一是能否快速、精確地對大量育種品系進(jìn)行全面的表型評價，從中篩選出具有期望表型的品系[2-3]。目前，基因標(biāo)記輔助選擇、基因組編輯[1]、高通量基因測序等技術(shù)已得到廣泛應(yīng)用，利用這些技術(shù)可以高效地挖掘調(diào)控作物重要性狀的基因信息[4]。然而，作物表型的研究技術(shù)和方法卻無法與高通量基因測序技術(shù)的發(fā)展相匹配[5-6]。

近年來興起的高通量獲取作物表型的方法或“表型組學(xué)”，被認(rèn)為可以顯著提高植物表型信息的獲取效率[7]。目前開發(fā)的表型平臺大多應(yīng)用于溫室等特定環(huán)境中[8-10]。田間表型平臺也逐漸開發(fā)出來，但往往體積龐大、價格昂貴[11]，且存在較大的應(yīng)用局限。如星載表型平臺較難應(yīng)用于面積較小的育種實(shí)驗(yàn)小區(qū)，不適合于高頻率的動態(tài)監(jiān)測[7,12]；機(jī)載平臺的花費(fèi)較大，需要較多的人員維護(hù)，且在低空獲取高分辨率圖像時對冠層的擾動較大[13]；地面表型平臺能夠同時獲取多光譜、高光譜及熱紅外等圖像信息，但通常為特定株高范圍的作物設(shè)計(jì)，其行走受限于田間地面狀態(tài)，靈活性和便攜性較差[7]。

目前，基于RGB相機(jī)獲取植株多視角圖像序列后，利用運(yùn)動恢復(fù)結(jié)構(gòu) (Structure from motion, SFM)算法可重建植株三維結(jié)構(gòu)。該方法已在溫室中廣泛應(yīng)用，如監(jiān)測溫室內(nèi)小麥、茄子、黃瓜、辣椒等生長過程[14-15]，對大田玉米也有初步的嘗試[16]，但重建范圍較小。無人機(jī)(Unmanned arial vehicle, UAV)平臺搭載RGB相機(jī)，以數(shù)十米或更高的飛行高度獲取冠層圖像序列，可以實(shí)現(xiàn)大田尺度的冠層點(diǎn)云重建[17]，并可進(jìn)一步計(jì)算株高[18]、植被指數(shù)(Normalized difference vegetation index, NDVI)[19-20]、葉面積指數(shù)(Leaf area index, LAI)[21]、植被覆蓋度[17]等群體冠層參數(shù)。然而，實(shí)施精準(zhǔn)農(nóng)業(yè)需要有作物冠層結(jié)構(gòu)精確信息的支持[22]，基于理想株型選育作物新品種[23]也需要精確的冠層結(jié)構(gòu)信息[24]。隨著無人機(jī)的微型化發(fā)展，其對冠層的擾動變小，可以通過超低空飛行獲取更清晰的冠層圖像，使快速獲取高精度的冠層結(jié)構(gòu)信息成為可能。

通過點(diǎn)云數(shù)據(jù)網(wǎng)格化可建立三維網(wǎng)格的冠層結(jié)構(gòu)模型。文獻(xiàn)[25]利用Geomagic Studio軟件實(shí)現(xiàn)了蘋果樹葉片的網(wǎng)格重建；文獻(xiàn)[26]對黃瓜、西瓜、甜瓜等植株的葉片利用Delaunay三角剖分的方法生成初始網(wǎng)格曲面，再對網(wǎng)格曲面進(jìn)行優(yōu)化處理，實(shí)現(xiàn)了冠層結(jié)構(gòu)建模；文獻(xiàn)[27]采用改進(jìn)的Crust方法實(shí)現(xiàn)了樹干的網(wǎng)格重建。這些研究大多在溫室內(nèi)進(jìn)行，且葉片形態(tài)較為簡單，受遮擋較少；在室外選擇的研究對象(如樹干)，其形態(tài)也不易受到風(fēng)的影響，重建效果較好。目前，基于點(diǎn)云數(shù)據(jù)對大田作物群體冠層結(jié)構(gòu)建模的研究較少。

本文利用超微小型無人機(jī)超低空航拍獲取大田不同生育期植株個體及群體的多視角圖像序列，基于偽極點(diǎn)-Crust算法構(gòu)建其冠層結(jié)構(gòu)模型，并基于大田實(shí)測數(shù)據(jù)對重建冠層結(jié)構(gòu)模型的精度進(jìn)行評估，以期為高效、高精度構(gòu)建大田作物冠層結(jié)構(gòu)模型、獲取特定的作物表型參數(shù)提供新途徑。

1 材料與方法

1.1 實(shí)驗(yàn)區(qū)概況

玉米大田實(shí)驗(yàn)于2019年在中國農(nóng)業(yè)大學(xué)吉林梨樹實(shí)驗(yàn)站(43°16′N，124°26′E)進(jìn)行。供試品種為先玉335，種植株距為20 cm，行距為50 cm。共3個小區(qū)，小區(qū)長、寬為24、9 m(圖1a)。種前施加基肥，用量為純N 80 kg/hm2、P2O5120 kg/hm2、K2O 100 kg/hm2，采用鏵式犁和旋耕機(jī)將肥料旋進(jìn)土壤耕層中。在玉米第6葉及第11葉展開時追施40 kg/hm2純N。大田出苗日為5月19日。

1.2 冠層航拍及田間測量

出苗后43 d(苗期)，選取一個先玉335小區(qū)的玉米群體為目標(biāo)，隨機(jī)選擇相鄰的4個條帶，并從中選取18棵植株(圖1b)，分別對每株上、中、下部位各葉片用標(biāo)簽紙標(biāo)記。為確定幾何換算系數(shù)，在目標(biāo)群體中心的地面上放置邊長為40 cm、啞光材質(zhì)的正方形校正板。在多云且風(fēng)速小的天氣條件下，采用搭載哈蘇相機(jī)的Mavic 2 Pro型(大疆創(chuàng)新科技有限公司，中國)超微小型無人機(jī)對目標(biāo)群體進(jìn)行航拍。航拍時采用圓周型航線飛行2圈，飛行高度及半徑均設(shè)置為5 m。相機(jī)鏡頭與水平面的夾角為45°，飛行速度為0.8 km/h，相鄰圖像的拍攝間隔為2 s，圖像重疊度為80%。每次航拍獲取160～180幅圖像(表1)。在出苗后112 d(成熟期)進(jìn)行了同樣方式的航拍，由于此時植株中下部相互遮蔽嚴(yán)重，因而選取一株、一行內(nèi)相鄰2、3、4株玉米為航拍目標(biāo)(圖1b)，航拍前去除目標(biāo)植株周邊的植株，并標(biāo)記葉片、放置校正板。每次航拍后，立即測量目標(biāo)植株標(biāo)記葉片的葉長、葉最大寬及株高(植株基部與其最高點(diǎn)間的垂直距離)。

表1 無人機(jī)平臺、飛行參數(shù)及圖像信息Tab.1 UAV platform, flight parameters and image information

1.3 數(shù)據(jù)預(yù)處理

使用Metashape Professional (Agisoft, RU)軟件對航拍獲取的多角度圖像序列(圖2a)進(jìn)行冠層點(diǎn)云重建。獲得點(diǎn)云后(圖2b、2c)，手動將選定的目標(biāo)植株分割出來，利用支持向量機(jī)(Support vector machine, SVM)分類器根據(jù)顏色信息進(jìn)行點(diǎn)云分類(圖2d)，刪除后得到僅包含植被顏色信息的點(diǎn)云(圖2e)。在成熟期時，由于地面有部分與植被顏色相近的落葉，采用隨機(jī)抽樣最大似然估計(jì)算法(Maximum likelihood estimation sample consensus, MLESAC)[28]對地面落葉點(diǎn)云進(jìn)行平面擬合(圖2f)，去除后得到植株點(diǎn)云(圖2g)。

由于點(diǎn)云數(shù)據(jù)來源于大田航拍圖像，因而不可避免地帶有噪聲點(diǎn)。為消除噪聲點(diǎn)對冠層結(jié)構(gòu)構(gòu)建的干擾，本文采用基于空間單元格搜索k-近鄰點(diǎn)的方法剔除噪聲點(diǎn)。設(shè)P={p1,p2, …，pn}是重建冠層上所有點(diǎn)的集合，P中與pi歐氏距離最小的k個點(diǎn)稱為pi的k-近鄰，記為N(pi)。該算法的核心思想是對點(diǎn)云數(shù)據(jù)的最大包圍盒進(jìn)行一定步長(本文設(shè)定為0.1 m)的空間劃分(圖2h)，通過計(jì)算一個點(diǎn)p所在的子包圍盒及其在空間上相鄰的26個子包圍盒(圖2i)中所有點(diǎn)的歐氏距離，由小到大進(jìn)行排序，選取前k個距離，即為該點(diǎn)與k個近鄰點(diǎn)間歐氏距離的集合。通過k-近鄰的歐氏距離平均值D(pi)來衡量此點(diǎn)的局部鄰域關(guān)系，根據(jù)點(diǎn)集P的整體鄰域關(guān)系確定距離閾值。如果D(pi)在設(shè)定閾值的范圍外，則認(rèn)為是噪聲點(diǎn)[25](圖2j)，刪除后得到最終的植株點(diǎn)云(圖2k)。

1.4 冠層結(jié)構(gòu)模型構(gòu)建

Crust算法為基于點(diǎn)云重建三維網(wǎng)格結(jié)構(gòu)的經(jīng)典算法之一。采用該算法時，先對點(diǎn)云中所有點(diǎn)進(jìn)行一次三維Delaunay剖分，得到樣本點(diǎn)集的Voronoi頂點(diǎn)集。在此基礎(chǔ)上計(jì)算各點(diǎn)的正負(fù)極點(diǎn)，將得到的極點(diǎn)集與已有點(diǎn)云點(diǎn)集的并集再次進(jìn)行三維Delaunay剖分[29], 采用傳統(tǒng)Crust算法計(jì)算量很大，效率較低。為了提高剖分的效率，本文在構(gòu)建冠層結(jié)構(gòu)模型時，在去噪后的點(diǎn)云中添加少量的偽極點(diǎn)，利用偽極點(diǎn)代替正負(fù)極點(diǎn)構(gòu)建原始點(diǎn)集，從而較大程度地減少剖分的計(jì)算量，將其稱之為偽極點(diǎn)-Crust算法[27]。

執(zhí)行偽極點(diǎn)-Crust算法時，首先確定樣本點(diǎn)云的最大包圍盒，并以包圍盒X、Y、Z方向上距離最大值的M倍(本文將M設(shè)置為0.2)作為步長向外拓展得到更大的四面體包絡(luò)面[27]。分別在包絡(luò)面上對各邊N等分處(本文N取4)均勻插值形成6N2個(即96個)偽極點(diǎn)集合(圖3a，以1行4株玉米點(diǎn)云群體偽極點(diǎn)的建立為例)。然后對偽極點(diǎn)與已有點(diǎn)云點(diǎn)集的并集進(jìn)行三維Delaunay剖分，得到初始四面體集合To。該集合中存在較多非冠層結(jié)構(gòu)的三角面片，需通過以下步驟進(jìn)行篩選。

(1)冠層結(jié)構(gòu)四面體的篩選遍歷集合To，判斷各四面體的頂點(diǎn)中是否含有偽極點(diǎn)。

如果包含偽極點(diǎn)，則將其標(biāo)記為“非冠層結(jié)構(gòu)”。然后，獲取To中所有的三角面片并剔除重復(fù)的三角面片得到集合T′o。遍歷T′o，判斷每個三角面片與To中各初始四面體的映射關(guān)系，標(biāo)記不在集合邊界的所有三角面片。為評估非邊界三角面片所在的2個相鄰四面體間的相似度，采用了由這2個相鄰四面體各自外接球半徑ri1、ri2之間夾角θ的余弦值表征的交集因子F(Intersection factor) (圖3b)

(1)

其中

F(i)∈[-1,1] (i=1, 2, …,j)

式中j——所有非邊界三角面片所在的四面體的個數(shù)

ci1、ci2——2個相鄰四面體外接球的球心坐標(biāo)

設(shè)置F的閾值Ft為0.95，F(xiàn)大于Ft表明θ小，即2個外接球重合的部分很大(圖3b)；F小于-Ft表明θ大，即2個外接球重合的部分很小(圖3b)。遍歷具有相同非邊界三角面片的相鄰四面體，通過標(biāo)記為“非冠層結(jié)構(gòu)”屬性的四面體及F與Ft的大小關(guān)系判斷另一個四面體的屬性。如果遍歷一次后未能標(biāo)記所有四面體的屬性，則將閾值逐步降低，直至每個四面體的屬性標(biāo)記完畢。將屬性為“非冠層結(jié)構(gòu)”的四面體剔除，得到屬于冠層結(jié)構(gòu)的四面體集合Tp。

(2)冠層三角面片的提取

從得到的冠層四面體集合Tp中取出所有三角面片，計(jì)算各三角面片3個頂點(diǎn)Z坐標(biāo)之和，以其和最大的三角面片為起始面片，確定其法向量的方向，然后計(jì)算其鄰接三角面片的法向量，將法向量變化最小的三角面片作為其相鄰的三角面片，直至遍歷所有三角面片。

由上述篩選過程得到的三角面片集構(gòu)建初始冠層結(jié)構(gòu)模型(圖3c)。由于在大田環(huán)境下獲取的點(diǎn)云密度不均一，利用該方法構(gòu)建的初始冠層結(jié)構(gòu)模型中存在大量的邊長畸長的異常面片(圖3c)。為此，計(jì)算了所構(gòu)建模型中所有三角面片最大邊的長度，統(tǒng)計(jì)確定邊長的閾值后，剔除最大邊長度大于閾值的三角面片(圖3d)。為消除所得到的器官曲面不光滑問題，將面片的頂點(diǎn)沿其法向量方向根據(jù)平均曲率進(jìn)行移動平滑[30]，得到最終的冠層結(jié)構(gòu)模型(圖3e)。

1.5 基于表型數(shù)據(jù)的模型精度評估

根據(jù)校正板點(diǎn)云的計(jì)算邊長及實(shí)際邊長確定長度換算系數(shù)后，基于所構(gòu)建模型計(jì)算目標(biāo)植株的株高；利用Geomagic Studio (Geomagic，美國)手動提取標(biāo)記葉片，計(jì)算其葉長及最大葉寬；對組成葉片的每個小三角面片面積進(jìn)行累加，獲得標(biāo)記葉片的葉面積。基于大田標(biāo)記葉的葉長和葉最大寬的測量值計(jì)算參比葉面積為

Am=αLmWm

(2)

式中Lm——標(biāo)記葉的葉長

Wm——最大葉寬

α——經(jīng)驗(yàn)系數(shù)，取0.75

通過比較實(shí)際測量的株高、葉長、葉最大寬、參比葉面積與由冠層結(jié)構(gòu)模型得到的計(jì)算值，以決定系數(shù)(R2)、均方根誤差(Root mean square error，RMSE)、相對均方根誤差(Root mean square relative error，rRMSE)、平均誤差(Mean error，ME)來評價基于UAV航拍圖像重建的冠層結(jié)構(gòu)模型的精度。

2 結(jié)果與分析

2.1 冠層點(diǎn)云重建、分類及去噪

在大田完成一個目標(biāo)冠層的航拍需要約5 min。利用Metashape Professional軟件對各圖像序列處理后，得到目標(biāo)群體的冠層點(diǎn)云(圖4a、4d，成熟期以1行4株群體為例，下同)。經(jīng)SVM分類器對點(diǎn)云分類并基于MLESAC算法對地面落葉點(diǎn)云進(jìn)行平面擬合后，去除非植被點(diǎn)云，得到目標(biāo)群體的植被點(diǎn)云(圖4b、4e)。對所得到的植被點(diǎn)云采用空間單元格搜索k-近鄰點(diǎn)的方法確定噪聲點(diǎn)。去噪后，不同生育期的冠層點(diǎn)云均能很好地保留冠層信息(圖4c、4f)。

2.2 冠層結(jié)構(gòu)模型構(gòu)建

基于偽極點(diǎn)-Crust法進(jìn)行了冠層結(jié)構(gòu)模型的初步構(gòu)建，所得到的冠層模型植株間及葉片間存在大量異常面片(圖5a、 5d，黑色矩形框內(nèi)為指定區(qū)域放大后的葉片形態(tài))。通過對每個三角面片最大邊長進(jìn)行篩選，剔除異常面片，得到相對精確的冠層結(jié)構(gòu)模型(圖5b、5e)。最后對冠層表面進(jìn)行平滑處理，獲得最終的冠層結(jié)構(gòu)模型(圖5c、 5f)。冠層面片顏色由構(gòu)成每個三角面片頂點(diǎn)的顏色插值得到。

2.3 冠層結(jié)構(gòu)模型精度評估

為評估基于無人機(jī)航拍圖像重建的冠層結(jié)構(gòu)模型的精度，比較了株高、葉長、葉最大寬、葉面積的模型計(jì)算值與大田實(shí)測值(圖6，圖中n為樣本數(shù))。結(jié)果表明，苗期冠層結(jié)構(gòu)模型的計(jì)算株高、葉長、葉最大寬、葉面積與大田實(shí)測值的R2均在0.90以上；RMSE、rRMSE和ME的值均很小(表2)；成熟期重建的冠層結(jié)構(gòu)模型的株高、葉長、葉最大寬計(jì)算值與大田實(shí)測值間的R2均在0.90以上；RMSE、rRMSE和ME的值均很小(表2)；計(jì)算葉面積與參比葉面積的一致性稍差，相應(yīng)評估參數(shù)R2、RMSE、rRMSE、ME分別為0.76、87.7 cm2、0.187、-27.1 cm2(表2)。

表2 冠層結(jié)構(gòu)模型精度評價參數(shù)Tab.2 Parameters for evaluating accuracy of canopy structure models

3 討論

農(nóng)田冠層輻射利用效率是影響作物產(chǎn)量的重要因素。將冠層結(jié)構(gòu)模型與冠層光分布模型、光合模型相耦合，可實(shí)現(xiàn)作物新株型的精確評估及大田栽培模式的優(yōu)化，從而有利于作物高產(chǎn)。進(jìn)行此研究的前提是獲取精確的大田冠層結(jié)構(gòu)信息。采用立體視角或多視角圖像方法[15,31]、三維數(shù)字化方法[32]等雖然可以較精確地獲取大田作物冠層結(jié)構(gòu)信息，但多費(fèi)時、費(fèi)力，通常只能獲取較小群體的結(jié)構(gòu)信息，使得所構(gòu)建的冠層結(jié)構(gòu)模型不具有很好的代表性。故迫切需要開發(fā)高通量、高精度獲取大田作物冠層結(jié)構(gòu)的方法。

采用無人機(jī)平臺獲取大田作物信息，具有高通量的優(yōu)勢。由于以往所采用的無人機(jī)在超低空飛行時，會明顯擾動冠層，從而難以獲得高精度的航拍圖像序列。本研究采用超微小型無人機(jī)，其攜帶的RGB相機(jī)分辨率高，且在距冠層數(shù)米高度處獲取冠層圖像時對冠層的擾動較小。獲取半徑為5 m的圓形區(qū)域內(nèi)數(shù)百株植株的冠層圖像僅需數(shù)分鐘，且操作簡單，采用該方法可無損、快速地獲取苗期玉米群體圖像。在生長中后期，植株間的遮蔽逐漸加大，使得無法獲得植株中下部的結(jié)構(gòu)信息。為此，在航拍前去除了目標(biāo)植株周圍的植株。結(jié)果表明，以單株或一行多株為目標(biāo)進(jìn)行航拍，均可獲取相對完整的冠層結(jié)構(gòu)，器官尺度參數(shù)的提取精度較高。本研究采用的偽極點(diǎn)-Crust算法，通過添加少量的偽極點(diǎn)代替?zhèn)鹘y(tǒng)Crust算法中原始點(diǎn)云對應(yīng)的正負(fù)極點(diǎn)，可大幅度減少極點(diǎn)的數(shù)量。以1行4株玉米為例，原有的正負(fù)極點(diǎn)數(shù)為1 392 490個，采用偽極點(diǎn)-Crust算法后降為96個，計(jì)算效率較高。實(shí)現(xiàn)5 m×5 m范圍內(nèi)的苗期植株群體結(jié)構(gòu)模型的構(gòu)建僅需10 min左右，對于成熟期4株玉米小群體的模型構(gòu)建僅需1 min左右。未來可以在生育中后期探索更大目標(biāo)群體的航拍方法，以提高獲取作物生長中后期冠層結(jié)構(gòu)信息的效率。該方法還可應(yīng)用于植株葉片相互遮蔽不是很嚴(yán)重的其他作物，如高粱、生長中前期的大豆等，以及部分果樹、蔬菜等，對這些作物的應(yīng)用還需要系統(tǒng)地評估。

在玉米生長的中后期，葉片相互遮蔽嚴(yán)重，使得葉片間的噪聲點(diǎn)較多，采用搜索k-近鄰方法去除噪聲點(diǎn)的效果不夠理想。需要進(jìn)一步優(yōu)化去噪算法，以減少噪聲點(diǎn)對構(gòu)建精確冠層結(jié)構(gòu)模型的干擾。目前所構(gòu)建的冠層模型，還存在因點(diǎn)云缺失導(dǎo)致所構(gòu)建的葉片形態(tài)不完整的問題，需要采用算法予以改進(jìn)。基于構(gòu)建的精確冠層結(jié)構(gòu)模型，可進(jìn)一步計(jì)算葉片角度[33]并進(jìn)行大田冠層空間光分布的模擬分析，為評估不同作物株型、篩選高產(chǎn)作物品種提供支撐。

4 結(jié)束語

采用超微小型無人機(jī)平臺獲取了大田玉米冠層前期和后期的航拍圖像，實(shí)現(xiàn)了冠層點(diǎn)云的快速重建?；赟VM分類器提取了目標(biāo)冠層的三維點(diǎn)云，并實(shí)現(xiàn)了冠層噪聲點(diǎn)的去除。采用偽極點(diǎn)-Crust算法初步構(gòu)建了冠層結(jié)構(gòu)模型，通過剔除異常的三角面片并進(jìn)行移動平滑，得到了精度較高的玉米冠層結(jié)構(gòu)模型?；谒鶚?gòu)建模型提取的株高、葉長、最大葉寬及葉面積與大田測量值有較好的一致性。本研究所采用的圖像獲取方法、點(diǎn)云預(yù)處理流程及冠層結(jié)構(gòu)建模的算法效率較高，并可提取較精確的表型參數(shù)，為高效、精確獲取大田作物冠層結(jié)構(gòu)信息提供了新的途徑。