GeographicLib News

Changes between 2.5 (released 2024-12-28) and 2.4 versions:

• New magnetic models provided

  • World Magnetic Model 2025, wmm2025; this covers the period
    2025-2030. This is now the model returned by
    MagneticModel::DefaultMagneticName and is default magnetic
    model for MagneticField (replacing wmm2020 which is only valid
    thru the end of 2024).
  • World Magnetic Model High Resolution 2025, wmmhr2025; this is a
    high-resolution version of wmm2025 including terms up to degree
    133 (as opposed to 12 for wmm2025).
  • International Geomagnetic Reference Field (14th generation),
    igrf14, which approximates the main magnetic field of the earth
    for the period 1900-2030.... read more

Changes between 2.4 (released 2024-07-14) and 2.3 versions:

• Math::sind, Math::sincosd, etc., now return accurate results for
  arguments = 30d, 45d, 60d, ...

• Add Math::hypot3.

• Add EllipticFunction::am.

• Make constructor for AuxAngle explicit.

• Include notes for plate carree areas.

• C++14 compliant compiler required (this requirement tracks the
  Boost Math library).

• CMake configuration changes

  • update cmake_minimum_required to 3.16.0.
  • exclude the tools and doc targets from the default build (make
    all) if these components won't be installed (BINDIR or DOCDIR
    set to OFF).
  • testprograms is now excluded from the default build and it now
    depends on tools. Before running the test suite, build
    testprograms, e.g., make testprograms test.
  • added a target libs to build only the static and/or shared
    libraries.
  • document the use of CMake's FetchContent as an alternative to
    find_package.
  • stop requiring that CMAKE_CROSSCOMPILING match with
    find_package(GeographicLib).
• Internal code changes
  • stop using reserved identifiers (_N in the DST class and _R in
    the Intersect class).
  • remove the workaround for the boost-quadmath bug with
    setprecision(0).
  • GEOGRAPHICLIB_PANIC(msg) macro now takes a message argument.

Changes between 2.3 (released 2023-07-25) and 2.2 versions:

• Add the Intersect class and the IntersectTool utility. The
  methods are described in "Geodesic intersections",
  https://arxiv.org/abs/2308.00495 (2023-08).

• Add typedefs Geodesic::LineClass, GeodesicExact::LineClass,
  Rhumb::LineClass, GeodesicLine::BaseClass,
  GeodesicLineExact::BaseClass, RhumbLine::BaseClass.

• Geodesic constructor accepts optional third argument exact, default
  false. If true, then the calculations are delegated to
  GeodesicExact. GeodSolve and Planimeter now use this simplified
  interface.... read more

Changes between 2.2 (released 2023-03-07) and 2.1.2 versions:

• Completely redo Rhumb and RhumbLine classes to yield accurate
  results for the course and the area even for very eccentric
  ellipsoids:
  • these methods are accurate for 1/100 < b/a < 100 (the previous
    methods for the course were limited to 1/4 < b/a < 4 for the
    course and the area calculation was only valid for |f| < 0.01);
  • the meaning of the exact flag in the constructor is now more
    straightforward: exact = false (the default) consistently uses
    series expansions valid for |f| < 0.01; exact = true uses exact
    equations for the course of the rhumb line and an accurate DFT
    fit for the area computation;
  • the RhumbSolve(1) utiity accepts the -u to unroll the
    longitudes for the results of direct calculations;
  • rhumb lines which include a pole as one of the endpoints are
    treated properly (previously point at the pole were moved
    slightly away from the pole);
  • this treatment of rhumb areas is described in polygons,
    https://arxiv.org/abs/2303.03219.... read more

Changes between 2.1.2 (released 2022-12-13) and 2.1.1 versions:

• Add MGRS::Decode to break an MGRS string into its components.

• Add definite integral overload for DST::integral. This is used in
  GeodesicExact.

• Add example code examples/AuxLatitude.[hc]pp implementating the
  AuxAngle and AuxLatitude classes. These classes implement the
  methods documented in the paper "On auxiliary latitudes",
  https://arxiv.org/abs/2212.05818. They are not part of
  GeographicLib. See Auxiliary latitudes for more details.... read more

I just posted this on arXiv:

C. F. F. Karney
Geodesics on an arbitrary ellipsoid of revolution
(Aug., 2022)
http://arxiv.org/abs/2208.00492

This documents the algorithms used for the GeodesicExact class and
the -E options for GeodSolve and Planimeter

Changes between 2.1.1 (released 2022-07-25) and 2.1 versions:

• The Planimeter utility used to accept polygon vertices as UTM/UPS
  or MGRS coordinates. However, since these are converted to
  latitude and longitude using the WGS84 ellipoid, this convention
  is incompatible with the -e option to specify the ellipsoid. So
  now, you have to provide the --geoconvert-input option to allow
  vertices to be specified in this way.... read more

There's a sign error in the DST class in GeographicLib 2.1. This affects the DST::refine member function which isn't used elsewhere in GeographicLib. Here's the fix

diff --git a/src/DST.cpp b/src/DST.cpp
index 6fc3ee90..1e70bfd1 100644
--- a/src/DST.cpp
+++ b/src/DST.cpp
@@ -67,7 +67,7 @@ namespace GeographicLib {
     fft_transform(data, F+_N, true);
     for (int i = 0; i < _N; ++i) data[i] = F[i+_N];
     for (int i = _N; i < 2*_N; ++i)
-      F[i] = (-data[2*_N-1-i] + F[2*_N-1-i])/2;
+      F[i] = (data[2*_N-1-i] - F[2*_N-1-i])/2;
     for (int i = 0; i < _N; ++i)
       F[i] = (data[i] + F[i])/2;
   }

Changes between 2.1 (released 2022-06-09) and 2.0 versions:

• Improve the accuracy of the area calculations in GeodesicExact
  and the Planimeter utility with the -E option. This now gives
  full double precision accuracy for b/a in [0.01, 100] or f in
  [-99, 0.99]. (Previously the range for full accuracy was b/a in
  [0.5, 2].) The method involves a direct computation of the
  Fourier series for the integrand of area integral which then
  allows the integral to be evaluated.... read more

Changes between 2.0 (released 2022-05-06) and 1.52 versions:

• Remove non C++ implementations from this package. These are now
  managed as separate packages. See
  https://geographiclib.sourceforge.io/doc/library.html#languages

• This allowed the cmake interface to be rationalized:

  • Replace GEOGRAPHICLIB_LIB_TYPE by BUILD_SHARED_LIBS (a standard
    cmake variable) and BUILD_BOTH_LIBS.
  • Remove COMMON_INSTALL_PATH and replace with variables to
    specify where various components are to be installed.
  • Remove BUILD_NETGEOGRAPHICLIB.
  • Remove legacy cmake config support (non-namespace).
  • make exampleprograms does a separate cmake configuration to
    mimic how a user's program might find GeographicLib.
  • find_package (GeographicLib) now sets and checks the value of
    GEOGRAPHICLIB_PRECISION.
  • Improve the separation of end-user and maintain cmake code.... read more

Changes between 1.52 (released 2021-06-22) and 1.51 versions:

• Add MagneticModel::FieldGeocentric and
  MagneticCircle::FieldGeocentric to return the field in geocentric
  coordinates (thanks to Marcelo Banik de Padua).

• Document realistic errors for PolygonAreaT and Planimeter.

• Geodesic routines: be more aggressive in preventing negative s12
  and m12 for short lines (all languages).... read more

Changes between 1.51 (released 2020-11-22) and 1.50.1 versions:

• C++11 compiler required for C++ library. As a consequence:
  • The workaround implementations for C++11 routines (Math::hypot,
    Math::expm1, Math::log1p, Math::asinh, Math::atanh,
    Math::copysign, Math::cbrt, Math::remainder, Math::remquo,
    Math::round, Math::lround, Math::fma, Math::isfinite, and
    Math::isnan) are now deprecated. Just use the versions in the
    std:: namespace instead.
  • SphericalEngine class, fix the namespace for using streamoff.
  • Some templated functions, e.g., Math::degree(), now have default
    template parameters, T = Math::real.... read more

I've uploaded the International Geomagnetic Reference Field (13th generation), igrf13, magnetic model to

https://sourceforge.net/projects/geographiclib/files/magnetic-distrib/

This can be used with the MagneticField utility and the MagneticModel class. To install it on a Unix system use

sudo geographiclib-get-magnetic -f igrf13

The -f is needed because igrf13 isn't yet known to geographiclib-get-magnetic.

The CMake command

find_package (GeographicLib ...)

(for compiling your software using GeographicLib) fails with CMake versions 3.4 and 3.5. It works with CMake version 3.7. So I recommend using CMake 3.7 or later with GeographicLib.

Changes between 1.50.1 (released 2019-12-13) and 1.50 versions:

• Add the World Magnetic Model 2020, wmm2020, covering the period
  2020-2025. This is now the model returned by
  MagneticModel::DefaultMagneticName and is default magnetic model
  for MagneticField (replacing wmm2015v2 which is only valid thru
  the end of 2019).

• Include float instantiations of those templated Math functions
  which migrated to Math.cpp in version 1.50.... read more

I've added the World Magnetic Model 2020, wmm2020. This can be
downloaded from the magnetic-distrib folder in the download area. No
changes are needed in GeographicLib to use this model. However, I will
shortly be releasing version 1.50.1 in which wmm2020 will replace
wmm2015v2 as the default magnetic model.

Changes between 1.50 (released 2019-09-24) and 1.49 versions:

• BUG fixes:
  • Java + JavaScript implementations of PolygonArea::TestEdge
    counted the pole encirclings wrong.
  • Fix typo in JavaScript implementation which affected unsigned
    areas.
  • Adding edges to a polygon counted pole encirclings inconsistent
    with the way the adding point counted them. This might have
    caused an incorrect result if a polygon vertex had longitude
    = 0. This affected all implementations except Fortran and
    MATLAB).
  • GARS::Forward: fix BUG in handling of longitude = +/-180d.
  • Fix bug in Rhumb class and RhumbSolve(1) utiity which caused
    incorrect area to be reported if an endpoint is at a pole.
    Thanks to Natalia Sabourova for reporting this.
  • Fix bug in MATLAB routine mgrs_inv which resulted in incorrect
    results for UPS zones with prec = -1.
  • In geodreckon.m geoddistance.m, suppress (innocuous) "warning:
    division by zero" messages from Octave.
  • In python implementation, work around problems caused by
    sin(inf) and fmod(inf) raising exceptions.
  • Geoid class, fix the use of streamoff.... read more

User Federico reported a problem with the declaration of streamoff when
compiling GeographicLib with gcc version 8.2 with -std=c++1z. The
following patch fixes this:

diff --git a/include/GeographicLib/Geoid.hpp b/include/GeographicLib/Geoid.hpp
index 3a55adc8..7f41b651 100644
--- a/include/GeographicLib/Geoid.hpp
+++ b/include/GeographicLib/Geoid.hpp
@@ -120,7 +120,7 @@ namespace GeographicLib {
     mutable real _v00, _v01, _v10, _v11;
     mutable real _t[nterms_];
     void filepos(int ix, int iy) const {
-      _file.seekg(std::ios::streamoff
+      _file.seekg(std::streamoff
                   (_datastart +
                    pixel_size_ * (unsigned(iy)*_swidth + unsigned(ix))));
     }... [read more](/p/geographiclib/news/2019/08/fix-to-declaration-of-streamoff-in-geoidhpp/)

There's a BUG in the PolygonArea::TestEdge routine in the Java
implementation of GeographicLib. The following patch fixes this:

diff --git a/java/src/main/java/net/sf/geographiclib/PolygonArea.java b/java/src/main/java/net/sf/geographiclib/PolygonArea.java
index 1e4ed1a1..b7379bcc 100644
--- a/java/src/main/java/net/sf/geographiclib/PolygonArea.java
+++ b/java/src/main/java/net/sf/geographiclib/PolygonArea.java
@@ -338,10 +338,10 @@ public class PolygonArea {
         _earth.Direct(_lat1, _lon1, azi, false, s, _mask);
       tempsum += g.S12;
       crossings += transitdirect(_lon1, g.lon2);
+      crossings += transit(g.lon2, _lon0);
       g = _earth.Inverse(g.lat2, g.lon2, _lat0, _lon0, _mask);
       perimeter += g.s12;
       tempsum += g.S12;
-      crossings += transit(g.lon2, _lon0);
     }... [read more](/p/geographiclib/news/2019/05/bug--fix-in-java-implementation-of-polygonareatestedge/)

The World Magnetic Model wmm2015v2 is now available for use with
GeographicLib. You can download it from

https://sourceforge.net/projects/geographiclib/files/magnetic-distrib/

This "out of cycle" release is described at

https://www.ncei.noaa.gov/news/world-magnetic-model-out-cycle-release
https://www.nytimes.com/2019/02/04/science/north-magnetic-pole-model.html

The original wmm2015 is still available, but this labeled "deprecated"
by NOAA.... read more

Changes between 1.49 (released 2017-10-05) and 1.48 versions:

• Add the Enhanced Magnetic Model 2017, emm2017. This is valid for
  2000 thru the end of 2021.
• Avoid potential problems with the order of initializations in DMS,
  GARS, Geohash, Georef, MGRS, OSGB, SphericalEngine; this only would
  have been an issue if GeographicLib objects were instantiated
  globally. Now no GeographicLib initialization code should be run
  prior to the entry of main().
• To support the previous fix, add an overload, Utility::lookup(const
  char* s, char c).
• NearestNeighbor::Search throws an error if pts is the wrong size
  (instead of merely returning no results).
• Use complex arithmetic for Clenshaw sums in TransverseMercator and
  tranmerc_{fwd,inv}.m.
• Changes in cmake support:
  • fix compiler flags for GEOGRAPHICLIB_PRECISION = 4;
  • add CONVERT_WARNINGS_TO_ERRORS option (default OFF), if ON then
    compiler warnings are treated as errors.
• Fix warnings about implicit conversions of doubles to bools in C++,
  C, and JavaScript packages.
• Binary installers for Windows now use Visual Studio 14 2015.

Changes between 1.48 (released 2017-04-09) and 1.47 versions:

• The "official" URL for GeographicLib is now
  https://geographiclib.sourceforge.io (instead of
  http://geographiclib.sourceforge.net).

• The default range for longitude and azimuth is now (-180d, 180d],
  instead of [-180d, 180d). This was already the case for the C++
  library; now the change has been made to the other implementations
  (C, Fortran, Java, JavaScript, Python, MATLAB, and Maxima).... read more

There's a bug in the MATLAB function geodreckon in calculating the area
when the distance argument is empty. The following patch fixes this.
This will be included in the next release of GeographicLib

diff --git a/matlab/geographiclib/geodreckon.m b/matlab/geographiclib/geodreckon.m
index 258e0e8..7751e7b 100644
--- a/matlab/geographiclib/geodreckon.m
+++ b/matlab/geographiclib/geodreckon.m
@@ -248,11 +248,14 @@ function [lat2, lon2, azi2, S12, m12, M12, M21, a12_s12] = geodreckon ...
                    ssig12 .* (csig1 .* ssig12 ./ (1 + csig12) + ssig1), ...
                    csig12 <= 0);
     calp12 = salp0.^2 + calp0.^2 .* csig1 .* csig2;
-    % Enlarge salp1, calp1 is case lat1 is an array and azi1 is a scalar
-    s = zeros(size(salp0)); salp1 = salp1 + s; calp1 = calp1 + s;
-    s = calp0 == 0 | salp0 == 0;
-    salp12(s) = salp2(s) .* calp1(s) - calp2(s) .* salp1(s);
-    calp12(s) = calp2(s) .* calp1(s) + salp2(s) .* salp1(s);
+    % Deal with geodreckon(10, 0, [], 0) which has calp2 = []
+    if length(Z) > 0
+      % Enlarge salp1, calp1 is case lat1 is an array and azi1 is a scalar
+      s = zeros(size(salp0)); salp1 = salp1 + s; calp1 = calp1 + s;
+      s = calp0 == 0 | salp0 == 0;
+      salp12(s) = salp2(s) .* calp1(s) - calp2(s) .* salp1(s);
+      calp12(s) = calp2(s) .* calp1(s) + salp2(s) .* salp1(s);
+    end
     if e2 ~= 0
       c2 = (a^2 + b^2 * eatanhe(1, e2) / e2) / 2;
     else

There's a bug in the code for reading a NearestNeighbor object. The
following patch fixes this. This will be included in the next release
of GeographicLib

diff --git a/include/GeographicLib/NearestNeighbor.hpp b/include/GeographicLib/NearestNeighbor.hpp
index ebefeb9..2aa91e1 100644
--- a/include/GeographicLib/NearestNeighbor.hpp
+++ b/include/GeographicLib/NearestNeighbor.hpp
@@ -569,7 +569,7 @@ namespace GeographicLib {

       // Sanity check on a Node
       void Check(int numpoints, int treesize, int bucket) const {
-        if (!( -1 <= index && index << numpoints ))
+        if (!( -1 <= index && index < numpoints ))
           throw GeographicLib::GeographicErr("Bad index");
         if (index >= 0) {
           if (!( -1 <= data.child[0] && data.child[0] < treesize &&

The following only affects users of the Python GeographicLib package who
use Python 2.x on Windows 32-bit machines...

On these systems, there's a bug in the handling of -0.0 by math.fmod.
The following patch works aroung this bug. This patch will be
incorporated in the next version of GeographicLib.

diff --git a/python/geographiclib/geomath.py b/python/geographiclib/geomath.py
index 74d2c39..d5a535b 100644
--- a/python/geographiclib/geomath.py
+++ b/python/geographiclib/geomath.py
@@ -129,9 +129,12 @@ class Math(object):
   def AngNormalize(x):
     """reduce angle to [-180,180)"""

-    x = math.fmod(x, 360)
-    return (x + 360 if x < -180 else
-            (x if x < 180 else x - 360))
+    y = math.fmod(x, 360)
+    # On Windows 32-bit with python 2.7, math.fmod(-0.0, 360) = +0.0
+    # This fixes this bug.  See also Math::sincosd in the C++ library.
+    y = x if x == 0 else y
+    return (y + 360 if y < -180 else
+            (y if y < 180 else y - 360))
   AngNormalize = staticmethod(AngNormalize)

   def LatFix(x):
@@ -162,6 +165,9 @@ class Math(object):
       s, c = 0.0-s, 0.0-c
     elif q == 3:
       s, c = 0.0-c,     s
+    # On Windows 32-bit with python 2.7, math.fmod(-0.0, 360) = +0.0
+    # This fixes this bug.  See also Math::sincosd in the C++ library.
+    s = x if x == 0 else s
     return s, c
   sincosd = staticmethod(sincosd)