A Lazy Sequence

SVG maps with Cairo and Mapnik

Mapping technology has become a commodity resource in the last few years, and supporting it is the wikipedia of mapping data, the Open Street Maps project. If you need something that google’s mapping technology doesn’t provide, it is very likely that the ecosystem surrounding OSM might have answers for you.

In this post I am going to demonstrate how to use the mapniki toolkit to generate SVG or PDF maps with custom overlays using the cairo library. Mapnik suits my needs well as it exports a straightforward Python API that integrates cleanly with cairo. The code in this post is based on these snippets. I followed this tutorialii to get a mapnik OSM environment running in a linux virtual machineiii on my mac. All the code for this example is in python. Prior to following this code, I suggest you run the generate_xml.py file from osm2psql as suggested in the tutorial to get an osm.xmliv file for rendering your maps.

What we are going to do here is render a map of the Central Business District of Christchurch, New Zealand and then use Cairo to overlay a vector graphic marker on the Cathedral. In a real world project you would be overlaying more sophisticated graphics; in my case I am rendering SVG objects onto the map. Note that if are just wanting to render raster graphics, or simple line art, you can offload all of that to mapnik to draw.

The first thing to do is to import the libraries we will be using, and define some constants:

import mapnik
import cairo

mapfile = "osm.xml"
svg_filename = "chch.svg"

# We are going to draw the CBD of Christchurch, New Zealand 
cathedral_ll =  (172.637014, -43.530909) # location of the cathedral, 
                                         # approximately the center of the CBD
cbd_lls = ((172.626758, -43.525931), (172.644954, -43.534954)) # boundary of the CBD

If you have done any GIS work in the past you might have noticed that these coordinates are lon, lat pairs rather than lat, lon; This is intentional and required for consistency with the rest of mapnik.

Next we create the map object, and then our projection. I’m just using a projection from an example file, if you need more info you should consult the mapnik wiki but you can safely ignore the details for now:

m = mapnik.Map(1000, 1300)
mapnik.load_map(m, mapfile)
prj = mapnik.Projection("+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +no_defs +over")

Next we tell mapnik what the minimum bounding view should be for our output. Mapnik will take care of fitting our box to the size of the image we have specified. To do this we create mapnik Coord objects and then use the projection we have created to project them from spherical coordinates into the same projection our map is working from and then tell the map to zoom to fit that box:

bbox = mapnik.Envelope(mapnik.Coord(*cbd_lls[0]), mapnik.Coord(*cbd_lls[1])).forward(prj)
m.zoom_to_box(bbox)

Now that we have have a bounding box that will show the CBD, we need the location of the cathedral in the same projection:

cathedral = mapnik.Coord(*cathedral_ll)
cathedral_projected = prj.forward(cathedral)

There is one last piece of view math we need to take care of before we get to rendering. Mapnik is going to take care of drawing our map and worrying about the view transformation from the projection’s abstract cartesian space to the drawing surface’s concrete cartesian space. However, because we are drawing the cathedral’s location ourselves we need to get the projected cathedral location transformed into view coordinates. Turns out this is really easy too:

cathedral_view_pos = m.view_transform().forward(cathedral_projected)

Lastly we have drawing. This is a two phase operation; first draw the map, second draw our custom graphics. One very important step here that you might miss if you follow some of the snippets on the mapnik wiki is that you want mapnik to render to the cairo surface’s context, rather than directly to the surface. Render is overloaded to handle either type. I found that if I passed the surface and then tried to render anything at all, my map would appear white.

f = open(svg_filename, 'w')
surface = cairo.SVGSurface(f.name, m.width, m.height)

ctx = cairo.Context(surface)

mapnik.render(m, ctx) # note that we render to the context, not the surface!

Lastly we need to render our custom graphic. This is going to be a simple black cross:

ctx.move_to (cathedral_view_pos.x - 5, cathedral_view_pos.y - 5)
ctx.line_to (cathedral_view_pos.x + 5, cathedral_view_pos.y + 5)
ctx.move_to (cathedral_view_pos.x - 5, cathedral_view_pos.y + 5)
ctx.line_to (cathedral_view_pos.x + 5, cathedral_view_pos.y - 5)
ctx.close_path()
ctx.stroke()

surface.finish()

And thats it. If you run this program with the required osm.xml file, it should generate chch.svg. As you can see, this turns out to be a very simple process. I strongly encourage you to read up on the the python apis to both mapnik and cairo.

Thanks to springmeyer in #mapnik on irc.freenode.org for helping me get this going.

  1. The other major contender for SVG generation is an XSLT package called OSMARender.
  2. For the username specified for the mapping data, instead of creating a new db user I just used my login to save some mucking about. I also opted to ignore the details about shared memory and related optimizations for my dev vm.
  3. Running in a headless Virtual Box environment.
  4. I used a small subset of the full world data that only includes the greater christchurch area. This is only ~30mb rather than ~9gb!
15 December 2010