9. Voronoi Diagrams
Introduction. In this chapter we will deal with the problem known as the ”post office
problem”. It is about dividing the territory of the city into areas around individual
post offices so that in each area the local post office is the closest. Another example is
the determination of the catchment areas of hospitals by distance. The mathematical
formulation of the problem is as follows:
Let P = {p1, p2, . . . , pn} be a set of points in the plane. The Voronoi diagram
(briefly V-diagram) is a plane subdivision which has faces
V (pi) = {q ∈ R2
; dist(q, pi) ≤ dist(q, pj) for all j ∈ P \ {i}},
where dist(q, p) denotes the classical Euclidean distance of points q and p.
FIGURE 9.1 An example of the Voronoi diagram.
The set of points in the plane having the distance to a point pi less than or equal
to the distance to the point pj = pi is a halfplane
h(pi, pj) = {q ∈ R2
; dist(q, pi) ≤ dist(q, pj)}.
The face V (pi), called the Voronoi cell, is the intersection of all such halfplanes for
j = i
V (pi) =
j=i
h(pi, pj).
According to Chapter 5 we can find the intersection of n − 1 halfplanes in time
O(n log n). To find all Voronoi cells by this way would need time at least O(n2
log n).
We will show an algorithm for finding the Voronoia diagram with running time O(n log n).
The fact that such running time is real is supported by the following statement:
Theorem 9.1. Any Voronoi diagram for the set of n ≥ 3 points (which do not lie on
one line) has at most 2n − 5 vertices and at most 3n − 6 edges.
Proof. Suppose that given points are not colinear. Denote m and h the number of
vertices and edges of the V-diagram, respectively. If we add a vertex v∞ ”in infinity”
to which go all the edges of the V-diagram which are halflines, we create a planar
graph with n faces, m + 1 vertices and h edges. See the figure.
FIGURE 9.2 The planar graph made of the Voronoi diagram.
The graph satisfies the Euler Formula
m + 1 − h + n = 2.
The degree of every vertex is at least 3, the sum of all degrees is 2h. Hence
2h ≥ 3(m + 1).
Substituting into this inequality h = m + n − 1 from the Euler Formula we get
2m + 2n − 2 ≥ 3(m + 1),
which gives m ≤ 2n − 5. Similarly after substituting m = h − n + 1 we obtain
h ≤ 3n − 6.
1
2
For a finite set of points P and each point q in the plane not lying in P define the
disc with the center q and the radius equal to the distance of the point q from the set
P
CP (q) = {r ∈ R2
; dist(q, r) ≤ dist(q, r)}.
The circle bounding this disc contains at least one point of the set P.
Obviously we have:
Lemma 9.2. A point q lies on an edge of the V-diagram if and only if there are at
least two points of the set P lying on the boundary of CP (q).
A point r is a vertex of the V-diagram if and only if there are at least three points
of the set P lying on the boundary of CP (r).
FIGURE 9.3 Illustration to Lemma 9.2.
Sweep line method. One of the possible algorithms for the construction of the Vdiagram
uses the sweep line method, but in a way that is somewhat more sophisticated
than those with which we have met so far. The sweep line l moves from the top to
the bottom and creates the V-diagram above it. This diagram does not occur in the
entire halfplane l+
above the line l, but only in its part, which is the union of areas
above certain parabols. How does this area look like?
A point q ∈ l+
is lying in this area if its distance from a point pi ∈ P ∩ l+
is less or
equal to its distance from the line l
(1) dist(q, pi) ≤ dist(q, l).
Then for all points pj ∈ P from the halfplane l−
under the line l
dist(q, pj) ≥ dist(q, l) ≥ dist(q, pi).
So the point q from the area determined by the enequality (1) does not lie in the
Voronoi cell V (pj) of the point pj ∈ P ∩ l−
.
Denote α(p, l) the parabola formed by the points which have the same distance from
the point p as from the line l. The point p is the focus of the parabola and the line l
is called the directrix of the parabola. Let α+
(p, l) is the set of points lying above or
on the parabola α(p, l). The area above the sweep line l where the V-diagram will be
formed, is the union
pi∈P∩l+
α+
(pi, l).
The boundary of this area, formed by arcs of parabolas, is called the beach line.
In this case the tree T connected with this sweep line method will be a balanced
binary tree that maintains the order of arcs of the beach line in its leaves. Breakpoints
of the beach line, i.e., the points where two consecutive arcs of the beach line with the
focuses p1 and p2 met, have the same distance from p1 and p2, as
dist(q, p1) = dist(q, l) = dist(q, p2),
and so they lie on edges of the Voronoi diagram. These breakpoints are internal nodes
of the tree T . Their description of the form p1 : p2 means that the node corresponds
to the intersection of the parabolas with the focuses p1 and p2 (both have the sweep
3
line l as the directrix) where the arc of the parabola with the focus p1 is on the left
and the arc of the parabola with focus p2 is on the right. Using such a tree, for a point
p on the sweep line l, you can find a beach line arc that is above the point p. The
following figure shows an example of a beach line and a corresponding tree.
FIGURE 9.4 A beach line.
FIGURE 9.5 A balanced binary tree corresponding the the previous beach line.
Events. The answers to the following two questions are of fundamental importance:
• In which position of the sweep line does a new arc appear on the beach line?
• In which position of the sweep line does an arc disappear from the beach line?
The following two lemmas answer these questions. Their proofs are technical, so we
will not do them, and rather limit ourselves to illustrating them with animations.
Lemma 9.3. New arcs in the beach line arise just as the sweep line passes through
some point of the set P.
Such a point is caled site event. The following animation illustrates this situation.
The sweep line l goes through the point p ∈ P. Let α be the arc of the beach line
above the point p with the focus at the point p1 ∈ P. After the sweep line passes
through the point p, the arc α splits into two arcs α1 and α2 and a new arc β, a part
of the parabola with focus p, appears between them in the beach line.
FIGURE 9.6 A site event p.
The intersections of the arcs α1 and β and of the arcs β and α2 are points of the
edge of the Voronoi diagram which separates the Voronoi cell V (p1) from V (p). The
change in the beach line corresponds to the following change of the tree T :
FIGURE 9.7 The change of T corresponding to the previous site event.
The internal nodes p1 : p2 , p1 : p a p : p1 of the tree T describe breakpoints
of the beach line using the focuses of the corresponding parabolas in the order from left
to right, as we have explained above. After changing the tree T we have to rebalance
it.
Now consider three consecutive arcs α, β and γ of the beach line. Let p1, p2,
and p3 be gradually the focuses of their respective parabolas. Consider the circle
curcumscribed to the tringle p1p2p3. Let us call s its center and denote q the lowest
point of this circle, that is, the point with the minimal coordinate y. If the point q
lies under the sweep line l, we call it a circle event for the arc β.
When the sweep line l passes this point, the center s of the curcumcircle of the
triangle p1p2p3 has the same distance from these points as from the sweep line l, so it
lies on all three arcs α, β and γ. When moving the sweep line l under q, the arc β
disappears from the beach line. The center s of the circle is a vertex of the Voronoi
diagram in which the edges separating the cell V (p1) from the cell V (p2), V (p2) from
V (p3) and V (p1) from V (p3) meet.
FIGURE 9.8 The sweep line passes a circle event q
FIGURE 9.9 The change of T corresponding to the previous circle event.
4
Lemma 9.4. An existing arc of the parabola disappears from the beach line if and only
if the sweep line passes through its circle event.
Algorithm. The start of the algorithm corresponds to the position of the sweep line
l above all points. The tree T is empty, and the event queue Q contains all the points
of the set P (site events) arranged lexicographically from top to bottom and from left
to right. When the sweep line passes through an event, depending on whether it is a
site or circle event, we do the following actions:
• We add an arc to the tree T and edges into the double-connected edge list for
the V-diagram,
• We remove an arc from the tree T and we add an edge and a vertex to the
V-diagram,
• From the queue Q we remove the event that we have just passed, eventually the
circle events that corresponded to the triples of consecutive arcs from which one
has disappeared after passing through the event. The order of arcs changes,
and the current circle events are counted from three consecutive arcs in the
new order. New circle events are added to the queue.
We do this procedure until the queue Q is empty. However, the tree T will not be
empty. Its current internal nodes determine the edges of the V-diagram that are
halflines. This allows us to find a rectangle R, in which all the points of the set P, all
the vertices of the V-diagram and the break points corresponding to these remaining
internal nodes of the tree T are located.
The course of the algorithm in the simplest possible situation of a set of three points
is demonstrated by the following animation.
ANIMATION
(1) The sweep line is above all the points. The beach line is empty. The queue is
Q = (p1, p2, p3).
(2) The beasch line is a single parabola. The queue is Q = (p2, p3).
(3) The sweep line passes through the site event p2. A new arc is arising in the
beach line.
(4) After the sweep line has passed p2, the beach line is formed by three arcs
α1, β.α2. The green segment is the edge of the V-diagram. The queue is
Q = (p3).
(5) The sweep line passes through the site event p3. There is a new arc γ in the
beach line, which divides the arc α2 in two parts α2 and α3. The circle event
q for the arc α2 is created The queue is Q = (q).
(6) After passing the point p3 the beach line consists of the arcs α1, β, α2, γ, α3.
The break points p1 : p3 and p3 : p1 determine the edge of the V-diagram
between the cells V (p1) and V (p3). The queue is Q = (q).
(7) The sweep line passes through the circle event q The arc α2 disappears from
the beach line. The edges of the V-diagram between the cells V (p2) and V (p1)
and between V (p1) and V (p3) meet at the center s of the circle curcumscribed
to the triangle p1p2p3.
(8) After passing the point q the queue is empty because a potential circle event
for arcs β and γ (it would be the point q) lies above the sweep line l. The break
5
point p2 : p3 determines a further edge of the V-diagram. The beach line is
formed by the arcs α1, β, γ, α3. The algorithm ends by creating a rectangle R
containing all points p1, p2, p3 and all break points. The resulting edges of the
V-diagram are ”touched” to the rectangle boundary.
The exact description of the algorithm is captured in the following pseudocodes.
PSEUDOCODE VoronoiDiagram(P) from pseudo.pdf, page 30.
PSEUDOCODE HandleSiteEvent(pi) from pseudo.pdf, page 31
PSEUDOCODE HandleCircleEvent(γ) from pseudo.pdf, page 32
9.1. Time complexity. Let the set P for which we implement the algorithm have n
points. After the first site event in the queue, the beach line will be a single parabola.
When passing through other site events, the number of arcs in the beach line increases
at most by two. Therefore, the total number of arcs that sometimes occur in the beach
line will not exceed 2n − 1. In a binary balanced tree with the maximum of 2n − 1
leaves, the operations of addition of two leaves, removing one leaf and rebalancing the
tree take the time
O(log(2n − 1)) = O(log n).
Since each arc in the beach line arises once and disappears at most once, the number of
exetuted events will be O(2n−1) = O(n). In each event, we execute a finite number of
steps which last a constant time or time O(log n). This gives the proof of the following
statement:
Theorem 9.5. The Voronoi diagram for n points can be computed using the above
sweep line algoritm in O(n log n) time, using O(n) storage.
An implementation of the algorithm using Matlab has been done by Jana Tomˇsov´a
in the bachelor thesis (in Czech) available at
https://is.muni.cz/auth/th/175345/prif b/