From ae5b95e88a9c40e802f6318a16a4614168e143f1 Mon Sep 17 00:00:00 2001 From: Lilleman auf Larv Date: Wed, 7 Oct 2026 02:14:12 +0200 Subject: [PATCH] Share the steps every country's geo import takes in geo.py, so geo-se.py and geo-us.py hold only their own rules --- data-import/geo-se.py | 58 +++++-------------------------- data-import/geo-us.py | 24 +++---------- data-import/geo.py | 80 +++++++++++++++++++++++++++++++++++++++++++ todo.md | 5 --- 4 files changed, 94 insertions(+), 73 deletions(-) create mode 100644 data-import/geo.py diff --git a/data-import/geo-se.py b/data-import/geo-se.py index 6a6c8b4..aa1619c 100644 --- a/data-import/geo-se.py +++ b/data-import/geo-se.py @@ -3,12 +3,10 @@ TRAFIKVERKET_API_KEY=… data-import/geo-se.py [--key-file FILE] [--cache DIR] [--streets-per-locality N] [--out DIR] """ -import argparse import collections import csv import io import json -import math import os import re import sys @@ -16,6 +14,7 @@ import urllib.request import zipfile from pathlib import Path +import geo import source import tsv import xlsx @@ -34,8 +33,6 @@ TATORTER = "https://geodata.scb.se/geoserver/stat/wfs?service=WFS&version=2.0.0& POSTAL_CODES = "https://download.geonames.org/export/zip/SE.zip" NVDB = "https://api.trafikinfo.trafikverket.se/v2/data.json" NVDB_PAGE = 50000 -OUT = Path(__file__).resolve().parent.parent / "data" / "geo" / "SE" -CACHE = Path(__file__).resolve().parent / "cache" TIMEZONE = "Europe/Stockholm" # This script's reading of the "exceptions in big cities" in docs/research/research-geo-se.md, Postal codes → Structure. ONE_POSITION = {"Stockholm", "Göteborg", "Malmö"} @@ -106,34 +103,6 @@ def nvdb_segments(cache, key): yield name, float(lat), float(lon) -class Nearest: - """Nearest point by an equirectangular distance, over a degree grid.""" - - def __init__(self, points, cell=0.05): - self.cell = cell - self.grid = collections.defaultdict(list) - for lat, lon, value in points: - self.grid[(int(lat // cell), int(lon // cell))].append((lat, lon, value)) - - def find(self, lat, lon): - ci, cj = int(lat // self.cell), int(lon // self.cell) - best, best_d = None, math.inf - ring = 0 - while ring < 400: - for i in range(ci - ring, ci + ring + 1): - for j in range(cj - ring, cj + ring + 1): - if max(abs(i - ci), abs(j - cj)) != ring: - continue - for plat, plon, value in self.grid.get((i, j), ()): - d = (plat - lat) ** 2 + ((plon - lon) * math.cos(math.radians(lat))) ** 2 - if d < best_d: - best, best_d = value, d - if best is not None and math.sqrt(best_d) < ring * self.cell * math.cos(math.radians(lat)): - return best - ring += 1 - return best - - def street_delivery(name, codes): """The codes delivered to a street. A postort's own prefix is three digits where one three-digit prefix starts at least half its codes, else two. Outside ONE_POSITION's cities, the digit after that prefix is 0 or 1 for boxes and 8 @@ -181,11 +150,10 @@ def localities(codes, tatorter, municipalities, population): municipality, method = municipality_of(name, rows, tatorter, municipalities) how[method] += 1 kept = street_delivery(name, [r["code"].replace(" ", "") for r in rows]) - with_point = [r for r in rows if r["lat"] is not None] - if municipality is None or not kept or not with_point or not well_cased(name): + point = geo.centroid(rows) + if municipality is None or not kept or point is None or not well_cased(name): continue - lat = sum(r["lat"] for r in with_point) / len(with_point) - lon = sum(r["lon"] for r in with_point) / len(with_point) + lat, lon = point out[name] = {"name": name, "municipality": municipality, "population": population_of(name, municipality, tatorter, municipalities, population), "lat": f"{lat:.4f}", "lon": f"{lon:.4f}", "codes": kept} print(f"municipality by {dict(how)}; {len(by_locality) - len(out)} postorter dropped", file=sys.stderr) return out @@ -193,37 +161,29 @@ def localities(codes, tatorter, municipalities, population): def streets(segments, codes, localities, per_locality): """The names with most segments per locality, each segment at its nearest code centroid.""" - nearest = Nearest((r["lat"], r["lon"], r["locality"]) for r in codes if r["lat"] is not None and r["locality"] in localities) + nearest = geo.Nearest((r["lat"], r["lon"], r["locality"]) for r in codes if r["lat"] is not None and r["locality"] in localities) count = collections.Counter() for name, lat, lon in segments: if name[0].isalpha(): count[(nearest.find(lat, lon), name)] += 1 - of = collections.defaultdict(list) - for (locality, name), n in count.items(): - of[locality].append((n, name)) - return {locality: [{"name": name, "locality": locality, "segments": n} for n, name in sorted(named, key=lambda s: (-s[0], s[1]))[:per_locality]] for locality, named in of.items()} + return geo.top_streets(count, per_locality, "segments") def main(): - p = argparse.ArgumentParser(description=__doc__.splitlines()[0]) - p.add_argument("--cache", default=str(CACHE)) + p = geo.parser(__doc__, "SE") p.add_argument("--key-file", help="file holding the Trafikverket API key; TRAFIKVERKET_API_KEY otherwise") - p.add_argument("--out", default=str(OUT)) - p.add_argument("--streets-per-locality", type=int, default=10) a = p.parse_args() key = Path(a.key_file).read_text().strip() if a.key_file else os.environ.get("TRAFIKVERKET_API_KEY") if not key: sys.exit("set TRAFIKVERKET_API_KEY or pass --key-file") - cache, out = Path(a.cache), Path(a.out) - cache.mkdir(parents=True, exist_ok=True) - out.mkdir(parents=True, exist_ok=True) + cache, out = geo.directories(a) regions, municipalities = scb_codes(cache) population = scb_population(cache) codes = geonames(cache) places = localities(codes, scb_tatorter(cache), municipalities, population) named = streets(nvdb_segments(cache, key), codes, places, a.streets_per_locality) - places = {name: l for name, l in places.items() if name in named} + places = geo.with_streets(places, named) empty = sorted(m for m in municipalities if not any(l["municipality"] == m for l in places.values())) if empty: sys.exit(f"municipalities without a locality: {empty}") diff --git a/data-import/geo-us.py b/data-import/geo-us.py index 4ce9c72..b2db9ab 100644 --- a/data-import/geo-us.py +++ b/data-import/geo-us.py @@ -3,7 +3,6 @@ data-import/geo-us.py [--cache DIR] [--min-population N] [--streets-per-locality N] [--out DIR] """ -import argparse import collections import concurrent.futures import csv @@ -12,8 +11,8 @@ import re import struct import sys import zipfile -from pathlib import Path +import geo import source import tsv @@ -24,8 +23,6 @@ COUNTIES = POPULATION.format("counties/totals/co-est2025-alldata.csv") PLACES = POPULATION.format("cities/totals/sub-est2025.csv") ZCTA_PLACE = "https://www2.census.gov/geo/docs/maps-data/data/rel2020/zcta520/tab20_zcta520_place20_natl.txt" TIGER = "https://www2.census.gov/geo/tiger/TIGER2025/{0}/tl_2025_{1}_{2}.zip" -OUT = Path(__file__).resolve().parent.parent / "data" / "geo" / "US" -CACHE = Path(__file__).resolve().parent / "cache" ESTIMATE = "POPESTIMATE2025" CDP = "57" HIGHWAY = re.compile(r"\b(I- |Hwy |Highway |Loop |Rte |Route |Rd )\d") @@ -153,23 +150,14 @@ def streets(cache, counties, locality_of_zcta, per_locality): if r["PAFLAG"] == "P" and r["FULLNAME"] and not HIGHWAY.search(r["FULLNAME"]): for z in zips.get(r["TLID"], ()): count[(locality_of_zcta[z], r["FULLNAME"])] += 1 - of = collections.defaultdict(list) - for (locality, name), n in count.items(): - of[locality].append((n, name)) - named = {locality: [{"name": name, "locality": locality, "addresses": n} for n, name in sorted(ranked, key=lambda s: (-s[0], s[1]))[:per_locality]] for locality, ranked in of.items()} - return named, addresses + return geo.top_streets(count, per_locality, "addresses"), addresses def main(): - p = argparse.ArgumentParser(description=__doc__.splitlines()[0]) - p.add_argument("--cache", default=str(CACHE)) + p = geo.parser(__doc__, "US") p.add_argument("--min-population", type=int, default=25000) - p.add_argument("--out", default=str(OUT)) - p.add_argument("--streets-per-locality", type=int, default=10) a = p.parse_args() - cache, out = Path(a.cache), Path(a.out) - cache.mkdir(parents=True, exist_ok=True) - out.mkdir(parents=True, exist_ok=True) + cache, out = geo.directories(a) state_population = {r["STATE"]: r[ESTIMATE] for r in csv_rows(cache, STATES, "nst-est2025.csv") if r["SUMLEV"] == "040"} regions = {r["USPS"]: {"abbr": r["USPS"], "code": r["GEOID"], "name": r["NAME"], "population": state_population[r["GEOID"]], "timezone": TIMEZONES[r["USPS"]]} for r in gazetteer(cache, "state")} @@ -185,9 +173,7 @@ def main(): places = localities(cache, a.min_population, counties) locality_of_zcta = postal_codes(cache, places) named, addresses = streets(cache, sorted({l["municipality"] for l in places.values()}), locality_of_zcta, a.streets_per_locality) - for geoid in [l for l in places if l not in named]: - print(f"{geoid} {places[geoid]['name']}: no streets, dropped", file=sys.stderr) - del places[geoid] + places = geo.with_streets(places, named) kept_counties = {l["municipality"] for l in places.values()} kept_regions = {counties[c]["region"] for c in kept_counties} diff --git a/data-import/geo.py b/data-import/geo.py new file mode 100644 index 0000000..7f743b5 --- /dev/null +++ b/data-import/geo.py @@ -0,0 +1,80 @@ +"""The steps every country's geo import takes.""" +import argparse +import collections +import math +import sys +from pathlib import Path + +DATA = Path(__file__).resolve().parent.parent / "data" / "geo" +CACHE = Path(__file__).resolve().parent / "cache" + + +def parser(doc, country): + """The options every country's import takes; a script adds its own before parse_args.""" + p = argparse.ArgumentParser(description=doc.splitlines()[0]) + p.add_argument("--cache", default=str(CACHE)) + p.add_argument("--out", default=str(DATA / country)) + p.add_argument("--streets-per-locality", type=int, default=10) + return p + + +def directories(args): + """The cache and out directories of the parsed options, both created.""" + cache, out = Path(args.cache), Path(args.out) + cache.mkdir(parents=True, exist_ok=True) + out.mkdir(parents=True, exist_ok=True) + return cache, out + + +def centroid(rows): + """The mean lat and lon of the rows with a point, or None where no row has one.""" + points = [r for r in rows if r["lat"] is not None] + if not points: + return None + return sum(r["lat"] for r in points) / len(points), sum(r["lon"] for r in points) / len(points) + + +class Nearest: + """Nearest point by an equirectangular distance, over a degree grid.""" + + def __init__(self, points, cell=0.05): + self.cell = cell + self.grid = collections.defaultdict(list) + for lat, lon, value in points: + self.grid[(int(lat // cell), int(lon // cell))].append((lat, lon, value)) + + def find(self, lat, lon): + ci, cj = int(lat // self.cell), int(lon // self.cell) + best, best_d = None, math.inf + ring = 0 + while ring < 400: + for i in range(ci - ring, ci + ring + 1): + for j in range(cj - ring, cj + ring + 1): + if max(abs(i - ci), abs(j - cj)) != ring: + continue + for plat, plon, value in self.grid.get((i, j), ()): + d = (plat - lat) ** 2 + ((plon - lon) * math.cos(math.radians(lat))) ** 2 + if d < best_d: + best, best_d = value, d + if best is not None and math.sqrt(best_d) < ring * self.cell * math.cos(math.radians(lat)): + return best + ring += 1 + return best + + +def top_streets(count, per_locality, column): + """Per locality, the per_locality street names counted most in count, keyed (locality, name), ties by name; each + row holds its count in column.""" + of = collections.defaultdict(list) + for (locality, name), n in count.items(): + of[locality].append((n, name)) + return {locality: [{"name": name, "locality": locality, column: n} for n, name in sorted(named, key=lambda s: (-s[0], s[1]))[:per_locality]] for locality, named in of.items()} + + +def with_streets(places, named): + """The places named holds streets for; each other one is logged as dropped.""" + for key, place in places.items(): + if key not in named: + label = key if key == place["name"] else f"{key} {place['name']}" + print(f"{label}: no streets, dropped", file=sys.stderr) + return {key: place for key, place in places.items() if key in named} diff --git a/todo.md b/todo.md index 0fd1d2e..c53dd5d 100644 --- a/todo.md +++ b/todo.md @@ -77,7 +77,6 @@ | 144 | 0.1.0 | decision | **Run the nine-seat comprehension panel and the scoring run again once items 130 to 143 have merged, and file what the panel names.** | 1 | 3 | 3 | 9 | 3 | 20.4 | | 91 | 0.1.0 | | **Accept a path through a name as a `calc` operand, so `{/misc.coordinate as c}{c}: {calc(c.lat * 60)}` computes from the coordinate `{c}` prints.** | 3 | 4 | 6 | 7 | 4.1, 5.3 | 19.9 | | 89 | 0.1.0 | | **Fail `New` without `WithSeed`, naming the option, and seed the CLI from the system.** | 3 | 3 | 6 | 7 | 10.1 | 16.0 | -| 138 | 0.1.0 | | **Split `geo-se.py` into the steps every country's import takes and the Swedish rules, and share the steps with `geo-us.py`, so a new country's script holds only its own rules.** | 2 | 3 | 3 | 6 | 2.1, 3.2 | 14.6 | | 127 | 0.1.0 | | **Give a Swedish street only postal codes whose delivery digit means street or rural delivery, dropping codes whose delivery digit is 9, or record why they stay.** | 3 | 3 | 4 | 6 | 4 | 14.4 | | 110 | 0.1.0 | | **Refuse an empty segment in a path `Fake`, `FakeRecord` and the CLI take, as a reference does: `sv_SE.person.` reports `path has an empty segment`.** | 1 | 1 | 3 | 6 | 6 | 13.8 | | 79 | 0.1.0 | decision | **Beside a `--data-path`, load in `New` the user's categories, the shipped ones they override or read, and the shipped ones whose reads reach an overridden one; load the rest on first reach.** | 6 | 6 | 7 | 8 | 13 | 12.9 | @@ -271,10 +270,6 @@ A calc operand parses as a bare identifier, so `{calc(c.lat * 60)}` fails with ` Goal 10.1 has the caller supply the seed, and today `New` without `WithSeed` seeds itself. This item ships in 0.1.0 beside item 77, so `New` breaks once before anyone depends on it. Mark `WithSeed` required in the Library section's options table. Its tests then stop swapping `randomBytes`, as `fejkdata_test.go` does today, which fails every `New` in the package while that test runs. -### 138. Split `geo-se.py` into the steps every country's import takes and the Swedish rules, and share the steps with `geo-us.py`, so a new country's script holds only its own rules. - -Needed by item 61. The extender seat of the 2026-10-07 panel, reading for item 61, could not tell in `localities` and its helpers which steps a Norwegian script keeps; `geo-us.py` copies the per-locality street ranking of `streets` in `geo-se.py`. - ### 127. Give a Swedish street only postal codes whose delivery digit means street or rural delivery, dropping codes whose delivery digit is 9, or record why they stay. `docs/research/research-geo-se.md` reads delivery digit 0 or 1 as boxes, 2–4 and 6–7 as street, 5 as rural, 8 as reply mail and 9 as competitions or temporary, with exceptions in big cities. `street_delivery` drops 0, 1 and 8 and keeps 9. The panel of 2026-10-05 found it; whether a shipped street carries a 9 code is unchecked.