This is part of Python for Geosciences notes.
=============
Matplotlib is a python 2D plotting library which produces publication quality figures in a variety of hardcopy formats and interactive environments across platforms.
Let's prepare some data:
x = linspace(0,10,20)
y = x ** 2
Plot is as easy as this:
plot(x,y)
Line style and labels are controlled in a way similar to Matlab:
plot(x, y, 'r--o')
xlabel('x')
ylabel('y')
title('title')
You can plot several individual lines at once:
plot(x, y, 'r--o', x, y ** 1.1, 'bs', x, y ** 1.2, 'g^-' )
One more example:
mu, sigma = 100, 15
x = mu + sigma * np.random.randn(10000)
# the histogram of the data
n, bins, patches = plt.hist(x, 50, normed=1, facecolor='g', alpha=0.75)
xlabel('Smarts')
ylabel('Probability')
title('Histogram of IQ')
text(60, .025, r'$\mu=100,\ \sigma=15$')
axis([40, 160, 0, 0.03])
grid(True)
If you feel a bit playful (only in matplotlib > 1.3):
with xkcd():
x = np.linspace(0, 1)
y = np.sin(4 * np.pi * x) * np.exp(-5 * x)
plt.fill(x, y, 'r')
plt.grid(True)
Following example is from matplotlib - 2D and 3D plotting in Python - great place to start for people interested in matplotlib.
n = array([0,1,2,3,4,5])
xx = np.linspace(-0.75, 1., 100)
x = linspace(0, 5, 10)
fig, axes = plt.subplots(1, 4, figsize=(12,3))
axes[0].scatter(xx, xx + 0.25*randn(len(xx)))
axes[1].step(n, n**2, lw=2)
axes[2].bar(n, n**2, align="center", width=0.5, alpha=0.5)
axes[3].fill_between(x, x**2, x**3, color="green", alpha=0.5);
When you going to plot something more or less complicated in Matplotlib, the first thing you do is open the Matplotlib example gallery and choose example closest to your case.
You can directly load python code (or basically any text file) to the notebook. This time we download code from the Matplotlib example gallery:
%load http://matplotlib.org/mpl_examples/pylab_examples/griddata_demo.py
from numpy.random import uniform, seed
from matplotlib.mlab import griddata
import matplotlib.pyplot as plt
import numpy as np
# make up data.
#npts = int(raw_input('enter # of random points to plot:'))
seed(0)
npts = 200
x = uniform(-2,2,npts)
y = uniform(-2,2,npts)
z = x*np.exp(-x**2-y**2)
# define grid.
xi = np.linspace(-2.1,2.1,100)
yi = np.linspace(-2.1,2.1,200)
# grid the data.
zi = griddata(x,y,z,xi,yi,interp='linear')
# contour the gridded data, plotting dots at the nonuniform data points.
CS = plt.contour(xi,yi,zi,15,linewidths=0.5,colors='k')
CS = plt.contourf(xi,yi,zi,15,cmap=plt.cm.rainbow,
vmax=abs(zi).max(), vmin=-abs(zi).max())
plt.colorbar() # draw colorbar
# plot data points.
plt.scatter(x,y,marker='o',c='b',s=5,zorder=10)
plt.xlim(-2,2)
plt.ylim(-2,2)
plt.title('griddata test (%d points)' % npts)
plt.show()
# test case that scikits.delaunay fails on, but natgrid passes..
#data = np.array([[-1, -1], [-1, 0], [-1, 1],
# [ 0, -1], [ 0, 0], [ 0, 1],
# [ 1, -1 - np.finfo(np.float_).eps], [ 1, 0], [ 1, 1],
# ])
#x = data[:,0]
#y = data[:,1]
#z = x*np.exp(-x**2-y**2)
## define grid.
#xi = np.linspace(-1.1,1.1,100)
#yi = np.linspace(-1.1,1.1,100)
## grid the data.
#zi = griddata(x,y,z,xi,yi)
Maps¶
In order to create a map we have to first import some data. We are going to use NCEP reanalysis file from previous section:
from scipy.io import netcdf
f =netcdf.netcdf_file('air.sig995.2012.nc')
Here we create netCDF variable objec for air (we would like to have acces to some of the attributes), but from lat and lon we import only data valies:
air = f.variables['air']
lat = f.variables['lat'][:]
lon = f.variables['lon'][:]
Easiest way to look at the array is imshow:
imshow(air[0,:,:])
colorbar()
But we want some real map :) First unpack and convert data from air:
air_c = ((air[:] * air.scale_factor) + air.add_offset) - 273.15
Our coordinate variables are vectors:
lat.shape
For the map we need 2d coordinate arrays. Convert lot lan to 2d:
lon2, lat2 = np.meshgrid(lon,lat)
Import Basemap - library for plotting 2D data on maps:
from mpl_toolkits.basemap import Basemap
Create Basemap instance (with certain characteristics) and convert lon lat to map coordinates
m = Basemap(projection='npstere',boundinglat=60,lon_0=0,resolution='l')
x, y = m(lon2, lat2)
Creating the map now is only two lines:
m.drawcoastlines()
m.contourf(x,y,air_c[0,:,:])
We can make the map look prettier by adding couple of lines:
fig = plt.figure(figsize=(15,7))
m.fillcontinents(color='gray',lake_color='gray')
m.drawcoastlines()
m.drawparallels(np.arange(-80.,81.,20.))
m.drawmeridians(np.arange(-180.,181.,20.))
m.drawmapboundary(fill_color='white')
m.contourf(x,y,air_c[0,:,:],40)
plt.title('Monthly mean SAT')
colorbar()
You can change map characteristics by changin the Basemap instance:
m = Basemap(projection='ortho',lat_0=45,lon_0=-100,resolution='l')
x, y = m(lon2, lat2)
While the rest of the code might be the same:
fig = plt.figure(figsize=(15,7))
#m.fillcontinents(color='gray',lake_color='gray')
m.drawcoastlines()
m.drawparallels(np.arange(-80.,81.,20.))
m.drawmeridians(np.arange(-180.,181.,20.))
m.drawmapboundary(fill_color='white')
cs = m.contourf(x,y,air_c[0,:,:],20)
plt.title('Monthly mean SAT')
One more map exampe:
m = Basemap(projection='cyl',llcrnrlat=-90,urcrnrlat=90,\
llcrnrlon=0,urcrnrlon=360,resolution='c')
x, y = m(lon2, lat2)
fig = plt.figure(figsize=(15,7))
#m.fillcontinents(color='gray',lake_color='gray')
m.drawcoastlines()
m.drawparallels(np.arange(-80.,81.,20.))
m.drawmeridians(np.arange(0.,360.,20.))
m.drawmapboundary(fill_color='white')
cs = m.contourf(x,y,air[0,:,:],20)
plt.title('Monthly mean SAT')
Comments !