distance_cross_track.hpp 9.4 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336
  1. // Boost.Geometry (aka GGL, Generic Geometry Library)
  2. // Copyright (c) 2007-2012 Barend Gehrels, Amsterdam, the Netherlands.
  3. // Use, modification and distribution is subject to the Boost Software License,
  4. // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
  5. // http://www.boost.org/LICENSE_1_0.txt)
  6. #ifndef BOOST_GEOMETRY_STRATEGIES_SPHERICAL_DISTANCE_CROSS_TRACK_HPP
  7. #define BOOST_GEOMETRY_STRATEGIES_SPHERICAL_DISTANCE_CROSS_TRACK_HPP
  8. #include <boost/config.hpp>
  9. #include <boost/concept_check.hpp>
  10. #include <boost/mpl/if.hpp>
  11. #include <boost/type_traits.hpp>
  12. #include <boost/geometry/core/cs.hpp>
  13. #include <boost/geometry/core/access.hpp>
  14. #include <boost/geometry/core/radian_access.hpp>
  15. #include <boost/geometry/strategies/distance.hpp>
  16. #include <boost/geometry/strategies/concepts/distance_concept.hpp>
  17. #include <boost/geometry/strategies/spherical/distance_haversine.hpp>
  18. #include <boost/geometry/util/promote_floating_point.hpp>
  19. #include <boost/geometry/util/math.hpp>
  20. #ifdef BOOST_GEOMETRY_DEBUG_CROSS_TRACK
  21. # include <boost/geometry/io/dsv/write.hpp>
  22. #endif
  23. namespace boost { namespace geometry
  24. {
  25. namespace strategy { namespace distance
  26. {
  27. /*!
  28. \brief Strategy functor for distance point to segment calculation
  29. \ingroup strategies
  30. \details Class which calculates the distance of a point to a segment, for points on a sphere or globe
  31. \see http://williams.best.vwh.net/avform.htm
  32. \tparam CalculationType \tparam_calculation
  33. \tparam Strategy underlying point-point distance strategy, defaults to haversine
  34. \qbk{
  35. [heading See also]
  36. [link geometry.reference.algorithms.distance.distance_3_with_strategy distance (with strategy)]
  37. }
  38. */
  39. template
  40. <
  41. typename CalculationType = void,
  42. typename Strategy = haversine<double, CalculationType>
  43. >
  44. class cross_track
  45. {
  46. public :
  47. template <typename Point, typename PointOfSegment>
  48. struct return_type
  49. : promote_floating_point
  50. <
  51. typename select_calculation_type
  52. <
  53. Point,
  54. PointOfSegment,
  55. CalculationType
  56. >::type
  57. >
  58. {};
  59. inline cross_track()
  60. {}
  61. explicit inline cross_track(typename Strategy::radius_type const& r)
  62. : m_strategy(r)
  63. {}
  64. inline cross_track(Strategy const& s)
  65. : m_strategy(s)
  66. {}
  67. // It might be useful in the future
  68. // to overload constructor with strategy info.
  69. // crosstrack(...) {}
  70. template <typename Point, typename PointOfSegment>
  71. inline typename return_type<Point, PointOfSegment>::type
  72. apply(Point const& p, PointOfSegment const& sp1, PointOfSegment const& sp2) const
  73. {
  74. #if !defined(BOOST_MSVC)
  75. BOOST_CONCEPT_ASSERT
  76. (
  77. (concept::PointDistanceStrategy<Strategy, Point, PointOfSegment>)
  78. );
  79. #endif
  80. typedef typename return_type<Point, PointOfSegment>::type return_type;
  81. // http://williams.best.vwh.net/avform.htm#XTE
  82. return_type d1 = m_strategy.apply(sp1, p);
  83. return_type d3 = m_strategy.apply(sp1, sp2);
  84. if (geometry::math::equals(d3, 0.0))
  85. {
  86. // "Degenerate" segment, return either d1 or d2
  87. return d1;
  88. }
  89. return_type d2 = m_strategy.apply(sp2, p);
  90. return_type crs_AD = course(sp1, p);
  91. return_type crs_AB = course(sp1, sp2);
  92. return_type crs_BA = crs_AB - geometry::math::pi<return_type>();
  93. return_type crs_BD = course(sp2, p);
  94. return_type d_crs1 = crs_AD - crs_AB;
  95. return_type d_crs2 = crs_BD - crs_BA;
  96. // d1, d2, d3 are in principle not needed, only the sign matters
  97. return_type projection1 = cos( d_crs1 ) * d1 / d3;
  98. return_type projection2 = cos( d_crs2 ) * d2 / d3;
  99. #ifdef BOOST_GEOMETRY_DEBUG_CROSS_TRACK
  100. std::cout << "Course " << dsv(sp1) << " to " << dsv(p) << " " << crs_AD * geometry::math::r2d << std::endl;
  101. std::cout << "Course " << dsv(sp1) << " to " << dsv(sp2) << " " << crs_AB * geometry::math::r2d << std::endl;
  102. std::cout << "Course " << dsv(sp2) << " to " << dsv(p) << " " << crs_BD * geometry::math::r2d << std::endl;
  103. std::cout << "Projection AD-AB " << projection1 << " : " << d_crs1 * geometry::math::r2d << std::endl;
  104. std::cout << "Projection BD-BA " << projection2 << " : " << d_crs2 * geometry::math::r2d << std::endl;
  105. #endif
  106. if(projection1 > 0.0 && projection2 > 0.0)
  107. {
  108. return_type XTD = radius() * geometry::math::abs( asin( sin( d1 / radius() ) * sin( d_crs1 ) ));
  109. #ifdef BOOST_GEOMETRY_DEBUG_CROSS_TRACK
  110. std::cout << "Projection ON the segment" << std::endl;
  111. std::cout << "XTD: " << XTD << " d1: " << d1 << " d2: " << d2 << std::endl;
  112. #endif
  113. // Return shortest distance, projected point on segment sp1-sp2
  114. return return_type(XTD);
  115. }
  116. else
  117. {
  118. #ifdef BOOST_GEOMETRY_DEBUG_CROSS_TRACK
  119. std::cout << "Projection OUTSIDE the segment" << std::endl;
  120. #endif
  121. // Return shortest distance, project either on point sp1 or sp2
  122. return return_type( (std::min)( d1 , d2 ) );
  123. }
  124. }
  125. inline typename Strategy::radius_type radius() const
  126. { return m_strategy.radius(); }
  127. private :
  128. Strategy m_strategy;
  129. /// Calculate course (bearing) between two points. Might be moved to a "course formula" ...
  130. template <typename Point1, typename Point2>
  131. inline typename return_type<Point1, Point2>::type
  132. course(Point1 const& p1, Point2 const& p2) const
  133. {
  134. typedef typename return_type<Point1, Point2>::type return_type;
  135. // http://williams.best.vwh.net/avform.htm#Crs
  136. return_type dlon = get_as_radian<0>(p2) - get_as_radian<0>(p1);
  137. return_type cos_p2lat = cos(get_as_radian<1>(p2));
  138. // "An alternative formula, not requiring the pre-computation of d"
  139. return atan2(sin(dlon) * cos_p2lat,
  140. cos(get_as_radian<1>(p1)) * sin(get_as_radian<1>(p2))
  141. - sin(get_as_radian<1>(p1)) * cos_p2lat * cos(dlon));
  142. }
  143. };
  144. #ifndef DOXYGEN_NO_STRATEGY_SPECIALIZATIONS
  145. namespace services
  146. {
  147. template <typename CalculationType, typename Strategy>
  148. struct tag<cross_track<CalculationType, Strategy> >
  149. {
  150. typedef strategy_tag_distance_point_segment type;
  151. };
  152. template <typename CalculationType, typename Strategy, typename P, typename PS>
  153. struct return_type<cross_track<CalculationType, Strategy>, P, PS>
  154. : cross_track<CalculationType, Strategy>::template return_type<P, PS>
  155. {};
  156. template <typename CalculationType, typename Strategy>
  157. struct comparable_type<cross_track<CalculationType, Strategy> >
  158. {
  159. // There is no shortcut, so the strategy itself is its comparable type
  160. typedef cross_track<CalculationType, Strategy> type;
  161. };
  162. template
  163. <
  164. typename CalculationType,
  165. typename Strategy
  166. >
  167. struct get_comparable<cross_track<CalculationType, Strategy> >
  168. {
  169. typedef typename comparable_type
  170. <
  171. cross_track<CalculationType, Strategy>
  172. >::type comparable_type;
  173. public :
  174. static inline comparable_type apply(cross_track<CalculationType, Strategy> const& strategy)
  175. {
  176. return comparable_type(strategy.radius());
  177. }
  178. };
  179. template
  180. <
  181. typename CalculationType,
  182. typename Strategy,
  183. typename P, typename PS
  184. >
  185. struct result_from_distance<cross_track<CalculationType, Strategy>, P, PS>
  186. {
  187. private :
  188. typedef typename cross_track<CalculationType, Strategy>::template return_type<P, PS> return_type;
  189. public :
  190. template <typename T>
  191. static inline return_type apply(cross_track<CalculationType, Strategy> const& , T const& distance)
  192. {
  193. return distance;
  194. }
  195. };
  196. template
  197. <
  198. typename CalculationType,
  199. typename Strategy
  200. >
  201. struct strategy_point_point<cross_track<CalculationType, Strategy> >
  202. {
  203. typedef Strategy type;
  204. };
  205. /*
  206. TODO: spherical polar coordinate system requires "get_as_radian_equatorial<>"
  207. template <typename Point, typename PointOfSegment, typename Strategy>
  208. struct default_strategy
  209. <
  210. segment_tag, Point, PointOfSegment,
  211. spherical_polar_tag, spherical_polar_tag,
  212. Strategy
  213. >
  214. {
  215. typedef cross_track
  216. <
  217. void,
  218. typename boost::mpl::if_
  219. <
  220. boost::is_void<Strategy>,
  221. typename default_strategy
  222. <
  223. point_tag, Point, PointOfSegment,
  224. spherical_polar_tag, spherical_polar_tag
  225. >::type,
  226. Strategy
  227. >::type
  228. > type;
  229. };
  230. */
  231. template <typename Point, typename PointOfSegment, typename Strategy>
  232. struct default_strategy
  233. <
  234. segment_tag, Point, PointOfSegment,
  235. spherical_equatorial_tag, spherical_equatorial_tag,
  236. Strategy
  237. >
  238. {
  239. typedef cross_track
  240. <
  241. void,
  242. typename boost::mpl::if_
  243. <
  244. boost::is_void<Strategy>,
  245. typename default_strategy
  246. <
  247. point_tag, Point, PointOfSegment,
  248. spherical_equatorial_tag, spherical_equatorial_tag
  249. >::type,
  250. Strategy
  251. >::type
  252. > type;
  253. };
  254. } // namespace services
  255. #endif // DOXYGEN_NO_STRATEGY_SPECIALIZATIONS
  256. }} // namespace strategy::distance
  257. #ifndef DOXYGEN_NO_STRATEGY_SPECIALIZATIONS
  258. #endif
  259. }} // namespace boost::geometry
  260. #endif // BOOST_GEOMETRY_STRATEGIES_SPHERICAL_DISTANCE_CROSS_TRACK_HPP