from matplotlib.pylab import *

# Primeiro gráfico.
plot([1, 2, 3], 'k--')
xlabel("1, 2, 3")
ylabel(u"Contando 1, 2, 3")
title(u"Nosso primeiro plot!")

plot?

plot([1, 2, 3], 'ro')
xlabel("a, b, c")
ylabel(u"?")
title(u"?")
axis([-0.1, 2.1, 0.8, 3.1])

import os
from glob import glob

import gsw
import numpy as np
import matplotlib.pyplot as plt

from netCDF4 import Dataset
from oceans.ff_tools import wrap_lon360
from mpl_toolkits.basemap import Basemap, shiftgrid, cm

def break_lines(line):
    return [float(num) for num in line.strip().split()]

def get_topo(url, m):
    etopodata = Dataset(url)
    try:
        topoin = etopodata.variables['ROSE'][:]
        lons = etopodata.variables['ETOPO05_X'][:]
        lats = etopodata.variables['ETOPO05_Y'][:]
        topoin, lons = shiftgrid(180., topoin, lons, start=False)
    except:
        topoin = etopodata.variables['topo'][:]
        lons = etopodata.variables['topo_lon'][:]
        lons = wrap_lon360(lons)
        lats = etopodata.variables['topo_lat'][:]        
        topoin, lons = shiftgrid(180., topoin, lons, start=False)
    nx = int((m.xmax - m.xmin) / 5000.) + 1
    ny = int((m.ymax - m.ymin) / 5000.) + 1
    return m.transform_scalar(topoin, lons, lats, nx, ny)


def bathymetry_angle(start=0, end=9):
    opposite = depth[end] - depth[start]
    adjacent = dist[end] - dist[start]
    return np.arctan2(-opposite, adjacent) * 180 / np.pi

lista = sorted(glob(os.path.join('./data/cruise', '*.dat')))

depth, lat, lon = [], [], []
for fname in lista:
    with open(fname, 'r') as f:
        lines = f.readlines()
        d, la, lo = break_lines(lines[0])[3:]
        lon.append(lo)
        lat.append(la)
        depth.append(d)

print('\n'.join(lista))

dist = gsw.distance(lon, lat)
dist = np.r_[0, np.cumsum(dist)] / 1e3
tangent = bathymetry_angle(0, 9)
print("Shelf angle: %.4f" % tangent)

fig, ax = plt.subplots(figsize=(7, 4), facecolor='w')
ax.fill_between(dist, max(depth), depth, color='0.75')
ax.plot(dist, depth, color='k')
ax.set_ylim(0, max(depth))
ax.set_xlim(-5, max(dist) + 10)
ax.invert_yaxis()
ax.set_ylabel("Depth [db]")
ax.set_xlabel("Distance [km]")
ax.plot(dist, len(dist) * [0], 'kv')
ax.xaxis.tick_bottom()

offset = 2.
llcrnrlon, urcrnrlon = min(lon) - offset, max(lon) + offset
llcrnrlat, urcrnrlat = min(lat) - offset, max(lat) + offset

m = Basemap(llcrnrlon=llcrnrlon, urcrnrlon=urcrnrlon,
            llcrnrlat=llcrnrlat, urcrnrlat=urcrnrlat,
            resolution='i', projection='merc')
fig, ax = plt.subplots(figsize=(7, 7))
m.ax = ax

try:
    topodata = get_topo('/home/filipe/00-NOBKP/OcFisData/etopo5.nc', m)
except:
    topodata = get_topo('http://ferret.pmel.noaa.gov/thredds/dodsC/data/PMEL/etopo5.nc', m)
    

im = m.imshow(topodata, cm.GMT_haxby)
m.drawstates()
m.drawcountries()
m.drawcoastlines()
m.fillcontinents()
meridians = np.arange(int(llcrnrlon), int(urcrnrlon), 1.)
parallels = np.arange(int(llcrnrlat), int(urcrnrlat), 1.)

m.plot(lon, lat, 'r*', latlon=True)
cb = m.colorbar(im, 'right', size='5%', pad='2%')
kw = dict(fontsize=12, fontweight='demibold')
m.drawparallels(parallels, labels=[1, 0, 0, 0], **kw)
_ = m.drawmeridians(meridians, labels=[0, 0, 0, 1], **kw)

import os
import numpy as np
import matplotlib.pyplot as plt
import mpl_toolkits.axisartist as AA
from mpl_toolkits.axes_grid1 import host_subplot

deg = u"\u00b0"

def plot_station(est=15):
    if est < 10:
        est = '0%s' % est
    fname = 'estacao%s.dat' % est
    fname = os.path.join('./data/cruise', fname)

    P, T, S = np.loadtxt(fname, unpack=True, usecols=(0, 3, 4),
                        skiprows=1)

    fig = plt.figure(figsize=(4.2, 6))

    ax0 = host_subplot(111, axes_class=AA.Axes)
    new = ax0.get_grid_helper()
    ax0.axis["bottom"] = new.new_fixed_axis(loc="bottom",
                        offset=(0, -40),
                        axes=ax0)
    p0, = ax0.plot(T, P, linewidth=2.0,
                color='green', label=r'Temperature')

    ax0.set_ylabel("Pressure [dbar]")
    ax0.set_xlabel(u"Temperature[%sC]" % deg)
    ax0.axis["bottom"].label.set_color(p0.get_color())
    ax0.set_title("CTD %s" % est)

    ax1 = ax0.twiny()
    new = ax1.get_grid_helper()
    ax1.axis["top"] = new.new_fixed_axis(loc="bottom",
                                                    offset=(0, 0),
                                                    axes=ax1)
    p1, = ax1.plot(S, P, linewidth=2.0,
                color='blue', label=r'Salinity')

    ax1.set_xlabel("Salinity [Kg g$^{-1}$]")
    ax1.axis["top"].label.set_color(p1.get_color())
    ax1.invert_yaxis()
    ax1.axis("tight")

    ax0.legend(shadow=True, fancybox=True, numpoints=1,
            loc=2, bbox_to_anchor=(.6, 0.2))

plot_station()

import os
from glob import glob

import gsw
import numpy as np
import matplotlib.pyplot as plt

from pandas import DataFrame

def basename(fname):
    return os.path.splitext(os.path.basename(fname))

def read_station(fname):
    P, T, S = np.loadtxt(fname, unpack=True, usecols=(0, 3, 4),
                         skiprows=1)
    column = basename(fname)[0]
    return (DataFrame(T, index=P, columns=[column]),
            DataFrame(S, index=P, columns=[column]))

lista = sorted(glob(os.path.join('./data/cruise', '*.dat')))

first = lista.pop(0)
tmp, sal = read_station(first)

for fname in lista:
    T, S = read_station(fname)
    tmp = tmp.join(T, how='right')
    sal = sal.join(S, how='right')

Te = np.arange(0, 32, 2)
Se = np.arange(32, 38.25, 0.25)

Sg, Tg = np.meshgrid(Se, Te)
cnt = np.arange(20, 34)

sigma_theta = gsw.sigma0_pt0_exact(Sg, Tg)

fig, ax = plt.subplots(figsize=(6, 6))
ax.plot(sal, tmp, 'k.')
ax.set_ylabel(u"Temperatura u\u00b0C")  # Unicode.
ax.set_xlabel(r"Salinidade [g kg$^{-1}$]")  # Latex.

cs = ax.contour(Se, Te, sigma_theta, colors='black', levels=cnt)
ax.clabel(cs, fontsize=9, inline=1, fmt='%2.1f')
_ = ax.axis([31.8, 37.2, 0.0, 30.0])

HTML(html)