Add PolygonShape.split_to_convex (polygon to convex tiles)

Ports geometry/planar/Polygon.java (corner de-duplication and collinear
removal, winding number) and PolygonShape.java, including the recursive
split_to_convex that decomposes a simple polygon into convex Simplex
tiles by dividing at concave corners along minimal axis-parallel lines.

Orientation and convexity tests are exact; the division-point search is
approximate (float line evaluation, split point rounded to an integer
corner) as upstream. The concave-corner search starts deterministically
at corner 0 rather than a seeded PRNG; this only affects which valid
decomposition is produced.

Supporting additions: Simplex.from_corners (convex polygon to simplex)
and Line.function_value_approx / function_in_y_value_approx.

Invariant tests over L, plus, staircase and square polygons: tile areas
sum to the polygon area (no gaps, no overlap), a point is in the polygon
iff in some tile, and no point is strictly inside more than one tile
(interiors disjoint). Also covers Polygon normalization and orientation.
This commit is contained in:
Ryan Malloy 2026-07-11 18:41:07 -06:00
parent 1599d181b5
commit 68ab16a7cf
6 changed files with 525 additions and 0 deletions

View File

@ -22,6 +22,8 @@ from .float_point import FloatPoint
from .limits import CRIT_DOUBLE, CRIT_INT from .limits import CRIT_DOUBLE, CRIT_INT
from .line import Line from .line import Line
from .point import IntPoint, Point, RationalPoint, point, rational_point from .point import IntPoint, Point, RationalPoint, point, rational_point
from .polygon import Polygon
from .polygon_shape import PolygonShape
from .side import Side, Signum from .side import Side, Signum
from .simplex import Simplex from .simplex import Simplex
from .tile import TileShape from .tile import TileShape
@ -47,4 +49,6 @@ __all__ = [
"IntBox", "IntBox",
"TileShape", "TileShape",
"Simplex", "Simplex",
"Polygon",
"PolygonShape",
] ]

View File

@ -198,6 +198,28 @@ class Line:
result = self.side_of(line_1.intersection(line_2)) result = self.side_of(line_1.intersection(line_2))
return result return result
def function_value_approx(self, x: float) -> float:
"""Approximate ``y`` on this (non-vertical) line at abscissa ``x``."""
p1 = self.a.to_float()
p2 = self.b.to_float()
dx = p2.x - p1.x
if dx == 0:
return 0.0
dy = p2.y - p1.y
det = p1.x * p2.y - p2.x * p1.y
return (dy * x - det) / dx
def function_in_y_value_approx(self, y: float) -> float:
"""Approximate ``x`` on this (non-horizontal) line at ordinate ``y``."""
p1 = self.a.to_float()
p2 = self.b.to_float()
dy = p2.y - p1.y
if dy == 0:
return 0.0
dx = p2.x - p1.x
det = p1.x * p2.y - p2.x * p1.y
return (dx * y + det) / dy
# --- intersection (exact) ---------------------------------------------- # --- intersection (exact) ----------------------------------------------
def intersection(self, other: Line) -> Point: def intersection(self, other: Line) -> Point:

View File

@ -0,0 +1,77 @@
"""``Polygon`` — a normalized corner list (no duplicate or collinear corners).
Ports ``geometry/planar/Polygon.java``. Used as the corner normalizer behind
:class:`~freeroute.geometry.polygon_shape.PolygonShape`.
**Arithmetic model:** corner de-duplication and collinearity removal are
**exact** (exact ``Point.side_of``); the winding-number test is approximate
(sum of ``angle_approx``), exactly as upstream.
"""
from __future__ import annotations
import math
from .point import Point
from .side import Side
__all__ = ["Polygon"]
class Polygon:
"""A list of corners with no two consecutive equal and no three collinear."""
__slots__ = ("corners",)
def __init__(self, points) -> None:
corners: list[Point] = list(points)
corner_removed = True
while corner_removed:
corner_removed = False
if not corners:
break
# remove consecutive duplicate points
deduped: list[Point] = [corners[0]]
for pt in corners[1:]:
if pt != deduped[-1]:
deduped.append(pt)
if len(deduped) != len(corners):
corner_removed = True
corners = deduped
# remove a point collinear with its previous and next point
if len(corners) >= 3:
for i in range(1, len(corners) - 1):
if corners[i].side_of(corners[i - 1], corners[i + 1]) == Side.COLLINEAR:
del corners[i]
corner_removed = True
break
self.corners: list[Point] = corners
def corner_array(self) -> list[Point]:
return list(self.corners)
def revert_corners(self) -> Polygon:
return Polygon(list(reversed(self.corners)))
def winding_number_after_closing(self) -> int:
"""Winding number of the closed polygon: >0 counter-clockwise, <0 clockwise."""
arr = self.corners
if len(arr) < 2:
return 0
first_side = arr[1].difference_by(arr[0])
prev_side = first_side
corner_count = len(arr)
if arr[0] == arr[corner_count - 1]:
corner_count -= 1
angle_sum = 0.0
for i in range(1, corner_count + 1):
if i == corner_count - 1:
next_side = arr[0].difference_by(arr[i])
elif i == corner_count:
next_side = first_side
else:
next_side = arr[i + 1].difference_by(arr[i])
angle_sum += prev_side.angle_approx(next_side)
prev_side = next_side
angle_sum /= 2.0 * math.pi
return round(angle_sum)

View File

@ -0,0 +1,262 @@
"""``PolygonShape`` — a simple (possibly concave) polygon and ``split_to_convex``.
Ports ``geometry/planar/PolygonShape.java``, in particular the recursive
``split_to_convex`` that decomposes a simple polygon into convex
:class:`~freeroute.geometry.simplex.Simplex` tiles the decomposition the
router's search tree relies on.
**Arithmetic model:** corner orientation / convexity tests are **exact**
(exact ``Point.side_of``); the division-point search is **approximate** (it uses
``float`` line evaluations and rounds the split point to an integer corner, as
upstream). Point containment is exact and delegated to the convex tiles.
The recursion starts the concave-corner search at corner 0 (deterministic);
upstream seeds a fixed PRNG to vary the start. The choice only affects *which*
valid decomposition is produced, never its correctness.
"""
from __future__ import annotations
from .float_point import FloatPoint
from .line import Line
from .point import Point
from .polygon import Polygon
from .side import Side
from .simplex import Simplex
__all__ = ["PolygonShape"]
class PolygonShape:
"""A simple polygon, normalized to a canonical counter-clockwise corner list."""
__slots__ = ("corners", "_convex_pieces")
def __init__(self, points) -> None:
poly = Polygon(points)
if poly.winding_number_after_closing() < 0:
poly = poly.revert_corners() # make counter-clockwise
curr = poly.corner_array()
last = len(curr) - 1
if last > 0 and curr[0] == curr[last]:
last -= 1
if last >= 2 and curr[last].side_of(curr[last - 1], curr[0]) == Side.COLLINEAR:
last -= 1
first = 0
if last - first >= 2 and curr[0].side_of(curr[1], curr[last]) == Side.COLLINEAR:
first += 1
# rotate so the corner with the lowest y (then lowest x) comes first
start = first
start_c = curr[start].to_float()
for i in range(start + 1, last + 1):
c = curr[i].to_float()
if c.y < start_c.y or (c.y == start_c.y and c.x < start_c.x):
start = i
start_c = c
self.corners: list[Point] = curr[start : last + 1] + curr[first:start]
self._convex_pieces: list[Simplex] | None = None
def __repr__(self) -> str:
return f"PolygonShape(<{len(self.corners)} corners>)"
# --- accessors ----------------------------------------------------------
def corner(self, no: int) -> Point:
return self.corners[no]
def border_line_count(self) -> int:
return len(self.corners)
def border_line(self, no: int) -> Line:
nxt = self.corners[0] if no == len(self.corners) - 1 else self.corners[no + 1]
return Line(self.corners[no], nxt)
def dimension(self) -> int:
n = len(self.corners)
if n == 0:
return -1
if n == 1:
return 0
if n == 2:
return 1
return 2
def is_empty(self) -> bool:
return len(self.corners) == 0
def area(self) -> float:
if self.dimension() <= 2 and len(self.corners) < 3:
return 0.0
result = 0.0
corners = [c.to_float() for c in self.corners]
n = len(corners)
prev_corner = corners[n - 2]
curr_corner = corners[n - 1]
for i in range(n):
next_corner = corners[i]
result += curr_corner.x * (next_corner.y - prev_corner.y)
prev_corner = curr_corner
curr_corner = next_corner
return 0.5 * abs(result)
# --- convexity ----------------------------------------------------------
def is_convex(self) -> bool:
corners = self.corners
n = len(corners)
if n <= 2:
return True
prev_point = corners[n - 1]
curr_point = corners[0]
next_point = corners[1]
for ind in range(n):
if next_point.side_of(prev_point, curr_point) == Side.ON_THE_RIGHT:
return False
prev_point = curr_point
curr_point = next_point
next_point = corners[0] if ind == n - 2 else corners[(ind + 2) % n]
return True
# --- containment (via convex decomposition) ----------------------------
def is_outside(self, point: Point) -> bool:
return all(piece.is_outside(point) for piece in self.split_to_convex())
def contains(self, point: Point) -> bool:
return not self.is_outside(point)
def contains_inside(self, point: Point) -> bool:
# strictly inside some tile and not merely on a shared division edge
pieces = self.split_to_convex()
return any(p.contains_inside(point) for p in pieces)
# --- the key operation: split into convex tiles ------------------------
def split_to_convex(self) -> list[Simplex]:
if self._convex_pieces is None:
pieces = _split_recu(self.corners)
if pieces is None:
self._convex_pieces = []
else:
self._convex_pieces = [Simplex.from_corners(p) for p in pieces]
return self._convex_pieces
def _split_recu(corners: list[Point]) -> list[list[Point]] | None:
"""Recursively divide a simple polygon into convex corner lists."""
n = len(corners)
if n < 3:
return [corners]
# find the first concave corner
concave_no = -1
for start in range(n):
prev_corner = corners[start - 1]
curr_corner = corners[start]
next_corner = corners[(start + 1) % n]
if next_corner.side_of(prev_corner, curr_corner) == Side.ON_THE_RIGHT:
concave_no = start
break
if concave_no < 0:
return [corners] # already convex
projection, corner_no_after = _division_point(corners, concave_no)
if projection is None:
return None # self-intersecting or degenerate
split_pt = projection.round()
# first piece: concave_no .. corner_no_after (exclusive) + split point
count = corner_no_after - concave_no
if count < 0:
count += n
count += 1
first_arr: list[Point] = []
idx = concave_no
for _ in range(count - 1):
first_arr.append(corners[idx])
idx = (idx + 1) % n
first_arr.append(split_pt)
# last piece: split point + corner_no_after .. concave_no
count = concave_no - corner_no_after
if count < 0:
count += n
count += 2
last_arr: list[Point] = [split_pt]
idx = corner_no_after
for _ in range(1, count):
last_arr.append(corners[idx])
idx = (idx + 1) % n
c1 = _split_recu(PolygonShape(first_arr).corners)
if c1 is None:
return None
c2 = _split_recu(PolygonShape(last_arr).corners)
if c2 is None:
return None
return c1 + c2
def _division_point(corners: list[Point], concave_no: int) -> tuple[FloatPoint | None, int]:
"""Minimal axis-parallel division line from a concave corner (approximate)."""
n = len(corners)
concave = corners[concave_no].to_float()
before = corners[concave_no - 1].to_float()
after = corners[(concave_no + 1) % n].to_float()
search_right = before.y > concave.y or concave.y > after.y
search_left = before.y < concave.y or concave.y < after.y
search_up = before.x < concave.x or concave.x < after.x
search_down = before.x > concave.x or concave.x > after.x
min_dist = float("inf")
min_projection: FloatPoint | None = None
corner_no_after_min = 0
idx_after = (concave_no + 2) % n
corner_before = corners[idx_after - 1] if idx_after != 0 else corners[n - 1]
before_approx = corner_before.to_float()
for _ in range(n - 2):
corner_after = corners[idx_after]
after_approx = corner_after.to_float()
if before_approx.y != after_approx.y: # horizontal division
lo_y = min(before_approx.y, after_approx.y)
hi_y = max(before_approx.y, after_approx.y)
if lo_y <= concave.y <= hi_y:
line = Line(corner_before, corner_after)
x_intersect = line.function_in_y_value_approx(concave.y)
dist = abs(x_intersect - concave.x)
ok = dist < min_dist and (
(search_right and x_intersect > concave.x and concave.y <= after_approx.y)
or (search_left and x_intersect < concave.x and concave.y >= after_approx.y)
)
if ok:
min_dist = dist
corner_no_after_min = idx_after
min_projection = FloatPoint(x_intersect, concave.y)
if before_approx.x != after_approx.x: # vertical division
lo_x = min(before_approx.x, after_approx.x)
hi_x = max(before_approx.x, after_approx.x)
if lo_x <= concave.x <= hi_x:
line = Line(corner_before, corner_after)
y_intersect = line.function_value_approx(concave.x)
dist = abs(y_intersect - concave.y)
ok = dist < min_dist and (
(search_up and y_intersect > concave.y and concave.x >= after_approx.x)
or (search_down and y_intersect < concave.y and concave.x <= after_approx.x)
)
if ok:
min_dist = dist
corner_no_after_min = idx_after
min_projection = FloatPoint(concave.x, y_intersect)
corner_before = corner_after
before_approx = after_approx
idx_after = 0 if idx_after == n - 1 else idx_after + 1
return min_projection, corner_no_after_min

View File

@ -52,6 +52,18 @@ class Simplex(TileShape):
arr.sort() arr.sort()
return Simplex(arr)._remove_redundant_lines() return Simplex(arr)._remove_redundant_lines()
@staticmethod
def from_corners(corners) -> Simplex:
"""Build a normalized simplex from the corners of a convex polygon.
Ports ``TileShape.get_instance(Point[])``: one directed border line per
edge (corner ``i`` to ``i+1``), then normalized.
"""
corners = list(corners)
n = len(corners)
lines = [Line(corners[j], corners[(j + 1) % n]) for j in range(n)]
return Simplex.get_instance(lines)
@staticmethod @staticmethod
def empty() -> Simplex: def empty() -> Simplex:
return Simplex(()) return Simplex(())

View File

@ -0,0 +1,148 @@
"""Invariant tests for Polygon / PolygonShape.split_to_convex.
No external oracle exists, so these assert the decomposition invariants from the
phase-2 brief for several concave polygons: the convex tiles' areas sum to the
polygon area (no gaps, no overlap), a point is contained in the polygon iff it
is contained in some tile, and no point is strictly inside more than one tile
(interiors disjoint).
Source: ``geometry/planar/{Polygon,PolygonShape}.java``.
"""
from __future__ import annotations
import pytest
from freeroute.geometry import IntPoint, Side
from freeroute.geometry.polygon import Polygon
from freeroute.geometry.polygon_shape import PolygonShape
def pts(*coords):
return [IntPoint(x, y) for x, y in coords]
# --- Polygon normalization --------------------------------------------------
def test_polygon_removes_consecutive_duplicates():
poly = Polygon(pts((0, 0), (0, 0), (10, 0), (10, 10), (10, 10), (0, 10)))
assert len(poly.corners) == 4
def test_polygon_removes_middle_collinear_corner():
# (5, 0) is collinear with its neighbours (0,0) and (10,0) -> dropped.
# Raw Polygon only removes *middle* collinear corners (endpoint-collinear
# removal is PolygonShape's job), so the trailing (0,5) is kept.
poly = Polygon(pts((0, 0), (5, 0), (10, 0), (10, 10), (0, 10), (0, 5)))
coords = [(c.x, c.y) for c in poly.corners]
assert (5, 0) not in coords
assert coords == [(0, 0), (10, 0), (10, 10), (0, 10), (0, 5)]
def test_polygon_shape_removes_endpoint_collinear_corner():
# PolygonShape additionally drops the endpoint-collinear (0,5).
ps = PolygonShape(pts((0, 0), (5, 0), (10, 0), (10, 10), (0, 10), (0, 5)))
coords = {(c.x, c.y) for c in ps.corners}
assert coords == {(0, 0), (10, 0), (10, 10), (0, 10)}
def test_winding_number_sign():
ccw = Polygon(pts((0, 0), (10, 0), (10, 10), (0, 10)))
cw = Polygon(pts((0, 0), (0, 10), (10, 10), (10, 0)))
assert ccw.winding_number_after_closing() > 0
assert cw.winding_number_after_closing() < 0
# --- convexity --------------------------------------------------------------
def test_is_convex():
square = PolygonShape(pts((0, 0), (10, 0), (10, 10), (0, 10)))
assert square.is_convex()
ell = PolygonShape(pts((0, 0), (60, 0), (60, 20), (20, 20), (20, 60), (0, 60)))
assert not ell.is_convex()
# --- split_to_convex invariants ---------------------------------------------
SHAPES = {
"L": pts((0, 0), (60, 0), (60, 20), (20, 20), (20, 60), (0, 60)),
"plus": pts(
(20, 0),
(40, 0),
(40, 20),
(60, 20),
(60, 40),
(40, 40),
(40, 60),
(20, 60),
(20, 40),
(0, 40),
(0, 20),
(20, 20),
),
"staircase": pts((0, 0), (30, 0), (30, 10), (20, 10), (20, 20), (10, 20), (10, 30), (0, 30)),
"square": pts((0, 0), (40, 0), (40, 40), (0, 40)),
}
def _bounds(corners):
xs = [c.x for c in corners]
ys = [c.y for c in corners]
return min(xs) - 3, max(xs) + 3, min(ys) - 3, max(ys) + 3
@pytest.mark.parametrize("name", list(SHAPES))
def test_split_area_conserved(name):
ps = PolygonShape(SHAPES[name])
tiles = ps.split_to_convex()
assert tiles, "expected at least one convex tile"
assert sum(t.area() for t in tiles) == pytest.approx(ps.area())
@pytest.mark.parametrize("name", list(SHAPES))
def test_split_union_equals_polygon(name):
ps = PolygonShape(SHAPES[name])
tiles = ps.split_to_convex()
lo_x, hi_x, lo_y, hi_y = _bounds(ps.corners)
for x in range(lo_x, hi_x + 1):
for y in range(lo_y, hi_y + 1):
p = IntPoint(x, y)
in_poly = ps.contains(p)
in_some_tile = any(t.contains(p) for t in tiles)
assert in_poly == in_some_tile, (x, y)
@pytest.mark.parametrize("name", list(SHAPES))
def test_split_interiors_are_disjoint(name):
ps = PolygonShape(SHAPES[name])
tiles = ps.split_to_convex()
lo_x, hi_x, lo_y, hi_y = _bounds(ps.corners)
for x in range(lo_x, hi_x + 1):
for y in range(lo_y, hi_y + 1):
p = IntPoint(x, y)
strict = sum(1 for t in tiles if t.contains_inside(p))
assert strict <= 1, (x, y, strict)
def test_convex_polygon_is_single_tile():
ps = PolygonShape(SHAPES["square"])
tiles = ps.split_to_convex()
assert len(tiles) == 1
assert tiles[0].dimension() == 2
def test_concave_corner_detected_by_side_of():
# the reflex corner of the L is (20, 20): its next corner is on the right
ell = PolygonShape(SHAPES["L"])
corners = ell.corners
reflex_found = False
n = len(corners)
for i in range(n):
prev_c = corners[i - 1]
curr = corners[i]
nxt = corners[(i + 1) % n]
if nxt.side_of(prev_c, curr) == Side.ON_THE_RIGHT:
reflex_found = True
assert reflex_found