1854 Cholera Outbreak

A quick overview of Fortune’s algorithm for computing Voronoi diagrams in \(O(n \log n)\) time for a input set of \(n\) points.

Fortune’s algorithm processes the points sequentially by first sorting the points by \(y\)-coordinate (\(O(n \log n)\) time). We then maintain a sweep line which processes events that modify the current Voronoi diagram at key locations of the sweep line. A separate beach line lags behind the sweep line and consists of a number of parabolic arcs. Voronoi edges will be generated at the intersection of adjacent parabolic arcs as the sweep line moves across the events.

We can maintain the current status of the beach line using a standard balanced binary search tree and maintain the event queue using a priority queue. A full analysis will show that there are never more then \(O(n)\) arcs in the tree or \(O(n)\) events in the queue. Since each event update, of which there are \(O(n)\), can be performed in \(O(\log n)\) time, we can get the overall \(O(n \log n)\) algorithm time bound.

One definition of a parabola is the set of points which are equidistant from a focus point and a directrix line. In our algorithm, we will associate a parabola with each input point above the sweep line and the directrix will be the current sweep line position. If two parabolas intersect on the beachline, that implies that intersection is equidistant from the two corresponding input points.

When the sweep line encounters a new point, we insert a new parabolic arc into the beach line.

When three parabolas intersect at a common point, this is typically called a circle or vertex event as it describes the point at which three input points lie on a circumcircle with the center being a vertex in the Voronoi diagram. These events always happen between three adjacent parabolic arcs in the beachline and we can detect when they will occur when we first insert a new point. These circle events can be added to the event queue for later processing, and when we encounter them later, we can update the beachline and output some Voronoi edges along with the newly created Voronoi vertex.

Today we’ll consider another computational geometry problem, perhaps a bit more removed from graphics, though we’ll make a brief connection back to Voronoi Diagrams. We will consider the problem of finding the closest pair of 2D points given a set of \(n\) points.

A naive upper bound would be to check all pairs of points and compute the minimum distance. How long would this take?

What if we had the Voronoi diagram for the set of points? A Voronoi Diagram is a planar graph with \(n\) faces, and \(O(n)\) vertices and edges, so its total size is linear. We can compute the Voronoi diagram in \(O(n \log n)\) time.

We’ll derive an alternative algorithm for computing the closest pair of points that does not require computing the Voronoi diagram. Consider a divide and conquer approach similar to that of merge sort. We will first split the set of points into a left half and right half and recursively answer the closest pair query on each half. When we get down to two or three points in a sub problem, we can compute the exact answer and return the result. Splitting the points in half can be done easily in \(O(n)\) time if we presort the points by \(x\) coordinate one time at the start of the algorithm.

The tricky part of this problem is merging the results of the left half and the right half. At first it might seem we just return the closest pair of the two halves. If the minimum distance between the pairs in the left and right halves are denoted \(d_l\) and \(d_r\), isn’t the minimum distance of the result \(d=\min(d_l, d_r)\)?

No, the closest pair could span the dividing line. But if such a pair \(p_i, p_j\) exists, both \(p_i\) and \(p_j\) must be within a distance \(d\) of the dividing line. Let \(S_y\) be a list of points within this slab of width \(2d\) of the dividing line, sorted by decreasing \(y\) coordinate. We can sort all the points once by \(y\) at the start of the algorithm and filter at each recursive merge step.

Checking all pairs in \(S_y\) could be too expensive, but since we have an estimate \(d\) on the minimum distance, we can do better.