基于機載激光雷達的杉木林分高生長變化分析

余 鑄，李 春，黃 媛

（1.廣西壯族自治區(qū)林業(yè)勘測設計院，廣西南寧 530011；2.中國林業(yè)科學研究院熱帶林業(yè)實驗中心，廣西憑祥 532600；3.自然資源部北部灣經(jīng)濟區(qū)自然資源監(jiān)測評價工程技術創(chuàng)新中心，廣西南寧 530011）

杉木（Cunninghamialanceolata）為我國南方主要造林樹種，具有生長快、產(chǎn)量高和材質(zhì)佳等特點，被廣泛應用于生產(chǎn)和生活中[1]。對杉木林分生長進行實時精準監(jiān)測，有助于森林管理者從宏觀和空間尺度掌握杉木生長狀況，提高杉木人工林經(jīng)營管理水平。隨著全球定位系統(tǒng)（Global Positioning System，GPS）和慣性導航系統(tǒng)（Inertial Navigation System，INS）的集成應用，激光掃描儀實現(xiàn)了精確定位定姿，點云位置精度得到保障；機載激光雷達成為當前應用在森林資源調(diào)查和生態(tài)監(jiān)測中的先進技術手段[2-4]。

近20年來，許多學者基于雷達變量與林分平均高的數(shù)據(jù)統(tǒng)計關系，對不同森林類型進行研究[5-6]。大多數(shù)林分平均高均可通過模型進行反演得到[7-9]，主要有參數(shù)回歸模型和非參數(shù)模型；這些模型在特定研究區(qū)、特定森林類型的林分平均高估測方面均取得了良好的反演效果。但這些模型的構建均需大量的樣地實測數(shù)據(jù)。如何更快速、更高效地獲取林木樹高，實現(xiàn)對更廣闊空間尺度林分高生長量的連續(xù)監(jiān)測，提供更全面的信息，對于維護生態(tài)平衡、保護生態(tài)系統(tǒng)和可持續(xù)利用森林資源至關重要。本研究運用機載激光雷達數(shù)據(jù)，基于冠層高度模型，直接獲取杉木人工林內(nèi)不同小班間和小班內(nèi)部樹高生長水平和變異情況，可實現(xiàn)區(qū)域尺度上精準、高效地監(jiān)測杉木林分年高生長狀況，為杉木人工林管理和森林資源管理與監(jiān)測提供參考。

1 材料與方法

1.1 試驗區(qū)概況

在廣西壯族自治區(qū)國有高峰林場內(nèi)設置試驗區(qū)（108°23′E，22°58′N），海拔為90 ～460 m，相對高度為50 ～100 m；主要包括東升分場、界牌分場和延河分場。林地主脈呈東北至西南走向，形狀近似矩形，寬度約為4.2 km，長度約為11.2 km。該試驗區(qū)幾乎全部為人工林，主要樹種為巨尾桉（Eucalyptus grandis×E.urohpylla）、馬尾松（Pinusmassoniana）、濕地松（Pinuselliottii）、杉木、八角（llliciumverum）、肉桂（Cinnamomumcassia）、紅錐（Castanopsishystrix）、米老排（Mytilarialaosenisis）、火力楠（Mytilaria macclurei）和厚莢相思（Acaciacrassicarpa）。

1.2 數(shù)據(jù)預處理

1.2.1 數(shù)據(jù)獲取

2016 和2017 年10 月，使用R44 直升機和RIEGL（美國）有限公司的VUX-1LR 激光雷達系統(tǒng)進行數(shù)據(jù)采集。該激光雷達系統(tǒng)的工作波長為1 550 nm，激光發(fā)散角為0.5 mrad，飛行時的高度為2 500 m，獲得的平均點云密度約為每平方米2 點。同時，獲取試驗區(qū)航空數(shù)字正射影像（Digital Orthophoto Map，DOM）。

1.2.2 冠層高度模型獲取

林分冠層高度模型（Canopy Height Model，CHM）是一種表面模型，能反映植被與地面高度的距離信息，常用于間接提取冠幅、樹高、郁閉度和森林生物量等重要森林參數(shù)。本研究中，CHM 主要通過機載激光雷達點云數(shù)據(jù)生成的數(shù)字表面模型（Digital Surface Model，DSM）和數(shù)字高程模型（Digital Elevation Model，DEM）獲取。具體方法為提取2016 和2017 年的DSM 和DEM，計算DSM 和DEM 的差，得到2016和2017年試驗區(qū)域的CHM（圖1）。

圖1 林分冠層高度模型Fig.1 Canopy height models of stands

圖中，不同顏色反映林木冠層的分布情況。其中，綠色代表低冠層林分，這些區(qū)域范圍較廣，分布連續(xù)；紅色代表高冠層林分，這些區(qū)域范圍較窄，零散分布。綠色部分中存在一些代表CHM 為負值的區(qū)域；實際情況中，CHM 不存在負值。這種情況出現(xiàn)的原因為在匹配算法過程中，出現(xiàn)凹坑或空洞現(xiàn)象，導致無法完全密集匹配像點，出現(xiàn)地面控制點誤差，使得DSM 低于對應區(qū)域的DEM。采用“Con”函數(shù)，將CHM中小于0的值替換成0。

1.2.3 重采樣

通過重采樣獲得不同分辨率的CHM，進行優(yōu)化。重采樣主要包括3 種方法，分別為最鄰近法（Nearest Neighbor）、雙線性插值法（Bilinear Interpolation）和三次卷積法（Cubic Convolution）。本研究采用雙線性插值法，該方法適用于連續(xù)和非連續(xù)數(shù)據(jù)類型。

對2016 和2017 年的CHM 圖進行重采樣，分別設置10 m×10 m、20 m×20 m 和30 m×30 m 3種分辨率，比較不同分辨率下冠層高生長圖提取差異，選出最優(yōu)的重采樣分辨率（圖2）。

圖2 不同重采樣下的林分冠層高度模型Fig.2 Canopy height models of stands under different resampling

1.2.4 冠層高生長圖獲取

利用試驗區(qū)2016 和2017 年的CHM 獲取2016 — 2017 年冠層年高生長量（d_CHM），計算公式為：

d_CHM=2017_CHM-2016_CHM

重采樣時，采用不同分辨率CHM 一一對應，分別獲取不同分辨率下的d_CHM（圖3）。為更直觀地顯示試驗區(qū)林木高生長狀況，將林木高生長量分為若干個級別，分別用不同顏色表示，直觀、形象地將林木高生長分布情況展現(xiàn)出來。采用ArcGIS 軟件中的歸一化分類進行處理并輸出圖片。

圖3 不同分辨率下冠層年高生長量Fig.3 Annual height growths of canopies under different resolutions

圖中，黃色部分表示林分高生長量逐漸增加的區(qū)域，其變化趨勢符合杉木生長規(guī)律；紅色和綠色部分分別代表林分高生長量過大和過低的區(qū)域，其林木高生長量異常，在現(xiàn)實情況中不存在。出現(xiàn)這種現(xiàn)象是因為林木樹冠有弧度且不平滑，隨時間推移，樹冠結構不斷發(fā)生變化，2017年CHM 與2016年CHM 中的林木位置不能一一對應，d_CHM 出現(xiàn)極低值和極高值。還可能存在其他原因，如地形起伏過大、陰影遮擋、部分地面信息缺失和密集匹配結果失衡等。

1.3 杉木林分年高生長量提取

利用小班矢量圖層、冠層年高生長量圖層和隨機點圖層進行杉木林分年高生長量提取。結合高峰林場小班資源數(shù)據(jù)和數(shù)字正射圖像，通過地物陰影、形狀、顏色、色調(diào)、大小、位置、紋理、空間關系和圖案等確定林分優(yōu)勢樹種和林齡[10]，最終勾繪出30個杉木林小班，按照“之”字編號原則對各小班進行編號。采用ArcGIS軟件中的創(chuàng)建隨機點工具創(chuàng)建1個隨機點圖層，設置不同隨機點數(shù)，得到隨機點圖層。將小班矢量圖層中的小班號字段和2017 年高峰林場小班圖層的林齡字段融合到隨機點圖層中。將冠層年高生長量圖層與隨機點圖層疊加，將林分年高生長量參數(shù)融合到隨機點圖層中。在杉木矢量圖層中，分別設置100、200、300 和400 個隨機點，并提取不同分辨下d_CHM，設計12個處理（表1）。

表1 不同處理Tab.1 Different treatments

將杉木林分年高生長量數(shù)據(jù)中的負值和過大值等不符合林木實際生長情況的異常數(shù)據(jù)剔除，剔除數(shù)據(jù)占比控制在5%以內(nèi)。統(tǒng)計各小班均值、標準差和變異系數(shù)，并進行方差分析，檢驗小班間年高生長量參數(shù)的差異顯著性。同時，計算不同林齡杉木林分在不同分辨率和提取隨機點數(shù)處理下的均值。上述過程均在Python 3.7軟件中進行。

2 結果與分析

2.1 不同處理下提取的杉木林分年高生長量分析

一般杉木人工林的樹高年生長量為0.5 m[11]。T7 處理的杉木林分年高生長量均值最?。?.41 m）；T11 處理的杉木林分年高生長量均值最大（0.62 m）（表2）。大部分處理獲得的杉木林分年高生長量均值在（0.50±0.10）m內(nèi)波動，符合實際杉木人工林樹高生長規(guī)律。少量杉木林分樹高生長量均值為4 ～5 m，數(shù)據(jù)離散程度和差異較大，這是由于不同年份的CHM 不能完全匹配，兩個年份數(shù)據(jù)相減出現(xiàn)正、負值。不同處理獲得的樹高生長量均呈正態(tài)分布，表明方差分析結果可靠（圖4）。

表2 不同處理下杉木林分年高年生長量統(tǒng)計Tab.2 Statistics on annual height growths of C.lanceolata stands under different treatments

圖4 不同處理下杉木林分樹高生長量正態(tài)分布Fig.4 Positive distributions of height growths of C.lanceolata stands under different treatments

異常率為均值出現(xiàn)負值的小班數(shù)除以小班總數(shù)的比例，用于反映使用該處理提取數(shù)據(jù)時出現(xiàn)異常情況的概率，其值越小，數(shù)據(jù)準確度越高。T9 處理的異常率最大（40.00%），T4 處理的異常率最?。?3.33%），即重采樣10 m×10 m 處理下取400 個隨機點數(shù)時獲取的數(shù)據(jù)更精確，重采樣30 m×30 m 處理下取100個隨機點數(shù)時獲取的數(shù)據(jù)精確度不高。

2.2 小班年高生長量分析

對T4 處理下不同杉木林分小班內(nèi)林木年高生長量進行統(tǒng)計；其中，1、3、6、9、10、13和18號小班的年高生長量為負數(shù)，與常規(guī)針葉人工林實際生長情況不符；剔除不符實際情況的數(shù)據(jù)后，計算試驗區(qū)杉木林分年高生長量（表3）。結果顯示，杉木林分年高生長量均值為0.80 m；說明試驗區(qū)杉木林分小班內(nèi)年高生長量總體處于較優(yōu)良范圍，樹高生長狀況良好，生長較旺盛。

表3 T4處理下杉木林分不同小班年高生長量統(tǒng)計Tab.3 Statistics on annual height growths of C.lanceolata stands in different subcompartments under T4 treatment

剔除異常數(shù)據(jù)后，試驗區(qū)杉木林分年高生長量標準差均值為4.44 m，各小班標準差為2.16 ～6.28 m；各小班年高生長量變異系數(shù)為0.57%～187.34%，即不同小班年高生長量離散程度較大；11 個小班年高生長量均值超過總體均值，小班生長狀況相差較大。小班內(nèi)杉木樹高生長不穩(wěn)定，可能是受到地形地勢、經(jīng)營管理活動等的影響；具體變動原因需進一步根據(jù)小班地類、立地類型和樹種等相關因子進行分析，開展外業(yè)核查。

2.3 不同林齡杉木林分年高生長量分析

計算不同林齡杉木林分年高生長量均值。排除明顯不符合生長規(guī)律的異常數(shù)值后，杉木林分年高生長量大致呈隨林齡增加先增加后減少的趨勢（表4，圖5）?？傮w上，在12 ～15 年，年高生長量隨林齡增加而增加；15 ～20 年，年高生長量隨林齡增加呈減少趨勢；20 ～27 年，年高生長量波動變化。說明12 ～15年為杉木樹高快速增長期，林木處于發(fā)育階段，生長空間大；隨林齡增加，林木間競爭加大，林木生長減緩；總體規(guī)律符合不同林齡下杉木林分高生長規(guī)律。

表4 不同林齡杉木林分年高生長量統(tǒng)計Tab.4 Statistics on annual height growths of C.lanceolata stands at different ages(m)

圖5 不同處理下杉木林分年高生長量Fig.5 Annual height growths of C.lanceolata stands under different treatments

3 討論與結論

目前，大多數(shù)林分高生長可通過模型進行反演，但模型在不同傳感器獲取的數(shù)據(jù)中普適性較低，其構建仍需大量實地調(diào)查樣本，耗費人力、物力，難以實時掌握林分動態(tài)變化，進行大面積林分高生長監(jiān)測。本研究通過直升飛機搭載的機載激光雷達系統(tǒng)獲取連續(xù)兩年的點云數(shù)據(jù)，運用雙線性插值法提取冠層高度模型，直接計算兩期冠層高度模型差，獲取杉木林分年高生長量，分析杉木生長趨勢，可實現(xiàn)高效率、低成本的實時林分監(jiān)測。

本研究設置不同分辨率重采樣與提取高度的隨機點數(shù)，比較它們對最終數(shù)據(jù)獲取的影響。結果顯示，不同設置下獲取的數(shù)據(jù)間差異不大，但相同設置下各小班的年高生長量變異系數(shù)普遍較大，數(shù)據(jù)存在一定誤差。Ben-arie 等[12]利用高效的半自動坑填充方法，統(tǒng)計的平均樹高差異百分比為0.75%。劉一成等[13]基于曼哈頓距離并結合形態(tài)學冠層控制的方法去除CHM 凹坑，提取平均樹高，其更接近實際值，差異百分比為0.59%。白明雄等[14]研究表明，直接使用CHM 表面模型，會丟失林木陰影、冠幅等細節(jié)，導致最終獲取的數(shù)據(jù)存在較大誤差。本研究采取將負值賦值為0 的方法去除CHM 負值，不能解決匹配算法在密集匹配像點時出現(xiàn)凹坑現(xiàn)象的問題，提取的林分年高生長量數(shù)據(jù)存在誤差。填補凹坑和空洞是否能真正提高提取數(shù)據(jù)精確度，減小數(shù)據(jù)差異，有待進一步研究。

汪垚等[15]研究表明，一般針葉人工林的年高生長量為0.5 m 左右。楊傳桃[16]研究顯示，樹高的年生長量一般規(guī)律為先增加后減少，最終趨于穩(wěn)定。本研究結果顯示，杉木人工林的年高生長量為0.80 m，其與林齡的關系大致符合先增后減的規(guī)律，但后期出現(xiàn)一定波動，可能是由于試驗數(shù)據(jù)提取較少，選取的林齡段有限，后續(xù)需進行更廣泛的研究。

總體而言，機載激光雷達技術在全域樹高的直接精度測量方面具有優(yōu)勢[17]，同時能快速獲取樹冠結構和高度信息[18]，實現(xiàn)數(shù)據(jù)化和圖像化的直觀呈現(xiàn)，操作便捷，精度高，受時間和空間影響小，可減少監(jiān)測工作的人力投入。此外，該技術還具有精準性好、獲取速度快等特點，能實現(xiàn)林分大面積數(shù)據(jù)采集，有利于森林資源管理與監(jiān)測。

采用機載激光雷達數(shù)據(jù)獲取冠層高生長量模型，并結合遙感數(shù)字正射影像，可間接提取杉木林分年高生長量數(shù)據(jù)，監(jiān)測杉木林分年高生長量；在重采樣10 m × 10 m 分辨率且隨機點數(shù)為400 的情況下，數(shù)據(jù)準性更高。利用機載激光雷達數(shù)據(jù)提取的杉木人工林林分年高生長量數(shù)據(jù)符合林木本身生長規(guī)律，具有現(xiàn)實意義。

利益沖突：所有作者聲明無利益沖突。

作者貢獻聲明：余鑄負責研究計劃制定、數(shù)據(jù)分析和論文撰寫與修改；李春負責數(shù)據(jù)分析、文獻檢索和論文修改；黃媛負責數(shù)據(jù)收集和圖片制作。