Chiebot-Mirror
/
opencv_contrib
8
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity

Merge pull request #1176 from vpisarev:anisodiff

Browse Source
pull/1050/merge
Vadim Pisarevsky 8 years ago
parent f78bb2c4e7 b885d28ca0
commit b74c25da7f
5 changed files with 483 additions and 0 deletions
Whitespace
Show all changes Ignore whitespace when comparing lines Ignore changes in amount of whitespace Ignore changes in whitespace at EOL
Split View
Diff Options
Show Stats Download Patch File Download Diff File
  1. 21
      modules/ximgproc/include/opencv2/ximgproc.hpp
  2. 105
      modules/ximgproc/samples/filterdemo.cpp
  3. 293
      modules/ximgproc/src/anisodiff.cpp
  4. 39
      modules/ximgproc/src/opencl/anisodiff.cl
  5. 25
      modules/ximgproc/test/test_anisodiff.cpp

21
modules/ximgproc/include/opencv2/ximgproc.hpp
Unescape View File

@ -122,6 +122,27 @@ The function transforms a binary blob image into a skeletized form using the tec
*/
CV_EXPORTS_W void thinning( InputArray src, OutputArray dst, int thinningType = THINNING_ZHANGSUEN);
/** @brief Performs anisotropic diffusian on an image.
The function applies Perona-Malik anisotropic diffusion to an image. This is the solution to the partial differential equation:
\f[{\frac {\partial I}{\partial t}}={\mathrm {div}}\left(c(x,y,t)\nabla I\right)=\nabla c\cdot \nabla I+c(x,y,t)\Delta I\f]
Suggested functions for c(x,y,t) are:
\f[c\left(\|\nabla I\|\right)=e^{{-\left(\|\nabla I\|/K\right)^{2}}}\f]
or
\f[ c\left(\|\nabla I\|\right)={\frac {1}{1+\left({\frac {\|\nabla I\|}{K}}\right)^{2}}} \f]
@param src Grayscale Source image.
@param dst Destination image of the same size and the same number of channels as src .
@param alpha The amount of time to step forward by on each iteration (normally, it's between 0 and 1).
@param K sensitivity to the edges
@param niters The number of iterations
*/
CV_EXPORTS_W void anisotropicDiffusion(InputArray src, OutputArray dst, float alpha, float K, int niters );
//! @}

105
modules/ximgproc/samples/filterdemo.cpp
Unescape View File

@ -0,0 +1,105 @@
/*M///////////////////////////////////////////////////////////////////////////////////////
//
// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
//
// By downloading, copying, installing or using the software you agree to this license.
// If you do not agree to this license, do not download, install,
// copy or use the software.
//
//
// License Agreement
// For Open Source Computer Vision Library
//
// Copyright (C) 2017, Intel Corporation, all rights reserved.
// Third party copyrights are property of their respective owners.
//
// Redistribution and use in source and binary forms, with or without modification,
// are permitted provided that the following conditions are met:
//
// * Redistribution's of source code must retain the above copyright notice,
// this list of conditions and the following disclaimer.
//
// * Redistribution's in binary form must reproduce the above copyright notice,
// this list of conditions and the following disclaimer in the documentation
// and/or other materials provided with the distribution.
//
// * The name of the copyright holders may not be used to endorse or promote products
// derived from this software without specific prior written permission.
//
// This software is provided by the copyright holders and contributors "as is" and
// any express or implied warranties, including, but not limited to, the implied
// warranties of merchantability and fitness for a particular purpose are disclaimed.
// In no event shall the Intel Corporation or contributors be liable for any direct,
// indirect, incidental, special, exemplary, or consequential damages
// (including, but not limited to, procurement of substitute goods or services;
// loss of use, data, or profits; or business interruption) however caused
// and on any theory of liability, whether in contract, strict liability,
// or tort (including negligence or otherwise) arising in any way out of
// the use of this software, even if advised of the possibility of such damage.
//
//M*/
#include "opencv2/core/utility.hpp"
#include "opencv2/imgproc.hpp"
#include "opencv2/imgcodecs.hpp"
#include "opencv2/highgui.hpp"
#include "opencv2/ximgproc.hpp"
#include <stdio.h>
using namespace cv;
using namespace std;
int main( int argc, const char** argv)
{
float alpha = 1.0f;
float sigma = 0.02f;
int rows0 = 480;
int niters = 10;
Mat frame, src, dst;
const char* window_name = "Anisodiff : Exponential Flux";
VideoCapture cap;
if( argc > 1 )
cap.open(argv[1]);
else
cap.open(0);
if (!cap.isOpened())
{
printf("Cannot initialize video capturing\n");
return 0;
}
// Create a window
namedWindow(window_name, 1);
// create a toolbar
createTrackbar("No. of time steps", window_name, &niters, 30, 0);
for(;;)
{
cap >> frame;
if( frame.empty() )
break;
if( frame.rows <= rows0 )
src = frame;
else
resize(frame, src, Size(cvRound(480.*frame.cols/frame.rows), 480));
float t = (float)getTickCount();
ximgproc::anisotropicDiffusion(src, dst, alpha, sigma, niters);
t = (float)getTickCount() - t;
printf("time: %.1fms\n", t*1000./getTickFrequency());
imshow(window_name, dst);
// Wait for a key stroke; the same function arranges events processing
char c = (char)waitKey(30);
if(c >= 0)
break;
}
return 0;
}

293
modules/ximgproc/src/anisodiff.cpp
Unescape View File

@ -0,0 +1,293 @@
/*M///////////////////////////////////////////////////////////////////////////////////////
//
// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
//
// By downloading, copying, installing or using the software you agree to this license.
// If you do not agree to this license, do not download, install,
// copy or use the software.
//
//
// License Agreement
// For Open Source Computer Vision Library
//
// Copyright (C) 2017, Intel Corporation, all rights reserved.
// Third party copyrights are property of their respective owners.
//
// Redistribution and use in source and binary forms, with or without modification,
// are permitted provided that the following conditions are met:
//
// * Redistribution's of source code must retain the above copyright notice,
// this list of conditions and the following disclaimer.
//
// * Redistribution's in binary form must reproduce the above copyright notice,
// this list of conditions and the following disclaimer in the documentation
// and/or other materials provided with the distribution.
//
// * The name of the copyright holders may not be used to endorse or promote products
// derived from this software without specific prior written permission.
//
// This software is provided by the copyright holders and contributors "as is" and
// any express or implied warranties, including, but not limited to, the implied
// warranties of merchantability and fitness for a particular purpose are disclaimed.
// In no event shall the Intel Corporation or contributors be liable for any direct,
// indirect, incidental, special, exemplary, or consequential damages
// (including, but not limited to, procurement of substitute goods or services;
// loss of use, data, or profits; or business interruption) however caused
// and on any theory of liability, whether in contract, strict liability,
// or tort (including negligence or otherwise) arising in any way out of
// the use of this software, even if advised of the possibility of such damage.
//
//M*/
/* the reference code has been contributed by Chris Sav */
#include "precomp.hpp"
#include "opencv2/core/hal/intrin.hpp"
#include "opencl_kernels_ximgproc.hpp"
namespace cv {
namespace ximgproc {
#if CV_SIMD128
inline void v_expand_s(const v_uint8x16& a, v_int16x8& b, v_int16x8& c)
{
v_uint16x8 t0, t1;
v_expand(a, t0, t1);
b = v_reinterpret_as_s16(t0);
c = v_reinterpret_as_s16(t1);
}
inline void v_expand_f32(const v_int16x8& a, v_float32x4& b, v_float32x4& c)
{
v_int32x4 t0, t1;
v_expand(a, t0, t1);
b = v_cvt_f32(t0);
c = v_cvt_f32(t1);
}
inline v_uint8x16 v_finalize_pix_ch(const v_int16x8& c0, const v_int16x8& c1,
const v_float32x4& s0, const v_float32x4& s1,
const v_float32x4& s2, const v_float32x4& s3,
const v_float32x4& alpha)
{
v_float32x4 f0, f1, f2, f3;
v_expand_f32(c0, f0, f1);
v_expand_f32(c1, f2, f3);
v_int16x8 d0 = v_pack(v_round(s0*alpha + f0), v_round(s1*alpha + f1));
v_int16x8 d1 = v_pack(v_round(s2*alpha + f2), v_round(s3*alpha + f3));
return v_pack_u(d0, d1);
}
#endif
class ADBody : public ParallelLoopBody
{
public:
ADBody(const Mat* src_, Mat* dst_, const float* exptab, float alpha)
{
src = src_;
dst = dst_;
exptab_ = exptab;
alpha_ = alpha;
}
void operator()(const Range& range) const
{
const int cn = 3;
int cols = src->cols;
int step = (int)src->step;
int tab[] = { -cn, cn, -step-cn, -step, -step+cn, step-cn, step, step+cn };
float alpha = alpha_;
const float* exptab = exptab_;
for( int i = range.start; i < range.end; i++ )
{
const uchar* psrc0 = src->ptr<uchar>(i);
uchar* pdst = dst->ptr<uchar>(i);
int j = 0;
#if CV_SIMD128
v_float32x4 v_alpha = v_setall_f32(alpha);
for( ; j <= cols - 16; j += 16 )
{
v_uint8x16 c0, c1, c2;
v_load_deinterleave(psrc0 + j*3, c0, c1, c2);
v_int16x8 c00, c01, c10, c11, c20, c21;
v_expand_s(c0, c00, c01);
v_expand_s(c1, c10, c11);
v_expand_s(c2, c20, c21);
v_float32x4 s00 = v_setzero_f32(), s01 = s00, s02 = s00, s03 = s00;
v_float32x4 s10 = v_setzero_f32(), s11 = s00, s12 = s00, s13 = s00;
v_float32x4 s20 = v_setzero_f32(), s21 = s00, s22 = s00, s23 = s00;
v_float32x4 fd0, fd1, fd2, fd3;
for( int k = 0; k < 8; k++ )
{
const uchar* psrc1 = psrc0 + j*3 + tab[k];
v_uint8x16 p0, p1, p2;
v_int16x8 p00, p01, p10, p11, p20, p21;
v_load_deinterleave(psrc1, p0, p1, p2);
v_expand_s(p0, p00, p01);
v_expand_s(p1, p10, p11);
v_expand_s(p2, p20, p21);
v_int16x8 d00 = p00 - c00, d01 = p01 - c01;
v_int16x8 d10 = p10 - c10, d11 = p11 - c11;
v_int16x8 d20 = p20 - c20, d21 = p21 - c21;
v_uint16x8 n0 = v_abs(d00) + v_abs(d10) + v_abs(d20);
v_uint16x8 n1 = v_abs(d01) + v_abs(d11) + v_abs(d21);
ushort CV_DECL_ALIGNED(16) nbuf[16];
v_store(nbuf, n0);
v_store(nbuf + 8, n1);
v_float32x4 w0(exptab[nbuf[0]], exptab[nbuf[1]], exptab[nbuf[2]], exptab[nbuf[3]]);
v_float32x4 w1(exptab[nbuf[4]], exptab[nbuf[5]], exptab[nbuf[6]], exptab[nbuf[7]]);
v_float32x4 w2(exptab[nbuf[8]], exptab[nbuf[9]], exptab[nbuf[10]], exptab[nbuf[11]]);
v_float32x4 w3(exptab[nbuf[12]], exptab[nbuf[13]], exptab[nbuf[14]], exptab[nbuf[15]]);
v_expand_f32(d00, fd0, fd1);
v_expand_f32(d01, fd2, fd3);
s00 += fd0*w0; s01 += fd1*w1; s02 += fd2*w2; s03 += fd3*w3;
v_expand_f32(d10, fd0, fd1);
v_expand_f32(d11, fd2, fd3);
s10 += fd0*w0; s11 += fd1*w1; s12 += fd2*w2; s13 += fd3*w3;
v_expand_f32(d20, fd0, fd1);
v_expand_f32(d21, fd2, fd3);
s20 += fd0*w0; s21 += fd1*w1; s22 += fd2*w2; s23 += fd3*w3;
}
c0 = v_finalize_pix_ch(c00, c01, s00, s01, s02, s03, v_alpha);
c1 = v_finalize_pix_ch(c10, c11, s10, s11, s12, s13, v_alpha);
c2 = v_finalize_pix_ch(c20, c21, s20, s21, s22, s23, v_alpha);
v_store_interleave(pdst + j*3, c0, c1, c2);
}
j *= 3;
#endif
for( ; j < cols*cn; j += cn )
{
int c0 = psrc0[j], c1 = psrc0[j+1], c2 = psrc0[j+2];
float s0 = 0.f, s1 = 0.f, s2 = 0.f;
for( int k = 0; k < 8; k++ )
{
const uchar* psrc1 = psrc0 + j + tab[k];
int delta0 = psrc1[0] - c0;
int delta1 = psrc1[1] - c1;
int delta2 = psrc1[2] - c2;
int nabla = std::abs(delta0) + std::abs(delta1) + std::abs(delta2);
float w = exptab[nabla];
s0 += delta0*w;
s1 += delta1*w;
s2 += delta2*w;
}
pdst[j] = saturate_cast<uchar>(c0 + alpha*s0);
pdst[j+1] = saturate_cast<uchar>(c1 + alpha*s1);
pdst[j+2] = saturate_cast<uchar>(c2 + alpha*s2);
}
}
}
const Mat* src;
Mat* dst;
const float* exptab_;
float alpha_;
};
#ifdef HAVE_OPENCL
static bool ocl_anisotropicDiffusion(InputArray src_, OutputArray dst_,
float alpha, int niters,
const std::vector<float>& exptab)
{
UMat src0 = src_.getUMat(), dst0 = dst_.getUMat();
int type = src0.type();
int rows = src0.rows, cols = src0.cols;
ocl::Kernel k("anisodiff", ocl::ximgproc::anisodiff_oclsrc, "");
if (k.empty())
return false;
UMat temp0x(rows + 2, cols + 2, type);
UMat temp1x(rows + 2, cols + 2, type);
UMat temp0(temp0x, Rect(1, 1, cols, rows));
UMat temp1(temp1x, Rect(1, 1, cols, rows));
int tabsz = (int)exptab.size();
UMat uexptab = Mat(1, tabsz, CV_32F, (void*)&exptab[0]).getUMat(ACCESS_READ);
for (int t = 0; t < niters; t++)
{
UMat src = temp0, dst = t == niters-1 ? dst0 : temp1;
copyMakeBorder(t == 0 ? src0 : src, temp0x, 1, 1, 1, 1, BORDER_REPLICATE);
k.args(ocl::KernelArg::ReadOnlyNoSize(src), ocl::KernelArg::WriteOnly(dst),
ocl::KernelArg::PtrReadOnly(uexptab), alpha);
size_t globalsize[] = { cols, rows };
if(!k.run(2, globalsize, NULL, true))
return false;
std::swap(temp0, temp1);
std::swap(temp0x, temp1x);
}
return true;
}
#endif
void anisotropicDiffusion(InputArray src_, OutputArray dst_, float alpha, float K, int niters )
{
if( niters == 0 )
{
src_.copyTo(dst_);
return;
}
int type = src_.type();
CV_Assert(src_.dims() == 2 && type == CV_8UC3);
CV_Assert(K != 0);
CV_Assert(alpha > 0);
CV_Assert(niters >= 0);
const int cn = 3;
float sigma = K * cn * 255.f;
float isigma2 = 1 / (sigma * sigma);
std::vector<float> exptab_(255*3);
float* exptab = &exptab_[0];
for( int k = 0; k < 255*3; k++ )
exptab[k] = std::exp(-k*k*isigma2);
dst_.create(src_.size(), type);
CV_OCL_RUN(dst_.isUMat(),
ocl_anisotropicDiffusion(src_, dst_, alpha, niters, exptab_));
Mat src0 = src_.getMat();
int rows = src0.rows, cols = src0.cols;
Mat dst0 = dst_.getMat();
Mat temp0x(rows + 2, cols + 2, type);
Mat temp1x(rows + 2, cols + 2, type);
Mat temp0(temp0x, Rect(1, 1, cols, rows));
Mat temp1(temp1x, Rect(1, 1, cols, rows));
for (int t = 0; t < niters; t++)
{
Mat src = temp0, dst = t == niters-1 ? dst0 : temp1;
copyMakeBorder(t == 0 ? src0 : src, temp0x, 1, 1, 1, 1, BORDER_REPLICATE);
ADBody body(&src, &dst, exptab, alpha);
parallel_for_(Range(0, rows), body, 8);
std::swap(temp0, temp1);
std::swap(temp0x, temp1x);
}
}
}
}

39
modules/ximgproc/src/opencl/anisodiff.cl
Unescape View File

@ -0,0 +1,39 @@
__kernel void anisodiff(__global const uchar * srcptr, int srcstep, int srcoffset,
__global uchar * dstptr, int dststep, int dstoffset,
int rows, int cols, __constant float* exptab, float alpha)
{
int x = get_global_id(0);
int y = get_global_id(1);
if( x < cols && y < rows )
{
int yofs = y*dststep + x*3;
int xofs = y*srcstep + x*3;
float4 s = 0.f;
float4 c = (float4)(srcptr[xofs], srcptr[xofs+1], srcptr[xofs+2], 0.f);
float4 delta, adelta;
float w;
#define UPDATE_SUM(xofs1) \
delta = (float4)(srcptr[xofs + xofs1], srcptr[xofs + xofs1 + 1], srcptr[xofs + xofs1 + 2], 0.f) - c; \
adelta = fabs(delta); \
w = exptab[convert_int(adelta.x + adelta.y + adelta.z)]; \
s += delta*w
UPDATE_SUM(3);
UPDATE_SUM(-3);
UPDATE_SUM(-srcstep-3);
UPDATE_SUM(-srcstep);
UPDATE_SUM(-srcstep+3);
UPDATE_SUM(srcstep-3);
UPDATE_SUM(srcstep);
UPDATE_SUM(srcstep+3);
s = s*alpha + c;
uchar4 d = convert_uchar4_sat(convert_int4_rte(s));
dstptr[yofs] = d.x;
dstptr[yofs+1] = d.y;
dstptr[yofs+2] = d.z;
}
}

25
modules/ximgproc/test/test_anisodiff.cpp
Unescape View File

@ -0,0 +1,25 @@
#include "test_precomp.hpp"
using namespace cv;
using namespace std;
TEST(ximgproc_AnisotropicDiffusion, regression)
{
string folder = string(cvtest::TS::ptr()->get_data_path()) + "cv/shared/";
string original_path = folder + "fruits.png";
Mat original = imread(original_path, IMREAD_COLOR);
ASSERT_FALSE(original.empty()) << "Could not load input image " << original_path;
ASSERT_EQ(3, original.channels()) << "Load color input image " << original_path;
Mat result;
float alpha = 1.0f;
float K = 0.02f;
int niters = 10;
ximgproc::anisotropicDiffusion(original, result, alpha, K, niters);
double adiff_psnr = cvtest::PSNR(original, result);
//printf("psnr=%.2f\n", adiff_psnr);
ASSERT_GT(adiff_psnr, 25.0);
}
Write Preview
Loading…
Cancel
Save