#!/opt/anaconda3/bin/python
"""pn-070/071 runner — see 2026-08-22-pn070-gex-regime.md (frozen design).
Rerunnable as the GEX archive accrues; judge at n >= 60 classified days."""
import gzip, json, glob, math
from collections import defaultdict
import pandas as pd

GEX = "<data-root>/fp-data/gex/2026"
TAPE = "<data-root>/fp-data/glbx/{sym}/1m/2026.csv"
MAP = {"NQ": "QQQ", "ES": "SPY"}

def f(x, d=0.0):
    try: return float(x)
    except (TypeError, ValueError): return d

def load_snapshots(ticker):
    """First snapshot per ET day with et_time <= 09:30, else earliest."""
    by_day = {}
    for path in sorted(glob.glob(f"{GEX}/{ticker}_*.json.gz")):
        try: d = json.load(gzip.open(path))
        except Exception: continue
        day, t = d.get("et_date"), d.get("et_time", "99:99")
        if not day: continue
        cur = by_day.get(day)
        if cur is None or (t <= "09:30" and not cur[0] <= "09:30") or (t <= "09:30" and t > cur[0]) or (cur[0] > "09:30" and t < cur[0]):
            by_day[day] = (t, d)
    return {day: d for day, (t, d) in by_day.items()}

def classify(d):
    """House computeGamma: tilt + nearest-significant call wall above spot."""
    spot = f(d["data"].get("stock_state", {}).get("close"))
    rows = d["data"].get("greek_exposure_strike") or []
    if not spot or not rows: return None
    sks = []
    for r in rows:
        k = f(r.get("strike"))
        if k <= 0: continue
        cg, pg = f(r.get("call_gex")), f(r.get("put_gex"))
        sks.append((k, cg, pg, cg + pg))
    if not sks: return None
    gross = sum(abs(c) + abs(p) for _, c, p, _ in sks)
    tilt = (sum(n for *_, n in sks) / gross) if gross > 0 else 0.0
    regime = "pos" if tilt > 0.08 else "neg" if tilt < -0.08 else None
    lo, hi = spot * 0.94, spot * 1.06
    win = [s for s in sks if lo <= s[0] <= hi] or sks
    above = [s for s in win if s[0] >= spot * 1.003]
    wall = None
    if above:
        peak = max(c for _, c, _, _ in above)
        sig = [s for s in above if s[1] >= 0.5 * peak] if peak > 0 else []
        pool = sig or above
        wall = min(pool, key=lambda s: abs(s[0] - spot))[0]
    return {"spot": spot, "tilt": tilt, "regime": regime, "wall": wall}

def load_tape(sym):
    df = pd.read_csv(TAPE.format(sym=sym), usecols=["ts_event", "open", "high", "low", "close"])
    ts = pd.to_datetime(df["ts_event"], utc=True).dt.tz_convert("America/New_York")
    df["day"] = ts.dt.strftime("%Y-%m-%d"); df["hm"] = ts.dt.strftime("%H:%M")
    return df[(df["hm"] >= "09:30") & (df["hm"] <= "15:59")]

def grid5(hm):  # 9:30-anchored 5m bar END minute is :34,:39,...
    h, m = int(hm[:2]), int(hm[3:])
    return (m % 5) == 4

def run_day(day_df, snap, fut_ratio_min="09:30"):
    r = {}
    bars = day_df.sort_values("hm")
    if len(bars) < 100: return None
    orr = bars[(bars["hm"] >= "09:30") & (bars["hm"] <= "09:59")]
    if orr.empty: return None
    or_hi, or_lo = orr["high"].max(), orr["low"].min()
    last = bars[bars["hm"] <= "15:55"].iloc[-1]["close"]
    # T1 breakout
    t1 = None
    for _, b in bars[(bars["hm"] >= "10:00") & (bars["hm"] <= "11:59")].iterrows():
        if not grid5(b["hm"]): continue
        if b["close"] > or_hi: t1 = last - b["close"]; break
        if b["close"] < or_lo: t1 = b["close"] - last; break
    r["t1"] = t1
    # T2 wall fade (short at scaled wall touch)
    t2 = None
    if snap.get("wall"):
        anchor = bars.iloc[0]["open"]
        wall_fut = snap["wall"] * (anchor / snap["spot"])
        touched = bars[(bars["hm"] >= "10:00") & (bars["high"] >= wall_fut)]
        if not touched.empty: t2 = wall_fut - last
    r["t2"] = t2
    return r

def welch(a, b):
    if len(a) < 2 or len(b) < 2: return float("nan")
    ma, mb = sum(a)/len(a), sum(b)/len(b)
    va = sum((x-ma)**2 for x in a)/(len(a)-1); vb = sum((x-mb)**2 for x in b)/(len(b)-1)
    se = math.sqrt(va/len(a) + vb/len(b))
    return (ma - mb) / se if se > 0 else float("nan")

out = defaultdict(list)
for sym, tick in MAP.items():
    snaps = load_snapshots(tick)
    tape = load_tape(sym)
    for day, sd in tape.groupby("day"):
        s = snaps.get(day)
        cls = classify(s) if s else None
        if not cls or not cls["regime"]: continue
        res = run_day(sd, cls)
        if not res: continue
        for test in ("t1", "t2"):
            if res[test] is not None:
                out[(sym, cls["regime"], test)].append(res[test])

print("pn-070/071 GEX regime — first read (archive begins 2026-07-26; HARNESS, not verdict)")
for (sym, reg, test), xs in sorted(out.items()):
    mu = sum(xs)/len(xs); wr = 100*sum(1 for x in xs if x > 0)/len(xs)
    print(f"  {sym} {reg:>3}GEX {test}: n={len(xs):2d}  mean {mu:+7.1f}pt  win {wr:3.0f}%")
for sym in MAP:
    a, b = out.get((sym, "neg", "t1"), []), out.get((sym, "pos", "t1"), [])
    print(f"H1 {sym} (breakout negGEX vs posGEX): t={welch(a,b):+.2f} (n {len(a)}/{len(b)})")
    c, d2 = out.get((sym, "pos", "t2"), []), out.get((sym, "neg", "t2"), [])
    print(f"H2 {sym} (fade posGEX vs negGEX):     t={welch(c,d2):+.2f} (n {len(c)}/{len(d2)})")