"""Coverage on the 30-row verified sample: Karak as printed vs as verified in the original papers (verify_A/B.jsonl).
Corrections below are transcribed from the verification JSON (see notes column there)."""
import csv, numpy as np
from score_karak import crps_gauss
OBS = 160.9
K = {int(r['row']): r for r in csv.DictReader(open('karak_table1.csv'))}
sample = [int(r['row']) for r in csv.DictReader(open('verify_sample.csv'))]
# verified (value, err_plus, err_minus, meaning); None = drop (not a Cycle 25 SN forecast); err None = no error stated
V = {2: None,                       # Cycle 24 forecast, not 25
     3: (119 * 1.43, 7.5 * 1.43, 7.5 * 1.43, 'range'),
     18: None,                      # no sunspot-number forecast in the paper
     26: (118, 21, 9, 'range'),
     28: (125, 32, 32, '2sigma'),
     41: (89, 29, 14, 'unstated'),
     56: (106, 19.75, 19.75, 'other'),
     68: (95, 20, 15, 'other'),
     100: [(131, 14, 14, 'ensemble sd'), (137, 2, 2, 'ensemble sd')],
     118: (127.59, 10.03, 10.03, 'ensemble sd'),
     127: (178.5, None, None, 'none')}
def printed(row):
    r = K[row]; s = 1.43 if int(r['year']) <= 2015 else 1.0
    return [(float(r['value']) * s, float(r['err_plus']) * s, float(r['err_minus']) * s, '')]
def verified(row):
    if row not in V: return printed(row)
    v = V[row]
    if v is None: return []
    return v if isinstance(v, list) else [v]
def stats(F, label, two_sigma_halved=False):
    zs, cr = [], []
    for v, ep, em, m in F:
        if ep is None:
            cr.append(abs(OBS - v)); continue
        if two_sigma_halved and m == '2sigma':
            ep, em = ep / 2, em / 2
        s = ep if OBS > v else em
        zs.append(abs(OBS - v) / s); cr.append(crps_gauss(v, s, OBS))
    zs = np.array(zs)
    print(f'{label:40s} n={len(F):2d} with err {len(zs):2d}: cover ±1 {np.mean(zs <= 1):.0%}, ±2 {np.mean(zs <= 2):.0%}; '
          f'median |z| {np.median(zs):.2f}; mean CRPS {np.mean(cr):.1f}; below outcome {np.mean([f[0] < OBS for f in F]):.0%}')
P = [x for r in sample for x in printed(r)]
Q = [x for r in sample for x in verified(r)]
stats(P, 'sample, Karak as printed')
stats(Q, 'sample, as verified (stated ±)')
stats(Q, 'sample, as verified (2sigma -> 1sigma)', True)
MISMATCH = [2, 18, 26, 41, 51, 56, 68, 100, 127]   # verifier 'matches_karak = no'; 3 and 118 are restatements/partial matches
print(f'rows differing substantively from the original paper: {len(MISMATCH)} of {len(sample)}; 1 unverifiable (row 96)')