diff options
Diffstat (limited to 'boost/geometry/strategies/geographic/distance_andoyer.hpp')
| -rw-r--r-- | boost/geometry/strategies/geographic/distance_andoyer.hpp | 71 |
1 files changed, 13 insertions, 58 deletions
diff --git a/boost/geometry/strategies/geographic/distance_andoyer.hpp b/boost/geometry/strategies/geographic/distance_andoyer.hpp index 64de8c1a41..1646727d09 100644 --- a/boost/geometry/strategies/geographic/distance_andoyer.hpp +++ b/boost/geometry/strategies/geographic/distance_andoyer.hpp @@ -1,9 +1,9 @@ // Boost.Geometry (aka GGL, Generic Geometry Library) -// Copyright (c) 2007-2012 Barend Gehrels, Amsterdam, the Netherlands. +// Copyright (c) 2007-2016 Barend Gehrels, Amsterdam, the Netherlands. -// This file was modified by Oracle on 2014. -// Modifications copyright (c) 2014 Oracle and/or its affiliates. +// This file was modified by Oracle on 2014, 2016. +// Modifications copyright (c) 2014-2016 Oracle and/or its affiliates. // Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle @@ -20,6 +20,7 @@ #include <boost/geometry/core/radius.hpp> #include <boost/geometry/core/srs.hpp> +#include <boost/geometry/algorithms/detail/andoyer_inverse.hpp> #include <boost/geometry/algorithms/detail/flattening.hpp> #include <boost/geometry/strategies/distance.hpp> @@ -50,6 +51,8 @@ are about the same as Vincenty. In my (Barend's) testcases the results didn't di \see http://www.codeguru.com/Cpp/Cpp/algorithms/article.php/c5115 (implementation) \see http://futureboy.homeip.net/frinksamp/navigation.frink (implementation) \see http://www.voidware.com/earthdist.htm (implementation) +\see http://www.dtic.mil/docs/citations/AD0627893 +\see http://www.dtic.mil/docs/citations/AD703541 */ template < @@ -82,16 +85,17 @@ public : : m_spheroid(spheroid) {} - template <typename Point1, typename Point2> inline typename calculation_type<Point1, Point2>::type apply(Point1 const& point1, Point2 const& point2) const { - return calc<typename calculation_type<Point1, Point2>::type> - ( - get_as_radian<0>(point1), get_as_radian<1>(point1), - get_as_radian<0>(point2), get_as_radian<1>(point2) - ); + return geometry::detail::andoyer_inverse + < + typename calculation_type<Point1, Point2>::type, + true, false + >::apply(get_as_radian<0>(point1), get_as_radian<1>(point1), + get_as_radian<0>(point2), get_as_radian<1>(point2), + m_spheroid).distance; } inline Spheroid const& model() const @@ -100,55 +104,6 @@ public : } private : - template <typename CT, typename T> - inline CT calc(T const& lon1, - T const& lat1, - T const& lon2, - T const& lat2) const - { - CT const G = (lat1 - lat2) / 2.0; - CT const lambda = (lon1 - lon2) / 2.0; - - if (geometry::math::equals(lambda, 0.0) - && geometry::math::equals(G, 0.0)) - { - return 0.0; - } - - CT const F = (lat1 + lat2) / 2.0; - - CT const sinG2 = math::sqr(sin(G)); - CT const cosG2 = math::sqr(cos(G)); - CT const sinF2 = math::sqr(sin(F)); - CT const cosF2 = math::sqr(cos(F)); - CT const sinL2 = math::sqr(sin(lambda)); - CT const cosL2 = math::sqr(cos(lambda)); - - CT const S = sinG2 * cosL2 + cosF2 * sinL2; - CT const C = cosG2 * cosL2 + sinF2 * sinL2; - - CT const c0 = 0; - CT const c1 = 1; - CT const c2 = 2; - CT const c3 = 3; - - if (geometry::math::equals(S, c0) || geometry::math::equals(C, c0)) - { - return c0; - } - - CT const radius_a = CT(get_radius<0>(m_spheroid)); - CT const flattening = geometry::detail::flattening<CT>(m_spheroid); - - CT const omega = atan(math::sqrt(S / C)); - CT const r3 = c3 * math::sqrt(S * C) / omega; // not sure if this is r or greek nu - CT const D = c2 * omega * radius_a; - CT const H1 = (r3 - c1) / (c2 * C); - CT const H2 = (r3 + c1) / (c2 * S); - - return D * (c1 + flattening * (H1 * sinF2 * cosG2 - H2 * cosF2 * sinG2) ); - } - Spheroid m_spheroid; }; |