Skip to main content

geometry_cs/
spheroidal.rs

1//! Canonical longitude/latitude intervals on spherical coordinate systems.
2//!
3//! Mirrors `boost::geometry::math::normalize_spheroidal_box_coordinates`
4//! from `util/normalize_spheroidal_box_coordinates.hpp`. The Rust coordinate
5//! systems expose only Boost's equatorial convention, so the polar-latitude
6//! conversion template parameter has no analogue here.
7
8use core::ops::{Add, Neg, Rem, Sub};
9
10use crate::{AngleUnit, Degree, Radian};
11
12/// Scalar operations needed to canonicalize spheroidal angles.
13///
14/// Mirrors the arithmetic requirements on `CoordinateType` in
15/// `util/normalize_spheroidal_coordinates.hpp:205-250`. Unlike trigonometric
16/// strategies, normalization also supports signed integer degree values.
17pub trait SpheroidalScalar:
18    Copy
19    + PartialOrd
20    + Add<Output = Self>
21    + Sub<Output = Self>
22    + Neg<Output = Self>
23    + Rem<Output = Self>
24{
25    /// Additive identity used for canonical pole longitudes.
26    const ZERO: Self;
27
28    /// Absolute value used to detect full longitude bands and poles.
29    #[must_use]
30    fn abs(self) -> Self;
31}
32
33macro_rules! impl_spheroidal_float {
34    ($($scalar:ty),* $(,)?) => {
35        $(
36            impl SpheroidalScalar for $scalar {
37                const ZERO: Self = 0.0;
38
39                #[inline]
40                fn abs(self) -> Self {
41                    self.abs()
42                }
43            }
44        )*
45    };
46}
47
48impl_spheroidal_float!(f32, f64);
49
50macro_rules! impl_spheroidal_integer {
51    ($($scalar:ty),* $(,)?) => {
52        $(
53            impl SpheroidalScalar for $scalar {
54                const ZERO: Self = 0;
55
56                #[inline]
57                fn abs(self) -> Self {
58                    self.abs()
59                }
60            }
61        )*
62    };
63}
64
65impl_spheroidal_integer!(i16, i32, i64, i128, isize);
66
67/// Canonical angular constants for a spheroidal unit and scalar pair.
68///
69/// Mirrors `constants_on_spheroid<CoordinateType, Units>` from
70/// `util/normalize_spheroidal_coordinates.hpp:34-145`.
71pub trait SpheroidalUnits<T: SpheroidalScalar>: AngleUnit {
72    /// One complete longitude turn.
73    const PERIOD: T;
74    /// The positive antimeridian.
75    const HALF_PERIOD: T;
76    /// The positive pole latitude in the equatorial convention.
77    const QUARTER_PERIOD: T;
78}
79
80macro_rules! impl_degree_float_constants {
81    ($($scalar:ty),* $(,)?) => {
82        $(
83            impl SpheroidalUnits<$scalar> for Degree {
84                const PERIOD: $scalar = 360.0;
85                const HALF_PERIOD: $scalar = 180.0;
86                const QUARTER_PERIOD: $scalar = 90.0;
87            }
88        )*
89    };
90}
91
92impl_degree_float_constants!(f32, f64);
93
94macro_rules! impl_degree_integer_constants {
95    ($($scalar:ty),* $(,)?) => {
96        $(
97            impl SpheroidalUnits<$scalar> for Degree {
98                const PERIOD: $scalar = 360;
99                const HALF_PERIOD: $scalar = 180;
100                const QUARTER_PERIOD: $scalar = 90;
101            }
102        )*
103    };
104}
105
106impl_degree_integer_constants!(i16, i32, i64, i128, isize);
107
108impl SpheroidalUnits<f32> for Radian {
109    const PERIOD: f32 = core::f32::consts::TAU;
110    const HALF_PERIOD: f32 = core::f32::consts::PI;
111    const QUARTER_PERIOD: f32 = core::f32::consts::FRAC_PI_2;
112}
113
114impl SpheroidalUnits<f64> for Radian {
115    const PERIOD: f64 = core::f64::consts::TAU;
116    const HALF_PERIOD: f64 = core::f64::consts::PI;
117    const QUARTER_PERIOD: f64 = core::f64::consts::FRAC_PI_2;
118}
119
120/// Normalize the two corners of an equatorial spheroidal box in place.
121///
122/// Mirrors `boost::geometry::math::normalize_spheroidal_box_coordinates`
123/// from `util/normalize_spheroidal_box_coordinates.hpp:113-145`.
124/// Longitudes are canonicalized to one increasing interval, an input spanning
125/// a complete turn becomes `[-half_period, half_period]`, and a box collapsed
126/// at either pole receives canonical zero longitudes.
127#[inline]
128pub fn normalize_spheroidal_box_coordinates<U, T>(
129    longitude1: &mut T,
130    latitude1: &mut T,
131    longitude2: &mut T,
132    latitude2: &mut T,
133) where
134    U: SpheroidalUnits<T>,
135    T: SpheroidalScalar,
136{
137    let band = (*longitude1 - *longitude2).abs() >= U::PERIOD;
138
139    normalize_longitude::<U, T>(longitude1);
140    normalize_longitude::<U, T>(longitude2);
141
142    let south = -U::QUARTER_PERIOD;
143    if (*latitude1 == south && *latitude2 == south)
144        || (*latitude1 == U::QUARTER_PERIOD && *latitude2 == U::QUARTER_PERIOD)
145    {
146        *longitude1 = T::ZERO;
147        *longitude2 = T::ZERO;
148    } else if band {
149        *longitude1 = -U::HALF_PERIOD;
150        *longitude2 = U::HALF_PERIOD;
151    } else if *longitude1 > *longitude2 {
152        *longitude2 = *longitude2 + U::PERIOD;
153    }
154}
155
156#[inline]
157fn normalize_longitude<U, T>(longitude: &mut T)
158where
159    U: SpheroidalUnits<T>,
160    T: SpheroidalScalar,
161{
162    if longitude.abs() == U::HALF_PERIOD {
163        *longitude = U::HALF_PERIOD;
164    } else if *longitude > U::HALF_PERIOD {
165        *longitude = (*longitude + U::HALF_PERIOD) % U::PERIOD - U::HALF_PERIOD;
166        if *longitude == -U::HALF_PERIOD {
167            *longitude = U::HALF_PERIOD;
168        }
169    } else if *longitude < -U::HALF_PERIOD {
170        *longitude = (*longitude - U::HALF_PERIOD) % U::PERIOD + U::HALF_PERIOD;
171    }
172}