feat(extract): Gebäudepolygone, Höhen aus Wahrzeichen und Vierteln, Palmenpunkte
buildings.py regularisiert Gebäude-Maskenpolygone (gebaeude/spekulativ) zu Rechtecken, wenn sie ihr gedrehtes Hüllrechteck zu mindestens 85 % ausfüllen, sonst zu vereinfachten Grundrissen; zu kleine oder entartete Polygone liefern None. heights.py weist Gebäuden eine Höhe zu (Wahrzeichen-Treffer schlägt Viertel-Bereich, deterministisch geseedet) und platziert Palmenpunkte entlang strandnaher Mittellinien: der Versatz vom Fahrbahnrand wird von innen nach außen gesucht (PALM_VERGE), bis eine Position sowohl PALM_KERB von der Fahrbahn als auch PALM_CLEARANCE von Gebäuden einhält; findet sich keine, entfällt die Palme. Start bei halber PALM_SPACING verhindert doppelte Palmen an Kreuzungen, wo Task 4 Mittellinien pro Ast trennt. pyproject.toml ergänzt einen gezielten filterwarnings-Eintrag: Shapelys oriented_envelope (minimum_rotated_rectangle) meldet für bereits achsen- parallele/rechtwinklige Polygone eine RuntimeWarning aus GEOS (Winkelterm wird 0/0), obwohl das Ergebnis nachweislich korrekt bleibt – u.a. löst der Briefs eigener L-Shape-Test das aus. Nur dieser eine Text+Kategorie- kombinierte Fall wird gefiltert, jede andere Warnung bleibt sichtbar. Auf dem echten Mosaik: 2141 Gebäudepolygone (258 spekulativ, 211 von regularize verworfen), Höhen 6.0-15.0 m (nur "district", da data/skyline.json und data/districts.json erst in Task 6 kommen), 2433 Palmenpunkte. Straßen- und Mittellinienzahlen (475 Stücke/40 Brücken, 2323 Linien/Median 27 m) decken sich exakt mit den in Task 4 unabhängig ermittelten Referenzwerten. Co-Authored-By: Claude Sonnet 5 <noreply@anthropic.com>
This commit is contained in:
knsoandClaude Sonnet 5 committed
1 parent
11f6aa76fe
commit
6de4b81942
5 files changed
+228
No files matched your search
@@ -0,0 +1,40 @@
|
||||
import warnings
|
||||
|
||||
import numpy as np
|
||||
from shapely.geometry import Polygon
|
||||
|
||||
from .classify import CLASSES
|
||||
from .vectors import mask_polygons
|
||||
|
||||
MIN_AREA = 12.0
|
||||
RECT_COVERAGE = 0.85
|
||||
|
||||
|
||||
def regularize(poly: Polygon) -> Polygon | None:
|
||||
"""Rechteck, wenn das Polygon sein gedrehtes Hüllrechteck gut ausfüllt; sonst vereinfacht."""
|
||||
poly = poly.buffer(0)
|
||||
if poly.geom_type != "Polygon" or poly.area < MIN_AREA:
|
||||
return None
|
||||
with warnings.catch_warnings():
|
||||
# Shapelys oriented_envelope (minimum_rotated_rectangle) wertet einen Winkelterm
|
||||
# aus, der für bereits achsenparallele/rechtwinklige Polygone exakt 0/0 wird -
|
||||
# GEOS meldet das als RuntimeWarning, obwohl das zurückgegebene Rechteck korrekt
|
||||
# bleibt (geprüft: Fläche stimmt, kein NaN). Nur diese eine Property-Abfrage ist
|
||||
# eingeschränkt, damit ein echter Warnhinweis an anderer Stelle sichtbar bliebe.
|
||||
warnings.filterwarnings("ignore", message=r".*oriented_envelope.*", category=RuntimeWarning)
|
||||
rect = poly.minimum_rotated_rectangle
|
||||
if rect.area > 0 and poly.area / rect.area >= RECT_COVERAGE:
|
||||
return rect
|
||||
simple = poly.simplify(1.0)
|
||||
return simple if simple.area >= MIN_AREA else None
|
||||
|
||||
|
||||
def building_polygons(labels: np.ndarray) -> list[tuple[Polygon, bool]]:
|
||||
out = []
|
||||
for name, speculative in (("gebaeude", False), ("spekulativ", True)):
|
||||
mask = labels == CLASSES[name]
|
||||
for poly in mask_polygons(mask, min_area=MIN_AREA, simplify=0.5):
|
||||
reg = regularize(poly)
|
||||
if reg is not None:
|
||||
out.append((reg, speculative))
|
||||
return out
|
||||
@@ -0,0 +1,102 @@
|
||||
import hashlib
|
||||
import math
|
||||
|
||||
from shapely.geometry import Point, Polygon
|
||||
|
||||
DISTRICT_RANGES = {
|
||||
"Downtown": (20.0, 60.0),
|
||||
"Ocean Beach": (8.0, 25.0),
|
||||
"default": (6.0, 15.0),
|
||||
}
|
||||
LANDMARK_SNAP = 20.0 # Meter: so weit darf ein Wahrzeichenpunkt vom Grundriss entfernt liegen
|
||||
|
||||
|
||||
def _seeded_unit(seed: int, key: str) -> float:
|
||||
digest = hashlib.sha256(f"{seed}:{key}".encode()).digest()
|
||||
return int.from_bytes(digest[:8], "big") / 2 ** 64
|
||||
|
||||
|
||||
def assign_heights(polys: list[Polygon], landmarks: list[dict],
|
||||
districts: list[dict], seed: int) -> list[dict]:
|
||||
out = []
|
||||
for poly in polys:
|
||||
centroid = poly.centroid
|
||||
key = f"{centroid.x:.1f},{centroid.y:.1f}"
|
||||
hit = None
|
||||
for lm in landmarks:
|
||||
p = Point(lm["xyz"][0], lm["xyz"][1])
|
||||
if poly.contains(p) or poly.distance(p) <= LANDMARK_SNAP:
|
||||
if hit is None or lm["height"] > hit["height"]:
|
||||
hit = lm
|
||||
if hit is not None:
|
||||
out.append({"height": float(hit["height"]), "source": "skyline"})
|
||||
continue
|
||||
name = "default"
|
||||
best = math.inf
|
||||
for d in districts:
|
||||
dist = math.dist((centroid.x, centroid.y), (d["x"], d["y"]))
|
||||
if dist < best:
|
||||
best, name = dist, d["name"] if d["name"] in DISTRICT_RANGES else "default"
|
||||
lo, hi = DISTRICT_RANGES[name]
|
||||
# Kleine Grundrisse bleiben niedriger: Fläche skaliert den oberen Rand.
|
||||
area_factor = min(1.0, poly.area / 400.0)
|
||||
height = lo + (hi - lo) * _seeded_unit(seed, key) * (0.4 + 0.6 * area_factor)
|
||||
out.append({"height": round(height, 1), "source": "district"})
|
||||
return out
|
||||
|
||||
|
||||
PALM_SPACING = 12.0 # Meter entlang der Mittellinie
|
||||
PALM_BEACH_RANGE = 40.0 # Meter Abstand zur nächsten Sandfläche
|
||||
PALM_CLEARANCE = 6.0 # Meter Abstand zu Gebäuden
|
||||
PALM_KERB = 1.5 # Meter Abstand zum Fahrbahnrand
|
||||
# Versatzkandidaten von innen nach außen; genommen wird der erste, der Fahrbahn und
|
||||
# Gebäude freihält. Ein fester Wert geht nicht: die Fahrbahn ist an den Mittellinien-
|
||||
# punkten p50 5,4 m, p90 11,7 m, p95 14,0 m halbbreit (am echten Mosaik gemessen).
|
||||
PALM_VERGE = (6.0, 8.0, 10.0, 12.0, 14.0, 16.0)
|
||||
|
||||
|
||||
def _verge_point(line, p, side: str, roads, buildings):
|
||||
"""Erster Versatz von innen nach außen, der Fahrbahn und Gebäude freihält."""
|
||||
for dist in PALM_VERGE:
|
||||
off = line.parallel_offset(dist, side)
|
||||
if off.is_empty:
|
||||
continue
|
||||
q = off.interpolate(off.project(p))
|
||||
if roads is not None and q.distance(roads) < PALM_KERB:
|
||||
continue
|
||||
if buildings is not None and q.distance(buildings) < PALM_CLEARANCE:
|
||||
continue
|
||||
return [(round(q.x, 2), round(q.y, 2))]
|
||||
return []
|
||||
|
||||
|
||||
def palm_points(sand_polys, park_polys, road_lines, building_polys, road_polys, seed: int):
|
||||
"""Palmen entlang strandnaher Mittellinien und verstreut in Parks."""
|
||||
from shapely.ops import unary_union
|
||||
sand = unary_union(sand_polys) if sand_polys else None
|
||||
roads = unary_union(list(road_polys)) if road_polys else None
|
||||
buildings = unary_union(list(building_polys)) if building_polys else None
|
||||
blocked = unary_union([g for g in (roads, buildings) if g is not None]) \
|
||||
if (roads is not None or buildings is not None) else None
|
||||
points = []
|
||||
for line in road_lines:
|
||||
if sand is None or line.distance(sand) > PALM_BEACH_RANGE:
|
||||
continue
|
||||
# Start bei halber Teilung: die Mittellinien sind an Kreuzungen getrennt, ein
|
||||
# Start bei 0 setzt je Ast eine Palme in denselben Punkt.
|
||||
pos = min(PALM_SPACING / 2, line.length / 2)
|
||||
while pos <= line.length:
|
||||
p = line.interpolate(pos)
|
||||
for side in ("left", "right"):
|
||||
points.extend(_verge_point(line, p, side, roads, buildings))
|
||||
pos += PALM_SPACING
|
||||
for park in park_polys:
|
||||
minx, miny, maxx, maxy = park.bounds
|
||||
n = int(park.area // 250)
|
||||
for k in range(n):
|
||||
u = _seeded_unit(seed, f"park:{minx:.0f}:{k}:x")
|
||||
v = _seeded_unit(seed, f"park:{miny:.0f}:{k}:y")
|
||||
q = Point(minx + u * (maxx - minx), miny + v * (maxy - miny))
|
||||
if park.contains(q) and (blocked is None or q.distance(blocked) > PALM_CLEARANCE):
|
||||
points.append((round(q.x, 2), round(q.y, 2)))
|
||||
return points
|
||||
@@ -13,3 +13,11 @@ dev = ["pytest>=8.3"]
|
||||
[tool.pytest.ini_options]
|
||||
testpaths = ["tests"]
|
||||
pythonpath = ["."]
|
||||
filterwarnings = [
|
||||
# Shapelys oriented_envelope (minimum_rotated_rectangle) meldet fuer bereits
|
||||
# achsenparallele/rechtwinklige Polygone eine RuntimeWarning aus GEOS (Winkelterm
|
||||
# wird 0/0), obwohl das zurueckgegebene Rechteck nachweislich korrekt bleibt. Nur
|
||||
# diese eine, an Text und Kategorie festgemachte Warnung wird gefiltert; jede andere
|
||||
# Warnung (auch jede echte DeprecationWarning) bleibt sichtbar.
|
||||
"ignore:.*oriented_envelope.*:RuntimeWarning",
|
||||
]
|
||||
@@ -0,0 +1,24 @@
|
||||
import numpy as np
|
||||
from shapely.geometry import Polygon
|
||||
from leonida.buildings import building_polygons, regularize
|
||||
from leonida.classify import CLASSES
|
||||
|
||||
def test_rough_rectangle_snaps_to_rectangle():
|
||||
poly = Polygon([(0, 0), (10.4, 0.3), (10.1, 20.2), (-0.2, 19.8)])
|
||||
out = regularize(poly)
|
||||
assert len(out.exterior.coords) == 5 # geschlossenes Rechteck
|
||||
assert abs(out.area - 205) < 25
|
||||
|
||||
def test_l_shape_is_kept_not_boxed():
|
||||
poly = Polygon([(0, 0), (30, 0), (30, 10), (10, 10), (10, 30), (0, 30)])
|
||||
out = regularize(poly)
|
||||
assert out.area < 0.8 * out.minimum_rotated_rectangle.area
|
||||
|
||||
def test_tiny_polygon_is_dropped():
|
||||
assert regularize(Polygon([(0, 0), (2, 0), (2, 2), (0, 2)])) is None
|
||||
|
||||
def test_two_buildings_are_separated():
|
||||
labels = np.full((60, 60), CLASSES["stadt"], dtype=np.uint8)
|
||||
labels[10:25, 10:25] = CLASSES["gebaeude"]
|
||||
labels[35:50, 35:50] = CLASSES["gebaeude"]
|
||||
assert len(building_polygons(labels)) == 2
|
||||
@@ -0,0 +1,54 @@
|
||||
from shapely.geometry import Polygon
|
||||
from leonida.heights import DISTRICT_RANGES, assign_heights
|
||||
|
||||
SQUARE = Polygon([(0, 0), (20, 0), (20, 20), (0, 20)])
|
||||
FAR = Polygon([(500, 500), (520, 500), (520, 520), (500, 520)])
|
||||
|
||||
def test_landmark_height_wins():
|
||||
res = assign_heights([SQUARE], landmarks=[{"xyz": [10, 10, 0], "height": 180.0}],
|
||||
districts=[], seed=1)
|
||||
assert res[0]["height"] == 180.0
|
||||
assert res[0]["source"] == "skyline"
|
||||
|
||||
def test_district_range_is_used_without_landmark():
|
||||
res = assign_heights([SQUARE], landmarks=[],
|
||||
districts=[{"name": "Downtown", "x": 5, "y": 5}], seed=1)
|
||||
lo, hi = DISTRICT_RANGES["Downtown"]
|
||||
assert lo <= res[0]["height"] <= hi
|
||||
assert res[0]["source"] == "district"
|
||||
|
||||
def test_same_seed_gives_same_height():
|
||||
args = dict(landmarks=[], districts=[{"name": "Downtown", "x": 5, "y": 5}], seed=42)
|
||||
assert assign_heights([SQUARE], **args)[0]["height"] == assign_heights([SQUARE], **args)[0]["height"]
|
||||
|
||||
def test_unknown_district_falls_back_to_default():
|
||||
res = assign_heights([FAR], landmarks=[], districts=[], seed=3)
|
||||
lo, hi = DISTRICT_RANGES["default"]
|
||||
assert lo <= res[0]["height"] <= hi
|
||||
|
||||
from shapely.geometry import LineString, Point
|
||||
from leonida.heights import PALM_KERB, palm_points
|
||||
|
||||
def test_no_palm_stands_on_the_carriageway():
|
||||
"""Der Versatz muss vom Fahrbahnrand her gedacht werden, nicht von einem festen Wert.
|
||||
|
||||
Am echten Mosaik hat die Fahrbahn an den Mittellinienpunkten p50 5,4 m, p90 11,7 m
|
||||
und p95 14,0 m Halbbreite. Ein fester Versatz von 4 m setzt 60 % der Palmen auf die
|
||||
Fahrbahn, einer von 11 m immer noch 12,5 % – und auf der Hälfte der Straßen stünde
|
||||
die Palme dann im Nachbargebäude.
|
||||
"""
|
||||
road = Polygon([(0, -9), (200, -9), (200, 9), (0, 9)]) # 18 m breite Fahrbahn
|
||||
line = LineString([(0, 0), (200, 0)])
|
||||
sand = Polygon([(0, 20), (200, 20), (200, 40), (0, 40)])
|
||||
pts = palm_points([sand], [], [line], [], [road], seed=1)
|
||||
assert pts
|
||||
assert all(road.distance(Point(x, y)) >= PALM_KERB for x, y in pts)
|
||||
|
||||
def test_palms_do_not_pile_up_where_centerlines_meet():
|
||||
"""Mittellinien sind an Kreuzungen getrennt (Task 4). Ein Start bei Position 0 setzt
|
||||
auf zwei gegenüberliegenden Ästen zwei Palmen in genau denselben Punkt."""
|
||||
sand = Polygon([(-80, 20), (80, 20), (80, 40), (-80, 40)])
|
||||
arms = [LineString([(0, 0), (60, 0)]), LineString([(0, 0), (-60, 0)])]
|
||||
pts = palm_points([sand], [], arms, [], [], seed=1)
|
||||
assert pts
|
||||
assert len(set(pts)) == len(pts)
|
||||
Reference in new issue
Block a user