Mi primer (jupyter) notebook: gráfica de la evolución de la temperatura en una ciudad determinada, a lo largo del tiempo

Os comparto mi primer notebook.

No se puede llamar un análisis de datos. Es simplemente una pequeña práctica que reúne lo que llevo aprendido en el tiempo que llevo cursando el Máster en Business Analytics y Big Data UCJC de IMF Business School. Python y su ecosistema (numpy, matplotlib pandas) forma parte de los contenidos del Módulo I (Fundamentos de Tratamiento de Datos para Data Science).

El script pinta una gráfica de la evolución de la temperatura en una ciudad determinada, a lo largo del tiempo (se puede cambiar la ciudad en la variable "ciudad" al comienzo del script).

Para ello utiliza Open Data del The KNMI Data Centre (KDC) (descarga el archivo de datos la primera vez que se ejecuta).

En este zip encontraréis dos archivos, el notebook (NetCDF.ipynb) y un archivo complementario (geolocalizacion.py) con funciones que obtienen la geolocalizacion de la ciudad de Google.

Esta es la previsualización, como adelanto:

In [1]:
import netCDF4 as nCDF4
import numpy as np


Fuente de datos: The KNMI Data Centre (KDC) Archivo: BeschikbaarPBL IMAGE SSP scenarios - Annual temperature (30 arcmin)

Nota: Recuperamos los datos de la URL solo la primera vez (si no tenemos el archivo)


In [2]:
from pathlib import Path
import urllib.request

fileName = "PBL_IMAGE_SSP2_SPA2_RCP_6.0_annual_temperature_30min.nc"
my_file = Path(fileName)

if not my_file.is_file():
    print("Descargando archivo...")
    remote_path = "https://data.knmi.nl/download/SSPs_GANNTMP/V17/noversion/1970/01/01/PBL_IMAGE_SSP2_SPA2_RCP_6.0_annual_temperature_30min.nc"
    urllib.request.urlretrieve(remote_path, fileName)

print("Comenzamos a procesar el archivo.")    
dataset = nCDF4.Dataset(fileName)

print(dataset.file_format)


Comenzamos a procesar el archivo.
NETCDF4


In [3]:
## Se puede cambiar la localizacion, funciona con otras ;-)

ciudad = 'Valencia, Spain'


In [4]:
## (antes de empezar:
##  averiguamos lo que tiene con ncdump en linea de comandos -desde el prompt del anaconda-)
##
## con esto vemos los atributos
for attr in dataset.ncattrs():
    print(attr, '=:=', getattr(dataset, attr))


title =:= IMAGE 3.0, project: SSPs, scenario: SSP2_SPA2_60, variable: GANNTMP_30MIN
conventions =:= CF-1.6
institution =:= PBL Netherlands Environmental Assessment Agency (http://www.pbl.nl/en )
source =:= IMAGE-LPJ; image25/MODEL/branches/imager1757#1833,interface_imagelpj/trunk#130,svn/LPJ/branches/LPJ_IMAGE#3702;image_WIP/SSPs#9643
references =:= http://www.pbl.nl/image;https://doi.org/10.1016/j.gloenvcha.2016.05.008
history =:=  
licence =:= CC BY. The IMAGE-team would appreciate to be involved in projects using the data. Here is how to cite: http://themasites.pbl.nl/models/image/index.php/contact


In [5]:
print(dataset.dimensions.keys())


odict_keys(['latitude', 'longitude', 'time'])


In [6]:
dataset.variables['latitude'].dimensions


Out[6]:
('latitude',)


In [7]:
dataset.variables['longitude'].dimensions


Out[7]:
('longitude',)


Averiguamos las coordenadas más cercanas a -c i u d a d- entre los valores del archivo (sabiendo que las coordenadas de Madrid -por ejemplo- son lat = 40.4893538; long = -3.6827461


In [8]:
def cantidadMasCercanaSuperior(arreglo, cantidad):
    arregloRestao = arreglo - cantidad
    arregloRestaoPositivo = arregloRestao[arregloRestao >= 0]
    miCantidad = arregloRestaoPositivo.min()
    goodvalues = [miCantidad]
    filtrador = np.isin(arregloRestao[:], goodvalues)
    return arreglo[filtrador][0]

def coordenadasMasCercanas(dataset, latitud, longitud):
  latitudAprox = cantidadMasCercanaSuperior(dataset.variables['latitude'][:], latitud)
  longitudAprox = cantidadMasCercanaSuperior(dataset.variables['longitude'][:], longitud)
  return latitudAprox, longitudAprox


Las funciones obtener_geocodificacion_google y obtener_coordenadas están en otro archivo geolocalizacion.py


In [9]:
from geolocalizacion import obtener_geocodificacion_google, obtener_coordenadas
geocodificacion = obtener_geocodificacion_google(ciudad)
geocodificacion = geocodificacion.decode("utf-8")
latitud, longitud = obtener_coordenadas(geocodificacion)


Recuperando  http://maps.googleapis.com/maps/api/geocode/json?sensor=false&address=Valencia%2C+Spain
lat= 39.4699075 lng= -0.3762881
location: Valencia, Spain


Si no queremos utilizarlo, podemos simplemente inicializar las variables latitud y longitud con sus coordenadas, así latitud = 40.4893538 longitud = -3.6827461 Esto es conveniente, por ejemplo, cuando google se "cansa", da el error

'error_message': 'You have exceeded your daily request quota for this API. We recommend registering for a key at the Google Developers Console

de todas formas, "se le pasa" pronto


In [10]:
latitudAprox, longitudAprox = coordenadasMasCercanas(dataset, latitud, longitud)
print(latitudAprox)
print(longitudAprox)


39.75
-0.25


In [11]:
dataset.variables['GANNTMP_30MIN'].dimensions


Out[11]:
('time', 'latitude', 'longitude')


De los valores de las variables 'latitude' y 'longitude' vamos a crear los indices:


In [12]:
tims = np.ones(len(dataset.variables['time']), dtype=bool)
tims = tuple(tims)

goodvalues = [latitudAprox] ## 40.75 en el caso de Madrid
lats = dataset.variables['latitude'][:]
lats = np.isin(lats, goodvalues)
lats = tuple(lats)

goodvalues = [longitudAprox] ## -3.25 en el caso de Madrid
lons = dataset.variables['longitude'][:]
lons = np.isin(lons, goodvalues)
lons = tuple(lons)


In [13]:
dataset.variables['GANNTMP_30MIN'].dimensions


Out[13]:
('time', 'latitude', 'longitude')


De la variable 'GANNTMP_30MIN' es de donde vamos a sacar la data.


In [14]:
miVar = (dataset.variables['GANNTMP_30MIN'])
annualTemperature = miVar[tims, lats, lons][:,0,0]
print(annualTemperature)
## esto obtiene el resultado como tupla
## resultado = tuple(annualTemperature)
## esto obtendría el resultado como una lista
## resultado = (annualTemperature).tolist()
## esto obtendría el resultado como un ndarray
## resultado = annualTemperature


[16.069872 16.126743 16.183615 16.328356 16.473095 16.625582 16.778063
 17.028795 17.279526 17.41638  17.608303 17.818867 18.03262  18.245716
 18.453545 18.659178 18.861135 19.060385 19.259811 19.44572  19.623772
 19.798883 19.952219 20.096075 20.223919 20.34253  20.460339]


In [15]:
dataset.variables['time'].dimensions


Out[15]:
('time',)


In [16]:
print(dataset.variables['time'][:])


[    0.  1826.  3652.  5479.  7305.  9131. 10957. 12784. 14610. 16436.
 18262. 20089. 21915. 23741. 25567. 27394. 29220. 31046. 32872. 34699.
 36525. 38351. 40177. 42004. 43830. 45656. 47482.]


Convertimos los scalar a datetimes de veras


In [17]:
times = dataset.variables['time']
idxTimes = nCDF4.num2date(times[:],times.units)


In [19]:
import pandas as pd
import matplotlib.pyplot as plt
pSerie = pd.Series(annualTemperature[:],index=idxTimes)
fig = plt.figure(figsize=(12,4))
ax = fig.add_subplot(111)
pSerie.plot(ax=ax,title='%s at %s' % (dataset.variables['GANNTMP_30MIN'].long_name,ciudad))
ax.set_ylabel(dataset.variables['GANNTMP_30MIN'].units)


Out[19]:
Text(0,0.5,'oC')






Os agradezco cualquier comentario al respecto, correcciones, agregados y en cualquier caso, espero que lo disfrutéis o por lo menos os sirva de algo, que no sea una pérdida de tiempo para vosotros.

Alberto Morales Morales

Software craftsman. Passion for developing quality code that can be proud of. Happily married.

Madrid, Spain.