use crate::ellipsoid::{self, Ellipsoid};
use crate::error::{Error, Result};
use crate::grid::GridDefinition;
use smallvec::SmallVec;
#[derive(Debug, Clone)]
pub struct Datum {
ellipsoid: Ellipsoid,
to_wgs84: DatumToWgs84,
}
impl Datum {
pub fn new(ellipsoid: Ellipsoid, to_wgs84: DatumToWgs84) -> Result<Self> {
to_wgs84.validate()?;
Ok(Self {
ellipsoid,
to_wgs84,
})
}
const fn new_unchecked(ellipsoid: Ellipsoid, to_wgs84: DatumToWgs84) -> Self {
Self {
ellipsoid,
to_wgs84,
}
}
pub const fn ellipsoid(&self) -> Ellipsoid {
self.ellipsoid
}
pub const fn to_wgs84(&self) -> &DatumToWgs84 {
&self.to_wgs84
}
pub fn is_wgs84_compatible(&self) -> bool {
matches!(self.to_wgs84, DatumToWgs84::Identity)
}
pub fn has_known_wgs84_transform(&self) -> bool {
!matches!(self.to_wgs84, DatumToWgs84::Unknown)
}
pub fn uses_grid_shift(&self) -> bool {
self.to_wgs84.uses_grid_shift()
}
pub fn helmert_to_wgs84(&self) -> Option<&HelmertParams> {
match &self.to_wgs84 {
DatumToWgs84::Helmert(params) => Some(params),
DatumToWgs84::Identity | DatumToWgs84::GridShift(_) | DatumToWgs84::Unknown => None,
}
}
pub fn approximate_helmert_to(&self, target: &Datum) -> Option<HelmertParams> {
if self.uses_grid_shift() || target.uses_grid_shift() {
return None;
}
if !self.has_known_wgs84_transform() || !target.has_known_wgs84_transform() {
return None;
}
let source_to_wgs84 = self
.helmert_to_wgs84()
.copied()
.unwrap_or_else(HelmertParams::identity);
let wgs84_to_target = target
.helmert_to_wgs84()
.map(HelmertParams::inverse)
.unwrap_or_else(HelmertParams::identity);
source_to_wgs84.compose_approx(&wgs84_to_target).ok()
}
pub fn same_datum(&self, other: &Datum) -> bool {
let same_ellipsoid =
(self.ellipsoid.semi_major_axis() - other.ellipsoid.semi_major_axis()).abs() < 1e-6
&& (self.ellipsoid.flattening() - other.ellipsoid.flattening()).abs() < 1e-12;
match (&self.to_wgs84, &other.to_wgs84) {
(DatumToWgs84::Identity, DatumToWgs84::Identity) => same_ellipsoid,
(DatumToWgs84::Helmert(a), DatumToWgs84::Helmert(b)) => {
same_ellipsoid && a.approx_eq(b)
}
(DatumToWgs84::GridShift(a), DatumToWgs84::GridShift(b)) => same_ellipsoid && a == b,
(DatumToWgs84::Unknown, DatumToWgs84::Unknown) => false,
_ => false,
}
}
}
#[derive(Debug, Clone, PartialEq)]
pub enum DatumToWgs84 {
Identity,
Helmert(HelmertParams),
GridShift(Box<DatumGridShift>),
Unknown,
}
impl DatumToWgs84 {
pub fn uses_grid_shift(&self) -> bool {
matches!(self, DatumToWgs84::GridShift(shift) if shift.uses_grid_shift())
}
pub fn validate(&self) -> Result<()> {
match self {
Self::Helmert(params) => params.validate(),
Self::Identity | Self::GridShift(_) | Self::Unknown => Ok(()),
}
}
}
#[derive(Debug, Clone, PartialEq)]
pub struct DatumGridShift {
entries: SmallVec<[DatumGridShiftEntry; 4]>,
}
impl DatumGridShift {
pub fn new(entries: SmallVec<[DatumGridShiftEntry; 4]>) -> Self {
Self { entries }
}
pub fn from_vec(entries: Vec<DatumGridShiftEntry>) -> Self {
Self {
entries: SmallVec::from_vec(entries),
}
}
pub fn entries(&self) -> &[DatumGridShiftEntry] {
&self.entries
}
pub fn uses_grid_shift(&self) -> bool {
self.entries
.iter()
.any(|entry| matches!(entry, DatumGridShiftEntry::Grid { .. }))
}
}
#[derive(Debug, Clone, PartialEq)]
pub enum DatumGridShiftEntry {
Grid {
definition: GridDefinition,
optional: bool,
},
Null,
}
#[derive(Debug, Clone, Copy, PartialEq)]
pub struct HelmertParams {
dx: f64,
dy: f64,
dz: f64,
rx: f64,
ry: f64,
rz: f64,
ds: f64,
}
impl HelmertParams {
pub fn new(dx: f64, dy: f64, dz: f64, rx: f64, ry: f64, rz: f64, ds: f64) -> Result<Self> {
let params = Self::new_unchecked(dx, dy, dz, rx, ry, rz, ds);
params.validate()?;
Ok(params)
}
pub fn translation(dx: f64, dy: f64, dz: f64) -> Result<Self> {
Self::new(dx, dy, dz, 0.0, 0.0, 0.0, 0.0)
}
pub const fn identity() -> Self {
Self::new_unchecked(0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0)
}
const fn translation_unchecked(dx: f64, dy: f64, dz: f64) -> Self {
Self::new_unchecked(dx, dy, dz, 0.0, 0.0, 0.0, 0.0)
}
const fn new_unchecked(dx: f64, dy: f64, dz: f64, rx: f64, ry: f64, rz: f64, ds: f64) -> Self {
Self {
dx,
dy,
dz,
rx,
ry,
rz,
ds,
}
}
pub const fn dx(&self) -> f64 {
self.dx
}
pub const fn dy(&self) -> f64 {
self.dy
}
pub const fn dz(&self) -> f64 {
self.dz
}
pub const fn rx(&self) -> f64 {
self.rx
}
pub const fn ry(&self) -> f64 {
self.ry
}
pub const fn rz(&self) -> f64 {
self.rz
}
pub const fn ds(&self) -> f64 {
self.ds
}
pub fn validate(&self) -> Result<()> {
if self.dx.is_finite()
&& self.dy.is_finite()
&& self.dz.is_finite()
&& self.rx.is_finite()
&& self.ry.is_finite()
&& self.rz.is_finite()
&& self.ds.is_finite()
{
return Ok(());
}
Err(Error::InvalidDefinition(
"Helmert parameters must be finite".into(),
))
}
pub fn inverse(&self) -> Self {
Self {
dx: -self.dx,
dy: -self.dy,
dz: -self.dz,
rx: -self.rx,
ry: -self.ry,
rz: -self.rz,
ds: -self.ds,
}
}
pub fn compose_approx(&self, next: &Self) -> Result<Self> {
Self::new(
self.dx + next.dx,
self.dy + next.dy,
self.dz + next.dz,
self.rx + next.rx,
self.ry + next.ry,
self.rz + next.rz,
self.ds + next.ds,
)
}
fn approx_eq(&self, other: &Self) -> bool {
(self.dx - other.dx).abs() < 1e-6
&& (self.dy - other.dy).abs() < 1e-6
&& (self.dz - other.dz).abs() < 1e-6
&& (self.rx - other.rx).abs() < 1e-9
&& (self.ry - other.ry).abs() < 1e-9
&& (self.rz - other.rz).abs() < 1e-9
&& (self.ds - other.ds).abs() < 1e-9
}
}
pub const WGS84: Datum = Datum::new_unchecked(ellipsoid::WGS84, DatumToWgs84::Identity);
pub const NAD83: Datum = Datum::new_unchecked(ellipsoid::GRS80, DatumToWgs84::Identity);
pub const NAD27: Datum = Datum::new_unchecked(
ellipsoid::CLARKE1866,
DatumToWgs84::Helmert(HelmertParams::translation_unchecked(-8.0, 160.0, 176.0)),
);
pub const ETRS89: Datum = Datum::new_unchecked(ellipsoid::GRS80, DatumToWgs84::Identity);
pub const OSGB36: Datum = Datum::new_unchecked(
ellipsoid::AIRY1830,
DatumToWgs84::Helmert(HelmertParams::new_unchecked(
446.448, -125.157, 542.060, 0.1502, 0.2470, 0.8421, -20.4894,
)),
);
pub const PULKOVO1942: Datum = Datum::new_unchecked(
ellipsoid::KRASSOWSKY,
DatumToWgs84::Helmert(HelmertParams::translation_unchecked(23.92, -141.27, -80.9)),
);
pub const ED50: Datum = Datum::new_unchecked(
ellipsoid::INTL1924,
DatumToWgs84::Helmert(HelmertParams::translation_unchecked(-87.0, -98.0, -121.0)),
);
pub const TOKYO: Datum = Datum::new_unchecked(
ellipsoid::BESSEL1841,
DatumToWgs84::Helmert(HelmertParams::translation_unchecked(
-146.414, 507.337, 680.507,
)),
);
#[cfg(test)]
mod tests {
use super::*;
#[test]
fn wgs84_is_wgs84_compatible() {
assert!(WGS84.is_wgs84_compatible());
assert!(NAD83.is_wgs84_compatible());
assert!(ETRS89.is_wgs84_compatible());
}
#[test]
fn nad27_is_not_wgs84_compatible() {
assert!(!NAD27.is_wgs84_compatible());
assert!(!OSGB36.is_wgs84_compatible());
}
#[test]
fn same_datum_identity() {
assert!(WGS84.same_datum(&WGS84));
assert!(NAD27.same_datum(&NAD27));
}
#[test]
fn different_datums() {
assert!(!WGS84.same_datum(&NAD27));
assert!(!NAD27.same_datum(&OSGB36));
}
#[test]
fn unknown_datums_are_not_collapsed_by_ellipsoid() {
let a = Datum::new(ellipsoid::WGS84, DatumToWgs84::Unknown).unwrap();
let b = Datum::new(ellipsoid::WGS84, DatumToWgs84::Unknown).unwrap();
assert!(!a.same_datum(&b));
}
#[test]
fn helmert_inverse_negates() {
let h = HelmertParams::new(1.0, 2.0, 3.0, 0.1, 0.2, 0.3, 0.5).unwrap();
let inv = h.inverse();
assert_eq!(inv.dx(), -1.0);
assert_eq!(inv.rx(), -0.1);
assert_eq!(inv.ds(), -0.5);
}
#[test]
fn helmert_params_reject_non_finite_values() {
let err = HelmertParams::new(f64::NAN, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0).unwrap_err();
assert!(matches!(err, Error::InvalidDefinition(_)), "got {err}");
let err = HelmertParams::new(0.0, 0.0, 0.0, 0.0, 0.0, f64::INFINITY, 0.0).unwrap_err();
assert!(err.to_string().contains("Helmert parameters"), "{err}");
}
}