% Chapter 9. Voronoi Diagrams
\documentclass[12pt,a4paper]{amsart}
%\usepackage[utf8]{inputenc}
%\usepackage[czech]{babel}

\textwidth=154mm
\hoffset=-14mm

\usepackage{amssymb}
\usepackage{amsmath}
\usepackage{amsthm}
\usepackage{amsfonts}
\usepackage{latexsym}
\usepackage{esvect}

\theoremstyle{plain}
\newtheorem{thm}{Theorem}[section]
\newtheorem{cor}[thm]{Corollary}
\newtheorem{lemma}[thm]{Lemma}
\newtheorem{prop}[thm]{Proposition}
\newtheorem{ap}[thm]{Appendix}

\theoremstyle{definition}
\newtheorem*{definition}{Definition}
\newtheorem{remark}[thm]{Remark}  
\newtheorem{example}[thm]{Example}
\newtheorem*{ack}{Acknowledgement}

\def\R{\mathbb{R}}
\def\d{\operatorname{dist}}

\begin{document}

\setcounter{section}{8}
\section{\large Voronoi Diagrams}

\subsection*{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=\{p_1,p_2,\dots,p_n\}$ be a set of points in the plane. The {\it Voronoi diagram} 
(briefly V-diagram) is a plane subdivision which has faces
$$V(p_i)=\{q\in \R^2;\ \d(q,p_i)\leq \d(q,p_j) \text{ for all }j\in P\setminus\{i\} \},$$
where $\d(q,p)$ denotes the classical Euclidean distance of points $q$ and $p$.

\medskip
FIGURE 9.1  An example of the Voronoi diagram. 
\medskip

The set of points in the plane having the distance to a point $ p_i $ less than or equal to the distance to the point $p_j \ne p_i$ is a halfplane
$$h(p_i,p_j)=\{q\in \R^2;\ \d(q,p_i)\leq \d(	q,p_j)\}.$$
The face $V(p_i)$, called the Voronoi cell, is the intersection of all such halfplanes for
$j\ne i$
$$V(p_i)=\bigcap_{j\ne i}h(p_i,p_j).$$
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(n^2\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:

\begin{thm} Any Voronoi diagram for the set of $n\geq 3$ points (which do not lie on one line) has at most $2n-5$ vertices and at most $3n-6$ edges.
\end{thm}

\begin{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_{\infty}$  "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. 

\medskip
FIGURE 9.2 The planar graph made of the  Voronoi diagram.
\medskip
 
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\geq 3(m+1).$$
Substituting into this inequality  $h=m+n-1$ from the Euler Formula we get
$$2m+2n-2\geq 3(m+1),$$
which gives $m\leq 2n-5$. Similarly after substituting $m=h-n+1$ we obtain
$$h\leq 3n-6.$$
\end{proof}

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$ 
$$ C_P (q) = \{r\in \mathbb R^2; \d(q, r)\leq \d(q, r)\}.$$ 
The circle bounding this disc contains at least one point of the set $P$.

Obviously we have:

\begin{lemma} 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 $C_P(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 $C_P(r)$. 
\end{lemma}

\medskip
FIGURE 9.3 Illustration to Lemma 9.2.
\medskip

\subsection*{Sweep line method}
One of the possible algorithms for the construction of the V-diagram 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\in l^+$ is lying in this area if its distance from a point $p_i\in P\cap l^+$ is less or equal to its distance from the line $l$
\begin{equation}
\d(q,p_i)\leq \d(q,l).
\end{equation}
Then for all points $p_j\in P$ from the halfplane $l^-$ under the line $l$ 
$$\d(q,p_j)\geq\d(q,l)\geq\d(q,p_i).$$
So the point $q$ from the area determined by the enequality (1) does not lie
in the Voronoi cell $V(p_j)$ of the point $p_j\in P\cap l^-$.

Denote $\alpha(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 $\alpha^+(p,l)$ is the set of points lying above or on the parabola $\alpha(p,l)$. The area above the sweep line $l$ where the V-diagram will be formed, is the union 
$$\bigcup_{p_i\in P\cap l^+}\alpha^+(p_i,l).$$
The boundary of this area, formed by arcs of parabolas, is called 
the  {\it beach line}.

In this case the tree $\mathcal 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 $p_1$ and $p_2$ met, have the same distance from $p_1$ and 
$p_2$, as
$$\d(q,p_1)=\d(q,l)=\d(q,p_2),$$
and so they lie on edges of the Voronoi diagram. These breakpoints are internal nodes of the tree $\mathcal T$. Their description of the form $\langle p_1\,:\,p_2\rangle$ 
means that the node corresponds to the intersection of the parabolas with
the focuses $p_1$ and $p_2$ (both have the sweep line  $l$ as the directrix) where the arc of the parabola with the focus $p_1$ is on the left and the arc of the parabola with focus $p_2$ 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.

\medskip
FIGURE 9.4 A beach line. 
\medskip

FIGURE 9.5 A balanced binary tree corresponding the the previous beach line.
\medskip

\subsection*{Events}
The answers to the following two questions are of fundamental importance:
\begin{itemize}
\item In which position of the sweep line does a new arc appear on the beach line?
\item In which position of the sweep line does an arc disappear from the beach line?
\end{itemize}
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.

\begin{lemma}
New arcs in the beach line arise just as the sweep line passes through some point of the set $P$.
\end{lemma}
Such a point is caled {\it site event}.  The following animation illustrates this situation. The sweep line $l$ goes through the point $p\in P$. Let $\alpha$ be the arc of the beach line above the point $p$  with the focus at the point $p_1\in P$. After the sweep line passes through the point $p$, the arc $\alpha$ splits into two arcs $\alpha_1$ and $\alpha_2$ and a new arc $\beta$, a part of the parabola with focus $p$, appears between them in the beach line.

\medskip
FIGURE 9.6 A site event $p$.
\medskip

The intersections of the arcs $\alpha_1$ and $\beta$ and of the arcs $\beta$ and 
$\alpha_2$ are points of  the edge of the Voronoi diagram which separates the 
Voronoi cell $V(p_1)$ from $V(p)$. The change in the beach line corresponds to the following change of the tree $\mathcal T$:

\medskip
FIGURE 9.7  The change of $\mathcal T$ corresponding to the previous site event. 
\medskip

The internal nodes $\langle p_1\, :\, p_2\rangle$,   $\langle p_1\,:\,p\rangle$ a  $\langle p\,:\,p_1\rangle$ of the tree  $\mathcal 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 $\mathcal T$ we have to rebalance it.

Now consider three consecutive arcs $\alpha$, $\beta$ and $\gamma$ of the beach line. Let $p_1$, $p_2$, and $p_3$ be gradually the focuses of their respective parabolas. Consider the circle curcumscribed to the tringle $p_1p_2p_3$. 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 {\it circle event} for the arc $\beta$.

When the sweep line $l$ passes this point, the center $s$ of the curcumcircle of the triangle $p_1p_2p_3$ has the same distance from these points as from the sweep line 
$l$, so it lies on all three arcs $\alpha$, $\beta$ and $\gamma$. When moving the sweep line $l$ under $q$, the arc $\beta$ 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(p_1)$ from the cell $V(p_2)$,  $V(p_2)$ from $V(p_3)$ and $V(p_1)$ from 
$V(p_3)$ meet.

\medskip
FIGURE 9.8 The sweep line passes a circle event $q$
\medskip

FIGURE 9.9  The change of $\mathcal T$ corresponding to the previous circle event. 
\medskip

\begin{lemma}
An existing arc of the parabola disappears from the beach line if and only if the sweep line passes through its circle event.
\end{lemma}

\subsection*{Algorithm}
The start of the algorithm corresponds to the position of the sweep line $l$ above all points. The tree $\mathcal T$ is empty, and the event queue $\mathcal 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:
\begin{itemize}
\item We add an arc to the tree $\mathcal T$ and  edges into the double-connected edge list for the V-diagram,
\item We remove an arc from the tree $\mathcal T$ and we add an edge and a vertex  to the V-diagram,
\item From the queue $\mathcal 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.
\end{itemize}  
We do this procedure until the queue $\mathcal Q$ is empty. However, the tree 
$\mathcal 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 $\mathcal 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.

\medskip
ANIMATION
\begin{enumerate}
\item The sweep line is above all the points. The beach line is empty. The queue is $\mathcal Q=(p_1,p_2,p_3)$.
\item The beasch line is a single parabola. The queue is $\mathcal Q=(p_2,p_3)$.
\item The sweep line passes through the site event $p_2$. A new arc is arising in the beach line.
\item After the sweep line has passed $p_2$, the beach line is formed by three arcs
$\alpha_1, \beta. \alpha_2$. The green segment is the edge of the V-diagram. The queue is $\mathcal Q=(p_3)$.
\item The sweep line passes through the site event $p_3$. There is a new arc $\gamma$  in the beach line, which divides the arc $\alpha_2$  in two parts
$\alpha_2$ and $\alpha_3$. The circle event $q$ for the arc $\alpha_2$ is created
The queue  is $\mathcal Q=(q)$.
\item After passing the point $p_3$ the beach line consists of the arcs $\alpha_1, \beta,\alpha_2,\gamma,\alpha_3$. The break points $\langle p_1\,:\,p_3\rangle$ and $\langle p_3\,:\,p_1\rangle$ determine the  edge of the V-diagram between the cells 
$V(p_1)$ and $V(p_3)$. The queue is $\mathcal Q=(q)$.
\item The sweep line passes through the  circle event $q$
The arc $\alpha_2$ disappears from the beach line. The edges of the V-diagram between the cells $V(p_2)$ and $V(p_1)$ and between $V(p_1)$ and $V(p_3)$ 
meet at the center $s$ of the circle curcumscribed to the triangle
$p_1p_2p_3$.
\item After passing the point $q$ the queue is empty because a potential circle event for arcs $\beta$ and $\gamma$ (it would be the point $q$)  lies above the sweep line $l$. 
The break point $\langle p_2\,:\, p_3\rangle$ determines a further edge of the V-diagram. 
The beach line is formed by the arcs $\alpha_1,\beta, \gamma, \alpha_3$. The algorithm ends by creating a rectangle $R$ containing all points $p_1,p_2,p_3$ and all break points. The resulting edges of the V-diagram are "touched" to the rectangle boundary.
\end{enumerate}
\medskip

The exact description of the algorithm is captured in the following pseudocodes.
\medskip

PSEUDOCODE VoronoiDiagram($P$) from pseudo.pdf, page 30.
\medskip

PSEUDOCODE HandleSiteEvent($p_i$) from pseudo.pdf, page 31
\medskip

 PSEUDOCODE HandleCircleEvent($\gamma$) from pseudo.pdf, page 32
\medskip

\subsection{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:

\begin{thm} 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. 
\end{thm} 

An implementation of the algorithm using Matlab has been done by Jana 
Tom\v sov\' a in the bachelor thesis (in Czech) available at

\centerline{https://is.muni.cz/auth/th/175345/prif$\_$b/}

\end{document} 





