proj_core/
datum.rs

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

proj_core/
datum.rs
