/**
* \file geodesic.h
* \brief API for the geodesic routines in C
*
* This an implementation in C of the geodesic algorithms described in
* - C. F. F. Karney,
* <a href="https://doi.org/10.1007/s00190-012-0578-z">
* Algorithms for geodesics</a>,
* J. Geodesy <b>87</b>, 43--55 (2013);
* DOI: <a href="https://doi.org/10.1007/s00190-012-0578-z">
* 10.1007/s00190-012-0578-z</a>;
* addenda: <a href="https://geographiclib.sourceforge.io/geod-addenda.html">
* geod-addenda.html</a>.
* .
* The principal advantages of these algorithms over previous ones (e.g.,
* Vincenty, 1975) are
* - accurate to round off for |<i>f</i>| < 1/50;
* - the solution of the inverse problem is always found;
* - differential and integral properties of geodesics are computed.
*
* The shortest path between two points on the ellipsoid at (\e lat1, \e
* lon1) and (\e lat2, \e lon2) is called the geodesic. Its length is
* \e s12 and the geodesic from point 1 to point 2 has forward azimuths
* \e azi1 and \e azi2 at the two end points.
*
* Traditionally two geodesic problems are considered:
* - the direct problem -- given \e lat1, \e lon1, \e s12, and \e azi1,
* determine \e lat2, \e lon2, and \e azi2. This is solved by the function
* geod_direct().
* - the inverse problem -- given \e lat1, \e lon1, and \e lat2, \e lon2,
* determine \e s12, \e azi1, and \e azi2. This is solved by the function
* geod_inverse().
*
* The ellipsoid is specified by its equatorial radius \e a (typically in
* meters) and flattening \e f. The routines are accurate to round off with
* double precision arithmetic provided that |<i>f</i>| < 1/50; for the
* WGS84 ellipsoid, the errors are less than 15 nanometers. (Reasonably
* accurate results are obtained for |<i>f</i>| < 1/5.) For a prolate
* ellipsoid, specify \e f < 0.
*
* The routines also calculate several other quantities of interest
* - \e S12 is the area between the geodesic from point 1 to point 2 and the
* equator; i.e., it is the area, measured counter-clockwise, of the
* quadrilateral with corners (\e lat1,\e lon1), (0,\e lon1), (0,\e lon2),
* and (\e lat2,\e lon2).
* - \e m12, the reduced length of the geodesic is defined such that if
* the initial azimuth is perturbed by \e dazi1 (radians) then the
* second point is displaced by \e m12 \e dazi1 in the direction
* perpendicular to the geodesic. On a curved surface the reduced
* length obeys a symmetry relation, \e m12 + \e m21 = 0. On a flat
* surface, we have \e m12 = \e s12.
* - \e M12 and \e M21 are geodesic scales. If two geodesics are
* parallel at point 1 and separated by a small distance \e dt, then
* they are separated by a distance \e M12 \e dt at point 2. \e M21
* is defined similarly (with the geodesics being parallel to one
* another at point 2). On a flat surface, we have \e M12 = \e M21
* = 1.
* - \e a12 is the arc length on the auxiliary sphere. This is a
* construct for converting the problem to one in spherical
* trigonometry. \e a12 is measured in degrees. The spherical arc
* length from one equator crossing to the next is always 180°.
*
* If points 1, 2, and 3 lie on a single geodesic, then the following
* addition rules hold:
* - \e s13 = \e s12 + \e s23
* - \e a13 = \e a12 + \e a23
* - \e S13 = \e S12 + \e S23
* - \e m13 = \e m12 \e M23 + \e m23 \e M21
* - \e M13 = \e M12 \e M23 − (1 − \e M12 \e M21) \e
* m23 / \e m12
* - \e M31 = \e M32 \e M21 − (1 − \e M23 \e M32) \e
* m12 / \e m23
*
* The shortest distance returned by the solution of the inverse problem is
* (obviously) uniquely defined. However, in a few special cases there are
* multiple azimuths which yield the same shortest distance. Here is a
* catalog of those cases:
* - \e lat1 = −\e lat2 (with neither point at a pole). If \e azi1 = \e
* azi2, the geodesic is unique. Otherwise there are two geodesics and the
* second one is obtained by setting [\e azi1, \e azi2] → [\e azi2, \e
* azi1], [\e M12, \e M21] → [\e M21, \e M12], \e S12 → −\e
* S12. (This occurs when the longitude difference is near ±180°
* for oblate ellipsoids.)
* - \e lon2 = \e lon1 ± 180° (with neither point at a pole). If \e
* azi1 = 0° or ±180°, the geodesic is unique. Otherwise
* there are two geodesics and the second one is obtained by setting [\e
* azi1, \e azi2] → [−\e azi1, −\e azi2], \e S12 →
* −\e S12. (This occurs when \e lat2 is near −\e lat1 for
* prolate ellipsoids.)
* - Points 1 and 2 at opposite poles. There are infinitely many geodesics
* which can be generated by setting [\e azi1, \e azi2] → [\e azi1, \e
* azi2] + [\e d, −\e d], for arbitrary \e d. (For spheres, this
* prescription applies when points 1 and 2 are antipodal.)
* - \e s12 = 0 (coincident points). There are infinitely many geodesics which
* can be generated by setting [\e azi1, \e azi2] → [\e azi1, \e azi2] +
* [\e d, \e d], for arbitrary \e d.
*
* These routines are a simple transcription of the corresponding C++ classes
* in <a href="https://geographiclib.sourceforge.io"> GeographicLib</a>. The
* "class data" is represented by the structs geod_geodesic, geod_geodesicline,
* geod_polygon and pointers to these objects are passed as initial arguments
* to the member functions. Most of the internal comments have been retained.
* However, in the process of transcription some documentation has been lost
* and the documentation for the C++ classes, GeographicLib::Geodesic,
* GeographicLib::GeodesicLine, and GeographicLib::PolygonAreaT, should be
* consulted. The C++ code remains the "reference implementation". Think
* twice about restructuring the internals of the C code since this may make
* porting fixes from the C++ code more difficult.
*
* Copyright (c) Charles Karney (2012-2018) <charles@karney.com> and licensed
* under the MIT/X11 License. For more information, see
* https://geographiclib.sourceforge.io/
*
* This library was distributed with
* <a href="../index.html">GeographicLib</a> 1.49.
**********************************************************************/
#if !defined(GEODESIC_H)
#define GEODESIC_H 1
/**
* The major version of the geodesic library. (This tracks the version of
* GeographicLib.)
**********************************************************************/
#define GEODESIC_VERSION_MAJOR 1
/**
* The minor version of the geodesic library. (This tracks the version of
* GeographicLib.)
**********************************************************************/
#define GEODESIC_VERSION_MINOR 49
/**
* The patch level of the geodesic library. (This tracks the version of
* GeographicLib.)
**********************************************************************/
#define GEODESIC_VERSION_PATCH 3
/**
* Pack the version components into a single integer. Users should not rely on
* this particular packing of the components of the version number; see the
* documentation for GEODESIC_VERSION, below.
**********************************************************************/
#define GEODESIC_VERSION_NUM(a,b,c) ((((a) * 10000 + (b)) * 100) + (c))
/**
* The version of the geodesic library as a single integer, packed as MMmmmmpp
* where MM is the major version, mmmm is the minor version, and pp is the
* patch level. Users should not rely on this particular packing of the
* components of the version number. Instead they should use a test such as
* @code{.c}
#if GEODESIC_VERSION >= GEODESIC_VERSION_NUM(1,40,0)
...
#endif
* @endcode
**********************************************************************/
#define GEODESIC_VERSION \
GEODESIC_VERSION_NUM(GEODESIC_VERSION_MAJOR, \
GEODESIC_VERSION_MINOR, \
GEODESIC_VERSION_PATCH)
#if defin