-
Notifications
You must be signed in to change notification settings - Fork 0
BCRS
Wolren edited this page Apr 30, 2026
·
4 revisions
BCRS computes the largest axis-aligned rectangle in an arbitrarily-oriented polygon using a two-phase approach:
- Vertex-coordinate grid solve (Stage 4)
- SDF-guided expansion (Stage 5)
- SDF-based certification (Stage 6)
flowchart TD
A[Stage 1: Geometry preparation] --> B[Stage 2: Heuristic candidates]
B --> C[Stage 3: Angle refinement]
C --> D[Stage 4: BCRS boundary-coordinate solve]
D --> E[Stage 5: SDF-guided expansion]
E --> F[Stage 6: SDF-based containment certification]
F --> G{Certification passes?}
G -- Yes --> H[Stage 7: Select best + output]
G -- No --> I{Fallback enabled?}
I -- Yes --> J[Best-effort symmetric shrink]
I -- No --> K[No rectangle]
J --> H
| Variant | Description | TOP_K |
|---|---|---|
| Standard | Full BCRS workflow | 3 |
| Fast | Trial ranking, reduced expensive runs | 3 |
Theorem [1]: The largest axis-aligned rectangle inscribed in a simple polygon always has at least two sides aligned with polygon vertex x- or y-coordinates.
Implication: Searching on the vertex-coordinate grid guarantees finding the global optimum.
Corollary: If the polygon has >300 unique coordinates on either axis, the vertex-coordinate approach becomes O(n²) in grid cells, making uniform grid seed preferable.
- Validate polygon: check
is_validandis_empty - If multipart: extract largest component by area
- Extract coordinate bounds:
minx, miny, maxx, maxy
- Extract all polygon edges from exterior and interiors
- Compute edge orientation:
atan2(dy, dx)for each edge - Map each orientation to
mod 90to handle symmetry - Bin orientations into histogram (unique edge angles)
- Compute convex hull upper bound on max area
- Prune angles with bound below current best × 0.90
- Run coarse uniform grid (40×40) at each remaining angle
- Sort candidates by area, keep top-K
For each candidate angle θ from Stage 2:
- Set search bounds:
[θ - 3°, θ + 3°] - Objective function: run grid solve at angle a, return area
- Apply bounded Brent minimization to find local maximum
- Refine grid at polished angle (120×120)
- Evaluate at polished angle + ±0.5° delta
flowchart TD
A[Extract vertex X coordinates] --> B[Extract vertex Y coordinates]
B --> C[Add min or max bounds]
C --> D[Sort unique to X, Y arrays]
D --> E{len X > 300 or len Y > 300?}
E -- Yes --> F[Abort BCRS, use uniform grid seed]
E -- No --> G[Build occupancy mask at cell centers]
G --> H[For each row: compute heights cumulative mask]
H --> I[Stack-based LRH on heights array]
I --> J[For each span compute width, height, depth]
J --> K[Return rectangle with maximum area]
This stage finds the maximum-area axis-aligned rectangle at vertex-coordinate precision.
Input: polygon P
Output: sorted unique X, Y arrays
1. Collect all vertex x-coordinates from exterior
2. Collect all vertex x-coordinates from each interior ring
3. Add minx, maxx to the set
4. Remove duplicates, sort → X = [x0, x1, ..., xm]
5. Repeat for y-coordinates → Y = [y0, y1, ..., yn]
If len(X) > 300 or len(Y) > 300: abort BCRS, use uniform grid seed instead.
For each grid cell [Xi, Xi+1] × [Yj, Yj+1]:
Compute cell center: ((Xi + Xi+1)/2, (Yj + Yj+1)/2)
Test if center is inside polygon using shapely prep.covers
Set mask[j][i] = 1 if inside, 0 otherwise
For each row j from 0 to n-1:
1. For each column i:
heights[i] += mask[j][i] // accumulate
heights[i] *= mask[j][i] // zero out if not occupied
2. Stack-based largest rectangle:
Initialize empty stack
For i from 0 to n_cols:
h = heights[i] (or 0 if i == n_cols)
start = i
While stack not empty and stack.top.height > h:
pop top → (start_idx, height)
width = X[i+1] - X[start_idx]
depth = Y[j+1] - Y[j - height + 1] // if height > 0
area = width × depth
update best rectangle if area larger
start = start_idx
push (start, h) onto stack
The key difference from standard LRH: width = X[i+1] - X[start_idx] uses variable pitch from actual vertex coordinates, not a uniform Δx.
best_rect = box(x0, y0, x1, y1) where bounds are from best area
best_area = (x1 - x0) × (y1 - y0)
flowchart TD
A[Query SDF at left side midpoint] --> B[Compute d_max]
B --> C{SDF less than 0?}
C -- Yes --> D[d_max equals min x0-minx or absolute SDF]
C -- No --> E[d_max equals x0 minus minx]
D --> F[Binary search 10 steps]
E --> F
F --> G{covers expanded rect?}
G -- Yes --> H[lo equals mid]
G -- No --> I[hi equals mid]
H --> J[Update left side x0]
I --> J
J --> K[Repeat for right bottom top]
K --> L[3 coordinate ascent iterations]
L --> M[Return expanded bounds]
This stage expands the rectangle outward while maintaining containment.
SDF(x, y) returns:
- positive distance if point outside polygon
- 0 if point on boundary
- negative distance (inside) to nearest boundary (exterior or hole)
Implementation:
polygon_sdf(point):
d = polygon.distance(point)
if d > 0: return d // outside
d_ext = polygon.exterior.distance(point)
if polygon.contains(point):
min_d = d_ext
for hole in polygon.interiors:
min_d = min(min_d, hole.distance(point))
return -min_d // inside, negative
if d_ext < 1e-12: return 0 // on exterior
for hole in polygon.interiors:
if hole.contains(point): return hole.exterior.distance(point)
return 0
For the left side (x0):
1. Query midpoint: s = SDF(x0, (y0 + y1)/2)
2. Compute upper bound:
If s < 0 (inside): d_max = min(x0 - minx, abs(s))
Else (outside/on boundary): d_max = x0 - minx
3. Binary search for maximum d:
lo = 0, hi = d_max
For step in 1 to 10:
mid = (lo + hi) / 2
if covers(P, box(x0 - mid, y0, x1, y1)):
lo = mid
else:
hi = mid
d_left = lo
Repeat for right side (x1), bottom (y0), top (y1).
For iteration in 1 to 3:
Expand left side (compute d_left)
Expand right side (compute d_right)
Expand bottom (compute d_bottom)
Expand top (compute d_top)
If no side expanded: break
The expansion is bounded by SDF, not by vertex coordinates. This provides a tighter bound because SDF uses actual distance to boundary.
flowchart TD
A[SDF at corner 1] --> B[SDF at corners 2,3,4]
B --> C[SDF at 4 edge midpoints]
C --> D[SDF_max equals maximum of 8 values]
D --> E{SDF_max less than or equal to epsilon?}
E -- Yes --> F[Rectangle certified - return]
E -- No --> G{shrink less than or equal to 20 percent of min side?}
G -- No --> H[Reject - return None]
G -- Yes --> I[Apply symmetric shrink]
I --> J{SDF_max after shrink less than or equal to epsilon?}
J -- Yes --> K[Return shrunk rectangle]
J -- No --> L[Reject - return None]
SDF_max = -infinity
For each corner:
SDF_max = max(SDF_max, SDF(corner))
For each edge midpoint:
SDF_max = max(SDF_max, SDF(midpoint))
Eight SDF evaluations total (4 corners + 4 midpoints).
If SDF_max <= epsilon (1e-7):
Rectangle is certified
Return rectangle as-is
Else:
// Attempt symmetric shrink
shrink = SDF_max + 2*epsilon
new_half_width = a - shrink
new_half_height = b - shrink
If new_half_width <= 0 or new_half_height <= 0:
Return None (cannot certify)
If shrink > min(a, b) × 0.20:
Return None (excessive shrink)
If max_ratio > 0 and new_half_width/new_half_height > max_ratio:
new_half_width = new_half_height × max_ratio
Rebuild rectangle with shrunk dimensions
Re-verify SDF_max <= epsilon
If passes: return shrunk rectangle
Else: return None
- For each angle candidate, keep best certified rectangle
- Sort by area descending
- Return highest-area rectangle with diagnostics:
- area, angle_deg, ratio
- cand_rank (0 = best)
- s2_gain (Stage 3 improvement)
- s4_gain (BCRS over coarse grid)
- s5_gain (SDF expansion)
- best_effort (1 if shrink fallback used)
| Parameter | Default | Description |
|---|---|---|
| GRID_COARSE | 40 | Coarse uniform grid size |
| GRID_FINE | 120 | Vertex-coordinate grid fallback size |
| TOP_K | 3 | Number of angle candidates for refinement |
| MAX_RATIO | 0 | Max aspect ratio (0 = unlimited) |
| ALWAYS_RETURN | True | Enable best-effort shrink fallback |
| ANGLE_STEP | 5.0 | Fallback sweep angle step |
| Variant | TOP_K | Time @290 (s) | Time @5406 (s) |
|---|---|---|---|
| Standard | 3 | 42.35 | TBD |
| Fast | 3 | 23.61 | 445.01 |
Best mode: N_WORKERS=1 (parallelization overhead exceeds per-feature time).
| Field | Description |
|---|---|
| area | Rectangle area in CRS units |
| angle_deg | Rotation angle in degrees (0-90) |
| ratio | Aspect ratio (long:short) |
| cand_rank | Rank of selected candidate (0 = best) |
| s2_gain | Area gain from Stage 3 (angle polishing) |
| s4_gain | Area gain from Stage 4 (BCRS over coarse) |
| s5_gain | Area gain from Stage 5 (SDF expansion) |
| best_effort | 1 if shrink fallback used |
[1] Daniels, K., Milenkovic, V., & Roel, G. (1997). Finding the largest axis-aligned rectangle in a polygon. Proceedings of the 13th Canadian Conference on Computational Geometry (CCCG). https://cgm.cs.mcgill.ca/~athena/workshop/2.html
[2] Brent, R.P. (1973). Algorithms for Finding Zeros and Extrema of Functions Without Calculating Derivatives. Prentice-Hall.
[3] Preparata, F.P., & Shamos, M.I. (1985). Computational Geometry: An Introduction. Springer-Verlag.
Navigation
- Home | Installation | Algorithms | Approximation | Skeleton | Axis-Aligned | Performance | Complexity | Foundations | Parameters | Folder Layout | Usage