等級(jí)Voronoi圖及加權(quán)Voronoi圖的ArcGIS矢量生成算法

田 松，崔 希 民，孫 云 華，宮 宇

(1.中國(guó)礦業(yè)大學(xué)(北京)地球科學(xué)與測(cè)繪工程學(xué)院，北京 100083；2.國(guó)家海洋局南海分局，廣東 廣州 510310)

?

等級(jí)Voronoi圖及加權(quán)Voronoi圖的ArcGIS矢量生成算法

田 松1，崔 希 民1，孫 云 華1，宮 宇2

(1.中國(guó)礦業(yè)大學(xué)(北京)地球科學(xué)與測(cè)繪工程學(xué)院，北京 100083；2.國(guó)家海洋局南海分局，廣東 廣州 510310)

等級(jí)Voronoi圖及加權(quán)Voronoi圖是以樹狀形式表達(dá)事物層次關(guān)系的方法，可作為一種空間數(shù)據(jù)模型應(yīng)用于城鎮(zhèn)等級(jí)體系、組織結(jié)構(gòu)關(guān)系表達(dá)等地學(xué)領(lǐng)域。該文首次提出一種等級(jí)Voronoi圖及加權(quán)Voronoi圖的ArcGIS矢量生成算法，以ArcGIS Engine為開發(fā)工具，利用增量法思想，通過區(qū)域分割和合并方法實(shí)現(xiàn)了對(duì)空間的多級(jí)劃分，為研究和發(fā)展GIS空間數(shù)據(jù)模型提供了重要的方法手段。

等級(jí)Voronoi圖；等級(jí)加權(quán)Voronoi圖；ArcGIS；區(qū)域分割；區(qū)域合并

0 引言

等級(jí)Voronoi圖(Hierarchical Voronoi Diagram)亦稱作分層Voronoi圖，簡(jiǎn)稱等級(jí)V圖或分層V圖,是對(duì)空間區(qū)域多級(jí)劃分，以樹狀形式表示區(qū)域?qū)哟侮P(guān)系的方法，是Voronoi圖理論和應(yīng)用的擴(kuò)展。等級(jí)V圖是一種數(shù)據(jù)結(jié)構(gòu)，用于提高數(shù)據(jù)組織、檢索及可視化效率，亦是一種空間數(shù)據(jù)模型，用于表達(dá)空間事物間的等級(jí)層次關(guān)系。Farin等[1]利用等級(jí)V圖進(jìn)行地形建模;Gold等[2]將等級(jí)V圖用于地形數(shù)據(jù)的組織與LOD顯示;Wallgrün[3]利用等級(jí)V圖解決了機(jī)器人路徑搜索問題；趙威、陳笑筑、薛瑩等[4-6]分別利用等級(jí)V圖表達(dá)了開封市、畢節(jié)市、蘇浙滬城市群的城市等級(jí)體系。

等級(jí)加權(quán)Voronoi圖(等級(jí)加權(quán)V圖)是等級(jí)Voronoi圖的擴(kuò)展，是在考慮母點(diǎn)權(quán)重的情況下，劃分多級(jí)空間的方法。Balzer等[7]研究了等級(jí)V圖及加法加權(quán)等級(jí)V圖的生成算法，并用于軟件度量可視化；Sud等[8]利用CPU加速技術(shù)實(shí)現(xiàn)了加法加權(quán)的等級(jí)V圖；Kang等[9]利用等級(jí)加權(quán)V圖進(jìn)行空間數(shù)據(jù)的組織和管理。目前，尚無結(jié)合GIS軟件或工具實(shí)現(xiàn)等級(jí)V圖及加權(quán)圖的方法。本文在文獻(xiàn)[10]方法的基礎(chǔ)上，利用ArcGIS Engine(簡(jiǎn)稱AE)實(shí)現(xiàn)了等級(jí)V圖及加權(quán)V圖的矢量生成方法，旨在建立一種適用于地學(xué)領(lǐng)域的空間數(shù)據(jù)模型，為表達(dá)空間事物的多級(jí)層次關(guān)系提供一種良好的可視化手段。

1 定義

1.1 等級(jí)V圖定義[9]

定義1：Vh(P)={Vh(p1),Vh(p2),…,Vh(pn)}是以pi(i=1,2,…,n)為母點(diǎn)的Voronoi圖，h表示等級(jí)，h=1時(shí)等級(jí)最高，Vh(pi)表示母點(diǎn)pi的Voronoi區(qū)域(簡(jiǎn)稱V區(qū)域)，令點(diǎn)集Qi={qj|qj∈Vh(pi),j≥1,1≤i≤n}，若Voronoi圖Vh+1(Qi)和Vh(pi)存有如下關(guān)系：

Vh+1(Qi)={Vh+1(qj)|Vh+1(qj)?Vh(pi),

qj∈Vh(pi),j≥1,1≤i≤n}

(1)

且有∑Vh+1(qj)≡Vh(pi)，則稱Vh(P)為等級(jí)Voronoi圖。公式為：

Vh(P)={Vh+1(Q1),Vh+1(Q2),…,Vh+1(Qn)}

(2)

1.2 等級(jí)加權(quán)V圖定義

定義2：Vh(P,λ)={Vh(p1,λ1),Vh(p2,λ2),…,Vh(pn,λn)}是以pi(i=1,2,…,n)為母點(diǎn)、λi為權(quán)重的加權(quán)Voronoi圖，Vh(pi,λi)表示母點(diǎn)pi的加權(quán)Voronoi區(qū)域(簡(jiǎn)稱加權(quán)V區(qū)域)，令點(diǎn)集Qi={qj|qj∈Vh(pi,λi),j≥1,1≤i≤n}，若加權(quán)Voronoi圖Vh+1(Qi,λ)和Vh(pi,λi)存在如下關(guān)系：

Vh+1(Qi,λ)={Vh+1(qj,λj)|Vh+1(qj,λj)?Vh(pi,λi),

qj∈Vh(pi,λi),j≥1,1≤i≤n}

(3)

且有∑Vh+1(qj,λj)≡Vh(pi,λi)，則稱Vh(P,λ)為等級(jí)加權(quán)Voronoi圖。公式為：Vh(P,λ)={Vh+1(Q1,λ),Vh+1(Q2,λ),…,Vh+1(Qn,λ)}

(4)

為方便后續(xù)描述，均用Vh表示第h級(jí)V圖或加權(quán)V圖，Vh(i)表示母點(diǎn)pi的V區(qū)域或加權(quán)V區(qū)域。

2 等級(jí)V圖、加權(quán)V圖ArcGIS矢量生成算法

2.1 算法基本思想

以等級(jí)加權(quán)V圖為例，設(shè)P={p1,p2,…,pn}是二維歐氏空間上的n個(gè)母點(diǎn)，λ為母點(diǎn)權(quán)重，h為母點(diǎn)等級(jí)，相同h值的母點(diǎn)處于同一等級(jí)。算法基本思想：生成點(diǎn)集最小外包矩形，用h=1的母點(diǎn)集合劃分最小外包矩形，生成加權(quán)V圖V1；判斷h=2的點(diǎn)集歸屬，并利用該點(diǎn)集劃分V1，生成第二層的加權(quán)V圖V2，以此類推，直至算法結(jié)束。

加權(quán)V圖生成方法為：利用增量法思想，每插入一個(gè)點(diǎn)，對(duì)原有加權(quán)區(qū)域進(jìn)行重新劃分，分割出屬于新插入點(diǎn)的加權(quán)V區(qū)域，當(dāng)所有點(diǎn)插入完成后，即得到最終的加權(quán)V圖。等級(jí)V圖是等級(jí)加權(quán)V圖權(quán)重相同時(shí)的特例，二者算法實(shí)現(xiàn)思想相同。

2.2 算法實(shí)現(xiàn)步驟

Input:P={p1,p2,…,pn}，pi∈P，包含屬性：ID，pi的唯一標(biāo)識(shí)(對(duì)應(yīng)圖1中的Belong字段)；Weight，pi權(quán)重(對(duì)應(yīng)圖1中的Weight字段)，與變量λ意義相同；Level，pi等級(jí)(對(duì)應(yīng)圖1中的Level字段)，與變量h意義相同。

Output:等級(jí)V圖、等級(jí)加權(quán)V圖。

圖1 ArcGIS點(diǎn)集示例

Fig.1 A set of points saved in ArcGIS attribute table

算法具體步驟如下：

(1)數(shù)據(jù)獲取及預(yù)處理。獲取母點(diǎn)集合P，并按其Level值劃分為不同等級(jí)的母點(diǎn)集合pL={pL1,pL2,…,pLm}，pLh∈pL，∑pLh≡P，pL1表示最高層次點(diǎn)集，pLm為最低層次點(diǎn)集。

(2)獲取P的最小外包矩形MBR。初始化加權(quán)V圖結(jié)構(gòu)體V0，并令V0=MBR。

(3)遞歸循環(huán)pLh集合，生成等級(jí)加權(quán)V圖，偽代碼如下：

for(inth=1;h≤m;h++){

//獲取當(dāng)前等級(jí)點(diǎn)集

pLh=GetLevelPoints(P,h)

//將pLh歸屬不同加權(quán)V區(qū)域

Vh-1=SetPnts2BelongRegion(pLh,Vh-1)

//初始化加權(quán)V圖結(jié)構(gòu)體Vh,_Vh

for(intj=1;j≤VorCount;j++){

//VorCount：當(dāng)前加權(quán)V區(qū)域個(gè)數(shù)

Vh←Vh-1(j)

for(intt=2;t≤VorPntsCount;t++){

//VorPntsCount：加權(quán)V區(qū)域內(nèi)包含的點(diǎn)個(gè)數(shù)

for(intk=1;k

if(pt.Weight=pk.Weight)

//垂直平分線劃分區(qū)域

Vh(tk)=DivideRegionWithPB(pt,pk,Vh(k))

else//圓弧劃分區(qū)域

Vh(tk)=DivideRegionWithArc(pt,pk,Vh(k)) }

//區(qū)域合并

Vh(tK)={Vh(t1),Vh(t2),…,Vh(tk),…,Vh(tt-1)}

Vh(t)=UnionRegion(Vh(tK))}

Vh←Vh

_Vh=null}

Vh=_Vh

//輸出h級(jí)加權(quán)V圖

ExportLevelVorRegion(Vh)}

等級(jí)V圖偽代碼和加權(quán)V圖類似，由于不需要判斷權(quán)重，所以只需保留垂直平分線劃分區(qū)域的方法，即代碼改變?nèi)缦拢?/p>

for(intk=1;k

Vh(tk)=DivideRegionWithPB(pt,pk,Vh(k))}

其他代碼不變。

2.3 算法說明

2.3.1 點(diǎn)集的區(qū)域歸屬判定 第一次循環(huán)，SetPnts2BelongRegion()方法認(rèn)定MBR為點(diǎn)集pL1的V區(qū)域或加權(quán)V區(qū)域V0，利用pL1劃分V0生成V1，自第二次循環(huán)起，將pLh中的點(diǎn)歸屬于Vh-1(j)(j=1,2,…,VorCount)。點(diǎn)集歸屬判別方法：如果pLh中的一個(gè)點(diǎn)p∩Vh-1(j)≠Ф，則認(rèn)為點(diǎn)p屬于Vh-1(j)。

2.3.2 區(qū)域分割和合并 假定pi(xi,yi)、pj(xj,yj)為空間中兩個(gè)不重疊的點(diǎn)，λi、λj是兩點(diǎn)權(quán)重。根據(jù)文獻(xiàn)[10,11]給出的V圖及加權(quán)V圖定義可知，兩點(diǎn)間Voronoi邊為一直線，其公式如下：

(5)

兩點(diǎn)間加權(quán)Voronoi邊為一圓弧，其公式如下：

(6)

圓心為：

(7)

半徑為：

(8)

式中:d(pi,pj)表示pi和pj間的歐氏距離。

以此數(shù)學(xué)理論為基礎(chǔ)，便可利用AE接口實(shí)現(xiàn)點(diǎn)集的加權(quán)V區(qū)域劃分，具體包括4個(gè)要點(diǎn)：直線生成、圓弧生成、V區(qū)域及加權(quán)V區(qū)域分割、V區(qū)域及加權(quán)V區(qū)域合并。

(1)直線生成。AE提供ILine接口，包含F(xiàn)romPoint和ToPoint屬性用于生成直線，表示通過設(shè)定起始點(diǎn)和終止點(diǎn)坐標(biāo)(通過式(5)和MBR邊界求交獲得)，生成一條兩點(diǎn)間的線段。該方法的具體實(shí)現(xiàn)包含在DivideRegionWithPB()方法中。

(2)弧段生成。IConstructCircularArc.ConstructCircle(IPointCenterPoint,doubleRadius)方法可用于生成弧段，其中圓心參數(shù)CenterPoint由式(7)獲得，半徑參數(shù)Radius由式(8)獲得。該方法的具體實(shí)現(xiàn)包含在DivideRegionWithArc()方法中。

(3)區(qū)域分割。如圖2所示，假設(shè)Vh(k)是點(diǎn)pk的加權(quán)V區(qū)域，pt是插入點(diǎn)，Lkt是pk和pt間的加權(quán)Voronoi邊，利用AE接口ITopologicalOperator.Cut(IPolylinecutter,refIGeometryleftGeom,refIGeometryrightGeom)方法，可將區(qū)域Vh(k)分割為leftGeom和rightGeom兩部分，其中，cutter即為L(zhǎng)kt，rightGeom和leftGeom對(duì)應(yīng)Vh(k) (分割后屬于pk的加權(quán)V區(qū)域)和Vh(tk) (從Vh(k)分割出來屬于pt的區(qū)域)。在AE提供的拓?fù)浞指罱涌谥校瑓?shù)Cutter即為由ILine和IConstructCircularArc接口生成的直線及弧段。該方法的具體實(shí)現(xiàn)包含在DivideRegionWithPB()和DivideRegionWithArc()中。

圖2 區(qū)域分割

Fig.2Regiondivisionbyanarc

加權(quán)V區(qū)域分割時(shí)，由于考慮權(quán)重因素，需判斷權(quán)重相等及不相等兩種情況，從而確定直線分割區(qū)域和弧段分割區(qū)域兩種方法。由于V區(qū)域劃分不涉及權(quán)重因素，所以只需考慮直線分割區(qū)域的情況。

(4)區(qū)域合并。UnionRegion()方法指定V區(qū)域和加權(quán)V區(qū)域合并方式。如圖3，設(shè)pk是h級(jí)加權(quán)V區(qū)域Vh(k) 所包含的點(diǎn)，k=1,2,…,c,c=VorPntsCount，若新插入一點(diǎn)pt(t=c+1)，其合并操作步驟如下：計(jì)算獲得pt與pk的加權(quán)Voronoi邊{L1t,L2t,…,Lkt,…,Lct}，利用Lkt對(duì)應(yīng)分割{Vh(1),Vh(2),…,Vh(c)}，將分割后屬于pt的區(qū)域{Vh(t1),Vh(t2),…,Vh(tc)}合并，生成pt的加權(quán)V區(qū)域Vh(t)。

圖3 區(qū)域合并

Fig.3 Regions union

2.3.3 等級(jí)加權(quán)V圖輸出 ExportLevelVorRegion()方法指定等級(jí)加權(quán)V圖的輸出方式：對(duì)于每一層級(jí)點(diǎn)集pLh，輸出對(duì)應(yīng)的加權(quán)V圖Vh，并以矢量圖形式存儲(chǔ)。將輸出圖層疊加，即可顯示最終的等級(jí)加權(quán)V圖效果。

3 應(yīng)用實(shí)例

圖4、圖5為等級(jí)V圖及加權(quán)V圖示例，其中母點(diǎn)集合為{pL1,pL2,pL3,pL4}，分別對(duì)應(yīng)于圖例中大圓、正方形、三角形、小圓所代表的點(diǎn)集，其中pL1包含4個(gè)母點(diǎn)，pL2包含10個(gè)母點(diǎn)，pL3包含11個(gè)母點(diǎn)，pL4包含125個(gè)母點(diǎn)。結(jié)果為四級(jí)等級(jí)V圖及加權(quán)圖{V1,V2,V3,V4}，邊界粗細(xì)程度代表等級(jí)程度差異，邊界越粗，等級(jí)越高。生成圖4所示的等級(jí)V圖平均用時(shí)5.21 s，生成圖5所示的等級(jí)加權(quán)V圖平均用時(shí)33.08 s。造成時(shí)間差異的主要原因是生成加權(quán)V圖時(shí)拓?fù)浞指詈秃喜⒉僮饕壬蒝圖復(fù)雜，所以用時(shí)較長(zhǎng)。

4 結(jié)論

圖4 等級(jí)V圖 圖5 加權(quán)等級(jí)V圖

Fig.4 HVD Fig.5 HMWVD

本文在文獻(xiàn)[10]的基礎(chǔ)上首次實(shí)現(xiàn)了等級(jí)V圖及加權(quán)V圖的ArcGIS矢量生成算法，該算法有如下特點(diǎn)：1)提供了一種適用于地學(xué)領(lǐng)域的空間數(shù)據(jù)模型，應(yīng)用該算法可以方便地建立等級(jí)V圖及加權(quán)V圖，用以表達(dá)城鎮(zhèn)等級(jí)體系、機(jī)構(gòu)布局設(shè)置、區(qū)域?qū)哟尾町惖葍?nèi)容，豐富了Voronoi圖在地學(xué)領(lǐng)域的應(yīng)用。2)利用ArcGIS提供接口，以拓?fù)浞指詈喜⒎绞綄?shí)現(xiàn)算法，精度高，實(shí)用性強(qiáng)。生成的等級(jí)V圖及加權(quán)V圖和原始生成元文件可以組合存儲(chǔ)于GIS空間數(shù)據(jù)庫中，便于管理，同時(shí)可以靈活地更改點(diǎn)、面樣式，豐富輸出結(jié)果的表達(dá)方式。等級(jí)V圖及加權(quán)V圖以矢量方式生成，以方便與其他矢量或柵格數(shù)據(jù)進(jìn)行疊加分析，挖掘更深層次的隱含信息。

[1] FARIN G,HAGEN H,HAMANN B.Hierarchical and Geometrical Methods in Scientific Visualization[M].Springer,2001.

[2] GOLD C,ANGEL P.Voronoi hierarchies[A].Proceedings of Geographic Information Science-4th International Conference,Geoscience[C].2006.99-111.

[3] WALLGRüN J O.Hierarchical Voronoi Graphs-Spatial Representation and Reasoning for Mobile Robots[D].Germany:University of Bremen,2008.

[4] 趙威，張仲元，秦耀辰,等.基于Voronoi圖的開封市域城鎮(zhèn)影像范圍和等級(jí)體系研究[J].河南大學(xué)學(xué)報(bào)(自然科學(xué)版)，2008，38(8)：271-275.

[5] 陳笑筑，王博，胡偉.基于Voronoi圖的鄉(xiāng)鎮(zhèn)等級(jí)體系劃分及影響范圍研究——以畢節(jié)市為例[J].安徽農(nóng)業(yè)科學(xué)，2011，39(35)：21947-21949.

[6] 薛瑩，俞路，韓貴峰.基于Voronoi圖剖分的區(qū)域旅游中心城市體系——以蘇浙滬地區(qū)為例[J].地域研究與開發(fā)，2008，27(4)：81-84.

[7] BALZER M,DEUSSEN O.Voronoi treemaps[A].Proceedings of the 2005 IEEE Symposium on Information Visualization[C].IEEE Press,2005.

[8] SUD A,FISHER D,LEE H P.Fast dynamic voronoi treemaps[A].Proceedings of the 2010 International Symposium on Voronoi Diagrams in Science and Engineering (ISVD′10)[C].2010.85-94.

[9] KANG S,LI J T,ZHANG F Q,et al.Hierarchical weighted Voronoi tessellations on visual application of spatial data[J].Journal of Computational Information Systems,2013,9(18):7341-7348.

[10] 范熙偉.加權(quán)Voronoi圖矢量生成算法研究及其實(shí)現(xiàn)[D].西安：西北大學(xué)，2011.21-35.

[11] 張有會(huì).加權(quán)Voronoi圖的畫法研究[J].計(jì)算機(jī)科學(xué)，2001，28(6)：126-130.

Vector-Based Realization of Hierarchical Voronoi Diagram and Multiplicatively Weighted Voronoi Diagram with ArcGIS

TIAN Song1,CUI Xi-min1,SUN Yun-hua1,GONG Yu2

(1.CollegeofGeoscienceandSurveyingEngineering,ChinaUniversityofMiningandTechnology,Beijing100083; 2.SouthChinaSeaBranch,StateOceanicAdministration,Guangzhou510310,China)

Hierarchical Voronoi diagram (HVD) and hierarchical multiplicatively weighted Voronoi diagram (HMWVD) are methods of expressing hierarchical relationships in the form of trees,which are widely applied in the field of geography,such as urban hierarchies exhibition,organization structures expression,etc.A vector-based approach is first proposed to generate HVD and HMWVD through the methods of regions division and regions union with ArcGIS Engine based on the thought of incremental method.It provides an important way to research and develop GIS spatial data model.

hierarchical Voronoi diagram;hierarchical multiplicatively weighted Voronoi diagram;ArcGIS;region division;regions union

2014-05-13；

2014-07-08

國(guó)家自然科學(xué)基金項(xiàng)目(71150001)

田松(1985-)，男，博士研究生，研究方向?yàn)镚IS算法研究。E-mail：tsisatc@163.com

10.3969/j.issn.1672-0504.2015.02.008

P208

A

1672-0504(2015)02-0034-04