toolbox_rs/
projection.rs

1pub mod mercator {
2    use crate::wgs84::EARTH_RADIUS_KM;
3
4    pub fn lon2x(lon: f64) -> f64 {
5        EARTH_RADIUS_KM * 1000. * lon.to_radians()
6    }
7
8    pub fn x2lon(x: f64) -> f64 {
9        (x / (EARTH_RADIUS_KM * 1000.)).to_degrees()
10    }
11
12    pub fn lat2y(lat: f64) -> f64 {
13        ((lat.to_radians() / 2. + std::f64::consts::PI / 4.).tan()).log(std::f64::consts::E)
14            * EARTH_RADIUS_KM
15            * 1000.
16    }
17
18    pub fn y2lat(y: f64) -> f64 {
19        (2. * ((y / (EARTH_RADIUS_KM * 1000.)).exp()).atan() - std::f64::consts::PI / 2.)
20            .to_degrees()
21    }
22}
23
24#[cfg(test)]
25mod tests {
26    use crate::projection::mercator::{lat2y, lon2x, x2lon, y2lat};
27
28    const ALLOWED_ERROR: f64 = 0.0000000000001;
29    // Allowed error in IEEE-754-based projection math.
30    // Note that this is way below a centimeter of error
31
32    #[test]
33    fn lon_conversion_roundtrip() {
34        // Roundtrip calculation of the projection with expected tiny errors
35
36        // longitude in [180. to -180.]
37        for i in -18_000..18_001 {
38            // off-by-one to be inclusive of 180.
39            let lon = f64::from(i) * 0.01;
40            let result = x2lon(lon2x(lon));
41            assert!((lon - result).abs() < ALLOWED_ERROR);
42        }
43    }
44
45    #[test]
46    fn lat_conversion_roundtrip() {
47        // Roundtrip calculation of the projection with expected tiny errors
48
49        // latitude in [90. to -90.]
50        for i in -90_00..90_01 {
51            // off-by-one to be inclusive of 90.
52            let lat = f64::from(i) * 0.01;
53            let result = y2lat(lat2y(lat));
54            assert!((lat - result).abs() < ALLOWED_ERROR);
55        }
56    }
57}
toolbox_rs/projection.rs

toolbox_rs/
projection.rs
