Commit 9c9edd41 authored by Rafael Pezzi's avatar Rafael Pezzi
Browse files

Adicionado algoritmo para analise de intervalo de tempo dentro da

série temporal; Licença AGPL V3.
 Código com mais comentários, teste da existência do arquivo de log
parent d3a3c1a7
# -*- coding:utf-8 -*- # -*- coding:utf-8 -*-
# Programa de visualização de dados meteorológicos coletados pela estação # Programa de visualização e estatísticas de dados meteorológicos
# Meteorolog do Centro de Tecnologia Acadêmica # coletados pela estação Meteorolog do Centro de Tecnologia Acadêmica
# Dados no formato DataHora Temperatura Umidade_do_Ar Solo Pressão Luminosidade # Dados no formato DataHora Temperatura Umidade_do_Ar Solo Pressão Luminosidade
# DataHora no formato AAAAMMDDHHmmSS onde # DataHora no formato AAAAMMDDHHmmSS onde
# AAAA = Ano # AAAA = Ano
...@@ -9,47 +9,106 @@ ...@@ -9,47 +9,106 @@
# HH = Hora # HH = Hora
# mm = Minuto # mm = Minuto
# ss = Segundo # ss = Segundo
#
# Mais informações: # Mais informações:
# http://cta.if.ufrgs.br/projects/estacao-meteorologica-modular/wiki/Meteorolog # http://cta.if.ufrgs.br/projects/estacao-meteorologica-modular/wiki/Meteorolog
#
# Autor: Rafael Pezzi - Centro de Tecnologia Acadêmica - IF/UFRGS
#
# Este programa é softare livre, disponível sob os termos da
# GNU Affero General Public License V3
# http://www.gnu.org/licenses/agpl-3.0.html
import numpy as np import numpy as np
import matplotlib.pyplot as plt import matplotlib.pyplot as plt
import matplotlib.dates as mdates import matplotlib.dates as mdates
import random as random import random as random
import sys import datetime as datetime
import sys , os
# Armazena os dados coletados em vetores unidimensionais
print "Estatísticas Meteorológicas Meteorolog - Centro de Tecnologia Acadêmica"
# Testa se o número de argumentos está correto
if len(sys.argv) < 2:
sys.stderr.write('Uso: python '+ sys.argv[0]+' arquivo_de_log\n' )
sys.exit(1)
# e se o arquivo de log existe
if not os.path.exists(sys.argv[1]):
sys.stderr.write('ERRO: Arquivo '+sys.argv[1]+' não encontrado!\n')
sys.exit(1)
# Ativa/desativa seleção de intervalo
# Caso negativo (False), todos dados do arquivo de log serão computados
seleciona_intervalo = False
#seleciona_intervalo = True
# Defina aqui o intervalo de tempo a se analisado
# TODO: Incluir seleção pela linha de comando de execução
sel_Inicio = datetime.datetime(2013,03,01,00,00,00)
sel_Fim = datetime.datetime(2013,03,31,23,59,59)
# Armazena os dados coletados em matrizes unidimensionais
data, temp, u_ar, u_solo, pressao, lum = np.loadtxt(sys.argv[1], unpack=True, data, temp, u_ar, u_solo, pressao, lum = np.loadtxt(sys.argv[1], unpack=True,
converters={ 0: mdates.strpdate2num('%Y%m%d%H%M%S')}) converters={ 0: mdates.strpdate2num('%Y%m%d%H%M%S')})
## Seria util incluir um algoritmo de filtragem para remover pontos expurios ## Seria útil incluir um algoritmo de filtragem para remover pontos expurios
## observados em algumas partes por mil medidas de temperatura ## observados em algumas partes por mil medidas de temperatura
# Determina o numero de amostras # Exibe informações dos dados
amostras = len(data) print "Dados de "+str(mdates.num2date(data[0]).strftime('%Y-%m-%d %H:%M:%S'))+u" até "+str(mdates.num2date(data[len(data)-1]).strftime('%Y-%m-%d %H:%M:%S'))
print "Série com "+str(amostras)+" amostras" print "Série com "+str(len(data))+" amostras\n"
# Seleciona dados em um intervalo de tempo
if (seleciona_intervalo == True):
# Cria listas temporárias para seleção
# Acho que pode ser otimizado, usando iteradores direto em arrays
# crítico para series de dados gigantes
selec_data=[]; selec_temp=[]; selec_u_ar=[]
selec_u_solo=[]; selec_pressao=[]; selec_lum=[]
for j in enumerate(data):
if mdates.date2num(sel_Inicio) <= j[1] <= mdates.date2num(sel_Fim):
selec_data.append(data[j[0]])
selec_temp.append(temp[j[0]])
selec_u_ar.append(u_ar[j[0]])
selec_u_solo.append(u_solo[j[0]])
selec_pressao.append(pressao[j[0]])
selec_lum.append(lum[j[0]])
# Copia dados selecionados para matrizes originais
data=np.asarray(selec_data)
temp=np.asarray(selec_temp)
u_ar=np.asarray(selec_u_ar)
u_solo=np.asarray(selec_u_solo)
pressao=np.asarray(selec_pressao)
lum=np.asarray(selec_lum)
# Exibe informações sobre o periodo selecionado
print "Selecionado período:"
print str(sel_Inicio)+" até "+str(sel_Fim)
print "Série com "+str(len(data))+" amostras\n"
## Análise da temperatura ## Análise da temperatura
T_min=min(temp) T_min=temp.min()
T_max=max(temp) T_max=temp.max()
T_med=sum(temp)/amostras T_med=temp.mean()
print "Temperatura (oC)\n \t min: "+str(T_min)+"\t max: "+str(T_max)+"\t média: "+ str(T_med) T_std=temp.std()
print "Temperatura: \nmin: "+str(T_min)+"°C\tmax: "+str(T_max)+"°C\tmédia: "+str("%.2f" % T_med)+"°C\tDesvio: "+str("%.1f" % T_std)+"°C"
## Análise da umidade relativa do Ar ## Análise da umidade relativa do Ar
U_ar_min=min(u_ar) U_ar_min=u_ar.min()
U_ar_max=max(u_ar) U_ar_max=u_ar.max()
U_ar_med=sum(u_ar)/amostras U_ar_med=u_ar.mean()
print "\nUmidade Relativa do Ar (%)\n \t min: "+str(U_ar_min)+"\t max: "+str(U_ar_max)+"\t média: "+ str(U_ar_med) U_ar_std=u_ar.std()
print "\nUmidade Relativa do Ar: \nmin: "+str(U_ar_min)+"%\tmax: "+str(U_ar_max)+"%\tmédia: "+str("%.1f" % U_ar_med)+"%\tDesvio: "+str("%.1f" % U_ar_std)+"%"
# Análise da pressão atomsférica # Análise da pressão atomsférica
print "\nPressão (Pa)\n \t min: "+str(min(pressao))+"\t max: "+str(max(pressao))+"\t média: "+ str(sum(pressao)/amostras)+"\n" P_min=pressao.min()
P_max=pressao.max()
#Outro atalho para calcular as maximas e minimas P_med=pressao.mean()
#print temp.max(), u_ar.max(), pressao.max() P_std=pressao.std()
#print temp.min(), u_ar.min(), pressao.min()
#print temp.mean(), u_ar.mean(), pressao.mean()
#print temp.std(), u_ar.std(), pressao.std()
print "\nPressão: \nmin: "+str(int(P_min))+" Pa\tmax: "+str(int(P_max))+" Pa\tmédia: "+ str("%.0f" % P_med)+" Pa\tDesvio: "+str("%.2f" % P_std)+" Pa"
# Apresentação de histrograma: # Apresentação de histrograma:
# Temperatura # Temperatura
...@@ -57,8 +116,9 @@ plt.figure() # Cria nova figura ...@@ -57,8 +116,9 @@ plt.figure() # Cria nova figura
plt.ylabel(u"Número de eventos") plt.ylabel(u"Número de eventos")
plt.xlabel(u"Temperatura (°C)") plt.xlabel(u"Temperatura (°C)")
plt.hist(temp,bins=25) plt.hist(temp,bins=25)
#Inclui texto com média e desvio padrão #Inclui título com média e desvio padrão
plt.text(plt.xlim()[0]+0.5,0.9*plt.ylim()[1],u"Média: "+str("%.2f" % temp.mean())+u"°C\nDesvio Padrão: "+str("%.2f" % temp.std())+u"°C") #plt.text(plt.xlim()[0]+0.5,0.9*plt.ylim()[1],u"Média: "+str("%.2f" % temp.mean())+u"°C\nDesvio Padrão: "+str("%.2f" % temp.std())+u"°C") # Este coloca dentro da figura
plt.title(u"Temperatura média: "+str("%.2f" % T_med)+u"°C | Desvio Padrão: "+str("%.2f" % temp.std())+u"°C")
plt.grid(True) plt.grid(True)
# Apresentação de histrograma: # Apresentação de histrograma:
...@@ -67,7 +127,8 @@ plt.figure() ...@@ -67,7 +127,8 @@ plt.figure()
plt.ylabel(u"Número de eventos") plt.ylabel(u"Número de eventos")
plt.xlabel("Umidade Relativa (%)") plt.xlabel("Umidade Relativa (%)")
plt.hist(u_ar,bins=range(26,100,2)) plt.hist(u_ar,bins=range(26,100,2))
plt.text(plt.xlim()[0]+1,0.9*plt.ylim()[1],u"Média: "+str("%.2f" % u_ar.mean())+u" %\nDesvio Padrão: "+str("%.2f" % u_ar.std())+" %") #plt.text(plt.xlim()[0]+1,0.9*plt.ylim()[1],u"Média: "+str("%.2f" % u_ar.mean())+u" %\nDesvio Padrão: "+str("%.2f" % u_ar.std())+" %")
plt.title(u"UR média: "+str("%.2f" % u_ar.mean())+u" % | Desvio Padrão: "+str("%.2f" % u_ar.std())+" %")
plt.grid(True) plt.grid(True)
# Apresentação de histrograma: # Apresentação de histrograma:
...@@ -76,22 +137,20 @@ plt.figure() ...@@ -76,22 +137,20 @@ plt.figure()
plt.ylabel(u"Número de eventos") plt.ylabel(u"Número de eventos")
plt.xlabel(u"Pressão atmosférica (Pa)") plt.xlabel(u"Pressão atmosférica (Pa)")
plt.hist(pressao,bins=25) plt.hist(pressao,bins=25)
plt.text(1.001*plt.xlim()[0],0.9*plt.ylim()[1],u"Média: "+str("%.2f" % pressao.mean())+u" Pa\nDesvio Padrão: "+str("%.2f" % pressao.std())+" Pa") #plt.text(1.001*plt.xlim()[0],0.9*plt.ylim()[1],u"Média: "+str("%.2f" % pressao.mean())+u" Pa\nDesvio Padrão: "+str("%.2f" % pressao.std())+" Pa")
plt.title(u"Pressão média: "+str("%.2f" % pressao.mean())+u" Pa | Desvio Padrão: "+str("%.2f" % pressao.std())+" Pa")
plt.grid(True) plt.grid(True)
# Apresentação de histrograma: # Apresentação de histrograma:
# Luminosidade # Luminosidade
plt.figure() plt.figure()
plt.ylabel(u"Número de eventos") plt.ylabel(u"Número de eventos")
plt.xlabel(u"Luminosidade (u.a.)") plt.xlabel(u"Luminosidade (u.a.)")
plt.hist(lum,bins=25) plt.hist(lum,bins=25)
plt.text(1*plt.xlim()[0]+60,0.9*plt.ylim()[1],u"Média: "+str("%.2f" % lum.mean())+u"\nDesvio Padrão: "+str("%.2f" % lum.std())+" ") #plt.text(1*plt.xlim()[0]+60,0.9*plt.ylim()[1],u"Média: "+str("%.2f" % lum.mean())+u"\nDesvio Padrão: "+str("%.2f" % lum.std())+" ")
plt.title(u"Luminosidade média: "+str("%.2f" % lum.mean())+u" | Desvio Padrão: "+str("%.2f" % lum.std())+" ")
plt.grid(True) plt.grid(True)
## Apresentação das series temporais de medidas ## Apresentação das series temporais de medidas
# Temperatura # Temperatura
plt.figure() plt.figure()
...@@ -102,7 +161,7 @@ plt.plot_date(x=data, y=temp, fmt="r.") ...@@ -102,7 +161,7 @@ plt.plot_date(x=data, y=temp, fmt="r.")
# Define os limites temporais do gráfico # Define os limites temporais do gráfico
# útil para evitar escalas exageradas devido a dadas incorretas - ocorreu na estação sem relógio # útil para evitar escalas exageradas devido a dadas incorretas - ocorreu na estação sem relógio
plt.xlim(data[0], data[len(data)-1]) plt.xlim(data[0], data[len(data)-1])
plt.ylabel("Temperatura (C)") plt.ylabel(u"Temperatura (°C)")
plt.grid(True) plt.grid(True)
# Altera as margens para a margem inferior comportar a data # Altera as margens para a margem inferior comportar a data
plt.subplots_adjust(left=0.1, right=0.9, top=0.9, bottom=0.2) plt.subplots_adjust(left=0.1, right=0.9, top=0.9, bottom=0.2)
...@@ -136,11 +195,10 @@ plt.ylabel(u"Pressão atmosférica (Pa)") ...@@ -136,11 +195,10 @@ plt.ylabel(u"Pressão atmosférica (Pa)")
plt.grid(True) plt.grid(True)
plt.subplots_adjust(left=0.13, right=0.9, top=0.9, bottom=0.2) plt.subplots_adjust(left=0.13, right=0.9, top=0.9, bottom=0.2)
# Luminosidade # Luminosidade
plt.figure() plt.figure()
plt.xticks(rotation=70) plt.xticks(rotation=70)
plt.plot_date(x=data, y=lum, fmt="ro") plt.plot_date(x=data, y=lum, fmt="r.")
plt.xlim(data[0], data[len(data)-1]) plt.xlim(data[0], data[len(data)-1])
plt.ylabel("Luminosidade (u.a.)") plt.ylabel("Luminosidade (u.a.)")
plt.grid(True) plt.grid(True)
...@@ -149,6 +207,8 @@ plt.subplots_adjust(left=0.1, right=0.9, top=0.9, bottom=0.2) ...@@ -149,6 +207,8 @@ plt.subplots_adjust(left=0.1, right=0.9, top=0.9, bottom=0.2)
''' '''
# Análise Monte Carlo para as médias # Análise Monte Carlo para as médias
# Para fins de ilustração didática do teorema do limite central
# Cria variaveis para calculo das médias # Cria variaveis para calculo das médias
temp_med = 0. temp_med = 0.
temp_med_vec=[] temp_med_vec=[]
...@@ -159,7 +219,7 @@ pressao_med_vec=[] ...@@ -159,7 +219,7 @@ pressao_med_vec=[]
for tentativas in [5,10,50]: for tentativas in [5,10,50]:
for j in range(5000): for j in range(5000):
for n in range(tentativas): for n in range(tentativas):
amostra=int(random.random()*amostras) amostra=int(random.random()*len(data))
temp_med+=temp[amostra] temp_med+=temp[amostra]
umid_med+=u_ar[amostra] umid_med+=u_ar[amostra]
pressao_med+=pressao[amostra] pressao_med+=pressao[amostra]
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment