「学习笔记」计算几何（一）——扫描线处理包含关系

西瓜nd · 2025-12-04 17:52:18 · 算法·理论

圆的扫描线算法

引入

平面上有 n 个圆，圆之间没有公共点，即圆之间要么包含，要么相离。此时圆要么都是一个包含一个，要么都是同级并列的，不难发现这很像树的结构，那我们不妨将圆的包含关系看做树的边然后建树，这样就会形成一个森林。我们可以给森林一个超级源点，将森林中所有树的根连起来，此时森林就变成了一棵树，相当于此时原点处有一个半径无限大的圆，所有圆都被包含其中。这样就能将平面上圆之间的问题转化为树上点之间的问题。

那么我们怎么找出圆之间的包含关系呢？此为问题的核心。

朴素的，我们可以枚举两个圆来判断，这样是 \mathcal{O}(n^2) 的，有没有更优的做法呢？

我们不妨先来思考一下如果是矩形该怎么处理。

如果按照扫描线正常的流程，遇到一个矩形的左边界，我们会把矩形上下边界围成的区间丢到线段树里维护，但这里我们希望要找出他们之间的包含关系，显然不能这么做，因为不能通过线段树找到一个区间的父区间。如果你仔细观察会发现，这些区间其实就像一个括号序列，想象一下矩形上边界是左括号，下边界是右括号，矩形之间的包含关系就像一个括号包含另一个括号，此时我们就相当于要维护一个括号序列并能够快速找到一个括号的父括号（直接包含此括号的括号）。此时我们把左右括号都丢到 set 里根据坐标维护相对位置，不难发现，插入一对括号时，我们想要找的父括号一定可以通过相邻括号（前驱或后继括号）找到，这样就能在 \mathcal{O}(n\log n) 的时间内对矩形建树。

了解了矩阵怎么处理，那么圆的做法就显而易见了。

算法过程

圆不同于矩形那样有规则的上下界，但我们可以类比上述的思想，实现圆的扫描线。

首先和传统扫描线的方法一样先把每个圆左右侧 x 坐标离散出来排个序标记好进出，然后开始扫描。在每个 x 处，扫描线从上到下会和一些圆相交，交点将扫描线划分为若干区间（如下图），因为圆之间没有公共点，那么从上到下扫描线穿过圆的顺序是不会变的，对应的这些区间的相对顺序也不会变，于是我们可以用一个 set （或其他平衡树，这里使用 set）来动态地维护这些区间。

具体的，我们将圆拆成上下半圆（相当于左右括号），set 中的区间需要维护的信息有上端点所交的半圆与下端点所交的半圆，以及区间所属的父圆。每次扫描到一个圆的左端点时，我们可以通过二分（lower_bound）在 set 中找到点所属的区间，通过该区间找到父圆并建边，然后将此区间分裂为两个区间并更新信息（如果你学过珂朵莉树的话，会发现这和珂朵莉树的 split 操作很像）。这样我们就能在 \mathcal{O}(n\log n) 的时间内对圆建树。

实现

以下代码是基于上述思想实现的，稍微进行了一些封装便于理解，如果想要代码简短一点其实定义一个结构体就够了，不需要将扫描线分成一个一个区间，直接把圆拆成上下半圆丢进 set 就行了，具体可以参考凸多边形扫描线的代码 （圆的懒得再写一遍了）。

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
const int maxn = 100050, mod = 998244353, INF = 1e9;
const double eps = 1e-9;
int sweep_line_x = -1e9;
vector<int> g[maxn];

struct Circle {
    int x, y, r;

    double get_y(bool up_side) const {
        int b = -2 * y;
        ll c = (ll)y * y + (sweep_line_x - x) * ll(sweep_line_x - x) - (ll)r * r;
        ll delta = (ll)b * b - 4 * c;
        return (-b + (up_side ? sqrt(delta) : -sqrt(delta))) / 2;
    }

    bool operator<(const Circle& rhs) const {
        return r < rhs.r;
    }
};

struct Curve {
    bool up_side;
    Circle c;

    bool operator<(const Curve& rhs) const {
        return c.get_y(up_side) + eps < rhs.c.get_y(rhs.up_side);
    }
};

struct Interval {
    int id;
    Curve up, down;

    bool operator<(const Interval& rhs) const {
        return down < rhs.down || up < rhs.up;
    }

    bool operator==(const Interval& rhs) const {
        return !(*this < rhs) && !(rhs < *this);
    }
};

int main() {
    ios::sync_with_stdio(0);
    int n;
    cin >> n;

    vector<Circle> c(n + 1);
    vector<pair<int, int>> events;
    for (int i = 1;i <= n;i++) {
        cin >> c[i].x >> c[i].y >> c[i].r;
        events.emplace_back(c[i].x - c[i].r, i);
        events.emplace_back(c[i].x + c[i].r, i);
    }

    c[0] = { 0,0,INF };
    sort(events.begin(), events.end());

    set<Interval> s;
    s.insert({ 0, {true, c[0]}, {false, c[0]} });

    for (auto &[x, id] : events) {
        sweep_line_x = x;
        Interval cur = { id, {true, c[id]}, {false, c[id]} };
        if (x == c[id].x - c[id].r) { //圆的左端点
            auto it = --s.lower_bound(cur);
            auto par = *it;
            s.erase(it);
            g[par.id].push_back(id);
            Interval up = par, down = par;
            up.down = { true, c[id] };
            down.up = { false, c[id] };
            s.insert(up);
            s.insert(down);
            s.insert(cur);
        }
        else {  //圆的右端点
            auto it = s.lower_bound(cur);
            auto joined = *prev(it);
            joined.up = next(it)->up;
            s.erase(prev(it));
            s.erase(next(it));
            s.erase(it);
            s.insert(joined);
        }
    }

    return 0;
}

例题

往日重现

题意：每次询问从 (x,y) 到 (p,q) 最少要穿过多少个圆的边界，n 个圆之间没有任何公共点，n\leq 10^5，值域 10^7，询问 10^5 次。

因为圆之间不相交，那么圆之间的包含关系可以看做树的边，穿过一个圆的边界就相当于走过树上的一条边，树上两点间距离即为答案。可以使用圆的扫描线算法建树，查询的时候通过 LCA 和树的深度就可求得两点所属圆在树上的最短距离。因为题目问的是两点的关系，所以要先求出每一个点属于哪一个圆（该圆是包含此点的圆中半径最小的），将点看做半径为 0 的圆一并扫描即可，需要先将询问点离线。

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
#define int ll
const int maxn = 100050, mod = 998244353, INF = 1e9;
const double eps = 1e-9;
int sweep_line_x = -1e9;
vector<int> g[maxn];
int dep[maxn], fa[maxn][22];

struct Circle {
    int x, y, r;

    double get_y(bool up_side) const {
        int b = -2 * y;
        ll c = (ll)y * y + (sweep_line_x - x) * ll(sweep_line_x - x) - (ll)r * r;
        ll delta = b * (ll)b - 4 * c;
        return (-b + (up_side ? sqrt(delta) : -sqrt(delta))) / 2;
    }

    bool operator<(const Circle& rhs) const {
        return r < rhs.r;
    }
};

struct Curve {
    bool up_side;
    Circle c;

    bool operator<(const Curve& rhs) const {
        return c.get_y(up_side) + eps < rhs.c.get_y(rhs.up_side);
    }
};

struct Interval {
    int id;
    Curve up, down;

    bool operator<(const Interval& rhs) const {
        return down < rhs.down || up < rhs.up;
    }

    bool operator==(const Interval& rhs) const {
        return !(*this < rhs) && !(rhs < *this);
    }
};

void dfs(int x, int fath) {
    fa[x][0] = fath;
    dep[x] = dep[fath] + 1;
    for (int i = 1;i <= 20;i++) {
        fa[x][i] = fa[fa[x][i - 1]][i - 1];
    }
    for (int to : g[x]) {
        if (to == fath)continue;
        dfs(to, x);
    }
}

int lca(int x, int y) {
    if (dep[x] < dep[y]) swap(x, y);

    for (int i = 20;i >= 0;i--) {
        if (dep[fa[x][i]] >= dep[y]) {
            x = fa[x][i];
        }
        if (x == y) return x;
    }
    for (int i = 20;i >= 0;i--) {
        if (fa[x][i] != fa[y][i]) {
            x = fa[x][i];
            y = fa[y][i];
        }
    }
    return fa[x][0];
}

struct Query {
    pair<int, int> a, b;
};

void solve() {
    int n, m;
    cin >> n >> m;

    vector<Circle> c(n + 1);
    vector<pair<int, int>> events, p;
    vector<Query> q(m + 1);
    for (int i = 1;i <= n;i++) {
        cin >> c[i].x >> c[i].y >> c[i].r ;
        events.emplace_back(c[i].x - c[i].r, i);
        events.emplace_back(c[i].x + c[i].r, i);
    }

    for (int i = 1;i <= m;i++) {
        pair<int, int> a, b;
        cin >> a.first >> a.second >> b.first >> b.second;
        q[i] = { a,b };
        p.push_back(a);
        p.push_back(b);
    }

    sort(p.begin(), p.end());
    p.erase(unique(p.begin(), p.end()), p.end());

    vector<int> bel(p.size() + 1, 0);

    for (int i = 0;i < p.size();i++) {
        events.emplace_back(p[i].first, i + n + 1);
    }

    c[0] = { 0,0,INF };
    sort(events.begin(), events.end());

    set<Interval> s;
    s.insert({ 0, {true, c[0]}, {false, c[0]} });

    for (auto &[x, id] : events) {
        sweep_line_x = x;
        if (id > n) {   //要询问的点
            id -= n + 1;
            Circle now = { p[id].first,p[id].second,0 };
            Interval cur = { id, {true, now}, {false, now} };
            auto it = --s.lower_bound(cur);
            auto par = *it;
            bel[id] = par.id;
            continue;
        }
        Interval cur = { id, {true, c[id]}, {false, c[id]} };
        if (x == c[id].x - c[id].r) { //圆的左端点
            auto it = --s.lower_bound(cur);
            auto par = *it;
            s.erase(it);
            g[par.id].push_back(id);
            Interval up = par, down = par;
            up.down = { true, c[id] };
            down.up = { false, c[id] };
            s.insert(up);
            s.insert(down);
            s.insert(cur);
        }
        else {  //圆的右端点
            auto it = s.lower_bound(cur);
            auto joined = *prev(it);
            joined.up = next(it)->up;
            s.erase(prev(it));
            s.erase(next(it));
            s.erase(it);
            s.insert(joined);
        }
    }

    dfs(0, 0);

    for (int i = 1;i <= m;i++) {
        pair<int,int> a=q[i].a,b=q[i].b;
        int x=lower_bound(p.begin(),p.end(),a)-p.begin();
        int y=lower_bound(p.begin(),p.end(),b)-p.begin();
        x=bel[x],y=bel[y];
        cout<<dep[x]+dep[y]-2*dep[lca(x,y)]<<'\n';
    }
}

signed main() {
    ios::sync_with_stdio(0);
    int t = 1;
    // cin>>t;
    while (t--) {
        solve();
    }

    return 0;
}

凸多边形的扫描线

回顾一下，由于圆没有规则的上下界，不能像矩形那样维护，但因为我们只是希望确定圆之间的包含关系，所以我们会想将圆拆成上下半圆，通过算半圆与扫描线的交点位置来确定相对位置。

我们发现，这一过程除了算交点外其实基本没用到圆本身的性质，那么如果不是圆能不能这么做呢？能不能推广一下呢？

不难发现，互不相交的凸多边形同样可做，因为我们依然能较快的求出与扫描线的交点，同时保持凸多边形的相对顺序不变，只是细节稍微多一点。

我们依然需要将凸多边形拆成上下两个部分，丢到 set 里面维护，不同的是，我们需要将每条边都加入扫描线，一条边一条边的扫，并在 set 里维护这些边的相对关系。

实现

实现方式相对于圆的进行了一些简化，直接把边丢到 set 里维护了。

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
#define int ll
const int maxn = 600050, mod = 998244353, INF = 1e9;
const double eps = 1e-9;
int sweep_line_x = -1e9;
vector<int> g[maxn];
int fa[maxn];

struct Segment {
    int id;
    pair<int, int> p, q;
    bool up_side;

    Segment() {}
    Segment(int id, pair<int, int> p, pair<int, int> q) : id(id), p(p), q(q)
    {
        if (p.first > q.first)up_side = true;
        else up_side = false;
    }
    double get_y(int x) const {
        if (p.first == q.first) return p.second;
        return p.second + 1.0 * (x - p.first) * (q.second - p.second) / (q.first - p.first);
    }

    bool operator<(const Segment& rhs) const {
        if (fabs(get_y(sweep_line_x) - rhs.get_y(sweep_line_x)) < eps) {
            if (up_side == rhs.up_side) return q.first < rhs.q.first;
            return up_side < rhs.up_side;
        }
        return get_y(sweep_line_x) + eps < rhs.get_y(sweep_line_x);
    }
};

struct Event {
    int x;
    Segment l;
    bool is_left;
    bool operator<(const Event& rhs) const {
        if (x == rhs.x) {
            if (is_left == rhs.is_left) return l.p.second > rhs.l.p.second;
            return is_left > rhs.is_left;
        }
        return x < rhs.x;
    }
};

signed main() {
    ios::sync_with_stdio(0);
    int n;
    cin >> n;

    vector<Segment> p;
    vector<Event> events;
    for (int i = 1;i <= n;i++) {
        int m;
        cin >> m;
        vector<pair<int, int>> tmp(m);
        for (int j = 0;j < m;j++) {
            cin >> tmp[j].first >> tmp[j].second;
        }
        for (int j = 0;j < m;j++) {
            int k = (j + 1) % m;
            if (tmp[j].first == tmp[k].first) continue;
            p.push_back({ i, tmp[j], tmp[k] });
        }
    }
    for (int i = 0; i < p.size(); i++) {
        if (p[i].p.first > p[i].q.first) swap(p[i].p, p[i].q);
        events.push_back({ p[i].p.first, p[i], true });
        events.push_back({ p[i].q.first, p[i], false });
    }

    sort(events.begin(), events.end());

    memset(fa,-1,sizeof fa);
    set<Segment> s;

    for (auto& [x, l, f] : events) {
        sweep_line_x = x;
        int id = l.id;
        if (l.up_side and fa[id] == -1) {
            auto it = s.upper_bound(l);
            if (it != s.end()) {
                Segment par = (*it);
                if (par.up_side) fa[id] = par.id;
                else fa[id] = fa[par.id];
            }
            else {
                fa[id] = 0;
            }
            g[fa[id]].push_back(id);
            g[id].push_back(fa[id]);
        }
        if (f == 1)s.insert(l);
        else s.erase(l);
    }

    return 0;
}

例题

Freedom from Prison

题意：给出了平面上 n 个没有公共点的凸多边形。平面上有两个点，你需要从其中一点走到另一点，再从该点走到一个不被任何凸多边形包含的地方，代价为路径上穿过的凸多边形边界数量。你可以随意安排两个点的位置，求代价的最小值最大是多少。n\leq 10^5，凸多边形总边数 6\times 10^5，值域 10^9。

这题就需要用到凸多边形的扫描线，将凸多边形之间的包含关系转化为一颗树。转化为树之后，就是求从树上一点走到另一点，再从该点走到超级源点的路径长度。显然，答案就是从直径的一端走到另一端，再走到树的根的路径长度，求个直径就做完了。

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;
#define int ll
const int maxn = 600050, mod = 998244353, INF = 1e9;
const double eps = 1e-9;
int sweep_line_x = -1e9;
vector<int> g[maxn];
int d, t;
int fa[maxn];

struct Segment {
    int id;
    pair<int, int> p, q;
    bool up_side;

    Segment() {}
    Segment(int id, pair<int, int> p, pair<int, int> q) : id(id), p(p), q(q)
    {
        if (p.first > q.first)up_side = true;
        else up_side = false;
    }
    double get_y(int x) const {
        if (p.first == q.first) return p.second;
        return p.second + 1.0 * (x - p.first) * (q.second - p.second) / (q.first - p.first);
    }

    bool operator<(const Segment& rhs) const {
        if (fabs(get_y(sweep_line_x) - rhs.get_y(sweep_line_x)) < eps) {
            if (up_side == rhs.up_side) return q.first < rhs.q.first;
            return up_side < rhs.up_side;
        }
        return get_y(sweep_line_x) + eps < rhs.get_y(sweep_line_x);
    }
};

struct Event {
    int x;
    Segment l;
    bool is_left;
    bool operator<(const Event& rhs) const {
        if (x == rhs.x) {
            if (is_left == rhs.is_left) return l.p.second > rhs.l.p.second;
            return is_left > rhs.is_left;
        }
        return x < rhs.x;
    }
};

void dfs(int x, int dep, int fath) {
    if (dep > d) {
        d = dep;
        t = x;
    }
    for (int to : g[x]) {
        if (to == fath)continue;
        dfs(to, dep+1, x);
    }
}

signed main() {
    ios::sync_with_stdio(0);
    int n;
    cin >> n;

    vector<Segment> p;
    vector<Event> events;
    for (int i = 1;i <= n;i++) {
        int m;
        cin >> m;
        vector<pair<int, int>> tmp(m);
        for (int j = 0;j < m;j++) {
            cin >> tmp[j].first >> tmp[j].second;
        }
        for (int j = 0;j < m;j++) {
            int k = (j + 1) % m;
            if (tmp[j].first == tmp[k].first) continue;
            p.push_back({ i, tmp[j], tmp[k] });
        }
    }
    for (int i = 0; i < p.size(); i++) {
        if (p[i].p.first > p[i].q.first) swap(p[i].p, p[i].q);
        events.push_back({ p[i].p.first, p[i], true });
        events.push_back({ p[i].q.first, p[i], false });
    }

    sort(events.begin(), events.end());

    memset(fa,-1,sizeof fa);
    set<Segment> s;

    for (auto& [x, l, f] : events) {
        sweep_line_x = x;
        int id = l.id;
        if (l.up_side and fa[id] == -1) {
            auto it = s.upper_bound(l);
            if (it != s.end()) {
                Segment par = (*it);
                if (par.up_side) fa[id] = par.id;
                else fa[id] = fa[par.id];
            }
            else {
                fa[id] = 0;
            }
            g[fa[id]].push_back(id);
            g[id].push_back(fa[id]);
        }
        if (f == 1)s.insert(l);
        else s.erase(l);
    }

    dfs(0, 0, -1);
    int st=t;
    dfs(t, 0, -1);
    int ans = d;

    vector<int> dis(maxn,INF);
    queue<int> q;
    dis[0] = 0;
    q.push(0);
    while (!q.empty()){
        int v = q.front();
        q.pop();
        for (auto const& to : g[v]){
            if (dis[to] > dis[v] + 1){
                dis[to] = dis[v] + 1;
                q.push(to);
            }
        }
    }
    ans += max(dis[t], dis[st]);
    cout << ans << '\n';

    return 0;
}

总结

本文主要介绍了如何使用扫描线将平面上圆或凸多边形之间的包含关系，高效地转化为树结构（森林，通过引入超级源点 0 转化为一棵树），从而将几何问题转化为树上问题，核心思想就是利用图形之间的相对关系，将图形拆成上下两部分，并在扫描线过程中使用 set 维护，时间复杂度为 \mathcal{O}(n \log n)。