diff --git a/doc/opencv.bib b/doc/opencv.bib index 6632271e4a..5531bb6dd5 100644 --- a/doc/opencv.bib +++ b/doc/opencv.bib @@ -1467,3 +1467,12 @@ volume = {60}, journal = {ISPRS Journal of Photogrammetry and Remote Sensing} } +@article{LowIlie2003, + author = {Kok-Lim Low, Adrian Ilie}, + year = {2003}, + pages = {3-15}, + title = {View Frustum Optimization to Maximize Object's Image Area}, + journal = {Journal of Graphics, (GPU, & Game) Tools (JGT)}, + volume = {8}, + url = {https://citeseerx.ist.psu.edu/document?repid=rep1&type=pdf&doi=1fbd43f3827fffeb76641a9c5ab5b625eb5a75ba} +} diff --git a/modules/imgproc/include/opencv2/imgproc.hpp b/modules/imgproc/include/opencv2/imgproc.hpp index a4852e80db..2f3c6f344f 100644 --- a/modules/imgproc/include/opencv2/imgproc.hpp +++ b/modules/imgproc/include/opencv2/imgproc.hpp @@ -4058,6 +4058,28 @@ CV_EXPORTS_W void approxPolyDP( InputArray curve, OutputArray approxCurve, double epsilon, bool closed ); +/** @brief Approximates a polygon with a convex hull with a specified accuracy and number of sides. + +The cv::approxPolyN function approximates a polygon with a convex hull +so that the difference between the contour area of the original contour and the new polygon is minimal. +It uses a greedy algorithm for contracting two vertices into one in such a way that the additional area is minimal. +Straight lines formed by each edge of the convex contour are drawn and the areas of the resulting triangles are considered. +Each vertex will lie either on the original contour or outside it. + +The algorithm based on the paper @cite LowIlie2003 . + +@param curve Input vector of a 2D points stored in std::vector or Mat, points must be float or integer. +@param approxCurve Result of the approximation. The type is vector of a 2D point (Point2f or Point) in std::vector or Mat. +@param nsides The parameter defines the number of sides of the result polygon. +@param epsilon_percentage defines the percentage of the maximum of additional area. +If it equals -1, it is not used. Otherwise algorighm stops if additional area is greater than contourArea(_curve) * percentage. +If additional area exceeds the limit, algorithm returns as many vertices as there were at the moment the limit was exceeded. +@param ensure_convex If it is true, algorithm creates a convex hull of input contour. Otherwise input vector should be convex. + */ +CV_EXPORTS_W void approxPolyN(InputArray curve, OutputArray approxCurve, + int nsides, float epsilon_percentage = -1.0, + bool ensure_convex = true); + /** @brief Calculates a contour perimeter or a curve length. The function computes a curve length or a closed contour perimeter. diff --git a/modules/imgproc/src/approx.cpp b/modules/imgproc/src/approx.cpp index f05a6bcf3c..6c89f1d879 100644 --- a/modules/imgproc/src/approx.cpp +++ b/modules/imgproc/src/approx.cpp @@ -39,6 +39,7 @@ // //M*/ #include "precomp.hpp" +#include /****************************************************************************************\ * Chain Approximation * @@ -860,4 +861,240 @@ cvApproxPoly( const void* array, int header_size, return dst_seq; } +enum class PointStatus : int8_t +{ + REMOVED = -1, + RECALCULATE = 0, + CALCULATED = 1 +}; + +struct neighbours +{ + PointStatus pointStatus; + cv::Point2f point; + int next; + int prev; + + explicit neighbours(int next_ = -1, int prev_ = -1, const cv::Point2f& point_ = { -1, -1 }) + { + next = next_; + prev = prev_; + point = point_; + pointStatus = PointStatus::CALCULATED; + } +}; + +struct changes +{ + float area; + int vertex; + cv::Point2f intersection; + + explicit changes(float area_, int vertex_, const cv::Point2f& intersection_) + { + area = area_; + vertex = vertex_; + intersection = intersection_; + } + + bool operator < (const changes& elem) const + { + return (area < elem.area) || ((area == elem.area) && (vertex < elem.vertex)); + } + bool operator > (const changes& elem) const + { + return (area > elem.area) || ((area == elem.area) && (vertex > elem.vertex)); + } +}; + +/* + returns intersection point and extra area +*/ +static void recalculation(std::vector& hull, int vertex_id, float& area_, float& x, float& y) +{ + cv::Point2f vertex = hull[vertex_id].point, + next_vertex = hull[hull[vertex_id].next].point, + extra_vertex_1 = hull[hull[vertex_id].prev].point, + extra_vertex_2 = hull[hull[hull[vertex_id].next].next].point; + + cv::Point2f curr_edge = next_vertex - vertex, + prev_edge = vertex - extra_vertex_1, + next_edge = extra_vertex_2 - next_vertex; + + float cross = prev_edge.x * next_edge.y - prev_edge.y * next_edge.x; + if (abs(cross) < 1e-8) + { + area_ = FLT_MAX; + x = -1; + y = -1; + return; + } + + float t = (curr_edge.x * next_edge.y - curr_edge.y * next_edge.x) / cross; + cv::Point2f intersection = vertex + cv::Point2f(prev_edge.x * t, prev_edge.y * t); + + float area = 0.5f * abs((next_vertex.x - vertex.x) * (intersection.y - vertex.y) + - (intersection.x - vertex.x) * (next_vertex.y - vertex.y)); + + area_ = area; + x = intersection.x; + y = intersection.y; +} + +static void update(std::vector& hull, int vertex_id) +{ + neighbours& v1 = hull[vertex_id], & removed = hull[v1.next], & v2 = hull[removed.next]; + + removed.pointStatus = PointStatus::REMOVED; + v1.pointStatus = PointStatus::RECALCULATE; + v2.pointStatus = PointStatus::RECALCULATE; + hull[v1.prev].pointStatus = PointStatus::RECALCULATE; + v1.next = removed.next; + v2.prev = removed.prev; +} + +/* + A greedy algorithm based on contraction of vertices for approximating a convex contour by a bounding polygon +*/ +void cv::approxPolyN(InputArray _curve, OutputArray _approxCurve, + int nsides, float epsilon_percentage, bool ensure_convex) +{ + CV_INSTRUMENT_REGION(); + + CV_Assert(epsilon_percentage > 0 || epsilon_percentage == -1); + CV_Assert(nsides > 2); + + if (_approxCurve.fixedType()) + { + CV_Assert(_approxCurve.type() == CV_32FC2 || _approxCurve.type() == CV_32SC2); + } + + Mat curve; + int depth = _curve.depth(); + + CV_Assert(depth == CV_32F || depth == CV_32S); + + if (ensure_convex) + { + cv::convexHull(_curve, curve); + } + else + { + CV_Assert(isContourConvex(_curve)); + curve = _curve.getMat(); + } + + CV_Assert((curve.cols == 1 && curve.rows >= nsides) + || (curve.rows == 1 && curve.cols >= nsides)); + + if (curve.rows == 1) + { + curve = curve.reshape(0, curve.cols); + } + + std::vector hull(curve.rows); + int size = curve.rows; + std::priority_queue, std::greater> areas; + float extra_area = 0, max_extra_area = epsilon_percentage * static_cast(contourArea(_curve)); + + if (curve.depth() == CV_32S) + { + for (int i = 0; i < size; ++i) + { + Point t = curve.at(i, 0); + hull[i] = neighbours(i + 1, i - 1, Point2f(static_cast(t.x), static_cast(t.y))); + } + } + else + { + for (int i = 0; i < size; ++i) + { + Point2f t = curve.at(i, 0); + hull[i] = neighbours(i + 1, i - 1, t); + } + } + hull[0].prev = size - 1; + hull[size - 1].next = 0; + + if (size > nsides) + { + for (int vertex_id = 0; vertex_id < size; ++vertex_id) + { + float area, new_x, new_y; + recalculation(hull, vertex_id, area, new_x, new_y); + + areas.push(changes(area, vertex_id, Point2f(new_x, new_y))); + } + } + + while (size > nsides) + { + changes base = areas.top(); + int vertex_id = base.vertex; + + if (hull[vertex_id].pointStatus == PointStatus::REMOVED) + { + areas.pop(); + } + else if (hull[vertex_id].pointStatus == PointStatus::RECALCULATE) + { + float area, new_x, new_y; + areas.pop(); + recalculation(hull, vertex_id, area, new_x, new_y); + + areas.push(changes(area, vertex_id, Point2f(new_x, new_y))); + hull[vertex_id].pointStatus = PointStatus::CALCULATED; + } + else + { + if (epsilon_percentage != -1) + { + extra_area += base.area; + if (extra_area > max_extra_area) + { + break; + } + } + + size--; + hull[vertex_id].point = base.intersection; + update(hull, vertex_id); + } + } + + if (_approxCurve.fixedType()) + { + depth = _approxCurve.depth(); + } + _approxCurve.create(1, size, CV_MAKETYPE(depth, 2)); + Mat buf = _approxCurve.getMat(); + int last_free = 0; + + if (depth == CV_32S) + { + for (int i = 0; i < curve.rows; ++i) + { + if (hull[i].pointStatus != PointStatus::REMOVED) + { + Point t = Point(static_cast(round(hull[i].point.x)), + static_cast(round(hull[i].point.y))); + + buf.at(0, last_free) = t; + last_free++; + } + } + } + else + { + for (int i = 0; i < curve.rows; ++i) + { + if (hull[i].pointStatus != PointStatus::REMOVED) + { + buf.at(0, last_free) = hull[i].point; + last_free++; + } + } + } +} + /* End of file. */ diff --git a/modules/imgproc/test/test_approxpoly.cpp b/modules/imgproc/test/test_approxpoly.cpp index f09475c9fc..2b07c1a7b5 100644 --- a/modules/imgproc/test/test_approxpoly.cpp +++ b/modules/imgproc/test/test_approxpoly.cpp @@ -377,4 +377,79 @@ TEST(Imgproc_ApproxPoly, bad_epsilon) ASSERT_ANY_THROW(approxPolyDP(inputPoints, outputPoints, eps, false)); } +struct ApproxPolyN: public testing::Test +{ + void SetUp() + { + vector> inputPoints = { + { {87, 103}, {100, 112}, {96, 138}, {80, 169}, {60, 183}, {38, 176}, {41, 145}, {56, 118}, {76, 104} }, + { {196, 102}, {205, 118}, {174, 196}, {152, 207}, {102, 194}, {100, 175}, {131, 109} }, + { {372, 101}, {377, 119}, {337, 238}, {324, 248}, {240, 229}, {199, 214}, {232, 123}, {245, 103} }, + { {463, 86}, {563, 112}, {574, 135}, {596, 221}, {518, 298}, {412, 266}, {385, 164}, {462, 86} } + }; + + Mat image(600, 600, CV_8UC1, Scalar(0)); + + for (vector& polygon : inputPoints) { + polylines(image, { polygon }, true, Scalar(255), 1); + } + + findContours(image, contours, RETR_LIST, CHAIN_APPROX_NONE); + } + + vector> contours; +}; + +TEST_F(ApproxPolyN, accuracyInt) +{ + vector> rightCorners = { + { {72, 187}, {37, 176}, {42, 127}, {133, 64} }, + { {168, 212}, {92, 192}, {131, 109}, {213, 100} }, + { {72, 187}, {37, 176}, {42, 127}, {133, 64} }, + { {384, 100}, {333, 251}, {197, 220}, {239, 103} }, + { {168, 212}, {92, 192}, {131, 109}, {213, 100} }, + { {333, 251}, {197, 220}, {239, 103}, {384, 100} }, + { {542, 6}, {596, 221}, {518, 299}, {312, 236} }, + { {596, 221}, {518, 299}, {312, 236}, {542, 6} } + }; + EXPECT_EQ(rightCorners.size(), contours.size()); + + for (size_t i = 0; i < contours.size(); ++i) { + std::vector corners; + approxPolyN(contours[i], corners, 4, -1, true); + ASSERT_EQ(rightCorners[i], corners ); + } +} + +TEST_F(ApproxPolyN, accuracyFloat) +{ + vector> rightCorners = { + { {72.f, 187.f}, {37.f, 176.f}, {42.f, 127.f}, {133.f, 64.f} }, + { {168.f, 212.f}, {92.f, 192.f}, {131.f, 109.f}, {213.f, 100.f} }, + { {72.f, 187.f}, {37.f, 176.f}, {42.f, 127.f}, {133.f, 64.f} }, + { {384.f, 100.f}, {333.f, 251.f}, {197.f, 220.f}, {239.f, 103.f} }, + { {168.f, 212.f}, {92.f, 192.f}, {131.f, 109.f}, {213.f, 100.f} }, + { {333.f, 251.f}, {197.f, 220.f}, {239.f, 103.f}, {384.f, 100.f} }, + { {542.f, 6.f}, {596.f, 221.f}, {518.f, 299.f}, {312.f, 236.f} }, + { {596.f, 221.f}, {518.f, 299.f}, {312.f, 236.f}, {542.f, 6.f} } + }; + EXPECT_EQ(rightCorners.size(), contours.size()); + + for (size_t i = 0; i < contours.size(); ++i) { + std::vector corners; + approxPolyN(contours[i], corners, 4, -1, true); + EXPECT_LT(cvtest::norm(rightCorners[i], corners, NORM_INF), .5f); + } +} + +TEST_F(ApproxPolyN, bad_args) +{ + Mat contour(10, 1, CV_32FC2); + vector> bad_contours; + vector corners; + ASSERT_ANY_THROW(approxPolyN(contour, corners, 0)); + ASSERT_ANY_THROW(approxPolyN(contour, corners, 3, 0)); + ASSERT_ANY_THROW(approxPolyN(bad_contours, corners, 4)); +} + }} // namespace