"""Is the initial-forecast over-forecast an artifact? (a) recount initial-forecast windows using the circle of each
mainshock's LAST forecast version (centroid-based, analyst-reviewed radius) and with a 2x radius; (b) compare the
initial magMain with the final ComCat magnitude."""
import sys; sys.path.insert(0, '/work/aftershock-forecasts')
from analyze import O, boot
from common import count
import csv, numpy as np, collections
V = collections.defaultdict(list)
for v in csv.DictReader(open('versions.csv')):
    if v['rkm']: V[v['eid']].append(v)
last = {e: max(vs, key=lambda v: int(v['creation_ms'])) for e, vs in V.items()}
vinfo = {v['vid']: v for vs in V.values() for v in vs}
C = [o for o in O if o['complete'] and o['window'] != '1 Year']
first = {}
for o in C:
    if o['eid'] not in first or int(o['creation_ms']) < first[o['eid']]: first[o['eid']] = int(o['creation_ms'])
F1 = [o for o in C if int(o['creation_ms']) == first[o['eid']] and (int(o['creation_ms']) - int(o['time_ms'])) / 6e4 <= 60]
def rep(name, T, yk):
    r, a, b = boot(T, lambda o: o[yk], lambda o: o['p'])
    print(f"{name:58s} obs={sum(o[yk] for o in T):5d} exp={sum(o['p'] for o in T):6.1f} obs/exp={r:.2f} [{a:.2f},{b:.2f}]")
for o in F1:
    L = last[o['eid']]; v = vinfo[o['vid']]
    o['y_last'] = int(count(o['eid'], float(L['rlat']), float(L['rlon']), float(L['rkm']), int(o['t0_ms']), int(o['t1_ms']), o['M']) >= 1)
    o['y_2r'] = int(count(o['eid'], float(v['rlat']), float(v['rlon']), 2 * float(v['rkm']), int(o['t0_ms']), int(o['t1_ms']), o['M']) >= 1)
    o['moved'] = (L['rlat'], L['rlon'], L['rkm']) != (v['rlat'], v['rlon'], v['rkm'])
rep('initial forecasts, own circle', F1, 'y')
rep('initial forecasts, counted in the last version\'s circle', F1, 'y_last')
rep('initial forecasts, own center, 2x radius', F1, 'y_2r')
print('share of mainshocks whose last circle differs from the initial one:', np.mean([o['moved'] for o in F1]).round(2))
# (b) magnitude revisions
ev = {f['id']: f for f in __import__('json').load(open('raw/eventlist.json'))['features']}
d = []
for e in {o['eid'] for o in F1}:
    vs = sorted(V[e], key=lambda v: int(v['creation_ms']))
    if vs[0]['magMain']: d.append(float(vs[0]['magMain']) - ev[e]['properties']['mag'])
d = np.array(d)
print(f'initial magMain - final ComCat mag: n={len(d)} median={np.median(d):+.2f} mean={d.mean():+.3f} share>+0.1={np.mean(d > 0.1):.2f} share<-0.1={np.mean(d < -0.1):.2f}')