This bugzilla service is closed. All entries have been migrated to https://gitlab.com/libeigen/eigen
Bug 1312 - AngleAxis from QuaternionBase , better precision
AngleAxis from QuaternionBase , better precision
 Status: RESOLVED FIXED None Eigen Unclassified Geometry (show other bugs) 3.2 All All Normal Accuracy Problem Nobody 1244 Show dependency tree / graph

 Reported: 2016-09-26 09:21 UTC by Michael Norel 2019-12-04 16:21 UTC (History) 4 users (show) chtz gael.guennebaud hauke.heibel jacob.benoit.1

Attachments

 Michael Norel 2016-09-26 09:21:27 UTC ```1. m_angle = Scalar(2)*acos((min)((max)(Scalar(-1),q.w()),Scalar(1))); Angle can be computed with better precision using atan2(sqrt(n2), q.w()), without check for "1" overflow. 2. Also, for better precision , it is good to switch sign of input quaternion to have positive q.w(), it keeps rotation angle in [0, PI]. 3. The check "if (n2 < NumTraits::dummy_precision()*NumTraits::dummy_precision())" can be changed to strong comparison with "zero", and provide same relative precision for small rotations without side effects.``` Michael Norel 2016-09-28 13:49:52 UTC `Sorry , this is dublicate, i posted this year ago. But still , not fixed in 3.xx.` Gael Guennebaud 2016-09-29 11:46:54 UTC ```hm, the current implementation is: using std::atan2; Scalar n = q.vec().norm(); if(n::epsilon()) n = q.vec().stableNorm(); if (n > Scalar(0)) { m_angle = Scalar(2)*atan2(n, q.w()); m_axis = q.vec() / n; } else { m_angle = Scalar(0); m_axis << Scalar(1), Scalar(0), Scalar(0); } return *this; So 1. and 3. were already "fixed". Regarding 2, I would prefer to keep it like this so that angle-axis -> quat -> angle-axis returns exactly the same axis.``` Michael Norel 2016-09-29 13:29:06 UTC ```Yes, my mistake , i checked "Eigen 3.2.9" instead of current "Eigen 3.3-rc1". Thank for correction, ill try to avoid such mistakes. About 2) "that angle-axis -> quat -> angle-axis returns exactly the same axis." If this is not important issue, a can demonstrate the loose of accuracy. double small_angle = 1.0e-16; // any value smaller epsilon double a1 = atan2(small_angle, -1.0);// == PI - small_angle; double a2 = atan2(small_angle, 1.0);// == small_angle; double a1Check = acos(-1.0) - a1;// a1Check is zero here :( For all angles less then epsilon * PI, rotation angle will be lost. So, current conversion: - loose absolute precision for small angles. - breaks relative precision for small angles. fix can be simple ... Scalar w = q.w(); if(w < 0) { w=-w; n=-n; } if (n != Scalar(0)) ... but introduce additional branch. Here, i prefer precision on full range of angles over performance. Anyway, it is not a big issue.``` Gael Guennebaud 2016-09-29 20:37:30 UTC `Sorry, you're right, and in my angleaxis-quat->angleaxis scenario, w will be positive anyway (unless the angle is outside [-2pi,2pi])` Gael Guennebaud 2016-09-29 20:44:20 UTC ```(In reply to Gael Guennebaud from comment #4) > (unless the angle is outside [-2pi,2pi]) I meant [-pi,pi].``` Gael Guennebaud 2016-09-29 21:24:22 UTC `Fixed: https://bitbucket.org/eigen/eigen/commits/5998f71df011` Michael Norel 2016-09-30 08:00:22 UTC ```As i see from patch, code looks template template AngleAxis& AngleAxis::operator=(const QuaternionBase& q) { using std::atan2; using std::abs; Scalar n = q.vec().norm(); if(n::epsilon()) n = q.vec().stableNorm(); if (n != Scalar(0)) { m_angle = Scalar(2)*atan2(n, std::abs(q.w())); if(q.w() < 0) n = -n; m_axis = q.vec() / n; } else { m_angle = Scalar(0); m_axis << Scalar(1), Scalar(0), Scalar(0); } return *this; } 1. "n" is positive here , and (n > Scalar(0)) still can be used 2. m_angle = Scalar(2)*atan2(n, std::abs(q.w())); returns value in [0,pi], so this comment should be changed "angle is in the [-pi,pi] range."``` Nobody 2019-12-04 16:21:26 UTC ```-- GitLab Migration Automatic Message -- This bug has been migrated to gitlab.com's GitLab instance and has been closed from further activity. You can subscribe and participate further through the new bug through this link to our GitLab instance: https://gitlab.com/libeigen/eigen/issues/1312.```

 Note You need to log in before you can comment on or make changes to this bug.