DGtal  0.9.3
cubical-complex-collapse.cpp
Go to the documentation of this file.
1 
67 #include <iostream>
69 #include <map>
70 #include "DGtal/base/Common.h"
71 #include "DGtal/helpers/StdDefs.h"
72 
73 #include "DGtal/io/DrawWithDisplay3DModifier.h"
74 #include "DGtal/io/viewers/Viewer3D.h"
75 #include "DGtal/topology/KhalimskySpaceND.h"
76 #include "DGtal/topology/CubicalComplex.h"
77 #include "DGtal/topology/CubicalComplexFunctions.h"
78 
80 
81 using namespace std;
82 using namespace DGtal;
83 using namespace DGtal::Z3i;
84 
85 
91 template <typename CubicalComplex>
92 struct DiagonalPriority {
93  typedef typename CubicalComplex::KSpace KSpace;
94  typedef typename CubicalComplex::Point Point;
95  typedef typename CubicalComplex::Cell Cell;
96  typedef typename CubicalComplex::CellMapIterator CellMapIterator;
97 
98  DiagonalPriority( const CubicalComplex& complex ) : myComplex( complex ) {}
99  bool operator()( const CellMapIterator& it1, const CellMapIterator& it2 ) const
100  {
101  Point k1 = myComplex.space().uKCoords( it1->first );
102  Point k2 = myComplex.space().uKCoords( it2->first );
103  double d1 = Point::diagonal( 1 ).dot( k1 ) / sqrt( (double) KSpace::dimension );
104  double d2 = Point::diagonal( 1 ).dot( k2 ) / sqrt( (double) KSpace::dimension );;
105  RealPoint v1( k1[ 0 ] - d1 * k1[ 0 ], k1[ 1 ] - d1 * k1[ 1 ], k1[ 2 ] - d1 * k1[ 2 ] );
106  RealPoint v2( k2[ 0 ] - d2 * k2[ 0 ], k2[ 1 ] - d2 * k2[ 1 ], k2[ 2 ] - d2 * k2[ 2 ] );
107  double n1 = v1.dot( v1 );
108  double n2 = v2.dot( v2 );
109  return ( n1 < n2 ) || ( ( n1 == n2 ) && ( it1->first < it2->first ) );
110  }
111 
112  const CubicalComplex& myComplex;
113 };
115 
116 int main( int argc, char** argv )
117 {
118  // JOL: unordered_map is approximately twice faster than map for
119  // collapsing.
120  typedef std::map<Cell, CubicalCellData> Map;
121  // typedef boost::unordered_map<Cell, CubicalCellData> Map;
123 
124  trace.beginBlock( "Creating Cubical Complex" );
125  KSpace K;
126  K.init( Point( 0,0,0 ), Point( 512,512,512 ), true );
127  CC complex( K );
128  Integer m = 40;
129  std::vector<Cell> S;
130  for ( Integer x = 0; x <= m; ++x )
131  for ( Integer y = 0; y <= m; ++y )
132  for ( Integer z = 0; z <= m; ++z )
133  {
134  Point k1 = Point( x, y, z );
135  S.push_back( K.uCell( k1 ) );
136  double d1 = Point::diagonal( 1 ).dot( k1 ) / (double) KSpace::dimension; // sqrt( (double) KSpace::dimension );
137  RealPoint v1( k1[ 0 ], k1[ 1 ], k1[ 2 ] );
138  v1 -= d1 * RealPoint::diagonal( 1.0 );
139  //RealPoint v1( k1[ 0 ] - d1 * k1[ 0 ], k1[ 1 ] - d1 * k1[ 1 ], k1[ 2 ] - d1 * k1[ 2 ] );
140  double n1 = v1.norm();
141  bool fixed = ( ( x == 0 ) && ( y == 0 ) && ( z == 0 ) )
142  || ( ( x == 0 ) && ( y == m ) && ( z == 0 ) )
143  || ( ( x == m ) && ( y == 0 ) && ( z == 0 ) )
144  || ( ( x == m ) && ( y == m ) && ( z == 0 ) )
145  || ( ( x == m/3 ) && ( y == 2*m/3 ) && ( z == 2*m/3 ) )
146  || ( ( x == 0 ) && ( y == 0 ) && ( z == m ) )
147  || ( ( x == 0 ) && ( y == m ) && ( z == m ) )
148  || ( ( x == m ) && ( y == 0 ) && ( z == m ) )
149  || ( ( x == m ) && ( y == m ) && ( z == m ) )
150  || ( ( x == 0 ) && ( y == m ) )
151  || ( ( x == m ) && ( y == m ) )
152  || ( ( z == 0 ) && ( y == m ) )
153  || ( ( z == m ) && ( y == m ) );
154  complex.insertCell( S.back(),
155  fixed ? CC::FIXED
156  : (DGtal::uint32_t) floor(64.0 * n1 ) // This is the priority for collapse
157  );
158  }
159  //complex.close();
160  trace.info() << "After close: " << complex << std::endl;
161  trace.endBlock();
162 
163  // for 3D display with Viewer3D
164  QApplication application(argc,argv);
166 
167  {
168  MyViewer viewer(K);
169  viewer.show();
171  for ( Dimension d = 0; d <= 3; ++d )
172  for ( CellMapConstIterator it = complex.begin( d ), itE = complex.end( d );
173  it != itE; ++it )
174  {
175  bool fixed = (it->second.data == CC::FIXED);
176  if ( fixed ) viewer.setFillColor( Color::Red );
177  else viewer.setFillColor( Color::White );
178  viewer << it->first;
179  }
180  viewer<< MyViewer::updateDisplay;
181  application.exec();
182  }
183 
184  trace.beginBlock( "Collapsing complex" );
186  DGtal::uint64_t removed
187  = functions::collapse( complex, S.begin(), S.end(), P, true, true, true );
188  trace.info() << "Collapse removed " << removed << " cells." << std::endl;
189  trace.info() << "After collapse: " << complex << std::endl;
190  trace.endBlock();
191 
192  {
193  MyViewer viewer(K);
194  viewer.show();
196  for ( Dimension d = 0; d <= 3; ++d )
197  for ( CellMapConstIterator it = complex.begin( d ), itE = complex.end( d );
198  it != itE; ++it )
199  {
200  bool fixed = (it->second.data == CC::FIXED);
201  if ( fixed ) viewer.setFillColor( Color::Red );
202  else viewer.setFillColor( Color::White );
203  viewer << it->first;
204  }
205  viewer<< MyViewer::updateDisplay;
206  return application.exec();
207  }
208 }
virtual void setFillColor(DGtal::Color aColor)
void beginBlock(const std::string &keyword="")
CellMap::iterator CellMapIterator
Iterator for visiting type CellMap.
std::unordered_map< Cell, CubicalCellData > Map
boost::uint32_t uint32_t
unsigned 32-bit integer.
Definition: BasicTypes.h:63
ConstIterator begin() const
Z3i this namespace gathers the standard of types for 3D imagery.
Trace trace
Definition: Common.h:137
virtual void show()
Overload QWidget method in order to add a call to updateList() method (to ensure that the lists are w...
DGtal::uint32_t Dimension
Definition: Common.h:120
STL namespace.
double endBlock()
TKSpace KSpace
Type of the cellular grid space.
CellMap::const_iterator CellMapConstIterator
Const iterator for visiting type CellMap.
int main(int argc, char **argv)
void insertCell(const Cell &aCell, const Data &data=Data())
CC::CellMapConstIterator CellMapConstIterator
bool init(const Point &lower, const Point &upper, bool isClosed)
double norm(const NormType type=L_2) const
Cell uCell(const PreCell &c) const
DGtal::int32_t Integer
Definition: StdDefs.h:143
Aim: This class represents an arbitrary cubical complex living in some Khalimsky space. Cubical complexes are sets of cells of different dimensions related together with incidence relations. Two cells in a cubical complex are incident if and only if they are incident in the surrounding Khalimsky space. In other words, cubical complexes are defined here as subsets of Khalimsky spaces.
DGtal is the top-level namespace which contains all DGtal functions and types.
MyPointD Point
Definition: testClone2.cpp:383
std::ostream & info()
CubicalComplex< KSpace, Map > CC
ConstIterator end() const
boost::uint64_t uint64_t
unsigned 64-bit integer.
Definition: BasicTypes.h:65
KSpace K
KSpace::Cell Cell
Aim: This class is a model of CCellularGridSpaceND. It represents the cubical grid as a cell complex...