map_engine/
affine.rs

1//Copyright (c) 2014, Sean C. Gillies
2//All rights reserved.
3//
4//Redistribution and use in source and binary forms, with or without
5//modification, are permitted provided that the following conditions are met:
6//
7//    * Redistributions of source code must retain the above copyright
8//      notice, this list of conditions and the following disclaimer.
9//    * Redistributions in binary form must reproduce the above copyright
10//      notice, this list of conditions and the following disclaimer in the
11//      documentation and/or other materials provided with the distribution.
12//    * Neither the name of Sean C. Gillies nor the names of
13//      its contributors may be used to endorse or promote products derived from
14//      this software without specific prior written permission.
15//
16//THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
17//AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
18//IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
19//ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE
20//LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
21//CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
22//SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
23//INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
24//CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
25//ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
26//POSSIBILITY OF SUCH DAMAGE.
27/*!Matrices describing affine transformation of the plane.
28
29Many of the functions are extracted from <https://github.com/rasterio/affine>.
30*/
31use crate::errors::MapEngineError;
32use num_traits::AsPrimitive;
33use serde::{Deserialize, Serialize};
34use std::ops::Mul;
35
36#[derive(Debug, PartialEq, Clone, Serialize, Deserialize)]
37pub struct GeoTransform {
38    pub geo: [f64; 6],
39}
40
41impl GeoTransform {
42    pub fn new(geo: &[f64; 6]) -> Self {
43        Self { geo: *geo }
44    }
45    #[allow(clippy::many_single_char_names)]
46    pub fn from_gdal(geo: &[f64; 6]) -> Self {
47        let c = geo[0];
48        let a = geo[1];
49        let b = geo[2];
50        let f = geo[3];
51        let d = geo[4];
52        let e = geo[5];
53        Self::new(&[a, b, c, d, e, f])
54    }
55
56    pub fn translation(tx: f64, ty: f64) -> Self {
57        let geo = [1.0, 0.0, tx, 0.0, 1.0, ty];
58        Self::new(&geo)
59    }
60
61    pub fn scale(sx: f64, sy: f64) -> Self {
62        let geo = [sx, 0.0, 0.0, 0.0, sy, 0.0];
63        Self::new(&geo)
64    }
65
66    pub fn shear(sx: f64, sy: f64) -> Self {
67        let geo = [
68            1.0,
69            sx.to_radians().tan(),
70            0.0,
71            sy.to_radians().tan(),
72            1.0,
73            0.0,
74        ];
75        Self::new(&geo)
76    }
77
78    #[allow(clippy::many_single_char_names)]
79    pub fn to_tuple(&self) -> (f64, f64, f64, f64, f64, f64) {
80        let a = self.geo[0];
81        let b = self.geo[1];
82        let c = self.geo[2];
83        let d = self.geo[3];
84        let e = self.geo[4];
85        let f = self.geo[5];
86        (a, b, c, d, e, f)
87    }
88
89    pub fn to_array(&self) -> [f64; 6] {
90        self.geo
91    }
92
93    pub fn inv(&self) -> Result<Self, MapEngineError> {
94        let (sa, sb, sc, sd, se, sf) = self.to_tuple();
95
96        let determinant = self.determinant();
97        if determinant == 0.0 {
98            return Err(MapEngineError::AffineError("Determinant is zero".into()));
99        }
100
101        let idet = 1.0 / determinant;
102
103        let ra = se * idet;
104        let rb = -sb * idet;
105        let rd = -sd * idet;
106        let re = sa * idet;
107
108        let geo = [ra, rb, -sc * ra - sf * rb, rd, re, -sc * rd - sf * re];
109        Ok(Self::new(&geo))
110    }
111
112    pub fn determinant(&self) -> f64 {
113        let (a, b, _, d, e, _) = self.to_tuple();
114        a * e - b * d
115    }
116
117    pub fn xy(&self, row: u32, col: u32) -> (f64, f64) {
118        let (coff, roff) = (0.5, 0.5);
119        let trans = Self::translation(coff, roff);
120        let tmp = self * &trans;
121        &tmp * (col as f64, row as f64)
122    }
123
124    pub fn rowcol(&self, x: f64, y: f64) -> Result<(i32, i32), MapEngineError> {
125        let eps = 2.220446049250313e-16;
126
127        let inv = self.inv()?;
128        let (fcol, frow) = &inv * (x + eps, y + eps);
129
130        Ok((frow.floor() as i32, fcol.floor() as i32))
131    }
132
133    #[allow(clippy::many_single_char_names)]
134    pub fn to_gdal(&self) -> Self {
135        let (a, b, c, d, e, f) = self.to_tuple();
136        let geo = [c, a, b, f, d, e];
137        Self::new(&geo)
138    }
139
140    pub fn xoff(&self) -> f64 {
141        self.geo[2]
142    }
143    pub fn yoff(&self) -> f64 {
144        self.geo[5]
145    }
146}
147
148impl Mul<&GeoTransform> for &GeoTransform {
149    type Output = GeoTransform;
150
151    fn mul(self, rhs: &GeoTransform) -> GeoTransform {
152        let (sa, sb, sc, sd, se, sf) = self.to_tuple();
153        let (oa, ob, oc, od, oe, of) = rhs.to_tuple();
154        let geo = [
155            sa * oa + sb * od,
156            sa * ob + sb * oe,
157            sa * oc + sb * of + sc,
158            sd * oa + se * od,
159            sd * ob + se * oe,
160            sd * oc + se * of + sf,
161        ];
162        GeoTransform::new(&geo)
163    }
164}
165
166impl<T> Mul<(T, T)> for &GeoTransform
167where
168    T: AsPrimitive<f64>,
169{
170    type Output = (f64, f64);
171
172    fn mul(self, rhs: (T, T)) -> Self::Output {
173        let (sa, sb, sc, sd, se, sf) = self.to_tuple();
174        let (vx, vy) = rhs;
175        let vx: f64 = vx.as_();
176        let vy: f64 = vy.as_();
177        (vx * sa + vy * sb + sc, vx * sd + vy * se + sf)
178    }
179}
180
181impl From<gdal::GeoTransform> for GeoTransform {
182    #[allow(clippy::many_single_char_names)]
183    fn from(geo: gdal::GeoTransform) -> Self {
184        let c = geo[0];
185        let a = geo[1];
186        let b = geo[2];
187        let f = geo[3];
188        let d = geo[4];
189        let e = geo[5];
190        Self::new(&[a, b, c, d, e, f])
191    }
192}
map_engine/affine.rs

map_engine/
affine.rs
