基于機載激光雷達冠層高度模型的小班區劃

熊 昊，龐 勇＊，莢 文，李春干

(1.中國林業科學研究院資源信息研究所，國家林業和草原局林業遙感與信息技術重點實驗室，北京 100091；2.廣西大學林學院，廣西 南寧 530004)

小班是內部特征基本一致，與相鄰地段有明顯區別而需要采取相同經營措施的森林地塊或小區，是森林資源統計、經營和管理的基本單位[1]。小班區劃通常是依據高分辨率遙感影像目視判讀或實地勾繪，工作量巨大，且會因調查員的知識水平、經驗、調查線路、觀察視角等主觀因素差異產生結果不一致的問題[2-3]。因此，對小班自動區劃進行研究十分必要。在小班的劃分條件[1]中，權屬、森林類別等因子是政策規定或造林時已經確定的，這些因子通常不因林分狀態的變化而改變，且可以從上一期林相圖中提取。而優勢樹種組成、齡級等因子可以通過遙感方法提取。

國內已有一些研究采用遙感影像來進行自動或半自動的小班區劃。這些研究通常采用Landsat、SPOT5、QuickBird、ALOS 等高分辨率衛星遙感影像作為輔助數據，使用多尺度分割和分類的方法獲得小班界線[3-6]。但衛星影像受到空間分辨率的約束，對小班區劃仍具有一定的局限性。國外有一些學者采用航空遙感影像進行小班區劃[7-8]，但多限于區分內部較均質的林分，沒有對我國森林調查規劃規程中規定的優勢樹種、小班面積等因子的嚴格界定。并且影像數據中不包含林分結構信息，可能會有樹種相同，但樹高、密度等不同的林分被劃分為同一小班，同時影像也可能有云和陰影等干擾因素。

近年來機載激光雷達（ALS）技術在林業中的應用日益廣泛[9]。ALS 數據是由搭載在有人機或無人機平臺上的激光雷達傳感器掃描獲取的，適用于森林結構觀測[10]。由ALS 數據生成的森林冠層高度模型可以為小班區劃提供準確的位置、林分平均高、林分密度、郁閉度等信息[11-12]。通過對CHM的分割可以將森林參數不同的地塊區分開，得到符合林區實際情況的小班邊界。Mustonen 等[13]比較了采用CHM 和RGB 航空遙感影像和綜合兩種數據的林分區劃，發現基于CHM 的區劃的林分比基于RGB 航空遙感影像和基于兩種數據區劃的林分更均勻。Koch 等[11]利用20 m 空間分辨率的低密度ALS 數據進行小班區劃，采用特征提取、創建和基于柵格的分類方法。Wu 等[14]從點云中提取樹的大小、林分密度和樹種特征，使用非監督均值漂移算法和區域生長法進行分割。Dechesne 等[15]從點云和多光譜影像中計算并篩選出了94 個特征，并使用隨機森林方法進行樹種監督分類得到林分區劃結果。Pukkala 等[16-18]利用16 m 的ALS 數據和森林資源調查數據，使用元胞自動機方法，通過調整柵格內部屬性參數來進行林分區劃。而國內現有研究中，僅有少數研究利用ALS 數據在ArboLiDAR軟件中進行了林分高度分割[19]，并用于反演云杉次生純林的林分平均高。Jia 等[20]討論了不同形狀和面積參數對元胞自動機林分區劃的影響。

已有的研究中體現了ALS 數據用于林分區劃的潛力，但這些研究沒有針對小班區劃需求，并且未充分考慮樹種信息。機載高光譜影像具有較高的空間和光譜分辨率，大量研究表明[21-25]，機載高光譜數據可以產生高精度的森林樹種分類結果，為小班區劃提供準確的樹種組成信息。因此，本研究利用機載激光雷達數據和高光譜影像樹種分類結果進行小班區劃，并與人工區劃小班進行對比，綜合探究在高光譜影像樹種信息的輔助下，ALS 數據生成的高精度CHM 在小班區劃中的應用潛力。

1 研究區概況及數據

1.1 研究區概況

研究區位于黑龍江省佳木斯市樺南縣東北部的孟家崗林場，地理坐標為130°32′～130°52′ E，46°20′～46°30′ N。地處完達山西麓余脈，以低山丘陵為主，坡度較為平緩，坡度在10°～20°之間。地勢東北高，西南低，平均海拔約為250 m。孟家崗林場人工林占林地面積的76.7%，其中落葉松（Larix olgensisHenry）、樟子松（Pinus sylvestrisvarmongolicaLitvin.）、紅松（Pinus koraiensisSieb.et Zucc.）約占80%。天然林主要分布在林場的東北方向，占全部森林面積27.5%，主要樹種包括白樺（Betula platyphyllaSuk.）、柞樹（Quercus mongolicaFisch.Ex Ledeb）和椴樹（Tilia tuanSzyszyl）。

1.2 數據

1.2.1 遙感數據 采用中國林業科學研究院的機載綜合遙感觀測平臺CAF-LiCHy[26]于2017 年5 月31 日至6 月15 日采集的LiDAR 數據，航攝平臺為運-5 多用途飛機，絕對飛行高度為1 000 m。CAF-LiCHy 機載系統包含全波形機載激光雷達（LMS-Q680i）、機載推掃式高光譜掃描儀（AISA EagleⅡ）和高分辨率CCD 相機，3 個傳感器共用一套POS 位置和姿態系統。飛行時航向和旁向交并比分別為70%和50%，獲取的數字正射影像空間分辨率為0.1 m。LMS-Q680i 激光傳感器是奧地利Riegl 公司的全波形激光雷達系統，工作波長為1 550 nm，激光束發散角為0.5 mrad，波形數據的記錄間隔為1 ns，最大脈沖重復頻率為400 kHz。激光雷達點云數據的平均密度為12 pts·m-2，航帶間點云誤差為水平0.1 m，垂直0.12 m。高光譜AISA Eagle II 傳感器光譜范圍400～970 nm，光譜分辨率為3.3 nm，空間分辨率為0.5 m。CAFLiCHy 機載系統同步獲取的激光雷達數據、數字正射影像和高光譜影像的空間位置的一致性精度優于1 m[27]。

1.2.2 樹種分類圖 樹種類別圖由CAF-LiCHy 機載系統同步獲取的高光譜影像監督分類獲得。對高光譜影像的每個波段使用灰度共生矩陣計算了均值、均勻性、異質性、相關性、方差、二階矩、熵、對比度8 種紋理特征，對紋理特征用最小噪聲分離法降維。使用提取的紋理信息和光譜信息進行SVM 監督分類，分為落葉松、紅松、樟子松、云杉、闊葉樹、耕地和其他7 類[24-25]。經2017、2019 和2020 年采集的679 個樣地調查數據驗證，總體分類精度為91.28%，Kappa 系數為0.88。

1.2.3 樣地實測數據 2017 年對孟家崗林場70 塊固定樣地進行每木檢尺，測量了樣地中每株樹的胸徑、樹高和冠幅，計算每個樣地的算術平均胸徑和算術平均高。采用差分GNSS 對樣地進行精確定位，定位精度優于1 m。70 個樣地中包括落葉松樣地40 個、樟子松樣地11 個和紅松樣地19 個。

1.2.4 人工區劃小班數據 2016 年以衛星影像為底圖，通過目視解譯得到林場林相圖，與航空遙感數據存在時間和空間位置上的差異。因此，以2017年機載高光譜影像為底圖，對2016 年林相圖進行修正，得到2017 年人工區劃小班數據。修正調整后的小班界線更符合實際情況，且保證了與機載高光譜和激光雷達數據空間位置的一致性。冠層高度模型、樹種類別圖、小班邊界和樣地分布見圖1。

圖1 孟家崗林場冠層高度模型、樹種類別圖及小班區劃圖Fig.1 Canopy height model,tree species map and related data for sub-compartments delineation in Mengjiagang forest farm

采用孟家崗林場2017—2019 年共有149 個伐區調查數據，面積為1～17 hm2。由于林分生長慢，對小班區劃而言，2017—2019 年的林分變化可以忽略不計，故采用3 個年度的主伐作業小班作為對照。

在森林資源調查監測工作中，參考對象主要有人工區劃的小班圖和比用于自動區劃的遙感數據分辨率更高的數據勾繪得到的小班圖，其中勾繪得到的小班圖誤差更少[4]。本研究基于0.1 m 空間分辨率DOM 屏幕勾繪了100 個小班作為參考小班。

2 研究方法

使用機載高光譜影像分類得到的樹種圖作為輔助數據，首先對1 m 空間分辨率的CHM 數據進行多尺度分割，得到過分割結果圖；然后對1 m 空間分辨率的CHM 數據均值下采樣到5 m，并使用對稱近鄰均值濾波（SNN）平滑，再對其進行多尺度分割，得到欠分割結果圖；最后將過分割和欠分割結果進行融合優化，得到最終小班區劃結果。總流程圖如圖2 所示。

圖2 基于機載激光雷達冠層高度模型和樹種類別的小班自動區劃流程Fig.2 The flowchart of automatic delineation of forest stands based on CHM from airborne LiDAR data and tree species

2.1 LiDAR 數據預處理

采用Lastools 軟件對LiDAR 點云數據進行預處理。去除噪聲點，進行點云分類，將分類后的地面點與非地面點分別采用TIN 內插算法與最大高程插值方法生成空間分辨率為1 m 的數字高程模型（DEM）和數字表面模型（DSM），用DSM與DEM 相減獲得CHM 數據，其空間分辨率為1 m[28]。

2.2 CHM 平滑濾波

CHM 空間分辨率小于研究區樹冠的平均冠幅，在圖像分割中將產生大量過分割的斑塊，為此，對1 m 空間分辨率的CHM 數據進行適當的重采樣和平滑濾波，對同一株樹的多個像元取平均值，以更好地顯示這株樹的信息。在林分尺度上，對CHM 數據適當的平滑不會對小班邊界產生大的影響。為減少像元重采樣和平滑對小班邊界的影響，確保邊界的突顯性，采用具有邊緣保持特性的平滑濾波器。

對KNN 均值、最小方差均值、Sigma 平滑和SNN4 種濾波器的5 m 空間分辨率CHM 數據平滑效果進行比較分析，綜合考慮異質區域間邊界的保持效果和同質區域內的平滑效果，選用SNN 濾波器進行5 m 分辨率CHM 的平滑處理。

SNN 濾波的核心思想是：比較窗口范圍內對稱像元的值，選取每對像元中與目標像元值更相近的值，對窗口內選取的所有像元取均值作為目標像元值，并保留與目標像元差異大的像元值不變。這樣在保持邊緣的同時，又可以以較低的計算量對圖像進行平滑[29]。根據地面調查數據，研究區成熟林分的平均冠幅約為5 m。比較3 × 3、5 × 5、7 ×7 和9 × 9 像素窗口的濾波效果，發現當窗口過小時會使平滑效果變差，最終區劃結果中產生過多面積較小的小班。窗口過大時不同小班的邊界更加模糊，難以區分不同小班，最終區劃結果中有很多差異較大的應區劃為不同小班的區域未分割開，窗口大小為5 × 5 時濾波結果最優。

2.3 多尺度分割

多尺度分割算法是從單像元對象開始的區域合并過程，合并規則基于局部光譜均勻性，即相鄰圖像對象之間的相似性。光譜均勻性標準定義為光譜值的標準偏差之和，每次合并時方差增加最小的對象被合并，當最小增量高于用戶定義的閾值(尺度參數)時，分割過程結束[30]。本研究中的過分割即分割出大量面積較小的對象，對象小于小班尺度。欠分割即分割出數量較少、面積較大的對象，對象接近小班尺度。其中面積很小的與周圍區域差異較大的區域沒有被分割出來，因此稱為欠分割。

將人工區劃小班疊加到樹種分類圖上，統計每個小班各樹種像元數占小班總像元數的比例，即樹種組成。選取樹種成數最大的優勢樹種作為該小班的樹種。融合相同樹種的區域形成樹種矢量圖。在eCognition 軟件中用CHM 數據作為底圖，加入樹種矢量圖作為輔助數據，在同類樹種區域內對CHM 進行多尺度分割得到林分高度和密度一致的小班斑塊。分割方案中尺度參數設定為10～300，形狀參數和緊致度參數設定為0.1～0.9，進行了多組試驗，最終確定了分割效果最優的參數。

2.4 分割結果優化

對CHM 數據進行多尺度分割的過程中，采用了多組參數進行試驗，結果顯示為1 m 和5 m 空間分辨率數據的分割結果最具互補性。5 m 空間分辨率CHM 的欠分割結果顯示了主要的林分區域，但其輪廓和林分邊緣的樹冠輪廓并不完全相符。1 m分辨率CHM 的過分割結果雖然形成了過多的斑塊，但其輪廓和林分邊緣的樹冠輪廓更相符。因此需要使用過分割結果來優化欠分割結果。

本研究使用ArcGIS 軟件中的update 工具實現多尺度分割結果的融合。其原理是對輸入圖層和用于更新的圖層進行幾何相交的計算，其中輸入圖層被更新圖層覆蓋部分的屬性用更新圖層的屬性代替。在更新時可以調整容差值，即輸出圖層邊界的X,Y坐標相對于輸入圖層邊界的X,Y坐標可移動的距離。將欠分割結果作為更新圖層，通過調整兩種尺度的分割結果間的容差值，經過多次試驗，確定容差值為3 m 時優化效果最佳。

2.5 區劃結果精度驗證

區劃結果采用形狀因子和定量因子2 種方式進行精度驗證。首先評估林分邊界的準確程度，采用人工區劃小班、主伐作業小班和0.1 m 空間分辨率DOM 屏幕勾繪的小班作為參考。根據自動區劃小班與參考小班的空間位置相對應的原則，采用最終測量精度準則，逐一測量自動區劃小班與參考小班的特征。最終測量精度準則評價的特征包括圓度（RO），緊致度（CO），形狀指數（SI），最小包絡圓短半徑（RE），橢圓度（EF）和形狀因子（P2A）[4]。

同時，使用交并比指標[31-32]，比較自動區劃小班和參考小班的重合情況。計算人工區劃小班與其對應的自動區劃小班的交集的面積與并集面積的比值[31]，即：

其次，采用樣地平均胸徑、平均樹高和5 m ×5 m CHM 的冠層平均高數據計算可解釋性方差精度，評價小班內部的均質性及不同小班之間的異質性[16-17,20]。即計算小班內部觀測值或冠層平均高的方差之和SSwithin，各樣地值與全研究區均值的方差之和SStotal的可解釋性方差R2，公式如下：

使用樣地數據計算時，k為小班數，ni為小班i中的樣地數，xij為小班i中樣地j的觀測值，為全研究區樣地觀測值的均值，為小班i中樣地觀測值的均值。使用CHM 數據計算時，ni為小班i中的柵格個數，xij為小班i中柵格j的冠層平均高度值，為全研究區冠層平均高的均值，為小班i中冠層平均高的均值。可解釋性方差越接近1，則小班內部一致性越高，不同小班間的差異性越大。

3 結果與分析

3.1 定性分析

經過多次試驗，對1 m 空間分辨率CHM 數據進行分割時，最優尺度參數為100，形狀參數為0.1，緊致度參數為0.5，得到過分割結果。接著對平滑后的5 m 空間分辨率CHM 數據進行分割，最優尺度參數為37，形狀參數為0.1，緊致度參數為0.5，得到欠分割結果。

最終區劃結果與過分割結果、欠分割結果和人工區劃小班的對比如圖3 所示。圖3-a 中顯示過分割結果劃分出了一些相對于小班尺度不必要分割開的區域。圖3-b 顯示欠分割結果沒有分割出相對于小班來說不必要的區域，但其邊界與小班邊緣的符合性較差。對比圖3-a 和3-b 可以看出欠分割結果對于小班區域的劃分效果比過分割結果好，但過分割結果與小班邊緣的符合性更好。優化后的最終區劃結果比欠分割結果與林分邊緣更相符，比過分割結果更好地區分了不同的林分區域。最終區劃結果與人工區劃小班的對照如圖3-c 所示，結果顯示最終區劃結果與人工區劃小班大體上近似，且最終區劃結果與林分邊界真實樹冠邊緣的一致性更高。

圖3 基于CHM 的最終區劃結果與欠分割結果、過分割結果、人工區劃小班對比Fig.3 Under-segmentation results,over-segmentation results,manual sub-compartments and automatic sub-compartments on CHM

3.2 定量分析

3.2.1 形狀因子一致性 自動區劃結果與參考小班邊界一致性的對比如表1 所示。表中各UMA特征從形狀和面積等多個角度反映了自動區劃結果邊界的準確程度，自動區劃結果與基于DOM 勾繪的小班、主伐作業小班和人工區劃的小班的各特征都較接近，相近程度依次遞增。當IoU＞0.7時，自動區劃結果與人工區劃小班相符性非常好，IoU＞0.5時，自動區劃結果與人工區劃結果相符性較好。

表1 自動區劃結果與參考小班的UMA 特征均值及lOU＞0.7、lOU＞0.5 的自動區劃小班占總數的百分比Table 1 The mean values of UMA features and the proportion of automatic sub-compartments with lOU＞0.7、lOU＞0.5compared with reference sub-compartments

由此可知自動區劃結果中46%、37%、43%與人工區劃小班、主伐作業小班、基于DOM 勾繪的小班相符性非常好，61%、54%、55%與人工區劃小班、主伐作業小班、基于DOM 勾繪的小班相符性較好。

交并比較低的小班主要對應以下3 種情況：①人工區劃小班邊界與地面真實小班邊界不符，而自動區劃結果邊界與真實邊界相符。②人工區劃小班人為地將同質的林分劃分為多個小班，而自動區劃時對這部分林分劃分的界線與人工區劃的不一致。③幼苗幼樹小班的樹高過低，接近耕地，CHM自動區劃方法沒有劃分開。

3.2.2 林分內部的均質性和外部異質性 自動區劃結果的平均胸徑可解釋性方差為0.97，平均樹高可解釋性方差為0.98，和人工區劃小班相同。可解釋性方差精度較高，說明自動區劃的小班和人工區劃的小班都具有內部一致性高，且不同小班之間區分度好的特點。自動區劃結果中小班與樣地個數的對應情況與人工區劃小班近似，并且二者的可解釋性方差一致，說明自動區劃方法產生了近似人工區劃結果的合理的小班。人工區劃小班的平均冠層高度的可解釋性方差為83.04%，自動區劃小班的平均冠層高度的可解釋性方差為84.81%。自動區劃小班比人工區劃小班的可解釋性方差提高了1.77%。說明自動區劃小班對不同林分的區分效果比人工區劃小班更優。

4 討論

本研究得到的全場小班區劃結果優于人工區劃結果，說明CHM 應用于小班區劃能夠提高邊界區劃精度，減少主觀性。高光譜圖像樹種分類結果圖使小班自動區劃結果包含了樹種屬性。此外，自動區劃結果的可解釋性方差較高，說明自動區劃結果符合小班內部特征基本一致，與相鄰小班有明顯區別的要求，而定性分析也顯示自動區劃的小班大部分符合實際情況。本研究方法在約15 000 hm2林場范圍的應用時間約為10 h，而人工區劃約1 周時間，工作效率顯著提高。

在交并比較低的3 種原因中，前2 種是由于人工區劃小班結果與地表真實小班的邊界不符造成的。這些交并比低的小班，不代表自動區劃不合適，而是表現出了自動區劃與真實值更加一致，不受主觀及其他因素影響的優勢。第3 種原因中，在人工區劃目視判讀時也會有部分幼苗幼樹和耕地難以區分開，這些區域需要進行實地調查才能確定。

不同于Dechesne 等[15]的研究需要計算大量的林分特征，本研究的方法可以直接利用高精度的CHM 和樹種圖進行自動區劃來得到和人工小班接近的區劃結果。Mustonen 等[13]計算的平均胸徑可解釋性方差為74%，平均樹高可解釋性方差為83%。Pukkala 等[16-18]計算的平均樹高、平均胸徑、斷面積可解釋性方差為66%～87%，最大樹高可解釋性方差為84.7%～94.2%。本研究自動區劃的小班的平均胸徑和平均樹高可解釋性方差精度相對于這些研究有所提高。本研究的方法在一個完整林場范圍內取得了內部一致性精度高且和人工區劃小班接近的結果。另外，其他相關區劃研究中沒有與人工區劃小班邊界對比，有的林分區劃結果雖然絕對精度較高，但圖斑較破碎，面積過小，與人工小班差異較大。本研究引入UMA 準則和交并比分析方法，定量計算了使用遙感數據自動區劃的小班與人工小班的相符性，并分析了其不同的原因，以及自動區劃小班的優勢和不足，有利于后續對自動分割小班更好的符合甚至優于人工區劃小班的進一步研究。

5 結論

本研究采用機載激光雷達數據提取的CHM 和高光譜影像產生的樹種類別圖進行了小班自動區劃的研究，主要結論如下：（1）本方法充分利用了樹高、林分密度、郁閉度等林分結構信息和高精度的樹種信息，區劃效率遠高于人工小班區劃；（2）自動區劃小班的內部一致性高，外部差異性大，在形狀、面積、交并比等多個指標上與人工區劃小班相符，自動區劃小班的邊界較人工區劃小班邊界更為精準；（3）本方法不需要計算大量特征，適用于林場及更大尺度的小班自動區劃。此外，采用機載遙感數據進行自動區劃能更好地應用于高精度的森林參數提取等工作，符合小班區劃精細化和森林質量精準提升的發展趨勢。