2025, Nov 03 03:00
Concave Polygon Ordering with a Guide Polyline: a Vectorized NumPy Method for Non-Intersecting Paths
Learn a robust Python and NumPy technique to order unordered boundary points into a non-self-intersecting polygon by leveraging an inner guide polyline.
Ordering an unordered set of boundary points into a single, non-self-intersecting polygon is easy when the shape is convex and hard when it isn’t. In this case you also have a guiding polyline drawn inside the point cloud. That line is the key: you can use it to define a local coordinate system for every point and derive a deterministic, reproducible traversal order that respects the inner constraint and visits every vertex once.
Problem setup
The following program sketches the scenario: a freehand stroke is smoothed into a polyline, normal samples are shown for inspection, and a separate set of draggable points represents the polygon vertices that need to be ordered. A naive angle sort doesn’t work for concave shapes, and a convex hull or an alpha shape may skip points. The goal is to find a robust ordering leveraging the inner line.
import pygame
import numpy as np
from scipy.interpolate import splprep, splev
import sys
import random
# Pygame setup
pygame.init()
surface = pygame.display.set_mode((1000, 1000))
pygame.display.set_caption("Smoothed Curve with Draggable Points")
ticker = pygame.time.Clock()
# Colors
COL_BG = (255, 255, 255)
COL_STROKE = (0, 0, 255)
COL_CURVE = (255, 0, 0)
COL_NORMAL = (0, 200, 0)
COL_POINT = (0, 0, 0)
is_tracing = False
stroke_pts = []
curve_pts = []
derivs = []
# Dot settings
PT_R = 6
PICK_R = 10
class Handle:
def __init__(self, pos):
self.x, self.y = pos
self._drag = False
def xy(self):
return (int(self.x), int(self.y))
def hit_test(self, mx, my):
return (self.x - mx)**2 + (self.y - my)**2 <= PICK_R**2
def begin_drag(self):
self._drag = True
def end_drag(self):
self._drag = False
def move_to(self, mx, my):
if self._drag:
self.x, self.y = mx, my
# crude test arrangement of dots.
anchors = [Handle(t) for t in [
(159, 459),(133, 193),(286, 481),(241, 345),(411, 404),
(280, 349),(352, 471),(395, 361),(85, 390),(203, 321),
(41, 281),(58, 348),(175, 275),(75, 185),(385, 443),
(44, 219),(148, 229),(215, 477),(338, 339),(122, 430)
]]
def take_every(seq, step=2):
return seq[::step]
def fit_curve(trace, smoothness=150):
if len(trace) < 4:
return [], []
trace = take_every(trace, step=2)
x, y = zip(*trace)
try:
tck, u = splprep([x, y], s=100.0, k=3)
u_new = np.linspace(0, 1, smoothness)
# Evaluate curve points
x_vals, y_vals = splev(u_new, tck)
# Derivatives: dx/du, dy/du
dx_du, dy_du = splev(u_new, tck, der=1)
slopes = np.array(dy_du) / np.array(dx_du)
return list(zip(x_vals, y_vals)), list(zip(dx_du, dy_du))
except:
return [], []
def render_normals(surf, curve, tangents, spacing=10, length=30):
for i in range(0, len(curve), spacing):
x, y = curve[i]
dx, dy = tangents[i]
nx, ny = -dy, dx
norm = np.hypot(nx, ny)
if norm == 0:
continue
nx /= norm
ny /= norm
x1 = x - nx * length / 2
y1 = y - ny * length / 2
x2 = x + nx * length / 2
y2 = y + ny * length / 2
pygame.draw.line(surf, COL_NORMAL, (x1, y1), (x2, y2), 2)
def spawn_anchors(n=20, w=800, h=600):
return [Handle((random.randint(0, w), random.randint(0, h))) for _ in range(n)]
# Main loop
active_handle = None
go = True
while go:
surface.fill(COL_BG)
mx, my = pygame.mouse.get_pos()
for event in pygame.event.get():
if event.type == pygame.QUIT:
go = False
elif event.type == pygame.MOUSEBUTTONDOWN:
is_tracing = True
stroke_pts = []
curve_pts = []
derivs = []
# Check for draggable handle
for h in anchors:
if h.hit_test(mx, my):
active_handle = h
h.begin_drag()
is_tracing = False
break
elif event.type == pygame.MOUSEBUTTONUP:
is_tracing = False
if active_handle:
active_handle.end_drag()
active_handle = None
else:
curve_pts, derivs = fit_curve(stroke_pts)
elif event.type == pygame.MOUSEMOTION:
if active_handle:
active_handle.move_to(mx, my)
elif is_tracing:
stroke_pts.append((mx, my))
elif event.type == pygame.KEYDOWN:
if event.key == pygame.K_p:
anchors = spawn_anchors()
if event.key == pygame.K_l:
for h in anchors: print(f"({h.x}, {h.y})")
# Draw anchors
for h in anchors:
pygame.draw.circle(surface, COL_POINT, h.xy(), PT_R)
# Draw raw stroke
if len(stroke_pts) > 1:
pygame.draw.lines(surface, COL_STROKE, False, stroke_pts, 1)
# Draw smoothed curve
if len(curve_pts) > 1:
pygame.draw.lines(surface, COL_CURVE, False, curve_pts, 3)
# Draw normals
if curve_pts and derivs:
render_normals(surface, curve_pts, derivs)
pygame.display.flip()
ticker.tick(60)
pygame.quit()
sys.exit()Why angle sorting and hulls fall short here
Angle-based strategies order points by their polar angle around a centroid or a chosen origin. That produces self-intersections for many concave shapes because the angular sweep can jump across the shape. Convex hulls, by definition, drop interior points. Concave hulls like alphashape depend on a parameter and may legally ignore points, which breaks the hard requirement that the polygon must pass through all given vertices. What changes the game is the presence of the inner line. If we know that polyline and the boundary must not cross it, we can linearize the problem along that line.
Core idea
Use the inner polyline to define a 2D local coordinate frame for each vertex. For every consecutive pair of points on the polyline, take the segment and its midpoint. For each boundary vertex, find the nearest segment midpoint and build a simple local frame aligned with that segment: one axis roughly parallel to the segment and the other roughly normal to it. Project the vertex into this frame to get two numbers: a longitudinal coordinate along the segment and a signed offset to indicate which side of the segment it lies on. Combine these into sort keys so that points on opposite sides of the polyline are traversed in opposite longitudinal directions, stitching everything into a single perimeter order. In practice, this benefits from decimating the inner polyline to a modest number of points.
Solution: vectorized ordering with NumPy
The implementation below takes two inputs: the unordered boundary vertices and the inner polyline. It returns the index order for the polygon together with diagnostic arrays that show which segment midpoint each vertex attached to.
import numpy as np
from matplotlib import pyplot as plt
from matplotlib.collections import LineCollection
def order_vertices(
verts: np.ndarray, guide: np.ndarray,
) -> tuple[
np.ndarray, # index order
np.ndarray, # segment midpoints
np.ndarray, # assigned midpoint per vertex
]:
# Consecutive segment pairs along the guide
seg_pairs = np.lib.stride_tricks.sliding_window_view(guide, window_shape=2, axis=0)
midpts = seg_pairs.mean(axis=-1)
# Local frames per segment: parallel and normal directions
frame = np.einsum(
'ijk,k,jmn->imn',
seg_pairs,
(-1e-4, +1e-4),
(
((1, 0),
(0, -1)),
((0, 1),
(1, 0)),
),
)
# Nearest segment midpoint for each vertex
delta = verts[:, np.newaxis] - midpts[np.newaxis]
near_idx = np.vecdot(delta, delta, axis=-1).argmin(axis=-1)
attach_midpts = midpts[near_idx]
# Projection into local (s, t) coordinates
s, t = np.einsum('ik,ijk->ji', verts - attach_midpts, frame[near_idx])
leftside = t > 0
# Sort by side, then by segment index, then by longitudinal coordinate.
# Flip directions on one side to make orientation consistent.
sort_keys = np.stack((s, near_idx, leftside))
sort_keys[:2, ~leftside] *= -1
return np.lexsort(sort_keys), midpts, attach_midpts
def visualize(
verts: np.ndarray, guide: np.ndarray, ordering: np.ndarray,
midpts: np.ndarray, attach_midpts: np.ndarray,
) -> plt.Figure:
fig, ax = plt.subplots()
ax.scatter(*verts.T)
ax.scatter(*midpts.T)
ax.plot(*guide.T)
ax.add_collection(LineCollection(
np.stack((attach_midpts, verts), axis=1)[ordering],
edgecolors='blue', alpha=0.2,
))
for o, (vx, vy) in enumerate(verts[ordering]):
ax.text(vx + 5, vy + 5, str(o))
return fig
def example_run() -> None:
# Sample data
verts = np.array((
(159, 459), (133, 193), (286, 481), (241, 345),
(411, 404), (280, 349), (352, 471), (395, 361),
( 85, 390), (203, 321), ( 41, 281), ( 58, 348),
(175, 275), ( 75, 185), (385, 443), ( 44, 219),
(148, 229), (215, 477), (338, 339), (122, 430),
))
guide = np.array((
( 94, 209),
(117, 329),
(200, 400),
(287, 419),
(375, 387),
))
ordering, midpts, attach_midpts = order_vertices(verts=verts, guide=guide)
visualize(
verts=verts, guide=guide, ordering=ordering,
midpts=midpts, attach_midpts=attach_midpts,
)
plt.show()
if __name__ == '__main__':
example_run()How this addresses the constraints
The approach does not rely on convexity and does not discard any vertices. By using the inner polyline as a scaffold, it removes the ambiguity that produces self-intersections in concave shapes. Each vertex is linked to a single guide segment via the nearest midpoint and is ordered consistently on either side of the guide. In practice, using a guide with a small number of points makes the segment assignment and sorting cleaner and more stable.
Why it’s worth knowing
When you have side information like an interior line that is correlated with the outer boundary, it can transform a seemingly combinatorial ordering problem into a straightforward coordinate transformation and sort. This pattern shows up in contour following, isocurve extraction, and path reconstruction tasks where you want reproducible, parameter-free behavior without dropping data.
Takeaways
If your data includes a polyline drawn inside the shape, align the ordering to that structure. Downsample the guide line to a manageable number of vertices, compute a per-segment local frame, project each boundary point to obtain longitudinal position and signed distance, and sort with a side-aware key. The result is a clean polygon ordering that stays consistent with the inner constraint and includes all points.
The article is based on a question from StackOverflow by Leo and an answer by Reinderien.