ogeom_math/
transform.rs

1//! Rigid and similarity transforms, and general affine transforms.
2//!
3//! [`Transform`] is a similarity: an orthonormal linear part, a uniform scale
4//! and a translation. That covers everything a solid modeller applies to a
5//! shape (placement, rotation, mirroring, uniform scaling) while preserving
6//! the two properties the geometry depends on: angles are unchanged, and an
7//! analytic surface stays the same *kind* of analytic surface. A cylinder
8//! remains a cylinder.
9//!
10//! [`GeneralTransform`] drops both guarantees, allowing non-uniform scaling and
11//! shear. It is a separate type on purpose: applying one turns a circle into an
12//! ellipse and a cylinder into something with no analytic form at all, so it
13//! cannot be used interchangeably.
14//!
15//! # Form classification
16//!
17//! Every [`Transform`] carries a [`TransformKind`], and applying one dispatches
18//! on it: a translation adds a vector, the identity does nothing at all. That
19//! matters because transforms are applied to every control point of every
20//! curve, every vertex of every tessellation, over an entire model: the
21//! difference between a branch and nine multiplies, repeated a hundred million
22//! times, is real.
23//!
24//! The kind is *derived from* the data rather than asserted alongside it, so it
25//! cannot drift out of agreement with the matrix it describes. Every
26//! constructor routes through one private classifier; there is no way to build a
27//! transform that claims more structure than it has.
28
29use core::ops::Mul;
30
31use ogeom_core::{OgeomResult, Tolerances, ogeom_bail};
32
33use crate::{
34    Axis, Direction, Direction2, Frame, Frame2, Matrix2, Matrix3, Point, Point2, Quaternion,
35    Vector, Vector2,
36};
37
38/// How much structure a [`Transform`] has, for dispatch.
39#[derive(Debug, Clone, Copy, PartialEq, Eq, Default)]
40pub enum TransformKind {
41    /// Does nothing.
42    #[default]
43    Identity,
44    /// Translation only.
45    Translation,
46    /// Rotation about an axis through the origin, possibly with a translation.
47    Rotation,
48    /// Reflection through a point (equivalently, a scale of `-1`).
49    PointMirror,
50    /// Reflection in a plane.
51    PlaneMirror,
52    /// Uniform scaling, possibly with a translation.
53    Scale,
54    /// Anything else: a combination with no simpler description.
55    Compound,
56}
57
58/// A similarity transform: orthonormal rotation or reflection, uniform scale,
59/// translation.
60///
61/// Applied as `p -> linear * (scale * p) + translation`.
62#[derive(Debug, Clone, Copy, PartialEq)]
63pub struct Transform {
64    linear: Matrix3,
65    scale: f64,
66    translation: Vector,
67    kind: TransformKind,
68}
69
70/// A similarity transform in the plane.
71#[derive(Debug, Clone, Copy, PartialEq)]
72pub struct Transform2 {
73    linear: Matrix2,
74    scale: f64,
75    translation: Vector2,
76    kind: TransformKind,
77}
78
79/// A general affine transform: any linear part, plus a translation.
80///
81/// Non-uniform scaling and shear are allowed, so angles are not preserved and
82/// analytic geometry does not survive intact. Kept distinct from [`Transform`]
83/// so that cannot happen by accident.
84#[derive(Debug, Clone, Copy, PartialEq)]
85pub struct GeneralTransform {
86    /// The linear part.
87    pub linear: Matrix3,
88    /// The translation.
89    pub translation: Vector,
90}
91
92/// How far from `1` a scale factor may sit and still count as unit, and how far
93/// a matrix may stray from a canonical form and still be classified as it.
94/// Dimensionless; classification is a fast-path hint, and a value that just
95/// misses the threshold is merely applied by the general path.
96const CLASSIFY_EPS: f64 = 1e-12;
97
98impl Default for Transform {
99    fn default() -> Self {
100        Self::IDENTITY
101    }
102}
103
104impl Transform {
105    /// The identity.
106    pub const IDENTITY: Self = Self {
107        linear: Matrix3::IDENTITY,
108        scale: 1.0,
109        translation: Vector::ZERO,
110        kind: TransformKind::Identity,
111    };
112
113    /// Derive the kind from the parts, and build the transform.
114    ///
115    /// The single constructor everything else routes through, so a transform's
116    /// kind can never disagree with what it actually does.
117    fn build(linear: Matrix3, scale: f64, translation: Vector) -> Self {
118        let kind = Self::classify(&linear, scale, translation);
119        Self {
120            linear,
121            scale,
122            translation,
123            kind,
124        }
125    }
126
127    /// A transform from the parts it is stored as.
128    ///
129    /// `linear` is the orthonormal part alone and `scale` the uniform factor
130    /// beside it, which is how a [`Transform`] holds them. The kind is
131    /// re-derived rather than taken on trust, since it is a function of the
132    /// other three.
133    ///
134    /// For reading a document back. Going the long way round (multiplying the
135    /// scale into the matrix and asking
136    /// [`GeneralTransform::to_similarity`](crate::GeneralTransform::to_similarity)
137    /// to factor it out again) recovers a transform that is *close*, not the
138    /// one that was written, and a round trip that drifts a little each time is
139    /// not a round trip.
140    ///
141    /// # Errors
142    ///
143    /// [`OgeomError::Construction`](ogeom_core::OgeomError::Construction) if `linear` is
144    /// not orthonormal within `eps`, or `scale` is not finite and non-zero; a
145    /// placement that squashes space is not a placement.
146    pub fn from_parts(
147        linear: Matrix3,
148        scale: f64,
149        translation: Vector,
150        eps: f64,
151    ) -> OgeomResult<Self> {
152        if !scale.is_finite() || scale == 0.0 {
153            ogeom_bail!(
154                Construction,
155                "a placement's scale must be finite and non-zero; got {scale}"
156            );
157        }
158        if !linear.is_orthonormal(eps) {
159            ogeom_bail!(
160                Construction,
161                "a placement's linear part must be orthonormal; this one shears                  or scales unevenly"
162            );
163        }
164        Ok(Self::build(linear, scale, translation))
165    }
166
167    /// Classify a similarity from its parts.
168    fn classify(linear: &Matrix3, scale: f64, translation: Vector) -> TransformKind {
169        let is_identity_linear = linear.is_equal(&Matrix3::IDENTITY, CLASSIFY_EPS);
170        let unit_scale = (scale - 1.0).abs() <= CLASSIFY_EPS;
171        let negative_unit_scale = (scale + 1.0).abs() <= CLASSIFY_EPS;
172        let no_translation = translation.square_magnitude() == 0.0;
173
174        if is_identity_linear && unit_scale {
175            return if no_translation {
176                TransformKind::Identity
177            } else {
178                TransformKind::Translation
179            };
180        }
181        if is_identity_linear && negative_unit_scale {
182            return TransformKind::PointMirror;
183        }
184        if is_identity_linear {
185            return TransformKind::Scale;
186        }
187        if unit_scale && linear.is_orthonormal(CLASSIFY_EPS) {
188            // Determinant separates a rotation from a reflection; both are
189            // orthonormal, and confusing them flips the sense of every face.
190            return if linear.determinant() > 0.0 {
191                TransformKind::Rotation
192            } else {
193                TransformKind::PlaneMirror
194            };
195        }
196        TransformKind::Compound
197    }
198
199    /// A translation.
200    #[must_use]
201    pub fn translation(v: Vector) -> Self {
202        Self::build(Matrix3::IDENTITY, 1.0, v)
203    }
204
205    /// A rotation of `angle` radians about `axis`.
206    #[must_use]
207    pub fn rotation(axis: Axis, angle: f64) -> Self {
208        let linear = Matrix3::rotation(axis.direction, angle);
209        // Rotating about an axis that misses the origin: move the axis point to
210        // the origin, rotate, move back. Written as one translation so the
211        // result stays a single transform.
212        let p = axis.location.to_vector();
213        Self::build(linear, 1.0, p - linear * p)
214    }
215
216    /// A rotation given as a quaternion, about an axis through the origin.
217    #[must_use]
218    pub fn from_quaternion(q: Quaternion) -> Self {
219        Self::build(q.to_matrix(), 1.0, Vector::ZERO)
220    }
221
222    /// A uniform scaling about `centre`.
223    ///
224    /// # Errors
225    ///
226    /// [`OgeomError::Construction`](ogeom_core::OgeomError::Construction) if `factor` is
227    /// zero or non-finite. A zero scale collapses every shape to a point and is
228    /// not invertible, so it is refused rather than allowed to produce
229    /// degenerate geometry.
230    pub fn scaling(centre: Point, factor: f64, tol: Tolerances) -> OgeomResult<Self> {
231        if !factor.is_finite() || factor.abs() <= tol.confusion() {
232            ogeom_bail!(Construction, "scale factor {factor} is degenerate");
233        }
234        let c = centre.to_vector();
235        Ok(Self::build(Matrix3::IDENTITY, factor, c - c * factor))
236    }
237
238    /// Reflection through a point.
239    #[must_use]
240    pub fn point_mirror(centre: Point) -> Self {
241        let c = centre.to_vector();
242        Self::build(Matrix3::IDENTITY, -1.0, c + c)
243    }
244
245    /// Reflection in the plane through `origin` with the given `normal`.
246    #[must_use]
247    pub fn plane_mirror(origin: Point, normal: Direction) -> Self {
248        let linear = Matrix3::reflection(normal);
249        let p = origin.to_vector();
250        Self::build(linear, 1.0, p - linear * p)
251    }
252
253    /// Reflection in a line: a half turn about it.
254    #[must_use]
255    pub fn axis_mirror(axis: Axis) -> Self {
256        Self::rotation(axis, core::f64::consts::PI)
257    }
258
259    /// The transform taking world coordinates into `frame`'s local coordinates.
260    #[must_use]
261    pub fn to_frame(frame: &Frame) -> Self {
262        let linear = frame.to_matrix().transposed();
263        Self::build(linear, 1.0, -(linear * frame.origin().to_vector()))
264    }
265
266    /// The transform taking `frame`'s local coordinates into world coordinates.
267    #[must_use]
268    pub fn from_frame(frame: &Frame) -> Self {
269        Self::build(frame.to_matrix(), 1.0, frame.origin().to_vector())
270    }
271
272    /// The transform taking `from`'s local coordinates into `to`'s.
273    #[must_use]
274    pub fn between_frames(from: &Frame, to: &Frame) -> Self {
275        Self::to_frame(to) * Self::from_frame(from)
276    }
277
278    /// This transform's classification.
279    #[must_use]
280    pub const fn kind(&self) -> TransformKind {
281        self.kind
282    }
283
284    /// The orthonormal part.
285    #[must_use]
286    pub const fn linear(&self) -> Matrix3 {
287        self.linear
288    }
289
290    /// The uniform scale factor. Negative for a point mirror.
291    #[must_use]
292    pub const fn scale_factor(&self) -> f64 {
293        self.scale
294    }
295
296    /// The translation.
297    #[must_use]
298    pub const fn translation_vector(&self) -> Vector {
299        self.translation
300    }
301
302    /// Whether this transform preserves handedness.
303    ///
304    /// A shape transformed by a transform that does not must have its
305    /// orientation flipped to stay consistent; otherwise a mirrored solid ends
306    /// up inside out.
307    #[must_use]
308    pub fn preserves_handedness(&self) -> bool {
309        self.linear.determinant() * self.scale.signum() > 0.0
310    }
311
312    /// Apply to a point.
313    #[must_use]
314    pub fn apply(&self, p: Point) -> Point {
315        match self.kind {
316            TransformKind::Identity => p,
317            TransformKind::Translation => p + self.translation,
318            TransformKind::PointMirror | TransformKind::Scale => {
319                Point::from_vector(p.to_vector() * self.scale + self.translation)
320            }
321            TransformKind::Rotation | TransformKind::PlaneMirror => {
322                Point::from_vector(self.linear * p.to_vector() + self.translation)
323            }
324            TransformKind::Compound => {
325                Point::from_vector(self.linear * (p.to_vector() * self.scale) + self.translation)
326            }
327        }
328    }
329
330    /// Apply to a free vector. Translation does not affect it.
331    #[must_use]
332    pub fn apply_vector(&self, v: Vector) -> Vector {
333        match self.kind {
334            TransformKind::Identity | TransformKind::Translation => v,
335            TransformKind::PointMirror | TransformKind::Scale => v * self.scale,
336            TransformKind::Rotation | TransformKind::PlaneMirror => self.linear * v,
337            TransformKind::Compound => self.linear * (v * self.scale),
338        }
339    }
340
341    /// Apply to a direction, renormalizing.
342    ///
343    /// A similarity maps unit vectors to vectors of length `|scale|`, so the
344    /// result is rescaled. Under a negative scale the direction reverses, which
345    /// is the correct behaviour and not a sign error.
346    ///
347    /// # Errors
348    ///
349    /// [`OgeomError::Construction`](ogeom_core::OgeomError::Construction) if the result
350    /// cannot be normalized, which a valid similarity never produces.
351    pub fn apply_direction(&self, d: Direction, tol: Tolerances) -> OgeomResult<Direction> {
352        match self.kind {
353            TransformKind::Identity | TransformKind::Translation | TransformKind::Scale => Ok(d),
354            TransformKind::PointMirror => Ok(d.reversed()),
355            _ => Direction::new(self.apply_vector(d.vector()), tol),
356        }
357    }
358
359    /// Apply to a frame, transforming origin and all three axes.
360    ///
361    /// # Errors
362    ///
363    /// [`OgeomError::Construction`](ogeom_core::OgeomError::Construction) if the
364    /// transformed axes cannot be renormalized.
365    pub fn apply_frame(&self, f: &Frame, tol: Tolerances) -> OgeomResult<Frame> {
366        Frame::from_axes(
367            self.apply(f.origin()),
368            self.apply_direction(f.x(), tol)?,
369            self.apply_direction(f.y(), tol)?,
370            self.apply_direction(f.z(), tol)?,
371            tol,
372        )
373    }
374
375    /// The inverse.
376    ///
377    /// # Errors
378    ///
379    /// [`OgeomError::Numeric`](ogeom_core::OgeomError::Numeric) if the linear part is
380    /// singular, which a valid similarity never is.
381    pub fn inverse(&self) -> OgeomResult<Self> {
382        if self.kind == TransformKind::Identity {
383            return Ok(Self::IDENTITY);
384        }
385        if self.kind == TransformKind::Translation {
386            return Ok(Self::translation(-self.translation));
387        }
388        if self.scale == 0.0 {
389            ogeom_bail!(Numeric, "transform has a zero scale and no inverse");
390        }
391        // The linear part is orthonormal, so its inverse is its transpose; no
392        // need to go through a general inversion, and no rounding beyond the
393        // transpose itself.
394        let inv_linear = self.linear.transposed();
395        let inv_scale = 1.0 / self.scale;
396        Ok(Self::build(
397            inv_linear,
398            inv_scale,
399            -(inv_linear * self.translation) * inv_scale,
400        ))
401    }
402
403    /// Whether two transforms agree in effect.
404    #[must_use]
405    pub fn is_equal(&self, other: &Self, tol: Tolerances) -> bool {
406        (self.scale - other.scale).abs() <= CLASSIFY_EPS
407            && self.linear.is_equal(&other.linear, CLASSIFY_EPS)
408            && self.translation.is_equal(other.translation, tol)
409    }
410
411    /// This transform as a general affine one.
412    #[must_use]
413    pub fn to_general(&self) -> GeneralTransform {
414        GeneralTransform {
415            linear: self.linear * self.scale,
416            translation: self.translation,
417        }
418    }
419}
420
421impl Mul for Transform {
422    type Output = Self;
423    /// Composition. `(a * b)` applies `b` first, then `a`.
424    fn mul(self, b: Self) -> Self {
425        if self.kind == TransformKind::Identity {
426            return b;
427        }
428        if b.kind == TransformKind::Identity {
429            return self;
430        }
431        Self::build(
432            self.linear * b.linear,
433            self.scale * b.scale,
434            self.linear * (b.translation * self.scale) + self.translation,
435        )
436    }
437}
438
439impl Default for Transform2 {
440    fn default() -> Self {
441        Self::IDENTITY
442    }
443}
444
445impl Transform2 {
446    /// The identity.
447    pub const IDENTITY: Self = Self {
448        linear: Matrix2::IDENTITY,
449        scale: 1.0,
450        translation: Vector2::ZERO,
451        kind: TransformKind::Identity,
452    };
453
454    fn build(linear: Matrix2, scale: f64, translation: Vector2) -> Self {
455        let is_identity_linear = linear.is_equal(&Matrix2::IDENTITY, CLASSIFY_EPS);
456        let unit_scale = (scale - 1.0).abs() <= CLASSIFY_EPS;
457        let kind = if is_identity_linear && unit_scale {
458            if translation.square_magnitude() == 0.0 {
459                TransformKind::Identity
460            } else {
461                TransformKind::Translation
462            }
463        } else if is_identity_linear && (scale + 1.0).abs() <= CLASSIFY_EPS {
464            TransformKind::PointMirror
465        } else if is_identity_linear {
466            TransformKind::Scale
467        } else if unit_scale && (linear.determinant().abs() - 1.0).abs() <= CLASSIFY_EPS {
468            if linear.determinant() > 0.0 {
469                TransformKind::Rotation
470            } else {
471                TransformKind::PlaneMirror
472            }
473        } else {
474            TransformKind::Compound
475        };
476        Self {
477            linear,
478            scale,
479            translation,
480            kind,
481        }
482    }
483
484    /// A translation.
485    #[must_use]
486    pub fn translation(v: Vector2) -> Self {
487        Self::build(Matrix2::IDENTITY, 1.0, v)
488    }
489
490    /// A rotation about `centre`.
491    #[must_use]
492    pub fn rotation(centre: Point2, angle: f64) -> Self {
493        let linear = Matrix2::rotation(angle);
494        let c = centre.to_vector();
495        Self::build(linear, 1.0, c - linear * c)
496    }
497
498    /// A uniform scaling about `centre`.
499    ///
500    /// # Errors
501    ///
502    /// [`OgeomError::Construction`](ogeom_core::OgeomError::Construction) if `factor` is
503    /// degenerate.
504    pub fn scaling(centre: Point2, factor: f64, tol: Tolerances) -> OgeomResult<Self> {
505        if !factor.is_finite() || factor.abs() <= tol.confusion() {
506            ogeom_bail!(Construction, "scale factor {factor} is degenerate");
507        }
508        let c = centre.to_vector();
509        Ok(Self::build(Matrix2::IDENTITY, factor, c - c * factor))
510    }
511
512    /// Reflection in the line through `origin` with the given `normal`.
513    #[must_use]
514    pub fn line_mirror(origin: Point2, normal: Direction2) -> Self {
515        let (x, y) = (normal.x(), normal.y());
516        let linear = Matrix2::new([
517            [(-2.0f64).mul_add(x * x, 1.0), -2.0 * x * y],
518            [-2.0 * x * y, (-2.0f64).mul_add(y * y, 1.0)],
519        ]);
520        let p = origin.to_vector();
521        Self::build(linear, 1.0, p - linear * p)
522    }
523
524    /// This transform's classification.
525    #[must_use]
526    pub const fn kind(&self) -> TransformKind {
527        self.kind
528    }
529
530    /// The orthonormal part.
531    #[must_use]
532    pub const fn linear(&self) -> Matrix2 {
533        self.linear
534    }
535
536    /// The uniform scale factor. Negative for a point mirror.
537    #[must_use]
538    pub const fn scale_factor(&self) -> f64 {
539        self.scale
540    }
541
542    /// The translation.
543    #[must_use]
544    pub const fn translation_vector(&self) -> Vector2 {
545        self.translation
546    }
547
548    /// Whether this transform preserves handedness.
549    #[must_use]
550    pub fn preserves_handedness(&self) -> bool {
551        self.linear.determinant() * self.scale.signum() > 0.0
552    }
553
554    /// Apply to a point.
555    #[must_use]
556    pub fn apply(&self, p: Point2) -> Point2 {
557        match self.kind {
558            TransformKind::Identity => p,
559            TransformKind::Translation => p + self.translation,
560            TransformKind::PointMirror | TransformKind::Scale => {
561                Point2::from_vector(p.to_vector() * self.scale + self.translation)
562            }
563            TransformKind::Rotation | TransformKind::PlaneMirror => {
564                Point2::from_vector(self.linear * p.to_vector() + self.translation)
565            }
566            TransformKind::Compound => {
567                Point2::from_vector(self.linear * (p.to_vector() * self.scale) + self.translation)
568            }
569        }
570    }
571
572    /// Apply to a direction, renormalizing.
573    ///
574    /// # Errors
575    ///
576    /// [`OgeomError::Construction`](ogeom_core::OgeomError::Construction) if the result
577    /// cannot be normalized, which a valid similarity never produces.
578    pub fn apply_direction(&self, d: Direction2, tol: Tolerances) -> OgeomResult<Direction2> {
579        match self.kind {
580            TransformKind::Identity | TransformKind::Translation | TransformKind::Scale => Ok(d),
581            TransformKind::PointMirror => Ok(d.reversed()),
582            _ => Direction2::new(self.apply_vector(d.vector()), tol),
583        }
584    }
585
586    /// Apply to a frame, transforming the origin and both axes.
587    ///
588    /// # Errors
589    ///
590    /// [`OgeomError::Construction`](ogeom_core::OgeomError::Construction) if the
591    /// transformed axes cannot be renormalized or are no longer perpendicular.
592    pub fn apply_frame(&self, f: &Frame2, tol: Tolerances) -> OgeomResult<Frame2> {
593        Frame2::from_axes(
594            self.apply(f.origin()),
595            self.apply_direction(f.x(), tol)?,
596            self.apply_direction(f.y(), tol)?,
597            tol,
598        )
599    }
600
601    /// Apply to a free vector.
602    #[must_use]
603    pub fn apply_vector(&self, v: Vector2) -> Vector2 {
604        match self.kind {
605            TransformKind::Identity | TransformKind::Translation => v,
606            TransformKind::PointMirror | TransformKind::Scale => v * self.scale,
607            TransformKind::Rotation | TransformKind::PlaneMirror => self.linear * v,
608            TransformKind::Compound => self.linear * (v * self.scale),
609        }
610    }
611
612    /// The inverse.
613    ///
614    /// # Errors
615    ///
616    /// [`OgeomError::Numeric`](ogeom_core::OgeomError::Numeric) if the transform is
617    /// degenerate.
618    pub fn inverse(&self) -> OgeomResult<Self> {
619        if self.kind == TransformKind::Identity {
620            return Ok(Self::IDENTITY);
621        }
622        if self.scale == 0.0 {
623            ogeom_bail!(Numeric, "transform has a zero scale and no inverse");
624        }
625        let inv_linear = self.linear.transposed();
626        let inv_scale = 1.0 / self.scale;
627        Ok(Self::build(
628            inv_linear,
629            inv_scale,
630            -(inv_linear * self.translation) * inv_scale,
631        ))
632    }
633
634    /// Whether two transforms agree in effect.
635    #[must_use]
636    pub fn is_equal(&self, other: &Self, tol: Tolerances) -> bool {
637        (self.scale - other.scale).abs() <= CLASSIFY_EPS
638            && self.linear.is_equal(&other.linear, CLASSIFY_EPS)
639            && self.translation.is_equal(other.translation, tol)
640    }
641}
642
643impl Mul for Transform2 {
644    type Output = Self;
645    fn mul(self, b: Self) -> Self {
646        if self.kind == TransformKind::Identity {
647            return b;
648        }
649        if b.kind == TransformKind::Identity {
650            return self;
651        }
652        Self::build(
653            self.linear * b.linear,
654            self.scale * b.scale,
655            self.linear * (b.translation * self.scale) + self.translation,
656        )
657    }
658}
659
660impl Default for GeneralTransform {
661    fn default() -> Self {
662        Self::IDENTITY
663    }
664}
665
666impl GeneralTransform {
667    /// The identity.
668    pub const IDENTITY: Self = Self {
669        linear: Matrix3::IDENTITY,
670        translation: Vector::ZERO,
671    };
672
673    /// From a linear part and a translation.
674    #[must_use]
675    pub const fn new(linear: Matrix3, translation: Vector) -> Self {
676        Self {
677            linear,
678            translation,
679        }
680    }
681
682    /// Non-uniform scaling about the origin.
683    #[must_use]
684    pub const fn scaling_xyz(x: f64, y: f64, z: f64) -> Self {
685        Self::new(Matrix3::scaling_xyz(x, y, z), Vector::ZERO)
686    }
687
688    /// Apply to a point.
689    #[must_use]
690    pub fn apply(&self, p: Point) -> Point {
691        Point::from_vector(self.linear * p.to_vector() + self.translation)
692    }
693
694    /// Apply to a free vector.
695    #[must_use]
696    pub fn apply_vector(&self, v: Vector) -> Vector {
697        self.linear * v
698    }
699
700    /// Apply to a normal vector.
701    ///
702    /// Normals transform by the inverse transpose, not by the linear part
703    /// itself. Using the linear part directly is only correct for a similarity;
704    /// under any shear or non-uniform scale it tilts normals off the surface
705    /// they belong to, which then breaks every orientation test downstream.
706    ///
707    /// The result is not renormalized; it is a direction, not a length.
708    ///
709    /// # Errors
710    ///
711    /// [`OgeomError::Numeric`](ogeom_core::OgeomError::Numeric) if the linear part is
712    /// singular.
713    pub fn apply_normal(&self, n: Vector) -> OgeomResult<Vector> {
714        Ok(self.linear.inverse()?.transposed() * n)
715    }
716
717    /// Whether this transform preserves handedness.
718    #[must_use]
719    pub fn preserves_handedness(&self) -> bool {
720        self.linear.determinant() > 0.0
721    }
722
723    /// The factor by which volumes are multiplied. Negative if handedness flips.
724    #[must_use]
725    pub fn volume_ratio(&self) -> f64 {
726        self.linear.determinant()
727    }
728
729    /// The inverse.
730    ///
731    /// # Errors
732    ///
733    /// [`OgeomError::Numeric`](ogeom_core::OgeomError::Numeric) if the linear part is
734    /// singular.
735    pub fn inverse(&self) -> OgeomResult<Self> {
736        let inv = self.linear.inverse()?;
737        Ok(Self::new(inv, -(inv * self.translation)))
738    }
739
740    /// Whether this is a similarity, and so can be narrowed to a [`Transform`].
741    #[must_use]
742    pub fn is_similarity(&self, eps: f64) -> bool {
743        self.to_similarity(eps).is_some()
744    }
745
746    /// This transform as a [`Transform`], if it is in fact a similarity.
747    ///
748    /// Returns `None` when the linear part contains shear or non-uniform
749    /// scaling, since no similarity describes it.
750    #[must_use]
751    pub fn to_similarity(&self, eps: f64) -> Option<Transform> {
752        // A similarity's linear part is `s * R` with `R` orthonormal, so its
753        // columns are mutually orthogonal and all of length |s|.
754        let det = self.linear.determinant();
755        if det == 0.0 {
756            return None;
757        }
758        let scale = det.abs().cbrt() * det.signum();
759        let rotation = self.linear * (1.0 / scale);
760        if !rotation.is_orthonormal(eps) {
761            return None;
762        }
763        Some(Transform::build(rotation, scale, self.translation))
764    }
765}
766
767impl Mul for GeneralTransform {
768    type Output = Self;
769    /// Composition. `(a * b)` applies `b` first, then `a`.
770    fn mul(self, b: Self) -> Self {
771        Self::new(
772            self.linear * b.linear,
773            self.linear * b.translation + self.translation,
774        )
775    }
776}
777
778impl From<Transform> for GeneralTransform {
779    fn from(t: Transform) -> Self {
780        t.to_general()
781    }
782}
783
784#[cfg(test)]
785#[allow(clippy::unwrap_used)]
786mod tests {
787    use super::*;
788    use approx::assert_relative_eq;
789
790    const T: Tolerances = Tolerances::millimetres();
791
792    fn sample_points() -> [Point; 4] {
793        [
794            Point::ORIGIN,
795            Point::new(1.0, 0.0, 0.0),
796            Point::new(-3.0, 7.5, 2.25),
797            Point::new(1e3, -1e3, 0.5),
798        ]
799    }
800
801    #[test]
802    fn classification_matches_what_the_transform_does() {
803        assert_eq!(Transform::IDENTITY.kind(), TransformKind::Identity);
804        assert_eq!(
805            Transform::translation(Vector::X).kind(),
806            TransformKind::Translation
807        );
808        assert_eq!(
809            Transform::rotation(Axis::Z, 0.5).kind(),
810            TransformKind::Rotation
811        );
812        assert_eq!(
813            Transform::point_mirror(Point::ORIGIN).kind(),
814            TransformKind::PointMirror
815        );
816        assert_eq!(
817            Transform::plane_mirror(Point::ORIGIN, Direction::Z).kind(),
818            TransformKind::PlaneMirror
819        );
820        assert_eq!(
821            Transform::scaling(Point::ORIGIN, 3.0, T).unwrap().kind(),
822            TransformKind::Scale
823        );
824        let compound =
825            Transform::rotation(Axis::Z, 0.5) * Transform::scaling(Point::ORIGIN, 3.0, T).unwrap();
826        assert_eq!(compound.kind(), TransformKind::Compound);
827    }
828
829    #[test]
830    fn a_zero_rotation_classifies_as_identity_not_rotation() {
831        // The classification is derived from the data, so it cannot claim more
832        // structure than the transform has, or less.
833        assert_eq!(
834            Transform::rotation(Axis::Z, 0.0).kind(),
835            TransformKind::Identity
836        );
837        assert_eq!(
838            Transform::translation(Vector::ZERO).kind(),
839            TransformKind::Identity
840        );
841        assert_eq!(
842            Transform::scaling(Point::ORIGIN, 1.0, T).unwrap().kind(),
843            TransformKind::Identity
844        );
845    }
846
847    #[test]
848    fn every_dispatch_path_gives_the_same_answer_as_the_general_one() {
849        // The whole point of classification is speed, so each fast path must
850        // agree exactly with the general formula it replaces.
851        let cases = [
852            Transform::IDENTITY,
853            Transform::translation(Vector::new(1.0, -2.0, 3.0)),
854            Transform::rotation(Axis::new(Point::new(1.0, 0.0, 0.0), Direction::Z), 0.7),
855            Transform::point_mirror(Point::new(2.0, 0.0, -1.0)),
856            Transform::plane_mirror(Point::new(0.0, 1.0, 0.0), Direction::Y),
857            Transform::scaling(Point::new(1.0, 1.0, 1.0), 2.5, T).unwrap(),
858        ];
859        for t in cases {
860            for p in sample_points() {
861                let general = Point::from_vector(
862                    t.linear() * (p.to_vector() * t.scale_factor()) + t.translation_vector(),
863                );
864                assert!(
865                    t.apply(p).is_equal(general, T),
866                    "fast path for {:?} disagrees",
867                    t.kind()
868                );
869            }
870        }
871    }
872
873    #[test]
874    fn rotation_about_an_off_origin_axis_leaves_the_axis_fixed() {
875        let axis = Axis::new(Point::new(5.0, 3.0, 0.0), Direction::Z);
876        let t = Transform::rotation(axis, 1.234);
877        assert!(t.apply(axis.location).is_equal(axis.location, T));
878        assert!(
879            t.apply(axis.point_at(10.0))
880                .is_equal(axis.point_at(10.0), T)
881        );
882        // A point off the axis moves, staying at the same radius.
883        let p = Point::new(6.0, 3.0, 0.0);
884        assert_relative_eq!(axis.distance_to(t.apply(p)), 1.0, epsilon = 1e-14);
885    }
886
887    #[test]
888    fn scaling_about_a_centre_leaves_the_centre_fixed() {
889        let c = Point::new(3.0, -1.0, 2.0);
890        let t = Transform::scaling(c, 4.0, T).unwrap();
891        assert!(t.apply(c).is_equal(c, T));
892        let p = c + Vector::new(1.0, 0.0, 0.0);
893        assert!(t.apply(p).is_equal(c + Vector::new(4.0, 0.0, 0.0), T));
894    }
895
896    #[test]
897    fn degenerate_scales_are_refused() {
898        assert!(Transform::scaling(Point::ORIGIN, 0.0, T).is_err());
899        assert!(Transform::scaling(Point::ORIGIN, f64::NAN, T).is_err());
900        assert!(Transform::scaling(Point::ORIGIN, f64::INFINITY, T).is_err());
901        assert!(
902            Transform::scaling(Point::ORIGIN, -2.0, T).is_ok(),
903            "negative is fine"
904        );
905    }
906
907    #[test]
908    fn handedness_tracks_mirroring() {
909        assert!(Transform::rotation(Axis::Z, 1.0).preserves_handedness());
910        assert!(Transform::translation(Vector::X).preserves_handedness());
911        assert!(
912            Transform::scaling(Point::ORIGIN, 3.0, T)
913                .unwrap()
914                .preserves_handedness()
915        );
916        assert!(!Transform::plane_mirror(Point::ORIGIN, Direction::Z).preserves_handedness());
917        assert!(!Transform::point_mirror(Point::ORIGIN).preserves_handedness());
918        // Two mirrors make a rotation.
919        let twice = Transform::plane_mirror(Point::ORIGIN, Direction::Z)
920            * Transform::plane_mirror(Point::ORIGIN, Direction::X);
921        assert!(twice.preserves_handedness());
922    }
923
924    #[test]
925    fn axis_mirror_is_a_half_turn() {
926        let t = Transform::axis_mirror(Axis::Z);
927        assert!(
928            t.apply(Point::new(1.0, 0.0, 5.0))
929                .is_equal(Point::new(-1.0, 0.0, 5.0), T)
930        );
931        assert!(t.preserves_handedness(), "a half turn is a rotation");
932    }
933
934    #[test]
935    fn inverse_round_trips_for_every_kind() {
936        let cases = [
937            Transform::IDENTITY,
938            Transform::translation(Vector::new(1.0, -2.0, 3.0)),
939            Transform::rotation(Axis::new(Point::new(1.0, 2.0, 3.0), Direction::Y), 2.1),
940            Transform::point_mirror(Point::new(1.0, 1.0, 1.0)),
941            Transform::plane_mirror(Point::new(0.0, 0.0, 4.0), Direction::Z),
942            Transform::scaling(Point::new(-1.0, 0.0, 0.0), 0.25, T).unwrap(),
943        ];
944        for t in cases {
945            let inv = t.inverse().unwrap();
946            for p in sample_points() {
947                assert!(inv.apply(t.apply(p)).is_equal(p, T), "{:?}", t.kind());
948                assert!(t.apply(inv.apply(p)).is_equal(p, T), "{:?}", t.kind());
949            }
950        }
951    }
952
953    #[test]
954    fn composition_applies_right_to_left() {
955        let a = Transform::translation(Vector::new(10.0, 0.0, 0.0));
956        let b = Transform::rotation(Axis::Z, core::f64::consts::FRAC_PI_2);
957        let p = Point::new(1.0, 0.0, 0.0);
958        assert!((a * b).apply(p).is_equal(a.apply(b.apply(p)), T));
959        assert!((b * a).apply(p).is_equal(b.apply(a.apply(p)), T));
960        // And the two orders genuinely differ.
961        assert!(!(a * b).is_equal(&(b * a), T));
962    }
963
964    #[test]
965    fn composition_is_associative() {
966        let a = Transform::rotation(Axis::X, 0.3);
967        let b = Transform::scaling(Point::new(1.0, 0.0, 0.0), 2.0, T).unwrap();
968        let c = Transform::translation(Vector::new(0.0, 5.0, 0.0));
969        assert!(((a * b) * c).is_equal(&(a * (b * c)), T));
970    }
971
972    #[test]
973    fn vectors_ignore_translation_and_directions_stay_unit() {
974        let t = Transform::translation(Vector::new(100.0, 0.0, 0.0))
975            * Transform::rotation(Axis::Z, 0.9);
976        let v = Vector::new(1.0, 2.0, 3.0);
977        assert!(
978            t.apply_vector(v)
979                .is_equal(Transform::rotation(Axis::Z, 0.9).apply_vector(v), T)
980        );
981        let d = t.apply_direction(Direction::X, T).unwrap();
982        assert_relative_eq!(d.vector().magnitude(), 1.0, epsilon = 1e-15);
983    }
984
985    #[test]
986    fn a_point_mirror_reverses_directions() {
987        let t = Transform::point_mirror(Point::new(5.0, 5.0, 5.0));
988        assert!(
989            t.apply_direction(Direction::X, T)
990                .unwrap()
991                .is_equal(-Direction::X, T)
992        );
993        // A positive scale does not.
994        let s = Transform::scaling(Point::ORIGIN, 3.0, T).unwrap();
995        assert!(
996            s.apply_direction(Direction::X, T)
997                .unwrap()
998                .is_equal(Direction::X, T)
999        );
1000    }
1001
1002    #[test]
1003    fn frame_transforms_round_trip_through_world() {
1004        let f = Frame::new(
1005            Point::new(1.0, 2.0, 3.0),
1006            Direction::from_coords(1.0, 1.0, 0.0, T).unwrap(),
1007            Direction::Z,
1008            T,
1009        )
1010        .unwrap();
1011        let to = Transform::to_frame(&f);
1012        let from = Transform::from_frame(&f);
1013        for p in sample_points() {
1014            assert!(from.apply(to.apply(p)).is_equal(p, T));
1015            // And they agree with the frame's own conversion.
1016            assert!(to.apply(p).is_equal(f.to_local(p), T));
1017            assert!(from.apply(f.to_local(p)).is_equal(p, T));
1018        }
1019    }
1020
1021    #[test]
1022    fn between_frames_composes_correctly() {
1023        let a = Frame::new(Point::new(1.0, 0.0, 0.0), Direction::Z, Direction::X, T).unwrap();
1024        let b = Frame::new(Point::new(0.0, 5.0, 0.0), Direction::X, Direction::Y, T).unwrap();
1025        let t = Transform::between_frames(&a, &b);
1026        // A point at local (1,2,3) in `a` must land at the same world position
1027        // when read back out of `b`.
1028        let local = Point::new(1.0, 2.0, 3.0);
1029        assert!(b.to_world(t.apply(local)).is_equal(a.to_world(local), T));
1030    }
1031
1032    #[test]
1033    fn general_transform_normals_use_the_inverse_transpose() {
1034        // Non-uniform scaling: the plane z = x has normal (1, 0, -1) up to
1035        // scale. Scale x by 2 and the plane becomes z = x/2, whose normal is
1036        // (1, 0, -2) up to scale, not (2, 0, -1), which is what applying the
1037        // linear part directly would give.
1038        let g = GeneralTransform::scaling_xyz(2.0, 1.0, 1.0);
1039        let n = Vector::new(1.0, 0.0, -1.0);
1040        let transformed = g.apply_normal(n).unwrap();
1041
1042        let on_plane = Vector::new(1.0, 0.0, 1.0);
1043        assert_relative_eq!(n.dot(on_plane), 0.0, epsilon = 1e-15);
1044        assert_relative_eq!(
1045            transformed.dot(g.apply_vector(on_plane)),
1046            0.0,
1047            epsilon = 1e-14,
1048            max_relative = 1e-14
1049        );
1050        // The naive answer does not stay perpendicular.
1051        assert!(g.apply_vector(n).dot(g.apply_vector(on_plane)).abs() > 1e-6);
1052    }
1053
1054    #[test]
1055    fn general_transform_recognizes_similarities() {
1056        let similar: GeneralTransform = Transform::rotation(Axis::Z, 0.4).into();
1057        assert!(similar.is_similarity(1e-12));
1058        let narrowed = similar.to_similarity(1e-12).unwrap();
1059        assert_eq!(narrowed.kind(), TransformKind::Rotation);
1060
1061        let scaled: GeneralTransform = Transform::scaling(Point::ORIGIN, 3.0, T).unwrap().into();
1062        assert!(scaled.is_similarity(1e-12));
1063        assert_relative_eq!(
1064            scaled.to_similarity(1e-12).unwrap().scale_factor(),
1065            3.0,
1066            epsilon = 1e-12
1067        );
1068
1069        assert!(!GeneralTransform::scaling_xyz(1.0, 2.0, 3.0).is_similarity(1e-12));
1070        let shear = GeneralTransform::new(
1071            Matrix3::new([[1.0, 0.5, 0.0], [0.0, 1.0, 0.0], [0.0, 0.0, 1.0]]),
1072            Vector::ZERO,
1073        );
1074        assert!(!shear.is_similarity(1e-12));
1075    }
1076
1077    #[test]
1078    fn general_transform_volume_ratio_and_inverse() {
1079        let g = GeneralTransform::scaling_xyz(2.0, 3.0, 4.0);
1080        assert_relative_eq!(g.volume_ratio(), 24.0);
1081        assert!(g.preserves_handedness());
1082        let inv = g.inverse().unwrap();
1083        for p in sample_points() {
1084            assert!(inv.apply(g.apply(p)).is_equal(p, T));
1085        }
1086        let flip = GeneralTransform::scaling_xyz(-1.0, 1.0, 1.0);
1087        assert!(!flip.preserves_handedness());
1088        assert!(
1089            GeneralTransform::scaling_xyz(0.0, 1.0, 1.0)
1090                .inverse()
1091                .is_err()
1092        );
1093    }
1094
1095    #[test]
1096    fn transform2_behaves_like_its_3d_counterpart() {
1097        let r = Transform2::rotation(Point2::new(1.0, 1.0), core::f64::consts::FRAC_PI_2);
1098        assert_eq!(r.kind(), TransformKind::Rotation);
1099        assert!(
1100            r.apply(Point2::new(1.0, 1.0))
1101                .is_equal(Point2::new(1.0, 1.0), T)
1102        );
1103        assert!(
1104            r.apply(Point2::new(2.0, 1.0))
1105                .is_equal(Point2::new(1.0, 2.0), T)
1106        );
1107        assert!(
1108            r.inverse()
1109                .unwrap()
1110                .apply(r.apply(Point2::ORIGIN))
1111                .is_equal(Point2::ORIGIN, T)
1112        );
1113
1114        let m = Transform2::line_mirror(Point2::ORIGIN, Direction2::Y);
1115        assert_eq!(m.kind(), TransformKind::PlaneMirror);
1116        assert!(!m.preserves_handedness());
1117        assert!(
1118            m.apply(Point2::new(3.0, 2.0))
1119                .is_equal(Point2::new(3.0, -2.0), T)
1120        );
1121
1122        assert!(Transform2::scaling(Point2::ORIGIN, 0.0, T).is_err());
1123    }
1124}
ogeom_math/transform.rs

ogeom_math/
transform.rs
