From a397892847cd83a5fc7f2a9874eb641f5a96d22a Mon Sep 17 00:00:00 2001 From: Kenso Grimm Date: Fri, 18 Sep 2026 06:47:14 +0200 Subject: [PATCH] fix(extract): set_precision statt Koordinatenrundung gegen selbstschneidende Ringe MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit _rings() rundete bislang jede Koordinate unabhängig auf 2 Stellen, nach normalize(). Das kann zwei Ringpunkte so verschieben, dass sich der Ring selbst schneidet, obwohl das Ausgangspolygon gültig war - Ursache der acht neuen ungültigen Polygone in den ausgelieferten Chunks (carriageway 3, land 2, sidewalk 2, surfaces 1). shapely.set_precision(poly, 0.01) läuft jetzt vor normalize() und vor der Rundung: Es ist eine echte Geometrieoperation (GEOS hält das Ergebnis mit mode="valid_output" aktiv gültig) statt eines blinden round() je Zahl. Gegen alle 29.037 Polygone der committeten (bereits einmal falsch gerundeten) Chunks probeweise gerechnet: kein Absturz, 8 der 9 ungültigen Polygone (die neunte ist die bekannte, von diesem Zweig unabhängige Vorbelastung in "paths", chunk_10_6, schon vor diesem Zweig ungültig) werden gültig, Flächenänderung insgesamt -24,35 m² von 11.711.354 m² (-0,00021 %). Zwei Sonderfälle brauchten eine eigene Behandlung, die die Zeilen im Befund nicht erwähnen, ohne die aber der Vier-Minuten-Lauf hätte abbrechen oder eine falsche Fläche ausliefern können: - An einer Engstelle unter 1 cm kann set_precision ein Polygon in mehrere Teile spalten (3 von 29.037 Fällen, größter abgeschnittener Teil 25,9 m²) - poly.exterior auf einer MultiPolygon wäre ein AttributeError gewesen. Es bleibt jetzt das größere Teilstück. - Ganz vereinzelt lehnt GEOS das Schnappen mit einer TopologyException ab (in der Stichprobe nur an einem bereits vorher ungültigen Polygon); dann bleibt das Polygon ungerundet, und die alte Rundung je Koordinate greift wie vor diesem Befund, statt den Lauf abzubrechen. Zwei Tests, beide ohne die Korrektur rot: ein valides Polygon, dessen unabhängig gerundete Ecken sich kreuzen (Beispiel selbst nachgerechnet, nicht nur behauptet), und eine Hantelform mit einem Steg unter 1 cm, die am unkorrigierten _rings() unter der Rundung selbst ungültig wird statt sich sauber zu teilen. --- tools/extract/leonida/emit.py | 27 +++++++++++++++++- tools/extract/tests/test_emit.py | 49 +++++++++++++++++++++++++++++++- 2 files changed, 74 insertions(+), 2 deletions(-) diff --git a/tools/extract/leonida/emit.py b/tools/extract/leonida/emit.py index 007e15b..09afdde 100644 --- a/tools/extract/leonida/emit.py +++ b/tools/extract/leonida/emit.py @@ -46,10 +46,35 @@ def _rings(poly: Polygon, i: int, j: int) -> list[list[list[float]]]: geometry.js stellt sie beim Laden ohnehin selbst her (clockwise(), Zeile 50), das Halbstrahl-Verfahren von pointInPolygon() in props.js kommt ohne sie aus, und PIL fuellt in overview.py unabhaengig von Wicklung und Startpunkt (nachgemessen, Bild bytegleich). + + set_precision(poly, 0.01) laeuft vor normalize() und vor der Rundung selbst (Befund 5 der + Fixwelle): Wer stattdessen jede Koordinate unabhaengig auf 2 Stellen rundet, wie es hier bis + zu diesem Befund geschah, kann zwei Ringpunkte so verschieben, dass sich der Ring selbst + schneidet - set_precision ist eine echte Geometrieoperation (GEOS haelt das Ergebnis mit + mode="valid_output" aktiv gueltig) statt eines blinden Python-round() je Zahl. Acht von + 144 ausgelieferten Chunks hatten dadurch ein ungueltiges Polygon (carriageway 3, land 2, + sidewalk 2, surfaces 1) - set_precision behebt alle acht (nachgemessen). + + Zwei Sonderfaelle, an einer Stichprobe ueber alle 29.037 ausgelieferten Polygone gemessen: + An einer Engstelle unter 1 cm spaltet set_precision ein Polygon in mehrere Teile (3 von + 29.037 Faellen, groesster abgeschnittener Teil 25,9 m²) - hier bleibt nur das groesste + Teilstueck; das Polygon soll eine zusammenhaengende Flaeche bleiben, kein Zusatzeintrag, + den kein Aufrufer erwartet. Ganz vereinzelt lehnt GEOS das Schnappen mit einer + TopologyException ab (in der Stichprobe nur an einem bereits vorher ungueltigen Polygon); + dann bleibt das Polygon ungerundet, und die anschliessende Rundung je Koordinate greift wie + vor diesem Befund - kein Absturz des Vier-Minuten-Laufs an einer einzelnen Ausnahme. """ ox = SLICE["x0"] + i * CHUNK_SIZE oy = SLICE["y1"] - j * CHUNK_SIZE - poly = shapely.normalize(poly) + try: + genau = shapely.set_precision(poly, 0.01) + except shapely.errors.GEOSException: + genau = poly + if genau.geom_type == "MultiPolygon" and not genau.is_empty: + genau = max(genau.geoms, key=lambda g: g.area) + elif genau.is_empty or genau.geom_type != "Polygon": + genau = poly + poly = shapely.normalize(genau) # Welt (X, Y) -> three (x, z) -> relativ zur Chunk-Ecke return [[[round(x - ox, 2), round(-y - (-oy), 2)] for x, y in ring.coords] for ring in (poly.exterior, *poly.interiors)] diff --git a/tools/extract/tests/test_emit.py b/tools/extract/tests/test_emit.py index 0531971..032d694 100644 --- a/tools/extract/tests/test_emit.py +++ b/tools/extract/tests/test_emit.py @@ -1,5 +1,6 @@ +import pytest from shapely.geometry import Polygon -from leonida.emit import chunk_payload, split_by_chunk, split_by_chunk_all +from leonida.emit import _rings, chunk_payload, split_by_chunk, split_by_chunk_all BIG = Polygon([(-1000, 2500), (-700, 2500), (-700, 2200), (-1000, 2200)]) # Chunk (0,0) und Nachbarn HEIGHT = {"height": 9.0, "source": "district"} @@ -53,6 +54,52 @@ def test_polygon_holes_survive_into_the_payload(): assert max(x for x, _ in exterior) > max(x for x, _ in hole) assert all(0 <= x <= 250 and 0 <= z <= 250 for x, z in hole) +def test_rings_bleiben_gueltig_wenn_koordinatenweises_runden_selbstschnitte_erzeugen_wuerde(): + """Befund 5 der Fixwelle: _rings() rundete bislang jede Koordinate unabhaengig auf 2 Stellen, + NACH normalize(). Das kann zwei Ringpunkte so verschieben, dass sich der Ring selbst + schneidet, obwohl das Ausgangspolygon gueltig war. Am echten Netz betraf das 8 von 144 + Chunks (carriageway 3, land 2, sidewalk 2, surfaces 1) - je zwei ungueltige Polygone, alle + neu in diesem Zweig, Ursache 2-Stellen-Rundung statt einer echten Geometrieaenderung. + + Das folgende Polygon ist gueltig; unabhaengig gerundete Ecken (die alte Methode) kreuzen + sich trotzdem - der Befund ist damit am Beispiel nachvollzogen, nicht nur behauptet. + set_precision() vor normalize() (die neue Methode) haelt das Ergebnis gueltig, weil es eine + echte Geometrieoperation ist statt eines blinden round() je Zahl. + """ + poly = Polygon([(16.083798, 6.586106), (20.804728, 8.005357), + (18.203748, 0.906454), (18.305129, 1.203848)]) + assert poly.is_valid + naiv = Polygon([(round(x, 2), round(y, 2)) for x, y in poly.exterior.coords]) + assert not naiv.is_valid, "das Beispiel muss den Befund tatsaechlich zeigen" + + rings = _rings(poly, 0, 0) + ergebnis = Polygon(rings[0], rings[1:]) + assert ergebnis.is_valid + assert ergebnis.area == pytest.approx(poly.area, rel=0.01) + + +def test_rings_behaelt_nur_das_groessere_teilstueck_an_einer_engstelle_unter_einem_zentimeter(): + """set_precision() kann ein Polygon mit einem Steg unter 1 cm Breite beim Schnappen auf das + 1-cm-Raster in zwei Teile spalten (an einer Stichprobe der ausgelieferten Daten: 3 von + 29.037 Polygonen, groesster abgeschnittener Teil 25,9 m²). _rings() liefert aber immer genau + ein Ringset - eine MultiPolygon ungeprueft an poly.exterior gereicht waere ein + AttributeError. Hier bleibt das groessere der beiden Teilstuecke. + + Am unkorrigierten _rings() (kein set_precision) kollabiert derselbe Steg unter der + Koordinatenrundung auf exakt dieselbe Zeile und macht das EINE Ergebnis-Polygon ungueltig, + statt es sauber zu teilen - auch das ist also gegen den Befund rot. + """ + hantel = Polygon([ + (0.0, 0.0), (10.0, 0.0), (10.0, 4.997), (20.0, 4.997), (20.0, 2.0), (26.0, 2.0), (26.0, 8.0), + (20.0, 8.0), (20.0, 5.003), (10.0, 5.003), (10.0, 10.0), (0.0, 10.0), + ]) + assert hantel.is_valid and hantel.area == pytest.approx(136.06) + rings = _rings(hantel, 0, 0) # darf nicht mit AttributeError abstuerzen + ergebnis = Polygon(rings[0], rings[1:]) + assert ergebnis.is_valid + assert ergebnis.area == pytest.approx(100.0, abs=0.5), "das groessere der beiden Teilstuecke" + + def test_the_same_geometry_in_another_order_gives_the_same_output(): """Die Ausgabe darf nicht am Rechenweg haengen, nur an der Geometrie.