邢艷秋 ,張錦繡 ,陳世培 ,高 立 ,關(guān) 雷 ,郭慧宇
(1.東北林業(yè)大學(xué) 森林作業(yè)與環(huán)境研究中心,黑龍江 哈爾濱 150040;2.黑龍江第三測繪工程院,黑龍江 哈爾濱 150025)
林分平均樹高,作為重要的森林測樹因子,既能反應(yīng)林分的生長狀況,又能為森林資源的經(jīng)營與管理提供基礎(chǔ)數(shù)據(jù)。隨著航空攝影測量技術(shù)與遙感技術(shù)的不斷革新,新的技術(shù)手段和遙感數(shù)據(jù)在森林資源監(jiān)測應(yīng)用中廣泛滲透[1-3]。
資源三號(ZY-3)衛(wèi)星作為我國首顆對地進行立體攝影測量的民用測繪衛(wèi)星,具有地面測高能力,目前主要用于1∶50 000地形測繪。ZY-3在無控制點情況下,定位精度優(yōu)于15 m,少量控制點的情況下高程精度優(yōu)于3 m,平面精度優(yōu)于4 m[4]。將ZY-3立體像對數(shù)據(jù)應(yīng)用于林業(yè)研究可以獲得大區(qū)域連續(xù)的森林冠層數(shù)字表面模型(digital surface model, DSM),能夠有效避免離散點因空間插值形成連續(xù)DSM時所引入的誤差。與此同時,ZY-3攜帶的多光譜信息可以用于輔助反映森林植被信息[5]。目前,已有研究人員利用攝影測量手段來提取樹高,該類樹高提取主要針對單木,且數(shù)據(jù)源多為航片。如:胥利紅利用航片進行解析空中三角測量,完成相對定向和絕對定向,配合影像判讀森林區(qū)域,所得樹高估測模型的決定系數(shù)為0.903的樹高估測模型,模型精度為87.22%[6]。王佳等以航片為數(shù)據(jù)源建立可量測的立體模型,進而從模型中提取單木樹高,最終所得模型估測精度為93.96%[7]。李明澤等同樣以航片為數(shù)據(jù)源,利用數(shù)字攝影測量平臺對立木高度進行量測,結(jié)果表明,此方法的估測精度為90.92%[8]。雖然上述研究結(jié)果表明,基于遙感航片數(shù)據(jù)能夠精確提取森林立木樹高,但目前基于ZY-3航空立體影像,且在無控制點的條件下大范圍獲取林分平均樹高的研究相對少見。這可能主要因為在森林冠層相對較密的林區(qū),被動光學(xué)信號較難穿透森林冠層到達地面,因此不能獲得林下高精度的地形數(shù)據(jù)[9]。
機載激光雷達(light detection and ranging,LiDAR)作為一種新型主動遙感技術(shù),既能獲取地物的水平結(jié)構(gòu)信息,同時又能獲得地物的垂直結(jié)構(gòu)信息,因而可以實現(xiàn)區(qū)域地形的高精度模擬[10]。如:吳芳等利用機載LiDAR點云數(shù)據(jù),通過對植被點數(shù)據(jù)進行分類移除,最終對植被覆蓋相對較高的礦山地面起伏情況進行模擬,生成數(shù)字高程模型(digital elevation model, DEM),結(jié)果DEM高程誤差為0.015 m[11]。尹艷豹等通對激光點云回波數(shù)據(jù)判別分析進而提取林地DEM,驗證得DEM高差總平均值為0.24 m,方差為0.05 m[12]。上述研究結(jié)果表明,機載LiDAR點云數(shù)據(jù)可以生成高精度DEM,彌補ZY-3立體影像在植被覆蓋區(qū)域提取DEM的不足。同時利用高精度的DEM作為高程基準(zhǔn)采集虛擬控制點進行平差能提高ZY-3立體影像數(shù)據(jù)的定位精度[13]。
因此,本研究通過聯(lián)合ZY-3數(shù)據(jù)與機載LiDAR點云數(shù)據(jù),以期實現(xiàn)無實測控制點條件下森林林分平均樹高的高精度估測,為我國測繪衛(wèi)星數(shù)據(jù)在林業(yè)結(jié)構(gòu)參數(shù)估測研究方面提供基礎(chǔ)數(shù)據(jù)及技術(shù)支撐。
本研究實驗區(qū)位于吉林省長春市凈月潭國家森林公園,其地處長白山余脈丘陵地帶(125°24′~125°28′E,43°45′~43°48′N),高程為175~406 m,森林由人工林與天然次生林組成,總面積達92.52 km2。樹種以落葉松Larix gmelinii、樟子松Pinus sylvestris和蒙古櫟Quercus mongolica為主。氣候為溫帶半濕潤季風(fēng)氣候,區(qū)域全年溫差約50 ℃,全年平均降水量約654.3 mm。
2012年10月對研究區(qū)進行了野外數(shù)據(jù)采集。根據(jù)林分的密度不同,在進行野外數(shù)據(jù)采集時布設(shè)了半徑分別為10 m和15 m的圓形樣地。利用Trimble GeoXT6000 GPS對樣地中心以及樣地內(nèi)所有活立木進行定位,定位精度約0.5 m。利用超聲波測高儀對樣地內(nèi)所有單木樹高使用進行3次測量,取3次量測均值作為單木樹高的實測值。同時利用胸徑尺對樣地內(nèi)的所有單木胸徑進行3次測量,同樣取測量的均值作為胸徑的實測值。最終在研究區(qū)共采集37個樣方(如圖1所示),其中包括34個針葉林(23個樟子松和11個落葉松)和3個闊葉林(2個蒙古櫟和1個楊樹)。
因遙感林分平均樹高估測研究中多采用胸高斷面積加權(quán)平均高[14],因此本研究同樣采用胸高斷面積加權(quán)平均高作為林分平均樹高,其具體計算公式如式(1)所示:
式中,HL為第L(L=1,2,…,n)個樣地的實測樹高,hi為第i株立木的樹高,gi為第i株立木的胸高斷面積,n為樣地內(nèi)立木株數(shù)。
ZY-3為太陽同步衛(wèi)星,高度約505km,回歸周期59 d,重訪周期5 d。本研究使用的ZY-3影像數(shù)據(jù)于2013年5月22日采集得到。
三線陣立體影像由3臺具有立體視角的光學(xué)相機沿衛(wèi)星飛行方向推掃成像,分別是前視22°和后視22°全色延時積分成像(TDI CCD)相機(分辨率3.5 m,幅寬52 km),下視全色TDI CCD相機(分辨率2.1 m,幅寬51 km)。譜段范圍0.5~0.8 μm,側(cè)擺能力±32°。
多光譜影像由一臺分辨率5.8 m,幅寬51 km, 擺 能 力 ±32°, 含有藍(0.45~0.52 μm)、綠(0.52~ 0.59 μm)、紅(0.63~0.69 μm)和近紅外(0.77~0.89 μm)四個通道的多光譜相機采集。
圖 1 研究區(qū)域Fig.1 Study area
本研究使用的機載LiDAR離散點云數(shù)據(jù)是由Leica ALS70型傳感器于2012年5月31日采集得到。傳感器的飛行高度為560 m,掃描頻率40.3 Hz,激光脈沖波長1 064 nm,光束發(fā)散角0.22 mrad,光斑直徑0.28 m,點云密度約10 個/m2。
設(shè)置點云數(shù)據(jù)的強度范圍為0~255,建立對應(yīng)的灰度圖像,提取與林區(qū)DEM“同源”的正射影像(分辨率1 m),作為ZY-3數(shù)據(jù)定向處理的平面控制基準(zhǔn)。
在林分平均樹高估測過程中,需要解決的幾個關(guān)鍵問題:①利用機載LiDAR數(shù)據(jù)提取DEM;②利用ZY-3數(shù)據(jù)提取DSM;③建立林分平均樹高估測模型并評價模型精度。具體的技術(shù)流程如圖2所示。
本研究獲取的原始LiDAR數(shù)據(jù)通常包含部分噪聲點,因此在使用之前需對其進行去噪處理。利用不規(guī)則三角格網(wǎng)(Triangulated irregular network,TIN)算法將去噪后的LiDAR點云數(shù)據(jù)分為地面點及非地面點,并利用TIN算法將離散地面點生成能夠有利描述地形細節(jié)的不規(guī)則三角網(wǎng),最后內(nèi)插生成數(shù)據(jù)結(jié)構(gòu)簡單,便于使用和管理的規(guī)則格網(wǎng)林區(qū)DEM(分辨率1 m)。
圖 2 技術(shù)流程圖Fig.2 Flow chart of the method
3.2.1 DSM提取方法
對立體影像進行絕對定向,相對定向,平差處理,最終利用多類像對組合提取包含地面林分樹高的林區(qū)DSM。
為完成絕對定向,以林區(qū)DEM為高程基準(zhǔn),“同源”正射影像為平面基準(zhǔn),采用快速傅里葉相位變換匹配法進行“虛擬控制點”采集,根據(jù)ZY-3的RPC參數(shù)文件建立有理函數(shù)模型,使用“虛擬控制點”將立體影像進行移動、轉(zhuǎn)動和放縮,建立像面點坐標(biāo)與地面點坐標(biāo)對應(yīng)關(guān)系。在控制點采集過程中,控制點搜索半徑不小于參考點如圖3a所示相對于三視影像上同名點各自偏移距離的最大值,本研究中參考點相對于前視,下視和后視影像上同名點偏移距離分別如圖3b,圖3c和圖3d所示,最終確定搜索半徑為200 m,同時設(shè)置搜索拒絕閾值為0.75。
圖 3 同名點偏移情況Fig.3 Homonymy point offsets
在絕對定向基礎(chǔ)上,以林區(qū)DEM為高程控制基準(zhǔn),采用快速傅里葉相位變換匹配法進行連接點采集,以光柵格網(wǎng)法進行采樣,還原立體像對中左右影像在成像時的相對位置,獲得所有成對的同名光線相交,完成相對定向。在連接點采集過程中,設(shè)置搜索半徑為100 m,拒絕閾值0.76。同名點采集完成后利用最小二乘原理進行平差,減小像點與實際地面點對應(yīng)的差異。
根據(jù)空中三角測量原理,利用“虛擬控制點”求出地面加密點的物方三維坐標(biāo),將三類立體像對進行組合,基于離散加密點并以“丘陵”地貌類型為約束,提取前下視DSM,下后視DSM和前后視DSM,通過自適應(yīng)方法將三類DSM最優(yōu)化合成三視DSM。
3.2.2 DSM高程精度評價方法
理想情況下,在無地物遮擋的地面點,DEM與DSM的高程應(yīng)保持一致。由于ZY-3數(shù)據(jù)與機載LiDAR數(shù)據(jù)的類型、成像時間、分辨率及處理方法等差異,加之地形起伏變化,實驗結(jié)果不可能達到理想情況。同時,林區(qū)DSM與高程基準(zhǔn)林區(qū)DEM的差值為圖上量測樹高值,該值的準(zhǔn)確性直接影像樹高估測模型精度,故需要對DSM數(shù)據(jù)進行高程精度評價,選擇高程精度最佳的DSM進行平均樹高估測。在無覆蓋的城區(qū)道路面提取同名點,同時,為保證研究區(qū)高程精度分析結(jié)果的準(zhǔn)確性,在林區(qū)寬闊無植被覆蓋的道路上也選取少量的同名點進行高程比較。求同名點DSM與DEM的差值,計算DSM數(shù)據(jù)高程均方根誤差(Root Mean Square Error, RMSE),高程平均正誤差值和高程平均負誤差值。RMSE和高程平均正、負誤差絕對值越小,說明DSM高程精度越好。
ZY-3多光譜影像以正射校正后的下視影像為參考影像進行幾何校正,將校正后的多光譜影像進行波段組合計算,提取10個植被指數(shù)如下所示:修正土壤調(diào)整植被指數(shù),比差植被指數(shù),修正單比植被指數(shù),修正三角植被指數(shù),非線性指數(shù),優(yōu)化土壤調(diào)整植被指數(shù),歸一化差值植被指數(shù),綠色歸一化植被指數(shù),三角植被指數(shù),紅外百分比植被指數(shù)。
3.4.1 基于回歸分析方法建立平均樹高估測模型
以圖上樹高量測值為自變量,實測樹高值為因變量進行一元回歸分析,建立初始樹高估測模型。圖上樹高量測值如式(2)所示:
式中,Hmea為圖上樹高量測值,Hdsm為DSM算數(shù)平均值,Hdem為DEM算數(shù)平均值。
在初始樹高估測模型的基礎(chǔ)上,對各項植被指數(shù)進行T檢驗,認為T值對應(yīng)概率小于0.05的植被指數(shù)引入回歸模型是有意義的,并將該指數(shù)作為回歸分析自變量,最后得到改進樹高估測模型如式(3)所示:
式中,Hpre為樹高改進預(yù)測值,Hmea為樹高圖上量測值,Vk為植被指數(shù)值,k為植被指數(shù)個數(shù)(1≤k≤10),a,bk和c為待求解的改進樹高預(yù)測模型的回歸系數(shù)。
3.4.2 利用DSM高程誤差修正樹高估測模型
ZY-3提取的DSM相對于高程基準(zhǔn)存在高程誤差。在DSM中地面高程高于基準(zhǔn)DEM的區(qū)域ZY-3量測高度值會高于實測值,而在DSM中地面高程低于基準(zhǔn)DEM的區(qū)域ZY-3量測高度值會低于實測值。為了減小DSM高程誤差對預(yù)測結(jié)果的影響,根據(jù)高程精度評價結(jié)果,利用DSM高程誤差對回歸樹高估測模型的預(yù)測值進行修正,在樹高預(yù)測值低于實測樹高的樣地加上DSM高程負誤差絕對值,在樹高預(yù)測值高于實測樹高的樣地減去DSM高程正誤差絕對值。進而得到修正初始樹高估測模型與修正改進樹高估測模型。
3.4.3 模型評價
26個樣方數(shù)據(jù)用于建立樹高估測模型,剩余的11個樣方數(shù)據(jù)用于模型精度評價,其中模型擬合優(yōu)度用AdjR2來評價,預(yù)測精度用P值評價。P值公式如式(4)所示:
式中,H為野外實測樹高,h?為模型預(yù)測樹高。
從機載LiDAR點云中提取了分辨率為1 m的DEM圖像,高程范圍為180.9~385.8 m。經(jīng)目視分析,DEM細節(jié)清晰,在森林覆蓋區(qū)域圖像上直接反應(yīng)林下地面起伏。各驗證樣地的DEM平均值如表1所示。
絕對定向過程中共采集控制點59個,控制點X,Y方向上的均方根誤差分別為1.27 m和1.29 m。相對定向過程中采集二度和三度連接點共275對,連接點X,Y方向上的均方根誤差分別為0.53 m和0.88 m。
最終從ZY-3數(shù)據(jù)提取了4個DSM數(shù)據(jù)(分辨率3.5 m),在全區(qū)共選取了168個同名點對DSM數(shù)據(jù)進行高程精度評價。研究表明,四類DSM數(shù)據(jù)的高程正負誤差分布均表現(xiàn)為無規(guī)律性,其中三視DSM具有最小高程均方根誤差1.433 m,最小高程差方差1.871,最小高程正誤差0.851 m和最小高程負誤差-1.175 m。經(jīng)目視分析,三視DSM中建筑物輪廓清晰,道路線條流暢,DSM過渡光滑,很少有異常劇烈的突起和凹陷,能夠從DSM中再現(xiàn)平面地物起伏的形態(tài)。認為三視DSM為林分平均樹高估測的最佳森林冠頂數(shù)字表面模型。各驗證樣地的DSM平均值如表1所示。
經(jīng)一元回歸建模后得到初始樹高估測模型,如式(5)所示:
經(jīng)多元回歸建模后得到改進樹高估測模型,如式(6)所示。10個植被指數(shù)經(jīng)t檢驗后,僅有優(yōu)化土壤調(diào)整植被指數(shù)(Optimized Soil Adjusted Vegetation Index, OSAVI)與修正土壤調(diào)整植被指 數(shù)(Modified Soil Adjusted Vegetation Index,MSAVI)對應(yīng)的T檢驗值概率顯著小于0.05,適合被引入改進樹高估測模型。這兩個植被指數(shù)均由對植被敏感的紅波段與近紅外波段組合計算得到,能增加植被信號的動態(tài)范圍,同時進一步最小化土壤背景的影響。
表1 實驗結(jié)果?Table 1 Experiment results m
利用11個驗證樣地來驗證模型精度。樹高實測值,圖上量測值,4個模型的預(yù)測值以及精度如表1所示。4個模型對應(yīng)的樹高預(yù)測趨勢及驗證樣地樹高誤差按模型精度由小到大排列,如圖4所示。
初始樹高估測模型精度為85.48%,其中,針葉林精度為86.46%,闊葉林精度為75.70%。改進樹高估測模型精度為87.47%,其中,針葉林精度為88.58%,闊葉林精度為76.41%。修正初始樹高估測模型精度提高到92.10%,其中,針葉林精度提高到93.32%,闊葉林精度提高到79.93%。修正改進樹高估測模型精度提高到93.29%,其中,針葉林精度提高到94.55%,闊葉林精度提高到80.63%。實驗結(jié)果表明,經(jīng)高程誤差修正的改進樹高估測模型為最佳樹高估測模型。
圖4 樹高預(yù)測趨勢及樹高誤差Fig.4 Predicting trends of mean canopy height and errors
最佳模型對針葉林樹高的預(yù)測精度較高,誤差約為1 m,而對闊葉林樹高的預(yù)測精度相對較小,誤差約為5 m??紤]到闊葉林樣地數(shù)量較少,樹高估測模型包含較少的闊葉林高度和光譜信息,所以導(dǎo)致闊葉林樹高的預(yù)測精度較低。
野外樣地實測數(shù)據(jù)是在10月采集的,而用于提取樹高的ZY-3數(shù)據(jù)是在5月收集的,林分在不同時間其生長狀態(tài)也不同,數(shù)據(jù)采集的時間差會導(dǎo)致預(yù)測樹高都要略低于實測樹高,這個現(xiàn)象是合理的。
本研究以吉林省長春市凈月潭國家森林公園林區(qū)為研究區(qū),基于ZY-3影像數(shù)據(jù)、機載LiDAR數(shù)據(jù)和野外實測數(shù)據(jù)進行林分平均樹高估測研究,所得結(jié)論主要如下:
(1)從機載LiDAR點云中提取的“同源”高精度DEM與正射影像數(shù)據(jù),作為ZY-3定向基準(zhǔn)可在無控制點條件下提高其定位測高能力。
(2)在ZY-3立體影像提取的多個DSM中,經(jīng)目視對比分析及高程精度評價分析,確定三視DSM的高程精度優(yōu)于二視DSM。
(3)經(jīng)高程誤差修正后的改進樹高估測模型精度最高,模型的調(diào)整決定系數(shù)AdjR2=0.913,預(yù)測精度為93.29%。
相比王佳,李明澤,胥利紅等學(xué)者利用航片數(shù)據(jù)進行對點的森林立木樹高提取研究,本研究利用ZY-3航空立體影像與機載LiDAR點云數(shù)據(jù),在無控制點條件下進行對面的林分平均樹高提取,既保證了較高的樹高模型估測精度,又實現(xiàn)了較大連續(xù)區(qū)域快速的林分平均樹高估測,亦成功實驗了將測繪數(shù)據(jù)應(yīng)用到林分測高。由于本研究中闊葉林樣本數(shù)量較少,樹高估測模型對該類林分建立的高度預(yù)測能力較針葉林林分弱,因此在下一步的研究中將加大闊葉林樣的比例,甚至增加針闊混交類型的林分,提高估測模型對各類林分平均樹高的估測精度以及模型自身的普適能力。