基于機載LiDAR的高郁閉度華北落葉松林單木識別

陳思宇，劉憲釗，王懿祥，梁 丹

（1.浙江農林大學 環(huán)境與資源學院，浙江 杭州 311300；2.浙江農林大學 省部共建亞熱帶森林培育國家重點實驗室，浙江 杭州 311300；3.浙江農林大學 浙江省森林生態(tài)系統(tǒng)碳循環(huán)與固碳減排重點實驗室，浙江 杭州 311300；4.中國林業(yè)科學研究院 資源信息研究所 國家林業(yè)和草原局森林經營與生長模擬重點實驗室，北京 100091）

華北落葉松Larix principis-rupprechtii是中國華北地區(qū)第一大樹種，具有較高的經濟和生態(tài)價值，探討華北落葉松單木參數(shù)的提取、林分蓄積量的估算具有重要意義[1]。高精度的單木識別是獲取單木和森林結構參數(shù)的重要前提，單木識別包括單木位置識別和樹冠分割[2-4]。相比高分辨率衛(wèi)星影像和光學影像，機載LiDAR技術不受區(qū)域信號影響，分辨率較高，是一種主動遙感技術，能獲取林木的水平和垂直結構信息[5-6]。利用機載LiDAR研究單木識別的方法主要有2類：基于冠層高度模型(CHM)識別方法和基于歸一化后的點云識別方法[3]。其中基于CHM的單木位置識別方法即在灰度圖像中探測單木樹頂，主要實現(xiàn)算法包括局部最大值法和多尺度分割法。局部最大值法假設樹冠頂點高于樹冠邊緣，通過設置固定尺寸[7-10]或動態(tài)尺寸[11-12]的搜索窗口從CHM中識別樹冠頂點特征，單木位置探測率可達85%。多尺度分割法以像元為單位，根據相似特征逐漸外擴分割，關鍵在于最優(yōu)分割尺度的選取[13-17]?；跉w一化后的點云單木位置識別方法以點云的結構信息為基礎，主流實現(xiàn)算法包括區(qū)域生長法和均值漂移法。區(qū)域生長法通過選擇合適的種子點，定義距離生長準則，判別不同屬類[18-20]，選擇的種子點和制定的距離生長規(guī)則不同，取得的單木識別精度也不同。均值漂移法是一種基于密度梯度上升的非參數(shù)方法，通過迭代運算求得同簇的數(shù)據點[21]，關鍵在于漂移均值向量的計算。單木識別時可利用林分樹冠的距離和多維特征空間為向量來完成[22-23]，其在中低郁閉度林分中的單木識別精度較高，但在高郁閉度林分條件下，由于樹冠重疊度增加，漂移均值向量計算會有偏差，識別結果不盡合理。為提高識別精度，結合樹冠冠型特點，本研究提出基于機載LiDAR點云空間特征的高斯核函數(shù)改進的均值漂移單木位置識別方法(MSP)。為驗證方法的有效性，對比了基于點云空間特征的區(qū)域生長圖像分割方法(RGP)、基于CHM的局部最大值單木探測方法(LMC)和基于CHM的多尺度分割方法(MSC)等3種目前應用較為廣泛的單木位置識別方法。最后，在單木位置識別結果上，進行單木樹冠的提取。

1 研究區(qū)概況

忻州市五寨縣羊圈溝林場 (38°44′～39°01′N，111°48′～112°03′E)位于山西省中北部，屬亞寒帶大陸性氣候，季節(jié)變化明顯，海拔為1 350～2 700 m，年均氣溫為5 ℃，年降水量為453.9 mm；華北落葉松-云杉Picea asperata天然次生林與華北落葉松人工林帶狀混交。本研究在華北落葉松人工林帶內設置8個20 m×20 m 的正方形樣地，植株密度為 1 875～2 925 株·hm-2，郁閉度為 0.8～1.0，是典型高郁閉度森林。

2 數(shù)據采集與研究方法

2.1 數(shù)據采集

標記并編號8個樣地內胸徑大于5 cm 的華北落葉松，調查記錄每木冠幅、樹高和胸徑；利用實時動態(tài)定位和全站型電子測距儀獲取每木的位置坐標。樣地林分基本情況如表1。

表1 樣地林分基本特征Table 1 Investigation results of basic characteristics of forest stand in sample plot

2018年6月2日，使用GV 1500六旋翼無人機搭載RIGEL VUX-1激光掃描儀獲取機載LiDAR點云數(shù)據；設置飛行高度為200 m，飛行速度為6 m·s-1，離散點云平均密度為70個·m-2，坐標系為WGS_1984_UTM_Zone_49N。

在同一天，使用GV 120六旋翼無人機搭載傳感器SONY Alpha7 RⅡ采集可見光遙感影像。通過Limapper V 2.0軟件生成數(shù)字正射影像作為輔助數(shù)據，用于目視解譯與可視化數(shù)據展示。

2.2 數(shù)據預處理

首先，對原始點云數(shù)據進行2 m緩沖帶裁剪，獲取8個樣地。然后，對樣地點云數(shù)據進行以下處理：①采用距離閾值法移除噪聲點。②采用漸進加密三角網濾波算法提取地面點云。為提高后期生成數(shù)字高程模型(DEM)的精度，設置迭代距離與迭代角度為1.4 m和8°，以便分離地面點云和非地面點云。③對去噪和分類后的點云數(shù)據進行克里金插值[23]。通過定義該點鄰域中采樣點的距離加權平均計算像元值。對地面點云和非地面點云分別插值，得到DEM和數(shù)字地表模型(DSM)，兩者相減獲得冠層高度模型(CHM) (圖1)。圖1C中，高程值較低的條帶為華北落葉松人工林帶，較高的為華北落葉松-云杉天然次生林。④基于LiDAR360 V 3.2軟件，借助DEM對去噪后的點云進行高度信息歸一化，避免地形起伏對點云高程的影響。

圖1 研究數(shù)據Figure 1 Study data

2.3 單木位置識別

2.3.1 基于點云空間特征的單木位置識別 均值漂移算法是一種基于非參數(shù)概率密度梯度估計的迭代聚類過程[24-25]。單木點云分布具有與高斯核函數(shù)相似的幾何形狀特征，可用高斯核函數(shù)改進的均值偏移算法對點云數(shù)據進行聚類，進而獲取單木位置。設置迭代運算的收斂閾值(ε)，基于點云空間特征的改進均值漂移單木位置識別方法，算法如下：(1) 在點云數(shù)據點中隨機選擇1個點作為初始中心點，找出離中心點距離在帶寬h之內的所有點，記做集合M，認為這些點屬于同一類簇；(2)選用點云三維空間特征，計算中心點到集合M中每個數(shù)據點的向量之和，得到偏移向量Mh(x)；(3)令中心點沿偏移向量方向移動，移動距離是偏移向量Mh(x)的模；(4)重復步驟(2)～(3)，直到偏移向量小于ε，記住此時的中心點；(5)重復步驟(1)～(4)，直到所有數(shù)據點都被歸類；(6)提取單木位置，每一個類簇的中心點即為單木位置。

2.3.2 基于CHM的單木位置識別 局部最大值法基于單木樹冠頂部高程大于樹冠邊緣高程的假設，設置搜索窗口進行局部最大值探測。本研究將搜索窗口的起始尺寸設置為林分的平均冠幅大小，對窗口尺寸按增加或減少1個像素進行動態(tài)調整，搜索出樹冠的局部最大值。基于CHM的圖像分割法通過圖像分割實現(xiàn)單木分割，根據區(qū)域中心確定單木位置。標記控制分水嶺分割法是最常用的分割方法，但在本研究中分割精度不理想，因此選用面向對象的多尺度分割。利用eCognition Developer V 8.9軟件的Estimation of Scale Parameter 2 (ESP2)工具包得到合適的分割尺度。本研究中，分割尺度為 66，形狀參數(shù)為0.8，緊致度參數(shù)為0.9，樣地分割效果最好。

2.4 單木樹冠提取

在單木位置識別的基礎上，劃分出同一類點云的投影范圍即為該株華北落葉松的樹冠范圍?？紤]到冠層提取過程中，林分樹冠間隙部分容易產生過分割現(xiàn)象，需要對其結果進行區(qū)域合并[26]，本研究將形狀異質性作為類別之間的合并準則。選擇eCognition軟件對點云初始分類結果進行區(qū)域合并，分割尺度為58，形狀參數(shù)和緊致度參數(shù)分別為0.8和0.8。

2.5 精度評價

2.5.1 單木位置識別精度評價 根據衛(wèi)星定位測量獲取坐標的定位精度，考慮樹頂與樹根的位置偏差，將華北落葉松位置誤差預計為1 m；以真值坐標點為圓心，創(chuàng)建半徑為1 m的緩沖區(qū)。單木位置識別結果落入此緩沖區(qū)內，則認為正確。測算正確識別的單木總數(shù)(NC)、實地調查單木實際總數(shù)(NF)、識別出的單木總數(shù)(ND)、重復識別的單木總數(shù)(NR)等數(shù)據。計算精度評價指標，并使用單因素方差分析檢驗不同單木位置識別方法在各個指標方面的效果。①總體精度(AO)即正確識別的單木總數(shù)占實地調查單木實際總數(shù)的比值，公式為：AO=NC/NF×100%。②識別精度(AD)即正確識別的單木總數(shù)占識別出的單木總數(shù)的比值，公式為：AD=NC/ND×100%。③漏分誤差(EO)即未識別到的單木總數(shù)占實地調查單木實際總數(shù)的比值，公式為：EO=(NF-NC)/NF×100%。④錯分誤差(EC)即被錯誤識別為單木的總數(shù)占識別出的單木總數(shù)的比值，公式為：EC=(ND-NC)/ND×100%。⑤重復性誤差(ER)即對同一株單木重復識別的總數(shù)占識別出的單木總數(shù)的比值，公式為：ER=NR/ND×100%。

2.5.2 單木樹冠提取精度評價 采用單株冠幅提取精度 (ω)來進行評價，公式為：ω=100%-∣(C1-C0)/C0∣×100%。其中：C0為實測冠幅，C1為單木對應的提取冠幅。

3 結果與分析

3.1 不同單木位置識別方法的識別效果

利用4種方法進行單木位置識別，空間分布圖(圖2)顯示：4種方法均存在漏分、錯分和重復性識別的現(xiàn)象，結合表2可知：4種方法的識別精度(AD)由高到低依次為 MSP (89.30%)、LMC (85.60%)、RGP (77.50%)、MSC (70.00%)；MSP 法的單木位置識別效果要好于其他3種方法。4種方法的識別精度均小于總體精度(AO)，說明存在單木錯誤識別現(xiàn)象；相比之下，MSP法錯分誤差最低(EC為8.0%),，MSC法錯分誤差值最高(EC為30.1%)，說明依靠鄰近相似特征識別單木的方法更容易產生將單木識別錯誤。相比RGP和MSC法，MSP法的漏分誤差(EO)、錯分誤差(EC)和重復性誤差(ER)均較低，說明高郁閉度針葉林中，按照鄰域分屬識別單木位置會受到樹冠重疊度的影響，并不適用華北落葉松的單木位置識別。盡管LMC法和MSP法的總體精度都高于90%，但后者錯分誤差(Ec)較低，說明MSP法更能避免將非樹對象錯誤識別為單木而增加識別總數(shù)的現(xiàn)象。

表2 4 種單木位置識別方法的精度統(tǒng)計Table 2 Accuracy statistics of four methods on individual tree detection

圖2 4 種單木位置識別的空間分布Figure 2 Spatial distribution of the result of individual tree position detection based on four methods

3.2 改進均值漂移法的樹冠提取效果

采用均值漂移方法提取華北落葉松單株樹冠(圖3)。以單株華北落葉松實測冠幅為真值，對利用MSP法提取的單株冠幅進行精度分析，結合表3可知：MSP法總體平均提取精度為90.18%。單株提取精度小于80%的占總體結果的10.20%，單株提取精度大于80%的占總體結果的89.80%。單株樹冠分割偏小的結果占總體的62.07%，要遠遠高于提取冠幅偏大的結果。

表3 MSP 法提取單株冠幅精度Table 3 Precision distribution of extracting individual canopy diameter by MSP

圖3 改進均值漂移法樹冠提取結果 (局部)Figure 3 Crown extraction results of improved mean shift method (local)

4 討論與結論

4.1 討論

本研究提出的改進均值漂移法對高郁閉度華北落葉松林的單木提取總體精度達91.3%，這與張怡卓等[25]對中低密度人工林的研究結果(總體精度達89.0%)相似，遠高于JING等[27]對加拿大高郁閉林分的總體精度(69.0%)。相比于其他3種提取方法，MSP法有2個優(yōu)勢。一是在林木樹冠重疊的情況下，MSP法直接利用了點云空間特征，避免了生成CHM過程中產生的柵格誤差；二是該方法利用統(tǒng)計迭代的方式確定聚類點位，根據訪問頻率分類，比直接比較鄰域特征而進行分類穩(wěn)健性更好。

本研究發(fā)現(xiàn)：MSP法提取單木冠幅的平均提取精度為90.18%，比同等郁閉度林分的樹冠提取精度要高出很多[28]，可滿足森林調查需求。約10%的單株冠幅提取精度低于80%，這些點主要位于林分下層。由于上層遮擋，下層點云密度遠低于上層，而點云聚類算法和點云密度息息相關，這也是本研究發(fā)現(xiàn)冠幅提取值偏小的原因。具體來說，由于上層冠層存在重疊，人工目視解譯可以考慮到冠層重疊的影響并準確提取冠幅，而在三維點云圖上，樹冠界限分類時存在距離閾值類似的情況，影響提取結果；實測中，研究者往往通過東—西、南—北方向來測量冠幅，但是這只是一種理想狀態(tài)下的測量，在現(xiàn)實中，華北落葉松林木樹冠是不規(guī)則的，測量因此會存在誤差；同時研究區(qū)域的坡度、地形等的限制，也會對人工測量造成影響。

4.2 結論

本研究基于點云單一空間特征，引入高斯核函數(shù)簡化并改進了均值漂移量的計算過程，同時所提方法不需要先驗知識，算法的穩(wěn)健性高，解決了高郁閉度的華北落葉松林受復雜冠層影響的問題，提高了在該林分內機載LiDAR單木識別的精度，為高郁閉度華北落葉松林單木識別提供了新的方法。