CassiniSoldner.hpp   CassiniSoldner.hpp 
/** /**
* \file CassiniSoldner.hpp * \file CassiniSoldner.hpp
* \brief Header for GeographicLib::CassiniSoldner class * \brief Header for GeographicLib::CassiniSoldner class
* *
* Copyright (c) Charles Karney (2009, 2010) <charles@karney.com> * Copyright (c) Charles Karney (2009, 2010) <charles@karney.com>
* and licensed under the LGPL. For more information, see * and licensed under the LGPL. For more information, see
* http://geographiclib.sourceforge.net/ * http://geographiclib.sourceforge.net/
**********************************************************************/ **********************************************************************/
#if !defined(GEOGRAPHICLIB_CASSINISOLDNER_HPP) #if !defined(GEOGRAPHICLIB_CASSINISOLDNER_HPP)
#define GEOGRAPHICLIB_CASSINISOLDNER_HPP "$Id: CassiniSoldner.hpp 6785 2010 -01-05 22:15:42Z karney $" #define GEOGRAPHICLIB_CASSINISOLDNER_HPP "$Id: CassiniSoldner.hpp 6827 2010 -05-20 19:56:18Z karney $"
#include "GeographicLib/Geodesic.hpp" #include "GeographicLib/Geodesic.hpp"
#include "GeographicLib/Constants.hpp" #include "GeographicLib/Constants.hpp"
namespace GeographicLib { namespace GeographicLib {
/** /**
* \brief Cassini-Soldner Projection. * \brief Cassini-Soldner Projection.
* *
* Cassini-Soldner projection centered at an arbitrary position, \e lat0, \e * Cassini-Soldner projection centered at an arbitrary position, \e lat0, \e
skipping to change at line 68 skipping to change at line 68
* geodesics between that point and, respectively, (\e lat1, \e lon1) and (\e * geodesics between that point and, respectively, (\e lat1, \e lon1) and (\e
* lat, \e lon). * lat, \e lon).
**********************************************************************/ **********************************************************************/
class CassiniSoldner { class CassiniSoldner {
private: private:
typedef Math::real real; typedef Math::real real;
const Geodesic _earth; const Geodesic _earth;
GeodesicLine _meridian; GeodesicLine _meridian;
real _sbet0, _cbet0; real _sbet0, _cbet0;
real Scale(const GeodesicLine& perp, real sig12) const throw();
static const real eps1, eps2; static const real eps1, eps2;
static const unsigned maxit = 10; static const unsigned maxit = 10;
static inline real sq(real x) throw() { return x * x; } static inline real sq(real x) throw() { return x * x; }
// The following private helper functions are copied from Geodesic.
static inline real AngNormalize(real x) throw() {
// Place angle in [-180, 180). Assumes x is in [-540, 540).
return x >= 180 ? x - 360 : x < -180 ? x + 360 : x;
}
static inline real AngRound(real x) throw() {
// The makes the smallest gap in x = 1/16 - nextafter(1/16, 0) = 1/2^
57
// for reals = 0.7 pm on the earth if x is an angle in degrees. (Thi
s
// is about 1000 times more resolution than we get with angles around
90
// degrees.) We use this to avoid having to deal with near singular
// cases when x is non-zero but tiny (e.g., 1.0e-200).
const real z = real(0.0625); // 1/16
volatile real y = std::abs(x);
// The compiler mustn't "simplify" z - (z - y) to y
y = y < z ? z - (z - y) : y;
return x < 0 ? -y : y;
}
static inline void SinCosNorm(real& sinx, real& cosx) throw() {
real r = Math::hypot(sinx, cosx);
sinx /= r;
cosx /= r;
}
public: public:
/** /**
* Constructor for CassiniSoldner setting the Geodesic object, \e earth , to * Constructor for CassiniSoldner setting the Geodesic object, \e earth , to
* use for geodesic calculations. By default this uses the WGS84 * use for geodesic calculations. By default this uses the WGS84
* ellipsoid. This constructor makes an "uninitialized" object. Call * ellipsoid. This constructor makes an "uninitialized" object. Call
* Reset to set the central latitude and longuitude, prior to calling * Reset to set the central latitude and longuitude, prior to calling
* Forward and Reverse. * Forward and Reverse.
********************************************************************** / ********************************************************************** /
explicit CassiniSoldner(const Geodesic& earth = Geodesic::WGS84) throw( ) explicit CassiniSoldner(const Geodesic& earth = Geodesic::WGS84) throw( )
 End of changes. 3 change blocks. 
2 lines changed or deleted 26 lines changed or added


 Constants.hpp   Constants.hpp 
/** /**
* \file Constants.hpp * \file Constants.hpp
* \brief Header for GeographicLib::Constants class * \brief Header for GeographicLib::Constants class
* *
* Copyright (c) Charles Karney (2008, 2009, 2010) <charles@karney.com> * Copyright (c) Charles Karney (2008, 2009, 2010) <charles@karney.com>
* and licensed under the LGPL. For more information, see * and licensed under the LGPL. For more information, see
* http://geographiclib.sourceforge.net/ * http://geographiclib.sourceforge.net/
**********************************************************************/ **********************************************************************/
#if !defined(GEOGRAPHICLIB_CONSTANTS_HPP) #if !defined(GEOGRAPHICLIB_CONSTANTS_HPP)
#define GEOGRAPHICLIB_CONSTANTS_HPP "$Id: Constants.hpp 6817 2010-02-08 14: 42:35Z karney $" #define GEOGRAPHICLIB_CONSTANTS_HPP "$Id: Constants.hpp 6827 2010-05-20 19: 56:18Z karney $"
/** /**
* A simple compile-time assert. This is designed to be compatible with th e * A simple compile-time assert. This is designed to be compatible with th e
* C++0X static_assert. * C++0X static_assert.
**********************************************************************/ **********************************************************************/
#if !defined(STATIC_ASSERT) #if !defined(STATIC_ASSERT)
#define STATIC_ASSERT(cond,reason) { enum{ STATIC_ASSERT_ENUM=1/int(cond) } ; } #define STATIC_ASSERT(cond,reason) { enum{ STATIC_ASSERT_ENUM=1/int(cond) } ; }
#endif #endif
#if defined(__GNUC__) #if defined(__GNUC__)
skipping to change at line 114 skipping to change at line 114
static inline double hypot(double x, double y) throw() static inline double hypot(double x, double y) throw()
{ return _hypot(x, y); } { return _hypot(x, y); }
static inline float hypot(float x, float y) throw() static inline float hypot(float x, float y) throw()
{ return _hypotf(x, y); } { return _hypotf(x, y); }
#else #else
// Use overloading to define generic versions // Use overloading to define generic versions
static inline double hypot(double x, double y) throw() static inline double hypot(double x, double y) throw()
{ return ::hypot(x, y); } { return ::hypot(x, y); }
static inline float hypot(float x, float y) throw() static inline float hypot(float x, float y) throw()
{ return ::hypotf(x, y); } { return ::hypotf(x, y); }
#if !defined(__NO_LONG_DOUBLE_MATH)
static inline long double hypot(long double x, long double y) throw() static inline long double hypot(long double x, long double y) throw()
{ return ::hypotl(x, y); } { return ::hypotl(x, y); }
#endif #endif
#endif
#if defined(DOXYGEN) || defined(_MSC_VER) #if defined(DOXYGEN) || defined(_MSC_VER)
/** /**
* Return exp(\e x) - 1 accurate near \e x = 0. This is taken from * Return exp(\e x) - 1 accurate near \e x = 0. This is taken from
* N. J. Higham, Accuracy and Stability of Numerical Algorithms, 2nd * N. J. Higham, Accuracy and Stability of Numerical Algorithms, 2nd
* Edition (SIAM, 2002), Sec 1.14.1, p 19. * Edition (SIAM, 2002), Sec 1.14.1, p 19.
********************************************************************** / ********************************************************************** /
static inline real expm1(real x) throw() { static inline real expm1(real x) throw() {
volatile real volatile real
y = std::exp(x), y = std::exp(x),
z = y - 1; z = y - 1;
// The reasoning here is similar to that for log1p. The expression // The reasoning here is similar to that for log1p. The expression
// mathematically reduces to exp(x) - 1, and the factor z/log(y) = (y - // mathematically reduces to exp(x) - 1, and the factor z/log(y) = (y -
// 1)/log(y) is a slowly varying quantity near y = 1 and is accuratel y // 1)/log(y) is a slowly varying quantity near y = 1 and is accuratel y
// computed. // computed.
return std::abs(x) > 1 ? z : z == 0 ? x : x * z / std::log(y); return std::abs(x) > 1 ? z : z == 0 ? x : x * z / std::log(y);
} }
#else #else
static inline double expm1(double x) throw() { return ::expm1(x); } static inline double expm1(double x) throw() { return ::expm1(x); }
static inline float expm1(float x) throw() { return ::expm1f(x); } static inline float expm1(float x) throw() { return ::expm1f(x); }
#if !defined(__NO_LONG_DOUBLE_MATH)
static inline long double expm1(long double x) throw() static inline long double expm1(long double x) throw()
{ return ::expm1l(x); } { return ::expm1l(x); }
#endif #endif
#endif
#if defined(DOXYGEN) || defined(_MSC_VER) #if defined(DOXYGEN) || defined(_MSC_VER)
/** /**
* Return log(\e x + 1) accurate near \e x = 0. This is taken See * Return log(\e x + 1) accurate near \e x = 0. This is taken See
* D. Goldberg, * D. Goldberg,
* <a href="http://docs.sun.com/source/806-3568/ncg_goldberg.html"> Wha t * <a href="http://docs.sun.com/source/806-3568/ncg_goldberg.html"> Wha t
* every computer scientist should know about floating-point arithmetic </a> * every computer scientist should know about floating-point arithmetic </a>
* (1991), Theorem 4. See also, Higham (op. cit.), Answer to Problem 1 .5, * (1991), Theorem 4. See also, Higham (op. cit.), Answer to Problem 1 .5,
* p 528. * p 528.
********************************************************************** / ********************************************************************** /
skipping to change at line 163 skipping to change at line 167
z = y - 1; z = y - 1;
// Here's the explanation for this magic: y = 1 + z, exactly, and z // Here's the explanation for this magic: y = 1 + z, exactly, and z
// approx x, thus log(y)/z (which is nearly constant near z = 0) retu rns // approx x, thus log(y)/z (which is nearly constant near z = 0) retu rns
// a good approximation to the true log(1 + x)/x. The multiplication x * // a good approximation to the true log(1 + x)/x. The multiplication x *
// (log(y)/z) introduces little additional error. // (log(y)/z) introduces little additional error.
return z == 0 ? x : x * std::log(y) / z; return z == 0 ? x : x * std::log(y) / z;
} }
#else #else
static inline double log1p(double x) throw() { return ::log1p(x); } static inline double log1p(double x) throw() { return ::log1p(x); }
static inline float log1p(float x) throw() { return ::log1pf(x); } static inline float log1p(float x) throw() { return ::log1pf(x); }
#if !defined(__NO_LONG_DOUBLE_MATH)
static inline long double log1p(long double x) throw() static inline long double log1p(long double x) throw()
{ return ::log1pl(x); } { return ::log1pl(x); }
#endif #endif
#endif
#if defined(DOXYGEN) || defined(_MSC_VER) #if defined(DOXYGEN) || defined(_MSC_VER)
/** /**
* Return asinh(\e x). This is defined in terms of log1p(\e x) in orde * Return asinh(\e x). This is defined in terms of Math::log1p(\e x) i
r to n
* maintain accuracy near \e x = 0. In addition, the odd parity of the * order to maintain accuracy near \e x = 0. In addition, the odd pari
* function is enforced. ty
* of the function is enforced.
********************************************************************** / ********************************************************************** /
static inline real asinh(real x) throw() { static inline real asinh(real x) throw() {
real y = std::abs(x); // Enforce odd parity real y = std::abs(x); // Enforce odd parity
y = log1p(y * (1 + y/(hypot(real(1), y) + 1))); y = log1p(y * (1 + y/(hypot(real(1), y) + 1)));
return x < 0 ? -y : y; return x < 0 ? -y : y;
} }
#else #else
static inline double asinh(double x) throw() { return ::asinh(x); } static inline double asinh(double x) throw() { return ::asinh(x); }
static inline float asinh(float x) throw() { return ::asinhf(x); } static inline float asinh(float x) throw() { return ::asinhf(x); }
#if !defined(__NO_LONG_DOUBLE_MATH)
static inline long double asinh(long double x) throw() static inline long double asinh(long double x) throw()
{ return ::asinhl(x); } { return ::asinhl(x); }
#endif #endif
#endif
#if defined(DOXYGEN) || defined(_MSC_VER) #if defined(DOXYGEN) || defined(_MSC_VER)
/** /**
* Return atanh(\e x). This is defined in terms of log1p(\e x) in orde * Return atanh(\e x). This is defined in terms of Math::log1p(\e x) i
r to n
* maintain accuracy near \e x = 0. In addition, the odd parity of the * order to maintain accuracy near \e x = 0. In addition, the odd pari
* function is enforced. ty
* of the function is enforced.
********************************************************************** / ********************************************************************** /
static inline real atanh(real x) throw() { static inline real atanh(real x) throw() {
real y = std::abs(x); // Enforce odd parity real y = std::abs(x); // Enforce odd parity
y = log1p(2 * y/(1 - y))/2; y = log1p(2 * y/(1 - y))/2;
return x < 0 ? -y : y; return x < 0 ? -y : y;
} }
#else #else
static inline double atanh(double x) throw() { return ::atanh(x); } static inline double atanh(double x) throw() { return ::atanh(x); }
static inline float atanh(float x) throw() { return ::atanhf(x); } static inline float atanh(float x) throw() { return ::atanhf(x); }
#if !defined(__NO_LONG_DOUBLE_MATH)
static inline long double atanh(long double x) throw() static inline long double atanh(long double x) throw()
{ return ::atanhl(x); } { return ::atanhl(x); }
#endif #endif
#endif
#if defined(DOXYGEN) || defined(_MSC_VER) #if defined(DOXYGEN) || defined(_MSC_VER)
/** /**
* Return the real cube root of \e x. * Return the real cube root of \e x.
********************************************************************** / ********************************************************************** /
static inline real cbrt(real x) throw() { static inline real cbrt(real x) throw() {
real y = std::pow(std::abs(x), 1/real(3)); // Return the real cube ro ot real y = std::pow(std::abs(x), 1/real(3)); // Return the real cube ro ot
return x < 0 ? -y : y; return x < 0 ? -y : y;
} }
#else #else
static inline double cbrt(double x) throw() { return ::cbrt(x); } static inline double cbrt(double x) throw() { return ::cbrt(x); }
static inline float cbrt(float x) throw() { return ::cbrtf(x); } static inline float cbrt(float x) throw() { return ::cbrtf(x); }
#if !defined(__NO_LONG_DOUBLE_MATH)
static inline long double cbrt(long double x) throw() { return ::cbrtl( x); } static inline long double cbrt(long double x) throw() { return ::cbrtl( x); }
#endif #endif
#endif
#if defined(DOXYGEN)
/** /**
* Return true if number is finite, false if NaN or infinite. * Return true if number is finite, false if NaN or infinite.
********************************************************************** / ********************************************************************** /
static inline bool isfinite(real x) throw() { static inline bool isfinite(real x) throw() {
return std::abs(x) < std::numeric_limits<real>::max(); #if defined(DOXYGEN)
} return std::abs(x) <= std::numeric_limits<real>::max();
#elif defined(_MSC_VER) #elif defined(_MSC_VER)
static inline real isfinite(real x) throw() { return _finite(x); } return _finite(x) != 0;
#else #else
static inline real isfinite(real x) throw() { return std::isfinite(x); } return std::isfinite(x) != 0;
#endif #endif
}
}; };
/** /**
* \brief %Constants needed by %GeographicLib * \brief %Constants needed by %GeographicLib
* *
* Define constants specifying the WGS84 ellipsoid, the UTM and UPS * Define constants specifying the WGS84 ellipsoid, the UTM and UPS
* projections, and various unit conversions. * projections, and various unit conversions.
**********************************************************************/ **********************************************************************/
class Constants { class Constants {
private: private:
 End of changes. 20 change blocks. 
14 lines changed or deleted 28 lines changed or added


 DMS.hpp   DMS.hpp 
/** /**
* \file DMS.hpp * \file DMS.hpp
* \brief Header for GeographicLib::DMS class * \brief Header for GeographicLib::DMS class
* *
* Copyright (c) Charles Karney (2008, 2009, 2010) <charles@karney.com> * Copyright (c) Charles Karney (2008, 2009, 2010) <charles@karney.com>
* and licensed under the LGPL. For more information, see * and licensed under the LGPL. For more information, see
* http://geographiclib.sourceforge.net/ * http://geographiclib.sourceforge.net/
**********************************************************************/ **********************************************************************/
#if !defined(GEOGRAPHICLIB_DMS_HPP) #if !defined(GEOGRAPHICLIB_DMS_HPP)
#define GEOGRAPHICLIB_DMS_HPP "$Id: DMS.hpp 6785 2010-01-05 22:15:42Z karne y $" #define GEOGRAPHICLIB_DMS_HPP "$Id: DMS.hpp 6827 2010-05-20 19:56:18Z karne y $"
#include "GeographicLib/Constants.hpp" #include "GeographicLib/Constants.hpp"
#include <sstream> #include <sstream>
#include <iomanip>
namespace GeographicLib { namespace GeographicLib {
/** /**
* \brief Convert between degrees and %DMS representation. * \brief Convert between degrees and %DMS representation.
* *
* Parse a string representing degree, minutes, and seconds and return th e * Parse a string representing degree, minutes, and seconds and return th e
* angle in degrees and format an angle in degrees as degree, minutes, an d * angle in degrees and format an angle in degrees as degree, minutes, an d
* seconds. * seconds.
**********************************************************************/ **********************************************************************/
skipping to change at line 46 skipping to change at line 47
static const std::string signs; static const std::string signs;
static const std::string digits; static const std::string digits;
static const std::string dmsindicators; static const std::string dmsindicators;
static const std::string components[3]; static const std::string components[3];
DMS(); // Disable constructor DMS(); // Disable constructor
public: public:
/** /**
* Indicator for presence of hemisphere indicator (N/S/E/W) on latitude s * Indicator for presence of hemisphere indicator (N/S/E/W) on latitude s
* and longitudes. AZIMUTH is used in Encode to indicate output in [00 * and longitudes. DMS::AZIMUTH is used in Encode to indicate output i
0, n
* 360) with no letter indicator. * [000, 360) with no letter indicator. DMS::NUMBER is used in Encode
to
* indicate the output of a plain number.
********************************************************************** / ********************************************************************** /
enum flag { NONE = 0, LATITUDE = 1, LONGITUDE = 2, AZIMUTH = 3 }; enum flag { NONE = 0, LATITUDE = 1, LONGITUDE = 2, AZIMUTH = 3,
NUMBER = 4, };
/** /**
* Indicator for trailing units on an angle. * Indicator for trailing units on an angle.
********************************************************************** / ********************************************************************** /
enum component { DEGREE = 0, MINUTE = 1, SECOND = 2 }; enum component { DEGREE = 0, MINUTE = 1, SECOND = 2 };
/** /**
* Read a string \e dms in DMS format and return the resulting angle in * Read a string \e dms in DMS format and return the resulting angle in
* degrees. Degrees, minutes, and seconds are indicated by the letters d, * degrees. Degrees, minutes, and seconds are indicated by the letters d,
* ', ", and these components may only be given in this order. Any (bu * ', &quot;, and these components may only be given in this order. An
t y
* not all) components may be omitted. The last component indicator ma * (but not all) components may be omitted. The last component indicat
y be or
* omitted and is assumed to be tbe next smallest unit (thus 33d10 is * may be omitted and is assumed to be tbe next smallest unit (thus 33d
* interpreted as 33d10'). The final component may be a decimal fracti 10
on * is interpreted as 33d10'). The final component may be a decimal
* but the non-final components must be integers. The integer parts of * fraction but the non-final components must be integers. The integer
the * parts of the minutes and seconds components must be less than 60. A
* minutes and seconds components must be less than 60. A single leadi * single leading sign is permitted. A hemisphere designator (N, E, W,
ng S)
* sign is permitted. A hemisphere designator (N, E, W, S) may be adde * may be added to tbe beginning or end of the string. The result is
d to * multiplied by the implied signed of the hemisphere designator (negat
* tbe beginning or end of the string. The result is multiplied by the ive
* implied signed of the hemisphere designator (negative for S and W). * for S and W). In addition \e flag is used to indicate whether such
In a
* addition \e flag is used to indicate whether such a designator was f * designator was found and whether it implies that the angle is a lati
ound tude
* and whether it implies that the angle is a latitude (N or S) or * (N or S) or longitude (E or W). Throws an error on a malformed stri
* longitude (E or W). Throws an error on a malformed string. No chec ng.
k is * No check is performed on the range of the result.
* performed on the range of the result.
********************************************************************** / ********************************************************************** /
static Math::real Decode(const std::string& dms, flag& ind); static Math::real Decode(const std::string& dms, flag& ind);
/** /**
* Convert real degrees \e d, minutes \e m, and seconds \e s, to degree s. * Convert real degrees \e d, minutes \e m, and seconds \e s, to degree s.
* This does not propagate the sign on \e d to the other components, so * This does not propagate the sign on \e d to the other components, so
* -3d20' would need to be represented as - DMS::Decode(3.0, 20.0) or * -3d20' would need to be represented as - DMS::Decode(3.0, 20.0) or
* DMS::Decode(-3.0, -20.0). * DMS::Decode(-3.0, -20.0).
********************************************************************** / ********************************************************************** /
static Math::real Decode(real d, real m = 0, real s = 0) throw() static Math::real Decode(real d, real m = 0, real s = 0) throw()
{ return d + (m + s/real(60))/real(60); } { return d + (m + s/real(60))/real(60); }
/** /**
* Convert a string \e str to a real number.
**********************************************************************
/
static Math::real Decode(const std::string& str);
/**
* Convert two strings \e dmsa and \e dmsb to a latitude, \e lat, and * Convert two strings \e dmsa and \e dmsb to a latitude, \e lat, and
* longitude, \e lon. By default, the \e lat (resp., \e lon) is assign ed * longitude, \e lon. By default, the \e lat (resp., \e lon) is assign ed
* to the results of decoding \e dmsa (resp., \e dmsb). However this i s * to the results of decoding \e dmsa (resp., \e dmsb). However this i s
* overridden if either \e dmsa or \e dmsb contain a latitude or longit ude * overridden if either \e dmsa or \e dmsb contain a latitude or longit ude
* hemisphere designator (N, S, E, W). Throws an error if the decoded * hemisphere designator (N, S, E, W). Throws an error if the decoded
* numbers are out of the ranges [-90<sup>o</sup>, 90<sup>o</sup>] for * numbers are out of the ranges [-90<sup>o</sup>, 90<sup>o</sup>] for
* latitude and [-180<sup>o</sup>, 360<sup>o</sup>] for longitude and, in * latitude and [-180<sup>o</sup>, 360<sup>o</sup>] for longitude and, in
* which case \e lat and \e lon are unchanged. Finally the longitude i s * which case \e lat and \e lon are unchanged. Finally the longitude i s
* reduced to the range [-180<sup>o</sup>, 180<sup>o</sup>). * reduced to the range [-180<sup>o</sup>, 180<sup>o</sup>).
********************************************************************** / ********************************************************************** /
static void DecodeLatLon(const std::string& dmsa, const std::string& dm sb, static void DecodeLatLon(const std::string& dmsa, const std::string& dm sb,
real& lat, real& lon); real& lat, real& lon);
/** /**
* Convert \e degree into a DMS string. \e trailing indicates the leas * Convert a string \e angstr to an angle in degrees. No hemisphere
t * designator is allowed and no check is done on the range of the resul
* significant component of the string (and this component is given as t.
a **********************************************************************
* decimal number if necessary). \e prec indicates the number of digit /
s static Math::real DecodeAngle(const std::string& angstr);
* after the decimal point for the trailing component. \e flag indicat
es /**
* additional formating as follows * Convert a string \e azistr to an azimuth in degrees. A hemisphere
* - flag == NONE, signed result no leading zeros on degrees except in * designator E/W can be used; the result is multiplied by -1 if W is
the * present. Throws an error if the result is out of the range
* units place, e.g., -8d03'. * [-180<sup>o</sup>, 360<sup>o</sup>]. Finally the azimuth is reduced
* - flag == LATITUDE, trailing N or S hemisphere designator, no sign, to
pad * the range [-180<sup>o</sup>, 180<sup>o</sup>).
* degress to 2 digits, e.g., 08d03'S. **********************************************************************
* - flag == LONGITUDE, trailing E or W hemisphere designator, no sign, /
pad static Math::real DecodeAzimuth(const std::string& azistr);
* degress to 3 digits, e.g., 008d03'W.
* - flag == AZIMUTH, convert to the range [0, 360<sup>o</sup>), no sig /**
n, * Convert \e angle (in degrees) into a DMS string. \e trailing indica
* pad degrees to 3 digits, , e.g., 351d57'. tes
* the least significant component of the string (and this component is
* given as a decimal number if necessary). \e prec indicates the numb
er
* of digits after the decimal point for the trailing component. \e fl
ag
* indicates additional formating as follows
* - flag == DMS::NONE, signed result no leading zeros on degrees excep
t in
* the units place, e.g., -8d03'.
* - flag == DMS::LATITUDE, trailing N or S hemisphere designator, no s
ign,
* pad degress to 2 digits, e.g., 08d03'S.
* - flag == DMS::LONGITUDE, trailing E or W hemisphere designator, no
* sign, pad degress to 3 digits, e.g., 008d03'W.
* - flag == DMS::AZIMUTH, convert to the range [0, 360<sup>o</sup>), n
o
* sign, pad degrees to 3 digits, , e.g., 351d57'.
* . * .
* The integer parts of the minutes and seconds components are always g iven * The integer parts of the minutes and seconds components are always g iven
* with 2 digits. * with 2 digits.
********************************************************************** / ********************************************************************** /
static std::string Encode(real degree, component trailing, unsigned pre c, static std::string Encode(real angle, component trailing, unsigned prec ,
flag ind = NONE); flag ind = NONE);
/** /**
* Convert \e degree into a DMS string selecting the trailing component * Convert \e angle into a DMS string selecting the trailing component
* based on \e prec. \e prec indicates the precision relative to 1 deg ree, * based on \e prec. \e prec indicates the precision relative to 1 deg ree,
* e.g., \e prec = 3 gives a result accurate to 0.1' and \e prec = 4 gi ves * e.g., \e prec = 3 gives a result accurate to 0.1' and \e prec = 4 gi ves
* a result accurate to 1". * a result accurate to 1&quot;. If \e ind is DMS::NUMBER, then merely
* format \e angle as a number in fixed format with precision \e prec.
********************************************************************** / ********************************************************************** /
static std::string Encode(real degree, unsigned prec, flag ind = NONE) static std::string Encode(real angle, unsigned prec, flag ind = NONE) {
{ if (ind == NUMBER) {
return Encode(degree, std::ostringstream s;
prec < 2 ? DEGREE : (prec < 4 ? MINUTE : SECOND), s << std::fixed << std::setprecision(prec) << angle;
prec < 2 ? prec : (prec < 4 ? prec - 2 : prec - 4), return s.str();
ind); } else
return Encode(angle,
prec < 2 ? DEGREE : (prec < 4 ? MINUTE : SECOND),
prec < 2 ? prec : (prec < 4 ? prec - 2 : prec - 4),
ind);
} }
/** /**
* Split angle, \e ang, into degrees, \e d, and minutes \e m. * Split angle, \e ang, into degrees, \e d, and minutes \e m.
********************************************************************** / ********************************************************************** /
static void Encode(real ang, real& d, real& m) throw() { static void Encode(real ang, real& d, real& m) throw() {
d = int(ang); m = 60 * (ang - d); d = int(ang); m = 60 * (ang - d);
} }
/** /**
 End of changes. 11 change blocks. 
57 lines changed or deleted 88 lines changed or added


 Geocentric.hpp   Geocentric.hpp 
/** /**
* \file Geocentric.hpp * \file Geocentric.hpp
* \brief Header for GeographicLib::Geocentric class * \brief Header for GeographicLib::Geocentric class
* *
* Copyright (c) Charles Karney (2008, 2009, 2010) <charles@karney.com> * Copyright (c) Charles Karney (2008, 2009, 2010) <charles@karney.com>
* and licensed under the LGPL. For more information, see * and licensed under the LGPL. For more information, see
* http://geographiclib.sourceforge.net/ * http://geographiclib.sourceforge.net/
**********************************************************************/ **********************************************************************/
#if !defined(GEOGRAPHICLIB_GEOCENTRIC_HPP) #if !defined(GEOGRAPHICLIB_GEOCENTRIC_HPP)
#define GEOGRAPHICLIB_GEOCENTRIC_HPP "$Id: Geocentric.hpp 6785 2010-01-05 2 2:15:42Z karney $" #define GEOGRAPHICLIB_GEOCENTRIC_HPP "$Id: Geocentric.hpp 6827 2010-05-20 1 9:56:18Z karney $"
#include "GeographicLib/Constants.hpp" #include "GeographicLib/Constants.hpp"
namespace GeographicLib { namespace GeographicLib {
/** /**
* \brief %Geocentric coordinates * \brief %Geocentric coordinates
* *
* Convert between geodetic coordinates latitude = \e lat, longitude = \e * Convert between geodetic coordinates latitude = \e lat, longitude = \e
* lon, height = \e h (measured vertically from the surface of the ellips oid) * lon, height = \e h (measured vertically from the surface of the ellips oid)
skipping to change at line 48 skipping to change at line 48
* *
* The errors in these routines are close to round-off. Specifically, fo r * The errors in these routines are close to round-off. Specifically, fo r
* points within 5000 km of the surface of the ellipsoid (either inside o r * points within 5000 km of the surface of the ellipsoid (either inside o r
* outside the ellipsoid), the error is bounded by 7 nm for the WGS84 * outside the ellipsoid), the error is bounded by 7 nm for the WGS84
* ellipsoid. See \ref geocentric for further information on the errors. * ellipsoid. See \ref geocentric for further information on the errors.
**********************************************************************/ **********************************************************************/
class Geocentric { class Geocentric {
private: private:
typedef Math::real real; typedef Math::real real;
const real _a, _r, _f, _e2, _e2m, _ax, _e2x, _e4x, _e2mx, _maxrad; const real _a, _r, _f, _e2, _e2m, _e2a, _e4a, _maxrad;
static inline real sq(real x) throw() { return x * x; } static inline real sq(real x) throw() { return x * x; }
public: public:
/** /**
* Constructor for a ellipsoid radius \e a (meters) and reciprocal * Constructor for a ellipsoid radius \e a (meters) and reciprocal
* flattening \e r. Setting \e r = 0 implies \e r = inf or flattening = 0 * flattening \e r. Setting \e r = 0 implies \e r = inf or flattening = 0
* (i.e., a sphere). Negative \e r indicates a prolate spheroid. An * (i.e., a sphere). Negative \e r indicates a prolate spheroid. An
* exception is thrown if \e a is not positive. * exception is thrown if \e a is not positive.
********************************************************************** / ********************************************************************** /
Geocentric(real a, real r); Geocentric(real a, real r);
 End of changes. 2 change blocks. 
2 lines changed or deleted 2 lines changed or added


 Geodesic.hpp   Geodesic.hpp 
/** /**
* \file Geodesic.hpp * \file Geodesic.hpp
* \brief Header for GeographicLib::Geodesic class * \brief Header for GeographicLib::Geodesic class
* *
* Copyright (c) Charles Karney (2009, 2010) <charles@karney.com> * Copyright (c) Charles Karney (2009, 2010) <charles@karney.com>
* and licensed under the LGPL. For more information, see * and licensed under the LGPL. For more information, see
* http://geographiclib.sourceforge.net/ * http://geographiclib.sourceforge.net/
**********************************************************************/ **********************************************************************/
#if !defined(GEOGRAPHICLIB_GEODESIC_HPP) #if !defined(GEOGRAPHICLIB_GEODESIC_HPP)
#define GEOGRAPHICLIB_GEODESIC_HPP "$Id: Geodesic.hpp 6785 2010-01-05 22:15 :42Z karney $" #define GEOGRAPHICLIB_GEODESIC_HPP "$Id: Geodesic.hpp 6827 2010-05-20 19:56 :18Z karney $"
#if !defined(GEOD_TAU_ORD) #include "GeographicLib/Constants.hpp"
/**
* The order of the series relating \e s and \e sigma when expanded in powe
rs
* of \e k1.
**********************************************************************/
#define GEOD_TAU_ORD \
(sizeof(real) == sizeof(double) ? 5 : sizeof(real) == sizeof(float) ? 4 : 6
)
#endif
#if !defined(GEOD_ETA_ORD) #if !defined(GEOD_ORD)
/** /**
* The order of the series relating \e lambda and \e eta when expanded in * The order of the expansions used by Geodesic.
* powers of \e fp.
**********************************************************************/ **********************************************************************/
#define GEOD_ETA_ORD \ #define GEOD_ORD (GEOGRAPHICLIB_PREC == 1 ? 6 : GEOGRAPHICLIB_PREC == 0 ? 3
(sizeof(real) == sizeof(double) ? 5 : sizeof(real) == sizeof(float) ? 4 : 6 : 7)
)
#endif #endif
#include "GeographicLib/Constants.hpp"
namespace GeographicLib { namespace GeographicLib {
class GeodesicLine; class GeodesicLine;
/** /**
* \brief %Geodesic calculations * \brief %Geodesic calculations
* *
* The shortest path between two points on a ellipsoid at (\e lat1, \e lo n1) * The shortest path between two points on a ellipsoid at (\e lat1, \e lo n1)
* and (\e lat2, \e lon2) is called the geodesic. Its length is \e s12 a nd * and (\e lat2, \e lon2) is called the geodesic. Its length is \e s12 a nd
* the geodesic from point 1 to point 2 has azimuths \e azi1 and \e azi2 at * the geodesic from point 1 to point 2 has azimuths \e azi1 and \e azi2 at
skipping to change at line 86 skipping to change at line 75
* The calculations are accurate to better than 15 nm. (See \ref geoderr ors * The calculations are accurate to better than 15 nm. (See \ref geoderr ors
* for details.) * for details.)
* *
* For more information on geodesics see \ref geodesic. * For more information on geodesics see \ref geodesic.
**********************************************************************/ **********************************************************************/
class Geodesic { class Geodesic {
private: private:
typedef Math::real real; typedef Math::real real;
friend class GeodesicLine; friend class GeodesicLine;
friend class CassiniSoldner; static const int nA1 = GEOD_ORD, nC1 = GEOD_ORD, nC1p = GEOD_ORD,
static const int tauord = GEOD_TAU_ORD; nA2 = GEOD_ORD, nC2 = GEOD_ORD, nA3 = GEOD_ORD, nC3 = GEOD_ORD;
static const int ntau = tauord;
static const int nsig = tauord;
static const int zetord = GEOD_TAU_ORD;
static const int nzet = zetord;
static const int etaord = GEOD_ETA_ORD;
// etaCoeff is multiplied by etaFactor which is O(f), so we reduce the
// order to which etaCoeff is computed by 1.
static const int neta = etaord > 0 ? etaord - 1 : 0;
static const unsigned maxit = 50; static const unsigned maxit = 50;
static inline real sq(real x) throw() { return x * x; } static inline real sq(real x) throw() { return x * x; }
void Lengths(real k1, real sig12, void Lengths(real eps, real sig12,
real ssig1, real csig1, real ssig2, real csig2, real ssig1, real csig1, real ssig2, real csig2,
real cbet1, real cbet2, real cbet1, real cbet2,
real& s12s, real& m12a, real& m0, real& s12s, real& m12a, real& m0,
real tc[], real zc[]) const throw(); real tc[], real zc[]) const throw();
static void Evolute(real R, real z, real& c, real& s) throw(); static real Astroid(real R, real z) throw();
void InverseStart(real sbet1, real cbet1, real sbet2, real cbet2, real InverseStart(real sbet1, real cbet1, real sbet2, real cbet2,
real lam12, real slam12, real clam12, real lam12,
real& salp1, real& calp1, real& salp1, real& calp1,
real tc[], real zc[]) const throw(); real& salp2, real& calp2,
real C1a[], real C2a[]) const throw();
real Lambda12(real sbet1, real cbet1, real sbet2, real cbet2, real Lambda12(real sbet1, real cbet1, real sbet2, real cbet2,
real salp1, real calp1, real salp1, real calp1,
real& salp2, real& calp2, real& sig12, real& salp2, real& calp2, real& sig12,
real& ssig1, real& csig1, real& ssig2, real& csig2, real& ssig1, real& csig1, real& ssig2, real& csig2,
real& k1, bool diffp, real& dlam12, real& eps, bool diffp, real& dlam12,
real tc[], real zc[], real ec[]) real C1a[], real C2a[], real C3a[])
const throw(); const throw();
static const real eps2, tol0, tol1, tol2, xthresh; static const real eps2, tol0, tol1, tol2, xthresh;
const real _a, _r, _f, _f1, _e2, _ep2, _n, _b; const real _a, _r, _f, _f1, _e2, _ep2, _n, _b, _etol2;
static real SinSeries(real sinx, real cosx, const real c[], int n) static real SinSeries(real sinx, real cosx, const real c[], int n)
throw(); throw();
static inline real AngNormalize(real x) throw() { static inline real AngNormalize(real x) throw() {
// Place angle in [-180, 180). Assumes x is in [-540, 540). // Place angle in [-180, 180). Assumes x is in [-540, 540).
return x >= 180 ? x - 360 : x < -180 ? x + 360 : x; return x >= 180 ? x - 360 : x < -180 ? x + 360 : x;
} }
static inline real AngRound(real x) throw() { static inline real AngRound(real x) throw() {
// The makes the smallest gap in x = 1/16 - nextafter(1/16, 0) = 1/2^ 57 // The makes the smallest gap in x = 1/16 - nextafter(1/16, 0) = 1/2^ 57
// for reals = 0.7 pm on the earth if x is an angle in degrees. (Thi s // for reals = 0.7 pm on the earth if x is an angle in degrees. (Thi s
// is about 1000 times more resolution than we get with angles around 90 // is about 1000 times more resolution than we get with angles around 90
skipping to change at line 142 skipping to change at line 124
volatile real y = std::abs(x); volatile real y = std::abs(x);
// The compiler mustn't "simplify" z - (z - y) to y // The compiler mustn't "simplify" z - (z - y) to y
y = y < z ? z - (z - y) : y; y = y < z ? z - (z - y) : y;
return x < 0 ? -y : y; return x < 0 ? -y : y;
} }
static inline void SinCosNorm(real& sinx, real& cosx) throw() { static inline void SinCosNorm(real& sinx, real& cosx) throw() {
real r = Math::hypot(sinx, cosx); real r = Math::hypot(sinx, cosx);
sinx /= r; sinx /= r;
cosx /= r; cosx /= r;
} }
// These are Maxima generated functions to provide series approximation s to // These are Maxima generated functions to provide series approximation s to
// the integrals for the ellipsoidal geodesic. // the integrals for the ellipsoidal geodesic.
static real tauFactorm1(real k1) throw(); static real A1m1f(real eps) throw();
static void tauCoeff(real k1, real t[]) throw(); static void C1f(real eps, real c[]) throw();
static void sigCoeff(real k1, real tp[]) throw(); static void C1pf(real eps, real c[]) throw();
static real zetFactorm1(real k1) throw(); static real A2m1f(real eps) throw();
static void zetCoeff(real k1, real t[]) throw(); static void C2f(real eps, real c[]) throw();
static real etaFactor(real f, real k1) throw(); static real A3f(real f, real eps) throw();
static void etaCoeff(real f, real k1, real h[]) throw(); static void C3f(real f, real eps, real c[]) throw();
public: public:
/** /**
* Constructor for a ellipsoid with equatorial radius \e a (meters) and * Constructor for a ellipsoid with equatorial radius \e a (meters) and
* reciprocal flattening \e r. Setting \e r = 0 implies \e r = inf or * reciprocal flattening \e r. Setting \e r = 0 implies \e r = inf or
* flattening = 0 (i.e., a sphere). Negative \e r indicates a prolate * flattening = 0 (i.e., a sphere). Negative \e r indicates a prolate
* ellipsoid. An exception is thrown if \e a is not positive. * ellipsoid. An exception is thrown if \e a is not positive.
********************************************************************** / ********************************************************************** /
Geodesic(real a, real r); Geodesic(real a, real r);
skipping to change at line 226 skipping to change at line 207
* used in the constructor. A value of 0 is returned for a sphere * used in the constructor. A value of 0 is returned for a sphere
* (infinite inverse flattening). * (infinite inverse flattening).
********************************************************************** / ********************************************************************** /
Math::real InverseFlattening() const throw() { return _r; } Math::real InverseFlattening() const throw() { return _r; }
/** /**
* A global instantiation of Geodesic with the parameters for the WGS84 * A global instantiation of Geodesic with the parameters for the WGS84
* ellipsoid. * ellipsoid.
********************************************************************** / ********************************************************************** /
const static Geodesic WGS84; const static Geodesic WGS84;
}; };
/** /**
* \brief A geodesic line. * \brief A geodesic line.
* *
* GeodesicLine facilitates the determination of a series of points on a * GeodesicLine facilitates the determination of a series of points on a
* single geodesic. Geodesic.Line returns a GeodesicLine object with the * single geodesic. Geodesic.Line returns a GeodesicLine object with the
* geodesic defined by by \e lat1, \e lon1, and \e azi1. * geodesic defined by by \e lat1, \e lon1, and \e azi1.
* GeodesicLine.Position returns the \e lat2, \e lon2, \e azi2, and \e m1 2 * GeodesicLine.Position returns the \e lat2, \e lon2, \e azi2, and \e m1 2
* given \e s12. An example of use of this class is: * given \e s12. An example of use of this class is:
skipping to change at line 268 skipping to change at line 250
* Similarly, a vector can be used to hold GeodesicLine objects. * Similarly, a vector can be used to hold GeodesicLine objects.
* *
* The calculations are accurate to better than 12 nm. (See \ref geoderr ors * The calculations are accurate to better than 12 nm. (See \ref geoderr ors
* for details.) * for details.)
**********************************************************************/ **********************************************************************/
class GeodesicLine { class GeodesicLine {
private: private:
typedef Math::real real; typedef Math::real real;
friend class Geodesic; friend class Geodesic;
friend class CassiniSoldner; static const int nC1 = Geodesic::nC1, nC1p = Geodesic::nC1p,
static const int ntau = Geodesic::ntau; nC2 = Geodesic::nC2, nC3 = Geodesic::nC3;
static const int nsig = Geodesic::nsig;
static const int nzet = Geodesic::nzet;
static const int neta = Geodesic::neta;
real _lat1, _lon1, _azi1; real _lat1, _lon1, _azi1;
real _a, _r, _b, _f1, _salp0, _calp0, _u2, real _a, _r, _b, _f1, _salp0, _calp0, _k2,
_ssig1, _csig1, _stau1, _ctau1, _somg1, _comg1, _ssig1, _csig1, _stau1, _ctau1, _somg1, _comg1,
_taufm1, _zetfm1, _etaf, _dtau1, _dzet1, _dlam1; _A1m1, _A2m1, _A3c, _B11, _B21, _B31;
real _tauCoeff[ntau ? ntau : 1], _sigCoeff[nsig ? nsig : 1], // index zero elements of these arrays are unused
_zetCoeff[nzet ? nzet : 1], _etaCoeff[neta ? neta : 1]; real _C1a[nC1 + 1], _C1pa[nC1p + 1], _C2a[nC2 + 1], _C3a[nC3];
static inline real sq(real x) throw() { return x * x; } static inline real sq(real x) throw() { return x * x; }
GeodesicLine(const Geodesic& g, real lat1, real lon1, real azi1) GeodesicLine(const Geodesic& g, real lat1, real lon1, real azi1)
throw(); throw();
public: public:
/** /**
* A default constructor. If GeodesicLine::Position is called on the * A default constructor. If GeodesicLine::Position is called on the
* resulting object, it returns immediately (without doing any * resulting object, it returns immediately (without doing any
* calculations). The object should be set with a call to Geodesic::Li ne. * calculations). The object should be set with a call to Geodesic::Li ne.
skipping to change at line 308 skipping to change at line 287
* \e s12 can be signed. If \e arcmode (default false) is set to true, \e * \e s12 can be signed. If \e arcmode (default false) is set to true, \e
* s12 is interpreted as the arc length \e a12 (in degrees) on the * s12 is interpreted as the arc length \e a12 (in degrees) on the
* auxiliary sphere. Returned value is the arc length \e a12 (degrees) if * auxiliary sphere. Returned value is the arc length \e a12 (degrees) if
* \e arcmode is false, otherwise it is the distance \e s12 (meters). * \e arcmode is false, otherwise it is the distance \e s12 (meters).
********************************************************************** / ********************************************************************** /
Math::real Position(real s12, real& lat2, real& lon2, Math::real Position(real s12, real& lat2, real& lon2,
real& azi2, real &m12, bool arcmode = false) real& azi2, real &m12, bool arcmode = false)
const throw(); const throw();
/** /**
* Return the scale of the geodesic line extending an arc length \e a12
* (degrees) from point 1 to point 2. \e M12 (a number) measures the
* convergence of initially parallel geodesics. It is defined by the
* following construction: starting at point 1 proceed at azimuth \e az
i1 +
* 90<sup>o</sup> a small distance \e dt; turn -90<sup>o</sup> and proc
eed
* a distance \e s12 (\e not the arc length \e a12); the distance to po
int
* 2 is given by \e M12 \e dt. \e M21 is defined analogously.
**********************************************************************
/
void Scale(real a12, real& M12, real& M21) const throw();
/**
* Has this object been initialized so that Position can be called? * Has this object been initialized so that Position can be called?
********************************************************************** / ********************************************************************** /
bool Init() const throw() { return _b > 0; } bool Init() const throw() { return _b > 0; }
/** /**
* Return the latitude of point 1 (degrees). * Return the latitude of point 1 (degrees).
********************************************************************** / ********************************************************************** /
Math::real Latitude() const throw() { return Init() ? _lat1 : 0; } Math::real Latitude() const throw() { return Init() ? _lat1 : 0; }
/** /**
skipping to change at line 329 skipping to change at line 319
********************************************************************** / ********************************************************************** /
Math::real Longitude() const throw() { return Init() ? _lon1 : 0; } Math::real Longitude() const throw() { return Init() ? _lon1 : 0; }
/** /**
* Return the azimuth (degrees) of the geodesic line as it passes throu gh * Return the azimuth (degrees) of the geodesic line as it passes throu gh
* point 1. * point 1.
********************************************************************** / ********************************************************************** /
Math::real Azimuth() const throw() { return Init() ? _azi1 : 0; } Math::real Azimuth() const throw() { return Init() ? _azi1 : 0; }
/** /**
* Return the azimuth (degrees) of the geodesic line as it crosses the
* equator in a northward direction.
**********************************************************************
/
Math::real EquatorialAzimuth() const throw() {
return Init() ? atan2(_salp0, _calp0) / Constants::degree() : 0;
}
/**
* Return the arc length (degrees) between the northward equatorial
* crossing and point 1.
**********************************************************************
/
Math::real EquatorialArc() const throw() {
return Init() ? atan2(_ssig1, _csig1) / Constants::degree() : 0;
}
/**
* The major radius of the ellipsoid (meters). This is that value of \ e a * The major radius of the ellipsoid (meters). This is that value of \ e a
* inherited from the Geodesic object used in the constructor. * inherited from the Geodesic object used in the constructor.
********************************************************************** / ********************************************************************** /
Math::real MajorRadius() const throw() { return Init() ? _a : 0; } Math::real MajorRadius() const throw() { return Init() ? _a : 0; }
/** /**
* The inverse flattening of the ellipsoid. This is that value of \e r * The inverse flattening of the ellipsoid. This is that value of \e r
* inherited from the Geodesic object used in the constructor. A value of * inherited from the Geodesic object used in the constructor. A value of
* 0 is returned for a sphere (infinite inverse flattening). * 0 is returned for a sphere (infinite inverse flattening).
********************************************************************** / ********************************************************************** /
 End of changes. 20 change blocks. 
54 lines changed or deleted 64 lines changed or added


 TransverseMercator.hpp   TransverseMercator.hpp 
/** /**
* \file TransverseMercator.hpp * \file TransverseMercator.hpp
* \brief Header for GeographicLib::TransverseMercator class * \brief Header for GeographicLib::TransverseMercator class
* *
* Copyright (c) Charles Karney (2008, 2009, 2010) <charles@karney.com> * Copyright (c) Charles Karney (2008, 2009, 2010) <charles@karney.com>
* and licensed under the LGPL. For more information, see * and licensed under the LGPL. For more information, see
* http://geographiclib.sourceforge.net/ * http://geographiclib.sourceforge.net/
**********************************************************************/ **********************************************************************/
#if !defined(GEOGRAPHICLIB_TRANSVERSEMERCATOR_HPP) #if !defined(GEOGRAPHICLIB_TRANSVERSEMERCATOR_HPP)
#define GEOGRAPHICLIB_TRANSVERSEMERCATOR_HPP "$Id: TransverseMercator.hpp 6 807 2010-02-01 11:26:34Z karney $" #define GEOGRAPHICLIB_TRANSVERSEMERCATOR_HPP "$Id: TransverseMercator.hpp 6 824 2010-04-19 14:25:10Z karney $"
#include "GeographicLib/Constants.hpp" #include "GeographicLib/Constants.hpp"
#if !defined(TM_TX_MAXPOW) #if !defined(TM_TX_MAXPOW)
/** /**
* The order of the series approximation used in * The order of the series approximation used in
* GeographicLib::TransverseMercator. TM_TX_MAXPOW can be set to any integ er * GeographicLib::TransverseMercator. TM_TX_MAXPOW can be set to any integ er
* in [4, 8]. * in [4, 8].
**********************************************************************/ **********************************************************************/
#define TM_TX_MAXPOW \ #define TM_TX_MAXPOW \
(sizeof(real) == sizeof(double) ? 6 : sizeof(real) == sizeof(float) ? 4 : 8 ) (GEOGRAPHICLIB_PREC == 1 ? 6 : GEOGRAPHICLIB_PREC == 0 ? 4 : 8)
#endif #endif
namespace GeographicLib { namespace GeographicLib {
/** /**
* \brief Transverse Mercator Projection * \brief Transverse Mercator Projection
* *
* This uses Kr&uuml;ger's method which evaluates the projection and its * This uses Kr&uuml;ger's method which evaluates the projection and its
* inverse in terms of a series. See * inverse in terms of a series. See
* - L. Kr&uuml;ger, * - L. Kr&uuml;ger,
* <a href="http://dx.doi.org/10.2312/GFZ.b103-krueger28"> Konforme * <a href="http://dx.doi.org/10.2312/GFZ.b103-krueger28"> Konforme
* Abbildung des Erdellipsoids in der Ebene</a> (Conformal mapping of the * Abbildung des Erdellipsoids in der Ebene</a> (Conformal mapping of the
* ellipsoidal earth to the plane), Royal Prussian Geodetic Institute, New * ellipsoidal earth to the plane), Royal Prussian Geodetic Institute, New
* Series 52, 172 pp. (1912). * Series 52, 172 pp. (1912).
* *
* Kr&uuml;ger's method has been extended from 4th to 6th order. The max imum * Kr&uuml;ger's method has been extended from 4th to 6th order. The max imum
* errors is 5 nm (ground distance) for all positions within 35 degrees o f * errors is 5 nm (ground distance) for all positions within 35 degrees o f
* the central meridian. The error in the convergence is 2e-15" and the * the central meridian. The error in the convergence is 2e-15&quot; and the
* relative error in the scale is 6e-12%%. (See \ref tmerrors for the we asel * relative error in the scale is 6e-12%%. (See \ref tmerrors for the we asel
* words.) The speed penalty in going to 6th order is only about 1%. * words.) The speed penalty in going to 6th order is only about 1%.
* GeographicLib::TransverseMercatorExact is an alternative implementatio n of * GeographicLib::TransverseMercatorExact is an alternative implementatio n of
* the projection using exact formulas which yield accurate (to 8 nm) * the projection using exact formulas which yield accurate (to 8 nm) res
* results over the entire ellipsoid. ults
* over the entire ellipsoid.
* *
* The ellipsoid parameters and the central scale are set in the construc tor. * The ellipsoid parameters and the central scale are set in the construc tor.
* The central meridian (which is a trivial shift of the longitude) is * The central meridian (which is a trivial shift of the longitude) is
* specified as the \e lon0 argument of the TransverseMercator::Forward a nd * specified as the \e lon0 argument of the TransverseMercator::Forward a nd
* TransverseMercator::Reverse functions. The latitude of origin is take n to * TransverseMercator::Reverse functions. The latitude of origin is take n to
* be the equator. There is no provision in this class for specifying a * be the equator. There is no provision in this class for specifying a
* false easting or false northing or a different latitude of origin. * false easting or false northing or a different latitude of origin.
* However these are can be simply included by the calling funtcion. For * However these are can be simply included by the calling funtcion. For
* example, the UTMUPS class applies the false easting and false northing for * example, the UTMUPS class applies the false easting and false northing for
* the UTM projections. A more complicated example is the British Nation al * the UTM projections. A more complicated example is the British Nation al
 End of changes. 4 change blocks. 
5 lines changed or deleted 6 lines changed or added


 TransverseMercatorExact.hpp   TransverseMercatorExact.hpp 
/** /**
* \file TransverseMercatorExact.hpp * \file TransverseMercatorExact.hpp
* \brief Header for GeographicLib::TransverseMercatorExact class * \brief Header for GeographicLib::TransverseMercatorExact class
* *
* Copyright (c) Charles Karney (2008, 2009, 2010) <charles@karney.com> * Copyright (c) Charles Karney (2008, 2009, 2010) <charles@karney.com>
* and licensed under the LGPL. For more information, see * and licensed under the LGPL. For more information, see
* http://geographiclib.sourceforge.net/ * http://geographiclib.sourceforge.net/
**********************************************************************/ **********************************************************************/
#if !defined(GEOGRAPHICLIB_TRANSVERSEMERCATOREXACT_HPP) #if !defined(GEOGRAPHICLIB_TRANSVERSEMERCATOREXACT_HPP)
#define GEOGRAPHICLIB_TRANSVERSEMERCATOREXACT_HPP "$Id: TransverseMercatorE xact.hpp 6807 2010-02-01 11:26:34Z karney $" #define GEOGRAPHICLIB_TRANSVERSEMERCATOREXACT_HPP "$Id: TransverseMercatorE xact.hpp 6824 2010-04-19 14:25:10Z karney $"
#include "GeographicLib/Constants.hpp" #include "GeographicLib/Constants.hpp"
#include "GeographicLib/EllipticFunction.hpp" #include "GeographicLib/EllipticFunction.hpp"
namespace GeographicLib { namespace GeographicLib {
/** /**
* \brief An exact implementation of the Transverse Mercator Projection * \brief An exact implementation of the Transverse Mercator Projection
* *
* Implementation of the Transverse Mercator Projection given in * Implementation of the Transverse Mercator Projection given in
skipping to change at line 33 skipping to change at line 33
* <a href="http://dx.doi.org/10.3138/X687-1574-4325-WM62"> Conformal * <a href="http://dx.doi.org/10.3138/X687-1574-4325-WM62"> Conformal
* Projections Based On Jacobian Elliptic Functions</a>, Part V of * Projections Based On Jacobian Elliptic Functions</a>, Part V of
* Conformal Projections Based on Elliptic Functions, * Conformal Projections Based on Elliptic Functions,
* (B. V. Gutsell, Toronto, 1976), 128pp., * (B. V. Gutsell, Toronto, 1976), 128pp.,
* ISBN: 0919870163 * ISBN: 0919870163
* (also appeared as: * (also appeared as:
* Monograph 16, Suppl. No. 1 to Canadian Cartographer, Vol 13). * Monograph 16, Suppl. No. 1 to Canadian Cartographer, Vol 13).
* *
* This method gives the correct results for forward and reverse * This method gives the correct results for forward and reverse
* transformations subject to the branch cut rules (see the description o f * transformations subject to the branch cut rules (see the description o f
* the \e extendp argument to the constructor).. The maximum error is ab * the \e extendp argument to the constructor). The maximum error is abo
out ut 8
* 8 nm (ground distance) for the forward and reverse transformations. T * nm (ground distance) for the forward and reverse transformations. The
he * error in the convergence is 2e-15&quot;, the relative error in the sca
* error in the convergence is 2e-15", the relative error in the scale is le
* 7e-12%%. (See \ref tmerrors for the weasel words.) The method is "ex * is 7e-12%%. (See \ref tmerrors for the weasel words.) The method is
act" * "exact" in the sense that the errors are close to the round-off limit
* in the sense that the errors are close to the round-off limit and that and
no * that no changes are needed in the algorithms for them to be used with
* changes are needed in the algorithms for them to be used with reals of * reals of a higher precision. Thus the errors using long double (with
a a
* higher precision. Thus the errors using long double (with a 64-bit * 64-bit fraction) are about 2000 times smaller than using double (with
* fraction) are about 2000 times smaller than using double (with a 53-bi a
t * 53-bit fraction).
* fraction).
* *
* This algorithm is about 4.5 times slower than the 6th-order Kr&uuml;ge r * This algorithm is about 4.5 times slower than the 6th-order Kr&uuml;ge r
* method, GeographicLib::TransverseMercator, taking about 11 us for a * method, GeographicLib::TransverseMercator, taking about 11 us for a
* combined forward and reverse projection on a 2.6 GHz Intel machine (g+ +, * combined forward and reverse projection on a 2.6 GHz Intel machine (g+ +,
* version 4.3.0, -O3). * version 4.3.0, -O3).
* *
* The ellipsoid parameters and the central scale are set in the construc tor. * The ellipsoid parameters and the central scale are set in the construc tor.
* The central meridian (which is a trivial shift of the longitude) is * The central meridian (which is a trivial shift of the longitude) is
* specified as the \e lon0 argument of the TransverseMercatorExact::Forw ard * specified as the \e lon0 argument of the TransverseMercatorExact::Forw ard
* and TransverseMercatorExact::Reverse functions. The latitude of origi n is * and TransverseMercatorExact::Reverse functions. The latitude of origi n is
 End of changes. 2 change blocks. 
16 lines changed or deleted 15 lines changed or added

This html diff was produced by rfcdiff 1.41. The latest version is available from http://tools.ietf.org/tools/rfcdiff/