geometry_strategy/
destination.rs

1//! Point-at-bearing-and-distance strategies.
2//!
3//! This is the strategy-facing wrapper around Boost's direct geodesic formulas.
4//! Bearings and direct-formula inputs are radians; output coordinates are
5//! converted back to the angular unit carried by the input point.
6
7use geometry_cs::{CoordinateSystem, GeographicFamily, SphericalFamily};
8use geometry_trait::Point;
9
10#[cfg(feature = "std")]
11use crate::normalise::{HasAngularUnits, lonlat_radians};
12#[cfg(feature = "std")]
13use geometry_cs::AngleUnit;
14#[cfg(feature = "std")]
15use geometry_model::Point2D;
16
17/// Strategy computing the endpoint reached from a point, bearing, and distance.
18pub trait DestinationStrategy<P: Point> {
19    /// Destination point type.
20    type Output: Point;
21
22    /// Compute the destination. `bearing` is measured clockwise from north in
23    /// radians; `distance` uses the strategy's radius or spheroid units.
24    fn destination(&self, origin: &P, bearing: f64, distance: f64) -> Self::Output;
25}
26
27/// Select the default destination strategy for a coordinate-system family.
28pub trait DefaultDestination<Family> {
29    /// Default strategy type.
30    type Strategy: Default;
31}
32
33impl DefaultDestination<SphericalFamily> for SphericalFamily {
34    type Strategy = crate::spherical::Haversine;
35}
36
37impl DefaultDestination<GeographicFamily> for GeographicFamily {
38    type Strategy = crate::geographic::KarneyDirect;
39}
40
41/// Default destination strategy for a point type.
42pub type DefaultDestinationStrategy<P> =
43    <<<P as Point>::Cs as CoordinateSystem>::Family as DefaultDestination<
44        <<P as Point>::Cs as CoordinateSystem>::Family,
45    >>::Strategy;
46
47#[cfg(feature = "std")]
48impl<P> DestinationStrategy<P> for crate::spherical::Haversine
49where
50    P: Point<Scalar = f64>,
51    P::Cs: HasAngularUnits,
52    <P::Cs as CoordinateSystem>::Family: geometry_tag::SameAs<SphericalFamily>,
53{
54    type Output = Point2D<f64, P::Cs>;
55
56    fn destination(&self, origin: &P, bearing: f64, distance: f64) -> Self::Output {
57        let (longitude, latitude) = lonlat_radians(origin);
58        let angular_distance = distance / self.radius;
59        let sin_latitude = latitude.sin();
60        let cos_latitude = latitude.cos();
61        let sin_distance = angular_distance.sin();
62        let cos_distance = angular_distance.cos();
63        let latitude2 =
64            (sin_latitude * cos_distance + cos_latitude * sin_distance * bearing.cos()).asin();
65        let longitude2 = longitude
66            + (bearing.sin() * sin_distance * cos_latitude)
67                .atan2(cos_distance - sin_latitude * latitude2.sin());
68        point_from_radians::<P>(normalize_longitude(longitude2), latitude2)
69    }
70}
71
72macro_rules! impl_geographic_destination {
73    ($strategy:ty) => {
74        #[cfg(feature = "std")]
75        impl<P> DestinationStrategy<P> for $strategy
76        where
77            P: Point<Scalar = f64>,
78            P::Cs: HasAngularUnits,
79            <P::Cs as CoordinateSystem>::Family: geometry_tag::SameAs<GeographicFamily>,
80        {
81            type Output = Point2D<f64, P::Cs>;
82
83            fn destination(&self, origin: &P, bearing: f64, distance: f64) -> Self::Output {
84                let (longitude, latitude) = lonlat_radians(origin);
85                let result = self.apply(longitude, latitude, distance, bearing);
86                point_from_radians::<P>(result.lon2, result.lat2)
87            }
88        }
89    };
90}
91
92impl_geographic_destination!(crate::geographic::KarneyDirect);
93impl_geographic_destination!(crate::geographic::ThomasDirect);
94impl_geographic_destination!(crate::geographic::VincentyDirect);
95
96#[cfg(feature = "std")]
97fn point_from_radians<P>(longitude: f64, latitude: f64) -> Point2D<f64, P::Cs>
98where
99    P: Point<Scalar = f64>,
100    P::Cs: HasAngularUnits,
101{
102    type Units<P> = <<P as Point>::Cs as HasAngularUnits>::Units;
103    Point2D::new(
104        Units::<P>::from_radians(longitude),
105        Units::<P>::from_radians(latitude),
106    )
107}
108
109#[cfg(feature = "std")]
110fn normalize_longitude(longitude: f64) -> f64 {
111    (longitude + core::f64::consts::PI).rem_euclid(core::f64::consts::TAU) - core::f64::consts::PI
112}
geometry_strategy/destination.rs

geometry_strategy/
destination.rs
