Skip to content
Draft
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
30 changes: 29 additions & 1 deletion polygon_maker.py
Original file line number Diff line number Diff line change
Expand Up @@ -22,6 +22,19 @@
# on where the ring happens to start
SIMPLIFY_TOLERANCE_CELLS = 1.0

# boundary smoothing (Chaikin corner cutting via QgsGeometry.smooth).
# The smoothing offset is a fraction of each segment's length, so the
# ring is densified to at most DENSIFY_INTERVAL_CELLS-long segments
# first: this caps the corner rounding radius at roughly
# DENSIFY_INTERVAL_CELLS * SMOOTH_OFFSET cells — about the scale of the
# pixel jaggedness — instead of letting long straight edges produce
# arbitrarily large rounded corners. Two iterations turn each chamfer
# into a reasonably smooth arc; more iterations only add vertices for
# little visual gain.
DENSIFY_INTERVAL_CELLS = 3.0
SMOOTH_ITERATIONS = 2
SMOOTH_OFFSET = 0.25


class PixelGrid:
"""Pixel<->map transform context for PolygonMaker.
Expand Down Expand Up @@ -135,7 +148,7 @@ def build_polygons(self, crs: QgsCoordinateReferenceSystem) -> list[QgsFeature]:
"OUTPUT": "memory:",
},
)["OUTPUT"]
return list(cleaned_layer.getFeatures())
return self.smooth_features(list(cleaned_layer.getFeatures()), cell_size)

def make_polygons(
self, crs: QgsCoordinateReferenceSystem, layer_id: str | None = None
Expand Down Expand Up @@ -218,6 +231,21 @@ def make_layer_by(
features_layer.updateExtents()
return features_layer

def smooth_features(
self, features: list[QgsFeature], cell_size: float
) -> list[QgsFeature]:
"""Round off the residual pixel jaggedness of the boundaries.

Works directly on the geometries (densify + Chaikin smoothing)
instead of going through processing algorithms, to keep the
per-click cost low."""
for feature in features:
geometry = feature.geometry().densifyByDistance(
cell_size * DENSIFY_INTERVAL_CELLS
)
feature.setGeometry(geometry.smooth(SMOOTH_ITERATIONS, SMOOTH_OFFSET))
return features

def noise_reduction(self, features, noise_multiply: float) -> list[QgsFeature]:
"""Drop features smaller than `noise_multiply` mask cells."""
return [
Expand Down
46 changes: 35 additions & 11 deletions tests/test_polygon_maker.py
Original file line number Diff line number Diff line change
Expand Up @@ -99,7 +99,8 @@ def test_returns_features_without_touching_project(
features = maker.build_polygons(crs=CRS)

assert len(features) == 1
assert features[0].geometry().area() == pytest.approx(5400)
# 54 cells x (10x10) map units, minus the slightly rounded corners
assert features[0].geometry().area() == pytest.approx(5400, rel=0.05)
# preview computation must not add layers to the project
assert len(QgsProject.instance().mapLayers()) == 0

Expand All @@ -110,12 +111,33 @@ def test_empty_mask_returns_no_features(self, canvas, polygon_maker_module):
assert maker.build_polygons(crs=CRS) == []


def max_turn_degrees(geometry) -> float:
"""Largest direction change (degrees) at any vertex of the
geometry's exterior ring."""
import math

ring = geometry.asPolygon()[0][:-1] # drop the closing point
worst = 0.0
for i in range(len(ring)):
a = ring[i - 1]
b = ring[i]
c = ring[(i + 1) % len(ring)]
heading_in = math.atan2(b.y() - a.y(), b.x() - a.x())
heading_out = math.atan2(c.y() - b.y(), c.x() - b.x())
turn = abs(math.degrees(heading_out - heading_in)) % 360
worst = max(worst, min(turn, 360 - turn))
return worst


@pytest.mark.usefixtures("native_processing", "qgis_new_project")
class TestSimplification:
def test_staircase_boundary_is_thinned(self, canvas, polygon_maker_module):
class TestSimplificationAndSmoothing:
def test_staircase_boundary_is_thinned_and_smoothed(
self, canvas, polygon_maker_module
):
# a pixel staircase (lower-left triangle of cells): the raw
# dissolved boundary has ~2 vertices per stair step; thinning
# should collapse it towards the diagonal without losing area
# dissolved boundary has ~2 vertices (a 90-degree zigzag) per
# stair step; thinning should collapse it towards the diagonal
# without losing area, and smoothing rounds what remains
bin_index = np.zeros((10, 20), dtype=bool)
for y in range(10):
bin_index[y, : y + 1] = True # 1+2+...+10 = 55 cells
Expand All @@ -126,9 +148,11 @@ def test_staircase_boundary_is_thinned(self, canvas, polygon_maker_module):
assert len(features) == 1
geometry = features[0].geometry()
assert geometry.area() == pytest.approx(5500, rel=0.1)
# raw staircase ring has ~23 vertices; the thinned ring must be
# substantially lighter
assert geometry.constGet().nCoordinates() <= 12
# no stair-step (90-degree) corners survive on the boundary
assert max_turn_degrees(geometry) < 60
# smoothing must stay cheap: the ring stays far lighter than
# Chaikin on the raw ~23-vertex staircase would produce
assert geometry.constGet().nCoordinates() <= 60


@pytest.mark.usefixtures("native_processing", "qgis_new_project")
Expand Down Expand Up @@ -198,9 +222,9 @@ def test_creates_new_layer_with_polygon(self, canvas, polygon_maker_module):
assert layer.name() == "magic_wand"
features = list(layer.getFeatures())
assert len(features) == 1
# 54 cells x (10x10) map units; the rectangle survives
# simplification exactly (area-based thinning keeps corners)
assert features[0].geometry().area() == pytest.approx(5400)
# 54 cells x (10x10) map units; area-based thinning keeps the
# rectangle exact, smoothing then rounds the corners slightly
assert features[0].geometry().area() == pytest.approx(5400, rel=0.05)

def test_appends_to_existing_layer(self, canvas, polygon_maker_module):
existing = QgsVectorLayer(f"Polygon?crs={CRS.authid()}", "existing", "memory")
Expand Down
8 changes: 4 additions & 4 deletions tests/test_preview_session.py
Original file line number Diff line number Diff line change
Expand Up @@ -57,7 +57,7 @@ def test_disjoint_seeds_produce_separate_features(self, modules):
)

areas = sorted(f.geometry().area() for f in features)
assert areas == [pytest.approx(20 * 25), pytest.approx(30 * 20)]
assert areas == [pytest.approx(20 * 25, rel=0.02), pytest.approx(30 * 20, rel=0.02)]

def test_seeds_in_the_same_region_do_not_duplicate(self, modules):
session_module, analyzer, grid = self._setup(modules, [(5, 5, 30, 20, RED)])
Expand All @@ -71,7 +71,7 @@ def test_seeds_in_the_same_region_do_not_duplicate(self, modules):
)

assert len(features) == 1
assert features[0].geometry().area() == pytest.approx(30 * 20)
assert features[0].geometry().area() == pytest.approx(30 * 20, rel=0.02)

def test_seed_colors_form_one_combined_model(self, modules):
# RED | BLUE | RED bands: seeding the left RED and the BLUE band
Expand All @@ -91,7 +91,7 @@ def test_seed_colors_form_one_combined_model(self, modules):
)

assert len(features) == 1
assert features[0].geometry().area() == pytest.approx(60 * 20)
assert features[0].geometry().area() == pytest.approx(60 * 20, rel=0.02)

def test_single_seed_matches_the_plain_flow(self, modules):
session_module, analyzer, grid = self._setup(modules, [(5, 5, 30, 20, RED)])
Expand All @@ -101,7 +101,7 @@ def test_single_seed_matches_the_plain_flow(self, modules):
)

assert len(features) == 1
assert features[0].geometry().area() == pytest.approx(30 * 20)
assert features[0].geometry().area() == pytest.approx(30 * 20, rel=0.02)

def test_no_seeds_returns_no_features(self, modules):
session_module, analyzer, grid = self._setup(modules, [])
Expand Down
10 changes: 5 additions & 5 deletions tests/test_processing_algorithm.py
Original file line number Diff line number Diff line change
Expand Up @@ -100,8 +100,8 @@ def test_one_selection_per_seed_feature(self, tmp_path):
output = result["OUTPUT"]
features = {f["seed_id"]: f for f in output.getFeatures()}
assert len(features) == 2
assert features[1].geometry().area() == pytest.approx(30 * 20)
assert features[2].geometry().area() == pytest.approx(20 * 25)
assert features[1].geometry().area() == pytest.approx(30 * 20, rel=0.02)
assert features[2].geometry().area() == pytest.approx(20 * 25, rel=0.02)
# each polygon contains its own seed
assert (
features[1].geometry().contains(QgsGeometry.fromPointXY(QgsPointXY(20, 25)))
Expand Down Expand Up @@ -137,7 +137,7 @@ def test_multipoint_seed_produces_one_merged_feature(self, tmp_path):
assert len(features) == 1 # one selection, one multipolygon feature
geometry = features[0].geometry()
assert geometry.isMultipart()
assert geometry.area() == pytest.approx(30 * 20 + 20 * 25)
assert geometry.area() == pytest.approx(30 * 20 + 20 * 25, rel=0.02)

def test_multipoint_seeds_in_the_same_region_do_not_duplicate(self, tmp_path):
from qgis import processing
Expand All @@ -157,7 +157,7 @@ def test_multipoint_seeds_in_the_same_region_do_not_duplicate(self, tmp_path):

features = list(result["OUTPUT"].getFeatures())
assert len(features) == 1
assert features[0].geometry().area() == pytest.approx(30 * 20)
assert features[0].geometry().area() == pytest.approx(30 * 20, rel=0.02)

def test_multipoint_colors_form_one_combined_model(self, tmp_path):
from qgis import processing
Expand Down Expand Up @@ -187,7 +187,7 @@ def test_multipoint_colors_form_one_combined_model(self, tmp_path):

features = list(result["OUTPUT"].getFeatures())
assert len(features) == 1
assert features[0].geometry().area() == pytest.approx(60 * 20)
assert features[0].geometry().area() == pytest.approx(60 * 20, rel=0.02)

def test_seed_outside_raster_is_skipped(self, tmp_path):
from qgis import processing
Expand Down
Loading