pdbtbx/
transformation.rs

1#![allow(dead_code)]
2
3/// A 3D affine transformation matrix
4#[cfg_attr(feature = "serde", derive(serde::Serialize, serde::Deserialize))]
5#[derive(Debug, Clone, PartialEq)]
6pub struct TransformationMatrix {
7    /// The matrix itself
8    matrix: [[f64; 4]; 3],
9}
10
11impl TransformationMatrix {
12    /// Get the raw matrix (row major order)
13    #[must_use]
14    pub const fn matrix(&self) -> [[f64; 4]; 3] {
15        self.matrix
16    }
17
18    /// Get the raw matrix (row major order)
19    #[must_use]
20    pub fn matrix_mut(&mut self) -> &mut [[f64; 4]; 3] {
21        &mut self.matrix
22    }
23
24    /// Set the raw matrix (row major order), the user needs to make sure the matrix is valid
25    pub fn set_matrix(&mut self, new_matrix: [[f64; 4]; 3]) {
26        self.matrix = new_matrix;
27    }
28
29    /// Create a matrix defining identity, so no transformation
30    #[must_use]
31    pub const fn identity() -> Self {
32        Self {
33            matrix: [
34                [1.0, 0.0, 0.0, 0.0],
35                [0.0, 1.0, 0.0, 0.0],
36                [0.0, 0.0, 1.0, 0.0],
37            ],
38        }
39    }
40
41    /// Create a matrix with the given matrix
42    #[must_use]
43    pub const fn from_matrix(matrix: [[f64; 4]; 3]) -> Self {
44        Self { matrix }
45    }
46
47    /// Create a matrix defining a rotation around the X axis
48    /// ## Arguments
49    /// * `deg` the rotation in degrees
50    /// ## Panics
51    /// It panics if `deg` is not finite (`f64.is_finite()`)
52    #[must_use]
53    pub fn rotation_x(deg: f64) -> Self {
54        assert!(deg.is_finite(), "The amount of degrees is not finite");
55        let (s, c) = deg.to_radians().sin_cos();
56        Self {
57            matrix: [[1.0, 0.0, 0.0, 0.0], [0.0, c, -s, 0.0], [0.0, s, c, 0.0]],
58        }
59    }
60
61    /// Create a matrix defining a rotation around the Y axis
62    /// ## Arguments
63    /// * `deg` the rotation in degrees
64    /// ## Panics
65    /// It panics if `deg` is not finite (`f64.is_finite()`)
66    #[must_use]
67    pub fn rotation_y(deg: f64) -> Self {
68        assert!(deg.is_finite(), "The amount of degrees is not finite");
69        let (s, c) = deg.to_radians().sin_cos();
70        Self {
71            matrix: [[c, 0.0, s, 0.0], [0.0, 1.0, 0.0, 0.0], [-s, 0.0, c, 0.0]],
72        }
73    }
74
75    /// Create a matrix defining a rotation around the Z axis
76    /// ## Arguments
77    /// * `deg` the rotation in degrees
78    /// ## Panics
79    /// It panics if `deg` is not finite (`f64.is_finite()`)
80    #[must_use]
81    pub fn rotation_z(deg: f64) -> Self {
82        assert!(deg.is_finite(), "The amount of degrees is not finite");
83        let c = deg.to_radians().cos();
84        let s = deg.to_radians().sin();
85        Self {
86            matrix: [[c, -s, 0.0, 0.0], [s, c, 0.0, 0.0], [0.0, 0.0, 1.0, 0.0]],
87        }
88    }
89
90    /// Create a matrix defining a translation
91    /// ## Panics
92    /// It panics if any of the arguments is not finite (`f64.is_finite()`)
93    #[must_use]
94    pub fn translation(x: f64, y: f64, z: f64) -> Self {
95        assert!(
96            x.is_finite() && y.is_finite() && z.is_finite(),
97            "One or more of the arguments is not finite"
98        );
99        Self {
100            matrix: [[1.0, 0.0, 0.0, x], [0.0, 1.0, 0.0, y], [0.0, 0.0, 1.0, z]],
101        }
102    }
103
104    /// Create a matrix defining a magnification
105    /// ## Arguments
106    /// * `f` the factor where 1.0 is the original size
107    /// ## Panics
108    /// It panics if `f` is not finite (`f64.is_finite()`)
109    #[must_use]
110    pub fn magnify(f: f64) -> Self {
111        assert!(f.is_finite(), "The factor is not finite");
112        Self {
113            matrix: [[f, 0.0, 0.0, 0.0], [0.0, f, 0.0, 0.0], [0.0, 0.0, f, 0.0]],
114        }
115    }
116
117    /// Create a matrix defining a magnification with three different factors
118    /// ## Arguments
119    /// * `x` the factor for the x dimension where 1.0 is the original size
120    /// * `y` the factor for the y dimension where 1.0 is the original size
121    /// * `z` the factor for the z dimension where 1.0 is the original size
122    /// ## Panics
123    /// It panics if any of the arguments is not finite (`f64.is_finite()`)
124    #[must_use]
125    pub fn scale(x: f64, y: f64, z: f64) -> Self {
126        assert!(
127            x.is_finite() && y.is_finite() && z.is_finite(),
128            "One or more of the arguments is not finite"
129        );
130        Self {
131            matrix: [[x, 0.0, 0.0, 0.0], [0.0, y, 0.0, 0.0], [0.0, 0.0, z, 0.0]],
132        }
133    }
134
135    /// This multiplies the translation with the given factors, this can be used to
136    /// convert fractional units into absolute units.
137    pub fn multiply_translation(&mut self, factors: (f64, f64, f64)) {
138        self.matrix[0][3] *= factors.0;
139        self.matrix[1][3] *= factors.1;
140        self.matrix[2][3] *= factors.2;
141    }
142
143    /// Apply this transformation to the given position.
144    /// It returns the new position.
145    /// ## Arguments
146    /// * `pos` the position (x, y, z)
147    #[must_use]
148    pub fn apply(&self, pos: (f64, f64, f64)) -> (f64, f64, f64) {
149        (
150            pos.2.mul_add(
151                self.matrix[0][2],
152                pos.0.mul_add(self.matrix[0][0], pos.1 * self.matrix[0][1]),
153            ) + self.matrix[0][3],
154            pos.2.mul_add(
155                self.matrix[1][2],
156                pos.0.mul_add(self.matrix[1][0], pos.1 * self.matrix[1][1]),
157            ) + self.matrix[1][3],
158            pos.2.mul_add(
159                self.matrix[2][2],
160                pos.0.mul_add(self.matrix[2][0], pos.1 * self.matrix[2][1]),
161            ) + self.matrix[2][3],
162        )
163    }
164
165    /// Combine this transformation with another transformation to deliver a new transformation.
166    /// This transformation is applied before the other transformation.
167    #[must_use]
168    pub fn combine(&self, other: &Self) -> Self {
169        Self {
170            matrix: [
171                [
172                    other.matrix[0][2].mul_add(
173                        self.matrix[2][0],
174                        other.matrix[0][0]
175                            .mul_add(self.matrix[0][0], other.matrix[0][1] * self.matrix[1][0]),
176                    ),
177                    other.matrix[0][2].mul_add(
178                        self.matrix[2][1],
179                        other.matrix[0][0]
180                            .mul_add(self.matrix[0][1], other.matrix[0][1] * self.matrix[1][1]),
181                    ),
182                    other.matrix[0][2].mul_add(
183                        self.matrix[2][2],
184                        other.matrix[0][0]
185                            .mul_add(self.matrix[0][2], other.matrix[0][1] * self.matrix[1][2]),
186                    ),
187                    other.matrix[0][2].mul_add(
188                        self.matrix[2][3],
189                        other.matrix[0][0]
190                            .mul_add(self.matrix[0][3], other.matrix[0][1] * self.matrix[1][3]),
191                    ) + other.matrix[0][3],
192                ],
193                [
194                    other.matrix[1][2].mul_add(
195                        self.matrix[2][0],
196                        other.matrix[1][0]
197                            .mul_add(self.matrix[0][0], other.matrix[1][1] * self.matrix[1][0]),
198                    ),
199                    other.matrix[1][2].mul_add(
200                        self.matrix[2][1],
201                        other.matrix[1][0]
202                            .mul_add(self.matrix[0][1], other.matrix[1][1] * self.matrix[1][1]),
203                    ),
204                    other.matrix[1][2].mul_add(
205                        self.matrix[2][2],
206                        other.matrix[1][0]
207                            .mul_add(self.matrix[0][2], other.matrix[1][1] * self.matrix[1][2]),
208                    ),
209                    other.matrix[1][2].mul_add(
210                        self.matrix[2][3],
211                        other.matrix[1][0]
212                            .mul_add(self.matrix[0][3], other.matrix[1][1] * self.matrix[1][3]),
213                    ) + other.matrix[1][3],
214                ],
215                [
216                    other.matrix[2][2].mul_add(
217                        self.matrix[2][0],
218                        other.matrix[2][0]
219                            .mul_add(self.matrix[0][0], other.matrix[2][1] * self.matrix[1][0]),
220                    ),
221                    other.matrix[2][2].mul_add(
222                        self.matrix[2][1],
223                        other.matrix[2][0]
224                            .mul_add(self.matrix[0][1], other.matrix[2][1] * self.matrix[1][1]),
225                    ),
226                    other.matrix[2][2].mul_add(
227                        self.matrix[2][2],
228                        other.matrix[2][0]
229                            .mul_add(self.matrix[0][2], other.matrix[2][1] * self.matrix[1][2]),
230                    ),
231                    other.matrix[2][2].mul_add(
232                        self.matrix[2][3],
233                        other.matrix[2][0]
234                            .mul_add(self.matrix[0][3], other.matrix[2][1] * self.matrix[1][3]),
235                    ) + other.matrix[2][3],
236                ],
237            ],
238        }
239    }
240}
241
242#[cfg(test)]
243#[allow(clippy::print_stdout, clippy::use_debug)]
244mod tests {
245    use super::TransformationMatrix;
246
247    #[test]
248    fn identity() {
249        let pos = (1.0, 1.0, 1.0);
250        let new_pos = TransformationMatrix::identity().apply(pos);
251        assert_eq!(pos, new_pos);
252    }
253
254    #[test]
255    fn combination() {
256        let pos = (10.0, 0.0, 0.0);
257        let new_pos = TransformationMatrix::rotation_y(90.0)
258            .combine(&TransformationMatrix::translation(0.0, 0.0, 10.0))
259            .apply(pos);
260        assert!(close_tuple(new_pos, (0.0, 0.0, 0.0)));
261    }
262    #[test]
263    fn rot_x() {
264        // 90 deg y
265        let mut pos = (0.0, 10.0, 0.0);
266        let mut new_pos = TransformationMatrix::rotation_x(90.0).apply(pos);
267        println!("{new_pos:?}");
268        assert!(close_tuple(new_pos, (0.0, 0.0, 10.0)));
269        // 90 deg z
270        pos = (0.0, 0.0, 10.0);
271        new_pos = TransformationMatrix::rotation_x(90.0).apply(pos);
272        println!("{new_pos:?}");
273        assert!(close_tuple(new_pos, (0.0, -10.0, 0.0)));
274        // -90 deg z
275        pos = (0.0, 0.0, 10.0);
276        new_pos = TransformationMatrix::rotation_x(-90.0).apply(pos);
277        println!("{new_pos:?}");
278        assert!(close_tuple(new_pos, (0.0, 10.0, 0.0)));
279        // 180 deg z
280        pos = (0.0, 0.0, 10.0);
281        new_pos = TransformationMatrix::rotation_x(180.0).apply(pos);
282        println!("{new_pos:?}");
283        assert!(close_tuple(new_pos, (0.0, 0.0, -10.0)));
284        // -180 deg z
285        pos = (0.0, 0.0, 10.0);
286        new_pos = TransformationMatrix::rotation_x(-180.0).apply(pos);
287        println!("{new_pos:?}");
288        assert!(close_tuple(new_pos, (0.0, 0.0, -10.0)));
289        // 360 deg z
290        pos = (0.0, 0.0, 10.0);
291        new_pos = TransformationMatrix::rotation_x(360.0).apply(pos);
292        println!("{new_pos:?}");
293        assert!(close_tuple(new_pos, (0.0, 0.0, 10.0)));
294        // 44.5 deg z
295        pos = (0.0, 0.0, -1.0);
296        new_pos = TransformationMatrix::rotation_x(44.5).apply(pos);
297        let end = (
298            0.0,
299            44.5_f64.to_radians().sin(),
300            (-44.5_f64).to_radians().cos(),
301        );
302        println!("{new_pos:?} vs {end:?}");
303        assert!(close_tuple(new_pos, end));
304        // 44.5 + 45.5 deg z
305        pos = (0.0, 0.0, -1.0);
306        new_pos = TransformationMatrix::rotation_x(44.5)
307            .combine(&TransformationMatrix::rotation_x(45.5))
308            .apply(pos);
309        println!("{new_pos:?}");
310        assert!(close_tuple(new_pos, (0.0, 1.0, 0.0)));
311    }
312    #[test]
313    fn rot_y() {
314        // 90 deg x
315        let mut pos = (10.0, 0.0, 0.0);
316        let mut new_pos = TransformationMatrix::rotation_y(90.0).apply(pos);
317        println!("{new_pos:?}");
318        assert!(close_tuple(new_pos, (0.0, 0.0, -10.0)));
319        // 90 deg z
320        pos = (0.0, 0.0, 10.0);
321        new_pos = TransformationMatrix::rotation_y(90.0).apply(pos);
322        println!("{new_pos:?}");
323        assert!(close_tuple(new_pos, (10.0, 0.0, 0.0)));
324        // -90 deg z
325        pos = (0.0, 0.0, 10.0);
326        new_pos = TransformationMatrix::rotation_y(-90.0).apply(pos);
327        println!("{new_pos:?}");
328        assert!(close_tuple(new_pos, (-10.0, 0.0, 0.0)));
329        // 180 deg z
330        pos = (0.0, 0.0, 10.0);
331        new_pos = TransformationMatrix::rotation_y(180.0).apply(pos);
332        println!("{new_pos:?}");
333        assert!(close_tuple(new_pos, (0.0, 0.0, -10.0)));
334        // -180 deg z
335        pos = (0.0, 0.0, 10.0);
336        new_pos = TransformationMatrix::rotation_y(-180.0).apply(pos);
337        println!("{new_pos:?}");
338        assert!(close_tuple(new_pos, (0.0, 0.0, -10.0)));
339        // 360 deg z
340        pos = (0.0, 0.0, 10.0);
341        new_pos = TransformationMatrix::rotation_y(360.0).apply(pos);
342        println!("{new_pos:?}");
343        assert!(close_tuple(new_pos, (0.0, 0.0, 10.0)));
344        // 44.5 deg z
345        pos = (0.0, 0.0, -1.0);
346        new_pos = TransformationMatrix::rotation_y(44.5).apply(pos);
347        let end = (
348            (-44.5_f64).to_radians().sin(),
349            0.0,
350            (-44.5_f64).to_radians().cos(),
351        );
352        println!("{new_pos:?} vs {end:?}");
353        assert!(close_tuple(new_pos, end));
354        // 44.5 + 45.5 deg z
355        pos = (0.0, 0.0, -1.0);
356        new_pos = TransformationMatrix::rotation_y(44.5)
357            .combine(&TransformationMatrix::rotation_y(45.5))
358            .apply(pos);
359        println!("{new_pos:?}");
360        assert!(close_tuple(new_pos, (-1.0, 0.0, 0.0)));
361    }
362
363    #[test]
364    fn rot_z() {
365        // 90 deg x
366        let mut pos = (10.0, 0.0, 0.0);
367        let mut new_pos = TransformationMatrix::rotation_z(90.0).apply(pos);
368        println!("{new_pos:?}");
369        assert!(close_tuple(new_pos, (0.0, 10.0, 0.0)));
370        // 90 deg y
371        pos = (0.0, 10.0, 0.0);
372        new_pos = TransformationMatrix::rotation_z(90.0).apply(pos);
373        println!("{new_pos:?}");
374        assert!(close_tuple(new_pos, (-10.0, 0.0, 0.0)));
375        // -90 deg y
376        pos = (0.0, 10.0, 0.0);
377        new_pos = TransformationMatrix::rotation_z(-90.0).apply(pos);
378        println!("{new_pos:?}");
379        assert!(close_tuple(new_pos, (10.0, 0.0, 0.0)));
380        // 180 deg y
381        pos = (0.0, 10.0, 0.0);
382        new_pos = TransformationMatrix::rotation_z(180.0).apply(pos);
383        println!("{new_pos:?}");
384        assert!(close_tuple(new_pos, (0.0, -10.0, 0.0)));
385        // -180 deg y
386        pos = (0.0, 10.0, 0.0);
387        new_pos = TransformationMatrix::rotation_z(-180.0).apply(pos);
388        println!("{new_pos:?}");
389        assert!(close_tuple(new_pos, (0.0, -10.0, 0.0)));
390        // 360 deg x
391        pos = (10.0, 0.0, 0.0);
392        new_pos = TransformationMatrix::rotation_z(360.0).apply(pos);
393        println!("{new_pos:?}");
394        assert!(close_tuple(new_pos, (10.0, 0.0, 0.0)));
395        // 44.5 deg y
396        pos = (0.0, -1.0, 0.0);
397        new_pos = TransformationMatrix::rotation_z(44.5).apply(pos);
398        let end = (
399            44.5_f64.to_radians().sin(),
400            (-44.5_f64).to_radians().cos(),
401            0.0,
402        );
403        println!("{new_pos:?} vs {end:?}");
404        assert!(close_tuple(new_pos, end));
405        // 44.5 + 45.5 deg y
406        pos = (0.0, -1.0, 0.0);
407        new_pos = TransformationMatrix::rotation_z(44.5)
408            .combine(&TransformationMatrix::rotation_z(45.5))
409            .apply(pos);
410        println!("{new_pos:?}");
411        assert!(close_tuple(new_pos, (1.0, 0.0, 0.0)));
412    }
413
414    #[test]
415    fn translation() {
416        let mut pos = (10.0, 0.0, 0.0);
417        let mut new_pos = TransformationMatrix::translation(-10.0, 0.0, 0.0).apply(pos);
418        assert!(close_tuple(new_pos, (0.0, 0.0, 0.0)));
419        pos = (-897.0, 0.0023, 1.0);
420        new_pos = TransformationMatrix::translation(897.0, -0.0023, -1.0).apply(pos);
421        assert!(close_tuple(new_pos, (0.0, 0.0, 0.0)));
422        pos = (0.0, 0.0, 0.0);
423        new_pos = TransformationMatrix::translation(0.0, 0.0, 0.0).apply(pos);
424        assert!(close_tuple(new_pos, (0.0, 0.0, 0.0)));
425    }
426
427    #[test]
428    fn magnification() {
429        let mut pos = (10.0, 0.0, 0.0);
430        let mut new_pos = TransformationMatrix::magnify(10.0).apply(pos);
431        assert!(close_tuple(new_pos, (100.0, 0.0, 0.0)));
432        pos = (-897.0, 0.0023, 1.0);
433        new_pos = TransformationMatrix::magnify(0.1).apply(pos);
434        assert!(close_tuple(new_pos, (-89.7, 0.00023, 0.1)));
435        pos = (0.0, 1.0, 0.0);
436        new_pos = TransformationMatrix::magnify(2.5).apply(pos);
437        assert!(close_tuple(new_pos, (0.0, 2.5, 0.0)));
438    }
439
440    #[test]
441    fn multiply_translation() {
442        let pos = (0.0, 0.0, 0.0);
443        let mut matrix = TransformationMatrix::translation(1.0, 2.0, -5.0);
444        matrix.multiply_translation((10.0, 5.0, -2.0));
445        let new_pos = matrix.apply(pos);
446        assert!(close_tuple(new_pos, (10.0, 10.0, 10.0)));
447        assert_eq!(
448            matrix.matrix(),
449            TransformationMatrix::translation(10.0, 10.0, 10.0).matrix()
450        );
451    }
452
453    #[test]
454    fn matrix() {
455        let normal = TransformationMatrix::rotation_x(45.0);
456        let raw = normal.matrix();
457        let from_matrix = TransformationMatrix::from_matrix(raw);
458        let mut set = TransformationMatrix::identity();
459        set.set_matrix(raw);
460        assert_eq!(normal, from_matrix);
461        assert_eq!(from_matrix, set);
462        assert_eq!(normal, set);
463    }
464
465    fn close_tuple(a: (f64, f64, f64), b: (f64, f64, f64)) -> bool {
466        close(a.0, b.0) && close(a.1, b.1) && close(a.2, b.2)
467    }
468
469    fn close(a: f64, b: f64) -> bool {
470        #[allow(clippy::float_cmp)]
471        if a == b {
472            true
473        } else if a == 0.0 || b == 0.0 {
474            (a - b) > -0.00000000000001 && (b - a) < 0.00000000000001
475        } else {
476            let dif = a / b;
477            (1.0 - dif) > -0.00000000000001 && (dif - 1.0) < 0.00000000000001
478        }
479    }
480}
pdbtbx/transformation.rs

pdbtbx/
transformation.rs
