Skip to content

Add fast grayscale morphology using sparse table data structure #3929

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
wants to merge 1 commit into
base: 4.x
Choose a base branch
from
Open
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 1 addition & 0 deletions modules/ximgproc/include/opencv2/ximgproc.hpp
Original file line number Diff line number Diff line change
@@ -63,6 +63,7 @@
#include "ximgproc/color_match.hpp"
#include "ximgproc/radon_transform.hpp"
#include "ximgproc/find_ellipses.hpp"
#include "ximgproc/sparse_table_morphology.hpp"


/**
167 changes: 167 additions & 0 deletions modules/ximgproc/include/opencv2/ximgproc/sparse_table_morphology.hpp
Original file line number Diff line number Diff line change
@@ -0,0 +1,167 @@
// 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_SPARSE_TABLE_MORPHOLOGY_HPP__
#define __OPENCV_SPARSE_TABLE_MORPHOLOGY_HPP__

#include <opencv2/core.hpp>
#include <vector>

namespace cv {
namespace ximgproc {
namespace stMorph {

//! @addtogroup imgproc_filter
//! @{

/**
* @struct kernelDecompInfo
* @brief struct to hold the results of decomposing the structuring element.
*/
struct CV_EXPORTS kernelDecompInfo
{
//! rows of the original kernel.
int rows;
//! cols of the original kernel.
int cols;
//!
//! set of rectangles to covers the kernel which height and width both are power of 2.
//! point stRects[rd][cd](c,r) means a rectangle left-top (c,r), width 2^rd and height 2^cd.
//!
std::vector<std::vector<std::vector<Point>>> stRects;
//!
//! Vec2b Mat which sotres the order to calculate sparse table.
//! The type of returned mat is Vec2b.
//! * if path[dr][dc][0] == 1 then st[dr+1][dc] will be calculated from st[dr][dc].
//! * if path[dr][dc][1] == 1 then st[dr][dc+1] will be calculated from st[dr][dc].
//!
Mat plan;
//! anchor position of the kernel.
Point anchor;
//! Number of times erosion and dilation are applied.
int iterations;
};

/**
* @brief Decompose the structuring element.
*
* @param kernel structuring element used for subsequent morphological operations.
* @param anchor position of the anchor within the element.
* default value (-1, -1) means that the anchor is at the element center.
* @param iterations number of times is applied.
*/
CV_EXPORTS kernelDecompInfo decompKernel(InputArray kernel,
Point anchor = Point(-1, -1), int iterations = 1);

/**
* @brief Erodes an image with a kernelDecompInfo using spase table method.
*
* @param src input image
* @param dst output image of the same size and type as src.
* @param kdi pre-computated kernelDecompInfo structure.
* @param borderType pixel extrapolation method, see #BorderTypes. #BORDER_WRAP is not supported.
* @param borderValue border value in case of a constant border
*/
CV_EXPORTS void erode( InputArray src, OutputArray dst, kernelDecompInfo kdi,
BorderTypes borderType = BORDER_CONSTANT,
const Scalar& borderValue = morphologyDefaultBorderValue() );

/**
* @brief Dilates an image with a kernelDecompInfo using spase table method.
*
* @param src input image;
* @param dst output image of the same size and type as src.
* @param kdi pre-computated kernelDecompInfo structure.
* @param borderType pixel extrapolation method, see #BorderTypes. #BORDER_WRAP is not supported.
* @param borderValue border value in case of a constant border
*/
CV_EXPORTS void dilate( InputArray src, OutputArray dst, kernelDecompInfo kdi,
BorderTypes borderType = BORDER_CONSTANT,
const Scalar& borderValue = morphologyDefaultBorderValue() );

/**
* @brief Performs advanced morphological transformations with a kernelDecompInfo.
*
* @param src input image;
* @param dst output image of the same size and type as src.
* @param op all operations supported by cv::morphologyEx (except cv::MORPH_HITMISS)
* @param kdi pre-computated kernelDecompInfo structure.
* @param borderType pixel extrapolation method, see #BorderTypes. #BORDER_WRAP is not supported.
* @param borderValue border value in case of a constant border
*/
CV_EXPORTS void morphologyEx( InputArray src, OutputArray dst, int op, kernelDecompInfo kdi,
BorderTypes borderType = BORDER_CONSTANT,
const Scalar& borderValue = morphologyDefaultBorderValue() );

/**
* @brief Faster implementation of cv::erode with sparse table concept.
*
* @param src input image; the number of channels can be arbitrary, but the depth should be one of
* CV_8U, CV_16U, CV_16S, CV_32F or CV_64F.
* @param dst output image of the same size and type as src.
* @param kernel structuring element used for erosion; if `element=Mat()`, a `3 x 3` rectangular
* structuring element is used. Kernel can be created using #getStructuringElement.
* @param anchor position of the anchor within the element; default value (-1, -1) means that the
* anchor is at the element center.
* @param iterations number of times erosion is applied.
* @param borderType pixel extrapolation method, see #BorderTypes. #BORDER_WRAP is not supported.
* @param borderValue border value in case of a constant border
*
* @see cv::erode
*/
CV_EXPORTS void erode( InputArray src, OutputArray dst, InputArray kernel,
Point anchor = Point(-1,-1), int iterations = 1,
BorderTypes borderType = BORDER_CONSTANT,
const Scalar& borderValue = morphologyDefaultBorderValue() );

/**
* @brief Faster implementation of cv::dilate with sparse table concept.
*
* @param src input image; the number of channels can be arbitrary, but the depth should be one of
* CV_8U, CV_16U, CV_16S, CV_32F or CV_64F.
* @param dst output image of the same size and type as src.
* @param kernel structuring element used for dilation; if element=Mat(), a 3 x 3 rectangular
* structuring element is used. Kernel can be created using #getStructuringElement
* @param anchor position of the anchor within the element; default value (-1, -1) means that the
* anchor is at the element center.
* @param iterations number of times dilation is applied.
* @param borderType pixel extrapolation method, see #BorderTypes. #BORDER_WRAP is not suported.
* @param borderValue border value in case of a constant border
*
* @see cv::dilate
*/
CV_EXPORTS void dilate( InputArray src, OutputArray dst, InputArray kernel,
Point anchor = Point(-1,-1), int iterations = 1,
BorderTypes borderType = BORDER_CONSTANT,
const Scalar& borderValue = morphologyDefaultBorderValue() );

/**
* @brief Faster implementation of cv::morphologyEx with sparse table concept.

* @param src Source image. The number of channels can be arbitrary. The depth should be one of
* CV_8U, CV_16U, CV_16S, CV_32F or CV_64F.
* @param dst Destination image of the same size and type as source image.
* @param op Type of a morphological operation, see #MorphTypes
* @param kernel Structuring element. It can be created using #getStructuringElement.
* @param anchor Anchor position with the kernel. Negative values mean that the anchor is at the
* kernel center.
* @param iterations Number of times erosion and dilation are applied.
* @param borderType Pixel extrapolation method, see #BorderTypes. #BORDER_WRAP is not supported.
* @param borderValue Border value in case of a constant border. The default value has a special
* meaning.
* @note The number of iterations is the number of times erosion or dilatation operation will be applied.
* For instance, an opening operation (#MORPH_OPEN) with two iterations is equivalent to apply
* successively: erode -> erode -> dilate -> dilate (and not erode -> dilate -> erode -> dilate).
*
* @see cv::morphologyEx
*/
CV_EXPORTS void morphologyEx( InputArray src, OutputArray dst,
int op, InputArray kernel,
Point anchor = Point(-1,-1), int iterations = 1,
BorderTypes borderType = BORDER_CONSTANT,
const Scalar& borderValue = morphologyDefaultBorderValue() );
//! @}

}}} // cv::ximgproc::stMorph::
#endif
40 changes: 40 additions & 0 deletions modules/ximgproc/perf/perf_sparse_table_morphology.cpp
Original file line number Diff line number Diff line change
@@ -0,0 +1,40 @@
// 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 "perf_precomp.hpp"

namespace opencv_test {
namespace {

typedef tuple<MorphTypes, MorphShapes> MorphTypes_MorphShapes_t;
typedef TestBaseWithParam<MorphTypes_MorphShapes_t> SparseTableMorphologyPerfTest;

PERF_TEST_P(SparseTableMorphologyPerfTest, perf,
testing::Combine(
testing::Values(
MORPH_ERODE, MORPH_DILATE, MORPH_OPEN, MORPH_CLOSE,
MORPH_GRADIENT, MORPH_TOPHAT, MORPH_BLACKHAT),
testing::Values(MORPH_RECT, MORPH_CROSS, MORPH_ELLIPSE)
) )
{
MorphTypes_MorphShapes_t params = GetParam();
int seSize = 51;
Size sz = sz1080p;
MorphTypes op = std::get<0>(params);
MorphShapes knType = std::get<1>(params);

Mat src(sz, CV_8UC3), dst(sz, CV_8UC3);
Mat kernel = getStructuringElement(knType, cv::Size(2 * seSize + 1, 2 * seSize + 1));

declare.in(src, WARMUP_RNG).out(dst);

TEST_CYCLE_N(5)
{
cv::ximgproc::stMorph::morphologyEx(src, dst, op, kernel);
}

SANITY_CHECK_NOTHING();
}

}} // opencv_test:: ::
429 changes: 429 additions & 0 deletions modules/ximgproc/src/sparse_table_morphology.cpp
Original file line number Diff line number Diff line change
@@ -0,0 +1,429 @@
// 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 <limits>
#include <utility>
#include <vector>

namespace cv {
namespace ximgproc {
namespace stMorph {

// normalizeAnchor; Copied from filterengine.hpp.
static inline Point normalizeAnchor(Point anchor, Size ksize)
{
if (anchor.x == -1)
anchor.x = ksize.width / 2;
if (anchor.y == -1)
anchor.y = ksize.height / 2;
CV_Assert(anchor.inside(Rect(0, 0, ksize.width, ksize.height)));
return anchor;
}

static int log2(int n)
{
int ans = -1;
while (n > 0)
{
n /= 2;
ans++;
}
return ans;
}

static int longestRowRunLength(const Mat& kernel)
{
int cnt = 0;
int maxLen = 0;
for (int c = 0; c < kernel.cols; c++)
{
cnt = 0;
for (int r = 0; r < kernel.rows; r++)
{
if (kernel.at<uchar>(r, c) == 0)
{
maxLen = std::max(maxLen, cnt);
cnt = 0;
}
else cnt++;
}
maxLen = std::max(maxLen, cnt);
}
return maxLen;
}

static int longestColRunLength(const Mat& kernel)
{
int cnt = 0;
int maxLen = 0;
for (int r = 0; r < kernel.rows; r++)
{
cnt = 0;
for (int c = 0; c < kernel.cols; c++)
{
if (kernel.at<uchar>(r, c) == 0)
{
maxLen = std::max(maxLen, cnt);
cnt = 0;
}
else cnt++;
}
maxLen = std::max(maxLen, cnt);
}
return maxLen;
}

static std::vector<Point> findP2RectCorners(const Mat& stNode, int rowDepth, int colDepth)
{
int rowOfst = 1 << rowDepth;
int colOfst = 1 << colDepth;
std::vector<Point> p2Rects;
for (int row = 0; row < stNode.rows; row++)
{
for (int col = 0; col < stNode.cols; col++)
{
// select white cells
if (stNode.at<uchar>(row, col) == 0) continue;

// select corner cells
if (col > 0 && stNode.at<uchar>(row, col - 1) == 1
&& col + 1 < stNode.cols && stNode.at<uchar>(row, col + 1) == 1) continue;
if (row > 0 && stNode.at<uchar>(row - 1, col) == 1
&& row + 1 < stNode.rows && stNode.at<uchar>(row + 1, col) == 1) continue;

// ignore if deeper cell is white
if (col + colOfst < stNode.cols && stNode.at<uchar>(row, col + colOfst) == 1) continue;
if (col - colOfst >= 0 && stNode.at<uchar>(row, col - colOfst) == 1) continue;
if (row + rowOfst < stNode.rows && stNode.at<uchar>(row + rowOfst, col) == 1) continue;
if (row - rowOfst >= 0 && stNode.at<uchar>(row - rowOfst, col) == 1) continue;

p2Rects.emplace_back(col, row);
}
}
return p2Rects;
}

/*
* Find a set of power-2-rectangles to cover the kernel.
* power-2-rectangles is a rectangle whose height and width are both power of 2.
*/
static std::vector<std::vector<std::vector<Point>>> genPow2RectsToCoverKernel(
const Mat& kernel, int rowDepthLim, int colDepthLim)
{
CV_Assert(kernel.type() == CV_8UC1);

std::vector<std::vector<std::vector<Point>>> p2Rects;
Mat stNodeCache = kernel;
for (int rowDepth = 0; rowDepth < rowDepthLim; rowDepth++)
{
Mat stNode = stNodeCache.clone();
p2Rects.emplace_back(std::vector<std::vector<Point>>());
for (int colDepth = 0; colDepth < colDepthLim; colDepth++)
{
p2Rects[rowDepth].emplace_back(findP2RectCorners(stNode, rowDepth, colDepth));
int colStep = 1 << colDepth;
if (stNode.cols - colStep < 0) break;
Rect s1(0, 0, stNode.cols - colStep, stNode.rows);
Rect s2(colStep, 0, stNode.cols - colStep, stNode.rows);
cv::min(stNode(s1), stNode(s2), stNode);
}
int rowStep = 1 << rowDepth;
if (stNodeCache.rows - rowStep < 0) break;
Rect s1(0, 0, stNodeCache.cols, stNodeCache.rows - rowStep);
Rect s2(0, rowStep, stNodeCache.cols, stNodeCache.rows - rowStep);
cv::min(stNodeCache(s1), stNodeCache(s2), stNodeCache);
}

return p2Rects;
}

/*
* Solves the rectilinear steiner arborescence problem greedy.
*/
static Mat SolveRSAPGreedy(const Mat& initialMap)
{
CV_Assert(initialMap.type() == CV_8UC1);
std::vector<Point> pos;
for (int r = 0; r < initialMap.rows; r++)
for (int c = 0; c < initialMap.cols; c++)
if (initialMap.at<uchar>(r, c) == 1) pos.emplace_back(c, r);
Mat resMap = Mat::zeros(initialMap.size(), CV_8UC2);

while (pos.size() > 1)
{
int maxCost = -1;
int maxI = 0;
int maxJ = 0;
int maxX = 0;
int maxY = 0;
for (uint i = 0; i < pos.size(); i++)
{
for (uint j = i + 1; j < pos.size(); j++)
{
int _x = std::min(pos[i].x, pos[j].x);
int _y = std::min(pos[i].y, pos[j].y);
int cost = _x + _y;
if (maxCost < cost)
{
maxCost = cost;
maxI = i;
maxJ = j;
maxX = _x;
maxY = _y;
}
}
}
for (int col = pos[maxI].x - 1; col >= maxX; col--)
resMap.at<Vec2b>(pos[maxI].y, col)[1] = 1;
for (int row = pos[maxI].y - 1; row >= maxY; row--)
resMap.at<Vec2b>(row, maxX)[0] = 1;
for (int col = pos[maxJ].x - 1; col >= maxX; col--)
resMap.at<Vec2b>(pos[maxJ].y, col)[1] = 1;
for (int row = pos[maxJ].y - 1; row >= maxY; row--)
resMap.at<Vec2b>(row, maxX)[0] = 1;

pos[maxI] = Point(maxX, maxY);
std::swap(pos[maxJ], pos[pos.size() - 1]);
pos.pop_back();
}
return resMap;
}

static Mat sparseTableFillPlanning(
std::vector<std::vector<std::vector<Point>>> pow2Rects, int rowDepthLim, int colDepthLim)
{
// list up required sparse table nodes.
Mat stMap = Mat::zeros(rowDepthLim, colDepthLim, CV_8UC1);
for (int rd = 0; rd < rowDepthLim; rd++)
for (int cd = 0; cd < colDepthLim; cd++)
if (pow2Rects[rd][cd].size() > 0)
stMap.at<uchar>(rd, cd) = 1;
stMap.at<uchar>(0, 0) = 1;
Mat path = SolveRSAPGreedy(stMap);
return path;
}

kernelDecompInfo decompKernel(InputArray kernel, Point anchor, int iterations)
{
Mat _kernel = kernel.getMat();
// Fix kernel in case of it is empty.
if (_kernel.empty())
{
_kernel = getStructuringElement(MORPH_RECT, Size(1 + iterations * 2, 1 + iterations * 2));
anchor = Point(iterations, iterations);
iterations = 1;
}
if (countNonZero(_kernel) == 0)
{
_kernel.at<uchar>(0, 0) = 1;
}

int rowDepthLim = log2(longestRowRunLength(_kernel)) + 1;
int colDepthLim = log2(longestColRunLength(_kernel)) + 1;
std::vector<std::vector<std::vector<Point>>> pow2Rects
= genPow2RectsToCoverKernel(_kernel, rowDepthLim, colDepthLim);

Mat stPlan
= sparseTableFillPlanning(pow2Rects, rowDepthLim, colDepthLim);

// Fix anchor to the center of the kernel.
anchor = stMorph::normalizeAnchor(anchor, _kernel.size());

return { _kernel.rows, _kernel.cols, pow2Rects, stPlan, anchor, iterations };
}

enum Op
{
Min, Max
};

static void morphDfs(int minmax, Mat& st, Mat& dst,
std::vector<std::vector<std::vector<Point>>> row2Rects, const Mat& stPlan,
int rowDepth, int colDepth)
{
for (Point p : row2Rects[rowDepth][colDepth])
{
Rect rect(p, dst.size());
if (minmax == Op::Min) cv::min(dst, st(rect), dst);
else cv::max(dst, st(rect), dst);
}

if (stPlan.at<Vec2b>(rowDepth, colDepth)[1] == 1)
{
// col direction
Mat st2 = st;
int ofs = 1 << colDepth;
Rect rect1(0, 0, st2.cols - ofs, st2.rows);
Rect rect2(ofs, 0, st2.cols - ofs, st2.rows);

if (minmax == Op::Min) cv::min(st2(rect1), st2(rect2), st2);
else cv::max(st2(rect1), st2(rect2), st2);
morphDfs(minmax, st2, dst, row2Rects, stPlan, rowDepth, colDepth + 1);
}
if (stPlan.at<Vec2b>(rowDepth, colDepth)[0] == 1)
{
// row direction
int ofs = 1 << rowDepth;
Rect rect1(0, 0, st.cols, st.rows - ofs);
Rect rect2(0, ofs, st.cols, st.rows - ofs);

if (minmax == Op::Min) cv::min(st(rect1), st(rect2), st);
else cv::max(st(rect1), st(rect2), st);
morphDfs(minmax, st, dst, row2Rects, stPlan, rowDepth + 1, colDepth);
}
}

template <typename T>
static void morphOp(Op minmax, InputArray _src, OutputArray _dst, kernelDecompInfo kdi,
BorderTypes borderType, const Scalar& borderVal)
{
T nil = (minmax == Op::Min) ? std::numeric_limits<T>::max() : std::numeric_limits<T>::min();

Mat src = _src.getMat();
_dst.create(_src.size(), _src.type());
Mat dst = _dst.getMat();

Scalar bV = borderVal;
if (borderType == BorderTypes::BORDER_CONSTANT && borderVal == morphologyDefaultBorderValue())
bV = Scalar::all(nil);

do
{
Mat expandedSrc;
copyMakeBorder(src, expandedSrc,
kdi.anchor.y, kdi.cols - 1 - kdi.anchor.y,
kdi.anchor.x, kdi.rows - 1 - kdi.anchor.x,
borderType, bV);
dst.setTo(nil);
morphDfs(minmax, expandedSrc, dst, kdi.stRects, kdi.plan, 0, 0);
src = dst;
} while (--kdi.iterations > 0);
}

static void morphOp(Op minmax, InputArray _src, OutputArray _dst, kernelDecompInfo kdi,
BorderTypes borderType, const Scalar& borderVal)
{
if (kdi.iterations == 0 || kdi.rows * kdi.cols == 1)
{
_src.copyTo(_dst);
return;
}

switch (_src.depth())
{
case CV_8U:
morphOp<uchar>(minmax, _src, _dst, kdi, borderType, borderVal);
return;
case CV_8S:
morphOp<char>(minmax, _src, _dst, kdi, borderType, borderVal);
return;
case CV_16U:
morphOp<ushort>(minmax, _src, _dst, kdi, borderType, borderVal);
return;
case CV_16S:
morphOp<short>(minmax, _src, _dst, kdi, borderType, borderVal);
return;
case CV_32S:
morphOp<int>(minmax, _src, _dst, kdi, borderType, borderVal);
return;
case CV_32F:
morphOp<float>(minmax, _src, _dst, kdi, borderType, borderVal);
return;
case CV_64F:
morphOp<double>(minmax, _src, _dst, kdi, borderType, borderVal);
return;
}
}

void erode(InputArray src, OutputArray dst, kernelDecompInfo kdi,
BorderTypes borderType, const Scalar& borderVal)
{
morphOp(Op::Min, src, dst, kdi, borderType, borderVal);
}

void dilate(InputArray src, OutputArray dst, kernelDecompInfo kdi,
BorderTypes borderType, const Scalar& borderVal)
{
morphOp(Op::Max, src, dst, kdi, borderType, borderVal);
}

void morphologyEx(InputArray src, OutputArray dst, int op, kernelDecompInfo kdi,
BorderTypes borderType, const Scalar& borderVal)
{
CV_INSTRUMENT_REGION();

CV_Assert(!src.empty());

Mat _src = src.getMat(), temp;
dst.create(_src.size(), _src.type());
Mat _dst = dst.getMat();

switch (op)
{
case MORPH_ERODE:
erode(src, dst, kdi, borderType, borderVal);
break;
case MORPH_DILATE:
dilate(src, dst, kdi, borderType, borderVal);
break;
case MORPH_OPEN:
stMorph::erode(src, dst, kdi, borderType, borderVal);
stMorph::dilate(dst, dst, kdi, borderType, borderVal);
break;
case MORPH_CLOSE:
stMorph::dilate(src, dst, kdi, borderType, borderVal);
stMorph::erode(dst, dst, kdi, borderType, borderVal);
break;
case MORPH_GRADIENT:
stMorph::erode(_src, temp, kdi, borderType, borderVal);
stMorph::dilate(_src, _dst, kdi, borderType, borderVal);
_dst -= temp;
break;
case MORPH_TOPHAT:
if (_src.data != _dst.data)
temp = _dst;
stMorph::erode(_src, temp, kdi, borderType, borderVal);
stMorph::dilate(temp, temp, kdi, borderType, borderVal);
_dst = _src - temp;
break;
case MORPH_BLACKHAT:
if (_src.data != _dst.data)
temp = _dst;
stMorph::dilate(_src, temp, kdi, borderType, borderVal);
stMorph::erode(temp, temp, kdi, borderType, borderVal);
_dst = temp - _src;
break;
case MORPH_HITMISS:
CV_Error(cv::Error::StsBadArg, "HIT-MISS operation is not supported.");
default:
CV_Error(cv::Error::StsBadArg, "Unknown morphological operation.");
}
}

void erode(InputArray src, OutputArray dst, InputArray kernel,
Point anchor, int iterations,
BorderTypes borderType, const Scalar& borderVal)
{
kernelDecompInfo kdi = decompKernel(kernel, anchor, iterations);
morphOp(Op::Min, src, dst, kdi, borderType, borderVal);
}

void dilate(InputArray src, OutputArray dst, InputArray kernel,
Point anchor, int iterations,
BorderTypes borderType, const Scalar& borderVal)
{
kernelDecompInfo kdi = decompKernel(kernel, anchor, iterations);
morphOp(Op::Max, src, dst, kdi, borderType, borderVal);
}

void morphologyEx(InputArray src, OutputArray dst, int op,
InputArray kernel, Point anchor, int iterations,
BorderTypes borderType, const Scalar& borderVal)
{
kernelDecompInfo kdi = decompKernel(kernel, anchor, iterations);
morphologyEx(src, dst, op, kdi, borderType, borderVal);
}

}}} // cv::ximgproc::stMorph::
496 changes: 496 additions & 0 deletions modules/ximgproc/test/test_sparse_table_morphology.cpp
Original file line number Diff line number Diff line change
@@ -0,0 +1,496 @@
// 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 "test_precomp.hpp"
#include "opencv2/ximgproc/sparse_table_morphology.hpp"
#include <vector>

namespace opencv_test {
namespace {

void assertArraysIdentical(InputArray ary1, InputArray ary2)
{
Mat xormat = ary1.getMat() ^ ary2.getMat();
CV_Assert(cv::countNonZero(xormat.reshape(1)) == 0);
}
Mat im(int type)
{
int depth = CV_MAT_DEPTH(type);
int ch = CV_MAT_CN(type);
Mat img = imread(cvtest::TS::ptr()->get_data_path() + "cv/shared/lena.png");
CV_Assert(img.type() == CV_8UC3);

if (ch == 1) cv::cvtColor(img, img, ColorConversionCodes::COLOR_BGR2GRAY, ch);
if (depth == CV_8S) img /= 2;
img.convertTo(img, depth);
if (depth == CV_16S) img *= (1 << 7);
if (depth == CV_16U) img *= (1 << 8);
if (depth == CV_32S) img *= (1 << 23);
if (depth == CV_32F) img /= (1 << 8);
if (depth == CV_64F) img /= (1 << 8);

return img;
}
Mat kn4() { return getStructuringElement(cv::MorphShapes::MORPH_ELLIPSE, Size(4, 4)); }
Mat kn5() { return getStructuringElement(cv::MorphShapes::MORPH_ELLIPSE, Size(5, 5)); }
Mat knBig() { return getStructuringElement(cv::MorphShapes::MORPH_RECT, Size(201, 201)); }
Mat kn1Zero() { return Mat::zeros(1, 1, CV_8UC1); }
Mat kn1One() { return Mat::ones(1, 1, CV_8UC1); }
Mat knEmpty() { return Mat(); }
Mat knZeros() { return Mat::zeros(5, 5, CV_8UC1); }
Mat knOnes() { return Mat::ones(5, 5, CV_8UC1); }
Mat knAsymm (){ return (Mat_<uchar>(5, 5) << 0,0,0,0,0,0,0,1,0,0,0,1,0,0,0,0,0,0,0,0,0,0,1,0,0); }
Mat knRnd(int size, int density)
{
Mat rndMat(size, size, CV_8UC1);
theRNG().state = getTickCount();
randu(rndMat, 2, 102);
density++;
rndMat.setTo(0, density < rndMat);
rndMat.setTo(1, 1 < rndMat);
return rndMat;
}

/*
* dilate regression tests.
*/
void dilate_rgr(InputArray src, InputArray kernel, Point anchor = Point(-1, -1),
int iterations = 1,
BorderTypes bdrType = BorderTypes::BORDER_CONSTANT,
const Scalar& bdrVal = morphologyDefaultBorderValue())
{
Mat expected, actual;
dilate(src, expected, kernel, anchor, iterations, bdrType, bdrVal);
stMorph::dilate(src, actual, kernel, anchor, iterations, bdrType, bdrVal);
assertArraysIdentical(expected, actual);
}
TEST(ximgproc_StMorph_dilate, regression_8UC1) { dilate_rgr(im(CV_8UC1), kn5()); }
TEST(ximgproc_StMorph_dilate, regression_8UC3) { dilate_rgr(im(CV_8UC3), kn5()); }
TEST(ximgproc_StMorph_dilate, regression_16UC1) { dilate_rgr(im(CV_16UC1), kn5()); }
TEST(ximgproc_StMorph_dilate, regression_16UC3) { dilate_rgr(im(CV_16UC3), kn5()); }
TEST(ximgproc_StMorph_dilate, regression_16SC1) { dilate_rgr(im(CV_16SC1), kn5()); }
TEST(ximgproc_StMorph_dilate, regression_16SC3) { dilate_rgr(im(CV_16SC3), kn5()); }
TEST(ximgproc_StMorph_dilate, regression_32FC1) { dilate_rgr(im(CV_32FC1), kn5()); }
TEST(ximgproc_StMorph_dilate, regression_32FC3) { dilate_rgr(im(CV_32FC3), kn5()); }
TEST(ximgproc_StMorph_dilate, regression_64FC1) { dilate_rgr(im(CV_64FC1), kn5()); }
TEST(ximgproc_StMorph_dilate, regression_64FC3) { dilate_rgr(im(CV_64FC3), kn5()); }
TEST(ximgproc_StMorph_dilate, regression_kn5) { dilate_rgr(im(CV_8UC3), kn5()); }
TEST(ximgproc_StMorph_dilate, regression_kn4) { dilate_rgr(im(CV_8UC3), kn4()); }
TEST(ximgproc_StMorph_dilate, wtf_regression_kn1Zero) { dilate_rgr(im(CV_8UC3), kn1Zero()); }
TEST(ximgproc_StMorph_dilate, regression_kn1One) { dilate_rgr(im(CV_8UC3), kn1One()); }
TEST(ximgproc_StMorph_dilate, wtf_regression_knEmpty) { dilate_rgr(im(CV_8UC3), knEmpty()); }
TEST(ximgproc_StMorph_dilate, wtf_regression_knZeros) { dilate_rgr(im(CV_8UC3), knZeros()); }
TEST(ximgproc_StMorph_dilate, regression_knOnes) { dilate_rgr(im(CV_8UC3), knOnes()); }
TEST(ximgproc_StMorph_dilate, regression_knBig) { dilate_rgr(im(CV_8UC3), knBig()); }
TEST(ximgproc_StMorph_dilate, regression_knAsymm) { dilate_rgr(im(CV_8UC3), knAsymm()); }
TEST(ximgproc_StMorph_dilate, regression_ancMid) { dilate_rgr(im(CV_8UC3), kn5(), Point(-1, -1)); }
TEST(ximgproc_StMorph_dilate, regression_ancEdge1) { dilate_rgr(im(CV_8UC3), kn5(), Point(0, 0)); }
TEST(ximgproc_StMorph_dilate, regression_ancEdge2) { dilate_rgr(im(CV_8UC3), kn5(), Point(4, 4)); }
TEST(ximgproc_StMorph_dilate, wtf_regression_it0) { dilate_rgr(im(CV_8UC3), kn5(), Point(-1, -1), 0); }
TEST(ximgproc_StMorph_dilate, regression_it1) { dilate_rgr(im(CV_8UC3), kn5(), Point(-1, -1), 1); }
TEST(ximgproc_StMorph_dilate, regression_it2) { dilate_rgr(im(CV_8UC3), kn5(), Point(-1, -1), 2); }
/*
* dilate feature tests.
*/
void dilate_ftr(InputArray src, InputArray kernel, Point anchor = Point(-1, -1),
int iterations = 1,
BorderTypes bdrType = BorderTypes::BORDER_CONSTANT,
const Scalar& bdrVal = morphologyDefaultBorderValue())
{
Mat expected, actual;
stMorph::dilate(src, actual, kernel, anchor, iterations, bdrType, bdrVal);
// todo: generate expected result.
// assertArraysIdentical(expected, actual);
}
/* CV_8S, CV_16F are not supported by morph.simd::getMorphologyFilter */
TEST(ximgproc_StMorph_dilate, feature_8SC1) { dilate_ftr(im(CV_8SC1), kn5()); }
TEST(ximgproc_StMorph_dilate, feature_8SC3) { dilate_ftr(im(CV_8SC3), kn5()); }
TEST(ximgproc_StMorph_dilate, feature_32SC1) { dilate_ftr(im(CV_32SC1), kn5()); }
TEST(ximgproc_StMorph_dilate, feature_32SC3) { dilate_ftr(im(CV_32SC3), kn5()); }

/*
* erode regression tests.
*/
void erode_rgr(InputArray src, InputArray kernel, Point anchor = Point(-1, -1),
int iterations = 1,
BorderTypes bdrType = BorderTypes::BORDER_CONSTANT,
const Scalar& bdrVal = morphologyDefaultBorderValue())
{
Mat expected, actual;
erode(src, expected, kernel, anchor, iterations, bdrType, bdrVal);
stMorph::erode(src, actual, kernel, anchor, iterations, bdrType, bdrVal);
assertArraysIdentical(expected, actual);
}
TEST(ximgproc_StMorph_erode, regression_8UC1) { erode_rgr(im(CV_8UC1), kn5()); }
TEST(ximgproc_StMorph_erode, regression_8UC3) { erode_rgr(im(CV_8UC3), kn5()); }
TEST(ximgproc_StMorph_erode, regression_16UC1) { erode_rgr(im(CV_16UC1), kn5()); }
TEST(ximgproc_StMorph_erode, regression_16UC3) { erode_rgr(im(CV_16UC3), kn5()); }
TEST(ximgproc_StMorph_erode, regression_16SC1) { erode_rgr(im(CV_16SC1), kn5()); }
TEST(ximgproc_StMorph_erode, regression_16SC3) { erode_rgr(im(CV_16SC3), kn5()); }
TEST(ximgproc_StMorph_erode, regression_32FC1) { erode_rgr(im(CV_32FC1), kn5()); }
TEST(ximgproc_StMorph_erode, regression_32FC3) { erode_rgr(im(CV_32FC3), kn5()); }
TEST(ximgproc_StMorph_erode, regression_64FC1) { erode_rgr(im(CV_64FC1), kn5()); }
TEST(ximgproc_StMorph_erode, regression_64FC3) { erode_rgr(im(CV_64FC3), kn5()); }
TEST(ximgproc_StMorph_erode, regression_kn5) { erode_rgr(im(CV_8UC3), kn5()); }
TEST(ximgproc_StMorph_erode, regression_kn4) { erode_rgr(im(CV_8UC3), kn4()); }
TEST(ximgproc_StMorph_erode, wtf_regression_kn1Zero) { erode_rgr(im(CV_8UC3), kn1Zero()); }
TEST(ximgproc_StMorph_erode, regression_kn1One) { erode_rgr(im(CV_8UC3), kn1One()); }
TEST(ximgproc_StMorph_erode, wtf_regression_knEmpty) { erode_rgr(im(CV_8UC3), knEmpty()); }
TEST(ximgproc_StMorph_erode, wtf_regression_knZeros) { erode_rgr(im(CV_8UC3), knZeros()); }
TEST(ximgproc_StMorph_erode, regression_knOnes) { erode_rgr(im(CV_8UC3), knOnes()); }
TEST(ximgproc_StMorph_erode, regression_knBig) { erode_rgr(im(CV_8UC3), knBig()); }
TEST(ximgproc_StMorph_erode, regression_knAsymm) { erode_rgr(im(CV_8UC3), knAsymm()); }
TEST(ximgproc_StMorph_erode, regression_ancMid) { erode_rgr(im(CV_8UC3), kn5(), Point(-1, -1)); }
TEST(ximgproc_StMorph_erode, regression_ancEdge1) { erode_rgr(im(CV_8UC3), kn5(), Point(0, 0)); }
TEST(ximgproc_StMorph_erode, regression_ancEdge2) { erode_rgr(im(CV_8UC3), kn5(), Point(4, 4)); }
TEST(ximgproc_StMorph_erode, wtf_regression_it0) { erode_rgr(im(CV_8UC3), kn5(), Point(-1, -1), 0); }
TEST(ximgproc_StMorph_erode, regression_it1) { erode_rgr(im(CV_8UC3), kn5(), Point(-1, -1), 1); }
TEST(ximgproc_StMorph_erode, regression_it2) { erode_rgr(im(CV_8UC3), kn5(), Point(-1, -1), 2); }
/*
* erode feature tests.
*/
void erode_ftr(InputArray src, InputArray kernel, Point anchor = Point(-1, -1),
int iterations = 1,
BorderTypes bdrType = BorderTypes::BORDER_CONSTANT,
const Scalar& bdrVal = morphologyDefaultBorderValue())
{
Mat expected, actual;
stMorph::erode(src, actual, kernel, anchor, iterations, bdrType, bdrVal);
// todo: generate expected result.
// assertArraysIdentical(expected, actual);
}
/* CV_8S, CV_16F are not supported by morph.simd::getMorphologyFilter */
TEST(ximgproc_StMorph_erode, feature_8SC1) { erode_ftr(im(CV_8SC1), kn5()); }
TEST(ximgproc_StMorph_erode, feature_8SC3) { erode_ftr(im(CV_8SC3), kn5()); }
TEST(ximgproc_StMorph_erode, feature_32SC1) { erode_ftr(im(CV_32SC1), kn5()); }
TEST(ximgproc_StMorph_erode, feature_32SC3) { erode_ftr(im(CV_32SC3), kn5()); }

/*
* morphologyEx regression tests.
*/
void ex_rgr(InputArray src, MorphTypes op, InputArray kernel, Point anchor = Point(-1, -1),
int iterations = 1,
BorderTypes bdrType = BorderTypes::BORDER_CONSTANT,
const Scalar& bdrVal = morphologyDefaultBorderValue())
{
Mat expected, actual;
morphologyEx(src, expected, op, kernel, anchor, iterations, bdrType, bdrVal);
stMorph::morphologyEx(src, actual, op, kernel, anchor, iterations, bdrType, bdrVal);
assertArraysIdentical(expected, actual);
}
TEST(ximgproc_StMorph_ex, regression_erode) { ex_rgr(im(CV_8UC3), MORPH_ERODE, kn5()); }
TEST(ximgproc_StMorph_ex, regression_dilage) { ex_rgr(im(CV_8UC3), MORPH_DILATE, kn5()); }
TEST(ximgproc_StMorph_ex, regression_open) { ex_rgr(im(CV_8UC3), MORPH_OPEN, kn5()); }
TEST(ximgproc_StMorph_ex, regression_close) { ex_rgr(im(CV_8UC3), MORPH_CLOSE, kn5()); }
TEST(ximgproc_StMorph_ex, regression_gradient) { ex_rgr(im(CV_8UC3), MORPH_GRADIENT, kn5()); }
TEST(ximgproc_StMorph_ex, regression_tophat) { ex_rgr(im(CV_8UC3), MORPH_TOPHAT, kn5()); }
TEST(ximgproc_StMorph_ex, regression_blackhat) { ex_rgr(im(CV_8UC3), MORPH_BLACKHAT, kn5()); }
TEST(ximgproc_StMorph_ex, regression_hitmiss)
{
EXPECT_THROW( { ex_rgr(im(CV_8UC1), MORPH_HITMISS, kn5()); }, cv::Exception);
}

stMorph::kernelDecompInfo ftr_decomp(InputArray kernel)
{
auto kdi = stMorph::decompKernel(kernel);
Mat expected = kernel.getMat();
Mat actual = Mat::zeros(kernel.size(), kernel.type());
for (uint r = 0; r < kdi.stRects.size(); r++)
{
for (uint c = 0; c < kdi.stRects[r].size(); c++)
{
for (Point p : kdi.stRects[r][c])
{
Rect rect(p.x, p.y, 1 << c, 1 << r);
actual(rect).setTo(1);
}
}
}
assertArraysIdentical(expected, actual);
return kdi;
}
Mat VisualizeCovering(Mat& kernel, const stMorph::kernelDecompInfo& kdi)
{
const int rate = 20;
const int fluct = 3;
const int colors = 20;
resize(kernel * 255, kernel, Size(), rate, rate, InterpolationFlags::INTER_NEAREST);
cvtColor(kernel, kernel, cv::COLOR_GRAY2BGR);
Scalar color[colors]{
Scalar(255, 127, 127), Scalar(255, 127, 191), Scalar(255, 127, 255), Scalar(191, 127, 255),
Scalar(127, 127, 255), Scalar(127, 191, 255), Scalar(127, 255, 255), Scalar(127, 255, 191),
Scalar(127, 255, 127), Scalar(191, 255, 127), Scalar(255, 255, 127), Scalar(255, 191, 127)
};
int i = 0;
for (int r = 0; r < kdi.rows; r++)
cv::line(kernel, Point(0, r * rate), Point(kdi.cols * rate, r * rate), Scalar(0));
for (int c = 0; c < kdi.cols; c++)
cv::line(kernel, Point(c * rate, 0), Point(c * rate, kdi.rows * rate), Scalar(0));
for (uint row = 0; row < kdi.stRects.size(); row++)
{
for (uint col = 0; col < kdi.stRects[row].size(); col++)
{
Size s(1 << col, 1 << row);
for (Point p : kdi.stRects[row][col])
{
Rect rect(p, s);
int l = (rect.x) * rate + i % fluct + 2;
int t = (rect.y) * rate + i % fluct + 2;
int r = (rect.x + rect.width) * rate - fluct + i % fluct - 1;
int b = (rect.y + rect.height) * rate - fluct + i % fluct - 1;
Point lt(l, t);
Point lb(l, b);
Point rb(r, b);
Point rt(r, t);
cv::line(kernel, lt, lb, color[i % colors], 1);
cv::line(kernel, lb, rb, color[i % colors], 1);
cv::line(kernel, rb, rt, color[i % colors], 1);
cv::line(kernel, rt, lt, color[i % colors], 1);
i++;
}
}
}
return kernel;
}
Mat VisualizePlanning(stMorph::kernelDecompInfo kdi)
{
int g = 30;
int rows = kdi.plan.rows;
int cols = kdi.plan.cols;
Scalar vCol = Scalar(20, 20, 255);
Scalar eCol = Scalar(100, 100, 100);
Mat m = Mat::zeros(rows * g, cols * g, CV_8UC3);
for (int row = 0; row < rows; row++)
{
for (int col = 0; col < cols; col++)
{
Rect nodeRect(col * g + g / 2 - 5, row * g + g / 2 - 5, 11, 11);
if (kdi.stRects[row][col].size() > 0)
cv::rectangle(m, nodeRect, vCol, -1);
else
cv::rectangle(m, nodeRect, vCol, 1);
}
}
for (int r = 0; r < rows; r++)
{
for (int c = 0; c < cols; c++)
{
Vec2b p = kdi.plan.at<Vec2b>(r, c);
Point sp(c * g + g / 2, r * g + g / 2);
if (p[0] == 1)
{
Point ep = Point(c * g + g / 2, (r + 1) * g + g / 2);
cv::line(m, sp, ep, eCol, 2);
}
if (p[1] == 1)
{
Point ep = Point((c + 1) * g + g / 2, r * g + g / 2);
cv::line(m, sp, ep, eCol, 2);
}
}
}
return m;
}
TEST(ximgproc_StMorph_decomp, feature_rnd1) { ftr_decomp(knRnd(1000, 1)); }
TEST(ximgproc_StMorph_decomp, feature_rnd10) { ftr_decomp(knRnd(1000, 10)); }
TEST(ximgproc_StMorph_decomp, feature_rnd30) { ftr_decomp(knRnd(1000, 30)); }
TEST(ximgproc_StMorph_decomp, feature_rnd50) { ftr_decomp(knRnd(1000, 50)); }
TEST(ximgproc_StMorph_decomp, feature_rnd80) { ftr_decomp(knRnd(1000, 80)); }
TEST(ximgproc_StMorph_decomp, feature_rnd90) { ftr_decomp(knRnd(1000, 90)); }
TEST(ximgproc_StMorph_decomp, feature_visualize) {
Mat kernel = getStructuringElement(MORPH_ELLIPSE, Size(5, 5));
auto kdi = ftr_decomp(kernel);
Mat covering = VisualizeCovering(kernel, kdi);
Mat plan = VisualizePlanning(kdi);
#if 0
imshow("Covering", covering);
imshow("Plan", plan);
waitKey();
destroyAllWindows();
#endif
}

TEST(ximgproc_StMorph_eval, pdi)
{
Mat img = im(CV_8UC3);
Mat dst;
int sizes[]{ 3, 5, 7, 9, 11, 13, 15, 17, 19, 21, 23, 25, 27, 29, 31, 35, 41, 45, 51, 55, 61, 71,
81, 91, 101, 121, 151, 171, 201, 221, 251, 301, 351, 401, 451, 501 };

std::ofstream ss("opencvlog_pdi.txt", std::ios_base::out);

for (int c = 0; c < 3; c++)
for (int i: sizes)
{
ss << i;
Size sz(i, i);
cv::TickMeter meter;
Mat kn;
stMorph::kernelDecompInfo kdi;

// cv-rect
kn = getStructuringElement(MORPH_RECT, sz);
if (i <= 401)
{
meter.start();
cv::erode(img, dst, kn);
meter.stop();
ss << "\t" << meter.getTimeMilli();
meter.reset();
}
else
{
ss << "\t";
}

// cv-cross
kn = getStructuringElement(MORPH_CROSS, sz);
if (i <= 401)
{
meter.start();
cv::erode(img, dst, kn);
meter.stop();
ss << "\t" << meter.getTimeMilli();
meter.reset();
}
else
{
ss << "\t";
}

// cv-ellipse
kn = getStructuringElement(MORPH_ELLIPSE, sz);
if (i <= 23)
{
meter.start();
cv::erode(img, dst, kn);
meter.stop();
ss << "\t" << meter.getTimeMilli();
meter.reset();
}
else
{
ss << "\t";
}

// st-rect
kn = getStructuringElement(MORPH_RECT, sz);
kdi = stMorph::decompKernel(kn);
meter.start();
stMorph::erode(img, dst, kdi);
meter.stop();
ss << "\t" << meter.getTimeMilli();
meter.reset();

// st-cross
kn = getStructuringElement(MORPH_CROSS, sz);
kdi = stMorph::decompKernel(kn);
meter.start();
stMorph::erode(img, dst, kdi);
meter.stop();
ss << "\t" << meter.getTimeMilli();
meter.reset();

// st-ellipse
kn = getStructuringElement(MORPH_ELLIPSE, sz);
kdi = stMorph::decompKernel(kn);
meter.start();
stMorph::erode(img, dst, kdi);
meter.stop();
ss << "\t" << meter.getTimeMilli() << "\n";
meter.reset();
}
ss.close();
}

TEST(ximgproc_StMorph_eval, integrated)
{
Mat img = im(CV_8UC3);
Mat dst;
int sizes[]{ 3, 5, 7, 9, 11, 13, 15, 17, 19, 21, 23, 25, 27, 29, 31, 35, 41, 45, 51, 55, 61, 71,
81, 91, 101, 121, 151, 171, 201, 221, 251, 301, 351, 401, 451, 501 };

std::ofstream ss("opencvlog_integrated.txt", std::ios_base::out);

for (int c = 0; c < 3; c++)
for (int i: sizes)
{
ss << i;
Size sz(i, i);
cv::TickMeter meter;
Mat kn;

// cv-rect
kn = getStructuringElement(MORPH_RECT, sz);
if (i <= 401)
{
meter.start();
cv::erode(img, dst, kn);
meter.stop();
ss << "\t" << meter.getTimeMilli();
meter.reset();
}
else
{
ss << "\t";
}

// cv-cross
kn = getStructuringElement(MORPH_CROSS, sz);
if (i <= 401)
{
meter.start();
cv::erode(img, dst, kn);
meter.stop();
ss << "\t" << meter.getTimeMilli();
meter.reset();
}
else
{
ss << "\t";
}

// cv-ellipse
kn = getStructuringElement(MORPH_ELLIPSE, sz);
if (i <= 23)
{
meter.start();
cv::erode(img, dst, kn);
meter.stop();
ss << "\t" << meter.getTimeMilli();
meter.reset();
}
else
{
ss << "\t";
}

// st-rect
kn = getStructuringElement(MORPH_RECT, sz);
meter.start();
stMorph::erode(img, dst, kn);
meter.stop();
ss << "\t" << meter.getTimeMilli();
meter.reset();

// st-cross
kn = getStructuringElement(MORPH_CROSS, sz);
meter.start();
stMorph::erode(img, dst, kn);
meter.stop();
ss << "\t" << meter.getTimeMilli();
meter.reset();

// st-ellipse
kn = getStructuringElement(MORPH_ELLIPSE, sz);
meter.start();
stMorph::erode(img, dst, kn);
meter.stop();
ss << "\t" << meter.getTimeMilli() << "\n";
meter.reset();
}
ss.close();
}

}} // opencv_test:: ::