Merge branch 'master' into master

chrisveness · Sep 29, 2022 · be35438 · be35438
2 parents 5477111 + df92e65
commit be35438
Show file tree

Hide file tree

Showing 24 changed files with 384 additions and 200 deletions.
diff --git a/.travis.yml b/.travis.yml
@@ -1,13 +1,16 @@
 language: node_js
 
 node_js:
-- "node"
-- "14"
+- node
+- lts/*
+- 14
 
 os:
 - linux
 - osx
-- windows
+# windows - removed 'cos travis have made non-backward-compatible move to nvs on windows (travis-ci.community/t/12393) - sigh!
+dist: jammy          # as of Sep 2022 travis defaults Linux to xenial, which is incompatible with node.js 18.x.x
+osx_image: xcode13.2 # as of Sep 2022 travis defaults macOS to 10.13, which is incompatible with node.js 18.x.x
 
 after_success:
 - c8 -r text-lcov npm test | coveralls
diff --git a/CHANGELOG.md b/CHANGELOG.md
@@ -4,12 +4,35 @@
 
 ### Fixed
 
-- Fix parsing of 'H' 500km squares (Scottish islands)
+- Truncate MGRS easting / northing values to max 1 metre resolution
+- Fix UTM constructor northing range check
+- Fix Mgrs.toUtm() edge case at zone boundaries (e.g. @ 64°S,0°E)
+- Fix rounding error in Utm.toMgrs() which caused UTM for 80°S,0°E to fail
+- Allow single-digit zone when parsing MGRS grid reference [#104]
+
+## [2.4.0] - 2022-03-16
+
+### Fixed
+
+- Fix check for coincident points (previously < ≈95mm got treated as coincident)
+- Add check for null arguments to LatLonEllipsoidal constructor
+
+### Added
+
+- LatLonNvectorSpherical.centreOf()
+
+## [2.3.0] - 2021-11-16
+
+### Fixed
+
+- Fix parsing of 'H' 500km squares (Scottish islands) [#96]
 - Fix Dms.wrap90(), Dms.wrap180() to work for all -ve degrees
+- LatLon_OsGridRef: Override super.convertDatum()
 
 ### Added
 
 - LatLonEllipsoidal_Vincenty.intermediatePointTo()
+- Extra type-checking (LatLonEllipsoidal_Vincenty.direct, LatLonNvectorSpherical.isEnclosedBy)
 
 ## [2.2.1] - 2020-04-22
 

diff --git a/README.md b/README.md
@@ -1,7 +1,7 @@
 Geodesy functions
 =================
 
-[![Build Status](https://travis-ci.org/chrisveness/geodesy.svg?branch=master)](https://travis-ci.org/chrisveness/geodesy)
+[![Build Status](https://travis-ci.com/chrisveness/geodesy.svg?branch=master)](https://app.travis-ci.com/github/chrisveness/geodesy)
 [![Coverage Status](https://coveralls.io/repos/github/chrisveness/geodesy/badge.svg)](https://coveralls.io/github/chrisveness/geodesy)
 [![Documentation](https://img.shields.io/badge/docs-www.movable--type.co.uk%2Fscripts%2Fgeodesy--library.html-lightgrey.svg)](https://www.movable-type.co.uk/scripts/geodesy-library.html)
 
@@ -96,7 +96,7 @@ The library can be used in the browser by taking a local copy, or loading it fro
 ```html
 <!doctype html><title>geodesy example</title><meta charset="utf-8">
 <script type="module">
-    import LatLon from 'https://cdn.jsdelivr.net/npm/geodesy@2.2.1/latlon-spherical.min.js';
+    import LatLon from 'https://cdn.jsdelivr.net/npm/geodesy@2.4.0/latlon-spherical.min.js';
 
     const p1 = new LatLon(50.06632, -5.71475);
     const p2 = new LatLon(58.64402, -3.07009);
@@ -112,13 +112,12 @@ The library can be used in the browser by taking a local copy, or loading it fro
 ### Usage in Node.js
 
 The library can be loaded from [npm](https://www.npmjs.com/package/geodesy) to be used in a Node.js app 
-(in Node.js v13.2.0+, or using the [esm](https://www.npmjs.com/package/esm) package in 
-v8.0.0–v12.15.0<sup title="v12.16.0+ is not compatible with esm@3.2.25">*</sup>):
+(in Node.js v13.2.0+, or Node.js v12.0.0+ using --experimental-modules, or v8.0.0–v12.15.0<sup title="v12.16.0+ is not compatible with esm@3.2.25">*</sup>) using the [esm](https://www.npmjs.com/package/esm) package:
 
 ```shell
-$ npm install geodesy esm
-$ node -r esm
-> import LatLon from 'geodesy/latlon-spherical.js';
+$ npm install geodesy
+$ node
+> const { default: LatLon } = await import('geodesy/latlon-spherical.js');
 > const p1 = new LatLon(50.06632, -5.71475);
 > const p2 = new LatLon(58.64402, -3.07009);
 > const d = p1.distanceTo(p2);
@@ -127,6 +126,22 @@ $ node -r esm
 > console.assert(mid.toString('dms') == '54° 21′ 44″ N, 004° 31′ 51″ W');
 ```
 
+To some extent, mixins can be used to add methods of a class to a different class, e.g.:
+
+```javascript
+import LatLon  from 'geodesy/latlon-nvector-ellipsoidal.js';
+import LatLonV from 'geodesy/latlon-ellipsoidal-vincenty.js';
+
+for (const method of Object.getOwnPropertyNames(LatLonV.prototype)) {
+    LatLon.prototype[method] = LatLonV.prototype[method];
+}
+
+const d = new LatLon(51, 0).distanceTo(new LatLon(52, 1)); // vincenty
+const δ = new LatLon(51, 0).deltaTo(new LatLon(52, 1));    // n-vector
+```
+
+More care is of course required if there are conflicting constructors or method names.
+
 For TypeScript users, type definitions are available from DefinitelyTyped: [www.npmjs.com/package/@types/geodesy](https://www.npmjs.com/package/@types/geodesy).
 
 ### Other examples

diff --git a/latlon-ellipsoidal-datum.js b/latlon-ellipsoidal-datum.js
@@ -1,5 +1,5 @@
 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  */
-/* Geodesy tools for conversions between (historical) datums          (c) Chris Veness 2005-2019  */
+/* Geodesy tools for conversions between (historical) datums          (c) Chris Veness 2005-2022  */
 /*                                                                                   MIT Licence  */
 /* www.movable-type.co.uk/scripts/latlong-convert-coords.html                                     */
 /* www.movable-type.co.uk/scripts/geodesy-library.html#latlon-ellipsoidal-datum                  */
@@ -39,12 +39,12 @@ const ellipsoids = {
     WGS84:         { a: 6378137,     b: 6356752.314245, f: 1/298.257223563 },
     Airy1830:      { a: 6377563.396, b: 6356256.909,    f: 1/299.3249646   },
     AiryModified:  { a: 6377340.189, b: 6356034.448,    f: 1/299.3249646   },
-    Bessel1841:    { a: 6377397.155, b: 6356078.962818, f: 1/299.1528128   },
+    Bessel1841:    { a: 6377397.155, b: 6356078.962822, f: 1/299.15281285  },
     Clarke1866:    { a: 6378206.4,   b: 6356583.8,      f: 1/294.978698214 },
     Clarke1880IGN: { a: 6378249.2,   b: 6356515.0,      f: 1/293.466021294 },
     GRS80:         { a: 6378137,     b: 6356752.314140, f: 1/298.257222101 },
-    Intl1924:      { a: 6378388,     b: 6356911.946,    f: 1/297           }, // aka Hayford
-    WGS72:         { a: 6378135,     b: 6356750.5,      f: 1/298.26        },
+    Intl1924:      { a: 6378388,     b: 6356911.946128, f: 1/297           }, // aka Hayford
+    WGS72:         { a: 6378135,     b: 6356750.52,     f: 1/298.26        },
 };
 
 

diff --git a/latlon-ellipsoidal-vincenty.js b/latlon-ellipsoidal-vincenty.js
@@ -1,8 +1,8 @@
 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  */
-/* Vincenty Direct and Inverse Solution of Geodesics on the Ellipsoid (c) Chris Veness 2002-2020  */
+/* Vincenty Direct and Inverse Solution of Geodesics on the Ellipsoid (c) Chris Veness 2002-2022  */
 /*                                                                                   MIT Licence  */
 /* www.ngs.noaa.gov/PUBS_LIB/inverse.pdf                                                          */
-/* www.movable-type.co.uk/scripts/latlong-ellipsoidal-vincenty.html                               */
+/* www.movable-type.co.uk/scripts/latlong-vincenty.html                                           */
 /* www.movable-type.co.uk/scripts/geodesy-library.html#latlon-ellipsoidal-vincenty                */
 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  */
 
@@ -69,8 +69,8 @@ class LatLonEllipsoidal_Vincenty extends LatLonEllipsoidal {
 
 
     /**
-     * Returns the initial bearing (forward azimuth) to travel along a geodesic from ‘this’ point to
-     * the given point, using Vincenty inverse solution.
+     * Returns the initial bearing to travel along a geodesic from ‘this’ point to the given point,
+     * using Vincenty inverse solution.
      *
      * @param   {LatLon} point - Latitude/longitude of destination point.
      * @returns {number} Initial bearing in degrees from north (0°..360°) or NaN if failed to converge.
@@ -92,8 +92,8 @@ class LatLonEllipsoidal_Vincenty extends LatLonEllipsoidal {
 
 
     /**
-     * Returns the final bearing (reverse azimuth) having travelled along a geodesic from ‘this’
-     * point to the given point, using Vincenty inverse solution.
+     * Returns the final bearing having travelled along a geodesic from ‘this’ point to the given
+     * point, using Vincenty inverse solution.
      *
      * @param   {LatLon} point - Latitude/longitude of destination point.
      * @returns {number} Final bearing in degrees from north (0°..360°) or NaN if failed to converge.
@@ -132,8 +132,8 @@ class LatLonEllipsoidal_Vincenty extends LatLonEllipsoidal {
 
 
     /**
-     * Returns the final bearing (reverse azimuth) having travelled along a geodesic given by initial
-     * bearing for a given distance from ‘this’ point, using Vincenty direct solution.
+     * Returns the final bearing having travelled along a geodesic given by initial bearing for a
+     * given distance from ‘this’ point, using Vincenty direct solution.
      * TODO: arg order? (this is consistent with destinationPoint, but perhaps less intuitive)
      *
      * @param   {number} distance - Distance travelled along the geodesic in metres.
@@ -164,6 +164,8 @@ class LatLonEllipsoidal_Vincenty extends LatLonEllipsoidal {
      */
     intermediatePointTo(point, fraction) {
         if (fraction == 0) return this;
+        if (fraction == 1) return point;
+
         const inverse = this.inverse(point);
         const dist = inverse.distance;
         const brng = inverse.initialBearing;
@@ -211,19 +213,18 @@ class LatLonEllipsoidal_Vincenty extends LatLonEllipsoidal {
         const A = 1 + uSq/16384*(4096+uSq*(-768+uSq*(320-175*uSq)));
         const B = uSq/1024 * (256+uSq*(-128+uSq*(74-47*uSq)));
 
-        let σ = s / (b*A), sinσ = null, cosσ = null, Δσ = null; // σ = angular distance P₁ P₂ on the sphere
+        let σ = s / (b*A), sinσ = null, cosσ = null; // σ = angular distance P₁ P₂ on the sphere
         let cos2σₘ = null; // σₘ = angular distance on the sphere from the equator to the midpoint of the line
 
         let σʹ = null, iterations = 0;
         do {
             cos2σₘ = Math.cos(2*σ1 + σ);
             sinσ = Math.sin(σ);
             cosσ = Math.cos(σ);
-            Δσ = B*sinσ*(cos2σₘ+B/4*(cosσ*(-1+2*cos2σₘ*cos2σₘ)-
-                B/6*cos2σₘ*(-3+4*sinσ*sinσ)*(-3+4*cos2σₘ*cos2σₘ)));
+            const Δσ = B*sinσ*(cos2σₘ+B/4*(cosσ*(-1+2*cos2σₘ*cos2σₘ)-B/6*cos2σₘ*(-3+4*sinσ*sinσ)*(-3+4*cos2σₘ*cos2σₘ)));
             σʹ = σ;
             σ = s / (b*A) + Δσ;
-        } while (Math.abs(σ-σʹ) > 1e-12 && ++iterations<100);
+        } while (Math.abs(σ-σʹ) > 1e-12 && ++iterations<100); // TV: 'iterate until negligible change in λ' (≈0.006mm)
         if (iterations >= 100) throw new EvalError('Vincenty formula failed to converge'); // not possible?
 
         const x = sinU1*sinσ - cosU1*cosσ*cosα1;
@@ -278,34 +279,32 @@ class LatLonEllipsoidal_Vincenty extends LatLonEllipsoidal {
         let λ = L, sinλ = null, cosλ = null; // λ = difference in longitude on an auxiliary sphere
         let σ = antipodal ? π : 0, sinσ = 0, cosσ = antipodal ? -1 : 1, sinSqσ = null; // σ = angular distance P₁ P₂ on the sphere
         let cos2σₘ = 1;                      // σₘ = angular distance on the sphere from the equator to the midpoint of the line
-        let sinα = null, cosSqα = 1;         // α = azimuth of the geodesic at the equator
-        let C = null;
+        let cosSqα = 1;                      // α = azimuth of the geodesic at the equator
 
         let λʹ = null, iterations = 0;
         do {
             sinλ = Math.sin(λ);
             cosλ = Math.cos(λ);
-            sinSqσ = (cosU2*sinλ) * (cosU2*sinλ) + (cosU1*sinU2-sinU1*cosU2*cosλ) * (cosU1*sinU2-sinU1*cosU2*cosλ);
-            if (Math.abs(sinSqσ) < ε) break;  // co-incident/antipodal points (falls back on λ/σ = L)
+            sinSqσ = (cosU2*sinλ)**2 + (cosU1*sinU2-sinU1*cosU2*cosλ)**2;
+            if (Math.abs(sinSqσ) < 1e-24) break;  // co-incident/antipodal points (σ < ≈0.006mm)
             sinσ = Math.sqrt(sinSqσ);
             cosσ = sinU1*sinU2 + cosU1*cosU2*cosλ;
             σ = Math.atan2(sinσ, cosσ);
-            sinα = cosU1 * cosU2 * sinλ / sinσ;
+            const sinα = cosU1 * cosU2 * sinλ / sinσ;
             cosSqα = 1 - sinα*sinα;
             cos2σₘ = (cosSqα != 0) ? (cosσ - 2*sinU1*sinU2/cosSqα) : 0; // on equatorial line cos²α = 0 (§6)
-            C = f/16*cosSqα*(4+f*(4-3*cosSqα));
+            const C = f/16*cosSqα*(4+f*(4-3*cosSqα));
             λʹ = λ;
             λ = L + (1-C) * f * sinα * (σ + C*sinσ*(cos2σₘ+C*cosσ*(-1+2*cos2σₘ*cos2σₘ)));
             const iterationCheck = antipodal ? Math.abs(λ)-π : Math.abs(λ);
             if (iterationCheck > π) throw new EvalError('λ > π');
-        } while (Math.abs(λ-λʹ) > 1e-12 && ++iterations<1000);
+        } while (Math.abs(λ-λʹ) > 1e-12 && ++iterations<1000); // TV: 'iterate until negligible change in λ' (≈0.006mm)
         if (iterations >= 1000) throw new EvalError('Vincenty formula failed to converge');
 
         const uSq = cosSqα * (a*a - b*b) / (b*b);
         const A = 1 + uSq/16384*(4096+uSq*(-768+uSq*(320-175*uSq)));
         const B = uSq/1024 * (256+uSq*(-128+uSq*(74-47*uSq)));
-        const Δσ = B*sinσ*(cos2σₘ+B/4*(cosσ*(-1+2*cos2σₘ*cos2σₘ)-
-            B/6*cos2σₘ*(-3+4*sinσ*sinσ)*(-3+4*cos2σₘ*cos2σₘ)));
+        const Δσ = B*sinσ*(cos2σₘ+B/4*(cosσ*(-1+2*cos2σₘ*cos2σₘ)-B/6*cos2σₘ*(-3+4*sinσ*sinσ)*(-3+4*cos2σₘ*cos2σₘ)));
 
         const s = b*A*(σ-Δσ); // s = length of the geodesic
 

diff --git a/latlon-ellipsoidal.js b/latlon-ellipsoidal.js
@@ -1,6 +1,8 @@
 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  */
-/* Geodesy tools for an ellipsoidal earth model                       (c) Chris Veness 2005-2019  */
+/* Geodesy tools for an ellipsoidal earth model                       (c) Chris Veness 2005-2022  */
 /*                                                                                   MIT Licence  */
+/* Core class for latlon-ellipsoidal-datum & latlon-ellipsoidal-referenceframe.                   */
+/*                                                                                                */
 /* www.movable-type.co.uk/scripts/latlong-convert-coords.html                                     */
 /* www.movable-type.co.uk/scripts/geodesy-library.html#latlon-ellipsoidal                         */
 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  */
@@ -77,9 +79,9 @@ class LatLonEllipsoidal {
      *   const p = new LatLon(51.47788, -0.00147, 17);
      */
     constructor(lat, lon, height=0) {
-        if (isNaN(lat)) throw new TypeError(`invalid lat ‘${lat}’`);
-        if (isNaN(lon)) throw new TypeError(`invalid lon ‘${lon}’`);
-        if (isNaN(height)) throw new TypeError(`invalid height ‘${height}’`);
+        if (isNaN(lat) || lat == null) throw new TypeError(`invalid lat ‘${lat}’`);
+        if (isNaN(lon) || lon == null) throw new TypeError(`invalid lon ‘${lon}’`);
+        if (isNaN(height) || height == null) throw new TypeError(`invalid height ‘${height}’`);
 
         this._lat = Dms.wrap90(Number(lat));
         this._lon = Dms.wrap180(Number(lon));

diff --git a/latlon-nvector-ellipsoidal.js b/latlon-nvector-ellipsoidal.js
@@ -1,11 +1,10 @@
 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  */
-/* Vector-based ellipsoidal geodetic (latitude/longitude) functions   (c) Chris Veness 2015-2020  */
+/* Vector-based ellipsoidal geodetic (latitude/longitude) functions   (c) Chris Veness 2015-2021  */
 /*                                                                                   MIT Licence  */
 /* www.movable-type.co.uk/scripts/latlong-vectors.html                                            */
 /* www.movable-type.co.uk/scripts/geodesy-library.html#latlon-nvector-ellipsoidal                 */
 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  */
 
-
 import LatLonEllipsoidal, { Cartesian, Vector3d, Dms } from './latlon-ellipsoidal.js';