31 #include <QGLViewer/qglviewer.h> 34 #include "DGtal/base/Common.h" 35 #include "DGtal/helpers/StdDefs.h" 36 #include "DGtal/io/viewers/Viewer3D.h" 37 #include "DGtal/io/DrawWithDisplay3DModifier.h" 38 #include "DGtal/io/readers/PointListReader.h" 39 #include "DGtal/io/readers/MeshReader.h" 40 #include "DGtal/io/colormaps/GradientColorMap.h" 42 #include "DGtal/topology/helpers/Surfaces.h" 43 #include "DGtal/topology/SurfelAdjacency.h" 44 #include "DGtal/topology/CanonicCellEmbedder.h" 45 #include "DGtal/math/Statistic.h" 47 #include "DGtal/io/Color.h" 48 #include "DGtal/io/colormaps/GradientColorMap.h" 49 #include "DGtal/io/readers/GenericReader.h" 51 #include <boost/program_options/options_description.hpp> 52 #include <boost/program_options/parsers.hpp> 53 #include <boost/program_options/variables_map.hpp> 58 using namespace DGtal;
62 namespace po = boost::program_options;
130 template <
typename Space = DGtal::Z3i::Space,
typename KSpace = DGtal::Z3i::KSpace>
131 struct ViewerSnap: DGtal::Viewer3D <Space, KSpace>
134 ViewerSnap(
const KSpace &KSEmb,
bool saveSnap): Viewer3D<Space, KSpace>(KSEmb), mySaveSnap(saveSnap){
139 DGtal::Viewer3D<>::init();
141 QObject::connect(
this, SIGNAL(drawFinished(
bool)),
this, SLOT(saveSnapshot(
bool)));
148 template <
typename Po
int>
150 getBoundingUpperAndLowerPoint(
const std::vector<Point> &vectorPt, Point &ptLower, Point &ptUpper){
151 for(
unsigned int i =1; i<vectorPt.size(); i++){
152 if(vectorPt.at(i)[0] < ptLower[0]){
153 ptLower[0] = vectorPt.at(i)[0];
155 if(vectorPt.at(i)[1] < ptLower[1]){
156 ptLower[1] = vectorPt.at(i)[1];
158 if(vectorPt.at(i)[2] < ptLower[2]){
159 ptLower[2] =vectorPt.at(i)[2];
161 if(vectorPt.at(i)[0] < ptLower[0]){
162 ptLower[0] = vectorPt.at(i)[0];
164 if(vectorPt.at(i)[1] < ptLower[1]){
165 ptLower[1] = vectorPt.at(i)[1];
167 if(vectorPt.at(i)[2] < ptLower[2]){
168 ptLower[2] =vectorPt.at(i)[2];
174 int main(
int argc,
char** argv )
177 typedef PointVector<4, double> Point4D;
178 typedef PointVector<1, int> Point1D;
181 po::options_description general_opt(
"Allowed options are: ");
182 general_opt.add_options()
183 (
"help,h",
"display this message")
184 (
"input,i", po::value<std::string>(),
"input file: sdpa (sequence of discrete points with attribute)" )
185 (
"reference,r", po::value<std::string>(),
"input reference file: sdpa (sequence of discrete points with attribute)" )
186 (
"compAccordingLabels,l",
"apply the comparisos only on points with same labels (by default fifth colomn)" )
187 (
"drawSurfelAssociations,a",
"Draw the surfel association." )
188 (
"fileMeasureOutput,o", po::value<std::string>(),
"specify the output file to store (append) the error stats else the result is given to std output. " )
189 (
"noWindows,n",
"Don't display Viewer windows." )
190 (
"doSnapShotAndExit,d", po::value<std::string>(),
"save display snapshot into file. Notes that the camera setting is set by default according the last saved configuration (use SHIFT+Key_M to save current camera setting in the Viewer3D)." )
191 (
"fixMaxColorValue", po::value<double>(),
"fix the maximal color value for the scale error display (else the scale is set from the maximal value)" )
192 (
"labelIndex", po::value<unsigned int>(),
"set the index of the label (by default set to 4) " )
193 (
"SDPindex", po::value<std::vector <unsigned int> >()->multitoken(),
"specify the sdp index (by default 0,1,2,3).");
197 bool cannotStart=
false;
198 po::variables_map vm;
202 po::store(po::parse_command_line(argc, argv, general_opt), vm);
203 }
catch(
const std::exception& ex){
205 trace.error()<<
"Error checking program options: "<< ex.what()<< endl;
208 if(parseOK && ! vm.count(
"input"))
210 trace.error() <<
" The input file name was not defined" << endl;
215 if( !parseOK || cannotStart || vm.count(
"help")||argc<=1)
217 trace.info() <<
"Usage: " << argv[0] <<
" [input]\n" 218 <<
"It computes generic scalar surfel data comparisons (squared error) ( given from an input data file and from a reference one). \n \n" 219 <<
"Each surfels are associated to the nearest one of the reference surfels (computed by a 'brut force' search) " 220 <<
"This association can also be limited to surfel of same label (if available in the data and by using the --compAccordingLabels option )." 221 <<
"The comparison and surfel association can be displayed and result statistics are saved on output file (--fileMeasureOutput)." 222 <<
"You can also remove the interactive 3d view by just doing a snapshot and exit with option --doSnapShotAndExit. \n \n" 223 <<
"Example of use: \n \n" 224 <<
"3dCompSurfelData -i surfelCurvatureInput.dat -r surfelCurvatureRef.dat --fixMaxColorValue 0.8 -d visuSEcurvature.png -o statMeasures.dat \n \n " 225 <<
"=> From the two compared files you should obtain the result of the comparison (statMeasures.dat) with the associated visualisation (visuSEcurvature.png). \n \n \n " 226 << general_opt <<
"\n";
231 string inputFilename = vm[
"input"].as<std::string>();
232 string referenceFilename = vm[
"reference"].as<std::string>();
234 std::vector<Point4D> surfelAndScalarInput;
235 std::vector<Point4D> surfelAndScalarReference;
236 std::vector<Point1D> vectLabelsInput;
237 std::vector<Point1D> vectLabelsReference;
238 bool useLabels = vm.count(
"compAccordingLabels");
241 if(vm.count(
"labelIndex")){
242 std::vector<unsigned int > vectIndex;
243 vectIndex.push_back(vm[
"lSDPindex"].as<unsigned int >());
244 vectLabelsInput = PointListReader<Point1D>::getPointsFromFile(inputFilename, vectIndex);
245 vectLabelsReference = PointListReader<Point1D>::getPointsFromFile(referenceFilename, vectIndex);
247 vectLabelsInput = PointListReader<Point1D>::getPointsFromFile(inputFilename);
248 vectLabelsReference = PointListReader<Point1D>::getPointsFromFile(referenceFilename);
254 if(vm.count(
"SDPindex")) {
255 std::vector<unsigned int > vectIndex = vm[
"SDPindex"].as<std::vector<unsigned int > >();
256 if(vectIndex.size()!=4){
257 trace.error() <<
"you need to specify the three indexes of vertex." << std::endl;
260 surfelAndScalarInput = PointListReader<Point4D>::getPointsFromFile(inputFilename, vectIndex);
261 surfelAndScalarReference = PointListReader<Point4D>::getPointsFromFile(referenceFilename, vectIndex);
263 surfelAndScalarInput = PointListReader<Point4D>::getPointsFromFile(inputFilename);
264 surfelAndScalarReference = PointListReader<Point4D>::getPointsFromFile(referenceFilename);
268 Point4D ptLower = surfelAndScalarReference.at(0);
269 Point4D ptUpper = surfelAndScalarReference.at(0);
270 getBoundingUpperAndLowerPoint(surfelAndScalarReference, ptLower, ptUpper);
271 getBoundingUpperAndLowerPoint(surfelAndScalarInput, ptLower, ptUpper);
273 K.init(Z3i::Point(2*ptLower[0]+1, 2*ptLower[1]+1, 2*ptLower[2]+1),
274 Z3i::Point(2*ptUpper[0]+1, 2*ptUpper[1]+1, 2*ptUpper[2]+1),
true);
277 std::vector<Cell> vectSurfelsReference;
278 std::vector<Cell> vectSurfelsInput;
281 for(
unsigned int i =0; i<surfelAndScalarReference.size(); i++){
282 Point4D pt4d = surfelAndScalarReference.at(i);
283 Cell c = K.uCell(Z3i::Point(pt4d[0], pt4d[1], pt4d[2]));
284 vectSurfelsReference.push_back(c);
287 for(
unsigned int i =0; i<surfelAndScalarInput.size(); i++){
288 Point4D pt4d = surfelAndScalarInput.at(i);
289 Cell c = K.uCell(Z3i::Point(pt4d[0], pt4d[1], pt4d[2]));
290 vectSurfelsInput.push_back(c);
299 CanonicCellEmbedder<KSpace> embeder(K);
300 std::vector<unsigned int> vectIndexMinToReference;
302 trace.info() <<
"Associating each input surfel to reference:" << std::endl;
303 trace.progressBar(0, surfelAndScalarInput.size());
305 for(
unsigned int i=0;i <vectSurfelsInput.size(); i++){
306 trace.progressBar(i, vectSurfelsInput.size());
307 Z3i::RealPoint ptCenterInput = embeder(vectSurfelsInput.at(i));
308 unsigned int indexDistanceMin = 0;
309 double distanceMin = std::numeric_limits<int>::max() ;
310 for(
unsigned int j=0; j <vectSurfelsReference.size(); j++){
311 if(useLabels && vectLabelsReference.at(j) != vectLabelsInput.at(i)){
314 Z3i::RealPoint ptCenterRef = embeder(vectSurfelsReference.at(j));
315 double distance = (ptCenterRef - ptCenterInput).norm();
316 if(distance < distanceMin){
317 distanceMin = distance;
318 indexDistanceMin = j;
321 vectIndexMinToReference.push_back(indexDistanceMin);
324 trace.progressBar(surfelAndScalarInput.size(), surfelAndScalarInput.size());
325 trace.info() << std::endl;
333 QApplication application(argc,argv);
334 typedef ViewerSnap<> Viewer;
337 Viewer viewer(K, vm.count(
"doSnapShotAndExit"));
338 if(vm.count(
"doSnapShotAndExit")){
339 viewer.setSnapshotFileName(QString(vm[
"doSnapShotAndExit"].as<std::string>().c_str()));
341 viewer.setWindowTitle(
"3dCompSurfel Viewer");
344 Statistic<double> statErrors(
true);
345 std::ofstream outputStatStream;
346 if(vm.count(
"fileMeasureOutput")){
347 outputStatStream.open(vm[
"fileMeasureOutput"].as<std::string>().c_str(), ios::app );
353 for(
unsigned int i=0;i <surfelAndScalarInput.size(); i++){
354 double scalarInput = surfelAndScalarInput.at(i)[3];
355 double scalarRef = surfelAndScalarReference.at(vectIndexMinToReference.at(i))[3];
356 double sqError = (scalarRef-scalarInput)*(scalarRef-scalarInput);
357 statErrors.addValue(sqError);
358 if(sqError> maxSqError){
362 double maxVal = vm.count(
"fixMaxColorValue")? vm[
"fixMaxColorValue"].as<
double>(): maxSqError;
364 GradientColorMap<double> gradientColorMap( 0, maxVal );
365 gradientColorMap.addColor( Color(255,255,255,100 ));
366 gradientColorMap.addColor( Color(255,0,0,100 ) );
367 gradientColorMap.addColor( Color(0,0,255,100 ) );
374 viewer << SetMode3D(vectSurfelsInput.at(0).className(),
"Basic");
375 for(
unsigned int i=0; i <surfelAndScalarInput.size(); i++){
376 double scalarInput = surfelAndScalarInput.at(i)[3];
377 double scalarRef = surfelAndScalarReference.at(vectIndexMinToReference.at(i))[3];
378 double sqError = (scalarRef-scalarInput)*(scalarRef-scalarInput);
379 viewer.setFillColor(gradientColorMap(sqError));
381 viewer << vectSurfelsInput.at(i);
382 if(vm.count(
"drawSurfelAssociations")){
383 viewer.addLine(embeder(vectSurfelsInput.at(i)),embeder(vectSurfelsReference.at(vectIndexMinToReference.at(i))));
387 statErrors.terminate();
388 if(vm.count(
"fileMeasureOutput")){
389 outputStatStream <<
"input= " << inputFilename <<
" reference=" << referenceFilename <<
" " ;
390 outputStatStream << statErrors << std::endl;
392 trace.info() << statErrors;
394 viewer << Viewer::updateDisplay;
396 if(vm.count(
"doSnapShotAndExit")){
398 viewer.restoreStateFromFile();
399 std::string name = vm[
"doSnapShotAndExit"].as<std::string>();
400 std::string extension = name.substr(name.find_last_of(
".") + 1);
401 std::string basename = name.substr(0, name.find_last_of(
"."));
402 for(
int i=0; i< viewer.snapshotCounter()-1; i++){
404 s << basename <<
"-"<< setfill(
'0') << setw(4)<< i <<
"." << extension;
405 trace.info() <<
"erase temp file: " << s.str() << std::endl;
406 remove(s.str().c_str());
409 s << basename <<
"-"<< setfill(
'0') << setw(4)<< viewer.snapshotCounter()-1 <<
"." << extension;
410 rename(s.str().c_str(), name.c_str());
414 if(vm.count(
"noWindows")){
417 return application.exec();