source: subversion/applications/rendering/mapnik/generate_tiles.py @ 4083

Last change on this file since 4083 was 3209, checked in by joerg, 13 years ago

start with rendering an overvirew

  • Property svn:executable set to *
File size: 5.0 KB
Line 
1#!/usr/bin/python
2from math import pi,cos,sin,log,exp,atan
3from subprocess import call
4import sys, os
5
6DEG_TO_RAD = pi/180
7RAD_TO_DEG = 180/pi
8
9def minmax (a,b,c):
10    a = max(a,b)
11    a = min(a,c)
12    return a
13
14class GoogleProjection:
15    def __init__(self,levels=18):
16        self.Bc = []
17        self.Cc = []
18        self.zc = []
19        self.Ac = []
20        c = 256
21        for d in range(0,levels):
22            e = c/2;
23            self.Bc.append(c/360.0)
24            self.Cc.append(c/(2 * pi))
25            self.zc.append((e,e))
26            self.Ac.append(c)
27            c *= 2
28               
29    def fromLLtoPixel(self,ll,zoom):
30         d = self.zc[zoom]
31         e = round(d[0] + ll[0] * self.Bc[zoom])
32         f = minmax(sin(DEG_TO_RAD * ll[1]),-0.9999,0.9999)
33         g = round(d[1] + 0.5*log((1+f)/(1-f))*-self.Cc[zoom])
34         return (e,g)
35     
36    def fromPixelToLL(self,px,zoom):
37         e = self.zc[zoom]
38         f = (px[0] - e[0])/self.Bc[zoom]
39         g = (px[1] - e[1])/-self.Cc[zoom]
40         h = RAD_TO_DEG * ( 2 * atan(exp(g)) - 0.5 * pi)
41         return (f,h)
42
43from mapnik import *
44
45def render_tiles(bbox, mapfile, tile_dir, minZoom=1,maxZoom=18, name="unknown"):
46    print "render_tiles(",bbox, mapfile, tile_dir, minZoom,maxZoom, name,")"
47
48    if not os.path.isdir(tile_dir):
49         os.mkdir(tile_dir)
50
51    gprj = GoogleProjection(maxZoom+1) 
52    m = Map(2 * 256,2 * 256)
53    load_map(m,mapfile)
54    prj = Projection("+proj=merc +datum=WGS84")
55   
56    ll0 = (bbox[0],bbox[3])
57    ll1 = (bbox[2],bbox[1])
58   
59    for z in range(minZoom,maxZoom + 1):
60        px0 = gprj.fromLLtoPixel(ll0,z)
61        px1 = gprj.fromLLtoPixel(ll1,z)
62        for x in range(int(px0[0]/256.0),int(px1[0]/256.0)+1):
63            for y in range(int(px0[1]/256.0),int(px1[1]/256.0)+1):
64                p0 = gprj.fromPixelToLL((x * 256.0, (y+1) * 256.0),z)
65                p1 = gprj.fromPixelToLL(((x+1) * 256.0, y * 256.0),z)
66
67                # render a new tile and store it on filesystem
68                c0 = prj.forward(Coord(p0[0],p0[1]))
69                c1 = prj.forward(Coord(p1[0],p1[1]))
70           
71                bbox = Envelope(c0.x,c0.y,c1.x,c1.y)
72                bbox.width(bbox.width() * 2)
73                bbox.height(bbox.height() * 2)
74                m.zoom_to_box(bbox)
75               
76                # check if we have directories in place
77                zoom = "%s" % z
78                str_x = "%s" % x
79                str_y = "%s" % y
80
81                if not os.path.isdir(tile_dir + zoom):
82                    os.mkdir(tile_dir + zoom)
83                if not os.path.isdir(tile_dir + zoom + '/' + str_x):
84                    os.mkdir(tile_dir + zoom + '/' + str_x)
85
86                tile_uri = tile_dir + zoom + '/' + str_x + '/' + str_y + '.png'
87
88                exists= ""
89                if os.path.isfile(tile_uri):
90                    exists= "exists"
91                else:
92                    im = Image(512, 512)
93                    render(m, im)
94                    view = im.view(128,128,256,256) # x,y,width,height
95                    save_to_file(tile_uri,'png',view)
96                    command = "convert  -colors 255 %s %s" % (tile_uri,tile_uri)
97                    call(command, shell=True)
98
99                bytes=os.stat(tile_uri)[6]
100                empty= ''
101                if bytes == 137:
102                    empty = " Empty Tile "
103
104                print name,"[",minZoom,"-",maxZoom,"]: " ,z,x,y,"p:",p0,p1,exists, empty
105
106if __name__ == "__main__":
107    home = os.environ['HOME']
108    mapfile = home + "/svn.openstreetmap.org/applications/rendering/mapnik/osm-local.xml"
109    tile_dir = home + "/osm/tiles/"
110
111    # Start with an overview
112    # World
113    bbox = (-180.0,-90.0, 180.0,90.0)
114    render_tiles(bbox, mapfile, tile_dir, 0, 5,"World")
115
116    minZoom = 10
117    maxZoom = 16
118    bbox = (-2, 50.0,1.0,52.0)
119    render_tiles(bbox, mapfile, tile_dir, minZoom, maxZoom)
120
121
122    # Muenchen
123    bbox = (11.4,48.07, 11.7,48.22)
124    render_tiles(bbox, mapfile, tile_dir, 1, 12 , "Muenchen")
125
126
127    # Muenchen+
128    bbox = (11.3,48.01, 12.15,48.44)
129    render_tiles(bbox, mapfile, tile_dir, 7, 12 , "Muenchen+")
130
131    # Muenchen++
132    bbox = (10.92,47.7, 12.24,48.61)
133    render_tiles(bbox, mapfile, tile_dir, 7, 12 , "Muenchen++")
134
135    # Nuernberg
136    bbox=(10.903198,49.560441,49.633534,11.038085)
137    render_tiles(bbox, mapfile, tile_dir, 10, 16,"Nuernberg")
138
139    # Karlsruhe
140    bbox=(8.179113,48.933617,8.489252,49.081707)
141    render_tiles(bbox, mapfile, tile_dir, 10, 16,"Karlsruhe")
142
143    # Karlsruhe+
144    bbox = (8.3,48.95,8.5,49.05)
145    render_tiles(bbox, mapfile, tile_dir, 1, 16, "Karlsruhe+")
146
147    # Augsburg
148    bbox = (8.3,48.95,8.5,49.05)
149    render_tiles(bbox, mapfile, tile_dir, 1, 16, "Augsburg")
150
151    # Augsburg+
152    bbox=(10.773251,48.369594,10.883834,48.438577)
153    render_tiles(bbox, mapfile, tile_dir, 10, 14, "Augsburg+")
154
155    # Europe+
156    bbox = (1.0,10.0, 20.6,50.0)
157    render_tiles(bbox, mapfile, tile_dir, 1, 11 , "Europe+")
158
159    # World
160    bbox = (-180.0,-90.0, 180.0,90.0)
161    render_tiles(bbox, mapfile, tile_dir, 1, 6,"World")
Note: See TracBrowser for help on using the repository browser.