cv_core/so3.rs
1use derive_more::{AsMut, AsRef, Deref, DerefMut, From, Into};
2use nalgebra::{Matrix3, Matrix4, Rotation3, Unit, Vector3};
3use num_traits::Float;
4
5/// Contains a member of the lie algebra so(3), a representation of the tangent space
6/// of 3d rotation. This is also known as the lie algebra of the 3d rotation group SO(3).
7///
8/// This is only intended to be used in optimization problems where it is desirable to
9/// have unconstranied variables representing the degrees of freedom of the rotation.
10/// In all other cases, a rotation matrix should be used to store rotations, since the
11/// conversion to and from a rotation matrix is non-trivial.
12#[derive(Debug, Clone, Copy, PartialEq, PartialOrd, AsMut, AsRef, Deref, DerefMut, From, Into)]
13pub struct Skew3(pub Vector3<f64>);
14
15impl Skew3 {
16 /// Converts the Skew3 to a Rotation3 matrix.
17 pub fn rotation(self) -> Rotation3<f64> {
18 self.into()
19 }
20
21 /// This converts a matrix in skew-symmetric form into a Skew3.
22 ///
23 /// Warning: Does no check to ensure matrix is actually skew-symmetric.
24 pub fn vee(mat: Matrix3<f64>) -> Self {
25 Self(Vector3::new(mat.m32, mat.m13, mat.m21))
26 }
27
28 /// This converts the Skew3 into its skew-symmetric matrix form.
29 #[rustfmt::skip]
30 pub fn hat(self) -> Matrix3<f64> {
31 self.0.cross_matrix()
32 }
33
34 /// This converts the Skew3 into its squared skew-symmetric matrix form efficiently.
35 #[rustfmt::skip]
36 pub fn hat2(self) -> Matrix3<f64> {
37 let w = self.0;
38 let w11 = w.x * w.x;
39 let w12 = w.x * w.y;
40 let w13 = w.x * w.z;
41 let w22 = w.y * w.y;
42 let w23 = w.y * w.z;
43 let w33 = w.z * w.z;
44 Matrix3::new(
45 -w22 - w33, w12, w13,
46 w12, -w11 - w33, w23,
47 w13, w23, -w11 - w22,
48 )
49 }
50
51 /// Computes the lie bracket [self, rhs].
52 pub fn bracket(self, rhs: Self) -> Self {
53 Self::vee(self.hat() * rhs.hat() - rhs.hat() * self.hat())
54 }
55
56 /// The jacobian of the output of a rotation in respect to the
57 /// input of a rotation.
58 ///
59 /// `y = R * x`
60 ///
61 /// `dy/dx = R`
62 ///
63 /// The formula is pretty simple and is just the rotation matrix created
64 /// from the exponential map of this so(3) element into SO(3). The result is converted
65 /// to homogeneous form (by adding a new dimension with a `1` in the diagonal) so
66 /// that it is compatible with homogeneous coordinates.
67 ///
68 /// If you have the rotation matrix already, please use the rotation matrix itself
69 /// rather than calling this method. Calling this method will waste time converting
70 /// the [`Skew3`] back into a [`Rotation3`], which is non-trivial.
71 pub fn jacobian_input(self) -> Matrix4<f64> {
72 let rotation: Rotation3<f64> = self.into();
73 let matrix: Matrix3<f64> = rotation.into();
74 matrix.to_homogeneous()
75 }
76
77 /// The jacobian of the output of a rotation in respect to the
78 /// rotation itself.
79 ///
80 /// `y = R * x`
81 ///
82 /// `dy/dR = -hat(y)`
83 ///
84 /// The derivative is purely based on the current output vector, and thus doesn't take `self`.
85 ///
86 /// Note that when working with homogeneous projective coordinates, only the first three components
87 /// (the bearing) are relevant, hence the resulting matrix is a [`Matrix3`].
88 pub fn jacobian_self(y: Vector3<f64>) -> Matrix3<f64> {
89 y.cross_matrix()
90 }
91}
92
93/// This is the exponential map.
94impl From<Skew3> for Rotation3<f64> {
95 fn from(w: Skew3) -> Self {
96 // This check is done to avoid the degenerate case where the angle is near zero.
97 let theta2 = w.0.norm_squared();
98 if theta2 <= f64::epsilon() {
99 Rotation3::from_matrix(&(Matrix3::identity() + w.hat()))
100 } else {
101 let theta = theta2.sqrt();
102 let axis = Unit::new_unchecked(w.0 / theta);
103 Self::from_axis_angle(&axis, theta)
104 }
105 }
106}
107
108/// This is the log map.
109impl From<Rotation3<f64>> for Skew3 {
110 fn from(r: Rotation3<f64>) -> Self {
111 Self(r.scaled_axis())
112 }
113}