"""ClearCast Space: EXPERIMENTAL underwater visibility estimate for Finistère (satellite Secchi depth -> diver visibility)."""
import joblib, numpy as np, pandas as pd, requests, gradio as gr
import plotly.graph_objects as go
from features import daily_weather_from_hourly
from spots import SPOTS

MODELS = joblib.load("zsd_models.joblib")
LATEST = pd.read_csv("latest.csv").set_index("spot")
HIND = pd.read_parquet("hindcast.parquet")
NAMES = list(SPOTS)
K, LINK = 0.76, {"0.1": -0.731, "0.5": 0.008, "0.9": 0.597}
SD_LINK = (LINK["0.9"] - LINK["0.1"]) / 2.563
LABEL = ("EXPERIMENTAL estimate. Visibility ≈ 0.76 × satellite Secchi depth, a relation checked against ~6,000 diver surveys "
         "(Reef Life Survey) worldwide, NOT in Brittany. No validation against Finistère visibility observations exists yet.")


def live_weather(lat, lon):
    base = dict(latitude=lat, longitude=lon, past_days=10, forecast_days=1, timezone="UTC")
    a = requests.get("https://api.open-meteo.com/v1/forecast", params=dict(base, hourly="wind_speed_10m,wind_direction_10m,precipitation"), timeout=30).json()
    m = requests.get("https://marine-api.open-meteo.com/v1/marine", params=dict(base, hourly="wave_height,wave_period,wave_direction"), timeout=30).json()
    w = pd.DataFrame(a["hourly"]).merge(pd.DataFrame(m["hourly"]), on="time")
    w["time"] = pd.to_datetime(w["time"])
    return w[w["time"] <= pd.Timestamp.utcnow().tz_localize(None)]


def predict(spot, horizon, last_zsd):
    lat, lon, _ = SPOTS[spot]
    h = int(horizon)
    d = daily_weather_from_hourly(live_weather(lat, lon))
    today = d.iloc[-1:].copy()
    today["last_kd"] = np.log(float(last_zsd))
    today["last_age"] = 0
    doy = pd.Timestamp.utcnow().dayofyear
    today["doy_sin"], today["doy_cos"] = np.sin(2 * np.pi * doy / 365.25), np.cos(2 * np.pi * doy / 365.25)
    today["lat"], today["lon"] = lat, lon
    m = MODELS[h]
    X = today[m["features"]]
    mg = m["conformal_margin"]
    l10, l50, l90 = (float(m["q"][k].predict(X)[0]) for k in ("0.1", "0.5", "0.9"))
    l10, l90 = l10 - mg, l90 + mg
    sd = np.sqrt(((l90 - l10) / 2.563) ** 2 + SD_LINK ** 2)
    mu = l50 + np.log(K) + LINK["0.5"]
    v10, v50, v90 = (float(np.exp(mu + z * sd)) for z in (-1.2816, 0, 1.2816))
    c = d.iloc[-1][["hs_max", "tp_mean", "wind", "rain_r3"]]
    return (f"### {spot}: estimated visibility in +{h} day(s)\n**≈ {v50:.0f} m** (80% range {v10:.0f} – {v90:.0f} m)\n\n"
            f"Forecast satellite Secchi depth ≈ {np.exp(l50):.1f} m. Today: Hs max {c.hs_max:.1f} m, period {c.tp_mean:.1f} s, "
            f"wind {c.wind:.1f} km/h, rain 3-day mean {c.rain_r3:.1f} mm.\n\n"
            f"Starting Secchi depth used: {float(last_zsd):.1f} m. Snapshot: {LATEST.loc[spot,'last_zsd']:.1f} m from **{LATEST.loc[spot,'last_date']}** "
            f"(stale: enter a recent value if you have one).\n\n> {LABEL}")


def hindcast_plot(spot):
    d = HIND[HIND.spot == spot].sort_values("date")
    f = go.Figure()
    f.add_traces([go.Scatter(x=d.date, y=d.p90, line=dict(width=0), showlegend=False),
                  go.Scatter(x=d.date, y=d.p10, fill="tonexty", line=dict(width=0), name="80% range", fillcolor="rgba(70,130,180,0.25)"),
                  go.Scatter(x=d.date, y=d.p50, name="Estimated visibility (24 h ahead)", line=dict(color="steelblue")),
                  go.Scatter(x=d.date, y=d.zsd_equiv_vis, mode="markers", name="0.76 × satellite Secchi (not a diver measurement)", marker=dict(color="black", size=4))])
    f.update_layout(title=f"2026 test period: {spot}", yaxis_title="Visibility (m)", height=380, margin=dict(l=10, r=10, t=40, b=10))
    return f


with gr.Blocks(title="ClearCast") as demo:
    gr.Markdown("# ClearCast — Finistère underwater visibility (experimental)\n" + LABEL +
                "\n\nPipeline: LightGBM quantile forecast of 4 km satellite Secchi depth (Copernicus Marine) from waves, wind and rain (Open-Meteo), "
                "then a physical scaling to visibility with uncertainty from the diver-survey calibration. Skill vs persistence is small (see model card).")
    with gr.Row():
        spot = gr.Dropdown(NAMES, value="Brest Iroise", label="Spot")
        hz = gr.Radio([1, 2, 3], value=1, label="Horizon (days)")
        lk = gr.Number(value=round(float(LATEST.loc["Brest Iroise", "median_zsd"]), 1), label="Latest satellite Secchi depth (m)")
    btn = gr.Button("Estimate visibility", variant="primary")
    out = gr.Markdown()
    plot = gr.Plot()
    btn.click(predict, [spot, hz, lk], out)
    spot.change(lambda s: (hindcast_plot(s), round(float(LATEST.loc[s, "median_zsd"]), 1)), spot, [plot, lk])
    demo.load(hindcast_plot, spot, plot)

if __name__ == "__main__":
    demo.launch()