嫦娥一號異源數(shù)據(jù)適應(yīng)精度加權(quán)插值融合生成DEM模型

王　麗

(宿州學(xué)院環(huán)境與測繪工程學(xué)院，安徽 宿州 234000)

?

嫦娥一號異源數(shù)據(jù)適應(yīng)精度加權(quán)插值融合生成DEM模型

王麗

(宿州學(xué)院環(huán)境與測繪工程學(xué)院，安徽 宿州 234000)

目的提出適應(yīng)精度加權(quán)插值融合的理論模型，實現(xiàn)嫦娥一號CCD立體相機(jī)正視影像圖制作的1∶250 萬全月分幅數(shù)字高程模型數(shù)據(jù)和激光高度計LAM獲取的高程數(shù)據(jù)兩種異源數(shù)據(jù)的融合，生成月球虹灣地帶30～40 °W、40～50 °N區(qū)域更精確的DEM模型。方法根據(jù)適應(yīng)精度加權(quán)插值融合的理論模型，選取月球虹灣地帶30～40 °W、40～50 °N區(qū)域作為研究對象，以異源數(shù)據(jù)為數(shù)據(jù)源，根據(jù)數(shù)據(jù)源不同精度分配權(quán)重，利用MATLAB編程，實現(xiàn)適應(yīng)精度加權(quán)插值融合的算法，生成了月球虹灣地帶30～40 °W、40～50 °N區(qū)域DEM模型。同時，根據(jù)克里金插值理論，生成相應(yīng)區(qū)域DEM模型。結(jié)果通過改變適應(yīng)精度加權(quán)插值融合的理論模型的權(quán)重參數(shù)，分別選取3個典型權(quán)重參數(shù)值，即分別分配激光高度計高程數(shù)據(jù)參數(shù)為0.9、0.5、0.1。結(jié)論克里金插值生成的DEM精度對比分析，說明適應(yīng)精度加權(quán)插值融合可以生成更精確的DEM模型，是切實可行、有效的方法。

DEM數(shù)據(jù)；LAM數(shù)據(jù)；克里金插值；異源數(shù)據(jù)適應(yīng)精度加權(quán)插值

1　引　言

2007年10月24日，中國發(fā)射了第一顆月球探測衛(wèi)星，并且在其上搭載了CCD三線陣立體相機(jī)(TLS)和激光高度計(LAM)。嫦娥一號衛(wèi)星利用CCD三線陣立體相機(jī)獲取了分辨率120 m的前視、正視和后視影像數(shù)據(jù)[1]；利用其上搭載的激光高度計獲取了包括月球極區(qū)在內(nèi)的2 500軌影像數(shù)據(jù)[2]。由于CCD三線陣立體相機(jī)獲取的影像數(shù)據(jù)平面精度高[3]、激光高度計獲取的數(shù)據(jù)高程精度高[4-6]，本文根據(jù)這2種數(shù)據(jù)特點，并結(jié)合二者的互補(bǔ)性，選取嫦娥一號CCD立體相機(jī)正視影像圖制作的1∶250萬全月分幅數(shù)字高程模型數(shù)據(jù)，以及激光高度計LAM獲取的高程數(shù)據(jù)作為兩種異源數(shù)據(jù)源，利用本文提出的異源數(shù)據(jù)適應(yīng)精度加權(quán)插值算法，結(jié)合MATLAB編程，實現(xiàn)以上2種異源數(shù)據(jù)的融合，生成了更為準(zhǔn)確的DEM模型。

2　數(shù)據(jù)源

本文選取的異源數(shù)據(jù)為激光高度計獲取的2B級數(shù)據(jù)以及由嫦娥一號CCD立體相機(jī)的正視影像圖制作的1∶250萬全月分幅數(shù)字高程模型數(shù)據(jù)。

本文研究區(qū)域：30～40 °W、40～50 °N，包括月球最為典型的虹灣地帶。

3　異源數(shù)據(jù)適應(yīng)精度加權(quán)插值

3.1原理

數(shù)據(jù)的空間插值方法有很多種，主要包括反距離權(quán)重法、最鄰近點法、徑向基函數(shù)法、趨勢面法、多元回歸法和克里金法等[7]。通過比較發(fā)現(xiàn)，克里金插值法綜合了數(shù)據(jù)的空間相關(guān)性和隨機(jī)性，比其他方法適用性更強(qiáng)，適用范圍更廣，且插值效果逼近程度更高[8，9]。因此，本文選用克里金插值法，在此基礎(chǔ)上提出適用于異源數(shù)據(jù)的插值方法，使其權(quán)重的選取能根據(jù)精度的大小進(jìn)行分配，對異源數(shù)據(jù)實現(xiàn)加權(quán)插值融合。

克里金插值法原理在此不再贅述。主要介紹在克里金插值方法的基礎(chǔ)上提出的異源數(shù)據(jù)適應(yīng)精度加權(quán)插值融合原理。

為了計算方便，將三線陣CCD影像生成的DEM數(shù)據(jù)的權(quán)重表示為λccd，激光高度計的數(shù)據(jù)權(quán)重表示為λlam，為了滿足所有權(quán)值之和等于1，即λlam+λccd=1，本文計算每類數(shù)據(jù)的權(quán)重時所運用的方程組為：

(1)

(2)

由公式(1)和(2)知，在一定的搜索半徑范圍內(nèi)，可計算出兩類異源數(shù)據(jù)的各自權(quán)重值。實現(xiàn)對每個數(shù)據(jù)分配權(quán)重，實現(xiàn)異源數(shù)據(jù)的融合，即：

(3)

(4)

Z*(x0)=Z*(xccd0)+Z*(xlam0)

(5)

由公式(3)和(4)得到三線陣CCD影像生成的DEM數(shù)據(jù)和激光高度計LAM加權(quán)插值后的高程值，再根據(jù)公式(5)即可得到搜索區(qū)域內(nèi)點插值融合后的高程值。

3.2步驟

異源數(shù)據(jù)適應(yīng)精度加權(quán)插值融合步驟如下：

(1)確定網(wǎng)格點數(shù)，計算出每個網(wǎng)格節(jié)點的坐標(biāo)；

(2)確定搜索半徑，利用近點距離搜索法，找出搜索半徑范圍內(nèi)的參估點，并對兩類參估點的屬性進(jìn)行判斷；

(3)計算變差函數(shù)值[10]，根據(jù)克里金插值方法變差函數(shù)的計算步驟，計算三線陣CCD和激光高度計的變差函數(shù)參數(shù)，求γ(xi,xj)和γ(xi,x0)；

(4)根據(jù)公式(1)和(2)，求2類數(shù)據(jù)各自在樣本數(shù)據(jù)中的權(quán)重；

(5)根據(jù)公式(3)、(4)和(5)計算出待估點插值融合后的高程值；

(6)對每個格網(wǎng)重復(fù)以上步驟，直到所有的網(wǎng)格節(jié)點坐標(biāo)全部求出，然后在MATLAB里顯示出異源數(shù)據(jù)插值融合后的高程模型DEM圖。

3.3實例分析

選取月球虹灣地帶30～40 °W、40～50 °N區(qū)域為例。

3.3.1數(shù)據(jù)源獲取與處理

激光高度計數(shù)據(jù)是2B級數(shù)據(jù)，無需處理，可以直接截取30～40 °W、40～50 °N區(qū)域的數(shù)據(jù)。

三線陣CCD正視影像圖需在Arcgis軟件里裁剪拼接出30～40 °W、40～50 °N區(qū)域(圖1)，并輸出其ASCII碼形式的數(shù)據(jù)(圖2)，便于提取經(jīng)緯度、高程信息數(shù)據(jù)。

ncols2688nrows2467xllcorner-41.204867401891yllcorner39.626919257964cellsize0.0044906302304878NODATA_Value-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999-9999

圖2ASCII碼形式的數(shù)據(jù)

由于ASCII碼數(shù)據(jù)不能直接運用，所以利用C++編程，將ASCII碼數(shù)據(jù)輸出成經(jīng)度、緯度、高程值排列的形式，如圖3。

-37.455250.0541-2692-37.450750.0541-2692-37.275650.0541-2689-37.271150.0541-2689-37.266650.0541-2689-37.095150.0541-2751-37.08750.0541-2751-36.934350.0541-2603-36.929850.0541-2603-36.925350.0541-2603-36.920850.0541-2603-36.916350.0541-2603

圖3三線陣CCD影像獲取的數(shù)據(jù)

3.3.2計算變差函數(shù)值

由原理可以看出，變差函數(shù)值的求取非常關(guān)鍵。變差函數(shù)值求取是否正確，決定了權(quán)值的求取。在求取變差函數(shù)值之前，需對數(shù)據(jù)源的空間特性作出判斷，看數(shù)據(jù)源是否符合正態(tài)分布特性。

利用MATLAB軟件編程，繪制出數(shù)據(jù)源的正態(tài)分布概率圖和直方圖，如圖4所示。

由圖4可以看出，數(shù)據(jù)源符合正態(tài)分布特征，可以求取變差函數(shù)值。變差函數(shù)值的求取過程，這里不再贅述。

分別求取了(0 °、45 °、90 °、135 °)4個方向上的變差函數(shù)，得到實驗變差函數(shù)模型，如圖5所示。

圖4　數(shù)據(jù)源的正態(tài)分布概率圖和直方圖

圖5　實驗變差函數(shù)模型圖

通過求取的實驗變差函數(shù)，利用球狀模型進(jìn)行擬合[11]，得到4個方向上的理論變差函數(shù)，再運用幾何各向異性套合公式，對4個方向上的理論變差函數(shù)套合，得到套合后的理論變差函數(shù)，如圖6所示。

圖6　套合后的理論變差函數(shù)

3.3.3異源數(shù)據(jù)適應(yīng)精度加權(quán)插值融合

利用異源數(shù)據(jù)適應(yīng)精度加權(quán)插值融合理論，結(jié)合求取的變差函數(shù)值，可以利用MATLAB程序?qū)崿F(xiàn)權(quán)重根據(jù)精度不同進(jìn)行分配，實現(xiàn)三線陣CCD的數(shù)字高程模型和激光高度計LAM的數(shù)字高程模型加權(quán)插值融合，生成高精度的DEM模型。同時根據(jù)普通克里金插值方法，對同一數(shù)據(jù)源進(jìn)行插值，生成常規(guī)的DEM模型。

由于三線陣CCD影像生成的DEM數(shù)據(jù)的高程精度比激光高度計LAM的數(shù)字高程精度低，因此，以λccd=0.4，λlam=0.6為例，對其計算分配權(quán)重，實現(xiàn)異源數(shù)據(jù)適應(yīng)精度加權(quán)插值融合。生成的DEM模型如圖7(a)所示，利用克里金插值原理生成的DEM模型如圖7(b)所示。

(a)異源數(shù)據(jù)加權(quán)插值融合后的DEM (b)克里金插值后的DEM

圖7對比圖

4　結(jié)論與討論

根據(jù)圖7，很難從直觀上判斷其精度的大小，所以從統(tǒng)計學(xué)的角度對其精度進(jìn)行統(tǒng)計分析。

實驗中因為沒有衛(wèi)星的軌道數(shù)據(jù)，無法利用三線陣CCD的影像數(shù)據(jù)生成DEM，只能通過對三線陣CCD的DEM模型數(shù)據(jù)加粗差降低精度，得到加粗差的三線陣CCD影像生成的DEM數(shù)據(jù)和激光高度計LAM數(shù)據(jù)兩種異源數(shù)據(jù)。并對這2種數(shù)據(jù)進(jìn)行插值，分析其精度。實驗中，分別將加粗差的三線陣數(shù)據(jù)權(quán)重λccd賦予0.1、0.5、0.9。故激光高度計LAM數(shù)據(jù)權(quán)重λlam值分別是0.9、0.5、0.1，加權(quán)插值融合后得到的高程差值如表1所示。

表1　加權(quán)插值融合后得到的高程差值

表1中的高程差值前3行值是由2種異源數(shù)據(jù)加權(quán)插值融合生成DEM高程值和三線陣CCD影像生成的DEM數(shù)據(jù)運用普通克里金插值生成DEM高程值相減所得值，第四行克里金插值所得到的高程差值是由異源數(shù)據(jù)利用普通克里金插值生成的DEM高程值和三線陣CCD影像生成的DEM數(shù)據(jù)生成的DEM高程值相減所得到的值。

從表1可以看出：一方面，隨著高程精度較高的激光高度計LAM權(quán)重的增大，中誤差值越來越小，說明增加精度較高數(shù)值的權(quán)重，適應(yīng)精度加權(quán)插值融合插值生成的DEM模型的精度越高；另一方面，利用適應(yīng)精度加權(quán)插值融合插值方法所得到的中誤差比利用普通克里金插值方法得到的中誤差小，說明異源數(shù)據(jù)適應(yīng)精度加權(quán)插值融合精度有所提高。綜上所述，本文提出的異源數(shù)據(jù)適應(yīng)精度加權(quán)插值融合方法是切實可行和有效的。

5　結(jié)束語

本文根據(jù)數(shù)據(jù)源精度的不同進(jìn)行分配權(quán)重，在普通克里金插值原理的基礎(chǔ)上，提出了一種能適應(yīng)數(shù)據(jù)源精度進(jìn)行加權(quán)插值融合方法，和普通克里金插值方法的主要區(qū)別就是在權(quán)重的分配上。本文中根據(jù)數(shù)據(jù)的精度將權(quán)重按比例分配，只是在初始分配時有很大的主觀因素，但是通過對實驗結(jié)果的分析，說明該方法比普通克里金插值的效果好、精度高，并且隨權(quán)重的不同插值的精度不同。與此同時，本研究也存在不足，由于沒有衛(wèi)星的軌道數(shù)據(jù)，因此不能利用三線陣CCD的影像數(shù)據(jù)生成DEM，實現(xiàn)其精度分析。

[1]李春來,劉建軍,任鑫,等.嫦娥一號圖像數(shù)據(jù)處理與全月球影像制圖[J].中國科學(xué):地球科學(xué),2010,40(03):294-306.

[2]蔡占川,鄭才目,唐澤圣，等.基于嫦娥一號衛(wèi)星激光測高數(shù)據(jù)的月球DEM及高程分布特征模型[J].中國科學(xué):技術(shù)科學(xué),2010,40(11):1300-1311.

[3]崔騰飛,陳圣波,王景然.基于“嫦娥”衛(wèi)星三線陣CCD立體相機(jī)的月球表面三維建模[J].國土資源遙感,2009,(04):31-34.

[4]鄒小端,劉建軍,任鑫,等.利用嫦娥一號激光高度計數(shù)據(jù)制作月球DEM的方法研究[J].測繪科學(xué),2011,36(01):13-16.

[5]李春來,任鑫,劉建軍,等.嫦娥一號激光測距數(shù)據(jù)及全月球DEM模型[J].中國科學(xué)：地球科學(xué),2010,40(03):281-293.

[6]李春來.嫦娥一號三線陣CCD數(shù)據(jù)攝影測量處理及全月球數(shù)字地形圖[J].測繪學(xué)報,2013,42(06):853-860.

[7]張靖.GIS空間內(nèi)插方法與應(yīng)用研究[D].呼和浩特：內(nèi)蒙古師范大學(xué),2008.

[8]靳國棟,劉衍聰,牛文杰.距離加權(quán)反比插值法和克里金插值法的比較[J].長春工業(yè)大學(xué)學(xué)報:自然科學(xué)版,2003,24(03):53-57.

[9]蘇姝,林愛文,劉慶華.普通Kriging法在空間內(nèi)插中的運用[J].江南大學(xué)學(xué)報:自然科學(xué)版,2004,3(01):18-21.

[10]陳天偉,盧獻(xiàn)健,楊子江.應(yīng)用遺傳算法解算克里金插值權(quán)系數(shù)[J].北京測繪,2015,(04):30-33.

[11]張小艷,譚勇.變差函數(shù)球狀模型的自動擬合與實現(xiàn)[J].物探與化探,2010,34(02):253-257.

[責(zé)任編輯：關(guān)金玉英文編輯：劉彥哲]

Generated DEM Model of Adaptive Weighted Interpolation Fusion Based on Chang’e-1 Heterologous Data Accuracy

WANG Li

(School of Environment and Surveying Engineering,Suzhou University,Suzhou,Anhui 234000,China)

ObjectiveTo propose the theoretical model of adaptive precision-weighted interpolation fusion,to implement heterologous data fusion of 1∶2500000 digital elevation model data produced by fore-sight images shot by Chang’e-1 CCD stereocamera and elevation data obtained by laser altimeter,and to generate more accurate digital elevation model of Sinus Iridum area at 30-40°W,40°-50°N.Methods With the theoretical model of weighted interpolation fusion,Sinus Iridum area at 30°-40°W,40°-50°N was selected as a research subject,and Heterologous Data as data resources.Weight was distributed weight based on the differerence of data accuracy.The MATLAB programming was used to complete the calculation process and the digital elevation model of Sinus Iridum area was generated at 30°-40°W,40°-50°N.Based on the Kriging interpolation theory meantime,the digital elevation model of related area was produced.ResultsVia revising the weighting parameter of the theoretical model of weighted interpolation fusion,three typical weighting parameters were selected to distribute the DEM of laser altimeter as 0.9,0.5 and 0.1.ConclusionCompared with the digital elevation model’s precision based on Kriging interpolation theory,the method of adaptive precision-weighted interpolation fusion could generate more precise DEM model,which showed that weighted interpolation fusion method was feasible and effective.

DEM data;LAM data;Kriging interpolation;heterologous data weighted interpolation fusion according to their accuracy

10.3969/j.issn.1673-1492.2016.05.001

安徽省精品資源共享課程(2013gxk097)；自然災(zāi)害過程與防控研究安徽省重點實驗室開放基金資助項目(ZH201307)；安徽省質(zhì)量工程項目(教高【2008】4號24);宿州學(xué)院科研平臺項目(2012YKF12)

王麗(1987-)，女，安徽安慶人，助教，碩士，主要從事攝影測量與遙感方面研究。

P 236

A

來稿日期：2015-12-24