多角形に関する計算

polygon_calc.py

import numpy as np
from numpy.linalg import norm

#------------------------------------------------------------------------------
# 凸多角形か判定
def isConvex(polygon):
    # n角形
    n = len(polygon)
    # 時計回りまたは反時計回りで、隣接する辺のベクトルのクロス積がすべて同符号なら凸多角形
    cp0 = 0
    for i in range(n):
        p1 = polygon[i]
        p2 = polygon[(i + 1) % n]
        p3 = polygon[(i + 2) % n]
        v21 = p2 - p1
        v32 = p3 - p2
        cp = v21[0] * v32[1] - v21[1] * v32[0]
        if i == 0: cp0 = cp
        elif cp0 * cp < 0:
            return False
    return True

print("Test: isConvex")
polygon = np.array([[  0,  0], [1, 1], [2, 0], [1, -1]])
print(isConvex(polygon))  # True

polygon = np.array([[  1, -1], [2, 0], [1, 1], [0,  0]])
print(isConvex(polygon))  # True

polygon = np.array([[1.5,  0], [1, 1], [2, 0], [1, -1]])
print(isConvex(polygon))  # False

#------------------------------------------------------------------------------
# 凸多角形の内部か判定
def isInConvex(polygon, point):
    # n角形
    n = len(polygon)
    # 時計回りまたは反時計回りで、辺のベクトルと頂点から点のベクトルのクロス積がすべて同符号なら内部の点
    cp0 = 0
    for i in range(n):
        p1 = polygon[i]
        p2 = polygon[(i + 1) % n]
        p3 = point
        v21 = p2 - p1
        v31 = p3 - p1
        cp = v21[0] * v31[1] - v21[1] * v31[0]
        if i == 0: cp0 = cp
        elif cp0 * cp < 0:
            return False
    return True

print("Test: isInConvex")
polygon = np.array([[0,  0], [1, 1], [2, 0], [1, -1]])
point = np.array([1,  0])
print(isInConvex(polygon, point))  # True

polygon = np.array([[1, -1], [2, 0], [1, 1], [0,  0]])
point = np.array([1,  0])
print(isInConvex(polygon, point))  # True

polygon = np.array([[0,  0], [1, 1], [2, 0], [1, -1]])
point = np.array([4,  0])
print(isInConvex(polygon, point))  # False

#------------------------------------------------------------------------------
# 点から多角形までの最短距離
def distance_to_poligon(polygon, point):
    # n角形
    n = len(polygon)
    # 各辺と点の距離の最小値
    for i in range(n):
        p1 = polygon[i]
        p2 = polygon[(i + 1) % n]
        p3 = point
        _d, _p4 = distance_to_segment(p1, p2, p3)
        if  i == 0 or _d < d:
            d = _d
            p4 = _p4
    return d, p4

# 点P3から線分V21までの最短距離dと最近傍点P4
def distance_to_segment(p1, p2, p3):
    v31 = p3 - p1
    v21 = p2 - p1
    v12 = p1 - p2
    v32 = p3 - p2
    # P1が最近傍点の場合
    if np.dot(v31, v21) < 0:
        d = norm(v31)
        p4 = p1
    # P2が最近傍点の場合
    elif np.dot(v32, v12) < 0:
        d = norm(v32)
        p4 = p2
    # P3からV21におろした垂点が最近傍点の場合
    else:
        p41_norm = np.dot(v31, v21) / norm(v21)
        p4 = p1 + v21 / norm(v21) * p41_norm
        d = norm(p3 - p4)
    return d, p4

print("Test: distance_to_segment")
p1 = np.array([1,  1])
p2 = np.array([2,  2])
p3 = np.array([2,  1])
print(distance_to_segment(p1, p2, p3))

print("Test: distance_to_poligon")
polygon = np.array([[0,  0], [1, 1], [2, 0], [1, -1]])
point = np.array([0,  1])
print(distance_to_poligon(polygon, point))
point = np.array([1,  2])
print(distance_to_poligon(polygon, point))