diff options
Diffstat (limited to 'boost/geometry/algorithms')
23 files changed, 1194 insertions, 89 deletions
diff --git a/boost/geometry/algorithms/area.hpp b/boost/geometry/algorithms/area.hpp index 919ab66fe1..8193200ab9 100644 --- a/boost/geometry/algorithms/area.hpp +++ b/boost/geometry/algorithms/area.hpp @@ -92,7 +92,7 @@ struct ring_area // An open ring has at least three points, // A closed ring has at least four points, // if not, there is no (zero) area - if (boost::size(ring) + if (int(boost::size(ring)) < core_detail::closure::minimum_ring_size<Closure>::value) { return type(); diff --git a/boost/geometry/algorithms/detail/for_each_range.hpp b/boost/geometry/algorithms/detail/for_each_range.hpp index b10017c209..7cb01fa9b4 100644 --- a/boost/geometry/algorithms/detail/for_each_range.hpp +++ b/boost/geometry/algorithms/detail/for_each_range.hpp @@ -22,6 +22,7 @@ #include <boost/geometry/core/tag_cast.hpp> #include <boost/geometry/util/add_const_if_c.hpp> +#include <boost/geometry/views/box_view.hpp> namespace boost { namespace geometry diff --git a/boost/geometry/algorithms/detail/get_left_turns.hpp b/boost/geometry/algorithms/detail/get_left_turns.hpp new file mode 100644 index 0000000000..62f0f7f0f4 --- /dev/null +++ b/boost/geometry/algorithms/detail/get_left_turns.hpp @@ -0,0 +1,367 @@ +// Boost.Geometry (aka GGL, Generic Geometry Library) + +// Copyright (c) 2012 Barend Gehrels, Amsterdam, the Netherlands. + +// 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) + +#ifndef BOOST_GEOMETRY_ALGORITHMS_DETAIL_GET_LEFT_TURNS_HPP +#define BOOST_GEOMETRY_ALGORITHMS_DETAIL_GET_LEFT_TURNS_HPP + +#include <boost/geometry/iterators/ever_circling_iterator.hpp> + +namespace boost { namespace geometry +{ + + +#ifndef DOXYGEN_NO_DETAIL +namespace detail +{ + +// TODO: move this to /util/ +template <typename T> +static inline std::pair<T, T> ordered_pair(T const& first, T const& second) +{ + return first < second ? std::make_pair(first, second) : std::make_pair(second, first); +} + +template <typename AngleInfo> +inline void debug_left_turn(AngleInfo const& ai, AngleInfo const& previous) +{ +#ifdef BOOST_GEOMETRY_DEBUG_BUFFER_OCCUPATION + std::cout << "Angle: " << (ai.incoming ? "i" : "o") + << " " << si(ai.seg_id) + << " " << (math::r2d * (ai.angle) ) + << " turn: " << ai.turn_index << "[" << ai.operation_index << "]" + ; + + if (ai.turn_index != previous.turn_index + || ai.operation_index != previous.operation_index) + { + std::cout << " diff: " << math::r2d * math::abs(previous.angle - ai.angle); + } + std::cout << std::endl; +#endif +} + +template <typename AngleInfo> +inline void debug_left_turn(std::string const& caption, AngleInfo const& ai, AngleInfo const& previous, + int code = -99, int code2 = -99, int code3 = -99, int code4 = -99) +{ +#ifdef BOOST_GEOMETRY_DEBUG_BUFFER_OCCUPATION + std::cout << " " << caption + << " turn: " << ai.turn_index << "[" << ai.operation_index << "]" + << " " << si(ai.seg_id) + << " " << (ai.incoming ? "i" : "o") + << " " << (math::r2d * (ai.angle) ) + << " turn: " << previous.turn_index << "[" << previous.operation_index << "]" + << " " << si(previous.seg_id) + << " " << (previous.incoming ? "i" : "o") + << " " << (math::r2d * (previous.angle) ) + ; + + if (code != -99) + { + std::cout << " code: " << code << " , " << code2 << " , " << code3 << " , " << code4; + } + std::cout << std::endl; +#endif +} + + +template <typename Operation> +inline bool include_operation(Operation const& op, + segment_identifier const& outgoing_seg_id, + segment_identifier const& incoming_seg_id) +{ + return op.seg_id == outgoing_seg_id + && op.other_id == incoming_seg_id + && (op.operation == detail::overlay::operation_union + ||op.operation == detail::overlay::operation_continue) + ; +} + +template <typename Turn> +inline bool process_include(segment_identifier const& outgoing_seg_id, segment_identifier const& incoming_seg_id, + int turn_index, Turn& turn, + std::set<int>& keep_indices, int priority) +{ + bool result = false; + for (int i = 0; i < 2; i++) + { + if (include_operation(turn.operations[i], outgoing_seg_id, incoming_seg_id)) + { + turn.operations[i].include_in_occupation_map = true; + if (priority > turn.priority) + { + turn.priority = priority; + } + keep_indices.insert(turn_index); + result = true; + } + } + return result; +} + +template <typename AngleInfo, typename Turns, typename TurnSegmentIndices> +inline bool include_left_turn_of_all( + AngleInfo const& outgoing, AngleInfo const& incoming, + Turns& turns, TurnSegmentIndices const& turn_segment_indices, + std::set<int>& keep_indices, int priority) +{ + segment_identifier const& outgoing_seg_id = turns[outgoing.turn_index].operations[outgoing.operation_index].seg_id; + segment_identifier const& incoming_seg_id = turns[incoming.turn_index].operations[incoming.operation_index].seg_id; + + if (outgoing.turn_index == incoming.turn_index) + { + return process_include(outgoing_seg_id, incoming_seg_id, outgoing.turn_index, turns[outgoing.turn_index], keep_indices, priority); + } + + bool result = false; + std::pair<segment_identifier, segment_identifier> pair = ordered_pair(outgoing_seg_id, incoming_seg_id); + auto it = turn_segment_indices.find(pair); + if (it != turn_segment_indices.end()) + { + for (auto sit = it->second.begin(); sit != it->second.end(); ++sit) + { + if (process_include(outgoing_seg_id, incoming_seg_id, *sit, turns[*sit], keep_indices, priority)) + { + result = true; + } + } + } + return result; +} + +template <std::size_t Index, typename Turn> +inline bool corresponds(Turn const& turn, segment_identifier const& seg_id) +{ + return turn.operations[Index].operation == detail::overlay::operation_union + && turn.operations[Index].seg_id == seg_id; +} + + +template <typename Turns, typename TurnSegmentIndices> +inline bool prefer_by_other(Turns const& turns, + TurnSegmentIndices const& turn_segment_indices, + std::set<int>& indices) +{ + std::map<segment_identifier, int> map; + for (auto sit = indices.begin(); sit != indices.end(); ++sit) + { + map[turns[*sit].operations[0].seg_id]++; + map[turns[*sit].operations[1].seg_id]++; + } + + std::set<segment_identifier> segment_occuring_once; + for (auto mit = map.begin(); mit != map.end(); ++mit) + { + if (mit->second == 1) + { + segment_occuring_once.insert(mit->first); + } +#ifdef BOOST_GEOMETRY_DEBUG_BUFFER_PREFER + std::cout << si(mit->first) << " " << mit->second << std::endl; +#endif + } + + if (segment_occuring_once.size() == 2) + { + // Try to find within all turns a turn with these two segments + + std::set<segment_identifier>::const_iterator soo_it = segment_occuring_once.begin(); + segment_identifier front = *soo_it; + soo_it++; + segment_identifier back = *soo_it; + + std::pair<segment_identifier, segment_identifier> pair = ordered_pair(front, back); + auto it = turn_segment_indices.find(pair); + if (it != turn_segment_indices.end()) + { + // debug_turns_by_indices("Found", it->second); + // Check which is the union/continue + segment_identifier good; + for (auto sit = it->second.begin(); sit != it->second.end(); ++sit) + { + if (turns[*sit].operations[0].operation == detail::overlay::operation_union) + { + good = turns[*sit].operations[0].seg_id; + } + else if (turns[*sit].operations[1].operation == detail::overlay::operation_union) + { + good = turns[*sit].operations[1].seg_id; + } + } + +#ifdef BOOST_GEOMETRY_DEBUG_BUFFER_PREFER + std::cout << "Good: " << si(good) << std::endl; +#endif + + // Find in indexes-to-keep this segment with the union. Discard the other one + std::set<int> ok_indices; + for (auto sit = indices.begin(); sit != indices.end(); ++sit) + { + if (corresponds<0>(turns[*sit], good) || corresponds<1>(turns[*sit], good)) + { + ok_indices.insert(*sit); + } + } + if (ok_indices.size() > 0 && ok_indices.size() < indices.size()) + { + indices = ok_indices; + std::cout << "^"; + return true; + } + } + } + return false; +} + +template <typename Turns> +inline void prefer_by_priority(Turns const& turns, std::set<int>& indices) +{ + // Find max prio + int min_prio = 1024, max_prio = 0; + for (auto sit = indices.begin(); sit != indices.end(); ++sit) + { + if (turns[*sit].priority > max_prio) + { + max_prio = turns[*sit].priority; + } + if (turns[*sit].priority < min_prio) + { + min_prio = turns[*sit].priority; + } + } + + if (min_prio == max_prio) + { + return; + } + + // Only keep indices with high prio + std::set<int> ok_indices; + for (auto sit = indices.begin(); sit != indices.end(); ++sit) + { + if (turns[*sit].priority >= max_prio) + { + ok_indices.insert(*sit); + } + } + if (ok_indices.size() > 0 && ok_indices.size() < indices.size()) + { + indices = ok_indices; + std::cout << "%"; + } +} + +template <typename AngleInfo, typename Angles, typename Turns, typename TurnSegmentIndices> +inline void calculate_left_turns(Angles const& angles, + Turns& turns, TurnSegmentIndices const& turn_segment_indices, + std::set<int>& keep_indices) +{ + bool debug_indicate_size = false; + + typedef typename strategy::side::services::default_strategy + < + typename cs_tag<typename AngleInfo::point_type>::type + >::type side_strategy; + + std::size_t i = 0; + std::size_t n = boost::size(angles); + + typedef geometry::ever_circling_range_iterator<Angles const> circling_iterator; + circling_iterator cit(angles); + + debug_left_turn(*cit, *cit); + for(circling_iterator prev = cit++; i < n; prev = cit++, i++) + { + debug_left_turn(*cit, *prev); + + bool const include = ! geometry::math::equals(prev->angle, cit->angle) + && ! prev->incoming + && cit->incoming; + + if (include) + { + // Go back to possibly include more same left outgoing angles: + // Because we check on side too we can take a large "epsilon" + circling_iterator back = prev - 1; + + typename AngleInfo::angle_type eps = 0.00001; + int b = 1; + for(std::size_t d = 0; + math::abs(prev->angle - back->angle) < eps + && ! back->incoming + && d < n; + d++) + { + --back; + ++b; + } + + // Same but forward to possibly include more incoming angles + int f = 1; + circling_iterator forward = cit + 1; + for(std::size_t d = 0; + math::abs(cit->angle - forward->angle) < eps + && forward->incoming + && d < n; + d++) + { + ++forward; + ++f; + } + +#ifdef BOOST_GEOMETRY_DEBUG_BUFFER_OCCUPATION + std::cout << "HANDLE " << b << "/" << f << " ANGLES" << std::endl; +#endif + for(circling_iterator bit = prev; bit != back; --bit) + { + int code = side_strategy::apply(bit->direction_point, prev->intersection_point, prev->direction_point); + int code2 = side_strategy::apply(prev->direction_point, bit->intersection_point, bit->direction_point); + for(circling_iterator fit = cit; fit != forward; ++fit) + { + int code3 = side_strategy::apply(fit->direction_point, cit->intersection_point, cit->direction_point); + int code4 = side_strategy::apply(cit->direction_point, fit->intersection_point, fit->direction_point); + + int priority = 2; + if (code == -1 && code2 == 1) + { + // This segment is lying right of the other one. + // Cannot ignore it (because of robustness, see a.o. #rt_p21 from unit test). + // But if we find more we can prefer later by priority + // (a.o. necessary for #rt_o2 from unit test) + priority = 1; + } + + bool included = include_left_turn_of_all(*bit, *fit, turns, turn_segment_indices, keep_indices, priority); + debug_left_turn(included ? "KEEP" : "SKIP", *fit, *bit, code, code2, code3, code4); + } + } + } + } + + if (debug_indicate_size) + { + std::cout << " size=" << keep_indices.size(); + } + + if (keep_indices.size() >= 2) + { + prefer_by_other(turns, turn_segment_indices, keep_indices); + } + if (keep_indices.size() >= 2) + { + prefer_by_priority(turns, keep_indices); + } +} + +} // namespace detail +#endif // DOXYGEN_NO_DETAIL + + +}} // namespace boost::geometry + +#endif // BOOST_GEOMETRY_ALGORITHMS_DETAIL_GET_LEFT_TURNS_HPP diff --git a/boost/geometry/algorithms/detail/occupation_info.hpp b/boost/geometry/algorithms/detail/occupation_info.hpp new file mode 100644 index 0000000000..e147ba12d8 --- /dev/null +++ b/boost/geometry/algorithms/detail/occupation_info.hpp @@ -0,0 +1,329 @@ +// Boost.Geometry (aka GGL, Generic Geometry Library) + +// Copyright (c) 2012 Barend Gehrels, Amsterdam, the Netherlands. + +// 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) + +#ifndef BOOST_GEOMETRY_ALGORITHMS_DETAIL_OCCUPATION_INFO_HPP +#define BOOST_GEOMETRY_ALGORITHMS_DETAIL_OCCUPATION_INFO_HPP + +#if ! defined(NDEBUG) + #define BOOST_GEOMETRY_DEBUG_BUFFER_OCCUPATION +#endif + +#include <algorithm> +#include <boost/range.hpp> + +#include <boost/geometry/core/coordinate_type.hpp> +#include <boost/geometry/core/point_type.hpp> + +#include <boost/geometry/algorithms/equals.hpp> +#include <boost/geometry/iterators/closing_iterator.hpp> + +#include <boost/geometry/algorithms/detail/get_left_turns.hpp> + + +namespace boost { namespace geometry +{ + + +#ifndef DOXYGEN_NO_DETAIL +namespace detail +{ + +template <typename P> +class relaxed_less +{ + typedef typename geometry::coordinate_type<P>::type coordinate_type; + + coordinate_type epsilon; + +public : + + inline relaxed_less() + { + // TODO: adapt for ttmath, and maybe build the map in another way + // (e.g. exact constellations of segment-id's), maybe adaptive. + epsilon = std::numeric_limits<double>::epsilon() * 100.0; + } + + inline bool operator()(P const& a, P const& b) const + { + coordinate_type const dx = math::abs(geometry::get<0>(a) - geometry::get<0>(b)); + coordinate_type const dy = math::abs(geometry::get<1>(a) - geometry::get<1>(b)); + + + if (dx < epsilon && dy < epsilon) + { + return false; + } + if (dx < epsilon) + { + return geometry::get<1>(a) < geometry::get<1>(b); + } + + return geometry::get<0>(a) < geometry::get<0>(b); + } + + inline bool equals(P const& a, P const& b) const + { + typedef typename geometry::coordinate_type<P>::type coordinate_type; + + coordinate_type const dx = math::abs(geometry::get<0>(a) - geometry::get<0>(b)); + coordinate_type const dy = math::abs(geometry::get<1>(a) - geometry::get<1>(b)); + + return dx < epsilon && dy < epsilon; + }; +}; + + +template <typename T, typename P1, typename P2> +inline T calculate_angle(P1 const& from_point, P2 const& to_point) +{ + typedef P1 vector_type; + vector_type v = from_point; + geometry::subtract_point(v, to_point); + return atan2(geometry::get<1>(v), geometry::get<0>(v)); +} + +template <typename Iterator, typename Vector> +inline Iterator advance_circular(Iterator it, Vector const& vector, segment_identifier& seg_id, bool forward = true) +{ + int const increment = forward ? 1 : -1; + if (it == boost::begin(vector) && increment < 0) + { + it = boost::end(vector); + seg_id.segment_index = boost::size(vector); + } + it += increment; + seg_id.segment_index += increment; + if (it == boost::end(vector)) + { + seg_id.segment_index = 0; + it = boost::begin(vector); + } + return it; +} + +template <typename Point, typename T> +struct angle_info +{ + typedef T angle_type; + typedef Point point_type; + + segment_identifier seg_id; + int turn_index; + int operation_index; + Point intersection_point; + Point direction_point; + T angle; + bool incoming; +}; + +template <typename AngleInfo> +class occupation_info +{ + typedef std::vector<AngleInfo> collection_type; + + struct angle_sort + { + inline bool operator()(AngleInfo const& left, AngleInfo const& right) const + { + // In this case we can compare even double using equals + // return geometry::math::equals(left.angle, right.angle) + return left.angle == right.angle + ? int(left.incoming) < int(right.incoming) + : left.angle < right.angle + ; + } + }; + +public : + collection_type angles; +private : + bool m_occupied; + bool m_calculated; + + inline bool is_occupied() + { + if (boost::size(angles) <= 1) + { + return false; + } + + std::sort(angles.begin(), angles.end(), angle_sort()); + + typedef geometry::closing_iterator<collection_type const> closing_iterator; + closing_iterator vit(angles); + closing_iterator end(angles, true); + + closing_iterator prev = vit++; + for( ; vit != end; prev = vit++) + { + if (! geometry::math::equals(prev->angle, vit->angle) + && ! prev->incoming + && vit->incoming) + { + return false; + } + } + return true; + } + +public : + inline occupation_info() + : m_occupied(false) + , m_calculated(false) + {} + + template <typename PointC, typename Point1, typename Point2> + inline void add(PointC const& map_point, Point1 const& direction_point, Point2 const& intersection_point, + int turn_index, int operation_index, + segment_identifier const& seg_id, bool incoming) + { + //std::cout << "-> adding angle " << geometry::wkt(direction_point) << " .. " << geometry::wkt(intersection_point) << " " << int(incoming) << std::endl; + if (geometry::equals(direction_point, intersection_point)) + { + //std::cout << "EQUAL! Skipping" << std::endl; + return; + } + + AngleInfo info; + info.incoming = incoming; + info.angle = calculate_angle<typename AngleInfo::angle_type>(direction_point, map_point); + info.seg_id = seg_id; + info.turn_index = turn_index; + info.operation_index = operation_index; + info.intersection_point = intersection_point; + info.direction_point = direction_point; + angles.push_back(info); + + m_calculated = false; + } + + inline bool occupied() + { + if (! m_calculated) + { + m_occupied = is_occupied(); + m_calculated = true; + } + return m_occupied; + } + + template <typename Turns, typename TurnSegmentIndices> + inline void get_left_turns( + Turns& turns, TurnSegmentIndices const& turn_segment_indices, + std::set<int>& keep_indices) + { + std::sort(angles.begin(), angles.end(), angle_sort()); + calculate_left_turns<AngleInfo>(angles, turns, turn_segment_indices, keep_indices); + } +}; + + +template <typename Point, typename Ring, typename Info> +inline void add_incoming_and_outgoing_angles(Point const& map_point, Point const& intersection_point, + Ring const& ring, + int turn_index, + int operation_index, + segment_identifier seg_id, + Info& info) +{ + typedef typename boost::range_iterator + < + Ring const + >::type iterator_type; + + int const n = boost::size(ring); + if (seg_id.segment_index >= n || seg_id.segment_index < 0) + { + return; + } + + segment_identifier real_seg_id = seg_id; + iterator_type it = boost::begin(ring) + seg_id.segment_index; + + // TODO: if we use turn-info ("to", "middle"), we know if to advance without resorting to equals + relaxed_less<Point> comparator; + + if (comparator.equals(intersection_point, *it)) + { + // It should be equal only once. But otherwise we skip it (in "add") + it = advance_circular(it, ring, seg_id, false); + } + + info.add(map_point, *it, intersection_point, turn_index, operation_index, real_seg_id, true); + + if (comparator.equals(intersection_point, *it)) + { + it = advance_circular(it, ring, real_seg_id); + } + else + { + // Don't upgrade the ID + it = advance_circular(it, ring, seg_id); + } + for (int defensive_check = 0; + comparator.equals(intersection_point, *it) && defensive_check < n; + defensive_check++) + { + it = advance_circular(it, ring, real_seg_id); + } + + info.add(map_point, *it, intersection_point, turn_index, operation_index, real_seg_id, false); +} + + +// Map in two senses of the word: it is a std::map where the key is a point. +// Per point an "occupation_info" record is kept +// Used for the buffer (but will also be used for intersections/unions having complex self-tangencies) +template <typename Point, typename OccupationInfo> +class occupation_map +{ +public : + typedef std::map<Point, OccupationInfo, relaxed_less<Point> > map_type; + + map_type map; + std::set<int> turn_indices; + + inline OccupationInfo& find_or_insert(Point const& point, Point& mapped_point) + { + typename map_type::iterator it = map.find(point); + if (it == boost::end(map)) + { + std::pair<typename map_type::iterator, bool> pair + = map.insert(std::make_pair(point, OccupationInfo())); + it = pair.first; + } + mapped_point = it->first; + return it->second; + } + + inline bool contains(Point const& point) const + { + typename map_type::const_iterator it = map.find(point); + return it != boost::end(map); + } + + inline bool contains_turn_index(int index) const + { + return turn_indices.count(index) > 0; + } + + inline void insert_turn_index(int index) + { + turn_indices.insert(index); + } +}; + + +} // namespace detail +#endif // DOXYGEN_NO_DETAIL + + +}} // namespace boost::geometry + +#endif // BOOST_GEOMETRY_ALGORITHMS_DETAIL_OCCUPATION_INFO_HPP diff --git a/boost/geometry/algorithms/detail/overlay/add_rings.hpp b/boost/geometry/algorithms/detail/overlay/add_rings.hpp index eb3e60e483..74595fedd0 100644 --- a/boost/geometry/algorithms/detail/overlay/add_rings.hpp +++ b/boost/geometry/algorithms/detail/overlay/add_rings.hpp @@ -9,6 +9,8 @@ #ifndef BOOST_GEOMETRY_ALGORITHMS_DETAIL_OVERLAY_ADD_RINGS_HPP #define BOOST_GEOMETRY_ALGORITHMS_DETAIL_OVERLAY_ADD_RINGS_HPP +#include <boost/geometry/core/closure.hpp> +#include <boost/geometry/algorithms/area.hpp> #include <boost/geometry/algorithms/detail/overlay/convert_ring.hpp> #include <boost/geometry/algorithms/detail/overlay/get_ring.hpp> @@ -73,6 +75,21 @@ inline OutputIterator add_rings(SelectionMap const& map, OutputIterator out) { typedef typename SelectionMap::const_iterator iterator; + typedef typename SelectionMap::mapped_type property_type; + typedef typename property_type::area_type area_type; + + area_type const zero = 0; + std::size_t const min_num_points = core_detail::closure::minimum_ring_size + < + geometry::closure + < + typename boost::range_value + < + RingCollection const + >::type + >::value + >::value; + for (iterator it = boost::begin(map); it != boost::end(map); @@ -99,7 +116,16 @@ inline OutputIterator add_rings(SelectionMap const& map, *child_it, mit->second.reversed, true); } } - *out++ = result; + + // Only add rings if they satisfy minimal requirements. + // This cannot be done earlier (during traversal), not + // everything is figured out yet (sum of positive/negative rings) + // TODO: individual rings can still contain less than 3 points. + if (geometry::num_points(result) >= min_num_points + && math::larger(geometry::area(result), zero)) + { + *out++ = result; + } } } return out; diff --git a/boost/geometry/algorithms/detail/overlay/assign_parents.hpp b/boost/geometry/algorithms/detail/overlay/assign_parents.hpp index 133530563e..5063f49eb4 100644 --- a/boost/geometry/algorithms/detail/overlay/assign_parents.hpp +++ b/boost/geometry/algorithms/detail/overlay/assign_parents.hpp @@ -130,7 +130,7 @@ struct assign_visitor return; } - if (outer.real_area > 0) + if (math::larger(outer.real_area, 0)) { if (inner.real_area < 0 || m_check_for_orientation) { @@ -317,13 +317,14 @@ template > inline void assign_parents(Geometry const& geometry, RingCollection const& collection, - RingMap& ring_map) + RingMap& ring_map, + bool check_for_orientation) { // Call it with an empty geometry // (ring_map should be empty for source_id==1) Geometry empty; - assign_parents(geometry, empty, collection, ring_map, true); + assign_parents(geometry, empty, collection, ring_map, check_for_orientation); } diff --git a/boost/geometry/algorithms/detail/overlay/calculate_distance_policy.hpp b/boost/geometry/algorithms/detail/overlay/calculate_distance_policy.hpp index 3e6a8897f5..2003d2350d 100644 --- a/boost/geometry/algorithms/detail/overlay/calculate_distance_policy.hpp +++ b/boost/geometry/algorithms/detail/overlay/calculate_distance_policy.hpp @@ -32,9 +32,18 @@ struct calculate_distance_policy { static bool const include_no_turn = false; static bool const include_degenerate = false; - - template <typename Point1, typename Point2, typename Info> - static inline void apply(Info& info, Point1 const& p1, Point2 const& p2) + static bool const include_opposite = false; + + template + < + typename Info, + typename Point1, + typename Point2, + typename IntersectionInfo, + typename DirInfo + > + static inline void apply(Info& info, Point1 const& p1, Point2 const& p2, + IntersectionInfo const&, DirInfo const&) { info.operations[0].enriched.distance = geometry::comparable_distance(info.point, p1); diff --git a/boost/geometry/algorithms/detail/overlay/copy_segment_point.hpp b/boost/geometry/algorithms/detail/overlay/copy_segment_point.hpp index 379444428f..5e18d0453a 100644 --- a/boost/geometry/algorithms/detail/overlay/copy_segment_point.hpp +++ b/boost/geometry/algorithms/detail/overlay/copy_segment_point.hpp @@ -55,14 +55,14 @@ struct copy_segment_point_range if (second) { index++; - if (index >= boost::size(range)) + if (index >= int(boost::size(range))) { index = 0; } } // Exception? - if (index >= boost::size(range)) + if (index >= int(boost::size(range))) { return false; } diff --git a/boost/geometry/algorithms/detail/overlay/copy_segments.hpp b/boost/geometry/algorithms/detail/overlay/copy_segments.hpp index b0183a3ac2..805f3923e3 100644 --- a/boost/geometry/algorithms/detail/overlay/copy_segments.hpp +++ b/boost/geometry/algorithms/detail/overlay/copy_segments.hpp @@ -78,7 +78,7 @@ struct copy_segments_ring int const from_index = seg_id.segment_index + 1; // Sanity check - BOOST_ASSERT(from_index < boost::size(view)); + BOOST_ASSERT(from_index < int(boost::size(view))); ec_iterator it(boost::begin(view), boost::end(view), boost::begin(view) + from_index); @@ -89,7 +89,7 @@ struct copy_segments_ring typedef typename boost::range_difference<Ring>::type size_type; size_type const count = from_index <= to_index ? to_index - from_index + 1 - : boost::size(view) - from_index + to_index + 1; + : int(boost::size(view)) - from_index + to_index + 1; for (size_type i = 0; i < count; ++i, ++it) { @@ -117,7 +117,7 @@ struct copy_segments_linestring int const from_index = seg_id.segment_index + 1; // Sanity check - if (from_index > to_index || from_index < 0 || to_index >= boost::size(ls)) + if (from_index > to_index || from_index < 0 || to_index >= int(boost::size(ls))) { return; } diff --git a/boost/geometry/algorithms/detail/overlay/debug_turn_info.hpp b/boost/geometry/algorithms/detail/overlay/debug_turn_info.hpp index 0bec816946..0cc34255ca 100644 --- a/boost/geometry/algorithms/detail/overlay/debug_turn_info.hpp +++ b/boost/geometry/algorithms/detail/overlay/debug_turn_info.hpp @@ -28,6 +28,7 @@ inline char method_char(detail::overlay::method_type const& method) case method_touch_interior : return 'm'; case method_collinear : return 'c'; case method_equal : return 'e'; + case method_error : return '!'; default : return '?'; } } @@ -42,6 +43,7 @@ inline char operation_char(detail::overlay::operation_type const& operation) case operation_intersection : return 'i'; case operation_blocked : return 'x'; case operation_continue : return 'c'; + case operation_opposite : return 'o'; default : return '?'; } } diff --git a/boost/geometry/algorithms/detail/overlay/follow.hpp b/boost/geometry/algorithms/detail/overlay/follow.hpp index 087092a5f2..b110cc9602 100644 --- a/boost/geometry/algorithms/detail/overlay/follow.hpp +++ b/boost/geometry/algorithms/detail/overlay/follow.hpp @@ -279,12 +279,12 @@ class follow { // In case of turn point at the same location, we want to have continue/blocked LAST // because that should be followed (intersection) or skipped (difference). - // By chance the enumeration is ordered like that but we keep the safe way here. inline int operation_order(Turn const& turn) const { operation_type const& operation = turn.operations[0].operation; switch(operation) { + case operation_opposite : return 0; case operation_none : return 0; case operation_union : return 1; case operation_intersection : return 2; diff --git a/boost/geometry/algorithms/detail/overlay/get_ring.hpp b/boost/geometry/algorithms/detail/overlay/get_ring.hpp index ad47665500..c2c6980577 100644 --- a/boost/geometry/algorithms/detail/overlay/get_ring.hpp +++ b/boost/geometry/algorithms/detail/overlay/get_ring.hpp @@ -83,7 +83,7 @@ struct get_ring<polygon_tag> BOOST_ASSERT ( id.ring_index >= -1 - && id.ring_index < boost::size(interior_rings(polygon)) + && id.ring_index < int(boost::size(interior_rings(polygon))) ); return id.ring_index < 0 ? exterior_ring(polygon) diff --git a/boost/geometry/algorithms/detail/overlay/get_turn_info.hpp b/boost/geometry/algorithms/detail/overlay/get_turn_info.hpp index 663d70d9af..b8320d9b7b 100644 --- a/boost/geometry/algorithms/detail/overlay/get_turn_info.hpp +++ b/boost/geometry/algorithms/detail/overlay/get_turn_info.hpp @@ -250,9 +250,15 @@ struct touch : public base_turn_handler int const side_pk_q2 = SideStrategy::apply(qj, qk, pk); int const side_pk_p = SideStrategy::apply(pi, pj, pk); int const side_qk_q = SideStrategy::apply(qi, qj, qk); + + bool const both_continue = side_pk_p == 0 && side_qk_q == 0; + bool const robustness_issue_in_continue = both_continue && side_pk_q2 != 0; + bool const q_turns_left = side_qk_q == 1; bool const block_q = side_qk_p1 == 0 - && ! same(side_qi_p1, side_qk_q); + && ! same(side_qi_p1, side_qk_q) + && ! robustness_issue_in_continue + ; // If Pk at same side as Qi/Qk // (the "or" is for collinear case) @@ -278,7 +284,7 @@ struct touch : public base_turn_handler if (side_pk_q1 == 0) { ti.operations[0].operation = operation_blocked; - // Q turns right -> union (both independant), + // Q turns right -> union (both independent), // Q turns left -> intersection ti.operations[1].operation = block_q ? operation_blocked : q_turns_left ? operation_intersection @@ -466,6 +472,45 @@ struct equal : public base_turn_handler template < typename TurnInfo, + typename AssignPolicy +> +struct equal_opposite : public base_turn_handler +{ + template + < + typename Point1, + typename Point2, + typename OutputIterator, + typename IntersectionInfo, + typename DirInfo + > + static inline void apply(Point1 const& pi, Point2 const& qi, + /* by value: */ TurnInfo tp, + OutputIterator& out, + IntersectionInfo const& intersection_info, + DirInfo const& dir_info) + { + // For equal-opposite segments, normally don't do anything. + if (AssignPolicy::include_opposite) + { + tp.method = method_equal; + for (int i = 0; i < 2; i++) + { + tp.operations[i].operation = operation_opposite; + } + for (unsigned int i = 0; i < intersection_info.count; i++) + { + geometry::convert(intersection_info.intersections[i], tp.point); + AssignPolicy::apply(tp, pi, qi, intersection_info, dir_info); + *out++ = tp; + } + } + } +}; + +template +< + typename TurnInfo, typename SideStrategy > struct collinear : public base_turn_handler @@ -494,6 +539,13 @@ struct collinear : public base_turn_handler - if P arrives and P turns right: intersection for P - if Q arrives and Q turns left: union for Q (=intersection for P) - if Q arrives and Q turns right: intersection for Q (=union for P) + + ROBUSTNESS: p and q are collinear, so you would expect + that side qk//p1 == pk//q1. But that is not always the case + in near-epsilon ranges. Then decision logic is different. + If p arrives, q is further, so the angle qk//p1 is (normally) + more precise than pk//p1 + */ template < @@ -516,12 +568,18 @@ struct collinear : public base_turn_handler // Should not be 0, this is checked before BOOST_ASSERT(arrival != 0); + int const side_p = SideStrategy::apply(pi, pj, pk); + int const side_q = SideStrategy::apply(qi, qj, qk); + // If p arrives, use p, else use q int const side_p_or_q = arrival == 1 - ? SideStrategy::apply(pi, pj, pk) - : SideStrategy::apply(qi, qj, qk) + ? side_p + : side_q ; + int const side_pk = SideStrategy::apply(qi, qj, pk); + int const side_qk = SideStrategy::apply(pi, pj, qk); + // See comments above, // resulting in a strange sort of mathematic rule here: // The arrival-info multiplied by the relevant side @@ -529,7 +587,15 @@ struct collinear : public base_turn_handler int const product = arrival * side_p_or_q; - if(product == 0) + // Robustness: side_p is supposed to be equal to side_pk (because p/q are collinear) + // and side_q to side_qk + bool const robustness_issue = side_pk != side_p || side_qk != side_q; + + if (robustness_issue) + { + handle_robustness(ti, arrival, side_p, side_q, side_pk, side_qk); + } + else if(product == 0) { both(ti, operation_continue); } @@ -538,6 +604,38 @@ struct collinear : public base_turn_handler ui_else_iu(product == 1, ti); } } + + static inline void handle_robustness(TurnInfo& ti, int arrival, + int side_p, int side_q, int side_pk, int side_qk) + { + // We take the longer one, i.e. if q arrives in p (arrival == -1), + // then p exceeds q and we should take p for a union... + + bool use_p_for_union = arrival == -1; + + // ... unless one of the sides consistently directs to the other side + int const consistent_side_p = side_p == side_pk ? side_p : 0; + int const consistent_side_q = side_q == side_qk ? side_q : 0; + if (arrival == -1 && (consistent_side_p == -1 || consistent_side_q == 1)) + { + use_p_for_union = false; + } + if (arrival == 1 && (consistent_side_p == 1 || consistent_side_q == -1)) + { + use_p_for_union = true; + } + + //std::cout << "ROBUSTNESS -> Collinear " + // << " arr: " << arrival + // << " dir: " << side_p << " " << side_q + // << " rev: " << side_pk << " " << side_qk + // << " cst: " << cside_p << " " << cside_q + // << std::boolalpha << " " << use_p_for_union + // << std::endl; + + ui_else_iu(use_p_for_union, ti); + } + }; template @@ -583,6 +681,7 @@ private : TurnInfo& tp, IntersectionInfo const& intersection_info) { int const side_rk_r = SideStrategy::apply(ri, rj, rk); + operation_type blocked = operation_blocked; switch(side_rk_r) { @@ -596,15 +695,24 @@ private : break; case 0 : // No turn on opposite collinear: block, do not traverse - // But this "xx" is ignored here, it is useless to include - // two operation blocked, so the whole point does not need + // But this "xx" is usually ignored, it is useless to include + // two operations blocked, so the whole point does not need // to be generated. // So return false to indicate nothing is to be done. - return false; + if (AssignPolicy::include_opposite) + { + tp.operations[Index].operation = operation_opposite; + blocked = operation_opposite; + } + else + { + return false; + } + break; } // The other direction is always blocked when collinear opposite - tp.operations[1 - Index].operation = operation_blocked; + tp.operations[1 - Index].operation = blocked; // If P arrives within Q, set info on P (which is done above, index=0), // this turn-info belongs to the second intersection point, index=1 @@ -633,32 +741,45 @@ public: IntersectionInfo const& intersection_info, DirInfo const& dir_info) { - /* - std::cout << "arrivals: " - << dir_info.arrival[0] - << "/" << dir_info.arrival[1] - << std::endl; - */ - TurnInfo tp = tp_model; tp.method = method_collinear; - // If P arrives within Q, there is a turn dependant on P + // If P arrives within Q, there is a turn dependent on P if (dir_info.arrival[0] == 1 && set_tp<0>(pi, pj, pk, tp, intersection_info)) { - AssignPolicy::apply(tp, pi, qi); + AssignPolicy::apply(tp, pi, qi, intersection_info, dir_info); *out++ = tp; } - // If Q arrives within P, there is a turn dependant on Q + // If Q arrives within P, there is a turn dependent on Q if (dir_info.arrival[1] == 1 && set_tp<1>(qi, qj, qk, tp, intersection_info)) { - AssignPolicy::apply(tp, pi, qi); + AssignPolicy::apply(tp, pi, qi, intersection_info, dir_info); *out++ = tp; } + + if (AssignPolicy::include_opposite) + { + // Handle cases not yet handled above + if ((dir_info.arrival[1] == -1 && dir_info.arrival[0] == 0) + || (dir_info.arrival[0] == -1 && dir_info.arrival[1] == 0)) + { + for (int i = 0; i < 2; i++) + { + tp.operations[i].operation = operation_opposite; + } + for (unsigned int i = 0; i < intersection_info.count; i++) + { + geometry::convert(intersection_info.intersections[i], tp.point); + AssignPolicy::apply(tp, pi, qi, intersection_info, dir_info); + *out++ = tp; + } + } + } + } }; @@ -722,9 +843,17 @@ struct assign_null_policy { static bool const include_no_turn = false; static bool const include_degenerate = false; - - template <typename Point1, typename Point2, typename Info> - static inline void apply(Info& , Point1 const& , Point2 const& ) + static bool const include_opposite = false; + + template + < + typename Info, + typename Point1, + typename Point2, + typename IntersectionInfo, + typename DirInfo + > + static inline void apply(Info& , Point1 const& , Point2 const&, IntersectionInfo const&, DirInfo const&) {} }; @@ -763,8 +892,9 @@ struct get_turn_info typedef typename si::segment_intersection_strategy_type strategy; - - + // Intersect pi-pj with qi-qj + // The points pk and qk are only used do determine more information + // about the turn. template <typename OutputIterator> static inline OutputIterator apply( Point1 const& pi, Point1 const& pj, Point1 const& pk, @@ -796,7 +926,7 @@ struct get_turn_info { only_convert<TurnInfo>::apply(tp, result.template get<0>()); - AssignPolicy::apply(tp, pi, qi); + AssignPolicy::apply(tp, pi, qi, result.template get<0>(), result.template get<1>()); *out++ = tp; } break; @@ -824,7 +954,7 @@ struct get_turn_info policy::template apply<1>(qi, qj, qk, pi, pj, pk, tp, result.template get<0>(), result.template get<1>()); } - AssignPolicy::apply(tp, pi, qi); + AssignPolicy::apply(tp, pi, qi, result.template get<0>(), result.template get<1>()); *out++ = tp; } break; @@ -838,7 +968,7 @@ struct get_turn_info policy::apply(pi, pj, pk, qi, qj, qk, tp, result.template get<0>(), result.template get<1>()); - AssignPolicy::apply(tp, pi, qi); + AssignPolicy::apply(tp, pi, qi, result.template get<0>(), result.template get<1>()); *out++ = tp; } break; @@ -853,7 +983,7 @@ struct get_turn_info policy::apply(pi, pj, pk, qi, qj, qk, tp, result.template get<0>(), result.template get<1>()); - AssignPolicy::apply(tp, pi, qi); + AssignPolicy::apply(tp, pi, qi, result.template get<0>(), result.template get<1>()); *out++ = tp; } break; @@ -871,10 +1001,18 @@ struct get_turn_info policy::apply(pi, pj, pk, qi, qj, qk, tp, result.template get<0>(), result.template get<1>()); - AssignPolicy::apply(tp, pi, qi); + AssignPolicy::apply(tp, pi, qi, result.template get<0>(), result.template get<1>()); *out++ = tp; } - // If they ARE opposite, don't do anything. + else + { + equal_opposite + < + TurnInfo, + AssignPolicy + >::apply(pi, qi, + tp, out, result.template get<0>(), result.template get<1>()); + } } break; case 'c' : @@ -906,7 +1044,7 @@ struct get_turn_info tp, result.template get<0>(), result.template get<1>()); } - AssignPolicy::apply(tp, pi, qi); + AssignPolicy::apply(tp, pi, qi, result.template get<0>(), result.template get<1>()); *out++ = tp; } else @@ -927,7 +1065,7 @@ struct get_turn_info if (AssignPolicy::include_degenerate) { only_convert<TurnInfo>::apply(tp, result.template get<0>()); - AssignPolicy::apply(tp, pi, qi); + AssignPolicy::apply(tp, pi, qi, result.template get<0>(), result.template get<1>()); *out++ = tp; } } diff --git a/boost/geometry/algorithms/detail/overlay/get_turns.hpp b/boost/geometry/algorithms/detail/overlay/get_turns.hpp index 5740a8b7ba..26629043cb 100644 --- a/boost/geometry/algorithms/detail/overlay/get_turns.hpp +++ b/boost/geometry/algorithms/detail/overlay/get_turns.hpp @@ -495,7 +495,7 @@ struct get_turns_cs int source_id1, Range const& range, int source_id2, Box const& box, Turns& turns, - InterruptPolicy& , + InterruptPolicy& interrupt_policy, int multi_index = -1, int ring_index = -1) { if (boost::size(range) <= 1) @@ -557,7 +557,7 @@ struct get_turns_cs get_turns_with_box(seg_id, source_id2, *prev, *it, *next, bp[0], bp[1], bp[2], bp[3], - turns); + turns, interrupt_policy); // Future performance enhancement: // return if told by the interrupt policy } @@ -596,7 +596,8 @@ private: box_point_type const& bp2, box_point_type const& bp3, // Output - Turns& turns) + Turns& turns, + InterruptPolicy& interrupt_policy) { // Depending on code some relations can be left out @@ -622,6 +623,12 @@ private: ti.operations[1].seg_id = segment_identifier(source_id2, -1, -1, 3); TurnPolicy::apply(rp0, rp1, rp2, bp3, bp0, bp1, ti, std::back_inserter(turns)); + + if (InterruptPolicy::enabled) + { + interrupt_policy.apply(turns); + } + } }; diff --git a/boost/geometry/algorithms/detail/overlay/overlay.hpp b/boost/geometry/algorithms/detail/overlay/overlay.hpp index ab5b6d123d..41665e0af1 100644 --- a/boost/geometry/algorithms/detail/overlay/overlay.hpp +++ b/boost/geometry/algorithms/detail/overlay/overlay.hpp @@ -233,7 +233,7 @@ std::cout << "traverse" << std::endl; std::cout << "map_turns: " << timer.elapsed() << std::endl; #endif - typedef ring_properties<typename geometry::point_type<Geometry1>::type> properties; + typedef ring_properties<typename geometry::point_type<GeometryOut>::type> properties; std::map<ring_identifier, properties> selected; select_rings<Direction>(geometry1, geometry2, map, selected, ! turn_points.empty()); diff --git a/boost/geometry/algorithms/detail/overlay/turn_info.hpp b/boost/geometry/algorithms/detail/overlay/turn_info.hpp index aa6b428f19..89a60b21ab 100644 --- a/boost/geometry/algorithms/detail/overlay/turn_info.hpp +++ b/boost/geometry/algorithms/detail/overlay/turn_info.hpp @@ -28,7 +28,8 @@ enum operation_type operation_union, operation_intersection, operation_blocked, - operation_continue + operation_continue, + operation_opposite }; @@ -102,6 +103,12 @@ struct turn_info { return has12(type, type); } + + inline bool has(operation_type type) const + { + return this->operations[0].operation == type + || this->operations[1].operation == type; + } inline bool combination(operation_type type1, operation_type type2) const { @@ -114,10 +121,13 @@ struct turn_info { return both(operation_blocked); } + inline bool opposite() const + { + return both(operation_opposite); + } inline bool any_blocked() const { - return this->operations[0].operation == operation_blocked - || this->operations[1].operation == operation_blocked; + return has(operation_blocked); } diff --git a/boost/geometry/algorithms/detail/partition.hpp b/boost/geometry/algorithms/detail/partition.hpp index 7a7de2cdd3..45ff52ccb1 100644 --- a/boost/geometry/algorithms/detail/partition.hpp +++ b/boost/geometry/algorithms/detail/partition.hpp @@ -153,12 +153,12 @@ class partition_one_collection static inline void next_level(Box const& box, InputCollection const& collection, index_vector_type const& input, - int level, int min_elements, + int level, std::size_t min_elements, Policy& policy, VisitBoxPolicy& box_policy) { if (boost::size(input) > 0) { - if (boost::size(input) > min_elements && level < 100) + if (std::size_t(boost::size(input)) > min_elements && level < 100) { sub_divide::apply(box, collection, input, level + 1, min_elements, policy, box_policy); @@ -176,7 +176,7 @@ public : InputCollection const& collection, index_vector_type const& input, int level, - int min_elements, + std::size_t min_elements, Policy& policy, VisitBoxPolicy& box_policy) { box_policy.apply(box, level); @@ -230,13 +230,13 @@ class partition_two_collections index_vector_type const& input1, InputCollection const& collection2, index_vector_type const& input2, - int level, int min_elements, + int level, std::size_t min_elements, Policy& policy, VisitBoxPolicy& box_policy) { if (boost::size(input1) > 0 && boost::size(input2) > 0) { - if (boost::size(input1) > min_elements - && boost::size(input2) > min_elements + if (std::size_t(boost::size(input1)) > min_elements + && std::size_t(boost::size(input2)) > min_elements && level < 100) { sub_divide::apply(box, collection1, input1, collection2, @@ -257,7 +257,7 @@ public : InputCollection const& collection1, index_vector_type const& input1, InputCollection const& collection2, index_vector_type const& input2, int level, - int min_elements, + std::size_t min_elements, Policy& policy, VisitBoxPolicy& box_policy) { box_policy.apply(box, level); @@ -335,11 +335,11 @@ public : template <typename InputCollection, typename VisitPolicy> static inline void apply(InputCollection const& collection, VisitPolicy& visitor, - int min_elements = 16, + std::size_t min_elements = 16, VisitBoxPolicy box_visitor = visit_no_policy() ) { - if (boost::size(collection) > min_elements) + if (std::size_t(boost::size(collection)) > min_elements) { index_vector_type index_vector; Box total; @@ -377,12 +377,12 @@ public : static inline void apply(InputCollection const& collection1, InputCollection const& collection2, VisitPolicy& visitor, - int min_elements = 16, + std::size_t min_elements = 16, VisitBoxPolicy box_visitor = visit_no_policy() ) { - if (boost::size(collection1) > min_elements - && boost::size(collection2) > min_elements) + if (std::size_t(boost::size(collection1)) > min_elements + && std::size_t(boost::size(collection2)) > min_elements) { index_vector_type index_vector1, index_vector2; Box total; diff --git a/boost/geometry/algorithms/detail/point_on_border.hpp b/boost/geometry/algorithms/detail/point_on_border.hpp index 33177924aa..b7e15ba3f9 100644 --- a/boost/geometry/algorithms/detail/point_on_border.hpp +++ b/boost/geometry/algorithms/detail/point_on_border.hpp @@ -25,6 +25,7 @@ #include <boost/geometry/geometries/concepts/check.hpp> #include <boost/geometry/algorithms/assign.hpp> +#include <boost/geometry/algorithms/detail/convert_point_to_point.hpp> #include <boost/geometry/algorithms/detail/disjoint.hpp> @@ -50,7 +51,8 @@ struct get_point template<typename Point, std::size_t Dimension, std::size_t DimensionCount> struct midpoint_helper { - static inline bool apply(Point& p, Point const& p1, Point const& p2) + template <typename InputPoint> + static inline bool apply(Point& p, InputPoint const& p1, InputPoint const& p2) { typename coordinate_type<Point>::type const two = 2; set<Dimension>(p, @@ -63,7 +65,8 @@ struct midpoint_helper template <typename Point, std::size_t DimensionCount> struct midpoint_helper<Point, DimensionCount, DimensionCount> { - static inline bool apply(Point& , Point const& , Point const& ) + template <typename InputPoint> + static inline bool apply(Point& , InputPoint const& , InputPoint const& ) { return true; } @@ -102,7 +105,7 @@ struct point_on_range if (n > 0) { - point = *boost::begin(range); + geometry::detail::conversion::convert_point_to_point(*boost::begin(range), point); return true; } return false; diff --git a/boost/geometry/algorithms/detail/sections/sectionalize.hpp b/boost/geometry/algorithms/detail/sections/sectionalize.hpp index 36bcbdd6e7..a6e6837fe4 100644 --- a/boost/geometry/algorithms/detail/sections/sectionalize.hpp +++ b/boost/geometry/algorithms/detail/sections/sectionalize.hpp @@ -254,7 +254,7 @@ struct sectionalize_part Range const& range, ring_identifier ring_id) { - if (boost::size(range) <= index) + if (int(boost::size(range)) <= index) { return; } diff --git a/boost/geometry/algorithms/disjoint.hpp b/boost/geometry/algorithms/disjoint.hpp index f36b8dbdd0..f986cc24af 100644 --- a/boost/geometry/algorithms/disjoint.hpp +++ b/boost/geometry/algorithms/disjoint.hpp @@ -27,6 +27,7 @@ #include <boost/geometry/core/reverse_dispatch.hpp> #include <boost/geometry/algorithms/detail/disjoint.hpp> +#include <boost/geometry/algorithms/detail/for_each_range.hpp> #include <boost/geometry/algorithms/detail/point_on_border.hpp> #include <boost/geometry/algorithms/detail/overlay/get_turns.hpp> #include <boost/geometry/algorithms/within.hpp> @@ -44,15 +45,57 @@ namespace boost { namespace geometry namespace detail { namespace disjoint { +template<typename Geometry> +struct check_each_ring_for_within +{ + bool has_within; + Geometry const& m_geometry; + + inline check_each_ring_for_within(Geometry const& g) + : has_within(false) + , m_geometry(g) + {} + + template <typename Range> + inline void apply(Range const& range) + { + typename geometry::point_type<Range>::type p; + geometry::point_on_border(p, range); + if (geometry::within(p, m_geometry)) + { + has_within = true; + } + } +}; + +template <typename FirstGeometry, typename SecondGeometry> +inline bool rings_containing(FirstGeometry const& geometry1, + SecondGeometry const& geometry2) +{ + check_each_ring_for_within<FirstGeometry> checker(geometry1); + geometry::detail::for_each_range(geometry2, checker); + return checker.has_within; +} + + struct assign_disjoint_policy { // We want to include all points: static bool const include_no_turn = true; static bool const include_degenerate = true; + static bool const include_opposite = true; // We don't assign extra info: - template <typename Point1, typename Point2, typename Info> - static inline void apply(Info& , Point1 const& , Point2 const& ) + template + < + typename Info, + typename Point1, + typename Point2, + typename IntersectionInfo, + typename DirInfo + > + static inline void apply(Info& , Point1 const& , Point2 const&, + IntersectionInfo const&, DirInfo const&) {} }; @@ -107,8 +150,6 @@ struct disjoint_segment } }; - - template <typename Geometry1, typename Geometry2> struct general_areal { @@ -119,20 +160,10 @@ struct general_areal return false; } - typedef typename geometry::point_type<Geometry1>::type point_type; - // If there is no intersection of segments, they might located // inside each other - point_type p1; - geometry::point_on_border(p1, geometry1); - if (geometry::within(p1, geometry2)) - { - return false; - } - - typename geometry::point_type<Geometry1>::type p2; - geometry::point_on_border(p2, geometry2); - if (geometry::within(p2, geometry1)) + if (rings_containing(geometry1, geometry2) + || rings_containing(geometry2, geometry1)) { return false; } diff --git a/boost/geometry/algorithms/distance.hpp b/boost/geometry/algorithms/distance.hpp index eca3b03c7b..11c2bc929b 100644 --- a/boost/geometry/algorithms/distance.hpp +++ b/boost/geometry/algorithms/distance.hpp @@ -130,7 +130,7 @@ struct point_to_range // check if other segments are closer for (++prev, ++it; it != boost::end(view); ++prev, ++it) { - return_type const ds = ps_strategy.apply(point, *prev, *it); + return_type const ds = eps_strategy.apply(point, *prev, *it); if (geometry::math::equals(ds, zero)) { return ds; diff --git a/boost/geometry/algorithms/touches.hpp b/boost/geometry/algorithms/touches.hpp new file mode 100644 index 0000000000..7d424af428 --- /dev/null +++ b/boost/geometry/algorithms/touches.hpp @@ -0,0 +1,181 @@ +// Boost.Geometry (aka GGL, Generic Geometry Library) + +// Copyright (c) 2007-2012 Barend Gehrels, Amsterdam, the Netherlands. +// Copyright (c) 2008-2012 Bruno Lalande, Paris, France. +// Copyright (c) 2009-2012 Mateusz Loskot, London, UK. + +// Parts of Boost.Geometry are redesigned from Geodan's Geographic Library +// (geolib/GGL), copyright (c) 1995-2010 Geodan, Amsterdam, the Netherlands. + +// 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) + +#ifndef BOOST_GEOMETRY_ALGORITHMS_TOUCHES_HPP +#define BOOST_GEOMETRY_ALGORITHMS_TOUCHES_HPP + + +#include <deque> + +#include <boost/geometry/geometries/concepts/check.hpp> +#include <boost/geometry/algorithms/detail/overlay/self_turn_points.hpp> +#include <boost/geometry/algorithms/detail/overlay/get_turns.hpp> +#include <boost/geometry/algorithms/disjoint.hpp> +#include <boost/geometry/algorithms/intersects.hpp> +#include <boost/geometry/algorithms/num_geometries.hpp> + + +namespace boost { namespace geometry +{ + +namespace detail { namespace touches +{ + +template <typename Turn> +inline bool ok_for_touch(Turn const& turn) +{ + return turn.both(detail::overlay::operation_union) + || turn.both(detail::overlay::operation_blocked) + || turn.combination(detail::overlay::operation_union, detail::overlay::operation_blocked) + ; +} + +template <typename Turns> +inline bool has_only_turns(Turns const& turns) +{ + bool has_touch = false; + typedef typename boost::range_iterator<Turns const>::type iterator_type; + for (iterator_type it = boost::begin(turns); it != boost::end(turns); ++it) + { + if (it->has(detail::overlay::operation_intersection)) + { + return false; + } + + switch(it->method) + { + case detail::overlay::method_crosses: + return false; + case detail::overlay::method_equal: + // Segment spatially equal means: at the right side + // the polygon internally overlaps. So return false. + return false; + case detail::overlay::method_touch: + case detail::overlay::method_touch_interior: + case detail::overlay::method_collinear: + if (ok_for_touch(*it)) + { + has_touch = true; + } + else + { + return false; + } + break; + case detail::overlay::method_none : + case detail::overlay::method_disjoint : + case detail::overlay::method_error : + break; + } + } + return has_touch; +} + +}} + +/*! +\brief \brief_check{has at least one touching point (self-tangency)} +\note This function can be called for one geometry (self-tangency) and + also for two geometries (touch) +\ingroup touches +\tparam Geometry \tparam_geometry +\param geometry \param_geometry +\return \return_check{is self-touching} + +\qbk{distinguish,one geometry} +\qbk{[def __one_parameter__]} +\qbk{[include reference/algorithms/touches.qbk]} +*/ +template <typename Geometry> +inline bool touches(Geometry const& geometry) +{ + concept::check<Geometry const>(); + + typedef detail::overlay::turn_info + < + typename geometry::point_type<Geometry>::type + > turn_info; + + typedef detail::overlay::get_turn_info + < + typename point_type<Geometry>::type, + typename point_type<Geometry>::type, + turn_info, + detail::overlay::assign_null_policy + > policy_type; + + std::deque<turn_info> turns; + detail::self_get_turn_points::no_interrupt_policy policy; + detail::self_get_turn_points::get_turns + < + Geometry, + std::deque<turn_info>, + policy_type, + detail::self_get_turn_points::no_interrupt_policy + >::apply(geometry, turns, policy); + + return detail::touches::has_only_turns(turns); +} + + +/*! +\brief \brief_check2{have at least one touching point (tangent - non overlapping)} +\ingroup touches +\tparam Geometry1 \tparam_geometry +\tparam Geometry2 \tparam_geometry +\param geometry1 \param_geometry +\param geometry2 \param_geometry +\return \return_check2{touch each other} + +\qbk{distinguish,two geometries} +\qbk{[include reference/algorithms/touches.qbk]} + */ +template <typename Geometry1, typename Geometry2> +inline bool touches(Geometry1 const& geometry1, Geometry2 const& geometry2) +{ + concept::check<Geometry1 const>(); + concept::check<Geometry2 const>(); + + + typedef detail::overlay::turn_info + < + typename geometry::point_type<Geometry1>::type + > turn_info; + + typedef detail::overlay::get_turn_info + < + typename point_type<Geometry1>::type, + typename point_type<Geometry2>::type, + turn_info, + detail::overlay::assign_null_policy + > policy_type; + + std::deque<turn_info> turns; + detail::get_turns::no_interrupt_policy policy; + boost::geometry::get_turns + < + false, false, + detail::overlay::assign_null_policy + >(geometry1, geometry2, turns, policy); + + return detail::touches::has_only_turns(turns) + && ! geometry::detail::disjoint::rings_containing(geometry1, geometry2) + && ! geometry::detail::disjoint::rings_containing(geometry2, geometry1) + ; +} + + + +}} // namespace boost::geometry + +#endif // BOOST_GEOMETRY_ALGORITHMS_TOUCHES_HPP diff --git a/boost/geometry/algorithms/within.hpp b/boost/geometry/algorithms/within.hpp index 0e0cb68d53..f1f0993d77 100644 --- a/boost/geometry/algorithms/within.hpp +++ b/boost/geometry/algorithms/within.hpp @@ -65,7 +65,7 @@ struct point_in_ring static inline int apply(Point const& point, Ring const& ring, Strategy const& strategy) { - if (boost::size(ring) + if (int(boost::size(ring)) < core_detail::closure::minimum_ring_size<Closure>::value) { return -1; |