A script to automate the download and generation of terrain tiles and contour shapes using GDAL for use in TileMill maps. Note - the contour generation doesn't currently work for merged tiles - trying to fix that!

get_terrain.rb

#!/usr/bin/env ruby
require 'fileutils'

# Based on instructions from:
# http://steveko.wordpress.com/2013/09/11/terrain-in-tilemill-a-walkthrough-for-non-gis-types/
#
# Further reading:
# https://www.mapbox.com/tilemill/docs/guides/terrain-data/
# https://www.mapbox.com/blog/tilemill-raster-colorizer/
# https://www.mapbox.com/tilemill/docs/guides/optimizing-shapefiles/


# Customise these next three lines
# Tile coords taken from: 
# http://srtm.csi.cgiar.org/SELECTION/inputCoord.asp
name        = "uk"
xtiles      = (35..37)
ytiles      = (1..2)

# name        = "alps"
# xtiles      = [38]   # (35..37)
# ytiles      = [3] # (1..2)



z_scale     = 5
contour     = 20
min_scale   = -40
max_scale   = 2000
url_prefix  = "http://srtm.csi.cgiar.org/SRT-ZIP/SRTM_V41/SRTM_Data_GeoTiff/"
cache_path  = "./data/cache/"
tmp_path    = "./data/tmp/"
output_path = "./data/#{name}/"
startTime   = Time.now.to_i

puts "--\nGenerating terrain files with prefix: #{name}"


FileUtils.mkdir_p cache_path
FileUtils.mkdir_p tmp_path
FileUtils.mkdir_p output_path


# Download elevation data
files = []
prefix = ''
xtiles.each do |x|
  x = "#{x}".rjust(2, '0')
  
  ytiles.each do |y|
    y = "#{y}".rjust(2, '0')
    
    prefix   = "srtm_#{x}_#{y}"
    filename = "#{prefix}.zip"
    filepath = File.join(cache_path, filename)
    print "Downloading tile - #{url_prefix}#{filename} : "
    
    if File.exists?(filepath)
      puts "already downloaded!"
    else
      t = Time.now.to_i
      `wget -q -O #{filepath} #{url_prefix}#{filename}`
      puts "done in #{Time.now.to_i - t}s"
    end
    `unzip -o '#{filepath}' -d '#{tmp_path}'`
    files << prefix
  end
end


# Process data for Tilemill
if files.count > 1
  print "Merging tiles: "

  if File.exists?("#{tmp_path}#{name}-src.tif")
    puts "already done!"
  else
    t = Time.now.to_i
    tifs = files.map{|f| "#{tmp_path}#{f}.tif"}.join(" ")
    `gdal_merge.py -o #{tmp_path}#{name}-src.tif #{tifs}`
    puts "done in #{Time.now.to_i - t}s"
  end
else
  `mv #{tmp_path}#{prefix}.tif #{tmp_path}#{name}-src.tif`
end


print "Re-projecting: "
if File.exists?("#{output_path}#{name}.tif")
  puts "already done!"
else
  t = Time.now.to_i
  `gdalwarp -s_srs EPSG:4326 -t_srs EPSG:3785 -r bilinear #{tmp_path}#{name}-src.tif #{output_path}#{name}.tif`
  puts "done in #{Time.now.to_i - t}s"
end


print "Generate hill shading (z-scaling=#{z_scale}): "
if File.exists?("#{output_path}#{name}-hillshade.tif")
  puts "already done!"
else
  t = Time.now.to_i
  `gdaldem hillshade -co compress=lzw -compute_edges -z #{z_scale} #{output_path}#{name}.tif #{output_path}#{name}-hillshade.tif`
  `gdaladdo -r average #{output_path}#{name}-hillshade.tif 2 4 8 16 32`
  puts "done in #{Time.now.to_i - t}s"
end


print "Generating slope files: "
if File.exists?("#{tmp_path}#{name}-slope.tif")
  puts "already done!"
else
  t = Time.now.to_i
  `gdaldem slope #{output_path}#{name}.tif #{tmp_path}#{name}-slope.tif`
  `gdal_translate -ot Byte -scale 0 90 #{tmp_path}#{name}-slope.tif #{output_path}#{name}-slope.tif`
  `gdaladdo -r average #{output_path}#{name}-slope.tif 2 4 8 16 32`
  puts "done in #{Time.now.to_i - t}s"
end


print "Generate contours at #{contour}m intervals (might take a while...): "
if File.exists?("#{output_path}#{name}-contour.shp")
  puts "already done!"
else
  t = Time.now.to_i
  `gdal_contour -a elev -snodata "-32768" -i #{contour} #{tmp_path}#{name}-src.tif #{tmp_path}#{name}-contour.shp`
  puts "done in #{Time.now.to_i - t}s"
  
  print "Reprojecting contours: "
  t = Time.now.to_i
  `ogr2ogr #{output_path}#{name}-contour.shp #{tmp_path}#{name}-contour.shp -t_srs EPSG:900913`
  puts "done in #{Time.now.to_i - t}s"
  
  print "Indexing shapefile: "
  `shapeindex #{output_path}#{name}-contour.shp`
  puts "done in #{Time.now.to_i - t}s"
end


print "Cleaning up tile data: "
`rm -rf #{tmp_path}`
puts "done"


puts "--\nTotal time: #{Time.now.to_i - startTime}s"