import json from itertools import groupby from pathlib import Path import shapely from shapely.geometry import Polygon, box from .coords import CHUNK_SIZE, SLICE N = (SLICE["x1"] - SLICE["x0"]) // CHUNK_SIZE # 12 def chunk_box(i: int, j: int) -> Polygon: x0 = SLICE["x0"] + i * CHUNK_SIZE y1 = SLICE["y1"] - j * CHUNK_SIZE return box(x0, y1 - CHUNK_SIZE, x0 + CHUNK_SIZE, y1) def split_by_chunk(poly: Polygon) -> dict[tuple[int, int], Polygon]: out = {} for j in range(N): for i in range(N): part = poly.intersection(chunk_box(i, j)) if not part.is_empty and part.area > 0.5: out[(i, j)] = part return out def _rings(poly: Polygon, i: int, j: int) -> list[list[list[float]]]: """Aussenring zuerst, danach jedes Loch, in derselben Transformation. Ruling 34: ein Loch (z.B. ein von Strassen umschlossener Block, ein Innenhof, ein Teich im Park) ist genauso echte Flaeche wie der Aussenring und darf beim Zuschneiden nicht verschwinden - sonst wird der Block zu Asphalt, der Innenhof zu Gebaeudemasse, der Teich zu Wiese. normalize() zuerst, weil dieselbe Flaeche sonst je nach Rechenweg anders in der Datei steht: GEOS waehlt Ringstartpunkt und Wicklung nach Belieben. unary_union() in apply_street() fasst dafuer das kartenweite Strassenpolygon neu zusammen und aenderte so 88 Chunkdateien, in denen sich geometrisch nichts bewegt hatte - symmetrische Differenz der Strassenvereinigung exakt 0,00 m2 in allen 88. normalize() legt beides fest (lexikographisch kleinster Punkt zuerst, Aussenring im Uhrzeigersinn, Loecher dagegen) und macht die Ausgabe damit vom Rechenweg unabhaengig. Die vereinheitlichte Wicklung ist unbedenklich, nachgeprueft an allen drei Verbrauchern: 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 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)] def _sorted_surfaces(entries: list[dict]) -> list[dict]: """Belaege nur innerhalb einer Art sortieren - zwischen den Arten traegt die Reihenfolge Bedeutung. overview.py malt die Belaege in Listenreihenfolge uebereinander, der spaetere gewinnt. An den ausgelieferten Chunks ueberlappen sich verschiedene Arten auf 36.873 m2: der Parkrasen liegt unter Feldern (19.573 m2), Plaetzen (11.872 m2), Parkplatz (2.878 m2), Pool (539 m2) und Spielplatzsand (270 m2), die Stadtflaeche ueber dem Strandsand (1.729 m2). Ein Sortieren ueber die Arten hinweg kehrte das um und liesse die Felder unter dem Rasen verschwinden - und eine feste Rangfolge je Art gaebe es dafuer nicht, weil "sand" als Strand vor und als Spielplatz hinter "park" steht. Die Abfolge der Arten selbst steht fest im Code (run.py, SURFACE_ORDER in park.py) und ist damit schon eindeutig; GEOS mischt nur innerhalb einer Art. """ return [entry for _, run in groupby(entries, key=lambda e: e["kind"]) for entry in sorted(run, key=lambda e: e["poly"])] def _line(line, i: int, j: int) -> list[list[float]]: """Ein Linienzug in Chunk-Koordinaten, richtungsunabhängig festgelegt. normalize() dreht eine Linie in eine feste Laufrichtung, damit dieselbe Linie unabhängig vom Rechenweg gleich in der Datei steht. Zulässig ist das, weil kein Verbraucher aus Richtung oder Listenposition etwas ableitet: Die Laternen (sidewalkLampPoints() in props.js) suchen ihre Straßenseite an der Gehwegfläche, und die Oberkante des Bordsteins (kerbGeometry() in geometry.js) findet ihre Seite an der Fahrbahnfläche. Bei den früheren centerlines war das nicht so, deshalb blieben sie unsortiert und unnormalisiert. """ ox = SLICE["x0"] + i * CHUNK_SIZE oy = SLICE["y1"] - j * CHUNK_SIZE return [[round(x - ox, 2), round(-y + oy, 2)] for x, y in shapely.normalize(line).coords] def chunk_payload(i, j, *, land, surfaces, carriageway, curb, sidewalk, markings, lots, bridges, edges, buildings, palms, paths=()) -> dict: """Ein Chunk als Ausgabestruktur, deterministisch sortiert. Sortiert wird nach der gerundeten Ausgabestruktur selbst, nicht nach Flaeche, Schwerpunkt oder Bounding-Box: nur sie ist total geordnet und ohne Sonderfall, denn zwei Eintraege mit gleichem Schluessel sind dann bytegleich und ihre Reihenfolge wieder egal. carriageway traegt keine Klasse (Befund 6 der Fixwelle): carriageway_area() vereinigt das ganze Strassennetz, ein Eintrag ist eine Zusammenhangskomponente und bekaeme damit nur die breiteste Klasse irgendeiner enthaltenen Strasse zugewiesen - am echten Netz 91,2 % der Flaeche "boulevard", unabhaengig davon, was dort tatsaechlich liegt. Die Klasse je Strasse steht in edges; kein Verbraucher las carriageway[].klasse. """ ox = SLICE["x0"] + i * CHUNK_SIZE oy = SLICE["y1"] - j * CHUNK_SIZE return { "chunk": [i, j], "origin": [ox, -oy], "land": sorted(_rings(p, i, j) for p in land), "surfaces": _sorted_surfaces([{"poly": _rings(p, i, j), "kind": k} for p, k in surfaces]), "carriageway": sorted(({"poly": _rings(p, i, j)} for p in carriageway), key=lambda c: c["poly"]), "curb": sorted(_line(l, i, j) for l in curb), "sidewalk": sorted(({"poly": _rings(p, i, j)} for p in sidewalk), key=lambda s: s["poly"]), "markings": sorted(({"poly": _rings(p, i, j), "art": a} for p, a in markings), key=lambda m: (m["poly"], m["art"])), "lots": sorted(({"poly": _rings(p, i, j)} for p in lots), key=lambda l: l["poly"]), "bridges": sorted(({"poly": _rings(p, i, j)} for p in bridges), key=lambda b: b["poly"]), # breite ist die gemessene Kantenbreite (Median aus der Distanztransformation) und kommt # ungerundet an (Befund 8 der Fixwelle: Werte wie 31.622776601683793), waehrend jede # andere Zahl in der Ausgabe gerundet ist; fahrbahn ist dagegen schon eine der vier festen # SECTION-Konstanten und braucht keine Rundung. "edges": sorted(({"line": _line(l, i, j), "klasse": k, "breite": round(w, 2), "fahrbahn": f} for l, k, w, f in edges), key=lambda e: (e["line"], e["klasse"], e["breite"], e["fahrbahn"])), "buildings": sorted(({"poly": _rings(p, i, j), "height": h["height"], "source": h["source"], "speculative": s} for p, h, s in buildings), key=lambda b: (b["poly"], b["height"], b["source"], b["speculative"])), "palms": sorted([round(x - ox, 2), round(-y + oy, 2)] for x, y in palms), "paths": sorted(({"poly": _rings(p, i, j)} for p in paths), key=lambda p: p["poly"]), } def chunk_box_of(i: int, j: int) -> Polygon: return chunk_box(i, j) def split_by_chunk_all(polys, i: int, j: int, min_area: float = 0.5): """Alle Anteile der übergebenen Polygone, die in Chunk (i, j) liegen. min_area verwirft die Splitter, die beim Zuschneiden an der Chunkgrenze entstehen. 0,5 m² passt für Flächen; Markierungen sind kleiner und brauchen eine eigene Grenze (run.py, chunk_markings). """ cell = chunk_box(i, j) for poly in polys: part = poly.intersection(cell) if not part.is_empty and part.area > min_area: if part.geom_type == "Polygon": yield part else: yield from (g for g in part.geoms if g.geom_type == "Polygon" and g.area > min_area) def in_chunk(point, i: int, j: int) -> bool: x0 = SLICE["x0"] + i * CHUNK_SIZE y1 = SLICE["y1"] - j * CHUNK_SIZE return x0 <= point[0] < x0 + CHUNK_SIZE and y1 - CHUNK_SIZE <= point[1] < y1 def write_debug_png(out_dir: Path, mosaic, payloads: list[dict]) -> None: """Ein Bild je Chunk: die extrahierten Vektoren über der Karte, zur Sichtprüfung.""" from PIL import Image, ImageDraw out_dir.mkdir(parents=True, exist_ok=True) for payload in payloads: i, j = payload["chunk"] x0, y0 = i * CHUNK_SIZE, j * CHUNK_SIZE crop = Image.fromarray(mosaic[y0:y0 + CHUNK_SIZE, x0:x0 + CHUNK_SIZE]).convert("RGB") draw = ImageDraw.Draw(crop) for rings in payload["land"]: ring = rings[0] # nur der Aussenring; Loecher werden nicht gezeichnet draw.line([tuple(pt) for pt in ring] + [tuple(ring[0])], fill=(0, 255, 255), width=1) for road in payload["carriageway"]: ring = road["poly"][0] draw.line([tuple(pt) for pt in ring], fill=(255, 0, 255), width=1) for b in payload["buildings"]: ring = b["poly"][0] if len(ring) >= 3: draw.polygon([tuple(pt) for pt in ring], outline=(255, 220, 0)) for x, z in payload["palms"]: draw.ellipse([x - 2, z - 2, x + 2, z + 2], outline=(0, 200, 0)) crop.save(out_dir / f"chunk_{i}_{j}.png") def write_chunks(out_dir: Path, payloads: list[dict], meta: dict) -> None: out_dir.mkdir(parents=True, exist_ok=True) for payload in payloads: i, j = payload["chunk"] (out_dir / f"chunk_{i}_{j}.json").write_text(json.dumps(payload, separators=(",", ":"))) (out_dir / "manifest.json").write_text(json.dumps(meta, indent=1))