凸包

• Spaghetti Source - 凸包 (Andrew's Monotone Chain)

上記サイト様の写経。

x,y座標は、複素数で与えることを前提とする。

def dot(a, b):
    return (a.conjugate() * b).real


def cross(a, b):
    return (a.conjugate() * b).imag


def ccw(a, b, c):
    b -= a
    c -= a
    cbc = cross(b, c)
    if cbc > 0:
        return 1
    if cbc < 0:
        return -1
    if dot(b, c) < 0:
        return 2
    if abs(b) < abs(c):
        return -2
    return 0


def convex_hull(n, ps):
    ps.sort(key=lambda p: (p.real, p.imag))
    hull = [0j] * (2 * n)
    i = 0
    k = 0
    while i < n:
        while k >= 2 and ccw(hull[k - 2], hull[k - 1], ps[i]) <= 0:
            k -= 1
        hull[k] = ps[i]
        i += 1
        k += 1
    i = n - 2
    t = k + 1
    while i >= 0:
        while k >= t and ccw(hull[k - 2], hull[k - 1], ps[i]) <= 0:
            k -= 1
        hull[k] = ps[i]
        i -= 1
        k += 1
    return hull[:k - 1]

凸包を構築し、点は時計回りに並んでいるという前提。

def convex_diameter(n, ps):
    iis, ijs = 0, 0
    for i in range(1, n):
        if ps[i].imag > ps[iis].imag:
            iis = i
        if ps[i].imag < ps[ijs].imag:
            ijs = i
    max_d = abs(ps[iis] - ps[ijs])
    i = max_i = iis
    j = max_j = ijs

    def diff(i):
        return ps[(i + 1) % n] - ps[i]

    def check(i, j):
        nonlocal max_d, max_i, max_j
        if cross(diff(i), diff(j)) >= 0:
            j = (j + 1) % n
        else:
            i = (i + 1) % n
        curr = abs(ps[i] - ps[j])
        if curr > max_d:
            max_d = curr
            max_i = i
            max_j = j
        return i, j

    i, j = check(i, j)
    while i != iis or j != ijs:
        i, j = check(i, j)

    return max_d