Skip to content

Instantly share code, notes, and snippets.

@mpiannucci
Last active March 4, 2022 21:20
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save mpiannucci/680c7c6e58367bb3603d8d3b73ca7f8f to your computer and use it in GitHub Desktop.
Save mpiannucci/680c7c6e58367bb3603d8d3b73ca7f8f to your computer and use it in GitHub Desktop.
D3 polygonize
import { contours } from 'd3-contour';
import * as GeoTIFF from 'geotiff';
import * as fs from 'fs';
const minLng = -95.9791062729049;
const minLat = 29.468568269199807;
const maxLng = -94.87582655116825;
const maxLat = 30.377703838800894;
let vals = undefined;
let width = undefined;
let height = undefined;
// Rotate a GeoTIFF’s longitude from [0, 360] to [-180, +180].
function rotate(values, width, height) {
var l = width >> 1;
for (var j = 0, k = 0; j < height; ++j, k += width) {
values.subarray(k, k + l).reverse();
values.subarray(k + l, k + width).reverse();
values.subarray(k, k + width).reverse();
}
return values;
}
function linspace(start, stop, num, endpoint = true) {
const div = endpoint ? (num - 1) : num;
const step = (stop - start) / div;
return Array.from({length: num}, (_, i) => start + step * i);
}
try {
const tiff = await GeoTIFF.fromFile('../CityofHoustonROM_100yr.tif');
const image = await tiff.getImage();
height = image.getHeight();
width = image.getWidth();
const raster = await image.readRasters();
//vals = rotate(raster[0], width, height);
vals = raster[0];
} catch (err) {
console.error(err)
}
const stops = linspace(-1, 30, 250);
console.log(stops);
const polygons = contours()
.size([height, width])
.thresholds(stops)
(vals)
.map(p => ({
...p,
coordinates: p.coordinates
.map(ring =>
ring.map(coords =>
coords
.map(point => ([
minLng + (maxLng - minLng) * (point[0] / height),
maxLat - (maxLat - minLat) * (point[1] / width)
]))
.map(point => ([
point[0] > 180 ? point[0] - 360 : point[0],
point[1] > 90 ? point[1] - 90 : point[1],
]))
)
)
}));
const features = polygons.map(p => ({
type: 'Feature',
properties: {
value: p.value,
},
geometry: {
type: 'MultiPolygon',
coordinates: p.coordinates,
},
}));
const featureCollection = {
type: 'FeatureCollection',
features: features,
};
const featureCollectionString = JSON.stringify(featureCollection);
fs.writeFileSync('../houston_max_depth_d3.json', featureCollectionString);
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment