DGtal 1.3.0
Loading...
Searching...
No Matches
testConvexHull2D.cpp
Go to the documentation of this file.
1
31#include <iostream>
32#include "DGtal/base/Common.h"
33#include "DGtal/base/Circulator.h"
34
35#include "DGtal/kernel/PointVector.h"
36#include "DGtal/geometry/tools/Hull2DHelpers.h"
37#include "DGtal/geometry/tools/PolarPointComparatorBy2x2DetComputer.h"
38#include "DGtal/geometry/tools/determinant/InHalfPlaneBySimple3x3Matrix.h"
39#include "DGtal/geometry/tools/determinant/PredicateFromOrientationFunctor2.h"
40#include "DGtal/io/boards/Board2D.h"
42
43using namespace std;
44using namespace DGtal;
45
46
47
48
50// Functions for testing the functions devoted to convex hull computations.
52
65template <typename ForwardIterator>
66bool circularlyEqual(const ForwardIterator& first1, const ForwardIterator& last1,
67 const ForwardIterator& first2, const ForwardIterator& last2 )
68{
69 ASSERT( first2 != last2 );
70
71 //find a element of the first range equal to the first element of the second range
72 ForwardIterator start1 = find( first1, last1, *first2 );
73 if (start1 == last1)
74 return false;
75 else
76 {
77 bool areEqual = true; //true if the two ranges are equal up to circular shifts
78 //check whether the two ranges are equal or not
79 Circulator<ForwardIterator> c1(start1, first1, last1);
80 Circulator<ForwardIterator> cEnd1 = c1;
81 Circulator<ForwardIterator> c2(first2, first2, last2);
82 do
83 {
84 if (*c1 != *c2)
85 areEqual = false;
86 else
87 {
88 ++c1;
89 ++c2;
90 }
91 } while ( (c1 != cEnd1) && (areEqual) );
92 return areEqual;
93
94 }
95}
96
97
102bool testConvexHull2D()
103{
104 unsigned int nbok = 0;
105 unsigned int nb = 0;
106
107 typedef PointVector<2,DGtal::int32_t> Point;
108
109 trace.beginBlock ( "One simple test..." );
110
111 vector<Point> data, g, res;
112 //data
113 data.push_back( Point(2,0) );
114 data.push_back( Point(4,0) );
115 data.push_back( Point(0,3) );
116 data.push_back( Point(0,-4) );
117 data.push_back( Point(3,4) );
118 data.push_back( Point(5,0) );
119 data.push_back( Point(4,3) );
120 data.push_back( Point(0,5) );
121 data.push_back( Point(-3,-4) );
122 data.push_back( Point(-5,0) );
123 data.push_back( Point(-4,-3) );
124 data.push_back( Point(0,-5) );
125 data.push_back( Point(3,-4) );
126 data.push_back( Point(4,-3) );
127 data.push_back( Point(-3,4) );
128 data.push_back( Point(-4,3) );
129 //ground truth
130 g.push_back( Point(5,0) );
131 g.push_back( Point(4,3) );
132 g.push_back( Point(3,4) );
133 g.push_back( Point(0,5) );
134 g.push_back( Point(-3,4) );
135 g.push_back( Point(-4,3) );
136 g.push_back( Point(-5,0) );
137 g.push_back( Point(-4,-3) );
138 g.push_back( Point(-3,-4) );
139 g.push_back( Point(0,-5) );
140 g.push_back( Point(3,-4) );
141 g.push_back( Point(4,-3) );
142
143 //geometric predicate
144 typedef InHalfPlaneBySimple3x3Matrix<Point, DGtal::int64_t> Functor;
145 Functor functor;
146 typedef PredicateFromOrientationFunctor2<Functor> Predicate;
147 Predicate predicate( functor );
148
149 //namespace
150 using namespace functions::Hull2D;
151
152 //andrew algorithm
153 trace.info() << " andrew algorithm " << std::endl;
154 andrewConvexHullAlgorithm( data.begin(), data.end(), back_inserter( res ), predicate );
155
156 copy(res.begin(), res.end(), ostream_iterator<Point>( cout, " " ) );
157 cout << endl;
158
159 if ( (res.size() == g.size()) &&
160 (circularlyEqual(res.begin(), res.end(), g.begin(), g.end())) )
161 nbok++;
162 nb++;
163 trace.info() << "(" << nbok << "/" << nb << ") " << endl;
164
165 //graham algorithm
166 res.clear();
167 trace.info() << " graham algorithm " << std::endl;
168 functors::PolarPointComparatorBy2x2DetComputer<Point> comparator;
169 grahamConvexHullAlgorithm( data.begin(), data.end(), back_inserter( res ), predicate, comparator );
170
171 copy(res.begin(), res.end(), ostream_iterator<Point>( cout, " " ) );
172 cout << endl;
173
174 if ( (res.size() == g.size()) &&
175 (circularlyEqual(res.begin(), res.end(), g.begin(), g.end())) )
176 nbok++;
177 nb++;
178 trace.info() << "(" << nbok << "/" << nb << ") " << endl;
179
180 //melkman algorithm
181 res.clear();
182 trace.info() << " melkman algorithm " << std::endl;
183 sort( data.begin(), data.end() );
184 melkmanConvexHullAlgorithm( data.begin(), data.end(), back_inserter( res ), functor );
185
186 copy(res.begin(), res.end(), ostream_iterator<Point>( cout, " " ) );
187 cout << endl;
188
189 if ( (res.size() == g.size()) &&
190 (circularlyEqual(res.begin(), res.end(), g.begin(), g.end())) )
191 nbok++;
192 nb++;
193 trace.info() << "(" << nbok << "/" << nb << ") " << endl;
194 // melkman on line construction of convex hull:
195 trace.info() << "on line convex hull construction" << std::endl;
196 DGtal::MelkmanConvexHull<Point, Functor> ch;
197 for(vector<Point>::const_iterator it = data.begin(); it != data.end(); it++)
198 {
199 ch.add(*it);
200 }
201
202
203 unsigned int cvSize = 0;
204 for(DGtal::MelkmanConvexHull<Point, Functor>::ConstIterator it = ch.begin(); it != ch.end(); it++, cvSize++)
205 {
206 trace.info() << *it ;
207 };
208 trace.info() << std::endl;
209
210 if(res.size() == cvSize && (cvSize == ch.size()))
211 nbok++;
212 nb++;
213
214 trace.info() << "(" << nbok << "/" << nb << ") " << endl;
215 // test copy and [] operator on convex hull:
216 trace.info() << "test copy and [] operator on convex hull:" << std::endl;
217 DGtal::MelkmanConvexHull<Point, Functor> ch2 = ch;
218 unsigned int cvSize2 = 0;
219 for(DGtal::MelkmanConvexHull<Point, Functor>::ConstIterator it = ch.begin(); it != ch.end(); it++, cvSize2++)
220 {
221 trace.info() << *it ;
222 };
223 if(res.size() == cvSize2 && ch[0] == ch2[0])
224 nbok++;
225 nb++;
226 trace.info() << "(" << nbok << "/" << nb << ") " << endl;
227 trace.endBlock();
228
229
230 trace.beginBlock ( "Random Tests..." );
231 vector<Point> randomData, res1, res2;
232 const int numberOfPoints = 1000;
233 const int numberOfTries = 50;
234
235 for (int i = 0; ( (i < numberOfTries)&&(nbok == nb) ); i++)
236 {
237 //new data
238 randomData.clear();
239 res1.clear();
240 res2.clear();
241 for (int j = 0; j < numberOfPoints; j++)
242 randomData.push_back( Point(rand()%256, rand()%256) );
243 //computation
244 andrewConvexHullAlgorithm( randomData.begin(), randomData.end(), back_inserter( res1 ), predicate );
245 grahamConvexHullAlgorithm( randomData.begin(), randomData.end(), back_inserter( res2 ), predicate, comparator );
246 //comparison
247 if ( (res1.size() == res2.size()) &&
248 (circularlyEqual(res1.begin(), res1.end(), res2.begin(), res2.end())) )
249 nbok++;
250 nb++;
251 trace.info() << "(" << nbok << "/" << nb << ") " << endl;
252 //another computation
253 res2.clear();
254 sort( randomData.begin(), randomData.end() );
255 melkmanConvexHullAlgorithm( randomData.begin(), randomData.end(), back_inserter( res2 ), functor );
256 //comparison
257 if ( (res1.size() == res2.size()) &&
258 (circularlyEqual(res1.begin(), res1.end(), res2.begin(), res2.end())) )
259 nbok++;
260 nb++;
261 trace.info() << "(" << nbok << "/" << nb << ") " << endl;
262 }
263
264 trace.endBlock();
265
266 return nbok == nb;
267}
268
269
274bool testConvexHullCompThickness()
275{
276 unsigned int nbok = 0;
277 unsigned int nb = 0;
278 typedef PointVector<2,DGtal::int32_t> Point;
279 typedef InHalfPlaneBySimple3x3Matrix<Point, DGtal::int32_t> Functor;
280
281 trace.beginBlock ( "One simple test..." );
282 DGtal::MelkmanConvexHull<Point, Functor> ch;
283 ch.add(Point(0,0));
284 ch.add(Point(11,1));
285 ch.add(Point(12,3));
286 ch.add(Point(8,3));
287 ch.add(Point(4,5));
288 ch.add(Point(2,6));
289 ch.add(Point(1,4));
290
291 Point antipodalP, antipodalQ, antipodalS;
292 double thicknessE = DGtal::functions::Hull2D::computeHullThickness(ch.begin(), ch.end(),
293 DGtal::functions::Hull2D::EuclideanThickness,
294 antipodalP,
295 antipodalQ,
296 antipodalS);
297 double thicknessHV = DGtal::functions::Hull2D::computeHullThickness(ch.begin(), ch.end(),
298 DGtal::functions::Hull2D::HorizontalVerticalThickness,
299 antipodalP, antipodalQ, antipodalS);
300 double thicknessHVb = DGtal::functions::Hull2D::computeHullThickness(ch.begin(), ch.end(),
301 DGtal::functions::Hull2D::HorizontalVerticalThickness);
302
303
304 Board2D aBoard;
305 for(DGtal::MelkmanConvexHull<Point, Functor>::ConstIterator it = ch.begin(); it != ch.end(); it++){
306 if(it != ch.end()-1)
307 aBoard.drawLine((*it)[0], (*it)[1], (*(it+1))[0], (*(it+1))[1]);
308 else{
309 aBoard.drawLine((*it)[0], (*it)[1], (*(ch.begin()))[0], (*(ch.begin()))[1]);
310 }
311 aBoard << *it;
312 }
313
314 aBoard.setPenColor(DGtal::Color::Red);
315 aBoard.drawCircle( antipodalS[0], antipodalS[1], 1.0) ;
316 aBoard.setPenColor(DGtal::Color::Blue);
317 aBoard.drawCircle(antipodalP[0], antipodalP[1], 1.0);
318 aBoard.drawCircle(antipodalQ[0], antipodalQ[1], 1.0);
319
320 aBoard.drawLine(antipodalP[0], antipodalP[1], antipodalQ[0], antipodalQ[1]);
321 double awaitedThHV = DGtal::functions::Hull2D::getThicknessAntipodalPair(Point(0,0), Point(11,1), Point(2,6), DGtal::functions::Hull2D::HorizontalVerticalThickness );
322 double awaitedThE = DGtal::functions::Hull2D::getThicknessAntipodalPair(Point(0,0), Point(11,1), Point(2,6), DGtal::functions::Hull2D::EuclideanThickness );
323 trace.info() << "Thickness HV = " << thicknessHV << std::endl;
324 trace.info() << "Expected Thickness HV = " << awaitedThHV << std::endl;
325 trace.info() << "Thickness Euclidean = " << thicknessE << std::endl;
326 trace.info() << "Expected Euclidean Thickness = " << awaitedThE << std::endl;
327 aBoard.saveEPS("testConvexHull2D_Thickness.eps");
328
329 // testing tickness after changing begin points.
330 std::vector<Z2i::RealPoint> hull {
331 {804.56832227024199, -68.471176393526704},
332 {804.96020257363512, -69.494933490400683},
333 {805.35232247974727, -70.519316530915930},
334 {825.15497274857830, -114.34711249674335},
335 {828.67425867587508, -121.69593677670120},
336 {829.05943325037651, -122.39546351563243},
337 {829.42436256177677, -122.73833140123513},
338 {827.82353937509288, -118.87280445059109},
339 {811.06291230576301, -82.124832029926566},
340 {806.83483124583609, -72.890457996792406},
341 {804.92531970702396, -68.803808853026638},
342 {804.55092650722997, -68.125792709291431},
343 };
344
345 double th = DGtal::functions::Hull2D::computeHullThickness(hull.begin(), hull.end(),
346 DGtal::functions::Hull2D::HorizontalVerticalThickness);
347 std::rotate(hull.begin(), hull.begin() + 1, hull.end());
348 double th2 = DGtal::functions::Hull2D::computeHullThickness(hull.begin(), hull.end(),
349 DGtal::functions::Hull2D::HorizontalVerticalThickness);
350 std::rotate(hull.begin(), hull.begin() + 5, hull.end());
351 double th3 = DGtal::functions::Hull2D::computeHullThickness(hull.begin(), hull.end(),
352 DGtal::functions::Hull2D::HorizontalVerticalThickness);
353 trace.info() << "Thickness (before change init point) = " << th << std::endl;
354 trace.info() << "Thickness (after change init point) = " << th2 << std::endl;
355 trace.info() << "Thickness (after change init point) = " << th3 << std::endl;
356
357 nbok += thicknessHV == awaitedThHV && thicknessE == awaitedThE &&
358 thicknessHVb == thicknessHV && abs(th - th2) < 0.000001 &&
359 th - 0.604414 < 0.000001 && abs(th - th3) < 0.000001 ;
360 nb++;
361 return nb==nbok;
362}
363
364
365
366
368// Standard services - public :
369
370int main( int argc, char** argv )
371{
372 trace.beginBlock ( "Testing convex hull function" );
373 trace.info() << "Args:";
374 for ( int i = 0; i < argc; ++i )
375 trace.info() << " " << argv[ i ];
376 trace.info() << endl;
377
378 bool res = testConvexHull2D() &&testConvexHullCompThickness();
379 trace.emphase() << ( res ? "Passed." : "Error." ) << endl;
380 trace.endBlock();
381 return res ? 0 : 1;
382}
383// //
DGtal::Board2D
Aim: This class specializes a 'Board' class so as to display DGtal objects more naturally (with <<)....
Definition: Board2D.h:71
DGtal::Circulator
Aim: Provides an adapter for classical iterators that can iterate through the underlying data structu...
Definition: Circulator.h:86
DGtal::Color::Red
static const Color Red
Definition: Color.h:416
DGtal::Color::Blue
static const Color Blue
Definition: Color.h:419
DGtal::InHalfPlaneBySimple3x3Matrix
Aim: Class that implements an orientation functor, ie. it provides a way to compute the orientation o...
Definition: InHalfPlaneBySimple3x3Matrix.h:92
DGtal::MelkmanConvexHull
Aim: This class implements the on-line algorithm of Melkman for the computation of the convex hull of...
Definition: MelkmanConvexHull.h:90
DGtal::MelkmanConvexHull::ConstIterator
std::deque< Point >::const_iterator ConstIterator
Definition: MelkmanConvexHull.h:123
DGtal::PointVector
Aim: Implements basic operations that will be used in Point and Vector classes.
Definition: PointVector.h:593
DGtal::PredicateFromOrientationFunctor2
Aim: Small adapter to models of COrientationFunctor2. It is a model of concepts::CPointPredicate....
Definition: PredicateFromOrientationFunctor2.h:82
DGtal::Trace::beginBlock
void beginBlock(const std::string &keyword="")
DGtal::Trace::emphase
std::ostream & emphase()
DGtal::Trace::info
std::ostream & info()
DGtal::Trace::endBlock
double endBlock()
DGtal::functors::PolarPointComparatorBy2x2DetComputer
Aim: Class that implements a binary point predicate, which is able to compare the position of two giv...
Definition: PolarPointComparatorBy2x2DetComputer.h:80
LibBoard::Board::setPenColor
Board & setPenColor(const DGtal::Color &color)
Definition: Board.cpp:298
LibBoard::Board::drawLine
void drawLine(double x1, double y1, double x2, double y2, int depthValue=-1)
Definition: Board.cpp:368
LibBoard::Board::drawCircle
void drawCircle(double x, double y, double radius, int depthValue=-1)
Definition: Board.cpp:451
LibBoard::Board::saveEPS
void saveEPS(const char *filename, PageSize size=Board::BoundingBox, double margin=10.0) const
Definition: Board.cpp:805
DGtal::functions::Hull2D::getThicknessAntipodalPair
double getThicknessAntipodalPair(const TInputPoint &p, const TInputPoint &q, const TInputPoint &r, const ThicknessDefinition &def)
DGtal::functions::Hull2D::computeHullThickness
double computeHullThickness(const ForwardIterator &itb, const ForwardIterator &ite, const ThicknessDefinition &def)
Procedure to compute the convex hull thickness given from different definitions (Horizontal/vertical ...
DGtal
DGtal is the top-level namespace which contains all DGtal functions and types.
Definition: ClosedIntegerHalfPlane.h:49
DGtal::trace
Trace trace
Definition: Common.h:154
std
STL namespace.
main
int main()
Definition: testBits.cpp:56
Point
MyPointD Point
Definition: testClone2.cpp:383
ch
DGtal::MelkmanConvexHull< Point, Functor > ch
Definition: testConvexHull2DReverse.cpp:54
Functor
InHalfPlaneBySimple3x3Matrix< Point, double > Functor
Definition: testConvexHull2DReverse.cpp:53
thicknessE
const double thicknessE
Definition: testConvexHull2DReverse.cpp:79
thicknessHV
const double thicknessHV
Definition: testConvexHull2DReverse.cpp:75
circularlyEqual
bool circularlyEqual(const ForwardIterator &first1, const ForwardIterator &last1, const ForwardIterator &first2, const ForwardIterator &last2)
Definition: testConvexHull2D.cpp:66
testConvexHullCompThickness
bool testConvexHullCompThickness()
Definition: testConvexHull2D.cpp:274
testConvexHull2D
bool testConvexHull2D()
Definition: testConvexHull2D.cpp:102
Generated on Sun Nov 27 2022 15:15:55 for DGtal by doxygen 1.9.5