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

#!/usr/bin/python
#

import mapnik2 as mapnik
import sys, os

def addShpLayer(m,shp,style):
   ''' Add a new layer to mapnik map m using shapefile datasource shp and the
   style named style.
   '''
   shp = shp.split('.')[0]  #strip off file extension
   lyr = mapnik.Layer('shp',"+proj=tmerc +lat_0=49 +lon_0=-2 +k=0.999601 +x_0=400000 +y_0=-100000 +ellps=airy +units=m +datum=OSGB36 +no_defs no_defs")
   lyr.datasource = mapnik.Shapefile(file='%s' % shp)
   lyr.styles.append(style)
   m.layers.append(lyr)

if __name__ == "__main__":
   mapfile = "styles_md2.xml"
   map_uri = "image.png"

    #---------------------------------------------------
    #  Change this to the bounding box you want
    #
    #ll = (-6.5, 49.5, 2.1, 59)
   Teesside = (-1.3,54.50,-1.17,54.75)
   Hartlepool = (-1.2,54.62,-1.17,54.72)
   ll = Hartlepool
    #---------------------------------------------------

   z = 10
   imgx = 500 * z
   imgy = 500 * z
   
   m = mapnik.Map(imgx,imgy)
   mapnik.load_map(m,mapfile)
   
   md2root = "/home/graham/OSM/data/Meridian2_Shape/data"

   print "md2map"
   print "Reading Shapefiles"
   # Add Tidal Boundaries
   shp = "%s/%s" % (md2root,"coast_ln_polyline.shp")
   addShpLayer(m,shp,'MHW')

    # Add Natural features (water, woodland and other features)
   shp = "%s/%s" % (md2root,"lake_region.shp")
   addShpLayer(m,shp,'water-poly')
   shp = "%s/%s" % (md2root,"river_polyline.shp")
   addShpLayer(m,shp,'water-line')
   shp = "%s/%s" % (md2root,"woodland_region.shp")
   addShpLayer(m,shp,'woodland-poly')

    # Add Railways
   shp = "%s/%s" % (md2root,"rail_ln_polyline.shp")
   addShpLayer(m,shp,'railway')
   shp = "%s/%s" % (md2root,"station_point.shp")
   addShpLayer(m,shp,'railway')

    # Add Roads
   shp = "%s/%s" % (md2root,"minor_rd_polyline.shp")
   addShpLayer(m,shp,'other-road') 
   shp = "%s/%s" % (md2root,"b_road_polyline.shp")
   addShpLayer(m,shp,'secondary')
   shp = "%s/%s" % (md2root,"a_road_polyline.shp")
   addShpLayer(m,shp,'primary')
   shp = "%s/%s" % (md2root,"motorway_polyline.shp")
   addShpLayer(m,shp,'motorway')

    # Add misc text
   shp = "%s/%s" % (md2root,"text.shp")
   addShpLayer(m,shp,'text')
   
   print "."
   print "rendering map..."

   prj = mapnik.Projection("+proj=tmerc +lat_0=49 +lon_0=-2 +k=0.999601 +x_0=400000 +y_0=-100000 +ellps=airy +units=m +datum=OSGB36 +no_defs no_defs")
   c0 = prj.forward(mapnik.Coord(ll[0],ll[1]))
   c1 = prj.forward(mapnik.Coord(ll[2],ll[3]))
   if hasattr(mapnik,'mapnik_version') and mapnik.mapnik_version() >= 800:
      bbox = mapnik.Box2d(c0.x,c0.y,c1.x,c1.y)
   else:
      bbox = mapnik.Envelope(c0.x,c0.y,c1.x,c1.y)
   m.zoom_to_box(bbox)
   im = mapnik.Image(imgx,imgy)
   mapnik.render(m, im)
   view = im.view(0,0,imgx,imgy) # x,y,width,height
   view.save(map_uri,'png')

   print "Done! - output in %s" % map_uri