In [ ]:
import iris
from datetime import datetime, timedelta

url = 'http://aviso-users:grid2010@opendap.aviso.oceanobs.com/thredds/dodsC/dataset-duacs-dt-upd-global-merged-msla-h'

time = lambda cell: datetime.today() <= cell <= datetime.today() - timedelta(days=365)
longitude = lambda cell: 300 <= cell <= 360
latitude = lambda cell: -60 <= cell <= 0
constraint = iris.Constraint(coord_values=dict(longitude=longitude, latitude=latitude, name='H'))
constraint = iris.Constraint(coord_values=dict(longitude=longitude, latitude=latitude))

cube = iris.load_cube(url, constraint)
print(cube)

In [ ]:
import os
if os.path.isfile('../../data/altimetry.nc'):
    cube = iris.load_cube('../../data/altimetry.nc')
else:
    cube.convert_units('meters')
    iris.save(cube, 'altimetry.nc')

time = cube.coord('time')
lon = cube.coord('longitude').points
lat = cube.coord('latitude').points
time_vec = time.units.num2date(time.points)
print(cube)

Vamos plotar o tempo mais próximo de hoje para checar os dados


In [ ]:
import iris.plot as iplt
import cartopy.crs as ccrs
import iris.quickplot as qplt
import matplotlib.pyplot as plt
from cartopy.feature import NaturalEarthFeature, LAND
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER

today = datetime.today()
itime = time.nearest_neighbour_index(time.units.date2num(today))

def plt_ssh(cube, itime):
    fig, ax = plt.subplots(figsize=(6, 6),
                           subplot_kw=dict(projection=ccrs.PlateCarree()))
    qplt.contourf(cube[itime, ...], cmap=plt.cm.GnBu, extend='both')
    ax.add_feature(LAND)
    ax.coastlines('50m', color='k')
    gl = ax.gridlines(crs=ccrs.PlateCarree(), draw_labels=True,
                      linewidth=1.5, color='gray', alpha=0.5, linestyle='--')
    gl.xlabels_top = False
    gl.ylabels_right = False
    gl.xformatter = LONGITUDE_FORMATTER
    gl.yformatter = LATITUDE_FORMATTER
    _ = iplt.citation(time_vec[itime].strftime('%Y-%m-%d'))
    return fig, ax

In [ ]:
fig, ax = plt_ssh(cube, itime)

In [ ]:
from JSAnimation import IPython_display
from matplotlib.animation import FuncAnimation
from matplotlib.colors import LinearSegmentedColormap

levels = [-0.8, -0.6, -0.4, -0.2, 0, 0.2, 0.4, 0.6, 0.8, 0.9]  # Must be even so zero is white.
cmap = LinearSegmentedColormap.from_list(name='brasil',
                                                 colors =[(1, 1, 0),  (1., 1., 1.),  (0, 1, 0)],
                                                 N=len(levels)-1)
def draw_coastline(ax):
    ax.coastlines('50m', color='k')
    gl = ax.gridlines(crs=ccrs.PlateCarree(), draw_labels=True,
                      linewidth=1.5, color='gray', alpha=0.5, linestyle='--')
    gl.xlabels_top = False
    gl.ylabels_right = False
    gl.xformatter = LONGITUDE_FORMATTER
    gl.yformatter = LATITUDE_FORMATTER
    ax.add_feature(LAND)
    states = NaturalEarthFeature(category='cultural', scale='50m', facecolor='none',
                                 name='admin_1_states_provinces_shp')
    ax.add_feature(states, edgecolor='gray')

In [ ]:
kw = dict(cmap=cmap, clim=(-1, 1), levels=levels)
fig, ax = plt.subplots(figsize=(6, 6), subplot_kw=dict(projection=ccrs.PlateCarree()))
ax.axis([-59, 0, -59, 0])
cs = ax.contourf(lon, lat, cube[0, ...].data.T, **kw)
fig.colorbar(cs, shrink=0.7)
text = ax.text(-7.5, -7.5, time_vec[0].strftime('%Y-%m-%d'), ha='center', va='center')
draw_coastline(ax)

def update_contourf(k):
    global ax, fig
    ax.cla()
    cs = ax.contourf(lon, lat, cube[k, ...].data.T, **kw)
    text = ax.text(-7.5, -7.5, time_vec[k].strftime('%Y-%m-%d'), ha='center', va='center')
    draw_coastline(ax)
    return cs

anim = FuncAnimation(fig, update_contourf, interval=100, frames=52)
# anim.save('ssh.mp4', fps=30, extra_args=['-vcodec', 'libx264'])
anim

In [ ]:
dist_degree = 15
degree2km = 111
time_weeks = 43
week2hour = 7 * 24

total_time = (dist_degree * degree2km) / (time_weeks * week2hour)
print("Tempo que leva para uma anomalia atravesar 15 graus de longitude é %0.2f km/h" % total_time)

In [ ]:
time_years = (time_weeks * 7 * 3) / 365.
print("Para chegar até São Paulo são 45 graus de longitude,\n"
      "o tempo total de viagem desde 0 graus foi de %1.2f anos" % time_years)