基于DEM的火山口提取的研究

■于华良 潘军 祝铭 刘兴春 张魁珂

（吉林大学地球探测科学与技术学院 吉林长春130026）

基于DEM的火山口提取的研究

■于华良 潘军 祝铭 刘兴春 张魁珂

（吉林大学地球探测科学与技术学院 吉林长春130026）

数字高程模型即DEM，是地表形态高程属性的数字化表达。本文基于二次高程异常的方法对DEM数据进行火山地形提取的实验。DEM中火山和平原区域的处理以及DEM分辨率的大小都直接影响着火山口提取的质量和效率。本文利用目视解译、坡度分析重分类以及Laplacian算子运算，并结合熔岩台地地貌分类提取出火山口。对长白山地区进行试验，分析结果证明了此方法的可行性。

DEM火山口个数Laplacian模板分析

1　引言

DTM即数字地形模型，是基于数学模型建立起来的地表地形起伏模型。DTM是以高等数学和计算机等学科为依托，随着这些学科和计算机成图的发展成长起来的。数字地形模型是将地形的真实的空间展布情况利用数字的形式表现出来。

数字高程模型DEM是构成DTM的基础，是区域地形表面简单的、离散的数字表示模型,它由规则水平间隔处地面点的抽样高程矩阵组成，一般以栅格数据的形式表示。DEM作为数字化的地形图，蕴含着大量的、各种各样的地形结构和特征信息，是定量描述地貌结构、水文过程、生物空间分布等空间变化的基础数据。对于地貌的相关研究工作，极为重要的一种手段就是地形分析。利用遥感图像配合DEM进行地形地貌的解译的方法，可以用来整合数据，提高遥感图像和DEM数据所隐含的地貌信息的提取。本文研究目的是提取长白山地区的火山口，将采用DEM数据的地形来提取。

在归纳遥感与DEM的特点运用于地形地貌方面，程维明等学者利用SRTMDEM与遥感数据探究提取基本地貌类型的方法，基本地貌类型被由SRTMDEM派生的各种模型自动提取出来，而且具有很大的实效性。利用高程模型及起伏度模型相结合的方式，初步实现了基本地貌类型的自动提取，对区域地貌宏观规律有较好的表达。2003年EndreD等利用遥感和DEM综合解译制作了区域土壤图，2004年刘学军等人对提取坡度、坡向的算法在DEM数据上进行的试验、研究并对得到的精度进行了分析,2009年周增波，程维明，周成虎等在长白山靖宇县地区运用DEM数据进行火山口地貌识别，并对提取和识别结果进行探。2010年王博，罗微进行了基于Aster G-DEM的海南岛地形地貌提，龙恩，张洪岩，程维明等进行了基于遥感与DEM的长白山区数字地貌制图方法的研。不同比例尺的DEM对提取地面坡度精度有较大影响，汤国安就这一问题进行研究并取得了一定结果。2003年刘泽慧提出利用DEM数据辅助山脊线、山谷线提取的方法。王培法等人利用不同比例尺和不同空间分辨率提取出了流域面积、河网密度、河道坡度、河道长度、平均高程、平均坡度和河道长度并做出了详细的分。

本研究以DEM地形因子的划分为主要技术手段对长白山地区进行研究，使用的是Landsat ETM+数据和30米分辨率DEM数据。采用基于二次微分的高程异常的方法，最后利用遥感图像对长白山地区的火山口个数进行研究和探讨。

2　提取方法

二次微分高程异常提取方法主要包括目视解译、Laplacian模版分析运算以及地形起伏度的计算。目视解译需要统计出火山口的大致半径。模板的如何选取是由火山口的形状大小确定的，需要将研究区域范围内火山口的相关信息进行统计。作为边缘增强滤波，Laplacian滤波在运用过程中可以忽略边缘的方向。经典的边缘检测算子的基本原理是考察图像上每个像素的相应邻域内的灰度值的变化情况，可以利用一阶或二阶方向导数变化规律，以此获得边缘信息。经典的边缘算子是通过各向同性的拉普拉斯算子演变出来的，最简单的各向同性微分算子是拉普拉斯算子。一个二维图像函数f(x,y)的拉普拉斯算子定义为：

因为对任意阶微分都是进行线性的操作，所以拉普拉斯变换也是一个线性算子。为了离散形式描述这一公式，可以用差分的形式来定义。在x方向上，

同理在y方向上，

由上述三个公式，可以推导出含两个变量的离散Laplacian算子的形式如下

Laplacian滤波强调图像中的最大值，它通过运用一个具有最高中心值的变换核来完成，所有的Laplacian滤波卷积核的维数都必须是奇数。下面是3*3的模版卷积核及5*5的卷积核。由这种规律性可以推导出更高阶的卷积核。

DEM地形起伏度可以在一定程度上反映出地貌的特征，地形起伏度是区域内最大高程与最小高程之差。寻找这个固定区域空间单位是关键。在研究区内以3*3的模版大小为起始依次增加到59*59来计算模版窗口的最大高差，再拟合成数学图像找到最大高差所对的模版大小，确定空间单位的大小。

其数学模型为：

熔岩台地在DEM上可以根据其地形起伏度来定义，长白山地区熔岩台地的地形起伏度在30m和70m之。

3　火山口提取实验及步骤

3.1研究区概况

本研究区位于吉林省东部长白山地区内,其地理坐标为东经127°～129°、北纬41°～43°。研究区面积4008.00km2。以长白山火山主峰为中心。

研究区主要地貌类型为火山地貌、构造剥蚀地貌和流水地貌。以长白山天池主峰为中心大致呈同心环状划分为巨型复式火山锥、山麓倾斜熔岩高原和熔岩台地三大地貌单元。

在火山活动时期，岩浆由火山口喷出。火山口的位置主要处于火山锥顶端，还有一部分位于其侧面，当岩浆不再喷发时，其流经的通道中的岩浆冷却，由此造成的环状凹陷的区域即为火山口。

在DEM上火山口通常是类似于山顶点但是又不同于山顶点这主要取决于其形状。长白山地区的火山口形态各异，一部分火山口为四周高中间低的形态特征，另一部分则表现出中间高四周低的形态特征，但是在数学模拟图像上其二阶导数不等于零，也就是说火山口的实际形态是上下起伏的，本文正是基于这一原理进行提取。

具体实验步骤如下：

①对DEM进行预处理后与遥感图像匹配融合进行目视解译；

②利用坡度重分类后与原图像相乘进行加强处理；

③进行Laplacian变换得到图像；

④统计分析出火山口的值的范围；

⑤通过波段运算进行进一步提取；

⑥对DEM进行地形起伏度的分析提取出熔岩台地；

图3-1　长白山地区的ETM+数据7、4、2波段合成影像

⑦提取后图像进行叠加从而精确结果，最后得出结论。

3.2数据源及预处理

研究采用的遥感数据为1999年9月2日的长白山地区ETM+遥感图像，ETM+遥感数据有7个波段，其中第6波段是热红外波段，故将其去掉，保留其余的6个波段。图3-1长白山地区的ETM+数据7、4、2波段合成影像。

图3-2　长白山地区DEM预处理图像

本研究采用的DEM数据为ASTER DEM数据分别是：41°N 127°E、41°N 128°E、42°N 127°E、42°N 128°E这四景图像，然后进行融合，镶嵌边界忽略零值。在运算中，如果单一考虑坡度因素，那么会出现许多噪声即零星的像元。为消除零星像元带来的干扰，可以在生成坡度之前，要对DEM进行预处理。如图3-2为长白山地区DEM预处理图像。

图3-3　目视解译效果图图像

3.2.2目视解译

在Arcgis中把landsat ETM+遥感图像与处理后的DEM图形进行校正和配准然后进行裁剪然后进行目视解译，在Arcgis中新建一个矢量文件以表示火山口。统计出火山口的个数，统计出大致的大小。

统计栅格的大小从而来确定模版的选择，把矢量数据转化成栅格文件，改为30米的分辨率和图像进行配准比较得到正确的统计值，经统计分析后得出最大的涵盖范围从而确定选择11*11的窗口进行试验。

3.3坡度分析加强处理

图3-4　坡度分析加强处理后图像

由于要突出火山口的地形坡度，所以要先分析出区域中DEM的坡度然后进行重分类，根据相关文献可知最佳坡度为3.5°[2]。重分类后的图像与原DEM相乘以突出平原和山地的差别，从而增强火山口的形状使其更容易识别。

3.4运用模版来提取

在envi下进行laplacian变换处理图像，通过试验后选择11*11的窗口。处理后进行训练样本的统计分析，从中分析出处理后火山口的最大值、最小值、均值。进行波段运算（B1 GT 10352 AND B1 LT 15923603）*B1+(B1 GT 15923603 OR B1 LT 10352)*0后设立阈值提取火山口。兴趣区统计量如下：

3.5用DEM地表起伏度因子来对的熔岩台地进行提取

在研究区内根据相关文献可知最佳的地形起伏度窗口选择为19*19，从而确定空间单位的大小。对DEM图像进行重分类，进而提取出熔岩台地来。

3.6叠加分析出火山口个数

在Arcgis中对上面两个结果图层进行叠加分析，然后分析提取出火山口，火山口的个数为184个。与目视解译的结果进行拓扑检查，确定提取后的点与原来目视解译的火山口重叠个数为106。

图3-5　处理以后的图像

图3-6　为熔岩台地的提取图

4　结果分析

本文在以长白山地区为背景研究提取火山口的方法，通过目视解译1999年的ETM+图像可以大体知道火山口个数为132个。二次微分的高程异常的方法经过Laplacian算子提取后与地貌分类相叠加提取火山口结果为184个，提取后与目视解译的火山口进行拓扑分析，得到重合的火山口为106个，其精度是在可以接受的范围内。某些分布于熔岩台地边沿或者分布在其他复杂的地理位置的火山口，是造成误差的主要原因，另一方面，有些火山口的规模较小，无法被提取是造成误差的另一原因。以上这些都给提取增加了难度，使得结果存在一定误差，但是总体上来说，此方法提取火山口的精度仍是较为准确的。

图3-7　叠加后提取火山口图像

5　结论

对于某种地貌的提取，如火山口，根据其不同于其他地质结构的特点，火山口的实际形态是上下起伏的，在数学模拟图像上其二阶导数不为零；以及火山口在遥感和DEM图像上的形态，最终确定相关参数的选取。利用模板提取和DEM地表起伏度因子提取熔岩台地的结果相叠加，最终确定火山口的数目和位置。本文采用二次高程异常方法进行提取，主要包括目视解译、坡度分析重分类以及Laplacian算子运算。在长白山地区这种方法的应用效果较好。证明二次高程异常方法具有较高的实用价值，在今后的工作中要继续深入研究，将此方法的准确性进一步提高。

[1]刘学军，王叶飞，等，基于DEM的坡度坡向误差空间分布特征研究 [B]，测绘通报，2004.12:22-28.

[2]周增波，程维明，周成虎等，基于DEM的火山口地貌识别方法探讨 [J]，地球信息科学报，2009.11(6):773-776.

[3]王博，罗微，基于Aster G-DEM的海南岛地形地貌提取。硕士学位论文，海南大学，2010.

[4]龙恩，张洪岩，程维明，基于遥感与DEM的长白山区数字地貌制图方法研究，硕士学位论文，东北师范大学，2005.

[5]汤国安，杨勤科，不同比例尺DEM提取地面坡度的精度研究-以在黄土丘陵沟壑区的试验为例 [J]。水土保持通报，2001.01:54-55.

[6]刘泽慧，黄培之，DEM数据辅助的山脊线和山谷线提取方法的研究测绘科学 [A]，山地学报，2003.12:34-36.

[7]汤国安，杨玮莹.一种基于DEM的明暗等高线制图方法 [B]，测绘通报，2001.07: 39-40.

[8]丁峰.高志海.魏怀东.DEM在旱耕地荒漠化评价中的应用-以武威、古浪二市 (县)为例 [J].水土保持学报，2000.14(3):31-35.

[9]王培法，栅格DEM的尺度与水平分辨率对流域特征提取的分析-以黄土岭流域为例 [J]，江西师范大学学报 (自然科学版)2004.06.

[10]张会平，张恒安，杨农等，基于GIS的岷江上游地貌形态初步分析 [J]，中国地质灾害与防治学报，2004.9:116-119.

[11]高杨，吕宁，薛重生，不同比例尺数字高程模型对土壤侵蚀强度分级的影响 [J]。中国水土保持，2007,(10):26-28.

[12]周启鸣，刘学军，数字地形分析，科学出版社，2008.

[13]基于遥感和DEM数据的柬埔寨洞里萨湖地区数字地貌制图研究，硕士学位论文，南京师范大学，2008.

P208[文献码]B

1000-405X（2016）-2-405-3

于华良（1990～），男，硕士研究生，研究方向为遥感解译与制图。

潘军（1971～），男，副教授，研究方向为遥感与地理信息系统教学和科研。

吉林大学研究生创新基金资助项目，项目号：2015137