Merge pull request #25146 from mshabunin:cpp-contours

Reworked findContours to reduce C-API usage #25146

What is done:
* rewritten `findContours` and `icvApproximateChainTC89` using C++ data structures
* extracted LINK_RUNS mode to separate new public functions - `findContoursLinkRuns` (it uses completely different algorithm)
* ~added new public `cv::approximateChainTC89`~ - **❌ decided to hide it**
* enabled chain code output (method = 0, no public enum value for this in C++ yet)
* kept old function as `findContours_old` (exported, but not exposed to user)
* added more tests for findContours (`test_contours_new.cpp`), some tests compare results of old function with new one. Following tests have been added:
  * contours of random rectangle
  * contours of many small (1-2px) blobs
  * contours of random noise
  * backport of old accuracy test
  * separate test for LINK RUNS variant

What is left to be done (can be done now or later):
* improve tests: 
  * some tests have limited verification (e.g. only verify contour sizes)
  * perhaps reference data can be collected and stored
  * maybe more test variants can be added (?)
* add enum value for chain code output and a method of returning starting points (e.g. first 8 elements of returned `vector<uchar>` can represent 2 int point coordinates)
* add documentation for new functions - **✔️ DONE**
* check and improve performance (my experiment showed 0.7x-1.1x some time ago)
* remove old functions completely (?)
* change contour return order (BFS) or allow to select it (?)
* return result tree as-is (?) (new data structures should be exposed, bindings should adapt)

modules/imgproc/include/opencv2/imgproc.hpp

@@ -4040,6 +4040,19 @@ A program using pyramid scaling, Canny, contours and contour simplification to f
 squares in the input image.
 */
 
+//! @brief Find contours using link runs algorithm
+//!
+//! This function implements an algorithm different from cv::findContours:
+//! - doesn't allocate temporary image internally, thus it has reduced memory consumption
+//! - supports CV_8UC1 images only
+//! - outputs 2-level hierarhy only (RETR_CCOMP mode)
+//! - doesn't support approximation change other than CHAIN_APPROX_SIMPLE
+//! In all other aspects this function is compatible with cv::findContours.
+CV_EXPORTS_W void findContoursLinkRuns(InputArray image, OutputArrayOfArrays contours, OutputArray hierarchy);
+
+//! @overload
+CV_EXPORTS_W void findContoursLinkRuns(InputArray image, OutputArrayOfArrays contours);
+
 /** @brief Approximates a polygonal curve(s) with the specified precision.
 
 The function cv::approxPolyDP approximates a curve or a polygon with another curve/polygon with less

modules/imgproc/include/opencv2/imgproc/detail/legacy.hpp

@@ -0,0 +1,30 @@
+// This file is part of OpenCV project.
+// It is subject to the license terms in the LICENSE file found in the top-level directory
+// of this distribution and at http://opencv.org/license.html
+
+#ifndef OPENCV_IMGPROC_DETAIL_LEGACY_HPP
+#define OPENCV_IMGPROC_DETAIL_LEGACY_HPP
+
+#include "opencv2/imgproc.hpp"
+
+namespace cv {
+
+#ifdef __OPENCV_BUILD
+
+CV_EXPORTS void findContours_legacy(InputArray _image,
+                                    OutputArrayOfArrays _contours,
+                                    OutputArray _hierarchy,
+                                    int mode,
+                                    int method,
+                                    Point offset = Point());
+CV_EXPORTS void findContours_legacy(InputArray image,
+                                    OutputArrayOfArrays contours,
+                                    int mode,
+                                    int method,
+                                    Point offset = Point());
+
+#endif
+
+}  // namespace cv
+
+#endif  // OPENCV_IMGPROC_DETAIL_LEGACY_HPP

modules/imgproc/src/contours.cpp

@@ -40,6 +40,7 @@
 //M*/
 #include "precomp.hpp"
 #include "opencv2/core/hal/intrin.hpp"
+#include "opencv2/imgproc/detail/legacy.hpp"
 
 using namespace cv;
 
@@ -1813,7 +1814,7 @@ cvFindContours( void*  img,  CvMemStorage*  storage,
     return cvFindContours_Impl(img, storage, firstContour, cntHeaderSize, mode, method, offset, 1);
 }
 
-void cv::findContours( InputArray _image, OutputArrayOfArrays _contours,
+void cv::findContours_legacy( InputArray _image, OutputArrayOfArrays _contours,
                    OutputArray _hierarchy, int mode, int method, Point offset )
 {
     CV_INSTRUMENT_REGION();
@@ -1878,7 +1879,7 @@ void cv::findContours( InputArray _image, OutputArrayOfArrays _contours,
     }
 }
 
-void cv::findContours( InputArray _image, OutputArrayOfArrays _contours,
+void cv::findContours_legacy( InputArray _image, OutputArrayOfArrays _contours,
                        int mode, int method, Point offset)
 {
     CV_INSTRUMENT_REGION();

modules/imgproc/src/contours_approx.cpp

@@ -0,0 +1,354 @@
+// This file is part of OpenCV project.
+// It is subject to the license terms in the LICENSE file found in the top-level directory
+// of this distribution and at http://opencv.org/license.html
+
+#include "opencv2/core/base.hpp"
+#include "opencv2/core/types.hpp"
+#include "opencv2/imgproc.hpp"
+#include "contours_common.hpp"
+#include <vector>
+
+using namespace std;
+using namespace cv;
+
+namespace {
+
+struct ApproxItem
+{
+    Point pt;
+    size_t k;  // support region
+    int s;  // 1-curvature
+    bool removed;
+    ApproxItem() : k(0), s(0), removed(false) {}
+    ApproxItem(const Point& pt_, int s_) : pt(pt_), k(0), s(s_), removed(false) {}
+};
+
+static const schar abs_diff[16] = {1, 2, 3, 4, 3, 2, 1, 0, 1, 2, 3, 4, 3, 2, 1};
+static const Point chainCodeDeltas[8] =
+    {{1, 0}, {1, -1}, {0, -1}, {-1, -1}, {-1, 0}, {-1, 1}, {0, 1}, {1, 1}};
+
+// Pass 0.
+// Restores all the digital curve points from the chain code.
+// Removes the points (from the resultant polygon)
+// that have zero 1-curvature
+static vector<ApproxItem> pass_0(const vector<schar>& chain, Point pt, bool isApprox, bool isFull)
+{
+    vector<ApproxItem> res;
+    const size_t len = chain.size();
+    res.reserve(len / 2);
+    for (size_t i = 0; i < len; ++i)
+    {
+        const schar prev = (i == 0) ? chain[len - 1] : chain[i - 1];
+        const schar cur = chain[i];
+        const schar s = abs_diff[cur - prev + 7];
+        if ((!isApprox && (isFull || s != 0)) || isApprox)
+        {
+            res.push_back(ApproxItem(pt, s));
+            if (s == 0)
+                (res.end() - 1)->removed = true;
+        }
+        pt += chainCodeDeltas[cur];
+    }
+    return res;
+}
+
+static vector<Point> gatherPoints(const vector<ApproxItem>& ares)
+{
+    vector<Point> res;
+    res.reserve(ares.size() / 2);
+    for (const ApproxItem& item : ares)
+    {
+        if (item.removed)
+            continue;
+        res.push_back(item.pt);
+    }
+    return res;
+}
+
+static size_t calc_support(const vector<ApproxItem>& ares, size_t i)
+{
+    const size_t len = ares.size();
+    /* determine support region */
+    int d_num = 0;
+    int l = 0;
+    size_t k = 1;
+    for (;; k++)
+    {
+        CV_Assert(k <= len);
+        /* calc indices */
+        const size_t i1 = (i >= k) ? (i - k) : (len - k + i);
+        const size_t i2 = (i + k < len) ? (i + k) : (i + k - len);
+
+        const int dx = ares[i2].pt.x - ares[i1].pt.x;
+        const int dy = ares[i2].pt.y - ares[i1].pt.y;
+
+        /* distance between p_(i - k) and p_(i + k) */
+        const int lk = dx * dx + dy * dy;
+
+        /* distance between p_i and the line (p_(i-k), p_(i+k)) */
+        const int dk_num =
+            (ares[i].pt.x - ares[i1].pt.x) * dy - (ares[i].pt.y - ares[i1].pt.y) * dx;
+
+        union
+        {
+            int i;
+            float f;
+        } d;
+        d.f = (float)(((double)d_num) * lk - ((double)dk_num) * l);
+
+        if (k > 1 && (l >= lk || ((d_num > 0 && d.i <= 0) || (d_num < 0 && d.i >= 0))))
+            break;
+
+        d_num = dk_num;
+        l = lk;
+    }
+    return k - 1;
+}
+
+static int calc_cosine(const vector<ApproxItem>& ares, size_t i)
+{
+    const size_t k = ares[i].k;
+    size_t j;
+    int s;
+    const size_t len = ares.size();
+    /* calc k-cosine curvature */
+    for (j = k, s = 0; j > 0; j--)
+    {
+        const size_t i1 = (i >= j) ? (i - j) : (len - j + i);
+        const size_t i2 = (i + j < len) ? (i + j) : (i + j - len);
+
+        const int dx1 = ares[i1].pt.x - ares[i].pt.x;
+        const int dy1 = ares[i1].pt.y - ares[i].pt.y;
+        const int dx2 = ares[i2].pt.x - ares[i].pt.x;
+        const int dy2 = ares[i2].pt.y - ares[i].pt.y;
+
+        if ((dx1 | dy1) == 0 || (dx2 | dy2) == 0)
+            break;
+
+        double temp_num = dx1 * dx2 + dy1 * dy2;
+        temp_num = (float)(temp_num / sqrt(((double)dx1 * dx1 + (double)dy1 * dy1) *
+                                           ((double)dx2 * dx2 + (double)dy2 * dy2)));
+        Cv32suf sk;
+        sk.f = (float)(temp_num + 1.1);
+
+        CV_Assert(0 <= sk.f && sk.f <= 2.2);
+        if (j < k && sk.i <= s)
+            break;
+
+        s = sk.i;
+    }
+    return s;
+}
+
+static bool calc_nms_cleanup(const vector<ApproxItem>& ares, size_t i)
+{
+    const size_t k2 = ares[i].k >> 1;
+    const int s = ares[i].s;
+    const size_t len = ares.size();
+    size_t j;
+    for (j = 1; j <= k2; j++)
+    {
+        const size_t i1 = (i >= j) ? (i - j) : (len - j + i);
+        const size_t i2 = (i + j < len) ? (i + j) : (i + j - len);
+        if (ares[i1].s > s || ares[i2].s > s)
+            break;
+    }
+    return j <= k2;
+}
+
+static bool calc_dominance(const vector<ApproxItem>& ares, size_t i)
+{
+    const size_t len = ares.size();
+    CV_Assert(len > 0);
+    const size_t i1 = (i >= 1) ? (i - 1) : (len - 1 + i);
+    const size_t i2 = (i + 1 < len) ? (i + 1) : (i + 1 - len);
+    return ares[i].s <= ares[i1].s || ares[i].s <= ares[i2].s;
+}
+
+inline size_t get_next_idx(const vector<ApproxItem>& ares, const size_t start)
+{
+    const size_t len = ares.size();
+    size_t res = start + 1;
+    for (; res < len; ++res)
+    {
+        if (!ares[res].removed)
+            break;
+    }
+    return res;
+}
+
+inline void clear_until(vector<ApproxItem>& ares, const size_t start, const size_t finish)
+{
+    const size_t len = ares.size();
+    for (size_t i = start + 1; i < finish && i < len; ++i)
+    {
+        ares[i].removed = true;
+    }
+}
+
+static bool calc_new_start(vector<ApproxItem>& ares, size_t& res)
+{
+    const size_t len = ares.size();
+    CV_Assert(len > 0);
+    size_t i1;
+    // remove all previous items from the beginning
+    for (i1 = 1; i1 < len && ares[i1].s != 0; i1++)
+    {
+        ares[i1 - 1].s = 0;
+    }
+    if (i1 == len)
+    {
+        // all points survived - skip to the end
+        return false;
+    }
+    i1--;
+
+    size_t i2;
+    // remove all following items from the end
+    for (i2 = len - 2; i2 > 0 && ares[i2].s != 0; i2--)
+    {
+        clear_until(ares, i2, len);
+        ares[i2 + 1].s = 0;
+    }
+    i2++;
+
+    // only two points left
+    if (i1 == 0 && i2 == len - 1)
+    {
+        // find first non-removed element from the start
+        i1 = get_next_idx(ares, 0);
+        // append first item to the end
+        ares.push_back(ares[0]);
+        (ares.end() - 1)->removed = false;
+    }
+    res = i1;
+    return true;
+}
+
+static void pass_cleanup(vector<ApproxItem>& ares, size_t start_idx)
+{
+    int count = 1;
+
+    const size_t len = ares.size();
+    size_t first = start_idx;
+    for (size_t i = start_idx, prev = i; i < len; ++i)
+    {
+        ApproxItem& item = ares[i];
+        if (item.removed)
+            continue;
+        size_t next_idx = get_next_idx(ares, i);
+        if (next_idx == len || next_idx - i != 1)
+        {
+            if (count >= 2)
+            {
+                if (count == 2)
+                {
+                    const int s1 = ares[prev].s;
+                    const int s2 = ares[i].s;
+
+                    if (s1 > s2 || (s1 == s2 && ares[prev].k <= ares[i].k))
+                        /* remove second */
+                        clear_until(ares, get_next_idx(ares, prev), get_next_idx(ares, i));
+                    else
+                        /* remove first */
+                        clear_until(ares, first, i);
+                }
+                else
+                {
+                    first = get_next_idx(ares, first);
+                    clear_until(ares, first, i);
+                }
+            }
+            clear_until(ares, first, i);
+            first = i;
+            count = 1;
+        }
+        else
+        {
+            ++count;
+        }
+        prev = i;
+    }
+}
+
+}  // namespace
+
+
+vector<Point> cv::approximateChainTC89(vector<schar> chain, const Point& origin, const int method)
+{
+    if (chain.size() == 0)
+    {
+        return vector<Point>({origin});
+    }
+
+    const bool isApprox = method == CHAIN_APPROX_TC89_L1 || method == CHAIN_APPROX_TC89_KCOS;
+
+    ApproxItem root;
+    vector<ApproxItem> ares = pass_0(chain, origin, isApprox, method == CHAIN_APPROX_NONE);
+
+    if (isApprox)
+    {
+        CV_DbgAssert(ares.size() < (size_t)numeric_limits<int>::max());
+
+        // Pass 1.
+        // Determines support region for all the remained points */
+        for (size_t i = 0; i < ares.size(); ++i)
+        {
+            ApproxItem& item = ares[i];
+            if (item.removed)
+                continue;
+            item.k = calc_support(ares, i);
+
+            if (method == CHAIN_APPROX_TC89_KCOS)
+                item.s = calc_cosine(ares, i);
+        }
+
+        // Pass 2.
+        // Performs non-maxima suppression
+        for (size_t i = 0; i < ares.size(); ++i)
+        {
+            ApproxItem& item = ares[i];
+            if (calc_nms_cleanup(ares, i))
+            {
+                item.s = 0;  // "clear"
+                item.removed = true;
+            }
+        }
+
+        // Pass 3.
+        // Removes non-dominant points with 1-length support region */
+        for (size_t i = 0; i < ares.size(); ++i)
+        {
+            ApproxItem& item = ares[i];
+            if (item.removed)
+                continue;
+            if (item.k == 1 && calc_dominance(ares, i))
+            {
+                item.s = 0;
+                item.removed = true;
+            }
+        }
+
+        if (method == cv::CHAIN_APPROX_TC89_L1)
+        {
+            // Pass 4.
+            // Cleans remained couples of points
+            bool skip = false;
+            size_t new_start_idx = 0;
+            const size_t len = ares.size();
+            if (ares[0].s != 0 && ares[len - 1].s != 0)
+            {
+                if (!calc_new_start(ares, new_start_idx))
+                {
+                    skip = true;
+                }
+            }
+            if (!skip)
+            {
+                pass_cleanup(ares, new_start_idx);
+            }
+        }
+    }
+
+    return gatherPoints(ares);
+}

modules/imgproc/src/contours_common.cpp

@@ -0,0 +1,75 @@
+// This file is part of OpenCV project.
+// It is subject to the license terms in the LICENSE file found in the top-level directory
+// of this distribution and at http://opencv.org/license.html
+
+#include "precomp.hpp"
+#include "contours_common.hpp"
+#include <map>
+#include <limits>
+
+using namespace std;
+using namespace cv;
+
+void cv::contourTreeToResults(CTree& tree,
+                              int res_type,
+                              OutputArrayOfArrays& _contours,
+                              OutputArray& _hierarchy)
+{
+    // check if there are no results
+    if (tree.isEmpty() || (tree.elem(0).body.isEmpty() && (tree.elem(0).first_child == -1)))
+    {
+        _contours.clear();
+        return;
+    }
+
+    // mapping for indexes (original -> resulting)
+    map<int, int> index_mapping;
+    index_mapping[-1] = -1;
+    index_mapping[0] = -1;
+
+    CV_Assert(tree.size() < (size_t)numeric_limits<int>::max());
+    const int total = (int)tree.size() - 1;
+    _contours.create(total, 1, 0, -1, true);
+    {
+        int i = 0;
+        CIterator it(tree);
+        while (!it.isDone())
+        {
+            const CNode& elem = it.getNext_s();
+            CV_Assert(elem.self() != -1);
+            if (elem.self() == 0)
+                continue;
+            index_mapping[elem.self()] = i;
+            CV_Assert(elem.body.size() < (size_t)numeric_limits<int>::max());
+            const int sz = (int)elem.body.size();
+            _contours.create(sz, 1, res_type, i, true);
+            if (sz > 0)
+            {
+                Mat cmat = _contours.getMat(i);
+                CV_Assert(cmat.isContinuous());
+                elem.body.copyTo(cmat.data);
+            }
+            ++i;
+        }
+    }
+
+    if (_hierarchy.needed())
+    {
+        _hierarchy.create(1, total, CV_32SC4, -1, true);
+        Mat h_mat = _hierarchy.getMat();
+        int i = 0;
+        CIterator it(tree);
+        while (!it.isDone())
+        {
+            const CNode& elem = it.getNext_s();
+            if (elem.self() == 0)
+                continue;
+            Vec4i& h_vec = h_mat.at<Vec4i>(i);
+            h_vec = Vec4i(index_mapping.at(elem.next),
+                          index_mapping.at(elem.prev),
+                          index_mapping.at(elem.first_child),
+                          index_mapping.at(elem.parent));
+            ++i;
+        }
+    }
+}

modules/imgproc/src/contours_common.hpp

@@ -0,0 +1,219 @@
+// This file is part of OpenCV project.
+// It is subject to the license terms in the LICENSE file found in the top-level directory
+// of this distribution and at http://opencv.org/license.html
+
+#ifndef OPENCV_CONTOURS_COMMON_HPP
+#define OPENCV_CONTOURS_COMMON_HPP
+
+#include "precomp.hpp"
+#include <stack>
+
+namespace cv {
+
+static const schar MAX_SIZE = 16;
+
+static const cv::Point chainCodeDeltas[8] =
+    {{1, 0}, {1, -1}, {0, -1}, {-1, -1}, {-1, 0}, {-1, 1}, {0, 1}, {1, 1}};
+
+static inline int getDelta(schar s, size_t step)
+{
+    CV_DbgAssert(s >= 0 && s < 16);
+    const cv::Point res = chainCodeDeltas[s % 8];
+    return res.x + res.y * (int)step;
+}
+
+inline schar clamp_direction(schar dir)
+{
+    return std::min(dir, (schar)15);
+}
+
+template <typename T>
+class TreeNode
+{
+private:
+    int self_;
+
+public:
+    // tree hierarchy (parent - children)
+    int parent;
+    int first_child;
+    // 1st linked list - bidirectional - sibling children
+    int prev;
+    int next;
+    // 2nd linked list - unidirectional - not related to 1st list
+    int ctable_next;
+    T body;
+
+public:
+    TreeNode(int self) :
+        self_(self), parent(-1), first_child(-1), prev(-1), next(-1), ctable_next(-1)
+    {
+        CV_Assert(self >= 0);
+    }
+    int self() const
+    {
+        return self_;
+    }
+};
+
+template <typename T>
+class Tree
+{
+private:
+    std::vector<TreeNode<T>> nodes;
+
+public:
+    TreeNode<T>& newElem()
+    {
+        const size_t idx = nodes.size();
+        CV_DbgAssert(idx < (size_t)std::numeric_limits<int>::max());
+        nodes.push_back(TreeNode<T>((int)idx));
+        return nodes[idx];
+    }
+    TreeNode<T>& elem(int idx)
+    {
+        CV_DbgAssert(idx >= 0 && (size_t)idx < nodes.size());
+        return nodes[(size_t)idx];
+    }
+    const TreeNode<T>& elem(int idx) const
+    {
+        CV_DbgAssert(idx >= 0 && (size_t)idx < nodes.size());
+        return nodes[(size_t)idx];
+    }
+    int lastSibling(int e) const
+    {
+        if (e != -1)
+        {
+            while (true)
+            {
+                const TreeNode<T>& cur_elem = elem(e);
+                if (cur_elem.next == -1)
+                    break;
+                e = cur_elem.next;
+            }
+        }
+        return e;
+    }
+    void addSiblingAfter(int prev, int idx)
+    {
+        TreeNode<T>& prev_item = nodes[prev];
+        TreeNode<T>& child = nodes[idx];
+        child.parent = prev_item.parent;
+        if (prev_item.next != -1)
+        {
+            nodes[prev_item.next].prev = idx;
+            child.next = prev_item.next;
+        }
+        child.prev = prev;
+        prev_item.next = idx;
+    }
+    void addChild(int parent_idx, int child_idx)
+    {
+        TreeNode<T>& parent = nodes[parent_idx];
+        TreeNode<T>& child = nodes[child_idx];
+        if (parent.first_child != -1)
+        {
+            TreeNode<T>& fchild_ = nodes[parent.first_child];
+            fchild_.prev = child_idx;
+            child.next = parent.first_child;
+        }
+        parent.first_child = child_idx;
+        child.parent = parent_idx;
+        child.prev = -1;
+    }
+    bool isEmpty() const
+    {
+        return nodes.size() == 0;
+    }
+    size_t size() const
+    {
+        return nodes.size();
+    }
+};
+
+template <typename T>
+class TreeIterator
+{
+public:
+    TreeIterator(Tree<T>& tree_) : tree(tree_)
+    {
+        CV_Assert(!tree.isEmpty());
+        levels.push(0);
+    }
+    bool isDone() const
+    {
+        return levels.empty();
+    }
+    const TreeNode<T>& getNext_s()
+    {
+        int idx = levels.top();
+        levels.pop();
+        const TreeNode<T>& res = tree.elem(idx);
+        int cur = tree.lastSibling(res.first_child);
+        while (cur != -1)
+        {
+            levels.push(cur);
+            cur = tree.elem(cur).prev;
+        }
+        return res;
+    }
+
+private:
+    std::stack<int> levels;
+    Tree<T>& tree;
+};
+
+//==============================================================================
+
+class Contour
+{
+public:
+    cv::Rect brect;
+    cv::Point origin;
+    std::vector<cv::Point> pts;
+    std::vector<schar> codes;
+    bool isHole;
+    bool isChain;
+
+    Contour() : isHole(false), isChain(false) {}
+    void updateBoundingRect() {}
+    bool isEmpty() const
+    {
+        return pts.size() == 0 && codes.size() == 0;
+    }
+    size_t size() const
+    {
+        return isChain ? codes.size() : pts.size();
+    }
+    void copyTo(void* data) const
+    {
+        // NOTE: Mat::copyTo doesn't work because it creates new Mat object
+        //       instead of reusing existing vector data
+        if (isChain)
+        {
+            memcpy(data, &codes[0], codes.size() * sizeof(codes[0]));
+        }
+        else
+        {
+            memcpy(data, &pts[0], pts.size() * sizeof(pts[0]));
+        }
+    }
+};
+
+typedef TreeNode<Contour> CNode;
+typedef Tree<Contour> CTree;
+typedef TreeIterator<Contour> CIterator;
+
+
+void contourTreeToResults(CTree& tree,
+                          int res_type,
+                          cv::OutputArrayOfArrays& _contours,
+                          cv::OutputArray& _hierarchy);
+
+
+std::vector<Point>
+    approximateChainTC89(std::vector<schar> chain, const Point& origin, const int method);
+
+}  // namespace cv
+
+#endif  // OPENCV_CONTOURS_COMMON_HPP

modules/imgproc/src/contours_link.cpp

@@ -0,0 +1,417 @@
+// This file is part of OpenCV project.
+// It is subject to the license terms in the LICENSE file found in the top-level directory
+// of this distribution and at http://opencv.org/license.html
+
+#include "precomp.hpp"
+#include "contours_common.hpp"
+#include "opencv2/core/hal/intrin.hpp"
+
+using namespace cv;
+using namespace std;
+
+//==============================================================================
+
+namespace {
+
+inline static int findStartContourPoint(uchar* src_data, Size img_size, int j)
+{
+#if (CV_SIMD || CV_SIMD_SCALABLE)
+    v_uint8 v_zero = vx_setzero_u8();
+    for (; j <= img_size.width - VTraits<v_uint8>::vlanes(); j += VTraits<v_uint8>::vlanes())
+    {
+        v_uint8 vmask = (v_ne(vx_load((uchar*)(src_data + j)), v_zero));
+        if (v_check_any(vmask))
+        {
+            j += v_scan_forward(vmask);
+            return j;
+        }
+    }
+#endif
+    for (; j < img_size.width && !src_data[j]; ++j)
+        ;
+    return j;
+}
+
+inline static int findEndContourPoint(uchar* src_data, Size img_size, int j)
+{
+#if (CV_SIMD || CV_SIMD_SCALABLE)
+    if (j < img_size.width && !src_data[j])
+    {
+        return j;
+    }
+    else
+    {
+        v_uint8 v_zero = vx_setzero_u8();
+        for (; j <= img_size.width - VTraits<v_uint8>::vlanes(); j += VTraits<v_uint8>::vlanes())
+        {
+            v_uint8 vmask = (v_eq(vx_load((uchar*)(src_data + j)), v_zero));
+            if (v_check_any(vmask))
+            {
+                j += v_scan_forward(vmask);
+                return j;
+            }
+        }
+    }
+#endif
+    for (; j < img_size.width && src_data[j]; ++j)
+        ;
+
+    return j;
+}
+
+//==============================================================================
+
+struct LinkRunPoint
+{
+    int link;
+    int next;
+    Point pt;
+    LinkRunPoint() : link(-1), next(-1) {}
+    LinkRunPoint(const Point& pt_) : link(-1), next(-1), pt(pt_) {}
+};
+
+typedef LinkRunPoint LRP;
+
+//==============================================================================
+
+class LinkRunner
+{
+public:
+    enum LinkConnectionDirection
+    {
+        ICV_SINGLE = 0,
+        ICV_CONNECTING_ABOVE = 1,
+        ICV_CONNECTING_BELOW = -1,
+    };
+
+    CTree tree;
+
+    vector<LRP> rns;
+    vector<int> ext_rns;
+    vector<int> int_rns;
+
+public:
+    LinkRunner()
+    {
+        tree.newElem();
+        rns.reserve(100);
+    }
+    void process(Mat& image);
+    void convertLinks(int& first, int& prev, bool isHole);
+    void establishLinks(int& prev_point,
+                        int upper_run,
+                        int lower_run,
+                        const int upper_total,
+                        const int lower_total);
+};
+
+void LinkRunner::convertLinks(int& first, int& prev, bool isHole)
+{
+    const vector<int>& contours = isHole ? int_rns : ext_rns;
+    int count = 0;
+    for (int j = 0; j < (int)contours.size(); j++, count++)
+    {
+        int start = contours[j];
+        int cur = start;
+
+        if (rns[cur].link == -1)
+            continue;
+
+        CNode& node = tree.newElem();
+        node.body.isHole = isHole;
+
+        do
+        {
+            node.body.pts.push_back(rns[cur].pt);
+            int p_temp = cur;
+            cur = rns[cur].link;
+            rns[p_temp].link = -1;
+        }
+        while (cur != start);
+
+        if (first == 0)
+        {
+            tree.addChild(0, node.self());
+            prev = first = node.self();
+        }
+        else
+        {
+            tree.addSiblingAfter(prev, node.self());
+            prev = node.self();
+        }
+    }
+}
+void LinkRunner::establishLinks(int& prev_point,
+                                int upper_run,
+                                int lower_run,
+                                const int upper_total,
+                                const int lower_total)
+{
+    int k, n;
+    int connect_flag = ICV_SINGLE;
+    for (k = 0, n = 0; k < upper_total / 2 && n < lower_total / 2;)
+    {
+        switch (connect_flag)
+        {
+            case ICV_SINGLE:
+                if (rns[rns[upper_run].next].pt.x < rns[rns[lower_run].next].pt.x)
+                {
+                    if (rns[rns[upper_run].next].pt.x >= rns[lower_run].pt.x - 1)
+                    {
+                        rns[lower_run].link = upper_run;
+                        connect_flag = ICV_CONNECTING_ABOVE;
+                        prev_point = rns[upper_run].next;
+                    }
+                    else
+                        rns[rns[upper_run].next].link = upper_run;
+                    k++;
+                    upper_run = rns[rns[upper_run].next].next;
+                }
+                else
+                {
+                    if (rns[upper_run].pt.x <= rns[rns[lower_run].next].pt.x + 1)
+                    {
+                        rns[lower_run].link = upper_run;
+                        connect_flag = ICV_CONNECTING_BELOW;
+                        prev_point = rns[lower_run].next;
+                    }
+                    else
+                    {
+                        rns[lower_run].link = rns[lower_run].next;
+                        // First point of contour
+                        ext_rns.push_back(lower_run);
+                    }
+                    n++;
+                    lower_run = rns[rns[lower_run].next].next;
+                }
+                break;
+            case ICV_CONNECTING_ABOVE:
+                if (rns[upper_run].pt.x > rns[rns[lower_run].next].pt.x + 1)
+                {
+                    rns[prev_point].link = rns[lower_run].next;
+                    connect_flag = ICV_SINGLE;
+                    n++;
+                    lower_run = rns[rns[lower_run].next].next;
+                }
+                else
+                {
+                    rns[prev_point].link = upper_run;
+                    if (rns[rns[upper_run].next].pt.x < rns[rns[lower_run].next].pt.x)
+                    {
+                        k++;
+                        prev_point = rns[upper_run].next;
+                        upper_run = rns[rns[upper_run].next].next;
+                    }
+                    else
+                    {
+                        connect_flag = ICV_CONNECTING_BELOW;
+                        prev_point = rns[lower_run].next;
+                        n++;
+                        lower_run = rns[rns[lower_run].next].next;
+                    }
+                }
+                break;
+            case ICV_CONNECTING_BELOW:
+                if (rns[lower_run].pt.x > rns[rns[upper_run].next].pt.x + 1)
+                {
+                    rns[rns[upper_run].next].link = prev_point;
+                    connect_flag = ICV_SINGLE;
+                    k++;
+                    upper_run = rns[rns[upper_run].next].next;
+                }
+                else
+                {
+                    // First point of contour
+                    int_rns.push_back(lower_run);
+
+                    rns[lower_run].link = prev_point;
+                    if (rns[rns[lower_run].next].pt.x < rns[rns[upper_run].next].pt.x)
+                    {
+                        n++;
+                        prev_point = rns[lower_run].next;
+                        lower_run = rns[rns[lower_run].next].next;
+                    }
+                    else
+                    {
+                        connect_flag = ICV_CONNECTING_ABOVE;
+                        k++;
+                        prev_point = rns[upper_run].next;
+                        upper_run = rns[rns[upper_run].next].next;
+                    }
+                }
+                break;
+        }
+    }  // k, n
+
+    for (; n < lower_total / 2; n++)
+    {
+        if (connect_flag != ICV_SINGLE)
+        {
+            rns[prev_point].link = rns[lower_run].next;
+            connect_flag = ICV_SINGLE;
+            lower_run = rns[rns[lower_run].next].next;
+            continue;
+        }
+        rns[rns[lower_run].next] = rns[rns[lower_run].next];
+        rns[lower_run].link = rns[lower_run].next;
+
+        // First point of contour
+        ext_rns.push_back(lower_run);
+        lower_run = rns[rns[lower_run].next].next;
+    }
+
+    for (; k < upper_total / 2; k++)
+    {
+        if (connect_flag != ICV_SINGLE)
+        {
+            rns[rns[upper_run].next].link = prev_point;
+            connect_flag = ICV_SINGLE;
+            upper_run = rns[rns[upper_run].next].next;
+            continue;
+        }
+        rns[rns[upper_run].next] = rns[rns[upper_run].next];
+        rns[rns[upper_run].next].link = upper_run;
+        upper_run = rns[rns[upper_run].next].next;
+    }
+}
+
+
+void LinkRunner::process(Mat& image)
+{
+    const Size sz = image.size();
+    int j;
+    int lower_total;
+    int upper_total;
+    int all_total;
+
+    Point cur_point;
+
+    rns.reserve(sz.height);  // optimization, assuming some contours exist
+
+    // First line. None of runs is binded
+    rns.push_back(LRP());
+    int upper_line = (int)rns.size() - 1;
+    int cur = upper_line;
+    for (j = 0; j < sz.width;)
+    {
+        j = findStartContourPoint(image.ptr<uchar>(), sz, j);
+
+        if (j == sz.width)
+            break;
+
+        cur_point.x = j;
+
+        rns.push_back(LRP(cur_point));
+        rns[cur].next = (int)rns.size() - 1;
+        cur = rns[cur].next;
+
+        j = findEndContourPoint(image.ptr<uchar>(), sz, j + 1);
+
+        cur_point.x = j - 1;
+
+        rns.push_back(LRP(cur_point));
+        rns[cur].next = (int)rns.size() - 1;
+        rns[cur].link = rns[cur].next;
+
+        // First point of contour
+        ext_rns.push_back(cur);
+        cur = rns[cur].next;
+    }
+    upper_line = rns[upper_line].next;
+    upper_total = (int)rns.size() - 1;  // runs->total - 1;
+
+    int last_elem = cur;
+    rns[cur].next = -1;
+    int prev_point = -1;
+    int lower_line = -1;
+    for (int i = 1; i < sz.height; i++)
+    {
+        // Find runs in next line
+        cur_point.y = i;
+        all_total = (int)rns.size();  // runs->total;
+        for (j = 0; j < sz.width;)
+        {
+            j = findStartContourPoint(image.ptr<uchar>(i), sz, j);
+
+            if (j == sz.width)
+                break;
+
+            cur_point.x = j;
+
+            rns.push_back(LRP(cur_point));
+            rns[cur].next = (int)rns.size() - 1;
+            cur = rns[cur].next;
+
+            j = findEndContourPoint(image.ptr<uchar>(i), sz, j + 1);
+
+            cur_point.x = j - 1;
+            rns.push_back(LRP(cur_point));
+            cur = rns[cur].next = (int)rns.size() - 1;
+        }  // j
+        lower_line = rns[last_elem].next;
+        lower_total = (int)rns.size() - all_total;  // runs->total - all_total;
+        last_elem = cur;
+        rns[cur].next = -1;
+
+        CV_DbgAssert(rns.size() < (size_t)numeric_limits<int>::max());
+
+        // Find links between runs of lower_line and upper_line
+        establishLinks(prev_point, upper_line, lower_line, upper_total, lower_total);
+
+        upper_line = lower_line;
+        upper_total = lower_total;
+    }  // i
+
+    // the last line of image
+    int upper_run = upper_line;
+    for (int k = 0; k < upper_total / 2; k++)
+    {
+        rns[rns[upper_run].next].link = upper_run;
+        upper_run = rns[rns[upper_run].next].next;
+    }
+
+    int first = 0;
+    int prev = 0;
+    convertLinks(first, prev, false);
+    convertLinks(first, prev, true);
+}
+
+}  // namespace
+
+//==============================================================================
+
+void cv::findContoursLinkRuns(InputArray _image,
+                              OutputArrayOfArrays _contours,
+                              OutputArray _hierarchy)
+{
+    CV_INSTRUMENT_REGION();
+
+    CV_CheckType(_image.type(),
+                 _image.type() == CV_8UC1 || _image.type() == CV_8SC1,
+                 "Bad input image type, must be CV_8UC1 or CV_8SC1");
+
+    // Sanity check: output must be of type vector<vector<Point>>
+    CV_Assert(_contours.kind() == _InputArray::STD_VECTOR_VECTOR ||
+              _contours.kind() == _InputArray::STD_VECTOR_MAT ||
+              _contours.kind() == _InputArray::STD_VECTOR_UMAT);
+
+    if (!_contours.empty())
+        CV_CheckTypeEQ(_contours.type(), CV_32SC2, "Contours must have type CV_32SC2");
+
+    if (_hierarchy.needed())
+        _hierarchy.clear();
+
+    Mat image = _image.getMat();
+
+    LinkRunner runner;
+    runner.process(image);
+
+    contourTreeToResults(runner.tree, CV_32SC2, _contours, _hierarchy);
+}
+
+
+void cv::findContoursLinkRuns(InputArray _image, OutputArrayOfArrays _contours)
+{
+    CV_INSTRUMENT_REGION();
+    findContoursLinkRuns(_image, _contours, noArray());
+}

modules/imgproc/src/contours_new.cpp

@@ -0,0 +1,697 @@
+// This file is part of OpenCV project.
+// It is subject to the license terms in the LICENSE file found in the top-level directory
+// of this distribution and at http://opencv.org/license.html
+
+#include "opencv2/imgproc.hpp"
+#include "precomp.hpp"
+#include "opencv2/core/hal/intrin.hpp"
+#include "opencv2/core/check.hpp"
+#include "opencv2/core/utils/logger.hpp"
+#include <iostream>
+#include <array>
+#include <limits>
+#include <map>
+
+#include "contours_common.hpp"
+
+using namespace std;
+using namespace cv;
+
+//==============================================================================
+
+namespace {
+
+template <typename T>
+struct Trait
+{
+};
+
+static const schar MASK8_RIGHT = '\x80';  // 1000 0000
+static const schar MASK8_NEW = '\x02';  // 0000 0010 (+2)
+static const schar MASK8_FLAGS = '\xFE';  // 1111 1110 (-2)
+static const schar MASK8_BLACK = '\x01';  // 0000 0001 - black pixel
+
+static const schar MASK8_LVAL = '\x7F';  // 0111 1111 (for table)
+
+template <>
+struct Trait<schar>
+{
+    static inline bool checkValue(const schar* elem, const schar*)
+    {
+        return *elem != 0;
+    }
+    static inline bool isVal(const schar* elem, const schar*)
+    {
+        return *elem == MASK8_BLACK;
+    }
+    static inline bool isRight(const schar* elem, const schar*)
+    {
+        return (*elem & MASK8_RIGHT) != 0;
+    }
+    static inline void setRightFlag(schar* elem, const schar*, schar nbd)
+    {
+        *elem = nbd | MASK8_RIGHT;
+    }
+    static inline void setNewFlag(schar* elem, const schar*, schar nbd)
+    {
+        *elem = nbd;
+    }
+};
+
+static const int MASK_RIGHT = 0x80000000;  // 100..000
+static const int MASK_NEW = 0x40000000;  // 010..000
+static const int MASK_FLAGS = 0xC0000000;  // right + new
+static const int MASK_VAL = 0x3FFFFFFF;  // ~flags - pixel label
+
+template <>
+struct Trait<int>
+{
+    static inline bool checkValue(const int* elem, const int* elem0)
+    {
+        return (*elem & MASK_VAL) == (*elem0 & MASK_VAL);
+    }
+    static inline bool isVal(const int* elem, const int* elem0)
+    {
+        return *elem == (*elem0 & MASK_VAL);
+    }
+    static inline bool isRight(const int* elem, const int* elem0)
+    {
+        return (*elem & MASK_RIGHT) == (*elem0 & MASK8_RIGHT);
+    }
+    static inline void setRightFlag(int* elem, const int* elem0, int)
+    {
+        *elem = (*elem0 & MASK_VAL) | MASK_NEW | MASK_RIGHT;
+    }
+    static inline void setNewFlag(int* elem, const int* elem0, int)
+    {
+        *elem = (*elem0 & MASK_VAL) | MASK_NEW;
+    }
+};
+
+}  // namespace
+
+
+//==============================================================================
+
+
+namespace {
+
+template <typename T>
+static bool icvTraceContour(Mat& image, const Point& start, const Point& end, bool isHole)
+{
+    const T* stop_ptr = image.ptr<T>(end.y, end.x);
+    const size_t step = image.step1();
+    const T *i0 = image.ptr<T>(start.y, start.x), *i1, *i3, *i4 = NULL;
+    const schar s_end = isHole ? 0 : 4;
+
+    schar s = s_end;
+    do
+    {
+        s = (s - 1) & 7;
+        i1 = i0 + getDelta(s, step);
+    }
+    while (!Trait<T>::checkValue(i1, i0) && s != s_end);
+
+    i3 = i0;
+
+    // check single pixel domain
+    if (s != s_end)
+    {
+        // follow border
+        for (;;)
+        {
+            CV_Assert(i3 != NULL);
+            s = clamp_direction(s);
+            while (s < MAX_SIZE - 1)
+            {
+                ++s;
+                i4 = i3 + getDelta(s, step);
+                CV_Assert(i4 != NULL);
+                if (Trait<T>::checkValue(i4, i0))
+                    break;
+            }
+
+            if (i3 == stop_ptr)
+            {
+                if (!Trait<T>::isRight(i3, i0))
+                {
+                    // it's the only contour
+                    return true;
+                }
+
+                // check if this is the last contour
+                // encountered during a raster scan
+                const T* i5;
+                schar t = s;
+                while (true)
+                {
+                    t = (t - 1) & 7;
+                    i5 = i3 + getDelta(t, step);
+                    if (*i5 != 0)
+                        break;
+                    if (t == 0)
+                        return true;
+                }
+            }
+
+            if ((i4 == i0 && i3 == i1))
+                break;
+
+            i3 = i4;
+            s = (s + 4) & 7;
+        }  // end of border following loop
+    }
+    else
+    {
+        return i3 == stop_ptr;
+    }
+
+    return false;
+}
+
+template <typename T>
+static void icvFetchContourEx(Mat& image,
+                              const Point& start,
+                              T nbd,
+                              Contour& res_contour,
+                              const bool isDirect)
+{
+    const size_t step = image.step1();
+    T *i0 = image.ptr<T>(start.y, start.x), *i1, *i3, *i4 = NULL;
+
+    Point pt = res_contour.origin;
+
+    cv::Rect rect(pt.x, pt.y, pt.x, pt.y);
+
+    schar s_end = res_contour.isHole ? 0 : 4;
+    schar s = s_end;
+    do
+    {
+        s = (s - 1) & 7;
+        i1 = i0 + getDelta(s, step);
+    }
+    while (!Trait<T>::checkValue(i1, i0) && s != s_end);
+
+    if (s == s_end)
+    {
+        Trait<T>::setRightFlag(i0, i0, nbd);
+        if (!res_contour.isChain)
+        {
+            res_contour.pts.push_back(pt);
+        }
+    }
+    else
+    {
+        i3 = i0;
+        schar prev_s = s ^ 4;
+
+        // follow border
+        for (;;)
+        {
+            s_end = s;
+            s = clamp_direction(s);
+            while (s < MAX_SIZE - 1)
+            {
+                ++s;
+                i4 = i3 + getDelta(s, step);
+                CV_Assert(i4 != NULL);
+                if (Trait<T>::checkValue(i4, i0))
+                    break;
+            }
+            s &= 7;
+
+            // check "right" bound
+            if ((unsigned)(s - 1) < (unsigned)s_end)
+            {
+                Trait<T>::setRightFlag(i3, i0, nbd);
+            }
+            else if (Trait<T>::isVal(i3, i0))
+            {
+                Trait<T>::setNewFlag(i3, i0, nbd);
+            }
+
+            if (res_contour.isChain)
+            {
+                res_contour.codes.push_back(s);
+            }
+            else if (s != prev_s || isDirect)
+            {
+                res_contour.pts.push_back(pt);
+            }
+
+            if (s != prev_s)
+            {
+                // update bounds
+                if (pt.x < rect.x)
+                    rect.x = pt.x;
+                else if (pt.x > rect.width)
+                    rect.width = pt.x;
+
+                if (pt.y < rect.y)
+                    rect.y = pt.y;
+                else if (pt.y > rect.height)
+                    rect.height = pt.y;
+            }
+
+            prev_s = s;
+            pt += chainCodeDeltas[s];
+
+            if (i4 == i0 && i3 == i1)
+                break;
+
+            i3 = i4;
+            s = (s + 4) & 7;
+        }
+    }
+    rect.width -= rect.x - 1;
+    rect.height -= rect.y - 1;
+    res_contour.brect = rect;
+}
+
+}  // namespace
+
+
+//==============================================================================
+
+//
+// Raster->Chain Tree (Suzuki algorithms)
+//
+
+// Structure that is used for sequential retrieving contours from the image.
+// It supports both hierarchical and plane variants of Suzuki algorithm.
+struct ContourScanner_
+{
+    Mat image;
+    Point offset;  // ROI offset: coordinates, added to each contour point
+    Point pt;  // current scanner position
+    Point lnbd;  // position of the last met contour
+    schar nbd;  // current mark val
+    int approx_method1;  // approx method when tracing
+    int approx_method2;  // final approx method
+    int mode;
+    CTree tree;
+    array<int, 128> ctable;
+
+public:
+    ContourScanner_() {}
+    ~ContourScanner_() {}
+    inline bool isInt() const
+    {
+        return (this->mode == CV_RETR_FLOODFILL);
+    }
+    inline bool isSimple() const
+    {
+        return (this->mode == RETR_EXTERNAL || this->mode == RETR_LIST);
+    }
+
+    CNode& makeContour(schar& nbd_, const bool is_hole, const int x, const int y);
+    bool contourScan(const int prev, int& p, Point& last_pos, const int x, const int y);
+    int findFirstBoundingContour(const Point& last_pos, const int y, const int lval, int par);
+    int findNextX(int x, int y, int& prev, int& p);
+    bool findNext();
+
+    static shared_ptr<ContourScanner_> create(Mat img, int mode, int method, Point offset);
+};  // class ContourScanner_
+
+typedef shared_ptr<ContourScanner_> ContourScanner;
+
+
+shared_ptr<ContourScanner_> ContourScanner_::create(Mat img, int mode, int method, Point offset)
+{
+    if (mode == RETR_CCOMP && img.type() == CV_32SC1)
+        mode = RETR_FLOODFILL;
+
+    if (mode == RETR_FLOODFILL)
+        CV_CheckTypeEQ(img.type(), CV_32SC1, "RETR_FLOODFILL mode supports only CV_32SC1 images");
+    else
+        CV_CheckTypeEQ(img.type(),
+                       CV_8UC1,
+                       "Modes other than RETR_FLOODFILL and RETR_CCOMP support only CV_8UC1 "
+                       "images");
+
+    CV_Check(mode,
+             mode == RETR_EXTERNAL || mode == RETR_LIST || mode == RETR_CCOMP ||
+                 mode == RETR_TREE || mode == RETR_FLOODFILL,
+             "Wrong extraction mode");
+
+    CV_Check(method,
+             method == 0 || method == CHAIN_APPROX_NONE || method == CHAIN_APPROX_SIMPLE ||
+                 method == CHAIN_APPROX_TC89_L1 || method == CHAIN_APPROX_TC89_KCOS,
+             "Wrong approximation method");
+
+    Size size = img.size();
+    CV_Assert(size.height >= 1);
+
+    shared_ptr<ContourScanner_> scanner = make_shared<ContourScanner_>();
+    scanner->image = img;
+    scanner->mode = mode;
+    scanner->offset = offset;
+    scanner->pt = Point(1, 1);
+    scanner->lnbd = Point(0, 1);
+    scanner->nbd = 2;
+    CNode& root = scanner->tree.newElem();
+    CV_Assert(root.self() == 0);
+    root.body.isHole = true;
+    root.body.brect = Rect(Point(0, 0), size);
+    scanner->ctable.fill(-1);
+    scanner->approx_method2 = scanner->approx_method1 = method;
+    if (method == CV_CHAIN_APPROX_TC89_L1 || method == CV_CHAIN_APPROX_TC89_KCOS)
+        scanner->approx_method1 = CV_CHAIN_CODE;
+    return scanner;
+}
+
+CNode& ContourScanner_::makeContour(schar& nbd_, const bool is_hole, const int x, const int y)
+{
+    const bool isChain = (this->approx_method1 == CV_CHAIN_CODE);  // TODO: get rid of old constant
+    const bool isDirect = (this->approx_method1 == CHAIN_APPROX_NONE);
+
+    const Point start_pt(x - (is_hole ? 1 : 0), y);
+
+    CNode& res = tree.newElem();
+    if (isChain)
+        res.body.codes.reserve(200);
+    else
+        res.body.pts.reserve(200);
+    res.body.isHole = is_hole;
+    res.body.isChain = isChain;
+    res.body.origin = start_pt + offset;
+    if (isSimple())
+    {
+        icvFetchContourEx<schar>(this->image, start_pt, MASK8_NEW, res.body, isDirect);
+    }
+    else
+    {
+        schar lval;
+        if (isInt())
+        {
+            const int start_val = this->image.at<int>(start_pt);
+            lval = start_val & MASK8_LVAL;
+            icvFetchContourEx<int>(this->image, start_pt, 0, res.body, isDirect);
+        }
+        else
+        {
+            lval = nbd_;
+            // change nbd
+            nbd_ = (nbd_ + 1) & MASK8_LVAL;
+            if (nbd_ == 0)
+                nbd_ = MASK8_BLACK | MASK8_NEW;
+            icvFetchContourEx<schar>(this->image, start_pt, lval, res.body, isDirect);
+        }
+        res.body.brect.x -= this->offset.x;
+        res.body.brect.y -= this->offset.y;
+        res.ctable_next = this->ctable[lval];
+        this->ctable[lval] = res.self();
+    }
+    const Point prev_origin = res.body.origin;
+    res.body.origin = start_pt;
+    if (this->approx_method1 != this->approx_method2)
+    {
+        CV_Assert(res.body.isChain);
+        res.body.pts = approximateChainTC89(res.body.codes, prev_origin, this->approx_method2);
+        res.body.isChain = false;
+    }
+    return res;
+}
+
+bool ContourScanner_::contourScan(const int prev, int& p, Point& last_pos, const int x, const int y)
+{
+    bool is_hole = false;
+
+    /* if not external contour */
+    if (isInt())
+    {
+        if (!(((prev & MASK_FLAGS) != 0 || prev == 0) && (p & MASK_FLAGS) == 0))
+        {
+            if ((prev & MASK_FLAGS) != 0 || ((p & MASK_FLAGS) != 0))
+                return false;
+
+            if (prev & MASK_FLAGS)
+            {
+                last_pos.x = x - 1;
+            }
+            is_hole = true;
+        }
+    }
+    else
+    {
+        if (!(prev == 0 && p == 1))
+        {
+            if (p != 0 || prev < 1)
+                return false;
+
+            if (prev & MASK8_FLAGS)
+            {
+                last_pos.x = x - 1;
+            }
+            is_hole = true;
+        }
+    }
+
+    if (mode == RETR_EXTERNAL && (is_hole || this->image.at<schar>(last_pos) > 0))
+    {
+        return false;
+    }
+
+    /* find contour parent */
+    int main_parent = -1;
+    if (isSimple() || (!is_hole && (mode == CV_RETR_CCOMP || mode == CV_RETR_FLOODFILL)) ||
+        last_pos.x <= 0)
+    {
+        main_parent = 0;
+    }
+    else
+    {
+        int lval;
+        if (isInt())
+            lval = this->image.at<int>(last_pos.y, last_pos.x) & MASK8_LVAL;
+        else
+            lval = this->image.at<schar>(last_pos.y, last_pos.x) & MASK8_LVAL;
+
+        main_parent = findFirstBoundingContour(last_pos, y, lval, main_parent);
+
+        // if current contour is a hole and previous contour is a hole or
+        // current contour is external and previous contour is external then
+        // the parent of the contour is the parent of the previous contour else
+        // the parent is the previous contour itself.
+        {
+            CNode& main_parent_elem = tree.elem(main_parent);
+            if (main_parent_elem.body.isHole == is_hole)
+            {
+                if (main_parent_elem.parent != -1)
+                {
+                    main_parent = main_parent_elem.parent;
+                }
+                else
+                {
+                    main_parent = 0;
+                }
+            }
+        }
+
+        // hole flag of the parent must differ from the flag of the contour
+        {
+            CNode& main_parent_elem = tree.elem(main_parent);
+            CV_Assert(main_parent_elem.body.isHole != is_hole);
+        }
+    }
+
+    last_pos.x = x - (is_hole ? 1 : 0);
+
+    schar nbd_ = this->nbd;
+    CNode& new_contour = makeContour(nbd_, is_hole, x, y);
+    if (new_contour.parent == -1)
+    {
+        tree.addChild(main_parent, new_contour.self());
+    }
+    this->pt.x = !isInt() ? (x + 1) : (x + 1 - (is_hole ? 1 : 0));
+    this->pt.y = y;
+    this->nbd = nbd_;
+    return true;
+}
+
+int ContourScanner_::findFirstBoundingContour(const Point& last_pos,
+                                              const int y,
+                                              const int lval,
+                                              int par)
+{
+    const Point end_point(last_pos.x, y);
+    int res = par;
+    int cur = ctable[lval];
+    while (cur != -1)
+    {
+        CNode& cur_elem = tree.elem(cur);
+        if (((last_pos.x - cur_elem.body.brect.x) < cur_elem.body.brect.width) &&
+            ((last_pos.y - cur_elem.body.brect.y) < cur_elem.body.brect.height))
+        {
+            if (res != -1)
+            {
+                CNode& res_elem = tree.elem(res);
+                const Point origin = res_elem.body.origin;
+                const bool isHole = res_elem.body.isHole;
+                if (isInt())
+                {
+                    if (icvTraceContour<int>(this->image, origin, end_point, isHole))
+                        break;
+                }
+                else
+                {
+                    if (icvTraceContour<schar>(this->image, origin, end_point, isHole))
+                        break;
+                }
+            }
+            res = cur;
+        }
+        cur = cur_elem.ctable_next;
+    }
+    return res;
+}
+
+int ContourScanner_::findNextX(int x, int y, int& prev, int& p)
+{
+    const int width = this->image.size().width - 1;
+    if (isInt())
+    {
+        for (; x < width &&
+               ((p = this->image.at<int>(y, x)) == prev || (p & MASK_VAL) == (prev & MASK_VAL));
+             x++)
+            prev = p;
+    }
+    else
+    {
+#if (CV_SIMD || CV_SIMD_SCALABLE)
+        if ((p = this->image.at<schar>(y, x)) != prev)
+        {
+            return x;
+        }
+        else
+        {
+            v_uint8 v_prev = vx_setall_u8((uchar)prev);
+            for (; x <= width - VTraits<v_uint8>::vlanes(); x += VTraits<v_uint8>::vlanes())
+            {
+                v_uint8 vmask = (v_ne(vx_load(this->image.ptr<uchar>(y, x)), v_prev));
+                if (v_check_any(vmask))
+                {
+                    x += v_scan_forward(vmask);
+                    p = this->image.at<schar>(y, x);
+                    return x;
+                }
+            }
+        }
+#endif
+        for (; x < width && (p = this->image.at<schar>(y, x)) == prev; x++)
+            ;
+    }
+    return x;
+}
+
+bool ContourScanner_::findNext()
+{
+    int x = this->pt.x;
+    int y = this->pt.y;
+    int width = this->image.size().width - 1;
+    int height = this->image.size().height - 1;
+    Point last_pos = this->lnbd;
+    int prev = isInt() ? this->image.at<int>(y, x - 1) : this->image.at<schar>(y, x - 1);
+
+    for (; y < height; y++)
+    {
+        int p = 0;
+        for (; x < width; x++)
+        {
+            x = findNextX(x, y, prev, p);
+            if (x >= width)
+                break;
+            if (contourScan(prev, p, last_pos, x, y))
+            {
+                this->lnbd = last_pos;
+                return true;
+            }
+            else
+            {
+                prev = p;
+                if ((isInt() && (prev & MASK_FLAGS)) || (!isInt() && (prev & MASK8_FLAGS)))
+                {
+                    last_pos.x = x;
+                }
+            }
+        }
+        last_pos = Point(0, y + 1);
+        x = 1;
+        prev = 0;
+    }
+
+    return false;
+}
+
+//==============================================================================
+
+void cv::findContours(InputArray _image,
+                      OutputArrayOfArrays _contours,
+                      OutputArray _hierarchy,
+                      int mode,
+                      int method,
+                      Point offset)
+{
+    CV_INSTRUMENT_REGION();
+
+    // TODO: remove this block in future
+    if (method == 5 /*CV_LINK_RUNS*/)
+    {
+        CV_LOG_ONCE_WARNING(NULL,
+                            "LINK_RUNS mode has been extracted to separate function: "
+                            "cv::findContoursLinkRuns. "
+                            "Calling through cv::findContours will be removed in future.");
+        CV_CheckTrue(!_hierarchy.needed() || mode == RETR_CCOMP,
+                     "LINK_RUNS mode supports only simplified hierarchy output (mode=RETR_CCOMP)");
+        findContoursLinkRuns(_image, _contours, _hierarchy);
+        return;
+    }
+
+    // TODO: need enum value, need way to return contour starting points with chain codes
+    if (method == 0 /*CV_CHAIN_CODE*/)
+    {
+        CV_LOG_ONCE_WARNING(NULL,
+                            "Chain code output is an experimental feature and might change in "
+                            "future!");
+    }
+
+    // Sanity check: output must be of type vector<vector<Point>>
+    CV_Assert((_contours.kind() == _InputArray::STD_VECTOR_VECTOR) ||
+              (_contours.kind() == _InputArray::STD_VECTOR_MAT) ||
+              (_contours.kind() == _InputArray::STD_VECTOR_UMAT));
+
+    const int res_type = (method == 0 /*CV_CHAIN_CODE*/) ? CV_8SC1 : CV_32SC2;
+    if (!_contours.empty())
+    {
+        CV_CheckTypeEQ(_contours.type(),
+                       res_type,
+                       "Contours must have type CV_8SC1 (chain code) or CV_32SC2 (other methods)");
+    }
+
+    if (_hierarchy.needed())
+        _hierarchy.clear();
+
+    // preprocess
+    Mat image;
+    copyMakeBorder(_image, image, 1, 1, 1, 1, BORDER_CONSTANT | BORDER_ISOLATED, Scalar(0));
+    if (image.type() != CV_32SC1)
+        threshold(image, image, 0, 1, THRESH_BINARY);
+
+    // find contours
+    ContourScanner scanner = ContourScanner_::create(image, mode, method, offset + Point(-1, -1));
+    while (scanner->findNext())
+    {
+    }
+
+    contourTreeToResults(scanner->tree, res_type, _contours, _hierarchy);
+}
+
+void cv::findContours(InputArray _image,
+                      OutputArrayOfArrays _contours,
+                      int mode,
+                      int method,
+                      Point offset)
+{
+    CV_INSTRUMENT_REGION();
+    findContours(_image, _contours, noArray(), mode, method, offset);
+}

modules/imgproc/test/test_contours.cpp

@@ -39,6 +39,7 @@
 //
 //M*/
 
+#include "opencv2/imgproc/types_c.h"
 #include "test_precomp.hpp"
 #include <opencv2/highgui.hpp>
 
@@ -459,7 +460,6 @@ TEST(Imgproc_FindContours, hilbert)
     dilate(img, img, Mat());
     vector<vector<Point> > contours;
     findContours(img, contours, noArray(), RETR_LIST, CHAIN_APPROX_SIMPLE);
-    printf("ncontours = %d, contour[0].npoints=%d\n", (int)contours.size(), (int)contours[0].size());
     img.setTo(Scalar::all(0));
 
     drawContours(img, contours, 0, Scalar::all(255), 1);
@@ -530,10 +530,12 @@ TEST(Imgproc_FindContours, regression_4363_shared_nbd)
 
     if (found)
     {
+        ASSERT_EQ(contours.size(), hierarchy.size());
         EXPECT_LT(hierarchy[index][3], 0) << "Desired result: (7,9) has no parent - Actual result: parent of (7,9) is another contour. index = " << index;
     }
 }
 
+
 TEST(Imgproc_PointPolygonTest, regression_10222)
 {
     vector<Point> contour;

modules/imgproc/test/test_contours_new.cpp

@@ -0,0 +1,606 @@
+// This file is part of OpenCV project.
+// It is subject to the license terms in the LICENSE file found in the top-level directory
+// of this distribution and at http://opencv.org/license.html
+
+#include "opencv2/imgproc/types_c.h"
+#include "test_precomp.hpp"
+#include "opencv2/ts/ocl_test.hpp"
+#include "opencv2/imgproc/detail/legacy.hpp"
+
+#define CHECK_OLD 1
+
+namespace opencv_test { namespace {
+
+// debug function
+template <typename T>
+inline static void print_pts(const T& c)
+{
+    for (const auto& one_pt : c)
+    {
+        cout << one_pt << " ";
+    }
+    cout << endl;
+}
+
+// debug function
+template <typename T>
+inline static void print_pts_2(vector<T>& cs)
+{
+    int cnt = 0;
+    cout << "Contours:" << endl;
+    for (const auto& one_c : cs)
+    {
+        cout << cnt++ << " : ";
+        print_pts(one_c);
+    }
+};
+
+// draw 1-2 px blob with orientation defined by 'kind'
+template <typename T>
+inline static void drawSmallContour(Mat& img, Point pt, int kind, int color_)
+{
+    const T color = static_cast<T>(color_);
+    img.at<T>(pt) = color;
+    switch (kind)
+    {
+        case 1: img.at<T>(pt + Point(1, 0)) = color; break;
+        case 2: img.at<T>(pt + Point(1, -1)) = color; break;
+        case 3: img.at<T>(pt + Point(0, -1)) = color; break;
+        case 4: img.at<T>(pt + Point(-1, -1)) = color; break;
+        case 5: img.at<T>(pt + Point(-1, 0)) = color; break;
+        case 6: img.at<T>(pt + Point(-1, 1)) = color; break;
+        case 7: img.at<T>(pt + Point(0, 1)) = color; break;
+        case 8: img.at<T>(pt + Point(1, 1)) = color; break;
+        default: break;
+    }
+}
+
+inline static void drawContours(Mat& img,
+                                const vector<vector<Point>>& contours,
+                                const Scalar& color = Scalar::all(255))
+{
+    for (const auto& contour : contours)
+    {
+        for (size_t n = 0, end = contour.size(); n < end; ++n)
+        {
+            size_t m = n + 1;
+            if (n == end - 1)
+                m = 0;
+            line(img, contour[m], contour[n], color, 1, LINE_8);
+        }
+    }
+}
+
+//==================================================================================================
+
+// Test parameters - mode + method
+typedef testing::TestWithParam<tuple<int, int>> Imgproc_FindContours_Modes1;
+
+
+// Draw random rectangle and find contours
+//
+TEST_P(Imgproc_FindContours_Modes1, rectangle)
+{
+    const int mode = get<0>(GetParam());
+    const int method = get<1>(GetParam());
+
+    const size_t ITER = 100;
+    RNG rng = TS::ptr()->get_rng();
+
+    for (size_t i = 0; i < ITER; ++i)
+    {
+        SCOPED_TRACE(cv::format("i=%zu", i));
+        const Size sz(rng.uniform(640, 1920), rng.uniform(480, 1080));
+        Mat img(sz, CV_8UC1, Scalar::all(0));
+        Mat img32s(sz, CV_32SC1, Scalar::all(0));
+        const Rect r(Point(rng.uniform(1, sz.width / 2 - 1), rng.uniform(1, sz.height / 2)),
+                     Point(rng.uniform(sz.width / 2 - 1, sz.width - 1),
+                           rng.uniform(sz.height / 2 - 1, sz.height - 1)));
+        rectangle(img, r, Scalar::all(255));
+        rectangle(img32s, r, Scalar::all(255), FILLED);
+
+        const vector<Point> ext_ref {r.tl(),
+                                     r.tl() + Point(0, r.height - 1),
+                                     r.br() + Point(-1, -1),
+                                     r.tl() + Point(r.width - 1, 0)};
+        const vector<Point> int_ref {ext_ref[0] + Point(0, 1),
+                                     ext_ref[0] + Point(1, 0),
+                                     ext_ref[3] + Point(-1, 0),
+                                     ext_ref[3] + Point(0, 1),
+                                     ext_ref[2] + Point(0, -1),
+                                     ext_ref[2] + Point(-1, 0),
+                                     ext_ref[1] + Point(1, 0),
+                                     ext_ref[1] + Point(0, -1)};
+        const size_t ext_perimeter = r.width * 2 + r.height * 2;
+        const size_t int_perimeter = ext_perimeter - 4;
+
+        vector<vector<Point>> contours;
+        vector<vector<schar>> chains;
+        vector<Vec4i> hierarchy;
+
+        // run functionn
+        if (mode == RETR_FLOODFILL)
+            if (method == 0)
+                findContours(img32s, chains, hierarchy, mode, method);
+            else
+                findContours(img32s, contours, hierarchy, mode, method);
+        else if (method == 0)
+            findContours(img, chains, hierarchy, mode, method);
+        else
+            findContours(img, contours, hierarchy, mode, method);
+
+        // verify results
+        if (mode == RETR_EXTERNAL)
+        {
+            if (method == 0)
+            {
+                ASSERT_EQ(1U, chains.size());
+            }
+            else
+            {
+                ASSERT_EQ(1U, contours.size());
+                if (method == CHAIN_APPROX_NONE)
+                {
+                    EXPECT_EQ(int_perimeter, contours[0].size());
+                }
+                else if (method == CHAIN_APPROX_SIMPLE)
+                {
+                    EXPECT_MAT_NEAR(Mat(ext_ref), Mat(contours[0]), 0);
+                }
+            }
+        }
+        else
+        {
+            if (method == 0)
+            {
+                ASSERT_EQ(2U, chains.size());
+            }
+            else
+            {
+                ASSERT_EQ(2U, contours.size());
+                if (mode == RETR_LIST)
+                {
+                    if (method == CHAIN_APPROX_NONE)
+                    {
+                        EXPECT_EQ(int_perimeter - 4, contours[0].size());
+                        EXPECT_EQ(int_perimeter, contours[1].size());
+                    }
+                    else if (method == CHAIN_APPROX_SIMPLE)
+                    {
+                        EXPECT_MAT_NEAR(Mat(int_ref), Mat(contours[0]), 0);
+                        EXPECT_MAT_NEAR(Mat(ext_ref), Mat(contours[1]), 0);
+                    }
+                }
+                else if (mode == RETR_CCOMP || mode == RETR_TREE)
+                {
+                    if (method == CHAIN_APPROX_NONE)
+                    {
+                        EXPECT_EQ(int_perimeter, contours[0].size());
+                        EXPECT_EQ(int_perimeter - 4, contours[1].size());
+                    }
+                    else if (method == CHAIN_APPROX_SIMPLE)
+                    {
+                        EXPECT_MAT_NEAR(Mat(ext_ref), Mat(contours[0]), 0);
+                        EXPECT_MAT_NEAR(Mat(int_ref), Mat(contours[1]), 0);
+                    }
+                }
+                else if (mode == RETR_FLOODFILL)
+                {
+                    if (method == CHAIN_APPROX_NONE)
+                    {
+                        EXPECT_EQ(int_perimeter + 4, contours[0].size());
+                    }
+                    else if (method == CHAIN_APPROX_SIMPLE)
+                    {
+                        EXPECT_EQ(int_ref.size(), contours[0].size());
+                        EXPECT_MAT_NEAR(Mat(ext_ref), Mat(contours[1]), 0);
+                    }
+                }
+            }
+        }
+
+#if CHECK_OLD
+        if (method != 0)  // old doesn't support chain codes
+        {
+            if (mode != RETR_FLOODFILL)
+            {
+                vector<vector<Point>> contours_o;
+                vector<Vec4i> hierarchy_o;
+                findContours_legacy(img, contours_o, hierarchy_o, mode, method);
+                ASSERT_EQ(contours.size(), contours_o.size());
+                for (size_t j = 0; j < contours.size(); ++j)
+                {
+                    SCOPED_TRACE(format("contour %zu", j));
+                    EXPECT_MAT_NEAR(Mat(contours[j]), Mat(contours_o[j]), 0);
+                }
+                EXPECT_MAT_NEAR(Mat(hierarchy), Mat(hierarchy_o), 0);
+            }
+            else
+            {
+                vector<vector<Point>> contours_o;
+                vector<Vec4i> hierarchy_o;
+                findContours_legacy(img32s, contours_o, hierarchy_o, mode, method);
+                ASSERT_EQ(contours.size(), contours_o.size());
+                for (size_t j = 0; j < contours.size(); ++j)
+                {
+                    SCOPED_TRACE(format("contour %zu", j));
+                    EXPECT_MAT_NEAR(Mat(contours[j]), Mat(contours_o[j]), 0);
+                }
+                EXPECT_MAT_NEAR(Mat(hierarchy), Mat(hierarchy_o), 0);
+            }
+        }
+#endif
+    }
+}
+
+
+// Draw many small 1-2px blobs and find contours
+//
+TEST_P(Imgproc_FindContours_Modes1, small)
+{
+    const int mode = get<0>(GetParam());
+    const int method = get<1>(GetParam());
+
+    const size_t DIM = 1000;
+    const Size sz(DIM, DIM);
+    const int num = (DIM / 10) * (DIM / 10);  // number of 10x10 squares
+
+    Mat img(sz, CV_8UC1, Scalar::all(0));
+    Mat img32s(sz, CV_32SC1, Scalar::all(0));
+    vector<Point> pts;
+    int extra_contours_32s = 0;
+    for (int j = 0; j < num; ++j)
+    {
+        const int kind = j % 9;
+        Point pt {(j % 100) * 10 + 4, (j / 100) * 10 + 4};
+        drawSmallContour<uchar>(img, pt, kind, 255);
+        drawSmallContour<int>(img32s, pt, kind, j + 1);
+        pts.push_back(pt);
+        // NOTE: for some reason these small diagonal contours (NW, SE)
+        //       result in 2 external contours for FLOODFILL mode
+        if (kind == 8 || kind == 4)
+            ++extra_contours_32s;
+    }
+    {
+        vector<vector<Point>> contours;
+        vector<vector<schar>> chains;
+        vector<Vec4i> hierarchy;
+
+        if (mode == RETR_FLOODFILL)
+        {
+            if (method == 0)
+            {
+                findContours(img32s, chains, hierarchy, mode, method);
+                ASSERT_EQ(pts.size() * 2 + extra_contours_32s, chains.size());
+            }
+            else
+            {
+                findContours(img32s, contours, hierarchy, mode, method);
+                ASSERT_EQ(pts.size() * 2 + extra_contours_32s, contours.size());
+#if CHECK_OLD
+                vector<vector<Point>> contours_o;
+                vector<Vec4i> hierarchy_o;
+                findContours_legacy(img32s, contours_o, hierarchy_o, mode, method);
+                ASSERT_EQ(contours.size(), contours_o.size());
+                for (size_t i = 0; i < contours.size(); ++i)
+                {
+                    SCOPED_TRACE(format("contour %zu", i));
+                    EXPECT_MAT_NEAR(Mat(contours[i]), Mat(contours_o[i]), 0);
+                }
+                EXPECT_MAT_NEAR(Mat(hierarchy), Mat(hierarchy_o), 0);
+#endif
+            }
+        }
+        else
+        {
+            if (method == 0)
+            {
+                findContours(img, chains, hierarchy, mode, method);
+                ASSERT_EQ(pts.size(), chains.size());
+            }
+            else
+            {
+                findContours(img, contours, hierarchy, mode, method);
+                ASSERT_EQ(pts.size(), contours.size());
+#if CHECK_OLD
+                vector<vector<Point>> contours_o;
+                vector<Vec4i> hierarchy_o;
+                findContours_legacy(img, contours_o, hierarchy_o, mode, method);
+                ASSERT_EQ(contours.size(), contours_o.size());
+                for (size_t i = 0; i < contours.size(); ++i)
+                {
+                    SCOPED_TRACE(format("contour %zu", i));
+                    EXPECT_MAT_NEAR(Mat(contours[i]), Mat(contours_o[i]), 0);
+                }
+                EXPECT_MAT_NEAR(Mat(hierarchy), Mat(hierarchy_o), 0);
+#endif
+            }
+        }
+    }
+}
+
+
+// Draw many nested rectangles and find contours
+//
+TEST_P(Imgproc_FindContours_Modes1, deep)
+{
+    const int mode = get<0>(GetParam());
+    const int method = get<1>(GetParam());
+
+    const size_t DIM = 1000;
+    const Size sz(DIM, DIM);
+    const size_t NUM = 249U;
+    Mat img(sz, CV_8UC1, Scalar::all(0));
+    Mat img32s(sz, CV_32SC1, Scalar::all(0));
+    Rect rect(1, 1, 998, 998);
+    for (size_t i = 0; i < NUM; ++i)
+    {
+        rectangle(img, rect, Scalar::all(255));
+        rectangle(img32s, rect, Scalar::all((double)i + 1), FILLED);
+        rect.x += 2;
+        rect.y += 2;
+        rect.width -= 4;
+        rect.height -= 4;
+    }
+    {
+        vector<vector<Point>> contours {{{0, 0}, {1, 1}}};
+        vector<vector<schar>> chains {{1, 2, 3}};
+        vector<Vec4i> hierarchy;
+
+        if (mode == RETR_FLOODFILL)
+        {
+            if (method == 0)
+            {
+                findContours(img32s, chains, hierarchy, mode, method);
+                ASSERT_EQ(2 * NUM, chains.size());
+            }
+            else
+            {
+                findContours(img32s, contours, hierarchy, mode, method);
+                ASSERT_EQ(2 * NUM, contours.size());
+#if CHECK_OLD
+                vector<vector<Point>> contours_o;
+                vector<Vec4i> hierarchy_o;
+                findContours_legacy(img32s, contours_o, hierarchy_o, mode, method);
+                ASSERT_EQ(contours.size(), contours_o.size());
+                for (size_t i = 0; i < contours.size(); ++i)
+                {
+                    SCOPED_TRACE(format("contour %zu", i));
+                    EXPECT_MAT_NEAR(Mat(contours[i]), Mat(contours_o[i]), 0);
+                }
+                EXPECT_MAT_NEAR(Mat(hierarchy), Mat(hierarchy_o), 0);
+#endif
+            }
+        }
+        else
+        {
+            const size_t expected_count = (mode == RETR_EXTERNAL) ? 1U : 2 * NUM;
+            if (method == 0)
+            {
+                findContours(img, chains, hierarchy, mode, method);
+                ASSERT_EQ(expected_count, chains.size());
+            }
+            else
+            {
+                findContours(img, contours, hierarchy, mode, method);
+                ASSERT_EQ(expected_count, contours.size());
+#if CHECK_OLD
+                vector<vector<Point>> contours_o;
+                vector<Vec4i> hierarchy_o;
+                findContours_legacy(img, contours_o, hierarchy_o, mode, method);
+                ASSERT_EQ(contours.size(), contours_o.size());
+                for (size_t i = 0; i < contours.size(); ++i)
+                {
+                    SCOPED_TRACE(format("contour %zu", i));
+                    EXPECT_MAT_NEAR(Mat(contours[i]), Mat(contours_o[i]), 0);
+                }
+                EXPECT_MAT_NEAR(Mat(hierarchy), Mat(hierarchy_o), 0);
+#endif
+            }
+        }
+    }
+}
+
+INSTANTIATE_TEST_CASE_P(
+    ,
+    Imgproc_FindContours_Modes1,
+    testing::Combine(
+        testing::Values(RETR_EXTERNAL, RETR_LIST, RETR_CCOMP, RETR_TREE, RETR_FLOODFILL),
+        testing::Values(0,
+                        CHAIN_APPROX_NONE,
+                        CHAIN_APPROX_SIMPLE,
+                        CHAIN_APPROX_TC89_L1,
+                        CHAIN_APPROX_TC89_KCOS)));
+
+//==================================================================================================
+
+typedef testing::TestWithParam<tuple<int, int>> Imgproc_FindContours_Modes2;
+
+// Very approximate backport of an old accuracy test
+//
+TEST_P(Imgproc_FindContours_Modes2, new_accuracy)
+{
+    const int mode = get<0>(GetParam());
+    const int method = get<1>(GetParam());
+
+    RNG& rng = TS::ptr()->get_rng();
+    const int blob_count = rng.uniform(1, 10);
+    const Size sz(rng.uniform(640, 1920), rng.uniform(480, 1080));
+    const int blob_sz = 50;
+
+    // prepare image
+    Mat img(sz, CV_8UC1, Scalar::all(0));
+    vector<RotatedRect> rects;
+    for (int i = 0; i < blob_count; ++i)
+    {
+        const Point2f center((float)rng.uniform(blob_sz, sz.width - blob_sz),
+                             (float)rng.uniform(blob_sz, sz.height - blob_sz));
+        const Size2f rsize((float)rng.uniform(1, blob_sz), (float)rng.uniform(1, blob_sz));
+        RotatedRect rect(center, rsize, rng.uniform(0.f, 180.f));
+        rects.push_back(rect);
+        ellipse(img, rect, Scalar::all(100), FILLED);
+    }
+
+    // draw contours manually
+    Mat cont_img(sz, CV_8UC1, Scalar::all(0));
+    for (int y = 1; y < sz.height - 1; ++y)
+    {
+        for (int x = 1; x < sz.width - 1; ++x)
+        {
+            if (img.at<uchar>(y, x) != 0 &&
+                ((img.at<uchar>(y - 1, x) == 0) || (img.at<uchar>(y + 1, x) == 0) ||
+                 (img.at<uchar>(y, x + 1) == 0) || (img.at<uchar>(y, x - 1) == 0)))
+            {
+                cont_img.at<uchar>(y, x) = 255;
+            }
+        }
+    }
+
+    // find contours
+    vector<vector<Point>> contours;
+    vector<Vec4i> hierarchy;
+    findContours(img, contours, hierarchy, mode, method);
+
+    // 0 < contours <= rects
+    EXPECT_GT(contours.size(), 0U);
+    EXPECT_GE(rects.size(), contours.size());
+
+    // draw contours
+    Mat res_img(sz, CV_8UC1, Scalar::all(0));
+    drawContours(res_img, contours);
+
+    // compare resulting drawn contours with manually drawn contours
+    const double diff1 = cvtest::norm(cont_img, res_img, NORM_L1) / 255;
+
+    if (method == CHAIN_APPROX_NONE || method == CHAIN_APPROX_SIMPLE)
+    {
+        EXPECT_EQ(0., diff1);
+    }
+#if CHECK_OLD
+    vector<vector<Point>> contours_o;
+    vector<Vec4i> hierarchy_o;
+    findContours(img, contours_o, hierarchy_o, mode, method);
+    ASSERT_EQ(contours_o.size(), contours.size());
+    for (size_t i = 0; i < contours_o.size(); ++i)
+    {
+        SCOPED_TRACE(format("contour = %zu", i));
+        EXPECT_MAT_NEAR(Mat(contours_o[i]), Mat(contours[i]), 0);
+    }
+    EXPECT_MAT_NEAR(Mat(hierarchy_o), Mat(hierarchy), 0);
+#endif
+}
+
+TEST_P(Imgproc_FindContours_Modes2, approx)
+{
+    const int mode = get<0>(GetParam());
+    const int method = get<1>(GetParam());
+
+    const Size sz {500, 500};
+    Mat img = Mat::zeros(sz, CV_8UC1);
+
+    for (int c = 0; c < 4; ++c)
+    {
+        if (c != 0)
+        {
+            // noise + filter + threshold
+            RNG& rng = TS::ptr()->get_rng();
+            cvtest::randUni(rng, img, 0, 255);
+
+            Mat fimg;
+            boxFilter(img, fimg, CV_8U, Size(5, 5));
+
+            Mat timg;
+            const int level = 44 + c * 42;
+            // 'level' goes through:
+            // 86 - some black speckles on white
+            // 128 - 50/50 black/white
+            // 170 - some white speckles on black
+            cv::threshold(fimg, timg, level, 255, THRESH_BINARY);
+        }
+        else
+        {
+            // circle with cut
+            const Point center {250, 250};
+            const int r {20};
+            const Point cut {r, r};
+            circle(img, center, r, Scalar(255), FILLED);
+            rectangle(img, center, center + cut, Scalar(0), FILLED);
+        }
+
+        vector<vector<Point>> contours;
+        vector<Vec4i> hierarchy;
+        findContours(img, contours, hierarchy, mode, method);
+
+#if CHECK_OLD
+        vector<vector<Point>> contours_o;
+        vector<Vec4i> hierarchy_o;
+        findContours_legacy(img, contours_o, hierarchy_o, mode, method);
+        ASSERT_EQ(contours_o.size(), contours.size());
+        for (size_t i = 0; i < contours_o.size(); ++i)
+        {
+            SCOPED_TRACE(format("c = %d, contour = %zu", c, i));
+            EXPECT_MAT_NEAR(Mat(contours_o[i]), Mat(contours[i]), 0);
+        }
+        EXPECT_MAT_NEAR(Mat(hierarchy_o), Mat(hierarchy), 0);
+#endif
+        // TODO: check something
+    }
+}
+
+// TODO: offset test
+
+// no RETR_FLOODFILL - no CV_32S input images
+INSTANTIATE_TEST_CASE_P(
+    ,
+    Imgproc_FindContours_Modes2,
+    testing::Combine(testing::Values(RETR_EXTERNAL, RETR_LIST, RETR_CCOMP, RETR_TREE),
+                     testing::Values(CHAIN_APPROX_NONE,
+                                     CHAIN_APPROX_SIMPLE,
+                                     CHAIN_APPROX_TC89_L1,
+                                     CHAIN_APPROX_TC89_KCOS)));
+
+TEST(Imgproc_FindContours, link_runs)
+{
+    const Size sz {500, 500};
+    Mat img = Mat::zeros(sz, CV_8UC1);
+
+    // noise + filter + threshold
+    RNG& rng = TS::ptr()->get_rng();
+    cvtest::randUni(rng, img, 0, 255);
+
+    Mat fimg;
+    boxFilter(img, fimg, CV_8U, Size(5, 5));
+
+    const int level = 135;
+    cv::threshold(fimg, img, level, 255, THRESH_BINARY);
+
+    vector<vector<Point>> contours;
+    vector<Vec4i> hierarchy;
+    findContoursLinkRuns(img, contours, hierarchy);
+
+    if (cvtest::debugLevel >= 10)
+    {
+        print_pts_2(contours);
+
+        Mat res = Mat::zeros(sz, CV_8UC1);
+        drawContours(res, contours);
+        imshow("res", res);
+        imshow("img", img);
+        waitKey(0);
+    }
+
+#if CHECK_OLD
+    vector<vector<Point>> contours_o;
+    vector<Vec4i> hierarchy_o;
+    findContours_legacy(img, contours_o, hierarchy_o, 0, 5);  // CV_LINK_RUNS method
+    ASSERT_EQ(contours_o.size(), contours.size());
+    for (size_t i = 0; i < contours_o.size(); ++i)
+    {
+        SCOPED_TRACE(format("contour = %zu", i));
+        EXPECT_MAT_NEAR(Mat(contours_o[i]), Mat(contours[i]), 0);
+    }
+    EXPECT_MAT_NEAR(Mat(hierarchy_o), Mat(hierarchy), 0);
+#endif
+}
+
+}}  // namespace opencv_test