DGtalTools  1.2.0
volSurfaceRegularization.cpp
1 
89 #include <DGtal/base/Common.h>
90 #include <DGtal/helpers/StdDefs.h>
91 
92 #include "CLI11.hpp"
93 
94 #include <DGtal/images/ImageSelector.h>
95 #include <DGtal/io/readers/GenericReader.h>
96 #include <DGtal/topology/SurfelNeighborhood.h>
97 #include <DGtal/topology/helpers/Surfaces.h>
98 #include <DGtal/geometry/volumes/distance/VoronoiMap.h>
99 #include <random>
100 #include <iomanip>
101 #include <regex>
102 #include "volSurfaceRegularization-details/surface_approx.h"
103 #include "volSurfaceRegularization-details/surface_extract.h"
104 #include "volSurfaceRegularization-details/shape.h"
105 
106 
107 struct Options
108 {
109  double noise_level {0};
110  std::string image_filename;
111  double normal_radius {4} ;
112  double regularization_position {1e-3};
113  double regularization_center {1e-2};
114  double align {1.0};
115  double fairness{0.0};
116  double barycenter{1e-1};
117  std::string regularized_obj_filename;
118  std::string cubical_obj_filename;
119 };
120 
121 
122 
123 
124 
125 
126 int main(int argc, char* argv[])
127 {
128  using DGtal::trace;
129  using std::endl;
130  using DGtal::PRIMAL;
131  using DGtal::DUAL;
132 
133  // parse command line using CLI ----------------------------------------------
134  CLI::App app;
135  Options options;
136 
137  app.description("Regularize a cubical complex into a smooth quadrangulated complex.");
138 
139  app.add_option("--image-filename,-i,1",options.image_filename, "input vol filename for image shape (object voxels have values > 0) or input cvs filename for surfels and normals") -> required();
140  app.add_option("--regularized-obj-filename,-o,1", options.regularized_obj_filename, "output regularized obj") -> required();
141  app.add_option("--cubical-obj-filename,-n", options.cubical_obj_filename, "output cubical obj");
142  app.add_option("--shape-noise,-k", options.noise_level,"noise shape parameter", true );
143  auto groupNormEst = app.add_option_group("Normal field estimator options");
144  groupNormEst->add_option("--normal-radius,-r", options.normal_radius, "radius of normal estimator", true);
145 
146  auto groupApprox = app.add_option_group("Surface approximation options");
147  groupApprox->add_option("--regularization-position,-p", options.regularization_position, "vertex position regularization coeff", true );
148  groupApprox->add_option("--regularization-center,-c",options.regularization_center, "face center regularization coeff", true);
149  groupApprox->add_option("--align,-a", options.align, "normal alignment coeff", true);
150  groupApprox->add_option("--fairness,-f",options.fairness,"face fairness coeff" ,true);
151  groupApprox->add_option("--barycenter,-b", options.barycenter, "barycenter fairness coeff", true);
152 
153  app.get_formatter()->column_width(40);
154  CLI11_PARSE(app, argc, argv);
155  // END parse command line using CLI ----------------------------------------------
156 
157  const KSpace kspace;
158 
159  const Calculus calculus;
160  const FlatVector original_face_normals;
161 
162  if (!ends_with(options.image_filename, ".csv"))
163  {
164  ASSERT( !options.image_filename.empty() );
165  typedef DGtal::Z3i::Domain Domain;
166  typedef DGtal::ImageSelector<Domain, unsigned char>::Type ImageUChar;
167  trace.info() << "image_filename=" << options.image_filename << endl;
168  ImageUChar image_uchar = DGtal::GenericReader<ImageUChar>::import(options.image_filename);
169  const Domain domain = image_uchar.domain();
170  trace.info() << "domain=" << domain << endl;
171  const Point center = (domain.upperBound()+domain.lowerBound())/2;
172  trace.info() << "center=" << center << endl;
173  const ImageShape<ImageUChar> shape(&image_uchar, center);
174 
175  const_cast<KSpace&>(kspace).init(domain.lowerBound()-center-Point::diagonal(1), domain.upperBound()-center+Point::diagonal(1), true);
176 
177  std::tie(const_cast<Calculus&>(calculus), const_cast<FlatVector&>(original_face_normals)) =
178  initCalculusAndNormalsWithNoise(kspace, shape, options.normal_radius, options.noise_level);
179  }
180 
181  if (ends_with(options.image_filename, ".csv"))
182  {
183  ASSERT( !options.image_filename.empty() );
184  trace.info() << "csv_filename=" << options.image_filename << endl;
185  std::tie(const_cast<Calculus&>(calculus), const_cast<FlatVector&>(original_face_normals)) =
186  initCalculusAndNormalsFromSurfelNormalsCSV(options.image_filename);
187  }
188 
189  const FlatVector original_vertex_normals = vertexNormals(calculus, original_face_normals);
190 
191  const FlatVector original_positions;
192  const FlatVector regularized_positions;
193  const FlatVector original_centers;
194  const FlatVector regularized_centers;
195  std::tie(const_cast<FlatVector&>(original_positions),
196  const_cast<FlatVector&>(regularized_positions),
197  const_cast<FlatVector&>(original_centers),
198  const_cast<FlatVector&>(regularized_centers)) =
199  approximateSurface(calculus, original_face_normals,
200  ApproxParams({options.regularization_position, options.regularization_center, options.align, options.fairness, options.barycenter}));
201  ASSERT( original_positions.size() == 3*calculus.kFormLength(0, PRIMAL) );
202  ASSERT( regularized_positions.size() == 3*calculus.kFormLength(0, PRIMAL) );
203  ASSERT( original_centers.size() == 3*calculus.kFormLength(2, PRIMAL) );
204  ASSERT( regularized_centers.size() == 3*calculus.kFormLength(2, PRIMAL) );
205 
206  {
207  trace.beginBlock( "computing energies" );
208 
209  {
210  double position_energy = 0;
211  double align_energy = 0;
212  std::tie( position_energy, align_energy ) = approximateSurfaceEnergies(
213  calculus, original_face_normals, original_positions );
214  align_energy *= options.align;
215  position_energy *= options.regularization_position;
216  trace.info() << "original_energies=" << position_energy << " "
217  << align_energy << " " << position_energy + align_energy
218  << endl;
219  }
220 
221  {
222  double position_energy = 0;
223  double align_energy = 0;
224  std::tie(position_energy, align_energy) = approximateSurfaceEnergies(calculus, original_face_normals, regularized_positions);
225  align_energy *= options.align;
226  position_energy *= options.regularization_position;
227  trace.info() << "regularized_energies=" << position_energy << " " << align_energy << " " << position_energy+align_energy << endl;
228  }
229 
230  trace.endBlock();
231  }
232 
233  {
234  ASSERT( !options.regularized_obj_filename.empty() );
235  trace.info() << "regularized_obj_filename=" << options.regularized_obj_filename << endl;
236  exportOBJ(calculus, regularized_positions, options.regularized_obj_filename);
237  }
238 
239  if (!options.cubical_obj_filename.empty())
240  {
241  ASSERT( !options.cubical_obj_filename.empty() );
242  trace.info() << "cubical_obj_filename=" << options.cubical_obj_filename << endl;
243  exportOBJ(calculus, original_positions, options.cubical_obj_filename);
244  }
245 
246  return 0;
247 }
248 
249