Commit e9f9b802 authored by Edwin Carlinet's avatar Edwin Carlinet
Browse files

Fix issue with 8-connexity in the merge procedure of the maxtree app.

    *  mln/morpho/maxtree_hqueue_parallel.hpp,
    *  mln/morpho/maxtree_pqueue.hpp,
    *  mln/morpho/maxtree_pqueue_parallel.hpp,
    *  mln/morpho/maxtree_ufind_parallel.hpp,
    *  mln/morpho/maxtree_ufindrank_parallel.hpp,
    *  mln/morpho/merge_tree.hpp: Add overload for c4/c8 computation.
parent 655bfaf4
......@@ -87,7 +87,7 @@ namespace mln
{
mln_precondition(m_has_previous);
merge_tree(m_ima, m_parent, this->m_current_domain, m_cmp);
merge_tree(m_ima, m_nbh, m_parent, this->m_current_domain, m_cmp);
if (joindomain)
{
m_current_domain.join(other.m_current_domain);
......
......@@ -147,91 +147,9 @@ namespace mln
if (!parallel)
*(--Send) = stack.top();
// std::cout << ima.domain().size() << std::endl;
// std::cout << cpt << std::endl;
// std::cout << cpt2 << std::endl;
// std::cout << "D:" << (cpt2-cpt) << std::endl;
// std::cout << cpt3 << std::endl;
// std::cout << cpt4 << std::endl;
if (!parallel)
assert((Send + ima.domain().size()) == Send_);
/* // Flood
while (!pqueue.empty())
{
size_type p = pqueue.top();
size_type repr = stack.top();
if (cmp(ima[p], ima[repr]))
{
// the component represented by repr is over
// we attach repr to its parent
stack.pop();
assert(!stack.empty());
size_type par = stack.top();
while (cmp(ima[p], ima[par]))
{
stack.pop();
parent[repr] = par;
if (!parallel) *(--Send) = repr;
repr = par;
par = stack.top();
}
if (cmp(ima[par], ima[p])) { // ima[par] < ima[p]
stack.push(p);
par = p;
}
parent[repr] = par;
if (!parallel) *(--Send) = repr;
repr = par;
}
assert(ima[repr] == ima[p]);
bool done = true;
mln_foreach(auto k, nbh_delta_indexes)
{
auto q = p + k;
bool processed;
if (use_dejavu)
processed = deja_vu[q];
else if (!parallel)
processed = (parent[q] != UNINITIALIZED);
else
processed = !(first_index <= q and q < last_index) or (parent[q] != UNINITIALIZED);
if (!processed) {
pqueue.push(q);
if (use_dejavu)
deja_vu[q] = true;
else
parent[q] = INQUEUE;
if (cmp(ima[p], ima[q])) {
stack.push(q);
done = false;
break;
}
}
}
if (done) {
parent[p] = repr;
if (!parallel and p != repr) *(--Send) = p;
pqueue.pop();
}
}
if (!parallel)
while (!stack.empty()) {
*(--Send) = stack.top();
stack.pop();
}
if (!parallel)
assert((Send + ima.domain().size()) == Send_);
*/
}
} // end of namespace mln::morpho::internal
......
......@@ -87,7 +87,7 @@ namespace mln
{
mln_precondition(m_has_previous);
merge_tree(m_ima, m_parent, this->m_current_domain, m_cmp);
merge_tree(m_ima, m_nbh, m_parent, this->m_current_domain, m_cmp);
if (joindomain)
{
m_current_domain.join(other.m_current_domain);
......
......@@ -227,7 +227,7 @@ namespace mln
if (m_has_previous)
{
merge_tree(m_ima, m_parent, this->m_current_domain, m_cmp);
merge_tree(m_ima, m_nbh, m_parent, this->m_current_domain, m_cmp);
m_current_domain.join(domain);
// const box2d& d = m_ima.domain();
......@@ -256,7 +256,7 @@ namespace mln
mln_precondition(m_has_previous);
// Merge trees
merge_tree(m_ima, m_parent, this->m_current_domain, m_cmp);
merge_tree(m_ima, m_nbh, m_parent, this->m_current_domain, m_cmp);
m_current_domain.join(other.m_current_domain);
m_nsplit += other.m_nsplit + 1;
}
......
......@@ -176,7 +176,7 @@ namespace mln
{
mln_precondition(m_has_previous);
merge_tree(m_ima, m_parent, this->m_current_domain, m_cmp);
merge_tree(m_ima, m_nbh, m_parent, this->m_current_domain, m_cmp);
if (joindomain)
{
......
......@@ -14,6 +14,7 @@ namespace mln
// tree merging
template <typename V, typename StrictWeakOrdering>
void merge_tree(const image2d<V>& ima,
const c4_t& _nbh,
image2d<typename image2d<V>::size_type>& parent,
box2d domain,
StrictWeakOrdering cmp)
......@@ -21,6 +22,8 @@ namespace mln
typedef typename image2d<V>::size_type size_type;
mln_precondition(!domain.empty());
(void) _nbh;
point2d p_ = domain.pmin;
point2d q_ = domain.pmin;
p_[0] = domain.pmax[0]-1;
......@@ -68,6 +71,148 @@ namespace mln
}
template <typename V, typename StrictWeakOrdering>
void merge_tree(const image2d<V>& ima,
const c8_t& _nbh,
image2d<typename image2d<V>::size_type>& parent,
box2d domain,
StrictWeakOrdering cmp)
{
typedef typename image2d<V>::size_type size_type;
mln_precondition(!domain.empty());
(void) _nbh;
point2d p_ = domain.pmin;
point2d q_ = domain.pmin;
p_[0] = domain.pmax[0]-1;
q_[0] = domain.pmax[0];
size_type p = ima.index_of_point(p_);
size_type q = ima.index_of_point(q_);
int ncols = ima.ncols();
for (int i = 0; i < ncols; ++i, ++p, ++q)
{
for (int j = -1; j <= 1; ++j)
if ( (i+j) >= 0 and (i+j) < ncols)
{
size_type x = internal::zfind_repr(ima, parent, p);
size_type y = internal::zfind_repr(ima, parent, q + j);
if (cmp(ima[x], ima[y]))
std::swap(x, y);
while (x != y)
{
//std::cout << "-- Merge: " << x << " @ " << y << std::endl;
// check that x and y are representative
mln_assertion(x == parent[x] or ima[parent[x]] != ima[x]);
mln_assertion(y == parent[y] or ima[parent[y]] != ima[y]);
mln_assertion(!cmp(ima[x], ima[y])); // ima(y) <= ima(x)
// we want to attach x to y
if (parent[x] == x)
{
parent[x] = y;
break;
}
else
{
size_type z = internal::zfind_parent(ima, parent, x);
if (!cmp(ima[z], ima[y])) // ima(y) <= ima(z)
x = z;
else
{
parent[x] = y;
x = y;
y = z;
}
}
}
}
}
}
// tree merging
template <typename V, typename StrictWeakOrdering>
void merge_tree_3(const image2d<V>& ima,
image2d<typename image2d<V>::size_type>& parent,
image2d<unsigned>& pos,
std::vector<unsigned>& S,
box2d domain,
StrictWeakOrdering cmp)
{
typedef typename image2d<V>::size_type size_type;
mln_precondition(!domain.empty());
point2d p_ = domain.pmin;
point2d q_ = domain.pmin;
p_[0] = domain.pmax[0]-1;
q_[0] = domain.pmax[0];
size_type p = ima.index_of_point(p_);
size_type q = ima.index_of_point(q_);
auto move_insert_at = [&pos, &S] (unsigned x, int i)
{
assert(pos[x] > (unsigned)i);
for (int j = pos[x]-1; j >= i; --j) {
S[j+1] = S[j];
pos[S[j+1]] = j+1;
}
S[i] = x;
pos[x] = i;
};
unsigned ncols = ima.ncols();
for (unsigned i = 0; i < ncols; ++i, ++p, ++q)
{
size_type x = internal::zfind_repr(ima, parent, p);
size_type y = internal::zfind_repr(ima, parent, q);
if (cmp(ima[x], ima[y]))
std::swap(x, y);
while (x != y)
{
mln_assertion(x == parent[x] or ima[parent[x]] != ima[x]);
mln_assertion(y == parent[y] or ima[parent[y]] != ima[y]);
mln_assertion(!cmp(ima[x], ima[y])); // ima(y) <= ima(x)
// we want to attach x to y
if (parent[x] == x)
{
parent[x] = y;
if (pos[x] < pos[y]) move_insert_at(y, pos[x]);
break;
}
else
{
size_type z = internal::zfind_parent(ima, parent, x);
if (!cmp(ima[z], ima[y])) // ima(y) <= ima(z)
{
if (pos[x] < pos[z]) move_insert_at(z, pos[x]);
x = z;
}
else
{
if (pos[x] < pos[y]) move_insert_at(y, pos[x]);
parent[x] = y;
x = y;
y = z;
}
}
}
}
}
template <typename size_type>
bool
check_S(const image2d<size_type>& parent, const size_type* begin, const size_type* end)
......@@ -85,73 +230,6 @@ namespace mln
return true;
}
/*
static
inline
void
addP(image2d<bool>& dejavu, const image2d<std::size_t>& parent, std::size_t p, (std::size_t*)& dst)
{
if (!deja_vu[p])
{
deja_vu[p] = true;
addAncestor(deja_vu, parent, parent[p], dst);
*(dst++) = p;
}
}
// tree merging
// note: dst may alias with src2 but not src1
void merge_S(const image2d<std::size_t>& parent,
image2d<bool>& deja_vu,
const std::size_t* src1, const std::size_t* end1,
const std::size_t* src2, const std::size_t* end2,
std::size_t* dst)
{
mln_precondition(dst != src1);
mln_precondition(src2 == end1 or (src2 == (dst + (end1-src1))));
std::size_t* buffer = dst;
std::size_t n = (end1-src1) + (end2-src2);
if (parent[*src1] == *src1)
*dst = *(src1++);
else
*dst = *(src2++);
mln_assertion(parent[*dst] == *dst);
deja_vu[*(dst++)] = true;
for (;src1 != end1; ++src1)
if (!deja_vu[*src1])
{
std::size_t q = parent[*src1];
while (!deja_vu[q]) {
*ds
}
*(dst++) = *src1;
}
if (deja_vu[parent[*src1]])
*dst = *(src1++);
else if (deja_vu[parent[*src2]])
*dst = *(src2++);
else
{
}
deja_vu[*(dst++)] = true;
}
mln_assertion((src1 == end1) != (src2 == end2));
if (src1 != end1)
std::copy(src1, end1, dst);
else if (src2 != end2 and src2 != dst)
std::copy(src2, end2, dst);
check_S(parent, buffer, buffer+n);
}
*/
}
}
......
Markdown is supported
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