王佳希, 鄧家勇, 張 巖?, 宋曉鵬, 馬 瑞

(1.北京林業(yè)大學(xué)水土保持學(xué)院，100083，北京；2.山西吉縣森林生態(tài)系統(tǒng)國(guó)家野外科學(xué)觀測(cè)研究站，100083，北京)

地形地貌是影響土壤侵蝕的重要因子。土壤侵蝕研究中，數(shù)字高程模型(digital elevation model, DEM)是重要的地形地貌基礎(chǔ)數(shù)據(jù)，可以用于提取坡度、坡向及水系特征；計(jì)算侵蝕溝形態(tài)指標(biāo)和土壤侵蝕模型中的地形因子；分析土壤侵蝕過(guò)程中的地貌變化[1-4]。DEM的制作早期只能通過(guò)實(shí)地勘測(cè)或者使用已有地形圖[5]。隨著測(cè)繪技術(shù)的不斷發(fā)展[6-7]，基于GPS[8]、衛(wèi)星和航空遙感[9]、攝影測(cè)量[10]、三維激光掃描[11]等都可制作DEM。三維激光掃描所產(chǎn)生的數(shù)據(jù)是具有空間坐標(biāo)及高程信息的點(diǎn)，稱為點(diǎn)云[11]。雖然三維激光掃描具有一定的穿透作用，但由于植被的反射作用等，點(diǎn)云中仍會(huì)存在植被等非地面點(diǎn)。所以點(diǎn)云仍需經(jīng)過(guò)濾波、插值[12-14]而最終生成高精度的DEM[13,15-19]。點(diǎn)云濾波方法有數(shù)學(xué)形態(tài)法濾波[20]、坡度濾波[21]、基于回波次數(shù)[13]、基于回波強(qiáng)度[22]等。常用的插值方法有反距離權(quán)重法(inverse distance weighted，IDW)、克里金法(Kriging)、自然鄰域法(natural neighbor，NN)[14]和不規(guī)則三角網(wǎng)(triangulated irregular network，TIN)[8]等。

點(diǎn)云濾波在平坦、植被較單一地區(qū)相對(duì)簡(jiǎn)單，在崎嶇、植被低矮復(fù)雜地區(qū)較為困難[22]。黃土高原丘陵溝壑區(qū)土壤侵蝕嚴(yán)重，地形破碎，植被多為次生演替的灌木和草本，屬于點(diǎn)云濾波困難的地區(qū)。潘少奇等[23]將點(diǎn)云用TIN法插值生成DEM，與1∶1萬(wàn)的地形圖比較，精度有所提升，但沒(méi)有植被濾波處理。馬鼎等[22]利用點(diǎn)云的激光回波強(qiáng)度衰減補(bǔ)償模型進(jìn)行了植被濾波，得到的DEM效果較好，但需要點(diǎn)云記錄有回波強(qiáng)度信息。

如上，植被對(duì)于DEM的影響研究多集中于植被濾波、插值方法的開(kāi)發(fā)[24-25]，而關(guān)于植被季相變化對(duì)DEM的影響研究較少。植被非生長(zhǎng)季時(shí)，枯萎凋落，激光腳點(diǎn)落在地面的可能性增大，DEM精度可能會(huì)提高；但是植被季相變化具體有何影響，缺乏數(shù)據(jù)支撐。筆者基于黃土高原丘陵溝壑區(qū)植物生長(zhǎng)季 (2018年9月)、非生長(zhǎng)季(2018年12月)2次實(shí)地測(cè)繪的點(diǎn)云，分析生長(zhǎng)季與非生長(zhǎng)季的點(diǎn)云以及DEM差異，討論植被季相變化是否對(duì)其具有影響，探究提高DEM精度的手段。

1 研究區(qū)概況

陜西省綏德縣橋溝小流域(E 110°17′24″～110°17′43″，N 37°29′33″～37°30′13″)位于無(wú)定河下游(圖1)，屬于黃土高原丘陵溝壑區(qū)第一副區(qū)。氣候類型為溫帶半干旱大陸性季風(fēng)氣候。土壤類型均為黃綿土亞類。植被類型以草地、灌木林地為主，有少量人工林地。植物群落多為次生演替，代表性群系有白羊草(Bothriochloaischaemum)群系、白草(Pennisetumflaccidum)-鐵桿蒿(Artemisiasacrorum)群系、茭蒿(Artemisiagiraldii)群系、茵陳蒿(Artemisiacapillaris)群系和興安胡枝子(Lespedezadavurica)群系。人工林造林主要有刺槐(Robiniapsdudoacacia)、油松(Pinustabulaeformis)、側(cè)柏(Platycladusorientalis)、檸條(Caraganakorshinskii)等[26]。在橋溝小流域中選取A、B和C樣區(qū)。樣區(qū)A(1 445 m2)切溝溝頭邊緣平滑，溝壁較陡，溝道比降在1.08～1.25之間；以禾本科，菊科等草本為主，溝底有檸條，冠層高度約20～60 cm；植被覆蓋度較低。樣區(qū)B(3 177 m2)溝壁較陡，溝頭邊緣鋒利，溝道比降在0.75～1.11之間，切溝較大；除垂直溝壁裸露外，在溝底、溝間有零星刺槐分布，其高度約2.5～3.5 m；灌草覆蓋度較高，灌木冠層高度約30～50 cm。樣區(qū)C(1 386 m2)切溝密集，溝道比降較大(1.2～2.0)，但切溝寬度較小；以低矮灌木和草本為主，灌木主要分布在溝底，冠層高度20～30 cm；植被覆蓋度中等。

A、B和C為三維激光測(cè)量樣區(qū)編號(hào)。底圖來(lái)自Google 影像。Sampling sites are named A，B and C for 3D laser measurement. Background images are from Google Earth images. 圖1 研究區(qū)位置圖Fig.1 Location of the study area

2 材料與方法

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

已有研究表明黃土高原植被覆蓋在7—9月達(dá)到峰值，12月達(dá)到最低值[27]。本研究分別于植物生長(zhǎng)季 (2018年9月)、非生長(zhǎng)季(2018年12月)，使用拓普康IS—IMAGING STATION型三維激光掃描儀對(duì)橋溝內(nèi)切溝特征差異較大的3個(gè)樣區(qū)(A、B和C)內(nèi)11條切溝進(jìn)行了掃描。為了保證測(cè)量精度，無(wú)障礙物遮擋，儀器架設(shè)點(diǎn)選擇在樣區(qū)對(duì)面山坡視野最佳且大致水平的位置，并用手持GPS確定架設(shè)點(diǎn)地理坐標(biāo)及高程。2次測(cè)量使用了相同的架設(shè)點(diǎn)，保證2次測(cè)量之間相對(duì)誤差較小(圖1)。點(diǎn)云平面間隔為0.15 m，密度約30～50 pts/m2。所得點(diǎn)云不包括回波次數(shù)、強(qiáng)度等。另外，在2次外業(yè)時(shí)，使用激光測(cè)距儀和卷尺抽樣測(cè)量部分優(yōu)勢(shì)種群的株高，用來(lái)輔助判斷植被季相變化。

2.2 數(shù)據(jù)處理與分析

2.2.1 點(diǎn)云處理流程 第1步，通過(guò)人工判讀除去包括飛鳥(niǎo)、電力線路、儀器偶然誤差等造成的離群值[13]。第2步，基于Terrascan，利用坡度濾波法[21]將點(diǎn)云分類為地面點(diǎn)和非地面點(diǎn)(大部分為植被點(diǎn))，隨后去除分類中的非地面點(diǎn)，僅保留地面點(diǎn)。第3步，將濾波前后點(diǎn)云用TIN法[8]、克里金法[28]分別插值構(gòu)建DEM(柵格分辨率0.15 m)。按照上述1～3步分別處理9月、12月數(shù)據(jù)。最終分別得到2季相的4幅DEM：TIN法、濾波-TIN法、克里金法、濾波-克里法。

2.2.2 點(diǎn)云和DEM精度評(píng)價(jià)方法 對(duì)點(diǎn)云高程進(jìn)行數(shù)理統(tǒng)計(jì)，計(jì)算點(diǎn)云平均高程、標(biāo)準(zhǔn)差、最大值和最小值等，對(duì)比2次點(diǎn)云整體差異。

采用中誤差和判定系數(shù)表征DEM和點(diǎn)云的符合程度。中誤差，即均方根誤差(root mean square error，RMSE)。它是最常用的體現(xiàn)DEM偏移情況的指標(biāo)，計(jì)算公式為：

(1)

式中：σ為中誤差，m；R為點(diǎn)云高程，m；Z為插值高程，m；n為測(cè)得的總點(diǎn)數(shù)；m為第m個(gè)點(diǎn)。中誤差越小，DEM越符合點(diǎn)云。

判定系數(shù)(R2)用以表征插值高程與實(shí)測(cè)高程總體吻合程度。計(jì)算公式為：

(2)

圖2 A、B和C樣區(qū)點(diǎn)云高程頻數(shù)分布圖Fig.2 Frequency distribution of point cloud elevation in the site A，B and C

本研究共選取7條切溝縱剖面，29條橫剖面來(lái)提取溝長(zhǎng)、溝寬和溝深，用來(lái)代表不同方法生成的DEM地貌形態(tài)差異。另外，通過(guò)三維可視化、坡度分析和切溝縱剖面圖對(duì)比不同季相，不同濾波插值方法對(duì)生成DEM的影響。三維可視化可以直觀判斷DEM質(zhì)量，地形地貌還原程度，植被覆蓋特征等。坡度圖可以直觀地顯示出DEM坡度變化，植被在其中多表現(xiàn)為“鼓包”，另外從坡度圖中也很容易識(shí)別出溝沿線等地貌特征。本研究區(qū)地面粗糙度不大，而三維可視化圖、坡度圖中粗糙度較大區(qū)域，也一定程度上反映了植被表面的影響。

3 結(jié)果與分析

3.1 植被季相變化

3.1.1 點(diǎn)云高程統(tǒng)計(jì)特征 點(diǎn)云高程統(tǒng)計(jì)分析表明，樣區(qū)A、B和C生長(zhǎng)季平均高程比非生長(zhǎng)季分別高0.22、0.47和0.07 m。實(shí)地抽樣方面，共測(cè)量?jī)?yōu)勢(shì)種灌木、多年生草本3種，9和12月分別測(cè)184和189株。其平均株高分別為0.73和0.68 m。測(cè)量的9株喬木高度無(wú)明顯變化。非生長(zhǎng)季點(diǎn)云標(biāo)準(zhǔn)差普遍比生長(zhǎng)季略低，且中位數(shù)、眾數(shù)也有不同程度減小，說(shuō)明非生長(zhǎng)季高程分布更趨于集中，且整體偏低。從圖2可看出3個(gè)樣區(qū)類似的偏態(tài)分布規(guī)律，且非生長(zhǎng)季比生長(zhǎng)季的高值部分略少，其他部分明顯增多。樣區(qū)B點(diǎn)云頻率分布特征差異最大。植被覆蓋度方面，樣區(qū)A約0.60，樣區(qū)B約0.85，樣區(qū)C約0.75，說(shuō)明點(diǎn)云高程變化和植被覆蓋度有一定關(guān)系，植被覆蓋度越高，點(diǎn)云高程季相變化越大。

3.1.2 植被季相變化與濾波效果 3個(gè)樣區(qū)12月所測(cè)點(diǎn)云高程平均標(biāo)準(zhǔn)差比9月有所降低，說(shuō)明非生長(zhǎng)季的點(diǎn)云高程離散趨勢(shì)減弱。12月所測(cè)數(shù)據(jù)濾波后保留的點(diǎn)云更多(表1)。從點(diǎn)云剖面圖可看出，9月所測(cè)點(diǎn)云中高大植被在12月幾乎消失，這說(shuō)明經(jīng)過(guò)坡度濾波后，高大喬灌被去除同時(shí)地面點(diǎn)盡可能地被保留(圖3)。樣區(qū)A、B和C點(diǎn)云平均高程分別降低了0.25、0.56和0.44 m(9月)。

表1 4種濾波插值方法所得DEM精度統(tǒng)計(jì)

深色點(diǎn)云為濾波前，淺色點(diǎn)云為濾波。 The dark point cloud is the one before filtering, and the light point cloud is the one after filtering圖3 濾波前后切溝橫剖面點(diǎn)云分布圖Fig.3 Distribution of point cloud of gully profiles before and after filtering

3.2 不同方法生成DEM對(duì)比

首先，9月所測(cè)點(diǎn)云經(jīng)過(guò)濾波之后，點(diǎn)云數(shù)目大幅度下降。12月所測(cè)點(diǎn)云經(jīng)過(guò)濾波之后，點(diǎn)云數(shù)目也有所下降，但幅度不及9月。尤其樣區(qū)A，可能是由于其以草本為主，且植被蓋度較低?？死锝鹋cTIN法所得點(diǎn)云高程均值幾乎相同，但濾波之后高程有所降低：9月所測(cè)點(diǎn)云大概0.4 m，12月所測(cè)點(diǎn)云在0.2 m左右。TIN法與克里金法交叉驗(yàn)證的R2都在0.9以上，且克里金和TIN法中誤差在0.30～0.50 m之間，濾波后插值中誤差降低至0.20～0.10 m(表1)。TIN法中誤差普遍比克里金法要小，說(shuō)明TIN法所得高程更接近于實(shí)測(cè)點(diǎn)云，地形地貌刻畫(huà)精度較高。

由圖4可見(jiàn)， b比a插值結(jié)果更平滑，粗糙度降低，但是c，d粗糙度差異不大。4種方法中溝沿線都十分清晰，a保留的地貌和植被細(xì)節(jié)更多。a、b中的高大植被在c、d中已被去除。a3中2次測(cè)得剖面相比，中上部基本吻合，下部12月比9月明顯更低，并且12月下部與中上部更為平滑，更符合實(shí)地特征。這是由于溝底的喬灌在12月落葉枯萎造成的。b3比a3相比，其不能去除9月植被且平滑了溝頭坡度。經(jīng)過(guò)濾波之后，c3，d3分別比a3、b3更為平滑，并且2次測(cè)得的結(jié)果基本一致，但是d3溝頭部分仍有明顯差異。

a為TIN法，b為克里金法，c為濾波-TIN法，d為濾波-克里金法；1為三維可視化疊加山體陰影效果圖，2為坡度圖，3為對(duì)應(yīng)三維可視化疊加山體陰影效果圖中的白色剖線。a refers to TIN, b refers to Kriging, c refers to filtering-TIN, and d refers to filtering-Kriging. 1 refers to 3D model with hill shade effect, 2 refers to slope map, and 3 refers to profiles of white lines in corresponding 3D model.圖4 4種方法生成的DEM比較(以樣區(qū)B為例)Fig.4 Comparison of DEMs by four processes (Taking site B as an sample)

在2次三維激光掃描之間，研究區(qū)內(nèi)沒(méi)有發(fā)生強(qiáng)烈的侵蝕性降雨，不存在劇烈侵蝕的情況。DEM中切溝形態(tài)參數(shù)變化主要是由植被變化引起的。切溝形態(tài)參數(shù)統(tǒng)計(jì)表明(表2)，克里金法所得的12月平均溝長(zhǎng)、溝寬和溝深都更大，相差0.05～0.35 m。而TIN法所得12月平均溝長(zhǎng)和溝寬比9月分別減小0.72和0.26 m，但溝深更大。濾波-TIN法所得12月平均溝長(zhǎng)、溝寬和溝深比9月分別增大0.03、0.25和0.10 m，說(shuō)明具有去除植被的效果。由于經(jīng)過(guò)濾波后，12月保留更多地貌特征點(diǎn)，即非生長(zhǎng)季比生長(zhǎng)季更接近于真實(shí)地貌。濾波-克里金法所得12月份的溝長(zhǎng)和溝深分別比9月份小0.36和0.34 m，說(shuō)明克里金法的平滑作用使得插值結(jié)果稍偏離實(shí)際地貌形態(tài)。

表2 4種處理方法提取的切溝形態(tài)參數(shù)比較

4 討論

4.1 植被季相變化對(duì)濾波效果的影響

點(diǎn)云平均高程和標(biāo)準(zhǔn)差等都表明12月點(diǎn)云高程總體上低于9月，高程分布更為集中。點(diǎn)云高程高值部分在9月數(shù)量比12月多，而中低值9月比12月數(shù)量要少(圖2)。這是由于非生長(zhǎng)季植被枯萎凋零，覆蓋度降低，植被平均高度降低，使激光腳點(diǎn)落在地面的概率增大。由于三維激光掃描儀無(wú)法測(cè)量到植株的最高點(diǎn)，所以點(diǎn)云中的植被點(diǎn)總小于最大植株高度[29]，即不能精確反映植被高度，加之抽樣株高不能完全反映植被總體情況，所以9月與12月點(diǎn)云平均高程差(0.10 m)，略高于抽樣平均株高差(0.05 m)?？傮w上看，植被變化與點(diǎn)云變化規(guī)律相符合，表明植被變化對(duì)于點(diǎn)云具有一定影響。

對(duì)濾波前后的點(diǎn)云和DEM進(jìn)行了對(duì)比分析，結(jié)果表明坡度濾波法對(duì)于去除植被具有一定效果。在植被濾波完全去除植被的理想狀態(tài)下，不同季相所得濾波結(jié)果應(yīng)完全一致。但本文9月與12月植被濾波后的點(diǎn)云并不一致，樣區(qū)A 2次點(diǎn)云平均高程相差在0.11 m左右，樣區(qū)B和C相差了0.3 m左右。所以即使濾波之后，植被季相變化對(duì)點(diǎn)云高程的影響仍然存在，進(jìn)而影響DEM精度。

4.2 插值方法對(duì)DEM精度的影響

前人比較克里金法與最鄰近法、反距離權(quán)重法等將點(diǎn)云插值生成DEM的精度，結(jié)果表明克里金法效果較好[28]。本研究將克里金法和TIN法進(jìn)行比較，發(fā)現(xiàn)TIN法比克里金法R2更高，σ更小，且σ均符合國(guó)家測(cè)繪局1∶500 DEM一級(jí)精度要求(坡度25°以上，柵格分辨率為0.5 m的，一級(jí)精度σ<0.70，坡度6°～25°，一級(jí)精度σ<0.50)。TIN法生成的DEM與點(diǎn)云擬合效果更好，精度更高，對(duì)濾波前后的地形變化更為敏感。而克里金法有一定的平滑效果，一方面有使平滑高大植被的效果，使其高程一定程度地降低(9月，σ=0.46)，但另一方面，也會(huì)出現(xiàn)“牛眼”現(xiàn)象[28]，造成與實(shí)際地貌的偏離。另外，未經(jīng)濾波處理直接生成DEM時(shí)，克里金法比TIN法得到的平均溝長(zhǎng)溝寬更大，但平均溝深更淺，而經(jīng)過(guò)濾波后，克里金法比TIN法生成的溝長(zhǎng)溝寬和溝深更小。可能是由于克里金法插值時(shí)，溝沿線附近、溝底的部分點(diǎn)云插值權(quán)重過(guò)低，導(dǎo)致溝沿線附近、溝底略緩，出現(xiàn)DEM平滑現(xiàn)象。

總之，TIN法生成的DEM溝長(zhǎng)溝寬等參數(shù)更接近于真實(shí)的切溝地貌形態(tài)。經(jīng)過(guò)濾波之后比未濾波生成DEM的切溝形態(tài)參數(shù)有所增大，說(shuō)明具有一定地植被去除效果。

5 結(jié)論

1)植被季相變化對(duì)于三維激光掃描具有一定影響，生長(zhǎng)季測(cè)量的點(diǎn)云高程均值高于非生長(zhǎng)季0.25 m，從而進(jìn)一步影響到DEM的精度和地形地貌刻畫(huà)。

2)基于濾波-TIN法生成的DEM平均溝長(zhǎng)、溝寬和溝深，2次相比，分別增大0.03、0.25和0.10 m，結(jié)合2次實(shí)地抽樣的植物株高，說(shuō)明坡度濾波能一定程度上去除植被。但濾波后9月份(生長(zhǎng)季)點(diǎn)云高程均值仍然比12月份(非生長(zhǎng)季)高出0.15 m，兩者并不完全一致，表明植被季相變化對(duì)DEM所造成的影響無(wú)法完全去除。

3)與克里金法相比，TIN法生成的DEM精度較高。前者存在平滑現(xiàn)象，植被生長(zhǎng)季時(shí)能粗略地去除高大植被，而在非生長(zhǎng)季時(shí)平滑溝道邊緣。TIN法比克里金法的插值結(jié)果更接近實(shí)際地貌。

黃土高原丘陵溝壑區(qū)植被季相變化對(duì)DEM精度具有影響，且僅通過(guò)植被濾波算法很難完全去除。選擇植被非生長(zhǎng)季時(shí)進(jìn)行測(cè)繪，并使用適合的濾波插值方法，才能有效地提高DEM的精度。