Add Equal Earth, Natural Earth and Lambert Conformal Conic projections (

#11091)

debug/projections.html

@@ -43,12 +43,14 @@
         <fieldset>
             <label>Projection:</label>
             <select id="projName">
-                <option value="mercator">Mercator</option>
-                <option value="winkel">Winkel Tripel</option>
-                <option value="sinusoidal">Sinusoidal</option>
-                <option value="wgs84">WGS84</option>
-                <option value="albers" selected>Albers USA</option>
                 <option value="alaska">Albers Alaska</option>
+                <option value="albers" selected>Albers USA</option>
+                <option value="equalEarth">Equal Earth</option>
+                <option value="equirectangular">Equirectangular</option>
+                <option value="lambertConformalConic">Lambert Conformal Conic</option>
+                <option value="mercator">Mercator</option>
+                <option value="naturalEarth">Natural Earth</option>
+                <option value="winkelTripel">Winkel Tripel</option>
             </select>
         </fieldset>
         <fieldset>
@@ -84,18 +86,24 @@
     if (map) map.remove();
     const el = document.getElementById('projName');
     const zooms = {
-        albers: 3,
         alaska: 4,
-        winkel: 1.2,
+        albers: 3,
+        equalEarth: 1.0,
+        equirectangular: 1.0,
+        lambertConformalConic: 1.0,
         mercator: 1.0,
-        wgs84: 1.0
+        naturalEarth: 1.0,
+        winkelTripel: 1.2
     };
     const centers = {
-        albers: [-122.414, 37.776],
         alaska: [-154, 63],
-        winkel: [0, 0],
+        albers: [-122.414, 37.776],
+        equalEarth: [0, 0],
+        equirectangular: [0, 0],
+        lambertConformalConic: [0, 0],
         mercator: [0, 0],
-        wgs84: [0, 0]
+        naturalEarth: [0, 0],
+        winkelTripel: [0, 0]
     };
     const projection = el.options[el.selectedIndex].value;
     const zoom = zooms[projection] || 2;

src/geo/projection/equalEarth.js

@@ -0,0 +1,56 @@
+// @flow
+import LngLat from '../lng_lat.js';
+import {clamp} from '../../util/util.js';
+
+const a1 = 1.340264;
+const a2 = -0.081106;
+const a3 = 0.000893;
+const a4 = 0.003796;
+const M = Math.sqrt(3) / 2;
+
+export default {
+    name: 'equalEarth',
+    center: [0, 0],
+    range: [3.5, 7],
+
+    project(lng: number, lat: number) {
+        // based on https://github.com/d3/d3-geo, MIT-licensed
+        lat = lat / 180 * Math.PI;
+        lng = lng / 180 * Math.PI;
+        const theta = Math.asin(M * Math.sin(lat));
+        const theta2 = theta * theta;
+        const theta6 = theta2 * theta2 * theta2;
+        const x = lng * Math.cos(theta) / (M * (a1 + 3 * a2 * theta2 + theta6 * (7 * a3 + 9 * a4 * theta2)));
+        const y = theta * (a1 + a2 * theta2 + theta6 * (a3 + a4 * theta2));
+
+        return {
+            x: (x / Math.PI + 0.5) * 0.5,
+            y: 1 - (y / Math.PI + 0.5) * 0.5
+        };
+    },
+
+    unproject(x: number, y: number) {
+        // based on https://github.com/d3/d3-geo, MIT-licensed
+        x = (2 * x - 0.5) * Math.PI;
+        y = (2 * (1 - y) - 0.5) * Math.PI;
+        let theta = y;
+        let theta2 = theta * theta;
+        let theta6 = theta2 * theta2 * theta2;
+
+        for (let i = 0, delta, fy, fpy; i < 12; ++i) {
+            fy = theta * (a1 + a2 * theta2 + theta6 * (a3 + a4 * theta2)) - y;
+            fpy = a1 + 3 * a2 * theta2 + theta6 * (7 * a3 + 9 * a4 * theta2);
+            theta -= delta = fy / fpy;
+            theta2 = theta * theta;
+            theta6 = theta2 * theta2 * theta2;
+            if (Math.abs(delta) < 1e-12) break;
+        }
+
+        const lambda = M * x * (a1 + 3 * a2 * theta2 + theta6 * (7 * a3 + 9 * a4 * theta2)) / Math.cos(theta);
+        const phi = Math.asin(Math.sin(theta) / M);
+        const lng = clamp(lambda * 180 / Math.PI, -180, 180);
+        const lat = clamp(phi * 180 / Math.PI, -90, 90);
+
+        return new LngLat(lng, lat);
+    }
+};

src/geo/projection/wgs84.js → src/geo/projection/equirectangular.js

@@ -3,7 +3,7 @@ import LngLat from '../lng_lat.js';
 import {clamp} from '../../util/util.js';
 
 export default {
-    name: 'wgs84',
+    name: 'equirectangular',
     center: [0, 0],
     project(lng: number, lat: number) {
         const x = 0.5 + lng / 360;

src/geo/projection/index.js

@@ -1,8 +1,11 @@
 // @flow
-import {albers, alaska} from './albers.js';
+import {alaska, albers} from './albers.js';
+import equalEarth from './equalEarth.js';
+import equirectangular from './equirectangular.js';
+import lambertConformalConic from './lambert.js';
 import mercator from './mercator.js';
-import wgs84 from './wgs84.js';
-import winkel from './winkelTripel.js';
+import naturalEarth from './naturalEarth.js';
+import winkelTripel from './winkelTripel.js';
 import LngLat from '../lng_lat.js';
 
 export type Projection = {
@@ -14,11 +17,14 @@ export type Projection = {
 };
 
 const projections = {
-    albers,
     alaska,
+    albers,
+    equalEarth,
+    equirectangular,
+    lambertConformalConic,
     mercator,
-    wgs84,
-    winkel
+    naturalEarth,
+    winkelTripel
 };
 
 export default function getProjection(config: {name: string} | string) {

src/geo/projection/lambert.js

@@ -0,0 +1,68 @@
+// @flow
+import LngLat from '../lng_lat.js';
+import {clamp} from '../../util/util.js';
+
+const halfPi = Math.PI / 2;
+
+function tany(y) {
+    return Math.tan((halfPi + y) / 2);
+}
+
+function getParams([lat0, lat1]) {
+    const y0 = lat0 * Math.PI / 180;
+    const y1 = lat1 * Math.PI / 180;
+    const cy0 = Math.cos(y0);
+    const n = y0 === y1 ? Math.sin(y0) : Math.log(cy0 / Math.cos(y1)) / Math.log(tany(y1) / tany(y0));
+    const f = cy0 * Math.pow(tany(y0), n) / n;
+
+    return {n, f};
+}
+
+export default {
+    name: 'lambertConformalConic',
+    range: [3.5, 7],
+
+    center: [0, 30],
+    parallels: [30, 30],
+
+    project(lng: number, lat: number) {
+        // based on https://github.com/d3/d3-geo, MIT-licensed
+        lat = lat / 180 * Math.PI;
+        lng = lng / 180 * Math.PI;
+
+        const epsilon = 1e-6;
+        const {n, f} = getParams(this.parallels);
+
+        if (f > 0) {
+            if (lat < -halfPi + epsilon) lat = -halfPi + epsilon;
+        } else {
+            if (lat > halfPi - epsilon) lat = halfPi - epsilon;
+        }
+
+        const r = f / Math.pow(tany(lat), n);
+        const x = r * Math.sin(n * lng);
+        const y = f - r * Math.cos(n * lng);
+
+        return {
+            x: (x / Math.PI + 0.5) * 0.5,
+            y: 1 - (y / Math.PI + 0.5) * 0.5
+        };
+    },
+
+    unproject(x: number, y: number) {
+        // based on https://github.com/d3/d3-geo, MIT-licensed
+        x = (2 * x - 0.5) * Math.PI;
+        y = (2 * (1 - y) - 0.5) * Math.PI;
+        const {n, f} = getParams(this.parallels);
+        const fy = f - y;
+        const r = Math.sign(n) * Math.sqrt(x * x + fy * fy);
+        let l = Math.atan2(x, Math.abs(fy)) * Math.sign(fy);
+
+        if (fy * n < 0) l -= Math.PI * Math.sign(x) * Math.sign(fy);
+
+        const lng = clamp((l / n)  * 180 / Math.PI, -180, 180);
+        const lat = clamp((2 * Math.atan(Math.pow(f / r, 1 / n)) - halfPi)  * 180 / Math.PI, -90, 90);
+
+        return new LngLat(lng, lat);
+    }
+};

src/geo/projection/naturalEarth.js

@@ -0,0 +1,51 @@
+// @flow
+import LngLat from '../lng_lat.js';
+import {clamp} from '../../util/util.js';
+
+export default {
+    name: 'naturalEarth',
+    center: [0, 0],
+    range: [3.5, 7],
+
+    project(lng: number, lat: number) {
+        // based on https://github.com/d3/d3-geo, MIT-licensed
+        lat = lat / 180 * Math.PI;
+        lng = lng / 180 * Math.PI;
+
+        const phi2 = lat * lat;
+        const phi4 = phi2 * phi2;
+        const x = lng * (0.8707 - 0.131979 * phi2 + phi4 * (-0.013791 + phi4 * (0.003971 * phi2 - 0.001529 * phi4)));
+        const y = lat * (1.007226 + phi2 * (0.015085 + phi4 * (-0.044475 + 0.028874 * phi2 - 0.005916 * phi4)));
+
+        return {
+            x: (x / Math.PI + 0.5) * 0.5,
+            y: 1 - (y / Math.PI + 0.5) * 0.5
+        };
+    },
+
+    unproject(x: number, y: number) {
+        // based on https://github.com/d3/d3-geo, MIT-licensed
+        x = (2 * x - 0.5) * Math.PI;
+        y = (2 * (1 - y) - 0.5) * Math.PI;
+        const epsilon = 1e-6;
+        let phi = y;
+        let i = 25;
+        let delta = 0;
+        let phi2 = phi * phi;
+
+        do {
+            phi2 = phi * phi;
+            const phi4 = phi2 * phi2;
+            phi -= delta = (phi * (1.007226 + phi2 * (0.015085 + phi4 * (-0.044475 + 0.028874 * phi2 - 0.005916 * phi4))) - y) /
+                (1.007226 + phi2 * (0.015085 * 3 + phi4 * (-0.044475 * 7 + 0.028874 * 9 * phi2 - 0.005916 * 11 * phi4)));
+        } while (Math.abs(delta) > epsilon && --i > 0);
+
+        phi2 = phi * phi;
+        const lambda = x / (0.8707 + phi2 * (-0.131979 + phi2 * (-0.013791 + phi2 * phi2 * phi2 * (0.003971 - 0.001529 * phi2))));
+
+        const lng = clamp(lambda * 180 / Math.PI, -180, 180);
+        const lat = clamp(phi * 180 / Math.PI, -90, 90);
+
+        return new LngLat(lng, lat);
+    }
+};

src/geo/projection/winkelTripel.js

@@ -3,7 +3,7 @@ import LngLat from '../lng_lat.js';
 import {clamp} from '../../util/util.js';
 
 export default {
-    name: 'winkel',
+    name: 'winkelTripel',
     center: [0, 0],
     range: [3.5, 7],
 

src/ui/map.js

@@ -268,7 +268,7 @@ const defaultOptions = {
  * @param {Object} [options.locale=null] A patch to apply to the default localization table for UI strings, e.g. control tooltips. The `locale` object maps namespaced UI string IDs to translated strings in the target language;
  *  see `src/ui/default_locale.js` for an example with all supported string IDs. The object may specify all UI strings (thereby adding support for a new translation) or only a subset of strings (thereby patching the default translation table).
  * @param {boolean} [options.testMode=false] Silences errors and warnings generated due to an invalid accessToken, useful when using the library to write unit tests.
- * @param {string} [options.projection='mercator'] The map projection to use when creating a map. Defaults to the Web Mercator ('mercator') projection. Other options are Winkel Tripel ('winkel'), Sinusoidal ('sinusoidal'), Albers ('albers'), Albers Alaska ('alaska'), and WGS84 ('wgs84').
+ * @param {string} [options.projection='mercator'] The map projection to use when creating a map. Defaults to the Web Mercator ('mercator') projection. Other options are Albers Alaska ('alaska'), Albers USA ('albers'), Equal Earth ('equalEarth'), Equirectangular/Plate Carrée/WGS84 ('equirectangular'), Lambert Conic Conformal ('lambertConformalConic'), Natural Earth ('naturalEarth') and Winkel Tripel ('winkelTripel').
  * @example
  * var map = new mapboxgl.Map({
  *   container: 'map',