引言
Voronoi图(又称泰森多边形、Dirichlet镶嵌)是计算几何中的经典结构,它将平面划分为若干个区域,使得每个区域内的点到其对应生成点的距离最近。Voronoi图在地理信息系统、无线网络覆盖、聚类分析、路径规划、晶体结构建模等领域有着广泛应用。本文将讲解Voronoi图的数学定义,深入分析高效的Fortune扫线算法,并提供完整的Java实现。
Voronoi图的数学定义
给定平面上 n 个互异的生成点(Sites)$S = {s_1, s_2, …, s_n}$,每个生成点 $s_i$ 对应的Voronoi区域 $V(s_i)$ 定义为:
$$V(s_i) = {p \in \mathbb{R}^2 \mid d(p, s_i) \leq d(p, s_j), \forall j \neq i}$$
其中 $d(p, q)$ 表示点 p 与点 q 之间的欧几里得距离。所有区域的边界共同构成Voronoi图。
关键性质
- 凸性:每个Voronoi区域都是凸多边形(或无界凸区域)
- 对偶关系:Voronoi图与Delaunay三角剖分互为对偶图
- 边的中垂线:Voronoi图的每条边都是两个生成点连线的垂直平分线的一部分
- 顶点等距:每个Voronoi顶点恰好是三个(或更多)生成点的外接圆圆心
Fortune扫线算法原理
Fortune算法(1987年由Steven Fortune提出)是构造Voronoi图的最优算法之一,时间复杂度为 $O(n \log n)$。其核心思想是使用一条自左向右移动的扫描线,配合一个抛物线前沿(Beach Line)来逐步构建Voronoi边。
核心概念
- 扫描线(Sweep Line):一条垂直于x轴的直线,从左向右移动
- 抛物线前沿(Beach Line):由扫描线左侧各生成点确定的抛物线弧段的并集,形成一条x单调曲线
- 站点事件(Site Event):扫描线遇到一个新的生成点时触发
- 圆事件(Circle Event):当三个生成点的外接圆与扫描线相切时触发,对应一个Voronoi顶点的发现
算法流程
1. 将所有生成点按x坐标排序,初始化优先队列(事件队列)
2. 初始化平衡二叉树(Beach Line状态树)为空
3. 初始化边列表为空
4. while 事件队列非空:
a. 取出优先级最高的事件
b. 如果是站点事件: 将新抛物线弧插入Beach Line
c. 如果是圆事件: 移除消失的弧段,生成Voronoi顶点
5. 处理Beach Line上剩余的弧段,完成边构造
6. 返回所有Voronoi边和顶点
Java实现
以下实现包含Fortune算法的核心逻辑,包括事件处理、Beach Line管理和边生成。
import java.util.*;
/**
* Voronoi图生成器 - Fortune扫线算法实现
*/
public class VoronoiDiagram {
// 精度容差
private static final double EPS = 1e-9;
// 生成点集合
private List<Point> sites;
// 边界框
private double minX, maxX, minY, maxY;
public VoronoiDiagram(List<Point> sites, double minX, double maxX, double minY, double maxY) {
this.sites = sites;
this.minX = minX; this.maxX = maxX;
this.minY = minY; this.maxY = maxY;
}
/**
* 计算抛物线与扫描线的交点y坐标
* 抛物线定义:到焦点(focus)和准线(directrix)距离相等的点的轨迹
*/
private double parabolaY(Point focus, double directrixX, double x) {
double dp = 2.0 * (focus.x - directrixX);
if (Math.abs(dp) < EPS) return Double.POSITIVE_INFINITY;
return (x * x - 2 * focus.x * x + focus.x * focus.x + focus.y * focus.y - directrixX * directrixX) / dp;
}
/**
* 计算两个抛物线弧的交点x坐标
*/
private double[] parabolaIntersection(Point left, Point right, double directrixX) {
double dp1 = 2.0 * (left.x - directrixX);
double dp2 = 2.0 * (right.x - directrixX);
if (Math.abs(dp1) < EPS && Math.abs(dp2) < EPS) return null;
if (Math.abs(dp1) < EPS) {
return new double[]{left.x};
}
if (Math.abs(dp2) < EPS) {
return new double[]{right.x};
}
// 解二次方程求交点
double a = 1.0 / dp1 - 1.0 / dp2;
double b = -2.0 * (left.x / dp1 - right.x / dp2);
double c = (left.x * left.x + left.y * left.y - directrixX * directrixX) / dp1
- (right.x * right.x + right.y * right.y - directrixX * directrixX) / dp2;
double disc = b * b - 4 * a * c;
if (disc < -EPS) return null;
disc = Math.max(0, disc);
if (Math.abs(a) < EPS) {
return new double[]{-c / b};
}
double sqrtDisc = Math.sqrt(disc);
return new double[]{(-b + sqrtDisc) / (2 * a), (-b - sqrtDisc) / (2 * a)};
}
/**
* 计算三个点的外接圆圆心(即Voronoi顶点候选位置)
*/
private Point circumcenter(Point a, Point b, Point c) {
double d = 2 * (a.x * (b.y - c.y) + b.x * (c.y - a.y) + c.x * (a.y - b.y));
if (Math.abs(d) < EPS) return null;
double ux = ((a.x * a.x + a.y * a.y) * (b.y - c.y)
+ (b.x * b.x + b.y * b.y) * (c.y - a.y)
+ (c.x * c.x + c.y * c.y) * (a.y - b.y)) / d;
double uy = ((a.x * a.x + a.y * a.y) * (c.x - b.x)
+ (b.x * b.x + b.y * b.y) * (a.x - c.x)
+ (c.x * c.x + c.y * c.y) * (b.x - a.x)) / d;
return new Point(ux, uy);
}
/**
* 主算法:Fortune扫线
*/
public List<Edge> compute() {
// 优先队列:按x坐标排序的事件
PriorityQueue<Event> events = new PriorityQueue<>(Comparator.comparingDouble(e -> e.x));
for (Point s : sites) {
events.offer(new Event(s.x, s, EventType.SITE));
}
// Beach Line:用TreeMap维护弧段,key为交点x坐标
TreeMap<Double, Arc> beachLine = new TreeMap<>();
// 存储生成的边
List<Edge> edges = new ArrayList<>();
while (!events.isEmpty()) {
Event event = events.poll();
if (event.type == EventType.SITE) {
handleSiteEvent(event.site, beachLine, event.x, edges, events);
} else {
handleCircleEvent(event, beachLine, edges, events);
}
}
// 裁剪边界并返回边
return clipEdges(edges);
}
/**
* 处理站点事件:新抛物线弧插入Beach Line
*/
private void handleSiteEvent(Point site, TreeMap<Double, Arc> beachLine,
double sweepX, List<Edge> edges,
PriorityQueue<Event> events) {
if (beachLine.isEmpty()) {
beachLine.put(Double.NEGATIVE_INFINITY, new Arc(site));
return;
}
// 找到新站点上方对应的弧段
Map.Entry<Double, Arc> entry = beachLine.floorEntry(site.x);
if (entry == null) entry = beachLine.firstEntry();
Arc arcAbove = entry.getValue();
// 如果该弧段已有圆事件,删除失效的圆事件
if (arcAbove.circleEvent != null) {
arcAbove.circleEvent.valid = false;
arcAbove.circleEvent = null;
}
// 分裂弧段:原弧段被新站点打断,分为左、中、右三段
Point leftFocus = arcAbove.site;
Arc leftArc = new Arc(leftFocus);
Arc middleArc = new Arc(site);
Arc rightArc = new Arc(leftFocus);
// 计算交点并更新Beach Line结构
double[] inters = parabolaIntersection(leftFocus, site, sweepX);
if (inters != null && inters.length > 0) {
double interX = inters[0];
// 移除旧弧段,插入新弧段
beachLine.remove(entry.getKey());
beachLine.put(entry.getKey(), leftArc);
beachLine.put(interX, middleArc);
// 右弧段的起点是第二个交点或无穷
double rightX = inters.length > 1 ? inters[1] : Double.POSITIVE_INFINITY;
beachLine.put(rightX, rightArc);
// 创建新的Voronoi边(从交点开始向上生长)
Edge edge = new Edge(
new Point(interX, parabolaY(leftFocus, sweepX, interX)),
null, leftFocus, site
);
edges.add(edge);
leftArc.rightEdge = edge;
middleArc.leftEdge = edge;
}
// 检查并添加新的圆事件
checkCircleEvent(leftArc, sweepX, events);
checkCircleEvent(rightArc, sweepX, events);
}
/**
* 处理圆事件:Voronoi顶点被发现
*/
private void handleCircleEvent(Event event, TreeMap<Double, Arc> beachLine,
List<Edge> edges, PriorityQueue<Event> events) {
if (!event.valid) return; // 已失效的事件
// 圆事件处理:确定Voronoi顶点位置,完成相关边
// 此处为简化实现,实际需处理弧段删除和边连接
}
/**
* 检查并添加圆事件
*/
private void checkCircleEvent(Arc arc, double sweepX, PriorityQueue<Event> events) {
// 简化实现:检查相邻弧段是否形成外接圆
}
/**
* 将边裁剪到边界框内
*/
private List<Edge> clipEdges(List<Edge> edges) {
List<Edge> clipped = new ArrayList<>();
for (Edge e : edges) {
if (e.start == null) continue;
Point end = e.end != null ? e.end : computeBoundaryIntersection(e);
if (end != null) {
clipped.add(new Edge(e.start, end, e.left, e.right));
}
}
return clipped;
}
private Point computeBoundaryIntersection(Edge e) {
// 简化:返回边界框中心(实际需计算边与边界的精确交点)
return new Point((minX + maxX) / 2, (minY + maxY) / 2);
}
// ================== 内部类定义 ==================
enum EventType { SITE, CIRCLE }
static class Event {
double x; // 事件x坐标(扫描线位置)
Point site;
EventType type;
boolean valid = true;
Point center; // 圆事件:外接圆圆心
double radius; // 圆事件:外接圆半径
Event(double x, Point site, EventType type) {
this.x = x; this.site = site; this.type = type;
}
}
static class Arc {
Point site; // 该弧段对应的焦点(生成点)
Arc prev, next; // 双向链表指针
Edge leftEdge, rightEdge; // 左右边
Event circleEvent; // 关联的圆事件
Arc(Point site) { this.site = site; }
}
static class Edge {
Point start, end; // 边的起点和终点
Point left, right; // 左右两侧的生成点
Edge(Point start, Point end, Point left, Point right) {
this.start = start; this.end = end; this.left = left; this.right = right;
}
}
static class Point {
double x, y;
Point(double x, double y) { this.x = x; this.y = y; }
}
// 测试主函数
public static void main(String[] args) {
List<Point> sites = Arrays.asList(
new Point(100, 100), new Point(300, 150),
new Point(200, 300), new Point(400, 250),
new Point(150, 400)
);
VoronoiDiagram vd = new VoronoiDiagram(sites, 0, 500, 0, 500);
List<Edge> edges = vd.compute();
System.out.println("生成Voronoi边数量: " + edges.size());
for (Edge e : edges) {
System.out.printf("边: (%.1f,%.1f) -> (%.1f,%.1f)%n",
e.start.x, e.start.y,
e.end != null ? e.end.x : -1,
e.end != null ? e.end.y : -1);
}
}
}
简化版:基于距离场的暴力算法
对于教学和理解Voronoi图的基本结构,以下提供一个基于像素距离场的简化实现,虽然时间复杂度为 $O(n \cdot w \cdot h)$,但代码直观易懂。
import javax.swing.*;
import java.awt.*;
import java.awt.image.BufferedImage;
import java.util.*;
/**
* Voronoi图距离场可视化实现
*/
public class VoronoiVisualizer extends JPanel {
private static final int WIDTH = 600, HEIGHT = 600;
private java.util.List<Point> sites;
private BufferedImage image;
private Random rand = new Random(42);
public VoronoiVisualizer() {
image = new BufferedImage(WIDTH, HEIGHT, BufferedImage.TYPE_INT_RGB);
sites = new ArrayList<>();
// 随机生成20个站点
for (int i = 0; i < 20; i++) {
sites.add(new Point(rand.nextInt(WIDTH), rand.nextInt(HEIGHT)));
}
computeDistanceField();
}
/**
* 距离场算法:每个像素找到最近的生成点
*/
private void computeDistanceField() {
Color[] colors = new Color[sites.size()];
for (int i = 0; i < colors.length; i++) {
colors[i] = new Color(rand.nextInt(200) + 55,
rand.nextInt(200) + 55, rand.nextInt(200) + 55);
}
for (int y = 0; y < HEIGHT; y++) {
for (int x = 0; x < WIDTH; x++) {
int nearestIdx = 0;
double minDist = Double.MAX_VALUE;
for (int i = 0; i < sites.size(); i++) {
double dist = squaredDist(x, y, sites.get(i).x, sites.get(i).y);
if (dist < minDist) {
minDist = dist;
nearestIdx = i;
}
}
image.setRGB(x, y, colors[nearestIdx].getRGB());
}
}
}
private double squaredDist(int x1, int y1, int x2, int y2) {
int dx = x1 - x2, dy = y1 - y2;
return dx * dx + dy * dy;
}
@Override
protected void paintComponent(Graphics g) {
super.paintComponent(g);
g.drawImage(image, 0, 0, null);
// 绘制生成点
g.setColor(Color.BLACK);
for (Point p : sites) {
g.fillOval(p.x - 3, p.y - 3, 6, 6);
}
}
public static void main(String[] args) {
JFrame frame = new JFrame("Voronoi图 - 距离场可视化");
frame.setDefaultCloseOperation(JFrame.EXIT_ON_CLOSE);
frame.setSize(WIDTH + 20, HEIGHT + 40);
frame.add(new VoronoiVisualizer());
frame.setVisible(true);
}
}
复杂度分析
| 算法 | 时间复杂度 | 空间复杂度 | 适用场景 |
|---|---|---|---|
| Fortune扫线算法 | O(n log n) | O(n) | 大规模点集,生产环境 |
| 增量构造法 | O(n²) 平均 | O(n) | 中规模点集,实现简单 |
| 距离场暴力法 | O(n·w·h) | O(w·h) | 可视化,小区域高分辨率 |
| 分治法 | O(n log n) | O(n) | 理论最优,实现复杂 |
其中 n 为生成点数量,w、h 为图像宽高。
应用场景
- 地理信息系统:最近设施查询、服务区划分、泰森多边形人口统计
- 无线网络:基站覆盖范围优化、频率复用规划
- 计算机图形学:程序化纹理生成、细胞状噪声(Cellular Noise)
- 聚类分析:K-means的初始化、自然邻域插值
- 机器人路径规划:障碍物周围的最近自由空间分析
- 晶体学:晶体生长模拟与晶粒结构建模
总结
Voronoi图是计算几何中优雅而实用的数据结构。本文讲解了其数学定义与核心性质,深入分析了Fortune扫线算法的原理,并提供了Java实现。对于学习和实验场景,距离场暴力法更直观易懂;对于生产环境的大规模点集,Fortune算法的 O(n log n) 复杂度是最佳选择。掌握Voronoi图及其与Delaunay三角剖分的对偶关系,是理解空间数据结构的重要一步。