summaryrefslogtreecommitdiff
path: root/boost/geometry/strategies/spherical/area_huiller.hpp
diff options
context:
space:
mode:
Diffstat (limited to 'boost/geometry/strategies/spherical/area_huiller.hpp')
-rw-r--r--boost/geometry/strategies/spherical/area_huiller.hpp58
1 files changed, 35 insertions, 23 deletions
diff --git a/boost/geometry/strategies/spherical/area_huiller.hpp b/boost/geometry/strategies/spherical/area_huiller.hpp
index e55fdb2b0e..37d8d20124 100644
--- a/boost/geometry/strategies/spherical/area_huiller.hpp
+++ b/boost/geometry/strategies/spherical/area_huiller.hpp
@@ -1,6 +1,12 @@
// Boost.Geometry (aka GGL, Generic Geometry Library)
-// Copyright (c) 2007-2012 Barend Gehrels, Amsterdam, the Netherlands.
+// Copyright (c) 2007-2015 Barend Gehrels, Amsterdam, the Netherlands.
+
+// This file was modified by Oracle on 2015.
+// Modifications copyright (c) 2015, Oracle and/or its affiliates.
+
+// Contributed and/or modified by Menelaos Karavelas, on behalf of Oracle
+// Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle
// Use, modification and distribution is subject to the Boost Software License,
// Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
@@ -31,18 +37,22 @@ namespace strategy { namespace area
\tparam PointOfSegment point type of segments of rings/polygons
\tparam CalculationType \tparam_calculation
\author Barend Gehrels. Adapted from:
-- http://www.soe.ucsc.edu/~pang/160/f98/Gems/GemsIV/sph_poly.c
+- http://webdocs.cs.ualberta.ca/~graphics/books/GraphicsGems/gemsiv/sph_poly.c
+- http://tog.acm.org/resources/GraphicsGems/gemsiv/sph_poly.c
- http://williams.best.vwh.net/avform.htm
-\note The version in Gems didn't account for polygons crossing the 180 meridian.
+\note The version in Graphics Gems IV (page 132-137) didn't account for
+polygons crossing the 0 and 180 meridians. The fix for this algorithm
+can be found in Graphics Gems V (pages 45-46). See:
+- http://kysmykseka.net/koti/wizardry/Game%20Development/Programming/Graphics%20Gems%204.pdf
+- http://kysmykseka.net/koti/wizardry/Game%20Development/Programming/Graphics%20Gems%205.pdf
\note This version works for convex and non-convex polygons, for 180 meridian
crossing polygons and for polygons with holes. However, some cases (especially
180 meridian cases) must still be checked.
\note The version which sums angles, which is often seen, doesn't handle non-convex
polygons correctly.
-\note The version which sums longitudes, see
-http://trs-new.jpl.nasa.gov/dspace/bitstream/2014/40409/1/07-03.pdf, is simple
-and works well in most cases but not in 180 meridian crossing cases. This probably
-could be solved.
+\note The version which sums longitudes, see http://hdl.handle.net/2014/40409,
+is simple and works well in most cases but not in 180 meridian crossing cases.
+This probably could be solved.
\note This version is made for spherical equatorial coordinate systems
@@ -113,8 +123,12 @@ public :
calculation_type const half = 0.5;
calculation_type const two = 2.0;
calculation_type const four = 4.0;
- calculation_type const two_pi = two * geometry::math::pi<calculation_type>();
- calculation_type const half_pi = half * geometry::math::pi<calculation_type>();
+ calculation_type const pi
+ = geometry::math::pi<calculation_type>();
+ calculation_type const two_pi
+ = geometry::math::two_pi<calculation_type>();
+ calculation_type const half_pi
+ = geometry::math::half_pi<calculation_type>();
// Distance p1 p2
calculation_type a = state.distance_over_unit_sphere.apply(p1, p2);
@@ -128,33 +142,31 @@ public :
// E: spherical excess, using l'Huiller's formula
// [tg(e / 4)]2 = tg[s / 2] tg[(s-a) / 2] tg[(s-b) / 2] tg[(s-c) / 2]
- calculation_type E = four
+ calculation_type excess = four
* atan(geometry::math::sqrt(geometry::math::abs(tan(s / two)
* tan((s - a) / two)
* tan((s - b) / two)
* tan((s - c) / two))));
- E = geometry::math::abs(E);
+ excess = geometry::math::abs(excess);
// In right direction: positive, add area. In left direction: negative, subtract area.
- // Longitude comparisons are not so obvious. If one is negative, other is positive,
+ // Longitude comparisons are not so obvious. If one is negative and other is positive,
// we have to take the dateline into account.
- // TODO: check this / enhance this, should be more robust. See also the "grow" for ll
- // TODO: use minmax or "smaller"/"compare" strategy for this
- calculation_type lon1 = geometry::get_as_radian<0>(p1) < 0
- ? geometry::get_as_radian<0>(p1) + two_pi
- : geometry::get_as_radian<0>(p1);
- calculation_type lon2 = geometry::get_as_radian<0>(p2) < 0
- ? geometry::get_as_radian<0>(p2) + two_pi
- : geometry::get_as_radian<0>(p2);
+ calculation_type lon_diff = geometry::get_as_radian<0>(p2)
+ - geometry::get_as_radian<0>(p1);
+ if (lon_diff <= 0)
+ {
+ lon_diff += two_pi;
+ }
- if (lon2 < lon1)
+ if (lon_diff > pi)
{
- E = -E;
+ excess = -excess;
}
- state.sum += E;
+ state.sum += excess;
}
}