地理信息系统下的空间分析-第七章_三维数据的空间分析方法

该方法基于地图上地形坡度越大等高线越密、坡度越 小等高线越稀这一地形地貌表示的基本逻辑，将所研究的 区域划分为m×n个矩形子区域(格网)，
然后计算各子区域内等高线的总长度，
再根据回归分析方法统计计算出单位面积内等高线长 度值与坡度值之间的回归模型， 然后将等高线的长度值转换成坡度值。 该算法的最大优点是可操作性强，且不受数据量的限 制，能够处理海量数据。
但这些算法有一个最大的缺点：当数据量较大时， 需要的计算机内存容量比较大，必要时必须将整个DEM 进行分块处理。
对于矢量数据来说，如果采用这类算法需要先将矢 量数据转换成格网数据，然后在此基础上进行曲面拟 合和坡度与坡向的计算。 在这一过程中，存在着两次信息损失：
(1)数据类型转换过程中的信息损失：将包含多个数 据点的矢量数据转换成规则化的格网点上的数据(栅格 数据)时的信息的丢失；
第七章 三维空间分析
主要内容：
1、表面积计算
2、体积计算
3、坡度、坡向计算
4、剖面分析
5、可视性分析
6、水文分析
7.1 三角形格网上的表面积计算
基于三角形格网的曲面插值一般使用一次多项式模型 (Z=a0+a1X+a2Y），计算三角格网上的曲面片的面积时， 首先将其转换成平面片，然后通过计算平面片的面积来计 算曲面片的面积。
该点的法线方程为：
，
f x1 ( x0 , y0 )(x x0 ) f y1 ( x0 , y0 )( y y0 ) ( z z0 )
其方向数为 f x ( x0 , y0 ) ，f y ( x0 , y0 ) 和－1，垂直方向Z的方 向数为0，0，1， 则有： cos a
b (b' 2 (h2 h3 ) 2 )1 / 2
c (c'2 (h3 h1 ) 2 )1/ 2
7.1.2 正方形格网上表面积的计算 1、曲面拟合重积分计算法
正方形格网上的曲面片表面积的计算问题要复杂得多， 正方形格网上最简单的曲面模型为双线性多项式，其拟合 面是一曲面，无法以简单的公式计算其曲面积。 根据数学分析，某定义域A上，空间单值曲面的面积由 以下重积分计算： S (1 f 2 f 2 )1 / 2 dxdy
2、基于正方形格网的体积计算
正方形格网的基本格网的体积计算方法如下：
V S A (h1 h2 h3 h4 ) / 4 (7.11)
其中，SA是基底格网正方形A的面积。如图(b)所示。
7.3 坡度计算
如果打算在山上建造一座房子，那就首先想知道哪里 是比较平坦的地方；如果要设计滑雪娱乐场，并且要把它 建在不同的坡度上。对于这些问题都需要考虑坡度。 坡度是地形描述中常用的参数，是一个具有方向与大 小的矢量。 作为地形的一个特征信息，坡度能间接表示地形的起 伏形态，在农业、交通、规划以及各类工程中有着很大的 用途，如，在农业土地开发中，坡度大于25°的土地一般 被认为是不宜开发的。
该方法的基本思想是： 设置一个小窗口，首先计算小窗口内单根矢量等高线 的坡向(等高线法线的倾角)，然后利用如下公式计算窗口 内的最终的坡向： ( (li i ) / li

i 1
/
i
式中， 为窗口内的最终坡向， l i 为窗口内单根等高线 的长度，l 为窗口内等高线的总长度，窗口内的坡向计 算以单根等高线的长度为权值。
以下分别给出基于三角形格网和正方形格网的体 积计算方法。 其基本思想均是：基底面积(三角形或正方形)乘以 格网点曲面高度的均值，区域总体积是这些基本格网 体积之和。
1、基于三角形格网的体积计算
三角形格网的基本格网的体积计算方法如下：
V S A (h1 h2 h3 ) / 3
其中，SA是基底格网三角形A的面积。如图(a)所示。
基准平面一般是水平面，基准平面的高度不同，计算出 的空间曲面的体积就不同， 当高度上升时，空间曲面的高度可能低于基准平面，此 时出现负的体积。 在对地形数据的处理中， 当体积为正时，工程中称 之为“挖方”，体积为负 时，称之为“填方”， 图中的阴影部分为“填 方”。
体积的计算通常也是近似方法。由于空间曲面的表 示方法不同，近似计算的方法也不一样。
当基于DEM计算坡向时，通常定义坡向为：
过格网单元所拟合的曲面片上某点的切平面的法线的 正方向在平面上的投影与正北方的夹角，即法方向水平投 影向量的方位角，坡向在下图中以β标识。
由数学分析知，设曲面 Z f ( x, y) 在点（x0，y0，z0）的切 平面方程为：Z Ax By C f x ( x0 , y0 ) x f y ( x0 , y0 ) y C 则该点的坡向为： arctg( A / B) 但根据此式计算的β在（（0，2π）中取值，因此还需要根据A，B的取值情况将β 转换成（0，2π）之间的值。
(2)曲面拟合时的信息损失：在用有限个格网点数据 来模拟真实地表起伏时，所拟合的曲面与实际地表间 存在误差。 因此，对于矢量数据来说，能否用已有的真实数据 直接(不丢失信息)比较精确地计算坡度和坡向是值得研 究的。
2、基于矢量数据的坡向计算
由于基于DEM的坡向计算方法的实现过程存在着两次 信息损失，因此，可以考虑直接利用矢量数据计算坡度。
7.3.1 单个点的坡度计算
空间曲面的坡度是点位的函数，除非曲面是一平面， 否则曲面上不同位置的坡度是不相等的， 给定点位的坡度是曲面上该点的法线方向N与垂直方 向Z之间的夹角a，如图所示。
由数学分析知，对于曲面上给定点（x0，y0，z0）的切 平面方程为：
f x ( x0 , y0 )(x x0 ) f y ( x0 , y0 )( y y0 ) ( z z0 ) 0
式中，α为坡度， z / x，z / y 一般采用2阶插分方法计 算。
对如图所示的网格中的（i，j）点，
z i , j 1 z i , j 1 z x 2 x z i 1, j z i 1, j z y 2 y
式中， x , y 为格网结点x，y方向的间隔。
1 f x2 ( x 0 , y 0 ) f y2 ( x 0 , y 0 ) 1
a arccos(f x2 ( x0 , y0 ) f y2 ( x0 , y0 ) 1) 1/ 2
由坡度的概念知0°≦a≦90°，故由上式比较容易确 定坡度值。 对于特殊的应用场合，例如对于Z=a0+ a1x+ a2y(三角形 格网上的曲面拟合)，其三角形格网单元上的曲面为一平 面，平面上的坡度处处相等， 可以直接计算如下：
如图所示，P1 P2 P3构成的三角形上的曲面片，P1’ P2’ P3’ 为使用一次多项式模型拟合得到的平面片，计算曲面片的 面积其实是计算拟合后的平面片的面积，利用海伦公式计 算面积，公式如下：
S [ P( P a)(P b)(P c)]1/ 2
P ( a b c) / 2

A
x
y
一般地说，上式无法直接计算，常用的方法是近似计算。 积分的近似计算方法很多(有关计算方法的著作对此都有详细 全面的讨论)，比较常用的方法是抛物线求积方法，亦称辛卜 生方法（Simpson）， 该方法的基本思想是：先用二次抛物面逼近面积计算函数， 进而将抛物面的表面积计算转换为函数值计算。
2、分解为三角形的计算方法
a arccos(a a 1)
2 1 2 2
1 / 2
7.3.2 基于格网数据单元的坡度计算 单独一个点上的坡度并无多大用处，通常总是计算基 本格网单元上的平均坡度。 平均坡度的计算可以通过计算若干点位上的坡度，然 后取其平均值。 但更常用的方法是在基本格网单元上用最小二乘逼近 的方法拟合一个平面，然后计算其平均坡度。 （1）三角形格网 在三角形格网上，最小二乘逼近的平面与插值方法的 平面是一致的，因为已知数据点与平面方程的待定系数个 数相等。 例如，如果用Z=a0+ a1x+ a2y来拟合三角形格网上的曲 面，坡度可以直接利用如下公式计算：
(1)首先根据DEM内插算法内插出A和B两点的高程值； (2)然后求出AB连线与DEM格网的所有交点，内插出各 交点的坐标和高程，并将交点按离起始点（A点或B点） 的距离进行排序；
最后，选择一定的垂直比例尺和水平比例尺，以距离 起始点的距离为横坐标，以各点的高程值为纵坐标绘制剖 面图。如图所示。
7.6 可视性分析
hl / p
式中，h为等高距， l 为测区等高线总长度，P为测区面积。
以上介绍的等高线计长法求出的是一个区域 内坡度的均值，且前提是量测区域内的等高距 相等。 该方法对于测区较大或等高距不等时，计算 出的坡度将有较大误差。
2、基于统计学的坡度计算方法
对于测区较大或等高距不等时，可以采用等高线计长方 法的一种变通方法，即一种基于统计学的方法。
a arccos( a 1) a
2 1 2 2
1 / 2
（2）正方形格网
对于正方形格网单元，平面拟合采用最小二乘逼近方 法。
在计算坡度时，如果在3*3的DEM网格窗口进行，窗 口在DEM数据矩阵中连续移动，然后完成整幅图的计算 工作。
z 2 z 2 1 / 2 tan( ) [( ) ( ) ] x y
2 2，）中取值，而坡向应为
无论是坡度还是坡向，在一个很小的范围内计算都只 有理论上的意义，但计算的原理是一样的。
以上的介绍都是从 Z f ( x, y) 出发，而不是直接从数据 点值计算，主要就是考虑了这一点。
基于DEM的坡度与坡向计算的算法到目前为止多是 基于格网点阵数据的，这类方法算法易于程序化，实施 起来比较容易。
7.3.3 基于矢量数据的坡度计算方法
坡度计算也可以直接根据矢量等高线来直接计算，具 体有如下几种方法：
1、等高线计长法
坡度计算的另一类方法是基于矢量数据的算法，它直 接根据数据点值来进行计算。 这类算法的原理是基于50年代原苏联著名的地图学家 伏尔科夫提出的等高线计长方法。