Milena (Olena)
User documentation 2.0a Id
|
00001 // Copyright (C) 2008, 2009 EPITA Research and Development Laboratory (LRDE) 00002 // 00003 // This file is part of Olena. 00004 // 00005 // Olena is free software: you can redistribute it and/or modify it under 00006 // the terms of the GNU General Public License as published by the Free 00007 // Software Foundation, version 2 of the License. 00008 // 00009 // Olena is distributed in the hope that it will be useful, 00010 // but WITHOUT ANY WARRANTY; without even the implied warranty of 00011 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 00012 // General Public License for more details. 00013 // 00014 // You should have received a copy of the GNU General Public License 00015 // along with Olena. If not, see <http://www.gnu.org/licenses/>. 00016 // 00017 // As a special exception, you may use this file as part of a free 00018 // software project without restriction. Specifically, if other files 00019 // instantiate templates or use macros or inline functions from this 00020 // file, or you compile this file and link it with other files to produce 00021 // an executable, this file does not by itself cause the resulting 00022 // executable to be covered by the GNU General Public License. This 00023 // exception does not however invalidate any other reasons why the 00024 // executable file might be covered by the GNU General Public License. 00025 00030 00047 #include <cstdlib> 00048 00049 #include <vector> 00050 00051 #include <mln/value/int_u8.hh> 00052 #include <mln/value/int_u16.hh> 00053 00054 #include <mln/core/routine/duplicate.hh> 00055 00056 #include <mln/core/image/image2d.hh> 00057 #include <mln/core/var.hh> 00058 00059 #include <mln/morpho/line_gradient.hh> 00060 #include <mln/morpho/closing/area_on_vertices.hh> 00061 #include <mln/labeling/regional_minima.hh> 00062 #include <mln/morpho/meyer_wst.hh> 00063 00064 #include <mln/io/pgm/load.hh> 00065 #include <mln/io/pgm/save.hh> 00066 00067 int main(int argc, char* argv[]) 00068 { 00069 if (argc != 4) 00070 { 00071 std::cerr 00072 << "usage: " << argv[0] << " max_nregions input.pgm output.pgm" 00073 << std::endl; 00074 std::exit(EXIT_FAILURE); 00075 } 00076 00077 /*--------. 00078 | Input. | 00079 `--------*/ 00080 00081 using namespace mln; 00082 using value::int_u8; 00083 using value::int_u16; 00084 00085 typedef int_u8 val_t; 00086 typedef image2d<val_t> orig_ima_t; 00087 00088 orig_ima_t input; 00089 io::pgm::load(input, argv[2]); 00090 if (!input.is_valid()) 00091 { 00092 std::cerr << "Error reading input " << argv[2] << std::endl; 00093 std::exit(2); 00094 } 00095 00096 /*----------------. 00097 | Line gradient. | 00098 `----------------*/ 00099 00100 // Type of the function mapping graph edges and image sites. 00101 typedef util::site_pair<point2d> P; 00102 typedef fun::i2v::array<P> fedge_site_t; 00103 00104 // Line graph image. 00105 typedef fun::i2v::array<val_t> fval_t; 00106 fval_t values; 00107 typedef edge_image<P,val_t> lg_ima_t; 00108 lg_ima_t lg_ima = morpho::line_gradient(input); 00109 00110 /*-----------. 00111 | Flooding. | 00112 `-----------*/ 00113 00114 /* FIXME: I'm not sure this is the way it should be done. Anyway, 00115 we should implement this as a canvas. */ 00116 00117 typedef lg_ima_t::nbh_t nbh_t; 00118 nbh_t nbh; 00119 00120 unsigned area = 0; 00121 unsigned max_area = input.nsites(); 00122 unsigned nregions = mln_max(unsigned); 00123 unsigned max_nregions = atoi(argv[1]); 00124 00125 lg_ima_t result = duplicate(lg_ima); 00126 while (area < max_area && nregions > max_nregions) 00127 { 00128 ++area; 00129 std::cerr << "area = " << area << " \t" 00130 << "nregions = " << nregions << std::endl; 00131 lg_ima_t work = duplicate(result); 00132 // Compute the closing. 00133 result = morpho::closing::area_on_vertices(work, nbh, area); 00134 // Compute the number of local minima (but get rid of the image, 00135 // as we don't need it). 00136 labeling::regional_minima(result, nbh, nregions); 00137 } 00138 00139 /*------. 00140 | WST. | 00141 `------*/ 00142 00143 // Perform a Watershed Transform. 00144 typedef int_u16 wst_val_t; 00145 wst_val_t nbasins; 00146 typedef edge_image<P,wst_val_t> wshed_t; 00147 00148 wshed_t wshed = morpho::meyer_wst(result, nbh, nbasins); 00149 std::cout << "nbasins = " << nbasins << std::endl; 00150 00151 /*---------. 00152 | Output. | 00153 `---------*/ 00154 00155 const wst_val_t wshed_label = 0; 00156 00157 // Create a 2D-equivalent of WSHED. 00158 image2d<wst_val_t> wshed2d(input.domain()); 00159 00160 /* FIXME: It'd better if we could iterate over the *vertices* of a 00161 line graph image. We could avoid all this lengthy code. */ 00162 // Iterate over each edge of the watershed image, and propagate the 00163 // label of an edge to its adjacent vertices when this edge is not 00164 // part of the watershed. 00165 mln_piter_(wshed_t) p(wshed.domain()); 00166 for_all(p) 00167 if (wshed(p) != wshed_label) 00168 { 00169 wshed2d(p.first()) = wshed(p); 00170 wshed2d(p.second()) = wshed(p); 00171 } 00172 00173 // For each basin, compute the average gray level. 00174 std::vector<mln_sum_(val_t)> sum(nbasins + 1, 0); 00175 std::vector<unsigned> nsites(nbasins + 1, 0); 00176 mln_piter_(orig_ima_t) q(input.domain()); 00177 for_all(q) 00178 { 00179 sum[wshed2d(q)] += input(q); 00180 ++nsites[wshed2d(q)]; 00181 } 00182 std::vector<float> average(nbasins + 1); 00183 for (unsigned i = 1; i <= nbasins; ++i) 00184 average[i] = float (sum[i]) / float (nsites[i]); 00185 00186 // Create an output image using the average gray levels of the basins. 00187 orig_ima_t output(input.domain()); 00188 for_all(q) 00189 output(q) = convert::to<mln_value_(orig_ima_t)>(average[wshed2d(q)]); 00190 00191 std::cout << "area = " << area << " \t" 00192 << "nregions = " << nregions << std::endl; 00193 io::pgm::save(output, argv[3]); 00194 }