霍朗寧 張曉麗

(北京林業(yè)大學(xué)精準(zhǔn)林業(yè)北京市重點(diǎn)實(shí)驗(yàn)室 北京 100083)

激光雷達(dá)(light detection and ranging，LiDAR)是一種集激光、全球定位系統(tǒng)(GPS)和慣性導(dǎo)航系統(tǒng)(inertial navigation systems，INS)于一身的用于獲取地面及地面目標(biāo)三維空間信息的主動(dòng)式雷達(dá)探測(cè)技術(shù)系統(tǒng)。將LiDAR應(yīng)用于林業(yè)領(lǐng)域，不僅能夠生成地面地形特征，還能夠得到林木樹冠表面和垂直結(jié)構(gòu)信息，如樹高、冠幅、胸徑、生物量等，可為森林生物量、生長(zhǎng)模型等研究提供數(shù)據(jù)支撐(N?sset，2002；Hyypp?etal.，2001)。近年來，隨著科學(xué)技術(shù)快速發(fā)展，LiDAR數(shù)據(jù)獲取成本不斷下降，應(yīng)用于林業(yè)領(lǐng)域的研究也越來越多，特別是機(jī)載LiDAR，探測(cè)森林垂直結(jié)構(gòu)的優(yōu)勢(shì)日益突出，在林分和單木尺度森林信息提取方面的探索逐步深入(Pedersenetal.，2012；Hyypp?etal.，2012；Packalenetal.，2015)。在單木尺度方面，如果能夠通過機(jī)載LiDAR數(shù)據(jù)獲取林分中單木位置、樹高、冠幅等信息，將極大提高森林經(jīng)營(yíng)管理效率，為林木生長(zhǎng)動(dòng)態(tài)評(píng)估、林分擇伐作業(yè)、森林可視化經(jīng)營(yíng)等提供數(shù)據(jù)支持。

在LiDAR單木提取研究中，樹冠2D分割和3D空間信息分割是2種主要的分割方式(Hyypp?etal.，1999；Solbergetal.，2006；Tangetal.，2007；Ferrazetal.，2012)。其中，常用的2D分割方法有分水嶺分割(Tangetal.，2007；Eneetal.，2012)、區(qū)域增長(zhǎng)法分割(Hyypp?etal.，1999；2001；Solbergetal.，2006)等，但這些方法存在一個(gè)共通的問題，即依賴CHM進(jìn)行分割操作，CHM的濾波效果會(huì)極大影響最終提取效果(Eneetal.，2012)，且很可能導(dǎo)致信息丟失。從原理上分析，基于CHM構(gòu)建并進(jìn)行分割的方法仍為基于表面模型信息的分割，往往是在對(duì)最上層點(diǎn)云進(jìn)行分割、劃分單木位置和冠幅范圍后將下層點(diǎn)云與之疊加，由于中下層點(diǎn)云利用率較低，數(shù)據(jù)信息沒有充分挖掘，復(fù)雜林分結(jié)構(gòu)條件下幾乎無法提取到中下層林木，因此研究結(jié)果應(yīng)用范圍受限。相比之下，3D空間信息分割則利用每一個(gè)點(diǎn)云的空間信息，單木提取更全面，信息也更完整(Lietal.，2012)。目前，常用的3D分割方法主要包括K均值聚類(Morsdorfetal.，2004；Guptaetal.，2010)、均值移動(dòng)聚類(Ferrazetal.，2012)、自適應(yīng)距離聚類(Leeetal.，2010；Guptaetal.，2010)等，然而這些方法并不一定都能獲得很好的提取效果，一些研究指出其中的單木提取方法誤差較大(Kaartinenetal.，2012；Eysnetal.，2015)，而對(duì)于提取單木正確率能夠達(dá)到多少的問題，很多研究避而不談(Guptaetal.，2010；Tangetal.，2013)。3D分割方法存在的困難主要有以下幾點(diǎn)：1)對(duì)高郁閉度林分下層林木提取能力較低，林分密度、郁閉度、林分結(jié)構(gòu)異質(zhì)性程度嚴(yán)重影響提取效果，提取效果較好的往往是林分密度較低的(Lietal.，2012)；2)提取林木有效性問題，由于3D分割方法對(duì)下層林木提取能力較低，導(dǎo)致能夠準(zhǔn)確提取的林木數(shù)量占總體比例偏低(Kandareetal.，2016)，如果沒有地面實(shí)測(cè)數(shù)據(jù)進(jìn)行匹配，將不能得知是哪部分林木被準(zhǔn)確提??；3)點(diǎn)云密度問題，盡管有研究指出點(diǎn)云密度對(duì)單木樹高、冠幅精度影響不大(Reitbergeretal.，2009；Kaartinenetal.，2012；Yaoetal.，2014)，但在復(fù)雜林分結(jié)構(gòu)條件下，點(diǎn)云密度過低將造成下層林木信息不完整，影響更新木的提取效果。

此外，以往研究在反演森林結(jié)構(gòu)方面沒有過多關(guān)注，LiDAR反演森林結(jié)構(gòu)的探索仍停留在林分尺度，且構(gòu)造的指標(biāo)僅適用于LiDAR數(shù)據(jù)，并沒有將其與實(shí)際生產(chǎn)實(shí)踐中常用的森林結(jié)構(gòu)指標(biāo)進(jìn)行轉(zhuǎn)換，從而導(dǎo)致應(yīng)用性較低。常用的森林結(jié)構(gòu)指標(biāo)如樹高基尼指數(shù)、樹高變異系數(shù)、聚集指數(shù)、競(jìng)爭(zhēng)指數(shù)等，在林分空間結(jié)構(gòu)優(yōu)化、撫育改造等森林經(jīng)營(yíng)管理工作中具有重要指導(dǎo)意義，如果能夠通過LiDAR數(shù)據(jù)反演出高精度的指標(biāo)，則可進(jìn)一步將遙感調(diào)查結(jié)果向生產(chǎn)實(shí)踐應(yīng)用轉(zhuǎn)化。

鑒于此，本研究聚焦上述問題，改進(jìn)單木提取策略和算法，主要特色和創(chuàng)新點(diǎn)如下：1)挑選1塊林分密度高、空間結(jié)構(gòu)異質(zhì)性大的復(fù)層林作為研究對(duì)象，重點(diǎn)關(guān)注中下層林木的提取效果；2)在聚類前將點(diǎn)云劃分為多個(gè)水平層，增加下層林木探測(cè)能力，同時(shí)制定合理的點(diǎn)云簇融合方法，重點(diǎn)分析水平層設(shè)置情況對(duì)信息提取精度的影響；3)構(gòu)建合理的精度評(píng)價(jià)指標(biāo)，重點(diǎn)關(guān)注提取林木的正確性、有效性，對(duì)提取效果進(jìn)行綜合評(píng)價(jià)，提升成果的可用性。

1 研究區(qū)概況與研究方法

1.1 研究區(qū)概況

研究區(qū)地處我國(guó)甘肅省祁連山國(guó)家級(jí)自然保護(hù)區(qū)肅南裕固族自治縣西水林場(chǎng)，為祁連山山區(qū)天然次生林、黑河流域水源涵養(yǎng)林，海拔2 700～3 000 m，屬溫帶高寒半干旱、半濕潤(rùn)山地森林草原氣候，陽(yáng)坡為山地草原，陰坡為森林景觀，林區(qū)主要森林類型為青海云杉(Piceacrassifolia)天然純林(成熟林)。研究區(qū)位置、樣地布置和原始點(diǎn)云示意如圖1所示。

圖1 研究區(qū)位置、樣地布置和原始點(diǎn)云示意

1.2 數(shù)據(jù)獲取與預(yù)處理

1.2.1 機(jī)載LiDAR數(shù)據(jù)與預(yù)處理 機(jī)載LiDAR數(shù)據(jù)獲取于2008年6月23日，無人機(jī)平均飛行海拔約3 560 m，距離地面飛行高度約760 m；搭載LiteMapper5600 LiDAR系統(tǒng)，激光掃描儀為RiegILMS-Q560，波長(zhǎng)1 550 nm，激光脈沖長(zhǎng)度3.5 ns，激光脈沖發(fā)散角小于等于0.5 m·rad，地面平均光斑直徑38 cm，可分辨目標(biāo)最小間隔0.6 m，脈沖重復(fù)頻率50 kHz。

從原始點(diǎn)云中裁出與地面實(shí)測(cè)數(shù)據(jù)相匹配的100 m×100 m區(qū)域，采用Terrasolid軟件進(jìn)行點(diǎn)云去噪、地面點(diǎn)與非地面點(diǎn)分離、植被標(biāo)準(zhǔn)化，并去掉高度小于2 m的點(diǎn)云，獲得樣地林層點(diǎn)云作為算法的試驗(yàn)數(shù)據(jù)。平均點(diǎn)云密度為2.86 m-2。

1.2.2 地面調(diào)查數(shù)據(jù) 地面調(diào)查時(shí)間為2008年6月1—13日，沿山坡走勢(shì)布設(shè)一塊100 m×100 m超級(jí)樣地，方位角約122°，按照25 m×25 m尺寸劃分為16塊子樣地(圖1)。對(duì)樣地中林木進(jìn)行每木檢尺，記錄林木胸徑、樹高、冠幅和枝下高，通過全站儀獲取單木位置，共得到1 435株單木信息。本研究重點(diǎn)討論LiDAR反演單木樹高、冠幅、位置以及林分結(jié)構(gòu)參數(shù)的精度，相關(guān)統(tǒng)計(jì)量如表1所示。

1.3 單木分割方法

聚類是一種無監(jiān)督學(xué)習(xí)，其將相似的對(duì)象歸到同一個(gè)簇中，簇內(nèi)對(duì)象越相似，聚類效果越好。本研究應(yīng)用K均值聚類算法，首先從數(shù)據(jù)集中隨機(jī)選取K個(gè)點(diǎn)作為初始聚類中心，然后計(jì)算各樣本到聚類中心的距離，將樣本歸到距離其最近的聚類中心所在的類。計(jì)算新形成的每一個(gè)聚類數(shù)據(jù)對(duì)象的平均值得到新的聚類中心，如果相鄰2次的聚類中心沒有任何變化，樣本調(diào)整結(jié)束，說明聚類準(zhǔn)則函數(shù)已收斂?；贛ATLAB中的聚類算法，完成單木分割、信息提取和精度評(píng)價(jià)。分層聚類方法提取單木工作流程如圖2所示。

圖2 分層聚類方法提取單木工作流程

1)點(diǎn)云水平切片 對(duì)于郁閉度高、株數(shù)密度大的復(fù)層林，直接對(duì)所有點(diǎn)云進(jìn)行聚類分析存在2個(gè)問題：一是數(shù)據(jù)量大、運(yùn)算速度慢，每加入一個(gè)新的聚類數(shù)據(jù)，則需要重新計(jì)算其聚類中心，將點(diǎn)云直接聚為2類的時(shí)間遠(yuǎn)大于聚為4類后對(duì)相近2類進(jìn)行融合的時(shí)間；二是初始聚類中心的選擇依賴于局部最大值位置，對(duì)于樹高異質(zhì)性大的林分，只對(duì)整體提取局部最大值會(huì)不可避免錯(cuò)過被壓木的樹冠頂點(diǎn)，從而無法分割優(yōu)勢(shì)層下的低矮木和更新木。通過對(duì)點(diǎn)云進(jìn)行水平切片(圖3)，可提取到各林層的局部最大值，尤其更利于下層林木頂點(diǎn)的偵測(cè)。切片片層高可根據(jù)林分分層情況人為設(shè)置，也可依照點(diǎn)云在垂直方向的分布規(guī)律自動(dòng)分層，保證每層點(diǎn)云數(shù)量相同。

圖3 點(diǎn)云水平切片和局部最大值位置提取

2)點(diǎn)云局部最大值位置提取 決定聚類算法運(yùn)行效率和聚類效果的主要因素之一是初始聚類中心選擇，點(diǎn)云數(shù)量龐大，初始聚類中心越接近最終聚類中心，運(yùn)算效率越高。本研究區(qū)主要樹種為青海云杉，其形態(tài)近似圓錐體，樹頂點(diǎn)大致位于該樹垂直投影的形狀中心，因此通過提取點(diǎn)云局部最大值位置確定樹頂(圖3)，并將其作為初始聚類中心。利用MATLAB中的imregionalmax函數(shù)提取局部最大值點(diǎn)，其原理是通過一系列腐蝕、膨脹運(yùn)算，從1個(gè)點(diǎn)的26個(gè)鄰近點(diǎn)判斷該點(diǎn)是否為局部最大值。由于點(diǎn)云數(shù)量龐大、異質(zhì)性高，因此需對(duì)點(diǎn)云進(jìn)行平滑濾波后再提取局部最大值。本研究采用高斯平滑濾波器對(duì)點(diǎn)云進(jìn)行多次濾波，再提取點(diǎn)云局部最大值，所得結(jié)果為單木樹頂點(diǎn)坐標(biāo)。

3)點(diǎn)云縱向壓縮K均值聚類算法適合于三維特征空間內(nèi)大體為球體、簇密度均勻的數(shù)據(jù)集聚類分析，本研究區(qū)青海云杉林形態(tài)近似圓錐體，因此需對(duì)其進(jìn)行縱向壓縮。采用3倍縮放比例，即點(diǎn)云x、y坐標(biāo)不變，z坐標(biāo)縮放為原來的1/3。

4)點(diǎn)云層聚類 計(jì)算每個(gè)縮放后的點(diǎn)云數(shù)據(jù)到每個(gè)聚類中心的歐氏距離，將點(diǎn)云歸到距離其最近的聚類中心所在的類，并重新計(jì)算聚類中心。重復(fù)該步驟至聚類中心不再變化，即完成點(diǎn)云聚類，可獲得標(biāo)注不同類別屬性的點(diǎn)云簇。

5)隸屬不同點(diǎn)云層的點(diǎn)云簇融合 對(duì)接近分割線且聚類中心距離滿足一定條件的點(diǎn)云簇進(jìn)行合并。合并條件具體設(shè)置和計(jì)算過程如下：經(jīng)凸包處理求得點(diǎn)云簇面積，將點(diǎn)云簇視為圓求圓形半徑作為點(diǎn)云簇半徑；融合條件設(shè)置為上下2層點(diǎn)云簇中心的水平距離不超過下層點(diǎn)云簇半徑。合并順序如下：從最上層開始，先合并1、2兩個(gè)點(diǎn)云層的點(diǎn)云簇，得到新的點(diǎn)云簇和聚類中心，再與第3層點(diǎn)云簇合并，以此類推。合并完所有符合條件的點(diǎn)云簇后，重新計(jì)算點(diǎn)云簇中心，作為單木位置信息。

1.4 參數(shù)反演

1.4.1 單木參數(shù)反演 單木分類工作完成后，可獲得標(biāo)注不同類別屬性的點(diǎn)云，具有同一屬性的點(diǎn)云為同一株單木的點(diǎn)云簇。對(duì)點(diǎn)云簇的x、y進(jìn)行平均，作為單木位置信息；提取點(diǎn)云簇中z的最大值，作為單木樹高。在冠幅提取方面，由于點(diǎn)云密度較低，尤其在提取最下層被壓木時(shí)，存在因點(diǎn)云數(shù)量過少、點(diǎn)云共線而無法形成凸包的情況，因此不能應(yīng)用凸包算法提取林木冠幅。本研究通過點(diǎn)云密度、數(shù)量計(jì)算冠幅面積的方式求得冠幅，即：

(1)

式中：c為冠幅；d為點(diǎn)云密度；n為點(diǎn)云簇中的點(diǎn)云個(gè)數(shù)。

1.4.2 林分結(jié)構(gòu)參數(shù)反演 獲得所有單木的位置、樹高和冠幅后，進(jìn)一步計(jì)算林分結(jié)構(gòu)參數(shù)。

1)林分樹高分布 在林分中，不同樹高林木按樹高組的分配狀態(tài)，以2 m為一級(jí)，分別統(tǒng)計(jì)各樹高級(jí)范圍內(nèi)的林木數(shù)量。為了量化LiDAR提取林分樹高分布與實(shí)際情況的差異性，計(jì)算各樹高級(jí)林木數(shù)量均方根誤差：

(2)

式中：n為樹高級(jí)數(shù)；N-test、N-ref分別為各樹高級(jí)中LiDAR提取數(shù)量和地面實(shí)測(cè)數(shù)量。

2)林分冠幅分布 在林分中，不同冠幅林木按樹高組的分配狀態(tài)，原理同樹高分布。

3)樹高基尼指數(shù)(Gc) 樹高基尼指數(shù)是目前成功應(yīng)用于樹高不均等性測(cè)度的極少數(shù)指標(biāo)之一(Gini，1912)，為樣地中所有個(gè)體間某一測(cè)定指標(biāo)成對(duì)比較的差值絕對(duì)值的算術(shù)平均，范圍在0～1之間，用于表示樹高指標(biāo)的不均勻性，即群落的樹高多樣性，指數(shù)值越接近1，則林分樹高越不均勻，指數(shù)值越接近0，則越均勻。計(jì)算公式為：

(3)

式中：Xi為第i株樹的樹高；i為研究區(qū)所有林木按高度從小到大排序后的順序；N為林木總數(shù)。

4)樹高變異系數(shù)(Cv) 樹高變異系數(shù)可反映林分垂直結(jié)構(gòu)的變化程度，為樣地中所有樹高標(biāo)準(zhǔn)差和平均值的比較。計(jì)算公式為：

(4)

5)聚集指數(shù)(Clark-Evans’s aggregation indes)聚集指數(shù)是常用的林木空間分布格局指數(shù)(Clorketal., 1954)，為樣地種群內(nèi)每個(gè)個(gè)體和其最近鄰體間的平均距離與其為隨機(jī)分布時(shí)的期望平均距離比較，所得比率為一隨機(jī)性偏離測(cè)度值，可以根據(jù)其統(tǒng)計(jì)顯著性判斷種群分布格局類型，取值范圍在0～2.149 6之間，若等于1為隨機(jī)分布，大于1為均勻分布，小于1為聚集分布。計(jì)算公式為：

(5)

式中：r0為樣地中每個(gè)個(gè)體和其最近鄰近體間的平均距離；ri為第i個(gè)個(gè)體和其最鄰近體間的距離；E(r)為樣地中所有個(gè)體隨機(jī)分布時(shí)r0的期望值；A為樣地面積(m2)。

1.5 精度評(píng)價(jià)

1.5.1 參數(shù)名稱設(shè)置 為便于結(jié)果呈現(xiàn)，本研究涉及的單木參數(shù)名稱及縮寫如表2所示。

表2 本研究涉及的單木參數(shù)名稱及縮寫

1.5.2 單木聚類效果評(píng)價(jià)指標(biāo) 聚類經(jīng)常會(huì)出現(xiàn)欠分類或過分類的情況，在聚類效果評(píng)價(jià)中，數(shù)量精度評(píng)價(jià)尤為重要。從匹配原理來說，由于匹配時(shí)較看重位置和樹高因素，而過分類對(duì)樹高影響不大，只增加了林木位置點(diǎn)，更容易將原本不能匹配的林木匹配在一起，當(dāng)聚類數(shù)量增加甚至超過實(shí)際數(shù)量時(shí)，可以成功匹配實(shí)測(cè)數(shù)據(jù)的林木數(shù)量往往隨之增加，因此片面提升成功匹配林木比例或單獨(dú)關(guān)注聚類數(shù)量與實(shí)際數(shù)量的差異，都不能全面反映聚類效果。本研究設(shè)置3個(gè)精度指標(biāo)，具體如下。

1)聚類數(shù)量精度(A-QT) 即單位林分中聚類出的單木數(shù)量與地面實(shí)測(cè)單木數(shù)量的差異性，代表聚類效果的數(shù)量精確程度，公式如下：

A-QT=1-|N-test-N-ref|/N-ref。

(6)

式中：N-test為聚類的單木數(shù)量；N-ref為地面實(shí)測(cè)單木數(shù)量。

2)聚類質(zhì)量精度(A-QL) 即聚類出的單木中可與地面實(shí)測(cè)單木相匹配數(shù)量占總體的比例，代表聚類質(zhì)量的好壞，公式如下：

A-QL=1-|N-test-N-match|/N-test。

(7)

式中：N-match為成功匹配的單木數(shù)量。

3)聚類有效精度(A-ED) 即單位林分中可匹配的聚類單木數(shù)量與地面實(shí)測(cè)單木數(shù)量的差異性，代表聚類效果相對(duì)實(shí)際情況而言有效部分所占比例，公式如下：

A-QL=1-|N-ref-N-match|/N-ref。

(8)

以數(shù)量精度、質(zhì)量精度和有效精度三者的平均值作為總體聚類效果的評(píng)價(jià)指標(biāo)，定義為聚類效果指數(shù)(A-k),代表聚類方法對(duì)單木的偵測(cè)和提取能力，公式如下：

A-k=(A-QT+ A-QL+ A-ED)/3。

(9)

2 結(jié)果與分析

2.1 不同分層方式對(duì)精度的影響

為了對(duì)比不同分層方式對(duì)精度的影響，暫將點(diǎn)云分為2層并設(shè)置分層界限高分別為3、4、5…12 m，統(tǒng)計(jì)聚類效果指數(shù)(A-k)、單木樹高精度(A-h)、單木冠幅精度(A-c)、林分平均高精度(A-H)、樹高基尼指數(shù)精度(A-Gc)如圖4a所示。可以看出，單木樹高精度、單木冠幅精度隨層高增大變動(dòng)幅度較小，規(guī)律性不大；而樹高基尼指數(shù)精度、林分平均高精度、聚類效果指數(shù)隨層高增大呈單峰曲線變化，尤其是樹高基尼指數(shù)精度在50%～85%范圍內(nèi)劇烈變化，三者的最佳點(diǎn)分別對(duì)應(yīng)層高5、7和9 m，且存在此消彼長(zhǎng)現(xiàn)象，無法選出最優(yōu)情況，故將點(diǎn)云分為2層不能同時(shí)滿足各方面要求，需進(jìn)行多層聚類。

圖4 不同分層方式對(duì)精度的影響

在多層聚類中，為了避免人為設(shè)置層高對(duì)精度的影響，統(tǒng)一采用均分百分位數(shù)方式進(jìn)行分層，即保證每層點(diǎn)云數(shù)量相同。由4b可知，分層數(shù)對(duì)林分水平參數(shù)精度的影響很大，尤其是樹高基尼指數(shù)精度變動(dòng)最大，由65%升至90%，林分平均高精度和聚類效果指數(shù)也有明顯提升；而單木樹高精度變化不明顯，單木冠幅精度顯著下降。

2.2 參數(shù)提取精度

由表3可知，無論分層數(shù)設(shè)置多少，單木樹高精度均在90%以上，可滿足林業(yè)調(diào)查對(duì)精度的要求，而單木冠幅精度較低，在60%左右，說明整體冠幅提取難度較大，需進(jìn)一步改進(jìn)；此外，聚類效果指數(shù)均在80%以上，最高可達(dá)88.28%，說明提取有效性較高，而林分平均高精度最高可達(dá)99.84%，樹高基尼指數(shù)精度最高可達(dá)89.65%，說明對(duì)林分整體信息和結(jié)構(gòu)反演能力較強(qiáng)。

表3 分層數(shù)設(shè)置對(duì)聚類效果和精度的影響

為了解具體影響情況，按樹高組的分配狀態(tài),以5 m為一級(jí)，分別對(duì)各樹高級(jí)的單木樹高精度和聚類效果指數(shù)進(jìn)行分析。由圖5a可知，無論分層數(shù)設(shè)置多少，精度曲線大致趨勢(shì)相同，即10～20 m主林層單木樹高精度較高，達(dá)94%左右，而少數(shù)樹高達(dá)20 m以上的林木，由于點(diǎn)云密度過低，可能錯(cuò)過樹頂點(diǎn)造成樹高精度較低。

對(duì)比不同分層數(shù)的單樹高精度和聚類效果(圖5b)可知，分層數(shù)對(duì)5 m以上主林層單木樹高精度影響不大，整體單木樹高精度差異主要是對(duì)5 m以下更新層林木的聚類。隨分層數(shù)不斷增加，對(duì)2～5 m更新層林木的提取能力顯著增強(qiáng)，由最低的34.86%升至74.43%，而其他林層聚類效果變動(dòng)不大，直接導(dǎo)致整體LiDAR提取能力較顯著提升，但也由此帶來了單木樹高、冠幅精度的波動(dòng)。更新層點(diǎn)云密度比主林層低，平均每株林木僅有4.9個(gè)點(diǎn)，聚類誤差引起1個(gè)點(diǎn)的聚類錯(cuò)誤，即會(huì)對(duì)單木整體造成較大影響。對(duì)更新層林木的提取也是所有LiDAR單木提取研究中較為困難的一部分。

圖5 不同分層數(shù)對(duì)單木樹高精度和聚類效果的影響

2.3 提取有效性

以分2層聚類結(jié)果為例，將LiDAR提取出的林木與地面實(shí)測(cè)林木以5 m為一級(jí)，分別統(tǒng)計(jì)各樹高級(jí)范圍內(nèi)的林木數(shù)量，計(jì)算二者比值代表該樹高級(jí)的提取率，超過100%的部分視為無效。不同層高對(duì)應(yīng)的各樹高級(jí)提取率如圖6所示。

圖6 層高3～12 m時(shí)各樹高級(jí)提取率

以7、10 m為界，各樹高級(jí)提取率大致分2種情況：1)分層高小于7 m或大于等于10 m時(shí)，10 m以下林木提取率較低，10 m以上林木提取率虛高，均超過100%，最高甚至達(dá)194%，即有近一半林木提取錯(cuò)誤，沒有實(shí)際意義；而決定最終聚類效果指數(shù)的往往是10 m以下林木提取的有效性；2)分層高在7～10 m之間，情況正相反。總體來說，提取率曲線在100%以下的面積可代表整體提取有效性高低，如以7 m為分層高，各樹高級(jí)提取率較為均衡，有效性最高。

2.4 林分結(jié)構(gòu)信息反演

對(duì)聚類結(jié)果和地面實(shí)測(cè)林木按樹高統(tǒng)計(jì)數(shù)量，結(jié)果如圖7所示。可以看出，不分層時(shí)樹高近似呈正態(tài)分布，而實(shí)際樹高分布為雙峰曲線，聚類結(jié)果與實(shí)際情況嚴(yán)重不符；相比之下，分層聚類對(duì)樹高分布的預(yù)測(cè)能力更為突出，尤其對(duì)10 m以上主林層，無論是數(shù)量還是比例均與實(shí)際情況十分吻合，而樹高2～8 m時(shí)還存在一定誤差，可通過提升對(duì)下層林木的提取能力來攻克這一難題。此外由上述可知，林分平均高精度最高達(dá)99.84%，樹高基尼指數(shù)精度最高達(dá)89.65%，樹高變異系數(shù)最高達(dá)89.23%，同時(shí)計(jì)算出最鄰近體間平均距離為1.53 m，與實(shí)際的1.31 m較為接近，均體現(xiàn)了分層聚類算法在林分結(jié)構(gòu)反演上的突出優(yōu)勢(shì)。

圖7 聚類結(jié)果與地面實(shí)測(cè)林木樹高分布情況對(duì)比

3 討論

3.1 LiDAR提取林木的有效性及評(píng)價(jià)指標(biāo)

以往LiDAR單木提取研究中，評(píng)判算法優(yōu)劣時(shí)較看重提取率、匹配精度和單木樹高精度等指標(biāo)。以分水嶺分割為代表的分割表面模型算法以及不分層、直接分割的相關(guān)算法由于對(duì)林下層點(diǎn)云利用率低，造成林木分割數(shù)量少于實(shí)際數(shù)量，而基于LiDAR提取林木只要不超過樹高、距離閾值均視為有效，因此LiDAR提取林木的有效性問題并不突出。但本研究采用的分層K均值聚類算法，很有可能因過分類出現(xiàn)LiDAR提取林木數(shù)量多于實(shí)際數(shù)量的情況，在實(shí)際應(yīng)用過程中，如果沒有地面實(shí)測(cè)數(shù)據(jù)匹配驗(yàn)證，有效林木和無效林木將無法區(qū)分，并可能進(jìn)一步影響整體結(jié)果的有效性。

以往研究表明，林木提取率(提取數(shù)量與實(shí)際數(shù)量比)在20%～60%之間不等(Morsdorfetal.，2004；Kandareetal., 2016；尹艷豹，2010)，大多數(shù)研究更加關(guān)注提取出并匹配的林木在三維信息方面的精度，忽視了提取有效性問題，即沒有驗(yàn)證提取出的林木究竟是主林層、優(yōu)勢(shì)木還是包括多少比例被壓木，因此盡管提取出的部分林木信息很精確，但由于不清楚其與整體的關(guān)系，不能夠確定這部分信息是否可以代表總體情況，在實(shí)際生產(chǎn)應(yīng)用中意義不大。在LiDAR單木提取研究中，追求能夠“更準(zhǔn)確”地提取出“更多”的林木，其中“更準(zhǔn)確”以單木信息精度來反映，而“更多”則以提取有效性來體現(xiàn)，即本研究中設(shè)置的“聚類效果指數(shù)”，用該指數(shù)綜合判斷“LiDAR提取了多少林木”“提取的林木中有多少比例是正確的”等問題，充分考慮了結(jié)果的可用性。聚類效果指數(shù)越高，提取有效性越高，只有達(dá)到一定標(biāo)準(zhǔn)，才能說明算法可以擺脫對(duì)地面實(shí)測(cè)數(shù)據(jù)的依賴，能夠直接應(yīng)用于生產(chǎn)實(shí)踐中。

3.2 LiDAR點(diǎn)云密度對(duì)單木及林分結(jié)構(gòu)信息的影響

點(diǎn)云密度對(duì)信息提取的影響是近年來比較熱門的研究?jī)?nèi)容之一。在單木信息精度方面，研究指出點(diǎn)云密度高于20 m-2帶來的精度提升微乎其微(Reitbergeretal.，2009；Kaartinenetal.，2012；Yaoetal.，2014)，點(diǎn)云密度2、4、8 m-2(Kaartinenetal.，2012)的單木信息提取精度差異不顯著，因此，很多后續(xù)研究采用先稀釋點(diǎn)云再進(jìn)行信息提取的方法(Magnussenetal.，2010；Hansenetal.，2015)。本研究中，林冠層平均點(diǎn)云密度為2.86 m-2，所得單木樹高精度達(dá)90%以上，單木冠幅精度達(dá)60%以上，與大多數(shù)研究結(jié)果一致，驗(yàn)證了該密度水平的點(diǎn)云數(shù)據(jù)可以很好完成單木信息提取工作。

然而在林分結(jié)構(gòu)方面，有研究對(duì)60 m-2點(diǎn)云先進(jìn)行不同程度稀釋后再反演林分結(jié)構(gòu)信息，包括林分密度、水平結(jié)構(gòu)和垂直結(jié)構(gòu)，結(jié)果發(fā)現(xiàn)60 m-2與6 m-2點(diǎn)云的效果存在顯著差異，極高密度點(diǎn)云在林分結(jié)構(gòu)反演方面更具有優(yōu)勢(shì)(Kandareetal., 2016)。本研究中，通過點(diǎn)云多次分層再聚類的方法反演林分結(jié)構(gòu)信息也得到很好效果，當(dāng)分層數(shù)大于8時(shí)，樹高基尼指數(shù)、樹高變異系數(shù)精度在85%以上，不以增加10倍信息量為代價(jià)也能夠較好反演森林結(jié)構(gòu)參數(shù)。

在提取有效性方面，目前還未有研究專門探討點(diǎn)云密度對(duì)提取有效性的影響。本研究中，通過點(diǎn)云分層可提升對(duì)下層林木的提取能力，但受困于下層林木的點(diǎn)云密度過低，平均每株林木僅有4.9個(gè)點(diǎn)云，在聚類時(shí)1個(gè)點(diǎn)分類錯(cuò)誤就有可能引起整體信息的劇烈變化，因此猜測(cè)如能提高點(diǎn)云密度，將會(huì)提升對(duì)下層林木的提取能力，從而提升整體的提取有效性。5 m以上的單木樹高精度達(dá)90%以上，提取有效性也在90%以上，可反映整體林分除更新層外其他林木的生長(zhǎng)情況；如果能夠通過提高點(diǎn)云密度從而攻克更新層林木的提取問題，將在森林精準(zhǔn)量測(cè)、林木生長(zhǎng)模型構(gòu)建、可視化經(jīng)營(yíng)等方面起到一定的推進(jìn)作用。

4 結(jié)論

本研究提出一種基于機(jī)載LiDAR點(diǎn)云多層聚類的單木信息提取方法，并將其應(yīng)用于1塊林木密度高、空間結(jié)構(gòu)異質(zhì)性大的復(fù)層林進(jìn)行驗(yàn)證，同時(shí)對(duì)比不同匹配方法對(duì)同一聚類結(jié)果的匹配精度。結(jié)果表明：1)分多層次K均值聚類進(jìn)行單木提取具有很好效果，尤其在下層林木的提取方面表現(xiàn)突出；2)可在低點(diǎn)云密度條件下完成工作，降低高密度點(diǎn)云數(shù)據(jù)信息的冗余；3)在單木信息精度方面，分層K均值聚類提取的單木樹高精度達(dá)90%以上，單木冠幅精度在60%左右；4)在林分結(jié)構(gòu)方面，采用多次分層方式可提升對(duì)林分結(jié)構(gòu)的反演能力，林分垂直結(jié)構(gòu)、水平結(jié)構(gòu)、樹高分布結(jié)構(gòu)等方面均有較好反映；5)通過構(gòu)建LiDAR林木提取有效性指標(biāo)，著重關(guān)注算法的適用性，提取林木有效性最高達(dá)88.70%，進(jìn)一步說明該方法可以擺脫對(duì)地面實(shí)測(cè)數(shù)據(jù)的依賴，能夠應(yīng)用于生產(chǎn)實(shí)踐中。