A continuación se analizarán los datos de las precipitaciones promedio anuales a escala regional para el norte de Perú. El conjunto de datos cuenta con 20 estaciones de la costa, los andes y la selva de Perú con las diferencias en la elevación y las tasas de precipitación.
El objetivo principal de esta es analizar los patrones de precipitación y sus tendencias utilizando Python.
Para esta tarea utilizamos IPython Notebook que usted puede instalar en su propio ordenador o ejecutarlo en la nube desde wakari.io. Además, puede descargar los datos de la estación para este ejercicio aquí.
La importación del archivo
Primero importamos los archivos en el Notebook IPython que se ejecuta en su navegador favorito.
In[]:%pylab inline
%cd C:\Users\Saul\Dropbox\Curso_19_Python_en_Hidrologia_v2\1_Doc\AuxData
import numpy as np
import matplotlib.pyplot as plt
datos = np.genfromtxt('DatosEstaciones.csv', delimiter=",")
datos
Out[]:
array([[ nan,nan,nan,nan,nan],
[-76.25, -10.68,4334.,nan,1052.04],
[-75.75,-9.9 ,1860.,nan, 378.],
[-76.78,-9.55,3430.,nan, 702.],
[-77.36, -10.06,3950.,nan, 755.04],
[-75.53, -10.03, 680.,nan,2424.96],
[-75.9 , -11.48,3632.,nan, 735.96],
[-77.43,-9.73,3561.,nan, 744.96],
[-76.21, -11.63,4193.,nan, 741.],
[-75.93, -11.55,3750.,nan, 567.],
[-75.13, -10.61,1050.,nan,3096.96],
[-75.3 , -11.13, 800.,nan,2007.96],
[-77.83, -10.68,46.,nan, 9.96],
[-76.38, -11.83,2379.,nan, 237.96],
[-75.95,-9.13, 669.,nan,3203.04],
[-77.51,-9.5 ,3050.,nan, 579.],
[-74.9 , -10.3 , 301.,nan,3309.],
[-77.6 ,-9.35,2760.,nan, 737.04],
[-78.05, -10.28, 100.,nan, 6.],
[-74.93,-9.86, 250.,nan,1701.96],
[-75.46, -11.78,3322.,nan, 714.]])La primera fila son los nombres de columna que no podían ser importados con un arreglo numpy. Nombres de columnas son: longitud, latitud, elevación, nombre y precipitación anual.
No hay necesidad de fijar la matriz para eliminar la variable no-numérica indeseada "nan", ya que vamos a hacer arreglo para cortar el script.
Gráfica de los valores por separado
Se hará un par de gráficos para ver la distribución espacial de la precipitación y la elevación de las estaciones.
Distribución de las precipitaciones
En la primera gráfica se muestran las estaciones como círculos escalados con su precipitación anual promedio. El color es uniforme para todas las estaciones. El script también guarda el gráfico resultante como un archivo png.
In[]:plt.scatter(datos[:,0],datos[:,1], marker="o", s=datos[:,4])
plt.savefig('puntosconmagnitud.png')
Out[]:
La precipitación aumenta hacia el este de Perú donde se encuentran los grandes círculos. En los andes las estaciones tienen círculos relativamente similares, mientras que en la costa de Perú la precipitación es mínima (dos pequeños puntos al lado derecho).
Distribución de Elevación
Hacemos un gráfico similar, pero esta vez con elevación de la estación trazó como hexágonos con su color y tamaño a escala de la elevación.
Podemos ver los hexágonos más grandes y rojos que corresponden a las estaciones de mayor altura en los Andes. En ambos lados de los Andes (la costa del desierto y la selva tropical) elevaciones son más bajos.
Gráficos 3D para mostrar las precipitaciones y la elevación
Se ha analizado la elevación y precipitación por separado, ahora se necesitan algunas herramientas como gráficos o diagramas para representarlos en el mismo gráfico. En esta parte vamos a navegar entre algunas opciones gráficas para representar la precipitación y elevación en un mismo gráfico.
Elevación y precipitación en forma de puntos
La siguiente gráfica muestra la altura como puntos rojos y las precipitaciones como puntos azules. En el fondo de la caja podemos ver una ubicación de la estación como puntos grises.
In[]:from mpl_toolkits.mplot3d import Axes3D import numpy as np import matplotlib.pyplot as plt fig = plt.figure() fig.set_size_inches(10,10) ax = fig.gca(projection='3d') este = datos[:,0] norte = datos[:,1] elev = datos[:,2] precip = datos[:,4] ax.scatter(este, norte, elev, marker='o', c='r',s=(datos[:,4]/20)) ax.scatter(este, norte, precip, marker='o', c='b',s=(datos[:,4]/20)) ax.scatter(este, norte, 0, marker='o', c='#4b4b4b',s=(datos[:,4]/20)) ax.set_zlim3d(0, 5000) plt.show() Out[]:
La representación de la elevación y la precipitación es posible en este gráfico 3D. En un enfoque "orientado al usuario" necesitamos una gráfica más "amigable" y más "habladora".
Se hace un gráfico similar, pero esta vez vamos a utilizar columnas en lugar de puntos. Tendremos dos columnas al lado del otro, una roja para la elevación y otra azul para la precipitación.
In[]:from mpl_toolkits.mplot3d import Axes3D
import numpy as np
import matplotlib.pyplot as plt
fig = plt.figure()
fig.set_size_inches(10,10)
ax = fig.gca(projection='3d')
este = datos[:,0]
norte = datos[:,1]
datum = np.zeros(este.size)
elev = datos[:,2]
precip = datos[:,4]
ax.bar3d(datos[:,0]-0.1,datos[:,1]-0.1,datum, .1, .1, datos[:,2], color='r', alpha=0.5)
ax.bar3d(datos[:,0],datos[:,1],datum, .1, .1, datos[:,4], color='b',alpha=0.5)
ax.set_xlabel('Este')
ax.set_ylabel('Norte')
ax.set_zlabel('Elevacion / Precipitacion')
ax.set_zlim3d(0, 5000)
plt.show()Out[]:
Esta gráfica es mejor para representar los valores de precipitación y elevación.
Interpolación y trazado de cuadrícula 3D
En lugar de utilizar valores discretos, vamos a interpolar la precipitación para tener una precipitación distribuida. Scipy tiene el módulo interpolate.griddata que es adecuado para este ejercicio. Al final de la secuencia de comandos se muestra la malla de interpolación.
In[]: from scipy import interpolate %pylab inline #generamos parametros para la grilla delta=0.25 este_r = np.arange(-78.5,-74.5,delta) norte_r = np.arange(-12.0,-9.0,delta) este_g,norte_r = np.meshgrid(este_r,norte_r) #generamos la grilla regular X,Y =np.meshgrid(np.linspace(-78.5,-74.5,100), np.linspace(-12,-8,100)) #realizamos la interpolación grid_datos = interpolate.griddata((datos[1:,0],datos[1:,1]), datos[1:,4],(X,Y) ,method='cubic') clf() contourf(X,Y,grid_datos,8,cmap='Blues') colorbar() scatter(datos[:,0],datos[:,1]) xlim(np.min(datos[1:,0]),np.max(datos[1:,0])) ylim(np.min(datos[1:,1]),np.max(datos[1:,1])) Out[]:
Es posible mostrar la precipitación interpolada como contornos con etiquetas.
In[]:clabel(contour(X,Y,grid_datos,8,cmap='Blues'),fmt='%r') scatter(datos[:,0],datos[:,1]) xlim(np.min(datos[1:,0]),np.max(datos[1:,0])) ylim(np.min(datos[1:,1]),np.max(datos[1:,1])) Out[]:
Finalmente se puede trazar una malla 3D de interpolación y tiene proyección cada eje.
In[]:from mpl_toolkits.mplot3d import Axes3D %pylab inline fig = figure() ax = Axes3D(fig) cset = ax.contourf(X, Y, grid_datos, zdir='z', offset=-500, cmap=cm.coolwarm) cset = ax.contourf(X, Y, grid_datos, zdir='x', offset=-78.5, cmap=cm.coolwarm) ax.plot_surface(X,Y,grid_datos,rstride=8, cstride=8,alpha=0.3, linewidth=0.1, antialiased=True) Out[]:
