"""OSM XML parsing and the local metric projection. Pure Python: no `bpy`, so this runs and tests outside Blender. """ import math import xml.etree.ElementTree as ET def _element_id(element): return element.attrib.get("id", "") def tags(element): return {t.attrib.get("k", ""): t.attrib.get("v", "") for t in element.findall("tag")} def _coords_from_refs(refs, nodes): return [nodes[r] for r in refs if r in nodes] def _open_ring(coords): if len(coords) > 1 and coords[0] == coords[-1]: return coords[:-1] return coords def _is_closed_ring(coords): return len(coords) >= 4 and coords[0] == coords[-1] and len(_open_ring(coords)) >= 3 def _join_member_rings(members): """Build closed rings from relation member coordinate runs. Handles the common OSM multipolygon cases: members are either already closed ways, or open way fragments whose endpoints can be stitched together. Malformed leftovers are dropped rather than aborting the whole import. """ rings = [] pending = [list(member) for member in members if len(member) >= 2] while pending: ring = pending.pop(0) changed = True while not _is_closed_ring(ring) and changed: changed = False for index, candidate in enumerate(pending): if ring[-1] == candidate[0]: ring.extend(candidate[1:]) elif ring[-1] == candidate[-1]: ring.extend(reversed(candidate[:-1])) elif ring[0] == candidate[-1]: ring = candidate[:-1] + ring elif ring[0] == candidate[0]: ring = list(reversed(candidate[1:])) + ring else: continue pending.pop(index) changed = True break if _is_closed_ring(ring): rings.append(ring) return rings def parse_osm(path): root = ET.parse(path).getroot() bounds_node = root.find("bounds") if bounds_node is None: raise RuntimeError("OSM file does not contain a bounds element") bounds = {"min_lon": float(bounds_node.attrib["minlon"]), "min_lat": float(bounds_node.attrib["minlat"]), "max_lon": float(bounds_node.attrib["maxlon"]), "max_lat": float(bounds_node.attrib["maxlat"])} nodes = {} point_features = [] for node in root.findall("node"): try: node_id = int(node.attrib["id"]) coord = (float(node.attrib["lon"]), float(node.attrib["lat"])) node_tags = tags(node) nodes[node_id] = coord if node_tags: point_features.append({"id": node.attrib.get("id", ""), "coord": coord, "tags": node_tags}) except (KeyError, ValueError): continue ways = [] way_coords_by_id = {} for way in root.findall("way"): if way.attrib.get("action") == "delete": continue refs = [] for ref in way.findall("nd"): try: refs.append(int(ref.attrib["ref"])) except (KeyError, ValueError): pass coords = _coords_from_refs(refs, nodes) way_id = _element_id(way) if len(coords) >= 2: way_coords_by_id[way_id] = coords way_tags = tags(way) if len(coords) >= 2: ways.append({"id": way_id, "coords": coords, "tags": way_tags}) for relation in root.findall("relation"): if relation.attrib.get("action") == "delete": continue relation_tags = tags(relation) if relation_tags.get("type") != "multipolygon" or "building" not in relation_tags: continue outer_members = [] inner_members = [] for member in relation.findall("member"): if member.attrib.get("type") != "way": continue coords = way_coords_by_id.get(member.attrib.get("ref", "")) if not coords: continue role = member.attrib.get("role", "") if role == "inner": inner_members.append(coords) elif role in ("", "outer"): outer_members.append(coords) outer_rings = _join_member_rings(outer_members) if not outer_rings: continue inner_rings = _join_member_rings(inner_members) relation_id = _element_id(relation) for index, outer in enumerate(outer_rings): synthetic_id = relation_id if len(outer_rings) == 1 else f"{relation_id}:{index + 1}" ways.append({ "id": synthetic_id, "coords": outer, "inner_coords": inner_rings, "tags": relation_tags, "source": "relation", }) return bounds, ways, point_features def parse_height(feature_tags, default): try: return max(0.5, float(feature_tags.get("height", default))) except (TypeError, ValueError): return default class Projector: """Equirectangular projection about the centre of the OSM bounds. Output is metres in a local ENU frame (X east, Y north), which is what both the Blender scene and the Cesium GLB are authored in. """ def __init__(self, bounds): self.bounds = bounds self.lon0 = (bounds["min_lon"] + bounds["max_lon"]) / 2 self.lat0 = (bounds["min_lat"] + bounds["max_lat"]) / 2 self.m_per_lat = 111320.0 self.m_per_lon = 111320.0 * math.cos(math.radians(self.lat0)) def xy(self, lon_lat): lon, lat = lon_lat return ((lon - self.lon0) * self.m_per_lon, (lat - self.lat0) * self.m_per_lat) def inside(self, lon_lat, pad=0.00035): lon, lat = lon_lat b = self.bounds return (b["min_lon"] - pad <= lon <= b["max_lon"] + pad and b["min_lat"] - pad <= lat <= b["max_lat"] + pad) def ring(self, coords): return [self.xy(c) for c in coords]