Skip to main content

proj_core/
datum.rs

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