diff --git a/doc/doxy/doxygen_input/groups/groups.hpp b/doc/doxy/doxygen_input/groups/groups.hpp index 9889ebd76..fabf1f376 100644 --- a/doc/doxy/doxygen_input/groups/groups.hpp +++ b/doc/doxy/doxygen_input/groups/groups.hpp @@ -2,6 +2,10 @@ // // Copyright (c) 2011-2012 Barend Gehrels, Amsterdam, the Netherlands. +// This file was modified by Oracle on 2018. +// Modifications copyright (c) 2018, Oracle and/or its affiliates. +// 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 // http://www.boost.org/LICENSE_1_0.txt) @@ -26,6 +30,7 @@ \defgroup covered_by covered_by: detect if a geometry is inside or on the border of another geometry, a.o. point-in-polygon (border included) \defgroup crosses crosses: detect if two geometries crosses each other \defgroup cs coordinate systems +\defgroup densify densify: add points to geometry, keeping shape \defgroup difference difference: difference of two geometries \defgroup disjoint disjoint: detect if geometries are not spatially related \defgroup distance distance: calculate distance between two geometries diff --git a/doc/imports.qbk b/doc/imports.qbk index 6c4863add..fd8aab313 100644 --- a/doc/imports.qbk +++ b/doc/imports.qbk @@ -5,6 +5,9 @@ Copyright (c) 2009-2012 Mateusz Loskot, London, UK. Copyright (c) 2009-2012 Bruno Lalande, Paris, France. + Copyright (c) 2018, Oracle and/or its affiliates. + 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 http://www.boost.org/LICENSE_1_0.txt) @@ -27,6 +30,8 @@ [import src/examples/algorithms/convert.cpp] [import src/examples/algorithms/convex_hull.cpp] [import src/examples/algorithms/correct.cpp] +[import src/examples/algorithms/densify.cpp] +[import src/examples/algorithms/densify_strategy.cpp] [import src/examples/algorithms/distance.cpp] [import src/examples/algorithms/difference.cpp] [import src/examples/algorithms/envelope.cpp] diff --git a/doc/make_qbk.py b/doc/make_qbk.py index cace58ff8..10824e120 100755 --- a/doc/make_qbk.py +++ b/doc/make_qbk.py @@ -5,6 +5,9 @@ # Copyright (c) 2008-2012 Bruno Lalande, Paris, France. # Copyright (c) 2009-2012 Mateusz Loskot (mateusz@loskot.net), London, UK # Copyright (c) 2017 Adam Wulkiewicz, Lodz, Poland +# +# Copyright (c) 2018, Oracle and/or its affiliates. +# 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 @@ -91,7 +94,7 @@ call_doxygen() algorithms = ["append", "assign", "make", "clear" , "area", "buffer", "centroid", "convert", "correct", "covered_by" - , "convex_hull", "crosses", "difference", "disjoint", "distance" + , "convex_hull", "crosses", "densify", "difference", "disjoint", "distance" , "envelope", "equals", "expand", "for_each", "is_empty" , "is_simple", "is_valid", "intersection", "intersects", "length" , "num_geometries", "num_interior_rings", "num_points" @@ -120,7 +123,8 @@ models = ["point", "linestring", "box" , "multi_linestring", "multi_point", "multi_polygon", "referring_segment"] -strategies = ["distance::pythagoras", "distance::pythagoras_box_box" +strategies = ["densify::cartesian", "densify::geographic", "densify::spherical" + , "distance::pythagoras", "distance::pythagoras_box_box" , "distance::pythagoras_point_box", "distance::haversine" , "distance::cross_track", "distance::cross_track_point_box" , "distance::projected_point" diff --git a/doc/quickref.xml b/doc/quickref.xml index 755fdd255..7c96b545f 100644 --- a/doc/quickref.xml +++ b/doc/quickref.xml @@ -10,8 +10,8 @@ Copyright (c) 2009-2015 Bruno Lalande, Paris, France. Copyright (c) 2013-2015 Adam Wulkiewicz, Lodz, Poland. - This file was modified by Oracle on 2014, 2015, 2017. - Modifications copyright (c) 2014-2017, Oracle and/or its affiliates. + This file was modified by Oracle on 2014, 2015, 2017, 2018. + Modifications copyright (c) 2014-2018, 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 @@ -389,6 +389,10 @@ + Densify + + densify + Distance distance @@ -410,17 +414,16 @@ for each (point, segment) + + Intersection intersection - - Length length - - + Num_ (counting) num_interior_rings @@ -428,8 +431,6 @@ num_points num_segments - - Perimeter perimeter @@ -439,12 +440,12 @@ relate relation + + Reverse reverse - - Simplify simplify @@ -554,6 +555,14 @@ strategy::convex_hull::graham_andrew + + Densify + + strategy::densify::cartesian + strategy::densify::geographic + strategy::densify::spherical + + Distance @@ -566,6 +575,8 @@ strategy::distance::haversine + + Side @@ -574,8 +585,6 @@ strategy::side::spherical_side_formula - - Simplify @@ -593,7 +602,9 @@ strategy::transform::rotate_transformer - + + + Within strategy::winding diff --git a/doc/reference.qbk b/doc/reference.qbk index 699712399..865989929 100644 --- a/doc/reference.qbk +++ b/doc/reference.qbk @@ -6,8 +6,8 @@ Copyright (c) 2009-2017 Bruno Lalande, Paris, France. Copyright (c) 2013-2017 Adam Wulkiewicz, Lodz, Poland. - This file was modified by Oracle on 2014, 2017. - Modifications copyright (c) 2014-2017, Oracle and/or its affiliates. + This file was modified by Oracle on 2014, 2017, 2018. + Modifications copyright (c) 2014-2018, 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 @@ -102,6 +102,10 @@ [include generated/crosses.qbk] [endsect] +[section:densify densify] +[include generated/densify.qbk] +[endsect] + [section:difference difference] [include generated/difference.qbk] [endsect] @@ -332,6 +336,9 @@ [section:strategies Strategies] +[include generated/densify_cartesian.qbk] +[include generated/densify_geographic.qbk] +[include generated/densify_spherical.qbk] [include generated/distance_pythagoras.qbk] [include generated/distance_pythagoras_box_box.qbk] [include generated/distance_pythagoras_point_box.qbk] diff --git a/doc/reference/algorithms/densify.qbk b/doc/reference/algorithms/densify.qbk new file mode 100644 index 000000000..676e59a13 --- /dev/null +++ b/doc/reference/algorithms/densify.qbk @@ -0,0 +1,25 @@ +[/============================================================================ + Boost.Geometry + + Copyright (c) 2018, Oracle and/or its affiliates. + 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 + http://www.boost.org/LICENSE_1_0.txt) +=============================================================================/] + +[def __this_function__ densify] + +[heading_conformance_no_ogc __this_function__] +[note PostGIS contains an algorithm ST_Segmentize with the same functionality. + See the [@http://www.postgis.org/docs/ST_Segmentize.html PostGIS documentation]. +] + +[heading Behavior] +The algorithm divides segments of a geometry if they are longer than passed +distance into smaller segments. + +[note The units of the distance depends on strategy. In order to change the +default behavior a user has to create a strategy and pass it explicitly into +the algorithm.] diff --git a/doc/src/examples/algorithms/Jamfile.v2 b/doc/src/examples/algorithms/Jamfile.v2 index 38c39dbd1..9b62521e1 100644 --- a/doc/src/examples/algorithms/Jamfile.v2 +++ b/doc/src/examples/algorithms/Jamfile.v2 @@ -4,10 +4,11 @@ # Copyright (c) 2008-2015 Bruno Lalande, Paris, France. # Copyright (c) 2009-2015 Mateusz Loskot, London, UK. -# This file was modified by Oracle on 2014, 2015. -# Modifications copyright (c) 2014-2015, Oracle and/or its affiliates. +# This file was modified by Oracle on 2014, 2015, 2018. +# Modifications copyright (c) 2014-2018, 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 @@ -40,6 +41,8 @@ exe convert : convert.cpp ; exe convex_hull : convex_hull.cpp ; exe correct : correct.cpp ; +exe densify : densify.cpp ; +exe densify_strategy : densify_strategy.cpp ; exe difference : difference.cpp ; exe distance : distance.cpp ; diff --git a/doc/src/examples/algorithms/densify.cpp b/doc/src/examples/algorithms/densify.cpp new file mode 100644 index 000000000..68f7426c6 --- /dev/null +++ b/doc/src/examples/algorithms/densify.cpp @@ -0,0 +1,47 @@ +// Boost.Geometry +// QuickBook Example + +// Copyright (c) 2018, Oracle and/or its affiliates +// 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 +// http://www.boost.org/LICENSE_1_0.txt) + +//[densify +//` Shows how to densify a polygon + +#include + +#include +#include +#include + +int main() +{ + typedef boost::geometry::model::d2::point_xy point_type; + typedef boost::geometry::model::polygon polygon_type; + + polygon_type poly; + boost::geometry::read_wkt( + "POLYGON((0 0,0 10,10 10,10 0,0 0),(1 1,4 1,4 4,1 4,1 1))", poly); + + polygon_type res; + + boost::geometry::densify(poly, res, 6.0); + + std::cout << "densified: " << boost::geometry::wkt(res) << std::endl; + + return 0; +} + +//] + +//[densify_output +/*` +Output: +[pre +densified: POLYGON((0 0,0 5,0 10,5 10,10 10,10 5,10 0,5 0,0 0),(1 1,4 1,4 4,1 4,1 1)) +] +*/ +//] diff --git a/doc/src/examples/algorithms/densify_strategy.cpp b/doc/src/examples/algorithms/densify_strategy.cpp new file mode 100644 index 000000000..eb9934995 --- /dev/null +++ b/doc/src/examples/algorithms/densify_strategy.cpp @@ -0,0 +1,50 @@ +// Boost.Geometry +// QuickBook Example + +// Copyright (c) 2018, Oracle and/or its affiliates +// 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 +// http://www.boost.org/LICENSE_1_0.txt) + +//[densify_strategy +//` Shows how to densify a linestring in geographic coordinate system + +#include + +#include +#include +#include + +int main() +{ + namespace bg = boost::geometry; + typedef bg::model::point > point_type; + typedef bg::model::linestring linestring_type; + + linestring_type ls; + bg::read_wkt("LINESTRING(0 0,1 1)", ls); + + linestring_type res; + + bg::srs::spheroid spheroid(6378137.0, 6356752.3142451793); + bg::strategy::densify::geographic<> strategy(spheroid); + + boost::geometry::densify(ls, res, 50000.0, strategy); + + std::cout << "densified: " << boost::geometry::wkt(res) << std::endl; + + return 0; +} + +//] + +//[densify_strategy_output +/*` +Output: +[pre +densified: LINESTRING(0 0,0.249972 0.250011,0.499954 0.500017,0.749954 0.750013,1 1) +] +*/ +//] diff --git a/include/boost/geometry/algorithms/densify.hpp b/include/boost/geometry/algorithms/densify.hpp new file mode 100644 index 000000000..ed948b346 --- /dev/null +++ b/include/boost/geometry/algorithms/densify.hpp @@ -0,0 +1,425 @@ +// Boost.Geometry + +// Copyright (c) 2017-2018, Oracle and/or its affiliates. + +// Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle + +// Licensed under the Boost Software License version 1.0. +// http://www.boost.org/users/license.html + +#ifndef BOOST_GEOMETRY_ALGORITHMS_DENSIFY_HPP +#define BOOST_GEOMETRY_ALGORITHMS_DENSIFY_HPP + + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +#include +#include + +#include + + +namespace boost { namespace geometry +{ + + +#ifndef DOXYGEN_NO_DETAIL +namespace detail { namespace densify +{ + +template +struct push_back_policy +{ + typedef typename boost::range_value::type point_type; + + inline explicit push_back_policy(Range & rng) + : m_rng(rng) + {} + + inline void apply(point_type const& p) + { + range::push_back(m_rng, p); + } + +private: + Range & m_rng; +}; + +template +inline void convert_and_push_back(Range & range, Point const& p) +{ + typename boost::range_value::type p2; + geometry::detail::conversion::convert_point_to_point(p, p2); + range::push_back(range, p2); +} + +template +struct densify_range +{ + template + static inline void apply(FwdRng const& rng, MutRng & rng_out, + T const& len, Strategy const& strategy) + { + typedef typename boost::range_iterator::type iterator_t; + typedef typename boost::range_value::type point_t; + + iterator_t it = boost::begin(rng); + iterator_t end = boost::end(rng); + + if (it == end) // empty(rng) + { + return; + } + + push_back_policy policy(rng_out); + + iterator_t prev = it; + for ( ++it ; it != end ; prev = it++) + { + point_t const& p0 = *prev; + point_t const& p1 = *it; + + convert_and_push_back(rng_out, p0); + + strategy.apply(p0, p1, policy, len); + } + + if (BOOST_GEOMETRY_CONDITION(AppendLastPoint)) + { + convert_and_push_back(rng_out, *prev); // back(rng) + } + } +}; + +template // false, X +struct densify_ring +{ + template + static inline void apply(Geometry const& ring, GeometryOut & ring_out, + T const& len, Strategy const& strategy) + { + geometry::detail::densify::densify_range + ::apply(ring, ring_out, len, strategy); + + if (boost::size(ring) <= 1) + return; + + typedef typename point_type::type point_t; + point_t const& p0 = range::back(ring); + point_t const& p1 = range::front(ring); + + push_back_policy policy(ring_out); + + strategy.apply(p0, p1, policy, len); + + if (BOOST_GEOMETRY_CONDITION(IsClosed2)) + { + convert_and_push_back(ring_out, p1); + } + } +}; + +template <> +struct densify_ring + : densify_range +{}; + +template <> +struct densify_ring + : densify_range +{}; + + +}} // namespace detail::densify +#endif // DOXYGEN_NO_DETAIL + + +#ifndef DOXYGEN_NO_DISPATCH +namespace dispatch +{ + + +template +< + typename Geometry, + typename GeometryOut, + typename Tag1 = typename tag::type, + typename Tag2 = typename tag::type +> +struct densify + : not_implemented +{}; + +template +struct densify + : geometry::detail::densify::densify_range<> +{}; + +template +struct densify +{ + template + static void apply(Geometry const& mls, GeometryOut & mls_out, + T const& len, Strategy const& strategy) + { + std::size_t count = boost::size(mls); + range::resize(mls_out, count); + + for (std::size_t i = 0 ; i < count ; ++i) + { + geometry::detail::densify::densify_range<> + ::apply(range::at(mls, i), range::at(mls_out, i), + len, strategy); + } + } +}; + +template +struct densify + : geometry::detail::densify::densify_ring + < + geometry::closure::value != geometry::open, + geometry::closure::value != geometry::open + > +{}; + +template +struct densify +{ + template + static void apply(Geometry const& poly, GeometryOut & poly_out, + T const& len, Strategy const& strategy) + { + apply_ring(exterior_ring(poly), exterior_ring(poly_out), + len, strategy); + + std::size_t count = boost::size(interior_rings(poly)); + range::resize(interior_rings(poly_out), count); + + for (std::size_t i = 0 ; i < count ; ++i) + { + apply_ring(range::at(interior_rings(poly), i), + range::at(interior_rings(poly_out), i), + len, strategy); + } + } + + template + static void apply_ring(Ring const& ring, RingOut & ring_out, + T const& len, Strategy const& strategy) + { + densify + ::apply(ring, ring_out, len, strategy); + } +}; + +template +struct densify +{ + template + static void apply(Geometry const& mpoly, GeometryOut & mpoly_out, + T const& len, Strategy const& strategy) + { + std::size_t count = boost::size(mpoly); + range::resize(mpoly_out, count); + + for (std::size_t i = 0 ; i < count ; ++i) + { + apply_poly(range::at(mpoly, i), + range::at(mpoly_out, i), + len, strategy); + } + } + + template + static void apply_poly(Poly const& poly, PolyOut & poly_out, + T const& len, Strategy const& strategy) + { + densify:: + apply(poly, poly_out, len, strategy); + } +}; + + +} // namespace dispatch +#endif // DOXYGEN_NO_DISPATCH + + +namespace resolve_strategy +{ + +struct densify +{ + template + static inline void apply(Geometry const& geometry, + Geometry& out, + Distance const& max_distance, + Strategy const& strategy) + { + dispatch::densify + ::apply(geometry, out, max_distance, strategy); + } + + template + static inline void apply(Geometry const& geometry, + Geometry& out, + Distance const& max_distance, + default_strategy) + { + typedef typename strategy::densify::services::default_strategy + < + typename cs_tag::type + >::type strategy_type; + + /*BOOST_CONCEPT_ASSERT( + (concepts::DensifyStrategy) + );*/ + + apply(geometry, out, max_distance, strategy_type()); + } +}; + +} // namespace resolve_strategy + + +namespace resolve_variant { + +template +struct densify +{ + template + static inline void apply(Geometry const& geometry, + Geometry& out, + Distance const& max_distance, + Strategy const& strategy) + { + resolve_strategy::densify::apply(geometry, out, max_distance, strategy); + } +}; + +template +struct densify > +{ + template + struct visitor: boost::static_visitor + { + Distance const& m_max_distance; + Strategy const& m_strategy; + + visitor(Distance const& max_distance, Strategy const& strategy) + : m_max_distance(max_distance) + , m_strategy(strategy) + {} + + template + void operator()(Geometry const& geometry, Geometry& out) const + { + densify::apply(geometry, out, m_max_distance, m_strategy); + } + }; + + template + static inline void + apply(boost::variant const& geometry, + boost::variant& out, + Distance const& max_distance, + Strategy const& strategy) + { + boost::apply_visitor( + visitor(max_distance, strategy), + geometry, + out + ); + } +}; + +} // namespace resolve_variant + + +/*! +\brief Densify a geometry using a specified strategy +\ingroup densify +\tparam Geometry \tparam_geometry +\tparam Distance A numerical distance measure +\tparam Strategy A type fulfilling a DensifyStrategy concept +\param geometry Input geometry, to be densified +\param out Output geometry, densified version of the input geometry +\param max_distance Distance threshold (in units depending on strategy) +\param strategy Densify strategy to be used for densification + +\qbk{distinguish,with strategy} +\qbk{[include reference/algorithms/densify.qbk]} + +\qbk{ +[heading Available Strategies] +\* [link geometry.reference.strategies.strategy_densify_cartesian Cartesian] +\* [link geometry.reference.strategies.strategy_densify_spherical Spherical] +\* [link geometry.reference.strategies.strategy_densify_geographic Geographic] + +[heading Example] +[densify_strategy] +[densify_strategy_output] +} +*/ +template +inline void densify(Geometry const& geometry, + Geometry& out, + Distance const& max_distance, + Strategy const& strategy) +{ + concepts::check(); + + if (max_distance <= Distance(0)) + { + BOOST_THROW_EXCEPTION(geometry::invalid_input_exception()); + } + + geometry::clear(out); + + resolve_variant::densify + < + Geometry + >::apply(geometry, out, max_distance, strategy); +} + + +/*! +\brief Densify a geometry +\ingroup densify +\tparam Geometry \tparam_geometry +\tparam Distance A numerical distance measure +\param geometry Input geometry, to be densified +\param out Output geometry, densified version of the input geometry +\param max_distance Distance threshold (in units depending on coordinate system) + +\qbk{[include reference/algorithms/densify.qbk]} + +\qbk{ +[heading Example] +[densify] +[densify_output] +} +*/ +template +inline void densify(Geometry const& geometry, + Geometry& out, + Distance const& max_distance) +{ + densify(geometry, out, max_distance, default_strategy()); +} + + +}} // namespace boost::geometry + +#endif // BOOST_GEOMETRY_ALGORITHMS_DENSIFY_HPP diff --git a/include/boost/geometry/geometry.hpp b/include/boost/geometry/geometry.hpp index aa8070110..f81a78c02 100644 --- a/include/boost/geometry/geometry.hpp +++ b/include/boost/geometry/geometry.hpp @@ -4,8 +4,8 @@ // Copyright (c) 2008-2015 Bruno Lalande, Paris, France. // Copyright (c) 2009-2015 Mateusz Loskot, London, UK. -// This file was modified by Oracle on 2014, 2015, 2016. -// Modifications copyright (c) 2014-2016 Oracle and/or its affiliates. +// This file was modified by Oracle on 2014-2018. +// Modifications copyright (c) 2014-2018 Oracle and/or its affiliates. // Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle // Contributed and/or modified by Menelaos Karavelas, on behalf of Oracle @@ -61,6 +61,7 @@ #include #include #include +#include #include #include #include diff --git a/include/boost/geometry/strategies/cartesian/densify.hpp b/include/boost/geometry/strategies/cartesian/densify.hpp new file mode 100644 index 000000000..23651637b --- /dev/null +++ b/include/boost/geometry/strategies/cartesian/densify.hpp @@ -0,0 +1,133 @@ +// Boost.Geometry + +// Copyright (c) 2017-2018, Oracle and/or its affiliates. + +// Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle + +// Licensed under the Boost Software License version 1.0. +// http://www.boost.org/users/license.html + +#ifndef BOOST_GEOMETRY_STRATEGIES_CARTESIAN_DENSIFY_HPP +#define BOOST_GEOMETRY_STRATEGIES_CARTESIAN_DENSIFY_HPP + + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + + +namespace boost { namespace geometry +{ + +namespace strategy { namespace densify +{ + + +/*! +\brief Densification of cartesian segment. +\ingroup strategies +\tparam CalculationType \tparam_calculation + +\qbk{ +[heading See also] +[link geometry.reference.algorithms.densify.densify_4_with_strategy densify (with strategy)] +} + */ +template +< + typename CalculationType = void +> +class cartesian +{ +public: + template + static inline void apply(Point const& p0, Point const& p1, AssignPolicy & policy, T const& length_threshold) + { + typedef typename AssignPolicy::point_type out_point_t; + typedef typename select_most_precise + < + typename coordinate_type::type, + typename coordinate_type::type, + CalculationType + >::type calc_t; + + typedef model::point::value, cs::cartesian> calc_point_t; + + calc_point_t cp0, cp1; + geometry::detail::conversion::convert_point_to_point(p0, cp0); + geometry::detail::conversion::convert_point_to_point(p1, cp1); + + // dir01 = xy1 - xy0 + calc_point_t dir01 = cp1; + geometry::subtract_point(dir01, cp0); + calc_t const dot01 = geometry::dot_product(dir01, dir01); + calc_t const len = math::sqrt(dot01); + + BOOST_GEOMETRY_ASSERT(length_threshold > T(0)); + + signed_size_type n = signed_size_type(len / length_threshold); + if (n <= 0) + { + return; + } + + // NOTE: Normalization will not work for integral coordinates + // normalize + //geometry::divide_value(dir01, len); + + calc_t step = len / (n + 1); + + calc_t d = step; + for (signed_size_type i = 0 ; i < n ; ++i, d += step) + { + // pd = xy0 + d * dir01 + calc_point_t pd = dir01; + + // without normalization + geometry::multiply_value(pd, calc_t(i + 1)); + geometry::divide_value(pd, calc_t(n + 1)); + // with normalization + //geometry::multiply_value(pd, d); + + geometry::add_point(pd, cp0); + + // NOTE: Only needed if types calc_point_t and out_point_t are different + // otherwise pd could simply be passed into policy + out_point_t p; + assert_dimension_equal(); + geometry::detail::conversion::convert_point_to_point(pd, p); + + policy.apply(p); + } + } +}; + + +#ifndef DOXYGEN_NO_STRATEGY_SPECIALIZATIONS +namespace services +{ + +template <> +struct default_strategy +{ + typedef strategy::densify::cartesian<> type; +}; + + +} // namespace services +#endif // DOXYGEN_NO_STRATEGY_SPECIALIZATIONS + + +}} // namespace strategy::densify + + +}} // namespace boost::geometry + +#endif // BOOST_GEOMETRY_STRATEGIES_CARTESIAN_DENSIFY_HPP diff --git a/include/boost/geometry/strategies/densify.hpp b/include/boost/geometry/strategies/densify.hpp new file mode 100644 index 000000000..cde061af9 --- /dev/null +++ b/include/boost/geometry/strategies/densify.hpp @@ -0,0 +1,42 @@ +// Boost.Geometry + +// Copyright (c) 2017, Oracle and/or its affiliates. + +// Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle + +// Licensed under the Boost Software License version 1.0. +// http://www.boost.org/users/license.html + +#ifndef BOOST_GEOMETRY_STRATEGIES_DENSIFY_HPP +#define BOOST_GEOMETRY_STRATEGIES_DENSIFY_HPP + + +#include + + +namespace boost { namespace geometry +{ + +namespace strategy { namespace densify +{ + +namespace services +{ + +template +struct default_strategy +{ + BOOST_MPL_ASSERT_MSG + ( + false, NOT_IMPLEMENTED_FOR_THIS_CS + , (types) + ); +}; + +} // namespace services + +}} // namespace strategy::densify + +}} // namespace boost::geometry + +#endif // BOOST_GEOMETRY_STRATEGIES_DENSIFY_HPP diff --git a/include/boost/geometry/strategies/geographic/densify.hpp b/include/boost/geometry/strategies/geographic/densify.hpp new file mode 100644 index 000000000..a2908ea4a --- /dev/null +++ b/include/boost/geometry/strategies/geographic/densify.hpp @@ -0,0 +1,135 @@ +// Boost.Geometry + +// Copyright (c) 2017-2018, Oracle and/or its affiliates. + +// Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle + +// Licensed under the Boost Software License version 1.0. +// http://www.boost.org/users/license.html + +#ifndef BOOST_GEOMETRY_STRATEGIES_GEOGRAPHIC_DENSIFY_HPP +#define BOOST_GEOMETRY_STRATEGIES_GEOGRAPHIC_DENSIFY_HPP + + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + + +namespace boost { namespace geometry +{ + +namespace strategy { namespace densify +{ + + +/*! +\brief Densification of geographic segment. +\ingroup strategies +\tparam FormulaPolicy The geodesic formulas used internally. +\tparam Spheroid The spheroid model. +\tparam CalculationType \tparam_calculation + +\qbk{ +[heading See also] +[link geometry.reference.algorithms.densify.densify_4_with_strategy densify (with strategy)] +} + */ +template +< + typename FormulaPolicy = strategy::andoyer, + typename Spheroid = srs::spheroid, + typename CalculationType = void +> +class geographic +{ +public: + geographic() + : m_spheroid() + {} + + explicit geographic(Spheroid const& spheroid) + : m_spheroid(spheroid) + {} + + template + inline void apply(Point const& p0, Point const& p1, AssignPolicy & policy, T const& length_threshold) const + { + typedef typename AssignPolicy::point_type out_point_t; + typedef typename select_most_precise + < + typename coordinate_type::type, + typename coordinate_type::type, + CalculationType + >::type calc_t; + + typedef typename FormulaPolicy::template direct direct_t; + typedef typename FormulaPolicy::template inverse inverse_t; + + typename inverse_t::result_type + inv_r = inverse_t::apply(get_as_radian<0>(p0), get_as_radian<1>(p0), + get_as_radian<0>(p1), get_as_radian<1>(p1), + m_spheroid); + + BOOST_GEOMETRY_ASSERT(length_threshold > T(0)); + + signed_size_type n = signed_size_type(inv_r.distance / length_threshold); + if (n <= 0) + return; + + calc_t step = inv_r.distance / (n + 1); + + calc_t current = step; + for (signed_size_type i = 0 ; i < n ; ++i, current += step) + { + typename direct_t::result_type + dir_r = direct_t::apply(get_as_radian<0>(p0), get_as_radian<1>(p0), + current, inv_r.azimuth, + m_spheroid); + + out_point_t p; + set_from_radian<0>(p, dir_r.lon2); + set_from_radian<1>(p, dir_r.lat2); + geometry::detail::conversion::point_to_point + < + Point, out_point_t, + 2, dimension::value + >::apply(p0, p); + + policy.apply(p); + } + } + +private: + Spheroid m_spheroid; +}; + + +#ifndef DOXYGEN_NO_STRATEGY_SPECIALIZATIONS +namespace services +{ + +template <> +struct default_strategy +{ + typedef strategy::densify::geographic<> type; +}; + + +} // namespace services +#endif // DOXYGEN_NO_STRATEGY_SPECIALIZATIONS + + +}} // namespace strategy::densify + + +}} // namespace boost::geometry + +#endif // BOOST_GEOMETRY_ALGORITHMS_DENSIFY_HPP diff --git a/include/boost/geometry/strategies/spherical/densify.hpp b/include/boost/geometry/strategies/spherical/densify.hpp new file mode 100644 index 000000000..a9b58dd41 --- /dev/null +++ b/include/boost/geometry/strategies/spherical/densify.hpp @@ -0,0 +1,186 @@ +// Boost.Geometry + +// Copyright (c) 2017-2018, Oracle and/or its affiliates. + +// Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle + +// Licensed under the Boost Software License version 1.0. +// http://www.boost.org/users/license.html + +#ifndef BOOST_GEOMETRY_STRATEGIES_SPHERICAL_DENSIFY_HPP +#define BOOST_GEOMETRY_STRATEGIES_SPHERICAL_DENSIFY_HPP + + +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + + +namespace boost { namespace geometry +{ + +namespace strategy { namespace densify +{ + + +/*! +\brief Densification of spherical segment. +\ingroup strategies +\tparam Sphere The sphere model. +\tparam CalculationType \tparam_calculation + +\qbk{ +[heading See also] +[link geometry.reference.algorithms.densify.densify_4_with_strategy densify (with strategy)] +} + */ +template +< + typename Sphere = srs::sphere, + typename CalculationType = void +> +class spherical +{ +public: + // For consistency with area strategy the radius is set to 1 + spherical() + : m_radius(1.0) + {} + + explicit spherical(Sphere const& sphere) + : m_radius(geometry::get_radius<0>(sphere)) + {} + + // TODO: use enable_if/disable_if to distinguish between Sphere and Radius? + template + explicit spherical(OtherSphere const& sphere) + : m_radius(geometry::get_radius<0>(sphere)) + {} + + template + inline void apply(Point const& p0, Point const& p1, AssignPolicy & policy, T const& length_threshold) const + { + typedef typename AssignPolicy::point_type out_point_t; + typedef typename select_most_precise + < + typename coordinate_type::type, + typename coordinate_type::type, + CalculationType + >::type calc_t; + + calc_t const c0 = 0; + calc_t const c1 = 1; + calc_t const pi = math::pi(); + + typedef model::point point3d_t; + point3d_t const xyz0 = formula::sph_to_cart3d(p0); + point3d_t const xyz1 = formula::sph_to_cart3d(p1); + calc_t const dot01 = geometry::dot_product(xyz0, xyz1); + calc_t const angle01 = acos(dot01); + + BOOST_GEOMETRY_ASSERT(length_threshold > T(0)); + + signed_size_type n = signed_size_type(angle01 * m_radius / length_threshold); + if (n <= 0) + return; + + point3d_t axis; + if (! math::equals(angle01, pi)) + { + axis = geometry::cross_product(xyz0, xyz1); + geometry::detail::vec_normalize(axis); + } + else // antipodal + { + calc_t const half_pi = math::half_pi(); + calc_t const lat = geometry::get_as_radian<1>(p0); + + if (math::equals(lat, half_pi)) + { + // pointing east, segment lies on prime meridian, going south + axis = point3d_t(c0, c1, c0); + } + else if (math::equals(lat, -half_pi)) + { + // pointing west, segment lies on prime meridian, going north + axis = point3d_t(c0, -c1, c0); + } + else + { + // lon rotated west by pi/2 at equator + calc_t const lon = geometry::get_as_radian<0>(p0); + axis = point3d_t(sin(lon), -cos(lon), c0); + } + } + + calc_t step = angle01 / (n + 1); + + calc_t a = step; + for (signed_size_type i = 0 ; i < n ; ++i, a += step) + { + // Axis-Angle rotation + // see: https://en.wikipedia.org/wiki/Axis-angle_representation + calc_t const cos_a = cos(a); + calc_t const sin_a = sin(a); + // cos_a * v + point3d_t s1 = xyz0; + geometry::multiply_value(s1, cos_a); + // sin_a * (n x v) + point3d_t s2 = geometry::cross_product(axis, xyz0); + geometry::multiply_value(s2, sin_a); + // (1 - cos_a)(n.v) * n + point3d_t s3 = axis; + geometry::multiply_value(s3, (c1 - cos_a) * geometry::dot_product(axis, xyz0)); + // v_rot = cos_a * v + sin_a * (n x v) + (1 - cos_a)(n.v) * e + point3d_t v_rot = s1; + geometry::add_point(v_rot, s2); + geometry::add_point(v_rot, s3); + + out_point_t p = formula::cart3d_to_sph(v_rot); + geometry::detail::conversion::point_to_point + < + Point, out_point_t, + 2, dimension::value + >::apply(p0, p); + + policy.apply(p); + } + } + +private: + typename geometry::radius_type::type m_radius; +}; + + +#ifndef DOXYGEN_NO_STRATEGY_SPECIALIZATIONS +namespace services +{ + +template <> +struct default_strategy +{ + typedef strategy::densify::spherical<> type; +}; + + +} // namespace services +#endif // DOXYGEN_NO_STRATEGY_SPECIALIZATIONS + + +}} // namespace strategy::densify + + +}} // namespace boost::geometry + +#endif // BOOST_GEOMETRY_ALGORITHMS_DENSIFY_HPP diff --git a/include/boost/geometry/strategies/strategies.hpp b/include/boost/geometry/strategies/strategies.hpp index e7f3604ab..77076a092 100644 --- a/include/boost/geometry/strategies/strategies.hpp +++ b/include/boost/geometry/strategies/strategies.hpp @@ -4,8 +4,8 @@ // Copyright (c) 2008-2012 Bruno Lalande, Paris, France. // Copyright (c) 2009-2012 Mateusz Loskot, London, UK. -// This file was modified by Oracle on 2014-2017. -// Modifications copyright (c) 2014-2017 Oracle and/or its affiliates. +// This file was modified by Oracle on 2014-2018. +// Modifications copyright (c) 2014-2018 Oracle and/or its affiliates. // Contributed and/or modified by Vissarion Fysikopoulos, on behalf of Oracle // Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle @@ -30,6 +30,7 @@ #include #include #include +#include #include #include #include @@ -54,6 +55,7 @@ #include #include #include +#include #include #include #include @@ -70,6 +72,7 @@ #include #include +#include #include #include #include @@ -82,6 +85,7 @@ #include #include +#include #include #include #include diff --git a/test/algorithms/Jamfile.v2 b/test/algorithms/Jamfile.v2 index 2806c1a10..c0ae853b7 100644 --- a/test/algorithms/Jamfile.v2 +++ b/test/algorithms/Jamfile.v2 @@ -29,6 +29,7 @@ test-suite boost-geometry-algorithms [ run correct.cpp : : : : algorithms_correct ] [ run correct_multi.cpp : : : : algorithms_correct_multi ] [ run correct_closure.cpp : : : : algorithms_correct_closure ] + [ run densify.cpp : : : : algorithms_densify ] [ run for_each.cpp : : : : algorithms_for_each ] [ run for_each_multi.cpp : : : : algorithms_for_each_multi ] [ run is_convex.cpp : : : : algorithms_is_convex ] diff --git a/test/algorithms/densify.cpp b/test/algorithms/densify.cpp new file mode 100644 index 000000000..bf723bbc1 --- /dev/null +++ b/test/algorithms/densify.cpp @@ -0,0 +1,227 @@ +// Boost.Geometry +// Unit Test + +// Copyright (c) 2017-2018, Oracle and/or its affiliates. + +// Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle + +// Licensed under the Boost Software License version 1.0. +// http://www.boost.org/users/license.html + + +#include + +#include + +#include +#include +#include +#include + +#include + +#include +#include +#include +#include +#include +#include + +#include + + +struct check_lengths +{ + template + void operator()(G const& g, G const& o, S const& s) const + { + double d1 = bg::length(g, s); + double d2 = bg::length(o, s); + + BOOST_CHECK_CLOSE(d1, d2, 0.0001); + } +}; + +struct check_perimeters +{ + template + void operator()(G const& g, G const& o, S const& s) const + { + double d1 = bg::perimeter(g, s); + double d2 = bg::perimeter(o, s); + + BOOST_CHECK_CLOSE(d1, d2, 0.0001); + } +}; + +template +double inline shortest_length(G const& g, DistS const& dist_s) +{ + double min_len = (std::numeric_limits::max)(); + for (bg::segment_iterator it = bg::segments_begin(g); + it != bg::segments_end(g); ++it) + { + double len = bg::length(*it, dist_s); + min_len = (std::min)(min_len, len); + } + return min_len; +} + +template +double inline greatest_length(G const& o, DistS const& dist_s) +{ + double max_len = 0.0; + for (bg::segment_iterator it = bg::segments_begin(o); + it != bg::segments_end(o); ++it) + { + double len = bg::length(*it, dist_s); + max_len = (std::max)(max_len, len); + } + return max_len; +} + +template ::type> +struct cs_data +{}; + +template +struct cs_data +{ + bg::strategy::densify::cartesian<> compl_s; + bg::strategy::distance::pythagoras<> dist_s; +}; + +template +struct cs_data +{ + cs_data() + : model(6378137.0) + , compl_s(model) + , dist_s(6378137.0) + {} + + bg::srs::sphere model; + bg::strategy::densify::spherical<> compl_s; + bg::strategy::distance::haversine dist_s; +}; + +template +struct cs_data +{ + cs_data() + : model(6378137.0, 6356752.3142451793) + , compl_s(model) + , dist_s(model) + {} + + bg::srs::spheroid model; + bg::strategy::densify::geographic<> compl_s; + bg::strategy::distance::geographic<> dist_s; +}; + +template +inline void check_result(G const& g, G const& o, double max_distance, + DistS const& dist_s, Check const& check) +{ + // geometry was indeed densified + std::size_t g_count = bg::num_points(g); + std::size_t o_count = bg::num_points(o); + BOOST_CHECK(g_count < o_count); + + // all segments have lengths smaller or equal to max_distance + double gr_len = greatest_length(o, dist_s); + // NOTE: Currently geographic strategies can generate segments that have + // lengths slightly greater than max_distance. In order to change + // this the generation of new points should e.g. be recursive with + // stop condition comparing the current distance calculated by + // inverse strategy. + // NOTE: Closeness value tweaked for Andoyer + bool is_close = (gr_len - max_distance) / (std::max)(gr_len, max_distance) < 0.0001; + BOOST_CHECK(gr_len <= max_distance || is_close); + + // the overall length or perimeter didn't change + check(g, o, dist_s); +} + +template +inline void test_geometry(std::string const& wkt, Check const& check) +{ + cs_data d; + + G g; + bg::read_wkt(wkt, g); + + { + bg::default_strategy def_s; + double max_distance = shortest_length(g, def_s) / 3.0; + + G o; + bg::densify(g, o, max_distance); + + check_result(g, o, max_distance, def_s, check); + } + + { + double max_distance = shortest_length(g, d.dist_s) / 3.0; + + G o; + bg::densify(g, o, max_distance, d.compl_s); + + check_result(g, o, max_distance, d.dist_s, check); + } +} + +template +inline void test_linear(std::string const& wkt) +{ + test_geometry(wkt, check_lengths()); +} + +template +inline void test_areal(std::string const& wkt) +{ + test_geometry(wkt, check_perimeters()); +} + +template +void test_all() +{ + typedef bg::model::linestring

ls_t; + typedef bg::model::multi_linestring mls_t; + + typedef bg::model::ring

ring_t; + typedef bg::model::polygon

poly_t; + typedef bg::model::multi_polygon mpoly_t; + + typedef bg::model::ring oring_t; + typedef bg::model::polygon opoly_t; + typedef bg::model::multi_polygon ompoly_t; + + test_linear("LINESTRING(4 -4, 4 -1)"); + test_linear("LINESTRING(4 4, 4 1)"); + test_linear("LINESTRING(0 0, 180 0)"); + test_linear("LINESTRING(1 1, -179 -1)"); + + test_linear("LINESTRING(1 1, 2 2, 4 2)"); + test_linear("MULTILINESTRING((1 1, 2 2),(2 2, 4 2))"); + + test_areal("POLYGON((1 1, 1 2, 2 2, 1 1))"); + test_areal("POLYGON((1 1, 1 4, 4 4, 4 1, 1 1),(1 1, 2 2, 2 3, 1 1))"); + test_areal("MULTIPOLYGON(((1 1, 1 4, 4 4, 4 1, 1 1),(1 1, 2 2, 2 3, 1 1)),((4 4, 5 5, 5 4, 4 4)))"); + + test_areal("POLYGON((1 1, 1 2, 2 2))"); + test_areal("POLYGON((1 1, 1 4, 4 4, 4 1),(1 1, 2 2, 2 3))"); + test_areal("MULTIPOLYGON(((1 1, 1 4, 4 4, 4 1),(1 1, 2 2, 2 3)),((4 4, 5 5, 5 4)))"); + + test_areal("POLYGON((0 0,0 40,40 40,40 0,0 0))"); + test_areal("POLYGON((0 0,0 40,40 40,40 0))"); +} + +int test_main(int, char* []) +{ + test_all< bg::model::point >(); + test_all< bg::model::point > >(); + test_all< bg::model::point > >(); + + return 0; +}