chart = map(projection)
import { zip } from "@fil/hello-zip-unzip"
// "https://prod-dcd-datasets-cache-zipfiles.s3.eu-west-1.amazonaws.com/x3zjcdpt7b-1.zip"
files = zip(await FileAttachment("x3zjcdpt7b.zip").arrayBuffer(), {
base64: false,
checkCRC32: true
}).files
rows = {
const text = files[
danseiji
.replace(/^D/, "d")
.replace(" ", "")
.concat(".csv")
].asText();
return d3.csvParseRows(text, d3.autoType);
}
row = rows[0] // get the header
N = row[0] // number of cells
A = row[1] // number of latitude rows
B = row[2] // number of longitude columns
E = row[3] // number of edge points on the map’s perimeter
C = row[4] // number of pixel rows
D = row[5] // number of pixel columns
W = row[6] // width of the map in raw coordinates
H = row[7] // height of the map in raw coordinates
vertices = rows.slice(1, 1 + row[0])
cells = rows.slice(1 + N, 1 + N + A * B) // each cell is described as [shape, ...vertices]
edge = rows.slice(1 + N + A * B, 1 + N + A * B + E).flat()
pixels = rows.slice(1 + N + A * B + E, 1 + N + A * B + E + C * D)
// l,p = lambda, phi (in radians)
danseijiRaw = {
function forward(l, p) {
const X = 0,
Y = 1;
let i = Math.floor(((Math.PI / 2 - p) / Math.PI) * A); // map it to the array
i = Math.max(0, Math.min(A - 1, i)); // coerce it into bounds
let j = Math.floor(((l + Math.PI) / (2 * Math.PI)) * B);
j = Math.max(0, Math.min(B - 1, j));
const [shape, ...v] = cells[i * B + j]; // now we know the planar Vertices with which we're working
if (!vertices) return [NaN, NaN];
const vP = v.map(d => vertices[d]);
const pN = Math.PI / 2 - (i * Math.PI) / A;
const pS = Math.PI / 2 - ((i + 1) * Math.PI) / A;
const yS = i + 1 - ((Math.PI / 2 - p) / Math.PI) * A; // do linear interpolation
let xS = ((l + Math.PI) / (2 * Math.PI)) * B - (j + .5);
xS *= yS * Math.cos(pN) + (1 - yS) * Math.cos(pS); // apply curvature (not strictly necessary, but helps near the poles)
const vSnw = [-.5 * Math.cos(pN), 1],
vSne = [.5 * Math.cos(pN), 1]; // compute the relative spherical vertex positions
const vSsw = [-.5 * Math.cos(pS), 0],
vSse = [.5 * Math.cos(pS), 0];
let triS; // the triangles from which we interpolate
let triP; // the triangles to which we interpolate
if (shape < 0) {
// for negative sloped cells,
triS = [[vSne, vSnw, vSse], [vSsw, vSse, vSnw]];
triP = [[vP[0], vP[1], vP[5]], [vP[3], vP[4], vP[2]]];
} else if (shape > 0) {
// for positive sloped cells,
triS = [[vSse, vSne, vSsw], [vSnw, vSsw, vSne]];
triP = [[vP[5], vP[0], vP[4]], [vP[2], vP[3], vP[1]]];
} else if (i < A / 2) {
// for the northern triangular cells,
triS = [[vSnw, vSsw, vSse]];
triP = [[vP[1], vP[2], vP[3]]];
} else {
// for the southern triangular cells,
triS = [[vSsw, vSne, vSnw]];
triP = [[vP[2], vP[0], vP[1]]];
}
for (let k = 0; k < triS.length; k++) {
const tS = triS[k],
tP = triP[k];
const detT =
(tS[1][Y] - tS[2][Y]) * (tS[2][X] - tS[0][X]) +
(tS[2][X] - tS[1][X]) * (tS[2][Y] - tS[0][Y]); // compute barycentric coordinates on sphere
const c0 =
((tS[1][Y] - tS[2][Y]) * (tS[2][X] - xS) +
(tS[2][X] - tS[1][X]) * (tS[2][Y] - yS)) /
detT;
if (c0 < 0) continue; // unless the point isn't in this triangle in which case you should skip to the next triangle
const c1 =
((tS[2][Y] - tS[0][Y]) * (tS[2][X] - xS) +
(tS[0][X] - tS[2][X]) * (tS[2][Y] - yS)) /
detT;
const c2 = 1 - c0 - c1;
return [
c0 * tP[0][X] + c1 * tP[1][X] + c2 * tP[2][X], // then interpolate into the plane!
c0 * tP[0][Y] + c1 * tP[1][Y] + c2 * tP[2][Y]
];
}
return [NaN, NaN]; // throw `[${xS},${yS}] doesn't seem to be in {${vSne},${vSnw},${vSsw},${vSse}}`;
}
forward.invert = invert;
return forward;
}
function invert(x, y) {
// this linear interpolation is much simpler
let inside = false;
for (let i = 0; i < edge.length; i++) {
const x0 = edge[i][0],
y0 = edge[i][1]; // for each segment of the edge
const x1 = edge[(i + 1) % edge.length][0],
y1 = edge[(i + 1) % edge.length][1];
if (y0 > y != y1 > y)
if (((y - y0) / (y1 - y0)) * (x1 - x0) + x0 > x)
// if the two points fall on either side of a rightward ray from (X,Y)
// and the line between them intersects our ray right of (X,Y)
inside = !inside; // toggle the boolean
}
const i = ((H / 2 - y) / H) * (C - 1);
const i0 = Math.min(Math.floor(i), C - 2);
const cy = i - i0;
const j = ((x + W / 2) / W) * (D - 1);
const j0 = Math.min(Math.floor(j), D - 2);
const cx = j - j0;
let X = 0,
Y = 0,
Z = 0;
for (let di = 0; di <= 1; di++) {
for (let dj = 0; dj <= 1; dj++) {
const weight = (di == 0 ? 1 - cy : cy) * (dj == 0 ? 1 - cx : cx);
const n = Math.max((i0 + di) * D + (j0 + dj), 0),
phiV = pixels[n][0],
lamV = pixels[n][1];
X += weight * Math.cos(phiV) * Math.cos(lamV);
Y += weight * Math.cos(phiV) * Math.sin(lamV);
Z += weight * Math.sin(phiV);
}
}
let phi = Math.atan2(Z, Math.hypot(X, Y)),
lam = Math.atan2(Y, X);
if (!inside) lam += 2 * Math.PI; // signal that this point is outside the normal map, if necessary
return [lam, phi];
}
height = width * (H / W)
projection = d3
.geoProjection(danseijiRaw)
.preclip(
["Danseji I", "Danseji II", "Danseji N"].includes(danseiji)
? d3.geoClipAntimeridian
: d3.geoClipPolygon(sphere)
)
/*
* The sphere serves as a clipping polygon. The difficulty is to take points
* that are a bit towards the inside of the surface: we do it by following
* the edge list of vertices, and taking a point just to the right-hand-side
* of each segment. — @fil
*/
sphere = {
function inv(i) {
return invert(...vertices[i]).map(d => d * (180 / Math.PI));
}
const coordinates = [];
let p0 = inv(edge[0]);
for (let i = edge.length - 1; i >= 0; i--) {
let p1 = inv(edge[i]);
const axis = attitude()
.arc(p1, p0)
.axis();
coordinates.push(
d3.geoInterpolate(d3.geoInterpolate(p0, p1)(0.5), axis)(0.01)
);
p0 = p1;
}
coordinates.push(coordinates[0]);
return {
type: "Polygon",
coordinates: [coordinates]
};
}
d3 = require("d3@7", "d3-geo-polygon@1")
attitude = require("attitude")
import { map } from "@fil/base-map"
Purpose-built for displays of data
Observable is your go-to platform for exploring data and creating expressive data visualizations. Use reactive JavaScript notebooks for prototyping and a collaborative canvas for visual data exploration and dashboard creation.
