[Olena-patches] cleanup-2008 2772: Inim: Add the code related to edge detection in the sandbox.

https://svn.lrde.epita.fr/svn/oln/branches/cleanup-2008/milena/sandbox Index: ChangeLog from Nicolas Ballas &lt;ballas(a)lrde.epita.fr&gt; Inim: Add the code related to edge detection in the sandbox. * ballas/color: New. * ballas/color/min_tree_volume_filter.cc: Min tree using volume filter. * ballas/color/reference.cc: Method using a gradient. * ballas/color/min_tree_area_filter.cc: Min tree using area filter. * ballas/color/min_tree_color.cc: Min tree using filter based on color. * ballas/color/reference2.cc: Method using the laplacien. * ballas/color/src: New. * ballas/color/src/graph.hh, * ballas/color/src/io.hh, * ballas/color/src/distance.hh, * ballas/color/src/convert.hh: New, factorize some code. * ballas/color/laplacien.cc: Laplacien tests. laplacien.cc | 129 ++++++++++ min_tree_area_filter.cc | 505 +++++++++++++++++++++++++++++++++++++++++ min_tree_color.cc | 525 +++++++++++++++++++++++++++++++++++++++++++ min_tree_volume_filter.cc | 514 ++++++++++++++++++++++++++++++++++++++++++ reference.cc | 556 ++++++++++++++++++++++++++++++++++++++++++++++ reference2.cc | 406 +++++++++++++++++++++++++++++++++ src/convert.hh | 37 +++ src/distance.hh | 52 ++++ src/graph.hh | 61 +++++ src/io.hh | 57 ++++ 10 files changed, 2842 insertions(+) Index: ballas/color/min_tree_volume_filter.cc --- ballas/color/min_tree_volume_filter.cc (revision 0) +++ ballas/color/min_tree_volume_filter.cc (revision 0) @@ -0,0 +1,514 @@ +# include <mln/core/var.hh> + +# include <mln/core/image/image2d.hh> +# include <mln/core/image/image_if.hh> +# include <mln/core/image/extended.hh> +# include <mln/core/routine/extend.hh> + +# include <mln/core/alias/window2d.hh> +# include <mln/core/alias/neighb2d.hh> +# include <mln/make/double_neighb2d.hh> +# include <mln/core/site_set/p_centered.hh> + +# include <mln/literal/origin.hh> +# include <mln/literal/black.hh> +# include <mln/literal/white.hh> + +# include <mln/value/int_u8.hh> +# include <mln/value/int_u16.hh> +# include <mln/io/pgm/load.hh> +# include <mln/io/pgm/save.hh> + +# include <mln/value/rgb8.hh> +# include <mln/io/ppm/load.hh> +# include <mln/io/ppm/save.hh> + +# include <mln/morpho/closing_area.hh> + +# include <mln/level/paste.hh> +# include <mln/level/fill.hh> +# include <mln/level/transform.hh> +# include <mln/extension/fill.hh> + +# include <mln/debug/println.hh> + +# include "src/distance.hh" + +namespace mln +{ + template <typename I, typename N, typename Ic, typename Nc> + struct min_tree_ + { + typedef mln_site(I) point; + typedef p_array<point> S; + + // in: + const I& f; + const N& nbh; + const Ic& ref; + const Nc& nbhc; + + // aux: + S s; + mln_ch_value(I, bool) deja_vu; + mln_ch_value(I, point) parent; + mln_ch_value(I, bool) resp; + mln_ch_value(I, point) zpar; + + // attached data: + int lambda; + mln_ch_value(I, int) volume; + //mln_ch_value(Ic, value::rgb8) values; + //initialize(values, ref); + //mln_ch_value(I, int) comp; + + min_tree_(const I& f, const N& nbh, const Ic& ref, const Nc& nbhc, + int lambda) + : f(f), + nbh(nbh), + ref(ref), + nbhc(nbhc), + lambda(lambda) + { + run(); + } + + void run() + { + // init + { + initialize(deja_vu, f); + initialize(parent, f); + initialize(resp, f); + initialize(zpar, f); + initialize(volume, f); + //initialize(comp, f); + + mln::level::fill(deja_vu, false); + //mln::level::fill(resp, false); + mln::level::fill(volume, 0); + + s = level::sort_psites_increasing(f); + } + + // first pass + { + mln_fwd_piter(S) p(s); + mln_niter(N) n(nbh, p); + for_all(p) + { + make_set(p); + for_all(n) + if (f.has(n) && deja_vu(n)) + do_union(n, p); + deja_vu(p) = true; + } + } + + // second pass: canonization + { + mln_bkd_piter(S) p(s); + for_all(p) + { + point q = parent(p); + if (f(parent(q)) == f(q)) + { + parent(p) = parent(q); + resp(q) = false; + } + } + } + + // third pass: Merging region with volume < lambda + { + mln_fwd_piter(S) p(s); + for_all(p) + { + if (resp(p) && (volume(p) < lambda)) + { + resp(p) = false; + update_data(parent(p), volume(p)); + } + } + } + + } // end of run() + + void make_set(const point& p) + { + parent(p) = p; + zpar(p) = p; + init_data(p); + } + + void set_parent(const point& r, const point& p) + { + parent(r) = p; + merge_data(r, p); + } + + bool is_root(const point& p) const + { + return parent(p) == p; + } + + bool is_node(const point& p) const + { + //return is_root(p) || f(parent(p)) != f(p); + return (is_root(p) || resp(p)); + } + + point find_root(const point& x) + { + if (zpar(x) == x) + return x; + else + return zpar(x) = find_root(zpar(x)); + } + + point find_representative(const point& x) + { + if (parent(x) == x || resp(x)) + return x; + else + return find_representative(parent(x)); + } + + void do_union(const point& n, const point& p) + { + point r = find_root(n); + if (r != p) + { + set_parent(r, p); + zpar(r) = p; + } + } + + void init_data(const point& p) + { + int red =0, green = 0, blue = 0; + + mln_niter(Nc) n(nbhc, p); + for_all(n) + { + red += ref(n).red(); + green += ref(n).green(); + blue += ref(n).blue(); + } + + red /= 2; + green /= 2; + blue /= 2; + + volume(p) = distance(value::rgb8(red, green, blue), + value::rgb8(0, 0, 0)); + resp(p) = true; + } + + void merge_data(const point& r, const point& p) + { + if (f(p) == f(r)) + { + resp(p) = false; + volume(r) += volume(p); + } + } + + void update_data(const point& p, int val) + { + volume(p) += val; + if (parent(p) != p && !resp(p)) + update_data(parent(p), val); + } + + }; +} + +namespace mln +{ + image2d<value::int_u16> convert_to_grey(const image2d<value::rgb8>& data) + { + image2d<value::int_u16> output(data.domain()); + mln_piter_(image2d<value::int_u16>) p(output.domain()); + for_all(p) + output(p) = (int) (data(p).red() * 0.3 + data(p).green() * 0.58 + data(p).blue()) * 0.12; + return output; + } +} // end of mln + +namespace mln +{ + + struct colorize : Function_v2v< colorize > + { + typedef value::rgb8 result; + colorize(unsigned max) + : lut(max + 1) + { + lut[0] = literal::black; + for (unsigned i = 1; i <= max; ++i) + lut[i] = result(100 + std::rand() % 150, + 100 + std::rand() % 150, + 100 + std::rand() % 150); + } + result operator()(unsigned i) const + { + return lut[i]; + } + std::vector<result> lut; + }; + + template <typename I> + I display_edge(const I& ima, mln_value(I) bg, unsigned zoom) + { + unsigned nrows = ima.nrows() / 2 + 1; + unsigned ncols = ima.ncols() / 2 + 1; + I output(nrows * (zoom + 1) - 1, + ncols * (zoom + 1) - 1); + level::fill(output, bg); + + mln_VAR(edge, ima | is_edge); + mln_piter(edge_t) p(edge.domain()); + for_all(p) + if (p.row() % 2) // horizontal edge + { + unsigned row = (p.row() / 2 + 1) * (zoom + 1) - 1; + unsigned col = (p.col() / 2) * (zoom + 1); + for (unsigned i = 0; i < zoom; ++i) + output.at(row, col + i) = ima(p); + } + else // vertical edge + { + unsigned row = (p.row() / 2) * (zoom + 1); + unsigned col = (p.col() / 2 + 1) * (zoom + 1) - 1; + for (unsigned i = 0; i < zoom; ++i) + output.at(row + i, col) = ima(p); + } + return output; + } + + + template <typename I> + I display_edge(const I& ima, unsigned zoom) + { + unsigned nrows = ima.nrows() / 2 + 1; + unsigned ncols = ima.ncols() / 2 + 1; + I output(nrows * (zoom + 1) - 1, + ncols * (zoom + 1) - 1); + + mln_VAR( cell, ima | is_cell ); + mln_piter(cell_t) q(cell.domain()); + for_all(q) + { + unsigned row = (q.row() / 2) * (zoom + 1); + unsigned col = (q.col() / 2) * (zoom + 1); + for (unsigned i = 0; i < zoom; ++i) + for (unsigned j = 0; j < zoom; ++j) + output.at(row + i, col + j) = ima(q); + } + + mln_VAR( edge, ima | is_edge ); + mln_piter(edge_t) p(edge.domain()); + for_all(p) + if (p.row() % 2) // horizontal edge + { + unsigned row = (p.row() / 2 + 1) * (zoom + 1) - 1; + unsigned col = (p.col() / 2) * (zoom + 1); + for (unsigned i = 0; i < zoom; ++i) + output.at(row, col + i) = ima(p); + } + else // vertical edge + { + unsigned row = (p.row() / 2) * (zoom + 1); + unsigned col = (p.col() / 2 + 1) * (zoom + 1) - 1; + for (unsigned i = 0; i < zoom; ++i) + output.at(row + i, col) = ima(p); + } + return output; + } + + + namespace morpho + { + + template <typename I, typename N> + mln_concrete(I) + dilation(const I& input, const N& nbh) + { + typedef mln_value(I) V; + + mln_concrete(I) output; + initialize(output, input); + + mln_piter(I) p(input.domain()); + mln_niter(N) n(nbh, p); + for_all(p) + { + for_all(n) + if (input.has(n) && input(n) != value::rgb8(0,0,0)) + output(p) = input(n); + } + return output; + } + } // mln::morpho + +} // mln + + + +template <typename T> +mln::image2d<T> +image2cells(const mln::image2d<T>& input) +{ + mln::image2d<T> output(2 * input.nrows() - 1, + 2 * input.ncols() - 1); + for (unsigned row = 0; row < input.nrows(); ++row) + for (unsigned col = 0; col < input.ncols(); ++col) + output.at(2 * row, 2 * col) = input.at(row, col); + return output; +} + + +template <typename T> +mln::image2d<T> +cells2image(const mln::image2d<T>& input) +{ + mln::image2d<T> output((input.nrows() + 1) / 2, + (input.ncols() + 1) / 2); + for (unsigned row = 0; row < input.nrows(); row += 2) + for (unsigned col = 0; col < input.ncols(); col += 2) + output.at(row / 2, col / 2) = input.at(row, col); + return output; +} + + +template <typename I, typename N, typename Ic, typename Nc> +unsigned min_tree(const I& f, const N& nbh, const Ic& ref, const Nc& nbhc, + int lambda) +{ + using namespace mln; + + min_tree_<I,N,Ic,Nc> run(f, nbh, ref, nbhc, lambda); + + mln_piter(I) p(f.domain()); + unsigned nnodes = 0; + for_all(p) + { + if (run.is_node(p)) + { + ++nnodes; + } + } + + colorize colors(nnodes); + image2d<value::rgb8> tmp(ref.domain()); + level::fill(tmp, ref); + + mln_piter(I) q(f.domain()); + unsigned int i = 0; + for_all(q) + { + if (run.is_node(q)) + { + tmp(q) = colors(i); + i++; + } + } + mln_piter(I) r(f.domain()); + for_all(r) + { + if (!run.is_node(r)) + { + tmp(r) = tmp(run.find_representative(r)); + } + } + + image2d<value::rgb8> to_display(tmp.domain()); + + level::fill(to_display, value::rgb8(255, 255, 255)); + level::paste((tmp | is_edge), to_display); + level::paste(morpho::dilation(to_display, c4()), to_display); + + io::ppm::save(display_edge(tmp, literal::black, 3), + "edge.ppm"); + io::ppm::save(tmp, "full.ppm"); + io::ppm::save(cells2image(to_display), "colorize.ppm"); + + + return nnodes; +} + + +template <typename I> +I +do_it(I& input, int lambda, unsigned& nbasins) +{ + using namespace mln; + + /// Graph creation + I graph; + create_graph(input, graph, value::rgb8(0, 0, 0)); + + // Initialization + image2d<value::int_u16> ima = convert_to_grey(graph); + + // Neigbhorhood + // e2c + bool e2c_h[] = { 0, 1, 0, + 0, 0, 0, + 0, 1, 0 }; + bool e2c_v[] = { 0, 0, 0, + 1, 0, 1, + 0, 0, 0 }; + + mln_VAR(e2c, make::double_neighb2d(is_row_odd, e2c_h, e2c_v)); + + bool e2e_h[] = { 0, 0, 1, 0, 0, + 0, 1, 0, 1, 0, + 0, 0, 0, 0, 0, + 0, 1, 0, 1, 0, + 0, 0, 1, 0, 0 }; + + bool e2e_v[] = { 0, 0, 0, 0, 0, + 0, 1, 0, 1, 0, + 1, 0, 0, 0, 1, + 0, 1, 0, 1, 0, + 0, 0, 0, 0, 0 }; + mln_VAR(e2e, make::double_neighb2d(is_row_odd, e2e_h, e2e_v)); + + // Algorithm + distance(extend((graph | is_edge).rw(), pw::value(graph)), e2c, ima); + + io::pgm::save(ima, "edge.pgm"); + + nbasins = min_tree((ima | is_edge), e2e, graph, e2c, lambda); + + return graph; +} + +void usage(char* argv[]) +{ + std::cerr << "usage: " << argv[0] << " input.pgm lambda" << std::endl; + std::cerr << " lambda >= 0" << std::endl; + abort(); +} + +int main(int argc, char* argv[]) +{ + using namespace mln; + + if (argc != 3) + usage(argv); + + int lambda = atoi(argv[2]); + if (lambda < 0) + usage(argv); + + image2d<value::rgb8> ima; + io::ppm::load(ima, argv[1]); + + unsigned nbasins; + image2d<value::rgb8> output = do_it(ima, lambda, nbasins); + + //io::ppm::save(output, argv[3]); +} Index: ballas/color/reference.cc --- ballas/color/reference.cc (revision 0) +++ ballas/color/reference.cc (revision 0) @@ -0,0 +1,556 @@ +# include <mln/core/var.hh> + +# include <mln/core/image/image2d.hh> +# include <mln/core/image/image_if.hh> +# include <mln/core/image/extended.hh> +# include <mln/core/routine/extend.hh> + +# include <mln/core/alias/window2d.hh> +# include <mln/core/alias/neighb2d.hh> +# include <mln/make/double_neighb2d.hh> +# include <mln/core/site_set/p_centered.hh> + +# include <mln/literal/origin.hh> +# include <mln/literal/black.hh> +# include <mln/literal/white.hh> + +# include <mln/value/int_u8.hh> +# include <mln/io/pgm/load.hh> +# include <mln/io/pgm/save.hh> + +# include <mln/value/rgb8.hh> +# include <mln/io/ppm/load.hh> +# include <mln/io/ppm/save.hh> + +# include <mln/accu/min_max.hh> +# include <mln/accu/mean.hh> + +# include <mln/fun/i2v/array.hh> +# include <mln/fun/p2v/iota.hh> + +# include <mln/level/paste.hh> +# include <mln/level/fill.hh> +# include <mln/level/transform.hh> +# include <mln/extension/fill.hh> +# include <mln/convert/to.hh> + +# include <mln/linear/gaussian.hh> + +# include <mln/morpho/meyer_wst.hh> +# include <mln/morpho/closing_volume.hh> + +# include <mln/make/w_window2d.hh> + +# include <mln/debug/println.hh> + +// Laplacian method +namespace mln +{ + namespace linear + { + // required to deal with a input image that differ from the output since I + //don't succeed in using a float image in entry. + template <class I, class O> + inline + void + gaussian_2nd_derivative(const Image<I>& input, float sigma, Image<O>& output) + { + mln_precondition(exact(input).has_data()); + + impl::recursivefilter_coef_ + coef(-1.331f, 3.661f, + 1.24f, 1.314f, + 0.3225f, -1.738f, + 0.748f, 2.166f, + sigma, impl::gaussian_2nd_deriv_coef_norm_); + impl::generic_filter_common_(mln_trait_value_nature(mln_value(I))(), + input, coef, sigma, output); + } + + template <class I, class O> + inline + void + gaussian_1st_derivative(const Image<I>& input, float sigma, Image<O>& output) + { + mln_precondition(exact(input).has_data()); + + impl::recursivefilter_coef_ + coef(-0.6472f, -4.531f, + 1.527f, 1.516f, + 0.6494f, 0.9557f, + 0.6719f, 2.072f, + sigma, impl::gaussian_1st_deriv_coef_norm_); + impl::generic_filter_common_(mln_trait_value_nature(mln_value(I))(), + input, coef, sigma, output); + } + } +} + +// Gradient + watershed method +namespace mln +{ + namespace morpho + { + template <typename I, typename N> + mln_concrete(I) + closing_volume(const I& input, const Neighborhood<N>& nbh, std::size_t lambda) + { + mln_concrete(I) output; + initialize(output, input); + closing_volume(input, nbh, lambda, output); + return output; + } + } +} // !mln + +namespace mln +{ + image2d<value::int_u8> convert_to_grey(const image2d<value::rgb8>& data) + { + image2d<value::int_u8> output(data.domain()); + mln_piter_(image2d<value::int_u8>) p(output.domain()); + for_all(p) + output(p) = (int) (data(p).red() * 0.3 + data(p).green() * 0.58 + data(p).blue()) * 0.12; + return output; + } +} // !mln + +// Functions + +inline +bool is_row_odd(const mln::point2d& p) +{ + return p.row() % 2; +} + +inline +bool is_cell(const mln::point2d& p) +{ + return p.row() % 2 == 0 && p.col() % 2 == 0; +} + +inline +bool is_edge(const mln::point2d& p) +{ + return p.row() % 2 + p.col() % 2 == 1; +} + +inline +bool is_point(const mln::point2d& p) +{ + return p.row() % 2 && p.col() % 2; +} + +inline +bool is_not_edge(const mln::point2d& p) +{ + return ! is_edge(p); +} + + + +namespace mln +{ + + namespace border + { + + template <typename I> + void + fill(I& ima, const mln_value(I)& v) + { + const int nrows = ima.nrows(); + const int ncols = ima.ncols(); + for (int r = -1; r <= nrows; ++r) + { + ima.at(r, -1) = v; + ima.at(r, ncols) = v; + } + for (int c = -1; c <= ncols; ++c) + { + ima.at(-1, c) = v; + ima.at(nrows, c) = v; + } + } + + } // mln::border + + namespace accu + { + + template <typename I, typename L, typename A, typename V> + inline + void + compute(const Image<I>& input_, + const Image<L>& label_, + const Accumulator<A>&, + V& v) + { + trace::entering("accu::compute"); + + const I& input = exact(input_); + const L& label = exact(label_); + + const unsigned n = v.size(); + std::vector<A> a(n); + + mln_piter(I) p(input.domain()); + for_all(p) + a[label(p)].take(input(p)); + + for (unsigned l = 1; l < n; ++l) + v(l) = a[l].to_result(); + + trace::exiting("accu::compute"); + } + + } // mln::accu + + namespace morpho + { + + template <typename I, typename N> + mln_concrete(I) + gradient(const I& input, const N& nbh) + { + mln_concrete(I) output; + initialize(output, input); + accu::min_max<mln_value(I)> mm; + + mln_piter(I) p(input.domain()); + mln_niter(N) n(nbh, p); + for_all(p) + { + mm.init(); + for_all(n) if (input.has(n)) + mm.take(input(n)); + output(p) = mm.second() - mm.first(); + } + return output; + } + + template <typename I, typename N> + mln_concrete(I) + dilation(const I& input, const N& nbh) + { + typedef mln_value(I) V; + + mln_concrete(I) output; + initialize(output, input); + + mln_piter(I) p(input.domain()); + mln_niter(N) n(nbh, p); + for_all(p) + { + for_all(n) + if (input.has(n) && input(n) != value::rgb8(0,0,0)) + output(p) = input(n); + } + return output; + } + } // mln::morpho + + + struct colorize : Function_v2v< colorize > + { + typedef value::rgb8 result; + colorize(unsigned max) + : lut(max + 1) + { + lut[0] = literal::black; + for (unsigned i = 1; i <= max; ++i) + lut[i] = result(100 + std::rand() % 150, + 100 + std::rand() % 150, + 100 + std::rand() % 150); + } + result operator()(unsigned i) const + { + return lut[i]; + } + std::vector<result> lut; + }; + + + template <typename I> + I display_edge_with_bg(const I& ima, unsigned zoom, mln_value(I) bg) + { + unsigned nrows = ima.nrows() / 2 + 1; + unsigned ncols = ima.ncols() / 2 + 1; + I output(nrows * (zoom + 1) - 1, + ncols * (zoom + 1) - 1); + level::fill(output, bg); + mln_VAR( edge, ima | is_edge ); + mln_piter(edge_t) p(edge.domain()); + for_all(p) + if (p.row() % 2) // horizontal edge + { + unsigned row = (p.row() / 2 + 1) * (zoom + 1) - 1; + unsigned col = (p.col() / 2) * (zoom + 1); + for (unsigned i = 0; i < zoom; ++i) + output.at(row, col + i) = ima(p); + } + else // vertical edge + { + unsigned row = (p.row() / 2) * (zoom + 1); + unsigned col = (p.col() / 2 + 1) * (zoom + 1) - 1; + for (unsigned i = 0; i < zoom; ++i) + output.at(row + i, col) = ima(p); + } + return output; + } + + template <typename I> + I display_edge(const I& ima, unsigned zoom) + { + unsigned nrows = ima.nrows() / 2 + 1; + unsigned ncols = ima.ncols() / 2 + 1; + I output(nrows * (zoom + 1) - 1, + ncols * (zoom + 1) - 1); + + mln_VAR( cell, ima | is_cell ); + mln_piter(cell_t) q(cell.domain()); + for_all(q) + { + unsigned row = (q.row() / 2) * (zoom + 1); + unsigned col = (q.col() / 2) * (zoom + 1); + for (unsigned i = 0; i < zoom; ++i) + for (unsigned j = 0; j < zoom; ++j) + output.at(row + i, col + j) = ima(q); + } + + mln_VAR( edge, ima | is_edge ); + mln_piter(edge_t) p(edge.domain()); + for_all(p) + if (p.row() % 2) // horizontal edge + { + unsigned row = (p.row() / 2 + 1) * (zoom + 1) - 1; + unsigned col = (p.col() / 2) * (zoom + 1); + for (unsigned i = 0; i < zoom; ++i) + output.at(row, col + i) = ima(p); + } + else // vertical edge + { + unsigned row = (p.row() / 2) * (zoom + 1); + unsigned col = (p.col() / 2 + 1) * (zoom + 1) - 1; + for (unsigned i = 0; i < zoom; ++i) + output.at(row + i, col) = ima(p); + } + return output; + } + +} // mln + + + +template <typename T> +mln::image2d<T> +image2cells(const mln::image2d<T>& input) +{ + mln::image2d<T> output(2 * input.nrows() - 1, + 2 * input.ncols() - 1); + for (int row = 0; row < input.nrows(); ++row) + for (int col = 0; col < input.ncols(); ++col) + output.at(2 * row, 2 * col) = input.at(row, col); + return output; +} + +namespace mln { + + template <typename I, typename N, typename M> + mln_concrete(I) + mean(const I& input, const N& nbh, const M& nbh2) + { + mln_concrete(I) output; + initialize(output, input); + + mln_piter(I) p(input.domain()); + mln_niter(N) n(nbh, p); + mln_niter(M) m(nbh2, p); + for_all(p) + { + if (is_edge(p)) + { + int nb = 0; + int r = 0, g = 0, b = 0; + for_all(n) + { + if (input.has(n)) + { + r += input(n).red(); + g += input(n).green(); + b += input(n).blue(); + ++nb; + } + } + output(p) = value::rgb8(r / nb, g / nb, b / nb); + } + if (is_point(p)) + { + int nb = 0; + int r = 0, g = 0, b = 0; + for_all(m) + { + if (input.has(m)) + { + r += input(m).red(); + g += input(m).green(); + b += input(m).blue(); + ++nb; + } + } + output(p) = value::rgb8(r / nb, g / nb, b / nb); + } + if (is_cell(p)) + output(p) = input(p); + } + return output; + } + +} + +template <typename T> +mln::image2d<T> +cells2image(const mln::image2d<T>& input) +{ + mln::image2d<T> output((input.nrows() + 1) / 2, + (input.ncols() + 1) / 2); + for (int row = 0; row < input.nrows(); row += 2) + for (int col = 0; col < input.ncols(); col += 2) + output.at(row / 2, col / 2) = input.at(row, col); + return output; +} + + + + +template <typename I> +mln_concrete(I) +do_it(I& input, float lambda, unsigned& nbasins) +{ + using namespace mln; + + /**************************/ + /* Neighborhood defintion */ + /**************************/ + + // e2c + bool e2c_h[] = { 0, 1, 0, + 0, 0, 0, + 0, 1, 0 }; + bool e2c_v[] = { 0, 0, 0, + 1, 0, 1, + 0, 0, 0 }; + mln_VAR( e2c, make::double_neighb2d(is_row_odd, e2c_h, e2c_v) ); + + // e2e + bool e2e_h[] = { 0, 0, 1, 0, 0, + 0, 1, 0, 1, 0, + 0, 0, 0, 0, 0, + 0, 1, 0, 1, 0, + 0, 0, 1, 0, 0 }; + bool e2e_v[] = { 0, 0, 0, 0, 0, + 0, 1, 0, 1, 0, + 1, 0, 0, 0, 1, + 0, 1, 0, 1, 0, + 0, 0, 0, 0, 0 }; + mln_VAR( e2e, make::double_neighb2d(is_row_odd, e2e_h, e2e_v) ); + + // e2p + bool e2p_h[] = { 1, 0, 1, + 0, 0, 0, + 1, 0, 1 }; + bool e2p_v[] = { 1, 0, 1, + 0, 0, 0, + 1, 0, 1 }; + mln_VAR( e2p, make::double_neighb2d(is_row_odd, e2p_h, e2p_v) ); + + /******************/ + /* Initialisation */ + /******************/ + + I output = mean(image2cells(input), e2c, e2p); + io::ppm::save(output, "tmp_input.ppm"); + //image2d<value::int_u8> ima = convert_to_grey(output); + image2d<value::int_u8> imau = convert_to_grey(output); + io::pgm::save(imau, "tmp_grey_input.pgm"); + + image2d<float> ima(exact(imau).domain()); + + // cell + mln_VAR(cell, imau | is_cell); + + // edge + mln_VAR(edge, extend((imau | is_edge).rw(), pw::value(imau))); + + // FIXME until laplacian is working use gradient / closing_area / wst + + linear::gaussian_2nd_derivative(imau, lambda, ima); + + mln_piter(edge_t) p(edge.domain()); + for_all(p) + { + if (is_row_odd(p)) + { + mln_value(image2d<float>) t = ima.at(p.row() - 1, p.col()); + mln_value(image2d<float>) b = ima.at(p.row() + 1, p.col()); + if ((t > 0 && b < 0) || (t < 0 && b > 0)) + output(p) = value::rgb8(255,0,0); + } + else + { + mln_value(image2d<float>) r = ima.at(p.row(), p.col() - 1); + mln_value(image2d<float>) d = ima.at(p.row(), p.col() + 1); + if ((r > 0 && d < 0) || (r < 0 && d > 0)) + output(p) = value::rgb8(255,0,0); + } + } + +#if 0 + level::paste(morpho::gradient(edge, e2c), edge); + level::paste(morpho::closing_volume(edge, e2e, lambda), edge); + level::fill(edge, morpho::meyer_wst(edge, e2e, nbasins)); + + // Fill regions (with colorize) (won't work with laplacian...) + + colorize colors(nbasins); + + image2d<value::rgb8> cells(ima.domain()); + level::fill(cells, literal::white); + level::paste(level::transform(edge, colors), cells); + io::ppm::save(display_edge_with_bg(cells, 3, literal::white), "tmp_edge.ppm"); + + // Move the color of an edge which is non black in the cell + level::paste(morpho::dilation(cells, c4()), cells); +#endif + + //cells = convert_to_rgb8(ima); + + return output; +} + +void usage(char* argv[]) +{ + std::cerr << "usage: " << argv[0] << " input.pgm lambda output.ppm" << std::endl; + std::cerr << " lambda >= 0" << std::endl; + abort(); +} + +int main(int argc, char* argv[]) +{ + using namespace mln; + + if (argc != 4) + usage(argv); + + float lambda = atof(argv[2]); + if (lambda < 0) + usage(argv); + + image2d<value::rgb8> ima; + io::ppm::load(ima, argv[1]); + + unsigned nbasins; + image2d<value::rgb8> output = do_it(ima, lambda, nbasins); + + io::ppm::save(display_edge(output, 3), argv[3]); +} Index: ballas/color/min_tree_area_filter.cc --- ballas/color/min_tree_area_filter.cc (revision 0) +++ ballas/color/min_tree_area_filter.cc (revision 0) @@ -0,0 +1,505 @@ +# include <mln/core/var.hh> + +# include <mln/core/image/image2d.hh> +# include <mln/core/image/image_if.hh> +# include <mln/core/image/extended.hh> +# include <mln/core/routine/extend.hh> + +# include <mln/core/alias/window2d.hh> +# include <mln/core/alias/neighb2d.hh> +# include <mln/make/double_neighb2d.hh> +# include <mln/core/site_set/p_centered.hh> + +# include <mln/literal/origin.hh> +# include <mln/literal/black.hh> +# include <mln/literal/white.hh> + +# include <mln/value/int_u8.hh> +# include <mln/value/int_u16.hh> +# include <mln/io/pgm/load.hh> +# include <mln/io/pgm/save.hh> + +# include <mln/value/rgb8.hh> +# include <mln/io/ppm/load.hh> +# include <mln/io/ppm/save.hh> + +# include <mln/accu/min_max.hh> + +# include <mln/fun/i2v/array.hh> +# include <mln/fun/p2v/iota.hh> + +# include <mln/level/paste.hh> +# include <mln/level/fill.hh> +# include <mln/level/transform.hh> +# include <mln/extension/fill.hh> + +# include <mln/morpho/meyer_wst.hh> +# include <mln/morpho/closing_area.hh> + +# include <mln/debug/println.hh> + +# include "src/distance.hh" + +namespace mln +{ + template <typename I, typename N, typename Ic, typename Nc> + struct min_tree_ + { + typedef mln_site(I) point; + typedef p_array<point> S; + + // in: + const I& f; + const N& nbh; + const Ic& ref; + const Nc& nbhc; + + // aux: + S s; + mln_ch_value(I, bool) deja_vu; + mln_ch_value(I, point) parent; + mln_ch_value(I, bool) resp; + mln_ch_value(I, point) zpar; + + // attached data: + int lambda; + mln_ch_value(I, int) area; + //mln_ch_value(Ic, value::rgb8) values; + //initialize(values, ref); + //mln_ch_value(I, int) comp; + + min_tree_(const I& f, const N& nbh, const Ic& ref, const Nc& nbhc, + int lambda) + : f(f), + nbh(nbh), + ref(ref), + nbhc(nbhc), + lambda(lambda) + { + run(); + } + + void run() + { + // init + { + initialize(deja_vu, f); + initialize(parent, f); + initialize(resp, f); + initialize(zpar, f); + initialize(area, f); + //initialize(comp, f); + + mln::level::fill(deja_vu, false); + //mln::level::fill(resp, false); + mln::level::fill(area, 0); + + s = level::sort_psites_increasing(f); + } + + // first pass + { + mln_fwd_piter(S) p(s); + mln_niter(N) n(nbh, p); + for_all(p) + { + make_set(p); + for_all(n) + if (f.has(n) && deja_vu(n)) + do_union(n, p); + deja_vu(p) = true; + } + } + + // second pass: canonization + { + mln_bkd_piter(S) p(s); + for_all(p) + { + point q = parent(p); + if (f(parent(q)) == f(q)) + { + parent(p) = parent(q); + resp(q) = false; + } + } + } + + // third pass: Merging region with area < lambda + { + mln_fwd_piter(S) p(s); + for_all(p) + { + if (resp(p) && area(p) < lambda) + { + resp(p) = false; + update_data(parent(p), area(p)); + } + } + } + } // end of run() + + void make_set(const point& p) + { + parent(p) = p; + zpar(p) = p; + init_data(p); + } + + void set_parent(const point& r, const point& p) + { + parent(r) = p; + merge_data(r, p); + } + + bool is_root(const point& p) const + { + return parent(p) == p; + } + + bool is_node(const point& p) const + { + //return is_root(p) || f(parent(p)) != f(p); + return (is_root(p) || resp(p)); + } + + point find_root(const point& x) + { + if (zpar(x) == x) + return x; + else + return zpar(x) = find_root(zpar(x)); + } + + point find_representative(const point& x) + { + if (parent(x) == x || resp(x)) + return x; + else + return find_representative(parent(x)); + } + + void do_union(const point& n, const point& p) + { + point r = find_root(n); + if (r != p) + { + set_parent(r, p); + zpar(r) = p; + } + } + + void init_data(const point& p) + { + area(p) = 1; + resp(p) = true; + } + + void merge_data(const point& r, const point& p) + { + if (f(p) == f(r)) + { + resp(p) = false; + area(p) += area(r); + } + } + + void update_data(const point& p, int val) + { + area(p) += val; + if (parent(p) != p && !resp(p)) + update_data(parent(p), val); + } + + }; +} + +namespace mln +{ + image2d<value::int_u16> convert_to_grey(const image2d<value::rgb8>& data) + { + image2d<value::int_u16> output(data.domain()); + mln_piter_(image2d<value::int_u16>) p(output.domain()); + for_all(p) + output(p) = (int) (data(p).red() * 0.3 + data(p).green() * 0.58 + data(p).blue()) * 0.12; + return output; + } +} // end of mln + +namespace mln +{ + + struct colorize : Function_v2v< colorize > + { + typedef value::rgb8 result; + colorize(unsigned max) + : lut(max + 1) + { + lut[0] = literal::black; + for (unsigned i = 1; i <= max; ++i) + lut[i] = result(100 + std::rand() % 150, + 100 + std::rand() % 150, + 100 + std::rand() % 150); + } + result operator()(unsigned i) const + { + return lut[i]; + } + std::vector<result> lut; + }; + + template <typename I> + I display_edge(const I& ima, mln_value(I) bg, unsigned zoom) + { + unsigned nrows = ima.nrows() / 2 + 1; + unsigned ncols = ima.ncols() / 2 + 1; + I output(nrows * (zoom + 1) - 1, + ncols * (zoom + 1) - 1); + level::fill(output, bg); + + mln_VAR(edge, ima | is_edge); + mln_piter(edge_t) p(edge.domain()); + for_all(p) + if (p.row() % 2) // horizontal edge + { + unsigned row = (p.row() / 2 + 1) * (zoom + 1) - 1; + unsigned col = (p.col() / 2) * (zoom + 1); + for (unsigned i = 0; i < zoom; ++i) + output.at(row, col + i) = ima(p); + } + else // vertical edge + { + unsigned row = (p.row() / 2) * (zoom + 1); + unsigned col = (p.col() / 2 + 1) * (zoom + 1) - 1; + for (unsigned i = 0; i < zoom; ++i) + output.at(row + i, col) = ima(p); + } + return output; + } + + template <typename I> + I display_edge(const I& ima, unsigned zoom) + { + unsigned nrows = ima.nrows() / 2 + 1; + unsigned ncols = ima.ncols() / 2 + 1; + I output(nrows * (zoom + 1) - 1, + ncols * (zoom + 1) - 1); + + mln_VAR( cell, ima | is_cell ); + mln_piter(cell_t) q(cell.domain()); + for_all(q) + { + unsigned row = (q.row() / 2) * (zoom + 1); + unsigned col = (q.col() / 2) * (zoom + 1); + for (unsigned i = 0; i < zoom; ++i) + for (unsigned j = 0; j < zoom; ++j) + output.at(row + i, col + j) = ima(q); + } + + mln_VAR( edge, ima | is_edge ); + mln_piter(edge_t) p(edge.domain()); + for_all(p) + if (p.row() % 2) // horizontal edge + { + unsigned row = (p.row() / 2 + 1) * (zoom + 1) - 1; + unsigned col = (p.col() / 2) * (zoom + 1); + for (unsigned i = 0; i < zoom; ++i) + output.at(row, col + i) = ima(p); + } + else // vertical edge + { + unsigned row = (p.row() / 2) * (zoom + 1); + unsigned col = (p.col() / 2 + 1) * (zoom + 1) - 1; + for (unsigned i = 0; i < zoom; ++i) + output.at(row + i, col) = ima(p); + } + return output; + } + + + namespace morpho + { + + template <typename I, typename N> + mln_concrete(I) + dilation(const I& input, const N& nbh) + { + typedef mln_value(I) V; + + mln_concrete(I) output; + initialize(output, input); + + mln_piter(I) p(input.domain()); + mln_niter(N) n(nbh, p); + for_all(p) + { + for_all(n) + if (input.has(n) && input(n) != value::rgb8(0,0,0)) + output(p) = input(n); + } + return output; + } + } // mln::morpho + + +} // mln + + + +template <typename T> +mln::image2d<T> +image2cells(const mln::image2d<T>& input) +{ + mln::image2d<T> output(2 * input.nrows() - 1, + 2 * input.ncols() - 1); + for (unsigned row = 0; row < input.nrows(); ++row) + for (unsigned col = 0; col < input.ncols(); ++col) + output.at(2 * row, 2 * col) = input.at(row, col); + return output; +} + + +template <typename T> +mln::image2d<T> +cells2image(const mln::image2d<T>& input) +{ + mln::image2d<T> output((input.nrows() + 1) / 2, + (input.ncols() + 1) / 2); + for (unsigned row = 0; row < input.nrows(); row += 2) + for (unsigned col = 0; col < input.ncols(); col += 2) + output.at(row / 2, col / 2) = input.at(row, col); + return output; +} + + +template <typename I, typename N, typename Ic, typename Nc> +unsigned min_tree(const I& f, const N& nbh, const Ic& ref, const Nc& nbhc, + int lambda) +{ + using namespace mln; + + min_tree_<I,N,Ic,Nc> run(f, nbh, ref, nbhc, lambda); + + + mln_piter(I) p(f.domain()); + unsigned nnodes = 0; + for_all(p) + { + if (run.is_node(p)) + { + std::cout << "nodes: " << p << std::endl; + ++nnodes; + } + } + + colorize colors(nnodes); + image2d<value::rgb8> tmp(ref.domain()); + level::fill(tmp, ref); + + mln_piter(I) q(f.domain()); + unsigned int i = 0; + for_all(q) + { + if (run.is_node(q)) + { + tmp(q) = colors(i); + i++; + } + } + mln_piter(I) r(f.domain()); + for_all(r) + { + if (!run.is_node(r)) + { + tmp(r) = tmp(run.find_representative(r)); + } + } + + + image2d<value::rgb8> to_display(tmp.domain()); + + level::fill(to_display, value::rgb8(255, 255, 255)); + level::paste((tmp | is_edge), to_display); + level::paste(morpho::dilation(to_display, c4()), to_display); + + io::ppm::save(display_edge(tmp, literal::white, 3), + "edge.ppm"); + io::ppm::save(tmp, "full.ppm"); + io::ppm::save(cells2image(to_display), "colorize.ppm"); + + return nnodes; +} + + +template <typename I> +I +do_it(I& input, int lambda, unsigned& nbasins) +{ + using namespace mln; + + /// Graph creation + I graph; + create_graph(input, graph, value::rgb8(0, 0, 0)); + + // Initialization + image2d<value::int_u16> ima = convert_to_grey(graph); + + // Neigbhorhood + // e2c + bool e2c_h[] = { 0, 1, 0, + 0, 0, 0, + 0, 1, 0 }; + bool e2c_v[] = { 0, 0, 0, + 1, 0, 1, + 0, 0, 0 }; + + mln_VAR(e2c, make::double_neighb2d(is_row_odd, e2c_h, e2c_v)); + + bool e2e_h[] = { 0, 0, 1, 0, 0, + 0, 1, 0, 1, 0, + 0, 0, 0, 0, 0, + 0, 1, 0, 1, 0, + 0, 0, 1, 0, 0 }; + + bool e2e_v[] = { 0, 0, 0, 0, 0, + 0, 1, 0, 1, 0, + 1, 0, 0, 0, 1, + 0, 1, 0, 1, 0, + 0, 0, 0, 0, 0 }; + mln_VAR(e2e, make::double_neighb2d(is_row_odd, e2e_h, e2e_v)); + + // Algorithm + distance(extend((graph | is_edge).rw(), pw::value(graph)), e2c, ima); + io::pgm::save(ima, "edge.pgm"); + + nbasins = min_tree((ima | is_edge), e2e, graph, e2c, lambda); + + return graph; +} + +void usage(char* argv[]) +{ + std::cerr << "usage: " << argv[0] << " input.pgm lambda" << std::endl; + std::cerr << " lambda >= 0" << std::endl; + abort(); +} + +int main(int argc, char* argv[]) +{ + using namespace mln; + + if (argc != 3) + usage(argv); + + int lambda = atoi(argv[2]); + if (lambda < 0) + usage(argv); + + image2d<value::rgb8> ima; + io::ppm::load(ima, argv[1]); + + unsigned nbasins; + image2d<value::rgb8> output = do_it(ima, lambda, nbasins); + + //io::ppm::save(output, argv[3]); +} Index: ballas/color/min_tree_color.cc --- ballas/color/min_tree_color.cc (revision 0) +++ ballas/color/min_tree_color.cc (revision 0) @@ -0,0 +1,525 @@ +# include <mln/core/var.hh> + +# include <mln/core/image/image2d.hh> +# include <mln/core/image/image_if.hh> +# include <mln/core/image/extended.hh> +# include <mln/core/routine/extend.hh> + +# include <mln/core/alias/window2d.hh> +# include <mln/core/alias/neighb2d.hh> +# include <mln/make/double_neighb2d.hh> +# include <mln/core/site_set/p_centered.hh> + +# include <mln/literal/origin.hh> +# include <mln/literal/black.hh> +# include <mln/literal/white.hh> + +# include <mln/value/int_u8.hh> +# include <mln/value/int_u16.hh> +# include <mln/io/pgm/load.hh> +# include <mln/io/pgm/save.hh> + +# include <mln/value/rgb8.hh> +# include <mln/io/ppm/load.hh> +# include <mln/io/ppm/save.hh> + +# include <mln/accu/min_max.hh> + +# include <mln/fun/i2v/array.hh> +# include <mln/fun/p2v/iota.hh> + +# include <mln/level/paste.hh> +# include <mln/level/fill.hh> +# include <mln/level/transform.hh> +# include <mln/extension/fill.hh> + +# include <mln/morpho/closing_area.hh> + + +# include <mln/debug/println.hh> + +# include "src/distance.hh" + +namespace mln +{ + template <typename I, typename N, typename Ic, typename Nc> + struct min_tree_ + { + typedef mln_site(I) point; + typedef p_array<point> S; + + // in: + const I& f; + const N& nbh; + const Ic& ref; + const Nc& nbhc; + + // aux: + S s; + mln_ch_value(I, bool) deja_vu; + mln_ch_value(I, point) parent; + mln_ch_value(I, bool) resp; + mln_ch_value(I, point) zpar; + + // attached data: + unsigned lambda; + mln_ch_value(I, value::rgb8) color; + //mln_ch_value(Ic, value::rgb8) values; + //initialize(values, ref); + //mln_ch_value(I, int) comp; + + min_tree_(const I& f, const N& nbh, const Ic& ref, const Nc& nbhc, + int lambda) + : f(f), + nbh(nbh), + ref(ref), + nbhc(nbhc), + lambda(lambda) + { + run(); + } + + void run() + { + // init + { + initialize(deja_vu, f); + initialize(parent, f); + initialize(resp, f); + initialize(zpar, f); + initialize(color, f); + //initialize(comp, f); + + mln::level::fill(deja_vu, false); + //mln::level::fill(resp, false); + mln::level::fill(color, value::rgb8(0, 0, 0)); + + s = level::sort_psites_increasing(f); + } + + // first pass + { + mln_fwd_piter(S) p(s); + mln_niter(N) n(nbh, p); + for_all(p) + { + make_set(p); + for_all(n) + if (f.has(n) && deja_vu(n)) + do_union(n, p); + deja_vu(p) = true; + } + } + + // second pass: canonization + { + mln_bkd_piter(S) p(s); + for_all(p) + { + point q = parent(p); + if (f(parent(q)) == f(q)) + { + parent(p) = parent(q); + resp(q) = false; + } + } + } + + // third pass: Merging region with distance(color) < lambda + { + mln_fwd_piter(S) p(s); + for_all(p) + { + point q = parent(p); + if (resp(p) && distance(color(p), color(q)) < lambda) + { + resp(p) = false; + update_data(q, color(p)); + } + } + } + + } // end of run() + + void make_set(const point& p) + { + parent(p) = p; + zpar(p) = p; + init_data(p); + } + + void set_parent(const point& r, const point& p) + { + parent(r) = p; + merge_data(r, p); + } + + bool is_root(const point& p) const + { + return parent(p) == p; + } + + bool is_node(const point& p) const + { + //return is_root(p) || f(parent(p)) != f(p); + return (is_root(p) || resp(p)); + } + + point find_root(const point& x) + { + if (zpar(x) == x) + return x; + else + return zpar(x) = find_root(zpar(x)); + } + + point find_representative(const point& x) + { + if (parent(x) == x || resp(x)) + return x; + else + return find_representative(parent(x)); + } + + void do_union(const point& n, const point& p) + { + point r = find_root(n); + if (r != p) + { + set_parent(r, p); + zpar(r) = p; + } + } + + void init_data(const point& p) + { + int red =0, green = 0, blue = 0; + + mln_niter(Nc) n(nbhc, p); + for_all(n) + { + red += ref(n).red(); + green += ref(n).green(); + blue += ref(n).blue(); + } + + red /= 2; + green /= 2; + blue /= 2; + + color(p).red() = red; + color(p).green() = green; + color(p).blue() = blue; + + resp(p) = true; + } + + void merge_data(const point& r, const point& p) + { + if (f(p) == f(r)) + { + resp(p) = false; + color(r) = (color(r) + color(p)) / 2; + } + } + + void update_data(const point& p, value::rgb8 val) + { + color(p) = (color(p) + val) / 2; + if (parent(p) != p && !resp(p)) + update_data(parent(p), color(p)); + } + + }; +} + +namespace mln +{ + image2d<value::int_u16> convert_to_grey(const image2d<value::rgb8>& data) + { + image2d<value::int_u16> output(data.domain()); + mln_piter_(image2d<value::int_u16>) p(output.domain()); + for_all(p) + output(p) = (int) (data(p).red() * 0.3 + data(p).green() * 0.58 + data(p).blue()) * 0.12; + return output; + } +} // end of mln + +namespace mln +{ + + struct colorize : Function_v2v< colorize > + { + typedef value::rgb8 result; + colorize(unsigned max) + : lut(max + 1) + { + lut[0] = literal::black; + for (unsigned i = 1; i <= max; ++i) + lut[i] = result(100 + std::rand() % 150, + 100 + std::rand() % 150, + 100 + std::rand() % 150); + } + result operator()(unsigned i) const + { + return lut[i]; + } + std::vector<result> lut; + }; + + template <typename I> + I display_edge(const I& ima, mln_value(I) bg, unsigned zoom) + { + unsigned nrows = ima.nrows() / 2 + 1; + unsigned ncols = ima.ncols() / 2 + 1; + I output(nrows * (zoom + 1) - 1, + ncols * (zoom + 1) - 1); + level::fill(output, bg); + + mln_VAR(edge, ima | is_edge); + mln_piter(edge_t) p(edge.domain()); + for_all(p) + if (p.row() % 2) // horizontal edge + { + unsigned row = (p.row() / 2 + 1) * (zoom + 1) - 1; + unsigned col = (p.col() / 2) * (zoom + 1); + for (unsigned i = 0; i < zoom; ++i) + output.at(row, col + i) = ima(p); + } + else // vertical edge + { + unsigned row = (p.row() / 2) * (zoom + 1); + unsigned col = (p.col() / 2 + 1) * (zoom + 1) - 1; + for (unsigned i = 0; i < zoom; ++i) + output.at(row + i, col) = ima(p); + } + return output; + } + + template <typename I> + I display_edge(const I& ima, unsigned zoom) + { + unsigned nrows = ima.nrows() / 2 + 1; + unsigned ncols = ima.ncols() / 2 + 1; + I output(nrows * (zoom + 1) - 1, + ncols * (zoom + 1) - 1); + + mln_VAR( cell, ima | is_cell ); + mln_piter(cell_t) q(cell.domain()); + for_all(q) + { + unsigned row = (q.row() / 2) * (zoom + 1); + unsigned col = (q.col() / 2) * (zoom + 1); + for (unsigned i = 0; i < zoom; ++i) + for (unsigned j = 0; j < zoom; ++j) + output.at(row + i, col + j) = ima(q); + } + + mln_VAR( edge, ima | is_edge ); + mln_piter(edge_t) p(edge.domain()); + for_all(p) + if (p.row() % 2) // horizontal edge + { + unsigned row = (p.row() / 2 + 1) * (zoom + 1) - 1; + unsigned col = (p.col() / 2) * (zoom + 1); + for (unsigned i = 0; i < zoom; ++i) + output.at(row, col + i) = ima(p); + } + else // vertical edge + { + unsigned row = (p.row() / 2) * (zoom + 1); + unsigned col = (p.col() / 2 + 1) * (zoom + 1) - 1; + for (unsigned i = 0; i < zoom; ++i) + output.at(row + i, col) = ima(p); + } + return output; + } + + + namespace morpho + { + + template <typename I, typename N> + mln_concrete(I) + dilation(const I& input, const N& nbh) + { + typedef mln_value(I) V; + + mln_concrete(I) output; + initialize(output, input); + + mln_piter(I) p(input.domain()); + mln_niter(N) n(nbh, p); + for_all(p) + { + for_all(n) + if (input.has(n) && input(n) != value::rgb8(0,0,0)) + output(p) = input(n); + } + return output; + } + } // mln::morpho + +} // mln + + + +template <typename T> +mln::image2d<T> +image2cells(const mln::image2d<T>& input) +{ + mln::image2d<T> output(2 * input.nrows() - 1, + 2 * input.ncols() - 1); + for (unsigned row = 0; row < input.nrows(); ++row) + for (unsigned col = 0; col < input.ncols(); ++col) + output.at(2 * row, 2 * col) = input.at(row, col); + return output; +} + + +template <typename T> +mln::image2d<T> +cells2image(const mln::image2d<T>& input) +{ + mln::image2d<T> output((input.nrows() + 1) / 2, + (input.ncols() + 1) / 2); + for (unsigned row = 0; row < input.nrows(); row += 2) + for (unsigned col = 0; col < input.ncols(); col += 2) + output.at(row / 2, col / 2) = input.at(row, col); + return output; +} + + +template <typename I, typename N, typename Ic, typename Nc> +unsigned min_tree(const I& f, const N& nbh, const Ic& ref, const Nc& nbhc, + int lambda) +{ + using namespace mln; + + min_tree_<I,N,Ic,Nc> run(f, nbh, ref, nbhc, lambda); + + + mln_piter(I) p(f.domain()); + unsigned nnodes = 0; + for_all(p) + { + if (run.is_node(p)) + ++nnodes; + } + +#if 1 + colorize colors(nnodes); + image2d<value::rgb8> tmp(ref.domain()); + level::fill(tmp, ref); + + + mln_piter(I) q(f.domain()); + unsigned int i = 0; + for_all(q) + { + if (run.is_node(q)) + { + tmp(q) = colors(i); + i++; + } + } + mln_piter(I) r(f.domain()); + for_all(r) + { + if (!run.is_node(r)) + tmp(r) = tmp(run.find_representative(r)); + } + + image2d<value::rgb8> to_display(tmp.domain()); + + level::fill(to_display, value::rgb8(255, 255, 255)); + level::paste((tmp | is_edge), to_display); + level::paste(morpho::dilation(to_display, c4()), to_display); + + io::ppm::save(display_edge(tmp, literal::white, 3), + "edge.ppm"); + io::ppm::save(tmp, "full.ppm"); + io::ppm::save(cells2image(to_display), "colorize.ppm"); +#endif + +#if 1 + // io::ppm::save(display_edge(run.values, literal::white, 3), "tmp_tree_colored.pgm"); +#endif + + return nnodes; +} + + +template <typename I> +I +do_it(I& input, int lambda, unsigned& nbasins) +{ + using namespace mln; + + /// Graph creation + I graph; + create_graph(input, graph, value::rgb8(0, 0, 0)); + + // Initialization + image2d<value::int_u16> ima = convert_to_grey(graph); + + // Neigbhorhood + // e2c + bool e2c_h[] = { 0, 1, 0, + 0, 0, 0, + 0, 1, 0 }; + bool e2c_v[] = { 0, 0, 0, + 1, 0, 1, + 0, 0, 0 }; + + mln_VAR(e2c, make::double_neighb2d(is_row_odd, e2c_h, e2c_v)); + + bool e2e_h[] = { 0, 0, 1, 0, 0, + 0, 1, 0, 1, 0, + 0, 0, 0, 0, 0, + 0, 1, 0, 1, 0, + 0, 0, 1, 0, 0 }; + + bool e2e_v[] = { 0, 0, 0, 0, 0, + 0, 1, 0, 1, 0, + 1, 0, 0, 0, 1, + 0, 1, 0, 1, 0, + 0, 0, 0, 0, 0 }; + mln_VAR(e2e, make::double_neighb2d(is_row_odd, e2e_h, e2e_v)); + + // Algorithm + distance(extend((graph | is_edge).rw(), pw::value(graph)), e2c, ima); + + io::pgm::save(ima, "edge.pgm"); + + nbasins = min_tree((ima | is_edge), e2e, graph, e2c, lambda); + + return graph; +} + +void usage(char* argv[]) +{ + std::cerr << "usage: " << argv[0] << " input.pgm lambda" << std::endl; + std::cerr << " lambda >= 0" << std::endl; + abort(); +} + +int main(int argc, char* argv[]) +{ + using namespace mln; + + if (argc != 3) + usage(argv); + + int lambda = atoi(argv[2]); + if (lambda < 0) + usage(argv); + + image2d<value::rgb8> ima; + io::ppm::load(ima, argv[1]); + + unsigned nbasins; + image2d<value::rgb8> output = do_it(ima, lambda, nbasins); + + //io::ppm::save(output, argv[3]); +} Index: ballas/color/reference2.cc --- ballas/color/reference2.cc (revision 0) +++ ballas/color/reference2.cc (revision 0) @@ -0,0 +1,406 @@ +# include <mln/core/var.hh> + +# include <mln/core/image/image2d.hh> +# include <mln/core/image/image_if.hh> +# include <mln/core/image/extended.hh> +# include <mln/core/routine/extend.hh> + +# include <mln/core/alias/window2d.hh> +# include <mln/core/alias/neighb2d.hh> +# include <mln/make/double_neighb2d.hh> +# include <mln/core/site_set/p_centered.hh> + +# include <mln/literal/origin.hh> +# include <mln/literal/black.hh> +# include <mln/literal/white.hh> + +# include <mln/value/int_u8.hh> +# include <mln/io/pgm/load.hh> +# include <mln/io/pgm/save.hh> + +# include <mln/value/rgb8.hh> +# include <mln/io/ppm/load.hh> +# include <mln/io/ppm/save.hh> + +# include <mln/accu/min_max.hh> +# include <mln/accu/mean.hh> + +# include <mln/fun/i2v/array.hh> +# include <mln/fun/p2v/iota.hh> + +# include <mln/level/paste.hh> +# include <mln/level/fill.hh> +# include <mln/level/transform.hh> +# include <mln/extension/fill.hh> + +# include <mln/morpho/meyer_wst.hh> +# include <mln/morpho/closing_volume.hh> + +# include <mln/linear/convolve.hh> +# include <mln/make/w_window2d.hh> + +# include <mln/debug/println.hh> + +namespace mln +{ + namespace morpho + { + template <typename I, typename N> + mln_concrete(I) + closing_volume(const I& input, const Neighborhood<N>& nbh, std::size_t lambda) + { + mln_concrete(I) output; + initialize(output, input); + closing_volume(input, nbh, lambda, output); + return output; + } + } +} + +namespace mln +{ + template <typename I, typename O> + void + LoG_17x17(const I& input, const O& output) + { + mln_precondition(exact(output).domain() == exact(input).domain()); + int ws[] = { +0, 0, 0, 0, 0, 0,-1,-1,-1,-1,-1, 0, 0, 0, 0, 0, 0, + +0, 0, 0, 0,-1,-1,-1,-1,-1,-1,-1,-1,-1, 0, 0, 0, 0, + +0, 0,-1,-1,-1,-2,-3,-3,-3,-3,-3,-2,-1,-1,-1, 0, 0, + +0, 0,-1,-1,-2,-3,-3,-3,-3,-3,-3,-3,-2,-1,-1, 0, 0, + +0,-1,-1,-2,-3,-3,-3,-2,-3,-2,-3,-3,-3,-2,-1,-1, 0, + +0,-1,-2,-3,-3,-3, 0, 2, 4, 2, 0,-3,-3,-3,-2,-1, 0, + -1,-1,-3,-3,-3, 0, 4,10,12,10, 4, 0,-3,-3,-3,-1,-1, + -1,-1,-3,-3,-2, 2,10,18,21,18,10, 2,-2,-3,-3,-1,-1, + -1,-1,-3,-3,-3, 4,12,21,24,21,12, 4,-3,-3,-3,-1,-1, + -1,-1,-3,-3,-2, 2,10,18,21,18,10, 2,-2,-3,-3,-1,-1, + -1,-1,-3,-3,-3, 0, 4,10,12,10, 4, 0,-3,-3,-3,-1,-1, + +0,-1,-2,-3,-3,-3, 0, 2, 4, 2, 0,-3,-3,-3,-2,-1, 0, + +0,-1,-1,-2,-3,-3,-3,-2,-3,-2,-3,-3,-3,-2,-1,-1, 0, + +0, 0,-1,-1,-2,-3,-3,-3,-3,-3,-3,-3,-2,-1,-1, 0, 0, + +0, 0,-1,-1,-1,-2,-3,-3,-3,-3,-3,-2,-1,-1,-1, 0, 0, + +0, 0, 0, 0,-1,-1,-1,-1,-1,-1,-1,-1,-1, 0, 0, 0, 0, + +0, 0, 0, 0, 0, 0,-1,-1,-1,-1,-1, 0, 0, 0, 0, 0, 0 }; + linear::convolve(input, make::w_window2d(ws), output); + } +} // !mln + +namespace mln +{ + image2d<value::int_u8> convert_to_grey(const image2d<value::rgb8>& data) + { + image2d<value::int_u8> output(data.domain()); + mln_piter_(image2d<value::int_u8>) p(output.domain()); + for_all(p) + output(p) = (int) (data(p).red() * 0.3 + data(p).green() * 0.58 + data(p).blue()) * 0.12; + return output; + } +} // !mln + +// Functions + +inline +bool is_row_odd(const mln::point2d& p) +{ + return p.row() % 2; +} + +inline +bool is_cell(const mln::point2d& p) +{ + return p.row() % 2 == 0 && p.col() % 2 == 0; +} + +inline +bool is_edge(const mln::point2d& p) +{ + return p.row() % 2 + p.col() % 2 == 1; +} + +inline +bool is_point(const mln::point2d& p) +{ + return p.row() % 2 && p.col() % 2; +} + +inline +bool is_not_edge(const mln::point2d& p) +{ + return ! is_edge(p); +} + + + +namespace mln +{ + + namespace border + { + + template <typename I> + void + fill(I& ima, const mln_value(I)& v) + { + const int nrows = ima.nrows(); + const int ncols = ima.ncols(); + for (int r = -1; r <= nrows; ++r) + { + ima.at(r, -1) = v; + ima.at(r, ncols) = v; + } + for (int c = -1; c <= ncols; ++c) + { + ima.at(-1, c) = v; + ima.at(nrows, c) = v; + } + } + + } // mln::border + + namespace accu + { + + template <typename I, typename L, typename A, typename V> + inline + void + compute(const Image<I>& input_, + const Image<L>& label_, + const Accumulator<A>&, + V& v) + { + trace::entering("accu::compute"); + + const I& input = exact(input_); + const L& label = exact(label_); + + const unsigned n = v.size(); + std::vector<A> a(n); + + mln_piter(I) p(input.domain()); + for_all(p) + a[label(p)].take(input(p)); + + for (unsigned l = 1; l < n; ++l) + v(l) = a[l].to_result(); + + trace::exiting("accu::compute"); + } + + } // mln::accu + + namespace morpho + { + + template <typename I, typename N> + mln_concrete(I) + gradient(const I& input, const N& nbh) + { + mln_concrete(I) output; + initialize(output, input); + accu::min_max<mln_value(I)> mm; + + mln_piter(I) p(input.domain()); + mln_niter(N) n(nbh, p); + for_all(p) + { + mm.init(); + for_all(n) if (input.has(n)) + mm.take(input(n)); + output(p) = mm.second() - mm.first(); + } + return output; + } + + template <typename I, typename N> + mln_concrete(I) + dilation(const I& input, const N& nbh) + { + typedef mln_value(I) V; + + mln_concrete(I) output; + initialize(output, input); + + mln_piter(I) p(input.domain()); + mln_niter(N) n(nbh, p); + for_all(p) + { + for_all(n) + if (input.has(n) && input(n) != value::rgb8(0,0,0)) + output(p) = input(n); + } + return output; + } + } // mln::morpho + + + struct colorize : Function_v2v< colorize > + { + typedef value::rgb8 result; + colorize(unsigned max) + : lut(max + 1) + { + lut[0] = literal::black; + for (unsigned i = 1; i <= max; ++i) + lut[i] = result(100 + std::rand() % 150, + 100 + std::rand() % 150, + 100 + std::rand() % 150); + } + result operator()(unsigned i) const + { + return lut[i]; + } + std::vector<result> lut; + }; + + + template <typename I> + I display_edge(const I& ima, mln_value(I) bg, unsigned zoom) + { + unsigned nrows = ima.nrows() / 2 + 1; + unsigned ncols = ima.ncols() / 2 + 1; + I output(nrows * (zoom + 1) - 1, + ncols * (zoom + 1) - 1); + level::fill(output, bg); + mln_VAR( edge, ima | is_edge ); + mln_piter(edge_t) p(edge.domain()); + for_all(p) + if (p.row() % 2) // horizontal edge + { + unsigned row = (p.row() / 2 + 1) * (zoom + 1) - 1; + unsigned col = (p.col() / 2) * (zoom + 1); + for (unsigned i = 0; i < zoom; ++i) + output.at(row, col + i) = ima(p); + } + else // vertical edge + { + unsigned row = (p.row() / 2) * (zoom + 1); + unsigned col = (p.col() / 2 + 1) * (zoom + 1) - 1; + for (unsigned i = 0; i < zoom; ++i) + output.at(row + i, col) = ima(p); + } + return output; + } + +} // mln + + + +template <typename T> +mln::image2d<T> +image2cells(const mln::image2d<T>& input) +{ + mln::image2d<T> output(2 * input.nrows() - 1, + 2 * input.ncols() - 1); + for (int row = 0; row < input.nrows(); ++row) + for (int col = 0; col < input.ncols(); ++col) + output.at(2 * row, 2 * col) = input.at(row, col); + return output; +} + + +template <typename T> +mln::image2d<T> +cells2image(const mln::image2d<T>& input) +{ + mln::image2d<T> output((input.nrows() + 1) / 2, + (input.ncols() + 1) / 2); + for (int row = 0; row < input.nrows(); row += 2) + for (int col = 0; col < input.ncols(); col += 2) + output.at(row / 2, col / 2) = input.at(row, col); + return output; +} + + + + +template <typename I> +mln_concrete(I) +do_it(I& input, int lambda, unsigned& nbasins) +{ + using namespace mln; + + /******************/ + /* Initialisation */ + /******************/ + + image2d<value::int_u8> ima = image2cells(convert_to_grey(input)); + + /**************************/ + /* Neighborhood defintion */ + /**************************/ + + // e2c + bool e2c_h[] = { 0, 1, 0, + 0, 0, 0, + 0, 1, 0 }; + bool e2c_v[] = { 0, 0, 0, + 1, 0, 1, + 0, 0, 0 }; + mln_VAR( e2c, make::double_neighb2d(is_row_odd, e2c_h, e2c_v) ); + + // e2e + bool e2e_h[] = { 0, 0, 1, 0, 0, + 0, 1, 0, 1, 0, + 0, 0, 0, 0, 0, + 0, 1, 0, 1, 0, + 0, 0, 1, 0, 0 }; + bool e2e_v[] = { 0, 0, 0, 0, 0, + 0, 1, 0, 1, 0, + 1, 0, 0, 0, 1, + 0, 1, 0, 1, 0, + 0, 0, 0, 0, 0 }; + mln_VAR( e2e, make::double_neighb2d(is_row_odd, e2e_h, e2e_v) ); + + // cell + mln_VAR(cell, ima | is_cell); + + // edge + mln_VAR(edge, extend((ima | is_edge).rw(), pw::value(ima))); + + // FIXME until laplacian is working use gradient / closing_area / wst + + level::paste(morpho::gradient(edge, e2c), edge); + level::paste(morpho::closing_volume(edge, e2e, lambda), edge); + level::fill(edge, morpho::meyer_wst(edge, e2e, nbasins)); + + // Fill regions (with colorize) (won't work with laplacian...) + + colorize colors(nbasins); + + image2d<value::rgb8> cells(ima.domain()); + level::fill(cells, literal::white); + level::paste(level::transform(edge, colors), cells); + io::ppm::save(display_edge(cells, literal::white, 3), "tmp_edge.ppm"); + + // Move the color of an edge which is non black in the cell + level::paste(morpho::dilation(cells, c4()), cells); + + return cells2image(cells); +} + +void usage(char* argv[]) +{ + std::cerr << "usage: " << argv[0] << " input.pgm lambda output.ppm" << std::endl; + std::cerr << " lambda >= 0" << std::endl; + abort(); +} + +int main(int argc, char* argv[]) +{ + using namespace mln; + + if (argc != 4) + usage(argv); + + int lambda = atoi(argv[2]); + if (lambda < 0) + usage(argv); + + image2d<value::rgb8> ima; + io::ppm::load(ima, argv[1]); + + unsigned nbasins; + image2d<value::rgb8> output = do_it(ima, lambda, nbasins); + + io::ppm::save(output, argv[3]); +} Index: ballas/color/src/graph.hh --- ballas/color/src/graph.hh (revision 0) +++ ballas/color/src/graph.hh (revision 0) @@ -0,0 +1,61 @@ +/*! \file src/graph.hh + * + */ + +#ifndef SRC_GRAPH_HH +# define SRC_GRAPH_HH + +# include <mln/value/int_u8.hh> +# include <mln/value/int_u16.hh> +# include <mln/value/rgb8.hh> + +# include <mln/level/fill.hh> + +# include <mln/core/image/image2d.hh> + + +// Neighborhood functions +inline +bool is_row_odd(const mln::point2d& p) +{ + return p.row() % 2; +} +inline +bool is_cell(const mln::point2d& p) +{ + return p.row() % 2 == 0 && p.col() % 2 == 0; +} +inline +bool is_edge(const mln::point2d& p) +{ + return (p.row() % 2 + p.col() % 2) == 1; +} +inline +bool is_point(const mln::point2d& p) +{ + return p.row() % 2 && p.col() % 2; +} +inline +bool is_not_edge(const mln::point2d& p) +{ + return ! is_edge(p); +} + +// Graph image creation function +// FIXME: add exact conversion.... +// FIXME: check that the input image is in 2 dimension +template <typename I> +void +create_graph(const I& ima, I& graph, mln_value(I) val) +{ + graph = I(ima.nrows() * 2 -1, ima.ncols() * 2 - 1); + + mln::level::fill(graph, val); + + mln_piter(I) p(ima.domain()); + for_all(p) + graph.at(p.row() * 2, p.col() * 2) = ima(p); +} + + +#endif // !SRC_GRAPH_HH Index: ballas/color/src/io.hh --- ballas/color/src/io.hh (revision 0) +++ ballas/color/src/io.hh (revision 0) @@ -0,0 +1,57 @@ +/*! \file src/io.hh + * + * Contains various method to load/save an image + */ + +#ifndef SRC_IO_HH +# define SRC_IO_HH + +#include <mln/core/image/image2d.hh> + +#include <mln/value/int_u16.hh> +#include <mln/value/int_s16.hh> + +# include <mln/io/pgm/load.hh> +# include <mln/io/pgm/save.hh> +# include <mln/io/ppm/load.hh> +# include <mln/io/ppm/save.hh> + +namespace IO +{ + + template <typename I> + void load(I& ima, const std::string& file) + { + mln::io::ppm::load(ima, file); + } + + template <> + void load(mln::image2d<mln::value::int_u16>& ima, + const std::string& file) + { + mln::io::pgm::load(ima, file); + } + + template <typename I> + void save(I& ima, const std::string& file) + { + mln::io::ppm::save(ima, file); + } + + template <> + void save(mln::image2d<mln::value::int_u16>& ima, + const std::string& file) + { + mln::io::pgm::save(ima, file); + } + + template <> + void save(mln::image2d<mln::value::int_s16>& ima, + const std::string& file) + { + mln::io::pgm::save(ima, file); + } + +} // !IO + +#endif // !SRC_IO_HH Index: ballas/color/src/distance.hh --- ballas/color/src/distance.hh (revision 0) +++ ballas/color/src/distance.hh (revision 0) @@ -0,0 +1,52 @@ +/*! \file src/distance.hh + * + */ + +#ifndef SRC_DISTANCE_HH +# define SRC_DISTANCE_HH + +# include "graph.hh" + +# include <cmath> + + +/// Manhatan distance +inline +unsigned distance(const mln::value::rgb8& lhs, + const mln::value::rgb8& rhs) +{ + return abs(lhs.red() - rhs.red()) + + abs(lhs.green() - rhs.green()) + + abs(lhs.blue() - rhs.blue()); +} + + +/// Store the distance between two points on edge +/// FIXME documentation +template <typename I, typename N, typename O> +void distance(const I& rgb_graph, + const N& nbh, + O& gl_graph) +{ + mln_piter(I) p(rgb_graph.domain()); + mln_niter(N) n(nbh, p); + + for_all(p) + { + mln::value::rgb8 v1; + mln::value::rgb8 v2; + + n.start(); + assert(n.is_valid() && rgb_graph.has(n)); + v1 = rgb_graph(n); + n.next(); + assert(n.is_valid() && rgb_graph.has(n)); + v2 = rgb_graph(n); + + gl_graph(p) = distance(v1, v2); + } +} + + + +#endif // !SRC_DISTANCE_HH Index: ballas/color/src/convert.hh --- ballas/color/src/convert.hh (revision 0) +++ ballas/color/src/convert.hh (revision 0) @@ -0,0 +1,37 @@ +/*! \file src/convert.hh + * + * Method that convert an rgb image 2d into gray level + */ + +#ifndef SRC_CONVERT_HH +# define SRC_CONVERT_HH + +# include <mln/value/int_u8.hh> +# include <mln/value/int_u16.hh> +# include <mln/value/int_s16.hh> +# include <mln/value/rgb8.hh> + +# include <mln/core/image/image2d.hh> + +// Convert function +void convert_to_gl(mln::image2d<mln::value::int_u16>& ima, + const mln::image2d<mln::value::rgb8>& data) +{ + mln_piter_(mln::image2d<mln::value::int_u16>) p(ima.domain()); + for_all(p) + ima(p) = (int) (data(p).red() * 0.3 + + data(p).green() * 0.58 + + data(p).blue() * 0.12); +} + +void convert_to_gl(mln::image2d<mln::value::int_s16>& ima, + const mln::image2d<mln::value::rgb8>& data) +{ + mln_piter_(mln::image2d<mln::value::int_s16>) p(ima.domain()); + for_all(p) + ima(p) = (int) (data(p).red() * 0.3 + + data(p).green() * 0.58 + + data(p).blue() * 0.12); +} + +#endif // !SRC_CONVERT_HH Index: ballas/color/laplacien.cc --- ballas/color/laplacien.cc (revision 0) +++ ballas/color/laplacien.cc (revision 0) @@ -0,0 +1,129 @@ +#include <mln/core/image/image2d.hh> +#include <mln/core/image/image_if.hh> + +# include <mln/core/alias/neighb2d.hh> +# include <mln/make/double_neighb2d.hh> + +#include <mln/value/rgb8.hh> +#include <mln/value/int_u16.hh> +#include <mln/value/int_s16.hh> + +#include <mln/linear/gaussian.hh> +#include <mln/morpho/erosion.hh> + +#include <mln/core/var.hh> +#include <mln/debug/println.hh> + +#include "src/io.hh" +#include "src/graph.hh" +#include "src/convert.hh" + + +void usage(char* argv[]) +{ + std::cerr << "usage: " << argv[0] << " input.ppm output.pgm" << std::endl; + abort(); +} + +/// FIXME Put these elsewhere +typedef mln::image2d<mln::value::rgb8> Icolor; +typedef mln::image2d<mln::value::int_u16> Igray; + +void process(Icolor& graph, const Icolor& input) +{ + using namespace mln; + + // Neighborhood definition + +// bool e2c_h[] = { 0, 1, 0, +// 0, 0, 0, +// 0, 1, 0 }; +// bool e2c_v[] = { 0, 0, 0, +// 1, 0, 1, +// 0, 0, 0 }; + +// bool e2e_h[] = { 0, 0, 1, 0, 0, +// 0, 1, 0, 1, 0, +// 0, 0, 0, 0, 0, +// 0, 1, 0, 1, 0, +// 0, 0, 1, 0, 0 }; +// bool e2e_v[] = { 0, 0, 0, 0, 0, +// 0, 1, 0, 1, 0, +// 1, 0, 0, 0, 1, +// 0, 1, 0, 1, 0, +// 0, 0, 0, 0, 0 }; + +// bool e2p_h[] = { 0, 0, 1, +// 0, 0, 0, +// 1, 0, 0 }; +// bool e2p_v[] = { 1, 0, 0, +// 0, 0, 0, +// 0, 0, 1 }; + //mln_VAR(e2c, make::double_neighb2d(is_row_odd, e2c_h, e2c_v)); + + + // convert the image into grey level + Igray gray_graph(graph.bbox()); + convert_to_gl(gray_graph, graph); + + Igray gray_input(input.bbox()); + convert_to_gl(gray_input, input); + + + //level::paste(morpho::dilation((gray_graph | is_cell), e2c.win()), gray_graph); +// Igray save(gray_graph.bbox()); +// level::paste(gray_graph | is_cell, save); +// IO::save(save, "tmp3.ppm"); + + // Create the laplacian image + //Igray laplacian(gray_input.bbox()); + image2d<int> laplacian(gray_input.bbox()); + linear::gaussian(gray_input, 2.0f, laplacian); + //IO::save(laplacian, "tmp1.ppm"); + linear::laplacian(gray_input, 2.0f, laplacian); + debug::println(laplacian); + //IO::save(laplacian, "tmp2.ppm"); + + // Display edge on the output images +// debug::println(laplacian); + +// mln_piter_(Igray) p(laplacian.domain()); +// for_all(p) +// if (laplacian(p) == 0u) +// { +// std::cout << "I'm here" << std::endl; +// graph.at(p.row() * 2, p.col() * 2) = value::rgb<8>(255, 0, 0); +// } + +} + + + +int main(int argc, char* argv[]) +{ + using namespace mln; + + // Initialisation + if (argc != 3) + usage(argv); + + std::string input_file(argv[1]); + std::string output_file(argv[2]); + + + Icolor input; + + // Load the image + IO::load(input, input_file); + + // create a graph image from the input + Icolor graph; + create_graph(input, graph); + + // Process + process(graph, input); + + + // Save + IO::save(graph, output_file); +}