4. LoRaWan – Funciones girni_lorawan_lib

El archivo contiene un grupo de  funciones usadas para procesar los datos de este prototipo, se separaron del bloque principal de instrucciones con  el objetivo de simplificar el desarrollo del las actividades principales.

Recuerde disponer de este archivo en la misma carpeta que el algoritmo que invoca a la librería.

archivo de libreria: girni_lorawan_libreria

Instrucciones Python

# Girni LoRaWan librerias 2021-10-28
# LoRaWan, lecturas de Rssi y SNR
# Propuesta: edelros@espol.edu.ec
# http://blog.espol.edu.ec/girni/
import numpy as np
import json
import pandas as pd
import datetime as dt
import os
import matplotlib.pyplot as plt
import matplotlib.dates as mdates
import matplotlib.units as munits

import geopandas as gpd
import fiona
import utm

def tablapunto(unarchivoDB,unsensor,gatewayDB,fechainicio,fechafin,zonaGMT=0):
    '''Lectura de archivo.json hacia una tabla0
       selecciona registros de "un sensor"
       graba unreporte.csv con pandas
    '''
    campos = ['domain','entity_id','state','attributes','created']

    # Lectura de archivo json, cambia a formato DataFrame
    tabla0 = pd.DataFrame()
    archivoexiste = os.path.exists(unarchivoDB)
    if archivoexiste:
        tabla0 = pd.read_json(unarchivoDB)
        tabla0 = pd.DataFrame(tabla0,columns=campos)
    else:
        print(' ERROR: NO se encuentra el archivo...')
        print('        revise el nombre de archivo. ')

    # Revisa cada registro 
    leidos = 0
    tabla  = {}
    if archivoexiste:
        # Intervalo de fecha como datetime
        hora_desplaza  = dt.timedelta(hours = abs(zonaGMT))
        fechaformatoDB = '%Y-%m-%dT%H:%M:%S.%f'
        fechaformato   = '%Y-%m-%d %H:%M:%S.%f'
        fechatxt     = fechainicio[0:np.min([26,len(fechainicio)])]
        fechainicio  = dt.datetime.strptime(fechatxt,fechaformato)
        fechatxt     = fechafin[0:np.min([26,len(fechafin)])]
        fechafin     = dt.datetime.strptime(fechatxt,fechaformato)

        # datos de trama
        for cadaregistro in tabla0.index:
            unatrama   = tabla0['attributes'][cadaregistro]
            trama_mqtt = json.loads(unatrama)
            
            # selecciona registro por sensor, en fecha y con datos
            cualsensor = (tabla0['entity_id'][cadaregistro] == unsensor)
            
            unafecha = tabla0['created'][cadaregistro]
            unafecha = unafecha[0:np.min([26,len(unafecha)])]
            unafecha = dt.datetime.strptime(unafecha,fechaformato)
            unafecha = unafecha - hora_desplaza
            
            enfecha  = (unafecha>=fechainicio) and (unafecha<=fechafin)
            condatos = 'applicationID' in trama_mqtt.keys()
            
            selecciona = cualsensor and condatos and enfecha
            if (selecciona == True):        
                # datos en la mensaje MQTT
                publishedAt = trama_mqtt["publishedAt"]
                publishedAt = publishedAt[0:np.min([26,len(publishedAt)])]
                publishedAt = dt.datetime.strptime(publishedAt,fechaformatoDB)
                publishedAt = publishedAt - hora_desplaza
                
                friendly_name = trama_mqtt["friendly_name"]
                deviceName    = trama_mqtt["deviceName"]
                dispositivo   = friendly_name.split('_')[2]
                
                dr    = trama_mqtt["dr"]
                fCnt  = trama_mqtt["fCnt"]
                fPort = trama_mqtt["fPort"]

                objectJSON  = trama_mqtt["objectJSON"]
                datosensor  = json.loads(objectJSON)
                
                freq_tx   = trama_mqtt["txInfo"]["frequency"]
                bandwidth   = trama_mqtt["txInfo"]["loRaModulationInfo"]["bandwidth"]
                spreadingFactor = trama_mqtt["txInfo"]["loRaModulationInfo"]["spreadingFactor"]
                     
                datostrama = {"publishedAt": publishedAt,
                              "dispositivo": dispositivo,
                              "fCnt": fCnt}
                
                # revisa gateway en la red LoRaWan
                tamano = len(trama_mqtt["rxInfo"])
                i = 0
                while i<tamano:
                    gatewayID = trama_mqtt["rxInfo"][i]["gatewayID"]
                    rssi      = trama_mqtt["rxInfo"][i]["rssi"]
                    loRaSNR   = trama_mqtt["rxInfo"][i]["loRaSNR"]
                    channel   = trama_mqtt["rxInfo"][i]["channel"]
                    rfChain   = trama_mqtt["rxInfo"][i]["rfChain"]
                    gtwNum    = gatewayDB[gatewayID]
                    
                    # registra en tabla, incluyendo tramas repetidas
                    datostrama['rssi_up']    = rssi
                    datostrama['snr_up']     = loRaSNR
                    datostrama['channel_up'] = channel
                    datostrama['rfChain_up'] = rfChain
                    datostrama['gtw_rx']     = gtwNum
                    i = i + 1

                # dato del sensor
                equivale = {'Down_datarate': 'dr_down',
                            'Down_rssi'    : 'rssi_down',
                            'Down_snr'     : 'snr_down',
                            'bateria_V'    : 'bateria_V'}
                for undato in datosensor:
                    if undato in equivale.keys():
                        datostrama[equivale[undato]] = datosensor[undato]
                    else:
                        datostrama[undato] = datosensor[undato]

                # datos restantes
                datostrama["frequency"] = freq_tx
                datostrama["bandwidth"] = bandwidth
                datostrama["spreadingFactor"] = spreadingFactor
                datostrama["fPort"]     = fPort
                datostrama["dr"]        = dr
                datostrama["created"]   = unafecha

                leidos = leidos + 1
                
                # revisa registro repetido
                repetido = False
                if leidos>1:
                    trama_num  = (fCnt == tabla[leidos-1]['fCnt'])
                    fecha_pub  = (publishedAt == tabla[leidos-1]['publishedAt'])
                    gtwNum_rep = (gtwNum == tabla[leidos-1]['gtw_rx'])
                    repetido   = (trama_num and fecha_pub and gtwNum_rep)
                if not(repetido):
                    tabla[leidos] = datostrama
                else:
                    leidos = leidos - 1
                
        # convierte diccionario en DataFrame
        if len(tabla)>0:
            tabla = pd.DataFrame(tabla)
            tabla = tabla.T
    return(tabla)

def describe_punto(tabla,codigopunto,carpeta_rsm,variables):
    ''' Descriptores estadisticos de los datos.csv
        de un dispositivo, revisa media, desviació estándar, pmf
        graba unreporte.csv con pandas y genera gráficas
    '''
    medida = variables['medida']
    medida_unidad = variables['medida_unidad']
    medida_normal = variables['medida_normal']
    movAvg_cual = variables['movAvg_cual']
    medida_grafica = variables['medida_grafica']
    movAvg_color = variables['movAvg_color']
    guarda = variables['guarda']
    precision = variables['precision']
                   
    fechaformato = "%Y-%m-%d %H:%M:%S.%f"
    # medida intervalo
    medida_min = np.min(medida_normal)
    medida_max = np.max(medida_normal)

    # fechas series a datetime
    tabla['publishedAt'] = pd.to_datetime(tabla['publishedAt'],
                                          format=fechaformato)
    tabla['created'] = pd.to_datetime(tabla['publishedAt'],
                                      format=fechaformato)

    # revisa errores de medida
    tabla["error_up"]   = 0
    tabla["error_down"] = 0
    for undato in tabla['publishedAt'].keys():   
        medida_up = tabla[medida+'_up'][undato]
        enrango = (medida_up>=np.min(medida_normal))
        enrango = (enrango and medida_up<=np.max(medida_normal))
        if not(enrango):
            tabla.at[undato,"error_up"] = 1
        
        medida_down = tabla[medida+'_down'][undato]
        enrango = (medida_down>=np.min(medida_normal))
        enrango = (enrango and medida_down<=np.max(medida_normal))
        if not(enrango):
            tabla.at[undato,"error_down"] = 1      

    # tasa error trama
    leidos = len(tabla)
    if leidos > 0:
        error_up   = np.sum(tabla['error_up'])
        error_up   = error_up/leidos
        error_down = np.sum(tabla['error_down'])
        error_down = error_down/leidos

    # descriptor estadístico, datos sin errores
    condicion_up = (tabla['error_up']==0)
    condicion_down = (tabla['error_down']==0)

    medida_up = tabla[condicion_up][medida+'_up']
    describe_up = medida_up.describe()
    describe_up['error_trama'] = error_up

    medida_down = tabla[condicion_down][medida+'_down']
    describe_down = medida_down.describe()
    describe_down['error_trama'] = error_down

    descriptor = describe_up.copy()
    descriptor = pd.concat([descriptor,describe_down],axis=1)
    descriptor['dispositivo'] = tabla['dispositivo'][0]

    pmf_up   = medida_pmf(medida_up,describe_up)
    pmf_down = medida_pmf(medida_down,describe_down)

    pmf_punto = {'pmf':{'pmf_up'   : pmf_up,
                        'pmf_down' : pmf_down}}
    pmf_punto = pd.DataFrame(pmf_punto)
    pmf_punto['dispositivo'] = tabla['dispositivo'][0]

    # Para gráficas
    # medias moviles en movAvg_cual[]
    serie_up  = pd.Series(medida_up)
    movAvg_up_mean = []
    movAvg_up_std = []
    m = len(movAvg_cual)
    for j in range(0,m,1):
        k = movAvg_cual[j]
        movAvg_up_mean.append(list(serie_up.rolling(k).mean()))
        movAvg_up_std.append(list(serie_up.rolling(k).std()))
        
    serie_down = pd.Series(medida_down)
    movAvg_down_mean = []
    movAvg_down_std = []
    for j in range(0,m,1):
        k = movAvg_cual[j]
        movAvg_down_mean.append(list(serie_down.rolling(k).mean()))
        movAvg_down_std.append(list(serie_down.rolling(k).std()))

    movAvgData ={'movAvg_cual'   : movAvg_cual,
                 'movAvg_color'  : movAvg_color,
                 'movAvg_up_mean'  : movAvg_up_mean,
                 'movAvg_down_mean': movAvg_down_mean,
                 'movAvg_up_std'   : movAvg_up_std,
                 'movAvg_down_std' : movAvg_down_std
                 }

    grafData ={'codigopunto' : codigopunto,
               'medida' : medida,
               'precision': precision,
               'medida_unidad' : medida_unidad,
               'medida_grafica': medida_grafica
               }
    resultado = [tabla,descriptor,pmf_punto,movAvgData,grafData]
    return(resultado)

# función de probabilidad de masa pmf
def medida_pmf(valores,undescriptor):
    pmin   = np.min(valores)
    pmax   = np.max(valores)
    tramo  = int(pmax-pmin)
    conteo = np.zeros(tramo+1,dtype=int)
    intervalo = np.arange(pmin,pmax+1,1)
    for valor in valores:
        donde = np.where(intervalo == valor)
        conteo[donde] = conteo[donde] + 1
    freq_relativa = np.array(conteo)/np.sum(conteo)
    unpmf = {'intervalo' : list(intervalo),
             'freq_relativa' : list(freq_relativa)}
    return(unpmf)

# GRAFICA -----
def graf_puntos_serie(tabla,descriptor,movAvgData,grafData):
    ''' grafica la serie de tiempo de cada punto
        añade medias móviles
    '''
    # ajuste de formato de fecha para eje x
    converter = mdates.ConciseDateConverter()
    munits.registry[np.datetime64] = converter
    munits.registry[dt.date] = converter
    munits.registry[dt.datetime] = converter

    # datos para grafica
    precision   = grafData['precision']
    medida = grafData['medida']
    medida_unidad  = grafData['medida_unidad']
    medida_grafica = grafData['medida_grafica']
    
    movAvg_cual  = movAvgData['movAvg_cual']
    movAvg_color = movAvgData['movAvg_color']
    
    media_up    = descriptor[medida+'_up']['mean']
    std_up      = descriptor[medida+'_up']['std']
    media_down  = descriptor[medida+'_down']['mean']
    std_down    = descriptor[medida+'_down']['std']

    # ajuste de intervalo eje y
    y_min = np.min([np.min(medida_grafica),
                    media_up - 2*std_up,
                    media_down - 2*std_down])
    y_max = np.max([np.max(medida_grafica),
                    media_up + 2*std_up,
                    media_down + 2*std_down])
    
    # selecciona sin error
    condicion_up = (tabla['error_up']==0)
    condicion_down = (tabla['error_down']==0)

    # grafica
    fig_serie,(graf_up,graf_down) = plt.subplots(2,1)
    
    # medida_up -----
    graf_up.plot(tabla[condicion_up]['publishedAt'],
                 tabla[condicion_up][medida+'_up'],
                 color='blue',marker ='.',
                 linestyle='')
    
    # medida_up, medias y std
    etiq_up = str(np.round(media_up,precision))+' +/- '
    etiq_up = etiq_up + str(np.round(std_up,precision))
    graf_up.axhline(media_up,
                    color='blue',label=etiq_up)
    graf_up.axhline(media_up-std_up,
                    color='blue',linestyle='dotted')
    graf_up.axhline(media_up+std_up,
                    color='blue',linestyle='dotted')
    
    # medida_up, medias móviles
    m = len(movAvg_cual)
    for j in range(0,m,1):
        k = str(movAvg_cual[j])
        graf_up.plot(tabla[condicion_up]['publishedAt'],
                     movAvgData['movAvg_up_mean'][j],
                     label='movAvg_'+k,
                     color=movAvg_color[j])
    
    graf_up.set_ylim(y_min,y_max)
    graf_up.set_ylabel(medida+'_up ('+medida_unidad+')',
                       color='blue')
    graf_up.legend()
    graf_up.grid(True,linestyle='dotted',
                 axis='x',which='both')

    # medida_down -------
    graf_down.plot(tabla[condicion_down]['publishedAt'],
                   tabla[condicion_down][medida+'_down'],
                   color='brown',marker ='.',
                   linestyle='')

    # medida_down, medias y std
    etiq_down = str(np.round(media_down,precision))+' +/- '
    etiq_down = etiq_down + str(np.round(std_down,precision))
    graf_down.axhline(media_down,
                      color='brown',label=etiq_down)
    graf_down.axhline(media_down+std_down,
                      color='brown',linestyle='dotted')
    graf_down.axhline(media_down-std_down,
                      color='brown',linestyle='dotted')
    
    # medida_down, medias moviles
    for j in range(0,m,1):
        k = str(movAvg_cual[j])
        graf_down.plot(tabla[condicion_down]['publishedAt'],
                       movAvgData['movAvg_down_mean'][j],
                       label='movAvg_'+k,
                       color=movAvg_color[j])
    
    graf_down.set_ylim(y_min,y_max)
    graf_down.set_xlabel('fecha')
    graf_down.set_ylabel(medida+'_down ('+medida_unidad+')',
                         color='brown')
    graf_down.legend()
    graf_down.grid(True,linestyle='dotted',
                   axis='x', which='both')
    graf_up.set_title('Serie: '+grafData['codigopunto']+' '+ medida)
    plt.tight_layout()
    return(fig_serie)

def graf_puntos_pmf(pmf_punto,descriptor,grafData):
    ''' grafica función de probabilida de masa
        para cada punto, media +/- std
    '''
    # datos para grafica
    x_pmfup   = pmf_punto['pmf']['pmf_up']['intervalo']
    y_pmfup   = pmf_punto['pmf']['pmf_up']['freq_relativa']
    x_pmfdown = pmf_punto['pmf']['pmf_down']['intervalo']
    y_pmfdown = pmf_punto['pmf']['pmf_down']['freq_relativa']
    
    precision = grafData['precision']
    medida = grafData['medida']
    medida_unidad  = grafData['medida_unidad']
    medida_grafica = grafData['medida_grafica']
    
    media_up   = descriptor[medida+'_up']['mean']
    std_up     = descriptor[medida+'_up']['std']
    media_down = descriptor[medida+'_down']['mean']
    std_down   = descriptor[medida+'_down']['std']

    prob_max = 0.40
    # ajuste de intervalo eje y
    y_min = np.min([np.min(medida_grafica),
                    media_up - 2*std_up,
                    media_down - 2*std_down])
    y_max = np.max([np.max(medida_grafica),
                    media_up + 2*std_up,
                    media_down + 2*std_down])
    # grafica
    fig_pmf,graf_pmf = plt.subplots()
    etiq_up = str(np.round(media_up,precision)) +' +/- '
    etiq_up = etiq_up + str(np.round(std_up,precision))
    graf_pmf.plot(x_pmfup,y_pmfup,
                  label='media_up '+etiq_up,
                  color='blue')
    graf_pmf.axvline(media_up,color='blue')
    graf_pmf.axvline(media_up+std_up,
                     linestyle='dotted',color='blue')
    graf_pmf.axvline(media_up-std_up,
                     linestyle='dotted',color='blue')

    etiq_down = str(np.round(media_down,precision))+' +/- '
    etiq_down = etiq_down + str(np.round(std_down,precision))
    graf_pmf.plot(x_pmfdown,y_pmfdown,
                  label='media_down '+etiq_down,
                  color='brown')
    graf_pmf.axvline(media_down,color='brown')
    graf_pmf.axvline(media_down+std_down,
                     linestyle='dotted',color='brown')
    graf_pmf.axvline(media_down-std_down,
                     linestyle='dotted',color='brown')

    graf_pmf.set_title('pmf: '+grafData['codigopunto']+' '+medida)
    graf_pmf.set_xlim(y_min,y_max)
    graf_pmf.set_ylim(0,prob_max)
    graf_pmf.set_xlabel(medida+' ('+medida_unidad+')')
    graf_pmf.set_ylabel('frecuencia relativa')
    graf_pmf.legend()
    graf_pmf.grid(True,linestyle='dotted',
                  axis='x', which='both')
    return(fig_pmf)

def graf_puntos_std(tabla,descriptor,movAvgData,grafData):
    ''' grafica serie de std usando medias moviles
        para cada punto, media_std
    '''
    # ajuste de formato de fecha para eje x
    converter = mdates.ConciseDateConverter()
    munits.registry[np.datetime64] = converter
    munits.registry[dt.date] = converter
    munits.registry[dt.datetime] = converter
    
    # datos para grafica
    precision = grafData['precision']
    medida    = grafData['medida']
    medida_unidad = grafData['medida_unidad']

    movAvg_cual = movAvgData['movAvg_cual']
    movAvg_color = movAvgData['movAvg_color']

    # selecciona sin error
    condicion_up   = (tabla['error_up']==0)
    condicion_down = (tabla['error_down']==0)
    
    # ajuste de intervalo eje y
    y_min = 0
    y_max = np.max([2, 2*descriptor[medida+'_up']['std'],
                    2*descriptor[medida+'_down']['std']])
    # grafica
    fig_std,(graf_stdUp,graf_stdDown) = plt.subplots(2,1)
    
    # std up
    std_up = np.round(descriptor[medida+'_up']['std'],precision)
    graf_stdUp.axhline(std_up,label='std '+str(std_up),
                       color='blue')
    m = len(movAvg_cual)
    for j in range(0,m,1):
        k = str(movAvg_cual[j])
        graf_stdUp.plot(tabla[condicion_up]['publishedAt'],
                        movAvgData['movAvg_up_std'][j],
                        label='movAvg_'+k,
                        color=movAvg_color[j])
    graf_stdUp.set_ylim(y_min,y_max)
    graf_stdUp.set_ylabel('std_up ('+medida_unidad+')',
                          color='blue')
    graf_stdUp.legend()
    graf_stdUp.grid(True,linestyle='dotted',
                    axis='x', which='both')
    graf_stdUp.set_title('std: '+grafData['codigopunto']+' '+ medida)

    # std down
    std_down = np.round(descriptor[medida+'_down']['std'],precision)
    graf_stdDown.axhline(std_down,label='std '+str(std_down),
                         color='brown')
    for j in range(0,m,1):
        k = str(movAvg_cual[j])
        graf_stdDown.plot(tabla[condicion_down]['publishedAt'],
                          movAvgData['movAvg_down_std'][j],
                          label='movAvg_'+k,color=movAvg_color[j])
    graf_stdDown.set_ylim(y_min,y_max)
    graf_stdDown.set_xlabel('fecha')
    graf_stdDown.set_ylabel('std_down ('+medida_unidad+')',
                            color='brown')
    graf_stdDown.legend()
    graf_stdDown.grid(True,linestyle='dotted',
                      axis='x', which='both')
    plt.tight_layout()
    return(fig_std)


def linealiza_lstsq(xi,yi,digitos = 3):
    ''' usa minimos cuadrados para entregar la ecuacion
        digitos: usados en expresion latex
    '''
    unaecuacion = {}
    # Eje x en log10()
    xilog = np.log10(xi)
    n = len(xi)
    
    # mínimos cuadrados (least square),
    # distancia vs medida
    A = np.vstack([xilog, np.ones(n)]).T
    [m0, b0] = np.linalg.lstsq(A, yi, rcond=None)[0]
    alpha = -m0/10
    beta  = b0

    # coeficiente de correlación
    coeficientes = np.corrcoef(xilog,yi)
    r = coeficientes[0,1]
    # coeficiente de determinación
    r2 = r**2

    # ecuaciones expresion rssi(d)
    fdist0 = lambda d: -10*alpha*(np.log10(d))+beta
    
    fdtxt0 = r'$ rssi = -10(' + str(np.round(alpha,digitos))
    fdtxt0 = fdtxt0 + ')log_{10}(d)' # +('
    texto = '+'
    if beta <0:
        texto = '-'
    fdtxt0 = fdtxt0 + texto + str(np.round(np.abs(beta),digitos))+' $'

    # Errores respecto a rssi(d) 
    yi0  = fdist0(xi)
    dyi0 = yi - yi0
    dyi0mean = np.mean(np.abs(dyi0))
    dyi0std  = np.std(dyi0, dtype=np.float64)

    # ecuaciones expresion d(rssi)
    grssi0 = lambda rssi: 10**((beta-rssi)/(10*alpha))
    grtxt0 = r"$ d = 10^{(" + str(np.round(beta,digitos)) + ' - '
    grtxt0 = grtxt0 + 'rssi)/' + '(10('+str(np.round(alpha,digitos))+'))} $'

    # Errores respecto a rssi(d) 
    xi0  = grssi0(yi)
    dxi0 = xi - xi0
    dxi0mean = np.mean(np.abs(dxi0))
    dxi0std  = np.std(dxi0, dtype=np.float64)
    
    yimin = np.around(np.min(yi),2)
    yimax = np.around(np.max(yi),2)
    
    unaecuacion = {'alpha'   : alpha,
                   'beta'    : beta,
                   'coef_correlacion'   : r,
                   'coef_determinacion' : r2,
                   'eq_latex': fdtxt0,
                   'intervalox' : [np.min(xi),np.max(xi)],
                   'error_medio': dyi0mean,
                   'error_std'  : dyi0std,
                   'eqg_latex'  : grtxt0,
                   'intervaloy' : [yimin,yimax],
                   'errorx_medio': dxi0mean,
                   'errorx_std'  : dxi0std,
                   }
    return(unaecuacion)

def ecuacion_gradiente(carpeta_rsm,arch_coord,var_gen):
    ''' desarrolla la ecuación para un gradiente
    '''
    var_gen['movAvg_cual'] = [2,4] #cada cuantas muestras
    partes = carpeta_rsm.strip('/').split('_')
    arch_nombre = partes[1]+'_'+partes[2]

    carp_coord = var_gen['carp_coord']
    medida = var_gen['medida']
    precision = var_gen['precision']
    
    # lee coordenadas y su distancia al gateway
    dist_Gtw = coordenadas_leer(arch_coord,carp_coord)

    # crea lista.txt de archivos a usar si no existe
    arch_lista = medida+"_"+arch_nombre+"_lista.txt"
    archivo_ruta  = carpeta_rsm + '/' + arch_lista
    archivoexiste = os.path.exists(archivo_ruta)
    if not(archivoexiste):
        puntoUsar = pd.DataFrame(columns=['punto', 'up', 'down'])
        for unarchivo in os.listdir(carpeta_rsm):
            verifica = unarchivo.startswith('describe_')
            verifica = verifica and unarchivo.endswith('.csv')
            if verifica:
                partes  = unarchivo.strip('.csv').split('_')
                unpunto = partes[2]
                puntoUsar = puntoUsar.append({'punto' : unpunto,
                                              'up' : 1,
                                              'down' : 1},
                                             ignore_index = True)
        puntoUsar = puntoUsar.set_index('punto')
        puntoUsar.to_csv(archivo_ruta)

    # usa lista.txt de archivos seleccionados 1 ó 0 para enlace up,down
    if (archivoexiste):
        puntoUsar = pd.read_csv(archivo_ruta)
        puntoUsar = puntoUsar.set_index('punto')
      
    # Datos para gráfica desde lista.txt
    punto_graf  = pd.DataFrame()
    puntoSinDist = []
    for unarchivo in os.listdir(carpeta_rsm):
        if unarchivo.startswith('describe_'):
            codigopunto = unarchivo.strip('.csv').split('_')[2]
            
            # lectura del archivo
            unresumen  = carpeta_rsm+'/'+unarchivo
            descriptor = pd.read_csv(unresumen,index_col='Unnamed: 0')

            if (codigopunto in dist_Gtw.index):
                undato = {'codigopunto': codigopunto,
                          'dist_Gtw'   : dist_Gtw[codigopunto],
                          medida+'_up' : descriptor['rssi_up']['mean'],
                          medida+'_up'+'_std': descriptor['rssi_up']['std'],
                          'usar_up' : puntoUsar['up'][codigopunto],
                          medida+'_down': descriptor['rssi_down']['mean'],
                          medida+'_down'+'_std' :descriptor['rssi_down']['std'],
                          'usar_down' : puntoUsar['down'][codigopunto],
                          'dispositivo': descriptor['dispositivo'][0]
                          }
                punto_graf = punto_graf.append(undato,ignore_index=True)
            else:
                puntoSinDist.append(codigopunto)
            
    punto_graf = punto_graf.set_index('codigopunto')
    punto_graf = punto_graf.sort_values('dist_Gtw' )

    # selecciona datos
    dist = punto_graf['dist_Gtw']
    usar_up = punto_graf['usar_up']
    usar_down = punto_graf['usar_down']

    # Eje x en log10()
    xi = np.array(dist)

    # enlace_up
    media_up = punto_graf[medida+'_up']
    std_up   = punto_graf[medida+'_up'+'_std']
    media_up_techo = media_up + std_up
    media_up_piso  = media_up - std_up

    # ecuacion enlace up
    yi_up = np.array(media_up)
    eq_up = linealiza_lstsq(xi[usar_up==1],
                            yi_up[usar_up==1],
                            digitos = precision)
    alpha_up = eq_up['alpha']
    beta_up  = eq_up['beta']
    fd_up = lambda d: -10*alpha_up*(np.log10(d))+beta_up

    # errores up de puntos vs ecuacion
    fi_up = fd_up(xi)
    dyi0std  = eq_up['error_std']
    punto_graf['fi_up'] = fi_up

    # enlace_down
    media_down = punto_graf[medida+'_down']
    std_down   = punto_graf[medida+'_down'+'_std']
    media_down_techo = media_down + std_down
    media_down_piso  = media_down - std_down

    # ecuación Enlace down
    yi_down = np.array(media_down)
    eq_down = linealiza_lstsq(xi[usar_down==1],
                              yi_down[usar_down==1],
                              digitos = precision)
    alpha_down = eq_down['alpha']
    beta_down  = eq_down['beta']
    fd_down = lambda d: -10*alpha_down*(np.log10(d))+beta_down

    # errores down de puntos vs ecuacion
    fi_down = fd_down(xi)
    dyi1std = eq_down['error_std']
    punto_graf['fi_down'] = fi_down

    resultado = {'tabla'         : punto_graf,
                 'ecuacion_up'   : eq_up,
                 'ecuacion_down' : eq_down,
                 'puntoSinDist'  : puntoSinDist}
    
    return (resultado)

def coordenadas_leer(arch_coord,carp_coord):
    ''' lista las coordenadas desde un archivo
        ubicado en la carpeta
    '''
    if len(carp_coord)>0:
        carp_coord = carp_coord+'/'
    archivo_ruta = carp_coord+arch_coord
    arch_coord_existe = os.path.exists(archivo_ruta)

    if arch_coord_existe:
        puntos = pd.read_csv(archivo_ruta,index_col='punto')
        # busca gateway y ubicación de punto observado
        Gw_nombre = ''
        for unpunto in puntos.index:
            if unpunto.startswith('Gw'):
                Gw_nombre = unpunto
        dist_Gtw = puntos['dist_'+Gw_nombre]   
    else:
        print(' ERROR: NO se encuentra el archivo...')
        print('        revise el ruta/nombre de archivo. ')
    return(dist_Gtw)

def graf_gradiente(arch_nombre,punto_graf,ecuacion,var_gen):
    ''' grafica de un gradiente
    '''
    escalabase = 10    # 10
    escala = 'log'
    medida = var_gen['medida']
    precision = var_gen['precision']

    # enlace_up
    media_up = punto_graf[medida+'_up']
    std_up   = punto_graf[medida+'_up'+'_std']
    media_up_techo = media_up + std_up
    media_up_piso  = media_up - std_up
    yi_up = np.array(media_up)
    fi_up = punto_graf['fi_up']

    # enlace_down
    media_down = punto_graf[medida+'_down']
    std_down   = punto_graf[medida+'_down'+'_std']
    media_down_techo = media_down + std_down
    media_down_piso  = media_down - std_down
    yi_down = np.array(media_down)
    fi_down = punto_graf['fi_down']

    dist = punto_graf['dist_Gtw']
    usar_up = punto_graf['usar_up']
    usar_down = punto_graf['usar_down']

    eq_up = ecuacion['ecuacion_up']
    eq_down = ecuacion['ecuacion_down']

    dyi0std = eq_up['error_std']
    dyi1std = eq_down['error_std']

    fig_gradnt,(graf_up,graf_down) = plt.subplots(2,1)
    if escala == 'log':
        graf_up.set_xscale(escala,base=escalabase)
        graf_down.set_xscale(escala,base=escalabase)

    # medida up +/- std
    graf_up.scatter(dist,media_up,color='blue',marker='o')
    graf_up.scatter(dist,media_up_techo,color='blue',marker='+')
    graf_up.scatter(dist,media_up_piso,color='blue',marker='+')

    # medida down +/- std
    graf_down.scatter(dist,media_down,color='orange',marker='o')
    graf_down.scatter(dist,media_down_techo,color='orange',marker='+')
    graf_down.scatter(dist,media_down_piso,color='orange',marker='+')

    # linealizado up
    graf_up.plot(dist,fi_up,label = '  up:'+eq_up['eq_latex'],
                 color='blue')
    graf_up.plot(dist,fi_up+dyi0std,
                 color='blue',linestyle='dotted')
    graf_up.plot(dist,fi_up-dyi0std,
                 color='blue',linestyle='dotted')
    # linealizado down
    graf_down.plot(dist,fi_down,label = 'down:'+eq_down['eq_latex'],
                   color='orange')
    graf_down.plot(dist,fi_down+dyi1std,
                   color='orange',linestyle='dotted')
    graf_down.plot(dist,fi_down-dyi1std,
                   color='orange',linestyle='dotted')

    # ajuste de eje y
    y_min = np.min([np.min(media_up_piso),np.min(media_down_piso)])
    y_max = np.max([np.max(yi_up+dyi0std),np.max(yi_down+dyi1std)])
    if y_min<-135:
        y_min = np.min([np.min([media_up]),np.min([media_down])])
    graf_up.set_ylim(y_min,y_max)
    graf_down.set_ylim(y_min,y_max)

    # etiquetas up
    for unpunto in punto_graf.index:
        media_etiq = np.round(punto_graf[medida+'_up'][unpunto],precision)
        usarpunto  = punto_graf['usar_up'][unpunto]
        
        etiqueta = unpunto #+str(media_etiq)
        xi = punto_graf['dist_Gtw'][unpunto]
        yi = punto_graf['rssi_up'][unpunto]
        if usarpunto:
            graf_up.annotate(etiqueta,(xi,yi), rotation=45)
        if not(usarpunto):
            graf_up.annotate(etiqueta,(xi,yi), rotation=45,color='red')
            graf_up.scatter(dist[unpunto],media_up[unpunto],color='red',marker='o')
            graf_up.scatter(dist[unpunto],media_up_techo[unpunto],
                            color='red',marker='+')

    # etiquetas down
    for unpunto in punto_graf.index:
        media_etiq = np.round(punto_graf[medida+'_down'][unpunto],precision)
        usarpunto = punto_graf['usar_down'][unpunto]
        etiqueta = unpunto #+str(media_etiq)
        xi = punto_graf['dist_Gtw'][unpunto]
        yi = punto_graf['rssi_down'][unpunto]
        if usarpunto:
            graf_down.annotate(etiqueta,(xi,yi), rotation=45)
        if not(usarpunto):
            graf_down.annotate(etiqueta,(xi,yi), rotation=45,color='red')
            graf_down.scatter(dist[unpunto],media_down[unpunto],color='red',marker='o')
            graf_down.scatter(dist[unpunto],media_down_techo[unpunto],
                              color='red',marker='+')

    graf_up.set_ylabel(medida+'_up (dBm)',color='blue')
    graf_up.set_title(arch_nombre+' '+medida)
    graf_up.legend()
    graf_up.grid(True,linestyle='dotted',
                 axis='x', which='both')

    graf_down.set_xlabel('distancia')
    graf_down.set_ylabel(medida+'_down (dBm)',color='brown')
    graf_down.legend()
    graf_down.grid(True,linestyle='dotted',
                   axis='x', which='both')
    return(fig_gradnt)

def graf_diferencia_updown(arch_nombre,punto_graf,var_gen):
    '''Grafica de diferencias up-down
    '''
    medida = var_gen['medida']
    medida_unidad = var_gen['medida_unidad']
    movAvg_cual = var_gen['movAvg_cual']
    movAvg_color = var_gen['movAvg_color']
    precision = var_gen['precision']

    # analiza diferencias de medias y std --------------
    punto_graf['usar_dif'] = punto_graf['usar_up']*punto_graf['usar_down']
    usar_dif = punto_graf['usar_dif']
    punto_graf['dif_mean'] = punto_graf['rssi_up']-punto_graf['rssi_down']
    punto_graf['dif_std']  = punto_graf[medida+'_up_std']-punto_graf[medida+'_down_std']
    punto_graf['dif_top']  = punto_graf['dif_mean']+punto_graf['dif_std']
    punto_graf['dif_bottom'] = punto_graf['dif_mean']-punto_graf['dif_std']

    dist = punto_graf['dist_Gtw']
    usar_dif = punto_graf['usar_dif']

    # linealiza diferencia
    xi_dif = np.array(punto_graf['dist_Gtw'][usar_dif==1])
    yi_dif = np.array(punto_graf['dif_mean'][usar_dif==1])
    eq_dif = linealiza_lstsq(xi_dif,yi_dif,digitos = precision)
    alpha_dif = eq_dif['alpha']
    beta_dif  = eq_dif['beta']
    fdistdif0 = lambda d: -10*alpha_dif*(np.log10(d))+beta_dif
    yi_dif0   = fdistdif0(xi_dif)

    ecua_dif = pd.DataFrame(eq_dif)

    # medias moviles en movAvg_cual[]
    serie_media  = pd.Series(punto_graf['dif_mean'][usar_dif==1])
    movAvg_dmedia = []
    m = len(movAvg_cual)
    for j in range(0,m,1):
        k = movAvg_cual[j]
        movAvg_dmedia.append(serie_media.rolling(k).mean())

    # Grafica diferencia
    escalabase = 10    # 10
    escala = 'log'
    
    fig_dif_updown,graf_dmean = plt.subplots()
    if escala == 'log':
        graf_dmean.set_xscale(escala,base=escalabase)

    # diferencia medida up - down
    graf_dmean.scatter(dist,punto_graf['dif_mean'],color='purple',marker='o')
    graf_dmean.plot(xi_dif,yi_dif0,label = '  dif:'+eq_dif['eq_latex'],color='blue')
    graf_dmean.set_ylabel('media Up-Down (dBm)')

    # medida_up, medias móviles
    for j in range(0,m,1):
        k = str(movAvg_cual[j])
        graf_dmean.plot(dist[usar_dif==1],movAvg_dmedia[j],
                      label='movAvg_dmedia_'+k,color=movAvg_color[j])

    # etiquetas
    i=0
    for unpunto in punto_graf.index:
        media_etiq = np.round(punto_graf[medida+'_up'][unpunto],precision)
        etiqueta = unpunto # + ' ' +str(media_etiq)
        usar_dif = punto_graf['usar_dif'][unpunto]
        xi = punto_graf['dist_Gtw'][unpunto]
        yi = punto_graf.iloc[i]['dif_mean']
        i = i+1
        if usar_dif:
            graf_dmean.annotate(etiqueta,(xi,yi), rotation=30)
        if not(usar_dif):
            graf_dmean.annotate(etiqueta,(xi,yi),color='red',rotation=30)
            graf_dmean.scatter(xi,yi,color='red',marker='o')

    graf_dmean.axhline(0,color='grey')
    graf_dmean.grid(True,linestyle='dotted',axis='x', which='both')

    graf_dmean.legend()
    graf_dmean.set_title(arch_nombre+': '+medida+' media Up-Down')
    return(fig_dif_updown)

def graf_gradiente_std(arch_nombre,punto_graf,medida,precision):
    ''' Grafica de std_up y std_down ------
    '''
    xi_std      = np.array(punto_graf['dist_Gtw'])
    yi_std_up   = np.array(punto_graf[medida+'_up_std'])
    yi_std_down = np.array(punto_graf[medida+'_down_std'])

    # linealiza std
    usar_up   = punto_graf['usar_up']
    usar_down = punto_graf['usar_down']

    eq_std_up = linealiza_lstsq(xi_std[usar_up==1],
                                      yi_std_up[usar_up==1],
                                      digitos = precision)
    alpha_std_up = eq_std_up['alpha']
    beta_std_up  = eq_std_up['beta']
    fdiststd_up0 = lambda d: -10*alpha_std_up*(np.log10(d))+beta_std_up
    yi_std_up0   = fdiststd_up0(xi_std)

    eq_std_down = linealiza_lstsq(xi_std[usar_down==1],
                                        yi_std_down[usar_down==1],
                                        digitos = precision)
    alpha_std_down = eq_std_down['alpha']
    beta_std_down  = eq_std_down['beta']
    fdiststd_down0 = lambda d: -10*alpha_std_down*(np.log10(d))+beta_std_down
    yi_std_down0   = fdiststd_down0(xi_std)

    # grafica std
    escalabase = 10    # 10
    escala = 'log'
    
    fig_grad_std,(graf_std_up,graf_std_down) = plt.subplots(2,1)
    if escala == 'log':
        graf_std_up.set_xscale(escala,base=escalabase)
        graf_std_down.set_xscale(escala,base=escalabase)

    graf_std_up.plot(xi_std,yi_std_up0,
                     label = '  std_up:'+eq_std_up['eq_latex'],
                     color='blue')
        
    # medida up +/- std
    graf_std_up.scatter(xi_std,yi_std_up,color='blue',marker='+')
    graf_std_down.scatter(xi_std,yi_std_down,color='orange',marker='+')

    # etiquetas
    i=0
    for unpunto in punto_graf.index:
        media_etiq = np.round(punto_graf[medida+'_up'][unpunto],precision)
        etiqueta   = unpunto # + ' ' +str(media_etiq)
        usar_up    = punto_graf['usar_up'][unpunto]
        usar_down  = punto_graf['usar_down'][unpunto]
        xi  = punto_graf['dist_Gtw'][unpunto]
        yiu = punto_graf.iloc[i][medida+'_up_std']
        yid = punto_graf.iloc[i][medida+'_down_std']
        i = i+1
        if usar_up:
            graf_std_up.annotate(etiqueta,(xi,yiu), rotation=30)
        if not(usar_up):
            graf_std_up.annotate(etiqueta,(xi,yiu),color='red',rotation=30)
            graf_std_up.scatter(xi,yiu,color='red',marker='o')
        if usar_down:
            graf_std_down.annotate(etiqueta,(xi,yid), rotation=30)
        if not(usar_down):
            graf_std_down.annotate(etiqueta,(xi,yid),color='red',rotation=30)
            graf_std_down.scatter(xi,yid,color='red',marker='o')

    graf_std_down.plot(xi_std,yi_std_down0,
                    label = 'std_down:'+eq_std_down['eq_latex'],
                    color='orange')
    # otras etiquetas
    graf_std_up.set_ylabel('std_up (dBm)', color='blue')
    graf_std_up.set_xlabel('distancia')
    graf_std_up.axhline(0,color='grey')
    graf_std_up.grid(True,linestyle='dotted',
                    axis='x', which='both')
    graf_std_up.set_title(arch_nombre+': '+medida+' std')
    graf_std_up.legend()

    graf_std_down.set_ylabel('std_down (dBm)', color='brown')
    graf_std_down.set_xlabel('distancia')
    graf_std_down.axhline(0,color='grey')
    graf_std_down.grid(True,linestyle='dotted',
                    axis='x', which='both')
    graf_std_down.legend()

    return(fig_grad_std)

def coordenadas_kml_utm(archivo_ruta,zonaNum,zonaLetra,precision):
    # Lectura de archivo
    gpd.io.file.fiona.drvsupport.supported_drivers['KML'] = 'rw'
    puntos  = gpd.read_file(archivo_ruta, driver='KML')
    puntos['dist'] = 0

    tabla = pd.DataFrame()
    for i in puntos.index:
        # en latitud, longitud y altura
        punto_nombre = puntos.loc[i]['Name']
        longitud = puntos.loc[i]['geometry'].x
        latitud  = puntos.loc[i]['geometry'].y
        altitud  = puntos.loc[i]['geometry'].z
        altitud  = np.round(altitud,precision)
        #en UTM
        coord_utm  = utm.from_latlon(latitud,longitud,
                                     zonaNum,zonaLetra)
        utm_este  = np.round(coord_utm[0],precision)
        utm_norte = np.round(coord_utm[1],precision)
        
        # distancias a Gateways
        dist = 0.0
        if punto_nombre.startswith('Gw'):
            n_Gtw = punto_nombre
            x1 = utm_este
            y1 = utm_norte
        else:
            x2 = utm_este
            y2 = utm_norte
            dist  = np.sqrt((x2-x1)**2 + (y2-y1)**2)
        dist = np.round(dist,precision)
        
        unpunto = {'punto'       : punto_nombre,
                   'utm_este'    : utm_este,
                   'utm_norte'   : utm_norte,
                   'utm_zNum'    : int(zonaNum),
                   'utm_zLetra'  : zonaLetra,
                   'altitud'     : altitud,
                   'altitud_gps' : altitud,
                   'dist_'+n_Gtw : dist,
                   'latitud'     : latitud,
                   'longitud'    : longitud
                   }
        tabla = tabla.append(unpunto,ignore_index=True)
    tabla = tabla.set_index('punto')
    return(tabla)

def revisa_perfil(archivo_ruta,alturaGw,alturapunto,
                  plantas,muestras=41,guarda=False):
    LineaObs = archivo_ruta.strip('.csv').split('_')[1]

    planta_desde  = plantas[0]
    planta_hasta  = plantas[1]
    planta_altura = plantas[2]

    # lectura de datos
    puntos = pd.read_csv(archivo_ruta,index_col='punto')

    # busca gateway
    n_Gtw = ''
    for unpunto in puntos.index:
        if unpunto.startswith('Gw'):
            n_Gtw = unpunto

    # perfil  
    xi = np.array(puntos['dist_'+n_Gtw])
    yi = np.array(puntos['altitud'])
    etiqueta = puntos.index

    # plantaciones
    conplantas = (xi > planta_desde) & (xi < planta_hasta)
    plantacion = yi[conplantas]+planta_altura

    # GRAFICA
    # plantacion
    plt.fill_between(xi[conplantas],
                     yi[conplantas],
                     plantacion,
                     color='lightgreen')
    # perfil
    xi = puntos['dist_'+n_Gtw]
    yi = puntos['altitud']
    plt.plot(xi,yi,color='brown')
        
    # antenas
    for unpunto in puntos.index:
        xip = xi[unpunto]
        yip = yi[unpunto]
        if unpunto.startswith('Gw'):
            xigtw = xip
            yigtw = yip+alturaGw
            plt.scatter(xip ,yigtw)
            plt.annotate(unpunto,(xigtw,yigtw))
        else:
            xid = xip
            yid = yip+alturapunto
            plt.plot((xigtw,xid),(yigtw,yid),
                     color = 'green',
                     linestyle='dotted')
            plt.scatter(xid ,yid)
            plt.annotate(unpunto,(xid,yid), rotation=45)
    plt.title('Linea :'+LineaObs)
    plt.xlabel('distancia')
    plt.ylabel('altura')
    plt.grid()

    # plt.axis('equal')   
    return(puntos)

def revisa_fresnel(analiza_punto, archivo_ruta,alturaGw,alturapunto,
                  plantas,guarda=False,nFres =1,freq=915,muestras=41):
    LineaObs = archivo_ruta.strip('.csv').split('_')[1]

    planta_desde  = plantas[0]
    planta_hasta  = plantas[1]
    planta_altura = plantas[2]

    # lambda Fresnel
    lamb = 300/freq #300e6(m/s)/(freq*1e6)
    puntos = pd.read_csv(archivo_ruta,index_col='punto')

    # busca gateway
    i_Gw    = 0
    i_punto = 0
    n_Gtw = ''
    i = 0
    for unpunto in puntos.index:
        if unpunto.startswith('Gw'):
            n_Gtw = unpunto
        if unpunto==analiza_punto:
            i_punto = i
        i = i+1

    # perfil   
    xi = np.array(puntos['dist_'+n_Gtw])
    yi = np.array(puntos['altitud'])
    etiqueta = puntos.index

    # plantaciones
    conplantas = (xi > planta_desde) & (xi < planta_hasta)
    plantacion = yi[conplantas]+planta_altura

    # añade alturas de antenas
    yi_ant = np.copy(yi) + alturapunto
    yi_ant[i_Gw] = yi[i_Gw] + alturaGw

    # Zona de Fresnel
    if i_punto>0:
        
        # linea directa
        dy = yi_ant[i_punto]-yi_ant[i_Gw]
        dx = xi[i_punto]-xi[i_Gw]
        m = dy/dx
        alpha = np.arctan(m)
        yi_f = lambda x: m*x + yi_ant[i_Gw]

        xi_D = np.linspace(xi[i_Gw],xi[i_punto],
                           muestras)
        yi_D = yi_f(xi_D)
        dist_D = np.sqrt(dx**2+dy**2)

        # zona Fresnel 1
        F1_up = np.zeros(len(xi_D))
        F1_down = np.zeros(len(xi_D))
        for i in range(0,len(xi_D),1):
            d1 = xi_D[i]/np.cos(alpha)
            d2 = dist_D-d1
            F1 = np.sqrt(np.abs(nFres*lamb*d1*d2)/(d1+d2))

            xiu = xi_D[i]-F1*np.sin(np.abs(alpha))
            d1_u = xiu
            d2_u = dist_D-d1_u

            xid = xi_D[i]+F1*np.sin(np.abs(alpha))
            d1_d = xid
            d2_d = dist_D-d1_d
            # Fresnel, formula
            Fup = np.sqrt(np.abs(nFres*lamb*d1_u*d2_u)/(d1_u+d2_u))
            Fdown = np.sqrt(np.abs(nFres*lamb*d1_d*d2_d)/(d1_d+d2_d))
            
            # Linea directa +/- Fresnel
            F1_up[i]   = yi_D[i] + Fup*np.cos(alpha)
            F1_down[i] = yi_D[i] - Fdown*np.cos(alpha)
    else:
        print('punto no encontrado',analiza_punto)

    # SALIDA
    # print('Radio Fresnel: ',np.max(F1u[i]))

    # GRAFICA
    plt.fill_between(xi[conplantas],yi[conplantas],
                     plantacion,color='lightgreen',
                     label = 'plantas')
    plt.plot(xi,yi,label='perfil', color = "brown")
    # antenas
    plt.scatter(xi,yi_ant,label=
                'Dispositivos',color='blue')
    if i_punto>0:
        # Fresnel
        plt.plot(xi_D,yi_D,
                 label='Directo',color='orange',
                 linestyle = 'dashed')
        plt.plot(xi_D,F1_down,
                 label='Fresnel1',color='orange')
        plt.plot(xi_D,F1_up,
                 color='orange')
        
    for i in range(0,len(xi),1):
        plt.annotate(etiqueta[i],
                     (xi[i],yi_ant[i]),rotation=45)
    plt.xlabel('distancia')
    plt.ylabel('altura')
    plt.title('Línea :'+LineaObs+', Enlace: '+n_Gtw+'-'+etiqueta[i_punto])
    plt.legend()
    plt.grid()
    # plt.axis('equal')
    # plt.show()
    return(puntos)