191 lines
6.3 KiB
Python
191 lines
6.3 KiB
Python
# Polygon helpers shared by the model builders: boundary tracing of a bitmap mask, simplification,
|
|
# offsetting, orientation and triangulation. Rings are Nx2 numpy arrays of (x, y) points.
|
|
import numpy as np
|
|
|
|
|
|
def traceBoundary(mask):
|
|
# Moore-neighbour tracing, clockwise on screen, starting from the top-left foreground pixel.
|
|
padded = np.pad(mask, 1)
|
|
ys, xs = np.nonzero(padded)
|
|
start = (ys[0], xs[0])
|
|
steps = [(0, 1), (1, 1), (1, 0), (1, -1), (0, -1), (-1, -1), (-1, 0), (-1, 1)]
|
|
boundary = [start]
|
|
current = start
|
|
came = 6
|
|
while True:
|
|
found = False
|
|
for k in range(8):
|
|
d = (came + 1 + k) % 8
|
|
ny = current[0] + steps[d][0]
|
|
nx = current[1] + steps[d][1]
|
|
if padded[ny, nx]:
|
|
came = (d + 4) % 8
|
|
# Next search starts just after the backtrack direction.
|
|
came = (d + 5) % 8
|
|
current = (ny, nx)
|
|
found = True
|
|
break
|
|
if not found or (current == start and len(boundary) > 1):
|
|
break
|
|
boundary.append(current)
|
|
if len(boundary) > 100000:
|
|
break
|
|
return np.array([(x - 1, y - 1) for y, x in boundary], dtype=float)
|
|
|
|
|
|
def simplify(points, tolerance):
|
|
# Douglas-Peucker on a closed ring, split at the two points farthest apart.
|
|
def dp(pts):
|
|
if len(pts) < 3:
|
|
return pts
|
|
a = pts[0]
|
|
b = pts[-1]
|
|
d = b - a
|
|
n = np.hypot(*d)
|
|
if n < 1e-9:
|
|
dist = np.hypot(*(pts - a).T)
|
|
else:
|
|
dist = np.abs(d[0] * (pts[:, 1] - a[1]) - d[1] * (pts[:, 0] - a[0])) / n
|
|
i = int(np.argmax(dist))
|
|
if dist[i] > tolerance:
|
|
return np.vstack([dp(pts[:i + 1])[:-1], dp(pts[i:])])
|
|
return np.array([a, b])
|
|
far = int(np.argmax(np.hypot(*(points - points[0]).T)))
|
|
first = dp(points[:far + 1])
|
|
second = dp(np.vstack([points[far:], points[:1]]))
|
|
ring = np.vstack([first[:-1], second[:-1]])
|
|
return ring
|
|
|
|
|
|
def signedArea(ring):
|
|
n = len(ring)
|
|
return 0.5 * sum(ring[i][0] * ring[(i + 1) % n][1] - ring[(i + 1) % n][0] * ring[i][1] for i in range(n))
|
|
|
|
|
|
def offsetRing(ring, amount):
|
|
# Grow (positive) or shrink a ring by a distance; mitres are clamped so acute tips do not spike.
|
|
n = len(ring)
|
|
sign = 1.0 if signedArea(ring) > 0 else -1.0
|
|
result = []
|
|
for i in range(n):
|
|
p = ring[i]
|
|
d0 = p - ring[i - 1]
|
|
d1 = ring[(i + 1) % n] - p
|
|
d0 = d0 / max(1e-9, np.hypot(*d0))
|
|
d1 = d1 / max(1e-9, np.hypot(*d1))
|
|
n0 = np.array([d0[1], -d0[0]]) * sign
|
|
n1 = np.array([d1[1], -d1[0]]) * sign
|
|
bis = n0 + n1
|
|
length = np.hypot(*bis)
|
|
if length < 1e-6:
|
|
result.append(p + n0 * amount)
|
|
continue
|
|
bis = bis / length
|
|
mitre = amount / max(0.35, bis @ n0)
|
|
result.append(p + bis * mitre)
|
|
return np.array(result)
|
|
|
|
|
|
def earClip(ring):
|
|
# Triangulates a simple polygon; returns index triples into ring.
|
|
n = len(ring)
|
|
indices = list(range(n))
|
|
if signedArea(ring) < 0:
|
|
indices.reverse()
|
|
triangles = []
|
|
|
|
def convex(a, b, c):
|
|
return (ring[b][0] - ring[a][0]) * (ring[c][1] - ring[a][1]) - (ring[b][1] - ring[a][1]) * (ring[c][0] - ring[a][0]) > 1e-9
|
|
|
|
def inside(p, a, b, c):
|
|
def side(u, v):
|
|
return (v[0] - u[0]) * (p[1] - u[1]) - (v[1] - u[1]) * (p[0] - u[0])
|
|
return side(ring[a], ring[b]) >= 0 and side(ring[b], ring[c]) >= 0 and side(ring[c], ring[a]) >= 0
|
|
|
|
guard = 0
|
|
while len(indices) > 3 and guard < 10 * n:
|
|
guard += 1
|
|
clipped = False
|
|
for k in range(len(indices)):
|
|
a = indices[k - 1]
|
|
b = indices[k]
|
|
c = indices[(k + 1) % len(indices)]
|
|
if not convex(a, b, c):
|
|
continue
|
|
if any(inside(ring[m], a, b, c) for m in indices if m not in (a, b, c)):
|
|
continue
|
|
triangles.append((a, b, c))
|
|
del indices[k]
|
|
clipped = True
|
|
break
|
|
if not clipped:
|
|
break
|
|
if len(indices) == 3:
|
|
triangles.append(tuple(indices))
|
|
return triangles
|
|
|
|
|
|
def isConvex(ring):
|
|
n = len(ring)
|
|
sign = 0.0
|
|
for i in range(n):
|
|
a = ring[i]
|
|
b = ring[(i + 1) % n]
|
|
c = ring[(i + 2) % n]
|
|
cross = (b[0] - a[0]) * (c[1] - b[1]) - (b[1] - a[1]) * (c[0] - b[0])
|
|
if abs(cross) < 1e-9:
|
|
continue
|
|
if sign == 0.0:
|
|
sign = np.sign(cross)
|
|
elif np.sign(cross) != sign:
|
|
return False
|
|
return sign != 0.0
|
|
|
|
|
|
def selfIntersects(ring):
|
|
# True when any two non-adjacent edges of the ring cross.
|
|
n = len(ring)
|
|
for i in range(n):
|
|
a0 = ring[i]
|
|
a1 = ring[(i + 1) % n]
|
|
for j in range(i + 2, n):
|
|
if i == 0 and j == n - 1:
|
|
continue
|
|
b0 = ring[j]
|
|
b1 = ring[(j + 1) % n]
|
|
d1 = (a1[0] - a0[0]) * (b0[1] - a0[1]) - (a1[1] - a0[1]) * (b0[0] - a0[0])
|
|
d2 = (a1[0] - a0[0]) * (b1[1] - a0[1]) - (a1[1] - a0[1]) * (b1[0] - a0[0])
|
|
d3 = (b1[0] - b0[0]) * (a0[1] - b0[1]) - (b1[1] - b0[1]) * (a0[0] - b0[0])
|
|
d4 = (b1[0] - b0[0]) * (a1[1] - b0[1]) - (b1[1] - b0[1]) * (a1[0] - b0[0])
|
|
if ((d1 > 0) != (d2 > 0)) and ((d3 > 0) != (d4 > 0)):
|
|
return True
|
|
return False
|
|
|
|
|
|
def removeSpikes(ring, degrees=8.0):
|
|
# Drops vertices where the outline turns back on itself, the hair-thin fingers a raster boundary
|
|
# leaves at one-pixel spurs.
|
|
limit = -np.cos(np.radians(degrees))
|
|
changed = True
|
|
while changed and len(ring) > 3:
|
|
changed = False
|
|
n = len(ring)
|
|
for i in range(n):
|
|
d0 = ring[i] - ring[i - 1]
|
|
d1 = ring[(i + 1) % n] - ring[i]
|
|
l0 = np.hypot(*d0)
|
|
l1 = np.hypot(*d1)
|
|
if l0 < 1e-9 or l1 < 1e-9 or (d0 @ d1) / (l0 * l1) < limit:
|
|
ring = np.delete(ring, i, axis=0)
|
|
changed = True
|
|
break
|
|
return ring
|
|
|
|
|
|
def simplifyClean(points, tolerance):
|
|
# Simplify, drop spikes, and back off the tolerance until the ring does not cross itself.
|
|
while True:
|
|
ring = removeSpikes(simplify(points, tolerance))
|
|
if not selfIntersects(ring) or tolerance < 0.2:
|
|
return ring
|
|
tolerance = tolerance / 2
|