1 // Boost.Geometry (aka GGL, Generic Geometry Library) 2 3 // Copyright (c) 2007-2014 Barend Gehrels, Amsterdam, the Netherlands. 4 // Copyright (c) 2013-2017 Adam Wulkiewicz, Lodz, Poland. 5 6 // This file was modified by Oracle on 2014, 2016, 2017, 2018, 2019. 7 // Modifications copyright (c) 2014-2019, Oracle and/or its affiliates. 8 9 // Contributed and/or modified by Menelaos Karavelas, on behalf of Oracle 10 // Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle 11 12 // Use, modification and distribution is subject to the Boost Software License, 13 // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at 14 // http://www.boost.org/LICENSE_1_0.txt) 15 16 #ifndef BOOST_GEOMETRY_STRATEGIES_CARTESIAN_INTERSECTION_HPP 17 #define BOOST_GEOMETRY_STRATEGIES_CARTESIAN_INTERSECTION_HPP 18 19 #include <algorithm> 20 21 #include <boost/geometry/core/exception.hpp> 22 23 #include <boost/geometry/geometries/concepts/point_concept.hpp> 24 #include <boost/geometry/geometries/concepts/segment_concept.hpp> 25 26 #include <boost/geometry/arithmetic/determinant.hpp> 27 #include <boost/geometry/algorithms/detail/assign_values.hpp> 28 #include <boost/geometry/algorithms/detail/assign_indexed_point.hpp> 29 #include <boost/geometry/algorithms/detail/equals/point_point.hpp> 30 #include <boost/geometry/algorithms/detail/recalculate.hpp> 31 32 #include <boost/geometry/util/math.hpp> 33 #include <boost/geometry/util/promote_integral.hpp> 34 #include <boost/geometry/util/select_calculation_type.hpp> 35 36 #include <boost/geometry/strategies/cartesian/area.hpp> 37 #include <boost/geometry/strategies/cartesian/disjoint_box_box.hpp> 38 #include <boost/geometry/strategies/cartesian/disjoint_segment_box.hpp> 39 #include <boost/geometry/strategies/cartesian/distance_pythagoras.hpp> 40 #include <boost/geometry/strategies/cartesian/envelope.hpp> 41 #include <boost/geometry/strategies/cartesian/expand_box.hpp> 42 #include <boost/geometry/strategies/cartesian/expand_segment.hpp> 43 #include <boost/geometry/strategies/cartesian/point_in_point.hpp> 44 #include <boost/geometry/strategies/cartesian/point_in_poly_winding.hpp> 45 #include <boost/geometry/strategies/cartesian/side_by_triangle.hpp> 46 #include <boost/geometry/strategies/covered_by.hpp> 47 #include <boost/geometry/strategies/intersection.hpp> 48 #include <boost/geometry/strategies/intersection_result.hpp> 49 #include <boost/geometry/strategies/side.hpp> 50 #include <boost/geometry/strategies/side_info.hpp> 51 #include <boost/geometry/strategies/within.hpp> 52 53 #include <boost/geometry/policies/robustness/rescale_policy_tags.hpp> 54 #include <boost/geometry/policies/robustness/robust_point_type.hpp> 55 56 #if defined(BOOST_GEOMETRY_DEBUG_ROBUSTNESS) 57 # include <boost/geometry/io/wkt/write.hpp> 58 #endif 59 60 61 namespace boost { namespace geometry 62 { 63 64 65 namespace strategy { namespace intersection 66 { 67 68 69 /*! 70 \see http://mathworld.wolfram.com/Line-LineIntersection.html 71 */ 72 template 73 < 74 typename CalculationType = void 75 > 76 struct cartesian_segments 77 { 78 typedef cartesian_tag cs_tag; 79 80 typedef side::side_by_triangle<CalculationType> side_strategy_type; 81 get_side_strategyboost::geometry::strategy::intersection::cartesian_segments82 static inline side_strategy_type get_side_strategy() 83 { 84 return side_strategy_type(); 85 } 86 87 template <typename Geometry1, typename Geometry2> 88 struct point_in_geometry_strategy 89 { 90 typedef strategy::within::cartesian_winding 91 < 92 typename point_type<Geometry1>::type, 93 typename point_type<Geometry2>::type, 94 CalculationType 95 > type; 96 }; 97 98 template <typename Geometry1, typename Geometry2> 99 static inline typename point_in_geometry_strategy<Geometry1, Geometry2>::type get_point_in_geometry_strategyboost::geometry::strategy::intersection::cartesian_segments100 get_point_in_geometry_strategy() 101 { 102 typedef typename point_in_geometry_strategy 103 < 104 Geometry1, Geometry2 105 >::type strategy_type; 106 return strategy_type(); 107 } 108 109 template <typename Geometry> 110 struct area_strategy 111 { 112 typedef area::cartesian 113 < 114 CalculationType 115 > type; 116 }; 117 118 template <typename Geometry> get_area_strategyboost::geometry::strategy::intersection::cartesian_segments119 static inline typename area_strategy<Geometry>::type get_area_strategy() 120 { 121 typedef typename area_strategy<Geometry>::type strategy_type; 122 return strategy_type(); 123 } 124 125 template <typename Geometry> 126 struct distance_strategy 127 { 128 typedef distance::pythagoras 129 < 130 CalculationType 131 > type; 132 }; 133 134 template <typename Geometry> get_distance_strategyboost::geometry::strategy::intersection::cartesian_segments135 static inline typename distance_strategy<Geometry>::type get_distance_strategy() 136 { 137 typedef typename distance_strategy<Geometry>::type strategy_type; 138 return strategy_type(); 139 } 140 141 typedef envelope::cartesian<CalculationType> envelope_strategy_type; 142 get_envelope_strategyboost::geometry::strategy::intersection::cartesian_segments143 static inline envelope_strategy_type get_envelope_strategy() 144 { 145 return envelope_strategy_type(); 146 } 147 148 typedef expand::cartesian_segment expand_strategy_type; 149 get_expand_strategyboost::geometry::strategy::intersection::cartesian_segments150 static inline expand_strategy_type get_expand_strategy() 151 { 152 return expand_strategy_type(); 153 } 154 155 typedef within::cartesian_point_point point_in_point_strategy_type; 156 get_point_in_point_strategyboost::geometry::strategy::intersection::cartesian_segments157 static inline point_in_point_strategy_type get_point_in_point_strategy() 158 { 159 return point_in_point_strategy_type(); 160 } 161 162 typedef within::cartesian_point_point equals_point_point_strategy_type; 163 get_equals_point_point_strategyboost::geometry::strategy::intersection::cartesian_segments164 static inline equals_point_point_strategy_type get_equals_point_point_strategy() 165 { 166 return equals_point_point_strategy_type(); 167 } 168 169 typedef disjoint::cartesian_box_box disjoint_box_box_strategy_type; 170 get_disjoint_box_box_strategyboost::geometry::strategy::intersection::cartesian_segments171 static inline disjoint_box_box_strategy_type get_disjoint_box_box_strategy() 172 { 173 return disjoint_box_box_strategy_type(); 174 } 175 176 typedef disjoint::segment_box disjoint_segment_box_strategy_type; 177 get_disjoint_segment_box_strategyboost::geometry::strategy::intersection::cartesian_segments178 static inline disjoint_segment_box_strategy_type get_disjoint_segment_box_strategy() 179 { 180 return disjoint_segment_box_strategy_type(); 181 } 182 183 typedef covered_by::cartesian_point_box disjoint_point_box_strategy_type; 184 typedef covered_by::cartesian_point_box covered_by_point_box_strategy_type; 185 typedef within::cartesian_point_box within_point_box_strategy_type; 186 typedef envelope::cartesian_box envelope_box_strategy_type; 187 typedef expand::cartesian_box expand_box_strategy_type; 188 189 template <typename CoordinateType, typename SegmentRatio> 190 struct segment_intersection_info 191 { 192 private : 193 typedef typename select_most_precise 194 < 195 CoordinateType, double 196 >::type promoted_type; 197 comparable_length_aboost::geometry::strategy::intersection::cartesian_segments::segment_intersection_info198 promoted_type comparable_length_a() const 199 { 200 return dx_a * dx_a + dy_a * dy_a; 201 } 202 comparable_length_bboost::geometry::strategy::intersection::cartesian_segments::segment_intersection_info203 promoted_type comparable_length_b() const 204 { 205 return dx_b * dx_b + dy_b * dy_b; 206 } 207 208 template <typename Point, typename Segment1, typename Segment2> assign_aboost::geometry::strategy::intersection::cartesian_segments::segment_intersection_info209 void assign_a(Point& point, Segment1 const& a, Segment2 const& ) const 210 { 211 assign(point, a, dx_a, dy_a, robust_ra); 212 } 213 template <typename Point, typename Segment1, typename Segment2> assign_bboost::geometry::strategy::intersection::cartesian_segments::segment_intersection_info214 void assign_b(Point& point, Segment1 const& , Segment2 const& b) const 215 { 216 assign(point, b, dx_b, dy_b, robust_rb); 217 } 218 219 template <typename Point, typename Segment> assignboost::geometry::strategy::intersection::cartesian_segments::segment_intersection_info220 void assign(Point& point, Segment const& segment, CoordinateType const& dx, CoordinateType const& dy, SegmentRatio const& ratio) const 221 { 222 // Calculate the intersection point based on segment_ratio 223 // Up to now, division was postponed. Here we divide using numerator/ 224 // denominator. In case of integer this results in an integer 225 // division. 226 BOOST_GEOMETRY_ASSERT(ratio.denominator() != 0); 227 228 typedef typename promote_integral<CoordinateType>::type calc_type; 229 230 calc_type const numerator 231 = boost::numeric_cast<calc_type>(ratio.numerator()); 232 calc_type const denominator 233 = boost::numeric_cast<calc_type>(ratio.denominator()); 234 calc_type const dx_calc = boost::numeric_cast<calc_type>(dx); 235 calc_type const dy_calc = boost::numeric_cast<calc_type>(dy); 236 237 set<0>(point, get<0, 0>(segment) 238 + boost::numeric_cast<CoordinateType>(numerator * dx_calc 239 / denominator)); 240 set<1>(point, get<0, 1>(segment) 241 + boost::numeric_cast<CoordinateType>(numerator * dy_calc 242 / denominator)); 243 } 244 245 public : 246 template <typename Point, typename Segment1, typename Segment2> calculateboost::geometry::strategy::intersection::cartesian_segments::segment_intersection_info247 void calculate(Point& point, Segment1 const& a, Segment2 const& b) const 248 { 249 bool use_a = true; 250 251 // Prefer one segment if one is on or near an endpoint 252 bool const a_near_end = robust_ra.near_end(); 253 bool const b_near_end = robust_rb.near_end(); 254 if (a_near_end && ! b_near_end) 255 { 256 use_a = true; 257 } 258 else if (b_near_end && ! a_near_end) 259 { 260 use_a = false; 261 } 262 else 263 { 264 // Prefer shorter segment 265 promoted_type const len_a = comparable_length_a(); 266 promoted_type const len_b = comparable_length_b(); 267 if (len_b < len_a) 268 { 269 use_a = false; 270 } 271 // else use_a is true but was already assigned like that 272 } 273 274 if (use_a) 275 { 276 assign_a(point, a, b); 277 } 278 else 279 { 280 assign_b(point, a, b); 281 } 282 } 283 284 CoordinateType dx_a, dy_a; 285 CoordinateType dx_b, dy_b; 286 SegmentRatio robust_ra; 287 SegmentRatio robust_rb; 288 }; 289 290 template <typename D, typename W, typename ResultType> cramers_ruleboost::geometry::strategy::intersection::cartesian_segments291 static inline void cramers_rule(D const& dx_a, D const& dy_a, 292 D const& dx_b, D const& dy_b, W const& wx, W const& wy, 293 // out: 294 ResultType& nominator, ResultType& denominator) 295 { 296 // Cramers rule 297 nominator = geometry::detail::determinant<ResultType>(dx_b, dy_b, wx, wy); 298 denominator = geometry::detail::determinant<ResultType>(dx_a, dy_a, dx_b, dy_b); 299 // Ratio r = nominator/denominator 300 // Collinear if denominator == 0, intersecting if 0 <= r <= 1 301 // IntersectionPoint = (x1 + r * dx_a, y1 + r * dy_a) 302 } 303 304 // Version for non-rescaled policies 305 template 306 < 307 typename UniqueSubRange1, 308 typename UniqueSubRange2, 309 typename Policy 310 > 311 static inline typename Policy::return_type applyboost::geometry::strategy::intersection::cartesian_segments312 apply(UniqueSubRange1 const& range_p, 313 UniqueSubRange2 const& range_q, 314 Policy const& policy) 315 { 316 // Pass the same ranges both as normal ranges and as modelled ranges 317 return apply(range_p, range_q, policy, range_p, range_q); 318 } 319 320 // Version for non rescaled versions. 321 // The "modelled" parameter might be rescaled (will be removed later) 322 template 323 < 324 typename UniqueSubRange1, 325 typename UniqueSubRange2, 326 typename Policy, 327 typename ModelledUniqueSubRange1, 328 typename ModelledUniqueSubRange2 329 > 330 static inline typename Policy::return_type applyboost::geometry::strategy::intersection::cartesian_segments331 apply(UniqueSubRange1 const& range_p, 332 UniqueSubRange2 const& range_q, 333 Policy const& policy, 334 ModelledUniqueSubRange1 const& modelled_range_p, 335 ModelledUniqueSubRange2 const& modelled_range_q) 336 { 337 typedef typename UniqueSubRange1::point_type point1_type; 338 typedef typename UniqueSubRange2::point_type point2_type; 339 340 BOOST_CONCEPT_ASSERT( (concepts::ConstPoint<point1_type>) ); 341 BOOST_CONCEPT_ASSERT( (concepts::ConstPoint<point2_type>) ); 342 343 point1_type const& p1 = range_p.at(0); 344 point1_type const& p2 = range_p.at(1); 345 point2_type const& q1 = range_q.at(0); 346 point2_type const& q2 = range_q.at(1); 347 348 // Declare segments, currently necessary for the policies 349 // (segment_crosses, segment_colinear, degenerate, one_degenerate, etc) 350 model::referring_segment<point1_type const> const p(p1, p2); 351 model::referring_segment<point2_type const> const q(q1, q2); 352 353 typedef typename select_most_precise 354 < 355 typename geometry::coordinate_type<typename ModelledUniqueSubRange1::point_type>::type, 356 typename geometry::coordinate_type<typename ModelledUniqueSubRange1::point_type>::type 357 >::type modelled_coordinate_type; 358 359 typedef segment_ratio<modelled_coordinate_type> ratio_type; 360 segment_intersection_info 361 < 362 typename select_calculation_type<point1_type, point2_type, CalculationType>::type, 363 ratio_type 364 > sinfo; 365 366 sinfo.dx_a = get<0>(p2) - get<0>(p1); // distance in x-dir 367 sinfo.dx_b = get<0>(q2) - get<0>(q1); 368 sinfo.dy_a = get<1>(p2) - get<1>(p1); // distance in y-dir 369 sinfo.dy_b = get<1>(q2) - get<1>(q1); 370 371 return unified<ratio_type>(sinfo, p, q, policy, modelled_range_p, modelled_range_q); 372 } 373 374 //! Returns true if two segments do not overlap. 375 //! If not, then no further calculations need to be done. 376 template 377 < 378 std::size_t Dimension, 379 typename CoordinateType, 380 typename PointP, 381 typename PointQ 382 > disjoint_by_rangeboost::geometry::strategy::intersection::cartesian_segments383 static inline bool disjoint_by_range(PointP const& p1, PointP const& p2, 384 PointQ const& q1, PointQ const& q2) 385 { 386 CoordinateType minp = get<Dimension>(p1); 387 CoordinateType maxp = get<Dimension>(p2); 388 CoordinateType minq = get<Dimension>(q1); 389 CoordinateType maxq = get<Dimension>(q2); 390 if (minp > maxp) 391 { 392 std::swap(minp, maxp); 393 } 394 if (minq > maxq) 395 { 396 std::swap(minq, maxq); 397 } 398 399 // In this case, max(p) < min(q) 400 // P Q 401 // <-------> <-------> 402 // (and the space in between is not extremely small) 403 return math::smaller(maxp, minq) || math::smaller(maxq, minp); 404 } 405 406 // Implementation for either rescaled or non rescaled versions. 407 template 408 < 409 typename RatioType, 410 typename SegmentInfo, 411 typename Segment1, 412 typename Segment2, 413 typename Policy, 414 typename UniqueSubRange1, 415 typename UniqueSubRange2 416 > 417 static inline typename Policy::return_type unifiedboost::geometry::strategy::intersection::cartesian_segments418 unified(SegmentInfo& sinfo, 419 Segment1 const& p, Segment2 const& q, Policy const&, 420 UniqueSubRange1 const& range_p, 421 UniqueSubRange2 const& range_q) 422 { 423 typedef typename UniqueSubRange1::point_type point1_type; 424 typedef typename UniqueSubRange2::point_type point2_type; 425 typedef typename select_most_precise 426 < 427 typename geometry::coordinate_type<point1_type>::type, 428 typename geometry::coordinate_type<point2_type>::type 429 >::type coordinate_type; 430 431 point1_type const& p1 = range_p.at(0); 432 point1_type const& p2 = range_p.at(1); 433 point2_type const& q1 = range_q.at(0); 434 point2_type const& q2 = range_q.at(1); 435 436 using geometry::detail::equals::equals_point_point; 437 bool const p_is_point = equals_point_point(p1, p2, point_in_point_strategy_type()); 438 bool const q_is_point = equals_point_point(q1, q2, point_in_point_strategy_type()); 439 440 if (p_is_point && q_is_point) 441 { 442 return equals_point_point(p1, q2, point_in_point_strategy_type()) 443 ? Policy::degenerate(p, true) 444 : Policy::disjoint() 445 ; 446 } 447 448 if (disjoint_by_range<0, coordinate_type>(p1, p2, q1, q2) 449 || disjoint_by_range<1, coordinate_type>(p1, p2, q1, q2)) 450 { 451 return Policy::disjoint(); 452 } 453 454 side_info sides; 455 sides.set<0>(side_strategy_type::apply(q1, q2, p1), 456 side_strategy_type::apply(q1, q2, p2)); 457 458 if (sides.same<0>()) 459 { 460 // Both points are at same side of other segment, we can leave 461 return Policy::disjoint(); 462 } 463 464 sides.set<1>(side_strategy_type::apply(p1, p2, q1), 465 side_strategy_type::apply(p1, p2, q2)); 466 467 if (sides.same<1>()) 468 { 469 // Both points are at same side of other segment, we can leave 470 return Policy::disjoint(); 471 } 472 473 bool collinear = sides.collinear(); 474 475 // Calculate the differences again 476 // (for rescaled version, this is different from dx_p etc) 477 coordinate_type const dx_p = get<0>(p2) - get<0>(p1); 478 coordinate_type const dx_q = get<0>(q2) - get<0>(q1); 479 coordinate_type const dy_p = get<1>(p2) - get<1>(p1); 480 coordinate_type const dy_q = get<1>(q2) - get<1>(q1); 481 482 // r: ratio 0-1 where intersection divides A/B 483 // (only calculated for non-collinear segments) 484 if (! collinear) 485 { 486 coordinate_type denominator_a, nominator_a; 487 coordinate_type denominator_b, nominator_b; 488 489 cramers_rule(dx_p, dy_p, dx_q, dy_q, 490 get<0>(p1) - get<0>(q1), 491 get<1>(p1) - get<1>(q1), 492 nominator_a, denominator_a); 493 494 cramers_rule(dx_q, dy_q, dx_p, dy_p, 495 get<0>(q1) - get<0>(p1), 496 get<1>(q1) - get<1>(p1), 497 nominator_b, denominator_b); 498 499 math::detail::equals_factor_policy<coordinate_type> 500 policy(dx_p, dy_p, dx_q, dy_q); 501 502 coordinate_type const zero = 0; 503 if (math::detail::equals_by_policy(denominator_a, zero, policy) 504 || math::detail::equals_by_policy(denominator_b, zero, policy)) 505 { 506 // If this is the case, no rescaling is done for FP precision. 507 // We set it to collinear, but it indicates a robustness issue. 508 sides.set<0>(0, 0); 509 sides.set<1>(0, 0); 510 collinear = true; 511 } 512 else 513 { 514 sinfo.robust_ra.assign(nominator_a, denominator_a); 515 sinfo.robust_rb.assign(nominator_b, denominator_b); 516 } 517 } 518 519 if (collinear) 520 { 521 std::pair<bool, bool> const collinear_use_first 522 = is_x_more_significant(geometry::math::abs(dx_p), 523 geometry::math::abs(dy_p), 524 geometry::math::abs(dx_q), 525 geometry::math::abs(dy_q), 526 p_is_point, q_is_point); 527 528 if (collinear_use_first.second) 529 { 530 // Degenerate cases: segments of single point, lying on other segment, are not disjoint 531 // This situation is collinear too 532 533 if (collinear_use_first.first) 534 { 535 return relate_collinear<0, Policy, RatioType>(p, q, 536 p1, p2, q1, q2, 537 p_is_point, q_is_point); 538 } 539 else 540 { 541 // Y direction contains larger segments (maybe dx is zero) 542 return relate_collinear<1, Policy, RatioType>(p, q, 543 p1, p2, q1, q2, 544 p_is_point, q_is_point); 545 } 546 } 547 } 548 549 return Policy::segments_crosses(sides, sinfo, p, q); 550 } 551 552 private: 553 // first is true if x is more significant 554 // second is true if the more significant difference is not 0 555 template <typename CoordinateType> 556 static inline std::pair<bool, bool> is_x_more_significantboost::geometry::strategy::intersection::cartesian_segments557 is_x_more_significant(CoordinateType const& abs_dx_a, 558 CoordinateType const& abs_dy_a, 559 CoordinateType const& abs_dx_b, 560 CoordinateType const& abs_dy_b, 561 bool const a_is_point, 562 bool const b_is_point) 563 { 564 //BOOST_GEOMETRY_ASSERT_MSG(!(a_is_point && b_is_point), "both segments shouldn't be degenerated"); 565 566 // for degenerated segments the second is always true because this function 567 // shouldn't be called if both segments were degenerated 568 569 if (a_is_point) 570 { 571 return std::make_pair(abs_dx_b >= abs_dy_b, true); 572 } 573 else if (b_is_point) 574 { 575 return std::make_pair(abs_dx_a >= abs_dy_a, true); 576 } 577 else 578 { 579 CoordinateType const min_dx = (std::min)(abs_dx_a, abs_dx_b); 580 CoordinateType const min_dy = (std::min)(abs_dy_a, abs_dy_b); 581 return min_dx == min_dy ? 582 std::make_pair(true, min_dx > CoordinateType(0)) : 583 std::make_pair(min_dx > min_dy, true); 584 } 585 } 586 587 template 588 < 589 std::size_t Dimension, 590 typename Policy, 591 typename RatioType, 592 typename Segment1, 593 typename Segment2, 594 typename RobustPoint1, 595 typename RobustPoint2 596 > 597 static inline typename Policy::return_type relate_collinearboost::geometry::strategy::intersection::cartesian_segments598 relate_collinear(Segment1 const& a, 599 Segment2 const& b, 600 RobustPoint1 const& robust_a1, RobustPoint1 const& robust_a2, 601 RobustPoint2 const& robust_b1, RobustPoint2 const& robust_b2, 602 bool a_is_point, bool b_is_point) 603 { 604 if (a_is_point) 605 { 606 return relate_one_degenerate<Policy, RatioType>(a, 607 get<Dimension>(robust_a1), 608 get<Dimension>(robust_b1), get<Dimension>(robust_b2), 609 true); 610 } 611 if (b_is_point) 612 { 613 return relate_one_degenerate<Policy, RatioType>(b, 614 get<Dimension>(robust_b1), 615 get<Dimension>(robust_a1), get<Dimension>(robust_a2), 616 false); 617 } 618 return relate_collinear<Policy, RatioType>(a, b, 619 get<Dimension>(robust_a1), 620 get<Dimension>(robust_a2), 621 get<Dimension>(robust_b1), 622 get<Dimension>(robust_b2)); 623 } 624 625 /// Relate segments known collinear 626 template 627 < 628 typename Policy, 629 typename RatioType, 630 typename Segment1, 631 typename Segment2, 632 typename Type1, 633 typename Type2 634 > 635 static inline typename Policy::return_type relate_collinearboost::geometry::strategy::intersection::cartesian_segments636 relate_collinear(Segment1 const& a, Segment2 const& b, 637 Type1 oa_1, Type1 oa_2, 638 Type2 ob_1, Type2 ob_2) 639 { 640 // Calculate the ratios where a starts in b, b starts in a 641 // a1--------->a2 (2..7) 642 // b1----->b2 (5..8) 643 // length_a: 7-2=5 644 // length_b: 8-5=3 645 // b1 is located w.r.t. a at ratio: (5-2)/5=3/5 (on a) 646 // b2 is located w.r.t. a at ratio: (8-2)/5=6/5 (right of a) 647 // a1 is located w.r.t. b at ratio: (2-5)/3=-3/3 (left of b) 648 // a2 is located w.r.t. b at ratio: (7-5)/3=2/3 (on b) 649 // A arrives (a2 on b), B departs (b1 on a) 650 651 // If both are reversed: 652 // a2<---------a1 (7..2) 653 // b2<-----b1 (8..5) 654 // length_a: 2-7=-5 655 // length_b: 5-8=-3 656 // b1 is located w.r.t. a at ratio: (8-7)/-5=-1/5 (before a starts) 657 // b2 is located w.r.t. a at ratio: (5-7)/-5=2/5 (on a) 658 // a1 is located w.r.t. b at ratio: (7-8)/-3=1/3 (on b) 659 // a2 is located w.r.t. b at ratio: (2-8)/-3=6/3 (after b ends) 660 661 // If both one is reversed: 662 // a1--------->a2 (2..7) 663 // b2<-----b1 (8..5) 664 // length_a: 7-2=+5 665 // length_b: 5-8=-3 666 // b1 is located w.r.t. a at ratio: (8-2)/5=6/5 (after a ends) 667 // b2 is located w.r.t. a at ratio: (5-2)/5=3/5 (on a) 668 // a1 is located w.r.t. b at ratio: (2-8)/-3=6/3 (after b ends) 669 // a2 is located w.r.t. b at ratio: (7-8)/-3=1/3 (on b) 670 Type1 const length_a = oa_2 - oa_1; // no abs, see above 671 Type2 const length_b = ob_2 - ob_1; 672 673 RatioType ra_from(oa_1 - ob_1, length_b); 674 RatioType ra_to(oa_2 - ob_1, length_b); 675 RatioType rb_from(ob_1 - oa_1, length_a); 676 RatioType rb_to(ob_2 - oa_1, length_a); 677 678 // use absolute measure to detect endpoints intersection 679 // NOTE: it'd be possible to calculate bx_wrt_a using ax_wrt_b values 680 int const a1_wrt_b = position_value(oa_1, ob_1, ob_2); 681 int const a2_wrt_b = position_value(oa_2, ob_1, ob_2); 682 int const b1_wrt_a = position_value(ob_1, oa_1, oa_2); 683 int const b2_wrt_a = position_value(ob_2, oa_1, oa_2); 684 685 // fix the ratios if necessary 686 // CONSIDER: fixing ratios also in other cases, if they're inconsistent 687 // e.g. if ratio == 1 or 0 (so IP at the endpoint) 688 // but position value indicates that the IP is in the middle of the segment 689 // because one of the segments is very long 690 // In such case the ratios could be moved into the middle direction 691 // by some small value (e.g. EPS+1ULP) 692 if (a1_wrt_b == 1) 693 { 694 ra_from.assign(0, 1); 695 rb_from.assign(0, 1); 696 } 697 else if (a1_wrt_b == 3) 698 { 699 ra_from.assign(1, 1); 700 rb_to.assign(0, 1); 701 } 702 703 if (a2_wrt_b == 1) 704 { 705 ra_to.assign(0, 1); 706 rb_from.assign(1, 1); 707 } 708 else if (a2_wrt_b == 3) 709 { 710 ra_to.assign(1, 1); 711 rb_to.assign(1, 1); 712 } 713 714 if ((a1_wrt_b < 1 && a2_wrt_b < 1) || (a1_wrt_b > 3 && a2_wrt_b > 3)) 715 //if ((ra_from.left() && ra_to.left()) || (ra_from.right() && ra_to.right())) 716 { 717 return Policy::disjoint(); 718 } 719 720 bool const opposite = math::sign(length_a) != math::sign(length_b); 721 722 return Policy::segments_collinear(a, b, opposite, 723 a1_wrt_b, a2_wrt_b, b1_wrt_a, b2_wrt_a, 724 ra_from, ra_to, rb_from, rb_to); 725 } 726 727 /// Relate segments where one is degenerate 728 template 729 < 730 typename Policy, 731 typename RatioType, 732 typename DegenerateSegment, 733 typename Type1, 734 typename Type2 735 > 736 static inline typename Policy::return_type relate_one_degenerateboost::geometry::strategy::intersection::cartesian_segments737 relate_one_degenerate(DegenerateSegment const& degenerate_segment, 738 Type1 d, Type2 s1, Type2 s2, 739 bool a_degenerate) 740 { 741 // Calculate the ratios where ds starts in s 742 // a1--------->a2 (2..6) 743 // b1/b2 (4..4) 744 // Ratio: (4-2)/(6-2) 745 RatioType const ratio(d - s1, s2 - s1); 746 747 if (!ratio.on_segment()) 748 { 749 return Policy::disjoint(); 750 } 751 752 return Policy::one_degenerate(degenerate_segment, ratio, a_degenerate); 753 } 754 755 template <typename ProjCoord1, typename ProjCoord2> position_valueboost::geometry::strategy::intersection::cartesian_segments756 static inline int position_value(ProjCoord1 const& ca1, 757 ProjCoord2 const& cb1, 758 ProjCoord2 const& cb2) 759 { 760 // S1x 0 1 2 3 4 761 // S2 |----------> 762 return math::equals(ca1, cb1) ? 1 763 : math::equals(ca1, cb2) ? 3 764 : cb1 < cb2 ? 765 ( ca1 < cb1 ? 0 766 : ca1 > cb2 ? 4 767 : 2 ) 768 : ( ca1 > cb1 ? 0 769 : ca1 < cb2 ? 4 770 : 2 ); 771 } 772 }; 773 774 775 #ifndef DOXYGEN_NO_STRATEGY_SPECIALIZATIONS 776 namespace services 777 { 778 779 template <typename CalculationType> 780 struct default_strategy<cartesian_tag, CalculationType> 781 { 782 typedef cartesian_segments<CalculationType> type; 783 }; 784 785 } // namespace services 786 #endif // DOXYGEN_NO_STRATEGY_SPECIALIZATIONS 787 788 789 }} // namespace strategy::intersection 790 791 namespace strategy 792 { 793 794 namespace within { namespace services 795 { 796 797 template <typename Geometry1, typename Geometry2, typename AnyTag1, typename AnyTag2> 798 struct default_strategy<Geometry1, Geometry2, AnyTag1, AnyTag2, linear_tag, linear_tag, cartesian_tag, cartesian_tag> 799 { 800 typedef strategy::intersection::cartesian_segments<> type; 801 }; 802 803 template <typename Geometry1, typename Geometry2, typename AnyTag1, typename AnyTag2> 804 struct default_strategy<Geometry1, Geometry2, AnyTag1, AnyTag2, linear_tag, polygonal_tag, cartesian_tag, cartesian_tag> 805 { 806 typedef strategy::intersection::cartesian_segments<> type; 807 }; 808 809 template <typename Geometry1, typename Geometry2, typename AnyTag1, typename AnyTag2> 810 struct default_strategy<Geometry1, Geometry2, AnyTag1, AnyTag2, polygonal_tag, linear_tag, cartesian_tag, cartesian_tag> 811 { 812 typedef strategy::intersection::cartesian_segments<> type; 813 }; 814 815 template <typename Geometry1, typename Geometry2, typename AnyTag1, typename AnyTag2> 816 struct default_strategy<Geometry1, Geometry2, AnyTag1, AnyTag2, polygonal_tag, polygonal_tag, cartesian_tag, cartesian_tag> 817 { 818 typedef strategy::intersection::cartesian_segments<> type; 819 }; 820 821 }} // within::services 822 823 namespace covered_by { namespace services 824 { 825 826 template <typename Geometry1, typename Geometry2, typename AnyTag1, typename AnyTag2> 827 struct default_strategy<Geometry1, Geometry2, AnyTag1, AnyTag2, linear_tag, linear_tag, cartesian_tag, cartesian_tag> 828 { 829 typedef strategy::intersection::cartesian_segments<> type; 830 }; 831 832 template <typename Geometry1, typename Geometry2, typename AnyTag1, typename AnyTag2> 833 struct default_strategy<Geometry1, Geometry2, AnyTag1, AnyTag2, linear_tag, polygonal_tag, cartesian_tag, cartesian_tag> 834 { 835 typedef strategy::intersection::cartesian_segments<> type; 836 }; 837 838 template <typename Geometry1, typename Geometry2, typename AnyTag1, typename AnyTag2> 839 struct default_strategy<Geometry1, Geometry2, AnyTag1, AnyTag2, polygonal_tag, linear_tag, cartesian_tag, cartesian_tag> 840 { 841 typedef strategy::intersection::cartesian_segments<> type; 842 }; 843 844 template <typename Geometry1, typename Geometry2, typename AnyTag1, typename AnyTag2> 845 struct default_strategy<Geometry1, Geometry2, AnyTag1, AnyTag2, polygonal_tag, polygonal_tag, cartesian_tag, cartesian_tag> 846 { 847 typedef strategy::intersection::cartesian_segments<> type; 848 }; 849 850 }} // within::services 851 852 } // strategy 853 854 }} // namespace boost::geometry 855 856 857 #endif // BOOST_GEOMETRY_STRATEGIES_CARTESIAN_INTERSECTION_HPP 858