1 // Boost.Geometry (aka GGL, Generic Geometry Library) 2 3 // Copyright (c) 2007-2012 Barend Gehrels, Amsterdam, the Netherlands. 4 5 // This file was modified by Oracle on 2013-2020. 6 // Modifications copyright (c) 2013-2020 Oracle and/or its affiliates. 7 8 // Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle 9 10 // Use, modification and distribution is subject to the Boost Software License, 11 // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at 12 // http://www.boost.org/LICENSE_1_0.txt) 13 14 #ifndef BOOST_GEOMETRY_ALGORITHMS_DETAIL_RELATE_LINEAR_AREAL_HPP 15 #define BOOST_GEOMETRY_ALGORITHMS_DETAIL_RELATE_LINEAR_AREAL_HPP 16 17 #include <boost/core/ignore_unused.hpp> 18 #include <boost/range/size.hpp> 19 20 #include <boost/geometry/core/assert.hpp> 21 #include <boost/geometry/core/topological_dimension.hpp> 22 23 #include <boost/geometry/util/condition.hpp> 24 #include <boost/geometry/util/range.hpp> 25 #include <boost/geometry/util/type_traits.hpp> 26 27 #include <boost/geometry/algorithms/num_interior_rings.hpp> 28 #include <boost/geometry/algorithms/detail/point_on_border.hpp> 29 #include <boost/geometry/algorithms/detail/sub_range.hpp> 30 #include <boost/geometry/algorithms/detail/single_geometry.hpp> 31 32 #include <boost/geometry/algorithms/detail/relate/point_geometry.hpp> 33 #include <boost/geometry/algorithms/detail/relate/turns.hpp> 34 #include <boost/geometry/algorithms/detail/relate/boundary_checker.hpp> 35 #include <boost/geometry/algorithms/detail/relate/follow_helpers.hpp> 36 37 #include <boost/geometry/views/detail/normalized_view.hpp> 38 39 namespace boost { namespace geometry 40 { 41 42 #ifndef DOXYGEN_NO_DETAIL 43 namespace detail { namespace relate { 44 45 // WARNING! 46 // TODO: In the worst case calling this Pred in a loop for MultiLinestring/MultiPolygon may take O(NM) 47 // Use the rtree in this case! 48 49 // may be used to set IE and BE for a Linear geometry for which no turns were generated 50 template 51 < 52 typename Geometry2, 53 typename Result, 54 typename PointInArealStrategy, 55 typename BoundaryChecker, 56 bool TransposeResult 57 > 58 class no_turns_la_linestring_pred 59 { 60 public: no_turns_la_linestring_pred(Geometry2 const & geometry2,Result & res,PointInArealStrategy const & point_in_areal_strategy,BoundaryChecker const & boundary_checker)61 no_turns_la_linestring_pred(Geometry2 const& geometry2, 62 Result & res, 63 PointInArealStrategy const& point_in_areal_strategy, 64 BoundaryChecker const& boundary_checker) 65 : m_geometry2(geometry2) 66 , m_result(res) 67 , m_point_in_areal_strategy(point_in_areal_strategy) 68 , m_boundary_checker(boundary_checker) 69 , m_interrupt_flags(0) 70 { 71 if ( ! may_update<interior, interior, '1', TransposeResult>(m_result) ) 72 { 73 m_interrupt_flags |= 1; 74 } 75 76 if ( ! may_update<interior, exterior, '1', TransposeResult>(m_result) ) 77 { 78 m_interrupt_flags |= 2; 79 } 80 81 if ( ! may_update<boundary, interior, '0', TransposeResult>(m_result) ) 82 { 83 m_interrupt_flags |= 4; 84 } 85 86 if ( ! may_update<boundary, exterior, '0', TransposeResult>(m_result) ) 87 { 88 m_interrupt_flags |= 8; 89 } 90 } 91 92 template <typename Linestring> operator ()(Linestring const & linestring)93 bool operator()(Linestring const& linestring) 94 { 95 std::size_t const count = boost::size(linestring); 96 97 // invalid input 98 if ( count < 2 ) 99 { 100 // ignore 101 // TODO: throw an exception? 102 return true; 103 } 104 105 // if those flags are set nothing will change 106 if ( m_interrupt_flags == 0xF ) 107 { 108 return false; 109 } 110 111 int const pig = detail::within::point_in_geometry(range::front(linestring), 112 m_geometry2, 113 m_point_in_areal_strategy); 114 //BOOST_GEOMETRY_ASSERT_MSG(pig != 0, "There should be no IPs"); 115 116 if ( pig > 0 ) 117 { 118 update<interior, interior, '1', TransposeResult>(m_result); 119 m_interrupt_flags |= 1; 120 } 121 else 122 { 123 update<interior, exterior, '1', TransposeResult>(m_result); 124 m_interrupt_flags |= 2; 125 } 126 127 // check if there is a boundary 128 if ( ( m_interrupt_flags & 0xC ) != 0xC // if wasn't already set 129 && ( m_boundary_checker.template 130 is_endpoint_boundary<boundary_front>(range::front(linestring)) 131 || m_boundary_checker.template 132 is_endpoint_boundary<boundary_back>(range::back(linestring)) ) ) 133 { 134 if ( pig > 0 ) 135 { 136 update<boundary, interior, '0', TransposeResult>(m_result); 137 m_interrupt_flags |= 4; 138 } 139 else 140 { 141 update<boundary, exterior, '0', TransposeResult>(m_result); 142 m_interrupt_flags |= 8; 143 } 144 } 145 146 return m_interrupt_flags != 0xF 147 && ! m_result.interrupt; 148 } 149 150 private: 151 Geometry2 const& m_geometry2; 152 Result & m_result; 153 PointInArealStrategy const& m_point_in_areal_strategy; 154 BoundaryChecker const& m_boundary_checker; 155 unsigned m_interrupt_flags; 156 }; 157 158 // may be used to set EI and EB for an Areal geometry for which no turns were generated 159 template <typename Result, bool TransposeResult> 160 class no_turns_la_areal_pred 161 { 162 public: no_turns_la_areal_pred(Result & res)163 no_turns_la_areal_pred(Result & res) 164 : m_result(res) 165 , interrupt(! may_update<interior, exterior, '2', TransposeResult>(m_result) 166 && ! may_update<boundary, exterior, '1', TransposeResult>(m_result) ) 167 {} 168 169 template <typename Areal> operator ()(Areal const & areal)170 bool operator()(Areal const& areal) 171 { 172 if ( interrupt ) 173 { 174 return false; 175 } 176 177 // TODO: 178 // handle empty/invalid geometries in a different way than below? 179 180 typedef typename geometry::point_type<Areal>::type point_type; 181 point_type dummy; 182 bool const ok = boost::geometry::point_on_border(dummy, areal); 183 184 // TODO: for now ignore, later throw an exception? 185 if ( !ok ) 186 { 187 return true; 188 } 189 190 update<interior, exterior, '2', TransposeResult>(m_result); 191 update<boundary, exterior, '1', TransposeResult>(m_result); 192 193 return false; 194 } 195 196 private: 197 Result & m_result; 198 bool const interrupt; 199 }; 200 201 // The implementation of an algorithm calculating relate() for L/A 202 template <typename Geometry1, typename Geometry2, bool TransposeResult = false> 203 struct linear_areal 204 { 205 // check Linear / Areal 206 BOOST_STATIC_ASSERT(topological_dimension<Geometry1>::value == 1 207 && topological_dimension<Geometry2>::value == 2); 208 209 static const bool interruption_enabled = true; 210 211 typedef typename geometry::point_type<Geometry1>::type point1_type; 212 typedef typename geometry::point_type<Geometry2>::type point2_type; 213 214 template <typename Geom1, typename Geom2, typename Strategy> 215 struct multi_turn_info 216 : turns::get_turns<Geom1, Geom2>::template turn_info_type<Strategy>::type 217 { multi_turn_infoboost::geometry::detail::relate::linear_areal::multi_turn_info218 multi_turn_info() : priority(0) {} 219 int priority; // single-geometry sorting priority 220 }; 221 222 template <typename Geom1, typename Geom2, typename Strategy> 223 struct turn_info_type 224 : std::conditional 225 < 226 util::is_multi<Geometry2>::value, 227 multi_turn_info<Geom1, Geom2, Strategy>, 228 typename turns::get_turns<Geom1, Geom2>::template turn_info_type<Strategy>::type 229 > 230 {}; 231 232 template <typename Result, typename IntersectionStrategy> applyboost::geometry::detail::relate::linear_areal233 static inline void apply(Geometry1 const& geometry1, Geometry2 const& geometry2, 234 Result & result, 235 IntersectionStrategy const& intersection_strategy) 236 { 237 // TODO: If Areal geometry may have infinite size, change the following line: 238 239 // The result should be FFFFFFFFF 240 relate::set<exterior, exterior, result_dimension<Geometry2>::value, TransposeResult>(result);// FFFFFFFFd, d in [1,9] or T 241 242 if ( BOOST_GEOMETRY_CONDITION( result.interrupt ) ) 243 return; 244 245 // get and analyse turns 246 typedef typename turn_info_type<Geometry1, Geometry2, IntersectionStrategy>::type turn_type; 247 std::vector<turn_type> turns; 248 249 interrupt_policy_linear_areal<Geometry2, Result> interrupt_policy(geometry2, result); 250 251 turns::get_turns<Geometry1, Geometry2>::apply(turns, geometry1, geometry2, interrupt_policy, intersection_strategy); 252 if ( BOOST_GEOMETRY_CONDITION( result.interrupt ) ) 253 return; 254 255 typedef typename IntersectionStrategy::template point_in_geometry_strategy<Geometry1, Geometry2>::type within_strategy_type; 256 within_strategy_type const within_strategy = intersection_strategy.template get_point_in_geometry_strategy<Geometry1, Geometry2>(); 257 258 typedef typename IntersectionStrategy::cs_tag cs_tag; 259 typedef typename within_strategy_type::equals_point_point_strategy_type eq_pp_strategy_type; 260 261 typedef boundary_checker 262 < 263 Geometry1, 264 eq_pp_strategy_type 265 > boundary_checker1_type; 266 boundary_checker1_type boundary_checker1(geometry1); 267 268 no_turns_la_linestring_pred 269 < 270 Geometry2, 271 Result, 272 within_strategy_type, 273 boundary_checker1_type, 274 TransposeResult 275 > pred1(geometry2, 276 result, 277 within_strategy, 278 boundary_checker1); 279 for_each_disjoint_geometry_if<0, Geometry1>::apply(turns.begin(), turns.end(), geometry1, pred1); 280 if ( BOOST_GEOMETRY_CONDITION( result.interrupt ) ) 281 return; 282 283 no_turns_la_areal_pred<Result, !TransposeResult> pred2(result); 284 for_each_disjoint_geometry_if<1, Geometry2>::apply(turns.begin(), turns.end(), geometry2, pred2); 285 if ( BOOST_GEOMETRY_CONDITION( result.interrupt ) ) 286 return; 287 288 if ( turns.empty() ) 289 return; 290 291 // This is set here because in the case if empty Areal geometry were passed 292 // those shouldn't be set 293 relate::set<exterior, interior, '2', TransposeResult>(result);// FFFFFF2Fd 294 if ( BOOST_GEOMETRY_CONDITION( result.interrupt ) ) 295 return; 296 297 { 298 sort_dispatch<cs_tag>(turns.begin(), turns.end(), util::is_multi<Geometry2>()); 299 300 turns_analyser<turn_type> analyser; 301 analyse_each_turn(result, analyser, 302 turns.begin(), turns.end(), 303 geometry1, geometry2, 304 boundary_checker1, 305 intersection_strategy.get_side_strategy()); 306 307 if ( BOOST_GEOMETRY_CONDITION( result.interrupt ) ) 308 return; 309 } 310 311 // If 'c' (insersection_boundary) was not found we know that any Ls isn't equal to one of the Rings 312 if ( !interrupt_policy.is_boundary_found ) 313 { 314 relate::set<exterior, boundary, '1', TransposeResult>(result); 315 } 316 // Don't calculate it if it's required 317 else if ( may_update<exterior, boundary, '1', TransposeResult>(result) ) 318 { 319 // TODO: REVISE THIS CODE AND PROBABLY REWRITE SOME PARTS TO BE MORE HUMAN-READABLE 320 // IN GENERAL IT ANALYSES THE RINGS OF AREAL GEOMETRY AND DETECTS THE ONES THAT 321 // MAY OVERLAP THE INTERIOR OF LINEAR GEOMETRY (NO IPs OR NON-FAKE 'u' OPERATION) 322 // NOTE: For one case std::sort may be called again to sort data by operations for data already sorted by ring index 323 // In the worst case scenario the complexity will be O( NlogN + R*(N/R)log(N/R) ) 324 // So always should remain O(NlogN) -> for R==1 <-> 1(N/1)log(N/1), for R==N <-> N(N/N)log(N/N) 325 // Some benchmarking should probably be done to check if only one std::sort should be used 326 327 // sort by multi_index and rind_index 328 std::sort(turns.begin(), turns.end(), less_ring()); 329 330 typedef typename std::vector<turn_type>::iterator turn_iterator; 331 332 turn_iterator it = turns.begin(); 333 segment_identifier * prev_seg_id_ptr = NULL; 334 // for each ring 335 for ( ; it != turns.end() ; ) 336 { 337 // it's the next single geometry 338 if ( prev_seg_id_ptr == NULL 339 || prev_seg_id_ptr->multi_index != it->operations[1].seg_id.multi_index ) 340 { 341 // if the first ring has no IPs 342 if ( it->operations[1].seg_id.ring_index > -1 ) 343 { 344 // we can be sure that the exterior overlaps the boundary 345 relate::set<exterior, boundary, '1', TransposeResult>(result); 346 break; 347 } 348 // if there was some previous ring 349 if ( prev_seg_id_ptr != NULL ) 350 { 351 signed_size_type const next_ring_index = prev_seg_id_ptr->ring_index + 1; 352 BOOST_GEOMETRY_ASSERT(next_ring_index >= 0); 353 354 // if one of the last rings of previous single geometry was ommited 355 if ( static_cast<std::size_t>(next_ring_index) 356 < geometry::num_interior_rings( 357 single_geometry(geometry2, *prev_seg_id_ptr)) ) 358 { 359 // we can be sure that the exterior overlaps the boundary 360 relate::set<exterior, boundary, '1', TransposeResult>(result); 361 break; 362 } 363 } 364 } 365 // if it's the same single geometry 366 else /*if ( previous_multi_index == it->operations[1].seg_id.multi_index )*/ 367 { 368 // and we jumped over one of the rings 369 if ( prev_seg_id_ptr != NULL // just in case 370 && prev_seg_id_ptr->ring_index + 1 < it->operations[1].seg_id.ring_index ) 371 { 372 // we can be sure that the exterior overlaps the boundary 373 relate::set<exterior, boundary, '1', TransposeResult>(result); 374 break; 375 } 376 } 377 378 prev_seg_id_ptr = boost::addressof(it->operations[1].seg_id); 379 380 // find the next ring first iterator and check if the analysis should be performed 381 has_boundary_intersection has_boundary_inters; 382 turn_iterator next = find_next_ring(it, turns.end(), has_boundary_inters); 383 384 // if there is no 1d overlap with the boundary 385 if ( !has_boundary_inters.result ) 386 { 387 // we can be sure that the exterior overlaps the boundary 388 relate::set<exterior, boundary, '1', TransposeResult>(result); 389 break; 390 } 391 // else there is 1d overlap with the boundary so we must analyse the boundary 392 else 393 { 394 // u, c 395 typedef turns::less<1, turns::less_op_areal_linear<1>, cs_tag> less; 396 std::sort(it, next, less()); 397 398 eq_pp_strategy_type const eq_pp_strategy = within_strategy.get_equals_point_point_strategy(); 399 400 // analyse 401 areal_boundary_analyser<turn_type> analyser; 402 for ( turn_iterator rit = it ; rit != next ; ++rit ) 403 { 404 // if the analyser requests, break the search 405 if ( !analyser.apply(it, rit, next, eq_pp_strategy) ) 406 break; 407 } 408 409 // if the boundary of Areal goes out of the Linear 410 if ( analyser.is_union_detected ) 411 { 412 // we can be sure that the boundary of Areal overlaps the exterior of Linear 413 relate::set<exterior, boundary, '1', TransposeResult>(result); 414 break; 415 } 416 } 417 418 it = next; 419 } 420 421 // if there was some previous ring 422 if ( prev_seg_id_ptr != NULL ) 423 { 424 signed_size_type const next_ring_index = prev_seg_id_ptr->ring_index + 1; 425 BOOST_GEOMETRY_ASSERT(next_ring_index >= 0); 426 427 // if one of the last rings of previous single geometry was ommited 428 if ( static_cast<std::size_t>(next_ring_index) 429 < geometry::num_interior_rings( 430 single_geometry(geometry2, *prev_seg_id_ptr)) ) 431 { 432 // we can be sure that the exterior overlaps the boundary 433 relate::set<exterior, boundary, '1', TransposeResult>(result); 434 } 435 } 436 } 437 } 438 439 template <typename It, typename Pred, typename Comp> for_each_equal_rangeboost::geometry::detail::relate::linear_areal440 static void for_each_equal_range(It first, It last, Pred pred, Comp comp) 441 { 442 if ( first == last ) 443 return; 444 445 It first_equal = first; 446 It prev = first; 447 for ( ++first ; ; ++first, ++prev ) 448 { 449 if ( first == last || !comp(*prev, *first) ) 450 { 451 pred(first_equal, first); 452 first_equal = first; 453 } 454 455 if ( first == last ) 456 break; 457 } 458 } 459 460 struct same_ip 461 { 462 template <typename Turn> operator ()boost::geometry::detail::relate::linear_areal::same_ip463 bool operator()(Turn const& left, Turn const& right) const 464 { 465 return left.operations[0].seg_id == right.operations[0].seg_id 466 && left.operations[0].fraction == right.operations[0].fraction; 467 } 468 }; 469 470 struct same_ip_and_multi_index 471 { 472 template <typename Turn> operator ()boost::geometry::detail::relate::linear_areal::same_ip_and_multi_index473 bool operator()(Turn const& left, Turn const& right) const 474 { 475 return same_ip()(left, right) 476 && left.operations[1].seg_id.multi_index == right.operations[1].seg_id.multi_index; 477 } 478 }; 479 480 template <typename OpToPriority> 481 struct set_turns_group_priority 482 { 483 template <typename TurnIt> operator ()boost::geometry::detail::relate::linear_areal::set_turns_group_priority484 void operator()(TurnIt first, TurnIt last) const 485 { 486 BOOST_GEOMETRY_ASSERT(first != last); 487 static OpToPriority op_to_priority; 488 // find the operation with the least priority 489 int least_priority = op_to_priority(first->operations[0]); 490 for ( TurnIt it = first + 1 ; it != last ; ++it ) 491 { 492 int priority = op_to_priority(it->operations[0]); 493 if ( priority < least_priority ) 494 least_priority = priority; 495 } 496 // set the least priority for all turns of the group 497 for ( TurnIt it = first ; it != last ; ++it ) 498 { 499 it->priority = least_priority; 500 } 501 } 502 }; 503 504 template <typename SingleLess> 505 struct sort_turns_group 506 { 507 struct less 508 { 509 template <typename Turn> operator ()boost::geometry::detail::relate::linear_areal::sort_turns_group::less510 bool operator()(Turn const& left, Turn const& right) const 511 { 512 return left.operations[1].seg_id.multi_index == right.operations[1].seg_id.multi_index ? 513 SingleLess()(left, right) : 514 left.priority < right.priority; 515 } 516 }; 517 518 template <typename TurnIt> operator ()boost::geometry::detail::relate::linear_areal::sort_turns_group519 void operator()(TurnIt first, TurnIt last) const 520 { 521 std::sort(first, last, less()); 522 } 523 }; 524 525 template <typename CSTag, typename TurnIt> sort_dispatchboost::geometry::detail::relate::linear_areal526 static void sort_dispatch(TurnIt first, TurnIt last, std::true_type const& /*is_multi*/) 527 { 528 // sort turns by Linear seg_id, then by fraction, then by other multi_index 529 typedef turns::less<0, turns::less_other_multi_index<0>, CSTag> less; 530 std::sort(first, last, less()); 531 532 // For the same IP and multi_index - the same other's single geometry 533 // set priorities as the least operation found for the whole single geometry 534 // so e.g. single geometries containing 'u' will always be before those only containing 'i' 535 typedef turns::op_to_int<0,2,3,1,4,0> op_to_int_xuic; 536 for_each_equal_range(first, last, 537 set_turns_group_priority<op_to_int_xuic>(), // least operation in xuic order 538 same_ip_and_multi_index()); // other's multi_index 539 540 // When priorities for single geometries are set now sort turns for the same IP 541 // if multi_index is the same sort them according to the single-less 542 // else use priority of the whole single-geometry set earlier 543 typedef turns::less<0, turns::less_op_linear_areal_single<0>, CSTag> single_less; 544 for_each_equal_range(first, last, 545 sort_turns_group<single_less>(), 546 same_ip()); 547 } 548 549 template <typename CSTag, typename TurnIt> sort_dispatchboost::geometry::detail::relate::linear_areal550 static void sort_dispatch(TurnIt first, TurnIt last, std::false_type const& /*is_multi*/) 551 { 552 // sort turns by Linear seg_id, then by fraction, then 553 // for same ring id: x, u, i, c 554 // for different ring id: c, i, u, x 555 typedef turns::less<0, turns::less_op_linear_areal_single<0>, CSTag> less; 556 std::sort(first, last, less()); 557 } 558 559 560 // interrupt policy which may be passed to get_turns to interrupt the analysis 561 // based on the info in the passed result/mask 562 template <typename Areal, typename Result> 563 class interrupt_policy_linear_areal 564 { 565 public: 566 static bool const enabled = true; 567 interrupt_policy_linear_areal(Areal const & areal,Result & result)568 interrupt_policy_linear_areal(Areal const& areal, Result & result) 569 : m_result(result), m_areal(areal) 570 , is_boundary_found(false) 571 {} 572 573 // TODO: since we update result for some operations here, we may not do it in the analyser! 574 575 template <typename Range> apply(Range const & turns)576 inline bool apply(Range const& turns) 577 { 578 typedef typename boost::range_iterator<Range const>::type iterator; 579 580 for (iterator it = boost::begin(turns) ; it != boost::end(turns) ; ++it) 581 { 582 if ( it->operations[0].operation == overlay::operation_intersection ) 583 { 584 bool const no_interior_rings 585 = geometry::num_interior_rings( 586 single_geometry(m_areal, it->operations[1].seg_id)) == 0; 587 588 // WARNING! THIS IS TRUE ONLY IF THE POLYGON IS SIMPLE! 589 // OR WITHOUT INTERIOR RINGS (AND OF COURSE VALID) 590 if ( no_interior_rings ) 591 update<interior, interior, '1', TransposeResult>(m_result); 592 } 593 else if ( it->operations[0].operation == overlay::operation_continue ) 594 { 595 update<interior, boundary, '1', TransposeResult>(m_result); 596 is_boundary_found = true; 597 } 598 else if ( ( it->operations[0].operation == overlay::operation_union 599 || it->operations[0].operation == overlay::operation_blocked ) 600 && it->operations[0].position == overlay::position_middle ) 601 { 602 // TODO: here we could also check the boundaries and set BB at this point 603 update<interior, boundary, '0', TransposeResult>(m_result); 604 } 605 } 606 607 return m_result.interrupt; 608 } 609 610 private: 611 Result & m_result; 612 Areal const& m_areal; 613 614 public: 615 bool is_boundary_found; 616 }; 617 618 // This analyser should be used like Input or SinglePass Iterator 619 // IMPORTANT! It should be called also for the end iterator - last 620 template <typename TurnInfo> 621 class turns_analyser 622 { 623 typedef typename TurnInfo::point_type turn_point_type; 624 625 static const std::size_t op_id = 0; 626 static const std::size_t other_op_id = 1; 627 628 template <typename TurnPointCSTag, typename PointP, typename PointQ, 629 typename SideStrategy, 630 typename Pi = PointP, typename Pj = PointP, typename Pk = PointP, 631 typename Qi = PointQ, typename Qj = PointQ, typename Qk = PointQ 632 > 633 struct la_side_calculator 634 { la_side_calculatorboost::geometry::detail::relate::linear_areal::turns_analyser::la_side_calculator635 inline la_side_calculator(Pi const& pi, Pj const& pj, Pk const& pk, 636 Qi const& qi, Qj const& qj, Qk const& qk, 637 SideStrategy const& side_strategy) 638 : m_pi(pi), m_pj(pj), m_pk(pk) 639 , m_qi(qi), m_qj(qj), m_qk(qk) 640 , m_side_strategy(side_strategy) 641 {} 642 pk_wrt_p1boost::geometry::detail::relate::linear_areal::turns_analyser::la_side_calculator643 inline int pk_wrt_p1() const { return m_side_strategy.apply(m_pi, m_pj, m_pk); } qk_wrt_p1boost::geometry::detail::relate::linear_areal::turns_analyser::la_side_calculator644 inline int qk_wrt_p1() const { return m_side_strategy.apply(m_pi, m_pj, m_qk); } pk_wrt_q2boost::geometry::detail::relate::linear_areal::turns_analyser::la_side_calculator645 inline int pk_wrt_q2() const { return m_side_strategy.apply(m_qj, m_qk, m_pk); } 646 647 private : 648 Pi const& m_pi; 649 Pj const& m_pj; 650 Pk const& m_pk; 651 Qi const& m_qi; 652 Qj const& m_qj; 653 Qk const& m_qk; 654 655 SideStrategy m_side_strategy; 656 }; 657 658 659 public: turns_analyser()660 turns_analyser() 661 : m_previous_turn_ptr(NULL) 662 , m_previous_operation(overlay::operation_none) 663 , m_boundary_counter(0) 664 , m_interior_detected(false) 665 , m_first_interior_other_id_ptr(NULL) 666 , m_first_from_unknown(false) 667 , m_first_from_unknown_boundary_detected(false) 668 {} 669 670 template <typename Result, 671 typename TurnIt, 672 typename Geometry, 673 typename OtherGeometry, 674 typename BoundaryChecker, 675 typename SideStrategy> apply(Result & res,TurnIt it,Geometry const & geometry,OtherGeometry const & other_geometry,BoundaryChecker const & boundary_checker,SideStrategy const & side_strategy)676 void apply(Result & res, TurnIt it, 677 Geometry const& geometry, 678 OtherGeometry const& other_geometry, 679 BoundaryChecker const& boundary_checker, 680 SideStrategy const& side_strategy) 681 { 682 overlay::operation_type op = it->operations[op_id].operation; 683 684 if ( op != overlay::operation_union 685 && op != overlay::operation_intersection 686 && op != overlay::operation_blocked 687 && op != overlay::operation_continue ) // operation_boundary / operation_boundary_intersection 688 { 689 return; 690 } 691 692 segment_identifier const& seg_id = it->operations[op_id].seg_id; 693 segment_identifier const& other_id = it->operations[other_op_id].seg_id; 694 695 const bool first_in_range = m_seg_watcher.update(seg_id); 696 697 // TODO: should apply() for the post-last ip be called if first_in_range ? 698 // this would unify how last points in ranges are handled 699 // possibly replacing parts of the code below 700 // e.g. for is_multi and m_interior_detected 701 702 // handle possible exit 703 bool fake_enter_detected = false; 704 if ( m_exit_watcher.get_exit_operation() == overlay::operation_union ) 705 { 706 // real exit point - may be multiple 707 // we know that we entered and now we exit 708 if ( ! turn_on_the_same_ip<op_id>(m_exit_watcher.get_exit_turn(), *it, 709 side_strategy.get_equals_point_point_strategy()) ) 710 { 711 m_exit_watcher.reset_detected_exit(); 712 713 update<interior, exterior, '1', TransposeResult>(res); 714 715 // next single geometry 716 if ( first_in_range && m_previous_turn_ptr ) 717 { 718 // NOTE: similar code is in the post-last-ip-apply() 719 segment_identifier const& prev_seg_id = m_previous_turn_ptr->operations[op_id].seg_id; 720 721 bool const prev_back_b = is_endpoint_on_boundary<boundary_back>( 722 range::back(sub_range(geometry, prev_seg_id)), 723 boundary_checker); 724 725 // if there is a boundary on the last point 726 if ( prev_back_b ) 727 { 728 update<boundary, exterior, '0', TransposeResult>(res); 729 } 730 } 731 } 732 // fake exit point, reset state 733 else if ( op == overlay::operation_intersection 734 || op == overlay::operation_continue ) // operation_boundary 735 { 736 m_exit_watcher.reset_detected_exit(); 737 fake_enter_detected = true; 738 } 739 } 740 else if ( m_exit_watcher.get_exit_operation() == overlay::operation_blocked ) 741 { 742 // ignore multiple BLOCKs for this same single geometry1 743 if ( op == overlay::operation_blocked 744 && seg_id.multi_index == m_previous_turn_ptr->operations[op_id].seg_id.multi_index ) 745 { 746 return; 747 } 748 749 if ( ( op == overlay::operation_intersection 750 || op == overlay::operation_continue ) 751 && turn_on_the_same_ip<op_id>(m_exit_watcher.get_exit_turn(), *it, 752 side_strategy.get_equals_point_point_strategy()) ) 753 { 754 fake_enter_detected = true; 755 } 756 757 m_exit_watcher.reset_detected_exit(); 758 } 759 760 if ( BOOST_GEOMETRY_CONDITION( util::is_multi<OtherGeometry>::value ) 761 && m_first_from_unknown ) 762 { 763 // For MultiPolygon many x/u operations may be generated as a first IP 764 // if for all turns x/u was generated and any of the Polygons doesn't contain the LineString 765 // then we know that the LineString is outside 766 // Similar with the u/u turns, if it was the first one it doesn't mean that the 767 // Linestring came from the exterior 768 if ( ( m_previous_operation == overlay::operation_blocked 769 && ( op != overlay::operation_blocked // operation different than block 770 || seg_id.multi_index != m_previous_turn_ptr->operations[op_id].seg_id.multi_index ) ) // or the next single-geometry 771 || ( m_previous_operation == overlay::operation_union 772 && ! turn_on_the_same_ip<op_id>(*m_previous_turn_ptr, *it, 773 side_strategy.get_equals_point_point_strategy()) ) 774 ) 775 { 776 update<interior, exterior, '1', TransposeResult>(res); 777 if ( m_first_from_unknown_boundary_detected ) 778 { 779 update<boundary, exterior, '0', TransposeResult>(res); 780 } 781 782 m_first_from_unknown = false; 783 m_first_from_unknown_boundary_detected = false; 784 } 785 } 786 787 // NOTE: THE WHOLE m_interior_detected HANDLING IS HERE BECAUSE WE CAN'T EFFICIENTLY SORT TURNS (CORRECTLY) 788 // BECAUSE THE SAME IP MAY BE REPRESENTED BY TWO SEGMENTS WITH DIFFERENT DISTANCES 789 // IT WOULD REQUIRE THE CALCULATION OF MAX DISTANCE 790 // TODO: WE COULD GET RID OF THE TEST IF THE DISTANCES WERE NORMALIZED 791 792 // UPDATE: THEY SHOULD BE NORMALIZED NOW 793 794 // TODO: THIS IS POTENTIALLY ERROREOUS! 795 // THIS ALGORITHM DEPENDS ON SOME SPECIFIC SEQUENCE OF OPERATIONS 796 // IT WOULD GIVE WRONG RESULTS E.G. 797 // IN THE CASE OF SELF-TOUCHING POINT WHEN 'i' WOULD BE BEFORE 'u' 798 799 // handle the interior overlap 800 if ( m_interior_detected ) 801 { 802 BOOST_GEOMETRY_ASSERT_MSG(m_previous_turn_ptr, "non-NULL ptr expected"); 803 804 // real interior overlap 805 if ( ! turn_on_the_same_ip<op_id>(*m_previous_turn_ptr, *it, 806 side_strategy.get_equals_point_point_strategy()) ) 807 { 808 update<interior, interior, '1', TransposeResult>(res); 809 m_interior_detected = false; 810 811 // new range detected - reset previous state and check the boundary 812 if ( first_in_range ) 813 { 814 segment_identifier const& prev_seg_id = m_previous_turn_ptr->operations[op_id].seg_id; 815 816 bool const prev_back_b = is_endpoint_on_boundary<boundary_back>( 817 range::back(sub_range(geometry, prev_seg_id)), 818 boundary_checker); 819 820 // if there is a boundary on the last point 821 if ( prev_back_b ) 822 { 823 update<boundary, interior, '0', TransposeResult>(res); 824 } 825 826 // The exit_watcher is reset below 827 // m_exit_watcher.reset(); 828 } 829 } 830 // fake interior overlap 831 else if ( op == overlay::operation_continue ) 832 { 833 m_interior_detected = false; 834 } 835 else if ( op == overlay::operation_union ) 836 { 837 // TODO: this probably is not a good way of handling the interiors/enters 838 // the solution similar to exit_watcher would be more robust 839 // all enters should be kept and handled. 840 // maybe integrate it with the exit_watcher -> enter_exit_watcher 841 if ( m_first_interior_other_id_ptr 842 && m_first_interior_other_id_ptr->multi_index == other_id.multi_index ) 843 { 844 m_interior_detected = false; 845 } 846 } 847 } 848 849 // NOTE: If post-last-ip apply() was called this wouldn't be needed 850 if ( first_in_range ) 851 { 852 m_exit_watcher.reset(); 853 m_boundary_counter = 0; 854 m_first_from_unknown = false; 855 m_first_from_unknown_boundary_detected = false; 856 } 857 858 // i/u, c/u 859 if ( op == overlay::operation_intersection 860 || op == overlay::operation_continue ) // operation_boundary/operation_boundary_intersection 861 { 862 bool const first_point = first_in_range || m_first_from_unknown; 863 bool no_enters_detected = m_exit_watcher.is_outside(); 864 m_exit_watcher.enter(*it); 865 866 if ( op == overlay::operation_intersection ) 867 { 868 if ( m_boundary_counter > 0 && it->operations[op_id].is_collinear ) 869 --m_boundary_counter; 870 871 if ( m_boundary_counter == 0 ) 872 { 873 // interiors overlaps 874 //update<interior, interior, '1', TransposeResult>(res); 875 876 // TODO: think about the implementation of the more robust version 877 // this way only the first enter will be handled 878 if ( !m_interior_detected ) 879 { 880 // don't update now 881 // we might enter a boundary of some other ring on the same IP 882 m_interior_detected = true; 883 m_first_interior_other_id_ptr = boost::addressof(other_id); 884 } 885 } 886 } 887 else // operation_boundary 888 { 889 // don't add to the count for all met boundaries 890 // only if this is the "new" boundary 891 if ( first_point || !it->operations[op_id].is_collinear ) 892 ++m_boundary_counter; 893 894 update<interior, boundary, '1', TransposeResult>(res); 895 } 896 897 bool const this_b 898 = is_ip_on_boundary<boundary_front>(it->point, 899 it->operations[op_id], 900 boundary_checker, 901 seg_id); 902 // going inside on boundary point 903 if ( this_b ) 904 { 905 update<boundary, boundary, '0', TransposeResult>(res); 906 } 907 // going inside on non-boundary point 908 else 909 { 910 update<interior, boundary, '0', TransposeResult>(res); 911 912 // if we didn't enter in the past, we were outside 913 if ( no_enters_detected 914 && ! fake_enter_detected 915 && it->operations[op_id].position != overlay::position_front ) 916 { 917 // TODO: calculate_from_inside() is only needed if the current Linestring is not closed 918 bool const from_inside = first_point 919 && calculate_from_inside(geometry, 920 other_geometry, 921 *it, 922 side_strategy); 923 924 if ( from_inside ) 925 update<interior, interior, '1', TransposeResult>(res); 926 else 927 update<interior, exterior, '1', TransposeResult>(res); 928 929 // if it's the first IP then the first point is outside 930 if ( first_point ) 931 { 932 bool const front_b = is_endpoint_on_boundary<boundary_front>( 933 range::front(sub_range(geometry, seg_id)), 934 boundary_checker); 935 936 // if there is a boundary on the first point 937 if ( front_b ) 938 { 939 if ( from_inside ) 940 update<boundary, interior, '0', TransposeResult>(res); 941 else 942 update<boundary, exterior, '0', TransposeResult>(res); 943 } 944 } 945 } 946 } 947 948 if ( BOOST_GEOMETRY_CONDITION( util::is_multi<OtherGeometry>::value ) ) 949 { 950 m_first_from_unknown = false; 951 m_first_from_unknown_boundary_detected = false; 952 } 953 } 954 // u/u, x/u 955 else if ( op == overlay::operation_union || op == overlay::operation_blocked ) 956 { 957 bool const op_blocked = op == overlay::operation_blocked; 958 bool const no_enters_detected = m_exit_watcher.is_outside() 959 // TODO: is this condition ok? 960 // TODO: move it into the exit_watcher? 961 && m_exit_watcher.get_exit_operation() == overlay::operation_none; 962 963 if ( op == overlay::operation_union ) 964 { 965 if ( m_boundary_counter > 0 && it->operations[op_id].is_collinear ) 966 --m_boundary_counter; 967 } 968 else // overlay::operation_blocked 969 { 970 m_boundary_counter = 0; 971 } 972 973 // we're inside, possibly going out right now 974 if ( ! no_enters_detected ) 975 { 976 if ( op_blocked 977 && it->operations[op_id].position == overlay::position_back ) // ignore spikes! 978 { 979 // check if this is indeed the boundary point 980 // NOTE: is_ip_on_boundary<>() should be called here but the result will be the same 981 if ( is_endpoint_on_boundary<boundary_back>(it->point, boundary_checker) ) 982 { 983 update<boundary, boundary, '0', TransposeResult>(res); 984 } 985 } 986 // union, inside, but no exit -> collinear on self-intersection point 987 // not needed since we're already inside the boundary 988 /*else if ( !exit_detected ) 989 { 990 update<interior, boundary, '0', TransposeResult>(res); 991 }*/ 992 } 993 // we're outside or inside and this is the first turn 994 else 995 { 996 bool const this_b = is_ip_on_boundary<boundary_any>(it->point, 997 it->operations[op_id], 998 boundary_checker, 999 seg_id); 1000 // if current IP is on boundary of the geometry 1001 if ( this_b ) 1002 { 1003 update<boundary, boundary, '0', TransposeResult>(res); 1004 } 1005 // if current IP is not on boundary of the geometry 1006 else 1007 { 1008 update<interior, boundary, '0', TransposeResult>(res); 1009 } 1010 1011 // TODO: very similar code is used in the handling of intersection 1012 if ( it->operations[op_id].position != overlay::position_front ) 1013 { 1014 // TODO: calculate_from_inside() is only needed if the current Linestring is not closed 1015 // NOTE: this is not enough for MultiPolygon and operation_blocked 1016 // For LS/MultiPolygon multiple x/u turns may be generated 1017 // the first checked Polygon may be the one which LS is outside for. 1018 bool const first_point = first_in_range || m_first_from_unknown; 1019 bool const first_from_inside = first_point 1020 && calculate_from_inside(geometry, 1021 other_geometry, 1022 *it, 1023 side_strategy); 1024 if ( first_from_inside ) 1025 { 1026 update<interior, interior, '1', TransposeResult>(res); 1027 1028 // notify the exit_watcher that we started inside 1029 m_exit_watcher.enter(*it); 1030 // and reset unknown flags since we know that we started inside 1031 m_first_from_unknown = false; 1032 m_first_from_unknown_boundary_detected = false; 1033 } 1034 else 1035 { 1036 if ( BOOST_GEOMETRY_CONDITION( util::is_multi<OtherGeometry>::value ) 1037 /*&& ( op == overlay::operation_blocked 1038 || op == overlay::operation_union )*/ ) // if we're here it's u or x 1039 { 1040 m_first_from_unknown = true; 1041 } 1042 else 1043 { 1044 update<interior, exterior, '1', TransposeResult>(res); 1045 } 1046 } 1047 1048 // first IP on the last segment point - this means that the first point is outside or inside 1049 if ( first_point && ( !this_b || op_blocked ) ) 1050 { 1051 bool const front_b = is_endpoint_on_boundary<boundary_front>( 1052 range::front(sub_range(geometry, seg_id)), 1053 boundary_checker); 1054 1055 // if there is a boundary on the first point 1056 if ( front_b ) 1057 { 1058 if ( first_from_inside ) 1059 { 1060 update<boundary, interior, '0', TransposeResult>(res); 1061 } 1062 else 1063 { 1064 if ( BOOST_GEOMETRY_CONDITION( util::is_multi<OtherGeometry>::value ) 1065 /*&& ( op == overlay::operation_blocked 1066 || op == overlay::operation_union )*/ ) // if we're here it's u or x 1067 { 1068 BOOST_GEOMETRY_ASSERT(m_first_from_unknown); 1069 m_first_from_unknown_boundary_detected = true; 1070 } 1071 else 1072 { 1073 update<boundary, exterior, '0', TransposeResult>(res); 1074 } 1075 } 1076 } 1077 } 1078 } 1079 } 1080 1081 // if we're going along a boundary, we exit only if the linestring was collinear 1082 if ( m_boundary_counter == 0 1083 || it->operations[op_id].is_collinear ) 1084 { 1085 // notify the exit watcher about the possible exit 1086 m_exit_watcher.exit(*it); 1087 } 1088 } 1089 1090 // store ref to previously analysed (valid) turn 1091 m_previous_turn_ptr = boost::addressof(*it); 1092 // and previously analysed (valid) operation 1093 m_previous_operation = op; 1094 } 1095 1096 // it == last 1097 template <typename Result, 1098 typename TurnIt, 1099 typename Geometry, 1100 typename OtherGeometry, 1101 typename BoundaryChecker> apply(Result & res,TurnIt first,TurnIt last,Geometry const & geometry,OtherGeometry const &,BoundaryChecker const & boundary_checker)1102 void apply(Result & res, 1103 TurnIt first, TurnIt last, 1104 Geometry const& geometry, 1105 OtherGeometry const& /*other_geometry*/, 1106 BoundaryChecker const& boundary_checker) 1107 { 1108 boost::ignore_unused(first, last); 1109 //BOOST_GEOMETRY_ASSERT( first != last ); 1110 1111 // For MultiPolygon many x/u operations may be generated as a first IP 1112 // if for all turns x/u was generated and any of the Polygons doesn't contain the LineString 1113 // then we know that the LineString is outside 1114 if ( BOOST_GEOMETRY_CONDITION( util::is_multi<OtherGeometry>::value ) 1115 && m_first_from_unknown ) 1116 { 1117 update<interior, exterior, '1', TransposeResult>(res); 1118 if ( m_first_from_unknown_boundary_detected ) 1119 { 1120 update<boundary, exterior, '0', TransposeResult>(res); 1121 } 1122 1123 // done below 1124 //m_first_from_unknown = false; 1125 //m_first_from_unknown_boundary_detected = false; 1126 } 1127 1128 // here, the possible exit is the real one 1129 // we know that we entered and now we exit 1130 if ( /*m_exit_watcher.get_exit_operation() == overlay::operation_union // THIS CHECK IS REDUNDANT 1131 ||*/ m_previous_operation == overlay::operation_union 1132 && !m_interior_detected ) 1133 { 1134 // for sure 1135 update<interior, exterior, '1', TransposeResult>(res); 1136 1137 BOOST_GEOMETRY_ASSERT(first != last); 1138 BOOST_GEOMETRY_ASSERT(m_previous_turn_ptr); 1139 1140 segment_identifier const& prev_seg_id = m_previous_turn_ptr->operations[op_id].seg_id; 1141 1142 bool const prev_back_b = is_endpoint_on_boundary<boundary_back>( 1143 range::back(sub_range(geometry, prev_seg_id)), 1144 boundary_checker); 1145 1146 // if there is a boundary on the last point 1147 if ( prev_back_b ) 1148 { 1149 update<boundary, exterior, '0', TransposeResult>(res); 1150 } 1151 } 1152 // we might enter some Areal and didn't go out, 1153 else if ( m_previous_operation == overlay::operation_intersection 1154 || m_interior_detected ) 1155 { 1156 // just in case 1157 update<interior, interior, '1', TransposeResult>(res); 1158 m_interior_detected = false; 1159 1160 BOOST_GEOMETRY_ASSERT(first != last); 1161 BOOST_GEOMETRY_ASSERT(m_previous_turn_ptr); 1162 1163 segment_identifier const& prev_seg_id = m_previous_turn_ptr->operations[op_id].seg_id; 1164 1165 bool const prev_back_b = is_endpoint_on_boundary<boundary_back>( 1166 range::back(sub_range(geometry, prev_seg_id)), 1167 boundary_checker); 1168 1169 // if there is a boundary on the last point 1170 if ( prev_back_b ) 1171 { 1172 update<boundary, interior, '0', TransposeResult>(res); 1173 } 1174 } 1175 1176 // This condition may be false if the Linestring is lying on the Polygon's collinear spike 1177 // if Polygon's spikes are not handled in get_turns() or relate() (they currently aren't) 1178 //BOOST_GEOMETRY_ASSERT_MSG(m_previous_operation != overlay::operation_continue, 1179 // "Unexpected operation! Probably the error in get_turns(L,A) or relate(L,A)"); 1180 // Currently one c/c turn is generated for the exit 1181 // when a Linestring is going out on a collinear spike 1182 // When a Linestring is going in on a collinear spike 1183 // the turn is not generated for the entry 1184 // So assume it's the former 1185 if ( m_previous_operation == overlay::operation_continue ) 1186 { 1187 update<interior, exterior, '1', TransposeResult>(res); 1188 1189 segment_identifier const& prev_seg_id = m_previous_turn_ptr->operations[op_id].seg_id; 1190 1191 bool const prev_back_b = is_endpoint_on_boundary<boundary_back>( 1192 range::back(sub_range(geometry, prev_seg_id)), 1193 boundary_checker); 1194 1195 // if there is a boundary on the last point 1196 if ( prev_back_b ) 1197 { 1198 update<boundary, exterior, '0', TransposeResult>(res); 1199 } 1200 } 1201 1202 // Reset exit watcher before the analysis of the next Linestring 1203 m_exit_watcher.reset(); 1204 m_boundary_counter = 0; 1205 m_first_from_unknown = false; 1206 m_first_from_unknown_boundary_detected = false; 1207 } 1208 1209 // check if the passed turn's segment of Linear geometry arrived 1210 // from the inside or the outside of the Areal geometry 1211 template <typename Turn, typename SideStrategy> calculate_from_inside(Geometry1 const & geometry1,Geometry2 const & geometry2,Turn const & turn,SideStrategy const & side_strategy)1212 static inline bool calculate_from_inside(Geometry1 const& geometry1, 1213 Geometry2 const& geometry2, 1214 Turn const& turn, 1215 SideStrategy const& side_strategy) 1216 { 1217 typedef typename cs_tag<typename Turn::point_type>::type cs_tag; 1218 1219 if ( turn.operations[op_id].position == overlay::position_front ) 1220 return false; 1221 1222 typename sub_range_return_type<Geometry1 const>::type 1223 range1 = sub_range(geometry1, turn.operations[op_id].seg_id); 1224 1225 typedef detail::normalized_view<Geometry2 const> const range2_type; 1226 typedef typename boost::range_iterator<range2_type>::type range2_iterator; 1227 range2_type range2(sub_range(geometry2, turn.operations[other_op_id].seg_id)); 1228 1229 BOOST_GEOMETRY_ASSERT(boost::size(range1)); 1230 std::size_t const s2 = boost::size(range2); 1231 BOOST_GEOMETRY_ASSERT(s2 > 2); 1232 std::size_t const seg_count2 = s2 - 1; 1233 1234 std::size_t const p_seg_ij = static_cast<std::size_t>(turn.operations[op_id].seg_id.segment_index); 1235 std::size_t const q_seg_ij = static_cast<std::size_t>(turn.operations[other_op_id].seg_id.segment_index); 1236 1237 BOOST_GEOMETRY_ASSERT(p_seg_ij + 1 < boost::size(range1)); 1238 BOOST_GEOMETRY_ASSERT(q_seg_ij + 1 < s2); 1239 1240 point1_type const& pi = range::at(range1, p_seg_ij); 1241 point2_type const& qi = range::at(range2, q_seg_ij); 1242 point2_type const& qj = range::at(range2, q_seg_ij + 1); 1243 point1_type qi_conv; 1244 geometry::convert(qi, qi_conv); 1245 bool const is_ip_qj = equals::equals_point_point(turn.point, qj, side_strategy.get_equals_point_point_strategy()); 1246 // TODO: test this! 1247 // BOOST_GEOMETRY_ASSERT(!equals::equals_point_point(turn.point, pi)); 1248 // BOOST_GEOMETRY_ASSERT(!equals::equals_point_point(turn.point, qi)); 1249 point1_type new_pj; 1250 geometry::convert(turn.point, new_pj); 1251 1252 if ( is_ip_qj ) 1253 { 1254 std::size_t const q_seg_jk = (q_seg_ij + 1) % seg_count2; 1255 // TODO: the following function should return immediately, however the worst case complexity is O(N) 1256 // It would be good to replace it with some O(1) mechanism 1257 range2_iterator qk_it = find_next_non_duplicated(boost::begin(range2), 1258 range::pos(range2, q_seg_jk), 1259 boost::end(range2), 1260 side_strategy.get_equals_point_point_strategy()); 1261 1262 // Will this sequence of points be always correct? 1263 la_side_calculator<cs_tag, point1_type, point2_type, SideStrategy> 1264 side_calc(qi_conv, new_pj, pi, qi, qj, *qk_it, side_strategy); 1265 1266 return calculate_from_inside_sides(side_calc); 1267 } 1268 else 1269 { 1270 point2_type new_qj; 1271 geometry::convert(turn.point, new_qj); 1272 1273 la_side_calculator<cs_tag, point1_type, point2_type, SideStrategy> 1274 side_calc(qi_conv, new_pj, pi, qi, new_qj, qj, side_strategy); 1275 1276 return calculate_from_inside_sides(side_calc); 1277 } 1278 } 1279 1280 template <typename It, typename EqPPStrategy> find_next_non_duplicated(It first,It current,It last,EqPPStrategy const & strategy)1281 static inline It find_next_non_duplicated(It first, It current, It last, 1282 EqPPStrategy const& strategy) 1283 { 1284 BOOST_GEOMETRY_ASSERT( current != last ); 1285 1286 It it = current; 1287 1288 for ( ++it ; it != last ; ++it ) 1289 { 1290 if ( !equals::equals_point_point(*current, *it, strategy) ) 1291 return it; 1292 } 1293 1294 // if not found start from the beginning 1295 for ( it = first ; it != current ; ++it ) 1296 { 1297 if ( !equals::equals_point_point(*current, *it, strategy) ) 1298 return it; 1299 } 1300 1301 return current; 1302 } 1303 1304 // calculate inside or outside based on side_calc 1305 // this is simplified version of a check from equal<> 1306 template <typename SideCalc> calculate_from_inside_sides(SideCalc const & side_calc)1307 static inline bool calculate_from_inside_sides(SideCalc const& side_calc) 1308 { 1309 int const side_pk_p = side_calc.pk_wrt_p1(); 1310 int const side_qk_p = side_calc.qk_wrt_p1(); 1311 // If they turn to same side (not opposite sides) 1312 if (! overlay::base_turn_handler::opposite(side_pk_p, side_qk_p)) 1313 { 1314 int const side_pk_q2 = side_calc.pk_wrt_q2(); 1315 return side_pk_q2 == -1; 1316 } 1317 else 1318 { 1319 return side_pk_p == -1; 1320 } 1321 } 1322 1323 private: 1324 exit_watcher<TurnInfo, op_id> m_exit_watcher; 1325 segment_watcher<same_single> m_seg_watcher; 1326 TurnInfo * m_previous_turn_ptr; 1327 overlay::operation_type m_previous_operation; 1328 unsigned m_boundary_counter; 1329 bool m_interior_detected; 1330 const segment_identifier * m_first_interior_other_id_ptr; 1331 bool m_first_from_unknown; 1332 bool m_first_from_unknown_boundary_detected; 1333 }; 1334 1335 // call analyser.apply() for each turn in range 1336 // IMPORTANT! The analyser is also called for the end iterator - last 1337 template <typename Result, 1338 typename TurnIt, 1339 typename Analyser, 1340 typename Geometry, 1341 typename OtherGeometry, 1342 typename BoundaryChecker, 1343 typename SideStrategy> analyse_each_turnboost::geometry::detail::relate::linear_areal1344 static inline void analyse_each_turn(Result & res, 1345 Analyser & analyser, 1346 TurnIt first, TurnIt last, 1347 Geometry const& geometry, 1348 OtherGeometry const& other_geometry, 1349 BoundaryChecker const& boundary_checker, 1350 SideStrategy const& side_strategy) 1351 { 1352 if ( first == last ) 1353 return; 1354 1355 for ( TurnIt it = first ; it != last ; ++it ) 1356 { 1357 analyser.apply(res, it, 1358 geometry, other_geometry, 1359 boundary_checker, 1360 side_strategy); 1361 1362 if ( BOOST_GEOMETRY_CONDITION( res.interrupt ) ) 1363 return; 1364 } 1365 1366 analyser.apply(res, first, last, 1367 geometry, other_geometry, 1368 boundary_checker); 1369 } 1370 1371 // less comparator comparing multi_index then ring_index 1372 // may be used to sort turns by ring 1373 struct less_ring 1374 { 1375 template <typename Turn> operator ()boost::geometry::detail::relate::linear_areal::less_ring1376 inline bool operator()(Turn const& left, Turn const& right) const 1377 { 1378 return left.operations[1].seg_id.multi_index < right.operations[1].seg_id.multi_index 1379 || ( left.operations[1].seg_id.multi_index == right.operations[1].seg_id.multi_index 1380 && left.operations[1].seg_id.ring_index < right.operations[1].seg_id.ring_index ); 1381 } 1382 }; 1383 1384 // policy/functor checking if a turn's operation is operation_continue 1385 struct has_boundary_intersection 1386 { has_boundary_intersectionboost::geometry::detail::relate::linear_areal::has_boundary_intersection1387 has_boundary_intersection() 1388 : result(false) {} 1389 1390 template <typename Turn> operator ()boost::geometry::detail::relate::linear_areal::has_boundary_intersection1391 inline void operator()(Turn const& turn) 1392 { 1393 if ( turn.operations[1].operation == overlay::operation_continue ) 1394 result = true; 1395 } 1396 1397 bool result; 1398 }; 1399 1400 // iterate through the range and search for the different multi_index or ring_index 1401 // also call fun for each turn in the current range 1402 template <typename It, typename Fun> find_next_ringboost::geometry::detail::relate::linear_areal1403 static inline It find_next_ring(It first, It last, Fun & fun) 1404 { 1405 if ( first == last ) 1406 return last; 1407 1408 signed_size_type const multi_index = first->operations[1].seg_id.multi_index; 1409 signed_size_type const ring_index = first->operations[1].seg_id.ring_index; 1410 1411 fun(*first); 1412 ++first; 1413 1414 for ( ; first != last ; ++first ) 1415 { 1416 if ( multi_index != first->operations[1].seg_id.multi_index 1417 || ring_index != first->operations[1].seg_id.ring_index ) 1418 { 1419 return first; 1420 } 1421 1422 fun(*first); 1423 } 1424 1425 return last; 1426 } 1427 1428 // analyser which called for turns sorted by seg/distance/operation 1429 // checks if the boundary of Areal geometry really got out 1430 // into the exterior of Linear geometry 1431 template <typename TurnInfo> 1432 class areal_boundary_analyser 1433 { 1434 public: areal_boundary_analyser()1435 areal_boundary_analyser() 1436 : is_union_detected(false) 1437 , m_previous_turn_ptr(NULL) 1438 {} 1439 1440 template <typename TurnIt, typename EqPPStrategy> apply(TurnIt,TurnIt it,TurnIt last,EqPPStrategy const & strategy)1441 bool apply(TurnIt /*first*/, TurnIt it, TurnIt last, 1442 EqPPStrategy const& strategy) 1443 { 1444 overlay::operation_type op = it->operations[1].operation; 1445 1446 if ( it != last ) 1447 { 1448 if ( op != overlay::operation_union 1449 && op != overlay::operation_continue ) 1450 { 1451 return true; 1452 } 1453 1454 if ( is_union_detected ) 1455 { 1456 BOOST_GEOMETRY_ASSERT(m_previous_turn_ptr != NULL); 1457 if ( !detail::equals::equals_point_point(it->point, m_previous_turn_ptr->point, strategy) ) 1458 { 1459 // break 1460 return false; 1461 } 1462 else if ( op == overlay::operation_continue ) // operation_boundary 1463 { 1464 is_union_detected = false; 1465 } 1466 } 1467 1468 if ( op == overlay::operation_union ) 1469 { 1470 is_union_detected = true; 1471 m_previous_turn_ptr = boost::addressof(*it); 1472 } 1473 1474 return true; 1475 } 1476 else 1477 { 1478 return false; 1479 } 1480 } 1481 1482 bool is_union_detected; 1483 1484 private: 1485 const TurnInfo * m_previous_turn_ptr; 1486 }; 1487 }; 1488 1489 template <typename Geometry1, typename Geometry2> 1490 struct areal_linear 1491 { 1492 typedef linear_areal<Geometry2, Geometry1, true> linear_areal_type; 1493 1494 static const bool interruption_enabled = linear_areal_type::interruption_enabled; 1495 1496 template <typename Result, typename IntersectionStrategy> applyboost::geometry::detail::relate::areal_linear1497 static inline void apply(Geometry1 const& geometry1, Geometry2 const& geometry2, 1498 Result & result, 1499 IntersectionStrategy const& intersection_strategy) 1500 { 1501 linear_areal_type::apply(geometry2, geometry1, result, intersection_strategy); 1502 } 1503 }; 1504 1505 }} // namespace detail::relate 1506 #endif // DOXYGEN_NO_DETAIL 1507 1508 }} // namespace boost::geometry 1509 1510 #endif // BOOST_GEOMETRY_ALGORITHMS_DETAIL_RELATE_LINEAR_AREAL_HPP 1511