rustial_engine/math/
equirectangular.rs

1//! Equirectangular (Plate Carree) projection (EPSG:4326 scaled to meters).
2//!
3//! # Overview
4//!
5//! The Equirectangular projection is the simplest cylindrical projection:
6//! longitude and latitude are mapped linearly to X and Y, then scaled by
7//! the WGS-84 semi-major axis so that output units are **meters** at the
8//! equator.
9//!
10//! ```text
11//! x = R * lon_rad
12//! y = R * lat_rad
13//! ```
14//!
15//! # Distortion
16//!
17//! - **Equator:** no distortion (scale factor = 1.0).
18//! - **Away from equator:** X distances are progressively stretched by
19//!   `sec(lat)`.  Y distances remain undistorted at all latitudes.
20//! - **Poles:** the projection is valid up to +/-90 degrees latitude,
21//!   but X distortion is infinite at the poles.
22//!
23//! # When to use
24//!
25//! Equirectangular is useful as a simple baseline, for data interchange
26//! (many raster datasets are stored in EPSG:4326 pixel grids), and for
27//! server-side tile pre-processing where rendering fidelity is not
28//! critical.  For map display, prefer [`WebMercator`](crate::WebMercator).
29//!
30//! # Altitude
31//!
32//! Altitude (`alt` / `z`) is passed through unchanged in both directions.
33
34use crate::coord::{GeoCoord, WorldCoord};
35use crate::ellipsoid::Ellipsoid;
36use crate::projection::Projection;
37
38/// Equirectangular projection: trivial linear mapping of lon/lat to meters.
39///
40/// Uses the WGS-84 semi-major axis (`R = 6 378 137 m`) as the scaling
41/// factor so that output units are meters at the equator.  This is a
42/// zero-size type -- construction is free.
43pub struct Equirectangular;
44
45impl Projection for Equirectangular {
46    /// Forward projection: geographic degrees to world meters.
47    ///
48    /// `x = R * lon_rad`, `y = R * lat_rad`, `z = alt` (passthrough).
49    fn project(&self, geo: &GeoCoord) -> WorldCoord {
50        let a = Ellipsoid::WGS84.a;
51        WorldCoord::new(a * geo.lon.to_radians(), a * geo.lat.to_radians(), geo.alt)
52    }
53
54    /// Inverse projection: world meters back to geographic degrees.
55    ///
56    /// `lat = degrees(y / R)`, `lon = degrees(x / R)`, `alt = z`.
57    fn unproject(&self, world: &WorldCoord) -> GeoCoord {
58        let a = Ellipsoid::WGS84.a;
59        let lat = (world.position.y / a).to_degrees().clamp(-90.0, 90.0);
60        let lon = ((world.position.x / a).to_degrees() + 180.0).rem_euclid(360.0) - 180.0;
61        GeoCoord::new(lat, lon, world.position.z)
62    }
63
64    /// Scale factor for Equirectangular: `sec(lat) = 1 / cos(lat)`.
65    ///
66    /// The projection maps longitude linearly, so X distances are
67    /// stretched relative to true distances by `sec(lat)`, the same
68    /// as the meridian convergence factor.  Along Y the scale is 1.0,
69    /// so this returns the *maximum* linear distortion at the point.
70    fn scale_factor(&self, geo: &GeoCoord) -> f64 {
71        1.0 / geo.lat.to_radians().cos()
72    }
73
74    // `projection_bounds` is not overridden: the default full-globe
75    // range (+/-90 lat, +/-180 lon) is correct for Equirectangular.
76}
77
78#[cfg(test)]
79mod tests {
80    use super::*;
81
82    // -- Roundtrip ---------------------------------------------------------
83
84    #[test]
85    fn roundtrip_origin() {
86        let geo = GeoCoord::from_lat_lon(0.0, 0.0);
87        let world = Equirectangular.project(&geo);
88        let back = Equirectangular.unproject(&world);
89        assert!((back.lat - geo.lat).abs() < 1e-10);
90        assert!((back.lon - geo.lon).abs() < 1e-10);
91    }
92
93    #[test]
94    fn roundtrip_nonzero() {
95        let geo = GeoCoord::from_lat_lon(45.0, 90.0);
96        let world = Equirectangular.project(&geo);
97        let back = Equirectangular.unproject(&world);
98        assert!((back.lat - geo.lat).abs() < 1e-8);
99        assert!((back.lon - geo.lon).abs() < 1e-8);
100    }
101
102    #[test]
103    fn poles_roundtrip() {
104        let north = GeoCoord::from_lat_lon(90.0, 0.0);
105        let north_w = Equirectangular.project(&north);
106        let north_back = Equirectangular.unproject(&north_w);
107        assert!((north_back.lat - 90.0).abs() < 1e-9);
108
109        let south = GeoCoord::from_lat_lon(-90.0, 0.0);
110        let south_w = Equirectangular.project(&south);
111        let south_back = Equirectangular.unproject(&south_w);
112        assert!((south_back.lat + 90.0).abs() < 1e-9);
113    }
114
115    #[test]
116    fn anti_meridian_roundtrip() {
117        // +180 and -180 represent the same meridian; the normalisation
118        // in `unproject` canonicalises both to the (-180, 180] range.
119        let east = GeoCoord::from_lat_lon(10.0, 180.0);
120        let east_back = Equirectangular.unproject(&Equirectangular.project(&east));
121        assert!((east_back.lon.abs() - 180.0).abs() < 1e-9);
122
123        let west = GeoCoord::from_lat_lon(10.0, -180.0);
124        let west_back = Equirectangular.unproject(&Equirectangular.project(&west));
125        assert!((west_back.lon.abs() - 180.0).abs() < 1e-9);
126    }
127
128    #[test]
129    fn extreme_latitudes_stable() {
130        for lat in [89.999, -89.999, 89.9, -89.9] {
131            let geo = GeoCoord::from_lat_lon(lat, 45.0);
132            let world = Equirectangular.project(&geo);
133            let back = Equirectangular.unproject(&world);
134            assert!(
135                (back.lat - lat).abs() < 1e-6,
136                "lat roundtrip failed for {lat}"
137            );
138        }
139    }
140
141    // -- Metric output at the equator -------------------------------------
142
143    #[test]
144    fn equator_scale() {
145        let a = GeoCoord::from_lat_lon(0.0, 0.0);
146        let b = GeoCoord::from_lat_lon(0.0, 1.0);
147        let wa = Equirectangular.project(&a);
148        let wb = Equirectangular.project(&b);
149        let dx = wb.position.x - wa.position.x;
150        // 1 degree of longitude at the equator ~ 111,319.49 m.
151        assert!((dx - 111_319.49).abs() < 1.0);
152    }
153
154    // -- Scale factor ------------------------------------------------------
155
156    #[test]
157    fn scale_factor_equator() {
158        let sf = Equirectangular.scale_factor(&GeoCoord::from_lat_lon(0.0, 0.0));
159        assert!((sf - 1.0).abs() < 1e-10);
160    }
161
162    #[test]
163    fn scale_factor_45_degrees() {
164        // sec(45 deg) = sqrt(2)
165        let sf = Equirectangular.scale_factor(&GeoCoord::from_lat_lon(45.0, 0.0));
166        assert!((sf - std::f64::consts::SQRT_2).abs() < 1e-10);
167    }
168
169    // -- Altitude passthrough ----------------------------------------------
170
171    #[test]
172    fn altitude_passthrough() {
173        let geo = GeoCoord::new(30.0, 60.0, 1234.5);
174        let world = Equirectangular.project(&geo);
175        let back = Equirectangular.unproject(&world);
176        assert!((world.position.z - 1234.5).abs() < 1e-12);
177        assert!((back.alt - 1234.5).abs() < 1e-12);
178    }
179
180    // -- Projection bounds default to full globe --------------------------
181
182    #[test]
183    fn projection_bounds_full_globe() {
184        let bounds = Equirectangular.projection_bounds();
185        assert!((bounds.sw().lat - (-90.0)).abs() < 1e-10);
186        assert!((bounds.ne().lat - 90.0).abs() < 1e-10);
187        assert!((bounds.sw().lon - (-180.0)).abs() < 1e-10);
188        assert!((bounds.ne().lon - 180.0).abs() < 1e-10);
189    }
190}
rustial_engine/math/equirectangular.rs

rustial_engine/math/
equirectangular.rs
