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, 2014, 2015, 2017, 2018, 2019. 6 // Modifications copyright (c) 2013-2019 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_AREAL_AREAL_HPP 15 #define BOOST_GEOMETRY_ALGORITHMS_DETAIL_RELATE_AREAL_AREAL_HPP 16 17 #include <boost/geometry/core/topological_dimension.hpp> 18 19 #include <boost/geometry/util/condition.hpp> 20 #include <boost/geometry/util/range.hpp> 21 22 #include <boost/geometry/algorithms/num_interior_rings.hpp> 23 #include <boost/geometry/algorithms/detail/point_on_border.hpp> 24 #include <boost/geometry/algorithms/detail/sub_range.hpp> 25 #include <boost/geometry/algorithms/detail/single_geometry.hpp> 26 27 #include <boost/geometry/algorithms/detail/relate/point_geometry.hpp> 28 #include <boost/geometry/algorithms/detail/relate/turns.hpp> 29 #include <boost/geometry/algorithms/detail/relate/boundary_checker.hpp> 30 #include <boost/geometry/algorithms/detail/relate/follow_helpers.hpp> 31 32 namespace boost { namespace geometry 33 { 34 35 #ifndef DOXYGEN_NO_DETAIL 36 namespace detail { namespace relate { 37 38 // WARNING! 39 // TODO: In the worst case calling this Pred in a loop for MultiPolygon/MultiPolygon may take O(NM) 40 // Use the rtree in this case! 41 42 // may be used to set EI and EB for an Areal geometry for which no turns were generated 43 template 44 < 45 typename OtherAreal, 46 typename Result, 47 typename PointInArealStrategy, 48 bool TransposeResult 49 > 50 class no_turns_aa_pred 51 { 52 public: no_turns_aa_pred(OtherAreal const & other_areal,Result & res,PointInArealStrategy const & point_in_areal_strategy)53 no_turns_aa_pred(OtherAreal const& other_areal, 54 Result & res, 55 PointInArealStrategy const& point_in_areal_strategy) 56 : m_result(res) 57 , m_point_in_areal_strategy(point_in_areal_strategy) 58 , m_other_areal(other_areal) 59 , m_flags(0) 60 { 61 // check which relations must be analysed 62 63 if ( ! may_update<interior, interior, '2', TransposeResult>(m_result) 64 && ! may_update<boundary, interior, '1', TransposeResult>(m_result) 65 && ! may_update<exterior, interior, '2', TransposeResult>(m_result) ) 66 { 67 m_flags |= 1; 68 } 69 70 if ( ! may_update<interior, exterior, '2', TransposeResult>(m_result) 71 && ! may_update<boundary, exterior, '1', TransposeResult>(m_result) ) 72 { 73 m_flags |= 2; 74 } 75 } 76 77 template <typename Areal> operator ()(Areal const & areal)78 bool operator()(Areal const& areal) 79 { 80 using detail::within::point_in_geometry; 81 82 // if those flags are set nothing will change 83 if ( m_flags == 3 ) 84 { 85 return false; 86 } 87 88 typedef typename geometry::point_type<Areal>::type point_type; 89 point_type pt; 90 bool const ok = boost::geometry::point_on_border(pt, areal); 91 92 // TODO: for now ignore, later throw an exception? 93 if ( !ok ) 94 { 95 return true; 96 } 97 98 // check if the areal is inside the other_areal 99 // TODO: This is O(N) 100 // Run in a loop O(NM) - optimize! 101 int const pig = point_in_geometry(pt, 102 m_other_areal, 103 m_point_in_areal_strategy); 104 //BOOST_GEOMETRY_ASSERT( pig != 0 ); 105 106 // inside 107 if ( pig > 0 ) 108 { 109 update<interior, interior, '2', TransposeResult>(m_result); 110 update<boundary, interior, '1', TransposeResult>(m_result); 111 update<exterior, interior, '2', TransposeResult>(m_result); 112 m_flags |= 1; 113 114 // TODO: OPTIMIZE! 115 // Only the interior rings of other ONE single geometry must be checked 116 // NOT all geometries 117 118 // Check if any interior ring is outside 119 ring_identifier ring_id(0, -1, 0); 120 std::size_t const irings_count = geometry::num_interior_rings(areal); 121 for ( ; static_cast<std::size_t>(ring_id.ring_index) < irings_count ; 122 ++ring_id.ring_index ) 123 { 124 typename detail::sub_range_return_type<Areal const>::type 125 range_ref = detail::sub_range(areal, ring_id); 126 127 if ( boost::empty(range_ref) ) 128 { 129 // TODO: throw exception? 130 continue; // ignore 131 } 132 133 // TODO: O(N) 134 // Optimize! 135 int const hpig = point_in_geometry(range::front(range_ref), 136 m_other_areal, 137 m_point_in_areal_strategy); 138 139 // hole outside 140 if ( hpig < 0 ) 141 { 142 update<interior, exterior, '2', TransposeResult>(m_result); 143 update<boundary, exterior, '1', TransposeResult>(m_result); 144 m_flags |= 2; 145 break; 146 } 147 } 148 } 149 // outside 150 else 151 { 152 update<interior, exterior, '2', TransposeResult>(m_result); 153 update<boundary, exterior, '1', TransposeResult>(m_result); 154 m_flags |= 2; 155 156 // Check if any interior ring is inside 157 ring_identifier ring_id(0, -1, 0); 158 std::size_t const irings_count = geometry::num_interior_rings(areal); 159 for ( ; static_cast<std::size_t>(ring_id.ring_index) < irings_count ; 160 ++ring_id.ring_index ) 161 { 162 typename detail::sub_range_return_type<Areal const>::type 163 range_ref = detail::sub_range(areal, ring_id); 164 165 if ( boost::empty(range_ref) ) 166 { 167 // TODO: throw exception? 168 continue; // ignore 169 } 170 171 // TODO: O(N) 172 // Optimize! 173 int const hpig = point_in_geometry(range::front(range_ref), 174 m_other_areal, 175 m_point_in_areal_strategy); 176 177 // hole inside 178 if ( hpig > 0 ) 179 { 180 update<interior, interior, '2', TransposeResult>(m_result); 181 update<boundary, interior, '1', TransposeResult>(m_result); 182 update<exterior, interior, '2', TransposeResult>(m_result); 183 m_flags |= 1; 184 break; 185 } 186 } 187 } 188 189 return m_flags != 3 && !m_result.interrupt; 190 } 191 192 private: 193 Result & m_result; 194 PointInArealStrategy const& m_point_in_areal_strategy; 195 OtherAreal const& m_other_areal; 196 int m_flags; 197 }; 198 199 // The implementation of an algorithm calculating relate() for A/A 200 template <typename Geometry1, typename Geometry2> 201 struct areal_areal 202 { 203 // check Linear / Areal 204 BOOST_STATIC_ASSERT(topological_dimension<Geometry1>::value == 2 205 && topological_dimension<Geometry2>::value == 2); 206 207 static const bool interruption_enabled = true; 208 209 typedef typename geometry::point_type<Geometry1>::type point1_type; 210 typedef typename geometry::point_type<Geometry2>::type point2_type; 211 212 template <typename Result, typename IntersectionStrategy> applyboost::geometry::detail::relate::areal_areal213 static inline void apply(Geometry1 const& geometry1, Geometry2 const& geometry2, 214 Result & result, 215 IntersectionStrategy const& intersection_strategy) 216 { 217 // TODO: If Areal geometry may have infinite size, change the following line: 218 219 // The result should be FFFFFFFFF 220 relate::set<exterior, exterior, result_dimension<Geometry2>::value>(result);// FFFFFFFFd, d in [1,9] or T 221 222 if ( BOOST_GEOMETRY_CONDITION(result.interrupt) ) 223 return; 224 225 // get and analyse turns 226 typedef typename turns::get_turns 227 < 228 Geometry1, Geometry2 229 >::template turn_info_type<IntersectionStrategy>::type turn_type; 230 std::vector<turn_type> turns; 231 232 interrupt_policy_areal_areal<Result> interrupt_policy(geometry1, geometry2, result); 233 234 turns::get_turns<Geometry1, Geometry2>::apply(turns, geometry1, geometry2, interrupt_policy, intersection_strategy); 235 if ( BOOST_GEOMETRY_CONDITION(result.interrupt) ) 236 return; 237 238 typedef typename IntersectionStrategy::cs_tag cs_tag; 239 240 typedef typename IntersectionStrategy::template point_in_geometry_strategy 241 < 242 Geometry1, Geometry2 243 >::type point_in_areal_strategy12_type; 244 point_in_areal_strategy12_type point_in_areal_strategy12 245 = intersection_strategy.template get_point_in_geometry_strategy<Geometry1, Geometry2>(); 246 typedef typename IntersectionStrategy::template point_in_geometry_strategy 247 < 248 Geometry2, Geometry1 249 >::type point_in_areal_strategy21_type; 250 point_in_areal_strategy21_type point_in_areal_strategy21 251 = intersection_strategy.template get_point_in_geometry_strategy<Geometry2, Geometry1>(); 252 253 no_turns_aa_pred<Geometry2, Result, point_in_areal_strategy12_type, false> 254 pred1(geometry2, result, point_in_areal_strategy12); 255 for_each_disjoint_geometry_if<0, Geometry1>::apply(turns.begin(), turns.end(), geometry1, pred1); 256 if ( BOOST_GEOMETRY_CONDITION(result.interrupt) ) 257 return; 258 259 no_turns_aa_pred<Geometry1, Result, point_in_areal_strategy21_type, true> 260 pred2(geometry1, result, point_in_areal_strategy21); 261 for_each_disjoint_geometry_if<1, Geometry2>::apply(turns.begin(), turns.end(), geometry2, pred2); 262 if ( BOOST_GEOMETRY_CONDITION(result.interrupt) ) 263 return; 264 265 if ( turns.empty() ) 266 return; 267 268 if ( may_update<interior, interior, '2'>(result) 269 || may_update<interior, exterior, '2'>(result) 270 || may_update<boundary, interior, '1'>(result) 271 || may_update<boundary, exterior, '1'>(result) 272 || may_update<exterior, interior, '2'>(result) ) 273 { 274 // sort turns 275 typedef turns::less<0, turns::less_op_areal_areal<0>, cs_tag> less; 276 std::sort(turns.begin(), turns.end(), less()); 277 278 /*if ( may_update<interior, exterior, '2'>(result) 279 || may_update<boundary, exterior, '1'>(result) 280 || may_update<boundary, interior, '1'>(result) 281 || may_update<exterior, interior, '2'>(result) )*/ 282 { 283 // analyse sorted turns 284 turns_analyser<turn_type, 0> analyser; 285 analyse_each_turn(result, analyser, turns.begin(), turns.end(), 286 point_in_areal_strategy12.get_equals_point_point_strategy()); 287 288 if ( BOOST_GEOMETRY_CONDITION(result.interrupt) ) 289 return; 290 } 291 292 if ( may_update<interior, interior, '2'>(result) 293 || may_update<interior, exterior, '2'>(result) 294 || may_update<boundary, interior, '1'>(result) 295 || may_update<boundary, exterior, '1'>(result) 296 || may_update<exterior, interior, '2'>(result) ) 297 { 298 // analyse rings for which turns were not generated 299 // or only i/i or u/u was generated 300 uncertain_rings_analyser<0, Result, Geometry1, Geometry2, point_in_areal_strategy12_type> 301 rings_analyser(result, geometry1, geometry2, point_in_areal_strategy12); 302 analyse_uncertain_rings<0>::apply(rings_analyser, turns.begin(), turns.end()); 303 304 if ( BOOST_GEOMETRY_CONDITION(result.interrupt) ) 305 return; 306 } 307 } 308 309 if ( may_update<interior, interior, '2', true>(result) 310 || may_update<interior, exterior, '2', true>(result) 311 || may_update<boundary, interior, '1', true>(result) 312 || may_update<boundary, exterior, '1', true>(result) 313 || may_update<exterior, interior, '2', true>(result) ) 314 { 315 // sort turns 316 typedef turns::less<1, turns::less_op_areal_areal<1>, cs_tag> less; 317 std::sort(turns.begin(), turns.end(), less()); 318 319 /*if ( may_update<interior, exterior, '2', true>(result) 320 || may_update<boundary, exterior, '1', true>(result) 321 || may_update<boundary, interior, '1', true>(result) 322 || may_update<exterior, interior, '2', true>(result) )*/ 323 { 324 // analyse sorted turns 325 turns_analyser<turn_type, 1> analyser; 326 analyse_each_turn(result, analyser, turns.begin(), turns.end(), 327 point_in_areal_strategy21.get_equals_point_point_strategy()); 328 329 if ( BOOST_GEOMETRY_CONDITION(result.interrupt) ) 330 return; 331 } 332 333 if ( may_update<interior, interior, '2', true>(result) 334 || may_update<interior, exterior, '2', true>(result) 335 || may_update<boundary, interior, '1', true>(result) 336 || may_update<boundary, exterior, '1', true>(result) 337 || may_update<exterior, interior, '2', true>(result) ) 338 { 339 // analyse rings for which turns were not generated 340 // or only i/i or u/u was generated 341 uncertain_rings_analyser<1, Result, Geometry2, Geometry1, point_in_areal_strategy21_type> 342 rings_analyser(result, geometry2, geometry1, point_in_areal_strategy21); 343 analyse_uncertain_rings<1>::apply(rings_analyser, turns.begin(), turns.end()); 344 345 //if ( result.interrupt ) 346 // return; 347 } 348 } 349 } 350 351 // interrupt policy which may be passed to get_turns to interrupt the analysis 352 // based on the info in the passed result/mask 353 template <typename Result> 354 class interrupt_policy_areal_areal 355 { 356 public: 357 static bool const enabled = true; 358 interrupt_policy_areal_areal(Geometry1 const & geometry1,Geometry2 const & geometry2,Result & result)359 interrupt_policy_areal_areal(Geometry1 const& geometry1, 360 Geometry2 const& geometry2, 361 Result & result) 362 : m_result(result) 363 , m_geometry1(geometry1) 364 , m_geometry2(geometry2) 365 {} 366 367 template <typename Range> apply(Range const & turns)368 inline bool apply(Range const& turns) 369 { 370 typedef typename boost::range_iterator<Range const>::type iterator; 371 372 for (iterator it = boost::begin(turns) ; it != boost::end(turns) ; ++it) 373 { 374 per_turn<0>(*it); 375 per_turn<1>(*it); 376 } 377 378 return m_result.interrupt; 379 } 380 381 private: 382 template <std::size_t OpId, typename Turn> per_turn(Turn const & turn)383 inline void per_turn(Turn const& turn) 384 { 385 //static const std::size_t other_op_id = (OpId + 1) % 2; 386 static const bool transpose_result = OpId != 0; 387 388 overlay::operation_type const op = turn.operations[OpId].operation; 389 390 if ( op == overlay::operation_union ) 391 { 392 // ignore u/u 393 /*if ( turn.operations[other_op_id].operation != overlay::operation_union ) 394 { 395 update<interior, exterior, '2', transpose_result>(m_result); 396 update<boundary, exterior, '1', transpose_result>(m_result); 397 }*/ 398 399 update<boundary, boundary, '0', transpose_result>(m_result); 400 } 401 else if ( op == overlay::operation_intersection ) 402 { 403 // ignore i/i 404 /*if ( turn.operations[other_op_id].operation != overlay::operation_intersection ) 405 { 406 // not correct e.g. for G1 touching G2 in a point where a hole is touching the exterior ring 407 // in this case 2 turns i/... and u/u will be generated for this IP 408 //update<interior, interior, '2', transpose_result>(m_result); 409 410 //update<boundary, interior, '1', transpose_result>(m_result); 411 }*/ 412 413 update<boundary, boundary, '0', transpose_result>(m_result); 414 } 415 else if ( op == overlay::operation_continue ) 416 { 417 update<boundary, boundary, '1', transpose_result>(m_result); 418 update<interior, interior, '2', transpose_result>(m_result); 419 } 420 else if ( op == overlay::operation_blocked ) 421 { 422 update<boundary, boundary, '1', transpose_result>(m_result); 423 update<interior, exterior, '2', transpose_result>(m_result); 424 } 425 } 426 427 Result & m_result; 428 Geometry1 const& m_geometry1; 429 Geometry2 const& m_geometry2; 430 }; 431 432 // This analyser should be used like Input or SinglePass Iterator 433 // IMPORTANT! It should be called also for the end iterator - last 434 template <typename TurnInfo, std::size_t OpId> 435 class turns_analyser 436 { 437 typedef typename TurnInfo::point_type turn_point_type; 438 439 static const std::size_t op_id = OpId; 440 static const std::size_t other_op_id = (OpId + 1) % 2; 441 static const bool transpose_result = OpId != 0; 442 443 public: turns_analyser()444 turns_analyser() 445 : m_previous_turn_ptr(0) 446 , m_previous_operation(overlay::operation_none) 447 , m_enter_detected(false) 448 , m_exit_detected(false) 449 {} 450 451 template <typename Result, typename TurnIt, typename EqPPStrategy> apply(Result & result,TurnIt it,EqPPStrategy const & strategy)452 void apply(Result & result, TurnIt it, EqPPStrategy const& strategy) 453 { 454 //BOOST_GEOMETRY_ASSERT( it != last ); 455 456 overlay::operation_type const op = it->operations[op_id].operation; 457 458 if ( op != overlay::operation_union 459 && op != overlay::operation_intersection 460 && op != overlay::operation_blocked 461 && op != overlay::operation_continue ) 462 { 463 return; 464 } 465 466 segment_identifier const& seg_id = it->operations[op_id].seg_id; 467 //segment_identifier const& other_id = it->operations[other_op_id].seg_id; 468 469 const bool first_in_range = m_seg_watcher.update(seg_id); 470 471 if ( m_previous_turn_ptr ) 472 { 473 if ( m_exit_detected /*m_previous_operation == overlay::operation_union*/ ) 474 { 475 // real exit point - may be multiple 476 if ( first_in_range 477 || ! turn_on_the_same_ip<op_id>(*m_previous_turn_ptr, *it, strategy) ) 478 { 479 update_exit(result); 480 m_exit_detected = false; 481 } 482 // fake exit point, reset state 483 else if ( op != overlay::operation_union ) 484 { 485 m_exit_detected = false; 486 } 487 } 488 /*else*/ 489 if ( m_enter_detected /*m_previous_operation == overlay::operation_intersection*/ ) 490 { 491 // real entry point 492 if ( first_in_range 493 || ! turn_on_the_same_ip<op_id>(*m_previous_turn_ptr, *it, strategy) ) 494 { 495 update_enter(result); 496 m_enter_detected = false; 497 } 498 // fake entry point, reset state 499 else if ( op != overlay::operation_intersection ) 500 { 501 m_enter_detected = false; 502 } 503 } 504 } 505 506 if ( op == overlay::operation_union ) 507 { 508 // already set in interrupt policy 509 //update<boundary, boundary, '0', transpose_result>(m_result); 510 511 // ignore u/u 512 //if ( it->operations[other_op_id].operation != overlay::operation_union ) 513 { 514 m_exit_detected = true; 515 } 516 } 517 else if ( op == overlay::operation_intersection ) 518 { 519 // ignore i/i 520 if ( it->operations[other_op_id].operation != overlay::operation_intersection ) 521 { 522 // this was set in the interrupt policy but it was wrong 523 // also here it's wrong since it may be a fake entry point 524 //update<interior, interior, '2', transpose_result>(result); 525 526 // already set in interrupt policy 527 //update<boundary, boundary, '0', transpose_result>(result); 528 m_enter_detected = true; 529 } 530 } 531 else if ( op == overlay::operation_blocked ) 532 { 533 // already set in interrupt policy 534 } 535 else // if ( op == overlay::operation_continue ) 536 { 537 // already set in interrupt policy 538 } 539 540 // store ref to previously analysed (valid) turn 541 m_previous_turn_ptr = boost::addressof(*it); 542 // and previously analysed (valid) operation 543 m_previous_operation = op; 544 } 545 546 // it == last 547 template <typename Result> apply(Result & result)548 void apply(Result & result) 549 { 550 //BOOST_GEOMETRY_ASSERT( first != last ); 551 552 if ( m_exit_detected /*m_previous_operation == overlay::operation_union*/ ) 553 { 554 update_exit(result); 555 m_exit_detected = false; 556 } 557 558 if ( m_enter_detected /*m_previous_operation == overlay::operation_intersection*/ ) 559 { 560 update_enter(result); 561 m_enter_detected = false; 562 } 563 } 564 565 template <typename Result> update_exit(Result & result)566 static inline void update_exit(Result & result) 567 { 568 update<interior, exterior, '2', transpose_result>(result); 569 update<boundary, exterior, '1', transpose_result>(result); 570 } 571 572 template <typename Result> update_enter(Result & result)573 static inline void update_enter(Result & result) 574 { 575 update<interior, interior, '2', transpose_result>(result); 576 update<boundary, interior, '1', transpose_result>(result); 577 update<exterior, interior, '2', transpose_result>(result); 578 } 579 580 private: 581 segment_watcher<same_ring> m_seg_watcher; 582 TurnInfo * m_previous_turn_ptr; 583 overlay::operation_type m_previous_operation; 584 bool m_enter_detected; 585 bool m_exit_detected; 586 }; 587 588 // call analyser.apply() for each turn in range 589 // IMPORTANT! The analyser is also called for the end iterator - last 590 template 591 < 592 typename Result, 593 typename Analyser, 594 typename TurnIt, 595 typename EqPPStrategy 596 > analyse_each_turnboost::geometry::detail::relate::areal_areal597 static inline void analyse_each_turn(Result & res, 598 Analyser & analyser, 599 TurnIt first, TurnIt last, 600 EqPPStrategy const& strategy) 601 { 602 if ( first == last ) 603 return; 604 605 for ( TurnIt it = first ; it != last ; ++it ) 606 { 607 analyser.apply(res, it, strategy); 608 609 if ( BOOST_GEOMETRY_CONDITION(res.interrupt) ) 610 return; 611 } 612 613 analyser.apply(res); 614 } 615 616 template 617 < 618 std::size_t OpId, 619 typename Result, 620 typename Geometry, 621 typename OtherGeometry, 622 typename PointInArealStrategy 623 > 624 class uncertain_rings_analyser 625 { 626 static const bool transpose_result = OpId != 0; 627 static const int other_id = (OpId + 1) % 2; 628 629 public: uncertain_rings_analyser(Result & result,Geometry const & geom,OtherGeometry const & other_geom,PointInArealStrategy const & point_in_areal_strategy)630 inline uncertain_rings_analyser(Result & result, 631 Geometry const& geom, 632 OtherGeometry const& other_geom, 633 PointInArealStrategy const& point_in_areal_strategy) 634 : geometry(geom) 635 , other_geometry(other_geom) 636 , interrupt(result.interrupt) // just in case, could be false as well 637 , m_result(result) 638 , m_point_in_areal_strategy(point_in_areal_strategy) 639 , m_flags(0) 640 { 641 // check which relations must be analysed 642 // NOTE: 1 and 4 could probably be connected 643 644 if ( ! may_update<interior, interior, '2', transpose_result>(m_result) ) 645 { 646 m_flags |= 1; 647 } 648 649 if ( ! may_update<interior, exterior, '2', transpose_result>(m_result) 650 && ! may_update<boundary, exterior, '1', transpose_result>(m_result) ) 651 { 652 m_flags |= 2; 653 } 654 655 if ( ! may_update<boundary, interior, '1', transpose_result>(m_result) 656 && ! may_update<exterior, interior, '2', transpose_result>(m_result) ) 657 { 658 m_flags |= 4; 659 } 660 } 661 no_turns(segment_identifier const & seg_id)662 inline void no_turns(segment_identifier const& seg_id) 663 { 664 // if those flags are set nothing will change 665 if ( m_flags == 7 ) 666 { 667 return; 668 } 669 670 typename detail::sub_range_return_type<Geometry const>::type 671 range_ref = detail::sub_range(geometry, seg_id); 672 673 if ( boost::empty(range_ref) ) 674 { 675 // TODO: throw an exception? 676 return; // ignore 677 } 678 679 // TODO: possible optimization 680 // if the range is an interior ring we may use other IPs generated for this single geometry 681 // to know which other single geometries should be checked 682 683 // TODO: optimize! e.g. use spatial index 684 // O(N) - running it in a loop gives O(NM) 685 using detail::within::point_in_geometry; 686 int const pig = point_in_geometry(range::front(range_ref), 687 other_geometry, 688 m_point_in_areal_strategy); 689 690 //BOOST_GEOMETRY_ASSERT(pig != 0); 691 if ( pig > 0 ) 692 { 693 update<interior, interior, '2', transpose_result>(m_result); 694 m_flags |= 1; 695 696 update<boundary, interior, '1', transpose_result>(m_result); 697 update<exterior, interior, '2', transpose_result>(m_result); 698 m_flags |= 4; 699 } 700 else 701 { 702 update<boundary, exterior, '1', transpose_result>(m_result); 703 update<interior, exterior, '2', transpose_result>(m_result); 704 m_flags |= 2; 705 } 706 707 // TODO: break if all things are set 708 // also some of them could be checked outside, before the analysis 709 // In this case we shouldn't relay just on dummy flags 710 // Flags should be initialized with proper values 711 // or the result should be checked directly 712 // THIS IS ALSO TRUE FOR OTHER ANALYSERS! in L/L and L/A 713 714 interrupt = m_flags == 7 || m_result.interrupt; 715 } 716 717 template <typename TurnIt> turns(TurnIt first,TurnIt last)718 inline void turns(TurnIt first, TurnIt last) 719 { 720 // if those flags are set nothing will change 721 if ( (m_flags & 6) == 6 ) 722 { 723 return; 724 } 725 726 bool found_ii = false; 727 bool found_uu = false; 728 729 for ( TurnIt it = first ; it != last ; ++it ) 730 { 731 if ( it->operations[0].operation == overlay::operation_intersection 732 && it->operations[1].operation == overlay::operation_intersection ) 733 { 734 found_ii = true; 735 } 736 else if ( it->operations[0].operation == overlay::operation_union 737 && it->operations[1].operation == overlay::operation_union ) 738 { 739 found_uu = true; 740 } 741 else // ignore 742 { 743 return; // don't interrupt 744 } 745 } 746 747 // only i/i was generated for this ring 748 if ( found_ii ) 749 { 750 update<interior, interior, '2', transpose_result>(m_result); 751 m_flags |= 1; 752 753 //update<boundary, boundary, '0', transpose_result>(m_result); 754 755 update<boundary, interior, '1', transpose_result>(m_result); 756 update<exterior, interior, '2', transpose_result>(m_result); 757 m_flags |= 4; 758 } 759 760 // only u/u was generated for this ring 761 if ( found_uu ) 762 { 763 update<boundary, exterior, '1', transpose_result>(m_result); 764 update<interior, exterior, '2', transpose_result>(m_result); 765 m_flags |= 2; 766 } 767 768 interrupt = m_flags == 7 || m_result.interrupt; // interrupt if the result won't be changed in the future 769 } 770 771 Geometry const& geometry; 772 OtherGeometry const& other_geometry; 773 bool interrupt; 774 775 private: 776 Result & m_result; 777 PointInArealStrategy const& m_point_in_areal_strategy; 778 int m_flags; 779 }; 780 781 template <std::size_t OpId> 782 class analyse_uncertain_rings 783 { 784 public: 785 template <typename Analyser, typename TurnIt> apply(Analyser & analyser,TurnIt first,TurnIt last)786 static inline void apply(Analyser & analyser, TurnIt first, TurnIt last) 787 { 788 if ( first == last ) 789 return; 790 791 for_preceding_rings(analyser, *first); 792 //analyser.per_turn(*first); 793 794 TurnIt prev = first; 795 for ( ++first ; first != last ; ++first, ++prev ) 796 { 797 // same multi 798 if ( prev->operations[OpId].seg_id.multi_index 799 == first->operations[OpId].seg_id.multi_index ) 800 { 801 // same ring 802 if ( prev->operations[OpId].seg_id.ring_index 803 == first->operations[OpId].seg_id.ring_index ) 804 { 805 //analyser.per_turn(*first); 806 } 807 // same multi, next ring 808 else 809 { 810 //analyser.end_ring(*prev); 811 analyser.turns(prev, first); 812 813 //if ( prev->operations[OpId].seg_id.ring_index + 1 814 // < first->operations[OpId].seg_id.ring_index) 815 { 816 for_no_turns_rings(analyser, 817 *first, 818 prev->operations[OpId].seg_id.ring_index + 1, 819 first->operations[OpId].seg_id.ring_index); 820 } 821 822 //analyser.per_turn(*first); 823 } 824 } 825 // next multi 826 else 827 { 828 //analyser.end_ring(*prev); 829 analyser.turns(prev, first); 830 for_following_rings(analyser, *prev); 831 for_preceding_rings(analyser, *first); 832 //analyser.per_turn(*first); 833 } 834 835 if ( analyser.interrupt ) 836 { 837 return; 838 } 839 } 840 841 //analyser.end_ring(*prev); 842 analyser.turns(prev, first); // first == last 843 for_following_rings(analyser, *prev); 844 } 845 846 private: 847 template <typename Analyser, typename Turn> for_preceding_rings(Analyser & analyser,Turn const & turn)848 static inline void for_preceding_rings(Analyser & analyser, Turn const& turn) 849 { 850 segment_identifier const& seg_id = turn.operations[OpId].seg_id; 851 852 for_no_turns_rings(analyser, turn, -1, seg_id.ring_index); 853 } 854 855 template <typename Analyser, typename Turn> for_following_rings(Analyser & analyser,Turn const & turn)856 static inline void for_following_rings(Analyser & analyser, Turn const& turn) 857 { 858 segment_identifier const& seg_id = turn.operations[OpId].seg_id; 859 860 signed_size_type 861 count = boost::numeric_cast<signed_size_type>( 862 geometry::num_interior_rings( 863 detail::single_geometry(analyser.geometry, seg_id))); 864 865 for_no_turns_rings(analyser, turn, seg_id.ring_index + 1, count); 866 } 867 868 template <typename Analyser, typename Turn> for_no_turns_rings(Analyser & analyser,Turn const & turn,signed_size_type first,signed_size_type last)869 static inline void for_no_turns_rings(Analyser & analyser, 870 Turn const& turn, 871 signed_size_type first, 872 signed_size_type last) 873 { 874 segment_identifier seg_id = turn.operations[OpId].seg_id; 875 876 for ( seg_id.ring_index = first ; seg_id.ring_index < last ; ++seg_id.ring_index ) 877 { 878 analyser.no_turns(seg_id); 879 } 880 } 881 }; 882 }; 883 884 }} // namespace detail::relate 885 #endif // DOXYGEN_NO_DETAIL 886 887 }} // namespace boost::geometry 888 889 #endif // BOOST_GEOMETRY_ALGORITHMS_DETAIL_RELATE_AREAL_AREAL_HPP 890