## Task 3 analysis: compare the FIFO (before) and fireTick-keyed (after) pairing
## for the three learned GF guns.
##
## Two analyses:
##   1. Per-fixture (per-run) paired comparison — the repo's convention. The 10
##      DrussGT fixtures are the independent runs; the paired delta is
##      after - before. An exact sign-flip permutation test (2^10 = 1024 sign
##      patterns) gives the p-value, and the per-run ranges give the overlap.
##   2. Bullet-level two-sample permutation test on the raw hit booleans dumped
##      by audit_wave_pairing.nim. ANTI-CONSERVATIVE: bullets within a fixture
##      share a trajectory and are correlated, so treat this as an upper bound on
##      significance, not the headline.
##
## Run: nim c -r common_libs/tests/compare_pairing.nim

import std/[math, strformat, random, os, strutils]

const
  FixtureNames = ["drussgt_vs_corners", "drussgt_vs_crazy", "drussgt_vs_drussgt",
                  "drussgt_vs_ramfire", "drussgt_vs_spinbot",
                  "tr_drussgt_vs_corners", "tr_drussgt_vs_crazy",
                  "tr_drussgt_vs_modularbot", "tr_drussgt_vs_modularbot_shield",
                  "tr_drussgt_vs_spinbot"]
  ShotsPerFixture = 400   # WindowSize(100) x 4 power bins
  GunNames = ["GuessFactor", "DecayGF", "KNN"]

  # Captured from `audit_wave_pairing.nim <tag> both` (fresh guns per fixture;
  # hits out of 400). Deterministic guns -> reproducible.
  BeforePath: array[3, array[10, int]] = [
    [36, 139, 23, 202, 233, 155, 9, 47, 42, 42],   # GuessFactor
    [61, 108, 8, 204, 232, 118, 9, 44, 108, 60],   # DecayGF
    [16, 94, 33, 188, 176, 123, 12, 33, 19, 37],   # KNN
  ]
  AfterPath: array[3, array[10, int]] = [
    [36, 139, 23, 202, 226, 155, 12, 44, 42, 42],
    [61, 105, 19, 213, 226, 140, 9, 47, 108, 42],
    [38, 94, 38, 193, 158, 129, 11, 29, 25, 37],
  ]
  BeforePoint: array[3, array[10, int]] = [
    [5, 33, 15, 59, 59, 36, 0, 26, 3, 21],
    [1, 8, 2, 37, 52, 48, 0, 0, 33, 6],
    [0, 19, 14, 53, 34, 39, 13, 15, 5, 25],
  ]
  AfterPoint: array[3, array[10, int]] = [
    [5, 33, 15, 59, 59, 36, 0, 28, 3, 21],
    [1, 9, 2, 37, 52, 48, 0, 0, 45, 6],
    [0, 20, 15, 48, 33, 38, 13, 14, 5, 25],
  ]

proc sum(a: array[10, int]): int =
  for x in a: result += x

proc meanPct(a: array[10, int]): float = sum(a).float / 10.0 / ShotsPerFixture.float * 100.0

proc minPct(a: array[10, int]): float =
  result = 1e9
  for x in a: result = min(result, x.float / ShotsPerFixture.float * 100.0)

proc maxPct(a: array[10, int]): float =
  result = -1e9
  for x in a: result = max(result, x.float / ShotsPerFixture.float * 100.0)

proc signFlipP(before, after: array[10, int]): tuple[p, obsMeanPp: float, nPos, nNeg, nZero: int] =
  ## Exact sign-flip permutation test on the paired per-fixture deltas.
  var deltas: array[10, float]
  for i in 0..<10:
    deltas[i] = (after[i] - before[i]).float / ShotsPerFixture.float * 100.0
    if deltas[i] > 1e-9: inc result.nPos
    elif deltas[i] < -1e-9: inc result.nNeg
    else: inc result.nZero
    result.obsMeanPp += deltas[i] / 10.0
  let obs = abs(result.obsMeanPp)
  var ge = 0
  for mask in 0..<(1 shl 10):
    var m = 0.0
    for i in 0..<10:
      let s = if ((mask shr i) and 1) == 1: -1.0 else: 1.0
      m += s * deltas[i] / 10.0
    if abs(m) >= obs - 1e-12: inc ge
  result.p = ge.float / 1024.0

proc loadDump(path: string): seq[bool] =
  if not fileExists(path):
    return @[]
  for line in lines(path):
    let s = line.strip()
    if s.len == 0: continue
    result.add (s == "1")

proc zTest(a, b: seq[bool]): tuple[p, diffPp, z: float] =
  ## Two-proportion z-test (analytic; the permutation equivalent is exact but
  ## 344k-element shuffles are needlessly slow). ANTI-CONSERVATIVE because the
  ## bullets are correlated within a fixture.
  if a.len == 0 or b.len == 0: return (1.0, 0.0, 0.0)
  var ha, hb: int
  for x in a: (if x: inc ha)
  for x in b: (if x: inc hb)
  let p1 = ha.float / a.len.float
  let p2 = hb.float / b.len.float
  result.diffPp = (p2 - p1) * 100.0
  let p = (ha + hb).float / (a.len + b.len).float
  let se = sqrt(max(1e-30, p * (1.0 - p) * (1.0/a.len.float + 1.0/b.len.float)))
  result.z = (p2 - p1) / se
  result.p = erfc(abs(result.z) / sqrt(2.0))

proc reportMetric(mname: string,
                  before, after: array[3, array[10, int]]) =
  echo "══════════════════════════════════════════════════════════════════"
  echo "  METRIC = ", mname
  echo "══════════════════════════════════════════════════════════════════"
  for gi in 0..<2:
    let b = before[gi]
    let a = after[gi]
    echo fmt"{GunNames[gi]}:"
    echo fmt"  before  {sum(b):>4}/{ShotsPerFixture*10} = {meanPct(b):5.2f}%  per-run {minPct(b):5.2f}..{maxPct(b):5.2f}%"
    echo fmt"  after   {sum(a):>4}/{ShotsPerFixture*10} = {meanPct(a):5.2f}%  per-run {minPct(a):5.2f}..{maxPct(a):5.2f}%"
    let sf = signFlipP(b, a)
    echo fmt"  delta   {meanPct(a)-meanPct(b):+5.2f}pp   paired sign-flip permutation p={sf.p:.3f}  (+{sf.nPos}/-{sf.nNeg}/0:{sf.nZero} of 10)"
  let b = before[2]
  let a = after[2]
  echo fmt"{GunNames[2]}:"
  echo fmt"  before  {sum(b):>4}/{ShotsPerFixture*10} = {meanPct(b):5.2f}%  per-run {minPct(b):5.2f}..{maxPct(b):5.2f}%"
  echo fmt"  after   {sum(a):>4}/{ShotsPerFixture*10} = {meanPct(a):5.2f}%  per-run {minPct(a):5.2f}..{maxPct(a):5.2f}%"
  let sf = signFlipP(b, a)
  echo fmt"  delta   {meanPct(a)-meanPct(b):+5.2f}pp   paired sign-flip permutation p={sf.p:.3f}  (+{sf.nPos}/-{sf.nNeg}/0:{sf.nZero} of 10)"
  echo ""

proc main() =
  randomize(12345)
  reportMetric("bmPath (shipped)", BeforePath, AfterPath)
  reportMetric("bmPoint", BeforePoint, AfterPoint)

  echo "══════════════════════════════════════════════════════════════════"
  echo "  bullet-level two-proportion z-test (bmPath dumps) — ANTI-CONSERVATIVE"
  echo "══════════════════════════════════════════════════════════════════"
  for gi in 0..<3:
    let bf = loadDump(fmt"/tmp/wavepair_before_path_{GunNames[gi]}.txt")
    let af = loadDump(fmt"/tmp/wavepair_after_path_{GunNames[gi]}.txt")
    let r = zTest(bf, af)
    var hb, ha: int
    for x in bf: (if x: inc hb)
    for x in af: (if x: inc ha)
    echo fmt"{GunNames[gi]:<12} before {hb:>6}/{bf.len:<6} {hb.float/bf.len.float*100:5.2f}%  after {ha:>6}/{af.len:<6} {ha.float/af.len.float*100:5.2f}%  diff {r.diffPp:+5.2f}pp  z={r.z:+5.2f}  p={r.p:.2g}"

when isMainModule:
  main()