DGtal  0.9.3beta
greedy-plane-segmentation.cpp
1 
30 #include <iostream>
32 #include <vector>
33 #include <set>
34 #include <map>
35 #include <iostream>
36 
37 #include "DGtal/base/Common.h"
38 #include "DGtal/io/readers/VolReader.h"
39 
40 #include "DGtal/io/Display3D.h"
41 #include "DGtal/io/viewers/Viewer3D.h"
42 #include "DGtal/io/DrawWithDisplay3DModifier.h"
43 #include "DGtal/images/ImageSelector.h"
44 #include "DGtal/images/imagesSetsUtils/SetFromImage.h"
45 #include "DGtal/topology/DigitalSurface.h"
46 #include "DGtal/topology/DigitalSetBoundary.h"
47 #include "DGtal/graph/BreadthFirstVisitor.h"
48 #include "DGtal/geometry/surfaces/COBANaivePlaneComputer.h"
49 #include "DGtal/helpers/StdDefs.h"
50 #include "ConfigExamples.h"
51 
53 
54 using namespace std;
55 using namespace DGtal;
56 
57 
59 using namespace Z3i;
60 typedef DGtal::int64_t InternalInteger;
62 // We choose the DigitalSetBoundary surface container in order to
63 // segment connected or unconnected surfaces.
66 typedef MyDigitalSurface::ConstIterator ConstIterator;
67 typedef MyDigitalSurface::Vertex Vertex;
68 typedef MyDigitalSurface::SurfelSet SurfelSet;
69 typedef SurfelSet::iterator SurfelSetIterator;
75 struct SegmentedPlane {
76  NaivePlaneComputer plane;
77  Color color;
78 };
80 
82 
83 int main( int argc, char** argv )
84 {
85 
87  trace.info() << "Segments the surface at given threshold within given volume into digital planes of rational width num/den." << std::endl;
88  // Setting default options: ----------------------------------------------
89  // input file used:
90  string inputFilename = examplesPath + "samples/Al.100.vol" ;
91  trace.info() << "input file used " << inputFilename << std::endl;
92  // parameter threshold
93  unsigned int threshold = 0;
94  trace.info() << "the value that defines the isosurface in the image (an integer between 0 and 255)= " << threshold<< std::endl;
95  // parameter widthNum
96  unsigned int widthNum = 1;
97  trace.info() << "the numerator of the rational width (a non-null integer) =" << widthNum<< std::endl;
98  // parameter widthDen
99  unsigned int widthDen = 1;
100  trace.info() << "the denominator of the rational width (a non-null integer)= " << widthDen<< std::endl;
101 
103 
105  QApplication application(argc,argv);
107  Image image = VolReader<Image>::importVol(inputFilename);
108  DigitalSet set3d (image.domain());
109  SetFromImage<DigitalSet>::append<Image>(set3d, image, threshold,255);
111 
113  trace.beginBlock( "Set up digital surface." );
114  // We initializes the cellular grid space used for defining the
115  // digital surface.
116  KSpace ks;
117  bool ok = ks.init( set3d.domain().lowerBound(),
118  set3d.domain().upperBound(), true );
119  if ( ! ok ) std::cerr << "[KSpace.init] Failed." << std::endl;
120  SurfelAdjacency<KSpace::dimension> surfAdj( true ); // interior in all directions.
121  MyDigitalSurfaceContainer* ptrSurfContainer =
122  new MyDigitalSurfaceContainer( ks, set3d, surfAdj );
123  MyDigitalSurface digSurf( ptrSurfContainer ); // acquired
124  trace.endBlock();
126 
128  trace.beginBlock( "Segment into planes." );
129  std::set<Vertex> processedVertices;
130  std::vector<SegmentedPlane*> segmentedPlanes;
131  std::map<Vertex,SegmentedPlane*> v2plane;
132  Point p;
133  Dimension axis;
134  unsigned int j = 0;
135  unsigned int nb = digSurf.size();
136  for ( ConstIterator it = digSurf.begin(), itE= digSurf.end(); it != itE; ++it )
137  {
138  if ( ( (++j) % 50 == 0 ) || ( j == nb ) ) trace.progressBar( j, nb );
139  Vertex v = *it;
140  if ( processedVertices.find( v ) != processedVertices.end() ) // already in set
141  continue; // process to next vertex
142 
143  SegmentedPlane* ptrSegment = new SegmentedPlane;
144  segmentedPlanes.push_back( ptrSegment ); // to delete them afterwards.
145  axis = ks.sOrthDir( v );
146  ptrSegment->plane.init( axis, 500, widthNum, widthDen );
147  // The visitor takes care of all the breadth-first traversal.
148  Visitor visitor( digSurf, v );
149  while ( ! visitor.finished() )
150  {
151  Visitor::Node node = visitor.current();
152  v = node.first;
153  if ( processedVertices.find( v ) == processedVertices.end() )
154  { // Vertex is not in processedVertices
155  axis = ks.sOrthDir( v );
156  p = ks.sCoords( ks.sDirectIncident( v, axis ) );
157  bool isExtended = ptrSegment->plane.extend( p );
158  if ( isExtended )
159  { // surfel is in plane.
160  processedVertices.insert( v );
161  v2plane[ v ] = ptrSegment;
162  visitor.expand();
163  }
164  else // surfel is not in plane and should not be used in the visit.
165  visitor.ignore();
166  }
167  else // surfel is already in some plane.
168  visitor.ignore();
169  }
170  // Assign random color for each plane.
171  ptrSegment->color = Color( rand() % 256, rand() % 256, rand() % 256, 255 );
172  }
173  trace.endBlock();
175 
177  Viewer3D<> viewer( ks );
178  viewer.show();
179  for ( std::map<Vertex,SegmentedPlane*>::const_iterator
180  it = v2plane.begin(), itE = v2plane.end();
181  it != itE; ++it )
182  {
183  viewer << CustomColors3D( it->second->color, it->second->color );
184  viewer << ks.unsigns( it->first );
185  }
186  viewer << Viewer3D<>::updateDisplay;
188 
190  for ( std::vector<SegmentedPlane*>::iterator
191  it = segmentedPlanes.begin(), itE = segmentedPlanes.end();
192  it != itE; ++it )
193  delete *it;
194  segmentedPlanes.clear();
195  v2plane.clear();
197 
198  return application.exec();
199 }
void beginBlock(const std::string &keyword="")
void progressBar(const double currentValue, const double maximalValue)
Aim: implements association bewteen points lying in a digital domain and values.
Definition: Image.h:69
Trace trace
Definition: Common.h:137
DigitalSurfaceContainer::SurfelConstIterator ConstIterator
DGtal::uint32_t Dimension
Definition: Common.h:120
KSpace::SurfelSet SurfelSet
STL namespace.
double endBlock()
Aim: Represents a set of n-1-cells in a nD space, together with adjacency relation between these cell...
Aim: A model of CDigitalSurfaceContainer which defines the digital surface as the boundary of a given...
SCell sDirectIncident(const SCell &p, Dimension k) const
Dimension sOrthDir(const SCell &s) const
bool init(const Point &lower, const Point &upper, bool isClosed)
Aim: implements methods to read a "Vol" file format.
Definition: VolReader.h:88
Surfel Vertex
Defines the type for a vertex.
Point sCoords(const SCell &c) const
DGtal is the top-level namespace which contains all DGtal functions and types.
Aim: A wrapper class around a STL associative container for storing sets of digital points within som...
Aim: Define utilities to convert a digital set into an image.
Definition: SetFromImage.h:63
static Dimension size()
std::pair< Vertex, Data > Node
FIXME.
std::ostream & info()
Aim: This class is useful to perform a breadth-first exploration of a graph given a starting point or...
Structure representing an RGB triple with alpha component.
Definition: Color.h:66
const Domain & domain() const
Definition: Image.h:192
Cell unsigns(const SCell &p) const
boost::int64_t int64_t
signed 94-bit integer.
Definition: BasicTypes.h:74
Aim: This class is a model of CCellularGridSpaceND. It represents the cubical grid as a cell complex...