Skip to main content

affine_rs/
lib.rs

1//! Two-dimensional affine transformations.
2//!
3//! An [`Affine`] stores the six effective coefficients of the augmented
4//! matrix below. The final row is fixed and cannot be put into an invalid
5//! state.
6//!
7//! ```text
8//! | a  b  c |
9//! | d  e  f |
10//! | 0  0  1 |
11//! ```
12//!
13//! # Example
14//!
15//! ```
16//! use affine_rs::Affine;
17//!
18//! let transform = Affine::translation(10.0, 20.0)
19//!     .compose(Affine::scale(2.0, 2.0));
20//! assert_eq!(transform.transform_point([3.0, 4.0]), [16.0, 28.0]);
21//! ```
22
23mod error;
24
25use core::fmt;
26use core::ops::Mul;
27
28pub use error::AffineError;
29
30/// Default tolerance for approximate comparisons.
31pub const DEFAULT_EPSILON: f64 = 1.0e-5;
32const DEFAULT_EPSILON_SQUARED: f64 = DEFAULT_EPSILON * DEFAULT_EPSILON;
33
34/// A two-dimensional affine transform.
35///
36/// The coefficients map an input point `[x, y]` to:
37///
38/// ```text
39/// x' = a*x + b*y + c
40/// y' = d*x + e*y + f
41/// ```
42#[derive(Clone, Copy, Debug, PartialEq)]
43pub struct Affine {
44    pub a: f64,
45    pub b: f64,
46    pub c: f64,
47    pub d: f64,
48    pub e: f64,
49    pub f: f64,
50}
51
52impl Affine {
53    /// The identity transform.
54    pub const IDENTITY: Self = Self::new(1.0, 0.0, 0.0, 0.0, 1.0, 0.0);
55
56    /// Creates a transform from its six effective coefficients.
57    #[must_use]
58    #[allow(clippy::many_single_char_names)]
59    pub const fn new(a: f64, b: f64, c: f64, d: f64, e: f64, f: f64) -> Self {
60        Self { a, b, c, d, e, f }
61    }
62
63    /// Creates a translation transform.
64    #[must_use]
65    pub const fn translation(x_offset: f64, y_offset: f64) -> Self {
66        Self::new(1.0, 0.0, x_offset, 0.0, 1.0, y_offset)
67    }
68
69    /// Creates independent x- and y-axis scaling.
70    #[must_use]
71    pub const fn scale(x_scale: f64, y_scale: f64) -> Self {
72        Self::new(x_scale, 0.0, 0.0, 0.0, y_scale, 0.0)
73    }
74
75    /// Creates uniform scaling on both axes.
76    #[must_use]
77    pub const fn uniform_scale(scale: f64) -> Self {
78        Self::scale(scale, scale)
79    }
80
81    /// Creates x- and y-axis shear from angles in degrees.
82    #[must_use]
83    pub fn shear(x_angle_degrees: f64, y_angle_degrees: f64) -> Self {
84        let x_shear = x_angle_degrees.to_radians().tan();
85        let y_shear = y_angle_degrees.to_radians().tan();
86        Self::new(1.0, x_shear, 0.0, y_shear, 1.0, 0.0)
87    }
88
89    /// Creates a counter-clockwise rotation around the origin.
90    #[must_use]
91    pub fn rotation(angle_degrees: f64) -> Self {
92        let (cosine, sine) = cos_sin_degrees(angle_degrees);
93        Self::new(cosine, -sine, 0.0, sine, cosine, 0.0)
94    }
95
96    /// Creates a counter-clockwise rotation around `pivot`.
97    #[must_use]
98    pub fn rotation_around(angle_degrees: f64, pivot: [f64; 2]) -> Self {
99        let (cosine, sine) = cos_sin_degrees(angle_degrees);
100        let [pivot_x, pivot_y] = pivot;
101        Self::new(
102            cosine,
103            -sine,
104            pivot_x - pivot_x * cosine + pivot_y * sine,
105            sine,
106            cosine,
107            pivot_y - pivot_x * sine - pivot_y * cosine,
108        )
109    }
110
111    /// Creates the non-identity permutation matrix for two dimensions.
112    #[must_use]
113    pub const fn permutation() -> Self {
114        Self::new(0.0, 1.0, 0.0, 1.0, 0.0, 0.0)
115    }
116
117    /// Creates a transform from GDAL `GeoTransform` coefficient order.
118    #[must_use]
119    #[allow(clippy::many_single_char_names)]
120    pub const fn from_gdal(coefficients: [f64; 6]) -> Self {
121        let [c, a, b, f, d, e] = coefficients;
122        Self::new(a, b, c, d, e, f)
123    }
124
125    /// Returns coefficients in GDAL `GeoTransform` order.
126    #[must_use]
127    pub const fn to_gdal(self) -> [f64; 6] {
128        [self.c, self.a, self.b, self.f, self.d, self.e]
129    }
130
131    /// Returns coefficients in Shapely affine-transform order.
132    #[must_use]
133    pub const fn to_shapely(self) -> [f64; 6] {
134        [self.a, self.b, self.d, self.e, self.c, self.f]
135    }
136
137    /// Returns the area scaling factor of the linear part.
138    #[must_use]
139    pub fn determinant(self) -> f64 {
140        self.a.mul_add(self.e, -(self.b * self.d))
141    }
142
143    /// Returns `true` when all stored coefficients are finite.
144    #[must_use]
145    pub fn is_finite(self) -> bool {
146        [self.a, self.b, self.c, self.d, self.e, self.f]
147            .into_iter()
148            .all(f64::is_finite)
149    }
150
151    /// Returns `true` when the determinant is exactly zero.
152    #[must_use]
153    pub fn is_degenerate(self) -> bool {
154        self.determinant() == 0.0
155    }
156
157    /// Returns `true` when the transform preserves orientation.
158    #[must_use]
159    pub fn is_proper(self) -> bool {
160        self.determinant() > 0.0
161    }
162
163    /// Returns the matrix as its three two-dimensional column vectors.
164    #[must_use]
165    pub const fn column_vectors(self) -> [[f64; 2]; 3] {
166        [[self.a, self.d], [self.b, self.e], [self.c, self.f]]
167    }
168
169    /// Returns `true` when the transform is approximately the identity.
170    #[must_use]
171    pub fn is_identity(self) -> bool {
172        self.approx_eq(Self::IDENTITY, DEFAULT_EPSILON)
173    }
174
175    /// Returns `true` when axis-aligned shapes remain axis-aligned.
176    #[must_use]
177    pub fn is_rectilinear(self) -> bool {
178        (self.a.abs() < DEFAULT_EPSILON && self.e.abs() < DEFAULT_EPSILON)
179            || (self.d.abs() < DEFAULT_EPSILON && self.b.abs() < DEFAULT_EPSILON)
180    }
181
182    /// Returns `true` when the linear transform preserves angles.
183    #[must_use]
184    pub fn is_conformal(self) -> bool {
185        self.a.mul_add(self.b, self.d * self.e).abs() < DEFAULT_EPSILON
186    }
187
188    /// Returns `true` when the linear transform is orthonormal.
189    #[must_use]
190    pub fn is_orthonormal(self) -> bool {
191        self.is_conformal()
192            && (1.0 - self.a.mul_add(self.a, self.d * self.d)).abs() < DEFAULT_EPSILON
193            && (1.0 - self.b.mul_add(self.b, self.e * self.e)).abs() < DEFAULT_EPSILON
194    }
195
196    /// Returns the two singular values, ordered from greatest to least.
197    #[must_use]
198    pub fn scaling(self) -> [f64; 2] {
199        let trace = self.a * self.a + self.b * self.b + self.d * self.d + self.e * self.e;
200        let determinant_squared = self.determinant().powi(2);
201        let mut delta = trace.mul_add(trace / 4.0, -determinant_squared);
202        if delta < DEFAULT_EPSILON_SQUARED {
203            delta = 0.0;
204        }
205        let root = delta.sqrt();
206        [(trace / 2.0 + root).sqrt(), (trace / 2.0 - root).sqrt()]
207    }
208
209    /// Returns the eccentricity induced by the linear transform.
210    #[must_use]
211    pub fn eccentricity(self) -> f64 {
212        let [major, minor] = self.scaling();
213        major.mul_add(major, -(minor * minor)).sqrt() / major
214    }
215
216    /// Returns the counter-clockwise rotation angle in degrees.
217    ///
218    /// # Errors
219    ///
220    /// Returns [`AffineError::UndefinedRotation`] for improper or degenerate
221    /// transforms, and [`AffineError::NonFiniteTransform`] for non-finite input.
222    pub fn rotation_angle(self) -> Result<f64, AffineError> {
223        if !self.is_finite() {
224            return Err(AffineError::NonFiniteTransform);
225        }
226        if !self.is_proper() || self.is_degenerate() {
227            return Err(AffineError::UndefinedRotation);
228        }
229        let [major, _] = self.scaling();
230        Ok((self.d / major).atan2(self.a / major).to_degrees())
231    }
232
233    /// Applies the transform to an `[x, y]` point.
234    #[must_use]
235    pub fn transform_point(self, [x, y]: [f64; 2]) -> [f64; 2] {
236        [
237            x.mul_add(self.a, y.mul_add(self.b, self.c)),
238            x.mul_add(self.d, y.mul_add(self.e, self.f)),
239        ]
240    }
241
242    /// Applies the transform to every point in a mutable slice.
243    pub fn transform_points_in_place(self, points: &mut [[f64; 2]]) {
244        for point in points {
245            *point = self.transform_point(*point);
246        }
247    }
248
249    /// Composes two transforms.
250    ///
251    /// `self.compose(rhs)` applies `rhs` first and `self` second.
252    #[must_use]
253    pub fn compose(self, rhs: Self) -> Self {
254        Self::new(
255            self.a.mul_add(rhs.a, self.b * rhs.d),
256            self.a.mul_add(rhs.b, self.b * rhs.e),
257            self.a.mul_add(rhs.c, self.b.mul_add(rhs.f, self.c)),
258            self.d.mul_add(rhs.a, self.e * rhs.d),
259            self.d.mul_add(rhs.b, self.e * rhs.e),
260            self.d.mul_add(rhs.c, self.e.mul_add(rhs.f, self.f)),
261        )
262    }
263
264    /// Returns the inverse transform.
265    ///
266    /// # Errors
267    ///
268    /// Returns [`AffineError::NonFiniteTransform`] for non-finite
269    /// coefficients and [`AffineError::NonInvertibleTransform`] when the
270    /// determinant is zero.
271    pub fn inverse(self) -> Result<Self, AffineError> {
272        if !self.is_finite() {
273            return Err(AffineError::NonFiniteTransform);
274        }
275
276        let determinant = self.determinant();
277        if determinant == 0.0 {
278            return Err(AffineError::NonInvertibleTransform);
279        }
280
281        let inverse_determinant = determinant.recip();
282        let a = self.e * inverse_determinant;
283        let b = -self.b * inverse_determinant;
284        let d = -self.d * inverse_determinant;
285        let e = self.a * inverse_determinant;
286
287        Ok(Self::new(
288            a,
289            b,
290            -self.c.mul_add(a, self.f * b),
291            d,
292            e,
293            -self.c.mul_add(d, self.f * e),
294        ))
295    }
296
297    /// Compares all six coefficients using an absolute tolerance.
298    #[must_use]
299    pub fn approx_eq(self, other: Self, epsilon: f64) -> bool {
300        let lhs = [self.a, self.b, self.c, self.d, self.e, self.f];
301        let rhs = [other.a, other.b, other.c, other.d, other.e, other.f];
302        lhs.into_iter()
303            .zip(rhs)
304            .all(|(left, right)| (left - right).abs() < epsilon)
305    }
306
307    /// Parses six coefficients in World File order.
308    ///
309    /// World Files store the center of the upper-left pixel. The returned
310    /// transform uses its upper-left corner, so a half-pixel translation is
311    /// applied automatically.
312    ///
313    /// # Errors
314    ///
315    /// Returns [`AffineError::InvalidWorldFile`] unless the input contains
316    /// exactly six finite floating-point values.
317    #[allow(clippy::many_single_char_names)]
318    pub fn from_world_file(text: &str) -> Result<Self, AffineError> {
319        let values = text
320            .split_whitespace()
321            .map(str::parse::<f64>)
322            .collect::<Result<Vec<_>, _>>()
323            .map_err(|_| AffineError::InvalidWorldFile)?;
324        let [a, d, b, e, c, f] =
325            <[f64; 6]>::try_from(values).map_err(|_| AffineError::InvalidWorldFile)?;
326        let center = Self::new(a, b, c, d, e, f);
327        if !center.is_finite() {
328            return Err(AffineError::InvalidWorldFile);
329        }
330        Ok(center.compose(Self::translation(-0.5, -0.5)))
331    }
332
333    /// Serializes this transform in six-line World File order.
334    #[must_use]
335    pub fn to_world_file(self) -> String {
336        let center = self.compose(Self::translation(0.5, 0.5));
337        [center.a, center.d, center.b, center.e, center.c, center.f]
338            .into_iter()
339            .map(|value| format!("{value:?}"))
340            .collect::<Vec<_>>()
341            .join("\n")
342            + "\n"
343    }
344}
345
346impl Default for Affine {
347    fn default() -> Self {
348        Self::IDENTITY
349    }
350}
351
352impl Mul for Affine {
353    type Output = Self;
354
355    fn mul(self, rhs: Self) -> Self::Output {
356        self.compose(rhs)
357    }
358}
359
360impl Mul<[f64; 2]> for Affine {
361    type Output = [f64; 2];
362
363    fn mul(self, rhs: [f64; 2]) -> Self::Output {
364        self.transform_point(rhs)
365    }
366}
367
368impl fmt::Display for Affine {
369    fn fmt(&self, formatter: &mut fmt::Formatter<'_>) -> fmt::Result {
370        write!(
371            formatter,
372            "|{:.2},{:.2},{:.2}|\n|{:.2},{:.2},{:.2}|\n|0.00,0.00,1.00|",
373            self.a, self.b, self.c, self.d, self.e, self.f
374        )
375    }
376}
377
378fn cos_sin_degrees(angle: f64) -> (f64, f64) {
379    let normalized = angle.rem_euclid(360.0);
380    match normalized {
381        90.0 => (0.0, 1.0),
382        180.0 => (-1.0, 0.0),
383        270.0 => (0.0, -1.0),
384        _ => {
385            let (sine, cosine) = normalized.to_radians().sin_cos();
386            (cosine, sine)
387        }
388    }
389}