boost/geometry/strategies/cartesian/distance_projected_point.hpp

   1 // Boost.Geometry (aka GGL, Generic Geometry Library)
   2
   3 // Copyright (c) 2008-2012 Bruno Lalande, Paris, France.
   4 // Copyright (c) 2008-2012 Barend Gehrels, Amsterdam, the Netherlands.
   5 // Copyright (c) 2009-2012 Mateusz Loskot, London, UK.
   6
   7 // Parts of Boost.Geometry are redesigned from Geodan's Geographic Library
   8 // (geolib/GGL), copyright (c) 1995-2010 Geodan, Amsterdam, the Netherlands.
   9
  10 // Use, modification and distribution is subject to the Boost Software License,
  11 // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at
  12 // http://www.boost.org/LICENSE_1_0.txt)
  13
  14 #ifndef BOOST_GEOMETRY_STRATEGIES_CARTESIAN_DISTANCE_PROJECTED_POINT_HPP
  15 #define BOOST_GEOMETRY_STRATEGIES_CARTESIAN_DISTANCE_PROJECTED_POINT_HPP
  16
  17
  18 #include <boost/concept_check.hpp>
  19 #include <boost/mpl/if.hpp>
  20 #include <boost/type_traits.hpp>
  21
  22 #include <boost/geometry/core/access.hpp>
  23 #include <boost/geometry/core/point_type.hpp>
  24
  25 #include <boost/geometry/algorithms/convert.hpp>
  26 #include <boost/geometry/arithmetic/arithmetic.hpp>
  27 #include <boost/geometry/arithmetic/dot_product.hpp>
  28
  29 #include <boost/geometry/strategies/tags.hpp>
  30 #include <boost/geometry/strategies/distance.hpp>
  31 #include <boost/geometry/strategies/default_distance_result.hpp>
  32 #include <boost/geometry/strategies/cartesian/distance_pythagoras.hpp>
  33
  34 #include <boost/geometry/util/select_coordinate_type.hpp>
  35
  36 // Helper geometry (projected point on line)
  37 #include <boost/geometry/geometries/point.hpp>
  38
  39
  40 namespace boost { namespace geometry
  41 {
  42
  43
  44 namespace strategy { namespace distance
  45 {
  46
  47
  48 /*!
  49 \brief Strategy for distance point to segment
  50 \ingroup strategies
  51 \details Calculates distance using projected-point method, and (optionally) Pythagoras
  52 \author Adapted from: http://geometryalgorithms.com/Archive/algorithm_0102/algorithm_0102.htm
  53 \tparam Point \tparam_point
  54 \tparam PointOfSegment \tparam_segment_point
  55 \tparam CalculationType \tparam_calculation
  56 \tparam Strategy underlying point-point distance strategy
  57 \par Concepts for Strategy:
  58 - cartesian_distance operator(Point,Point)
  59 \note If the Strategy is a "comparable::pythagoras", this strategy
  60     automatically is a comparable projected_point strategy (so without sqrt)
  61
  62 \qbk{
  63 [heading See also]
  64 [link geometry.reference.algorithms.distance.distance_3_with_strategy distance (with strategy)]
  65 }
  66
  67 */
  68 template
  69 <
  70     typename Point,
  71     typename PointOfSegment = Point,
  72     typename CalculationType = void,
  73     typename Strategy = pythagoras<Point, PointOfSegment, CalculationType>
  74 >
  75 class projected_point
  76 {
  77 public :
  78     // The three typedefs below are necessary to calculate distances
  79     // from segments defined in integer coordinates.
  80
  81     // Integer coordinates can still result in FP distances.
  82     // There is a division, which must be represented in FP.
  83     // So promote.
  84     typedef typename promote_floating_point
  85         <
  86             typename strategy::distance::services::return_type
  87                 <
  88                     Strategy
  89                 >::type
  90         >::type calculation_type;
  91
  92 private :
  93
  94     // A projected point of points in Integer coordinates must be able to be
  95     // represented in FP.
  96     typedef model::point
  97         <
  98             calculation_type,
  99             dimension<PointOfSegment>::value,
 100             typename coordinate_system<PointOfSegment>::type
 101         > fp_point_type;
 102
 103     // For convenience
 104     typedef fp_point_type fp_vector_type;
 105
 106     // We have to use a strategy using FP coordinates (fp-type) which is
 107     // not always the same as Strategy (defined as point_strategy_type)
 108     // So we create a "similar" one
 109     typedef typename strategy::distance::services::similar_type
 110         <
 111             Strategy,
 112             Point,
 113             fp_point_type
 114         >::type fp_strategy_type;
 115
 116 public :
 117
 118     inline calculation_type apply(Point const& p,
 119                     PointOfSegment const& p1, PointOfSegment const& p2) const
 120     {
 121         assert_dimension_equal<Point, PointOfSegment>();
 122
 123         /*
 124             Algorithm [p1: (x1,y1), p2: (x2,y2), p: (px,py)]
 125             VECTOR v(x2 - x1, y2 - y1)
 126             VECTOR w(px - x1, py - y1)
 127             c1 = w . v
 128             c2 = v . v
 129             b = c1 / c2
 130             RETURN POINT(x1 + b * vx, y1 + b * vy)
 131         */
 132
 133         // v is multiplied below with a (possibly) FP-value, so should be in FP
 134         // For consistency we define w also in FP
 135         fp_vector_type v, w;
 136
 137         geometry::convert(p2, v);
 138         geometry::convert(p, w);
 139         subtract_point(v, p1);
 140         subtract_point(w, p1);
 141
 142         Strategy strategy;
 143         boost::ignore_unused_variable_warning(strategy);
 144
 145         calculation_type const zero = calculation_type();
 146         calculation_type const c1 = dot_product(w, v);
 147         if (c1 <= zero)
 148         {
 149             return strategy.apply(p, p1);
 150         }
 151         calculation_type const c2 = dot_product(v, v);
 152         if (c2 <= c1)
 153         {
 154             return strategy.apply(p, p2);
 155         }
 156
 157         // See above, c1 > 0 AND c2 > c1 so: c2 != 0
 158         calculation_type const b = c1 / c2;
 159
 160         fp_strategy_type fp_strategy
 161             = strategy::distance::services::get_similar
 162                 <
 163                     Strategy, Point, fp_point_type
 164                 >::apply(strategy);
 165
 166         fp_point_type projected;
 167         geometry::convert(p1, projected);
 168         multiply_value(v, b);
 169         add_point(projected, v);
 170
 171         //std::cout << "distance " << dsv(p) << " .. " << dsv(projected) << std::endl;
 172
 173         return fp_strategy.apply(p, projected);
 174     }
 175 };
 176
 177 #ifndef DOXYGEN_NO_STRATEGY_SPECIALIZATIONS
 178 namespace services
 179 {
 180
 181 template <typename Point, typename PointOfSegment, typename CalculationType, typename Strategy>
 182 struct tag<projected_point<Point, PointOfSegment, CalculationType, Strategy> >
 183 {
 184     typedef strategy_tag_distance_point_segment type;
 185 };
 186
 187
 188 template <typename Point, typename PointOfSegment, typename CalculationType, typename Strategy>
 189 struct return_type<projected_point<Point, PointOfSegment, CalculationType, Strategy> >
 190 {
 191     typedef typename projected_point<Point, PointOfSegment, CalculationType, Strategy>::calculation_type type;
 192 };
 193
 194 template <typename Point, typename PointOfSegment, typename CalculationType, typename Strategy>
 195 struct strategy_point_point<projected_point<Point, PointOfSegment, CalculationType, Strategy> >
 196 {
 197     typedef Strategy type;
 198 };
 199
 200
 201 template
 202 <
 203     typename Point,
 204     typename PointOfSegment,
 205     typename CalculationType,
 206     typename Strategy,
 207     typename P1,
 208     typename P2
 209 >
 210 struct similar_type<projected_point<Point, PointOfSegment, CalculationType, Strategy>, P1, P2>
 211 {
 212     typedef projected_point<P1, P2, CalculationType, Strategy> type;
 213 };
 214
 215
 216 template
 217 <
 218     typename Point,
 219     typename PointOfSegment,
 220     typename CalculationType,
 221     typename Strategy,
 222     typename P1,
 223     typename P2
 224 >
 225 struct get_similar<projected_point<Point, PointOfSegment, CalculationType, Strategy>, P1, P2>
 226 {
 227     static inline typename similar_type
 228         <
 229             projected_point<Point, PointOfSegment, CalculationType, Strategy>, P1, P2
 230         >::type apply(projected_point<Point, PointOfSegment, CalculationType, Strategy> const& )
 231     {
 232         return projected_point<P1, P2, CalculationType, Strategy>();
 233     }
 234 };
 235
 236
 237 template <typename Point, typename PointOfSegment, typename CalculationType, typename Strategy>
 238 struct comparable_type<projected_point<Point, PointOfSegment, CalculationType, Strategy> >
 239 {
 240     // Define a projected_point strategy with its underlying point-point-strategy
 241     // being comparable
 242     typedef projected_point
 243         <
 244             Point,
 245             PointOfSegment,
 246             CalculationType,
 247             typename comparable_type<Strategy>::type
 248         > type;
 249 };
 250
 251
 252 template <typename Point, typename PointOfSegment, typename CalculationType, typename Strategy>
 253 struct get_comparable<projected_point<Point, PointOfSegment, CalculationType, Strategy> >
 254 {
 255     typedef typename comparable_type
 256         <
 257             projected_point<Point, PointOfSegment, CalculationType, Strategy>
 258         >::type comparable_type;
 259 public :
 260     static inline comparable_type apply(projected_point<Point, PointOfSegment, CalculationType, Strategy> const& )
 261     {
 262         return comparable_type();
 263     }
 264 };
 265
 266
 267 template <typename Point, typename PointOfSegment, typename CalculationType, typename Strategy>
 268 struct result_from_distance<projected_point<Point, PointOfSegment, CalculationType, Strategy> >
 269 {
 270 private :
 271     typedef typename return_type<projected_point<Point, PointOfSegment, CalculationType, Strategy> >::type return_type;
 272 public :
 273     template <typename T>
 274     static inline return_type apply(projected_point<Point, PointOfSegment, CalculationType, Strategy> const& , T const& value)
 275     {
 276         Strategy s;
 277         return result_from_distance<Strategy>::apply(s, value);
 278     }
 279 };
 280
 281
 282 // Get default-strategy for point-segment distance calculation
 283 // while still have the possibility to specify point-point distance strategy (PPS)
 284 // It is used in algorithms/distance.hpp where users specify PPS for distance
 285 // of point-to-segment or point-to-linestring.
 286 // Convenient for geographic coordinate systems especially.
 287 template <typename Point, typename PointOfSegment, typename Strategy>
 288 struct default_strategy<segment_tag, Point, PointOfSegment, cartesian_tag, cartesian_tag, Strategy>
 289 {
 290     typedef strategy::distance::projected_point
 291     <
 292         Point,
 293         PointOfSegment,
 294         void,
 295         typename boost::mpl::if_
 296             <
 297                 boost::is_void<Strategy>,
 298                 typename default_strategy
 299                     <
 300                         point_tag, Point, PointOfSegment,
 301                         cartesian_tag, cartesian_tag
 302                     >::type,
 303                 Strategy
 304             >::type
 305     > type;
 306 };
 307
 308
 309 } // namespace services
 310 #endif // DOXYGEN_NO_STRATEGY_SPECIALIZATIONS
 311
 312
 313 }} // namespace strategy::distance
 314
 315
 316 }} // namespace boost::geometry
 317
 318
 319 #endif // BOOST_GEOMETRY_STRATEGIES_CARTESIAN_DISTANCE_PROJECTED_POINT_HPP