I have a Python script, png_from_osm.py, that I’ve used a few times to make paper maps from osm data. I think I first got the script from here, but I’ve modified it since then (I’ve included it below). Usage looks like this:
png_from_osm.py
Usage: png_from_osm.py xmin xmax ymin ymax osm_zoom_level output_file
xmin xmax ymin ymax are the osm tile indices corresponding to osm_zoom_level. It can be useful to use this website to find the tiles of the area that you’re interested in. Below is an example url with the xmin xmax ymin ymax info in it. “Scale” doesn’t seem to do anything, but zoom level and tile indices are important (both at this website and for the png_from_osm.py script).
png_from_osm.py uses osm’s pre-rendered tiles to make the map. These tiles are fully rendered at various zoom scales. For zoom parameter this website may be useful:
http://wiki.openstreetmap.org/wiki/Zoom_levels
and some of my observations:
- 10m contours exist starting at zoom 13 and go till at least zoom 18.
- Village roads (at least the yayla roads in the Kaçkar) start at zoom 13.
- Contour labels unfortunately don’t start until zoom level 15.
- I printed out zoom level 14 at 1:40K, and it’s good for a walking map.
Once you know the tile indices and zoom level you’d like, you can use that script to make a map (png). From there if you know what scale you’d like, you can get your meter_per_pixel_calc.py by using the script of the same name. Then you can calculate dpi or pixels/cm for your printer. Mulitply that by the number of pixels to figure out how big your map will be.
Here’s an example:
png_from_osm.py 10050 10077 6149 6164 14 map14-160809-0753.png
By looking at the German url example (in a web browser; it’s the same tile indices as used above), we learn that that map is centered around 40.764 latitude. The meter_per_pixel_calc.py script shows us that there are approximately 7.237 meters on the ground for every pixel in the png (at the given latitude and osm zoom level):
meter_per_pixel_calc.py 40.764 14
7.23675372446
If we want to print the map at 1:40000 (i.e. 1 cm on the map will be 400 meters on ground), then we need to print the file at 55 pixels/cm:
400 / 7.236
55.2791597567717
Here are the scripts:
png_from_osm.py
#!/usr/bin/python
import io, urllib2, datetime, time, re, random
import sys
import os
from PIL import Image, ImageDraw
def Usage():
print "Usage: png_from_osm.py xmin xmax ymin ymax osm_zoom_level output_file"
print " This routine will create a png file from osm data. The"
print " opencyclemap tiles are hard-coded. xmin, xmax, ymin, ymax are"
print " osm tile indices for the given osm_zoom_level"
return
if len(sys.argv) != 7:
Usage()
sys.exit()
#parse command line arguments
xmin = None
xmax = None
ymin = None
ymax = None
zoom = None
outputFileName = None
i = 1
while (i < len(sys.argv)):
arg = sys.argv[i]
if (xmin is None):
try:
xmin = float(arg)
except:
print "xmin must be a number."
sys.exit()
xmin = int(arg)
if (xmin != float(arg)):
print "xmin must be an integer."
sys.exit()
elif (xmax is None):
try:
xmax = float(arg)
except:
print "xmax must be a number."
sys.exit()
xmax = int(arg)
if (xmax != float(arg)):
print "xmax must be an integer."
sys.exit()
elif (ymin is None):
try:
ymin = float(arg)
except:
print "ymin must be a number."
sys.exit()
ymin = int(arg)
if (ymin != float(arg)):
print "ymin must be an integer."
sys.exit()
elif (ymax is None):
try:
ymax = float(arg)
except:
print "ymax must be a number."
sys.exit()
ymax = int(arg)
if (ymax != float(arg)):
print "ymax must be an integer."
sys.exit()
elif (zoom is None):
try:
zoom = float(arg)
except:
print "osm_zoom_level must be a number."
sys.exit()
zoom = int(arg)
if (zoom != float(arg)):
print "osm_zoon_level must be an integer."
sys.exit()
elif (outputFileName is None):
outputFileName = arg
if (os.path.isfile(outputFileName)):
print "%s exists." % outputFileName
sys.exit()
i += 1
if (outputFileName is None):
Usage()
sys.exit()
#I now (June 2017) need an API key to get opencycltmap tiles without a
#watermark
#I got this key from thunderforest.com
APIkey = 'you need to get your own API key'
sURLAPI = '?apikey=%s' % APIkey
#for testing
#sURLAPI = ''
#layers = ["http://{abc}.tile.opencyclemap.org/cycle/!z/!x/!y.png%s" % sURLAPI]
#There's a new server to pull these tiles from to use the api key
layers = ["https://{abc}.tile.thunderforest.com/cycle/!z/!x/!y.png%s" % sURLAPI]
attribution = 'Map data (c) OpenStreetMap, Tiles (c) Andy Allan, compilation Bryan Keith'
xsize = xmax - xmin + 1
ysize = ymax - ymin + 1
resultImage = Image.new("RGBA", (xsize * 256, ysize * 256), (0,0,0,0))
counter = 0
for x in range(xmin, xmax+1):
for y in range(ymin, ymax+1):
for layer in layers:
url = layer.replace("!x", str(x)).replace("!y", str(y)).replace("!z", str(zoom))
match = re.search("{([a-z0-9]+)}", url)
if match:
url = url.replace(match.group(0), random.choice(match.group(1)))
print url, "... ";
try:
#I tried different combinations to get the api key to work
#req = urllib2.Request(url, headers={'User-Agent': 'BigMap/2.0'})
#req = urllib2.Request(url, headers={'apikey': APIkey})
req = urllib2.Request(url)
#req.add_header("Authorization", "Basic %s" % APIkey)
#req.add_header("auth_appkey", APIkey)
#req.add_header("apikey", APIkey)
tile = urllib2.urlopen(req).read()
except Exception, e:
print "Error", e
continue;
image = Image.open(io.BytesIO(tile))
resultImage.paste(image, ((x-xmin)*256, (y-ymin)*256), image.convert("RGBA"))
counter += 1
if counter == 10:
time.sleep(2);
counter = 0
draw = ImageDraw.Draw(resultImage)
draw.text((5, ysize*256-15), attribution, (0,0,0))
del draw
now = datetime.datetime.now()
resultImage.save(outputFileName)meter_per_degree_latitude.py
#!/usr/bin/python
#approximately calculate the meters/degree of longitude at the specified
#latitude
#see https://knowledge.safe.com/articles/725/calculating-accurate-length-in-meters-for-lat-long.html
import math
import sys
def Usage():
print "Usage: meter_per_degree_latiude.py latitude_in_degrees"
print " This routine will write out the approximate number of meters per"
print " degree of longitude for the given latitude."
print " meters are hard-coded"
return
if len(sys.argv) != 2:
Usage()
sys.exit()
#parse command line arguments
lat = None
i = 1
while (i < len(sys.argv)):
arg = sys.argv[i]
if (lat is None):
try:
lat = float(arg)
except:
print "latitude_in_degrees must be a number."
sys.exit()
i += 1
#convert latitude degrees to radians
pi = math.pi
rlat = lat * pi / 180
#calculate meters/degree of longitude
mlon = 111412.84 * math.cos(rlat) - 93.5 * math.cos(3 * rlat)
#calculate meters/degree of latitude
mlat = 111132.92 - 559.82 * math.cos(2 * rlat) + 1.175 * math.cos(4 * rlat)
print "%f meters/degree of longitude at %f degrees latitude" % (mlon, lat)
print "%f meters/degree of latitude at %f degrees latitude" % (mlat, lat)
meter_per_pixel_calc.py
#!/usr/bin/python
#calculate the distance/pixel of osm tiles at your desired tile zoom and
#latitude
#see http://wiki.openstreetmap.org/wiki/Zoom_levels
import math
import sys
def Usage():
print "Usage: meter_per_pixel_calc.py latitude_in_degrees osm_zoom_level"
print " This routine will write out the approximate numbers of meters per"
print " pixel for the given latitude and openstreetmap tiled zoom level."
print " meters are hard-coded"
return
if len(sys.argv) != 3:
Usage()
sys.exit()
#parse command line arguments
lat = None
z = None
i = 1
while (i < len(sys.argv)):
arg = sys.argv[i]
if (lat is None):
try:
lat = float(arg)
except:
print "latitude_in_degrees must be a number."
sys.exit()
elif (z is None):
try:
z = float(arg)
except:
print "osm_zoom_level must be an integer."
sys.exit()
z = int(arg)
if (z != float(arg)):
print "osm_zoom_level must be an integer."
sys.exit()
i += 1
if (z is None):
Usage()
sys.exit()
#set the openstreetmap tile zoom level you're interested in
#z = 14
#set the circumference of the earth at the equator in meters (or whatever
#unit you want per pixel)
C = 40075.16 * 1000
#set the latitude that we're interested in (in degrees)
#lat = 0
#convert latitude degrees to radians
pi = math.pi
y = lat * pi / 180
#calculate distance/pixel
S = C * math.cos(y) / 2 ** (z + 8)
print S