- Olena-patches

last-svn-commit-188-gf36081f Actually use the maximal curvature in mesh-complex-max-curv-segm.
by Roland Levillain 14 Jun '10

14 Jun '10

* apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc: Here. --- milena/ChangeLog | 6 ++++++ .../mesh-segm-skel/mesh-complex-max-curv-segm.cc | 11 +++++++---- 2 files changed, 13 insertions(+), 4 deletions(-) diff --git a/milena/ChangeLog b/milena/ChangeLog index 00afa02..923f16b 100644 --- a/milena/ChangeLog +++ b/milena/ChangeLog @@ -1,5 +1,11 @@ 2010-06-14 Roland Levillain <roland(a)lrde.epita.fr> + Actually use the maximal curvature in mesh-complex-max-curv-segm. + + * apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc: Here. + +2010-06-14 Roland Levillain <roland(a)lrde.epita.fr> + Improve the precision of some mesh-segm-skel apps. * apps/mesh-segm-skel/mesh-complex-segm.cc, diff --git a/milena/apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc b/milena/apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc index 1f43381..a7a741a 100644 --- a/milena/apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc +++ b/milena/apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc @@ -91,14 +91,17 @@ int main(int argc, char* argv[]) mln::p_n_faces_fwd_piter<D, G> v(float_ima.domain(), 0); for_all(v) { + // FIXME: The program should allow the user to choose the kind + // of measure (curvature). +#if 0 float h = (curv.first(v) + curv.second(v)) / 2; // Pseudo-inverse curvature. float h_inv = 1 / pi * (atan(-h) + pi / 2); float_ima(v) = h_inv; - // FIXME: The program should allow the user to choose the kind - // of measure. -// input(v) = mln::math::max(mln::math::sqr(curv.first(v)), -// mln::math::sqr(curv.second(v))); +#endif + // Max curvature. + float_ima(v) = mln::math::max(mln::math::sqr(curv.first(v)), + mln::math::sqr(curv.second(v))); } // Values on edges. -- 1.5.6.5

1 0

last-svn-commit-188-gc68caaf Actually use the maximal curvature in mesh-complex-max-curv-segm.
by Roland Levillain 14 Jun '10

14 Jun '10

* apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc: Here. --- milena/ChangeLog | 6 ++++++ .../mesh-segm-skel/mesh-complex-max-curv-segm.cc | 11 +++++++---- 2 files changed, 13 insertions(+), 4 deletions(-) diff --git a/milena/ChangeLog b/milena/ChangeLog index 8d0f1be..1fe8399 100644 --- a/milena/ChangeLog +++ b/milena/ChangeLog @@ -1,5 +1,11 @@ 2010-06-14 Roland Levillain <roland(a)lrde.epita.fr> + Actually use the maximal curvature in mesh-complex-max-curv-segm. + + * apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc: Here. + +2010-06-14 Roland Levillain <roland(a)lrde.epita.fr> + Improve the precision of some mesh-segm-skel apps. * apps/mesh-segm-skel/mesh-complex-segm.cc, diff --git a/milena/apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc b/milena/apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc index 1f43381..a7a741a 100644 --- a/milena/apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc +++ b/milena/apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc @@ -91,14 +91,17 @@ int main(int argc, char* argv[]) mln::p_n_faces_fwd_piter<D, G> v(float_ima.domain(), 0); for_all(v) { + // FIXME: The program should allow the user to choose the kind + // of measure (curvature). +#if 0 float h = (curv.first(v) + curv.second(v)) / 2; // Pseudo-inverse curvature. float h_inv = 1 / pi * (atan(-h) + pi / 2); float_ima(v) = h_inv; - // FIXME: The program should allow the user to choose the kind - // of measure. -// input(v) = mln::math::max(mln::math::sqr(curv.first(v)), -// mln::math::sqr(curv.second(v))); +#endif + // Max curvature. + float_ima(v) = mln::math::max(mln::math::sqr(curv.first(v)), + mln::math::sqr(curv.second(v))); } // Values on edges. -- 1.5.6.5

1 0

last-svn-commit-183-gc9b58ad Improve the precision of some mesh-segm-skel apps.
by Roland Levillain 14 Jun '10

14 Jun '10

* apps/mesh-segm-skel/mesh-complex-segm.cc, * apps/mesh-segm-skel/mesh-complex-skel.cc, * apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc, * apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc: Convert the floating point value input to an unsigned integer image to avoid precision errors due to float comparisons. --- milena/ChangeLog | 11 ++++ .../mesh-segm-skel/mesh-complex-max-curv-segm.cc | 39 ++++++++---- .../mesh-segm-skel/mesh-complex-max-curv-skel.cc | 61 +++++++++++++------- milena/apps/mesh-segm-skel/mesh-complex-segm.cc | 24 ++++++-- milena/apps/mesh-segm-skel/mesh-complex-skel.cc | 41 ++++++++++--- 5 files changed, 127 insertions(+), 49 deletions(-) diff --git a/milena/ChangeLog b/milena/ChangeLog index 3de4207..2ce43e0 100644 --- a/milena/ChangeLog +++ b/milena/ChangeLog @@ -1,5 +1,16 @@ 2010-06-14 Roland Levillain <roland(a)lrde.epita.fr> + Improve the precision of some mesh-segm-skel apps. + + * apps/mesh-segm-skel/mesh-complex-segm.cc, + * apps/mesh-segm-skel/mesh-complex-skel.cc, + * apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc, + * apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc: + Convert the floating point value input to an unsigned integer + image to avoid precision errors due to float comparisons. + +2010-06-14 Roland Levillain <roland(a)lrde.epita.fr> + New application: mesh-complex-max-curv-skel. * apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc: New. diff --git a/milena/apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc b/milena/apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc index 477ccdf..1f43381 100644 --- a/milena/apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc +++ b/milena/apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc @@ -75,25 +75,26 @@ int main(int argc, char* argv[]) `----------------*/ // Image type. - typedef mln::float_2complex_image3df ima_t; + typedef mln::float_2complex_image3df float_ima_t; // Dimension of the image (and therefore of the complex). - static const unsigned D = ima_t::dim; + static const unsigned D = float_ima_t::dim; // Geometry of the image. - typedef mln_geom_(ima_t) G; + typedef mln_geom_(float_ima_t) G; mln::bin_2complex_image3df bin_input; mln::io::off::load(bin_input, input_filename); - std::pair<ima_t, ima_t> curv = mln::geom::mesh_curvature(bin_input.domain()); + std::pair<float_ima_t, float_ima_t> curv = + mln::geom::mesh_curvature(bin_input.domain()); // Compute the pseudo_inverse curvature at each vertex. - ima_t input(bin_input.domain()); - mln::p_n_faces_fwd_piter<D, G> v(input.domain(), 0); + float_ima_t float_ima(bin_input.domain()); + mln::p_n_faces_fwd_piter<D, G> v(float_ima.domain(), 0); for_all(v) { float h = (curv.first(v) + curv.second(v)) / 2; // Pseudo-inverse curvature. float h_inv = 1 / pi * (atan(-h) + pi / 2); - input(v) = h_inv; + float_ima(v) = h_inv; // FIXME: The program should allow the user to choose the kind // of measure. // input(v) = mln::math::max(mln::math::sqr(curv.first(v)), @@ -101,7 +102,7 @@ int main(int argc, char* argv[]) } // Values on edges. - mln::p_n_faces_fwd_piter<D, G> e(input.domain(), 1); + mln::p_n_faces_fwd_piter<D, G> e(float_ima.domain(), 1); typedef mln::complex_lower_neighborhood<D, G> adj_vertices_nbh_t; adj_vertices_nbh_t adj_vertices_nbh; mln_niter_(adj_vertices_nbh_t) adj_v(adj_vertices_nbh, e); @@ -113,14 +114,26 @@ int main(int argc, char* argv[]) // Iterate on vertices (0-faces). for_all(adj_v) { - s += input(adj_v); + s += float_ima(adj_v); ++n; } - input(e) = s / n; + float_ima(e) = s / n; // An edge should be adjacent to exactly two vertices. mln_invariant(n <= 2); } + // Convert the float image into an unsigned image because some + // algorithms do not handle float images well. + /* FIXME: This is bad: float images should be handled as-is. At + least, we should use a decent conversion, using an optimal affine + transformation (stretching/shrinking) or any other kind of + interpolation. */ + typedef mln::unsigned_2complex_image3df ima_t; + ima_t ima(float_ima.domain()); + mln_piter_(ima_t) p(float_ima.domain()); + for_all(p) + ima(p) = 1000 * float_ima(p); + /*-----------------. | Simplification. | `-----------------*/ @@ -131,7 +144,7 @@ int main(int argc, char* argv[]) adj_edges_nbh_t; adj_edges_nbh_t adj_edges_nbh; - ima_t closed_input = mln::morpho::closing::area(input, adj_edges_nbh, lambda); + ima_t closed_ima = mln::morpho::closing::area(ima, adj_edges_nbh, lambda); /*------. | WST. | @@ -141,7 +154,7 @@ int main(int argc, char* argv[]) not constrained the domain of the image to 1-faces only. It would be great if we could use something like this: - closed_input | mln::p_faces<1, D, G>(closed_input.domain()) + closed_ima | mln::p_faces<1, D, G>(closed_ima.domain()) as input of the WST. */ @@ -150,7 +163,7 @@ int main(int argc, char* argv[]) wst_val_t nbasins; typedef mln::unsigned_2complex_image3df wst_ima_t; wst_ima_t wshed = - mln::morpho::meyer_wst(closed_input, adj_edges_nbh, nbasins); + mln::morpho::meyer_wst(closed_ima, adj_edges_nbh, nbasins); std::cout << "nbasins = " << nbasins << std::endl; // Label polygons (i.e., propagate labels from edges to polygons). diff --git a/milena/apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc b/milena/apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc index 6265a25..a174dd6 100644 --- a/milena/apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc +++ b/milena/apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc @@ -75,37 +75,38 @@ main(int argc, char* argv[]) | Complex image. | `----------------*/ - // Image type. - typedef mln::float_2complex_image3df ima_t; + // Curvature image type. + typedef mln::float_2complex_image3df float_ima_t; // Dimension of the image (and therefore of the complex). - static const unsigned D = ima_t::dim; + static const unsigned D = float_ima_t::dim; // Geometry of the image. - typedef mln_geom_(ima_t) G; + typedef mln_geom_(float_ima_t) G; mln::bin_2complex_image3df bin_input; mln::io::off::load(bin_input, input_filename); - std::pair<ima_t, ima_t> curv = mln::geom::mesh_curvature(bin_input.domain()); + std::pair<float_ima_t, float_ima_t> curv = + mln::geom::mesh_curvature(bin_input.domain()); // Compute the pseudo_inverse curvature at each vertex. - ima_t input(bin_input.domain()); - mln::p_n_faces_fwd_piter<D, G> v(input.domain(), 0); + float_ima_t float_ima(bin_input.domain()); + mln::p_n_faces_fwd_piter<D, G> v(float_ima.domain(), 0); for_all(v) { // FIXME: The program should allow the user to choose the kind - // of measure. + // of measure (curvature). #if 0 float h = (curv.first(v) + curv.second(v)) / 2; // Pseudo-inverse curvature. float h_inv = 1 / pi * (atan(-h) + pi / 2); - input(v) = h_inv; + float_ima(v) = h_inv; #endif // Max curvature. - input(v) = mln::math::max(mln::math::sqr(curv.first(v)), - mln::math::sqr(curv.second(v))); + float_ima(v) = mln::math::max(mln::math::sqr(curv.first(v)), + mln::math::sqr(curv.second(v))); } // Values on triangles. - mln::p_n_faces_fwd_piter<D, G> t(input.domain(), 2); + mln::p_n_faces_fwd_piter<D, G> t(float_ima.domain(), 2); // For each triangle (2-face) T, iterate on the the set of vertices // (0-faces) transitively adjacent to T. typedef mln::complex_m_face_neighborhood<D, G> adj_vertices_nbh_t; @@ -123,21 +124,34 @@ main(int argc, char* argv[]) // Iterate on vertices (0-faces). for_all(adj_v) { - s += input(adj_v); + s += float_ima(adj_v); ++n; } - input(t) = s / n; + float_ima(t) = s / n; // A triangle should be adjacent to exactly two vertices. mln_invariant(n <= 3); } + // Convert the float image into an unsigned image because some + // algorithms do not handle float images well. + /* FIXME: This is bad: float images should be handled as-is. At + least, we should use a decent conversion, using an optimal affine + transformation (stretching/shrinking) or any other kind of + interpolation. */ + typedef mln::unsigned_2complex_image3df ima_t; + ima_t ima(float_ima.domain()); + // Process only triangles, as edges and vertices are set afterwards + // (see FIXME below). + for_all(t) + ima(t) = 1000 * float_ima(t); + /* FIXME: Workaround: Set maximal values on vertices and edges to exclude them from the set of minimal values. */ for_all(v) - input(v) = mln_max(float); - mln::p_n_faces_fwd_piter<D, G> e(input.domain(), 1); + ima(v) = mln_max(mln_value_(ima_t)); + mln::p_n_faces_fwd_piter<D, G> e(float_ima.domain(), 1); for_all(e) - input(e) = mln_max(float); + ima(e) = mln_max(mln_value_(ima_t)); /*-----------------. | Simplification. | @@ -147,7 +161,7 @@ main(int argc, char* argv[]) typedef mln::complex_lower_dim_connected_n_face_neighborhood<D, G> nbh_t; nbh_t nbh; - ima_t closed_input = mln::morpho::closing::area(input, nbh, lambda); + ima_t closed_ima = mln::morpho::closing::area(ima, nbh, lambda); /*---------------. | Local minima. | @@ -161,7 +175,7 @@ main(int argc, char* argv[]) to 2-faces. */ typedef mln_ch_value_(ima_t, label_t) label_ima_t; label_ima_t minima = - mln::labeling::regional_minima(closed_input, nbh, nminima); + mln::labeling::regional_minima(closed_ima, nbh, nminima); typedef mln::complex_higher_neighborhood<D, G> higher_nbh_t; higher_nbh_t higher_nbh; @@ -201,7 +215,7 @@ main(int argc, char* argv[]) for_all(adj_e) if (minima(adj_e) == mln::literal::zero) { - // If V is adjcent to a non-minimal triangle, then it must + // If V is adjacent to a non-minimal triangle, then it must // not belong to a minima. ref_adj_minimum = mln::literal::zero; break; @@ -223,6 +237,13 @@ main(int argc, char* argv[]) | Initial binary image. | `-----------------------*/ + /* Careful: creating ``holes'' in the surface obviously changes its + topology, but it may also split a single connected component in + two or more components, resulting in a disconnected skeleton. We + may want to improve this step either by forbidding any splitting, + or by incrementally ``digging'' a regional minima as long as no + splitting occurs. */ + typedef mln_ch_value_(ima_t, bool) bin_ima_t; bin_ima_t surface(minima.domain()); mln::data::fill(surface, true); diff --git a/milena/apps/mesh-segm-skel/mesh-complex-segm.cc b/milena/apps/mesh-segm-skel/mesh-complex-segm.cc index 8d349fe..655a46c 100644 --- a/milena/apps/mesh-segm-skel/mesh-complex-segm.cc +++ b/milena/apps/mesh-segm-skel/mesh-complex-segm.cc @@ -66,15 +66,27 @@ int main(int argc, char* argv[]) | Complex image. | `----------------*/ - // Image type. - typedef mln::float_2complex_image3df ima_t; + // Input type. + typedef mln::float_2complex_image3df float_input_t; // Dimension of the image (and therefore of the complex). - static const unsigned D = ima_t::dim; + static const unsigned D = float_input_t::dim; // Geometry of the image. - typedef mln_geom_(ima_t) G; - - ima_t input; - mln::io::off::load(input, input_filename); + typedef mln_geom_(float_input_t) G; + + float_input_t float_input; + mln::io::off::load(float_input, input_filename); + + // Convert the float image into an unsigned image because some + // algorithms do not handle float images well. + /* FIXME: This is bad: float images should be handled as-is. At + least, we should use a decent conversion, using an optimal affine + transformation (stretching/shrinking) or any other kind of + interpolation. */ + typedef mln::unsigned_2complex_image3df ima_t; + ima_t input(float_input.domain()); + mln_piter_(ima_t) p(input.domain()); + for_all(p) + input(p) = 1000 * float_input(p); // Values on edges. mln::p_n_faces_fwd_piter<D, G> e(input.domain(), 1); diff --git a/milena/apps/mesh-segm-skel/mesh-complex-skel.cc b/milena/apps/mesh-segm-skel/mesh-complex-skel.cc index b7e1f1f..5cd97b4 100644 --- a/milena/apps/mesh-segm-skel/mesh-complex-skel.cc +++ b/milena/apps/mesh-segm-skel/mesh-complex-skel.cc @@ -72,23 +72,37 @@ main(int argc, char* argv[]) `----------------*/ // Image type. - typedef mln::float_2complex_image3df ima_t; + typedef mln::float_2complex_image3df float_input_t; // Dimension of the image (and therefore of the complex). - static const unsigned D = ima_t::dim; + static const unsigned D = float_input_t::dim; // Geometry of the image. - typedef mln_geom_(ima_t) G; - - ima_t input; - mln::io::off::load(input, input_filename); + typedef mln_geom_(float_input_t) G; + + float_input_t float_input; + mln::io::off::load(float_input, input_filename); + + // Convert the float image into an unsigned image because some + // algorithms do not handle float images well. + /* FIXME: This is bad: float images should be handled as-is. At + least, we should use a decent conversion, using an optimal affine + transformation (stretching/shrinking) or any other kind of + interpolation. */ + typedef mln::unsigned_2complex_image3df ima_t; + ima_t input(float_input.domain()); + // Process only triangles, as edges and vertices are set afterwards + // (see FIXME below). + mln::p_n_faces_fwd_piter<D, G> t(input.domain(), 2); + for_all(t) + input(t) = 1000 * float_input(t); /* FIXME: Workaround: Set maximal values on vertices and edges to exclude them from the set of minimal values. */ mln::p_n_faces_fwd_piter<D, G> v(input.domain(), 0); for_all(v) - input(v) = mln_max(float); + input(v) = mln_max(mln_value_(ima_t)); mln::p_n_faces_fwd_piter<D, G> e(input.domain(), 1); for_all(e) - input(e) = mln_max(float); + input(e) = mln_max(mln_value_(ima_t)); /*-----------------. | Simplification. | @@ -126,7 +140,7 @@ main(int argc, char* argv[]) for_all(adj_t) if (minima(adj_t) == mln::literal::zero) { - // If E is adjcent to a non-minimal triangle, then it must + // If E is adjacent to a non-minimal triangle, then it must // not belong to a minima. ref_adj_minimum = mln::literal::zero; break; @@ -152,7 +166,7 @@ main(int argc, char* argv[]) for_all(adj_e) if (minima(adj_e) == mln::literal::zero) { - // If V is adjcent to a non-minimal triangle, then it must + // If V is adjacent to a non-minimal triangle, then it must // not belong to a minima. ref_adj_minimum = mln::literal::zero; break; @@ -174,6 +188,13 @@ main(int argc, char* argv[]) | Initial binary image. | `-----------------------*/ + /* Careful: creating ``holes'' in the surface obviously changes its + topology, but it may also split a single connected component in + two or more components, resulting in a disconnected skeleton. We + may want to improve this step either by forbidding any splitting, + or by incrementally ``digging'' a regional minima as long as no + splitting occurs. */ + typedef mln_ch_value_(ima_t, bool) bin_ima_t; bin_ima_t surface(minima.domain()); mln::data::fill(surface, true); -- 1.5.6.5

1 0

last-svn-commit-187-g3f6e1c5 Improve the precision of some mesh-segm-skel apps.
by Roland Levillain 14 Jun '10

14 Jun '10

* apps/mesh-segm-skel/mesh-complex-segm.cc, * apps/mesh-segm-skel/mesh-complex-skel.cc, * apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc, * apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc: Convert the floating point value input to an unsigned integer image to avoid precision errors due to float comparisons. --- milena/ChangeLog | 11 ++++ .../mesh-segm-skel/mesh-complex-max-curv-segm.cc | 39 ++++++++---- .../mesh-segm-skel/mesh-complex-max-curv-skel.cc | 61 +++++++++++++------- milena/apps/mesh-segm-skel/mesh-complex-segm.cc | 24 ++++++-- milena/apps/mesh-segm-skel/mesh-complex-skel.cc | 41 ++++++++++--- 5 files changed, 127 insertions(+), 49 deletions(-) diff --git a/milena/ChangeLog b/milena/ChangeLog index bb52bf8..00afa02 100644 --- a/milena/ChangeLog +++ b/milena/ChangeLog @@ -1,5 +1,16 @@ 2010-06-14 Roland Levillain <roland(a)lrde.epita.fr> + Improve the precision of some mesh-segm-skel apps. + + * apps/mesh-segm-skel/mesh-complex-segm.cc, + * apps/mesh-segm-skel/mesh-complex-skel.cc, + * apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc, + * apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc: + Convert the floating point value input to an unsigned integer + image to avoid precision errors due to float comparisons. + +2010-06-14 Roland Levillain <roland(a)lrde.epita.fr> + New application: mesh-complex-max-curv-skel. * apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc: New. diff --git a/milena/apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc b/milena/apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc index 477ccdf..1f43381 100644 --- a/milena/apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc +++ b/milena/apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc @@ -75,25 +75,26 @@ int main(int argc, char* argv[]) `----------------*/ // Image type. - typedef mln::float_2complex_image3df ima_t; + typedef mln::float_2complex_image3df float_ima_t; // Dimension of the image (and therefore of the complex). - static const unsigned D = ima_t::dim; + static const unsigned D = float_ima_t::dim; // Geometry of the image. - typedef mln_geom_(ima_t) G; + typedef mln_geom_(float_ima_t) G; mln::bin_2complex_image3df bin_input; mln::io::off::load(bin_input, input_filename); - std::pair<ima_t, ima_t> curv = mln::geom::mesh_curvature(bin_input.domain()); + std::pair<float_ima_t, float_ima_t> curv = + mln::geom::mesh_curvature(bin_input.domain()); // Compute the pseudo_inverse curvature at each vertex. - ima_t input(bin_input.domain()); - mln::p_n_faces_fwd_piter<D, G> v(input.domain(), 0); + float_ima_t float_ima(bin_input.domain()); + mln::p_n_faces_fwd_piter<D, G> v(float_ima.domain(), 0); for_all(v) { float h = (curv.first(v) + curv.second(v)) / 2; // Pseudo-inverse curvature. float h_inv = 1 / pi * (atan(-h) + pi / 2); - input(v) = h_inv; + float_ima(v) = h_inv; // FIXME: The program should allow the user to choose the kind // of measure. // input(v) = mln::math::max(mln::math::sqr(curv.first(v)), @@ -101,7 +102,7 @@ int main(int argc, char* argv[]) } // Values on edges. - mln::p_n_faces_fwd_piter<D, G> e(input.domain(), 1); + mln::p_n_faces_fwd_piter<D, G> e(float_ima.domain(), 1); typedef mln::complex_lower_neighborhood<D, G> adj_vertices_nbh_t; adj_vertices_nbh_t adj_vertices_nbh; mln_niter_(adj_vertices_nbh_t) adj_v(adj_vertices_nbh, e); @@ -113,14 +114,26 @@ int main(int argc, char* argv[]) // Iterate on vertices (0-faces). for_all(adj_v) { - s += input(adj_v); + s += float_ima(adj_v); ++n; } - input(e) = s / n; + float_ima(e) = s / n; // An edge should be adjacent to exactly two vertices. mln_invariant(n <= 2); } + // Convert the float image into an unsigned image because some + // algorithms do not handle float images well. + /* FIXME: This is bad: float images should be handled as-is. At + least, we should use a decent conversion, using an optimal affine + transformation (stretching/shrinking) or any other kind of + interpolation. */ + typedef mln::unsigned_2complex_image3df ima_t; + ima_t ima(float_ima.domain()); + mln_piter_(ima_t) p(float_ima.domain()); + for_all(p) + ima(p) = 1000 * float_ima(p); + /*-----------------. | Simplification. | `-----------------*/ @@ -131,7 +144,7 @@ int main(int argc, char* argv[]) adj_edges_nbh_t; adj_edges_nbh_t adj_edges_nbh; - ima_t closed_input = mln::morpho::closing::area(input, adj_edges_nbh, lambda); + ima_t closed_ima = mln::morpho::closing::area(ima, adj_edges_nbh, lambda); /*------. | WST. | @@ -141,7 +154,7 @@ int main(int argc, char* argv[]) not constrained the domain of the image to 1-faces only. It would be great if we could use something like this: - closed_input | mln::p_faces<1, D, G>(closed_input.domain()) + closed_ima | mln::p_faces<1, D, G>(closed_ima.domain()) as input of the WST. */ @@ -150,7 +163,7 @@ int main(int argc, char* argv[]) wst_val_t nbasins; typedef mln::unsigned_2complex_image3df wst_ima_t; wst_ima_t wshed = - mln::morpho::meyer_wst(closed_input, adj_edges_nbh, nbasins); + mln::morpho::meyer_wst(closed_ima, adj_edges_nbh, nbasins); std::cout << "nbasins = " << nbasins << std::endl; // Label polygons (i.e., propagate labels from edges to polygons). diff --git a/milena/apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc b/milena/apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc index 6265a25..a174dd6 100644 --- a/milena/apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc +++ b/milena/apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc @@ -75,37 +75,38 @@ main(int argc, char* argv[]) | Complex image. | `----------------*/ - // Image type. - typedef mln::float_2complex_image3df ima_t; + // Curvature image type. + typedef mln::float_2complex_image3df float_ima_t; // Dimension of the image (and therefore of the complex). - static const unsigned D = ima_t::dim; + static const unsigned D = float_ima_t::dim; // Geometry of the image. - typedef mln_geom_(ima_t) G; + typedef mln_geom_(float_ima_t) G; mln::bin_2complex_image3df bin_input; mln::io::off::load(bin_input, input_filename); - std::pair<ima_t, ima_t> curv = mln::geom::mesh_curvature(bin_input.domain()); + std::pair<float_ima_t, float_ima_t> curv = + mln::geom::mesh_curvature(bin_input.domain()); // Compute the pseudo_inverse curvature at each vertex. - ima_t input(bin_input.domain()); - mln::p_n_faces_fwd_piter<D, G> v(input.domain(), 0); + float_ima_t float_ima(bin_input.domain()); + mln::p_n_faces_fwd_piter<D, G> v(float_ima.domain(), 0); for_all(v) { // FIXME: The program should allow the user to choose the kind - // of measure. + // of measure (curvature). #if 0 float h = (curv.first(v) + curv.second(v)) / 2; // Pseudo-inverse curvature. float h_inv = 1 / pi * (atan(-h) + pi / 2); - input(v) = h_inv; + float_ima(v) = h_inv; #endif // Max curvature. - input(v) = mln::math::max(mln::math::sqr(curv.first(v)), - mln::math::sqr(curv.second(v))); + float_ima(v) = mln::math::max(mln::math::sqr(curv.first(v)), + mln::math::sqr(curv.second(v))); } // Values on triangles. - mln::p_n_faces_fwd_piter<D, G> t(input.domain(), 2); + mln::p_n_faces_fwd_piter<D, G> t(float_ima.domain(), 2); // For each triangle (2-face) T, iterate on the the set of vertices // (0-faces) transitively adjacent to T. typedef mln::complex_m_face_neighborhood<D, G> adj_vertices_nbh_t; @@ -123,21 +124,34 @@ main(int argc, char* argv[]) // Iterate on vertices (0-faces). for_all(adj_v) { - s += input(adj_v); + s += float_ima(adj_v); ++n; } - input(t) = s / n; + float_ima(t) = s / n; // A triangle should be adjacent to exactly two vertices. mln_invariant(n <= 3); } + // Convert the float image into an unsigned image because some + // algorithms do not handle float images well. + /* FIXME: This is bad: float images should be handled as-is. At + least, we should use a decent conversion, using an optimal affine + transformation (stretching/shrinking) or any other kind of + interpolation. */ + typedef mln::unsigned_2complex_image3df ima_t; + ima_t ima(float_ima.domain()); + // Process only triangles, as edges and vertices are set afterwards + // (see FIXME below). + for_all(t) + ima(t) = 1000 * float_ima(t); + /* FIXME: Workaround: Set maximal values on vertices and edges to exclude them from the set of minimal values. */ for_all(v) - input(v) = mln_max(float); - mln::p_n_faces_fwd_piter<D, G> e(input.domain(), 1); + ima(v) = mln_max(mln_value_(ima_t)); + mln::p_n_faces_fwd_piter<D, G> e(float_ima.domain(), 1); for_all(e) - input(e) = mln_max(float); + ima(e) = mln_max(mln_value_(ima_t)); /*-----------------. | Simplification. | @@ -147,7 +161,7 @@ main(int argc, char* argv[]) typedef mln::complex_lower_dim_connected_n_face_neighborhood<D, G> nbh_t; nbh_t nbh; - ima_t closed_input = mln::morpho::closing::area(input, nbh, lambda); + ima_t closed_ima = mln::morpho::closing::area(ima, nbh, lambda); /*---------------. | Local minima. | @@ -161,7 +175,7 @@ main(int argc, char* argv[]) to 2-faces. */ typedef mln_ch_value_(ima_t, label_t) label_ima_t; label_ima_t minima = - mln::labeling::regional_minima(closed_input, nbh, nminima); + mln::labeling::regional_minima(closed_ima, nbh, nminima); typedef mln::complex_higher_neighborhood<D, G> higher_nbh_t; higher_nbh_t higher_nbh; @@ -201,7 +215,7 @@ main(int argc, char* argv[]) for_all(adj_e) if (minima(adj_e) == mln::literal::zero) { - // If V is adjcent to a non-minimal triangle, then it must + // If V is adjacent to a non-minimal triangle, then it must // not belong to a minima. ref_adj_minimum = mln::literal::zero; break; @@ -223,6 +237,13 @@ main(int argc, char* argv[]) | Initial binary image. | `-----------------------*/ + /* Careful: creating ``holes'' in the surface obviously changes its + topology, but it may also split a single connected component in + two or more components, resulting in a disconnected skeleton. We + may want to improve this step either by forbidding any splitting, + or by incrementally ``digging'' a regional minima as long as no + splitting occurs. */ + typedef mln_ch_value_(ima_t, bool) bin_ima_t; bin_ima_t surface(minima.domain()); mln::data::fill(surface, true); diff --git a/milena/apps/mesh-segm-skel/mesh-complex-segm.cc b/milena/apps/mesh-segm-skel/mesh-complex-segm.cc index 8d349fe..655a46c 100644 --- a/milena/apps/mesh-segm-skel/mesh-complex-segm.cc +++ b/milena/apps/mesh-segm-skel/mesh-complex-segm.cc @@ -66,15 +66,27 @@ int main(int argc, char* argv[]) | Complex image. | `----------------*/ - // Image type. - typedef mln::float_2complex_image3df ima_t; + // Input type. + typedef mln::float_2complex_image3df float_input_t; // Dimension of the image (and therefore of the complex). - static const unsigned D = ima_t::dim; + static const unsigned D = float_input_t::dim; // Geometry of the image. - typedef mln_geom_(ima_t) G; - - ima_t input; - mln::io::off::load(input, input_filename); + typedef mln_geom_(float_input_t) G; + + float_input_t float_input; + mln::io::off::load(float_input, input_filename); + + // Convert the float image into an unsigned image because some + // algorithms do not handle float images well. + /* FIXME: This is bad: float images should be handled as-is. At + least, we should use a decent conversion, using an optimal affine + transformation (stretching/shrinking) or any other kind of + interpolation. */ + typedef mln::unsigned_2complex_image3df ima_t; + ima_t input(float_input.domain()); + mln_piter_(ima_t) p(input.domain()); + for_all(p) + input(p) = 1000 * float_input(p); // Values on edges. mln::p_n_faces_fwd_piter<D, G> e(input.domain(), 1); diff --git a/milena/apps/mesh-segm-skel/mesh-complex-skel.cc b/milena/apps/mesh-segm-skel/mesh-complex-skel.cc index b7e1f1f..5cd97b4 100644 --- a/milena/apps/mesh-segm-skel/mesh-complex-skel.cc +++ b/milena/apps/mesh-segm-skel/mesh-complex-skel.cc @@ -72,23 +72,37 @@ main(int argc, char* argv[]) `----------------*/ // Image type. - typedef mln::float_2complex_image3df ima_t; + typedef mln::float_2complex_image3df float_input_t; // Dimension of the image (and therefore of the complex). - static const unsigned D = ima_t::dim; + static const unsigned D = float_input_t::dim; // Geometry of the image. - typedef mln_geom_(ima_t) G; - - ima_t input; - mln::io::off::load(input, input_filename); + typedef mln_geom_(float_input_t) G; + + float_input_t float_input; + mln::io::off::load(float_input, input_filename); + + // Convert the float image into an unsigned image because some + // algorithms do not handle float images well. + /* FIXME: This is bad: float images should be handled as-is. At + least, we should use a decent conversion, using an optimal affine + transformation (stretching/shrinking) or any other kind of + interpolation. */ + typedef mln::unsigned_2complex_image3df ima_t; + ima_t input(float_input.domain()); + // Process only triangles, as edges and vertices are set afterwards + // (see FIXME below). + mln::p_n_faces_fwd_piter<D, G> t(input.domain(), 2); + for_all(t) + input(t) = 1000 * float_input(t); /* FIXME: Workaround: Set maximal values on vertices and edges to exclude them from the set of minimal values. */ mln::p_n_faces_fwd_piter<D, G> v(input.domain(), 0); for_all(v) - input(v) = mln_max(float); + input(v) = mln_max(mln_value_(ima_t)); mln::p_n_faces_fwd_piter<D, G> e(input.domain(), 1); for_all(e) - input(e) = mln_max(float); + input(e) = mln_max(mln_value_(ima_t)); /*-----------------. | Simplification. | @@ -126,7 +140,7 @@ main(int argc, char* argv[]) for_all(adj_t) if (minima(adj_t) == mln::literal::zero) { - // If E is adjcent to a non-minimal triangle, then it must + // If E is adjacent to a non-minimal triangle, then it must // not belong to a minima. ref_adj_minimum = mln::literal::zero; break; @@ -152,7 +166,7 @@ main(int argc, char* argv[]) for_all(adj_e) if (minima(adj_e) == mln::literal::zero) { - // If V is adjcent to a non-minimal triangle, then it must + // If V is adjacent to a non-minimal triangle, then it must // not belong to a minima. ref_adj_minimum = mln::literal::zero; break; @@ -174,6 +188,13 @@ main(int argc, char* argv[]) | Initial binary image. | `-----------------------*/ + /* Careful: creating ``holes'' in the surface obviously changes its + topology, but it may also split a single connected component in + two or more components, resulting in a disconnected skeleton. We + may want to improve this step either by forbidding any splitting, + or by incrementally ``digging'' a regional minima as long as no + splitting occurs. */ + typedef mln_ch_value_(ima_t, bool) bin_ima_t; bin_ima_t surface(minima.domain()); mln::data::fill(surface, true); -- 1.5.6.5

1 0

last-svn-commit-187-gca0ba34 Improve the precision of some mesh-segm-skel apps.
by Roland Levillain 14 Jun '10

14 Jun '10

* apps/mesh-segm-skel/mesh-complex-segm.cc, * apps/mesh-segm-skel/mesh-complex-skel.cc, * apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc, * apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc: Convert the floating point value input to an unsigned integer image to avoid precision errors due to float comparisons. --- milena/ChangeLog | 11 ++++ .../mesh-segm-skel/mesh-complex-max-curv-segm.cc | 39 ++++++++---- .../mesh-segm-skel/mesh-complex-max-curv-skel.cc | 61 +++++++++++++------- milena/apps/mesh-segm-skel/mesh-complex-segm.cc | 24 ++++++-- milena/apps/mesh-segm-skel/mesh-complex-skel.cc | 41 ++++++++++--- 5 files changed, 127 insertions(+), 49 deletions(-) diff --git a/milena/ChangeLog b/milena/ChangeLog index 009db1b..8d0f1be 100644 --- a/milena/ChangeLog +++ b/milena/ChangeLog @@ -1,5 +1,16 @@ 2010-06-14 Roland Levillain <roland(a)lrde.epita.fr> + Improve the precision of some mesh-segm-skel apps. + + * apps/mesh-segm-skel/mesh-complex-segm.cc, + * apps/mesh-segm-skel/mesh-complex-skel.cc, + * apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc, + * apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc: + Convert the floating point value input to an unsigned integer + image to avoid precision errors due to float comparisons. + +2010-06-14 Roland Levillain <roland(a)lrde.epita.fr> + New application: mesh-complex-max-curv-skel. * apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc: New. diff --git a/milena/apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc b/milena/apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc index 477ccdf..1f43381 100644 --- a/milena/apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc +++ b/milena/apps/mesh-segm-skel/mesh-complex-max-curv-segm.cc @@ -75,25 +75,26 @@ int main(int argc, char* argv[]) `----------------*/ // Image type. - typedef mln::float_2complex_image3df ima_t; + typedef mln::float_2complex_image3df float_ima_t; // Dimension of the image (and therefore of the complex). - static const unsigned D = ima_t::dim; + static const unsigned D = float_ima_t::dim; // Geometry of the image. - typedef mln_geom_(ima_t) G; + typedef mln_geom_(float_ima_t) G; mln::bin_2complex_image3df bin_input; mln::io::off::load(bin_input, input_filename); - std::pair<ima_t, ima_t> curv = mln::geom::mesh_curvature(bin_input.domain()); + std::pair<float_ima_t, float_ima_t> curv = + mln::geom::mesh_curvature(bin_input.domain()); // Compute the pseudo_inverse curvature at each vertex. - ima_t input(bin_input.domain()); - mln::p_n_faces_fwd_piter<D, G> v(input.domain(), 0); + float_ima_t float_ima(bin_input.domain()); + mln::p_n_faces_fwd_piter<D, G> v(float_ima.domain(), 0); for_all(v) { float h = (curv.first(v) + curv.second(v)) / 2; // Pseudo-inverse curvature. float h_inv = 1 / pi * (atan(-h) + pi / 2); - input(v) = h_inv; + float_ima(v) = h_inv; // FIXME: The program should allow the user to choose the kind // of measure. // input(v) = mln::math::max(mln::math::sqr(curv.first(v)), @@ -101,7 +102,7 @@ int main(int argc, char* argv[]) } // Values on edges. - mln::p_n_faces_fwd_piter<D, G> e(input.domain(), 1); + mln::p_n_faces_fwd_piter<D, G> e(float_ima.domain(), 1); typedef mln::complex_lower_neighborhood<D, G> adj_vertices_nbh_t; adj_vertices_nbh_t adj_vertices_nbh; mln_niter_(adj_vertices_nbh_t) adj_v(adj_vertices_nbh, e); @@ -113,14 +114,26 @@ int main(int argc, char* argv[]) // Iterate on vertices (0-faces). for_all(adj_v) { - s += input(adj_v); + s += float_ima(adj_v); ++n; } - input(e) = s / n; + float_ima(e) = s / n; // An edge should be adjacent to exactly two vertices. mln_invariant(n <= 2); } + // Convert the float image into an unsigned image because some + // algorithms do not handle float images well. + /* FIXME: This is bad: float images should be handled as-is. At + least, we should use a decent conversion, using an optimal affine + transformation (stretching/shrinking) or any other kind of + interpolation. */ + typedef mln::unsigned_2complex_image3df ima_t; + ima_t ima(float_ima.domain()); + mln_piter_(ima_t) p(float_ima.domain()); + for_all(p) + ima(p) = 1000 * float_ima(p); + /*-----------------. | Simplification. | `-----------------*/ @@ -131,7 +144,7 @@ int main(int argc, char* argv[]) adj_edges_nbh_t; adj_edges_nbh_t adj_edges_nbh; - ima_t closed_input = mln::morpho::closing::area(input, adj_edges_nbh, lambda); + ima_t closed_ima = mln::morpho::closing::area(ima, adj_edges_nbh, lambda); /*------. | WST. | @@ -141,7 +154,7 @@ int main(int argc, char* argv[]) not constrained the domain of the image to 1-faces only. It would be great if we could use something like this: - closed_input | mln::p_faces<1, D, G>(closed_input.domain()) + closed_ima | mln::p_faces<1, D, G>(closed_ima.domain()) as input of the WST. */ @@ -150,7 +163,7 @@ int main(int argc, char* argv[]) wst_val_t nbasins; typedef mln::unsigned_2complex_image3df wst_ima_t; wst_ima_t wshed = - mln::morpho::meyer_wst(closed_input, adj_edges_nbh, nbasins); + mln::morpho::meyer_wst(closed_ima, adj_edges_nbh, nbasins); std::cout << "nbasins = " << nbasins << std::endl; // Label polygons (i.e., propagate labels from edges to polygons). diff --git a/milena/apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc b/milena/apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc index 6265a25..a174dd6 100644 --- a/milena/apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc +++ b/milena/apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc @@ -75,37 +75,38 @@ main(int argc, char* argv[]) | Complex image. | `----------------*/ - // Image type. - typedef mln::float_2complex_image3df ima_t; + // Curvature image type. + typedef mln::float_2complex_image3df float_ima_t; // Dimension of the image (and therefore of the complex). - static const unsigned D = ima_t::dim; + static const unsigned D = float_ima_t::dim; // Geometry of the image. - typedef mln_geom_(ima_t) G; + typedef mln_geom_(float_ima_t) G; mln::bin_2complex_image3df bin_input; mln::io::off::load(bin_input, input_filename); - std::pair<ima_t, ima_t> curv = mln::geom::mesh_curvature(bin_input.domain()); + std::pair<float_ima_t, float_ima_t> curv = + mln::geom::mesh_curvature(bin_input.domain()); // Compute the pseudo_inverse curvature at each vertex. - ima_t input(bin_input.domain()); - mln::p_n_faces_fwd_piter<D, G> v(input.domain(), 0); + float_ima_t float_ima(bin_input.domain()); + mln::p_n_faces_fwd_piter<D, G> v(float_ima.domain(), 0); for_all(v) { // FIXME: The program should allow the user to choose the kind - // of measure. + // of measure (curvature). #if 0 float h = (curv.first(v) + curv.second(v)) / 2; // Pseudo-inverse curvature. float h_inv = 1 / pi * (atan(-h) + pi / 2); - input(v) = h_inv; + float_ima(v) = h_inv; #endif // Max curvature. - input(v) = mln::math::max(mln::math::sqr(curv.first(v)), - mln::math::sqr(curv.second(v))); + float_ima(v) = mln::math::max(mln::math::sqr(curv.first(v)), + mln::math::sqr(curv.second(v))); } // Values on triangles. - mln::p_n_faces_fwd_piter<D, G> t(input.domain(), 2); + mln::p_n_faces_fwd_piter<D, G> t(float_ima.domain(), 2); // For each triangle (2-face) T, iterate on the the set of vertices // (0-faces) transitively adjacent to T. typedef mln::complex_m_face_neighborhood<D, G> adj_vertices_nbh_t; @@ -123,21 +124,34 @@ main(int argc, char* argv[]) // Iterate on vertices (0-faces). for_all(adj_v) { - s += input(adj_v); + s += float_ima(adj_v); ++n; } - input(t) = s / n; + float_ima(t) = s / n; // A triangle should be adjacent to exactly two vertices. mln_invariant(n <= 3); } + // Convert the float image into an unsigned image because some + // algorithms do not handle float images well. + /* FIXME: This is bad: float images should be handled as-is. At + least, we should use a decent conversion, using an optimal affine + transformation (stretching/shrinking) or any other kind of + interpolation. */ + typedef mln::unsigned_2complex_image3df ima_t; + ima_t ima(float_ima.domain()); + // Process only triangles, as edges and vertices are set afterwards + // (see FIXME below). + for_all(t) + ima(t) = 1000 * float_ima(t); + /* FIXME: Workaround: Set maximal values on vertices and edges to exclude them from the set of minimal values. */ for_all(v) - input(v) = mln_max(float); - mln::p_n_faces_fwd_piter<D, G> e(input.domain(), 1); + ima(v) = mln_max(mln_value_(ima_t)); + mln::p_n_faces_fwd_piter<D, G> e(float_ima.domain(), 1); for_all(e) - input(e) = mln_max(float); + ima(e) = mln_max(mln_value_(ima_t)); /*-----------------. | Simplification. | @@ -147,7 +161,7 @@ main(int argc, char* argv[]) typedef mln::complex_lower_dim_connected_n_face_neighborhood<D, G> nbh_t; nbh_t nbh; - ima_t closed_input = mln::morpho::closing::area(input, nbh, lambda); + ima_t closed_ima = mln::morpho::closing::area(ima, nbh, lambda); /*---------------. | Local minima. | @@ -161,7 +175,7 @@ main(int argc, char* argv[]) to 2-faces. */ typedef mln_ch_value_(ima_t, label_t) label_ima_t; label_ima_t minima = - mln::labeling::regional_minima(closed_input, nbh, nminima); + mln::labeling::regional_minima(closed_ima, nbh, nminima); typedef mln::complex_higher_neighborhood<D, G> higher_nbh_t; higher_nbh_t higher_nbh; @@ -201,7 +215,7 @@ main(int argc, char* argv[]) for_all(adj_e) if (minima(adj_e) == mln::literal::zero) { - // If V is adjcent to a non-minimal triangle, then it must + // If V is adjacent to a non-minimal triangle, then it must // not belong to a minima. ref_adj_minimum = mln::literal::zero; break; @@ -223,6 +237,13 @@ main(int argc, char* argv[]) | Initial binary image. | `-----------------------*/ + /* Careful: creating ``holes'' in the surface obviously changes its + topology, but it may also split a single connected component in + two or more components, resulting in a disconnected skeleton. We + may want to improve this step either by forbidding any splitting, + or by incrementally ``digging'' a regional minima as long as no + splitting occurs. */ + typedef mln_ch_value_(ima_t, bool) bin_ima_t; bin_ima_t surface(minima.domain()); mln::data::fill(surface, true); diff --git a/milena/apps/mesh-segm-skel/mesh-complex-segm.cc b/milena/apps/mesh-segm-skel/mesh-complex-segm.cc index 8d349fe..655a46c 100644 --- a/milena/apps/mesh-segm-skel/mesh-complex-segm.cc +++ b/milena/apps/mesh-segm-skel/mesh-complex-segm.cc @@ -66,15 +66,27 @@ int main(int argc, char* argv[]) | Complex image. | `----------------*/ - // Image type. - typedef mln::float_2complex_image3df ima_t; + // Input type. + typedef mln::float_2complex_image3df float_input_t; // Dimension of the image (and therefore of the complex). - static const unsigned D = ima_t::dim; + static const unsigned D = float_input_t::dim; // Geometry of the image. - typedef mln_geom_(ima_t) G; - - ima_t input; - mln::io::off::load(input, input_filename); + typedef mln_geom_(float_input_t) G; + + float_input_t float_input; + mln::io::off::load(float_input, input_filename); + + // Convert the float image into an unsigned image because some + // algorithms do not handle float images well. + /* FIXME: This is bad: float images should be handled as-is. At + least, we should use a decent conversion, using an optimal affine + transformation (stretching/shrinking) or any other kind of + interpolation. */ + typedef mln::unsigned_2complex_image3df ima_t; + ima_t input(float_input.domain()); + mln_piter_(ima_t) p(input.domain()); + for_all(p) + input(p) = 1000 * float_input(p); // Values on edges. mln::p_n_faces_fwd_piter<D, G> e(input.domain(), 1); diff --git a/milena/apps/mesh-segm-skel/mesh-complex-skel.cc b/milena/apps/mesh-segm-skel/mesh-complex-skel.cc index b7e1f1f..5cd97b4 100644 --- a/milena/apps/mesh-segm-skel/mesh-complex-skel.cc +++ b/milena/apps/mesh-segm-skel/mesh-complex-skel.cc @@ -72,23 +72,37 @@ main(int argc, char* argv[]) `----------------*/ // Image type. - typedef mln::float_2complex_image3df ima_t; + typedef mln::float_2complex_image3df float_input_t; // Dimension of the image (and therefore of the complex). - static const unsigned D = ima_t::dim; + static const unsigned D = float_input_t::dim; // Geometry of the image. - typedef mln_geom_(ima_t) G; - - ima_t input; - mln::io::off::load(input, input_filename); + typedef mln_geom_(float_input_t) G; + + float_input_t float_input; + mln::io::off::load(float_input, input_filename); + + // Convert the float image into an unsigned image because some + // algorithms do not handle float images well. + /* FIXME: This is bad: float images should be handled as-is. At + least, we should use a decent conversion, using an optimal affine + transformation (stretching/shrinking) or any other kind of + interpolation. */ + typedef mln::unsigned_2complex_image3df ima_t; + ima_t input(float_input.domain()); + // Process only triangles, as edges and vertices are set afterwards + // (see FIXME below). + mln::p_n_faces_fwd_piter<D, G> t(input.domain(), 2); + for_all(t) + input(t) = 1000 * float_input(t); /* FIXME: Workaround: Set maximal values on vertices and edges to exclude them from the set of minimal values. */ mln::p_n_faces_fwd_piter<D, G> v(input.domain(), 0); for_all(v) - input(v) = mln_max(float); + input(v) = mln_max(mln_value_(ima_t)); mln::p_n_faces_fwd_piter<D, G> e(input.domain(), 1); for_all(e) - input(e) = mln_max(float); + input(e) = mln_max(mln_value_(ima_t)); /*-----------------. | Simplification. | @@ -126,7 +140,7 @@ main(int argc, char* argv[]) for_all(adj_t) if (minima(adj_t) == mln::literal::zero) { - // If E is adjcent to a non-minimal triangle, then it must + // If E is adjacent to a non-minimal triangle, then it must // not belong to a minima. ref_adj_minimum = mln::literal::zero; break; @@ -152,7 +166,7 @@ main(int argc, char* argv[]) for_all(adj_e) if (minima(adj_e) == mln::literal::zero) { - // If V is adjcent to a non-minimal triangle, then it must + // If V is adjacent to a non-minimal triangle, then it must // not belong to a minima. ref_adj_minimum = mln::literal::zero; break; @@ -174,6 +188,13 @@ main(int argc, char* argv[]) | Initial binary image. | `-----------------------*/ + /* Careful: creating ``holes'' in the surface obviously changes its + topology, but it may also split a single connected component in + two or more components, resulting in a disconnected skeleton. We + may want to improve this step either by forbidding any splitting, + or by incrementally ``digging'' a regional minima as long as no + splitting occurs. */ + typedef mln_ch_value_(ima_t, bool) bin_ima_t; bin_ima_t surface(minima.domain()); mln::data::fill(surface, true); -- 1.5.6.5

1 0

last-svn-commit-182-ge4978ab Configure milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-skel.
by Roland Levillain 14 Jun '10

14 Jun '10

* configure.ac: Here. --- ChangeLog | 6 ++++++ configure.ac | 2 ++ 2 files changed, 8 insertions(+), 0 deletions(-) diff --git a/ChangeLog b/ChangeLog index 6c4d3ff..3799e38 100644 --- a/ChangeLog +++ b/ChangeLog @@ -1,3 +1,9 @@ +2010-06-14 Roland Levillain <roland(a)lrde.epita.fr> + + Configure milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-skel. + + * configure.ac: Here. + 2010-05-11 Roland Levillain <roland(a)lrde.epita.fr> Configure milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-extrema. diff --git a/configure.ac b/configure.ac index 4d6bfaf..8c96e04 100644 --- a/configure.ac +++ b/configure.ac @@ -515,6 +515,8 @@ AC_CONFIG_FILES([milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-segm], [chmod +x milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-segm]) AC_CONFIG_FILES([milena/apps/mesh-segm-skel/test-mesh-complex-skel], [chmod +x milena/apps/mesh-segm-skel/test-mesh-complex-skel]) +AC_CONFIG_FILES([milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-skel], + [chmod +x milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-skel]) AC_CONFIG_FILES([milena/apps/constrained-connectivity/test-constrained-connectivity], [chmod +x milena/apps/constrained-connectivity/test-constrained-connectivity]) -- 1.5.6.5

1 0

last-svn-commit-186-gb5e8235 Configure milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-skel.
by Roland Levillain 14 Jun '10

14 Jun '10

* configure.ac: Here. --- ChangeLog | 6 ++++++ configure.ac | 2 ++ 2 files changed, 8 insertions(+), 0 deletions(-) diff --git a/ChangeLog b/ChangeLog index bfce0d3..6e0e5b5 100644 --- a/ChangeLog +++ b/ChangeLog @@ -1,3 +1,9 @@ +2010-06-14 Roland Levillain <roland(a)lrde.epita.fr> + + Configure milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-skel. + + * configure.ac: Here. + 2010-05-11 Roland Levillain <roland(a)lrde.epita.fr> Configure milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-extrema. diff --git a/configure.ac b/configure.ac index 4d6bfaf..8c96e04 100644 --- a/configure.ac +++ b/configure.ac @@ -515,6 +515,8 @@ AC_CONFIG_FILES([milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-segm], [chmod +x milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-segm]) AC_CONFIG_FILES([milena/apps/mesh-segm-skel/test-mesh-complex-skel], [chmod +x milena/apps/mesh-segm-skel/test-mesh-complex-skel]) +AC_CONFIG_FILES([milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-skel], + [chmod +x milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-skel]) AC_CONFIG_FILES([milena/apps/constrained-connectivity/test-constrained-connectivity], [chmod +x milena/apps/constrained-connectivity/test-constrained-connectivity]) -- 1.5.6.5

1 0

last-svn-commit-186-gf452dc6 Configure milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-skel.
by Roland Levillain 14 Jun '10

14 Jun '10

* configure.ac: Here. --- ChangeLog | 6 ++++++ configure.ac | 2 ++ 2 files changed, 8 insertions(+), 0 deletions(-) diff --git a/ChangeLog b/ChangeLog index bfce0d3..6e0e5b5 100644 --- a/ChangeLog +++ b/ChangeLog @@ -1,3 +1,9 @@ +2010-06-14 Roland Levillain <roland(a)lrde.epita.fr> + + Configure milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-skel. + + * configure.ac: Here. + 2010-05-11 Roland Levillain <roland(a)lrde.epita.fr> Configure milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-extrema. diff --git a/configure.ac b/configure.ac index 4d6bfaf..8c96e04 100644 --- a/configure.ac +++ b/configure.ac @@ -515,6 +515,8 @@ AC_CONFIG_FILES([milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-segm], [chmod +x milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-segm]) AC_CONFIG_FILES([milena/apps/mesh-segm-skel/test-mesh-complex-skel], [chmod +x milena/apps/mesh-segm-skel/test-mesh-complex-skel]) +AC_CONFIG_FILES([milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-skel], + [chmod +x milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-skel]) AC_CONFIG_FILES([milena/apps/constrained-connectivity/test-constrained-connectivity], [chmod +x milena/apps/constrained-connectivity/test-constrained-connectivity]) -- 1.5.6.5

1 0

last-svn-commit-181-g47481b4 New application: mesh-complex-max-curv-skel.
by Roland Levillain 14 Jun '10

14 Jun '10

* apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc: New. * apps/mesh-segm-skel/test-mesh-complex-max-curv-skel.in: New. * apps/mesh-segm-skel/Makefile.am (bin_PROGRAMS): Add mesh-complex-max-curv-skel. (mesh_complex_max_curv_skel_SOURCES): New. (TESTS): Add test-mesh-complex-max-curv-skel. (CLEANFILES): Add socket-complex-max-curv-skel.off and teapot-complex-max-curv-skel.off. --- milena/ChangeLog | 13 ++++ milena/apps/mesh-segm-skel/Makefile.am | 7 ++ ...mplex-skel.cc => mesh-complex-max-curv-skel.cc} | 65 +++++++++++++++++--- ...-segm.in => test-mesh-complex-max-curv-skel.in} | 8 +- 4 files changed, 81 insertions(+), 12 deletions(-) copy milena/apps/mesh-segm-skel/{mesh-complex-skel.cc => mesh-complex-max-curv-skel.cc} (76%) copy milena/apps/mesh-segm-skel/{test-mesh-segm.in => test-mesh-complex-max-curv-skel.in} (66%) mode change 100755 => 100644 diff --git a/milena/ChangeLog b/milena/ChangeLog index 8b45871..3de4207 100644 --- a/milena/ChangeLog +++ b/milena/ChangeLog @@ -1,3 +1,16 @@ +2010-06-14 Roland Levillain <roland(a)lrde.epita.fr> + + New application: mesh-complex-max-curv-skel. + + * apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc: New. + * apps/mesh-segm-skel/test-mesh-complex-max-curv-skel.in: New. + * apps/mesh-segm-skel/Makefile.am (bin_PROGRAMS): + Add mesh-complex-max-curv-skel. + (mesh_complex_max_curv_skel_SOURCES): New. + (TESTS): Add test-mesh-complex-max-curv-skel. + (CLEANFILES): Add socket-complex-max-curv-skel.off and + teapot-complex-max-curv-skel.off. + 2010-05-11 Roland Levillain <roland(a)lrde.epita.fr> Clean up and new FIXMES in in apps/mesh-segm-skel. diff --git a/milena/apps/mesh-segm-skel/Makefile.am b/milena/apps/mesh-segm-skel/Makefile.am index b702ccb..debed78 100644 --- a/milena/apps/mesh-segm-skel/Makefile.am +++ b/milena/apps/mesh-segm-skel/Makefile.am @@ -189,3 +189,10 @@ CLEANFILES += \ # FIXME: Implement `mesh-complex-{max,pinv}-curv-skel'. # ... +bin_PROGRAMS += mesh-complex-max-curv-skel +mesh_complex_max_curv_skel_SOURCES = \ + mesh-complex-max-curv-skel.cc save_bin_alt.hh +TESTS += test-mesh-complex-max-curv-skel +CLEANFILES += \ + socket-complex-max-curv-skel.off \ + teapot-complex-max-curv-skel.off diff --git a/milena/apps/mesh-segm-skel/mesh-complex-skel.cc b/milena/apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc similarity index 76% copy from milena/apps/mesh-segm-skel/mesh-complex-skel.cc copy to milena/apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc index b7e1f1f..6265a25 100644 --- a/milena/apps/mesh-segm-skel/mesh-complex-skel.cc +++ b/milena/apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc @@ -1,4 +1,5 @@ -// Copyright (C) 2008, 2009 EPITA Research and Development Laboratory (LRDE) +// Copyright (C) 2008, 2009, 2010 EPITA Research and Development +// Laboratory (LRDE) // // This file is part of the Milena Library. This library is free // software; you can redistribute it and/or modify it under the terms @@ -25,9 +26,10 @@ // reasons why the executable file might be covered by the GNU General // Public License. -/// \file apps/mesh-segm-skel/mesh-complex-skel.cc -/// \brief A program computing a skeleton of the surface of the -/// (triangle) mesh of a statue (by thinning), using a complex-based +/// \file apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc +/// \brief A program computing the maximal curvature values from the +/// surface of the (triangle) mesh of a statue, then computing a +/// skeleton of this surface (by/ thinning), using a complex-based /// image. #include <iostream> @@ -52,6 +54,8 @@ morphed mln::complex_image (i.e., seen through image_if). */ #include "save_bin_alt.hh" +#include "trimesh/misc.hh" + int main(int argc, char* argv[]) @@ -78,12 +82,57 @@ main(int argc, char* argv[]) // Geometry of the image. typedef mln_geom_(ima_t) G; - ima_t input; - mln::io::off::load(input, input_filename); + mln::bin_2complex_image3df bin_input; + mln::io::off::load(bin_input, input_filename); + std::pair<ima_t, ima_t> curv = mln::geom::mesh_curvature(bin_input.domain()); + + // Compute the pseudo_inverse curvature at each vertex. + ima_t input(bin_input.domain()); + mln::p_n_faces_fwd_piter<D, G> v(input.domain(), 0); + for_all(v) + { + // FIXME: The program should allow the user to choose the kind + // of measure. +#if 0 + float h = (curv.first(v) + curv.second(v)) / 2; + // Pseudo-inverse curvature. + float h_inv = 1 / pi * (atan(-h) + pi / 2); + input(v) = h_inv; +#endif + // Max curvature. + input(v) = mln::math::max(mln::math::sqr(curv.first(v)), + mln::math::sqr(curv.second(v))); + } + + // Values on triangles. + mln::p_n_faces_fwd_piter<D, G> t(input.domain(), 2); + // For each triangle (2-face) T, iterate on the the set of vertices + // (0-faces) transitively adjacent to T. + typedef mln::complex_m_face_neighborhood<D, G> adj_vertices_nbh_t; + adj_vertices_nbh_t adj_vertices_nbh; + mln_niter_(adj_vertices_nbh_t) adj_v(adj_vertices_nbh, t); + /* FIXME: We should be able to pass this value (m) either at the + construction of the neighborhood or at the construction of the + iterator. */ + adj_v.iter().set_m(0); + // Iterate on triangles (2-faces). + for_all(t) + { + float s = 0.0f; + unsigned n = 0; + // Iterate on vertices (0-faces). + for_all(adj_v) + { + s += input(adj_v); + ++n; + } + input(t) = s / n; + // A triangle should be adjacent to exactly two vertices. + mln_invariant(n <= 3); + } /* FIXME: Workaround: Set maximal values on vertices and edges to exclude them from the set of minimal values. */ - mln::p_n_faces_fwd_piter<D, G> v(input.domain(), 0); for_all(v) input(v) = mln_max(float); mln::p_n_faces_fwd_piter<D, G> e(input.domain(), 1); @@ -126,7 +175,7 @@ main(int argc, char* argv[]) for_all(adj_t) if (minima(adj_t) == mln::literal::zero) { - // If E is adjcent to a non-minimal triangle, then it must + // If E is adjacent to a non-minimal triangle, then it must // not belong to a minima. ref_adj_minimum = mln::literal::zero; break; diff --git a/milena/apps/mesh-segm-skel/test-mesh-segm.in b/milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-skel.in old mode 100755 new mode 100644 similarity index 66% copy from milena/apps/mesh-segm-skel/test-mesh-segm.in copy to milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-skel.in index e70e63a..ea7f604 --- a/milena/apps/mesh-segm-skel/test-mesh-segm.in +++ b/milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-skel.in @@ -1,6 +1,6 @@ #! /bin/sh -# Copyright (C) 2008, 2009 EPITA Research and Development Laboratory (LRDE). +# Copyright (C) 2010 EPITA Research and Development Laboratory (LRDE). # # This file is part of Olena. # @@ -20,7 +20,7 @@ set -ex mesh_dir=@top_srcdir@/milena/mesh -time ./mesh-segm $mesh_dir/socket.off 25 socket-segm.off -time ./mesh-segm $mesh_dir/teapot.off 50 teapot-segm.off +time ./mesh-complex-max-curv-skel $mesh_dir/socket.off 10 socket-complex-max-curv-skel.off +time ./mesh-complex-max-curv-skel $mesh_dir/teapot.off 20 teapot-complex-max-curv-skel.off # FIXME: Too long. Should be run with (future) target `check-full'. -#time ./mesh-segm $mesh_dir/bunny-holefilled.off 500 bunny-holefilled-segm.off +#time ./mesh-complex-max-curv-skel $mesh_dir/bunny-holefilled.off 50 bunny-holefilled-complex-max-curv-skel.off -- 1.5.6.5

1 0

last-svn-commit-185-gf1c3f6a New application: mesh-complex-max-curv-skel.
by Roland Levillain 14 Jun '10

14 Jun '10

* apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc: New. * apps/mesh-segm-skel/test-mesh-complex-max-curv-skel.in: New. * apps/mesh-segm-skel/Makefile.am (bin_PROGRAMS): Add mesh-complex-max-curv-skel. (mesh_complex_max_curv_skel_SOURCES): New. (TESTS): Add test-mesh-complex-max-curv-skel. (CLEANFILES): Add socket-complex-max-curv-skel.off and teapot-complex-max-curv-skel.off. --- milena/ChangeLog | 13 ++++ milena/apps/mesh-segm-skel/Makefile.am | 7 ++ ...mplex-skel.cc => mesh-complex-max-curv-skel.cc} | 65 +++++++++++++++++--- ...-segm.in => test-mesh-complex-max-curv-skel.in} | 8 +- 4 files changed, 81 insertions(+), 12 deletions(-) copy milena/apps/mesh-segm-skel/{mesh-complex-skel.cc => mesh-complex-max-curv-skel.cc} (76%) copy milena/apps/mesh-segm-skel/{test-mesh-segm.in => test-mesh-complex-max-curv-skel.in} (66%) mode change 100755 => 100644 diff --git a/milena/ChangeLog b/milena/ChangeLog index 64370e9..bb52bf8 100644 --- a/milena/ChangeLog +++ b/milena/ChangeLog @@ -1,3 +1,16 @@ +2010-06-14 Roland Levillain <roland(a)lrde.epita.fr> + + New application: mesh-complex-max-curv-skel. + + * apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc: New. + * apps/mesh-segm-skel/test-mesh-complex-max-curv-skel.in: New. + * apps/mesh-segm-skel/Makefile.am (bin_PROGRAMS): + Add mesh-complex-max-curv-skel. + (mesh_complex_max_curv_skel_SOURCES): New. + (TESTS): Add test-mesh-complex-max-curv-skel. + (CLEANFILES): Add socket-complex-max-curv-skel.off and + teapot-complex-max-curv-skel.off. + 2010-05-11 Roland Levillain <roland(a)lrde.epita.fr> Clean up and new FIXMES in in apps/mesh-segm-skel. diff --git a/milena/apps/mesh-segm-skel/Makefile.am b/milena/apps/mesh-segm-skel/Makefile.am index 6853f76..e41aa3e 100644 --- a/milena/apps/mesh-segm-skel/Makefile.am +++ b/milena/apps/mesh-segm-skel/Makefile.am @@ -189,3 +189,10 @@ MOSTLYCLEANFILES += \ # FIXME: Implement `mesh-complex-{max,pinv}-curv-skel'. # ... +bin_PROGRAMS += mesh-complex-max-curv-skel +mesh_complex_max_curv_skel_SOURCES = \ + mesh-complex-max-curv-skel.cc save_bin_alt.hh +TESTS += test-mesh-complex-max-curv-skel +CLEANFILES += \ + socket-complex-max-curv-skel.off \ + teapot-complex-max-curv-skel.off diff --git a/milena/apps/mesh-segm-skel/mesh-complex-skel.cc b/milena/apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc similarity index 76% copy from milena/apps/mesh-segm-skel/mesh-complex-skel.cc copy to milena/apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc index b7e1f1f..6265a25 100644 --- a/milena/apps/mesh-segm-skel/mesh-complex-skel.cc +++ b/milena/apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc @@ -1,4 +1,5 @@ -// Copyright (C) 2008, 2009 EPITA Research and Development Laboratory (LRDE) +// Copyright (C) 2008, 2009, 2010 EPITA Research and Development +// Laboratory (LRDE) // // This file is part of the Milena Library. This library is free // software; you can redistribute it and/or modify it under the terms @@ -25,9 +26,10 @@ // reasons why the executable file might be covered by the GNU General // Public License. -/// \file apps/mesh-segm-skel/mesh-complex-skel.cc -/// \brief A program computing a skeleton of the surface of the -/// (triangle) mesh of a statue (by thinning), using a complex-based +/// \file apps/mesh-segm-skel/mesh-complex-max-curv-skel.cc +/// \brief A program computing the maximal curvature values from the +/// surface of the (triangle) mesh of a statue, then computing a +/// skeleton of this surface (by/ thinning), using a complex-based /// image. #include <iostream> @@ -52,6 +54,8 @@ morphed mln::complex_image (i.e., seen through image_if). */ #include "save_bin_alt.hh" +#include "trimesh/misc.hh" + int main(int argc, char* argv[]) @@ -78,12 +82,57 @@ main(int argc, char* argv[]) // Geometry of the image. typedef mln_geom_(ima_t) G; - ima_t input; - mln::io::off::load(input, input_filename); + mln::bin_2complex_image3df bin_input; + mln::io::off::load(bin_input, input_filename); + std::pair<ima_t, ima_t> curv = mln::geom::mesh_curvature(bin_input.domain()); + + // Compute the pseudo_inverse curvature at each vertex. + ima_t input(bin_input.domain()); + mln::p_n_faces_fwd_piter<D, G> v(input.domain(), 0); + for_all(v) + { + // FIXME: The program should allow the user to choose the kind + // of measure. +#if 0 + float h = (curv.first(v) + curv.second(v)) / 2; + // Pseudo-inverse curvature. + float h_inv = 1 / pi * (atan(-h) + pi / 2); + input(v) = h_inv; +#endif + // Max curvature. + input(v) = mln::math::max(mln::math::sqr(curv.first(v)), + mln::math::sqr(curv.second(v))); + } + + // Values on triangles. + mln::p_n_faces_fwd_piter<D, G> t(input.domain(), 2); + // For each triangle (2-face) T, iterate on the the set of vertices + // (0-faces) transitively adjacent to T. + typedef mln::complex_m_face_neighborhood<D, G> adj_vertices_nbh_t; + adj_vertices_nbh_t adj_vertices_nbh; + mln_niter_(adj_vertices_nbh_t) adj_v(adj_vertices_nbh, t); + /* FIXME: We should be able to pass this value (m) either at the + construction of the neighborhood or at the construction of the + iterator. */ + adj_v.iter().set_m(0); + // Iterate on triangles (2-faces). + for_all(t) + { + float s = 0.0f; + unsigned n = 0; + // Iterate on vertices (0-faces). + for_all(adj_v) + { + s += input(adj_v); + ++n; + } + input(t) = s / n; + // A triangle should be adjacent to exactly two vertices. + mln_invariant(n <= 3); + } /* FIXME: Workaround: Set maximal values on vertices and edges to exclude them from the set of minimal values. */ - mln::p_n_faces_fwd_piter<D, G> v(input.domain(), 0); for_all(v) input(v) = mln_max(float); mln::p_n_faces_fwd_piter<D, G> e(input.domain(), 1); @@ -126,7 +175,7 @@ main(int argc, char* argv[]) for_all(adj_t) if (minima(adj_t) == mln::literal::zero) { - // If E is adjcent to a non-minimal triangle, then it must + // If E is adjacent to a non-minimal triangle, then it must // not belong to a minima. ref_adj_minimum = mln::literal::zero; break; diff --git a/milena/apps/mesh-segm-skel/test-mesh-segm.in b/milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-skel.in old mode 100755 new mode 100644 similarity index 66% copy from milena/apps/mesh-segm-skel/test-mesh-segm.in copy to milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-skel.in index e70e63a..ea7f604 --- a/milena/apps/mesh-segm-skel/test-mesh-segm.in +++ b/milena/apps/mesh-segm-skel/test-mesh-complex-max-curv-skel.in @@ -1,6 +1,6 @@ #! /bin/sh -# Copyright (C) 2008, 2009 EPITA Research and Development Laboratory (LRDE). +# Copyright (C) 2010 EPITA Research and Development Laboratory (LRDE). # # This file is part of Olena. # @@ -20,7 +20,7 @@ set -ex mesh_dir=@top_srcdir@/milena/mesh -time ./mesh-segm $mesh_dir/socket.off 25 socket-segm.off -time ./mesh-segm $mesh_dir/teapot.off 50 teapot-segm.off +time ./mesh-complex-max-curv-skel $mesh_dir/socket.off 10 socket-complex-max-curv-skel.off +time ./mesh-complex-max-curv-skel $mesh_dir/teapot.off 20 teapot-complex-max-curv-skel.off # FIXME: Too long. Should be run with (future) target `check-full'. -#time ./mesh-segm $mesh_dir/bunny-holefilled.off 500 bunny-holefilled-segm.off +#time ./mesh-complex-max-curv-skel $mesh_dir/bunny-holefilled.off 50 bunny-holefilled-complex-max-curv-skel.off -- 1.5.6.5

1 0