[Olena-patches] 3212: Fix regsitration.

https://svn.lrde.epita.fr/svn/oln/trunk/milena Index: ChangeLog from Ugo Jardonnet &lt;ugo.jardonnet(a)lrde.epita.fr&gt; Fix regsitration. * mln/registration/multiscale.hh: Fix error. * mln/registration/icp.hh: Fix call. * mln/registration/get_rtransf.hh: Fix minor error. fun/x2x/composed.hh | 4 ++-- registration/get_rtransf.hh | 14 ++++++++++++-- registration/icp.hh | 20 ++++++++++++++++++-- registration/multiscale.hh | 16 ++++++++++++++-- 4 files changed, 46 insertions(+), 8 deletions(-) Index: mln/registration/multiscale.hh --- mln/registration/multiscale.hh (revision 3212) +++ mln/registration/multiscale.hh (working copy) @@ -40,6 +40,17 @@ namespace registration { + /*! Register point in \p c using a multiscale icp algorithm + * + * \param[in] cloud Image to register + * \param[in] surface Image using to register onto + * \param[in] q Subsampling quotient + * \param[in] nb_it Number of registrations + * \param[out] The rigid transformation obtained. + * + * \pre \p cloud has to be initialized. + * \pre \p surface has to be initialized. + */ template <typename I, typename J> inline composed< rotation<I::psite::dim, float>, translation<I::Psite::dim, float> > @@ -113,7 +124,7 @@ //working box const box<mln_psite(I)> working_box = - enlarge(bigger(geom::bbox(c), geom::bbox(x)), 100); + enlarge(bigger(geom::bbox(c), geom::bbox(x)), 1000); //warning 100 could be insufficient //make a lazy_image map via function closest_point fun::x2p::closest_point<mln_psite(I)> fun(x, working_box); @@ -121,12 +132,13 @@ map(fun, fun.domain()); //init rigid transform qk + // FIXME: check qk initialisation. composed< rotation<I::psite::dim, float>, translation<I::psite::dim, float> > qk; //run registration for (int e = nb_it-1; e >= 0; e--) { - unsigned int l = cloud.nsites() / std::pow(q, e); + unsigned int l = c.nsites() / std::pow(q, e); l = (l < 1) ? 1 : l; icp_subset(c, l, map, qk); } Index: mln/registration/icp.hh --- mln/registration/icp.hh (revision 3212) +++ mln/registration/icp.hh (working copy) @@ -169,8 +169,18 @@ //work on a reduced version of c_ p_array<P> c = c_; + + std::cout << "c_length" << c_length << std::endl; + std::cout << "-nsites" << c_.nsites() << std::endl; + + std::cout << "bad point" << c[9100] << std::endl; + + //FIXME: Use c_ and specify c_length to function using c_ c.hook_std_vector_().resize(c_length); + std::cout << "+nsites" << c.nsites() << std::endl; + std::cout << "bad point" << c[9100] << std::endl; + p_array<P> ck(c); algebra::vec<P::dim,float> mu_c = center(c); @@ -180,6 +190,7 @@ unsigned int k = 0; do { + std::cout << "*------" << std::endl; buf_dk.store(d_k); //compute qk @@ -196,14 +207,16 @@ k++; //FIXME: Add matrix norm //} while ((qk - buf_qk[1]).sqr_norm() / qk.sqr_norm() > epsilon); - } while (d_k > epsilon); + + std::cout << d_k << std::endl; + + } while (d_k > 9.2); trace::exiting("registration::impl::icp_"); } } // end of namespace mln::registration::impl - template <typename P, typename M> inline composed<rotation<P::dim, float>, translation<P::dim, float> > @@ -211,6 +224,9 @@ const float epsilon = 1e-3) { composed<rotation<P::dim, float>, translation<P::dim, float> > qk; + + std::cout << " icp " << std::endl; + impl::icp_(c, c.nsites(), map, qk, epsilon); return qk; } Index: mln/registration/get_rtransf.hh --- mln/registration/get_rtransf.hh (revision 3212) +++ mln/registration/get_rtransf.hh (working copy) @@ -63,9 +63,19 @@ { //mu_xk = center map(Ck) algebra::vec<P::dim,float> mu_xk(literal::zero); + + std::cout << c.nsites() << std::endl; for (unsigned i = 0; i < c.nsites(); ++i) - mu_xk += convert::to< algebra::vec<P::dim,float> >(map( - convert::to< algebra::vec<P::dim,float> >(ck[i]))); + { + if (not map.has(ck[i])) + { + std::cout << ck[i] << " " << i << std::endl; + } + else + mu_xk += convert::to< algebra::vec<P::dim,float> > + (map(convert::to< algebra::vec<P::dim,float> >(ck[i]))); + + } mu_xk /= c.nsites(); // qR Index: mln/fun/x2x/composed.hh --- mln/fun/x2x/composed.hh (revision 3212) +++ mln/fun/x2x/composed.hh (working copy) @@ -87,7 +87,7 @@ /// Return the inverse function. invert inv() const; - protected: + //protected: F f_; G g_; @@ -115,7 +115,7 @@ /// Set the new second transformation. void set_second(const G& g); - protected: + //protected: F f_; G g_;