https://svn.lrde.epita.fr/svn/oln/branches/cleanup-2008/milena/sandbox
Index: ChangeLog
from Alexandre Abraham <abraham(a)lrde.epita.fr>
Add bijective functors, 2-way functors and meta-fun.
* abraham/tests/morpho/test_component_tree.cc: Remove
This test was outdated.
* abraham/tests/morpho/test_watershed.cc:
Clean.
* abraham/tests/morpho/test_watershed_topo.cc:
Clean.
* abraham/tests/fun/meta: New.
* abraham/tests/fun/meta/red.cc: New
(mln::meta::red) : meta function to access red in rgb.
* abraham/tests/core/image/thru_norm.cc: New
test for norm functor.
* abraham/tests/core/image/thru_v2v.cc: New
test for abs functor.
* abraham/tests/core/image/thru.cc: Remove.
* abraham/tests/core/image/thru_const.cc: New
test of writing in a const image.
* abraham/tests/core/image/thru_v2w2v.cc: New
test for cosinus functor.
* abraham/tests/core/concept: New.
* abraham/tests/core/concept/test.cc: New.
* abraham/mln/morpho/najman_wst.hh: New
M-W watershed.
* abraham/mln/morpho/basic_najman.hh:
Obsolete, will eventually be removed.
* abraham/mln/morpho/topo_wst.hh: New
1 pass topological watershed.
* abraham/mln/core/image/thru.hh:
(mln::thru): see an image throught a function.
* abraham/mln/core/image/shell.hh: New
(mln::shell): value shell.
* abraham/mln/core/concept/function.hh:
Add bijective and 2-way function.
* abraham/mln/core/concept/meta_fun.hh: New
Add concept for meta function.
* abraham/mln/fun/meta: New.
* abraham/mln/fun/meta/red.hh: New
red fun object.
* abraham/mln/fun/v2w_w2v/norm.hh:
2-way functor.
* vigouroux/tests.cc: Rename image_if_value > image_if.
* vigouroux/color/rgb_to_cmy.hh: Rename image_if_value > image_if.
* vigouroux/color/rgb_to_xyz.hh: Rename image_if_value > image_if.
* vigouroux/color/rgb_to_hsi.hh: Rename image_if_value > image_if.
* vigouroux/color/rgb_to_yuv.hh: Rename image_if_value > image_if.
abraham/mln/core/concept/function.hh | 2
abraham/mln/core/concept/meta_fun.hh | 83 ++
abraham/mln/core/image/shell.hh | 118 ++++
abraham/mln/core/image/thru.hh | 40 -
abraham/mln/fun/meta/red.hh | 45 +
abraham/mln/fun/v2w_w2v/norm.hh | 4
abraham/mln/morpho/basic_najman.hh | 196 +++---
abraham/mln/morpho/najman_wst.hh | 791 ++++++++++++++++++++++++++++
abraham/mln/morpho/topo_wst.hh | 744 ++++++++++++++++++++++++++
abraham/tests/core/concept/test.cc | 75 ++
abraham/tests/core/image/thru.cc | 56 -
abraham/tests/core/image/thru_const.cc | 64 ++
abraham/tests/core/image/thru_norm.cc | 58 ++
abraham/tests/core/image/thru_v2v.cc | 56 +
abraham/tests/core/image/thru_v2w2v.cc | 64 ++
abraham/tests/fun/meta/red.cc | 60 ++
abraham/tests/morpho/test_component_tree.cc | 230 --------
abraham/tests/morpho/test_watershed.cc | 60 --
abraham/tests/morpho/test_watershed_topo.cc | 8
vigouroux/color/rgb_to_cmy.hh | 8
vigouroux/color/rgb_to_hsi.hh | 2
vigouroux/color/rgb_to_xyz.hh | 14
vigouroux/color/rgb_to_yuv.hh | 8
vigouroux/tests.cc | 6
24 files changed, 2311 insertions(+), 481 deletions(-)
Index: abraham/tests/morpho/test_watershed.cc
--- abraham/tests/morpho/test_watershed.cc (revision 2672)
+++ abraham/tests/morpho/test_watershed.cc (working copy)
@@ -7,7 +7,7 @@
#include <mln/io/pgm/load.hh>
#include <mln/io/pgm/save.hh>
-#include <mln/morpho/basic_najman.hh>
+#include <mln/morpho/najman_wst.hh>
#include <string>
#include <iostream>
@@ -22,15 +22,15 @@
// #define TEST
- io::pgm::load(input, "./images/test_watershed.pgm");
+ // io::pgm::load(input, "./images/test_watershed.pgm");
// io::pgm::load(input, "./images/little_test.pgm");
- // io::pgm::load(input, "./images/test_2.pgm");
+ io::pgm::load(input, "./images/test.pgm");
// io::pgm::load(input, "./images/lena_light.pgm");
// io::pgm::load(mverif, "./images/result_m_watershed.pgm");
// io::pgm::load(wverif, "./images/result_topo_watershed.pgm");
- morpho::basic_najman< image2d<int_u8>, neighb2d> n(input, c4());
+ morpho::najman_wst< image2d<int_u8>, neighb2d> n(input, c4());
/*
image2dint::fwd_piter it(input.domain());
@@ -43,60 +43,8 @@
// io::tikz::save(input, "start.tex");
- std::cout << "Building Component tree..." << std::endl;
-
n.go();
- std::cout << "M-Watershed" << std::endl;
-
- n.m_watershed();
-
- io::pgm::save(n.pima, "int.pgm");
-
- // io::tikz::save(n.pima, "step.tex");
-
-#ifdef TEST
-
- 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;
-
-#endif
-
- std::cout << "W-watershed" << std::endl;
-
- n.w_watershed();
- // io::pgm::save(n.pima, "int.pgm");
- // std::cout << "~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~" <<
std::endl;
- // n.w_watershed();
-
- // for_all(it)
- // n.pima(it) = n.pima(it) * 17;
-
-#ifdef TEST
- equal = n.pima == wverif;
-
- if (!equal)
- {
- std::cout << "W-Watershed broken" << std::endl;
- return 1;
- }
- else
- std::cout << "W-Watershed good as chocolate ice cream" <<
std::endl;
-#endif
-
- // io::tikz::save(n.pima, "end.tex");
-
- // n.image_output(n.pima);
-
io::pgm::save(n.pima, "out.pgm");
return 0;
Index: abraham/tests/morpho/test_watershed_topo.cc
--- abraham/tests/morpho/test_watershed_topo.cc (revision 2672)
+++ abraham/tests/morpho/test_watershed_topo.cc (working copy)
@@ -9,7 +9,7 @@
#include <mln/io/pgm/load.hh>
#include <mln/io/pgm/save.hh>
-#include <mln/morpho/basic_najman.hh>
+#include <mln/morpho/topo_wst.hh>
#include <string>
#include <iostream>
@@ -32,7 +32,7 @@
// io::pgm::load(input, "./images/test_watershed.pgm");
// io::pgm::load(input, "./images/little_test.pgm");
- io::pgm::load(input, "./images/test_4.pgm");
+ io::pgm::load(input, "./images/test.pgm");
// io::pgm::load(input, "../../img/dots.pgm");
//io::pgm::load(input, "./images/+irm6.pgm");
@@ -40,7 +40,7 @@
// io::pgm::load(mverif, "./images/result_m_watershed.pgm");
// io::pgm::load(wverif, "./images/result_topo_watershed.pgm");
- morpho::basic_najman< image2d<int_u8>, neighb2d> n(input, c4());
+ morpho::topo_wst< image2d<int_u8>, neighb2d> n(input, c4());
/*
image2dint::fwd_piter it(input.domain());
@@ -53,7 +53,7 @@
// io::tikz::save(input, "start.tex");
- n.gotopo();
+ n.go();
io::pgm::save(n.pima, "out.pgm");
Index: abraham/tests/fun/meta/red.cc
--- abraham/tests/fun/meta/red.cc (revision 0)
+++ abraham/tests/fun/meta/red.cc (revision 0)
@@ -0,0 +1,60 @@
+#include <mln/fun/meta/red.hh>
+#include <mln/core/image/image2d.hh>
+#include <mln/core/image/thru.hh>
+
+namespace mln
+{
+
+ template <class T>
+ struct rgb
+ {
+ T red() const { return r; }
+ T& red() { return r; }
+ T r;
+ };
+
+ template <class T>
+ struct function< meta::red, rgb<T> > : public
Function_v2w_w2v<function< meta::red, rgb<T> > >
+ {
+ typedef T result;
+ T read(const rgb<T>& c)
+ {
+ std::cout << "read red rgb<T>" << std::endl;
+ return c.red();
+ }
+
+ typedef T& lresult;
+ T& write(rgb<T>& c)
+ {
+ std::cout << "write red rgb<T>" << std::endl;
+ return c.red();
+ }
+ };
+}
+
+int main ()
+{
+ typedef mln::rgb<int> C;
+ mln::image2d<C> i(3, 2, 0);
+ C c;
+ c.r = 51;
+ i(mln::point2d(0,0)) = c;
+ c.r = 23;
+ i(mln::point2d(0,1)) = c;
+ c.r = 43;
+ i(mln::point2d(1,0)) = c;
+ c.r = 0;
+ i(mln::point2d(1,1)) = c;
+ c.r = 65;
+ i(mln::point2d(2,0)) = c;
+ c.r = 1;
+ i(mln::point2d(2,1)) = c;
+
+ mln::thru<mln::meta::red, mln::image2d<C> > out(i);
+
+ for_all (p)
+ std::cout << out(p) << std::endl;
+
+
+}
+
Index: abraham/tests/core/image/thru_norm.cc
--- abraham/tests/core/image/thru_norm.cc (revision 0)
+++ abraham/tests/core/image/thru_norm.cc (revision 0)
@@ -0,0 +1,58 @@
+// 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.
+
+# include <mln/core/image/image2d.hh>
+# include <mln/core/image/thru.hh>
+# include <mln/fun/v2w_w2v/norm.hh>
+# include <mln/algebra/vec.hh>
+# include <mln/level/fill.hh>
+# include <mln/core/image/violent_cast_image.hh>
+
+
+int main ()
+{
+
+ using namespace mln;
+ typedef image2d<algebra::vec<2, double> > I;
+
+ I ima(3, 2, 0);
+
+ ima(point2d(0,0)).set (1, 1);
+ ima(point2d(0,1)).set (1, 3);
+ ima(point2d(1,0)).set (4, 4);
+ ima(point2d(1,1)).set (-1, 3);
+ ima(point2d(2,0)).set (23, 23);
+ ima(point2d(2,1)).set (3, 1);
+
+ thru<mln::fun::v2w_w2v::l1_norm<algebra::vec<2, double>, double>, I >
out(ima);
+ level::fill(out, 1);
+
+ box_fwd_piter_<point2d> p(ima.domain());
+
+ for_all (p)
+ std::cout << ima(p) << std::endl;
+}
Index: abraham/tests/core/image/thru_v2v.cc
--- abraham/tests/core/image/thru_v2v.cc (revision 0)
+++ abraham/tests/core/image/thru_v2v.cc (revision 0)
@@ -0,0 +1,56 @@
+// 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.
+
+# include <mln/core/image/image2d.hh>
+# include <mln/core/image/thru.hh>
+# include <mln/fun/v2v/abs.hh>
+
+int main ()
+{
+
+ using namespace mln;
+
+ typedef image2d<int> I;
+
+ int vs[6][5] = {
+
+ { -3, -3, -4, -4, -4 },
+ { -2, -1, -1, -1, -1 },
+ { -1, -4, -4, -4, -1 },
+ { -1, -4, -3, -4, -1 },
+ { -1, -4, -5, -3, -1 },
+ { -1, -1, -1, -1, -1 }
+
+ };
+
+ image2d<int> ima(make::image2d(vs));
+ thru<mln::fun::v2v::abs<int>, image2d<int> > out(ima);
+
+ box_fwd_piter_<point2d> p(ima.domain());
+ for_all (p)
+ mln_assertion (out(p) >= 0);
+}
Index: abraham/tests/core/image/thru_const.cc
--- abraham/tests/core/image/thru_const.cc (revision 0)
+++ abraham/tests/core/image/thru_const.cc (revision 0)
@@ -0,0 +1,64 @@
+// 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.
+
+# include <mln/core/image/image2d.hh>
+# include <mln/core/image/thru.hh>
+# include <mln/fun/v2w2v/cos.hh>
+
+int main ()
+{
+
+ using namespace mln;
+
+ typedef const image2d<double> I;
+
+ double vs[6][5] = {
+
+ { 12, -3, 8, -4, 6 },
+ { -2, 22, -1, 45, -1 },
+ { -1, -4, -4, -4, -1 },
+ { -1, -4, -3, -4, -1 },
+ { -1, -4, -5, -3, -1 },
+ { -1, -1, -1, -1, -1 }
+
+ };
+
+ I ima(make::image2d(vs));
+ thru<mln::fun::v2w2v::cos<double>, I > out(ima);
+
+ double i = 0;
+
+ box_fwd_piter_<point2d> p(ima.domain());
+ for_all (p)
+ {
+ out(p) = i;
+ i += 1./40.;
+ }
+
+ for_all (p)
+ std::cout << out(p) << std::endl;
+}
Index: abraham/tests/core/image/thru_v2w2v.cc
--- abraham/tests/core/image/thru_v2w2v.cc (revision 0)
+++ abraham/tests/core/image/thru_v2w2v.cc (revision 0)
@@ -0,0 +1,64 @@
+// 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.
+
+# include <mln/core/image/image2d.hh>
+# include <mln/core/image/thru.hh>
+# include <mln/fun/v2w2v/cos.hh>
+
+int main ()
+{
+
+ using namespace mln;
+
+ typedef image2d<double> I;
+
+ double vs[6][5] = {
+
+ { 12, -3, 8, -4, 6 },
+ { -2, 22, -1, 45, -1 },
+ { -1, -4, -4, -4, -1 },
+ { -1, -4, -3, -4, -1 },
+ { -1, -4, -5, -3, -1 },
+ { -1, -1, -1, -1, -1 }
+
+ };
+
+ image2d<double> ima(make::image2d(vs));
+ thru<mln::fun::v2w2v::cos<double>, image2d<double> > out(ima);
+
+ double i = 0;
+
+ box_fwd_piter_<point2d> p(ima.domain());
+ for_all (p)
+ {
+ out(p) = i;
+ i += 1./40.;
+ }
+
+ for_all (p)
+ std::cout << out(p) << std::endl;
+}
Index: abraham/tests/core/concept/test.cc
--- abraham/tests/core/concept/test.cc (revision 0)
+++ abraham/tests/core/concept/test.cc (revision 0)
@@ -0,0 +1,75 @@
+#include <iostream>
+#include <mln/core/concept/meta_fun.hh>
+
+using namespace mln;
+
+meta::red red; // fun-obj
+
+
+
+template <class T>
+struct rgb
+{
+ T red() const { return r; }
+ T& red() { return r; }
+ T r;
+};
+
+
+
+template <class T>
+struct fun< meta::red, rgb<T> > : Function_v2v< fun< meta::red,
rgb<T> > >
+{
+ typedef T result;
+ T read(const rgb<T>& c)
+ {
+ std::cout << "read red rgb<T>" << std::endl;
+ return c.red();
+ }
+
+ typedef T& lresult;
+ T& write(rgb<T>& c)
+ {
+ std::cout << "write red rgb<T>" << std::endl;
+ return c.red();
+ }
+};
+
+
+/*
+
+ morpher apply(M m, I ima)
+ {
+ M et value(I) -> F -> nature
+ }
+
+ */
+
+
+
+int main()
+{
+ typedef rgb<int> C;
+ C c;
+ c.r = 51;
+
+ int r = red(c);
+ std::cout << r << std::endl;
+
+
+ // norm(v) = 1;
+
+
+ // norm( ima(p) ) = n;
+
+ // norm( v& ) = n;
+ // v = n * v / norm(v);
+
+ // norm( v& ) = 1;
+
+
+
+ // ima' = apply(norm, ima)
+ // ima'(p) = 1;
+
+}
Index: abraham/mln/morpho/najman_wst.hh
--- abraham/mln/morpho/najman_wst.hh (revision 0)
+++ abraham/mln/morpho/najman_wst.hh (revision 0)
@@ -0,0 +1,791 @@
+// 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.
+
+#ifndef MLN_MORPHO_NAJMAN_WST_HH
+# define MLN_MORPHO_NAJMAN_WST_HH
+
+
+
+#include <mln/level/sort_psites.hh>
+#include <mln/level/fill.hh>
+#include <mln/core/image/image2d.hh>
+#include <mln/core/site_set/p_set.hh>
+#include <mln/estim/min_max.hh>
+#include <mln/math/sqr.hh>
+#include <mln/util/greater_psite.hh>
+#include <mln/util/ord.hh>
+
+#include <mln/core/site_set/p_priority.hh>
+#include <mln/core/site_set/p_queue_fast.hh>
+
+#include <queue>
+#include <vector>
+#include <map>
+
+namespace mln
+{
+ namespace morpho
+ {
+
+ template <class I, class N>
+ struct najman_wst
+ {
+
+ // Typedefs
+ // --------
+ typedef mln_site(I) site;
+
+ // Component tree management
+ // -------------------------
+ private:
+ struct node {
+ mln_value(I) level;
+ int area;
+ int highest;
+ typedef mln_site(I) site;
+ // Can't modify the sites in a p_array
+ // p_array<mln_site(I)> children;
+ std::vector<site> children;
+
+ void addChildren(const node& n)
+ {
+ // typename p_array<mln_site(I)>::fwd_piter it(n.children);
+ // for (it.start();
+ // it.is_valid();
+ // it.next())
+ // children.append(it.to_site());
+ for (unsigned i=0; i < n.children.size(); ++i)
+ children.push_back(n.children[i]);
+ }
+
+ void addChild(const site p)
+ {
+ // children.append(n);
+ children.push_back(p);
+ }
+ }; // struct node
+
+ const mln_exact(N)& nbh;
+ mln_ch_value(I, site) Par_node;
+ mln_ch_value(I, site) Par_tree;
+
+ mln_ch_value(I, int) Rnk_tree;
+ mln_ch_value(I, int) Rnk_node;
+ mln_ch_value(I, site) subtreeRoot;
+ mln_ch_value(I, node) nodes;
+ site Root;
+
+ void MakeSet_tree(site x);
+ void MakeSet_node(site x);
+ site Find_tree(site x);
+ void swap(site& x, site& y);
+ site Find_node(site x);
+ void BuildComponentTree();
+ site MergeNode(site& node1, site& node2);
+ site Link_tree(site x, site y);
+ site Link_node(site x, site y);
+ node MakeNode(int level);
+
+
+ // Watershed algorithm
+ // -------------------
+
+ private:
+ mln_ch_value(I, bool) isproc;
+
+ // Ctor
+ public:
+ najman_wst(const Image<I>& i,
+ const Neighborhood<N>& n);
+
+
+ // Result
+ public:
+ mln_exact(I) pima;
+
+ public:
+ void go();
+
+ private:
+ site min (p_set<site>& components);
+ site max (p_set<site>& components);
+ bool highest_fork (p_set<site>& components);
+ bool highest_fork (p_set<site>& components, site &r);
+ bool m_destructible(site p);
+ bool w_destructible(site p);
+ bool m_destructible(site p, site &r);
+ bool w_destructible(site p, site &r);
+ void m_watershed ();
+ void w_watershed ();
+
+ // Optimized LCA Algorithm
+ private:
+ int *euler;
+ int *depth;
+ int ctree_size;
+ std::map<site, int, mln::util::ord<site> > pos;
+ site *repr;
+
+ int *minim;
+ int **Minim;
+
+ void compute_ctree_size (site p);
+ void build_euler_tour_rec(site p, int &position, int d);
+ void build_euler_tour();
+ void build_minim ();
+ site lca (site a, site b);
+ void removeOneSonNodes(site *p, mln_ch_value(I, site) &newPar_node);
+ void compressTree();
+ }; // struct najman_wst
+
+# ifndef MLN_INCLUDE_ONLY
+
+ // Ctor
+ template <class I, class N>
+ najman_wst<I, N>::najman_wst(const Image<I>& i,
+ const Neighborhood<N>& n)
+ : pima(exact(i)),
+ nbh(exact(n)),
+ Par_node(exact(i).domain(), exact(i).border()),
+ Par_tree(exact(i).domain(), exact(i).border()),
+ Rnk_tree(exact(i).domain(), exact(i).border()),
+ Rnk_node(exact(i).domain(), exact(i).border()),
+ subtreeRoot(exact(i).domain(), exact(i).border()),
+ nodes(exact(i).domain(), exact(i).border()),
+ isproc(exact(i).domain(), exact(i).border())
+ {
+ }
+
+ template <class I, class N>
+ void najman_wst<I, N>::MakeSet_tree(site x)
+ {
+ Par_tree(x) = x;
+ Rnk_tree(x) = 0;
+ }
+
+ template <class I, class N>
+ void najman_wst<I, N>::MakeSet_node(site x)
+ {
+ Par_node(x) = x;
+ Rnk_node(x) = 0;
+ }
+
+ template <class I, class N>
+ typename najman_wst<I, N>::site najman_wst<I, N>::Find_tree(site x)
+ {
+ if (Par_tree(x) != x)
+ Par_tree(x) = Find_tree(Par_tree(x));
+ return Par_tree(x);
+ }
+
+ template <class I, class N>
+ void najman_wst<I, N>::swap(site& x, site& y)
+ {
+ site memo = x;
+ x = y;
+ y = memo;
+ }
+
+ template <class I, class N>
+ typename najman_wst<I, N>::site najman_wst<I, N>::Find_node(site x)
+ {
+ if (Par_node(x) != x)
+ Par_node(x) = Find_node(Par_node(x));
+ return Par_node(x);
+ }
+
+ template <class I, class N>
+ void najman_wst<I, N>::go()
+ {
+ BuildComponentTree();
+ compressTree();
+ build_euler_tour();
+ build_minim();
+ m_watershed();
+ w_watershed();
+ }
+
+ template <class I, class N>
+ void najman_wst<I, N>::BuildComponentTree()
+ {
+ // Init:
+
+ // Sort the sites in increasing order
+ p_array<mln_site(I)> S;
+ S = level::sort_psites_increasing(pima);
+
+ // Clear the marker map
+ level::fill(isproc, false);
+ for (int ip = 0; ip < int(S.nsites()); ++ip)
+ {
+ site p = S[ip];
+ MakeSet_node(p);
+ MakeSet_tree(p);
+ // if (subtreeRoot.hold(p))
+ subtreeRoot(p) = p;
+ // if (nodes.hold(p))
+ nodes(p) = MakeNode(pima(p));
+ }
+
+
+
+ typename p_array<mln_site(I)>::fwd_piter ip(S);
+ for_all(ip)
+ {
+ site p = ip.to_site();
+
+ site curCanonicalElt = Find_tree(p);
+ site curNode = Find_node(subtreeRoot(curCanonicalElt));
+
+ mln_niter(N) q(nbh, ip);
+ for_all(q)
+ if (pima.has(q) and isproc(q) and pima(q) <= pima(p))
+ {
+ site adjCanonicalElt = Find_tree(q);
+ site adjNode = Find_node(subtreeRoot(adjCanonicalElt));
+ if (curNode != adjNode)
+ {
+ if (nodes(curNode).level == nodes(adjNode).level)
+ curNode = MergeNode(adjNode, curNode);
+ else
+ {
+ nodes(curNode).addChild(adjNode);
+ nodes(curNode).area += nodes(adjNode).area;
+ nodes(curNode).highest += nodes(adjNode).highest;
+ }
+ }
+
+ curCanonicalElt = Link_tree(adjCanonicalElt, curCanonicalElt);
+ subtreeRoot(curCanonicalElt) = curNode;
+ }
+ isproc(p) = true;
+ }
+ // Pour garder une map de correspondance site <-> local_root
+ // for (int ip = 0; ip < int(S.size()); ++ip)
+ // {
+ // site p = S[ip];
+ // M(p) = Find_node(p);
+ // }
+
+ mln_piter(I) r(Par_node.domain());
+ for_all(r)
+ Par_node(r) = Find_node(r);
+
+ Root = subtreeRoot(Find_tree(Find_node(site(0, 0))));
+ }
+
+
+ template <class I, class N>
+ typename najman_wst<I, N>::site najman_wst<I, N>::MergeNode(site&
node1, site& node2)
+ {
+ site tmpNode = Link_node(node1, node2);
+ site tmpNode2;
+ if (tmpNode == node2)
+ {
+ nodes(node2).addChildren(nodes(node1));
+ tmpNode2 = node1;
+ }
+ else
+ {
+ nodes(node1).addChildren(nodes(node2));
+ tmpNode2 = node2;
+ }
+ nodes(tmpNode).area += nodes(tmpNode2).area;
+ nodes(tmpNode).highest += nodes(tmpNode2).highest;
+ return tmpNode;
+ }
+
+ template <class I, class N>
+ typename najman_wst<I, N>::site najman_wst<I, N>::Link_tree(site x, site
y)
+ {
+ if (Rnk_tree(x) > Rnk_tree(y))
+ swap(x, y);
+ else
+ if (Rnk_tree(x) == Rnk_tree(y))
+ Rnk_tree(y) += 1;
+ Par_tree(x) = y;
+ return y;
+ }
+
+ template <class I, class N>
+ typename najman_wst<I, N>::site najman_wst<I, N>::Link_node(site x, site
y)
+ {
+ if (Rnk_node(x) > Rnk_node(y))
+ swap(x, y);
+ else
+ if (Rnk_node(x) == Rnk_node(y))
+ Rnk_node(y) += 1;
+ Par_node(x) = y;
+ return y;
+ }
+
+ template <class I, class N>
+ typename najman_wst<I, N>::node najman_wst<I, N>::MakeNode(int level)
+ {
+ node n;
+ n.level = level;
+ n.area = 1;
+ n.highest = level;
+ return n;
+ }
+
+ template <class I, class N>
+ typename najman_wst<I, N>::site najman_wst<I, N>::min
(p_set<site>& components)
+ {
+ if (components.nsites() == 0)
+ return site(-1, -1);
+
+ typename p_set<site>::fwd_piter it(components);
+ site min = components[0];
+
+ for_all(it)
+ if (pima(it.to_site()) < pima(min))
+ min = it.to_site();
+
+ return min;
+ }
+
+ template <class I, class N>
+ typename najman_wst<I, N>::site najman_wst<I, N>::max
(p_set<site>& components)
+ {
+ if (components.nsites() == 0)
+ return site(-1, -1);
+
+ typename p_set<site>::fwd_piter it(components);
+ site max = components[0];
+
+ for_all(it)
+ if (pima(it.to_site()) > pima(max))
+ max = it.to_site();
+
+ return max;
+ }
+
+
+ template <class I, class N>
+ bool najman_wst<I, N>::highest_fork (p_set<site>& components, site
&r)
+ {
+ if (components.nsites() == 0)
+ {
+ std::cerr << "highest fork : empty set" << std::endl;
+ return false;
+ }
+
+ site
+ m = min(components),
+ hfork = m;
+
+ typename p_set<site>::fwd_piter it(components);
+ for_all(it)
+ {
+ // Can't remove the site from the set
+ if (it.to_site() == m)
+ continue;
+
+ site c = lca(hfork, it.to_site());
+ if (c != it.to_site())
+ hfork = c;
+ }
+
+ if (nodes(m).level == nodes(hfork).level)
+ return false;
+
+ r = hfork;
+ return true;
+ }
+
+ template <class I, class N>
+ bool najman_wst<I, N>::w_destructible(site p, site &r)
+ {
+ mln_niter(N) q(nbh, p);
+ p_set<site> v;
+
+ for_all(q)
+ if (pima.has(q) && pima(q) < pima(p))
+ v.insert(Par_node(q));
+
+ if (v.nsites() == 0)
+ return false;
+
+ site hf;
+ bool is_hf = highest_fork(v, hf);
+
+ if (!is_hf) {
+ r = min(v);
+ return true;
+ }
+
+ if (nodes(hf).level <= pima(p)) {
+ r = hf;
+ return true;
+ }
+ return false;
+ }
+
+ template <class I, class N>
+ bool najman_wst<I, N>::highest_fork (p_set<site>& components) {
+ site i;
+ return highest_fork(components, i);
+ }
+
+ template <class I, class N>
+ bool najman_wst<I, N>::m_destructible(site p) {
+ site i;
+ return m_destructible(p, i);
+ }
+
+ template <class I, class N>
+ bool najman_wst<I, N>::w_destructible(site p) {
+ site i;
+ return w_destructible(p, i);
+ }
+
+ template <class I, class N>
+ bool najman_wst<I, N>::m_destructible(site p, site &r)
+ {
+ mln_niter(N) q(nbh, p);
+ p_set<site> v = p_set<site>();
+
+ for_all(q)
+ if (pima.has(q) && pima(q) < pima(p))
+ v.insert(Par_node(q));
+
+ if (v.nsites() == 0)
+ return false;
+
+ // if (nodes(min(v)).children.nsites() != 0)
+ if (nodes(min(v)).children.size() != 0)
+ return false;
+
+ bool is_hf = highest_fork(v);
+
+ if (!is_hf) {
+ r = min(v);
+ return true;
+ }
+
+ return false;
+ }
+
+
+ template <class I, class N>
+ void najman_wst<I, N>::m_watershed ()
+ {
+ typedef
+ // p_priority_queue<site, mln_value(I) >
+ std::priority_queue<site, std::vector<site>, util::greater_psite<I> >
+ ordered_queue_type;
+
+
+ ordered_queue_type l(util::make_greater_psite(pima));
+
+ // Clear the marker map
+ level::fill(isproc, false);
+ mln_piter(I) it(pima.domain());
+
+ for_all(it)
+ if (m_destructible(it.to_site()))
+ {
+ l.push(it.to_site());
+ isproc(it.to_site()) = true;
+ }
+
+ site p, i;
+ while (!l.empty())
+ {
+ // Extract priority queue
+ p = l.top();
+ l.pop();
+
+ // unmark p
+ isproc(p) = false;
+
+ bool is_m = m_destructible(p, i);
+
+ if (is_m)
+ {
+ pima(p) = nodes(i).level;
+ Par_node(p) = i;
+
+ mln_niter(N) q(nbh, p);
+ for_all(q)
+ if (pima.has(q) && !isproc(q))
+ if (m_destructible(q))
+ {
+ // Add priority queue
+ l.push(q);
+
+ // mark q
+ isproc(q) = true;
+ }
+ }
+ }
+ }
+
+
+ template <class I, class N>
+ void najman_wst<I, N>::w_watershed ()
+ {
+ mln::p_priority< mln_value(I), p_queue_fast<site> > L;
+
+ mln_value(I) max = level::compute(accu::meta::max(), pima);
+
+ // I K(pima.domain(), pima.border());
+ mln_ch_value(I, unsigned) K(pima.domain(), pima.border());
+ mln_ch_value(I, site) H(pima.domain(), pima.border());
+
+ // For All p of E do
+ mln_piter(I) it(pima.domain());
+ for_all(it)
+ {
+ site p = it.to_site();
+ site i;
+
+ // i <- W-Destructible(p)
+ bool is_w = w_destructible(p, i);
+
+ // If i != infinity then
+ if (is_w)
+ {
+ L.push(max - nodes(i).level, p);
+ K(p) = nodes(i).level;
+ H(p) = i;
+ }
+ }
+
+ while (!L.is_empty())
+ {
+ mln_value(I) k = max - L.highest_priority();
+
+ site p = L.front();
+ L.pop();
+
+ if (K(p) == k)
+ {
+ pima(p) = k;
+ Par_node(p) = H(p);
+
+ mln_niter(N) q(nbh, p);
+ for_all(q)
+ if (pima.has(q) && k < pima(q))
+ {
+ site i;
+ bool is_w = w_destructible(q, i);
+ if (!is_w)
+ K(q) = 10000; // FIXME : supposed to be infinity...
+ else
+ if (K(q) != nodes(i).level)
+ {
+ L.push(max - nodes(i).level, q);
+ K(q) = nodes(i).level;
+ H(q) = i;
+ }
+ }
+ }
+ }
+ }
+
+
+ template <class I, class N>
+ void najman_wst<I, N>::compute_ctree_size (site p)
+ {
+ ctree_size += 1;
+ node& n = nodes(p);
+
+ // typename p_array<mln_site(I)>::fwd_piter it(n.children);
+ // for_all(it)
+ // compute_ctree_size(it.to_site());
+
+ for (unsigned i=0; i < n.children.size(); ++i)
+ compute_ctree_size(n.children[i]);
+ }
+
+
+ template <class I, class N>
+ void najman_wst<I, N>::build_euler_tour_rec(site p, int &position, int d)
+ {
+ if (pos.find(p) == pos.end())
+ pos[p] = position;
+
+ repr[position] = p;
+ depth[position] = d;
+ euler[position] = pos[p];
+ ++position;
+ node& n = nodes(p);
+
+ // typename p_array<mln_site(I)>::fwd_piter it(n.children);
+ // for_all(it)
+ // {
+ // build_euler_tour_rec(it.to_site(), position, d+1);
+ // depth[position] = d; // Not optimized
+ // euler[position] = pos[p];
+ // repr[position] = p; // Pas necessaire?
+ // ++position;
+ // }
+
+ for(unsigned i=0; i < n.children.size(); ++i)
+ {
+ build_euler_tour_rec(n.children[i], position, d+1);
+ depth[position] = d; // Not optimized
+ euler[position] = pos[p];
+ repr[position] = p; // Pas necessaire?
+ ++position;
+ }
+ }
+
+
+ template <class I, class N>
+ void najman_wst<I, N>::build_euler_tour ()
+ {
+ ctree_size = 0;
+ compute_ctree_size(Root);
+
+ int size = 2 * ctree_size - 1;
+
+ // FIXME : free this
+ euler = new int[size];
+ depth = new int[size];
+ repr = new site[size];
+
+ int position = 0;
+ build_euler_tour_rec(Root, position, 0);
+ }
+
+
+ template <class I, class N>
+ void najman_wst<I, N>::build_minim ()
+ {
+ // compute_tree_size(); // Already done in build_euler_tour
+ int size = 2 * ctree_size - 1;
+ int logn = (int)(ceil(log((double)(size))/log(2.0)));
+ // minim.reserve(size);
+ minim = new int [logn * size]; // FIXME : Think about freeing this
+ Minim = new int* [logn];
+
+ Minim[0] = minim;
+
+ // Init
+ for (int i = 0; i < size - 1; ++i)
+ if (depth[euler[i]] < depth[euler[i+1]]) {
+ Minim[0][i] = i;
+ } else {
+ Minim[0][i] = i+1;
+ }
+
+ // Minim[0][size - 1] = size - 1;
+
+ int k;
+ for (int j = 1; j < logn; j++) {
+ k = 1 << (j - 1);
+ Minim[j] = &minim[j * size];
+ for (int i = 0; i < size; i++) {
+ if ((i + (k << 1)) >= size) {
+ //Minim[j][i] = size - 1;
+ }
+ else {
+ if (depth[euler[Minim[j - 1][i]]] <= depth[euler[Minim[j - 1][i + k]]])
+ Minim[j][i] = Minim[j - 1][i];
+ else
+ Minim[j][i] = Minim[j - 1][i + k];
+ }
+ }
+ }
+
+ } // void build_minim ()
+
+
+ template <class I, class N>
+ typename najman_wst<I, N>::site najman_wst<I, N>::lca (site a, site b)
+ {
+ int
+ m = pos[a],
+ n = pos[b],
+ k;
+
+ if (m == n)
+ return repr[m];
+
+ if (m > n)
+ {
+ k = n;
+ n = m;
+ m = k;
+ }
+
+ k = (int)(log((double)(n - m))/log(2.));
+
+ if (depth[euler[Minim[k][m]]] < depth[euler[Minim[k][n - (1 << k)]]]) {
+ return repr[euler[Minim[k][m]]];
+ } else {
+ return repr[euler[Minim[k][n - (1 << k)]]];
+ }
+ }
+
+
+ template <class I, class N>
+ void najman_wst<I, N>::removeOneSonNodes(site *p, mln_ch_value(I, site)
&newPar_node)
+ {
+ node &n = nodes(*p);
+
+ if (n.children.size() == 1) // this node has 1 son, delete it
+ {
+ n.area = -1;
+ newPar_node(*p) = n.children[0];
+ *p = n.children[0];
+ removeOneSonNodes(p, newPar_node);
+ }
+ else // there is more than one son, recursive call
+ {
+ for (unsigned i = 0; i < n.children.size(); ++i)
+ removeOneSonNodes(&(n.children[i]), newPar_node);
+ }
+ }
+
+
+ template <class I, class N>
+ void najman_wst<I, N>::compressTree()
+ {
+ mln_ch_value(I, site) newPar_node(Par_node.domain(), Par_node.border());
+
+ // Remove the nodes with one son
+ removeOneSonNodes(&Root, newPar_node);
+
+ // Update the references on deleted nodes
+ mln_piter(I) p(Par_node.domain());
+ for_all(p)
+ while (nodes(Par_node(p)).area == -1)
+ Par_node(p) = newPar_node(Par_node(p));
+ }
+
+# endif // MLN_INCLUDE_ONLY
+
+
+ }; // namespace morpho
+
+}; // namespace mln
+
+#endif // MLN_MORPHO_NAJMAN_WST_HH
Index: abraham/mln/morpho/basic_najman.hh
--- abraham/mln/morpho/basic_najman.hh (revision 2672)
+++ abraham/mln/morpho/basic_najman.hh (working copy)
@@ -33,8 +33,10 @@
#include <mln/estim/min_max.hh>
#include <mln/math/sqr.hh>
#include <mln/util/greater_psite.hh>
+#include <mln/util/ord.hh>
-#include <mln/core/site_set/p_priority_queue.hh>
+#include <mln/core/site_set/p_priority.hh>
+#include <mln/core/site_set/p_queue_fast.hh>
#include <queue>
#include <vector>
@@ -78,8 +80,8 @@
}
}; // struct node
- I pima;
- const Neighborhood<N>& nbh;
+ mln_exact(I) pima;
+ const mln_exact(N)& nbh;
mln_ch_value(I, site) Par_node;
mln_ch_value(I, site) Par_tree;
@@ -148,6 +150,8 @@
compressTree();
build_euler_tour();
build_minim();
+ m_watershed();
+ w_watershed();
}
@@ -269,22 +273,6 @@
return n;
}
- template <typename C>
- void image_output(image2d<C>& img)
- {
- for(unsigned int i = 0; i < img.domain().len(0); ++i)
- {
- for(unsigned int j = 0; j < img.domain().len(1); ++j)
- {
- C l = (img(site(i, j)));//.row() * img.domain().len(1) + (img(site(i, j))).col();
- std::cout << " " << l << " ";
- }
- std::cout << std::endl;
- }
- }
-
-
-
void print_tree(site p)
{
node& n = nodes(p);
@@ -335,12 +323,12 @@
}
- site highest_fork (p_set<site>& components)
+ bool highest_fork (p_set<site>& components, site &r)
{
if (components.nsites() == 0)
{
std::cerr << "highest fork : empty set" << std::endl;
- return site(-1, -1);
+ return false;
}
site
@@ -360,12 +348,13 @@
}
if (nodes(m).level == nodes(hfork).level)
- return site(-1, -1);
+ return false;
- return hfork;
+ r = hfork;
+ return true;
}
- site w_destructible(site p)
+ bool w_destructible(site p, site &r)
{
mln_niter(N) q(nbh, p);
p_set<site> v;
@@ -375,20 +364,45 @@
v.insert(Par_node(q));
if (v.nsites() == 0)
- return site(-1, -1);
+ return false;
- site hf = highest_fork(v);
+ site hf;
+ bool is_hf = highest_fork(v, hf);
- if (hf == site(-1, -1))
- return min(v);
+ if (!is_hf) {
+ r = min(v);
+ return true;
+ }
- if (nodes(hf).level <= pima(p))
- return hf;
+ if (nodes(hf).level <= pima(p)) {
+ r = hf;
+ return true;
+ }
+ return false;
+ }
- return site(-1, -1);
+ bool highest_fork (p_set<site>& components) {
+ site i;
+ return highest_fork(components, i);
}
- site m_destructible(site p)
+ bool m_destructible(site p) {
+ site i;
+ return m_destructible(p, i);
+ }
+
+ bool w_destructible(site p) {
+ site i;
+ return w_destructible(p, i);
+ }
+
+ bool w_constructible(site p) {
+ site i;
+ return w_constructible(p, i);
+ }
+
+
+ bool m_destructible(site p, site &r)
{
mln_niter(N) q(nbh, p);
p_set<site> v = p_set<site>();
@@ -398,21 +412,23 @@
v.insert(Par_node(q));
if (v.nsites() == 0)
- return site(-1, -1);
+ return false;
// if (nodes(min(v)).children.nsites() != 0)
if (nodes(min(v)).children.size() != 0)
- return (site(-1, -1));
+ return false;
- site hf = highest_fork(v);
+ bool is_hf = highest_fork(v);
- if (hf == site(-1, -1))
- return min(v);
+ if (!is_hf) {
+ r = min(v);
+ return true;
+ }
- return site(-1, -1);
+ return false;
}
- site w_constructible(site p)
+ bool w_constructible(site p, site &r)
{
mln_niter(N) q(nbh, p);
p_set<site> v;
@@ -422,10 +438,12 @@
v.insert(Par_node(q));
if (v.nsites() == 0)
- return site(-1, -1);
+ return false;
- if (v.nsites() == 1)
- return v[0];
+ if (v.nsites() == 1) {
+ r = v[0];
+ return true;
+ }
site
c = max(v),
@@ -445,9 +463,10 @@
}
if (nodes(c).level <= pima(p))
- return site(-1, -1);
+ return false;
- return c;
+ r = c;
+ return true;
}
void m_watershed ()
@@ -465,7 +484,7 @@
mln_piter(I) it(pima.domain());
for_all(it)
- if (m_destructible(it.to_site()) != site(-1, -1))
+ if (m_destructible(it.to_site()))
{
l.push(it.to_site());
isproc(it.to_site()) = true;
@@ -481,9 +500,9 @@
// unmark p
isproc(p) = false;
- i = m_destructible(p);
+ bool is_m = m_destructible(p, i);
- if (i != site(-1, -1))
+ if (is_m)
{
pima(p) = nodes(i).level;
Par_node(p) = i;
@@ -491,7 +510,7 @@
mln_niter(N) q(nbh, p);
for_all(q)
if (pima.has(q) && !isproc(q))
- if (m_destructible(q) != site(-1, -1))
+ if (m_destructible(q))
{
// Add priority queue
l.push(q);
@@ -505,15 +524,9 @@
void w_watershed()
{
- std::map< mln_value(I), std::set<site> > L;
+ mln::p_priority< mln_value(I), p_queue_fast<site> > L;
- // Setting the min and the max of the image
- mln_value(I) k, kmin, kmax;
- mln::estim::min_max(pima, kmin, kmax);
-
- // For k From kmin to kmax - 1 Do Lk <- <empty set>
- for (k = kmin; k < kmax; k++)
- L[k] = std::set<site>();
+ mln_value(I) max = level::compute(accu::meta::max(), pima);
// I K(pima.domain(), pima.border());
mln_ch_value(I, unsigned) K(pima.domain(), pima.border());
@@ -524,27 +537,26 @@
for_all(it)
{
site p = it.to_site();
+ site i;
// i <- W-Destructible(p)
- site i = w_destructible(p);
+ bool is_w = w_destructible(p, i);
// If i != infinity then
- if (i != site(-1, -1))
+ if (is_w)
{
- L[nodes(i).level].insert(p);
+ L.push(max - nodes(i).level, p);
K(p) = nodes(i).level;
H(p) = i;
}
}
- for (k = kmin; k < kmax; k++)
+ while (!L.is_empty())
{
- std::set<site>& Lk = L[k];
+ mln_value(I) k = max - L.highest_priority();
- while (!Lk.empty())
- {
- site p = *(Lk.begin());
- Lk.erase(p);
+ site p = L.front();
+ L.pop();
if (K(p) == k)
{
@@ -555,13 +567,14 @@
for_all(q)
if (pima.has(q) && k < pima(q))
{
- site i = w_destructible(q);
- if (i == site(-1, -1))
+ site i;
+ bool is_w = w_destructible(q, i);
+ if (!is_w)
K(q) = 10000; // FIXME : supposed to be infinity...
else
if (K(q) != nodes(i).level)
{
- L[nodes(i).level].insert(q);
+ L.push(max - nodes(i).level, q);
K(q) = nodes(i).level;
H(q) = i;
}
@@ -569,8 +582,6 @@
}
}
}
- }
-
void revert_tree (site p)
{
@@ -628,7 +639,7 @@
// Mark enqueued sites
mln_ch_value(I, bool) enqueued(pima.domain(), pima.border());
- p_priority_queue < site, mln_value(I) > l;
+ p_priority< mln_value(I), p_queue_fast<site> > l;
// p_queue < site > m;
std::queue<site> m;
mln_value(I) max = mln_max(mln_value(I));
@@ -655,7 +666,7 @@
if (cmax.has(q) && !cmax(q) && !enqueued(q))
{
enqueued(q) = true;
- l.push(q, pima(q));
+ l.push(pima(q), q);
}
m.pop();
}
@@ -670,9 +681,10 @@
l.pop();
enqueued(x) = false;
- site c = w_constructible(x);
+ site c;
+ bool is_w = w_constructible(x, c);
- if (c != site(-1, -1))
+ if (is_w)
{
pima(x) = nodes(c).level;
Par_node(x) = c;
@@ -692,7 +704,7 @@
{
enqueued(q) = true;
- l.push(q, pima(q));
+ l.push(pima(q), q);
}
} // if (c != site(-1, -1))
} // while(!l.empty())
@@ -707,7 +719,7 @@
int *euler;
int *depth;
int ctree_size;
- std::map<site, int> pos;
+ std::map<site, int, mln::util::ord<site> > pos;
site *repr;
int *minim;
Index: abraham/mln/morpho/topo_wst.hh
--- abraham/mln/morpho/topo_wst.hh (revision 0)
+++ abraham/mln/morpho/topo_wst.hh (revision 0)
@@ -0,0 +1,744 @@
+// 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.
+
+#ifndef MLN_MORPHO_TOPO_WST_HH
+# define MLN_MORPHO_TOPO_WST_HH
+
+
+
+#include <mln/level/sort_psites.hh>
+#include <mln/level/fill.hh>
+#include <mln/core/image/image2d.hh>
+#include <mln/core/site_set/p_set.hh>
+#include <mln/estim/min_max.hh>
+#include <mln/math/sqr.hh>
+#include <mln/util/greater_psite.hh>
+#include <mln/util/ord.hh>
+#include <mln/arith/revert.hh>
+
+#include <mln/core/site_set/p_priority.hh>
+#include <mln/core/site_set/p_queue_fast.hh>
+
+#include <queue>
+#include <vector>
+#include <map>
+
+namespace mln
+{
+ namespace morpho
+ {
+
+ template <class I, class N>
+ struct topo_wst
+ {
+
+ // Typedefs
+ // --------
+ typedef mln_site(I) site;
+
+ // Component tree management
+ // -------------------------
+ private:
+ struct node {
+ mln_value(I) level;
+ int area;
+ int highest;
+ typedef mln_site(I) site;
+ // Can't modify the sites in a p_array
+ // p_array<mln_site(I)> children;
+ std::vector<site> children;
+
+ void addChildren(const node& n)
+ {
+ // typename p_array<mln_site(I)>::fwd_piter it(n.children);
+ // for (it.start();
+ // it.is_valid();
+ // it.next())
+ // children.append(it.to_site());
+ for (unsigned i=0; i < n.children.size(); ++i)
+ children.push_back(n.children[i]);
+ }
+
+ void addChild(const site p)
+ {
+ // children.append(n);
+ children.push_back(p);
+ }
+
+ mln_value(I) min() { return mln_min(mln_value(I)); }
+ mln_value(I) max() { return mln_max(mln_value(I)); }
+
+ }; // struct node
+
+ const mln_exact(N)& nbh;
+ mln_ch_value(I, site) Par_node;
+ mln_ch_value(I, site) Par_tree;
+
+ mln_ch_value(I, int) Rnk_tree;
+ mln_ch_value(I, int) Rnk_node;
+ mln_ch_value(I, site) subtreeRoot;
+ mln_ch_value(I, node) nodes;
+ site Root;
+
+ void MakeSet_tree(site x);
+ void MakeSet_node(site x);
+ site Find_tree(site x);
+ void swap(site& x, site& y);
+ site Find_node(site x);
+ void BuildComponentTree();
+ site MergeNode(site& node1, site& node2);
+ site Link_tree(site x, site y);
+ site Link_node(site x, site y);
+ node MakeNode(int level);
+
+
+ // Watershed algorithm
+ // -------------------
+
+ private:
+ mln_ch_value(I, bool) isproc;
+
+ // Ctor
+ public:
+ topo_wst(const Image<I>& i,
+ const Neighborhood<N>& n);
+
+
+ // Result
+ public:
+ mln_exact(I) pima;
+
+ public:
+ void go();
+
+ private:
+ site min (p_set<site>& components);
+ site max (p_set<site>& components);
+ bool highest_fork (p_set<site>& components);
+ bool highest_fork (p_set<site>& components, site &r);
+ bool w_constructible(site p, site &r);
+ bool w_constructible(site p);
+ void topo_watershed();
+
+
+ // Optimized LCA Algorithm
+ private:
+ int *euler;
+ int *depth;
+ int ctree_size;
+ std::map<site, int, mln::util::ord<site> > pos;
+ site *repr;
+
+ int *minim;
+ int **Minim;
+
+ void compute_ctree_size (site p);
+ void build_euler_tour_rec(site p, int &position, int d);
+ void build_euler_tour();
+ void build_minim ();
+ site lca (site a, site b);
+ void removeOneSonNodes(site *p, mln_ch_value(I, site) &newPar_node);
+ void compressTree();
+ }; // struct topo_wst
+
+# ifndef MLN_INCLUDE_ONLY
+
+ // Ctor
+ template <class I, class N>
+ topo_wst<I, N>::topo_wst(const Image<I>& i,
+ const Neighborhood<N>& n)
+ : pima(exact(i)),
+ nbh(exact(n)),
+ Par_node(exact(i).domain(), exact(i).border()),
+ Par_tree(exact(i).domain(), exact(i).border()),
+ Rnk_tree(exact(i).domain(), exact(i).border()),
+ Rnk_node(exact(i).domain(), exact(i).border()),
+ subtreeRoot(exact(i).domain(), exact(i).border()),
+ nodes(exact(i).domain(), exact(i).border()),
+ isproc(exact(i).domain(), exact(i).border())
+ {
+ }
+
+ template <class I, class N>
+ void topo_wst<I, N>::MakeSet_tree(site x)
+ {
+ Par_tree(x) = x;
+ Rnk_tree(x) = 0;
+ }
+
+ template <class I, class N>
+ void topo_wst<I, N>::MakeSet_node(site x)
+ {
+ Par_node(x) = x;
+ Rnk_node(x) = 0;
+ }
+
+ template <class I, class N>
+ typename topo_wst<I, N>::site topo_wst<I, N>::Find_tree(site x)
+ {
+ if (Par_tree(x) != x)
+ Par_tree(x) = Find_tree(Par_tree(x));
+ return Par_tree(x);
+ }
+
+ template <class I, class N>
+ void topo_wst<I, N>::swap(site& x, site& y)
+ {
+ site memo = x;
+ x = y;
+ y = memo;
+ }
+
+ template <class I, class N>
+ typename topo_wst<I, N>::site topo_wst<I, N>::Find_node(site x)
+ {
+ if (Par_node(x) != x)
+ Par_node(x) = Find_node(Par_node(x));
+ return Par_node(x);
+ }
+
+ template <class I, class N>
+ void topo_wst<I, N>::go()
+ {
+ BuildComponentTree();
+ compressTree();
+ arith::revert_inplace(pima);
+ arith::revert_inplace(nodes);
+ build_euler_tour();
+ build_minim();
+ topo_watershed();
+ arith::revert_inplace(pima);
+ }
+
+ template <class I, class N>
+ void topo_wst<I, N>::BuildComponentTree()
+ {
+ // Init:
+
+ // Sort the sites in increasing order
+ p_array<mln_site(I)> S;
+ S = level::sort_psites_increasing(pima);
+
+ // Clear the marker map
+ level::fill(isproc, false);
+ for (int ip = 0; ip < int(S.nsites()); ++ip)
+ {
+ site p = S[ip];
+ MakeSet_node(p);
+ MakeSet_tree(p);
+ // if (subtreeRoot.hold(p))
+ subtreeRoot(p) = p;
+ // if (nodes.hold(p))
+ nodes(p) = MakeNode(pima(p));
+ }
+
+
+
+ typename p_array<mln_site(I)>::fwd_piter ip(S);
+ for_all(ip)
+ {
+ site p = ip.to_site();
+
+ site curCanonicalElt = Find_tree(p);
+ site curNode = Find_node(subtreeRoot(curCanonicalElt));
+
+ mln_niter(N) q(nbh, ip);
+ for_all(q)
+ if (pima.has(q) and isproc(q) and pima(q) <= pima(p))
+ {
+ site adjCanonicalElt = Find_tree(q);
+ site adjNode = Find_node(subtreeRoot(adjCanonicalElt));
+ if (curNode != adjNode)
+ {
+ if (nodes(curNode).level == nodes(adjNode).level)
+ curNode = MergeNode(adjNode, curNode);
+ else
+ {
+ nodes(curNode).addChild(adjNode);
+ nodes(curNode).area += nodes(adjNode).area;
+ nodes(curNode).highest += nodes(adjNode).highest;
+ }
+ }
+
+ curCanonicalElt = Link_tree(adjCanonicalElt, curCanonicalElt);
+ subtreeRoot(curCanonicalElt) = curNode;
+ }
+ isproc(p) = true;
+ }
+ // Pour garder une map de correspondance site <-> local_root
+ // for (int ip = 0; ip < int(S.size()); ++ip)
+ // {
+ // site p = S[ip];
+ // M(p) = Find_node(p);
+ // }
+
+ mln_piter(I) r(Par_node.domain());
+ for_all(r)
+ Par_node(r) = Find_node(r);
+
+ Root = subtreeRoot(Find_tree(Find_node(site(0, 0))));
+ }
+
+
+ template <class I, class N>
+ typename topo_wst<I, N>::site topo_wst<I, N>::MergeNode(site& node1,
site& node2)
+ {
+ site tmpNode = Link_node(node1, node2);
+ site tmpNode2;
+ if (tmpNode == node2)
+ {
+ nodes(node2).addChildren(nodes(node1));
+ tmpNode2 = node1;
+ }
+ else
+ {
+ nodes(node1).addChildren(nodes(node2));
+ tmpNode2 = node2;
+ }
+ nodes(tmpNode).area += nodes(tmpNode2).area;
+ nodes(tmpNode).highest += nodes(tmpNode2).highest;
+ return tmpNode;
+ }
+
+ template <class I, class N>
+ typename topo_wst<I, N>::site topo_wst<I, N>::Link_tree(site x, site y)
+ {
+ if (Rnk_tree(x) > Rnk_tree(y))
+ swap(x, y);
+ else
+ if (Rnk_tree(x) == Rnk_tree(y))
+ Rnk_tree(y) += 1;
+ Par_tree(x) = y;
+ return y;
+ }
+
+ template <class I, class N>
+ typename topo_wst<I, N>::site topo_wst<I, N>::Link_node(site x, site y)
+ {
+ if (Rnk_node(x) > Rnk_node(y))
+ swap(x, y);
+ else
+ if (Rnk_node(x) == Rnk_node(y))
+ Rnk_node(y) += 1;
+ Par_node(x) = y;
+ return y;
+ }
+
+ template <class I, class N>
+ typename topo_wst<I, N>::node topo_wst<I, N>::MakeNode(int level)
+ {
+ node n;
+ n.level = level;
+ n.area = 1;
+ n.highest = level;
+ return n;
+ }
+
+ template <class I, class N>
+ typename topo_wst<I, N>::site topo_wst<I, N>::min (p_set<site>&
components)
+ {
+ if (components.nsites() == 0)
+ return site(-1, -1);
+
+ typename p_set<site>::fwd_piter it(components);
+ site min = components[0];
+
+ for_all(it)
+ if (pima(it.to_site()) < pima(min))
+ min = it.to_site();
+
+ return min;
+ }
+
+ template <class I, class N>
+ typename topo_wst<I, N>::site topo_wst<I, N>::max (p_set<site>&
components)
+ {
+ if (components.nsites() == 0)
+ return site(-1, -1);
+
+ typename p_set<site>::fwd_piter it(components);
+ site max = components[0];
+
+ for_all(it)
+ if (pima(it.to_site()) > pima(max))
+ max = it.to_site();
+
+ return max;
+ }
+
+
+ template <class I, class N>
+ bool topo_wst<I, N>::highest_fork (p_set<site>& components, site
&r)
+ {
+ if (components.nsites() == 0)
+ {
+ std::cerr << "highest fork : empty set" << std::endl;
+ return false;
+ }
+
+ site
+ m = min(components),
+ hfork = m;
+
+ typename p_set<site>::fwd_piter it(components);
+ for_all(it)
+ {
+ // Can't remove the site from the set
+ if (it.to_site() == m)
+ continue;
+
+ site c = lca(hfork, it.to_site());
+ if (c != it.to_site())
+ hfork = c;
+ }
+
+ if (nodes(m).level == nodes(hfork).level)
+ return false;
+
+ r = hfork;
+ return true;
+ }
+
+ template <class I, class N>
+ bool topo_wst<I, N>::highest_fork (p_set<site>& components) {
+ site i;
+ return highest_fork(components, i);
+ }
+
+ template <class I, class N>
+ bool topo_wst<I, N>::w_constructible(site p) {
+ site i;
+ return w_constructible(p, i);
+ }
+
+ template <class I, class N>
+ bool topo_wst<I, N>::w_constructible(site p, site &r)
+ {
+ mln_niter(N) q(nbh, p);
+ p_set<site> v;
+
+ for_all(q)
+ if (pima.has(q) && pima(q) > pima(p))
+ v.insert(Par_node(q));
+
+ if (v.nsites() == 0)
+ return false;
+
+ if (v.nsites() == 1) {
+ r = v[0];
+ return true;
+ }
+
+ site
+ c = max(v),
+ cmax = c;
+
+ typename p_set<site>::fwd_piter it(v);
+ for_all(it)
+ {
+ // Can't remove the site from the set
+ if (it.to_site() == cmax)
+ continue;
+
+ site c1 = lca(c, it.to_site());
+
+ if (c1 != it.to_site())
+ c = c1;
+ }
+
+ if (nodes(c).level <= pima(p))
+ return false;
+
+ r = c;
+ return true;
+ }
+
+
+ template <class I, class N>
+ void topo_wst<I, N>::topo_watershed()
+ {
+ // Maxima components
+ mln_ch_value(I, bool) cmax(pima.domain(), pima.border());
+
+ // Mark enqueued sites
+ mln_ch_value(I, bool) enqueued(pima.domain(), pima.border());
+
+ p_priority< mln_value(I), p_queue_fast<site> > l;
+ // p_queue < site > m;
+ std::queue<site> m;
+ mln_value(I) max = mln_max(mln_value(I));
+
+ std::cout << "Init" << std::endl;
+
+ // Flag C-maxima
+ level::fill(cmax, false);
+ level::fill(enqueued, false);
+
+ mln_piter(I) it(Par_node.domain());
+ for_all(it)
+ // if (nodes(Par_node(it.to_site())).children.nsites() == 0)
+ if (nodes(Par_node(it.to_site())).children.size() == 0)
+ {
+ cmax(it.to_site()) = true;
+ m.push(it.to_site());
+ }
+
+ while (!m.empty())
+ {
+ mln_niter(N) q(nbh, m.front());
+ for_all(q)
+ if (cmax.has(q) && !cmax(q) && !enqueued(q))
+ {
+ enqueued(q) = true;
+ l.push(pima(q), q);
+ }
+ m.pop();
+ }
+
+
+ std::cout << "end" << std::endl;
+
+ // Main loop
+ while(!l.is_empty())
+ {
+ site x = l.front();
+ l.pop();
+ enqueued(x) = false;
+
+ site c;
+ bool is_w = w_constructible(x, c);
+
+ if (is_w)
+ {
+ pima(x) = nodes(c).level;
+ Par_node(x) = c;
+
+ // if (nodes(c).children.nsites() == 0)
+ if (nodes(c).children.size() == 0)
+ cmax(x) = true;
+ else
+ // if (nodes(c).children.nsites() > 1)
+ if (nodes(c).children.size() == 1)
+ std::cerr << "ERREUR COMPOSANTE BRANCHE " <<
nodes(c).children.size() << std::endl;
+
+
+ mln_niter(N) q(nbh, x);
+ for_all(q)
+ if (pima.has(q) && !cmax(q) && !enqueued(q))
+ {
+ enqueued(q) = true;
+
+ l.push(pima(q), q);
+ }
+ } // if (c != site(-1, -1))
+ } // while(!l.empty())
+
+ for_all(it)
+ pima(it.to_site()) = nodes(Par_node(it.to_site())).level;
+ }
+
+
+
+ template <class I, class N>
+ void topo_wst<I, N>::compute_ctree_size (site p)
+ {
+ ctree_size += 1;
+ node& n = nodes(p);
+
+ // typename p_array<mln_site(I)>::fwd_piter it(n.children);
+ // for_all(it)
+ // compute_ctree_size(it.to_site());
+
+ for (unsigned i=0; i < n.children.size(); ++i)
+ compute_ctree_size(n.children[i]);
+ }
+
+
+ template <class I, class N>
+ void topo_wst<I, N>::build_euler_tour_rec(site p, int &position, int d)
+ {
+ if (pos.find(p) == pos.end())
+ pos[p] = position;
+
+ repr[position] = p;
+ depth[position] = d;
+ euler[position] = pos[p];
+ ++position;
+ node& n = nodes(p);
+
+ // typename p_array<mln_site(I)>::fwd_piter it(n.children);
+ // for_all(it)
+ // {
+ // build_euler_tour_rec(it.to_site(), position, d+1);
+ // depth[position] = d; // Not optimized
+ // euler[position] = pos[p];
+ // repr[position] = p; // Pas necessaire?
+ // ++position;
+ // }
+
+ for(unsigned i=0; i < n.children.size(); ++i)
+ {
+ build_euler_tour_rec(n.children[i], position, d+1);
+ depth[position] = d; // Not optimized
+ euler[position] = pos[p];
+ repr[position] = p; // Pas necessaire?
+ ++position;
+ }
+ }
+
+
+ template <class I, class N>
+ void topo_wst<I, N>::build_euler_tour ()
+ {
+ ctree_size = 0;
+ compute_ctree_size(Root);
+
+ int size = 2 * ctree_size - 1;
+
+ // FIXME : free this
+ euler = new int[size];
+ depth = new int[size];
+ repr = new site[size];
+
+ int position = 0;
+ build_euler_tour_rec(Root, position, 0);
+ }
+
+
+ template <class I, class N>
+ void topo_wst<I, N>::build_minim ()
+ {
+ // compute_tree_size(); // Already done in build_euler_tour
+ int size = 2 * ctree_size - 1;
+ int logn = (int)(ceil(log((double)(size))/log(2.0)));
+ // minim.reserve(size);
+ minim = new int [logn * size]; // FIXME : Think about freeing this
+ Minim = new int* [logn];
+
+ Minim[0] = minim;
+
+ // Init
+ for (int i = 0; i < size - 1; ++i)
+ if (depth[euler[i]] < depth[euler[i+1]]) {
+ Minim[0][i] = i;
+ } else {
+ Minim[0][i] = i+1;
+ }
+
+ // Minim[0][size - 1] = size - 1;
+
+ int k;
+ for (int j = 1; j < logn; j++) {
+ k = 1 << (j - 1);
+ Minim[j] = &minim[j * size];
+ for (int i = 0; i < size; i++) {
+ if ((i + (k << 1)) >= size) {
+ //Minim[j][i] = size - 1;
+ }
+ else {
+ if (depth[euler[Minim[j - 1][i]]] <= depth[euler[Minim[j - 1][i + k]]])
+ Minim[j][i] = Minim[j - 1][i];
+ else
+ Minim[j][i] = Minim[j - 1][i + k];
+ }
+ }
+ }
+
+ } // void build_minim ()
+
+
+ template <class I, class N>
+ typename topo_wst<I, N>::site topo_wst<I, N>::lca (site a, site b)
+ {
+ int
+ m = pos[a],
+ n = pos[b],
+ k;
+
+ if (m == n)
+ return repr[m];
+
+ if (m > n)
+ {
+ k = n;
+ n = m;
+ m = k;
+ }
+
+ k = (int)(log((double)(n - m))/log(2.));
+
+ if (depth[euler[Minim[k][m]]] < depth[euler[Minim[k][n - (1 << k)]]]) {
+ return repr[euler[Minim[k][m]]];
+ } else {
+ return repr[euler[Minim[k][n - (1 << k)]]];
+ }
+ }
+
+
+ template <class I, class N>
+ void topo_wst<I, N>::removeOneSonNodes(site *p, mln_ch_value(I, site)
&newPar_node)
+ {
+ node &n = nodes(*p);
+
+ if (n.children.size() == 1) // this node has 1 son, delete it
+ {
+ n.area = -1;
+ newPar_node(*p) = n.children[0];
+ *p = n.children[0];
+ removeOneSonNodes(p, newPar_node);
+ }
+ else // there is more than one son, recursive call
+ {
+ for (unsigned i = 0; i < n.children.size(); ++i)
+ removeOneSonNodes(&(n.children[i]), newPar_node);
+ }
+ }
+
+
+ template <class I, class N>
+ void topo_wst<I, N>::compressTree()
+ {
+ mln_ch_value(I, site) newPar_node(Par_node.domain(), Par_node.border());
+
+ // Remove the nodes with one son
+ removeOneSonNodes(&Root, newPar_node);
+
+ // Update the references on deleted nodes
+ mln_piter(I) p(Par_node.domain());
+ for_all(p)
+ while (nodes(Par_node(p)).area == -1)
+ Par_node(p) = newPar_node(Par_node(p));
+ }
+
+# endif // MLN_INCLUDE_ONLY
+
+
+ }; // namespace morpho
+
+}; // namespace mln
+
+#endif // MLN_MORPHO_TOPO_WST_HH
Index: abraham/mln/core/image/thru.hh
--- abraham/mln/core/image/thru.hh (revision 2672)
+++ abraham/mln/core/image/thru.hh (working copy)
@@ -36,7 +36,9 @@
*/
# include <mln/core/internal/image_value_morpher.hh>
+# include <mln/trait/images.hh>
# include <mln/value/set.hh>
+# include <mln/core/image/shell.hh>
namespace mln
{
@@ -50,8 +52,8 @@
template <typename F, typename I>
struct data< thru<F,I> >
{
- data(const I& ima);
- const I& ima_;
+ data( I& ima);
+ I& ima_;
};
} // end of namespace mln::internal
@@ -64,7 +66,9 @@
template <typename F, typename I>
struct image_< thru<F,I> > : default_image_morpher< I, mln_result(F),
thru<F,I> >
{
- typedef trait::image::value_io::read_only value_io;
+ typedef trait::image::value_io::read_write value_io;
+ typedef trait::image::value_access::computed value_access;
+ typedef trait::image::value_storage::disrupted value_storage;
};
} // end of namespace mln::trait
@@ -85,28 +89,22 @@
typedef mln_result(F) rvalue;
/// Return type of read-write access.
- typedef mln_result(F) lvalue;
+ typedef shell<F,I> lvalue;
/// Skeleton.
typedef thru< tag::value_<mln_result(F)>, tag::image_<I> >
skeleton;
/// Constructor.
- thru(const Image<I>& ima);
+ thru(Image<I>& ima);
/// Initialize an empty image.
- void init_(const Image<I>& ima);
+ void init_(Image<I>& ima);
/// Read-only access of pixel value at point site \p p.
mln_result(F) operator()(const mln_psite(I)& p) const;
- /// Mutable access is only OK for reading (not writing).
- // T operator()(const mln_psite(I)& p);
- };
-
- template <typename F, typename I>
- struct thru <F, const I>
- {
- mln_result(F) operator()(const mln_psite(I)& p);
+ /// Mutable access is for reading and writing.
+ lvalue operator()(const mln_psite(I)& p);
};
# ifndef MLN_INCLUDE_ONLY
@@ -119,7 +117,7 @@
template <typename F, typename I>
inline
- data< thru<F,I> >::data(const I& ima)
+ data< thru<F,I> >::data(I& ima)
: ima_(ima)
{
}
@@ -130,22 +128,22 @@
template <typename F, typename I>
inline
- thru<F,I>::thru(const Image<I>& ima)
+ thru<F,I>::thru(Image<I>& ima)
{
mln_precondition(exact(ima).has_data());
this->data_ = new internal::data< thru<F,I> >(exact(ima));
}
- /*
+
template <typename F, typename I>
inline
void
- thru<F,I>::init_(const Image<I>& ima)
+ thru<F,I>::init_(Image<I>& ima)
{
mln_precondition(exact(ima).has_data());
this->data_ = new internal::data<thru<F,I> >(exact(ima));
}
- */
+
template <typename F, typename I>
inline
@@ -158,10 +156,10 @@
template <typename F, typename I>
inline
- mln_result(F)
+ shell<F, I>
thru<F,I>::operator()(const mln_psite(I)& p)
{
- return *(F*)(void*)&( this->data_->ima_(p) );
+ return shell<F, I>( this->data_->ima_, p );
}
# endif // ! MLN_INCLUDE_ONLY
Index: abraham/mln/core/image/shell.hh
--- abraham/mln/core/image/shell.hh (revision 0)
+++ abraham/mln/core/image/shell.hh (revision 0)
@@ -0,0 +1,118 @@
+// Copyright (C) 2007 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_CORE_VALUE_SHELL_HH
+# define MLN_CORE_VALUE_SHELL_HH
+
+# include <mln/core/concept/image.hh>
+
+namespace mln
+{
+
+ // FIXME : traits
+
+ template <typename F, typename I, class C>
+ struct shell_write;
+
+ template <typename F, typename I>
+ struct shell_write<F, I, Function_v2w2v<void> >
+ {
+ mln_value(I) set_ (I &ima, const mln_site(I) &s, const typename F::result
&v);
+ };
+
+ template <typename F, typename I>
+ struct shell_write<F, I, Function_v2w_w2v<void> >
+ {
+ mln_value(I) set_ (I &ima, const mln_site(I) &s, const typename F::result
&v);
+ };
+
+
+ template <typename F, typename I>
+ struct shell;
+
+ template <typename F, typename I>
+ struct shell : shell_write<F, I, typename F::category> // FIXME : inherit of
value_base or sth?
+ {
+ typedef mln_value(I) value;
+
+ // Ctor
+ shell(Image<I> &ima, const mln_site(I) &s);
+
+ // Read
+ operator value ();
+
+ // Write
+ mln_value(I) operator= (typename F::result);
+
+ protected :
+ I &ima;
+ mln_site(I) s;
+ };
+
+# ifndef MLN_INCLUDE_ONLY
+
+ // Ctor
+
+ template <typename F, typename I>
+ shell<F, I>::shell(Image<I> &ima, const mln_site(I) &s)
+ : ima(exact(ima)), s(s)
+ { }
+
+
+ // Read for everyone
+ template <typename F, typename I>
+ shell<F, I>::operator value ()
+ {
+ return F()(ima(s));
+ }
+
+ // Write for everyone
+ template <typename F, typename I>
+ mln_value(I) shell<F, I>::operator= (typename F::result v)
+ {
+ set_(ima, s, v);
+ }
+
+ template <typename F, typename I>
+ mln_value(I) shell_write<F, I, Function_v2w2v<void> >::set_ (I &ima,
const mln_site(I) &s, const typename F::result &v)
+ {
+ ima(s) = F().f_1(v);
+ return ima(s);
+ }
+
+ template <typename F, typename I>
+ mln_value(I) shell_write<F, I, Function_v2w_w2v<void> >::set_ (I &ima,
const mln_site(I) &s, const typename F::result &v)
+ {
+ ima(s) = F().f_1(ima(s), v);
+ return ima(s);
+ }
+
+# endif // MLN_INCLUDE_ONLY
+
+}; // end of namespace mln
+
+#endif // MLN_CORE_VALUE_SHELL_HH
Index: abraham/mln/core/concept/function.hh
--- abraham/mln/core/concept/function.hh (revision 2672)
+++ abraham/mln/core/concept/function.hh (working copy)
@@ -160,7 +160,7 @@
template <typename E>
struct Function_v2w_w2v : public Function<E>
{
- typedef Function_v2w2v<void> category;
+ typedef Function_v2w_w2v<void> category;
/*
result operator() (value);
Index: abraham/mln/core/concept/meta_fun.hh
--- abraham/mln/core/concept/meta_fun.hh (revision 0)
+++ abraham/mln/core/concept/meta_fun.hh (revision 0)
@@ -0,0 +1,83 @@
+// Copyright (C) 2007 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_CORE_CONCEPT_META_FUN_HH
+# define MLN_CORE_CONCEPT_META_FUN_HH
+
+/*! \file mln/core/concept/meta_fun.hh
+ *
+ * \brief Definition of concept of meta function.
+ */
+
+# include <mln/core/concept/object.hh>
+# include <mln/core/concept/function.hh>
+
+namespace mln
+{
+
+ template <class M, class T>
+ struct function;
+
+ namespace meta
+ {
+
+ template <class M>
+ struct impl
+ {
+
+ template <class T>
+ struct info
+ {
+ typedef function<M, T> F;
+ typedef typename F::result result;
+ typedef typename F::lresult lresult;
+ };
+
+ template <class T>
+ typename info<T>::result
+ operator()(const T& t) const
+ {
+ function<M,T> f;
+ return f.read(t);
+ }
+
+ template <class T>
+ T&
+ f_1(typename info<T>::result v, T& t)
+ {
+ function<M,T> f;
+ f.write(t) = v;
+ return t;
+ }
+
+ };
+
+ } // end of namespace meta
+
+} // end of namespace mln
+
+#endif // MLN_CORE_CONCEPT_META_FUN_HH
Index: abraham/mln/fun/meta/red.hh
--- abraham/mln/fun/meta/red.hh (revision 0)
+++ abraham/mln/fun/meta/red.hh (revision 0)
@@ -0,0 +1,45 @@
+// Copyright (C) 2007 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_FUN_META_RED_HH
+# define MLN_FUN_META_RED_HH
+
+# include <mln/core/concept/meta_fun.hh>
+
+namespace mln {
+
+ namespace meta {
+
+ struct red : impl<red> { };
+
+ };
+
+ meta::red red; // fun obj
+
+};
+
+#endif // MLN_FUN_META_RED_HH
Index: abraham/mln/fun/v2w_w2v/norm.hh
--- abraham/mln/fun/v2w_w2v/norm.hh (revision 2672)
+++ abraham/mln/fun/v2w_w2v/norm.hh (working copy)
@@ -71,7 +71,7 @@
* \see mln::norm::l2.
*/
template <typename V, typename R>
- struct l2_norm : public Function_v2v< l2_norm<V, R> >
+ struct l2_norm : public Function_v2w_w2v< l2_norm<V, R> >
{
typedef R result;
R operator()(const V& v) const;
@@ -85,7 +85,7 @@
* \see mln::norm::linfty.
*/
template <typename V, typename R>
- struct linfty_norm : public Function_v2v< linfty_norm<V, R> >
+ struct linfty_norm : public Function_v2w_w2v< linfty_norm<V, R> >
{
typedef R result;
R operator()(const V& v) const;
Index: vigouroux/tests.cc
--- vigouroux/tests.cc (revision 2672)
+++ vigouroux/tests.cc (working copy)
@@ -1,4 +1,4 @@
-#include <mln/core/image_if_value.hh>
+#include <mln/core/image/image_if.hh>
#include <mln/core/alias/w_window2d_int.hh>
#include <mln/display/show.hh>
Index: vigouroux/color/rgb_to_cmy.hh
--- vigouroux/color/rgb_to_cmy.hh (revision 2672)
+++ vigouroux/color/rgb_to_cmy.hh (working copy)
@@ -1,4 +1,4 @@
-#include <mln/core/image_if_value.hh>
+#include <mln/core/image/image_if.hh>
#include <mln/core/alias/w_window2d_int.hh>
#include <mln/display/show.hh>
Index: vigouroux/color/rgb_to_xyz.hh
--- vigouroux/color/rgb_to_xyz.hh (revision 2672)
+++ vigouroux/color/rgb_to_xyz.hh (working copy)
@@ -1,7 +1,7 @@
#ifndef OLENA_CONVERT_RGBXYZ_HH
# define OLENA_CONVERT_RGBXYZ_HH
-# include <mln/core/image_if_value.hh>
+# include <mln/core/image/image_if.hh>
# include <mln/core/alias/w_window2d_int.hh>
# include <mln/display/show.hh>
Index: vigouroux/color/rgb_to_hsi.hh
--- vigouroux/color/rgb_to_hsi.hh (revision 2672)
+++ vigouroux/color/rgb_to_hsi.hh (working copy)
@@ -1,7 +1,7 @@
#include <cmath>
-#include <mln/core/image_if_value.hh>
+#include <mln/core/image/image_if.hh>
#include <mln/core/alias/w_window2d_int.hh>
#include <mln/display/show.hh>
Index: vigouroux/color/rgb_to_yuv.hh
--- vigouroux/color/rgb_to_yuv.hh (revision 2672)
+++ vigouroux/color/rgb_to_yuv.hh (working copy)
@@ -1,4 +1,4 @@
-#include <mln/core/image_if_value.hh>
+#include <mln/core/image/image_if.hh>
#include <mln/core/alias/w_window2d_int.hh>
#include <mln/display/show.hh>