Protoyping a Matplotlib/Agg Engine for PCL

Matplotlib (http://matplotlib.sourceforge.net), a 2D plotting library, is attracting well deserved attention from Python users. It is a readily installed and rich environment for 2D plotting and visualization with a fat users guide and many nice examples. There even exists a cartographic module to generate background basemaps for plots of geophysical data. What intrigues me the most is its Agg backend and API that provides an easy way to learn about anti-grain geometry. It seemed like it wouldn't take long to find out if matplotlib and Agg were a good fit to PCL's mapping classes.

I began by revisiting the two slightly overlapping triangle polygons from a previous entry, then creating a polygon symbolizer, and defining map image output parameters. PCL polygon symbolizers have a default opacity of 50% -- we'll see the overlap effect in the output image:

from cartography import mapping, spatial
from cartography.spatial import Point, LinearRing, Polygon

t1wkt = 'POLYGON ((-1.0 50.5, -0.5 51.2, 0.3 50.9, -1 50.5))'
triangle1 = spatial.Polygon(wkt=t1wkt)

t2wkt = 'POLYGON((-0.7 50.3, 0.1 51.0, 0.6 50.1, -0.7 50.3))'
triangle2 = spatial.Polygon(wkt=t2wkt)

# Symbolizer for triangles
symb = mapping.PolygonSymbolizer(fill={'fill': '#FF8080'})

# Define output map parameters
epsg4326 = spatial.SpatialReference(epsg=4326)
view = mapping.View(epsg4326, spatial.Point(-1.2,51.25),
                    spatial.Point(1,50.05))

width = 550
height = 300

Define two utility functions: hex_to_rgb() converts PCL colors to matplotlib's normalized tuple, and transform_point() maps PCL points in world units to matplotlib image units. Image coordinates in matplotlib have their origin at the lower left corner, contrary to convention:

def hex_to_rgb(color):
    """convert a hex color such as #FFFFFF to a red, green, blue tuple"""
    c = eval('0x' + color[1:])
    r = (c >> 16) & 0xFF
    g = (c >> 8) & 0xFF
    b = c & 0xFF
    return (r, g, b)

def transform_point(point, view, width, height):
    """transform from world to image coordinates"""
    dx = (view.lr.x - view.ul.x)/width
    dy = (view.lr.y - view.ul.y)/height
    p = (point.x - view.ul.x)/dx
    l = (point.y - view.lr.y)/dy
    return (p, l)

Import the matplotlib module and render the two triangles. The matplotlib graphic context 'gc' is configured from properties of the PCL polygon symbolizer:

import matplotlib
matplotlib.use('Agg')
from matplotlib.backends.backend_agg import RendererAgg
from matplotlib.transforms import Value

dpi = Value(72.0)
o = RendererAgg(w, h, dpi)
gc = o.new_gc()

for polygon in [triangle1, triangle2]:
    points = [transform_point(p, view, width, height) for p in polygon[0]]

    # Fill
    gc.set_alpha(symb.fill['fill-opacity'])
    rgb = hex_to_rgb(symb.fill['fill'])
    face = (rgb[0]/255.0, rgb[1]/255.0, rgb[2]/255.0)
    o.draw_polygon(gc, face, points)

    # Stroke
    gc.set_alpha(symb.stroke['stroke-opacity'])
    gc.set_foreground(symb.stroke['stroke'])
    gc.set_linewidth(symb.stroke['stroke-width'])
    o.draw_polygon(gc, None, points)

o._renderer.write_png('agg_rendered.png')

the result, agg_rendered.png:

My impressions? So far, I like what I see. The matplotlib API is better suited to PCL's SLD-ish mapping objects than MapServer's mapscript module. That's a big plus for me. I don't yet appreciate why matplotlib excludes shape fill colors from the graphic context, or the quirky choice of lower left for image coordinate origin.

Clearly, the code within the loop over triangles needs to be written in C++ (like Agg) rather than Python to get decent performance when rendering hundreds or thousands of polygons at a time. I haven't written any C++ extensions for Python yet (C only), but this might be my first opportunity.