计算几何算法合集

//#include <algorithm>
namespace Geometry {
#define eps (1e-8)
    class point {
    public:
        double x, y;
        point() {}
        point(const point &p): x(p.x), y(p.y) {}
        point(double a, double b): x(a), y(b) {}
        point operator + (const point & p)const {
            point ret;
            ret.x = x + p.x, ret.y = y + p.y;
            return ret;
        }
        point operator - (const point & p)const {
            point ret;
            ret.x = x - p.x, ret.y = y - p.y;
            return ret;
        }
        //dot product
        double operator * (const point & p)const {
            return x * p.x + y * p.y;
        }
        //cross product
        double operator ^ (const point & p)const {
            return x * p.y - p.x * y;
        }
        bool operator < (const point & p)const {
            if (fabs(x - p.x) < eps) {
                return y < p.y;
            }
            return x < p.x;
        }
        double mold() {
            return sqrt(x * x + y * y);
        }
    };
    double cp(point a, point b, point o) {
        return (a - o) ^ (b - o);
    }
    double dp(point a, point b, point o) {
        return (a - o) * (b - o);
    }
    class line {
    public:
        point A, B;
        line() {}
        line(point a, point b): A(a), B(b) {}
        bool IsLineCrossed(const line &l)const {
            point v1, v2;
            double c1, c2;
            v1 = B - A, v2 = l.A - A;
            c1 = v1 ^ v2;
            v2 = l.B - A;
            c2 = v1 ^ v2;
            if (c1 * c2 >= 0) {
                return false;
            }
            v1 = l.B - l.A, v2 = A - l.A;
            c1 = v1 ^ v2;
            v2 = B - l.A;
            c2 = v1 ^ v2;
            if (c1 * c2 >= 0) {
                return false;
            }
            return true;
        }
    };
    /*
    **  get the convex closure of dot set,store in array s.
    **  return the amount of the dot in the convex closure
    */
    int Graham(point * p, point * s, int n) {
        std::sort(p, p + n);
        int top, m;
        s[0] = p[0];
        s[1] = p[1];
        top = 1;
        for (int i = 2; i < n; i++) {
            while (top > 0 && cp(p[i], s[top], s[top - 1]) >= 0) {
                top--;
            }
            s[++top] = p[i];
        }
        m = top;
        s[++top] = p[n - 2];
        for (int i = n - 3; i >= 0; i--) {
            while (top > m && cp(p[i], s[top], s[top - 1]) >= 0) {
                top--;
            }
            s[++top] = p[i];
        }
        return top;
    }
    int dcmp(double x) {
        if (x < -eps) {
            return -1;
        } else {
            return (x > eps);
        }
    }
    //if the point p0 on the segment consists of point p1 and p2
    int PointOnSegment(point p0, point p1, point p2) {
        return dcmp(cp(p1, p2, p0)) == 0 && dcmp(dp(p1, p2, p0)) <= 0;
    }
    /*
    **  if the point pt in polygon consists of the dots in array p
    **  0:outside
    **  1:inside
    **  2:on the border
    */
    int PointInPolygon(point pt, point * p, int n) {
        int i, k, d1, d2, wn = 0;
        p[n] = p[0];
        for (i = 0; i < n; i++) {
            if (PointOnSegment(pt, p[i], p[i + 1])) {
                return 2;
            }
            k = dcmp(cp(p[i + 1], pt, p[i]));
            d1 = dcmp(p[i + 0].y - pt.y);
            d2 = dcmp(p[i + 1].y - pt.y);
            if (k > 0 && d1 <= 0 && d2 > 0) {
                wn++;
            }
            if (k < 0 && d2 <= 0 && d1 > 0) {
                wn--;
            }
        }
        return wn != 0 ? 1 : 0;
    }
}
//using namespace Geometry;