feat(s2): New module @math.gl/s2

docs/whats-new.md

@@ -31,6 +31,10 @@
 
 - Add `isTypedArray()` and `isNumericArray()` utilities that both check values and return properly restricted types to help write strictly typed code (avoids the `DataView` issue with `ArrayBuffer.isView()`).
 
+**`@math.gl/s2`** (NEW MODULE)
+
+- New module that contains a lightweight implementation of the S2 DGGS (Discrete Global Grid System).
+
 ## v3.6
 
 Target Release Date: Q4, 2021

modules/s2/README.md

@@ -0,0 +1,7 @@
+# @math.gl/s2
+
+[math.gl](https://math.gl/docs) is a suite of math modules for 3D and geospatial applications.
+
+This module contains math to work with the S2 DGGS (Discrete Global Grid System).
+
+For documentation please visit the [website](https://math.gl).

modules/s2/docs/README.md

@@ -0,0 +1,21 @@
+# Overview
+
+`@math.gl/s2` is a small JavaScript library for working with the S2 DGGS (Discrete Global Grid System).
+
+## Installation
+
+```bash
+npm install @math.gl/s2
+```
+
+## Usage
+
+```js
+import {} from '@math.gl/s2';
+```
+
+## Attribution
+
+This module is a fork of a subset of the s2-geometry module under ISC License (ISC)
+Copyright (c) 2012-2016, Jon Atkins <[email protected]>
+Copyright (c) 2016, AJ ONeal <[email protected]>

modules/s2/docs/api-reference/s2.md

@@ -0,0 +1,5 @@
+# S2
+
+See [s2geometry.io](https://s2geometry.io/) for more information.
+
+> TBA

modules/s2/package.json

@@ -0,0 +1,35 @@
+{
+  "name": "@math.gl/s2",
+  "description": "A minimal S2 DGGS (Discrete Global Grid System) implementation",
+  "license": "MIT",
+  "publishConfig": {
+    "access": "public"
+  },
+  "version": "3.6.3",
+  "keywords": [
+    "javascript",
+    "math",
+    "geospatial",
+    "sun",
+    "suncalc",
+    "solar"
+  ],
+  "repository": {
+    "type": "git",
+    "url": "https://github.com/uber-web/math.gl.git"
+  },
+  "types": "dist/index.d.ts",
+  "main": "dist/es5/index.js",
+  "module": "dist/esm/index.js",
+  "files": [
+    "dist",
+    "src"
+  ],
+  "dependencies": {
+    "@babel/runtime": "^7.12.0",
+    "long": "^5.2.1"
+  },
+  "devDependencies": {
+    "s2-geometry": "^1.2.10"
+  }
+}

modules/s2/src/index.ts

@@ -0,0 +1,13 @@
+// math.gl MIT license
+
+export {getIdFromToken, getGeoBounds, getS2QuadKey, getS2Polygon} from './s2-utils';
+
+export {
+  getTokenFromId,
+  getS2ChildCellId,
+  getS2CellFromToken,
+  get2dRegionFromS2Cell
+} from './s2-utils-ext';
+
+export type {S2HeightInfo} from './s2-to-obb-points';
+export {getOrientedBoundingBoxCornerPoints} from './s2-to-obb-points';

modules/s2/src/s2-geometry.ts

@@ -0,0 +1,204 @@
+// math.gl, MIT license
+/*
+Adapted from s2-geometry under ISC License (ISC)
+Copyright (c) 2012-2016, Jon Atkins <[email protected]>
+Copyright (c) 2016, AJ ONeal <[email protected]>
+Permission to use, copy, modify, and/or distribute this software for any purpose with or without fee is hereby granted, provided that the above copyright notice and this permission notice appear in all copies.
+THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.
+*/
+
+import Long from 'long';
+
+//
+// Functional Style
+//
+const FACE_BITS = 3;
+const MAX_LEVEL = 30;
+const POS_BITS = 2 * MAX_LEVEL + 1; // 61 (60 bits of data, 1 bit lsb marker)
+const RADIAN_TO_DEGREE = 180 / Math.PI;
+
+export function IJToST(
+  ij: [number, number],
+  order: number,
+  offsets: [number, number]
+): [number, number] {
+  const maxSize = 1 << order;
+
+  return [(ij[0] + offsets[0]) / maxSize, (ij[1] + offsets[1]) / maxSize];
+}
+
+function singleSTtoUV(st: number): number {
+  if (st >= 0.5) {
+    return (1 / 3.0) * (4 * st * st - 1);
+  }
+  return (1 / 3.0) * (1 - 4 * (1 - st) * (1 - st));
+}
+
+export function STToUV(st: [number, number]): [number, number] {
+  return [singleSTtoUV(st[0]), singleSTtoUV(st[1])];
+}
+
+export function FaceUVToXYZ(face: number, [u, v]: [number, number]): [number, number, number] {
+  switch (face) {
+    case 0:
+      return [1, u, v];
+    case 1:
+      return [-u, 1, v];
+    case 2:
+      return [-u, -v, 1];
+    case 3:
+      return [-1, -v, -u];
+    case 4:
+      return [v, -1, -u];
+    case 5:
+      return [v, u, -1];
+    default:
+      throw new Error('Invalid face');
+  }
+}
+
+export function XYZToLngLat([x, y, z]: [number, number, number]): [number, number] {
+  const lat = Math.atan2(z, Math.sqrt(x * x + y * y));
+  const lng = Math.atan2(y, x);
+
+  return [lng * RADIAN_TO_DEGREE, lat * RADIAN_TO_DEGREE];
+}
+
+/*
+  Here is the version of toHilbertQuadkey taken from deck.gl
+  We replace it with the function that takes Long instead of string.
+  The original function does not support the case of level==0
+
+export function toHilbertQuadkey(idS: string): string {
+  let bin = Long.fromString(idS, true, 10).toString(2);
+
+  while (bin.length < FACE_BITS + POS_BITS) {
+    // eslint-disable-next-line prefer-template
+    bin = '0' + bin;
+  }
+
+  // MUST come AFTER binstr has been left-padded with '0's
+  const lsbIndex = bin.lastIndexOf('1');
+  // substr(start, len)
+  // substring(start, end) // includes start, does not include end
+  const faceB = bin.substring(0, 3);
+  // posB will always be a multiple of 2 (or it's invalid)
+  const posB = bin.substring(3, lsbIndex);
+  const levelN = posB.length / 2;
+
+  const faceS = Long.fromString(faceB, true, 2).toString(10);
+  let posS = Long.fromString(posB, true, 2).toString(4);
+
+  while (posS.length < levelN) {
+    // eslint-disable-next-line prefer-template
+    posS = '0' + posS;
+  }
+
+  return `${faceS}/${posS}`;
+}
+*/
+
+export function toHilbertQuadkey(id: Long): string {
+  let bin = id.toString(2);
+
+  while (bin.length < FACE_BITS + POS_BITS) {
+    // eslint-disable-next-line prefer-template
+    bin = '0' + bin;
+  }
+
+  // MUST come AFTER binstr has been left-padded with '0's
+  const lsbIndex = bin.lastIndexOf('1');
+  // substr(start, len)
+  // substring(start, end) // includes start, does not include end
+  const faceB = bin.substring(0, 3);
+  // posB will always be a multiple of 2 (or it's invalid)
+  const posB = bin.substring(3, lsbIndex);
+  const levelN = posB.length / 2;
+
+  const faceS = Long.fromString(faceB, true, 2).toString(10);
+
+  /*
+    Here is a fix for the case when posB is an empty string that causes an exception in Long.fromString
+
+    let posS = Long.fromString(posB, true, 2).toString(4);
+  */
+  let posS = '0';
+  if (levelN !== 0) {
+    // posB is not an empty string< because levelN!==0
+    posS = Long.fromString(posB, true, 2).toString(4);
+
+    while (posS.length < levelN) {
+      // eslint-disable-next-line prefer-template
+      posS = '0' + posS;
+    }
+  }
+  // Note, posS will be "0" for the level==0
+  // TODO: Is it ok?
+
+  return `${faceS}/${posS}`;
+}
+
+function rotateAndFlipQuadrant(n: number, point: [number, number], rx: number, ry: number): void {
+  if (ry === 0) {
+    if (rx === 1) {
+      point[0] = n - 1 - point[0];
+      point[1] = n - 1 - point[1];
+    }
+
+    const x = point[0];
+    point[0] = point[1];
+    point[1] = x;
+  }
+}
+
+/*
+  Original function taken from deck.gl doesn't support the case of (face <= 5)
+  It's fixed here.
+*/
+export function FromHilbertQuadKey(hilbertQuadkey: string): {
+  face: number;
+  ij: [number, number];
+  level: number;
+} {
+  const parts = hilbertQuadkey.split('/');
+  const face = parseInt(parts[0], 10);
+  const position = parts[1];
+  /*
+    Fix for the case of level==0 that corresponds to (face <= 5)
+    const maxLevel = position.length;
+    let level;
+  */
+  const maxLevel = face > 5 ? position.length : 0;
+  let level = 0;
+
+  const point = [0, 0] as [number, number];
+
+  for (let i = maxLevel - 1; i >= 0; i--) {
+    level = maxLevel - i;
+    const bit = position[i];
+    let rx = 0;
+    let ry = 0;
+    if (bit === '1') {
+      ry = 1;
+    } else if (bit === '2') {
+      rx = 1;
+      ry = 1;
+    } else if (bit === '3') {
+      rx = 1;
+    }
+
+    const val = Math.pow(2, level - 1);
+    rotateAndFlipQuadrant(val, point, rx, ry);
+
+    point[0] += val * rx;
+    point[1] += val * ry;
+  }
+
+  if (face % 2 === 1) {
+    const t = point[0];
+    point[0] = point[1];
+    point[1] = t;
+  }
+
+  return {face, ij: point, level};
+}

modules/s2/src/s2-to-obb-points.ts

@@ -0,0 +1,72 @@
+import {getS2CellFromToken, get2dRegionFromS2Cell} from './s2-utils-ext';
+import {Vector3} from '@math.gl/core';
+
+// import {OrientedBoundingBox, makeOrientedBoundingBoxFromPoints} from '@math.gl/culling';
+// import {Ellipsoid} from '@math.gl/geospatial';
+
+export type S2HeightInfo = {
+  minimumHeight: number;
+  maximumHeight: number;
+};
+
+/**
+ * Converts S2HeightInfo to corner points of an oriented bounding box
+ * Can be used to constuct an OrientedBoundingBox instance
+ *
+ * @param s2bv - s2 bounding volume to convert
+ * @param [result] Optional object onto which to store the result.
+ * @returns The modified result parameter or a new `OrientedBoundingBox` instance if not provided.
+ */
+export function getOrientedBoundingBoxCornerPoints(
+  token: string, // This can be an S2 key or token
+  heighInfo?: S2HeightInfo
+): Vector3[] {
+  const min: number = heighInfo?.minimumHeight || 0;
+  const max: number = heighInfo?.maximumHeight || 0;
+
+  const s2cell = getS2CellFromToken(token);
+  const region = get2dRegionFromS2Cell(s2cell);
+
+  // region is {lngWest, latSouth, lngEast, latNorth} in degrees
+  const W = region[0];
+  const S = region[1];
+  const E = region[2];
+  const N = region[3];
+
+  const points: Vector3[] = [];
+
+  points.push(new Vector3(W, N, min));
+  points.push(new Vector3(E, N, min));
+  points.push(new Vector3(E, S, min));
+  points.push(new Vector3(W, S, min));
+
+  points.push(new Vector3(W, N, max));
+  points.push(new Vector3(E, N, max));
+  points.push(new Vector3(E, S, max));
+  points.push(new Vector3(W, S, max));
+
+  return points;
+}
+
+/*
+function convertCartToXYZ(longitude, latitude, height, result?) {
+  const point = Ellipsoid.WGS84.cartographicToCartesian([longitude, latitude, height]);
+  return new Vector3(point[0], point[1], point[2]);
+}
+
+// Add a point that doesn't allow the box dive under the Earth
+  // This point is actually a center of a face that could be a tangent to the Earth surface if max==0
+  points.push(convertCartToXYZ((W + E) / 2.0, (S + N) / 2.0, max));
+
+  // points should be an array of Vector3 (XYZ)
+  // Passing result===null throws an exception from makeOrientedBoundingBoxFromPoints
+  const obb: OrientedBoundingBox = makeOrientedBoundingBoxFromPoints(
+    points,
+    result !== null ? result : undefined
+  );
+
+  const box: number[] = [...obb.center, ...obb.halfAxes];
+
+  return box;
+}
+*/