《p13第十三章-数字高程模型(简版)课件.ppt》由会员分享，可在线阅读，更多相关《p13第十三章-数字高程模型(简版)课件.ppt（87页珍藏版）》请在三一办公上搜索。
1、第十三章、数字高程模型DEM,表面分析,DEM的表示和生成方法DEM基本地形因子的计算方法DEM复杂地形分析,第一节、数字高程模型DEM的表示和生成方法,概念区分,数字地面模型（Digital Terrain Model,DTM）是地表形态等多种信息的一个数字表示，其中的属性包括地形（x,y,z）、地貌、地物、自然资源、环境、社会经济等等信息的定量或定性描述；因而DTM是GIS中一种非常重要的产品；DTM也称为数字地形模型；数字高程模型（Digital Elevation Model,DEM）是表示区域上地形三维向量的有限序列Vi=(Xi,Yi,Zi)Xi,YiD是平面坐标，Zi是（Xi,Yi
2、）对应的高程，对于规则格网Vi=Zi；DHM(Digital Height Model)是一个与DEM等价的概念；4D产品：数字高程模型(DEM)、数字正射影像图(Digital Ortho image Map，DOM)、数字线划图(Digital Line Graphic,DLG)和数字栅格地图(Digital Raster Graphic,DRG)等前3D为国家空间数据基础设施（NSDI）的框架数据,一、DEM的表示方法,数字地形的表达可分为：矢量点等高线三角网栅格（连续栅格）,二、DEMGRID的数据来源和生成方法,DEM的来源DEM的生成方法人工网格法立体像对分析等值线插值三角网转换法
3、高程点插值方法,DEM的数据来源,直接定位测量：GPS接收仪、普通测量设备（如经纬仪）；用于小范围内各种大比例尺高精度的地形建模，这种数据获取方法的工作量很大，效率不高，费用高昂；地形图：数字化或由等高线自动生成；对于经济发达地区，由于土地开发利用使得地形地貌变化剧烈而且迅速，既有地图往往也不宜作为DEM的数据源；但对于其他经济落后地区如山区，因地形变化小，既有地图无疑是物美价廉的数据源。航空影像：航空影像测量方法可直接得到DEM；航空摄影测量一直是地形图测绘和更新最有效也是最主要的手段，其获取的影像是高精度大范围DEM生产最有价值的数据源。利用该数据源，可以快速获取或更新大面积的DEM数据，
4、从而满足对数据现势性的要求。,DEM的数据来源,遥感卫星：卫星扫描系统（如SPOT卫星上的立体扫描仪）获取的图像也能提供DEM；SAR：雷达干涉波可直接得到DEM；机载激光扫描仪也可直接获取DEM；,1、人工网格法,将地形图蒙上格网，逐格读取中心或交叉点的高程值、构成数字高程模型。,2、摄影测量法,利用遥感立体像对，根据视差模型，通过选配左右影像的同名点，可建立数字高程模型。,3、等值线插值,根据各局部等值线上的高程点，通过插值公式计算各点的高程，得到DEM。等值线插值法是比较常用的方法，输入等值线后，可在矢量格式的等值线数据基础上进行，插值效果较好。,等高线插值算法生成DEM的方法,采用移动
5、拟合加权平均插值方法。设A点为待内插的点，从A点按45的方位间隔引出八条搜索射线，八条射线与A点相邻等高线的交点为C1，C2Ci，其高程分别为Z1，Z2Zi，它们到P点的距离设为d1，d2di，则P点的插值高程Zp为,4、高程点插值方法,以不规则点图元组织的Z变量的数据，并不适合于图形显示，也不适于进行分析。高级曲面分析要求将Z值转换成一个规则间距空间格网，或者转换成不规则三角形网。栅格法可用来将不规则的空间数据转化为规则格网的空间模式。栅格法以一系列不规则分布的点图元为基础，它涉及某个有效范围内任意点Z值的内插，这个生成的点通常在在格网垂直线与水平线的交叉处。反距离加权法样条函数法趋势面法克
6、利金法,5、三角网转换法,对有限个离散点，每三个邻近点联接成三角形，构成三角网。每个三角形代表一个局部平面，再根据每个平面的拟合方程，可计算各格网点高程，生成DEM。,（1）、TIN的生成方法,首先取其中任一点P，在其余各点中寻找与此点距离最近的点P1，连接P和P1构成三角形的第一边，然后在其余所有点中寻找与这条边最近的点P2，找到后即构成第一个三角形，再以这个三角形新生成的两边为底边分别寻找距它们最近的点构成第二个、第三个三角形，依此类推，直到把所有的点全部连入三角网中。,（2）三角网的插值：双线性方法,求落在各个三角形内的网络点高程值（包括落在三角形边上的点）待求点落在三角形ABC内，先用
7、线性插值的方法，求D，E两点的值。设A，B，C，D，E，P处的值分别为VA、VB、VC、VD、VE，其中VA、VB、VC为已知，在DEM中实质上为高程值，则D、E两点处的插值为,则P点的插值为：,ArcGIS的命令,TINLATTICE LINEAR|QUINTIC z_factor FLOAT|INT-待转换的三角网；-要生成的DEM；LINEAR 插值方法选择线性插值；QUINTIC-插值方法采用双变量五次插值方程z_factor Z坐标的调正因子；FLOAT|INT GRID的输出数据类型。,生成方法小结,DEM数据采集的几点结论摄影测量是DEM的重要数据源，是进行数据库更新的最有效方法
8、之一；现有地形图是DEM的另一重要数据源，从等高线生产DEM的方法已经完全成熟，被广泛地用于生产；使用GPS、激光扫描、干涉雷达等新型技术进行DEM数据采集是很有发展前景的方法；利用基于不规则三角网TIN的方法进行数据建模和栅格转换，是快速可靠地生产高精度DEM的切实可行方案；无论哪种方法，最大限度地采集重要的地形特征点是保证DEM质量和提高作业效率的基本前提；,第二节、DEM基本地形因子计算,地形分析分类,基本地形因子计算,复杂地形分析,地形分析,坡度计算,。,地形水文分析,地形特征提取,可视性分析,坡向计算,。,粗糙度计算,可从DEM数据中通过地形分析方法得到的一些主要地形属性（据I.D.
9、MOORE 等）,1、高程分级,高程影响地表物质和能量的分布；等间距或不等间距划分为若干高程等级，如用来区分丘陵、低山、中山、高山等，划分地貌类型。,ArcGIS命令,LATTICEPOLY lookup_table z_factor转换栅格中的坡度、坡向、高程范围和边界成多边形SLOPE 根据查找表中的坡度等级生成多边形；ASPECT 根据查找表中的坡向等级生成多边形；RANGE 根据查找表中的高程等级生成多边形；NODATA 根据栅格单元是否有值或空值来生成多边形；BOX 生成最小矩形区域包含输入格网中的所有栅格EXTENT 按栅格的extent生成单个矩形框；,查找表,lookup_ta
10、ble 包含分类范围和相关输出多边形代码lookup item-PERCENT_SLOPE,DEGREE_SLOPE,ASPECT or RANGE作为lookup_table的索引字段.类型可以为I,B,F,N例如PERCENT_SLOPESLOPE-CODEData ranges31 0.0 slope=352 3 slope=5103 5 slope=10254 10 slope=25405 25 slope=40606 40 slope=60807 60 slope=80,2、平均高程,式中n的计算单元内栅格个数；h（Pk）为第k点的高程。ArcGIS中使用Describe查找或者使用
11、ZonalMean()获取,3、极值高程和高差,4、区域相对高程,设参考高程为hm，则各栅格点上相对高程为：k=1，NArcGIS下使用newdemdem-heght命令，其中Height可以为常量,5、局部相对高差,格点面元是在格网DEM 的水平投影面上，以四个相邻格点(i,j),(i,j+1),(i+1,j+1),(i+1,j)为顶点的面积范围；格点面元的相对高差指在格点面元的四个格点中，最高点与最低点之差，按下式求算：,h=MAX(h00,h01,h02,h03)MIN(h00,h01,h02,h03),式中hij=(i,j=0,1)为四个格点的高程;ArcGIS-Grid的实现方法？F
12、ocalmax-Focalmin,6、剖面（Profile）,用于检查和量测沿某条线的高度变化SURFACEPROFILE profile_info_table sample_distance,ArcView下的图示,剖面积计算,根据某穿越地形表面的线路，求该线路垂直剖面积；具体计算时先求得该线路与DEM 格网边的所有交点Pi(Xi,Yi,Zi)，再按下面公式计算得到:,n 为交点数;Di,i+1 为Pi 与Pi+1 间距离。,计算曲面沿线长度,SURFACELENGTH z_factor sample_distance surface_length_item计算沿弧段的曲面长度,7、坡度和坡
13、向,坡度和坡向是表示地表面在地面某一点处的倾斜程度和倾斜方位的一个量；它是一个矢量，既有大小又有方向。某个点的坡度是其数值等于地表曲面在该点的切平面与水平面夹角，某个点的坡向是指该切平面上沿最大倾斜方向矢量在水平面上的投影的方位角；,坡度的表示,输出的坡度网格单位为度或百分数。,注意:坡度表示为0-90，0为代表水平面的情况，当坡度为90百分数为无贫大,坡度坡向的计算,拟合曲面法坡度的计算方法有很多，经证明，拟合曲面法是解求坡度的最佳方法；,一般采用二次曲面，使用一个3*3 窗口，每个窗口单元的中心有一高程点，,坡度坡向的计算公式,坡度坡向的计算公式,例子,(dz/dx)=(a+2d+g)-(
14、c+2f+i)/(8*x_mesh_spacing)(dz/dy)=(a+2b+c)-(g+2h+i)/(8*y_mesh_spacing)rise_run=SQRT(SQR(dz/dx)+SQR(dz/dy)degree_slope=ATAN(rise_run)*57.29578,ESRI,ArcGIS中的SLOPE函数,SLOPE(,DEGREE|PERCENTRISE),ArcGIS中的Aspect函数,ESRI,ASPECT(),注意坡向的规定,8、Curvature地表曲率,profile curvature：坡度方向的曲率，区别凸坡和凹坡，剖面曲率影响水流的加速和减速，从而影响侵蚀
15、和沉积过程。在凸坡区域，水流加速，水流可聚集更多的能量，输送量也增加，发生侵蚀过程，相反，在凹坡区域，水流减速，水流失去能量，发生沉积过程）Plan curvature：表示等高线方向的曲率，影响水流的聚集和分散（Converging/diverging flow）,以及土壤水分的含量（soil water content），用于根据坡形的水平凹凸变化情况，确定沟谷线、山脊和鞍部的位置，划分流域范围；Tangential curvature切线曲率：水平曲率和坡度的乘积（Plan curvature multiplied by slope），Provides alternative measu
16、re of local flow convergence and divergence，影响局部水流的聚集和分散。,剖面曲率的算法原理,待算格点面元的四个格点中，最高点与其对角点的连线称为格点主轴，主轴两端点高程的平均值与格点面元平均高程的比，称为格点面元凹凸系数（CD）：,式中hmax 为最高格点高程，hmax 为最高格点的对角格点高程，h 为格点面元高程的平均值，当CD 为正时，格点面元的实际表面为凸形坡，为负时为凹形坡；,剖面曲率,表示坡面是凸坡或凹坡，可导致水流的加速或减速，分别导致侵蚀和沉积。,ArcGIS计算方法,CURVATURE(,out_profile_curve,out_p
17、lan_curve,out_slope,out_aspect)正的曲率值表示凸坡，负值表示凹坡，0表示平地,9、表面积计算,即求算格网表面的面积对于格网分解为两个三角形，然后求三角形的表面积，整个DEM的表面积为所有格网单元面积之和；对于TIN 则直接求三角形面积；求三角形面积使用每点的(x,y,z)值;,表面积计算公式,具体的计算公式为:,式中，Di 表示第i 对三角形两顶点间的表面距离，S 表示三角形的表面积，P 表面三角形的周长的一半,10、投影面积计算,投影面积指的是任意多边形在水平面上投影的面积，可直接采用海伦公式进行计算，一种简单的方法是根据梯形法则，如一个多边形由顺序排列的N 个
18、点(xi,yi,i=1,N)组成，并且第N 点与第1 点相同，则水平投影面积计算公式如下，如果多边形顶点按顺时针方向排列，则计算的面积值为负，否则为正。,11、地表粗糙度,反映某一面积单元内地势伏变化的复杂程度格点面元的粗糙度（roughness）指格点面元所对应的DEM 上表面积与其水平投影面积之比,记为CZ:CZ=S 表面积/S 投影面积当CZ=1 时,粗糙度最小，格点面元的实际表面即为水平面;,12、体积计算,DEM 的体积可由四棱柱和三棱柱的体积进行累加得到，四棱柱上表面可用抛物双曲面拟合，三棱柱上表面可用斜平面拟合，下表面均为水平面或参考平面，计算公式分别为:,其中S3 与S4 分别
19、是三棱柱与四棱柱的底面积,体积计算的ARCGIS命令,VOLUMEbase-valueout-info-filez-factor输出结果：NAMEWIDTHOUTPUTTYPEN.DECTIN3232C-ZMIN412F3ZMAX412F3DATUM412F3AREA818F5 planimetric平面面积VOLUME818F5UNITS3232C-ZUNITS3232C-ZFACTOR412F3,14、挖方填方计算,通过对两个输入曲面进行对比，计算出挖方区域和填方区域CUTFILL z-factor,计算流程,15、山体阴影,HILLSHADE azimuth altitude ALL|S
20、HADE|SHADOW z_factorAzimuth：太阳方位角Altitude：太阳高度角 SHADE 仅考虑 local illumination angles;不考虑遮盖的影响。输出范围在0 and 255,0 表示黑，255表示亮.SHADOW 输出为0表示遮盖区，1表示非遮盖区,Shade选项,Shadow选项,地表辐照度,计算辐照度需考虑日照条件（太阳赤纬、高度角、时角及大气状况）与坡面几何条件的相互关系。辐照度由下式决定：式中，大气透过率，与太阳高度和大气状况有关；Sc为太阳常数；Sa为太阳高度角，可由球面三角公式求出；t是时角；a、b为坡面方程系数；为坡度。,16、生成等高线
21、,TINCONTOURLATTICECONTOUR,二、地表水再分配分析,利用数字地形模型作为输入，通过确定栅格任意点上坡区域的贡献和下坡区的水流方向，勾画出水系并且对其相关特征定量化。利用DEM生成的流域和水系，是大多数水文分析的数据源，可用于决定洪水淹没范围分析，进行流域渗透水文的模拟研究。,1、地表水再分配的基本过程,水系和流域的概念水系连接的类型径流过程洼地和峰的问题,1）、水系和流域的概念,水系（drainage system）为水在地表流动并且流向出口的网络。流域（drainage basin,watershed,basin,catchment,）是地表水和其他物质流向统一出口点(
22、outlet,pour point)的区域。出口点为流域边界的最低点。两个流域的边界被定义为分水岭（drainage）,2）、水系连接的类型,水系可以理解为树，出口点为树根，树干为河沟，河沟的交点被称为节点或交点（node or junction.）。连接两个内流的节点，或连接节点和出口点的河沟叫做内流连接（interior links.）。外连接（Exterior links）定义为树最外面的分支。,3）、径流过程,水流的方向取决于每个点的坡向；水流的能量取决于地表坡度，陡的坡度产生大的能量，可以输送更多和更大的固体物质，并具有较大的侵蚀潜力；水平曲率用于区别斜坡是山脊或山沟，山脊和山沟分别
23、导致水流的聚集和分散；,4）、洼地（sinks）和 峰（peaks）,DEMS中的错误通常划分为洼地（sinks）和 峰（peaks）。洼地（sinks）是比较高的邻域所包围的区域，洼地可以是自然形成的。峰（peaks）是被比较低的邻域所包围的地区，通常是自然地貌。洼地的数量在低分辨率时较多。另外当珊格的数据格式采用整形时，在低高差的地区洼地数量也会增多。当计算流向时，洼地会产生不合适的结果，在提取流向时剔除。DEM可能包含明显的水平条带，这来源于生成dem时的系统取样误差，这在平坦地表的整形珊格中非常明显。,2、流域河网、流域边界提取算法,基于栅格的地形分析技术在水文学中的应用，一般采用了O
24、Callaghan和Mark的坡面流模拟方法，Jense和Domingue（1988）、Martz和De Jong（1988）、Garbrecht（1997）在此基础上做了改进，这种算法一般称之为D8（Deterministic eight-neighbours）算法；,D8算法示例,D8 模型,高程格网示例,格网流向模型,水流聚集模型,算法流程图,1）洼地的填充,洼地的出现可导致流向的计算错误。在某些情况下，可能有合理的洼地存在，这就需要对当地的地形情况有较多了解，以辨别真的洼地和错误。洼地和流向不确定的情况：（1）周围单元格的值高于处理的单元格的值（2）两个单元格间的水流互相流动（3）单元
25、格在多个方向上有同样的高差。要精确的表示流向，应该对洼地进行填充；洼地可以用Sink函数定位，并进行填充。,填充洼地的过程计算洼地深度,1、利用SINK函数识别洼地Grid:sinks=sink(flowdir)2、利用WATERSHED函数生成每个洼地代表的流域区域Grid:sink_areas=watershed(flowdir,sinks)3、生成每个洼地流域范围内的最低高程值栅格Grid:sink_min=zonalmin(sink_areas,elevation)4、将洼地填充到沿洼地边界的最低高程值，使水流能流出洼地。表示为Sink_maxGrid:sink_max=zonalfi
26、ll(sink_areas,elevation)5、最低高程值网格减去最低高程值网格得到洼地深度栅格Grid:sink_depth=sink_max-sink_min洼地深度的计算便于下述fill指令确定z_limit的取值,填充洼地的过程填充,FILL SINK|PEAK z_limit out_dir_grid：表示连续曲面的网格：剔除洼地或山顶的输出网格SINK|PEAK：表示洼地将被填充或山顶将被削去z_limit 洼地和其出口点的最大高差或山顶与其相邻最高单元格之间的高差。洼地深度大于这个值，则该洼地不被填充；如果山顶与其相邻最高单元格之间的高差大于这个值，则山顶不被移去。默认情况下
27、时填充所有洼地或移去所有山顶，而不考虑深度。out_dir_grid 可选输出流向栅格,2）Flow Direction,按流向选最大坡度的原理，计算每个单元格流向其邻域的方向（流出方向）；如果到相邻栅格的坡度相同，则放大邻域以选择最陡的方向,流向的规定,A.判断指标：距离权落差,距离权落差按下式计算，对角方向距离D为sqrt(2)格网尺寸，其它方向D为格网尺寸,B.确定水流方向步骤,1)对所有DEM 边缘的格网，赋以指向边缘的方向值；(可增加判断流向内部格网）2)计算中心格网对8 个邻域网的距离权落差值；3)确定具有最大落差值的格网；（分以下几种情况分别处理）如果最大落差0，则赋负值；如果最
28、大落差0，赋最大值方向编码；如果最大落差0，有一个以上最大值，按逻辑判断作为流向编码；如果最大落差0，则0值对应的方向值相加，作为其方向代码；4)对没有赋以负值或水流方向值为1,2，4，128，的格网，检查最大落差邻域格网，如果此格网正确编码，且方向没有指向中心格网，则以此格网的方向值作为中心格网的方向值；5)重复第4 步，直至没有任何格网能被赋以方向值，对方向值不为1,2，4，128 的格网赋以负值。,ArcGIS流向命令,Grid:flow_dir=flowdirection(elevation),3）Flow Accumulation,累积流表示流入该单元格的栅格数目；高的累积流可用于生
29、成水系网；低的累积流表示局部的地形高区，可用于识别山脊；,ArcGIS中累积流的计算,FLOWACCUMULATION函数计算流入每一个上坡方向的栅格累积值作为这个栅格的累积流。也可以考虑权重因子，比如降雨量；Grid:flow_acc=flowaccumulation(flow_dir,rainfall),4）流域河网的提取,确定上游集水区面积阈值；标注集水量累计值大于阈值的格网；将得到的图形进行矢量化处理；通过阈值调整控制生成河网密度；,水系网提取的ArcGIS,水系网可从FLOWACCUMULATION函数中得到的结果中勾画出来。通过使用地图代数表达式对累积流网格设置阈值，可以勾画出水系
30、网。例如,要生成一个在无值背景上值1代表的水系网，可以使用：Grid:stream_net=con(flowacc 100,1)Grid:stream_net=setnull(flowacc 100,1)上述单元格取值大于100的所有单元格被赋值为1，其它单元格被赋值为无值NODATA。,5）、水系网级别的划分,水系定级是一种给水系网中的连接赋级别代码。目的是按水系的支流数识别和分类河谷类别。水系的某些特征可以通过其级别推断出来。例如，一级支流仅受坡面漫流的控制，它没有其它支流的水流贡献,它们极易受非点源污染的影响.Usage:STREAMORDER(,STRAHLER|SHREVE),STR
31、AHLER/SHREVE,STRAHLER外连接总被赋值为1。但当同样级别的水系交叉时，水系级别增加。所以，两个第二级的水系交叉时，会生成第三级的水系。但当两个不同级别的水系交叉时，则不会产生水系别的增加，例如，第一和第二级的水系交叉不会产生第三级的水系，而保持较高级别的水系等级。Shreve法计算水系网中的所有连接数，外连接被赋值为1，但对内连接其级别数是增加的。例如，两个一级支流交叉时会产生二级支流，而一级支流和二级支流交叉时则生成三级支流。,6)、水系矢量化,STREAMLINE(,out_item,weed),7）ArcGIS流域的划分,watersheds=watershed(flo
32、w_dir,pour_points)BASIN(),ArcGIS中的流域命令,将流向网格作为输入，可以利用WATERSHED或 BASIN 函数生成流域。Grid:watersheds=watershed(flow_dir,pour_points)Grid:wsheds=watershed(flowdirection(elev),selectpoint(elev,*)对任意单元格可以勾画流域，在WATERSHED函数中可以利用SELECTPOINT函数标记流域出口点(pour points)。在交互选择流域出口点时，显示累积流网格，有助于选择主要的水系分段，否则，可能选出接近主水系的网格，生成
33、的流域要比期望的小。,小结：ArcGIS中的处理流程,3、其它基于DEM的高级指数,地表湿度指数：假定在恒定饱和条件下，径流是上坡集水面积、土壤渗透性和坡度的函数,水流侵蚀能量指数：假定流量是流域面积的函数，则水流侵蚀能表示为。侵蚀发生在剖面凸坡和切线凸坡（水流加速和聚集），沉积发生在剖面凹坡处,4、流域生态水文的次生参数（1）,流域生态水文的次生参数（2）,二、视见分析,ARC/INFO中提供3种方法做视见分析：ARC中的VISIBILITYARCPLOT中的SURFACEVIEWSHEDARCPLOT中的SURFACEDRAPE,1、VISIBILITY,ARC：VISIBILITY PO
34、LY|GRID FREQUENCY|OBSERVERS-输入栅格，做视见分析-输入的图层，为点或线，定义观测点.POINT-中的点定义观测位置LINE-线中的节点和拐点定义观测位置 POLY 输出结果为多边形 GRID 输出结果为栅格，每个栅格记录了被观测的次数(FREQUENCY option),或者对应的观测点编码(OBSERVERS option).,视见分析的结果,2、SURFACEVIEWSHED,SURFACEVIEWSHED 在定义的视见分析环境下，该命令高亮显示可视范围 SURFACEVIEWSHED 使用以下命令视见分析环境:SURFACEOBSERVER 设置观测点位置SU
35、RFACETARGET 设置目标点位置SURFACEEXTENT 设置分析的范围SURFACEVIEWFIELD 设置观测场（方位）（可选AUTO）SURFACERANGE 设置从观测点到曲面的距离限制SURFACERESOLUTION 设置临时结果栅格的分辨率,3、SURFACEDRAPE,SURFACEDRAPE MESH FISHNET|DIAGONAL|ALONGX|ALONGY distance lookup_table|linesymbol_grid item SURFACEDRAPE LOCATOR distance*SURFACEDRAPE XYZ z_valueSURFACEDRAPE 通过将对象叠置在曲面上，生成一个曲面视见范围，对象可以为线网孔、点或图形文件 MESH.叠置网孔线.LOCATOR.用于交互定位网孔线 XYZ.叠置点GRAPHICSFILE.叠置图形distance 网孔线之间的距离 在执行SURFACEDRAPE 前，必须使用SURFACE 命令设置当前处理的曲面，并且建立视见分析环境，然后再使用SURFACEDRAPE或SURFACEVIEWSHED命令,
