"""A1 · IMPLEMENTACIÓN 1 (DuckDB + Python) del cálculo sellado en DEFINICIONES_PREVIAS.md §3-§6. Sin red, sin modelos.
Entradas: priv/a1_base.parquet, priv/lectura_final.csv (puede no existir o estar vacío), salidas/15_ref_mensual.csv,
salidas/15_factor_trimestral.csv, et/lotes_cartel_v2.csv, et/banco_entidades.parquet, /d identidad (organos.parquet).
Salidas: priv/30_contratos_calculo.csv, salidas/30_razones_producto.csv, salidas/31_hechos_organo.csv (interno),
salidas/publicable_organo_anio.csv, salidas/32_validaciones.json, salidas/33_mensual.csv."""
import duckdb, json, math, random, csv, os, re, statistics as st, sys
P = '/p'; SEED_BOOT = 20261021; SEED_VAL = 20261022
PRODS = ('gasoleo_a', 'gasoleo_b', 'gasoleo_c', 'gasolina95')
con = duckdb.connect(); q = lambda s: con.execute(s).fetchall()
con.execute(f"CREATE TABLE ref AS SELECT * FROM read_csv('{P}/salidas/15_ref_mensual.csv', header=true, columns={{'mes':'VARCHAR','producto':'VARCHAR','semanas_wob':'INT','ref_nacional_sin_iva_eur_l':'DOUBLE'}})")
con.execute(f"CREATE TABLE fac AS SELECT * FROM read_csv('{P}/salidas/15_factor_trimestral.csv', header=true, columns={{'anio':'INT','trimestre':'INT','dia':'VARCHAR','prov_ine':'VARCHAR','producto':'VARCHAR','estaciones':'INT','factor':'DOUBLE'}})")
LF = f'{P}/priv/lectura_final.csv'
if os.path.exists(LF) and os.path.getsize(LF) > 0:
    con.execute(f"CREATE TABLE lf AS SELECT * FROM read_csv('{LF}', header=true, all_varchar=true)")
else:
    con.execute("CREATE TABLE lf (source VARCHAR, id_canonical VARCHAR, licitacion_id VARCHAR, caso VARCHAR, lote VARCHAR, via VARCHAR, producto VARCHAR, cantidad VARCHAR, unidad VARCHAR, periodo VARCHAR, duracion_leida VARCHAR, importe_tipo VARCHAR, duda VARCHAR)")
# §3.1 Universo y §3.2 cantidad
con.execute(f"""CREATE TABLE u AS SELECT a.*, coalesce(a.fecha_adj, a.fecha_formalizacion, a.fecha_pub) f0
 FROM read_parquet('{P}/priv/a1_base.parquet') a
 WHERE mercado = 'carburantes' AND NOT es_copia AND NOT excluido_v2 AND imp > 0 AND coalesce(rol_am, '') NOT IN ('acuerdo', 'sda')
   AND coalesce(prov_ine, '') NOT IN ('35', '38', '51', '52') AND NOT regexp_matches(coalesce(nuts, ''), '^(ES70|ES63|ES64)')""")
con.execute("""CREATE TABLE c AS SELECT u.source, u.id_canonical, u.organo_clave, u.organo_nombre, u.tipo_organo, u.anio_ref, u.prov_ine, u.imp,
   u.menor, u.licitacion_id,
   CASE WHEN f0 BETWEEN DATE '2022-01-01' AND DATE '2025-12-31' THEN strftime(f0, '%Y-%m') ELSE anio_ref::VARCHAR || '-07' END mes,
   CASE WHEN l.id_canonical IS NOT NULL AND TRY_CAST(l.cantidad AS DOUBLE) IS NOT NULL THEN 'L'
        WHEN u.cant_v2 IS NOT NULL THEN 'R' END via,
   CASE WHEN l.id_canonical IS NOT NULL AND TRY_CAST(l.cantidad AS DOUBLE) IS NOT NULL THEN l.producto ELSE u.producto_v2 END producto,
   l.unidad l_unidad, l.periodo l_periodo, l.importe_tipo, l.duda l_duda, TRY_CAST(l.cantidad AS DOUBLE) l_cant,
   coalesce(u.duracion_meses, TRY_CAST(l.duracion_leida AS DOUBLE)) dur, u.cant_v2, u.n_cant_v2
 FROM u LEFT JOIN lf l ON l.source = u.source AND l.id_canonical = u.id_canonical""")
con.execute("""CREATE TABLE c2 AS SELECT *,
   CASE WHEN via = 'R' THEN cant_v2
        WHEN via = 'L' AND l_unidad = 'litros' AND l_periodo = 'total_contrato' THEN l_cant
        WHEN via = 'L' AND l_unidad = 'litros' AND l_periodo = 'anual' AND dur > 0 THEN l_cant * dur / 12
        WHEN via = 'L' AND l_unidad = 'litros' AND l_periodo = 'mensual' AND dur > 0 THEN l_cant * dur END cantidad_total,
   CASE WHEN via IS NULL THEN 'sin_cantidad'
        WHEN producto NOT IN ('gasoleo_a','gasoleo_b','gasoleo_c','gasolina95') THEN 'sin_producto_comparable'
        WHEN via = 'L' AND coalesce(l_unidad, '') <> 'litros' THEN 'unidad_no_litros'
        WHEN via = 'L' AND importe_tipo = 'descuento_surtidor' THEN 'descuento_surtidor'
        WHEN via = 'L' AND l_periodo NOT IN ('total_contrato', 'anual', 'mensual') THEN 'periodo_desconocido'
        WHEN via = 'L' AND l_periodo IN ('anual', 'mensual') AND NOT (dur > 0) THEN 'sin_duracion'
        ELSE 'comparable' END estado
 FROM c""")
con.execute("""CREATE TABLE c3 AS SELECT c2.*, r.ref_nacional_sin_iva_eur_l *
   CASE WHEN c2.producto IN ('gasoleo_a', 'gasolina95') THEN coalesce(f.factor, 1.0) ELSE 1.0 END ref_eur_l,
   (c2.producto IN ('gasoleo_a', 'gasolina95') AND f.factor IS NOT NULL) con_factor
 FROM c2 LEFT JOIN ref r ON r.mes = c2.mes AND r.producto = c2.producto
 LEFT JOIN fac f ON f.prov_ine = c2.prov_ine AND f.producto = c2.producto AND f.anio = TRY_CAST(substr(c2.mes, 1, 4) AS INT)
      AND f.trimestre = (TRY_CAST(substr(c2.mes, 6, 2) AS INT) - 1) // 3 + 1""")
con.execute("""CREATE TABLE k AS SELECT *, imp / cantidad_total pu, imp / cantidad_total / ref_eur_l razon,
   estado = 'comparable' AND cantidad_total >= 200 AND (via = 'L' OR n_cant_v2 = 1) nucleo,
   estado = 'comparable' AND cantidad_total >= 200 AND (via = 'L' OR n_cant_v2 = 1)
     AND imp / cantidad_total / ref_eur_l BETWEEN 0.5 AND 2.0 valido
 FROM c3""")
con.execute(f"COPY (SELECT source, id_canonical, organo_clave, anio_ref, mes, prov_ine, via, producto, estado, cantidad_total, imp, ref_eur_l, razon, nucleo, valido FROM k ORDER BY source, id_canonical) TO '{P}/priv/30_contratos_calculo.csv' (HEADER)")
# §3.4 Intervalo esperado (dejando fuera al órgano) y K de N
con.execute("CREATE TABLE v AS SELECT * FROM k WHERE valido AND organo_clave IS NOT NULL")
con.execute("CREATE TABLE un AS SELECT organo_clave, producto, count(*) n FROM v GROUP BY 1, 2")
con.execute("""CREATE TABLE iv AS SELECT un.organo_clave, un.producto, un.n,
   quantile_cont(v.razon, 0.10) q10, quantile_cont(v.razon, 0.90) q90
 FROM un JOIN v ON v.producto = un.producto AND v.organo_clave <> un.organo_clave GROUP BY 1, 2, 3""")
con.execute("""CREATE TABLE kn AS SELECT iv.*, sum((v.razon > iv.q90)::int) k_sup, sum((v.razon < iv.q10)::int) k_inf, median(v.razon) mediana
 FROM iv JOIN v USING (organo_clave, producto) GROUP BY ALL""")
def cola(n, k, p=0.10):
    return sum(math.comb(n, i) * p ** i * (1 - p) ** (n - i) for i in range(k, n + 1))
filas = [dict(zip([d[0] for d in con.execute("SELECT * FROM kn").description], r)) for r in q("SELECT * FROM kn ORDER BY organo_clave, producto")]
ev = [f for f in filas if f['n'] >= 10]
pruebas = []
for f in ev:
    f['p_sup'] = cola(f['n'], f['k_sup']); f['p_inf'] = cola(f['n'], f['k_inf'])
    pruebas += [(f['p_sup'], id(f), 'sup'), (f['p_inf'], id(f), 'inf')]
pruebas.sort(); m = len(pruebas); rech = set(); seguir = True
for i, (p, fid, lado) in enumerate(pruebas):           # Holm-Bonferroni, alfa 0,05 en todas las pruebas órgano × producto × lado
    if seguir and p <= 0.05 / (m - i): rech.add((fid, lado))
    else: seguir = False
rnd = random.Random(SEED_BOOT)
vals = {}
for o, p, r in q("SELECT organo_clave, producto, razon FROM v ORDER BY organo_clave, producto, source, id_canonical"): vals.setdefault((o, p), []).append(r)
for f in filas:
    f['evaluable'] = f['n'] >= 10
    f['hecho'] = ('por_encima' if (id(f), 'sup') in rech else 'por_debajo' if (id(f), 'inf') in rech else 'ninguno') if f['evaluable'] else 'no_evaluable'
    if f['evaluable']:
        xs = vals[(f['organo_clave'], f['producto'])]; b = sorted(st.median(rnd.choices(xs, k=len(xs))) for _ in range(10000))
        f['ic_lo'], f['ic_hi'] = b[249], b[9749]
    else: f['p_sup'] = f['p_inf'] = f['ic_lo'] = f['ic_hi'] = None
PF = re.compile(r'^([0-9]{8}|[XYZKLM][0-9]{7})[A-Z]$')
nombre = dict(q("SELECT organo_clave, any_value(organo_nombre) FROM k WHERE organo_clave IS NOT NULL GROUP BY 1"))
for f in filas: f['persona_fisica'] = bool(PF.match((f['organo_clave'] or '').split('|')[0]))
cols = ['organo_clave', 'producto', 'n', 'q10', 'q90', 'k_sup', 'k_inf', 'p_sup', 'p_inf', 'mediana', 'ic_lo', 'ic_hi', 'evaluable', 'hecho', 'persona_fisica']
with open(f'{P}/salidas/31_hechos_organo.csv', 'w', newline='') as fh:
    w = csv.writer(fh); w.writerow(cols + ['organo_nombre'])
    for f in filas: w.writerow([round(f[c], 6) if isinstance(f[c], float) else f[c] for c in cols] + [nombre.get(f['organo_clave'])])
hecho = {(f['organo_clave'], f['producto']): f['hecho'] for f in filas}
with open(f'{P}/salidas/publicable_organo_anio.csv', 'w', newline='') as fh:
    w = csv.writer(fh); w.writerow(['organo', 'tipo_organo', 'anio', 'producto', 'contratos_validos', 'por_encima_p90', 'por_debajo_p10', 'hecho_2022_2025'])
    for o, t, a, p, n, ks, ki in q("""SELECT v.organo_clave, any_value(v.tipo_organo), v.anio_ref, v.producto, count(*),
           sum((v.razon > kn.q90)::int), sum((v.razon < kn.q10)::int) FROM v JOIN kn USING (organo_clave, producto) GROUP BY 1, 3, 4 ORDER BY 1, 3, 4"""):
        if PF.match((o or '').split('|')[0]): continue
        w.writerow([nombre.get(o), t, a, p, n, ks, ki, hecho.get((o, p))])
# Distribuciones
def dist(xs):
    xs = sorted(xs); n = len(xs)
    if not n: return dict(n=0)
    qq = lambda p: xs[min(n - 1, int(p * n))]
    return dict(n=n, mediana=round(st.median(xs), 4), p10=round(qq(.1), 4), p25=round(qq(.25), 4), p75=round(qq(.75), 4), p90=round(qq(.9), 4),
                dentro_15=round(100 * sum(.85 <= x <= 1.15 for x in xs) / n, 1))
tab = []
for grupo, cond in [('valido', 'valido'), ('nucleo', 'nucleo AND razon IS NOT NULL')]:
    for p in PRODS + ('todos',):
        for via in ('R', 'L', 'todas'):
            xs = [r[0] for r in q(f"SELECT razon FROM k WHERE {cond} {'' if p == 'todos' else f"AND producto='{p}'"} {'' if via == 'todas' else f"AND via='{via}'"}")]
            tab.append(dict(grupo=grupo, producto=p, via=via, **dist(xs)))
with open(f'{P}/salidas/30_razones_producto.csv', 'w', newline='') as fh:
    w = csv.DictWriter(fh, fieldnames=['grupo', 'producto', 'via', 'n', 'mediana', 'p10', 'p25', 'p75', 'p90', 'dentro_15']); w.writeheader(); w.writerows(tab)
with open(f'{P}/salidas/33_mensual.csv', 'w', newline='') as fh:
    w = csv.writer(fh); w.writerow(['mes', 'producto', 'n', 'mediana_razon', 'p25', 'p75'])
    for r in q("SELECT mes, producto, count(*), median(razon), quantile_cont(razon, .25), quantile_cont(razon, .75) FROM k WHERE valido GROUP BY 1, 2 ORDER BY 1, 2"): w.writerow([r[0], r[1], r[2]] + [round(x, 4) for x in r[3:]])
# §5 Validaciones
V = {}
nuc = [r[0] for r in q("SELECT razon FROM k WHERE nucleo AND razon IS NOT NULL")]
V['mediana_mercado_nucleo'] = dict(valor=round(st.median(nuc), 4) if nuc else None, n=len(nuc), pasa=bool(nuc) and 0.85 <= st.median(nuc) <= 1.15)
cc = [r[0] for r in q("SELECT razon FROM k WHERE nucleo AND producto='gasoleo_c' AND razon IS NOT NULL")]
V['control_gasoleo_c_dentro_15'] = dict(valor=round(100 * sum(.85 <= x <= 1.15 for x in cc) / len(cc), 1) if cc else None, n=len(cc), umbral=75)
V['control_gasoleo_c_dentro_15']['pasa'] = (V['control_gasoleo_c_dentro_15']['valor'] or 0) >= 75
V['control_nucleo_dentro_15'] = dict(valor=round(100 * sum(.85 <= x <= 1.15 for x in nuc) / len(nuc), 1) if nuc else None, n=len(nuc), umbral=50)
V['control_nucleo_dentro_15']['pasa'] = (V['control_nucleo_dentro_15']['valor'] or 0) >= 50
# Estabilidad (Spearman de la mediana por órgano × producto, 2022-2023 frente a 2024-2025, ≥ 5 válidos en cada uno)
par = q("""SELECT organo_clave, producto, median(razon) FILTER (WHERE anio_ref <= 2023), median(razon) FILTER (WHERE anio_ref >= 2024)
  FROM v GROUP BY 1, 2 HAVING count(*) FILTER (WHERE anio_ref <= 2023) >= 5 AND count(*) FILTER (WHERE anio_ref >= 2024) >= 5 ORDER BY 1, 2""")
def rangos(xs):
    o = sorted(range(len(xs)), key=lambda i: xs[i]); r = [0.0] * len(xs); i = 0
    while i < len(o):
        j = i
        while j + 1 < len(o) and xs[o[j + 1]] == xs[o[i]]: j += 1
        for t in range(i, j + 1): r[o[t]] = (i + j) / 2 + 1
        i = j + 1
    return r
def spearman(a, b):
    if len(a) < 3: return None
    ra, rb = rangos(a), rangos(b); ma, mb = st.mean(ra), st.mean(rb)
    num = sum((x - ma) * (y - mb) for x, y in zip(ra, rb)); den = math.sqrt(sum((x - ma) ** 2 for x in ra) * sum((y - mb) ** 2 for y in rb))
    return num / den if den else None
a = [r[2] for r in par]; b = [r[3] for r in par]; rho = spearman(a, b); rv = random.Random(SEED_VAL); bs = []
for _ in range(2000):
    idx = [rv.randrange(len(par)) for _ in par] if par else []
    s = spearman([a[i] for i in idx], [b[i] for i in idx]) if idx else None
    if s is not None: bs.append(s)
bs.sort()
ic = [round(bs[int(.025 * len(bs))], 3), round(bs[int(.975 * len(bs)) - 1], 3)] if len(bs) > 40 else None
V['estabilidad'] = dict(unidades=len(par), spearman=round(rho, 3) if rho is not None else None, ic95=ic,
                        resultado=('no_evaluable' if len(par) < 15 else 'estable' if (rho or 0) >= 0.30 and ic and ic[0] > 0 else 'no_estable'))
# Etiquetas (solo si evaluables: ≥ 10 positivos y ≥ 10 negativos)
def auc(pos, neg):
    if not pos or not neg: return None
    s = sum((x > y) + 0.5 * (x == y) for x in pos for y in neg); return s / len(pos) / len(neg)
cart = q(f"""SELECT v.razon, c.id_md5 IS NOT NULL, v.organo_clave FROM v LEFT JOIN (SELECT DISTINCT id_md5 FROM read_csv_auto('{P}/et/lotes_cartel_v2.csv') WHERE var IN ('V2','V3') AND source='PLACSP') c
  ON c.id_md5 = md5(v.source || '|' || v.id_canonical) WHERE v.source='PLACSP'""")
pos = [r[0] for r in cart if r[1]]; neg = [r[0] for r in cart if not r[1]]
V['carteles_lote'] = dict(positivos=len(pos), negativos=len(neg), auc=None, ic95=None, resultado='no_evaluable')
if len(pos) >= 10 and len(neg) >= 10:
    V['carteles_lote']['auc'] = round(auc(pos, neg), 3); V['carteles_lote']['resultado'] = 'calcular_ic'   # (no se alcanza: 0 positivos medidos antes de sellar)
p6 = q(f"""SELECT b.entidad_clave, b.anio, b.etiqueta, median(v.razon), count(*) FROM read_parquet('{P}/et/banco_entidades.parquet') b
  JOIN (SELECT organo_clave, entidad_clave FROM read_parquet('/d/unified_20261018.identidad/organos.parquet')) o ON o.entidad_clave = b.entidad_clave
  JOIN v ON v.organo_clave = o.organo_clave AND v.anio_ref = b.anio WHERE b.etiqueta IN ('positivo','negativo') GROUP BY 1, 2, 3 HAVING count(*) >= 3""")
pp = [r[3] for r in p6 if r[2] == 'positivo']; nn = [r[3] for r in p6 if r[2] == 'negativo']
V['fiscalizacion_p6'] = dict(positivos=len(pp), negativos=len(nn), auc=None, ic95=None, resultado='no_evaluable')
if len(pp) >= 10 and len(nn) >= 10:
    V['fiscalizacion_p6']['auc'] = round(auc(pp, nn), 3); V['fiscalizacion_p6']['resultado'] = 'calcular_ic'
# Cobertura
V['cobertura'] = dict(zip(['universo', 'con_cantidad_R', 'con_cantidad_L', 'comparables', 'nucleo', 'validos', 'importe_universo_meur', 'importe_validos_meur'],
    q("""SELECT count(*), sum((via='R')::int), sum((via='L')::int), sum((estado='comparable')::int), sum(nucleo::int), sum(valido::int),
         round(sum(imp)/1e6, 2), round(sum(imp) FILTER (WHERE valido)/1e6, 2) FROM k""")[0]))
V['estados'] = dict(q("SELECT estado, count(*) FROM k GROUP BY 1 ORDER BY 2 DESC"))
V['unidades_organo_producto'] = dict(total=len(filas), evaluables=len(ev), pruebas_holm=m,
    por_encima=sum(f['hecho'] == 'por_encima' for f in filas), por_debajo=sum(f['hecho'] == 'por_debajo' for f in filas),
    persona_fisica_excluidas=sum(f['persona_fisica'] for f in filas))
json.dump(V, open(f'{P}/salidas/32_validaciones.json', 'w'), indent=1, ensure_ascii=False, default=str)
print(json.dumps(V, ensure_ascii=False, default=str))