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
参考
-
O’rourke, Joseph. Computational geometry in C. Cambridge university press, 1998.