distance_projected_point_ax.hpp 9.6 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315
  1. // Boost.Geometry (aka GGL, Generic Geometry Library)
  2. // Copyright (c) 2008-2014 Bruno Lalande, Paris, France.
  3. // Copyright (c) 2008-2014 Barend Gehrels, Amsterdam, the Netherlands.
  4. // Copyright (c) 2009-2014 Mateusz Loskot, London, UK.
  5. // This file was modified by Oracle on 2014.
  6. // Modifications copyright (c) 2014, 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_DISTANCE_PROJECTED_POINT_AX_HPP
  15. #define BOOST_GEOMETRY_STRATEGIES_CARTESIAN_DISTANCE_PROJECTED_POINT_AX_HPP
  16. #include <algorithm>
  17. #include <boost/concept_check.hpp>
  18. #include <boost/core/ignore_unused.hpp>
  19. #include <boost/geometry/core/access.hpp>
  20. #include <boost/geometry/core/point_type.hpp>
  21. #include <boost/geometry/algorithms/convert.hpp>
  22. #include <boost/geometry/arithmetic/arithmetic.hpp>
  23. #include <boost/geometry/arithmetic/dot_product.hpp>
  24. #include <boost/geometry/strategies/tags.hpp>
  25. #include <boost/geometry/strategies/distance.hpp>
  26. #include <boost/geometry/strategies/default_distance_result.hpp>
  27. #include <boost/geometry/strategies/cartesian/distance_pythagoras.hpp>
  28. #include <boost/geometry/strategies/cartesian/distance_projected_point.hpp>
  29. #include <boost/geometry/util/select_coordinate_type.hpp>
  30. // Helper geometry (projected point on line)
  31. #include <boost/geometry/geometries/point.hpp>
  32. namespace boost { namespace geometry
  33. {
  34. namespace strategy { namespace distance
  35. {
  36. #ifndef DOXYGEN_NO_DETAIL
  37. namespace detail
  38. {
  39. template <typename T>
  40. struct projected_point_ax_result
  41. {
  42. typedef T value_type;
  43. projected_point_ax_result(T const& c = T(0))
  44. : atd(c), xtd(c)
  45. {}
  46. projected_point_ax_result(T const& a, T const& x)
  47. : atd(a), xtd(x)
  48. {}
  49. friend inline bool operator<(projected_point_ax_result const& left,
  50. projected_point_ax_result const& right)
  51. {
  52. return left.xtd < right.xtd || left.atd < right.atd;
  53. }
  54. T atd, xtd;
  55. };
  56. // This less-comparator may be used as a parameter of detail::douglas_peucker.
  57. // In this simplify strategy distances are compared in 2 places
  58. // 1. to choose the furthest candidate (md < dist)
  59. // 2. to check if the candidate is further than max_distance (max_distance < md)
  60. template <typename Distance>
  61. class projected_point_ax_less
  62. {
  63. public:
  64. projected_point_ax_less(Distance const& max_distance)
  65. : m_max_distance(max_distance)
  66. {}
  67. inline bool operator()(Distance const& left, Distance const& right) const
  68. {
  69. //return left.xtd < right.xtd && right.atd < m_max_distance.atd;
  70. typedef typename Distance::value_type value_type;
  71. value_type const lx = left.xtd > m_max_distance.xtd ? left.xtd - m_max_distance.xtd : 0;
  72. value_type const rx = right.xtd > m_max_distance.xtd ? right.xtd - m_max_distance.xtd : 0;
  73. value_type const la = left.atd > m_max_distance.atd ? left.atd - m_max_distance.atd : 0;
  74. value_type const ra = right.atd > m_max_distance.atd ? right.atd - m_max_distance.atd : 0;
  75. value_type const l = (std::max)(lx, la);
  76. value_type const r = (std::max)(rx, ra);
  77. return l < r;
  78. }
  79. private:
  80. Distance const& m_max_distance;
  81. };
  82. // This strategy returns 2-component Point/Segment distance.
  83. // The ATD (along track distance) is parallel to the Segment
  84. // and is a distance between Point projected into a line defined by a Segment and the nearest Segment's endpoint.
  85. // If the projected Point intersects the Segment the ATD is equal to 0.
  86. // The XTD (cross track distance) is perpendicular to the Segment
  87. // and is a distance between input Point and its projection.
  88. // If the Segment has length equal to 0, ATD and XTD has value equal
  89. // to the distance between the input Point and one of the Segment's endpoints.
  90. //
  91. // p3 p4
  92. // ^ 7
  93. // | /
  94. // p1<-----e========e----->p2
  95. //
  96. // p1: atd=D, xtd=0
  97. // p2: atd=D, xtd=0
  98. // p3: atd=0, xtd=D
  99. // p4: atd=D/2, xtd=D
  100. template
  101. <
  102. typename CalculationType = void,
  103. typename Strategy = pythagoras<CalculationType>
  104. >
  105. class projected_point_ax
  106. {
  107. public :
  108. template <typename Point, typename PointOfSegment>
  109. struct calculation_type
  110. : public projected_point<CalculationType, Strategy>
  111. ::template calculation_type<Point, PointOfSegment>
  112. {};
  113. template <typename Point, typename PointOfSegment>
  114. struct result_type
  115. {
  116. typedef projected_point_ax_result
  117. <
  118. typename calculation_type<Point, PointOfSegment>::type
  119. > type;
  120. };
  121. public :
  122. template <typename Point, typename PointOfSegment>
  123. inline typename result_type<Point, PointOfSegment>::type
  124. apply(Point const& p, PointOfSegment const& p1, PointOfSegment const& p2) const
  125. {
  126. assert_dimension_equal<Point, PointOfSegment>();
  127. typedef typename calculation_type<Point, PointOfSegment>::type calculation_type;
  128. // A projected point of points in Integer coordinates must be able to be
  129. // represented in FP.
  130. typedef model::point
  131. <
  132. calculation_type,
  133. dimension<PointOfSegment>::value,
  134. typename coordinate_system<PointOfSegment>::type
  135. > fp_point_type;
  136. // For convenience
  137. typedef fp_point_type fp_vector_type;
  138. /*
  139. Algorithm [p: (px,py), p1: (x1,y1), p2: (x2,y2)]
  140. VECTOR v(x2 - x1, y2 - y1)
  141. VECTOR w(px - x1, py - y1)
  142. c1 = w . v
  143. c2 = v . v
  144. b = c1 / c2
  145. RETURN POINT(x1 + b * vx, y1 + b * vy)
  146. */
  147. // v is multiplied below with a (possibly) FP-value, so should be in FP
  148. // For consistency we define w also in FP
  149. fp_vector_type v, w, projected;
  150. geometry::convert(p2, v);
  151. geometry::convert(p, w);
  152. geometry::convert(p1, projected);
  153. subtract_point(v, projected);
  154. subtract_point(w, projected);
  155. Strategy strategy;
  156. boost::ignore_unused(strategy);
  157. typename result_type<Point, PointOfSegment>::type result;
  158. calculation_type const zero = calculation_type();
  159. calculation_type const c2 = dot_product(v, v);
  160. if ( math::equals(c2, zero) )
  161. {
  162. result.xtd = strategy.apply(p, projected);
  163. // assume that the 0-length segment is perpendicular to the Pt->ProjPt vector
  164. result.atd = 0;
  165. return result;
  166. }
  167. calculation_type const c1 = dot_product(w, v);
  168. calculation_type const b = c1 / c2;
  169. multiply_value(v, b);
  170. add_point(projected, v);
  171. result.xtd = strategy.apply(p, projected);
  172. if (c1 <= zero)
  173. {
  174. result.atd = strategy.apply(p1, projected);
  175. }
  176. else if (c2 <= c1)
  177. {
  178. result.atd = strategy.apply(p2, projected);
  179. }
  180. else
  181. {
  182. result.atd = 0;
  183. }
  184. return result;
  185. }
  186. };
  187. } // namespace detail
  188. #endif // DOXYGEN_NO_DETAIL
  189. #ifndef DOXYGEN_NO_STRATEGY_SPECIALIZATIONS
  190. namespace services
  191. {
  192. template <typename CalculationType, typename Strategy>
  193. struct tag<detail::projected_point_ax<CalculationType, Strategy> >
  194. {
  195. typedef strategy_tag_distance_point_segment type;
  196. };
  197. template <typename CalculationType, typename Strategy, typename P, typename PS>
  198. struct return_type<detail::projected_point_ax<CalculationType, Strategy>, P, PS>
  199. {
  200. typedef typename detail::projected_point_ax<CalculationType, Strategy>
  201. ::template result_type<P, PS>::type type;
  202. };
  203. template <typename CalculationType, typename Strategy>
  204. struct comparable_type<detail::projected_point_ax<CalculationType, Strategy> >
  205. {
  206. // Define a projected_point strategy with its underlying point-point-strategy
  207. // being comparable
  208. typedef detail::projected_point_ax
  209. <
  210. CalculationType,
  211. typename comparable_type<Strategy>::type
  212. > type;
  213. };
  214. template <typename CalculationType, typename Strategy>
  215. struct get_comparable<detail::projected_point_ax<CalculationType, Strategy> >
  216. {
  217. typedef typename comparable_type
  218. <
  219. detail::projected_point_ax<CalculationType, Strategy>
  220. >::type comparable_type;
  221. public :
  222. static inline comparable_type apply(detail::projected_point_ax<CalculationType, Strategy> const& )
  223. {
  224. return comparable_type();
  225. }
  226. };
  227. template <typename CalculationType, typename Strategy, typename P, typename PS>
  228. struct result_from_distance<detail::projected_point_ax<CalculationType, Strategy>, P, PS>
  229. {
  230. private :
  231. typedef typename return_type<detail::projected_point_ax<CalculationType, Strategy>, P, PS>::type return_type;
  232. public :
  233. template <typename T>
  234. static inline return_type apply(detail::projected_point_ax<CalculationType, Strategy> const& , T const& value)
  235. {
  236. Strategy s;
  237. return_type ret;
  238. ret.atd = result_from_distance<Strategy, P, PS>::apply(s, value.atd);
  239. ret.xtd = result_from_distance<Strategy, P, PS>::apply(s, value.xtd);
  240. return ret;
  241. }
  242. };
  243. } // namespace services
  244. #endif // DOXYGEN_NO_STRATEGY_SPECIALIZATIONS
  245. }} // namespace strategy::distance
  246. }} // namespace boost::geometry
  247. #endif // BOOST_GEOMETRY_STRATEGIES_CARTESIAN_DISTANCE_PROJECTED_POINT_AX_HPP