1 // Boost.Geometry (aka GGL, Generic Geometry Library) 2 3 // Copyright (c) 2007-2015 Barend Gehrels, Amsterdam, the Netherlands. 4 // Copyright (c) 2008-2015 Bruno Lalande, Paris, France. 5 // Copyright (c) 2009-2015 Mateusz Loskot, London, UK. 6 7 // This file was modified by Oracle on 2015. 8 // Modifications copyright (c) 2015 Oracle and/or its affiliates. 9 10 // Contributed and/or modified by Menelaos Karavelas, on behalf of Oracle 11 12 // Parts of Boost.Geometry are redesigned from Geodan's Geographic Library 13 // (geolib/GGL), copyright (c) 1995-2010 Geodan, Amsterdam, the Netherlands. 14 15 // Use, modification and distribution is subject to the Boost Software License, 16 // Version 1.0. (See accompanying file LICENSE_1_0.txt or copy at 17 // http://www.boost.org/LICENSE_1_0.txt) 18 19 #ifndef BOOST_GEOMETRY_STRATEGIES_STRATEGY_TRANSFORM_HPP 20 #define BOOST_GEOMETRY_STRATEGIES_STRATEGY_TRANSFORM_HPP 21 22 #include <cstddef> 23 #include <cmath> 24 #include <functional> 25 26 #include <boost/numeric/conversion/cast.hpp> 27 28 #include <boost/geometry/algorithms/convert.hpp> 29 #include <boost/geometry/arithmetic/arithmetic.hpp> 30 #include <boost/geometry/core/access.hpp> 31 #include <boost/geometry/core/radian_access.hpp> 32 #include <boost/geometry/core/coordinate_dimension.hpp> 33 #include <boost/geometry/strategies/transform.hpp> 34 35 #include <boost/geometry/util/math.hpp> 36 #include <boost/geometry/util/promote_floating_point.hpp> 37 #include <boost/geometry/util/select_coordinate_type.hpp> 38 39 namespace boost { namespace geometry 40 { 41 42 namespace strategy { namespace transform 43 { 44 45 #ifndef DOXYGEN_NO_DETAIL 46 namespace detail 47 { 48 49 template 50 < 51 typename Src, typename Dst, 52 std::size_t D, std::size_t N, 53 template <typename> class F 54 > 55 struct transform_coordinates 56 { 57 template <typename T> transformboost::geometry::strategy::transform::detail::transform_coordinates58 static inline void transform(Src const& source, Dst& dest, T value) 59 { 60 typedef typename select_coordinate_type<Src, Dst>::type coordinate_type; 61 62 F<coordinate_type> function; 63 set<D>(dest, boost::numeric_cast<coordinate_type>(function(get<D>(source), value))); 64 transform_coordinates<Src, Dst, D + 1, N, F>::transform(source, dest, value); 65 } 66 }; 67 68 template 69 < 70 typename Src, typename Dst, 71 std::size_t N, 72 template <typename> class F 73 > 74 struct transform_coordinates<Src, Dst, N, N, F> 75 { 76 template <typename T> transformboost::geometry::strategy::transform::detail::transform_coordinates77 static inline void transform(Src const& , Dst& , T ) 78 { 79 } 80 }; 81 82 } // namespace detail 83 #endif // DOXYGEN_NO_DETAIL 84 85 86 /*! 87 \brief Transformation strategy to copy one point to another using assignment operator 88 \ingroup transform 89 \tparam P point type 90 */ 91 template <typename P> 92 struct copy_direct 93 { applyboost::geometry::strategy::transform::copy_direct94 inline bool apply(P const& p1, P& p2) const 95 { 96 p2 = p1; 97 return true; 98 } 99 }; 100 101 /*! 102 \brief Transformation strategy to do copy a point, copying per coordinate. 103 \ingroup transform 104 \tparam P1 first point type 105 \tparam P2 second point type 106 */ 107 template <typename P1, typename P2> 108 struct copy_per_coordinate 109 { applyboost::geometry::strategy::transform::copy_per_coordinate110 inline bool apply(P1 const& p1, P2& p2) const 111 { 112 // Defensive check, dimensions are equal, selected by specialization 113 assert_dimension_equal<P1, P2>(); 114 115 geometry::convert(p1, p2); 116 return true; 117 } 118 }; 119 120 121 /*! 122 \brief Transformation strategy to go from degree to radian and back 123 \ingroup transform 124 \tparam P1 first point type 125 \tparam P2 second point type 126 \tparam F additional functor to divide or multiply with d2r 127 */ 128 template <typename P1, typename P2, template <typename> class F> 129 struct degree_radian_vv 130 { applyboost::geometry::strategy::transform::degree_radian_vv131 inline bool apply(P1 const& p1, P2& p2) const 132 { 133 // Spherical coordinates always have 2 coordinates measured in angles 134 // The optional third one is distance/height, provided in another strategy 135 // Polar coordinates having one angle, will be also in another strategy 136 assert_dimension<P1, 2>(); 137 assert_dimension<P2, 2>(); 138 139 typedef typename promote_floating_point 140 < 141 typename select_coordinate_type<P1, P2>::type 142 >::type calculation_type; 143 144 detail::transform_coordinates 145 < 146 P1, P2, 0, 2, F 147 >::transform(p1, p2, math::d2r<calculation_type>()); 148 return true; 149 } 150 }; 151 152 template <typename P1, typename P2, template <typename> class F> 153 struct degree_radian_vv_3 154 { applyboost::geometry::strategy::transform::degree_radian_vv_3155 inline bool apply(P1 const& p1, P2& p2) const 156 { 157 assert_dimension<P1, 3>(); 158 assert_dimension<P2, 3>(); 159 160 typedef typename promote_floating_point 161 < 162 typename select_coordinate_type<P1, P2>::type 163 >::type calculation_type; 164 165 detail::transform_coordinates 166 < 167 P1, P2, 0, 2, F 168 >::transform(p1, p2, math::d2r<calculation_type>()); 169 170 // Copy height or other third dimension 171 set<2>(p2, get<2>(p1)); 172 return true; 173 } 174 }; 175 176 177 #ifndef DOXYGEN_NO_DETAIL 178 namespace detail 179 { 180 181 /// Helper function for conversion, phi/theta are in radians 182 template <typename P, typename T, typename R> spherical_polar_to_cartesian(T phi,T theta,R r,P & p)183 inline void spherical_polar_to_cartesian(T phi, T theta, R r, P& p) 184 { 185 assert_dimension<P, 3>(); 186 187 // http://en.wikipedia.org/wiki/List_of_canonical_coordinate_transformations#From_spherical_coordinates 188 // http://www.vias.org/comp_geometry/math_coord_convert_3d.htm 189 // https://moodle.polymtl.ca/file.php/1183/Autres_Documents/Derivation_for_Spherical_Co-ordinates.pdf 190 // http://en.citizendium.org/wiki/Spherical_polar_coordinates 191 192 // Phi = first, theta is second, r is third, see documentation on cs::spherical 193 194 // (calculations are splitted to implement ttmath) 195 196 T r_sin_theta = r; 197 T r_cos_theta = r; 198 r_sin_theta *= sin(theta); 199 r_cos_theta *= cos(theta); 200 201 set<0>(p, r_sin_theta * cos(phi)); 202 set<1>(p, r_sin_theta * sin(phi)); 203 set<2>(p, r_cos_theta); 204 } 205 206 /// Helper function for conversion, lambda/delta (lon lat) are in radians 207 template <typename P, typename T, typename R> spherical_equatorial_to_cartesian(T lambda,T delta,R r,P & p)208 inline void spherical_equatorial_to_cartesian(T lambda, T delta, R r, P& p) 209 { 210 assert_dimension<P, 3>(); 211 212 // http://mathworld.wolfram.com/GreatCircle.html 213 // http://www.spenvis.oma.be/help/background/coortran/coortran.html WRONG 214 215 T r_cos_delta = r; 216 T r_sin_delta = r; 217 r_cos_delta *= cos(delta); 218 r_sin_delta *= sin(delta); 219 220 set<0>(p, r_cos_delta * cos(lambda)); 221 set<1>(p, r_cos_delta * sin(lambda)); 222 set<2>(p, r_sin_delta); 223 } 224 225 226 /// Helper function for conversion 227 template <typename P, typename T> cartesian_to_spherical2(T x,T y,T z,P & p)228 inline bool cartesian_to_spherical2(T x, T y, T z, P& p) 229 { 230 assert_dimension<P, 2>(); 231 232 // http://en.wikipedia.org/wiki/List_of_canonical_coordinate_transformations#From_Cartesian_coordinates 233 234 #if defined(BOOST_GEOMETRY_TRANSFORM_CHECK_UNIT_SPHERE) 235 // TODO: MAYBE ONLY IF TO BE CHECKED? 236 T const r = /*sqrt not necessary, sqrt(1)=1*/ (x * x + y * y + z * z); 237 238 // Unit sphere, so r should be 1 239 if (geometry::math::abs(r - 1.0) > T(1e-6)) 240 { 241 return false; 242 } 243 // end todo 244 #endif 245 246 set_from_radian<0>(p, atan2(y, x)); 247 set_from_radian<1>(p, acos(z)); 248 return true; 249 } 250 251 template <typename P, typename T> cartesian_to_spherical_equatorial2(T x,T y,T z,P & p)252 inline bool cartesian_to_spherical_equatorial2(T x, T y, T z, P& p) 253 { 254 assert_dimension<P, 2>(); 255 256 set_from_radian<0>(p, atan2(y, x)); 257 set_from_radian<1>(p, asin(z)); 258 return true; 259 } 260 261 262 template <typename P, typename T> cartesian_to_spherical3(T x,T y,T z,P & p)263 inline bool cartesian_to_spherical3(T x, T y, T z, P& p) 264 { 265 assert_dimension<P, 3>(); 266 267 // http://en.wikipedia.org/wiki/List_of_canonical_coordinate_transformations#From_Cartesian_coordinates 268 T const r = math::sqrt(x * x + y * y + z * z); 269 set<2>(p, r); 270 set_from_radian<0>(p, atan2(y, x)); 271 if (r > 0.0) 272 { 273 set_from_radian<1>(p, acos(z / r)); 274 return true; 275 } 276 return false; 277 } 278 279 template <typename P, typename T> cartesian_to_spherical_equatorial3(T x,T y,T z,P & p)280 inline bool cartesian_to_spherical_equatorial3(T x, T y, T z, P& p) 281 { 282 assert_dimension<P, 3>(); 283 284 // http://en.wikipedia.org/wiki/List_of_canonical_coordinate_transformations#From_Cartesian_coordinates 285 T const r = math::sqrt(x * x + y * y + z * z); 286 set<2>(p, r); 287 set_from_radian<0>(p, atan2(y, x)); 288 if (r > 0.0) 289 { 290 set_from_radian<1>(p, asin(z / r)); 291 return true; 292 } 293 return false; 294 } 295 296 } // namespace detail 297 #endif // DOXYGEN_NO_DETAIL 298 299 300 /*! 301 \brief Transformation strategy for 2D spherical (phi,theta) to 3D cartesian (x,y,z) 302 \details on Unit sphere 303 \ingroup transform 304 \tparam P1 first point type 305 \tparam P2 second point type 306 */ 307 template <typename P1, typename P2> 308 struct from_spherical_polar_2_to_cartesian_3 309 { applyboost::geometry::strategy::transform::from_spherical_polar_2_to_cartesian_3310 inline bool apply(P1 const& p1, P2& p2) const 311 { 312 assert_dimension<P1, 2>(); 313 detail::spherical_polar_to_cartesian(get_as_radian<0>(p1), get_as_radian<1>(p1), 1.0, p2); 314 return true; 315 } 316 }; 317 318 template <typename P1, typename P2> 319 struct from_spherical_equatorial_2_to_cartesian_3 320 { applyboost::geometry::strategy::transform::from_spherical_equatorial_2_to_cartesian_3321 inline bool apply(P1 const& p1, P2& p2) const 322 { 323 assert_dimension<P1, 2>(); 324 detail::spherical_equatorial_to_cartesian(get_as_radian<0>(p1), get_as_radian<1>(p1), 1.0, p2); 325 return true; 326 } 327 }; 328 329 330 /*! 331 \brief Transformation strategy for 3D spherical (phi,theta,r) to 3D cartesian (x,y,z) 332 \ingroup transform 333 \tparam P1 first point type 334 \tparam P2 second point type 335 */ 336 template <typename P1, typename P2> 337 struct from_spherical_polar_3_to_cartesian_3 338 { applyboost::geometry::strategy::transform::from_spherical_polar_3_to_cartesian_3339 inline bool apply(P1 const& p1, P2& p2) const 340 { 341 assert_dimension<P1, 3>(); 342 detail::spherical_polar_to_cartesian( 343 get_as_radian<0>(p1), get_as_radian<1>(p1), get<2>(p1), p2); 344 return true; 345 } 346 }; 347 348 template <typename P1, typename P2> 349 struct from_spherical_equatorial_3_to_cartesian_3 350 { applyboost::geometry::strategy::transform::from_spherical_equatorial_3_to_cartesian_3351 inline bool apply(P1 const& p1, P2& p2) const 352 { 353 assert_dimension<P1, 3>(); 354 detail::spherical_equatorial_to_cartesian( 355 get_as_radian<0>(p1), get_as_radian<1>(p1), get<2>(p1), p2); 356 return true; 357 } 358 }; 359 360 361 /*! 362 \brief Transformation strategy for 3D cartesian (x,y,z) to 2D spherical (phi,theta) 363 \details on Unit sphere 364 \ingroup transform 365 \tparam P1 first point type 366 \tparam P2 second point type 367 \note If x,y,z point is not lying on unit sphere, transformation will return false 368 */ 369 template <typename P1, typename P2> 370 struct from_cartesian_3_to_spherical_polar_2 371 { applyboost::geometry::strategy::transform::from_cartesian_3_to_spherical_polar_2372 inline bool apply(P1 const& p1, P2& p2) const 373 { 374 assert_dimension<P1, 3>(); 375 return detail::cartesian_to_spherical2(get<0>(p1), get<1>(p1), get<2>(p1), p2); 376 } 377 }; 378 379 template <typename P1, typename P2> 380 struct from_cartesian_3_to_spherical_equatorial_2 381 { applyboost::geometry::strategy::transform::from_cartesian_3_to_spherical_equatorial_2382 inline bool apply(P1 const& p1, P2& p2) const 383 { 384 assert_dimension<P1, 3>(); 385 return detail::cartesian_to_spherical_equatorial2(get<0>(p1), get<1>(p1), get<2>(p1), p2); 386 } 387 }; 388 389 390 /*! 391 \brief Transformation strategy for 3D cartesian (x,y,z) to 3D spherical (phi,theta,r) 392 \ingroup transform 393 \tparam P1 first point type 394 \tparam P2 second point type 395 */ 396 template <typename P1, typename P2> 397 struct from_cartesian_3_to_spherical_polar_3 398 { applyboost::geometry::strategy::transform::from_cartesian_3_to_spherical_polar_3399 inline bool apply(P1 const& p1, P2& p2) const 400 { 401 assert_dimension<P1, 3>(); 402 return detail::cartesian_to_spherical3(get<0>(p1), get<1>(p1), get<2>(p1), p2); 403 } 404 }; 405 406 template <typename P1, typename P2> 407 struct from_cartesian_3_to_spherical_equatorial_3 408 { applyboost::geometry::strategy::transform::from_cartesian_3_to_spherical_equatorial_3409 inline bool apply(P1 const& p1, P2& p2) const 410 { 411 assert_dimension<P1, 3>(); 412 return detail::cartesian_to_spherical_equatorial3(get<0>(p1), get<1>(p1), get<2>(p1), p2); 413 } 414 }; 415 416 #ifndef DOXYGEN_NO_STRATEGY_SPECIALIZATIONS 417 418 namespace services 419 { 420 421 /// Specialization for same coordinate system family, same system, same dimension, same point type, can be copied 422 template <typename CoordSysTag, typename CoordSys, std::size_t D, typename P> 423 struct default_strategy<CoordSysTag, CoordSysTag, CoordSys, CoordSys, D, D, P, P> 424 { 425 typedef copy_direct<P> type; 426 }; 427 428 /// Specialization for same coordinate system family and system, same dimension, different point type, copy per coordinate 429 template <typename CoordSysTag, typename CoordSys, std::size_t D, typename P1, typename P2> 430 struct default_strategy<CoordSysTag, CoordSysTag, CoordSys, CoordSys, D, D, P1, P2> 431 { 432 typedef copy_per_coordinate<P1, P2> type; 433 }; 434 435 /// Specialization to transform from degree to radian for any coordinate system / point type combination 436 template <typename CoordSysTag, template<typename> class CoordSys, typename P1, typename P2> 437 struct default_strategy<CoordSysTag, CoordSysTag, CoordSys<degree>, CoordSys<radian>, 2, 2, P1, P2> 438 { 439 typedef degree_radian_vv<P1, P2, std::multiplies> type; 440 }; 441 442 /// Specialization to transform from radian to degree for any coordinate system / point type combination 443 template <typename CoordSysTag, template<typename> class CoordSys, typename P1, typename P2> 444 struct default_strategy<CoordSysTag, CoordSysTag, CoordSys<radian>, CoordSys<degree>, 2, 2, P1, P2> 445 { 446 typedef degree_radian_vv<P1, P2, std::divides> type; 447 }; 448 449 450 /// Specialization degree->radian in 3D 451 template <typename CoordSysTag, template<typename> class CoordSys, typename P1, typename P2> 452 struct default_strategy<CoordSysTag, CoordSysTag, CoordSys<degree>, CoordSys<radian>, 3, 3, P1, P2> 453 { 454 typedef degree_radian_vv_3<P1, P2, std::multiplies> type; 455 }; 456 457 /// Specialization radian->degree in 3D 458 template <typename CoordSysTag, template<typename> class CoordSys, typename P1, typename P2> 459 struct default_strategy<CoordSysTag, CoordSysTag, CoordSys<radian>, CoordSys<degree>, 3, 3, P1, P2> 460 { 461 typedef degree_radian_vv_3<P1, P2, std::divides> type; 462 }; 463 464 /// Specialization to transform from unit sphere(phi,theta) to XYZ 465 template <typename CoordSys1, typename CoordSys2, typename P1, typename P2> 466 struct default_strategy<spherical_polar_tag, cartesian_tag, CoordSys1, CoordSys2, 2, 3, P1, P2> 467 { 468 typedef from_spherical_polar_2_to_cartesian_3<P1, P2> type; 469 }; 470 471 /// Specialization to transform from sphere(phi,theta,r) to XYZ 472 template <typename CoordSys1, typename CoordSys2, typename P1, typename P2> 473 struct default_strategy<spherical_polar_tag, cartesian_tag, CoordSys1, CoordSys2, 3, 3, P1, P2> 474 { 475 typedef from_spherical_polar_3_to_cartesian_3<P1, P2> type; 476 }; 477 478 template <typename CoordSys1, typename CoordSys2, typename P1, typename P2> 479 struct default_strategy<spherical_equatorial_tag, cartesian_tag, CoordSys1, CoordSys2, 2, 3, P1, P2> 480 { 481 typedef from_spherical_equatorial_2_to_cartesian_3<P1, P2> type; 482 }; 483 484 template <typename CoordSys1, typename CoordSys2, typename P1, typename P2> 485 struct default_strategy<spherical_equatorial_tag, cartesian_tag, CoordSys1, CoordSys2, 3, 3, P1, P2> 486 { 487 typedef from_spherical_equatorial_3_to_cartesian_3<P1, P2> type; 488 }; 489 490 /// Specialization to transform from XYZ to unit sphere(phi,theta) 491 template <typename CoordSys1, typename CoordSys2, typename P1, typename P2> 492 struct default_strategy<cartesian_tag, spherical_polar_tag, CoordSys1, CoordSys2, 3, 2, P1, P2> 493 { 494 typedef from_cartesian_3_to_spherical_polar_2<P1, P2> type; 495 }; 496 497 template <typename CoordSys1, typename CoordSys2, typename P1, typename P2> 498 struct default_strategy<cartesian_tag, spherical_equatorial_tag, CoordSys1, CoordSys2, 3, 2, P1, P2> 499 { 500 typedef from_cartesian_3_to_spherical_equatorial_2<P1, P2> type; 501 }; 502 503 /// Specialization to transform from XYZ to sphere(phi,theta,r) 504 template <typename CoordSys1, typename CoordSys2, typename P1, typename P2> 505 struct default_strategy<cartesian_tag, spherical_polar_tag, CoordSys1, CoordSys2, 3, 3, P1, P2> 506 { 507 typedef from_cartesian_3_to_spherical_polar_3<P1, P2> type; 508 }; 509 template <typename CoordSys1, typename CoordSys2, typename P1, typename P2> 510 struct default_strategy<cartesian_tag, spherical_equatorial_tag, CoordSys1, CoordSys2, 3, 3, P1, P2> 511 { 512 typedef from_cartesian_3_to_spherical_equatorial_3<P1, P2> type; 513 }; 514 515 516 } // namespace services 517 518 519 #endif // DOXYGEN_NO_STRATEGY_SPECIALIZATIONS 520 521 522 }} // namespace strategy::transform 523 524 525 }} // namespace boost::geometry 526 527 #endif // BOOST_GEOMETRY_STRATEGIES_STRATEGY_TRANSFORM_HPP 528