diff options
Diffstat (limited to 'boost/geometry/algorithms/detail/envelope/segment.hpp')
| -rw-r--r-- | boost/geometry/algorithms/detail/envelope/segment.hpp | 209 |
1 files changed, 93 insertions, 116 deletions
diff --git a/boost/geometry/algorithms/detail/envelope/segment.hpp b/boost/geometry/algorithms/detail/envelope/segment.hpp index 6186e72a3a..7631e84883 100644 --- a/boost/geometry/algorithms/detail/envelope/segment.hpp +++ b/boost/geometry/algorithms/detail/envelope/segment.hpp @@ -4,9 +4,10 @@ // Copyright (c) 2008-2015 Bruno Lalande, Paris, France. // Copyright (c) 2009-2015 Mateusz Loskot, London, UK. -// This file was modified by Oracle on 2015, 2016. -// Modifications copyright (c) 2015-2016, Oracle and/or its affiliates. +// This file was modified by Oracle on 2015-2017. +// Modifications copyright (c) 2015-2017, Oracle and/or its affiliates. +// Contributed and/or modified by Vissarion Fysikopoulos, on behalf of Oracle // Contributed and/or modified by Menelaos Karavelas, on behalf of Oracle // Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle @@ -26,6 +27,7 @@ #include <boost/geometry/core/coordinate_system.hpp> #include <boost/geometry/core/coordinate_type.hpp> #include <boost/geometry/core/cs.hpp> +#include <boost/geometry/core/srs.hpp> #include <boost/geometry/core/point_type.hpp> #include <boost/geometry/core/radian_access.hpp> #include <boost/geometry/core/tags.hpp> @@ -34,10 +36,9 @@ #include <boost/geometry/geometries/helper_geometry.hpp> -#include <boost/geometry/strategies/compare.hpp> +#include <boost/geometry/formulas/vertex_latitude.hpp> #include <boost/geometry/algorithms/detail/assign_indexed_point.hpp> -#include <boost/geometry/algorithms/detail/normalize.hpp> #include <boost/geometry/algorithms/detail/envelope/point.hpp> #include <boost/geometry/algorithms/detail/envelope/transform_units.hpp> @@ -46,7 +47,6 @@ #include <boost/geometry/algorithms/dispatch/envelope.hpp> - namespace boost { namespace geometry { @@ -54,63 +54,36 @@ namespace boost { namespace geometry namespace detail { namespace envelope { - -template <std::size_t Dimension, std::size_t DimensionCount> -struct envelope_one_segment +template <typename CalculationType, typename CS_Tag> +struct envelope_segment_call_vertex_latitude { - template<typename Point, typename Box> - static inline void apply(Point const& p1, Point const& p2, Box& mbr) + template <typename T1, typename T2, typename Strategy> + static inline CalculationType apply(T1 const& lat1, + T2 const& alp1, + Strategy const& ) { - envelope_one_point<Dimension, DimensionCount>::apply(p1, mbr); - detail::expand::point_loop - < - strategy::compare::default_strategy, - strategy::compare::default_strategy, - Dimension, - DimensionCount - >::apply(mbr, p2); + return geometry::formula::vertex_latitude<CalculationType, CS_Tag> + ::apply(lat1, alp1); } }; - -// Computes the MBR of a segment given by (lon1, lat1) and (lon2, -// lat2), and with azimuths a1 and a2 at the two endpoints of the -// segment. -// It is assumed that the spherical coordinates of the segment are -// normalized and in radians. -// The longitudes and latitudes of the endpoints are overridden by -// those of the box. -class compute_mbr_of_segment +template <typename CalculationType> +struct envelope_segment_call_vertex_latitude<CalculationType, geographic_tag> { -private: - // computes the azimuths of the segment with endpoints (lon1, lat1) - // and (lon2, lat2) - // radians - template <typename CalculationType> - static inline void azimuths(CalculationType const& lon1, - CalculationType const& lat1, - CalculationType const& lon2, - CalculationType const& lat2, - CalculationType& a1, - CalculationType& a2) + template <typename T1, typename T2, typename Strategy> + static inline CalculationType apply(T1 const& lat1, + T2 const& alp1, + Strategy const& strategy) { - BOOST_GEOMETRY_ASSERT(lon1 <= lon2); - - CalculationType dlon = lon2 - lon1; - CalculationType sin_dlon = sin(dlon); - CalculationType cos_dlon = cos(dlon); - CalculationType cos_lat1 = cos(lat1); - CalculationType cos_lat2 = cos(lat2); - CalculationType sin_lat1 = sin(lat1); - CalculationType sin_lat2 = sin(lat2); - - a1 = atan2(sin_dlon * cos_lat2, - cos_lat1 * sin_lat2 - sin_lat1 * cos_lat2 * cos_dlon); - - a2 = atan2(-sin_dlon * cos_lat1, - cos_lat2 * sin_lat1 - sin_lat2 * cos_lat1 * cos_dlon); - a2 += math::pi<CalculationType>(); + return geometry::formula::vertex_latitude<CalculationType, geographic_tag> + ::apply(lat1, alp1, strategy.model()); } +}; + +template <typename CS_Tag> +class envelope_segment_impl +{ +private: // degrees or radians template <typename CalculationType> @@ -134,8 +107,8 @@ private: static CalculationType const pi_half = math::half_pi<CalculationType>(); return (a1 < a2) - ? (a1 < pi_half && pi_half < a2) - : (a1 > pi_half && pi_half > a2); + ? (a1 < pi_half && pi_half < a2) + : (a1 > pi_half && pi_half > a2); } // radians or degrees @@ -151,21 +124,13 @@ private: return math::abs(lon1 - lon2) > constants::half_period(); // > pi } - // radians - template <typename CalculationType> - static inline CalculationType max_latitude(CalculationType const& azimuth, - CalculationType const& latitude) - { - // azimuth and latitude are assumed to be in radians - return acos( math::abs(cos(latitude) * sin(azimuth)) ); - } - // degrees or radians - template <typename Units, typename CalculationType> + template <typename Units, typename CalculationType, typename Strategy> static inline void compute_box_corners(CalculationType& lon1, CalculationType& lat1, CalculationType& lon2, - CalculationType& lat2) + CalculationType& lat2, + Strategy const& strategy) { // coordinates are assumed to be in radians BOOST_GEOMETRY_ASSERT(lon1 <= lon2); @@ -175,13 +140,14 @@ private: CalculationType lon2_rad = math::as_radian<Units>(lon2); CalculationType lat2_rad = math::as_radian<Units>(lat2); - CalculationType a1 = 0, a2 = 0; - azimuths(lon1_rad, lat1_rad, lon2_rad, lat2_rad, a1, a2); + CalculationType a1, a2; + strategy.apply(lon1_rad, lat1_rad, lon2_rad, lat2_rad, a1, a2); if (lat1 > lat2) { std::swap(lat1, lat2); std::swap(lat1_rad, lat2_rad); + std::swap(a1, a2); } if (math::equals(a1, a2)) @@ -192,12 +158,16 @@ private: if (contains_pi_half(a1, a2)) { + CalculationType p_max = envelope_segment_call_vertex_latitude + <CalculationType, CS_Tag>::apply(lat1_rad, a1, strategy); + CalculationType const mid_lat = lat1 + lat2; if (mid_lat < 0) { // update using min latitude - CalculationType const lat_min_rad = -max_latitude(a1, lat1_rad); - CalculationType const lat_min = math::from_radian<Units>(lat_min_rad); + CalculationType const lat_min_rad = -p_max; + CalculationType const lat_min + = math::from_radian<Units>(lat_min_rad); if (lat1 > lat_min) { @@ -207,8 +177,9 @@ private: else if (mid_lat > 0) { // update using max latitude - CalculationType const lat_max_rad = max_latitude(a1, lat1_rad); - CalculationType const lat_max = math::from_radian<Units>(lat_max_rad); + CalculationType const lat_max_rad = p_max; + CalculationType const lat_max + = math::from_radian<Units>(lat_max_rad); if (lat2 < lat_max) { @@ -218,11 +189,12 @@ private: } } - template <typename Units, typename CalculationType> + template <typename Units, typename CalculationType, typename Strategy> static inline void apply(CalculationType& lon1, CalculationType& lat1, CalculationType& lon2, - CalculationType& lat2) + CalculationType& lat2, + Strategy const& strategy) { typedef math::detail::constants_on_spheroid < @@ -278,16 +250,22 @@ private: swap(lon1, lat1, lon2, lat2); } - compute_box_corners<Units>(lon1, lat1, lon2, lat2); + compute_box_corners<Units>(lon1, lat1, lon2, lat2, strategy); } public: - template <typename Units, typename CalculationType, typename Box> + template < + typename Units, + typename CalculationType, + typename Box, + typename Strategy + > static inline void apply(CalculationType lon1, CalculationType lat1, CalculationType lon2, CalculationType lat2, - Box& mbr) + Box& mbr, + Strategy const& strategy) { typedef typename coordinate_type<Box>::type box_coordinate_type; @@ -298,7 +276,7 @@ public: helper_box_type radian_mbr; - apply<Units>(lon1, lat1, lon2, lat2); + apply<Units>(lon1, lat1, lon2, lat2, strategy); geometry::set < @@ -324,29 +302,42 @@ public: } }; +template <std::size_t Dimension, std::size_t DimensionCount> +struct envelope_one_segment +{ + template<typename Point, typename Box, typename Strategy> + static inline void apply(Point const& p1, + Point const& p2, + Box& mbr, + Strategy const& strategy) + { + envelope_one_point<Dimension, DimensionCount>::apply(p1, mbr, strategy); + detail::expand::point_loop + < + strategy::compare::default_strategy, + strategy::compare::default_strategy, + Dimension, + DimensionCount + >::apply(mbr, p2, strategy); + } +}; + template <std::size_t DimensionCount> -struct envelope_segment_on_sphere +struct envelope_segment { - template <typename Point, typename Box> - static inline void apply(Point const& p1, Point const& p2, Box& mbr) + template <typename Point, typename Box, typename Strategy> + static inline void apply(Point const& p1, + Point const& p2, + Box& mbr, + Strategy const& strategy) { // first compute the envelope range for the first two coordinates - Point p1_normalized = detail::return_normalized<Point>(p1); - Point p2_normalized = detail::return_normalized<Point>(p2); - - typedef typename coordinate_system<Point>::type::units units_type; - - compute_mbr_of_segment::template apply<units_type>( - geometry::get<0>(p1_normalized), - geometry::get<1>(p1_normalized), - geometry::get<0>(p2_normalized), - geometry::get<1>(p2_normalized), - mbr); + strategy.apply(p1, p2, mbr); // now compute the envelope range for coordinates of // dimension 2 and higher - envelope_one_segment<2, DimensionCount>::apply(p1, p2, mbr); + envelope_one_segment<2, DimensionCount>::apply(p1, p2, mbr, strategy); } template <typename Segment, typename Box> @@ -359,21 +350,6 @@ struct envelope_segment_on_sphere } }; - - -template <std::size_t DimensionCount, typename CS_Tag> -struct envelope_segment - : envelope_one_segment<0, DimensionCount> -{}; - - -template <std::size_t DimensionCount> -struct envelope_segment<DimensionCount, spherical_equatorial_tag> - : envelope_segment_on_sphere<DimensionCount> -{}; - - - }} // namespace detail::envelope #endif // DOXYGEN_NO_DETAIL @@ -383,23 +359,24 @@ namespace dispatch { -template <typename Segment, typename CS_Tag> -struct envelope<Segment, segment_tag, CS_Tag> +template <typename Segment> +struct envelope<Segment, segment_tag> { - template <typename Box> - static inline void apply(Segment const& segment, Box& mbr) + template <typename Box, typename Strategy> + static inline void apply(Segment const& segment, + Box& mbr, + Strategy const& strategy) { typename point_type<Segment>::type p[2]; detail::assign_point_from_index<0>(segment, p[0]); detail::assign_point_from_index<1>(segment, p[1]); detail::envelope::envelope_segment < - dimension<Segment>::value, CS_Tag - >::apply(p[0], p[1], mbr); + dimension<Segment>::value + >::apply(p[0], p[1], mbr, strategy); } }; - } // namespace dispatch #endif // DOXYGEN_NO_DISPATCH |