#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Simulação reprodutível da reforma do FPM-Interior (Programa SOMA / Doutrina 01).

Dados reais via API pública do IBGE:
  - População: Censo 2022 (agregado 4709, variável 93).
  - Privação D: rendimento domiciliar per capita MÉDIO, Censo 2022
    (agregado 10295, variável 13431, R$/mês). D = inverso normalizado do log
    da renda → renda baixa = privação alta. (Medida oficial e melhor que o PIB
    per capita; o CadÚnico/MDS seria equivalente, mas não tem API estável.)
  - Bolo do FPM-Interior 2024: ~R$ 147,9 bi (86,4% de R$ 171,2 bi; Tesouro/CNM).

Mecânica: tabela vigente (DL 1.881/1981) -> coeficiente contínuo c(P)
(rampa-piso + interpolação dos nós) -> fator de equidade (1+α·D). Quota sobre
bolo fixo (massa preservada). Imprime sensibilidade a α em R$ e % de proteção
dos municípios pequenos-e-pobres.

Uso:  python3 simulacao_fpm.py
"""
import json, urllib.request, gzip, math, time

POP ="https://servicodados.ibge.gov.br/api/v3/agregados/4709/periodos/2022/variaveis/93?localidades=N6[all]"
REND="https://servicodados.ibge.gov.br/api/v3/agregados/10295/periodos/2022/variaveis/13431?localidades=N6[all]&classificacao=2[6794]|86[95251]|58[95253]"

CAPITAIS={'1200401','2704302','1600303','1302603','2927408','2304400','5300108','3205309','5208707',
'2111300','5103403','5002704','3106200','1501402','2507507','4106902','2611606','2211001','3304557',
'2408102','4314902','1100205','1400100','4205407','3550308','2800308','1721000'}
TAB=[(10188,0.6),(13584,0.8),(16980,1.0),(23772,1.2),(30564,1.4),(37356,1.6),(44148,1.8),(50940,2.0),
(61128,2.2),(71316,2.4),(81504,2.6),(91692,2.8),(101880,3.0),(115464,3.2),(129048,3.4),(142632,3.6),(156216,3.8)]
CMIN=0.40; POOL=147.9e9

def get(url):
    for att in range(4):
        try:
            r=urllib.request.Request(url,headers={'Accept-Encoding':'identity','User-Agent':'musoa/1.0'})
            raw=urllib.request.urlopen(r,timeout=180).read()
            if raw[:2]==b'\x1f\x8b': raw=gzip.decompress(raw)
            return json.loads(raw)
        except Exception:
            if att==3: raise
            time.sleep(3)
def serie(url):
    d=get(url); o={}
    for s in d[0]['resultados'][0]['series']:
        try: o[s['localidade']['id']]=(s['localidade']['nome'],float(list(s['serie'].values())[0]))
        except: pass
    return o
def k(p):
    for lim,c in TAB:
        if p<=lim: return c
    return 4.0
def cc(p):
    if p<=10188: return CMIN+(0.6-CMIN)*(p/10188)
    if p>156216: return 4.0
    for i in range(len(TAB)-1):
        l0,c0=TAB[i]; l1,c1=TAB[i+1]
        if l0<p<=l1: return c0+(c1-c0)*(p-l0)/(l1-l0)
    return 3.8

def main():
    pop=serie(POP); rend=serie(REND)
    I=[]
    for cod,(nome,p) in pop.items():
        if cod in CAPITAIS: continue
        p=int(p); inc=rend.get(cod,(None,0))[1]
        I.append([cod,nome,p,inc])
    incs=[r[3] for r in I if r[3]>0]; lo,hi=math.log(min(incs)),math.log(max(incs))
    for r in I: r.append((hi-math.log(r[3]))/(hi-lo) if r[3]>0 else 1.0)  # D
    Kt=sum(k(p) for _,_,p,_,_ in I)
    ra=lambda p: POOL*k(p)/Kt
    print(f"FPM-Interior: {len(I)} municípios | bolo ~ R$ {POOL/1e9:.1f} bi")
    print(f"piso 0,6 = {sum(1 for _,_,p,_,_ in I if k(p)==0.6)} municípios; cada um recebe ~ R$ {ra(900)/1e6:.1f} mi/ano")
    print("\nsensibilidade ao peso de equidade α (renda do Censo 2022 como privação):")
    print(" α  | realocado (R$/%)        | pequenos-pobres protegidos")
    for ALPHA in (0.5,1.0,1.5,2.0):
        We=sum(cc(p)*(1+ALPHA*D) for _,_,p,_,D in I)
        re=lambda p,D: POOL*cc(p)*(1+ALPHA*D)/We
        pob=[r for r in I if r[2]<=10188 and r[4]>=0.6]
        gp=sum(1 for c,n,p,i,D in pob if re(p,D)>=ra(p))
        mov=sum(abs(re(p,D)-ra(p)) for _,_,p,_,D in I)/2
        print(f" {ALPHA:<3}| R$ {mov/1e9:4.2f} bi ({100*mov/POOL:4.1f}%) | {gp}/{len(pob)} ({100*gp/len(pob):.0f}%)")

if __name__=='__main__': main()