"""FILONE 1 — LA CURVA FRONTIERA (r0703_frontier_curve.py)

MISURA (non ricerca segnali) della frontiera rischio/rendimento dei DUE book che GIA' esistono:
  - book RESEARCH 5-sleeve: TP01 33 / XS01 15 / VRP01 12 / SKH01 20 / GTAA01 20 (combine_outer)
  - book LIVE 2-sleeve:      TP01 75 / SKH01 25   (cio' che gira su soldi reali su Deribit)

Per decidere a quale TARGET-VOL conviene girare il book a capitale piccolo (2k / 5k).

LEVA STUDIATA = moltiplicatore lambda sui rendimenti giornalieri del book combinato:
    r_scaled(t) = clip(lambda * r_book(t), -1, None)     (clip: non si perde >100% in un giorno)
Sweep di lambda per coprire target-vol annuo ~5% -> ~30%.  Sharpe INVARIANTE per costruzione
(lo verifichiamo: cambia solo con l'attivazione del clip a leva alta).

CONFIDENZA PER-SLEEVE (anti-alllucinazione #2): 3 livelli. La de-luck e l'haircut sono un
HAIRCUT SULLA MEDIA (edge), applicato DE-MEANando lo sleeve: s' = s - h*mean(s).
=> std invariata (stessa vol/coda), drift ridotto. Quindi tra i 3 livelli cambia SOLO il drift
   (CAGR / P-rovina via minor deriva), NON la forma della coda -> il messaggio e' la BANDA.
  (a) FULL      : book as-is.
  (b) DE-LUCK   : rimuove l'anchor timing-luck (audit 4/4):
                    SKH01 mean x0.60 (canonica = 93-98 pctl di 23 offset, ~+0.5 HOLD di fortuna),
                    TP01  mean x0.85 (hold-out 0.31 = migliore di 24 ancore, mediana 0.04).
                  VRP01/XS01 NON haircuttati qui: le loro ancore canoniche sono CONSERVATIVE
                    (VRP 7 pctl = peggiore; XS 15 pctl di DD).  GTAA: nessuna ancora-luck.
  (c) HAIRCUT   : de-luck + taglio -50% della MEDIA degli sleeve a FIDUCIA BASSA sull'edge
                    (VRP01 premio MODELLATO, f di stress mai catturato; XS01 STAT-MODE, 2.5y).
                  = la "frontiera fidata". Il GAP full<->haircut e' il margine di onesta'.
  Variante extra: EXCL = escludi del tutto VRP01+XS01 (ultra-conservativo, riferimento).
  Per il book LIVE 2-sleeve (niente VRP/XS): haircut = de-luck + SKH01 mean a x0.50 dall'originale
    (SKH e' la gamba a fiducia MEDIA: research, ETH DD margine sottile vs 30%).

TAIL (anti-allucinazione #1): la coda NON scala lineare-Gaussiano. Quantificato con:
  - BLOCK BOOTSTRAP (blocchi 10 e 20g, B=5000) sui rendimenti REALI del book (fat-tail+autocorr):
      P(rovina-50%|5y) = P(equity tocca <=50% del capitale iniziale in un path di 5y)
      P(DD>30%|5y)     = P(max drawdown peak-to-trough > 30% in 5y)
  - BENCHMARK GAUSSIANO iid (stessa media/vol) -> il GAP misura il contributo fat-tail/autocorr.
  - INIEZIONE DI CODA SINTETICA: un giorno-crash = 1.5x il peggior giorno storico, ~1 volta per
      path di 5y -> sensibilita' al "il campione puo' non contenere il crash peggiore".

Nessuna selezione su hold-out: FULL e HOLD affiancati solo per PROVARE che la frontiera non e'
un artefatto in-sample.  NON tocca src/config/live/tests: importa solo (sola lettura) i builder.

Output CSV+report testuale in scratchpad.  Girare:  uv run python scripts/research/r0703_frontier_curve.py
"""
from __future__ import annotations
import sys, os, time
sys.path.insert(0, '/opt/docker/PythagorasGoal')
sys.path.insert(0, '/opt/docker/PythagorasGoal/scripts/research/alt')
import numpy as np
import pandas as pd

from src.portfolio.sleeves import (_tp01_returns, _skyhook_returns, _xsec_returns,
                                    _vrp_combo_returns, _gtaa_daily_returns)
from src.portfolio.portfolio import combine_outer, to_daily, HOLDOUT

DPY = 365.25
OUT = '/tmp/claude-1001/-opt-docker-PythagorasGoal/e00896d3-d4bb-4f2a-b471-55a1d88a12ba/scratchpad'
os.makedirs(OUT, exist_ok=True)
RNG = np.random.default_rng(20260703)

# ---------------------------------------------------------------- sleeves
def build_sleeves() -> dict[str, pd.Series]:
    d = {}
    for nm, fn in [("TP01", _tp01_returns), ("XS01", _xsec_returns), ("VRP01", _vrp_combo_returns),
                   ("SKH01", _skyhook_returns), ("GTAA01", _gtaa_daily_returns)]:
        d[nm] = to_daily(fn())
    return d


def demean_haircut(s: pd.Series, h: float) -> pd.Series:
    """s' = s - h*mean(s): riduce il drift del fattore (1-h), std INVARIATA."""
    if h == 0.0:
        return s
    return s - h * float(s.mean())


# per-livello: fattore di haircut sulla MEDIA per sleeve  (h)
LEVELS_5 = {
    "FULL":    dict(TP01=0.0,  XS01=0.0, VRP01=0.0, SKH01=0.0,  GTAA01=0.0),
    "DELUCK":  dict(TP01=0.15, XS01=0.0, VRP01=0.0, SKH01=0.40, GTAA01=0.0),
    "HAIRCUT": dict(TP01=0.15, XS01=0.5, VRP01=0.5, SKH01=0.40, GTAA01=0.0),
    "EXCL":    dict(TP01=0.15, XS01=1.0, VRP01=1.0, SKH01=0.40, GTAA01=0.0),  # escludi VRP+XS (mean->0)
}
W5 = {"TP01": 0.33, "XS01": 0.15, "VRP01": 0.12, "SKH01": 0.20, "GTAA01": 0.20}

LEVELS_2 = {
    "FULL":    dict(TP01=0.0,  SKH01=0.0),
    "DELUCK":  dict(TP01=0.15, SKH01=0.40),
    "HAIRCUT": dict(TP01=0.15, SKH01=0.50),   # SKH = gamba fiducia MEDIA -> taglio piu' profondo
}
W2 = {"TP01": 0.75, "SKH01": 0.25}


def make_book(sleeves: dict, weights: dict, hair: dict) -> pd.Series:
    cols = {nm: demean_haircut(sleeves[nm], hair.get(nm, 0.0)) for nm in weights}
    return combine_outer(cols, weights)


# ---------------------------------------------------------------- metrics (full sample, scaled)
def full_metrics(book: pd.Series, lam: float) -> dict:
    r = np.clip(lam * book.values.astype(float), -1.0, None)
    s = pd.Series(r, index=book.index)
    eq = np.cumprod(1.0 + r)
    pk = np.maximum.accumulate(eq)
    yrs = len(r) / DPY
    dd = float(np.max((pk - eq) / pk))
    cagr = float(eq[-1] ** (1 / yrs) - 1) if eq[-1] > 0 else -1.0
    shp = float(r.mean() / r.std() * np.sqrt(DPY)) if r.std() > 0 else 0.0
    ww = float(((1.0 + s).resample('168h').prod() - 1.0).min())   # worst-week (168h, non '7D')
    vol = float(r.std() * np.sqrt(DPY))
    return dict(sharpe=shp, cagr=cagr, maxdd=dd, worstweek=ww, vol=vol)


# ---------------------------------------------------------------- block bootstrap engine
def block_index(n: int, L: int, npath: int, B: int) -> np.ndarray:
    """Matrice (B, npath) di indici in [0,n) via block-bootstrap circolare-troncato (blocchi L)."""
    nblk = int(np.ceil(npath / L))
    starts = RNG.integers(0, n - L + 1, size=(B, nblk))       # start di ogni blocco
    offs = np.arange(L)
    idx = (starts[:, :, None] + offs[None, None, :]).reshape(B, nblk * L)[:, :npath]
    return idx


def ruin_dd(Rmat: np.ndarray, lam: float) -> tuple[float, float]:
    """P(rovina-50%|5y), P(DD>30%|5y) da una matrice (B,npath) di rendimenti NATIVI, scalata per lam.
    Versione veloce: la matrice bootstrap (o la coda iniettata, o il gauss iid) e' PRE-costruita una
    volta per (book,level); qui si applica solo lo scaling di leva + cumprod. Numerica identica a
    boot_probs/gauss_probs, senza ri-gather/ri-generazione per ogni target-vol."""
    R = np.clip(lam * Rmat, -1.0, None)
    eq = np.cumprod(1.0 + R, axis=1)
    ruin = float((eq.min(axis=1) <= 0.5).mean())
    pk = np.maximum.accumulate(eq, axis=1)
    dd = ((pk - eq) / pk).max(axis=1)
    return ruin, float((dd > 0.30).mean())


def boot_probs(r_native: np.ndarray, lam: float, idx: np.ndarray,
               inject_crash: float | None = None) -> tuple[float, float]:
    """P(rovina-50%|5y), P(DD>30%|5y) su path bootstrap.
    inject_crash: se non None, in OGNI path sostituisce 1 giorno casuale con questo rendimento
                  NATIVO (pre-scala) -> stress "una coda sintetica per 5 anni"."""
    R = r_native[idx].copy()                                  # (B, npath) nativi
    if inject_crash is not None:
        B, npath = R.shape
        pos = RNG.integers(0, npath, size=B)
        R[np.arange(B), pos] = inject_crash
    R = np.clip(lam * R, -1.0, None)
    eq = np.cumprod(1.0 + R, axis=1)
    ruin = float((eq.min(axis=1) <= 0.5).mean())              # tocca 50% del capitale iniziale
    pk = np.maximum.accumulate(eq, axis=1)
    dd = ((pk - eq) / pk).max(axis=1)
    dd30 = float((dd > 0.30).mean())
    return ruin, dd30


def gauss_probs(mu: float, sd: float, npath: int, B: int) -> tuple[float, float]:
    """Benchmark iid-Gaussiano con stessa media/std giornaliera (post-scala)."""
    R = RNG.normal(mu, sd, size=(B, npath))
    R = np.clip(R, -1.0, None)
    eq = np.cumprod(1.0 + R, axis=1)
    ruin = float((eq.min(axis=1) <= 0.5).mean())
    pk = np.maximum.accumulate(eq, axis=1)
    dd = ((pk - eq) / pk).max(axis=1)
    return ruin, float((dd > 0.30).mean())


# ---------------------------------------------------------------- main
def run():
    t0 = time.time()
    print("building sleeves...", flush=True)
    sl = build_sleeves()

    NPATH = int(round(5 * DPY))     # 5 anni ~ 1826 giorni-calendario
    B = 5000
    TV_GRID = np.round(np.arange(0.05, 0.305, 0.01), 3)    # 5%..30% target-vol (fine)
    TV_FINE = np.round(np.arange(0.04, 0.305, 0.005), 4)   # per la ricerca del punto-decisione

    books = {
        "5SLEEVE": dict(weights=W5, levels=LEVELS_5),
        "2LIVE":   dict(weights=W2, levels=LEVELS_2),
    }

    # ---- costruzione book per livello + sanity + banda di confidenza
    print("\n" + "=" * 78)
    print("SANITY + BANDA DI CONFIDENZA (FULL vs HOLD affiancati)")
    print("=" * 78)
    built = {}   # (book,level) -> pd.Series (native, lam=1)
    for bk, cfg in books.items():
        for lv, hair in cfg["levels"].items():
            b = make_book(sl, cfg["weights"], hair)
            built[(bk, lv)] = b
        # print sanity for this book, all levels
        print(f"\n[{bk}]  weights={cfg['weights']}")
        print(f"  {'level':9s} {'nat.vol':>8s} {'FULL Sh':>8s} {'HOLD Sh':>8s} "
              f"{'FULL CAGR':>10s} {'HOLD CAGR':>10s} {'FULL DD':>8s}")
        for lv in cfg["levels"]:
            b = built[(bk, lv)]
            mF = full_metrics(b, 1.0)
            bh = b[b.index >= HOLDOUT]
            mH = full_metrics(bh, 1.0)
            print(f"  {lv:9s} {mF['vol']*100:7.2f}% {mF['sharpe']:8.3f} {mH['sharpe']:8.3f} "
                  f"{mF['cagr']*100:9.2f}% {mH['cagr']*100:9.2f}% {mF['maxdd']*100:7.2f}%")

    # ---- anchor band (citata dagli audit; non ricalcolata qui)
    print("\n" + "-" * 78)
    print("BANDA D'ANCORA (dove il DE-LUCK e' ancorato — audit 4/4, diari 2026-07-02/03):")
    print("  TP01 : hold-out 0.31 = MIGLIORE di 24 ancore orarie (mediana 0.04, banda [-0.13,+0.30],")
    print("         P~0.86 che un'ancora mostri tale spike per caso) -> de-luck mean x0.85")
    print("  SKH01: canonica = 93-98 pctl di 23 offset (griglia 230/690m), ~+0.5 HOLD di fortuna sul")
    print("         book; gate DD<30% fallisce in 15/23 offset -> de-luck mean x0.60")
    print("  XS01 : canonica = 15 pctl di DD (10 fasi) -> CONSERVATIVA sul DD, nessun de-luck edge")
    print("  VRP01: canonica = PEGGIORE (7 pctl, 7 giorni) -> CONSERVATIVA, nessun de-luck edge")
    print("  GTAA : validato OOS 30y indipendente, nessuna ancora-luck")
    print("  => il DE-LUCK haircutta SOLO gli sleeve con ancora-FORTUNATA (TP01/SKH01).")

    # ---- frontiera + bootstrap
    rows = []          # per CSV
    print("\n" + "=" * 78)
    print("FRONTIERA + BLOCK-BOOTSTRAP  (B=%d, 5y=%d gg, blocchi 10 & 20)" % (B, NPATH))
    print("=" * 78)

    # crash sintetico per book: 1.5x il peggior giorno storico del book FULL (nativo)
    for bk, cfg in books.items():
        native_min = float(built[(bk, "FULL")].values.min())
        crash = 1.5 * native_min
        n = len(built[(bk, "FULL")])
        idx10 = block_index(n, 10, NPATH, B)
        idx20 = block_index(n, 20, NPATH, B)
        print(f"\n########## BOOK {bk}  (n={n} gg, worst-day nativo {native_min*100:.2f}%, "
              f"crash sintetico {crash*100:.2f}%) ##########")
        for lv in cfg["levels"]:
            b = built[(bk, lv)]
            r_nat = b.values.astype(float)
            bh = b[b.index >= HOLDOUT]
            r_hold = bh.values.astype(float)
            nat_vol = float(r_nat.std() * np.sqrt(DPY))
            # PRE-costruzione (una volta per level) delle matrici bootstrap NATIVE + coda + gauss:
            R20 = r_nat[idx20]                                        # (B, npath) block-20
            R10 = r_nat[idx10]                                        # (B, npath) block-10
            R20c = R20.copy()                                         # coda sintetica iniettata 1x/path
            _pos = RNG.integers(0, R20c.shape[1], size=R20c.shape[0])
            R20c[np.arange(R20c.shape[0]), _pos] = crash
            Z = RNG.standard_normal(R20.shape)                        # normali std fisse (gauss iid)
            mu0 = float(r_nat.mean()); sd0 = float(r_nat.std())
            print(f"\n--- {bk} / {lv}  (native vol {nat_vol*100:.2f}%, "
                  f"FULL Sh {full_metrics(b,1.0)['sharpe']:.3f}, "
                  f"HOLD Sh {full_metrics(bh,1.0)['sharpe']:.3f}) ---")
            hdr = (f"  {'tvol':>5s} {'lam':>5s} {'Sh':>6s} {'CAGRf':>7s} {'CAGRh':>7s} "
                   f"{'maxDDf':>7s} {'wWeek':>7s} | {'Pruin20':>8s} {'Pdd30_20':>9s} "
                   f"{'Pruin10':>8s} {'Pdd30_10':>9s} | {'Gruin':>7s} {'Gdd30':>7s} | "
                   f"{'Pruin+cr':>9s} {'Pdd30+cr':>9s}")
            print(hdr)
            for tv in TV_GRID:
                lam = tv / nat_vol
                mF = full_metrics(b, lam)
                mH = full_metrics(bh, lam)
                r20, d20 = ruin_dd(R20, lam)
                r10, d10 = ruin_dd(R10, lam)
                # gaussiano iid con media/std della serie SCALATA (Rg = lam*mu0 + lam*sd0*Z)
                gr, gd = ruin_dd((lam * mu0) + (lam * sd0) * Z, 1.0)
                # tail injection (blocco 20)
                rc, dc = ruin_dd(R20c, lam)
                print(f"  {tv*100:4.0f}% {lam:5.2f} {mF['sharpe']:6.2f} "
                      f"{mF['cagr']*100:6.1f}% {mH['cagr']*100:6.1f}% "
                      f"{mF['maxdd']*100:6.1f}% {mF['worstweek']*100:6.1f}% | "
                      f"{r20*100:7.2f}% {d20*100:8.2f}% {r10*100:7.2f}% {d10*100:8.2f}% | "
                      f"{gr*100:6.2f}% {gd*100:6.2f}% | {rc*100:8.2f}% {dc*100:8.2f}%")
                rows.append(dict(book=bk, level=lv, target_vol=tv, lam=round(lam, 3),
                                 sharpe=round(mF['sharpe'], 3), cagr_full=round(mF['cagr'], 4),
                                 cagr_hold=round(mH['cagr'], 4), maxdd_full=round(mF['maxdd'], 4),
                                 worstweek=round(mF['worstweek'], 4),
                                 p_ruin50_b20=round(r20, 4), p_dd30_b20=round(d20, 4),
                                 p_ruin50_b10=round(r10, 4), p_dd30_b10=round(d10, 4),
                                 p_ruin50_gauss=round(gr, 4), p_dd30_gauss=round(gd, 4),
                                 p_ruin50_crash=round(rc, 4), p_dd30_crash=round(dc, 4)))

    df = pd.DataFrame(rows)
    csv = os.path.join(OUT, 'frontier_curve.csv')
    df.to_csv(csv, index=False)
    print(f"\nCSV completo -> {csv}")

    # ---- punto OGGI (lam=1) marcato
    print("\n" + "=" * 78)
    print("PUNTO OPERATIVO OGGI (lambda=1, book as-run):")
    for bk, cfg in books.items():
        b = built[(bk, "FULL")]
        m = full_metrics(b, 1.0)
        print(f"  {bk}: native vol {m['vol']*100:.2f}%  =>  target-vol as-run ~{m['vol']*100:.1f}%, "
              f"lambda=1, FULL maxDD {m['maxdd']*100:.2f}%")

    # ---- NUMERO-DECISIONE: max CAGR con P(rovina-50%|5y) <= 1%
    print("\n" + "=" * 78)
    print("NUMERO-DECISIONE: max target-vol/CAGR con P(rovina-50%|5y) <= 1.0%")
    print("  (headline = block-20 baseline; verifica block-10 e coda-iniettata)")
    print("=" * 78)
    dec = []
    for bk, cfg in books.items():
        n = len(built[(bk, "FULL")])
        idx20 = block_index(n, 20, NPATH, B)
        idx10 = block_index(n, 10, NPATH, B)
        native_min = float(built[(bk, "FULL")].values.min())
        crash = 1.5 * native_min
        for lv in cfg["levels"]:
            b = built[(bk, lv)]
            r_nat = b.values.astype(float)
            bh = b[b.index >= HOLDOUT]
            nat_vol = float(r_nat.std() * np.sqrt(DPY))
            R20 = r_nat[idx20]; R10 = r_nat[idx10]
            R20c = R20.copy()
            _pos = RNG.integers(0, R20c.shape[1], size=R20c.shape[0])
            R20c[np.arange(R20c.shape[0]), _pos] = crash
            best = None
            for tv in TV_FINE:
                lam = tv / nat_vol
                r20, _ = ruin_dd(R20, lam)
                if r20 <= 0.01:
                    best = tv
                else:
                    break
            if best is None:
                dec.append(dict(book=bk, level=lv, tvmax=np.nan))
                print(f"  {bk:8s}/{lv:8s}: nessun punto con P(ruin)<=1% (troppo grasso)")
                continue
            lam = best / nat_vol
            mF = full_metrics(b, lam); mH = full_metrics(bh, lam)
            r20, d20 = ruin_dd(R20, lam)
            r10, d10 = ruin_dd(R10, lam)
            rc, dc = ruin_dd(R20c, lam)
            dec.append(dict(book=bk, level=lv, tvmax=best, lam=round(lam, 3),
                            cagr_full=mF['cagr'], cagr_hold=mH['cagr'], maxdd=mF['maxdd'],
                            pruin20=r20, pruin10=r10, pruin_crash=rc, pdd30=d20))
            print(f"  {bk:8s}/{lv:8s}: tvol_max {best*100:4.1f}% (lam {lam:.2f}) -> "
                  f"CAGR_full {mF['cagr']*100:5.1f}% / CAGR_hold {mH['cagr']*100:5.1f}% ; "
                  f"maxDDfull {mF['maxdd']*100:4.1f}% ; P(ruin) b20 {r20*100:.2f}% / "
                  f"b10 {r10*100:.2f}% / +crash {rc*100:.2f}% ; P(DD>30) {d20*100:.2f}%")
    pd.DataFrame(dec).to_csv(os.path.join(OUT, 'decision_pruin1pct.csv'), index=False)

    # ---- €/giorno a 2k e 5k dal CAGR (haircut = frontiera fidata)
    print("\n" + "-" * 78)
    print("TRADUZIONE €/giorno (CAGR haircut = frontiera fidata):")
    for bk in books:
        row = next((x for x in dec if x['book'] == bk and x['level'] == 'HAIRCUT'), None)
        if row and not np.isnan(row.get('tvmax', np.nan)):
            for cap in (2000, 5000):
                eur = cap * row['cagr_full'] / 365.25
                print(f"  {bk}/HAIRCUT @ tvol {row['tvmax']*100:.1f}%: {cap}EUR -> "
                      f"~{eur:.2f} EUR/giorno (CAGR {row['cagr_full']*100:.1f}%)")

    print(f"\ndone in {time.time()-t0:.0f}s")


if __name__ == "__main__":
    run()