Commit 2a645e50 authored by Baptiste Esteban's avatar Baptiste Esteban
Browse files

Integration of the directional hqueue in the alphatree algorithm

parent a7ad24fc
Pipeline #27855 failed with stages
in 10 minutes and 48 seconds
......@@ -4,13 +4,15 @@
#include <mln/core/concepts/image.hpp>
#include <mln/core/concepts/neighborhood.hpp>
#include <mln/core/functional_ops.hpp>
#include <mln/core/algorithm/for_each.hpp>
#include <mln/morpho/component_tree.hpp>
#include <mln/core/functional_ops.hpp>
#include <mln/morpho/canvas/unionfind.hpp>
#include <mln/morpho/component_tree.hpp>
#include <mln/morpho/private/directional_hqueue.hpp>
#include <type_traits>
#include <algorithm>
#include <type_traits>
#include <range/v3/functional/concepts.hpp>
......@@ -33,26 +35,78 @@ namespace mln::morpho
/******************************************/
namespace internal
{
template <typename I, typename N, typename W>
class alphatree_edges;
template <class P, class W>
struct alphatree_edge_t
template <typename I, typename N, typename W>
requires(std::is_integral_v<W>&& std::is_unsigned_v<W> && sizeof(W) <= 2) class alphatree_edges<I, N, W>
{
P a;
P b;
W w;
static_assert(is_a_v<I, mln::details::Image>);
using P = image_point_t<I>;
public:
alphatree_edges(){};
void push(std::size_t dir, W w, P p) { m_cont.insert(dir, w, p); }
std::tuple<P, P, W> pop() { return m_cont.pop(); }
W top() const { return m_cont.current_level(); }
bool empty() const { return m_cont.empty(); }
void on_finish_insert() {}
private:
details::directional_hqueue<I, N, W> m_cont;
};
template <typename I, typename N, typename W>
class alphatree_edges
{
static_assert(is_a_v<I, mln::details::Image>);
using P = image_point_t<I>;
public:
alphatree_edges()
: m_current(0)
{
}
void push(std::size_t dir, W w, P p) { m_cont.push_back({p, p + offsets[dir], w}); }
std::tuple<P, P, W> pop()
{
assert(m_current < m_cont.size());
const auto e = m_cont[m_current++];
return {std::move(e.p), std::move(e.q), std::move(e.w)};
}
W top() const
{
assert(m_current < m_cont.size());
return m_cont[m_current].w;
}
bool empty() const { return m_cont.size() == m_current; }
void on_finish_insert()
{
std::sort(m_cont.begin(), m_cont.end(), [](const edge_t& a, const edge_t& b) { return a.w < b.w; });
}
private:
static constexpr auto offsets = N::after_offsets();
struct edge_t
{
P p;
P q;
W w;
};
std::vector<edge_t> m_cont;
std::size_t m_current;
};
// Compute a node_id for each flat zone
template <class I, class J>
[[gnu::noinline]]
std::size_t alphatree_create_nodemap(I node_map, J zpar)
[[gnu::noinline]] std::size_t alphatree_create_nodemap(I node_map, J zpar)
{
std::size_t node_count = 0;
mln_foreach(auto px, node_map.pixels())
mln_foreach (auto px, node_map.pixels())
{
auto p = px.point();
auto rp = canvas::impl::zfindroot(zpar, p);
......@@ -66,21 +120,18 @@ namespace mln::morpho
// Compute flat zone of the image using union-find structures
// It returns the number of edges processed
template <class E, class J>
std::size_t alphatree_compute_flatzones(const E* edges, std::size_t n, J zpar)
void alphatree_compute_flatzones(E& edges, J zpar)
{
canvas::impl::union_find_init_par(zpar);
std::size_t i;
for (i = 0; i < n && edges[i].w == 0; ++i)
while (!edges.empty() && edges.top() == 0)
{
auto p = edges[i].a;
auto q = edges[i].b;
auto rp = canvas::impl::zfindroot(zpar, p);
auto rq = canvas::impl::zfindroot(zpar, q);
const auto [p, q, w] = edges.pop();
auto rp = canvas::impl::zfindroot(zpar, p);
auto rq = canvas::impl::zfindroot(zpar, q);
if (rp != rq) // Merge set
zpar(rq) = rp;
}
return i;
}
// Inverse the node array, updating the parent relation
......@@ -104,39 +155,47 @@ namespace mln::morpho
}
template <class I, class N, class F, class E>
void alphatree_compute_edges(I input, N nbh, F distance, std::vector<E>& edges)
template <class I, class N, class F>
auto alphatree_compute_edges(I input, N nbh, F distance)
{
auto dom = input.domain();
static_assert(is_a_v<I, mln::details::Image>);
using V = image_value_t<I>;
auto edges = alphatree_edges<I, N, std::invoke_result_t<F, V, V>>();
auto dom = input.domain();
mln_foreach (auto p, dom)
{
std::size_t i = 0;
for (auto n : nbh.after(p))
{
if (dom.has(n))
edges.push_back({p, n, distance(input(p), input(n))});
edges.push(i, distance(input(p), input(n)), p);
++i;
}
}
edges.on_finish_insert();
return edges;
}
//
//
template <class E, class I, class W>
std::size_t alphatree_compute_hierarchy(const E* edges, std::size_t n, I node_map, //
std::size_t node_count, //
std::vector<int>& par, //
std::size_t alphatree_compute_hierarchy(E& edges, I node_map, //
std::size_t node_count, //
std::vector<int>& par, //
std::vector<W>& levels)
{
static_assert(mln::is_a<I, mln::details::Image>());
std::vector<int> zpar(node_count); // Parent of ufind structure
std::vector<int> links(node_count); // Node index in the alpha tree
std::vector<int> zpar(node_count); // Parent of ufind structure
std::vector<int> links(node_count); // Node index in the alpha tree
std::iota(std::begin(zpar), std::end(zpar), 0);
std::iota(std::begin(links), std::end(links), 0);
for (std::size_t i = 0; i < n; ++i)
while (!edges.empty())
{
auto [p, q, w] = edges[i];
auto [p, q, w] = edges.pop();
auto rp = canvas::impl::zfindroot(zpar.data(), node_map(p));
auto rq = canvas::impl::zfindroot(zpar.data(), node_map(q));
......@@ -189,24 +248,15 @@ namespace mln::morpho
static_assert(::concepts::totally_ordered<W>);
using edge_t = internal::alphatree_edge_t<P, W>;
// 1. Get the list of edges
std::vector<edge_t> edges;
internal::alphatree_compute_edges(std::move(input), std::move(nbh), std::move(distance), edges);
// Sort the edges by increasing weights
std::sort(std::begin(edges), std::end(edges), [](const auto& a, const auto& b) { return a.w < b.w; });
auto edges = internal::alphatree_compute_edges(std::move(input), std::move(nbh), std::move(distance));
std::size_t edges_null_count;
std::size_t flatzones_count;
image_ch_value_t<I, int> node_map = imchvalue<int>(input).set_init_value(-1);
{
image_ch_value_t<I, P> zpar = imchvalue<P>(input);
image_ch_value_t<I, P> zpar = imchvalue<P>(input);
// 2. Compute flat zone of the image
edges_null_count = internal::alphatree_compute_flatzones(edges.data(), edges.size(), zpar);
internal::alphatree_compute_flatzones(edges, zpar);
// 3. Compute a node_id for each flat zone
flatzones_count = internal::alphatree_create_nodemap(node_map, zpar);
......@@ -215,12 +265,11 @@ namespace mln::morpho
std::size_t node_count = flatzones_count;
std::vector<int> par(node_count);
std::vector<W> levels(node_count, 0);
std::vector<W> levels(node_count, 0);
// 4. Compute the hierarchy
{
std::iota(std::begin(par), std::end(par), 0);
node_count = internal::alphatree_compute_hierarchy(edges.data() + edges_null_count, edges.size() - edges_null_count,
node_map, node_count, par, levels);
node_count = internal::alphatree_compute_hierarchy(edges, node_map, node_count, par, levels);
}
// 5. Parent / levels are ordered from leaves to root, we need to reverse
......@@ -236,4 +285,4 @@ namespace mln::morpho
return {std::move(t), std::move(node_map)};
}
}
} // namespace mln::morpho
......@@ -87,14 +87,14 @@ namespace mln::morpho::details
assert(m_size > 0);
const auto p = m_queues[m_current_dir]->pop_front(m_current_level);
const auto q = p + N::after_offsets()[m_current_dir];
const int w = m_current_level;
const int w = m_current_level;
m_size--;
// Update the current level and the current dir if needed
if (m_size > 0 && m_queues[m_current_dir]->empty(m_current_level))
{
int lvl = m_queues[0]->lower_bound(m_current_level);
int dir = 0;
std::size_t dir = 0;
for (std::size_t i = 1; i < m_ndir; i++)
{
int tmp = m_queues[i]->lower_bound(m_current_level);
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment