agent_tk/
projection.rs

1//! Handle conversion between WGS84 and a carthesian representation
2
3/// Projection from UTM to WGS84
4#[derive(Clone)]
5pub struct Projection
6{
7  projection_utm: proj4rs::proj::Proj,
8  projection_wgs84: proj4rs::proj::Proj,
9  offset_x: f64,
10  offset_y: f64,
11}
12
13fn utm_zone(lon: f64, lat: f64) -> i32
14{
15  if (56.0..64.0).contains(&lat) && (3.0..12.0).contains(&lon)
16  {
17    return 32;
18  }
19  // Special zones for Svalbard
20  else if (72.0..84.0).contains(&lat)
21  {
22    if (0.0..9.0).contains(&lon)
23    {
24      return 31;
25    }
26    else if (9.0..21.0).contains(&lon)
27    {
28      return 33;
29    }
30    else if (21.0..33.0).contains(&lon)
31    {
32      return 35;
33    }
34    else if (33.0..42.0).contains(&lon)
35    {
36      return 37;
37    }
38  }
39  ((lon + 180.0) / 6.0) as i32 + 1
40}
41
42impl Projection
43{
44  /// Create a new projection centered around (lon, lat).
45  pub fn new(lon: f64, lat: f64) -> Projection
46  {
47    let utm_z: i32 = utm_zone(lon, lat);
48    let proj_str = format!(
49      "+proj=utm +zone={} +datum=WGS84 +units=m +no_defs +type=crs",
50      utm_z
51    );
52    let proj_from = proj4rs::proj::Proj::from_proj_string(&proj_str[..]).unwrap();
53    let proj_to = proj4rs::proj::Proj::from_proj_string(concat!(
54      "+proj=longlat +ellps=WGS84",
55      " +datum=WGS84 +no_defs"
56    ))
57    .unwrap();
58
59    let mut point_3d = (lon.to_radians(), lat.to_radians(), 0.0);
60    proj4rs::transform::transform(&proj_to, &proj_from, &mut point_3d).unwrap();
61
62    Projection {
63      projection_utm: proj_from,
64      projection_wgs84: proj_to,
65      offset_x: point_3d.0,
66      offset_y: point_3d.1,
67    }
68  }
69  /// Convert from (lon, lat) to (x, y)
70  pub fn from_wgs84(&self, lon: f64, lat: f64) -> (f64, f64)
71  {
72    let mut point_3d = (lon.to_radians(), lat.to_radians(), 0.0);
73    proj4rs::transform::transform(&self.projection_wgs84, &self.projection_utm, &mut point_3d)
74      .unwrap();
75    (point_3d.0 - self.offset_x, point_3d.1 - self.offset_y)
76  }
77  /// Convert from (x, y) to (lon, lat)
78  pub fn to_wgs84(&self, x: f64, y: f64) -> (f64, f64)
79  {
80    let mut point_3d = (x + self.offset_x, y + self.offset_y, 0.0);
81    proj4rs::transform::transform(&self.projection_utm, &self.projection_wgs84, &mut point_3d)
82      .unwrap();
83    (point_3d.0.to_degrees(), point_3d.1.to_degrees())
84  }
85}
agent_tk/projection.rs

agent_tk/
projection.rs
