import json
import math
import os
from pathlib import Path

ROOT = Path(".")

INPUT_FIRST = ROOT / "results" / "row20_h0_attribution" / "row20_NZ_CFGL_bridge_operator.json"

OUTPUT_JSON = ROOT / "results" / "row20_h0_attribution" / "row20_BAO_RT_coupling_operator.json"
OUTPUT_REPORT = ROOT / "results" / "row20_h0_attribution" / "row20_BAO_RT_coupling_operator.md"


def load(path):
    with open(path, "r", encoding="utf-8") as f:
        return json.load(f)


def rss(vals):
    return math.sqrt(sum(v * v for v in vals))


def chi2(vals):
    return sum(v * v for v in vals)


def main():
    data = load(INPUT_FIRST)

    r = data["transformed_residual_vector"]

    r_H = r["r_H_BAO"]
    r_DM = r["r_D_M_BAO"]

    # Derive theta from BAO structure (NOT free)
    numerator = (r_H - r_DM)
    denominator = (r_H + r_DM) if (r_H + r_DM) != 0 else 1e-9

    theta = 0.5 * math.atan(numerator / denominator)

    cos_t = math.cos(theta)
    sin_t = math.sin(theta)

    # Apply coupling rotation
    r_H_new = cos_t * r_H - sin_t * r_DM
    r_DM_new = sin_t * r_H + cos_t * r_DM

    transformed = {
        "r_ell_A": r["r_ell_A"],
        "r_H_BAO": r_H_new,
        "r_D_M_BAO": r_DM_new,
        "r_SN": r["r_SN"],
        "r_SH0ES": r["r_SH0ES"],
    }

    original_vals = list(r.values())
    transformed_vals = list(transformed.values())

    original_rss = rss(original_vals)
    transformed_rss = rss(transformed_vals)

    original_chi2 = chi2(original_vals)
    transformed_chi2 = chi2(transformed_vals)

    rss_improvement = original_rss - transformed_rss
    chi2_improvement = original_chi2 - transformed_chi2

    result = {
        "signature": "ROW20_BAO_RT_COUPLING_OPERATOR_v1",
        "theta": theta,
        "cos_theta": cos_t,
        "sin_theta": sin_t,
        "original_residuals": r,
        "transformed_residuals": transformed,
        "original_rss": original_rss,
        "transformed_rss": transformed_rss,
        "rss_improvement": rss_improvement,
        "original_chi2": original_chi2,
        "transformed_chi2": transformed_chi2,
        "chi2_improvement": chi2_improvement,
        "verdict": "diagnostic only; no physics or parameter changes"
    }

    OUTPUT_JSON.parent.mkdir(parents=True, exist_ok=True)

    with open(OUTPUT_JSON, "w", encoding="utf-8") as f:
        json.dump(result, f, indent=2)

    with open(OUTPUT_REPORT, "w", encoding="utf-8") as f:
        f.write("# BAO Radial/Transverse Coupling Operator\n\n")
        f.write(f"theta = {theta:.6f}\n\n")
        f.write(f"RSS before = {original_rss:.6f}\n")
        f.write(f"RSS after  = {transformed_rss:.6f}\n")
        f.write(f"improvement = {rss_improvement:.6f}\n\n")

    print("BAO RT coupling:")
    print(f"theta={theta:.6f}")
    print(f"RSS before={original_rss:.6f}")
    print(f"RSS after={transformed_rss:.6f}")
    print(f"improvement={rss_improvement:.6f}")


if __name__ == "__main__":
    main()