繆志修 齊 華,2 王國(guó)昌 程 昂

(1.西南交通大學(xué)信息科學(xué)與技術(shù)學(xué)院,四川成都 610031；2.西南交通大學(xué)土木工程學(xué)院,四川成都 610031；3.中鐵二院工程集團(tuán)有限責(zé)任公司測(cè)繪分院,四川成都 610031)

機(jī)載激光雷達(dá)(LiDAR)是取高精度和高密度地形數(shù)據(jù)的一種新的測(cè)量手段。LiDAR數(shù)據(jù)擁有的精度和密度使其在生成DEM數(shù)據(jù)方面的應(yīng)用得到進(jìn)一步的發(fā)展。然而,高密度的數(shù)據(jù)導(dǎo)致數(shù)據(jù)量的大幅增加,給數(shù)據(jù)的組織、處理、后期應(yīng)用都帶來(lái)了諸多不便。

在大多數(shù)情況下,對(duì)一個(gè)具體DEM項(xiàng)目來(lái)說(shuō),并不需要得到DEM表面所表達(dá)出來(lái)的全部信息,只需量測(cè)表達(dá)相應(yīng)地表所需要的數(shù)據(jù)點(diǎn),以達(dá)到一定的地形表面精度和可信度即可[1]。特別是在地勢(shì)較平坦地區(qū),不需很高密度的激光腳點(diǎn)。一般說(shuō)來(lái),采樣的地形數(shù)據(jù)越密集,生成的DEM精度越高。但是,在LiDAR數(shù)據(jù)采集過(guò)程中,沒有將地形類型和采樣密度相匹配,過(guò)度采樣是不可避免的,結(jié)果增加了數(shù)據(jù)存儲(chǔ)量和處理時(shí)間[9]。因此,如何有根據(jù)地進(jìn)行地形數(shù)據(jù)抽稀(如信息含量與數(shù)據(jù)集容量的比例),得到一個(gè)更易于處理和操作,用于DEM生成的地形數(shù)據(jù)集,具有重要的實(shí)用價(jià)值。

數(shù)據(jù)抽稀的主要目的是在采樣密度和數(shù)據(jù)精度之間達(dá)到一種最優(yōu)平衡,以此來(lái)得到一個(gè)有效的高精度DEM?？紤]到不同數(shù)據(jù)元素對(duì)生成的DEM精度有不同的貢獻(xiàn),數(shù)據(jù)抽稀應(yīng)該遵循保留關(guān)鍵要素,去掉次要元素的原則[8]。

1 國(guó)內(nèi)外研究情況

目前國(guó)內(nèi)外關(guān)于LiDAR數(shù)據(jù)構(gòu)建DEM抽稀的算法可以歸結(jié)為以下幾種方法：基于非選擇性的數(shù)據(jù)抽稀、基于TIN的數(shù)據(jù)抽稀、基于距離及角度的數(shù)據(jù)抽稀等。

1.1 基于非選擇性的抽稀算法

基于非選擇性的數(shù)據(jù)抽稀,主要包括兩種：第一種是基于格網(wǎng)的數(shù)據(jù)抽稀。基于格網(wǎng)的數(shù)據(jù)抽稀是將輸入的所有數(shù)據(jù)建立虛擬格網(wǎng),然后只保留虛擬格網(wǎng)內(nèi)一個(gè)數(shù)據(jù)點(diǎn)。其中虛擬格網(wǎng)的大小要依據(jù)數(shù)據(jù)點(diǎn)的密度來(lái)設(shè)定,如果數(shù)據(jù)點(diǎn)的密度達(dá)到每平方米幾個(gè)點(diǎn),則可以設(shè)定格網(wǎng)大小為1 m。這種算法雖然比較簡(jiǎn)單,但無(wú)法顧及到地形特征,這將一定程度上降低數(shù)據(jù)的精度。該方法隨壓縮程度的增加,特征點(diǎn)被刪除的可能性也增加,可能丟失關(guān)鍵數(shù)據(jù)信息。目前,這種方法一般只應(yīng)用到對(duì)精度要求較低的工程項(xiàng)目上。

另外一種是基于系統(tǒng)的抽稀。在系統(tǒng)抽稀中,對(duì)于一個(gè)大樣本N,首先確定抽樣間隔,比如每隔10個(gè)點(diǎn)作為抽樣間隔,然后在這10個(gè)樣本點(diǎn)中隨機(jī)選擇一個(gè)樣本數(shù)據(jù)點(diǎn)。該方法和基于虛擬格網(wǎng)的方法一樣,不能很好地保留地形特征。這種方法主要用來(lái)做快速顯示,比如在TerraScan軟件中就采用該種方法來(lái)快速查看點(diǎn)云數(shù)據(jù)。

1.2 基于Tin的抽稀算法

基于非選擇性抽稀不能很好的保留地形特征,在一定程度上降低了DEM精度。為了能夠更好地保留地形特征,有學(xué)者提出了基于Tin的抽稀方法。劉春等提出的基于Tin的抽稀方法的算法原理是利用已經(jīng)建立好的不規(guī)則三角網(wǎng)模型,然后依據(jù)平坦區(qū)域上的點(diǎn)構(gòu)成三角網(wǎng)包含的信息量最少,同時(shí)這些三角形構(gòu)成面的法線向量也就趨近于平行,也就是說(shuō)法線向量的夾角越小。利用三角網(wǎng)的這個(gè)性質(zhì)可以實(shí)現(xiàn)對(duì)DEM數(shù)據(jù)的抽稀[3]。

該算法能夠較好的保留地形特征點(diǎn)。但是,該算法在判斷一點(diǎn)是否去除時(shí),需查找包含該點(diǎn)的三角面,判斷不同三角面之間的夾角,這在一定程度上降低了數(shù)據(jù)處理的速度。

1.3 基于距離及高差的抽稀算法

美國(guó)PAMELA S.LAFONTAINE女士提出的DDR(Data Density Reduction)算法主要依靠?jī)蓚€(gè)參數(shù)：搜索半徑及高差閾值。搜索半徑也就是搜索的水平距離,用來(lái)判定滿足指定半徑鄰域內(nèi)點(diǎn)云數(shù)據(jù)。這個(gè)參數(shù)主要控制抽稀的程度；高差閾值也就是點(diǎn)云數(shù)據(jù)之間的垂直距離,主要用來(lái)控制是否接收點(diǎn),也就是控制抽稀的精度。如果在給定鄰域內(nèi)的一個(gè)點(diǎn)的高程值和中心指定點(diǎn)的高程小于給定的高差閾值,則可以認(rèn)為這些點(diǎn)中存在冗余數(shù)據(jù)[10]。

該算法能夠較好地保留地形的特征,其保留地形特征點(diǎn)的能力介于基于非選擇性的抽稀方法和基于TIN的抽稀方法之間。由于該算法不斷的搜索一定半徑范圍內(nèi)的數(shù)據(jù)點(diǎn),如果數(shù)據(jù)量過(guò)大,會(huì)在很大程度上降低數(shù)據(jù)處理速度。因此需采用一種策略,即將原始數(shù)據(jù)分成一系列的小塊,則在搜索一定距離范圍內(nèi)的數(shù)據(jù)點(diǎn)時(shí),只需搜尋該塊以及與該塊相鄰的其他數(shù)據(jù)塊,并不用搜索所有的數(shù)據(jù)塊。這樣將會(huì)提高數(shù)據(jù)的處理速度。

2 基于地形坡度的抽稀算法

以上的幾種算法中,基于非選擇性的抽稀算法不能夠很好的保留地形特征信息,即地形變化幅度大的地形與地形變化幅度較小的地形同等地受到削弱。而基于TIN的抽稀算法和基于角度及距離限制的抽稀算法,對(duì)于不同的地形在運(yùn)行時(shí)采用相同的參數(shù),一般較難得到理想的結(jié)果。

為了考慮能夠?qū)Σ煌匦蔚腄EM進(jìn)行數(shù)據(jù)的抽稀,本文提出了基于地形坡度的抽稀算法。該算法思想如下：首先利用所有的數(shù)據(jù)構(gòu)建TIN,為了判斷一個(gè)數(shù)據(jù)點(diǎn)是否保留(如圖1中的P點(diǎn)),則需判斷以該點(diǎn)相鄰所有三角面(如圖1中與P點(diǎn)相鄰的三角形有5個(gè))的坡度。如果這些三角形中的最大坡度值與最小坡度值之差大于給定的閾值,則認(rèn)為這個(gè)點(diǎn)對(duì)地形的影響較大,需保留。

地形是比較復(fù)雜的,可以根據(jù)坡度與高差的不同將地形分成平原、丘陵、山地以及高山地[7]。如表1所示,可將地形的坡度分成四個(gè)段：第Ⅰ坡度段為0～2°；第Ⅱ坡度段為2°～6°；第Ⅲ坡度段為6°～25°；第Ⅳ坡度段為大于25°。因此,本文提出根據(jù)不同地形的坡度段來(lái)設(shè)置抽稀閾值,即每個(gè)坡度段設(shè)置相應(yīng)的抽稀閾值。LiDAR數(shù)據(jù)構(gòu)建的DEM數(shù)據(jù)點(diǎn)的平均間距約為1 m,則與P點(diǎn)相鄰的三角形所構(gòu)成的面積一般在4 m2左右。在這個(gè)小區(qū)域面積的范圍內(nèi),可以通過(guò)求得該點(diǎn)相鄰三角面的平均坡度,來(lái)判斷該局部區(qū)域坡度屬于哪個(gè)坡度段。比如該點(diǎn)相鄰的所有三角面的平均坡度小于2°,則可以認(rèn)為該局部區(qū)域是屬于較平坦的地形,則最大坡度差將和第Ⅰ坡度段設(shè)置的抽稀閾值進(jìn)行比較。同理,如果與該點(diǎn)相鄰的所有三角面的平均坡度在2°～6°之間,則最大坡度差將和第Ⅱ坡度段設(shè)置的抽稀閾值進(jìn)行比較。

表1 根據(jù)坡度和高差的地形分類

圖1 與P點(diǎn)相鄰的所有三角形

算法步驟如下：

(1)讀入所有的數(shù)據(jù)點(diǎn),生成不規(guī)則三角網(wǎng)。

(2)在生成的三角網(wǎng)數(shù)據(jù)點(diǎn)中,讀入第一個(gè)數(shù)據(jù)點(diǎn),并查找與該點(diǎn)相鄰的三角形。

(3)在找到的所有三角形中,計(jì)算每個(gè)三角面的坡度。

(4)求得這些三角面的平均坡度。

(5)計(jì)算這些三角面的最大坡度和最小坡度,并求得最大坡度和最小坡度之差,即最大坡度差。

(6)如果最大坡度差為零(如圖2所示),則需判斷該點(diǎn)到所在投影三角面的距離(如圖3中E點(diǎn)到△ABC的距離),如果這個(gè)距離大于給定的閾值,則保留這個(gè)點(diǎn)。

圖2 與E點(diǎn)相鄰三角面最大坡度差為零的局部地形

(7)如果最大坡度差不為零,判斷平均坡度的大小,以此來(lái)判斷該局部區(qū)域坡度屬于哪個(gè)坡度段。如果屬于第Ⅰ坡度段,則最大坡度差將和第Ⅰ坡度段設(shè)置的閾值進(jìn)行比較；同理,如果屬于第Ⅱ坡度段,則最大坡度差將和第Ⅱ坡度段設(shè)置的閾值進(jìn)行比較。

(8)繼續(xù)讀入三角網(wǎng)數(shù)據(jù)中的下一個(gè)點(diǎn),直到所有的數(shù)據(jù)點(diǎn)都進(jìn)行了判斷。

算法流程如圖3所示。

圖3 基于坡度抽稀算法流程

3 多種抽稀算法實(shí)驗(yàn)分析

3.1 實(shí)驗(yàn)區(qū)域數(shù)據(jù)介紹

為了能更好的說(shuō)明不同的算法對(duì)抽稀后DEM精度的影響,本節(jié)的實(shí)驗(yàn)都采用以下數(shù)據(jù),如圖4所示。該試驗(yàn)區(qū)數(shù)據(jù)位于中國(guó)西部某省,數(shù)據(jù)點(diǎn)個(gè)數(shù)為6 923個(gè),面積大小為29 640 m2,激光點(diǎn)間距為1.2 m,海拔高度在608.71～624.09 m之間。該實(shí)驗(yàn)區(qū)域數(shù)據(jù)的地形暈渲圖如圖5所示。

圖4 點(diǎn)云數(shù)據(jù)

圖5 該區(qū)域地形暈渲

3.2 精度評(píng)定方法

本次實(shí)驗(yàn)的精度評(píng)價(jià)方法如下：

(1)采用ArcGis9.2地統(tǒng)計(jì)分析模塊,利用該模塊從原始數(shù)據(jù)7 287個(gè)點(diǎn)中抽取5%的數(shù)據(jù)點(diǎn),即364個(gè)數(shù)據(jù)點(diǎn)作為樣本點(diǎn)。剩下的6 923個(gè)點(diǎn)作為測(cè)試點(diǎn)。

(2)對(duì)測(cè)試點(diǎn)進(jìn)行抽稀,利用抽稀后的數(shù)據(jù)點(diǎn)內(nèi)插出這364個(gè)樣本點(diǎn)的高程。不同的內(nèi)插方法對(duì)DEM的精度影響是不同的,本實(shí)驗(yàn)的內(nèi)插方法為基于TIN的線性內(nèi)插方法。

(3)通過(guò)計(jì)算,內(nèi)插出樣本點(diǎn)高程和實(shí)際高程的標(biāo)準(zhǔn)差和中誤差來(lái)實(shí)現(xiàn)精度的評(píng)定。

3.3 基于非選擇性抽稀方法的實(shí)驗(yàn)分析

基于非選擇性的抽稀方法中,主要有基于系統(tǒng)抽稀和基于格網(wǎng)抽稀兩種。

基于系統(tǒng)抽稀方法,本次實(shí)驗(yàn)分別將原始數(shù)據(jù)抽稀到9.995%,24.599%,33.323%和49.992%的程度。然后利用這些抽稀后的數(shù)據(jù)內(nèi)插出這364個(gè)樣本點(diǎn)的高程,通過(guò)計(jì)算實(shí)際樣本點(diǎn)高程和內(nèi)插樣本點(diǎn)高程的中誤差和標(biāo)準(zhǔn)差進(jìn)行精度評(píng)定。表2顯示了抽稀程度與DEM精度關(guān)系。

表2 DEM精度與抽稀率

從表2可以看出,當(dāng)抽稀程度大于49.992%時(shí),其中誤差小于30 cm。圖6顯示DEM精度與抽稀程度之間的關(guān)系,從圖6可以更好的看出采用該種方法時(shí),當(dāng)抽稀程度大于50%時(shí),其DEM精度變化幅度是較小的。

圖6 DEM精度與抽稀率

基于格網(wǎng)的DEM數(shù)據(jù)抽稀,本次實(shí)驗(yàn)將原始數(shù)據(jù)抽稀到9.54%,24.931%,41.167%,54.817%,74.707%的程度。表3顯示了抽稀程度與DEM精度關(guān)系。

表3 DEM精度與抽稀率

從表3可以看出,當(dāng)抽稀率大于41.167%時(shí),其DEM精度小于25 cm。而從圖7中可以看出,當(dāng)抽稀率大于40%時(shí),其精度變化幅度也是較小的。

圖7 DEM精度與抽稀率

3.4 基于Tin抽稀方法的實(shí)驗(yàn)分析

基于Tin抽稀方法,本次實(shí)驗(yàn)分別將原始數(shù)據(jù)抽稀到10.169%,24.353%,40.647%,54.759%和74.534%的程度。表4顯示基于Tin方法得到的數(shù)據(jù)抽稀率與DEM的精度關(guān)系。

表4 DEM精度與抽稀率

從表4可以看出,當(dāng)抽稀率大于40.647%時(shí),其中誤差小于20 cm。而從圖8中可以看出,當(dāng)抽稀率大于28%時(shí),其DEM精度的變化幅度較小。

圖8 DEM精度與抽稀率

3.5 基于距離及高差抽稀方法的實(shí)驗(yàn)分析

基于距離及高差抽稀方法,本次實(shí)驗(yàn)分別將原始數(shù)據(jù)抽稀到9.518%,23.154%,43.232%,55.120%和75.501%的程度。表5顯示基于距離及高差抽稀方法得到的數(shù)據(jù)抽稀率與DEM的精度關(guān)系。

表5 DEM精度與抽稀率

從表5可以看出,采用該算法時(shí),抽稀率大于43.232%,其DEM精度小于28 cm。而從圖9中可以看出,當(dāng)抽稀率大于50%時(shí),其DEM精度的變化幅度較小。

3.6 基于地形坡度抽稀方法的實(shí)驗(yàn)分析

基于地形坡度抽稀方法,本次實(shí)驗(yàn)分別將原始數(shù)據(jù)抽稀到13.606%,24.151%,39.000%,53.242%和71.053%的程度。表6顯示基于地形坡度抽稀方法得到的數(shù)據(jù)抽稀率與DEM的精度關(guān)系。

圖9 DEM精度與抽稀率

表6 DEM精度與抽稀率

從表6可以看出,采用該算法時(shí),抽稀率大于24.151%,其DEM精度小于23 cm。而從圖10中可以看出,當(dāng)抽稀率大于24%時(shí),其DEM精度的變化幅度較小。

圖10 DEM精度與抽稀率

3.7 各種不同抽稀方法綜合比較

抽稀結(jié)果表明,無(wú)論數(shù)據(jù)抽稀到何種程度,丟失一些重要的信息是不可避免的,抽稀的目標(biāo)是最大的降低數(shù)據(jù)的密度而損失最小的精度。采用一種合適的抽稀方法可以使本實(shí)驗(yàn)區(qū)域的數(shù)據(jù)在抽稀率達(dá)到25%時(shí),其精度仍能好于23 cm。

為了能夠更好的說(shuō)明各種抽稀算法對(duì)DEM精度的影響,將以上的五種抽稀算法的抽稀率與DEM精度整合到一個(gè)圖中,如圖11所示。

從圖11可以看出,當(dāng)抽稀程度大于24%時(shí),基于Tin的抽稀方法和基于地形坡度的抽稀算法得到的DEM精度的幅度變化很小,而這時(shí)的生產(chǎn)效率能提高到4倍。這是由于采用該方法能夠較好的保留地形特征信息。而采用格網(wǎng)或是基于距離及高差的抽稀算法當(dāng)抽稀率大于50%時(shí),其DEM精度變化幅度也較小。而基于系統(tǒng)的抽稀方法,整體上得到的DEM精度較其他方法低。同時(shí)可以從圖中看出,當(dāng)壓縮率大于90%時(shí),各種抽稀算法所得到的DEM精度基本相同。

圖11 各種抽稀算法的DEM精度與抽稀率

4 基于坡度抽稀算法應(yīng)用于委內(nèi)瑞拉數(shù)據(jù)實(shí)驗(yàn)

該實(shí)驗(yàn)區(qū)數(shù)據(jù)位于委內(nèi)瑞拉中部,屬于地勢(shì)平坦地區(qū),數(shù)據(jù)點(diǎn)個(gè)數(shù)為1 669 465個(gè)(如圖12所示),面積大小為420 m×8200 m,激光點(diǎn)平均間距為1 m,海拔高度在147.56～176.39 m之間。該實(shí)驗(yàn)區(qū)域數(shù)據(jù)的地形暈渲圖如圖13所示。

圖12 點(diǎn)云數(shù)據(jù)

圖13 該區(qū)域地形暈渲

為了能夠更快的處理數(shù)據(jù),在做實(shí)驗(yàn)之前采用TerraSolid軟件將該實(shí)驗(yàn)數(shù)據(jù)分成7塊,如圖14所示。其中A1塊的實(shí)驗(yàn)數(shù)據(jù)為176 199個(gè)點(diǎn),A2塊實(shí)驗(yàn)數(shù)據(jù)為207 201個(gè)點(diǎn),A3塊實(shí)驗(yàn)數(shù)據(jù)為215 372個(gè)點(diǎn),A4塊實(shí)驗(yàn)數(shù)據(jù)為216912個(gè)點(diǎn),A5塊實(shí)驗(yàn)數(shù)據(jù)為258 952個(gè)點(diǎn),A6塊實(shí)驗(yàn)數(shù)據(jù)為284 807個(gè)點(diǎn),A7塊實(shí)驗(yàn)數(shù)據(jù)為310022個(gè)點(diǎn)。

圖14 數(shù)據(jù)分塊

實(shí)驗(yàn)結(jié)果及精度分析見表7。

表7 各數(shù)據(jù)塊DEM精度與抽稀率

經(jīng)過(guò)抽稀后該實(shí)驗(yàn)區(qū)域的數(shù)據(jù)變成269 284個(gè)點(diǎn),壓縮率達(dá)到16.13%,而從表7的中誤差可以看出,各數(shù)據(jù)塊中誤差都低于10 cm。實(shí)驗(yàn)區(qū)域的工程精度要求是不低于10 cm,因此,對(duì)此數(shù)據(jù)壓縮到16.13%時(shí)仍能滿足精度要求。

抽稀后的實(shí)驗(yàn)數(shù)據(jù)如圖15與圖16所示。

圖15 壓縮后點(diǎn)云數(shù)據(jù)

圖16 壓縮后該區(qū)域地形暈渲

5 結(jié)論

機(jī)載激光雷達(dá)技術(shù)為生成高精度的DEM提供了一種有效的方法。然而,由于采用該技術(shù)生成的DEM數(shù)據(jù)量十分龐大,對(duì)數(shù)據(jù)的存儲(chǔ)和處理帶來(lái)不便。數(shù)據(jù)抽稀的目的是去除DEM數(shù)據(jù)中的冗余數(shù)據(jù),從而提高數(shù)據(jù)存儲(chǔ)和處理的效率。通過(guò)對(duì)幾種抽稀算法的實(shí)驗(yàn),表明基于坡度的抽稀算法能夠得到較好的結(jié)果。同時(shí)通過(guò)實(shí)驗(yàn)驗(yàn)證得出,由機(jī)載激光雷達(dá)數(shù)據(jù)構(gòu)建的DEM能夠承受一定的抽稀,而仍保持一定的精度。采用基于坡度的抽稀方法對(duì)委內(nèi)瑞拉數(shù)據(jù)進(jìn)行抽稀,實(shí)驗(yàn)結(jié)果表明,對(duì)此數(shù)據(jù)抽稀到16.13%時(shí),仍能滿足精度的要求。在DEM數(shù)據(jù)抽稀中,地形特征信息對(duì)描述地形表面特征具有重要的意義。因此,在DEM抽稀中如果能考慮地形特征信息,將會(huì)在減少數(shù)據(jù)量的同時(shí)而仍能保持較高的精度。

[1]李志林,朱 慶.數(shù)字高程模型[M].武漢：武漢大學(xué)出版社,2003

[2]劉經(jīng)南,張小紅.激光掃描測(cè)高技術(shù)的發(fā)展及現(xiàn)狀[J].武漢大學(xué)學(xué)報(bào)：信息科學(xué)版,2003,28(2)：132-137

[3]劉 春,吳杭彬.基于平面不規(guī)則三角網(wǎng)的DEM數(shù)據(jù)壓縮與質(zhì)量分析[J].中國(guó)圖像圖形學(xué)報(bào),2007,12(5)：836-840

[4]徐景中,萬(wàn)幼川,張圣望.LIDAR地面點(diǎn)云的簡(jiǎn)化方法研究[J].測(cè)繪信息與工程,2008,33(1)：32-34

[5]周興華,姚藝強(qiáng),趙吉先.DEM內(nèi)插方法與精度評(píng)定[J].測(cè)繪科學(xué),2005,30(5)：86-88

[6]張小紅.機(jī)載激光雷達(dá)測(cè)量技術(shù)理論與方法[M].武漢：武漢大學(xué)出版社,2007

[7]周啟鳴,劉學(xué)軍.數(shù)字地形分析[M].北京：科學(xué)出版社,2006

[8]Chou Y H,Liu PS,Dezzani RJ.Terrain complexity and reduction of topographic data.Geographical Systems,1999,1(2)：179-197

[9]Liu X Z,Peterson J,Chandra S.The effect of LiDARdatadensity on DEM accuracy.Proceedingsof International congress on modeling and simulation,Christchurch,New Zealand,2007：1363-1369

[10]Pamela S.A data density reduction algorithm for post-processed airborne LiDAR bathymetric survey data.University of Florida,2000