geo/algorithm/
affine_ops.rs

1use crate::{Coord, CoordFloat, CoordNum, MapCoords, MapCoordsInPlace};
2use std::fmt;
3
4/// Apply an [`AffineTransform`] like [`scale`](AffineTransform::scale),
5/// [`skew`](AffineTransform::skew), or [`rotate`](AffineTransform::rotate) to a
6/// [`Geometry`](crate::geometry::Geometry).
7///
8/// Multiple transformations can be composed in order to be efficiently applied in a single
9/// operation. See [`AffineTransform`] for more on how to build up a transformation.
10///
11/// If you are not composing operations, traits that leverage this same machinery exist which might
12/// be more readable. See: [`Scale`](crate::algorithm::Scale),
13/// [`Translate`](crate::algorithm::Translate), [`Rotate`](crate::algorithm::Rotate),
14/// and [`Skew`](crate::algorithm::Skew).
15///
16/// # Examples
17/// ## Build up transforms by beginning with a constructor, then chaining mutation operations
18/// ```
19/// use geo::{AffineOps, AffineTransform};
20/// use geo::{line_string, BoundingRect, Point, LineString};
21/// use approx::assert_relative_eq;
22///
23/// let ls: LineString = line_string![
24///     (x: 0.0f64, y: 0.0f64),
25///     (x: 0.0f64, y: 10.0f64),
26/// ];
27/// let center = ls.bounding_rect().unwrap().center();
28///
29/// let transform = AffineTransform::skew(40.0, 40.0, center).rotated(45.0, center);
30///
31/// let skewed_rotated = ls.affine_transform(&transform);
32///
33/// assert_relative_eq!(skewed_rotated, line_string![
34///     (x: 0.5688687f64, y: 4.4311312),
35///     (x: -0.5688687, y: 5.5688687)
36/// ], max_relative = 1.0);
37/// ```
38pub trait AffineOps<T: CoordNum> {
39    /// Apply `transform` immutably, outputting a new geometry.
40    #[must_use]
41    fn affine_transform(&self, transform: &AffineTransform<T>) -> Self;
42
43    /// Apply `transform` to mutate `self`.
44    fn affine_transform_mut(&mut self, transform: &AffineTransform<T>);
45}
46
47impl<T: CoordNum, M: MapCoordsInPlace<T> + MapCoords<T, T, Output = Self>> AffineOps<T> for M {
48    fn affine_transform(&self, transform: &AffineTransform<T>) -> Self {
49        self.map_coords(|c| transform.apply(c))
50    }
51
52    fn affine_transform_mut(&mut self, transform: &AffineTransform<T>) {
53        self.map_coords_in_place(|c| transform.apply(c))
54    }
55}
56
57/// A general affine transformation matrix, and associated operations.
58///
59/// Note that affine ops are **already implemented** on most `geo-types` primitives, using this module.
60///
61/// Affine transforms using the same numeric type (e.g. [`CoordFloat`](crate::CoordFloat)) can be **composed**,
62/// and the result can be applied to geometries using e.g. [`MapCoords`](crate::MapCoords). This allows the
63/// efficient application of transforms: an arbitrary number of operations can be chained.
64/// These are then composed, producing a final transformation matrix which is applied to the geometry coordinates.
65///
66/// `AffineTransform` is a row-major matrix.
67/// 2D affine transforms require six matrix parameters:
68///
69/// `[a, b, xoff, d, e, yoff]`
70///
71/// these map onto the `AffineTransform` rows as follows:
72/// ```ignore
73/// [[a, b, xoff],
74/// [d, e, yoff],
75/// [0, 0, 1]]
76/// ```
77/// The equations for transforming coordinates `(x, y) -> (x', y')` are given as follows:
78///
79/// `x' = ax + by + xoff`
80///
81/// `y' = dx + ey + yoff`
82///
83/// # Usage
84///
85/// Two types of operation are provided: construction and mutation. **Construction** functions create a *new* transform
86/// and are denoted by the use of the **present tense**: `scale()`, `translate()`, `rotate()`, and `skew()`.
87///
88/// **Mutation** methods *add* a transform to the existing `AffineTransform`, and are denoted by the use of the past participle:
89/// `scaled()`, `translated()`, `rotated()`, and `skewed()`.
90///
91/// # Examples
92/// ## Build up transforms by beginning with a constructor, then chaining mutation operations
93/// ```
94/// use geo::{AffineOps, AffineTransform};
95/// use geo::{line_string, BoundingRect, Point, LineString};
96/// use approx::assert_relative_eq;
97///
98/// let ls: LineString = line_string![
99///     (x: 0.0f64, y: 0.0f64),
100///     (x: 0.0f64, y: 10.0f64),
101/// ];
102/// let center = ls.bounding_rect().unwrap().center();
103///
104/// let transform = AffineTransform::skew(40.0, 40.0, center).rotated(45.0, center);
105///
106/// let skewed_rotated = ls.affine_transform(&transform);
107///
108/// assert_relative_eq!(skewed_rotated, line_string![
109///     (x: 0.5688687f64, y: 4.4311312),
110///     (x: -0.5688687, y: 5.5688687)
111/// ], max_relative = 1.0);
112/// ```
113#[derive(Copy, Clone, PartialEq, Eq)]
114pub struct AffineTransform<T: CoordNum = f64>([[T; 3]; 3]);
115
116impl<T: CoordNum> Default for AffineTransform<T> {
117    fn default() -> Self {
118        // identity matrix
119        Self::identity()
120    }
121}
122
123impl<T: CoordNum> AffineTransform<T> {
124    /// Create a new affine transformation by composing two `AffineTransform`s.
125    ///
126    /// This is a **cumulative** operation; the new transform is *added* to the existing transform.
127    #[must_use]
128    pub fn compose(&self, other: &Self) -> Self {
129        // lol
130        Self([
131            [
132                (self.0[0][0] * other.0[0][0])
133                    + (self.0[0][1] * other.0[1][0])
134                    + (self.0[0][2] * other.0[2][0]),
135                (self.0[0][0] * other.0[0][1])
136                    + (self.0[0][1] * other.0[1][1])
137                    + (self.0[0][2] * other.0[2][1]),
138                (self.0[0][0] * other.0[0][2])
139                    + (self.0[0][1] * other.0[1][2])
140                    + (self.0[0][2] * other.0[2][2]),
141            ],
142            [
143                (self.0[1][0] * other.0[0][0])
144                    + (self.0[1][1] * other.0[1][0])
145                    + (self.0[1][2] * other.0[2][0]),
146                (self.0[1][0] * other.0[0][1])
147                    + (self.0[1][1] * other.0[1][1])
148                    + (self.0[1][2] * other.0[2][1]),
149                (self.0[1][0] * other.0[0][2])
150                    + (self.0[1][1] * other.0[1][2])
151                    + (self.0[1][2] * other.0[2][2]),
152            ],
153            [
154                // this section isn't technically necessary since the last row is invariant: [0, 0, 1]
155                (self.0[2][0] * other.0[0][0])
156                    + (self.0[2][1] * other.0[1][0])
157                    + (self.0[2][2] * other.0[2][0]),
158                (self.0[2][0] * other.0[0][1])
159                    + (self.0[2][1] * other.0[1][1])
160                    + (self.0[2][2] * other.0[2][1]),
161                (self.0[2][0] * other.0[0][2])
162                    + (self.0[2][1] * other.0[1][2])
163                    + (self.0[2][2] * other.0[2][2]),
164            ],
165        ])
166    }
167    /// Create the identity matrix
168    ///
169    /// The matrix is:
170    /// ```ignore
171    /// [[1, 0, 0],
172    /// [0, 1, 0],
173    /// [0, 0, 1]]
174    /// ```
175    pub fn identity() -> Self {
176        Self::new(
177            T::one(),
178            T::zero(),
179            T::zero(),
180            T::zero(),
181            T::one(),
182            T::zero(),
183        )
184    }
185
186    /// Whether the transformation is equivalent to the [identity matrix](Self::identity),
187    /// that is, whether it's application will be a a no-op.
188    ///
189    /// ```
190    /// use geo::AffineTransform;
191    /// let mut transform = AffineTransform::identity();
192    /// assert!(transform.is_identity());
193    ///
194    /// // mutate the transform a bit
195    /// transform = transform.translated(1.0, 2.0);
196    /// assert!(!transform.is_identity());
197    ///
198    /// // put it back
199    /// transform = transform.translated(-1.0, -2.0);
200    /// assert!(transform.is_identity());
201    /// ```
202    pub fn is_identity(&self) -> bool {
203        self == &Self::identity()
204    }
205
206    /// **Create** a new affine transform for scaling, scaled by factors along the `x` and `y` dimensions.
207    /// The point of origin is *usually* given as the 2D bounding box centre of the geometry, but
208    /// any coordinate may be specified.
209    /// Negative scale factors will mirror or reflect coordinates.
210    ///
211    /// The matrix is:
212    /// ```ignore
213    /// [[xfact, 0, xoff],
214    /// [0, yfact, yoff],
215    /// [0, 0, 1]]
216    ///
217    /// xoff = origin.x - (origin.x * xfact)
218    /// yoff = origin.y - (origin.y * yfact)
219    /// ```
220    pub fn scale(xfact: T, yfact: T, origin: impl Into<Coord<T>>) -> Self {
221        let (x0, y0) = origin.into().x_y();
222        let xoff = x0 - (x0 * xfact);
223        let yoff = y0 - (y0 * yfact);
224        Self::new(xfact, T::zero(), xoff, T::zero(), yfact, yoff)
225    }
226
227    /// **Add** an affine transform for scaling, scaled by factors along the `x` and `y` dimensions.
228    /// The point of origin is *usually* given as the 2D bounding box centre of the geometry, but
229    /// any coordinate may be specified.
230    /// Negative scale factors will mirror or reflect coordinates.
231    /// This is a **cumulative** operation; the new transform is *added* to the existing transform.
232    #[must_use]
233    pub fn scaled(mut self, xfact: T, yfact: T, origin: impl Into<Coord<T>>) -> Self {
234        self.0 = self.compose(&Self::scale(xfact, yfact, origin)).0;
235        self
236    }
237
238    /// **Create** an affine transform for translation, shifted by offsets along the `x` and `y` dimensions.
239    ///
240    /// The matrix is:
241    /// ```ignore
242    /// [[1, 0, xoff],
243    /// [0, 1, yoff],
244    /// [0, 0, 1]]
245    /// ```
246    pub fn translate(xoff: T, yoff: T) -> Self {
247        Self::new(T::one(), T::zero(), xoff, T::zero(), T::one(), yoff)
248    }
249
250    /// **Add** an affine transform for translation, shifted by offsets along the `x` and `y` dimensions
251    ///
252    /// This is a **cumulative** operation; the new transform is *added* to the existing transform.
253    #[must_use]
254    pub fn translated(mut self, xoff: T, yoff: T) -> Self {
255        self.0 = self.compose(&Self::translate(xoff, yoff)).0;
256        self
257    }
258
259    /// Apply the current transform to a coordinate
260    pub fn apply(&self, coord: Coord<T>) -> Coord<T> {
261        Coord {
262            x: (self.0[0][0] * coord.x + self.0[0][1] * coord.y + self.0[0][2]),
263            y: (self.0[1][0] * coord.x + self.0[1][1] * coord.y + self.0[1][2]),
264        }
265    }
266
267    /// Create a new custom transform matrix
268    ///
269    /// The argument order matches that of the affine transform matrix:
270    ///```ignore
271    /// [[a, b, xoff],
272    /// [d, e, yoff],
273    /// [0, 0, 1]] <-- not part of the input arguments
274    /// ```
275    pub fn new(a: T, b: T, xoff: T, d: T, e: T, yoff: T) -> Self {
276        Self([[a, b, xoff], [d, e, yoff], [T::zero(), T::zero(), T::one()]])
277    }
278}
279
280impl<T: CoordNum> fmt::Debug for AffineTransform<T> {
281    fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
282        f.debug_struct("AffineTransform")
283            .field("a", &self.0[0][0])
284            .field("b", &self.0[0][1])
285            .field("xoff", &self.0[1][2])
286            .field("d", &self.0[1][0])
287            .field("e", &self.0[1][1])
288            .field("yoff", &self.0[1][2])
289            .finish()
290    }
291}
292
293impl<T: CoordNum> From<[T; 6]> for AffineTransform<T> {
294    fn from(arr: [T; 6]) -> Self {
295        Self::new(arr[0], arr[1], arr[2], arr[3], arr[4], arr[5])
296    }
297}
298
299impl<T: CoordNum> From<(T, T, T, T, T, T)> for AffineTransform<T> {
300    fn from(tup: (T, T, T, T, T, T)) -> Self {
301        Self::new(tup.0, tup.1, tup.2, tup.3, tup.4, tup.5)
302    }
303}
304
305impl<U: CoordFloat> AffineTransform<U> {
306    /// **Create** an affine transform for rotation, using an arbitrary point as its centre.
307    ///
308    /// Note that this operation is only available for geometries with floating point coordinates.
309    ///
310    /// `angle` is given in **degrees**.
311    ///
312    /// The matrix (angle denoted as theta) is:
313    /// ```ignore
314    /// [[cos_theta, -sin_theta, xoff],
315    /// [sin_theta, cos_theta, yoff],
316    /// [0, 0, 1]]
317    ///
318    /// xoff = origin.x - (origin.x * cos(theta)) + (origin.y * sin(theta))
319    /// yoff = origin.y - (origin.x * sin(theta)) + (origin.y * cos(theta))
320    /// ```
321    pub fn rotate(degrees: U, origin: impl Into<Coord<U>>) -> Self {
322        let (sin_theta, cos_theta) = degrees.to_radians().sin_cos();
323        let (x0, y0) = origin.into().x_y();
324        let xoff = x0 - (x0 * cos_theta) + (y0 * sin_theta);
325        let yoff = y0 - (x0 * sin_theta) - (y0 * cos_theta);
326        Self::new(cos_theta, -sin_theta, xoff, sin_theta, cos_theta, yoff)
327    }
328
329    /// **Add** an affine transform for rotation, using an arbitrary point as its centre.
330    ///
331    /// Note that this operation is only available for geometries with floating point coordinates.
332    ///
333    /// `angle` is given in **degrees**.
334    ///
335    /// This is a **cumulative** operation; the new transform is *added* to the existing transform.
336    #[must_use]
337    pub fn rotated(mut self, angle: U, origin: impl Into<Coord<U>>) -> Self {
338        self.0 = self.compose(&Self::rotate(angle, origin)).0;
339        self
340    }
341
342    /// **Create** an affine transform for skewing.
343    ///
344    /// Note that this operation is only available for geometries with floating point coordinates.
345    ///
346    /// Geometries are sheared by angles along x (`xs`) and y (`ys`) dimensions.
347    /// The point of origin is *usually* given as the 2D bounding box centre of the geometry, but
348    /// any coordinate may be specified. Angles are given in **degrees**.
349    /// The matrix is:
350    /// ```ignore
351    /// [[1, tan(x), xoff],
352    /// [tan(y), 1, yoff],
353    /// [0, 0, 1]]
354    ///
355    /// xoff = -origin.y * tan(xs)
356    /// yoff = -origin.x * tan(ys)
357    /// ```
358    pub fn skew(xs: U, ys: U, origin: impl Into<Coord<U>>) -> Self {
359        let Coord { x: x0, y: y0 } = origin.into();
360        let mut tanx = xs.to_radians().tan();
361        let mut tany = ys.to_radians().tan();
362        // These checks are stolen from Shapely's implementation -- may not be necessary
363        if tanx.abs() < U::from::<f64>(2.5e-16).unwrap() {
364            tanx = U::zero();
365        }
366        if tany.abs() < U::from::<f64>(2.5e-16).unwrap() {
367            tany = U::zero();
368        }
369        let xoff = -y0 * tanx;
370        let yoff = -x0 * tany;
371        Self::new(U::one(), tanx, xoff, tany, U::one(), yoff)
372    }
373
374    /// **Add** an affine transform for skewing.
375    ///
376    /// Note that this operation is only available for geometries with floating point coordinates.
377    ///
378    /// Geometries are sheared by angles along x (`xs`) and y (`ys`) dimensions.
379    /// The point of origin is *usually* given as the 2D bounding box centre of the geometry, but
380    /// any coordinate may be specified. Angles are given in **degrees**.
381    ///
382    /// This is a **cumulative** operation; the new transform is *added* to the existing transform.
383    #[must_use]
384    pub fn skewed(mut self, xs: U, ys: U, origin: impl Into<Coord<U>>) -> Self {
385        self.0 = self.compose(&Self::skew(xs, ys, origin)).0;
386        self
387    }
388}
389
390#[cfg(test)]
391mod tests {
392    use super::*;
393    use crate::{polygon, Point};
394
395    // given a matrix with the shape
396    // [[a, b, xoff],
397    // [d, e, yoff],
398    // [0, 0, 1]]
399    #[test]
400    fn matrix_multiply() {
401        let a = AffineTransform::new(1, 2, 5, 3, 4, 6);
402        let b = AffineTransform::new(7, 8, 11, 9, 10, 12);
403        let composed = a.compose(&b);
404        assert_eq!(composed.0[0][0], 25);
405        assert_eq!(composed.0[0][1], 28);
406        assert_eq!(composed.0[0][2], 40);
407        assert_eq!(composed.0[1][0], 57);
408        assert_eq!(composed.0[1][1], 64);
409        assert_eq!(composed.0[1][2], 87);
410    }
411    #[test]
412    fn test_transform_composition() {
413        let p0 = Point::new(0.0f64, 0.0);
414        // scale once
415        let mut scale_a = AffineTransform::default().scaled(2.0, 2.0, p0);
416        // rotate
417        scale_a = scale_a.rotated(45.0, p0);
418        // rotate back
419        scale_a = scale_a.rotated(-45.0, p0);
420        // scale up again, doubling
421        scale_a = scale_a.scaled(2.0, 2.0, p0);
422        // scaled once
423        let scale_b = AffineTransform::default().scaled(2.0, 2.0, p0);
424        // scaled once, but equal to 2 + 2
425        let scale_c = AffineTransform::default().scaled(4.0, 4.0, p0);
426        assert_ne!(&scale_a.0, &scale_b.0);
427        assert_eq!(&scale_a.0, &scale_c.0);
428    }
429
430    #[test]
431    fn affine_transformed() {
432        let transform = AffineTransform::translate(1.0, 1.0).scaled(2.0, 2.0, (0.0, 0.0));
433        let mut poly = polygon![(x: 0.0, y: 0.0), (x: 0.0, y: 2.0), (x: 1.0, y: 2.0)];
434        poly.affine_transform_mut(&transform);
435
436        let expected = polygon![(x: 1.0, y: 1.0), (x: 1.0, y: 5.0), (x: 3.0, y: 5.0)];
437        assert_eq!(expected, poly);
438    }
439}
geo/algorithm/affine_ops.rs

geo/algorithm/
affine_ops.rs
