root/modules/imgproc/src/min_enclosing_triangle.cpp

/* [<][>][^][v][top][bottom][index][help] */

DEFINITIONS

This source file includes following definitions.
  1. minEnclosingTriangle
  2. findMinEnclosingTriangle
  3. createConvexHull
  4. findMinEnclosingTriangle
  5. copyResultingTriangle
  6. initialise
  7. findMinimumAreaEnclosingTriangle
  8. returnMinimumAreaEnclosingTriangle
  9. advanceBToRightChain
  10. moveAIfLowAndBIfHigh
  11. searchForBTangency
  12. isNotBTangency
  13. updateSidesCA
  14. updateSidesBA
  15. updateSideB
  16. isLocalMinimalTriangle
  17. isValidMinimalTriangle
  18. updateMinimumAreaEnclosingTriangle
  19. middlePointOfSideB
  20. intersectsBelow
  21. intersectsAbove
  22. intersects
  23. intersectsAboveOrBelow
  24. gamma
  25. findVertexCOnSideB
  26. findGammaIntersectionPoints
  27. areIdenticalLines
  28. areIntersectingLines
  29. lineEquationParameters
  30. height
  31. height
  32. advance
  33. successor
  34. predecessor
  35. isFlushAngleBtwPredAndSucc
  36. isGammaAngleEqualTo
  37. isGammaAngleBtw
  38. angleOfLineWrtOxAxis
  39. isAngleBetweenNonReflex
  40. isOppositeAngleBetweenNonReflex
  41. isAngleBetween
  42. oppositeAngle
  43. distanceFromPointToLine
  44. distanceBtwPoints
  45. areaOfTriangle
  46. middlePoint
  47. lineIntersection
  48. lineIntersection
  49. lineEquationDeterminedByPoints
  50. areOnTheSameSideOfLine
  51. isPointOnLineSegment
  52. areIdenticalLines
  53. areEqualPoints
  54. sign
  55. maximum
  56. almostEqual
  57. greaterOrEqual
  58. lessOrEqual

/*M///////////////////////////////////////////////////////////////////////////////////////
//
//  IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
//
//  By downloading, copying, installing or using the software you agree to this license.
//  If you do not agree to this license, do not download, install,
//  copy or use the software.
//
//  INFORMATION REGARDING THE CONTRIBUTION:
//
//  Author: Ovidiu Parvu
//  Affiliation: Brunel University
//  Created: 11.09.2013
//  E-mail: <ovidiu.parvu[AT]gmail.com>
//  Web: http://people.brunel.ac.uk/~cspgoop
//
//  These functions were implemented during Ovidiu Parvu's first year as a PhD student at
//  Brunel University, London, UK. The PhD project is supervised by prof. David Gilbert (principal)
//  and prof. Nigel Saunders (second).
//
//  THE IMPLEMENTATION OF THE MODULES IS BASED ON THE FOLLOWING PAPERS:
//
//  [1] V. Klee and M. C. Laskowski, "Finding the smallest triangles containing a given convex
//  polygon", Journal of Algorithms, vol. 6, no. 3, pp. 359-375, Sep. 1985.
//  [2] J. O'Rourke, A. Aggarwal, S. Maddila, and M. Baldwin, "An optimal algorithm for finding
//  minimal enclosing triangles", Journal of Algorithms, vol. 7, no. 2, pp. 258-269, Jun. 1986.
//
//  The overall complexity of the algorithm is theta(n) where "n" represents the number
//  of vertices in the convex polygon.
//
//
//
//                           License Agreement
//                For Open Source Computer Vision Library
//
// Copyright (C) 2000, Intel Corporation, all rights reserved.
// Copyright (C) 2013, OpenCV Foundation, all rights reserved.
// Copyright (C) 2013, Ovidiu Parvu, all rights reserved.
// Third party copyrights are property of their respective owners.
//
// Redistribution and use in source and binary forms, with or without modification,
// are permitted provided that the following conditions are met:
//
//   * Redistribution's of source code must retain the above copyright notice,
//     this list of conditions and the following disclaimer.
//
//   * Redistribution's in binary form must reproduce the above copyright notice,
//     this list of conditions and the following disclaimer in the documentation
//     and/or other materials provided with the distribution.
//
//   * The name of the copyright holders may not be used to endorse or promote products
//     derived from this software without specific prior written permission.
//
// This software is provided by the copyright holders and contributors "as is" and
// any express or implied warranties, including, but not limited to, the implied
// warranties of merchantability and fitness for a particular purpose are disclaimed.
// In no event shall the Intel Corporation or contributors be liable for any direct,
// indirect, incidental, special, exemplary, or consequential damages
// (including, but not limited to, procurement of substitute goods or services;
// loss of use, data, or profits; or business interruption) however caused
// and on any theory of liability, whether in contract, strict liability,
// or tort (including negligence or otherwise) arising in any way out of
// the use of this software, even if advised of the possibility of such damage.
//
//M*/

#include "precomp.hpp"

#include <algorithm>
#include <cmath>
#include <limits>
#include <vector>


///////////////////////////////// Constants definitions //////////////////////////////////


// Intersection of line and polygon

#define INTERSECTS_BELOW        1
#define INTERSECTS_ABOVE        2
#define INTERSECTS_CRITICAL     3

// Error messages

#define ERR_SIDE_B_GAMMA        "The position of side B could not be determined, because gamma(b) could not be computed."
#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."

// Possible values for validation flag

#define VALIDATION_SIDE_A_TANGENT   0
#define VALIDATION_SIDE_B_TANGENT   1
#define VALIDATION_SIDES_FLUSH      2

// Threshold value for comparisons

#define EPSILON 1E-5


////////////////////////////// Helper functions declarations /////////////////////////////


namespace minEnclosingTriangle {

static void advance(unsigned int &index, unsigned int nrOfPoints);

static void advanceBToRightChain(const std::vector<cv::Point2f> &polygon,
                                 unsigned int nrOfPoints, unsigned int &b,
                                 unsigned int c);

static bool almostEqual(double number1, double number2);

static double angleOfLineWrtOxAxis(const cv::Point2f &a, const cv::Point2f &b);

static bool areEqualPoints(const cv::Point2f &point1, const cv::Point2f &point2);

static bool areIdenticalLines(const std::vector<double> &side1Params,
                              const std::vector<double> &side2Params, double sideCExtraParam);

static bool areIdenticalLines(double a1, double b1, double c1, double a2, double b2, double c2);

static bool areIntersectingLines(const std::vector<double> &side1Params,
                                 const std::vector<double> &side2Params,
                                 double sideCExtraParam, cv::Point2f &intersectionPoint1,
                                 cv::Point2f &intersectionPoint2);

static bool areOnTheSameSideOfLine(const cv::Point2f &p1, const cv::Point2f &p2,
                                   const cv::Point2f &a, const cv::Point2f &b);

static double areaOfTriangle(const cv::Point2f &a, const cv::Point2f &b, const cv::Point2f &c);

static void copyResultingTriangle(const std::vector<cv::Point2f> &resultingTriangle, cv::OutputArray triangle);

static void createConvexHull(cv::InputArray points, std::vector<cv::Point2f> &polygon);

static double distanceBtwPoints(const cv::Point2f &a, const cv::Point2f &b);

static double distanceFromPointToLine(const cv::Point2f &a, const cv::Point2f &linePointB,
                                      const cv::Point2f &linePointC);

static bool findGammaIntersectionPoints(const std::vector<cv::Point2f> &polygon, unsigned int nrOfPoints,
                                        unsigned int c, unsigned int polygonPointIndex,
                                        const cv::Point2f &side1StartVertex, const cv::Point2f &side1EndVertex,
                                        const cv::Point2f &side2StartVertex, const cv::Point2f &side2EndVertex,
                                        cv::Point2f &intersectionPoint1, cv::Point2f &intersectionPoint2);

static void findMinEnclosingTriangle(cv::InputArray points,
                                     CV_OUT cv::OutputArray triangle, CV_OUT double &area);

static void findMinEnclosingTriangle(const std::vector<cv::Point2f> &polygon,
                                     std::vector<cv::Point2f> &triangle, double &area);

static void findMinimumAreaEnclosingTriangle(const std::vector<cv::Point2f> &polygon,
                                             std::vector<cv::Point2f> &triangle, double &area);

static cv::Point2f findVertexCOnSideB(const std::vector<cv::Point2f> &polygon, unsigned int nrOfPoints,
                                      unsigned int a, unsigned int c,
                                      const cv::Point2f &sideBStartVertex,
                                      const cv::Point2f &sideBEndVertex,
                                      const cv::Point2f &sideCStartVertex,
                                      const cv::Point2f &sideCEndVertex);

static bool gamma(unsigned int polygonPointIndex, cv::Point2f &gammaPoint,
                  const std::vector<cv::Point2f> &polygon, unsigned int nrOfPoints,
                  unsigned int a, unsigned int c);

static bool greaterOrEqual(double number1, double number2);

static double height(const cv::Point2f &polygonPoint, const std::vector<cv::Point2f> &polygon,
                     unsigned int nrOfPoints, unsigned int c);

static double height(unsigned int polygonPointIndex, const std::vector<cv::Point2f> &polygon,
                     unsigned int nrOfPoints, unsigned int c);

static void initialise(std::vector<cv::Point2f> &triangle, double &area);

static unsigned int intersects(double angleGammaAndPoint, unsigned int polygonPointIndex,
                               const std::vector<cv::Point2f> &polygon, unsigned int nrOfPoints,
                               unsigned int c);

static bool intersectsAbove(const cv::Point2f &gammaPoint, unsigned int polygonPointIndex,
                            const std::vector<cv::Point2f> &polygon, unsigned int nrOfPoints,
                            unsigned int c);

static unsigned int intersectsAboveOrBelow(unsigned int succPredIndex, unsigned int pointIndex,
                                           const std::vector<cv::Point2f> &polygon,
                                           unsigned int nrOfPoints, unsigned int c);

static bool intersectsBelow(const cv::Point2f &gammaPoint, unsigned int polygonPointIndex,
                            const std::vector<cv::Point2f> &polygon, unsigned int nrOfPoints,
                            unsigned int c);

static bool isAngleBetween(double angle1, double angle2, double angle3);

static bool isAngleBetweenNonReflex(double angle1, double angle2, double angle3);

static bool isFlushAngleBtwPredAndSucc(double &angleFlushEdge, double anglePred, double angleSucc);

static bool isGammaAngleBtw(double &gammaAngle, double angle1, double angle2);

static bool isGammaAngleEqualTo(double &gammaAngle, double angle);

static bool isLocalMinimalTriangle(cv::Point2f &vertexA, cv::Point2f &vertexB,
                                   cv::Point2f &vertexC, const std::vector<cv::Point2f> &polygon,
                                   unsigned int nrOfPoints, unsigned int a, unsigned int b,
                                   unsigned int validationFlag, const cv::Point2f &sideAStartVertex,
                                   const cv::Point2f &sideAEndVertex, const cv::Point2f &sideBStartVertex,
                                   const cv::Point2f &sideBEndVertex, const cv::Point2f &sideCStartVertex,
                                   const cv::Point2f &sideCEndVertex);

static bool isNotBTangency(const std::vector<cv::Point2f> &polygon,
                           unsigned int nrOfPoints, unsigned int a, unsigned int b,
                           unsigned int c);

static bool isOppositeAngleBetweenNonReflex(double angle1, double angle2, double angle3);

static bool isPointOnLineSegment(const cv::Point2f &point, const cv::Point2f &lineSegmentStart,
                                 const cv::Point2f &lineSegmentEnd);

static bool isValidMinimalTriangle(const cv::Point2f &vertexA, const cv::Point2f &vertexB,
                                   const cv::Point2f &vertexC, const std::vector<cv::Point2f> &polygon,
                                   unsigned int nrOfPoints, unsigned int a, unsigned int b,
                                   unsigned int validationFlag, const cv::Point2f &sideAStartVertex,
                                   const cv::Point2f &sideAEndVertex, const cv::Point2f &sideBStartVertex,
                                   const cv::Point2f &sideBEndVertex, const cv::Point2f &sideCStartVertex,
                                   const cv::Point2f &sideCEndVertex);

static bool lessOrEqual(double number1, double number2);

static void lineEquationDeterminedByPoints(const cv::Point2f &p, const cv::Point2f &q,
                                           double &a, double &b, double &c);

static std::vector<double> lineEquationParameters(const cv::Point2f& p, const cv::Point2f &q);

static bool lineIntersection(const cv::Point2f &a1, const cv::Point2f &b1, const cv::Point2f &a2,
                             const cv::Point2f &b2, cv::Point2f &intersection);

static bool lineIntersection(double a1, double b1, double c1, double a2, double b2, double c2,
                             cv::Point2f &intersection);

static double maximum(double number1, double number2, double number3);

static cv::Point2f middlePoint(const cv::Point2f &a, const cv::Point2f &b);

static bool middlePointOfSideB(cv::Point2f &middlePointOfSideB, const cv::Point2f &sideAStartVertex,
                               const cv::Point2f &sideAEndVertex, const cv::Point2f &sideBStartVertex,
                               const cv::Point2f &sideBEndVertex, const cv::Point2f &sideCStartVertex,
                               const cv::Point2f &sideCEndVertex);

static void moveAIfLowAndBIfHigh(const std::vector<cv::Point2f> &polygon,
                                 unsigned int nrOfPoints, unsigned int &a, unsigned int &b,
                                 unsigned int c);

static double oppositeAngle(double angle);

static unsigned int predecessor(unsigned int index, unsigned int nrOfPoints);

static void returnMinimumAreaEnclosingTriangle(const std::vector<cv::Point2f> &polygon,
                                               std::vector<cv::Point2f> &triangle, double &area);

static void searchForBTangency(const std::vector<cv::Point2f> &polygon,
                               unsigned int nrOfPoints, unsigned int a, unsigned int &b,
                               unsigned int c);

static int sign(double number);

static unsigned int successor(unsigned int index, unsigned int nrOfPoints);

static void updateMinimumAreaEnclosingTriangle(std::vector<cv::Point2f> &triangle, double &area,
                                               const cv::Point2f &vertexA, const cv::Point2f &vertexB,
                                               const cv::Point2f &vertexC);

static void updateSideB(const std::vector<cv::Point2f> &polygon,
                        unsigned int nrOfPoints, unsigned int a, unsigned int b,
                        unsigned int c, unsigned int &validationFlag,
                        cv::Point2f &sideBStartVertex, cv::Point2f &sideBEndVertex);

static void updateSidesBA(const std::vector<cv::Point2f> &polygon,
                          unsigned int nrOfPoints, unsigned int a, unsigned int b,
                          unsigned int c, unsigned int &validationFlag,
                          cv::Point2f &sideAStartVertex, cv::Point2f &sideAEndVertex,
                          cv::Point2f &sideBStartVertex, cv::Point2f &sideBEndVertex,
                          const cv::Point2f &sideCStartVertex, const cv::Point2f &sideCEndVertex);

static void updateSidesCA(const std::vector<cv::Point2f> &polygon,
                          unsigned int nrOfPoints, unsigned int a, unsigned int c,
                          cv::Point2f &sideAStartVertex, cv::Point2f &sideAEndVertex,
                          cv::Point2f &sideCStartVertex, cv::Point2f &sideCEndVertex);

}


///////////////////////////////////// Main functions /////////////////////////////////////


//! Find the minimum enclosing triangle for the given set of points and return its area
/*!
* @param points         Set of points
* @param triangle       Minimum area triangle enclosing the given set of points
*/
double cv::minEnclosingTriangle(cv::InputArray points, CV_OUT cv::OutputArray triangle) {
    double area;

    minEnclosingTriangle::findMinEnclosingTriangle(points, triangle, area);

    return area;
}


/////////////////////////////// Helper functions definition //////////////////////////////


namespace minEnclosingTriangle {

//! Find the minimum enclosing triangle and its area
/*!
* @param points         Set of points
* @param triangle       Minimum area triangle enclosing the given set of points
* @param area           Area of the minimum area enclosing triangle
*/
static void findMinEnclosingTriangle(cv::InputArray points,
                                     CV_OUT cv::OutputArray triangle, CV_OUT double &area) {
    std::vector<cv::Point2f> resultingTriangle, polygon;

    createConvexHull(points, polygon);
    findMinEnclosingTriangle(polygon, resultingTriangle, area);
    copyResultingTriangle(resultingTriangle, triangle);
}

//! Create the convex hull of the given set of points
/*!
* @param points     The provided set of points
* @param polygon    The polygon representing the convex hull of the points
*/
static void createConvexHull(cv::InputArray points, std::vector<cv::Point2f> &polygon) {
    cv::Mat pointsMat = points.getMat();
    std::vector<cv::Point2f> pointsVector;

    CV_Assert((pointsMat.checkVector(2) > 0) &&
              ((pointsMat.depth() == CV_32F) || (pointsMat.depth() == CV_32S)));

    pointsMat.convertTo(pointsVector, CV_32F);

    convexHull(pointsVector, polygon, true, true);
}

//! Find the minimum enclosing triangle and its area
/*!
* The overall complexity of the algorithm is theta(n) where "n" represents the number
* of vertices in the convex polygon
*
* @param polygon    The polygon representing the convex hull of the points
* @param triangle   Minimum area triangle enclosing the given polygon
* @param area       Area of the minimum area enclosing triangle
*/
static void findMinEnclosingTriangle(const std::vector<cv::Point2f> &polygon,
                                     std::vector<cv::Point2f> &triangle, double &area) {
    initialise(triangle, area);

    if (polygon.size() > 3) {
        findMinimumAreaEnclosingTriangle(polygon, triangle, area);
    } else {
        returnMinimumAreaEnclosingTriangle(polygon, triangle, area);
    }
}

//! Copy resultingTriangle to the OutputArray triangle
/*!
* @param resultingTriangle  Minimum area triangle enclosing the given polygon found by the algorithm
* @param triangle           Minimum area triangle enclosing the given polygon returned to the user
*/
static void copyResultingTriangle(const std::vector<cv::Point2f> &resultingTriangle,
                                  cv::OutputArray triangle) {
    cv::Mat(resultingTriangle).copyTo(triangle);
}

//! Initialisation function
/*!
* @param triangle       Minimum area triangle enclosing the given polygon
* @param area           Area of the minimum area enclosing triangle
*/
static void initialise(std::vector<cv::Point2f> &triangle, double &area) {
    area = std::numeric_limits<double>::max();

    // Clear all points previously stored in the vector
    triangle.clear();
}

//! Find the minimum area enclosing triangle for the given polygon
/*!
* @param polygon    The polygon representing the convex hull of the points
* @param triangle   Minimum area triangle enclosing the given polygon
* @param area       Area of the minimum area enclosing triangle
*/
static void findMinimumAreaEnclosingTriangle(const std::vector<cv::Point2f> &polygon,
                                             std::vector<cv::Point2f> &triangle, double &area) {
    // Algorithm specific variables

    unsigned int validationFlag;

    cv::Point2f vertexA, vertexB, vertexC;

    cv::Point2f sideAStartVertex, sideAEndVertex;
    cv::Point2f sideBStartVertex, sideBEndVertex;
    cv::Point2f sideCStartVertex, sideCEndVertex;

    unsigned int a, b, c;
    unsigned int nrOfPoints;

    // Variables initialisation

    nrOfPoints = static_cast<unsigned int>(polygon.size());

    a = 1;
    b = 2;
    c = 0;

    // Main algorithm steps

    for (c = 0; c < nrOfPoints; c++) {
        advanceBToRightChain(polygon, nrOfPoints, b, c);
        moveAIfLowAndBIfHigh(polygon, nrOfPoints, a, b, c);
        searchForBTangency(polygon, nrOfPoints, a ,b, c);

        updateSidesCA(polygon, nrOfPoints, a, c, sideAStartVertex, sideAEndVertex,
                      sideCStartVertex, sideCEndVertex);

        if (isNotBTangency(polygon, nrOfPoints, a, b, c)) {
            updateSidesBA(polygon, nrOfPoints, a, b, c, validationFlag, sideAStartVertex,
                          sideAEndVertex, sideBStartVertex, sideBEndVertex,
                          sideCStartVertex, sideCEndVertex);
        } else {
            updateSideB(polygon, nrOfPoints, a, b, c, validationFlag,
                        sideBStartVertex,  sideBEndVertex);
        }

        if (isLocalMinimalTriangle(vertexA, vertexB, vertexC, polygon, nrOfPoints, a, b,
                                   validationFlag, sideAStartVertex, sideAEndVertex,
                                   sideBStartVertex, sideBEndVertex, sideCStartVertex,
                                   sideCEndVertex)) {
            updateMinimumAreaEnclosingTriangle(triangle, area, vertexA, vertexB, vertexC);
        }
    }
}

//! Return the minimum area enclosing (pseudo-)triangle in case the convex polygon has at most three points
/*!
* @param polygon    The polygon representing the convex hull of the points
* @param triangle   Minimum area triangle enclosing the given polygon
* @param area       Area of the minimum area enclosing triangle
*/
static void returnMinimumAreaEnclosingTriangle(const std::vector<cv::Point2f> &polygon,
                                               std::vector<cv::Point2f> &triangle, double &area) {
    unsigned int nrOfPoints = static_cast<unsigned int>(polygon.size());

    for (int i = 0; i < 3; i++) {
        triangle.push_back(polygon[i % nrOfPoints]);
    }

    area = areaOfTriangle(triangle[0], triangle[1], triangle[2]);
}

//! Advance b to the right chain
/*!
* See paper [2] for more details
*
* @param polygon            The polygon representing the convex hull of the points
* @param nrOfPoints         Number of points defining the convex polygon
* @param b                  Index b
* @param c                  Index c
*/
static void advanceBToRightChain(const std::vector<cv::Point2f> &polygon,
                                 unsigned int nrOfPoints, unsigned int &b,
                                 unsigned int c) {
    while (greaterOrEqual(height(successor(b, nrOfPoints), polygon, nrOfPoints, c),
                          height(b, polygon, nrOfPoints, c))) {
        advance(b, nrOfPoints);
    }
}

//! Move "a" if it is low and "b" if it is high
/*!
* See paper [2] for more details
*
* @param polygon            The polygon representing the convex hull of the points
* @param nrOfPoints         Number of points defining the convex polygon
* @param a                  Index a
* @param b                  Index b
* @param c                  Index c
*/
static void moveAIfLowAndBIfHigh(const std::vector<cv::Point2f> &polygon,
                                 unsigned int nrOfPoints, unsigned int &a, unsigned int &b,
                                 unsigned int c) {
    cv::Point2f gammaOfA;

    while(height(b, polygon, nrOfPoints, c) > height(a, polygon, nrOfPoints, c)) {
        if ((gamma(a, gammaOfA, polygon, nrOfPoints, a, c)) && (intersectsBelow(gammaOfA, b, polygon, nrOfPoints, c))) {
            advance(b, nrOfPoints);
        } else {
            advance(a, nrOfPoints);
        }
    }
}

//! Search for the tangency of side B
/*!
* See paper [2] for more details
*
* @param polygon            The polygon representing the convex hull of the points
* @param nrOfPoints         Number of points defining the convex polygon
* @param a                  Index a
* @param b                  Index b
* @param c                  Index c
*/
static void searchForBTangency(const std::vector<cv::Point2f> &polygon,
                               unsigned int nrOfPoints, unsigned int a, unsigned int &b,
                               unsigned int c) {
    cv::Point2f gammaOfB;

    while (((gamma(b, gammaOfB, polygon, nrOfPoints, a, c)) &&
            (intersectsBelow(gammaOfB, b, polygon, nrOfPoints, c))) &&
           (greaterOrEqual(height(b, polygon, nrOfPoints, c),
                           height(predecessor(a, nrOfPoints), polygon, nrOfPoints, c)))
          ) {
        advance(b, nrOfPoints);
    }
}

//! Check if tangency for side B was not obtained
/*!
* See paper [2] for more details
*
* @param polygon            The polygon representing the convex hull of the points
* @param nrOfPoints         Number of points defining the convex polygon
* @param a                  Index a
* @param b                  Index b
* @param c                  Index c
*/
static bool isNotBTangency(const std::vector<cv::Point2f> &polygon,
                           unsigned int nrOfPoints, unsigned int a, unsigned int b,
                           unsigned int c) {
    cv::Point2f gammaOfB;

    if (((gamma(b, gammaOfB, polygon, nrOfPoints, a, c)) &&
         (intersectsAbove(gammaOfB, b, polygon, nrOfPoints, c))) ||
        (height(b, polygon, nrOfPoints, c) < height(predecessor(a, nrOfPoints), polygon, nrOfPoints, c))) {
        return true;
    }

    return false;
}

//! Update sides A and C
/*!
* Side C will have as start and end vertices the polygon points "c" and "c-1"
* Side A will have as start and end vertices the polygon points "a" and "a-1"
*
* @param polygon            The polygon representing the convex hull of the points
* @param nrOfPoints         Number of points defining the convex polygon
* @param a                  Index a
* @param c                  Index c
* @param sideAStartVertex   Start vertex for defining side A
* @param sideAEndVertex     End vertex for defining side A
* @param sideCStartVertex   Start vertex for defining side C
* @param sideCEndVertex     End vertex for defining side C
*/
static void updateSidesCA(const std::vector<cv::Point2f> &polygon,
                          unsigned int nrOfPoints, unsigned int a, unsigned int c,
                          cv::Point2f &sideAStartVertex, cv::Point2f &sideAEndVertex,
                          cv::Point2f &sideCStartVertex, cv::Point2f &sideCEndVertex) {
    sideCStartVertex = polygon[predecessor(c, nrOfPoints)];
    sideCEndVertex = polygon[c];

    sideAStartVertex = polygon[predecessor(a, nrOfPoints)];
    sideAEndVertex = polygon[a];
}

//! Update sides B and possibly A if tangency for side B was not obtained
/*!
* See paper [2] for more details
*
* @param polygon            The polygon representing the convex hull of the points
* @param nrOfPoints         Number of points defining the convex polygon
* @param a                  Index a
* @param b                  Index b
* @param c                  Index c
* @param validationFlag     Flag used for validation
* @param sideAStartVertex   Start vertex for defining side A
* @param sideAEndVertex     End vertex for defining side A
* @param sideBStartVertex   Start vertex for defining side B
* @param sideBEndVertex     End vertex for defining side B
* @param sideCStartVertex   Start vertex for defining side C
* @param sideCEndVertex     End vertex for defining side C
*/
static void updateSidesBA(const std::vector<cv::Point2f> &polygon,
                          unsigned int nrOfPoints, unsigned int a, unsigned int b,
                          unsigned int c, unsigned int &validationFlag,
                          cv::Point2f &sideAStartVertex, cv::Point2f &sideAEndVertex,
                          cv::Point2f &sideBStartVertex, cv::Point2f &sideBEndVertex,
                          const cv::Point2f &sideCStartVertex, const cv::Point2f &sideCEndVertex) {
    // Side B is flush with edge [b, b-1]
    sideBStartVertex = polygon[predecessor(b, nrOfPoints)];
    sideBEndVertex = polygon[b];

    // Find middle point of side B
    cv::Point2f sideBMiddlePoint;

    if ((middlePointOfSideB(sideBMiddlePoint, sideAStartVertex, sideAEndVertex, sideBStartVertex,
                            sideBEndVertex, sideCStartVertex, sideCEndVertex)) &&
        (height(sideBMiddlePoint, polygon, nrOfPoints, c) <
         height(predecessor(a, nrOfPoints), polygon, nrOfPoints, c))) {
        sideAStartVertex = polygon[predecessor(a, nrOfPoints)];
        sideAEndVertex = findVertexCOnSideB(polygon, nrOfPoints, a, c,
                                            sideBStartVertex, sideBEndVertex,
                                            sideCStartVertex, sideCEndVertex);

        validationFlag = VALIDATION_SIDE_A_TANGENT;
    } else {
        validationFlag = VALIDATION_SIDES_FLUSH;
    }
}

//! Set side B if tangency for side B was obtained
/*!
* See paper [2] for more details
*
* @param polygon            The polygon representing the convex hull of the points
* @param nrOfPoints         Number of points defining the convex polygon
* @param a                  Index a
* @param b                  Index b
* @param c                  Index c
* @param validationFlag     Flag used for validation
* @param sideBStartVertex   Start vertex for defining side B
* @param sideBEndVertex     End vertex for defining side B
*/
static void updateSideB(const std::vector<cv::Point2f> &polygon,
                        unsigned int nrOfPoints, unsigned int a, unsigned int b,
                        unsigned int c, unsigned int &validationFlag,
                        cv::Point2f &sideBStartVertex, cv::Point2f &sideBEndVertex) {
    if (!gamma(b, sideBStartVertex, polygon, nrOfPoints, a, c)) {
        CV_Error(cv::Error::StsInternal, ERR_SIDE_B_GAMMA);
    }

    sideBEndVertex = polygon[b];

    validationFlag = VALIDATION_SIDE_B_TANGENT;
}

//! Update the triangle vertices after all sides were set and check if a local minimal triangle was found or not
/*!
* See paper [2] for more details
*
* @param vertexA            Vertex A of the enclosing triangle
* @param vertexB            Vertex B of the enclosing triangle
* @param vertexC            Vertex C of the enclosing triangle
* @param polygon            The polygon representing the convex hull of the points
* @param nrOfPoints         Number of points defining the convex polygon
* @param a                  Index a
* @param b                  Index b
* @param validationFlag     Flag used for validation
* @param sideAStartVertex   Start vertex for defining side A
* @param sideAEndVertex     End vertex for defining side A
* @param sideBStartVertex   Start vertex for defining side B
* @param sideBEndVertex     End vertex for defining side B
* @param sideCStartVertex   Start vertex for defining side C
* @param sideCEndVertex     End vertex for defining side C
*/
static bool isLocalMinimalTriangle(cv::Point2f &vertexA, cv::Point2f &vertexB,
                                   cv::Point2f &vertexC, const std::vector<cv::Point2f> &polygon,
                                   unsigned int nrOfPoints, unsigned int a, unsigned int b,
                                   unsigned int validationFlag, const cv::Point2f &sideAStartVertex,
                                   const cv::Point2f &sideAEndVertex, const cv::Point2f &sideBStartVertex,
                                   const cv::Point2f &sideBEndVertex, const cv::Point2f &sideCStartVertex,
                                   const cv::Point2f &sideCEndVertex) {
    if ((!lineIntersection(sideAStartVertex, sideAEndVertex,
                           sideBStartVertex, sideBEndVertex, vertexC)) ||
        (!lineIntersection(sideAStartVertex, sideAEndVertex,
                           sideCStartVertex, sideCEndVertex, vertexB)) ||
        (!lineIntersection(sideBStartVertex, sideBEndVertex,
                           sideCStartVertex, sideCEndVertex, vertexA))) {
        return false;
    }

    return isValidMinimalTriangle(vertexA, vertexB, vertexC, polygon, nrOfPoints, a, b,
                                  validationFlag, sideAStartVertex, sideAEndVertex,
                                  sideBStartVertex, sideBEndVertex, sideCStartVertex,
                                  sideCEndVertex);
}

//! Check if the found minimal triangle is valid
/*!
* This means that all midpoints of the triangle should touch the polygon
*
* See paper [2] for more details
*
* @param vertexA            Vertex A of the enclosing triangle
* @param vertexB            Vertex B of the enclosing triangle
* @param vertexC            Vertex C of the enclosing triangle
* @param polygon            The polygon representing the convex hull of the points
* @param nrOfPoints         Number of points defining the convex polygon
* @param a                  Index a
* @param b                  Index b
* @param validationFlag     Flag used for validation
* @param sideAStartVertex   Start vertex for defining side A
* @param sideAEndVertex     End vertex for defining side A
* @param sideBStartVertex   Start vertex for defining side B
* @param sideBEndVertex     End vertex for defining side B
* @param sideCStartVertex   Start vertex for defining side C
* @param sideCEndVertex     End vertex for defining side C
*/
static bool isValidMinimalTriangle(const cv::Point2f &vertexA, const cv::Point2f &vertexB,
                                   const cv::Point2f &vertexC, const std::vector<cv::Point2f> &polygon,
                                   unsigned int nrOfPoints, unsigned int a, unsigned int b,
                                   unsigned int validationFlag, const cv::Point2f &sideAStartVertex,
                                   const cv::Point2f &sideAEndVertex, const cv::Point2f &sideBStartVertex,
                                   const cv::Point2f &sideBEndVertex, const cv::Point2f &sideCStartVertex,
                                   const cv::Point2f &sideCEndVertex) {
    cv::Point2f midpointSideA = middlePoint(vertexB, vertexC);
    cv::Point2f midpointSideB = middlePoint(vertexA, vertexC);
    cv::Point2f midpointSideC = middlePoint(vertexA, vertexB);

    bool sideAValid = (validationFlag == VALIDATION_SIDE_A_TANGENT)
                        ? (areEqualPoints(midpointSideA, polygon[predecessor(a, nrOfPoints)]))
                        : (isPointOnLineSegment(midpointSideA, sideAStartVertex, sideAEndVertex));

    bool sideBValid = (validationFlag == VALIDATION_SIDE_B_TANGENT)
                          ? (areEqualPoints(midpointSideB, polygon[b]))
                          : (isPointOnLineSegment(midpointSideB, sideBStartVertex, sideBEndVertex));

    bool sideCValid = isPointOnLineSegment(midpointSideC, sideCStartVertex, sideCEndVertex);

    return (sideAValid && sideBValid && sideCValid);
}

//! Update the current minimum area enclosing triangle if the newly obtained one has a smaller area
/*!
* @param triangle   Minimum area triangle enclosing the given polygon
* @param area       Area of the minimum area triangle enclosing the given polygon
* @param vertexA    Vertex A of the enclosing triangle
* @param vertexB    Vertex B of the enclosing triangle
* @param vertexC    Vertex C of the enclosing triangle
*/
static void updateMinimumAreaEnclosingTriangle(std::vector<cv::Point2f> &triangle, double &area,
                                               const cv::Point2f &vertexA, const cv::Point2f &vertexB,
                                               const cv::Point2f &vertexC) {
    double triangleArea = areaOfTriangle(vertexA, vertexB, vertexC);

    if (triangleArea < area) {
        triangle.clear();

        triangle.push_back(vertexA);
        triangle.push_back(vertexB);
        triangle.push_back(vertexC);

        area = triangleArea;
    }
}

//! Return the middle point of side B
/*!
* @param middlePointOfSideB Middle point of side B
* @param sideAStartVertex   Start vertex for defining side A
* @param sideAEndVertex     End vertex for defining side A
* @param sideBStartVertex   Start vertex for defining side B
* @param sideBEndVertex     End vertex for defining side B
* @param sideCStartVertex   Start vertex for defining side C
* @param sideCEndVertex     End vertex for defining side C
*/
static bool middlePointOfSideB(cv::Point2f &middlePointOfSideB, const cv::Point2f &sideAStartVertex,
                               const cv::Point2f &sideAEndVertex, const cv::Point2f &sideBStartVertex,
                               const cv::Point2f &sideBEndVertex, const cv::Point2f &sideCStartVertex,
                               const cv::Point2f &sideCEndVertex) {
    cv::Point2f vertexA, vertexC;

    if ((!lineIntersection(sideBStartVertex, sideBEndVertex, sideCStartVertex, sideCEndVertex, vertexA)) ||
        (!lineIntersection(sideBStartVertex, sideBEndVertex, sideAStartVertex, sideAEndVertex, vertexC))) {
        return false;
    }

    middlePointOfSideB = middlePoint(vertexA, vertexC);

    return true;
}

//! Check if the line intersects below
/*!
* Check if the line determined by gammaPoint and polygon[polygonPointIndex] intersects
* the polygon below the point polygon[polygonPointIndex]
*
* @param gammaPoint         Gamma(p)
* @param polygonPointIndex  Index of the polygon point which is considered when determining the line
* @param polygon            The polygon representing the convex hull of the points
* @param nrOfPoints         Number of points defining the convex polygon
* @param c                  Index c
*/
static bool intersectsBelow(const cv::Point2f &gammaPoint, unsigned int polygonPointIndex,
                            const std::vector<cv::Point2f> &polygon, unsigned int nrOfPoints,
                            unsigned int c) {
    double angleOfGammaAndPoint = angleOfLineWrtOxAxis(polygon[polygonPointIndex], gammaPoint);

    return (intersects(angleOfGammaAndPoint, polygonPointIndex, polygon, nrOfPoints, c) == INTERSECTS_BELOW);
}

//! Check if the line intersects above
/*!
* Check if the line determined by gammaPoint and polygon[polygonPointIndex] intersects
* the polygon above the point polygon[polygonPointIndex]
*
* @param gammaPoint         Gamma(p)
* @param polygonPointIndex  Index of the polygon point which is considered when determining the line
* @param polygon            The polygon representing the convex hull of the points
* @param nrOfPoints         Number of points defining the convex polygon
* @param c                  Index c
*/
static bool intersectsAbove(const cv::Point2f &gammaPoint, unsigned int polygonPointIndex,
                            const std::vector<cv::Point2f> &polygon, unsigned int nrOfPoints,
                            unsigned int c) {
    double angleOfGammaAndPoint = angleOfLineWrtOxAxis(gammaPoint, polygon[polygonPointIndex]);

    return (intersects(angleOfGammaAndPoint, polygonPointIndex, polygon, nrOfPoints, c) == INTERSECTS_ABOVE);
}

//! Check if/where the line determined by gammaPoint and polygon[polygonPointIndex] intersects the polygon
/*!
* @param angleGammaAndPoint     Angle determined by gammaPoint and polygon[polygonPointIndex] wrt Ox axis
* @param polygonPointIndex      Index of the polygon point which is considered when determining the line
* @param polygon                The polygon representing the convex hull of the points
* @param nrOfPoints             Number of points defining the convex polygon
* @param c                      Index c
*/
static unsigned int intersects(double angleGammaAndPoint, unsigned int polygonPointIndex,
                               const std::vector<cv::Point2f> &polygon, unsigned int nrOfPoints,
                               unsigned int c) {
    double anglePointPredecessor = angleOfLineWrtOxAxis(polygon[predecessor(polygonPointIndex, nrOfPoints)],
                                                        polygon[polygonPointIndex]);
    double anglePointSuccessor   = angleOfLineWrtOxAxis(polygon[successor(polygonPointIndex, nrOfPoints)],
                                                        polygon[polygonPointIndex]);
    double angleFlushEdge        = angleOfLineWrtOxAxis(polygon[predecessor(c, nrOfPoints)],
                                                        polygon[c]);

    if (isFlushAngleBtwPredAndSucc(angleFlushEdge, anglePointPredecessor, anglePointSuccessor)) {
        if ((isGammaAngleBtw(angleGammaAndPoint, anglePointPredecessor, angleFlushEdge)) ||
            (almostEqual(angleGammaAndPoint, anglePointPredecessor))) {
            return intersectsAboveOrBelow(predecessor(polygonPointIndex, nrOfPoints),
                                          polygonPointIndex, polygon, nrOfPoints, c);
        } else if ((isGammaAngleBtw(angleGammaAndPoint, anglePointSuccessor, angleFlushEdge)) ||
                  (almostEqual(angleGammaAndPoint, anglePointSuccessor))) {
            return intersectsAboveOrBelow(successor(polygonPointIndex, nrOfPoints),
                                          polygonPointIndex, polygon, nrOfPoints, c);
        }
    } else {
        if (
            (isGammaAngleBtw(angleGammaAndPoint, anglePointPredecessor, anglePointSuccessor)) ||
            (
                (isGammaAngleEqualTo(angleGammaAndPoint, anglePointPredecessor)) &&
                (!isGammaAngleEqualTo(angleGammaAndPoint, angleFlushEdge))
            ) ||
            (
                (isGammaAngleEqualTo(angleGammaAndPoint, anglePointSuccessor)) &&
                (!isGammaAngleEqualTo(angleGammaAndPoint, angleFlushEdge))
            )
           ) {
            return INTERSECTS_BELOW;
        }
    }

    return INTERSECTS_CRITICAL;
}

//! If (gamma(x) x) intersects P between successorOrPredecessorIndex and pointIntex is it above/below?
/*!
* @param succPredIndex  Index of the successor or predecessor
* @param pointIndex     Index of the point x in the polygon
* @param polygon        The polygon representing the convex hull of the points
* @param nrOfPoints     Number of points defining the convex polygon
* @param c              Index c
*/
static unsigned int intersectsAboveOrBelow(unsigned int succPredIndex, unsigned int pointIndex,
                                           const std::vector<cv::Point2f> &polygon,
                                           unsigned int nrOfPoints, unsigned int c) {
    if (height(succPredIndex, polygon, nrOfPoints, c) > height(pointIndex, polygon, nrOfPoints, c)) {
        return INTERSECTS_ABOVE;
    } else {
        return INTERSECTS_BELOW;
    }
}

//! Find gamma for a given point "p" specified by its index
/*!
* The function returns true if gamma exists i.e. if lines (a a-1) and (x y) intersect
* and false otherwise. In case the two lines intersect in point intersectionPoint, gamma is computed.
*
* Considering that line (x y) is a line parallel to (c c-1) and that the distance between the lines is equal
* to 2 * height(p), we can have two possible (x y) lines.
*
* Therefore, we will compute two intersection points between the lines (x y) and (a a-1) and take the
* point which is on the same side of line (c c-1) as the polygon.
*
* See paper [2] and formula for distance from point to a line for more details
*
* @param polygonPointIndex Index of the polygon point
* @param gammaPoint        Point gamma(polygon[polygonPointIndex])
* @param polygon           The polygon representing the convex hull of the points
* @param nrOfPoints        Number of points defining the convex polygon
* @param a                 Index a
* @param c                 Index c
*/
static bool gamma(unsigned int polygonPointIndex, cv::Point2f &gammaPoint,
                  const std::vector<cv::Point2f> &polygon, unsigned int nrOfPoints,
                  unsigned int a, unsigned int c) {
    cv::Point2f intersectionPoint1, intersectionPoint2;

    // Get intersection points if they exist
    if (!findGammaIntersectionPoints(polygon, nrOfPoints, c, polygonPointIndex,
                                     polygon[a], polygon[predecessor(a, nrOfPoints)],
                                     polygon[c], polygon[predecessor(c, nrOfPoints)],
                                     intersectionPoint1, intersectionPoint2)) {
        return false;
    }

    // Select the point which is on the same side of line C as the polygon
    if (areOnTheSameSideOfLine(intersectionPoint1, polygon[successor(c, nrOfPoints)],
                               polygon[c], polygon[predecessor(c, nrOfPoints)])) {
        gammaPoint = intersectionPoint1;
    } else {
        gammaPoint = intersectionPoint2;
    }

    return true;
}

//! Find vertex C which lies on side B at a distance = 2 * height(a-1) from side C
/*!
* Considering that line (x y) is a line parallel to (c c-1) and that the distance between the lines is equal
* to 2 * height(a-1), we can have two possible (x y) lines.
*
* Therefore, we will compute two intersection points between the lines (x y) and (b b-1) and take the
* point which is on the same side of line (c c-1) as the polygon.
*
* See paper [2] and formula for distance from point to a line for more details
*
* @param polygon            The polygon representing the convex hull of the points
* @param nrOfPoints         Number of points defining the convex polygon
* @param a                  Index a
* @param c                  Index c
* @param sideBStartVertex   Start vertex for defining side B
* @param sideBEndVertex     End vertex for defining side B
* @param sideCStartVertex   Start vertex for defining side C
* @param sideCEndVertex     End vertex for defining side C
*/
static cv::Point2f findVertexCOnSideB(const std::vector<cv::Point2f> &polygon, unsigned int nrOfPoints,
                                      unsigned int a, unsigned int c,
                                      const cv::Point2f &sideBStartVertex,
                                      const cv::Point2f &sideBEndVertex,
                                      const cv::Point2f &sideCStartVertex,
                                      const cv::Point2f &sideCEndVertex) {
    cv::Point2f intersectionPoint1, intersectionPoint2;

    // Get intersection points if they exist
    if (!findGammaIntersectionPoints(polygon, nrOfPoints, c, predecessor(a, nrOfPoints),
                                     sideBStartVertex, sideBEndVertex,
                                     sideCStartVertex, sideCEndVertex,
                                     intersectionPoint1, intersectionPoint2)) {
        CV_Error(cv::Error::StsInternal, ERR_VERTEX_C_ON_SIDE_B);
    }

    // Select the point which is on the same side of line C as the polygon
    if (areOnTheSameSideOfLine(intersectionPoint1, polygon[successor(c, nrOfPoints)],
                               polygon[c], polygon[predecessor(c, nrOfPoints)])) {
        return intersectionPoint1;
    } else {
        return intersectionPoint2;
    }
}

//! Find the intersection points to compute gamma(point)
/*!
* @param polygon                The polygon representing the convex hull of the points
* @param nrOfPoints             Number of points defining the convex polygon
* @param c                      Index c
* @param polygonPointIndex      Index of the polygon point for which the distance is known
* @param side1StartVertex       Start vertex for side 1
* @param side1EndVertex         End vertex for side 1
* @param side2StartVertex       Start vertex for side 2
* @param side2EndVertex         End vertex for side 2
* @param intersectionPoint1     First intersection point between one pair of lines
* @param intersectionPoint2     Second intersection point between other pair of lines
*/
static bool findGammaIntersectionPoints(const std::vector<cv::Point2f> &polygon, unsigned int nrOfPoints,
                                        unsigned int c, unsigned int polygonPointIndex,
                                        const cv::Point2f &side1StartVertex, const cv::Point2f &side1EndVertex,
                                        const cv::Point2f &side2StartVertex, const cv::Point2f &side2EndVertex,
                                        cv::Point2f &intersectionPoint1, cv::Point2f &intersectionPoint2) {
    std::vector<double> side1Params = lineEquationParameters(side1StartVertex, side1EndVertex);
    std::vector<double> side2Params = lineEquationParameters(side2StartVertex, side2EndVertex);

    // Compute side C extra parameter using the formula for distance from a point to a line
    double polygonPointHeight = height(polygonPointIndex, polygon, nrOfPoints, c);
    double distFormulaDenom = sqrt((side2Params[0] * side2Params[0]) + (side2Params[1] * side2Params[1]));
    double sideCExtraParam = 2 * polygonPointHeight * distFormulaDenom;

    // Get intersection points if they exist or if lines are identical
    if (!areIntersectingLines(side1Params, side2Params, sideCExtraParam, intersectionPoint1, intersectionPoint2)) {
        return false;
    } else if (areIdenticalLines(side1Params, side2Params, sideCExtraParam)) {
        intersectionPoint1 = side1StartVertex;
        intersectionPoint2 = side1EndVertex;
    }

    return true;
}

//! Check if the given lines are identical or not
/*!
* The lines are specified as:
*      ax + by + c = 0
*  OR
*      ax + by + c (+/-) sideCExtraParam = 0
*
* @param side1Params       Vector containing the values of a, b and c for side 1
* @param side2Params       Vector containing the values of a, b and c for side 2
* @param sideCExtraParam   Extra parameter for the flush edge C
*/
static bool areIdenticalLines(const std::vector<double> &side1Params,
                              const std::vector<double> &side2Params, double sideCExtraParam) {
    return (
        (areIdenticalLines(side1Params[0], side1Params[1], -(side1Params[2]),
                           side2Params[0], side2Params[1], -(side2Params[2]) - sideCExtraParam)) ||
        (areIdenticalLines(side1Params[0], side1Params[1], -(side1Params[2]),
                           side2Params[0], side2Params[1], -(side2Params[2]) + sideCExtraParam))
    );
}

//! Check if the given lines intersect or not. If the lines intersect find their intersection points.
/*!
* The lines are specified as:
*      ax + by + c = 0
*  OR
*      ax + by + c (+/-) sideCExtraParam = 0
*
* @param side1Params           Vector containing the values of a, b and c for side 1
* @param side2Params           Vector containing the values of a, b and c for side 2
* @param sideCExtraParam       Extra parameter for the flush edge C
* @param intersectionPoint1    The first intersection point, if it exists
* @param intersectionPoint2    The second intersection point, if it exists
*/
static bool areIntersectingLines(const std::vector<double> &side1Params,
                                 const std::vector<double> &side2Params,
                                 double sideCExtraParam, cv::Point2f &intersectionPoint1,
                                 cv::Point2f &intersectionPoint2) {
    return (
        (lineIntersection(side1Params[0], side1Params[1], -(side1Params[2]),
                          side2Params[0], side2Params[1], -(side2Params[2]) - sideCExtraParam,
                          intersectionPoint1)) &&
        (lineIntersection(side1Params[0], side1Params[1], -(side1Params[2]),
                          side2Params[0], side2Params[1], -(side2Params[2]) + sideCExtraParam,
                          intersectionPoint2))
    );
}

//! Get the line equation parameters "a", "b" and "c" for the line determined by points "p" and "q"
/*!
* The equation of the line is considered in the general form:
* ax + by + c = 0
*
* @param p One point for defining the equation of the line
* @param q Second point for defining the equation of the line
*/
static std::vector<double> lineEquationParameters(const cv::Point2f& p, const cv::Point2f &q) {
    std::vector<double> lineEquationParameters;
    double a, b, c;

    lineEquationDeterminedByPoints(p, q, a, b, c);

    lineEquationParameters.push_back(a);
    lineEquationParameters.push_back(b);
    lineEquationParameters.push_back(c);

    return lineEquationParameters;
}

//! Compute the height of the point
/*!
* See paper [2] for more details
*
* @param polygonPoint       Polygon point
* @param polygon            The polygon representing the convex hull of the points
* @param nrOfPoints         Number of points defining the convex polygon
* @param c                  Index c
*/
static double height(const cv::Point2f &polygonPoint, const std::vector<cv::Point2f> &polygon,
                     unsigned int nrOfPoints, unsigned int c) {
    cv::Point2f pointC = polygon[c];
    cv::Point2f pointCPredecessor = polygon[predecessor(c, nrOfPoints)];

    return distanceFromPointToLine(polygonPoint, pointC, pointCPredecessor);
}

//! Compute the height of the point specified by the given index
/*!
* See paper [2] for more details
*
* @param polygonPointIndex  Index of the polygon point
* @param polygon            The polygon representing the convex hull of the points
* @param nrOfPoints         Number of points defining the convex polygon
* @param c                  Index c
*/
static double height(unsigned int polygonPointIndex, const std::vector<cv::Point2f> &polygon,
                     unsigned int nrOfPoints, unsigned int c) {
    cv::Point2f pointC = polygon[c];
    cv::Point2f pointCPredecessor = polygon[predecessor(c, nrOfPoints)];

    cv::Point2f polygonPoint = polygon[polygonPointIndex];

    return distanceFromPointToLine(polygonPoint, pointC, pointCPredecessor);
}

//! Advance the given index with one position
/*!
* @param index          Index of the point
* @param nrOfPoints     Number of points defining the convex polygon
*/
static void advance(unsigned int &index, unsigned int nrOfPoints) {
    index = successor(index, nrOfPoints);
}

//! Return the succesor of the provided point index
/*!
* The succesor of the last polygon point is the first polygon point
* (circular referencing)
*
* @param index          Index of the point
* @param nrOfPoints     Number of points defining the convex polygon
*/
static unsigned int successor(unsigned int index, unsigned int nrOfPoints) {
    return ((index + 1) % nrOfPoints);
}

//! Return the predecessor of the provided point index
/*!
* The predecessor of the first polygon point is the last polygon point
* (circular referencing)
*
* @param index          Index of the point
* @param nrOfPoints     Number of points defining the convex polygon
*/
static unsigned int predecessor(unsigned int index, unsigned int nrOfPoints) {
    return (index == 0) ? (nrOfPoints - 1)
                        : (index - 1);
}

//! Check if the flush edge angle/opposite angle lie between the predecessor and successor angle
/*!
* Check if the angle of the flush edge or its opposite angle lie between the angle of
* the predecessor and successor
*
* @param angleFlushEdge    Angle of the flush edge
* @param anglePred         Angle of the predecessor
* @param angleSucc         Angle of the successor
*/
static bool isFlushAngleBtwPredAndSucc(double &angleFlushEdge, double anglePred, double angleSucc) {
    if (isAngleBetweenNonReflex(angleFlushEdge, anglePred, angleSucc)) {
        return true;
    } else if (isOppositeAngleBetweenNonReflex(angleFlushEdge, anglePred, angleSucc)) {
        angleFlushEdge = oppositeAngle(angleFlushEdge);

        return true;
    }

    return false;
}

//! Check if the angle of the line (gamma(p) p) or its opposite angle is equal to the given angle
/*!
* @param gammaAngle    Angle of the line (gamma(p) p)
* @param angle         Angle to compare against
*/
static bool isGammaAngleEqualTo(double &gammaAngle, double angle) {
    return (almostEqual(gammaAngle, angle));
}

//! Check if the angle of the line (gamma(p) p) or its opposite angle lie between angle1 and angle2
/*!
* @param gammaAngle    Angle of the line (gamma(p) p)
* @param angle1        One of the boundary angles
* @param angle2        Another boundary angle
*/
static bool isGammaAngleBtw(double &gammaAngle, double angle1, double angle2) {
    return (isAngleBetweenNonReflex(gammaAngle, angle1, angle2));
}

//! Get the angle of the line measured from the Ox axis in counterclockwise direction
/*!
* The line is specified by points "a" and "b". The value of the angle is expressed in degrees.
*
* @param a Point a
* @param b Point b
*/
static double angleOfLineWrtOxAxis(const cv::Point2f &a, const cv::Point2f &b) {
    double y = b.y - a.y;
    double x = b.x - a.x;

    double angle = (std::atan2(y, x) * 180 / CV_PI);

    return (angle < 0) ? (angle + 360)
                       : angle;
}

//! Check if angle1 lies between non reflex angle determined by angles 2 and 3
/*!
* @param angle1 The angle which lies between angle2 and angle3 or not
* @param angle2 One of the boundary angles
* @param angle3 The other boundary angle
*/
static bool isAngleBetweenNonReflex(double angle1, double angle2, double angle3) {
    if (std::abs(angle2 - angle3) > 180) {
        if (angle2 > angle3) {
            return (((angle2 < angle1) && (lessOrEqual(angle1, 360))) ||
                    ((lessOrEqual(0, angle1)) && (angle1 < angle3)));
        } else {
            return (((angle3 < angle1) && (lessOrEqual(angle1, 360))) ||
                    ((lessOrEqual(0, angle1)) && (angle1 < angle2)));
        }
    } else {
        return isAngleBetween(angle1, angle2, angle3);
    }
}

//! Check if the opposite of angle1, ((angle1 + 180) % 360), lies between non reflex angle determined by angles 2 and 3
/*!
* @param angle1 The angle which lies between angle2 and angle3 or not
* @param angle2 One of the boundary angles
* @param angle3 The other boundary angle
*/
static bool isOppositeAngleBetweenNonReflex(double angle1, double angle2, double angle3) {
    double angle1Opposite = oppositeAngle(angle1);

    return (isAngleBetweenNonReflex(angle1Opposite, angle2, angle3));
}

//! Check if angle1 lies between angles 2 and 3
/*!
* @param angle1 The angle which lies between angle2 and angle3 or not
* @param angle2 One of the boundary angles
* @param angle3 The other boundary angle
*/
static bool isAngleBetween(double angle1, double angle2, double angle3) {
    if ((((int)(angle2 - angle3)) % 180) > 0) {
        return ((angle3 < angle1) && (angle1 < angle2));
    } else {
        return ((angle2 < angle1) && (angle1 < angle3));
    }
}

//! Return the angle opposite to the given angle
/*!
* if (angle < 180) then
*      return (angle + 180);
* else
*      return (angle - 180);
* endif
*
* @param angle Angle
*/
static double oppositeAngle(double angle) {
    return (angle > 180) ? (angle - 180)
                         : (angle + 180);
}

//! Compute the distance from a point "a" to a line specified by two points "B" and "C"
/*!
* Formula used:
*
*     |(x_c - x_b) * (y_b - y_a) - (x_b - x_a) * (y_c - y_b)|
* d = -------------------------------------------------------
*            sqrt(((x_c - x_b)^2) + ((y_c - y_b)^2))
*
* Reference: http://mathworld.wolfram.com/Point-LineDistance2-Dimensional.html
* (Last access: 15.09.2013)
*
* @param a             Point from which the distance is measures
* @param linePointB    One of the points determining the line
* @param linePointC    One of the points determining the line
*/
static double distanceFromPointToLine(const cv::Point2f &a, const cv::Point2f &linePointB,
                                      const cv::Point2f &linePointC) {
    double term1 = linePointC.x - linePointB.x;
    double term2 = linePointB.y - a.y;
    double term3 = linePointB.x - a.x;
    double term4 = linePointC.y - linePointB.y;

    double nominator = std::abs((term1 * term2) - (term3 * term4));
    double denominator = std::sqrt((term1 * term1) + (term4 * term4));

    return (denominator != 0) ? (nominator / denominator)
                              : 0;
}

//! Compute the distance between two points
/*! Compute the Euclidean distance between two points
*
* @param a Point a
* @param b Point b
*/
static double distanceBtwPoints(const cv::Point2f &a, const cv::Point2f &b) {
    double xDiff = a.x - b.x;
    double yDiff = a.y - b.y;

    return std::sqrt((xDiff * xDiff) + (yDiff * yDiff));
}

//! Compute the area of a triangle defined by three points
/*!
* The area is computed using the determinant method.
* An example is depicted at http://demonstrations.wolfram.com/TheAreaOfATriangleUsingADeterminant/
* (Last access: 15.09.2013)
*
* @param a Point a
* @param b Point b
* @param c Point c
*/
static double areaOfTriangle(const cv::Point2f &a, const cv::Point2f &b, const cv::Point2f &c) {
    double posTerm = (a.x * b.y) + (a.y * c.x) + (b.x * c.y);
    double negTerm = (b.y * c.x) + (a.x * c.y) + (a.y * b.x);

    double determinant = posTerm - negTerm;

    return std::abs(determinant) / 2;
}

//! Get the point in the middle of the segment determined by points "a" and "b"
/*!
* @param a Point a
* @param b Point b
*/
static cv::Point2f middlePoint(const cv::Point2f &a, const cv::Point2f &b) {
    double middleX = static_cast<double>((a.x + b.x) / 2);
    double middleY = static_cast<double>((a.y + b.y) / 2);

    return cv::Point2f(static_cast<float>(middleX), static_cast<float>(middleY));
}

//! Determine the intersection point of two lines, if this point exists
/*! Two lines intersect if they are not parallel (Parallel lines intersect at
* +/- infinity, but we do not consider this case here).
*
* The lines are specified in the following form:
*      A1x + B1x = C1
*      A2x + B2x = C2
*
* If det (= A1*B2 - A2*B1) == 0, then lines are parallel
*                                else they intersect
*
* If they intersect, then let us denote the intersection point with P(x, y) where:
*      x = (C1*B2 - C2*B1) / (det)
*      y = (C2*A1 - C1*A2) / (det)
*
* @param a1 A1
* @param b1 B1
* @param c1 C1
* @param a2 A2
* @param b2 B2
* @param c2 C2
* @param intersection The intersection point, if this point exists
*/
static bool lineIntersection(double a1, double b1, double c1, double a2, double b2, double c2,
                             cv::Point2f &intersection) {
    double det = (a1 * b2) - (a2 * b1);

    if (!(almostEqual(det, 0))) {
        intersection.x = static_cast<float>(((c1 * b2) - (c2 * b1)) / (det));
        intersection.y = static_cast<float>(((c2 * a1) - (c1 * a2)) / (det));

        return true;
    }

    return false;
}

//! Determine the intersection point of two lines, if this point exists
/*! Two lines intersect if they are not parallel (Parallel lines intersect at
* +/- infinity, but we do not consider this case here).
*
* The lines are specified by a pair of points each. If they intersect, then
* the function returns true, else it returns false.
*
* Lines can be specified in the following form:
*      A1x + B1x = C1
*      A2x + B2x = C2
*
* If det (= A1*B2 - A2*B1) == 0, then lines are parallel
*                                else they intersect
*
* If they intersect, then let us denote the intersection point with P(x, y) where:
*      x = (C1*B2 - C2*B1) / (det)
*      y = (C2*A1 - C1*A2) / (det)
*
* @param a1 First point for determining the first line
* @param b1 Second point for determining the first line
* @param a2 First point for determining the second line
* @param b2 Second point for determining the second line
* @param intersection The intersection point, if this point exists
*/
static bool lineIntersection(const cv::Point2f &a1, const cv::Point2f &b1, const cv::Point2f &a2,
                             const cv::Point2f &b2, cv::Point2f &intersection) {
    double A1 = b1.y - a1.y;
    double B1 = a1.x - b1.x;
    double C1 = (a1.x * A1) + (a1.y * B1);

    double A2 = b2.y - a2.y;
    double B2 = a2.x - b2.x;
    double C2 = (a2.x * A2) + (a2.y * B2);

    double det = (A1 * B2) - (A2 * B1);

    if (!almostEqual(det, 0)) {
        intersection.x = static_cast<float>(((C1 * B2) - (C2 * B1)) / (det));
        intersection.y = static_cast<float>(((C2 * A1) - (C1 * A2)) / (det));

        return true;
    }

    return false;
}

//! 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
/*!
* a = q.y - p.y
* b = p.x - q.x
* c = - (p.x * a) - (p.y * b)
*
* @param p Point p
* @param q Point q
* @param a Parameter "a" from the line equation
* @param b Parameter "b" from the line equation
* @param c Parameter "c" from the line equation
*/
static void lineEquationDeterminedByPoints(const cv::Point2f &p, const cv::Point2f &q,
                                           double &a, double &b, double &c) {
    CV_Assert(areEqualPoints(p, q) == false);

    a = q.y - p.y;
    b = p.x - q.x;
    c = ((-p.y) * b) - (p.x * a);
}

//! Check if p1 and p2 are on the same side of the line determined by points a and b
/*!
* @param p1    Point p1
* @param p2    Point p2
* @param a     First point for determining line
* @param b     Second point for determining line
*/
static bool areOnTheSameSideOfLine(const cv::Point2f &p1, const cv::Point2f &p2,
                                   const cv::Point2f &a, const cv::Point2f &b) {
    double a1, b1, c1;

    lineEquationDeterminedByPoints(a, b, a1, b1, c1);

    double p1OnLine = (a1 * p1.x) + (b1 * p1.y) + c1;
    double p2OnLine = (a1 * p2.x) + (b1 * p2.y) + c1;

    return (sign(p1OnLine) == sign(p2OnLine));
}

//! Check if one point lies between two other points
/*!
* @param point             Point lying possibly outside the line segment
* @param lineSegmentStart  First point determining the line segment
* @param lineSegmentEnd    Second point determining the line segment
*/
static bool isPointOnLineSegment(const cv::Point2f &point, const cv::Point2f &lineSegmentStart,
                                 const cv::Point2f &lineSegmentEnd) {
    double d1 = distanceBtwPoints(point, lineSegmentStart);
    double d2 = distanceBtwPoints(point, lineSegmentEnd);
    double lineSegmentLength = distanceBtwPoints(lineSegmentStart, lineSegmentEnd);

    return (almostEqual(d1 + d2, lineSegmentLength));
}

//! Check if two lines are identical
/*!
* Lines are be specified in the following form:
*      A1x + B1x = C1
*      A2x + B2x = C2
*
* If (A1/A2) == (B1/B2) == (C1/C2), then the lines are identical
*                                   else they are not
*
* @param a1 A1
* @param b1 B1
* @param c1 C1
* @param a2 A2
* @param b2 B2
* @param c2 C2
*/
static bool areIdenticalLines(double a1, double b1, double c1, double a2, double b2, double c2) {
    double a1B2 = a1 * b2;
    double a2B1 = a2 * b1;
    double a1C2 = a1 * c2;
    double a2C1 = a2 * c1;
    double b1C2 = b1 * c2;
    double b2C1 = b2 * c1;

    return ((almostEqual(a1B2, a2B1)) && (almostEqual(b1C2, b2C1)) && (almostEqual(a1C2, a2C1)));
}

//! Check if points point1 and point2 are equal or not
/*!
* @param point1 One point
* @param point2 The other point
*/
static bool areEqualPoints(const cv::Point2f &point1, const cv::Point2f &point2) {
    return (almostEqual(point1.x, point2.x) && almostEqual(point1.y, point2.y));
}

//! Return the sign of the number
/*!
* The sign function returns:
*  -1, if number < 0
*  +1, if number > 0
*  0, otherwise
*/
static int sign(double number) {
    return (number > 0) ? 1 : ((number < 0) ? -1 : 0);
}

//! Return the maximum of the provided numbers
static double maximum(double number1, double number2, double number3) {
    return std::max(std::max(number1, number2), number3);
}

//! Check if the two numbers are equal (almost)
/*!
* The expression for determining if two real numbers are equal is:
* if (Abs(x - y) <= EPSILON * Max(1.0f, Abs(x), Abs(y))).
*
* @param number1 First number
* @param number2 Second number
*/
static bool almostEqual(double number1, double number2) {
    return (std::abs(number1 - number2) <= (EPSILON * maximum(1.0, std::abs(number1), std::abs(number2))));
}

//! Check if the first number is greater than or equal to the second number
/*!
* @param number1 First number
* @param number2 Second number
*/
static bool greaterOrEqual(double number1, double number2) {
    return ((number1 > number2) || (almostEqual(number1, number2)));
}

//! Check if the first number is less than or equal to the second number
/*!
* @param number1 First number
* @param number2 Second number
*/
static bool lessOrEqual(double number1, double number2) {
    return ((number1 < number2) || (almostEqual(number1, number2)));
}

}

/* [<][>][^][v][top][bottom][index][help] */