DGtal 1.3.0
Loading...
Searching...
No Matches
testDistanceTransformation.cpp
Go to the documentation of this file.
1
31#include <iostream>
32#include <iomanip>
33#include "DGtal/base/Common.h"
34
35#include "DGtal/kernel/SpaceND.h"
36#include "DGtal/kernel/domains/HyperRectDomain.h"
37#include "DGtal/images/ImageSelector.h"
38#include "DGtal/geometry/volumes/distance/DistanceTransformation.h"
39#include "DGtal/geometry/volumes/distance/VoronoiMap.h"
40#include "DGtal/geometry/volumes/distance/ExactPredicateLpSeparableMetric.h"
41#include "DGtal/geometry/volumes/distance/InexactPredicateLpSeparableMetric.h"
42#include "DGtal/io/colormaps/HueShadeColorMap.h"
43#include "DGtal/io/colormaps/GrayscaleColorMap.h"
44#include "DGtal/shapes/Shapes.h"
45#include "DGtal/helpers/StdDefs.h"
46#include "DGtal/shapes/ShapeFactory.h"
47#include "DGtal/io/boards/Board2D.h"
48#include "DGtal/images/SimpleThresholdForegroundPredicate.h"
50
51using namespace std;
52using namespace DGtal;
53using namespace DGtal::functors;
54
56// Functions for testing class DistanceTransformation.
58
59template<typename Image>
60void randomSeeds(Image &input, const unsigned int nb, const int value)
61{
62 typename Image::Point p, low = input.domain().lowerBound(), up = input.domain().upperBound();
63 typename Image::Vector ext;
64
65 for (Dimension i = 0; i < Image::Domain::dimension; i++)
66 ext[i] = up[i] - low[i] + 1;
67
68
69 for (unsigned int k = 0 ; k < nb; k++)
70 {
71 for (unsigned int dim = 0; dim < Image::dimension; dim++)
72 {
73 p[dim] = rand() % (ext[dim]) + low[dim];
74 }
75 input.setValue(p, value);
76 }
77}
78
79
80
86{
87 bool allfine;
88
89 trace.beginBlock ( "Testing the whole DT computation" );
90
91 typedef SpaceND<2> TSpace;
92 typedef TSpace::Point Point;
95 Point a ( 2, 2 );
96 Point b ( 15, 15 );
98 Image image ( Domain(a, b ));
99
100 for ( unsigned k = 0; k < 49; k++ )
101 {
102 a[0] = ( k / 7 ) + 5;
103 a[1] = ( k % 7 ) + 5;
104 image.setValue ( a, 128 );
105 }
106 a= Point(2,2);
107
109 Predicate aPredicate(image,0);
110
112 L2Metric l2;
113 Domain dom(a,b);
115 VoronoiMap<Z2i::Space, Predicate, L2Metric> voro(&dom,&aPredicate,&l2);
116
117 Board2D board;
119 Display2DFactory::drawImage<Gray>(board, image, (unsigned int)0, (unsigned int)255);
120 board.saveSVG ( "image-preDT.svg" );
121 //We just iterate on the Domain points and print out the point coordinates.
122 std::copy ( image.begin(),
123 image.end(),
124 std::ostream_iterator<unsigned int> ( std::cout, " " ) );
125
126 trace.info()<<std::endl;
127 for(int i=2;i<=15;++i)
128 {
129 for(int j=2;j<=15;++j)
130 trace.info()<< image(Point(i,j))<<" ";
131 trace.info()<<std::endl;
132 }
133
134
135 trace.warning() << dt << endl;
136 trace.info() <<std::endl;
137
140
141 for (; it != itend; ++it)
142 {
143 std::cout << (*it) << " ";
144 }
145 std::cout << std::endl;
146
147 trace.info()<<std::endl;
148 for(int i=2;i<=15;++i)
149 {
150 for(int j=2;j<=15;++j)
151 trace.info()<< dt(Point(i,j))<<" ";
152 trace.info()<<std::endl;
153 }
154
155
156 trace.info()<<std::endl;
157 for(int i=2;i<=15;++i)
158 {
159 for(int j=2;j<=15;++j)
160 {
161 Point p= dt.getVoronoiSite(Point(i,j));
162 if (p==Point(i,j))
163 trace.info()<<"-,- ";
164 else
165 trace.info()<< p[0]<<","<<p[1]<<" ";
166 }
167 trace.info()<<std::endl;
168 }
169
170
171 allfine = true;
172
173 trace.info()<<std::endl;
174 for(int i=2;i<=15;++i)
175 {
176 for(int j=2;j<=15;++j)
177 {
178 Point p= voro(Point(i,j));
179 if (p != dt.getVoronoiSite(Point(i,j)))
180 allfine = false;
181 if (p==Point(i,j))
182 trace.info()<<"-,- ";
183 else
184 trace.info()<< p[0]<<","<<p[1]<<" ";
185 }
186 trace.info()<<std::endl;
187 }
188
189 board.clear();
190 Display2DFactory::drawImage<Gray>(board, dt, (DGtal::int64_t)0, (DGtal::int64_t)16);
191 board.saveSVG ( "image-postDT.svg" );
192 trace.info() << dt << endl;
193 trace.endBlock();
194
195 return allfine;
196}
202{
203 unsigned int nbok = 0;
204 unsigned int nb = 0;
205
206 trace.beginBlock ( "Testing the Neg DT computation" );
207
208 typedef Z2i::Space TSpace;
209 typedef TSpace::Point Point;
210 typedef Z2i::Domain Domain;
212 Point a ( -10, -10 );
213 Point b ( 10, 10 );
215 Image image ( Domain( a, b ));
216
217 for(int y=-10; y<=10;y++)
218 for(int x=-10; x<=10;x++)
219 {
220 if ((abs(x)<7) && (abs(y)<5))
221 image.setValue(Point(x,y),1);
222 else
223 image.setValue(Point(x,y),0);
224 }
225
227 Predicate aPredicate(image,0);
228
230 L2Metric l2;
231 Board2D board;
233 Display2DFactory::drawImage<Gray>(board, image, (unsigned int)0, (unsigned int)1);
234 board.saveSVG ( "image-preDT-neg.svg" );
235
236
237 for(int y=-10; y<=10;y++)
238 {
239 for(int x=-10; x<=10;x++)
240 {
241 std::cout<<image(Point(x,y))<<" ";
242 }
243 std::cout<<std::endl;
244 }
245
246 trace.info()<<"Domain "<<Domain(a,b)<<std::endl;
247 Domain dom(a,b);
249
252 itend = dt.constRange().end();
253 it != itend ; ++it)
254 if ((*it) > maxv)
255 maxv = (*it);
256
257 for(int y=-10; y<=10;y++)
258 {
259 for(int x=-10; x<=10;x++)
260 {
261 std::cout<<dt(Point(x,y))<<" ";
262 }
263 std::cout<<std::endl;
264 }
265
266
267
268 trace.warning() << dt << endl;
269 trace.warning() << dt.domain() << endl;
270
271
272 trace.info() << "Exporting..." << endl;
273 board.clear();
274 Display2DFactory::drawImage<Gray>(board, dt, 0, maxv);
275 board.saveSVG ( "image-postDT-neg.svg" );
276
277 trace.info() << "Done..." << endl;
278
279 trace.info() << dt << endl;
280
281 trace.endBlock();
282
283 return nbok == nb;
284}
285
286
288{
289 unsigned int nbok = 0;
290 unsigned int nb = 0;
291
292 trace.beginBlock ( "Testing the whole DT computation from a Set" );
293
294 typedef SpaceND<2> TSpace;
296
297
298 Board2D board;
299
300 AccFlower2D<Z2i::Space> flower(Z2i::Point(0,0), 30, 5,2,0);
301 Z2i::Domain domain(flower.getLowerBound(), flower.getUpperBound());
303
305
306 // Since 0.6, models of CDigitalSet are models of CPointPredicate.
307 // SetPredicate<Z2i::DigitalSet> aPredicate(aSet);
310 L2Metric l2;
311 L2DT dt(&domain,&aSet, &l2);
314 L1Metric l1;
315 L1DT dt1(&domain,&aSet, &l1);
316
317 L2DT::Value maxv = 0;
318 for ( L2DT::ConstRange::ConstIterator it = dt.constRange().begin(), itend = dt.constRange().end();
319 it != itend; ++it)
320 if ( (*it) > maxv)
321 maxv = (*it);
322 trace.error() << "MaxV="<<maxv<<std::endl;
323 Display2DFactory::drawImage<Hue>(board, dt, 0, maxv+1);
324 board.saveSVG ( "image-DTSet.svg" );
325
326 board.clear();
327 maxv = 0;
328 for ( L1DT::ConstRange::ConstIterator it = dt1.constRange().begin(), itend = dt1.constRange().end();
329 it != itend; ++it)
330 if ( (*it) > maxv)
331 maxv = (*it);
332 trace.error() << "MaxV="<<maxv<<std::endl;
333 Display2DFactory::drawImage<Hue>(board, dt1, 0, maxv+1);
334 board.saveSVG ( "image-DTSet-l1.svg" );
335 trace.endBlock();
336
337 return nbok == nb;
338}
339
345{
346 unsigned int nbok = 0;
347 unsigned int nb = 0;
348
349 trace.beginBlock ( "Testing DT computation with Infinity values at the first step" );
350
351 typedef SpaceND<2> TSpace;
352 typedef TSpace::Point Point;
355
356 Point a (0, 0 );
357 Point b ( 128, 128 );
358
360 Image image ( Domain(a, b ));
361
362 for ( Image::Iterator it = image.begin(), itend = image.end();it != itend; ++it)
363 *it = 128 ;
364
365
366 randomSeeds(image, 19, 0);
367
369 Predicate aPredicate(image,0);
370
372 L2Metric l2;
373 Domain dom(a,b);
375
376 Board2D board;
378 Display2DFactory::drawImage<Hue>(board, image, (unsigned int)0, (unsigned int)150);
379 board.saveSVG ( "image-preDT-border.svg" );
380
382 for ( DistanceTransformation<TSpace, Predicate, L2Metric>::ConstRange::ConstIterator it = dt.constRange().begin(), itend = dt.constRange().end();it != itend; ++it)
383 if ( (*it) > maxv)
384 maxv = (*it);
385
387 for (unsigned int y = 0; y < 33; y++)
388 {
389 for (unsigned int x = 0; x < 33; x++)
390 {
391 std::cout << std::setw(4) << (*it) << " ";
392 ++it;
393 }
394 std::cout << std::endl;
395 }
396
397 trace.warning() << dt << "MaxV = " << maxv << endl;
398
399
400 board.clear();
401 Display2DFactory::drawImage<Hue>(board, dt, (DGtal::int64_t)0, (DGtal::int64_t)maxv+1);
402 board.saveSVG ( "image-postDT-border.svg" );
403
404
405 trace.info() << dt << endl;
406
407 trace.endBlock();
408
409 return nbok == nb;
410}
411
412
418{
419 unsigned int nbok = 0;
420 unsigned int nb = 0;
421
422 trace.beginBlock ( "Testing 3D DT computation" );
423
424 typedef SpaceND<3> TSpace;
425 typedef TSpace::Point Point;
427 Point a ( 0, 0, 0 );
428 Point b ( 15, 15, 15 );
430 Image image ( Domain(a, b ));
431 Point c(8, 8, 8);
432 Domain dom(a, b);
433
434 for (Domain::ConstIterator it = dom.begin(),
435 itend = dom.end(); it != itend; ++it)
436 {
437 if ( ((*it) - c).norm() < 7)
438 image.setValue ( *it, 128 );
439 }
440
442 Predicate aPredicate(image,0);
443
445 L2Metric l2;
447
448 //We display the values on a 2D slice
449 for (unsigned int y = 0; y < 16; y++)
450 {
451 for (unsigned int x = 0; x < 16; x++)
452 {
453 Point p(x, y, 8);
454 std::cout << dt(p) << " ";
455 }
456 std::cout << std::endl;
457 }
458
459
460 trace.warning() << dt << endl;
461
462 trace.endBlock();
463
464 return nbok == nb;
465}
466
467
469{
470 unsigned int nbok = 0;
471 unsigned int nb = 0;
472
473 trace.beginBlock ( "Testing DT computation with Infinity values at the first step" );
474
475 typedef SpaceND<2> TSpace;
476 typedef TSpace::Point Point;
479
480 Point a (0, 0 );
481 Point b ( 128, 128 );
482 Domain dom(a,b);
483
485 Image image ( dom );
486
487 for ( Image::Iterator it = image.begin(), itend = image.end();it != itend; ++it)
488 (*it) = 128;
489
490
491 randomSeeds(image, 19, 0);
492
494 Predicate aPredicate(image,0);
495
496
497 // L_euc metric
499 L2Metric l2;
501 DT2 dt2(&dom, &aPredicate, &l2);
502 DT2 dt2_periodic( &dom, &aPredicate, &l2, {true, true} );
503
504 //L_infinity metric
505 //typedef DistanceTransformation<TSpace,Predicate, 0> DT;
506 //DT dt(Domain(a,b), aPredicate);;
507
508 //L_1 metric
510 L1Metric l1;
512 DT1 dt1(&dom,&aPredicate,&l1);
513 DT1 dt1_periodic( &dom, &aPredicate, &l1, {true, true} );
514
515 DGtal::int64_t maxv = 0;
516 for ( DistanceTransformation<TSpace,Predicate, L2Metric>::ConstRange::ConstIterator it = dt2.constRange().begin(), itend = dt2.constRange().end();it != itend; ++it)
517 if ( (*it) > maxv)
518 maxv = (*it);
519
520 trace.warning() << dt2 << "MaxV = " << maxv << endl;
521 //We display the values on a 2D slice
522 for (unsigned int y = 0; y < 16; y++)
523 {
524 for (unsigned int x = 0; x < 16; x++)
525 {
526 Point p(x, y);
527 std::cout << std::setw(4) << dt2(p) << " ";
528 }
529 std::cout << std::endl;
530 }
531
532 trace.info()<< "Exporting to SVG"<<endl;
533
534 Board2D board;
536 Display2DFactory::drawImage<Hue>(board, dt2, (DGtal::int64_t)0, (DGtal::int64_t)maxv+1);
537 board.saveSVG ( "image-DT-l2.svg" );
538
539 board.clear();
540 Display2DFactory::drawImage<Hue>(board, dt2_periodic, (DGtal::int64_t)0, (DGtal::int64_t)maxv+1);
541 board.saveSVG ( "image-DT-l2-periodic.svg" );
542
543 trace.info()<< "done"<<endl;
544
545
546
547 trace.info()<< "max L1"<<endl;
548 maxv = 0;
550 itend = dt1.constRange().end();
551 it2 != itend; ++it2)
552 {
553 if ( *it2 > maxv)
554 maxv = (*it2);
555 }
556
557 trace.info()<< "Exporting to SVG L1"<<endl;
558 board.clear();
559 Display2DFactory::drawImage<Hue>(board, dt1, (DGtal::int64_t)0, (DGtal::int64_t)maxv+1);
560 board.saveSVG ( "image-DT-l1.svg" );
561
562 board.clear();
563 Display2DFactory::drawImage<Hue>(board, dt1_periodic, (DGtal::int64_t)0, (DGtal::int64_t)maxv+1);
564 board.saveSVG ( "image-DT-l1-periodic.svg" );
565
566 trace.info()<< "done"<<endl;
567
568
569 trace.endBlock();
570
571 return nbok == nb;
572}
573
574template <typename Space, int norm>
575bool testCompareExactInexact(unsigned int size, unsigned int nb)
576{
577 trace.beginBlock("Checking Exact/Inexct predicate metrics");
581 typedef typename Space::Point Point;
582 typedef DigitalSetBySTLSet<Domain> Set;
583 // typedef NotPointPredicate<SetPredicate<Set> > NegPredicate;
584 typedef functors::NotPointPredicate<Set> NegPredicate;
585
586 Point low=Point::diagonal(0),
587 up=Point::diagonal(size);
588
589 Domain domain(low,up);
590 Set set(domain);
591
592 for(unsigned int i = 0; i<nb; ++i)
593 {
594 Point p;
595 for(unsigned int dim=0; dim<Space::dimension;++dim)
596 p[dim] = rand() % size;
597 set.insert(p);
598 }
599
600 trace.info()<< "Testing metrics "<<MetricEx()<<" "<<MetricInex(norm)<<std::endl;
601 trace.info()<< "Testing space dimension "<<Space::dimension<<std::endl;
602 trace.info()<< "Inserting "<<set.size() << " points."<<std::endl;
603
604 // SetPredicate<Set> setPred(set);
605 NegPredicate negPred(set);
606
609 MetricEx metricEx;
610 MetricInex metricInex(norm);
611 DTEx dtex(&domain, &negPred, &metricEx);
612 DTIn dtinex(&domain, &negPred, &metricInex);
613
614 double MSE=0.0;
615 typename DTEx::ConstRange::ConstIterator it=dtex.constRange().begin(), itend=dtex.constRange().end();
616 typename DTIn::ConstRange::ConstIterator it2 = dtinex.constRange().begin();
617 for( ; it != itend; ++it, ++it2)
618 MSE += ((*it) - (*it2))*((*it) - (*it2));
619
620 trace.warning()<<"Resulting MSE= "<<MSE;
621 trace.endBlock();
622 return true;
623}
624
626// Standard services - public :
627
628int main ( int argc, char** argv ){
629
630 trace.beginBlock ( "Testing class DistanceTransformation" );
631 trace.info() << "Args:";
632 for ( int i = 0; i < argc; ++i )
633 trace.info() << " " << argv[ i ];
634 trace.info() << endl;
635
637 && testDTFromSet()
640 && testChessboard()
641 && testDTFromSet()
642 && testCompareExactInexact<Z2i::Space, 2>(50, 50)
643 && testCompareExactInexact<Z3i::Space, 2>(50, 50)
644 && testCompareExactInexact<Z2i::Space, 4>(50, 50)
645 && testCompareExactInexact<Z3i::Space, 4>(50, 50)
646 ;
647 //&& ... other tests
648 trace.emphase() << ( res ? "Passed." : "Error." ) << endl;
649 trace.endBlock();
650 return res ? 0 : 1;
651}
652// //
Aim: Model of the concept StarShaped represents any accelerated flower in the plane.
Definition: AccFlower2D.h:65
RealPoint getLowerBound() const
Definition: AccFlower2D.h:134
RealPoint getUpperBound() const
Definition: AccFlower2D.h:143
Aim: This class specializes a 'Board' class so as to display DGtal objects more naturally (with <<)....
Definition: Board2D.h:71
This class adapts any iterator so that operator* returns another element than the one pointed to by t...
Aim: A wrapper class around a STL associative container for storing sets of digital points within som...
Aim: A container class for storing sets of digital points within some given domain.
Aim: Implementation of the linear in time distance transformation for separable metrics.
SeparableMetric::Value Value
Definition of the image value type.
Aim: implements separable l_p metrics with exact predicates.
Aim: This class template may be used to (linearly) convert scalar values in a given range into gray l...
Aim: This class template may be used to (linearly) convert scalar values in a given range into a colo...
Iterator for HyperRectDomain.
Aim: Parallelepidec region of a digital space, model of a 'CDomain'.
const ConstIterator & begin() const
const ConstIterator & end() const
Aim: implements association bewteen points lying in a digital domain and values.
Definition: Image.h:70
Aim: implements separable l_p metrics with approximated predicates.
static void euclideanShaper(TDigitalSet &aSet, const TShapeFunctor &aFunctor, const double h=1.0)
static const Dimension dimension
static constants to store the dimension.
Definition: SpaceND.h:132
void beginBlock(const std::string &keyword="")
std::ostream & warning()
std::ostream & emphase()
std::ostream & error()
std::ostream & info()
double endBlock()
Aim: Implementation of the linear in time Voronoi map construction.
Definition: VoronoiMap.h:127
Aim: Define a simple Foreground predicate thresholding image values given a single thresold....
void clear(const DGtal::Color &color=DGtal::Color::None)
Definition: Board.cpp:152
void saveSVG(const char *filename, PageSize size=Board::BoundingBox, double margin=10.0) const
Definition: Board.cpp:1012
void setUnit(Unit unit)
Definition: Board.cpp:240
functors namespace gathers all DGtal functors.
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 predicate returns true when the point predicate given at construction return false....
int main()
Definition: testBits.cpp:56
MyPointD Point
Definition: testClone2.cpp:383
bool testChessboard()
bool testDistanceTransformationBorder()
bool testCompareExactInexact(unsigned int size, unsigned int nb)
bool testDistanceTransformation()
void randomSeeds(Image &input, const unsigned int nb, const int value)
bool testDistanceTransformationNeg()
bool testDTFromSet()
bool testDistanceTransformation3D()
Domain domain
HyperRectDomain< Space > Domain