基于约束Delaunay三角网的道路密度算法

金 澄, 安晓亚，徐道柱

1.西安测绘研究所，陕西 西安，710054;2. 地理信息工程国家重点实验室，陕西 西安，710054



基于约束Delaunay三角网的道路密度算法

金 澄1,2, 安晓亚1,2，徐道柱1,2

1.西安测绘研究所，陕西 西安，710054;2. 地理信息工程国家重点实验室，陕西 西安，710054

本文针对道路密度求解过程中道路覆盖区域面积难以确定的问题，提出一种基于约束Delaunay三角网的道路覆盖区面积计算方法。采取化整为零的思想，先对道路网分布区域进行约束Delaunay三角网(D-TIN)剖分，并根据道路所关联三角形对覆盖区面积的贡献度，将三角形分为四类，给出每类三角形面积贡献度的计算方法，通过累加所有关联三角形贡献的面积，得到道路覆盖区的面积；然后计算道路的密度。实验结果表明，此算法简化了道路密度计算过程中覆盖区域面积的求解过程，算法可靠和高效。

地图综合；道路密度；Delaunay三角网

1 引 言

线状道路是一种重要的人文地理要素，道路网密度是体现城市道路建设数量、水平以及城市路网布局是否合理和均衡的重要评价指标之一[1]，也是地图综合选取操作的重要综合指标之一[2]。在进行道路网制图综合选取时，一般根据地图比例尺和地图综合区域的特点制定不同的密度指标，然后根据密度选取[3]。由于目前计算机的软硬件技术还无法达到人脑的智能化识别能力，无法做到像人眼一样、一眼看出路网的疏密程度，因此，如何自动计算道路的密度是实施路网自动综合选取的关键所在。本文提出的道路网密度计算方法主要用于道路网制图综合与密度合理性评价、与道路关系紧密的居民地制图综合，并可推广到水系网密度的计算与制图综合。目前的研究中，文献[1]利用GIS软件计算道路网的密度，但没有给出具体的实现方法；胡云岗等提出一种基于道路网眼的路网密度计算方法[4]，先将道路网纵横交错形成的最小闭合区块称为“道路网眼”，以网眼密度代替道路的真实密度实施道路综合选取。这种方法一般只适合存在网络回路的小比例尺地图道路或路网比较发达的城市区域路网密度的计算，对于道路比较稀疏的地区，由于道路没有形成网络回路，所以此方法不适用。文献[5]提出了一种基于道路网架的计算方法，先构建路网的骨架线；然后将骨架线连成封闭的多边形，构成道路的覆盖区域，并计算多边形面积；最后根据道路的长度和覆盖区域的面积计算道路的密度。该方法能够实现每条道路的密度计算，但由于数据的复杂性，很难将骨架线提取完整，导致构建道路覆盖区域多边形失败。本文提出一种基于约束Delaunay三角网的道路密度分析方法，采取化整为零的思想，先对道路网分布区域进行约束Delaunay三角网(D-TIN)剖分，同时根据每条道路所关联的三角形对覆盖区域面积的贡献度，将三角形分为四类，给出每类三角形面积贡献度的计算方法，通过累加所有关联三角形贡献的面积，得到道路覆盖区的面积，然后计算道路的密度。

2 基于D-TIN计算道路密度方法

2.1 道路密度计算公式

道路密度一般是指一定区域内道路长度与该区域面积的比值，可理解为一种线密度，其计算公式为[4]：

D=L/A

(1)

其中，D表示道路密度，A表示区域面积，L表示区域内道路的总长度。

2.2 用D-TIN计算道路密度的方法

根据公式(1)可知，要计算某条道路的密度，就要分别计算其长度L以及对应生长区域的面积A。道路的长度很容易计算，但其辐射区域的面积比较难确定。该区域可以看作是以道路网中每条道路为中心，以相同的速率向外扩张，直到彼此相遇为止而在平面上形成的图形；除最外围的道路形成的是开放区域外，其余每条道路将被一个多边形所包围，这个多边形就是道路所辐射的区域，称之为道路的“生长区”。这一几何构造与Voronoi图类似(Voronoi图又称泰森多边形)，但按照计算几何的严格定义[6]，它并不是Voronoi图。Voronoi图是针对点群目标的，具有以下性质：①剖分元多边形是凸的；②相邻多边形中点的连接得到其对偶Delaunay三角网[7]。显然，道路的“生长区”几何构造不满足这2个性质。从空间相关性上分析，生长区域多边形可以认为是其包含线状道路的生存范围，是与相邻道路竞争的结果，这一特征与二维平面进行D-TIN剖分后形成的结构骨架线表达的区域等分性相似。因此，对道路网进行D-TIN剖分，然后用其结构骨架线来表达道路的生长区域是比较合适的方法。当确定了道路生长区域的边界后，就可以依次连接这些边界构成一个封闭的多边形，然后计算其面积A，从而用公式(1)即可计算每条道路的密度。

根据这一思想，用D-TIN计算道路网密度的方法如下：

(1)以道路的弧段集{Li}作为约束边，对道路网构建D-TIN；

(2)对D-TIN剖分的区域提取骨架线；

(3)对所有骨架线进行线拓扑成区(成区方法参考文献[11])，构成道路弧段生长的多边形集{Pi}；

(4)将{Li}和{Pi}进行空间位置匹配，依次计算出每条弧段Li的长度和对应生长区域Pi的面积，然后用公式(1)计算该弧段的密度。

笔者根据以上步骤进行初步试验，在对道路网数据构建D-TIN后，发现第2步提取剖分区域的骨架线非常复杂，容易漏提和重复提取骨架线；并且所提取的骨架线通常不封闭，导致第3步由线自动生成区时失败，从而无法进行第4步处理。因此，笔者经过仔细研究，提出了一种可以回避第2步和第3步的处理过程，直接计算道路生长区面积的方法。

3 道路生长区面积计算优化

由于道路的密度计算只与道路的长度和面积相关，通过分析道路所关联三角形的面积与要求解生长区面积的关系，得出一种直接计算道路生长区面积的方法。其基本思想是：将道路生长区看作是所关联的三角形切割部分后的集合，采取化整为零的思想，根据每条道路所关联的三角形对生长区面积贡献度的差异，将其分为4类，计算每类三角形对生长区面积的贡献度，然后累加所有关联的三角形的贡献面积，从而直接得到道路生长区的面积。下面详细介绍4类三角形的分法和每类三角形对生长区面积贡献度的计算方法。

3.1 三角形的分类

图1是道路网进行D-TIN剖分后截取的部分区域，P为道路L1的生长区，其通过连接道路L1所关联的三角形的边界的中点所得。根据L1所关联三角形对生长区面积贡献度的差异，可以将其分为以下4类，如图2所示：

Ⅰ型：三角形的三个顶点都落在弧段上，如△abc的三个顶点都落在L1上；

Ⅱ型：三角形有两个顶点落在弧段上，如△cde的两个顶点c和d，△ace的两个顶点a和c落在弧段L1上；

Ⅲ型：三角形只有一个顶点落在弧段上，如△dfe只有一个顶点d落在弧段L1上；

Ⅳ型：三角形有两个顶点落在弧段上，但其中有一个顶点被其它弧段共享为路网的结点，如△aeg的两个顶点a和g都在弧段L1上，g为结点。

图1 道路剖分后形成的生长区及其所关联的4类三角形

图2 各类三角形对道路生长区面积的贡献度

这4类三角形包含了道路网进行D-TIN剖分后形成的所有三角形的类型，根据每类三角形对生长区面积的贡献量，可以计算每类三角形对生长区面积的贡献度，通过累加道路所关联的所有三角形的面积贡献度，即可得到该道路生长区的面积。下面介绍其具体计算方法。

3.2 道路生长区域面积计算

用S表示面积，则L1生长区P的面积可以表示为：

Sp=∑SⅠ+∑SⅡ+∑SⅢ+∑SⅣ

(2)

其中，Sp表示道路生长区的面积，SⅠ、SⅡ、SⅢ、SⅣ分别表示Ⅰ型、Ⅱ型、Ⅲ型、Ⅳ型三角形对生长区P的面积贡献度。

图2中，阴影表示每类三角对生长区面积的贡献度，e和f分别为三角形对应边上的中点。根据三角形的中线特征和面积割补法求解4类三角形的阴影面积如下：

SⅠ=SΔabc；

SⅡ=SΔabc-SΔaef=3/4SΔabc;

SⅢ=SΔaef=1/4SΔabc；

SⅣ=1/2SΔabc;

计算出4类三角形的面积贡献度后，通过道路L1上的坐标点搜索出所有关联的三角形，并将其分类，然后用公式(2)即可得到生长区p的面积。

4 算法步骤

综上所述，道路网密度计算的优化方法如下：

(1)道路网自动相交打断，清除道路重叠坐标及重叠弧段，形成路网弧段集{Li}；

(2)以{Li}的弧段边界为约束，构建约束D-TIN；

(3)依次取{Li}中每条弧段Li，根据Li上的坐标点搜索该弧段关联的所有三角形{Ti}；

(4)根据2.1节的分类规则对{Ti}进行分类，根据2.2节的4类三角形面积计算方法计算{Ti}中每个三角形面积的贡献度，然后用公式(2)计算弧段Li生长区域面积Si；

(5)计算弧段Li的长度，然后根据公式(1)计算该弧段的密度。

遍历完{Li}中所有弧段后即可计算出每条道路的密度。步骤(1)是对道路网数据进行预处理，以防止重叠的弧段和重叠坐标对构建约束D-TIN产生干扰，若数据质量符合要求，可以直接进行第(2)步计算。

5 算法应用实例(实验与分析)

本算法采用在MapGIS K9平台上进行C++二次开发的方式进行试验，以便利用平台提供数据处理和可视化等基本功能。实验数据为某地区1∶10万公路网，如图3所示。图4是对道路网数据构建的D-TIN，利用本文提出的密度计算方法得到结果如图5所示。图中道路颜色越深则表明密度越大，越浅则表明密度越小。从路网的整体分布情况来看，本文提出的道路密度分析算法所计算的结果与路网的稀疏分布程度一致。

图3 实验道路网数据

图4 对道路网构建D-TIN

图5 道路密度计算结果

6 结 论

道路密度分析是制图综合及城市道路规划重要参考指标之一。本文在对道路网区域进行D-TIN剖分的基础上，提出一种直接计算道路生长区面积的方法。该方法避免了构建道路生长区域的中间过程，采取化整为零的思想，将道路生长区域面积转化为对关联的三角形面积贡献度的求解，简化了面积计算过程，提高了算法效率，增强了算法的可靠度。该方法为道路网的选取和城市路网结构合理性的评价，提供了有力的道路密度分析工具。

[1]任慧，周振红，周鑫鑫. 基于RS与GIS的城市道路网密度计算[J]. 计算机辅助工程. 2009, 18(2):85-87.

[2]黄远林. 地图图形综合指标体系框架与图形结构识别研究[D]. 武汉：武汉大学, 2012.

[3]陈波，武芳，钱海忠. 道路网自动选取方法研究[J]. 中国图象图形学报，2008(12): 2388-2393.

[4]胡云岗，陈军，李志林等. 基于网眼密度的道路选取方法[J]. 测绘学报，2007(3): 351-357.

[5]Liu Xing-jian, Ai Ting-hua, Liu Yao-lin. Road Density Analysis Based on Skeleton Partitioning for Road Generalization[J]. Geo-spatial Information Science，2009, 12(2): 110-116.

[6]周培德. 计算几何——算法分析与设计[M]. 北京: 清华大学出版社, 2000.

[7]艾廷华. 基于场论分析的建筑物群的移位[J]. 测绘学报，2004(1): 89-94.

[8]武晓波，王世新，肖春生. Delaunay三角网的生成算法研究[J]. 测绘学报，1999(1): 30-37.

[9]刘学军，龚健雅. 约束数据域的Delaunay三角剖分与修改算法[J]. 测绘学报，2001(1): 82-88.

[10]艾廷华，刘耀林，黄亚锋. 河网汇水区域的层次化剖分与地图综合[J]. 测绘学报，2007(2): 231-236.

[11]齐华. 自动建立多边形拓扑关系算法步骤的优化与改进[J]. 测绘学报，1997, 26(3): 254-260.

Road Density Algorithm Based on Constrained Delaunay Triangulated Irregular Network

Jin Cheng1,2, Aa Xiaoya1,2, Xu Daozhu1,2

1. Xi’an Research Institute of Surveying and Mapping, Xi’an 710054, China;2.State Key Laboratory of Geo-information Engineering, Xi’an 710054, China

Aiming at the problem that it is difficult to calculate the area of road covered region in road density analysis, this paper proposes an algorithm for calculating the area of road covered region based on the constrained Delaunay triangulation irregular network (D-TIN). Taking the idea of breaking up the whole into parts, the triangulations are classified into four types according to the contribution of triangulations associated with road to the area of road covered region, and the method for calculating the area of each type triangulation is presented. Meanwhile the area of the road covered region is calculated by accumulating all the area of triangulations related to the road. Then the road network density is calculated. The experiment results show that the algorithm simplifies the road covered area calculation in road density analysis, also it is reliable and efficient.

map generalization; road density; Delaunay triangulation network

2015-02-12。

国家自然科学基金资助项目(41201469)，地理信息工程国家重点实验室开放基金资助项目(NO.SKLGIE2014-M-4-3)。

金澄(1976—)，男，高级工程师，主要从事地理信息处理与服务方面的研究。

P283.1

A