콘텐츠로 이동

점이 볼록 다각형에 속하는지 $O(\log N)$에 확인하기

다음 문제를 생각해 봅시다. 정수 좌표의 정점을 가진 볼록 다각형과 많은 수의 쿼리가 주어집니다. 각 쿼리는 하나의 점이며, 그 점이 다각형의 내부 또는 경계 위에 있는지 여부를 판별해야 합니다. 다각형의 꼭짓점이 반시계 방향으로 정렬되어 있다고 가정합니다. 각 쿼리에 $O(\log n)$에 온라인으로 답하겠습니다.

알고리즘

x좌표가 가장 작은 점을 고릅시다. 그러한 점이 여러 개라면 그중 y좌표가 가장 작은 점을 고릅니다. 이 점을 $p_0$로 표기합시다. 이제 다각형의 나머지 모든 점 $p_1,\dots,p_n$은 선택한 점을 기준으로 한 편각 순서로 정렬되어 있습니다(다각형이 반시계 방향으로 정렬되어 있기 때문입니다).

점이 다각형에 속한다면 그 점은 어떤 삼각형 $p_0, p_i, p_{i + 1}$에 속합니다(삼각형들의 경계 위에 있다면 둘 이상의 삼각형에 속할 수도 있습니다). $p$가 속하는 삼각형 $p_0, p_i, p_{i + 1}$ 중에서 $i$가 최대인 삼각형을 생각해 봅시다.

한 가지 특수한 경우가 있습니다. $p$가 선분 $(p_0, p_1)$ 위에 있는 경우입니다. 이 경우는 따로 확인하겠습니다. 그렇지 않다면 $j \le i$인 모든 점 $p_j$는 $p_0$를 기준으로 $p$에서 반시계 방향에 있고, 나머지 모든 점은 $p$에서 반시계 방향에 있지 않습니다. 따라서 $p_0$을 기준으로 $p_i$가 $p$로부터 반시계 방향에 있지 않으면서 그러한 모든 점 중 $i$가 최대인 점 $p_i$를 이진 탐색으로 찾을 수 있습니다. 그다음 그 점이 실제로 결정된 삼각형 안에 있는지 확인합니다.

$(a - c) \times (b - c)$의 부호를 보면 점 $c$를 기준으로 점 $a$가 점 $b$에서 시계 방향에 있는지 반시계 방향에 있는지 알 수 있습니다. $(a - c) \times (b - c) > 0$이면 점 $a$는 $c$에서 $b$로 향하는 벡터의 오른쪽에 있으며, 이는 $c$를 기준으로 $b$에서 시계 방향에 있다는 뜻입니다. 그리고 $(a - c) \times (b - c) < 0$이면 그 점은 왼쪽, 즉 반시계 방향에 있습니다. 0이면 점 $a$는 정확히 점 $b$와 $c$를 지나는 직선 위에 있습니다.

다시 알고리즘으로 돌아갑시다. 쿼리 점 $p$를 생각해 봅시다. 먼저, 점이 $p_1$과 $p_n$ 사이에 있는지 확인해야 합니다. 그렇지 않다면 그 점은 다각형에 속할 수 없다는 것을 이미 알 수 있습니다. 이는 외적 $(p_1 - p_0)\times(p - p_0)$이 0이거나 $(p_1 - p_0)\times(p_n - p_0)$과 부호가 같은지, 그리고 $(p_n - p_0)\times(p - p_0)$이 0이거나 $(p_n - p_0)\times(p_1 - p_0)$과 부호가 같은지를 확인하여 수행할 수 있습니다. 그다음 $p$가 직선 $(p_0, p_1)$ 위에 있는 특수한 경우를 처리합니다. 그다음 $p_1,\dots p_n$ 중에서 $p_0$을 기준으로 $p$로부터 반시계 방향에 있지 않은 마지막 점을 이진 탐색으로 찾을 수 있습니다. 하나의 점 $p_i$에 대해서는 $(p_i - p_0)\times(p - p_0) \le 0$인지 확인하여 이 조건을 검사할 수 있습니다. 이러한 점 $p_i$를 찾은 후에는 $p$가 삼각형 $p_0, p_i, p_{i + 1}$ 내부에 있는지 검사해야 합니다. 삼각형에 속하는지 검사하려면 간단히 $|(p_i - p_0)\times(p_{i + 1} - p_0)| = |(p_0 - p)\times(p_i - p)| + |(p_i - p)\times(p_{i + 1} - p)| + |(p_{i + 1} - p)\times(p_0 - p)|$인지 확인하면 됩니다. 이는 삼각형 $p_0, p_i, p_{i+1}$의 넓이가 삼각형 $p_0, p_i, p$, 삼각형 $p_0, p, p_{i+1}$, 삼각형 $p_i, p_{i+1}, p$의 넓이의 합과 정확히 같은지를 확인하는 것입니다. $p$가 외부에 있다면 세 삼각형의 넓이의 합은 원래 삼각형의 넓이보다 커집니다. 내부에 있다면 두 값은 같아집니다.

구현

함수 prepare는 사전순으로 가장 작은 점(x 값이 가장 작고, 같을 경우 y 값이 가장 작은 점)이 $p_0$이 되도록 하고, 벡터 $p_i - p_0$을 계산합니다. 그 후 함수 pointInConvexPolygon이 쿼리의 결과를 계산합니다. 벡터에는 시작점이 없으므로, 올바른 거리를 계산하기 위해 점 $p_0$을 추가로 기억해 두고 모든 질의 점을 이 점에 대해 상대적으로 평행 이동합니다. 쿼리 점들을 평행 이동하면 모든 벡터가 원점 $(0, 0)$에서 시작한다고 가정할 수 있으며, 거리와 길이 계산을 단순화할 수 있습니다.

struct pt {
    long long x, y;
    pt() {}
    pt(long long _x, long long _y) : x(_x), y(_y) {}
    pt operator+(const pt &p) const { return pt(x + p.x, y + p.y); }
    pt operator-(const pt &p) const { return pt(x - p.x, y - p.y); }
    long long cross(const pt &p) const { return x * p.y - y * p.x; }
    long long dot(const pt &p) const { return x * p.x + y * p.y; }
    long long cross(const pt &a, const pt &b) const { return (a - *this).cross(b - *this); }
    long long dot(const pt &a, const pt &b) const { return (a - *this).dot(b - *this); }
    long long sqrLen() const { return this->dot(*this); }
};

bool lexComp(const pt &l, const pt &r) {
    return l.x < r.x || (l.x == r.x && l.y < r.y);
}

int sgn(long long val) { return val > 0 ? 1 : (val == 0 ? 0 : -1); }

vector<pt> seq;
pt translation;
int n;

bool pointInTriangle(pt a, pt b, pt c, pt point) {
    long long s1 = abs(a.cross(b, c));
    long long s2 = abs(point.cross(a, b)) + abs(point.cross(b, c)) + abs(point.cross(c, a));
    return s1 == s2;
}

void prepare(vector<pt> &points) {
    n = points.size();
    int pos = 0;
    for (int i = 1; i < n; i++) {
        if (lexComp(points[i], points[pos]))
            pos = i;
    }
    rotate(points.begin(), points.begin() + pos, points.end());

    n--;
    seq.resize(n);
    for (int i = 0; i < n; i++)
        seq[i] = points[i + 1] - points[0];
    translation = points[0];
}

bool pointInConvexPolygon(pt point) {
    point = point - translation;
    if (seq[0].cross(point) != 0 &&
            sgn(seq[0].cross(point)) != sgn(seq[0].cross(seq[n - 1])))
        return false;
    if (seq[n - 1].cross(point) != 0 &&
            sgn(seq[n - 1].cross(point)) != sgn(seq[n - 1].cross(seq[0])))
        return false;

    if (seq[0].cross(point) == 0)
        return seq[0].sqrLen() >= point.sqrLen();

    int l = 0, r = n - 1;
    while (r - l > 1) {
        int mid = (l + r) / 2;
        int pos = mid;
        if (seq[pos].cross(point) >= 0)
            l = mid;
        else
            r = mid;
    }
    int pos = l;
    return pointInTriangle(seq[pos], seq[pos + 1], pt(0, 0), point);
}

문제