Reputation: 769
The following code creates a map of continental US states that is shaded by population density. I want to create a similar map (my data is not actually pop density, but this is an easy example), except that it also includes the states of Alaska and Hawaii.
Specifically I would like to have Alaska/Hawaii show up in the figure, but be moved so that they are below the part of the figure showing the continental US. Or something along those lines.
Any idea how I would create such a map using Cartopy?
import matplotlib.pyplot as plt
import cartopy.crs as ccrs
import cartopy.io.shapereader as shpreader
fig = plt.figure()
ax = fig.add_axes([0, 0, 1, 1], projection=ccrs.LambertConformal())
ax.set_extent([-125, -66.5, 20, 50], ccrs.Geodetic())
shapename = 'admin_1_states_provinces_lakes_shp'
states_shp = shpreader.natural_earth(resolution='110m',
category='cultural', name=shapename)
popdensity = {
'New Jersey': 438.00,
'Rhode Island': 387.35,
'Massachusetts': 312.68,
'Connecticut': 271.40,
'Maryland': 209.23,
'New York': 155.18,
'Delaware': 154.87,
'Florida': 114.43,
'Ohio': 107.05,
'Pennsylvania': 105.80,
'Illinois': 86.27,
'California': 83.85,
'Virginia': 69.03,
'Michigan': 67.55,
'Indiana': 65.46,
'North Carolina': 63.80,
'Georgia': 54.59,
'Tennessee': 53.29,
'New Hampshire': 53.20,
'South Carolina': 51.45,
'Louisiana': 39.61,
'Kentucky': 39.28,
'Wisconsin': 38.13,
'Washington': 34.20,
'Alabama': 33.84,
'Missouri': 31.36,
'Texas': 30.75,
'West Virginia': 29.00,
'Vermont': 25.41,
'Minnesota': 23.86,
'Mississippi': 23.42,
'Iowa': 20.22,
'Arkansas': 19.82,
'Oklahoma': 19.40,
'Arizona': 17.43,
'Colorado': 16.01,
'Maine': 15.95,
'Oregon': 13.76,
'Kansas': 12.69,
'Utah': 10.50,
'Nebraska': 8.60,
'Nevada': 7.03,
'Idaho': 6.04,
'New Mexico': 5.79,
'South Dakota': 3.84,
'North Dakota': 3.59,
'Montana': 2.39,
'Wyoming': 1.96}
ax.background_patch.set_visible(False)
ax.outline_patch.set_visible(False)
ax.set_title('State Population Density')
for state in shpreader.Reader(states_shp).records():
edgecolor = 'black'
try:
# use the name of this state to get pop_density
state_dens = popdensity[ state.attributes['name'] ]
except:
state_dens = 0
# simple scheme to assign color to each state
if state_dens < 40:
facecolor = "lightyellow"
elif state_dens > 200:
facecolor = "red"
else:
facecolor = "pink"
# `state.geometry` is the polygon to plot
ax.add_geometries([state.geometry], ccrs.PlateCarree(),
facecolor=facecolor, edgecolor=edgecolor)
plt.show()
The figure this (currently) creates is as follows:
Upvotes: 0
Views: 1962
Reputation: 18812
To plot inset maps as parts of a main map is challenging. You will need to create an axes
for plotting each inset map and place it on the figure
at proper location and relative scale. Here is a working code that you can experiment with.
import matplotlib.pyplot as plt
import cartopy
import cartopy.crs as ccrs
import cartopy.io.shapereader as shpreader
import shapely.geometry as sgeom
# A function that draws inset map, ++
# ===================================
def add_insetmap(axes_extent, map_extent, state_name, facecolor, edgecolor, geometry):
# create new axes, set its projection
use_projection = ccrs.Mercator() # preserve shape well
#use_projection = ccrs.PlateCarree() # large distortion in E-W for Alaska
geodetic = ccrs.Geodetic(globe=ccrs.Globe(datum='WGS84'))
sub_ax = plt.axes(axes_extent, projection=use_projection) # normal units
sub_ax.set_extent(map_extent, geodetic) # map extents
# add basic land, coastlines of the map
# you may comment out if you don't need them
sub_ax.add_feature(cartopy.feature.LAND)
sub_ax.coastlines()
sub_ax.set_title(state_name)
# add map `geometry` here
sub_ax.add_geometries([geometry], ccrs.PlateCarree(), \
facecolor=facecolor, edgecolor=edgecolor)
# +++ more features can be added here +++
# plot box around the map
extent_box = sgeom.box(map_extent[0], map_extent[2], map_extent[1], map_extent[3])
sub_ax.add_geometries([extent_box], ccrs.PlateCarree(), color='none', linewidth=0.05)
fig = plt.figure()
ax = fig.add_axes([0, 0, 1, 1], projection=ccrs.LambertConformal())
ax.set_extent([-125, -66.5, 20, 50], ccrs.Geodetic())
shapename = 'admin_1_states_provinces_lakes_shp'
states_shp = shpreader.natural_earth(resolution='110m',
category='cultural', name=shapename)
popdensity = {
'New Jersey': 438.00,
'Rhode Island': 387.35,
'Massachusetts': 312.68,
'Connecticut': 271.40,
'Maryland': 209.23,
'New York': 155.18,
'Delaware': 154.87,
'Florida': 114.43,
'Ohio': 107.05,
'Pennsylvania': 105.80,
'Illinois': 86.27,
'California': 83.85,
'Virginia': 69.03,
'Michigan': 67.55,
'Indiana': 65.46,
'North Carolina': 63.80,
'Georgia': 54.59,
'Tennessee': 53.29,
'New Hampshire': 53.20,
'South Carolina': 51.45,
'Louisiana': 39.61,
'Kentucky': 39.28,
'Wisconsin': 38.13,
'Washington': 34.20,
'Alabama': 33.84,
'Missouri': 31.36,
'Texas': 30.75,
'West Virginia': 29.00,
'Vermont': 25.41,
'Minnesota': 23.86,
'Mississippi': 23.42,
'Iowa': 20.22,
'Arkansas': 19.82,
'Oklahoma': 19.40,
'Arizona': 17.43,
'Colorado': 16.01,
'Maine': 15.95,
'Oregon': 13.76,
'Kansas': 12.69,
'Utah': 10.50,
'Nebraska': 8.60,
'Nevada': 7.03,
'Idaho': 6.04,
'New Mexico': 5.79,
'South Dakota': 3.84,
'North Dakota': 3.59,
'Montana': 2.39,
'Wyoming': 1.96}
ax.background_patch.set_visible(False)
ax.outline_patch.set_visible(False)
ax.set_title('State Population Density')
for state in shpreader.Reader(states_shp).records():
edgecolor = 'black'
try:
# use the name of this state to get pop_density
state_dens = popdensity[ state.attributes['name'] ]
except:
state_dens = 0
# simple scheme to assign color to each state
if state_dens < 40:
facecolor = "lightyellow"
elif state_dens > 200:
facecolor = "red"
else:
facecolor = "pink"
# special handling for the 2 states
# ---------------------------------
if state.attributes['name'] in ("Alaska", "Hawaii"):
# print("state.attributes['name']:", state.attributes['name'])
state_name = state.attributes['name']
# prep map settings
# experiment with the numbers in both `_extents` for your best results
if state_name == "Alaska":
# (1) Alaska
map_extent = (-178, -135, 46, 73) # degrees: (lonmin,lonmax,latmin,latmax)
axes_extent = (0.04, 0.06, 0.29, 0.275) # axes units: 0 to 1, (LLx,LLy,width,height)
if state_name == "Hawaii":
# (2) Hawii
map_extent = (-162, -152, 15, 25)
axes_extent = (0.27, 0.06, 0.15, 0.15)
# add inset maps
add_insetmap(axes_extent, map_extent, state_name, \
facecolor, \
edgecolor, \
state.geometry)
# the other (conterminous) states go here
else:
# `state.geometry` is the polygon to plot
ax.add_geometries([state.geometry], ccrs.PlateCarree(),
facecolor=facecolor, edgecolor=edgecolor)
plt.show()
The output plot will be:
Upvotes: 5