Check if a point is inside a polygon

概要

多角形 \(P = (v_1, v_2, ..., v_N)\) に点 \(w\) が含まれているかを判定する。

Ray Casting Algorithm

点\(w\) から半直線を出し、多角形\(P\)と何回交差するかで包含判定する。

多角形\(P\) と奇数回交差している場合に包含されていると判定する。

計算量: \(O(N)\)

Winding Number Algorithm

点\(w\) から見て、多角形\(P\) 上を移動する時に何回転するかによって包含判定を行う。この時、回転しなければ包含されず、1回転以上すれば包含していることを意味する。

半直線 \(wv_i\) と 半直線 \(wv_{i+1}\) のなす角 \(\theta_i\) を求め、その総和 \(\displaystyle \sum_{1 \le i \le N} \theta_i\) を確認する。

この値が 0 以外 (この時は \(2 \pi\) またはその倍数) になれば包含され、0の場合は包含されないことが分かる。

計算量: \(O(N)\)

実装

多角形の頂点は反時計回りのリストとして持つとする。

from math import atan2, pi

# Ray casting algorithm
def inside_polygon(p0, qs):
    cnt = 0
    L = len(qs)
    x, y = p0
    for i in range(L):
        x0, y0 = qs[i-1]; x1, y1 = qs[i]
        x0 -= x; y0 -= y
        x1 -= x; y1 -= y

        cv = x0*x1 + y0*y1
        sv = x0*y1 - x1*y0
        if sv == 0 and cv <= 0:
            # a point is on a segment
            return True

        if not y0 < y1:
            x0, x1 = x1, x0
            y0, y1 = y1, y0

        if y0 <= 0 < y1 and x0*(y1 - y0) > y0*(x1 - x0):
            cnt += 1
    return (cnt % 2 == 1)

# Winding number algorithm
def inside_polygon1(p0, qs):
    x, y = p0
    L = len(qs)
    theta = 0
    for i in range(L):
        x0, y0 = qs[i-1]; x1, y1 = qs[i]
        x0 -= x; y0 -= y
        x1 -= x; y1 -= y

        cv = x0*x1 + y0*y1
        sv = x0*y1 - x1*y0
        if sv == 0 and cv <= 0:
            # a point is on a segment
            return True

        theta += atan2(sv, cv)
    return abs(theta) > 1

Verified

Ray Casting Algorithm

  • AOJ: "CGL_3_C: Polygon-Point Containment": source (Python3, 0.07sec)

Winding Number Algorithm

  • AOJ: "CGL_3_C: Polygon-Point Containment": source (Python3, 0.11sec)

参考