In [ ]:
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.basemap import Basemap
def make_map(lonStart=-48, lonEnd=-32, latStart=-30, latEnd=-18, img=None):
m = Basemap(projection='merc', llcrnrlon=lonStart, urcrnrlon=lonEnd,
llcrnrlat=latStart, urcrnrlat=latEnd, resolution='c',
lat_ts=(latStart + latEnd) / 2.)
fig, ax = plt.subplots(figsize=(6, 6), facecolor='w')
m.ax = ax
image = plt.imread(img)
m.imshow(image, origin='upper', alpha=0.75)
lonStart, latStart = -42, -26 # Crop the image.
lon_lim, lat_lim = m([lonStart, lonEnd], [latStart, latEnd])
m.ax.axis([lon_lim[0], lon_lim[1], lat_lim[0], lat_lim[1]])
dx = dy = 2
meridians = np.arange(lonStart, lonEnd + dy, dy)
parallels = np.arange(latStart, latEnd + dx, dx)
xoffset = -lon_lim[0] + 1e4
yoffset = -lat_lim[0] + 1e4
kw = dict(linewidth=0)
m.drawparallels(parallels, xoffset=xoffset, labels=[1, 0, 0, 0], **kw)
m.drawmeridians(meridians, yoffset=yoffset, labels=[0, 0, 0, 1], **kw)
return fig, m
In [ ]:
lon = [-40.77, -40.51, -40.30, -40.23, -40.13, -40.06, -39.99,
-39.87, -39.72, -39.52, -39.32, -39.11, -38.91, -38.71]
lat = [-21.29, -21.39, -21.48, -21.51, -21.56, -21.58, -21.62,
-21.69, -21.76, -21.86, -21.96, -22.08, -22.15, -22.25]
fig, m = make_map(img='./data/AVHRR.png')
kw = dict(marker='o', markerfacecolor='k', markeredgecolor='w', markersize=6, linestyle='none')
_ = m.plot(*m(lon, lat), **kw)
In [ ]:
import seaborn
import numpy as np
import matplotlib.pyplot as plt
from io import BytesIO
from pandas import read_csv
kw = dict(na_values='NaN', sep=',', encoding='utf-8',
skipinitialspace=True, index_col=False)
df = read_csv("./data/fish.csv", **kw)
df.head()
In [ ]:
kw = {'axes.edgecolor': '0', 'text.color': '0', 'ytick.color': '0', 'xtick.color': '0',
'ytick.major.size': 5, 'xtick.major.size': 5, 'axes.labelcolor': '0'}
seaborn.set_style("whitegrid", kw)
ax = seaborn.corrplot(df, annot=False, diag_names=False)
In [ ]:
g = df.groupby('Days')
mean_df = g.mean()
g.describe().head(16)
In [ ]:
ax = seaborn.jointplot("Days", "BDE 99 (ng/g)", df, kind="reg")
In [ ]:
ax = seaborn.jointplot("Days", "BDE 47 (ng/g)", df, kind="reg")
In [ ]:
ax = seaborn.residplot("Days", "BDE 47 (ng/g)", df)
In [ ]:
import iris
url = './data/wodStandardLevels.nc'
sal = iris.load_cube(url, 'sea_water_salinity')
print(sal)
In [ ]:
import cartopy.crs as ccrs
import matplotlib.pyplot as plt
from cartopy.feature import LAND, COASTLINE
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER
x = sal.coord('longitude').points
y = sal.coord('latitude').points
fig, ax = plt.subplots(figsize=(11, 13),
subplot_kw=dict(projection=ccrs.PlateCarree()))
ax.set_extent([-180, 180, -90, 90])
ax.stock_img()
ax.plot(x, y, 'go')
ax.add_feature(COASTLINE)
kw = dict(linewidth=1.5, color='gray', alpha=0.5, linestyle='--')
gl = ax.gridlines(crs=ccrs.PlateCarree(), draw_labels=True, **kw)
gl.xlabels_top = gl.ylabels_right = False
gl.xformatter = LONGITUDE_FORMATTER
gl.yformatter = LATITUDE_FORMATTER
from shapely.geometry.polygon import LinearRing
lons = [105, 105, 125, 125]
lats = [-40, -20, -20, -40]
ring = LinearRing(list(zip(lons, lats)))
sq = ax.add_geometries([ring], ccrs.PlateCarree(), facecolor='none', edgecolor='red')
In [ ]:
lon = iris.Constraint(longitude=lambda l:min(lons) < l < max(lons))
lat = iris.Constraint(latitude=lambda l:min(lats) < l < max(lats))
alt = iris.Constraint(altitude=lambda a: a <= 30)
sal_profile = sal.extract(alt & lon & lat)
print(sal_profile)
In [ ]:
import iris.plot as iplt
lon = sal_profile.coord(axis='X').points.squeeze()
lat = sal_profile.coord(axis='Y').points.squeeze()
depth = sal_profile.coord(axis='Z').points.max()
fig, ax = plt.subplots(figsize=(5, 6))
kw = dict(linewidth=2, color='k', marker='o')
iplt.plot(sal_profile, sal_profile.coord('altitude'), **kw)
ax.grid()
ax.set_ylabel('Depth (m)')
ax.set_xlabel('Salinity (unknown)')
t = ax.set_title('lon: %s\nlat: %s\nMax depth = %s' % (lon, lat, depth))
In [ ]:
import numpy as np
from pandas import DataFrame
x = np.array([2.5, 0.5, 2.2, 1.9, 3.1, 2.3, 2.0, 1.0, 1.5, 1.1])
y = np.array([2.4, 0.7, 2.9, 2.2, 3.0, 2.7, 1.6, 1.1, 1.6, 0.9])
Z = np.c_[x - x.mean(), y - y.mean()]
df = DataFrame(Z, columns=['x', 'y'])
df.index.name = 'Medidas'
df.T
In [ ]:
%matplotlib inline
import matplotlib.pyplot as plt
line = dict(linewidth=1, linestyle='--', color='k')
marker = dict(linestyle='none', marker='o', markersize=7, color='blue', alpha=0.5)
fig, ax = plt.subplots(figsize=(5, 5))
ax.plot(x, y, **marker)
ax.axhline(**line)
ax.axvline(**line)
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_title(u'Dados Originais (com a média)')
_ = ax.axis([-1, 4, -1, 4])
In [ ]:
cov = np.cov(Z.T)
cov
In [ ]:
eigenvalues, eigenvectors = np.linalg.eig(cov)
eigenvalues
In [ ]:
eigenvectors
In [ ]:
pc = eigenvectors[1] / eigenvectors[0]
pc
In [ ]:
fig, ax = plt.subplots(figsize=(5, 5))
ax.plot(Z[:, 0], Z[:, 1], **marker)
ax.set_xlim([-2, 2])
ax.set_ylim([-2, 2])
ax.axhline(**line)
ax.axvline(**line)
ax.set_xlabel('x')
ax.set_ylabel('y')
arrowprops = dict(width=0.01, head_width=0.05, alpha = 0.5,
length_includes_head=False)
a1 = ax.arrow(0, 0, eigenvalues[0], pc[0] * eigenvalues[0],
color='k', **arrowprops)
a2 = ax.arrow(0, 0, eigenvalues[1], pc[1] * eigenvalues[1],
color='k', **arrowprops)
a3 = ax.arrow(0, 0, -eigenvalues[0], -pc[0] * eigenvalues[0],
color='r', **arrowprops)
a4 = ax.arrow(0, 0, -eigenvalues[1], -pc[1] * eigenvalues[1],
color='r', **arrowprops)
In [ ]:
from sklearn.decomposition import PCA
pca = PCA(n_components=2, copy=True)
X = pca.fit_transform(Z)
fig, ax = plt.subplots()
ax.plot(X[:, 0], X[:, 1], **marker)
ax.set_xlim([-2, 2])
ax.set_ylim([-2, 2])
ax.axhline(**line)
ax.axvline(**line)
_ = ax.set_title("PCA")
In [ ]:
print(eigenvectors)
print('')
print(pca.components_)
Quanto da variância é explicada por cada componente.
In [ ]:
pca.explained_variance_ratio_
In [ ]:
from scipy.stats.stats import pearsonr
r, p = pearsonr(x-x.mean(), y-y.mean())
print('Pearson r: {}\nPC 1: {}:'.format(r, pc[0]))
In [ ]:
import statsmodels.api as sm
def OLSreg(y, Xmat):
return sm.OLS(y, sm.add_constant(Xmat, prepend=True)).fit()
scale = lambda x: (x - x.mean()) / x.std()
ols_fit = OLSreg(Z[:, 1], Z[:, 0])
fig, ax = plt.subplots(figsize=(5, 5))
ax.plot(Z[:, 0], Z[:, 1], **marker)
ax.set_xlim([-2, 2])
ax.set_ylim([-2, 2])
ax.axhline(**line)
ax.axvline(**line)
ax.set_xlabel('x')
ax.set_xlabel('y')
ax.plot(Z[:, 0], ols_fit.fittedvalues, 'r', alpha=0.5)
ax.text(-1, 1, r'R$^2$ = %4.3f' % round(ols_fit.rsquared, 3))
a1 = ax.arrow(0, 0, eigenvalues[1], pc[1] * eigenvalues[1],
**arrowprops)
a2 = ax.arrow(0, 0, -eigenvalues[1], -pc[1] * eigenvalues[1],
**arrowprops)
In [ ]:
%%file data.csv
Food,England,Wales,Scotland,N Ireland
Cheese,105,103,103,66
Carcass meat,245,227,242,267
Other meat,685,803,750,586
Fish,147,160,122,93
Fats and oils,193,235,184,209
Sugars,156,175,147,139
Fresh potatoes,720,874,566,1033
Fresh Veg,253,265,171,143
Other Veg,488,570,418,355
Processed potatoes,198,203,220,187
Processed Veg,360,365,337,334
Fresh fruit,1102,1137,957,674
Cereals,1472,1582,1462,1494
Beverages,57,73,53,47
Soft drinks,1374,1256,1572,1506
Alcoholic drinks,375,475,458,135
Confectionery,54,64,62,41
In [ ]:
from pandas import read_csv
df = read_csv('data.csv', index_col='Food')
df
In [ ]:
def z_score(x):
"""Remove a média e normaliza os pelo desvio padrão"""
return (x - x.mean()) / x.std()
In [ ]:
from sklearn.decomposition import PCA
pca = PCA(n_components=None)
pca.fit(df.apply(z_score).T)
In [ ]:
from pandas import DataFrame
loadings = DataFrame(pca.components_.T)
loadings.index = ['PC %s' % pc for pc in loadings.index + 1]
loadings.columns = ['TS %s' % pc for pc in loadings.columns + 1]
loadings
In [ ]:
PCs = np.dot(loadings.values.T, df)
marker = dict(linestyle='none', marker='o', markersize=7, color='blue', alpha=0.5)
fig, ax = plt.subplots(figsize=(7, 2.75))
ax.plot(PCs[0], np.zeros_like(PCs[0]),
label="Scores", **marker)
[ax.text(x, y, t) for x, y, t in zip(PCs[0], loadings.values[1, :], df.columns)]
ax.set_xlabel("PC1")
_ = ax.set_ylim(-1, 1)
marker = dict(linestyle='none', marker='o', markersize=7, color='blue', alpha=0.5)
In [ ]:
fig, ax = plt.subplots(figsize=(7, 2.75))
ax.plot(PCs[0], PCs[1], label="Scores", **marker)
ax.set_xlabel("PC1")
ax.set_ylabel("PC2")
text = [ax.text(x, y, t) for x, y, t in
zip(PCs[0], PCs[1]+0.5, df.columns)]
In [ ]:
perc = pca.explained_variance_ratio_ * 100
perc = DataFrame(perc, columns=['Percentage explained ratio'], index=['PC %s' % pc for pc in np.arange(len(perc)) + 1])
ax = perc.plot(kind='bar')
In [ ]:
TS1 = loadings['TS 1']
TS1.index = df.index
ax = TS1.plot(kind='bar')
In [ ]:
marker = dict(linestyle='none', marker='o', markersize=7, color='blue', alpha=0.5)
fig, ax = plt.subplots(figsize=(7, 7))
ax.plot(loadings.icol(0), loadings.icol(1), label="Loadings", **marker)
ax.set_xlabel("non-projected PC1")
ax.set_ylabel("non-projected PC2")
ax.axis([-1, 1, -1, 1])
text = [ax.text(x, y, t) for
x, y, t in zip(loadings.icol(0), loadings.icol(1), df.index)]
In [ ]:
from sympy.interactive import printing
from sympy import Function, symbols, sqrt, S
printing.init_printing()
x, y = symbols("x y")
In [ ]:
from IPython.display import Latex
red = S("(x / 7) ** 2 * sqrt(abs(abs(x) - 3) / (abs(x) - 3)) + (y / 3) ** 2 * sqrt(abs(y + 3 / 7 * sqrt(33)) / (y + 3 / 7 * sqrt(33))) - 1")
red
In [ ]:
orange = S("abs(x / 2) - ((3 * sqrt(33) - 7) / 112) * x ** 2 -"
"3 + sqrt(1 - (abs(abs(x) - 2) - 1) ** 2) - y")
orange
In [ ]:
green = S("9 * sqrt(abs((abs(x) - 1) * (abs(x) - 0.75)) /"
"((1 - abs(x)) * (abs(x) - 0.75))) - 8 * abs(x) - y")
green
In [ ]:
blue = S("3 * abs(x) + 0.75 * sqrt(abs((abs(x) - 0.75) * (abs(x) - 0.5)) /"
"((0.75 - abs(x)) * (abs(x) - 0.5))) - y")
blue
In [ ]:
pink = S("2.25 * sqrt(abs((x - 0.5) * (x + 0.5)) / ((0.5 - x) * (0.5 + x))) - y")
pink
In [ ]:
brown = S("6 * sqrt(10) / 7 + (1.5 - 0.5 * abs(x)) * sqrt(abs(abs(x) - 1) /"
"(abs(x) - 1)) - (6 * sqrt(10) / 14) * sqrt(4 - (abs(x) - 1) ** 2) - y")
brown
In [ ]:
from sympy import lambdify
def lamb(func):
return lambdify((x, y), func, modules='numpy')
red, orange, green, blue, pink, brown = map(lamb, (red, orange, green, blue, pink, brown))
In [ ]:
import numpy as np
spacing = 0.01
x, y = np.meshgrid(np.arange(-7.25, 7.25, spacing),
np.arange(-5, 5, spacing))
red = red(x, y)
orange = orange(x, y)
green = green(x, y)
blue = blue(x, y)
pink = pink(x, y)
brown = brown(x, y)
In [ ]:
%matplotlib inline
import matplotlib.pyplot as plt
colors = dict(red='#FF0000', orange='#FFA500', green='#008000',
blue='#003399', pink='#CC0033', brown='#800000')
equations = dict(red=red, orange=orange, green=green,
blue=blue, pink=pink, brown=brown)
fig, ax = plt.subplots(figsize=(8, 6))
for key, color in colors.items():
ax.contour(x, y, equations[key], [0], colors=color, linewidths=3)
_ = ax.set_xticks([])
_ = ax.set_yticks([])
In [ ]:
import numpy as np
lon = np.loadtxt('./data/copano_bay/lon.txt')
lat = np.loadtxt('./data/copano_bay/lat.txt')
depth = np.loadtxt('./data/copano_bay/depth.txt')
lon.shape, lat.shape, depth.shape
In [ ]:
from oceans import wrap_lon180
lon = wrap_lon180(lon)
In [ ]:
%matplotlib inline
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
from cartopy.io import shapereader
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER
def make_map(projection=ccrs.PlateCarree()):
fig, ax = plt.subplots(figsize=(9, 13),
subplot_kw=dict(projection=projection))
gl = ax.gridlines(draw_labels=True)
gl.xlabels_top = gl.ylabels_right = False
gl.xformatter = LONGITUDE_FORMATTER
gl.yformatter = LATITUDE_FORMATTER
return fig, ax
In [ ]:
import cartopy.io.img_tiles as cimgt
extent = [lon.min(), lon.max(), lat.min(), lat.max()]
request = cimgt.GoogleTiles()
fig, ax = make_map(projection=request.crs)
ax.set_extent(extent)
ax.add_image(request, 10)
In [ ]:
%matplotlib tk
fig, ax = make_map()
ax.set_extent(extent)
cs = ax.pcolormesh(lon, lat, depth, cmap=plt.cm.ocean,
transform=ccrs.PlateCarree(), zorder=1)
poly = plt.ginput()
In [ ]:
from matplotlib.path import Path
polygon = Path(poly)
def in_polygon(polygon, xp, yp, transform=None, radius=0.0):
return polygon.contains_points(np.atleast_2d([xp, yp]).T,
transform=None, radius=0.0)
In [ ]:
x_poly = np.array(poly)[:, 0]
y_poly = np.array(poly)[:, 1]
xp = lon.ravel()
yp = lat.ravel()
inside = in_polygon(polygon, xp, yp)
mask = inside.reshape(lon.shape)
In [ ]:
%matplotlib inline
fig, ax = make_map()
ax.set_extent(extent)
cs = ax.pcolormesh(lon, lat, depth, cmap=plt.cm.ocean,
transform=ccrs.PlateCarree(), zorder=0)
ax.set_boundary(polygon, transform=ccrs.PlateCarree())
In [ ]:
import numpy.ma as ma
masked_depth = ma.masked_array(depth, ~mask)
In [ ]:
fig, ax = make_map()
ax.set_extent(extent)
cs = ax.pcolormesh(lon, lat, masked_depth, cmap=plt.cm.ocean,
transform=ccrs.PlateCarree(), zorder=1)
In [ ]:
%matplotlib inline
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER
def make_map(projection=ccrs.PlateCarree(), figsize=(11, 11)):
fig, ax = plt.subplots(figsize=figsize,
subplot_kw=dict(projection=projection))
gl = ax.gridlines(draw_labels=True)
gl.xlabels_top = gl.ylabels_right = False
gl.xformatter = LONGITUDE_FORMATTER
gl.yformatter = LATITUDE_FORMATTER
return fig, ax
In [ ]:
import geopandas
basins = geopandas.read_file('./data/hidrografia_dupla_cut/')
basins = basins['geometry'].values.tolist()
In [ ]:
from shapely import ops
basin = ops.cascaded_union(basins)
In [ ]:
bra = geopandas.read_file('./data/BRA_adm2/')
salvador = bra[bra['NAME_2'] == 'Salvador'].geometry.values[0]
camacari = bra[bra['NAME_2'] == u'Camaçari'].geometry.values[0]
salvador = salvador.geoms[0] # Throw the islands away.
In [ ]:
projection = ccrs.PlateCarree()
fig, ax = make_map(projection=projection)
extent = [-39.15, -38.10, -13.30, -12.4]
ax.set_extent(extent)
ax.add_geometries(basin, projection, facecolor='0.85', edgecolor='black')
ax.add_geometries([salvador], projection, facecolor='0.65')
ax.add_geometries([camacari], projection, facecolor='0.65')
kw = dict(fontweight='bold', color='white', verticalalignment='top')
lon, lat = camacari.centroid.bounds[0], camacari.centroid.bounds[1]
ax.text(lon, lat, u'Camaçari', horizontalalignment='center', **kw)
lon, lat = salvador.centroid.bounds[0], salvador.centroid.bounds[1]
ax.text(lon, lat, 'Salvador', **kw);
In [ ]:
from ctd import DataFrame
cast = DataFrame.from_cnv('./data/CTD_046.cnv.bz2', compression='bz2')
In [ ]:
lon, lat = cast.longitude.mean(), cast.latitude.mean()
hours = cast['timeS'].max() / 60. / 60.
print("CTD cast info\nlon: {:.2f}\nlat: "
"{:.2f}\nduration: {:.2f} hours".format(lon, lat, hours))
In [ ]:
cast = cast[cast['pumps']]
In [ ]:
downcast, upcast = cast.split()
In [ ]:
downcast[['scan', 'c0S/m', 't090C', 'timeS', 'altM', 'longitude', 'latitude']].tail()
In [ ]:
downcast[['t090C', 't190C', 'c0S/m', 'c1S/m']].describe()
In [ ]:
downcast.columns
In [ ]:
keep = set(['t090C', 'c0S/m', 'sbeox0Mm/Kg', 'dz/dtM'])
downcast = downcast.drop(keep.symmetric_difference(cast.columns), axis=1)
In [ ]:
downcast.head()
In [ ]:
from ctd import movingaverage
cast['dz/dtM'] = movingaverage(cast['dz/dtM'], window_size=48)
cast['sbeox0Mm/Kg'] = movingaverage(cast['sbeox0Mm/Kg'], window_size=48)
In [ ]:
from ctd import lp_filter
kw = dict(sample_rate=24.0, time_constant=0.15)
cast.index = lp_filter(cast.index, **kw)
In [ ]:
downcast = downcast.press_check() # Remove pressure reversals.
downcast = downcast[downcast['dz/dtM'] >= 0.25] # Threshold velocity.
In [ ]:
from ctd import Series
kw = dict(n1=2, n2=20, block=150)
downcast = downcast.apply(Series.despike, **kw)
In [ ]:
downcast = downcast.apply(Series.bindata, **dict(delta=1.))
downcast = downcast.apply(Series.interpolate)
In [ ]:
pmax = max(cast.index)
if pmax >= 500:
window_len = 21
elif pmax >= 100:
window_len = 11
else:
window_len = 5
kw = dict(window_len=window_len, window='hanning')
downcast = downcast.apply(Series.smooth, **kw)
In [ ]:
import gsw
from oceans.sw_extras import gamma_GP_from_SP_pt
p = downcast.index.values.astype(float)
downcast.index.name = "Pressure"
downcast['depth'] = -gsw.z_from_p(p, lat)
downcast['SP'] = gsw.SP_from_C(downcast['c0S/m'].values * 10, downcast['t090C'].values, p)
downcast['SA'] = gsw.SA_from_SP(downcast['SP'].values, p, lon, lat)
downcast['CT'] = gsw.CT_from_t(downcast['SA'].values, downcast['t090C'].values, p)
downcast['sigma0_CT'] = gsw.sigma0_CT_exact(downcast['SA'].values, downcast['CT'].values)
downcast['sound'] = gsw.sound_speed_CT_exact(downcast['SA'].values, downcast['CT'].values, p)
downcast['gamma'] = gamma_GP_from_SP_pt(downcast['SA'].values, downcast['CT'].values, p, lon, lat)
# Does't fit the DataFrame because the index in at the middle (p_mid) point of the original data.
N2, p_mid = gsw.Nsquared(downcast['SA'].values, downcast['CT'].values, p, lat=lat)
N2 = DataFrame(N2, index=p_mid, columns=['N2'])
In [ ]:
from matplotlib import style
style.use('ggplot')
import matplotlib.pyplot as plt
import mpl_toolkits.axisartist as AA
from mpl_toolkits.axes_grid1 import host_subplot
def adjust_xlim(ax, offset):
x1, x2 = ax.get_xlim()
ax.set_xlim(x1 - offset, x2 + offset)
def plot_vars(self, title, **kw):
fig = plt.figure(figsize=(6, 9))
fig.subplots_adjust(bottom=0.1, top=0.85)
if title:
fig.suptitle(title)
ax0 = host_subplot(111, axes_class=AA.Axes)
ax0.invert_yaxis()
ax1, ax2, ax3 = ax0.twiny(), ax0.twiny(), ax0.twiny()
new_axis0 = ax0.get_grid_helper().new_fixed_axis
new_axis1 = ax1.get_grid_helper().new_fixed_axis
new_axis2 = ax2.get_grid_helper().new_fixed_axis
new_axis3 = ax3.get_grid_helper().new_fixed_axis
ax0.axis["bottom"] = new_axis0(loc="bottom", axes=ax0, offset=(0, 0))
ax1.axis["bottom"] = new_axis1(loc="bottom", axes=ax1, offset=(0, -30))
ax2.axis["top"] = new_axis2(loc="top", axes=ax2, offset=(0, 0))
ax3.axis["top"] = new_axis3(loc="top", axes=ax3, offset=(0, 30))
ax0.axis["top"].toggle(all=False)
ax1.axis["top"].toggle(all=False)
ax2.axis["bottom"].toggle(all=False)
ax3.axis["bottom"].toggle(all=False)
ax0.set_ylabel("Pressure [dbar]")
ax0.set_xlabel(r"Oxygen [$\mu$mol kg$^{-1}$]")
ax1.set_xlabel(r"Absolute Salinity [g kg$^{-1}$]")
ax2.set_xlabel(u"Conservative Temperature [\xb0C]")
ax3.set_xlabel(r"$\gamma^{GP}$ [kg m$^{-3}$]")
l0, = ax0.plot(self['sbeox0Mm/Kg'], self.index, 'orchid')
l1, = ax1.plot(self['SA'], self.index, 'darkorange')
l2, = ax2.plot(self['CT'], self.index, 'darkgreen')
l3, = ax3.plot(self['gamma'], self.index, 'crimson')
adjust_xlim(ax0, offset=1)
adjust_xlim(ax1, offset=0.15)
adjust_xlim(ax2, offset=0.05)
adjust_xlim(ax3, offset=0.05)
ax0.axis["bottom"].label.set_color(l0.get_color())
ax1.axis["bottom"].label.set_color(l1.get_color())
ax2.axis["top"].label.set_color(l2.get_color())
ax3.axis["top"].label.set_color(l3.get_color())
return fig, (ax0, ax1, ax2, ax3), (l0, l1, l2, l3)
In [ ]:
fig, axes, lines = plot_vars(downcast, title=None, **kw)
In [ ]:
!cat ./data/nt_20150326_f17_nrt_n.bin.txt
In [ ]:
import numpy as np
import numpy.ma as ma
infile = './data/nt_20150326_f17_nrt_n.bin'
with open(infile, 'rb') as fr:
hdr = fr.read(300)
ice = np.fromfile(fr, dtype=np.uint8)
ice = ice.reshape(448, 304)
ice = ice / 250.
ice = ma.masked_greater(ice, 1.0)
In [ ]:
dx = dy = 25000
x = np.arange(-3850000, +3750000, +dx)
y = np.arange(+5850000, -5350000, -dy)
x.shape, y.shape, ice.shape
In [ ]:
%matplotlib inline
import cartopy.crs as ccrs
import matplotlib.pyplot as plt
fig = plt.figure(figsize=(9, 9))
ax = plt.axes(projection=ccrs.NorthPolarStereo(central_longitude=0))
cs = ax.coastlines(resolution='110m', linewidth=0.5)
ax.gridlines()
ax.set_extent([-180, 180, 40, 90], crs=ccrs.PlateCarree())
kw = dict(central_latitude=90, central_longitude=-45, true_scale_latitude=70)
cs = ax.pcolormesh(x, y, ice, cmap=plt.cm.Blues,
transform=ccrs.Stereographic(**kw))
conda install oct2pyVeja também: https://ocefpaf.github.io/python4oceanographers/blog/2013/11/04/t_tide/
In [ ]:
%load_ext oct2py.ipython
In [ ]:
%%file fahr_to_kelvin.m
function ktemp = fahr_to_kelvin(ftemp)
ktemp = ((ftemp - 32) * (5/9)) + 273.15;
end
In [ ]:
%%octave
fahr_to_kelvin(32)
In [ ]:
import numpy as np
angle = np.linspace(-2*np.pi, 2*np.pi, 100)
In [ ]:
%%octave -i angle
plot(cos(angle))
In [ ]:
%%octave -o angle
angle = linspace(-2*pi, 2*pi, 100);
In [ ]:
%matplotlib inline
import matplotlib.pyplot as plt
plt.plot(np.cos(angle.squeeze()));
Veja também: https://ocefpaf.github.io/python4oceanographers/blog/2015/01/26/fortran/
In [ ]:
%load_ext fortranmagic
In [ ]:
%%fortran
subroutine calc_sin(x, y)
real, intent(in) :: x
real, intent(out) :: y
y = sin(x)
end subroutine calc_sin
In [ ]:
calc_sin(np.pi/2)
In [ ]:
calc_sin?
In [ ]:
%%fortran --link lapack -vv
subroutine solve(A, b, x, n)
! solve the matrix equation A*x=b using LAPACK
implicit none
real*8, dimension(n,n), intent(in) :: A
real*8, dimension(n), intent(in) :: b
real*8, dimension(n), intent(out) :: x
integer :: i, j, pivot(n), ok
integer, intent(in) :: n
x = b
! find the solution using the LAPACK routine SGESV
call DGESV(n, 1, A, n, pivot, x, n, ok)
end subroutine
In [ ]:
A = [[1, 2.5], [-3, 4]]
b = [1, 2.5]
solve(A, b)
In [ ]:
np.linalg.solve(A, b)
In [ ]:
solve?
conda install -c r rconda install -c r r-irkernelObs: options(menu.graphics=FALSE) and install.packages("oce")
In [ ]:
sum(c(1, 2), na.r=F)
In [ ]:
sum(c(1, 2), na.r=T)
In [ ]:
library(oce)
daylength <- function(t, lon=-38.476667, lat=-12.974722)
{
t <- as.numeric(t)
alt <- function(t)
sunAngle(t, longitude=lon, latitude=lat)$altitude
rise <- uniroot(alt, lower=t-86400/2, upper=t)$root
set <- uniroot(alt, lower=t, upper=t+86400/2)$root
set - rise
}
t0 <- as.POSIXct("2015-01-01 12:00:00", tz="UTC")
t <- seq.POSIXt(t0, by="1 day", length.out=1*356)
dayLength <- unlist(lapply(t, daylength))
par(mfrow=c(2,1), mar=c(3, 3, 1, 1), mgp=c(2, 0.7, 0))
plot(t, dayLength/3600, type='o', pch=20,
xlab="", ylab="Day length (hours)")
grid()
In [ ]:
Sa <- 30
Ta <- 10
Sb <- 40
library(oce)
# Should not need to edit below this line
rho0 <- swRho(Sa, Ta, 0)
Tb <- uniroot(function(T) rho0-swRho(Sb,T,0), lower=0, upper=100)$root
Sc <- (Sa + Sb) /2
Tc <- (Ta + Tb) /2
## density change, and equiv temp change
drho <- swRho(Sc, Tc, 0) - rho0
dT <- drho / rho0 / swAlpha(Sc, Tc, 0)
plotTS(as.ctd(c(Sa, Sb, Sc), c(Ta, Tb, Tc), 0), pch=20, cex=2)
drawIsopycnals(levels=rho0, col="red", cex=0)
segments(Sa, Ta, Sb, Tb, col="blue")
text(Sb, Tb, "b", pos=4)
text(Sa, Ta, "a", pos=4)
text(Sc, Tc, "c", pos=4)
legend("topleft",
legend=sprintf("Sa=%.1f, Ta=%.1f, Sb=%.1f -> Tb=%.1f, drho=%.2f, dT=%.2f",
Sa, Ta, Sb, Tb, drho, dT),
bg="white")