Projections (or AI)

challenges/06-raster/README.md

@@ -0,0 +1,16 @@
+# Raster data: rendering 2d raster map in 3d
+
+Weather models cannot 'see' the city as we do. The grid cells are too coarse to
+represent individual buildings explicitly. Instead, they see the city at block
+level.
+
+Local climate zones are archetypes of city blocks. They are used to describe the
+important characteristic of a city succinctly. LCZs can be displayed on 2d maps;
+however, that doesn't really appeal to ones imagination.
+
+I figured, why not render these nice prototype tiles in 3d? To this end, I tried
+to work with a file format that was new to me, called CityJSON. This is how far
+I got in limited time. There is still lots of room for improvements, but I think
+the concept is already clear and quite nice.
+
+![Image for mapchallenge](LCZ_CityJSON.png)

challenges/26-projections/gingery.js

@@ -0,0 +1,182 @@
+/** Modified from https://github.com/d3/d3-geo-projection/blob/main/src/gingery.js */
+
+import * as d3 from "https://cdn.jsdelivr.net/npm/d3@7/+esm";
+import {
+  abs,
+  asin,
+  atan2,
+  cos,
+  degrees,
+  epsilon,
+  epsilon2,
+  halfPi,
+  pi,
+  radians,
+  round,
+  sin,
+  sqrt,
+} from "./math.js";
+
+export function gingeryRaw(rho, n) {
+  var k = (2 * pi) / n,
+    rho2 = rho * rho;
+
+  function forward(lambda, phi) {
+    const lambdaTwisted = lambda + 0.3 * pi * phi;
+    var p = d3.geoAzimuthalEquidistantRaw(lambdaTwisted, phi),
+      x = p[0],
+      y = p[1],
+      r2 = x * x + y * y;
+
+    if (r2 > rho2) {
+      var r = sqrt(r2),
+        theta = atan2(y, x),
+        theta0 = k * round(theta / k),
+        alpha = theta - theta0,
+        rhoCosAlpha = rho * cos(alpha),
+        k_ =
+          (rho * sin(alpha) - alpha * sin(rhoCosAlpha)) /
+          (halfPi - rhoCosAlpha),
+        s_ = gingeryLength(alpha, k_),
+        e = (pi - rho) / gingeryIntegrate(s_, rhoCosAlpha, pi);
+
+      x = r;
+      var i = 50,
+        delta;
+      do {
+        x -= delta =
+          (rho + gingeryIntegrate(s_, rhoCosAlpha, x) * e - r) / (s_(x) * e);
+      } while (abs(delta) > epsilon && --i > 0);
+
+      y = alpha * sin(x);
+      if (x < halfPi) y -= k_ * (x - halfPi);
+
+      var s = sin(theta0),
+        c = cos(theta0);
+      p[0] = x * c - y * s;
+      p[1] = x * s + y * c;
+    }
+    return p;
+  }
+
+  //   forward.invert = function (x, y) {
+  //     var r2 = x * x + y * y;
+  //     if (r2 > rho2) {
+  //       var r = sqrt(r2),
+  //         theta = atan2(y, x),
+  //         theta0 = k * round(theta / k),
+  //         dTheta = theta - theta0;
+
+  //       x = r * cos(dTheta);
+  //       y = r * sin(dTheta);
+
+  //       var x_halfPi = x - halfPi,
+  //         sinx = sin(x),
+  //         alpha = y / sinx,
+  //         delta = x < halfPi ? Infinity : 0,
+  //         i = 10;
+
+  //       while (true) {
+  //         var rhosinAlpha = rho * sin(alpha),
+  //           rhoCosAlpha = rho * cos(alpha),
+  //           sinRhoCosAlpha = sin(rhoCosAlpha),
+  //           halfPi_RhoCosAlpha = halfPi - rhoCosAlpha,
+  //           k_ = (rhosinAlpha - alpha * sinRhoCosAlpha) / halfPi_RhoCosAlpha,
+  //           s_ = gingeryLength(alpha, k_);
+
+  //         if (abs(delta) < epsilon2 || !--i) break;
+
+  //         alpha -= delta =
+  //           (alpha * sinx - k_ * x_halfPi - y) /
+  //           (sinx -
+  //             (x_halfPi *
+  //               2 *
+  //               (halfPi_RhoCosAlpha *
+  //                 (rhoCosAlpha +
+  //                   alpha * rhosinAlpha * cos(rhoCosAlpha) -
+  //                   sinRhoCosAlpha) -
+  //                 rhosinAlpha * (rhosinAlpha - alpha * sinRhoCosAlpha))) /
+  //               (halfPi_RhoCosAlpha * halfPi_RhoCosAlpha));
+  //       }
+  //       r =
+  //         rho +
+  //         (gingeryIntegrate(s_, rhoCosAlpha, x) * (pi - rho)) /
+  //           gingeryIntegrate(s_, rhoCosAlpha, pi);
+  //       theta = theta0 + alpha;
+  //       x = r * cos(theta);
+  //       y = r * sin(theta);
+  //     }
+  //     return azimuthalEquidistantRaw.invert(x, y);
+  //   };
+
+  return forward;
+}
+
+function gingeryLength(alpha, k) {
+  return function (x) {
+    var y_ = alpha * cos(x);
+    if (x < halfPi) y_ -= k;
+    return sqrt(1 + y_ * y_);
+  };
+}
+
+// Numerical integration: trapezoidal rule.
+function gingeryIntegrate(f, a, b) {
+  var n = 50,
+    h = (b - a) / n,
+    s = f(a) + f(b);
+  for (var i = 1, x = a; i < n; ++i) s += 2 * f((x += h));
+  return s * 0.5 * h;
+}
+
+export default function () {
+  var n = 6,
+    rho = 30 * radians,
+    cRho = cos(rho),
+    sRho = sin(rho),
+    m = d3.geoProjectionMutator(gingeryRaw),
+    p = m(rho, n),
+    stream_ = p.stream,
+    epsilon = 1e-2,
+    cr = -cos(epsilon * radians),
+    sr = sin(epsilon * radians);
+
+  //   p.radius = function (_) {
+  //     if (!arguments.length) return rho * degrees;
+  //     cRho = cos((rho = _ * radians));
+  //     sRho = sin(rho);
+  //     return m(rho, n);
+  //   };
+
+  //   p.lobes = function (_) {
+  //     if (!arguments.length) return n;
+  //     return m(rho, (n = +_));
+  //   };
+
+  p.stream = function (stream) {
+    var rotate = p.rotate(),
+      rotateStream = stream_(stream),
+      sphereStream = (p.rotate([0, 0]), stream_(stream));
+    p.rotate(rotate);
+    rotateStream.sphere = function () {
+      sphereStream.polygonStart(), sphereStream.lineStart();
+      for (var i = 0, delta = (2 * pi) / n, phi = 0; i < n; ++i, phi -= delta) {
+        sphereStream.point(
+          atan2(sr * cos(phi), cr) * degrees,
+          asin(sr * sin(phi)) * degrees
+        );
+        sphereStream.point(
+          atan2(sRho * cos(phi - delta / 2), cRho) * degrees,
+          asin(sRho * sin(phi - delta / 2)) * degrees
+        );
+      }
+      sphereStream.lineEnd(), sphereStream.polygonEnd();
+    };
+    return rotateStream;
+  };
+
+  return p
+    .rotate([0, 90])
+    .scale(91.7095)
+    .clipAngle(180 - 1e-3);
+}

challenges/26-projections/gingeryspiral.html

@@ -0,0 +1,92 @@
+<!DOCTYPE html>
+<html lang="en">
+  <head>
+    <meta charset="UTF-8" />
+    <meta name="viewport" content="width=device-width, initial-scale=1.0" />
+    <script src="//d3js.org/topojson.v1.min.js"></script>
+    <title>Custom Gingery-Spiral Projection Map</title>
+    <style>
+      .stroke {
+        fill: none;
+        stroke: #000;
+        stroke-width: 3px;
+      }
+
+      .fill {
+        fill: #fff;
+      }
+
+      .graticule {
+        fill: none;
+        stroke: #777;
+        stroke-width: 0.5px;
+        stroke-opacity: 0.5;
+      }
+
+      .land {
+        fill: #222;
+      }
+
+      .boundary {
+        fill: none;
+        stroke: #fff;
+        stroke-width: 0.5px;
+      }
+    </style>
+  </head>
+
+  <body>
+    <svg id="map" width="960" height="600"></svg>
+    <script type="module">
+      import * as d3 from "https://cdn.jsdelivr.net/npm/d3@7/+esm";
+      //   import { geoGingery } from "https://cdn.skypack.dev/d3-geo-projection@4";
+      import gingery from "./gingery.js";
+      //   const projection = geoGingery().rotate([0, -90]);
+
+      // Set up the SVG container
+      var width = 960;
+      var height = 600;
+
+      var svg = d3.select("svg");
+      var projection = gingery();
+
+      // Create a path generator with the custom projection
+      var graticule = d3.geoGraticule();
+      var path = d3.geoPath().projection(projection);
+
+      // Add sphere to clip the map
+      var defs = svg.append("defs");
+      defs
+        .append("path")
+        .datum({ type: "Sphere" })
+        .attr("id", "sphere")
+        .attr("d", path);
+      defs
+        .append("clipPath")
+        .attr("id", "clip")
+        .append("use")
+        .attr("xlink:href", "#sphere");
+
+      svg.append("use").attr("class", "stroke").attr("xlink:href", "#sphere");
+      svg.append("use").attr("class", "fill").attr("xlink:href", "#sphere");
+
+      svg
+        .append("path")
+        .datum(graticule)
+        .attr("class", "graticule")
+        .attr("clip-path", "url(#clip)")
+        .attr("d", path);
+
+      // Load and display the world map
+      d3.json("./world-110m.json").then(function (world) {
+        svg
+          .append("path")
+          .datum(topojson.feature(world, world.objects.countries))
+          .attr("d", path)
+          .attr("fill", "#cccccc")
+          .attr("stroke", "#000000")
+          .attr("clip-path", "url(#clip)");
+      });
+    </script>
+  </body>
+</html>

challenges/26-projections/math.js

@@ -0,0 +1,68 @@
+export var abs = Math.abs;
+export var atan = Math.atan;
+export var atan2 = Math.atan2;
+export var ceil = Math.ceil;
+export var cos = Math.cos;
+export var exp = Math.exp;
+export var floor = Math.floor;
+export var log = Math.log;
+export var max = Math.max;
+export var min = Math.min;
+export var pow = Math.pow;
+export var round = Math.round;
+export var sign =
+  Math.sign ||
+  function (x) {
+    return x > 0 ? 1 : x < 0 ? -1 : 0;
+  };
+export var sin = Math.sin;
+export var tan = Math.tan;
+
+export var epsilon = 1e-6;
+export var epsilon2 = 1e-12;
+export var pi = Math.PI;
+export var halfPi = pi / 2;
+export var quarterPi = pi / 4;
+export var sqrt1_2 = Math.SQRT1_2;
+export var sqrt2 = sqrt(2);
+export var sqrtPi = sqrt(pi);
+export var tau = pi * 2;
+export var degrees = 180 / pi;
+export var radians = pi / 180;
+
+export function sinci(x) {
+  return x ? x / Math.sin(x) : 1;
+}
+
+export function asin(x) {
+  return x > 1 ? halfPi : x < -1 ? -halfPi : Math.asin(x);
+}
+
+export function acos(x) {
+  return x > 1 ? 0 : x < -1 ? pi : Math.acos(x);
+}
+
+export function sqrt(x) {
+  return x > 0 ? Math.sqrt(x) : 0;
+}
+
+export function tanh(x) {
+  x = exp(2 * x);
+  return (x - 1) / (x + 1);
+}
+
+export function sinh(x) {
+  return (exp(x) - exp(-x)) / 2;
+}
+
+export function cosh(x) {
+  return (exp(x) + exp(-x)) / 2;
+}
+
+export function arsinh(x) {
+  return log(x + sqrt(x * x + 1));
+}
+
+export function arcosh(x) {
+  return log(x + sqrt(x * x - 1));
+}