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 00029 00030 // FIXME: Factor commons parts between mesh-segm and mesh-skel. 00031 00032 #include <cstdlib> 00033 #include <cmath> 00034 00035 #include <utility> 00036 #include <iostream> 00037 00038 #include <TriMesh.h> 00039 00040 #include <mln/core/alias/point3d.hh> 00041 #include <mln/core/alias/point3d.hh> 00042 00043 #include <mln/util/graph.hh> 00044 #include <mln/core/image/graph_image.hh> 00045 #include <mln/core/image/graph_elt_neighborhood.hh> 00046 00047 #include <mln/morpho/closing/area.hh> 00048 #include <mln/labeling/regional_minima.hh> 00049 00050 #include "io.hh" 00051 00052 00053 // Doesn't C++ have a better way to express Pi ? 00054 const float pi = 4 * atanf(1); 00055 00056 00057 int main(int argc, char* argv[]) 00058 { 00059 if (argc != 4) 00060 { 00061 std::cerr << "usage: " << argv[0] << " input.off lambda output.off" 00062 << std::endl; 00063 std::exit(1); 00064 } 00065 00066 std::string input_filename = argv[1]; 00067 unsigned lambda = atoi(argv[2]); 00068 std::string output_filename = argv[3]; 00069 00070 00071 /*-------. 00072 | Mesh. | 00073 `-------*/ 00074 00075 // TriMesh is a pain: it systematically allocates on the heap. 00076 // Introduce another name to manipulate the mesh as a (non-pointer) 00077 // object. 00078 TriMesh* mesh_ptr = TriMesh::read(input_filename.c_str()); 00079 if (!mesh_ptr) 00080 std::exit(2); 00081 TriMesh& mesh = *mesh_ptr; 00082 00083 // Computes faces (triangles). 00084 mesh.need_faces(); 00085 // Computation of the mean curvature on each vertex of the mesh. 00086 mesh.need_curvatures(); 00087 std::vector<float> vertex_h_inv(mesh.faces.size(), 0.f); 00088 for (unsigned v = 0; v < mesh.vertices.size(); ++v) 00089 { 00090 float h = (mesh.curv1[v] + mesh.curv2[v]) / 2; 00091 float h_inv = 1 / pi * atan(-h) + pi / 2; 00092 vertex_h_inv[v] = h_inv; 00093 } 00094 // Computation of the mean curvature on each face of the mesh. 00095 /* FIXME: Our implementation of the computation of the regional 00096 minima doesn't work well with floats (yet). Convert floating 00097 point values to a proportional integer value for the moment. */ 00098 typedef int curv_t; 00099 std::vector<curv_t> face_h_inv(mesh.faces.size(), 0.f); 00100 for (unsigned f = 0; f < mesh.faces.size(); ++f) 00101 { 00102 float h_inv = 00103 (vertex_h_inv[mesh.faces[f][0]] + 00104 vertex_h_inv[mesh.faces[f][1]] + 00105 vertex_h_inv[mesh.faces[f][2]]) 00106 / 3; 00107 /* FIXME: This coefficient is used to distinguish small 00108 curvature values. We should get rid of it as soon as 00109 labeling::regional_minima works correctly on images of float 00110 values. */ 00111 face_h_inv[f] = 1000 * h_inv; 00112 } 00113 00114 /*--------. 00115 | Graph. | 00116 `--------*/ 00117 00118 /* Build a graph whose vertices correspond to the faces of the mesh, 00119 whose edges (between two vertices) correspond to edges (between 00120 two faces) of the mesh. */ 00121 00122 /* FIXME: We don't have the required site type yet. Simulate with a 00123 dummy type (point3d). */ 00124 mln::util::graph<mln::point3d> g; 00125 // Populate the graph with vertices. 00126 for (unsigned i = 0; i < mesh.faces.size(); ++i) 00127 g.add_vertex (mln::point3d(i, i, i)); 00128 00129 // Populate the graph with edges. 00130 mesh.need_across_edge(); 00131 for (unsigned f = 0; f < mesh.faces.size(); ++f) 00132 for (unsigned e = 0; e < 3; ++e) 00133 { 00134 int f_adj = mesh.across_edge[f][e]; 00135 if (f_adj != -1) 00136 // Add an edge into the graph. 00137 g.add_edge(f, f_adj); 00138 } 00139 00140 /*--------------. 00141 | Graph image. | 00142 `--------------*/ 00143 00144 mln::p_graph<mln::point3d> pg(g); 00145 00146 typedef mln::graph_image<mln::point3d, curv_t> ima_t; 00147 ima_t g_ima(pg, face_h_inv); 00148 00149 /*-----------------. 00150 | Simplification. | 00151 `-----------------*/ 00152 00153 typedef mln::graph_elt_neighborhood<mln::point3d> nbh_t; 00154 nbh_t nbh; 00155 00156 ima_t closed_g_ima = mln::morpho::closing::area(g_ima, nbh, lambda); 00157 00158 /*------------------. 00159 | Regional minima. | 00160 `------------------*/ 00161 00162 typedef unsigned label_t; 00163 label_t nlabels; 00164 typedef mln::graph_image<mln::point3d, label_t> label_ima_t; 00165 label_ima_t minima = 00166 mln::labeling::regional_minima(closed_g_ima, nbh, nlabels); 00167 std::cout << "nlabels = " << nlabels << std::endl; 00168 00169 /*-----------. 00170 | Skeleton. | 00171 `-----------*/ 00172 00173 // FIXME: To do. 00174 00175 /*---------. 00176 | Output. | 00177 `---------*/ 00178 00179 /* FIXME We should created a boolean graph_image instead. But as 00180 we cannot directly save a graph_image as an OFF file now, just 00181 store the values of this would-be image into an array. */ 00182 // Assign a boolean value to graph vertices (mesh faces). 00183 std::vector<bool> face_value (minima.domain().nvertices(), true); 00184 mln_piter_(label_ima_t) pm(minima.domain()); 00185 for_all(pm) 00186 // FIXME: Use literal::zero. 00187 if (minima(pm) != 0) 00188 { 00189 // The face belongs to a regional minima: ``remove'' it from 00190 // the mesh by tagging it as false. 00191 mln_psite_(label_ima_t) pp(pm); 00192 face_value[pp.id().to_equiv()] = false; 00193 } 00194 00195 // Taken and adapted from TriMesh_io.cc 00196 FILE* f_out = fopen(output_filename.c_str(), "wb"); 00197 if (!f_out) 00198 { 00199 std::cerr << "Error opening " << output_filename.c_str() 00200 << " for writing." << std::endl; 00201 std::exit(2); 00202 } 00203 write_off_binary(mesh_ptr, face_value, f_out); 00204 fclose(f_out); 00205 00206 delete mesh_ptr; 00207 }