利用單次無人機影像的果樹精準識別方法※

李福根，段玉林，史 云，吳文斌，黃 平

（1. 中國農業科學院農業資源與農業區劃研究所/農業農村部農業信息技術重點實驗室，北京100081；2. 中國科學院遙感與數字地球研究所，北京100101；3. 四川省農業科學院遙感應用研究所，成都610066）

0 引言

果樹的種植和管理是現代農業發展的重要組成部分，據聯合國糧農組織統計，全世界果樹種植面積約占總耕地面積的3%～4%［1］。準確識別果樹并統計果樹數量對監測果樹長勢、果園產量估算以及種植管理至關重要。相較于傳統的人工果樹統計費時費力且主觀性強的特點，遙感技術以其經濟性強、持續性好、可靠性高的特點，可以提供大面積、長時期必要的果樹信息數據，成為保證果樹數量可持續統計的重要方式［2］。

利用遙感技術對果樹進行識別并統計的理論基礎是提取果樹在影像中的光譜特征。Wulder等［3-5］對高分辨率衛星影像進行了系列研究，假設波段組合后的影像局部峰值即為樹的位置，利用局部最大濾波器即可識別樹的位置。Gebreslasie等［6］研究發現，利用局部峰值識別并定位單株樹木的方法其實質是基于高斯濾波器平滑衛星影像以消除噪聲，同時還提出采用半變異函數來確定局部峰值檢測窗口的大小。

然而，考慮到果樹樹冠的大小和形狀，要對果樹進行精準識別并計數需要較高空間分辨率的遙感影像，這使得空間分辨率較低的衛星影像在此方面的應用得到很大限制。隨著無人機遙感技術的發展，使得利用無人機獲取果園尺度高空間、高時間和多光譜分辨率影像成為可能；同時，無人機影像處理過程中利用SFM方法［7-8］產生的密集點云可以生成質量較高的數字表面模型（Digital Surface Model，DSM）影像［9-10］，也為利用遙感技術進行植被識別與計數提供了新的可靠數據源。Wan等［11］提出了一種利用無人機可見光和多光譜影像對油菜花的開花數量進行統計的方法，但該方法需要分類數據做配合提供油菜花叢的位置信息。Wu等［12］將無人機影像與深度學習結合對小麥抽穗數量進行識別，得到了較好識別精度，但該方法使用的前提是需要人工事先勾畫小麥生長區，而且應用區域受限于訓練樣本區域。Harris公司［13］發布了普適性較強的Count Crops Tools工具，該工具僅需輸入作物的最小和最大直徑就能對作物進行自動識別和計數，考慮到果樹之間通常有較明顯的間隙，可以單獨區分，因此該方法被認為在果樹識別與計數方面有著廣闊的應用前景。然而，由于遙感影像通常只能表達作物冠層頂部信息，該方法在果樹識別中不可避免地會將果園中的雜草和果園旁的其他樹木誤識別為果樹，從而影響統計精度，使其在實際果園管理中的應用受到很大限制。

目前還沒有利用DSM影像進行作物識別的研究，但已有部分研究利用無人機DSM求取作物高度。Ziliani等［14］利用無耕種時期的無人機數字表面模型作為基底數據，利用獲取的植被生長期無人機數字表面模型與基底數據做差求取植被高度，但該方法需要在無耕種時利用無人機獲取DSM，無法利用單次無人機數據計算植被高度，對于已種植果園來說是不可能實現的。

文章通過對目前研究現狀的分析，首次將無人機密集點云產生的DSM影像作為輔助數據加入到果樹識別與計數研究中，結合利用果樹冠層直徑對果樹的識別算法，可以在不用人為對影像研究區進行勾畫的情況下利用單次無人機拍攝的一組數據（DSM和正射影像）實現對果樹的識別與計數。

1 研究區與數據

1.1 研究區概況

該文選擇位于美國加利福尼亞州弗雷斯諾縣里德利市郊區的一個果園，中心坐標為北緯36°34′33″，西經119°26′06″。該果園旁邊有高大喬木和低矮灌木，果園中間也有成簇的雜草叢，是進行果樹精準識別理想的研究區域（圖1）。

1.2 數據獲取及預處理

該研究所用數據由旋翼無人機搭載MicaSense公司的RedEdge相機采集。拍攝時間為2017年10月19日，無人機飛行高度30 m。

研究選用Agisoft PhotoScan Profession軟件來處理無人機影像，處理過程包括特征點提取、特征點匹配、密集點云生成、DSM和正射影像生成。最終獲得的影像空間分辨率為0.01 m。

圖1 研究區果園無人機影像Fig.1 Image of unmanned aerial vehicle in the orchard of study area

2 果樹識別及計數方法

該研究目的是利用高分辨率無人機遙感影像對果樹進行識別。目前，無人機影像主要是普通數碼相機提供的可見光影像和RedEdge相機提供的多光譜影像，這2種數據在生成正射影像的同時都會依據影像處理過程中產生的密集點云生成DSM影像。因此，將分別對這2類數據進行研究，利用植被指數對植被區域進行提取，將非植被區域設置為背景，消除其對果樹識別時產生的影響。

將消除非植被區域的植被指數影像作為輸入影像，果樹直徑范圍作為輸入參數，利用ENVI5.5軟件的Count Crops工具，獲得初步識別的果樹影像。由于該方法會將果樹直徑范圍內的其他樹種和草叢誤識別為果樹，因此需要對產生的結果進行二次處理，從而提高果樹的識別精度。

研究將初步識別的果樹影像與無人機處理后同步生成的DSM影像疊加，以每一棵識別到的果樹范圍作為一個搜索窗口，搜索其在DSM影像中的相應位置，然后向外擴展一個像素形成新的計算窗口。假設搜索窗口中最大值為植被最高點的海拔高度，計算窗口中最小值為植被臨近地表的海拔高度，求取兩個高度的差值即可計算出植被的高度。利用果樹的高度范圍和求取的植被高度數據，對初步識別出的果樹進行再次識別，最終確定果樹的位置和數量。具體流程如圖2所示。

圖2 果樹精準識別流程Fig.2 Flowchart of fruit trees accurate detection

2.1 植被指數計算

利用圖像光譜信息進行目標識別的第一步是利用圖像光譜特征減少信息量［15］。由于該研究是針對果樹的識別，因此選擇能增強植被特性的植被指數來對研究區內的植被進行提取，同時去除非植被區域對研究結果的干擾，形成植被影像，提高植被識別精度。

目前，基于可見光影像和多光譜影像的植被指數有40多種［16］，它們的提出都對應不同的作物性質反演，例如植被健康程度、植被葉綠素含量和植被過火面積等。為了使該方法有更好的收斂性，后續流程能平穩運行，提高方法的普適性，該研究計算了5種基于可見光（2種）和多光譜影像（3種）的歸一化植被指數（詳見表1），并對植被指數在研究區植被的識別精度進行了分析討論。

表1 5種植被指數及其表達形式Table 1 Five vegetation indices and their band-specific formulations and associated principal reference

續表1

2.2 基于果樹樹冠直徑的果樹初步識別

ENVI5.5軟件以及配套的IDL87編程軟件提供了一種基于作物冠層直徑范圍的作物識別計數方法（函數）［13］。該方法以遙感影像的像元分辨率為計量依據，以一定容差范圍內的聚類像元為對象對作物進行識別并統計。該文基于該方法以及編程函數對植被指數影像中的果園進行初步識別，編程函數表達形式如下：

該函數中，ENVIAgCropCount為調用的函數名；outCrops代表函數的輸出對象，該對象包括兩部分，一部分為作物位置影像，另一部分為包括作物中心坐標以及冠層直徑的數據庫；Input代表輸入的柵格影像，該研究采用閾值分割后的植被指數影像；Minimum_Crop_Diameter為果樹冠層的最小直徑（單位：m）；Maximum_Crop_Diameter為果樹冠層的最大直徑（單位：m）；GAUSSIAN_FACTOR為高斯平滑因子；/INCLUDE_EDGES為可選擇輸入，代表作物識別時是否包括影像邊緣的作物；INTENSITY_THRESHOLD為作物強度閾值，用于影像像元的聚類分析；NUMBER_OF_INCREMENTS為搜索步長；OUTPUT_NCROPS為輸出的作物數量，并將結果儲存到outputNumCrops變量中；PERCENT_OVERLAP為作物重疊度。

2.3 利用初步識別的果樹位置求取果樹高度

該研究利用IDL87的ENVIAgCropCount函數生成了作物位置影像和包括作物中心坐標以及冠層直徑的數據庫，并將該位置信息影像疊加到無人機同步生成的DSM影像中，根據DSM影像分辨率和輸出數據庫中的果樹中心坐標以及冠層直徑信息，可以得到研究區果樹海拔數據集。

式（1）～（2）中，D（u）代表果樹的像素集合，d為識別到的果樹的直徑，PPI為DSM影像的空間分辨率，R為識別到的果樹的半徑所對應的像素數量，o（x，y）為識別到的果樹中心坐標。

根據果樹的鄰近空間關系［24-25］，將果樹的半徑向外擴展一個像素，形成新集合。該集合表示識別到的單棵果樹和其最鄰近的地表的海拔數據集。

此時，果樹高度即可表示為：

式（4）中，H代表果樹高度，hmax代表數據集D（u）中的最大值，hmin代表數據集D′（u）的最小值。

2.4 利用果樹高度信息對果樹進行再識別

在縱向空間中，果樹的高度通常介于高大喬木與低矮灌木之間，形成其獨特的縱向特征。依據此特性，以果樹高度范圍為閾值，利用DSM和位置信息求出的高度數據進行二次識別，剔除不符合果樹高度范圍的對象，生成全新的果樹位置影像并統計果樹數量。

式（5）中，Hmin，Hmax分別為果樹高度的最小值和最大值。

3 結果與分析

為了評價該文提出的方法對果樹進行識別的精度，該研究利用無人機拍攝生成的可見光、多光譜和DSM影像進行了果樹識別實驗。實驗包括植被指數篩選、基于果樹冠層直徑的果樹初識別、基于果樹初識別結果的果樹高度計算和果樹二次識別等分析，并基于實驗誤差對方法可行性進行了討論。

3.1 植被指數選擇

該文計算了5種植被指數，并依據植被指數的特性來提取植被區域，從而消除非植被區域對果樹識別產生的影響。為了篩選出識別植被最佳的植被指數，通過目視解譯和定量評價2個角度分別對5種植被指數的結果進行了比較。

圖3為從研究區選取的示例區域原始圖與5種植被指數在該區域的植被識別結果。結果發現，基于多光譜影像的MRENDVI對植被的識別結果最好，植被區域識別較為完整，幾乎沒有誤識別區域，可以作為植被識別的最佳特征指數。對于可見光影像，GLI也有很高的識別精度，可以作為可見光影像對植被識別的特征指數。但GLI會將部分行壟間的非植被區域誤識別為植被，這可能是因為地物散射造成的可見光影像噪聲引起的。而NDI、NDVI和GNDVI這3個指數的識別結果都摻雜著大量的植被陰影，無法滿足研究對植被識別的需求。

圖3 原始影像與5種植被指數提取植被區域結果：（a）原始影像；（b）NDI提取的植被區域；（c）GLI提取的植被區域；（d）NDVI提取的植被區域；（e）GNDVI提取的植被區域；（f）MRENDVI提取的植被區域Fig.3 Original Image and vegetation area extracted from five vegetation indices

目視解譯的結果在一定程度上存在主觀性。為了更加客觀地比較5種植被指數對植被識別的優劣，該研究選用Matusita距離法［15，26］來定量化描述植被區域和非植被區域植被指數頻率直方圖的差異性，進而比較植被區域和非植被區域的非相似性，從而評價5種植被指數的識別精度。計算公式如下：

式（6）中，DistM（v，s）代表Matusita距離，v（x）代表植被區域頻率直方圖，s（x）代表非植被區域頻率直方圖。

表2為5種植被指數對植被區域進行識別后，植被區域和非植被區域頻率直方圖的Matusita距離值。研究發現該結果與人工目視解譯結果高度吻合。MRENDVI的Matusita距離值最大，為0.43，說明識別到的植被區域和非植被區域有高度的非相似性，可以作為植被識別的一種重要指數。GLI的Matusita距離為0.40，說明GLI也具有植被識別能力，但效果不如MRENDVI，這是因為GLI本質是一種將可見光波段中的綠波段增強的指數，該指數雖然很大程度上突出綠葉在影像識別中的貢獻，但不可避免地會將影像中其他綠波段反射率較高的像素點識別為植被。NDI、NDVI和GNDVI的Matusita距離值相似，分別為0.29、0.23和0.22，說明這3個指數識別的植被影像和非植被影像的非相似性不高，在該研究中不適合作為植被識別的植被指數來使用。

表2 植被區域與非植被區域頻率直方圖的Matusita距離Table 2 Matusita distance of the histograms of vegetation area and non-vegetation area

根據上述比較分析，該研究最終選擇GLI為可見光影像的植被識別指數，MRENDVI為多光譜影像的植被識別指數。為了對該文提出的研究方法進行更加科學地分析，研究選用了MRENDVI提取的植被區域進行后續操作。

3.2 果樹初識別及高度

研究將已知的果樹冠層直徑范圍（1.5～2.5 m）輸入到ENVIAgCropCount函數中，得到初識別的果樹影像。

如圖4所示，果園中果樹識別較為完整可靠。對識別到的果樹對象計數，共識別到3 191棵果樹。該方法仍將一些雜草和果園周圍的其他樹木識別為果樹。同時，果園中的果樹有部分漏識別，主要原因是部分果樹已經枯萎（圖5（a）），無明顯冠層，研究區內這種情況的果樹共有86棵，該研究中不再把此類果樹加入到實際果樹計數中。還有部分漏識別是因為冠層長勢較差，直徑太小，未在輸入的已知果樹直徑范圍內（圖5（b）），該情景漏識別的果樹共有7棵（占總果樹數量的0.2%），這也是ENVIAgCropCount函數存在的局限性之一。對于該情景還需要進一步論證，尋求擴大果樹直徑范圍和識別精度之間的平衡。

圖4 果樹初識別結果以及局部特征Fig.4 Fruit trees preliminary detection results and certain area

排除掉已枯萎的果樹，通過目視解譯，研究區內共有果樹3 034棵，識別精度為94.8%。

圖5 兩種果樹漏識別示例：（a）枯萎果樹漏識別；（b）未在已知冠層直徑范圍內果樹漏識別Fig.5 Two kinds of examples of non-detection in fruit trees

為了進一步提高果樹識別精度，該文利用無人機同步生成的DSM影像和果樹高度求取方法計算了研究區內識別出的果樹高度，如圖6所示。研究發現，果園內果樹的高度普遍在1.8～2.5 m區間內，這與該果園果樹實際高度吻合。果園西側誤識別為果樹的喬木高度普遍大于2.8 m，果園東側和南側誤識別為果樹的灌木的高度普遍高于0.5 m且低于1 m，果園中果樹間誤識別為果樹的草叢高度則接近0 m，它們的高度與果樹高度形成了鮮明對比。

圖6 果樹高度示意圖Fig.6 Height of fruit trees

3.3 果樹再識別

研究利用計算出的果樹高度和果樹再識別方法，將1.8～2.5 m作為果樹的高度范圍，對果樹初識別的結果進行二次識別，并根據識別后對象的數量進行計數統計，識別結果如圖7所示。從圖7中發現，識別到的果樹均在果園范圍，且均勻整齊排列，這與現實情況相符。果園以外其他植被的誤識別和果園內部果樹間草叢誤識別基本消除，這說明該方法對于提高果樹識別精度有較明顯的效果。但對果園中的果樹有部分漏識別，這主要是因為ENVIAgCropCount函數進行初識別時對部分果樹未識別（詳見3.2節），說明該方法對漏識別果樹沒有補救措施，需要在日后研究中進一步改進。通過果樹高度對果樹再次識別后，共識別3 027棵果樹，與實際3 034棵果樹相比，識別精度為99.8%。相較于初識別結果，精度提高了5%。由于研究區中混雜非果園植被區域面積有限，造成該方法識別精度提高受到限制，需要在后續研究中對該方法進行改進，進一步提高該方法對果樹的識別精度和普適性。

圖7 果樹識別結果Fig.7 Fruit trees detection results

4 結論

為了能更加方便、準確地在宏觀尺度對果樹進行識別，該文基于現階段作物識別和無人機應用的研究現狀，首次將無人機影像拍攝的正射影像和DSM相結合，提出了一種利用單次無人機影像的果樹精準識別方法。該方法通過對無人機生成的正射影像進行植被指數計算，提取研究區內的植被區域；根據植被區域和果樹冠層直徑范圍對果樹進行初步識別；基于果樹位置和無人機生成的DSM影像計算果樹高度；最后依據果樹高度范圍實現對果樹的精準識別。

該方法在美國加利福尼亞州弗雷斯諾縣里德利市郊區的一個果園中進行了應用測試。研究發現，利用該方法對果樹進行精準識別，精度達到了99.8%，較果樹初識別結果精度提高了5%。識別誤差來源主要是果樹初識別過程中漏識別的果樹。該方法基本消除了果樹初識別過程中，果園周圍植被和果園內部草叢被誤識別為果樹的誤差，但對于初識別過程中漏識別的果樹還無法進行補救，需要進一步提高該方法的識別精度。

該文提出了一種利用果樹基本表型參量實現果樹精準識別的方法，具有較高的普適性和可推廣性。該方法的提出促進了利用無人機遙感影像對植被表型參量的反演從二維尺度走向了三維空間。隨著人工智能技術的發展，在未來可以將該方法與人工智能技術相結合，實現更加快速、準確的果樹識別。