Commit 6790b231 authored by Guillaume Lazzara's avatar Guillaume Lazzara
Browse files

Add fastest versions for skeleton constrained related algorithms.

	* mln/morpho/skeleton_constrained.hh,
	* mln/topo/skeleton/crest.hh: Add fastest versions.

	* mln/topo/skeleton/is_simple_point.hh: Rewrite as functor and add
	fastest versions.

	* tests/topo/skeleton/is_simple_point.cc: Fix test.
parent c287ea3f
2011-05-17 Guillaume Lazzara <z@lrde.epita.fr>
Add fastest versions for skeleton constrained related algorithms.
* mln/morpho/skeleton_constrained.hh,
* mln/topo/skeleton/crest.hh: Add fastest versions.
* mln/topo/skeleton/is_simple_point.hh: Rewrite as functor and add
fastest versions.
* tests/topo/skeleton/is_simple_point.cc: Fix test.
2011-05-05 Guillaume Lazzara <lazzara@fidji.lrde.epita.fr>
 
* mln/io/magick/save.hh: Add support for opacity.
// Copyright (C) 2008, 2009 EPITA Research and Development Laboratory (LRDE)
// Copyright (C) 2008, 2009, 2011 EPITA Research and Development
// Laboratory (LRDE)
//
// This file is part of Olena.
//
......@@ -39,6 +40,7 @@
# include <mln/extension/adjust_duplicate.hh>
# include <mln/data/fill.hh>
# include <mln/util/timer.hh>
namespace mln
{
......@@ -58,71 +60,242 @@ namespace mln
# ifndef MLN_INCLUDE_ONLY
template <typename I,
typename N, typename F,
typename K, typename R>
inline
mln_ch_value(I, bool)
skeleton_constrained(const Image<I>& input_,
const Neighborhood<N>& nbh_, const F& is_simple,
const Image<K>& constraint_, const Image<R>& priority_)
namespace impl
{
trace::entering("morpho::skeleton_constrained");
const I& input = exact(input_);
const N& nbh = exact(nbh_);
const K& constraint = exact(constraint_);
const R& priority = exact(priority_);
namespace generic
{
mln_precondition(input.is_valid());
mln_precondition(nbh.is_valid());
mln_precondition(constraint.is_valid());
mln_precondition(priority.is_valid());
template <typename I,
typename N, typename F,
typename K, typename R>
inline
mln_ch_value(I, bool)
skeleton_constrained(const Image<I>& input_,
const Neighborhood<N>& nbh_, const F& is_simple,
const Image<K>& constraint_, const Image<R>& priority_)
{
trace::entering("morpho::skeleton_constrained");
extension::adjust_duplicate(input, nbh);
const I& input = exact(input_);
const N& nbh = exact(nbh_);
const K& constraint = exact(constraint_);
const R& priority = exact(priority_);
// FIXME: Tests!
mln_precondition(input.is_valid());
mln_precondition(nbh.is_valid());
mln_precondition(constraint.is_valid());
mln_precondition(priority.is_valid());
typedef mln_psite(I) P;
typedef p_queue_fast<P> Q;
p_priority<mln_value(R), Q> q;
typedef mln_value(I) V;
mlc_is(V, bool)::check();
mln_ch_value(I, bool) output;
extension::adjust_duplicate(input, nbh);
// Initialization.
{
initialize(output, input);
data::fill(output, input);
extension::adjust_duplicate(output, nbh);
mln_piter(I) p(input.domain());
for_all(p)
if (input(p) == false &&
is_simple(input, nbh, p)) // p is a simple point of the background.
{
q.push(priority(p), p);
}
}
// FIXME: Tests!
// Propagation.
{
P p;
mln_niter(N) n(nbh, p);
while (! q.is_empty())
typedef mln_psite(I) P;
typedef p_queue_fast<P> Q;
p_priority<mln_value(R), Q> q;
mln_concrete(I) output;
// Initialization.
{
output = duplicate(input);
extension::adjust_duplicate(output, nbh);
mln_piter(I) p(output.domain());
for_all(p)
if (output(p) == false &&
is_simple.check(input, p)) // <-- is_simple.check
// p is a simple point of the background.
{
q.push(priority(p), p);
}
}
// Propagation.
{
p = q.pop_front();
for_all(n)
if (output.has(n) &&
output(n) == true &&
constraint(n) == false &&
is_simple(output, nbh, n))
P p;
mln_niter(N) n(nbh, p);
while (! q.is_empty())
{
p = q.pop_front();
for_all(n)
if (output.has(n) &&
output(n) == true &&
constraint(n) == false &&
is_simple.check(output, n)) // <-- is_simple.check
{
output(n) = false; // Remove n from object.
q.push(priority(n), n);
}
}
}
trace::exiting("morpho::skeleton_constrained");
return output;
}
} // end of namespace mln::morpho::impl::generic
template <typename I,
typename N, typename F,
typename K, typename R>
inline
mln_ch_value(I, bool)
skeleton_constrained_fast(const Image<I>& input_,
const Neighborhood<N>& nbh_,
const F& is_simple,
const Image<K>& constraint_,
const Image<R>& priority_)
{
trace::entering("morpho::skeleton_constrained_fast");
const I& input = exact(input_);
const N& nbh = exact(nbh_);
const K& constraint = exact(constraint_);
const R& priority = exact(priority_);
mln_precondition(input.is_valid());
mln_precondition(nbh.is_valid());
mln_precondition(constraint.is_valid());
mln_precondition(priority.is_valid());
typedef mln_value(I) V;
mlc_is(V, bool)::check();
// Whatever the value of the extension, results of this fast
// version may differ from the generic one. Indeed, we
// cannot check whether a site belongs to the image or not,
// so we try to not take border site in consideration by
// setting their value to false.
extension::adjust_fill(input, nbh, false);
extension::adjust_fill(constraint, nbh, false);
// FIXME: Tests!
typedef p_queue_fast<unsigned> Q;
p_priority<mln_value(R), Q> q;
mln_concrete(I) output;
// Initialization.
{
output = duplicate(input);
extension::adjust_fill(output, nbh, false);
mln_pixter(I) p_out(output);
for_all(p_out)
if (p_out.val() == false &&
is_simple.check__(input, p_out)) // <-- is_simple.check
// p is a simple point of the background.
{
q.push(priority.element(p_out.offset()), p_out);
}
}
// Propagation.
{
mln_pixter(I) p(output);
mln_nixter(I, N) n(p, nbh);
util::array<int> dp = offsets_wrt(input, nbh);
const unsigned n_nbhs = dp.nelements();
while (! q.is_empty())
{
unsigned p = q.pop_front();
for (unsigned i = 0; i < n_nbhs; ++i)
{
unsigned n = p + dp[i];
if (output.element(n) == true &&
constraint.element(n) == false &&
is_simple.check__(output, n)) // <-- is_simple.check
{
output.element(n) = false; // Remove n from object.
q.push(priority.element(n), n);
}
}
}
}
trace::exiting("morpho::skeleton_constrained_fast");
return output;
}
} // end of namespace mln::morpho::impl
namespace internal
{
template <typename I,
typename N, typename F,
typename K, typename R>
inline
mln_ch_value(I, bool)
skeleton_constrained_dispatch(
mln::trait::image::value_access::any,
mln::trait::image::value_storage::any,
const Image<I>& input,
const Neighborhood<N>& nbh,
const F& is_simple,
const Image<K>& constraint,
const Image<R>& priority)
{
return impl::generic::skeleton_constrained(input, nbh,
is_simple,
constraint,
priority);
}
template <typename I,
typename N, typename F,
typename K, typename R>
inline
mln_ch_value(I, bool)
skeleton_constrained_dispatch(
mln::trait::image::value_access::direct,
mln::trait::image::value_storage::one_block,
const Image<I>& input,
const Neighborhood<N>& nbh,
const F& is_simple,
const Image<K>& constraint,
const Image<R>& priority)
{
return impl::skeleton_constrained_fast(input, nbh,
is_simple,
constraint,
priority);
}
} // end of namespace mln::morpho::internal
template <typename I,
typename N, typename F,
typename K, typename R>
inline
mln_ch_value(I, bool)
skeleton_constrained(const Image<I>& input,
const Neighborhood<N>& nbh, const F& is_simple,
const Image<K>& constraint, const Image<R>& priority)
{
trace::entering("morpho::skeleton_constrained");
mln_ch_value(I, bool)
output = internal::skeleton_constrained_dispatch(
mln_trait_image_value_access(I)(),
mln_trait_image_value_storage(I)(),
input, nbh, is_simple, constraint, priority);
trace::exiting("morpho::skeleton_constrained");
return output;
}
......
// Copyright (C) 2009, 2010 EPITA Research and Development Laboratory
// (LRDE)
// Copyright (C) 2009, 2010, 2011 EPITA Research and Development
// Laboratory (LRDE)
//
// This file is part of Olena.
//
......@@ -34,6 +34,8 @@
# include <mln/core/concept/image.hh>
# include <mln/core/concept/neighborhood.hh>
# include <mln/data/fill.hh>
# include <mln/extension/adjust.hh>
# include <mln/border/equalize.hh>
namespace mln
......@@ -103,55 +105,205 @@ namespace mln
# ifndef MLN_INCLUDE_ONLY
// IMPLEMENTATIONS
template <typename I, typename D, typename N>
mln_concrete(I)
crest(const Image<I>& input_, const Image<D>& dist_map_,
const Neighborhood<N>& nbh_, unsigned psi_threshold)
namespace impl
{
trace::entering("topo::skeleton::crest");
const I& input = exact(input_);
const D& dist_map = exact(dist_map_);
const N& nbh = exact(nbh_);
mlc_equal(mln_value(I), bool)::check();
mln_precondition(input.is_valid());
mln_precondition(dist_map.is_valid());
mln_precondition(nbh.is_valid());
namespace generic
{
template <typename I, typename D, typename N>
mln_concrete(I)
crest(const Image<I>& input_, const Image<D>& dist_map_,
const Neighborhood<N>& nbh_, unsigned psi_threshold)
{
trace::entering("topo::skeleton::impl::generic::crest");
const I& input = exact(input_);
const D& dist_map = exact(dist_map_);
const N& nbh = exact(nbh_);
mlc_equal(mln_value(I), bool)::check();
mln_precondition(input.is_valid());
mln_precondition(dist_map.is_valid());
mln_precondition(nbh.is_valid());
mln_concrete(I) is_crest;
initialize(is_crest, input);
data::fill(is_crest, false);
mln_piter(I) p(input.domain());
mln_niter(N) n(nbh, p);
for_all(p)
{
if (!input(p) || dist_map(p) < static_cast<mln_value(D)>(0))
continue;
unsigned nb_eq = 0;
unsigned nb_lt = 0;
for_all(n)
if (input.domain().has(n)
// We want to only consider sites which are part of
// the skeleton. If this test is removed, sometimes
// edge sites are considered as sites with a high PSI
// which is wrong.
&& dist_map(n) > static_cast<mln_value(D)>(0))
{
if (dist_map(n) == dist_map(p))
++nb_eq;
else if (dist_map(n) < dist_map(p))
++nb_lt;
}
if ((nb_lt + nb_eq) >= psi_threshold) // Pixel Superiority index
is_crest(p) = true;
}
trace::exiting("topo::skeleton::impl::generic::crest");
return is_crest;
}
} // end of namespace mln::topo::skeleton::impl::generic
mln_concrete(I) is_crest;
initialize(is_crest, input);
data::fill(is_crest, false);
mln_piter(I) p(input.domain());
mln_niter(N) n(nbh, p);
for_all(p)
template <typename I, typename D, typename N>
mln_concrete(I)
crest_fastest_2d(const Image<I>& input_, const Image<D>& dist_map_,
const Neighborhood<N>& nbh_, unsigned psi_threshold)
{
if (!input(p) || dist_map(p) < static_cast<mln_value(D)>(0))
continue;
unsigned nb_eq = 0;
unsigned nb_lt = 0;
for_all(n)
if (input.domain().has(n)
// We want to only consider sites which are part of
// the skeleton. If this test is removed, sometimes
// edge sites are considered as sites with a high PSI
// which is wrong.
&& dist_map(n) > static_cast<mln_value(D)>(0))
trace::entering("topo::skeleton::impl::crest_fastest_2d");
const I& input = exact(input_);
const D& dist_map = exact(dist_map_);
const N& nbh = exact(nbh_);
mlc_equal(mln_value(I), bool)::check();
mln_precondition(input.is_valid());
mln_precondition(dist_map.is_valid());
mln_precondition(nbh.is_valid());
mln_precondition(input.domain() == dist_map.domain());
extension::adjust(input, nbh);
border::equalize(input, dist_map, input.border());
mln_concrete(I) is_crest;
initialize(is_crest, input);
data::fill(is_crest, false);
mln_pixter(const I) p(input);
util::array<int> dp = mln::offsets_wrt(input, nbh);
unsigned n_nbhs = dp.nelements();
for_all(p)
{
if (!p.val()
|| dist_map.element(p) < static_cast<mln_value(D)>(0))
continue;
unsigned nb_eq = 0;
unsigned nb_lt = 0;
for (unsigned i = 0; i < n_nbhs; ++i)
{
if (dist_map(n) == dist_map(p))
++nb_eq;
else if (dist_map(n) < dist_map(p))
++nb_lt;
unsigned n = p.offset() + dp[i];
if (// We want to only consider sites which are part of
// the skeleton. If this test is removed, sometimes
// edge sites are considered as sites with a high PSI
// which is wrong.
dist_map.element(n) > static_cast<mln_value(D)>(0))
{
if (dist_map.element(n) == dist_map.element(p))
++nb_eq;
else if (dist_map.element(n) < dist_map.element(p))
++nb_lt;
}
if ((nb_lt + nb_eq) >= psi_threshold) // Pixel Superiority index
is_crest.element(p) = true;
}
if ((nb_lt + nb_eq) >= psi_threshold) // Pixel Superiority index
is_crest(p) = true;
}
trace::exiting("topo::skeleton::impl::crest_fastest_2d");
return is_crest;
}
} // end of namespace mln::topo::skeleton::impl
// DISPATCH
namespace internal
{
template <typename I, typename D, typename N>
mln_concrete(I)
crest_dispatch_2d(mln::trait::image::value_storage::any,
mln::trait::image::value_access::any,
mln::trait::image::ext_domain::any,
const Image<I>& input, const Image<D>& dist_map,
const Neighborhood<N>& nbh, unsigned psi_threshold)
{
return skeleton::impl::generic::crest(input, dist_map,
nbh, psi_threshold);
}
template <typename I, typename D, typename N>
mln_concrete(I)
crest_dispatch_2d(mln::trait::image::value_storage::one_block,
mln::trait::image::value_access::direct,
mln::trait::image::ext_domain::some,
const Image<I>& input, const Image<D>& dist_map,
const Neighborhood<N>& nbh, unsigned psi_threshold)
{
return skeleton::impl::crest_fastest_2d(input, dist_map,
nbh, psi_threshold);
}
template <typename I, typename D, typename N>
mln_concrete(I)
crest_dispatch(const Image<I>& input, const Image<D>& dist_map,
const Neighborhood<N>& nbh, unsigned psi_threshold)
{
unsigned dim = mln_site_(I)::dim;
if (dim == 2)
return
crest_dispatch_2d(mln_trait_image_value_storage(I)(),
mln_trait_image_value_access(I)(),
mln_trait_image_ext_domain(I)(),
input, dist_map, nbh, psi_threshold);
return impl::generic::crest(input, dist_map, nbh, psi_threshold);
}
} // end of namespace mln::topo::skeleton::internal
// FACADES
template <typename I, typename D, typename N>
mln_concrete(I)
crest(const Image<I>& input, const Image<D>& dist_map,
const Neighborhood<N>& nbh, unsigned psi_threshold)
{
trace::entering("topo::skeleton::crest");
mlc_equal(mln_value(I), bool)::check();
mln_precondition(exact(input).is_valid());
mln_precondition(exact(dist_map).is_valid());
mln_precondition(exact(nbh).is_valid());
mln_concrete(I)
output = internal::crest_dispatch(input, dist_map,
nbh, psi_threshold);
trace::exiting("topo::skeleton::crest");
return is_crest;
return output;
}
......
// Copyright (C) 2009 EPITA Research and Development Laboratory (LRDE)
// Copyright (C) 2009, 2011 EPITA Research and Development Laboratory
// (LRDE)
//
// This file is part of Olena.
//
......@@ -60,20 +61,35 @@ namespace mln
\endverbatim
*/
template <typename I, typename N>
bool
is_simple_point(const Image<I>& ima,
const Neighborhood<N>& nbh,
const mln_site(I)& p);
template <typename N>
struct is_simple_point
{
is_simple_point(const Neighborhood<N>& nbh);
template <typename I>
bool check(const I& ima, const mln_psite(I)& p) const;
template <typename I>