Differences between revisions 12 and 13

Deletions are marked like this. Additions are marked like this.
Line 17: Line 17:
map = Basemap(projection='ortho',lat_0=50,lon_0=-100,
              resolution='l',area_thresh=1000.)
map = Basemap(projection='ortho', lat_0 = 50, lon_0 = -100,
              resolution = 'l', area_thresh = 1000.)
Line 22: Line 22:
map.fillcontinents(color='coral') map.fillcontinents(color = 'coral')
Line 26: Line 26:
map.drawmeridians(np.arange(0,360,30))
map.drawparallels(np.arange(-90,90,30))
map.drawmeridians(np.arange(0, 360, 30))
map.drawparallels(np.arange(-90, 90, 30))
Line 34: Line 34:
Now, suppose you would like to plot the locations of five cities on this map. Add the following just before the {{{p.show()}}} in the above script: Now, suppose you would like to plot the locations of five cities on this map. Add the following just before the {{{plt.show()}}} in the above script:
Line 38: Line 38:
lats=[40.02,32.73,38.55,48.25,17.29]
lons=[-105.16,-117.16,-77.00,-114.21,-88.10]
lats = [40.02, 32.73, 38.55, 48.25, 17.29]
lons = [-105.16, -117.16, -77.00, -114.21, -88.10]

Plotting data on map projections is easy with the basemap toolkit. Toolkits are collections of application-specific functions that extend matplotlib.

The basemap toolkit is not in the default matplotlib install - you can download it from the matplotlib sourceforge download page.

Suppose you'd like to make a map of the world using an orthographic, or satellite projection and plot some data on it. Here's how to make the map (using matplotlib >= 0.98.0 and basemap >= 0.99):

   1 from mpl_toolkits.basemap import Basemap
   2 import matplotlib.pyplot as plt
   3 import numpy as np
   4 # set up orthographic map projection with
   5 # perspective of satellite looking down at 50N, 100W.
   6 # use low resolution coastlines.
   7 # don't plot features that are smaller than 1000 square km.
   8 map = Basemap(projection='ortho', lat_0 = 50, lon_0 = -100,
   9               resolution = 'l', area_thresh = 1000.)
  10 # draw coastlines, country boundaries, fill continents.
  11 map.drawcoastlines()
  12 map.drawcountries()
  13 map.fillcontinents(color = 'coral')
  14 # draw the edge of the map projection region (the projection limb)
  15 map.drawmapboundary()
  16 # draw lat/lon grid lines every 30 degrees.
  17 map.drawmeridians(np.arange(0, 360, 30))
  18 map.drawparallels(np.arange(-90, 90, 30))
  19 plt.show()

There are many other map projections available, probably more than you've even heard of before. A complete list is available in the basemap docstrings. Coastlines, political boundaries and rivers are available in four resolutions, crude, low, intermediate and high. Here's what the low resolution coastlines look like.

Now, suppose you would like to plot the locations of five cities on this map. Add the following just before the plt.show() in the above script:

   1 # lat/lon coordinates of five cities.
   2 lats = [40.02, 32.73, 38.55, 48.25, 17.29]
   3 lons = [-105.16, -117.16, -77.00, -114.21, -88.10]
   4 cities=['Boulder, CO','San Diego, CA',
   5         'Washington, DC','Whitefish, MT','Belize City, Belize']
   6 # compute the native map projection coordinates for cities.
   7 x,y = map(lons,lats)
   8 # plot filled circles at the locations of the cities.
   9 map.plot(x,y,'bo')
  10 # plot the names of those five cities.
  11 for name,xpt,ypt in zip(cities,x,y):
  12     plt.text(xpt+50000,ypt+50000,name)

Calling a basemap class instance with arrays of longitudes and latitudes returns those locations in native map projection coordinates using the proj4 library. Now suppose you have some data on a regular latitude/longitude grid and you would like to plot contours of that data over the map. Try adding the following lines just before p.show()

   1 # make up some data on a regular lat/lon grid.
   2 nlats = 73; nlons = 145; delta = 2.*np.pi/(nlons-1)
   3 lats = (0.5*np.pi-delta*np.indices((nlats,nlons))[0,:,:])
   4 lons = (delta*np.indices((nlats,nlons))[1,:,:])
   5 wave = 0.75*(np.sin(2.*lats)**8*np.cos(4.*lons))
   6 mean = 0.5*np.cos(2.*lats)*((np.sin(2.*lats))**2 + 2.)
   7 # compute native map projection coordinates of lat/lon grid.
   8 x, y = map(lons*180./np.pi, lats*180./np.pi)
   9 # contour data over the map.
  10 CS = map.contour(x,y,wave+mean,15,linewidths=1.5)

Instead of plotting the continents and coastlines, you can use an image as a map background using the warpimage method. The default background image is the NASA 'blue marble' image, which you can apply by using

   1 map.bluemarble()

in place of

   1 map.drawcoastlines()
   2 map.drawcountries()
   3 map.fillcontinents(color='coral')

Here's what the resulting plot looks like (using white text instead of black, and white dots instead of blue)

You can also plot images, pcolor plots and vectors over map projections. Examples that illustrate this and more can be found in the examples directory of the basemap source distribution.


CategoryCookbookMatplotlib

Cookbook/Matplotlib/Maps (last edited 2010-12-28 07:29:08 by JoonRo)