"""E8 · SEGUNDA IMPLEMENTACIÓN del cálculo sellado (Python puro + pyarrow solo para leer el parquet; sin DuckDB, sin reglas_e8, sin reglas_v2 ni comun_a1 de A1). Escrita aparte a partir de DEFINICIONES_PREVIAS_E8.md: rehace el universo, las reglas (impl2_reglas), las REFERENCIAS desde los ficheros crudos (geoportal, WOB, Eurostat, ofertas de los AM), la razón, el estado, «por encima» y los hechos por órgano; y compara con la implementación 1. Uso: impl2_calculo.py [base.parquet] [carpeta_salida] [crudos_ref] [contratos_impl1.csv] [organos_impl1.csv] [salida.json]""" import sys, os, csv, json, gzip, glob, math, re, zipfile, collections, datetime as dt import xml.etree.ElementTree as ET sys.path.insert(0, os.path.dirname(os.path.abspath(__file__)) if os.path.isdir(os.path.dirname(os.path.abspath(__file__))) else '/p/scripts') sys.path.insert(0, '/p/scripts') import pyarrow.parquet as pq import impl2_reglas as R2 A = sys.argv[1:] + [None] * 6 BASE = A[0] or '/p/priv/e8_base.parquet'; OUT = A[1] or '/p'; CRU = A[2] or '/ref' C1 = A[3] or f'{OUT}/priv/e8_contratos_impl1.csv'; O1 = A[4] or f'{OUT}/salidas/e8_organos_impl1.csv' SAL = A[5] or f'{OUT}/salidas/e8_impl2.json' PAPEL_CSV = os.environ.get('E8_REF_PAPEL', '/p/salidas/ref_papel_am.csv') # ---------------------------------------------------------------- provincias (tabla propia: NUTS3 → INE y nombre → INE) NUTS = dict(ES111='15', ES112='27', ES113='32', ES114='36', ES120='33', ES130='39', ES211='01', ES212='20', ES213='48', ES220='31', ES230='26', ES241='22', ES242='44', ES243='50', ES300='28', ES411='05', ES412='09', ES413='24', ES414='34', ES415='37', ES416='40', ES417='42', ES418='47', ES419='49', ES421='02', ES422='13', ES423='16', ES424='19', ES425='45', ES431='06', ES432='10', ES511='08', ES512='17', ES513='25', ES514='43', ES521='03', ES522='12', ES523='46', ES531='07', ES532='07', ES533='07', ES611='04', ES612='11', ES613='14', ES614='18', ES615='21', ES616='23', ES617='29', ES618='41', ES620='30', ES630='51', ES640='52', ES703='38', ES704='35', ES705='35', ES706='38', ES707='38', ES708='35', ES709='38') NOMBRES = {'01': 'alava araba araba/alava', '02': 'albacete', '03': 'alicante alacant alicante/alacant', '04': 'almeria', '05': 'avila', '06': 'badajoz', '07': 'baleares|illes balears|balears|islas baleares', '08': 'barcelona', '09': 'burgos', '10': 'caceres', '11': 'cadiz', '12': 'castellon castello castellon/castello', '13': 'ciudad real', '14': 'cordoba', '15': 'a coruna|coruna|la coruna|coruna, a', '16': 'cuenca', '17': 'girona gerona', '18': 'granada', '19': 'guadalajara', '20': 'gipuzkoa guipuzcoa', '21': 'huelva', '22': 'huesca', '23': 'jaen', '24': 'leon', '25': 'lleida lerida', '26': 'la rioja|rioja|rioja, la', '27': 'lugo', '28': 'madrid', '29': 'malaga', '30': 'murcia', '31': 'navarra nafarroa', '32': 'ourense orense', '33': 'asturias', '34': 'palencia', '35': 'las palmas|palmas, las', '36': 'pontevedra', '37': 'salamanca', '38': 'santa cruz de tenerife', '39': 'cantabria', '40': 'segovia', '41': 'sevilla', '42': 'soria', '43': 'tarragona', '44': 'teruel', '45': 'toledo', '46': 'valencia valencia/valencia', '47': 'valladolid', '48': 'bizkaia vizcaya', '49': 'zamora', '50': 'zaragoza', '51': 'ceuta', '52': 'melilla'} NOMBRE_A_INE = {} for k, v in NOMBRES.items(): for nm in (v.split('|') if '|' in v else v.split(' ') if k not in ('13', '38') else [v]): NOMBRE_A_INE[nm] = k def prov(nuts, nombre): if nuts and nuts[:5] in NUTS: return NUTS[nuts[:5]] return NOMBRE_A_INE.get(R2._base(nombre or '').strip()) if nombre else None # ---------------------------------------------------------------- referencias desde crudos def ref_carburantes(): campos = {'gasoleo_a': 'Precio Gasoleo A', 'gasoleo_b': 'Precio Gasoleo B', 'gasolina95': 'Precio Gasolina 95 E5'} por_mes = collections.defaultdict(list) # (mes, prov, prod) → valores diarios for f in sorted(glob.glob(f'{CRU}/carburantes/geoportal/geo_hist_*.json.gz')): doc = json.loads(gzip.decompress(open(f, 'rb').read()).decode('utf-8')) dd, mm, aa = doc['Fecha'][:10].split('/') mes = f'{aa}-{mm}' suma = collections.defaultdict(float); cuenta = collections.Counter() for est in doc['ListaEESSPrecio']: p_ = est.get('IDProvincia') for prod, campo in campos.items(): t = (est.get(campo) or '').strip() if t: v = float(t.replace(',', '.')) if v > 0: suma[(p_, prod)] += v; cuenta[(p_, prod)] += 1 for prod in campos: claves = [(p_, pr) for (p_, pr) in cuenta if pr == prod and p_ not in ('35', '38', '51', '52')] nac = sum(suma[c] for c in claves) / sum(cuenta[c] for c in claves) por_mes[(mes, '00', prod)].append(nac) for (p_, pr) in claves: por_mes[(mes, p_, prod)].append(suma[(p_, pr)] / cuenta[(p_, pr)] if cuenta[(p_, pr)] >= 5 else nac) ref = {k: (sum(v) / len(v)) / 1.21 for k, v in por_mes.items()} # WOB: hoja 1 («with taxes»), columna ES_price_with_tax_heating_oil, €/1000 l z = zipfile.ZipFile(f'{CRU}/carburantes/wob/wob_history.xlsx') NS = '{http://schemas.openxmlformats.org/spreadsheetml/2006/main}' comp = [''.join(x.itertext()) for x in ET.fromstring(z.read('xl/sharedStrings.xml')).iter(NS + 'si')] hoja = ET.fromstring(z.read('xl/worksheets/sheet1.xml')) col_hoil = col_die = col_95 = None; sem = collections.defaultdict(lambda: collections.defaultdict(list)) for fila in hoja.iter(NS + 'row'): celdas = {} for c in fila.iter(NS + 'c'): v = c.find(NS + 'v') if v is None: continue x = comp[int(v.text)] if c.get('t') == 's' else v.text celdas[re.sub(r'\d', '', c.get('r'))] = x nfila = int(fila.get('r')) if nfila == 1: for k, v in celdas.items(): if v == 'ES_price_with_tax_heating_oil': col_hoil = k if v == 'ES_price_with_tax_diesel': col_die = k if v == 'ES_price_with_tax_euro95': col_95 = k continue if nfila < 4: continue try: fecha = dt.date(1899, 12, 30) + dt.timedelta(days=int(float(celdas.get('A')))) except (TypeError, ValueError): continue for nombre, col in (('c', col_hoil), ('a', col_die), ('95', col_95)): if col and celdas.get(col) not in (None, ''): sem[f'{fecha:%Y-%m}'][nombre].append(float(celdas[col]) / 1000) for mes, d in sem.items(): if d['c']: ref[(mes, '00', 'gasoleo_c')] = (sum(d['c']) / len(d['c'])) / 1.21 wob = {mes: {k: sum(v) / len(v) for k, v in d.items() if v} for mes, d in sem.items()} return ref, wob def ref_electricidad(): e = json.load(open(f'{CRU}/electricidad/eurostat/nrg_pc_205_ES.json')) dims, tam = e['id'], e['size'] pos = {d: {c: i for c, i in e['dimension'][d]['category']['index'].items()} for d in dims} salida = {} for t, it in pos['time'].items(): for b, ib in pos['nrg_cons'].items(): idx = [] for d in dims: idx.append({'time': it, 'nrg_cons': ib, 'tax': pos['tax'].get('X_VAT')}.get(d, 0)) flat = 0 for i, d in enumerate(dims): flat = flat * tam[i] + idx[i] v = e['value'].get(str(flat)) if v is not None: salida[(t, b)] = v return salida def ref_papel(): out = [] for r in csv.DictReader(open(PAPEL_CSV)): of = sorted(float(x) for x in r['ofertas_adjudicatarios_eur_caja'].split()) n = len(of); med = of[n // 2] if n % 2 else (of[n // 2 - 1] + of[n // 2]) / 2 out.append((int(r['lote']), int(r['gramaje']), r['vigencia_desde'][:7], r['vigencia_hasta'][:7], med)) return out # ---------------------------------------------------------------- cálculo def calcular(T, RC, RE, RP): res = {} for r in T: m = r['mercado'] o = r['objeto'] or '' f0 = next((x for x in (r['fecha_adj'], r['fecha_formalizacion'], r['fecha_pub']) if x is not None), None) mes = f0.strftime('%Y-%m') if (f0 and dt.date(2015, 1, 1) <= f0 <= dt.date(2026, 10, 31)) else (f"{r['anio_ref']}-07" if r['anio_ref'] else None) pv = prov(r['nuts'], r['provincia']) gref = m == 'papel' and (re.search(r"(acuerdo marco|a\.?\s?m\.?)\s*(n.?\s*)?20/(20)?(19|21|23|25)(?![0-9a-z_])", o.lower()) is not None or 'racionalizaci' in (r['organo_nombre'] or '').lower()) k = (r['source'], r['id_canonical'], m) d = dict(mes=mes, grupo=gref, razon=None, valido=False, encima=None) def fin(est): d['estado'] = est; res[k] = d if r['copia_de'] is not None or r['correccion_imp_adj'] == 'copia_otra_fuente' or r['cuenta'] is False: fin('copia'); continue imp = r['imp'] if imp is None or imp <= 0: fin('sin_importe'); continue if r['rol_am'] in ('acuerdo', 'sda') and not gref: fin('acuerdo_marco'); continue if pv in ('35', '38', '51', '52') or (r['nuts'] or '')[:4] in ('ES70', 'ES63', 'ES64'): fin('fuera_territorio'); continue lim = {'carburantes': ('2015-01', '2026-10'), 'electricidad': ('2015-01', '2026-06'), 'papel': ('2020-05', '2027-03')}[m] if mes is None or mes < lim[0] or mes > lim[1]: fin('fuera_periodo'); continue if m == 'carburantes': if re.search(r"limpiador|aditivo|hidrolimpiadora|congelador", R2._base(o)): fin('excluido_v2'); continue c = R2.cantidad_carb(o) if not c or not c[0] > 0: fin('sin_cantidad'); continue # desvío 1: cantidad 0 = sin cantidad p = R2.producto_carb(o) if p not in ('gasoleo_a', 'gasoleo_b', 'gasoleo_c', 'gasolina95'): fin('sin_producto_comparable'); continue ref = RC.get((mes, '00', p)) if p == 'gasoleo_c' else RC.get((mes, pv, p), RC.get((mes, '00', p))) if ref is None: fin('sin_referencia'); continue q = c[0]; nucleo = q >= 200 and c[2] == 1 elif m == 'electricidad': e = R2.energia_elec(o) if e['motivo'] != 'ok': fin(e['motivo']); continue if not e['kwh'] > 0: fin('sin_cantidad'); continue dur = r['duracion_meses'] if r['duracion_meses'] and r['duracion_meses'] > 0 else 12.0 q = e['kwh'] * dur / 12.0 if e['anual'] else e['kwh'] anual_mwh = (e['kwh'] if e['anual'] else e['kwh'] * 12.0 / dur) / 1000.0 semestre = mes[:4] + ('-S1' if int(mes[5:7]) < 7 else '-S2') ref = RE.get((semestre, R2.banda_eurostat(anual_mwh))) if ref is None: fin('sin_referencia'); continue nucleo = q >= 1000 else: pp = R2.papel(o) if pp['motivo'] != 'ok': fin(pp['motivo']); continue lote = 2 if pp['fibra'] == 'reciclado' else 1 cand = [x[4] for x in RP if x[0] == lote and x[1] == pp['gramaje'] and x[2] <= mes <= x[3]] if not cand: fin('sin_referencia'); continue ref = cand[0]; q = pp['cajas_eq']; nucleo = q >= 2 razon = imp / q / ref d['razon'] = razon if not nucleo: fin('fuera_nucleo'); continue alto = 2.5 if m == 'papel' else 2.0 if razon < 0.5 or razon > alto: fin('dudoso'); continue d['valido'] = True; d['encima'] = razon > 1.15 d['organo'] = r['organo_clave']; d['producto'] = p if m == 'carburantes' else 'todos' fin('valido') return res def organos(res, nif_de): out = {} for m, por_prod in (('carburantes', True), ('papel', False)): V = [(k, d) for k, d in res.items() if k[2] == m and d['estado'] == 'valido' and not d['grupo']] tot = collections.Counter(); enc = collections.Counter(); g = collections.defaultdict(lambda: [0, 0]) for k, d in V: tot[d['producto']] += 1; enc[d['producto']] += d['encima'] if d['organo']: g[(d['organo'], d['producto'])][0] += 1; g[(d['organo'], d['producto'])][1] += d['encima'] lista = [] for (oc, pr), (n, kk) in g.items(): if n < 10: continue n0 = tot[pr] - n; p0 = (enc[pr] - kk) / n0 if n0 else None if p0 is None or p0 == 0: pval = 0.0 if kk > 0 else 1.0 else: pval = 0.0 for i in range(kk, n + 1): pval += math.comb(n, i) * (p0 ** i) * ((1 - p0) ** (n - i)) pval = min(1.0, pval) if kk > 0 else 1.0 lista.append([pval, oc, pr, n, kk]) lista.sort() M = len(lista); sigue = True for i, u in enumerate(lista): sigue = sigue and u[0] <= 0.05 / (M - i) out[(m, u[1], u[2])] = (u[3], u[4], 'por_encima_mas_frecuente' if sigue else 'ninguno') return out def main(): T = pq.read_table(BASE).to_pylist() RC, WOB = ref_carburantes(); RE = ref_electricidad(); RP = ref_papel() res = calcular(T, RC, RE, RP) org2 = organos(res, None) # validación de la referencia de carburantes (§7.4): geoportal nacional vs WOB cmp = [] for (mes, p_, prod), v in RC.items(): if p_ == '00' and prod in ('gasoleo_a', 'gasolina95') and mes in WOB: w = WOB[mes].get('a' if prod == 'gasoleo_a' else '95') if w: cmp.append((mes, prod, v * 1.21 / w - 1)) dentro = sum(1 for x in cmp if abs(x[2]) <= 0.03) salida = dict(referencia_carb_vs_wob=dict(meses_producto=len(cmp), dentro_3pct=dentro, pct=100 * dentro / len(cmp) if cmp else None, peor=sorted(cmp, key=lambda x: -abs(x[2]))[:5])) # comparación con la implementación 1 if os.path.exists(C1): I1 = {(r['source'], r['id_canonical'], r['mercado']): r for r in csv.DictReader(open(C1))} UNIV_FUERA = {'copia', 'sin_importe', 'acuerdo_marco', 'fuera_territorio', 'fuera_periodo'} for m in ('carburantes', 'electricidad', 'papel'): ks = [k for k in set(res) | set(k for k in I1 if k[2] == m) if k[2] == m] dif = collections.Counter(); ej = [] univ = sum(1 for k in ks if k in res and res[k]['estado'] not in UNIV_FUERA and not res[k]['grupo']) val = sum(1 for k in ks if (k in res and res[k]['valido']) or (k in I1 and I1[k]['valido'] == 'True')) for k in ks: a, b = I1.get(k), res.get(k) if a is None or b is None: dif['falta'] += 1; continue fallos = [] if a['estado'] != b['estado']: fallos.append('estado') if (a['valido'] == 'True') != b['valido']: fallos.append('valido') ra = float(a['razon']) if a['razon'] else None if (ra is None) != (b['razon'] is None) or (ra is not None and abs(ra - b['razon']) > 1e-6 * max(1, abs(ra))): fallos.append('razon') if b['valido'] and (a['por_encima'] == 'True') != b['encima']: fallos.append('por_encima') if fallos: dif['contratos'] += 1; dif.update(fallos) dif['en_validos'] += (a['valido'] == 'True' or b['valido']) if len(ej) < 15: ej.append(dict(k=list(k), i1=a['estado'], i2=b['estado'], r1=a['razon'], r2=b['razon'], fallos=fallos)) salida[m] = dict(filas=len(ks), universo=univ, validos_union=val, diferencias=dict(dif), pct_universo=100 * dif['contratos'] / univ if univ else 0, pct_validos=100 * dif['en_validos'] / val if val else 0, ejemplos=ej, validos_impl2=sum(1 for k in ks if k in res and res[k]['valido'] and not res[k]['grupo'])) if os.path.exists(O1): O = {(r['mercado'], r['organo_clave'], r['producto']): (int(r['n']), int(r['k_por_encima']), r['hecho']) for r in csv.DictReader(open(O1))} claves = set(O) | set(org2) iguales = sum(1 for c in claves if O.get(c) == org2.get(c)) salida['organos'] = dict(unidades=len(claves), iguales=iguales, pct_iguales=100 * iguales / len(claves) if claves else 100.0, distintas=[[list(c), O.get(c), org2.get(c)] for c in claves if O.get(c) != org2.get(c)][:10], hechos_impl2=sum(1 for v in org2.values() if v[2] != 'ninguno')) json.dump(salida, open(SAL, 'w'), indent=1, ensure_ascii=False, default=str) print(json.dumps({k: (v if k in ('organos', 'referencia_carb_vs_wob') else {kk: v[kk] for kk in ('universo', 'validos_impl2', 'diferencias', 'pct_universo', 'pct_validos')}) for k, v in salida.items()}, ensure_ascii=False, default=str)[:3000]) if __name__ == '__main__': main()