一種多發生元Voronoi圖的柵格生成方法

江 威,吳 艷 蘭,譚 樹 東,馬 藝 文

(1.武漢市測繪研究院，湖北 武漢 430022；2.安徽大學資源與環境工程學院，安徽 合肥 230601；3.國家海洋信息中心，天津300171；4.武漢市國土資源和規劃信息中心，湖北 武漢 430014)

?

一種多發生元Voronoi圖的柵格生成方法

江 威1,吳 艷 蘭2*,譚 樹 東3,馬 藝 文4

(1.武漢市測繪研究院，湖北 武漢 430022；2.安徽大學資源與環境工程學院，安徽 合肥 230601；3.國家海洋信息中心，天津300171；4.武漢市國土資源和規劃信息中心，湖北 武漢 430014)

利用距離變換和柵格疊加分析，提出一種實現任意距離定義的2-site Voronoi圖生成方法。首先進行距離變換得到距離圖，然后通過鄰近關系對邊界進行劃分得到2-site Voronoi圖，最后將生成的距離圖和2-site Voronoi圖疊加。實驗表明，該文提出的2-site Voronoi圖生成方法可以快速構建多種距離類型和不同鄰近關系下的Voronoi圖，共生成了21種距離函數下的最遠、最鄰近和次鄰近Voronoi圖，解決了Voronoi圖的多樣性問題。該方法并不局限于點狀發生元，可以生成任意形態發生元Voronoi圖，并可以擴展生成N-site Voronoi圖，生成的廣義距離圖可用來模擬成組的發生元在諸多約束條件下的區域增長過程。

距離函數；Voronoi圖；距離變換；生長模型

0 引言

Voronoi圖(簡稱V圖)因其在空間鄰接和空間分割中的獨特性質，被廣泛應用于計算機輔助設計、機器人路徑規劃、模式識別等領域，也是計算幾何的一個重要研究內容。影響V圖外觀和屬性的要素有發生元、距離類型和空間類型,其中，距離類型由空間類型決定。常規V圖是由平面上n個發生元通過其歐氏距離最近空間點把平面劃分成n個獨立區域，其發生元是一組點集，每個Voronoi區域只包含一個點發生元，但這種V圖并不適用于發生元具有預先分組信息的場合。本文引用文獻[1]中的術語將常規V圖稱為1-site V圖；類似地，將發生元具有預先分組信息的V圖稱為N-site V圖，其主要外觀特征是每個Voronoi區域對應著一組點發生元而不是某個發生元。與單發生元V圖相比，多發生元V圖除在發生元的形態和數量上不同之外，其距離定義也不同。單發生元V圖中的距離值是每個像元到其最近單發生元的距離，而多發生元V圖的距離含義由距離函數決定，如距離和函數對應于空間每個像元到多個發生元的距離和，此時的距離是由多個發生元共同決定的。

目前對于單發生元V圖已有廣泛研究[2,3],對于多發生元V圖的研究則較少，并且主要研究了2-site V圖，即成對發生元V圖。點發生元成組之后會產生許多距離類型，因此與常規V圖相比，2-site V圖有更多的變形體。目前對于2-site V圖的研究，主要是從理論上探討不同距離類型對應的2-site V圖的復雜度。Barequet定義的2-site的距離函數包括距離和、距離差、距離乘積、三角形面積等，并利用矢量生成方法中的分治法有效生成了對應的V圖。但對于更復雜的距離函數，如這3個點組成的三角形的周長，其2-site V圖有待進一步研究。對此，Hanniel證明了距離函數為三角形周長的2-site V圖的計算復雜度為O(n2+ε)(ε>0)[4]。Dickerson等提出一個統一的模型形象地描述各種2-site V圖，如距離和函數和三角形周長對應的V圖，但沒有介紹構建每個V圖的方法[5]。最近，人們提出一些新的2-site距離函數的變形。例如，文獻[6，7]提出“視角”函數，文獻[7]還提出一些基于圓的距離函數(如3個點外接圓的半徑)。Vyatkina等進一步研究了平面上曼哈頓距離和切比雪夫距離下的2-site V的復雜度，并可以擴展到閔可夫斯基距離度量和更高維的空間[8]。Dickerson等利用2-site V圖，根據鄰近關系，給地理網絡圖中的點加了標注[9]。

2-site V圖具有很強的趣味性和特殊性。在生成方法上，Barequet等利用矢量方法中的分治法生成2-site V圖，但存在算法復雜等矢量方法固有的缺點，此外，其距離函數主要是算術距離和幾何距離(即三角形面積和三角形周長)，不具有靈活性。對于任意距離函數的2-site V圖，這種基于矢量的生成方法并不實用，目前缺乏適用于任意距離函數的2-site V圖的生成方法。本文利用距離變換和柵格疊加分析，提出一種實現任意距離定義的2-site V圖生成方法，并將其擴展為適用于任意形態發生元、任意預先分組的多發生元V圖。

1 方法原理

生成V圖，最重要的就是進行距離變換，也就是模擬發生元的區域增長過程。類似的，生成2-site V圖，首先也是進行距離變換得到距離圖，由于此時控制距離尺度的距離函數有多種，就會得到多種擴展的距離圖；然后，通過鄰近關系對邊界進行劃分，從而得到2-site V圖；最后將生成的距離圖和2-site V圖疊加，就得到了最終用來模擬區域增長過程的V圖。本文利用距離變換和柵格疊加，提出一種2-site V圖的生成方法，關鍵步驟為：1)由每個發生元分別生成距離圖，稱為單發生元距離圖(記為1-site距離圖)；2)對每對發生元，將其對應的兩個1-site 距離圖通過柵格疊加生成新的距離圖，稱為單對發生元距離圖(記為1-pair-site距離圖)，代表該距離圖中只有一對發生元；3)對于包含多對的一組發生元，再次通過柵格疊加運算，將其對應的單對距離圖合成為多對發生元距離圖(記為N-pair-site距離圖)；4)由N-pair-site 距離圖得到對應的2-site V圖。

值得說明的是：步驟2中柵格疊加對應的疊加算子為距離函數，步驟3中柵格疊加對應的疊加算子與步驟2不同，可定義為最小值(或最大值)函數，由此得到多對發生元的最近(或最遠)距離圖。為便于描述，將步驟2和步驟3中的柵格疊加分別稱為第一次柵格疊加和第二次柵格疊加。本文利用ArcObjects開發了一個構建2-site V圖和N-site V圖的實驗程序(圖1),可快速實現歐氏距離變換、柵格疊加(兩次)和V圖可視化。

圖1 2-site V圖生成流程

Fig.1 Flow chart of generation of 2-site V diagram

1.1 任意距離類型的單對發生元距離圖

單個發生元在平面上的增長過程可用圍繞發生元的同心等距線的1-site距離圖直觀表示。對于2-site V圖，由于區域增長過程中有兩個中心，因而其等距線圖的形狀更復雜。Barequet等基于矢量方法，利用膨脹的同心橢圓和同心帶等形狀的等距線，生成了前文提到的幾種距離函數對應的V圖[2],但在其后來的文章里提到的更加復雜的距離函數，其等距線圖很難用矢量的方法生成[5]。通過對1-site距離圖的分析，并結合距離變換和柵格疊加的特點，可在圖上直觀地表示平面上的一對發生元的區域增長過程。與常規距離圖相比，1-pair-site距離圖是由一對發生元生成，因此每個像元值為像元到這對發生元的距離，而不是像元到一個發生元的距離。通過第一次柵格疊加可生成1-pair-site距離圖。

如圖2所示，p和q為空間上任意兩點，其分別生成了1-site距離圖。因此任意空間點v到點p和q的距離可通過查詢這兩個1-site距離圖中的像元值獲取。在此基礎上，用不同的距離函數作為疊加算子，將兩個1-site距離圖中的像元值疊加，可生成一個新的距離圖，即1-pair-site距離圖。由于疊加算子由距離函數控制，疊加算子不同會直接影響1-pair-site距離圖的外觀和性質。因此，構建不同距離函數下的2-site V圖，在第一次柵格疊加時，要使用對應的距離函數作為疊加算子。例如，加、減、乘算子可用來生成距離函數分別為和、差、乘積的1-pair-site距離圖，并最終生成2-site V圖。圖2a為1-site 距離圖的生成說明，疊加算子使用的是3×3窗口，并用不同距離函數指定的疊加算子生成各類1-pair-site距離圖；圖2b為一對點發生元p和q利用不同的疊加算子生成1-pair-site距離圖的實例。

圖2 1-pair-site距離圖的生成

Fig.2 Generation of 1-pair-site DM

1.2 任意鄰近關系的多對發生元距離圖

與1-pair-site距離圖類似，N-pair-site距離圖是一種新的距離圖，用來模擬平面上N個點對的區域增長過程。本文中，N-pair-site距離圖的生成過程包含兩個關鍵步驟：1)計算每個像元到所有成對發生元的距離值；2)給每個像元分配最小(或最大)距離值，可以生成N-pair-site最近(或最遠)距離圖。

步驟1可通過查詢1-pair-site距離圖獲得,通過查詢圖中的像元值，可獲得每個像元到這對發生元的距離。因此，從N組發生元生成的一組1-pair-site距離圖，在每個1-pair-site距離圖中，通過查詢在相同位置的每個像元的值，可獲得每個像元分別到這N組發生元的距離。

步驟2的詳細說明見圖3a。假設有4組發生元，ID分別為p1、p2、p3、p4，利用距離乘函數分別生成這4組發生元的1-pair-site距離圖，并用ID值代表生成的4個1-pair-site距離圖。對于每個像元，都可以查詢到其在這4個1-pair-site距離圖中的像元值，即為每個像素到這些在特殊的距離函數下生成的1-pair-site距離圖的距離值。例如，這4個1-pair-site距離圖對應的3×3像素窗口的左上角的值分別為0、2、2、2。執行了第二次疊加之后，每個1-pair-site距離圖中的像素值通過一個特定的算子合并成一個值，并生成了4-pair-site距離圖和2-site V圖。同樣以3×3窗口的左上角像素值為例，在使用最小算子的情況下，4-pair-site距離圖的左上角像素值為0。圖3b為生成4-pair-site距離圖的實例，兩次疊加所用的疊加算子分別為距離乘函數和最近的鄰近關系。步驟2通過圖4所示的第二次疊加完成，其與第一次疊加類似，但是添加的數據和疊加算子不同。第二次柵格疊加，將一組1-pair-site距離圖作為添加的數據，在合成的新距離圖中，1-pair-site距離圖中的每個像元會生成一個新的像素值。這個合成圖即為N-pair-site距離圖。此處使用的疊加算子根據鄰近關系的不同而不同。例如，最小(或最大)算子用來尋找最近(或最遠)領域，并生成包含任意鄰近關系的N-pair-site距離圖。

圖3 N-pair-site距離圖和2-site V圖的生成

Fig.3 Generation ofN-pair-site DM and 2-site V diagram

1.3 成對發生元V圖的生成和可視化

生成距離圖后，對應的V圖可通過連接距離圖中的最大距離值或提取等距線相遇處的像素生成。使用生成4-pair-site距離圖時的相關信息，可生成2-site V圖。如圖3a所示，N-pair-site距離圖中的每個像素的距離值均來源于某個1-pair-site距離圖,這些像素值的來源決定了V圖的空間像素分配。因此，在生成4-pair-site距離圖的過程中，直接把像素值來源的1-pair-site距離圖的ID值賦給該像素，即可同時生成2-site V 圖?？紤]3×3窗口左上角的像素值，由于4-pair-site中該像素的值0來源的1-pair-site距離圖的ID值為p1，則2-site V圖中該像素的值為p1。實際上，N-pair-site距離圖中的距離值來源決定了Voronoi區域每個像元的值。改變疊加算子，如使用最遠和次鄰近的鄰域關系，可生成對應的4-pair-site距離圖和2-site V圖(圖3b)。

V圖可將平面劃分為若干區域或多邊形，將每個點對的ID值對應的顏色渲染給對應的區域，可將V圖可視化。渲染之后的Voronoi區域圖中，能夠迅速查詢出每個區域所歸屬的點對和鄰近關系，然而發生元的傳播方式和每個發生元到各個區域的距離等信息無法從V圖中獲取。在一幅圖中展現V圖的發生元、傳播距離和Voronoi區域這3個要素，就不會破壞V圖的完整性。尤其是2-site V圖，相比1-site V圖，其與距離定義和鄰近關系間的關系更密切。對于柵格疊加，關鍵步驟是把發生元和V圖疊加在距離圖之上，并將V圖調整至一定的透明度。

2 實驗結果

為驗證該方法的靈活性，本文使用文獻[2，5]中討論的算術距離函數和幾何距離函數及新定義的一些距離函數，共21種(表1)。利用這些距離函數，執行第一次柵格疊加，可生成對應的1-pair-site距離圖(圖4)。由表1和圖4可見，21種距離函數對應的距離圖存在形態相似情況，因此可將表1中的21種距離函數和對應的圖6中的1-pair-site距離圖分為幾組，每組距離函數本質上相同。如圖4(1)和圖4(10)，其分別對應距離函數(f1)和2-site的三角形周長公式(f10)，因這兩個函數的差為一常數，所以圖4(1)和圖4(10)的等距線形狀相同。類似的還有3組：(f2，f20)，(f7，f9)，(f11，f14，f15，f16)。對這4組距離函數，每組只選擇一個函數進行后續實驗,則后續實驗所用距離函數為16種。

表1 21種距離函數的定義

Table 1 The definition of twenty-one kinds of distance functions

距離函數說明f1d(v，p)+d(v，q)距離和f2|d(v，p)-d(v，q)|距離差的絕對值f3d(v，p)×d(v，q)距離乘積f4Max(d(v，p)，d(v，q))Min(d(v．p)，d(v，q))距離相除f5Max(d(v，p)，d(v，q))取較大距離f6Min(d(v，p)，d(v，q))取較小距離f7d(v，lpq)空間點v到直線lpq的距離f8d(v，pq)空間點v到線段pq的距離f9SΔvpqΔvpq的面積f10PcΔvpq帶參數的Δvpq周長：PcΔvpq=d(v，p)+d(v，q)+c·d(p，q)且c≥-1。當c=1，PcΔvpq為Δvpq的周長f11R1c(v，p，q)Δvpq的外接圓半徑距離函數說明f12R2cv，p，qΔvpq的內切圓半徑f13R3Ovpq包含v、p、q的最小圓的半徑f14d(o，pq)o到線段pq的距離，o為Δvpq內切圓圓心f15SΔopqΔopq的面積，o為Δvpq內切圓圓心f16pΔopqΔopq的周長，o為Δvpq內切圓圓心f17∠pvq夾角∠pvqf18w1·d(v，p)+w2·d(v，q)加權距離，w1和w2分別為d(v，p)和d(v，q)的權重f19d2(v，p)+d2(v，q)距離平方和f20Var(d(v，(p，q)))V到點p和q距離的方差：Var(d(v，p)，d(v，q))=[(d(v，p)-d)2+(d(v，q)-d)2]/2且d=(d(v，p)+d(v，q))/2f21Var(d(v，p，q))v，p，q3個點之間的3個距離的方差：Var(d(v，p，q))=[(d(v，p)-d)2+(d(v，q)-d)2+(d(p，q)-d)2]/3且d=(d(v，p)+d(v，q)+d(p，q))/3

圖4 21種距離函數下的1-pair-site距離圖

Fig.4 1-pair-site DMs under 21 kinds of distance functions

3 討論

對本文方法的特性作如下討論：1)由圖4和圖5可見，本文方法可以生成多種類型的2-site V圖，從而實現2-site V圖的多樣性。通過改變距離函數和鄰近關系，可以控制生成的2-site V圖的外觀和性質。本方法的主要特點是兩次柵格疊加過程：第一次柵格疊加的疊加算子是定義的21種距離函數，第二次柵格疊加的疊加算子為鄰近關系。通過改變第一次疊加的算子，不同的距離函數會生成不同結構的2-site距離圖(圖4)。類似地，改變第二次疊加的算子，可以生成不同鄰近關系的N-pair-site距離圖和2-site V圖(圖5)。由于對多種距離函數和鄰近關系的一般控制，本文的方法具有靈活性和可控性。2)本文方法并不局限于點狀發生元，其可以用于生成任意形態發生元的V圖。這種基于柵格的處理方式，由于任何發生元(點、弧段和面)都按照相同的方法計算距離，因而極大地簡化了不同形態發生元2-site V圖的生成過程。3)本文方法可以生成N-site V圖，此時的距離函數定義為每個點到多個發生元的距離。例如，3-site V圖的一個距離和函數定義的距離為點到3個發生元的距離。只需改變第一次疊加，該方法可生成3-site V圖。4)本文方法可以利用生成的距離圖模擬多種距離函數和鄰近關系下的成組發生元的區域增長過程(圖5),且本文方法簡單、健全且容易實現，其主要步驟是距離變換和柵格疊加,計算復雜度小，穩定性高。5)作為一種基于柵格的方法，當發生元數量比較多時，該方法比較費時。為了提高效率，可以結合并行計算等方法，從而節省運算時間。

4 結語

本文提出一種基于距離變換和柵格疊加的多發生元生成方法，生成了21種距離函數下的最遠、次鄰近和最鄰近3種鄰近關系下的2-site V圖，并將其擴展生成為多發生元V圖。實驗表明：與現有矢量方法相比，本文方法具有普遍性和靈活性，具體體現在：1)可根據任意距離函數和鄰近關系，生成多種2-site V圖；2)方法并不局限于點發生元，可以應用于任意形態發生元(如弧段和面)；3)由于采用柵格運算方式，發生元由點類型擴展到任意形態類型，并沒有增加算法的復雜度；4)可以在成對發生元V圖2-site V圖的基礎上擴展生成N-site V圖；5)除V圖外，本文方法生成的廣義距離圖，可用來模擬成組的發生元在諸多約束條件下的區域增長過程。因此，它很可能成為分析生長模型在相關應用領域(如材料科學和植物生態學)的一個有力工具。

[1] BAREQUET G,DICKERSON M T,SCOT DRYSDALE R L.2-Point site Voronoi diagrams[J].Discrete Applied Mathematics,2002,122(1):37-54.

[2] CHEN Z G,XIAO Y Y,CAO J.Approximation by piecewise polynomials on Voronoi tessellation[J].Graphical Models,2014,76(5):522-531.

[3] EMIRIS I Z,MANTZAFLARIS A,MOURRAIN B.Voronoi diagrams of algebraic distance fields[J].Computer-Aided Design,2013,45(2)：511-516.

[4] HANNIEL I,BAREQUET G.On the triangle-perimeter two-site Voronoi diagram[A].6th International Symposium on Voronoi Diagrams in Science and Engineering[C].2009.

[5] DICKERSON M T,EPPSTEIN D.Animating a continuous family of two-site Voronoi diagrams (and a proof of a bound on the number of regions)[A].Proceedings of the 25th Annual Symposium on Computational Geometry[C].2009.

[6] ASANO T,TAMAKI H.Angular Voronoi diagram with applications[A].3rd International Symposium on Voronoi Diagrams in Science and Engineering[C].2006.

[7] BAREQUET G,DICKERSON M T,EPPSTEIN D,et al.On 2-site Voronoi diagrams under geometric distance functions[A].8th International Symposium on Voronoi Diagrams in Science and Engineering[C].2011.

[8] VYATKINA K,BAREQUET G.On 2-site Voronoi diagrams under arithmetic combinations of point-to-point distances[A].7th International Symposium on Voronoi Diagrams in Science and Engineering[C].2010.

[9] DICKERSON M T,GOODRICH M T.Two-site Voronoi diagrams in geographic networks[A].16th ACM SIGSPATIAL International Conference on Advances in Geopraphic Information Systems[C].2008.

[10] VACAVANT A.Fast distance transformation on irregular two-dimensional grids[J].Pattern Recognition,2010,43(10)：3348-3358.

[11] CHEN J.A raster-based method for computing Voronoi diagrams of spatial objects using dynamic distance transformation[J].International Journal of Geographical Information Science,1999,13(3)：209-225.

A Raster-Based Method to Generate Two-Site Voronoi Diagrams

JIANG Wei1,WU Yan-lan2,TAN Shu-dong3,MA Yi-wen4

(1.WuhanGeomaticInstitute,Wuhan430022;2.SchoolofResoucesandEnvironmentEngineering,AnhuiUniversity,Hefei230601;3.NationalMarineDataandInformationService,Tianjin300171;4.WuhanLandResourcesandPlanningBureau,Wuhan430014,China)

In contrast with the regular (1-site) Voronoi Diagram,the 2-site VD has more variants due to the change of the distance function.Generating the 2-site VD with respect to arbitrary distance functions and constraints remains a challenge in the field.This paper proposes a flexible and general approach to generate 2-site VD by combining the distance transform and raster overlay into a unified framework.This framework is characterized by the two procedures of raster overlay：the first overlay procedure whose operator is used to control the distance function,and the second overlay procedure whose operator is specified by the neighbor relationship.By manipulating the two overlay operators,the nearest-(furthest-,etc.) neighbor 2-site VDs with respect to the various distance functions can be obtained.The proposed approach was implemented and tested with 21 kinds of distance functions.The results show improved flexibility and robustness over existing vector-based approaches and emphasize the convenience of extending to general sites andN-site VD.The proposed approach also producesN-pair-site distance map：a new type of distance map,providing an easy and convenient way to simulate the region-growing process of the sites previously grouped into different pairs.

distance function;Voronoi diagram;distance transform;growth model

2014-12-22；

2015-03-11

國家自然科學基金項目(41271443);安徽省自然科學基金項目(1308085MD52)

江威(1988-)，男，碩士研究生，研究方向為地圖學與地理信息系統。*通訊作者E-mail：wylmq@sina.com

10.3969/j.issn.1672-0504.2015.05.009

P208

A

1672-0504(2015)05-0039-05