import os
from zipfile import ZipFile
import datetime as dt
import numpy as np
from matplotlib import pyplot as plt
import netCDF4 as nc
from lxml import etree
from pykml.factory import nsmap
from pykml.factory import KML_ElementMaker as KML
from pykml.factory import GX_ElementMaker as GX
# Desde la versión 5.0, google ha añadido una serie de extensiones
# a kml que lo hacen más potente e interactivo
# https://developers.google.com/kml/documentation/kmlreference?hl=es#kmlextensions
# Definimos una variable para el espacio de nombres de las
# extensiones de google que vamos a usar.
gxns = '{' + nsmap['gx'] + '}'
# Ponemos la información referente al huracán:
# La posición en longitud cada 6h desde que es tormenta tropical
lon_hur = [-30.3, -32.1, -33.6, -35, -36.5, -38.2, -39.9,
-41.4, -43.4, -45.1, -46.8, -48.5, -50.5, -52.5,
-54.4, -56.1, -57.8, -59.4, -61.1, -62.6, -64.1,
-65.5, -67, -68.3, -69.5, -70.8, -71.9, -72.8,
-73.8, -74.7, -75.8, -76.5, -77.6, -78.4, -79,
-79.6, -80.4, -81.2, -82.1, -82.8, -83.5, -84.1,
-84.7, -85.1, -85.6, -86, -86.5, -87, -87.4,
-87.9, -88.2, -88.2, -87.9, -87.7, -87.4, -86.5,
-85.7, -84, -82.3, -80.5, -78.5, -76.7, -75.5,
-74, -74, -74.5, -75.8, -77.5, -78.5, -78.7, -79.1,
-79.7, -80.6, -81.7, -82.8, -84.1, -86.1, -87.3,
-88.6, -89.5, -91, -92.2, -92.7, -93.2, -94.2]
# La posición en latitud cada 6h desde que es tormenta tropical
lat_hur = [9.7, 9.5, 9.3, 9.1, 8.9, 8.9, 9, 9.3, 9.5, 9.8,
10.2, 10.6, 10.8, 11, 11.3, 11.2, 11.3, 11.6,
11.8, 12, 12.3, 12.6, 13, 13.3, 13.7, 14.2, 14.7,
15.2, 15.7, 16.2, 16.8, 17.3, 17.4, 17.7, 18,
18.2, 18.4, 18.8, 19.1, 19.5, 19.9, 20.4, 20.9,
21.6, 22.4, 23, 23.7, 24.7, 25.6, 26.7, 27.9,
28.9, 30, 31.4, 32.5, 33.8, 34.7, 35.4, 36.2,
37, 37.7, 38.4, 38, 37.5, 36, 34.5, 32.8, 31,
29, 27.5, 26.4, 26.1, 25.9, 25.8, 25.2, 24.8,
25.1, 26, 26.5, 27.1, 27.9, 28.9, 29.2, 29.6, 30.1]
# El estado del huracan cada 6h
tipo = ['TS','TS','TS','TS','TS','TS','TS',
'TS','H1','H2','H3','H4','H3','H3',
'H2','H2','H2','H3','H3','H4','H4',
'H4','H4','H4','H5','H5','H4','H4',
'H4','H4','H4','H4','H4','H4','H5',
'H5','H4','H4','H4','H5','H5','H5',
'H5','H5','H5','H4','H4','H4','H4',
'H4','H4','H3','H3','H1','TS','TD',
'TD','TD','TD','TD','TD','TD','TD',
'TS','TS','TS','TS','TS','TS','TD',
'TD','TD','TD','TD','TD','TD','TD',
'TD','TD','TS','TS','TS','TS','TD','TD']
# La fecha cada 6h desde el 21/09/1998 a las 18Z
fecha_inicial = dt.datetime(2004, 9, 3, 6, 0, 0)
deltat = dt.timedelta(hours = 6)
fechas = [fecha_inicial + deltat * i for i in range(len(tipo))]
# El icono que vamos a asignar a cada 6h, en función del estado de la tormenta.
iconos = {'TD':'http://www.srh.noaa.gov/gis/kml/hurricanetrack/hurrimages/td.gif',
'TS':'http://www.srh.noaa.gov/gis/kml/hurricanetrack/hurrimages/ts.gif',
'H1':'http://www.srh.noaa.gov/gis/kml/hurricanetrack/hurrimages/h1.gif',
'H2':'http://www.srh.noaa.gov/gis/kml/hurricanetrack/hurrimages/h2.gif',
'H3':'http://www.srh.noaa.gov/gis/kml/hurricanetrack/hurrimages/h3.gif',
'H4':'http://www.srh.noaa.gov/gis/kml/hurricanetrack/hurrimages/h4.gif',
'H5':'http://www.srh.noaa.gov/gis/kml/hurricanetrack/hurrimages/h5.gif'}
# Función que crea las imágenes
def pinta_mapas(nombre, lon, lat, mslp, u, v, loni, lati):
fig = plt.figure()
fig.set_size_inches(8, 8. * lat.shape[0] / lon.shape[0])
ax = plt.Axes(fig, [0., 0., 1., 1.])
ax.set_axis_off()
fig.add_axes(ax)
ax.contour(lon, lat, mslp / 100., np.arange(900,1100.,2.),
colors='y',linewidths=2, aspect='normal')
wspd = np.sqrt(u * u + v * v) * 3.6
ax.contourf(lon, lat, wspd, np.arange(50, 250, 25),
colors = plt.hot(), aspect='normal')
plt.savefig(nombre, dpi = 80, transparent=True)
# Leemos los datos de presión, velocidad, latitud y longitud.
data = nc.Dataset('Usando google earth con ayuda de python y pykml (II)/data.nc')
lon = data.variables['longitude'][:]
lat = data.variables['latitude'][:]
# presión = mslp por mean sea level pressure,
# presión media al nivel del mar
mslp = data.variables['msl'][1:-2,:,:]
lon[lon > 180] = lon - 360
u = data.variables['u10'][1:-2,:,:]
v = data.variables['v10'][1:-2,:,:]
# Como vamos a crear un fichero kmz que incluirá todas las imágenes
# abrimos el fichero kmz donde iremos guardándolo todo
fich_kmz = ZipFile('HuracanIvan.kmz', 'w')
# Desde la versión 5.0, google ha añadido una serie de extensiones
# a kml que lo hacen más potente e interactivo
# https://developers.google.com/kml/documentation/kmlreference?hl=es#kmlextensions
# Definimos una variable para el espacio de nombres de las
# extensiones de google que vamos a usar.
gxns = '{' + nsmap['gx'] + '}'
# start with a base KML tour and playlist
fich_kml = KML.kml(
KML.Document(
GX.Tour(
KML.name(u"¡Reprodúceme!"),
GX.Playlist(),
),
KML.Folder(
KML.name('Ruta del huracan Ivan'),
id='lugares',
),
)
)
# Hacemos un bucle para recorrer todos los 'momentos'
# del huracán definidos en el bloque de código anterior
for i in range(len(tipo)):
# Primero volamos hasta cada posición del huracán (cada 6h)
# y nos quedamos observándolo desde el espacio
# A unos 2000 km de altura.
# Para ello usamos FlyTo donde especificamos
# La duración en segundos, si queremos que vaya más rápido o más
# lento, el modo de vuelo, que puede ser más suave (smooth) o
# más rápido (bounce) y la posición.
# La posición la determinan longitude, latitude y altitude,
# tilt es el ángulo desde la vertical (en este caso lo vemos
# verticalmente (0º)) y range es la distancia desde la posición
# determinada por longitude, latitude y altitude, teniendo en
# cuenta el ángulo tilt (en este caso hemos puesto 2.000 km)
# Las unidades para latitude y longitude son grados, para altitude
# y range son metros y para tilt son grados desde la vertical.
# En el gráfico de este enlace se verá mejor:
# https://developers.google.com/kml/documentation/images/lookAt.gif
fich_kml.Document[gxns+"Tour"].Playlist.append(
GX.FlyTo(
GX.duration(1),
GX.flyToMode("bounce"),
KML.LookAt(
KML.longitude(lon_hur[i]),
KML.latitude(lat_hur[i]),
KML.altitude(0),
KML.heading(0),
KML.tilt(0),
KML.range(10000000.0),
KML.altitudeMode("relativeToGround"),
)
),
)
fich_kml.Document[gxns+"Tour"].Playlist.append(GX.Wait(GX.duration(0.1)))
# Añadimos información de cada lugar en la carpeta
# como 'placemarks' o marcas de posición. En el nombre y en la descripción
# del placemark se puede meter HTML y CSS, por si queréis ser un poco más
# creativos y dejar la información de una forma visible y bonita.
# extrude == 1 indica que el nombre del lugar se hará visible. Si queremos
# hacerlo visible solo cuando lo seleccionemos pues le ponemos el valor 0
desc = 'lon: {1}
lat: {2}
Estado del huracan: {3}]]>'
desc = desc.format(fechas[i].isoformat()[0:13], lon_hur[i], lat_hur[i], tipo[i])
fich_kml.Document.Folder.append(
KML.Placemark(
KML.name(fechas[i].isoformat()[0:13]),
KML.description(
"