叠置分析(2):数据结构:topology模型/DCEL

这一节我计划介绍叠置分析算法的基础-数据结构,即拓扑(topology)数据模型,在不同的领域/软件中,对拓扑数据模型的实现有不一样的叫法,比如geos中称之为planargraph(平面图)或者geomgraph,计算几何(Computational Geometry)的书籍以及CGAL中将其称为DCEL(Double Connected Edge List),ArcGIS早期的数据模型称之为Coverage。拓扑数据模型是一种用GeoRelational模型,即显式的表达了空间关系,叠置算法就是建立在这种数据模型之上的。

1)发展趋势的不完全总结

早期的GIS用户一定对拓扑数据模型不陌生,因为拓扑数据模型曾经是GIS软件系统的主流数据模型,比如ArcGIS就曾采用了(现在也支持)converage/E00作为主要的数据模型。topo模型建立的初衷主要有两点:

1.强调空间要素的拓扑关系,将拓扑关系(邻接关系)进行了保存,方便系统进行空间关系分析(GIS软件的基本功能)。

2.通过共享几何达到节约物理存储的目的。因此受到上古系统的偏爱。海量的地理信息数据的存储和分发本来就是早期GIS系统的难题。这里可以提一下E00文件有意识的命名规则。为什么叫E00?这是因为早期的软盘只有1.44MB,因此ARC/INFO设计,可以通过多张软盘(或其他介质)保存一个大文件,名字依次是x.e00,x.e01,x.e02…

这种模型后来被强大的第三代地理数据库模型所取代,今天我们熟知的ESRI FileGDB,ArcSDE,Oracle spatial,PostgreSQL(with PostGIS)都是面向对象的地理信息数据库。这不是说拓扑数据模型不重要了,而是存储模型的变革。拓扑数据模型由持久化保存转向了通过计算生成,在内存中进行保存,如geos与CGAL。拓扑数据模型仍然主宰着多边形布尔运算、拓扑关系查询、overlay等主要的应用。

目前,部分软件仍然定义了采用关系型数据库存储表达拓扑模型。比如postGIS Topology,及Oracle Spatial Topology 就定义了一套拓扑数据模型和对应的数据库函数。这些软件提供了一种持久化保存拓扑模型以及在其上执行分析操作的方法。

2)定义

Topology 是一个数学名词,在GIS领域中,用来指节点(node)、弧段(arc)、面(face)相互关联构成的一种图数据结构,用来表达地理要素的联通与邻接关系。拓扑关系的表达主要借助四张表:

1)弧段-节点关系:表达弧段的起始、结束节点。弧段的几何是折线(polyline)

2)节点-弧段关系

表达节点上关联的弧段的列表。其中节点关联的弧段的个数成为节点的度数(degree)。度数为1的节点也成为悬挂节点(dangle node),在上图中,共有10,11,…,17共7个节点,10,12,15,16,17为悬挂节点。

节点保存的弧段列表的顺序不是无序的,而是按照时针顺序保存的。假设按照顺时针顺序保存数据,节点11的边表为1,2,3,节点13的边表为3,5,4。… 在后面具体实现的博客中,我会给大家介绍一下geos对采用弧段末端象限与std::sort生成边表的算法,非常优雅。

3)弧段-面关系

表达弧段左右两侧的面的信息。

细心的你可能发现了一个非常特殊的面,就是面A。我们假设所有的面充斥了整个平面,因此除了B,C,D,E,F之外,还有他们的补集面A,面A比较特殊,没有外轮廓。在CGAL中,面A被称为unbounded face ,可以通过 ArrangementDcelFace::is_unbounded()方法判断面是否为unbounded face。

4)面-弧段关系

表达面由哪些弧段构成。

面的边表也要按照一定的顺序,首先要按照邻接边的顺序记录。其次,要确定边的顺序为你是逆时针(CCW)还是顺时针(CW)。逆时针顺序表示边的左侧为面的内部。顺时针则相反。一般的惯例为CCW,如CGAL,geos的设定。上图是ArcGIS的设定,为CW顺时针顺序。ArcGIS在其软件中一贯按照边的右侧为多边形的内部计算,这点与开源软件不同。

另外,F含有一个岛,一个面可能含有0个或若干岛屿,需要在记录时加以区别。如上图,在岛屿和外环之间加入了“0”作为间隔符,这样做也是为了保证每个面的边表只含有一个List。

我们可以看到,在弧段-面关系表中,弧段是有方向的,但是在面-弧段关系示意图中,弧段的方向被刻意隐去不表达。其实,为了解决这个问题,可以在边的序号上加上正负号表达方向。正号表示同向,负号表示反方向。这种表达方式见于ESRI E00与Oracle spatial Topology。实际上,在逻辑上产生了有向边的概念。对于每一条边,其都含有两个方向不同的有向边。goes中的类型为Directed Edge,CGAL中的概念为HalfEdge(半边)。如下图所示:

3)实现

根据定义,一共有节点(Node),弧段(Edge),有向弧段(DirectedEdge),面(Face)四种实体。当然,可以直接按照定义的方式定义上述几种数据结构,实现拓扑模型。比如geos的实现就与上述定义的方式相同。但是,从拓扑模型的构建与查询算法上考虑,CGAL, PostGIS,Oracle Spatial等软件提出了一种更合理的变体,将边表转换为边的邻接关系表,不仅节约了内存空间,使用起来还非常方便。

以ORACLE的实现为例:

1.边结构

  • START_NODE_ID,END_NODE_ID:起止node,与定义相同。
  • LEFT_FACE_ID,RIGHT_FACE_ID:左右多边形,与定义相同。
  • NEXT_LEFT_EDGE_ID:按正方向逆时针转动的下一条EDGE。
  • NEXT_RIGHT_EDGE_ID:按反方向逆时针转动的下一条EDGE。
  • PREV_LEFT_EDGE_ID:按正方向上一条EDGE。
  • PREV_RIGHT_EDGE_ID:按反方向上一条EDGE。

看上图,边E4正方向逆时针转动的下一条Edge为-E5,因此其NEXT_LEFT_EDGE_ID为-E5。同理,其PREV_LEFT_EDGE_ID为E3。NEXT_RIGHT_EDGE_ID为E2,NEXT_PREV_EDGE_ID为-E6。

2.节点结构

  • EDGE_ID:表示节点所关联的任意一个边的ID。

按照定义,节点记录的是边表,但是这里为啥一个边就够了?这正是这个设计的精妙之处。如果已知一个边,求在以相同的START_NODE上的其他边,只需要调用 this->opposite()->next_left_edge()方法就OK了。其中opposite()采用的是CGAL的叫法,指的是有向边的对向边(twin)。

3.面结构

  • BOUNDARY_EDGE_ID: 同理,面的外环不需要一个边表,只需要任意一个EDGE的ID就好了,可以通过循环调用 next_left_edge() 实现对所有边的遍历。
  • ISLAND_EDGE_ID_LIST:因为面有可能有0到多个岛屿,因此采取每个岛屿记录一个EDGE的方式构成一个ID_LIST。

这种改良的拓扑模型表达方式确实具有简洁、应用方便的优点。


需要说明的是,PostGIS,Oracle,CGAL均对拓扑模型进行了扩展,支持面中包含孤立点(isolated node),因此,会出现有的node的EDGE_ID为空,并且增加FACE_ID表达孤立点的情况。因为没有在定义中涉及到,我这里暂时不做表述。


参考文献

PostGIS Topology

Oracle Spatial : Spatial Topology and Network Data Models

ArcGIS Desktop Help: Coverage Topology

叠置分析(2):数据结构:topology模型/DCEL》上有1条评论

发表评论

电子邮件地址不会被公开。