Cartopy is a Python package designed for geospatial data processing in order to produce maps and other geospatial data analyses. Making a map with Python and Cartopy is actually easier than you might think, but still it has its own tricks. Once you get the idea how it is done, hopefully it becomes a bit easier (I am still waiting for this moment myself!). As I was really fascinated with what the library can do, I have decided to write it down. Hopefully, it will help you as well.
These are the topics in the article:
- Higlight a single country
- What data cartopy has per country?
- How can we make a heat map of the “richest countries”?
- Countries with multiple territories
- How to put a label over a country?
- Where are the rivers and the lakes?
- Coloring a country with two colors, dividing it into North and South
- Coloring a small country / city
- GitHub link to the Jupyter Notebook with all examples
Higlight a single country
Pretty much, every country has its shape, written on the Cartopy library and could be referred by name. For example, if you want to show Bulgaria in green on the map, this is ok:
|
1 2 3 4 5 6 |
import cartopy import cartopy.crs as ccrs import cartopy.io.shapereader as shpreader import matplotlib.pyplot as plt import matplotlib.patheffects as PathEffects from shapely.geometry import LineString, MultiLineString |
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 |
ax = plt.axes(projection=ccrs.PlateCarree()) ax.add_feature(cartopy.feature.LAND) ax.add_feature(cartopy.feature.OCEAN) ax.add_feature(cartopy.feature.COASTLINE) ax.add_feature(cartopy.feature.BORDERS, linestyle='-', alpha=.5) ax.add_feature(cartopy.feature.LAKES, alpha=0.95) ax.add_feature(cartopy.feature.RIVERS) ax.set_extent([-150, 60, -25, 60]) shpfilename = shpreader.natural_earth(resolution='110m', category='cultural', name='admin_0_countries') reader = shpreader.Reader(shpfilename) countries = reader.records() for country in countries: if country.attributes['SOVEREIGNT'] == "Bulgaria": ax.add_geometries(country.geometry, ccrs.PlateCarree(), facecolor=(0, 1, 0)) else: ax.add_geometries(country.geometry, ccrs.PlateCarree(), facecolor=(1, 1, 1)) plt.rcParams["figure.figsize"] = (20,20) plt.show() |
What data cartopy has per country? How can we make a heat map of the “richest countries”?
How does it do it? If you add print(country) somewhere in the code, you would see this information, per country (of course, not formatted in such a nice way):
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 |
MULTIPOLYGON (( (22.65714969248299 44.23492300066128, 22.94483239105185 43.82378530534713,23.33230228037633 43.89701080990471, etc... { 'featurecla':'Admin-0 country', 'scalerank':1, 'LABELRANK':4, 'SOVEREIGNT':'Bulgaria', 'SOV_A3':'BGR', 'ADM0_DIF':0, 'LEVEL':2, 'TYPE':'Sovereign country', 'ADMIN':'Bulgaria', 'ADM0_A3':'BGR', 'GEOU_DIF':0, 'GEOUNIT':'Bulgaria', 'GU_A3':'BGR', 'SU_DIF':0, 'SUBUNIT':'Bulgaria', 'SU_A3':'BGR', 'BRK_DIFF':0, 'NAME':'Bulgaria', 'NAME_LONG':'Bulgaria', 'BRK_A3':'BGR', 'BRK_NAME':'Bulgaria', 'BRK_GROUP':'', 'ABBREV':'Bulg.', 'POSTAL':'BG', 'FORMAL_EN':'Republic of Bulgaria', 'FORMAL_FR':'', 'NAME_CIAWF':'Bulgaria', 'NOTE_ADM0':'', 'NOTE_BRK':'', 'NAME_SORT':'Bulgaria', 'NAME_ALT':'', 'MAPCOLOR7':4, 'MAPCOLOR8':5, 'MAPCOLOR9':1, 'MAPCOLOR13':8, 'POP_EST':7101510, 'POP_RANK':13, 'GDP_MD_EST':143100.0, 'POP_YEAR':2017, 'LASTCENSUS':2011, 'GDP_YEAR':2016, 'ECONOMY':'2. Developed region: nonG7', 'INCOME_GRP':'3. Upper middle income', 'WIKIPEDIA':-99, 'FIPS_10_':'BU', 'ISO_A2':'BG', 'ISO_A3':'BGR', 'ISO_A3_EH':'BGR', 'ISO_N3':'100', 'UN_A3':'100', 'WB_A2':'BG', 'WB_A3':'BGR', 'WOE_ID':23424771, 'WOE_ID_EH':23424771, 'WOE_NOTE':'Exact WOE match as country', 'ADM0_A3_IS':'BGR', 'ADM0_A3_US':'BGR', 'ADM0_A3_UN':-99, 'ADM0_A3_WB':-99, 'CONTINENT':'Europe', 'REGION_UN':'Europe', 'SUBREGION':'Eastern Europe', 'REGION_WB':'Europe & Central Asia', 'NAME_LEN':8, 'LONG_LEN':8, 'ABBREV_LEN':5, 'TINY':-99, 'HOMEPART':1, 'MIN_ZOOM':0.0, 'MIN_LABEL':4.0, 'MAX_LABEL':9.0, 'NE_ID':1159320409, 'WIKIDATAID':'Q219', 'NAME_AR':'بلغاريا', 'NAME_BN':'বুলগেরিয়া', 'NAME_DE':'Bulgarien', 'NAME_EN':'Bulgaria', 'NAME_ES':'Bulgaria', 'NAME_FR':'Bulgarie', 'NAME_EL':'Βουλγαρία', 'NAME_HI':'बुल्गारिया', 'NAME_HU':'Bulgária', 'NAME_ID':'Bulgaria', 'NAME_IT':'Bulgaria', 'NAME_JA':'ブルガリア', 'NAME_KO':'불가리아', 'NAME_NL':'Bulgarije', 'NAME_PL':'Bułgaria', 'NAME_PT':'Bulgária', 'NAME_RU':'Болгария', 'NAME_SV':'Bulgarien', 'NAME_TR':'Bulgaristan', 'NAME_VI':'Bulgaria', 'NAME_ZH':'保加利亚' |
As you see, there are the names of the country in 20+ languages, information for the subregion, economy, formal name, income gdp, wikipedia data aid (Q219) and other stuff. Thus, making a heat map of “rich countries” should be easy – there are 5 types of countries in the library, based on how much their GRP is:
|
1 2 3 4 5 |
1. High income: OECD 2. High income: nonOECD 3. Upper middle income 4. Lower middle income 5. Low income |
Thus a loop with 5 condition results pretty much in a map, split by the GRP:
|
1 2 3 4 5 6 7 8 9 10 11 |
for country in countries: if country.attributes['INCOME_GRP'].startswith("1"): ax.add_geometries(country.geometry, ccrs.PlateCarree(), facecolor=(0, 1, 1)) elif country.attributes['INCOME_GRP'].startswith("2"): ax.add_geometries(country.geometry, ccrs.PlateCarree(), facecolor=(0, 0, 1)) elif country.attributes['INCOME_GRP'].startswith("3"): ax.add_geometries(country.geometry, ccrs.PlateCarree(), facecolor=(0, 0, 0.5)) elif country.attributes['INCOME_GRP'].startswith("4"): ax.add_geometries(country.geometry, ccrs.PlateCarree(), facecolor=(1, 0, 0.5)) elif country.attributes['INCOME_GRP'].startswith("5"): ax.add_geometries(country.geometry, ccrs.PlateCarree(), facecolor=(1, 0, 0)) |
Not going to speak about the code, it is rather self-explanatory (or take a look at the github link of the Jupyter notebook at the end and run it).
Countries with multiple territories
The story is a bit different, if you put one of these countries – ["Denmark", "France", "Israel", "United Kingdom", "United States of America"], as these have some territories, which are under other jurisdiction. Thus, highlighting the “main country” in blue and the territory under other jurisdiction in green would result in this:
The code to achieve the map above is a bit different, and in it we differentiate between country.attributes['NAME'] and country.attributes['SOVEREIGNT']:
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 |
multiple_country_list = ["Denmark", "France", "Israel", "United Kingdom", "United States of America"] shpfilename = shpreader.natural_earth(resolution='110m', category='cultural', name='admin_0_countries') reader = shpreader.Reader(shpfilename) countries = reader.records() for country in countries: if country.attributes['SOVEREIGNT'] in multiple_country_list: if country.attributes['NAME'] != country.attributes['SOVEREIGNT']: ax.add_geometries(country.geometry, ccrs.PlateCarree(), facecolor=(0, 1, 0)) else: ax.add_geometries(country.geometry, ccrs.PlateCarree(), facecolor=(0, 0, 1)) else: ax.add_geometries(country.geometry, ccrs.PlateCarree(), facecolor=(0, 0, 0), label = country.attributes['SOVEREIGNT']) plt.rcParams["figure.figsize"] = (30,120) plt.show() |
How to put a label over a country?
Writing over a country was something that took me quite some time to achieve. At the end, I decided it would be easier, if I asked in StackOverflow and I got quite a good answer – the trick in writing over a country is to retrieve the centroid of the country.geometry and to plot the text at that location. If a good map size and font is used, you may write the name of every country in the library. The colors are a bit repeated, because of the code I have created only allows about 10 different colors. Anyway, the borders are visible.
|
1 2 3 4 5 6 7 8 9 |
for country in countries: my_color = my_color + .1 if my_color > 1: my_color = .01 g = ax.add_geometries(country.geometry, ccrs.PlateCarree(), facecolor=(my_color, 1, 0)) x = country.geometry.centroid.x y = country.geometry.centroid.y ax.text(x, y, country.attributes['NAME'], color='red', size=15, ha='center', va='center', transform=ccrs.PlateCarree(), path_effects=[PathEffects.withStroke(linewidth=5, foreground="k", alpha=.6)]) |
Where are the rivers and the lakes?
Probably you have seen the ax.add_feature(cartopy.feature.LAKES, alpha=0.95) , ax.add_feature(cartopy.feature.RIVERS) on the top of the first code and you are wondering, that you could not find any river or a lake on the maps so far. Well, true. As far as you are coloring the whole country, somehow the lakes and the rivers are covered as well. Anyway, it is always a possibility to show them, running this:
|
1 2 3 4 5 6 7 |
ax = plt.axes(projection=ccrs.PlateCarree()) ax.add_feature(cartopy.feature.LAND) ax.add_feature(cartopy.feature.OCEAN) ax.add_feature(cartopy.feature.COASTLINE) ax.add_feature(cartopy.feature.BORDERS, linestyle='-', alpha=.5) ax.add_feature(cartopy.feature.LAKES, alpha=0.3) ax.add_feature(cartopy.feature.RIVERS) |
Coloring a country with two colors, dividing it into North and South
If you want to split a country in two, the centroid could be helpful again. The code below divides Brazil into south and north, based on its centre. Somehow useful, if you need it for some reasons:
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 |
for country in countries: if country.attributes['NAME'] in 'Brazil': l = LineString([(-180, country.geometry.centroid.y), (180, country.geometry.centroid.y)]) north_line = [(-180, country.geometry.centroid.x), (180, country.geometry.centroid.x)] south_line = [(country.geometry.centroid.y, -180), (180, country.geometry.centroid.x)] north_poly = MultiLineString([l, north_line]).convex_hull south_poly = MultiLineString([l, south_line]).convex_hull g = ax.add_geometries([country.geometry.intersection(north_poly)], ccrs.PlateCarree(), facecolor=(.8, 0, 0, 0.8), zorder=99) g = ax.add_geometries([country.geometry.intersection(south_poly)], ccrs.PlateCarree(), facecolor=(0, 0.2, 0.9, 0.2), zorder=99) x = country.geometry.centroid.x y = country.geometry.centroid.y ax.text(x, y, country.attributes['NAME'], color='k', size=16, ha='center', va='center', transform=ccrs.PlateCarree(), path_effects=[PathEffects.withStroke(linewidth=5, foreground="w", alpha=1)], zorder=100) else: ax.add_geometries(country.geometry, ccrs.PlateCarree(), facecolor=(1, 1, 1), label=country.attributes['NAME']) |
Coloring a small country / city
There are some small countries, which are not visible initially on the “admin_0_countries” map. That’s why, there are other maps, in which these are visible. Make sure to adjust the resolution correctly and you would see green points for Singapore, Monaco, Malta and Lichtenstein. Scotland is also visible on the map, and as far as it does not fall in the condition country.geometry.area < 2, it does not have a green point, but is colored completely:
|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 |
ax.set_extent([-150, 60, -25, 60]) shpfilename = shpreader.natural_earth(resolution='10m', category='cultural', name='admin_0_map_units') highlight = ['Singapore', 'Liechtenstein','Monaco', 'Scotland', 'Malta'] reader = shpreader.Reader(shpfilename) countries = reader.records() for country in countries: if country.attributes['NAME'] in highlight: if country.geometry.area < 2: geom = [country.geometry.buffer(2)] else: geom = [country.geometry] g = ax.add_geometries(geom, ccrs.PlateCarree(), facecolor=(0, 0.5, 0, 0.6)) x = country.geometry.centroid.x y = country.geometry.centroid.y ax.text(x, y+5, country.attributes['NAME'], color='red', size=14, ha='center', va='center', transform=ccrs.PlateCarree(), path_effects=[PathEffects.withStroke(linewidth=3, foreground="k", alpha=.8)]) else: ax.add_geometries(country.geometry, ccrs.PlateCarree(), facecolor=(1, 1, 1), label=country.attributes['NAME']) plt.rcParams["figure.figsize"] = (50,50) |
The whole code is available in GitHub as a JupyterNotebook. I hope it works: