"""F4-conf-03: harmonise five printed effect sizes to target definitions.
Run from workspace root: ~/wp29/venv/bin/python out/analysis.py
All inputs are parsed from data/ extracts; the source lines used are printed."""
import json, math, os, re

D = "data"
S = os.path.join(D, "sources")

def lines(fn):
    with open(fn, encoding="utf-8") as f:
        return f.read().split("\n")

def find(fn, pattern):
    """Return (line_no, line) of the unique-first line matching regex pattern."""
    for i, l in enumerate(lines(fn), 1):
        if re.search(pattern, l):
            print(f"[{os.path.basename(fn)} L{i}] {l.strip()}")
            return i, l
    raise ValueError(f"pattern {pattern!r} not found in {fn}")

def cells(line):
    return [c.strip() for c in line.split("|")]

def num(s):
    return float(re.sub(r"[^0-9.\-]", "", s))

effects = json.load(open(os.path.join(D, "effects.json"), encoding="utf-8"))
ids = [e["id"] for e in effects]
assert ids == ["E1", "E2", "E3", "E4", "E5"]

H = {}
issues = []

# ---------- E1: DCR overall ATT, 5 km, log points -> percent ----------
f1 = os.path.join(S, "arxiv-2609.27732v1.txt")
_, hdr = find(f1, r"^\| 1 km \| 3 km \| 5 km$")
ln, row = find(f1, r"^ATT \| -0\.399")
c = cells(row); att5_str = c[3]; att5 = num(att5_str)        # "-0.431**"
_, pct_row = find(f1, r"^% change \| -32\.9%")
printed_pct = num(cells(pct_row)[3])
v1 = 100 * (math.exp(att5) - 1)
H["E1"] = dict(
    value=round(v1, 2), unit="percent change in the drug-related hospitalisation rate (relative to counterfactual)",
    transformation=f"100*(exp(b)-1) with b = {att5:.3f} (printed '{att5_str}') -> 100*(exp(-0.431)-1) = {v1:.4f}%; "
                   f"matches the printed '% change' row ({printed_pct}%). Exact under the log-point convention.",
    harmonisable=True, missing_quantity=None,
    locator=f"arxiv-2609.27732v1 | Table 6 | row 'ATT' (line {ln}) | column '5 km' ; cross-check: Table 6 | row '% change' | column '5 km'")

# ---------- E2: RESEA reemployment ATE pp -> percent of control mean ----------
f2 = os.path.join(S, "arxiv-2606.14621v1.txt")
ln2, row2 = find(f2, r"^Reemployment \| 75")
c = cells(row2)
ctrl = num(c[2])                        # "75.0%"
ate = num(c[4].split(";")[0])           # "1.47 pp"
v2 = 100 * ate / ctrl
H["E2"] = dict(
    value=round(v2, 2), unit="percent change in the reemployment probability (relative to control-group mean)",
    transformation=f"100*Delta/p with Delta = {ate:.2f} pp (ATE) and p = {ctrl:.1f}% (control mean) -> 100*1.47/75.0 = {v2:.4f}%",
    harmonisable=True, missing_quantity=None,
    locator=f"arxiv-2606.14621v1 | Table 2 | row 'Reemployment' (line {ln2}) | column 'Average treatment effect estimate' (1.47 pp) and column 'Control' (75.0%)")

# ---------- E3: MW elasticity x printed % MW change ----------
f3 = os.path.join(S, "arxiv-2505.04555v3.txt")
ln3a, r3a = find(f3, r"^Employment elasticity w\.r\.t MW")
eta = num(cells(r3a)[1])
ln3b, r3b = find(f3, r"^% MW changes")
g = num(cells(r3b)[1])
_, note3 = find(f3, r"about 4% minimum wage increase")
_, r_miss = find(f3, r"^Missing jobs below new MW")
_, r_exc = find(f3, r"^Excess jobs above new MW")
miss, exc = num(cells(r_miss)[1]), num(cells(r_exc)[1])
v3 = 100 * eta * g
H["E3"] = dict(
    value=round(v3, 2), unit="percent change in employment (for the observed minimum-wage increase)",
    transformation=f"100*eta*g with eta = {eta:.3f} and g = {g:.3f} (proportional MW change printed in the same table) -> 100*(-0.388)*0.047 = {v3:.4f}%. "
                   "Linear (first-order) approximation; exact only if the elasticity is interpreted as a ratio of proportional changes.",
    harmonisable=True, missing_quantity=None,
    locator=f"arxiv-2505.04555v3 | Table 2 | row 'Employment elasticity w.r.t MW' (line {ln3a}) | column 'Estimate (SE)' ; Table 2 | row '% MW changes' (line {ln3b}) | column 'Estimate (SE)'")
cs_check = (miss + exc) / g
print(f"Check E3: (missing+excess)/%MW = ({miss}+{exc})/{g} = {cs_check:.4f} vs printed elasticity {eta}")

# ---------- E4: fatal hospitalisations, log points -> absolute per 100,000 ----------
ln4, r4 = find(f1, r"^Drug-related \(fatal\)")
b4_str = cells(r4)[3]; b4 = num(b4_str)
pct4 = 100 * (math.exp(b4) - 1)
_, abs4 = find(f1, r"fatal outcomes, which fall by 35%")
_, mech4 = find(f1, r"by approximately 30% at 5 km")
print(f"E4 relative-change equivalent (not the target): 100*(exp({b4})-1) = {pct4:.4f}%")
H["E4"] = dict(
    value=None, unit="hospitalisations per 100,000 inhabitants per year (absolute change)",
    transformation=f"Would require Delta_abs = Y0*(exp(b)-1) with b = {b4:.3f} (printed '{b4_str}'); only the relative change "
                   f"100*(exp(-0.354)-1) = {pct4:.2f}% can be computed from the extract. Y0 is not printed, so the target cannot be reached.",
    harmonisable=False,
    missing_quantity="baseline (counterfactual / pre-treatment mean) level of the fatal drug-related hospitalisation rate per 100,000 inhabitants per year for the treated 5 km regions (not printed anywhere in the extract; also the time unit of the rate is not stated)",
    locator=f"arxiv-2609.27732v1 | Table 7 | row 'Drug-related (fatal)' (line {ln4}) | column '5 km'")

# ---------- E5: LPM coefficient (0-1) -> percentage points ----------
f5 = os.path.join(S, "arxiv-2510.03658v1.txt")
_, hdr5 = find(f5, r"^\| \(1\) \| \(2\) \| \(3\) \| \(4\)")
ln5, r5 = find(f5, r"^Treatment \| -0\.0911")
b5 = num(cells(r5)[4])
_, mean5 = find(f5, r"^Mean Dep\. Var\.")
_, pc5 = find(f5, r"^% Change \|")
m5 = num(cells(mean5)[4]); pcp5 = num(cells(pc5)[4])
v5 = 100 * b5
print(f"Check E5: b/mean = {b5}/{m5} = {b5/m5:.4f} vs printed % Change {pcp5}")
H["E5"] = dict(
    value=round(v5, 2), unit="percentage points (change in probability of infant death or stillbirth)",
    transformation=f"100*b with b = {b5:.4f} (0-1 probability scale) -> {v5:.2f} pp. Exact rescaling.",
    harmonisable=True, missing_quantity=None,
    locator=f"arxiv-2510.03658v1 | Table 3 | row 'Treatment' (line {ln5}) | column '(4)'")

# ---------- data issues detected ----------
_, se_a19 = find(f1, r"^\| \(0\.218\)")
_, n640 = find(f1, r"never-treated comparison regions \(640 regions\)")
_, t2 = find(f5, r"result in table 2 suggests")
issues = [
    {"issue": "arxiv-2609.27732v1: the abstract states fatal drug-related hospitalisations 'fall by 35%%', but Table 7 (5 km) prints -0.354 log points, i.e. 100*(exp(-0.354)-1) = %.1f%%; Section 6.1 says 'approximately 30%%'. The 35%% figure matches the overall outcome (Table 6), not the fatal sub-outcome." % pct4,
     "handling": "Used the Table 7 coefficient as specified in effects.json; flagged the abstract figure as internally inconsistent. Does not affect E4 (not harmonisable anyway)."},
    {"issue": "arxiv-2609.27732v1: the same CS-DR 5 km ATT (-0.431) has SE (0.210) in Table 6 but (0.218) in Table A19; the A19 note states 640 never-treated comparison regions while its table body and Table 6 state 608.",
     "handling": "Point estimates agree, so E1 is unaffected; no uncertainty propagation is required. Recorded only."},
    {"issue": "arxiv-2505.04555v3: Table 2 note describes the 2023 increase as 'about 4%' while the table row '% MW changes' prints 0.047 (4.7%).",
     "handling": "Used the printed table value 0.047, as instructed by the target definition ('as printed in the same table'). With the note's 'about 4%%' (g = 0.04) the result would be %.2f percent; reported as sensitivity only." % (100*eta*0.04)},
    {"issue": "arxiv-2505.04555v3: Table 2 rows are on different scales (e.g. 'Missing jobs below new MW' -0.03 printed with two decimals); the employment elasticity is approximately (missing+excess)/(MW change) = %.3f, consistent with -0.388 up to rounding." % cs_check,
     "handling": "Consistency check only; no change."},
    {"issue": "arxiv-2510.03658v1: the abstract/Section 4.1 headline (-0.0855, mean 0.455, '18.79 percentage change') refers to a Table 2 not included in the extract; it is also labelled inconsistently ('percentage change' vs a relative change). Table 3 col (4) relative change printed as -0.3685 vs -0.0962/0.261 = %.4f (rounding of the printed mean)." % (b5/m5),
     "handling": "Used Table 3 col (4) as specified; the Table 2 numbers were not used."},
    {"issue": "data/README.md mentions a re-set 'London Table 3 header', but no extract in this task concerns a London study (apparent copy-over from another task's README).",
     "handling": "Ignored; no numbers depend on it."},
    {"issue": "arxiv-2606.14621v1: the Table 2 'Control' mean (75.0%) is the raw aggregate control-group mean, not a block(week)-adjusted counterfactual level; the authors stress that raw group comparisons are biased by week of randomisation.",
     "handling": "Used 75.0% because the target definition explicitly asks for percent of the control-group mean; noted that a block-adjusted counterfactual mean is not printed."},
]

# ---------- input row counts ----------
def wc_l(fn):
    with open(fn, "rb") as f:
        return f.read().count(b"\n")
n_rows = {"effects.json": len(effects), "README.md": wc_l(os.path.join(D, "README.md"))}
for fn in sorted(os.listdir(S)):
    n_rows[fn] = wc_l(os.path.join(S, fn))

n_h = sum(1 for k in H if H[k]["value"] is not None)
n_nh = sum(1 for k in H if not H[k]["harmonisable"])
res = {
    "harmonised": H,
    "n_harmonised": n_h,
    "n_not_harmonisable": n_nh,
    "task_id": "F4-conf-03",
    "feasibility": "complete",
    "feasibility_reason": None,
    "data_issues": issues,
    "n_input_rows": n_rows,
    "n_analysed": len(effects),
    "random_seed": None,
}
os.makedirs("out", exist_ok=True)
with open(os.path.join("out", "results.json"), "w", encoding="utf-8") as f:
    json.dump(res, f, indent=2, ensure_ascii=False)
print(json.dumps({k: H[k]["value"] for k in H}, indent=1))