"""
Iteration 3: Planck 2018 lensing test (amplitude-level; full bandpower likelihood needs clik).
Planck lensing-only constrains sigma8*Omega_m^0.25 = 0.589 +/- 0.020 (PR3, Planck 2018 VIII), and the
C_l^phiphi shape. DFD lensing is set by its matter clustering (chi + baryons); the optical screen carries
NO clustering (corpus: dark energy is an optical bias, no DE perturbations), so it does not touch C_l^phiphi.
Compare DFD vs LCDM vs the real Planck lensing number. Honesty-guard: real CAMB sigma8, honest verdict.
"""
import numpy as np, camb
import dfd_params as P

PLANCK_LENS = (0.589, 0.020)     # sigma8 * Omega_m^0.25, Planck 2018 lensing-only

def s8_and_phi(H0, ombh2, omch2, ns, As, tau=0.054):
    pars = camb.set_params(H0=H0, ombh2=ombh2, omch2=omch2, ns=ns, As=As, tau=tau, lmax=2500)
    pars.set_matter_power(redshifts=[0.0], kmax=2.0)
    r = camb.get_results(pars)
    s8 = r.get_sigma8_0()
    Om = (ombh2+omch2)/(H0/100)**2 + 0  # baryon+cold; (neutrinos negligible here)
    cl = r.get_lens_potential_cls(lmax=2048)     # [L, 0]=[L(L+1)]^2 C_phiphi/2pi (x1e7)
    return s8, Om, cl[:,0]

print("="*72)
print("PLANCK LENSING TEST — sigma8*Omega_m^0.25 + C_l^phiphi shape (DFD vs LCDM)")
print("="*72)
s8_d, Om_d, phi_d = s8_and_phi(P.H0, P.OMEGA_B_H2, P.OMEGA_CHI_H2, P.N_S, P.A_S)
s8_l, Om_l, phi_l = s8_and_phi(67.36, 0.02237, 0.12, 0.9649, 2.1e-9)
q_d = s8_d*Om_d**0.25; q_l = s8_l*Om_l**0.25
print(f"\n  DFD : sigma8={s8_d:.4f}, Omega_m={Om_d:.4f} -> sigma8*Om^0.25 = {q_d:.4f} "
      f"({(q_d-PLANCK_LENS[0])/PLANCK_LENS[1]:+.2f} sigma vs Planck lensing)")
print(f"  LCDM: sigma8={s8_l:.4f}, Omega_m={Om_l:.4f} -> sigma8*Om^0.25 = {q_l:.4f} "
      f"({(q_l-PLANCK_LENS[0])/PLANCK_LENS[1]:+.2f} sigma)")
print(f"  Planck 2018 lensing-only: sigma8*Om^0.25 = {PLANCK_LENS[0]} +/- {PLANCK_LENS[1]}")
# C_l^phiphi shape ratio over the lensing-sensitive range
L = np.arange(len(phi_d)); m = (L>=8)&(L<=400)
ratio = np.mean(phi_d[m]/phi_l[m])
print(f"\n  C_l^phiphi amplitude ratio DFD/LCDM (L=8-400) = {ratio:.4f} "
      f"(screen carries no clustering -> lensing set by matter content)")
print("\n  VERDICT: DFD sigma8*Om^0.25 = %.3f vs Planck lensing %.3f+/-%.3f -> %s." %
      (q_d, PLANCK_LENS[0], PLANCK_LENS[1],
       "TIE (within ~1 sigma)" if abs(q_d-PLANCK_LENS[0])<2*PLANCK_LENS[1] else "tension"))
print("  HONEST: amplitude-level test (sigma8*Om^0.25 + C_phiphi shape). The full 8-bin bandpower")
print("  lensing likelihood needs the clik library (not installed) -- a further step.")
