張 瓊 劉 芳 范文義 李 典 楊樹文 陳 成

(東北林業(yè)大學(xué)，哈爾濱，150040)

激光雷達(dá)LIDAR(Light detection and ranging)用于估測森林參數(shù)的研究始于上世紀(jì)80年代中期［1－3］，它最大的優(yōu)勢之一就是能直接獲取目標(biāo)物體的三維坐標(biāo)信息，這對于定量估測森林參數(shù)，尤其在估測森林冠層高度及林木空間結(jié)構(gòu)信息方面具有得天獨厚的優(yōu)勢，因此與傳統(tǒng)的采用測高器量取樹高和利用樹冠垂直投影量測冠幅的方法相比，可信度大大增加。高精度的森林冠層高度信息的獲取一直是林業(yè)遙感中難以解決的問題，而LIDAR憑借其極高的角分辨率能力、距離分辨率能力、抗干擾能力，能夠快速準(zhǔn)確地獲取地表物體的高度信息［4］，進(jìn)而獲取林地數(shù)字高程模型。本研究旨在探討利用高密度小光斑激光雷達(dá)和同步獲取的大比例尺航片數(shù)據(jù)提取我國小興安嶺地區(qū)林木信息的可用性，具體包括:(1)研究單株木樹高和樹冠的提取技術(shù);(2)將提取出的林木參數(shù)與樣地實測林木參數(shù)進(jìn)行相關(guān)性分析，建立回歸模型反演林木參數(shù)，并分析影響林木參數(shù)提取精度的因素。

1 研究區(qū)域概況

研究區(qū)選在黑龍江省伊春市帶嶺區(qū)的涼水國家級自然保護(hù)區(qū)，地理坐標(biāo)為北緯47°6＇49″～47°16＇10″，東經(jīng)128°47＇8″～128°57＇19″。保護(hù)區(qū)總面積 12 133 hm2，森林覆被率 98% ，屬小興安嶺南部達(dá)里帶嶺支脈東坡，海拔高度280～707 m。本區(qū)的主要保護(hù)對象是以紅松(Pinus koraiensis)為主的溫帶針闊葉混交林生態(tài)系統(tǒng)，區(qū)內(nèi)既有處于演替頂極狀態(tài)的原始闊葉紅松林、冷杉(Abies ziyuanensis)林、云杉(P.kor-aiensis)林和興安落葉松(Larix gmelini Rupr.)林，又有處于不同演替階段的次生林，幾乎囊括了小興安嶺山脈的所有森林植被類型。

2 研究方法

2.1 樣地數(shù)據(jù)獲取

涼水自然保護(hù)區(qū)的外業(yè)數(shù)據(jù)采集時間為2010年7月28日—8月6日，在研究區(qū)內(nèi)共選擇了24塊10 m×10 m的正方形樣地，并使用差分GPS定位樣地中心地理坐標(biāo)。根據(jù)機載數(shù)據(jù)覆蓋情況，選取試驗區(qū)的代表樹種，選擇被測單株木的標(biāo)準(zhǔn)是該樹冠的大部分沒有被其他樹冠遮擋。采用胸徑尺測量胸徑，手持超聲波測高儀測量單株木樹高，皮尺測量東西和南北的樹冠直徑，利用差分GPS對樣地中的每株樹進(jìn)行定位。選擇的大多數(shù)樣地為天然的針闊葉混交林。

2.2 LIDAR數(shù)據(jù)獲取及預(yù)處理

研究中采用的激光雷達(dá)數(shù)據(jù)由LiteMapper—5600系統(tǒng)獲取，激光掃描儀為Riegl LMS—Q560，數(shù)據(jù)獲取時間為2009年9月。航飛相對高度為800 m，激光脈沖長度是3.5 ns，回波寬度分辨率0.15 m，激光離散角度是0.5 mrad，垂直精度可達(dá)0.15 m，采樣間隔為1 ns，采樣的點云密度大于2點/m2，平均對地飛行速度是180 km/h。

LIDAR數(shù)據(jù)每個激光點包含了激光點三維坐標(biāo)值、強度值、回波類型等信息。原始點云數(shù)據(jù)會包含大量粗差、錯誤或無關(guān)信息。在前期數(shù)據(jù)處理過程中，采用人工交互編輯的方法，利用比較成熟的商業(yè)LIDAR數(shù)據(jù)處理軟件Terrasolid進(jìn)行預(yù)處理，剔除錯誤信息并補足遺漏信息，以保證后期生成冠層高度模型的精確性。首先，用分類軟件Terrasolid中的Terrascan分類模塊對LIDAR數(shù)據(jù)進(jìn)行粗分類，將LIDAR數(shù)據(jù)分成純地面激光回波點和非地面激光回波點，再將分類后的純地面激光回波點使用TIN插值方法生成數(shù)字高程模型(Digital Elevation Model，簡稱 DEM)，同樣，將非地面激光回波點內(nèi)插生成數(shù)字表面模型(digital surface model，簡稱DSM)。運用ERDAS IMAGINEG軟件將插值生成的DSM和DEM模型作相減運算，直接獲得研究區(qū)的正規(guī)化數(shù)字表面模型nDSM。對于森林地區(qū)而言，nDSM也等同于冠層高度模型CHM(Canopy Heigh Model，簡稱CHM)。CHM消除了地形起伏變化對DSM中樹木高度及形狀的干擾，可以獲得相對準(zhǔn)確的樹木形狀信息，被定義為樹冠在水平表面上的水平和垂直分布模式［5］。以1塊1 000 m×1 000 m的試驗區(qū)數(shù)據(jù)(空間分辨率為0.2 m)為例，如圖1為由LIDAR數(shù)據(jù)處理得到的nDSM，圖像上由暗到亮表示高度值由小到大。

圖1 試驗區(qū)的nDSM

此外，在飛行時還同步獲取了空間分辨率為0.2 m的CCD影像共189幅，比例尺為1∶2 000，總覆蓋面積約為1.2 km×1.6 km。利用ERDAS Imagine 9.1軟件中的LPS(數(shù)字?jǐn)z影測量系統(tǒng))模塊，導(dǎo)入相機文件及內(nèi)外方位元素等航攝參數(shù)，同時，基于雷達(dá)點云數(shù)據(jù)獲取的DEM對該大比例尺航片進(jìn)行正射糾正，得到空間分辨率為0.2 m的正射影像(Digital Orthophoto Map，即DOM)，圖2為覆蓋23號樣地的正射影像。

2.3 單株木樹高提取

從圖1可以看出，樹冠頂部比樹冠的其他區(qū)域要亮一些，樹冠之間的間隙區(qū)域最暗。因此，假設(shè)樹冠頂部高于樹冠的其他區(qū)域，通過計算局部最大值，可以得到樹冠的頂點［6］。本文利用eCognition軟件對空間分辨率為0.2 m的正射影像進(jìn)行分割，識別出單株木樹冠區(qū)域。將輸出的樹冠多邊形矢量文件與CHM疊加，利用ERDAS的GIS空間區(qū)域分析功能，確定出單株木樹冠區(qū)域的局部最大值，得到樹冠頂點的高度信息，此即為利用雷達(dá)數(shù)據(jù)提取的單株木樹高。單株木樹高提取流程圖見圖3。

圖2 覆蓋23號樣地的正射影像

圖3 單株木樹高提取流程圖

2.4 單株木樹冠識別

采用黑龍江省涼水國家自然保護(hù)區(qū)部分試驗區(qū)的大比例尺航片作為試驗對象進(jìn)行樹冠提取。利用面向?qū)ο蟮膱D像分割方法提取單株木水平信息。首先將植被信息完整地提取出來，并將其作為一個整體多邊形為下一步的樹冠信息提取做準(zhǔn)備。因此，先對影像進(jìn)行大尺度分割，然后用基于樣本的分類方法進(jìn)行植被信息提取，并將植被區(qū)對象多邊形進(jìn)行基于類的融合，即將所有被分為植被的影像對象重新合并為一個對象，再對合并后的植被區(qū)影像對象進(jìn)行小尺度分割，提取出單株木樹冠?？紤]到單木分割精度，本文對試驗區(qū)10個樣地對應(yīng)的正射影像進(jìn)行研究。以23號樣地為例，首先將影像分割形成“粗”的影像對象，即將影像分割結(jié)果形成多邊形;再對影像對象分類，建立類別，然后分別對各個類別選取2～3個訓(xùn)練樣本，利用最大似然分類法分離出植被信息，并對漏分或錯分的樹冠多邊形進(jìn)行手動編輯，得到研究區(qū)的分類結(jié)果。最后，對植被多邊形進(jìn)行基于類的合并(merge)，將所分出的全部植被的影像重新合并為植被信息多邊形，再對合并后的植被多邊形進(jìn)行二次分割，最終得到精度較高的單株木樹冠分割結(jié)果，如圖4。圖4為單株樹冠分割出的多邊形，從圖上可見分割出的單株木樹冠大多比較完整和準(zhǔn)確。選擇單株木樹冠，輸出其矢量格式文件。通過對各項參數(shù)進(jìn)行組合試驗、篩選，最終選出適合本研究區(qū)的最佳分割參數(shù)組合，見表1。

表1 多尺度分割最佳參數(shù)組合

圖4 單株木樹冠分割圖

3 結(jié)果與分析

3.1 單株木冠幅反演

把分割出的單株木樹冠以多邊形矢量格式輸出，并計算單株木冠幅的大小。冠幅計算的方法有幾種，如:用與樹冠面積相同的直徑作為樹冠冠幅［7－9］;用樹冠4個主方向的半徑的平方的2倍作為樹冠的冠幅［10］等。本研究采用與面積相同的圓的直徑作為樹冠冠幅的方法。將利用差分GPS定位的實地樣木坐標(biāo)值生成矢量文件，并與研究區(qū)的正射影像疊加，之后進(jìn)行目視對照，精準(zhǔn)定位每株樹，并輸出單株木樹冠的多邊形面積，再求出其直徑，即為單株木樹冠冠幅。

由于單木分割時并不能把每株樹的冠幅分割出來，且航片拍攝期間太陽照射產(chǎn)生的陰影、研究區(qū)的地形以及樹冠與樹冠之間的疊合等因素都能影響分割的精度，造成一些地區(qū)的冠幅被錯分為其他類別。因此，研究中從10塊樣地中選取了96株能代表該地區(qū)林分情況且分割效果較佳的樹木，建立估測樹冠與實測樹冠的相關(guān)關(guān)系。圖5為估測樹冠與實測樹冠的散點圖，二者相關(guān)系數(shù)R2為0.8925，單株木冠幅估測精度都高于81%，平均精度高達(dá)90%，由圖5可見，單株樹冠估測的精度較高。

圖5 實測冠幅與估測冠幅的散點圖

3.2 單株木樹高反演

3.2.1 提取單株木樹高

樹高即樹冠的最高點到地面的距離。將輸出的單株木樹冠多邊形矢量文件與植被區(qū)的冠層高度模型(CHM)疊加(圖6)，然后利用ERDAS IMAGING中的GIS區(qū)域空間分析功能，計算出單株木樹冠多邊形內(nèi)CHM模型上的最大高程差值，即得到利用LIDAR數(shù)據(jù)提取出的單株木樹高。

圖6 單株木樹冠多邊形與CHM疊加圖

3.2.2 反演單株木樹高

以外業(yè)得到的實測樹高為因變量、以通過LIDAR數(shù)據(jù)提取出的單株木樹高為自變量進(jìn)行相關(guān)分析，繪制散點圖(圖7)，計算得出二者的相關(guān)系數(shù)R2為0.936 1，可以看出利用LIDAR測得的單株木樹高與實測樹高線性相關(guān)。

圖7 LIDAR單株木樹高與實測樹高的比較

根據(jù)LIDAR數(shù)據(jù)的離散特點和外業(yè)測量的實際條件，為避免真實結(jié)果受到影響，本研究結(jié)果分析中不考慮利用LIDAR提取出的單株木樹高與實測樹高的差值大于3.5 m的單株木(共14株)。從10個樣地中選擇82株樹進(jìn)行相關(guān)分析，如圖7。

以62株樹作為訓(xùn)練樣本、另外20株作為檢驗樣本進(jìn)行樹高反演，62株訓(xùn)練樣本與實測樹高的相關(guān)系數(shù)R2為0.930 5，建立的線性回歸方程為:

利用方程(1)對20株檢驗樣本進(jìn)行樹高反演，反演的樹高與實測的樹高比較的散點圖如圖8所示。

圖8 LIDAR反演樹高與實測樹高的比較

3.2.3 樹高反演精度評價

利用建立的LIDAR樹高反演方程(1)，對20株檢驗樣本進(jìn)行樹高反演，結(jié)果見表2。

表2 單株木樹高反演結(jié)果與精度評價

表2中，精度計算公式為:

由表2可以看出:估測樹高低于實測樹高的單株木占多數(shù)，估測精度最低的為86%，最高達(dá)到98%。20株檢驗樣本的平均精度為96%。其中:估測精度低于90%的有5株，占總數(shù)的25%;估測誤差大于3 m的樹木僅有1株;誤差小于1 m的樹木有4株;平均估測誤差為0.85 m。70%的估測樹高均低于實測樹高，這可能是由于機載激光雷達(dá)系統(tǒng)采樣時有一定間隔，容易錯失樹頂，從而造成樹高的低估。

分析表2的精度結(jié)果:20株單株木中闊葉樹8株，精度均高于95%，平均精度為96%;針葉樹共12株，精度都高于86%，平均精度91%。由此可見，利用LIDAR數(shù)據(jù)反演單木樹高時，闊葉樹精度要高于針葉樹，其主要原因是大部分激光點都不會直接來自樹頂?shù)姆瓷?，闊葉樹的樹冠形狀一般呈橢球形，而針葉樹的樹冠形狀一般呈圓錐形或鉛筆形，因此，闊葉樹的反射率更高，其樹冠的垂直結(jié)構(gòu)上升得更快，在LIDAR點云采樣中碰到樹頂或上層樹冠的概率更大。

4 結(jié)論與展望

利用TerraSolid專業(yè)激光雷達(dá)處理軟件對涼水國家級自然保護(hù)區(qū)的高密度機載LIDAR數(shù)據(jù)進(jìn)行處理，結(jié)合同步獲取的大比例尺航片反演單株木冠幅和樹高，結(jié)果表明:利用高密度機載小光斑激光雷達(dá)數(shù)據(jù)結(jié)合同步獲取的大比例尺航片提取我國小興安嶺地區(qū)的林木信息是可行的。研究采用面向?qū)ο蟮亩喑叨确指罘椒▽崿F(xiàn)了單株木樹冠的提取，結(jié)合高密度LIDAR數(shù)據(jù)提取出單株木樹高，并進(jìn)行單株木的樹高和樹冠的反演試驗，雖然精度相對較高，但目前這一技術(shù)的應(yīng)用仍具有局限性。主要表現(xiàn)在本研究所選樣地郁閉度均較低(＜0.8)，森林結(jié)構(gòu)比較簡單，因而樹冠提取的效果較好，而對于郁閉度更高、森林結(jié)構(gòu)更復(fù)雜的林分，可望提出一種更新、更精確的分割算法來自動識別單株木。激光雷達(dá)的高度信息能較準(zhǔn)確地估測樹高，今后可以對其它地區(qū)進(jìn)行此類樹高估測研究，如果估測精度較高且結(jié)果穩(wěn)定，那么此激光雷達(dá)變量提取方法用于其它地區(qū)的樹高估計是相當(dāng)有前景的。因此，隨著LIDAR技術(shù)的不斷發(fā)展和數(shù)據(jù)處理算法的不斷完善，對于高密度、多回波的激光雷達(dá)點云數(shù)據(jù)的進(jìn)一步研究和應(yīng)用，將獲得更為理想的結(jié)果。

［1］Arp H，Griesbach J，Burns J.Mapping in tropical forests:A new approach using the laser APR［J］.Photogrammetric Engineering ＆Remote Sensing，1982，48(1):91 －100.

［2］Nelson R，Krabill W，Maclean G.Determining forest canopy characteristics using airborne laser data［J］.Remote Sensing of Environment，1984，15(3):201 －212.

［3］MacLean G A，Krabill W B.Gross－merchantable timber volume estimation using an airborne LiDAR system［J］.Canadian Journal of Remote Sensing，1986，12:7 －18.

［4］楊伯鋼，馮仲科，羅旭，等.LIDAR技術(shù)在樹高測量上的應(yīng)用于精度分析［J］.北京林業(yè)大學(xué)學(xué)報，2007，29(增刊2):78 －81.

［5］Koukoulas S，Blackbum G A.Quantifying the spatial properties of forest canopy gaps using LIDAR imagery and GIS［J］.International Journal of Remote Sensing，2004，25(15):3049 －3071.

［6］劉清旺，李增元，陳爾學(xué)，等.利用機載激光雷達(dá)數(shù)據(jù)提取單株木樹高和冠幅［J］.北京林業(yè)大學(xué)學(xué)報，2008，30(6):83 －88.

［7］Persson˙A，Holmgren J，S?dermann U.Detecting and measuring individual trees using an airborne laser scanner［J］.Photogr Engineering and Remote Sensing，2002，68(9):925 －932.

［8］Hyyppa J，Kelle O，Lehikoinen M，et al.A segmentation－based method to retrieve stem volume estimates from 3－D tree height models produced by Laser scanners［J］.IEEE Transactions on Geoscience and Remote Sensing，2001，39(5):969 －975.

［9］F Morsdor F，Meier E，Kotz B.LIDAR－based geometric reconstruction of boreal type forest stands at single tree level for forest and wildland fire management［J］.Remote Sensing of Environment，2004，92(3):353 －362.

［10］Solberg S，Naesset E，Bollandsas O M.Single Tree Segmentation Using Airborne Laser Scanner Data in a Structurally Heterogeneous Spruce Forest［J］.Photogrammetric Engineering ＆ Remote Sensing，2006，72(12):1369 －1378.