00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00048 #include <map>
00049 #include <vector>
00050
00051 #include <mln/util/ord.hh>
00052
00053 #include <mln/core/image/image2d.hh>
00054 #include <mln/core/alias/point2d.hh>
00055 #include <mln/core/alias/window2d.hh>
00056 #include <mln/core/alias/neighb2d.hh>
00057
00058 #include <mln/convert/to_window.hh>
00059
00061 #include <mln/core/image/edge_image.hh>
00062 #include <mln/core/var.hh>
00063 #include <mln/fun/i2v/array.hh>
00064 #include <mln/util/graph.hh>
00065
00066 #include <mln/morpho/gradient.hh>
00067 #include <mln/morpho/closing/area.hh>
00068 #include <mln/morpho/meyer_wst.hh>
00069
00070 #include <mln/value/int_u8.hh>
00071 #include <mln/value/rgb8.hh>
00072 #include <mln/literal/black.hh>
00073 #include <mln/literal/colors.hh>
00074
00075 #include <mln/io/pgm/load.hh>
00076 #include <mln/io/ppm/save.hh>
00077
00078 #include <mln/math/max.hh>
00079 #include <mln/math/abs.hh>
00080
00081 #include <mln/util/site_pair.hh>
00082
00083 #include "tests/data.hh"
00084
00085
00086
00087 int main()
00088 {
00089 using namespace mln;
00090 using value::int_u8;
00091 using value::rgb8;
00092
00093
00094
00095
00096
00097 typedef int_u8 input_val_t;
00098 image2d<input_val_t> input;
00099 io::pgm::load(input, MLN_IMG_DIR "/tiny.pgm");
00100
00101
00102
00103 image2d<input_val_t> gradient =
00104 morpho::gradient (input, convert::to_window(c4()));
00105
00106
00107 image2d<input_val_t> work(input.domain());
00108 work = morpho::closing::area(gradient, c4(), 10);
00109
00110
00111
00112
00113
00114
00115
00116 util::graph g;
00117
00118
00119 image2d<unsigned> equiv_vertex;
00120 initialize(equiv_vertex, work);
00121
00122
00123 mln_fwd_piter_(image2d<input_val_t>) p(work.domain());
00124 for_all(p)
00125 equiv_vertex(p) = g.add_vertex();
00126
00127
00128 window2d next_c4_win;
00129 next_c4_win.insert(0, 1).insert(1, 0);
00130 typedef fun::i2v::array<int> edge_values_t;
00131 typedef fun::i2v::array< util::site_pair<point2d> > edge_sites_t;
00132 edge_values_t edge_values;
00133 edge_sites_t edge_sites;
00134 mln_fwd_qiter_(window2d) q(next_c4_win, p);
00135 for_all (p)
00136 for_all(q)
00137 if (work.domain().has(q))
00138 {
00139 g.add_edge(equiv_vertex(p), equiv_vertex(q));
00140 edge_values.append(math::max(work(p), work(q)));
00141 edge_sites.append(util::site_pair<point2d>(p, q));
00142 }
00143
00144
00145 typedef edge_image<util::site_pair<point2d>,int,util::graph> lg_ima_t;
00146 lg_ima_t lg_ima(g, edge_sites, edge_values);
00147
00148
00149
00150
00151
00152 typedef lg_ima_t::nbh_t nbh_t;
00153 nbh_t nbh;
00154
00155
00156 int_u8 nbasins;
00157 typedef edge_image<util::site_pair<point2d>,int_u8,util::graph> wshed_t;
00158 wshed_t wshed = morpho::meyer_wst(lg_ima, nbh, nbasins);
00159 mln_assertion(nbasins == 5u);
00160
00161
00162
00163
00164
00165
00166
00167
00168
00169
00170 typedef rgb8 output_val_t;
00171 typedef image2d<output_val_t> output_t;
00172 point2d output_pmin = input.domain().pmin();
00173 point2d output_pmax(input.domain().pmax()[0] * 2,
00174 input.domain().pmax()[1] * 2);
00175 output_t output(box2d(output_pmin, output_pmax));
00176 data::fill(output, literal::black);
00177
00178 for_all(p)
00179 {
00180
00181 point2d q(p[0] * 2, p[1] * 2);
00182 input_val_t v = input(p);
00183
00184
00185 output(q) = output_val_t(v, v, v);
00186 }
00187
00188 mln_piter_(output_t) p_out(output.domain());
00189 for_all(p_out)
00190 {
00191
00192 if (p_out[0] % 2 == 0 && p_out[1] % 2 == 1)
00193 output(p_out) = (output(p_out + left) + output(p_out + right)) / 2;
00194
00195 if (p_out[0] % 2 == 1 && p_out[1] % 2 == 0)
00196 output(p_out) = (output(p_out + up) + output(p_out + down)) / 2;
00197
00198 if (p_out[0] % 2 == 1 && p_out[1] % 2 == 1)
00199 output(p_out) =
00200 (output(p_out + dpoint2d(-1, -1)) +
00201 output(p_out + dpoint2d(-1, +1)) +
00202 output(p_out + dpoint2d(+1, -1)) +
00203 output(p_out + dpoint2d(+1, +1))) / 4;
00204 }
00205
00206
00207
00208
00209 mln_piter_(wshed_t) pw(wshed.domain());
00210 for_all(pw)
00211 {
00212 if (wshed(pw) == 0u)
00213 {
00214 mln_psite_(lg_ima_t) pp(pw);
00215
00216 int row1 = pp.first()[0] * 2;
00217 int col1 = pp.first()[1] * 2;
00218 int row2 = pp.second()[0] * 2;
00219 int col2 = pp.second()[1] * 2;
00220 point2d q((row1 + row2) / 2, (col1 + col2) / 2);
00221
00222 output(q) = literal::red;
00223 }
00224 }
00225
00226
00227
00228
00229
00230
00231
00232
00233
00234 for_all (p_out)
00235
00236 if (p_out[0] % 2 == 1 && p_out[1] % 2 == 1)
00237 {
00238
00239
00240
00241
00242 unsigned nwsheds =
00243 (output.has(p_out + up ) && output(p_out + up ) == literal::red) +
00244 (output.has(p_out + down ) && output(p_out + down ) == literal::red) +
00245 (output.has(p_out + left ) && output(p_out + right) == literal::red) +
00246 (output.has(p_out + right) && output(p_out + left ) == literal::red);
00247 if (nwsheds >= 2)
00248 output(p_out) = literal::red;
00249 }
00250 io::ppm::save(output, "lena_line_graph_image_wst1-out.ppm");
00251 }