一種顧及體元鄰域關系的行道樹提取方法

徐嘉淼，劉如飛，魏曉東，柴永寧

(1.山東科技大學 測繪科學與工程學院，山東 青島 266590；2.青島秀山移動測量有限公司，山東 青島 266590)

0 引言

道路兩側的行道樹是城市建設管理中的重要組成部分，其在景觀功能、生態功能和意向功能等方面起著積極的作用[1]。因此，行道樹信息的自動化快速提取是智慧城市建設中不可或缺的環節。車載移動測量作為測繪領域內一種前沿的測量技術，可以在高速運動狀態下迅速獲取道路兩側地物的高精度三維信息，對城市行道樹信息化管理意義重大。

針對車載LiDAR點云行道樹目標自動化提取，國內外學者已經展開了廣泛的研究。一些方法對城市典型地物(行道樹、路面、建筑物等)統一進行分類與提取，從而間接提取出行道樹點云；另一些方法則單獨對行道樹點云進行提取。上述方法可以分為4類。①掃描線法。這類方法主要對樹干點云在每條掃描線上的分布特征進行分析，可以快速簡單地自動化提取出樹干點云[2]。但是雜亂無章的樹冠點云無法通過掃描線直接反映出幾何形態，不利于樹冠點云的進一步提取，適用性不強。②語義模型法。直接將離散的點云分割成獨立的點云塊，以行道樹的幾何特征為約束建立相應的語義模型，從點云塊中分類出行道樹目標[3-4]。這類方法可以直接對原始點云進行處理，但是在點云分割與聚類過程中容易受周圍鄰近地物影響，且計算量偏大。③點云投影法。將點云投影至平面內并建立二維格網，統計格網內點云的密度和高差等信息，從而對行道樹點云進行提??；或者生成規則的點云特征圖像，然后借助于圖像處理算法來提取出行道樹點云[5-6]。這類方法采用了降維的思想，減少了數據量，提高數據處理效率，但是在點云投影和生成規則圖像過程中，存在精度損失，且格網尺寸選擇對最終提取結果有影響。④體元法。將行道樹實體分割成一系列非空體元，既保留了行道樹的三維信息，也方便對其局部細節進行展示和分析，充分利用了行道樹的整體特征，減少數據冗余，提高計算效率[7-8]。Wu等[9]提出了一種基于分層體元點密度的行道樹提取方法，實現單株行道樹的完整提取，僅考慮單一的體元密度特征無法有效區分人造桿目標與行道樹樹干；Yang等[10]基于多尺度超體素進行區域聚類分割，根據先驗知識設定各類目標對象幾何特征分類閾值，實現多類目標層次化提取，實驗結果表明，該算法適用于復雜城市環境下的行道樹點云提取；Ordonez等[11]將點云進行三維空間體元簡化，對單個體元進行二維平面分析，利用桿狀地物的孤立性識別出目標數據，并分組、立體還原，最后根據支持向量機算法分離出行道樹目標；Kang等[12]建立規則體元索引，通過分析桿狀地物的孤立性與豎直連續性特征，采用自適應半徑的平面圓模型和垂直區域生長算法對桿狀地物體元進行檢測與分割，建立多類桿目標語義規則，進一步區分樹干、路燈與電線桿；羅海峰等[13]將體元與深度學習相結合，對地物點云建立八叉樹體元索引，歸一化分割重疊點云獲得獨立目標對象，基于多個方向生成目標對象的體元化全局特征，輸入深度置信網絡(deep belief networks，DBN)進行訓練，實現行道樹點云的提取。

上述算法無法對體元內點云進行細致的形態與屬性分析，并且沒有充分挖掘體元之間的鄰域關系。此外，當行道樹與其他地物距離太近，特別是路燈穿插在樹冠中時，極易造成誤分類。針對以上不足，本文提出了一種以體元為研究對象且顧及自身空間鄰域關系的行道樹點云自動提取算法，在城市部件采集、大規模園林普查等方面具有重要意義。

1 行道樹點云提取算法

行道樹廣泛種植于道路兩側或分車帶內，其自身結構復雜且常與其他地物混合，一定程度上造成了提取困難。以體元為研究對象且顧及自身空間鄰域關系的行道樹點云自動提取算法主要包括數據預處理、樹干提取、樹冠提取3個步驟，如圖1所示。

圖1 行道樹自動化提取流程圖

1.1 原始點云預處理

首先對原始點云數據進行雙邊濾波去噪和高程歸一化處理，然后采用文獻[14]中方法濾除地面點，保留非地面點，減少了數據量，提高數據處理效率。

車載移動測量記錄了激光點云的反射強度、三維坐標等信息。結合點云自身的三維坐標與體元尺寸參數實現三維虛擬規則格網的創建，方便對點云進行快速索引和鄰域搜索。假設體元在三軸方向上的尺寸分別為dx、dy、dz，存在點云P(Xp，Yp，Zp)落入對應的體元Voxel內，則Voxel的行號row、列號col、層號layer的計算方法如式(1)所示。

(1)

1.2 體元特征計算與鄰域分析

點云曲率是衡量地物表面平坦程度的重要特征。灌木和樹冠點云呈散亂團狀，樹干點云呈規則圓柱狀，導致前者的曲率大于后者，在體元內部表現為點云平均曲率的差異。點云回波反射強度是衡量地物表面反射率的重要特征。人造桿目標表面涂抹了特殊的反光材料且十分光滑，樹干表皮粗糙，導致前者的回波反射強度大于后者，在體元內部表現為點云平均回波反射強度的差異。

為提高樹干點云分類與識別的準確性，本文選擇了體元的曲率與回波反射強度特征，結合樹干幾何形態，充分挖掘體元空間鄰域關系，最終建立了表1所示的體元空間鄰域關系描述規則。

表1 體元空間鄰域分析

SVI和SVc表示種子體元的平均反射強度值和平均曲率值；VI和Vc表示非種子體元的平均反射強度值和平均曲率值；OL和RL表示第L層體元模板內點云擬合平面圓的圓心和半徑；Shape和Direction表示點云集合的形狀系數和主方向。

1.3 樹干與樹冠點云提取

1)體元水平鄰域一致性。行道樹周圍存在的人造桿目標(主要是路燈與交通標志牌)與灌木是行道樹點云提取過程中的主要干擾因素。人造桿與樹干點云之間的回波反射強度存在差異；同理，灌木與樹干點云之間的曲率存在差異。以體元的平均回波反射強度與平均曲率特征為約束，結合表1中的體元水平鄰域內特征描述，在同層格網內進行體元水平鄰域分析，可有效濾除人造桿目標與灌木體元。具體步驟如下。

①在距離地面1.2～1.4 m的那層格網內隨機選擇一個非空體元作為種子體元SV。

②若SV的平均回波反射強度與平均曲率特征同時滿足閾值條件，則將其標記屬性為Treek；反之，標記屬性為False，并重新選擇種子體元，直到滿足閾值條件為止。

③如圖2所示，以種子體元SV為中心創建5×5體元模板，搜索模板內其他非空體元V，計算二者的平均回波反射強度與平均曲率特征在二維歐式空間下的距離Distance，根據距離遠近判斷體元相似度。若滿足相似度閾值條件，則將體元V標記屬性Treek；反之，標記屬性為False。

圖2 體元水平鄰域一致性分析

④在同層格網內未被標記任何屬性的非空體元集合中重新選擇種子體元SV，并重復步驟②、步驟③，直到不存在未被標記屬性的非空體元為止。

2)體元豎直鄰域一致性。行道樹包括樹干和樹冠2個部分，可利用體元的豎直鄰域一致性來確定第一分枝點所在位置。樹干點云在豎直方向上連續分布，因此，在相鄰的2層網格中一定存在具有相同行列號的非空體元，將多層格網疊加處理，使得樹干在各層格網中的豎直鄰域關系得以體現。

對三維格網逐層進行體元水平鄰域分析，最終保留一定數量的5×5體元模板。搜索存在相同行列號體元的體元模板并進行聚類，則聚類體元內點云集合為疑似樹干點云。某層聚類體元內點云如圖3所示，利用隨機抽樣一致算法(random sample consensus，RANSAC)進行平面圓擬合，計算圓心O與半徑R，同時利用主成分分析法(principle components analysis，PCA)計算點云集合的形狀系數Shape與主方向Direction。利用表1中的體元豎直鄰域特征描述對疑似樹干點云作進一步判斷，判斷條件如式(2)所示。

圖3 體元豎直鄰域一致性分析

|OL-OL+1|

(2)

樹干呈豎直圓柱狀，即相鄰2層聚類體元內點云的擬合圓圓心坐標偏移量與半徑之差同時小于閾值，且點云集合的形狀系數為1(線狀)，主方向近似平行于Z軸，據此提取行道樹樹干點云。

3)分層點云投影面積計算。經過上述體元鄰域分析之后，可濾除大部分路燈與交通標志牌噪點，但是樹干與燈桿外形接近，且后者長期暴露在戶外，難免會發生損壞老化，其表面反射率不達標，無法通過體元平均回波反射強度和平均曲率特征進行濾除，造成個別路燈燈桿被誤分類為樹干。根據樹冠與路燈的形態差異，設計了一種分層點云投影面積判別法來進一步區分出行道樹與路燈，同時提取完整的樹冠點云，實現測區內行道樹的單體化輸出。

在原始三維格網索引下，從第一分枝點開始，自下而上逐層聚類搜索疑似樹冠點云。實際道路環境中很多樹木樹冠相連，且個別路燈穿插在樹冠中，在聚類搜索時需逐層設置約束半徑，防止樹冠點云過度生長。算法以樹干點云重心坐標為圓心，根據垂直于行駛方向的道路內側點云計算當前分層的約束半徑。

在數據采集過程中，道路外側的樹冠點云必然存在缺失，為提高分層點云投影面積計算的準確性，采用文獻[15]算法引入樹冠點云方位殘缺率，加權計算各分層點云投影面積值。

4)樹冠點云判別。行道樹樹冠點云的外形近似于“卵形”，若沿豎直方向等距切分，可分解為一個圓錐與若干圓臺模型。如圖4所示，樹冠點云被等距切分為19個冠層，按照上文所述方法計算各冠層點云投影面積，對19個面積值進行3次樣條插值，獲得光滑曲線。當冠層序號layer=7時，投影面積S取得最大值，投影面積曲線在宏觀上滿足先單調遞增后單調遞減的變化趨勢，于是建立冠層點云投影面積理論模型，如式(3)所示。

圖4 單個樹冠分層點云投影面積曲線

Si-1Si>Si+1

(3)

式中：Si表示第i個冠層點云投影面積。

隨機選取了10棵行道樹樣本，對其樹冠點云進行分層，分別計算各冠層內點云投影面積S，并繪制投影面積曲線(圖5)。不難發現，雖然不同行道樹的樹冠體積存在較大差異，但自下而上各冠層點云投影面積總體上呈先遞增后遞減的變化趨勢，仍然滿足式(3)的投影面積理論模型。

圖5 10棵行道樹樹冠分層點云投影面積曲線

如圖6所示，路燈貫穿樹冠，且下部燈桿被誤分類為樹干。這類路燈的上部聚類點云包含樹冠與燈頭，分層計算聚類點云投影面積并繪制投影面積曲線(圖6)。從圖6中可知，在第10冠層(layer=10)處點云的投影面積S取得最大值；當layer<10時，冠層點云的S值單調遞增，與理論模型相符；當layer>10時，S值理論上單調遞減，但是受上部燈桿與路燈燈臂影響，導致局部冠層點云(layer=13、14、15)投影面積反而增大；當layer>20時，僅存在燈頭點云，多個相鄰冠層點云的S值近似相等，曲線趨于平緩。故存在一定數量的分層點云投影面積不符合理論模型。

圖6 樹冠與路燈分層點云投影面積曲線

結合樹冠分層點云投影面積理論模型來對疑似樹冠點云進行判別，進一步剔除被誤分類為樹干的路燈。為量化樹冠點云的判別指標，引入道格拉斯-普克算法[16]提取投影面積曲線輪廓特征點，根據特征點計算樹冠判別值δ。

道格拉斯-普克算法作為一種全局的曲線化簡算法，能在一定程度上保持化簡后曲線的整體形態。對圖4樹冠分層點云投影面積曲線使用道格拉斯-普克算法進行化簡，結果如圖7(a)所示，最終保留7個特征點，且曲線在特征點3處取得最大值。依次計算相鄰特征點i和i+1之間的連續斜率ki，同時記錄2點之間所包含的若干冠層(如特征點4、5之間包含了第10、11、12冠層，連線斜率為k4)。在樹冠分層點云投影面積理論模型約束下，特征點1、2、3位于單調遞增區間內，故相鄰特征點連線斜率ki>0(i=1，2)；特征點4、5、6、7位于單調遞減區間內，故相鄰特征點連線斜率ki<0(i=3，4，5，6)。且相鄰多個冠層點云投影面積變化明顯，斜率絕對值不小于閾值Kmin。

同理，圖7(b)為圖6中投影面積曲線的化簡結果，保留了9個特征點，曲線在特征點5處取得最大值。位于單調遞減區間內的特征點6、7的連線斜率k6>0；特征點8、9的連線過于平緩，其斜率k8的絕對值小于閾值Kmin。顯然，上述2組相鄰特征點內的冠層點云投影面積與理論模型相悖。

圖7 道格拉斯-普克算法曲線簡化

一次迭代完成之后，統計滿足斜率條件的相鄰特征點之間的冠層數量，用于計算樹冠判別值δ(滿足理論模型的冠層數量與冠層總數之比)。若δ>80%，則認為當前聚類點云為樹冠，否則作為偽樹冠剔除。

2 實驗

2.1 實驗數據與參數設置

為驗證算法的有效性，采用車載移動測量系統所采集的2組實際路段點云進行實驗，數據點位精度優于3 cm。如圖8(a)所示，道路Ⅰ全長約0.547 km，為雙向4車道，地形平坦，道路兩側的行道樹主要為白蠟樹，數據共有10 811 215個點。如圖8(b)所示，道路Ⅱ全長0.835 km，為雙向2車道，地形起伏較大，道路兩側的行道樹主要為法國梧桐，數據共有16 280 769個點。2條道路均存在低矮綠化，主要的桿目標有行道樹、路燈和交通標志牌。

圖8 實驗數據

根據城市行道樹種植規范，道路兩側的行道樹規格接近，且樹干位于所在區域地面1.3 m處，為林業部門定義胸徑的高度[17]。經過實地量測，道路Ⅰ的白蠟樹樹干平均胸徑為25 cm，道路Ⅱ的法國梧桐樹干平均胸徑為38 cm。以行道樹樹干的平均胸徑作為體元的劃分尺度，樹干點云提取過程中需人工設定的參數見表2。對于樹冠點云的提取，分層點云投影面積曲線，經過道格拉斯-普克算法簡化后，保留若干特征點，需人工設定相鄰特征點連線斜率絕對值的最小閾值Kmin。若Kmin設置過大，則無法判別樹冠點云；設置過小，則無法有效去除被誤分類為樹干的人造桿目標。經過大量樣本實驗，確定Kmin的取值為0.6。

2.2 結果與討論

本文利用C++實現上述算法，設置好參數對原始點云進實驗，2條道路上的行道樹點云提取結果分別如圖9和圖10所示。圖中樹干點云被統一賦成紅色，樹冠點云按照行道樹ID號隨機賦色。道路兩側的行道樹被準確提取，不受周圍人造桿目標和灌木的影響。

圖9 道路Ⅰ提取結果

圖10 道路Ⅱ提取結果

圖11(a)和圖11(b)分別是圖9中區域A和區域B的局部放大圖，將提取得到的行道樹目標與原始點云疊加顯示。路燈、監控探頭、交通標志牌和行道樹這4類桿目標相互交錯且周圍存在大量低矮灌木，一定程度上對行道樹目標的提取造成了困難。算法依據體元的平均回波反射強度與平均曲率特征進行空間鄰域分析，可以準確提取樹干點云，不受鄰近地物噪點的干擾。

圖11 行道樹提取結果與原始點云疊加

圖12為圖10中區域A的局部放大圖，樹木之間距離較近，部分樹冠相連。算法以道路內側點云為參考，逐層計算冠幅半徑Rmin，然后以Rmin為聚類生長約束條件，逐層提取樹冠點云。根據最后提取結果可以看出，樹冠分割效果較好。

圖12 相連行道樹提取結果

樹干是行道樹的重要組成部分，本文算法首先提取樹干，然后從樹干頂端向上聚類獲取樹冠點云。圖13(a)為圖10中區域B的原始點云，路燈貫穿樹冠。圖13(b)為樹干提取結果，其中，路燈的下部燈桿因表面老化導致反射率不達標，從而被誤分類為樹干。算法建立了樹冠分層點云投影面積理論模型，搜索疑似樹冠點云，并逐層計算投影面積，在理論模型約束下計算樹冠判別閾值。最終提取結果如圖13(c)所示，有效濾除了被誤分類為樹干的燈桿。

圖13 樹冠點云判別

但2組實驗數據中仍然存在一些錯誤提取及遺漏的情況。個別樹木距離掃描儀較遠且受到車輛遮擋，樹干點云缺失，不滿足桿狀地物特征，無法通過體元鄰域分析進行識別。此外，個別樹冠的枝葉過于稀疏，其外形不符合常規樹冠的“卵形”，無法使用冠層投影面積變化理論模型對其進行判別。可見，點云數據的完整程度與密度直接影響整體的提取效果。對原始數據和實驗結果進行了統計分析，根據式(4)評價行道樹點云提取的準確率與召回率。

(4)

式中：TP為正確提取的數量；FN未被提取的數量；FP為錯誤提取的數量。

統計結果如表3所示，算法對2組數據中行道樹目標提取的準確率分別達到94.77%和92.61%，召回率分別達到91.78%和89.13%。將本文算法和同樣采用體元法的文獻[9]進行了對比。編程實現文獻[9]算法，再次對2組實驗數據進行處理，對比分析正確率與召回率。如表4所示，本文方法對行道樹提取的正確率與召回率均有提高。由于部分人造桿目標穿插在樹冠中，文獻[9]直接對非空體元聚類會將其誤分類為樹木，降低正確率。行道樹底部存在大量低矮灌木，與樹干距離較近的上部灌木噪點會破壞樹干的形態，造成漏提，降低召回率。

表3 行道樹提取結果定量分析

表4 方法對比

3 結束語

針對車載激光掃描數據，本文提出了一種顧及體元鄰域關系的行道樹提取方法。首先建立規則格網索引，然后以樹干幾何特征為基礎，結合點云的回波反射強度與曲率屬性，充分挖掘體元的空間鄰域關系，成功提取出樹干點云。在此基礎上，建立樹冠分層點云投影面積理論模型，聚類搜索疑似樹冠點云并分層計算投影面積，將計算結果與理論模型進行比較，通過設定分類閾值來識別樹冠，同時濾除被誤分類為樹干的燈桿，最終實現行道樹點云的自動化提取。實驗結果表明，本文算法可以避免鄰近地物噪點影響，準確提取復雜城市道路環境下的行道樹點云，提取結果的正確率和召回率分別達到了93.41%和90.12%，具有較強的適應性與穩健性。在后續的研究中，可嘗試引入新的體元特征，以增強行道樹提取的穩健性。此外，也可以根據樹冠分層點云投影面積進行精細化建模，用于行道樹種類的識別。