基于ＭａｐＩｎｆｏ的线面拓扑分析算法与实现

摘　要:由于不具备拓扑关系的数据结构,ＭａｐＩｎｆｏ的拓扑分析能力受到限制。为了进行线面拓扑关系判断,通过算法分析提出了两种高效的建立拓扑关系的方案,并加以实现。
　引言
在ＧＩＳ的应用中要求能表达复杂的图形关系,有较多的图形运算及操作功能、ＤＥＭ分析功能等,而在土地资源管理中以二维图形处理为主。ＭａｐＩｎｆｏ开发平台具备了所有ＧＩＳ开发系统的基本功能,目前许多软、硬件厂家对ＭａｐＩｎｆｏ提供了专门支持,解决了开发技术的后顾之忧。ＭａｐＩｎｆｏ已经成为地理信息系统领域应用最广泛的的工具之一,具有良好的发展前景。
在基于ＧＩＳ技术的土地资源管理信息系统中,常常要求计算图斑的净面积(即图斑面积减去其中包含的线状、零星地物的面积)。于是怎样获得拓扑信息,判断一个线状物与几个图斑相邻接,以及与哪几个图斑相邻接成为解决问题的关键。本文涉及的系统,采用在Ｄｅｌｐｈｉ程序中链接嵌入ＭａｐＩｎｆｏ的方式进行开发。
2　ＭａｐＩｎｆｏ中的拓扑关系拓扑关系是指图形保持连续状态下变形,但图形关系不变的性质,描述了图形实体之间的邻接、关联和包含关系。ＧＩＳ中空间数据的拓扑关系,对数据处理和空间分析具有重要的意义,因为:拓扑关系能清楚地反映实体之间的逻辑结构关系;利用拓扑关系有利于空间要素的查询;可以根据拓扑关系重建地理实体。ＭａｐＩｎｆｏ作为面向应用的ＧＩＳ产品,通过采用"空间实体模型+空间索引"的描述方法,来建立空间对象的拓扑关系。
这种拓扑结构的基础是"空间实体"。空间实体是地理实体的抽象,主要包括点、线、面三种基本类型。任一空间实体,都是由一个或多个"部分"组成的,而"部分"则由"点集"组成,"点集"是若干具有固定ｘ、ｙ坐标的节点集合。这样,在一个实体对象内部,存储了其全部节点的坐标信息和构成关系(点集、部分的排列顺序和特性),也就完成了对其自身的拓扑模型的建立。
在对空间信息进行各种操作的过程中需要经常提取和检索空间数据,"空间索引"的目的就是为了加快空间定位数据的存取速度,尤其对大型空间数据库而言更是如此。为此ＭａｐＩｎｆｏ采用了Ｒ-Ｔｒｅｅ树技术。简单地说,Ｒ-Ｔｒｅｅ将各种空间实体的ＭＢＲ(最小外接矩形)存储在索引中,并按从大到小的顺序进行空间索引搜索。这样,能根据给定的坐标,快速定位到某一空间对象。该技术实现了一种动态的"拓扑关系"机制。只有在需要时,系统才根据空间索引,建立并使用实体间的拓扑关系,而不具备表达拓扑关系的数据结构。
在ＭａｐＩｎｆｏ中进行空间分析有多种方法,比较常用的,也是最有特色的,是通过ＳＱＬ语言来进行。因为ＳＱＬ语言是关系数据库的标准查询语言,而ＭａｐＩｎｆｏＰｒｏ也是以关系数据库为中心的,因此整个系统的功能能统一用ＳＱＬ来实现。ＭａｐＩｎｆｏ将拓扑实现的细节隐藏在系统中,通过在ＳＱＬ语言中调用其定义的空间分析函数、运算符,"隐式"地利用拓扑关系进行空间分析,可以说是实现了一种"隐式"的拓扑关系。
3　问题描述土地资源管理的图形运算中主要涉及图斑、线状地物、零星地物三个图层。线状地物层中的每个线状物被处理为仅有以下两种情况:悬于一个图斑中或位于两个图斑的交界线上,即一个线状物仅与一个或两个图斑邻接,称为长线断链。应用系统中要求实现对图形属性数据的提取和入库,具体如下:在提取图斑的属性数据入库的同时,统计并输入图斑的净面积;零星地物面积的扣除规则为:图斑面积减去该图斑包含的所有零星地物的面积;线状地物面积的扣除规则为:若线状物与两个图斑邻接,这两个图斑在计算净面积时,各扣除该线状物一半的面积;若线状物与一个图斑邻接,图斑净面积等于图斑面积减去线状物的全面积。
4　算法分析要达到上述要求,需判断一个线状地物与几个及与哪几个图斑相邻接。由于ＭａｐＩｎｆｏ中并没有"显式"的拓扑信息,要获得线-面拓扑关系,只有通过调用相关的空间分析函数和运算符来实现。此处需在ＭａｐＢａｓｉｃ的ＳＱＬ语句中用到以下几个图形对象关系运算算子和函数:
1)ＩｎｔｅｒＳｅｃｔｓ:两个对象在某些点上相交,它们至少有一个公共点;
2)ＯｖｅｒＬａｐ:计算两个对象相交部分,并以对象形式返回相交部分;如果一个对象是线而另一个对象是封闭的面,ＯｖｅｒＬａｐ()函数返回被封闭区域所覆盖的线性对象;
3)ＯｂｊｅｃｔＬｅｎ:返回线或折线对象的长度。则线面有邻接关系的条件可描述为:一条线与一个或两个面相交,且线面相叠部分的长度大于零。若采用每次判断时临时生成拓扑信息的方法,则需在入库过程中即针对图斑对象的循环中,先找出与当前图斑相邻接的线状地物,再对这些线状地物进行循环,求出这些线状地物各与几个图斑相邻接,根据其邻接的图斑个数进行处理。其程序流程如图1。图1
　动态生成拓扑信息的程序流程图
分析该算法的时间复杂度：
设图斑个数为Ｎ;设位于两个图斑交界线上的线状地物个数为Ｄ,则对此类线状物执行模块Ｂ的次数为2 Ｄ;设悬于一个图斑中的线状地物个数为Ｓ,则对此类线状物执行模块Ｂ的次数为Ｓ;线状地物个数Ｘ=Ｄ+Ｓ。图1中的模块Ａ和Ｂ的结构相同,只是操作对象不同:模块Ａ根据面对象在线状层表中查找线对象,模块Ｂ根据线对象在图斑层表中查找面对象。若将Ａ、Ｂ中的"查找"抽象为基本操作,则模块Ａ的时间复杂度为Ｘ,模块Ｂ的时间复杂度为Ｎ。整个算法的时间复杂度
Ｔ(Ｎ,Ｘ)=Ｎ*Ｘ+2Ｄ*Ｎ+Ｓ*Ｎ
=2Ｎ*Ｘ+Ｄ*Ｎ
=Ｎ (2Ｘ+Ｄ)≈3Ｎ*Ｘ　　　　　　(Ｄ≈Ｘ)
从以上分析可看出,对位于两个图斑交界线上的线状地物重复执行了模块Ｂ。由于此类线状物所占比例较高,再加上以ＯＬＥ方式调用ＭａｐＩｎｆｏ的效率较低,该算法浪费了较多的系统资源,需构造效率更高的算法。于是,可以考虑先将拓扑信息"显示"化和"固定"化,再进行处理。即以空间换时间,把拓扑信息以某种形式保存下来。鉴于拓扑信息不同的生成和保存方式,本文提出以下两个实现方案。
5　算法的实现
5.1　方案1在线状层的表结构中扩展必需的字段,将线面拓扑信息保存在线状层中。图斑层需有标识字段(如:Ｕｉｄ)。线状层扩展的字段有:
ＰｌａｎｅＮｕｍ---记录线状对象邻接的图斑个数;
Ｐｌａｎｅ1Ｉｄ---记录第1个图斑的标识码;
Ｐｌａｎｅ2Ｉｄ---记录第2个图斑的标识码。
把图1中的模块Ｂ转移到"预处理过程"(对线状物的循环)中,将每个线状物邻接的图斑数及图斑的标识码填入相应字段。这样,避免了在入库过程中对位于两个图斑交界线上的线状物进行重复的拓扑判断操作。
关键步骤的程序实现如下
:预处理过程
ＲｅｃＮｕｍ：=Ｍａｐｉｎｆｏ．ＥＶＡＬ（′ｔａｂｌｅＩｎｆｏ（′+ＸＺＬ-Ｔａｂｌｅ+′，8）′）；
／／取线状层表中记录个数，ＸＺＬ-Ｔａｂｌｅ为线状层表名
ＭａｐＩｎｆｏ．Ｄｏ（′ＦｅｔｃｈＦｉｒｓｔＦｒｏｍ′+ＸＺＬ-Ｔａｂｌｅ）；ｆｏｒｋ：=1ｔｏＲｅｃＮｕｍｄｏ
／／循环线状记录，选择当前线状物相交且相叠部分
／／长度大于零的图斑到临时表Ｓｅｌｔｂｐ2ｂｅｇｉｎＭａｐＩｎｆｏ．Ｄｏ（′Ｂｏｂｊ=′+ＸＺＬ-Ｔａｂｌｅ+′．Ｏｂｊ′）；
／／Ｂｏｂｊ为对象变量
ＭａｐＩｎｆｏ．Ｄｏ（′Ｓｅｌｅｃｔ ｆｒｏｍ′+ＴＢＰ-Ｔａｂｌｅ+′ＷｈｅｒｅＯｂｊＩｎｔｅｒｓｅｃｔｓＢｏｂｊｉｎｔｏＳｅｌｔｂｐ1′）；
ＩｆＭａｐｉｎｆｏ．ＥＶＡＬ（′ＳｅｌｅｃｔｉｏｎＩｎｆｏ（3）′）>0Ｔｈｅｎｂｅｇｉｎ
ＭａｐＩｎｆｏ．ｄｏ（′Ｓｅｌｅｃｔ Ｆｒｏｍ Ｓｅｌｔｂｐ1 Ｗｈｅｒｅ ＯｂｊｅｃｔＬｅｎ（Ｏｖｅｒｌａｐ（Ｂｏｂｊ，ｏｂｊ），″ｍ″）>0ｉｎｔｏＳｅｌｔｂｐ2′）；
ＴｂｐＳｅｌＮｕｍ：=ＭａｐＩｎｆｏ．ＥＶＡＬ（′ＳｅｌｅｃｔｉｏｎＩｎｆｏ（3）′）；
ＭａｐＩｎｆｏ．Ｄｏ（′ＦｅｔｃｈＦｉｒｓｔＦｒｏｍＳｅｌｔｂｐ2′）；ｆｏｒＬ：=1ｔｏＴｂｐＳｅｌＮｕｍｄｏ／／循环图斑记录
ｂｅｇｉｎ
Ｔｂｉｄ：=ＭａｐＩｎｆｏ．ＥＶＡＬ（′Ｓｅｌｔｂｐ2．Ｕｉｄ′）；／／取图斑的标识码
Ｍａｐｉｎｆｏ．ｄｏ（′ｕｐｄａｔｅ′+ＸＺＬ-Ｔａｂｌｅ+′ｓｅｔＰｌａｎｅ′+ＩｎｔＴｏＳｔｒ（Ｌ）+′=′+Ｔｂｉｄ+′，ＰｌａｎｅＮｕｍ=′+ＩｎｔＴｏＳｔｒ（ＴｂｐＳｅｌＮｕｍ）+′ＷｈｅｒｅＲｏｗｉｄ=′+ＩｎｔＴｏＳｔｒ（ｋ））；
／／在线状记录中填写拓扑信息
ＭａｐＩｎｆｏ．Ｄｏ（′ＦｅｔｃｈｎｅｘｔＦｒｏｍＳｅｌｔｂｐ2′）；ｅｎｄ；ｅｎｄ；
ＭａｐＩｎｆｏ．Ｄｏ（′ＦｅｔｃｈＮｅｘｔＦｒｏｍ′+ＸＺＬ-Ｔａｂｌｅ）；ｅｎｄ；
入库过程的修改结果
Ｔｂｉｄ：=Ｍａｐｉｎｆｏ．ＥＶＡＬ（ＴＢＰ-Ｔａｂｌｅ+′．Ｕｉｄ′）；
／／提取当前图斑的标识码
ＭａｐＩｎｆｏ．Ｄｏ（′Ｓｅｌｅｃｔ Ｆｒｏｍ′+ＸＺＬ-Ｔａｂｌｅ+′Ｗｈｅｒｅ（Ｐｌａｎｅ1Ｉｄ=″′+Ｔｂｉｄ+′″ｏｒＰｌａｎｅ2Ｉｄ=″′+Ｔｂｉｄ+′″）ＩｎｔｏＳｅｌＸｚｌ′）；
／／选择与当前图斑邻接的线状物到临时表ＳｅｌＸｚｌＰｓ：=Ｍａｐｉｎｆｏ．ＥＶＡＬ（′ＳｅｌｅｃｔｉｏｎＩｎｆｏ（3）′）；ＩｆＰｓ>0ｔｈｅｎｂｅｇｉｎＦｏｒＫ：=1ｔｏＰｓｄｏ
／／循环线状记录
ｂｅｇｉｎ
ＭａｐＩｎｆｏ．Ｄｏ（′ＦｅｔｃｈＲｅｃ′+ＩｎｔＴｏＳｔｒ（Ｋ）+′ＦｒｏｍＳｅｌＸｚｌ′）；
Ｉ：=ＭａｐＩｎｆｏ．ＥＶＡＬ（′ＳｅｌＸｚｌ-2．ＰｌａｎｅＮｕｍ′）；
ＩｆＩ=1Ｔｈｅｎ ；
／／如果与1个图斑对象邻接，就扣除其全部面积并入库
ＩｆＩ=2Ｔｈｅｎ ；
／／如果与2个图斑对象邻接，就扣除一半面积并入库
ｅｎｄ；
ｅｎｄ；
该算法除了避免对线状物进行重复的拓扑判断操作,由于拓扑信息的事先生成,还使得图1中的模块Ａ由复杂的空间分析查询过程变成了简单的属性查询过程。若将进行属性查询生成临时表的时间复杂度忽略,仅就"查找"操作而言,整个算法的时间复杂度降到了Ｎ Ｘ,显著地提高了程序的执行效率。
5.2　方案2通过两表关联生成拓扑信息,并将关联结果保存为.Ｔａｂ文件,供以后查询。数据库中可以通过关联将两个表连接在一起,结果表的列是两个表的列组合,行是匹配的行。ＭａｐＩｎｆｏ扩展了连接的功能,支持地理条件的连接。例如,不用数字进行匹配,而是利用一个表中的对象含有另一个表中的对象连接成一个新表,尤其是当没有匹配字段时,这一特性非常有用。该方案不需要在已有表中扩展字段,只需线状物和图斑拥有各自的标识码,如线状层表的标识码字段为"Ｘｚｂｈ",图斑层表的标识码字段为"Ｕｉｄ"。以线面邻接条件连接图斑和线状层表,连接查询表的每条记录代表了一个线状物与一个图斑的邻接关系,而同一线状物的标识码在连接查询表中出现的次数则代表了该线状物与几个图斑相邻接。该算法中,预处理过程和入库过程的实现代码如下。
预处理过程
ＭａｐＩｎｆｏ．ｄｏ（′ｓｅｌｅｃｔＸｚｂｈ，Ｕｉｄｆｒｏｍ′+ＸＺＬ-Ｔａｂｌｅ+′，′+ＴＢＰ-Ｔａｂｌｅ+′ｗｈｅｒｅ′+ＸＺＬ-Ｔａｂｌｅ+′．ｏｂｊＩｎｔｅｒｓｅｃｔｓ′+ＴＢＰ-Ｔａｂｌｅ+′．ｏｂｊａｎｄｏｂｊｅｃｔｌｅｎ（ｏｖｅｒｌａｐ（′+ＸＺＬ-Ｔａｂｌｅ+′．ｏｂｊ，′+ＴＢＰ-Ｔａｂｌｅ+′．ｏｂｊ），″Ｍ″）>0ｉｎｔｏＸＭＴｅｍｐ′）；
／／连接线状层表和图斑层表，产生查询表ＸＭＴｅｍｐＭａｐＩｎｆｏ．ｄｏ（′ＣｏｍｍｉｔｔａｂｌｅＸＭＴｅｍｐａｓ″′+ｓＰａｔｈ+′ＸＭＴｅｍｐ．ｔａｂ″′）；
／／将ＸＭＴｅｍｐ保存到磁盘
ＭａｐＩｎｆｏ．ｄｏ（′ＣｌｏｓｅｔａｂｌｅＸＭＴｅｍｐ′）；
ＭａｐＩｎｆｏ．ｄｏ（′Ｏｐｅｎｔａｂｌｅ″′+ｓＰａｔｈ+′ＸＭＴｅｍｐ．ｔａｂ″′）；
入库过程
Ｔｂｉｄ：=ＭａｐＩｎｆｏ．Ｅｖａｌ（ＴＢＰ-Ｔａｂｌｅ+′．Ｕｉｄ′）；
ＭａｐＩｎｆｏ．ｄｏ（′Ｓｅｌｅｃｔ ｆｒｏｍＸＭＴｅｍｐＷｈｅｒｅＵｉｄ=″′+Ｔｂｉｄ+′″ｉｎｔｏＳｅｌＸｚｌ′）；
／／选择与当前图斑邻接的线状物到临时表ＳｅｌＸｚｌ
Ｐｓ：=Ｍａｐｉｎｆｏ．ＥＶＡＬ（′ＳｅｌｅｃｔｉｏｎＩｎｆｏ（3）′）；
ＩｆＰｓ>0 ｔｈｅｎ
ｂｅｇｉｎ
ＦｏｒＫ：=1ｔｏＰｓｄｏ
ｂｅｇｉｎ
ＭａｐＩｎｆｏ．ｄｏ（′Ｆｅｔｃｈｒｅｃ′+ＩｎｔＴｏＳｔｒ（Ｋ）+′ＦｒｏｍＳｅｌＸｚｌ′）；
Ｘｚｉｄ：=ＭａｐＩｎｆｏ．Ｅｖａｌ（′ＳｅｌＸｚｌ．ｘｚｂｈ′）；
ＭａｐＩｎｆｏ．ｄｏ（′ＳｅｌｅｃｔＣｏｕｎｔ（ ）″Ｃｏｕｎｔ″ｆｒｏｍＸＭＴｅｍｐｗｈｅｒｅｘｚｂｈ=″′+Ｘｚｉｄ+′″ｉｎｔｏｔｅｍｐ′）；
／／统计该线状物出现的次数
Ｉ：=ＭａｐＩｎｆｏ．Ｅｖａｌ（′Ｔｅｍｐ．Ｃｏｕｎｔ′）；ＩｆＩ=1Ｔｈｅｎ ；ＩｆＩ=2Ｔｈｅｎ ；
ｅｎｄ；
ｅｎｄ；
对于该方案要特别注意,在预处理过程中,要将查询表ＸＭＴｅｍｐ保存到磁盘形成基表,以后的查询在这张基表上进行。这是因为,对ＭａｐＩｎｆｏ的查询表的操作,与相关基表紧密相联,在关联查询表上进行查询操作,同时要在多个基表中进行检索。这样,查询效率非常低。而将查询表保存到磁盘,则切断了它同原基表的联系,查询效率大大提高。
6　结论
本文提出的两个方案均能在图形属性入库过程中有效地进行线-面拓扑关系判断。以包含1328个图斑、933个线状物的图幅为样本,在ＣＰＵ为Ｋ6-2300、内存128Ｍ的计算机上进行对比实验(入库过程中还包括许多其它操作),如表1所示,方案1、2较图1所示算法效率明显提高;方案2的入库效率较方案1略胜一筹,但其预处理效率是方案1的8.5倍。读者可以根据不同需求选择不同的方案。目前,方案2已成功应用到某土地资源管理信息系统中。表1　各解决方案的执行效果比较算法预处理时间入库时间总的时间图1所示算法---22′30″22′30″方案12′48″12′49″15′37″方案220″11′23″11′43″虽然本文仅讨论了线-面拓扑分析,但其算法思想,也可以运用到其它拓扑分析中去。