基于高精度DEM的地形监测方法研究

，， ，

(1.中国地质大学(北京)，北京 100083； 2.北京市国土资源信息研究开发重点实验室，北京 100083)

0 引 言

随着社会的发展，人类对能源的需求量也越来越大，对地球能源的开发与开采等，势必造成地球原有地形的变化，地形的变化又会引发滑坡、崩塌等众多不良连锁反应，给人类带来灾害。因此，对于地形的监测也逐渐受到人们的关注与重视。

近些年，遥感技术的发展不仅表现在光谱信息成像化,雷达成像多极化,光学探测多向化等光学方面(郭凯等,2012)，通过遥感手段获取DEM(数字高程模型)的非光学遥感技术也发展迅速。

1 高精度DEM的获取与处理

DEM数据包括平面位置和高程数据2种信息。为了建立DEM,需量测一些点的三维坐标,这就是DEM数据采集(周启鸣等,2006)。随着科学技术的发展，计算机技术、测量技术和传感器质量都得到很大的提高，促使遥感技术上升到了一个新的台阶。DEM相关技术及应用也迎来了新的契机，主要表现在精度更高、空间分辨率更好、测量效率更高、经济快速等方面。目前，高精度DEM的主要获取方式有激光雷达(LiDAR)、合成孔径雷达(InSAR)、数字摄影测量3种方式。以下简单介绍高精度DEM获取的方法及数据处理流程。

1.1 激光雷达(LiDAR)

LiDAR(Light Detection And Ranging)，即光探测与测量。自激光发出，到达建筑物、植被等顶端，返回第一次值，直至最后一次到达地面返回最后一次值，从而可以计算出每一次激光返回时的高度，而最后一次激光返回时的高度即为裸露地表的高程，还可以通过第一次和最后一次返回值所得的高度之差得到植被和地面建筑等的高度。

机载LiDAR是一种主动式对地观测系统，是20世纪90年代初首先由西方国家发展起来并投入商业化应用的一门新兴技术。它集成激光测距技术、计算机技术、惯性测量单元(IMU)/DGPS差分定位技术于一体，该技术在三维空间信息的实时获取方面取得了重大突破，为获取高时空分辨率地球空间信息提供了一种全新的技术手段。LiDAR具有高精度、高密度、高效率、高分辨率的优势，一般平面精度可以达到0.1～0.5 m。

LiDAR获取的是点云数据，在对点云数据进行处理前，首先应当对数据进行检验，剔除错误的点和高程异常的点，例如特别低的点(地面以下的点)或特别高的点(云或飞行中的鸟)。然后建立地面数据类型，分离出非地面数据。再根据要求，对非地面数据进行分类，如建筑物类数据、输电线路数据、铁路类数据等。另外，可以根据高出地面的高度对非地面点数据进行分类，如可以对不同的植被进行分类。然后输出分类的数据，通过滤波处理，得到DEM等(图1)。

图1 LiDAR数据处理流程图

1.2 合成孔径雷达(InSAR)

合成孔径雷达(InSAR)技术是依靠SAR成像处理后的单视复数影像，将其相位信息作为信息源，以提取地表三维信息或地表变化信息(王超等，2002)。InSAR 通过雷达的2副天线同时对地观测或雷达2次近似平行的对地观测来获取地面同一地区的干涉对，利用干涉对间微小的差别，在复图像上产生相位差，通过图像处理，即可形成干涉纹图。斜距向上的点与两天线位置之差等信息在干涉纹图上也有所体现，利用传感器高度、波束视向、天线基线距及雷达波长之间的几何关系，就可准确测出图像上每一点的三维信息和变化信息(Tarolli et al,2009)。

为了从干涉图中提取高度信息，InSAR技术需要经过选择包含共同观测区的具有一定相关性的2幅SAR单视复影像，通过图像的精确匹配将两图像匹配到亚像元级，然后通过复数据共轭相乘得到干涉图，由于得到的干涉图的相位值不仅包含高程的信息，还包含由平地效应引起的斜距信息，所以必须进行去平地效应，然后再进行干涉图的噪声滤除。又因为干涉图的相位差为缠绕相位，因此还要进行相位解缠，这样才能得到非模糊相位。然后进行相高转换将相位转化为高程信息，最后再进行地理编码，将图像从斜距坐标转化为地面坐标，即可得到观测区的地形图,其数据处理流程如图2。

图2 InSAR数据处理流程

1.3 数字摄影测量

数字摄影测量是DEM数据采集常用的方法之一。利用附有的自动记录装置(接口)的立体测图仪或立体坐标仪、解析测图仪及数字摄影测量系统，进行人工、半自动或全自动的量测来获取数据。立体像对是由不同摄站获取的同一地区的具有一定重叠度的2张影像，其实质是传感器模拟人的双眼，构成人体立体视觉，即将相对上的视差反应为人眼的生理视差后得到立体视觉，一般立体相对的重叠度至少为60%。立体相对影像存在高差引起的投影差，同一高程的地物在2张影像上投影差不同，当用两像片进行立体观察时将引起视差，产生立体感。

数字摄影测量法提取DEM的精度很大程度上取决于控制点和同名点的精度及分布情况，要保证和提高遥感影像成像模型的精度，需要更多的控制点，而同名点要求在整个工作区域均匀分布(图3)。

图3 立体像对处理流程

2 基于高精度DEM的地形监测

DEM是通过数值阵列所表示的一组实地模型，为数字地形模型中的一个分支(李志林等, 2000)，主要描述各种地形因子，例如地面的坡度、坡向等非线性或线性组合而成的空间分布情况。DEM有着定量化以及定位的特点(任莉莉，2014)，更高精度的DEM使得地形定量化分析及监测更加的精确、准确，高精度DEM在地形监测方面有着传统低分辨率DEM无法比拟的优势。以下选取某地不同分辨率DEM进行对比说明。

但实际上，王维不仅是一位诗人，还是一位画家。尽管目前公认出自他手的画作还没有一幅能确定是他的真迹，但这对于我们这并不重要，因为我们真正关注的，是那画中如诗般的意境。

图4中，可以依次把1,10,30 m 3种不同分辨率的DEM称为高精度、较高精度、一般精度。从图中不难看出高精度DEM显示更为细腻，地形细节显示更为清晰明了；对于较高精度DEM，已可以较清晰地看出栅格像元，地形细节也不够凸显；一般精度DEM，可以明显看到栅格像元，已分辨不出地形细节。

图4 不同分辨率DEM对比图

2.1 山体阴影

山体阴影是模拟太阳光照射地形所引起的明暗对比，对地形图进行渲染，使之看起来更有立体效果。山体阴影通过为栅格中的每个像元确定照明度，来获取表面的假定照明度。通过设置假定光源的位置和计算与相邻像元相关的每个像元的照明度值，即可得出假定照明度。进行分析或图形显示时，特别是使用透明度时，“山体阴影”工具可大大增强表面的可视化。如图5中，通过对3种不同分辨率的DEM进行山体阴影渲染后，增强了其立体效果，地形特征更为突出。同时可以看出，高精度DEM地形特征的凸显效果更好，较高精度次之，一般精度分辨不出任何地形特征。

图5 不同分辨率山体栅格阴影图

2.2 坡度

坡度表示的地面的陡缓程度，坡度为每个像元高程值在从该像元到与其相邻的像元方向上的最大变化率，即是DEM的一阶导数。坡度值是地形起伏度的定量化表现，坡度值越大表示地形起伏越大，其地形走势越陡峭，如山区地势陡峭，坡度值则比较大；反之，坡度值越小则表示地形起伏越小，地形平缓，如平原地区地势平缓，坡度值则比较小。

坡度是滑坡发育的内在因素之一(郭果等,2013)，同时，滑坡的发生也伴随着坡度值的改变，并且滑坡区域一般坡度值都比较高，图6是基于1 m DEM计算所得坡度图，图中红色为高值，绿色为低值。可以看出，该区北边为高值区，山地地形；中间区域总体为低值，地势较平缓，然而图中蓝色框体中又有着比较高的坡度值，结合影像数据，对蓝框中感兴趣的区进行对比分析，初步判定为滑坡体，进一步的确认有待野外考察验证。

综合上述分析可知坡度对地形监测的重要性，不仅能表示地形的陡缓程度，还在滑坡的人机交互提取中起着重要作用。

图6 滑坡区坡度图

2.3 地形曲率

地形曲率是表达地形曲面结构的主要参数之一，也是地表过程模拟、土壤侵蚀模型、土地利用分类等环境模型的基本变量。在地貌、地理、土壤、资源与环境等领域有着重要的应用(刘学军等,2006)。地形曲率用于描述流域盆地的物理特征，从而便于理解侵蚀过程和径流形成过程。坡度会影响下坡时的总体移动速率。坡向将决定流向。剖面曲率将影响流动的加速和减速，进而将影响到侵蚀和沉积。

地形表面曲率是地形曲面在各个截面方向上的形状、凹凸变化的反映，他们是平面点位的函数。地形表面曲率反映了地形结构和形态，地形曲面曲率的计算与地形曲面的二阶导数有关，在DEM上的实现通常要通过曲面拟合或差分的方式实现。

图7是基于1 m DEM计算所得的最小曲率和最大曲率图，使用窗口为5×5。最小曲率图中，蓝色表示低值，红色表示高值，不难发现低值区域成线状分布，结合影像可以发现低值区所表示的是凹陷型线状地形，初步判定为水沟等。最大曲率图中，蓝色表示低值，红色表示高值，而不同于最小曲率的是此处高值区域成线状分布，结合影像高值区所表示的是突出型线状地形，初步判定为山脊线。对于最小曲率和最大曲率所展示出的线状地形的类型，只能结合影像作初步判断，进一步的验证需要野外调查。

图7 最小曲率与最大曲率图

地形曲率中，表面曲率反映的是地形在各个截面方向上的形状、凹凸变化。最小、最大曲率则与地形中线状地形相关，在线状地形的自动识别与提取方面起着重要作用。

3 结 论

通过计算山体阴影、坡度和曲率，概述了基于高精度DEM的地形监测方法：通过坡度的计算，可以提取滑坡、崩塌等；通过地形曲率的计算，可以提取水沟、山脊线等线状地形，从而达到地形监测的目的。

研究所采用的是同一时间、不同空间分辨率的DEM数据，所做的地形监测方法研究均是基于单期DEM所进行的，并没有实现对地形的动态监测。在今后的工作中，可以收集多期不同时相的高精度DEM数据，计算不同期次DEM的坡度、地形曲率。通过不同期次坡度、地形曲率差分等方法的综合分析，从而实现地形的动态监测。

白贵霞,程传录,郭春喜,等.2010.汶川地震地形形变监测与分析[J].地球科学与环境学报,32(2):205-210.

陈强.2006.基于永久散射体雷达差分干涉探测区域地表形变的研究[D].成都:西南交通大学.

郭凯,吕勇.2012.微波遥感技术在雷达中的应用探究[J].信息通信,(3):53-53.

郭果,陈钧,李明惠,等.2013.土质滑坡发育概率与坡度间关系研究术[J].工程地质学报,21(4):607-612.

匡薇.2014.Radarsat-2 数据处理方法研究及其在地面沉降中的应用:以阜阳地区为例[D].北京:中国地质大学(北京).

李志林,朱庆.2000.数字高程模型[M].武汉:武汉大学出版社.

刘学军,王叶飞,曹志东,等.2006.基于DEM的地形曲率计算模型误差分析[J].测绘科学,31(5):50-53

任莉莉.2014.DEM上地学分析方法研究与实现[J].黑龙江科技信息,(11):32-32.

王超,张红,刘智.2002.星载合成孔径雷达干涉测量[M].北京:科学出版社.

巫谢香.2014.我国遥感技术的运用与发展趋势[J].读写算：教育教学研究,(2):87-88.

薛强,张茂省,朱立峰,等.2013.基于多期DEM数据的滑坡变形定量分析[J].地质通报,32(6):935-942.

周启鸣,刘学军.2006.数字地形分析[M].北京:科学出版社.

张丽丽,王小平,张璞.2009.基于无人机影像生产高精度DEM的实践[J].测绘技术装备,11(1):33-34.

EVANS I S.1980.An integrated system of terrain analysis and slopemapping[J].Zeitschrift für Geomorphologic Suppl-Bd,36:274-295.

TAROLLI P,FONTANA G D.2009.Hillslope to valley transition morphology: new opportunities from high resolution DTMs[J].Geomorphology,113(1/2): 47-56.

TAROLLI P,SOFIA G,FONTANA G D.2012.Geomorphic features extraction from high resolution topography:landslide crowns and bank erosion [J].Natural Hazards,61(1):65-83.