* mln/fun/x2x/rotation.hh: (get_rot_h_mat 3d): fix matrix construction. (rotation(algebra::quat)): update attributes in order to produce a valid inverted rotation. * tests/fun/x2x/rotation.cc: improve test. --- milena/ChangeLog | 19 +++++++++++++ milena/mln/fun/x2x/rotation.hh | 36 ++++++++++++++++++------- milena/tests/fun/x2x/rotation.cc | 55 ++++++++++++++----------------------- 3 files changed, 66 insertions(+), 44 deletions(-) diff --git a/milena/ChangeLog b/milena/ChangeLog index 62c0696..93b5737 100644 --- a/milena/ChangeLog +++ b/milena/ChangeLog @@ -1,5 +1,24 @@ 2009-02-18 Guillaume Lazzara <z@lrde.epita.fr> + Fix fun::x2v::rotation. + + * mln/fun/x2x/rotation.hh: + (get_rot_h_mat 3d): fix matrix construction. + (rotation(algebra::quat)): update attributes in order to produce a + valid inverted rotation. + + * tests/fun/x2x/rotation.cc: improve test. + +2009-02-18 Guillaume Lazzara <z@lrde.epita.fr> + + Add fun::x2v::trilinear. + + * mln/fun/x2v/trilinear.hh: new file. New interpolation algorithm. + + * mln/fun/x2v/all.hh: include new file. + +2009-02-18 Guillaume Lazzara <z@lrde.epita.fr> + Add aliases for 2D algebra::vec. * mln/core/alias/vec2d.hh: new. new aliases for float and double 2d diff --git a/milena/mln/fun/x2x/rotation.hh b/milena/mln/fun/x2x/rotation.hh index dce37ed..9628c26 100644 --- a/milena/mln/fun/x2x/rotation.hh +++ b/milena/mln/fun/x2x/rotation.hh @@ -32,6 +32,9 @@ /// \file mln/fun/x2x/rotation.hh /// /// Define a rotation function. +/// +/// \todo store the quaternion instead of (axis, alpha) +/// => better precision while composing two rotation matrices. # include <cmath> # include <mln/core/concept/function.hh> @@ -70,6 +73,10 @@ namespace mln //FIXME: cannot use '!=' operator. mln_precondition(!(axis_ == vec_t(literal::zero))); + /// Be sure that axis is normalized. + algebra::vec<3,C> axis = axis_; + axis.normalize(); + const float cos_a = cos(alpha_); const float sin_a = sin(alpha_); const float u = axis_[0]; @@ -78,23 +85,22 @@ namespace mln const float u2 = u * u; const float v2 = v * v; const float w2 = w * w; - const float uvw2 = u2 + v2 + w2; algebra::h_mat<3, C> m_; - m_(0,0) = (u2 + (v2 + w2) * cos_a) / uvw2; - m_(0,1) = (u*v * (1 - cos_a) - u * std::sqrt(uvw2) * sin_a) / uvw2; - m_(0,2) = (u*w * (1 - cos_a) + v * std::sqrt(uvw2) * sin_a) / uvw2; + m_(0,0) = u2 + (1 - u2) * cos_a; + m_(0,1) = u*v * (1 - cos_a) - w * sin_a; + m_(0,2) = u*w * (1 - cos_a) + v * sin_a; m_(0,3) = 0; - m_(1,0) = (u*v * (1 - cos_a) + w * std::sqrt(uvw2) * sin_a) / uvw2; - m_(1,1) = (v2 + (u2 + w2) * cos_a) / uvw2; - m_(1,2) = (v*w * (1 - cos_a) - u * std::sqrt(uvw2) * sin_a) / uvw2; + m_(1,0) = u*v * (1 - cos_a) + w * sin_a; + m_(1,1) = v2 + (1 - v2) * cos_a; + m_(1,2) = v * w * (1 - cos_a) - u * sin_a; m_(1,3) = 0; - m_(2,0) = (u*w * (1 - cos_a) - v * std::sqrt(uvw2) * sin_a) / uvw2; - m_(2,1) = (v*w * (1 - cos_a) + u * std::sqrt(uvw2) * sin_a) / uvw2; - m_(2,2) = (u2 + (u2 + v2) * cos_a) / uvw2; + m_(2,0) = u * w * (1 - cos_a) - v * sin_a; + m_(2,1) = u * w * (1 - cos_a) + u * sin_a; + m_(2,2) = w2 + (1 - w2) * cos_a; m_(2,3) = 0; m_(3,0) = 0; @@ -189,19 +195,29 @@ namespace mln rotation<n,C>::rotation(const algebra::quat& q) { mln_precondition(q.is_unit()); + // FIXME: Should also work for 2d. mln_precondition(n == 3); + float w = q.to_vec()[0], x = q.to_vec()[1], x2 = 2*x*x, xw = 2*x*w, y = q.to_vec()[2], y2 = 2*y*y, xy = 2*x*y, yw = 2*y*w, z = q.to_vec()[3], z2 = 2*z*z, xz = 2*x*z, yz = 2*y*z, zw = 2*z*w; + float t[9] = {1.f - y2 - z2, xy - zw, xz + yw, xy + zw, 1.f - x2 - z2, yz - xw, xz - yw, yz + xw, 1.f - x2 - y2}; this->m_ = mln::make::h_mat(t); mln_assertion(check_rotation(q)); + + /// Update attributes + + alpha_ = acos(w) * 2; + axis_[0] = x; + axis_[1] = y; + axis_[2] = z; } diff --git a/milena/tests/fun/x2x/rotation.cc b/milena/tests/fun/x2x/rotation.cc index b96c551..fc36a17 100644 --- a/milena/tests/fun/x2x/rotation.cc +++ b/milena/tests/fun/x2x/rotation.cc @@ -33,14 +33,9 @@ #include <iostream> #include <mln/fun/x2x/rotation.hh> -#include <mln/core/image/image2d.hh> -#include <mln/value/int_u8.hh> -#include <mln/io/pgm/load.hh> -#include <mln/io/pgm/save.hh> -#include <mln/core/image/interpolated.hh> +#include <mln/core/alias/vec2d.hh> +#include <mln/core/alias/vec3d.hh> #include <mln/make/vec.hh> -#include <mln/fun/x2v/bilinear.hh> -#include <mln/extension/adjust.hh> #include "tests/data.hh" @@ -48,36 +43,28 @@ int main() { using namespace mln; - using value::int_u8; - algebra::vec<2,float> axis; - axis[0] = 0; - axis[1] = 1; - typedef image2d<int_u8> ima_t; - ima_t lena; - io::pgm::load(lena, MLN_IMG_DIR "/lena.pgm"); + { + algebra::vec<2,float> axis; + axis[0] = 0; + axis[1] = 1; + fun::x2x::rotation<2,float> rot(0.1f, axis); - ima_t out; - initialize(out, lena); + vec2d_f v = make::vec(4, 4); + vec2d_f res = rot(v); + mln_assertion(rot.inv()(res) == v); + } - interpolated<ima_t, fun::x2v::bilinear> inter(lena); + { + algebra::vec<3,float> axis; + axis[0] = 1; + axis[1] = 0; + axis[2] = 0; + fun::x2x::rotation<3,float> rot(0.1f, axis); - fun::x2x::rotation<2,float> rot1(0.1, axis); + vec3d_f v = make::vec(4, 1, 2); + vec3d_f res = rot(v); + mln_assertion(rot.inv()(res) == v); + } - mln_piter_(ima_t) p(out.domain()); - for_all(p) - { - algebra::vec<2,float> v = rot1.inv()(p.to_site().to_vec()); - if (inter.domain().has(v)) - out(p) = inter(v); - else - out(p) = 255; - } - io::pgm::save(out, "out.pgm"); - - fun::x2x::rotation<2,float> rot2(3.14116, axis); - mln_assertion(fabs(rot2(make::vec(0.0, 1.0))[0] - - make::vec(0.0, -1.0)[0]) <= 0.125); - mln_assertion(fabs(rot2(make::vec(0.0, 1.0))[1] - - make::vec(0.0, -1.0)[1]) <= 0.125); } -- 1.5.6.5