diff --git a/wpimath/src/main/java/edu/wpi/first/math/geometry/Pose3d.java b/wpimath/src/main/java/edu/wpi/first/math/geometry/Pose3d.java index 38d7513128..9926b33465 100644 --- a/wpimath/src/main/java/edu/wpi/first/math/geometry/Pose3d.java +++ b/wpimath/src/main/java/edu/wpi/first/math/geometry/Pose3d.java @@ -205,18 +205,25 @@ public class Pose3d implements Interpolatable { final var u = VecBuilder.fill(twist.dx, twist.dy, twist.dz); final var rvec = VecBuilder.fill(twist.rx, twist.ry, twist.rz); final var omega = rotationVectorToMatrix(rvec); - final var omegaSq = omega.pow(2); + final var omegaSq = omega.times(omega); double theta = rvec.norm(); double thetaSq = theta * theta; double A; double B; double C; - if (theta < 1E-9) { + if (Math.abs(theta) < 1E-9) { // Taylor Expansions around θ = 0 // A = 1/1! - θ²/3! + θ⁴/5! // B = 1/2! - θ²/4! + θ⁴/6! // C = 1/3! - θ²/5! + θ⁴/7! + // sources: + // A: + // https://www.wolframalpha.com/input?i2d=true&i=series+expansion+of+Divide%5Bsin%5C%2840%29x%5C%2841%29%2Cx%5D+at+x%3D0 + // B: + // https://www.wolframalpha.com/input?i2d=true&i=series+expansion+of+Divide%5B1-cos%5C%2840%29x%5C%2841%29%2CPower%5Bx%2C2%5D%5D+at+x%3D0 + // C: + // https://www.wolframalpha.com/input?i2d=true&i=series+expansion+of+Divide%5B1-Divide%5Bsin%5C%2840%29x%5C%2841%29%2Cx%5D%2CPower%5Bx%2C2%5D%5D+at+x%3D0 A = 1 - thetaSq / 6 + thetaSq * thetaSq / 120; B = 1 / 2.0 - thetaSq / 24 + thetaSq * thetaSq / 720; C = 1 / 6.0 - thetaSq / 120 + thetaSq * thetaSq / 5040; @@ -257,16 +264,23 @@ public class Pose3d implements Interpolatable { final var rvec = transform.getRotation().getQuaternion().toRotationVector(); final var omega = rotationVectorToMatrix(rvec); - final var theta = transform.getRotation().getAngle(); + final var theta = rvec.norm(); final var thetaSq = theta * theta; double C; - if (theta < 1E-9) { + if (Math.abs(theta) < 1E-9) { // Taylor Expansions around θ = 0 // A = 1/1! - θ²/3! + θ⁴/5! // B = 1/2! - θ²/4! + θ⁴/6! // C = 1/6 * (1/2 + θ²/5! + θ⁴/7!) - C = 1 / 6.0 - thetaSq / 120 + thetaSq * thetaSq / 5040; + // sources: + // A: + // https://www.wolframalpha.com/input?i2d=true&i=series+expansion+of+Divide%5Bsin%5C%2840%29x%5C%2841%29%2Cx%5D+at+x%3D0 + // B: + // https://www.wolframalpha.com/input?i2d=true&i=series+expansion+of+Divide%5B1-cos%5C%2840%29x%5C%2841%29%2CPower%5Bx%2C2%5D%5D+at+x%3D0 + // C: + // https://www.wolframalpha.com/input?i2d=true&i=series+expansion+of+Divide%5B1-Divide%5BDivide%5Bsin%5C%2840%29x%5C%2841%29%2Cx%5D%2C2Divide%5B1-cos%5C%2840%29x%5C%2841%29%2CPower%5Bx%2C2%5D%5D%5D%2CPower%5Bx%2C2%5D%5D+at+x%3D0 + C = 1 / 12.0 + thetaSq / 720 + thetaSq * thetaSq / 30240; } else { // A = sin(θ)/θ // B = (1 - cos(θ)) / θ² @@ -276,7 +290,8 @@ public class Pose3d implements Interpolatable { C = (1 - A / (2 * B)) / thetaSq; } - final var V_inv = Matrix.eye(Nat.N3()).minus(omega.times(0.5)).plus(omega.pow(2).times(C)); + final var V_inv = + Matrix.eye(Nat.N3()).minus(omega.times(0.5)).plus(omega.times(omega).times(C)); final var twist_translation = V_inv.times(VecBuilder.fill(transform.getX(), transform.getY(), transform.getZ())); diff --git a/wpimath/src/main/java/edu/wpi/first/math/geometry/Rotation3d.java b/wpimath/src/main/java/edu/wpi/first/math/geometry/Rotation3d.java index a0f46f6626..4eb518d23e 100644 --- a/wpimath/src/main/java/edu/wpi/first/math/geometry/Rotation3d.java +++ b/wpimath/src/main/java/edu/wpi/first/math/geometry/Rotation3d.java @@ -72,6 +72,17 @@ public class Rotation3d implements Interpolatable { cr * cp * sy - sr * sp * cy); } + /** + * Constructs a Rotation3d with the given rotation vector representation. This representation is + * equivalent to axis-angle, where the normalized axis is multiplied by the rotation around the + * axis in radians. + * + * @param rvec The rotation vector. + */ + public Rotation3d(Vector rvec) { + this(rvec, rvec.norm()); + } + /** * Constructs a Rotation3d with the given axis-angle representation. The axis doesn't have to be * normalized. diff --git a/wpimath/src/main/native/cpp/geometry/Pose3d.cpp b/wpimath/src/main/native/cpp/geometry/Pose3d.cpp index 65afa58cb4..7eef6d5c9a 100644 --- a/wpimath/src/main/native/cpp/geometry/Pose3d.cpp +++ b/wpimath/src/main/native/cpp/geometry/Pose3d.cpp @@ -72,21 +72,28 @@ Pose3d Pose3d::RelativeTo(const Pose3d& other) const { Pose3d Pose3d::Exp(const Twist3d& twist) const { // Implementation from Section 3.2 of https://ethaneade.org/lie.pdf - auto u = Vectord<3>{twist.dx.value(), twist.dy.value(), twist.dz.value()}; - auto rvec = Vectord<3>{twist.rx.value(), twist.ry.value(), twist.rz.value()}; - auto omega = RotationVectorToMatrix(rvec); - auto omegaSq = omega * omega; - auto theta = rvec.norm(); - auto thetaSq = theta * theta; + Eigen::Vector3d u{twist.dx.value(), twist.dy.value(), twist.dz.value()}; + Eigen::Vector3d rvec{twist.rx.value(), twist.ry.value(), twist.rz.value()}; + Eigen::Matrix3d omega = RotationVectorToMatrix(rvec); + Eigen::Matrix3d omegaSq = omega * omega; + double theta = rvec.norm(); + double thetaSq = theta * theta; double A; double B; double C; - if (theta < 1E-9) { + if (std::abs(theta) < 1E-9) { // Taylor Expansions around θ = 0 // A = 1/1! - θ²/3! + θ⁴/5! // B = 1/2! - θ²/4! + θ⁴/6! // C = 1/3! - θ²/5! + θ⁴/7! + // sources: + // A: + // https://www.wolframalpha.com/input?i2d=true&i=series+expansion+of+Divide%5Bsin%5C%2840%29x%5C%2841%29%2Cx%5D+at+x%3D0 + // B: + // https://www.wolframalpha.com/input?i2d=true&i=series+expansion+of+Divide%5B1-cos%5C%2840%29x%5C%2841%29%2CPower%5Bx%2C2%5D%5D+at+x%3D0 + // C: + // https://www.wolframalpha.com/input?i2d=true&i=series+expansion+of+Divide%5B1-Divide%5Bsin%5C%2840%29x%5C%2841%29%2Cx%5D%2CPower%5Bx%2C2%5D%5D+at+x%3D0 A = 1 - thetaSq / 6 + thetaSq * thetaSq / 120; B = 1 / 2.0 - thetaSq / 24 + thetaSq * thetaSq / 720; C = 1 / 6.0 - thetaSq / 120 + thetaSq * thetaSq / 5040; @@ -99,8 +106,8 @@ Pose3d Pose3d::Exp(const Twist3d& twist) const { C = (1 - A) / thetaSq; } - auto R = Matrixd<3, 3>::Identity() + A * omega + B * omegaSq; - auto V = Matrixd<3, 3>::Identity() + B * omega + C * omegaSq; + Eigen::Matrix3d R = Eigen::Matrix3d::Identity() + A * omega + B * omegaSq; + Eigen::Matrix3d V = Eigen::Matrix3d::Identity() + B * omega + C * omegaSq; auto translation_component = V * u; const Transform3d transform{ @@ -116,22 +123,30 @@ Twist3d Pose3d::Log(const Pose3d& end) const { // Implementation from Section 3.2 of https://ethaneade.org/lie.pdf const auto transform = end.RelativeTo(*this); - auto u = Vectord<3>{transform.X().value(), transform.Y().value(), - transform.Z().value()}; - auto rvec = transform.Rotation().GetQuaternion().ToRotationVector(); + Eigen::Vector3d u{transform.X().value(), transform.Y().value(), + transform.Z().value()}; + Eigen::Vector3d rvec = + transform.Rotation().GetQuaternion().ToRotationVector(); - auto omega = RotationVectorToMatrix(rvec); - auto omegaSq = omega * omega; - auto theta = rvec.norm(); - auto thetaSq = theta * theta; + Eigen::Matrix3d omega = RotationVectorToMatrix(rvec); + Eigen::Matrix3d omegaSq = omega * omega; + double theta = rvec.norm(); + double thetaSq = theta * theta; double C; - if (theta < 1E-9) { + if (std::abs(theta) < 1E-9) { // Taylor Expansions around θ = 0 // A = 1/1! - θ²/3! + θ⁴/5! // B = 1/2! - θ²/4! + θ⁴/6! // C = 1/6 * (1/2 + θ²/5! + θ⁴/7!) - C = 1 / 6.0 - thetaSq / 120 + thetaSq * thetaSq / 5040; + // sources: + // A: + // https://www.wolframalpha.com/input?i2d=true&i=series+expansion+of+Divide%5Bsin%5C%2840%29x%5C%2841%29%2Cx%5D+at+x%3D0 + // B: + // https://www.wolframalpha.com/input?i2d=true&i=series+expansion+of+Divide%5B1-cos%5C%2840%29x%5C%2841%29%2CPower%5Bx%2C2%5D%5D+at+x%3D0 + // C: + // https://www.wolframalpha.com/input?i2d=true&i=series+expansion+of+Divide%5B1-Divide%5BDivide%5Bsin%5C%2840%29x%5C%2841%29%2Cx%5D%2C2Divide%5B1-cos%5C%2840%29x%5C%2841%29%2CPower%5Bx%2C2%5D%5D%5D%2CPower%5Bx%2C2%5D%5D+at+x%3D0 + C = 1 / 12.0 + thetaSq / 720 + thetaSq * thetaSq / 30240; } else { // A = std::sin(θ)/θ // B = (1 - std::cos(θ)) / θ² @@ -141,9 +156,10 @@ Twist3d Pose3d::Log(const Pose3d& end) const { C = (1 - A / (2 * B)) / thetaSq; } - auto V_inv = Matrixd<3, 3>::Identity() - 0.5 * omega + C * omegaSq; + Eigen::Matrix3d V_inv = + Eigen::Matrix3d::Identity() - 0.5 * omega + C * omegaSq; - auto translation_component = V_inv * u; + Eigen::Vector3d translation_component = V_inv * u; return Twist3d{units::meter_t{translation_component(0)}, units::meter_t{translation_component(1)}, diff --git a/wpimath/src/main/native/cpp/geometry/Rotation3d.cpp b/wpimath/src/main/native/cpp/geometry/Rotation3d.cpp index b109682f70..a0bc176de7 100644 --- a/wpimath/src/main/native/cpp/geometry/Rotation3d.cpp +++ b/wpimath/src/main/native/cpp/geometry/Rotation3d.cpp @@ -38,6 +38,9 @@ Rotation3d::Rotation3d(units::radian_t roll, units::radian_t pitch, cr * sp * cy + sr * cp * sy, cr * cp * sy - sr * sp * cy}; } +Rotation3d::Rotation3d(const Eigen::Vector3d& rvec) + : Rotation3d{rvec, units::radian_t{rvec.norm()}} {} + Rotation3d::Rotation3d(const Vectord<3>& axis, units::radian_t angle) { double norm = axis.norm(); if (norm == 0.0) { diff --git a/wpimath/src/main/native/include/frc/geometry/Rotation3d.h b/wpimath/src/main/native/include/frc/geometry/Rotation3d.h index 7c1a60d65a..df01f3b689 100644 --- a/wpimath/src/main/native/include/frc/geometry/Rotation3d.h +++ b/wpimath/src/main/native/include/frc/geometry/Rotation3d.h @@ -59,13 +59,22 @@ class WPILIB_DLLEXPORT Rotation3d { */ Rotation3d(const Vectord<3>& axis, units::radian_t angle); + /** + * Constructs a Rotation3d with the given rotation vector representation. This + * representation is equivalent to axis-angle, where the normalized axis is + * multiplied by the rotation around the axis in radians. + * + * @param rvec The rotation vector. + */ + explicit Rotation3d(const Eigen::Vector3d& rvec); + /** * Constructs a Rotation3d from a rotation matrix. * * @param rotationMatrix The rotation matrix. * @throws std::domain_error if the rotation matrix isn't special orthogonal. */ - explicit Rotation3d(const Matrixd<3, 3>& rotationMatrix); + explicit Rotation3d(const Eigen::Matrix3d& rotationMatrix); /** * Constructs a Rotation3d that rotates the initial vector onto the final diff --git a/wpimath/src/test/java/edu/wpi/first/math/geometry/Rotation3dTest.java b/wpimath/src/test/java/edu/wpi/first/math/geometry/Rotation3dTest.java index 072a2e64bd..95ec904ff1 100644 --- a/wpimath/src/test/java/edu/wpi/first/math/geometry/Rotation3dTest.java +++ b/wpimath/src/test/java/edu/wpi/first/math/geometry/Rotation3dTest.java @@ -24,17 +24,23 @@ class Rotation3dTest { final var xAxis = VecBuilder.fill(1.0, 0.0, 0.0); final var rot1 = new Rotation3d(xAxis, Math.PI / 3); final var rot2 = new Rotation3d(Math.PI / 3, 0.0, 0.0); + final var rvec1 = new Rotation3d(xAxis.times(Math.PI / 3)); assertEquals(rot1, rot2); + assertEquals(rot1, rvec1); final var yAxis = VecBuilder.fill(0.0, 1.0, 0.0); final var rot3 = new Rotation3d(yAxis, Math.PI / 3); final var rot4 = new Rotation3d(0.0, Math.PI / 3, 0.0); + final var rvec2 = new Rotation3d(yAxis.times(Math.PI / 3)); assertEquals(rot3, rot4); + assertEquals(rot3, rvec2); final var zAxis = VecBuilder.fill(0.0, 0.0, 1.0); final var rot5 = new Rotation3d(zAxis, Math.PI / 3); final var rot6 = new Rotation3d(0.0, 0.0, Math.PI / 3); + final var rvec3 = new Rotation3d(zAxis.times(Math.PI / 3)); assertEquals(rot5, rot6); + assertEquals(rot5, rvec3); } @Test diff --git a/wpimath/src/test/native/cpp/geometry/Rotation3dTest.cpp b/wpimath/src/test/native/cpp/geometry/Rotation3dTest.cpp index 4709ed0d0b..e41e170cb8 100644 --- a/wpimath/src/test/native/cpp/geometry/Rotation3dTest.cpp +++ b/wpimath/src/test/native/cpp/geometry/Rotation3dTest.cpp @@ -17,17 +17,23 @@ TEST(Rotation3dTest, InitAxisAngleAndRollPitchYaw) { const Eigen::Vector3d xAxis{1.0, 0.0, 0.0}; const Rotation3d rot1{xAxis, units::radian_t{std::numbers::pi / 3}}; const Rotation3d rot2{units::radian_t{std::numbers::pi / 3}, 0_rad, 0_rad}; + const Rotation3d rvec1{Eigen::Vector3d{xAxis * std::numbers::pi / 3}}; EXPECT_EQ(rot1, rot2); + EXPECT_EQ(rot1, rvec1); const Eigen::Vector3d yAxis{0.0, 1.0, 0.0}; const Rotation3d rot3{yAxis, units::radian_t{std::numbers::pi / 3}}; const Rotation3d rot4{0_rad, units::radian_t{std::numbers::pi / 3}, 0_rad}; + const Rotation3d rvec2{Eigen::Vector3d{yAxis * std::numbers::pi / 3}}; EXPECT_EQ(rot3, rot4); + EXPECT_EQ(rot3, rvec2); const Eigen::Vector3d zAxis{0.0, 0.0, 1.0}; const Rotation3d rot5{zAxis, units::radian_t{std::numbers::pi / 3}}; const Rotation3d rot6{0_rad, 0_rad, units::radian_t{std::numbers::pi / 3}}; + const Rotation3d rvec3{Eigen::Vector3d{zAxis * std::numbers::pi / 3}}; EXPECT_EQ(rot5, rot6); + EXPECT_EQ(rot5, rvec3); } TEST(Rotation3dTest, InitRotationMatrix) {