王　丽

(宿州学院环境与测绘工程学院，安徽 宿州 234000)



嫦娥一号异源数据适应精度加权插值融合生成DEM模型

王丽

(宿州学院环境与测绘工程学院，安徽 宿州 234000)

目的提出适应精度加权插值融合的理论模型，实现嫦娥一号CCD立体相机正视影像图制作的1∶250 万全月分幅数字高程模型数据和激光高度计LAM获取的高程数据两种异源数据的融合，生成月球虹湾地带30～40 °W、40～50 °N区域更精确的DEM模型。方法根据适应精度加权插值融合的理论模型，选取月球虹湾地带30～40 °W、40～50 °N区域作为研究对象，以异源数据为数据源，根据数据源不同精度分配权重，利用MATLAB编程，实现适应精度加权插值融合的算法，生成了月球虹湾地带30～40 °W、40～50 °N区域DEM模型。同时，根据克里金插值理论，生成相应区域DEM模型。结果通过改变适应精度加权插值融合的理论模型的权重参数，分别选取3个典型权重参数值，即分别分配激光高度计高程数据参数为0.9、0.5、0.1。结论克里金插值生成的DEM精度对比分析，说明适应精度加权插值融合可以生成更精确的DEM模型，是切实可行、有效的方法。

DEM数据；LAM数据；克里金插值；异源数据适应精度加权插值

1　引　言

2007年10月24日，中国发射了第一颗月球探测卫星，并且在其上搭载了CCD三线阵立体相机(TLS)和激光高度计(LAM)。嫦娥一号卫星利用CCD三线阵立体相机获取了分辨率120 m的前视、正视和后视影像数据[1]；利用其上搭载的激光高度计获取了包括月球极区在内的2 500轨影像数据[2]。由于CCD三线阵立体相机获取的影像数据平面精度高[3]、激光高度计获取的数据高程精度高[4-6]，本文根据这2种数据特点，并结合二者的互补性，选取嫦娥一号CCD立体相机正视影像图制作的1∶250万全月分幅数字高程模型数据，以及激光高度计LAM获取的高程数据作为两种异源数据源，利用本文提出的异源数据适应精度加权插值算法，结合MATLAB编程，实现以上2种异源数据的融合，生成了更为准确的DEM模型。

2　数据源

本文选取的异源数据为激光高度计获取的2B级数据以及由嫦娥一号CCD立体相机的正视影像图制作的1∶250万全月分幅数字高程模型数据。

本文研究区域：30～40 °W、40～50 °N，包括月球最为典型的虹湾地带。

3　异源数据适应精度加权插值

3.1原理

数据的空间插值方法有很多种，主要包括反距离权重法、最邻近点法、径向基函数法、趋势面法、多元回归法和克里金法等[7]。通过比较发现，克里金插值法综合了数据的空间相关性和随机性，比其他方法适用性更强，适用范围更广，且插值效果逼近程度更高[8，9]。因此，本文选用克里金插值法，在此基础上提出适用于异源数据的插值方法，使其权重的选取能根据精度的大小进行分配，对异源数据实现加权插值融合。

克里金插值法原理在此不再赘述。主要介绍在克里金插值方法的基础上提出的异源数据适应精度加权插值融合原理。

为了计算方便，将三线阵CCD影像生成的DEM数据的权重表示为λccd，激光高度计的数据权重表示为λlam，为了满足所有权值之和等于1，即λlam+λccd=1，本文计算每类数据的权重时所运用的方程组为：

(1)

(2)

由公式(1)和(2)知，在一定的搜索半径范围内，可计算出两类异源数据的各自权重值。实现对每个数据分配权重，实现异源数据的融合，即：

(3)

(4)

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

(5)

由公式(3)和(4)得到三线阵CCD影像生成的DEM数据和激光高度计LAM加权插值后的高程值，再根据公式(5)即可得到搜索区域内点插值融合后的高程值。

3.2步骤

异源数据适应精度加权插值融合步骤如下：

(1)确定网格点数，计算出每个网格节点的坐标；

(2)确定搜索半径，利用近点距离搜索法，找出搜索半径范围内的参估点，并对两类参估点的属性进行判断；

(3)计算变差函数值[10]，根据克里金插值方法变差函数的计算步骤，计算三线阵CCD和激光高度计的变差函数参数，求γ(xi,xj)和γ(xi,x0)；

(4)根据公式(1)和(2)，求2类数据各自在样本数据中的权重；

(5)根据公式(3)、(4)和(5)计算出待估点插值融合后的高程值；

(6)对每个格网重复以上步骤，直到所有的网格节点坐标全部求出，然后在MATLAB里显示出异源数据插值融合后的高程模型DEM图。

3.3实例分析

选取月球虹湾地带30～40 °W、40～50 °N区域为例。

3.3.1数据源获取与处理

激光高度计数据是2B级数据，无需处理，可以直接截取30～40 °W、40～50 °N区域的数据。

三线阵CCD正视影像图需在Arcgis软件里裁剪拼接出30～40 °W、40～50 °N区域(图1)，并输出其ASCII码形式的数据(图2)，便于提取经纬度、高程信息数据。

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码形式的数据

由于ASCII码数据不能直接运用，所以利用C++编程，将ASCII码数据输出成经度、纬度、高程值排列的形式，如图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影像获取的数据

3.3.2计算变差函数值

由原理可以看出，变差函数值的求取非常关键。变差函数值求取是否正确，决定了权值的求取。在求取变差函数值之前，需对数据源的空间特性作出判断，看数据源是否符合正态分布特性。

利用MATLAB软件编程，绘制出数据源的正态分布概率图和直方图，如图4所示。

由图4可以看出，数据源符合正态分布特征，可以求取变差函数值。变差函数值的求取过程，这里不再赘述。

分别求取了(0 °、45 °、90 °、135 °)4个方向上的变差函数，得到实验变差函数模型，如图5所示。

图4　数据源的正态分布概率图和直方图

图5　实验变差函数模型图

通过求取的实验变差函数，利用球状模型进行拟合[11]，得到4个方向上的理论变差函数，再运用几何各向异性套合公式，对4个方向上的理论变差函数套合，得到套合后的理论变差函数，如图6所示。

图6　套合后的理论变差函数

3.3.3异源数据适应精度加权插值融合

利用异源数据适应精度加权插值融合理论，结合求取的变差函数值，可以利用MATLAB程序实现权重根据精度不同进行分配，实现三线阵CCD的数字高程模型和激光高度计LAM的数字高程模型加权插值融合，生成高精度的DEM模型。同时根据普通克里金插值方法，对同一数据源进行插值，生成常规的DEM模型。

由于三线阵CCD影像生成的DEM数据的高程精度比激光高度计LAM的数字高程精度低，因此，以λccd=0.4，λlam=0.6为例，对其计算分配权重，实现异源数据适应精度加权插值融合。生成的DEM模型如图7(a)所示，利用克里金插值原理生成的DEM模型如图7(b)所示。

(a)异源数据加权插值融合后的DEM (b)克里金插值后的DEM

图7对比图

4　结论与讨论

根据图7，很难从直观上判断其精度的大小，所以从统计学的角度对其精度进行统计分析。

实验中因为没有卫星的轨道数据，无法利用三线阵CCD的影像数据生成DEM，只能通过对三线阵CCD的DEM模型数据加粗差降低精度，得到加粗差的三线阵CCD影像生成的DEM数据和激光高度计LAM数据两种异源数据。并对这2种数据进行插值，分析其精度。实验中，分别将加粗差的三线阵数据权重λccd赋予0.1、0.5、0.9。故激光高度计LAM数据权重λlam值分别是0.9、0.5、0.1，加权插值融合后得到的高程差值如表1所示。

表1　加权插值融合后得到的高程差值

表1中的高程差值前3行值是由2种异源数据加权插值融合生成DEM高程值和三线阵CCD影像生成的DEM数据运用普通克里金插值生成DEM高程值相减所得值，第四行克里金插值所得到的高程差值是由异源数据利用普通克里金插值生成的DEM高程值和三线阵CCD影像生成的DEM数据生成的DEM高程值相减所得到的值。

从表1可以看出：一方面，随着高程精度较高的激光高度计LAM权重的增大，中误差值越来越小，说明增加精度较高数值的权重，适应精度加权插值融合插值生成的DEM模型的精度越高；另一方面，利用适应精度加权插值融合插值方法所得到的中误差比利用普通克里金插值方法得到的中误差小，说明异源数据适应精度加权插值融合精度有所提高。综上所述，本文提出的异源数据适应精度加权插值融合方法是切实可行和有效的。

5　结束语

本文根据数据源精度的不同进行分配权重，在普通克里金插值原理的基础上，提出了一种能适应数据源精度进行加权插值融合方法，和普通克里金插值方法的主要区别就是在权重的分配上。本文中根据数据的精度将权重按比例分配，只是在初始分配时有很大的主观因素，但是通过对实验结果的分析，说明该方法比普通克里金插值的效果好、精度高，并且随权重的不同插值的精度不同。与此同时，本研究也存在不足，由于没有卫星的轨道数据，因此不能利用三线阵CCD的影像数据生成DEM，实现其精度分析。

[1]李春来,刘建军,任鑫,等.嫦娥一号图像数据处理与全月球影像制图[J].中国科学:地球科学,2010,40(03):294-306.

[2]蔡占川,郑才目,唐泽圣，等.基于嫦娥一号卫星激光测高数据的月球DEM及高程分布特征模型[J].中国科学:技术科学,2010,40(11):1300-1311.

[3]崔腾飞,陈圣波,王景然.基于“嫦娥”卫星三线阵CCD立体相机的月球表面三维建模[J].国土资源遥感,2009,(04):31-34.

[4]邹小端,刘建军,任鑫,等.利用嫦娥一号激光高度计数据制作月球DEM的方法研究[J].测绘科学,2011,36(01):13-16.

[5]李春来,任鑫,刘建军,等.嫦娥一号激光测距数据及全月球DEM模型[J].中国科学：地球科学,2010,40(03):281-293.

[6]李春来.嫦娥一号三线阵CCD数据摄影测量处理及全月球数字地形图[J].测绘学报,2013,42(06):853-860.

[7]张靖.GIS空间内插方法与应用研究[D].呼和浩特：内蒙古师范大学,2008.

[8]靳国栋,刘衍聪,牛文杰.距离加权反比插值法和克里金插值法的比较[J].长春工业大学学报:自然科学版,2003,24(03):53-57.

[9]苏姝,林爱文,刘庆华.普通Kriging法在空间内插中的运用[J].江南大学学报:自然科学版,2004,3(01):18-21.

[10]陈天伟,卢献健,杨子江.应用遗传算法解算克里金插值权系数[J].北京测绘,2015,(04):30-33.

[11]张小艳,谭勇.变差函数球状模型的自动拟合与实现[J].物探与化探,2010,34(02):253-257.

[责任编辑：关金玉英文编辑：刘彦哲]

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)；自然灾害过程与防控研究安徽省重点实验室开放基金资助项目(ZH201307)；安徽省质量工程项目(教高【2008】4号24);宿州学院科研平台项目(2012YKF12)

王丽(1987-)，女，安徽安庆人，助教，硕士，主要从事摄影测量与遥感方面研究。

P 236

A

来稿日期：2015-12-24