DGtal  1.4.beta
testChordNaivePlaneComputer.cpp
Go to the documentation of this file.
1 
31 #include <cstdlib>
32 #include <iostream>
33 #include "DGtal/base/Common.h"
34 #include "DGtal/helpers/StdDefs.h"
35 #include "DGtal/kernel/CPointPredicate.h"
36 #include "DGtal/geometry/surfaces/CAdditivePrimitiveComputer.h"
37 #include "DGtal/geometry/surfaces/ChordNaivePlaneComputer.h"
38 #include "DGtal/geometry/surfaces/ChordGenericNaivePlaneComputer.h"
40 
41 using namespace std;
42 using namespace DGtal;
43 using namespace DGtal::concepts;
44 
46 // Functions for testing class ChordNaivePlaneComputer.
48 
49 template <typename Integer>
50 Integer getRandomInteger( Integer first, Integer after_last )
51 {
52  Integer r = (Integer) rand();
53  return ( r % (after_last - first) ) + first;
54 }
55 
59 template <typename Integer, typename NaivePlaneComputer>
60 bool
61 checkPlane( Integer a, Integer b, Integer c, Integer d,
62  int diameter, unsigned int nbtries )
63 {
64  typedef typename NaivePlaneComputer::Point Point;
65  typedef typename Point::Component PointInteger;
66  IntegerComputer<Integer> ic;
67  Integer absA = ic.abs( a );
68  Integer absB = ic.abs( b );
69  Integer absC = ic.abs( c );
70  Integer x, y, z;
71  Dimension axis;
72  if ( ( absA >= absB ) && ( absA >= absC ) )
73  axis = 0;
74  else if ( ( absB >= absA ) && ( absB >= absC ) )
75  axis = 1;
76  else
77  axis = 2;
78  Point p;
79  NaivePlaneComputer plane;
80  plane.init( axis, 1, 1 );
81  // Checks that points within the naive plane are correctly recognized.
82  unsigned int nb = 0;
83  unsigned int nbok = 0;
84  while ( nb != nbtries )
85  {
86  p[ 0 ] = getRandomInteger<PointInteger>( -diameter+1, diameter );
87  p[ 1 ] = getRandomInteger<PointInteger>( -diameter+1, diameter );
88  p[ 2 ] = getRandomInteger<PointInteger>( -diameter+1, diameter );
89  x = (Integer) p[ 0 ];
90  y = (Integer) p[ 1 ];
91  z = (Integer) p[ 2 ];
92  switch ( axis ) {
93  case 0: p[ 0 ] = (PointInteger)NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - b * y - c * z, a ) ); break;
94  case 1: p[ 1 ] = (PointInteger)NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - a * x - c * z, b ) ); break;
95  case 2: p[ 2 ] = (PointInteger)NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - a * x - b * y, c ) ); break;
96  }
97  bool ok_ext = plane.isExtendable( p ); // should be ok
98  bool ok = plane.extend( p ); // should be ok
99  ++nb; nbok += ok_ext ? 1 : 0;
100  ++nb; nbok += ok ? 1 : 0;
101  if ( ! ok )
102  {
103  std::cerr << "[ERROR] p=" << p << " NOT IN plane=" << plane << std::endl;
104  for ( typename NaivePlaneComputer::ConstIterator it = plane.begin(), itE = plane.end();
105  it != itE; ++it )
106  std::cerr << " " << *it;
107  std::cerr << endl;
108  std::cerr << "d <= a*x+b*y+c*z <= d+max(a,b,c)"
109  << d << " <= " << a << "*" << p[0]
110  << "+" << b << "*" << p[1]
111  << "+" << c << "*" << p[2]
112  << " = " << (a*p[0]+b*p[1]+c*p[2])
113  << std::endl;
114  break;
115  }
116  if ( ! ok_ext )
117  {
118  std::cerr << "[ERROR] p=" << p << " was NOT extendable IN plane=" << plane << std::endl;
119  break;
120  }
121  // else
122  // std::cerr << "[OK] p=" << p << " IN plane=" << plane << std::endl;
123  }
124 
125  // Checks that points outside the naive plane are correctly recognized as outliers.
126  while ( nb != (nbtries * 11 ) / 10 )
127  {
128  p[ 0 ] = getRandomInteger<PointInteger>( -diameter+1, diameter );
129  p[ 1 ] = getRandomInteger<PointInteger>( -diameter+1, diameter );
130  p[ 2 ] = getRandomInteger<PointInteger>( -diameter+1, diameter );
131  x = (Integer) p[ 0 ];
132  y = (Integer) p[ 1 ];
133  z = (Integer) p[ 2 ];
134  switch ( axis ) {
135  case 0: p[ 0 ] = (PointInteger)NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - b * y - c * z, a ) ); break;
136  case 1: p[ 1 ] = (PointInteger)NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - a * x - c * z, b ) ); break;
137  case 2: p[ 2 ] = (PointInteger)NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - a * x - b * y, c ) ); break;
138  }
139  PointInteger tmp = getRandomInteger<PointInteger>( 2, 5 )
140  * (2*getRandomInteger<PointInteger>( 0, 2 ) - 1 );
141  p[ axis ] += tmp;
142  bool ok_ext = ! plane.isExtendable( p ); // should *not* be ok
143  bool ok = ! plane.extend( p ); // should *not* be ok
144  ++nb; nbok += ok ? 1 : 0;
145  ++nb; nbok += ok_ext ? 1 : 0;
146  if ( ! ok )
147  {
148  std::cerr << "[ERROR] p=" << p << " IN plane=" << plane << std::endl;
149  break;
150  }
151  if ( ! ok_ext )
152  {
153  std::cerr << "[ERROR] p=" << p << " was extendable IN plane=" << plane << std::endl;
154  break;
155  }
156  // else
157  // std::cerr << "[OK] p=" << p << " IN plane=" << plane << std::endl;
158  }
159  return nb == nbok;
160 }
161 
165 template <typename Integer, typename NaivePlaneComputer>
166 bool
167 checkPlaneGroupExtension( Integer a, Integer b, Integer c, Integer d,
168  int diameter, unsigned int nbtries )
169 {
170  typedef typename NaivePlaneComputer::Point Point;
171  typedef typename Point::Component PointInteger;
172  IntegerComputer<Integer> ic;
173  Integer absA = ic.abs( a );
174  Integer absB = ic.abs( b );
175  Integer absC = ic.abs( c );
176  Integer x, y, z;
177  Dimension axis;
178  if ( ( absA >= absB ) && ( absA >= absC ) )
179  axis = 0;
180  else if ( ( absB >= absA ) && ( absB >= absC ) )
181  axis = 1;
182  else
183  axis = 2;
184  Point p;
185  NaivePlaneComputer plane;
186  plane.init( axis, 1, 1 );
187  // Checks that points within the naive plane are correctly recognized.
188  unsigned int nb = 0;
189  unsigned int nbok = 0;
190  while ( nb < nbtries )
191  {
192  std::vector<Point> points( 5 );
193  for ( unsigned int i = 0; i < 5; ++i )
194  {
195  Point & pp = points[ i ];
196  pp[ 0 ] = getRandomInteger<PointInteger>( -diameter+1, diameter );
197  pp[ 1 ] = getRandomInteger<PointInteger>( -diameter+1, diameter );
198  pp[ 2 ] = getRandomInteger<PointInteger>( -diameter+1, diameter );
199  x = (Integer) pp[ 0 ];
200  y = (Integer) pp[ 1 ];
201  z = (Integer) pp[ 2 ];
202  switch ( axis ) {
203  case 0: pp[ 0 ] = (PointInteger)NumberTraits<Integer>::castToInt64_t
204  ( ic.ceilDiv( d - b * y - c * z, a ) ); break;
205  case 1: pp[ 1 ] = (PointInteger)NumberTraits<Integer>::castToInt64_t
206  ( ic.ceilDiv( d - a * x - c * z, b ) ); break;
207  case 2: pp[ 2 ] = (PointInteger)NumberTraits<Integer>::castToInt64_t
208  ( ic.ceilDiv( d - a * x - b * y, c ) ); break;
209  }
210  }
211  bool ok_ext = plane.isExtendable( points.begin(), points.end() ); // should be ok
212  bool ok = plane.extend( points.begin(), points.end() ); // should be ok
213  ++nb; nbok += ok_ext ? 1 : 0;
214  ++nb; nbok += ok ? 1 : 0;
215  if ( ! ok )
216  {
217  std::cerr << "[ERROR] p=" << points[ 0 ] << " NOT IN plane=" << plane << std::endl;
218  for ( typename NaivePlaneComputer::ConstIterator it = plane.begin(), itE = plane.end();
219  it != itE; ++it )
220  std::cerr << " " << *it;
221  std::cerr << endl;
222  std::cerr << "d <= a*x+b*y+c*z <= d+max(a,b,c)"
223  << d << " <= " << a << "*" << p[0]
224  << "+" << b << "*" << p[1]
225  << "+" << c << "*" << p[2]
226  << " = " << (a*p[0]+b*p[1]+c*p[2])
227  << std::endl;
228  break;
229  }
230  if ( ! ok_ext )
231  {
232  std::cerr << "[ERROR] p=" << p << " was NOT extendable IN plane=" << plane << std::endl;
233  break;
234  }
235  // else
236  // std::cerr << "[OK] p=" << p << " IN plane=" << plane << std::endl;
237  }
238 
239  // Checks that points outside the naive plane are correctly recognized as outliers.
240  while ( nb < (nbtries * 11 ) / 10 )
241  {
242  p[ 0 ] = getRandomInteger<PointInteger>( -diameter+1, diameter );
243  p[ 1 ] = getRandomInteger<PointInteger>( -diameter+1, diameter );
244  p[ 2 ] = getRandomInteger<PointInteger>( -diameter+1, diameter );
245  x = (Integer) p[ 0 ];
246  y = (Integer) p[ 1 ];
247  z = (Integer) p[ 2 ];
248  switch ( axis ) {
249  case 0: p[ 0 ] = (PointInteger)NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - b * y - c * z, a ) ); break;
250  case 1: p[ 1 ] = (PointInteger)NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - a * x - c * z, b ) ); break;
251  case 2: p[ 2 ] = (PointInteger)NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - a * x - b * y, c ) ); break;
252  }
253  PointInteger tmp = getRandomInteger<PointInteger>( 2, 5 )
254  * (2*getRandomInteger<PointInteger>( 0, 2 ) - 1 );
255  p[ axis ] += tmp;
256  bool ok_ext = ! plane.isExtendable( p ); // should *not* be ok
257  bool ok = ! plane.extend( p ); // should *not* be ok
258  ++nb; nbok += ok ? 1 : 0;
259  ++nb; nbok += ok_ext ? 1 : 0;
260  if ( ! ok )
261  {
262  std::cerr << "[ERROR] p=" << p << " IN plane=" << plane << std::endl;
263  break;
264  }
265  if ( ! ok_ext )
266  {
267  std::cerr << "[ERROR] p=" << p << " was extendable IN plane=" << plane << std::endl;
268  break;
269  }
270  // else
271  // std::cerr << "[OK] p=" << p << " IN plane=" << plane << std::endl;
272  }
273  return nb == nbok;
274 }
275 
276 
277 
281 template <typename Integer, typename GenericNaivePlaneComputer>
282 bool
283 checkGenericPlane( Integer a, Integer b, Integer c, Integer d,
284  int diameter, unsigned int nbtries )
285 {
286  typedef typename GenericNaivePlaneComputer::Point Point;
287  typedef typename Point::Component PointInteger;
288  IntegerComputer<Integer> ic;
289  Integer absA = ic.abs( a );
290  Integer absB = ic.abs( b );
291  Integer absC = ic.abs( c );
292  Integer x, y, z;
293  Dimension axis;
294  if ( ( absA >= absB ) && ( absA >= absC ) )
295  axis = 0;
296  else if ( ( absB >= absA ) && ( absB >= absC ) )
297  axis = 1;
298  else
299  axis = 2;
300  Point p;
301  GenericNaivePlaneComputer plane;
302  plane.init( 1, 1 );
303  // Checks that points within the naive plane are correctly recognized.
304  unsigned int nb = 0;
305  unsigned int nbok = 0;
306  while ( nb != nbtries )
307  {
308  p[ 0 ] = getRandomInteger<PointInteger>( -diameter+1, diameter );
309  p[ 1 ] = getRandomInteger<PointInteger>( -diameter+1, diameter );
310  p[ 2 ] = getRandomInteger<PointInteger>( -diameter+1, diameter );
311  x = (Integer) p[ 0 ];
312  y = (Integer) p[ 1 ];
313  z = (Integer) p[ 2 ];
314  switch ( axis ) {
315  case 0: p[ 0 ] = (PointInteger)NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - b * y - c * z, a ) ); break;
316  case 1: p[ 1 ] = (PointInteger)NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - a * x - c * z, b ) ); break;
317  case 2: p[ 2 ] = (PointInteger)NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - a * x - b * y, c ) ); break;
318  }
319  bool ok_ext = plane.isExtendable( p ); // should be ok
320  bool ok = plane.extend( p ); // should be ok
321  ++nb; nbok += ok_ext ? 1 : 0;
322  ++nb; nbok += ok ? 1 : 0;
323  if ( ! ok )
324  {
325  std::cerr << "[ERROR] p=" << p << " NOT IN plane=" << plane << std::endl;
326  for ( typename GenericNaivePlaneComputer::ConstIterator it = plane.begin(), itE = plane.end();
327  it != itE; ++it )
328  std::cerr << " " << *it;
329  std::cerr << endl;
330  std::cerr << "d <= a*x+b*y+c*z <= d+max(a,b,c)"
331  << d << " <= " << a << "*" << p[0]
332  << "+" << b << "*" << p[1]
333  << "+" << c << "*" << p[2]
334  << " = " << (a*p[0]+b*p[1]+c*p[2])
335  << std::endl;
336  break;
337  }
338  if ( ! ok_ext )
339  {
340  std::cerr << "[ERROR] p=" << p << " was NOT extendable IN plane=" << plane << std::endl;
341  break;
342  }
343  // else
344  // std::cerr << "[OK] p=" << p << " IN plane=" << plane << std::endl;
345  }
346 
347  // Checks that points outside the naive plane are correctly recognized as outliers.
348  while ( nb != (nbtries * 11 ) / 10 )
349  {
350  p[ 0 ] = getRandomInteger<PointInteger>( -diameter+1, diameter );
351  p[ 1 ] = getRandomInteger<PointInteger>( -diameter+1, diameter );
352  p[ 2 ] = getRandomInteger<PointInteger>( -diameter+1, diameter );
353  x = (Integer) p[ 0 ];
354  y = (Integer) p[ 1 ];
355  z = (Integer) p[ 2 ];
356  switch ( axis ) {
357  case 0: p[ 0 ] = (PointInteger)NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - b * y - c * z, a ) ); break;
358  case 1: p[ 1 ] = (PointInteger)NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - a * x - c * z, b ) ); break;
359  case 2: p[ 2 ] = (PointInteger)NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - a * x - b * y, c ) ); break;
360  }
361  PointInteger tmp = getRandomInteger<PointInteger>( 2, 5 )
362  * (2*getRandomInteger<PointInteger>( 0, 2 ) - 1 );
363  p[ axis ] += tmp;
364  bool ok_ext = ! plane.isExtendable( p ); // should *not* be ok
365  bool ok = ! plane.extend( p ); // should *not* be ok
366  ++nb; nbok += ok ? 1 : 0;
367  ++nb; nbok += ok_ext ? 1 : 0;
368  if ( ! ok )
369  {
370  std::cerr << "[ERROR] p=" << p << " IN plane=" << plane << std::endl;
371  break;
372  }
373  if ( ! ok_ext )
374  {
375  std::cerr << "[ERROR] p=" << p << " was extendable IN plane=" << plane << std::endl;
376  break;
377  }
378  // else
379  // std::cerr << "[OK] p=" << p << " IN plane=" << plane << std::endl;
380  }
381  std::cerr << "plane = " << plane << std::endl;
382  return nb == nbok;
383 }
384 
385 
386 template <typename Integer, typename NaivePlaneComputer>
387 bool
388 checkPlanes( unsigned int nbplanes, int diameter, unsigned int nbtries )
389 {
390  //using namespace Z3i;
391  //typedef ChordNaivePlaneComputer<Z3, Integer> NaivePlaneComputer;
392  unsigned int nb = 0;
393  unsigned int nbok = 0;
394  for ( unsigned int nbp = 0; nbp < nbplanes; ++nbp )
395  {
396  Integer a = getRandomInteger<Integer>( (Integer) 0, (Integer) diameter / 2 );
397  Integer b = getRandomInteger<Integer>( (Integer) 0, (Integer) diameter / 2 );
398  Integer c = getRandomInteger<Integer>( (Integer) 0, (Integer) diameter / 2 );
399  Integer d = getRandomInteger<Integer>( (Integer) 0, (Integer) diameter / 2 );
400  if ( ( a != 0 ) || ( b != 0 ) || ( c != 0 ) )
401  {
402  ++nb; nbok += checkPlane<Integer, NaivePlaneComputer>( a, b, c, d, diameter, nbtries ) ? 1 : 0;
403  if ( nb != nbok )
404  {
405  std::cerr << "[ERROR] (Simple extension) for plane " << a << " * x + "
406  << b << " * y + " << c << " * z = " << d << std::endl;
407  break;
408  }
409  ++nb; nbok += checkPlaneGroupExtension<Integer, NaivePlaneComputer>( a, b, c, d, diameter, nbtries ) ? 1 : 0;
410  if ( nb != nbok )
411  {
412  std::cerr << "[ERROR] (Group extension) for plane " << a << " * x + "
413  << b << " * y + " << c << " * z = " << d << std::endl;
414  break;
415  }
416  }
417  }
418  return nb == nbok;
419 }
420 
424 template <typename Integer, typename NaivePlaneComputer>
425 bool
426 checkWidth( Integer a, Integer b, Integer c, Integer d,
427  int diameter, unsigned int nbtries )
428 {
429  typedef typename NaivePlaneComputer::Point Point;
430  typedef typename NaivePlaneComputer::InternalScalar InternalScalar;
431  typedef typename Point::Component PointInteger;
432  IntegerComputer<Integer> ic;
433  Integer absA = ic.abs( a );
434  Integer absB = ic.abs( b );
435  Integer absC = ic.abs( c );
436  Integer x, y, z;
437  Dimension axis;
438  if ( ( absA >= absB ) && ( absA >= absC ) )
439  axis = 0;
440  else if ( ( absB >= absA ) && ( absB >= absC ) )
441  axis = 1;
442  else
443  axis = 2;
444  // Checks that points within the naive plane are correctly recognized.
445  unsigned int nb = 0;
446  unsigned int nbok = 0;
447  std::vector<Point> points( nbtries );
448  for ( unsigned int i = 0; i != nbtries; ++i )
449  {
450  Point & p = points[ i ];
451  p[ 0 ] = getRandomInteger<Integer>( -diameter+1, diameter );
452  p[ 1 ] = getRandomInteger<Integer>( -diameter+1, diameter );
453  p[ 2 ] = getRandomInteger<Integer>( -diameter+1, diameter );
454  x = (Integer) p[ 0 ];
455  y = (Integer) p[ 1 ];
456  z = (Integer) p[ 2 ];
457  switch ( axis ) {
458  case 0: p[ 0 ] = (PointInteger)NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - b * y - c * z, a ) ); break;
459  case 1: p[ 1 ] = (PointInteger)NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - a * x - c * z, b ) ); break;
460  case 2: p[ 2 ] = (PointInteger)NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - a * x - b * y, c ) ); break;
461  }
462  }
463  trace.beginBlock( "Computing axis width." );
464  trace.info() << "- plane is "
465  << d << " <= " << a << "*x"
466  << "+" << b << "*y"
467  << "+" << c << "*z"
468  << " <= d + max(|a|,|b|,|c|)"
469  << std::endl;
470  trace.info() << "- " << points.size() << " points tested in diameter " << diameter
471  << std::endl;
472  double min = -1.0;
473  for ( unsigned int i = 0; i < 3; ++i )
474  {
475  std::pair<InternalScalar, InternalScalar> width
476  = NaivePlaneComputer::computeAxisWidth( i, points.begin(), points.end() );
477  double wn = NumberTraits<InternalScalar>::castToDouble( width.first );
478  double wd = NumberTraits<InternalScalar>::castToDouble( width.second );
479  trace.info() << " (" << i << ") width=" << (wn/wd) << std::endl;
480  if ( min < 0.0 ) min = wn/wd;
481  else if ( wn/wd < min ) min = wn/wd;
482  }
483  ++nb; nbok += (min < 1.0 ) ? 1 : 0;
484  trace.info() << "(" << nbok << "/" << nb << ") min width = " << min
485  << " < 1.0" << std::endl;
486  ++nb; nbok += (0.9 < min ) ? 1 : 0;
487  trace.info() << "(" << nbok << "/" << nb << ") min width = " << min
488  << " > 0.9" << std::endl;
489  trace.endBlock();
490  return nb == nbok;
491 }
492 
493 template <typename Integer, typename NaivePlaneComputer>
494 bool
495 checkWidths( unsigned int nbplanes, int diameter, unsigned int nbtries )
496 {
497  //using namespace Z3i;
498  //typedef ChordNaivePlaneComputer<Z3, Integer> NaivePlaneComputer;
499  unsigned int nb = 0;
500  unsigned int nbok = 0;
501  for ( unsigned int nbp = 0; nbp < nbplanes; ++nbp )
502  {
503  Integer a = getRandomInteger<Integer>( (Integer) 0, (Integer) diameter / 2 );
504  Integer b = getRandomInteger<Integer>( (Integer) 0, (Integer) diameter / 2 );
505  Integer c = getRandomInteger<Integer>( (Integer) 0, (Integer) diameter / 2 );
506  Integer d = getRandomInteger<Integer>( (Integer) 0, (Integer) diameter / 2 );
507  if ( ( a != 0 ) || ( b != 0 ) || ( c != 0 ) )
508  {
509  ++nb; nbok += checkWidth<Integer, NaivePlaneComputer>( a, b, c, d, diameter, nbtries ) ? 1 : 0;
510  if ( nb != nbok )
511  {
512  std::cerr << "[ERROR] (checkWidth) for plane " << a << " * x + "
513  << b << " * y + " << c << " * z = " << d << std::endl;
514  break;
515  }
516  }
517  }
518  return nb == nbok;
519 }
520 
521 
526 bool testChordNaivePlaneComputer()
527 {
528  unsigned int nbok = 0;
529  unsigned int nb = 0;
530  typedef DGtal::int64_t Integer;
531  typedef DGtal::Z3i::Z3 Space;
532  typedef DGtal::Z3i::Point Point;
533  typedef ChordNaivePlaneComputer<Space, Point, Integer> NaivePlaneComputer;
534  typedef ChordGenericNaivePlaneComputer<Space, Point, Integer> GenericNaivePlaneComputer;
535 
536  BOOST_CONCEPT_ASSERT(( CAdditivePrimitiveComputer< NaivePlaneComputer > ));
537  BOOST_CONCEPT_ASSERT(( CAdditivePrimitiveComputer< GenericNaivePlaneComputer > ));
538  BOOST_CONCEPT_ASSERT(( boost::ForwardContainer< NaivePlaneComputer > ));
539  BOOST_CONCEPT_ASSERT(( boost::ForwardContainer< GenericNaivePlaneComputer > ));
540  BOOST_CONCEPT_ASSERT(( CPointPredicate< NaivePlaneComputer::Primitive > ));
541  BOOST_CONCEPT_ASSERT(( CPointPredicate< GenericNaivePlaneComputer::Primitive > ));
542 
543  trace.beginBlock ( "Testing block: ChordNaivePlaneComputer instantiation." );
544  NaivePlaneComputer plane;
545  Point pt0( 0, 0, 0 );
546  plane.init( 2, 1, 1 );
547  bool pt0_inside = plane.extend( pt0 );
548  ++nb; nbok += pt0_inside == true ? 1 : 0;
549  trace.info() << "(" << nbok << "/" << nb << ") Plane=" << plane
550  << std::endl;
551  Point pt1( Point( 2, 0, 0 ) );
552  bool pt1_inside = plane.extend( pt1 );
553  ++nb; nbok += pt1_inside == true ? 1 : 0;
554  trace.info() << "(" << nbok << "/" << nb << ") add " << pt1
555  << " Plane=" << plane << std::endl;
556  Point pt2( Point( 0, 2, 2 ) );
557  bool pt2_inside = plane.extend( pt2 );
558  ++nb; nbok += pt2_inside == true ? 1 : 0;
559  trace.info() << "(" << nbok << "/" << nb << ") add " << pt2
560  << " Plane=" << plane << std::endl;
561 
562  Point pt3( Point( 1, 1, 1 ) );
563  bool pt3_inside = plane.extend( pt3 );
564  ++nb; nbok += pt3_inside == true ? 1 : 0;
565  trace.info() << "(" << nbok << "/" << nb << ") add " << pt3
566  << " Plane=" << plane << std::endl;
567 
568  Point pt4( Point( -10, -10, 10 ) );
569  bool pt4_inside = plane.extend( pt4 );
570  ++nb; nbok += pt4_inside == false ? 1 : 0;
571  trace.info() << "(" << nbok << "/" << nb << ") impossible add " << pt4
572  << " Plane=" << plane << std::endl;
573 
574  Point pt5 = pt2 + Point( 1, 0, 1 );
575  bool pt5_inside = plane.extend( pt5 );
576  ++nb; nbok += pt5_inside == true ? 1 : 0;
577  trace.info() << "(" << nbok << "/" << nb << ") add " << pt5
578  << " Plane=" << plane << std::endl;
579 
580  Point pt6 = pt5 + Point( 6, 0, 2 );
581  bool pt6_inside = plane.extend( pt6 );
582  ++nb; nbok += pt6_inside == true ? 1 : 0;
583  trace.info() << "(" << nbok << "/" << nb << ") add " << pt6
584  << " Plane=" << plane << std::endl;
585 
586  NaivePlaneComputer plane2;
587  plane2.init( 2, 1, 1 );
588  plane2.extend( Point( 10, 0, 0 ) );
589  plane2.extend( Point( 0, 8, 0 ) );
590  plane2.extend( Point( 0, 0, 6 ) );
591  trace.info() << "(" << nbok << "/" << nb << ") "
592  << " Plane2=" << plane2 << std::endl;
593 
594  ++nb; nbok += checkPlane<Integer,NaivePlaneComputer>( 11, 5, 19, 20, 100, 100 ) ? 1 : 0;
595  trace.info() << "(" << nbok << "/" << nb
596  << ") checkPlane<Integer,NaivePlaneComputer>( 11, 5, 19, 20, 100, 100 )"
597  << std::endl;
598 
599  ++nb; nbok += checkGenericPlane<Integer,GenericNaivePlaneComputer>( 11, 5, 19, 20, 100, 100 ) ? 1 : 0;
600  trace.info() << "(" << nbok << "/" << nb
601  << ") checkGenericPlane<Integer,GenericNaivePlaneComputer>( 11, 5, 19, 20, 100, 100 )"
602  << std::endl;
603  ++nb; nbok += checkGenericPlane<Integer,GenericNaivePlaneComputer>( 17, 33, 7, 10, 100, 100 ) ? 1 : 0;
604  trace.info() << "(" << nbok << "/" << nb
605  << ") checkGenericPlane<Integer,GenericNaivePlaneComputer>( 17, 33, 7, 10, 100, 100 )"
606  << std::endl;
607  ++nb; nbok += checkPlane<Integer,NaivePlaneComputer>( 15, 8, 13, 15, 100, 100 ) ? 1 : 0;
608  trace.info() << "(" << nbok << "/" << nb
609  << ") checkPlane<Integer,NaivePlaneComputer>( 15, 8, 13, 15, 100, 100 )"
610  << std::endl;
611  ++nb; nbok += checkGenericPlane<Integer,GenericNaivePlaneComputer>( 15, 8, 13, 15, 100, 100 ) ? 1 : 0;
612  trace.info() << "(" << nbok << "/" << nb
613  << ") checkGenericPlane<Integer,GenericNaivePlaneComputer>( 15, 8, 13, 15, 100, 100 )"
614  << std::endl;
615  trace.endBlock();
616 
617  {
618  trace.beginBlock ( "Testing block: ChordNaivePlaneComputer vertical instantiation." );
619  NaivePlaneComputer ppplane;
620  Point pppt0( 0, 0, 0 );
621  ppplane.init( 2, 5, 2 );
622  bool pppt0_inside = ppplane.extend( pppt0 );
623  ++nb; nbok += pppt0_inside == true ? 1 : 0;
624  trace.info() << "(" << nbok << "/" << nb << ") Plane=" << ppplane
625  << std::endl;
626  Point pppt1( 3, 2, 2 );
627  bool pppt1_inside = ppplane.extend( pppt1 );
628  ++nb; nbok += pppt1_inside == true ? 1 : 0;
629  trace.info() << "(" << nbok << "/" << nb << ") Plane=" << ppplane
630  << std::endl;
631  Point pppt2( 0, 0, 1 );
632  bool pppt2_inside = ppplane.extend( pppt2 );
633  ++nb; nbok += pppt2_inside == true ? 1 : 0;
634  trace.info() << "(" << nbok << "/" << nb << ") Plane=" << ppplane
635  << std::endl;
636  Point pppt3 = pppt1 + Point( 0, 0, 1 );
637  bool pppt3_inside = ppplane.extend( pppt3 );
638  ++nb; nbok += pppt3_inside == true ? 1 : 0;
639  trace.info() << "(" << nbok << "/" << nb << ") Plane=" << ppplane
640  << std::endl;
641  Point pppt4 = pppt3 + Point( 0, 0, 1 );
642  bool pppt4_inside = ppplane.extend( pt4 );
643  ++nb; nbok += pppt4_inside == true ? 1 : 0;
644  trace.info() << "(" << nbok << "/" << nb << ") Plane=" << ppplane
645  << std::endl;
646  trace.endBlock();
647  }
648 
649  {
650  trace.beginBlock ( "Testing block: ChordNaivePlaneComputer vertical instantiation 2." );
651  NaivePlaneComputer pplane;
652  pplane.init( 1, 1, 1 );
653  Point ppt0( -6, -3, 5 );
654  bool ppt0_inside = pplane.extend( ppt0 );
655  ++nb; nbok += ppt0_inside == true ? 1 : 0;
656  trace.info() << "(" << nbok << "/" << nb << ") Plane=" << pplane
657  << std::endl;
658  Point ppt1( 4, 4, -5 );
659  bool ppt1_inside = pplane.extend( ppt1 );
660  ++nb; nbok += ppt1_inside == true ? 1 : 0;
661  trace.info() << "(" << nbok << "/" << nb << ") Plane=" << pplane
662  << std::endl;
663  Point ppt2( -5, -2, 4 );
664  bool ppt2_inside = pplane.extend( ppt2 );
665  ++nb; nbok += ppt2_inside == true ? 1 : 0;
666  trace.info() << "(" << nbok << "/" << nb << ") Plane=" << pplane
667  << std::endl;
668  trace.endBlock();
669  }
670 
671  return nbok == nb;
672 }
673 
674 template <typename NaivePlaneComputer>
675 bool
676 checkManyPlanes( unsigned int diameter,
677  unsigned int nbplanes,
678  unsigned int nbpoints )
679 {
680  unsigned int nbok = 0;
681  unsigned int nb = 0;
682  typedef typename NaivePlaneComputer::InternalScalar Scalar;
683  stringstream ss (stringstream::out);
684  ss << "Testing block: Diameter is " << diameter << ". Check " << nbplanes << " planes with " << nbpoints << " points each.";
685  trace.beginBlock ( ss.str() );
686  ++nb; nbok += checkPlanes<Scalar,NaivePlaneComputer>( nbplanes, diameter, nbpoints ) ? 1 : 0;
687  trace.info() << "(" << nbok << "/" << nb
688  << ") checkPlanes<Scalar,NaivePlaneComputer>()"
689  << std::endl;
690  trace.endBlock();
691  return nbok == nb;
692 }
693 
694 template <typename GenericNaivePlaneComputer>
695 bool
696 checkExtendWithManyPoints( unsigned int diameter,
697  unsigned int nbplanes,
698  unsigned int nbpoints )
699 {
700  unsigned int nbok = 0;
701  unsigned int nb = 0;
702  typedef typename GenericNaivePlaneComputer::InternalScalar Integer;
703  typedef typename GenericNaivePlaneComputer::Point Point;
704  typedef typename Point::Coordinate PointInteger;
705  IntegerComputer<Integer> ic;
706 
707  trace.beginBlock( "checkExtendWithManyPoints" );
708  for ( unsigned int j = 0; j < nbplanes; ++j )
709  {
710  Integer a = getRandomInteger<Integer>( (Integer) 0, (Integer) diameter / 2 );
711  Integer b = getRandomInteger<Integer>( (Integer) 0, (Integer) diameter / 2 );
712  Integer c = getRandomInteger<Integer>( (Integer) 1, (Integer) diameter / 2 );
713  Integer d = getRandomInteger<Integer>( (Integer) 0, (Integer) diameter / 2 );
714  GenericNaivePlaneComputer plane;
715  Dimension axis;
716  if ( ( a >= b ) && ( a >= c ) ) axis = 0;
717  else if ( ( b >= a ) && ( b >= c ) ) axis = 1;
718  else axis = 2;
719  plane.init( 1, 1 );
720 
721  std::vector<Point> pts;
722  for ( unsigned int i = 0; i < nbpoints; ++i )
723  {
724  Point p;
725  p[ 0 ] = getRandomInteger<PointInteger>( -diameter+1, diameter );
726  p[ 1 ] = getRandomInteger<PointInteger>( -diameter+1, diameter );
727  p[ 2 ] = getRandomInteger<PointInteger>( -diameter+1, diameter );
728  Integer x = (Integer) p[ 0 ];
729  Integer y = (Integer) p[ 1 ];
730  Integer z = (Integer) p[ 2 ];
731  switch( axis ) {
732  case 0: p[ 0 ] = (Integer)NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - b * y - c * z, a ) ); break;
733  case 1: p[ 1 ] = (Integer)NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - a * x - c * z, b ) ); break;
734  case 2: p[ 2 ] = (Integer)NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - a * x - b * y, c ) ); break;
735  }
736  pts.push_back( p );
737  }
738  ++nb; nbok += plane.isExtendable( pts.begin(), pts.end() ); // should be ok
739  trace.info() << "(" << nbok << "/" << nb
740  << ") plane.isExtendable( pts.begin(), pts.end() )"
741  << std::endl;
742  Point & any0 = pts[ getRandomInteger<int>( 0, pts.size() ) ];
743  pts.push_back( any0 + Point(1,0,0) );
744  Point & any1 = pts[ getRandomInteger<int>( 0, pts.size() ) ];
745  pts.push_back( any1 + Point(0,1,0) );
746  Point & any2 = pts[ getRandomInteger<int>( 0, pts.size() ) ];
747  pts.push_back( any2 + Point(0,0,1) );
748  bool check = ! plane.isExtendable( pts.begin(), pts.end() ); // should not be ok
749  ++nb; nbok += check ? 1 : 0;
750  trace.info() << "(" << nbok << "/" << nb
751  << ") ! plane.isExtendable( pts.begin(), pts.end() )"
752  << std::endl;
753  if ( ! check )
754  trace.warning() << plane << " last=" << pts.back() << std::endl
755  << "a=" << a << " b=" << b << " c=" << c << " d=" << d << std::endl;
756  ++nb; nbok += plane.extend( pts.begin(), pts.end() - 3 ); // should be ok
757  trace.info() << "(" << nbok << "/" << nb
758  << ") plane.extend( pts.begin(), pts.end() - 3)"
759  << std::endl;
760  ++nb; nbok += ! plane.extend( pts.end() - 3, pts.end() ); // should not be ok
761  trace.info() << "(" << nbok << "/" << nb
762  << ") ! plane.extend( pts.end() - 3, pts.end() )"
763  << std::endl;
764  }
765  trace.endBlock();
766  return nb == nbok;
767 }
768 
769 
770 
772 // Standard services - public :
773 
774 int main( int /*argc*/, char**/* argv */)
775 {
776  using namespace Z3i;
777 
778  trace.beginBlock ( "Testing class ChordNaivePlaneComputer" );
779  bool res = true
780  && testChordNaivePlaneComputer()
781  && checkManyPlanes<ChordNaivePlaneComputer<Z3i::Space, Z3i::Point, DGtal::int32_t> >( 4, 100, 200 )
782  && checkManyPlanes<ChordNaivePlaneComputer<Z3i::Space, Z3i::Point, DGtal::int32_t> >( 8, 100, 200 )
783  && checkManyPlanes<ChordNaivePlaneComputer<Z3i::Space, Z3i::Point, DGtal::int32_t> >( 20, 100, 200 )
784  && checkManyPlanes<ChordNaivePlaneComputer<Z3i::Space, Z3i::Point, DGtal::int32_t> >( 100, 100, 200 )
785  && checkManyPlanes<ChordNaivePlaneComputer<Z3i::Space, Z3i::Point, DGtal::int64_t> >( 2000, 100, 200 )
786  && checkWidths<DGtal::int64_t, ChordNaivePlaneComputer<Z3i::Space, Z3i::Point, DGtal::int64_t> >( 100, 1000000, 1000 )
787  && checkExtendWithManyPoints<ChordGenericNaivePlaneComputer<Z3i::Space, Z3i::Point, DGtal::int64_t> >( 100, 100, 200 );
788 
789  trace.emphase() << ( res ? "Passed." : "Error." ) << endl;
790  trace.endBlock();
791  return res ? 0 : 1;
792 }
793 // //
DGtal::COBANaivePlaneComputer::init
void init(Dimension axis, InternalInteger diameter, InternalInteger widthNumerator=NumberTraits< InternalInteger >::ONE, InternalInteger widthDenominator=NumberTraits< InternalInteger >::ONE)
DGtal::COBANaivePlaneComputer::end
ConstIterator end() const
DGtal::COBANaivePlaneComputer::begin
ConstIterator begin() const
DGtal::COBANaivePlaneComputer::isExtendable
bool isExtendable(const Point &p) const
DGtal::COBANaivePlaneComputer::extend
bool extend(const Point &p)
DGtal::ChordGenericNaivePlaneComputer
Aim: A class that recognizes pieces of digital planes of given axis width. When the width is 1,...
Definition: ChordGenericNaivePlaneComputer.h:119
DGtal::ChordNaivePlaneComputer
Aim: A class that contains the chord-based algorithm for recognizing pieces of digital planes of give...
Definition: ChordNaivePlaneComputer.h:158
DGtal::IntegerComputer::abs
static Integer abs(IntegerParamType a)
DGtal::IntegerComputer::ceilDiv
Integer ceilDiv(IntegerParamType na, IntegerParamType nb) const
DGtal::Trace::beginBlock
void beginBlock(const std::string &keyword="")
DGtal::Trace::emphase
std::ostream & emphase()
DGtal::Trace::info
std::ostream & info()
DGtal::Trace::warning
std::ostream & warning()
DGtal::Trace::endBlock
double endBlock()
NaivePlaneComputer
COBANaivePlaneComputer< Z3, InternalInteger > NaivePlaneComputer
Definition: greedy-plane-segmentation-ex2.cpp:88
ConstIterator
MyDigitalSurface::ConstIterator ConstIterator
Definition: greedy-plane-segmentation-ex2.cpp:93
DGtal::concepts
Aim: Gathers several functions useful for concept checks.
Definition: CPositiveIrreducibleFraction.h:54
DGtal
DGtal is the top-level namespace which contains all DGtal functions and types.
DGtal::int64_t
boost::int64_t int64_t
signed 94-bit integer.
Definition: BasicTypes.h:74
DGtal::Dimension
DGtal::uint32_t Dimension
Definition: Common.h:136
DGtal::trace
Trace trace
Definition: Common.h:153
DGtal::NumberTraits
Aim: The traits class for all models of Cinteger.
Definition: NumberTraits.h:564
DGtal::concepts::CAdditivePrimitiveComputer
Aim: Defines the concept describing an object that computes some primitive from input points given gr...
Definition: CAdditivePrimitiveComputer.h:100
DGtal::concepts::CPointPredicate
Aim: Defines a predicate on a point.
Definition: CPointPredicate.h:81
boost::ForwardContainer
Go to http://www.sgi.com/tech/stl/ForwardContainer.html.
Definition: Boost.dox:110
checkWidth
bool checkWidth(Integer a, Integer b, Integer c, Integer d, int diameter, unsigned int nbtries)
Definition: testChordNaivePlaneComputer.cpp:426
testChordNaivePlaneComputer
bool testChordNaivePlaneComputer()
Definition: testChordNaivePlaneComputer.cpp:526
main
int main(int, char **)
Definition: testChordNaivePlaneComputer.cpp:774
checkPlane
bool checkPlane(Integer a, Integer b, Integer c, Integer d, int diameter, unsigned int nbtries)
Definition: testChordNaivePlaneComputer.cpp:61
checkPlanes
bool checkPlanes(unsigned int nbplanes, int diameter, unsigned int nbtries)
Definition: testChordNaivePlaneComputer.cpp:388
checkGenericPlane
bool checkGenericPlane(Integer a, Integer b, Integer c, Integer d, int diameter, unsigned int nbtries)
Definition: testChordNaivePlaneComputer.cpp:283
checkExtendWithManyPoints
bool checkExtendWithManyPoints(unsigned int diameter, unsigned int nbplanes, unsigned int nbpoints)
Definition: testChordNaivePlaneComputer.cpp:696
getRandomInteger
Integer getRandomInteger(Integer first, Integer after_last)
Definition: testChordNaivePlaneComputer.cpp:50
checkWidths
bool checkWidths(unsigned int nbplanes, int diameter, unsigned int nbtries)
Definition: testChordNaivePlaneComputer.cpp:495
checkPlaneGroupExtension
bool checkPlaneGroupExtension(Integer a, Integer b, Integer c, Integer d, int diameter, unsigned int nbtries)
Definition: testChordNaivePlaneComputer.cpp:167
checkManyPlanes
bool checkManyPlanes(unsigned int diameter, unsigned int nbplanes, unsigned int nbpoints)
Definition: testChordNaivePlaneComputer.cpp:676
Point
MyPointD Point
Definition: testClone2.cpp:383
Generated on Thu Mar 9 2023 13:45:01 for DGtal by doxygen 1.9.1