"""Does tick-derived order flow predict the NEXT bars, or only explain its own?

This is the test the whole tick pipeline exists for, and it is the first configuration in
this project with enough independent trades to settle the question rather than fail to
reject it. H1 with a 128-bar horizon gives ~300 independent trades in 18 years - enough to
resolve a +6pp edge when real edges are 1-3pp. M5 with a 12-48 bar horizon over 23 years
gives tens of thousands, which resolves ~1pp.

THE PRIOR, STATED BEFORE LOOKING
--------------------------------
Order-flow imbalance is well established as a CONTEMPORANEOUS explainer of price change
(Cont/Kukanov/Stoikov), and its predictive power is reported to decay within seconds. On
the build sample OFI correlated +0.56 with the SAME-bar return, reproducing that. So the
honest expectation is that it explains the bar it is measured in and says nothing about the
next one, and the contemporaneous correlation must never be quoted as evidence of edge.

What is being tested is the gap between that literature (equities, sub-second, size-weighted
book data) and this setting (retail FX CFD feed, 5-minute bars, event-count OFI without
sizes). That gap is worth one honest measurement.

DISCIPLINE
----------
Signals are computed on bar i and entered at the OPEN of bar i+1. Barriers are scanned from
the entry bar forward only. Trades are sequential and NON-OVERLAPPING, so each is
independent and the confidence interval means something - the pseudo-replication that
produced a fake +2.66pp at 2.9 sigma earlier in this project came from exactly this being
skipped. Break-even == chance by the gambler's-ruin identity, so "beats a coin" and "makes
money" are one question, and the spread is charged inside every barrier.
"""
import numpy as np, sys, os, time
sys.stdout.reconfigure(encoding='utf-8', errors='replace')

BARS = 'c:/Users/admin/Documents/Workspaces/Market Data/bars/'


#--- bars per day, used to keep the trailing z-score window at ~1 day on every timeframe
PER_DAY = {'M5': 288, 'M15': 96, 'H1': 24, 'H4': 6}


def load(sym, tf='M5'):
    z = np.load(f"{BARS}{sym}_{tf}_ticks.npz", allow_pickle=True)
    arr = z['bars']
    cols = [str(c) for c in z['columns']]
    return arr, {c: k for k, c in enumerate(cols)}


def atr_bars(h, l, c, n=14):
    pc = np.roll(c, 1); pc[0] = c[0]
    tr = np.maximum(h - l, np.maximum(np.abs(h - pc), np.abs(l - pc)))
    out = np.convolve(tr, np.ones(n) / n, mode='full')[:len(tr)]
    out[:n] = tr[:n].mean()
    return out


def zscore(x, n=288):
    """Trailing z-score over n bars (288 M5 bars = one day), causal and shifted by one so a
    bar never contributes to its own baseline."""
    x = np.asarray(x, dtype=float)
    c1 = np.concatenate([[0.0], np.cumsum(x)])
    c2 = np.concatenate([[0.0], np.cumsum(x * x)])
    out = np.zeros(len(x))
    for i in range(n + 1, len(x)):
        s = c1[i] - c1[i - n]
        s2 = c2[i] - c2[i - n]
        mu = s / n
        var = max(s2 / n - mu * mu, 1e-18)
        out[i] = (x[i] - mu) / np.sqrt(var)
    return np.clip(out, -8, 8)


def build_signals(a, I, zwin=288):
    g = lambda c: a[:, I[c]]
    up, dn = g('upticks'), g('downticks')
    bu, bd, au, ad = g('bid_up'), g('bid_dn'), g('ask_up'), g('ask_dn')
    tot = np.maximum(up + dn, 1.0)
    # --- the two SIGNED candidates. Everything else this pipeline produces is unsigned and
    # cannot point a trade, however well it measures.
    tick_imb = (up - dn) / tot
    ofi_raw = (bu + au) - (bd + ad)
    ofi = ofi_raw / np.maximum(bu + bd + au + ad, 1.0)
    sig = {
        'tick_imbalance': tick_imb,
        'ofi': ofi,
        'ofi_z': zscore(ofi_raw, zwin),
        'tick_imb_z': zscore(tick_imb, zwin),
    }
    return sig


def barrier_outcomes(o, h, l, a_sig, sl_m, tp_m, H, sp):
    """Entry at the OPEN of bar i, both directions. Stop tested before target within a bar,
    so a bar spanning both books the loss.

    `sp` is PER BAR, not a global median: this data spans 2003-2026 and FX spreads narrowed
    by roughly an order of magnitude over it. A single median charges the modern cost to the
    2000s and the 2000s cost to today, which flatters exactly the era with the most bars.

    Timeouts (neither barrier touched within H) are returned separately rather than folded
    silently into the loss column - a timeout is not a stop-out, and if they are common the
    quoted break-even no longer describes the experiment."""
    n = len(o)
    INF = np.iinfo(np.int32).max
    winL = np.zeros(n, bool); winS = np.zeros(n, bool)
    toL = np.zeros(n, bool); toS = np.zeros(n, bool)
    risk, rew = sl_m * a_sig, tp_m * a_sig
    lTp, lSl = o + rew + sp, o - risk + sp
    sTp, sSl = o - rew - sp, o + risk - sp
    CH = max(4000000 // max(H, 1), 1)
    for s in range(0, n, CH):
        e2 = min(s + CH, n - H)
        if e2 <= s:
            break
        wi = np.arange(0, H)[None, :] + np.arange(s, e2)[:, None]
        wh, wl = h[wi], l[wi]

        def first(mask):
            any_ = mask.any(axis=1)
            return np.where(any_, mask.argmax(axis=1), INF)
        lsl = first(wl <= lSl[s:e2, None]); ltp = first(wh >= lTp[s:e2, None])
        ssl = first(wh >= sSl[s:e2, None]); stp = first(wl <= sTp[s:e2, None])
        winL[s:e2] = ltp < lsl
        winS[s:e2] = stp < ssl
        toL[s:e2] = (ltp == INF) & (lsl == INF)
        toS[s:e2] = (stp == INF) & (ssl == INF)
    return winL, winS, toL, toS


def sequential(fire, dirs, winL, winS, H, n):
    out_i, out_w, out_d = [], [], []
    busy = -1
    for j, d in zip(fire, dirs):
        if j <= busy or j + 1 + H >= n:
            continue
        e = j + 1
        out_i.append(e)
        out_w.append(bool(winL[e] if d > 0 else winS[e]))
        out_d.append(int(d))
        busy = e + H
    return np.array(out_i, int), np.array(out_w, bool), np.array(out_d, int)


def run(sym, sl_m, tp_m, H, thresholds=(0.5, 1.0, 1.5, 2.0), nperm=2000, seed=3, tf='M5',
        cost_q=None):
    a, I = load(sym, tf)
    g = lambda c: a[:, I[c]]
    o, h, l, c = g('open'), g('high'), g('low'), g('close')
    n = len(c)
    atr = atr_bars(h, l, c, 14)
    a_sig = np.concatenate([[atr[0]], atr[:-1]])
    #--- per-bar spread, shifted one bar so entry cost is known before entering
    spb = g('spread_mean')
    spb = np.concatenate([[spb[0]], spb[:-1]])
    winL, winS, toL, toS = barrier_outcomes(o, h, l, a_sig, sl_m, tp_m, H, spb)
    be = sl_m / (sl_m + tp_m)
    sigs = build_signals(a, I, PER_DAY.get(tf, 288))
    rng = np.random.default_rng(seed)

    #--- COST CONDITIONING. spread/ATR varies by an order of magnitude within a timeframe
    #--- (thin Asian hours and news spikes vs the London/NY overlap), so the average cost
    #--- is not the cost you must pay - you can choose to trade only the cheap bars. This
    #--- is the one honest use of an UNSIGNED feature: it cannot point a direction, but it
    #--- can decline to trade. Both values are known at the entry decision (shifted one
    #--- bar), so this is a filter, not hindsight.
    cheap = None
    if cost_q is not None:
        ratio = spb / np.maximum(a_sig, 1e-12)
        thr_c = np.nanquantile(ratio, cost_q)
        cheap = ratio <= thr_c
        print(f"  [cost filter] spread/ATR <= {thr_c:.4f} (lowest {100*cost_q:.0f}%), "
              f"{cheap.mean():.1%} of bars eligible")

    print(f"\n=== {sym} {tf}  SL{sl_m}:TP{tp_m}  H={H} bars  n={n:,}  "
          f"spread {np.nanmedian(spb):.6f} ({np.nanmedian(spb)/np.nanmedian(atr):.3f} ATR)"
          f"  break-even {100*be:.2f}% ===")
    print(f"{'signal':<16}{'thr':>6}{'trades':>9}{'win%':>8}{'null%':>8}{'vs null':>9}"
          f"{'z':>7}{'exp R':>8}{'expR rev':>10}{'t/o%':>7}")
    acc = []
    rows = []
    for name, s in sigs.items():
        scale = 1.0 if name.endswith('_z') else 1.0
        for thr in thresholds:
            t = thr if name.endswith('_z') else thr * 0.25
            m = np.abs(s) >= t
            m[:PER_DAY.get(tf, 288) + 12] = False      # z-score burn-in
            if cheap is not None:
                m &= cheap
            fire = np.nonzero(m)[0]
            if len(fire) < 50:
                continue
            d = np.sign(s[fire]).astype(int)
            ti, tw, td = sequential(fire, d, winL, winS, H, n)
            if len(ti) < 100:
                continue
            wr = tw.mean(); nT = len(ti)
            to = np.where(td > 0, toL[ti], toS[ti]).mean()
            #--- the REVERSED rule on the same bars: the money question if the signal turns
            #--- out to be anti-predictive. Not a second hypothesis - it is the same test
            #--- read backwards, so it does not enlarge the family.
            wrev = np.where(td > 0, winS[ti], winL[ti]).mean()
            expR = wr * tp_m - (1 - wr) * sl_m
            expRrev = wrev * tp_m - (1 - wrev) * sl_m
            # NULL: keep the firing bars and the long/short MIX, shuffle which trade gets
            # which direction. A 50/50 coin flip would be the wrong null for a directionally
            # skewed signal on a trending instrument - it would let drift alone look like
            # timing skill. Permuting the observed directions holds the mix fixed and tests
            # only the pairing of direction to bar, which is the actual claim.
            wl_, ws_ = winL[ti], winS[ti]
            long_ = td > 0
            pw = np.empty(nperm)
            # batched: a (nperm x nT) array is ~2 GB at these trade counts
            B = max(1, 4000000 // max(nT, 1))
            for b0 in range(0, nperm, B):
                b1 = min(b0 + B, nperm)
                pl = rng.permuted(np.broadcast_to(long_, (b1 - b0, nT)), axis=1)
                pw[b0:b1] = np.where(pl, wl_, ws_).mean(axis=1)
            #--- Standardise against the EMPIRICAL null, not the textbook break-even.
            #--- sl/(sl+tp) is the break-even of a costless coin. These barriers charge the
            #--- spread and book a loss when one bar spans both levels, so random entry
            #--- sits WELL below it - about 39.5% where the textbook says 50%. Measuring
            #--- against 50% reports that fixed cost as if it were signal, which produced a
            #--- -75 sigma "result" that was almost entirely the cost of trading.
            nmu, nsd = pw.mean(), max(pw.std(ddof=1), 1e-12)
            z = (wr - nmu) / nsd
            rows.append((f"{name}", t, nT, 100 * wr, 100 * nmu, 100 * (wr - nmu), z,
                         expR, expRrev, 100 * to))
            acc.append(np.abs((pw - nmu) / nsd))
    if not rows:
        print("  no signal fired often enough")
        return
    crit = float(np.quantile(np.maximum.reduce(acc), 0.95))
    for r in sorted(rows, key=lambda x: -abs(x[6])):
        star = ' *' if abs(r[6]) > crit else ''
        print(f"{r[0]:<16}{r[1]:>6.2f}{r[2]:>9d}{r[3]:>8.2f}{r[4]:>8.2f}{r[5]:>+9.2f}"
              f"{r[6]:>+7.2f}{r[7]:>+8.3f}{r[8]:>+10.3f}{r[9]:>7.1f}{star}")
    print(f"  family-wise 5% bar over {len(rows)} tests: |z| > {crit:.2f}   (* = clears it)")
    print(f"  NOTE break-even {100*be:.1f}% is the COSTLESS coin; random entry at these "
          f"barriers wins ~{np.mean([r[4] for r in rows]):.1f}% after spread.")
    print(f"  A rule only makes money if win% > {100*be:.1f}%, i.e. exp R > 0 - "
          f"beating the null is necessary but NOT sufficient.")


#--- Narrow barriers are where the signal lives and where the spread is fatal: at 1 ATR the
#--- spread is ~10% of the target and costs 13pp of win rate against a ~1pp effect. Cost
#--- falls roughly as 1/width, so 8 ATR should cost under 2pp. The honest expectation is
#--- that the effect decays with horizon faster than the cost does - order flow is a
#--- seconds-to-minutes phenomenon and these wide barriers run 8 to 24 hours - but that is
#--- exactly the trade-off worth measuring rather than assuming.
NARROW = [(1, 1, 12), (1, 2, 24), (2, 3, 48)]
WIDE = [(4, 4, 96), (4, 8, 192), (8, 8, 288)]

if __name__ == '__main__':
    args = [a for a in sys.argv[1:] if not a.startswith('-')]
    geos = WIDE if '--wide' in sys.argv else NARROW if '--narrow' in sys.argv else NARROW + WIDE
    tf, cost_q = 'M5', None
    for a_ in sys.argv[1:]:
        if a_.startswith('--tf='):
            tf = a_.split('=', 1)[1]
        if a_.startswith('--cheap='):
            cost_q = float(a_.split('=', 1)[1])
    syms = args or ['EURUSD']
    for sym in syms:
        if not os.path.exists(f"{BARS}{sym}_{tf}_ticks.npz"):
            print(f"{sym}: {tf} bars not built yet")
            continue
        for (s, p, H) in geos:
            run(sym, s, p, H, tf=tf, cost_q=cost_q)