DGtal 1.3.0
Searching...
No Matches
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"
37#include "DGtal/geometry/surfaces/ChordNaivePlaneComputer.h"
38#include "DGtal/geometry/surfaces/ChordGenericNaivePlaneComputer.h"
40
41using namespace std;
42using namespace DGtal;
43using namespace DGtal::concepts;
44
46// Functions for testing class ChordNaivePlaneComputer.
48
49template <typename Integer>
51{
52 Integer r = (Integer) rand();
53 return ( r % (after_last - first) ) + first;
54}
55
59template <typename Integer, typename NaivePlaneComputer>
60bool
62 int diameter, unsigned int nbtries )
63{
64 typedef typename NaivePlaneComputer::Point Point;
65 typedef typename Point::Component PointInteger;
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;
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 ] = NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - b * y - c * z, a ) ); break;
94 case 1: p[ 1 ] = NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - a * x - c * z, b ) ); break;
95 case 2: p[ 2 ] = 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 ] = NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - b * y - c * z, a ) ); break;
136 case 1: p[ 1 ] = NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - a * x - c * z, b ) ); break;
137 case 2: p[ 2 ] = 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
165template <typename Integer, typename NaivePlaneComputer>
166bool
168 int diameter, unsigned int nbtries )
169{
170 typedef typename NaivePlaneComputer::Point Point;
171 typedef typename Point::Component PointInteger;
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 ) {
204 ( ic.ceilDiv( d - b * y - c * z, a ) ); break;
206 ( ic.ceilDiv( d - a * x - c * z, b ) ); break;
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 ] = NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - b * y - c * z, a ) ); break;
250 case 1: p[ 1 ] = NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - a * x - c * z, b ) ); break;
251 case 2: p[ 2 ] = 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
281template <typename Integer, typename GenericNaivePlaneComputer>
282bool
284 int diameter, unsigned int nbtries )
285{
286 typedef typename GenericNaivePlaneComputer::Point Point;
287 typedef typename Point::Component PointInteger;
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 ] = NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - b * y - c * z, a ) ); break;
316 case 1: p[ 1 ] = NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - a * x - c * z, b ) ); break;
317 case 2: p[ 2 ] = 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 ] = NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - b * y - c * z, a ) ); break;
358 case 1: p[ 1 ] = NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - a * x - c * z, b ) ); break;
359 case 2: p[ 2 ] = 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
386template <typename Integer, typename NaivePlaneComputer>
387bool
388checkPlanes( 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
424template <typename Integer, typename NaivePlaneComputer>
425bool
427 int diameter, unsigned int nbtries )
428{
429 typedef typename NaivePlaneComputer::Point Point;
430 typedef typename NaivePlaneComputer::InternalScalar InternalScalar;
432 Integer absA = ic.abs( a );
433 Integer absB = ic.abs( b );
434 Integer absC = ic.abs( c );
435 Integer x, y, z;
436 Dimension axis;
437 if ( ( absA >= absB ) && ( absA >= absC ) )
438 axis = 0;
439 else if ( ( absB >= absA ) && ( absB >= absC ) )
440 axis = 1;
441 else
442 axis = 2;
443 // Checks that points within the naive plane are correctly recognized.
444 unsigned int nb = 0;
445 unsigned int nbok = 0;
446 std::vector<Point> points( nbtries );
447 for ( unsigned int i = 0; i != nbtries; ++i )
448 {
449 Point & p = points[ i ];
450 p[ 0 ] = getRandomInteger<Integer>( -diameter+1, diameter );
451 p[ 1 ] = getRandomInteger<Integer>( -diameter+1, diameter );
452 p[ 2 ] = getRandomInteger<Integer>( -diameter+1, diameter );
453 x = (Integer) p[ 0 ];
454 y = (Integer) p[ 1 ];
455 z = (Integer) p[ 2 ];
456 switch ( axis ) {
457 case 0: p[ 0 ] = NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - b * y - c * z, a ) ); break;
458 case 1: p[ 1 ] = NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - a * x - c * z, b ) ); break;
459 case 2: p[ 2 ] = NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - a * x - b * y, c ) ); break;
460 }
461 }
462 trace.beginBlock( "Computing axis width." );
463 trace.info() << "- plane is "
464 << d << " <= " << a << "*x"
465 << "+" << b << "*y"
466 << "+" << c << "*z"
467 << " <= d + max(|a|,|b|,|c|)"
468 << std::endl;
469 trace.info() << "- " << points.size() << " points tested in diameter " << diameter
470 << std::endl;
471 double min = -1.0;
472 for ( unsigned int i = 0; i < 3; ++i )
473 {
474 std::pair<InternalScalar, InternalScalar> width
475 = NaivePlaneComputer::computeAxisWidth( i, points.begin(), points.end() );
476 double wn = NumberTraits<InternalScalar>::castToDouble( width.first );
477 double wd = NumberTraits<InternalScalar>::castToDouble( width.second );
478 trace.info() << " (" << i << ") width=" << (wn/wd) << std::endl;
479 if ( min < 0.0 ) min = wn/wd;
480 else if ( wn/wd < min ) min = wn/wd;
481 }
482 ++nb; nbok += (min < 1.0 ) ? 1 : 0;
483 trace.info() << "(" << nbok << "/" << nb << ") min width = " << min
484 << " < 1.0" << std::endl;
485 ++nb; nbok += (0.9 < min ) ? 1 : 0;
486 trace.info() << "(" << nbok << "/" << nb << ") min width = " << min
487 << " > 0.9" << std::endl;
488 trace.endBlock();
489 return nb == nbok;
490}
491
492template <typename Integer, typename NaivePlaneComputer>
493bool
494checkWidths( unsigned int nbplanes, int diameter, unsigned int nbtries )
495{
496 //using namespace Z3i;
497 //typedef ChordNaivePlaneComputer<Z3, Integer> NaivePlaneComputer;
498 unsigned int nb = 0;
499 unsigned int nbok = 0;
500 for ( unsigned int nbp = 0; nbp < nbplanes; ++nbp )
501 {
502 Integer a = getRandomInteger<Integer>( (Integer) 0, (Integer) diameter / 2 );
503 Integer b = getRandomInteger<Integer>( (Integer) 0, (Integer) diameter / 2 );
504 Integer c = getRandomInteger<Integer>( (Integer) 0, (Integer) diameter / 2 );
505 Integer d = getRandomInteger<Integer>( (Integer) 0, (Integer) diameter / 2 );
506 if ( ( a != 0 ) || ( b != 0 ) || ( c != 0 ) )
507 {
508 ++nb; nbok += checkWidth<Integer, NaivePlaneComputer>( a, b, c, d, diameter, nbtries ) ? 1 : 0;
509 if ( nb != nbok )
510 {
511 std::cerr << "[ERROR] (checkWidth) for plane " << a << " * x + "
512 << b << " * y + " << c << " * z = " << d << std::endl;
513 break;
514 }
515 }
516 }
517 return nb == nbok;
518}
519
520
526{
527 unsigned int nbok = 0;
528 unsigned int nb = 0;
529 typedef DGtal::int64_t Integer;
530 typedef DGtal::Z3i::Z3 Space;
531 typedef DGtal::Z3i::Point Point;
533 typedef ChordGenericNaivePlaneComputer<Space, Point, Integer> GenericNaivePlaneComputer;
534
535 BOOST_CONCEPT_ASSERT(( CAdditivePrimitiveComputer< NaivePlaneComputer > ));
537 BOOST_CONCEPT_ASSERT(( boost::ForwardContainer< NaivePlaneComputer > ));
539 BOOST_CONCEPT_ASSERT(( CPointPredicate< NaivePlaneComputer::Primitive > ));
541
542 trace.beginBlock ( "Testing block: ChordNaivePlaneComputer instantiation." );
543 NaivePlaneComputer plane;
544 Point pt0( 0, 0, 0 );
545 plane.init( 2, 1, 1 );
546 bool pt0_inside = plane.extend( pt0 );
547 ++nb; nbok += pt0_inside == true ? 1 : 0;
548 trace.info() << "(" << nbok << "/" << nb << ") Plane=" << plane
549 << std::endl;
550 Point pt1( Point( 2, 0, 0 ) );
551 bool pt1_inside = plane.extend( pt1 );
552 ++nb; nbok += pt1_inside == true ? 1 : 0;
553 trace.info() << "(" << nbok << "/" << nb << ") add " << pt1
554 << " Plane=" << plane << std::endl;
555 Point pt2( Point( 0, 2, 2 ) );
556 bool pt2_inside = plane.extend( pt2 );
557 ++nb; nbok += pt2_inside == true ? 1 : 0;
558 trace.info() << "(" << nbok << "/" << nb << ") add " << pt2
559 << " Plane=" << plane << std::endl;
560
561 Point pt3( Point( 1, 1, 1 ) );
562 bool pt3_inside = plane.extend( pt3 );
563 ++nb; nbok += pt3_inside == true ? 1 : 0;
564 trace.info() << "(" << nbok << "/" << nb << ") add " << pt3
565 << " Plane=" << plane << std::endl;
566
567 Point pt4( Point( -10, -10, 10 ) );
568 bool pt4_inside = plane.extend( pt4 );
569 ++nb; nbok += pt4_inside == false ? 1 : 0;
570 trace.info() << "(" << nbok << "/" << nb << ") impossible add " << pt4
571 << " Plane=" << plane << std::endl;
572
573 Point pt5 = pt2 + Point( 1, 0, 1 );
574 bool pt5_inside = plane.extend( pt5 );
575 ++nb; nbok += pt5_inside == true ? 1 : 0;
576 trace.info() << "(" << nbok << "/" << nb << ") add " << pt5
577 << " Plane=" << plane << std::endl;
578
579 Point pt6 = pt5 + Point( 6, 0, 2 );
580 bool pt6_inside = plane.extend( pt6 );
581 ++nb; nbok += pt6_inside == true ? 1 : 0;
582 trace.info() << "(" << nbok << "/" << nb << ") add " << pt6
583 << " Plane=" << plane << std::endl;
584
585 NaivePlaneComputer plane2;
586 plane2.init( 2, 1, 1 );
587 plane2.extend( Point( 10, 0, 0 ) );
588 plane2.extend( Point( 0, 8, 0 ) );
589 plane2.extend( Point( 0, 0, 6 ) );
590 trace.info() << "(" << nbok << "/" << nb << ") "
591 << " Plane2=" << plane2 << std::endl;
592
593 ++nb; nbok += checkPlane<Integer,NaivePlaneComputer>( 11, 5, 19, 20, 100, 100 ) ? 1 : 0;
594 trace.info() << "(" << nbok << "/" << nb
595 << ") checkPlane<Integer,NaivePlaneComputer>( 11, 5, 19, 20, 100, 100 )"
596 << std::endl;
597
598 ++nb; nbok += checkGenericPlane<Integer,GenericNaivePlaneComputer>( 11, 5, 19, 20, 100, 100 ) ? 1 : 0;
599 trace.info() << "(" << nbok << "/" << nb
600 << ") checkGenericPlane<Integer,GenericNaivePlaneComputer>( 11, 5, 19, 20, 100, 100 )"
601 << std::endl;
602 ++nb; nbok += checkGenericPlane<Integer,GenericNaivePlaneComputer>( 17, 33, 7, 10, 100, 100 ) ? 1 : 0;
603 trace.info() << "(" << nbok << "/" << nb
604 << ") checkGenericPlane<Integer,GenericNaivePlaneComputer>( 17, 33, 7, 10, 100, 100 )"
605 << std::endl;
606 ++nb; nbok += checkPlane<Integer,NaivePlaneComputer>( 15, 8, 13, 15, 100, 100 ) ? 1 : 0;
607 trace.info() << "(" << nbok << "/" << nb
608 << ") checkPlane<Integer,NaivePlaneComputer>( 15, 8, 13, 15, 100, 100 )"
609 << std::endl;
610 ++nb; nbok += checkGenericPlane<Integer,GenericNaivePlaneComputer>( 15, 8, 13, 15, 100, 100 ) ? 1 : 0;
611 trace.info() << "(" << nbok << "/" << nb
612 << ") checkGenericPlane<Integer,GenericNaivePlaneComputer>( 15, 8, 13, 15, 100, 100 )"
613 << std::endl;
614 trace.endBlock();
615
616 {
617 trace.beginBlock ( "Testing block: ChordNaivePlaneComputer vertical instantiation." );
618 NaivePlaneComputer ppplane;
619 Point pppt0( 0, 0, 0 );
620 ppplane.init( 2, 5, 2 );
621 bool pppt0_inside = ppplane.extend( pppt0 );
622 ++nb; nbok += pppt0_inside == true ? 1 : 0;
623 trace.info() << "(" << nbok << "/" << nb << ") Plane=" << ppplane
624 << std::endl;
625 Point pppt1( 3, 2, 2 );
626 bool pppt1_inside = ppplane.extend( pppt1 );
627 ++nb; nbok += pppt1_inside == true ? 1 : 0;
628 trace.info() << "(" << nbok << "/" << nb << ") Plane=" << ppplane
629 << std::endl;
630 Point pppt2( 0, 0, 1 );
631 bool pppt2_inside = ppplane.extend( pppt2 );
632 ++nb; nbok += pppt2_inside == true ? 1 : 0;
633 trace.info() << "(" << nbok << "/" << nb << ") Plane=" << ppplane
634 << std::endl;
635 Point pppt3 = pppt1 + Point( 0, 0, 1 );
636 bool pppt3_inside = ppplane.extend( pppt3 );
637 ++nb; nbok += pppt3_inside == true ? 1 : 0;
638 trace.info() << "(" << nbok << "/" << nb << ") Plane=" << ppplane
639 << std::endl;
640 Point pppt4 = pppt3 + Point( 0, 0, 1 );
641 bool pppt4_inside = ppplane.extend( pt4 );
642 ++nb; nbok += pppt4_inside == true ? 1 : 0;
643 trace.info() << "(" << nbok << "/" << nb << ") Plane=" << ppplane
644 << std::endl;
645 trace.endBlock();
646 }
647
648 {
649 trace.beginBlock ( "Testing block: ChordNaivePlaneComputer vertical instantiation 2." );
650 NaivePlaneComputer pplane;
651 pplane.init( 1, 1, 1 );
652 Point ppt0( -6, -3, 5 );
653 bool ppt0_inside = pplane.extend( ppt0 );
654 ++nb; nbok += ppt0_inside == true ? 1 : 0;
655 trace.info() << "(" << nbok << "/" << nb << ") Plane=" << pplane
656 << std::endl;
657 Point ppt1( 4, 4, -5 );
658 bool ppt1_inside = pplane.extend( ppt1 );
659 ++nb; nbok += ppt1_inside == true ? 1 : 0;
660 trace.info() << "(" << nbok << "/" << nb << ") Plane=" << pplane
661 << std::endl;
662 Point ppt2( -5, -2, 4 );
663 bool ppt2_inside = pplane.extend( ppt2 );
664 ++nb; nbok += ppt2_inside == true ? 1 : 0;
665 trace.info() << "(" << nbok << "/" << nb << ") Plane=" << pplane
666 << std::endl;
667 trace.endBlock();
668 }
669
670 return nbok == nb;
671}
672
673template <typename NaivePlaneComputer>
674bool
675checkManyPlanes( unsigned int diameter,
676 unsigned int nbplanes,
677 unsigned int nbpoints )
678{
679 unsigned int nbok = 0;
680 unsigned int nb = 0;
681 typedef typename NaivePlaneComputer::InternalScalar Scalar;
682 stringstream ss (stringstream::out);
683 ss << "Testing block: Diameter is " << diameter << ". Check " << nbplanes << " planes with " << nbpoints << " points each.";
684 trace.beginBlock ( ss.str() );
685 ++nb; nbok += checkPlanes<Scalar,NaivePlaneComputer>( nbplanes, diameter, nbpoints ) ? 1 : 0;
686 trace.info() << "(" << nbok << "/" << nb
687 << ") checkPlanes<Scalar,NaivePlaneComputer>()"
688 << std::endl;
689 trace.endBlock();
690 return nbok == nb;
691}
692
693template <typename GenericNaivePlaneComputer>
694bool
695checkExtendWithManyPoints( unsigned int diameter,
696 unsigned int nbplanes,
697 unsigned int nbpoints )
698{
699 unsigned int nbok = 0;
700 unsigned int nb = 0;
701 typedef typename GenericNaivePlaneComputer::InternalScalar Integer;
702 typedef typename GenericNaivePlaneComputer::Point Point;
703 typedef typename Point::Coordinate PointInteger;
705
706 trace.beginBlock( "checkExtendWithManyPoints" );
707 for ( unsigned int j = 0; j < nbplanes; ++j )
708 {
709 Integer a = getRandomInteger<Integer>( (Integer) 0, (Integer) diameter / 2 );
710 Integer b = getRandomInteger<Integer>( (Integer) 0, (Integer) diameter / 2 );
711 Integer c = getRandomInteger<Integer>( (Integer) 1, (Integer) diameter / 2 );
712 Integer d = getRandomInteger<Integer>( (Integer) 0, (Integer) diameter / 2 );
713 GenericNaivePlaneComputer plane;
714 Dimension axis;
715 if ( ( a >= b ) && ( a >= c ) ) axis = 0;
716 else if ( ( b >= a ) && ( b >= c ) ) axis = 1;
717 else axis = 2;
718 plane.init( 1, 1 );
719
720 std::vector<Point> pts;
721 for ( unsigned int i = 0; i < nbpoints; ++i )
722 {
723 Point p;
724 p[ 0 ] = getRandomInteger<PointInteger>( -diameter+1, diameter );
725 p[ 1 ] = getRandomInteger<PointInteger>( -diameter+1, diameter );
726 p[ 2 ] = getRandomInteger<PointInteger>( -diameter+1, diameter );
727 Integer x = (Integer) p[ 0 ];
728 Integer y = (Integer) p[ 1 ];
729 Integer z = (Integer) p[ 2 ];
730 switch( axis ) {
731 case 0: p[ 0 ] = NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - b * y - c * z, a ) ); break;
732 case 1: p[ 1 ] = NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - a * x - c * z, b ) ); break;
733 case 2: p[ 2 ] = NumberTraits<Integer>::castToInt64_t( ic.ceilDiv( d - a * x - b * y, c ) ); break;
734 }
735 pts.push_back( p );
736 }
737 ++nb; nbok += plane.isExtendable( pts.begin(), pts.end() ); // should be ok
738 trace.info() << "(" << nbok << "/" << nb
739 << ") plane.isExtendable( pts.begin(), pts.end() )"
740 << std::endl;
741 Point & any0 = pts[ getRandomInteger<int>( 0, pts.size() ) ];
742 pts.push_back( any0 + Point(1,0,0) );
743 Point & any1 = pts[ getRandomInteger<int>( 0, pts.size() ) ];
744 pts.push_back( any1 + Point(0,1,0) );
745 Point & any2 = pts[ getRandomInteger<int>( 0, pts.size() ) ];
746 pts.push_back( any2 + Point(0,0,1) );
747 bool check = ! plane.isExtendable( pts.begin(), pts.end() ); // should not be ok
748 ++nb; nbok += check ? 1 : 0;
749 trace.info() << "(" << nbok << "/" << nb
750 << ") ! plane.isExtendable( pts.begin(), pts.end() )"
751 << std::endl;
752 if ( ! check )
753 trace.warning() << plane << " last=" << pts.back() << std::endl
754 << "a=" << a << " b=" << b << " c=" << c << " d=" << d << std::endl;
755 ++nb; nbok += plane.extend( pts.begin(), pts.end() - 3 ); // should be ok
756 trace.info() << "(" << nbok << "/" << nb
757 << ") plane.extend( pts.begin(), pts.end() - 3)"
758 << std::endl;
759 ++nb; nbok += ! plane.extend( pts.end() - 3, pts.end() ); // should not be ok
760 trace.info() << "(" << nbok << "/" << nb
761 << ") ! plane.extend( pts.end() - 3, pts.end() )"
762 << std::endl;
763 }
764 trace.endBlock();
765 return nb == nbok;
766}
767
768
769
771// Standard services - public :
772
773int main( int /*argc*/, char**/* argv */)
774{
775 using namespace Z3i;
776
777 trace.beginBlock ( "Testing class ChordNaivePlaneComputer" );
778 bool res = true
780 && checkManyPlanes<ChordNaivePlaneComputer<Z3i::Space, Z3i::Point, DGtal::int32_t> >( 4, 100, 200 )
782 && checkManyPlanes<ChordNaivePlaneComputer<Z3i::Space, Z3i::Point, DGtal::int32_t> >( 20, 100, 200 )
784 && checkManyPlanes<ChordNaivePlaneComputer<Z3i::Space, Z3i::Point, DGtal::int64_t> >( 2000, 100, 200 )
786 && checkExtendWithManyPoints<ChordGenericNaivePlaneComputer<Z3i::Space, Z3i::Point, DGtal::int64_t> >( 100, 100, 200 );
787
788 trace.emphase() << ( res ? "Passed." : "Error." ) << endl;
789 trace.endBlock();
790 return res ? 0 : 1;
791}
792// //
void init(Dimension axis, InternalInteger diameter, InternalInteger widthNumerator=NumberTraits< InternalInteger >::ONE, InternalInteger widthDenominator=NumberTraits< InternalInteger >::ONE)
ConstIterator end() const
ConstIterator begin() const
bool isExtendable(const Point &p) const
bool extend(const Point &p)
Aim: A class that recognizes pieces of digital planes of given axis width. When the width is 1,...
Aim: A class that contains the chord-based algorithm for recognizing pieces of digital planes of give...
Aim: This class gathers several types and methods to make computation with integers.
static Integer abs(IntegerParamType a)
Integer ceilDiv(IntegerParamType na, IntegerParamType nb) const
void beginBlock(const std::string &keyword="")
std::ostream & warning()
std::ostream & emphase()
std::ostream & info()
double endBlock()
Point::Coordinate Integer
COBANaivePlaneComputer< Z3, InternalInteger > NaivePlaneComputer
Aim: Gathers several functions useful for concept checks.
DGtal is the top-level namespace which contains all DGtal functions and types.
boost::int64_t int64_t
signed 94-bit integer.
Definition: BasicTypes.h:74
DGtal::uint32_t Dimension
Definition: Common.h:137
Trace trace
Definition: Common.h:154
STL namespace.
Aim: The traits class for all models of Cinteger.
Definition: NumberTraits.h:564
Aim: Defines the concept describing an object that computes some primitive from input points given gr...
Aim: Defines a predicate on a point.
Go to http://www.sgi.com/tech/stl/ForwardContainer.html.
Definition: Boost.dox:110
int main()
Definition: testBits.cpp:56
bool checkWidth(Integer a, Integer b, Integer c, Integer d, int diameter, unsigned int nbtries)
bool testChordNaivePlaneComputer()
bool checkPlane(Integer a, Integer b, Integer c, Integer d, int diameter, unsigned int nbtries)
bool checkPlanes(unsigned int nbplanes, int diameter, unsigned int nbtries)
bool checkGenericPlane(Integer a, Integer b, Integer c, Integer d, int diameter, unsigned int nbtries)
bool checkExtendWithManyPoints(unsigned int diameter, unsigned int nbplanes, unsigned int nbpoints)
Integer getRandomInteger(Integer first, Integer after_last)
bool checkWidths(unsigned int nbplanes, int diameter, unsigned int nbtries)
bool checkPlaneGroupExtension(Integer a, Integer b, Integer c, Integer d, int diameter, unsigned int nbtries)
bool checkManyPlanes(unsigned int diameter, unsigned int nbplanes, unsigned int nbpoints)
MyPointD Point
Definition: testClone2.cpp:383