fetch_map_data.py

#!/usr/bin/env python3
"""
Fetch Kyiv metro lines/stations and building footprints for bird.rent listings
from the Overpass API (OpenStreetMap data).

Outputs map_data.js with JS constants for the HTML map:
  METRO_LINES        - GeoJSON FeatureCollection of LineString features
  METRO_STATIONS     - GeoJSON FeatureCollection of Point features (station entrances)
  BUILDING_FOOTPRINTS - GeoJSON FeatureCollection of Polygon features (matched buildings)
  LISTING_POINTS     - GeoJSON FeatureCollection of Point features (one per building, for low-zoom dots)

Usage:
  pip install requests
  python fetch_map_data.py

Takes ~5-15 minutes on first run (batched Overpass queries for ~3700 buildings).
Re-run any time the listing data or OSM data changes.
"""

import json
import math
import time
import sys
from collections import defaultdict

import requests

# Try multiple mirrors in order; the main server sometimes blocks default UAs
OVERPASS_MIRRORS = [
    "https://overpass-api.de/api/interpreter",
    "https://overpass.kumi.systems/api/interpreter",
    "https://maps.mail.ru/osm/tools/overpass/api/interpreter",
]
LISTINGS_FILE = "bird_listings.json"
OUTPUT_JS = "map_data.js"

HEADERS = {
    "User-Agent": "birdrent-map-builder/1.0 (educational project; contact eleazar.levchenko@gmail.com)",
    "Accept": "application/json",
}

# Kyiv bounding box  (south, west, north, east)
BBOX = (50.20, 30.15, 50.65, 30.95)

# Colour map for Kyiv metro lines (by OSM ref tag)
METRO_LINE_COLORS = {
    "1": "#d52b1e",  # M1 Sviatoshynsko-Brovarska (red)
    "2": "#0066b3",  # M2 Obolonsko-Teremkivska   (blue)
    "3": "#00a651",  # M3 Syretsko-Pecherska       (green)
}


# ---------------------------------------------------------------------------
# Overpass helpers
# ---------------------------------------------------------------------------

def overpass_query(query: str, retries: int = 3) -> dict:
    """POST query to Overpass API, rotating mirrors on failure."""
    last_exc = None
    for mirror in OVERPASS_MIRRORS:
        for attempt in range(retries):
            try:
                resp = requests.post(
                    mirror,
                    data={"data": query},
                    headers=HEADERS,
                    timeout=180,
                )
                resp.raise_for_status()
                return resp.json()
            except Exception as exc:
                last_exc = exc
                if attempt < retries - 1:
                    wait = 8 * (attempt + 1)
                    print(f"    [retry {attempt+1}/{retries-1} on {mirror} in {wait}s: {exc}]")
                    time.sleep(wait)
                else:
                    print(f"    [failed on {mirror}: {exc}, trying next mirror]")
    raise RuntimeError(f"All Overpass mirrors failed. Last error: {last_exc}")


# ---------------------------------------------------------------------------
# Metro
# ---------------------------------------------------------------------------

def _dist2d(a, b) -> float:
    return math.sqrt((a[0] - b[0]) ** 2 + (a[1] - b[1]) ** 2)


def _stitch_segments(segments: list) -> list:
    """Greedily stitch a list of [lon,lat] coordinate arrays into continuous lines."""
    if not segments:
        return []

    THRESHOLD = 0.0015  # ~150 m tolerance for segment join
    result = []
    current = list(segments[0])
    remaining = list(segments[1:])

    while remaining:
        best_i, best_rev, best_prepend = None, False, False
        best_d = float("inf")

        for i, seg in enumerate(remaining):
            for rev in (False, True):
                s = seg[::-1] if rev else seg
                for prepend in (False, True):
                    anchor = current[0] if prepend else current[-1]
                    d = _dist2d(anchor, s[0] if not prepend else s[-1])
                    if d < best_d:
                        best_d, best_i, best_rev, best_prepend = d, i, rev, prepend

        seg = remaining.pop(best_i)
        if best_rev:
            seg = seg[::-1]

        if best_d > THRESHOLD:
            # Gap too large — save current line and start fresh
            result.append(current)
            current = list(seg)
        elif best_prepend:
            current = list(seg[:-1]) + current
        else:
            current = current + list(seg[1:])

    result.append(current)
    return result


def _point_to_line_dist(px: float, py: float, line_coords: list) -> float:
    """Minimum distance from point to a polyline (in degrees, for comparison only)."""
    min_d = float("inf")
    for i in range(len(line_coords) - 1):
        ax, ay = line_coords[i]
        bx, by = line_coords[i + 1]
        # Project point onto segment
        dx, dy = bx - ax, by - ay
        t = max(0.0, min(1.0, ((px - ax) * dx + (py - ay) * dy) / (dx*dx + dy*dy + 1e-18)))
        cx, cy = ax + t * dx, ay + t * dy
        d = _dist2d([px, py], [cx, cy])
        if d < min_d:
            min_d = d
    return min_d


def fetch_metro() -> tuple:
    """Return (lines_fc, stations_fc) as GeoJSON FeatureCollections."""
    print("Fetching Kyiv metro from Overpass …")
    s, w, n, e = BBOX

    # Two separate queries:
    #   1. Route relations  → line geometry (ways carry geometry via out geom)
    #   2. Station nodes    → named stops with tags (out body gives tags, out geom gives coords)
    lines_query = f"""
[out:json][timeout:120];
relation["route"="subway"]({s},{w},{n},{e});
out geom;
"""
    stations_query = f"""
[out:json][timeout:120];
(
  node["railway"="station"]["station"="subway"]({s},{w},{n},{e});
  node["railway"="stop_position"]["subway"="yes"]({s},{w},{n},{e});
);
out body;
out skel qt;
"""
    lines_data    = overpass_query(lines_query)
    stations_data = overpass_query(stations_query)

    # ---- Parse route relations → line features + per-line colour map ----
    line_features = []
    line_meta = []  # list of (line_name, colour, stitched_coords_list)

    for rel in lines_data["elements"]:
        if rel["type"] != "relation":
            continue
        tags = rel.get("tags", {})
        ref  = tags.get("ref", "")
        colour = tags.get("colour", METRO_LINE_COLORS.get(ref, "#888888"))
        if colour and not colour.startswith("#"):
            colour = f"#{colour}"
        line_name = f"M{ref}" if ref else tags.get("name:en", tags.get("name", "?"))

        segments = []
        for member in rel.get("members", []):
            if member["type"] == "way" and "geometry" in member:
                coords = [[round(nd["lon"], 7), round(nd["lat"], 7)]
                          for nd in member["geometry"]]
                if len(coords) >= 2:
                    segments.append(coords)

        stitched_segs = _stitch_segments(segments)
        for seg in stitched_segs:
            if len(seg) >= 2:
                line_features.append({
                    "type": "Feature",
                    "properties": {"line": line_name, "color": colour},
                    "geometry": {"type": "LineString", "coordinates": seg},
                })
        # Flatten all coords for this line (used for nearest-line assignment)
        all_coords = [c for seg in stitched_segs for c in seg]
        if all_coords:
            line_meta.append((line_name, colour, all_coords))

    # ---- Parse station nodes and assign each to nearest metro line ----
    station_features = []
    seen_stations: set = set()

    for el in stations_data["elements"]:
        if el["type"] != "node":
            continue
        lat = el.get("lat")
        lon = el.get("lon")
        if lat is None or lon is None:
            continue
        node_tags = el.get("tags", {})
        name = node_tags.get("name:en") or node_tags.get("name", "")
        if not name:
            continue

        # Find the closest metro line
        best_line, best_colour, best_d = "?", "#888888", float("inf")
        for lname, lcolour, lcoords in line_meta:
            d = _point_to_line_dist(lon, lat, lcoords)
            if d < best_d:
                best_d, best_line, best_colour = d, lname, lcolour

        key = (name, best_line)
        if key in seen_stations:
            continue
        seen_stations.add(key)

        station_features.append({
            "type": "Feature",
            "properties": {"name": name, "line": best_line, "color": best_colour},
            "geometry": {"type": "Point", "coordinates": [round(lon, 7), round(lat, 7)]},
        })

    print(f"  {len(line_features)} line segments, {len(station_features)} station stops")
    return (
        {"type": "FeatureCollection", "features": line_features},
        {"type": "FeatureCollection", "features": station_features},
    )


# ---------------------------------------------------------------------------
# Building footprints
# ---------------------------------------------------------------------------

def _point_in_polygon(px: float, py: float, ring: list) -> bool:
    """Ray-casting point-in-polygon test."""
    inside = False
    n = len(ring)
    j = n - 1
    for i in range(n):
        xi, yi = ring[i]
        xj, yj = ring[j]
        if ((yi > py) != (yj > py)) and (
            px < (xj - xi) * (py - yi) / (yj - yi) + xi
        ):
            inside = not inside
        j = i
    return inside


def _polygon_area(ring: list) -> float:
    """Shoelace formula — returns area in degrees² (used only for relative comparison)."""
    n = len(ring)
    a = 0.0
    for i in range(n):
        j = (i + 1) % n
        a += ring[i][0] * ring[j][1] - ring[j][0] * ring[i][1]
    return abs(a) / 2.0


def _centroid(ring: list) -> tuple:
    return (
        sum(c[0] for c in ring) / len(ring),
        sum(c[1] for c in ring) / len(ring),
    )


def fetch_building_footprints(all_listings: list) -> dict:
    """
    For each unique building (by building_id), query Overpass for OSM building
    polygons within 70 m, then match each listing coordinate to the smallest
    containing polygon (or nearest centroid as fallback).

    Returns a GeoJSON FeatureCollection where each feature's properties include
    the aggregated listing data for that building.
    """
    # Aggregate listings by building_id
    by_building: dict = defaultdict(list)
    for lst in all_listings:
        by_building[lst["building_id"]].append(lst)

    unique_buildings = [v[0] for v in by_building.values()]  # one representative per building
    total = len(unique_buildings)
    print(f"Fetching building footprints for {total} unique buildings …")

    BATCH_SIZE = 50
    batches = [unique_buildings[i: i + BATCH_SIZE]
               for i in range(0, total, BATCH_SIZE)]

    all_osm_ways: dict = {}  # osm_id -> element

    for idx, batch in enumerate(batches):
        print(f"  batch {idx+1}/{len(batches)} … ", end="", flush=True)
        union = "".join(
            f'way["building"](around:70,{b["latitude"]},{b["longitude"]});'
            for b in batch
        )
        query = f"[out:json][timeout:120];({union});out geom;"
        result = overpass_query(query)

        added = 0
        for el in result["elements"]:
            if el["type"] == "way" and "geometry" in el and el["id"] not in all_osm_ways:
                all_osm_ways[el["id"]] = el
                added += 1

        print(f"+{added} (total {len(all_osm_ways)})")
        if idx < len(batches) - 1:
            time.sleep(3)

    print(f"  Downloaded {len(all_osm_ways)} unique OSM building polygons")
    print("  Matching listings to polygons …")

    # Pre-parse OSM way geometries
    parsed_ways = {}  # osm_id -> {"ring": [[lon,lat],...], "tags": {...}}
    for osm_id, el in all_osm_ways.items():
        ring = [[round(n["lon"], 7), round(n["lat"], 7)] for n in el["geometry"]]
        if len(ring) >= 3:
            parsed_ways[osm_id] = {"ring": ring, "tags": el.get("tags", {})}

    # Match each building's coordinate to an OSM polygon
    building_to_osm: dict = {}   # building_id -> osm_id
    unmatched = 0

    for bld in unique_buildings:
        bid = bld["building_id"]
        px, py = bld["longitude"], bld["latitude"]

        containing = []
        for osm_id, pw in parsed_ways.items():
            if _point_in_polygon(px, py, pw["ring"]):
                containing.append((osm_id, _polygon_area(pw["ring"])))

        if containing:
            # Smallest containing polygon = most specific building
            containing.sort(key=lambda x: x[1])
            building_to_osm[bid] = containing[0][0]
        else:
            # Nearest centroid fallback (within 60 m ≈ 0.0005°)
            best_d, best_id = float("inf"), None
            for osm_id, pw in parsed_ways.items():
                cx, cy = _centroid(pw["ring"])
                d = _dist2d([px, py], [cx, cy])
                if d < best_d:
                    best_d, best_id = d, osm_id
            if best_id and best_d < 0.0006:
                building_to_osm[bid] = best_id
            else:
                unmatched += 1

    print(f"  Matched {len(building_to_osm)}/{total} buildings ({unmatched} unmatched)")

    # Build GeoJSON features — one per OSM building, aggregating all listings
    osm_to_listings: dict = defaultdict(list)
    for bid, osm_id in building_to_osm.items():
        osm_to_listings[osm_id].extend(by_building[bid])

    features = []
    for osm_id, listings in osm_to_listings.items():
        pw = parsed_ways[osm_id]
        ring = pw["ring"]
        if ring[0] != ring[-1]:
            ring = ring + [ring[0]]

        tags = pw["tags"]
        osm_levels = 0
        for key in ("building:levels", "levels"):
            try:
                osm_levels = int(float(tags.get(key, 0)))
                if osm_levels:
                    break
            except (ValueError, TypeError):
                pass

        max_floor = max((l.get("floor_count") or 0) for l in listings)
        height_m = (osm_levels * 3) if osm_levels else (max_floor * 3 if max_floor else 30)

        has_premium = any(l.get("tier") == "premium" for l in listings)
        first = listings[0]
        label = f"{first['street']} {first['number']}"
        # Compact listing summary for popup
        listing_summaries = [
            {
                "id": l["listing_id"],
                "price": l.get("price"),
                "currency": l.get("currency", "USD"),
                "rooms": l.get("rooms"),
                "area": l.get("area"),
                "floor": l.get("floor"),
                "tier": l.get("tier"),
                "url": l.get("source_url"),
            }
            for l in listings
        ]

        features.append({
            "type": "Feature",
            "properties": {
                "osm_id": osm_id,
                "label": label,
                "count": len(listings),
                "height": height_m,
                "is_premium": 1 if has_premium else 0,
                "listings": listing_summaries,
            },
            "geometry": {"type": "Polygon", "coordinates": [ring]},
        })

    return {"type": "FeatureCollection", "features": features}


def build_listing_points(building_footprints_fc: dict, all_listings: list) -> dict:
    """
    One Point feature per building (for low-zoom circle markers).
    Uses the centroid of the matched footprint polygon if available,
    otherwise the raw lat/lon from the listing.
    """
    by_building: dict = defaultdict(list)
    for lst in all_listings:
        by_building[lst["building_id"]].append(lst)

    # Map label → centroid from footprints
    label_to_centroid = {}
    for feat in building_footprints_fc["features"]:
        ring = feat["geometry"]["coordinates"][0]
        cx, cy = _centroid(ring)
        label_to_centroid[feat["properties"]["label"]] = [round(cx, 7), round(cy, 7)]

    features = []
    seen_labels: set = set()

    for feat in building_footprints_fc["features"]:
        props = feat["properties"]
        label = props["label"]
        if label in seen_labels:
            continue
        seen_labels.add(label)

        coord = label_to_centroid.get(label)
        if not coord:
            continue

        features.append({
            "type": "Feature",
            "properties": {
                "label": label,
                "count": props["count"],
                "is_premium": props["is_premium"],
            },
            "geometry": {"type": "Point", "coordinates": coord},
        })

    return {"type": "FeatureCollection", "features": features}


# ---------------------------------------------------------------------------
# Output
# ---------------------------------------------------------------------------

def write_js(metro_lines, metro_stations, building_footprints, listing_points):
    def compact(obj):
        return json.dumps(obj, ensure_ascii=False, separators=(",", ":"))

    js = f"""// Auto-generated by fetch_map_data.py — do not edit manually.
// Re-run fetch_map_data.py to refresh from OpenStreetMap / bird_listings.json.
/* eslint-disable */
const METRO_LINES = {compact(metro_lines)};
const METRO_STATIONS = {compact(metro_stations)};
const BUILDING_FOOTPRINTS = {compact(building_footprints)};
const LISTING_POINTS = {compact(listing_points)};
"""
    with open(OUTPUT_JS, "w", encoding="utf-8") as f:
        f.write(js)

    size_kb = len(js.encode()) / 1024
    print(f"\nWrote {OUTPUT_JS} ({size_kb:.0f} KB)")
    print(f"  Metro lines:        {len(metro_lines['features'])} features")
    print(f"  Metro stations:     {len(metro_stations['features'])} features")
    print(f"  Building footprints:{len(building_footprints['features'])} features")
    print(f"  Listing points:     {len(listing_points['features'])} features")


# ---------------------------------------------------------------------------
# Main
# ---------------------------------------------------------------------------

def main():
    print(f"Loading {LISTINGS_FILE} …")
    with open(LISTINGS_FILE, encoding="utf-8") as f:
        all_listings = json.load(f)
    by_building = defaultdict(list)
    for l in all_listings:
        by_building[l["building_id"]].append(l)
    print(f"  {len(all_listings)} listings across {len(by_building)} unique buildings")

    metro_lines, metro_stations = fetch_metro()
    building_footprints = fetch_building_footprints(all_listings)
    listing_points = build_listing_points(building_footprints, all_listings)

    write_js(metro_lines, metro_stations, building_footprints, listing_points)
    print("\nDone. Open bird_rent_3d_map_with_metro_and_csv.html in a browser.")


if __name__ == "__main__":
    main()