Add orphan bus stops (platforms not in route relations)

Some bus stop platform nodes in OSM are not members of any route
relation (e.g. Borgergade). These were missed by the route-membership
download approach. Now download_osm.py:
- tracks bus platform node ids captured by route queries
- runs a supplementary query for all highway=bus_stop + platform nodes
- appends orphans (with empty route_refs) to the bus stops layer

render.py: orphan stops (empty route_refs) fall back to distance
clustering instead of route-membership grouping.

Results: 1164 route-member stops + 707 orphans = 1871 total
-> 1370 after clip -> 1049 markers after grouping

.gitignore

@@ -0,0 +1,18 @@
+# Generated/regenerable data — keep directory structure via .gitkeep only
+copenhagen/data/raw/**
+copenhagen/data/processed/**
+!copenhagen/data/raw/.gitkeep
+!copenhagen/data/processed/.gitkeep
+
+# Generated map outputs (regenerable via render.py)
+copenhagen/output/*.png
+copenhagen/output/*.svg
+copenhagen/output/*.pdf
+
+# Python
+.venv/
+__pycache__/
+*.pyc
+
+# OS
+.DS_Store

copenhagen/config/area.json

@@ -0,0 +1,9 @@
+{
+  "description": "City Pass area approximation: Copenhagen Kommune + Frederiksberg Kommune + Amager island (dissolved).",
+  "source": "OpenStreetMap API relation/{id}/full.json",
+  "relations": [
+    {"name": "Kobenhavns Kommune", "id": 2192363},
+    {"name": "Frederiksberg Kommune", "id": 2186660},
+    {"name": "Amager", "id": 5175924}
+  ]
+}

copenhagen/config/modes.yaml

@@ -0,0 +1,40 @@
+# Mapping of map modes to OSM route filters, GTFS route types, and style keys.
+# `style` references a key in styling.yaml.
+#
+# OSM filter keys map directly to Overpass tag filters. A key suffixed with
+# `_neq` emits a negated tag filter (["k"!="v"]); any other key emits an
+# equality filter (["k"="v"]). The bbox is appended by download_osm.py.
+#
+# S-tog are tagged route=light_rail in OSM (NOT route=train); they are scoped
+# by network="Takst Sjælland". The separate `light_rail` mode catches other
+# light rail (e.g. Hovedstadens Letbane, opened Aug 2026) via network!=Takst.
+# GTFS route_type: 0=tram,1=subway,2=rail,3=bus,4=ferry.
+modes:
+  subway:
+    osm:
+      route: subway
+    gtfs_route_type: 1
+    style: metro
+  s_tog:
+    osm:
+      route: light_rail
+      network: "Takst Sjælland"
+    gtfs_route_type: 2
+    style: s_tog
+  light_rail:
+    osm:
+      route: light_rail
+      network_neq: "Takst Sjælland"
+    gtfs_route_type: 2
+    style: light_rail
+  regional:
+    osm:
+      route: train
+    gtfs_route_type: 2
+    style: regional
+  bus:
+    osm:
+      route: bus
+      network: Movia
+    gtfs_route_type: 3
+    style: bus

copenhagen/config/styling.yaml

@@ -0,0 +1,96 @@
+# Styling for the Copenhagen transit map. Referenced by render.py and prepare.py.
+# Colours keyed by `style` (from modes.yaml) then route ref.
+# zorder is draw order bottom->top (higher drawn later / on top).
+
+zorder:
+  bus: 1
+  s_tog: 2
+  light_rail: 3
+  regional: 4
+  metro: 5
+
+line_width:           # body width in points (outline = body + outline_width)
+  bus: 0.8
+  s_tog: 2.0
+  light_rail: 2.0
+  regional: 1.8
+  metro: 2.6
+
+outline_width: 1.2    # added to line_width for the dark casing
+outline_color: "#2b2b2b"
+
+alpha:
+  bus: 0.45
+  s_tog: 0.95
+  light_rail: 0.95
+  regional: 0.9
+  metro: 1.0
+
+# Fallback palette when no GTFS route_color and no OSM colour tag are present.
+# Metro / S-tog values mirror the OSM colour tags (verified). Regional trains
+# have no OSM colour, so the palette supplies one.
+palette:
+  metro:
+    M1: "#008d41"
+    M2: "#ffc600"
+    M3: "#ff0a0a"
+    M4: "#009cd3"
+  s_tog:
+    A: "#00a4eb"
+    B: "#50ae30"
+    Bx: "#adce6d"
+    C: "#f68b1f"
+    E: "#7670b3"
+    F: "#fcc019"
+    H: "#e63511"
+  regional:
+    default: "#7a1b3e"
+  light_rail:
+    default: "#e89b3a"
+  bus:
+    A: "#d6264a"
+    C: "#16a085"
+    S: "#2a6fb5"
+    default: "#c8a80e"
+
+bus_filters:
+  # ref-based bus categorisation (regex applied to route ref)
+  categories:
+    A: "^\\d+A$"
+    C: "^\\d+C$"
+    S: "^\\d+S$"
+    regular: "^\\d+$"
+
+mask:
+  fade_distance_m: 2500   # shapeburst fade width outside the area
+  max_alpha: 0.85
+  color: white
+
+basemap:
+  provider: Esri.WorldGrayCanvas
+  zoom: 14
+
+stations:
+  show: true              # draw station markers (symbols without text)
+  styles: [metro, s_tog, regional, bus]
+  stop_cluster_m:          # merge same-name stops within this distance (per mode)
+    metro: 200             # covers platform/entrance spread (~170m max)
+    s_tog: 400             # covers platform/entrance spread (~380m max)
+    regional: 500          # covers platform spread (~475m max)
+    bus: 100               # opposite sides of road, wide boulevards
+  marker:
+    shape: circle         # circle | square
+    fill: white
+    edge: "#2b2b2b"
+    size:                 # diameter in points
+      metro: 3.5
+      s_tog: 3.0
+      regional: 2.5
+      bus: 1.5
+    linewidth: 0.8
+
+labels:
+  show: false             # off by default; enable via --labels
+  styles: [metro, s_tog]   # which station types to label
+  min_fontsize: 5
+  max_fontsize: 9

copenhagen/scripts/_common.py

@@ -0,0 +1,14 @@
+"""Shared paths and constants for the Copenhagen map scripts."""
+from pathlib import Path
+
+HERE = Path(__file__).resolve().parent
+CPH = HERE.parent
+CONFIG = CPH / "config"
+DATA = CPH / "data"
+RAW = DATA / "raw"
+OSM_RAW = RAW / "osm"
+GTFS_RAW = RAW / "gtfs"
+PROCESSED = DATA / "processed"
+OUTPUT = CPH / "output"
+
+UA = "jetlag-maps/0.1 (https://github.com/marvin/jetlag-maps)"

copenhagen/scripts/build_area.py

@@ -0,0 +1,149 @@
+#!/usr/bin/env python3
+"""Build the City Pass area polygon from three OSM boundary relations.
+
+Fetches relation/{id}/full.json from the OSM API for each relation in
+config/area.json, assembles member ways into rings (respecting outer/inner
+roles), builds polygons, and dissolves the union into a single multipolygon.
+
+Outputs (in data/processed):
+    area.geojson  (EPSG:4326, human-readable + portable)
+    area.gpkg     (EPSG:4326, for geopandas/scripts)
+"""
+import json
+import sys
+from pathlib import Path
+
+import geopandas as gpd
+import requests
+from shapely.geometry import MultiPolygon, Polygon
+from shapely.ops import unary_union
+
+from _common import CONFIG, PROCESSED, UA
+
+OSM_API = "https://api.openstreetmap.org/api/0.6"
+
+
+def fetch_relation(rel_id):
+    url = f"{OSM_API}/relation/{rel_id}/full.json"
+    r = requests.get(url, headers={"User-Agent": UA}, timeout=120)
+    r.raise_for_status()
+    return r.json()["elements"]
+
+
+def stitch_rings(ways):
+    """Stitch a list of node-coordinate polylines into closed rings.
+
+    `ways` is a list of lists of (lon, lat). Returns a list of closed rings
+    (each starting point == ending point).
+    """
+    eps = 1e-9
+
+    def close(a, b):
+        return abs(a[0] - b[0]) < eps and abs(a[1] - b[1]) < eps
+
+    pool = [w[:] for w in ways if w]
+    rings = []
+    while pool:
+        cur = pool.pop(0)
+        changed = True
+        while changed and not close(cur[0], cur[-1]):
+            changed = False
+            for i, w in enumerate(pool):
+                if close(cur[-1], w[0]):
+                    cur = cur + w[1:]
+                    pool.pop(i)
+                    changed = True
+                    break
+                if close(cur[-1], w[-1]):
+                    cur = cur + w[::-1][1:]
+                    pool.pop(i)
+                    changed = True
+                    break
+        rings.append(cur)
+    return [r for r in rings if close(r[0], r[-1]) and len(r) >= 4]
+
+
+def build_polygon(elements):
+    nodes = {}
+    ways = {}
+    relation = None
+    for e in elements:
+        t = e["type"]
+        if t == "node":
+            nodes[e["id"]] = (e["lon"], e["lat"])
+        elif t == "way":
+            ways[e["id"]] = e["nodes"]
+        elif t == "relation":
+            relation = e
+    if relation is None:
+        raise ValueError("no relation in element set")
+
+    outer, inner = [], []
+    for m in relation["members"]:
+        if m["type"] != "way" or m["ref"] not in ways:
+            continue
+        seq = [nodes[n] for n in ways[m["ref"]] if n in nodes]
+        if len(seq) < 2:
+            continue
+        (inner if m.get("role") == "inner" else outer).append(seq)
+
+    outer_rings = stitch_rings(outer)
+    inner_rings = stitch_rings(inner)
+
+    polys = []
+    for oring in outer_rings:
+        op = Polygon(oring)
+        if not op.is_valid:
+            op = op.buffer(0)
+        holes = []
+        for iring in inner_rings:
+            ip = Polygon(iring)
+            if not ip.is_valid:
+                ip = ip.buffer(0)
+            if op.contains(ip.representative_point()):
+                holes.append(list(ip.exterior.coords))
+        try:
+            polys.append(Polygon(oring, holes=holes))
+        except Exception:
+            polys.append(Polygon(oring))
+    if not polys:
+        raise ValueError("could not assemble any outer rings")
+    if len(polys) == 1:
+        return polys[0]
+    return MultiPolygon(polys)
+
+
+def main():
+    area_cfg = json.loads((CONFIG / "area.json").read_text())
+    geoms = []
+    for rel in area_cfg["relations"]:
+        rid = rel["id"]
+        print(f"fetching relation {rid} ({rel['name']})...", flush=True)
+        elements = fetch_relation(rid)
+        geom = build_polygon(elements)
+        if not geom.is_valid:
+            geom = geom.buffer(0)
+        geoms.append(geom)
+        print(f"  -> {geom.geom_type}, area={geom.area:.5f} deg^2", flush=True)
+
+    dissolved = unary_union(geoms)
+    if not dissolved.is_valid:
+        dissolved = dissolved.buffer(0)
+
+    PROCESSED.mkdir(parents=True, exist_ok=True)
+    gdf = gpd.GeoDataFrame(
+        {"name": ["City Pass area"]}, geometry=[dissolved], crs="EPSG:4326"
+    )
+    gdf.to_file(PROCESSED / "area.geojson", driver="GeoJSON")
+    gdf.to_file(PROCESSED / "area.gpkg", driver="GPKG", layer="area")
+    bounds = dissolved.bounds
+    print(
+        f"dissolved: {dissolved.geom_type} | bounds "
+        f"({bounds[0]:.4f},{bounds[1]:.4f})-({bounds[2]:.4f},{bounds[3]:.4f})",
+        flush=True,
+    )
+    print(f"wrote {PROCESSED/'area.geojson'} and {PROCESSED/'area.gpkg'}")
+
+
+if __name__ == "__main__":
+    main()

copenhagen/scripts/download_gtfs.py

@@ -0,0 +1,45 @@
+#!/usr/bin/env python3
+"""Download the Rejseplanen GTFS feed for the Copenhagen area.
+
+STATUS: BLOCKED. Rejseplanen (journey planner) does not currently expose a
+public GTFS download. This script documents the intended flow and exits
+non-zero with guidance.
+
+When a feed URL becomes available, set it here and the script will download
+and unzip into data/raw/gtfs/.
+"""
+import sys
+import urllib.request
+import zipfile
+from pathlib import Path
+
+from _common import GTFS_RAW
+
+# Fill in once access is granted. Likely candidates:
+#   - Rejseplanen / DOT open-data portal
+#   - a static GTFS zip provided on request
+GTFS_URL = None  # e.g. "https://.../rejseplanen.zip"
+
+
+def main():
+    if not GTFS_URL:
+        print(
+            "download_gtfs: BLOCKED — no GTFS feed URL configured.\n"
+            "Set GTFS_URL in this script once Rejseplanen grants access, then re-run.\n"
+            "The OSM-only pipeline (build_area -> download_osm -> prepare -> render)\n"
+            "is fully functional without GTFS; GTFS only enriches colours/stops.",
+            file=sys.stderr,
+        )
+        sys.exit(2)
+    GTFS_RAW.mkdir(parents=True, exist_ok=True)
+    zip_path = GTFS_RAW / "feed.zip"
+    print(f"downloading {GTFS_URL} -> {zip_path} ...", flush=True)
+    urllib.request.urlretrieve(GTFS_URL, zip_path)
+    with zipfile.ZipFile(zip_path) as z:
+        z.extractall(GTFS_RAW)
+    zip_path.unlink(missing_ok=True)
+    print(f"extracted GTFS txt files into {GTFS_RAW}")
+
+
+if __name__ == "__main__":
+    main()

copenhagen/scripts/download_osm.py

@@ -0,0 +1,302 @@
+#!/usr/bin/env python3
+"""Download transit route lines + stops from Overpass, per mode.
+
+For each mode in config/modes.yaml, runs an Overpass query (with mirror
+fallback + exponential backoff) and writes:
+    data/raw/osm/{mode}.gpkg          (layers: "lines", "stops", EPSG:4326)
+    data/raw/osm/{mode}.overpassql    (the exact query text)
+
+Lines: one MultiLineString per route relation (concatenation of member way
+geometries), carrying ref/name/network/colour/operator/route/mode tags.
+Stops: one Point per platform node member, carrying ref/name/mode/route_refs
+(route refs of routes using this stop, for disambiguation).
+"""
+import sys
+import time
+from pathlib import Path
+
+import geopandas as gpd
+import requests
+import yaml
+from shapely.geometry import LineString, MultiLineString, Point
+
+from _common import CONFIG, OSM_RAW, PROCESSED, UA
+
+MIRRORS = [
+    "https://overpass-api.de/api/interpreter",
+    "https://overpass.kumi.systems/api/interpreter",
+    "https://overpass.private.coffee/api/interpreter",
+]
+
+PLATFORM_ROLES = {"platform", "platform_entry", "platform_exit"}
+STOP_ROLES = {"stop", "stop_entry", "stop_exit", "stop_entry_only", "stop_exit_only"}
+# bus stops: use platform nodes (where passengers wait)
+# rail stops: use stop nodes (platforms are often areas, stop nodes are reliable points)
+BUS_STOP_ROLES = PLATFORM_ROLES
+RAIL_STOP_ROLES = STOP_ROLES | PLATFORM_ROLES
+
+
+def build_query(mode_cfg, bbox):
+    s, n, w, e = bbox
+    filters = []
+    for k, v in mode_cfg["osm"].items():
+        if k.endswith("_neq"):
+            filters.append(f'["{k[:-4]}"!="{v}"]')
+        else:
+            filters.append(f'["{k}"="{v}"]')
+    filt = "".join(filters)
+    return f"""[out:json][timeout:180];
+relation{filt}({s},{w},{n},{e});
+out geom;
+>;
+out body qt;
+"""
+
+
+def overpass(ql):
+    last_err = None
+    for mi, mirror in enumerate(MIRRORS):
+        for attempt in range(4):
+            try:
+                r = requests.post(
+                    mirror,
+                    data={"data": ql},
+                    headers={"User-Agent": UA},
+                    timeout=300,
+                )
+                if r.status_code == 429 or r.status_code >= 500:
+                    raise RuntimeError(f"HTTP {r.status_code}")
+                r.raise_for_status()
+                return r.json()
+            except Exception as e:
+                last_err = e
+                wait = 2 ** (attempt + mi)
+                print(f"  [{mirror}] attempt {attempt+1} failed: {e}; retry in {wait}s", flush=True)
+                time.sleep(wait)
+    raise RuntimeError(f"all mirrors failed: {last_err}")
+
+
+def parse(data, mode):
+    is_bus = mode == "bus"
+    stop_roles = BUS_STOP_ROLES if is_bus else RAIL_STOP_ROLES
+
+    nodes = {}
+    for e in data["elements"]:
+        if e["type"] == "node":
+            nodes[e["id"]] = e
+
+    lines = []
+    # collect route_refs per stop node: {node_id: set(route_refs)}
+    node_routes = {}
+    for e in data["elements"]:
+        if e["type"] != "relation" or "tags" not in e:
+            continue
+        tags = e["tags"]
+        if tags.get("route") is None:
+            continue
+        route_ref = tags.get("ref")
+        # collect line geometries
+        segs = []
+        for m in e["members"]:
+            if m["type"] == "way" and "geometry" in m and m.get("role") not in PLATFORM_ROLES:
+                coords = [(g["lon"], g["lat"]) for g in m["geometry"]]
+                if len(coords) >= 2:
+                    segs.append(LineString(coords))
+        geom = MultiLineString(segs) if len(segs) > 1 else (segs[0] if segs else None)
+        if geom is None:
+            continue
+        lines.append({
+            "mode": mode,
+            "ref": route_ref,
+            "name": tags.get("name"),
+            "network": tags.get("network"),
+            "colour": tags.get("colour"),
+            "operator": tags.get("operator"),
+            "route": tags.get("route"),
+            "geometry": geom,
+        })
+        # record route_ref on each stop member node
+        for m in e["members"]:
+            if m["type"] != "node" or m.get("role") not in stop_roles:
+                continue
+            nid = m["ref"]
+            if nid not in nodes:
+                continue
+            node_routes.setdefault(nid, set())
+            if route_ref:
+                node_routes[nid].add(route_ref)
+
+    # build stops from collected nodes
+    stops = []
+    stop_seen = set()
+    for nid, routes in node_routes.items():
+        if nid in stop_seen:
+            continue
+        nd = nodes[nid]
+        if "lat" not in nd or "lon" not in nd:
+            continue
+        stop_seen.add(nid)
+        ntags = nd.get("tags", {})
+        stops.append({
+            "mode": mode,
+            "name": ntags.get("name") or ntags.get("public_transport") or "stop",
+            "ref": ntags.get("ref"),
+            "public_transport": ntags.get("public_transport"),
+            "railway": ntags.get("railway"),
+            "route_refs": ";".join(sorted(routes)),
+            "geometry": Point(nd["lon"], nd["lat"]),
+        })
+    return lines, stops
+
+
+def stations_query(bbox):
+    s, n, w, e = bbox
+    return f"""[out:json][timeout:180];
+(
+  node["railway"="station"]({s},{w},{n},{e});
+  node["public_transport"="station"]({s},{w},{n},{e});
+  way["railway"="station"]({s},{w},{n},{e});
+  way["public_transport"="station"]({s},{w},{n},{e});
+);
+out center;
+"""
+
+
+def parse_stations(data):
+    """Stations (nodes or area centroids) with proper human-readable names."""
+    rows = []
+    seen = set()
+    for e in data["elements"]:
+        if e["type"] not in ("node", "way"):
+            continue
+        tags = e.get("tags", {})
+        name = tags.get("name")
+        if not name:
+            continue
+        if e["type"] == "node":
+            geom = Point(e["lon"], e["lat"])
+        else:
+            c = e.get("center")
+            if not c:
+                continue
+            geom = Point(c["lon"], c["lat"])
+        key = (name, round(geom.x, 5), round(geom.y, 5))
+        if key in seen:
+            continue
+        seen.add(key)
+        rows.append({
+            "name": name,
+            "railway": tags.get("railway"),
+            "public_transport": tags.get("public_transport"),
+            "station": tags.get("station"),
+            "subway": tags.get("subway"),
+            "light_rail": tags.get("light_rail"),
+            "train": tags.get("train"),
+            "uic_ref": tags.get("uic_ref"),
+            "geometry": geom,
+        })
+    return rows
+
+
+def orphan_bus_stops_query(bbox):
+    s, n, w, e = bbox
+    return f"""[out:json][timeout:180];
+(
+  node["highway"="bus_stop"]["public_transport"="platform"]({s},{w},{n},{e});
+);
+out body;
+"""
+
+
+def parse_orphan_bus_stops(data, known_node_ids):
+    """Parse standalone bus platform nodes not already captured by routes."""
+    stops = []
+    for e in data["elements"]:
+        if e["type"] != "node" or e["id"] in known_node_ids:
+            continue
+        tags = e.get("tags", {})
+        if tags.get("public_transport") != "platform":
+            continue
+        name = tags.get("name")
+        if not name:
+            continue
+        stops.append({
+            "mode": "bus",
+            "name": name,
+            "ref": tags.get("ref"),
+            "public_transport": tags.get("public_transport"),
+            "railway": tags.get("railway"),
+            "route_refs": "",  # orphan — no route membership
+            "geometry": Point(e["lon"], e["lat"]),
+        })
+    return stops
+
+
+def main():
+    modes = yaml.safe_load((CONFIG / "modes.yaml").read_text())["modes"]
+    area = gpd.read_file(PROCESSED / "area.gpkg")
+    w, s, e, n = area.total_bounds  # (minx, miny, maxx, maxy) = (west, south, east, north)
+    bbox = (s, n, w, e)
+
+    OSM_RAW.mkdir(parents=True, exist_ok=True)
+    bus_stop_node_ids = set()
+    for mode, cfg in modes.items():
+        ql = build_query(cfg, bbox)
+        (OSM_RAW / f"{mode}.overpassql").write_text(ql)
+        print(f"[{mode}] querying Overpass bbox=({s:.4f},{w:.4f},{n:.4f},{e:.4f})...", flush=True)
+        data = overpass(ql)
+        lines, stops = parse(data, mode)
+        # track bus stop node ids to avoid dupes with orphans
+        if mode == "bus":
+            for e in data["elements"]:
+                if e["type"] == "node" and e.get("tags", {}).get("public_transport") == "platform":
+                    bus_stop_node_ids.add(e["id"])
+        print(f"  -> {len(lines)} lines, {len(stops)} stops", flush=True)
+        if not lines:
+            continue
+        lgdf = gpd.GeoDataFrame(lines, crs="EPSG:4326")
+        lgdf.to_file(OSM_RAW / f"{mode}.gpkg", driver="GPKG", layer="lines")
+        if stops:
+            sgdf = gpd.GeoDataFrame(stops, crs="EPSG:4326")
+            sgdf.to_file(OSM_RAW / f"{mode}.gpkg", driver="GPKG", layer="stops")
+
+    # orphan bus stops (platform nodes not in any route relation)
+    oql = orphan_bus_stops_query(bbox)
+    (OSM_RAW / "orphan_bus_stops.overpassql").write_text(oql)
+    print("[orphan_bus_stops] querying Overpass...", flush=True)
+    odata = overpass(oql)
+    orphans = parse_orphan_bus_stops(odata, bus_stop_node_ids)
+    print(f"  -> {len(orphans)} orphan bus stops (not in route relations)", flush=True)
+    if orphans:
+        # append to bus.gpkg stops layer
+        existing = gpd.read_file(OSM_RAW / "bus.gpkg", layer="stops")
+        combined = gpd.GeoDataFrame(
+            __import__("pandas").concat([existing, gpd.GeoDataFrame(orphans, crs="EPSG:4326")],
+                                        ignore_index=True),
+            crs="EPSG:4326"
+        )
+        # rewrite stops layer
+        import tempfile, shutil
+        tmp = OSM_RAW / "bus_tmp.gpkg"
+        lines_gdf = gpd.read_file(OSM_RAW / "bus.gpkg", layer="lines")
+        lines_gdf.to_file(tmp, driver="GPKG", layer="lines")
+        combined.to_file(tmp, driver="GPKG", layer="stops")
+        shutil.move(str(tmp), str(OSM_RAW / "bus.gpkg"))
+        print(f"  -> bus stops layer now: {len(combined)} total", flush=True)
+
+    # stations (named station nodes/areas) for labelling
+    sql = stations_query(bbox)
+    (OSM_RAW / "stations.overpassql").write_text(sql)
+    print("[stations] querying Overpass...", flush=True)
+    sdata = overpass(sql)
+    st = parse_stations(sdata)
+    print(f"  -> {len(st)} named stations", flush=True)
+    if st:
+        gpd.GeoDataFrame(st, crs="EPSG:4326").to_file(
+            OSM_RAW / "stations.gpkg", driver="GPKG", layer="stations"
+        )
+    print("done.")
+
+
+if __name__ == "__main__":
+    main()

copenhagen/scripts/gtfs_to_geopackage.py

@@ -0,0 +1,42 @@
+#!/usr/bin/env python3
+"""Convert raw GTFS txt files into portable GeoPackage layers + a colour table.
+
+STATUS: BLOCKED on download_gtfs.py (no feed URL yet).
+
+Intended outputs (data/processed):
+    gtfs_shapes.gpkg   (shapes.txt -> LineString per shape_id, EPSG:4326)
+    gtfs_stops.gpkg    (stops.txt -> Point per stop, EPSG:4326)
+    route_colors.csv   (route_id, route_short_name, route_type, route_color)
+
+Uses partridge for fast GTFS parsing. prepare.py reads these when present and
+falls back to OSM-only data when absent.
+"""
+import sys
+from pathlib import Path
+
+from _common import GTFS_RAW, PROCESSED
+
+
+def main():
+    if not (GTFS_RAW / "routes.txt").exists():
+        print(
+            "gtfs_to_geopackage: BLOCKED — no GTFS data found in "
+            f"{GTFS_RAW}. Run download_gtfs.py first once a feed URL is set.",
+            file=sys.stderr,
+        )
+        sys.exit(2)
+
+    # TODO (when GTFS available):
+    #   import partridge as pt
+    #   import geopandas as gpd
+    #   from shapely.geometry import LineString, Point
+    #   feed = pt.load_geo_feed(str(GTFS_RAW))
+    #   shapes -> gtfs_shapes.gpkg (LineString per shape_id)
+    #   stops  -> gtfs_stops.gpkg
+    #   routes -> route_colors.csv (route_id, route_short_name, route_type, route_color)
+    print("gtfs_to_geopackage: not yet implemented (GTFS feed unavailable).")
+    sys.exit(2)
+
+
+if __name__ == "__main__":
+    main()

copenhagen/scripts/prepare.py

@@ -0,0 +1,205 @@
+#!/usr/bin/env python3
+"""Merge OSM layers into a single master GeoPackage, enriched for rendering.
+
+Reads data/raw/osm/{mode}.gpkg (lines + stops) and data/processed/area.gpkg,
+and writes data/processed/master.gpkg with layers "lines" and "stops".
+
+Enrichment:
+  - style      : from config/modes.yaml (mode -> style key)
+  - colour     : GTFS route_color (if present) -> OSM colour tag -> palette
+  - bus_category : A / C / S / regular (regex from styling.yaml)
+Lines and stops are clipped to the City Pass area polygon.
+"""
+import re
+from pathlib import Path
+
+import geopandas as gpd
+import pandas as pd
+import yaml
+
+from _common import CONFIG, OSM_RAW, PROCESSED
+
+MIXED_CRS_WARN = "mixed CRS"
+
+
+def load_modes():
+    return yaml.safe_load((CONFIG / "modes.yaml").read_text())["modes"]
+
+
+def load_styling():
+    return yaml.safe_load((CONFIG / "styling.yaml").read_text())
+
+
+def collect_lines(modes):
+    frames = []
+    for mode, cfg in modes.items():
+        p = OSM_RAW / f"{mode}.gpkg"
+        if not p.exists():
+            continue
+        g = gpd.read_file(p, layer="lines")
+        if g.empty:
+            continue
+        g["mode"] = mode
+        g["style"] = cfg["style"]
+        g["gtfs_route_type"] = cfg.get("gtfs_route_type")
+        frames.append(g)
+    if not frames:
+        return gpd.GeoDataFrame(columns=["mode", "ref", "style", "geometry"], crs="EPSG:4326")
+    return pd.concat(frames, ignore_index=True)
+
+
+def collect_stops(modes):
+    frames = []
+    for mode, cfg in modes.items():
+        p = OSM_RAW / f"{mode}.gpkg"
+        if not p.exists():
+            continue
+        try:
+            g = gpd.read_file(p, layer="stops")
+        except Exception:
+            continue
+        if g.empty:
+            continue
+        g["mode"] = mode
+        g["style"] = cfg["style"]
+        frames.append(g)
+    if not frames:
+        return None
+    return pd.concat(frames, ignore_index=True)
+
+
+def gtfs_colour_map():
+    """Return {(route_short_name, route_type): route_color} if GTFS data exists."""
+    csv = PROCESSED / "route_colors.csv"
+    if not csv.exists():
+        return {}
+    df = pd.read_csv(csv)
+    out = {}
+    for _, r in df.iterrows():
+        name = str(r.get("route_short_name") or "").strip()
+        if name:
+            out[(name, r.get("route_type"))] = r.get("route_color")
+    return out
+
+
+def categorise_bus(ref, patterns):
+    if ref is None:
+        return "regular"
+    s = str(ref)
+    for cat, pat in patterns.items():
+        if re.search(pat, s):
+            return cat
+    return "regular"
+
+
+def resolve_colour(row, palette, gtfs):
+    style = row["style"]
+    ref = row.get("ref")
+    colour = row.get("colour")
+    # 1. GTFS (if available)
+    if gtfs and isinstance(ref, str):
+        c = gtfs.get((ref, row.get("gtfs_route_type")))
+        if c:
+            return str(c)
+    # 2. OSM colour tag (only accept real hex/named strings, not NA/None)
+    if isinstance(colour, str) and colour.strip():
+        return colour
+    # 3. palette
+    p = palette.get(style, {})
+    if style == "bus":
+        cat = row.get("bus_category")
+        if not isinstance(cat, str):
+            cat = "regular"
+        return p.get(cat) or p.get("default")
+    key = ref if isinstance(ref, str) else None
+    return p.get(key) or p.get("default")
+
+
+def main():
+    modes = load_modes()
+    styling = load_styling()
+    palette = styling["palette"]
+    bus_patterns = styling["bus_filters"]["categories"]
+
+    area = gpd.read_file(PROCESSED / "area.gpkg")
+    area_geom = area.geometry.union_all()
+
+    lines = collect_lines(modes)
+    print(f"collected {len(lines)} raw lines across modes", flush=True)
+
+    # bus categorisation (only meaningful for bus style)
+    lines["bus_category"] = [
+        categorise_bus(r, bus_patterns) if s == "bus" else None
+        for r, s in zip(lines.get("ref"), lines["style"])
+    ]
+
+    gtfs = gtfs_colour_map()
+    if gtfs:
+        print(f"GTFS colours loaded: {len(gtfs)} routes", flush=True)
+    else:
+        print("no GTFS colour data; using OSM colour -> palette", flush=True)
+
+    lines["colour_final"] = [
+        resolve_colour(r, palette, gtfs) for _, r in lines.iterrows()
+    ]
+
+    # clip to area
+    lines = lines[~lines.geometry.isna() & lines.geometry.is_valid]
+    clipped = gpd.clip(lines, area_geom)
+    print(f"clipped to area: {len(clipped)} line features remain", flush=True)
+
+    out_cols = [
+        "mode", "style", "ref", "name", "network", "operator", "route",
+        "colour", "bus_category", "colour_final", "geometry",
+    ]
+    for c in out_cols:
+        if c not in clipped.columns:
+            clipped[c] = None
+    clipped = clipped.set_geometry("geometry")
+    clipped = gpd.GeoDataFrame(clipped[out_cols], crs="EPSG:4326")
+
+    PROCESSED.mkdir(parents=True, exist_ok=True)
+    # drop existing master.gpkg so layer overwrite is clean
+    out = PROCESSED / "master.gpkg"
+    if out.exists():
+        out.unlink()
+    clipped.to_file(out, driver="GPKG", layer="lines")
+    print(f"wrote lines layer: {len(clipped)} features", flush=True)
+
+    # stops (stop_positions along routes)
+    stops = collect_stops(modes)
+    if stops is not None and not stops.empty:
+        stops = stops[~stops.geometry.isna() & stops.geometry.is_valid]
+        stops_clipped = gpd.clip(stops, area_geom)
+        keep = ["mode", "style", "name", "ref", "public_transport", "railway", "route_refs", "geometry"]
+        for c in keep:
+            if c not in stops_clipped.columns:
+                stops_clipped[c] = None
+        stops_clipped = gpd.GeoDataFrame(stops_clipped[keep], crs="EPSG:4326")
+        stops_clipped.to_file(out, driver="GPKG", layer="stops")
+        print(f"wrote stops layer: {len(stops_clipped)} features", flush=True)
+    else:
+        print("no stops to write", flush=True)
+
+    # stations (named station nodes/areas) for labelling
+    st_path = OSM_RAW / "stations.gpkg"
+    if st_path.exists():
+        st = gpd.read_file(st_path, layer="stations")
+        if not st.empty:
+            st = st[~st.geometry.isna() & st.geometry.is_valid]
+            st_clipped = gpd.clip(st, area_geom)
+            st_clipped = gpd.GeoDataFrame(st_clipped, crs="EPSG:4326")
+            st_clipped.to_file(out, driver="GPKG", layer="stations")
+            print(f"wrote stations layer: {len(st_clipped)} named stations", flush=True)
+    else:
+        print("no stations file; labelling will be limited", flush=True)
+
+    # summary by style
+    print("\nsummary by style:")
+    for style, grp in clipped.groupby("style"):
+        refs = grp["ref"].dropna().unique()
+        print(f"  {style:10s}: {len(grp):4d} feats, {len(refs):3d} refs")
+
+
+if __name__ == "__main__":
+    main()

copenhagen/scripts/render.py

@@ -0,0 +1,509 @@
+#!/usr/bin/env python3
+"""Render the Copenhagen transit map to PNG/SVG/PDF.
+
+Reads data/processed/master.gpkg (layers: lines, stations) + area.gpkg +
+config/styling.yaml, applies CLI filters, and composes a printable map:
+  - Carto Positron (no labels) basemap via contextily (EPSG:3857)
+  - two-tone lines: dark casing + colour body, per route
+  - z-order: bus (bottom) -> s_tog -> light_rail -> regional -> metro (top)
+  - shapeburst fade mask outside the City Pass area
+  - station labels (metro + S-tog) via adjustText
+
+Usage:
+  uv run python render.py                       # default: all modes, PNG, markers only
+  uv run python render.py --labels              # add station name labels
+  uv run python render.py --no-stops            # no station markers
+  uv run python render.py --no-buses --format svg --out map.svg
+  uv run python render.py --bus-subset A,C,S     # only trunk buses
+  uv run python render.py --modes metro,s_tog     # rail-only
+"""
+import argparse
+from pathlib import Path
+
+import geopandas as gpd
+import matplotlib
+matplotlib.use("Agg")
+import matplotlib.pyplot as plt
+import numpy as np
+import yaml
+from affine import Affine
+from adjustText import adjust_text
+from scipy.ndimage import distance_transform_edt
+
+from _common import CONFIG, OUTPUT, PROCESSED, UA
+
+TARGET_CRS = "EPSG:3857"
+
+
+def load_args():
+    p = argparse.ArgumentParser(description="Render the Copenhagen transit map.")
+    p.add_argument("--modes", default=None,
+                   help="comma-separated styles to include (metro,s_tog,light_rail,regional,bus)")
+    p.add_argument("--no-buses", action="store_true", help="omit bus layer")
+    p.add_argument("--bus-subset", default=None,
+                   help="comma-separated bus categories to keep (A,C,S,regular)")
+    p.add_argument("--max-bus-routes", type=int, default=None,
+                   help="cap on number of distinct bus refs drawn")
+    p.add_argument("--dpi", type=int, default=200)
+    p.add_argument("--size", default="16",
+                   help="figure size in inches: 'W' (height auto) or 'WxH'")
+    p.add_argument("--format", default="png", help="output format (png/svg/pdf)")
+    p.add_argument("--out", default=None, help="output path (default: output/copenhagen.<fmt>)")
+    p.add_argument("--no-basemap", action="store_true", help="skip basemap tiles (offline/faster)")
+    p.add_argument("--labels", action="store_true", help="draw station name labels (off by default)")
+    p.add_argument("--no-stops", action="store_true", help="skip station markers")
+    return p.parse_args()
+
+
+def parse_size(spec, aspect):
+    if "x" in spec:
+        w, h = spec.split("x")
+        return float(w), float(h)
+    w = float(spec)
+    return w, w / aspect
+
+
+def filter_lines(lines, args):
+    keep = set(lines["style"].unique())
+    if args.modes:
+        keep &= {m.strip() for m in args.modes.split(",")}
+    if args.no_buses:
+        keep.discard("bus")
+    out = lines[lines["style"].isin(keep)].copy()
+
+    if "bus" in keep and out["style"].eq("bus").any():
+        bus = out[out["style"].eq("bus")]
+        other = out[~out["style"].eq("bus")]
+        if args.bus_subset:
+            cats = {c.strip() for c in args.bus_subset.split(",")}
+            bus = bus[bus["bus_category"].isin(cats)]
+        if args.max_bus_routes is not None:
+            refs = list(dict.fromkeys(bus["ref"].dropna()))
+            keep_refs = set(refs[:args.max_bus_routes])
+            bus = bus[bus["ref"].isin(keep_refs)]
+        out = gpd.GeoDataFrame(pd_concat([other, bus]), crs=out.crs)
+    return out
+
+
+def pd_concat(frames):
+    import pandas as pd
+    return pd.concat(frames, ignore_index=True)
+
+
+def add_basemap(ax, styling, zoom=None):
+    if zoom is None:
+        zoom = styling["basemap"]["zoom"]
+    try:
+        import contextily as cx
+        provider = styling["basemap"]["provider"]
+        src = cx.providers
+        for part in provider.split("."):
+            src = getattr(src, part)
+        cx.add_basemap(ax, crs=TARGET_CRS, source=src, zoom=zoom,
+                       attribution=False, zorder=0,
+                       headers={"User-Agent": UA})
+        return True
+    except Exception as e:
+        print(f"  [basemap] unavailable ({e}); using plain background", flush=True)
+        ax.set_facecolor("#f2f2f0")
+        return False
+
+
+def shapeburst_mask(ax, area_3857, extent, styling, res=1600):
+    fade_m = styling["mask"]["fade_distance_m"]
+    max_alpha = styling["mask"]["max_alpha"]
+    color = styling["mask"]["color"]
+    x0, x1, y0, y1 = extent
+    nx = res
+    ny = max(1, int(round(res * (y1 - y0) / (x1 - x0))))
+    pxw = (x1 - x0) / nx
+    pxh = (y1 - y0) / ny
+    transform = Affine.translation(x0, y1) * Affine.scale(pxw, -pxh)
+
+    from rasterio.features import rasterize
+    geom = area_3857.geometry.union_all()
+    mask = rasterize([(geom, 1)], out_shape=(ny, nx), transform=transform,
+                     fill=0, dtype=np.uint8, all_touched=False).astype(bool)
+    dist_px = distance_transform_edt(~mask)
+    dist_m = dist_px * pxw
+    alpha = np.where(mask, 0.0, np.clip(dist_m / fade_m, 0.0, 1.0)) * max_alpha
+
+    cmap = matplotlib.colors.to_rgba(color)
+    rgba = np.zeros((ny, nx, 4))
+    rgba[..., :3] = cmap[:3]
+    rgba[..., 3] = alpha
+    ax.imshow(rgba, extent=(x0, x1, y0, y1), origin="upper",
+              aspect="equal", interpolation="bilinear", zorder=0.5)
+
+
+def plot_lines(ax, lines, styling):
+    zord = styling["zorder"]
+    widths = styling["line_width"]
+    alphas = styling["alpha"]
+    outline_w = styling["outline_width"]
+    outline_c = styling["outline_color"]
+    for style in sorted(zord, key=lambda k: zord[k]):
+        grp = lines[lines["style"] == style]
+        if grp.empty:
+            continue
+        w = widths.get(style, 1.0)
+        z = zord[style]
+        # dark casing
+        grp.plot(ax=ax, color=outline_c, linewidth=w + outline_w,
+                 zorder=z, alpha=1.0, capstyle="round", joinstyle="round")
+        # colour body
+        body_colors = grp["colour_final"].tolist()
+        grp.plot(ax=ax, color=body_colors, linewidth=w,
+                 zorder=z + 0.1, alpha=alphas.get(style, 1.0),
+                 capstyle="round", joinstyle="round")
+
+
+def classify_station(r):
+    """Return the style key for a station row, or None."""
+    if r.get("subway") == "yes" or r.get("station") == "subway":
+        return "metro"
+    if r.get("light_rail") == "yes" or r.get("station") == "light_rail":
+        return "s_tog"
+    if r.get("railway") == "station" or r.get("train") == "yes":
+        return "regional"
+    return None
+
+
+DISTANCE_CRS = "EPSG:25832"  # UTM 32N — accurate meters for Copenhagen
+
+
+def cluster_stops(coords, threshold_m):
+    """Cluster (x, y) coordinates in EPSG:3857 within threshold_m (real meters).
+
+    Reprojects to UTM 32N for accurate distance computation, clusters, then
+    returns centroids in EPSG:3857 for plotting.
+    """
+    from pyproj import Transformer
+    from scipy.cluster.hierarchy import fcluster, linkage
+    from scipy.spatial.distance import pdist
+
+    coords = np.asarray(coords)
+    if len(coords) <= 1:
+        return coords.tolist() if len(coords) else []
+
+    # reproject to UTM for accurate distances
+    to_utm = Transformer.from_crs(TARGET_CRS, DISTANCE_CRS, always_xy=True)
+    back = Transformer.from_crs(DISTANCE_CRS, TARGET_CRS, always_xy=True)
+    utm = np.array([to_utm.transform(x, y) for x, y in coords])
+
+    dists = pdist(utm)
+    links = linkage(dists, method="single")
+    labels = fcluster(links, t=threshold_m, criterion="distance")
+
+    centroids = []
+    for label in set(labels):
+        members = utm[labels == label]
+        ux, uy = members.mean(axis=0)
+        cx, cy = back.transform(ux, uy)
+        centroids.append((cx, cy))
+    return centroids
+
+
+def plot_stations(ax, stations, styling, active_styles):
+    if stations is None or stations.empty:
+        return
+    cfg = styling.get("stations", {})
+    if not cfg.get("show", True):
+        return
+    marker_styles = set(cfg.get("styles", ["metro", "s_tog", "regional"]))
+    cluster_cfg = cfg.get("stop_cluster_m", {})
+    mk = cfg.get("marker", {})
+    sizes = mk.get("size", {})
+    shape = mk.get("shape", "circle")
+    fill = mk.get("fill", "white")
+    edge = mk.get("edge", "#2b2b2b")
+    lw = mk.get("linewidth", 0.8)
+    marker = "o" if shape == "circle" else "s"
+
+    # group coordinates by style, then by name for clustering
+    pts_by_style = {}
+    names_by_style = {}
+    for _, r in stations.iterrows():
+        s = classify_station(r)
+        if s is None or s not in marker_styles or s not in active_styles:
+            continue
+        coord = (r.geometry.x, r.geometry.y)
+        name = r.get("name") if isinstance(r.get("name"), str) else None
+        pts_by_style.setdefault(s, []).append(coord)
+        names_by_style.setdefault(s, {}).setdefault(name, []).append(coord)
+
+    zord = styling["zorder"]
+    top_z = max(zord.values()) + 1  # all markers above all lines
+    for s in sorted(pts_by_style, key=lambda k: zord.get(k, 0)):
+        threshold = cluster_cfg.get(s, 50) if isinstance(cluster_cfg, dict) else cluster_cfg
+        # cluster same-name stations to one marker; keep unnamed as-is
+        centroids = []
+        for name, coords in names_by_style[s].items():
+            if name and len(coords) > 1:
+                centroids.extend(cluster_stops(coords, threshold))
+            else:
+                centroids.extend(coords)
+        sz = sizes.get(s, 3.0)
+        ax.scatter([p[0] for p in centroids], [p[1] for p in centroids],
+                   s=sz ** 2, marker=marker, c=fill, edgecolors=edge,
+                   linewidths=lw, zorder=top_z, alpha=1.0)
+
+
+def rail_station_names(stations):
+    """Return the set of names for rail stations (metro/s_tog/regional)."""
+    if stations is None or stations.empty:
+        return set()
+    names = set()
+    for _, r in stations.iterrows():
+        if classify_station(r) is not None:
+            name = r.get("name")
+            if isinstance(name, str) and name.strip():
+                names.add(name)
+    return names
+
+
+def group_by_route_membership(rows, cluster_m):
+    """Group stop rows by (name, shared route_ref) into distinct stops.
+
+    Two same-name stops that share at least one route_ref are the same stop.
+    Stops on disjoint routes are different stops. Within each group,
+    platform duplicates are merged by distance clustering.
+    Returns a list of (x, y) centroids.
+    """
+    if rows is None or rows.empty:
+        return []
+
+    coords = [(r.geometry.x, r.geometry.y) for _, r in rows.iterrows()]
+    route_sets = []
+    for _, r in rows.iterrows():
+        rr = r.get("route_refs")
+        if isinstance(rr, str) and rr:
+            route_sets.append(set(rr.split(";")))
+        else:
+            route_sets.append(set())
+
+    n = len(coords)
+    # union-find: stops sharing a route are connected
+    parent = list(range(n))
+
+    def find(x):
+        while parent[x] != x:
+            parent[x] = parent[parent[x]]
+            x = parent[x]
+        return x
+
+    def union(a, b):
+        ra, rb = find(a), find(b)
+        if ra != rb:
+            parent[ra] = rb
+
+    for i in range(n):
+        for j in range(i + 1, n):
+            if route_sets[i] & route_sets[j]:
+                union(i, j)
+
+    # group indices by connected component
+    components = {}
+    for i in range(n):
+        root = find(i)
+        components.setdefault(root, []).append(i)
+
+    centroids = []
+    for indices in components.values():
+        comp_coords = [coords[i] for i in indices]
+        comp_routes = [route_sets[i] for i in indices]
+        has_routes = any(rs for rs in comp_routes)
+        if len(comp_coords) == 1:
+            centroids.append(comp_coords[0])
+        elif has_routes:
+            # stops sharing a route are the same stop; merge to centroid
+            centroids.append(tuple(np.mean(comp_coords, axis=0)))
+        else:
+            # no route info (orphans); fall back to distance clustering
+            centroids.extend(cluster_stops(comp_coords, cluster_m))
+    return centroids
+
+
+def plot_bus_stops(ax, stops, styling, active_styles, rail_names=None):
+    """Draw small markers for bus stops from the stops layer.
+
+    Bus stops whose name matches a rail station are skipped — the rail
+    station marker represents that stop. Stops are grouped by shared
+    route membership: same-name stops on disjoint routes are distinct.
+    """
+    if stops is None or stops.empty or "bus" not in active_styles:
+        return
+    cfg = styling.get("stations", {})
+    if not cfg.get("show", True) or "bus" not in set(cfg.get("styles", [])):
+        return
+    bus = stops[stops["style"] == "bus"]
+    if bus.empty:
+        return
+
+    if rail_names is None:
+        rail_names = set()
+
+    cluster_cfg = cfg.get("stop_cluster_m", {})
+    cluster_m = cluster_cfg.get("bus", 50) if isinstance(cluster_cfg, dict) else cluster_cfg
+    mk = cfg.get("marker", {})
+    sizes = mk.get("size", {})
+    sz = sizes.get("bus", 1.5)
+    fill = mk.get("fill", "white")
+    edge = mk.get("edge", "#2b2b2b")
+    lw = mk.get("linewidth", 0.8)
+    shape = mk.get("shape", "circle")
+    marker = "o" if shape == "circle" else "s"
+    top_z = max(styling["zorder"].values()) + 1
+
+    # group by name, then by route membership; skip rail station names
+    centroids = []
+    skipped = 0
+    for name, grp in bus.groupby("name"):
+        if name and name in rail_names:
+            skipped += len(grp)
+            continue
+        centroids.extend(group_by_route_membership(grp, cluster_m))
+
+    ax.scatter([p[0] for p in centroids], [p[1] for p in centroids],
+               s=sz ** 2, marker=marker, c=fill, edgecolors=edge,
+               linewidths=lw, zorder=top_z, alpha=0.8)
+
+
+def label_stations(ax, stations, styling, active_styles):
+    if stations is None or stations.empty:
+        return
+    cfg = styling.get("labels", {})
+    if not cfg.get("show", False):
+        return
+    label_styles = set(cfg.get("styles", ["metro", "s_tog"])) & active_styles
+    fmin = cfg.get("min_fontsize", 5)
+    fmax = cfg.get("max_fontsize", 9)
+
+    def is_metro(r):
+        return classify_station(r) == "metro"
+
+    def is_stog(r):
+        return classify_station(r) == "s_tog"
+
+    pts = []
+    for _, r in stations.iterrows():
+        s = classify_station(r)
+        if s not in label_styles:
+            continue
+        name = r.get("name")
+        if not isinstance(name, str) or not name.strip():
+            continue
+        pts.append((r.geometry.x, r.geometry.y, name, s == "metro"))
+
+    if not pts:
+        return
+    # dedupe by name (keep first location)
+    seen = {}
+    uniq = []
+    for x, y, name, is_m in pts:
+        if name in seen:
+            continue
+        seen[name] = True
+        uniq.append((x, y, name, is_m))
+
+    texts = []
+    for x, y, name, is_m in uniq:
+        fs = fmax if is_m else fmin
+        bbox = dict(boxstyle="round,pad=0.15", fc="white", ec="none", alpha=0.7)
+        t = ax.text(x, y, name, fontsize=fs, color="#111111",
+                    zorder=10, ha="center", va="center", bbox=bbox)
+        texts.append(t)
+    try:
+        adjust_text(texts, ax=ax, expand_points=(1.4, 1.6),
+                    force_text=(0.4, 0.6), lim=300,
+                    arrowprops=dict(arrowstyle="-", color="#888", lw=0.4))
+    except Exception as e:
+        print(f"  [labels] adjustText failed ({e}); labels may overlap", flush=True)
+
+
+def main():
+    args = load_args()
+    styling = yaml.safe_load((CONFIG / "styling.yaml").read_text())
+
+    area = gpd.read_file(PROCESSED / "area.gpkg").to_crs(TARGET_CRS)
+    lines = gpd.read_file(PROCESSED / "master.gpkg", layer="lines").to_crs(TARGET_CRS)
+    try:
+        stations = gpd.read_file(PROCESSED / "master.gpkg", layer="stations").to_crs(TARGET_CRS)
+    except Exception:
+        stations = None
+    try:
+        stops = gpd.read_file(PROCESSED / "master.gpkg", layer="stops").to_crs(TARGET_CRS)
+    except Exception:
+        stops = None
+
+    lines = filter_lines(lines, args)
+    print(f"rendering {len(lines)} line features", flush=True)
+
+    bounds = area.total_bounds  # (minx, miny, maxx, maxy)
+    dx = bounds[2] - bounds[0]
+    dy = bounds[3] - bounds[1]
+    aspect = dx / dy
+    pad = 0.12
+    extent = (bounds[0] - pad * dx, bounds[2] + pad * dx,
+              bounds[1] - pad * dy, bounds[3] + pad * dy)
+
+    fig_w, fig_h = parse_size(args.size, aspect)
+    fig, ax = plt.subplots(figsize=(fig_w, fig_h), dpi=args.dpi)
+    ax.set_xlim(extent[0], extent[1])
+    ax.set_ylim(extent[2], extent[3])
+    ax.set_aspect("equal")
+    ax.axis("off")
+
+    # basemap
+    if not args.no_basemap:
+        add_basemap(ax, styling)
+
+    # fade mask outside area
+    shapeburst_mask(ax, area, extent, styling)
+
+    # lines
+    plot_lines(ax, lines, styling)
+
+    # active styles (which mode layers are actually drawn)
+    active_styles = set(lines["style"].unique())
+
+    # station markers (symbols, no text)
+    if not args.no_stops:
+        plot_stations(ax, stations, styling, active_styles)
+        rail_names = rail_station_names(stations)
+        plot_bus_stops(ax, stops, styling, active_styles, rail_names=rail_names)
+
+    # station labels (opt-in via --labels)
+    if args.labels:
+        styling["labels"]["show"] = True
+    if styling.get("labels", {}).get("show", False):
+        label_stations(ax, stations, styling, active_styles)
+
+    # title + attribution
+    ax.text(0.5, 0.985, "Kurragömma Köpenhamn",
+            transform=ax.transAxes, ha="center", va="top",
+            fontsize=fig_w * 1.1, fontweight="bold", color="#222",
+            bbox=dict(boxstyle="round,pad=0.3", fc="white", ec="none", alpha=0.7))
+    ax.text(0.01, 0.01,
+            "Data: © OpenStreetMap contributors (ODbL) · Base: Esri, HERE",
+            transform=ax.transAxes, ha="left", va="bottom",
+            fontsize=max(5, fig_w * 0.5), color="#666")
+
+    out_dir = OUTPUT
+    out_dir.mkdir(parents=True, exist_ok=True)
+    if args.out:
+        out_path = Path(args.out)
+        if not out_path.is_absolute():
+            out_path = out_dir / out_path
+    else:
+        out_path = out_dir / f"copenhagen.{args.format}"
+    out_path.parent.mkdir(parents=True, exist_ok=True)
+
+    plt.savefig(out_path, dpi=args.dpi, bbox_inches="tight",
+                pad_inches=0.2, facecolor="white")
+    plt.close(fig)
+    print(f"wrote {out_path} ({args.format}, {fig_w}x{fig_h}in @ {args.dpi}dpi)")
+
+
+if __name__ == "__main__":
+    main()

pyproject.toml

@@ -0,0 +1,21 @@
+[project]
+name = "jetlag-maps"
+version = "0.1.0"
+description = "Scriptable transit maps (Copenhagen) built from OSM + GTFS."
+requires-python = ">=3.12"
+dependencies = [
+    "geopandas>=1.1",
+    "shapely>=2.0",
+    "matplotlib>=3.9",
+    "contextily>=1.6",
+    "pyogrio>=0.10",
+    "partridge>=1.1",
+    "requests>=2.31",
+    "pyyaml>=6.0",
+    "adjusttext>=1.2",
+    "rtree>=1.0",
+    "scipy>=1.13",
+]
+
+[tool.uv]
+package = false