matplotlib basemap toolkit

mpl_toolkits.basemap

class mpl_toolkits.basemap.Basemap(llcrnrlon=None, llcrnrlat=None, urcrnrlon=None, urcrnrlat=None, llcrnrx=None, llcrnry=None, urcrnrx=None, urcrnry=None, width=None, height=None, projection='cyl', resolution='c', area_thresh=None, rsphere=6370997.0, ellps=None, lat_ts=None, lat_1=None, lat_2=None, lat_0=None, lon_0=None, lon_1=None, lon_2=None, o_lon_p=None, o_lat_p=None, k_0=None, no_rot=False, suppress_ticks=True, satellite_height=35786000, boundinglat=None, fix_aspect=True, anchor='C', celestial=False, round=False, epsg=None, ax=None)

Sets up a basemap with specified map projection. and creates the coastline data structures in map projection coordinates.

Calling a Basemap class instance with the arguments lon, lat will convert lon/lat (in degrees) to x/y map projection coordinates (in meters). The inverse transformation is done if the optional keyword inverse is set to True.

The desired projection is set with the projection keyword. Default is cyl. Supported values for the projection keyword are:

Value

Description

cyl

Cylindrical Equidistant

merc

Mercator

tmerc

Transverse Mercator

omerc

Oblique Mercator

mill

Miller Cylindrical

gall

Gall Stereographic Cylindrical

cea

Cylindrical Equal Area

lcc

Lambert Conformal

laea

Lambert Azimuthal Equal Area

nplaea

North-Polar Lambert Azimuthal

splaea

South-Polar Lambert Azimuthal

eqdc

Equidistant Conic

aeqd

Azimuthal Equidistant

npaeqd

North-Polar Azimuthal Equidistant

spaeqd

South-Polar Azimuthal Equidistant

aea

Albers Equal Area

stere

Stereographic

npstere

North-Polar Stereographic

spstere

South-Polar Stereographic

cass

Cassini-Soldner

poly

Polyconic

ortho

Orthographic

geos

Geostationary

nsper

Near-Sided Perspective

sinu

Sinusoidal

moll

Mollweide

hammer

Hammer

robin

Robinson

kav7

Kavrayskiy VII

eck4

Eckert IV

vandg

van der Grinten

mbtfpq

McBryde-Thomas Flat-Polar Quartic

gnom

Gnomonic

rotpole

Rotated Pole

For most map projections, the map projection region can either be specified by setting these keywords:

Keyword

Description

llcrnrlon

longitude of lower left hand corner of the desired map domain (degrees).

llcrnrlat

latitude of lower left hand corner of the desired map domain (degrees).

urcrnrlon

longitude of upper right hand corner of the desired map domain (degrees).

urcrnrlat

latitude of upper right hand corner of the desired map domain (degrees).

or these

Keyword

Description

width

width of desired map domain in projection coordinates (meters).

height

height of desired map domain in projection coordinates (meters).

lon_0

center of desired map domain (in degrees).

lat_0

center of desired map domain (in degrees).

For sinu, moll, hammer, npstere, spstere, nplaea, splaea, npaeqd, spaeqd, robin, eck4, kav7, or mbtfpq, the values of llcrnrlon, llcrnrlat, urcrnrlon, urcrnrlat, width and height are ignored (because either they are computed internally, or entire globe is always plotted).

For the cylindrical projections (cyl, merc, mill, cea and gall), the default is to use llcrnrlon=-180,llcrnrlat=-90, urcrnrlon=180 and urcrnrlat=90). For all other projections except ortho, geos and nsper, either the lat/lon values of the corners or width and height must be specified by the user.

For ortho, geos and nsper, the lat/lon values of the corners may be specified, or the x/y values of the corners (llcrnrx,llcrnry,urcrnrx,urcrnry) in the coordinate system of the global projection (with x=0,y=0 at the center of the global projection). If the corners are not specified, the entire globe is plotted.

For rotpole, the lat/lon values of the corners on the unrotated sphere may be provided as llcrnrlon,llcrnrlat,urcrnrlon,urcrnrlat, or the lat/lon values of the corners on the rotated sphere can be given as llcrnrx,llcrnry,urcrnrx,urcrnry.

Other keyword arguments:

Keyword

Description

resolution

resolution of boundary database to use. Can be c (crude), l (low), i (intermediate), h (high), f (full) or None. If None, no boundary data will be read in (and class methods such as drawcoastlines will raise an if invoked). Resolution drops off by roughly 80% between datasets. Higher res datasets are much slower to draw. Default c. Coastline data is from the GSHHS (http://www.soest.hawaii.edu/wessel/gshhs/gshhs.html). State, country and river datasets from the Generic Mapping Tools (http://gmt.soest.hawaii.edu).

area_thresh

coastline or lake with an area smaller than area_thresh in km^2 will not be plotted. Default 10000,1000,100,10,1 for resolution c, l, i, h, f.

rsphere

radius of the sphere used to define map projection (default 6370997 meters, close to the arithmetic mean radius of the earth). If given as a sequence, the first two elements are interpreted as the radii of the major and minor axes of an ellipsoid. Note: sometimes an ellipsoid is specified by the major axis and an inverse flattening parameter (if). The minor axis (b) can be computed from the major axis (a) and the inverse flattening parameter using the formula if = a/(a-b).

ellps

string describing ellipsoid (‘GRS80’ or ‘WGS84’, for example). If both rsphere and ellps are given, rsphere is ignored. Default None. See pyproj.pj_ellps for allowed values.

suppress_ticks

suppress automatic drawing of axis ticks and labels in map projection coordinates. Default True, so parallels and meridians can be labelled instead. If parallel or meridian labelling is requested (using drawparallels and drawmeridians methods), automatic tick labelling will be supressed even if suppress_ticks=False. suppress_ticks=False is useful if you want to use your own custom tick formatter, or if you want to let matplotlib label the axes in meters using map projection coordinates.

fix_aspect

fix aspect ratio of plot to match aspect ratio of map projection region (default True).

anchor

determines how map is placed in axes rectangle (passed to axes.set_aspect). Default is C, which means map is centered. Allowed values are C, SW, S, SE, E, NE, N, NW, and W.

celestial

use astronomical conventions for longitude (i.e. negative longitudes to the east of 0). Default False. Implies resolution=None.

ax

set default axes instance (default None - matplotlib.pyplot.gca() may be used to get the current axes instance). If you do not want matplotlib.pyplot to be imported, you can either set this to a pre-defined axes instance, or use the ax keyword in each Basemap method call that does drawing. In the first case, all Basemap method calls will draw to the same axes instance. In the second case, you can draw to different axes with the same Basemap instance. You can also use the ax keyword in individual method calls to selectively override the default axes instance.

The following keywords are map projection parameters which all default to None. Not all parameters are used by all projections, some are ignored. The module variable projection_params is a dictionary which lists which parameters apply to which projections.

Keyword

Description

lat_ts

latitude of true scale. Optional for stereographic, cylindrical equal area and mercator projections. default is lat_0 for stereographic projection. default is 0 for mercator and cylindrical equal area projections.

lat_1

first standard parallel for lambert conformal, albers equal area and equidistant conic. Latitude of one of the two points on the projection centerline for oblique mercator. If lat_1 is not given, but lat_0 is, lat_1 is set to lat_0 for lambert conformal, albers equal area and equidistant conic.

lat_2

second standard parallel for lambert conformal, albers equal area and equidistant conic. Latitude of one of the two points on the projection centerline for oblique mercator. If lat_2 is not given it is set to lat_1 for lambert conformal, albers equal area and equidistant conic.

lon_1

Longitude of one of the two points on the projection centerline for oblique mercator.

lon_2

Longitude of one of the two points on the projection centerline for oblique mercator.

k_0

Scale factor at natural origin (used by ‘tmerc’, ‘omerc’, ‘stere’ and ‘lcc’).

no_rot

only used by oblique mercator. If set to True, the map projection coordinates will not be rotated to true North. Default is False (projection coordinates are automatically rotated).

lat_0

central latitude (y-axis origin) - used by all projections.

lon_0

central meridian (x-axis origin) - used by all projections.

o_lat_p

latitude of rotated pole (only used by ‘rotpole’)

o_lon_p

longitude of rotated pole (only used by ‘rotpole’)

boundinglat

bounding latitude for pole-centered projections (npstere,spstere,nplaea,splaea,npaeqd,spaeqd). These projections are square regions centered on the north or south pole. The longitude lon_0 is at 6-o’clock, and the latitude circle boundinglat is tangent to the edge of the map at lon_0.

round

cut off pole-centered projection at boundinglat (so plot is a circle instead of a square). Only relevant for npstere,spstere,nplaea,splaea,npaeqd or spaeqd projections. Default False.

satellite_height

height of satellite (in m) above equator - only relevant for geostationary and near-sided perspective (geos or nsper) projections. Default 35,786 km.

Useful instance variables:

Variable Name

Description

projection

map projection. Print the module variable supported_projections to see a list of allowed values.

epsg

EPSG code defining projection (see http://spatialreference.org for a list of EPSG codes and their definitions).

aspect

map aspect ratio (size of y dimension / size of x dimension).

llcrnrlon

longitude of lower left hand corner of the selected map domain.

llcrnrlat

latitude of lower left hand corner of the selected map domain.

urcrnrlon

longitude of upper right hand corner of the selected map domain.

urcrnrlat

latitude of upper right hand corner of the selected map domain.

llcrnrx

x value of lower left hand corner of the selected map domain in map projection coordinates.

llcrnry

y value of lower left hand corner of the selected map domain in map projection coordinates.

urcrnrx

x value of upper right hand corner of the selected map domain in map projection coordinates.

urcrnry

y value of upper right hand corner of the selected map domain in map projection coordinates.

rmajor

equatorial radius of ellipsoid used (in meters).

rminor

polar radius of ellipsoid used (in meters).

resolution

resolution of boundary dataset being used (c for crude, l for low, etc.). If None, no boundary dataset is associated with the Basemap instance.

proj4string

the string describing the map projection that is used by PROJ.4.

Converting from Geographic (lon/lat) to Map Projection (x/y) Coordinates

Calling a Basemap class instance with the arguments lon, lat will convert lon/lat (in degrees) to x/y map projection coordinates (in meters). If optional keyword inverse is True (default is False), the inverse transformation from x/y to lon/lat is performed.

For cylindrical equidistant projection (cyl), this does nothing (i.e. x,y == lon,lat).

For non-cylindrical projections, the inverse transformation always returns longitudes between -180 and 180 degrees. For cylindrical projections (self.projection == cyl, mill, cea, gall or merc) the inverse transformation will return longitudes between self.llcrnrlon and self.llcrnrlat.

Input arguments lon, lat can be either scalar floats, sequences or numpy arrays.

Example Usage:

>>> from mpl_toolkits.basemap import Basemap
>>> import numpy as np
>>> import matplotlib.pyplot as plt
>>> # read in topo data (on a regular lat/lon grid)
>>> etopo = np.loadtxt('etopo20data.gz')
>>> lons  = np.loadtxt('etopo20lons.gz')
>>> lats  = np.loadtxt('etopo20lats.gz')
>>> # create Basemap instance for Robinson projection.
>>> m = Basemap(projection='robin',lon_0=0.5*(lons[0]+lons[-1]))
>>> # compute map projection coordinates for lat/lon grid.
>>> x, y = m(*np.meshgrid(lons,lats))
>>> # make filled contour plot.
>>> cs = m.contourf(x,y,etopo,30,cmap=plt.cm.jet)
>>> m.drawcoastlines() # draw coastlines
>>> m.drawmapboundary() # draw a line around the map region
>>> m.drawparallels(np.arange(-90.,120.,30.),labels=[1,0,0,0]) # draw parallels
>>> m.drawmeridians(np.arange(0.,420.,60.),labels=[0,0,0,1]) # draw meridians
>>> plt.title('Robinson Projection') # add a title
>>> plt.show()

[this example (simpletest.py) plus many others can be found in the examples directory of source distribution. The “OO” version of this example (which does not use matplotlib.pyplot) is called “simpletest_oo.py”.]

arcgisimage(server='http://server.arcgisonline.com/ArcGIS', service='ESRI_Imagery_World_2D', xpixels=400, ypixels=None, dpi=96, verbose=False, **kwargs)

Retrieve an image using the ArcGIS Server REST API and display it on the map. In order to use this method, the Basemap instance must be created using the epsg keyword to define the map projection, unless the cyl projection is used (in which case the epsg code 4326 is assumed).

Keywords

Description

server

web map server URL (default http://server.arcgisonline.com/ArcGIS).

service

service (image type) hosted on server (default ESRI_Imagery_World_2D, which is NASA ‘Blue Marble’ image).

xpixels

requested number of image pixels in x-direction (default 400).

ypixels

requested number of image pixels in y-direction. Default (None) is to infer the number from from xpixels and the aspect ratio of the map projection region.

dpi

The device resolution of the exported image (dots per inch, default 96).

verbose

if True, print URL used to retrieve image (default False).

Extra keyword ax can be used to override the default axis instance.

returns a matplotlib.image.AxesImage instance.

barbs(x, y, u, v, *args, **kwargs)

Make a wind barb plot (u, v) with on the map. (see matplotlib.pyplot.barbs documentation).

If latlon keyword is set to True, x,y are intrepreted as longitude and latitude in degrees. Data and longitudes are automatically shifted to match map projection region for cylindrical and pseudocylindrical projections, and x,y are transformed to map projection coordinates. If latlon is False (default), x and y are assumed to be map projection coordinates.

Extra keyword ax can be used to override the default axis instance.

Other *args and **kwargs passed on to matplotlib.pyplot.barbs

Returns two matplotlib.axes.Barbs instances, one for the Northern Hemisphere and one for the Southern Hemisphere.

bluemarble(ax=None, scale=None, **kwargs)

display blue marble image (from http://visibleearth.nasa.gov) as map background. Default image size is 5400x2700, which can be quite slow and use quite a bit of memory. The scale keyword can be used to downsample the image (scale=0.5 downsamples to 2700x1350).

**kwargs passed on to imshow().

returns a matplotlib.image.AxesImage instance.

colorbar(mappable=None, location='right', size='5%', pad='2%', fig=None, ax=None, **kwargs)

Add colorbar to axes associated with a map. The colorbar axes instance is created using the axes_grid toolkit.

Keywords

Description

mappable

the Image, ContourSet, etc. to which the colorbar applies. Default None, matplotlib.pyplot.gci() is used to retrieve the current image mappable.

location

where to put colorbar (‘top’,’bottom’,’left’,’right’) Default ‘right’.

size

width of colorbar axes (string ‘N%’, where N is an integer describing the fractional width of the parent axes). Default ‘5%’.

pad

Padding between parent axes and colorbar axes in same units as size. Default ‘2%’.

fig

Figure instance the map axes instance is associated with. Default None, and matplotlib.pyplot.gcf() is used to retrieve the current active figure instance.

ax

The axes instance which the colorbar will be associated with. Default None, searches for self.ax, and if None uses matplotlib.pyplot.gca().

**kwargs

extra keyword arguments passed on to colorbar method of the figure instance.

Returns a matplotlib colorbar instance.

contour(x, y, data, *args, **kwargs)

Make a contour plot over the map (see matplotlib.pyplot.contour documentation).

If latlon keyword is set to True, x,y are intrepreted as longitude and latitude in degrees. Data and longitudes are automatically shifted to match map projection region for cylindrical and pseudocylindrical projections, and x,y are transformed to map projection coordinates. If latlon is False (default), x and y are assumed to be map projection coordinates.

Extra keyword ax can be used to override the default axis instance.

If tri is set to True, an unstructured grid is assumed (x,y,data must be 1-d) and matplotlib.pyplot.tricontour is used.

Other *args and **kwargs passed on to matplotlib.pyplot.contour (or tricontour if tri=True).

contourf(x, y, data, *args, **kwargs)

Make a filled contour plot over the map (see matplotlib.pyplot.contourf documentation).

If latlon keyword is set to True, x,y are intrepreted as longitude and latitude in degrees. Data and longitudes are automatically shifted to match map projection region for cylindrical and pseudocylindrical projections, and x,y are transformed to map projection coordinates. If latlon is False (default), x and y are assumed to be map projection coordinates.

If x or y are outside projection limb (i.e. they have values > 1.e20), the corresponing data elements will be masked.

Extra keyword ‘ax’ can be used to override the default axis instance.

If tri is set to True, an unstructured grid is assumed (x,y,data must be 1-d) and matplotlib.pyplot.tricontourf is used.

Other *args and **kwargs passed on to matplotlib.pyplot.contourf (or tricontourf if tri=True).

drawcoastlines(linewidth=1.0, linestyle='solid', color='k', antialiased=1, ax=None, zorder=None)

Draw coastlines.

Keyword

Description

linewidth

coastline width (default 1.)

linestyle

coastline linestyle (default solid)

color

coastline color (default black)

antialiased

antialiasing switch for coastlines (default True).

ax

axes instance (overrides default axes instance)

zorder

sets the zorder for the coastlines (if not specified, uses default zorder for matplotlib.patches.LineCollections).

returns a matplotlib.patches.LineCollection object.

drawcounties(linewidth=0.1, linestyle='solid', color='k', antialiased=1, facecolor='none', ax=None, zorder=None, drawbounds=False)

Draw county boundaries in US. The county boundary shapefile originates with the NOAA Coastal Geospatial Data Project (http://coastalgeospatial.noaa.gov/data_gis.html).

Keyword

Description

linewidth

county boundary line width (default 0.1)

linestyle

coastline linestyle (default solid)

color

county boundary line color (default black)

facecolor

fill color of county (default is no fill)

antialiased

antialiasing switch for county boundaries (default True).

ax

axes instance (overrides default axes instance)

zorder

sets the zorder for the county boundaries (if not specified, uses default zorder for matplotlib.patches.LineCollections).

returns a matplotlib.patches.LineCollection object.

drawcountries(linewidth=0.5, linestyle='solid', color='k', antialiased=1, ax=None, zorder=None)

Draw country boundaries.

Keyword

Description

linewidth

country boundary line width (default 0.5)

linestyle

coastline linestyle (default solid)

color

country boundary line color (default black)

antialiased

antialiasing switch for country boundaries (default True).

ax

axes instance (overrides default axes instance)

zorder

sets the zorder for the country boundaries (if not specified uses default zorder for matplotlib.patches.LineCollections).

returns a matplotlib.patches.LineCollection object.

drawgreatcircle(lon1, lat1, lon2, lat2, del_s=100.0, **kwargs)

Draw a great circle on the map from the longitude-latitude pair lon1,lat1 to lon2,lat2

Keyword

Description

del_s

points on great circle computed every del_s kilometers (default 100).

**kwargs

other keyword arguments are passed on to plot() method of Basemap instance.

Returns a list with a single matplotlib.lines.Line2D object like a call to pyplot.plot().

drawlsmask(land_color='0.8', ocean_color='w', lsmask=None, lsmask_lons=None, lsmask_lats=None, lakes=True, resolution='l', grid=5, **kwargs)

Draw land-sea mask image.

Note

The land-sea mask image cannot be overlaid on top of other images, due to limitations in matplotlib image handling (you can’t specify the zorder of an image).

Keywords

Description

land_color

desired land color (color name or rgba tuple). Default gray (“0.8”).

ocean_color

desired water color (color name or rgba tuple). Default white.

lsmask

An array of 0’s for ocean pixels, 1’s for land pixels and 2’s for lake/pond pixels. Default is None (default 5-minute resolution land-sea mask is used).

lakes

Plot lakes and ponds (Default True)

lsmask_lons

1d array of longitudes for lsmask (ignored if lsmask is None). Longitudes must be ordered from -180 W eastward.

lsmask_lats

1d array of latitudes for lsmask (ignored if lsmask is None). Latitudes must be ordered from -90 S northward.

resolution

gshhs coastline resolution used to define land/sea mask (default ‘l’, available ‘c’,’l’,’i’,’h’ or ‘f’)

grid

land/sea mask grid spacing in minutes (Default 5; 10, 2.5 and 1.25 are also available).

**kwargs

extra keyword arguments passed on to imshow()

If any of the lsmask, lsmask_lons or lsmask_lats keywords are not set, the built in GSHHS land-sea mask datasets are used.

Extra keyword ax can be used to override the default axis instance.

returns a matplotlib.image.AxesImage instance.

drawmapboundary(color='k', linewidth=1.0, fill_color=None, zorder=None, ax=None)

draw boundary around map projection region, optionally filling interior of region.

Keyword

Description

linewidth

line width for boundary (default 1.)

color

color of boundary line (default black)

fill_color

fill the map region background with this color (default is to fill with axis background color). If set to the string ‘none’, no filling is done.

zorder

sets the zorder for filling map background (default 0).

ax

axes instance to use (default None, use default axes instance).

returns matplotlib.collections.PatchCollection representing map boundary.

drawmapscale(lon, lat, lon0, lat0, length, barstyle='simple', units='km', fontsize=9, yoffset=None, labelstyle='simple', fontcolor='k', fillcolor1='w', fillcolor2='k', ax=None, format='%d', zorder=None, linecolor=None, linewidth=None)

Draw a map scale at lon,lat of length length representing distance in the map projection coordinates at lon0,lat0.

Keywords

Description

units

the units of the length argument (Default km).

barstyle

simple or fancy (roughly corresponding to the styles provided by Generic Mapping Tools). Default simple.

fontsize

for map scale annotations, default 9.

fontcolor

for map scale annotations, default black.

labelstyle

simple (default) or fancy. For fancy the map scale factor (ratio betwee the actual distance and map projection distance at lon0,lat0) and the value of lon0,lat0 are also displayed on the top of the scale bar. For simple, just the units are display on top and the distance below the scale bar. If equal to False, plot an empty label.

format

a string formatter to format numeric values

yoffset

yoffset controls how tall the scale bar is, and how far the annotations are offset from the scale bar. Default is 0.02 times the height of the map (0.02*(self.ymax-self.ymin)).

fillcolor1(2)

colors of the alternating filled regions (default white and black). Only relevant for ‘fancy’ barstyle.

zorder

sets the zorder for the map scale.

linecolor

sets the color of the scale, by default, fontcolor is used

linewidth

linewidth for scale and ticks

Extra keyword ax can be used to override the default axis instance.

drawmeridians(meridians, color='k', textcolor='k', linewidth=1.0, zorder=None, dashes=[1, 1], labels=[0, 0, 0, 0], labelstyle=None, fmt='%g', xoffset=None, yoffset=None, ax=None, latmax=None, **text_kwargs)

Draw and label meridians (longitude lines) for values (in degrees) given in the sequence meridians.

Keyword

Description

color

color to draw meridians (default black).

textcolor

color to draw labels (default black).

linewidth

line width for meridians (default 1.)

zorder

sets the zorder for meridians (if not specified, uses default zorder for matplotlib.lines.Line2D objects).

dashes

dash pattern for meridians (default [1,1], i.e. 1 pixel on, 1 pixel off).

labels

list of 4 values (default [0,0,0,0]) that control whether meridians are labelled where they intersect the left, right, top or bottom of the plot. For example labels=[1,0,0,1] will cause meridians to be labelled where they intersect the left and and bottom of the plot, but not the right and top.

labelstyle

if set to “+/-“, east and west longitudes are labelled with “+” and “-“, otherwise they are labelled with “E” and “W”.

fmt

a format string to format the meridian labels (default ‘%g’) or a function that takes a longitude value in degrees as it’s only argument and returns a formatted string.

xoffset

label offset from edge of map in x-direction (default is 0.01 times width of map in map projection coordinates).

yoffset

label offset from edge of map in y-direction (default is 0.01 times height of map in map projection coordinates).

ax

axes instance (overrides default axes instance)

latmax

absolute value of latitude to which meridians are drawn (default is 80).

**text_kwargs

additional keyword arguments controlling text for labels that are passed on to the text method of the axes instance (see matplotlib.pyplot.text documentation).

returns a dictionary whose keys are the meridian values, and whose values are tuples containing lists of the matplotlib.lines.Line2D and matplotlib.text.Text instances associated with each meridian. Deleting an item from the dictionary removes the correpsonding meridian from the plot.

drawparallels(circles, color='k', textcolor='k', linewidth=1.0, zorder=None, dashes=[1, 1], labels=[0, 0, 0, 0], labelstyle=None, fmt='%g', xoffset=None, yoffset=None, ax=None, latmax=None, **text_kwargs)

Draw and label parallels (latitude lines) for values (in degrees) given in the sequence circles.

Keyword

Description

color

color to draw parallels (default black).

textcolor

color to draw labels (default black).

linewidth

line width for parallels (default 1.)

zorder

sets the zorder for parallels (if not specified, uses default zorder for matplotlib.lines.Line2D objects).

dashes

dash pattern for parallels (default [1,1], i.e. 1 pixel on, 1 pixel off).

labels

list of 4 values (default [0,0,0,0]) that control whether parallels are labelled where they intersect the left, right, top or bottom of the plot. For example labels=[1,0,0,1] will cause parallels to be labelled where they intersect the left and and bottom of the plot, but not the right and top.

labelstyle

if set to “+/-“, north and south latitudes are labelled with “+” and “-“, otherwise they are labelled with “N” and “S”.

fmt

a format string to format the parallel labels (default ‘%g’) or a function that takes a latitude value in degrees as it’s only argument and returns a formatted string.

xoffset

label offset from edge of map in x-direction (default is 0.01 times width of map in map projection coordinates).

yoffset

label offset from edge of map in y-direction (default is 0.01 times height of map in map projection coordinates).

ax

axes instance (overrides default axes instance)

latmax

absolute value of latitude to which meridians are drawn (default is 80).

**text_kwargs

additional keyword arguments controlling text for labels that are passed on to the text method of the axes instance (see matplotlib.pyplot.text documentation).

returns a dictionary whose keys are the parallel values, and whose values are tuples containing lists of the matplotlib.lines.Line2D and matplotlib.text.Text instances associated with each parallel. Deleting an item from the dictionary removes the corresponding parallel from the plot.

drawrivers(linewidth=0.5, linestyle='solid', color='k', antialiased=1, ax=None, zorder=None)

Draw major rivers.

Keyword

Description

linewidth

river boundary line width (default 0.5)

linestyle

coastline linestyle (default solid)

color

river boundary line color (default black)

antialiased

antialiasing switch for river boundaries (default True).

ax

axes instance (overrides default axes instance)

zorder

sets the zorder for the rivers (if not specified uses default zorder for matplotlib.patches.LineCollections).

returns a matplotlib.patches.LineCollection object.

drawstates(linewidth=0.5, linestyle='solid', color='k', antialiased=1, ax=None, zorder=None)

Draw state boundaries in Americas.

Keyword

Description

linewidth

state boundary line width (default 0.5)

linestyle

coastline linestyle (default solid)

color

state boundary line color (default black)

antialiased

antialiasing switch for state boundaries (default True).

ax

axes instance (overrides default axes instance)

zorder

sets the zorder for the state boundaries (if not specified, uses default zorder for matplotlib.patches.LineCollections).

returns a matplotlib.patches.LineCollection object.

etopo(ax=None, scale=None, **kwargs)

display etopo relief image (from http://www.ngdc.noaa.gov/mgg/global/global.html) as map background. Default image size is 5400x2700, which can be quite slow and use quite a bit of memory. The scale keyword can be used to downsample the image (scale=0.5 downsamples to 5400x2700).

**kwargs passed on to imshow().

returns a matplotlib.image.AxesImage instance.

fillcontinents(color='0.8', lake_color=None, ax=None, zorder=None, alpha=None)

Fill continents.

Keyword

Description

color

color to fill continents (default gray).

lake_color

color to fill inland lakes (default axes background).

ax

axes instance (overrides default axes instance).

zorder

sets the zorder for the continent polygons (if not specified, uses default zorder for a Polygon patch). Set to zero if you want to paint over the filled continents).

alpha

sets alpha transparency for continent polygons

After filling continents, lakes are re-filled with axis background color.

returns a list of matplotlib.patches.Polygon objects.

gcpoints(lon1, lat1, lon2, lat2, npoints)

compute points points along a great circle with endpoints (lon1,lat1) and (lon2,lat2).

Returns arrays x,y with map projection coordinates.

hexbin(x, y, **kwargs)

Make a hexagonal binning plot of x versus y, where x, y are 1-D sequences of the same length, N. If C is None (the default), this is a histogram of the number of occurences of the observations at (x[i],y[i]).

If C is specified, it specifies values at the coordinate (x[i],y[i]). These values are accumulated for each hexagonal bin and then reduced according to reduce_C_function, which defaults to the numpy mean function (np.mean). (If C is specified, it must also be a 1-D sequence of the same length as x and y.)

x, y and/or C may be masked arrays, in which case only unmasked points will be plotted.

(see matplotlib.pyplot.hexbin documentation).

Extra keyword ax can be used to override the default axis instance.

Other **kwargs passed on to matplotlib.pyplot.hexbin

imshow(*args, **kwargs)

Display an image over the map (see matplotlib.pyplot.imshow documentation).

extent and origin keywords set automatically so image will be drawn over map region.

Extra keyword ax can be used to override the default axis instance.

Other **kwargs passed on to matplotlib.pyplot.plot.

returns an matplotlib.image.AxesImage instance.

is_land(xpt, ypt)

Returns True if the given x,y point (in projection coordinates) is over land, False otherwise. The definition of land is based upon the GSHHS coastline polygons associated with the class instance. Points over lakes inside land regions are not counted as land points.

makegrid(nx, ny, returnxy=False)

return arrays of shape (ny,nx) containing lon,lat coordinates of an equally spaced native projection grid.

If returnxy = True, the x,y values of the grid are returned also.

nightshade(date, color='k', delta=0.25, alpha=0.5, ax=None, zorder=2)

Shade the regions of the map that are in darkness at the time specifed by date. date is a datetime instance, assumed to be UTC.

Keywords

Description

color

color to shade night regions (default black).

delta

day/night terminator is computed with a a resolution of delta degrees (default 0.25).

alpha

alpha transparency for shading (default 0.5, so map background shows through).

zorder

zorder for shading (default 2).

Extra keyword ax can be used to override the default axis instance.

returns a matplotlib.contour.ContourSet instance.

pcolor(x, y, data, **kwargs)

Make a pseudo-color plot over the map (see matplotlib.pyplot.pcolor documentation).

If latlon keyword is set to True, x,y are intrepreted as longitude and latitude in degrees. Data and longitudes are automatically shifted to match map projection region for cylindrical and pseudocylindrical projections, and x,y are transformed to map projection coordinates. If latlon is False (default), x and y are assumed to be map projection coordinates.

If x or y are outside projection limb (i.e. they have values > 1.e20) they will be convert to masked arrays with those values masked. As a result, those values will not be plotted.

If tri is set to True, an unstructured grid is assumed (x,y,data must be 1-d) and matplotlib.pyplot.tripcolor is used.

Extra keyword ax can be used to override the default axis instance.

Other **kwargs passed on to matplotlib.pyplot.pcolor (or tripcolor if tri=True).

Note: (taken from matplotlib.pyplot.pcolor documentation) Ideally the dimensions of x and y should be one greater than those of data; if the dimensions are the same, then the last row and column of data will be ignored.

pcolormesh(x, y, data, **kwargs)

Make a pseudo-color plot over the map (see matplotlib.pyplot.pcolormesh documentation).

If latlon keyword is set to True, x,y are intrepreted as longitude and latitude in degrees. Data and longitudes are automatically shifted to match map projection region for cylindrical and pseudocylindrical projections, and x,y are transformed to map projection coordinates. If latlon is False (default), x and y are assumed to be map projection coordinates.

Extra keyword ax can be used to override the default axis instance.

Other **kwargs passed on to matplotlib.pyplot.pcolormesh.

Note: (taken from matplotlib.pyplot.pcolor documentation) Ideally the dimensions of x and y should be one greater than those of data; if the dimensions are the same, then the last row and column of data will be ignored.

plot(*args, **kwargs)

Draw lines and/or markers on the map (see matplotlib.pyplot.plot documentation).

If latlon keyword is set to True, x,y are intrepreted as longitude and latitude in degrees. Data and longitudes are automatically shifted to match map projection region for cylindrical and pseudocylindrical projections, and x,y are transformed to map projection coordinates. If latlon is False (default), x and y are assumed to be map projection coordinates.

Extra keyword ax can be used to override the default axis instance.

Other **kwargs passed on to matplotlib.pyplot.plot.

quiver(x, y, u, v, *args, **kwargs)

Make a vector plot (u, v) with arrows on the map.

Arguments may be 1-D or 2-D arrays or sequences (see matplotlib.pyplot.quiver documentation for details).

If latlon keyword is set to True, x,y are intrepreted as longitude and latitude in degrees. Data and longitudes are automatically shifted to match map projection region for cylindrical and pseudocylindrical projections, and x,y are transformed to map projection coordinates. If latlon is False (default), x and y are assumed to be map projection coordinates.

Extra keyword ax can be used to override the default axis instance.

Other *args and **kwargs passed on to matplotlib.pyplot.quiver.

readshapefile(shapefile, name, drawbounds=True, zorder=None, linewidth=0.5, color='k', antialiased=1, ax=None, default_encoding='utf-8')

Read in shape file, optionally draw boundaries on map.

Note

  • Assumes shapes are 2D

  • only works for Point, MultiPoint, Polyline and Polygon shapes.

  • vertices/points must be in geographic (lat/lon) coordinates.

Mandatory Arguments:

Argument

Description

shapefile

path to shapefile components. Example: shapefile=’/home/jeff/esri/world_borders’ assumes that world_borders.shp, world_borders.shx and world_borders.dbf live in /home/jeff/esri.

name

name for Basemap attribute to hold the shapefile vertices or points in map projection coordinates. Class attribute name+’_info’ is a list of dictionaries, one for each shape, containing attributes of each shape from dbf file, For example, if name=’counties’, self.counties will be a list of x,y vertices for each shape in map projection coordinates and self.counties_info will be a list of dictionaries with shape attributes. Rings in individual Polygon shapes are split out into separate polygons, and additional keys ‘RINGNUM’ and ‘SHAPENUM’ are added to the shape attribute dictionary.

The following optional keyword arguments are only relevant for Polyline and Polygon shape types, for Point and MultiPoint shapes they are ignored.

Keyword

Description

drawbounds

draw boundaries of shapes (default True).

zorder

shape boundary zorder (if not specified, default for mathplotlib.lines.LineCollection is used).

linewidth

shape boundary line width (default 0.5)

color

shape boundary line color (default black)

antialiased

antialiasing switch for shape boundaries (default True).

ax

axes instance (overrides default axes instance)

A tuple (num_shapes, type, min, max) containing shape file info is returned. num_shapes is the number of shapes, type is the type code (one of the SHPT* constants defined in the shapelib module, see http://shapelib.maptools.org/shp_api.html) and min and max are 4-element lists with the minimum and maximum values of the vertices. If drawbounds=True a matplotlib.patches.LineCollection object is appended to the tuple.

rotate_vector(uin, vin, lons, lats, returnxy=False)

Rotate a vector field (uin,vin) on a rectilinear grid with longitudes = lons and latitudes = lats from geographical (lat/lon) into map projection (x/y) coordinates.

Differs from transform_vector in that no interpolation is done. The vector is returned on the same grid, but rotated into x,y coordinates.

The input vector field is defined in spherical coordinates (it has eastward and northward components) while the output vector field is rotated to map projection coordinates (relative to x and y). The magnitude of the vector is preserved.

Arguments

Description

uin, vin

input vector field on a lat/lon grid.

lons, lats

Arrays containing longitudes and latitudes (in degrees) of input data in increasing order. For non-cylindrical projections (those other than cyl, merc, cyl, gall and mill) lons must fit within range -180 to 180.

Returns uout, vout (rotated vector field). If the optional keyword argument returnxy is True (default is False), returns uout,vout,x,y (where x,y are the map projection coordinates of the grid defined by lons,lats).

scatter(*args, **kwargs)

Plot points with markers on the map (see matplotlib.pyplot.scatter documentation).

If latlon keyword is set to True, x,y are intrepreted as longitude and latitude in degrees. Data and longitudes are automatically shifted to match map projection region for cylindrical and pseudocylindrical projections, and x,y are transformed to map projection coordinates. If latlon is False (default), x and y are assumed to be map projection coordinates.

Extra keyword ax can be used to override the default axes instance.

Other **kwargs passed on to matplotlib.pyplot.scatter.

set_axes_limits(ax=None)

Final step in Basemap method wrappers of Axes plotting methods:

Set axis limits, fix aspect ratio for map domain using current or specified axes instance. This is done only once per axes instance.

In interactive mode, this method always calls draw_if_interactive before returning.

shadedrelief(ax=None, scale=None, **kwargs)

display shaded relief image (from http://www.shadedrelief.com) as map background. Default image size is 10800x5400, which can be quite slow and use quite a bit of memory. The scale keyword can be used to downsample the image (scale=0.5 downsamples to 5400x2700).

**kwargs passed on to imshow().

returns a matplotlib.image.AxesImage instance.

shiftdata(lonsin, datain=None, lon_0=None, fix_wrap_around=True)

Shift longitudes (and optionally data) so that they match map projection region. Only valid for cylindrical/pseudo-cylindrical global projections and data on regular lat/lon grids. longitudes and data can be 1-d or 2-d, if 2-d it is assumed longitudes are 2nd (rightmost) dimension.

Arguments

Description

lonsin

original 1-d or 2-d longitudes.

if datain given, returns dataout,lonsout (data and longitudes shifted to fit in interval [lon_0-180,lon_0+180]), otherwise just returns longitudes. If transformed longitudes lie outside map projection region, data is masked and longitudes are set to 1.e30.

streamplot(x, y, u, v, *args, **kwargs)

Draws streamlines of a vector flow. (see matplotlib.pyplot.streamplot documentation).

If latlon keyword is set to True, x,y are intrepreted as longitude and latitude in degrees. Data and longitudes are automatically shifted to match map projection region for cylindrical and pseudocylindrical projections, and x,y are transformed to map projection coordinates. If latlon is False (default), x and y are assumed to be map projection coordinates.

Extra keyword ax can be used to override the default axis instance.

Other *args and **kwargs passed on to matplotlib.pyplot.streamplot.

tissot(lon_0, lat_0, radius_deg, npts, ax=None, **kwargs)

Draw a polygon centered at lon_0,lat_0. The polygon approximates a circle on the surface of the earth with radius radius_deg degrees latitude along longitude lon_0, made up of npts vertices. The polygon represents a Tissot’s indicatrix (http://en.wikipedia.org/wiki/Tissot’s_Indicatrix), which when drawn on a map shows the distortion inherent in the map projection.

Note

Cannot handle situations in which the polygon intersects the edge of the map projection domain, and then re-enters the domain.

Extra keyword ax can be used to override the default axis instance.

Other **kwargs passed on to matplotlib.patches.Polygon.

returns a matplotlib.patches.Polygon object.

transform_scalar(datin, lons, lats, nx, ny, returnxy=False, checkbounds=False, order=1, masked=False)

Interpolate a scalar field (datin) from a lat/lon grid with longitudes = lons and latitudes = lats to a ny by nx map projection grid. Typically used to transform data to map projection coordinates for plotting on a map with the imshow().

Argument

Description

datin

input data on a lat/lon grid.

lons, lats

rank-1 arrays containing longitudes and latitudes (in degrees) of input data in increasing order. For non-cylindrical projections (those other than cyl, merc, cea, gall and mill) lons must fit within range -180 to 180.

nx, ny

The size of the output regular grid in map projection coordinates

Keyword

Description

returnxy

If True, the x and y values of the map projection grid are also returned (Default False).

checkbounds

If True, values of lons and lats are checked to see that they lie within the map projection region. Default is False, and data outside map projection region is clipped to values on boundary.

masked

If True, interpolated data is returned as a masked array with values outside map projection region masked (Default False).

order

0 for nearest-neighbor interpolation, 1 for bilinear, 3 for cubic spline (Default 1). Cubic spline interpolation requires scipy.ndimage.

Returns datout (data on map projection grid). If returnxy=True, returns data,x,y.

transform_vector(uin, vin, lons, lats, nx, ny, returnxy=False, checkbounds=False, order=1, masked=False)

Rotate and interpolate a vector field (uin,vin) from a lat/lon grid with longitudes = lons and latitudes = lats to a ny by nx map projection grid.

The input vector field is defined in spherical coordinates (it has eastward and northward components) while the output vector field is rotated to map projection coordinates (relative to x and y). The magnitude of the vector is preserved.

Arguments

Description

uin, vin

input vector field on a lat/lon grid.

lons, lats

rank-1 arrays containing longitudes and latitudes (in degrees) of input data in increasing order. For non-cylindrical projections (those other than cyl, merc, cea, gall and mill) lons must fit within range -180 to 180.

nx, ny

The size of the output regular grid in map projection coordinates

Keyword

Description

returnxy

If True, the x and y values of the map projection grid are also returned (Default False).

checkbounds

If True, values of lons and lats are checked to see that they lie within the map projection region. Default is False, and data outside map projection region is clipped to values on boundary.

masked

If True, interpolated data is returned as a masked array with values outside map projection region masked (Default False).

order

0 for nearest-neighbor interpolation, 1 for bilinear, 3 for cubic spline (Default 1). Cubic spline interpolation requires scipy.ndimage.

Returns uout, vout (vector field on map projection grid). If returnxy=True, returns uout,vout,x,y.

warpimage(image='bluemarble', scale=None, **kwargs)

Display an image (filename given by image keyword) as a map background. If image is a URL (starts with ‘http’), it is downloaded to a temp file using urllib.urlretrieve.

Default (if image not specified) is to display ‘blue marble next generation’ image from http://visibleearth.nasa.gov/.

Specified image must have pixels covering the whole globe in a regular lat/lon grid, starting and -180W and the South Pole. Works with the global images from http://earthobservatory.nasa.gov/Features/BlueMarble/BlueMarble_monthlies.php.

The scale keyword can be used to downsample (rescale) the image. Values less than 1.0 will speed things up at the expense of image resolution.

Extra keyword ax can be used to override the default axis instance.

**kwargs passed on to imshow().

returns a matplotlib.image.AxesImage instance.

wmsimage(server, xpixels=400, ypixels=None, format='png', alpha=None, verbose=False, **kwargs)

Retrieve an image using from a WMS server using the Open Geospatial Consortium (OGC) standard interface and display on the map. Requires OWSLib (http://pypi.python.org/pypi/OWSLib). In order to use this method, the Basemap instance must be created using the epsg keyword to define the map projection, unless the cyl projection is used (in which case the epsg code 4326 is assumed).

Keywords

Description

server

WMS server URL.

xpixels

requested number of image pixels in x-direction (default 400).

ypixels

requested number of image pixels in y-direction. Default (None) is to infer the number from from xpixels and the aspect ratio of the map projection region.

format

desired image format (default ‘png’)

alpha

The alpha blending value, between 0 (transparent) and 1 (opaque) (default None)

verbose

if True, print WMS server info (default False).

**kwargs

extra keyword arguments passed on to OWSLib.wms.WebMapService.getmap.

Extra keyword ax can be used to override the default axis instance.

returns a matplotlib.image.AxesImage instance.

mpl_toolkits.basemap.addcyclic(*arr, **kwargs)

Adds cyclic (wraparound) points in longitude to one or several arrays, the last array being longitudes in degrees. e.g.

data1out, data2out, lonsout = addcyclic(data1,data2,lons)

Keywords

Description

axis

the dimension representing longitude (default -1, or right-most)

cyclic

width of periodic domain (default 360)

mpl_toolkits.basemap.interp(datain, xin, yin, xout, yout, checkbounds=False, masked=False, order=1)

Interpolate data (datain) on a rectilinear grid (with x = xin y = yin) to a grid with x = xout, y= yout.

Arguments

Description

datain

a rank-2 array with 1st dimension corresponding to y, 2nd dimension x.

xin, yin

rank-1 arrays containing x and y of datain grid in increasing order.

xout, yout

rank-2 arrays containing x and y of desired output grid.

Keywords

Description

checkbounds

If True, values of xout and yout are checked to see that they lie within the range specified by xin and xin. If False, and xout,yout are outside xin,yin, interpolated values will be clipped to values on boundary of input grid (xin,yin) Default is False.

masked

If True, points outside the range of xin and yin are masked (in a masked array). If masked is set to a number, then points outside the range of xin and yin will be set to that number. Default False.

order

0 for nearest-neighbor interpolation, 1 for bilinear interpolation, 3 for cublic spline (default 1). order=3 requires scipy.ndimage.

Note

If datain is a masked array and order=1 (bilinear interpolation) is used, elements of dataout will be masked if any of the four surrounding points in datain are masked. To avoid this, do the interpolation in two passes, first with order=1 (producing dataout1), then with order=0 (producing dataout2). Then replace all the masked values in dataout1 with the corresponding elements in dataout2 (using numpy.where). This effectively uses nearest neighbor interpolation if any of the four surrounding points in datain are masked, and bilinear interpolation otherwise.

Returns dataout, the interpolated data on the grid xout, yout.

mpl_toolkits.basemap.maskoceans(lonsin, latsin, datain, inlands=True, resolution='l', grid=5)

mask data (datain), defined on a grid with latitudes latsin longitudes lonsin so that points over water will not be plotted.

Arguments

Description

lonsin, latsin

rank-2 arrays containing longitudes and latitudes of grid.

datain

rank-2 input array on grid defined by lonsin and latsin.

inlands

if False, masked only ocean points and not inland lakes (Default True).

resolution

gshhs coastline resolution used to define land/sea mask (default ‘l’, available ‘c’,’l’,’i’,’h’ or ‘f’)

grid

land/sea mask grid spacing in minutes (Default 5; 10, 2.5 and 1.25 are also available).

returns a masked array the same shape as datain with “wet” points masked.

mpl_toolkits.basemap.shiftgrid(lon0, datain, lonsin, start=True, cyclic=360.0)

Shift global lat/lon grid east or west.

Arguments

Description

lon0

starting longitude for shifted grid (ending longitude if start=False). lon0 must be on input grid (within the range of lonsin).

datain

original data with longitude the right-most dimension.

lonsin

original longitudes.

Keywords

Description

start

if True, lon0 represents the starting longitude of the new grid. if False, lon0 is the ending longitude. Default True.

cyclic

width of periodic domain (default 360)

returns dataout,lonsout (data and longitudes on shifted grid).