1 /*M///////////////////////////////////////////////////////////////////////////////////////
2 //
3 // IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
4 //
5 // By downloading, copying, installing or using the software you agree to this license.
6 // If you do not agree to this license, do not download, install,
7 // copy or use the software.
8 //
9 //
10 // Intel License Agreement
11 // For Open Source Computer Vision Library
12 //
13 // Copyright (C) 2000, Intel Corporation, all rights reserved.
14 // Third party copyrights are property of their respective owners.
15 //
16 // Redistribution and use in source and binary forms, with or without modification,
17 // are permitted provided that the following conditions are met:
18 //
19 // * Redistribution's of source code must retain the above copyright notice,
20 // this list of conditions and the following disclaimer.
21 //
22 // * Redistribution's in binary form must reproduce the above copyright notice,
23 // this list of conditions and the following disclaimer in the documentation
24 // and/or other materials provided with the distribution.
25 //
26 // * The name of Intel Corporation may not be used to endorse or promote products
27 // derived from this software without specific prior written permission.
28 //
29 // This software is provided by the copyright holders and contributors "as is" and
30 // any express or implied warranties, including, but not limited to, the implied
31 // warranties of merchantability and fitness for a particular purpose are disclaimed.
32 // In no event shall the Intel Corporation or contributors be liable for any direct,
33 // indirect, incidental, special, exemplary, or consequential damages
34 // (including, but not limited to, procurement of substitute goods or services;
35 // loss of use, data, or profits; or business interruption) however caused
36 // and on any theory of liability, whether in contract, strict liability,
37 // or tort (including negligence or otherwise) arising in any way out of
38 // the use of this software, even if advised of the possibility of such damage.
39 //
40 //M*/
41 #include "precomp.hpp"
42
43 namespace cv
44 {
45
nextEdge(int edge) const46 int Subdiv2D::nextEdge(int edge) const
47 {
48 CV_DbgAssert((size_t)(edge >> 2) < qedges.size());
49 return qedges[edge >> 2].next[edge & 3];
50 }
51
rotateEdge(int edge,int rotate) const52 int Subdiv2D::rotateEdge(int edge, int rotate) const
53 {
54 return (edge & ~3) + ((edge + rotate) & 3);
55 }
56
symEdge(int edge) const57 int Subdiv2D::symEdge(int edge) const
58 {
59 return edge ^ 2;
60 }
61
getEdge(int edge,int nextEdgeType) const62 int Subdiv2D::getEdge(int edge, int nextEdgeType) const
63 {
64 CV_DbgAssert((size_t)(edge >> 2) < qedges.size());
65 edge = qedges[edge >> 2].next[(edge + nextEdgeType) & 3];
66 return (edge & ~3) + ((edge + (nextEdgeType >> 4)) & 3);
67 }
68
edgeOrg(int edge,CV_OUT Point2f * orgpt) const69 int Subdiv2D::edgeOrg(int edge, CV_OUT Point2f* orgpt) const
70 {
71 CV_DbgAssert((size_t)(edge >> 2) < qedges.size());
72 int vidx = qedges[edge >> 2].pt[edge & 3];
73 if( orgpt )
74 {
75 CV_DbgAssert((size_t)vidx < vtx.size());
76 *orgpt = vtx[vidx].pt;
77 }
78 return vidx;
79 }
80
edgeDst(int edge,CV_OUT Point2f * dstpt) const81 int Subdiv2D::edgeDst(int edge, CV_OUT Point2f* dstpt) const
82 {
83 CV_DbgAssert((size_t)(edge >> 2) < qedges.size());
84 int vidx = qedges[edge >> 2].pt[(edge + 2) & 3];
85 if( dstpt )
86 {
87 CV_DbgAssert((size_t)vidx < vtx.size());
88 *dstpt = vtx[vidx].pt;
89 }
90 return vidx;
91 }
92
93
getVertex(int vertex,CV_OUT int * firstEdge) const94 Point2f Subdiv2D::getVertex(int vertex, CV_OUT int* firstEdge) const
95 {
96 CV_DbgAssert((size_t)vertex < vtx.size());
97 if( firstEdge )
98 *firstEdge = vtx[vertex].firstEdge;
99 return vtx[vertex].pt;
100 }
101
102
Subdiv2D()103 Subdiv2D::Subdiv2D()
104 {
105 validGeometry = false;
106 freeQEdge = 0;
107 freePoint = 0;
108 recentEdge = 0;
109 }
110
Subdiv2D(Rect rect)111 Subdiv2D::Subdiv2D(Rect rect)
112 {
113 validGeometry = false;
114 freeQEdge = 0;
115 freePoint = 0;
116 recentEdge = 0;
117
118 initDelaunay(rect);
119 }
120
121
QuadEdge()122 Subdiv2D::QuadEdge::QuadEdge()
123 {
124 next[0] = next[1] = next[2] = next[3] = 0;
125 pt[0] = pt[1] = pt[2] = pt[3] = 0;
126 }
127
QuadEdge(int edgeidx)128 Subdiv2D::QuadEdge::QuadEdge(int edgeidx)
129 {
130 CV_DbgAssert((edgeidx & 3) == 0);
131 next[0] = edgeidx;
132 next[1] = edgeidx+3;
133 next[2] = edgeidx+2;
134 next[3] = edgeidx+1;
135
136 pt[0] = pt[1] = pt[2] = pt[3] = 0;
137 }
138
isfree() const139 bool Subdiv2D::QuadEdge::isfree() const
140 {
141 return next[0] <= 0;
142 }
143
Vertex()144 Subdiv2D::Vertex::Vertex()
145 {
146 firstEdge = 0;
147 type = -1;
148 }
149
Vertex(Point2f _pt,bool _isvirtual,int _firstEdge)150 Subdiv2D::Vertex::Vertex(Point2f _pt, bool _isvirtual, int _firstEdge)
151 {
152 firstEdge = _firstEdge;
153 type = (int)_isvirtual;
154 pt = _pt;
155 }
156
isvirtual() const157 bool Subdiv2D::Vertex::isvirtual() const
158 {
159 return type > 0;
160 }
161
isfree() const162 bool Subdiv2D::Vertex::isfree() const
163 {
164 return type < 0;
165 }
166
splice(int edgeA,int edgeB)167 void Subdiv2D::splice( int edgeA, int edgeB )
168 {
169 int& a_next = qedges[edgeA >> 2].next[edgeA & 3];
170 int& b_next = qedges[edgeB >> 2].next[edgeB & 3];
171 int a_rot = rotateEdge(a_next, 1);
172 int b_rot = rotateEdge(b_next, 1);
173 int& a_rot_next = qedges[a_rot >> 2].next[a_rot & 3];
174 int& b_rot_next = qedges[b_rot >> 2].next[b_rot & 3];
175 std::swap(a_next, b_next);
176 std::swap(a_rot_next, b_rot_next);
177 }
178
setEdgePoints(int edge,int orgPt,int dstPt)179 void Subdiv2D::setEdgePoints(int edge, int orgPt, int dstPt)
180 {
181 qedges[edge >> 2].pt[edge & 3] = orgPt;
182 qedges[edge >> 2].pt[(edge + 2) & 3] = dstPt;
183 vtx[orgPt].firstEdge = edge;
184 vtx[dstPt].firstEdge = edge ^ 2;
185 }
186
connectEdges(int edgeA,int edgeB)187 int Subdiv2D::connectEdges( int edgeA, int edgeB )
188 {
189 int edge = newEdge();
190
191 splice(edge, getEdge(edgeA, NEXT_AROUND_LEFT));
192 splice(symEdge(edge), edgeB);
193
194 setEdgePoints(edge, edgeDst(edgeA), edgeOrg(edgeB));
195 return edge;
196 }
197
swapEdges(int edge)198 void Subdiv2D::swapEdges( int edge )
199 {
200 int sedge = symEdge(edge);
201 int a = getEdge(edge, PREV_AROUND_ORG);
202 int b = getEdge(sedge, PREV_AROUND_ORG);
203
204 splice(edge, a);
205 splice(sedge, b);
206
207 setEdgePoints(edge, edgeDst(a), edgeDst(b));
208
209 splice(edge, getEdge(a, NEXT_AROUND_LEFT));
210 splice(sedge, getEdge(b, NEXT_AROUND_LEFT));
211 }
212
triangleArea(Point2f a,Point2f b,Point2f c)213 static double triangleArea( Point2f a, Point2f b, Point2f c )
214 {
215 return ((double)b.x - a.x) * ((double)c.y - a.y) - ((double)b.y - a.y) * ((double)c.x - a.x);
216 }
217
isRightOf(Point2f pt,int edge) const218 int Subdiv2D::isRightOf(Point2f pt, int edge) const
219 {
220 Point2f org, dst;
221 edgeOrg(edge, &org);
222 edgeDst(edge, &dst);
223 double cw_area = triangleArea( pt, dst, org );
224
225 return (cw_area > 0) - (cw_area < 0);
226 }
227
newEdge()228 int Subdiv2D::newEdge()
229 {
230 if( freeQEdge <= 0 )
231 {
232 qedges.push_back(QuadEdge());
233 freeQEdge = (int)(qedges.size()-1);
234 }
235 int edge = freeQEdge*4;
236 freeQEdge = qedges[edge >> 2].next[1];
237 qedges[edge >> 2] = QuadEdge(edge);
238 return edge;
239 }
240
deleteEdge(int edge)241 void Subdiv2D::deleteEdge(int edge)
242 {
243 CV_DbgAssert((size_t)(edge >> 2) < (size_t)qedges.size());
244 splice( edge, getEdge(edge, PREV_AROUND_ORG) );
245 int sedge = symEdge(edge);
246 splice(sedge, getEdge(sedge, PREV_AROUND_ORG) );
247
248 edge >>= 2;
249 qedges[edge].next[0] = 0;
250 qedges[edge].next[1] = freeQEdge;
251 freeQEdge = edge;
252 }
253
newPoint(Point2f pt,bool isvirtual,int firstEdge)254 int Subdiv2D::newPoint(Point2f pt, bool isvirtual, int firstEdge)
255 {
256 if( freePoint == 0 )
257 {
258 vtx.push_back(Vertex());
259 freePoint = (int)(vtx.size()-1);
260 }
261 int vidx = freePoint;
262 freePoint = vtx[vidx].firstEdge;
263 vtx[vidx] = Vertex(pt, isvirtual, firstEdge);
264
265 return vidx;
266 }
267
deletePoint(int vidx)268 void Subdiv2D::deletePoint(int vidx)
269 {
270 CV_DbgAssert( (size_t)vidx < vtx.size() );
271 vtx[vidx].firstEdge = freePoint;
272 vtx[vidx].type = -1;
273 freePoint = vidx;
274 }
275
locate(Point2f pt,int & _edge,int & _vertex)276 int Subdiv2D::locate(Point2f pt, int& _edge, int& _vertex)
277 {
278 int vertex = 0;
279
280 int i, maxEdges = (int)(qedges.size() * 4);
281
282 if( qedges.size() < (size_t)4 )
283 CV_Error( CV_StsError, "Subdivision is empty" );
284
285 if( pt.x < topLeft.x || pt.y < topLeft.y || pt.x >= bottomRight.x || pt.y >= bottomRight.y )
286 CV_Error( CV_StsOutOfRange, "" );
287
288 int edge = recentEdge;
289 CV_Assert(edge > 0);
290
291 int location = PTLOC_ERROR;
292
293 int right_of_curr = isRightOf(pt, edge);
294 if( right_of_curr > 0 )
295 {
296 edge = symEdge(edge);
297 right_of_curr = -right_of_curr;
298 }
299
300 for( i = 0; i < maxEdges; i++ )
301 {
302 int onext_edge = nextEdge( edge );
303 int dprev_edge = getEdge( edge, PREV_AROUND_DST );
304
305 int right_of_onext = isRightOf( pt, onext_edge );
306 int right_of_dprev = isRightOf( pt, dprev_edge );
307
308 if( right_of_dprev > 0 )
309 {
310 if( right_of_onext > 0 || (right_of_onext == 0 && right_of_curr == 0) )
311 {
312 location = PTLOC_INSIDE;
313 break;
314 }
315 else
316 {
317 right_of_curr = right_of_onext;
318 edge = onext_edge;
319 }
320 }
321 else
322 {
323 if( right_of_onext > 0 )
324 {
325 if( right_of_dprev == 0 && right_of_curr == 0 )
326 {
327 location = PTLOC_INSIDE;
328 break;
329 }
330 else
331 {
332 right_of_curr = right_of_dprev;
333 edge = dprev_edge;
334 }
335 }
336 else if( right_of_curr == 0 &&
337 isRightOf( vtx[edgeDst(onext_edge)].pt, edge ) >= 0 )
338 {
339 edge = symEdge( edge );
340 }
341 else
342 {
343 right_of_curr = right_of_onext;
344 edge = onext_edge;
345 }
346 }
347 }
348
349 recentEdge = edge;
350
351 if( location == PTLOC_INSIDE )
352 {
353 Point2f org_pt, dst_pt;
354 edgeOrg(edge, &org_pt);
355 edgeDst(edge, &dst_pt);
356
357 double t1 = fabs( pt.x - org_pt.x );
358 t1 += fabs( pt.y - org_pt.y );
359 double t2 = fabs( pt.x - dst_pt.x );
360 t2 += fabs( pt.y - dst_pt.y );
361 double t3 = fabs( org_pt.x - dst_pt.x );
362 t3 += fabs( org_pt.y - dst_pt.y );
363
364 if( t1 < FLT_EPSILON )
365 {
366 location = PTLOC_VERTEX;
367 vertex = edgeOrg( edge );
368 edge = 0;
369 }
370 else if( t2 < FLT_EPSILON )
371 {
372 location = PTLOC_VERTEX;
373 vertex = edgeDst( edge );
374 edge = 0;
375 }
376 else if( (t1 < t3 || t2 < t3) &&
377 fabs( triangleArea( pt, org_pt, dst_pt )) < FLT_EPSILON )
378 {
379 location = PTLOC_ON_EDGE;
380 vertex = 0;
381 }
382 }
383
384 if( location == PTLOC_ERROR )
385 {
386 edge = 0;
387 vertex = 0;
388 }
389
390 _edge = edge;
391 _vertex = vertex;
392
393 return location;
394 }
395
396
397 inline int
isPtInCircle3(Point2f pt,Point2f a,Point2f b,Point2f c)398 isPtInCircle3( Point2f pt, Point2f a, Point2f b, Point2f c)
399 {
400 const double eps = FLT_EPSILON*0.125;
401 double val = ((double)a.x * a.x + (double)a.y * a.y) * triangleArea( b, c, pt );
402 val -= ((double)b.x * b.x + (double)b.y * b.y) * triangleArea( a, c, pt );
403 val += ((double)c.x * c.x + (double)c.y * c.y) * triangleArea( a, b, pt );
404 val -= ((double)pt.x * pt.x + (double)pt.y * pt.y) * triangleArea( a, b, c );
405
406 return val > eps ? 1 : val < -eps ? -1 : 0;
407 }
408
409
insert(Point2f pt)410 int Subdiv2D::insert(Point2f pt)
411 {
412 int curr_point = 0, curr_edge = 0, deleted_edge = 0;
413 int location = locate( pt, curr_edge, curr_point );
414
415 if( location == PTLOC_ERROR )
416 CV_Error( CV_StsBadSize, "" );
417
418 if( location == PTLOC_OUTSIDE_RECT )
419 CV_Error( CV_StsOutOfRange, "" );
420
421 if( location == PTLOC_VERTEX )
422 return curr_point;
423
424 if( location == PTLOC_ON_EDGE )
425 {
426 deleted_edge = curr_edge;
427 recentEdge = curr_edge = getEdge( curr_edge, PREV_AROUND_ORG );
428 deleteEdge(deleted_edge);
429 }
430 else if( location == PTLOC_INSIDE )
431 ;
432 else
433 CV_Error_(CV_StsError, ("Subdiv2D::locate returned invalid location = %d", location) );
434
435 assert( curr_edge != 0 );
436 validGeometry = false;
437
438 curr_point = newPoint(pt, false);
439 int base_edge = newEdge();
440 int first_point = edgeOrg(curr_edge);
441 setEdgePoints(base_edge, first_point, curr_point);
442 splice(base_edge, curr_edge);
443
444 do
445 {
446 base_edge = connectEdges( curr_edge, symEdge(base_edge) );
447 curr_edge = getEdge(base_edge, PREV_AROUND_ORG);
448 }
449 while( edgeDst(curr_edge) != first_point );
450
451 curr_edge = getEdge( base_edge, PREV_AROUND_ORG );
452
453 int i, max_edges = (int)(qedges.size()*4);
454
455 for( i = 0; i < max_edges; i++ )
456 {
457 int temp_dst = 0, curr_org = 0, curr_dst = 0;
458 int temp_edge = getEdge( curr_edge, PREV_AROUND_ORG );
459
460 temp_dst = edgeDst( temp_edge );
461 curr_org = edgeOrg( curr_edge );
462 curr_dst = edgeDst( curr_edge );
463
464 if( isRightOf( vtx[temp_dst].pt, curr_edge ) > 0 &&
465 isPtInCircle3( vtx[curr_org].pt, vtx[temp_dst].pt,
466 vtx[curr_dst].pt, vtx[curr_point].pt ) < 0 )
467 {
468 swapEdges( curr_edge );
469 curr_edge = getEdge( curr_edge, PREV_AROUND_ORG );
470 }
471 else if( curr_org == first_point )
472 break;
473 else
474 curr_edge = getEdge( nextEdge( curr_edge ), PREV_AROUND_LEFT );
475 }
476
477 return curr_point;
478 }
479
insert(const std::vector<Point2f> & ptvec)480 void Subdiv2D::insert(const std::vector<Point2f>& ptvec)
481 {
482 for( size_t i = 0; i < ptvec.size(); i++ )
483 insert(ptvec[i]);
484 }
485
initDelaunay(Rect rect)486 void Subdiv2D::initDelaunay( Rect rect )
487 {
488 float big_coord = 3.f * MAX( rect.width, rect.height );
489 float rx = (float)rect.x;
490 float ry = (float)rect.y;
491
492 vtx.clear();
493 qedges.clear();
494
495 recentEdge = 0;
496 validGeometry = false;
497
498 topLeft = Point2f( rx, ry );
499 bottomRight = Point2f( rx + rect.width, ry + rect.height );
500
501 Point2f ppA( rx + big_coord, ry );
502 Point2f ppB( rx, ry + big_coord );
503 Point2f ppC( rx - big_coord, ry - big_coord );
504
505 vtx.push_back(Vertex());
506 qedges.push_back(QuadEdge());
507
508 freeQEdge = 0;
509 freePoint = 0;
510
511 int pA = newPoint(ppA, false);
512 int pB = newPoint(ppB, false);
513 int pC = newPoint(ppC, false);
514
515 int edge_AB = newEdge();
516 int edge_BC = newEdge();
517 int edge_CA = newEdge();
518
519 setEdgePoints( edge_AB, pA, pB );
520 setEdgePoints( edge_BC, pB, pC );
521 setEdgePoints( edge_CA, pC, pA );
522
523 splice( edge_AB, symEdge( edge_CA ));
524 splice( edge_BC, symEdge( edge_AB ));
525 splice( edge_CA, symEdge( edge_BC ));
526
527 recentEdge = edge_AB;
528 }
529
530
clearVoronoi()531 void Subdiv2D::clearVoronoi()
532 {
533 size_t i, total = qedges.size();
534
535 for( i = 0; i < total; i++ )
536 qedges[i].pt[1] = qedges[i].pt[3] = 0;
537
538 total = vtx.size();
539 for( i = 0; i < total; i++ )
540 {
541 if( vtx[i].isvirtual() )
542 deletePoint((int)i);
543 }
544
545 validGeometry = false;
546 }
547
548
computeVoronoiPoint(Point2f org0,Point2f dst0,Point2f org1,Point2f dst1)549 static Point2f computeVoronoiPoint(Point2f org0, Point2f dst0, Point2f org1, Point2f dst1)
550 {
551 double a0 = dst0.x - org0.x;
552 double b0 = dst0.y - org0.y;
553 double c0 = -0.5*(a0 * (dst0.x + org0.x) + b0 * (dst0.y + org0.y));
554
555 double a1 = dst1.x - org1.x;
556 double b1 = dst1.y - org1.y;
557 double c1 = -0.5*(a1 * (dst1.x + org1.x) + b1 * (dst1.y + org1.y));
558
559 double det = a0 * b1 - a1 * b0;
560
561 if( det != 0 )
562 {
563 det = 1. / det;
564 return Point2f((float) ((b0 * c1 - b1 * c0) * det),
565 (float) ((a1 * c0 - a0 * c1) * det));
566 }
567
568 return Point2f(FLT_MAX, FLT_MAX);
569 }
570
571
calcVoronoi()572 void Subdiv2D::calcVoronoi()
573 {
574 // check if it is already calculated
575 if( validGeometry )
576 return;
577
578 clearVoronoi();
579 int i, total = (int)qedges.size();
580
581 // loop through all quad-edges, except for the first 3 (#1, #2, #3 - 0 is reserved for "NULL" pointer)
582 for( i = 4; i < total; i++ )
583 {
584 QuadEdge& quadedge = qedges[i];
585
586 if( quadedge.isfree() )
587 continue;
588
589 int edge0 = (int)(i*4);
590 Point2f org0, dst0, org1, dst1;
591
592 if( !quadedge.pt[3] )
593 {
594 int edge1 = getEdge( edge0, NEXT_AROUND_LEFT );
595 int edge2 = getEdge( edge1, NEXT_AROUND_LEFT );
596
597 edgeOrg(edge0, &org0);
598 edgeDst(edge0, &dst0);
599 edgeOrg(edge1, &org1);
600 edgeDst(edge1, &dst1);
601
602 Point2f virt_point = computeVoronoiPoint(org0, dst0, org1, dst1);
603
604 if( fabs( virt_point.x ) < FLT_MAX * 0.5 &&
605 fabs( virt_point.y ) < FLT_MAX * 0.5 )
606 {
607 quadedge.pt[3] = qedges[edge1 >> 2].pt[3 - (edge1 & 2)] =
608 qedges[edge2 >> 2].pt[3 - (edge2 & 2)] = newPoint(virt_point, true);
609 }
610 }
611
612 if( !quadedge.pt[1] )
613 {
614 int edge1 = getEdge( edge0, NEXT_AROUND_RIGHT );
615 int edge2 = getEdge( edge1, NEXT_AROUND_RIGHT );
616
617 edgeOrg(edge0, &org0);
618 edgeDst(edge0, &dst0);
619 edgeOrg(edge1, &org1);
620 edgeDst(edge1, &dst1);
621
622 Point2f virt_point = computeVoronoiPoint(org0, dst0, org1, dst1);
623
624 if( fabs( virt_point.x ) < FLT_MAX * 0.5 &&
625 fabs( virt_point.y ) < FLT_MAX * 0.5 )
626 {
627 quadedge.pt[1] = qedges[edge1 >> 2].pt[1 + (edge1 & 2)] =
628 qedges[edge2 >> 2].pt[1 + (edge2 & 2)] = newPoint(virt_point, true);
629 }
630 }
631 }
632
633 validGeometry = true;
634 }
635
636
637 static int
isRightOf2(const Point2f & pt,const Point2f & org,const Point2f & diff)638 isRightOf2( const Point2f& pt, const Point2f& org, const Point2f& diff )
639 {
640 double cw_area = ((double)org.x - pt.x)*diff.y - ((double)org.y - pt.y)*diff.x;
641 return (cw_area > 0) - (cw_area < 0);
642 }
643
644
findNearest(Point2f pt,Point2f * nearestPt)645 int Subdiv2D::findNearest(Point2f pt, Point2f* nearestPt)
646 {
647 if( !validGeometry )
648 calcVoronoi();
649
650 int vertex = 0, edge = 0;
651 int loc = locate( pt, edge, vertex );
652
653 if( loc != PTLOC_ON_EDGE && loc != PTLOC_INSIDE )
654 return vertex;
655
656 vertex = 0;
657
658 Point2f start;
659 edgeOrg(edge, &start);
660 Point2f diff = pt - start;
661
662 edge = rotateEdge(edge, 1);
663
664 int i, total = (int)vtx.size();
665
666 for( i = 0; i < total; i++ )
667 {
668 Point2f t;
669
670 for(;;)
671 {
672 CV_Assert( edgeDst(edge, &t) > 0 );
673 if( isRightOf2( t, start, diff ) >= 0 )
674 break;
675
676 edge = getEdge( edge, NEXT_AROUND_LEFT );
677 }
678
679 for(;;)
680 {
681 CV_Assert( edgeOrg( edge, &t ) > 0 );
682
683 if( isRightOf2( t, start, diff ) < 0 )
684 break;
685
686 edge = getEdge( edge, PREV_AROUND_LEFT );
687 }
688
689 Point2f tempDiff;
690 edgeDst(edge, &tempDiff);
691 edgeOrg(edge, &t);
692 tempDiff -= t;
693
694 if( isRightOf2( pt, t, tempDiff ) >= 0 )
695 {
696 vertex = edgeOrg(rotateEdge( edge, 3 ));
697 break;
698 }
699
700 edge = symEdge( edge );
701 }
702
703 if( nearestPt && vertex > 0 )
704 *nearestPt = vtx[vertex].pt;
705
706 return vertex;
707 }
708
getEdgeList(std::vector<Vec4f> & edgeList) const709 void Subdiv2D::getEdgeList(std::vector<Vec4f>& edgeList) const
710 {
711 edgeList.clear();
712
713 for( size_t i = 4; i < qedges.size(); i++ )
714 {
715 if( qedges[i].isfree() )
716 continue;
717 if( qedges[i].pt[0] > 0 && qedges[i].pt[2] > 0 )
718 {
719 Point2f org = vtx[qedges[i].pt[0]].pt;
720 Point2f dst = vtx[qedges[i].pt[2]].pt;
721 edgeList.push_back(Vec4f(org.x, org.y, dst.x, dst.y));
722 }
723 }
724 }
725
getTriangleList(std::vector<Vec6f> & triangleList) const726 void Subdiv2D::getTriangleList(std::vector<Vec6f>& triangleList) const
727 {
728 triangleList.clear();
729 int i, total = (int)(qedges.size()*4);
730 std::vector<bool> edgemask(total, false);
731
732 for( i = 4; i < total; i += 2 )
733 {
734 if( edgemask[i] )
735 continue;
736 Point2f a, b, c;
737 int edge = i;
738 edgeOrg(edge, &a);
739 edgemask[edge] = true;
740 edge = getEdge(edge, NEXT_AROUND_LEFT);
741 edgeOrg(edge, &b);
742 edgemask[edge] = true;
743 edge = getEdge(edge, NEXT_AROUND_LEFT);
744 edgeOrg(edge, &c);
745 edgemask[edge] = true;
746 triangleList.push_back(Vec6f(a.x, a.y, b.x, b.y, c.x, c.y));
747 }
748 }
749
getVoronoiFacetList(const std::vector<int> & idx,CV_OUT std::vector<std::vector<Point2f>> & facetList,CV_OUT std::vector<Point2f> & facetCenters)750 void Subdiv2D::getVoronoiFacetList(const std::vector<int>& idx,
751 CV_OUT std::vector<std::vector<Point2f> >& facetList,
752 CV_OUT std::vector<Point2f>& facetCenters)
753 {
754 calcVoronoi();
755 facetList.clear();
756 facetCenters.clear();
757
758 std::vector<Point2f> buf;
759
760 size_t i, total;
761 if( idx.empty() )
762 i = 4, total = vtx.size();
763 else
764 i = 0, total = idx.size();
765
766 for( ; i < total; i++ )
767 {
768 int k = idx.empty() ? (int)i : idx[i];
769
770 if( vtx[k].isfree() || vtx[k].isvirtual() )
771 continue;
772 int edge = rotateEdge(vtx[k].firstEdge, 1), t = edge;
773
774 // gather points
775 buf.clear();
776 do
777 {
778 buf.push_back(vtx[edgeOrg(t)].pt);
779 t = getEdge( t, NEXT_AROUND_LEFT );
780 }
781 while( t != edge );
782
783 facetList.push_back(buf);
784 facetCenters.push_back(vtx[k].pt);
785 }
786 }
787
788
checkSubdiv() const789 void Subdiv2D::checkSubdiv() const
790 {
791 int i, j, total = (int)qedges.size();
792
793 for( i = 0; i < total; i++ )
794 {
795 const QuadEdge& qe = qedges[i];
796
797 if( qe.isfree() )
798 continue;
799
800 for( j = 0; j < 4; j++ )
801 {
802 int e = (int)(i*4 + j);
803 int o_next = nextEdge(e);
804 int o_prev = getEdge(e, PREV_AROUND_ORG );
805 int d_prev = getEdge(e, PREV_AROUND_DST );
806 int d_next = getEdge(e, NEXT_AROUND_DST );
807
808 // check points
809 CV_Assert( edgeOrg(e) == edgeOrg(o_next));
810 CV_Assert( edgeOrg(e) == edgeOrg(o_prev));
811 CV_Assert( edgeDst(e) == edgeDst(d_next));
812 CV_Assert( edgeDst(e) == edgeDst(d_prev));
813
814 if( j % 2 == 0 )
815 {
816 CV_Assert( edgeDst(o_next) == edgeOrg(d_prev));
817 CV_Assert( edgeDst(o_prev) == edgeOrg(d_next));
818 CV_Assert( getEdge(getEdge(getEdge(e,NEXT_AROUND_LEFT),NEXT_AROUND_LEFT),NEXT_AROUND_LEFT) == e );
819 CV_Assert( getEdge(getEdge(getEdge(e,NEXT_AROUND_RIGHT),NEXT_AROUND_RIGHT),NEXT_AROUND_RIGHT) == e);
820 }
821 }
822 }
823 }
824
825 }
826
827 /* End of file. */
828