Files
settled-reach/tooling/planet-gen/generate_atlas.py
T
jpmschweitzer 5619218e31 feat(tooling): generate_atlas.py city + infrastructure pipeline and atlas DB index (#832)
Implements the Phase 3 atlas content generator per D-191 §3, §8, and §9.

Pipeline per body (terrain-aware, deterministic per seed + body):
  1. Simulate terrain via planet_simulation.simulate().
  2. Analyse continents (flood-fill), habitability (temp/moisture/slope +
     coastal bonus), river mouths, and a terrain A* cost grid.
  3. Place cities sequentially — capital first (habitability + river-mouth
     bias), then corridor growth via multi-source Dijkstra, quadrant-spread
     penalty after 2 cities in a quadrant, port-on-new-continent bonus at
     cities 3–4. ±25% noise for seed variation.
  4. Generate roads and railroads as an MST over city positions, with
     A* paths on the terrain cost grid (rail follows roads where possible).
  5. Place a transit POI at the capital (15% chance to scatter to a
     secondary city).

Output (canonical markers.json schema, pixel space per D-191 §8):
  - cities:    {id, name, kind, center:[r,c], population}
  - roads:     {id, name, kind, path:[[r,c],...]}
  - railroads: {id, name, kind, path:[[r,c],...]}
  - pois:      {id, name, kind, center:[r,c]}
  - existing rivers/oceans/mountain_ranges preserved untouched.

City names are left empty for gemma_naming.py (#833). Body population is
split across cities with geometric decay (capital ~50%, each subsequent
city half the previous). The 6 hand-authored bodies (Lendel, Edict,
Vuurkloof, Røros, Cairnside, Estrade) are detected by existing
`cities` and skipped for regeneration; their markers are still synced
to the DB index below.

Atlas index in systems.db (new):
  - atlas_body_grids, atlas_cities, atlas_roads, atlas_railroads,
    atlas_pois, atlas_rivers, atlas_oceans, atlas_mountain_ranges
  - Scalar metadata mirror of every markers.json — the implant atlas app
    and development queries can lookup cities/POIs/features without
    scanning 267 JSON files. Polyline geometry stays in the markers.json
    files next to the heightmaps (used by the renderer); the DB only
    stores filterable scalar fields plus `point_count` as a length proxy.
  - Schema lives in server/data/systems-schema.sql; generate_atlas.py
    mirrors the CREATE TABLE IF NOT EXISTS block so it runs against any
    DB state (matches the economy-db importer pattern).
  - Populated and refreshed on every run. Each body's rows are deleted
    and reinserted deterministically — no stale state.

Also fixes a pre-existing WIP bug in the quadrant-saturation penalty
loop (a stray outer `for r in range(GRID_H)` with unreachable breaks
meant only the NW quadrant was ever checked).

Runtime: 280s for all 267 inhabited bodies on a single core. 265 bodies
updated this run, 6 hand-authored bodies synced to DB without
regeneration.

Atlas index after run:
  atlas_cities             329    (15 hand-authored + 314 awaiting #833)
  atlas_roads               46
  atlas_railroads           44
  atlas_pois               287
  atlas_rivers            2034
  atlas_oceans             696
  atlas_mountain_ranges   1953
  atlas_body_grids         267
2026-04-15 08:46:55 +02:00

1334 lines
46 KiB
Python
Raw Blame History

This file contains ambiguous Unicode characters
This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.
#!/usr/bin/env python3
"""
generate_atlas.py — Terrain-aware sequential city placement and infrastructure
generation for The Settled Reach Atlas (Phase 3, D-191 §3, §8, §9, #832).
Pipeline per body:
1. Load body definition from body index.md frontmatter
2. Simulate terrain via planet_simulation.simulate()
3. Analyse terrain: continents, habitability, river mouths, cost grid
4. Place cities sequentially (capital first, corridor growth, quadrant spread)
5. Generate infrastructure: A* roads + rail MST connecting all cities
6. Place gate terminal POI at largest population centre (occasionally scatter)
7. Write updated markers.json (preserves existing rivers/oceans/mountains)
City names are left as empty strings; gemma_naming.py (#833) fills them.
Usage:
python3 tooling/planet-gen/generate_atlas.py
python3 tooling/planet-gen/generate_atlas.py --body GJ380c
python3 tooling/planet-gen/generate_atlas.py --force
python3 tooling/planet-gen/generate_atlas.py --dry-run
python3 tooling/planet-gen/generate_atlas.py --seed 12345
Decisions: D-191 (atlas pipeline), D-188 (planet_class field name)
"""
import argparse
import hashlib
import heapq
import json
import math
import os
import sys
import time
from pathlib import Path
# ---------------------------------------------------------------------------
# Venv bootstrap
# ---------------------------------------------------------------------------
TOOLING_DIR = Path(__file__).resolve().parent
REPO_ROOT = (TOOLING_DIR / ".." / "..").resolve()
_venv_python = REPO_ROOT / ".venv" / "bin" / "python"
if _venv_python.exists() and Path(sys.executable).resolve() != _venv_python.resolve():
os.execv(str(_venv_python), [str(_venv_python)] + sys.argv)
import numpy as np
import sqlite3
import yaml
from planet_simulation import simulate
# ---------------------------------------------------------------------------
# Constants
# ---------------------------------------------------------------------------
GRID_W = 512
GRID_H = 256
# Terrain cost grid (D-191 §3: A* on terrain cost grid)
COST_WATER = 1e9 # impassable
COST_MOUNTAIN = 10.0 # expensive
COST_RIVER = 0.5 # cheap corridor
COST_FLAT = 1.0 # baseline land
# City placement
MOUNTAIN_ELEV_THRESHOLD = 0.55 # normalised elevation above sea_level → mountain
SETTLEMENT_PATTERN_MODIFIERS = {
"urban_concentrated": -1,
"dispersed_rural": +1,
"orbital_only": None, # no surface cities
"domed": None, # special: forced to 1
"cave": None, # special: forced to 1
}
# Quadrant distribution (D-191 §3: after 2 cities in same quadrant, prefer others)
QUADRANT_SATURATION = 2
DB_PATH = REPO_ROOT / "server" / "data" / "systems.db"
WIKI_SYSTEMS = REPO_ROOT / "wiki" / "star-systems"
# ---------------------------------------------------------------------------
# Atlas index tables (D-191 §8, #832) — scalar mirror of markers.json.
# Duplicated from server/data/systems-schema.sql so generate_atlas.py works
# against any DB state without needing a wiki_sync first. Keep in sync with
# the canonical schema when adding columns.
# ---------------------------------------------------------------------------
ATLAS_MIGRATION_SQL = """
CREATE TABLE IF NOT EXISTS atlas_body_grids (
body_id TEXT PRIMARY KEY REFERENCES bodies(body_id),
grid_w INTEGER NOT NULL,
grid_h INTEGER NOT NULL,
updated_at TEXT DEFAULT (datetime('now'))
);
CREATE TABLE IF NOT EXISTS atlas_cities (
city_id TEXT PRIMARY KEY,
body_id TEXT NOT NULL REFERENCES bodies(body_id),
local_id TEXT NOT NULL,
name TEXT NOT NULL DEFAULT '',
kind TEXT NOT NULL,
center_row INTEGER NOT NULL,
center_col INTEGER NOT NULL,
population INTEGER NOT NULL DEFAULT 0,
updated_at TEXT DEFAULT (datetime('now'))
);
CREATE TABLE IF NOT EXISTS atlas_roads (
road_id TEXT PRIMARY KEY,
body_id TEXT NOT NULL REFERENCES bodies(body_id),
local_id TEXT NOT NULL,
name TEXT NOT NULL DEFAULT '',
kind TEXT NOT NULL,
point_count INTEGER NOT NULL,
updated_at TEXT DEFAULT (datetime('now'))
);
CREATE TABLE IF NOT EXISTS atlas_railroads (
railroad_id TEXT PRIMARY KEY,
body_id TEXT NOT NULL REFERENCES bodies(body_id),
local_id TEXT NOT NULL,
name TEXT NOT NULL DEFAULT '',
kind TEXT NOT NULL,
point_count INTEGER NOT NULL,
updated_at TEXT DEFAULT (datetime('now'))
);
CREATE TABLE IF NOT EXISTS atlas_pois (
poi_id TEXT PRIMARY KEY,
body_id TEXT NOT NULL REFERENCES bodies(body_id),
local_id TEXT NOT NULL,
name TEXT NOT NULL DEFAULT '',
kind TEXT NOT NULL,
center_row INTEGER NOT NULL,
center_col INTEGER NOT NULL,
updated_at TEXT DEFAULT (datetime('now'))
);
CREATE TABLE IF NOT EXISTS atlas_rivers (
river_id TEXT PRIMARY KEY,
body_id TEXT NOT NULL REFERENCES bodies(body_id),
local_id TEXT NOT NULL,
name TEXT NOT NULL DEFAULT '',
point_count INTEGER NOT NULL,
updated_at TEXT DEFAULT (datetime('now'))
);
CREATE TABLE IF NOT EXISTS atlas_oceans (
water_id TEXT PRIMARY KEY,
body_id TEXT NOT NULL REFERENCES bodies(body_id),
local_id TEXT NOT NULL,
name TEXT NOT NULL DEFAULT '',
kind TEXT NOT NULL,
center_row INTEGER NOT NULL,
center_col INTEGER NOT NULL,
area_fraction REAL NOT NULL DEFAULT 0.0,
updated_at TEXT DEFAULT (datetime('now'))
);
CREATE TABLE IF NOT EXISTS atlas_mountain_ranges (
range_id TEXT PRIMARY KEY,
body_id TEXT NOT NULL REFERENCES bodies(body_id),
local_id TEXT NOT NULL,
name TEXT NOT NULL DEFAULT '',
center_row INTEGER NOT NULL,
center_col INTEGER NOT NULL,
peak_row INTEGER NOT NULL,
peak_col INTEGER NOT NULL,
area_cells INTEGER NOT NULL DEFAULT 0,
updated_at TEXT DEFAULT (datetime('now'))
);
CREATE INDEX IF NOT EXISTS idx_atlas_cities_body ON atlas_cities(body_id);
CREATE INDEX IF NOT EXISTS idx_atlas_cities_kind ON atlas_cities(kind);
CREATE INDEX IF NOT EXISTS idx_atlas_cities_population ON atlas_cities(population);
CREATE INDEX IF NOT EXISTS idx_atlas_cities_name ON atlas_cities(name);
CREATE INDEX IF NOT EXISTS idx_atlas_roads_body ON atlas_roads(body_id);
CREATE INDEX IF NOT EXISTS idx_atlas_railroads_body ON atlas_railroads(body_id);
CREATE INDEX IF NOT EXISTS idx_atlas_pois_body ON atlas_pois(body_id);
CREATE INDEX IF NOT EXISTS idx_atlas_pois_kind ON atlas_pois(kind);
CREATE INDEX IF NOT EXISTS idx_atlas_rivers_body ON atlas_rivers(body_id);
CREATE INDEX IF NOT EXISTS idx_atlas_oceans_body ON atlas_oceans(body_id);
CREATE INDEX IF NOT EXISTS idx_atlas_mountain_ranges_body ON atlas_mountain_ranges(body_id);
"""
def ensure_atlas_schema(conn: sqlite3.Connection) -> None:
"""Apply CREATE TABLE / INDEX IF NOT EXISTS for all atlas_* tables."""
conn.executescript(ATLAS_MIGRATION_SQL)
def _first_int(values, default: int = 0) -> int:
"""Coerce the first numeric element of an iterable (e.g. [r, c]) to int."""
try:
return int(values[0])
except (TypeError, ValueError, IndexError):
return default
def sync_markers_to_db(
conn: sqlite3.Connection,
body_id: str,
markers: dict,
) -> dict:
"""Refresh atlas_* rows for a single body from its markers.json contents.
Deletes any existing rows for the body (cascading across all seven atlas
tables) and re-inserts from the markers dict. Returns a count dict for
reporting.
"""
# Wipe existing rows for this body — fully deterministic rebuild.
for table in (
"atlas_cities",
"atlas_roads",
"atlas_railroads",
"atlas_pois",
"atlas_rivers",
"atlas_oceans",
"atlas_mountain_ranges",
):
conn.execute(f"DELETE FROM {table} WHERE body_id = ?", (body_id,))
grid = markers.get("grid") or {}
grid_w = int(grid.get("w", GRID_W))
grid_h = int(grid.get("h", GRID_H))
conn.execute(
"""INSERT INTO atlas_body_grids (body_id, grid_w, grid_h, updated_at)
VALUES (?, ?, ?, datetime('now'))
ON CONFLICT(body_id) DO UPDATE SET
grid_w = excluded.grid_w,
grid_h = excluded.grid_h,
updated_at = excluded.updated_at""",
(body_id, grid_w, grid_h),
)
counts = {
"cities": 0, "roads": 0, "railroads": 0, "pois": 0,
"rivers": 0, "oceans": 0, "mountain_ranges": 0,
}
for city in markers.get("cities") or []:
local = city.get("id")
if not local:
continue
center = city.get("center") or [0, 0]
conn.execute(
"""INSERT INTO atlas_cities
(city_id, body_id, local_id, name, kind, center_row, center_col, population)
VALUES (?, ?, ?, ?, ?, ?, ?, ?)""",
(
f"{body_id}/{local}",
body_id,
local,
city.get("name") or "",
city.get("kind") or "city",
_first_int(center, 0),
_first_int(center[1:] if len(center) > 1 else [0], 0),
int(city.get("population") or 0),
),
)
counts["cities"] += 1
for road in markers.get("roads") or []:
local = road.get("id")
if not local:
continue
path = road.get("path") or []
conn.execute(
"""INSERT INTO atlas_roads
(road_id, body_id, local_id, name, kind, point_count)
VALUES (?, ?, ?, ?, ?, ?)""",
(
f"{body_id}/{local}",
body_id,
local,
road.get("name") or "",
road.get("kind") or "commercial",
len(path),
),
)
counts["roads"] += 1
for rail in markers.get("railroads") or []:
local = rail.get("id")
if not local:
continue
path = rail.get("path") or []
conn.execute(
"""INSERT INTO atlas_railroads
(railroad_id, body_id, local_id, name, kind, point_count)
VALUES (?, ?, ?, ?, ?, ?)""",
(
f"{body_id}/{local}",
body_id,
local,
rail.get("name") or "",
rail.get("kind") or "passenger_freight",
len(path),
),
)
counts["railroads"] += 1
for poi in markers.get("pois") or []:
local = poi.get("id")
if not local:
continue
center = poi.get("center") or [0, 0]
conn.execute(
"""INSERT INTO atlas_pois
(poi_id, body_id, local_id, name, kind, center_row, center_col)
VALUES (?, ?, ?, ?, ?, ?, ?)""",
(
f"{body_id}/{local}",
body_id,
local,
poi.get("name") or "",
poi.get("kind") or "transit",
_first_int(center, 0),
_first_int(center[1:] if len(center) > 1 else [0], 0),
),
)
counts["pois"] += 1
for river in markers.get("rivers") or []:
local = river.get("id")
if not local:
continue
path = river.get("path") or []
conn.execute(
"""INSERT INTO atlas_rivers
(river_id, body_id, local_id, name, point_count)
VALUES (?, ?, ?, ?, ?)""",
(
f"{body_id}/{local}",
body_id,
local,
river.get("name") or "",
len(path),
),
)
counts["rivers"] += 1
for water in markers.get("oceans") or []:
local = water.get("id")
if not local:
continue
center = water.get("center") or [0, 0]
conn.execute(
"""INSERT INTO atlas_oceans
(water_id, body_id, local_id, name, kind, center_row, center_col, area_fraction)
VALUES (?, ?, ?, ?, ?, ?, ?, ?)""",
(
f"{body_id}/{local}",
body_id,
local,
water.get("name") or "",
water.get("kind") or "ocean",
_first_int(center, 0),
_first_int(center[1:] if len(center) > 1 else [0], 0),
float(water.get("area_fraction") or 0.0),
),
)
counts["oceans"] += 1
for rng_feat in markers.get("mountain_ranges") or []:
local = rng_feat.get("id")
if not local:
continue
center = rng_feat.get("center") or [0, 0]
peak = rng_feat.get("peak") or center
conn.execute(
"""INSERT INTO atlas_mountain_ranges
(range_id, body_id, local_id, name, center_row, center_col,
peak_row, peak_col, area_cells)
VALUES (?, ?, ?, ?, ?, ?, ?, ?, ?)""",
(
f"{body_id}/{local}",
body_id,
local,
rng_feat.get("name") or "",
_first_int(center, 0),
_first_int(center[1:] if len(center) > 1 else [0], 0),
_first_int(peak, 0),
_first_int(peak[1:] if len(peak) > 1 else [0], 0),
int(rng_feat.get("area_cells") or 0),
),
)
counts["mountain_ranges"] += 1
return counts
# ---------------------------------------------------------------------------
# Grid helpers
# ---------------------------------------------------------------------------
# markers.json is stored in pixel space per D-191 §8 — the canonical format is
# `center: [row, col]` and `path: [[row, col], ...]`. Lat/lon is a display-time
# derivation in the atlas UI, not a storage format.
def quadrant(row: int, col: int) -> int:
"""Return quadrant index 0-3: NW=0, NE=1, SW=2, SE=3."""
h = 0 if row < GRID_H // 2 else 1
w = 0 if col < GRID_W // 2 else 1
return h * 2 + w
def body_id_salt(body_id: str) -> int:
"""Deterministic per-body salt from body_id."""
return int(hashlib.sha256(body_id.encode()).hexdigest()[:8], 16)
# ---------------------------------------------------------------------------
# Body definition loading
# ---------------------------------------------------------------------------
def load_body_def(body_dir: Path) -> dict | None:
"""Load body definition from body index.md frontmatter."""
index_md = body_dir / "index.md"
if not index_md.exists():
return None
content = index_md.read_text()
if not content.startswith("---"):
return None
try:
end = content.index("---", 3)
bd = yaml.safe_load(content[3:end])
except (ValueError, yaml.YAMLError):
return None
if not bd or "id" not in bd or "planet_class" not in bd:
return None
return bd
# ---------------------------------------------------------------------------
# Terrain analysis
# ---------------------------------------------------------------------------
def _analyse_terrain(terrain: dict, seed_rng: np.random.Generator) -> dict:
"""Analyse terrain for atlas city placement.
Returns:
continents: list of continent dicts {cells, area, id}
habitability: float32 (H, W) — 0=uninhabitable, 1=ideal
river_mouths: list of (row, col) — land cells where rivers meet water
cost_grid: float32 (H, W) — A* traversal cost per cell
continent_map: int32 (H, W) — continent label per land cell (-1=water)
land_mask: bool (H, W)
"""
from scipy.ndimage import label
elevation = terrain["elevation"]
temperature = terrain["temperature"]
moisture = terrain["moisture"]
surface_water = terrain["surface_water"]
river_grid = terrain.get("river_grid", np.zeros((GRID_H, GRID_W), bool))
sea_level = terrain["sea_level"]
land = ~surface_water
# Normalised elevation above sea level
elev_norm = np.where(land,
(elevation - sea_level) / (1.0 - sea_level + 1e-9),
0.0)
mountains = land & (elev_norm > MOUNTAIN_ELEV_THRESHOLD)
# --- Continent detection (flood-fill connected land) ---
continent_labels, n_continents = label(land)
continents = []
for lbl in range(1, n_continents + 1):
mask = continent_labels == lbl
area = int(mask.sum())
if area < 50:
continue # ignore tiny rock outcroppings
ys, xs = np.where(mask)
continents.append({
"id": lbl,
"area": area,
"cells": list(zip(ys.tolist(), xs.tolist())),
"_center": (int(ys.mean()), int(xs.mean())),
})
continents.sort(key=lambda c: c["area"], reverse=True)
continent_map = np.where(land, continent_labels, -1).astype(np.int32)
# --- Habitability scoring ---
# temperature: 0.3–0.7 is ideal (mid-range)
temp_score = 1.0 - np.abs(temperature - 0.5) * 2.0
temp_score = np.clip(temp_score, 0.0, 1.0)
# moisture: 0.3–0.6 is ideal
moist_score = np.where(
moisture < 0.3, moisture / 0.3,
np.where(moisture < 0.6, 1.0, 1.0 - (moisture - 0.6) / 0.4)
)
moist_score = np.clip(moist_score, 0.0, 1.0)
# slope: low slope = ideal (flat land)
dy = np.gradient(elevation, axis=0)
dx = np.gradient(elevation, axis=1)
slope = np.sqrt(dy**2 + dx**2)
slope_max = float(slope[land].max()) if land.any() else 1.0
slope_score = 1.0 - np.clip(slope / (slope_max + 1e-9), 0.0, 1.0)
habitability = (temp_score * 0.4 + moist_score * 0.3 + slope_score * 0.3)
habitability = np.where(land, habitability, 0.0).astype(np.float32)
# Coastal access bonus: land cells adjacent to water
water_f = surface_water.astype(np.float32)
from scipy.ndimage import uniform_filter
water_proximity = uniform_filter(water_f, size=7)
coastal_bonus = np.clip(water_proximity * 2.0, 0.0, 0.3)
habitability = np.where(land, habitability + coastal_bonus.astype(np.float32), 0.0)
habitability = np.clip(habitability, 0.0, 1.0).astype(np.float32)
# --- River mouth detection ---
# A river mouth is a river-grid cell that is land but adjacent to water
river_mouths = []
if river_grid.any() and surface_water.any():
# Dilate surface_water by 1 cell
from scipy.ndimage import binary_dilation
water_dilated = binary_dilation(surface_water)
# River cells on land that touch water
mouth_mask = river_grid & land & water_dilated
rys, rxs = np.where(mouth_mask)
# Deduplicate by proximity (cluster to strongest mouth)
if len(rys) > 0:
used = np.zeros(len(rys), bool)
for i in range(len(rys)):
if used[i]:
continue
used[i] = True
r, c = int(rys[i]), int(rxs[i])
# Suppress nearby duplicates within radius 10
dist = np.sqrt((rys - r)**2 + (rxs - c)**2)
used[dist < 10] = True
river_mouths.append((r, c))
# --- Terrain cost grid (D-191 §3) ---
cost_grid = np.full((GRID_H, GRID_W), COST_FLAT, dtype=np.float32)
cost_grid[surface_water] = COST_WATER
cost_grid[mountains] = COST_MOUNTAIN
cost_grid[river_grid & land] = COST_RIVER
return {
"continents": continents,
"habitability": habitability,
"river_mouths": river_mouths,
"cost_grid": cost_grid,
"continent_map": continent_map,
"land_mask": land,
}
# ---------------------------------------------------------------------------
# City count computation (D-191 §8)
# ---------------------------------------------------------------------------
def distribute_population(total_population: int, n_cities: int) -> list[int]:
"""Split a body's population across n_cities using a simple geometric decay.
Capital gets ~50%, each subsequent city gets half of what the previous one
did, and the remaining tail is absorbed by the capital so totals still match.
Returns a list of int populations (rounded to the nearest 10,000 for realism).
"""
if n_cities <= 0 or total_population <= 0:
return []
if n_cities == 1:
return [int(total_population)]
weights = [0.5 ** i for i in range(n_cities)]
weight_sum = sum(weights)
shares = [w / weight_sum for w in weights]
populations = []
remaining = total_population
for i, share in enumerate(shares):
if i == len(shares) - 1:
pop = remaining
else:
pop = int(round(total_population * share / 10_000) * 10_000)
pop = min(pop, remaining - (len(shares) - i - 1) * 10_000)
pop = max(pop, 10_000)
populations.append(pop)
remaining -= pop
return populations
def compute_city_count(population: int, settlement_pattern: str | None) -> int:
"""Compute number of cities for a body.
Base count: floor(log10(pop / 1_000_000)), minimum 1 if inhabited.
Settlement pattern modifiers from D-191 §8.
"""
if population <= 0:
return 0
if settlement_pattern in ("orbital_only",):
return 0
if settlement_pattern in ("domed", "cave"):
return 1
base = int(math.floor(math.log10(max(population, 1_000_001) / 1_000_000)))
base = max(1, base) # at least 1 city if inhabited
mod = SETTLEMENT_PATTERN_MODIFIERS.get(settlement_pattern, 0) or 0
result = max(1, base + mod)
return min(result, 12) # cap at 12 to avoid overcrowding small maps
# ---------------------------------------------------------------------------
# City placement (D-191 §3)
# ---------------------------------------------------------------------------
def _score_capital_sites(
analysis: dict,
rng: np.random.Generator,
noise_factor: float = 0.25,
) -> np.ndarray:
"""Score all land cells for capital placement.
Capital scoring factors:
- Habitability score (temperature, moisture, slope)
- River mouth bonus (~50% of capitals at river mouths per D-191 §3)
- ±25% noise variation per seed
"""
habitability = analysis["habitability"]
land = analysis["land_mask"]
river_mouths = analysis["river_mouths"]
# Base score
score = habitability.copy()
# River mouth bonus (makes ~50% of capital choices go there)
river_bonus = np.zeros((GRID_H, GRID_W), dtype=np.float32)
for r, c in river_mouths:
# Smooth bonus around the river mouth
from scipy.ndimage import gaussian_filter
tmp = np.zeros((GRID_H, GRID_W), np.float32)
tmp[r, c] = 1.0
river_bonus += gaussian_filter(tmp, sigma=8.0)
if river_mouths:
river_bonus = river_bonus / (river_bonus.max() + 1e-9)
score = score * 0.5 + river_bonus * 0.5
# ±25% noise for seed variation (D-191 §3)
noise = rng.uniform(1.0 - noise_factor, 1.0 + noise_factor,
size=(GRID_H, GRID_W)).astype(np.float32)
score = score * noise
score = np.where(land, score, 0.0)
return score
def _quadrant_counts(placed: list[tuple[int, int]]) -> np.ndarray:
"""Return count of placed cities per quadrant (4 quadrants)."""
counts = np.zeros(4, dtype=int)
for r, c in placed:
counts[quadrant(r, c)] += 1
return counts
def _dijkstra_from_sources(
cost_grid: np.ndarray,
sources: list[tuple[int, int]],
) -> np.ndarray:
"""Multi-source Dijkstra on cost_grid from all source cells.
Returns distance array (H, W). Uses 8-connectivity.
"""
dist = np.full((GRID_H, GRID_W), np.inf, dtype=np.float64)
heap = []
for r, c in sources:
if cost_grid[r, c] < 1e8:
dist[r, c] = 0.0
heapq.heappush(heap, (0.0, r, c))
while heap:
d, r, c = heapq.heappop(heap)
if d > dist[r, c]:
continue
for dr in (-1, 0, 1):
for dc in (-1, 0, 1):
if dr == 0 and dc == 0:
continue
nr, nc = r + dr, c + dc
if nr < 0 or nr >= GRID_H or nc < 0 or nc >= GRID_W:
continue
edge = cost_grid[nr, nc]
if edge >= 1e8:
continue
step = math.sqrt(2) if dr != 0 and dc != 0 else 1.0
nd = d + edge * step
if nd < dist[nr, nc]:
dist[nr, nc] = nd
heapq.heappush(heap, (nd, nr, nc))
return dist
def place_cities(
analysis: dict,
n_cities: int,
rng: np.random.Generator,
noise_factor: float = 0.25,
body_population: int = 0,
) -> list[dict]:
"""Place cities sequentially per D-191 §3.
Returns a list of city dicts in the canonical markers.json schema
(D-191 §8 — pixel space):
{id, name, kind, center: [row, col], population}
with internal `_row`, `_col`, `_continent_id` fields used by downstream
infrastructure generation. The underscore-prefixed fields are stripped
before writing to disk.
"""
if n_cities <= 0:
return []
continents = analysis["continents"]
habitability = analysis["habitability"]
cost_grid = analysis["cost_grid"]
continent_map = analysis["continent_map"]
land = analysis["land_mask"]
if not continents:
return []
placed_coords: list[tuple[int, int]] = []
placed_continents: set[int] = set()
cities = []
# ── 1. Capital ──────────────────────────────────────────────────────────
capital_score = _score_capital_sites(analysis, rng, noise_factor)
# Suppress edge pixels (poles tend to be degenerate)
capital_score[:5, :] = 0
capital_score[-5:, :] = 0
best = int(np.argmax(capital_score))
cap_r, cap_c = divmod(best, GRID_W)
continent_id = int(continent_map[cap_r, cap_c]) if continent_map[cap_r, cap_c] > 0 else 1
placed_coords.append((cap_r, cap_c))
placed_continents.add(continent_id)
cities.append({
"id": "city_0",
"name": "",
"kind": "capital",
"center": [cap_r, cap_c],
"population": 0, # filled after all cities are placed
"_row": cap_r,
"_col": cap_c,
"_continent_id": continent_id,
})
if n_cities == 1:
_assign_populations(cities, body_population)
return cities
# Minimum exclusion radius: varies by city count (more cities → tighter spacing)
excl_radius = max(20, GRID_H // (n_cities + 1))
# ── 2. Subsequent cities (corridor growth + quadrant spread) ────────────
for i in range(1, n_cities):
q_counts = _quadrant_counts(placed_coords)
# Multi-source Dijkstra from all placed cities along cost grid
corridor_dist = _dijkstra_from_sources(cost_grid, placed_coords)
# Score for new city: habitability weighted by corridor distance
# Prefer connected but not too-close (distance 50–300 from existing cities)
ideal_min = excl_radius
ideal_max = min(300, excl_radius * 6)
dist_score = np.where(
(corridor_dist >= ideal_min) & (corridor_dist < ideal_max),
1.0 - (corridor_dist - ideal_min) / (ideal_max - ideal_min + 1e-9),
np.where(corridor_dist < ideal_min, 0.0, 0.1),
).astype(np.float32)
site_score = habitability * 0.6 + dist_score * 0.4
# Quadrant penalty: if a quadrant is already at QUADRANT_SATURATION,
# reduce its score so subsequent cities prefer emptier quadrants.
q_penalty = np.ones((GRID_H, GRID_W), np.float32)
for qr in (0, 1):
for qc in (0, 1):
q_idx = qr * 2 + qc
if q_counts[q_idx] >= QUADRANT_SATURATION:
rlo, rhi = (0, GRID_H // 2) if qr == 0 else (GRID_H // 2, GRID_H)
clo, chi = (0, GRID_W // 2) if qc == 0 else (GRID_W // 2, GRID_W)
q_penalty[rlo:rhi, clo:chi] = 0.3
site_score = site_score * q_penalty
# New-continent bonus at cities 3-4 (D-191 §3: port on new continent)
if i in (2, 3) and len(continents) > 1:
unexplored_conts = [
c for c in continents
if c["id"] not in placed_continents and c["area"] > 200
]
if unexplored_conts:
cont = unexplored_conts[0]
# Find coastline cells on that continent
coast_bonus = np.zeros((GRID_H, GRID_W), np.float32)
cont_mask = continent_map == cont["id"]
# Coastal cells: continent land cells adjacent to water
from scipy.ndimage import binary_dilation
water_dilated = binary_dilation(analysis["land_mask"] == False) # noqa: E712
coast_mask = cont_mask & water_dilated
coast_bonus[coast_mask] = 0.4
site_score = site_score + coast_bonus
# ±25% noise
noise = rng.uniform(1.0 - noise_factor, 1.0 + noise_factor,
size=(GRID_H, GRID_W)).astype(np.float32)
site_score = site_score * noise
site_score = np.where(land, site_score, 0.0)
# Suppress already-placed city zones
for pr, pc in placed_coords:
rlo = max(0, pr - excl_radius)
rhi = min(GRID_H, pr + excl_radius)
clo = max(0, pc - excl_radius)
chi = min(GRID_W, pc + excl_radius)
site_score[rlo:rhi, clo:chi] = 0.0
if site_score.max() < 1e-6:
break # no more habitable land
best = int(np.argmax(site_score))
r, c = divmod(best, GRID_W)
cont_id = int(continent_map[r, c]) if continent_map[r, c] > 0 else continent_id
placed_coords.append((r, c))
placed_continents.add(cont_id)
cities.append({
"id": f"city_{i}",
"name": "",
"kind": "city",
"center": [r, c],
"population": 0,
"_row": r,
"_col": c,
"_continent_id": cont_id,
})
_assign_populations(cities, body_population)
return cities
def _assign_populations(cities: list[dict], body_population: int) -> None:
"""Fill the `population` field on each city from the body's total."""
pops = distribute_population(body_population, len(cities))
for city, pop in zip(cities, pops):
city["population"] = int(pop)
# ---------------------------------------------------------------------------
# Infrastructure generation (D-191 §3)
# ---------------------------------------------------------------------------
def _astar_path(
cost_grid: np.ndarray,
start: tuple[int, int],
goal: tuple[int, int],
) -> list[tuple[int, int]]:
"""A* pathfinding on cost_grid (8-connectivity). Returns path or []."""
sr, sc = start
gr, gc = goal
def h(r, c):
return math.sqrt((r - gr)**2 + (c - gc)**2)
dist = {(sr, sc): 0.0}
prev = {}
heap = [(h(sr, sc), 0.0, sr, sc)]
while heap:
_, d, r, c = heapq.heappop(heap)
if r == gr and c == gc:
# Reconstruct path
path = []
cur = (gr, gc)
while cur in prev:
path.append(cur)
cur = prev[cur]
path.append((sr, sc))
path.reverse()
return path
if d > dist.get((r, c), float("inf")) + 1e-9:
continue
for dr in (-1, 0, 1):
for dc in (-1, 0, 1):
if dr == 0 and dc == 0:
continue
nr, nc = r + dr, c + dc
if nr < 0 or nr >= GRID_H or nc < 0 or nc >= GRID_W:
continue
edge = cost_grid[nr, nc]
if edge >= 1e8:
continue
step = math.sqrt(2) if dr != 0 and dc != 0 else 1.0
nd = d + edge * step
if nd < dist.get((nr, nc), float("inf")):
dist[(nr, nc)] = nd
prev[(nr, nc)] = (r, c)
heapq.heappush(heap, (nd + h(nr, nc), nd, nr, nc))
return []
def _mst_edges(n_cities: int, city_coords: list[tuple[int, int]]) -> list[tuple[int, int]]:
"""Prim's MST on Euclidean distances between cities. Returns list of (i, j) edges."""
if n_cities <= 1:
return []
in_tree = {0}
edges = []
while len(in_tree) < n_cities:
best_cost = float("inf")
best_edge = None
for i in in_tree:
for j in range(n_cities):
if j in in_tree:
continue
r1, c1 = city_coords[i]
r2, c2 = city_coords[j]
cost = math.sqrt((r1 - r2)**2 + (c1 - c2)**2)
if cost < best_cost:
best_cost = cost
best_edge = (i, j)
if best_edge is None:
break
edges.append(best_edge)
in_tree.add(best_edge[1])
return edges
def _subsample_path(path: list[tuple[int, int]], max_points: int = 64) -> list[list[int]]:
"""Subsample a path to at most max_points for compact JSON storage."""
if len(path) <= max_points:
return [[r, c] for r, c in path]
step = len(path) / max_points
result = []
for i in range(max_points):
idx = int(i * step)
r, c = path[idx]
result.append([r, c])
# Always include endpoint
r, c = path[-1]
if result[-1] != [r, c]:
result.append([r, c])
return result
def generate_infrastructure(
analysis: dict,
cities: list[dict],
) -> tuple[list[dict], list[dict]]:
"""Generate roads and railroads connecting cities.
Roads: A* paths on terrain cost grid connecting each MST edge.
Railroads: same MST edges but via a cost grid favouring corridors.
Returns (roads, railroads) — each is a list of path dicts.
"""
if len(cities) < 2:
return [], []
cost_grid = analysis["cost_grid"]
city_coords = [((c["_row"], c["_col"])) for c in cities]
n = len(city_coords)
mst = _mst_edges(n, city_coords)
roads = []
railroads = []
road_cost = cost_grid.copy()
# Rail cost: slightly prefer following roads (lower cost after first pass)
rail_cost = cost_grid.copy()
for edge_idx, (i, j) in enumerate(mst):
start = city_coords[i]
goal = city_coords[j]
road_path = _astar_path(road_cost, start, goal)
if road_path:
roads.append({
"id": f"road_{edge_idx}",
"name": "",
"kind": "commercial",
"path": _subsample_path(road_path),
})
# After placing a road, reduce cost along it for rail (rail follows roads)
for r, c in road_path:
if rail_cost[r, c] < 1e8:
rail_cost[r, c] = max(0.3, rail_cost[r, c] * 0.5)
rail_path = _astar_path(rail_cost, start, goal)
if rail_path:
railroads.append({
"id": f"railroad_{edge_idx}",
"name": "",
"kind": "passenger_freight",
"path": _subsample_path(rail_path),
})
return roads, railroads
# ---------------------------------------------------------------------------
# Gate terminal POI placement (D-191 §8)
# ---------------------------------------------------------------------------
def place_gate_terminal(
cities: list[dict],
rng: np.random.Generator,
) -> list[dict]:
"""Place a gate terminal POI at the largest population centre (sometimes scatter).
Primary placement: city 0 (capital, largest population center).
Occasional scatter: 15% chance to a secondary city (per D-191 §8).
Returns a list with a single POI dict in the canonical template schema:
`{id, name, kind: "transit", center: [row, col]}`. The name is left empty
for gemma_naming.py (#833) to fill.
"""
if not cities:
return []
scatter_prob = 0.15
gate_city = cities[0]
if len(cities) > 1 and rng.random() < scatter_prob:
gate_city = cities[int(rng.integers(1, len(cities)))]
return [{
"id": "poi_0",
"name": "",
"kind": "transit",
"center": [gate_city["_row"], gate_city["_col"]],
}]
# ---------------------------------------------------------------------------
# Markers.json update
# ---------------------------------------------------------------------------
def load_markers(body_dir: Path) -> dict:
"""Load existing markers.json, return empty structure if missing."""
markers_path = body_dir / "markers.json"
if markers_path.exists():
with open(markers_path) as f:
return json.load(f)
return {
"grid": {"w": GRID_W, "h": GRID_H},
"rivers": [],
"oceans": [],
"mountain_ranges": [],
"roads": [],
"cities": [],
"railroads": [],
"pois": [],
}
# ---------------------------------------------------------------------------
# DB queries
# ---------------------------------------------------------------------------
def query_inhabited_bodies(conn: sqlite3.Connection) -> list[dict]:
"""Query systems.db for all inhabited bodies with terrain_reference set."""
rows = conn.execute("""
SELECT
b.body_id,
b.system_id,
b.terrain_reference,
b.population,
b.settlement_pattern,
b.planet_class,
b.economic_role,
COALESCE(b.cultural_corridor, s.cultural_corridor) AS cultural_corridor
FROM bodies b
JOIN star_systems s ON b.system_id = s.system_id
WHERE b.inhabited = 1
AND b.terrain_reference IS NOT NULL
AND b.population > 0
ORDER BY b.system_id, b.body_id
""").fetchall()
return [
{
"body_id": row[0],
"system_id": row[1],
"terrain_reference": row[2],
"population": row[3] or 0,
"settlement_pattern": row[4],
"planet_class": row[5],
"economic_role": row[6],
"cultural_corridor": row[7],
}
for row in rows
]
# ---------------------------------------------------------------------------
# Main pipeline
# ---------------------------------------------------------------------------
def process_body(
body_info: dict,
seed: int,
noise_factor: float,
dry_run: bool,
force: bool,
verbose: bool,
) -> dict:
"""Process one body.
Returns a dict with:
status: 'generated' | 'already_populated' | 'gas_giant' | 'error'
message: human-readable detail (on error)
markers: the markers.json dict if one was produced or loaded (for DB sync)
"""
body_id = body_info["body_id"]
terrain_ref = body_info["terrain_reference"]
population = body_info["population"]
settlement_pattern = body_info["settlement_pattern"]
# Locate body directory from terrain_reference (repo-root relative).
body_dir = REPO_ROOT / Path(terrain_ref).parent
if not body_dir.exists():
return {"status": "error", "message": f"body_dir not found: {body_dir}"}
# Incremental check: if the markers file already has cities, we keep it as
# the source of truth and return it for DB sync instead of regenerating.
if not force:
markers_path = body_dir / "markers.json"
if markers_path.exists():
try:
with open(markers_path) as f:
existing = json.load(f)
if existing.get("cities"):
return {
"status": "already_populated",
"markers": existing,
}
except json.JSONDecodeError:
pass # corrupted markers.json → regenerate below
bd = load_body_def(body_dir)
if not bd:
return {"status": "error", "message": f"no body definition for {body_id}"}
# Seeded RNG: deterministic per (world_seed, body_id)
body_salt = body_id_salt(body_id)
rng = np.random.default_rng(seed ^ body_salt)
try:
terrain = simulate(bd)
except Exception as e:
return {"status": "error", "message": f"simulate failed: {e}"}
if not terrain:
return {"status": "gas_giant"}
analysis = _analyse_terrain(terrain, rng)
n_cities = compute_city_count(population, settlement_pattern)
if verbose:
print(f" {body_id}: pop={population:,} pattern={settlement_pattern} "
f"→ {n_cities} cities, {len(analysis['river_mouths'])} river mouths, "
f"{len(analysis['continents'])} continents")
cities = place_cities(
analysis, n_cities, rng, noise_factor, body_population=population
)
roads, railroads = generate_infrastructure(analysis, cities)
pois = place_gate_terminal(cities, rng)
# Load existing markers (preserves rivers/oceans/mountain_ranges) and
# overlay the new generator output.
markers = load_markers(body_dir)
output_cities = [
{k: v for k, v in city.items() if not k.startswith("_")}
for city in cities
]
markers["cities"] = output_cities
markers["roads"] = roads
markers["railroads"] = railroads
markers["pois"] = pois
if not dry_run:
with open(body_dir / "markers.json", "w") as f:
json.dump(markers, f, indent=2)
return {"status": "generated", "markers": markers}
def main():
parser = argparse.ArgumentParser(
description="Atlas generation — terrain-aware city placement and infrastructure (#832)"
)
parser.add_argument("--db", default=str(DB_PATH), help="Path to systems.db")
parser.add_argument("--body", help="Process only this body_id")
parser.add_argument("--force", action="store_true",
help="Regenerate even if cities already populated")
parser.add_argument("--dry-run", action="store_true",
help="Run analysis but do not write markers.json")
parser.add_argument("--seed", type=int, default=42,
help="World seed for deterministic placement (default: 42)")
parser.add_argument("--noise", type=float, default=0.25,
help="City placement noise factor ±N (default: 0.25 = ±25%%)")
parser.add_argument("--verbose", action="store_true",
help="Print per-body detail")
args = parser.parse_args()
db_path = Path(args.db)
if not db_path.exists():
print(f"error: {db_path} not found", file=sys.stderr)
sys.exit(1)
print(f"\n Atlas Generation Pipeline (#832)")
print(f" DB: {db_path}")
print(f" Seed: {args.seed}")
if args.dry_run:
print(f" Mode: DRY RUN (no markers.json written, no DB sync)")
if args.force:
print(f" Force: enabled (will overwrite existing city placements)")
print()
conn = sqlite3.connect(str(db_path))
conn.execute("PRAGMA foreign_keys=ON")
ensure_atlas_schema(conn)
bodies = query_inhabited_bodies(conn)
if args.body:
bodies = [b for b in bodies if b["body_id"] == args.body]
if not bodies:
print(f"error: body '{args.body}' not found or has no terrain_reference",
file=sys.stderr)
conn.close()
sys.exit(1)
print(f" {len(bodies)} inhabited bodies with terrain_reference\n")
t_total = time.time()
n_generated = 0
n_already = 0
n_gas = 0
n_errors = 0
total_counts = {
"cities": 0, "roads": 0, "railroads": 0, "pois": 0,
"rivers": 0, "oceans": 0, "mountain_ranges": 0,
}
for i, body_info in enumerate(bodies):
body_id = body_info["body_id"]
t0 = time.time()
result = process_body(
body_info,
seed=args.seed,
noise_factor=args.noise,
dry_run=args.dry_run,
force=args.force,
verbose=args.verbose,
)
elapsed = time.time() - t0
status = result["status"]
if status == "generated":
n_generated += 1
if not args.dry_run:
counts = sync_markers_to_db(conn, body_id, result["markers"])
for k, v in counts.items():
total_counts[k] += v
print(f" [{i+1}/{len(bodies)}] {body_id:20s} generated ({elapsed:.1f}s)")
elif status == "already_populated":
n_already += 1
if not args.dry_run:
counts = sync_markers_to_db(conn, body_id, result["markers"])
for k, v in counts.items():
total_counts[k] += v
if args.verbose:
print(f" [{i+1}/{len(bodies)}] {body_id:20s} already populated — DB sync only")
elif status == "gas_giant":
n_gas += 1
if args.verbose:
print(f" [{i+1}/{len(bodies)}] {body_id:20s} gas giant — no surface")
elif status == "error":
n_errors += 1
print(f" [{i+1}/{len(bodies)}] {body_id:20s} ERROR: {result.get('message', '')}")
if not args.dry_run:
conn.commit()
conn.close()
elapsed_total = time.time() - t_total
print(f"\n Done: {elapsed_total:.0f}s")
print(f" generated: {n_generated}")
print(f" already populated: {n_already}")
print(f" gas giants: {n_gas}")
print(f" errors: {n_errors}")
if args.dry_run:
print(f"\n Dry run — no markers.json files were written and no DB rows were touched.")
else:
print(f"\n Atlas index refreshed:")
for k in ("cities", "roads", "railroads", "pois", "rivers", "oceans", "mountain_ranges"):
print(f" atlas_{k:16s} {total_counts[k]:6d} rows")
if n_generated > 0:
print(f"\n {n_generated} markers.json files updated with cities, roads, railroads, pois.")
print(f" Run gemma_naming.py (#833) to fill city and feature names.")
print()
if n_errors > 0:
sys.exit(1)
if __name__ == "__main__":
main()