]>
Commit | Line | Data |
---|---|---|
7c673cae FG |
1 | // Boost.Geometry (aka GGL, Generic Geometry Library) |
2 | ||
3 | // Copyright (c) 2007-2014 Barend Gehrels, Amsterdam, the Netherlands. | |
4 | // Copyright (c) 2013-2014 Adam Wulkiewicz, Lodz, Poland. | |
5 | ||
b32b8144 FG |
6 | // This file was modified by Oracle on 2014, 2016, 2017. |
7 | // Modifications copyright (c) 2014-2017, Oracle and/or its affiliates. | |
7c673cae FG |
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 | ||
b32b8144 FG |
36 | #include <boost/geometry/strategies/cartesian/area_surveyor.hpp> |
37 | #include <boost/geometry/strategies/cartesian/distance_pythagoras.hpp> | |
38 | #include <boost/geometry/strategies/cartesian/envelope_segment.hpp> | |
39 | #include <boost/geometry/strategies/cartesian/point_in_poly_winding.hpp> | |
7c673cae | 40 | #include <boost/geometry/strategies/cartesian/side_by_triangle.hpp> |
b32b8144 | 41 | #include <boost/geometry/strategies/covered_by.hpp> |
7c673cae FG |
42 | #include <boost/geometry/strategies/intersection.hpp> |
43 | #include <boost/geometry/strategies/intersection_result.hpp> | |
b32b8144 FG |
44 | #include <boost/geometry/strategies/side.hpp> |
45 | #include <boost/geometry/strategies/side_info.hpp> | |
46 | #include <boost/geometry/strategies/within.hpp> | |
7c673cae FG |
47 | |
48 | #include <boost/geometry/policies/robustness/robust_point_type.hpp> | |
49 | #include <boost/geometry/policies/robustness/segment_ratio_type.hpp> | |
50 | ||
51 | ||
52 | #if defined(BOOST_GEOMETRY_DEBUG_ROBUSTNESS) | |
53 | # include <boost/geometry/io/wkt/write.hpp> | |
54 | #endif | |
55 | ||
56 | ||
57 | namespace boost { namespace geometry | |
58 | { | |
59 | ||
60 | ||
61 | namespace strategy { namespace intersection | |
62 | { | |
63 | ||
64 | ||
65 | /*! | |
66 | \see http://mathworld.wolfram.com/Line-LineIntersection.html | |
67 | */ | |
b32b8144 FG |
68 | template |
69 | < | |
70 | typename CalculationType = void | |
71 | > | |
72 | struct cartesian_segments | |
7c673cae | 73 | { |
b32b8144 FG |
74 | typedef side::side_by_triangle<CalculationType> side_strategy_type; |
75 | ||
76 | static inline side_strategy_type get_side_strategy() | |
77 | { | |
78 | return side_strategy_type(); | |
79 | } | |
80 | ||
81 | template <typename Geometry1, typename Geometry2> | |
82 | struct point_in_geometry_strategy | |
83 | { | |
84 | typedef strategy::within::cartesian_winding | |
85 | < | |
86 | typename point_type<Geometry1>::type, | |
87 | typename point_type<Geometry2>::type, | |
88 | CalculationType | |
89 | > type; | |
90 | }; | |
91 | ||
92 | template <typename Geometry1, typename Geometry2> | |
93 | static inline typename point_in_geometry_strategy<Geometry1, Geometry2>::type | |
94 | get_point_in_geometry_strategy() | |
95 | { | |
96 | typedef typename point_in_geometry_strategy | |
97 | < | |
98 | Geometry1, Geometry2 | |
99 | >::type strategy_type; | |
100 | return strategy_type(); | |
101 | } | |
102 | ||
103 | template <typename Geometry> | |
104 | struct area_strategy | |
105 | { | |
106 | typedef area::surveyor | |
107 | < | |
108 | typename point_type<Geometry>::type, | |
109 | CalculationType | |
110 | > type; | |
111 | }; | |
112 | ||
113 | template <typename Geometry> | |
114 | static inline typename area_strategy<Geometry>::type get_area_strategy() | |
115 | { | |
116 | typedef typename area_strategy<Geometry>::type strategy_type; | |
117 | return strategy_type(); | |
118 | } | |
119 | ||
120 | template <typename Geometry> | |
121 | struct distance_strategy | |
122 | { | |
123 | typedef distance::pythagoras | |
124 | < | |
125 | CalculationType | |
126 | > type; | |
127 | }; | |
128 | ||
129 | template <typename Geometry> | |
130 | static inline typename distance_strategy<Geometry>::type get_distance_strategy() | |
131 | { | |
132 | typedef typename distance_strategy<Geometry>::type strategy_type; | |
133 | return strategy_type(); | |
134 | } | |
135 | ||
136 | typedef envelope::cartesian_segment<CalculationType> | |
137 | envelope_strategy_type; | |
138 | ||
139 | static inline envelope_strategy_type get_envelope_strategy() | |
140 | { | |
141 | return envelope_strategy_type(); | |
142 | } | |
7c673cae FG |
143 | |
144 | template <typename CoordinateType, typename SegmentRatio> | |
145 | struct segment_intersection_info | |
146 | { | |
147 | typedef typename select_most_precise | |
148 | < | |
149 | CoordinateType, double | |
150 | >::type promoted_type; | |
151 | ||
152 | promoted_type comparable_length_a() const | |
153 | { | |
154 | return dx_a * dx_a + dy_a * dy_a; | |
155 | } | |
156 | ||
157 | promoted_type comparable_length_b() const | |
158 | { | |
159 | return dx_b * dx_b + dy_b * dy_b; | |
160 | } | |
161 | ||
162 | template <typename Point, typename Segment1, typename Segment2> | |
163 | void assign_a(Point& point, Segment1 const& a, Segment2 const& ) const | |
164 | { | |
165 | assign(point, a, dx_a, dy_a, robust_ra); | |
166 | } | |
167 | template <typename Point, typename Segment1, typename Segment2> | |
168 | void assign_b(Point& point, Segment1 const& , Segment2 const& b) const | |
169 | { | |
170 | assign(point, b, dx_b, dy_b, robust_rb); | |
171 | } | |
172 | ||
173 | template <typename Point, typename Segment> | |
174 | void assign(Point& point, Segment const& segment, CoordinateType const& dx, CoordinateType const& dy, SegmentRatio const& ratio) const | |
175 | { | |
176 | // Calculate the intersection point based on segment_ratio | |
177 | // Up to now, division was postponed. Here we divide using numerator/ | |
178 | // denominator. In case of integer this results in an integer | |
179 | // division. | |
180 | BOOST_GEOMETRY_ASSERT(ratio.denominator() != 0); | |
181 | ||
182 | typedef typename promote_integral<CoordinateType>::type promoted_type; | |
183 | ||
184 | promoted_type const numerator | |
185 | = boost::numeric_cast<promoted_type>(ratio.numerator()); | |
186 | promoted_type const denominator | |
187 | = boost::numeric_cast<promoted_type>(ratio.denominator()); | |
188 | promoted_type const dx_promoted = boost::numeric_cast<promoted_type>(dx); | |
189 | promoted_type const dy_promoted = boost::numeric_cast<promoted_type>(dy); | |
190 | ||
191 | set<0>(point, get<0, 0>(segment) + boost::numeric_cast | |
192 | < | |
193 | CoordinateType | |
194 | >(numerator * dx_promoted / denominator)); | |
195 | set<1>(point, get<0, 1>(segment) + boost::numeric_cast | |
196 | < | |
197 | CoordinateType | |
198 | >(numerator * dy_promoted / denominator)); | |
199 | } | |
200 | ||
201 | CoordinateType dx_a, dy_a; | |
202 | CoordinateType dx_b, dy_b; | |
203 | SegmentRatio robust_ra; | |
204 | SegmentRatio robust_rb; | |
205 | }; | |
206 | ||
207 | template <typename D, typename W, typename ResultType> | |
208 | static inline void cramers_rule(D const& dx_a, D const& dy_a, | |
209 | D const& dx_b, D const& dy_b, W const& wx, W const& wy, | |
210 | // out: | |
211 | ResultType& d, ResultType& da) | |
212 | { | |
213 | // Cramers rule | |
214 | d = geometry::detail::determinant<ResultType>(dx_a, dy_a, dx_b, dy_b); | |
215 | da = geometry::detail::determinant<ResultType>(dx_b, dy_b, wx, wy); | |
216 | // Ratio is da/d , collinear if d == 0, intersecting if 0 <= r <= 1 | |
217 | // IntersectionPoint = (x1 + r * dx_a, y1 + r * dy_a) | |
218 | } | |
219 | ||
220 | ||
221 | // Relate segments a and b | |
b32b8144 FG |
222 | template |
223 | < | |
224 | typename Segment1, | |
225 | typename Segment2, | |
226 | typename Policy, | |
227 | typename RobustPolicy | |
228 | > | |
229 | static inline typename Policy::return_type | |
230 | apply(Segment1 const& a, Segment2 const& b, | |
231 | Policy const& policy, RobustPolicy const& robust_policy) | |
7c673cae FG |
232 | { |
233 | // type them all as in Segment1 - TODO reconsider this, most precise? | |
234 | typedef typename geometry::point_type<Segment1>::type point_type; | |
235 | ||
236 | typedef typename geometry::robust_point_type | |
237 | < | |
238 | point_type, RobustPolicy | |
239 | >::type robust_point_type; | |
240 | ||
241 | point_type a0, a1, b0, b1; | |
242 | robust_point_type a0_rob, a1_rob, b0_rob, b1_rob; | |
243 | ||
244 | detail::assign_point_from_index<0>(a, a0); | |
245 | detail::assign_point_from_index<1>(a, a1); | |
246 | detail::assign_point_from_index<0>(b, b0); | |
247 | detail::assign_point_from_index<1>(b, b1); | |
248 | ||
249 | geometry::recalculate(a0_rob, a0, robust_policy); | |
250 | geometry::recalculate(a1_rob, a1, robust_policy); | |
251 | geometry::recalculate(b0_rob, b0, robust_policy); | |
252 | geometry::recalculate(b1_rob, b1, robust_policy); | |
253 | ||
b32b8144 | 254 | return apply(a, b, policy, robust_policy, a0_rob, a1_rob, b0_rob, b1_rob); |
7c673cae FG |
255 | } |
256 | ||
257 | // The main entry-routine, calculating intersections of segments a / b | |
258 | // NOTE: Robust* types may be the same as Segments' point types | |
b32b8144 FG |
259 | template |
260 | < | |
261 | typename Segment1, | |
262 | typename Segment2, | |
263 | typename Policy, | |
264 | typename RobustPolicy, | |
265 | typename RobustPoint1, | |
266 | typename RobustPoint2 | |
267 | > | |
268 | static inline typename Policy::return_type | |
269 | apply(Segment1 const& a, Segment2 const& b, | |
270 | Policy const&, RobustPolicy const& /*robust_policy*/, | |
271 | RobustPoint1 const& robust_a1, RobustPoint1 const& robust_a2, | |
272 | RobustPoint2 const& robust_b1, RobustPoint2 const& robust_b2) | |
7c673cae FG |
273 | { |
274 | BOOST_CONCEPT_ASSERT( (concepts::ConstSegment<Segment1>) ); | |
275 | BOOST_CONCEPT_ASSERT( (concepts::ConstSegment<Segment2>) ); | |
276 | ||
277 | using geometry::detail::equals::equals_point_point; | |
278 | bool const a_is_point = equals_point_point(robust_a1, robust_a2); | |
279 | bool const b_is_point = equals_point_point(robust_b1, robust_b2); | |
280 | ||
281 | if(a_is_point && b_is_point) | |
282 | { | |
283 | return equals_point_point(robust_a1, robust_b2) | |
284 | ? Policy::degenerate(a, true) | |
285 | : Policy::disjoint() | |
286 | ; | |
287 | } | |
288 | ||
7c673cae | 289 | side_info sides; |
b32b8144 FG |
290 | sides.set<0>(side_strategy_type::apply(robust_b1, robust_b2, robust_a1), |
291 | side_strategy_type::apply(robust_b1, robust_b2, robust_a2)); | |
7c673cae FG |
292 | |
293 | if (sides.same<0>()) | |
294 | { | |
295 | // Both points are at same side of other segment, we can leave | |
296 | return Policy::disjoint(); | |
297 | } | |
298 | ||
b32b8144 FG |
299 | sides.set<1>(side_strategy_type::apply(robust_a1, robust_a2, robust_b1), |
300 | side_strategy_type::apply(robust_a1, robust_a2, robust_b2)); | |
7c673cae FG |
301 | |
302 | if (sides.same<1>()) | |
303 | { | |
304 | // Both points are at same side of other segment, we can leave | |
305 | return Policy::disjoint(); | |
306 | } | |
307 | ||
308 | bool collinear = sides.collinear(); | |
309 | ||
310 | typedef typename select_most_precise | |
311 | < | |
312 | typename geometry::coordinate_type<RobustPoint1>::type, | |
313 | typename geometry::coordinate_type<RobustPoint2>::type | |
314 | >::type robust_coordinate_type; | |
315 | ||
316 | typedef typename segment_ratio_type | |
b32b8144 FG |
317 | < |
318 | typename geometry::point_type<Segment1>::type, // TODO: most precise point? | |
319 | RobustPolicy | |
320 | >::type ratio_type; | |
7c673cae FG |
321 | |
322 | segment_intersection_info | |
b32b8144 FG |
323 | < |
324 | typename select_calculation_type<Segment1, Segment2, CalculationType>::type, | |
325 | ratio_type | |
326 | > sinfo; | |
7c673cae FG |
327 | |
328 | sinfo.dx_a = get<1, 0>(a) - get<0, 0>(a); // distance in x-dir | |
329 | sinfo.dx_b = get<1, 0>(b) - get<0, 0>(b); | |
330 | sinfo.dy_a = get<1, 1>(a) - get<0, 1>(a); // distance in y-dir | |
331 | sinfo.dy_b = get<1, 1>(b) - get<0, 1>(b); | |
332 | ||
333 | robust_coordinate_type const robust_dx_a = get<0>(robust_a2) - get<0>(robust_a1); | |
334 | robust_coordinate_type const robust_dx_b = get<0>(robust_b2) - get<0>(robust_b1); | |
335 | robust_coordinate_type const robust_dy_a = get<1>(robust_a2) - get<1>(robust_a1); | |
336 | robust_coordinate_type const robust_dy_b = get<1>(robust_b2) - get<1>(robust_b1); | |
337 | ||
338 | // r: ratio 0-1 where intersection divides A/B | |
339 | // (only calculated for non-collinear segments) | |
340 | if (! collinear) | |
341 | { | |
342 | robust_coordinate_type robust_da0, robust_da; | |
343 | robust_coordinate_type robust_db0, robust_db; | |
344 | ||
345 | cramers_rule(robust_dx_a, robust_dy_a, robust_dx_b, robust_dy_b, | |
346 | get<0>(robust_a1) - get<0>(robust_b1), | |
347 | get<1>(robust_a1) - get<1>(robust_b1), | |
348 | robust_da0, robust_da); | |
349 | ||
350 | cramers_rule(robust_dx_b, robust_dy_b, robust_dx_a, robust_dy_a, | |
351 | get<0>(robust_b1) - get<0>(robust_a1), | |
352 | get<1>(robust_b1) - get<1>(robust_a1), | |
353 | robust_db0, robust_db); | |
354 | ||
355 | math::detail::equals_factor_policy<robust_coordinate_type> | |
356 | policy(robust_dx_a, robust_dy_a, robust_dx_b, robust_dy_b); | |
357 | robust_coordinate_type const zero = 0; | |
358 | if (math::detail::equals_by_policy(robust_da0, zero, policy) | |
359 | || math::detail::equals_by_policy(robust_db0, zero, policy)) | |
360 | { | |
361 | // If this is the case, no rescaling is done for FP precision. | |
362 | // We set it to collinear, but it indicates a robustness issue. | |
363 | sides.set<0>(0,0); | |
364 | sides.set<1>(0,0); | |
365 | collinear = true; | |
366 | } | |
367 | else | |
368 | { | |
369 | sinfo.robust_ra.assign(robust_da, robust_da0); | |
370 | sinfo.robust_rb.assign(robust_db, robust_db0); | |
371 | } | |
372 | } | |
373 | ||
374 | if (collinear) | |
375 | { | |
376 | std::pair<bool, bool> const collinear_use_first | |
377 | = is_x_more_significant(geometry::math::abs(robust_dx_a), | |
378 | geometry::math::abs(robust_dy_a), | |
379 | geometry::math::abs(robust_dx_b), | |
380 | geometry::math::abs(robust_dy_b), | |
381 | a_is_point, b_is_point); | |
382 | ||
383 | if (collinear_use_first.second) | |
384 | { | |
385 | // Degenerate cases: segments of single point, lying on other segment, are not disjoint | |
386 | // This situation is collinear too | |
387 | ||
388 | if (collinear_use_first.first) | |
389 | { | |
b32b8144 | 390 | return relate_collinear<0, Policy, ratio_type>(a, b, |
7c673cae FG |
391 | robust_a1, robust_a2, robust_b1, robust_b2, |
392 | a_is_point, b_is_point); | |
393 | } | |
394 | else | |
395 | { | |
396 | // Y direction contains larger segments (maybe dx is zero) | |
b32b8144 | 397 | return relate_collinear<1, Policy, ratio_type>(a, b, |
7c673cae FG |
398 | robust_a1, robust_a2, robust_b1, robust_b2, |
399 | a_is_point, b_is_point); | |
400 | } | |
401 | } | |
402 | } | |
403 | ||
404 | return Policy::segments_crosses(sides, sinfo, a, b); | |
405 | } | |
406 | ||
407 | private: | |
408 | // first is true if x is more significant | |
409 | // second is true if the more significant difference is not 0 | |
410 | template <typename RobustCoordinateType> | |
411 | static inline std::pair<bool, bool> | |
412 | is_x_more_significant(RobustCoordinateType const& abs_robust_dx_a, | |
413 | RobustCoordinateType const& abs_robust_dy_a, | |
414 | RobustCoordinateType const& abs_robust_dx_b, | |
415 | RobustCoordinateType const& abs_robust_dy_b, | |
416 | bool const a_is_point, | |
417 | bool const b_is_point) | |
418 | { | |
419 | //BOOST_GEOMETRY_ASSERT_MSG(!(a_is_point && b_is_point), "both segments shouldn't be degenerated"); | |
420 | ||
421 | // for degenerated segments the second is always true because this function | |
422 | // shouldn't be called if both segments were degenerated | |
423 | ||
424 | if (a_is_point) | |
425 | { | |
426 | return std::make_pair(abs_robust_dx_b >= abs_robust_dy_b, true); | |
427 | } | |
428 | else if (b_is_point) | |
429 | { | |
430 | return std::make_pair(abs_robust_dx_a >= abs_robust_dy_a, true); | |
431 | } | |
432 | else | |
433 | { | |
434 | RobustCoordinateType const min_dx = (std::min)(abs_robust_dx_a, abs_robust_dx_b); | |
435 | RobustCoordinateType const min_dy = (std::min)(abs_robust_dy_a, abs_robust_dy_b); | |
436 | return min_dx == min_dy ? | |
437 | std::make_pair(true, min_dx > RobustCoordinateType(0)) : | |
438 | std::make_pair(min_dx > min_dy, true); | |
439 | } | |
440 | } | |
441 | ||
442 | template | |
443 | < | |
444 | std::size_t Dimension, | |
b32b8144 | 445 | typename Policy, |
7c673cae FG |
446 | typename RatioType, |
447 | typename Segment1, | |
448 | typename Segment2, | |
449 | typename RobustPoint1, | |
450 | typename RobustPoint2 | |
451 | > | |
b32b8144 FG |
452 | static inline typename Policy::return_type |
453 | relate_collinear(Segment1 const& a, | |
454 | Segment2 const& b, | |
455 | RobustPoint1 const& robust_a1, RobustPoint1 const& robust_a2, | |
456 | RobustPoint2 const& robust_b1, RobustPoint2 const& robust_b2, | |
457 | bool a_is_point, bool b_is_point) | |
7c673cae FG |
458 | { |
459 | if (a_is_point) | |
460 | { | |
b32b8144 | 461 | return relate_one_degenerate<Policy, RatioType>(a, |
7c673cae FG |
462 | get<Dimension>(robust_a1), |
463 | get<Dimension>(robust_b1), get<Dimension>(robust_b2), | |
464 | true); | |
465 | } | |
466 | if (b_is_point) | |
467 | { | |
b32b8144 | 468 | return relate_one_degenerate<Policy, RatioType>(b, |
7c673cae FG |
469 | get<Dimension>(robust_b1), |
470 | get<Dimension>(robust_a1), get<Dimension>(robust_a2), | |
471 | false); | |
472 | } | |
b32b8144 | 473 | return relate_collinear<Policy, RatioType>(a, b, |
7c673cae FG |
474 | get<Dimension>(robust_a1), |
475 | get<Dimension>(robust_a2), | |
476 | get<Dimension>(robust_b1), | |
477 | get<Dimension>(robust_b2)); | |
478 | } | |
479 | ||
480 | /// Relate segments known collinear | |
481 | template | |
482 | < | |
b32b8144 | 483 | typename Policy, |
7c673cae FG |
484 | typename RatioType, |
485 | typename Segment1, | |
486 | typename Segment2, | |
487 | typename RobustType1, | |
488 | typename RobustType2 | |
489 | > | |
b32b8144 FG |
490 | static inline typename Policy::return_type |
491 | relate_collinear(Segment1 const& a, Segment2 const& b, | |
492 | RobustType1 oa_1, RobustType1 oa_2, | |
493 | RobustType2 ob_1, RobustType2 ob_2) | |
7c673cae FG |
494 | { |
495 | // Calculate the ratios where a starts in b, b starts in a | |
496 | // a1--------->a2 (2..7) | |
497 | // b1----->b2 (5..8) | |
498 | // length_a: 7-2=5 | |
499 | // length_b: 8-5=3 | |
500 | // b1 is located w.r.t. a at ratio: (5-2)/5=3/5 (on a) | |
501 | // b2 is located w.r.t. a at ratio: (8-2)/5=6/5 (right of a) | |
502 | // a1 is located w.r.t. b at ratio: (2-5)/3=-3/3 (left of b) | |
503 | // a2 is located w.r.t. b at ratio: (7-5)/3=2/3 (on b) | |
504 | // A arrives (a2 on b), B departs (b1 on a) | |
505 | ||
506 | // If both are reversed: | |
507 | // a2<---------a1 (7..2) | |
508 | // b2<-----b1 (8..5) | |
509 | // length_a: 2-7=-5 | |
510 | // length_b: 5-8=-3 | |
511 | // b1 is located w.r.t. a at ratio: (8-7)/-5=-1/5 (before a starts) | |
512 | // b2 is located w.r.t. a at ratio: (5-7)/-5=2/5 (on a) | |
513 | // a1 is located w.r.t. b at ratio: (7-8)/-3=1/3 (on b) | |
514 | // a2 is located w.r.t. b at ratio: (2-8)/-3=6/3 (after b ends) | |
515 | ||
516 | // If both one is reversed: | |
517 | // a1--------->a2 (2..7) | |
518 | // b2<-----b1 (8..5) | |
519 | // length_a: 7-2=+5 | |
520 | // length_b: 5-8=-3 | |
521 | // b1 is located w.r.t. a at ratio: (8-2)/5=6/5 (after a ends) | |
522 | // b2 is located w.r.t. a at ratio: (5-2)/5=3/5 (on a) | |
523 | // a1 is located w.r.t. b at ratio: (2-8)/-3=6/3 (after b ends) | |
524 | // a2 is located w.r.t. b at ratio: (7-8)/-3=1/3 (on b) | |
525 | RobustType1 const length_a = oa_2 - oa_1; // no abs, see above | |
526 | RobustType2 const length_b = ob_2 - ob_1; | |
527 | ||
528 | RatioType ra_from(oa_1 - ob_1, length_b); | |
529 | RatioType ra_to(oa_2 - ob_1, length_b); | |
530 | RatioType rb_from(ob_1 - oa_1, length_a); | |
531 | RatioType rb_to(ob_2 - oa_1, length_a); | |
532 | ||
533 | // use absolute measure to detect endpoints intersection | |
534 | // NOTE: it'd be possible to calculate bx_wrt_a using ax_wrt_b values | |
535 | int const a1_wrt_b = position_value(oa_1, ob_1, ob_2); | |
536 | int const a2_wrt_b = position_value(oa_2, ob_1, ob_2); | |
537 | int const b1_wrt_a = position_value(ob_1, oa_1, oa_2); | |
538 | int const b2_wrt_a = position_value(ob_2, oa_1, oa_2); | |
539 | ||
540 | // fix the ratios if necessary | |
541 | // CONSIDER: fixing ratios also in other cases, if they're inconsistent | |
542 | // e.g. if ratio == 1 or 0 (so IP at the endpoint) | |
543 | // but position value indicates that the IP is in the middle of the segment | |
544 | // because one of the segments is very long | |
545 | // In such case the ratios could be moved into the middle direction | |
546 | // by some small value (e.g. EPS+1ULP) | |
547 | if (a1_wrt_b == 1) | |
548 | { | |
549 | ra_from.assign(0, 1); | |
550 | rb_from.assign(0, 1); | |
551 | } | |
552 | else if (a1_wrt_b == 3) | |
553 | { | |
554 | ra_from.assign(1, 1); | |
555 | rb_to.assign(0, 1); | |
556 | } | |
557 | ||
558 | if (a2_wrt_b == 1) | |
559 | { | |
560 | ra_to.assign(0, 1); | |
561 | rb_from.assign(1, 1); | |
562 | } | |
563 | else if (a2_wrt_b == 3) | |
564 | { | |
565 | ra_to.assign(1, 1); | |
566 | rb_to.assign(1, 1); | |
567 | } | |
568 | ||
569 | if ((a1_wrt_b < 1 && a2_wrt_b < 1) || (a1_wrt_b > 3 && a2_wrt_b > 3)) | |
570 | //if ((ra_from.left() && ra_to.left()) || (ra_from.right() && ra_to.right())) | |
571 | { | |
572 | return Policy::disjoint(); | |
573 | } | |
574 | ||
575 | bool const opposite = math::sign(length_a) != math::sign(length_b); | |
576 | ||
577 | return Policy::segments_collinear(a, b, opposite, | |
578 | a1_wrt_b, a2_wrt_b, b1_wrt_a, b2_wrt_a, | |
579 | ra_from, ra_to, rb_from, rb_to); | |
580 | } | |
581 | ||
582 | /// Relate segments where one is degenerate | |
583 | template | |
584 | < | |
b32b8144 | 585 | typename Policy, |
7c673cae FG |
586 | typename RatioType, |
587 | typename DegenerateSegment, | |
588 | typename RobustType1, | |
589 | typename RobustType2 | |
590 | > | |
b32b8144 FG |
591 | static inline typename Policy::return_type |
592 | relate_one_degenerate(DegenerateSegment const& degenerate_segment, | |
593 | RobustType1 d, RobustType2 s1, RobustType2 s2, | |
594 | bool a_degenerate) | |
7c673cae FG |
595 | { |
596 | // Calculate the ratios where ds starts in s | |
597 | // a1--------->a2 (2..6) | |
598 | // b1/b2 (4..4) | |
599 | // Ratio: (4-2)/(6-2) | |
600 | RatioType const ratio(d - s1, s2 - s1); | |
601 | ||
602 | if (!ratio.on_segment()) | |
603 | { | |
604 | return Policy::disjoint(); | |
605 | } | |
606 | ||
607 | return Policy::one_degenerate(degenerate_segment, ratio, a_degenerate); | |
608 | } | |
609 | ||
610 | template <typename ProjCoord1, typename ProjCoord2> | |
611 | static inline int position_value(ProjCoord1 const& ca1, | |
612 | ProjCoord2 const& cb1, | |
613 | ProjCoord2 const& cb2) | |
614 | { | |
615 | // S1x 0 1 2 3 4 | |
616 | // S2 |----------> | |
617 | return math::equals(ca1, cb1) ? 1 | |
618 | : math::equals(ca1, cb2) ? 3 | |
619 | : cb1 < cb2 ? | |
620 | ( ca1 < cb1 ? 0 | |
621 | : ca1 > cb2 ? 4 | |
622 | : 2 ) | |
623 | : ( ca1 > cb1 ? 0 | |
624 | : ca1 < cb2 ? 4 | |
625 | : 2 ); | |
626 | } | |
627 | }; | |
628 | ||
629 | ||
630 | #ifndef DOXYGEN_NO_STRATEGY_SPECIALIZATIONS | |
631 | namespace services | |
632 | { | |
633 | ||
b32b8144 FG |
634 | template <typename CalculationType> |
635 | struct default_strategy<cartesian_tag, CalculationType> | |
7c673cae | 636 | { |
b32b8144 | 637 | typedef cartesian_segments<CalculationType> type; |
7c673cae FG |
638 | }; |
639 | ||
640 | } // namespace services | |
641 | #endif // DOXYGEN_NO_STRATEGY_SPECIALIZATIONS | |
642 | ||
643 | ||
644 | }} // namespace strategy::intersection | |
645 | ||
b32b8144 FG |
646 | namespace strategy |
647 | { | |
648 | ||
649 | namespace within { namespace services | |
650 | { | |
651 | ||
652 | template <typename Geometry1, typename Geometry2, typename AnyTag1, typename AnyTag2> | |
653 | struct default_strategy<Geometry1, Geometry2, AnyTag1, AnyTag2, linear_tag, linear_tag, cartesian_tag, cartesian_tag> | |
654 | { | |
655 | typedef strategy::intersection::cartesian_segments<> type; | |
656 | }; | |
657 | ||
658 | template <typename Geometry1, typename Geometry2, typename AnyTag1, typename AnyTag2> | |
659 | struct default_strategy<Geometry1, Geometry2, AnyTag1, AnyTag2, linear_tag, polygonal_tag, cartesian_tag, cartesian_tag> | |
660 | { | |
661 | typedef strategy::intersection::cartesian_segments<> type; | |
662 | }; | |
663 | ||
664 | template <typename Geometry1, typename Geometry2, typename AnyTag1, typename AnyTag2> | |
665 | struct default_strategy<Geometry1, Geometry2, AnyTag1, AnyTag2, polygonal_tag, linear_tag, cartesian_tag, cartesian_tag> | |
666 | { | |
667 | typedef strategy::intersection::cartesian_segments<> type; | |
668 | }; | |
669 | ||
670 | template <typename Geometry1, typename Geometry2, typename AnyTag1, typename AnyTag2> | |
671 | struct default_strategy<Geometry1, Geometry2, AnyTag1, AnyTag2, polygonal_tag, polygonal_tag, cartesian_tag, cartesian_tag> | |
672 | { | |
673 | typedef strategy::intersection::cartesian_segments<> type; | |
674 | }; | |
675 | ||
676 | }} // within::services | |
677 | ||
678 | namespace covered_by { namespace services | |
679 | { | |
680 | ||
681 | template <typename Geometry1, typename Geometry2, typename AnyTag1, typename AnyTag2> | |
682 | struct default_strategy<Geometry1, Geometry2, AnyTag1, AnyTag2, linear_tag, linear_tag, cartesian_tag, cartesian_tag> | |
683 | { | |
684 | typedef strategy::intersection::cartesian_segments<> type; | |
685 | }; | |
686 | ||
687 | template <typename Geometry1, typename Geometry2, typename AnyTag1, typename AnyTag2> | |
688 | struct default_strategy<Geometry1, Geometry2, AnyTag1, AnyTag2, linear_tag, polygonal_tag, cartesian_tag, cartesian_tag> | |
689 | { | |
690 | typedef strategy::intersection::cartesian_segments<> type; | |
691 | }; | |
692 | ||
693 | template <typename Geometry1, typename Geometry2, typename AnyTag1, typename AnyTag2> | |
694 | struct default_strategy<Geometry1, Geometry2, AnyTag1, AnyTag2, polygonal_tag, linear_tag, cartesian_tag, cartesian_tag> | |
695 | { | |
696 | typedef strategy::intersection::cartesian_segments<> type; | |
697 | }; | |
698 | ||
699 | template <typename Geometry1, typename Geometry2, typename AnyTag1, typename AnyTag2> | |
700 | struct default_strategy<Geometry1, Geometry2, AnyTag1, AnyTag2, polygonal_tag, polygonal_tag, cartesian_tag, cartesian_tag> | |
701 | { | |
702 | typedef strategy::intersection::cartesian_segments<> type; | |
703 | }; | |
704 | ||
705 | }} // within::services | |
706 | ||
707 | } // strategy | |
708 | ||
7c673cae FG |
709 | }} // namespace boost::geometry |
710 | ||
711 | ||
712 | #endif // BOOST_GEOMETRY_STRATEGIES_CARTESIAN_INTERSECTION_HPP |