GeoTIFF Contours II

An example of d3-contour, reprojecting contours computed from a GeoTIFF.

forked from mbostock's block: GeoTIFF Contours II

<!DOCTYPE html>
<svg width="960" height="500"></svg>
<script src="https://unpkg.com/d3-array@1"></script>
<script src="https://unpkg.com/d3-contour@1"></script>
<script src="https://unpkg.com/d3-collection@1"></script>
<script src="https://unpkg.com/d3-color@1"></script>
<script src="https://unpkg.com/d3-dispatch@1"></script>
<script src="https://unpkg.com/d3-geo@1"></script>
<script src="https://unpkg.com/d3-geo-projection@2"></script>
<script src="https://unpkg.com/d3-interpolate@1"></script>
<script src="https://unpkg.com/d3-request@1"></script>
<script src="https://unpkg.com/d3-selection@1"></script>
<script src="https://unpkg.com/d3-scale@1"></script>
<script src="https://unpkg.com/geotiff@0.4/dist/geotiff.browserify.min.js"></script>
<script>

d3.request("sfctmp.tiff").responseType("arraybuffer").get(function(error, request) {
  if (error) throw error;

var tiff = GeoTIFF.parse(request.response),
      image = tiff.getImage(),
      values = image.readRasters()[0],
      m = image.getHeight(),
      n = image.getWidth();

var color = d3.scaleSequential(d3.interpolateMagma)
      .domain(d3.extent(values));

var projection = d3.geoNaturalEarth()
      .rotate([-180, 0])
      .precision(0.1);

var path = d3.geoPath(projection);

var contours = d3.contours()
      .size([n, m]);

d3.select("svg")
      .attr("stroke", "#000")
      .attr("stroke-width", 0.5)
      .attr("stroke-linejoin", "round")
    .selectAll("path")
    .data(contours(values).map(invert))
    .enter().append("path")
      .attr("fill", function(d) { return color(d.value); })
      .attr("d", path);

// Invert the pixel coordinates to [longitude, latitude]. This assumes the
  // source GeoTIFF is in equirectangular coordinates. This particular GeoTIFF
  // also appears to be rotated by 180° longitude, so the antimeridian cut is
  // actually at the meridian (0°).
  //
  // Note that inverting the projection breaks the polygon ring associations:
  // holes are no longer inside their exterior rings. Fortunately, since the
  // winding order of the rings is consistent and we’re now in spherical
  // coordinates, we can just merge everything into a single polygon!
  function invert(d) {
    var shared = {};

// Record the y-intersections with the antimeridian.
    p.coordinates.forEach(function(ring) {
      ring.forEach(function(p) {
        if (p[0] === -180 || p[0] === 180) {
          shared[p[1]] |= p[0] === -180 ? 1 : 2;
        }
      });
    });

// Offset any unshared antimeridian points to prevent their stitching.
    p.coordinates.forEach(function(ring) {
      ring.forEach(function(p) {
        if ((p[0] === -180 || p[0] === 180) && shared[p[1]] !== 3) {
          p[0] = p[0] === -180 ? -179.9995 : 179.9995;
        }
      });
    });

p = d3.geoStitch(p);
    if (!p.coordinates.length) p = {type: "Sphere"}; // TODO fix d3.geoStitch
    p.value = d.value;
    return p;
  }
});

</script>

https://unpkg.com/d3-array@1

https://unpkg.com/d3-contour@1

https://unpkg.com/d3-collection@1

https://unpkg.com/d3-color@1

https://unpkg.com/d3-dispatch@1

https://unpkg.com/d3-geo@1

https://unpkg.com/d3-geo-projection@2

https://unpkg.com/d3-interpolate@1

https://unpkg.com/d3-request@1

https://unpkg.com/d3-selection@1

https://unpkg.com/d3-scale@1

https://unpkg.com/geotiff@0.4/dist/geotiff.browserify.min.js