邢鵬威,唐詩(shī)華,張 曦,張 躍,何廣煥
(1.廣西空間信息與測(cè)繪重點(diǎn)實(shí)驗(yàn)室,廣西 桂林 541006;2.桂林理工大學(xué)測(cè)繪地理信息學(xué)院,廣西 桂林 541006;3.河北工程大學(xué)礦業(yè)與測(cè)繪工程學(xué)院,河北 邯鄲 056038)
隨著無(wú)人機(jī)行業(yè)的快速崛起,無(wú)人機(jī)在我國(guó)智慧城市、智慧工廠、智慧礦山等諸多領(lǐng)域發(fā)揮出必不可少的優(yōu)勢(shì)。基于水利、建筑等行業(yè)對(duì)高精度DEM(Digital Elevation Model)的需求[1],DEM的數(shù)據(jù)獲取方式也隨之改變。目前,常采用機(jī)載Lidar(Light detection and ranging)、機(jī)載相機(jī)等方法采集數(shù)據(jù)。機(jī)載Lidar數(shù)據(jù)雖然光斑小可以穿透植被,但是雷達(dá)點(diǎn)云密度較無(wú)人機(jī)點(diǎn)云稀疏,且雷達(dá)點(diǎn)云無(wú)顏色信息[2];而無(wú)人機(jī)影像匹配點(diǎn)云數(shù)據(jù)具有顏色(RGB)信息和數(shù)據(jù)量大的特點(diǎn),但是匹配點(diǎn)云無(wú)強(qiáng)度與回波信息,點(diǎn)云會(huì)受到地表植被的影響,經(jīng)分類(lèi)后的點(diǎn)云易產(chǎn)生漏洞,從而影響DEM的精度。為此,本文提出一種基于LS-SVM的無(wú)人機(jī)匹配點(diǎn)云修補(bǔ)方法,對(duì)漏洞周邊點(diǎn)云進(jìn)行提取,通過(guò)LS-SVM方法對(duì)濾波后點(diǎn)云數(shù)據(jù)漏洞進(jìn)行預(yù)測(cè),通過(guò)擬合預(yù)測(cè)的點(diǎn)云數(shù)據(jù)生成DEM,從而解決去噪后無(wú)人機(jī)匹配點(diǎn)云影響DEM模型精度的問(wèn)題。
本文采用不規(guī)則三角網(wǎng)(TIN)對(duì)高山區(qū)點(diǎn)云數(shù)據(jù)進(jìn)行處理。首先,采用人工目譯的方法去除異常噪聲點(diǎn),再以測(cè)區(qū)內(nèi)最大建筑物為依據(jù),將測(cè)區(qū)內(nèi)點(diǎn)云進(jìn)行網(wǎng)格劃分,以每個(gè)網(wǎng)格中最低點(diǎn)的種子點(diǎn)云構(gòu)成三角網(wǎng)[3]。設(shè)定其他點(diǎn)到三角網(wǎng)的角度和長(zhǎng)度閾值,判別是否為地面點(diǎn),并將判斷好的地面點(diǎn)不斷計(jì)入三角網(wǎng)進(jìn)行迭代計(jì)算,直到將所有點(diǎn)云判別完畢。TIN降噪后的點(diǎn)云見(jiàn)圖1。
圖1 TIN降噪匹配點(diǎn)云數(shù)據(jù)
再通過(guò)植被指數(shù)剔除未過(guò)濾掉的植被點(diǎn),獲得更精細(xì)的地面點(diǎn)點(diǎn)云。濾除植被后點(diǎn)云數(shù)據(jù)見(jiàn)圖2。
圖2 去植被匹配點(diǎn)云數(shù)據(jù)
Kriging插值算法是以函數(shù)加區(qū)域變量為基礎(chǔ),應(yīng)用于空間數(shù)據(jù)分布的自然問(wèn)題[4]。Kriging的多種優(yōu)化算法都是基于微小變量,公式如下
(1)
Kriging插值方法對(duì)求解權(quán)重λi極為重要,權(quán)重應(yīng)滿足無(wú)偏估計(jì),使估計(jì)方差達(dá)到最優(yōu)。該算法無(wú)偏估計(jì)表達(dá)式為
(2)
估計(jì)方差表達(dá)式為
(3)
式中,E為數(shù)學(xué)期望;C為協(xié)方差函數(shù),運(yùn)用協(xié)方差函數(shù)和變異函數(shù)之間的關(guān)系,以拉格朗日乘法為衡量標(biāo)準(zhǔn),構(gòu)建函數(shù),建立n+1維線性方程組,即
(4)
將已知點(diǎn)數(shù)據(jù)代入式(4)即可求得λi。采用Kriging局部插值算法獲取的成果見(jiàn)圖3。
圖3 基于Kriging插值修補(bǔ)成果
反距離加權(quán)算法(Inverse Distance Weighted,IDW)以權(quán)重為核心,插值點(diǎn)的大小會(huì)被采樣點(diǎn)影響,針對(duì)未知點(diǎn)的距離權(quán)重,以權(quán)函數(shù)為依據(jù)[5],其權(quán)函數(shù)公式如下
(5)
式中,λj為j點(diǎn)權(quán)重;dj為j點(diǎn)和待預(yù)測(cè)點(diǎn)間的距離;u為權(quán)指數(shù)。權(quán)重是以距離為自變量的衰減函數(shù),采樣點(diǎn)與待預(yù)測(cè)點(diǎn)之間距離越大,權(quán)重越小,當(dāng)預(yù)測(cè)點(diǎn)距離采樣點(diǎn)一定距離以外,可忽略不計(jì)權(quán)重。各個(gè)采樣點(diǎn)權(quán)重的和便是每1個(gè)預(yù)測(cè)點(diǎn)的數(shù)值,公式為
(6)
式中,zj為j點(diǎn)的高程值;zp為預(yù)測(cè)點(diǎn)的高程值;d-u為距離衰減函數(shù)。
由于預(yù)測(cè)點(diǎn)的值取決于不同距離采樣點(diǎn)的權(quán)重,是一種精確的插值算法,因此,IDW易受采樣點(diǎn)精度的影響,使采樣點(diǎn)附近模型出現(xiàn)凸凹現(xiàn)象。采用IDW插值獲取的成果見(jiàn)圖4。
圖4 基于IDW插值修補(bǔ)成果
改進(jìn)Shepare插值算法是一種標(biāo)準(zhǔn)導(dǎo)數(shù)距離加權(quán)過(guò)程[6],加權(quán)公式如下
(7)
式中,wi為權(quán)重;di為i點(diǎn)到插值點(diǎn)的距離;r為調(diào)整距離。改進(jìn)Shepard算法調(diào)整權(quán)重的形式有2種:
(1)在給定半徑范圍內(nèi),權(quán)重wi根據(jù)最遠(yuǎn)點(diǎn)的距離來(lái)調(diào)整,若范圍內(nèi)距預(yù)測(cè)點(diǎn)最遠(yuǎn)的距離為r,則修正的距離倒數(shù)加權(quán)公式為
(8)
(9)
式中,u為權(quán)指數(shù);δ為平滑因子。當(dāng)δ=0時(shí)為精確性插值算法,δ≠0時(shí)為非精確插值算法。
(2)權(quán)重wi結(jié)合了局部二次多項(xiàng)式,以擬合二次多項(xiàng)式修正高程值為參與倒數(shù)加權(quán)函數(shù)的高程值,公式為
(10)
式中,zj為預(yù)測(cè)點(diǎn);Qi為二次多項(xiàng)式函數(shù)。采用改進(jìn)Shepard插值算法獲取的成果見(jiàn)圖5。
圖5 改進(jìn)Shepare插值修補(bǔ)成果
徑向基函數(shù)插值(Radial Basis Function,RBF)算法是用多個(gè)曲面的線性組合去逼近真實(shí)曲面。RBF算法由2部分組成,公式為
(11)
式中,zp為預(yù)測(cè)點(diǎn)高程;λi為第i點(diǎn)權(quán)重值;di為預(yù)測(cè)點(diǎn)到第i采樣點(diǎn)的距離;φ(di)為徑向基函數(shù);fj(x)為趨勢(shì)函數(shù),是次數(shù)小于m的基本多項(xiàng)式。
采用徑向基函數(shù)插值算法可以選取不同種類(lèi)的徑向基函數(shù)φ(di),如多重二次曲面、倒數(shù)多重二次曲面等對(duì)離散數(shù)據(jù)進(jìn)行插值[7]。與IDW插值算法相比,IDW插值以采樣點(diǎn)作為依據(jù)進(jìn)行插值,不能預(yù)測(cè)超出采樣點(diǎn)范圍的數(shù)值,而徑向基函數(shù)算法通過(guò)趨勢(shì)函數(shù)可以計(jì)算超出采樣點(diǎn)范圍的數(shù)值。此次試驗(yàn)以薄板樣條函數(shù)為徑向基函數(shù),對(duì)去噪后匹配點(diǎn)云數(shù)據(jù)進(jìn)行插值,公式為
(12)
式中,c為光滑因子。采用RBF算法獲取的成果見(jiàn)圖6。
圖6 基于RBF插值修補(bǔ)成果
LS-SVM是一種等式約束代替不等式約束,避免了在對(duì)偶空間求解二次規(guī)劃問(wèn)題的算法,在計(jì)算速度上有著極大優(yōu)越性,與其他方法相比更適用于動(dòng)態(tài)模型的構(gòu)建。同時(shí),也繼承了支持向量機(jī)模型(SVM)解決非線性、局部極值等問(wèn)題的功能,這使其在函數(shù)估計(jì)和預(yù)測(cè)逼近中得到了廣泛應(yīng)用[8-9]。運(yùn)用已知數(shù)據(jù)非線性映射到特征空間構(gòu)造函數(shù),即
g(x,w)=wTφ(x)+b
(13)
式中,g為線性回歸函數(shù);w為權(quán)值向量;φ(x)為非線性函數(shù)映射;b為閾值。LS-SVM的最小化問(wèn)題見(jiàn)下式
(14)
其約束條件為
s.t.yi=wTφ(xi)+b+ξi
(15)
式中,w為法向量;C為正則化參數(shù);ξi為誤差項(xiàng)。對(duì)于式(14)、(15)同時(shí)構(gòu)造拉格朗日函數(shù),得
(16)
式中,ai為拉格朗日乘子。根據(jù)ai與KKT(Karush Kuhn Tucker)條件進(jìn)行計(jì)算,同時(shí)引入k(xi,xj)核函數(shù)得到優(yōu)化的線性方程,即
(17)
(18)
采用LS-SVM算法插值獲取數(shù)據(jù)結(jié)果見(jiàn)圖7。
圖7 基于LS-SVM算法修補(bǔ)成果
在分析5組插值后的DEM精度時(shí),以測(cè)區(qū)實(shí)測(cè)CORS實(shí)時(shí)動(dòng)態(tài)測(cè)量數(shù)據(jù)為檢查點(diǎn),檢查點(diǎn)點(diǎn)位分布見(jiàn)圖8。將同時(shí)在測(cè)區(qū)測(cè)量的CORS實(shí)時(shí)動(dòng)態(tài)測(cè)量高程數(shù)據(jù)作為真值Zi,取插值試驗(yàn)數(shù)據(jù)為預(yù)測(cè)值Ri,以Z方向偏差、平均誤差MAE、均方差RMSE、擬合優(yōu)越度R2等交叉驗(yàn)證指標(biāo)評(píng)判插值算法精度。公式如下
(19)
圖8 檢查點(diǎn)分布
表1 5種插值方法插值精度統(tǒng)計(jì)
為了更加清晰分辨數(shù)據(jù)精度,將Kriging算法、IDW算法、RBF算法與LS-SVM算法的檢查點(diǎn)殘差4組數(shù)據(jù)殘差生成折線圖,由于改進(jìn)Shepard插值數(shù)據(jù)結(jié)果RMSE與R2較差,不能滿足模型精度需求,故不繪制其殘差折線圖。4種插值算法殘差見(jiàn)圖9。從圖9可知,以改進(jìn)Shepard插值算法的數(shù)據(jù)誤差較大,最大誤差值達(dá)到-1.151 m,將大量誤差較大數(shù)據(jù)剔除后,DEM模型產(chǎn)生豎向空洞與豎向凸起,其擬合優(yōu)越度較低(R2=0.992 4);IDW算法插值數(shù)據(jù)的采樣點(diǎn)誤差最大為-0.735 m,遠(yuǎn)大于CH/T 9008.2—2010《基礎(chǔ)地理信息數(shù)字成果1∶500、1∶1 000、1∶2 000數(shù)字搞成模型》規(guī)定的采樣點(diǎn)最大誤差應(yīng)為中誤差2倍(0.735 m>2×0.306 m);RBF與Kriging插值精度均劣于IDW算法,同時(shí)2組結(jié)果的MAE、RMSE與R2基本一致;LS-SVM插值數(shù)據(jù)相比較前4組數(shù)據(jù)而言精度明顯提高,LS-SVM插值的MAE=-0.148 m、RMSE=0.250 m、最大殘差0.447 m。
圖9 不同修補(bǔ)方法檢查點(diǎn)殘差折線
針對(duì)無(wú)人機(jī)匹配點(diǎn)云數(shù)據(jù)降噪后的漏洞修補(bǔ)問(wèn)題,本文提出了基于LS-SVM的無(wú)人機(jī)匹配點(diǎn)云數(shù)據(jù)修補(bǔ)法。與其他幾種常規(guī)插值算法相比,預(yù)測(cè)數(shù)據(jù)與實(shí)際測(cè)區(qū)高程更接近,能較好地還原實(shí)際地貌特征。經(jīng)此方法修補(bǔ)后數(shù)據(jù)的R2=0.999 5,擬合優(yōu)越度較高,RMSE=0.250 m,最大殘差0.447 m,滿足山地RMSE≤0.5 m、高山地RMSE≤1.0 m和采樣點(diǎn)數(shù)據(jù)最大誤差小于2倍RMSE的規(guī)范要求,滿足山地、高山地1∶500A級(jí)數(shù)字高程數(shù)字成果生產(chǎn)要求。同時(shí),還能減少外業(yè)人員對(duì)漏洞的補(bǔ)測(cè)工作量以及在行業(yè)設(shè)計(jì)過(guò)程中土方工程量與水庫(kù)蓄水量等相關(guān)設(shè)計(jì)誤差。