side_by_triangle.hpp 9.4 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297
  1. // Boost.Geometry (aka GGL, Generic Geometry Library)
  2. // Copyright (c) 2007-2015 Barend Gehrels, Amsterdam, the Netherlands.
  3. // Copyright (c) 2008-2015 Bruno Lalande, Paris, France.
  4. // Copyright (c) 2009-2015 Mateusz Loskot, London, UK.
  5. // This file was modified by Oracle on 2015, 2017, 2018.
  6. // Modifications copyright (c) 2015-2018, Oracle and/or its affiliates.
  7. // Contributed and/or modified by Menelaos Karavelas, on behalf of Oracle
  8. // Contributed and/or modified by Adam Wulkiewicz, on behalf of Oracle
  9. // Parts of Boost.Geometry are redesigned from Geodan's Geographic Library
  10. // (geolib/GGL), copyright (c) 1995-2010 Geodan, Amsterdam, the Netherlands.
  11. // Use, modification and distribution is subject to the Boost Software License,
  12. // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
  13. // http://www.boost.org/LICENSE_1_0.txt)
  14. #ifndef BOOST_GEOMETRY_STRATEGIES_CARTESIAN_SIDE_BY_TRIANGLE_HPP
  15. #define BOOST_GEOMETRY_STRATEGIES_CARTESIAN_SIDE_BY_TRIANGLE_HPP
  16. #include <boost/mpl/if.hpp>
  17. #include <boost/type_traits/is_integral.hpp>
  18. #include <boost/type_traits/is_void.hpp>
  19. #include <boost/geometry/arithmetic/determinant.hpp>
  20. #include <boost/geometry/core/access.hpp>
  21. #include <boost/geometry/util/select_coordinate_type.hpp>
  22. #include <boost/geometry/strategies/cartesian/disjoint_segment_box.hpp>
  23. #include <boost/geometry/strategies/cartesian/envelope.hpp>
  24. #include <boost/geometry/strategies/cartesian/point_in_point.hpp>
  25. #include <boost/geometry/strategies/compare.hpp>
  26. #include <boost/geometry/strategies/side.hpp>
  27. #include <boost/geometry/algorithms/detail/equals/point_point.hpp>
  28. namespace boost { namespace geometry
  29. {
  30. namespace strategy { namespace side
  31. {
  32. /*!
  33. \brief Check at which side of a segment a point lies:
  34. left of segment (> 0), right of segment (< 0), on segment (0)
  35. \ingroup strategies
  36. \tparam CalculationType \tparam_calculation
  37. */
  38. template <typename CalculationType = void>
  39. class side_by_triangle
  40. {
  41. template <typename Policy>
  42. struct eps_policy
  43. {
  44. eps_policy() {}
  45. template <typename Type>
  46. eps_policy(Type const& a, Type const& b, Type const& c, Type const& d)
  47. : policy(a, b, c, d)
  48. {}
  49. Policy policy;
  50. };
  51. struct eps_empty
  52. {
  53. eps_empty() {}
  54. template <typename Type>
  55. eps_empty(Type const&, Type const&, Type const&, Type const&) {}
  56. };
  57. public :
  58. typedef strategy::envelope::cartesian<CalculationType> envelope_strategy_type;
  59. static inline envelope_strategy_type get_envelope_strategy()
  60. {
  61. return envelope_strategy_type();
  62. }
  63. typedef strategy::disjoint::segment_box disjoint_strategy_type;
  64. static inline disjoint_strategy_type get_disjoint_strategy()
  65. {
  66. return disjoint_strategy_type();
  67. }
  68. typedef strategy::within::cartesian_point_point equals_point_point_strategy_type;
  69. static inline equals_point_point_strategy_type get_equals_point_point_strategy()
  70. {
  71. return equals_point_point_strategy_type();
  72. }
  73. // Template member function, because it is not always trivial
  74. // or convenient to explicitly mention the typenames in the
  75. // strategy-struct itself.
  76. // Types can be all three different. Therefore it is
  77. // not implemented (anymore) as "segment"
  78. template
  79. <
  80. typename CoordinateType,
  81. typename PromotedType,
  82. typename P1,
  83. typename P2,
  84. typename P,
  85. typename EpsPolicy
  86. >
  87. static inline
  88. PromotedType side_value(P1 const& p1, P2 const& p2, P const& p, EpsPolicy & eps_policy)
  89. {
  90. CoordinateType const x = get<0>(p);
  91. CoordinateType const y = get<1>(p);
  92. CoordinateType const sx1 = get<0>(p1);
  93. CoordinateType const sy1 = get<1>(p1);
  94. CoordinateType const sx2 = get<0>(p2);
  95. CoordinateType const sy2 = get<1>(p2);
  96. PromotedType const dx = sx2 - sx1;
  97. PromotedType const dy = sy2 - sy1;
  98. PromotedType const dpx = x - sx1;
  99. PromotedType const dpy = y - sy1;
  100. eps_policy = EpsPolicy(dx, dy, dpx, dpy);
  101. return geometry::detail::determinant<PromotedType>
  102. (
  103. dx, dy,
  104. dpx, dpy
  105. );
  106. }
  107. template
  108. <
  109. typename CoordinateType,
  110. typename PromotedType,
  111. typename P1,
  112. typename P2,
  113. typename P
  114. >
  115. static inline
  116. PromotedType side_value(P1 const& p1, P2 const& p2, P const& p)
  117. {
  118. eps_empty dummy;
  119. return side_value<CoordinateType, PromotedType>(p1, p2, p, dummy);
  120. }
  121. template
  122. <
  123. typename CoordinateType,
  124. typename PromotedType,
  125. bool AreAllIntegralCoordinates
  126. >
  127. struct compute_side_value
  128. {
  129. template <typename P1, typename P2, typename P, typename EpsPolicy>
  130. static inline PromotedType apply(P1 const& p1, P2 const& p2, P const& p, EpsPolicy & epsp)
  131. {
  132. return side_value<CoordinateType, PromotedType>(p1, p2, p, epsp);
  133. }
  134. };
  135. template <typename CoordinateType, typename PromotedType>
  136. struct compute_side_value<CoordinateType, PromotedType, false>
  137. {
  138. template <typename P1, typename P2, typename P, typename EpsPolicy>
  139. static inline PromotedType apply(P1 const& p1, P2 const& p2, P const& p, EpsPolicy & epsp)
  140. {
  141. // For robustness purposes, first check if any two points are
  142. // the same; in this case simply return that the points are
  143. // collinear
  144. if (equals_point_point(p1, p2)
  145. || equals_point_point(p1, p)
  146. || equals_point_point(p2, p))
  147. {
  148. return PromotedType(0);
  149. }
  150. // The side_by_triangle strategy computes the signed area of
  151. // the point triplet (p1, p2, p); as such it is (in theory)
  152. // invariant under cyclic permutations of its three arguments.
  153. //
  154. // In the context of numerical errors that arise in
  155. // floating-point computations, and in order to make the strategy
  156. // consistent with respect to cyclic permutations of its three
  157. // arguments, we cyclically permute them so that the first
  158. // argument is always the lexicographically smallest point.
  159. typedef compare::cartesian<compare::less> less;
  160. if (less::apply(p, p1))
  161. {
  162. if (less::apply(p, p2))
  163. {
  164. // p is the lexicographically smallest
  165. return side_value<CoordinateType, PromotedType>(p, p1, p2, epsp);
  166. }
  167. else
  168. {
  169. // p2 is the lexicographically smallest
  170. return side_value<CoordinateType, PromotedType>(p2, p, p1, epsp);
  171. }
  172. }
  173. if (less::apply(p1, p2))
  174. {
  175. // p1 is the lexicographically smallest
  176. return side_value<CoordinateType, PromotedType>(p1, p2, p, epsp);
  177. }
  178. else
  179. {
  180. // p2 is the lexicographically smallest
  181. return side_value<CoordinateType, PromotedType>(p2, p, p1, epsp);
  182. }
  183. }
  184. };
  185. template <typename P1, typename P2, typename P>
  186. static inline int apply(P1 const& p1, P2 const& p2, P const& p)
  187. {
  188. typedef typename coordinate_type<P1>::type coordinate_type1;
  189. typedef typename coordinate_type<P2>::type coordinate_type2;
  190. typedef typename coordinate_type<P>::type coordinate_type3;
  191. typedef typename boost::mpl::if_c
  192. <
  193. boost::is_void<CalculationType>::type::value,
  194. typename select_most_precise
  195. <
  196. typename select_most_precise
  197. <
  198. coordinate_type1, coordinate_type2
  199. >::type,
  200. coordinate_type3
  201. >::type,
  202. CalculationType
  203. >::type coordinate_type;
  204. // Promote float->double, small int->int
  205. typedef typename select_most_precise
  206. <
  207. coordinate_type,
  208. double
  209. >::type promoted_type;
  210. bool const are_all_integral_coordinates =
  211. boost::is_integral<coordinate_type1>::value
  212. && boost::is_integral<coordinate_type2>::value
  213. && boost::is_integral<coordinate_type3>::value;
  214. eps_policy< math::detail::equals_factor_policy<promoted_type> > epsp;
  215. promoted_type s = compute_side_value
  216. <
  217. coordinate_type, promoted_type, are_all_integral_coordinates
  218. >::apply(p1, p2, p, epsp);
  219. promoted_type const zero = promoted_type();
  220. return math::detail::equals_by_policy(s, zero, epsp.policy) ? 0
  221. : s > zero ? 1
  222. : -1;
  223. }
  224. private:
  225. template <typename P1, typename P2>
  226. static inline bool equals_point_point(P1 const& p1, P2 const& p2)
  227. {
  228. typedef equals_point_point_strategy_type strategy_t;
  229. return geometry::detail::equals::equals_point_point(p1, p2, strategy_t());
  230. }
  231. };
  232. #ifndef DOXYGEN_NO_STRATEGY_SPECIALIZATIONS
  233. namespace services
  234. {
  235. template <typename CalculationType>
  236. struct default_strategy<cartesian_tag, CalculationType>
  237. {
  238. typedef side_by_triangle<CalculationType> type;
  239. };
  240. }
  241. #endif
  242. }} // namespace strategy::side
  243. }} // namespace boost::geometry
  244. #endif // BOOST_GEOMETRY_STRATEGIES_CARTESIAN_SIDE_BY_TRIANGLE_HPP