proj_core/
ellipsoid.rs

1//! Reference ellipsoid definitions.
2
3use std::f64::consts::PI;
4
5/// A reference ellipsoid.
6#[derive(Debug, Clone, Copy)]
7pub struct Ellipsoid {
8    /// Semi-major axis (equatorial radius) in meters.
9    pub a: f64,
10    /// Flattening (f = (a - b) / a).
11    pub f: f64,
12}
13
14impl Ellipsoid {
15    /// Create an ellipsoid from semi-major axis and inverse flattening.
16    pub const fn from_a_rf(a: f64, rf: f64) -> Self {
17        Self { a, f: 1.0 / rf }
18    }
19
20    /// Create a sphere (flattening = 0).
21    pub const fn sphere(radius: f64) -> Self {
22        Self { a: radius, f: 0.0 }
23    }
24
25    /// Semi-minor axis (polar radius) in meters.
26    pub fn b(&self) -> f64 {
27        self.a * (1.0 - self.f)
28    }
29
30    /// First eccentricity squared: e² = 2f - f².
31    pub fn e2(&self) -> f64 {
32        2.0 * self.f - self.f * self.f
33    }
34
35    /// First eccentricity.
36    pub fn e(&self) -> f64 {
37        self.e2().sqrt()
38    }
39
40    /// Third flattening: n = f / (2 - f) = (a - b) / (a + b).
41    pub fn n(&self) -> f64 {
42        self.f / (2.0 - self.f)
43    }
44
45    /// Second eccentricity squared: e'² = e² / (1 - e²).
46    pub fn ep2(&self) -> f64 {
47        let e2 = self.e2();
48        e2 / (1.0 - e2)
49    }
50
51    /// Radius of curvature in the prime vertical at latitude phi (radians).
52    pub fn n_radius(&self, phi: f64) -> f64 {
53        let sin_phi = phi.sin();
54        self.a / (1.0 - self.e2() * sin_phi * sin_phi).sqrt()
55    }
56
57    /// Radius of curvature in the meridian at latitude phi (radians).
58    pub fn m_radius(&self, phi: f64) -> f64 {
59        let e2 = self.e2();
60        let sin_phi = phi.sin();
61        let denom = 1.0 - e2 * sin_phi * sin_phi;
62        self.a * (1.0 - e2) / denom.powf(1.5)
63    }
64}
65
66// Well-known ellipsoids.
67
68/// WGS 84 ellipsoid.
69pub const WGS84: Ellipsoid = Ellipsoid::from_a_rf(6378137.0, 298.257223563);
70
71/// GRS 1980 ellipsoid (used by NAD83).
72pub const GRS80: Ellipsoid = Ellipsoid::from_a_rf(6378137.0, 298.257222101);
73
74/// Clarke 1866 ellipsoid (used by NAD27).
75pub const CLARKE1866: Ellipsoid = Ellipsoid::from_a_rf(6378206.4, 294.978698214);
76
77/// International 1924 (Hayford) ellipsoid.
78pub const INTL1924: Ellipsoid = Ellipsoid::from_a_rf(6378388.0, 297.0);
79
80/// Bessel 1841 ellipsoid.
81pub const BESSEL1841: Ellipsoid = Ellipsoid::from_a_rf(6377397.155, 299.1528128);
82
83/// Krassowsky 1940 ellipsoid.
84pub const KRASSOWSKY: Ellipsoid = Ellipsoid::from_a_rf(6378245.0, 298.3);
85
86/// Airy 1830 ellipsoid (used by OSGB36).
87pub const AIRY1830: Ellipsoid = Ellipsoid::from_a_rf(6377563.396, 299.3249646);
88
89/// Unit sphere (radius = 1).
90pub const UNIT_SPHERE: Ellipsoid = Ellipsoid::sphere(1.0);
91
92/// Authalic sphere matching WGS84 surface area.
93pub const WGS84_SPHERE: Ellipsoid = Ellipsoid::sphere(6371007.181);
94
95/// Convert degrees to radians.
96pub fn deg_to_rad(deg: f64) -> f64 {
97    deg * PI / 180.0
98}
99
100/// Convert radians to degrees.
101pub fn rad_to_deg(rad: f64) -> f64 {
102    rad * 180.0 / PI
103}
104
105#[cfg(test)]
106mod tests {
107    use super::*;
108
109    #[test]
110    fn wgs84_basics() {
111        let e = WGS84;
112        assert!((e.a - 6378137.0).abs() < 1e-6);
113        assert!((e.b() - 6356752.314245).abs() < 1.0);
114        assert!((e.e2() - 0.00669437999014).abs() < 1e-11);
115    }
116
117    #[test]
118    fn sphere_has_zero_eccentricity() {
119        let s = UNIT_SPHERE;
120        assert_eq!(s.e2(), 0.0);
121        assert_eq!(s.b(), 1.0);
122    }
123
124    #[test]
125    fn deg_rad_roundtrip() {
126        let deg = 45.0;
127        assert!((rad_to_deg(deg_to_rad(deg)) - deg).abs() < 1e-12);
128    }
129}
proj_core/ellipsoid.rs

proj_core/
ellipsoid.rs
