block by vlandham 8075631

8075631

Full Screen

Shaded Relief

This is an example of how to create a shaded relief raster with a vector map overlay (using SVG and d3.js).

Step 1 was to create the raster. I used tiled GeoTiffs from the SRTM project, downloading four tiles that completed a map of Costa Rica. To combine the tiff files into a single raster with the correct projection and dimensions, I used gdalwarp:

gdalwarp \
  -r lanczos \
  -te -250000 -156250 250000 156250 \
  -t_srs "+proj=aea +lat_1=8 +lat_2=11.5 +lat_0=9.7 +lon_0=-84.2 +x_0=0 +y_0=0" \
  -ts 960 0 \
  srtm_19_10.tif \
  srtm_20_10.tif \
  srtm_19_11.tif \
  srtm_20_11.tif \
  relief.tiff

The t_srs option sets an albers equal area projection that will center on Costa Rica. The te option defines the extent of the map, using SRS coordinates. I don’t fully understand how this works and used some trial and error. Note that the x/y has a ratio of 1.6, the same as the intended output resolution (960x600).

Note that the projection here mirrors the projection set in index.html.

Step 2 is to create, from this GeoTiff file, two images: one, grayscale, that represents “shade” — given a certain direction of sunlight, it simulates the effect of light on the relief map:

gdaldem \
  hillshade \
  relief.tiff \
  hill-relief-shaded.tiff \
  -z 4 -az 20

The second image is a “color relief” that maps certain colors to certain elevations. The color_relief.txt file provides this information in the format: elevation r g b.

gdaldem \
  color-relief \
  relief.tiff \
  color_relief.txt \
  hill-relief-c.tiff \

These files are combined using the program hsv_merge.py:

hsv_merge.py \
  hill-relief-c.tiff \
  hill-relief-shaded.tiff \
  hill-relief-merged.tiff

index.html

<!DOCTYPE html>
<meta charset='utf-8'>
<style>
  /* Style here... */
  path#CRI {
    fill: none;
    stroke: #000;
  }
  image.bg {
    opacity: 0.2;
  }
</style>

<div></div>

<script src="//d3js.org/d3.v3.min.js"></script>
<script src="//d3js.org/topojson.v0.min.js"></script>
<script>
  var width = 960,
      height = 600;

  var projection = d3.geo.albers()
    .center([0, 9.7])
    .rotate([84.2,0])
    .parallels([8, 11.5])
    //.scale(10200)
    .scale(12240)
    .translate([width / 2, height / 2]);

  var path = d3.geo.path()
    .projection(projection);

  var svg = d3.select("body").append("svg")
    .attr("class", "map")
    .attr("width", width)
    .attr("height", height);

  d3.json("costarica_min_topo.json", function(error, data) {
    var costarica = topojson.object(data, data.objects.costarica);

    svg.append("image")
      //.attr("clip-path", "url(#clip)")
      .attr("xlink:href", "hill-relief.jpg")
      .attr("width", width)
      .attr("height", height)
      .attr("class", "bg");

    svg.selectAll(".cr-subunit")
      .data(costarica.geometries)
    .enter().append("path")
      .attr("id", function(d) { return "CRI"; })
      .attr("d", path);

    svg.append("clipPath")
      .attr("id", "clip")
    .append("use")
      .attr("xlink:href", "#CRI");

    svg.append("image")
      .attr("clip-path", "url(#clip)")
      .attr("xlink:href", "hill-relief.jpg")
      .attr("width", width)
      .attr("height", height)
      .attr("class", "fg");
  });

  // Script here...
</script>

Makefile

all: hill-relief.jpg costarica_min_topo.json

# -------------
# Relief raster
# -------------
#
# Notice the `zip` file requirements here have no download.
# You will need to search for them online. They are from the
# SRTM project: http://www2.jpl.nasa.gov/srtm/
# (which appears to have multiple versions of files).

srtm_20_11.tif: srtm_20_11.zip
	unzip srtm_20_11.zip
	rm readme.txt
	touch srtm_20_11.tif

srtm_20_10.tif: srtm_20_10.zip
	unzip srtm_20_10.zip
	rm readme.txt
	touch srtm_20_10.tif

srtm_19_11.tif: srtm_19_11.zip
	unzip srtm_19_11.zip
	rm readme.txt
	touch srtm_19_11.tif

srtm_19_10.tif: srtm_19_10.zip
	unzip srtm_19_10.zip
	rm readme.txt
	touch srtm_19_10.tif

# Generate, from the four tiles, a relief GeoTiff file with the correct
# projection. Right now I'm finding the correct `-te` values by trial and error.
# -te -300000 -250000 300000 250000
relief.tiff: srtm_19_10.tif srtm_20_10.tif srtm_20_11.tif srtm_19_11.tif
	gdalwarp \
		-r lanczos \
		-te -250000 -156250 250000 156250 \
		-t_srs "+proj=aea +lat_1=8 +lat_2=11.5 +lat_0=9.7 +lon_0=-84.2 +x_0=0 +y_0=0" \
		-ts 960 0 \
		srtm_19_10.tif \
		srtm_20_10.tif \
		srtm_19_11.tif \
		srtm_20_11.tif \
		relief.tiff

# Shaded relief.
hill-relief-shaded.tiff: relief.tiff
	gdaldem \
		hillshade \
		relief.tiff \
		hill-relief-shaded.tiff \
		-z 4 -az 20

# Colorized elevation.
hill-relief-c.tiff: relief.tiff
	gdaldem \
		color-relief \
		relief.tiff \
		color_relief.txt \
		hill-relief-c.tiff \

# Merge the colorized and shaded rasters. Requires `hsv_merge.py`:
# http://fwarmerdam.blogspot.com/2010/01/hsvmergepy.html
hill-relief-merged.tiff: hill-relief-c.tiff hill-relief-shaded.tiff
	~/bin/hsv_merge.py \
		hill-relief-c.tiff \
		hill-relief-shaded.tiff \
		hill-relief-merged.tiff

hill-relief.jpg: hill-relief-merged.tiff
	convert hill-relief-merged.tiff hill-relief.jpg

# --------------------
# Costa Rica Geography
# --------------------

# Data from http://gadm.org/
CRI_adm.zip:
	curl -o CRI_adm.zip http://gadm.org/data/shp/CRI_adm.zip

CRI_adm0.shp: CRI_adm.zip
	unzip CRI_adm.zip
	touch CRI_adm0.shp

costarica.json: CRI_adm0.shp
	ogr2ogr \
		-f GeoJSON \
		costarica.json \
		CRI_adm0.shp

# Require topojson: https://github.com/mbostock/topojson
# (this minifies/simplifies the data)
costarica_min_topo.json: costarica.json
	topojson \
		-p name=NAME \
		-p name \
		-q 1e4 \
		-o costarica_min_topo.json \
		costarica.json

.PHONY: clean

# Clean does not remove downloaded files.
clean:
	rm -f hill-relief*
	rm -f relief*
	rm -f srtm_19_10.hdr
	rm -f srtm_19_10.t*
	rm -f srtm_19_11.hdr
	rm -f srtm_19_11.t*
	rm -f srtm_20_10.hdr
	rm -f srtm_20_10.t*
	rm -f srtm_20_11.hdr
	rm -f srtm_20_11.t*
	rm -f CRI_adm0*
	rm -f CRI_adm1*
	rm -f CRI_adm2*
	rm -f CRI_readme.txt
	rm -f costarica.json
	rm -f costarica_min_topo.json

color_relief.txt

65535 255 255 255
5800 254 254 254
3000 121 117 10
1500 151 106 47
800 127 166 122
500 213 213 149
1 201 213 166