DGtal  1.1.0
viewMarchingCubes.cpp
Go to the documentation of this file.
1 
14 #include <iostream>
17 #include <queue>
18 #include "DGtal/base/Common.h"
19 #include "DGtal/kernel/CanonicEmbedder.h"
20 #include "DGtal/io/readers/VolReader.h"
21 #include "DGtal/images/ImageSelector.h"
22 #include "DGtal/images/imagesSetsUtils/SetFromImage.h"
23 #include "DGtal/images/ImageLinearCellEmbedder.h"
24 #include "DGtal/shapes/Shapes.h"
25 #include "DGtal/helpers/StdDefs.h"
26 #include "DGtal/topology/helpers/Surfaces.h"
27 #include "DGtal/topology/DigitalSurface.h"
28 #include "DGtal/topology/SetOfSurfels.h"
29 #include "DGtal/shapes/Mesh.h"
30 #include "DGtal/shapes/TriangulatedSurface.h"
31 #include "DGtal/shapes/MeshHelpers.h"
32 #include "DGtal/io/viewers/Viewer3D.h"
34 
36 
37 using namespace std;
38 using namespace DGtal;
39 using namespace Z3i;
40 
42 
43 void usage( int, char** argv )
44 {
45  std::cerr << "Usage: " << argv[ 0 ] << " <fileName.vol> <minT> [<maxT>=255] [<Adj>=0]" << std::endl;
46  std::cerr << "\t - displays the boundary of the shape stored in vol file <fileName.vol>" << std::endl;
47  std::cerr << "\t as a Marching-Cube triangulated surface (more precisely a dual" << std::endl;
48  std::cerr << "\t surface to the digital boundary)." << std::endl;
49  std::cerr << "\t - voxel v belongs to the shape iff its value I(v) follows minT < I(v) <= maxT." << std::endl;
50  std::cerr << "\t - minT is the iso-surface level." << std::endl;
51  std::cerr << "\t - maxT should be equal to the maximum possible value in the image." << std::endl;
52  std::cerr << "\t - 0: interior adjacency, 1: exterior adjacency (rules used to connect surface elements unambiguously)." << std::endl;
53 }
54 
55 int main( int argc, char** argv )
56 {
57  if ( argc < 5 )
58  {
59  usage( argc, argv );
60  return 1;
61  }
62  std::string inputFilename = argv[ 1 ];
63  unsigned int minThreshold = atoi( argv[ 2 ] );
64  unsigned int maxThreshold = argc > 3 ? atoi( argv[ 3 ] ) : 255;
65  bool intAdjacency = argc > 4 ? (atoi( argv[ 4 ] ) == 0) : true;
66 
68 
70  trace.beginBlock( "Reading vol file into an image." );
71  Image image = VolReader<Image>::importVol(inputFilename);
72  DigitalSet set3d (image.domain());
74  minThreshold, maxThreshold);
75  trace.endBlock();
77 
78 
80  trace.beginBlock( "Construct the Khalimsky space from the image domain." );
81  KSpace ks;
82  bool space_ok = ks.init( image.domain().lowerBound(),
83  image.domain().upperBound(), true );
84  if (!space_ok)
85  {
86  trace.error() << "Error in the Khamisky space construction."<<std::endl;
87  return 2;
88  }
89  trace.endBlock();
91 
93  typedef SurfelAdjacency<KSpace::dimension> MySurfelAdjacency;
94  MySurfelAdjacency surfAdj( intAdjacency ); // interior in all directions.
96 
98  trace.beginBlock( "Extracting boundary by scanning the space. " );
100  typedef SetOfSurfels< KSpace, SurfelSet > MySetOfSurfels;
102  MySetOfSurfels theSetOfSurfels( ks, surfAdj );
103  Surfaces<KSpace>::sMakeBoundary( theSetOfSurfels.surfelSet(),
104  ks, set3d,
105  image.domain().lowerBound(),
106  image.domain().upperBound() );
107  MyDigitalSurface digSurf( theSetOfSurfels );
108  trace.info() << "Digital surface has " << digSurf.size() << " surfels."
109  << std::endl;
110  trace.endBlock();
112 
114  trace.beginBlock( "Making triangulated surface. " );
115  typedef CanonicEmbedder< Space > TrivialEmbedder;
119  typedef Mesh< RealPoint > ViewMesh;
120  typedef std::map< MyDigitalSurface::Vertex, TriMesh::Index > VertexMap;
121  TriMesh trimesh;
122  ViewMesh viewmesh;
123  TrivialEmbedder trivialEmbedder;
124  CellEmbedder cellEmbedder;
125  // The +0.5 is to avoid isosurface going exactly through a voxel
126  // center, especially for binary volumes.
127  cellEmbedder.init( ks, image, trivialEmbedder,
128  ( (double) minThreshold ) + 0.5 );
129  VertexMap vmap; // stores the map Vertex -> Index
130  MeshHelpers::digitalSurface2DualTriangulatedSurface
131  ( digSurf, cellEmbedder, trimesh, vmap );
132  trace.info() << "Triangulated surface is " << trimesh << std::endl;
133  MeshHelpers::triangulatedSurface2Mesh( trimesh, viewmesh );
134  trace.info() << "Mesh has " << viewmesh.nbVertex()
135  << " vertices and " << viewmesh.nbFaces() << " faces." << std::endl;
136  trace.endBlock();
138 
139  QApplication application(argc,argv);
140  Viewer3D<> viewer;
141  viewer.show();
142  viewer.setLineColor(Color(150,0,0,254));
143  viewer << viewmesh;
144  viewer << Viewer3D<>::updateDisplay;
145  application.exec();
146 
147 }
148 
DGtal::Surfaces
Aim: A utility class for constructing surfaces (i.e. set of (n-1)-cells).
Definition: Surfaces.h:79
TriMesh
TriangulatedSurface< RealPoint > TriMesh
Definition: testTriangulatedSurface.cpp:52
DGtal::SetFromImage
Aim: Define utilities to convert a digital set into an image.
Definition: SetFromImage.h:64
DGtal::Trace::endBlock
double endBlock()
DGtal::KhalimskySpaceND::SurfelSet
std::set< SCell > SurfelSet
Preferred type for defining a set of surfels (always signed cells).
Definition: KhalimskySpaceND.h:450
DGtal::ImageContainerBySTLVector
Definition: ImageContainerBySTLVector.h:127
DGtal::DigitalSurface
Aim: Represents a set of n-1-cells in a nD space, together with adjacency relation between these cell...
Definition: DigitalSurface.h:140
DGtal::SurfelAdjacency< KSpace::dimension >
DGtal::Trace::error
std::ostream & error()
DGtal::Color
Structure representing an RGB triple with alpha component.
Definition: Color.h:67
DGtal::KhalimskySpaceND::init
bool init(const Point &lower, const Point &upper, bool isClosed)
Specifies the upper and lower bounds for the maximal cells in this space.
DGtal::DigitalSurface::size
Size size() const
DGtal::trace
Trace trace
Definition: Common.h:150
DGtal::Trace::beginBlock
void beginBlock(const std::string &keyword="")
DGtal::KhalimskySpaceND::lowerBound
const Point & lowerBound() const
Return the lower bound for digital points in this space.
DGtal::Trace::info
std::ostream & info()
DGtal::Viewer3D
Definition: Viewer3D.h:133
DGtal::Mesh
Aim: This class is defined to represent a surface mesh through a set of vertices and faces....
Definition: Mesh.h:92
Image
ImageContainerBySTLVector< Domain, Value > Image
Definition: testSimpleRandomAccessRangeFromPoint.cpp:45
DGtal
DGtal is the top-level namespace which contains all DGtal functions and types.
Definition: ClosedIntegerHalfPlane.h:49
DGtal::CanonicEmbedder
Aim: A trivial embedder for digital points, which corresponds to the canonic injection of Zn into Rn.
Definition: CanonicEmbedder.h:65
usage
void usage(int, char **argv)
Definition: viewMarchingCubes.cpp:43
DGtal::TriangulatedSurface
Aim: Represents a triangulated surface. The topology is stored with a half-edge data structure....
Definition: TriangulatedSurface.h:86
DGtal::Image
Aim: implements association bewteen points lying in a digital domain and values.
Definition: Image.h:70
DGtal::Display3D< SpaceND< 3 >, KhalimskySpaceND< 3 > >::setLineColor
virtual void setLineColor(DGtal::Color aColor)
DGtal::VolReader
Aim: implements methods to read a "Vol" file format.
Definition: VolReader.h:90
DGtal::SetOfSurfels
Aim: A model of CDigitalSurfaceContainer which defines the digital surface as connected surfels....
Definition: SetOfSurfels.h:74
DGtal::Viewer3D::show
virtual void show()
Overload QWidget method in order to add a call to updateList() method (to ensure that the lists are w...
DGtal::ImageLinearCellEmbedder
Aim: a cellular embedder for images. (default constructible, copy constructible, assignable)....
Definition: ImageLinearCellEmbedder.h:70
RealPoint
Z2i::RealPoint RealPoint
Definition: testAstroid2D.cpp:46
Value
double Value
Definition: testSimpleRandomAccessRangeFromPoint.cpp:38
MyDigitalSurface
DigitalSurface< MyDigitalSurfaceContainer > MyDigitalSurface
Definition: greedy-plane-segmentation-ex2.cpp:92
SurfelSet
MyDigitalSurface::SurfelSet SurfelSet
Definition: greedy-plane-segmentation-ex2.cpp:95
DGtal::KhalimskySpaceND
Aim: This class is a model of CCellularGridSpaceND. It represents the cubical grid as a cell complex,...
Definition: KhalimskySpaceND.h:394
main
int main(int argc, char **argv)
Definition: viewMarchingCubes.cpp:55