[Olena-patches] 1818: Add personal folder in sandbox with prototype of najman component tree.

1 Apr 2008

https://svn.lrde.epita.fr/svn/oln/trunk/milena


Index: ChangeLog
from  Alexandre Abraham  <abraham@lrde.epita.fr>

	Add personal folder in sandbox with prototype of najman component tree.

	* sandbox/abraham: New.
	* sandbox/abraham/morpho: New.
	* sandbox/abraham/morpho/basic_najman.hh: New, contains Najman component tree.
	* sandbox/abraham/morpho/test.cc: New.
	* sandbox/abraham/morpho/images: New.
	* sandbox/abraham/morpho/images/test_component_tree.pgm: New.
	* sandbox/abraham/morpho/topo_wst.hh: New, draft of topological watershed.
	* sandbox/abraham/morpho/Makefile: New.

 Makefile        |    8 +
 basic_najman.hh |  287 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 test.cc         |   43 ++++++++
 topo_wst.hh     |  130 +++++++++++++++++++++++++
 4 files changed, 468 insertions(+)

Index: sandbox/abraham/morpho/basic_najman.hh

--- sandbox/abraham/morpho/basic_najman.hh	(revision 0)
+++ sandbox/abraham/morpho/basic_najman.hh	(revision 0)
@@ -0,0 +1,287 @@
+#include <mln/level/sort_psites.hh>
+#include <mln/level/fill.hh>
+#include <mln/core/image2d.hh>
+
+namespace mln
+{
+
+  namespace morpho
+  {
+
+    template <class I, class N>
+    struct basic_najman
+    {
+
+      typedef mln_psite(I) psite;
+
+      struct node {
+	int level;
+	int area;
+	int highest;
+	p_array<mln_psite(I)> children;
+
+	void addChildren(const node& n)
+	{
+	  typename p_array<mln_psite(I)>::fwd_piter it(n.children);
+	  for (it.start();
+	       it.is_valid();
+	       it.next())
+	    children.append(it.to_psite());
+	}
+	void addChild(const psite n)
+	{
+	  children.append(n);
+	}
+      }; // struct node
+
+
+      const Image<I>& ima;
+      const Neighborhood<N>& nbh;
+      mln_ch_value(I, psite) par;
+
+      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) lowestNode;
+      mln_ch_value(I, node) nodes;
+      mln_ch_value(I, bool) isproc;
+      psite Root;
+      p_array<mln_psite(I)> S;
+      // std::map<psite,psite> M; // to keep the root of any point.
+
+      basic_najman(const Image<I>& i,
+		   const Neighborhood<N>& n)
+	: ima(i),
+	  nbh(n),
+	  par(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()),
+	  lowestNode(exact(i).domain(), exact(i).border()),
+	  nodes(exact(i).domain(), exact(i).border()),
+	  isproc(exact(i).domain(), exact(i).border())
+      {
+      }
+
+
+      void MakeSet_tree(psite x)
+      {
+	Par_tree(x) = x;
+	Rnk_tree(x) = 0;
+      }
+
+      void MakeSet_node(psite x)
+      {
+	par(x) = x; // was Par_node(x) = x;
+	Rnk_node(x) = 0;
+      }
+
+      psite Find_tree(psite x)
+      {
+	if (Par_tree(x) != x)
+	  Par_tree(x) = Find_tree(Par_tree(x));
+	return Par_tree(x);
+      }
+
+      void swap(psite& x, psite& y)
+      {
+	psite memo = x;
+	x = y;
+	y = memo;
+      }
+
+      bool is_root(const psite& x) const
+      {
+	return par(x) == x;
+      }
+
+      bool is_level_root(const psite& x) const
+      {
+	return is_root(x) || exact(ima)(par(x)) < exact(ima)(x);
+      }
+
+      psite find_level_root(const psite& x)
+      {
+	if (is_level_root(x))
+	  return x;
+	else
+	  return par(x) = find_level_root(par(x));
+      }
+
+      psite Find_node(psite x)
+      {
+	if (par(x) != x)
+	  par(x) = Find_node(par(x));
+	return par(x);
+      }
+
+      void go()
+      {
+	init();
+	BuildComponentTree();
+      }
+
+
+      void init()
+      {
+	// Sort the points in increasing order
+	S = level::sort_psites_increasing(ima);
+
+	// Clear the marker map
+	level::fill(isproc, false);
+	for (int ip = 0; ip < int(S.npoints()); ++ip)
+	  {
+	    psite p = S[ip];
+	    MakeSet_node(p);
+	    MakeSet_tree(p);
+	    // 	      if (lowestNode.hold(p))
+	    lowestNode(p) = p;
+	    // 	      if (nodes.hold(p))
+	    nodes(p) = MakeNode(exact(ima)(p));
+	  }
+      }
+
+
+      void BuildComponentTree()
+      {
+	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(lowestNode(curTree));
+
+	    mln_niter(N) q(nbh, p);
+
+	    for_all(q)
+	      if (exact(ima).has(q) and isproc(q) and exact(ima)(q) >= exact(ima)(p)) // f[q] >= f[p] parce qu'on prend l'ordre de l'histogramme
+		{
+		  psite adjTree = Find_tree(q);
+		  psite adjNode = Find_node(lowestNode(adjTree));
+		  if (curNode != adjNode)
+		    {
+		      if (nodes(curNode).level == nodes(adjNode).level)
+			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);
+		  lowestNode(curTree) = curNode;
+		}
+	  }
+	Root = lowestNode(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)
+	//              {
+	//                psite p = S[ip];
+	// 	       M(p) = Find_node(p);
+	// 	     }
+
+	mln_piter(I) r(par.domain());
+	for_all(r)
+	  par(r) = Find_node(r);
+      }
+
+
+      psite MergeNode(psite& node1, psite& node2)
+      {
+	psite tmpNode = Link_node(node1, node2);
+	psite 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;
+      }
+
+      psite Link_tree(psite x, psite 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;
+      }
+
+      psite Link_node(psite x, psite y)
+      {
+	if (Rnk_node(x) > Rnk_node(y))
+	  swap(x, y);
+	else
+	  if (Rnk_node(x) == Rnk_node(y))
+	    Rnk_node(y) += 1;
+	par(x) = y; // was Par_node(x) = y;
+	return y;
+      }
+
+      node MakeNode(int level)
+      {
+	node n;
+	n.level = level;
+	n.area = 1;
+	n.highest = level;
+	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(psite(i, j)));//.row() * img.domain().len(1) + (img(psite(i, j))).col();
+		std::cout << " " << l << " ";
+	      }
+	    std::cout << std::endl;
+	  }
+      }
+
+
+
+      void print_tree(psite p)
+      {
+	node& n = nodes(p);
+	std::cout << "psite " << p << "=" << (p.row() * exact(lowestNode).domain().len(1) + p.col()) << " : ";
+
+	typename p_array<mln_psite(I)>::fwd_piter it(n.children);
+	for (it.start();
+	     it.is_valid();
+	     it.next())
+	  {
+	    psite q = it.to_psite();
+	    std::cout << q << "=" << (q.row() * lowestNode.domain().len(2) + q.col()) << " ";
+	  }
+	std::cout << std::endl;
+
+	for (it.start();
+	     it.is_valid();
+	     it.next())
+	  {
+	    print_tree(it.to_psite());
+	  }
+      }
+
+    }; // struct basic_najman
+
+  }; // namespace morpho
+
+}; // namespace mln
Index: sandbox/abraham/morpho/test.cc
--- sandbox/abraham/morpho/test.cc	(revision 0)
+++ sandbox/abraham/morpho/test.cc	(revision 0)
@@ -0,0 +1,43 @@
+#include <mln/core/image2d.hh>
+#include <mln/core/window2d.hh>
+#include <mln/core/neighb2d.hh>
+
+#include <mln/value/int_u8.hh>
+
+//#include <mln/morpho/meyer_wst.hh>
+
+#include <mln/io/pgm/load.hh>
+#include <mln/io/pgm/save.hh>
+
+# include "basic_najman.hh"
+// # include "topo_wst.hh"
+
+int main ()
+{
+  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");
+
+  morpho::basic_najman< image2d<int_u8>, neighb2d> n(input, c8());
+  n.image_output(input);
+  n.go();
+  n.image_output(n.Par_tree);
+  n.print_tree(n.Root);
+  n.image_output(n.par);
+  std::cout << " --------------------------- " << std::endl;
+  n.image_output(n.Rnk_tree);
+  std::cout << " --------------------------- " << std::endl;
+  n.image_output(n.Rnk_node);
+
+  // image2d<int_u8> output = morpho::topo_wst(input, c4());
+  //  io::pgm::save(output, "out.pgm");
+
+  return 0;
+}
Index: sandbox/abraham/morpho/images/test_component_tree.pgm
Cannot display: file marked as a binary type.
svn:mime-type = application/octet-stream

Property changes on: sandbox/abraham/morpho/images/test_component_tree.pgm
___________________________________________________________________
Name: svn:mime-type
   + application/octet-stream

Index: sandbox/abraham/morpho/topo_wst.hh
--- sandbox/abraham/morpho/topo_wst.hh	(revision 0)
+++ sandbox/abraham/morpho/topo_wst.hh	(revision 0)
@@ -0,0 +1,130 @@
+// 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
+
+/** \file mln/morpho/topo_wst.hh
+    \brief Topological Watershed Transform (WST) algorithm.
+
+    Reference:
+      Quasi-linear algorithms for the topological watershed. Michel
+      Courie and Laurent Najman and Gilles Bertrand.
+*/
+
+# include <mln/trait/ch_value.hh>
+
+//# include <mln/util/greater_psite.hh>
+# include <mln/morpho/includes.hh>
+//# include <mln/morpho/extrema_components.hh>
+
+# include "component_tree.hh"
+
+namespace mln
+{
+
+  namespace morpho
+  {
+
+    /** \brief Topological Watershed Transform (WST) algorithm.
+
+	\param[in]  input   The input image.
+	\param[in]  nbh     The connexity of markers.
+	\param[out] nbasins The number of basins.
+
+	\li \p L is the type of labels, used to number the watershed
+	itself (with the minimal value), and the basins.
+	\li \p I is the exact type of the input image.
+	\li \p N is the exact type of the neighborhood used to express
+	\a input's connexity.  */
+    template <typename L, typename I, typename N>
+    mln_ch_value(I, L)
+    topo_wst(const Image<I>& input, const Neighborhood<N>& nbh);
+
+    template <typename I, typename P>
+    bool is_m_destructible(const Image<I>& ima, const Point<P>& point, component_tree c);
+
+    template <typename L, typename I, typename N>
+    mln_ch_value(I, L)
+    m_wst(const Image<I>& input, const Neighborhood<N>& nbh);
+
+
+# ifndef MLN_INCLUDE_ONLY
+
+    template <typename I, typename P>
+    int is_m_destructible(const Image<I>& ima, const mln_psite(I) p, const Neighborhood<N>& nbh, component_tree c)
+    {
+      p_set<P> v;
+
+      mln_niter(N) n(nbh, p);
+      for_all(n)
+	if (n < p)
+	  v = uni(V, c(n));
+      if (v.is_empty())
+	return -1;
+
+      p_set<P>::fwd_piter p = v.begin();
+      if (c(*p).children.size() != 0)
+	return -1;
+
+      component_tree highest_fork = c.highest_fork(v);
+      if (highest_fork.level == -1)
+	return -1;
+
+      return highest_fork.level;
+    }
+
+    template <typename L, typename I, typename N>
+    mln_ch_value(I, L)
+    m_wst(const Image<I>& input, const Neighborhood<N>& nbh);
+    {
+      typedef
+	std::priority_queue< psite, std::vector<psite>, util::greater_psite<I> >
+	ordered_queue_type;
+
+      ordered_queue_type queue();
+      component_tree c(input);
+
+      mln_piter(I) p(input.domain());
+
+      for_all(p)
+      {
+	if (is_m_destructible(input, p, c4(), c
+      }
+
+    }
+
+
+
+# endif // ! MLN_INCLUDE_ONLY
+
+  } // end of namespace mln::morpho
+
+} // end of namespace mln
+
+
+#endif // ! MLN_MORPHO_TOPOLOGICAL_WST_HH
Index: sandbox/abraham/morpho/Makefile
--- sandbox/abraham/morpho/Makefile	(revision 0)
+++ sandbox/abraham/morpho/Makefile	(revision 0)
@@ -0,0 +1,8 @@
+SOURCE=test.cc
+OBJECTS=test.o
+CXXFLAGS=-Wall -W -I ../../.. -g
+
+all:${OBJECTS}
+	g++ ${FLAGS} ${OBJECTS} -o test
+
+test.o: basic_najman.hh
\ No newline at end of file

    