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 // INFORMATION REGARDING THE CONTRIBUTION: 10 // 11 // Author: Ovidiu Parvu 12 // Affiliation: Brunel University 13 // Created: 11.09.2013 14 // E-mail: <ovidiu.parvu[AT]gmail.com> 15 // Web: http://people.brunel.ac.uk/~cspgoop 16 // 17 // These functions were implemented during Ovidiu Parvu's first year as a PhD student at 18 // Brunel University, London, UK. The PhD project is supervised by prof. David Gilbert (principal) 19 // and prof. Nigel Saunders (second). 20 // 21 // THE IMPLEMENTATION OF THE MODULES IS BASED ON THE FOLLOWING PAPERS: 22 // 23 // [1] V. Klee and M. C. Laskowski, "Finding the smallest triangles containing a given convex 24 // polygon", Journal of Algorithms, vol. 6, no. 3, pp. 359-375, Sep. 1985. 25 // [2] J. O'Rourke, A. Aggarwal, S. Maddila, and M. Baldwin, "An optimal algorithm for finding 26 // minimal enclosing triangles", Journal of Algorithms, vol. 7, no. 2, pp. 258-269, Jun. 1986. 27 // 28 // The overall complexity of the algorithm is theta(n) where "n" represents the number 29 // of vertices in the convex polygon. 30 // 31 // 32 // 33 // License Agreement 34 // For Open Source Computer Vision Library 35 // 36 // Copyright (C) 2000, Intel Corporation, all rights reserved. 37 // Copyright (C) 2013, OpenCV Foundation, all rights reserved. 38 // Copyright (C) 2013, Ovidiu Parvu, all rights reserved. 39 // Third party copyrights are property of their respective owners. 40 // 41 // Redistribution and use in source and binary forms, with or without modification, 42 // are permitted provided that the following conditions are met: 43 // 44 // * Redistribution's of source code must retain the above copyright notice, 45 // this list of conditions and the following disclaimer. 46 // 47 // * Redistribution's in binary form must reproduce the above copyright notice, 48 // this list of conditions and the following disclaimer in the documentation 49 // and/or other materials provided with the distribution. 50 // 51 // * The name of the copyright holders may not be used to endorse or promote products 52 // derived from this software without specific prior written permission. 53 // 54 // This software is provided by the copyright holders and contributors "as is" and 55 // any express or implied warranties, including, but not limited to, the implied 56 // warranties of merchantability and fitness for a particular purpose are disclaimed. 57 // In no event shall the Intel Corporation or contributors be liable for any direct, 58 // indirect, incidental, special, exemplary, or consequential damages 59 // (including, but not limited to, procurement of substitute goods or services; 60 // loss of use, data, or profits; or business interruption) however caused 61 // and on any theory of liability, whether in contract, strict liability, 62 // or tort (including negligence or otherwise) arising in any way out of 63 // the use of this software, even if advised of the possibility of such damage. 64 // 65 //M*/ 66 67 #include "precomp.hpp" 68 69 #include <algorithm> 70 #include <cmath> 71 #include <limits> 72 #include <vector> 73 74 75 ///////////////////////////////// Constants definitions ////////////////////////////////// 76 77 78 // Intersection of line and polygon 79 80 #define INTERSECTS_BELOW 1 81 #define INTERSECTS_ABOVE 2 82 #define INTERSECTS_CRITICAL 3 83 84 // Error messages 85 86 #define ERR_SIDE_B_GAMMA "The position of side B could not be determined, because gamma(b) could not be computed." 87 #define ERR_VERTEX_C_ON_SIDE_B "The position of the vertex C on side B could not be determined, because the considered lines do not intersect." 88 89 // Possible values for validation flag 90 91 #define VALIDATION_SIDE_A_TANGENT 0 92 #define VALIDATION_SIDE_B_TANGENT 1 93 #define VALIDATION_SIDES_FLUSH 2 94 95 // Threshold value for comparisons 96 97 #define EPSILON 1E-5 98 99 100 ////////////////////////////// Helper functions declarations ///////////////////////////// 101 102 103 namespace minEnclosingTriangle { 104 105 static void advance(unsigned int &index, unsigned int nrOfPoints); 106 107 static void advanceBToRightChain(const std::vector<cv::Point2f> &polygon, 108 unsigned int nrOfPoints, unsigned int &b, 109 unsigned int c); 110 111 static bool almostEqual(double number1, double number2); 112 113 static double angleOfLineWrtOxAxis(const cv::Point2f &a, const cv::Point2f &b); 114 115 static bool areEqualPoints(const cv::Point2f &point1, const cv::Point2f &point2); 116 117 static bool areIdenticalLines(const std::vector<double> &side1Params, 118 const std::vector<double> &side2Params, double sideCExtraParam); 119 120 static bool areIdenticalLines(double a1, double b1, double c1, double a2, double b2, double c2); 121 122 static bool areIntersectingLines(const std::vector<double> &side1Params, 123 const std::vector<double> &side2Params, 124 double sideCExtraParam, cv::Point2f &intersectionPoint1, 125 cv::Point2f &intersectionPoint2); 126 127 static bool areOnTheSameSideOfLine(const cv::Point2f &p1, const cv::Point2f &p2, 128 const cv::Point2f &a, const cv::Point2f &b); 129 130 static double areaOfTriangle(const cv::Point2f &a, const cv::Point2f &b, const cv::Point2f &c); 131 132 static void copyResultingTriangle(const std::vector<cv::Point2f> &resultingTriangle, cv::OutputArray triangle); 133 134 static void createConvexHull(cv::InputArray points, std::vector<cv::Point2f> &polygon); 135 136 static double distanceBtwPoints(const cv::Point2f &a, const cv::Point2f &b); 137 138 static double distanceFromPointToLine(const cv::Point2f &a, const cv::Point2f &linePointB, 139 const cv::Point2f &linePointC); 140 141 static bool findGammaIntersectionPoints(const std::vector<cv::Point2f> &polygon, unsigned int nrOfPoints, 142 unsigned int c, unsigned int polygonPointIndex, 143 const cv::Point2f &side1StartVertex, const cv::Point2f &side1EndVertex, 144 const cv::Point2f &side2StartVertex, const cv::Point2f &side2EndVertex, 145 cv::Point2f &intersectionPoint1, cv::Point2f &intersectionPoint2); 146 147 static void findMinEnclosingTriangle(cv::InputArray points, 148 CV_OUT cv::OutputArray triangle, CV_OUT double &area); 149 150 static void findMinEnclosingTriangle(const std::vector<cv::Point2f> &polygon, 151 std::vector<cv::Point2f> &triangle, double &area); 152 153 static void findMinimumAreaEnclosingTriangle(const std::vector<cv::Point2f> &polygon, 154 std::vector<cv::Point2f> &triangle, double &area); 155 156 static cv::Point2f findVertexCOnSideB(const std::vector<cv::Point2f> &polygon, unsigned int nrOfPoints, 157 unsigned int a, unsigned int c, 158 const cv::Point2f &sideBStartVertex, 159 const cv::Point2f &sideBEndVertex, 160 const cv::Point2f &sideCStartVertex, 161 const cv::Point2f &sideCEndVertex); 162 163 static bool gamma(unsigned int polygonPointIndex, cv::Point2f &gammaPoint, 164 const std::vector<cv::Point2f> &polygon, unsigned int nrOfPoints, 165 unsigned int a, unsigned int c); 166 167 static bool greaterOrEqual(double number1, double number2); 168 169 static double height(const cv::Point2f &polygonPoint, const std::vector<cv::Point2f> &polygon, 170 unsigned int nrOfPoints, unsigned int c); 171 172 static double height(unsigned int polygonPointIndex, const std::vector<cv::Point2f> &polygon, 173 unsigned int nrOfPoints, unsigned int c); 174 175 static void initialise(std::vector<cv::Point2f> &triangle, double &area); 176 177 static unsigned int intersects(double angleGammaAndPoint, unsigned int polygonPointIndex, 178 const std::vector<cv::Point2f> &polygon, unsigned int nrOfPoints, 179 unsigned int c); 180 181 static bool intersectsAbove(const cv::Point2f &gammaPoint, unsigned int polygonPointIndex, 182 const std::vector<cv::Point2f> &polygon, unsigned int nrOfPoints, 183 unsigned int c); 184 185 static unsigned int intersectsAboveOrBelow(unsigned int succPredIndex, unsigned int pointIndex, 186 const std::vector<cv::Point2f> &polygon, 187 unsigned int nrOfPoints, unsigned int c); 188 189 static bool intersectsBelow(const cv::Point2f &gammaPoint, unsigned int polygonPointIndex, 190 const std::vector<cv::Point2f> &polygon, unsigned int nrOfPoints, 191 unsigned int c); 192 193 static bool isAngleBetween(double angle1, double angle2, double angle3); 194 195 static bool isAngleBetweenNonReflex(double angle1, double angle2, double angle3); 196 197 static bool isFlushAngleBtwPredAndSucc(double &angleFlushEdge, double anglePred, double angleSucc); 198 199 static bool isGammaAngleBtw(double &gammaAngle, double angle1, double angle2); 200 201 static bool isGammaAngleEqualTo(double &gammaAngle, double angle); 202 203 static bool isLocalMinimalTriangle(cv::Point2f &vertexA, cv::Point2f &vertexB, 204 cv::Point2f &vertexC, const std::vector<cv::Point2f> &polygon, 205 unsigned int nrOfPoints, unsigned int a, unsigned int b, 206 unsigned int validationFlag, const cv::Point2f &sideAStartVertex, 207 const cv::Point2f &sideAEndVertex, const cv::Point2f &sideBStartVertex, 208 const cv::Point2f &sideBEndVertex, const cv::Point2f &sideCStartVertex, 209 const cv::Point2f &sideCEndVertex); 210 211 static bool isNotBTangency(const std::vector<cv::Point2f> &polygon, 212 unsigned int nrOfPoints, unsigned int a, unsigned int b, 213 unsigned int c); 214 215 static bool isOppositeAngleBetweenNonReflex(double angle1, double angle2, double angle3); 216 217 static bool isPointOnLineSegment(const cv::Point2f &point, const cv::Point2f &lineSegmentStart, 218 const cv::Point2f &lineSegmentEnd); 219 220 static bool isValidMinimalTriangle(const cv::Point2f &vertexA, const cv::Point2f &vertexB, 221 const cv::Point2f &vertexC, const std::vector<cv::Point2f> &polygon, 222 unsigned int nrOfPoints, unsigned int a, unsigned int b, 223 unsigned int validationFlag, const cv::Point2f &sideAStartVertex, 224 const cv::Point2f &sideAEndVertex, const cv::Point2f &sideBStartVertex, 225 const cv::Point2f &sideBEndVertex, const cv::Point2f &sideCStartVertex, 226 const cv::Point2f &sideCEndVertex); 227 228 static bool lessOrEqual(double number1, double number2); 229 230 static void lineEquationDeterminedByPoints(const cv::Point2f &p, const cv::Point2f &q, 231 double &a, double &b, double &c); 232 233 static std::vector<double> lineEquationParameters(const cv::Point2f& p, const cv::Point2f &q); 234 235 static bool lineIntersection(const cv::Point2f &a1, const cv::Point2f &b1, const cv::Point2f &a2, 236 const cv::Point2f &b2, cv::Point2f &intersection); 237 238 static bool lineIntersection(double a1, double b1, double c1, double a2, double b2, double c2, 239 cv::Point2f &intersection); 240 241 static double maximum(double number1, double number2, double number3); 242 243 static cv::Point2f middlePoint(const cv::Point2f &a, const cv::Point2f &b); 244 245 static bool middlePointOfSideB(cv::Point2f &middlePointOfSideB, const cv::Point2f &sideAStartVertex, 246 const cv::Point2f &sideAEndVertex, const cv::Point2f &sideBStartVertex, 247 const cv::Point2f &sideBEndVertex, const cv::Point2f &sideCStartVertex, 248 const cv::Point2f &sideCEndVertex); 249 250 static void moveAIfLowAndBIfHigh(const std::vector<cv::Point2f> &polygon, 251 unsigned int nrOfPoints, unsigned int &a, unsigned int &b, 252 unsigned int c); 253 254 static double oppositeAngle(double angle); 255 256 static unsigned int predecessor(unsigned int index, unsigned int nrOfPoints); 257 258 static void returnMinimumAreaEnclosingTriangle(const std::vector<cv::Point2f> &polygon, 259 std::vector<cv::Point2f> &triangle, double &area); 260 261 static void searchForBTangency(const std::vector<cv::Point2f> &polygon, 262 unsigned int nrOfPoints, unsigned int a, unsigned int &b, 263 unsigned int c); 264 265 static int sign(double number); 266 267 static unsigned int successor(unsigned int index, unsigned int nrOfPoints); 268 269 static void updateMinimumAreaEnclosingTriangle(std::vector<cv::Point2f> &triangle, double &area, 270 const cv::Point2f &vertexA, const cv::Point2f &vertexB, 271 const cv::Point2f &vertexC); 272 273 static void updateSideB(const std::vector<cv::Point2f> &polygon, 274 unsigned int nrOfPoints, unsigned int a, unsigned int b, 275 unsigned int c, unsigned int &validationFlag, 276 cv::Point2f &sideBStartVertex, cv::Point2f &sideBEndVertex); 277 278 static void updateSidesBA(const std::vector<cv::Point2f> &polygon, 279 unsigned int nrOfPoints, unsigned int a, unsigned int b, 280 unsigned int c, unsigned int &validationFlag, 281 cv::Point2f &sideAStartVertex, cv::Point2f &sideAEndVertex, 282 cv::Point2f &sideBStartVertex, cv::Point2f &sideBEndVertex, 283 const cv::Point2f &sideCStartVertex, const cv::Point2f &sideCEndVertex); 284 285 static void updateSidesCA(const std::vector<cv::Point2f> &polygon, 286 unsigned int nrOfPoints, unsigned int a, unsigned int c, 287 cv::Point2f &sideAStartVertex, cv::Point2f &sideAEndVertex, 288 cv::Point2f &sideCStartVertex, cv::Point2f &sideCEndVertex); 289 290 } 291 292 293 ///////////////////////////////////// Main functions ///////////////////////////////////// 294 295 296 //! Find the minimum enclosing triangle for the given set of points and return its area 297 /*! 298 * @param points Set of points 299 * @param triangle Minimum area triangle enclosing the given set of points 300 */ 301 double cv::minEnclosingTriangle(cv::InputArray points, CV_OUT cv::OutputArray triangle) { 302 double area; 303 304 minEnclosingTriangle::findMinEnclosingTriangle(points, triangle, area); 305 306 return area; 307 } 308 309 310 /////////////////////////////// Helper functions definition ////////////////////////////// 311 312 313 namespace minEnclosingTriangle { 314 315 //! Find the minimum enclosing triangle and its area 316 /*! 317 * @param points Set of points 318 * @param triangle Minimum area triangle enclosing the given set of points 319 * @param area Area of the minimum area enclosing triangle 320 */ 321 static void findMinEnclosingTriangle(cv::InputArray points, 322 CV_OUT cv::OutputArray triangle, CV_OUT double &area) { 323 std::vector<cv::Point2f> resultingTriangle, polygon; 324 325 createConvexHull(points, polygon); 326 findMinEnclosingTriangle(polygon, resultingTriangle, area); 327 copyResultingTriangle(resultingTriangle, triangle); 328 } 329 330 //! Create the convex hull of the given set of points 331 /*! 332 * @param points The provided set of points 333 * @param polygon The polygon representing the convex hull of the points 334 */ 335 static void createConvexHull(cv::InputArray points, std::vector<cv::Point2f> &polygon) { 336 cv::Mat pointsMat = points.getMat(); 337 std::vector<cv::Point2f> pointsVector; 338 339 CV_Assert((pointsMat.checkVector(2) > 0) && 340 ((pointsMat.depth() == CV_32F) || (pointsMat.depth() == CV_32S))); 341 342 pointsMat.convertTo(pointsVector, CV_32F); 343 344 convexHull(pointsVector, polygon, true, true); 345 } 346 347 //! Find the minimum enclosing triangle and its area 348 /*! 349 * The overall complexity of the algorithm is theta(n) where "n" represents the number 350 * of vertices in the convex polygon 351 * 352 * @param polygon The polygon representing the convex hull of the points 353 * @param triangle Minimum area triangle enclosing the given polygon 354 * @param area Area of the minimum area enclosing triangle 355 */ 356 static void findMinEnclosingTriangle(const std::vector<cv::Point2f> &polygon, 357 std::vector<cv::Point2f> &triangle, double &area) { 358 initialise(triangle, area); 359 360 if (polygon.size() > 3) { 361 findMinimumAreaEnclosingTriangle(polygon, triangle, area); 362 } else { 363 returnMinimumAreaEnclosingTriangle(polygon, triangle, area); 364 } 365 } 366 367 //! Copy resultingTriangle to the OutputArray triangle 368 /*! 369 * @param resultingTriangle Minimum area triangle enclosing the given polygon found by the algorithm 370 * @param triangle Minimum area triangle enclosing the given polygon returned to the user 371 */ 372 static void copyResultingTriangle(const std::vector<cv::Point2f> &resultingTriangle, 373 cv::OutputArray triangle) { 374 cv::Mat(resultingTriangle).copyTo(triangle); 375 } 376 377 //! Initialisation function 378 /*! 379 * @param triangle Minimum area triangle enclosing the given polygon 380 * @param area Area of the minimum area enclosing triangle 381 */ 382 static void initialise(std::vector<cv::Point2f> &triangle, double &area) { 383 area = std::numeric_limits<double>::max(); 384 385 // Clear all points previously stored in the vector 386 triangle.clear(); 387 } 388 389 //! Find the minimum area enclosing triangle for the given polygon 390 /*! 391 * @param polygon The polygon representing the convex hull of the points 392 * @param triangle Minimum area triangle enclosing the given polygon 393 * @param area Area of the minimum area enclosing triangle 394 */ 395 static void findMinimumAreaEnclosingTriangle(const std::vector<cv::Point2f> &polygon, 396 std::vector<cv::Point2f> &triangle, double &area) { 397 // Algorithm specific variables 398 399 unsigned int validationFlag; 400 401 cv::Point2f vertexA, vertexB, vertexC; 402 403 cv::Point2f sideAStartVertex, sideAEndVertex; 404 cv::Point2f sideBStartVertex, sideBEndVertex; 405 cv::Point2f sideCStartVertex, sideCEndVertex; 406 407 unsigned int a, b, c; 408 unsigned int nrOfPoints; 409 410 // Variables initialisation 411 412 nrOfPoints = static_cast<unsigned int>(polygon.size()); 413 414 a = 1; 415 b = 2; 416 c = 0; 417 418 // Main algorithm steps 419 420 for (c = 0; c < nrOfPoints; c++) { 421 advanceBToRightChain(polygon, nrOfPoints, b, c); 422 moveAIfLowAndBIfHigh(polygon, nrOfPoints, a, b, c); 423 searchForBTangency(polygon, nrOfPoints, a ,b, c); 424 425 updateSidesCA(polygon, nrOfPoints, a, c, sideAStartVertex, sideAEndVertex, 426 sideCStartVertex, sideCEndVertex); 427 428 if (isNotBTangency(polygon, nrOfPoints, a, b, c)) { 429 updateSidesBA(polygon, nrOfPoints, a, b, c, validationFlag, sideAStartVertex, 430 sideAEndVertex, sideBStartVertex, sideBEndVertex, 431 sideCStartVertex, sideCEndVertex); 432 } else { 433 updateSideB(polygon, nrOfPoints, a, b, c, validationFlag, 434 sideBStartVertex, sideBEndVertex); 435 } 436 437 if (isLocalMinimalTriangle(vertexA, vertexB, vertexC, polygon, nrOfPoints, a, b, 438 validationFlag, sideAStartVertex, sideAEndVertex, 439 sideBStartVertex, sideBEndVertex, sideCStartVertex, 440 sideCEndVertex)) { 441 updateMinimumAreaEnclosingTriangle(triangle, area, vertexA, vertexB, vertexC); 442 } 443 } 444 } 445 446 //! Return the minimum area enclosing (pseudo-)triangle in case the convex polygon has at most three points 447 /*! 448 * @param polygon The polygon representing the convex hull of the points 449 * @param triangle Minimum area triangle enclosing the given polygon 450 * @param area Area of the minimum area enclosing triangle 451 */ 452 static void returnMinimumAreaEnclosingTriangle(const std::vector<cv::Point2f> &polygon, 453 std::vector<cv::Point2f> &triangle, double &area) { 454 unsigned int nrOfPoints = static_cast<unsigned int>(polygon.size()); 455 456 for (int i = 0; i < 3; i++) { 457 triangle.push_back(polygon[i % nrOfPoints]); 458 } 459 460 area = areaOfTriangle(triangle[0], triangle[1], triangle[2]); 461 } 462 463 //! Advance b to the right chain 464 /*! 465 * See paper [2] for more details 466 * 467 * @param polygon The polygon representing the convex hull of the points 468 * @param nrOfPoints Number of points defining the convex polygon 469 * @param b Index b 470 * @param c Index c 471 */ 472 static void advanceBToRightChain(const std::vector<cv::Point2f> &polygon, 473 unsigned int nrOfPoints, unsigned int &b, 474 unsigned int c) { 475 while (greaterOrEqual(height(successor(b, nrOfPoints), polygon, nrOfPoints, c), 476 height(b, polygon, nrOfPoints, c))) { 477 advance(b, nrOfPoints); 478 } 479 } 480 481 //! Move "a" if it is low and "b" if it is high 482 /*! 483 * See paper [2] for more details 484 * 485 * @param polygon The polygon representing the convex hull of the points 486 * @param nrOfPoints Number of points defining the convex polygon 487 * @param a Index a 488 * @param b Index b 489 * @param c Index c 490 */ 491 static void moveAIfLowAndBIfHigh(const std::vector<cv::Point2f> &polygon, 492 unsigned int nrOfPoints, unsigned int &a, unsigned int &b, 493 unsigned int c) { 494 cv::Point2f gammaOfA; 495 496 while(height(b, polygon, nrOfPoints, c) > height(a, polygon, nrOfPoints, c)) { 497 if ((gamma(a, gammaOfA, polygon, nrOfPoints, a, c)) && (intersectsBelow(gammaOfA, b, polygon, nrOfPoints, c))) { 498 advance(b, nrOfPoints); 499 } else { 500 advance(a, nrOfPoints); 501 } 502 } 503 } 504 505 //! Search for the tangency of side B 506 /*! 507 * See paper [2] for more details 508 * 509 * @param polygon The polygon representing the convex hull of the points 510 * @param nrOfPoints Number of points defining the convex polygon 511 * @param a Index a 512 * @param b Index b 513 * @param c Index c 514 */ 515 static void searchForBTangency(const std::vector<cv::Point2f> &polygon, 516 unsigned int nrOfPoints, unsigned int a, unsigned int &b, 517 unsigned int c) { 518 cv::Point2f gammaOfB; 519 520 while (((gamma(b, gammaOfB, polygon, nrOfPoints, a, c)) && 521 (intersectsBelow(gammaOfB, b, polygon, nrOfPoints, c))) && 522 (greaterOrEqual(height(b, polygon, nrOfPoints, c), 523 height(predecessor(a, nrOfPoints), polygon, nrOfPoints, c))) 524 ) { 525 advance(b, nrOfPoints); 526 } 527 } 528 529 //! Check if tangency for side B was not obtained 530 /*! 531 * See paper [2] for more details 532 * 533 * @param polygon The polygon representing the convex hull of the points 534 * @param nrOfPoints Number of points defining the convex polygon 535 * @param a Index a 536 * @param b Index b 537 * @param c Index c 538 */ 539 static bool isNotBTangency(const std::vector<cv::Point2f> &polygon, 540 unsigned int nrOfPoints, unsigned int a, unsigned int b, 541 unsigned int c) { 542 cv::Point2f gammaOfB; 543 544 if (((gamma(b, gammaOfB, polygon, nrOfPoints, a, c)) && 545 (intersectsAbove(gammaOfB, b, polygon, nrOfPoints, c))) || 546 (height(b, polygon, nrOfPoints, c) < height(predecessor(a, nrOfPoints), polygon, nrOfPoints, c))) { 547 return true; 548 } 549 550 return false; 551 } 552 553 //! Update sides A and C 554 /*! 555 * Side C will have as start and end vertices the polygon points "c" and "c-1" 556 * Side A will have as start and end vertices the polygon points "a" and "a-1" 557 * 558 * @param polygon The polygon representing the convex hull of the points 559 * @param nrOfPoints Number of points defining the convex polygon 560 * @param a Index a 561 * @param c Index c 562 * @param sideAStartVertex Start vertex for defining side A 563 * @param sideAEndVertex End vertex for defining side A 564 * @param sideCStartVertex Start vertex for defining side C 565 * @param sideCEndVertex End vertex for defining side C 566 */ 567 static void updateSidesCA(const std::vector<cv::Point2f> &polygon, 568 unsigned int nrOfPoints, unsigned int a, unsigned int c, 569 cv::Point2f &sideAStartVertex, cv::Point2f &sideAEndVertex, 570 cv::Point2f &sideCStartVertex, cv::Point2f &sideCEndVertex) { 571 sideCStartVertex = polygon[predecessor(c, nrOfPoints)]; 572 sideCEndVertex = polygon[c]; 573 574 sideAStartVertex = polygon[predecessor(a, nrOfPoints)]; 575 sideAEndVertex = polygon[a]; 576 } 577 578 //! Update sides B and possibly A if tangency for side B was not obtained 579 /*! 580 * See paper [2] for more details 581 * 582 * @param polygon The polygon representing the convex hull of the points 583 * @param nrOfPoints Number of points defining the convex polygon 584 * @param a Index a 585 * @param b Index b 586 * @param c Index c 587 * @param validationFlag Flag used for validation 588 * @param sideAStartVertex Start vertex for defining side A 589 * @param sideAEndVertex End vertex for defining side A 590 * @param sideBStartVertex Start vertex for defining side B 591 * @param sideBEndVertex End vertex for defining side B 592 * @param sideCStartVertex Start vertex for defining side C 593 * @param sideCEndVertex End vertex for defining side C 594 */ 595 static void updateSidesBA(const std::vector<cv::Point2f> &polygon, 596 unsigned int nrOfPoints, unsigned int a, unsigned int b, 597 unsigned int c, unsigned int &validationFlag, 598 cv::Point2f &sideAStartVertex, cv::Point2f &sideAEndVertex, 599 cv::Point2f &sideBStartVertex, cv::Point2f &sideBEndVertex, 600 const cv::Point2f &sideCStartVertex, const cv::Point2f &sideCEndVertex) { 601 // Side B is flush with edge [b, b-1] 602 sideBStartVertex = polygon[predecessor(b, nrOfPoints)]; 603 sideBEndVertex = polygon[b]; 604 605 // Find middle point of side B 606 cv::Point2f sideBMiddlePoint; 607 608 if ((middlePointOfSideB(sideBMiddlePoint, sideAStartVertex, sideAEndVertex, sideBStartVertex, 609 sideBEndVertex, sideCStartVertex, sideCEndVertex)) && 610 (height(sideBMiddlePoint, polygon, nrOfPoints, c) < 611 height(predecessor(a, nrOfPoints), polygon, nrOfPoints, c))) { 612 sideAStartVertex = polygon[predecessor(a, nrOfPoints)]; 613 sideAEndVertex = findVertexCOnSideB(polygon, nrOfPoints, a, c, 614 sideBStartVertex, sideBEndVertex, 615 sideCStartVertex, sideCEndVertex); 616 617 validationFlag = VALIDATION_SIDE_A_TANGENT; 618 } else { 619 validationFlag = VALIDATION_SIDES_FLUSH; 620 } 621 } 622 623 //! Set side B if tangency for side B was obtained 624 /*! 625 * See paper [2] for more details 626 * 627 * @param polygon The polygon representing the convex hull of the points 628 * @param nrOfPoints Number of points defining the convex polygon 629 * @param a Index a 630 * @param b Index b 631 * @param c Index c 632 * @param validationFlag Flag used for validation 633 * @param sideBStartVertex Start vertex for defining side B 634 * @param sideBEndVertex End vertex for defining side B 635 */ 636 static void updateSideB(const std::vector<cv::Point2f> &polygon, 637 unsigned int nrOfPoints, unsigned int a, unsigned int b, 638 unsigned int c, unsigned int &validationFlag, 639 cv::Point2f &sideBStartVertex, cv::Point2f &sideBEndVertex) { 640 if (!gamma(b, sideBStartVertex, polygon, nrOfPoints, a, c)) { 641 CV_Error(cv::Error::StsInternal, ERR_SIDE_B_GAMMA); 642 } 643 644 sideBEndVertex = polygon[b]; 645 646 validationFlag = VALIDATION_SIDE_B_TANGENT; 647 } 648 649 //! Update the triangle vertices after all sides were set and check if a local minimal triangle was found or not 650 /*! 651 * See paper [2] for more details 652 * 653 * @param vertexA Vertex A of the enclosing triangle 654 * @param vertexB Vertex B of the enclosing triangle 655 * @param vertexC Vertex C of the enclosing triangle 656 * @param polygon The polygon representing the convex hull of the points 657 * @param nrOfPoints Number of points defining the convex polygon 658 * @param a Index a 659 * @param b Index b 660 * @param validationFlag Flag used for validation 661 * @param sideAStartVertex Start vertex for defining side A 662 * @param sideAEndVertex End vertex for defining side A 663 * @param sideBStartVertex Start vertex for defining side B 664 * @param sideBEndVertex End vertex for defining side B 665 * @param sideCStartVertex Start vertex for defining side C 666 * @param sideCEndVertex End vertex for defining side C 667 */ 668 static bool isLocalMinimalTriangle(cv::Point2f &vertexA, cv::Point2f &vertexB, 669 cv::Point2f &vertexC, const std::vector<cv::Point2f> &polygon, 670 unsigned int nrOfPoints, unsigned int a, unsigned int b, 671 unsigned int validationFlag, const cv::Point2f &sideAStartVertex, 672 const cv::Point2f &sideAEndVertex, const cv::Point2f &sideBStartVertex, 673 const cv::Point2f &sideBEndVertex, const cv::Point2f &sideCStartVertex, 674 const cv::Point2f &sideCEndVertex) { 675 if ((!lineIntersection(sideAStartVertex, sideAEndVertex, 676 sideBStartVertex, sideBEndVertex, vertexC)) || 677 (!lineIntersection(sideAStartVertex, sideAEndVertex, 678 sideCStartVertex, sideCEndVertex, vertexB)) || 679 (!lineIntersection(sideBStartVertex, sideBEndVertex, 680 sideCStartVertex, sideCEndVertex, vertexA))) { 681 return false; 682 } 683 684 return isValidMinimalTriangle(vertexA, vertexB, vertexC, polygon, nrOfPoints, a, b, 685 validationFlag, sideAStartVertex, sideAEndVertex, 686 sideBStartVertex, sideBEndVertex, sideCStartVertex, 687 sideCEndVertex); 688 } 689 690 //! Check if the found minimal triangle is valid 691 /*! 692 * This means that all midpoints of the triangle should touch the polygon 693 * 694 * See paper [2] for more details 695 * 696 * @param vertexA Vertex A of the enclosing triangle 697 * @param vertexB Vertex B of the enclosing triangle 698 * @param vertexC Vertex C of the enclosing triangle 699 * @param polygon The polygon representing the convex hull of the points 700 * @param nrOfPoints Number of points defining the convex polygon 701 * @param a Index a 702 * @param b Index b 703 * @param validationFlag Flag used for validation 704 * @param sideAStartVertex Start vertex for defining side A 705 * @param sideAEndVertex End vertex for defining side A 706 * @param sideBStartVertex Start vertex for defining side B 707 * @param sideBEndVertex End vertex for defining side B 708 * @param sideCStartVertex Start vertex for defining side C 709 * @param sideCEndVertex End vertex for defining side C 710 */ 711 static bool isValidMinimalTriangle(const cv::Point2f &vertexA, const cv::Point2f &vertexB, 712 const cv::Point2f &vertexC, const std::vector<cv::Point2f> &polygon, 713 unsigned int nrOfPoints, unsigned int a, unsigned int b, 714 unsigned int validationFlag, const cv::Point2f &sideAStartVertex, 715 const cv::Point2f &sideAEndVertex, const cv::Point2f &sideBStartVertex, 716 const cv::Point2f &sideBEndVertex, const cv::Point2f &sideCStartVertex, 717 const cv::Point2f &sideCEndVertex) { 718 cv::Point2f midpointSideA = middlePoint(vertexB, vertexC); 719 cv::Point2f midpointSideB = middlePoint(vertexA, vertexC); 720 cv::Point2f midpointSideC = middlePoint(vertexA, vertexB); 721 722 bool sideAValid = (validationFlag == VALIDATION_SIDE_A_TANGENT) 723 ? (areEqualPoints(midpointSideA, polygon[predecessor(a, nrOfPoints)])) 724 : (isPointOnLineSegment(midpointSideA, sideAStartVertex, sideAEndVertex)); 725 726 bool sideBValid = (validationFlag == VALIDATION_SIDE_B_TANGENT) 727 ? (areEqualPoints(midpointSideB, polygon[b])) 728 : (isPointOnLineSegment(midpointSideB, sideBStartVertex, sideBEndVertex)); 729 730 bool sideCValid = isPointOnLineSegment(midpointSideC, sideCStartVertex, sideCEndVertex); 731 732 return (sideAValid && sideBValid && sideCValid); 733 } 734 735 //! Update the current minimum area enclosing triangle if the newly obtained one has a smaller area 736 /*! 737 * @param triangle Minimum area triangle enclosing the given polygon 738 * @param area Area of the minimum area triangle enclosing the given polygon 739 * @param vertexA Vertex A of the enclosing triangle 740 * @param vertexB Vertex B of the enclosing triangle 741 * @param vertexC Vertex C of the enclosing triangle 742 */ 743 static void updateMinimumAreaEnclosingTriangle(std::vector<cv::Point2f> &triangle, double &area, 744 const cv::Point2f &vertexA, const cv::Point2f &vertexB, 745 const cv::Point2f &vertexC) { 746 double triangleArea = areaOfTriangle(vertexA, vertexB, vertexC); 747 748 if (triangleArea < area) { 749 triangle.clear(); 750 751 triangle.push_back(vertexA); 752 triangle.push_back(vertexB); 753 triangle.push_back(vertexC); 754 755 area = triangleArea; 756 } 757 } 758 759 //! Return the middle point of side B 760 /*! 761 * @param middlePointOfSideB Middle point of side B 762 * @param sideAStartVertex Start vertex for defining side A 763 * @param sideAEndVertex End vertex for defining side A 764 * @param sideBStartVertex Start vertex for defining side B 765 * @param sideBEndVertex End vertex for defining side B 766 * @param sideCStartVertex Start vertex for defining side C 767 * @param sideCEndVertex End vertex for defining side C 768 */ 769 static bool middlePointOfSideB(cv::Point2f &middlePointOfSideB, const cv::Point2f &sideAStartVertex, 770 const cv::Point2f &sideAEndVertex, const cv::Point2f &sideBStartVertex, 771 const cv::Point2f &sideBEndVertex, const cv::Point2f &sideCStartVertex, 772 const cv::Point2f &sideCEndVertex) { 773 cv::Point2f vertexA, vertexC; 774 775 if ((!lineIntersection(sideBStartVertex, sideBEndVertex, sideCStartVertex, sideCEndVertex, vertexA)) || 776 (!lineIntersection(sideBStartVertex, sideBEndVertex, sideAStartVertex, sideAEndVertex, vertexC))) { 777 return false; 778 } 779 780 middlePointOfSideB = middlePoint(vertexA, vertexC); 781 782 return true; 783 } 784 785 //! Check if the line intersects below 786 /*! 787 * Check if the line determined by gammaPoint and polygon[polygonPointIndex] intersects 788 * the polygon below the point polygon[polygonPointIndex] 789 * 790 * @param gammaPoint Gamma(p) 791 * @param polygonPointIndex Index of the polygon point which is considered when determining the line 792 * @param polygon The polygon representing the convex hull of the points 793 * @param nrOfPoints Number of points defining the convex polygon 794 * @param c Index c 795 */ 796 static bool intersectsBelow(const cv::Point2f &gammaPoint, unsigned int polygonPointIndex, 797 const std::vector<cv::Point2f> &polygon, unsigned int nrOfPoints, 798 unsigned int c) { 799 double angleOfGammaAndPoint = angleOfLineWrtOxAxis(polygon[polygonPointIndex], gammaPoint); 800 801 return (intersects(angleOfGammaAndPoint, polygonPointIndex, polygon, nrOfPoints, c) == INTERSECTS_BELOW); 802 } 803 804 //! Check if the line intersects above 805 /*! 806 * Check if the line determined by gammaPoint and polygon[polygonPointIndex] intersects 807 * the polygon above the point polygon[polygonPointIndex] 808 * 809 * @param gammaPoint Gamma(p) 810 * @param polygonPointIndex Index of the polygon point which is considered when determining the line 811 * @param polygon The polygon representing the convex hull of the points 812 * @param nrOfPoints Number of points defining the convex polygon 813 * @param c Index c 814 */ 815 static bool intersectsAbove(const cv::Point2f &gammaPoint, unsigned int polygonPointIndex, 816 const std::vector<cv::Point2f> &polygon, unsigned int nrOfPoints, 817 unsigned int c) { 818 double angleOfGammaAndPoint = angleOfLineWrtOxAxis(gammaPoint, polygon[polygonPointIndex]); 819 820 return (intersects(angleOfGammaAndPoint, polygonPointIndex, polygon, nrOfPoints, c) == INTERSECTS_ABOVE); 821 } 822 823 //! Check if/where the line determined by gammaPoint and polygon[polygonPointIndex] intersects the polygon 824 /*! 825 * @param angleGammaAndPoint Angle determined by gammaPoint and polygon[polygonPointIndex] wrt Ox axis 826 * @param polygonPointIndex Index of the polygon point which is considered when determining the line 827 * @param polygon The polygon representing the convex hull of the points 828 * @param nrOfPoints Number of points defining the convex polygon 829 * @param c Index c 830 */ 831 static unsigned int intersects(double angleGammaAndPoint, unsigned int polygonPointIndex, 832 const std::vector<cv::Point2f> &polygon, unsigned int nrOfPoints, 833 unsigned int c) { 834 double anglePointPredecessor = angleOfLineWrtOxAxis(polygon[predecessor(polygonPointIndex, nrOfPoints)], 835 polygon[polygonPointIndex]); 836 double anglePointSuccessor = angleOfLineWrtOxAxis(polygon[successor(polygonPointIndex, nrOfPoints)], 837 polygon[polygonPointIndex]); 838 double angleFlushEdge = angleOfLineWrtOxAxis(polygon[predecessor(c, nrOfPoints)], 839 polygon[c]); 840 841 if (isFlushAngleBtwPredAndSucc(angleFlushEdge, anglePointPredecessor, anglePointSuccessor)) { 842 if ((isGammaAngleBtw(angleGammaAndPoint, anglePointPredecessor, angleFlushEdge)) || 843 (almostEqual(angleGammaAndPoint, anglePointPredecessor))) { 844 return intersectsAboveOrBelow(predecessor(polygonPointIndex, nrOfPoints), 845 polygonPointIndex, polygon, nrOfPoints, c); 846 } else if ((isGammaAngleBtw(angleGammaAndPoint, anglePointSuccessor, angleFlushEdge)) || 847 (almostEqual(angleGammaAndPoint, anglePointSuccessor))) { 848 return intersectsAboveOrBelow(successor(polygonPointIndex, nrOfPoints), 849 polygonPointIndex, polygon, nrOfPoints, c); 850 } 851 } else { 852 if ( 853 (isGammaAngleBtw(angleGammaAndPoint, anglePointPredecessor, anglePointSuccessor)) || 854 ( 855 (isGammaAngleEqualTo(angleGammaAndPoint, anglePointPredecessor)) && 856 (!isGammaAngleEqualTo(angleGammaAndPoint, angleFlushEdge)) 857 ) || 858 ( 859 (isGammaAngleEqualTo(angleGammaAndPoint, anglePointSuccessor)) && 860 (!isGammaAngleEqualTo(angleGammaAndPoint, angleFlushEdge)) 861 ) 862 ) { 863 return INTERSECTS_BELOW; 864 } 865 } 866 867 return INTERSECTS_CRITICAL; 868 } 869 870 //! If (gamma(x) x) intersects P between successorOrPredecessorIndex and pointIntex is it above/below? 871 /*! 872 * @param succPredIndex Index of the successor or predecessor 873 * @param pointIndex Index of the point x in the polygon 874 * @param polygon The polygon representing the convex hull of the points 875 * @param nrOfPoints Number of points defining the convex polygon 876 * @param c Index c 877 */ 878 static unsigned int intersectsAboveOrBelow(unsigned int succPredIndex, unsigned int pointIndex, 879 const std::vector<cv::Point2f> &polygon, 880 unsigned int nrOfPoints, unsigned int c) { 881 if (height(succPredIndex, polygon, nrOfPoints, c) > height(pointIndex, polygon, nrOfPoints, c)) { 882 return INTERSECTS_ABOVE; 883 } else { 884 return INTERSECTS_BELOW; 885 } 886 } 887 888 //! Find gamma for a given point "p" specified by its index 889 /*! 890 * The function returns true if gamma exists i.e. if lines (a a-1) and (x y) intersect 891 * and false otherwise. In case the two lines intersect in point intersectionPoint, gamma is computed. 892 * 893 * Considering that line (x y) is a line parallel to (c c-1) and that the distance between the lines is equal 894 * to 2 * height(p), we can have two possible (x y) lines. 895 * 896 * Therefore, we will compute two intersection points between the lines (x y) and (a a-1) and take the 897 * point which is on the same side of line (c c-1) as the polygon. 898 * 899 * See paper [2] and formula for distance from point to a line for more details 900 * 901 * @param polygonPointIndex Index of the polygon point 902 * @param gammaPoint Point gamma(polygon[polygonPointIndex]) 903 * @param polygon The polygon representing the convex hull of the points 904 * @param nrOfPoints Number of points defining the convex polygon 905 * @param a Index a 906 * @param c Index c 907 */ 908 static bool gamma(unsigned int polygonPointIndex, cv::Point2f &gammaPoint, 909 const std::vector<cv::Point2f> &polygon, unsigned int nrOfPoints, 910 unsigned int a, unsigned int c) { 911 cv::Point2f intersectionPoint1, intersectionPoint2; 912 913 // Get intersection points if they exist 914 if (!findGammaIntersectionPoints(polygon, nrOfPoints, c, polygonPointIndex, 915 polygon[a], polygon[predecessor(a, nrOfPoints)], 916 polygon[c], polygon[predecessor(c, nrOfPoints)], 917 intersectionPoint1, intersectionPoint2)) { 918 return false; 919 } 920 921 // Select the point which is on the same side of line C as the polygon 922 if (areOnTheSameSideOfLine(intersectionPoint1, polygon[successor(c, nrOfPoints)], 923 polygon[c], polygon[predecessor(c, nrOfPoints)])) { 924 gammaPoint = intersectionPoint1; 925 } else { 926 gammaPoint = intersectionPoint2; 927 } 928 929 return true; 930 } 931 932 //! Find vertex C which lies on side B at a distance = 2 * height(a-1) from side C 933 /*! 934 * Considering that line (x y) is a line parallel to (c c-1) and that the distance between the lines is equal 935 * to 2 * height(a-1), we can have two possible (x y) lines. 936 * 937 * Therefore, we will compute two intersection points between the lines (x y) and (b b-1) and take the 938 * point which is on the same side of line (c c-1) as the polygon. 939 * 940 * See paper [2] and formula for distance from point to a line for more details 941 * 942 * @param polygon The polygon representing the convex hull of the points 943 * @param nrOfPoints Number of points defining the convex polygon 944 * @param a Index a 945 * @param c Index c 946 * @param sideBStartVertex Start vertex for defining side B 947 * @param sideBEndVertex End vertex for defining side B 948 * @param sideCStartVertex Start vertex for defining side C 949 * @param sideCEndVertex End vertex for defining side C 950 */ 951 static cv::Point2f findVertexCOnSideB(const std::vector<cv::Point2f> &polygon, unsigned int nrOfPoints, 952 unsigned int a, unsigned int c, 953 const cv::Point2f &sideBStartVertex, 954 const cv::Point2f &sideBEndVertex, 955 const cv::Point2f &sideCStartVertex, 956 const cv::Point2f &sideCEndVertex) { 957 cv::Point2f intersectionPoint1, intersectionPoint2; 958 959 // Get intersection points if they exist 960 if (!findGammaIntersectionPoints(polygon, nrOfPoints, c, predecessor(a, nrOfPoints), 961 sideBStartVertex, sideBEndVertex, 962 sideCStartVertex, sideCEndVertex, 963 intersectionPoint1, intersectionPoint2)) { 964 CV_Error(cv::Error::StsInternal, ERR_VERTEX_C_ON_SIDE_B); 965 } 966 967 // Select the point which is on the same side of line C as the polygon 968 if (areOnTheSameSideOfLine(intersectionPoint1, polygon[successor(c, nrOfPoints)], 969 polygon[c], polygon[predecessor(c, nrOfPoints)])) { 970 return intersectionPoint1; 971 } else { 972 return intersectionPoint2; 973 } 974 } 975 976 //! Find the intersection points to compute gamma(point) 977 /*! 978 * @param polygon The polygon representing the convex hull of the points 979 * @param nrOfPoints Number of points defining the convex polygon 980 * @param c Index c 981 * @param polygonPointIndex Index of the polygon point for which the distance is known 982 * @param side1StartVertex Start vertex for side 1 983 * @param side1EndVertex End vertex for side 1 984 * @param side2StartVertex Start vertex for side 2 985 * @param side2EndVertex End vertex for side 2 986 * @param intersectionPoint1 First intersection point between one pair of lines 987 * @param intersectionPoint2 Second intersection point between other pair of lines 988 */ 989 static bool findGammaIntersectionPoints(const std::vector<cv::Point2f> &polygon, unsigned int nrOfPoints, 990 unsigned int c, unsigned int polygonPointIndex, 991 const cv::Point2f &side1StartVertex, const cv::Point2f &side1EndVertex, 992 const cv::Point2f &side2StartVertex, const cv::Point2f &side2EndVertex, 993 cv::Point2f &intersectionPoint1, cv::Point2f &intersectionPoint2) { 994 std::vector<double> side1Params = lineEquationParameters(side1StartVertex, side1EndVertex); 995 std::vector<double> side2Params = lineEquationParameters(side2StartVertex, side2EndVertex); 996 997 // Compute side C extra parameter using the formula for distance from a point to a line 998 double polygonPointHeight = height(polygonPointIndex, polygon, nrOfPoints, c); 999 double distFormulaDenom = sqrt((side2Params[0] * side2Params[0]) + (side2Params[1] * side2Params[1])); 1000 double sideCExtraParam = 2 * polygonPointHeight * distFormulaDenom; 1001 1002 // Get intersection points if they exist or if lines are identical 1003 if (!areIntersectingLines(side1Params, side2Params, sideCExtraParam, intersectionPoint1, intersectionPoint2)) { 1004 return false; 1005 } else if (areIdenticalLines(side1Params, side2Params, sideCExtraParam)) { 1006 intersectionPoint1 = side1StartVertex; 1007 intersectionPoint2 = side1EndVertex; 1008 } 1009 1010 return true; 1011 } 1012 1013 //! Check if the given lines are identical or not 1014 /*! 1015 * The lines are specified as: 1016 * ax + by + c = 0 1017 * OR 1018 * ax + by + c (+/-) sideCExtraParam = 0 1019 * 1020 * @param side1Params Vector containing the values of a, b and c for side 1 1021 * @param side2Params Vector containing the values of a, b and c for side 2 1022 * @param sideCExtraParam Extra parameter for the flush edge C 1023 */ 1024 static bool areIdenticalLines(const std::vector<double> &side1Params, 1025 const std::vector<double> &side2Params, double sideCExtraParam) { 1026 return ( 1027 (areIdenticalLines(side1Params[0], side1Params[1], -(side1Params[2]), 1028 side2Params[0], side2Params[1], -(side2Params[2]) - sideCExtraParam)) || 1029 (areIdenticalLines(side1Params[0], side1Params[1], -(side1Params[2]), 1030 side2Params[0], side2Params[1], -(side2Params[2]) + sideCExtraParam)) 1031 ); 1032 } 1033 1034 //! Check if the given lines intersect or not. If the lines intersect find their intersection points. 1035 /*! 1036 * The lines are specified as: 1037 * ax + by + c = 0 1038 * OR 1039 * ax + by + c (+/-) sideCExtraParam = 0 1040 * 1041 * @param side1Params Vector containing the values of a, b and c for side 1 1042 * @param side2Params Vector containing the values of a, b and c for side 2 1043 * @param sideCExtraParam Extra parameter for the flush edge C 1044 * @param intersectionPoint1 The first intersection point, if it exists 1045 * @param intersectionPoint2 The second intersection point, if it exists 1046 */ 1047 static bool areIntersectingLines(const std::vector<double> &side1Params, 1048 const std::vector<double> &side2Params, 1049 double sideCExtraParam, cv::Point2f &intersectionPoint1, 1050 cv::Point2f &intersectionPoint2) { 1051 return ( 1052 (lineIntersection(side1Params[0], side1Params[1], -(side1Params[2]), 1053 side2Params[0], side2Params[1], -(side2Params[2]) - sideCExtraParam, 1054 intersectionPoint1)) && 1055 (lineIntersection(side1Params[0], side1Params[1], -(side1Params[2]), 1056 side2Params[0], side2Params[1], -(side2Params[2]) + sideCExtraParam, 1057 intersectionPoint2)) 1058 ); 1059 } 1060 1061 //! Get the line equation parameters "a", "b" and "c" for the line determined by points "p" and "q" 1062 /*! 1063 * The equation of the line is considered in the general form: 1064 * ax + by + c = 0 1065 * 1066 * @param p One point for defining the equation of the line 1067 * @param q Second point for defining the equation of the line 1068 */ 1069 static std::vector<double> lineEquationParameters(const cv::Point2f& p, const cv::Point2f &q) { 1070 std::vector<double> lineEquationParameters; 1071 double a, b, c; 1072 1073 lineEquationDeterminedByPoints(p, q, a, b, c); 1074 1075 lineEquationParameters.push_back(a); 1076 lineEquationParameters.push_back(b); 1077 lineEquationParameters.push_back(c); 1078 1079 return lineEquationParameters; 1080 } 1081 1082 //! Compute the height of the point 1083 /*! 1084 * See paper [2] for more details 1085 * 1086 * @param polygonPoint Polygon point 1087 * @param polygon The polygon representing the convex hull of the points 1088 * @param nrOfPoints Number of points defining the convex polygon 1089 * @param c Index c 1090 */ 1091 static double height(const cv::Point2f &polygonPoint, const std::vector<cv::Point2f> &polygon, 1092 unsigned int nrOfPoints, unsigned int c) { 1093 cv::Point2f pointC = polygon[c]; 1094 cv::Point2f pointCPredecessor = polygon[predecessor(c, nrOfPoints)]; 1095 1096 return distanceFromPointToLine(polygonPoint, pointC, pointCPredecessor); 1097 } 1098 1099 //! Compute the height of the point specified by the given index 1100 /*! 1101 * See paper [2] for more details 1102 * 1103 * @param polygonPointIndex Index of the polygon point 1104 * @param polygon The polygon representing the convex hull of the points 1105 * @param nrOfPoints Number of points defining the convex polygon 1106 * @param c Index c 1107 */ 1108 static double height(unsigned int polygonPointIndex, const std::vector<cv::Point2f> &polygon, 1109 unsigned int nrOfPoints, unsigned int c) { 1110 cv::Point2f pointC = polygon[c]; 1111 cv::Point2f pointCPredecessor = polygon[predecessor(c, nrOfPoints)]; 1112 1113 cv::Point2f polygonPoint = polygon[polygonPointIndex]; 1114 1115 return distanceFromPointToLine(polygonPoint, pointC, pointCPredecessor); 1116 } 1117 1118 //! Advance the given index with one position 1119 /*! 1120 * @param index Index of the point 1121 * @param nrOfPoints Number of points defining the convex polygon 1122 */ 1123 static void advance(unsigned int &index, unsigned int nrOfPoints) { 1124 index = successor(index, nrOfPoints); 1125 } 1126 1127 //! Return the succesor of the provided point index 1128 /*! 1129 * The succesor of the last polygon point is the first polygon point 1130 * (circular referencing) 1131 * 1132 * @param index Index of the point 1133 * @param nrOfPoints Number of points defining the convex polygon 1134 */ 1135 static unsigned int successor(unsigned int index, unsigned int nrOfPoints) { 1136 return ((index + 1) % nrOfPoints); 1137 } 1138 1139 //! Return the predecessor of the provided point index 1140 /*! 1141 * The predecessor of the first polygon point is the last polygon point 1142 * (circular referencing) 1143 * 1144 * @param index Index of the point 1145 * @param nrOfPoints Number of points defining the convex polygon 1146 */ 1147 static unsigned int predecessor(unsigned int index, unsigned int nrOfPoints) { 1148 return (index == 0) ? (nrOfPoints - 1) 1149 : (index - 1); 1150 } 1151 1152 //! Check if the flush edge angle/opposite angle lie between the predecessor and successor angle 1153 /*! 1154 * Check if the angle of the flush edge or its opposite angle lie between the angle of 1155 * the predecessor and successor 1156 * 1157 * @param angleFlushEdge Angle of the flush edge 1158 * @param anglePred Angle of the predecessor 1159 * @param angleSucc Angle of the successor 1160 */ 1161 static bool isFlushAngleBtwPredAndSucc(double &angleFlushEdge, double anglePred, double angleSucc) { 1162 if (isAngleBetweenNonReflex(angleFlushEdge, anglePred, angleSucc)) { 1163 return true; 1164 } else if (isOppositeAngleBetweenNonReflex(angleFlushEdge, anglePred, angleSucc)) { 1165 angleFlushEdge = oppositeAngle(angleFlushEdge); 1166 1167 return true; 1168 } 1169 1170 return false; 1171 } 1172 1173 //! Check if the angle of the line (gamma(p) p) or its opposite angle is equal to the given angle 1174 /*! 1175 * @param gammaAngle Angle of the line (gamma(p) p) 1176 * @param angle Angle to compare against 1177 */ 1178 static bool isGammaAngleEqualTo(double &gammaAngle, double angle) { 1179 return (almostEqual(gammaAngle, angle)); 1180 } 1181 1182 //! Check if the angle of the line (gamma(p) p) or its opposite angle lie between angle1 and angle2 1183 /*! 1184 * @param gammaAngle Angle of the line (gamma(p) p) 1185 * @param angle1 One of the boundary angles 1186 * @param angle2 Another boundary angle 1187 */ 1188 static bool isGammaAngleBtw(double &gammaAngle, double angle1, double angle2) { 1189 return (isAngleBetweenNonReflex(gammaAngle, angle1, angle2)); 1190 } 1191 1192 //! Get the angle of the line measured from the Ox axis in counterclockwise direction 1193 /*! 1194 * The line is specified by points "a" and "b". The value of the angle is expressed in degrees. 1195 * 1196 * @param a Point a 1197 * @param b Point b 1198 */ 1199 static double angleOfLineWrtOxAxis(const cv::Point2f &a, const cv::Point2f &b) { 1200 double y = b.y - a.y; 1201 double x = b.x - a.x; 1202 1203 double angle = (std::atan2(y, x) * 180 / CV_PI); 1204 1205 return (angle < 0) ? (angle + 360) 1206 : angle; 1207 } 1208 1209 //! Check if angle1 lies between non reflex angle determined by angles 2 and 3 1210 /*! 1211 * @param angle1 The angle which lies between angle2 and angle3 or not 1212 * @param angle2 One of the boundary angles 1213 * @param angle3 The other boundary angle 1214 */ 1215 static bool isAngleBetweenNonReflex(double angle1, double angle2, double angle3) { 1216 if (std::abs(angle2 - angle3) > 180) { 1217 if (angle2 > angle3) { 1218 return (((angle2 < angle1) && (lessOrEqual(angle1, 360))) || 1219 ((lessOrEqual(0, angle1)) && (angle1 < angle3))); 1220 } else { 1221 return (((angle3 < angle1) && (lessOrEqual(angle1, 360))) || 1222 ((lessOrEqual(0, angle1)) && (angle1 < angle2))); 1223 } 1224 } else { 1225 return isAngleBetween(angle1, angle2, angle3); 1226 } 1227 } 1228 1229 //! Check if the opposite of angle1, ((angle1 + 180) % 360), lies between non reflex angle determined by angles 2 and 3 1230 /*! 1231 * @param angle1 The angle which lies between angle2 and angle3 or not 1232 * @param angle2 One of the boundary angles 1233 * @param angle3 The other boundary angle 1234 */ 1235 static bool isOppositeAngleBetweenNonReflex(double angle1, double angle2, double angle3) { 1236 double angle1Opposite = oppositeAngle(angle1); 1237 1238 return (isAngleBetweenNonReflex(angle1Opposite, angle2, angle3)); 1239 } 1240 1241 //! Check if angle1 lies between angles 2 and 3 1242 /*! 1243 * @param angle1 The angle which lies between angle2 and angle3 or not 1244 * @param angle2 One of the boundary angles 1245 * @param angle3 The other boundary angle 1246 */ 1247 static bool isAngleBetween(double angle1, double angle2, double angle3) { 1248 if ((((int)(angle2 - angle3)) % 180) > 0) { 1249 return ((angle3 < angle1) && (angle1 < angle2)); 1250 } else { 1251 return ((angle2 < angle1) && (angle1 < angle3)); 1252 } 1253 } 1254 1255 //! Return the angle opposite to the given angle 1256 /*! 1257 * if (angle < 180) then 1258 * return (angle + 180); 1259 * else 1260 * return (angle - 180); 1261 * endif 1262 * 1263 * @param angle Angle 1264 */ 1265 static double oppositeAngle(double angle) { 1266 return (angle > 180) ? (angle - 180) 1267 : (angle + 180); 1268 } 1269 1270 //! Compute the distance from a point "a" to a line specified by two points "B" and "C" 1271 /*! 1272 * Formula used: 1273 * 1274 * |(x_c - x_b) * (y_b - y_a) - (x_b - x_a) * (y_c - y_b)| 1275 * d = ------------------------------------------------------- 1276 * sqrt(((x_c - x_b)^2) + ((y_c - y_b)^2)) 1277 * 1278 * Reference: http://mathworld.wolfram.com/Point-LineDistance2-Dimensional.html 1279 * (Last access: 15.09.2013) 1280 * 1281 * @param a Point from which the distance is measures 1282 * @param linePointB One of the points determining the line 1283 * @param linePointC One of the points determining the line 1284 */ 1285 static double distanceFromPointToLine(const cv::Point2f &a, const cv::Point2f &linePointB, 1286 const cv::Point2f &linePointC) { 1287 double term1 = linePointC.x - linePointB.x; 1288 double term2 = linePointB.y - a.y; 1289 double term3 = linePointB.x - a.x; 1290 double term4 = linePointC.y - linePointB.y; 1291 1292 double nominator = std::abs((term1 * term2) - (term3 * term4)); 1293 double denominator = std::sqrt((term1 * term1) + (term4 * term4)); 1294 1295 return (denominator != 0) ? (nominator / denominator) 1296 : 0; 1297 } 1298 1299 //! Compute the distance between two points 1300 /*! Compute the Euclidean distance between two points 1301 * 1302 * @param a Point a 1303 * @param b Point b 1304 */ 1305 static double distanceBtwPoints(const cv::Point2f &a, const cv::Point2f &b) { 1306 double xDiff = a.x - b.x; 1307 double yDiff = a.y - b.y; 1308 1309 return std::sqrt((xDiff * xDiff) + (yDiff * yDiff)); 1310 } 1311 1312 //! Compute the area of a triangle defined by three points 1313 /*! 1314 * The area is computed using the determinant method. 1315 * An example is depicted at http://demonstrations.wolfram.com/TheAreaOfATriangleUsingADeterminant/ 1316 * (Last access: 15.09.2013) 1317 * 1318 * @param a Point a 1319 * @param b Point b 1320 * @param c Point c 1321 */ 1322 static double areaOfTriangle(const cv::Point2f &a, const cv::Point2f &b, const cv::Point2f &c) { 1323 double posTerm = (a.x * b.y) + (a.y * c.x) + (b.x * c.y); 1324 double negTerm = (b.y * c.x) + (a.x * c.y) + (a.y * b.x); 1325 1326 double determinant = posTerm - negTerm; 1327 1328 return std::abs(determinant) / 2; 1329 } 1330 1331 //! Get the point in the middle of the segment determined by points "a" and "b" 1332 /*! 1333 * @param a Point a 1334 * @param b Point b 1335 */ 1336 static cv::Point2f middlePoint(const cv::Point2f &a, const cv::Point2f &b) { 1337 double middleX = static_cast<double>((a.x + b.x) / 2); 1338 double middleY = static_cast<double>((a.y + b.y) / 2); 1339 1340 return cv::Point2f(static_cast<float>(middleX), static_cast<float>(middleY)); 1341 } 1342 1343 //! Determine the intersection point of two lines, if this point exists 1344 /*! Two lines intersect if they are not parallel (Parallel lines intersect at 1345 * +/- infinity, but we do not consider this case here). 1346 * 1347 * The lines are specified in the following form: 1348 * A1x + B1x = C1 1349 * A2x + B2x = C2 1350 * 1351 * If det (= A1*B2 - A2*B1) == 0, then lines are parallel 1352 * else they intersect 1353 * 1354 * If they intersect, then let us denote the intersection point with P(x, y) where: 1355 * x = (C1*B2 - C2*B1) / (det) 1356 * y = (C2*A1 - C1*A2) / (det) 1357 * 1358 * @param a1 A1 1359 * @param b1 B1 1360 * @param c1 C1 1361 * @param a2 A2 1362 * @param b2 B2 1363 * @param c2 C2 1364 * @param intersection The intersection point, if this point exists 1365 */ 1366 static bool lineIntersection(double a1, double b1, double c1, double a2, double b2, double c2, 1367 cv::Point2f &intersection) { 1368 double det = (a1 * b2) - (a2 * b1); 1369 1370 if (!(almostEqual(det, 0))) { 1371 intersection.x = static_cast<float>(((c1 * b2) - (c2 * b1)) / (det)); 1372 intersection.y = static_cast<float>(((c2 * a1) - (c1 * a2)) / (det)); 1373 1374 return true; 1375 } 1376 1377 return false; 1378 } 1379 1380 //! Determine the intersection point of two lines, if this point exists 1381 /*! Two lines intersect if they are not parallel (Parallel lines intersect at 1382 * +/- infinity, but we do not consider this case here). 1383 * 1384 * The lines are specified by a pair of points each. If they intersect, then 1385 * the function returns true, else it returns false. 1386 * 1387 * Lines can be specified in the following form: 1388 * A1x + B1x = C1 1389 * A2x + B2x = C2 1390 * 1391 * If det (= A1*B2 - A2*B1) == 0, then lines are parallel 1392 * else they intersect 1393 * 1394 * If they intersect, then let us denote the intersection point with P(x, y) where: 1395 * x = (C1*B2 - C2*B1) / (det) 1396 * y = (C2*A1 - C1*A2) / (det) 1397 * 1398 * @param a1 First point for determining the first line 1399 * @param b1 Second point for determining the first line 1400 * @param a2 First point for determining the second line 1401 * @param b2 Second point for determining the second line 1402 * @param intersection The intersection point, if this point exists 1403 */ 1404 static bool lineIntersection(const cv::Point2f &a1, const cv::Point2f &b1, const cv::Point2f &a2, 1405 const cv::Point2f &b2, cv::Point2f &intersection) { 1406 double A1 = b1.y - a1.y; 1407 double B1 = a1.x - b1.x; 1408 double C1 = (a1.x * A1) + (a1.y * B1); 1409 1410 double A2 = b2.y - a2.y; 1411 double B2 = a2.x - b2.x; 1412 double C2 = (a2.x * A2) + (a2.y * B2); 1413 1414 double det = (A1 * B2) - (A2 * B1); 1415 1416 if (!almostEqual(det, 0)) { 1417 intersection.x = static_cast<float>(((C1 * B2) - (C2 * B1)) / (det)); 1418 intersection.y = static_cast<float>(((C2 * A1) - (C1 * A2)) / (det)); 1419 1420 return true; 1421 } 1422 1423 return false; 1424 } 1425 1426 //! Get the values of "a", "b" and "c" of the line equation ax + by + c = 0 knowing that point "p" and "q" are on the line 1427 /*! 1428 * a = q.y - p.y 1429 * b = p.x - q.x 1430 * c = - (p.x * a) - (p.y * b) 1431 * 1432 * @param p Point p 1433 * @param q Point q 1434 * @param a Parameter "a" from the line equation 1435 * @param b Parameter "b" from the line equation 1436 * @param c Parameter "c" from the line equation 1437 */ 1438 static void lineEquationDeterminedByPoints(const cv::Point2f &p, const cv::Point2f &q, 1439 double &a, double &b, double &c) { 1440 CV_Assert(areEqualPoints(p, q) == false); 1441 1442 a = q.y - p.y; 1443 b = p.x - q.x; 1444 c = ((-p.y) * b) - (p.x * a); 1445 } 1446 1447 //! Check if p1 and p2 are on the same side of the line determined by points a and b 1448 /*! 1449 * @param p1 Point p1 1450 * @param p2 Point p2 1451 * @param a First point for determining line 1452 * @param b Second point for determining line 1453 */ 1454 static bool areOnTheSameSideOfLine(const cv::Point2f &p1, const cv::Point2f &p2, 1455 const cv::Point2f &a, const cv::Point2f &b) { 1456 double a1, b1, c1; 1457 1458 lineEquationDeterminedByPoints(a, b, a1, b1, c1); 1459 1460 double p1OnLine = (a1 * p1.x) + (b1 * p1.y) + c1; 1461 double p2OnLine = (a1 * p2.x) + (b1 * p2.y) + c1; 1462 1463 return (sign(p1OnLine) == sign(p2OnLine)); 1464 } 1465 1466 //! Check if one point lies between two other points 1467 /*! 1468 * @param point Point lying possibly outside the line segment 1469 * @param lineSegmentStart First point determining the line segment 1470 * @param lineSegmentEnd Second point determining the line segment 1471 */ 1472 static bool isPointOnLineSegment(const cv::Point2f &point, const cv::Point2f &lineSegmentStart, 1473 const cv::Point2f &lineSegmentEnd) { 1474 double d1 = distanceBtwPoints(point, lineSegmentStart); 1475 double d2 = distanceBtwPoints(point, lineSegmentEnd); 1476 double lineSegmentLength = distanceBtwPoints(lineSegmentStart, lineSegmentEnd); 1477 1478 return (almostEqual(d1 + d2, lineSegmentLength)); 1479 } 1480 1481 //! Check if two lines are identical 1482 /*! 1483 * Lines are be specified in the following form: 1484 * A1x + B1x = C1 1485 * A2x + B2x = C2 1486 * 1487 * If (A1/A2) == (B1/B2) == (C1/C2), then the lines are identical 1488 * else they are not 1489 * 1490 * @param a1 A1 1491 * @param b1 B1 1492 * @param c1 C1 1493 * @param a2 A2 1494 * @param b2 B2 1495 * @param c2 C2 1496 */ 1497 static bool areIdenticalLines(double a1, double b1, double c1, double a2, double b2, double c2) { 1498 double a1B2 = a1 * b2; 1499 double a2B1 = a2 * b1; 1500 double a1C2 = a1 * c2; 1501 double a2C1 = a2 * c1; 1502 double b1C2 = b1 * c2; 1503 double b2C1 = b2 * c1; 1504 1505 return ((almostEqual(a1B2, a2B1)) && (almostEqual(b1C2, b2C1)) && (almostEqual(a1C2, a2C1))); 1506 } 1507 1508 //! Check if points point1 and point2 are equal or not 1509 /*! 1510 * @param point1 One point 1511 * @param point2 The other point 1512 */ 1513 static bool areEqualPoints(const cv::Point2f &point1, const cv::Point2f &point2) { 1514 return (almostEqual(point1.x, point2.x) && almostEqual(point1.y, point2.y)); 1515 } 1516 1517 //! Return the sign of the number 1518 /*! 1519 * The sign function returns: 1520 * -1, if number < 0 1521 * +1, if number > 0 1522 * 0, otherwise 1523 */ 1524 static int sign(double number) { 1525 return (number > 0) ? 1 : ((number < 0) ? -1 : 0); 1526 } 1527 1528 //! Return the maximum of the provided numbers 1529 static double maximum(double number1, double number2, double number3) { 1530 return std::max(std::max(number1, number2), number3); 1531 } 1532 1533 //! Check if the two numbers are equal (almost) 1534 /*! 1535 * The expression for determining if two real numbers are equal is: 1536 * if (Abs(x - y) <= EPSILON * Max(1.0f, Abs(x), Abs(y))). 1537 * 1538 * @param number1 First number 1539 * @param number2 Second number 1540 */ 1541 static bool almostEqual(double number1, double number2) { 1542 return (std::abs(number1 - number2) <= (EPSILON * maximum(1.0, std::abs(number1), std::abs(number2)))); 1543 } 1544 1545 //! Check if the first number is greater than or equal to the second number 1546 /*! 1547 * @param number1 First number 1548 * @param number2 Second number 1549 */ 1550 static bool greaterOrEqual(double number1, double number2) { 1551 return ((number1 > number2) || (almostEqual(number1, number2))); 1552 } 1553 1554 //! Check if the first number is less than or equal to the second number 1555 /*! 1556 * @param number1 First number 1557 * @param number2 Second number 1558 */ 1559 static bool lessOrEqual(double number1, double number2) { 1560 return ((number1 < number2) || (almostEqual(number1, number2))); 1561 } 1562 1563 } 1564
