DGtal  0.6.devel
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
testLatticePolytope2D.cpp
1 
30 
31 #include <iostream>
32 
33 //#define DEBUG_LatticePolytope2D
34 
35 #include "DGtal/base/Common.h"
36 #include "DGtal/arithmetic/LatticePolytope2D.h"
37 #include "DGtal/shapes/Shapes.h"
38 #include "DGtal/io/boards/Board2D.h"
39 
41 
42 using namespace std;
43 using namespace DGtal;
44 
46 // Functions for testing class LatticePolytope2D.
48 
49 template <typename LatticePolytope2D>
50 bool
51 checkCut( LatticePolytope2D & cip,
52  typename LatticePolytope2D::HalfSpace hs )
53 {
54  trace.beginBlock ( "Check cut, see <cip.eps> and <cip2.eps>" );
55  typedef typename LatticePolytope2D::Space Space;
56  typedef typename LatticePolytope2D::Domain Domain;
58  typedef typename DigitalSet::ConstIterator ConstIterator;
59 
60  Domain d = cip.boundingBoxDomain();
61 #ifdef DEBUG_LatticePolytope2D
62  Board2D board;
63  board << SetMode( d.className(), "Grid" ) << d;
64  DigitalSet aSet( d );
66  Color col1( 100, 100, 255 );
67  Color col2( 180, 180, 255 );
68  board << CustomStyle( aSet.className(), new CustomColors( col1, col2 ) )
69  << aSet;
70  board << SetMode( cip.className(), "Transparent" ) << cip;
71  board.saveEPS( "cip.eps" );
72  board.clear();
73 #endif
74 
75  DigitalSet cipSet( d );
76  DigitalSet cipSet2( d );
77  cip.getIncludedDigitalPoints( cipSet );
78  cip.cut( hs );
79  cip.getIncludedDigitalPoints( cipSet2 );
80 #ifdef DEBUG_LatticePolytope2D
81  board << SetMode( d.className(), "Grid" ) << d;
82  board << CustomStyle( aSet.className(), new CustomColors( col1, col2 ) )
83  << cipSet;
84  board << CustomStyle( aSet.className(), new CustomColors( Color( 255, 180, 20 ), Color( 200, 170, 0 ) ) )
85  << cipSet2;
86  board << SetMode( cip.className(), "Transparent" ) << cip;
87  board.saveEPS( "cip2.eps" );
88 #endif
89 
90  unsigned int nbok = 0;
91  unsigned int nb = 0;
92  for ( ConstIterator it = cipSet2.begin(), it_end = cipSet2.end();
93  it != it_end; ++it )
94  {
95  nbok += ( cipSet.find( *it ) != cipSet.end() ) ? 1 : 0;
96  ++nb;
97  nbok += hs( *it );
98  ++nb;
99  }
100  for ( ConstIterator it = cipSet.begin(), it_end = cipSet.end();
101  it != it_end; ++it )
102  {
103  if ( cipSet2.find( *it ) == cipSet2.end() )
104  nbok += ! hs( *it );
105  else
106  nbok += hs( *it );
107  ++nb;
108  }
109  trace.info() << "(" << nbok << "/" << nb << ")"
110  << " cip.size()=" << cip.size()
111  << " #before=" << cipSet.size()
112  << " #after=" << cipSet2.size()
113  << std::endl;
114  trace.endBlock();
115  return nbok == nb;
116 }
117 
122 template <typename Space>
123 bool testLatticePolytope2D()
124 {
125  unsigned int nbok = 0;
126  unsigned int nb = 0;
127 
128  trace.beginBlock ( "Testing block LatticePolytope2D area and centroid" );
129  typedef typename Space::Point Point;
130  typedef typename Space::Vector Vector;
131  typedef typename Space::Integer Integer;
132  typedef LatticePolytope2D<Space> CIP;
133  typedef typename CIP::Point3I Point3I;
134  typedef typename CIP::Domain Domain;
135  typedef typename CIP::HalfSpace HalfSpace;
136  typedef typename CIP::Iterator Iterator;
137  typedef typename CIP::SizeCouple SizeCouple;
138  typedef typename DigitalSetSelector< Domain, BIG_DS+HIGH_BEL_DS >::Type DigitalSet;
139 
140  CIP cip;
141  cip.pushBack( Point( 0, 0 ) );
142  cip.pushBack( Point( 5, 0 ) );
143  cip.pushBack( Point( 0, 3 ) );
144  Integer area2 = cip.twiceArea();
145  trace.info() << "- 2*area = " << area2 << std::endl;
146  ++nb, nbok += ( area2 == 15 ) ? 1 : 0;
147  trace.info() << "(" << nbok << "/" << nb << ") "
148  << "2*area == 15" << std::endl;
149  Point3I c = cip.centroid( area2 );
150  trace.info() << "- centroid = " << c << std::endl;
151  ++nb, nbok += ( c == Point3I( 75, 45, 45 ) ) ? 1 : 0;
152  trace.info() << "(" << nbok << "/" << nb << ") "
153  << "centroid == [75,45,45]" << std::endl;
154  Domain d = cip.boundingBoxDomain();
155  trace.info() << "- domain = " << d << std::endl;
156  trace.endBlock();
157 
158  trace.beginBlock ( "Output LatticePolytope2D in <cip.eps>" );
159  cip.pushBack( Point( -4, 4 ) );
160  cip.pushBack( Point( -7, 2 ) );
161  cip.pushBack( Point( -5, 1 ) );
162  Board2D board;
163  d = cip.boundingBoxDomain();
164  board << SetMode( d.className(), "Grid" ) << d;
165  DigitalSet aSet( d );
166  HalfSpace h( Vector( 1, 3 ), 8 );
167  //HalfSpace h = cip.halfSpace( ++cip.begin() );
169  Color col1( 100, 100, 255 );
170  Color col2( 180, 180, 255 );
171  board << CustomStyle( aSet.className(), new CustomColors( col1, col2 ) )
172  << aSet;
173  cip.getIncludedDigitalPoints( aSet );
174  board << CustomStyle( aSet.className(), new CustomColors( Color( 255, 180, 20 ), Color( 200, 130, 0 ) ) )
175  << aSet;
176  board << SetMode( cip.className(), "Transparent" ) << cip;
177  Iterator itA1, itB2;
178  SizeCouple nbs = cip.findCut( itA1, itB2, h );
179  Iterator itB1 = itA1; ++itB1;
180  if ( itB1 == cip.end() ) itB1 = cip.begin();
181  Iterator itA2 = itB2; ++itA2;
182  if ( itA2 == cip.end() ) itA2 = cip.begin();
183  Color col3( 0, 255, 0 );
184  Color col4( 255, 0, 0 );
185  board << CustomStyle( Point().className(), new CustomColors( col3, col3 ) )
186  << *itA1 << *itA2;
187  board << CustomStyle( Point().className(), new CustomColors( col4, col4 ) )
188  << *itB1 << *itB2;
189  board.saveEPS( "cip.eps" );
190  board.saveSVG( "cip.svg" );
191  trace.endBlock();
192 
193  trace.beginBlock ( "Output cut LatticePolytope2D in <cip2.eps>" );
194  board.clear();
195  board << SetMode( d.className(), "Grid" ) << d;
196  board << SetMode( cip.className(), "Transparent" ) << cip;
197  bool wasCut = cip.cut( h );
198  board << SetMode( cip.className(), "Filled" ) << cip;
199  board.saveEPS( "cip2.eps" );
200  board.saveSVG( "cip2.svg" );
201  trace.endBlock();
202 
203  checkCut( cip, HalfSpace( Vector( -2, 3 ), 4 ) );
204 
205  return nbok == nb;
206 }
207 
208 int myRandom( int nb )
209 {
210  return random() % nb;
211 }
212 
217 template <typename Space>
218 bool exhaustiveTestLatticePolytope2D()
219 {
220  unsigned int nbok = 0;
221  unsigned int nb = 0;
222 
223  trace.beginBlock ( "Testing block LatticePolytope2D area and centroid" );
224  typedef typename Space::Point Point;
225  typedef typename Space::Vector Vector;
226  typedef typename Space::Integer Integer;
227  typedef LatticePolytope2D<Space> CIP;
228  typedef typename CIP::Point3I Point3I;
229  typedef typename CIP::Domain Domain;
230  typedef typename CIP::HalfSpace HalfSpace;
231  typedef typename CIP::ConstIterator ConstIterator;
232  typedef typename CIP::SizeCouple SizeCouple;
233  typedef typename DigitalSetSelector< Domain, BIG_DS+HIGH_BEL_DS >::Type DigitalSet;
234 
235  CIP cip;
236  cip.pushBack( Point( 0, 0 ) );
237  cip.pushBack( Point( 8, -3 ) );
238  cip.pushBack( Point( 17, 2 ) );
239  cip.pushBack( Point( 21, 13 ) );
240  cip.pushBack( Point( 13, 19 ) );
241  cip.pushBack( Point( 6, 17 ) );
242  cip.pushBack( Point( -3, 6 ) );
243  Integer area2 = cip.twiceArea();
244  trace.info() << "- 2*area = " << area2 << std::endl;
245  Point3I c = cip.centroid( area2 );
246  trace.info() << "- centroid = " << c << std::endl;
247  Domain d = cip.boundingBoxDomain();
248  trace.info() << "- domain = " << d << std::endl;
249  trace.endBlock();
250 
251  trace.beginBlock ( "Check cuts..." );
252  for ( unsigned int j = 0; j < 100; ++j )
253  {
254  CIP cip2 = cip;
255  int x = 0;
256  int y = 0;
257  while ( ( x == 0 ) && ( y == 0 ) )
258  {
259  x = myRandom( 63 ) - 31;
260  y = myRandom( 63 ) - 31;
261  }
262  int g = IntegerComputer<int>::staticGcd( x , y );
263  x /= g; y /= g;
264  int c = myRandom( 4 ) *x + myRandom( 4 ) * y + myRandom( 40 ) + 40;
265  for ( unsigned int i = 0; i < 10; ++i, c -= myRandom( 40 ) )
266  {
267  HalfSpace h( Vector( x, y ), c );
268  trace.info() << "[" << j << " size=" << cip2.size() << "]"
269  << " cut by (" << x << "," << y << ")," << c << std::endl;
270  ++nb, nbok += checkCut( cip2, h ) ? 1 : 0;
271  trace.info() << "(" << nbok << "/" << nb << ") cuts" << std::endl;
272  //std::cerr << " " << cip2.size() << flush;
273  if ( nb != nbok ) break;
274  }
275  std::cerr << std::endl;
276  if ( nb != nbok ) break;
277  }
278  trace.endBlock();
279  return nbok == nb;
280 }
281 
286 template <typename Space>
287 bool specificTestLatticePolytope2D()
288 {
289  unsigned int nbok = 0;
290  unsigned int nb = 0;
291 
292  typedef typename Space::Point Point;
293  typedef typename Space::Vector Vector;
294  typedef typename Space::Integer Integer;
295  typedef LatticePolytope2D<Space> CIP;
296  typedef typename CIP::Point3I Point3I;
297  typedef typename CIP::Domain Domain;
298  typedef typename CIP::HalfSpace HalfSpace;
299  typedef typename CIP::Iterator Iterator;
300  typedef typename CIP::ConstIterator ConstIterator;
301  typedef typename CIP::SizeCouple SizeCouple;
302  typedef typename DigitalSetSelector< Domain, BIG_DS+HIGH_BEL_DS >::Type DigitalSet;
303 
304  CIP cip;
305  cip.pushBack( Point( 0, 0 ) );
306  cip.pushBack( Point( 8, -3 ) );
307  cip.pushBack( Point( 17, 2 ) );
308  cip.pushBack( Point( 21, 13 ) );
309  cip.pushBack( Point( 13, 19 ) );
310  cip.pushBack( Point( 6, 17 ) );
311  cip.pushBack( Point( -3, 6 ) );
312  HalfSpace hs1( Vector( -5, 11 ), 8 ); //12
313  HalfSpace hs2( Vector( -3, -11 ), 8 ); //-37
314  ++nb, nbok += checkCut( cip, hs1 ) ? 1 : 0;
315 
316  Domain d = cip.boundingBoxDomain();
317  Board2D board;
318  board << SetMode( d.className(), "Grid" ) << d;
319  DigitalSet aSet( d );
320  ConstIterator itv = cip.begin();
321  ++itv; ++itv; ++itv; ++itv; ++itv; ++itv; ++itv;
322  HalfSpace hs( cip.halfSpace( itv ) );
324  Color col1( 100, 180, 100 );
325  Color col2( 130, 200, 130 );
326  board << CustomStyle( aSet.className(), new CustomColors( col1, col2 ) )
327  << aSet;
328  board << SetMode( cip.className(), "Transparent" ) << cip;
329  Point p( 0, 0 );
330  board << CustomStyle( p.className(), new CustomColors( Color::Red, Color::Red ) )
331  << p;
332  board.saveEPS( "cip3.eps" );
333  board.clear();
334 
335  ++nb, nbok += checkCut( cip, hs2 ) ? 1 : 0;
336  return nbok == nb;
337 }
338 
339 template <typename Space>
340 bool
341 checkOutputConvexHullBorder()
342 {
343  typedef typename Space::Point Point;
344  typedef typename Space::Vector Vector;
345  typedef typename Space::Integer Integer;
346  typedef LatticePolytope2D<Space> CIP;
347  typedef typename CIP::Point3I Point3I;
348  typedef typename CIP::Domain Domain;
349  typedef typename CIP::HalfSpace HalfSpace;
350  typedef typename DigitalSetSelector< Domain, BIG_DS+HIGH_BEL_DS >::Type DigitalSet;
352 
353  BOOST_CONCEPT_ASSERT(( boost::Container< CIP > ));
354 
355  Point c1( 0, 0 );
356  Point c3( 30, 0 );
357  Point A( -10, -2 );
358  Point B( 40, 27 );
359  Point In( 1, 0 );
360  HalfSpace hs1( c1, c1 + Vector( 0, 1 ), In, ic );
361  HalfSpace hs2( A, B, In, ic );
362  HalfSpace hs3( c3, c3 + Vector( 0, 1 ), In, ic );
363  CIP cip;
364  cip.pushBack( c1 );
365  cip.computeConvexHullBorder( std::back_inserter( cip ), c1, c3, hs1, hs2, hs3 );
366  cip.pushBack( c3 );
367 
368  Domain d = cip.boundingBoxDomain();
369  Board2D board;
370  board << SetMode( d.className(), "Grid" ) << d;
371  DigitalSet aSet( d );
373  Color col1( 100, 180, 100 );
374  Color col2( 220, 250, 220 );
375  board << CustomStyle( aSet.className(), new CustomColors( col1, col2 ) )
376  << aSet;
377  board << SetMode( cip.className(), "Transparent" ) << cip;
378  // Point p( 0, 0 );
379  // board << CustomStyle( p.className(), new CustomColors( Color::Red, Color::Red ) )
380  // << p;
381  board.saveEPS( "cip4.eps" );
382  board.clear();
383  return true;
384 }
385 
387 // Standard services - public :
388 
389 int main( int, char** )
390 {
391  trace.beginBlock ( "Testing class LatticePolytope2D" );
392 
394  typedef SpaceND<2, DGtal::BigInteger> Z2I;
395  bool res = testLatticePolytope2D<Z2>()
396  && testLatticePolytope2D<Z2I>()
397  && exhaustiveTestLatticePolytope2D<Z2>()
398  && checkOutputConvexHullBorder<Z2>();
399  //&& specificTestLatticePolytope2D<Z2>();
400  //&& exhaustiveTestLatticePolytope2D<Z2I>();
401  trace.emphase() << ( res ? "Passed." : "Error." ) << endl;
402  trace.endBlock();
403  return res ? 0 : 1;
404 }
405 // //