Replace the per-cell parquet cache with TimescaleDB hypertables as the production backend for the raw daily climate record, and drop pg_duckdb. Parquet stays the backend whenever THERMOGRAPH_DATABASE_URL is not a Postgres URL (dev, tests, offline tooling), the same dialect switch the accounts DB and derived store already use, so CI stays Postgres-free. - data/climate_store.py: psycopg + polars bridge over climate_history (a hypertable), climate_recent, and climate_sync (per-cell freshness). Reads via pl.read_database, writes via COPY + ON CONFLICT upsert; fail-soft to a cache miss so a DB hiccup degrades to upstream refetch. - data/climate.py: route every cache/mtime touchpoint through a backend dispatch. recent_stamp becomes int(recent_synced_at) on Postgres; the stale-serve path still avoids bumping it, so derived-payload tokens invalidate on exactly the same events as before. - alembic 0002: CREATE EXTENSION timescaledb plus the hypertable schema (compression policy on year-old chunks), guarded to no-op off Postgres. - migrate_cache_to_pg.py (make migrate-cache): idempotent backfill of the parquet cache into the hypertables, preserving file mtimes as sync timestamps so recent_stamp is unchanged across cutover. - db image -> stock timescale/timescaledb:latest-pg18; drop the custom pg_duckdb Dockerfile, the read-only /parquet mount, and the duckdb tuning GUC. Docs updated for the new backend and cutover. Co-authored-by: Claude <noreply@anthropic.com>
905 lines
41 KiB
Python
905 lines
41 KiB
Python
"""Fetch historical + recent daily weather from Open-Meteo and cache to parquet.
|
||
|
||
The full multi-decade daily record for a grid cell is fetched once and stored as
|
||
a zstd-compressed parquet file keyed by cell id. Percentiles are computed on the
|
||
fly from this raw record (see grading.py), which keeps the cache small and lets
|
||
us handle ties (e.g. many zero-precip days) correctly.
|
||
"""
|
||
import datetime
|
||
import os
|
||
import threading
|
||
import time
|
||
|
||
import httpx
|
||
import numpy as np
|
||
import polars as pl
|
||
|
||
from core import audit
|
||
from core import metrics
|
||
import paths
|
||
from data import climate_store
|
||
from data import store
|
||
|
||
CACHE_DIR = os.path.join(paths.DATA_DIR, "cache")
|
||
MAX_ATTEMPTS = 3 # per upstream call, before giving up
|
||
START_DATE = "1980-01-01" # ERA5 reaches back to 1940; 1980 = 45 yrs, fast + robust
|
||
ARCHIVE_LATENCY_DAYS = 6 # ERA5 archive lags real time by a few days
|
||
# A real archive spans decades (~16k daily rows from START_DATE). A self-hosted
|
||
# Open-Meteo instance that isn't backfilled yet answers with all-null values or only
|
||
# its recent sync window, which _finalize_frame reduces to a near-empty frame. Below
|
||
# this many days we treat the archive as "not ready" — fall through to the NASA backup
|
||
# and do NOT cache the short result as a complete, indefinitely-held record.
|
||
MIN_ARCHIVE_DAYS = 3650 # ~10 yrs: far above any sync window, far below a full archive
|
||
# History rarely changes, so the full multi-decade archive is fetched once and
|
||
# cached indefinitely (refetched only to add new metric columns). Only the recent
|
||
# tail is refreshed — a small incremental fetch, at most hourly.
|
||
HISTORY_TOPUP_INTERVAL = 3600 # seconds between recent-tail refresh attempts
|
||
FORECAST_TTL_HOURS = 1 # refetch the forward forecast hourly to track updates
|
||
|
||
# The archive (historical) endpoint. Overridable so it can point at a self-hosted
|
||
# Open-Meteo instance (ERA5 in object storage) instead of the rate-limited public
|
||
# API; the response shape is identical either way.
|
||
ARCHIVE_URL = os.environ.get("THERMOGRAPH_ARCHIVE_URL",
|
||
"https://archive-api.open-meteo.com/v1/archive")
|
||
# The historical model: the ERA5 "seamless" blend — 0.1° ERA5-Land for
|
||
# temperature/precip/humidity/wind, 0.25° ERA5 for wind gusts (which ERA5-Land
|
||
# lacks). This is the public API's default; sent explicitly so a self-hosted
|
||
# instance serves the same 0.1° resolution. Archive requests only, not forecast.
|
||
ARCHIVE_MODEL = "era5_seamless"
|
||
FORECAST_URL = "https://api.open-meteo.com/v1/forecast"
|
||
# Backup archive when Open-Meteo is unavailable (e.g. its daily rate limit). NASA
|
||
# POWER is free + keyless, global, daily from 1981. It lacks gusts + apparent temp,
|
||
# so gusts read as unavailable and "feels like" is computed from heat index/chill.
|
||
NASA_POWER_URL = "https://power.larc.nasa.gov/api/temporal/daily/point"
|
||
NASA_START = "19810101"
|
||
NASA_FILL = -900.0 # POWER's missing-value sentinel is ~-999
|
||
|
||
# Backup forward forecast when Open-Meteo's forecast API is unavailable. MET Norway
|
||
# (yr.no) is free + keyless + global, mirroring NASA POWER's role for history. It
|
||
# returns a sub-daily timeseries in metric units with no gusts or apparent temp, so
|
||
# we aggregate to daily, convert units, and derive feels-like like the NASA path.
|
||
# It is forecast-only — no recent past days — so it's a degraded-but-working fallback.
|
||
METNO_URL = "https://api.met.no/weatherapi/locationforecast/2.0/complete"
|
||
# MET Norway's ToS requires an identifying User-Agent (a missing/generic one is
|
||
# 403'd); include the app and a contact URL so they can reach us if usage misbehaves.
|
||
METNO_UA = "Thermograph/0.2 (+https://thermograph.org)"
|
||
|
||
DAILY_VARS = (
|
||
"temperature_2m_max,temperature_2m_min,precipitation_sum,"
|
||
"wind_speed_10m_max,wind_gusts_10m_max,"
|
||
"apparent_temperature_max,apparent_temperature_min,"
|
||
"relative_humidity_2m_mean"
|
||
)
|
||
|
||
# Columns added after the original tmax/tmin/precip schema. A cached parquet
|
||
# missing any of these predates the wind/feels/humidity features (or the split
|
||
# apparent high/low) and is refetched so the new metrics show up immediately
|
||
# instead of after the 30-day cache TTL.
|
||
NEW_COLS = ("wind", "gust", "feels", "humid", "fmax", "fmin")
|
||
|
||
# Thermoneutral baseline (°F). The combined "feels like" metric reports whichever
|
||
# apparent-temperature extreme — the heat-index-driven daily max or the
|
||
# wind-chill-driven daily min — sits further from this comfort point, so a single
|
||
# daily value captures both heat index (hot side) and wind chill (cold side).
|
||
COMFORT_F = 65.0
|
||
|
||
|
||
def _cache_path(cell_id: str) -> str:
|
||
return os.path.join(CACHE_DIR, f"{cell_id}.parquet")
|
||
|
||
|
||
# One lock per cell so concurrent requests/refreshes don't each pull the full
|
||
# multi-decade archive (which quickly trips the upstream rate limit).
|
||
_LOCKS_GUARD = threading.Lock()
|
||
_CELL_LOCKS: dict[str, threading.Lock] = {}
|
||
|
||
# When the archive rate-limits us (429), back off globally for a bit rather than
|
||
# re-hitting it on every refresh — that only prolongs the limit.
|
||
ARCHIVE_COOLDOWN = 120 # seconds
|
||
_archive_cooldown_until = 0.0
|
||
|
||
|
||
_archive_limit_daily = False # is the current cooldown a daily-quota exhaustion?
|
||
|
||
|
||
class WeatherUnavailable(RuntimeError):
|
||
"""Weather data can't be fetched right now (rate limit, upstream outage with
|
||
no cached fallback). Carries user-facing text; ``daily`` marks a daily-quota
|
||
exhaustion (resets after UTC midnight) rather than a transient burst limit.
|
||
The API layer maps this to a retryable 503 — it never needs to parse the
|
||
message, and new upstream sources classify their own failures here."""
|
||
|
||
def __init__(self, message: str, daily: bool = False):
|
||
super().__init__(message)
|
||
self.daily = daily
|
||
|
||
|
||
def limit_message(daily: bool) -> str:
|
||
"""The user-facing rate-limit copy, in one place."""
|
||
if daily:
|
||
return ("Open-Meteo's daily request limit is exhausted, so new locations will work again "
|
||
"tomorrow. Places you've already viewed still work.")
|
||
return "The weather service is rate-limited right now. Please try again in a minute."
|
||
|
||
|
||
def is_rate_limit(e) -> bool:
|
||
return getattr(getattr(e, "response", None), "status_code", None) == 429
|
||
|
||
|
||
def _rate_limit_reason(e) -> str:
|
||
"""Upstream's human-readable 429 reason, if present in the response body."""
|
||
resp = getattr(e, "response", None)
|
||
if resp is not None:
|
||
try:
|
||
j = resp.json()
|
||
if isinstance(j, dict) and j.get("reason"):
|
||
return str(j["reason"])
|
||
except Exception: # noqa: BLE001
|
||
pass
|
||
return ""
|
||
|
||
|
||
def _seconds_to_utc_reset() -> float:
|
||
"""Seconds until just after the next UTC midnight (when daily quotas reset)."""
|
||
now = datetime.datetime.now(datetime.timezone.utc)
|
||
reset = (now + datetime.timedelta(days=1)).replace(hour=0, minute=10, second=0, microsecond=0)
|
||
return max((reset - now).total_seconds(), 600)
|
||
|
||
|
||
def _note_rate_limit(e) -> None:
|
||
"""Record a rate limit: back off ~2 min, or until tomorrow for a daily quota."""
|
||
global _archive_cooldown_until, _archive_limit_daily
|
||
reason = _rate_limit_reason(e).lower()
|
||
_archive_limit_daily = "daily" in reason or "tomorrow" in reason
|
||
_archive_cooldown_until = time.time() + (
|
||
_seconds_to_utc_reset() if _archive_limit_daily else ARCHIVE_COOLDOWN)
|
||
|
||
|
||
def _cell_lock(cell_id: str) -> threading.Lock:
|
||
with _LOCKS_GUARD:
|
||
lk = _CELL_LOCKS.get(cell_id)
|
||
if lk is None:
|
||
lk = _CELL_LOCKS[cell_id] = threading.Lock()
|
||
return lk
|
||
|
||
|
||
def _request(url, params, timeout, *, phase, headers=None, attempts=MAX_ATTEMPTS):
|
||
"""GET with bounded retries; every retry and the final failure are logged to
|
||
the errors folder (tagged ``retry`` / ``error``). A 429 (rate limit) fails fast
|
||
without retrying, so we don't hammer the limit and make it worse."""
|
||
last = None
|
||
for attempt in range(1, attempts + 1):
|
||
try:
|
||
r = httpx.get(url, params=params, timeout=timeout, headers=headers)
|
||
r.raise_for_status()
|
||
metrics.record_outbound(phase, "ok")
|
||
return r
|
||
except Exception as e: # noqa: BLE001 - upstream/network failures are expected
|
||
last = e
|
||
status = getattr(getattr(e, "response", None), "status_code", None)
|
||
rate_limited = status == 429
|
||
final = rate_limited or attempt == attempts
|
||
metrics.record_outbound(
|
||
phase, "rate_limited" if rate_limited else ("error" if final else "retry"))
|
||
audit.log_event(
|
||
"error" if final else "retry",
|
||
{"phase": phase, "attempt": attempt, "max_attempts": attempts,
|
||
"url": url, "status": status, "error": repr(e)},
|
||
)
|
||
if final:
|
||
break
|
||
time.sleep(min(0.5 * 2 ** (attempt - 1), 4.0))
|
||
raise last
|
||
|
||
|
||
def _normalize_read(df: pl.DataFrame) -> pl.DataFrame:
|
||
"""Normalize a frame read from the parquet cache: the daily `date` column is a
|
||
calendar date, so coerce it to ``pl.Date`` (older files were written by pandas
|
||
as a nanosecond ``Datetime``). Downstream date math and comparisons all key on
|
||
``pl.Date`` / stdlib ``datetime.date``."""
|
||
if df.schema["date"] != pl.Date:
|
||
df = df.with_columns(pl.col("date").cast(pl.Date))
|
||
return df
|
||
|
||
|
||
def _combined_feels_expr(hi: str = "fmax", lo: str = "fmin") -> pl.Expr:
|
||
"""Expression for one daily "feels like" value: the apparent-temperature extreme
|
||
furthest from the comfort baseline — the heat-index high on warm days, the
|
||
wind-chill low on cold ones. Falls back to whichever side is present if one is
|
||
missing (the coalesce picks the non-null side when a comparison is null)."""
|
||
amax, amin = pl.col(hi), pl.col(lo)
|
||
hot, cold = amax - COMFORT_F, COMFORT_F - amin
|
||
chosen = pl.when(hot >= cold).then(amax).otherwise(amin)
|
||
return pl.coalesce([chosen, amax, amin])
|
||
|
||
|
||
def _derive_humidity(df: pl.DataFrame) -> pl.DataFrame:
|
||
"""Replace the raw mean *relative* humidity column with *absolute* humidity
|
||
(grams of water vapor per m³).
|
||
|
||
Absolute humidity is derived from the day's mean RH and its mean temperature
|
||
((tmax+tmin)/2) via the Magnus saturation-vapor-pressure formula. It's a far
|
||
more informative "how muggy was it" signal than relative humidity, which mostly
|
||
tracks the day/night temperature swing (cold nights read ~100% RH regardless of
|
||
actual moisture). Applied at the read boundary so the parquet cache keeps the
|
||
raw RH the archive returns — no refetch needed for existing cells."""
|
||
if "humid" not in df.columns:
|
||
return df
|
||
tmean_c = ((pl.col("tmax") + pl.col("tmin")) / 2.0 - 32.0) * 5.0 / 9.0
|
||
rh = pl.col("humid").cast(pl.Float64, strict=False)
|
||
es = 6.112 * (17.67 * tmean_c / (tmean_c + 243.5)).exp() # sat. vapor pressure, hPa
|
||
return df.with_columns(
|
||
(es * rh * 2.1674 / (273.15 + tmean_c)).round(1).alias("humid")) # abs. humidity, g/m³
|
||
|
||
|
||
def _derive_wetbulb(df: pl.DataFrame) -> pl.DataFrame:
|
||
"""Add a `wetbulb` column (°F): the daytime-peak wet-bulb temperature via the
|
||
Stull (2011) single-value approximation from the day's high (`tmax`) and mean
|
||
relative humidity.
|
||
|
||
Wet bulb is the temperature a parcel reaches by evaporative cooling to
|
||
saturation — the ceiling on how much the body can shed heat by sweating, so it
|
||
is the sharper heat-stress signal than dry-bulb temperature or humidity alone.
|
||
Must run BEFORE ``_derive_humidity`` (which replaces the raw RH column with
|
||
absolute humidity, the input this needs). Read-time only, like the humidity
|
||
derivation, so the parquet cache is untouched. Stull's fit is valid for
|
||
RH 5-99% and −20..50 °C; days outside that range grade as null."""
|
||
if "humid" not in df.columns:
|
||
return df
|
||
t = (pl.col("tmax") - 32.0) * 5.0 / 9.0 # daily high, °C
|
||
rh = pl.col("humid").cast(pl.Float64, strict=False)
|
||
tw = (t * (0.151977 * (rh + 8.313659).sqrt()).arctan()
|
||
+ (t + rh).arctan()
|
||
- (rh - 1.676331).arctan()
|
||
+ 0.00391838 * rh.pow(1.5) * (0.023101 * rh).arctan()
|
||
- 4.686035) # wet bulb, °C
|
||
tw_f = (tw * 9.0 / 5.0 + 32.0).round(1)
|
||
valid = (rh >= 5) & (rh <= 99) & (t >= -20) & (t <= 50)
|
||
return df.with_columns(pl.when(valid).then(tw_f).otherwise(None).alias("wetbulb"))
|
||
|
||
|
||
def _derive_metrics(df: pl.DataFrame) -> pl.DataFrame:
|
||
"""Read-boundary derivations that depend on the raw relative-humidity column.
|
||
Wet bulb must be computed before ``_derive_humidity`` replaces raw RH with
|
||
absolute humidity, so both live behind this single wrapper to keep the order
|
||
right at every read site."""
|
||
return _derive_humidity(_derive_wetbulb(df))
|
||
|
||
|
||
def _with_doy(df: pl.DataFrame) -> pl.DataFrame:
|
||
"""(Re)attach the int16 day-of-year column the grading windows key on."""
|
||
return df.with_columns(pl.col("date").dt.ordinal_day().cast(pl.Int16).alias("doy"))
|
||
|
||
|
||
def _write_cache(df: pl.DataFrame, path: str) -> None:
|
||
"""Persist a daily record to the parquet cache — the raw record only, never
|
||
the derived doy column (it's recomputed at read time)."""
|
||
os.makedirs(CACHE_DIR, exist_ok=True)
|
||
df.drop("doy", strict=False).write_parquet(path, compression="zstd")
|
||
|
||
|
||
# --- cache backend dispatch --------------------------------------------------
|
||
# The raw daily record lives in TimescaleDB hypertables on Postgres (prod) and in
|
||
# per-cell parquet files otherwise (dev / tests / offline). These helpers hide the
|
||
# choice so the loaders below stay backend-agnostic: they return the same shapes a
|
||
# parquet read produced — a frame WITHOUT the derived `doy` column, plus an age in
|
||
# seconds since the last write (the mtime age the topup/TTL logic keyed on). The
|
||
# active backend is climate_store.is_postgres() (same THERMOGRAPH_DATABASE_URL
|
||
# switch as accounts/db.py and data/store.py).
|
||
|
||
|
||
def _read_history_backed(cell_id: str):
|
||
"""(schema-complete history frame without doy, age_s) or None."""
|
||
if climate_store.is_postgres():
|
||
hit = climate_store.read_history(cell_id)
|
||
if hit is None:
|
||
return None
|
||
df, synced_at = hit
|
||
return df, time.time() - synced_at
|
||
return _read_history_cache(_cache_path(cell_id))
|
||
|
||
|
||
def _stale_history_backed(cell_id: str):
|
||
"""Any cached history rows regardless of schema completeness (for the
|
||
all-fetches-failed stale serve), or None."""
|
||
if climate_store.is_postgres():
|
||
return climate_store.read_history_raw(cell_id)
|
||
path = _cache_path(cell_id)
|
||
return _normalize_read(pl.read_parquet(path)) if os.path.exists(path) else None
|
||
|
||
|
||
def _write_history_backed(cell_id: str, full_df: pl.DataFrame,
|
||
delta_df: "pl.DataFrame | None" = None) -> None:
|
||
"""Persist a cell's history. Parquet rewrites the whole file from ``full_df``;
|
||
Postgres upserts ``delta_df`` when given (just the topped-up tail) else the full
|
||
frame — the upsert makes a delta write equivalent to a full rewrite."""
|
||
if climate_store.is_postgres():
|
||
climate_store.write_history(
|
||
cell_id, delta_df if delta_df is not None else full_df, time.time())
|
||
else:
|
||
_write_cache(full_df, _cache_path(cell_id))
|
||
|
||
|
||
def _touch_history_backed(cell_id: str) -> None:
|
||
"""Reset the hourly topup timer when the tail is already current, without a
|
||
rewrite — bump history_synced_at on Postgres, or the file mtime on parquet."""
|
||
if climate_store.is_postgres():
|
||
climate_store.touch_history(cell_id, time.time())
|
||
else:
|
||
try:
|
||
os.utime(_cache_path(cell_id), None)
|
||
except OSError:
|
||
pass
|
||
|
||
|
||
def _read_recent_backed(cell_id: str):
|
||
"""(recent+forecast frame without doy, age_s) or None."""
|
||
if climate_store.is_postgres():
|
||
hit = climate_store.read_recent(cell_id)
|
||
if hit is None:
|
||
return None
|
||
df, synced_at = hit
|
||
return df, time.time() - synced_at
|
||
path = _rf_cache_path(cell_id)
|
||
if not os.path.exists(path):
|
||
return None
|
||
try:
|
||
return _normalize_read(pl.read_parquet(path)), time.time() - os.path.getmtime(path)
|
||
except Exception: # noqa: BLE001 - a truncated/corrupt cache file is a miss
|
||
return None
|
||
|
||
|
||
def _write_recent_backed(cell_id: str, df: pl.DataFrame) -> None:
|
||
"""Persist a cell's recent+forecast bundle (full replace on both backends)."""
|
||
if climate_store.is_postgres():
|
||
climate_store.write_recent(cell_id, df, time.time())
|
||
else:
|
||
_write_cache(df, _rf_cache_path(cell_id))
|
||
|
||
|
||
def _om_daily_params(cell: dict, **window) -> dict:
|
||
"""The Open-Meteo daily-request params every fetch shares — location, the
|
||
variable set, imperial units. The date/window selectors (start_date/end_date
|
||
or past_days/forecast_days) come in as kwargs."""
|
||
return {
|
||
"latitude": cell["center_lat"],
|
||
"longitude": cell["center_lon"],
|
||
"daily": DAILY_VARS,
|
||
"timezone": "auto",
|
||
"temperature_unit": "fahrenheit",
|
||
"precipitation_unit": "inch",
|
||
"wind_speed_unit": "mph",
|
||
**window,
|
||
}
|
||
|
||
|
||
def _finalize_frame(df: pl.DataFrame) -> pl.DataFrame:
|
||
"""Shared tail of every source→frame mapping: unify missing values as null, add
|
||
the combined feels-like, apply the valid-day filter, and attach day-of-year. A
|
||
usable climate day needs a real high/low; other columns may be missing and
|
||
simply grade as None. Float NaN (from numpy-derived columns) is folded into null
|
||
so the whole pipeline models "missing" one way — the grading boundary drops it."""
|
||
df = df.with_columns(pl.col(pl.Float32, pl.Float64).fill_nan(None))
|
||
df = df.with_columns(_combined_feels_expr().alias("feels"))
|
||
df = df.drop_nulls(subset=["tmax", "tmin"])
|
||
return _with_doy(df)
|
||
|
||
|
||
def _to_frame(daily: dict) -> pl.DataFrame:
|
||
n = len(daily["time"])
|
||
# Older upstream responses (or a narrowed variable set) may omit a series; fall
|
||
# back to an all-null column of the right length so the frame shape is stable.
|
||
def col(key):
|
||
vals = daily.get(key)
|
||
return vals if vals else [None] * n
|
||
|
||
df = pl.DataFrame(
|
||
{
|
||
"date": daily["time"],
|
||
"tmax": daily["temperature_2m_max"],
|
||
"tmin": daily["temperature_2m_min"],
|
||
"precip": daily["precipitation_sum"],
|
||
"wind": col("wind_speed_10m_max"),
|
||
"gust": col("wind_gusts_10m_max"),
|
||
"humid": col("relative_humidity_2m_mean"),
|
||
# The apparent (felt) high and low are kept as their own columns — graded
|
||
# independently — alongside the combined `feels` (whichever side is further
|
||
# from the fixed comfort baseline), which the weekly/day views still use.
|
||
"fmax": col("apparent_temperature_max"),
|
||
"fmin": col("apparent_temperature_min"),
|
||
}
|
||
)
|
||
df = df.with_columns(
|
||
pl.col("date").str.to_date(),
|
||
*[pl.col(c).cast(pl.Float64, strict=False)
|
||
for c in ("tmax", "tmin", "precip", "wind", "gust", "fmax", "fmin")],
|
||
)
|
||
return _finalize_frame(df)
|
||
|
||
|
||
def _fetch_history(cell: dict) -> pl.DataFrame:
|
||
end = (datetime.date.today() - datetime.timedelta(days=ARCHIVE_LATENCY_DAYS)).isoformat()
|
||
params = _om_daily_params(cell, start_date=START_DATE, end_date=end, models=ARCHIVE_MODEL)
|
||
r = _request(ARCHIVE_URL, params, 180, phase="history_fetch")
|
||
return _to_frame(r.json()["daily"])
|
||
|
||
|
||
def _heat_index(t_f, rh):
|
||
"""NWS heat index (°F); below ~80°F it's just the air temperature."""
|
||
T, R = t_f, rh
|
||
hi = (-42.379 + 2.04901523 * T + 10.14333127 * R - 0.22475541 * T * R
|
||
- 0.00683783 * T * T - 0.05481717 * R * R + 0.00122874 * T * T * R
|
||
+ 0.00085282 * T * R * R - 0.00000199 * T * T * R * R)
|
||
return np.where(np.asarray(T, dtype="float64") >= 80, hi, T)
|
||
|
||
|
||
def _wind_chill(t_f, v_mph):
|
||
"""NWS wind chill (°F); applies only when cold + breezy, else the air temp."""
|
||
T = np.asarray(t_f, dtype="float64")
|
||
V = np.clip(np.asarray(v_mph, dtype="float64"), 0, None)
|
||
Vp = np.power(V, 0.16)
|
||
wc = 35.74 + 0.6215 * T - 35.75 * Vp + 0.4275 * T * Vp
|
||
return np.where((T <= 50) & (V >= 3), wc, T)
|
||
|
||
|
||
def _nasa_to_frame(param: dict) -> pl.DataFrame:
|
||
"""Map a NASA POWER daily response to our (tmax/tmin/precip/wind/gust/humid/feels)
|
||
schema, converting units (°C→°F, mm→in, m/s→mph) and computing feels-like."""
|
||
dates = sorted(param.get("T2M_MAX", {}).keys())
|
||
|
||
def col(key, transform):
|
||
d = param.get(key, {})
|
||
return [None if (v is None or v <= NASA_FILL) else transform(v)
|
||
for v in (d.get(k) for k in dates)]
|
||
|
||
c2f = lambda c: c * 9.0 / 5.0 + 32.0
|
||
df = pl.DataFrame({
|
||
"date": dates,
|
||
"tmax": col("T2M_MAX", c2f),
|
||
"tmin": col("T2M_MIN", c2f),
|
||
"precip": col("PRECTOTCORR", lambda mm: mm / 25.4),
|
||
"wind": col("WS10M_MAX", lambda ms: ms * 2.2369362920544),
|
||
"gust": [None] * len(dates), # POWER has no gusts
|
||
"humid": col("RH2M", lambda v: v),
|
||
})
|
||
df = df.with_columns(
|
||
pl.col("date").str.to_date("%Y%m%d"),
|
||
*[pl.col(c).cast(pl.Float64, strict=False)
|
||
for c in ("tmax", "tmin", "precip", "wind", "gust")],
|
||
)
|
||
# POWER has no apparent temperature; approximate the felt high with the NWS
|
||
# heat index and the felt low with NWS wind chill (each falls back to the air
|
||
# temperature outside its regime), mirroring the Open-Meteo columns. numpy NaN
|
||
# from these lands in float columns and is folded to null in _finalize_frame.
|
||
df = df.with_columns(
|
||
pl.Series("fmax", _heat_index(df["tmax"].to_numpy(), df["humid"].to_numpy())),
|
||
pl.Series("fmin", _wind_chill(df["tmin"].to_numpy(), df["wind"].to_numpy())),
|
||
)
|
||
return _finalize_frame(df)
|
||
|
||
|
||
def _fetch_history_nasa(cell: dict) -> pl.DataFrame:
|
||
"""Backup history fetch from NASA POWER (used when Open-Meteo is unavailable)."""
|
||
end = (datetime.date.today() - datetime.timedelta(days=ARCHIVE_LATENCY_DAYS)).strftime("%Y%m%d")
|
||
params = {
|
||
"parameters": "T2M_MAX,T2M_MIN,PRECTOTCORR,WS10M_MAX,RH2M",
|
||
"community": "RE",
|
||
"latitude": cell["center_lat"],
|
||
"longitude": cell["center_lon"],
|
||
"start": NASA_START,
|
||
"end": end,
|
||
"format": "JSON",
|
||
}
|
||
r = _request(NASA_POWER_URL, params, 180, phase="history_nasa")
|
||
return _nasa_to_frame(r.json()["properties"]["parameter"])
|
||
|
||
|
||
def _metno_to_frame(props: dict) -> pl.DataFrame:
|
||
"""Map a MET Norway Locationforecast (yr.no) response to our daily schema.
|
||
|
||
MET Norway returns a per-timestep timeseries (hourly near-term, then 6-hourly)
|
||
in metric units, with no gusts or apparent temperature, so we aggregate to daily
|
||
extremes/means, convert units (°C→°F, mm→in, m/s→mph), and approximate feels-like
|
||
from the NWS heat index / wind chill — the same treatment as the NASA POWER path."""
|
||
# Aggregate the sub-daily steps into per-day values, keyed by (UTC) calendar day.
|
||
agg: dict[str, dict] = {}
|
||
for step in props.get("timeseries", []):
|
||
day = step["time"][:10]
|
||
data = step.get("data", {})
|
||
inst = data.get("instant", {}).get("details", {})
|
||
a = agg.setdefault(day, {"t": [], "rh": [], "wind": [], "precip": None})
|
||
if (t := inst.get("air_temperature")) is not None:
|
||
a["t"].append(t)
|
||
if (rh := inst.get("relative_humidity")) is not None:
|
||
a["rh"].append(rh)
|
||
if (w := inst.get("wind_speed")) is not None:
|
||
a["wind"].append(w)
|
||
# Precip: prefer the 1-hour block (hourly near-term), else the 6-hour block
|
||
# (6-hourly far-term). Never add both — they overlap — so summing each step's
|
||
# chosen block avoids double counting across the resolution switch.
|
||
p = (data.get("next_1_hours") or {}).get("details", {}).get("precipitation_amount")
|
||
if p is None:
|
||
p = (data.get("next_6_hours") or {}).get("details", {}).get("precipitation_amount")
|
||
if p is not None:
|
||
a["precip"] = (a["precip"] or 0.0) + p
|
||
|
||
dates = sorted(agg)
|
||
c2f = lambda c: c * 9.0 / 5.0 + 32.0
|
||
df = pl.DataFrame({
|
||
"date": dates,
|
||
"tmax": [c2f(max(agg[d]["t"])) if agg[d]["t"] else None for d in dates],
|
||
"tmin": [c2f(min(agg[d]["t"])) if agg[d]["t"] else None for d in dates],
|
||
"precip": [agg[d]["precip"] / 25.4 if agg[d]["precip"] is not None else None
|
||
for d in dates],
|
||
"wind": [max(agg[d]["wind"]) * 2.2369362920544 if agg[d]["wind"] else None
|
||
for d in dates],
|
||
"gust": [None] * len(dates), # MET Norway has no gusts
|
||
"humid": [sum(agg[d]["rh"]) / len(agg[d]["rh"]) if agg[d]["rh"] else None
|
||
for d in dates],
|
||
})
|
||
df = df.with_columns(
|
||
pl.col("date").str.to_date(),
|
||
*[pl.col(c).cast(pl.Float64, strict=False)
|
||
for c in ("tmax", "tmin", "precip", "wind", "gust")],
|
||
)
|
||
df = df.with_columns(
|
||
pl.Series("fmax", _heat_index(df["tmax"].to_numpy(), df["humid"].to_numpy())),
|
||
pl.Series("fmin", _wind_chill(df["tmin"].to_numpy(), df["wind"].to_numpy())),
|
||
)
|
||
return _finalize_frame(df)
|
||
|
||
|
||
def _fetch_forecast_metno(cell: dict) -> pl.DataFrame:
|
||
"""Backup forward-forecast fetch from MET Norway / yr.no (used when Open-Meteo's
|
||
forecast API is unavailable). Coordinates are rounded to 4 decimals per MET's
|
||
ToS (improves their cache hit rate)."""
|
||
r = _request(
|
||
METNO_URL,
|
||
{"lat": round(cell["center_lat"], 4), "lon": round(cell["center_lon"], 4)},
|
||
60,
|
||
phase="forecast_metno",
|
||
headers={"User-Agent": METNO_UA},
|
||
)
|
||
return _metno_to_frame(r.json().get("properties", {}))
|
||
|
||
|
||
def _fetch_history_range(cell: dict, start_date: str, end_date: str) -> pl.DataFrame:
|
||
"""Fetch just a date range of archive history (used to top up the recent tail)."""
|
||
params = _om_daily_params(cell, start_date=start_date, end_date=end_date, models=ARCHIVE_MODEL)
|
||
r = _request(ARCHIVE_URL, params, 60, phase="history_topup")
|
||
return _to_frame(r.json()["daily"])
|
||
|
||
|
||
def _read_history_cache(path):
|
||
"""Schema-complete cached frame (no doy) + its file age in seconds, or None.
|
||
|
||
No age expiry — history is cached indefinitely; a pre-wind/humidity schema is
|
||
the only reason to refetch (to add the new metric columns)."""
|
||
if not os.path.exists(path):
|
||
return None
|
||
df = _normalize_read(pl.read_parquet(path))
|
||
if not all(c in df.columns for c in NEW_COLS):
|
||
return None
|
||
return df, time.time() - os.path.getmtime(path)
|
||
|
||
|
||
def _topup_tail(cell: dict, df: pl.DataFrame) -> pl.DataFrame:
|
||
"""Append newly-available archive days to a cached record. Best-effort and
|
||
serialized per cell; a small incremental fetch, not the full multi-decade pull."""
|
||
global _archive_cooldown_until
|
||
cell_id = cell["id"]
|
||
with _cell_lock(cell_id):
|
||
fresh = _read_history_backed(cell_id) # another thread may have just refreshed it
|
||
if fresh is not None:
|
||
df, age_s = fresh
|
||
if age_s < HISTORY_TOPUP_INTERVAL:
|
||
return df
|
||
if time.time() < _archive_cooldown_until:
|
||
return df
|
||
expected = datetime.date.today() - datetime.timedelta(days=ARCHIVE_LATENCY_DAYS)
|
||
cached_max = df["date"].max()
|
||
if cached_max >= expected:
|
||
_touch_history_backed(cell_id) # tail already current — reset the hourly timer
|
||
return df
|
||
try:
|
||
recent = _fetch_history_range(
|
||
cell, (cached_max + datetime.timedelta(days=1)).isoformat(), expected.isoformat())
|
||
except Exception as e: # noqa: BLE001 - keep the cached record on failure
|
||
if is_rate_limit(e):
|
||
_note_rate_limit(e)
|
||
return df
|
||
# Concatenate archive-first, recent-last so `keep="last"` prefers a freshly
|
||
# fetched day over its cached duplicate; maintain_order keeps that precedence
|
||
# before the final chronological sort.
|
||
tail = recent.drop("doy", strict=False)
|
||
merged = (pl.concat([df, tail], how="diagonal_relaxed")
|
||
.unique(subset="date", keep="last", maintain_order=True)
|
||
.sort("date"))
|
||
# Postgres only needs the new/refetched tail rows upserted (ON CONFLICT gives
|
||
# the same keep="last" precedence); parquet rewrites the whole file.
|
||
_write_history_backed(cell_id, merged, delta_df=tail)
|
||
return merged
|
||
|
||
|
||
def get_history(cell: dict) -> tuple[pl.DataFrame, dict]:
|
||
"""Return (daily history frame, cache metadata) for a cell, with humidity as
|
||
absolute humidity (g/m³). Thin wrapper over the raw loader (see below)."""
|
||
df, meta = _load_history(cell)
|
||
return _derive_metrics(df), meta
|
||
|
||
|
||
def load_cached_history(cell: dict) -> pl.DataFrame | None:
|
||
"""History from the cache ONLY — never fetches upstream and never tops up the
|
||
tail. Powers the warm-only prefetch path (which must not spend upstream quota)
|
||
and the offline migrate script. None when the cell has no (schema-complete)
|
||
cached record."""
|
||
hit = _read_history_backed(cell["id"])
|
||
if hit is None:
|
||
return None
|
||
return _derive_metrics(_with_doy(hit[0]))
|
||
|
||
|
||
def _load_history(cell: dict) -> tuple[pl.DataFrame, dict]:
|
||
"""Return (daily history frame, cache metadata) for a cell.
|
||
|
||
The full archive is cached indefinitely (fetched once); only the recent tail is
|
||
topped up, at most hourly. Concurrent callers are serialized so only one archive
|
||
fetch happens; on an upstream failure (e.g. 429) any existing cache is served."""
|
||
cell_id = cell["id"]
|
||
|
||
hit = _read_history_backed(cell_id)
|
||
if hit is not None:
|
||
df, age_s = hit
|
||
if age_s > HISTORY_TOPUP_INTERVAL:
|
||
df = _topup_tail(cell, df) # refresh just the recent days
|
||
return _with_doy(df), {"cached": True, "cache_age_days": round(age_s / 86400.0, 1)}
|
||
|
||
global _archive_cooldown_until
|
||
|
||
with _cell_lock(cell_id):
|
||
hit = _read_history_backed(cell_id) # another thread may have populated it while we waited
|
||
if hit is not None:
|
||
df, age_s = hit
|
||
return _with_doy(df), {"cached": True, "cache_age_days": round(age_s / 86400.0, 1)}
|
||
|
||
stale = _stale_history_backed(cell_id)
|
||
|
||
def _serve_stale():
|
||
return _with_doy(stale), {"cached": True, "stale": True, "cache_age_days": None}
|
||
|
||
# Fetch. Open-Meteo is primary (richer: gusts + apparent temp); NASA POWER is
|
||
# the backup when Open-Meteo is unavailable (network error or its daily limit).
|
||
# Skip Open-Meteo entirely while it's in a rate-limit cooldown.
|
||
df = None
|
||
source = None
|
||
if time.time() >= _archive_cooldown_until:
|
||
try:
|
||
fetched = _fetch_history(cell)
|
||
# Guard a self-hosted archive that isn't backfilled yet: an all-null or
|
||
# recent-only response reduces to a near-empty frame, which would
|
||
# otherwise be cached indefinitely as a complete record and never
|
||
# fall to NASA. Require a plausibly-full span before accepting it.
|
||
if fetched.height >= MIN_ARCHIVE_DAYS:
|
||
df = fetched
|
||
source = "open-meteo"
|
||
except Exception as e: # noqa: BLE001
|
||
if is_rate_limit(e):
|
||
_note_rate_limit(e)
|
||
if df is None:
|
||
try:
|
||
df = _fetch_history_nasa(cell)
|
||
source = "nasa-power"
|
||
except Exception: # noqa: BLE001 - backup unavailable too
|
||
df = None
|
||
if df is None:
|
||
if stale is not None:
|
||
return _serve_stale()
|
||
raise WeatherUnavailable(limit_message(_archive_limit_daily),
|
||
daily=_archive_limit_daily)
|
||
|
||
# Store the raw record; percentiles are derived at request time.
|
||
_write_history_backed(cell_id, df)
|
||
return df, {"cached": False, "cache_age_days": 0, "source": source}
|
||
|
||
|
||
RECENT_PAST_DAYS = 25 # recent observations window (covers the ~2-week graded view)
|
||
FORECAST_DAYS = 8 # today + 7 days ahead
|
||
|
||
|
||
def _rf_cache_path(cell_id: str) -> str:
|
||
return os.path.join(CACHE_DIR, f"{cell_id}_rf.parquet")
|
||
|
||
|
||
def recent_stamp(cell_id: str) -> int:
|
||
"""Identity stamp (whole seconds) of the cell's cached recent+forecast bundle;
|
||
0 when absent. On Postgres it's int(recent_synced_at) from climate_sync; on
|
||
parquet it's the file mtime. Either way it changes exactly when the
|
||
recent/forecast data is rewritten (and NOT on a stale serve), so it's the
|
||
freshness token for every payload that grades recent or future days."""
|
||
if climate_store.is_postgres():
|
||
return int(climate_store.recent_synced_at(cell_id))
|
||
try:
|
||
return int(os.path.getmtime(_rf_cache_path(cell_id)))
|
||
except OSError:
|
||
return 0
|
||
|
||
|
||
def get_recent_forecast(cell: dict) -> pl.DataFrame:
|
||
"""Recent observations + forward forecast, with humidity as absolute humidity
|
||
(g/m³). Thin wrapper over the raw loader (see below)."""
|
||
return _derive_metrics(_load_recent_forecast(cell))
|
||
|
||
|
||
def load_cached_recent_forecast(cell: dict) -> "pl.DataFrame | None":
|
||
"""Recent+forecast from the cache ONLY — never fetches upstream, and unlike the
|
||
loader below it does not care how stale the record is. The sibling of
|
||
load_cached_history, for the same reason: the homepage precompute sweeps every
|
||
cached city and must not spend a single upstream request doing it. None when
|
||
the cell has no cached recent/forecast record."""
|
||
hit = _read_recent_backed(cell["id"])
|
||
if hit is None:
|
||
return None
|
||
try:
|
||
return _derive_metrics(_with_doy(hit[0]))
|
||
except Exception: # noqa: BLE001 - a truncated/corrupt cache record is a miss, not a crash
|
||
return None
|
||
|
||
|
||
def _load_recent_forecast(cell: dict) -> pl.DataFrame:
|
||
"""Recent observations AND the forward forecast in ONE forecast-API call.
|
||
|
||
Both the recent (past) view and the forecast (future) view slice from this, so
|
||
a cell needs just one upstream forecast request per hour (plus the ~monthly
|
||
archive fetch). Cached per cell for one hour so it still follows model updates.
|
||
"""
|
||
cell_id = cell["id"]
|
||
hit = _read_recent_backed(cell_id)
|
||
if hit is not None:
|
||
cached, age_s = hit
|
||
if age_s / 3600.0 < FORECAST_TTL_HOURS:
|
||
return _with_doy(cached)
|
||
|
||
params = {
|
||
"latitude": cell["center_lat"],
|
||
"longitude": cell["center_lon"],
|
||
"daily": DAILY_VARS,
|
||
"timezone": "auto",
|
||
"temperature_unit": "fahrenheit",
|
||
"precipitation_unit": "inch",
|
||
"wind_speed_unit": "mph",
|
||
"past_days": RECENT_PAST_DAYS,
|
||
"forecast_days": FORECAST_DAYS,
|
||
}
|
||
try:
|
||
r = _request(FORECAST_URL, params, 60, phase="recent_forecast_fetch")
|
||
df = _to_frame(r.json()["daily"])
|
||
except Exception as e: # noqa: BLE001
|
||
# Open-Meteo forecast unavailable (rate limit or outage). Fall back to MET
|
||
# Norway (yr.no) — global + keyless — mirroring the archive's NASA POWER
|
||
# backup. MET Norway is forecast-only (no recent past days), so this keeps
|
||
# the forecast / day-ahead views working in a degraded form.
|
||
try:
|
||
df = _fetch_forecast_metno(cell)
|
||
except Exception: # noqa: BLE001 - backup unavailable too
|
||
# Last resort: serve the stale cache if we have one. An hours-old bundle
|
||
# (which still carries the recent observed days MET Norway lacks) beats a
|
||
# hard failure, mirroring _load_history's stale-serve. Return WITHOUT
|
||
# rewriting it, so recent_stamp stays old and the next request still
|
||
# retries upstream first rather than serving this as if it were fresh.
|
||
if hit is not None:
|
||
return _with_doy(hit[0])
|
||
# No cache either: surface the original error, classifying an Open-Meteo
|
||
# rate limit as the typed, daily-aware WeatherUnavailable (the archive
|
||
# path classifies its own inside _load_history).
|
||
if is_rate_limit(e):
|
||
daily = "daily" in _rate_limit_reason(e).lower()
|
||
raise WeatherUnavailable(limit_message(daily), daily=daily) from e
|
||
raise
|
||
_write_recent_backed(cell_id, df)
|
||
return df
|
||
|
||
|
||
_REVGEO_CACHE: dict[str, str | None] = {}
|
||
|
||
# Nominatim's usage policy allows ~1 request/second and rejects bursts. The compare
|
||
# page loads several locations at once, so serialize the reverse lookups here and
|
||
# space them out — otherwise the burst is rate-limited, a null label gets cached, and
|
||
# those locations are left showing bare coordinates.
|
||
_REVGEO_LOCK = threading.Lock()
|
||
_REVGEO_MIN_INTERVAL = 1.1 # seconds between successive Nominatim reverse calls
|
||
_revgeo_last = 0.0
|
||
|
||
|
||
def reverse_geocode_cached(lat: float, lon: float) -> tuple[bool, str | None]:
|
||
"""(found, label) from the in-memory + SQLite revgeo caches only — never calls
|
||
Nominatim. Used by paths that must not add upstream traffic (prefetch)."""
|
||
key = store.revgeo_key(lat, lon)
|
||
if key in _REVGEO_CACHE:
|
||
return (True, _REVGEO_CACHE[key])
|
||
found, label = store.get_revgeo(key)
|
||
if found:
|
||
_REVGEO_CACHE[key] = label
|
||
return (found, label)
|
||
|
||
|
||
def reverse_geocode(lat: float, lon: float) -> str | None:
|
||
"""Best-effort neighbourhood/city label for a point (OpenStreetMap Nominatim).
|
||
|
||
Cached per ~cell — in memory for the process and in SQLite across restarts —
|
||
so panning around (or redeploying) doesn't re-hit the service, and failures
|
||
return None so the caller can fall back to bare coordinates.
|
||
"""
|
||
key = store.revgeo_key(lat, lon)
|
||
found, label = reverse_geocode_cached(lat, lon)
|
||
if found:
|
||
return label
|
||
# Serialize + rate-limit the upstream call. Concurrent callers (the compare page
|
||
# loads several locations at once) would otherwise burst past Nominatim's ~1/sec
|
||
# limit and get rate-limited, caching a null label. Under the lock we re-check the
|
||
# cache (a peer may have just resolved this cell) and space successive calls out.
|
||
global _revgeo_last
|
||
with _REVGEO_LOCK:
|
||
found, label = reverse_geocode_cached(lat, lon)
|
||
if found:
|
||
return label
|
||
wait = _REVGEO_MIN_INTERVAL - (time.monotonic() - _revgeo_last)
|
||
if wait > 0:
|
||
time.sleep(wait)
|
||
label = None
|
||
try:
|
||
r = _request(
|
||
"https://nominatim.openstreetmap.org/reverse",
|
||
# zoom 14 resolves to the suburb/neighbourhood level so we can lead
|
||
# with it when OSM has one (zoom 10 only ever returns the city).
|
||
# accept-language=en keeps labels in one script worldwide (matches
|
||
# the forward geocoder's language=en).
|
||
{"lat": lat, "lon": lon, "format": "jsonv2", "zoom": 14,
|
||
"addressdetails": 1, "accept-language": "en"},
|
||
15,
|
||
phase="reverse_geocode",
|
||
headers={"User-Agent": "Thermograph/0.1 (local weather grading app)"},
|
||
)
|
||
a = r.json().get("address", {}) or {}
|
||
# Lead with the neighbourhood when available, then the city, then region.
|
||
neighborhood = (a.get("neighbourhood") or a.get("suburb")
|
||
or a.get("quarter") or a.get("city_district")
|
||
or a.get("borough"))
|
||
city = (a.get("city") or a.get("town") or a.get("village")
|
||
or a.get("hamlet") or a.get("county"))
|
||
region = a.get("state") or a.get("province") or a.get("region")
|
||
country = a.get("country")
|
||
# dict.fromkeys drops duplicates (e.g. neighbourhood == city) in order.
|
||
# Country trails the neighbourhood/city/region so the label reads as a full
|
||
# hierarchy ("West Seattle, Seattle, Washington, United States"); the
|
||
# frontend leads with the first part and mutes the rest.
|
||
parts = dict.fromkeys(p for p in (neighborhood, city, region, country) if p)
|
||
label = ", ".join(parts) or None
|
||
except Exception: # noqa: BLE001 - reverse geocoding is a nicety, never fatal
|
||
label = None
|
||
_revgeo_last = time.monotonic()
|
||
_REVGEO_CACHE[key] = label
|
||
store.put_revgeo(key, label) # survive restarts (a None label retries after its TTL)
|
||
return label
|
||
|
||
|
||
def geocode(name: str, count: int = 5) -> list[dict]:
|
||
"""Look up places by name worldwide via Open-Meteo's geocoder."""
|
||
r = _request(
|
||
"https://geocoding-api.open-meteo.com/v1/search",
|
||
{"name": name, "count": count, "language": "en", "format": "json"},
|
||
30,
|
||
phase="geocode",
|
||
)
|
||
results = r.json().get("results", []) or []
|
||
return [
|
||
{
|
||
"name": g.get("name"),
|
||
"admin1": g.get("admin1"),
|
||
"country": g.get("country"),
|
||
"country_code": g.get("country_code"),
|
||
"lat": g.get("latitude"),
|
||
"lon": g.get("longitude"),
|
||
"population": g.get("population"),
|
||
}
|
||
for g in results
|
||
]
|