/* Utilities for generating and validating 2D polygons.
 *
 * Contains the following algorithms:
 *
 * bool simple_polygon(polygon<coord> P)
 *   Check if a polygon is simple (no line segment intersections,
 *   except that consecutive segments intersect in their common
 *   vertex).  The polygon can be given in either clockwise or
 *   counter-clockwise order.  When using a coordinate type with exact
 *   arithmetic (such as integers or rationals), the check is exact
 *   (assuming the type is large enough, see Assumptions below).
 *   Using a floating point coordinate type also works in principle
 *   but the results may not be exact.
 *   Assumptions: - the coordinate type is signed
 *                - area of bounding box of the polygon points union
 *                  origin fits in the coordinate type.
 *   Complexity: O(N*log(N)) arithmetic operations on coordinates
 *   Special cases: returns false if P has fewer than 3 points.
 *
 *
 * bool ccw(polygon<coord> P)
 *   Check if a polygon is given in counter-clockwise order.
 *   Assumptions: - P is simple
 *                - the coordinate type is signed
 *                - the area of the polygon fits in the coordinate type
 *   Complexity: O(N) arithmetic operations
 *
 *
 * bool line_line_isect(const point<Coord> &p1, const point<Coord> &p2,
 *   	                const point<Coord> &q1, const point<Coord> &q2,
 *                      line_type p1_end = CLOSED, line_type p2_end = CLOSED,
 *                      line_type q1_end = CLOSED, line_type q2_end = CLOSED)
 *   Check if two lines/line segments intersect.
 *      p1_end/p2_end and q1_end/q2_end indicate whether the
 *      lines/line end at the respective endpoint (and if so, if the
 *      endpoint is included or not), or if they extend infinitely.
 *      Possible values are OPEN, CLOSED, and INFINITE.
 *
 *
 * Author: Per Austrin
 */

#include <cassert>
#include <vector>
#include <set>
#include <algorithm>
#include <utility>
#include <tuple>

template <typename T>
int sgn(T x) { return x < 1e-9 ? x < -1e-9 ? -1 : 0 : 1; }

template <typename Coord>
struct point {
    Coord x, y;
    point(Coord x=Coord(), Coord y=Coord()): x(x), y(y) {}
    point operator+(const point &p) const { return point(x+p.x, y+p.y); }
    point operator-(const point &p) const { return point(x-p.x, y-p.y); }
    bool operator<(const point &p) const {
        return x != p.x ? x < p.x : y < p.y;
    }
    bool operator==(const point &p) const {
        return x == p.x && y == p.y;
    }
    Coord dot(const point &p) const {
        return x*p.x + y*p.y;
    }
    Coord cross(const point &p) const {
        return x*p.y - y*p.x;
    }
    Coord dist2() const {
        return dot(*this);
    }
};

template <typename Coord>
using polygon = std::vector<point<Coord>>;

typedef enum {INFINITE=-1, OPEN=0, CLOSED=1} line_type;

template <typename Coord>
bool __inside_ival(const Coord &t, const Coord &B,
                   line_type left_end = CLOSED, line_type right_end = CLOSED) {
    if (t < 0) return left_end == INFINITE;
    if (t > B) return right_end == INFINITE;
    if (t == 0) return left_end != OPEN;
    if (t == B) return right_end != OPEN;
    return true;
}

template <typename Coord>
bool line_line_isect(const point<Coord> &p1, const point<Coord> &p2,
                     const point<Coord> &q1, const point<Coord> &q2,
                     line_type p1_end = CLOSED, line_type p2_end = CLOSED,
                     line_type q1_end = CLOSED, line_type q2_end = CLOSED) {
    const point<Coord> dp = p2-p1, dq = q2-q1, dpq = p1-q1;
    Coord det = dp.cross(dq);
    if (det == 0) {
        if (dp.cross(dpq) != 0) return false;
        Coord t1 = dp.dot(q1-p1), t2 = dp.dot(q2-p1);
        Coord s1 = dq.dot(p1-q1), s2 = dq.dot(p2-q1);
        Coord pL = dp.dot(dp), qL = dq.dot(dq);
        return (__inside_ival(t1, pL, std::min(p1_end, q1_end), std::min(p2_end, q1_end)) ||
                __inside_ival(t2, pL, std::min(p1_end, q2_end), std::min(p2_end, q2_end)) ||
                __inside_ival(s1, qL, std::min(q1_end, p1_end), std::min(q2_end, p1_end)) ||
                __inside_ival(s2, qL, std::min(q1_end, p2_end), std::min(q2_end, p2_end)));
    }
    Coord s = dq.cross(dpq);
    Coord t = dp.cross(dpq);
    if (det < 0) {
        det = -det;
        s = -s;
        t = -t;
    }
    return __inside_ival(s, det, p1_end, p2_end) && __inside_ival(t, det, q1_end, q2_end);
}

template <typename Coord>
bool simple_polygon(const polygon<Coord>& poly) {
    typedef std::pair<point<Coord>, int> event;

    struct line {
        point<Coord> p1, p2, dp;
        bool reversed;

        line(const point<Coord> &_p1, const point<Coord> &_p2): p1(_p1), p2(_p2) {
            reversed = p2 < p1;
            if (reversed) std::swap(p1, p2);
            dp = p2-p1;
        }

        Coord cmp(const line &L) const {
            if (L.p1.x < p1.x) return -L.cmp(*this);
            Coord res = dp.cross(L.p1-p1);
            return res ? res : dp.cross(L.p2-p1);
        }

        bool operator<(const line &L) const { return cmp(L) < 0; }

        // NB this function will not give the right answer for two colinear
        // lines (it will say that they overlap, regardless of whether or
        // not they do).  However in the polygon sweep algorithm, this
        // function will never be called for two non-overlapping colinear
        // lines
        bool intersects(const line &L) const {
            return line_line_isect(p1, p2, L.p1, L.p2,
                                   reversed ? CLOSED : OPEN, reversed ? OPEN : CLOSED,
                                   L.reversed ? CLOSED : OPEN, L.reversed ? OPEN : CLOSED);
        }
    };

    int n = poly.size();
    if (n < 3) return false;

    std::vector<event> events;
    for (int j = 0, i = n-1; j < n; i = j++) {
        auto A = poly[j], B = poly[i];
        if (B < A) std::swap(A, B);
        events.push_back(event(A, +(i+1)));
        events.push_back(event(B, -(i+1)));
    }
    sort(events.begin(), events.end());

    std::set<line> frontier;
    for (auto ev: events) {
        point<Coord> X;
        int type;
        std::tie(X, type) = ev;
        line L(poly[abs(type)-1], poly[abs(type)%n]);
        auto it = frontier.end();
        if (type > 0) {
            auto res = frontier.insert(L);
            if (!res.second) return false;
            it = res.first;
        } else {
            it = frontier.find(L);
            assert(it != frontier.end());
            frontier.erase(it++);
        }
        if (it != frontier.end()) {
            auto jt = it;
            if (it != frontier.begin() && it->intersects(*--jt)) {
                return false;
            }
            jt = it;
            if (++jt != frontier.end() && it->intersects(*jt)) {
                return false;
            }
        }
    }
    return true;
}

template <typename Coord>
bool ccw(polygon<Coord> poly) {
    Coord A = 0;
    int n = poly.size();
    for (int j = 0, i = n-1; j < n; i = j++)
        A += poly[i].cross(poly[j]);
    return A > 0;
}