1 // Boost.Geometry (aka GGL, Generic Geometry Library) 2 3 // Copyright (c) 2007-2014 Barend Gehrels, Amsterdam, the Netherlands. 4 5 // This file was modified by Oracle on 2014. 6 // Modifications copyright (c) 2014, Oracle and/or its affiliates. 7 8 // Contributed and/or modified by Menelaos Karavelas, on behalf of Oracle 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_SPHERICAL_DISTANCE_CROSS_TRACK_HPP 15 #define BOOST_GEOMETRY_STRATEGIES_SPHERICAL_DISTANCE_CROSS_TRACK_HPP 16 17 #include <algorithm> 18 19 #include <boost/config.hpp> 20 #include <boost/concept_check.hpp> 21 #include <boost/mpl/if.hpp> 22 #include <boost/type_traits/is_void.hpp> 23 24 #include <boost/geometry/core/cs.hpp> 25 #include <boost/geometry/core/access.hpp> 26 #include <boost/geometry/core/radian_access.hpp> 27 #include <boost/geometry/core/tags.hpp> 28 29 #include <boost/geometry/algorithms/detail/course.hpp> 30 31 #include <boost/geometry/strategies/distance.hpp> 32 #include <boost/geometry/strategies/concepts/distance_concept.hpp> 33 #include <boost/geometry/strategies/spherical/distance_haversine.hpp> 34 35 #include <boost/geometry/util/math.hpp> 36 #include <boost/geometry/util/promote_floating_point.hpp> 37 #include <boost/geometry/util/select_calculation_type.hpp> 38 39 #ifdef BOOST_GEOMETRY_DEBUG_CROSS_TRACK 40 # include <boost/geometry/io/dsv/write.hpp> 41 #endif 42 43 44 namespace boost { namespace geometry 45 { 46 47 namespace strategy { namespace distance 48 { 49 50 51 namespace comparable 52 { 53 54 /* 55 Given a spherical segment AB and a point D, we are interested in 56 computing the distance of D from AB. This is usually known as the 57 cross track distance. 58 59 If the projection (along great circles) of the point D lies inside 60 the segment AB, then the distance (cross track error) XTD is given 61 by the formula (see http://williams.best.vwh.net/avform.htm#XTE): 62 63 XTD = asin( sin(dist_AD) * sin(crs_AD-crs_AB) ) 64 65 where dist_AD is the great circle distance between the points A and 66 B, and crs_AD, crs_AB is the course (bearing) between the points A, 67 D and A, B, respectively. 68 69 If the point D does not project inside the arc AB, then the distance 70 of D from AB is the minimum of the two distances dist_AD and dist_BD. 71 72 Our reference implementation for this procedure is listed below 73 (this was the old Boost.Geometry implementation of the cross track distance), 74 where: 75 * The member variable m_strategy is the underlying haversine strategy. 76 * p stands for the point D. 77 * sp1 stands for the segment endpoint A. 78 * sp2 stands for the segment endpoint B. 79 80 ================= reference implementation -- start ================= 81 82 return_type d1 = m_strategy.apply(sp1, p); 83 return_type d3 = m_strategy.apply(sp1, sp2); 84 85 if (geometry::math::equals(d3, 0.0)) 86 { 87 // "Degenerate" segment, return either d1 or d2 88 return d1; 89 } 90 91 return_type d2 = m_strategy.apply(sp2, p); 92 93 return_type crs_AD = geometry::detail::course<return_type>(sp1, p); 94 return_type crs_AB = geometry::detail::course<return_type>(sp1, sp2); 95 return_type crs_BA = crs_AB - geometry::math::pi<return_type>(); 96 return_type crs_BD = geometry::detail::course<return_type>(sp2, p); 97 return_type d_crs1 = crs_AD - crs_AB; 98 return_type d_crs2 = crs_BD - crs_BA; 99 100 // d1, d2, d3 are in principle not needed, only the sign matters 101 return_type projection1 = cos( d_crs1 ) * d1 / d3; 102 return_type projection2 = cos( d_crs2 ) * d2 / d3; 103 104 if (projection1 > 0.0 && projection2 > 0.0) 105 { 106 return_type XTD 107 = radius() * math::abs( asin( sin( d1 / radius() ) * sin( d_crs1 ) )); 108 109 // Return shortest distance, projected point on segment sp1-sp2 110 return return_type(XTD); 111 } 112 else 113 { 114 // Return shortest distance, project either on point sp1 or sp2 115 return return_type( (std::min)( d1 , d2 ) ); 116 } 117 118 ================= reference implementation -- end ================= 119 120 121 Motivation 122 ---------- 123 In what follows we develop a comparable version of the cross track 124 distance strategy, that meets the following goals: 125 * It is more efficient than the original cross track strategy (less 126 operations and less calls to mathematical functions). 127 * Distances using the comparable cross track strategy can not only 128 be compared with other distances using the same strategy, but also with 129 distances computed with the comparable version of the haversine strategy. 130 * It can serve as the basis for the computation of the cross track distance, 131 as it is more efficient to compute its comparable version and 132 transform that to the actual cross track distance, rather than 133 follow/use the reference implementation listed above. 134 135 Major idea 136 ---------- 137 The idea here is to use the comparable haversine strategy to compute 138 the distances d1, d2 and d3 in the above listing. Once we have done 139 that we need also to make sure that instead of returning XTD (as 140 computed above) that we return a distance CXTD that is compatible 141 with the comparable haversine distance. To achieve this CXTD must satisfy 142 the relation: 143 XTD = 2 * R * asin( sqrt(XTD) ) 144 where R is the sphere's radius. 145 146 Below we perform the mathematical analysis that show how to compute CXTD. 147 148 149 Mathematical analysis 150 --------------------- 151 Below we use the following trigonometric identities: 152 sin(2 * x) = 2 * sin(x) * cos(x) 153 cos(asin(x)) = sqrt(1 - x^2) 154 155 Observation: 156 The distance d1 needed when the projection of the point D is within the 157 segment must be the true distance. However, comparable::haversine<> 158 returns a comparable distance instead of the one needed. 159 To remedy this, we implicitly compute what is needed. 160 More precisely, we need to compute sin(true_d1): 161 162 sin(true_d1) = sin(2 * asin(sqrt(d1))) 163 = 2 * sin(asin(sqrt(d1)) * cos(asin(sqrt(d1))) 164 = 2 * sqrt(d1) * sqrt(1-(sqrt(d1))^2) 165 = 2 * sqrt(d1 - d1 * d1) 166 This relation is used below. 167 168 As we mentioned above the goal is to find CXTD (named "a" below for 169 brevity) such that ("b" below stands for "d1", and "c" for "d_crs1"): 170 171 2 * R * asin(sqrt(a)) == R * asin(2 * sqrt(b-b^2) * sin(c)) 172 173 Analysis: 174 2 * R * asin(sqrt(a)) == R * asin(2 * sqrt(b-b^2) * sin(c)) 175 <=> 2 * asin(sqrt(a)) == asin(sqrt(b-b^2) * sin(c)) 176 <=> sin(2 * asin(sqrt(a))) == 2 * sqrt(b-b^2) * sin(c) 177 <=> 2 * sin(asin(sqrt(a))) * cos(asin(sqrt(a))) == 2 * sqrt(b-b^2) * sin(c) 178 <=> 2 * sqrt(a) * sqrt(1-a) == 2 * sqrt(b-b^2) * sin(c) 179 <=> sqrt(a) * sqrt(1-a) == sqrt(b-b^2) * sin(c) 180 <=> sqrt(a-a^2) == sqrt(b-b^2) * sin(c) 181 <=> a-a^2 == (b-b^2) * (sin(c))^2 182 183 Consider the quadratic equation: x^2-x+p^2 == 0, 184 where p = sqrt(b-b^2) * sin(c); its discriminant is: 185 d = 1 - 4 * p^2 = 1 - 4 * (b-b^2) * (sin(c))^2 186 187 The two solutions are: 188 a_1 = (1 - sqrt(d)) / 2 189 a_2 = (1 + sqrt(d)) / 2 190 191 Which one to choose? 192 "a" refers to the distance (on the unit sphere) of D from the 193 supporting great circle Circ(A,B) of the segment AB. 194 The two different values for "a" correspond to the lengths of the two 195 arcs delimited D and the points of intersection of Circ(A,B) and the 196 great circle perperdicular to Circ(A,B) passing through D. 197 Clearly, the value we want is the smallest among these two distances, 198 hence the root we must choose is the smallest root among the two. 199 200 So the answer is: 201 CXTD = ( 1 - sqrt(1 - 4 * (b-b^2) * (sin(c))^2) ) / 2 202 203 Therefore, in order to implement the comparable version of the cross 204 track strategy we need to: 205 (1) Use the comparable version of the haversine strategy instead of 206 the non-comparable one. 207 (2) Instead of return XTD when D projects inside the segment AB, we 208 need to return CXTD, given by the following formula: 209 CXTD = ( 1 - sqrt(1 - 4 * (d1-d1^2) * (sin(d_crs1))^2) ) / 2; 210 211 212 Complexity Analysis 213 ------------------- 214 In the analysis that follows we refer to the actual implementation below. 215 In particular, instead of computing CXTD as above, we use the more 216 efficient (operation-wise) computation of CXTD shown here: 217 218 return_type sin_d_crs1 = sin(d_crs1); 219 return_type d1_x_sin = d1 * sin_d_crs1; 220 return_type d = d1_x_sin * (sin_d_crs1 - d1_x_sin); 221 return d / (0.5 + math::sqrt(0.25 - d)); 222 223 Notice that instead of computing: 224 0.5 - 0.5 * sqrt(1 - 4 * d) = 0.5 - sqrt(0.25 - d) 225 we use the following formula instead: 226 d / (0.5 + sqrt(0.25 - d)). 227 This is done for numerical robustness. The expression 0.5 - sqrt(0.25 - x) 228 has large numerical errors for values of x close to 0 (if using doubles 229 the error start to become large even when d is as large as 0.001). 230 To remedy that, we re-write 0.5 - sqrt(0.25 - x) as: 231 0.5 - sqrt(0.25 - d) 232 = (0.5 - sqrt(0.25 - d) * (0.5 - sqrt(0.25 - d)) / (0.5 + sqrt(0.25 - d)). 233 The numerator is the difference of two squares: 234 (0.5 - sqrt(0.25 - d) * (0.5 - sqrt(0.25 - d)) 235 = 0.5^2 - (sqrt(0.25 - d))^ = 0.25 - (0.25 - d) = d, 236 which gives the expression we use. 237 238 For the complexity analysis, we distinguish between two cases: 239 (A) The distance is realized between the point D and an 240 endpoint of the segment AB 241 242 Gains: 243 Since we are using comparable::haversine<> which is called 244 3 times, we gain: 245 -> 3 calls to sqrt 246 -> 3 calls to asin 247 -> 6 multiplications 248 249 Loses: None 250 251 So the net gain is: 252 -> 6 function calls (sqrt/asin) 253 -> 6 arithmetic operations 254 255 If we use comparable::cross_track<> to compute 256 cross_track<> we need to account for a call to sqrt, a call 257 to asin and 2 multiplications. In this case the net gain is: 258 -> 4 function calls (sqrt/asin) 259 -> 4 arithmetic operations 260 261 262 (B) The distance is realized between the point D and an 263 interior point of the segment AB 264 265 Gains: 266 Since we are using comparable::haversine<> which is called 267 3 times, we gain: 268 -> 3 calls to sqrt 269 -> 3 calls to asin 270 -> 6 multiplications 271 Also we gain the operations used to compute XTD: 272 -> 2 calls to sin 273 -> 1 call to asin 274 -> 1 call to abs 275 -> 2 multiplications 276 -> 1 division 277 So the total gains are: 278 -> 9 calls to sqrt/sin/asin 279 -> 1 call to abs 280 -> 8 multiplications 281 -> 1 division 282 283 Loses: 284 To compute a distance compatible with comparable::haversine<> 285 we need to perform a few more operations, namely: 286 -> 1 call to sin 287 -> 1 call to sqrt 288 -> 2 multiplications 289 -> 1 division 290 -> 1 addition 291 -> 2 subtractions 292 293 So roughly speaking the net gain is: 294 -> 8 fewer function calls and 3 fewer arithmetic operations 295 296 If we were to implement cross_track directly from the 297 comparable version (much like what haversine<> does using 298 comparable::haversine<>) we need additionally 299 -> 2 function calls (asin/sqrt) 300 -> 2 multiplications 301 302 So it pays off to re-implement cross_track<> to use 303 comparable::cross_track<>; in this case the net gain would be: 304 -> 6 function calls 305 -> 1 arithmetic operation 306 307 Summary/Conclusion 308 ------------------ 309 Following the mathematical and complexity analysis above, the 310 comparable cross track strategy (as implemented below) satisfies 311 all the goal mentioned in the beginning: 312 * It is more efficient than its non-comparable counter-part. 313 * Comparable distances using this new strategy can also be compared 314 with comparable distances computed with the comparable haversine 315 strategy. 316 * It turns out to be more efficient to compute the actual cross 317 track distance XTD by first computing CXTD, and then computing 318 XTD by means of the formula: 319 XTD = 2 * R * asin( sqrt(CXTD) ) 320 */ 321 322 template 323 < 324 typename CalculationType = void, 325 typename Strategy = comparable::haversine<double, CalculationType> 326 > 327 class cross_track 328 { 329 public : 330 template <typename Point, typename PointOfSegment> 331 struct return_type 332 : promote_floating_point 333 < 334 typename select_calculation_type 335 < 336 Point, 337 PointOfSegment, 338 CalculationType 339 >::type 340 > 341 {}; 342 343 typedef typename Strategy::radius_type radius_type; 344 cross_track()345 inline cross_track() 346 {} 347 cross_track(typename Strategy::radius_type const & r)348 explicit inline cross_track(typename Strategy::radius_type const& r) 349 : m_strategy(r) 350 {} 351 cross_track(Strategy const & s)352 inline cross_track(Strategy const& s) 353 : m_strategy(s) 354 {} 355 356 357 // It might be useful in the future 358 // to overload constructor with strategy info. 359 // crosstrack(...) {} 360 361 362 template <typename Point, typename PointOfSegment> 363 inline typename return_type<Point, PointOfSegment>::type apply(Point const & p,PointOfSegment const & sp1,PointOfSegment const & sp2) const364 apply(Point const& p, PointOfSegment const& sp1, PointOfSegment const& sp2) const 365 { 366 367 #if !defined(BOOST_MSVC) 368 BOOST_CONCEPT_ASSERT 369 ( 370 (concepts::PointDistanceStrategy<Strategy, Point, PointOfSegment>) 371 ); 372 #endif 373 374 typedef typename return_type<Point, PointOfSegment>::type return_type; 375 376 #ifdef BOOST_GEOMETRY_DEBUG_CROSS_TRACK 377 std::cout << "Course " << dsv(sp1) << " to " << dsv(p) << " " 378 << crs_AD * geometry::math::r2d<return_type>() << std::endl; 379 std::cout << "Course " << dsv(sp1) << " to " << dsv(sp2) << " " 380 << crs_AB * geometry::math::r2d<return_type>() << std::endl; 381 std::cout << "Course " << dsv(sp2) << " to " << dsv(p) << " " 382 << crs_BD * geometry::math::r2d << std::endl; 383 std::cout << "Projection AD-AB " << projection1 << " : " 384 << d_crs1 * geometry::math::r2d<return_type>() << std::endl; 385 std::cout << "Projection BD-BA " << projection2 << " : " 386 << d_crs2 * geometry::math::r2d<return_type>() << std::endl; 387 #endif 388 389 // http://williams.best.vwh.net/avform.htm#XTE 390 return_type d1 = m_strategy.apply(sp1, p); 391 return_type d3 = m_strategy.apply(sp1, sp2); 392 393 if (geometry::math::equals(d3, 0.0)) 394 { 395 // "Degenerate" segment, return either d1 or d2 396 return d1; 397 } 398 399 return_type d2 = m_strategy.apply(sp2, p); 400 401 return_type crs_AD = geometry::detail::course<return_type>(sp1, p); 402 return_type crs_AB = geometry::detail::course<return_type>(sp1, sp2); 403 return_type crs_BA = crs_AB - geometry::math::pi<return_type>(); 404 return_type crs_BD = geometry::detail::course<return_type>(sp2, p); 405 return_type d_crs1 = crs_AD - crs_AB; 406 return_type d_crs2 = crs_BD - crs_BA; 407 408 // d1, d2, d3 are in principle not needed, only the sign matters 409 return_type projection1 = cos( d_crs1 ) * d1 / d3; 410 return_type projection2 = cos( d_crs2 ) * d2 / d3; 411 412 if (projection1 > 0.0 && projection2 > 0.0) 413 { 414 #ifdef BOOST_GEOMETRY_DEBUG_CROSS_TRACK 415 return_type XTD = radius() * geometry::math::abs( asin( sin( d1 ) * sin( d_crs1 ) )); 416 417 std::cout << "Projection ON the segment" << std::endl; 418 std::cout << "XTD: " << XTD 419 << " d1: " << (d1 * radius()) 420 << " d2: " << (d2 * radius()) 421 << std::endl; 422 #endif 423 return_type const half(0.5); 424 return_type const quarter(0.25); 425 426 return_type sin_d_crs1 = sin(d_crs1); 427 /* 428 This is the straightforward obvious way to continue: 429 430 return_type discriminant 431 = 1.0 - 4.0 * (d1 - d1 * d1) * sin_d_crs1 * sin_d_crs1; 432 return 0.5 - 0.5 * math::sqrt(discriminant); 433 434 Below we optimize the number of arithmetic operations 435 and account for numerical robustness: 436 */ 437 return_type d1_x_sin = d1 * sin_d_crs1; 438 return_type d = d1_x_sin * (sin_d_crs1 - d1_x_sin); 439 return d / (half + math::sqrt(quarter - d)); 440 } 441 else 442 { 443 #ifdef BOOST_GEOMETRY_DEBUG_CROSS_TRACK 444 std::cout << "Projection OUTSIDE the segment" << std::endl; 445 #endif 446 447 // Return shortest distance, project either on point sp1 or sp2 448 return return_type( (std::min)( d1 , d2 ) ); 449 } 450 } 451 radius() const452 inline typename Strategy::radius_type radius() const 453 { return m_strategy.radius(); } 454 455 private : 456 Strategy m_strategy; 457 }; 458 459 } // namespace comparable 460 461 462 /*! 463 \brief Strategy functor for distance point to segment calculation 464 \ingroup strategies 465 \details Class which calculates the distance of a point to a segment, for points on a sphere or globe 466 \see http://williams.best.vwh.net/avform.htm 467 \tparam CalculationType \tparam_calculation 468 \tparam Strategy underlying point-point distance strategy, defaults to haversine 469 470 \qbk{ 471 [heading See also] 472 [link geometry.reference.algorithms.distance.distance_3_with_strategy distance (with strategy)] 473 } 474 475 */ 476 template 477 < 478 typename CalculationType = void, 479 typename Strategy = haversine<double, CalculationType> 480 > 481 class cross_track 482 { 483 public : 484 template <typename Point, typename PointOfSegment> 485 struct return_type 486 : promote_floating_point 487 < 488 typename select_calculation_type 489 < 490 Point, 491 PointOfSegment, 492 CalculationType 493 >::type 494 > 495 {}; 496 497 typedef typename Strategy::radius_type radius_type; 498 cross_track()499 inline cross_track() 500 {} 501 cross_track(typename Strategy::radius_type const & r)502 explicit inline cross_track(typename Strategy::radius_type const& r) 503 : m_strategy(r) 504 {} 505 cross_track(Strategy const & s)506 inline cross_track(Strategy const& s) 507 : m_strategy(s) 508 {} 509 510 511 // It might be useful in the future 512 // to overload constructor with strategy info. 513 // crosstrack(...) {} 514 515 516 template <typename Point, typename PointOfSegment> 517 inline typename return_type<Point, PointOfSegment>::type apply(Point const & p,PointOfSegment const & sp1,PointOfSegment const & sp2) const518 apply(Point const& p, PointOfSegment const& sp1, PointOfSegment const& sp2) const 519 { 520 521 #if !defined(BOOST_MSVC) 522 BOOST_CONCEPT_ASSERT 523 ( 524 (concepts::PointDistanceStrategy<Strategy, Point, PointOfSegment>) 525 ); 526 #endif 527 typedef typename return_type<Point, PointOfSegment>::type return_type; 528 typedef cross_track<CalculationType, Strategy> this_type; 529 530 typedef typename services::comparable_type 531 < 532 this_type 533 >::type comparable_type; 534 535 comparable_type cstrategy 536 = services::get_comparable<this_type>::apply(m_strategy); 537 538 return_type const a = cstrategy.apply(p, sp1, sp2); 539 return_type const c = return_type(2.0) * asin(math::sqrt(a)); 540 return c * radius(); 541 } 542 radius() const543 inline typename Strategy::radius_type radius() const 544 { return m_strategy.radius(); } 545 546 private : 547 548 Strategy m_strategy; 549 }; 550 551 552 553 #ifndef DOXYGEN_NO_STRATEGY_SPECIALIZATIONS 554 namespace services 555 { 556 557 template <typename CalculationType, typename Strategy> 558 struct tag<cross_track<CalculationType, Strategy> > 559 { 560 typedef strategy_tag_distance_point_segment type; 561 }; 562 563 564 template <typename CalculationType, typename Strategy, typename P, typename PS> 565 struct return_type<cross_track<CalculationType, Strategy>, P, PS> 566 : cross_track<CalculationType, Strategy>::template return_type<P, PS> 567 {}; 568 569 570 template <typename CalculationType, typename Strategy> 571 struct comparable_type<cross_track<CalculationType, Strategy> > 572 { 573 typedef comparable::cross_track 574 < 575 CalculationType, typename comparable_type<Strategy>::type 576 > type; 577 }; 578 579 580 template 581 < 582 typename CalculationType, 583 typename Strategy 584 > 585 struct get_comparable<cross_track<CalculationType, Strategy> > 586 { 587 typedef typename comparable_type 588 < 589 cross_track<CalculationType, Strategy> 590 >::type comparable_type; 591 public : 592 static inline comparable_type applyboost::geometry::strategy::distance::services::get_comparable593 apply(cross_track<CalculationType, Strategy> const& strategy) 594 { 595 return comparable_type(strategy.radius()); 596 } 597 }; 598 599 600 template 601 < 602 typename CalculationType, 603 typename Strategy, 604 typename P, 605 typename PS 606 > 607 struct result_from_distance<cross_track<CalculationType, Strategy>, P, PS> 608 { 609 private : 610 typedef typename cross_track 611 < 612 CalculationType, Strategy 613 >::template return_type<P, PS>::type return_type; 614 public : 615 template <typename T> 616 static inline return_type applyboost::geometry::strategy::distance::services::result_from_distance617 apply(cross_track<CalculationType, Strategy> const& , T const& distance) 618 { 619 return distance; 620 } 621 }; 622 623 624 // Specializations for comparable::cross_track 625 template <typename RadiusType, typename CalculationType> 626 struct tag<comparable::cross_track<RadiusType, CalculationType> > 627 { 628 typedef strategy_tag_distance_point_segment type; 629 }; 630 631 632 template 633 < 634 typename RadiusType, 635 typename CalculationType, 636 typename P, 637 typename PS 638 > 639 struct return_type<comparable::cross_track<RadiusType, CalculationType>, P, PS> 640 : comparable::cross_track 641 < 642 RadiusType, CalculationType 643 >::template return_type<P, PS> 644 {}; 645 646 647 template <typename RadiusType, typename CalculationType> 648 struct comparable_type<comparable::cross_track<RadiusType, CalculationType> > 649 { 650 typedef comparable::cross_track<RadiusType, CalculationType> type; 651 }; 652 653 654 template <typename RadiusType, typename CalculationType> 655 struct get_comparable<comparable::cross_track<RadiusType, CalculationType> > 656 { 657 private : 658 typedef comparable::cross_track<RadiusType, CalculationType> this_type; 659 public : applyboost::geometry::strategy::distance::services::get_comparable660 static inline this_type apply(this_type const& input) 661 { 662 return input; 663 } 664 }; 665 666 667 template 668 < 669 typename RadiusType, 670 typename CalculationType, 671 typename P, 672 typename PS 673 > 674 struct result_from_distance 675 < 676 comparable::cross_track<RadiusType, CalculationType>, P, PS 677 > 678 { 679 private : 680 typedef comparable::cross_track<RadiusType, CalculationType> strategy_type; 681 typedef typename return_type<strategy_type, P, PS>::type return_type; 682 public : 683 template <typename T> applyboost::geometry::strategy::distance::services::result_from_distance684 static inline return_type apply(strategy_type const& strategy, 685 T const& distance) 686 { 687 return_type const s 688 = sin( (distance / strategy.radius()) / return_type(2.0) ); 689 return s * s; 690 } 691 }; 692 693 694 695 /* 696 697 TODO: spherical polar coordinate system requires "get_as_radian_equatorial<>" 698 699 template <typename Point, typename PointOfSegment, typename Strategy> 700 struct default_strategy 701 < 702 segment_tag, Point, PointOfSegment, 703 spherical_polar_tag, spherical_polar_tag, 704 Strategy 705 > 706 { 707 typedef cross_track 708 < 709 void, 710 typename boost::mpl::if_ 711 < 712 boost::is_void<Strategy>, 713 typename default_strategy 714 < 715 point_tag, Point, PointOfSegment, 716 spherical_polar_tag, spherical_polar_tag 717 >::type, 718 Strategy 719 >::type 720 > type; 721 }; 722 */ 723 724 template <typename Point, typename PointOfSegment, typename Strategy> 725 struct default_strategy 726 < 727 point_tag, segment_tag, Point, PointOfSegment, 728 spherical_equatorial_tag, spherical_equatorial_tag, 729 Strategy 730 > 731 { 732 typedef cross_track 733 < 734 void, 735 typename boost::mpl::if_ 736 < 737 boost::is_void<Strategy>, 738 typename default_strategy 739 < 740 point_tag, point_tag, Point, PointOfSegment, 741 spherical_equatorial_tag, spherical_equatorial_tag 742 >::type, 743 Strategy 744 >::type 745 > type; 746 }; 747 748 749 template <typename PointOfSegment, typename Point, typename Strategy> 750 struct default_strategy 751 < 752 segment_tag, point_tag, PointOfSegment, Point, 753 spherical_equatorial_tag, spherical_equatorial_tag, 754 Strategy 755 > 756 { 757 typedef typename default_strategy 758 < 759 point_tag, segment_tag, Point, PointOfSegment, 760 spherical_equatorial_tag, spherical_equatorial_tag, 761 Strategy 762 >::type type; 763 }; 764 765 766 } // namespace services 767 #endif // DOXYGEN_NO_STRATEGY_SPECIALIZATIONS 768 769 }} // namespace strategy::distance 770 771 }} // namespace boost::geometry 772 773 #endif // BOOST_GEOMETRY_STRATEGIES_SPHERICAL_DISTANCE_CROSS_TRACK_HPP 774