Skip to main content

grib/
projection.rs

1//! Map projection functionality.
2
3pub use lcc::{Params as LccParams, Projection as Lcc};
4pub use merc::{Params as MercParams, Projection as Merc};
5
6/// Map projection functionality.
7pub trait Project {
8    fn forward(&self, xy: &(f64, f64)) -> Result<(f64, f64), &'static str>;
9    fn inverse(&self, xy: &(f64, f64)) -> Result<(f64, f64), &'static str>;
10    fn a(&self) -> &f64;
11    fn lam0(&self) -> &f64;
12
13    /// Performs a projection.
14    ///
15    /// For forward transformation (`inverse = false`), `xy` is treated as
16    /// `(lambda, phi)`, where `lambda` represents longitude and `phi`
17    /// represents latitude, both expressed in degrees. The return value
18    /// should be considered to be `(x, y)`, and both `x` and `y` are values in
19    /// meters.
20    ///
21    /// For inverse transformation (`inverse = true`), `xy` is treated as `(x,
22    /// y)` and both `x` and `y` are values in meters.
23    /// The return value should be considered to be `(lambda, phi)`, where
24    /// `lambda` represents longitude and `phi` represents latitude, both
25    /// expressed in degrees.
26    fn project(&self, xy: &(f64, f64), inverse: bool) -> Result<(f64, f64), &'static str> {
27        if inverse {
28            let &(x, y) = xy;
29            let x = x / self.a();
30            let y = y / self.a();
31            let (lambda, phi) = self.inverse(&(x, y))?;
32            let lambda = helpers::normalize_longitude(lambda + self.lam0());
33            Ok((lambda, phi))
34        } else {
35            let &(lambda, phi) = xy;
36            let lambda = helpers::normalize_longitude(lambda - self.lam0());
37            let (x, y) = self.forward(&(lambda, phi))?;
38            let x = x * self.a();
39            let y = y * self.a();
40            Ok((x, y))
41        }
42    }
43}
44
45#[cfg(feature = "gridpoints-proj")]
46pub(crate) trait OsgeoProj {
47    fn proj_args(&self) -> String;
48}
49
50/// Parameters for Stereographic projection.
51#[derive(Debug, PartialEq, Clone)]
52pub struct StereParams {
53    /// Ellipsoid definition.
54    pub ellipsoid: Ellipsoid,
55    /// Latitude where scale is not distorted (in degrees).
56    pub lat_ts: f64,
57    /// Latitude of origin (in degrees).
58    pub lat_0: f64,
59    /// Central meridian (in degrees).
60    pub lon_0: f64,
61}
62
63#[cfg(feature = "gridpoints-proj")]
64impl OsgeoProj for StereParams {
65    fn proj_args(&self) -> String {
66        let Self {
67            ellipsoid: Ellipsoid { a, b, .. },
68            lat_ts,
69            lat_0,
70            lon_0,
71        } = self;
72        format!("+a={a} +b={b} +proj=stere +lat_ts={lat_ts} +lat_0={lat_0} +lon_0={lon_0}")
73    }
74}
75
76/// Ellipsoid definition.
77#[derive(Debug, PartialEq, Clone)]
78pub struct Ellipsoid {
79    /// Semimajor radius of the ellipsoid axis (in meters).
80    pub a: f64,
81    /// Semiminor radius of the ellipsoid axis (in meters).
82    pub b: f64,
83    /// Eccentricity.
84    pub e: f64,
85    /// Eccentricity squared.
86    pub e_sq: f64,
87}
88
89impl Ellipsoid {
90    pub fn from_a_and_b(a: f64, b: f64) -> Self {
91        let f = (a - b) / a;
92        let e_sq = 2. * f - f * f;
93        let e = e_sq.sqrt();
94        Self { a, b, e, e_sq }
95    }
96}
97
98mod helpers;
99mod lcc;
100mod merc;