Skip to main content

proj_core/
datum.rs

1use crate::ellipsoid::{self, Ellipsoid};
2
3/// A geodetic datum, defined by a reference ellipsoid and its relationship to WGS84.
4#[derive(Debug, Clone, Copy)]
5pub struct Datum {
6    /// The reference ellipsoid.
7    pub ellipsoid: Ellipsoid,
8    /// Explicit relationship from this datum to WGS84.
9    pub to_wgs84: DatumToWgs84,
10}
11
12impl Datum {
13    /// Returns true if this datum is WGS84 or functionally identical (no Helmert shift needed).
14    pub fn is_wgs84_compatible(&self) -> bool {
15        matches!(self.to_wgs84, DatumToWgs84::Identity)
16    }
17
18    /// Returns true if this datum has a known path to WGS84.
19    pub fn has_known_wgs84_transform(&self) -> bool {
20        !matches!(self.to_wgs84, DatumToWgs84::Unknown)
21    }
22
23    /// Return the Helmert parameters for this datum's path to WGS84, when available.
24    pub fn helmert_to_wgs84(&self) -> Option<&HelmertParams> {
25        match &self.to_wgs84 {
26            DatumToWgs84::Helmert(params) => Some(params),
27            DatumToWgs84::Identity | DatumToWgs84::Unknown => None,
28        }
29    }
30
31    pub fn approximate_helmert_to(&self, target: &Datum) -> Option<HelmertParams> {
32        if !self.has_known_wgs84_transform() || !target.has_known_wgs84_transform() {
33            return None;
34        }
35
36        let source_to_wgs84 = self
37            .helmert_to_wgs84()
38            .copied()
39            .unwrap_or_else(|| HelmertParams::translation(0.0, 0.0, 0.0));
40        let wgs84_to_target = target
41            .helmert_to_wgs84()
42            .map(HelmertParams::inverse)
43            .unwrap_or_else(|| HelmertParams::translation(0.0, 0.0, 0.0));
44
45        Some(source_to_wgs84.compose_approx(&wgs84_to_target))
46    }
47
48    /// Returns true if two datums are the same (same ellipsoid, same Helmert parameters).
49    pub fn same_datum(&self, other: &Datum) -> bool {
50        let same_ellipsoid = (self.ellipsoid.a - other.ellipsoid.a).abs() < 1e-6
51            && (self.ellipsoid.f - other.ellipsoid.f).abs() < 1e-12;
52
53        match (&self.to_wgs84, &other.to_wgs84) {
54            (DatumToWgs84::Identity, DatumToWgs84::Identity) => same_ellipsoid,
55            (DatumToWgs84::Helmert(a), DatumToWgs84::Helmert(b)) => {
56                same_ellipsoid && a.approx_eq(b)
57            }
58            (DatumToWgs84::Unknown, DatumToWgs84::Unknown) => false,
59            _ => false,
60        }
61    }
62}
63
64/// Explicit WGS84 relationship for a datum.
65#[derive(Debug, Clone, Copy)]
66pub enum DatumToWgs84 {
67    /// The datum can be treated as WGS84-compatible in the current model.
68    Identity,
69    /// The datum requires the provided Helmert transform to reach WGS84.
70    Helmert(HelmertParams),
71    /// The datum's path to WGS84 is not known.
72    Unknown,
73}
74
75/// 7-parameter Helmert (Bursa-Wolf) transformation parameters.
76///
77/// Defines the transformation from one datum to WGS84 geocentric coordinates:
78/// ```text
79/// [X']   [dx]         [  1  -rz   ry] [X]
80/// [Y'] = [dy] + (1+ds)[  rz   1  -rx] [Y]
81/// [Z']   [dz]         [ -ry  rx   1 ] [Z]
82/// ```
83#[derive(Debug, Clone, Copy, PartialEq)]
84pub struct HelmertParams {
85    /// X-axis translation in meters.
86    pub dx: f64,
87    /// Y-axis translation in meters.
88    pub dy: f64,
89    /// Z-axis translation in meters.
90    pub dz: f64,
91    /// X-axis rotation in arc-seconds.
92    pub rx: f64,
93    /// Y-axis rotation in arc-seconds.
94    pub ry: f64,
95    /// Z-axis rotation in arc-seconds.
96    pub rz: f64,
97    /// Scale difference in parts-per-million (ppm).
98    pub ds: f64,
99}
100
101impl HelmertParams {
102    /// Create a translation-only (3-parameter) transformation.
103    pub const fn translation(dx: f64, dy: f64, dz: f64) -> Self {
104        Self {
105            dx,
106            dy,
107            dz,
108            rx: 0.0,
109            ry: 0.0,
110            rz: 0.0,
111            ds: 0.0,
112        }
113    }
114
115    /// Return the inverse parameters (WGS84 → this datum).
116    pub fn inverse(&self) -> Self {
117        Self {
118            dx: -self.dx,
119            dy: -self.dy,
120            dz: -self.dz,
121            rx: -self.rx,
122            ry: -self.ry,
123            rz: -self.rz,
124            ds: -self.ds,
125        }
126    }
127
128    pub fn compose_approx(&self, next: &Self) -> Self {
129        Self {
130            dx: self.dx + next.dx,
131            dy: self.dy + next.dy,
132            dz: self.dz + next.dz,
133            rx: self.rx + next.rx,
134            ry: self.ry + next.ry,
135            rz: self.rz + next.rz,
136            ds: self.ds + next.ds,
137        }
138    }
139
140    fn approx_eq(&self, other: &Self) -> bool {
141        (self.dx - other.dx).abs() < 1e-6
142            && (self.dy - other.dy).abs() < 1e-6
143            && (self.dz - other.dz).abs() < 1e-6
144            && (self.rx - other.rx).abs() < 1e-9
145            && (self.ry - other.ry).abs() < 1e-9
146            && (self.rz - other.rz).abs() < 1e-9
147            && (self.ds - other.ds).abs() < 1e-9
148    }
149}
150
151// ---------------------------------------------------------------------------
152// Well-known datums
153// ---------------------------------------------------------------------------
154
155/// WGS 84 datum.
156pub const WGS84: Datum = Datum {
157    ellipsoid: ellipsoid::WGS84,
158    to_wgs84: DatumToWgs84::Identity,
159};
160
161/// NAD83 datum (functionally identical to WGS84 for sub-meter work).
162pub const NAD83: Datum = Datum {
163    ellipsoid: ellipsoid::GRS80,
164    to_wgs84: DatumToWgs84::Identity,
165};
166
167/// NAD27 datum (Clarke 1866 ellipsoid).
168/// Helmert parameters from EPSG dataset (approximate continental US average).
169pub const NAD27: Datum = Datum {
170    ellipsoid: ellipsoid::CLARKE1866,
171    to_wgs84: DatumToWgs84::Helmert(HelmertParams::translation(-8.0, 160.0, 176.0)),
172};
173
174/// ETRS89 datum (European Terrestrial Reference System 1989).
175/// Functionally identical to WGS84 for most purposes.
176pub const ETRS89: Datum = Datum {
177    ellipsoid: ellipsoid::GRS80,
178    to_wgs84: DatumToWgs84::Identity,
179};
180
181/// OSGB36 datum (Ordnance Survey Great Britain 1936).
182pub const OSGB36: Datum = Datum {
183    ellipsoid: ellipsoid::AIRY1830,
184    to_wgs84: DatumToWgs84::Helmert(HelmertParams {
185        dx: 446.448,
186        dy: -125.157,
187        dz: 542.060,
188        rx: 0.1502,
189        ry: 0.2470,
190        rz: 0.8421,
191        ds: -20.4894,
192    }),
193};
194
195/// Pulkovo 1942 datum (used in Russia and former Soviet states).
196pub const PULKOVO1942: Datum = Datum {
197    ellipsoid: ellipsoid::KRASSOWSKY,
198    to_wgs84: DatumToWgs84::Helmert(HelmertParams::translation(23.92, -141.27, -80.9)),
199};
200
201/// ED50 datum (European Datum 1950).
202pub const ED50: Datum = Datum {
203    ellipsoid: ellipsoid::INTL1924,
204    to_wgs84: DatumToWgs84::Helmert(HelmertParams::translation(-87.0, -98.0, -121.0)),
205};
206
207/// Tokyo datum (used in Japan).
208pub const TOKYO: Datum = Datum {
209    ellipsoid: ellipsoid::BESSEL1841,
210    to_wgs84: DatumToWgs84::Helmert(HelmertParams::translation(-146.414, 507.337, 680.507)),
211};
212
213#[cfg(test)]
214mod tests {
215    use super::*;
216
217    #[test]
218    fn wgs84_is_wgs84_compatible() {
219        assert!(WGS84.is_wgs84_compatible());
220        assert!(NAD83.is_wgs84_compatible());
221        assert!(ETRS89.is_wgs84_compatible());
222    }
223
224    #[test]
225    fn nad27_is_not_wgs84_compatible() {
226        assert!(!NAD27.is_wgs84_compatible());
227        assert!(!OSGB36.is_wgs84_compatible());
228    }
229
230    #[test]
231    fn same_datum_identity() {
232        assert!(WGS84.same_datum(&WGS84));
233        assert!(NAD27.same_datum(&NAD27));
234    }
235
236    #[test]
237    fn different_datums() {
238        assert!(!WGS84.same_datum(&NAD27));
239        assert!(!NAD27.same_datum(&OSGB36));
240    }
241
242    #[test]
243    fn unknown_datums_are_not_collapsed_by_ellipsoid() {
244        let a = Datum {
245            ellipsoid: ellipsoid::WGS84,
246            to_wgs84: DatumToWgs84::Unknown,
247        };
248        let b = Datum {
249            ellipsoid: ellipsoid::WGS84,
250            to_wgs84: DatumToWgs84::Unknown,
251        };
252
253        assert!(!a.same_datum(&b));
254    }
255
256    #[test]
257    fn helmert_inverse_negates() {
258        let h = HelmertParams {
259            dx: 1.0,
260            dy: 2.0,
261            dz: 3.0,
262            rx: 0.1,
263            ry: 0.2,
264            rz: 0.3,
265            ds: 0.5,
266        };
267        let inv = h.inverse();
268        assert_eq!(inv.dx, -1.0);
269        assert_eq!(inv.rx, -0.1);
270        assert_eq!(inv.ds, -0.5);
271    }
272}