gistools/proj/project/
tcc.rs

1use crate::proj::{CoordinateStep, EPS10, Proj, ProjectCoordinates, TransformCoordinates};
2use alloc::rc::Rc;
3use core::cell::RefCell;
4use libm::{atan2, cos, sin, sqrt, tan};
5
6/// Transverse Central Cylindrical Projection
7#[derive(Debug, Clone, PartialEq)]
8pub struct TransverseCentralCylindricalProjection {
9    proj: Rc<RefCell<Proj>>,
10}
11impl ProjectCoordinates for TransverseCentralCylindricalProjection {
12    fn code(&self) -> i64 {
13        -1
14    }
15    fn name(&self) -> &'static str {
16        "Transverse Central Cylindrical"
17    }
18    fn names() -> &'static [&'static str] {
19        &["Transverse Central Cylindrical", "tcc"]
20    }
21}
22impl CoordinateStep for TransverseCentralCylindricalProjection {
23    fn new(proj: Rc<RefCell<Proj>>) -> Self {
24        proj.borrow_mut().es = 0.;
25        TransverseCentralCylindricalProjection { proj }
26    }
27    fn forward<P: TransformCoordinates>(&self, p: &mut P) {
28        tcc_s_forward(p);
29    }
30    fn inverse<P: TransformCoordinates>(&self, p: &mut P) {
31        tcc_s_inverse(p);
32    }
33}
34
35/// Transverse Central Cylindrical forward project
36/// let b = cos(φ) * sin(λ);
37/// x = b / sqrt(1 - b²);
38/// y = atan2(tan(φ), cos(λ));
39pub fn tcc_s_forward<P: TransformCoordinates>(p: &mut P) {
40    let b = cos(p.phi()) * sin(p.lam());
41    let bt = 1. - b * b;
42    if bt < EPS10 {
43        panic!("Coordinate outside projection domain");
44    }
45    p.set_x(b / sqrt(bt));
46    p.set_y(atan2(tan(p.phi()), cos(p.lam())));
47}
48
49/// Transverse Central Cylindrical inverse project
50pub fn tcc_s_inverse<P: TransformCoordinates>(p: &mut P) {
51    let x = p.x();
52    let y = p.y();
53
54    let denom = sqrt(1. + x * x);
55    let phi = (y.sin() / denom).atan();
56    let lam = (x / y.cos()).atan();
57
58    p.set_phi(phi);
59    p.set_lam(lam);
60}