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
44 CV_IMPL CvRect
cvMaxRect(const CvRect * rect1,const CvRect * rect2)45 cvMaxRect( const CvRect* rect1, const CvRect* rect2 )
46 {
47 if( rect1 && rect2 )
48 {
49 CvRect max_rect;
50 int a, b;
51
52 max_rect.x = a = rect1->x;
53 b = rect2->x;
54 if( max_rect.x > b )
55 max_rect.x = b;
56
57 max_rect.width = a += rect1->width;
58 b += rect2->width;
59
60 if( max_rect.width < b )
61 max_rect.width = b;
62 max_rect.width -= max_rect.x;
63
64 max_rect.y = a = rect1->y;
65 b = rect2->y;
66 if( max_rect.y > b )
67 max_rect.y = b;
68
69 max_rect.height = a += rect1->height;
70 b += rect2->height;
71
72 if( max_rect.height < b )
73 max_rect.height = b;
74 max_rect.height -= max_rect.y;
75 return max_rect;
76 }
77 else if( rect1 )
78 return *rect1;
79 else if( rect2 )
80 return *rect2;
81 else
82 return cvRect(0,0,0,0);
83 }
84
85
86 CV_IMPL void
cvBoxPoints(CvBox2D box,CvPoint2D32f pt[4])87 cvBoxPoints( CvBox2D box, CvPoint2D32f pt[4] )
88 {
89 if( !pt )
90 CV_Error( CV_StsNullPtr, "NULL vertex array pointer" );
91 cv::RotatedRect(box).points((cv::Point2f*)pt);
92 }
93
94
pointPolygonTest(InputArray _contour,Point2f pt,bool measureDist)95 double cv::pointPolygonTest( InputArray _contour, Point2f pt, bool measureDist )
96 {
97 double result = 0;
98 Mat contour = _contour.getMat();
99 int i, total = contour.checkVector(2), counter = 0;
100 int depth = contour.depth();
101 CV_Assert( total >= 0 && (depth == CV_32S || depth == CV_32F));
102
103 bool is_float = depth == CV_32F;
104 double min_dist_num = FLT_MAX, min_dist_denom = 1;
105 Point ip(cvRound(pt.x), cvRound(pt.y));
106
107 if( total == 0 )
108 return measureDist ? -DBL_MAX : -1;
109
110 const Point* cnt = contour.ptr<Point>();
111 const Point2f* cntf = (const Point2f*)cnt;
112
113 if( !is_float && !measureDist && ip.x == pt.x && ip.y == pt.y )
114 {
115 // the fastest "purely integer" branch
116 Point v0, v = cnt[total-1];
117
118 for( i = 0; i < total; i++ )
119 {
120 int dist;
121 v0 = v;
122 v = cnt[i];
123
124 if( (v0.y <= ip.y && v.y <= ip.y) ||
125 (v0.y > ip.y && v.y > ip.y) ||
126 (v0.x < ip.x && v.x < ip.x) )
127 {
128 if( ip.y == v.y && (ip.x == v.x || (ip.y == v0.y &&
129 ((v0.x <= ip.x && ip.x <= v.x) || (v.x <= ip.x && ip.x <= v0.x)))) )
130 return 0;
131 continue;
132 }
133
134 dist = (ip.y - v0.y)*(v.x - v0.x) - (ip.x - v0.x)*(v.y - v0.y);
135 if( dist == 0 )
136 return 0;
137 if( v.y < v0.y )
138 dist = -dist;
139 counter += dist > 0;
140 }
141
142 result = counter % 2 == 0 ? -1 : 1;
143 }
144 else
145 {
146 Point2f v0, v;
147 Point iv;
148
149 if( is_float )
150 {
151 v = cntf[total-1];
152 }
153 else
154 {
155 v = cnt[total-1];
156 }
157
158 if( !measureDist )
159 {
160 for( i = 0; i < total; i++ )
161 {
162 double dist;
163 v0 = v;
164 if( is_float )
165 v = cntf[i];
166 else
167 v = cnt[i];
168
169 if( (v0.y <= pt.y && v.y <= pt.y) ||
170 (v0.y > pt.y && v.y > pt.y) ||
171 (v0.x < pt.x && v.x < pt.x) )
172 {
173 if( pt.y == v.y && (pt.x == v.x || (pt.y == v0.y &&
174 ((v0.x <= pt.x && pt.x <= v.x) || (v.x <= pt.x && pt.x <= v0.x)))) )
175 return 0;
176 continue;
177 }
178
179 dist = (double)(pt.y - v0.y)*(v.x - v0.x) - (double)(pt.x - v0.x)*(v.y - v0.y);
180 if( dist == 0 )
181 return 0;
182 if( v.y < v0.y )
183 dist = -dist;
184 counter += dist > 0;
185 }
186
187 result = counter % 2 == 0 ? -1 : 1;
188 }
189 else
190 {
191 for( i = 0; i < total; i++ )
192 {
193 double dx, dy, dx1, dy1, dx2, dy2, dist_num, dist_denom = 1;
194
195 v0 = v;
196 if( is_float )
197 v = cntf[i];
198 else
199 v = cnt[i];
200
201 dx = v.x - v0.x; dy = v.y - v0.y;
202 dx1 = pt.x - v0.x; dy1 = pt.y - v0.y;
203 dx2 = pt.x - v.x; dy2 = pt.y - v.y;
204
205 if( dx1*dx + dy1*dy <= 0 )
206 dist_num = dx1*dx1 + dy1*dy1;
207 else if( dx2*dx + dy2*dy >= 0 )
208 dist_num = dx2*dx2 + dy2*dy2;
209 else
210 {
211 dist_num = (dy1*dx - dx1*dy);
212 dist_num *= dist_num;
213 dist_denom = dx*dx + dy*dy;
214 }
215
216 if( dist_num*min_dist_denom < min_dist_num*dist_denom )
217 {
218 min_dist_num = dist_num;
219 min_dist_denom = dist_denom;
220 if( min_dist_num == 0 )
221 break;
222 }
223
224 if( (v0.y <= pt.y && v.y <= pt.y) ||
225 (v0.y > pt.y && v.y > pt.y) ||
226 (v0.x < pt.x && v.x < pt.x) )
227 continue;
228
229 dist_num = dy1*dx - dx1*dy;
230 if( dy < 0 )
231 dist_num = -dist_num;
232 counter += dist_num > 0;
233 }
234
235 result = std::sqrt(min_dist_num/min_dist_denom);
236 if( counter % 2 == 0 )
237 result = -result;
238 }
239 }
240
241 return result;
242 }
243
244
245 CV_IMPL double
cvPointPolygonTest(const CvArr * _contour,CvPoint2D32f pt,int measure_dist)246 cvPointPolygonTest( const CvArr* _contour, CvPoint2D32f pt, int measure_dist )
247 {
248 cv::AutoBuffer<double> abuf;
249 cv::Mat contour = cv::cvarrToMat(_contour, false, false, 0, &abuf);
250 return cv::pointPolygonTest(contour, pt, measure_dist != 0);
251 }
252
253 /*
254 This code is described in "Computational Geometry in C" (Second Edition),
255 Chapter 7. It is not written to be comprehensible without the
256 explanation in that book.
257
258 Written by Joseph O'Rourke.
259 Last modified: December 1997
260 Questions to orourke@cs.smith.edu.
261 --------------------------------------------------------------------
262 This code is Copyright 1997 by Joseph O'Rourke. It may be freely
263 redistributed in its entirety provided that this copyright notice is
264 not removed.
265 --------------------------------------------------------------------
266 */
267
268 namespace cv
269 {
270 typedef enum { Pin, Qin, Unknown } tInFlag;
271
areaSign(Point2f a,Point2f b,Point2f c)272 static int areaSign( Point2f a, Point2f b, Point2f c )
273 {
274 static const double eps = 1e-5;
275 double area2 = (b.x - a.x) * (double)(c.y - a.y) - (c.x - a.x ) * (double)(b.y - a.y);
276 return area2 > eps ? 1 : area2 < -eps ? -1 : 0;
277 }
278
279 //---------------------------------------------------------------------
280 // Returns true iff point c lies on the closed segement ab.
281 // Assumes it is already known that abc are collinear.
282 //---------------------------------------------------------------------
between(Point2f a,Point2f b,Point2f c)283 static bool between( Point2f a, Point2f b, Point2f c )
284 {
285 Point2f ba, ca;
286
287 // If ab not vertical, check betweenness on x; else on y.
288 if ( a.x != b.x )
289 return ((a.x <= c.x) && (c.x <= b.x)) ||
290 ((a.x >= c.x) && (c.x >= b.x));
291 else
292 return ((a.y <= c.y) && (c.y <= b.y)) ||
293 ((a.y >= c.y) && (c.y >= b.y));
294 }
295
parallelInt(Point2f a,Point2f b,Point2f c,Point2f d,Point2f & p,Point2f & q)296 static char parallelInt( Point2f a, Point2f b, Point2f c, Point2f d, Point2f& p, Point2f& q )
297 {
298 char code = 'e';
299 if( areaSign(a, b, c) != 0 )
300 code = '0';
301 else if( between(a, b, c) && between(a, b, d))
302 p = c, q = d;
303 else if( between(c, d, a) && between(c, d, b))
304 p = a, q = b;
305 else if( between(a, b, c) && between(c, d, b))
306 p = c, q = b;
307 else if( between(a, b, c) && between(c, d, a))
308 p = c, q = a;
309 else if( between(a, b, d) && between(c, d, b))
310 p = d, q = b;
311 else if( between(a, b, d) && between(c, d, a))
312 p = d, q = a;
313 else
314 code = '0';
315 return code;
316 }
317
318 //---------------------------------------------------------------------
319 // segSegInt: Finds the point of intersection p between two closed
320 // segments ab and cd. Returns p and a char with the following meaning:
321 // 'e': The segments collinearly overlap, sharing a point.
322 // 'v': An endpoint (vertex) of one segment is on the other segment,
323 // but 'e' doesn't hold.
324 // '1': The segments intersect properly (i.e., they share a point and
325 // neither 'v' nor 'e' holds).
326 // '0': The segments do not intersect (i.e., they share no points).
327 // Note that two collinear segments that share just one point, an endpoint
328 // of each, returns 'e' rather than 'v' as one might expect.
329 //---------------------------------------------------------------------
segSegInt(Point2f a,Point2f b,Point2f c,Point2f d,Point2f & p,Point2f & q)330 static char segSegInt( Point2f a, Point2f b, Point2f c, Point2f d, Point2f& p, Point2f& q )
331 {
332 double s, t; // The two parameters of the parametric eqns.
333 double num, denom; // Numerator and denoninator of equations.
334 char code = '?'; // Return char characterizing intersection.
335
336 denom = a.x * (double)( d.y - c.y ) +
337 b.x * (double)( c.y - d.y ) +
338 d.x * (double)( b.y - a.y ) +
339 c.x * (double)( a.y - b.y );
340
341 // If denom is zero, then segments are parallel: handle separately.
342 if (denom == 0.0)
343 return parallelInt(a, b, c, d, p, q);
344
345 num = a.x * (double)( d.y - c.y ) +
346 c.x * (double)( a.y - d.y ) +
347 d.x * (double)( c.y - a.y );
348 if ( (num == 0.0) || (num == denom) ) code = 'v';
349 s = num / denom;
350
351 num = -( a.x * (double)( c.y - b.y ) +
352 b.x * (double)( a.y - c.y ) +
353 c.x * (double)( b.y - a.y ) );
354 if ( (num == 0.0) || (num == denom) ) code = 'v';
355 t = num / denom;
356
357 if ( (0.0 < s) && (s < 1.0) &&
358 (0.0 < t) && (t < 1.0) )
359 code = '1';
360 else if ( (0.0 > s) || (s > 1.0) ||
361 (0.0 > t) || (t > 1.0) )
362 code = '0';
363
364 p.x = (float)(a.x + s*(b.x - a.x));
365 p.y = (float)(a.y + s*(b.y - a.y));
366
367 return code;
368 }
369
inOut(Point2f p,tInFlag inflag,int aHB,int bHA,Point2f * & result)370 static tInFlag inOut( Point2f p, tInFlag inflag, int aHB, int bHA, Point2f*& result )
371 {
372 if( p != result[-1] )
373 *result++ = p;
374 // Update inflag.
375 return aHB > 0 ? Pin : bHA > 0 ? Qin : inflag;
376 }
377
378 //---------------------------------------------------------------------
379 // Advances and prints out an inside vertex if appropriate.
380 //---------------------------------------------------------------------
advance(int a,int * aa,int n,bool inside,Point2f v,Point2f * & result)381 static int advance( int a, int *aa, int n, bool inside, Point2f v, Point2f*& result )
382 {
383 if( inside && v != result[-1] )
384 *result++ = v;
385 (*aa)++;
386 return (a+1) % n;
387 }
388
addSharedSeg(Point2f p,Point2f q,Point2f * & result)389 static void addSharedSeg( Point2f p, Point2f q, Point2f*& result )
390 {
391 if( p != result[-1] )
392 *result++ = p;
393 if( q != result[-1] )
394 *result++ = q;
395 }
396
397
intersectConvexConvex_(const Point2f * P,int n,const Point2f * Q,int m,Point2f * result,float * _area)398 static int intersectConvexConvex_( const Point2f* P, int n, const Point2f* Q, int m,
399 Point2f* result, float* _area )
400 {
401 Point2f* result0 = result;
402 // P has n vertices, Q has m vertices.
403 int a=0, b=0; // indices on P and Q (resp.)
404 Point2f Origin(0,0);
405 tInFlag inflag=Unknown; // {Pin, Qin, Unknown}: which inside
406 int aa=0, ba=0; // # advances on a & b indices (after 1st inter.)
407 bool FirstPoint=true;// Is this the first point? (used to initialize).
408 Point2f p0; // The first point.
409 *result++ = Point2f(FLT_MAX, FLT_MAX);
410
411 do
412 {
413 // Computations of key variables.
414 int a1 = (a + n - 1) % n; // a-1, b-1 (resp.)
415 int b1 = (b + m - 1) % m;
416
417 Point2f A = P[a] - P[a1], B = Q[b] - Q[b1]; // directed edges on P and Q (resp.)
418
419 int cross = areaSign( Origin, A, B ); // sign of z-component of A x B
420 int aHB = areaSign( Q[b1], Q[b], P[a] ); // a in H(b).
421 int bHA = areaSign( P[a1], P[a], Q[b] ); // b in H(A);
422
423 // If A & B intersect, update inflag.
424 Point2f p, q;
425 int code = segSegInt( P[a1], P[a], Q[b1], Q[b], p, q );
426 if( code == '1' || code == 'v' )
427 {
428 if( inflag == Unknown && FirstPoint )
429 {
430 aa = ba = 0;
431 FirstPoint = false;
432 p0 = p;
433 *result++ = p;
434 }
435 inflag = inOut( p, inflag, aHB, bHA, result );
436 }
437
438 //-----Advance rules-----
439
440 // Special case: A & B overlap and oppositely oriented.
441 if( code == 'e' && A.ddot(B) < 0 )
442 {
443 addSharedSeg( p, q, result );
444 return (int)(result - result0);
445 }
446
447 // Special case: A & B parallel and separated.
448 if( cross == 0 && aHB < 0 && bHA < 0 )
449 return (int)(result - result0);
450
451 // Special case: A & B collinear.
452 else if ( cross == 0 && aHB == 0 && bHA == 0 ) {
453 // Advance but do not output point.
454 if ( inflag == Pin )
455 b = advance( b, &ba, m, inflag == Qin, Q[b], result );
456 else
457 a = advance( a, &aa, n, inflag == Pin, P[a], result );
458 }
459
460 // Generic cases.
461 else if( cross >= 0 )
462 {
463 if( bHA > 0)
464 a = advance( a, &aa, n, inflag == Pin, P[a], result );
465 else
466 b = advance( b, &ba, m, inflag == Qin, Q[b], result );
467 }
468 else
469 {
470 if( aHB > 0)
471 b = advance( b, &ba, m, inflag == Qin, Q[b], result );
472 else
473 a = advance( a, &aa, n, inflag == Pin, P[a], result );
474 }
475 // Quit when both adv. indices have cycled, or one has cycled twice.
476 }
477 while ( ((aa < n) || (ba < m)) && (aa < 2*n) && (ba < 2*m) );
478
479 // Deal with special cases: not implemented.
480 if( inflag == Unknown )
481 {
482 // The boundaries of P and Q do not cross.
483 // ...
484 }
485
486 int i, nr = (int)(result - result0);
487 double area = 0;
488 Point2f prev = result0[nr-1];
489 for( i = 1; i < nr; i++ )
490 {
491 result0[i-1] = result0[i];
492 area += (double)prev.x*result0[i].y - (double)prev.y*result0[i].x;
493 prev = result0[i];
494 }
495
496 *_area = (float)(area*0.5);
497
498 if( result0[nr-2] == result0[0] && nr > 1 )
499 nr--;
500 return nr-1;
501 }
502
503 }
504
intersectConvexConvex(InputArray _p1,InputArray _p2,OutputArray _p12,bool handleNested)505 float cv::intersectConvexConvex( InputArray _p1, InputArray _p2, OutputArray _p12, bool handleNested )
506 {
507 Mat p1 = _p1.getMat(), p2 = _p2.getMat();
508 CV_Assert( p1.depth() == CV_32S || p1.depth() == CV_32F );
509 CV_Assert( p2.depth() == CV_32S || p2.depth() == CV_32F );
510
511 int n = p1.checkVector(2, p1.depth(), true);
512 int m = p2.checkVector(2, p2.depth(), true);
513
514 CV_Assert( n >= 0 && m >= 0 );
515
516 if( n < 2 || m < 2 )
517 {
518 _p12.release();
519 return 0.f;
520 }
521
522 AutoBuffer<Point2f> _result(n*2 + m*2 + 1);
523 Point2f *fp1 = _result, *fp2 = fp1 + n;
524 Point2f* result = fp2 + m;
525 int orientation = 0;
526
527 for( int k = 1; k <= 2; k++ )
528 {
529 Mat& p = k == 1 ? p1 : p2;
530 int len = k == 1 ? n : m;
531 Point2f* dst = k == 1 ? fp1 : fp2;
532
533 Mat temp(p.size(), CV_MAKETYPE(CV_32F, p.channels()), dst);
534 p.convertTo(temp, CV_32F);
535 CV_Assert( temp.ptr<Point2f>() == dst );
536 Point2f diff0 = dst[0] - dst[len-1];
537 for( int i = 1; i < len; i++ )
538 {
539 double s = diff0.cross(dst[i] - dst[i-1]);
540 if( s != 0 )
541 {
542 if( s < 0 )
543 {
544 orientation++;
545 flip( temp, temp, temp.rows > 1 ? 0 : 1 );
546 }
547 break;
548 }
549 }
550 }
551
552 float area = 0.f;
553 int nr = intersectConvexConvex_(fp1, n, fp2, m, result, &area);
554 if( nr == 0 )
555 {
556 if( !handleNested )
557 {
558 _p12.release();
559 return 0.f;
560 }
561
562 if( pointPolygonTest(_InputArray(fp1, n), fp2[0], false) >= 0 )
563 {
564 result = fp2;
565 nr = m;
566 }
567 else if( pointPolygonTest(_InputArray(fp2, n), fp1[0], false) >= 0 )
568 {
569 result = fp1;
570 nr = n;
571 }
572 else
573 {
574 _p12.release();
575 return 0.f;
576 }
577 area = (float)contourArea(_InputArray(result, nr), false);
578 }
579
580 if( _p12.needed() )
581 {
582 Mat temp(nr, 1, CV_32FC2, result);
583 // if both input contours were reflected,
584 // let's orient the result as the input vectors
585 if( orientation == 2 )
586 flip(temp, temp, 0);
587
588 temp.copyTo(_p12);
589 }
590 return (float)fabs(area);
591 }
592