compro_library

This documentation is automatically generated by online-judge-tools/verification-helper

View the Project on GitHub siro53/compro_library

:heavy_check_mark: geometry/convex-cut.hpp

Depends on

Required by

Verified with

Code

#pragma once

#include <vector>

#include "ccw.hpp"
#include "cross-point.hpp"

namespace geometry {
    // 凸多角形pを直線lで切断し、その左側を返す
    inline std::vector<Point> ConvexCut(std::vector<Point> p, Line l) {
        std::vector<Point> ret;
        int sz = (int)p.size();
        for(int i = 0; i < sz; i++) {
            Point now = p[i];
            Point nxt = p[i == sz - 1 ? 0 : i + 1];
            if(ccw(l.a, l.b, now) != -1) ret.emplace_back(now);
            if(ccw(l.a, l.b, now) * ccw(l.a, l.b, nxt) < 0) {
                ret.emplace_back(crossPoint(Line(now, nxt), l));
            }
        }
        return ret;
    }
} // namespace geometry
#line 2 "geometry/convex-cut.hpp"

#include <vector>

#line 2 "geometry/ccw.hpp"

#line 2 "geometry/cross.hpp"

#line 2 "geometry/base.hpp"

#include <cmath>
#include <complex>

namespace geometry {
    // Point : 複素数型を位置ベクトルとして扱う
    // 実軸(real)をx軸、挙軸(imag)をy軸として見る
    using D = long double;
    using Point = std::complex<D>;
    const D EPS = 1e-7;
    const D PI = std::acos(D(-1));

    inline bool equal(const D &a, const D &b) { return std::fabs(a - b) < EPS; }
} // namespace geometry
#line 4 "geometry/cross.hpp"

namespace geometry {
    // 外積(cross product) : a×b = |a||b|sinΘ
    inline D cross(const Point &a, const Point &b) {
        return (a.real() * b.imag() - a.imag() * b.real());
    }
} // namespace geometry
#line 2 "geometry/dot.hpp"

#line 4 "geometry/dot.hpp"

namespace geometry {
    // 内積(dot product) : a・b = |a||b|cosΘ
    inline D dot(const Point &a, const Point &b) {
        return (a.real() * b.real() + a.imag() * b.imag());
    }
} // namespace geometry
#line 5 "geometry/ccw.hpp"

namespace geometry {
    // 点の回転方向
    // 点a, b, cの位置関係について(aが基準点)
    inline int ccw(const Point &a, Point b, Point c) {
        b -= a, c -= a;
        // 点a, b, c が
        // 反時計回りの時、
        if(cross(b, c) > EPS) return 1;
        // 時計回りの時、
        if(cross(b, c) < -EPS) return -1;
        // c, a, bがこの順番で同一直線上にある時、
        if(dot(b, c) < 0) return 2;
        // a, b, cがこの順番で同一直線上にある場合、
        if(std::norm(b) < std::norm(c)) return -2;
        // cが線分ab上にある場合、
        return 0;
    }
} // namespace geometry
#line 2 "geometry/cross-point.hpp"

#line 5 "geometry/cross-point.hpp"

#line 2 "geometry/circle.hpp"

#line 4 "geometry/circle.hpp"

namespace geometry {
    // Circle : 円を表す構造体
    // pが中心の位置ベクトル、rは半径
    struct Circle {
        Point p;
        D r;

        Circle() = default;

        Circle(Point p, D r) : p(p), r(r) {}
    };
} // namespace geometry
#line 2 "geometry/distance-between-line-and-point.hpp"

#line 2 "geometry/line.hpp"

#line 4 "geometry/line.hpp"

namespace geometry {
    // Line : 直線を表す構造体
    // b - a で直線・線分を表せる
    struct Line {
        Point a, b;
        Line() = default;
        Line(Point a, Point b) : a(a), b(b) {}
        // Ax+By=C
        Line(D A, D B, D C) {
            if(equal(A, 0)) {
                a = Point(0, C / B), b = Point(1, C / B);
            } else if(equal(B, 0)) {
                a = Point(C / A, 0), b = Point(C / A, 1);
            } else if(equal(C, 0)) {
                a = Point(0, C / B), b = Point(1, (C - A) / B);
            } else {
                a = Point(0, C / B), b = Point(C / A, 0);
            }
        }
    };
} // namespace geometry
#line 5 "geometry/distance-between-line-and-point.hpp"

namespace geometry {
    // 直線lと点pの距離を求める
    inline D distanceBetweenLineAndPoint(const Line &l, const Point &p) {
        return std::abs(cross(l.b - l.a, p - l.a)) / std::abs(l.b - l.a);
    }
} // namespace geometry
#line 2 "geometry/is-intersect.hpp"

#line 2 "geometry/segment.hpp"

#line 4 "geometry/segment.hpp"

namespace geometry {
    // Segment : 線分を表す構造体
    // Lineと同じ
    struct Segment : Line {
        Segment() = default;

        Segment(Point a, Point b) : Line(a, b) {}
        D get_dist() { return std::abs(a - b); }
    };
} // namespace geometry
#line 6 "geometry/is-intersect.hpp"

namespace geometry {
    // 線分sと線分tが交差しているかどうか
    // bound:線分の端点を含むか
    inline bool isIntersect(const Segment &s, const Segment &t, bool bound) {
        return ccw(s.a, s.b, t.a) * ccw(s.a, s.b, t.b) < bound &&
               ccw(t.a, t.b, s.a) * ccw(t.a, t.b, s.b) < bound;
    }

    // 2つの円の交差判定
    // 返り値は共通接線の数
    inline int isIntersect(const Circle &c1, const Circle &c2) {
        D d = std::abs(c1.p - c2.p);
        // 2つの円が離れている場合
        if(d > c1.r + c2.r + EPS) return 4;
        // 外接している場合
        if(equal(d, c1.r + c2.r)) return 3;
        // 内接している場合
        if(equal(d, std::abs(c1.r - c2.r))) return 1;
        // 内包している場合
        if(d < std::abs(c1.r - c2.r) - EPS) return 0;
        return 2;
    }
} // namespace geometry
#line 2 "geometry/projection.hpp"

#line 6 "geometry/projection.hpp"

namespace geometry {
    // 射影(projection)
    // 直線(線分)lに点pから引いた垂線の足を求める
    inline Point projection(const Line &l, const Point &p) {
        D t = dot(p - l.a, l.a - l.b) / std::norm(l.a - l.b);
        return l.a + (l.a - l.b) * t;
    }

    inline Point projection(const Segment &l, const Point &p) {
        D t = dot(p - l.a, l.a - l.b) / std::norm(l.a - l.b);
        return l.a + (l.a - l.b) * t;
    }
} // namespace geometry
#line 2 "geometry/unit-vector.hpp"

#line 4 "geometry/unit-vector.hpp"

namespace geometry {
    // 単位ベクトル(unit vector)を求める
    inline Point unitVector(const Point &a) { return a / std::abs(a); }
} // namespace geometry
#line 12 "geometry/cross-point.hpp"

namespace geometry {
    // 直線s, tの交点の計算
    inline Point crossPoint(const Line &s, const Line &t) {
        D d1 = cross(s.b - s.a, t.b - t.a);
        D d2 = cross(s.b - s.a, s.b - t.a);
        if(equal(std::abs(d1), 0) && equal(std::abs(d2), 0)) return t.a;
        return t.a + (t.b - t.a) * (d2 / d1);
    }

    // 線分s, tの交点の計算
    inline Point crossPoint(const Segment &s, const Segment &t) {
        return crossPoint(Line(s), Line(t));
    }

    // 2つの円の交点
    inline std::vector<Point> crossPoint(const Circle &c1, const Circle &c2) {
        std::vector<Point> res;
        int mode = isIntersect(c1, c2);
        // 2つの中心の距離
        D d = std::abs(c1.p - c2.p);
        // 2円が離れている場合
        if(mode == 4) return res;
        // 1つの円がもう1つの円に内包されている場合
        if(mode == 0) return res;
        // 2円が外接する場合
        if(mode == 3) {
            D t = c1.r / (c1.r + c2.r);
            res.emplace_back(c1.p + (c2.p - c1.p) * t);
            return res;
        }
        // 内接している場合
        if(mode == 1) {
            if(c2.r < c1.r - EPS) {
                res.emplace_back(c1.p + (c2.p - c1.p) * (c1.r / d));
            } else {
                res.emplace_back(c2.p + (c1.p - c2.p) * (c2.r / d));
            }
            return res;
        }
        // 2円が重なる場合
        D rc1 = (c1.r * c1.r + d * d - c2.r * c2.r) / (2 * d);
        D rs1 = std::sqrt(c1.r * c1.r - rc1 * rc1);
        if(c1.r - std::abs(rc1) < EPS) rs1 = 0;
        Point e12 = (c2.p - c1.p) / std::abs(c2.p - c1.p);
        res.emplace_back(c1.p + rc1 * e12 + rs1 * e12 * Point(0, 1));
        res.emplace_back(c1.p + rc1 * e12 + rs1 * e12 * Point(0, -1));
        return res;
    }

    // 円cと直線lの交点
    inline std::vector<Point> crossPoint(const Circle &c, const Line &l) {
        std::vector<Point> res;
        D d = distanceBetweenLineAndPoint(l, c.p);
        // 交点を持たない
        if(d > c.r + EPS) return res;
        // 接する
        Point h = projection(l, c.p);
        if(equal(d, c.r)) {
            res.emplace_back(h);
            return res;
        }
        Point e = unitVector(l.b - l.a);
        D ph = std::sqrt(c.r * c.r - d * d);
        res.emplace_back(h - e * ph);
        res.emplace_back(h + e * ph);
        return res;
    }
} // namespace geometry
#line 7 "geometry/convex-cut.hpp"

namespace geometry {
    // 凸多角形pを直線lで切断し、その左側を返す
    inline std::vector<Point> ConvexCut(std::vector<Point> p, Line l) {
        std::vector<Point> ret;
        int sz = (int)p.size();
        for(int i = 0; i < sz; i++) {
            Point now = p[i];
            Point nxt = p[i == sz - 1 ? 0 : i + 1];
            if(ccw(l.a, l.b, now) != -1) ret.emplace_back(now);
            if(ccw(l.a, l.b, now) * ccw(l.a, l.b, nxt) < 0) {
                ret.emplace_back(crossPoint(Line(now, nxt), l));
            }
        }
        return ret;
    }
} // namespace geometry
Back to top page