"""GOES-16 rescaling test (lead-1 M-class). Periods: C23 1996-08..2008-12, C24 2009-01..2019-12-08 (legacy GOES, 0.7 factor),
C25a 2019-12-09..2024-12 and C25b 2025-01..2026-10 (GOES-16+ unscaled). Activity covariate scale-independent of GOES:
F10.7 90-day mean stated in the RSGA (f107_90d). Outputs results/scale_bins.csv and results/scale_glm.csv.
GLM: logistic regression of (a) operational outcome, (b) legacy-scale outcome, (c) SWPC probability (fractional logit
via weighted duplicate rows) on log F10.7_90d + period dummies (C24 = reference); exp(coef) = odds ratio vs C24."""
import numpy as np, pandas as pd
from sklearn.linear_model import LogisticRegression
import os; R = os.environ.get('FLARE_ROOT', '/work/flare-forecasts')
rs = pd.read_csv(f'{R}/data/rsga.csv', parse_dates=['issue_date'])[['issue_date', 'f107_90d']].rename(columns={'issue_date': 'issue'})
a = pd.read_csv(f'{R}/data/flare_pairs.csv', parse_dates=['issue']); b = pd.read_csv(f'{R}/data/flare_pairs_oldscale.csv', parse_dates=['issue'])
out = []
for ev in ('M', 'X'):
    d = a[(a.event == ev) & (a.lead == 1)][['issue', 'p', 'y']].merge(b[(b.event == ev) & (b.lead == 1)][['issue', 'y']].rename(columns={'y': 'y_old'}), on='issue').merge(rs, on='issue').dropna()
    d['period'] = pd.cut(d.issue, pd.to_datetime(['1996-01-01', '2009-01-01', '2019-12-09', '2025-01-01', '2027-01-01']), labels=['C23', 'C24', 'C25a', 'C25b'], right=False)
    d['fbin'] = pd.cut(d.f107_90d, [60, 80, 100, 125, 150, 175, 200, 300])
    t = d.groupby(['period', 'fbin'], observed=True).agg(n=('p', 'size'), mean_p=('p', 'mean'), freq_oper=('y', 'mean'), freq_old=('y_old', 'mean')).reset_index()
    t.insert(0, 'event', ev); t['fbin'] = t.fbin.astype(str); out.append(t)
    Xd = pd.get_dummies(d.period).drop(columns='C24').astype(float); Xd['logF'] = np.log(d.f107_90d)
    res = {'event': ev}
    for nm, y in (('oper', d.y), ('old', d.y_old)):
        m = LogisticRegression(C=1e4, max_iter=5000).fit(Xd.values, y.values)
        for k, c in zip(Xd.columns, m.coef_[0]): res[f'OR_{nm}_{k}'] = np.exp(c)
    XX = np.vstack([Xd.values, Xd.values]); yy = np.r_[np.ones(len(d)), np.zeros(len(d))]; ww = np.r_[d.p.values, 1 - d.p.values]
    m = LogisticRegression(C=1e4, max_iter=5000).fit(XX, yy, sample_weight=ww)
    for k, c in zip(Xd.columns, m.coef_[0]): res[f'OR_swpc_{k}'] = np.exp(c)
    out.append(pd.DataFrame([res]))
    print(ev, {k: round(v, 2) for k, v in res.items() if k != 'event'})
bins = pd.concat([o for o in out if 'fbin' in o]); bins.to_csv(f'{R}/results/scale_bins.csv', index=False)
pd.concat([o for o in out if 'fbin' not in o]).to_csv(f'{R}/results/scale_glm.csv', index=False)
pd.set_option('display.width', 200); print(bins[bins.event == 'M'].round(3).to_string(index=False))