dem制作流程(dem建立流程)
4D(Dsm,Dom,Dem,Dlg)的产品技术流程。
简单洁说,4D就是DOM(数字正射影像图)、DEM(数字高程模型)、DRG(数字栅格地图)、DLG(数字线划地图),其产品技术流程过于复杂,福州三维前沿。恰好最近有一个学习会。
DEM地形操作(geotools方式与nga方式)
由于项目中需要用到大范围tiff的图像(全中国30米分辨率的dem影像),并且需要单点获取高程,以及实现部分范围的dem裁切与获取趋于范围极值,当时在网上查找的部分,很多都不满足预期,或者计算结果与实际并不够契合,因此单独开一篇专门讲这块内容。
引入基础的依赖包,使用中间组件完成所需功能。
由于组件之间存在交互关系,因此将其全部放到一起。
由于接下去的操作都需要依赖同一份数据源,因此在最开始的时候,需要先封装关于数据源的相关操作。
通过获取gridCoverage2D 的的投影信息,加载投影并使用经纬度进行计算位置点高程
总体来说,整个的实现流程基本上可以概括为[引入依赖]-[加载数据源]-[业务数据读取]。读取DEM的高程流程相对简单,获取范围内的高程极值操作也不复杂,主要是需要算出范围对应的格网,以及获取每个格网像元的高程值;较为复杂的是DEM的图像导出,由于使用geotools导出的dem图像存在一定的偏移与投影问题,导致其他平台软件无法正常读取dem的一些属性信息,因此在操作上,先使用geotools自定义的导出,先输出一份DEM,再使用正常平台导出的dem文件,通过读取它的属性信息,赋值给空的dem影像,再将像元的高程值写入,最后输出成一份符合规范的dem。此操作较为繁琐,但基础满足功能所需,此为抛砖引玉,如有更好的办法实现,望诸位不吝赐教,将感激不尽!
遥感与DEM集成处理
遥感与地质地理信息的集成可制作三维影像地质图、构建虚拟地质环境。传统的地质图是由地质专业信息和地理信息组成的平面图,图注、图例纷繁复杂,专业性很强,制约了包括政府部门在内的其他各行业非地质专业人员对地质成果的积极和有效利用。三维影像地质图能够将地质调查成果及认识直观、鲜明、通俗易懂地表达出来,是使用者所喜闻乐见的一种表达方式,它不仅能服务于地质专业,以直观地分析地质构造规律,而且能被社会各个方面的用户所接受和利用(李建星等,2003),对遥感信息产业将起到有力的促进作用。
制作三维影像地质图的技术流程大致要经过数据源的准备、纹理的生成和DEM与纹理的套合整饰等步骤(图8.3)。数据源的准备工作见第二章相关内容。纹理的生成是在MapGIS6.5平台上实现的,即先将相山地区地质图的所有区都设置为透明输出,再利用MapGIS的光栅输出功能,通过依次添加已进行融合的遥感影像图(如SIFM543.img)和地质构造图,生成新的JPEG格式的平面影像地质图,并在图像处理模块中转换成MSI图像,这即是所需要的纹理。DEM与纹理的套合是在电子沙盘模块中进行的。
图8.3 三维影像地质图制作流程
图8.4即为相山地区形象直观、通俗易懂的三维影像地质图,相山火山-侵入杂岩的地形地貌景观跃入眼中。
arcgis地形预处理工具在吗
本章导读:ArcGIS的水文分析工具是基于DEM进行地表水流动的模拟,其本身不涉及到精确数值的水流流量。在形成径流的过程中考虑的全是地形因素,D8单流向算法决定了其必须针对无凹陷的DEM数据才能正确的分析出结果。在上一章中笔者已经就如何制作适合于ArcGIS水文分析的DEM数据。在本章中笔者将会介绍如何针对已有的DEM进行预处理,制作出适合水文分析的无凹陷DEM数据。
在ArcGIS原生的水文分析工具中包含两个常用的地形分析工具:汇、填洼。
###汇
汇是指流向栅格中流向无法被赋予八个有效值之一的一个或一组空间连接像元。汇被视为具有未定义的流向,并被赋予等于其可能方向总和的值。例如,如果最陡下落及其产生的流向都是向右 (1) 和向左 (16),则会分配值 17 作为该像元的流向。
在介绍ArcGIS流向分析算法的时候那个流向九宫格其定义是2的n次方,所以,通过汇的计算之后形成的栅格数据,可以通过其值反演出汇所在的方位。对于分析人员来说没什么作用,但对于后续编写程序的人员去改进水文分析是非常重要的。
这一段话中包含了非常大的信息量:
1、产生汇的主要原因之一就是DEM数据制作的时候采样效果和高程值设置为整数时导致的。这就是笔者在上一章中详细介绍制作适合于水文分析的DEM的原因。有些参数的设置在测绘行业上生成DEM基本上是无所谓的,但对水文分析来说就非常重要的。如果还有一些疑问,可以翻看上一章《ArcGIS水文分析实战教程(3)DEM数据准备》。如果能在DEM制作一环上已经对DEM数据进行控制,那是效果最好的。
2、冰川和喀斯特地形不适合使用ArcGIS水文分析。这主要是D8算法的原因。
3、像元大小小于10米的DEM数据很少自然产生汇。如果这种精度的数据都产生了汇,一定要检查这些汇是不是在现实中存在,通过叠加一些地形数据可以对比查看。否则就是数据中存在一些致命的错误。
4、像元大小增大,汇也随之增加。所以精度太差的数据其实是不适合于做小河流或小流域的分析,因为汇的大量存在基本上上会导致一部分的径流断流,由于汇入大江大河的径流会比较多,所以,针对一定级别的河流,还是具有参考价值。
以下是汇的剖面图
由于汇的地势都低于周边8个区域,所以水流都会汇入其中,导致径流最终断流。但现实中就算出现一些小的洼地,只要降水充足,这些洼地都会被填平,填满后径流将会继续往外流出。所以,D8算法就是假设有无限的降水,雨水不断的在地表形成径流。
###填洼
填平汇的过程就是填洼。ArcGIS中的填洼工具,使用与焦点流、流向、汇、分水岭和区域填充等工具等效的功能来定位和填充汇。该工具的执行过程会进行迭代,直到指定 z 限制内的所有汇均填充完毕。在填充汇的同时,可能会在填充区域的边界处创建其他汇,这些汇将在下个迭代中移除。
接下来看看填洼工具的一些参数,如下图
填洼工具必填的两个参数,其一是DEM,这个DEM有可能是原始的DEM,还有可能是经过填洼之后还发现汇的DEM(也就是n次填洼后的DEM),其二是输出的路径。关于Z值限制,是个可选项,估计大部分人都看不懂这个可选项到底代表什么意思,到底怎么设置。这个Z值限制在下面再做介绍。
先来看一张非常经典的水分分析的流程图,这张图似乎被引用过无数次,但大部分人只是做了一个单向的操作。如下图
这张流程图左边部分是制作无凹陷点DEM数据的流程。这个过程其实是一个循环的过程,不断的填洼并重新判断是否存在汇,如果还存在汇,则继续填洼,直到出现无凹陷点DEM。从图上来看,这最起码是个多次填洼的过程,而大部分人的做法就是不管汇是否存在,做一次填洼后再进行后续的水文分析。
精确填洼
直接使用填洼工具对DEM数据进行填洼,不管汇的情况,确实可以形成无凹陷点的DEM。但前面也提及到,如果只是负责填,但不是精确的填,大量的凹陷点被填平(这些凹陷点不一定就都是汇),会导致径流的线形形状与实际不符,在弯曲细节等会存在一定的损失。所以,如果水文分析人员要研究非常准确的地形与河流的关系,而不是单纯的河网拓扑关系,精确填洼就非常重要了。
不管是网上关于ArcGIS水文分析的教程还是ArcGIS官方的帮助,都没有具体说明使用【汇】工具之后找到汇后作什么样的处理。大部分人的处理更是简单粗暴,有汇即填洼,也不管具体哪里需要怎么填洼,填洼完之后直接进入分析环节。
笔者花了不少时间来研究其汇和填洼的原理,总结出以下步骤:
1、寻找汇
寻找汇相对来说比较简单,先通过DEM计算流向,然后直接使用【汇】工具进行查找,最终汇得出汇的栅格数据,如下图,使用的是ArcGIS spatial中自带的DEM数据进行分析。
图上叠加了第一次查找汇的栅格,发亮部分的点就是凹陷点,如果放大一点显示,那么可以清楚看到其像元形状,如下图
计算出Z值限制
z 限制指定凹陷点深度和倾泻点间的最大允许差值并确定要填充的凹陷点和保持不变的凹陷点。z 限制并非要填充的最大深度。
从描述来看,Z值限制要做的事情就是确保填洼工具不要一刀切,不要把所有的低洼点都填充,因为有些不是真的凹陷点;另一作用就是符合条件的洼地都会被指定填充。
例如,假设一个凹陷点区域中倾泻点的高程为 210 英尺,凹陷点的最深点为 204 英尺(相差 6 英尺)。如果将 z 限制设置为 8,则会填充该特殊凹陷点。但是,如果将 z 限制设置为 4,则不会填充该凹陷点,因为该凹陷点的深度超过该限制值,将其视为有效凹陷点。
填洼
填洼是地形处理中必不可少的一环。如果按照工具论,直接使用ArcGIS的填洼工具的默认配置,可以对DEM数据进行无限量的填洼。但笔者不推荐这种方式,而且Esri也不推荐。
正确的填洼流程应该是先要计算出Z值限制,然后根据一个或者一组Z值限制进行精确填洼。关于Z值限制前面已经说得很清楚,在这个阈值范围内汇才会被填充。那么Z值限制怎么计算出来呢?在【汇】的原理里面,帮助有具体提及到,但都是使用arcpy脚本的形式进行计算。以下引用ArcGIS帮助的Z值计算流程:
使用汇创建通过深度进行编码的汇的栅格。
输入流向栅格数据:flowdir
输出栅格:汇点
使用分水岭为每个汇创建汇流区域栅格。
输入流向栅格数据:flowdir
输入栅格数据或要素类倾泻点数据:汇点
输出栅格:sink_areas
将分区统计与最小值统计数据结合使用,以在每个汇的分水岭中创建最小高程的栅格。
输入栅格数据或要素区域数据:sink_areas
区域字段:值
输入赋值栅格:高程
输出栅格:sink_min
统计类型:MINIMUM
使用区域填充在每个汇的分水岭中创建最大高程的栅格。
输入区域栅格数据:sink_areas
输入权重栅格数据:高程
输出栅格:sink_max
使用减将最大值减去最小值以查找深度。
输入栅格数据 1:sink_max
输入栅格数据 2:sink_min
输出栅格:sink_depth
ArcPy的代码如下
为了与默认的填洼工具的结果进行对比,笔者特意针对该流程,使用modelbuilder进行了建模,以求出Z值限制。建模流程如下图
最后求出的是一组Z值限制,最大的Z值为32,最小是1,共有20处限制值。如下图所示
笔者做了几组测试
测试1:
对原始DEM数据进行无Z值限制填洼,直接执行填洼工具,使用默认设置,不设置Z值;
利用填洼完之后的DEM数据计算流向,并执行【汇】工具,发现所有的汇已经被填充,再也没有汇的存在了。如下图
测试2
对原始DEM数据进行有Z值限制的填洼,将Z值设置为33(之前已经计算出Z值最大值为32)。同样,对填洼过的DEM进行流向计算,并查找汇,发现存在汇,如下图
如果按照帮助里面的说明,只要设置了最大的Z值,应该会将所有的汇都填充掉,但事实上没有。唯一的解释了填充汇的同时又产生了新的汇。
测试3
使用笔者按流程做好的Z值计算工具,针对填洼限制值为33的DEM数据重新计算Z值限制,如下图
计算出来新的Z值为60,如下图
再利用61作为Z值限制重新对填洼为Z值为33的DEM数据进行再次填洼,并重新查找一次汇。这一次的结果没有再发现任何的汇,跟测试1同样的结果,如下图
测试4
流程与测试1一样,直接使用原始DEM数据填洼,但填洼的Z值限制一步到位设置为61,然后查找汇。结果跟预测的一样,所有的汇都被填充,如下图
测试5
该测试才是最后的测试结论,分别都这个产生无汇结果的填洼过的DEM数据进行对比,用减法工具执行。如果这些DEM是一样的话,减法的结果将会没有任何的数据。论证的目的在于设置Z值和不设置Z值到底有多大的差别。
第一组对比为测试1的无Z值填洼的DEM数据VS测试3的第二次填洼的DEM数据(Z值为61的那次),减法结果如下
减法的结果为0,也就是说这两个填洼后的栅格无凹陷点DEM是一模一样的。如下图
第二组对比为测试1的无Z值填洼的DEM数据VS测试4的直接使用Z值为61填洼的DEM数据。同样做减法,结果如下
结果同样是0,不用做第三组测试也知道测试2与测试3的填洼后的DEM也是一样的。
测试6
对原始DEM数据进行Z值为300的填洼,然后查找汇,结果显然也是完全被填充,不再存在汇;然后用减法对填洼后的数据与无Z值填洼的DEM数据做减法对比,结果为0。也就是说只要设定一个无限大的Z值,所有的汇都会被填充,然后构成无凹陷DEM数据前一步的填洼数据,与无Z值填洼结果是一样的。
###总结
通过6个测试和3次对比,基本上可以确定,精确计算Z值多次填洼与不设置Z值填洼效果是等价的。基本上可以确定如果要生成无凹陷DEM,直接使用填洼工具做一次无Z值填洼即可。ArcGIS的创建无凹陷DEM方法是从8.x时代就传下来的,笔者的测试版本是10.5,估计填洼工具早已经实现了自动化计算Z值以及自动迭代操作,只是工具帮助上没有作更新的说明罢了。因此,在使用ArcGIS水文分析之前可以大胆的对原始DEM数据做无Z值限制的填洼操作。
————————————————
版权声明:本文为CSDN博主「李远祥」的原创文章,遵循CC 4.0 BY-SA版权协议,转载请附上原文出处链接及本声明。
原文链接:ArcGIS水文分析实战教程(4)地形预处理_李远祥的博客-CSDN博客_arcgis憛急智
DEM及数字地理底图制作
(一)1:5万调查区的DEM
调查区的DEM是由17幅1:5万图幅的分幅DEM数据拼接而成的。将该17幅地形图进行扫描,在ENVI图像处理软件中进行校正、配准和拼接,形成整幅1:5万调查区地形图,而后进行地形线矢量化,再结合日本卫星ASTER立体像对生成的15m栅格的DEM及国家地理信息中心提供的境内部分地区的DEM共三部分数据,在MAPGIS软件平台生成1:5万调查区的DEM。
(二)1:5万调查区的数字地理底图
首先,在矢量化地形等高线时,也将河流、道路、山峰、高程点、居民地等要素矢量化;将已完成的1:5万调查区DEM转换成Surfer格式的网格数据,再根据需要在MAP?GIS中绘制出高程间隔为100m、50m或20m的高程等值线图;最终编辑形成调查区数字地理底图。本图的投影方式为高斯投影,中央经线为东经81°,采用以克拉索夫斯基椭球为基准的北京54坐标系。
(三)1:1万调查区的DEM
1.技术难点
高精度DEM是1:1万灾害与地质环境定量遥感调查与监测工作的基础,在山岭起伏地区制作高精度DEM是当今国内外的技术难点。其主要技术难点有两方面:一是当今只有很少的建立高精度立体模型的卫星数据;二是缺少在高差起伏较大地区生成高精度DEM的技术方法。
2.技术难点攻关及作业过程
(1)寻求高分辨率卫星立体像对
本项目要求建立1~5m栅格DEM,目前广泛使用的SPOT-5卫星的2.5m立体像对不能满足精度要求。经过调研,除了SAR以外,目前只有美国OrbView卫星立体像对可能制作这样高精度的DEM。经过一年多的努力,直到2006年11月份才获得该卫星数据。OrbView-3卫星是世界上最早提供高分辨率影像的商业卫星之一。卫星轨道高度470km,回访周期<3天,全色波段的波谱范围为450-900nm,空间分辨率1m。本项目采用了12幅共6个像对的1m分辨率的OrbView卫星影像数据建立立体模型,生成DEM。
(2)软件平台
开始试采用VirtuoZo作业,但普通的VirtuoZo全数字测图系统软件不支持OrbView卫星影像,经向VirtuoZo供应商要求提供技术援助后,获得了为西部测图新开发的可以支持OrbView卫星影像的VirtuoZoSeri软件的有限使用权。
该项工作还使用了ERDAS、ENVI和PHOTOSHOP等辅助。
(3)三种作业流程方案及对比
高精度DEM是在调查区1:5万工作DEM和数字地理底图完成后进行的。由于制作大起伏山区的高精度DEM是一项探索性工作,所以我们设计了三套方案的工作流程:①从1:5万地形图上选择平面控制点及从1:5万DEM上确定的高程来校正用RSAT模块定向OrbView卫星立体像对形成的DEM;②通过自由网平差来校正用RSAT模块定向Orb?View卫星立体像对建立的DEM,而后再用地形图上的控制点校正;③无控制点,根据卫星轨道参数,通过自由网平差用RSAT模块定向OrbView卫星立体像对建立DEM,如图1?2所示。
图1?2 建立1:1万DEM工作流程的三种方案
在执行“方案一”的作业过程中,定向中误差非常大,最大定向中误差达17.852m。究其原因是控制点本身误差太大,所以在参与定向时也不能控制住。分析影响控制点精度的主要因素有以下几点:①栅格地形图误差,控制点是在纠正后的1:5万栅格地图上读取的,1:5万栅格图的一个像素尺度为约4m,现要制作1m栅格的DEM,所以其精度相对较低;尽管已经对1:5万地图采取逐格网纠正,也会有较大误差;作为地理控制的地图资料与影像资料的时间间隔超过20年,在该强风化地区,地形地貌会有一定变化,不容易选择同名点。②地形变化误差,调查区属于高山峡谷地形,难以找到比较固定的参考地形,基本上都是通过河流来选择控制点,由于水面季节性变动及强烈冲刷等原因,20年来河流的边线或形状发生了较大变化。③两种坐标系统转换误差及DEM误差,虽然每幅都有自己的转换参数,但仍存在不同椭球系统之间的转换差,从国家地理信息中心提供的DEM读取控制点高程,该DEM格网间隔为25m,相对1:1万工作,误差太大。
后执行方案二,先用立体像对,通过数字摄影测量的自由网平差方法,制作一套正射影像(DOM),利用影像本身的经纬度,通过坐标转换和移位,使地形图和生成的DOM的位置相关,并参照该地区的ASTER影像图寻找栅格图和影像的同名点,读取所选控制点的54平面坐标。再将控制点的54坐标转换为80坐标,把80坐标的控制点与已制作完成的1:5万80坐标的DEM进行套合,读取控制点的高程数据。这样虽然确定了控制点,但由于上述地形图与影像资料时间差太大和特殊地形,获取的成果精度仍不合格。对控制点分析结果表明,控制点参与定向后,残差比没有控制点参与的要大得多,引入控制点作业会加大作业区的内部误差。
因此,最终采用方案3-主要使用卫星的轨道参数来控制。
(4)提高DEM精度的方法
本项目采取以下解决办法:①在纠正地形图时采取逐点(每个格网点都参与)二次多项式纠正法,尽量减少纠正误差;②该高山峡谷地区在地形图和影像图上选取控制点,难度均很大,后来以该地区的ASTER彩色影像辅助参照选点,并在控制点套合DEM读取控制点高程信息时,尽量将所有控制点对应的DEM处放到最大,以减少人为选择平面控制点误差;③创建完立体模型后在显示立体工具栏下可以看见生成的立体影像,但由于地形高差太大,在测图模块下不能显示立体;此外,创建的立体模型不能编辑DEM,但可以自动匹配DEM,也可以生成正射影像。对这些问题,均与协作方联合攻关,最后所有软、硬件问题都一一得到解决。
(5)图像处理
ETM、SPOT、ASTER、CBERS-2各类卫星数据的图像处理,包括多光谱合成、数据融合、镶嵌、几何校正与图像配准工作,主要在ENVI、PCI和PHOTOSHOP平台上进行。
在获取高精度DEM以前,地面分辨率≤1m的高分辨率图像的校正是基于1:5万DEM的,所以其绝对精度只有1:5万。1:1万高精度正射影像及各时相影像之间的精确配准是滑坡及地质环境定量解译与监测的基础与保证。在建立合格的1:1万DEM后,将已获取的2004-2007年度QUICKBIRD、ALOS共8个时相的多光谱数据重新进行3、4、2波段合成及与全色波段融合,并全部与OrbView DOM(1个时相)进行图像对图像校正、配准,并统一重采样成1m分辨率的图像,至此完成调查区1:1万9个时相的多光谱正射图像制作。
(6)人机交互解译及验证
人机交互遥感解译,就是基于滑坡地学原理,在处理合格的解译基础上,采用人机交互方法进行解译,获取滑坡及地质环境基本信息。解译主要在MAPGIS、ENVI和PHOTO?SHOP平台上进行。
1:5万灾害与地质环境解译以5m分辨率的SPOT-5多光谱正射影像为基础,同时参照ASTER、ETM及ALOS影像。本区的地质工作程度较低,区内唯一详细的资料是1:25万扎达幅和斯诺乌山幅区域地质图。但据访问,由于地形复杂及气候恶劣等原因,填图工作未能到达帕里河流域。本项目遥感解译,首先参照该图及文字说明,结合影像特征建立解译标志,然后据解译标志逐片解译。初步解译完成后曾去西藏现场验证,虽已是6月,但由扎达通往帕里河调查区需翻越的多座5000m高程以上的垭口,积雪覆盖太厚,虽雇了当地民工及马匹,还是未能到达帕里河流域。由于喜马拉雅山脉东西两端气候虽有较大差别,但地形是基本对称相似的,所以我们便辗转到了东端的南迦巴瓦峰山脉,考察了那里的冰川与泥石流地形与环境。此外又通过访问当地曾去过帕里河的水利及地质环境监测站人员了解实地情况,收集了帕里河的野外照片,并通过附近卫星影像对比解译来验证调查区的灾害与地质环境情况。野外验证返回后,再次对全区灾害与地质环境进一步解译分析。
(7)GIS和空间分析
将以上解译获取的基本信息在GIS系统中进行空间分析及计算,包括重点调查区的灾害类型、性质及环境分析,灾害体位置、形态及规模估算;1:5万调查区重力侵蚀类型与位置确定、规模计算、危险性评价及与环境关系分析。该项工作主要在MAPGIS、ARC?VIEW和ENVI平台上进行。
(8)成果精度
1)1:1万遥感调查。本项目调查区总体地形困难程度应属最高的三级高山地,但对于局部滑坡而言也有相对较平缓的地形,对多时相滑坡监测,要求有更严格的几何校正及各时相图像的配准,所以要求中误差达到1m以内。需要说明的是,这只是重点区范围内部的相对精度,如表1?2所示。
表1-2 本项目重点区内部1:1万DEM精度
另需说明的是,项目工作的前一阶段,由于未能获得建立用于1:1万调查的高精度DEM的数据源,所以只能先建立1:5万DEM,相应的重点工作区虽然购买了0.6m分辨率的卫星数据,但校正及配准精度还是1:5万的,解译基础(正射影像、DEM和数字地形图)也只能是1:5万精度的。直至2006年12月才重新建立了重点区的高精度DEM及解译基础。
2)1:5万遥感调查。本项目采用的1:5万DEM由前述三部分组成,境内部分满足国家测绘标准,境外部分精度难以统计。
1:5万灾害与地质环境解译以5m分辨率的SPOT-5多光谱正射影像为基础,同时参照ASTER、ETM及ALOS影像。就地面分辨率而言,足以满足1:5万调查的要求。
在图像处理过程中,主要用满足国家测绘标准的境内DEM作校正及与地理坐标配准,调查区的SPOT图像各景季节不同,PAN数据与多光谱时相也不同,加之在高山峡谷地区,故校正及融合难度都很大。经多种方法比较,最终采用了有限元计算处理,最终融合数据校正误差不超过10个像元。ASTER、ETM及ALOS则与已融合校正的SPOT图像采用图像对图像校正,误差控制在2个像元内。