# -*- coding: utf-8 -*-
'''
#Autores
Fidel Serrano, Víctor Hernández
Afiliación
##############
Área de Planeación Colaborativa, Laboratorio Nacional de Ciencias de la Sostenibilidad, Instituto de Ecología, Universidad Nacional Autónoma de México
Reconocimiento
#################
Este trabajo fue realizado con financiamiento del Programa de Apoyo a Proyectos de Investigación e Innovación Tecnológica (PAPIIT IV100118)
del proyecto Análisis integrado de sistemas socio-ambientales acoplados: desarrollo de capacidades para la evaluación de la vulnerabilidad costera.
'''
import copy
import pprint
import string
import qgis
import qgis.core
from osgeo import gdal
import gdal_calc
import os
import processing as pr
[documentos]def raster_min_max(path_raster):
    '''
    Esta funcion regresa los valores maximos y minimos de una capa raster
    :param path_raster: ruta de la capa raster
    :type path_raster: str
    '''
    rlayer = QgsRasterLayer(path_raster,"raster")
    extent = rlayer.extent()
    provider = rlayer.dataProvider()
    stats = provider.bandStatistics(1,
                                    QgsRasterBandStats.All,
                                    extent,
                                    0)
    v_min = stats.minimumValue
    v_max = stats.maximumValue
    return v_min, v_max 
[documentos]def get_region(path_layer):
    '''
    Esta función regresa en forma de cadena de texto
    las coordenadas de la extensión de una capa raster
    param path_layer: ruta de la capa raster
    type path_layer: str
    '''
    layer = QgsRasterLayer(path_layer,"")
    ext = layer.extent()
    xmin = ext.xMinimum()
    xmax = ext.xMaximum()
    ymin = ext.yMinimum()
    ymax = ext.yMaximum()
    region = "%f,%f,%f,%f" % (xmin, xmax, ymin, ymax)
    return region 
[documentos]def nombre_capa(path_capa):
    '''
    Esta función regresa el nombre de una capa sin extensión
    :param path_capa: ruta de la capa
    :type path_capa: str
    '''
    nombre_capa=(path_capa.split("/")[-1:])[0].split(".")[0]
    return nombre_capa 
def pesos_factores_vulnerabilidad(dicc):
    lista = []
    for k,v in dicc.items():
        lista.append(dicc[k]['w'])
    return lista
[documentos]def ecuacion_vulnerabilidad(n):
    '''
    Esta función expresa la ecuación para el cálculo de la vulnerabilidad
    .. math::
        vulnerabilidad = \exp^{( 1 - sus)^{(1 + ca)}}
        | exp = Exposición
        | sus = Susceptibilidad
        | ca = Capacidad adaptativa
    :returns: str ecuacion
    '''
    if n==1:
        ecuacion = 'pow(A,(1-B))'
    if n==2:
        ecuacion = 'pow(pow(A,(1-B)),(1+C))'
    return ecuacion 
[documentos]def lista_criterios(dicc):
    '''
    Esta función regresa una lista de los criterios de un diccionario
    :param dicc: Diccionario que contiene nombres, rutas y pesos para el
    análisis de vulnerabilidad / sensibilidad
    :type dicc: diccionario python
    '''
    criterios = []
    for k1,v1 in dicc.items():
        for k2,v2, in v1['criterios'].items():
            elementos = list(v2.keys())
            if 'ruta' in elementos:
                criterios.append(k2)
            if 'criterios' in elementos:
                for k3,v3, in v2['criterios'].items():
                    elementos = list(v3.keys())
                    if 'ruta' in elementos:
                        criterios.append(k3)
                    if 'criterios' in elementos:
                        for k4,v4, in v3['criterios'].items():
                            elementos = list(v4.keys())
                            if 'ruta' in elementos:
                                criterios.append(k4)
                            if 'criterios' in elementos:
                                for k5,v5, in v4['criterios'].items():
                                    elementos = list(v5.keys())
                                    if 'ruta' in elementos:
                                        criterios.append(k5)
    return criterios 
[documentos]def crea_capa(ecuacion,rasters_input,salida):
    '''
    Esta función crea una capa mediante la calculadora raster
    de GDAL, esta función esta limitada hasta 14 variables en la ecuación.
    :param ecuacion: ecuación expresada en formato gdal,\
                    es este caso es la salida de la funcion *ecuacion_clp*
    :type ecuacion: String
    :param rasters_input: lista de los paths de los archivos rasters, salida de la función *separa_ruta_pesos*
    :type rasters_input: lista
    :param salida: ruta con extensión tiff de la salida
    :type salida: String
    '''
    path_A=''
    path_B=''
    path_C=''
    path_D=''
    path_E=''
    path_F=''
    path_G=''
    path_H=''
    path_I=''
    path_J=''
    path_K=''
    path_L=''
    path_M=''
    path_N=''
    total_raster = len(rasters_input)
    for a,b in zip(range(total_raster), rasters_input):
        if a == 0:
            path_A=b
        elif a == 1:
            path_B=b
        elif a == 2:
            path_C=b
        elif a == 3:
            path_D=b
        elif a == 4:
            path_E=b
        elif a == 5:
            path_F=b
        elif a == 6:
            path_G=b
        elif a == 7:
            path_H=b
        elif a == 8:
            path_I=b
        elif a == 9:
            path_J=b
        elif a == 10:
            path_K=b
        elif a == 11:
            path_L=b
        elif a == 12:
            path_M=b
        elif a == 13:
            path_N=b
    if total_raster == 1:
        gdal_calc.Calc(calc=ecuacion,
                            A=path_A,
                            outfile=salida,
                            NoDataValue=-9999.0,
                            quiet=True)
    if total_raster == 2:
        gdal_calc.Calc(calc=ecuacion,
                        A=path_A,
                        B=path_B,
                        outfile=salida,
                        NoDataValue=-9999.0,
                        quiet=True)
    if total_raster == 3:
            gdal_calc.Calc(calc=ecuacion,
                            A=path_A,
                            B=path_B,
                            C=path_C,
                            outfile=salida,
                            NoDataValue=-9999.0,
                            quiet=True)
    if total_raster == 4:
        gdal_calc.Calc(calc=ecuacion,
                        A=path_A,
                        B=path_B,
                        C=path_C,
                        D=path_D,
                        outfile=salida,
                        NoDataValue=-9999.0,
                        quiet=True)
    if total_raster == 5:
            gdal_calc.Calc(calc=ecuacion,
                            A=path_A,
                            B=path_B,
                            C=path_C,
                            D=path_D,
                            E=path_E,
                            outfile=salida,
                            NoDataValue=-9999.0,
                            quiet=True)
    if total_raster == 6:
        gdal_calc.Calc(calc=ecuacion,
                        A=path_A,
                        B=path_B,
                        C=path_C,
                        D=path_D,
                        E=path_E,
                        F=path_F,
                        outfile=salida,
                        NoDataValue=-9999.0,
                        quiet=True)
    if total_raster == 7:
            gdal_calc.Calc(calc=ecuacion,
                            A=path_A,
                            B=path_B,
                            C=path_C,
                            D=path_D,
                            E=path_E,
                            F=path_F,
                            G=path_G,
                            outfile=salida,
                            NoDataValue=-9999.0,
                            quiet=True)
    if total_raster == 8:
        gdal_calc.Calc(calc=ecuacion,
                        A=path_A,
                        B=path_B,
                        C=path_C,
                        D=path_D,
                        E=path_E,
                        F=path_F,
                        G=path_G,
                        H=path_H,
                        outfile=salida,
                        NoDataValue=-9999.0,
                        quiet=True)
    if total_raster == 9:
        gdal_calc.Calc(calc=ecuacion,
                        A=path_A,
                        B=path_B,
                        C=path_C,
                        D=path_D,
                        E=path_E,
                        F=path_F,
                        G=path_G,
                        H=path_H,
                        I=path_I,
                        outfile=salida,
                        NoDataValue=-9999.0,
                        quiet=True)
    if total_raster == 10:
        gdal_calc.Calc(calc=ecuacion,
                        A=path_A,
                        B=path_B,
                        C=path_C,
                        D=path_D,
                        E=path_E,
                        F=path_F,
                        G=path_G,
                        H=path_H,
                        I=path_I,
                        J=path_J,
                        outfile=salida,
                        NoDataValue=-9999.0,
                        quiet=True)
    if total_raster == 11:
        gdal_calc.Calc(calc=ecuacion,
                        A=path_A,
                        B=path_B,
                        C=path_C,
                        D=path_D,
                        E=path_E,
                        F=path_F,
                        G=path_G,
                        H=path_H,
                        I=path_I,
                        J=path_J,
                        K=path_K,
                        outfile=salida,
                        NoDataValue=-9999.0,
                        quiet=True)
    if total_raster == 12:
        gdal_calc.Calc(calc=ecuacion,
                        A=path_A,
                        B=path_B,
                        C=path_C,
                        D=path_D,
                        E=path_E,
                        F=path_F,
                        G=path_G,
                        H=path_H,
                        I=path_I,
                        J=path_J,
                        K=path_K,
                        L=path_L,
                        outfile=salida,
                        NoDataValue=-9999.0,
                        quiet=True)
    if total_raster == 13:
        gdal_calc.Calc(calc=ecuacion,
                        A=path_A,
                        B=path_B,
                        C=path_C,
                        D=path_D,
                        E=path_E,
                        F=path_F,
                        G=path_G,
                        H=path_H,
                        I=path_I,
                        J=path_J,
                        K=path_K,
                        L=path_L,
                        M=path_M,
                        outfile=salida,
                        NoDataValue=-9999.0,
                        quiet=True)
    if total_raster == 14:
        gdal_calc.Calc(calc=ecuacion,
                        A=path_A,
                        B=path_B,
                        C=path_C,
                        D=path_D,
                        E=path_E,
                        F=path_F,
                        G=path_G,
                        H=path_H,
                        I=path_I,
                        J=path_J,
                        K=path_K,
                        L=path_L,
                        M=path_M,
                        N=path_N,
                        outfile=salida,
                        NoDataValue=-9999.0,
                        quiet=True) 
[documentos]def ecuacion_clp(pesos):
    '''
    Esta función recibe una lista de pesos para regresar la ecuación
    en la estructura requerida por gdal para la combinación lineal ponderada.
    :param pesos: lista de los pesos de las capas, salida de la función *separa_ruta_pesos*
    :type pesos: lista
    '''
    n_variables=len(pesos)
    abc = list(string.ascii_uppercase)
    ecuacion = ''
    for a,b in zip(range(n_variables),pesos):
        if a < n_variables-1:
            ecuacion+= (str(b)+str(' * ')+str(abc[a])+' + ' )
        else:
            ecuacion+= (str(b)+str(' * ')+str(abc[a]))
    return ecuacion 
[documentos]def lista_pesos_ruta(dicc):
    '''
    Funcion para sacar listas por subcriterio
    '''
    exp_fis=[]
    exp_bio=[]
    sus_fis=[]
    sus_bio=[]
    pesos_1n = []
    for k1, v1 in dicc.items():
        for k2, v2 in v1['criterios'].items():
            pesos_1n.append(k1+"|"+str(v1['w'])+"|"+k2+"|"+str(v2['w']))
            for k3, v3 in v2['criterios'].items():
                if k1=='exposicion' and k2=='fisico':
                    exp_fis.append(v3['ruta'] + "|" + str(v3['w']))
                    #print (k1,v1['w'],k2,v2['w'],k3,v3['w'],v3['ruta'])
                elif k1=='exposicion' and k2=='biologico':
                    exp_bio.append(v3['ruta'] + "|" + str(v3['w']))
                    #print (k1,v1['w'],k2,v2['w'],k3,v3['w'],v3['ruta'])
                elif k1=='susceptibilidad' and k2=='fisico':
                    sus_fis.append(v3['ruta'] + "|" + str(v3['w']))
                    #print (k1,v1['w'],k2,v2['w'],k3,v3['w'],v3['ruta']
                elif k1=='susceptibilidad' and k2=='biologico':
                    sus_bio.append(v3['ruta'] + "|" + str(v3['w']))
                    #print (k1,v1['w'],k2,v2['w'],k3,v3['w'],v3['ruta']
    return exp_fis, exp_bio, sus_fis, sus_bio,pesos_1n 
def separa_ruta_pesos(lista):
    rutas =[]
    pesos= []
    for capa in lista:
        rutas.append(capa.split("|")[0])
        pesos.append(capa.split("|")[1])
    return rutas,pesos
def rutas_pesos_globales(factor,dicc):
    lista = []   # se tiene que pasar a nivel recursiva
    for k1, v1 in dicc[factor].items():
        if type(v1) is dict:
            for k2, v2 in v1.items(): # nivel de biologicos o fisico
                if type(v2) is dict:
                    for k3,v3 in v2.items():
                        if type(v3) is dict:
                             for k4,v4 in v3.items(): #k4 e sel nombre de los criterios  como v_acuatica_exp
                                elementos = list(v4.keys())
                                if 'ruta' in elementos:
                                    lista.append(v4['ruta']+"|"+str(round(v2['w']*v4['w'],3)))
                                elif 'criterios' in elementos and  type(v4) is dict:
                                    for k5,v5, in v4.items():
                                        if type(v5) is dict:
                                            for k6,v6 in v5.items():
                                                elementos = list(v6.keys())
                                                if 'ruta' in elementos:
                                                    lista.append(v6['ruta']+"|"+str(round(v2['w']*v4['w']*v6['w'],3)))
                                                elif 'criterios' in elementos and  type(v6) is dict:
                                                    for k7,v7, in v6.items():
                                                        if type(v7) is dict:
                                                            for k8,v8 in v7.items():
                                                                elementos = list(v8.keys())
                                                                if 'ruta' in elementos:
                                                                    lista.append(v8['ruta']+"|"+str(round(v2['w']*v4['w']*v6['w']*v8['w'],3)))
                                                                elif 'criterios' in elementos and  type(v8) is dict:
                                                                        for k9,v9, in v8.items():
                                                                            if type(v9) is dict:
                                                                                for k10,v10 in v9.items():
                                                                                    elementos = list(v10.keys())
                                                                                    if 'ruta' in elementos:
                                                                                        lista.append(v10['ruta']+"|"+str(round(v2['w']*v4['w']*v6['w']*v8['w']*v10['w'],3)))
    return lista
def find_key(d, key):
    for k,v in d.items():
        if isinstance(v, dict):
            p = find_key(v, key)
            if p:
                return [k] + p
        if k == key:
            return [k]
[documentos]def quita(dicc,key):
    '''
    Esta función retira un elemento del diccionario y regresa un nuevo diccionario
    sin dicho elemento <<dicc_q>>.
    :param dicc: Diccionario con la estructura requerida
    :type dicc: diccionario
    :param key: nombre de la variable a quitar
    :type key: String
    '''
    dicc_q = copy.deepcopy(dicc)
    niveles = find_key(dicc_q,key)
    if len(niveles)==1:
        dicc_q.pop(niveles[0])
    if len(niveles)==3:
        dicc_q[niveles[0]]['criterios'].pop(niveles[2])
        if len(dicc_q[niveles[0]]['criterios'])==0:
            dicc_q.pop(niveles[0])
    if len(niveles)==5:
        dicc_q[niveles[0]]['criterios'][niveles[2]]['criterios'].pop(niveles[4])
        if len(dicc_q[niveles[0]]['criterios'][niveles[2]]['criterios'])==0:
            dicc_q[niveles[0]]['criterios'].pop(niveles[2])
        if len(dicc_q[niveles[0]]['criterios'])==0:
            dicc_q.pop(niveles[0])
    if len(niveles)==7:
        dicc_q[niveles[0]]['criterios'][niveles[2]]['criterios'][niveles[4]]['criterios'].pop(niveles[6])
        if len(dicc_q[niveles[0]]['criterios'][niveles[2]]['criterios'][niveles[4]]['criterios'])==0:
            dicc_q[niveles[0]]['criterios'][niveles[2]]['criterios'].pop(niveles[4])
        if len(dicc_q[niveles[0]]['criterios'][niveles[2]]['criterios'])==0:
            dicc_q[niveles[0]]['criterios'].pop(niveles[2])
        if len(dicc_q[niveles[0]]['criterios'])==0:
            dicc_q.pop(niveles[0])
    if len(niveles)==9:
        dicc_q[niveles[0]]['criterios'][niveles[2]]['criterios'][niveles[4]]['criterios'][niveles[6]]['criterios'].pop(niveles[8])
        if len(dicc_q[niveles[0]]['criterios'][niveles[2]]['criterios'][niveles[4]]['criterios'][niveles[6]]['criterios'])==0:
            dicc_q[niveles[0]]['criterios'][niveles[2]]['criterios'][niveles[4]]['criterios'].pop(niveles[6])
        if len(dicc_q[niveles[0]]['criterios'][niveles[2]]['criterios'][niveles[4]]['criterios'])==0:
            dicc_q[niveles[0]]['criterios'][niveles[2]]['criterios'].pop(niveles[4])
        if len(dicc_q[niveles[0]]['criterios'][niveles[2]]['criterios'])==0:
            dicc_q[niveles[0]]['criterios'].pop(niveles[2])
        if len(dicc_q[niveles[0]]['criterios'])==0:
            dicc_q.pop(niveles[0])
    return dicc_q 
[documentos]def reescala(dicc_q):
    '''
    Esta función rescala un diccionario que se le a quitado un criterio
    y regresa el diccionario con los pesos rescalados
    :param dicc_q: salida de la función *quita*
    '''
    dicc_r = copy.deepcopy(dicc_q)
    suma = 0
    for k1, v1 in dicc_r.items():
        suma += v1['w']
    for k1, v1 in dicc_r.items():
        v1['w'] = v1['w']/suma
    for k1, v1 in dicc_r.items():
        suma = 0
        elementos = list(v1.keys())
        if 'criterios' in elementos:
            for k2, v2 in v1['criterios'].items():
                suma += v2['w']
            for k2, v2 in v1['criterios'].items():
                v2['w'] = v2['w'] / suma
    for k1, v1 in dicc_r.items():
        elementos = list(v1.keys())
        if 'criterios' in elementos:
            for k2, v2 in v1['criterios'].items():
                suma = 0
                elementos = list(v2.keys())
                if 'criterios' in elementos:
                    for k3, v3 in v2['criterios'].items():
                        suma += v3['w']
                    for k3, v3 in v2['criterios'].items():
                        v3['w'] = v3['w'] / suma
    for k1, v1 in dicc_r.items():
        elementos = list(v1.keys())
        if 'criterios' in elementos:
            for k2, v2 in v1['criterios'].items():
                elementos = list(v2.keys())
                if 'criterios' in elementos:
                    for k3, v3 in v2['criterios'].items():
                        suma = 0
                        elementos = list(v3.keys())
                        if 'criterios' in elementos:
                            for k4, v4 in v3['criterios'].items():
                                suma += v4['w']
                            for k4, v4 in v3['criterios'].items():
                                v4['w'] = v4['w'] / suma
    for k1, v1 in dicc_r.items():
        elementos = list(v1.keys())
        if 'criterios' in elementos:
            for k2, v2 in v1['criterios'].items():
                elementos = list(v2.keys())
                if 'criterios' in elementos:
                    for k3, v3 in v2['criterios'].items():
                        elementos = list(v3.keys())
                        if 'criterios' in elementos:
                            for k4, v4 in v3['criterios'].items():
                                suma = 0
                                elementos = list(v4.keys())
                                if 'criterios' in elementos:
                                    for k5, v5 in v4['criterios'].items():
                                        suma += v5['w']
                                    for k5, v5 in v4['criterios'].items():
                                        v5['w'] = v5['w'] / suma
    return dicc_r 
[documentos]def quita_reescala(dicc,key):
    '''
    Función que integra las funciones quita y reescala y regresa
    un diccionario sin la variable y con los pesos reescalados.
    :param dicc: Diccionario con la estructura requerida
    :type dicc: diccionario
    :param key: nombre de la variable a quitar
    :type key: String
    '''
    dicc_q = quita(dicc,key)
    dicc_r = reescala(dicc_q)
    return dicc_r 
def pesos_superiores(dicc):
    pesos=[]
    for k1,v1 in dicc.items():
        #print (k1,v1['w'])
        pesos.append([k1,str(v1['w'])])
    return pesos
def ideales(path_raster, path_raster_n):
    min,max = raster_min_max(path_raster)
    no_data =raster_nodata(path_raster)
    ec_norm ='(A' + ') / (' + str(max) +')'  # llevar a ideal
    #ec_norm ='(A - '+str(min) + ') / (' + str(max)+'-'+str(min) +')'  # normalizar
    dicc ={
        'INPUT_A':path_raster,
        'BAND_A':1,
        'FORMULA':ec_norm,
        'NO_DATA': no_data,
        'RTYPE':5,
        'OUTPUT':path_raster_n}
    pr.run("gdal:rastercalculator",dicc)
def raster_nodata(path_raster):
    rlayer = QgsRasterLayer(path_raster,"raster")
    extent = rlayer.extent()
    provider = rlayer.dataProvider()
    rows = rlayer.rasterUnitsPerPixelY()
    cols = rlayer.rasterUnitsPerPixelX()
    block = provider.block(1, extent,  rows, cols)
    no_data = block.noDataValue()
    return no_data
def lineal_decreciente(path_raster, path_raster_n):
    min,max = raster_min_max(path_raster)
    no_data =raster_nodata(path_raster)
    xmax_menos_xmin = max-min
    xmax_menos_xmin_xmax = xmax_menos_xmin / max
    xmax_mas_xmin = max + min
    xmax_mas_xmin_xmax = xmax_mas_xmin / max
    ec_norm ='((-A * '+str(xmax_menos_xmin_xmax)+') / '+str(xmax_menos_xmin)+') + '+str(xmax_mas_xmin_xmax)  # llevar a ideal
    dicc ={
        'INPUT_A':path_raster,
        'BAND_A':1,
        'FORMULA':ec_norm,
        'NO_DATA': no_data,
        'RTYPE':5,
        'OUTPUT':path_raster_n}
    pr.run("gdal:rastercalculator",dicc)
''' # <- Quitar comillas para usar código
## rutas de entrada y salida
p_sig_exp= 'C:/Dropbox (LANCIS)/SIG/desarrollo/sig_papiit/entregables/\
exposicion/'
p_sig_sens= 'C:/Dropbox (LANCIS)/SIG/desarrollo/sig_papiit/entregables/\
susceptibilidad/'
p_sig_res= 'C:/Dropbox (LANCIS)/SIG/desarrollo/sig_papiit/entregables/\
resiliencia/'
p_procesamiento = 'C:/vulnerabilidad13/'
zonas_ambiente_construido  = 'C:/Dropbox (LANCIS)/SIG/desarrollo/sig_papiit/entregables/mascaras/zonas_ambiente_construido.tif'
dicc = {
    'exposicion': {'w':0.333,
                        'criterios':{'biologico':{'w':0.50,
                                                        'criterios':{'v_acuatica_exp':{'w':0.16,'ruta':p_sig_exp+'biologica/v_acuatica_yuc/fv_v_acuatica_yuc_100m.tif'},
                                                                        'v_costera_exp':{ 'w':0.84,
                                                                                            'criterios':{'distancia_manglar_exp':{'w':0.75,'ruta':p_sig_exp + 'biologica/v_costera_yuc/fv_distancia_manglar_yuc_100m.tif'},
                                                                                                        'distancia_dunas_exp':{'w':0.25,'ruta':p_sig_exp + 'biologica/v_costera_yuc/fv_distancia_dunas_yuc_100m.tif'}}}}},
                                        'fisico':{'w':0.50,
                                                    'criterios':{'elevacion_exp':{ 'w':0.87,'ruta':p_sig_exp+ 'fisica/elev_yuc/fv_elevacion_yuc_100m.tif'},
                                                                    'ancho_playa_exp':{'w':0.13,'ruta':p_sig_exp+ 'fisica/ancho_playa_yuc/fv_distancia_playa_yuc_100m.tif'}
                                                                }
                                                }
                                    }
                    },
    'susceptibilidad': {'w':0.333,
                            'criterios':{'biologico':{'w':0.50 ,
                                                            'criterios':{'v_costera_sus':{ 'w':1.0,'ruta':p_sig_sens + 'biologica/v_costera_yuc/fv_v_costera_presencia_yuc_100m.tif'}}},
                                            'fisico':{'w':0.50,
                                                        'criterios':{'elevacion_sus':{ 'w':0.26,'ruta':p_sig_sens + 'fisica/elev_yuc/fv_elevacion_yuc_100m.tif' },
                                                                        'duna_costera_sus':{'w':0.10,'ruta':p_sig_sens + 'fisica/duna_yuc/fv_duna_yuc_100m.tif'},
                                                                        'tipo_litoral_sus':{'w':0.64,'ruta':p_sig_sens + 'fisica/t_litoral_yuc/fv_tipo_litoral_yuc_100m_corregida.tif'}
                                                                    }
                                                    }
                                        }
                        },
    'resiliencia': {'w':0.333,
                            'criterios':{'biologico':{'w':0.50 ,
                                                            'criterios':{'biodiversidad_res':{'w':0.50,'ruta':p_sig_res + 'biologica/biodiversidad/fv_biodiversidad_yuc_100m.tif'},
                                                                        'servicios_ambientales_res':{'w':0.50,
                                                                            'criterios':{'proteccion_costera_res':{'w':0.75,
                                                                                                                'criterios':{'entidades_protectoras_res':{'w':0.84,'ruta':p_sig_res + 'biologica/serv_ambientales/prot_costera_yuc/fv_entidades_protectoras.tif'},
                                                                                                                            'v_acuatica_res':{'w':0.16,'ruta':p_sig_res + 'biologica/serv_ambientales/prot_costera_yuc/fv_v_acuatica_yuc_100m.tif'}
                                                                                                                            }
                                                                                                                    },
                                                                                        'provision_res':{'w':0.25,'ruta':p_sig_res + 'biologica/serv_ambientales/provision_yuc/fv_provision_yuc_100m.tif'}}
                                                                                                    }
                                                                        }
                                                    },
                                        'fisico':{'w':0.50,
                                                        'criterios':{'elevacion_res':{'w':0.60,'ruta':p_sig_res + 'fisica/elev_yuc/fv_elev_yuc.tif'},
                                                                    'tipo_litoral_res':{'w':0.40,'ruta':p_sig_res + 'fisica/t_litoral_yuc/fv_tipo_litoral_yuc.tif'}
                                                                    }
                                                }
                                        }
                    }
    }
# ----- GENERACIÓN DE LAS CAPAS CONSIDERANDO TODOS LOS ELEMENTOS ----- #
exposicion_total,susceptibilidad_total,resiliencia_total = rutas_pesos_globales('exposicion',dicc),rutas_pesos_globales('susceptibilidad',dicc),rutas_pesos_globales('resiliencia',dicc)
path_exp_t,w_exp_t = separa_ruta_pesos(exposicion_total)
path_sus_t,w_sus_t = separa_ruta_pesos(susceptibilidad_total)
path_res_t,w_res_t = separa_ruta_pesos(resiliencia_total)
ecuacion_exp_t = ecuacion_clp(w_exp_t)
ecuacion_sus_t = ecuacion_clp(w_sus_t)
ecuacion_res_t = ecuacion_clp(w_res_t)
dir_tmp = p_procesamiento+"tmp"
if "tmp" not in os.listdir(p_procesamiento):
    os.mkdir(dir_tmp)
exposicion = dir_tmp+"/tp_exposicion.tif"
susceptibilidad = dir_tmp+"/tp_susceptibilidad.tif"
tp_resiliencia = dir_tmp+"/tp_resiliencia.tif"
exposicion_ze = dir_tmp+"/tp_exposicion_ze.tif"
susceptibilidad_ze = dir_tmp+"/tp_susceptibilidad_ze.tif"
tp_resiliencia_ze = dir_tmp+"/tp_resiliencia_ze.tif"
crea_capa(ecuacion_exp_t,path_exp_t,exposicion)
crea_capa(ecuacion_sus_t,path_sus_t,susceptibilidad)
crea_capa(ecuacion_res_t,path_res_t,tp_resiliencia)
crea_capa('A*B',[exposicion,zonas_ambiente_construido],exposicion_ze)
crea_capa('A*B',[susceptibilidad, zonas_ambiente_construido],susceptibilidad_ze)
crea_capa('A*B',[tp_resiliencia, zonas_ambiente_construido],tp_resiliencia_ze)
resiliencia = dir_tmp+"/tp_resiliencia_dec.tif"
# --- Creación de la capa de vulnerabilidad ---##
vulnerabilidad= dir_tmp +"/tp_vulnerabilidad.tif"
lista_exp_sus_res = [exposicion,susceptibilidad,resiliencia]
#lista_exp_sus_res = [exposicion_ze,susceptibilidad_ze,resiliencia]  ## CAPAS CON ZONAS DE EXCLUSION
w_vulnerabilidad = pesos_factores_vulnerabilidad(dicc)
ecuacion_vul_lineal = ecuacion_clp(w_vulnerabilidad)
crea_capa(ecuacion_vul_lineal,lista_exp_sus_res,vulnerabilidad)
### Datos para el analisis de sensibilidad ####
vulnerabilidad_total_media = media_raster(vulnerabilidad)
exp_media_total = media_raster(exposicion)
sus_media_total = media_raster(susceptibilidad)
res_media_total = media_raster(tp_resiliencia)
## PROMEDIOS DE CAPAS CON ZONAS DE EXCLUSION
# exp_media_total = media_raster(exposicion_ze)
# sus_media_total = media_raster(susceptibilidad_ze)
# res_media_total = media_raster(tp_resiliencia_ze)
#---- TERMINA ---#
#--- ANALISIS DE SENSIBILIDAD ----#
archivo = open(p_procesamiento+"sensibilidad_por_criterio_ze.csv","w")
archivo.write("criterio,exposicion,sensibilidad_exp,suceptibilidad,sensibilidad_sus,resiliencia,sensibilidad_res,vulnerabilidad,sensibilidad_vul\n")
archivo.write("total,"+str(round(exp_media_total,3))+",,"+str(round(sus_media_total,3))+",,"+str(round(res_media_total,3))+",,"+str(round(vulnerabilidad_total_media,3))+","+"\n")
criterios = lista_criterios(dicc) #obtiene criterios del diccionario principal
cont=0
for criterio in criterios:
    cont +=1
    print ("procensado criterio: ",criterio,"  ",cont,"de",len(criterios))
    dicc2  = quita_reescala(dicc,criterio)
    exposicion,susceptibilidad,resiliencia = rutas_pesos_globales('exposicion',dicc2),rutas_pesos_globales('susceptibilidad',dicc2),rutas_pesos_globales('resiliencia',dicc2) #separa los subcriterios por criterios
    path_exp,w_exp = separa_ruta_pesos(exposicion)
    path_sus,w_sus = separa_ruta_pesos(susceptibilidad)
    path_res,w_res = separa_ruta_pesos(resiliencia)
    ecuacion_exp = ecuacion_clp(w_exp)
    ecuacion_sus = ecuacion_clp(w_sus)
    ecuacion_res = ecuacion_clp(w_res)
    salida_exposicion = dir_tmp+"/tp_exposicion_sin_"+criterio+".tif"
    salida_susceptibilidad = dir_tmp+"/tp_suscep_sin_"+criterio+".tif"
    salida_resiliencia = dir_tmp+"/tp_resilencia_sin_"+criterio+".tif"
    salida_exposicion_ze = dir_tmp+"/tp_exposicion_sin_"+criterio+"_ze.tif"
    salida_susceptibilidad_ze = dir_tmp+"/tp_suscep_sin_"+criterio+"_ze.tif"
    salida_resiliencia_ze = dir_tmp+"/tp_resilencia_sin_"+criterio+"_ze.tif"
    crea_capa(ecuacion_exp,path_exp,salida_exposicion)
    crea_capa(ecuacion_sus,path_sus,salida_susceptibilidad)
    crea_capa(ecuacion_res,path_res,salida_resiliencia)
    crea_capa('A*B',[salida_exposicion,zonas_ambiente_construido],salida_exposicion_ze)
    crea_capa('A*B',[salida_susceptibilidad, zonas_ambiente_construido],salida_susceptibilidad_ze)
    crea_capa('A*B',[salida_resiliencia, zonas_ambiente_construido],salida_resiliencia_ze)
    salida_resiliencia_lineal_dec = dir_tmp+"/tp_resilencia_sin_"+criterio+"_ld.tif"
    #lineal_decreciente(salida_resiliencia_ze,salida_resiliencia_lineal_dec) # LINEAL DECRECIENTE DE RESILIENCIA A ZONAS DE EXCLUSION
    lineal_decreciente(salida_resiliencia,salida_resiliencia_lineal_dec)
    salida_vulnerabilidad = dir_tmp+"/tp_vulnerabilidad_sin_"+criterio+".tif"
    # --- Creación de la capa de vulnerabilidad --- ##
    # lista_c = [salida_exposicion_ze,salida_susceptibilidad_ze,salida_resiliencia_lineal_dec] #VULNERABILIDAD CON ZONAS DE EXCLUSION
    lista_c = [salida_exposicion,salida_susceptibilidad,salida_resiliencia_lineal_dec] # VULNERABILIDAD SIN ZONAS DE EXCLUSION
    w_vulnerabilidad = pesos_factores_vulnerabilidad(dicc2)
    ecuacion_vul_lineal = ecuacion_clp(w_vulnerabilidad)
    crea_capa(ecuacion_vul_lineal,lista_c,salida_vulnerabilidad)
    print('pesos vulnerabilidad',w_vulnerabilidad)
    # --- Promedios de los criterios ---- ##
    vulnerabilidad_media = media_raster(salida_vulnerabilidad)
    exp_media = media_raster(salida_exposicion)
    sus_media = media_raster(salida_susceptibilidad)
    res_media = media_raster(salida_resiliencia)
    ### PROEMDIOS DE CRITERIOS CON ZONAS DE EXCLUSION ##
    # exp_media = media_raster(salida_exposicion_ze)
    # sus_media = media_raster(salida_susceptibilidad_ze)
    # res_media = media_raster(salida_resiliencia_ze)
    #---  Cálculo de la sensibilidad ---##
    sensibilidad_exp_calculada = round((abs((exp_media_total-exp_media))/exp_media_total),3)
    sensibilidad_sus_calculada = round((abs((sus_media_total-sus_media))/sus_media_total),3)
    sensibilidad_res_calculada = round((abs((res_media_total-res_media))/res_media_total),3)
    sensibilidad_calculada = round((abs((vulnerabilidad_total_media-vulnerabilidad_media))/vulnerabilidad_total_media),3)
    archivo.write(criterio+","+\
                    str(exp_media)+","+\
                    str(sensibilidad_exp_calculada)+","+\
                    str(round(sus_media,3))+","+\
                    str(sensibilidad_sus_calculada)+","+\
                    str(round(res_media,3))+","+\
                    str(sensibilidad_res_calculada)+","+\
                    str(round(vulnerabilidad_media,3))+","+\
                    str(sensibilidad_calculada)+"\n")
archivo.close()
''' # <- Quitar comillas para usar código