Introducción a la programación con Python

Jupyter Notebook 2 para Temas Selectos de Modelación Numérica 2021-2

Facultad de Ciencias, UNAM

El material de este "notebook" está basado en la lección Programming with Python de los Software Carpentry Workshops, publicada bajo la licencia Creative Commons - Attribution License. Traducción y modificaciones por Karina Ramos Musalem.

Antes de empezar

Al terminar este notebook (intro02_python):

  • Reforzarán los conocimientos adquiridos en el notebook pasado.
  • Aprenderán los conceptos de programación de ciclos, condicionales y funciones.
  • Aprenderán algunas funciones de matplotlib para visualización de datos 2D.
  • Aprenderán qué es un archivo netCDF y cómo manipularlo.

4. Repetir acciones con ciclos o loops

Preguntas a responder

  1. ¿Cómo repetir la misma operación en muchos valores distintos?

Objetivos

  • Explicar qué es un ciclo for (for loop).
  • Escribir correctamente un ciclo for para repetir cálculos simples.
  • Monitorear cambios en la variable del ciclo conforme corre el ciclo.
  • Monitorear cambios en las otras variables conforme se actualizan durante el ciclo for.

En la clase pasada escribimos código de Python que grafica valores de interés de nuestro archivo de datos (precip_mensual_estacion01.txt). Tenemos otros 6 archivos de datos de otras variables meteorológicas y podríamos tener muchos más. Queremos crear gráficas para todos nuestros archivos y para eso tenemos que enseñarle a la computadora a repetir acciones.

Por ejemplo, podemos repetir la impresión de cada caracter en un palabra:

In [1]:
palabra = 'plomo'

En Python, un string o cadena es básicamente una colección ordenada de caracteres y cada caracter tiene un número asociado - su índice (similar al índice que estudiamos en el notebook 1 para arreglos). Así que podemos acceder a los caracteres de una cadena usando sus índices. Por ejemplo, podemos obtener el primer caracter de la palabra 'plomo'usando palabra[0]. Una forma de imprimir cada caracter es usando 5 comandos print:

In [2]:
print(palabra[0])
print(palabra[1])
print(palabra[2])
print(palabra[3])
print(palabra[4])
p
l
o
m
o

¡Guácala! Esta forma de hacerlo no es buena por tres razones:

  1. No es escalable. Imagínate que tienes que imprimir algo con cientos de caracteres. Es fácil hacerlo con cinco, pero hacelo con cien es muy impráctico.
  2. Es dificil de mantener y modificar. Si queremos decorar cada caracter con un asterisco o algún otro caracter, tendríamos que cambiar línea por línea del código. De nuevo, parece fácil cuado hay 5 líneas pero no lo es tanto para cientos de líneas.
  3. Es frágil. Si usamos este código con otra palabra que tenga un número distinto de caracteres, el código no va a servir. Si la palabra es más larga, solo aparecerán parte de las letras, y si es más corta se producirá un error porque no habrá nada que imprimir.

Por ejemplo:

In [3]:
palabra = 'uno'
print(palabra[0])
print(palabra[1])
print(palabra[2])
print(palabra[3])
print(palabra[4])
u
n
o
---------------------------------------------------------------------------
IndexError                                Traceback (most recent call last)
<ipython-input-3-c8a801f7cd23> in <module>
      3 print(palabra[1])
      4 print(palabra[2])
----> 5 print(palabra[3])
      6 print(palabra[4])

IndexError: string index out of range

Pero hay una mejor solución:

In [4]:
palabra = 'plomo'

for char in palabra:
    print(char)
p
l
o
m
o

Esto es más corto - al menos más corto que 100 prints - y es más robusto (podemos usar el mismo código con cualquier palabra):

In [5]:
palabra = 'oxygen'

for char in palabra:
    print(char)
o
x
y
g
e
n

Esta versión usa un loop o ciclo for para repetir una operación - en este caso imprimir un caracter - una vez por cada 'elemento' en un secuencia. La forma general de un ciclo for es:

for variable in coleccion:

# haz algo usando variable, por ejemplo imprimir

Usando el ejemplo anterior (solo cambia word por palabra):

En donde cada caracter (char) en la variable word o en nuestro caso palabra es pasada por el ciclo e impresa, un caracter tras otro. Los números en el diagrama denotan en qué ciclo se imprimió el caracter (1 es el primer ciclo y 6 es es último).

Podemos llamar a la variable del ciclo como queramos, pero debe haber un símbolo :al final de la línea inicial de ciclo y siempre debemos respetar la alineación: Lo que queremos correr dentro del loop debe estar "indentado" hacia la derecha. En Python no hay un símbolo para indicar el final del ciclo; en cambio, lo que está indentado despues del for será parte del ciclo.

Aquí hay otro ciclo que actualiza el valor de una variable cada vuelta:

In [6]:
largo = 0
for vocal in 'aeiou':
    largo = largo + 1
print('Hay', largo, 'vocales')
Hay 5 vocales

Vale la pena desmenuzar la ejecución de este programa paso por paso. Como hay cinco caracteres en 'aeiou', el código en la tercera línea se ejecutará cinco veces. La primera vez, largoes cero (el valor que le asinamos en la primera línea del código) y vocal es 'a'. La línea 3 agrega 1 al valor de largo, produciendo 1, y actualiza largopara referir al nuevo valor. En la segunda vuelta, vocal es 'e' y largoes 1, así que largo se actualiza al valor 2. Después de las tres vueltas siguentes, largoes 5; como no hay nungún otro elemento en 'aeiou'para que Python procese, el ciclo se termina y el comando print de la cuarta línea se ejecuta y nos dice la respuesta final.

Nota que la variable del ciclo se está usando para llevar un registro del progreso del ciclo. Sigue existiendo cuando termina el ciclo y también podemos reutilizar varables que se usaron como variables de ciclo anterirormente:

In [7]:
letra = 'z'
for letra in 'abc':
    print(letra)
print('Después del ciclo, letra es', letra)
a
b
c
Después del ciclo, letra es c

Nota también que encontrar la longitud de un cadena de caracteres es una operación tan común que Python tiene una función integrada llamada lenpara hacerlo:

In [8]:
print(len('aeiou'))
5

len es mucho más rápida que cualquier función que podríamos escribir y mucho más fácil de leer que el ciclo que escribimos; nos puede dar la longitud de muchas otras cosas además de cadenas, así que la seguiremos usando cuando podamos.

Pon en práctica lo que aprendiste

  1. Python tiene una función llamada rangeque genera una secuencia de números. range acepta 1, 2, o 3 parámetros.

    • Si le damos un parámetro, range genera una secuencia de esa longitud que comienza en 0 e incrementa en 1. Por ejemplo range(3) produce números 0,1,2.
    • Si le damos 2 parámetros, rangecomienza en el primero y va hasta el segundo en incrementos de 1. Por ejemplo range(2,5) produce 2, 3, 4.
    • Si le damos 3 parámetros, empieza en el primero, va hasta el segundo e incrementa en pasos dados por el tercero. Por ejemplo range(3, 10, 2) produce 3, 5, 7, 9.

    Usando range, escribe un ciclo que imprima los primeros 10 números naturles (1,2,3,...,10)

In [9]:
# Respuesta

for ii in range(1,11):
    print(ii)
1
2
3
4
5
6
7
8
9
10
  1. Podemos exponenciar algo usando la notación **. Por ejemplo $2^3$ es 2**3. Escribe un ciclo que calcule el mismo resultado que 2 ** 3 utilizando únicamente la multiplicación * (sin usar exponenciación).
In [16]:
# Respuesta
a = 1
for ii in range(3):
    a = a*2

print(a)

#Comprobación
print(2**3)
8
8

4.1 En resumen

  • Usa for variable in secuencia: para procesar uno por uno los elementos de una secuencia.
  • El cuerpo de un ciclo for debe estar indentado.
  • Usa len(cosa)para determinar la longitud de alguna cosa que contiene valores.

5. Listas

Preguntas a responder

  • ¿Cómo puedo guardar muchos valores juntos?

Objetivos

  • Explicar qué es una lista.
  • Crear e indizar listas de valores simples.
  • Cambiar el valor de elementos individuales.
  • Reordenar y rebanar elementos de una lista.

De manera similar a un cadena que contiene muchos caracteres, una lista en un contenedor que puede guardar muchos valores. A diferencia de los arreglos de NumPy, las listas son intrínsecas de Python, por lo cual no necesitamos importar ninguna bibliotca para usarlas. Podemos crear listas poniendo valores entre [] y separandolos con ,:

In [17]:
impares = [1, 3, 5, 7]
print('los impares son:', impares)
los impares son: [1, 3, 5, 7]

Podemos acceder a los elementos de un alista usando índices:

In [18]:
print('El primer elemento es ', impares[0])
print('El último elemento es ', impares[3])
print('El elemento -1 es ', impares[-1])
El primer elemento es  1
El último elemento es  7
El elemento -1 es  7

Sí, podemos usar índices negativos en Pyhton. El elemento -1 nos da el último de la lista, el -2el penúltimo, etc.

Si hacemos un ciclo sobre un a lista, la variable del ciclo se asigna los elementos de la lista uno por uno:

In [19]:
for numero in impares:
    print(numero)
1
3
5
7

Hay una diferencia importante entre las cadenas de caracteres y las listas: podemos cambiar los valores de una lista, pero no podemos cambiar caracteres individuales en una cadena. Por ejemplo:

In [20]:
nombres = ['Curie', 'Franklyn', 'Turing'] # Franklin mal escrito
print('nombres era originalmente:', nombres)
nombres[1] = 'Franklin' # Nombre correcto
print('nombres es ahora:', nombres)
nombres era originalmente: ['Curie', 'Franklyn', 'Turing']
nombres es ahora: ['Curie', 'Franklin', 'Turing']

pero

In [21]:
nombre = 'Franklin'
nombre[0] = 'f'
---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
<ipython-input-21-36e713ca9a8b> in <module>
      1 nombre = 'Franklin'
----> 2 nombre[0] = 'f'

TypeError: 'str' object does not support item assignment

no funciona.

Los números y las cadenas son inmutables, es decir, no pueden modificarse. Las listas y los arreglos son mutables y pueden modificarse después de ser creadas.

Las listas en Python pueden contener distintos tipos de elementos. Por ejemplo:

In [22]:
edades = [10, 20.4, 'desconocido'] # enteros, decimales y cadenas

Hay otras maneras de cambiar el contenido de una lista además de reasignando valores. Podemos agregar (append):

In [23]:
impares.append(11)
print('impares después de agregar un valor:', impares)
impares después de agregar un valor: [1, 3, 5, 7, 11]
In [24]:
elemento_eliminado = impares.pop(0)
print('impares después de quitar el primer elemento:', impares)
print('elemento eliminado:', elemento_eliminado)
impares después de quitar el primer elemento: [3, 5, 7, 11]
elemento eliminado: 1
In [25]:
impares.reverse()
print('impares al revés:', impares)
impares al revés: [11, 7, 5, 3]

Hay que tener cuidado cuando modificamos una lista. Si queremos hacer una copia de la lista y modificar la copia, es posible que queramos hacer lo siguiente:

In [26]:
impares = [1, 3, 5, 7]
primos = impares # Este es el paso mortal
primos.append(2)
print('lista primos', primos)
print('lista impares', impares)
lista primos [1, 3, 5, 7, 2]
lista impares [1, 3, 5, 7, 2]

¿Qué paso aquí? Solo queríamos modificar la lista primos pero también cambió la lista impares. Esto es porque Python guarda una lista en memoria y luego puede usar distintos nombres para referirse a la misma lista. Si solo queremos copiar una lista simple, podemos usar la función list de manera que no modifiquemos una lista que no queríamos modificar:

In [27]:
impares = [1, 3, 5, 7]
primos = list(impares) # Este es el paso mortal
primos.append(2)
print('lista primos', primos)
print('lista impares', impares)
lista primos [1, 3, 5, 7, 2]
lista impares [1, 3, 5, 7]

Se puede acceder a rebanadas o subconjuntos de listas y cadenas si especificamos los intervalos de valores entre [] como hicimos con los arreglos:

In [28]:
binomial_nombre = 'Drosophila melanogaster'
grupo = binomial_nombre[0:10]
print('grupo:', grupo)

especie = binomial_nombre[11:23]
print('especie:', especie)

chromosomas = ['X', 'Y', '2', '3', '4']
autosomas = chromosomas[2:5]
print('autosoma:', autosomas)

ultimo = chromosomas[-1]
print('último:', ultimo)
grupo: Drosophila
especie: melanogaster
autosoma: ['2', '3', '4']
último: 4

Pon en práctica lo que aprendiste

  1. Usa las rebanadas para seleccionar solo los últimos 4 caracteres de una cadena y de una lista:
In [29]:
cadena_para_rebanar = 'Observaciones del 27-abr-2013'
lista_para_rebanar = ['fluor', 'cloro', 'bromo', 
                      'neodimio', 'disprosio', 'einstenio']

¿Tu solución funciona sin importar si conoces de antemano o no la longitud de la cadena y la lista? Si la respuesta es no, intenta cambiar tu estrategia para que sea más robusta.

HINT: Recuerda que los índices pueden ser negativos también.

In [30]:
# Respuesta

print(cadena_para_rebanar[-4:])
print(lista_para_rebanar[-4:])
2013
['bromo', 'neodimio', 'disprosio', 'einstenio']

5.1 En Resumen

  • [valor1, valor2, valor3, ...] crea una lista.
  • Las listas pueden contener cualquier objeto de Python, incluyendo otras listas.
  • Las listas están indizadas y se rebanan con [] como los arreglos y las cadenas.
  • Las listas son mutables (los valores de sus elementos puede cambiarse en el mismo lugar).
  • Las cadenas son inmutables (no les puedes cambiar un caracter).

6. Analizar datos de muchos archivos

Preguntas a responder

  • ¿Cómo puedo hacer las mismas operaciones en varios archivos distintos?

Objetivos

  • Usar una función de una biblioteca para obtener una lista de nombres de archivos que concuerdan con un patrón.
  • Escribir un ciclo for para procesar múltiples archivos.

Ya tenemos casi todo para procesar todos los archivos de datos que tenemos. Solo falta importar una biblioteca especial:

In [31]:
import glob

La biblioteca glob contiene una función, también llamada glob, que encuentra archivos y directorios (carpetas) cuyos nombres concuerdan con un patrón dado. Le proporcionamos esos partrones como cadenas: el caracter * representa cero o más caracteres, mientras que ?representa solo un caracter. Podemos usarlos para obtener los nombres de todos los archivos CSV en la carpeta, o directorio, en la que se encuentra este notebook:

In [32]:
print(glob.glob('*_mensual_estacion01.csv'))
['temp_mensual_estacion01.csv', 'indiceUV_mensual_estacion01.csv', 'humedad_mensual_estacion01.csv', 'presion_mensual_estacion01.csv', 'viento_mensual_estacion01.csv', 'radiacion_mensual_estacion01.csv', 'precip_mensual_estacion01.csv']

Como muestran estos ejemplos, el resultado de glob.glob es una lista de archivos y directorios en orden arbitrario. Así que podemos hacer un ciclo sobre esta lista para hacer algo con cada nombre de archivo. En este caso, ese algo es generar gráficas para todos los datos. Si queremos empezar por analizar los primeros cinco archivos por orden alfabético, podemos usar la función sorted(ordenado) para generar una lista ordenada alfabéticamente del output de glob.glob:

In [33]:
import glob
import numpy
import matplotlib.pyplot

filenames = sorted(glob.glob('*_mensual_estacion01.csv'))
filenames = filenames[0:5]
for filename in filenames:
    print(filename)

    data = numpy.loadtxt(fname=filename, delimiter=',')

    fig = matplotlib.pyplot.figure(figsize=(10.0, 3.0))

    axes1 = fig.add_subplot(1, 3, 1)
    axes2 = fig.add_subplot(1, 3, 2)
    axes3 = fig.add_subplot(1, 3, 3)

    axes1.set_ylabel('promedio')
    axes1.plot(numpy.mean(data, axis=0))

    axes2.set_ylabel('max')
    axes2.plot(numpy.max(data, axis=0))

    axes3.set_ylabel('min')
    axes3.plot(numpy.min(data, axis=0))

    fig.tight_layout()
    matplotlib.pyplot.show()
humedad_mensual_estacion01.csv
indiceUV_mensual_estacion01.csv
precip_mensual_estacion01.csv
presion_mensual_estacion01.csv
radiacion_mensual_estacion01.csv

Si nos fjamos en los mínimos, vemos que hay muchos ceros en todas las variables. Es un poco sospechoso que la radiación solar mínima en el promeido mensual sea cero, o peor aún, ¡que la humedad mínima sea cero! Hay algo raro con estos datos que debemos averiguar...

6.1 En resumen

  • Usa glob.glob(patrón) para crear una lista de nombres de archivos cuyos nombres concuerdan con el patrón.
  • Usa * en un patrón como comodín para representar cero o más caracteres, y ? como comodín para representar un solo caracter.

7. Condicionales

Preguntas a responder:

  • ¿Cómo puedo hacer que mis programas hagan distintas cosas con base en los valores de datos?

Objetivos:

  • Utilizar los condicionales if, elif y else.
  • Evaluar correctamente expresiones que contienen and y or.

Anteriormente descubrimos que había algo raro con los datos. ¿Cómo podemos usar Python para descubrir automáticamente los datos sospechosos y operar de distinta manera con ellos? En esta parte aprenderemos cómo escribir código que corra sólo cuando ciertas condiciones son ciertas.

Podemos pedirle a Python que ejecute ciertas acciones, dependiendo de alguna condición, utilizando el comando if (si). Si (if) se cumple la condición, entonces haz ... :

In [34]:
num = 37
if num > 100:
    print('mayor que')
else:
    print('menor que')
print('fin')
menor que
fin

La segunda línea usa la palabra if para decirle a Python que vamos a hacer una elección. Si la condición que sigue despuésd de if es verdadera, los comandos después del if (indentados dentro del if) se ejecutan y se imprime "mayor que". Si la prueba es falsa, se ejeculan los comandos después de elsey entonces se imprime "menor que". Sólo se ejecuta una u otra opción antes de continuar con el programa y que se imprima "fin":

No es absolutamente necesario tener un else. Si no hay else, Python simplemente no hace nada si la prueba es falsa:

In [35]:
num = 53
print('Antes del condicional...')
if num > 100:
    print(num, ' es más grande que 100')
print('... después del condicional')
Antes del condicional...
... después del condicional

También podemos encadenar varias pruebas usando elif, que es una abreviación de else if (si no). El siguiente programa usa elif para imprimir el signo de un número:

In [36]:
num = -3

if num >0:
    print(num, 'es positivo')
elif num == 0:
    print(num, 'es cero')
else:
    print(num, 'es negativo')
-3 es negativo

Nota que para probar igualdad usamos == porque =es solamente para asignar valores.

También podemos combinar pruebas usando and (y) y or (o). La prueba and es verdadera solo si ambas partes son verdaderas:

In [37]:
if (1 > 0) and (-1 > 0):
    print('ambas partes son verdaderas')
else:
    print('al menos una parte es falsa')
al menos una parte es falsa

Mientras que or es verdadero si al menos una parte es verdadera:

In [38]:
if (1 < 0) or (-1 < 0):
    print('al menos una parte es verdadera')
al menos una parte es verdadera

Verdadero y Falso

Las palabras True (verdadero) y False (falso) son especiales en Python. Se llaman booleanos y representan la "verdad" en los valores. Una expresión como 1<0 regresa el valor False, mientras que -1<0regresa el valor True.

Revisemos nuestros datos

Anteriormente vimos que el mínimo mensual de los datos para todos los años era cero. Esto es sospechoso porque ninguna variable de las que estamos analizando debería ser cero, salvo la precipitación (temperatura, presión, radiación, índice UV y humedad). ¿Puedes pensar por qué no son cero? ¿Qué crees que esté pasando aquí?

Probablemente hay meses en donde no hubo mediciones y se reportaron como ceros. Cuando no hay mediciones o no son confiables, regularmente se elige una bandera o etiqueta especial para indicar estos valores faltantes o missing values. Se pueden poner valores fuera del rango esperado de mediciones como 999999, o NaN (Not a Number). NaN es una manera de decirle a la computadora que no hay número asignado.

Revisemos en qué años hay valores cero en la radiación:

In [39]:
# Abrimos el archivo de mediciones de radiación solar
radiacion = numpy.loadtxt(fname='radiacion_mensual_estacion01.csv', delimiter=',')

# Ciclo sobre los 8 años de datos (filas)
for fila in range(8):
    # Prueba si el valor mínimo en la fila #fila es cero
    if numpy.min(radiacion[fila,:]) == 0:
        # Si la pueba es verdadera imprime esto:
        print('Hay ceros en el año ', fila)
    else:
        # Si la prueba es falsa imprime esto:
        print('No hay ceros en el año ', fila)
Hay ceros en el año  0
No hay ceros en el año  1
Hay ceros en el año  2
No hay ceros en el año  3
Hay ceros en el año  4
No hay ceros en el año  5
No hay ceros en el año  6
No hay ceros en el año  7

Por lo tanto, podemos descartar los datos para los años 0 (2011), 2 (2013), 4 (2015).

Pon en práctica lo que aprendiste

  1. Escribe un ciclo que cuente el número de vocales en una cadena.
  2. Comprueba que sirve probando distintas palabras y oraciones.
In [40]:
cadena = "Ciencias De La Tierra"
vocales = ["a", "e", "i", "o", "u", "á", "é", "í", "ó", "ú"]
contador = 0

for char in cadena:
    if char in vocales:
        contador = contador + 1

print(contador)
9

7.1 En resumen

  • Usa if condición: para empezar un condicional, elif condición: para agregar pruebas adicionales y else: para agregar una acción default.
  • Las acciones a realizar dentro de un condicional deben estar indentados.
  • Usa == para probar que dos cosas son iguales.
  • X and Y es verdadero solo si Xy Y son ambos verdaderos.
  • X or Y es verdadero si alguno o ambos X o Y son verdaderos.
  • El cero, una lista vacía, o una cadena vacía se consideran falsos; todos los demás números, cadenas y listas son verdaderos.
  • True y False representan valores de verdad verdadero y falso.

8. Funciones

Preguntas a resolver:

  • ¿Cómo puedo definir funciones nuevas?
  • ¿Cuál es la diferencia entre definir una función y llamar a un a función?
  • ¿Qué pasa cuando llamo a una función?

Objetivos:

  • Definir una función que toma parámetros.
  • Regresar un valor de una función.
  • Probar y depurar (debug) una función.
  • Definir valores por omisión para una función.
  • Explicar por qué debemos dividir los programas en unidades pequeñas con propósitos únicos (funciones).

Hasta ahora, hemos escrito código para graficar algunas estadísticas de nuestros datos, hemos aprendido a abrir y analizar todos nuenstros archivos en un solo ciclo, a hacer que Python haga cosas con los datos dependiendo de alguna condición, pero el código se está volviendo largo y complicado. Imagina que tuviéramos cientos de archivos y que no quisieramos generar una gráfica de cada uno. Comentar las líneas de código que grafican los datos es tedioso y poco robusto. Cortar y pegar tampoco es una buena solución... Lo que necesitamos es una forma de empaquetar nuestro código para que sea fácil de reutilizar y esa solución es hacer "funciones". Esta es una manera rápida de ejecutar una y otra vez pedazos más largos de código.

Empecemos definiendo una función fahr_a_celsius que convierte temperaturas en escala Fahrenheit a Celsius:

In [41]:
def fahr_a_celsius(temp):
    return((temp - 32) * (5/9))

La definición comienza con def seguido del nombre que queremos darle a la función (fahr_a_celsius) y entre paréntesis una lista de nombres de parámetros (temp). El cuerpo de la función - los comandos que se ejecutan cuando corre la función - deben estar indentados. El cuerpo concluye con return seguido del valor que regresará la función.

Cuando llamamos a la función, los valores que le pasamos de le asignan a las variables que se usarán dentro de la función. Dentro de la función, usamos la expresión return regresar el resultado a quien lo pidió.

Corramos nuestra función:

In [42]:
fahr_a_celsius(32)
Out[42]:
0.0

Este comando llama a la función, usando "32" como input y regresa el valor de la función.

Llamar a nuestra función es lo mismo que llamar a cualquier función que hemos usado antes:

In [43]:
print('El papel se quema a:', fahr_a_celsius(451), 'C')
print('El punto de ebullición del agua es:', fahr_a_celsius(212), 'C')
El papel se quema a: 232.7777777777778 C
El punto de ebullición del agua es: 100.0 C

8.1 Composición de funciones

Ahora escibamos una función que transforme de Celsius a Kelvin:

In [44]:
def celsius_a_kelvin(temp_c):
    return temp_c + 273.15

print('El punto de fusión del agua en Kelvin:', celsius_a_kelvin(0.))
El punto de fusión del agua en Kelvin: 273.15

¿Qué hay de convertir de Fahrenheit a Kelvin? Podríamos escribir otra función pero no es necesario. Podemos componer dos funciones (como en su clase de matemáticas) que ya creamos:

In [45]:
def fahr_a_kelvin(temp_f):
    temp_c = fahr_a_celsius(temp_f)
    temp_k = celsius_a_kelvin(temp_c)
    return temp_k

print('El punto de ebullición  del agua en Kelvin:', fahr_a_kelvin(212.0))
El punto de ebullición  del agua en Kelvin: 373.15

Esta es una probadita de cómo se contruyen programas más grandes: definimos funciones básicas y las combinamos en pedazos más grandes. En la vida real, las funciones son un poco mas largas que las que definimos aquí - unas pocas decenas de líneas - pero no deben ser mucho más largas que eso.

8.2 Una limpiada...

Ahora que sabemos como empacar pedazos de código, podemos hacer que nuestro análisis de los datos PEMBU sea más fácil de leer y de reusar. Hagámos una función visualiza que genere nuestras gráficas:

In [46]:
def visualiza(filename):

    data = numpy.loadtxt(fname=filename, delimiter=',')

    fig = matplotlib.pyplot.figure(figsize=(10.0, 3.0))

    axes1 = fig.add_subplot(1, 3, 1)
    axes2 = fig.add_subplot(1, 3, 2)
    axes3 = fig.add_subplot(1, 3, 3)

    axes1.set_ylabel('average')
    axes1.plot(numpy.mean(data, axis=0))

    axes2.set_ylabel('max')
    axes2.plot(numpy.max(data, axis=0))

    axes3.set_ylabel('min')
    axes3.plot(numpy.min(data, axis=0))

    fig.tight_layout()
    matplotlib.pyplot.show()

y otra función detecta_problemas que revisa si tenemos esos ceros extraños cuando faltan mediciones en los archivos:

In [47]:
def detecta_problemas(filename):
    
    data = numpy.loadtxt(fname=filename, delimiter=',')
    dimensiones = data.shape
    for fila in range(dimensiones[0]):
        if numpy.min(data[fila,:]) == 0:
            print('Hay ceros en el año ', fila, 'archivo', filename)
        

Pero, ¿no olvidamos el return? Ah, en Python no es neceario incluir un return. La funciones pueden usarse con el simple propósito de agrupar código con una función específica.

Ahora, en vez de juntar el código en un ciclo gigante, podemos leer y reusar ambas ideas separadas. Podemos reproducir lo que hemos hecho en este notebook en un código mucho más simple:

In [48]:
archivos = sorted(glob.glob('*_mensual_estacion01.csv'))

for archivo in archivos[:3]:
    print(archivo)
    visualiza(archivo)
    detecta_problemas(archivo)
humedad_mensual_estacion01.csv
Hay ceros en el año  0 archivo humedad_mensual_estacion01.csv
Hay ceros en el año  2 archivo humedad_mensual_estacion01.csv
Hay ceros en el año  4 archivo humedad_mensual_estacion01.csv
indiceUV_mensual_estacion01.csv
Hay ceros en el año  0 archivo indiceUV_mensual_estacion01.csv
Hay ceros en el año  2 archivo indiceUV_mensual_estacion01.csv
Hay ceros en el año  4 archivo indiceUV_mensual_estacion01.csv
precip_mensual_estacion01.csv
Hay ceros en el año  0 archivo precip_mensual_estacion01.csv
Hay ceros en el año  1 archivo precip_mensual_estacion01.csv
Hay ceros en el año  2 archivo precip_mensual_estacion01.csv
Hay ceros en el año  4 archivo precip_mensual_estacion01.csv
Hay ceros en el año  5 archivo precip_mensual_estacion01.csv
Hay ceros en el año  6 archivo precip_mensual_estacion01.csv
Hay ceros en el año  7 archivo precip_mensual_estacion01.csv

Definición de parámetros por omisión (default)

Hasta ahora hemos pasado parámetros de dos formas: directamente como en type(data), y por nombre, como en numpy.loadtxt(fname=algo.csv', delimiter=','). De hecho, podemos pasar el nombre del archivo a loadtxt sin usar fname=:

In [49]:
numpy.loadtxt('precip_mensual_estacion01.csv', delimiter=',')
Out[49]:
array([[0.   , 0.   , 0.   , 0.   , 0.   , 0.   , 0.   , 0.173, 0.021,
        0.035, 0.014, 0.   ],
       [0.005, 0.021, 0.026, 0.001, 0.003, 0.079, 0.113, 0.084, 0.064,
        0.012, 0.007, 0.   ],
       [0.   , 0.001, 0.001, 0.014, 0.032, 0.069, 0.069, 0.115, 0.135,
        0.071, 0.033, 0.   ],
       [0.001, 0.006, 0.005, 0.015, 0.068, 0.085, 0.061, 0.078, 0.07 ,
        0.067, 0.006, 0.015],
       [0.   , 0.003, 0.   , 0.009, 0.103, 0.063, 0.03 , 0.068, 0.098,
        0.009, 0.001, 0.   ],
       [0.005, 0.   , 0.013, 0.009, 0.028, 0.06 , 0.092, 0.076, 0.081,
        0.01 , 0.042, 0.   ],
       [0.   , 0.   , 0.012, 0.016, 0.048, 0.085, 0.165, 0.126, 0.126,
        0.   , 0.   , 0.   ],
       [0.008, 0.001, 0.008, 0.02 , 0.024, 0.059, 0.063, 0.107, 0.065,
        0.046, 0.012, 0.   ]])

pero debemos decir delimiter=:

In [50]:
numpy.loadtxt('precip_mensual_estacion01.csv', ',')
Traceback (most recent call last):

  File "/Users/Karina/anaconda3/lib/python3.7/site-packages/IPython/core/interactiveshell.py", line 3331, in run_code
    exec(code_obj, self.user_global_ns, self.user_ns)

  File "<ipython-input-50-c12d5d880ab3>", line 1, in <module>
    numpy.loadtxt('precip_mensual_estacion01.csv', ',')

  File "/Users/Karina/anaconda3/lib/python3.7/site-packages/numpy/lib/npyio.py", line 1080, in loadtxt
    dtype = np.dtype(dtype)

  File "/Users/Karina/anaconda3/lib/python3.7/site-packages/numpy/core/_internal.py", line 199, in _commastring
    newitem = (dtype, ast.literal_eval(repeats))

  File "/Users/Karina/anaconda3/lib/python3.7/ast.py", line 46, in literal_eval
    node_or_string = parse(node_or_string, mode='eval')

  File "/Users/Karina/anaconda3/lib/python3.7/ast.py", line 35, in parse
    return compile(source, filename, mode, PyCF_ONLY_AST)

  File "<unknown>", line 1
    ,
    ^
SyntaxError: unexpected EOF while parsing

Para enteder qué está pasando veamos la ayuda de la función loadtxt:

In [51]:
help(numpy.loadtxt)
Help on function loadtxt in module numpy:

loadtxt(fname, dtype=<class 'float'>, comments='#', delimiter=None, converters=None, skiprows=0, usecols=None, unpack=False, ndmin=0, encoding='bytes', max_rows=None)
    Load data from a text file.
    
    Each row in the text file must have the same number of values.
    
    Parameters
    ----------
    fname : file, str, or pathlib.Path
        File, filename, or generator to read.  If the filename extension is
        ``.gz`` or ``.bz2``, the file is first decompressed. Note that
        generators should return byte strings.
    dtype : data-type, optional
        Data-type of the resulting array; default: float.  If this is a
        structured data-type, the resulting array will be 1-dimensional, and
        each row will be interpreted as an element of the array.  In this
        case, the number of columns used must match the number of fields in
        the data-type.
    comments : str or sequence of str, optional
        The characters or list of characters used to indicate the start of a
        comment. None implies no comments. For backwards compatibility, byte
        strings will be decoded as 'latin1'. The default is '#'.
    delimiter : str, optional
        The string used to separate values. For backwards compatibility, byte
        strings will be decoded as 'latin1'. The default is whitespace.
    converters : dict, optional
        A dictionary mapping column number to a function that will parse the
        column string into the desired value.  E.g., if column 0 is a date
        string: ``converters = {0: datestr2num}``.  Converters can also be
        used to provide a default value for missing data (but see also
        `genfromtxt`): ``converters = {3: lambda s: float(s.strip() or 0)}``.
        Default: None.
    skiprows : int, optional
        Skip the first `skiprows` lines, including comments; default: 0.
    usecols : int or sequence, optional
        Which columns to read, with 0 being the first. For example,
        ``usecols = (1,4,5)`` will extract the 2nd, 5th and 6th columns.
        The default, None, results in all columns being read.
    
        .. versionchanged:: 1.11.0
            When a single column has to be read it is possible to use
            an integer instead of a tuple. E.g ``usecols = 3`` reads the
            fourth column the same way as ``usecols = (3,)`` would.
    unpack : bool, optional
        If True, the returned array is transposed, so that arguments may be
        unpacked using ``x, y, z = loadtxt(...)``.  When used with a structured
        data-type, arrays are returned for each field.  Default is False.
    ndmin : int, optional
        The returned array will have at least `ndmin` dimensions.
        Otherwise mono-dimensional axes will be squeezed.
        Legal values: 0 (default), 1 or 2.
    
        .. versionadded:: 1.6.0
    encoding : str, optional
        Encoding used to decode the inputfile. Does not apply to input streams.
        The special value 'bytes' enables backward compatibility workarounds
        that ensures you receive byte arrays as results if possible and passes
        'latin1' encoded strings to converters. Override this value to receive
        unicode arrays and pass strings as input to converters.  If set to None
        the system default is used. The default value is 'bytes'.
    
        .. versionadded:: 1.14.0
    max_rows : int, optional
        Read `max_rows` lines of content after `skiprows` lines. The default
        is to read all the lines.
    
        .. versionadded:: 1.16.0
    
    Returns
    -------
    out : ndarray
        Data read from the text file.
    
    See Also
    --------
    load, fromstring, fromregex
    genfromtxt : Load data with missing values handled as specified.
    scipy.io.loadmat : reads MATLAB data files
    
    Notes
    -----
    This function aims to be a fast reader for simply formatted files.  The
    `genfromtxt` function provides more sophisticated handling of, e.g.,
    lines with missing values.
    
    .. versionadded:: 1.10.0
    
    The strings produced by the Python float.hex method can be used as
    input for floats.
    
    Examples
    --------
    >>> from io import StringIO   # StringIO behaves like a file object
    >>> c = StringIO("0 1\n2 3")
    >>> np.loadtxt(c)
    array([[0., 1.],
           [2., 3.]])
    
    >>> d = StringIO("M 21 72\nF 35 58")
    >>> np.loadtxt(d, dtype={'names': ('gender', 'age', 'weight'),
    ...                      'formats': ('S1', 'i4', 'f4')})
    array([(b'M', 21, 72.), (b'F', 35, 58.)],
          dtype=[('gender', 'S1'), ('age', '<i4'), ('weight', '<f4')])
    
    >>> c = StringIO("1,0,2\n3,0,4")
    >>> x, y = np.loadtxt(c, delimiter=',', usecols=(0, 2), unpack=True)
    >>> x
    array([1., 3.])
    >>> y
    array([2., 4.])
    
    This example shows how `converters` can be used to convert a field
    with a trailing minus sign into a negative number.
    
    >>> s = StringIO('10.01 31.25-\n19.22 64.31\n17.57- 63.94')
    >>> def conv(fld):
    ...     return -float(fld[:-1]) if fld.endswith(b'-') else float(fld)
    ...
    >>> np.loadtxt(s, converters={0: conv, 1: conv})
    array([[ 10.01, -31.25],
           [ 19.22,  64.31],
           [-17.57,  63.94]])

Hay mucha información pero en las dos primeras lineas vemos:

loadtxt(fname, dtype=<class 'float'>, comments='#', delimiter=None, converters=None, skiprows=0, usecols=None, unpack=False, ndmin=0, encoding='bytes', max_rows=None)

Esto nos dice que loadtxttiene un parámetro llamado fnameque no tiene un valor default y ocho parámetros que sí tienen valor default. Esto lo sabemos porque el valor default está después del =. fname no tiene esto por lo tanto la función siempre espera que la llames pasándole este parámetro. Si no especificas los otros ocho, entonces éstos tomarán el valor default.

DOCUMENTACIÓN: Podemos agregar información de ayuda a la función de la siguente forma:

In [52]:
def fahr_a_kelvin(temp_f):
    '''Esta función convierte el valor temp_f de Fahrenheit a Kelvin'''
    temp_c = fahr_a_celsius(temp_f)
    temp_k = celsius_a_kelvin(temp_c)
    return temp_k

print('El punto de ebullición  del agua en Kelvin:', fahr_a_kelvin(212.0))
El punto de ebullición  del agua en Kelvin: 373.15
In [53]:
help(fahr_a_kelvin)
Help on function fahr_a_kelvin in module __main__:

fahr_a_kelvin(temp_f)
    Esta función convierte el valor temp_f de Fahrenheit a Kelvin

Pon en práctica lo que aprendiste

  1. Dado el siguiente código:
In [54]:
def func(a, b=3, c=6):
    print('a: ', a, 'b: ', b, 'c:', c)

¿Cuál será el resultado de correr?

func(-1, 2)

  1. Si la variable s es una cadena, entonces s[0] es el primer caracter de la cadena y s[-1] es el último. Escribe una función que se llame extremos que regrese una cadena constituida por el primer y úlimo elementos del input de la función. La llamada a tu función deberá verse como:

print(extremos('helio'))

Output: hm

8.3 En resumen

  • Define un función usando def nombre_de_función(parámetros)
  • El cuerpo de la función debe estar indentado.
  • Llama a una función usando nombre_de_función(valor)
  • Las variables definidas dentro de un a función solo pueden ser vistas y usadas dentro de esa función.
  • Si una variable no está definida dentro de la función, Python buscará la definición dentro del código que se encuentra antes de la llamada a la función.
  • Usa help(algo) para ver la ayuda de algo.
  • Agrega texto para documentar tu función usando ''' texto '''.
  • Especifica valores por default para los parámetros usando parámetro=valor default en la lista de parámetros dentro de la definición de la función.
  • Los parámetros se pueden pasar a una función basados en el nombre, la posición u omitiéndolos, en cuyo caso se usa el valor default.

9. Archivos NetCDF y más de matplotlib

Finalmente aprenderemos qué son los archivos en formato netCDF y aprovecharemos para explorar más funciones de visualización de la biblioteca matplotlib. Esta breve introducción a netCDF está basada en los tutoriales READING NETCDF4 DATA IN PYTHON y WRITING NETCDF DATA IN PYTHON de Uptal Kumar en el portal Earth Inversion y Python - NetCDF reading and writing example with plotting de Chris Slocum. Les recomiendo que los revisen completos porque tienen temas adicionales como el manejo de fechas, datos nulos y cómo usar la biblioteca xarray, que es súper útil y rápida.

En Ciencias de la Tierra es muy común lidiar con estructuras de datos de muchas dimensiones (generalmente 3 espaciales +1 temporal) como el estado del océano, la atmósfera, el interior de un planeta, etc. Es impráctico guardar estos datos en archivos de texto (como los .csv que hemos estado usando) porque necesitarimos mucha capacidad de memoria para guardarlos, leerlos y procesarlos. Una de las mejores herramientas para manipular datos multidimensionales son los netCDF. Estos archivos almacenan datos en formato HDF (Hierarchical Data Format).

En esta sección aprenderemos a leer y escribir datos netCDF pero antes de eso necesitamos instalar el paquete netcdf4 de python.

Ve a la terminal y escribe el siguente comando:

conda install -c conda-forge netcdf4

Cuando haya terminado la instalación continúa con la sección.

In [55]:
# Importa la biblioteca netCDF4
import netCDF4

Ahora leamos un archivo netCDF con datos de temperatura del aire cerca de la superficie en todo el planeta a lo largo de 1 año que provienen de la salida de un modelo numérico de NOAA. Los archivos NetCDF tienen extensión .nc.

In [56]:
f = netCDF4.Dataset('air.sig995.2012.nc')
print(f)
<class 'netCDF4._netCDF4.Dataset'>
root group (NETCDF3_CLASSIC data model, file format NETCDF3):
    Conventions: COARDS
    title: mean daily NMC reanalysis (2012)
    history: created 2011/12 by Hoop (netCDF2.3)
    description: Data is from NMC initialized reanalysis
(4x/day).  These are the 0.9950 sigma level values.
    platform: Model
    references: http://www.esrl.noaa.gov/psd/data/gridded/data.ncep.reanalysis.html
    dimensions(sizes): lat(73), lon(144), time(366)
    variables(dimensions): float32 lat(lat), float32 lon(lon), float64 time(time), int16 air(time,lat,lon)
    groups: 
/Users/Karina/anaconda3/lib/python3.7/site-packages/ipykernel_launcher.py:2: DeprecationWarning: tostring() is deprecated. Use tobytes() instead.
  

En este caso leímos el archivo “air.sig995.2012.nc”. Cuando imprimimos el objeto f podemos ver que el archivo es la clase NETCDF3. Hay más información como título, intitución, etc. Esto se conoce como metadatos.

Al final vemos las dimensiones y variables del conjunto de datos. Este conjunto de datos tiene 3 dimensiones: lat (size:73), lon(size:144), time(size:366). Luego tenemos las variables, que están definidas con base en las dimensiones anteriores. En este caso tenemos lon, lat, time y aitr basadas en las dimensiones lat, lon y time. Podemos tener variables basadas en solo una dimensión, como time o en tres como air.

Podemos acceder a la información de este objeto como leeríamos un diccionario en Python:

In [57]:
print(f.variables.keys()) # Imprime los nombres de todas las variables en f
dict_keys(['lat', 'lon', 'time', 'air'])

Esto imprime los nombres de todas las variables que se leyeron del archivo netCDF referidas por el objeto "f".

Podemos acceder a cada variable individualemente:

In [58]:
dep = f.variables['air'] 
print(dep)
<class 'netCDF4._netCDF4.Variable'>
int16 air(time, lat, lon)
    long_name: mean Daily Air temperature at sigma level 995
    unpacked_valid_range: [185.16 331.16]
    actual_range: [193.80002 317.52002]
    units: degK
    add_offset: 512.81
    scale_factor: 0.01
    missing_value: 32766
    precision: 2
    least_significant_digit: 1
    GRIB_id: 11
    GRIB_name: TMP
    var_desc: Air temperature
    dataset: NCEP Reanalysis Daily Averages
    level_desc: Surface
    statistic: Mean
M
    parent_stat: Individual Obs
I
    valid_range: [-32765 -18165]
unlimited dimensions: time
current shape = (366, 73, 144)
filling on, default _FillValue of -32767 used
/Users/Karina/anaconda3/lib/python3.7/site-packages/ipykernel_launcher.py:2: DeprecationWarning: tostring() is deprecated. Use tobytes() instead.
  

La variable "air" tiene 3 dimensiones - lon,lat, time. También podemos obtener más información (metadatos) como las coordenadas, el tiempo o las unidades de la variable. Las coordenadas son las variables 1D que tienen el mismo nombre que las dimensiones. La unidad de la variable air es Kelvin (degK). La forma o dimensiones nos dan imformación acerca del tamaño de la variable. Aqui la forma (shape) es (366, 73, 144).

También podemos revisar el tamaño de las dimensiones individualmente:

In [59]:
for d in f.dimensions.items():
  print(d)
('lat', <class 'netCDF4._netCDF4.Dimension'>: name = 'lat', size = 73)
('lon', <class 'netCDF4._netCDF4.Dimension'>: name = 'lon', size = 144)
('time', <class 'netCDF4._netCDF4.Dimension'> (unlimited): name = 'time', size = 366)
/Users/Karina/anaconda3/lib/python3.7/site-packages/ipykernel_launcher.py:2: DeprecationWarning: tostring() is deprecated. Use tobytes() instead.
  

Cuando tenemos "tiempo" como dimensión es común que ésta sea tipo unlimited (ilimitadas). Esto significa que el tamaño de esa dimensión puede incrementar indefinidamente. En este caso todas las dimsnsiones están fijas.

Si solo queremos encontrar las dimensiones de una variable particular podemos hacer:

In [60]:
print(dep.dimensions)
print(dep.shape)
('time', 'lat', 'lon')
(366, 73, 144)

Podemos inspeccionar las variables asociadas a las dimensiones de la siguiente forma:

In [61]:
time = f.variables['time']
air = f.variables['air']
lon,lat = f.variables['lon'], f.variables['lat']
print(time)
print(lon)
print(lat)
<class 'netCDF4._netCDF4.Variable'>
float64 time(time)
    units: hours since 1-1-1 00:00:0.0
    long_name: Time
    actual_range: [17628096. 17636856.]
    delta_t: 0000-00-01 00:00:00
    standard_name: time
    axis: T
    avg_period: 0000-00-01 00:00:00
unlimited dimensions: time
current shape = (366,)
filling on, default _FillValue of 9.969209968386869e+36 used
<class 'netCDF4._netCDF4.Variable'>
float32 lon(lon)
    units: degrees_east
    long_name: Longitude
    actual_range: [  0.  357.5]
    standard_name: longitude
    axis: X
unlimited dimensions: 
current shape = (144,)
filling on, default _FillValue of 9.969209968386869e+36 used
<class 'netCDF4._netCDF4.Variable'>
float32 lat(lat)
    units: degrees_north
    actual_range: [ 90. -90.]
    long_name: Latitude
    standard_name: latitude
    axis: Y
unlimited dimensions: 
current shape = (73,)
filling on, default _FillValue of 9.969209968386869e+36 used
/Users/Karina/anaconda3/lib/python3.7/site-packages/ipykernel_launcher.py:4: DeprecationWarning: tostring() is deprecated. Use tobytes() instead.
  after removing the cwd from sys.path.
/Users/Karina/anaconda3/lib/python3.7/site-packages/ipykernel_launcher.py:5: DeprecationWarning: tostring() is deprecated. Use tobytes() instead.
  """
/Users/Karina/anaconda3/lib/python3.7/site-packages/ipykernel_launcher.py:6: DeprecationWarning: tostring() is deprecated. Use tobytes() instead.
  

Ahora,¿cómo podemos acceder a los datos guaradados en cada variable? Las variables NetCDF se comportan de forma similar a arreglos de Numpy: podemos rebanarlas, enmascararlas y hacer cálculos con ellas.

Por ejemplo, leámos los datos de la variable lon:

In [62]:
lon_data = lon[:]
print(lon_data)

lat_data = lat[:]
print(lat_data)
[  0.    2.5   5.    7.5  10.   12.5  15.   17.5  20.   22.5  25.   27.5
  30.   32.5  35.   37.5  40.   42.5  45.   47.5  50.   52.5  55.   57.5
  60.   62.5  65.   67.5  70.   72.5  75.   77.5  80.   82.5  85.   87.5
  90.   92.5  95.   97.5 100.  102.5 105.  107.5 110.  112.5 115.  117.5
 120.  122.5 125.  127.5 130.  132.5 135.  137.5 140.  142.5 145.  147.5
 150.  152.5 155.  157.5 160.  162.5 165.  167.5 170.  172.5 175.  177.5
 180.  182.5 185.  187.5 190.  192.5 195.  197.5 200.  202.5 205.  207.5
 210.  212.5 215.  217.5 220.  222.5 225.  227.5 230.  232.5 235.  237.5
 240.  242.5 245.  247.5 250.  252.5 255.  257.5 260.  262.5 265.  267.5
 270.  272.5 275.  277.5 280.  282.5 285.  287.5 290.  292.5 295.  297.5
 300.  302.5 305.  307.5 310.  312.5 315.  317.5 320.  322.5 325.  327.5
 330.  332.5 335.  337.5 340.  342.5 345.  347.5 350.  352.5 355.  357.5]
[ 90.   87.5  85.   82.5  80.   77.5  75.   72.5  70.   67.5  65.   62.5
  60.   57.5  55.   52.5  50.   47.5  45.   42.5  40.   37.5  35.   32.5
  30.   27.5  25.   22.5  20.   17.5  15.   12.5  10.    7.5   5.    2.5
   0.   -2.5  -5.   -7.5 -10.  -12.5 -15.  -17.5 -20.  -22.5 -25.  -27.5
 -30.  -32.5 -35.  -37.5 -40.  -42.5 -45.  -47.5 -50.  -52.5 -55.  -57.5
 -60.  -62.5 -65.  -67.5 -70.  -72.5 -75.  -77.5 -80.  -82.5 -85.  -87.5
 -90. ]

Ahora, grafiquemos los datos de temperatura del aire guardados en la variable air. Esta variable tiene 3 dimensiones (time, lat, lon), por lo que podemos hacer representaciones de los datos en 2D rebanando el arreglo. Hay varias formas de representar estas rebanadas gráficamente. Probemos las funciones contourf y pcolormesh de matplotlib. Para ello rebanemos los datos de air sleccionando todas las latitudes y longitudes para un solo día, es decir, pasaremos de un arreglo tamaño (366, 73, 144) a uno tamaño (73, 144), que sí podemos visualizar:

In [63]:
time_index = 100 # elegimos un índice temporal al azar para visualizar 
lat_data = lat[:]    
air_data = air[time_index,:,:] # rebanamos los datos de air, eligiendo solo un día 
                               # y todas las lats y lons
In [64]:
# contourf y pcolormesh
fig, (ax1,ax2) = matplotlib.pyplot.subplots(1,2, figsize=(13,4))
ax1.pcolormesh(lon_data, lat_data,  air_data)
ax1.set_xlabel('lon')
ax1.set_ylabel('lat')
ax1.set_title('Pcolormesh')
ax2.contourf(lon_data, lat_data,  air_data, 10) # dibuja 10 contornos
ax2.set_xlabel('lon')
ax2.set_ylabel('lat')
ax2.set_title('Contourf')
Out[64]:
Text(0.5, 1.0, 'Contourf')

Ambos graficos tienen similitudes y diferencias. Pcolormesh genera una función que asigna un color, dentro de un mapa de colores, al rango de datos dentro de los cuales se encuentran los valores del arreglo. Dependiendo del valor del elemento se le asigna un color. En cambio, contourf genera contornos, como curvas de nivel, para los valores del arreglo. Tú puedes elegir el número de contornos a graficar y, para ambas funcones, el mapa de color que quieras utilizar.

Ahora, necesitamos una barra de color para poder traducir de colores a temperaturas. Por ejemplo ¿qué significa que una región sea amarilla o azul?

In [65]:
# contourf y pcolormesh
fig, (ax1,ax2) = matplotlib.pyplot.subplots(1,2, figsize=(13,4))
pc = ax1.pcolormesh(lon_data, lat_data,  air_data)
matplotlib.pyplot.colorbar(pc, ax=ax1, label='Temperatura del aire (K)')
ax1.set_xlabel('lon')
ax1.set_ylabel('lat')
ax1.set_title('Pcolormesh')

pc2 = ax2.contourf(lon_data, lat_data,  air_data)
matplotlib.pyplot.colorbar(pc2, ax=ax2, label='Temperatura del aire (K)')
ax2.set_xlabel('lon')
ax2.set_ylabel('lat')
ax2.set_title('Contourf')
Out[65]:
Text(0.5, 1.0, 'Contourf')

Nota como las barras de color son distintas. La barra de pcolormesh es continua mientras que la de contourf es discreta, mostrando todos los intervalos representados en el gráfico.

Pon en práctica lo que aprendiste

  1. Reproduce la figura anterior cambiando los mapas de color con el argumento opcional cmap=así como los límites del mapa de color usando los argumentos vmax= y vmin=. Recuerda que puedes consultar la documentación para más información.
In [70]:
# contourf y pcolormesh
fig, (ax1,ax2) = matplotlib.pyplot.subplots(1,2, figsize=(13,4))
pc = ax1.pcolormesh(lon_data, lat_data,  air_data,cmap='gist_heat', vmin=200, vmax=315,)
matplotlib.pyplot.colorbar(pc, ax=ax1, label='Temperatura del aire (K)')
ax1.set_xlabel('lon')
ax1.set_ylabel('lat')
ax1.set_title('Pcolormesh')

pc2 = ax2.contourf(lon_data, lat_data,  air_data, cmap='gist_heat', vmin=200, vmax=315,)
matplotlib.pyplot.colorbar(pc2, ax=ax2, label='Temperatura del aire (K)')
ax2.set_xlabel('lon')
ax2.set_ylabel('lat')
ax2.set_title('Contourf')
Out[70]:
Text(0.5, 1.0, 'Contourf')

Pon en práctica lo que aprendiste

Siguiendo el tutorial WRITING NETCDF DATA IN PYTHON crea tu propio archivo netcdf. Puede tener los datos aleatorios que generan en el tutorial o tu propio set de datos (puedes usar los del ejemplo anterior).

  1. Lee el archivo que generaste
  2. Imprime las variables del archivo.
  3. Grafica los datos, ya sea en 2D (contourf o pcolormesh) o como gráficos de líneas.
In [ ]: