__Voronoi Diagrams__

__Voronoi Diagrams__

Recall the *post-oﬃce problem *mentioned in the introduction. Here we want to preprocess a set *S *of *n *points, referred to as *sites*, in the plane such that we can answer the query: which site in *S *is closest to a given query point *q*? In order to solve this problem we can divide the plane into regions according to the nearest-neighbor rule: each site *s *gets assigned the region which is closest to *s*. This subdivision, which compactly captures the distance information inherent in a given conﬁguration, is called the *Voronoi diagram *of *S*—see Fig. 62.7(a).

More formally, the Voronoi diagram of a set of sites *S *= *{**s*1*,…**, s**n**} *in R*d*, which we refer to as Vor(*S*), partitions space into *n *regions—one for each site—such that the region for a site *s**i *consists of all points that are closer to *s**i *than to any other site *s**j **∈ **S*. The set of points that are closest to a particular site *s**i *forms the so-called *Voronoi cell *of *s**i*, and is denoted by *V *(*s**i*). Thus, when *S *is a set of sites in the plane we have

where dist(*.**, .*) denotes the Euclidean distance.

Now consider the dual graph of the Voronoi diagram, that is, the graph that has a node for every Voronoi cell and an arc between any two Voronoi cells that share a common edge—see Fig. 62.7(b). (Observe that the concept of dual graph used here has nothing to do with the duality transform discussed in Section 62.2.3.) Suppose we embed this graph in the plane, by using the site *s**i *to represent the node corresponding to the cell *V *(*s**i*) and by drawing the edges as straight line segments, as in Fig. 62.7(c). Somewhat surprisingly perhaps, this graph is always planar. Moreover, it is actually a triangulation of the point set *S*, assuming that no four points in *S *are co-circular. More details on this special triangulation, which is called the *Delaunay triangulation*, will be given in Section 62.5.1.

There exists a fascinating connection between Voronoi diagrams in R*d *and half-space intersections in R*d*+1. Assume for simplicity that *d *= 2, and consider the transformation that maps a site *s *= (*s**x**, s**y *) in R2 to the non-vertical plane *h*(*s*) : *z *= 2*s**x**x*+2*s**y **y **−*(*s*2 +*s*2 ) in . Geometrically, *h*(*s*) is the plane tangent to the unit paraboloid *z *= *x *+*y *at the point vertically above (*s**x**, s**y **, *0). Let *H*(*S*) be the set of planes that are the image of a set of point sites *S *in the plane. Let *S *denote the convex polyhedron that is formed by the intersection of the positive half-spaces deﬁned by the planes in ( ), that is, = n *h**∈**H*(*S*) *h*+, where *h*+ denotes the half-space above h. Surprisingly, the projection of the edges and vertices of *S *vertically downward on the *x**y*-plane is exactly the Voronoi diagram of *S*.

The Voronoi diagram can be deﬁned for various sets of sites, for example points, line segments, circles, or circular arcs, and for various metrics. Sections 62.4.1 and 62.4.2 discuss the complexity and construction algorithms of Voronoi diagrams for the usual case of point sites in R*d*, while Section 62.4.3 describes some of the possible variations. Additional details and proofs of the material presented in this section can be found in [4, 5, 22].

Let *S *be a set of *n *points in R*d*. The Voronoi diagram of *S *is a cell complex in R*d*. If *d *= 2 then the Voronoi cell of a site is the interior of a convex, possibly inﬁnite polygon. Its boundary consists of *Voronoi edges*, which are equidistant from two sites, and *Voronoi **ver**t**i**c**e**s*, which are equidistant from at least three sites. The Voronoi diagram of *n **≥ *3 sites has at most 2*n **− *5 vertices and at most 3*n **− *6 edges, which implies that the Delaunay triangulation has at most 2*n **− *5 triangles and 3*n **− *6 edges. What is the complexity of the Voronoi diagram in *d **≥ *3? Here the connection between Voronoi diagrams in R*d *and intersections of half-spaces in R*d*+1 comes in handy: we know from the Upper Bound Theorem that the intersection of *n *half-spaces in R*d*+1 has complexity *O*(*n**l*(*d*+1)*/*2*J *) = *O*(*n**l**d/*2*1*). This bound is tight, so the maximum complexity of the Voronoi diagram—and of the Delaunay triangulation, for that matter—in R*d *is Θ(*n**l**d/*2*1*).

**Construction**

In this section we present several algorithms to compute the Voronoi diagram or its dual, the Delaunay triangulation, for point sites in R*d*. Several of these algorithms can be generalized to work with metrics other than the Euclidean, and to sites other than points. Since the Voronoi diagram in the plane has only linear complexity one might be tempted to search for a linear time construction algorithm. However the problem of sorting *n *real numbers is reducible to the problem of computing Voronoi diagrams and, hence, any algorithm for computing the Voronoi diagram must take Ω(*n *log *n*) time in the worst case.

Two data structures that are particularly well suited to working with planar subdivi- sions like the Voronoi diagram are the *doubly-connected edge list *(DCEL) by Muller and Preparata [38] and the *quad-edge structure *by Guibas and Stolﬁ [29]. Both structures re- quire *O*(*n*) space, where *n *is the complexity of the subdivision, and allow to eﬃciently traverse the edges adjacent to a vertex and the edges bounding a face. In both structures, one can easily obtain in *O*(*n*) time a structure for the Voronoi diagram from a structure for the Delaunay triangulation, and vice versa. In fact, the quad-edge structure, as well as a variant of the DCEL [7], are representations of a planar subdivision and its dual at the same time, so there is nothing to convert (except, perhaps, that one may wish to explicitly compute the coordinates of the vertices in the Voronoi diagram, if the quad-edge structure or DCEL describes the Delaunay triangulation).

**D****ivide-and-conquer**. The ﬁrst deterministic worst-case optimal algorithm to compute the Voronoi diagram was given by Shamos and Hoey [52]. Their divide-and-conquer algorithm splits the set *S *of point sites by a dividing line into subsets *L *and *R *of approximately the same size. The Voronoi diagrams Vor(*L*) and Vor(*R*) are computed recursively and then merged into Vor(*S*) in linear time. This algorithm constructs the Voronoi diagram of a set of *n *points in the plane in *O*(*n *log *n*) time and linear space.

**P****l****a****ne Sweep**. The strategy of a sweep line algorithm is to move a line, called the *sweep** **li**ne*, from top to bottom over the plane. While the sweep is performed, information is maintained regarding the structure one wants to compute. It is tempting to apply the same approach to Voronoi diagrams, by keeping track of the Voronoi edges that are currently intersected by the sweep line. It is problematic, however, to discover a new Voronoi region in time: when the sweep line reaches a new site, then the line has actually been intersecting the Voronoi edges of its region for a while. Fortune [23] was the ﬁrst to ﬁnd a way around this problem. *Fortune’s algorithm *applies the plane sweep paradigm in a slightly diﬀerent fashion: instead of maintaining the intersection of the sweep line with the Voronoi diagram, it maintains information of the part of the Voronoi diagram of the sites above the line that can not be changed by sites below it. This algorithm provides an alternative way of computing the Voronoi diagram of *n *points in the plane in *O*(*n *log *n*) time and linear space.

**Randomized incremental construction**. A natural idea is to construct the Voronoi diagram by *incremental insertion*, that is, to obtain Vor(*S*) from Vor(*S**{**s**}*) by inserting the site *s*. Insertion of a site means integrating its Voronoi region into the diagram constructed so far. Unfortunately the region of *s *can have up to *n **− *1 edges, for *|**S**| *= *n*, which may lead to a running time of *O*(*n*2). The insertion process is probably better described and implemented in the dual environment, for the Delaunay triangulation *D**T *: construct

*D**T**i *= *D**T *(*{**s*1*,…**, s**i**}*) by inserting *s**i *into *D**T**i**−*1. The advantage of this approach over a direct construction of Vor(*S*) is that Voronoi vertices that appear only in intermediate diagrams but not in the ﬁnal one need not be computed or stored. *D**T**i *is constructed by exchanging edges, using edge ﬂips [36], until all edges invalidated by *s**i *have been removed.

Still, the worst-case running time of this algorithm can be quadratic. However, if we insert the sites in *random order*, and the algorithm is implemented carefully, then one can prove that the expected running time is *O*(*n *log *n*), and that the expected amount of storage is *O*(*n*).

**Other approaches. **Finally, recall the connection between Delaunay triangulations and convex hulls. Since there exist *O*(*n *log *n*) algorithms to compute the convex hull of points in R3 (see Section 62.3) we therefore have yet another optimal algorithm for the computation of Voronoi diagrams.

**Variations**

In this section we present some of the common variations on Voronoi diagrams. The ﬁrst is the *order-**k **V**o**r**o**noi diagram *of a set *S *of *n *sites, which partitions R*d *on the basis of the ﬁrst *k *closest point sites. In other words, each cell in the order-*k *Voronoi diagram of a set *S *of sites in the plane corresponds to a *k*-tuple of sites from *S *and consists of those points in the plane for which that *k*-tuple are the *k *closest sites. One might fear that the order-*k** *Voronoi diagram has Θ(*n**k*) cells, but this is not the case. In two dimensions, for example, its complexity is *O*(*k*(*n**−**k*)), and it can be computed in *O*(*k*(*n**−**k*) log *n*+*n *log3 *n*) expected time [2].

The *furthest-site Voronoi diagram *partitions R*d *according to the furthest site, or equivalently, according to the closest *n **− *1 of *n *sites. The furthest-site Voronoi diagram can be computed in *O*(*n *log *n*) time in two dimensions, and in *O*(*n**l**d/*2*1 *) in dimension *d **≥ *3.

One can also consider diﬀerent distance functions than the Euclidean distance. For ex- ample, one can alter the distance function by the addition of *additive *or *multiplicative** **weigh**t**s*. In this case every point site *s**i *is associated with a weight *w**i *and the distance function *d*(*s**i**, x*) between a point *x *and a site *s**i *becomes *d*(*s**i**, x*) = *w**i *+ dist(*s**i**, x*) (additive weights) or *d*(*s**i**, x*) = *w**i **· *dist(*s**i**, x*) where dist(*s**i**, x*) denotes the Euclidean distance between *s**i *and *x*. The Voronoi diagram for point sites in 2 dimensions with additive weights can be computed in *O*(*n *log *n*) time, for multiplicative weights the time increases to *O*(*n*2) time.

Finally the *power diagram*, or *Laguerre diagram*, is another Voronoi diagram for point sites *s**i *that are associated with weights *w**i*. Here the distance function is the *power distance *intro- duced by Aurenhammer [3], where the distance from a point *x *to a site *s**i *i__s measured alon g__a line tangent to the sphere of radius *√**w**i *centered at *s**i*, i.e., *d*(*s**i**, x*) = jdist(*s**i**, **x*)2 *− **w**i*.

The power diagram can be computed in *O*(*n *log *n*) time in two dimensions.