Source Copenhagen transit data from Rejseplanen GTFS

- Replace OSM Overpass transit data with the Rejseplanen GTFS feed
  (routes keyed by (agency, short name); styles from modes.yaml)
- Draw one shape per (style, ref, direction), choosing the shape that
  serves the most in-area stops so lines pass the stops we show
- Collapse stops by names (verified unambiguous); prune stops whose
  serving refs have no drawn line within 300 m
- Patch Københavns Havn + Nordhavn into the area polygon so ferry
  routes and sub-harbour metro tunnels survive clipping; 100 m buffer
  closes relation boundary slivers
- Drop OSM download pipeline; keep tiled basemap

.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,26 @@
+{
+  "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}
+  ],
+  "waterways": [
+    {
+      "name": "Kobenhavns Havn + Nordhavn",
+      "comment": "The harbour channel and Nordhavn basin belong to the City Pass zone in any practical sense: the metro tunnels under it (M1/M2/M4), harbour ferries 991/992 sail it, and the fade-out should not start mid-harbour. The kommune boundary relations exclude all water, so this polygon is patched in.",
+      "ring": [
+        [12.532, 55.640],
+        [12.545, 55.666],
+        [12.565, 55.680],
+        [12.585, 55.715],
+        [12.612, 55.717],
+        [12.632, 55.685],
+        [12.632, 55.660],
+        [12.600, 55.636],
+        [12.560, 55.634]
+      ]
+    }
+  ]
+}

copenhagen/config/modes.yaml

@@ -0,0 +1,34 @@
+# Mapping of map styles to GTFS agency / route_type filters.
+#
+# Used by prepare.py to classify GTFS routes (lines layer) and the routes
+# serving each stop pole (stops/stations layers). Style keys are referenced
+# by styling.yaml (palette, zorder, widths).
+#
+# `agencies` matches GTFS agency_name from agency.txt.
+# `route_types` optionally restricts GTFS route_type values, including the
+# extended codes used in this feed (109 = S-tog suburban rail, 700/715 =
+# Movia bus service types). Without `route_types`, all of the agency's
+# routes match (e.g. Movia including harbour ferries 991/992).
+#
+# ORDER MATTERS: the styles listed here are also the priority order when a
+# pole is served by several modes (e.g. S-tog and regional trains share
+# platforms at hub stations — the higher style wins).
+modes:
+  metro:
+    agencies: [Metroselskabet]
+  s_tog:
+    agencies: [DSB S-tog]
+    route_types: [109, 2]
+  light_rail:
+    agencies: [Hovedstadens Letbane]
+  regional:
+    agencies: [DSB, Lokaltog A/S, Skånetrafiken, Snälltåget AB, DSB Vores Tog]
+    route_types: [2]
+  ferry:
+    # harbour buses 991/992 (route_type 4); their shapes cross the harbour
+    # water and are intentionally NOT clipped to the land area polygon
+    agencies: [Movia]
+    route_types: [4]
+  bus:
+    agencies: [Movia]
+    route_types: [3, 700, 715]

copenhagen/config/styling.yaml

@@ -0,0 +1,99 @@
+# 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
+  ferry: 2
+  s_tog: 3
+  light_rail: 4
+  regional: 5
+  metro: 6
+
+line_width:           # body width in points (outline = body + outline_width)
+  bus: 0.8
+  ferry: 1.4
+  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
+  ferry: 1.0
+  s_tog: 0.95
+  light_rail: 0.95
+  regional: 0.9
+  metro: 1.0
+
+# Line colours. The Rejseplanen GTFS feed leaves route_color empty for the
+# Copenhagen operators, so in practice everything resolves to this palette.
+# Metro / S-tog values mirror the official line colours; regional trains and
+# the light rail have no operator 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"
+  ferry:
+    default: "#2aa8e0"   # harbour-bus blue
+
+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]
+  # NOTE: stops/stations arrive pre-collapsed from prepare.py (name-keyed,
+  # one point per stop/station) — no clustering configuration lives here.
+  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,13 @@
+"""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"
+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,168 @@
+#!/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, dissolves the union, and patches in any
+config/area.json "waterways" polygons (the harbour is excluded from the
+kommune boundaries but belongs to the City Pass zone in practice).
+
+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)
+
+    # patch waterways into the zone: kommune boundaries exclude water, but
+    # the harbour is practically part of the City Pass area (metro tunnels
+    # beneath it, harbour ferries sail it). See area.json -> waterways.
+    for w in area_cfg.get("waterways", []):
+        wp = Polygon(w["ring"])
+        print(f"adding waterway: {w['name']}", flush=True)
+        dissolved = unary_union([dissolved, wp])
+
+    # Small outward buffer (100 m) to close boundary slivers: the three
+    # relation outlines don't abut perfectly, leaving metre-wide cracks that
+    # otherwise fragment lines clipped to the area (bridge nicks included).
+    dissolved = (
+        gpd.GeoSeries([dissolved], crs="EPSG:4326")
+        .to_crs("EPSG:25832").buffer(100)
+        .to_crs("EPSG:4326").iloc[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,32 @@
+#!/usr/bin/env python3
+"""Download the Rejseplanen GTFS feed (all of Denmark) and unzip into data/raw/gtfs/.
+
+Source: https://www.rejseplanen.info/labs — a static GTFS zip published by
+Rejseplanen covering DSB, DSB S-tog, Metroselskabet, Movia, Lokaltog,
+Hovedstadens Letbane, Skånetrafiken, etc. Nationwide feed; gtfs_to_geopackage.py
+filters it down to the Copenhagen area.
+"""
+import urllib.request
+import zipfile
+
+from _common import GTFS_RAW, UA
+
+GTFS_URL = "https://www.rejseplanen.info/labs/GTFS.zip"
+
+
+def main():
+    GTFS_RAW.mkdir(parents=True, exist_ok=True)
+    zip_path = GTFS_RAW / "feed.zip"
+    print(f"downloading {GTFS_URL} -> {zip_path} ...", flush=True)
+    req = urllib.request.Request(GTFS_URL, headers={"User-Agent": UA})
+    with urllib.request.urlopen(req) as r, zip_path.open("wb") as f:
+        while chunk := r.read(1 << 20):
+            f.write(chunk)
+    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/gtfs_to_geopackage.py

@@ -0,0 +1,151 @@
+#!/usr/bin/env python3
+"""Convert raw GTFS txt files into portable GeoPackage layers + a colour table.
+
+Inputs (data/raw/gtfs): routes.txt, trips.txt, stops.txt, shapes.txt
+Outputs (data/processed):
+    gtfs_shapes.gpkg   (shapes.txt -> LineString per shape_id, EPSG:4326,
+                        each shape tagged with its route_id)
+    gtfs_stops.gpkg    (stops.txt -> Point per stop, EPSG:4326)
+    route_colors.csv   (route_id, agency_id, route_short_name, route_long_name,
+                        route_type, route_color, route_text_color)
+
+The feed is nationwide, so everything is pre-filtered to the Copenhagen area:
+stops by location, shapes to those intersecting the area polygon (or a
+fallback bbox when data/processed/area.gpkg does not exist yet), and routes
+to those having at least one surviving shape. prepare.py clips precisely to
+the area later.
+
+Note: this feed leaves route_color empty for the Copenhagen operators, so
+prepare.py will mostly fall through to OSM colour tags / the styling palette.
+stop_times.txt is not needed here and is intentionally not parsed (220 MB).
+"""
+import sys
+
+import geopandas as gpd
+import pandas as pd
+from shapely.geometry import LineString, box
+from shapely.geometry.base import BaseGeometry
+
+from _common import GTFS_RAW, PROCESSED
+
+# Generous Copenhagen bbox (covers København, Frederiksberg, Amager, and
+# immediate surroundings: airport, Hellerup, Lyngby, Brøndby, ...).
+CPH_BBOX = box(12.30, 55.55, 12.75, 55.85)
+
+
+def load_area() -> BaseGeometry:
+    """City Pass area polygon if build_area has run, else the fallback bbox."""
+    area_path = PROCESSED / "area.gpkg"
+    if area_path.exists():
+        return gpd.read_file(area_path).geometry.union_all()
+    print(f"no {area_path}; using fallback Copenhagen bbox", flush=True)
+    return CPH_BBOX
+
+
+def build_stops(area: BaseGeometry) -> gpd.GeoDataFrame:
+    stops = pd.read_csv(
+        GTFS_RAW / "stops.txt",
+        usecols=["stop_id", "stop_name", "stop_lat", "stop_lon",
+                 "location_type", "parent_station"],
+        dtype={"stop_id": str, "parent_station": str},
+    )
+    minx, miny, maxx, maxy = area.bounds
+    in_bbox = (
+        stops["stop_lat"].between(miny, maxy)
+        & stops["stop_lon"].between(minx, maxx)
+    )
+    stops = stops[in_bbox].copy()
+    g = gpd.GeoDataFrame(
+        stops,
+        geometry=gpd.points_from_xy(stops["stop_lon"], stops["stop_lat"]),
+        crs="EPSG:4326",
+    )
+    return g[g.intersects(area)].drop(columns=["stop_lat", "stop_lon"])
+
+
+def build_shapes(area: BaseGeometry) -> gpd.GeoDataFrame:
+    pts = pd.read_csv(GTFS_RAW / "shapes.txt", dtype={"shape_id": str})
+    minx, miny, maxx, maxy = area.bounds
+    hits = pts[
+        pts["shape_pt_lat"].between(miny, maxy)
+        & pts["shape_pt_lon"].between(minx, maxx)
+    ]
+    pts = pts[pts["shape_id"].isin(set(hits["shape_id"].unique()))]
+    pts = pts.sort_values(["shape_id", "shape_pt_sequence"])
+    lines = (
+        pts.groupby("shape_id")
+        .apply(
+            lambda g: LineString(zip(g["shape_pt_lon"], g["shape_pt_lat"])),
+            include_groups=False,
+        )
+        .rename("geometry")
+        .reset_index()
+    )
+    g = gpd.GeoDataFrame(lines, geometry="geometry", crs="EPSG:4326")
+    return g[g.intersects(area)]
+
+
+def build_routes(
+    shapes: gpd.GeoDataFrame,
+) -> tuple[gpd.GeoDataFrame, pd.DataFrame]:
+    """Tag shapes with route_id and keep routes having >= 1 surviving shape."""
+    trips = pd.read_csv(
+        GTFS_RAW / "trips.txt",
+        usecols=["route_id", "shape_id"],
+        dtype=str,
+    ).dropna(subset=["shape_id"])
+    shape_to_route = (
+        trips[trips["shape_id"].isin(set(shapes["shape_id"]))]
+        .drop_duplicates("shape_id")
+        .set_index("shape_id")["route_id"]
+    )
+    shapes = shapes.copy()
+    shapes["route_id"] = shapes["shape_id"].map(shape_to_route)
+
+    routes = pd.read_csv(GTFS_RAW / "routes.txt", dtype=str).fillna("")
+    routes = routes[
+        routes["route_id"].isin(set(shape_to_route.unique()))
+    ].copy()
+    return shapes, routes[
+        [
+            "route_id", "agency_id", "route_short_name", "route_long_name",
+            "route_type", "route_color", "route_text_color",
+        ]
+    ]
+
+
+def main():
+    if not (GTFS_RAW / "routes.txt").exists():
+        print(
+            "gtfs_to_geopackage: no GTFS data found in "
+            f"{GTFS_RAW}. Run download_gtfs.py first.",
+            file=sys.stderr,
+        )
+        sys.exit(2)
+
+    PROCESSED.mkdir(parents=True, exist_ok=True)
+    area = load_area()
+
+    stops = build_stops(area)
+    print(f"stops in area: {len(stops)}", flush=True)
+
+    shapes = build_shapes(area)
+    print(f"shapes in area: {len(shapes)}", flush=True)
+
+    shapes, routes = build_routes(shapes)
+    print(f"routes in area: {len(routes)}", flush=True)
+
+    stops.to_file(PROCESSED / "gtfs_stops.gpkg", driver="GPKG", layer="stops")
+    shapes.to_file(
+        PROCESSED / "gtfs_shapes.gpkg", driver="GPKG", layer="shapes"
+    )
+    routes.to_csv(PROCESSED / "route_colors.csv", index=False)
+    print(
+        "wrote gtfs_stops.gpkg, gtfs_shapes.gpkg, route_colors.csv "
+        f"into {PROCESSED}",
+        flush=True,
+    )
+
+
+if __name__ == "__main__":
+    main()

copenhagen/scripts/prepare.py

@@ -0,0 +1,320 @@
+#!/usr/bin/env python3
+"""Merge GTFS-derived layers into a single master GeoPackage, enriched for rendering.
+
+Reads (from data/processed):
+    gtfs_shapes.gpkg   route-tagged line geometries (gtfs_to_geopackage.py)
+    gtfs_stops.gpkg    raw stop poles (one row per physical pole)
+    route_colors.csv   route metadata (agency_id, short/long name, type, colour)
+    area.gpkg          City Pass area polygon (build_area.py, incl. waterways)
+and from data/raw/gtfs: agency.txt, trips.txt, stop_times.txt.
+
+Writes data/processed/master.gpkg with layers:
+    lines     one feature per (style, ref, direction) — the shape that serves
+              the most in-area stops (so drawn lines pass the drawn stops)
+    stops     bus/ferry stops, one point per stop name, pruned to stops
+              actually served by a drawn line (see STOP_LINE_MARGIN_M)
+    stations  rail-family stations, one point per (name, style)
+
+Stop names in the Rejseplanen feed are unambiguous per location (verified:
+max spread within a name is ~350 m), so every merge is key-based on names —
+no distance-based clustering anywhere.
+
+Style classification comes from config/modes.yaml (agency + route_type);
+file order is the priority when a pole is served by several modes.
+Colours: GTFS route_color -> styling.yaml palette (the feed leaves
+route_color empty for the Copenhagen operators, so the palette wins).
+"""
+import re
+import sys
+
+import geopandas as gpd
+import pandas as pd
+import yaml
+
+from _common import CONFIG, GTFS_RAW, PROCESSED
+
+LENGTH_CRS = "EPSG:25832"  # UTM 32N — metres, for length/centroid computations
+STATION_STYLES = {"metro", "s_tog", "light_rail", "regional"}
+# A stop is drawn only if a drawn line of one of its serving refs passes
+# within this distance. Routes have variants; without this check a stop can
+# end up further from the map than the variant we chose not to draw.
+STOP_LINE_MARGIN_M = 300
+
+
+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 load_area():
+    area = gpd.read_file(PROCESSED / "area.gpkg")
+    return area.geometry.union_all()
+
+
+def classify(agency, route_type, modes):
+    """Map a GTFS (agency_name, route_type) pair to a style key, or None."""
+    for style, cfg in modes.items():
+        if agency not in cfg.get("agencies", []):
+            continue
+        types = cfg.get("route_types")
+        if types and route_type not in types:
+            continue
+        return style
+    return None
+
+
+def load_routes(modes):
+    """route_colors.csv joined with agency names and classified by style."""
+    agencies = pd.read_csv(GTFS_RAW / "agency.txt", dtype=str)
+    agency_names = dict(zip(agencies["agency_id"], agencies["agency_name"]))
+    routes = pd.read_csv(PROCESSED / "route_colors.csv", dtype=str)
+    routes["route_type"] = pd.to_numeric(routes["route_type"], errors="coerce")
+    routes["agency_name"] = routes["agency_id"].map(agency_names)
+    routes["style"] = [
+        classify(a, t, modes)
+        for a, t in zip(routes["agency_name"], routes["route_type"])
+    ]
+    return routes
+
+
+def categorise_bus(ref, patterns):
+    s = "" if ref is None else str(ref)
+    for cat, pat in patterns.items():
+        if re.search(pat, s):
+            return cat
+    return "regular"
+
+
+def resolve_colour(style, ref, bus_category, gtfs_colour, palette):
+    # 1. GTFS route_color (mostly empty in this feed)
+    if isinstance(gtfs_colour, str) and gtfs_colour.strip():
+        return gtfs_colour
+    # 2. palette
+    p = palette.get(style, {})
+    if style == "bus":
+        return p.get(bus_category) or p.get("default")
+    return p.get(ref) or p.get("default")
+
+
+def build_lines(routes, styling, area_geom, trips, st, area_stop_ids):
+    """One feature per (style, ref, direction).
+
+    Geometry per group: the shape that serves the most in-area stops
+    (tie-break: longest). Picking by length alone can draw a variant that
+    skips stops shown on the map (terminal stubs, short-turn branches).
+    """
+    shapes = gpd.read_file(PROCESSED / "gtfs_shapes.gpkg")
+
+    shape_rows = trips.dropna(subset=["shape_id"]).drop_duplicates("shape_id")
+    shape_direction = dict(zip(shape_rows["shape_id"], shape_rows["direction_id"]))
+
+    # coverage: distinct in-area stops visited per shape
+    trip_shape = dict(zip(trips["trip_id"], trips["shape_id"]))
+    sv = st.assign(shape_id=st["trip_id"].map(trip_shape))
+    sv = sv[sv["stop_id"].isin(area_stop_ids)]
+    coverage = sv.groupby("shape_id")["stop_id"].nunique()
+
+    use = routes.dropna(subset=["style"])[
+        ["route_id", "route_short_name", "route_long_name",
+         "route_color", "style"]
+    ]
+    g = shapes.merge(use, on="route_id", how="inner")
+    g["direction_id"] = g["shape_id"].map(shape_direction)
+
+    g = g.to_crs(LENGTH_CRS)
+    g["_len"] = g.geometry.length
+    g["_cov"] = g["shape_id"].map(coverage).fillna(0)
+    g = (
+        g.sort_values(["_cov", "_len"], ascending=False)
+        .drop_duplicates(["style", "route_short_name", "direction_id"])
+        .to_crs("EPSG:4326")
+    )
+
+    g["ref"] = g["route_short_name"]
+    g["name"] = [
+        ln if isinstance(ln, str) and ln.strip() else ref
+        for ln, ref in zip(g["route_long_name"], g["ref"])
+    ]
+    patterns = styling["bus_filters"]["categories"]
+    g["bus_category"] = [
+        categorise_bus(ref, patterns) if style == "bus" else None
+        for ref, style in zip(g["ref"], g["style"])
+    ]
+    palette = styling["palette"]
+    g["colour_final"] = [
+        resolve_colour(style, ref, cat, gtfs_c, palette)
+        for style, ref, cat, gtfs_c
+        in zip(g["style"], g["ref"], g["bus_category"], g["route_color"])
+    ]
+
+    g = g[["ref", "name", "style", "bus_category", "colour_final", "geometry"]]
+    g = g[~g.geometry.isna() & g.geometry.is_valid]
+    # the area polygon includes patched-in waterways (area.json), so ferry
+    # shapes and sub-harbour metro tunnels survive the clip unfragmented
+    return gpd.clip(g, area_geom)
+
+
+def build_pole_classes(routes, modes, trips, st, area_stop_ids):
+    """Per-pole classification and serving refs, from stop_times.
+
+    Returns (pole_styles, pole_refs):
+      pole_styles[stop_id] = highest-priority style among classified serving
+                             routes (modes.yaml order)
+      pole_refs[stop_id]   = set of serving route short names (classified)
+    """
+    classified = routes.dropna(subset=["style"])
+    route_style = dict(zip(classified["route_id"], classified["style"]))
+    route_ref = dict(
+        zip(classified["route_id"], classified["route_short_name"])
+    )
+    trip_route = dict(zip(trips["trip_id"], trips["route_id"]))
+
+    sv = st[st["stop_id"].isin(area_stop_ids)].copy()
+    sv["route_id"] = sv["trip_id"].map(trip_route)
+    sv = sv.dropna(subset=["route_id"])
+    sv["style"] = sv["route_id"].map(route_style)
+    sv = sv.dropna(subset=["style"])
+
+    order = {s: i for i, s in enumerate(modes)}
+    sv["_p"] = sv["style"].map(order)
+    best = sv.loc[sv.groupby("stop_id")["_p"].idxmin(), ["stop_id", "style"]]
+    pole_styles = dict(zip(best["stop_id"], best["style"]))
+
+    pole_refs = sv.groupby("stop_id")["route_id"].apply(
+        lambda ids: {route_ref[i] for i in ids}
+    ).to_dict()
+    return pole_styles, pole_refs
+
+
+def build_stops_and_stations(pole_styles, pole_refs, area_geom, lines):
+    poles = gpd.read_file(PROCESSED / "gtfs_stops.gpkg")
+    poles["style"] = poles["stop_id"].map(pole_styles)
+    poles["refs"] = poles["stop_id"].map(lambda s: pole_refs.get(s, set()))
+    poles = poles.dropna(subset=["style"])
+
+    # project for accurate centroids and distances
+    poles = poles.to_crs(LENGTH_CRS)
+
+    # stops (everything not rail-family: buses + harbour ferries), one row
+    # per name, union of serving refs across its poles
+    not_rail = poles[~poles["style"].isin(STATION_STYLES)]
+    stop_rows = []
+    for (name, style), grp in not_rail.groupby(["stop_name", "style"]):
+        refs = set().union(*grp["refs"]) if len(grp) else set()
+        stop_rows.append({
+            "name": name,
+            "style": style,
+            "n_poles": len(grp),
+            "refs": refs,
+            "geometry": grp.geometry.union_all().centroid,
+        })
+
+    # stations: one row per (name, style)
+    rail = poles[poles["style"].isin(STATION_STYLES)]
+    station_rows = []
+    for (name, style), grp in rail.groupby(["stop_name", "style"]):
+        station_rows.append({
+            "name": name,
+            "style": style,
+            "n_poles": len(grp),
+            "geometry": grp.geometry.union_all().centroid,
+        })
+
+    area_25832 = gpd.GeoSeries([area_geom], crs="EPSG:4326").to_crs(LENGTH_CRS).union_all()
+
+    stops = gpd.GeoDataFrame(stop_rows, crs=LENGTH_CRS)
+
+    # prune stops not served by any drawn line within STOP_LINE_MARGIN_M
+    if len(stops) and len(lines):
+        line_geom = (
+            lines.to_crs(LENGTH_CRS)
+            .groupby("ref")["geometry"]
+            .agg(lambda g: g.union_all())
+            .to_dict()
+        )
+        keep = []
+        for _, r in stops.iterrows():
+            dists = [
+                geom.distance(r.geometry)
+                for rf in r["refs"]
+                if (geom := line_geom.get(rf)) is not None
+            ]
+            keep.append(bool(dists) and min(dists) <= STOP_LINE_MARGIN_M)
+        n_dropped = (~pd.Series(keep, index=stops.index)).sum()
+        if n_dropped:
+            print(f"pruned {n_dropped} stops not served by a drawn line "
+                  f"(>{STOP_LINE_MARGIN_M} m from nearest)", flush=True)
+        stops = stops[keep].drop(columns=["refs"])
+    else:
+        stops = stops.drop(columns=["refs"])
+
+    stops = gpd.clip(stops, area_25832)
+    stations = gpd.clip(
+        gpd.GeoDataFrame(station_rows, crs=LENGTH_CRS), area_25832
+    )
+    return stops.to_crs("EPSG:4326"), stations.to_crs("EPSG:4326")
+
+
+def main():
+    for dep in ("gtfs_shapes.gpkg", "gtfs_stops.gpkg", "route_colors.csv"):
+        if not (PROCESSED / dep).exists():
+            print(f"prepare: missing {PROCESSED / dep} — "
+                  "run gtfs_to_geopackage.py first.", file=sys.stderr)
+            sys.exit(2)
+
+    modes = load_modes()
+    styling = load_styling()
+    area_geom = load_area()
+
+    routes = load_routes(modes)
+    print(f"routes classified: {routes['style'].notna().sum()} of "
+          f"{len(routes)} map to a style", flush=True)
+
+    area_stop_ids = set(
+        gpd.read_file(PROCESSED / "gtfs_stops.gpkg")["stop_id"]
+    )
+    trips = pd.read_csv(
+        GTFS_RAW / "trips.txt",
+        usecols=["trip_id", "route_id", "shape_id", "direction_id"],
+        dtype=str,
+    )
+    st = pd.read_csv(
+        GTFS_RAW / "stop_times.txt",
+        usecols=["trip_id", "stop_id"],
+        dtype=str,
+    )
+
+    lines = build_lines(routes, styling, area_geom, trips, st, area_stop_ids)
+    print(f"lines: {len(lines)} features (one per style/ref/direction)",
+          flush=True)
+
+    pole_styles, pole_refs = build_pole_classes(
+        routes, modes, trips, st, area_stop_ids
+    )
+    stops, stations = build_stops_and_stations(
+        pole_styles, pole_refs, area_geom, lines
+    )
+    print(f"stops: {len(stops)} (bus + ferry); stations: {len(stations)}",
+          flush=True)
+
+    out = PROCESSED / "master.gpkg"
+    if out.exists():
+        out.unlink()
+    lines.to_file(out, driver="GPKG", layer="lines")
+    print(f"wrote lines layer: {len(lines)} features", flush=True)
+    stops.to_file(out, driver="GPKG", layer="stops")
+    print(f"wrote stops layer: {len(stops)} features", flush=True)
+    stations.to_file(out, driver="GPKG", layer="stations")
+    print(f"wrote stations layer: {len(stations)} features", flush=True)
+
+    print("\nsummary by style:")
+    for style, grp in lines.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,375 @@
+#!/usr/bin/env python3
+"""Render the Copenhagen transit map to PNG/SVG/PDF.
+
+Reads data/processed/master.gpkg (layers: lines, stops, stations — all
+GTFS-derived by prepare.py) + area.gpkg + config/styling.yaml, applies CLI
+filters, and composes a printable map:
+  - Carto/Esri grey 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
+
+All transit data comes from the Rejseplanen GTFS feed. Stop/station layers
+are pre-collapsed (one point per stop name / station) by prepare.py, so
+rendering is pure plotting — no clustering or merging here.
+
+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 pandas as pd
+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], ignore_index=True),
+                               crs=out.crs)
+    return out
+
+
+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")
+
+
+RAIL_STYLES = ("metro", "s_tog", "light_rail", "regional")
+
+
+def classify_station(r):
+    """Style key for a station row, or None. Stations are pre-classified."""
+    s = r.get("style")
+    return s if s in RAIL_STYLES else None
+
+
+def plot_stations(ax, stations, styling, active_styles):
+    """Draw station markers. Input is pre-collapsed: one row per (name, style)."""
+    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"]))
+    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"
+
+    pts_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
+        pts_by_style.setdefault(s, []).append((r.geometry.x, r.geometry.y))
+
+    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)):
+        pts = pts_by_style[s]
+        sz = sizes.get(s, 3.0)
+        ax.scatter([p[0] for p in pts], [p[1] for p in pts],
+                   s=sz ** 2, marker=marker, c=fill, edgecolors=edge,
+                   linewidths=lw, zorder=top_z, alpha=1.0)
+
+
+def rail_station_names(stations):
+    """Set of names of rail stations (metro/s_tog/light_rail/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 plot_bus_stops(ax, stops, styling, active_styles, rail_names=None):
+    """Draw small markers for bus stops (pre-collapsed: one point per name).
+
+    A bus stop whose name exactly matches a rail station is skipped — the
+    station marker represents it.
+    """
+    if stops is None or stops.empty or not ({"bus", "ferry"} & 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"].isin(["bus", "ferry"])]
+    if bus.empty:
+        return
+
+    rail_names = rail_names or set()
+    mk = cfg.get("marker", {})
+    sz = mk.get("size", {}).get("bus", 1.5)
+    fill = mk.get("fill", "white")
+    edge = mk.get("edge", "#2b2b2b")
+    lw = mk.get("linewidth", 0.8)
+    marker = "o" if mk.get("shape", "circle") == "circle" else "s"
+    top_z = max(styling["zorder"].values()) + 1
+
+    pts = [(r.geometry.x, r.geometry.y)
+           for _, r in bus.iterrows()
+           if not (isinstance(r.get("name"), str) and r.get("name") in rail_names)]
+    ax.scatter([p[0] for p in pts], [p[1] for p in pts],
+               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)
+
+    # stations are unique per (name, style); dedupe by name for labelling
+    seen = set()
+    uniq = []
+    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() or name in seen:
+            continue
+        seen.add(name)
+        uniq.append((r.geometry.x, r.geometry.y, name, s == "metro"))
+
+    if not uniq:
+        return
+    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,
+            "Transit: Rejseplanen GTFS · Area: © 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,20 @@
+[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",
+    "requests>=2.31",
+    "pyyaml>=6.0",
+    "adjusttext>=1.2",
+    "rtree>=1.0",
+    "scipy>=1.13",
+]
+
+[tool.uv]
+package = false