"""F4-conf-03: harmonise five effect sizes from printed quantities only.
Run from workspace root: ~/wp29/venv/bin/python out/analysis.py
"""
import json, math, re, os

ROOT = os.path.dirname(os.path.dirname(os.path.abspath(__file__)))
D = os.path.join(ROOT, "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(fname, pattern):
    """Return (line_no, line) of the unique line starting with pattern; print it."""
    L = lines(os.path.join(S, fname))
    hits = [(i + 1, l) for i, l in enumerate(L) if l.startswith(pattern)]
    assert len(hits) == 1, (fname, pattern, hits)
    print(f"[{fname} L{hits[0][0]}] {hits[0][1]}")
    return hits[0]

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

def num(s):
    m = re.search(r"-?\d+(?:,\d{3})*(?:\.\d+)?", s.replace("−", "-"))
    return float(m.group(0).replace(",", ""))

effects = json.load(open(os.path.join(D, "effects.json")))
H = {}

# ---- E1: log-point ATT -> percent change
DCR = "arxiv-2609.27732v1.txt"
ln_hdr, hdr = find(DCR, "| 1 km | 3 km | 5 km")
ln_att, att = find(DCR, "ATT | -0.399")
ln_pct, pct = find(DCR, "% change | -32.9%")
c = cells(att); b = num(c[3])
e1 = 100 * (math.exp(b) - 1)
printed_pct = num(cells(pct)[3])
# consistency check of printed % change row for all columns
for k in (1, 2, 3):
    bb = num(cells(att)[k]); pp = num(cells(pct)[k])
    print(f"  check col {k}: 100(exp({bb})-1) = {100*(math.exp(bb)-1):.2f}% vs printed {pp}%")
H["E1"] = dict(value=round(e1, 2), unit="percent change in the drug-related hospitalisation rate relative to the counterfactual",
    transformation=f"100*(exp(b)-1) with b = {c[3]} (log points) -> 100*(exp(-0.431)-1) = {e1:.4f}%; agrees with printed '% change' {cells(pct)[3]}",
    harmonisable=True, missing_quantity=None,
    locator=f"arxiv-2609.27732v1 | Table 6 | row 'ATT' (line {ln_att}) | column '5 km' (value -0.431**); cross-check row '% change' (line {ln_pct}) column '5 km' (-35.0%)")

# ---- E2: pp effect -> relative percent of control mean
RE = "arxiv-2606.14621v1.txt"
find(RE, "Outcome | Treatment | Control |")
ln_r, r = find(RE, "Reemployment | 75.2%")
c = cells(r)
p_ctrl = num(c[2]); delta = num(c[4])
e2 = 100 * delta / p_ctrl
H["E2"] = dict(value=round(e2, 2), unit="percent of the control-group mean reemployment rate",
    transformation=f"100*Delta/p with Delta = {delta} pp (ATE) and p = {p_ctrl}% (control mean) -> 100*1.47/75.0 = {e2:.4f}%",
    harmonisable=True, missing_quantity=None,
    locator=f"arxiv-2606.14621v1 | Table 2 | row 'Reemployment' (line {ln_r}) | column 'Average treatment effect estimate' (1.47 pp) and column 'Control' (75.0%)")

# ---- E3: elasticity x printed % MW change
MW = "arxiv-2505.04555v3.txt"
ln_el, el = find(MW, "Employment elasticity w.r.t MW")
ln_g, gl = find(MW, "% MW changes")
find(MW, "Note: The table reports the effects of a minimum wage increase")
eta = num(cells(el)[1]); g = num(cells(gl)[1])
e3 = 100 * eta * g
H["E3"] = dict(value=round(e3, 2), unit="percent change in employment",
    transformation=f"100*eta*g with eta = {eta} and g = {g} (printed '% MW changes', a proportion) -> 100*(-0.388)*(0.047) = {e3:.4f}% (first-order/linear approximation)",
    harmonisable=True, missing_quantity=None,
    locator=f"arxiv-2505.04555v3 | Table 2 | row 'Employment elasticity w.r.t MW' (line {ln_el}) | column 'Estimate (SE)' (-0.388); row '% MW changes' (line {ln_g}) | column 'Estimate (SE)' (0.047)")

# ---- E4: log-point ATT -> absolute change per 100,000: needs baseline level
ln_f, fl = find(DCR, "Drug-related (fatal) |")
bf = num(cells(fl)[3])
rel4 = 100 * (math.exp(bf) - 1)
print(f"  E4 relative change (not the target): 100(exp({bf})-1) = {rel4:.2f}%")
H["E4"] = dict(value=None, unit="hospitalisations per 100,000 inhabitants per year",
    transformation=f"Absolute change = Y0*(exp(b)-1) with b = {cells(fl)[3]}; the relative part (exp(-0.354)-1 = {rel4/100:.4f}) is computable, but the counterfactual/baseline level Y0 of the fatal drug-related hospitalisation rate is not printed in the extract, so the absolute change cannot be computed",
    harmonisable=False,
    missing_quantity="baseline (pre-treatment or counterfactual) mean fatal drug-related hospitalisation rate per 100,000 inhabitants per year for the treated regions within 5 km (Y0), which the extract does not print",
    locator=f"arxiv-2609.27732v1 | Table 7 | row 'Drug-related (fatal)' (line {ln_f}) | column '5 km' (-0.354**)")

# ---- E5: LPM coefficient (0-1 scale) -> pp
BE = "arxiv-2510.03658v1.txt"
ln_t, tl = find(BE, "Treatment | -0.0911")
ln_m, ml = find(BE, "Mean Dep. Var.")
beta = num(cells(tl)[4])
e5 = 100 * beta
H["E5"] = dict(value=round(e5, 2), unit="percentage points",
    transformation=f"100*beta with beta = {cells(tl)[4]} (probability on 0-1 scale) -> {e5:.2f} pp (exact rescaling)",
    harmonisable=True, missing_quantity=None,
    locator=f"arxiv-2510.03658v1 | Table 3 | row 'Treatment' (line {ln_t}) | column '(4)' (-0.0962)")
print(f"  E5 cross-check: -0.0962/0.261 = {beta/num(cells(ml)[4]):.4f} vs printed % Change -0.3685")

data_issues = [
 {"issue": "arxiv-2505.04555v3 Table 2 note describes the 2023 minimum-wage changes as 'about 4% minimum wage increase', whereas the table row '% MW changes' prints 0.047 (4.7%).",
  "handling": f"Used the table value 0.047 as instructed by the target definition ('size of that increase as printed in the same table'); with 4% the implied change would be {100*eta*0.04:.2f}% instead of {e3:.2f}%."},
 {"issue": "arxiv-2505.04555v3 abstract reports 'a 2% employment decline in affected bins', which is a different quantity (affected-bin employment) and not identical to eta*g = -1.82%.",
  "handling": "Not used as an input; noted only for context."},
 {"issue": f"arxiv-2609.27732v1: the abstract says fatal outcomes 'fall by 35%', Section 6.1 says 'approximately 30% at 5 km', while the Table 7 coefficient -0.354 implies 100(exp(-0.354)-1) = {rel4:.1f}%. The abstract's 35% matches the overall ATT (Table 6), not the fatal sub-outcome.",
  "handling": "Relied only on the Table 7 coefficient; the text percentages were not used. E4 is not harmonisable regardless (absolute target)."},
 {"issue": "arxiv-2609.27732v1 Table A19 prints SE (0.218) for the CS-DR ATT -0.431 whereas Table 6 prints (0.210) for the same estimate; the A19 note mentions 640 never-treated comparison regions while the A19 table and Tables 6-7 state 608 control regions.",
  "handling": "Point estimates agree (-0.431, -35.0%); E1 uses Table 6 as specified. SEs are not needed (no uncertainty propagation)."},
 {"issue": "arxiv-2609.27732v1 Table 7 note lists only four sub-outcome definitions; the 'Drug-related (fatal)' category is not defined in the note.",
  "handling": "Recorded; does not affect the computation."},
 {"issue": "arxiv-2606.14621v1 Table 2: the 75.0% control mean is the raw aggregate control-group mean, whereas the 1.47 pp ATE adjusts for randomisation week (blocks); the raw treatment-control difference is only 0.26 pp. The control mean is therefore only an approximate counterfactual level for the block-adjusted ATE.",
  "handling": "Used the printed control-group mean as the target definition requires; flagged as approximation."},
 {"issue": "arxiv-2510.03658v1 text refers to 'table 2' (estimate 0.0855, mean 0.455), which is omitted from the extract; footnote 11 is garbled in rendering (three percentages listed against five wealth categories); Table 3 prints mean dep. var. 0.261 and % change -0.3685 consistent with -0.0962/0.261.",
  "handling": "E5 uses only Table 3 column (4); the Table-2-based numbers in the abstract were not used."},
 {"issue": "arxiv-2606.14621v1 text contains rendering artefact 'p ¡ 0.001'.",
  "handling": "Cosmetic; no effect."},
]

def nlines(p):
    with open(p, "rb") as f:
        return f.read().count(b"\n")

n_input = {"effects.json": len(effects), "README.md": nlines(os.path.join(D, "README.md"))}
for fn in sorted(os.listdir(S)):
    n_input[fn] = nlines(os.path.join(S, fn))

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