"""How far the sequence-updated productivity a moves from the generic prior: a(version nearest day 7, Bayesian) - a(generic),
per mainshock, by region. Negative = sequence less productive than the generic prior."""
import sys; sys.path.insert(0, '/work/aftershock-forecasts')
import csv, numpy as np, collections
from region import regions
REG = regions()
V = collections.defaultdict(list)
for v in csv.DictReader(open('/work/aftershock-forecasts/versions.csv')):
    if v['a']: V[v['eid']].append(v)
rows = []
for e, vs in V.items():
    g = [v for v in vs if 'Generic' in v['model']]
    b = [v for v in vs if 'Bayesian' in v['model'] and 'Reasenberg' in v['model']]
    if not g or not b: continue
    lagd = lambda v: (int(v['creation_ms']) - int(v['time_ms'])) / 864e5
    v7 = min(b, key=lambda v: abs(lagd(v) - 7))
    if abs(lagd(v7) - 7) > 1: continue
    rows.append(dict(eid=e, reg=REG[e], mag=float(v7['mag']), a_generic=float(g[0]['a']), a_day7=float(v7['a']), da=float(v7['a']) - float(g[0]['a'])))
with open('/work/aftershock-forecasts/results/productivity_shift.csv', 'w', newline='') as f:
    w = csv.DictWriter(f, fieldnames=list(rows[0])); w.writeheader(); w.writerows(rows)
for reg in ['ALL', 'CONUS+HI', 'Alaska', 'Marianas/Guam', 'Non-US', 'Puerto Rico']:
    d = np.array([r['da'] for r in rows if reg == 'ALL' or r['reg'] == reg])
    if len(d): print(f'{reg:14s} n={len(d):3d} median da={np.median(d):+.2f} mean={d.mean():+.2f} frac<0={np.mean(d < 0):.2f} -> median factor 10^da={10**np.median(d):.2f}')
# generic parameters and day-30 update for deep (>=70 km) vs shallow mainshocks
print('--- generic a by depth, and median a(day 30) - a(generic)')
for reg in ['Alaska', 'Marianas/Guam']:
    for deep in (False, True):
        G, D = [], []
        for e, vs in V.items():
            if REG[e] != reg or (float(vs[0]['depth']) >= 70) != deep: continue
            g = [v for v in vs if 'Generic' in v['model']]
            b = [v for v in vs if 'Bayesian' in v['model'] and 'Reasenberg' in v['model']]
            if not g or not b: continue
            lagd = lambda v: (int(v['creation_ms']) - int(v['time_ms'])) / 864e5
            v30 = min(b, key=lambda v: abs(lagd(v) - 30))
            G.append(float(g[0]['a'])); D.append(float(v30['a']) - float(g[0]['a']))
        if G: print(f"{reg:14s} {'deep' if deep else 'shallow':7s} n={len(G):3d} generic a values {sorted(set(G))} median shift by day 30 {np.median(D):+.2f}")