test various exotic projections with mapnik and osm data

mapnik_proj_test.py

#!/usr/bin/python
import mapnik2 as mapnik

#mapfile = "tests/data/good_maps/bounds_clipping.xml"
mapfile = "/Users/dane/src/osm_mapnik3/osm2.xml"

projections = {
   "latlon": "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs", #EPSG:4326
   "google": "+proj=merc +lon_0=0 +lat_ts=0 +x_0=0 +y_0=0 +a=6378137 +b=6378137 +units=m +no_defs", # SR-ORG:95
   "mercator_world": "+proj=merc +lon_0=0 +k=1 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs", # SR-ORG:16
   "lambert": "+proj=lcc +lat_1=49.83333333333334 +lat_2=51.16666666666666 +lat_0=90 +lon_0=0 +x_0=150000 +y_0=5400000 +ellps=intl +pm=brussels +units=m +no_defs",
   "albert": "+proj=aea +lat_1=29.83333333333334 +lat_2=45.83333333333334 +lat_0=37.5 +lon_0=-96 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs ",
   "pytest": "+proj=lcc +lat_1=45.89893890000052 +lat_2=47.69601440000037 +lat_0=46.8 +lon_0=2.33722917 +x_0=600000 +y_0=200000 +a=6378249.145 +b=6356514.96582849 +pm=2.337229167 +units=m +no_defs",
   '21500':'+proj=lcc +lat_1=49.83333333333334 +lat_2=51.16666666666666 +lat_0=90 +lon_0=0 +x_0=150000 +y_0=5400000 +ellps=intl +pm=brussels +units=m +no_defs'
}

# brussels
box = mapnik.Box2d(4.184933,50.785753,4.582157,50.919701)

imgx = 600
imgy = 400

# Render an image for every projection
for projName, projDef in sorted(projections.items()):

   #print "Projection: " + projName

   # Load the map data
   m = mapnik.Map(imgx, imgy)
   mapnik.load_map(m, mapfile)

   # Override projection defined in osm.xml
   m.srs = projDef

   # Calculate projected boundaries
   prj = mapnik.Projection(projDef)
   wgs84 = mapnik.Projection('+init=epsg:4326')
   tr = mapnik.ProjTransform(wgs84,prj)
   pbox = tr.forward(box)

   # Apply bounding box
   m.zoom_to_box(pbox)

   # Render image
   im = mapnik.Image(imgx, imgy)
   mapnik.render(m, im)
   im.save("img/" + projName + ".png", "png256")