April 2008 - Olena-patches

1871: Add functional M-Watershed algorithm and tests.
by Alexandre Abraham 16 Apr '08

16 Apr '08

https://svn.lrde.epita.fr/svn/oln/trunk/milena Index: ChangeLog from Alexandre Abraham <abraham(a)lrde.epita.fr> Add functional M-Watershed algorithm and tests. * sandbox/abraham/morpho/basic_najman.hh (mln::morpho::basic_najman::lca) (mln::morpho::basic_najman::highest_fork) (mln::morpho::basic_najman::m-destructible) (mln::morpho::basic_najman::w-destructible) (mln::morpho::basic_najman::m-watershed) (mln::morpho::basic_najman::w-watershed): New methods. * sandbox/abraham/morpho/test.cc: Add debug tests. * sandbox/abraham/morpho/test_component_tree.cc: New, Add tests for component tree algorithms. * sandbox/abraham/morpho/images/result_m_watershed.pgm: New, Reference for m-watershed algorithm. * sandbox/abraham/morpho/images/result_topo_watershed.pgm: New, Reference for w-watershed algorithm. * sandbox/abraham/morpho/images/test_watershed.pgm: New, Test image for watershed algorithms. * sandbox/abraham/morpho/images/test_component_mapping.pgm: New, Reference image for component mapping. * sandbox/abraham/morpho/test_watershed.cc: New, Add tests for watershed algorithms. * sandbox/abraham/morpho/Makefile: Add rules to compile new tests. Makefile | 16 +- basic_najman.hh | 295 +++++++++++++++++++++++++++++++------- images/test_component_mapping.pgm | 5 test.cc | 15 - test_component_tree.cc | 224 ++++++++++++++++++++++++++++ test_watershed.cc | 55 +++++++ 6 files changed, 542 insertions(+), 68 deletions(-) Index: sandbox/abraham/morpho/basic_najman.hh --- sandbox/abraham/morpho/basic_najman.hh (revision 1870) +++ sandbox/abraham/morpho/basic_najman.hh (working copy) @@ -1,6 +1,10 @@ #include <mln/level/sort_psites.hh> #include <mln/level/fill.hh> #include <mln/core/image2d.hh> +#include <mln/core/p_set.hh> +//#include <mln/util/greater_psite.hh> +#include <queue> +#include <set> namespace mln { @@ -28,19 +32,19 @@ it.next()) children.append(it.to_psite()); } + void addChild(const psite n) { children.append(n); } }; // struct node - + I pima; const Image& ima; const Neighborhood<N>& nbh; mln_ch_value(I, psite) Par_node; - mln_ch_value(I, psite) Par_tree; - // Par_node is handled by par (from ap_maxtree) + mln_ch_value(I, int) Rnk_tree; mln_ch_value(I, int) Rnk_node; mln_ch_value(I, psite) subtreeRoot; @@ -52,7 +56,8 @@ basic_najman(const Image& i, const Neighborhood<N>& n) - : ima(i), + : pima(exact(i)), + ima(i), nbh(n), Par_node(exact(i).domain(), exact(i).border()), Par_tree(exact(i).domain(), exact(i).border()), @@ -73,7 +78,7 @@ void MakeSet_node(psite x) { - Par_node(x) = x; // was Par_node(x) = x; + Par_node(x) = x; Rnk_node(x) = 0; } @@ -91,24 +96,6 @@ y = memo; } -// bool is_root(const psite& x) const -// { -// return Par_node(x) == x; -// } - -// bool is_level_root(const psite& x) const -// { -// return is_root(x) || exact(ima)(Par_node(x)) < exact(ima)(x); -// } - -// psite find_level_root(const psite& x) -// { -// if (is_level_root(x)) -// return x; -// else -// return Par_node(x) = find_level_root(Par_node(x)); -// } - psite Find_node(psite x) { if (Par_node(x) != x) @@ -148,7 +135,6 @@ for (int ip = 0; ip < int(S.npoints()); ++ip) { psite p = S[ip]; - // isproc(p) = true; psite curTree = Find_tree(p); psite curNode = Find_node(subtreeRoot(curTree)); @@ -162,36 +148,22 @@ psite adjNode = Find_node(subtreeRoot(adjTree)); if (curNode != adjNode) { - std::cout << "curNode != adjNode; " << nodes(curNode).level << " vs. " << nodes(adjNode).level << std::endl; if (nodes(curNode).level == nodes(adjNode).level) { - // curNode = MergeNode(adjNode, curNode); - psite tmpNode = MergeNode(adjNode, curNode); - /*if (tmpNode == curNode) - nodes(curNode).addChildren(nodes(adjNode)); - else - nodes(adjNode).addChildren(nodes(curNode)); - nodes(adjNode) = nodes(curNode);*/ - curNode = tmpNode; + curNode = MergeNode(adjNode, curNode); } else { - // we have nodes[curNode].level < nodes[adjNode].level nodes(curNode).addChild(adjNode); - //apparemment NON :Par_node[adjNode] = curNode; // car adjNode devient fils de curNode nodes(curNode).area += nodes(adjNode).area; nodes(curNode).highest += nodes(adjNode).highest; } - // curTree = Link_tree(adjTree, curTree); - // subtreeRoot(curTree) = curNode; - } curTree = Link_tree(adjTree, curTree); subtreeRoot(curTree) = curNode; } isproc(p) = true; } - Root = subtreeRoot(Find_tree(Find_node(psite(0, 0)))); // Pour garder une map de correspondance point <-> local_root // for (int ip = 0; ip < int(S.size()); ++ip) // { @@ -203,6 +175,7 @@ for_all(r) Par_node(r) = Find_node(r); + Root = subtreeRoot(Find_tree(Find_node(psite(0, 0)))); } @@ -243,7 +216,7 @@ else if (Rnk_node(x) == Rnk_node(y)) Rnk_node(y) += 1; - Par_node(x) = y; // was Par_node(x) = y; + Par_node(x) = y; return y; } @@ -275,40 +248,254 @@ void print_tree(psite p) { node& n = nodes(p); - std::cout << "psite " << p << "=" << (p.row() * exact(subtreeRoot).domain().len(1) + p.col()) << " : "; + std::cout << "psite " <::fwd_piter it(n.children); - for (it.start(); - it.is_valid(); - it.next()) + for_all(it) { psite q = it.to_psite(); std::cout << q << "=" << (q.row() * subtreeRoot.domain().len(1) + q.col()) << " "; } std::cout << std::endl; - for (it.start(); - it.is_valid(); - it.next()) + for_all(it) { print_tree(it.to_psite()); } } - void print_sup_tree() + psite lca (psite a, psite b) { - psite max = psite(0, 0); + return lca_rec(Root, Par_node(a), Par_node(b)); + } - for(unsigned int i = 0; i < nodes.domain().len(0); ++i) - for(unsigned int j = 0; j < nodes.domain().len(1); ++j) + psite lca_rec (psite cur, psite a, psite b) + { + // FIXME : naive implementation, make something better + + if (cur == a) + return a; + + if (cur == b) + return b; + + if (nodes(cur).children.npoints() == 0) + return psite (-1, -1); + + psite tmp, acc = psite(-1, -1); + int n = 0; + typename p_array<mln_psite(I)>::fwd_piter it(nodes(cur).children); + for_all(it) + { + tmp = lca_rec(it.to_psite(), a, b); + if (tmp != psite(-1, -1)) + { + if (tmp == a) { - if (nodes(psite(i, j)).area > nodes(max).area) - max = psite(i, j); + acc = tmp; + n++; + continue; } + if (tmp == b) + { + acc = tmp; + n++; + continue; + } + return tmp; + } + } + if (n == 2) + return cur; + + return acc; + } + + psite min (p_set<psite>& components) + { + if (components.npoints() == 0) + return psite(-1, -1); + + typename p_set<psite>::fwd_piter it(components); + psite min = components[0]; + + for_all(it) + if (nodes(it.to_point()).level > nodes(min).level) + min = it.to_point(); + + return min; + } + + psite highest_fork (p_set<psite>& components) + { + if (components.npoints() == 0) + return psite(-1, -1); + + psite hfork = components[0], min = hfork; + typename p_set<psite>::fwd_piter it(components); - print_tree(max); + for_all(it) + hfork = lca(hfork, it.to_point()); + + if (nodes(min).level == nodes(hfork).level) + return psite(-1, -1); + + return hfork; } + psite w_destructible(psite p) + { + mln_niter(N) q(nbh, p); + p_set<psite> v; + + for_all(q) + if (exact(ima)(q) < exact(ima)(p)) + v.append(Par_node(q)); + + if (v.npoints() == 0) + return psite(-1, -1); + + psite hf = highest_fork(v); + + if (hf == psite(-1, -1)) + return min(v); + + if (nodes(hf).level <= exact(ima)(p)) + return hf; + + return psite(-1, -1); + } + + psite m_destructible(psite p) + { + mln_niter(N) q(nbh, p); + p_set<psite> v; + + for_all(q) + if (exact(ima)(q) < exact(ima)(p)) + v.insert(Par_node(q)); + + if (v.npoints() == 0) + return psite(-1, -1); + + if (nodes(min(v)).children.npoints() != 0) + return (psite(-1, -1)); + + psite hf = highest_fork(v); + + if (hf == psite(-1, -1)) + return min(v); + + return psite(-1, -1); + } + + void m_watershed () + { + + // FIXME : make a better priority function + // typedef + // std::priority_queue< psite, std::vector<psite>, util::greater_psite > + // ordered_queue_type; + + + std::priority_queue<psite> l; + + // Clear the marker map + level::fill(isproc, false); + mln_piter(I) it(exact(ima).domain()); + + for_all(it) + if (m_destructible(it.to_point()) != psite(-1, -1)) + { + l.push(it.to_point()); + isproc(it.to_point()) = true; + } + + psite p, i; + while (!l.empty()) + { + p = l.top(); + l.pop(); + isproc(p) = false; + + i = m_destructible(p); + + if (i != psite(-1, -1)) + { + pima(p) = nodes(i).level ; + Par_node(p) = i; + mln_niter(N) q(nbh, p); + + for_all(q) + if (!isproc(q)) + if (m_destructible(q) != psite(-1, -1)) + { + l.push(q); + isproc(q) = true; + } + } + } + } + + void w_watershed() + { + p_array< p_set<psite> > L; + // TODO : replace int with the type of value + mln_ch_value(I, int) K; + mln_ch_value(I, psite) H; + + mln_piter(I) it(exact(ima).domain()); + + psite i; + for_all(it) + { + psite p = it.to_point(); + + i = w_destructible(p); + if (i != psite(-1, -1)) + { + L[nodes(i).level].insert(p); + K(p) = nodes(i).level; + H(p) = i; + + typename p_array< p_set<psite> >::fwd_piter n(L); + for_all(n) + { + while (!n.empty()) + { + psite p = *n.begin(); + n.erase(p); + // TODO : replace int with the type of value + int k = nodes(p).level + 1; + if (K(p) == k) + { + pima(p) = k; + Par_node(p) = H(p); + + mln_niter(N) q(nbh, p); + for_all(q) + if (k < pima(q)) + { + psite j = w_destructible(q); + if (j != psite(-1, -1)) + K(q) = -1; + else + if (K(q) != nodes(j).level) + { + n.insert(q); + K(q) = i-1; + H(q) = j; + } + } + } + } + } + } + } + + } + + }; // struct basic_najman }; // namespace morpho Index: sandbox/abraham/morpho/test.cc --- sandbox/abraham/morpho/test.cc (revision 1870) +++ sandbox/abraham/morpho/test.cc (working copy) @@ -17,21 +17,13 @@ using namespace mln; using value::int_u8; - // component_tree< image2d<int_u8> > c; - using namespace mln; using value::int_u8; image2d<int_u8> input; io::pgm::load(input, "./images/test_component_tree.pgm"); - - for(unsigned int i = 0; i < input.domain().len(0); ++i) - for(unsigned int j = 0; j < input.domain().len(1); ++j) - input.at(i, j) = 255 - input.at(i, j); - - - morpho::basic_najman< image2d<int_u8>, neighb2d> n(input, c8()); + morpho::basic_najman< image2d<int_u8>, neighb2d> n(input, c4()); n.image_output(input); n.go(); n.image_output(n.Par_tree); @@ -42,7 +34,10 @@ std::cout << " --------------------------- " << std::endl; n.image_output(n.Rnk_node); std::cout << " --------------------------- " << std::endl; - n.print_sup_tree(); + n.print_tree(n.Root); + std::cout << " --------------------------- " << std::endl; + + n.m_watershed(); // image2d<int_u8> output = morpho::topo_wst(input, c4()); // io::pgm::save(output, "out.pgm"); Index: sandbox/abraham/morpho/test_component_tree.cc --- sandbox/abraham/morpho/test_component_tree.cc (revision 0) +++ sandbox/abraham/morpho/test_component_tree.cc (revision 0) @@ -0,0 +1,224 @@ +#include <mln/core/image2d.hh> +#include <mln/core/window2d.hh> +#include <mln/core/neighb2d.hh> +#include <mln/core/p_set.hh> + +#include <mln/value/int_u8.hh> + +#include <mln/io/pgm/load.hh> +#include <mln/io/pgm/save.hh> + +#include "basic_najman.hh" +#include <string> +#include <iostream> + +int print_and_exit (std::string s) +{ + std::cerr << s << std::endl; + return 1; +} + +template <typename T> +void swap (T& x, T& y) +{ + T memo = x; + x = y; + y = memo; +} + +template <typename T> +class test_tree; + +template <typename T> +class test_tree +{ +public: + + test_tree(T pdata) : data(pdata), children(NULL), brother(NULL) { } + + T data; + test_tree<T> *children; + test_tree<T> *brother; + + test_tree<T>* insert_child(test_tree<T> *n) + { + test_tree<T> **p = &children; + while (*p && (**p).data < n->data) + p = &((**p).brother); + n->brother = *p; + *p = n; + return n; + } + + test_tree<T>* insert_child(T pdata) + { + test_tree<T> *n = new test_tree<T>(pdata); + return insert_child(n); + } +}; + +template <typename T> +bool compare_tree(test_tree<T>* t1, test_tree<T>* t2) +{ + if (!t1 && !t2) + return true; + + if ((!t1 || !t2) || (t1->data != t2->data)) + return false; + + return (compare_tree(t1->brother, t2->brother) + && compare_tree(t1->children, t2->children)); +} + +template <typename T> +void show_tree (test_tree<T> *t) +{ + if (!t) + return ; + + if (t->brother) + { + std::cout << t->data << " - "; + show_tree(t->brother); + } + else + std::cout << t->data << std::endl; + + if (t->children) + { + std::cout << "sons of " << t->data << " : "; + show_tree(t->children); + } +} + + +using namespace mln; +using value::int_u8; + + +test_tree<image2d<int_u8>::psite>* test_convert (morpho::basic_najman<image2d<int_u8>, neighb2d>& ima, image2d<int_u8>::psite root); + +test_tree<image2d<int_u8>::psite>* test_convert (morpho::basic_najman<image2d<int_u8>, neighb2d>& ima, image2d<int_u8>::psite root) +{ + morpho::basic_najman<image2d<int_u8>, neighb2d>::node& n = ima.nodes(ima.Par_node(root)); + + test_tree<image2d<int_u8>::psite> *res = new test_tree<image2d<int_u8>::psite>(ima.Par_node(root)); + + p_array<image2d<int_u8>::psite>::fwd_piter it (n.children); + for_all(it) + res->insert_child(test_convert(ima, it.to_psite())); + + return res; +} + +int main () +{ + typedef image2d<int_u8> image2dint; + + image2dint input, verif; + io::pgm::load(input, "./images/test_component_tree.pgm"); + io::pgm::load(verif, "./images/test_component_mapping.pgm"); + + morpho::basic_najman< image2d<int_u8>, neighb2d> n(input, c4()); + n.go(); + + /* + Component mapping : + + 3 3 3 3 9 2 2 2 2 + 1 1 1 9 7 9 2 2 2 + 1 1 9 7 7 7 9 2 2 + 1 9 6 9 7 9 5 9 2 + 1 1 9 1 7 2 7 2 2 + 1 1 1 1 7 2 2 2 2 + + + Component tree + + 9 + / \ + 7 6 + /|\ + 3 2 5 + | + 1 + */ + + // Component mapping test + + image2dint::fwd_piter it(input.domain()); + image2dint::psite c[10]; + c[1] = n.Par_node(image2dint::psite(1, 0)); + c[2] = n.Par_node(image2dint::psite(0, 5)); + c[3] = n.Par_node(image2dint::psite(0, 0)); + c[5] = n.Par_node(image2dint::psite(3, 6)); + c[6] = n.Par_node(image2dint::psite(3, 2)); + c[7] = n.Par_node(image2dint::psite(1, 4)); + c[9] = n.Par_node(image2dint::psite(0, 4)); + + bool id = true; + + for_all(it) + id = (c[verif(it.to_point())] == n.Par_node(it.to_point())); + + if (!id) + { + std::cerr << "Component mapping is fucked up!" << std::endl; + return 1; + } + + // Component tree test + + test_tree<image2dint::psite> *gen, ref(c[9]), *l1, *l2; + + // Create the ref + ref.insert_child(c[6]); + l1 = ref.insert_child(c[7]); + l1->insert_child(c[2]); + l1->insert_child(c[5]); + l2 = l1->insert_child(c[3]); + l2->insert_child(c[1]); + + // Create the generated tree + gen = test_convert(n, n.Par_node(n.Root)); + + // Compare the trees + if (!compare_tree(gen, &ref)) + { + std::cout << "Your tree :" << std::endl; + show_tree(gen); + + std::cout << "Reference tree :" << std::endl; + show_tree(&ref); + + std::cerr << "Component tree is fucked up!" << std::endl; + return 1; + } + + // Test the LCA alorithm + image2dint::psite p(0, 0), q(0, 6), r(2, 4); + if (n.lca(p, q) != image2dint::psite(2,4)) + return 1; + + if (n.lca(r, q) != image2dint::psite(2,4)) + return 1; + + // Test the Highest Fork + p_set<image2dint::psite> comp; + comp.insert(image2dint::psite(0,0)); + comp.insert(image2dint::psite(2,4)); + comp.insert(image2dint::psite(0,6)); + comp.insert(image2dint::psite(1,0)); + + if (n.highest_fork(comp) != image2dint::psite(2, 4)) + return 1; + + comp.insert(image2dint::psite(3,2)); + + if (n.highest_fork(comp) != image2dint::psite(1, 3)) + return 1; + + std::cout << "Test is successfull" << std::endl; + + return 0; +} Index: sandbox/abraham/morpho/images/result_m_watershed.pgm Cannot display: file marked as a binary type. svn:mime-type = application/octet-stream Property changes on: sandbox/abraham/morpho/images/result_m_watershed.pgm ___________________________________________________________________ Name: svn:mime-type + application/octet-stream Index: sandbox/abraham/morpho/images/result_topo_watershed.pgm Cannot display: file marked as a binary type. svn:mime-type = application/octet-stream Property changes on: sandbox/abraham/morpho/images/result_topo_watershed.pgm ___________________________________________________________________ Name: svn:mime-type + application/octet-stream Index: sandbox/abraham/morpho/images/test_watershed.pgm Cannot display: file marked as a binary type. svn:mime-type = application/octet-stream Property changes on: sandbox/abraham/morpho/images/test_watershed.pgm ___________________________________________________________________ Name: svn:mime-type + application/octet-stream Index: sandbox/abraham/morpho/images/test_component_mapping.pgm --- sandbox/abraham/morpho/images/test_component_mapping.pgm (revision 0) +++ sandbox/abraham/morpho/images/test_component_mapping.pgm (revision 0) @@ -0,0 +1,5 @@ +P5 +# CREATOR: GIMP PNM Filter Version 1.1 +9 6 +255 + \ No newline at end of file Index: sandbox/abraham/morpho/test_watershed.cc --- sandbox/abraham/morpho/test_watershed.cc (revision 0) +++ sandbox/abraham/morpho/test_watershed.cc (revision 0) @@ -0,0 +1,55 @@ +#include <mln/core/image2d.hh> +#include <mln/core/window2d.hh> +#include <mln/core/neighb2d.hh> +#include <mln/core/p_set.hh> + +#include <mln/value/int_u8.hh> +#include <mln/level/compare.hh> + +#include <mln/io/pgm/load.hh> +#include <mln/io/pgm/save.hh> + +#include "basic_najman.hh" +#include <string> +#include <iostream> + +int print_and_exit (std::string s) +{ + std::cerr << s << std::endl; + return 1; +} + +int main () +{ + using namespace mln; + using value::int_u8; + + typedef image2d<int_u8> image2dint; + + image2dint input, mverif; + io::pgm::load(input, "./images/test_watershed.pgm"); + // io::pgm::load(input, "./images/lena.pgm"); + io::pgm::load(mverif, "./images/result_m_watershed.pgm"); + + morpho::basic_najman< image2d<int_u8>, neighb2d> n(input, c4()); + n.go(); + n.m_watershed(); + + bool equal = true; + + equal = n.pima == mverif; + + if (!equal) + { + std::cout << "M-Watershed broken" << std::endl; + return 1; + } + else + std::cout << "M-Watershed good as chocolate ice cream" << std::endl; + + // n.w_watershed(); + + io::pgm::save(n.pima, "out.pgm"); + + return 0; +} Index: sandbox/abraham/morpho/Makefile --- sandbox/abraham/morpho/Makefile (revision 1870) +++ sandbox/abraham/morpho/Makefile (working copy) @@ -1,8 +1,16 @@ -SOURCE=test.cc -OBJECTS=test.o CXXFLAGS=-Wall -W -I ../../.. -g -all:${OBJECTS} - g++ ${FLAGS} ${OBJECTS} -o test +test: test.o + g++ ${CXXFLAGS} test.o -o test +test_component_tree: test_component_tree.o + g++ ${CXXFLAGS} test_component_tree.o -o test_component_tree + +test_watershed: test_watershed.o + g++ ${CXXFLAGS} test_watershed.o -o test_watershed + +all: tree + +test_watershed.o: basic_najman.hh +test_component_tree.o: basic_najman.hh test.o: basic_najman.hh \ No newline at end of file

1 0

milena r1870: Commit my last version of fllt
by Matthieu Garrigues 16 Apr '08

16 Apr '08

URL: https://svn.lrde.epita.fr/svn/oln/trunk/milena ChangeLog: 2008-04-16 Matthieu Garrigues <garrigues(a)lrde.epita.fr> Commit my last version of fllt. * sandbox/garrigues/fllt/compute_level_set_fast.hh: New. * sandbox/garrigues/fllt/compute_level_set_fast2.hh: New. * sandbox/garrigues/fllt/debug.hh: . * sandbox/garrigues/fllt/doc.hh: . * sandbox/garrigues/fllt/fllt.hh: . * sandbox/garrigues/fllt/give_confs.cc: New. * sandbox/garrigues/fllt/local_configurations.hh: New. * sandbox/garrigues/fllt/lower.hh: . * sandbox/garrigues/fllt/test.cc: New. * sandbox/garrigues/fllt/test_fllt2.cc: . * sandbox/garrigues/fllt/types.hh: . * sandbox/garrigues/fllt/upper.hh: . --- compute_level_set_fast.hh | 486 +++++++++++++++++++++++++++++++++++++++++++++ compute_level_set_fast2.hh | 471 +++++++++++++++++++++++++++++++++++++++++++ debug.hh | 83 +++++++ doc.hh | 1 fllt.hh | 49 ---- give_confs.cc | 56 +++++ local_configurations.hh | 144 +++++++++++++ lower.hh | 4 test.cc | 62 +++++ test_fllt2.cc | 4 types.hh | 150 +++++++++++++ upper.hh | 4 12 files changed, 1460 insertions(+), 54 deletions(-) Index: trunk/milena/sandbox/garrigues/fllt/lower.hh =================================================================== --- trunk/milena/sandbox/garrigues/fllt/lower.hh (revision 1869) +++ trunk/milena/sandbox/garrigues/fllt/lower.hh (revision 1870) @@ -36,6 +36,7 @@ */ # include <mln/core/neighb2d.hh> +# include <mln/core/clock_neighb2d.hh> # include <mln/labeling/regional_minima.hh> # include <mln/accu/min.hh> @@ -75,6 +76,9 @@ static const neighb2d& bdr_nbh() { return c8(); } static const neighb2d& reg_nbh() { return c4(); } + static const clock_neighb2d bdr_c_nbh(dpoint2d& dp) { return cc8(dp); } + static const clock_neighb2d reg_c_nbh(dpoint2d& dp) { return cc4(dp); } + }; } // end of namespace mln::fllt Index: trunk/milena/sandbox/garrigues/fllt/upper.hh =================================================================== --- trunk/milena/sandbox/garrigues/fllt/upper.hh (revision 1869) +++ trunk/milena/sandbox/garrigues/fllt/upper.hh (revision 1870) @@ -36,6 +36,7 @@ */ # include <mln/core/neighb2d.hh> +# include <mln/core/clock_neighb2d.hh> # include <mln/labeling/regional_maxima.hh> # include <mln/accu/max.hh> @@ -74,6 +75,9 @@ static const neighb2d& bdr_nbh() { return c4(); } static const neighb2d& reg_nbh() { return c8(); } + + static const clock_neighb2d bdr_c_nbh(dpoint2d& dp) { return cc4(dp); } + static const clock_neighb2d reg_c_nbh(dpoint2d& dp) { return cc8(dp); } }; } // end of namespace mln::fllt Index: trunk/milena/sandbox/garrigues/fllt/types.hh =================================================================== --- trunk/milena/sandbox/garrigues/fllt/types.hh (revision 1869) +++ trunk/milena/sandbox/garrigues/fllt/types.hh (revision 1870) @@ -39,6 +39,12 @@ # include <mln/util/tree.hh> # include <mln/util/branch_iter_ind.hh> + +# define fllt_tree(P, V) mln::util::tree< mln::fllt::fllt_node_elt<P, V> > +# define fllt_node(P, V) mln::util::tree_node< mln::fllt::fllt_node_elt<P, V> > +# define fllt_branch(P, V) mln::util::branch< mln::fllt::fllt_node_elt<P, V> > +# define fllt_branch_iter_ind(P, V) mln::util::branch_iter_ind< mln::fllt::fllt_node_elt<P, V> > + namespace mln { namespace fllt @@ -55,13 +61,147 @@ bool brighter; }; -# define fllt_tree(P, V) util::tree< fllt_node_elt<P, V> > -# define fllt_node(P, V) util::tree_node< fllt_node_elt<P, V> > -# define fllt_branch(P, V) util::branch< fllt_node_elt<P, V> > -# define fllt_branch_iter_ind(P, V) util::branch_iter_ind< fllt_node_elt<P, V> > - // # define fllt_node(P, V) typename fllt_tree(P, V)::node_t + + class ran_domains + { + public: + + /// Forward Point_Iterator associated type. + typedef mln_fwd_piter_(box2d) fwd_piter; + + /// Backward Point_Iterator associated type. + typedef mln_bkd_piter_(box2d) bkd_piter; + + /// Constructor. + ran_domains(const box2d& b); + + /// Add a point to a domain. + template <char domain> + ran_domains& add_to(const point2d& p); + + /// Test if a point belong to a domain. + template <char domain> + bool belongs_to(const point2d& p) const; + + /// Move a point from a domain to another domain. + template <char src, char dest> + ran_domains& move_to(const point2d& p); + + /// Clear the image. + void clear(); + + /// Give the exact bounding box. + const box2d& bbox() const; + + /// Give the number of points. + unsigned npoints() const; + + /// Hook to the image2d containing the points. + sub_image<image2d<value::int_u8>, box2d> image(); + + private: + image2d<value::int_u8> ima_; + unsigned npoints_; + accu::bbox<point2d> bb_; + }; + + +# ifndef MLN_INCLUDE_ONLY + + inline + ran_domains::ran_domains(const box2d& b) + : ima_(b) + { + bb_.init(); + npoints_ = 0; + + level::fill(ima_, false); + } + + template <char domain> + inline + ran_domains& + ran_domains::add_to(const point2d& p) + { + bb_.take(p); + ima_(p) = domain; + npoints_++; + return *this; + } + + + template <char src, char dest> + ran_domains& + ran_domains::move_to(const point2d& p) + { + mln_assertion(ima_(p) == src); + ima_(p) += dest - src; + return *this; + } + + template <char domain> + inline + bool + ran_domains::belongs_to(const point2d& p) const + { + return ima_(p) == domain; + } + + + void + ran_domains::clear() + { + if (npoints_ == 0) + return; + + // unsigned bb_nrows = geom::nrows(bb_.to_result()); + // unsigned ima_nrows = geom::nrows(points_); + + // if (bb_nrows * 3 < ima_nrows * 2) + // { + // unsigned bb_ncols = geom::ncols(bb_.to_result()); + // mln_line_piter_(image2d<value::int_u8>) p(bb_.to_result()); + // for_all(p) + // { + // level::memset_(ima_, p, false, bb_ncols); + // } + // } + // else + level::fill(ima_, false); + + npoints_ = 0; + bb_.init(); + } + + + inline + const box2d& + ran_domains::bbox() const + { + mln_precondition(npoints_ != 0); + return bb_.to_result(); + } + + inline + unsigned + ran_domains::npoints() const + { + return npoints_; + } + + inline + sub_image<image2d<value::int_u8>, box2d> + ran_domains::image() + { + mln_precondition(npoints_ > 0); + mln_assertion(ima_.has_data()); + return ima_ | bb_.to_result(); + } + +# endif // ! MLN_INCLUDE_ONLY + } // end of namespace mln::fllt } // end of namespace mln Index: trunk/milena/sandbox/garrigues/fllt/doc.hh =================================================================== --- trunk/milena/sandbox/garrigues/fllt/doc.hh (revision 1869) +++ trunk/milena/sandbox/garrigues/fllt/doc.hh (revision 1870) @@ -58,6 +58,7 @@ // gn <- min u(x) x belongs to N. // R <- R union A // tag the pixels of A. + // Update the number of holes of R. // 4) // IF g < gn Index: trunk/milena/sandbox/garrigues/fllt/compute_level_set_fast2.hh =================================================================== --- trunk/milena/sandbox/garrigues/fllt/compute_level_set_fast2.hh (revision 0) +++ trunk/milena/sandbox/garrigues/fllt/compute_level_set_fast2.hh (revision 1870) @@ -0,0 +1,471 @@ +// Copyright (C) 2007, 2008 EPITA Research and Development Laboratory +// +// This file is part of the Olena Library. This library is free +// software; you can redistribute it and/or modify it under the terms +// of the GNU General Public License version 2 as published by the +// Free Software Foundation. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// General Public License for more details. +// +// You should have received a copy of the GNU General Public License +// along with this library; see the file COPYING. If not, write to +// the Free Software Foundation, 51 Franklin Street, Fifth Floor, +// Boston, MA 02111-1307, USA. +// +// As a special exception, you may use this file as part of a free +// software library without restriction. Specifically, if other files +// instantiate templates or use macros or inline functions from this +// file, or you compile this file and link it with other files to +// produce an executable, this file does not by itself cause the +// resulting executable to be covered by the GNU General Public +// License. This exception does not however invalidate any other +// reasons why the executable file might be covered by the GNU General +// Public License. + + +#ifndef MLN_FIXME_FLLT_COMPUTE_LEVEL_SET_FAST_HH +# define MLN_FIXME_FLLT_COMPUTE_LEVEL_SET_FAST_HH + +/*! \file compute_level_set_fast2.hh + * + * \brief Compute Level Set for Fast level line transform. + * Version with one image for A, R, N. + * + */ + +# include <mln/core/p_image2d.hh> +# include <mln/core/image_if_value.hh> +# include "local_configurations.hh" + +# define SET_R 3 +# define SET_A 2 +# define SET_N 1 + +# define image_sub_if_value image_if_value<const sub_image<image2d<value::int_u8>, box2d> > + +namespace mln +{ + namespace fllt_fast + { + using namespace fllt; + + class fllt::ran_domains; + + template <typename V> + void step_fast1 (const image2d<V>& ima, + point2d p, + V& g, + point2d& x0) + { + //std::cout << "entering step_fast 1" << std::endl; + // x0 <- a not tagged local mininum of ima. + //std::cout << std::endl << "x0 = " << p << std::endl; + x0 = p; + // g <- u(x0) + g = ima(x0); + //std::cout << "g = " << g << std::endl; + //std::cout << "exiting step_fast 1" << std::endl; + } + + void step_fast2 (ran_domains& r_a_n, + point2d& x0) + { + //std::cout << "entering step_fast 2" << std::endl; + // R <- {} + // N <- {} + // A <- {x0} + r_a_n.clear(); + r_a_n.add_to<SET_A>(x0); + + + //std::cout << "exiting step_fast 2" << std::endl; + } + + + template <typename V, typename F> + void step_fast3 (const image2d<V>& u, + image2d<bool>& tagged, + ran_domains& r_a_n, + V& gn, + fllt_node(point2d, V)*& current_region, + image2d<fllt_node(point2d, V)*>& regions) + { + static bool finished = false; + //std::cout << "entering step_fast 3" << std::endl; + + // Stop the algorithm. + if (finished) + { finished = false; gn -= 2 * F::inc; return; } + + // N <- N union {x neighbor of a pixel in a\R} + //mln_piter(p_image2d) qa(A); + + // p_image2d_fwd_pixter<point2d> qa(A); + image_sub_if_value ima(r_a_n.image() | SET_A); + mln_assertion(ima.has_data()); + image_sub_if_value::fwd_piter qa(ima.domain()); + mln_assertion(ima.has_data()); + for_all(qa) + { + mln_assertion(ima.has_data()); + mln_niter(neighb2d) n(F::reg_nbh(), qa); + for_all (n) + if (u.has(n) && !r_a_n.belongs_to<SET_R>(n) && + !r_a_n.belongs_to<SET_A>(n)) + r_a_n.add_to<SET_N>(n); + } + + // debug::println(u); + +// std::cout << "A :" << A << std::endl; +// debug::println(A.image()); +// std::cout << "N :" << N << std::endl; +// debug::println(N.image()); +// std::cout << "R :" << R << std::endl; +// debug::println(R.image()); + + // gn <- min u(x) x belongs to N. + + image_sub_if_value tmp(r_a_n.image() | SET_N); + finished = tmp.npoints() == 0; + + if (!finished) + gn = level::compute< typename F::accu_for_gn >(u | tmp.domain()); + + if (finished) + gn += F::inc; + + //std::cout << std::endl << "gN = " << gn << std::endl; + // R <- R union A + // tag the pixels of A. + + for_all(qa) + { + if (!r_a_n.belongs_to<SET_A>(qa)) + continue; +// std::cout << "Added to R" << tmp << std::endl; + mln_assertion(!r_a_n.belongs_to<SET_R>(qa)); + r_a_n.move_to<SET_A, SET_R>(qa); + tagged(qa) = true; + // Update the number of holes of R. + //n_added_holes<point2d, F>(R, qa); + } + + + //std::cout << "exiting step_fast 3" << std::endl; + } + + + template <typename V> + void update_set(const image2d<V>& u, + ran_domains& r_a_n, + V& g) + { + // A <- {x belongs to N / u(x) == g} + + // This is slow. + // A.insert(set::inter(N, u.domain()) | pw::value(u) == pw::cst(g)); + + // N <- N\{x belongs to N / u(x) == g} + // N.remove(set::inter(N, u.domain()) | pw::value(u) == pw::cst(g)); + // mln_assertion(check_step_4_2(u, N, g)); + + // This is faster. + +// if (N.npoints() >= 0) +// { + image_sub_if_value tmp(r_a_n.image() | SET_N); + image_sub_if_value::fwd_piter qn(tmp.domain()); + for_all(qn) + { + if (u(qn) == g) + r_a_n.move_to<SET_N, SET_A>(qn); + } +// } +// else { +// mln_fwd_piter(p_image2d) p(N); +// for_all(p) +// { +// if (u(p) == g) +// { +// A.insert(p); +// N.remove(p); +// } +// } +// mln_assertion(check_step_4_2(u, N, g)); +// } + } + /// IF g < gn. + template <typename V, typename F> + void step_fast4_1 (image2d<V>& u, + ran_domains& r_a_n, + V& g, + V& gn, + fllt_node(point2d, V)*& current_region, + image2d<fllt_node(point2d, V)*>& regions + ) + { +// std::cout << "entering step_fast 4_1" << std::endl; + + // If the region is bounded + // Create a new conected component. + // FIXME : we can make it faster. + + // FIXME : rewrite this test. + // if ((R.bbox() < u.domain()) || (R.npoints() == u.domain().npoints())) +// { + image_sub_if_value tmp(r_a_n.image() | SET_R); + image_sub_if_value::fwd_piter p(tmp.domain()); + current_region = new fllt_node(point2d, V)(); + current_region->elt().brighter = F::parent_is_brighter; + current_region->elt().value = g; + for_all(p) + { + current_region->elt().points.insert(p); + + if (regions(p) == 0) + { + //current_region->elt().points.insert(p); + regions(p) = current_region; + } + else + { + if (regions(p)->parent() == 0) + regions(p)->set_parent(current_region); + } + } + + +// // Count the number of conected components of the border of R. +// static image2d<unsigned> tmp(u.domain().to_larger(1)); +// static image2d<bool> border_ima(tmp.domain()); +// level::fill(border_ima, false); + +// // level::fill(inplace(border_ima | N), true); +// // std::cout << "tmp border = " << tmp.border () << std::endl; +// // std::cout << "ima border = " << border_ima.border () << std::endl; +// mln_piter(p_image2d) z(N); +// for_all(z) +// { +// mln_assertion(border_ima.owns_(z)); +// border_ima(z) = true; +// } +// unsigned n; +// tmp = labeling::level(border_ima, true, F::bdr_nbh(), n); + +// // debug::println(border_ima); +// //std::cout << "nb composantes :" << n << std::endl; +// // debug::println(tmp); +// if (n > 1) +// { + +// // IF number of conected components of the border > 1 +// for (int i = 2; i <= n; i++) +// { +// // follow each border to find which is the exterior border +// // and which are the holes. Keep one pixel of each holes. + +// // WARNING : We trust labeling::level to label the exterior border with 1. +// current_region->elt().holes.insert(a_point_of(tmp | pw::value(tmp) == pw::cst(i))); + +// // FIXME : [optimisation] Remove from N border of holes???. +// // Recompute gn <- min u(x) x belongs to A +// } +// } + + +// } + + g = gn; + + update_set(u, r_a_n, g); + + // std::cout << "A :" << std::endl; + // if (A.npoints()) + // debug::println(u | A); + // std::cout << "N :" << std::endl; + // if (N.npoints()) + // debug::println(u | N); + +// std::cout << "exiting step_fast 4_1" << std::endl; + } + + /// IF g == gn. + template <typename V, typename P> + void step_fast4_2 (const image2d<V>& u, + ran_domains& r_a_n, + V& g, + fllt_node(P, V)* current_region, + image2d<fllt_node(P, V)*>& regions + ) + { +// std::cout << "entering step_fast 4_2" << std::endl; + + update_set(u,r_a_n,g); + + // std::cout << "A :" << std::endl; + // if (A.npoints()) + // debug::println(u | A); + // std::cout << "N :" << std::endl; + // if (N.npoints()) + // debug::println(u | N); + +// std::cout << "exiting step_fast 4_2" << std::endl; + } + + /// IF g > gn. + template <typename V> + void step_fast4_3 (image2d<V>& u, + const image2d<bool>& tagged, + ran_domains& r_a_n, + const V& g) + { +// std::cout << "entering step_fast 4_3" << std::endl; + + // set the gray-level of the pixels of R to g. + image_sub_if_value tmp(r_a_n.image() | SET_R); + image_sub_if_value::fwd_piter p(tmp.domain()); + for_all(p) + { + mln_assertion (tagged(p)); + u (p) = g; + } +// std::cout << "exiting step_fast 4_3" << std::endl; + + } + + + template <typename V, typename F> + fllt_tree(point2d, V)& + compute_level_set_fast(const image2d<V>& ima, + image2d< fllt_node(point2d, V)* >& regions) + { + typedef point2d P; + typedef image2d<V> I; + + // FIXME: not nice. + typedef mln::image_if< + mln::image2d<unsigned>, + mln::fun::greater_p2b_expr_<mln::pw::value_<mln::image2d<unsigned> >, + mln::pw::cst_<int> > + > I_IF; + + // Check + mln_assertion(ima.domain() == regions.domain()); + + // Declarations. + ran_domains r_a_n(ima.domain()); + V g, gn; + point2d x0; + image2d<unsigned> min_locals(ima.domain()); + image2d<V> u = clone(ima); + border::fill(u, 0); + + //std::cout << "image U:" << std::endl; + // debug::println_with_border(u); + image2d<bool> tagged(ima.domain()); + fllt_node(P, V)* current_region; + + // INIT + g= 0; + gn = 0; + current_region = 0; + + level::fill(regions, 0); + level::fill(tagged, false); + + // Get the locals extremums + unsigned nlabels; + min_locals = F::regional_extremum(ima, F::reg_nbh(), nlabels); + + // debug::println(min_locals); + // debug::println(min_locals | (pw::value(min_locals) > pw::cst(0))); + + /// Algorithm. + { + // For all locals extremums + I_IF min_locals_list(min_locals | (pw::value(min_locals) > pw::cst(0))); + mln_piter(I_IF) p(min_locals_list.domain()); + int cpt = 0, snapshop = 0; + for_all(p) + { + if (tagged(p)) + continue; + + //std::cout << "boucle no " << cpt++ << std::endl; + step_fast1(ima, p, g, x0); + step_fast2(r_a_n, x0); + while (1) + { + // To have a view of the progression of the algotithm. +#ifdef FLLT_TRACE + save_state(ima, r_a_n); +#endif + //std::cout << "g = " << g << std::endl; + step_fast3<V, F>(u, tagged, r_a_n, gn, current_region, regions); + /// step_fast4. + if (F::compare(g, gn)) + { + step_fast4_1<V, F>(u, r_a_n, g, gn, current_region, regions); + // GO TO 3) + continue; + } + + + if (g == gn) + { + step_fast4_2(u, r_a_n, g, current_region, regions); + // GO TO 3) + continue; + } + + + if (!F::compare(g, gn)) + { + step_fast4_3(u, tagged, r_a_n, g); + // GO TO 1) + break; + } + } + //std::cout << "current_region = " << current_region << std::endl; + } + } // end of Algorithm + + image2d<value::int_u8> output (ima.domain ()); + fllt_tree(P, V)& tree = *new fllt_tree(P, V)(current_region); + util::tree_to_image (tree, output); + + // util::display_tree(ima, tree); + + // debug::println(output); + // std::cout << std::endl; + // debug::println(ima); + + // if (output != ima) + // { + // std::cerr << "BUG!!!" << std::endl; + // abort(); + // } + + // io::pgm::save(output, "out.pgm"); + // std::cout << "out.pgm generate" + // << std::endl; + + + // debug::println(regions); + //debug::println(ima | regions(make:defined reference to `mln::fllt::lower<mln::value::int_u<8u> >::inc':point2d(-4,-1))->elt().points); + + return (tree); + + } // end of compute_level_set_fast + + } // end of namespace mln::fllt + +} // end of namespace mln + + + +#endif // ! MLN_FIXME_FLLT_COMPUTE_LEVEL_SET_FAST_HH Index: trunk/milena/sandbox/garrigues/fllt/debug.hh =================================================================== --- trunk/milena/sandbox/garrigues/fllt/debug.hh (revision 1869) +++ trunk/milena/sandbox/garrigues/fllt/debug.hh (revision 1870) @@ -35,6 +35,12 @@ * */ +# include <iomanip> +# include <iostream> +# include <sstream> +# include "mln/core/p_image2d.hh" +# include "mln/literal/all.hh" +# include "mln/io/ppm/save.hh" # include "types.hh" namespace mln @@ -113,6 +119,83 @@ } } + template <typename P, typename V> + void + debug(const image2d<V>& ima, + fllt_tree(P, V)& result_tree) + { + std::cout << "4/ Generate outputs." << std::endl; + + image2d<value::int_u8> output (ima.domain ()); + util::tree_to_image (result_tree, output); + + + // io::pgm::save(output, "out_final.pgm"); + // std::cout << "out_final.pgm generate" + // << std::endl; + + + // util::display_tree(ima, lower_tree); + draw_tree(ima, result_tree); + + // debug::println(ima); + // debug::println(output); + + // if (output != ima) + // { + // std::cerr << "BUG!!!" << std::endl; + // abort(); + // } + + image2d<value::int_u8> viz(ima.domain()); + // image2d<value::int_u8> viz2(ima.domain()); + + // visualize_deepness(viz, lower_tree); + // level::stretch(viz, viz2); + // debug::println(viz); + // debug::println(viz2); + // io::pgm::save(viz2, "fllt.pgm"); + + visualize_bounds(viz, result_tree, 200); + //debug::println(viz); + io::pgm::save(viz, "fllt_bounds_200.pgm"); + + visualize_bounds(viz, result_tree, 100); + io::pgm::save(viz, "fllt_bounds_100.pgm"); + + visualize_bounds(viz, result_tree, 50); + io::pgm::save(viz, "fllt_bounds_50.pgm"); + + visualize_bounds(viz, result_tree, 20); + io::pgm::save(viz, "fllt_bounds_20.pgm"); + + } + + template <typename P, typename V> + void save_state(const image2d<V>& ima, + const p_image2d& R, + const p_image2d& A, + const p_image2d& N) + { + static int id = 0; + std::stringstream filename; + filename << "fllt_trace_" << std::setw(5) << std::setfill('0') + << std::right << id++ << ".ppm"; + + image2d<value::rgb8> out(ima.domain()); + + level::fill(out, literal::white); + + if (R.npoints() != 0) + level::fill(inplace(out | R), literal::green); + if (A.npoints() != 0) + level::fill(inplace(out | A), literal::blue); + if (N.npoints() != 0) + level::fill(inplace(out | N), literal::red); + + io::ppm::save(out, filename.str()); + } + } // end of namespace mln::fllt } // end of namespace mln Index: trunk/milena/sandbox/garrigues/fllt/compute_level_set_fast.hh =================================================================== --- trunk/milena/sandbox/garrigues/fllt/compute_level_set_fast.hh (revision 0) +++ trunk/milena/sandbox/garrigues/fllt/compute_level_set_fast.hh (revision 1870) @@ -0,0 +1,486 @@ +// Copyright (C) 2007, 2008 EPITA Research and Development Laboratory +// +// This file is part of the Olena Library. This library is free +// software; you can redistribute it and/or modify it under the terms +// of the GNU General Public License version 2 as published by the +// Free Software Foundation. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// General Public License for more details. +// +// You should have received a copy of the GNU General Public License +// along with this library; see the file COPYING. If not, write to +// the Free Software Foundation, 51 Franklin Street, Fifth Floor, +// Boston, MA 02111-1307, USA. +// +// As a special exception, you may use this file as part of a free +// software library without restriction. Specifically, if other files +// instantiate templates or use macros or inline functions from this +// file, or you compile this file and link it with other files to +// produce an executable, this file does not by itself cause the +// resulting executable to be covered by the GNU General Public +// License. This exception does not however invalidate any other +// reasons why the executable file might be covered by the GNU General +// Public License. + + +#ifndef MLN_FIXME_FLLT_COMPUTE_LEVEL_SET_FAST_HH +# define MLN_FIXME_FLLT_COMPUTE_LEVEL_SET_FAST_HH + +/*! \file compute_level_set_fast.hh + * + * \brief Compute Level Set for Fast level line transform. + * Version with three p_image2d for A, R, N. + * + */ + +# include <mln/core/p_image2d.hh> +# include "local_configurations.hh" + +namespace mln +{ + namespace fllt_fast + { + using namespace fllt; + + template <typename V> + void step_fast1 (const image2d<V>& ima, + point2d p, + V& g, + point2d& x0) + { + //std::cout << "entering step_fast 1" << std::endl; + // x0 <- a not tagged local mininum of ima. + //std::cout << std::endl << "x0 = " < + void step_fast2 (p_image2d& A, + p_image2d& R, + p_image2d& N, + point2d& x0) + { + //std::cout << "entering step_fast 2" << std::endl; + // A <- {x0} + A.clear(); + A.insert(x0); + // R <- {} + R.clear(); + // N <- {} + N.clear(); + //std::cout << "exiting step_fast 2" << std::endl; + } + + + template <typename V, typename P, typename F> + void step_fast3 (const image2d<V>& u, + image2d<bool>& tagged, + p_image2d& A, + p_image2d& R, + p_image2d& N, + V& gn, + fllt_node(P, V)*& current_region, + image2d<fllt_node(P, V)*>& regions) + { + static bool finished = false; + //std::cout << "entering step_fast 3" << std::endl; + + // Stop the algorithm. + if (finished) + { finished = false; gn -= 2 * F::inc; return; } + + // N <- N union {x neighbor of a pixel in a\R} + //mln_piter(p_image2d) qa(A); + + // p_image2d_fwd_pixter<point2d> qa(A); + p_image2d_fwd_piter_<point2d> qa(A); + for_all(qa) + { + mln_niter(neighb2d) n(F::reg_nbh(), qa); + for_all (n) + if (u.has(n) && !R.has (n) && !A.has (n)) + N.insert (n); + } + + // debug::println(u); + +// std::cout << "A :" << A << std::endl; +// debug::println(A.image()); +// std::cout << "N :" << N << std::endl; +// debug::println(N.image()); +// std::cout << "R :" << R << std::endl; +// debug::println(R.image()); + + // gn <- min u(x) x belongs to N. + + finished = N.npoints() == 0; + + if (!finished) + gn = level::compute< typename F::accu_for_gn >(u | N); + else + gn += F::inc; + + //std::cout << std::endl << "gN = " << gn << std::endl; + // R <- R union A + // tag the pixels of A. + + for_all(qa) + { + point2d tmp = qa; +// std::cout << "Added to R" << tmp << std::endl; + mln_assertion(!R.has(tmp)); + R.insert(tmp); + tagged(tmp) = true; + // Update the number of holes of R. + //n_added_holes<point2d, F>(R, qa); + } + + + //std::cout << "exiting step_fast 3" << std::endl; + } + + + template <typename V, typename P> + void update_set(const image2d<V>& u, + p_image2d& A, + p_image2d& N, + V& g) + { + // A <- {x belongs to N / u(x) == g} + A.clear(); + mln_assertion(A.is_empty()); + + // This is slow. + // A.insert(set::inter(N, u.domain()) | pw::value(u) == pw::cst(g)); + + // N <- N\{x belongs to N / u(x) == g} + // N.remove(set::inter(N, u.domain()) | pw::value(u) == pw::cst(g)); + // mln_assertion(check_step_4_2(u, N, g)); + + // This is faster. + +// if (N.npoints() >= 0) +// { + p_image2d_fwd_pixter<point2d> pxl(N); + for_all(pxl) + { + point2d tmp = pxl.to_point(); + if (u(tmp) == g) + { + A.insert(tmp); + N.remove(tmp); + } + } +// } +// else { +// mln_fwd_piter(p_image2d) p(N); +// for_all(p) +// { +// if (u(p) == g) +// { +// A.insert(p); +// N.remove(p); +// } +// } +// mln_assertion(check_step_4_2(u, N, g)); +// } + } + /// IF g < gn. + template <typename V, typename P, typename F> + void step_fast4_1 (image2d<V>& u, + p_image2d& A, + p_image2d& R, + p_image2d& N, + V& g, + V& gn, + fllt_node(P, V)*& current_region, + image2d<fllt_node(P, V)*>& regions + ) + { +// std::cout << "entering step_fast 4_1" << std::endl; + + // If the region is bounded + // Create a new conected component. + // FIXME : we can make it faster. + + if ((R.bbox() < u.domain()) || (R.npoints() == u.domain().npoints())) + { + mln_piter(p_image2d) p(R); + current_region = new fllt_node(P, V)(); + current_region->elt().brighter = F::parent_is_brighter; + current_region->elt().value = g; + for_all(p) + { + current_region->elt().points.insert(p); + + if (regions(p) == 0) + { + //current_region->elt().points.insert(p); + regions(p) = current_region; + } + else + { + if (regions(p)->parent() == 0) + regions(p)->set_parent(current_region); + } + } + + +// // Count the number of conected components of the border of R. +// static image2d<unsigned> tmp(u.domain().to_larger(1)); +// static image2d<bool> border_ima(tmp.domain()); +// level::fill(border_ima, false); + +// // level::fill(inplace(border_ima | N), true); +// // std::cout << "tmp border = " << tmp.border () << std::endl; +// // std::cout << "ima border = " << border_ima.border () << std::endl; +// mln_piter(p_image2d) z(N); +// for_all(z) +// { +// mln_assertion(border_ima.owns_(z)); +// border_ima(z) = true; +// } +// unsigned n; +// tmp = labeling::level(border_ima, true, F::bdr_nbh(), n); + +// // debug::println(border_ima); +// //std::cout << "nb composantes :" << n << std::endl; +// // debug::println(tmp); +// if (n > 1) +// { + +// // IF number of conected components of the border > 1 +// for (int i = 2; i <= n; i++) +// { +// // follow each border to find which is the exterior border +// // and which are the holes. Keep one pixel of each holes. + +// // WARNING : We trust labeling::level to label the exterior border with 1. +// current_region->elt().holes.insert(a_point_of(tmp | pw::value(tmp) == pw::cst(i))); + +// // FIXME : [optimisation] Remove from N border of holes???. +// // Recompute gn <- min u(x) x belongs to A +// } +// } + + } + + g = gn; + + update_set(u,A,N,g); + + // std::cout << "A :" << std::endl; + // if (A.npoints()) + // debug::println(u | A); + // std::cout << "N :" << std::endl; + // if (N.npoints()) + // debug::println(u | N); + +// std::cout << "exiting step_fast 4_1" << std::endl; + } + + + template <typename V, typename P> + bool check_step_4_2(const image2d<V>& u, + p_image2d& N, + V& g) + { + mln_fwd_piter(p_image2d) p(N); + for_all(p) + { + mln_assertion(u(p) != g); + } + return true; + } + + + /// IF g == gn. + template <typename V, typename P> + void step_fast4_2 (const image2d<V>& u, + p_image2d& A, + p_image2d& N, + V& g, + fllt_node(P, V)* current_region, + image2d<fllt_node(P, V)*>& regions + ) + { +// std::cout << "entering step_fast 4_2" << std::endl; + + update_set(u,A,N,g); + + // std::cout << "A :" << std::endl; + // if (A.npoints()) + // debug::println(u | A); + // std::cout << "N :" << std::endl; + // if (N.npoints()) + // debug::println(u | N); + +// std::cout << "exiting step_fast 4_2" << std::endl; + } + + /// IF g > gn. + template <typename V, typename P> + void step_fast4_3 (image2d<V>& u, + const image2d<bool>& tagged, + const p_image2d& R, + const V& g) + { +// std::cout << "entering step_fast 4_3" << std::endl; + + // set the gray-level of the pixels of R to g. + mln_piter(p_image2d) p(R); + for_all(p) + { + mln_assertion (tagged(p)); + u (p) = g; + } +// std::cout << "exiting step_fast 4_3" << std::endl; + + } + + + template <typename V, typename F> + fllt_tree(point2d, V)& + compute_level_set_fast(const image2d<V>& ima, + image2d< fllt_node(point2d, V)* >& regions) + { + typedef point2d P; + typedef image2d<V> I; + + // FIXME: not nice. + typedef mln::image_if< + mln::image2d<unsigned>, + mln::fun::greater_p2b_expr_<mln::pw::value_<mln::image2d<unsigned> >, + mln::pw::cst_<int> > + > I_IF; + + // Check + mln_assertion(ima.domain() == regions.domain()); + + // Declarations. + p_image2d R(ima.domain()), + N(ima.domain()), + A(ima.domain()); + + V g, gn; + point2d x0; + image2d<unsigned> min_locals(ima.domain()); + image2d<V> u = clone(ima); + border::fill(u, 0); + + //std::cout << "image U:" << std::endl; + // debug::println_with_border(u); + image2d<bool> tagged(ima.domain()); + fllt_node(P, V)* current_region; + + // INIT + R.clear(); + N.clear(); + A.clear(); + g= 0; + gn = 0; + current_region = 0; + + level::fill(regions, 0); + level::fill(tagged, false); + + // Get the locals extremums + unsigned nlabels; + min_locals = F::regional_extremum(ima, F::reg_nbh(), nlabels); + + // debug::println(min_locals); + // debug::println(min_locals | (pw::value(min_locals) > pw::cst(0))); + + /// Algorithm. + { + // For all locals extremums + I_IF min_locals_list(min_locals | (pw::value(min_locals) > pw::cst(0))); + mln_piter(I_IF) p(min_locals_list.domain()); + int cpt = 0, snapshop = 0; + for_all(p) + { + if (tagged(p)) + continue; + + //std::cout << "boucle no " << cpt++ << std::endl; + step_fast1(ima, p, g, x0); + step_fast2(A, R, N, x0); + mln_assertion(R.is_empty() && N.is_empty()); + while (1) + { + // To have a view of the progression of the algotithm. +#ifdef FLLT_TRACE + save_state(ima, R, A, N); +#endif + //std::cout << "g = " << g << std::endl; + step_fast3<V, P, F>(u, tagged, A, R, N, gn, current_region, regions); + /// step_fast4. + if (F::compare(g, gn)) + { + step_fast4_1<V, P, F>(u, A, R, N, g, gn, current_region, regions); + // GO TO 3) + continue; + } + + + if (g == gn) + { + step_fast4_2(u, A, N, g, current_region, regions); + // GO TO 3) + continue; + } + + + if (!F::compare(g, gn)) + { + step_fast4_3(u, tagged, R, g); + // GO TO 1) + break; + } + } + //std::cout << "current_region = " << current_region << std::endl; + } + } // end of Algorithm + + image2d<value::int_u8> output (ima.domain ()); + fllt_tree(P, V)& tree = *new fllt_tree(P, V)(current_region); + util::tree_to_image (tree, output); + + // util::display_tree(ima, tree); + + // debug::println(output); + // std::cout << std::endl; + // debug::println(ima); + + // if (output != ima) + // { + // std::cerr << "BUG!!!" << std::endl; + // abort(); + // } + + // io::pgm::save(output, "out.pgm"); + // std::cout << "out.pgm generate" + // << std::endl; + + + // debug::println(regions); + //debug::println(ima | regions(make:defined reference to `mln::fllt::lower<mln::value::int_u<8u> >::inc':point2d(-4,-1))->elt().points); + + return (tree); + + } // end of compute_level_set_fast + + } // end of namespace mln::fllt + +} // end of namespace mln + + + +#endif // ! MLN_FIXME_FLLT_COMPUTE_LEVEL_SET_FAST_HH Index: trunk/milena/sandbox/garrigues/fllt/fllt.hh =================================================================== --- trunk/milena/sandbox/garrigues/fllt/fllt.hh (revision 1869) +++ trunk/milena/sandbox/garrigues/fllt/fllt.hh (revision 1870) @@ -94,7 +94,7 @@ template <typename V> // Fixme : return type - void + fllt_tree(point2d, V) fllt(const image2d<V>& ima) { typedef point2d P; @@ -124,52 +124,7 @@ //fllt_tree(P, V) result_tree = merge_trees(upper_tree, lower_tree, upp_reg, low_reg, ima); fllt_tree(P, V) result_tree = merge_trees(lower_tree, upper_tree, low_reg, upp_reg, ima); - - std::cout << "4/ Generate outputs." << std::endl; - - image2d<value::int_u8> output (ima.domain ()); - util::tree_to_image (result_tree, output); - - - // io::pgm::save(output, "out_final.pgm"); - // std::cout << "out_final.pgm generate" - // << std::endl; - - - // util::display_tree(ima, lower_tree); - draw_tree(ima, result_tree); - - // debug::println(ima); - // debug::println(output); - - // if (output != ima) - // { - // std::cerr << "BUG!!!" << std::endl; - // abort(); - // } - - image2d<value::int_u8> viz(ima.domain()); - // image2d<value::int_u8> viz2(ima.domain()); - - // visualize_deepness(viz, lower_tree); - // level::stretch(viz, viz2); - // debug::println(viz); - // debug::println(viz2); - // io::pgm::save(viz2, "fllt.pgm"); - - visualize_bounds(viz, result_tree, 200); - //debug::println(viz); - io::pgm::save(viz, "fllt_bounds_200.pgm"); - - visualize_bounds(viz, result_tree, 100); - io::pgm::save(viz, "fllt_bounds_100.pgm"); - - visualize_bounds(viz, result_tree, 50); - io::pgm::save(viz, "fllt_bounds_50.pgm"); - - visualize_bounds(viz, result_tree, 20); - io::pgm::save(viz, "fllt_bounds_20.pgm"); - + return result_tree; } } // end of namespace mln::fllt Index: trunk/milena/sandbox/garrigues/fllt/local_configurations.hh =================================================================== --- trunk/milena/sandbox/garrigues/fllt/local_configurations.hh (revision 0) +++ trunk/milena/sandbox/garrigues/fllt/local_configurations.hh (revision 1870) @@ -0,0 +1,144 @@ +// Copyright (C) 2007, 2008 EPITA Research and Development Laboratory +// +// This file is part of the Olena Library. This library is free +// software; you can redistribute it and/or modify it under the terms +// of the GNU General Public License version 2 as published by the +// Free Software Foundation. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// General Public License for more details. +// +// You should have received a copy of the GNU General Public License +// along with this library; see the file COPYING. If not, write to +// the Free Software Foundation, 51 Franklin Street, Fifth Floor, +// Boston, MA 02111-1307, USA. +// +// As a special exception, you may use this file as part of a free +// software library without restriction. Specifically, if other files +// instantiate templates or use macros or inline functions from this +// file, or you compile this file and link it with other files to +// produce an executable, this file does not by itself cause the +// resulting executable to be covered by the GNU General Public +// License. This exception does not however invalidate any other +// reasons why the executable file might be covered by the GNU General +// Public License. + + +#ifndef MLN_FIXME_FLLT_LOCAL_CONFIGURATIONS_HH +# define MLN_FIXME_FLLT_LOCAL_CONFIGURATIONS_HH + +/*! \file local_configurations.hh + * + * \brief Stuffs related to local configurations for the FLLT. + * + */ + +#include <mln/core/clock_neighb2d.hh> +#include "upper.hh" +#include "lower.hh" + +namespace mln +{ + namespace fllt + { + + template <typename V> + struct lower; + + template <typename V> + struct upper; + + template <typename F> + struct added_holes; + + template <typename V> + struct added_holes<upper<V> > + { + // For each possible configuration (of the c8-neighbourghood) of the + // added pixel, stock the number of added holes. + // ==> region with c8, holes with c4. + static char ret[256]; + }; + + template <typename V> + char added_holes<upper<V> >::ret[256] = + { + 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, + 0, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, + 0, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 1, 1, 1, 1, + 0, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, + + 0, 0, 1, 0, 1, 0, 1, 0, 1, 1, 2, 1, 1, 0, 1, 0, + 0, 0, 1, 0, 0, -1, 0, -1, 0, 0, 1, 0, 0, -1, 0, -1, + 0, 0, 1, 0, 1, 0, 1, 0, 1, 1, 2, 1, 1, 0, 1, 0, + 0, 0, 1, 0, 0, -1, 0, -1, 0, 0, 1, 0, 0, -1, 0, -1, + + 0, 0, 1, 0, 1, 0, 1, 0, 1, 1, 2, 1, 1, 0, 1, 0, + 1, 1, 2, 1, 1, 0, 1, 0, 1, 1, 2, 1, 1, 0, 1, 0, + 1, 1, 2, 1, 2, 1, 2, 1, 2, 2, 3, 2, 2, 1, 2, 1, + 1, 1, 2, 1, 1, 0, 1, 0, 1, 1, 2, 1, 1, 0, 1, 0, + + 0, 0, 1, 0, 1, 0, 1, 0, 1, 1, 2, 1, 1, 0, 1, 0, + 0, 0, 1, 0, 0, -1, 0, -1, 0, 0, 1, 0, 0, -1, 0, -1, + 0, 0, 1, 0, 1, 0, 1, 0, 1, 1, 2, 1, 1, 0, 1, 0, + 0, 0, 1, 0, 0, -1, 0, -1, 0, 0, 1, 0, 0, -1, 0, -1 + }; + + template <typename V> + struct added_holes<lower<V> > + { + // For each possible configuration (of the c4-neighbourghood) of the + // added pixel, stock the number of added holes. + // ==> region with c4, holes with c8. + static char ret[256]; + }; + + template <typename V> + char added_holes<lower<V> >::ret[256] = + { + 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, + 0, 1, 0, 1, 1, 2, 1, 1, 0, 1, 0, 1, 0, 1, 0, 0, + 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, + 0, 1, 0, 1, 1, 2, 1, 1, 0, 1, 0, 1, 0, 1, 0, 0, + + 0, 1, 0, 1, 1, 2, 1, 1, 0, 1, 0, 1, 1, 2, 1, 1, + 1, 2, 1, 2, 2, 3, 2, 2, 1, 2, 1, 2, 1, 2, 1, 1, + 0, 1, 0, 1, 1, 2, 1, 1, 0, 1, 0, 1, 1, 2, 1, 1, + 0, 1, 0, 1, 1, 2, 1, 1, 0, 1, 0, 1, 0, 1, 0, 0, + + 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, + 0, 1, 0, 1, 1, 2, 1, 1, 0, 1, 0, 1, 0, 1, 0, 0, + 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, + 0, 1, 0, 1, 1, 2, 1, 1, 0, 1, 0, 1, 0, 1, 0, 0, + + 0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 0, + 1, 1, 1, 1, 2, 2, 2, 1, 1, 1, 1, 1, 1, 1, 1, 0, + 0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 0, + 0, 0, 0, 0, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, -1 + }; + + template <typename P, typename F> + char n_added_holes(const p_image2d& set, const point2d& p) + { + unsigned char res = 0; + dpoint2d dp(1,1); + clock_neighb2d nbh = F::reg_c_nbh(dp); + mln_fwd_niter_(clock_neighb2d) n(nbh , p); + for_all(n) + { + res = (res << 1); + if (set.has(n)) + res = res | 1; + } +// std::cout << int(res) << ":" << int(added_holes<F>::ret[res]) << " hole added." << std::endl; + return added_holes<F>::ret[res]; + } + } // end of namespace mln::fllt + +} // end of namespace mln + + + +#endif // ! MLN_FIXME_FLLT_HH Index: trunk/milena/sandbox/garrigues/fllt/test.cc =================================================================== --- trunk/milena/sandbox/garrigues/fllt/test.cc (revision 0) +++ trunk/milena/sandbox/garrigues/fllt/test.cc (revision 1870) @@ -0,0 +1,62 @@ +# include "fllt.hh" +//# include "compute_level_set.hh" + +# include "compute_level_set_fast2.hh" + +# include <mln/core/image2d.hh> +# include <mln/core/cast_image.hh> +# include <mln/core/clone.hh> +# include <mln/value/int_u8.hh> +# include <mln/debug/println.hh> +# include <mln/convert/to_w_window.hh> +# include <mln/core/w_window2d_int.hh> +# include <mln/convert/to_image.hh> +# include <mln/level/fill.hh> +# include <mln/io/pgm/load.hh> +# include <mln/io/pbm/load.hh> + +int main() +{ + + using namespace mln; + + typedef point2d P; + typedef int V; + + //trace::quiet = false; + +// V ws[81] = {5,5,5,5,5,5,5,5,5, +// 5,5,5,5,5,5,5,5,5, +// 5,5,5,5,5,5,8,8,5, +// 5,5,5,9,5,5,8,8,5, +// 5,5,5,9,5,5,8,8,5, +// 5,5,5,5,5,5,8,8,5, +// 5,5,5,5,5,5,8,8,5, +// 5,5,5,5,5,5,5,5,5, +// 5,5,5,5,5,5,5,5,5}; + +// w_window2d_int w_win = make::w_window2d(ws); +// image2d<V> ima = convert::to_image(w_win); + + image2d<value::int_u8> ima_ = io::pgm::load<value::int_u8>("small.pgm"); + image2d<V> ima(ima_.domain()); + + level::fill(ima, ima_); + +// { +// image2d<fllt_node(P, V)*> low_reg(ima.domain()); +// fllt_tree(P, V) lower_tree; +// lower_tree +// = fllt::compute_level_set<V, fllt::lower<V> >(ima, low_reg); +// fllt::draw_tree(ima, lower_tree); +// } + + + { + image2d<fllt_node(P, V)*> low_reg(ima.domain()); + fllt_tree(P, V) lower_tree; + lower_tree + = fllt_fast::compute_level_set_fast<V, fllt::upper<V> >(ima, low_reg); + //fllt::draw_tree(ima, lower_tree); + } +} Index: trunk/milena/sandbox/garrigues/fllt/give_confs.cc =================================================================== --- trunk/milena/sandbox/garrigues/fllt/give_confs.cc (revision 0) +++ trunk/milena/sandbox/garrigues/fllt/give_confs.cc (revision 1870) @@ -0,0 +1,56 @@ +#include <mln/core/image2d.hh> +#include <mln/core/point2d.hh> +#include <mln/core/p_set.hh> +#include <mln/core/clock_neighb2d.hh> +#include <mln/core/neighb2d.hh> + +#include <mln/labeling/level.hh> +#include <mln/level/fill.hh> +#include <mln/debug/println.hh> +#include <iomanip> + +int main() +{ + using namespace mln; + using namespace mln::value; + + image2d<bool> ima(3,3); + point2d p(1,1); + dpoint2d dp(-1,0); + clock_neighb2d nbh = cc8(dp); + int n_holes = 0; + + bool tab[8][8]; + for (int i = 0; i < 256; i++) + { + level::fill(ima, false); + int_u8 tmp = i; + + mln_fwd_niter_(clock_neighb2d) n(nbh , p); + for_all(n) + { + if (tmp % 2) + ima(n) = true; + tmp = tmp >> 1; + } + + + int x, y; + labeling::level(ima, false, c8(), x); + + + ima(p) = true; + labeling::level(ima, false, c8(), y); + +// std::cout << "----- conf no " << i << "------" << std::endl; +// debug::println(ima); +// std::cout << "----- " << x << " holes before------" << std::endl; +// std::cout << "----- " << y << " holes after------" << std::endl; +// std::cout << "----- " << y - x << " holes added------" << std::endl << std::endl; + + std::cout << std::setw(2) << y - x << ", "; + if (! ((i + 1) % 4)) std::cout << " "; + if (! ((i + 1) % 16)) std::cout << std::endl; + if (! ((i + 1) % 64)) std::cout << std::endl; + } +} Index: trunk/milena/sandbox/garrigues/fllt/test_fllt2.cc =================================================================== --- trunk/milena/sandbox/garrigues/fllt/test_fllt2.cc (revision 1869) +++ trunk/milena/sandbox/garrigues/fllt/test_fllt2.cc (revision 1870) @@ -1,4 +1,4 @@ -# include "fllt_optimized.hh" +# include "fllt.hh" # include <mln/core/image2d.hh> # include <mln/core/clone.hh> # include <mln/value/int_u8.hh> @@ -29,5 +29,5 @@ image2d<int> ima = convert::to_image(w_win); fllt_tree(point2d, int) t = fllt::fllt(ima); - fllt::debug(ima, t); + fllt::draw_tree(ima, t); }

1 0

milena r1869: Problems on image_if_value<sub_image>
by Matthieu Garrigues 16 Apr '08

16 Apr '08

URL: https://svn.lrde.epita.fr/svn/oln/trunk/milena ChangeLog: 2008-04-16 Matthieu Garrigues <garrigues(a)lrde.epita.fr> * sandbox/garrigues/fllt/essai.cc: New, Problems on image_if_value<sub_image>. --- essai.cc | 104 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 104 insertions(+) Index: trunk/milena/sandbox/garrigues/fllt/essai.cc =================================================================== --- trunk/milena/sandbox/garrigues/fllt/essai.cc (revision 0) +++ trunk/milena/sandbox/garrigues/fllt/essai.cc (revision 1869) @@ -0,0 +1,104 @@ +#include <mln/core/image2d.hh> + +#include <mln/core/sub_image.hh> + +#include <mln/core/image_if_value.hh> + +#include <mln/accu/bbox.hh> +#include <mln/level/fill.hh> +#include <mln/debug/println.hh> + +#include "types.hh" + +using namespace mln; +using value::int_u8; + +// class ran_domains +// { +// public: + +// ran_domains(const box2d& b) +// : ima_(b) +// { +// bb_.init(); +// level::fill(ima_, false); +// bb_.take(make::point2d(2,2)); +// bb_.take(make::point2d(2,3)); +// } + +// sub_image<image2d<value::int_u8>, box2d> image() +// { return ima_ | bb_.to_result(); } + +// image2d<value::int_u8>& full_image() { return ima_; } + +// const box2d& bbox() const { return bb_.to_result(); } + +// private: + +// image2d<int_u8> ima_; +// accu::bbox<point2d> bb_; +// }; + +# define image_if_value_sub image_if_value<const sub_image<image2d<value::int_u8>, box2d> > +# define sub_if_value_image sub_image<const image_if_value<image2d<int_u8 > >, box2d> +# define SET_R 3 +# define SET_A 2 +# define SET_N 1 + +int main() +{ + + { + int_u8 vs[5][5] = { {1, 2, 3, 4, 9}, + {5, 6, 7, 8, 9}, + {9, 10, 6,12, 9}, + {13,14,15, 6, 9}, + {13,14,15,16, 9} }; + + image2d<int_u8> ima = make::image2d(vs); + + box2d b(make::box2d(1,1,3,3)); + + // Browse a image_if_value. + { + image_if_value<image2d<value::int_u8> > ima2 = ima | 6; + image_if_value<image2d<value::int_u8> >::fwd_piter qa(ima2.domain()); + for_all(qa) + {} + } + + // Browse a sub_image. + { + sub_image<image2d<value::int_u8>, box2d> ima2 = ima | b; + sub_image<image2d<value::int_u8>, box2d>::fwd_piter qa(ima2.domain()); + for_all(qa) + {} + } + + // Sub_image, image_if_value conversion. + { + image_if_value<const sub_image<image2d<value::int_u8>, box2d> > s = ima | b | 6; + sub_image<image2d<value::int_u8>, box2d> ima2 = s; + } + + // Problem 1 : Browse a sub_image of a image_if_value. + { + sub_if_value_image ima2 = ima | 6 | b; + debug::println(ima2); + } + + // Problem 2 : Browse a image_if_value of a sub_image. + { + image_if_value_sub s = ima | b | 6; + + image_if_value_sub::fwd_piter qa(s.domain()); + for_all(qa) + {} + } + + } +} + + +// mln::image_if_value<const mln::sub_image<mln::image2d<mln::value::int_u<8u> >, mln::box_<mln::point_<mln::grid::square, int> > > > +// mln::image_if_value<mln::sub_image<mln::image2d<mln::value::int_u<8u> >, mln::box_<mln::point_<mln::grid::square, int> > > >

1 0

1868: Sandbox: ICP: Fix Qk driving.
by Ugo Jardonnet 16 Apr '08

16 Apr '08

https://svn.lrde.epita.fr/svn/oln/trunk/milena Index: ChangeLog from Ugo Jardonnet <ugo.jardonnet(a)lrde.epita.fr> Sandbox: ICP: Fix Qk driving. * sandbox/jardonnet/test/bench: Add Qk driving log. * sandbox/jardonnet/registration/quat7.hh (to_vec): Add function. * sandbox/jardonnet/registration/icp.hh: Add qk/dk buffers. * sandbox/jardonnet/registration/update_qk.hh: Fix. registration/icp.hh | 16 +++++++----- registration/interpolation.hh | 2 - registration/quat7.hh | 54 +++++++++++++++++------------------------- registration/update_qk.hh | 54 ++++++++++++++++++++++-------------------- test/bench | 54 ++++++++++++++++++++++++++++++++++++++++++ test/icp.cc | 3 +- 6 files changed, 118 insertions(+), 65 deletions(-) Index: sandbox/jardonnet/test/bench --- sandbox/jardonnet/test/bench (revision 1867) +++ sandbox/jardonnet/test/bench (working copy) @@ -32,3 +32,57 @@ gprof icp_map 1 call chamfer 0.33s icp_lazy 16 call memo 0.21 + +**icp 01.pbm 02.pbm 1 1** +Register 943 points +k error dqk +0 18.9244 41.2757 +1 16.6149 21.4739 +2 15.2707 12.8804 +3 14.1741 10.3871 +4 13.4248 7.97125 +5 12.8224 7.72872 +6 12.5292 6.52063 +7 12.4356 4.50726 +8 12.3936 3.98019 +9 12.3728 3.42463 +10 12.3798 2.88974 +11 12.3838 2.42067 +12 12.3916 1.56065 +13 12.3627 1.42373 +14 12.3857 1.32764 +15 12.3562 1.2033 +16 12.3713 1.01804 +17 12.3657 0.672616 +18 12.3829 0.599551 +19 12.3859 0.423572 +20 12.3897 0.284606 +21 12.3894 0.183171 +22 12.3888 0.0807438 +closest_point(Ck[i]) = 12212 +Pts processed = 21689 +./icpD 01.pbm 02.pbm 1 1 18,19s user 0,05s system 97% cpu 18,619 total + + +**icp_drived 01.pbm 02.pbm 1 1** +Register 943 points +k error dqk +0 18.9244 41.2757 +1 16.6149 21.4739 +2 15.2707 12.8804 +3 14.1741 10.3871 +4 12.8711 15.9425 +5 12.3806 13.1241 +6 12.3601 3.38092 +7 12.3878 3.1961 +8 12.3975 5.48096 +9 12.3887 3.08797 +10 12.3611 2.27939 +11 12.3779 0.969556 +12 12.3955 0.581234 +13 12.3795 0.435324 +14 12.3931 0.0650237 +closest_point(Ck[i]) = 9847 +Pts processed = 14145 +./icpD 01.pbm 02.pbm 1 1 13,50s user 0,03s system 98% cpu 13,750 total + Index: sandbox/jardonnet/test/icp.cc --- sandbox/jardonnet/test/icp.cc (revision 1867) +++ sandbox/jardonnet/test/icp.cc (working copy) @@ -63,7 +63,7 @@ qk.apply_on(c, c, c.npoints()); - const box_<point2d> box2d(600,600); + const box_<point2d> box2d(400,700); image2d<bool> output(box2d, 1); //to 2d : projection (FIXME:if 3d) @@ -75,5 +75,6 @@ } io::pbm::save(output, "registred.pbm"); + } Index: sandbox/jardonnet/registration/interpolation.hh --- sandbox/jardonnet/registration/interpolation.hh (revision 1867) +++ sandbox/jardonnet/registration/interpolation.hh (working copy) @@ -10,7 +10,7 @@ { int k,j,i; float phi,ff,b; - std::vector<float> s(n,0); + std::vector<float> s(n); for (i = 0; i <= n; i++) s[i] = cof[i] = 0.0; Index: sandbox/jardonnet/registration/quat7.hh --- sandbox/jardonnet/registration/quat7.hh (revision 1867) +++ sandbox/jardonnet/registration/quat7.hh (working copy) @@ -79,6 +79,28 @@ return os; } + float operator[](int i) + { + if (i < n) + return _qT[i]; + else + return _qR.to_vec()[i - n]; + } + + algebra::vec<7,float> to_vec() + { + algebra::vec<7,float> v; + + for (int i = 0; i < 7; i++) + { + if (i < n) + v[i] = _qT[i]; + else + v[i] = _qR.to_vec()[i - n]; + } + return v; + } + quat7 operator*(float f) { return quat7(_qR * f, _qT * f); @@ -100,21 +122,6 @@ } }; - - // very usefull routine - - template <unsigned p, unsigned q, unsigned n, unsigned m> - void put(const algebra::mat<p,q,float>& in, // a matrix to put into... - unsigned row, unsigned col, // top-left location - algebra::mat<n,m,float>& inout) // ...a matrix to modify - { - assert(row + p <= n && col + q <= m); - for (unsigned i = 0; i < p; ++i) - for (unsigned j = 0; j < q; ++j) - inout(row + i, col + j) = in(i,j); - } - - template <typename P, typename M> quat7<P::dim> match(const p_array& C, const algebra::vec<P::dim,float>& mu_C, @@ -124,7 +131,6 @@ { //FIXME: mix mu_Xk and qk loop - //mu_Xk = center map(Ck) algebra::vec<P::dim,float> mu_Xk(literal::zero); for (size_t i = 0; i < c_length; ++i) @@ -144,22 +150,6 @@ } Mk /= c_length; - - /* - const algebra::mat<P::dim,P::dim,float> Ak = Mk - trans(Mk); - - const float v[3] = {Ak(1,2), Ak(2,0), Ak(0,1)}; - const algebra::mat<3,1,float> D = make::mat<3,1,3,float>(v); - - - algebra::mat<4,4,float> Qk(literal::zero); - Qk(0,0) = tr(Mk); - put(trans(D), 0,1, Qk); - put(D, 1,0, Qk); - - put(Mk + trans(Mk) - algebra::mat<P::dim,P::dim,float>::identity() * tr(Mk), 1,1, Qk); - */ - algebra::vec<3,float> a; a[0] = Mk(1,2) - Mk(2,1); a[1] = Mk(2,0) - Mk(0,2); Index: sandbox/jardonnet/registration/icp.hh --- sandbox/jardonnet/registration/icp.hh (revision 1867) +++ sandbox/jardonnet/registration/icp.hh (working copy) @@ -92,10 +92,10 @@ std::cout << "k\terror\tdqk" << std::endl; #endif - quat7<P::dim> buf_qk[3]; + quat7<P::dim> buf_qk[4]; float buf_dk[3]; - float err; + float err = 100; //float err_bis; p_array Ck(C); //FIXME: Xk copy C @@ -106,7 +106,7 @@ unsigned int k = 0; do { - //update buff dk + //update buff dk //FIXME: rewrite it buf_dk[2] = buf_dk[1]; buf_dk[1] = buf_dk[0]; buf_dk[0] = err; @@ -114,13 +114,17 @@ //compute qk qk = match(C, mu_C, Ck, map, c_length); - //update buff qk + //update buff qk //FIXME: rewrite it + buf_qk[3] = buf_qk[2]; buf_qk[2] = buf_qk[1]; buf_qk[1] = buf_qk[0]; buf_qk[0] = qk; //update qk + if (k > 3) qk = update_qk(buf_qk, buf_dk); + qk._qR.set_unit(); + buf_qk[0] = qk; //Ck+1 = qk(C) qk.apply_on(C, Ck, c_length); @@ -151,13 +155,13 @@ //plot file std::cout << k << '\t' << err << '\t' - << (qk - old_qk).sqr_norm() << '\t' + << (qk - buf_qk[1]).sqr_norm() << '\t' << std::endl; //count the number of points processed pts += c_length; #endif k++; - } while (k < 3 || (qk - buf_qk[0]).sqr_norm() > epsilon); + } while (k < 3 || (qk - buf_qk[1]).sqr_norm() > epsilon); trace::exiting("registration::impl::icp_"); return Ck; Index: sandbox/jardonnet/registration/update_qk.hh --- sandbox/jardonnet/registration/update_qk.hh (revision 1867) +++ sandbox/jardonnet/registration/update_qk.hh (working copy) @@ -9,23 +9,23 @@ float arc(quat7<3> dqk, quat7<3> dqk_1) { - //algebra::mat - - return 0; + float res = 0.0; + for (int i = 0; i < 7; i++) + res += dqk[i] * dqk_1[i]; + res = res / (dqk.sqr_norm() * dqk_1.sqr_norm()); + return acos(res) * 180.0 / 3.14159265; } - quat7<3> update_qk(quat7<3> qk[3], + quat7<3> update_qk(quat7<3> qk[4], float dk[3]) { - static float tetak = 11; - static quat7<3> dqk; - - quat7<3> dqk_1 = dqk; - dqk = qk[0] - qk[1]; + quat7<3> dqk_2 = qk[2] - qk[3]; + quat7<3> dqk_1 = qk[1] - qk[2]; + quat7<3> dqk = qk[0] - qk[1]; - float tetak_1 = tetak; - tetak = arc(dqk, dqk_1); + float tetak_1 = arc(dqk_1, dqk_2); + float tetak = arc(dqk, dqk_1); if (tetak < 10 and tetak_1 < 10) //10degrees { @@ -34,33 +34,37 @@ float coef[3]; vk[0] = 0; - vk[1] = - dqk.norm(); - vk[2] = - dqk_1.norm() + vk[1]; + vk[1] = - dqk.sqr_norm(); + vk[2] = - dqk_1.sqr_norm() + vk[1]; - float a1 = (dk[0] + dk[2]) / (vk[2] - vk[0]); - float b1 = dk[1] - ((dk[0] + dk[2]) / (vk[2] - vk[0])) * vk[1]; - float v1 = - b1 / a1; // > 0 + float a1 = (dk[2] + dk[0]) / (vk[2] - vk[0]); + float b1 = dk[1] - ((dk[2] + dk[0]) / (vk[2] - vk[0])) * vk[1]; + float v1 = math::abs(- b1 / a1); // > 0 - polcoe(vk, dk, 3, coef); + polcoe(vk, dk, 2, coef); float a2 = coef[0]; float b2 = coef[1]; float v2 = -b2 / (2. * a2); - float vmax = 25 * dqk.norm(); - - if (0 < v2 and v2 < v1 and v1 < vmax) - return qk[0] + (dqk / dqk.norm()) * v2; + float vmax = 25 * dqk.sqr_norm(); - if (0 < v1 and v1 < v2 and v2 < vmax) - return qk[0] + (dqk / dqk.norm()) * v1; + if ((0 < v2 and v2 < v1 and v1 < vmax) or + (0 < v2 and v2 < vmax and vmax < v1)) + return qk[0] + (dqk / dqk.sqr_norm()) * v2; + + if ((0 < v1 and v1 < v2 and v2 < vmax) or + (0 < v1 and v1 < vmax and vmax < v2) or + (v2 < 0 and 0 < v1 and v1 < vmax)) + return qk[0] + (dqk / dqk.sqr_norm()) * v1; if (v1 > vmax and v2 > vmax) - return qk[0] + (dqk / dqk.norm()) * vmax; + return qk[0] + (dqk / dqk.sqr_norm()) * vmax; + std::cout << "echec" << std::endl; + } return qk[0]; } - } }

1 0

1867: Sandbox: ICP: Qk driving.
by Ugo Jardonnet 15 Apr '08

15 Apr '08

https://svn.lrde.epita.fr/svn/oln/trunk/milena Index: ChangeLog from Ugo Jardonnet <ugo.jardonnet(a)lrde.epita.fr> Sandbox: ICP: Qk driving. * sandbox/jardonnet/test/Makefile: Add rules. * sandbox/jardonnet/TODO: Update. * sandbox/jardonnet/registration/interpolation.hh: give coef of the interpolated polynome. * sandbox/jardonnet/registration/quat7.hh: Update interface */-+. * sandbox/jardonnet/registration/cloud.hh: Update center. * sandbox/jardonnet/registration/icp.hh: Add Optimization on qk. * sandbox/jardonnet/registration/update_qk.hh: Optimization on qk. TODO | 8 ++++ registration/cloud.hh | 12 +------ registration/icp.hh | 28 ++++++++++++----- registration/interpolation.hh | 41 +++++++++++++++++++++++++ registration/quat7.hh | 27 ++++++++++++++-- registration/update_qk.hh | 68 ++++++++++++++++++++++++++++++++++++++++++ test/Makefile | 19 ++++++++--- test/plotscript | 2 - 8 files changed, 180 insertions(+), 25 deletions(-) Index: sandbox/jardonnet/test/plotscript --- sandbox/jardonnet/test/plotscript (revision 1866) +++ sandbox/jardonnet/test/plotscript (working copy) @@ -2,4 +2,4 @@ set ylabel "e" set zlabel "time" splot "log.dat" -pause 10 +pause -1 Index: sandbox/jardonnet/test/Makefile --- sandbox/jardonnet/test/Makefile (revision 1866) +++ sandbox/jardonnet/test/Makefile (working copy) @@ -13,17 +13,27 @@ gau: g++ gaussian.cc $(FLAG) -o 'gau' -icpD: icp.cc +icpD: icp.cc +depend g++ icp.cc -I../../.. -g -o 'icpD' -icp: icp.cc - g++ icp.cc -I../../.. -O3 -DNDEBUG -o 'icp' + +icp: icp.cc +depend icp.o + g++ icp.o -I../../.. -O3 -DNDEBUG -o 'icp' + +icp.o: + g++ -c icp.cc -I../../.. -O3 -DNDEBUG bench: - ./check.rb + ruby bench.rb clean: rm -f -- ./bin/* rm -f sub gsub gau icpD icp rm -f log.dat registred.pbm + rm -f i_* + +.PHONY: check bench icpD + ++depend: + g++ -M icp.cc $(FLAG) > +depend -.PHONY: check \ No newline at end of file +include +depend \ No newline at end of file Index: sandbox/jardonnet/TODO --- sandbox/jardonnet/TODO (revision 1866) +++ sandbox/jardonnet/TODO (working copy) @@ -11,3 +11,11 @@ adapt test for threshold (old thresholding) - - + +La translation est mauvaise car certain point de Xk sont compt� plusieurs fois. /bof + +- - +write +trans : vec --> quat + v |-> m = trans(t) + Index: sandbox/jardonnet/registration/interpolation.hh --- sandbox/jardonnet/registration/interpolation.hh (revision 0) +++ sandbox/jardonnet/registration/interpolation.hh (revision 0) @@ -0,0 +1,41 @@ +#ifndef MLN_INTERPOLATION_HH +# define MLN_INTERPOLATION_HH + + +namespace mln +{ + + void polcoe(float x[], float y[], int n, + float cof[]) + { + int k,j,i; + float phi,ff,b; + std::vector<float> s(n,0); + + for (i = 0; i <= n; i++) + s[i] = cof[i] = 0.0; + s[n] = -x[0]; + for (i = 1; i <= n; i++) + { + for (j = n - i; j <= n - 1; j++) + s[j] -= x[i] * s[j+1]; + s[n] -= x[i]; + } + for (j = 0; j <= n; j++) + { + phi = n + 1; + for (k = n; k >= 1; k--) + phi = k * s[k] + x[j] * phi; + ff = y[j] / phi; + b = 1.0; + for (k = n; k >= 0; k--) + { + cof[k] += b * ff; + b = s[k] + x[j] * b; + } + } + } + +} + +#endif // ! MLN_INTERPOLATION_HH Index: sandbox/jardonnet/registration/quat7.hh --- sandbox/jardonnet/registration/quat7.hh (revision 1866) +++ sandbox/jardonnet/registration/quat7.hh (working copy) @@ -13,7 +13,7 @@ # include "power_it.hh" # include "frankel_young.hh" -# include <mln/norm/l2.hh> +# include <mln/norm/all.hh> # include <mln/trait/all.hh> @@ -46,11 +46,12 @@ return norm::l2(_qR.to_vec()) + norm::l2(_qT); } - quat7 operator-(const quat7& rhs) const + float norm() const { - return quat7(_qR - rhs._qR, _qT - rhs._qT); + return norm::l1(_qR.to_vec()) + norm::l1(_qT); } + // quat7 is an object-function algebra::vec<n,float> operator()(const algebra::vec<n,float>& v) const @@ -77,6 +78,26 @@ std::cout << q._qT << std::endl; return os; } + + quat7 operator*(float f) + { + return quat7(_qR * f, _qT * f); + } + + quat7 operator/(float f) + { + return quat7(_qR / f, _qT / f); + } + + quat7 operator-(const quat7& rhs) const + { + return quat7(_qR - rhs._qR, _qT - rhs._qT); + } + + quat7 operator+(const quat7& rhs) const + { + return quat7(_qR + rhs._qR, _qT + rhs._qT); + } }; Index: sandbox/jardonnet/registration/cloud.hh --- sandbox/jardonnet/registration/cloud.hh (revision 1866) +++ sandbox/jardonnet/registration/cloud.hh (working copy) @@ -41,24 +41,18 @@ return tmp; } - template <typename P> + template <typename P, typename M> float rms(const p_array& a1, - const p_array& a2, + M& map, const size_t length) { assert(length <= a1.npoints()); - assert(length <= a2.npoints()); - /* - float f = 0.f; - for (size_t i = 0; i < a1.npoints(); ++i) - f += sqr_norm(a1[i] - a2[i]); - return f / a1.npoints();*/ float f = 0.f; for (size_t i = 0; i < length; ++i) { algebra::vec<3,float> a1f = a1[i]; - algebra::vec<3,float> a2f = a2[i]; + algebra::vec<3,float> a2f = map(a1[i]); f += norm::l2(a1f - a2f); } return f / a1.npoints(); Index: sandbox/jardonnet/registration/icp.hh --- sandbox/jardonnet/registration/icp.hh (revision 1866) +++ sandbox/jardonnet/registration/icp.hh (working copy) @@ -46,7 +46,7 @@ # include "cloud.hh" # include "quat7.hh" -# include "projection.hh" +# include "update_qk.hh" # include "chamfer.hh" namespace mln @@ -91,10 +91,13 @@ << c_length << " points" << std::endl; std::cout << "k\terror\tdqk" << std::endl; #endif - quat7<P::dim> old_qk; + + quat7<P::dim> buf_qk[3]; + float buf_dk[3]; + float err; //float err_bis; - p_array Ck(C), Xk(C); //FIXME: Xk copy C + p_array Ck(C); //FIXME: Xk copy C algebra::vec<P::dim,float> mu_C = center(C, c_length), mu_Xk; @@ -103,15 +106,27 @@ unsigned int k = 0; do { + //update buff dk + buf_dk[2] = buf_dk[1]; + buf_dk[1] = buf_dk[0]; + buf_dk[0] = err; + //compute qk - old_qk = qk; qk = match(C, mu_C, Ck, map, c_length); + //update buff qk + buf_qk[2] = buf_qk[1]; + buf_qk[1] = buf_qk[0]; + buf_qk[0] = qk; + + //update qk + qk = update_qk(buf_qk, buf_dk); + //Ck+1 = qk(C) qk.apply_on(C, Ck, c_length); //err = d(Ck+1,Xk) - err = rms(Ck, Xk, c_length); + err = rms(Ck, map, c_length); #ifndef NDEBUG { @@ -141,9 +156,8 @@ //count the number of points processed pts += c_length; #endif - k++; - } while (k < 3 || (qk - old_qk).sqr_norm() > epsilon); + } while (k < 3 || (qk - buf_qk[0]).sqr_norm() > epsilon); trace::exiting("registration::impl::icp_"); return Ck; Index: sandbox/jardonnet/registration/update_qk.hh --- sandbox/jardonnet/registration/update_qk.hh (revision 0) +++ sandbox/jardonnet/registration/update_qk.hh (revision 0) @@ -0,0 +1,68 @@ +#ifndef MLN_UPDATE_QK_HH +# define MLN_UPDATE_QK_HH + +#include "interpolation.hh" +#include "quat7.hh" + +namespace mln +{ + + float arc(quat7<3> dqk, quat7<3> dqk_1) + { + //algebra::mat + + return 0; + } + + + quat7<3> update_qk(quat7<3> qk[3], + float dk[3]) + { + static float tetak = 11; + static quat7<3> dqk; + + quat7<3> dqk_1 = dqk; + dqk = qk[0] - qk[1]; + + float tetak_1 = tetak; + tetak = arc(dqk, dqk_1); + + if (tetak < 10 and tetak_1 < 10) //10degrees + { + float vk[3]; + float dk[3]; + float coef[3]; + + vk[0] = 0; + vk[1] = - dqk.norm(); + vk[2] = - dqk_1.norm() + vk[1]; + + float a1 = (dk[0] + dk[2]) / (vk[2] - vk[0]); + float b1 = dk[1] - ((dk[0] + dk[2]) / (vk[2] - vk[0])) * vk[1]; + float v1 = - b1 / a1; // > 0 + + polcoe(vk, dk, 3, coef); + float a2 = coef[0]; + float b2 = coef[1]; + float v2 = -b2 / (2. * a2); + + float vmax = 25 * dqk.norm(); + + if (0 < v2 and v2 < v1 and v1 < vmax) + return qk[0] + (dqk / dqk.norm()) * v2; + + if (0 < v1 and v1 < v2 and v2 < vmax) + return qk[0] + (dqk / dqk.norm()) * v1; + + if (v1 > vmax and v2 > vmax) + return qk[0] + (dqk / dqk.norm()) * vmax; + + + return qk[0]; + } + } + +} + + +#endif // ! MLN_UPDATE_QK_HH

1 0

1865: Update documentation about image types and properties.
by Nicolas Ballas 15 Apr '08

15 Apr '08

https://svn.lrde.epita.fr/svn/oln/trunk/milena Index: ChangeLog from Nicolas Ballas <ballas(a)lrde.epita.fr> Update documentation about image types and properties. * sandbox/ballas/doc/image_tours.txt: add the morpher types. * sandbox/ballas/doc/draft.txt: update. draft.txt | 10 ++++++ image_tours.txt | 85 ++++++++++++++++++++++++++++++++++++++++++++++++++------ 2 files changed, 86 insertions(+), 9 deletions(-) Index: sandbox/ballas/doc/image_tours.txt --- sandbox/ballas/doc/image_tours.txt (revision 1864) +++ sandbox/ballas/doc/image_tours.txt (working copy) @@ -432,20 +432,87 @@ No change in the property. +FIXME: Is there any property define for these types? + +**** decorated_image + templated image: I (image), D (decorator). +Image that can have additional features. + + +**** plain + Prevent an image from sharing data. + +**** safe + Make an image accessible at undefined location. + + +**** interpolated + Make an image readable with float coordinates. + +**** neigh::image + Add a neighbohood to an image. + +Changed properties: +category -> morpher +neighb -> some + + *** Domain morpher -*** Value morpher +**** sub_image + +Make image restricted to a given pset. + +Changed Properties: + category -> domain_morpher + border -> none -*** Image_if +**** hexa -??? + Make an hexagonal Mesh of the image. +(work only with 2D image). + +Changed Properties: + category -> domain_morpher + grid -> regular + +**** Image_if + Parameter I, F + Make an image which has its domain restricted to a function. + +Changed Properties: + category -> domain_morpher + io -> read_only if I is cont, trait::image_::io otherwise + data -> if the data property of I is linear then stored (FIXME why?) + else like I + + + +**** Image_if_interval + Image which domain is restricted by an interval + +Changed Property: + Same as image_if + + + +**** translate + Translate an Image with a delta point. +Changed Property: + category -> domain_morpher_category. + +*** Value morpher +**** cast_image +Make the user see the same image but with another data type. +Changed properties: + io -> read_only +**** value::stack_image + Represent a Stack of image. -Note: -The different grid (retrieve from point_ class): - - tick, - - square, - - hexa, - - cube +Changed properties: + category -> value_morpher + value -> value::vectorial + speed -> fast Index: sandbox/ballas/doc/draft.txt --- sandbox/ballas/doc/draft.txt (revision 1864) +++ sandbox/ballas/doc/draft.txt (working copy) @@ -315,6 +315,16 @@ the value. We can access directly to the value of an image image.value + values_access: none, + random, + browsing + + values_io: read, + write, + read_only < read,// This is dummy + write_only < write, + read_write < both read'n write +

1 0

1864: Update tests/tests.mk.
by Roland Levillain 14 Apr '08

14 Apr '08

https://svn.lrde.epita.fr/svn/oln/trunk/milena Index: ChangeLog from Roland Levillain <roland(a)lrde.epita.fr> Update tests/tests.mk. * tests/tests.mk (LONG_TESTS_CXXFLAGS): Remove. (TESTS_CXXFLAGS_SPEED): New. Set to @TESTS_CXXFLAGS_SPEED@. tests.mk | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) Index: tests/tests.mk --- tests/tests.mk (revision 1863) +++ tests/tests.mk (working copy) @@ -13,6 +13,4 @@ TESTS_CXXFLAGS = @TESTS_CXXFLAGS@ AM_CXXFLAGS = $(TESTS_CXXFLAGS) -# FIXME: Likewise, we should compute these values at configure time. -# Hard-code them for the moment. -LONG_TESTS_CXXFLAGS = -O3 +TESTS_CXXFLAGS_SPEED = @TESTS_CXXFLAGS_SPEED@

1 0

1863: Clean up tests/morpho/dilation a bit.
by Roland Levillain 14 Apr '08

14 Apr '08

https://svn.lrde.epita.fr/svn/oln/trunk/milena Index: ChangeLog from Roland Levillain <roland(a)lrde.epita.fr> Clean up tests/morpho/dilation a bit. * tests/morpho/dilation.cc: Reorganize. Disable neigh::image-related tests. Remove outdated tests. * tests/morpho/Makefile.am (CXXFLAGS_SPEED): Remove variable. s/CXXFLAGS_SPEED/TESTS_CXXFLAGS_SPEED/. (dilation_CXXFLAGS): New. Set to $(TESTS_CXXFLAGS_SPEED). Makefile.am | 25 ++++++++++++++++++------- dilation.cc | 55 +++++++++++++++++++++---------------------------------- 2 files changed, 39 insertions(+), 41 deletions(-) Index: tests/morpho/dilation.cc --- tests/morpho/dilation.cc (revision 1862) +++ tests/morpho/dilation.cc (working copy) @@ -48,9 +48,6 @@ #include <mln/pw/cst.hh> #include <mln/fun/ops.hh> -#include <mln/convert/to_p_array.hh> -#include <mln/convert/to_window.hh> - #include <mln/core/neighb2d.hh> #include <mln/neighb/image.hh> @@ -62,18 +59,14 @@ using namespace mln; using value::int_u8; + // FIXME: Do we really need such a high value? border::thickness = 66; - // FIXME: Maybe we should split all these tests cases into severals - // files. Sure it's a pain to create, but we want to be able to - // choose the granularity of the executed test suite (fast, - // comprehensive, etc.). - image2d<int_u8> lena; io::pgm::load(lena, MLN_IMG_DIR "/lena.pgm"); + win::rectangle2d rec(5, 3); { - win::rectangle2d rec(5, 3); image2d<int_u8> out = morpho::dilation(lena, rec); io::pgm::save(out, "out1.pgm"); } @@ -85,35 +78,29 @@ } { - image2d<int_u8> out = morpho::dilation(lena + c4()); - io::pgm::save(out, "out3.pgm"); + image2d<bool> bin(lena.domain()), out(lena.domain()); + level::fill(bin, pw::value(lena) > pw::cst(127u)); + morpho::dilation(bin, rec, out); + + image2d<int_u8> test(lena.domain()); + image2d<int_u8>::fwd_piter p(lena.domain()); + for_all(p) + test(p) = out(p) ? 255 : 0; + io::pgm::save(test, "out3.pgm"); } + /* FIXME: Re-enable these tests for Olena 1.1, when associated + neighborhoods (i.e., the neighb::image morpher) are supported and + shipped. */ +#if 0 { - image2d<int_u8> out = morpho::dilation(lena + c8()); + image2d<int_u8> out = morpho::dilation(lena + c4()); io::pgm::save(out, "out4.pgm"); } -// { -// p_array<point2d> vec = convert::to_p_array(rec, point2d::zero); -// window2d win = convert::to_window(vec); - -// image2d<int_u8> out(lena.domain()); -// level::ero(lena, win, out); -// morpho::dilation(lena, win, out); -// io::pgm::save(out, "out5.pgm"); -// } - -// { -// image2d<bool> bin(lena.domain()), out(lena.domain()); -// level::fill(bin, pw::value(lena) > pw::cst(127)); -// morpho::dilation(bin, rec, out); - -// image2d<int_u8> test(lena.domain()); -// image2d<int_u8>::fwd_piter p(lena.domain()); -// for_all(p) -// test(p) = out(p) ? 255 : 0; -// io::pgm::save(test, "out6.pgm"); -// } - + { + image2d<int_u8> out = morpho::dilation(lena + c8()); + io::pgm::save(out, "out5.pgm"); + } +#endif } Index: tests/morpho/Makefile.am --- tests/morpho/Makefile.am (revision 1862) +++ tests/morpho/Makefile.am (working copy) @@ -19,7 +19,11 @@ opening_area \ thinning -dilation_SOURCES = dilation.cc +# -------------- # +# Normal tests. # +# -------------- # + +# FIXME: Have erosion and dilation perform symmetric tests. erosion_SOURCES = erosion.cc dilation_max_h_SOURCES = dilation_max_h.cc @@ -36,14 +40,21 @@ meyer_wst_SOURCES = meyer_wst.cc -# FIXME: We should isolate those tests, as they take a long time to -# execute. -CXXFLAGS_SPEED = -O3 -ggdb -DNDEBUG +# --------------- # +# Complex tests. # +# --------------- # + +dilation_SOURCES = dilation.cc +dilation_CXXFLAGS = $(TESTS_CXXFLAGS_SPEED) + lena_line_graph_image_wst1_SOURCES = lena_line_graph_image_wst1.cc -lena_line_graph_image_wst1_CXXFLAGS = $(CXXFLAGS_SPEED) +lena_line_graph_image_wst1_CXXFLAGS = $(TESTS_CXXFLAGS_SPEED) + lena_line_graph_image_wst2_SOURCES = lena_line_graph_image_wst2.cc -lena_line_graph_image_wst2_CXXFLAGS = $(CXXFLAGS_SPEED) +lena_line_graph_image_wst2_CXXFLAGS = $(TESTSCXXFLAGS_SPEED) + meyer_wst_long_SOURCES = meyer_wst_long.cc -meyer_wst_long_CXXFLAGS = $(CXXFLAGS_SPEED) +meyer_wst_long_CXXFLAGS = $(TESTS_CXXFLAGS_SPEED) + TESTS = $(check_PROGRAMS)

1 0

1862: More on CXXFLAGS in configure.ac.
by Roland Levillain 14 Apr '08

14 Apr '08

https://svn.lrde.epita.fr/svn/oln/trunk Index: ChangeLog from Roland Levillain <roland(a)lrde.epita.fr> More on CXXFLAGS in configure.ac. * configure.ac (TESTS_CXXFLAGS_SPEED): New variable. (TOOLS_CXXFLAGS): Append `-Wall -W'. configure.ac | 15 ++++++++++++--- 1 file changed, 12 insertions(+), 3 deletions(-) Index: configure.ac --- configure.ac (revision 1861) +++ configure.ac (working copy) @@ -33,23 +33,32 @@ AC_LANG([C++]) AC_PROG_CXX -# Speed the compilation up. +# Speed up compiling times. if test "$GXX" = yes; then CXXFLAGS="$CXXFLAGS -pipe" fi # C++ compiler flags for tests. AC_ARG_VAR([TESTS_CXXFLAGS]) -# We want no optimization for the tests (too slow), and a lot of debugging. +# We want no optimization for the tests (it slows down compiling +# times), and a lot of debugging. if test "$GXX" = yes && test -z "$TESTS_CXXFLAGS"; then TESTS_CXXFLAGS="-O0 -ggdb -Wall -W" fi +# C++ compiler flags for complex tests. +AC_ARG_VAR([TESTS_CXXFLAGS_SPEED]) +# We want optimization for complex tests, and keep debugging flags +# (still useful). +if test "$GXX" = yes && test -z "$TESTS_CXXFLAGS_SPEED"; then + TESTS_CXXFLAGS_SPEED="-O3 -DNDEBUG -ggdb -Wall -W" +fi + # C++ compiler flags for tools. AC_ARG_VAR([TOOLS_CXXFLAGS]) # On the contrary, we want fast binaries for tools. if test "$GXX" = yes && test -z "$TOOLS_CXXFLAGS"; then - TOOLS_CXXFLAGS="-O3 -DNDEBUG -ggdb" + TOOLS_CXXFLAGS="-O3 -DNDEBUG -ggdb -Wall -W" fi

1 0

1861: Add a test on regular neighborhoods.
by Roland Levillain 14 Apr '08

14 Apr '08

https://svn.lrde.epita.fr/svn/oln/trunk Index: ChangeLog from Roland Levillain <roland(a)lrde.epita.fr> Add a test on regular neighborhoods. * milena/tests/core/neighb.cc: New. * milena/tests/core/Makefile.am (check_PROGRAMS): Add neighb. (neighb_SOURCES): New. Makefile.am | 2 ++ neighb.cc | 58 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 2 files changed, 60 insertions(+) Index: milena/tests/core/neighb.cc --- milena/tests/core/neighb.cc (revision 0) +++ milena/tests/core/neighb.cc (revision 0) @@ -0,0 +1,58 @@ +// Copyright (C) 2008 EPITA Research and Development Laboratory (LRDE) +// +// This file is part of the Olena Library. This library is free +// software; you can redistribute it and/or modify it under the terms +// of the GNU General Public License version 2 as published by the +// Free Software Foundation. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// General Public License for more details. +// +// You should have received a copy of the GNU General Public License +// along with this library; see the file COPYING. If not, write to +// the Free Software Foundation, 51 Franklin Street, Fifth Floor, +// Boston, MA 02111-1307, USA. +// +// As a special exception, you may use this file as part of a free +// software library without restriction. Specifically, if other files +// instantiate templates or use macros or inline functions from this +// file, or you compile this file and link it with other files to +// produce an executable, this file does not by itself cause the +// resulting executable to be covered by the GNU General Public +// License. This exception does not however invalidate any other +// reasons why the executable file might be covered by the GNU General +// Public License. + +/// \file tests/core/neighb.cc +/// \brief Tests on mln::neighb<D> specializations. + +#include <mln/core/neighb1d.hh> +#include <mln/core/neighb2d.hh> +#include <mln/core/neighb3d.hh> + +using namespace mln; + +template <typename N> +void test(const Neighborhood<N>& nbh, const mln_point(N)& p_ref, + unsigned size) +{ + mln_niter(N) n(nbh, p_ref); + unsigned i = 0; + for_all(n) + ++i; + mln_assertion(i == size); +} + +int main() +{ + test(c2(), make::point1d(0), 2); + + test(c4(), make::point2d(0, 0), 4); + test(c8(), make::point2d(0, 0), 8); + + test( c6(), make::point3d(0, 0, 0), 6); + test(c18(), make::point3d(0, 0, 0), 18); + test(c26(), make::point3d(0, 0, 0), 26); +} Index: milena/tests/core/Makefile.am --- milena/tests/core/Makefile.am (revision 1860) +++ milena/tests/core/Makefile.am (working copy) @@ -17,6 +17,7 @@ line_graph_image_wst \ mono_obased_rle_image \ mono_rle_image \ + neighb \ obased_rle_image \ p_priority_queue \ p_priority_queue_fast \ @@ -44,6 +45,7 @@ line_graph_image_wst_SOURCES = line_graph_image_wst.cc mono_obased_rle_image_SOURCES = mono_obased_rle_image.cc mono_rle_image_SOURCES = mono_rle_image.cc +neighb_SOURCES = neighb.cc obased_rle_image_SOURCES = obased_rle_image.cc p_priority_queue_SOURCES = p_priority_queue.cc p_priority_queue_fast_SOURCES = p_priority_queue_fast.cc

1 0