Matplotlib Basemap Tutorial: Making a simple map

Update:

This tutorial has been updated and moved to http://introtopython.org/visualization_earthquakes.html.

Now that you have matplotlib and basemap installed, let’s try to make a simple map of the world.  This tutorial is adapted from the Scipy Cookbook.  Open your favorite text editor, and enter the following:

from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
import numpy as np

# make sure the value of resolution is a lowercase L,
#  for 'low', not a numeral 1
map = Basemap(projection='ortho', lat_0=50, lon_0=-100,
              resolution='l', area_thresh=1000.0)

map.drawcoastlines()

plt.show()

If you run this program you should see a plot of the world, with the coastlines mostly filled in:

globe coastlines

A world map, with coastlines shown.

Let’s add some more detail to this map, starting with country borders.  Add the following line after map.drawcoastlines():

from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
import numpy as np

# make sure the value of resolution is a lowercase L,
#  for 'low', not a numeral 1
map = Basemap(projection='ortho', lat_0=50, lon_0=-100,
              resolution='l', area_thresh=1000.0)

map.drawcoastlines()
map.drawcountries()

plt.show()

Run the program again, and you will see countries outlined:

globe with countries

A world map, with country outlines drawn.

Now let’s fill in the continents:

from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
import numpy as np

# make sure the value of resolution is a lowercase L,
#  for 'low', not a numeral 1
map = Basemap(projection='ortho', lat_0=50, lon_0=-100,
              resolution='l', area_thresh=1000.0)

map.drawcoastlines()
map.drawcountries()
map.fillcontinents(color='coral')

plt.show()

You should see the continents filled in.  If the entire globe is filled in, you may need to update your version of matplotlib, the basemap toolkit, or both.  That said, I watched one student follow this tutorial.  The entire globe ended up colored in, but on subsequent runs the colors came out right.  I’m pretty sure this is a glitch in either matplotlib or Basemap.  If you are having issues with this, you may want to go back and install basemap from source.

globe with continents

A world map, with the continents filled in.

world map with color error

If your entire map ends up colored, you may need to install a more recent version of matplotlib, the Baseemap toolkit, or both.

Let’s clean up the edge of the globe:

from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
import numpy as np

# make sure the value of resolution is a lowercase L,
#  for 'low', not a numeral 1
map = Basemap(projection='ortho', lat_0=50, lon_0=-100,
              resolution='l', area_thresh=1000.0)

map.drawcoastlines()
map.drawcountries()
map.fillcontinents(color='coral')
map.drawmapboundary()

plt.show()

You should see a cleaner circle outlining the globe.

globe with clean boundary

A world map, with a cleaner edge.

Now let’s draw latitude and longitude lines:

from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
import numpy as np

# make sure the value of resolution is a lowercase L,
#  for 'low', not a numeral 1
map = Basemap(projection='ortho', lat_0=50, lon_0=-100,
              resolution='l', area_thresh=1000.0)

map.drawcoastlines()
map.drawcountries()
map.fillcontinents(color='coral')
map.drawmapboundary()

map.drawmeridians(np.arange(0, 360, 30))
map.drawparallels(np.arange(-90, 90, 30))

plt.show()

The np.arange arguments tell where your latitude and longitude lines should begin and end, and how far apart they should be spaced.

globe with lat and lon lines

A world map, complete with latitude and longitude lines.

Let’s play with two of the map settings, and then we will move on to plotting data on this globe.  Let’s start by adjusting the perspective.  Change the latitude and longitude parameters in the original Basemap definition to 0 and 100:

from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
import numpy as np

# make sure the value of resolution is a lowercase L,
#  for 'low', not a numeral 1
map = Basemap(projection='ortho', lat_0=0, lon_0=-100,
              resolution='l', area_thresh=1000.0)

map.drawcoastlines()
map.drawcountries()
map.fillcontinents(color='coral')
map.drawmapboundary()

map.drawmeridians(np.arange(0, 360, 30))
map.drawparallels(np.arange(-90, 90, 30))

plt.show()

Run the program, and you should see our map centered along the equator:

globe centered on equator

A world map, centered on the equator.

Now let’s change the kind of map we are producing.  Change the projection type to ‘robin’:

from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
import numpy as np

# make sure the value of resolution is a lowercase L,
#  for 'low', not a numeral 1
map = Basemap(projection='robin', lat_0=0, lon_0=-100,
              resolution='l', area_thresh=1000.0)

map.drawcoastlines()
map.drawcountries()
map.fillcontinents(color='coral')
map.drawmapboundary()

map.drawmeridians(np.arange(0, 360, 30))
map.drawparallels(np.arange(-90, 90, 30))

plt.show()

Now we have a Robinson projection instead of a globe:

Robinson projection map

A Robinson projection of the world.

Before we move on to plotting points on the map, let’s see how to zoom in on a region.  This is good to know because there are many data sets specific to one region of the world, which would get lost when plotted on a map of the whole world.  Some projections can not be zoomed in at all, so if things are not working well, make sure to look at the documentation.

I live on Baranof Island in southeast Alaska, so let’s zoom in on that region.  One way to zoom in is to specify the latitude and longitude of the lower left and upper right corners of the region you want to show.  Let’s use a mercator projection, which supports this method of zooming.  The notation for “lower left corner at 136.25 degrees west and 56 degrees north” is:

llcrnrlon = -136.25, llcrnrlat = 56.0

So, the full set of parameters we will try is:

from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
import numpy as np

# make sure the value of resolution is a lowercase L,
#  for 'low', not a numeral 1
map = Basemap(projection='merc', lat_0=57, lon_0=-135,
    resolution = 'l', area_thresh = 1000.0,
    llcrnrlon=-136.25, llcrnrlat=56,
    urcrnrlon=-134.25, urcrnrlat=57.75)

map.drawcoastlines()
map.drawcountries()
map.fillcontinents(color='coral')
map.drawmapboundary()

map.drawmeridians(np.arange(0, 360, 30))
map.drawparallels(np.arange(-90, 90, 30))

plt.show()

Note that the center of the map, given by lat_0 and lon_0, must be within the region you are zoomed in on.

This works, but the map is pretty ugly:

A map of Baranof Island, at low resolution.

We are missing an entire island to the west of here!  Let’s change the resolution to ‘h’ for ‘high’, and see what we get:

from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
import numpy as np

# make sure the value of resolution is a lowercase L,
#  for 'low', not a numeral 1
map = Basemap(projection='merc', lat_0=57, lon_0=-135,
    resolution = 'h', area_thresh = 1000.0,
    llcrnrlon=-136.25, llcrnrlat=56,
    urcrnrlon=-134.25, urcrnrlat=57.75)

map.drawcoastlines()
map.drawcountries()
map.fillcontinents(color='coral')
map.drawmapboundary()

map.drawmeridians(np.arange(0, 360, 30))
map.drawparallels(np.arange(-90, 90, 30))

plt.show()

A map of Baranof Island at high resolution, but still with a low level of detail.

This is much better, but we are still missing an entire island to the west.  This is because of the “area_thresh” setting.  This setting specifies how large a feature must be in order to appear on the map.  The current setting will only show features larger than 1000 square kilometers.  This is a reasonable setting for a low-resolution map of the world, but it is a really bad choice for a small-scale map.  Let’s change that setting to 0.1, and see how much detail we get:

from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt
import numpy as np

# make sure the value of resolution is a lowercase L,
#  for 'low', not a numeral 1
map = Basemap(projection='merc', lat_0=57, lon_0=-135,
    resolution = 'h', area_thresh = 0.1,
    llcrnrlon=-136.25, llcrnrlat=56,
    urcrnrlon=-134.25, urcrnrlat=57.75)

map.drawcoastlines()
map.drawcountries()
map.fillcontinents(color='coral')
map.drawmapboundary()

map.drawmeridians(np.arange(0, 360, 30))
map.drawparallels(np.arange(-90, 90, 30))

plt.show()

A map of Baranof Island, at high resolution, showing details larger than 0.1 square kilometers.

This is now a meaningful map.  We can see Kruzof island, the large island to the west of Baranof, and many other islands in the area.  Settings lower than area_thresh=0.1 do not add any new details at this level of zoom.

Basemap is an incredibly flexible package.  If you are curious to play around with other settings, take a look at the Basemap documentation.  Next, we will learn how to plot points on our maps.

Next:  Plotting points on a simple map

Back to Tutorial contents

Advertisements

About ehmatthes

Teacher, hacker, new dad, outdoor guy
This entry was posted in programming and tagged , , , . Bookmark the permalink.

One Response to Matplotlib Basemap Tutorial: Making a simple map

  1. Hugh says:

    Hey, This tutorial (I mean the ilearnpython book) is great… it should be linked to from the documentation of the module itself.

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s