gistools/proj/project/
eqc.rs

1use crate::proj::{
2    CoordinateStep, EQUIDISTANT_CYLINDRICAL, LATITUDE_STD_PARALLEL, Proj, ProjValue,
3    ProjectCoordinates, TransformCoordinates,
4};
5use alloc::rc::Rc;
6use core::cell::RefCell;
7use libm::cos;
8
9/// Equidistant Cylindrical variables
10#[derive(Debug, Default, Clone, Copy, PartialEq)]
11pub struct Eqc {
12    rc: f64,
13}
14
15/// # Equidistant Cylindrical (Plate Carrée)
16///
17/// **Classification**: Conformal cylindrical
18///
19/// **Available forms**: Forward and inverse
20///
21/// **Defined area**: Global, but best used near the equator
22///
23/// **Alias**: eqc, plate_carrée, simple_cylindrical
24///
25/// **Domain**: 2D
26///
27/// **Input type**: Geodetic coordinates
28///
29/// **Output type**: Projected coordinates
30///
31/// ## Projection String
32/// ```ini
33/// +proj=eqc
34/// ```
35///
36/// ## Usage
37///
38/// Because of the distortions introduced by this projection, it has little use in navigation or
39/// cadastral mapping and finds its main use in thematic mapping.
40/// In particular, the Plate Carrée has become a standard for global raster datasets, such as
41/// Celestia and NASA World Wind, because of the particularly simple relationship between the
42/// position of an image pixel on the map and its corresponding geographic location on Earth.
43///
44/// ### Special Cases of Cylindrical Equidistant Projection:
45///
46/// - Plain/Plane Chart: $0°$
47/// - Simple Cylindrical: $0°$
48/// - Plate Carrée: $0°$
49/// - Ronald Miller—minimum overall scale distortion: $37°30'$
50/// - E. Grafarend and A. Niermann: $42°$
51/// - Ronald Miller—minimum continental scale distortion: $43°30'$
52/// - Gall Isographic: $45°$
53/// - Ronald Miller Equirectangular: $50°30'$
54/// - E. Grafarend and A. Niermann minimum linear distortion: $61°7'$
55///
56/// ## Example
57///
58/// Example using EPSG 32662 (WGS84 Plate Carrée):
59/// ```bash
60/// echo 2 47 | proj +proj=eqc +ellps=WGS84
61/// ```
62/// Output: 222638.98 5232016.07
63///
64/// Example using Plate Carrée projection with true scale at latitude 30° and central meridian 90°W:
65/// ```bash
66/// echo -88 30 | proj +proj=eqc +lat_ts=30 +lon_0=90w
67/// ```
68/// Output: 192811.01 3339584.72
69///
70/// ## Parameters
71///
72/// - `+lon_0` (Central meridian)
73/// - `+lat_0` (Latitude of origin)
74/// - `+lat_ts` (Latitude of true scale)
75/// - `+x_0` (False easting)
76/// - `+y_0` (False northing)
77/// - `+ellps` (Ellipsoid name)
78/// - `+R` (Radius of the sphere)
79///
80/// ## Mathematical Definition
81///
82/// ### Forward projection:
83/// $$x = \lambda \cos(\phi_{ts})$$
84/// $$y = \phi - \phi_0$$
85///
86/// ### Inverse projection:
87/// $$\lambda = x / \cos(\phi_{ts})$$
88/// $$\phi = y + \phi_0$$
89///
90/// ## Further Reading
91///
92/// - [Wikipedia](https://en.wikipedia.org/wiki/Equirectangular_projection)
93/// - [Wolfram Mathworld](http://mathworld.wolfram.com/CylindricalEquidistantProjection.html)
94///
95/// ![Equidistant Cylindrical](https://github.com/Open-S2/gis-tools/blob/master/assets/proj4/projections/images/eqc.png?raw=true)
96#[derive(Debug, Clone, PartialEq)]
97pub struct EquidistantCylindricalProjection {
98    proj: Rc<RefCell<Proj>>,
99    store: RefCell<Eqc>,
100}
101impl ProjectCoordinates for EquidistantCylindricalProjection {
102    fn code(&self) -> i64 {
103        EQUIDISTANT_CYLINDRICAL
104    }
105    fn name(&self) -> &'static str {
106        "Equidistant Cylindrical"
107    }
108    fn names() -> &'static [&'static str] {
109        &[
110            "Equidistant Cylindrical",
111            "EquidistantCylindrical",
112            "Equidistant Cylindrical (Plate Carree)",
113            "eqc",
114        ]
115    }
116}
117impl CoordinateStep for EquidistantCylindricalProjection {
118    fn new(proj: Rc<RefCell<Proj>>) -> Self {
119        let mut store = Eqc::default();
120        {
121            let proj = &mut proj.borrow_mut();
122            let lat_ts = proj
123                .params
124                .get(&LATITUDE_STD_PARALLEL) // (lat_ts)
125                .unwrap_or(&ProjValue::default())
126                .f64()
127                .to_radians();
128            if cos(lat_ts) <= 0. {
129                panic!("Invalid value for lat_ts: |lat_ts| should be <= 90°");
130            }
131            store.rc = lat_ts;
132            proj.es = 0.;
133        }
134        EquidistantCylindricalProjection { proj, store: store.into() }
135    }
136    fn forward<P: TransformCoordinates>(&self, p: &mut P) {
137        eqc_s_forward(&self.store.borrow(), &self.proj.borrow(), p);
138    }
139    fn inverse<P: TransformCoordinates>(&self, p: &mut P) {
140        eqc_s_inverse(&self.store.borrow(), &self.proj.borrow(), p);
141    }
142}
143
144/// Equidistant Cylindrical Spheroidal forward project
145pub fn eqc_s_forward<P: TransformCoordinates>(eqc: &Eqc, proj: &Proj, p: &mut P) {
146    p.set_x(eqc.rc * p.lam());
147    p.set_y(p.phi() - proj.phi0);
148}
149
150/// Equidistant Cylindrical Spheroidal inverse project
151pub fn eqc_s_inverse<P: TransformCoordinates>(eqc: &Eqc, proj: &Proj, p: &mut P) {
152    p.set_lam(p.x() / eqc.rc);
153    p.set_phi(p.y() + proj.phi0);
154}
gistools/proj/project/eqc.rs

gistools/proj/project/
eqc.rs
