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 6891 2010-11-16 16: 10:53Z karney $" #define GEOGRAPHICLIB_CONSTANTS_HPP "$Id: Constants.hpp 6918 2010-12-21 12: 56:07Z karney $"
/** /**
* A simple compile-time assert. This is designed to be compatible with th * Are C++0X math functions available?
e **********************************************************************/
* C++0X static_assert. #if !defined(GEOGRAPHICLIB_CPLUSPLUS0X_MATH)
#if defined(__GXX_EXPERIMENTAL_CXX0X__)
#define GEOGRAPHICLIB_CPLUSPLUS0X_MATH 1
#else
#define GEOGRAPHICLIB_CPLUSPLUS0X_MATH 0
#endif
#endif
/**
* A compile-time assert. Use C++0X static_assert, if available.
**********************************************************************/ **********************************************************************/
#if !defined(STATIC_ASSERT) #if !defined(STATIC_ASSERT)
#if defined(__GXX_EXPERIMENTAL_CXX0X__)
#define STATIC_ASSERT static_assert
#elif defined(_MSC_VER) && _MSC_VER >= 1600
#define STATIC_ASSERT static_assert
#else
#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
#endif
#if defined(__GNUC__) #if defined(__GNUC__)
// Suppress "defined but not used" warnings // Suppress "defined but not used" warnings
#define RCSID_DECL(x) namespace \ #define RCSID_DECL(x) namespace \
{ char VAR_ ## x [] __attribute__((used)) = x; } { char VAR_ ## x [] __attribute__((used)) = x; }
#else #else
/** /**
* Insertion of RCS Id strings into the object file. * Insertion of RCS Id strings into the object file.
**********************************************************************/ **********************************************************************/
#define RCSID_DECL(x) namespace { char VAR_ ## x [] = x; } #define RCSID_DECL(x) namespace { char VAR_ ## x [] = x; }
skipping to change at line 133 skipping to change at line 149
* @return sqrt(\e x<sup>2</sup> + \e y<sup>2</sup>). * @return sqrt(\e x<sup>2</sup> + \e y<sup>2</sup>).
********************************************************************** / ********************************************************************** /
static inline real hypot(real x, real y) throw() { static inline real hypot(real x, real y) throw() {
x = std::abs(x); x = std::abs(x);
y = std::abs(y); y = std::abs(y);
real real
a = std::max(x, y), a = std::max(x, y),
b = std::min(x, y) / a; b = std::min(x, y) / a;
return a * std::sqrt(1 + b * b); return a * std::sqrt(1 + b * b);
} }
#elif GEOGRAPHICLIB_CPLUSPLUS0X_MATH
static inline real hypot(real x, real y) throw() {return std::hypot(x,
y);}
#elif defined(_MSC_VER) #elif defined(_MSC_VER)
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) #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 _hypot(x, y); } { return _hypot(x, y); }
#endif #endif
#else #else
skipping to change at line 154 skipping to change at line 172
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) #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 #endif
#if defined(DOXYGEN) || defined(_MSC_VER) #if defined(DOXYGEN) || (defined(_MSC_VER) && !GEOGRAPHICLIB_CPLUSPLUS0X_MA TH)
/** /**
* exp(\e x) - 1 accurate near \e x = 0. This is taken from * 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.
* *
* @param[in] x * @param[in] x
* @return exp(\e x) - 1. * @return exp(\e x) - 1.
********************************************************************** / ********************************************************************** /
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);
} }
#elif GEOGRAPHICLIB_CPLUSPLUS0X_MATH
static inline real expm1(real x) throw() { return std::expm1(x); }
#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) #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 #endif
#if defined(DOXYGEN) || defined(_MSC_VER) #if defined(DOXYGEN) || (defined(_MSC_VER) && !GEOGRAPHICLIB_CPLUSPLUS0X_MA TH)
/** /**
* log(\e x + 1) accurate near \e x = 0. This is taken See * 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.
* *
* @param[in] x * @param[in] x
* @return log(\e x + 1). * @return log(\e x + 1).
skipping to change at line 204 skipping to change at line 224
static inline real log1p(real x) throw() { static inline real log1p(real x) throw() {
volatile real volatile real
y = 1 + x, y = 1 + x,
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;
} }
#elif GEOGRAPHICLIB_CPLUSPLUS0X_MATH
static inline real log1p(real x) throw() { return std::log1p(x); }
#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) #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 #endif
#if defined(DOXYGEN) || defined(_MSC_VER) #if defined(DOXYGEN) || (defined(_MSC_VER) && !GEOGRAPHICLIB_CPLUSPLUS0X_MA TH)
/** /**
* The inverse hyperbolic sine function. This is defined in terms of * The inverse hyperbolic sine function. This is defined in terms of
* Math::log1p(\e x) in order to maintain accuracy near \e x = 0. In * Math::log1p(\e x) in order to maintain accuracy near \e x = 0. In
* addition, the odd parity of the function is enforced. * addition, the odd parity of the function is enforced.
* *
* @param[in] x * @param[in] x
* @return asinh(\e x). * @return asinh(\e x).
********************************************************************** / ********************************************************************** /
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;
} }
#elif GEOGRAPHICLIB_CPLUSPLUS0X_MATH
static inline real asinh(real x) throw() { return std::asinh(x); }
#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) #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 #endif
#if defined(DOXYGEN) || defined(_MSC_VER) #if defined(DOXYGEN) || (defined(_MSC_VER) && !GEOGRAPHICLIB_CPLUSPLUS0X_MA TH)
/** /**
* The inverse hyperbolic tangent function. This is defined in terms o f * The inverse hyperbolic tangent function. This is defined in terms o f
* Math::log1p(\e x) in order to maintain accuracy near \e x = 0. In * Math::log1p(\e x) in order to maintain accuracy near \e x = 0. In
* addition, the odd parity of the function is enforced. * addition, the odd parity of the function is enforced.
* *
* @param[in] x * @param[in] x
* @return atanh(\e x). * @return atanh(\e x).
********************************************************************** / ********************************************************************** /
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;
} }
#elif GEOGRAPHICLIB_CPLUSPLUS0X_MATH
static inline real atanh(real x) throw() { return std::atanh(x); }
#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) #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 #endif
#if defined(DOXYGEN) || defined(_MSC_VER) #if defined(DOXYGEN) || (defined(_MSC_VER) && !GEOGRAPHICLIB_CPLUSPLUS0X_MA TH)
/** /**
* The cube root function. * The cube root function.
* *
* @param[in] x * @param[in] x
* @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;
} }
#elif GEOGRAPHICLIB_CPLUSPLUS0X_MATH
static inline real cbrt(real x) throw() { return std::cbrt(x); }
#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) #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 #endif
/** /**
* Test for finiteness. * Test for finiteness.
* *
* @param[in] x * @param[in] x
* @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() {
#if defined(DOXYGEN) #if defined(DOXYGEN)
return std::abs(x) <= std::numeric_limits<real>::max(); return std::abs(x) <= std::numeric_limits<real>::max();
#elif defined(_MSC_VER) #elif (defined(_MSC_VER) && !GEOGRAPHICLIB_CPLUSPLUS0X_MATH)
return _finite(x) != 0; return _finite(x) != 0;
#else #else
return std::isfinite(x) != 0; return std::isfinite(x);
#endif #endif
} }
/** /**
* The NaN (not a number) * The NaN (not a number)
* *
* @return NaN if available, otherwise return the max real. * @return NaN if available, otherwise return the max real.
********************************************************************** / ********************************************************************** /
static inline real NaN() throw() { static inline real NaN() throw() {
return std::numeric_limits<real>::has_quiet_NaN ? return std::numeric_limits<real>::has_quiet_NaN ?
std::numeric_limits<real>::quiet_NaN() : std::numeric_limits<real>::quiet_NaN() :
std::numeric_limits<real>::max(); std::numeric_limits<real>::max();
} }
/** /**
* Test for NaN.
*
* @param[in] x
* @return true if argument is a NaN.
**********************************************************************
/
static inline bool isnan(real x) throw() {
#if defined(DOXYGEN) || (defined(_MSC_VER) && !GEOGRAPHICLIB_CPLUSPLUS0X_MA
TH)
return x != x;
#else
return std::isnan(x);
#endif
}
/**
* Infinity
*
* @return infinity if available, otherwise return the max real.
**********************************************************************
/
static inline real infinity() throw() {
return std::numeric_limits<real>::has_infinity ?
std::numeric_limits<real>::infinity() :
std::numeric_limits<real>::max();
}
/**
* @return \e pi in extended precision * @return \e pi in extended precision
********************************************************************** / ********************************************************************** /
static inline extended epi() throw() static inline extended epi() throw()
// good for about 168-bit accuracy // good for about 168-bit accuracy
{ return extended(3.1415926535897932384626433832795028841971693993751L) ; } { return extended(3.1415926535897932384626433832795028841971693993751L) ; }
/** /**
* @return the number of radians in a degree in extended precision. * @return the number of radians in a degree in extended precision.
********************************************************************** / ********************************************************************** /
static inline extended edegree() throw() { return epi() / 180; } static inline extended edegree() throw() { return epi() / 180; }
 End of changes. 18 change blocks. 
11 lines changed or deleted 67 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 6866 2010-09-11 02:15:29Z karne y $" #define GEOGRAPHICLIB_DMS_HPP "$Id: DMS.hpp 6911 2010-12-09 23:13:55Z karne y $"
#include "GeographicLib/Constants.hpp" #include "GeographicLib/Constants.hpp"
#include <sstream> #include <sstream>
#include <iomanip> #include <iomanip>
namespace GeographicLib { namespace GeographicLib {
/** /**
* \brief Convert between degrees and %DMS representation. * \brief Convert between degrees and %DMS representation.
* *
skipping to change at line 41 skipping to change at line 41
return r == std::string::npos ? -1 : int(r); return r == std::string::npos ? -1 : int(r);
} }
template<typename T> static std::string str(T x) { template<typename T> static std::string str(T x) {
std::ostringstream s; s << x; return s.str(); std::ostringstream s; s << x; return s.str();
} }
static const std::string hemispheres; static const std::string hemispheres;
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];
static Math::real NumMatch(const std::string& s);
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. * and longitudes.
********************************************************************** / ********************************************************************** /
enum flag { enum flag {
/** /**
skipping to change at line 234 skipping to change at line 235
* @return formatting string * @return formatting string
* *
* \e prec indicates the precision relative to 1 degree, e.g., \e prec = 3 * \e prec indicates the precision relative to 1 degree, e.g., \e prec = 3
* gives a result accurate to 0.1' and \e prec = 4 gives a result accur ate * gives a result accurate to 0.1' and \e prec = 4 gives a result accur ate
* to 1&quot;. \e ind is interpreted as in DMS::Encode with the additi onal * to 1&quot;. \e ind is interpreted as in DMS::Encode with the additi onal
* facility at DMS::NUMBER treats \e angle a number in fixed format wit h * facility at DMS::NUMBER treats \e angle a number in fixed format wit h
* precision \e prec. * precision \e prec.
********************************************************************** / ********************************************************************** /
static std::string Encode(real angle, unsigned prec, flag ind = NONE) { static std::string Encode(real angle, unsigned prec, flag ind = NONE) {
if (ind == NUMBER) { if (ind == NUMBER) {
if (!Math::isfinite(angle))
return angle < 0 ? std::string("-inf") :
(angle > 0 ? std::string("inf") : std::string("nan"));
std::ostringstream s; std::ostringstream s;
s << std::fixed << std::setprecision(prec) << angle; s << std::fixed << std::setprecision(prec) << angle;
return s.str(); return s.str();
} else } else
return Encode(angle, return Encode(angle,
prec < 2 ? DEGREE : (prec < 4 ? MINUTE : SECOND), prec < 2 ? DEGREE : (prec < 4 ? MINUTE : SECOND),
prec < 2 ? prec : (prec < 4 ? prec - 2 : prec - 4), prec < 2 ? prec : (prec < 4 ? prec - 2 : prec - 4),
ind); ind);
} }
 End of changes. 3 change blocks. 
1 lines changed or deleted 5 lines changed or added


 Geoid.hpp   Geoid.hpp 
/** /**
* \file Geoid.hpp * \file Geoid.hpp
* \brief Header for GeographicLib::Geoid class * \brief Header for GeographicLib::Geoid 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_GEOID_HPP) #if !defined(GEOGRAPHICLIB_GEOID_HPP)
#define GEOGRAPHICLIB_GEOID_HPP "$Id: Geoid.hpp 6888 2010-11-13 15:19:42Z k arney $" #define GEOGRAPHICLIB_GEOID_HPP "$Id: Geoid.hpp 6910 2010-12-08 22:01:30Z k arney $"
#include "GeographicLib/Constants.hpp" #include "GeographicLib/Constants.hpp"
#include <string> #include <string>
#include <vector> #include <vector>
#include <fstream> #include <fstream>
namespace GeographicLib { namespace GeographicLib {
/** /**
* \brief Computing the height of the geoid * \brief Looking up the height of the geoid
* *
* This class evaluated the height of one of the standard geoids, EGM84, * This class evaluated the height of one of the standard geoids, EGM84,
* EGM96, or EGM2008 by bilinear or cubic interpolation into a rectangula r * EGM96, or EGM2008 by bilinear or cubic interpolation into a rectangula r
* grid of data. These geoid models are documented in * grid of data. These geoid models are documented in
* - EGM84: * - EGM84:
* http://earth-info.nga.mil/GandG/wgs84/gravitymod/wgs84_180/wgs84_180 .html * http://earth-info.nga.mil/GandG/wgs84/gravitymod/wgs84_180/wgs84_180 .html
* - EGM96: * - EGM96:
* http://earth-info.nga.mil/GandG/wgs84/gravitymod/egm96/egm96.html * http://earth-info.nga.mil/GandG/wgs84/gravitymod/egm96/egm96.html
* - EGM2008: * - EGM2008:
* http://earth-info.nga.mil/GandG/wgs84/gravitymod/egm2008 * http://earth-info.nga.mil/GandG/wgs84/gravitymod/egm2008
 End of changes. 2 change blocks. 
2 lines changed or deleted 2 lines changed or added


 Gnomonic.hpp   Gnomonic.hpp 
/** /**
* \file Gnomonic.hpp * \file Gnomonic.hpp
* \brief Header for GeographicLib::Gnomonic class * \brief Header for GeographicLib::Gnomonic class
* *
* Copyright (c) Charles Karney (2009, 2010) <charles@karney.com> * Copyright (c) Charles Karney (2010) <charles@karney.com> and licensed un
* and licensed under the LGPL. For more information, see der
* http://geographiclib.sourceforge.net/ * the LGPL. For more information, see http://geographiclib.sourceforge.ne
t/
**********************************************************************/ **********************************************************************/
#if !defined(GEOGRAPHICLIB_GNOMONIC_HPP) #if !defined(GEOGRAPHICLIB_GNOMONIC_HPP)
#define GEOGRAPHICLIB_GNOMONIC_HPP "$Id: Gnomonic.hpp 6867 2010-09-11 13:04 :26Z karney $" #define GEOGRAPHICLIB_GNOMONIC_HPP "$Id: Gnomonic.hpp 6911 2010-12-09 23:13 :55Z karney $"
#include "GeographicLib/Geodesic.hpp" #include "GeographicLib/Geodesic.hpp"
#include "GeographicLib/GeodesicLine.hpp" #include "GeographicLib/GeodesicLine.hpp"
#include "GeographicLib/Constants.hpp" #include "GeographicLib/Constants.hpp"
namespace GeographicLib { namespace GeographicLib {
/** /**
* \brief %Gnomonic Projection. * \brief %Gnomonic Projection.
* *
 End of changes. 2 change blocks. 
4 lines changed or deleted 5 lines changed or added


 LambertConformalConic.hpp   LambertConformalConic.hpp 
/** /**
* \file LambertConformalConic.hpp * \file LambertConformalConic.hpp
* \brief Header for GeographicLib::LambertConformalConic class * \brief Header for GeographicLib::LambertConformalConic class
* *
* Copyright (c) Charles Karney (2009, 2010) <charles@karney.com> * Copyright (c) Charles Karney (2010) <charles@karney.com> and licensed un
* and licensed under the LGPL. For more information, see der
* http://geographiclib.sourceforge.net/ * the LGPL. For more information, see http://geographiclib.sourceforge.ne
t/
**********************************************************************/ **********************************************************************/
#if !defined(GEOGRAPHICLIB_LAMBERTCONFORMALCONIC_HPP) #if !defined(GEOGRAPHICLIB_LAMBERTCONFORMALCONIC_HPP)
#define GEOGRAPHICLIB_LAMBERTCONFORMALCONIC_HPP "$Id: LambertConformalConic .hpp 6898 2010-11-19 19:21:49Z karney $" #define GEOGRAPHICLIB_LAMBERTCONFORMALCONIC_HPP "$Id: LambertConformalConic .hpp 6916 2010-12-20 23:03:47Z karney $"
#include "GeographicLib/Constants.hpp" #include "GeographicLib/Constants.hpp"
#include <algorithm> #include <algorithm>
namespace GeographicLib { namespace GeographicLib {
/** /**
* \brief Lambert Conformal Conic Projection * \brief Lambert Conformal Conic Projection
* *
* Implementation taken from the report, * Implementation taken from the report,
skipping to change at line 38 skipping to change at line 37
* differences have been used to transform the expressions into ones whic h * differences have been used to transform the expressions into ones whic h
* may be evaluated accurately and that Newton's method is used to invert the * may be evaluated accurately and that Newton's method is used to invert the
* projection. In this implementation, the projection correctly becomes the * projection. In this implementation, the projection correctly becomes the
* Mercator projection or the polar sterographic projection when the stan dard * Mercator projection or the polar sterographic projection when the stan dard
* latitude is the equator or a pole. The accuracy of the projections is * latitude is the equator or a pole. The accuracy of the projections is
* about 10 nm. * about 10 nm.
* *
* The ellipsoid parameters, the standard parallels, and the scale on the * The ellipsoid parameters, the standard parallels, and the scale on the
* standard parallels are set in the constructor. Internally, the case w ith * standard parallels are set in the constructor. Internally, the case w ith
* two standard parallels is converted into a single standard parallel, t he * two standard parallels is converted into a single standard parallel, t he
* latitude of tangency (also the latitude of minimum scale). This latit * latitude of tangency (also the latitude of minimum scale), with a scal
ude e
* is also used as the latitude of origin which is returned by * specified on this parallel. This latitude is also used as the latitud
* LambertConformalConic::OriginLatitude. The scale on the latitude of e of
* origin is given by LambertConformalConic::CentralScale. The case with * origin which is returned by LambertConformalConic::OriginLatitude. Th
two e
* distinct standard parallels where one is a pole is singular and is * scale on the latitude of origin is given by
* disallowed. The central meridian (which is a trivial shift * LambertConformalConic::CentralScale. The case with two distinct stand
* of the longitude) is specified as the \e lon0 argument of the ard
* LambertConformalConic::Forward and LambertConformalConic::Reverse * parallels where one is a pole is singular and is disallowed. The cent
* functions. There is no provision in this class for specifying a false ral
* easting or false northing or a different latitude of origin. However * meridian (which is a trivial shift of the longitude) is specified as t
* these are can be simply included by the calling function. For example he
the * \e lon0 argument of the LambertConformalConic::Forward and
* Pennsylvania South state coordinate system * LambertConformalConic::Reverse functions. There is no provision in th
* (<a href="http://www.spatialreference.org/ref/epsg/3364/"> EPSG:3364</ is
a>) * class for specifying a false easting or false northing or a different
* is obtained by: * latitude of origin. However these are can be simply included by the
\code * calling function. For example the Pennsylvania South state coordinate
const double * system (<a href="http://www.spatialreference.org/ref/epsg/3364/">
a = GeographicLib::Constants::WGS84_a(), r = 298.257222101, // GRS80 * EPSG:3364</a>) is obtained by: \code const double a =
lat1 = 39 + 56/60.0, lat1 = 40 + 58/60.0, // standard parallels * GeographicLib::Constants::WGS84_a(), r = 298.257222101, // GRS80 lat1
k1 = 1, // scale = 39
lat0 = 39 + 20/60.0, lon0 = 75 + 45/60.0, // origin * + 56/60.0, lat1 = 40 + 58/60.0, // standard parallels k1 = 1, // scale
fe = 600000, fn = 0; // false easting and northin * lat0 = 39 + 20/60.0, lon0 = 75 + 45/60.0, // origin fe = 600000, fn =
g * 0; // false easting and northing
// Set up basic projection // Set up basic projection
const GeographicLib::LambertConformalConic PASouth(a, r, lat1, lat2, k1) ; const GeographicLib::LambertConformalConic PASouth(a, r, lat1, lat2, k1) ;
double x0, y0; double x0, y0;
{ {
// Transform origin point // Transform origin point
PASouth.Forward(lon0, lat0, lon0, x0, y0); PASouth.Forward(lon0, lat0, lon0, x0, y0);
x0 -= fe; y0 -= fn; // Combine result with false origin x0 -= fe; y0 -= fn; // Combine result with false origin
} }
double lat, lon, x, y; double lat, lon, x, y;
// Sample conversion from geodetic to PASouth grid // Sample conversion from geodetic to PASouth grid
skipping to change at line 106 skipping to change at line 102
// Definition: Df(x,y) = (f(x)-f(y))/(x-y) // Definition: Df(x,y) = (f(x)-f(y))/(x-y)
// See: W. M. Kahan and R. J. Fateman, // See: W. M. Kahan and R. J. Fateman,
// Symbolic computation of divided differences, // Symbolic computation of divided differences,
// SIGSAM Bull. 33(3), 7-28 (1999) // SIGSAM Bull. 33(3), 7-28 (1999)
// http://doi.acm.org/10.1145/334714.334716 // http://doi.acm.org/10.1145/334714.334716
// http://www.cs.berkeley.edu/~fateman/papers/divdiff.pdf // http://www.cs.berkeley.edu/~fateman/papers/divdiff.pdf
// //
// General rules // General rules
// h(x) = f(g(x)): Dh(x,y) = Df(g(x),g(y))*Dg(x,y) // h(x) = f(g(x)): Dh(x,y) = Df(g(x),g(y))*Dg(x,y)
// h(x) = f(x)*g(x): // h(x) = f(x)*g(x):
// Dh(x,y) = Df(x,y)*(g(x)+g(y))/2 + Dg(x,y)*(f(x)+f(y))/2 // Dh(x,y) = Df(x,y)*g(x) + Dg(x,y)*f(y)
// = Df(x,y)*g(y) + Dg(x,y)*f(x)
// = Df(x,y)*(g(x)+g(y))/2 + Dg(x,y)*(f(x)+f(y))/2
// //
// hyp(x) = sqrt(1+x^2): Dhyp(x,y) = (x+y)/(hyp(x)+hyp(y)) // hyp(x) = sqrt(1+x^2): Dhyp(x,y) = (x+y)/(hyp(x)+hyp(y))
static inline real Dhyp(real x, real y, real hx, real hy) throw() static inline real Dhyp(real x, real y, real hx, real hy) throw()
// hx = hyp(x) // hx = hyp(x)
{ return (x + y) / (hx + hy); } { return (x + y) / (hx + hy); }
// sn(x) = x/sqrt(1+x^2): Dsn(x,y) = (x+y)/((sn(x)+sn(y))*(1+x^2)*(1+y^ 2)) // sn(x) = x/sqrt(1+x^2): Dsn(x,y) = (x+y)/((sn(x)+sn(y))*(1+x^2)*(1+y^ 2))
static inline real Dsn(real x, real y, real sx, real sy) throw() { static inline real Dsn(real x, real y, real sx, real sy) throw() {
// sx = x/hyp(x) // sx = x/hyp(x)
real t = x * y; real t = x * y;
return t > 0 ? (x+y) * sq( (sx*sy)/t ) / (sx+sy) : return t > 0 ? (x + y) * sq( (sx * sy)/t ) / (sx + sy) :
(x-y != 0 ? (sx-sy) / (x-y) : 1); (x - y != 0 ? (sx - sy) / (x - y) : 1);
} }
// Dlog1p(x,y) = log1p((x-y)/(1+y)/(x-y) // Dlog1p(x,y) = log1p((x-y)/(1+y)/(x-y)
static inline real Dlog1p(real x, real y) throw() { static inline real Dlog1p(real x, real y) throw() {
real t = x - y; if (t < 0) { t = -t; y = x; } real t = x - y; if (t < 0) { t = -t; y = x; }
return t != 0 ? Math::log1p(t/(1+y)) / t : 1/(1+x); return t != 0 ? Math::log1p(t / (1 + y)) / t : 1 / (1 + x);
} }
// Dexp(x,y) = exp((x+y)/2) * 2*sinh((x-y)/2)/(x-y) // Dexp(x,y) = exp((x+y)/2) * 2*sinh((x-y)/2)/(x-y)
static inline real Dexp(real x, real y) throw() { static inline real Dexp(real x, real y) throw() {
real t = (x - y)/2; real t = (x - y)/2;
return (t != 0 ? sinh(t)/t : real(1)) * exp((x + y)/2); return (t != 0 ? sinh(t)/t : real(1)) * exp((x + y)/2);
} }
// Dsinh(x,y) = 2*sinh((x-y)/2)/(x-y) * cosh((x+y)/2) // Dsinh(x,y) = 2*sinh((x-y)/2)/(x-y) * cosh((x+y)/2)
// cosh((x+y)/2) = (c+sinh(x)*sinh(y)/c)/2 // cosh((x+y)/2) = (c+sinh(x)*sinh(y)/c)/2
// c=sqrt((1+cosh(x))*(1+cosh(y))) // c=sqrt((1+cosh(x))*(1+cosh(y)))
// cosh((x+y)/2) = sqrt( (sinh(x)*sinh(y) + cosh(x)*cosh(y) + 1)/2 ) // cosh((x+y)/2) = sqrt( (sinh(x)*sinh(y) + cosh(x)*cosh(y) + 1)/2 )
skipping to change at line 152 skipping to change at line 150
// = asinh((x*sqrt(1+y^2)-y*sqrt(1+x^2)))/(x-y) // = asinh((x*sqrt(1+y^2)-y*sqrt(1+x^2)))/(x-y)
static inline real Dasinh(real x, real y, real hx, real hy) throw() { static inline real Dasinh(real x, real y, real hx, real hy) throw() {
// hx = hyp(x) // hx = hyp(x)
real t = x - y; real t = x - y;
return t != 0 ? return t != 0 ?
Math::asinh(x*y > 0 ? t * (x+y) / (x*hy + y*hx) : x*hy - y*hx) / t : Math::asinh(x*y > 0 ? t * (x+y) / (x*hy + y*hx) : x*hy - y*hx) / t :
1/hx; 1/hx;
} }
// Deatanhe(x,y) = eatanhe((x-y)/(1-e^2*x*y))/(x-y) // Deatanhe(x,y) = eatanhe((x-y)/(1-e^2*x*y))/(x-y)
inline real Deatanhe(real x, real y) const throw() { inline real Deatanhe(real x, real y) const throw() {
real t = x - y, d = (1 - _e2 * x * y); real t = x - y, d = 1 - _e2 * x * y;
return t != 0 ? eatanhe(t / d) / t : _e2 / d; return t != 0 ? eatanhe(t / d) / t : _e2 / d;
} }
void Init(real sphi1, real cphi1, real sphi2, real cphi2, real k1) thro w(); void Init(real sphi1, real cphi1, real sphi2, real cphi2, real k1) thro w();
public: public:
/** /**
* Constructor with a single standard parallel. * Constructor with a single standard parallel.
* *
* @param[in] a equatorial radius of ellipsoid (meters) * @param[in] a equatorial radius of ellipsoid (meters)
* @param[in] r reciprocal flattening of ellipsoid. Setting \e r = 0 * @param[in] r reciprocal flattening of ellipsoid. Setting \e r = 0
* implies \e r = inf or flattening = 0 (i.e., a sphere). Negative \ e r * implies \e r = inf or flattening = 0 (i.e., a sphere). Negative \ e r
* indicates a prolate ellipsoid. * indicates a prolate ellipsoid.
* @param[in] stdlat standard parallel (degrees), the circle of tangenc y. * @param[in] stdlat standard parallel (degrees), the circle of tangenc y.
* @param[in] k0 scale on standard parallel. * @param[in] k0 scale on the standard parallel.
* *
* An exception is thrown if \e a or \e k0 is not positive or if \e std lat * An exception is thrown if \e a or \e k0 is not positive or if \e std lat
* is not in the range [-90, 90]. * is not in the range [-90, 90].
********************************************************************** / ********************************************************************** /
LambertConformalConic(real a, real r, real stdlat, real k0); LambertConformalConic(real a, real r, real stdlat, real k0);
/** /**
* Constructor with two standard parallels. * Constructor with two standard parallels.
* *
* @param[in] a equatorial radius of ellipsoid (meters) * @param[in] a equatorial radius of ellipsoid (meters)
* @param[in] r reciprocal flattening of ellipsoid. Setting \e r = 0 * @param[in] r reciprocal flattening of ellipsoid. Setting \e r = 0
* implies \e r = inf or flattening = 0 (i.e., a sphere). Negative \ e r * implies \e r = inf or flattening = 0 (i.e., a sphere). Negative \ e r
* indicates a prolate ellipsoid. * indicates a prolate ellipsoid.
* @param[in] stdlat1 first standard parallel (degrees). * @param[in] stdlat1 first standard parallel (degrees).
* @param[in] stdlat2 second standard parallel (degrees). * @param[in] stdlat2 second standard parallel (degrees).
* @param[in] k1 scale on the standard parallels. * @param[in] k1 scale on the standard parallels.
* *
* An exception is thrown if \e a or \e k0 is not positive or if \e std lat1 * An exception is thrown if \e a or \e k0 is not positive or if \e std lat1
* or \e stdlat2 is not in the range [-90, 90]. In addition, if either \e * or \e stdlat2 is not in the range [-90, 90]. In addition, if either \e
* stdlat1 or \e stdlat2 is a pole, then an exception is thrown if \e * stdlat1 or \e stdlat2 is a pole, then an exception is thrown if \e
* stdlat1 is not equal \e stdlat2 * stdlat1 is not equal \e stdlat2.
********************************************************************** / ********************************************************************** /
LambertConformalConic(real a, real r, real stdlat1, real stdlat2, real k1); LambertConformalConic(real a, real r, real stdlat1, real stdlat2, real k1);
/** /**
* Constructor with two standard parallels specified by sines and cosin es. * Constructor with two standard parallels specified by sines and cosin es.
* *
* @param[in] a equatorial radius of ellipsoid (meters) * @param[in] a equatorial radius of ellipsoid (meters)
* @param[in] r reciprocal flattening of ellipsoid. Setting \e r = 0 * @param[in] r reciprocal flattening of ellipsoid. Setting \e r = 0
* implies \e r = inf or flattening = 0 (i.e., a sphere). Negative \ e r * implies \e r = inf or flattening = 0 (i.e., a sphere). Negative \ e r
* indicates a prolate ellipsoid. * indicates a prolate ellipsoid.
* @param[in] sinlat1 sine of first standard parallel. * @param[in] sinlat1 sine of first standard parallel.
* @param[in] coslat1 cosine of first standard parallel. * @param[in] coslat1 cosine of first standard parallel.
* @param[in] sinlat2 sine of second standard parallel. * @param[in] sinlat2 sine of second standard parallel.
* @param[in] coslat2 cosine of second standard parallel. * @param[in] coslat2 cosine of second standard parallel.
* @param[in] k1 scale on the standard parallels. * @param[in] k1 scale on the standard parallels.
* *
* This allows parallels close to the poles to be specified accurately. * This allows parallels close to the poles to be specified accurately.
* This routine computes the latitude of origin and the scale at this * This routine computes the latitude of origin and the scale at this
* latitude. In the case where \e lat1 and \e lat2 are different, the * latitude. In the case where \e lat1 and \e lat2 are different, the
* errors in this routines are as follows: if abs(\e lat2 - \e lat1) <= * errors in this routines are as follows: if \e dlat = abs(\e lat2 - \
* 160<sup>o</sup> and max(abs(\e lat1), abs(\e lat2)) <= e
* 89.9999<sup>o</sup>, then the error in the latitude of origin is les * lat1) <= 160<sup>o</sup> and max(abs(\e lat1), abs(\e lat2)) <= 90 -
s * min(0.0002, 2.2e-6(180 - \e dlat), 6e-8\e dlat<sup>2</sup>) (in
* than 4e-14d and the relative error in the scale is less than 7e-15. * degrees), then the error in the latitude of origin is less than
* 4.5e-14<sup>o</sup> and the relative error in the scale is less than
* 7e-15.
********************************************************************** / ********************************************************************** /
LambertConformalConic(real a, real r, LambertConformalConic(real a, real r,
real sinlat1, real coslat1, real sinlat1, real coslat1,
real sinlat2, real coslat2, real sinlat2, real coslat2,
real k1); real k1);
/** /**
* Set the scale for the projection. * Set the scale for the projection.
* *
* @param[in] lat (degrees). * @param[in] lat (degrees).
* @param[in] k scale at latitude \e lat (default 1). * @param[in] k scale at latitude \e lat (default 1).
* *
* This allows a "latitude of true scale" to be specified. An exceptio n is * This allows a "latitude of true scale" to be specified. An exceptio n is
* thrown if \e k is not positive. * thrown if \e k is not positive or if \e stdlat is not in the range [
-90,
* 90]
********************************************************************** / ********************************************************************** /
void SetScale(real lat, real k = real(1)); void SetScale(real lat, real k = real(1));
/** /**
* Forward projection, from geographic to Lambert conformal conic. * Forward projection, from geographic to Lambert conformal conic.
* *
* @param[in] lon0 central meridian longitude (degrees). * @param[in] lon0 central meridian longitude (degrees).
* @param[in] lat latitude of point (degrees). * @param[in] lat latitude of point (degrees).
* @param[in] lon longitude of point (degrees). * @param[in] lon longitude of point (degrees).
* @param[out] x easting of point (meters). * @param[out] x easting of point (meters).
skipping to change at line 318 skipping to change at line 319
* @return latitude of the origin for the projection (degrees). * @return latitude of the origin for the projection (degrees).
* *
* This is the latitude of minimum scale and equals the \e stdlat in th e * This is the latitude of minimum scale and equals the \e stdlat in th e
* 1-parallel constructor and lies between \e stdlat1 and \e stdlat2 in the * 1-parallel constructor and lies between \e stdlat1 and \e stdlat2 in the
* 2-parallel constructors. * 2-parallel constructors.
********************************************************************** / ********************************************************************** /
Math::real OriginLatitude() const throw() { return _lat0; } Math::real OriginLatitude() const throw() { return _lat0; }
/** /**
* @return central scale for the projection. This is the scale on the * @return central scale for the projection. This is the scale on the
* latitude of origin.. * latitude of origin.
********************************************************************** / ********************************************************************** /
Math::real CentralScale() const throw() { return _k0; } Math::real CentralScale() const throw() { return _k0; }
///@} ///@}
/** /**
* A global instantiation of LambertConformalConic with the WGS84 ellip * A global instantiation of LambertConformalConic with the WGS84
soid * ellipsoid, \e stdlat = 0, and \e k0 = 1. This degenerates to the
* and the UPS scale factor and \e stdlat = 0. This degenerates to the
* Mercator projection. * Mercator projection.
********************************************************************** / ********************************************************************** /
static const LambertConformalConic Mercator; static const LambertConformalConic Mercator;
}; };
} // namespace GeographicLib } // namespace GeographicLib
#endif #endif
 End of changes. 13 change blocks. 
47 lines changed or deleted 53 lines changed or added


 MGRS.hpp   MGRS.hpp 
/** /**
* \file MGRS.hpp * \file MGRS.hpp
* \brief Header for GeographicLib::MGRS class * \brief Header for GeographicLib::MGRS 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_MGRS_HPP) #if !defined(GEOGRAPHICLIB_MGRS_HPP)
#define GEOGRAPHICLIB_MGRS_HPP "$Id: MGRS.hpp 6876 2010-10-18 13:47:55Z kar ney $" #define GEOGRAPHICLIB_MGRS_HPP "$Id: MGRS.hpp 6911 2010-12-09 23:13:55Z kar ney $"
#include "GeographicLib/Constants.hpp" #include "GeographicLib/Constants.hpp"
#include "GeographicLib/UTMUPS.hpp" #include "GeographicLib/UTMUPS.hpp"
#include <sstream> #include <sstream>
namespace GeographicLib { namespace GeographicLib {
/** /**
* \brief Convert between UTM/UPS and %MGRS * \brief Convert between UTM/UPS and %MGRS
* *
skipping to change at line 186 skipping to change at line 186
* the zone number followed by one of [C&ndash;M] for the southern * the zone number followed by one of [C&ndash;M] for the southern
* hemisphere and [N&ndash;X] for the northern hemisphere. For \e zone = * hemisphere and [N&ndash;X] for the northern hemisphere. For \e zone =
* 0, the MGRS coordinates begins with one of [AB] for the southern * 0, the MGRS coordinates begins with one of [AB] for the southern
* hemisphere and [XY] for the northern hemisphere. * hemisphere and [XY] for the northern hemisphere.
* *
* The conversion to the MGRS is exact for prec in [0, 5] except that a * The conversion to the MGRS is exact for prec in [0, 5] except that a
* neighboring latitude band letter may be given if the point is within 5nm * neighboring latitude band letter may be given if the point is within 5nm
* of a band boundary. For prec in [6, 11], the conversion is accurate to * of a band boundary. For prec in [6, 11], the conversion is accurate to
* roundoff. * roundoff.
* *
* If \e x or \e y is NaN or if \e zone is UTMUPS::INVALID, the returne
d
* MGRS string is "INVALID".
*
* Return the result via a reference argument to avoid the overhead of * Return the result via a reference argument to avoid the overhead of
* allocating a potentially large number of small strings. If an error is * allocating a potentially large number of small strings. If an error is
* thrown, then \e mgrs is unchanged. * thrown, then \e mgrs is unchanged.
********************************************************************** / ********************************************************************** /
static void Forward(int zone, bool northp, real x, real y, static void Forward(int zone, bool northp, real x, real y,
int prec, std::string& mgrs); int prec, std::string& mgrs);
/** /**
* Convert UTM or UPS coordinate to an MGRS coordinate when the latitud e is * Convert UTM or UPS coordinate to an MGRS coordinate when the latitud e is
* known. * known.
 End of changes. 2 change blocks. 
1 lines changed or deleted 5 lines changed or added


 OSGB.hpp   OSGB.hpp 
/** /**
* \file OSGB.hpp * \file OSGB.hpp
* \brief Header for GeographicLib::OSGB class * \brief Header for GeographicLib::OSGB class
* *
* Copyright (c) Charles Karney (2008, 2009, 2010) <charles@karney.com> * Copyright (c) Charles Karney (2010) <charles@karney.com> and licensed un
* and licensed under the LGPL. For more information, see der
* http://geographiclib.sourceforge.net/ * the LGPL. For more information, see http://geographiclib.sourceforge.ne
t/
**********************************************************************/ **********************************************************************/
#if !defined(GEOGRAPHICLIB_OSGB_HPP) #if !defined(GEOGRAPHICLIB_OSGB_HPP)
#define GEOGRAPHICLIB_OSGB_HPP "$Id: OSGB.hpp 6881 2010-10-29 22:13:52Z kar ney $" #define GEOGRAPHICLIB_OSGB_HPP "$Id: OSGB.hpp 6911 2010-12-09 23:13:55Z kar ney $"
#include "GeographicLib/Constants.hpp" #include "GeographicLib/Constants.hpp"
#include "GeographicLib/TransverseMercator.hpp" #include "GeographicLib/TransverseMercator.hpp"
#include <string> #include <string>
#include <sstream> #include <sstream>
namespace GeographicLib { namespace GeographicLib {
/** /**
* \brief Ordnance Survery grid system for Great Britain * \brief Ordnance Survey grid system for Great Britain
* *
* The class implements the coordinate system used by the Ordnance Survey for * The class implements the coordinate system used by the Ordnance Survey for
* maps of Great Britain and conversions to the grid reference system. * maps of Great Britain and conversions to the grid reference system.
* *
* See * See
* - <a href="http://www.ordnancesurvey.co.uk/oswebsite/gps/docs/A_Guide_ to_Coordinate_Systems_in_Great_Britain.pdf"> * - <a href="http://www.ordnancesurvey.co.uk/oswebsite/gps/docs/A_Guide_ to_Coordinate_Systems_in_Great_Britain.pdf">
* A guide to coordinate systems in Great Britain</a> * A guide to coordinate systems in Great Britain</a>
* - <a href="http://www.ordnancesurvey.co.uk/oswebsite/gps/information/c oordinatesystemsinfo/guidetonationalgrid/page1.html"> * - <a href="http://www.ordnancesurvey.co.uk/oswebsite/gps/information/c oordinatesystemsinfo/guidetonationalgrid/page1.html">
* Guide to National Grid</a> * Guide to National Grid</a>
* *
* \b WARNING: the latitudes and longitudes for the Ordnance Survery grid * \b WARNING: the latitudes and longitudes for the Ordnance Survey grid
* system do not use the WGS84 datum. Do not use the values returned by this * system do not use the WGS84 datum. Do not use the values returned by this
* class in the UTMUPS, MGRS, or Geoid classes without first converting t he * class in the UTMUPS, MGRS, or Geoid classes without first converting t he
* datum (and vice versa). * datum (and vice versa).
**********************************************************************/ **********************************************************************/
class OSGB { class OSGB {
private: private:
typedef Math::real real; typedef Math::real real;
static const std::string letters; static const std::string letters;
static const std::string digits; static const std::string digits;
static const TransverseMercator OSGBTM; static const TransverseMercator OSGBTM;
 End of changes. 4 change blocks. 
6 lines changed or deleted 7 lines changed or added


 PolarStereographic.hpp   PolarStereographic.hpp 
/** /**
* \file PolarStereographic.hpp * \file PolarStereographic.hpp
* \brief Header for GeographicLib::PolarStereographic class * \brief Header for GeographicLib::PolarStereographic 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_POLARSTEREOGRAPHIC_HPP) #if !defined(GEOGRAPHICLIB_POLARSTEREOGRAPHIC_HPP)
#define GEOGRAPHICLIB_POLARSTEREOGRAPHIC_HPP "$Id: PolarStereographic.hpp 6 876 2010-10-18 13:47:55Z karney $" #define GEOGRAPHICLIB_POLARSTEREOGRAPHIC_HPP "$Id: PolarStereographic.hpp 6 906 2010-12-02 22:10:56Z karney $"
#include "GeographicLib/Constants.hpp" #include "GeographicLib/Constants.hpp"
namespace GeographicLib { namespace GeographicLib {
/** /**
* \brief Polar Stereographic Projection * \brief Polar Stereographic Projection
* *
* Implementation taken from the report, * Implementation taken from the report,
* - J. P. Snyder, * - J. P. Snyder,
skipping to change at line 62 skipping to change at line 62
* @param[in] k0 central scale factor. * @param[in] k0 central scale factor.
* *
* An exception is thrown if either of the axes of the ellipsoid is * An exception is thrown if either of the axes of the ellipsoid is
* not positive \e a or if \e k0 is not positive. * not positive \e a or if \e k0 is not positive.
********************************************************************** / ********************************************************************** /
PolarStereographic(real a, real r, real k0); PolarStereographic(real a, real r, real k0);
/** /**
* Set the scale for the projection. * Set the scale for the projection.
* *
* @param[in] lat (degrees). * @param[in] lat (degrees) assuming \e northp = true.
* @param[in] k scale at latitude \e lat (default 1). * @param[in] k scale at latitude \e lat (default 1).
* *
* This allows a "latitude of true scale" to be specified. An exceptio n is * This allows a "latitude of true scale" to be specified. An exceptio n is
* thrown if \e k is not positive. * thrown if \e k is not positive or if \e lat is not in the range (-90
,
* 90].
********************************************************************** / ********************************************************************** /
void SetScale(real lat, real k = real(1)); void SetScale(real lat, real k = real(1));
/** /**
* Forward projection, from geographic to polar stereographic. * Forward projection, from geographic to polar stereographic.
* *
* @param[in] northp the pole which is the center of projection (true m eans * @param[in] northp the pole which is the center of projection (true m eans
* north, false means south). * north, false means south).
* @param[in] lat latitude of point (degrees). * @param[in] lat latitude of point (degrees).
* @param[in] lon longitude of point (degrees). * @param[in] lon longitude of point (degrees).
 End of changes. 3 change blocks. 
3 lines changed or deleted 5 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 876 2010-10-18 13:47:55Z karney $" #define GEOGRAPHICLIB_TRANSVERSEMERCATOR_HPP "$Id: TransverseMercator.hpp 6 905 2010-12-01 21:28:56Z 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 TransverseMercator. * The order of the series approximation used in TransverseMercator.
* TM_TX_MAXPOW can be set to any integer in [4, 8]. * TM_TX_MAXPOW can be set to any integer in [4, 8].
**********************************************************************/ **********************************************************************/
#define TM_TX_MAXPOW \ #define TM_TX_MAXPOW \
(GEOGRAPHICLIB_PREC == 1 ? 6 : GEOGRAPHICLIB_PREC == 0 ? 4 : 8) (GEOGRAPHICLIB_PREC == 1 ? 6 : GEOGRAPHICLIB_PREC == 0 ? 4 : 8)
skipping to change at line 77 skipping to change at line 77
static const int maxpow = TM_TX_MAXPOW; static const int maxpow = TM_TX_MAXPOW;
static const real tol, overflow; static const real tol, overflow;
static const int numit = 5; static const int numit = 5;
const real _a, _r, _f, _k0, _e2, _e, _e2m, _c, _n; const real _a, _r, _f, _k0, _e2, _e, _e2m, _c, _n;
// _alp[0] and _bet[0] unused // _alp[0] and _bet[0] unused
real _a1, _b1, _alp[maxpow + 1], _bet[maxpow + 1]; real _a1, _b1, _alp[maxpow + 1], _bet[maxpow + 1];
static inline real sq(real x) throw() { return x * x; } static inline real sq(real x) throw() { return x * x; }
// tan(x) for x in [-pi/2, pi/2] ensuring that the sign is right // tan(x) for x in [-pi/2, pi/2] ensuring that the sign is right
static inline real tanx(real x) throw() { static inline real tanx(real x) throw() {
real t = std::tan(x); real t = std::tan(x);
return x >= 0 ? (t >= 0 ? t : overflow) : (t < 0 ? t : -overflow); // Write the tests this way to ensure that tanx(NaN()) is NaN()
return x >= 0 ? (!(t < 0) ? t : overflow) : (!(t >= 0) ? t : -overflo
w);
} }
// Return e * atanh(e * x) for f >= 0, else return // Return e * atanh(e * x) for f >= 0, else return
// - sqrt(-e2) * atan( sqrt(-e2) * x) for f < 0 // - sqrt(-e2) * atan( sqrt(-e2) * x) for f < 0
inline real eatanhe(real x) const throw() { inline real eatanhe(real x) const throw() {
return _f >= 0 ? _e * Math::atanh(_e * x) : - _e * std::atan(_e * x); return _f >= 0 ? _e * Math::atanh(_e * x) : - _e * std::atan(_e * x);
} }
public: public:
/** /**
* Constructor for a ellipsoid with * Constructor for a ellipsoid with
 End of changes. 2 change blocks. 
2 lines changed or deleted 4 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 6876 2010-10-18 13:47:55Z karney $" #define GEOGRAPHICLIB_TRANSVERSEMERCATOREXACT_HPP "$Id: TransverseMercatorE xact.hpp 6905 2010-12-01 21:28:56Z 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 72 skipping to change at line 72
typedef Math::real real; typedef Math::real real;
static const real tol, tol1, tol2, taytol, overflow; static const real tol, tol1, tol2, taytol, overflow;
static const int numit = 10; static const int numit = 10;
const real _a, _r, _f, _k0, _mu, _mv, _e, _ep2; const real _a, _r, _f, _k0, _mu, _mv, _e, _ep2;
const bool _extendp; const bool _extendp;
const EllipticFunction _Eu, _Ev; const EllipticFunction _Eu, _Ev;
static inline real sq(real x) throw() { return x * x; } static inline real sq(real x) throw() { return x * x; }
// tan(x) for x in [-pi/2, pi/2] ensuring that the sign is right // tan(x) for x in [-pi/2, pi/2] ensuring that the sign is right
static inline real tanx(real x) throw() { static inline real tanx(real x) throw() {
real t = std::tan(x); real t = std::tan(x);
return x >= 0 ? (t >= 0 ? t : overflow) : (t < 0 ? t : -overflow); // Write the tests this way to ensure that tanx(NaN()) is NaN()
return x >= 0 ? (!(t < 0) ? t : overflow) : (!(t >= 0) ? t : -overflo
w);
} }
real taup(real tau) const throw(); real taup(real tau) const throw();
real taupinv(real taup) const throw(); real taupinv(real taup) const throw();
void zeta(real u, real snu, real cnu, real dnu, void zeta(real u, real snu, real cnu, real dnu,
real v, real snv, real cnv, real dnv, real v, real snv, real cnv, real dnv,
real& taup, real& lam) const throw(); real& taup, real& lam) const throw();
void dwdzeta(real u, real snu, real cnu, real dnu, void dwdzeta(real u, real snu, real cnu, real dnu,
 End of changes. 2 change blocks. 
2 lines changed or deleted 4 lines changed or added


 UTMUPS.hpp   UTMUPS.hpp 
/** /**
* \file UTMUPS.hpp * \file UTMUPS.hpp
* \brief Header for GeographicLib::UTMUPS class * \brief Header for GeographicLib::UTMUPS 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_UTMUPS_HPP) #if !defined(GEOGRAPHICLIB_UTMUPS_HPP)
#define GEOGRAPHICLIB_UTMUPS_HPP "$Id: UTMUPS.hpp 6888 2010-11-13 15:19:42Z karney $" #define GEOGRAPHICLIB_UTMUPS_HPP "$Id: UTMUPS.hpp 6905 2010-12-01 21:28:56Z karney $"
#include "GeographicLib/Constants.hpp" #include "GeographicLib/Constants.hpp"
#include <sstream> #include <sstream>
namespace GeographicLib { namespace GeographicLib {
/** /**
* \brief Convert between Geographic coordinates and UTM/UPS * \brief Convert between Geographic coordinates and UTM/UPS
* *
* UTM and UPS are defined * UTM and UPS are defined
skipping to change at line 94 skipping to change at line 94
* The UTM divides the earth between latitudes -80 and 84 into 60 zones * The UTM divides the earth between latitudes -80 and 84 into 60 zones
* numbered 1 thru 60. Zone assign zone number 0 to the UPS regions, * numbered 1 thru 60. Zone assign zone number 0 to the UPS regions,
* covering the two poles. Within UTMUPS, non-negative zone numbers re fer * covering the two poles. Within UTMUPS, non-negative zone numbers re fer
* to one of the "physical" zones, 0 for UPS and [1, 60] for UTM. Nega tive * to one of the "physical" zones, 0 for UPS and [1, 60] for UTM. Nega tive
* "pseudo-zone" numbers are used to select one of the physical zones. * "pseudo-zone" numbers are used to select one of the physical zones.
********************************************************************** / ********************************************************************** /
enum zonespec { enum zonespec {
/** /**
* The smallest pseudo-zone number. * The smallest pseudo-zone number.
******************************************************************** **/ ******************************************************************** **/
MINPSEUDOZONE = -3, MINPSEUDOZONE = -4,
/**
* A marker for an undefined or invalid zone. Equivalent to NaN.
********************************************************************
**/
INVALID = -4,
/** /**
* If a coordinate already include zone information (e.g., it is an M GRS * If a coordinate already include zone information (e.g., it is an M GRS
* coordinate), use that, otherwise apply the UTMUPS::STANDARD rules. * coordinate), use that, otherwise apply the UTMUPS::STANDARD rules.
******************************************************************** **/ ******************************************************************** **/
MATCH = -3, MATCH = -3,
/** /**
* Apply the standard rules for UTM zone assigment extending the UTM zone * Apply the standard rules for UTM zone assigment extending the UTM zone
* to each pole to give a zone number in [1, 60]. For example, use U TM * to each pole to give a zone number in [1, 60]. For example, use U TM
* zone 38 for longitude in [42, 48). The rules include the Norway a nd * zone 38 for longitude in [42, 48). The rules include the Norway a nd
* Svalbard exceptions. * Svalbard exceptions.
skipping to change at line 152 skipping to change at line 156
/** /**
* The standard zone. * The standard zone.
* *
* @param[in] lat latitude (degrees). * @param[in] lat latitude (degrees).
* @param[in] lon longitude (degrees). * @param[in] lon longitude (degrees).
* @param[in] setzone zone override (optional) * @param[in] setzone zone override (optional)
* *
* This is exact. If the optional argument \e setzone is given then us e * This is exact. If the optional argument \e setzone is given then us e
* that zone if it is non-negative, otherwise apply the rules given in * that zone if it is non-negative, otherwise apply the rules given in
* UTMUPS::zonespec. Throws an error if \e setzone is outsize the rang e * UTMUPS::zonespec. Throws an error if \e setzone is outsize the rang e
* [UTMUPS::MINPSEUDOZONE, UTMUPS::MAXZONE] = [-3, 60]. * [UTMUPS::MINPSEUDOZONE, UTMUPS::MAXZONE] = [-4, 60].
********************************************************************** / ********************************************************************** /
static int StandardZone(real lat, real lon, int setzone = STANDARD); static int StandardZone(real lat, real lon, int setzone = STANDARD);
/** /**
* Forward projection, from geographic to UTM/UPS. * Forward projection, from geographic to UTM/UPS.
* *
* @param[in] lat latitude of point (degrees). * @param[in] lat latitude of point (degrees).
* @param[in] lon longitude of point (degrees). * @param[in] lon longitude of point (degrees).
* @param[out] zone the UTM zone (zero means UPS). * @param[out] zone the UTM zone (zero means UPS).
* @param[out] northp hemisphere of location (true means northern, fals e * @param[out] northp hemisphere of location (true means northern, fals e
skipping to change at line 258 skipping to change at line 262
* @param[out] zone the UTM zone (zero means UPS). * @param[out] zone the UTM zone (zero means UPS).
* @param[out] northp the hemisphere (true means northern, false * @param[out] northp the hemisphere (true means northern, false
* means southern). * means southern).
* *
* For UTM, \e zonestr has the form of a zone number in the range * For UTM, \e zonestr has the form of a zone number in the range
* [UTMUPS::MINUTMZONE, UTMUPS::MAXUTMZONE] = [1, 60] followed by a * [UTMUPS::MINUTMZONE, UTMUPS::MAXUTMZONE] = [1, 60] followed by a
* hemisphere letter, N or S. For UPS, it consists just of the hemisph ere * hemisphere letter, N or S. For UPS, it consists just of the hemisph ere
* letter. The returned value of \e zone is UTMUPS::UPS = 0 for UPS. Note * letter. The returned value of \e zone is UTMUPS::UPS = 0 for UPS. Note
* well that "38S" indicates the southern hemisphere of zone 38 and not * well that "38S" indicates the southern hemisphere of zone 38 and not
* latitude band S, [32, 40]. N, 01S, 2N, 38S are legal. 0N, 001S, 61 N, * latitude band S, [32, 40]. N, 01S, 2N, 38S are legal. 0N, 001S, 61 N,
* 38P are illegal. Throws an error is the zone string is malformed. * 38P are illegal. INV is a special value for which the returned valu
e of
* \e is UTMUPS::INVALID. Throws an error is the zone string is malfor
med.
********************************************************************** / ********************************************************************** /
static void DecodeZone(const std::string& zonestr, int& zone, bool& nor thp); static void DecodeZone(const std::string& zonestr, int& zone, bool& nor thp);
/** /**
* Encode a UTM/UPS zone string. * Encode a UTM/UPS zone string.
* *
* @param[out] zone the UTM zone (zero means UPS). * @param[out] zone the UTM zone (zero means UPS).
* @param[out] northp the hemisphere (true means northern, false * @param[out] northp the hemisphere (true means northern, false
* means southern). * means southern).
* @return string represention of zone and hemisphere. * @return string represention of zone and hemisphere.
* *
* \e zone must be in the range [UTMUPS::MINZONE, UTMUPS::MAXZONE] = [0 , * \e zone must be in the range [UTMUPS::MINZONE, UTMUPS::MAXZONE] = [0 ,
* 60] with \e zone = UTMUPS::UPS, 0, indicating UPS (but the resulting * 60] with \e zone = UTMUPS::UPS, 0, indicating UPS (but the resulting
* string does not contain "0"). This reverses UTMUPS::DecodeZone. * string does not contain "0"). \e zone may also be UTMUPS::INVALID,
in
* which case the returned string is "INV". This reverses
* UTMUPS::DecodeZone.
********************************************************************** / ********************************************************************** /
static std::string EncodeZone(int zone, bool northp); static std::string EncodeZone(int zone, bool northp);
/** /**
* @return shift (meters) necessary to align N and S halves of a UTM zo ne * @return shift (meters) necessary to align N and S halves of a UTM zo ne
* (10<sup>7</sup>). * (10<sup>7</sup>).
********************************************************************** / ********************************************************************** /
static Math::real UTMShift() throw(); static Math::real UTMShift() throw();
/** \name Inspector functions /** \name Inspector functions
 End of changes. 5 change blocks. 
5 lines changed or deleted 16 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/