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)
