convex_hull
(geometry/convex_hull.hpp)

• View this file on GitHub
• Last update: 2023-06-08 15:34:10+09:00
• Include: #include "geometry/convex_hull.hpp"

説明

base_ld を基礎とした凸包ライブラリ。

vector<vec> convex_hull(vector<vec> a)

重複のない vec の配列 a に対し、それらの点が成す凸包の頂点を反時計回りに返す。(重複があってもうまく動き得るが、保証しない)

凸包が $N$ 点からなるとき、長さ $N$ の vector が返される。

使用例

int main(){
    int n; cin >> n;
    vector<vec> a(n);
    rep(i,0,n){
        ld x, y; cin >> x >> y;
        a[i] = vec(x,y);
    }
    auto hull = convex_hull(a);
    for (vec v : hull){
        cout << v << endl;
    }
}

入力例

6
0 0
2 0
1 1
2 1
3 1
2 2

出力例

0 0
2 0
3 1
2 2

この出力例は正確ではなく、本当は (0.00000,0.00000) などと出力される。

Depends on

• base_ld (geometry/base_ld.hpp)
• テンプレート (template/template.hpp)

Verified with

• test/geometry/Convex_Hull.test.cpp

Code

#pragma once

#include "../geometry/base_ld.hpp"
#include "../template/template.hpp"

namespace lib {

vector<vec> convex_hull(vector<vec> a) {
    int n = a.size();
    if (n <= 2) return a;
    auto comp = [&](vec lhs, vec rhs) {
        if (lhs.real() == rhs.real()) return lhs.imag() < rhs.imag();
        return lhs.real() < rhs.real();
    };
    sort(all(a), comp);
    stack<int> uid, did;
    vec ri = a[n - 1];
    auto make_half = [&](bool isupper) {
        auto &id = (isupper ? uid : did);
        id.push(0);
        rep(i, 1, n - 1) {
            vec le = a[id.top()];
            auto cr = cross(ri - le, a[i] - le);
            if ((cr > 0 && isupper) || (cr < 0 && !isupper)) {
                while (id.size() >= 2) {
                    int test = id.top();
                    id.pop();
                    vec from = a[id.top()];
                    auto cr2 = cross(a[i] - from, a[test] - from);
                    if ((cr2 > 0 && isupper) || (cr2 < 0 && !isupper)) {
                        id.push(test);
                        break;
                    }
                }
                id.push(i);
            }
        }
    };
    make_half(true);
    make_half(false);
    vector<int> ids(1, n - 1);
    while (!did.empty()) ids.emplace_back(did.top()), did.pop();
    reverse(all(ids));
    while (!uid.empty()) ids.emplace_back(uid.top()), uid.pop();
    ids.pop_back();
    vector<vec> ans(ids.size());
    rep(i, 0, ids.size()) ans[i] = a[ids[i]];
    return ans;
}

}  // namespace lib

#line 2 "geometry/convex_hull.hpp"

#line 2 "geometry/base_ld.hpp"

#line 2 "template/template.hpp"

#include <bits/stdc++.h>

#define rep(i, s, n) for (int i = (int)(s); i < (int)(n); i++)
#define rrep(i, s, n) for (int i = (int)(n)-1; i >= (int)(s); i--)
#define all(v) v.begin(), v.end()

using ll = long long;
using ld = long double;
using ull = unsigned long long;

template <typename T> bool chmin(T &a, const T &b) {
    if (a <= b) return false;
    a = b;
    return true;
}
template <typename T> bool chmax(T &a, const T &b) {
    if (a >= b) return false;
    a = b;
    return true;
}

namespace lib {

using namespace std;

}  // namespace lib

// using namespace lib;
#line 4 "geometry/base_ld.hpp"

namespace lib {

using vec = complex<ld>;

const ld PI = acos(-1);

void ldout(int len = 20) {
    cout << fixed << setprecision(len);
}

int sgn(ld a, const ld eps = 1e-7) {
    return (a < -eps) ? -1 : (a > eps) ? 1 : 0;
}

bool same_vec(vec a, vec b) {
    a -= b;
    return sgn(a.real()) == 0 && sgn(a.imag()) == 0;
}

ld dot(const vec &a, const vec &b) {
    return (conj(a) * b).real();
}

ld cross(const vec &a, const vec &b) {
    return (conj(a) * b).imag();
}

int isp(const vec &a, const vec &b, const vec &c) {
    int cross_sgn = sgn(cross(b - a, c - a));
    if (cross_sgn == 0) {
        if (sgn(dot(b - a, c - a)) < 0) return -2;
        if (sgn(dot(a - b, c - b)) < 0) return 2;
    }
    return cross_sgn;
}

vec rot90(const vec &a) {
    return {-a.imag(), a.real()};
}

vec rot(const vec &a, ld rad) {
    return a * vec(cosl(rad), sinl(rad));
}

bool comp_for_argument_sort(const vec &lhs, const vec &rhs) {
    // if (abs(arg(lhs)-arg(rhs)) < eps) return false; // need ?
    return arg(lhs) < arg(rhs);
}

}  // namespace lib
#line 5 "geometry/convex_hull.hpp"

namespace lib {

vector<vec> convex_hull(vector<vec> a) {
    int n = a.size();
    if (n <= 2) return a;
    auto comp = [&](vec lhs, vec rhs) {
        if (lhs.real() == rhs.real()) return lhs.imag() < rhs.imag();
        return lhs.real() < rhs.real();
    };
    sort(all(a), comp);
    stack<int> uid, did;
    vec ri = a[n - 1];
    auto make_half = [&](bool isupper) {
        auto &id = (isupper ? uid : did);
        id.push(0);
        rep(i, 1, n - 1) {
            vec le = a[id.top()];
            auto cr = cross(ri - le, a[i] - le);
            if ((cr > 0 && isupper) || (cr < 0 && !isupper)) {
                while (id.size() >= 2) {
                    int test = id.top();
                    id.pop();
                    vec from = a[id.top()];
                    auto cr2 = cross(a[i] - from, a[test] - from);
                    if ((cr2 > 0 && isupper) || (cr2 < 0 && !isupper)) {
                        id.push(test);
                        break;
                    }
                }
                id.push(i);
            }
        }
    };
    make_half(true);
    make_half(false);
    vector<int> ids(1, n - 1);
    while (!did.empty()) ids.emplace_back(did.top()), did.pop();
    reverse(all(ids));
    while (!uid.empty()) ids.emplace_back(uid.top()), uid.pop();
    ids.pop_back();
    vector<vec> ans(ids.size());
    rep(i, 0, ids.size()) ans[i] = a[ids[i]];
    return ans;
}

}  // namespace lib

Back to top page