diff --git a/modules/imgproc/perf/opencl/perf_imgwarp.cpp b/modules/imgproc/perf/opencl/perf_imgwarp.cpp index c31aba67a0..00420eb706 100644 --- a/modules/imgproc/perf/opencl/perf_imgwarp.cpp +++ b/modules/imgproc/perf/opencl/perf_imgwarp.cpp @@ -73,10 +73,6 @@ OCL_PERF_TEST_P(WarpAffineFixture, WarpAffine, const Size srcSize = get<0>(params); const int type = get<1>(params), interpolation = get<2>(params); - // BUG: OpenCL and CPU version diverges a bit - // Ticket: https://github.com/opencv/opencv/issues/26235 - const double eps = CV_MAT_DEPTH(type) <= CV_32S ? 2 : interpolation == INTER_CUBIC ? 2e-3 : 3e-2; - checkDeviceMaxMemoryAllocSize(srcSize, type); UMat src(srcSize, type), dst(srcSize, type); @@ -84,7 +80,16 @@ OCL_PERF_TEST_P(WarpAffineFixture, WarpAffine, OCL_TEST_CYCLE() cv::warpAffine(src, dst, M, srcSize, interpolation); - SANITY_CHECK(dst, eps); + if (interpolation == INTER_CUBIC) + { + SANITY_CHECK_NOTHING(); + } + else + { + // Ticket: https://github.com/opencv/opencv/issues/26235 + const double eps = CV_MAT_DEPTH(type) <= CV_32S ? 2 : interpolation == INTER_CUBIC ? 2e-3 : 3e-2; + SANITY_CHECK(dst, eps); + } } ///////////// WarpPerspective //////////////////////// @@ -107,7 +112,6 @@ OCL_PERF_TEST_P(WarpPerspectiveFixture, WarpPerspective, const WarpPerspectiveParams params = GetParam(); const Size srcSize = get<0>(params); const int type = get<1>(params), interpolation = get<2>(params); - const double eps = CV_MAT_DEPTH(type) <= CV_32S ? 1 : 1e-4; checkDeviceMaxMemoryAllocSize(srcSize, type); @@ -116,7 +120,15 @@ OCL_PERF_TEST_P(WarpPerspectiveFixture, WarpPerspective, OCL_TEST_CYCLE() cv::warpPerspective(src, dst, M, srcSize, interpolation); - SANITY_CHECK(dst, eps); + if (interpolation == INTER_CUBIC) + { + SANITY_CHECK_NOTHING(); + } + else + { + const double eps = CV_MAT_DEPTH(type) <= CV_32S ? 1 : 1e-4; + SANITY_CHECK(dst, eps); + } } ///////////// Resize //////////////////////// diff --git a/modules/imgproc/perf/perf_remap.cpp b/modules/imgproc/perf/perf_remap.cpp index bd2129d233..693650ebd8 100644 --- a/modules/imgproc/perf/perf_remap.cpp +++ b/modules/imgproc/perf/perf_remap.cpp @@ -64,7 +64,14 @@ PERF_TEST_P( TestRemap, Remap, int runs = (sz.width <= 640) ? 3 : 1; TEST_CYCLE_MULTIRUN(runs) remap(src, dst, map1, map2, inter_type); - SANITY_CHECK(dst); + if (inter_type == INTER_CUBIC) + { + SANITY_CHECK_NOTHING(); + } + else + { + SANITY_CHECK(dst); + } } } // namespace diff --git a/modules/imgproc/src/imgwarp.cpp b/modules/imgproc/src/imgwarp.cpp index 78ca11db34..b98ee04223 100644 --- a/modules/imgproc/src/imgwarp.cpp +++ b/modules/imgproc/src/imgwarp.cpp @@ -55,14 +55,162 @@ #include "opencv2/core/softfloat.hpp" #include "imgwarp.hpp" +#include "warp_common.hpp" #include "warp_kernels.simd.hpp" #include "warp_kernels.simd_declarations.hpp" -using namespace cv; - namespace cv { +///////////////////////////// new-style kernel-base image warping without tables //////////////////////// + +ImgWarpFunc getImgWarpFunc(int type, int interpolation) +{ + if (interpolation == INTER_CUBIC) { + CV_CPU_DISPATCH(getBicubicWarpFunc_, (type), CV_CPU_DISPATCH_MODES_ALL); + } + return (ImgWarpFunc)nullptr; +} + +static bool genericWarp(const Mat& src, const Mat& M_, const Mat& mapx, const Mat& mapy, Mat& dst, + int interpolation, int borderType_, const Scalar& borderValue, bool relativeMap_) +{ + constexpr int MAX_CHANNELS = 4; + int srctype = src.type(); + + ImgWarpFunc func = getImgWarpFunc(srctype, interpolation); + if (!func) { + return false; // not implemented; + } + + CV_Assert(src.type() == dst.type()); + CV_Assert(src.channels() <= MAX_CHANNELS); + Matx33f Mdata(1.f, 0.f, 0.f, 0.f, 1.f, 0.f, 0.f, 0.f, 1.f); + + if (!M_.empty()) { + if (!mapx.empty() || !mapy.empty()) { + CV_Error(Error::StsBadArg, "Either M or both mapx and mapy must be empty"); + } + CV_Assert(M_.size() == Size(3, 2) || M_.size() == Size(3, 3)); + CV_Assert(M_.type() == CV_32F || M_.type() == CV_64F); + Mat Mcopy(M_.size(), CV_32F, Mdata.val); + M_.convertTo(Mcopy, CV_32F); + } else { + if (mapx.empty()) { + CV_Error(Error::StsBadArg, "When M is empty, mapx and/or mapy must be non-empty"); + } + CV_Assert(mapx.size() == dst.size()); + if (mapy.empty()) { + CV_Assert(mapx.type() == CV_32FC2); + } else { + if (mapx.type() == mapy.type()) { + CV_Assert(mapx.type() == CV_32FC1); + } else { + CV_Assert(mapx.type() == CV_16SC2 && (mapy.type() == CV_16UC1 || mapy.type() == CV_16SC1)); + } + CV_Assert(mapx.size() == mapy.size()); + } + } + + int64_t borderValBuf_[MAX_CHANNELS]; + scalarToRawData(borderValue, borderValBuf_, src.type()); + + parallel_for_(Range(0, dst.rows), [&](const Range& range) { + constexpr float FIXPT_SCALE = 1.f/32; + constexpr int BLOCK_SIZE = 128; + const uint8_t* srcdata = src.data; + const float* fparams = nullptr; + size_t srcstep = src.step; + Size srcsize = src.size(); + int dstcols = dst.cols; + int borderType = borderType_; + size_t bpp = src.elemSize(); + + bool warping = !M_.empty(); + bool warpAffine = warping && M_.rows == 2; + bool warpPerspective = warping && M_.rows == 3; + + bool relativeMap = relativeMap_ && !warping; + float relscale = relativeMap ? 1.f : 0.f; + bool interleavedmap = !warping && (mapx.type() == CV_32FC2); + bool planarmap = !warping && (mapx.type() == mapy.type()); + bool fixedmap = !warping && (mapx.type() == CV_16SC2); + + CV_Assert(warping || interleavedmap || planarmap || fixedmap); + + Matx33f M = Mdata; + const uint8_t* borderValBuf = (const uint8_t*)borderValBuf_; + float M_xx = M(0, 0), M_yx = M(1, 0), M_zx = M(2, 0); + float xbuf[BLOCK_SIZE], ybuf[BLOCK_SIZE]; + const float* xbufptr = xbuf; + const float* ybufptr = ybuf; + + for (int y = range.start; y < range.end; y++) { + uint8_t* dstptr = dst.ptr(y); + float M_x = float(y)*M(0, 1) + M(0, 2); + float M_y = float(y)*M(1, 1) + M(1, 2); + float M_z = float(y)*M(2, 1) + M(2, 2); + const Vec2s* xysptr = fixedmap ? mapx.ptr(y) : nullptr; + const ushort* idxptr = fixedmap ? mapy.ptr(y) : nullptr; + const Vec2f* xyfptr = interleavedmap ? mapx.ptr(y) : nullptr; + const float* xfptr = planarmap ? mapx.ptr(y) : nullptr; + const float* yfptr = planarmap ? mapy.ptr(y) : nullptr; + float y0 = relativeMap ? float(y) : 0.f; + + for (int x = 0; x < dstcols; x += BLOCK_SIZE) { + int blocksize = std::min(BLOCK_SIZE, dstcols - x); + if (warpAffine) { + for (int dx = 0; dx < blocksize; dx++) { + float xf = float(x + dx); + xbuf[dx] = M_x + M_xx*xf; + ybuf[dx] = M_y + M_yx*xf; + } + } else if (warpPerspective) { + for (int dx = 0; dx < blocksize; dx++) { + double xf = double(x + dx); + double invz = 1./(M_z + M_zx*xf); + xbuf[dx] = float((M_x + M_xx*xf)*invz); + ybuf[dx] = float((M_y + M_yx*xf)*invz); + } + } else if (fixedmap) { + for (int dx = 0; dx < blocksize; dx++) { + Vec2s xy = xysptr[x + dx]; + ushort idx = idxptr[x + dx]; + float xf = float(xy[0]) + (idx & 31)*FIXPT_SCALE + (x + dx)*relscale; + float yf = float(xy[1]) + ((idx >> 5) & 31)*FIXPT_SCALE + y0; + xbuf[dx] = xf; + ybuf[dx] = yf; + } + } else if (interleavedmap) { + for (int dx = 0; dx < blocksize; dx++) { + Vec2f xy = xyfptr[x + dx]; + float xf = xy[0] + (x + dx)*relscale; + float yf = xy[1] + y0; + xbuf[dx] = xf; + ybuf[dx] = yf; + } + } else if (relativeMap) { + for (int dx = 0; dx < blocksize; dx++) { + float xf = xfptr[x + dx] + (x + dx)*relscale; + float yf = yfptr[x + dx] + y0; + xbuf[dx] = xf; + ybuf[dx] = yf; + } + } else { + // planar absolute map — just use it as-is without copying + xbufptr = xfptr + x; + ybufptr = yfptr + x; + } + + func(xbufptr, ybufptr, blocksize, srcdata, srcstep, srcsize, + dstptr + x*bpp, fparams, borderType, borderValBuf); + } + } + }); + + return true; +} + /************** interpolation formulas and tables ***************/ const int INTER_REMAP_COEF_BITS=15; @@ -78,9 +226,6 @@ static short BilinearTab_iC4_buf[INTER_TAB_SIZE2+2][2][8]; static short (*BilinearTab_iC4)[2][8] = (short (*)[2][8])alignPtr(BilinearTab_iC4_buf, 16); #endif -static float BicubicTab_f[INTER_TAB_SIZE2][4][4]; -static short BicubicTab_i[INTER_TAB_SIZE2][4][4]; - static float Lanczos4Tab_f[INTER_TAB_SIZE2][8][8]; static short Lanczos4Tab_i[INTER_TAB_SIZE2][8][8]; @@ -90,16 +235,6 @@ static inline void interpolateLinear( float x, float* coeffs ) coeffs[1] = x; } -static inline void interpolateCubic( float x, float* coeffs ) -{ - const float A = -0.75f; - - coeffs[0] = ((A*(x + 1) - 5*A)*(x + 1) + 8*A)*(x + 1) - 4*A; - coeffs[1] = ((A + 2)*x - (A + 3))*x*x + 1; - coeffs[2] = ((A + 2)*(1 - x) - (A + 3))*(1 - x)*(1 - x) + 1; - coeffs[3] = 1.f - coeffs[0] - coeffs[1] - coeffs[2]; -} - static inline void interpolateLanczos4( float x, float* coeffs ) { static const double s45 = 0.70710678118654752440084436210485; @@ -138,8 +273,7 @@ static void initInterTab1D(int method, float* tab, int tabsz) } else if( method == INTER_CUBIC ) { - for( int i = 0; i < tabsz; i++, tab += 4 ) - interpolateCubic( i*scale, tab ); + ; // pass } else if( method == INTER_LANCZOS4 ) { @@ -154,13 +288,15 @@ static void initInterTab1D(int method, float* tab, int tabsz) static const void* initInterTab2D( int method, bool fixpt ) { static bool inittab[INTER_MAX+1] = {false}; + + if( method == INTER_CUBIC ) + return nullptr; + float* tab = 0; short* itab = 0; int ksize = 0; if( method == INTER_LINEAR ) tab = BilinearTab_f[0][0], itab = BilinearTab_i[0][0], ksize=2; - else if( method == INTER_CUBIC ) - tab = BicubicTab_f[0][0], itab = BicubicTab_i[0][0], ksize=4; else if( method == INTER_LANCZOS4 ) tab = Lanczos4Tab_f[0][0], itab = Lanczos4Tab_i[0][0], ksize=8; else @@ -232,8 +368,6 @@ static bool initAllInterTab2D() { return initInterTab2D( INTER_LINEAR, false ) && initInterTab2D( INTER_LINEAR, true ) && - initInterTab2D( INTER_CUBIC, false ) && - initInterTab2D( INTER_CUBIC, true ) && initInterTab2D( INTER_LANCZOS4, false ) && initInterTab2D( INTER_LANCZOS4, true ); } @@ -845,111 +979,6 @@ static void remapBilinear( const Mat& _src, Mat& _dst, const Mat& _xy, } -template -static void remapBicubic( const Mat& _src, Mat& _dst, const Mat& _xy, - const Mat& _fxy, const void* _wtab, - int borderType, const Scalar& _borderValue, const Point& _offset ) -{ - typedef typename CastOp::rtype T; - typedef typename CastOp::type1 WT; - Size ssize = _src.size(), dsize = _dst.size(); - const int cn = _src.channels(); - const AT* wtab = (const AT*)_wtab; - const T* S0 = _src.ptr(); - size_t sstep = _src.step/sizeof(S0[0]); - T cval[CV_CN_MAX]; - CastOp castOp; - - for(int k = 0; k < cn; k++ ) - cval[k] = saturate_cast(_borderValue[k & 3]); - - int borderType1 = borderType != BORDER_TRANSPARENT ? borderType : BORDER_REFLECT_101; - - unsigned width1 = std::max(ssize.width-3, 0), height1 = std::max(ssize.height-3, 0); - - if( _dst.isContinuous() && _xy.isContinuous() && _fxy.isContinuous() && !isRelative ) - { - dsize.width *= dsize.height; - dsize.height = 1; - } - - for(int dy = 0; dy < dsize.height; dy++ ) - { - T* D = _dst.ptr(dy); - const short* XY = _xy.ptr(dy); - const ushort* FXY = _fxy.ptr(dy); - const int off_y = isRelative ? (_offset.y+dy) : 0; - for(int dx = 0; dx < dsize.width; dx++, D += cn ) - { - const int off_x = isRelative ? (_offset.x+dx) : 0; - int sx = XY[dx*2]-1+off_x, sy = XY[dx*2+1]-1+off_y; - const AT* w = wtab + FXY[dx]*16; - if( (unsigned)sx < width1 && (unsigned)sy < height1 ) - { - const T* S = S0 + sy*sstep + sx*cn; - for(int k = 0; k < cn; k++ ) - { - WT sum = S[0]*w[0] + S[cn]*w[1] + S[cn*2]*w[2] + S[cn*3]*w[3]; - S += sstep; - sum += S[0]*w[4] + S[cn]*w[5] + S[cn*2]*w[6] + S[cn*3]*w[7]; - S += sstep; - sum += S[0]*w[8] + S[cn]*w[9] + S[cn*2]*w[10] + S[cn*3]*w[11]; - S += sstep; - sum += S[0]*w[12] + S[cn]*w[13] + S[cn*2]*w[14] + S[cn*3]*w[15]; - S -= sstep * 3 - 1; - D[k] = castOp(sum); - } - } - else - { - int x[4], y[4]; - if( borderType == BORDER_TRANSPARENT && - ((unsigned)(sx+1) >= (unsigned)ssize.width || - (unsigned)(sy+1) >= (unsigned)ssize.height) ) - continue; - - if( borderType1 == BORDER_CONSTANT && - (sx >= ssize.width || sx+4 <= 0 || - sy >= ssize.height || sy+4 <= 0)) - { - for(int k = 0; k < cn; k++ ) - D[k] = cval[k]; - continue; - } - - for(int i = 0; i < 4; i++ ) - { - x[i] = borderInterpolate(sx + i, ssize.width, borderType1)*cn; - y[i] = borderInterpolate(sy + i, ssize.height, borderType1); - } - - for(int k = 0; k < cn; k++, S0++, w -= 16 ) - { - WT cv = cval[k], sum = cv*ONE; - for(int i = 0; i < 4; i++, w += 4 ) - { - int yi = y[i]; - if( yi < 0 ) - continue; - const T* S = S0 + yi*sstep; - if( x[0] >= 0 ) - sum += (S[x[0]] - cv)*w[0]; - if( x[1] >= 0 ) - sum += (S[x[1]] - cv)*w[1]; - if( x[2] >= 0 ) - sum += (S[x[2]] - cv)*w[2]; - if( x[3] >= 0 ) - sum += (S[x[3]] - cv)*w[3]; - } - D[k] = castOp(sum); - } - S0 -= cn; - } - } - } -} - - template static void remapLanczos4( const Mat& _src, Mat& _dst, const Mat& _xy, const Mat& _fxy, const void* _wtab, @@ -1395,7 +1424,6 @@ void cv::remap( InputArray _src, OutputArray _dst, _dst.create( map1.size(), src.type() ); Mat dst = _dst.getMat(); - CV_Assert( dst.cols < SHRT_MAX && dst.rows < SHRT_MAX && src.cols < SHRT_MAX && src.rows < SHRT_MAX ); if( dst.data == src.data ) @@ -1423,6 +1451,10 @@ void cv::remap( InputArray _src, OutputArray _dst, if( interpolation == INTER_AREA ) interpolation = INTER_LINEAR; + if (genericWarp(src, Mat(), map1, map2, dst, interpolation, borderType, borderValue, hasRelativeFlag)) { + return; + } + int type = src.type(), depth = CV_MAT_DEPTH(type); if (interpolation == INTER_NEAREST && map1.depth() == CV_32F) { @@ -1575,24 +1607,6 @@ void cv::remap( InputArray _src, OutputArray _dst, } }; - static RemapFunc cubic_tab[2][CV_DEPTH_MAX] = - { - { - remapBicubic, short, INTER_REMAP_COEF_SCALE, false>, 0, - remapBicubic, float, 1, false>, - remapBicubic, float, 1, false>, 0, - remapBicubic, float, 1, false>, - remapBicubic, float, 1, false>, 0 - }, - { - remapBicubic, short, INTER_REMAP_COEF_SCALE, true>, 0, - remapBicubic, float, 1, true>, - remapBicubic, float, 1, true>, 0, - remapBicubic, float, 1, true>, - remapBicubic, float, 1, true>, 0 - } - }; - static RemapFunc lanczos4_tab[2][8] = { { @@ -1621,10 +1635,6 @@ void cv::remap( InputArray _src, OutputArray _dst, { if( interpolation == INTER_LINEAR ) ifunc = linear_tab[relativeOptionIndex][depth]; - else if( interpolation == INTER_CUBIC ){ - ifunc = cubic_tab[relativeOptionIndex][depth]; - CV_Assert( _src.channels() <= 4 ); - } else if( interpolation == INTER_LANCZOS4 ){ ifunc = lanczos4_tab[relativeOptionIndex][depth]; CV_Assert( _src.channels() <= 4 ); @@ -2052,7 +2062,7 @@ static bool ocl_warpTransform_cols4(InputArray _src, OutputArray _dst, InputArra if ( !dev.isIntel() || !(type == CV_8UC1) || !(dtype == CV_8UC1) || !(_dst.cols() % 4 == 0) || - (op_type == OCL_OP_PERSPECTIVE && interpolation == INTER_LINEAR && (cn == 1 || cn == 3 || cn == 4)) || + (op_type == OCL_OP_PERSPECTIVE && interpolation != INTER_NEAREST) || !(borderType == cv::BORDER_CONSTANT && (interpolation == cv::INTER_NEAREST || interpolation == cv::INTER_LINEAR || interpolation == cv::INTER_CUBIC))) return false; @@ -2062,7 +2072,7 @@ static bool ocl_warpTransform_cols4(InputArray _src, OutputArray _dst, InputArra ocl::ProgramSource program = ocl::imgproc::warp_transform_oclsrc; String kernelName = format("warp%s_%s_8u", warp_op[op_type], interpolationMap[interpolation]); - bool is32f = (interpolation == INTER_CUBIC || interpolation == INTER_LINEAR) && op_type == OCL_OP_AFFINE; + bool is32f = interpolation == INTER_CUBIC || interpolation == INTER_LINEAR; int wdepth = interpolation == INTER_NEAREST ? depth : std::max(is32f ? CV_32F : CV_32S, depth); int sctype = CV_MAKETYPE(wdepth, cn); @@ -2145,9 +2155,7 @@ static bool ocl_warpTransform(InputArray _src, OutputArray _dst, InputArray _M0, const char * const kernelName = op_type == OCL_OP_AFFINE ? "warpAffine" : "warpPerspective"; int scalarcn = cn == 3 ? 4 : cn; - bool is32f = op_type == OCL_OP_AFFINE ? - /* Affine*/ !dev.isAMD() && (interpolation == INTER_CUBIC || interpolation == INTER_LINEAR) : - /* Perspective*/ interpolation == INTER_LINEAR; + bool is32f = interpolation == INTER_CUBIC || interpolation == INTER_LINEAR; int wdepth = interpolation == INTER_NEAREST ? depth : std::max(is32f ? CV_32F : CV_32S, depth); int sctype = CV_MAKETYPE(wdepth, scalarcn); @@ -2490,6 +2498,10 @@ void cv::warpAffine( InputArray _src, OutputArray _dst, M[2] = b1; M[5] = b2; } + if (genericWarp(src, matM, Mat(), Mat(), dst, interpolation, borderType, borderValue, false)) { + return; + } + hal::warpAffine(src.type(), src.data, src.step, src.cols, src.rows, dst.data, dst.step, dst.cols, dst.rows, M, interpolation, borderType, borderValue.val, hint); } @@ -3009,6 +3021,8 @@ void cv::warpPerspective( InputArray _src, OutputArray _dst, InputArray _M0, CV_Assert( _src.total() > 0 ); + int interpolation = flags & INTER_MAX; + CV_OCL_RUN(_src.dims() <= 2 && _dst.isUMat() && _src.cols() <= SHRT_MAX && _src.rows() <= SHRT_MAX, ocl_warpTransform_cols4(_src, _dst, _M0, dsize, flags, borderType, borderValue, @@ -3027,7 +3041,7 @@ void cv::warpPerspective( InputArray _src, OutputArray _dst, InputArray _M0, double M[9]; Mat matM(3, 3, CV_64F, M); - int interpolation = flags & INTER_MAX; + if( interpolation == INTER_AREA ) interpolation = INTER_LINEAR; @@ -3037,6 +3051,10 @@ void cv::warpPerspective( InputArray _src, OutputArray _dst, InputArray _M0, if( !(flags & WARP_INVERSE_MAP) ) invert(matM, matM); + if (genericWarp(src, matM, Mat(), Mat(), dst, interpolation, borderType, borderValue, false)) { + return; + } + hal::warpPerspective(src.type(), src.data, src.step, src.cols, src.rows, dst.data, dst.step, dst.cols, dst.rows, matM.ptr(), interpolation, borderType, borderValue.val, hint); } diff --git a/modules/imgproc/src/opencl/remap.cl b/modules/imgproc/src/opencl/remap.cl index 9a5b30d193..106b80460e 100644 --- a/modules/imgproc/src/opencl/remap.cl +++ b/modules/imgproc/src/opencl/remap.cl @@ -51,6 +51,13 @@ #endif #endif +#define CV_64F 6 +#if defined SRC_DEPTH && SRC_DEPTH == CV_64F +#define WT1 double +#else +#define WT1 float +#endif + #define noconvert #if CN != 3 @@ -424,8 +431,8 @@ __kernel void remap_2_32FC1(__global const uchar * srcptr, int src_step, int src int sx = convert_int_rtn(X0); int sy = convert_int_rtn(Y0); - float ax = X0 - (float) sx; - float ay = Y0 - (float) sy; + WT1 ax = (WT1)(X0 - (float) sx); + WT1 ay = (WT1)(Y0 - (float) sy); int2 map_data0 = (int2)(sx, sy); int2 map_data1 = (int2)(sx+1, sy); @@ -493,8 +500,8 @@ __kernel void remap_32FC2(__global const uchar * srcptr, int src_step, int src_o int sx = convert_int_rtn(X0); int sy = convert_int_rtn(Y0); - float ax = X0 - (float) sx; - float ay = Y0 - (float) sy; + WT1 ax = (WT1)(X0 - (float) sx); + WT1 ay = (WT1)(Y0 - (float) sy); int2 map_data0 = (int2)(sx, sy); int2 map_data1 = (int2)(sx+1, sy); diff --git a/modules/imgproc/src/opencl/warp_affine.cl b/modules/imgproc/src/opencl/warp_affine.cl index f83aa17b0b..f385cf8cd3 100644 --- a/modules/imgproc/src/opencl/warp_affine.cl +++ b/modules/imgproc/src/opencl/warp_affine.cl @@ -51,6 +51,13 @@ #endif #endif +#define CV_64F 6 +#if defined SRC_DEPTH && SRC_DEPTH == CV_64F +#define WT1 double +#else +#define WT1 float +#endif + #define INTER_BITS 5 #define INTER_TAB_SIZE (1 << INTER_BITS) #define INTER_SCALE 1.f/INTER_TAB_SIZE @@ -137,8 +144,8 @@ __kernel void warpAffine(__global const uchar * srcptr, int src_step, int src_of int sx = convert_short_rtn(X0); int sy = convert_short_rtn(Y0); - float ax = X0 - (CT)sx; - float ay = Y0 - (CT)sy; + WT1 ax = (WT1)(X0 - (CT)sx); + WT1 ay = (WT1)(Y0 - (CT)sy); WT v0 = scalar, v1 = scalar, v2 = scalar, v3 = scalar; if (sx >= 0 && sx < src_cols) @@ -168,35 +175,16 @@ __kernel void warpAffine(__global const uchar * srcptr, int src_step, int src_of #elif defined INTER_CUBIC -#ifdef AMD_DEVICE - -inline void interpolateCubic( float x, float* coeffs ) +inline void interpolateCubic( float x, WT1* coeffs ) { const float A = -0.75f; - coeffs[0] = fma(fma(fma(A, (x + 1.f), - 5.0f*A), (x + 1.f), 8.0f*A), x + 1.f, - 4.0f*A); - coeffs[1] = fma(fma(A + 2.f, x, - (A + 3.f)), x*x, 1.f); - coeffs[2] = fma(fma(A + 2.f, 1.f - x, - (A + 3.f)), (1.f - x)*(1.f - x), 1.f); - coeffs[3] = 1.f - coeffs[0] - coeffs[1] - coeffs[2]; + coeffs[0] = (WT1)fma(fma(fma(A, (x + 1.f), - 5.0f*A), (x + 1.f), 8.0f*A), x + 1.f, - 4.0f*A); + coeffs[1] = (WT1)fma(fma(A + 2.f, x, - (A + 3.f)), x*x, 1.f); + coeffs[2] = (WT1)fma(fma(A + 2.f, 1.f - x, - (A + 3.f)), (1.f - x)*(1.f - x), 1.f); + coeffs[3] = (WT1)(1. - coeffs[0] - coeffs[1] - coeffs[2]); } -#else - -__constant float coeffs[128] = - { 0.000000f, 1.000000f, 0.000000f, 0.000000f, -0.021996f, 0.997841f, 0.024864f, -0.000710f, -0.041199f, 0.991516f, 0.052429f, -0.002747f, - -0.057747f, 0.981255f, 0.082466f, -0.005974f, -0.071777f, 0.967285f, 0.114746f, -0.010254f, -0.083427f, 0.949837f, 0.149040f, -0.015450f, - -0.092834f, 0.929138f, 0.185120f, -0.021423f, -0.100136f, 0.905418f, 0.222755f, -0.028038f, -0.105469f, 0.878906f, 0.261719f, -0.035156f, - -0.108971f, 0.849831f, 0.301781f, -0.042641f, -0.110779f, 0.818420f, 0.342712f, -0.050354f, -0.111031f, 0.784904f, 0.384285f, -0.058159f, - -0.109863f, 0.749512f, 0.426270f, -0.065918f, -0.107414f, 0.712471f, 0.468437f, -0.073494f, -0.103821f, 0.674011f, 0.510559f, -0.080750f, - -0.099220f, 0.634361f, 0.552406f, -0.087547f, -0.093750f, 0.593750f, 0.593750f, -0.093750f, -0.087547f, 0.552406f, 0.634361f, -0.099220f, - -0.080750f, 0.510559f, 0.674011f, -0.103821f, -0.073494f, 0.468437f, 0.712471f, -0.107414f, -0.065918f, 0.426270f, 0.749512f, -0.109863f, - -0.058159f, 0.384285f, 0.784904f, -0.111031f, -0.050354f, 0.342712f, 0.818420f, -0.110779f, -0.042641f, 0.301781f, 0.849831f, -0.108971f, - -0.035156f, 0.261719f, 0.878906f, -0.105469f, -0.028038f, 0.222755f, 0.905418f, -0.100136f, -0.021423f, 0.185120f, 0.929138f, -0.092834f, - -0.015450f, 0.149040f, 0.949837f, -0.083427f, -0.010254f, 0.114746f, 0.967285f, -0.071777f, -0.005974f, 0.082466f, 0.981255f, -0.057747f, - -0.002747f, 0.052429f, 0.991516f, -0.041199f, -0.000710f, 0.024864f, 0.997841f, -0.021996f }; - -#endif - __kernel void warpAffine(__global const uchar * srcptr, int src_step, int src_offset, int src_rows, int src_cols, __global uchar * dstptr, int dst_step, int dst_offset, int dst_rows, int dst_cols, __constant CT * M, ST scalar_) @@ -206,103 +194,60 @@ __kernel void warpAffine(__global const uchar * srcptr, int src_step, int src_of if (dx < dst_cols && dy < dst_rows) { - int tmp = (dx << AB_BITS); - int X0 = rint(M[0] * tmp) + rint(fma(M[1], (CT)dy, M[2]) * AB_SCALE) + ROUND_DELTA; - int Y0 = rint(M[3] * tmp) + rint(fma(M[4], (CT)dy, M[5]) * AB_SCALE) + ROUND_DELTA; + float X0 = (float)fma(M[1], (CT)dy, fma(M[0], (CT)dx, M[2])); + float Y0 = (float)fma(M[4], (CT)dy, fma(M[3], (CT)dx, M[5])); - X0 = X0 >> (AB_BITS - INTER_BITS); - Y0 = Y0 >> (AB_BITS - INTER_BITS); + int sx = convert_short_rtn(X0); + int sy = convert_short_rtn(Y0); + float ax = X0 - (float)sx; + float ay = Y0 - (float)sy; - int sx = (short)(X0 >> INTER_BITS) - 1, sy = (short)(Y0 >> INTER_BITS) - 1; - int ay = (short)(Y0 & (INTER_TAB_SIZE - 1)), ax = (short)(X0 & (INTER_TAB_SIZE - 1)); + sx--; + sy--; -#ifdef AMD_DEVICE - WT v[16]; - #pragma unroll - for (int y = 0; y < 4; y++) - { - if (sy+y >= 0 && sy+y < src_rows) - { - #pragma unroll - for (int x = 0; x < 4; x++) - v[mad24(y, 4, x)] = sx+x >= 0 && sx+x < src_cols ? - CONVERT_TO_WT(loadpix(srcptr + mad24(sy+y, src_step, mad24(sx+x, pixsize, src_offset)))) : scalar; - } - else - { - #pragma unroll - for (int x = 0; x < 4; x++) - v[mad24(y, 4, x)] = scalar; - } - } - - float tab1y[4], tab1x[4]; - - float ayy = INTER_SCALE * ay; - float axx = INTER_SCALE * ax; - interpolateCubic(ayy, tab1y); - interpolateCubic(axx, tab1x); - - int dst_index = mad24(dy, dst_step, mad24(dx, pixsize, dst_offset)); + WT1 taby[4], tabx[4]; + interpolateCubic(ay, taby); + interpolateCubic(ax, tabx); WT sum = (WT)(0); -#if SRC_DEPTH <= 4 - int itab[16]; - #pragma unroll - for (int i = 0; i < 16; i++) - itab[i] = rint(tab1y[(i>>2)] * tab1x[(i&3)] * INTER_REMAP_COEF_SCALE); - - #pragma unroll - for (int i = 0; i < 16; i++) - sum = mad24(v[i], itab[i], sum); - storepix(CONVERT_TO_T( (sum + (1 << (INTER_REMAP_COEF_BITS-1))) >> INTER_REMAP_COEF_BITS ), dstptr + dst_index); -#else - #pragma unroll - for (int i = 0; i < 16; i++) - sum = fma(v[i], tab1y[(i>>2)] * tab1x[(i&3)], sum); - storepix(CONVERT_TO_T( sum ), dstptr + dst_index); -#endif -#else // INTEL_DEVICE - __constant float * coeffs_y = coeffs + (ay << 2), * coeffs_x = coeffs + (ax << 2); - - int src_index0 = mad24(sy, src_step, mad24(sx, pixsize, src_offset)), src_index; - int dst_index = mad24(dy, dst_step, mad24(dx, pixsize, dst_offset)); - - WT sum = (WT)(0), xsum; - #pragma unroll - for (int y = 0; y < 4; y++) - { - src_index = mad24(y, src_step, src_index0); - if (sy + y >= 0 && sy + y < src_rows) - { - xsum = (WT)(0); - if (sx >= 0 && sx + 4 < src_cols) - { -#if SRC_DEPTH == 0 && CN == 1 - uchar4 value = vload4(0, srcptr + src_index); - xsum = dot(convert_float4(value), (float4)(coeffs_x[0], coeffs_x[1], coeffs_x[2], coeffs_x[3])); -#else - #pragma unroll - for (int x = 0; x < 4; x++) - xsum = fma(CONVERT_TO_WT(loadpix(srcptr + mad24(x, pixsize, src_index))), coeffs_x[x], xsum); -#endif - } - else - { - #pragma unroll - for (int x = 0; x < 4; x++) - xsum = fma(sx + x >= 0 && sx + x < src_cols ? - CONVERT_TO_WT(loadpix(srcptr + mad24(x, pixsize, src_index))) : scalar, coeffs_x[x], xsum); - } - sum = fma(xsum, coeffs_y[y], sum); + if (0 <= sy && sy + 4 <= src_rows && + 0 <= sx && sx + 4 <= src_cols) { + #pragma unroll + for (int y = 0; y < 4; y++) { + int row_offset = mad24(sy+y, src_step, src_offset); + WT v0 = CONVERT_TO_WT(loadpix(srcptr + mad24(sx, pixsize, row_offset))); + WT v1 = CONVERT_TO_WT(loadpix(srcptr + mad24(sx + 1, pixsize, row_offset))); + WT v2 = CONVERT_TO_WT(loadpix(srcptr + mad24(sx + 2, pixsize, row_offset))); + WT v3 = CONVERT_TO_WT(loadpix(srcptr + mad24(sx + 3, pixsize, row_offset))); + WT wsum = (WT)(0); + wsum = fma(v0, tabx[0], wsum); + wsum = fma(v1, tabx[1], wsum); + wsum = fma(v2, tabx[2], wsum); + wsum = fma(v3, tabx[3], wsum); + sum = fma(wsum, taby[y], sum); + } + } + else { + #pragma unroll + for (int y = 0; y < 4; y++) { + if (sy+y >= 0 && sy+y < src_rows) { + int row_offset = mad24(sy+y, src_step, src_offset); + #pragma unroll + for (int x = 0; x < 4; x++) { + WT v = sx+x >= 0 && sx+x < src_cols ? + CONVERT_TO_WT(loadpix(srcptr + mad24(sx + x, pixsize, row_offset))) : scalar; + sum = fma(v, taby[y] * tabx[x], sum); + } + } + else { + sum = fma(scalar, taby[y], sum); + } } - else - sum = fma(scalar, coeffs_y[y], sum); } - storepix(CONVERT_TO_T(sum), dstptr + dst_index); -#endif + int dst_index = mad24(dy, dst_step, mad24(dx, pixsize, dst_offset)); + storepix(CONVERT_TO_T( sum ), dstptr + dst_index); } } diff --git a/modules/imgproc/src/opencl/warp_perspective.cl b/modules/imgproc/src/opencl/warp_perspective.cl index b28eca85b8..2050660d3f 100644 --- a/modules/imgproc/src/opencl/warp_perspective.cl +++ b/modules/imgproc/src/opencl/warp_perspective.cl @@ -51,6 +51,13 @@ #endif #endif +#define CV_64F 6 +#if defined SRC_DEPTH && SRC_DEPTH == CV_64F +#define WT1 double +#else +#define WT1 float +#endif + #define INTER_BITS 5 #define INTER_TAB_SIZE (1 << INTER_BITS) #define INTER_SCALE 1.f / INTER_TAB_SIZE @@ -150,14 +157,14 @@ __kernel void warpPerspective(__global const uchar * srcptr, int src_step, int s #elif defined INTER_CUBIC -inline void interpolateCubic( float x, float* coeffs ) +inline void interpolateCubic( float x, WT1* coeffs ) { const float A = -0.75f; - coeffs[0] = ((A*(x + 1.f) - 5.0f*A)*(x + 1.f) + 8.0f*A)*(x + 1.f) - 4.0f*A; - coeffs[1] = ((A + 2.f)*x - (A + 3.f))*x*x + 1.f; - coeffs[2] = ((A + 2.f)*(1.f - x) - (A + 3.f))*(1.f - x)*(1.f - x) + 1.f; - coeffs[3] = 1.f - coeffs[0] - coeffs[1] - coeffs[2]; + coeffs[0] = (WT1)fma(fma(fma(A, (x + 1.f), - 5.0f*A), (x + 1.f), 8.0f*A), x + 1.f, - 4.0f*A); + coeffs[1] = (WT1)fma(fma(A + 2.f, x, - (A + 3.f)), x*x, 1.f); + coeffs[2] = (WT1)fma(fma(A + 2.f, 1.f - x, - (A + 3.f)), (1.f - x)*(1.f - x), 1.f); + coeffs[3] = (WT1)(1. - coeffs[0] - coeffs[1] - coeffs[2]); } __kernel void warpPerspective(__global const uchar * srcptr, int src_step, int src_offset, int src_rows, int src_cols, @@ -169,52 +176,61 @@ __kernel void warpPerspective(__global const uchar * srcptr, int src_step, int s if (dx < dst_cols && dy < dst_rows) { - CT X0 = M[0] * dx + M[1] * dy + M[2]; - CT Y0 = M[3] * dx + M[4] * dy + M[5]; - CT W = M[6] * dx + M[7] * dy + M[8]; - W = W != 0.0f ? INTER_TAB_SIZE / W : 0.0f; - int X = rint(X0 * W), Y = rint(Y0 * W); + float W = (float)fma(M[7], (CT)dy, fma(M[6], (CT)dx, M[8])); + float X0 = (float)fma(M[1], (CT)dy, fma(M[0], (CT)dx, M[2])) / W; + float Y0 = (float)fma(M[4], (CT)dy, fma(M[3], (CT)dx, M[5])) / W; - short sx = convert_short_sat(X >> INTER_BITS) - 1; - short sy = convert_short_sat(Y >> INTER_BITS) - 1; - short ay = (short)(Y & (INTER_TAB_SIZE-1)); - short ax = (short)(X & (INTER_TAB_SIZE-1)); + int sx = convert_short_rtn(X0); + int sy = convert_short_rtn(Y0); + float ax = X0 - (float)sx; + float ay = Y0 - (float)sy; - WT v[16]; - #pragma unroll - for (int y = 0; y < 4; y++) - #pragma unroll - for (int x = 0; x < 4; x++) - v[mad24(y, 4, x)] = (sx+x >= 0 && sx+x < src_cols && sy+y >= 0 && sy+y < src_rows) ? - CONVERT_TO_WT(loadpix(srcptr + mad24(sy+y, src_step, src_offset + (sx+x) * pixsize))) : scalar; + sx--; + sy--; - float tab1y[4], tab1x[4]; - - float ayy = INTER_SCALE * ay; - float axx = INTER_SCALE * ax; - interpolateCubic(ayy, tab1y); - interpolateCubic(axx, tab1x); - - int dst_index = mad24(dy, dst_step, dst_offset + dx * pixsize); + WT1 taby[4], tabx[4]; + interpolateCubic(ay, taby); + interpolateCubic(ax, tabx); WT sum = (WT)(0); -#if SRC_DEPTH <= 4 - int itab[16]; - #pragma unroll - for (int i = 0; i < 16; i++) - itab[i] = rint(tab1y[(i>>2)] * tab1x[(i&3)] * INTER_REMAP_COEF_SCALE); + if (0 <= sy && sy + 4 <= src_rows && + 0 <= sx && sx + 4 <= src_cols) { + #pragma unroll + for (int y = 0; y < 4; y++) { + int row_offset = mad24(sy+y, src_step, src_offset); + WT v0 = CONVERT_TO_WT(loadpix(srcptr + mad24(sx, pixsize, row_offset))); + WT v1 = CONVERT_TO_WT(loadpix(srcptr + mad24(sx + 1, pixsize, row_offset))); + WT v2 = CONVERT_TO_WT(loadpix(srcptr + mad24(sx + 2, pixsize, row_offset))); + WT v3 = CONVERT_TO_WT(loadpix(srcptr + mad24(sx + 3, pixsize, row_offset))); + WT wsum = (WT)(0); + wsum = fma(v0, tabx[0], wsum); + wsum = fma(v1, tabx[1], wsum); + wsum = fma(v2, tabx[2], wsum); + wsum = fma(v3, tabx[3], wsum); + sum = fma(wsum, taby[y], sum); + } + } + else { + #pragma unroll + for (int y = 0; y < 4; y++) { + if (sy+y >= 0 && sy+y < src_rows) { + int row_offset = mad24(sy+y, src_step, src_offset); + #pragma unroll + for (int x = 0; x < 4; x++) { + WT v = sx+x >= 0 && sx+x < src_cols ? + CONVERT_TO_WT(loadpix(srcptr + mad24(sx + x, pixsize, row_offset))) : scalar; + sum = fma(v, taby[y] * tabx[x], sum); + } + } + else { + sum = fma(scalar, taby[y], sum); + } + } + } - #pragma unroll - for (int i = 0; i < 16; i++) - sum += v[i] * itab[i]; - storepix(CONVERT_TO_T( (sum + (1 << (INTER_REMAP_COEF_BITS-1))) >> INTER_REMAP_COEF_BITS ), dstptr + dst_index); -#else - #pragma unroll - for (int i = 0; i < 16; i++) - sum += v[i] * tab1y[(i>>2)] * tab1x[(i&3)]; + int dst_index = mad24(dy, dst_step, mad24(dx, pixsize, dst_offset)); storepix(CONVERT_TO_T( sum ), dstptr + dst_index); -#endif } } diff --git a/modules/imgproc/src/warp_common.hpp b/modules/imgproc/src/warp_common.hpp index 4af7006abd..fe1fd697d1 100644 --- a/modules/imgproc/src/warp_common.hpp +++ b/modules/imgproc/src/warp_common.hpp @@ -8,4 +8,15 @@ #include "warp_common.vector.hpp" #include "warp_common.scalar.hpp" +namespace cv { + +typedef void (*ImgWarpFunc)(const float* x, const float* y, int len, + const void* src, size_t srcstep, Size size, + void* dst, const float* params, + int borderType, const void* borderVal); + +ImgWarpFunc getImgWarpFunc(int type, int interpolation); + +} + #endif // __OPENCV_IMGPROC_WARP_COMMON_HPP__ diff --git a/modules/imgproc/src/warp_kernels.simd.hpp b/modules/imgproc/src/warp_kernels.simd.hpp index d48609656d..3b819beafc 100644 --- a/modules/imgproc/src/warp_kernels.simd.hpp +++ b/modules/imgproc/src/warp_kernels.simd.hpp @@ -106,6 +106,8 @@ namespace cv{ CV_CPU_OPTIMIZATION_NAMESPACE_BEGIN +ImgWarpFunc getBicubicWarpFunc_(int type); + void warpAffineNearestInvoker_8UC1(const uint8_t *src_data, size_t src_step, int src_rows, int src_cols, uint8_t *dst_data, size_t dst_step, int dst_rows, int dst_cols, const double M[6], int border_type, const double border_value[4]); @@ -6991,6 +6993,1693 @@ void remapLinearApproxInvoker_8UC4(const uint8_t *src_data, size_t src_step, int #endif } + +namespace { + +static CV_ALWAYS_INLINE void +bicubicWeights(float alpha, float A, float& w0, float& w1, float& w2, float& w3) +{ + float a2 = alpha*alpha; + float b = 1.f - alpha; + float b2 = b*b; + + w0 = A*alpha*b2; + w3 = A*a2*b; + w1 = a2*((A + 2.f)*alpha - (A + 3.f)) + 1.f; + w2 = 1.f - w0 - w1 - w3; +} + +static CV_ALWAYS_INLINE int +bicubicCoeffs(float xs0, float ys0, + size_t srcstep, Size size, int bpp, float A, + int bordertype, int& tl_x, int& tl_y, int& tl_ofs, + float& wx0, float& wx1, float& wx2, float& wx3, + float& wy0, float& wy1, float& wy2, float& wy3) +{ + constexpr int MIN_SIZE = 16; + int width = size.width, height = size.height; + int bigwidth = std::max(width, MIN_SIZE); + int bigheight = std::max(height, MIN_SIZE); + float minx = float(-bigwidth), maxx = float(bigwidth*2); + float miny = float(-bigheight), maxy = float(bigheight*2); + + // clamp coordinates in floating-point to try to avoid unpredictable + // behavior if the floating-point coordinates are huge + float vx0 = std::clamp(xs0, minx, maxx); + float vy0 = std::clamp(ys0, miny, maxy); + + int ix0 = cvFloor(vx0); + int iy0 = cvFloor(vy0); + float alpha = vx0 - ix0; + float beta = vy0 - iy0; + + ix0--; + iy0--; + + int all_outliers = int((unsigned)(ix0 + 4) >= (unsigned)(width + 4)) | + int((unsigned)(iy0 + 4) >= (unsigned)(height + 4)); + if (all_outliers && (bordertype == BORDER_CONSTANT || bordertype == BORDER_TRANSPARENT)) + return -1; + + int all_inliers = int((unsigned)ix0 < (unsigned)std::max(width - 3, 0)) & + int((unsigned)iy0 < (unsigned)std::max(height - 3, 0)); + + tl_x = ix0; + tl_y = iy0; + tl_ofs = iy0*(int)srcstep + ix0*bpp; + + bicubicWeights(alpha, A, wx0, wx1, wx2, wx3); + bicubicWeights(beta, A, wy0, wy1, wy2, wy3); + + return all_inliers; +} + +template +static void bicubicFetchPixels(const _Tp* src, size_t srcstep, Size size, int cn, + const int32_t* tl_x, const int32_t* tl_y, + int32_t* goodx, int row, _Fp* pixbuf, int len, + int borderType, const _Tp* defVal) +{ + int width = size.width, height = size.height; + srcstep /= sizeof(_Tp); + if (borderType == BORDER_CONSTANT || borderType == BORDER_TRANSPARENT) { + _Tp defR = defVal[0], defG = cn > 1 ? defVal[1] : _Tp(); + _Tp defB = cn > 2 ? defVal[2] : _Tp(), defA = cn > 3 ? defVal[3] : _Tp(); + for (int i = 0; i < len; i++) { + int x0 = tl_x[i], y0 = tl_y[i] + row; + int x1 = x0 + 1, x2 = x0 + 2, x3 = x0 + 3; + int my = int((unsigned)y0 < (unsigned)height); + int mx0 = int((unsigned)x0 < (unsigned)width) & my; + int mx1 = int((unsigned)x1 < (unsigned)width) & my; + int mx2 = int((unsigned)x2 < (unsigned)width) & my; + int mx3 = int((unsigned)x3 < (unsigned)width) & my; + y0 = std::clamp(y0, 0, height-1); + x0 = std::clamp(x0, 0, width-1)*cn; + x1 = std::clamp(x1, 0, width-1)*cn; + x2 = std::clamp(x2, 0, width-1)*cn; + x3 = std::clamp(x3, 0, width-1)*cn; + const _Tp* srcrow = src + srcstep*y0; + pixbuf[i] = (_Fp)(srcrow[x0]*mx0 + defR*(1 - mx0)); + pixbuf[i + len] = (_Fp)(srcrow[x1]*mx1 + defR*(1 - mx1)); + pixbuf[i + len*2] = (_Fp)(srcrow[x2]*mx2 + defR*(1 - mx2)); + pixbuf[i + len*3] = (_Fp)(srcrow[x3]*mx3 + defR*(1 - mx3)); + if (cn > 1) { + pixbuf[i + len*4] = (_Fp)(srcrow[x0 + 1]*mx0 + defG*(1 - mx0)); + pixbuf[i + len*5] = (_Fp)(srcrow[x1 + 1]*mx1 + defG*(1 - mx1)); + pixbuf[i + len*6] = (_Fp)(srcrow[x2 + 1]*mx2 + defG*(1 - mx2)); + pixbuf[i + len*7] = (_Fp)(srcrow[x3 + 1]*mx3 + defG*(1 - mx3)); + if (cn > 2) { + pixbuf[i + len*8] = (_Fp)(srcrow[x0 + 2]*mx0 + defB*(1 - mx0)); + pixbuf[i + len*9] = (_Fp)(srcrow[x1 + 2]*mx1 + defB*(1 - mx1)); + pixbuf[i + len*10] = (_Fp)(srcrow[x2 + 2]*mx2 + defB*(1 - mx2)); + pixbuf[i + len*11] = (_Fp)(srcrow[x3 + 2]*mx3 + defB*(1 - mx3)); + if (cn > 3) { + pixbuf[i + len*12] = (_Fp)(srcrow[x0 + 3]*mx0 + defA*(1 - mx0)); + pixbuf[i + len*13] = (_Fp)(srcrow[x1 + 3]*mx1 + defA*(1 - mx1)); + pixbuf[i + len*14] = (_Fp)(srcrow[x2 + 3]*mx2 + defA*(1 - mx2)); + pixbuf[i + len*15] = (_Fp)(srcrow[x3 + 3]*mx3 + defA*(1 - mx3)); + } + } + } + } + } else { + for (int i = 0; i < len; i++) { + int y0 = borderInterpolate_fast(tl_y[i] + row, height, borderType); + int x0, x1, x2, x3; + if (row == 0) { + int x0_ = tl_x[i]; + if ((unsigned)x0_ < (unsigned)std::max(width - 3, 0)) { + x0 = x0_*cn; + x1 = (x0_ + 1)*cn; + x2 = (x0_ + 2)*cn; + x3 = (x0_ + 3)*cn; + } else { + x0 = borderInterpolate_fast(x0_, width, borderType)*cn; + x1 = borderInterpolate_fast(x0_ + 1, width, borderType)*cn; + x2 = borderInterpolate_fast(x0_ + 2, width, borderType)*cn; + x3 = borderInterpolate_fast(x0_ + 3, width, borderType)*cn; + } + goodx[i*4] = x0; + goodx[i*4 + 1] = x1; + goodx[i*4 + 2] = x2; + goodx[i*4 + 3] = x3; + } else { + x0 = goodx[i*4]; + x1 = goodx[i*4 + 1]; + x2 = goodx[i*4 + 2]; + x3 = goodx[i*4 + 3]; + } + const _Tp* srcrow = src + srcstep*y0; + pixbuf[i] = (_Fp)srcrow[x0]; + pixbuf[i + len] = (_Fp)srcrow[x1]; + pixbuf[i + len*2] = (_Fp)srcrow[x2]; + pixbuf[i + len*3] = (_Fp)srcrow[x3]; + if (cn > 1) { + pixbuf[i + len*4] = (_Fp)srcrow[x0 + 1]; + pixbuf[i + len*5] = (_Fp)srcrow[x1 + 1]; + pixbuf[i + len*6] = (_Fp)srcrow[x2 + 1]; + pixbuf[i + len*7] = (_Fp)srcrow[x3 + 1]; + if (cn > 2) { + pixbuf[i + len*8] = (_Fp)srcrow[x0 + 2]; + pixbuf[i + len*9] = (_Fp)srcrow[x1 + 2]; + pixbuf[i + len*10] = (_Fp)srcrow[x2 + 2]; + pixbuf[i + len*11] = (_Fp)srcrow[x3 + 2]; + if (cn > 3) { + pixbuf[i + len*12] = (_Fp)srcrow[x0 + 3]; + pixbuf[i + len*13] = (_Fp)srcrow[x1 + 3]; + pixbuf[i + len*14] = (_Fp)srcrow[x2 + 3]; + pixbuf[i + len*15] = (_Fp)srcrow[x3 + 3]; + } + } + } + } + } +} + +#undef BICUBIC_UPDATE_ACC +#define BICUBIC_UPDATE_ACC(acc, v0, v1, v2, v3, wy) \ + acc += (v0*wx0 + v1*wx1 + v2*wx2 + v3*wx3)*wy + +#undef BICUBIC_UPDATE_ACC_VEC +#define BICUBIC_UPDATE_ACC_VEC(acc, v0, v1, v2, v3, wy) \ + sumwx = v_fma(v1, wx1, v_mul(v0, wx0)); \ + sumwx = v_fma(v2, wx2, sumwx); \ + sumwx = v_fma(v3, wx3, sumwx); \ + acc = v_fma(sumwx, wy, acc) + +#undef BICUBIC_PROCESS_INLIERS_C1_SCALAR +#define BICUBIC_PROCESS_INLIERS_C1_SCALAR(row) \ + srcrow = (const chtype*)((const uint8_t*)src + row*srcstep + tl_ofs); \ + V0 = float(srcrow[0]); \ + V1 = float(srcrow[1]); \ + V2 = float(srcrow[2]); \ + V3 = float(srcrow[3]); \ + BICUBIC_UPDATE_ACC(acc_r, V0, V1, V2, V3, wy##row) + +#undef BICUBIC_PROCESS_INLIERS_C2_SCALAR +#define BICUBIC_PROCESS_INLIERS_C2_SCALAR(row) \ + srcrow = (const chtype*)((const uint8_t*)src + row*srcstep + tl_ofs); \ + V0 = float(srcrow[0]); \ + V1 = float(srcrow[2]); \ + V2 = float(srcrow[4]); \ + V3 = float(srcrow[6]); \ + BICUBIC_UPDATE_ACC(acc_r, V0, V1, V2, V3, wy##row); \ + V0 = float(srcrow[1]); \ + V1 = float(srcrow[3]); \ + V2 = float(srcrow[5]); \ + V3 = float(srcrow[7]); \ + BICUBIC_UPDATE_ACC(acc_g, V0, V1, V2, V3, wy##row) + +#undef BICUBIC_PROCESS_INLIERS_C3_SCALAR +#define BICUBIC_PROCESS_INLIERS_C3_SCALAR(row) \ + srcrow = (const chtype*)((const uint8_t*)src + row*srcstep + tl_ofs); \ + V0 = float(srcrow[0]); \ + V1 = float(srcrow[3]); \ + V2 = float(srcrow[6]); \ + V3 = float(srcrow[9]); \ + BICUBIC_UPDATE_ACC(acc_r, V0, V1, V2, V3, wy##row); \ + V0 = float(srcrow[1]); \ + V1 = float(srcrow[4]); \ + V2 = float(srcrow[7]); \ + V3 = float(srcrow[10]); \ + BICUBIC_UPDATE_ACC(acc_g, V0, V1, V2, V3, wy##row); \ + V0 = float(srcrow[2]); \ + V1 = float(srcrow[5]); \ + V2 = float(srcrow[8]); \ + V3 = float(srcrow[11]); \ + BICUBIC_UPDATE_ACC(acc_b, V0, V1, V2, V3, wy##row) + +#undef BICUBIC_PROCESS_INLIERS_C4_SCALAR +#define BICUBIC_PROCESS_INLIERS_C4_SCALAR(row) \ + srcrow = (const chtype*)((const uint8_t*)src + row*srcstep + tl_ofs); \ + V0 = float(srcrow[0]); \ + V1 = float(srcrow[4]); \ + V2 = float(srcrow[8]); \ + V3 = float(srcrow[12]); \ + BICUBIC_UPDATE_ACC(acc_r, V0, V1, V2, V3, wy##row); \ + V0 = float(srcrow[1]); \ + V1 = float(srcrow[5]); \ + V2 = float(srcrow[9]); \ + V3 = float(srcrow[13]); \ + BICUBIC_UPDATE_ACC(acc_g, V0, V1, V2, V3, wy##row); \ + V0 = float(srcrow[2]); \ + V1 = float(srcrow[6]); \ + V2 = float(srcrow[10]); \ + V3 = float(srcrow[14]); \ + BICUBIC_UPDATE_ACC(acc_b, V0, V1, V2, V3, wy##row); \ + V0 = float(srcrow[3]); \ + V1 = float(srcrow[7]); \ + V2 = float(srcrow[11]); \ + V3 = float(srcrow[15]); \ + BICUBIC_UPDATE_ACC(acc_a, V0, V1, V2, V3, wy##row) + +template +static void bicubicRef(const float* srcx, const float* srcy, int len, + const void* src, size_t srcstep, Size size, + chtype* dst, const float* params, int borderType, chtype* borderVal) +{ + constexpr float defaultA = -0.75f; + float A = params ? *params : defaultA; + + constexpr int BPP = int(NCHANNELS*sizeof(dst[0])); + using buftype = std::conditional_t, float, int>; + using pixtype = std::conditional_t>; + + pixtype bval = borderVal ? *(pixtype*)borderVal : pixtype(); + + for (int i = 0; i < len; i++, dst += NCHANNELS) { + buftype pixbuf[NCHANNELS][4]; + int tl_x, tl_y, tl_ofs, goodx[4]; + float V0, V1, V2, V3; + float wx0, wx1, wx2, wx3, wy0, wy1, wy2, wy3; + float xs = srcx[i], ys = srcy[i]; + + int code = bicubicCoeffs(xs, ys, srcstep, size, BPP, A, + borderType, tl_x, tl_y, tl_ofs, + wx0, wx1, wx2, wx3, wy0, wy1, wy2, wy3); + if (code > 0) { + const chtype* srcrow; + if constexpr (NCHANNELS == 1) { + float acc_r = 0.f; + BICUBIC_PROCESS_INLIERS_C1_SCALAR(0); + BICUBIC_PROCESS_INLIERS_C1_SCALAR(1); + BICUBIC_PROCESS_INLIERS_C1_SCALAR(2); + BICUBIC_PROCESS_INLIERS_C1_SCALAR(3); + dst[0] = saturate_cast(acc_r); + } else if constexpr (NCHANNELS == 2) { + float acc_r = 0.f, acc_g = 0.f; + BICUBIC_PROCESS_INLIERS_C2_SCALAR(0); + BICUBIC_PROCESS_INLIERS_C2_SCALAR(1); + BICUBIC_PROCESS_INLIERS_C2_SCALAR(2); + BICUBIC_PROCESS_INLIERS_C2_SCALAR(3); + dst[0] = saturate_cast(acc_r); + dst[1] = saturate_cast(acc_g); + } else if constexpr (NCHANNELS == 3) { + float acc_r = 0.f, acc_g = 0.f, acc_b = 0.f; + BICUBIC_PROCESS_INLIERS_C3_SCALAR(0); + BICUBIC_PROCESS_INLIERS_C3_SCALAR(1); + BICUBIC_PROCESS_INLIERS_C3_SCALAR(2); + BICUBIC_PROCESS_INLIERS_C3_SCALAR(3); + dst[0] = saturate_cast(acc_r); + dst[1] = saturate_cast(acc_g); + dst[2] = saturate_cast(acc_b); + } else if constexpr (NCHANNELS == 4) { + float acc_r = 0.f, acc_g = 0.f, acc_b = 0.f, acc_a = 0.f; + BICUBIC_PROCESS_INLIERS_C4_SCALAR(0); + BICUBIC_PROCESS_INLIERS_C4_SCALAR(1); + BICUBIC_PROCESS_INLIERS_C4_SCALAR(2); + BICUBIC_PROCESS_INLIERS_C4_SCALAR(3); + dst[0] = saturate_cast(acc_r); + dst[1] = saturate_cast(acc_g); + dst[2] = saturate_cast(acc_b); + dst[3] = saturate_cast(acc_a); + } + } else if (code < 0) { + if (borderType == BORDER_CONSTANT) { + *((pixtype*)dst) = bval; + } + continue; + } else { + float acc[NCHANNELS] = {}; + float wys[] = {wy0, wy1, wy2, wy3}; + const chtype* defVal = borderType == BORDER_TRANSPARENT ? dst : borderVal; + for (int row = 0; row < 4; row++) { + bicubicFetchPixels((const chtype*)src, srcstep, size, NCHANNELS, &tl_x, &tl_y, + goodx, row, &pixbuf[0][0], 1, borderType, defVal); + float wy = wys[row]; + for (int c = 0; c < NCHANNELS; c++) { + V0 = float(pixbuf[c][0]); + V1 = float(pixbuf[c][1]); + V2 = float(pixbuf[c][2]); + V3 = float(pixbuf[c][3]); + BICUBIC_UPDATE_ACC(acc[c], V0, V1, V2, V3, wy); + } + } + for (int c = 0; c < NCHANNELS; c++) { + dst[c] = saturate_cast(acc[c]); + } + } + } +} + +#if CV_SIMD + +CV_ALWAYS_INLINE void +bicubicWeights(const v_float32& alpha, float A, + v_float32& w0, v_float32& w1, + v_float32& w2, v_float32& w3) +{ + const v_float32 vA = vx_setall_f32(A); + const v_float32 vAp2 = vx_setall_f32(A + 2.0f); + const v_float32 vAp3 = vx_setall_f32(-(A + 3.0f)); + const v_float32 v1 = vx_setall_f32(1.0f); + + const v_float32 a2 = v_mul(alpha, alpha); // α² + const v_float32 b = v_sub(v1, alpha); // b = 1-α + const v_float32 b2 = v_mul(b, b); // b² + + w0 = v_mul(vA, v_mul(alpha, b2)); // A·α·b² + w3 = v_mul(vA, v_mul(a2, b)); // A·α²·b + w1 = v_fma(a2, v_fma(vAp2, alpha, vAp3), v1); // a²·((A+2)α-(A+3))+1 + w2 = v_sub(v_sub(v_sub(v1, w0), w1), w3); +} + +static CV_ALWAYS_INLINE int +bicubicCoeffs(const float* srcx, const float* srcy, + size_t srcstep, Size size, int bpp, float A, + int bordertype, int32_t* tl_x, int32_t* tl_y, int32_t* tl_ofs, + v_float32& wx0, v_float32& wx1, v_float32& wx2, v_float32& wx3, + v_float32& wy0, v_float32& wy1, v_float32& wy2, v_float32& wy3) +{ + constexpr int MIN_SIZE = 16; + int width = size.width, height = size.height; + int bigwidth = std::max(width, MIN_SIZE); + int bigheight = std::max(height, MIN_SIZE); + v_float32 minx = vx_setall_f32(float(-bigwidth)), maxx = vx_setall_f32(float(bigwidth*2)); + v_float32 miny = vx_setall_f32(float(-bigheight)), maxy = vx_setall_f32(float(bigheight*2)); + + // clamp coordinates in floating-point to try to avoid unpredictable + // behavior if the floating-point coordinates are huge + v_float32 xs0 = vx_load(srcx), ys0 = vx_load(srcy); + v_float32 vx0 = v_min(v_max(xs0, minx), maxx); + v_float32 vy0 = v_min(v_max(ys0, miny), maxy); + + v_int32 ix0 = v_floor(vx0); + v_int32 iy0 = v_floor(vy0); + v_float32 alpha = v_sub(vx0, v_cvt_f32(ix0)); + v_float32 beta = v_sub(vy0, v_cvt_f32(iy0)); + + v_int32 one_i = vx_setall_s32(1), four_i = vx_setall_s32(4); + ix0 = v_sub(ix0, one_i); + iy0 = v_sub(iy0, one_i); + + v_uint32 width_outer = vx_setall_u32((uint32_t)(width + 4)); + v_uint32 height_outer = vx_setall_u32((uint32_t)(height + 4)); + + v_uint32 outliers_mask = v_or(v_ge(v_reinterpret_as_u32(v_add(ix0, four_i)), width_outer), + v_ge(v_reinterpret_as_u32(v_add(iy0, four_i)), height_outer)); + + bool all_outliers = v_check_all(outliers_mask); + if (all_outliers && (bordertype == BORDER_CONSTANT || bordertype == BORDER_TRANSPARENT)) + return -1; + + v_uint32 width_inner = vx_setall_u32((uint32_t)std::max(width - 3, 0)); + v_uint32 height_inner = vx_setall_u32((uint32_t)std::max(height - 5, 0)); + + v_uint32 inliers_mask = v_and(v_lt(v_reinterpret_as_u32(ix0), width_inner), + v_lt(v_reinterpret_as_u32(iy0), height_inner)); + + bool all_inliers = v_check_all(inliers_mask); + v_int32 tl_ofs0 = v_add(v_mul(iy0, vx_setall_s32((int)srcstep)), v_mul(ix0, vx_setall_s32(bpp))); + + v_store(tl_x, ix0); + v_store(tl_y, iy0); + v_store(tl_ofs, tl_ofs0); + + bicubicWeights(alpha, A, wx0, wx1, wx2, wx3); + bicubicWeights(beta, A, wy0, wy1, wy2, wy3); + + return int(all_inliers); +} + +#if CV_SIMD_FP16 +CV_ALWAYS_INLINE void +bicubicWeights(const v_float16& alpha, float A, + v_float16& w0, v_float16& w1, + v_float16& w2, v_float16& w3) +{ + const v_float16 vA = vx_setall_f16(hfloat(A)); + const v_float16 vAp2 = vx_setall_f16(hfloat(A + 2.0f)); + const v_float16 vAp3 = vx_setall_f16(hfloat(-(A + 3.0f))); + const v_float16 v1 = vx_setall_f16(hfloat(1.0f)); + + const v_float16 a2 = v_mul(alpha, alpha); // α² + const v_float16 b = v_sub(v1, alpha); // b = 1-α + const v_float16 b2 = v_mul(b, b); // b² + + w0 = v_mul(vA, v_mul(alpha, b2)); // A·α·b² + w3 = v_mul(vA, v_mul(a2, b)); // A·α²·b + w1 = v_fma(a2, v_fma(vAp2, alpha, vAp3), v1); // a²·((A+2)α-(A+3))+1 + w2 = v_sub(v_sub(v_sub(v1, w0), w1), w3); +} + +static CV_ALWAYS_INLINE int +bicubicCoeffs(const float* srcx, const float* srcy, + size_t srcstep, Size size, int bpp, float A, + int bordertype, int32_t* tl_x, int32_t* tl_y, int32_t* tl_ofs, + v_float16& wx0, v_float16& wx1, v_float16& wx2, v_float16& wx3, + v_float16& wy0, v_float16& wy1, v_float16& wy2, v_float16& wy3) +{ + constexpr int nlanes32 = VTraits::nlanes; + constexpr int MIN_SIZE = 16; + int width = size.width, height = size.height; + int bigwidth = std::max(width, MIN_SIZE); + int bigheight = std::max(height, MIN_SIZE); + v_float32 minx = vx_setall_f32(float(-bigwidth)), maxx = vx_setall_f32(float(bigwidth*2)); + v_float32 miny = vx_setall_f32(float(-bigheight)), maxy = vx_setall_f32(float(bigheight*2)); + + // clamp coordinates in floating-point to try to avoid unpredictable + // behavior if the floating-point coordinates are huge + v_float32 xs0 = vx_load(srcx), ys0 = vx_load(srcy); + v_float32 xs1 = vx_load(srcx + nlanes32), ys1 = vx_load(srcy + nlanes32); + v_float32 vx0 = v_min(v_max(xs0, minx), maxx); + v_float32 vy0 = v_min(v_max(ys0, miny), maxy); + v_float32 vx1 = v_min(v_max(xs1, minx), maxx); + v_float32 vy1 = v_min(v_max(ys1, miny), maxy); + + v_int32 ix0 = v_floor(vx0); + v_int32 iy0 = v_floor(vy0); + v_int32 ix1 = v_floor(vx1); + v_int32 iy1 = v_floor(vy1); + v_float32 alpha0 = v_sub(vx0, v_cvt_f32(ix0)); + v_float32 beta0 = v_sub(vy0, v_cvt_f32(iy0)); + v_float32 alpha1 = v_sub(vx1, v_cvt_f32(ix1)); + v_float32 beta1 = v_sub(vy1, v_cvt_f32(iy1)); + + hfloat abuf[nlanes32*2], bbuf[nlanes32*2]; + + v_pack_store(abuf, alpha0); + v_pack_store(abuf + nlanes32, alpha1); + v_pack_store(bbuf, beta0); + v_pack_store(bbuf + nlanes32, beta1); + + v_float16 alpha = vx_load(abuf); + v_float16 beta = vx_load(bbuf); + + v_int32 one_i = vx_setall_s32(1), four_i = vx_setall_s32(4); + ix0 = v_sub(ix0, one_i); + iy0 = v_sub(iy0, one_i); + ix1 = v_sub(ix1, one_i); + iy1 = v_sub(iy1, one_i); + + v_uint32 width_outer = vx_setall_u32((uint32_t)(width + 4)); + v_uint32 height_outer = vx_setall_u32((uint32_t)(height + 4)); + + v_uint32 outliers_mask0 = v_or(v_ge(v_reinterpret_as_u32(v_add(ix0, four_i)), width_outer), + v_ge(v_reinterpret_as_u32(v_add(iy0, four_i)), height_outer)); + v_uint32 outliers_mask1 = v_or(v_ge(v_reinterpret_as_u32(v_add(ix1, four_i)), width_outer), + v_ge(v_reinterpret_as_u32(v_add(iy1, four_i)), height_outer)); + + bool all_outliers = v_check_all(v_and(outliers_mask0, outliers_mask1)); + if (all_outliers && (bordertype == BORDER_CONSTANT || bordertype == BORDER_TRANSPARENT)) + return -1; + + v_uint32 width_inner = vx_setall_u32((uint32_t)std::max(width - 3, 0)); + v_uint32 height_inner = vx_setall_u32((uint32_t)std::max(height - 5, 0)); + + v_uint32 inliers_mask0 = v_and(v_lt(v_reinterpret_as_u32(ix0), width_inner), + v_lt(v_reinterpret_as_u32(iy0), height_inner)); + v_uint32 inliers_mask1 = v_and(v_lt(v_reinterpret_as_u32(ix1), width_inner), + v_lt(v_reinterpret_as_u32(iy1), height_inner)); + + bool all_inliers = v_check_all(v_and(inliers_mask0, inliers_mask1)); + v_int32 tl_ofs0 = v_add(v_mul(iy0, vx_setall_s32((int)srcstep)), v_mul(ix0, vx_setall_s32(bpp))); + v_int32 tl_ofs1 = v_add(v_mul(iy1, vx_setall_s32((int)srcstep)), v_mul(ix1, vx_setall_s32(bpp))); + + v_store(tl_x, ix0); + v_store(tl_y, iy0); + v_store(tl_ofs, tl_ofs0); + v_store(tl_x + nlanes32, ix1); + v_store(tl_y + nlanes32, iy1); + v_store(tl_ofs + nlanes32, tl_ofs1); + + bicubicWeights(alpha, A, wx0, wx1, wx2, wx3); + bicubicWeights(beta, A, wy0, wy1, wy2, wy3); + + return int(all_inliers); +} + +CV_ALWAYS_INLINE v_float16 v_mul(const v_int16& a, const v_float16& b) +{ + return v_mul(v_cvt_f16(a), b); +} + +CV_ALWAYS_INLINE v_float16 v_fma(const v_int16& a, const v_float16& b, const v_float16& c) +{ + return v_fma(v_cvt_f16(a), b, c); +} + +#endif + +CV_ALWAYS_INLINE v_float32 v_mul(const v_int32& a, const v_float32& b) +{ + return v_mul(v_cvt_f32(a), b); +} + +CV_ALWAYS_INLINE v_float32 v_fma(const v_int32& a, const v_float32& b, const v_float32& c) +{ + return v_fma(v_cvt_f32(a), b, c); +} + +#undef FETCH_INLIERS_C1_DEFAULT +#define FETCH_INLIERS_C1_DEFAULT(row) \ + for (int j = 0; j < BATCH; j++) { \ + const chtype* srcj = (const chtype*)((const uint8_t*)src + row*srcstep + tl_ofs[j]); \ + pixbuf[0][j] = buftype(srcj[1]); \ + pixbuf[0][j + BATCH] = buftype(srcj[2]); \ + pixbuf[0][j + BATCH*2] = buftype(srcj[3]); \ + pixbuf[0][j + BATCH*3] = buftype(srcj[4]); \ + } \ + R0 = vx_load(&pixbuf[0][0]); \ + R1 = vx_load(&pixbuf[0][BATCH]); \ + R2 = vx_load(&pixbuf[0][BATCH*2]); \ + R3 = vx_load(&pixbuf[0][BATCH*3]) + +#undef FETCH_INLIERS_C3_DEFAULT +#define FETCH_INLIERS_C3_DEFAULT(row) \ + for (int j = 0; j < BATCH; j++) { \ + const chtype* srcj = (const chtype*)((const uint8_t*)src + row*srcstep + tl_ofs[j]); \ + pixbuf[0][j] = buftype(srcj[0]); \ + pixbuf[0][j + BATCH] = buftype(srcj[3]); \ + pixbuf[0][j + BATCH*2] = buftype(srcj[6]); \ + pixbuf[0][j + BATCH*3] = buftype(srcj[9]); \ + pixbuf[1][j] = buftype(srcj[1]); \ + pixbuf[1][j + BATCH] = buftype(srcj[4]); \ + pixbuf[1][j + BATCH*2] = buftype(srcj[7]); \ + pixbuf[1][j + BATCH*3] = buftype(srcj[10]); \ + pixbuf[2][j] = buftype(srcj[2]); \ + pixbuf[2][j + BATCH] = buftype(srcj[5]); \ + pixbuf[2][j + BATCH*2] = buftype(srcj[8]); \ + pixbuf[2][j + BATCH*3] = buftype(srcj[11]); \ + } \ + R0 = vx_load(&pixbuf[0][0]); \ + R1 = vx_load(&pixbuf[0][BATCH]); \ + R2 = vx_load(&pixbuf[0][BATCH*2]); \ + R3 = vx_load(&pixbuf[0][BATCH*3]); \ + G0 = vx_load(&pixbuf[1][0]); \ + G1 = vx_load(&pixbuf[1][BATCH]); \ + G2 = vx_load(&pixbuf[1][BATCH*2]); \ + G3 = vx_load(&pixbuf[1][BATCH*3]); \ + B0 = vx_load(&pixbuf[2][0]); \ + B1 = vx_load(&pixbuf[2][BATCH]); \ + B2 = vx_load(&pixbuf[2][BATCH*2]); \ + B3 = vx_load(&pixbuf[2][BATCH*3]) + +#undef FETCH_INLIERS_C4_DEFAULT +#define FETCH_INLIERS_C4_DEFAULT(row) \ + for (int j = 0; j < BATCH; j++) { \ + const chtype* srcj = (const chtype*)((const uint8_t*)src + row*srcstep + tl_ofs[j]); \ + pixbuf[0][j] = buftype(srcj[0]); \ + pixbuf[0][j + BATCH] = buftype(srcj[4]); \ + pixbuf[0][j + BATCH*2] = buftype(srcj[8]); \ + pixbuf[0][j + BATCH*3] = buftype(srcj[12]); \ + pixbuf[1][j] = buftype(srcj[1]); \ + pixbuf[1][j + BATCH] = buftype(srcj[5]); \ + pixbuf[1][j + BATCH*2] = buftype(srcj[9]); \ + pixbuf[1][j + BATCH*3] = buftype(srcj[13]); \ + pixbuf[2][j] = buftype(srcj[2]); \ + pixbuf[2][j + BATCH] = buftype(srcj[6]); \ + pixbuf[2][j + BATCH*2] = buftype(srcj[10]); \ + pixbuf[2][j + BATCH*3] = buftype(srcj[14]); \ + pixbuf[3][j] = buftype(srcj[3]); \ + pixbuf[3][j + BATCH] = buftype(srcj[7]); \ + pixbuf[3][j + BATCH*2] = buftype(srcj[11]); \ + pixbuf[3][j + BATCH*3] = buftype(srcj[15]); \ + } \ + R0 = vx_load(&pixbuf[0][0]); \ + R1 = vx_load(&pixbuf[0][BATCH]); \ + R2 = vx_load(&pixbuf[0][BATCH*2]); \ + R3 = vx_load(&pixbuf[0][BATCH*3]); \ + G0 = vx_load(&pixbuf[1][0]); \ + G1 = vx_load(&pixbuf[1][BATCH]); \ + G2 = vx_load(&pixbuf[1][BATCH*2]); \ + G3 = vx_load(&pixbuf[1][BATCH*3]); \ + B0 = vx_load(&pixbuf[2][0]); \ + B1 = vx_load(&pixbuf[2][BATCH]); \ + B2 = vx_load(&pixbuf[2][BATCH*2]); \ + B3 = vx_load(&pixbuf[2][BATCH*3]); \ + A0 = vx_load(&pixbuf[3][0]); \ + A1 = vx_load(&pixbuf[3][BATCH]); \ + A2 = vx_load(&pixbuf[3][BATCH*2]); \ + A3 = vx_load(&pixbuf[3][BATCH*3]) + +#undef FETCH_INLIERS_8UC1 +#define FETCH_INLIERS_8UC1(row) \ + FETCH_INLIERS_C1_DEFAULT(row) + +#undef FETCH_INLIERS_8UC3 +#define FETCH_INLIERS_8UC3(row) \ + FETCH_INLIERS_C3_DEFAULT(row) + +#undef FETCH_INLIERS_8UC4 +#define FETCH_INLIERS_8UC4(row) \ + FETCH_INLIERS_C4_DEFAULT(row) + +#undef FETCH_INLIERS_16UC1 +#define FETCH_INLIERS_16UC1(row) \ + FETCH_INLIERS_C1_DEFAULT(row) + +#undef FETCH_INLIERS_16UC3 +#define FETCH_INLIERS_16UC3(row) \ + FETCH_INLIERS_C3_DEFAULT(row) + +#undef FETCH_INLIERS_16UC4 +#define FETCH_INLIERS_16UC4(row) \ + FETCH_INLIERS_C4_DEFAULT(row) + +#undef FETCH_INLIERS_32FC1 +#define FETCH_INLIERS_32FC1(row) \ + FETCH_INLIERS_C1_DEFAULT(row) + +#undef FETCH_INLIERS_32FC3 +#define FETCH_INLIERS_32FC3(row) \ + FETCH_INLIERS_C3_DEFAULT(row) + +#undef FETCH_INLIERS_32FC4 +#define FETCH_INLIERS_32FC4(row) \ + FETCH_INLIERS_C4_DEFAULT(row) + +// NEON-optimized macros for super-fast pixels retrieval and reordering +#if defined __ARM_NEON && defined __aarch64__ + +#if CV_SIMD_FP16 + +#undef FETCH_INLIERS_8UC1 +#define FETCH_INLIERS_8UC1(row) { \ + const uint8_t* srcrow = (const uint8_t*)src + row*srcstep; \ + uint8x8_t _r0 = vld1_u8(srcrow + tl_ofs[0]); \ + uint8x8_t _r1 = vld1_u8(srcrow + tl_ofs[1]); \ + uint8x8_t _r2 = vld1_u8(srcrow + tl_ofs[2]); \ + uint8x8_t _r3 = vld1_u8(srcrow + tl_ofs[3]); \ + uint8x8_t _r4 = vld1_u8(srcrow + tl_ofs[4]); \ + uint8x8_t _r5 = vld1_u8(srcrow + tl_ofs[5]); \ + uint8x8_t _r6 = vld1_u8(srcrow + tl_ofs[6]); \ + uint8x8_t _r7 = vld1_u8(srcrow + tl_ofs[7]); \ + \ + uint8x16_t _r0415 = vcombine_u8(vzip1_u8(_r0, _r4), vzip1_u8(_r1, _r5)); \ + uint8x16_t _r2637 = vcombine_u8(vzip1_u8(_r2, _r6), vzip1_u8(_r3, _r7)); \ + uint8x16_t _r0246 = vzip1q_u8(_r0415, _r2637); \ + uint8x16_t _r1357 = vzip2q_u8(_r0415, _r2637); \ + uint8x16_t _r_c01 = vzip1q_u8(_r0246, _r1357); \ + uint8x16_t _r_c23 = vzip2q_u8(_r0246, _r1357); \ + \ + R0.val = vreinterpretq_s16_u16(vmovl_u8(vget_low_u8(_r_c01))); \ + R1.val = vreinterpretq_s16_u16(vmovl_high_u8(_r_c01)); \ + R2.val = vreinterpretq_s16_u16(vmovl_u8(vget_low_u8(_r_c23))); \ + R3.val = vreinterpretq_s16_u16(vmovl_high_u8(_r_c23)); } + +#undef FETCH_INLIERS_8UC3 +#define FETCH_INLIERS_8UC3(row) { \ + const uint8x16_t mask_rgb = { 0, 3, 6, 9, 1, 4, 7, 10, 2, 5, 8, 11, 255, 255, 255, 255 }; \ + const uint8_t* srcrow = (const uint8_t*)src + row*srcstep; \ + uint8x16_t _r0 = vqtbl1q_u8(vld1q_u8((const chtype*)(srcrow + tl_ofs[0])), mask_rgb); \ + uint8x16_t _r1 = vqtbl1q_u8(vld1q_u8((const chtype*)(srcrow + tl_ofs[1])), mask_rgb); \ + uint8x16_t _r2 = vqtbl1q_u8(vld1q_u8((const chtype*)(srcrow + tl_ofs[2])), mask_rgb); \ + uint8x16_t _r3 = vqtbl1q_u8(vld1q_u8((const chtype*)(srcrow + tl_ofs[3])), mask_rgb); \ + uint8x16_t _r4 = vqtbl1q_u8(vld1q_u8((const chtype*)(srcrow + tl_ofs[4])), mask_rgb); \ + uint8x16_t _r5 = vqtbl1q_u8(vld1q_u8((const chtype*)(srcrow + tl_ofs[5])), mask_rgb); \ + uint8x16_t _r6 = vqtbl1q_u8(vld1q_u8((const chtype*)(srcrow + tl_ofs[6])), mask_rgb); \ + uint8x16_t _r7 = vqtbl1q_u8(vld1q_u8((const chtype*)(srcrow + tl_ofs[7])), mask_rgb); \ + \ + uint8x16_t _z04a = vzip1q_u8(_r0, _r4); \ + uint8x16_t _z15a = vzip1q_u8(_r1, _r5); \ + uint8x16_t _z26a = vzip1q_u8(_r2, _r6); \ + uint8x16_t _z37a = vzip1q_u8(_r3, _r7); \ + uint8x16_t _z04b = vzip2q_u8(_r0, _r4); \ + uint8x16_t _z15b = vzip2q_u8(_r1, _r5); \ + uint8x16_t _z26b = vzip2q_u8(_r2, _r6); \ + uint8x16_t _z37b = vzip2q_u8(_r3, _r7); \ + \ + uint8x16_t _r0246 = vzip1q_u8(_z04a, _z26a); \ + uint8x16_t _g0246 = vzip2q_u8(_z04a, _z26a); \ + uint8x16_t _r1357 = vzip1q_u8(_z15a, _z37a); \ + uint8x16_t _g1357 = vzip2q_u8(_z15a, _z37a); \ + uint8x16_t _b0246 = vzip1q_u8(_z04b, _z26b); \ + uint8x16_t _b1357 = vzip1q_u8(_z15b, _z37b); \ + \ + uint8x16_t _r_c01 = vzip1q_u8(_r0246, _r1357); \ + uint8x16_t _r_c23 = vzip2q_u8(_r0246, _r1357); \ + uint8x16_t _g_c01 = vzip1q_u8(_g0246, _g1357); \ + uint8x16_t _g_c23 = vzip2q_u8(_g0246, _g1357); \ + uint8x16_t _b_c01 = vzip1q_u8(_b0246, _b1357); \ + uint8x16_t _b_c23 = vzip2q_u8(_b0246, _b1357); \ + \ + R0.val = vreinterpretq_s16_u16(vmovl_u8(vget_low_u8(_r_c01))); \ + R1.val = vreinterpretq_s16_u16(vmovl_high_u8(_r_c01)); \ + R2.val = vreinterpretq_s16_u16(vmovl_u8(vget_low_u8(_r_c23))); \ + R3.val = vreinterpretq_s16_u16(vmovl_high_u8(_r_c23)); \ + G0.val = vreinterpretq_s16_u16(vmovl_u8(vget_low_u8(_g_c01))); \ + G1.val = vreinterpretq_s16_u16(vmovl_high_u8(_g_c01)); \ + G2.val = vreinterpretq_s16_u16(vmovl_u8(vget_low_u8(_g_c23))); \ + G3.val = vreinterpretq_s16_u16(vmovl_high_u8(_g_c23)); \ + B0.val = vreinterpretq_s16_u16(vmovl_u8(vget_low_u8(_b_c01))); \ + B1.val = vreinterpretq_s16_u16(vmovl_high_u8(_b_c01)); \ + B2.val = vreinterpretq_s16_u16(vmovl_u8(vget_low_u8(_b_c23))); \ + B3.val = vreinterpretq_s16_u16(vmovl_high_u8(_b_c23)); } + +#undef FETCH_INLIERS_8UC4 +#define FETCH_INLIERS_8UC4(row) { \ + const uint8x16_t mask_rgb = { 0, 4, 8, 12, 1, 5, 9, 13, 2, 6, 10, 14, 3, 7, 11, 15 }; \ + const uint8_t* srcrow = (const uint8_t*)src + row*srcstep; \ + uint8x16_t _r0 = vqtbl1q_u8(vld1q_u8((const chtype*)(srcrow + tl_ofs[0])), mask_rgb); \ + uint8x16_t _r1 = vqtbl1q_u8(vld1q_u8((const chtype*)(srcrow + tl_ofs[1])), mask_rgb); \ + uint8x16_t _r2 = vqtbl1q_u8(vld1q_u8((const chtype*)(srcrow + tl_ofs[2])), mask_rgb); \ + uint8x16_t _r3 = vqtbl1q_u8(vld1q_u8((const chtype*)(srcrow + tl_ofs[3])), mask_rgb); \ + uint8x16_t _r4 = vqtbl1q_u8(vld1q_u8((const chtype*)(srcrow + tl_ofs[4])), mask_rgb); \ + uint8x16_t _r5 = vqtbl1q_u8(vld1q_u8((const chtype*)(srcrow + tl_ofs[5])), mask_rgb); \ + uint8x16_t _r6 = vqtbl1q_u8(vld1q_u8((const chtype*)(srcrow + tl_ofs[6])), mask_rgb); \ + uint8x16_t _r7 = vqtbl1q_u8(vld1q_u8((const chtype*)(srcrow + tl_ofs[7])), mask_rgb); \ + \ + uint8x16_t _z04a = vzip1q_u8(_r0, _r4); \ + uint8x16_t _z15a = vzip1q_u8(_r1, _r5); \ + uint8x16_t _z26a = vzip1q_u8(_r2, _r6); \ + uint8x16_t _z37a = vzip1q_u8(_r3, _r7); \ + uint8x16_t _z04b = vzip2q_u8(_r0, _r4); \ + uint8x16_t _z15b = vzip2q_u8(_r1, _r5); \ + uint8x16_t _z26b = vzip2q_u8(_r2, _r6); \ + uint8x16_t _z37b = vzip2q_u8(_r3, _r7); \ + \ + uint8x16_t _r0246 = vzip1q_u8(_z04a, _z26a); \ + uint8x16_t _g0246 = vzip2q_u8(_z04a, _z26a); \ + uint8x16_t _r1357 = vzip1q_u8(_z15a, _z37a); \ + uint8x16_t _g1357 = vzip2q_u8(_z15a, _z37a); \ + uint8x16_t _b0246 = vzip1q_u8(_z04b, _z26b); \ + uint8x16_t _a0246 = vzip2q_u8(_z04b, _z26b); \ + uint8x16_t _b1357 = vzip1q_u8(_z15b, _z37b); \ + uint8x16_t _a1357 = vzip2q_u8(_z15b, _z37b); \ + \ + uint8x16_t _r_c01 = vzip1q_u8(_r0246, _r1357); \ + uint8x16_t _r_c23 = vzip2q_u8(_r0246, _r1357); \ + uint8x16_t _g_c01 = vzip1q_u8(_g0246, _g1357); \ + uint8x16_t _g_c23 = vzip2q_u8(_g0246, _g1357); \ + uint8x16_t _b_c01 = vzip1q_u8(_b0246, _b1357); \ + uint8x16_t _b_c23 = vzip2q_u8(_b0246, _b1357); \ + uint8x16_t _a_c01 = vzip1q_u8(_a0246, _a1357); \ + uint8x16_t _a_c23 = vzip2q_u8(_a0246, _a1357); \ + \ + R0.val = vreinterpretq_s16_u16(vmovl_u8(vget_low_u8(_r_c01))); \ + R1.val = vreinterpretq_s16_u16(vmovl_high_u8(_r_c01)); \ + R2.val = vreinterpretq_s16_u16(vmovl_u8(vget_low_u8(_r_c23))); \ + R3.val = vreinterpretq_s16_u16(vmovl_high_u8(_r_c23)); \ + G0.val = vreinterpretq_s16_u16(vmovl_u8(vget_low_u8(_g_c01))); \ + G1.val = vreinterpretq_s16_u16(vmovl_high_u8(_g_c01)); \ + G2.val = vreinterpretq_s16_u16(vmovl_u8(vget_low_u8(_g_c23))); \ + G3.val = vreinterpretq_s16_u16(vmovl_high_u8(_g_c23)); \ + B0.val = vreinterpretq_s16_u16(vmovl_u8(vget_low_u8(_b_c01))); \ + B1.val = vreinterpretq_s16_u16(vmovl_high_u8(_b_c01)); \ + B2.val = vreinterpretq_s16_u16(vmovl_u8(vget_low_u8(_b_c23))); \ + B3.val = vreinterpretq_s16_u16(vmovl_high_u8(_b_c23)); \ + A0.val = vreinterpretq_s16_u16(vmovl_u8(vget_low_u8(_a_c01))); \ + A1.val = vreinterpretq_s16_u16(vmovl_high_u8(_a_c01)); \ + A2.val = vreinterpretq_s16_u16(vmovl_u8(vget_low_u8(_a_c23))); \ + A3.val = vreinterpretq_s16_u16(vmovl_high_u8(_a_c23)); } + +#else +#undef FETCH_INLIERS_8UC1 +#define FETCH_INLIERS_8UC1(row) { \ + const uint8_t* srcrow = (const uint8_t*)src + row*srcstep; \ + uint8x8_t _r0 = vld1_u8(srcrow + tl_ofs[0]); \ + uint8x8_t _r1 = vld1_u8(srcrow + tl_ofs[1]); \ + uint8x8_t _r2 = vld1_u8(srcrow + tl_ofs[2]); \ + uint8x8_t _r3 = vld1_u8(srcrow + tl_ofs[3]); \ + uint16x4_t _rw0 = vget_low_u16(vmovl_u8(_r0)); \ + uint16x4_t _rw1 = vget_low_u16(vmovl_u8(_r1)); \ + uint16x4_t _rw2 = vget_low_u16(vmovl_u8(_r2)); \ + uint16x4_t _rw3 = vget_low_u16(vmovl_u8(_r3)); \ + \ + uint16x4_t _r02a = vzip1_u16(_rw0, _rw2); \ + uint16x4_t _r02b = vzip2_u16(_rw0, _rw2); \ + uint16x4_t _r13a = vzip1_u16(_rw1, _rw3); \ + uint16x4_t _r13b = vzip2_u16(_rw1, _rw3); \ + R0.val = vreinterpretq_s32_u32(vmovl_u16(vzip1_u16(_r02a, _r13a))); \ + R1.val = vreinterpretq_s32_u32(vmovl_u16(vzip2_u16(_r02a, _r13a))); \ + R2.val = vreinterpretq_s32_u32(vmovl_u16(vzip1_u16(_r02b, _r13b))); \ + R3.val = vreinterpretq_s32_u32(vmovl_u16(vzip2_u16(_r02b, _r13b))); } + +#undef FETCH_INLIERS_8UC3 +#define FETCH_INLIERS_8UC3(row) { \ + const uint8x16_t mask_rgb = { 0, 3, 6, 9, 1, 4, 7, 10, 2, 5, 8, 11, 255, 255, 255, 255 }; \ + const uint8_t* srcrow = (const uint8_t*)src + row*srcstep; \ + uint8x16_t _r0 = vqtbl1q_u8(vld1q_u8(srcrow + tl_ofs[0]), mask_rgb); \ + uint8x16_t _r1 = vqtbl1q_u8(vld1q_u8(srcrow + tl_ofs[1]), mask_rgb); \ + uint8x16_t _r2 = vqtbl1q_u8(vld1q_u8(srcrow + tl_ofs[2]), mask_rgb); \ + uint8x16_t _r3 = vqtbl1q_u8(vld1q_u8(srcrow + tl_ofs[3]), mask_rgb); \ + \ + uint8x16_t _z02a = vzip1q_u8(_r0, _r2); \ + uint8x16_t _z02b = vzip2q_u8(_r0, _r2); \ + uint8x16_t _z13a = vzip1q_u8(_r1, _r3); \ + uint8x16_t _z13b = vzip2q_u8(_r1, _r3); \ + \ + uint8x16_t _r_c0123 = vzip1q_u8(_z02a, _z13a); \ + uint8x16_t _g_c0123 = vzip2q_u8(_z02a, _z13a); \ + uint8x16_t _b_c0123 = vzip1q_u8(_z02b, _z13b); \ + \ + uint16x8_t _rl = vmovl_u8(vget_low_u8(_r_c0123)); \ + uint16x8_t _rh = vmovl_high_u8(_r_c0123); \ + uint16x8_t _gl = vmovl_u8(vget_low_u8(_g_c0123)); \ + uint16x8_t _gh = vmovl_high_u8(_g_c0123); \ + uint16x8_t _bl = vmovl_u8(vget_low_u8(_b_c0123)); \ + uint16x8_t _bh = vmovl_high_u8(_b_c0123); \ + \ + R0.val = vreinterpretq_s32_u32(vmovl_u16(vget_low_u16(_rl))); \ + R1.val = vreinterpretq_s32_u32(vmovl_high_u16(_rl)); \ + R2.val = vreinterpretq_s32_u32(vmovl_u16(vget_low_u16(_rh))); \ + R3.val = vreinterpretq_s32_u32(vmovl_high_u16(_rh)); \ + G0.val = vreinterpretq_s32_u32(vmovl_u16(vget_low_u16(_gl))); \ + G1.val = vreinterpretq_s32_u32(vmovl_high_u16(_gl)); \ + G2.val = vreinterpretq_s32_u32(vmovl_u16(vget_low_u16(_gh))); \ + G3.val = vreinterpretq_s32_u32(vmovl_high_u16(_gh)); \ + B0.val = vreinterpretq_s32_u32(vmovl_u16(vget_low_u16(_bl))); \ + B1.val = vreinterpretq_s32_u32(vmovl_high_u16(_bl)); \ + B2.val = vreinterpretq_s32_u32(vmovl_u16(vget_low_u16(_bh))); \ + B3.val = vreinterpretq_s32_u32(vmovl_high_u16(_bh)); } + +#undef FETCH_INLIERS_8UC4 +#define FETCH_INLIERS_8UC4(row) { \ + const uint8x16_t mask_rgb = { 0, 4, 8, 12, 1, 5, 9, 13, 2, 6, 10, 14, 3, 7, 11, 15 }; \ + const uint8_t* srcrow = (const uint8_t*)src + row*srcstep; \ + uint8x16_t _r0 = vqtbl1q_u8(vld1q_u8(srcrow + tl_ofs[0]), mask_rgb); \ + uint8x16_t _r1 = vqtbl1q_u8(vld1q_u8(srcrow + tl_ofs[1]), mask_rgb); \ + uint8x16_t _r2 = vqtbl1q_u8(vld1q_u8(srcrow + tl_ofs[2]), mask_rgb); \ + uint8x16_t _r3 = vqtbl1q_u8(vld1q_u8(srcrow + tl_ofs[3]), mask_rgb); \ + \ + uint8x16_t _z02a = vzip1q_u8(_r0, _r2); \ + uint8x16_t _z02b = vzip2q_u8(_r0, _r2); \ + uint8x16_t _z13a = vzip1q_u8(_r1, _r3); \ + uint8x16_t _z13b = vzip2q_u8(_r1, _r3); \ + \ + uint8x16_t _r_c0123 = vzip1q_u8(_z02a, _z13a); \ + uint8x16_t _g_c0123 = vzip2q_u8(_z02a, _z13a); \ + uint8x16_t _b_c0123 = vzip1q_u8(_z02b, _z13b); \ + uint8x16_t _a_c0123 = vzip2q_u8(_z02b, _z13b); \ + \ + uint16x8_t _rl = vmovl_u8(vget_low_u8(_r_c0123)); \ + uint16x8_t _rh = vmovl_high_u8(_r_c0123); \ + uint16x8_t _gl = vmovl_u8(vget_low_u8(_g_c0123)); \ + uint16x8_t _gh = vmovl_high_u8(_g_c0123); \ + uint16x8_t _bl = vmovl_u8(vget_low_u8(_b_c0123)); \ + uint16x8_t _bh = vmovl_high_u8(_b_c0123); \ + uint16x8_t _al = vmovl_u8(vget_low_u8(_a_c0123)); \ + uint16x8_t _ah = vmovl_high_u8(_a_c0123); \ + \ + R0.val = vreinterpretq_s32_u32(vmovl_u16(vget_low_u16(_rl))); \ + R1.val = vreinterpretq_s32_u32(vmovl_high_u16(_rl)); \ + R2.val = vreinterpretq_s32_u32(vmovl_u16(vget_low_u16(_rh))); \ + R3.val = vreinterpretq_s32_u32(vmovl_high_u16(_rh)); \ + G0.val = vreinterpretq_s32_u32(vmovl_u16(vget_low_u16(_gl))); \ + G1.val = vreinterpretq_s32_u32(vmovl_high_u16(_gl)); \ + G2.val = vreinterpretq_s32_u32(vmovl_u16(vget_low_u16(_gh))); \ + G3.val = vreinterpretq_s32_u32(vmovl_high_u16(_gh)); \ + B0.val = vreinterpretq_s32_u32(vmovl_u16(vget_low_u16(_bl))); \ + B1.val = vreinterpretq_s32_u32(vmovl_high_u16(_bl)); \ + B2.val = vreinterpretq_s32_u32(vmovl_u16(vget_low_u16(_bh))); \ + B3.val = vreinterpretq_s32_u32(vmovl_high_u16(_bh)); \ + A0.val = vreinterpretq_s32_u32(vmovl_u16(vget_low_u16(_al))); \ + A1.val = vreinterpretq_s32_u32(vmovl_high_u16(_al)); \ + A2.val = vreinterpretq_s32_u32(vmovl_u16(vget_low_u16(_ah))); \ + A3.val = vreinterpretq_s32_u32(vmovl_high_u16(_ah)); } +#endif + +#undef FETCH_INLIERS_16UC1 +#define FETCH_INLIERS_16UC1(row) { \ + const uint8_t* srcrow = (const uint8_t*)src + row*srcstep; \ + uint16x4_t _r0 = vld1_u16((const chtype*)(srcrow + tl_ofs[0])); \ + uint16x4_t _r1 = vld1_u16((const chtype*)(srcrow + tl_ofs[1])); \ + uint16x4_t _r2 = vld1_u16((const chtype*)(srcrow + tl_ofs[2])); \ + uint16x4_t _r3 = vld1_u16((const chtype*)(srcrow + tl_ofs[3])); \ + \ + uint16x4_t _r02a = vzip1_u16(_r0, _r2); \ + uint16x4_t _r02b = vzip2_u16(_r0, _r2); \ + uint16x4_t _r13a = vzip1_u16(_r1, _r3); \ + uint16x4_t _r13b = vzip2_u16(_r1, _r3); \ + R0.val = vreinterpretq_s32_u32(vmovl_u16(vzip1_u16(_r02a, _r13a))); \ + R1.val = vreinterpretq_s32_u32(vmovl_u16(vzip2_u16(_r02a, _r13a))); \ + R2.val = vreinterpretq_s32_u32(vmovl_u16(vzip1_u16(_r02b, _r13b))); \ + R3.val = vreinterpretq_s32_u32(vmovl_u16(vzip2_u16(_r02b, _r13b))); } + +#undef FETCH_INLIERS_16UC3 +#define FETCH_INLIERS_16UC3(row) { \ + const uint8_t* srcrow = (const uint8_t*)src + row*srcstep; \ + uint16x4x3_t _r0 = vld3_u16((const chtype*)(srcrow + tl_ofs[0])); \ + uint16x4x3_t _r1 = vld3_u16((const chtype*)(srcrow + tl_ofs[1])); \ + uint16x4x3_t _r2 = vld3_u16((const chtype*)(srcrow + tl_ofs[2])); \ + uint16x4x3_t _r3 = vld3_u16((const chtype*)(srcrow + tl_ofs[3])); \ + uint16x4_t _r02a = vzip1_u16(_r0.val[0], _r2.val[0]); \ + uint16x4_t _r02b = vzip2_u16(_r0.val[0], _r2.val[0]); \ + uint16x4_t _r13a = vzip1_u16(_r1.val[0], _r3.val[0]); \ + uint16x4_t _r13b = vzip2_u16(_r1.val[0], _r3.val[0]); \ + uint16x4_t _g02a = vzip1_u16(_r0.val[1], _r2.val[1]); \ + uint16x4_t _g02b = vzip2_u16(_r0.val[1], _r2.val[1]); \ + uint16x4_t _g13a = vzip1_u16(_r1.val[1], _r3.val[1]); \ + uint16x4_t _g13b = vzip2_u16(_r1.val[1], _r3.val[1]); \ + uint16x4_t _b02a = vzip1_u16(_r0.val[2], _r2.val[2]); \ + uint16x4_t _b02b = vzip2_u16(_r0.val[2], _r2.val[2]); \ + uint16x4_t _b13a = vzip1_u16(_r1.val[2], _r3.val[2]); \ + uint16x4_t _b13b = vzip2_u16(_r1.val[2], _r3.val[2]); \ + R0.val = vreinterpretq_s32_u32(vmovl_u16(vzip1_u16(_r02a, _r13a))); \ + R1.val = vreinterpretq_s32_u32(vmovl_u16(vzip2_u16(_r02a, _r13a))); \ + R2.val = vreinterpretq_s32_u32(vmovl_u16(vzip1_u16(_r02b, _r13b))); \ + R3.val = vreinterpretq_s32_u32(vmovl_u16(vzip2_u16(_r02b, _r13b))); \ + G0.val = vreinterpretq_s32_u32(vmovl_u16(vzip1_u16(_g02a, _g13a))); \ + G1.val = vreinterpretq_s32_u32(vmovl_u16(vzip2_u16(_g02a, _g13a))); \ + G2.val = vreinterpretq_s32_u32(vmovl_u16(vzip1_u16(_g02b, _g13b))); \ + G3.val = vreinterpretq_s32_u32(vmovl_u16(vzip2_u16(_g02b, _g13b))); \ + B0.val = vreinterpretq_s32_u32(vmovl_u16(vzip1_u16(_b02a, _b13a))); \ + B1.val = vreinterpretq_s32_u32(vmovl_u16(vzip2_u16(_b02a, _b13a))); \ + B2.val = vreinterpretq_s32_u32(vmovl_u16(vzip1_u16(_b02b, _b13b))); \ + B3.val = vreinterpretq_s32_u32(vmovl_u16(vzip2_u16(_b02b, _b13b))); } + +#undef FETCH_INLIERS_16UC4 +#define FETCH_INLIERS_16UC4(row) { \ + const uint8_t* srcrow = (const uint8_t*)src + row*srcstep; \ + uint16x4x4_t _r0 = vld4_u16((const chtype*)(srcrow + tl_ofs[0])); \ + uint16x4x4_t _r1 = vld4_u16((const chtype*)(srcrow + tl_ofs[1])); \ + uint16x4x4_t _r2 = vld4_u16((const chtype*)(srcrow + tl_ofs[2])); \ + uint16x4x4_t _r3 = vld4_u16((const chtype*)(srcrow + tl_ofs[3])); \ + uint16x4_t _r02a = vzip1_u16(_r0.val[0], _r2.val[0]); \ + uint16x4_t _r02b = vzip2_u16(_r0.val[0], _r2.val[0]); \ + uint16x4_t _r13a = vzip1_u16(_r1.val[0], _r3.val[0]); \ + uint16x4_t _r13b = vzip2_u16(_r1.val[0], _r3.val[0]); \ + uint16x4_t _g02a = vzip1_u16(_r0.val[1], _r2.val[1]); \ + uint16x4_t _g02b = vzip2_u16(_r0.val[1], _r2.val[1]); \ + uint16x4_t _g13a = vzip1_u16(_r1.val[1], _r3.val[1]); \ + uint16x4_t _g13b = vzip2_u16(_r1.val[1], _r3.val[1]); \ + uint16x4_t _b02a = vzip1_u16(_r0.val[2], _r2.val[2]); \ + uint16x4_t _b02b = vzip2_u16(_r0.val[2], _r2.val[2]); \ + uint16x4_t _b13a = vzip1_u16(_r1.val[2], _r3.val[2]); \ + uint16x4_t _b13b = vzip2_u16(_r1.val[2], _r3.val[2]); \ + uint16x4_t _a02a = vzip1_u16(_r0.val[3], _r2.val[3]); \ + uint16x4_t _a02b = vzip2_u16(_r0.val[3], _r2.val[3]); \ + uint16x4_t _a13a = vzip1_u16(_r1.val[3], _r3.val[3]); \ + uint16x4_t _a13b = vzip2_u16(_r1.val[3], _r3.val[3]); \ + R0.val = vreinterpretq_s32_u32(vmovl_u16(vzip1_u16(_r02a, _r13a))); \ + R1.val = vreinterpretq_s32_u32(vmovl_u16(vzip2_u16(_r02a, _r13a))); \ + R2.val = vreinterpretq_s32_u32(vmovl_u16(vzip1_u16(_r02b, _r13b))); \ + R3.val = vreinterpretq_s32_u32(vmovl_u16(vzip2_u16(_r02b, _r13b))); \ + G0.val = vreinterpretq_s32_u32(vmovl_u16(vzip1_u16(_g02a, _g13a))); \ + G1.val = vreinterpretq_s32_u32(vmovl_u16(vzip2_u16(_g02a, _g13a))); \ + G2.val = vreinterpretq_s32_u32(vmovl_u16(vzip1_u16(_g02b, _g13b))); \ + G3.val = vreinterpretq_s32_u32(vmovl_u16(vzip2_u16(_g02b, _g13b))); \ + B0.val = vreinterpretq_s32_u32(vmovl_u16(vzip1_u16(_b02a, _b13a))); \ + B1.val = vreinterpretq_s32_u32(vmovl_u16(vzip2_u16(_b02a, _b13a))); \ + B2.val = vreinterpretq_s32_u32(vmovl_u16(vzip1_u16(_b02b, _b13b))); \ + B3.val = vreinterpretq_s32_u32(vmovl_u16(vzip2_u16(_b02b, _b13b))); \ + A0.val = vreinterpretq_s32_u32(vmovl_u16(vzip1_u16(_a02a, _a13a))); \ + A1.val = vreinterpretq_s32_u32(vmovl_u16(vzip2_u16(_a02a, _a13a))); \ + A2.val = vreinterpretq_s32_u32(vmovl_u16(vzip1_u16(_a02b, _a13b))); \ + A3.val = vreinterpretq_s32_u32(vmovl_u16(vzip2_u16(_a02b, _a13b))); } + +#undef FETCH_INLIERS_32FC1 +#define FETCH_INLIERS_32FC1(row) { \ + const uint8_t* srcrow = (const uint8_t*)src + row*srcstep; \ + float32x4_t _r0 = vld1q_f32((const chtype*)(srcrow + tl_ofs[0])); \ + float32x4_t _r1 = vld1q_f32((const chtype*)(srcrow + tl_ofs[1])); \ + float32x4_t _r2 = vld1q_f32((const chtype*)(srcrow + tl_ofs[2])); \ + float32x4_t _r3 = vld1q_f32((const chtype*)(srcrow + tl_ofs[3])); \ + \ + float32x4_t _r02a = vzip1q_f32(_r0, _r2); \ + float32x4_t _r02b = vzip2q_f32(_r0, _r2); \ + float32x4_t _r13a = vzip1q_f32(_r1, _r3); \ + float32x4_t _r13b = vzip2q_f32(_r1, _r3); \ + R0.val = vzip1q_f32(_r02a, _r13a); \ + R1.val = vzip2q_f32(_r02a, _r13a); \ + R2.val = vzip1q_f32(_r02b, _r13b); \ + R3.val = vzip2q_f32(_r02b, _r13b); } + +#undef FETCH_INLIERS_32FC3 +#define FETCH_INLIERS_32FC3(row) { \ + const uint8_t* srcrow = (const uint8_t*)src + row*srcstep; \ + float32x4x3_t _r0 = vld3q_f32((const chtype*)(srcrow + tl_ofs[0])); \ + float32x4x3_t _r1 = vld3q_f32((const chtype*)(srcrow + tl_ofs[1])); \ + float32x4x3_t _r2 = vld3q_f32((const chtype*)(srcrow + tl_ofs[2])); \ + float32x4x3_t _r3 = vld3q_f32((const chtype*)(srcrow + tl_ofs[3])); \ + float32x4_t _r02a = vzip1q_f32(_r0.val[0], _r2.val[0]); \ + float32x4_t _r02b = vzip2q_f32(_r0.val[0], _r2.val[0]); \ + float32x4_t _r13a = vzip1q_f32(_r1.val[0], _r3.val[0]); \ + float32x4_t _r13b = vzip2q_f32(_r1.val[0], _r3.val[0]); \ + float32x4_t _g02a = vzip1q_f32(_r0.val[1], _r2.val[1]); \ + float32x4_t _g02b = vzip2q_f32(_r0.val[1], _r2.val[1]); \ + float32x4_t _g13a = vzip1q_f32(_r1.val[1], _r3.val[1]); \ + float32x4_t _g13b = vzip2q_f32(_r1.val[1], _r3.val[1]); \ + float32x4_t _b02a = vzip1q_f32(_r0.val[2], _r2.val[2]); \ + float32x4_t _b02b = vzip2q_f32(_r0.val[2], _r2.val[2]); \ + float32x4_t _b13a = vzip1q_f32(_r1.val[2], _r3.val[2]); \ + float32x4_t _b13b = vzip2q_f32(_r1.val[2], _r3.val[2]); \ + R0.val = vzip1q_f32(_r02a, _r13a); \ + R1.val = vzip2q_f32(_r02a, _r13a); \ + R2.val = vzip1q_f32(_r02b, _r13b); \ + R3.val = vzip2q_f32(_r02b, _r13b); \ + G0.val = vzip1q_f32(_g02a, _g13a); \ + G1.val = vzip2q_f32(_g02a, _g13a); \ + G2.val = vzip1q_f32(_g02b, _g13b); \ + G3.val = vzip2q_f32(_g02b, _g13b); \ + B0.val = vzip1q_f32(_b02a, _b13a); \ + B1.val = vzip2q_f32(_b02a, _b13a); \ + B2.val = vzip1q_f32(_b02b, _b13b); \ + B3.val = vzip2q_f32(_b02b, _b13b); } + +#undef FETCH_INLIERS_32FC4 +#define FETCH_INLIERS_32FC4(row) { \ + const uint8_t* srcrow = (const uint8_t*)src + row*srcstep; \ + float32x4x4_t _r0 = vld4q_f32((const chtype*)(srcrow + tl_ofs[0])); \ + float32x4x4_t _r1 = vld4q_f32((const chtype*)(srcrow + tl_ofs[1])); \ + float32x4x4_t _r2 = vld4q_f32((const chtype*)(srcrow + tl_ofs[2])); \ + float32x4x4_t _r3 = vld4q_f32((const chtype*)(srcrow + tl_ofs[3])); \ + float32x4_t _r02a = vzip1q_f32(_r0.val[0], _r2.val[0]); \ + float32x4_t _r02b = vzip2q_f32(_r0.val[0], _r2.val[0]); \ + float32x4_t _r13a = vzip1q_f32(_r1.val[0], _r3.val[0]); \ + float32x4_t _r13b = vzip2q_f32(_r1.val[0], _r3.val[0]); \ + float32x4_t _g02a = vzip1q_f32(_r0.val[1], _r2.val[1]); \ + float32x4_t _g02b = vzip2q_f32(_r0.val[1], _r2.val[1]); \ + float32x4_t _g13a = vzip1q_f32(_r1.val[1], _r3.val[1]); \ + float32x4_t _g13b = vzip2q_f32(_r1.val[1], _r3.val[1]); \ + float32x4_t _b02a = vzip1q_f32(_r0.val[2], _r2.val[2]); \ + float32x4_t _b02b = vzip2q_f32(_r0.val[2], _r2.val[2]); \ + float32x4_t _b13a = vzip1q_f32(_r1.val[2], _r3.val[2]); \ + float32x4_t _b13b = vzip2q_f32(_r1.val[2], _r3.val[2]); \ + float32x4_t _a02a = vzip1q_f32(_r0.val[3], _r2.val[3]); \ + float32x4_t _a02b = vzip2q_f32(_r0.val[3], _r2.val[3]); \ + float32x4_t _a13a = vzip1q_f32(_r1.val[3], _r3.val[3]); \ + float32x4_t _a13b = vzip2q_f32(_r1.val[3], _r3.val[3]); \ + R0.val = vzip1q_f32(_r02a, _r13a); \ + R1.val = vzip2q_f32(_r02a, _r13a); \ + R2.val = vzip1q_f32(_r02b, _r13b); \ + R3.val = vzip2q_f32(_r02b, _r13b); \ + G0.val = vzip1q_f32(_g02a, _g13a); \ + G1.val = vzip2q_f32(_g02a, _g13a); \ + G2.val = vzip1q_f32(_g02b, _g13b); \ + G3.val = vzip2q_f32(_g02b, _g13b); \ + B0.val = vzip1q_f32(_b02a, _b13a); \ + B1.val = vzip2q_f32(_b02a, _b13a); \ + B2.val = vzip1q_f32(_b02b, _b13b); \ + B3.val = vzip2q_f32(_b02b, _b13b); \ + A0.val = vzip1q_f32(_a02a, _a13a); \ + A1.val = vzip2q_f32(_a02a, _a13a); \ + A2.val = vzip1q_f32(_a02b, _a13b); \ + A3.val = vzip2q_f32(_a02b, _a13b); } + +#endif + +#undef BICUBIC_C1_PROCESS_ROW +#define BICUBIC_C1_PROCESS_ROW(row, fetch_inliers) \ + fetch_inliers(row); \ + BICUBIC_UPDATE_ACC_VEC(acc_r, R0, R1, R2, R3, wy##row) + +#undef BICUBIC_C3_PROCESS_ROW +#define BICUBIC_C3_PROCESS_ROW(row, fetch_inliers) \ + fetch_inliers(row); \ + BICUBIC_UPDATE_ACC_VEC(acc_r, R0, R1, R2, R3, wy##row); \ + BICUBIC_UPDATE_ACC_VEC(acc_g, G0, G1, G2, G3, wy##row); \ + BICUBIC_UPDATE_ACC_VEC(acc_b, B0, B1, B2, B3, wy##row) + +#undef BICUBIC_C4_PROCESS_ROW +#define BICUBIC_C4_PROCESS_ROW(row, fetch_inliers) \ + fetch_inliers(row); \ + BICUBIC_UPDATE_ACC_VEC(acc_r, R0, R1, R2, R3, wy##row); \ + BICUBIC_UPDATE_ACC_VEC(acc_g, G0, G1, G2, G3, wy##row); \ + BICUBIC_UPDATE_ACC_VEC(acc_b, B0, B1, B2, B3, wy##row); \ + BICUBIC_UPDATE_ACC_VEC(acc_a, A0, A1, A2, A3, wy##row) + +template +static void +bicubicVec(const float* srcx, const float* srcy, int len, + const void* src, size_t srcstep, Size size, + chtype* dst, const float* params, + int borderType, chtype* borderVal) +{ + constexpr float defaultA = -0.75f; + float A = params ? *params : defaultA; + + using pixtype = std::conditional_t>; + using buftype = typename VTraits::lane_type; + + constexpr int BPP = int(NCHANNELS*sizeof(dst[0])); + constexpr int BATCH = VTraits::nlanes; + float xbuf[BATCH] = {}, ybuf[BATCH] = {}; + pixtype bval = borderVal ? *(pixtype*)borderVal : pixtype(); + + int savestorelen = borderType == BORDER_TRANSPARENT ? 0 : len; + + vecfptype acc_r = v_setzero_(), acc_g = v_setzero_(), + acc_b = v_setzero_(), acc_a = v_setzero_(); + + for (int i = 0; i < len; i += BATCH, dst += BATCH*NCHANNELS) { + if (i + BATCH > len && i > 0 && borderType != BORDER_TRANSPARENT) { + int i1 = len - BATCH; + dst -= (i - i1)*NCHANNELS; + i = i1; + } + int dlen = std::min(BATCH, len - i); + int32_t tl_x[BATCH], tl_y[BATCH], tl_ofs[BATCH], goodx[BATCH*4]; + buftype pixbuf[NCHANNELS][BATCH*4]; + vecfptype wx0, wx1, wx2, wx3, wy0, wy1, wy2, wy3, sumwx; + const float *xptr = srcx + i, *yptr = srcy + i; + if (dlen < BATCH) { + memcpy(xbuf, srcx + i, dlen*sizeof(srcx[0])); + memcpy(ybuf, srcy + i, dlen*sizeof(srcy[0])); + xptr = xbuf; + yptr = ybuf; + } + int code = bicubicCoeffs(xptr, yptr, srcstep, size, BPP, A, + borderType, tl_x, tl_y, tl_ofs, + wx0, wx1, wx2, wx3, wy0, wy1, wy2, wy3); + + if (code > 0) { + if constexpr (NCHANNELS == 1) { + vecbuftype R0, R1, R2, R3; + acc_r = v_setzero_(); + if constexpr (std::is_same_v) { + BICUBIC_C1_PROCESS_ROW(0, FETCH_INLIERS_8UC1); + BICUBIC_C1_PROCESS_ROW(1, FETCH_INLIERS_8UC1); + BICUBIC_C1_PROCESS_ROW(2, FETCH_INLIERS_8UC1); + BICUBIC_C1_PROCESS_ROW(3, FETCH_INLIERS_8UC1); + } else if constexpr (std::is_same_v) { + BICUBIC_C1_PROCESS_ROW(0, FETCH_INLIERS_16UC1); + BICUBIC_C1_PROCESS_ROW(1, FETCH_INLIERS_16UC1); + BICUBIC_C1_PROCESS_ROW(2, FETCH_INLIERS_16UC1); + BICUBIC_C1_PROCESS_ROW(3, FETCH_INLIERS_16UC1); + } else if constexpr (std::is_same_v) { + BICUBIC_C1_PROCESS_ROW(0, FETCH_INLIERS_32FC1); + BICUBIC_C1_PROCESS_ROW(1, FETCH_INLIERS_32FC1); + BICUBIC_C1_PROCESS_ROW(2, FETCH_INLIERS_32FC1); + BICUBIC_C1_PROCESS_ROW(3, FETCH_INLIERS_32FC1); + } + } else if constexpr (NCHANNELS == 3) { + vecbuftype R0, R1, R2, R3, G0, G1, G2, G3, B0, B1, B2, B3; + acc_r = acc_g = acc_b = v_setzero_(); + if constexpr (std::is_same_v) { + BICUBIC_C3_PROCESS_ROW(0, FETCH_INLIERS_8UC3); + BICUBIC_C3_PROCESS_ROW(1, FETCH_INLIERS_8UC3); + BICUBIC_C3_PROCESS_ROW(2, FETCH_INLIERS_8UC3); + BICUBIC_C3_PROCESS_ROW(3, FETCH_INLIERS_8UC3); + } else if constexpr (std::is_same_v) { + BICUBIC_C3_PROCESS_ROW(0, FETCH_INLIERS_16UC3); + BICUBIC_C3_PROCESS_ROW(1, FETCH_INLIERS_16UC3); + BICUBIC_C3_PROCESS_ROW(2, FETCH_INLIERS_16UC3); + BICUBIC_C3_PROCESS_ROW(3, FETCH_INLIERS_16UC3); + } else if constexpr (std::is_same_v) { + BICUBIC_C3_PROCESS_ROW(0, FETCH_INLIERS_32FC3); + BICUBIC_C3_PROCESS_ROW(1, FETCH_INLIERS_32FC3); + BICUBIC_C3_PROCESS_ROW(2, FETCH_INLIERS_32FC3); + BICUBIC_C3_PROCESS_ROW(3, FETCH_INLIERS_32FC3); + } + } else if constexpr (NCHANNELS == 4) { + vecbuftype R0, R1, R2, R3, G0, G1, G2, G3, B0, B1, B2, B3, A0, A1, A2, A3; + acc_r = acc_g = acc_b = acc_a = v_setzero_(); + if constexpr (std::is_same_v) { + BICUBIC_C4_PROCESS_ROW(0, FETCH_INLIERS_8UC4); + BICUBIC_C4_PROCESS_ROW(1, FETCH_INLIERS_8UC4); + BICUBIC_C4_PROCESS_ROW(2, FETCH_INLIERS_8UC4); + BICUBIC_C4_PROCESS_ROW(3, FETCH_INLIERS_8UC4); + } else if constexpr (std::is_same_v) { + BICUBIC_C4_PROCESS_ROW(0, FETCH_INLIERS_16UC4); + BICUBIC_C4_PROCESS_ROW(1, FETCH_INLIERS_16UC4); + BICUBIC_C4_PROCESS_ROW(2, FETCH_INLIERS_16UC4); + BICUBIC_C4_PROCESS_ROW(3, FETCH_INLIERS_16UC4); + } else if constexpr (std::is_same_v) { + BICUBIC_C4_PROCESS_ROW(0, FETCH_INLIERS_32FC4); + BICUBIC_C4_PROCESS_ROW(1, FETCH_INLIERS_32FC4); + BICUBIC_C4_PROCESS_ROW(2, FETCH_INLIERS_32FC4); + BICUBIC_C4_PROCESS_ROW(3, FETCH_INLIERS_32FC4); + } + } + } else if (code < 0) { + if (borderType == BORDER_CONSTANT) { + for (int j = 0; j < dlen; j++) { + ((pixtype*)dst)[j] = bval; + } + } + continue; + } else { + const chtype* defVal = borderType == BORDER_TRANSPARENT ? dst : borderVal; + vecfptype wys[4] = {wy0, wy1, wy2, wy3}; + vecbuftype V0, V1, V2, V3; + acc_r = acc_g = acc_b = acc_a = v_setzero_(); + + for (int row = 0; row < 4; row++) { + bicubicFetchPixels((const chtype*)src, srcstep, size, NCHANNELS, tl_x, tl_y, + goodx, row, &pixbuf[0][0], BATCH, borderType, defVal); + vecfptype wy = wys[row]; + V0 = vx_load(&pixbuf[0][0]); + V1 = vx_load(&pixbuf[0][BATCH]); + V2 = vx_load(&pixbuf[0][BATCH*2]); + V3 = vx_load(&pixbuf[0][BATCH*3]); + BICUBIC_UPDATE_ACC_VEC(acc_r, V0, V1, V2, V3, wy); + if constexpr (NCHANNELS > 1) { + V0 = vx_load(&pixbuf[1][0]); + V1 = vx_load(&pixbuf[1][BATCH]); + V2 = vx_load(&pixbuf[1][BATCH*2]); + V3 = vx_load(&pixbuf[1][BATCH*3]); + BICUBIC_UPDATE_ACC_VEC(acc_g, V0, V1, V2, V3, wy); + V0 = vx_load(&pixbuf[2][0]); + V1 = vx_load(&pixbuf[2][BATCH]); + V2 = vx_load(&pixbuf[2][BATCH*2]); + V3 = vx_load(&pixbuf[2][BATCH*3]); + BICUBIC_UPDATE_ACC_VEC(acc_b, V0, V1, V2, V3, wy); + if constexpr (NCHANNELS > 3) { + V0 = vx_load(&pixbuf[3][0]); + V1 = vx_load(&pixbuf[3][BATCH]); + V2 = vx_load(&pixbuf[3][BATCH*2]); + V3 = vx_load(&pixbuf[3][BATCH*3]); + BICUBIC_UPDATE_ACC_VEC(acc_a, V0, V1, V2, V3, wy); + } + } + } + } + + // store the result + chtype outbuf[BATCH*NCHANNELS*4]; + if constexpr (NCHANNELS == 1) { + if constexpr (std::is_same_v) { + #if CV_SIMD_FP16 + v_int16 wacc_r = v_round(acc_r); + #else + v_int32 iacc_r = v_round(acc_r); + v_int16 wacc_r = v_pack(iacc_r, iacc_r); + #endif + v_uint8 bacc_r = v_pack_u(wacc_r, wacc_r); + if (i + BATCH*4 <= savestorelen) { + v_store(dst, bacc_r); + continue; + } + v_store(outbuf, bacc_r); + } else if constexpr (std::is_same_v) { + v_int32 iacc_r = v_round(acc_r); + v_uint16 wacc_r = v_pack_u(iacc_r, iacc_r); + if (i + BATCH*2 <= savestorelen) { + v_store(dst, wacc_r); + continue; + } + v_store(outbuf, wacc_r); + } else if constexpr (std::is_same_v) { + if (i + BATCH <= savestorelen) { + v_store(dst, acc_r); + continue; + } + v_store(outbuf, acc_r); + } + } else if constexpr (NCHANNELS == 3) { + if constexpr (std::is_same_v) { + #if CV_SIMD_FP16 + v_int16 wacc_r = v_round(acc_r); + v_int16 wacc_g = v_round(acc_g); + v_int16 wacc_b = v_round(acc_b); + #else + v_int32 iacc_r = v_round(acc_r); + v_int32 iacc_g = v_round(acc_g); + v_int32 iacc_b = v_round(acc_b); + v_int16 wacc_r = v_pack(iacc_r, iacc_r); + v_int16 wacc_g = v_pack(iacc_g, iacc_g); + v_int16 wacc_b = v_pack(iacc_b, iacc_b); + #endif + v_uint8 bacc_r = v_pack_u(wacc_r, wacc_r); + v_uint8 bacc_g = v_pack_u(wacc_g, wacc_g); + v_uint8 bacc_b = v_pack_u(wacc_b, wacc_b); + if (i + BATCH*4 <= savestorelen) { + v_store_interleave(dst, bacc_r, bacc_g, bacc_b); + continue; + } + v_store_interleave(outbuf, bacc_r, bacc_g, bacc_b); + } else if constexpr (std::is_same_v) { + v_int32 iacc_r = v_round(acc_r), iacc_g = v_round(acc_g), iacc_b = v_round(acc_b); + v_uint16 wacc_r = v_pack_u(iacc_r, iacc_r); + v_uint16 wacc_g = v_pack_u(iacc_g, iacc_g); + v_uint16 wacc_b = v_pack_u(iacc_b, iacc_b); + if (i + BATCH*2 <= savestorelen) { + v_store_interleave(dst, wacc_r, wacc_g, wacc_b); + continue; + } + v_store_interleave(outbuf, wacc_r, wacc_g, wacc_b); + } else if constexpr (std::is_same_v) { + if (i + BATCH <= savestorelen) { + v_store_interleave(dst, acc_r, acc_g, acc_b); + continue; + } + v_store_interleave(outbuf, acc_r, acc_g, acc_b); + } + } else if constexpr (NCHANNELS == 4) { + if constexpr (std::is_same_v) { + #if CV_SIMD_FP16 + v_int16 wacc_r = v_round(acc_r); + v_int16 wacc_g = v_round(acc_g); + v_int16 wacc_b = v_round(acc_b); + v_int16 wacc_a = v_round(acc_a); + #else + v_int32 iacc_r = v_round(acc_r); + v_int32 iacc_g = v_round(acc_g); + v_int32 iacc_b = v_round(acc_b); + v_int32 iacc_a = v_round(acc_a); + v_int16 wacc_r = v_pack(iacc_r, iacc_r); + v_int16 wacc_g = v_pack(iacc_g, iacc_g); + v_int16 wacc_b = v_pack(iacc_b, iacc_b); + v_int16 wacc_a = v_pack(iacc_a, iacc_a); + #endif + v_uint8 bacc_r = v_pack_u(wacc_r, wacc_r); + v_uint8 bacc_g = v_pack_u(wacc_g, wacc_g); + v_uint8 bacc_b = v_pack_u(wacc_b, wacc_b); + v_uint8 bacc_a = v_pack_u(wacc_a, wacc_a); + if (i + BATCH*4 <= savestorelen) { + v_store_interleave(dst, bacc_r, bacc_g, bacc_b, bacc_a); + continue; + } + v_store_interleave(outbuf, bacc_r, bacc_g, bacc_b, bacc_a); + } else if constexpr (std::is_same_v) { + v_int32 iacc_r = v_round(acc_r), iacc_g = v_round(acc_g); + v_int32 iacc_b = v_round(acc_b), iacc_a = v_round(acc_a); + v_uint16 wacc_r = v_pack_u(iacc_r, iacc_r); + v_uint16 wacc_g = v_pack_u(iacc_g, iacc_g); + v_uint16 wacc_b = v_pack_u(iacc_b, iacc_b); + v_uint16 wacc_a = v_pack_u(iacc_a, iacc_a); + if (i + BATCH*2 <= savestorelen) { + v_store_interleave(dst, wacc_r, wacc_g, wacc_b, wacc_a); + continue; + } + v_store_interleave(outbuf, wacc_r, wacc_g, wacc_b, wacc_a); + } else if constexpr (std::is_same_v) { + if (i + BATCH <= savestorelen) { + v_store_interleave(dst, acc_r, acc_g, acc_b, acc_a); + continue; + } + v_store_interleave(outbuf, acc_r, acc_g, acc_b, acc_a); + } + } + + memcpy(dst, outbuf, dlen*NCHANNELS*sizeof(dst[0])); + } +} +#endif + +static void +bicubic8uC1(const float* srcx, const float* srcy, int len, + const void* src, size_t srcstep, Size size, + uint8_t* dst, const float* params, + int borderType, uint8_t* borderVal) +{ +#if CV_SIMD_FP16 + bicubicVec(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +#elif CV_SIMD + bicubicVec(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +#else + bicubicRef(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +#endif +} + +static void +bicubic8uC2(const float* srcx, const float* srcy, int len, + const void* src, size_t srcstep, Size size, + uint8_t* dst, const float* params, + int borderType, uint8_t* borderVal) +{ + bicubicRef(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +} + +static void +bicubic8uC3(const float* srcx, const float* srcy, int len, + const void* src, size_t srcstep, Size size, + uint8_t* dst, const float* params, + int borderType, uint8_t* borderVal) +{ +#if CV_SIMD_FP16 + bicubicVec(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +#elif CV_SIMD + bicubicVec(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +#else + bicubicRef(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +#endif +} + +static void +bicubic8uC4(const float* srcx, const float* srcy, int len, + const void* src, size_t srcstep, Size size, + uint8_t* dst, const float* params, + int borderType, uint8_t* borderVal) +{ +#if CV_SIMD_FP16 + bicubicVec(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +#elif CV_SIMD + bicubicVec(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +#else + bicubicRef(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +#endif +} + +static void +bicubic16uC1(const float* srcx, const float* srcy, int len, + const void* src, size_t srcstep, Size size, + uint16_t* dst, const float* params, + int borderType, uint16_t* borderVal) +{ +#if CV_SIMD + bicubicVec(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +#else + bicubicRef(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +#endif +} + +static void +bicubic16uC2(const float* srcx, const float* srcy, int len, + const void* src, size_t srcstep, Size size, + uint16_t* dst, const float* params, + int borderType, uint16_t* borderVal) +{ + bicubicRef(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +} + +static void +bicubic16uC3(const float* srcx, const float* srcy, int len, + const void* src, size_t srcstep, Size size, + uint16_t* dst, const float* params, + int borderType, uint16_t* borderVal) +{ +#if CV_SIMD + bicubicVec(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +#else + bicubicRef(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +#endif +} + +static void +bicubic16uC4(const float* srcx, const float* srcy, int len, + const void* src, size_t srcstep, Size size, + uint16_t* dst, const float* params, + int borderType, uint16_t* borderVal) +{ +#if CV_SIMD + bicubicVec(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +#else + bicubicRef(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +#endif +} + +static void +bicubic16sC1(const float* srcx, const float* srcy, int len, + const void* src, size_t srcstep, Size size, + int16_t* dst, const float* params, + int borderType, int16_t* borderVal) +{ + bicubicRef(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +} + +static void +bicubic16sC2(const float* srcx, const float* srcy, int len, + const void* src, size_t srcstep, Size size, + int16_t* dst, const float* params, + int borderType, int16_t* borderVal) +{ + bicubicRef(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +} + +static void +bicubic16sC3(const float* srcx, const float* srcy, int len, + const void* src, size_t srcstep, Size size, + int16_t* dst, const float* params, + int borderType, int16_t* borderVal) +{ + bicubicRef(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +} + +static void +bicubic16sC4(const float* srcx, const float* srcy, int len, + const void* src, size_t srcstep, Size size, + int16_t* dst, const float* params, + int borderType, int16_t* borderVal) +{ + bicubicRef(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +} + +static void +bicubic32fC1(const float* srcx, const float* srcy, int len, + const void* src, size_t srcstep, Size size, + float* dst, const float* params, + int borderType, float* borderVal) +{ +#if CV_SIMD + bicubicVec(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +#else + bicubicRef(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +#endif +} + +static void +bicubic32fC2(const float* srcx, const float* srcy, int len, + const void* src, size_t srcstep, Size size, + float* dst, const float* params, + int borderType, float* borderVal) +{ + bicubicRef(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +} + +static void +bicubic32fC3(const float* srcx, const float* srcy, int len, + const void* src, size_t srcstep, Size size, + float* dst, const float* params, + int borderType, float* borderVal) +{ +#if CV_SIMD + bicubicVec(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +#else + bicubicRef(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +#endif +} + +static void +bicubic32fC4(const float* srcx, const float* srcy, int len, + const void* src, size_t srcstep, Size size, + float* dst, const float* params, + int borderType, float* borderVal) +{ +#if CV_SIMD + bicubicVec(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +#else + bicubicRef(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +#endif +} + +static void +bicubic64fC1(const float* srcx, const float* srcy, int len, + const void* src, size_t srcstep, Size size, + double* dst, const float* params, + int borderType, double* borderVal) +{ + bicubicRef(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +} + +static void +bicubic64fC2(const float* srcx, const float* srcy, int len, + const void* src, size_t srcstep, Size size, + double* dst, const float* params, + int borderType, double* borderVal) +{ + bicubicRef(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +} + +static void +bicubic64fC3(const float* srcx, const float* srcy, int len, + const void* src, size_t srcstep, Size size, + double* dst, const float* params, + int borderType, double* borderVal) +{ + bicubicRef(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +} + +static void +bicubic64fC4(const float* srcx, const float* srcy, int len, + const void* src, size_t srcstep, Size size, + double* dst, const float* params, + int borderType, double* borderVal) +{ + bicubicRef(srcx, srcy, len, src, srcstep, size, + dst, params, borderType, borderVal); +} + +} + +ImgWarpFunc getBicubicWarpFunc_(int type) +{ + if (type == CV_8UC1) { + return (ImgWarpFunc)bicubic8uC1; + } + if (type == CV_8UC2) { + return (ImgWarpFunc)bicubic8uC2; + } + if (type == CV_8UC3) { + return (ImgWarpFunc)bicubic8uC3; + } + if (type == CV_8UC4) { + return (ImgWarpFunc)bicubic8uC4; + } + if (type == CV_16UC1) { + return (ImgWarpFunc)bicubic16uC1; + } + if (type == CV_16UC2) { + return (ImgWarpFunc)bicubic16uC2; + } + if (type == CV_16UC3) { + return (ImgWarpFunc)bicubic16uC3; + } + if (type == CV_16UC4) { + return (ImgWarpFunc)bicubic16uC4; + } + if (type == CV_16SC1) { + return (ImgWarpFunc)bicubic16sC1; + } + if (type == CV_16SC2) { + return (ImgWarpFunc)bicubic16sC2; + } + if (type == CV_16SC3) { + return (ImgWarpFunc)bicubic16sC3; + } + if (type == CV_16SC4) { + return (ImgWarpFunc)bicubic16sC4; + } + if (type == CV_32FC1) { + return (ImgWarpFunc)bicubic32fC1; + } + if (type == CV_32FC2) { + return (ImgWarpFunc)bicubic32fC2; + } + if (type == CV_32FC3) { + return (ImgWarpFunc)bicubic32fC3; + } + if (type == CV_32FC4) { + return (ImgWarpFunc)bicubic32fC4; + } + if (type == CV_64FC1) { + return (ImgWarpFunc)bicubic64fC1; + } + if (type == CV_64FC2) { + return (ImgWarpFunc)bicubic64fC2; + } + if (type == CV_64FC3) { + return (ImgWarpFunc)bicubic64fC3; + } + if (type == CV_64FC4) { + return (ImgWarpFunc)bicubic64fC4; + } + return (ImgWarpFunc)nullptr; +} + #endif // CV_CPU_OPTIMIZATION_DECLARATIONS_ONLY diff --git a/modules/imgproc/test/ocl/test_warp.cpp b/modules/imgproc/test/ocl/test_warp.cpp index 30f1df68d3..a5add1a69c 100644 --- a/modules/imgproc/test/ocl/test_warp.cpp +++ b/modules/imgproc/test/ocl/test_warp.cpp @@ -107,7 +107,7 @@ PARAM_TEST_CASE(WarpTestBase, MatType, Interpolation, bool, bool) void Near(double threshold = 0.0) { if (depth < CV_32F) - EXPECT_MAT_N_DIFF(dst_roi, udst_roi, cvRound(dst_roi.total()*threshold)); + EXPECT_MAT_N_DIFF_EPS(dst_roi, udst_roi, 1, cvRound(dst_roi.total()*threshold)); else OCL_EXPECT_MATS_NEAR_RELATIVE(dst, threshold); } diff --git a/modules/imgproc/test/test_imgwarp.cpp b/modules/imgproc/test/test_imgwarp.cpp index cd70bf3581..6975c31a3b 100644 --- a/modules/imgproc/test/test_imgwarp.cpp +++ b/modules/imgproc/test/test_imgwarp.cpp @@ -1114,5 +1114,90 @@ TEST(Imgproc_getPerspectiveTransform, issue_26916) EXPECT_MAT_NEAR(obtained_homogeneous_dst_points, expected_homogeneous_dst_points, 1e-10); } +static void rotation2affine(float scale, float angle, float cx, float cy, float* M) +{ + // (x - cx)*cos(a) + (y - cy)*sin(a) + cx + // -(x - cx)*sin(a) + (y - cy)*cos(a) + cy + float ca = cosf(angle), sa = sinf(angle); + M[0] = scale*ca; M[1] = scale*sa; M[2] = scale*(-cx*ca - cy*sa) + cx; + M[3] = -scale*sa; M[4] = scale*ca; M[5] = scale*(cx*sa - cy*ca) + cy; +} + +TEST(Imgproc_Warping, DISABLED_playground) +{ + int imgtype = CV_32F; + int imgcn = 3; + bool useOpenCL = true; + + auto ts = cvtest::TS::ptr(); + Mat img0 = imread(string(ts->get_data_path()) + "stereomatching/datasets/tsukuba/im2.png", 1), img1, img; + int iangle = -1; + int borderType = BORDER_CONSTANT; + Scalar borderValue(0, 128, 0); + + double cvtscale = imgtype == CV_16U ? 256. : imgtype == CV_32F ? 1./255 : 1.; + if (imgcn == 1) { + cvtColor(img0, img1, COLOR_BGR2GRAY); + } else if (imgcn == 4) { + cvtColor(img0, img1, COLOR_BGR2BGRA); + } else if (imgcn == 3) { + img1 = img0; + } else { + CV_Assert(imgcn == 2); + std::vector ch; + split(img0, ch); + ch.pop_back(); + merge(ch, img1); + } + img1.convertTo(img, imgtype, cvtscale); + Mat canvas0(img.size(), imgtype), canvas8; + float cx = img.cols*0.5f, cy = img.rows*0.5f; + if (img.depth() == CV_32F) { + borderValue = Scalar(100*cvtscale, 0*cvtscale, 100*cvtscale); + } + UMat uimg, ucanvas; + if (useOpenCL) { + img.copyTo(uimg); + } + + for(;;) { + iangle = (iangle + 1) % (360*4); + float angle = float(iangle*CV_PI/180.f*0.25f); + float scale = float(1 + 0.2f*sin(angle)); + float Mdata[9] = {1, 0, 0, 0, 1, 0, 0, 0, 1}; + rotation2affine(scale, angle, cx, cy, Mdata); + Mat M(2, 3, CV_32F, Mdata); + + double t0 = (double)getTickCount(); + if (!useOpenCL) { + warpAffine(img, canvas0, M, canvas0.size(), INTER_CUBIC, borderType, borderValue); + } else { + warpAffine(uimg, ucanvas, M, uimg.size(), INTER_CUBIC, borderType, borderValue); + ocl::finish(); + } + t0 = (double)getTickCount() - t0; + + if (useOpenCL) { + ucanvas.copyTo(canvas0); + } + canvas0.convertTo(canvas8, CV_8U, 1./cvtscale); + if (canvas8.channels() == 2) { + std::vector ch; + split(canvas8, ch); + ch.push_back(Mat::zeros(canvas8.size(), CV_8U)); + merge(ch, canvas8); + } + printf("opencv time = %.1fms\n", t0*1000./getTickFrequency()); + imshow("result (opencv)", canvas8); + int c = waitKey(1); + if (c < 0) + continue; + if ((c & 255) == 27) + break; + if ((waitKey() & 255) == 27) + break; + } +} + }} // namespace /* End of file. */ diff --git a/modules/imgproc/test/test_imgwarp_strict.cpp b/modules/imgproc/test/test_imgwarp_strict.cpp index 92c76e254f..58cc1aa9d4 100644 --- a/modules/imgproc/test/test_imgwarp_strict.cpp +++ b/modules/imgproc/test/test_imgwarp_strict.cpp @@ -231,13 +231,13 @@ void CV_ImageWarpBaseTest::run(int) float CV_ImageWarpBaseTest::get_success_error_level(int _interpolation, int) const { if (_interpolation == INTER_CUBIC) - return 1.0f; + return 2.0f; else if (_interpolation == INTER_LANCZOS4) return 1.0f; else if (_interpolation == INTER_NEAREST) return 255.0f; // FIXIT: check is not reliable for Black/White (0/255) images else if (_interpolation == INTER_AREA) - return 2.0f; + return 1.0f; else return 1.0f; } @@ -304,7 +304,7 @@ void CV_ImageWarpBaseTest::validate_results() const const int radius = 3; int rmin = MAX(dy - radius, 0), rmax = MIN(dy + radius, dsize.height); - int cmin = MAX(dx / cn - radius, 0), cmax = MIN(dx / cn + radius, dsize.width); + int cmin = MAX(dx / cn - radius, 0), cmax = MIN(dx / cn + radius, dsize.width/cn); std::cout << "opencv result:\n" << dst(Range(rmin, rmax), Range(cmin, cmax)) << std::endl; std::cout << "reference result:\n" << reference_dst(Range(rmin, rmax), Range(cmin, cmax)) << std::endl; @@ -851,16 +851,23 @@ void CV_Remap_Test::run_func() void CV_Remap_Test::convert_maps() { - if (mapx.type() != CV_16SC2) - convertMaps(mapx.clone(), mapy.clone(), mapx, mapy, CV_16SC2, interpolation == INTER_NEAREST); - else if (interpolation != INTER_NEAREST) - if (mapy.type() != CV_16UC1) - mapy.clone().convertTo(mapy, CV_16UC1); + if (interpolation == INTER_CUBIC) { + if (mapx.type() != CV_32FC1) + convertMaps(mapx.clone(), mapy.clone(), mapx, mapy, CV_32FC1, interpolation == INTER_NEAREST); + } + else { + if (mapx.type() != CV_16SC2) + convertMaps(mapx.clone(), mapy.clone(), mapx, mapy, CV_16SC2, interpolation == INTER_NEAREST); + else if (interpolation != INTER_NEAREST) { + if (mapy.type() != CV_16UC1) + mapy.clone().convertTo(mapy, CV_16UC1); + } - if (interpolation == INTER_NEAREST) - mapy = Mat(); - CV_Assert(((interpolation == INTER_NEAREST && mapy.empty()) || mapy.type() == CV_16UC1 || - mapy.type() == CV_16SC1) && mapx.type() == CV_16SC2); + if (interpolation == INTER_NEAREST) + mapy = Mat(); + CV_Assert(((interpolation == INTER_NEAREST && mapy.empty()) || mapy.type() == CV_16UC1 || + mapy.type() == CV_16SC1) && mapx.type() == CV_16SC2); + } } const char* CV_Remap_Test::borderType_to_string() const @@ -1164,8 +1171,10 @@ void CV_Remap_Test::remap_nearest(const Mat& _src, Mat& _dst) void CV_Remap_Test::remap_generic(const Mat& _src, Mat& _dst) { - CV_Assert(mapx.type() == CV_16SC2 && mapy.type() == CV_16UC1); + CV_Assert((mapx.type() == CV_16SC2 && mapy.type() == CV_16UC1) || + (mapx.type() == CV_32FC1 && mapy.type() == CV_32FC1)); + bool fixedpt = mapx.type() == CV_16SC2; int ksize = 2; if (interpolation == INTER_CUBIC) ksize = 4; @@ -1185,19 +1194,42 @@ void CV_Remap_Test::remap_generic(const Mat& _src, Mat& _dst) for (int dy = 0; dy < dsize.height; ++dy) { - const short* yMx = mapx.ptr(dy); - const ushort* yMy = mapy.ptr(dy); + const short* Mxy = nullptr; + const ushort* Mab = nullptr; + const float* Mx = nullptr; + const float* My = nullptr; + + if (fixedpt) { + Mxy = mapx.ptr(dy); + Mab = mapy.ptr(dy); + } else { + Mx = mapx.ptr(dy); + My = mapy.ptr(dy); + } float* yD = _dst.ptr(dy); for (int dx = 0; dx < dsize.width; ++dx) { float* xyD = yD + dx * cn; - float sx = yMx[dx * 2], sy = yMx[dx * 2 + 1]; - int isx = cvFloor(sx), isy = cvFloor(sy); + int isx, isy; + float alpha, beta; + if (fixedpt) { + isx = Mxy[dx * 2]; + isy = Mxy[dx * 2 + 1]; + alpha = (Mab[dx] & (INTER_TAB_SIZE - 1)) / static_cast(INTER_TAB_SIZE); + beta = ((Mab[dx] >> INTER_BITS) & (INTER_TAB_SIZE - 1)) / static_cast(INTER_TAB_SIZE); + } else { + float sx = Mx[dx]; + float sy = My[dx]; + isx = cvFloor(sx); + isy = cvFloor(sy); + alpha = sx - isx; + beta = sy - isy; + } - inter_func((yMy[dx] & (INTER_TAB_SIZE - 1)) / static_cast(INTER_TAB_SIZE), w); - inter_func(((yMy[dx] >> INTER_BITS) & (INTER_TAB_SIZE - 1)) / static_cast(INTER_TAB_SIZE), w + ksize); + inter_func(alpha, w); + inter_func(beta, w + ksize); isx -= ofs; isy -= ofs; @@ -1332,7 +1364,7 @@ void CV_WarpAffine_Test::run_func() float CV_WarpAffine_Test::get_success_error_level(int _interpolation, int _depth) const { - return _depth == CV_8U ? 0.f : CV_ImageWarpBaseTest::get_success_error_level(_interpolation, _depth); + return _depth == CV_8U ? 2.f : CV_ImageWarpBaseTest::get_success_error_level(_interpolation, _depth); } void CV_WarpAffine_Test::run_reference_func() @@ -1390,6 +1422,44 @@ void CV_WarpAffine_Test::new_warpAffine(const Mat &_src, Mat &_dst, const Mat &t } } +static void computeWarpMaps(const Mat& M, Mat& mapx, Mat& mapy) +{ + CV_Assert(M.type() == CV_64F); + CV_Assert(M.size() == Size(3, 2) || + M.size() == Size(3, 3)); + CV_Assert(mapx.size() == mapy.size()); + CV_Assert(mapx.type() == CV_32FC1 && mapy.type() == CV_32FC1); + + bool perspective = M.rows == 3; + const double* dataM = M.ptr(); + Size dsize = mapx.size(); + + for (int dy = 0; dy < dsize.height; ++dy) + { + float* mapxrow = mapx.ptr(dy); + float* mapyrow = mapy.ptr(dy); + for (int dx = 0; dx < dsize.width; ++dx) + { + double x = dataM[0] * dx + dataM[1] * dy + dataM[2]; + double y = dataM[3] * dx + dataM[4] * dy + dataM[5]; + + if (perspective) { + double z = dataM[6] * dx + dataM[7] * dy + dataM[8]; + if (z != 0.0) { + x /= z; + y /= z; + } else { + x = 0.0; + y = 0.0; + } + } + + mapxrow[dx] = float(x); + mapyrow[dx] = float(y); + } + } +} + void CV_WarpAffine_Test::warpAffine(const Mat& _src, Mat& _dst) { Size dsize = _dst.size(); @@ -1405,12 +1475,6 @@ void CV_WarpAffine_Test::warpAffine(const Mat& _src, Mat& _dst) if (inter == INTER_AREA) inter = INTER_LINEAR; - mapx.create(dsize, CV_16SC2); - if (inter != INTER_NEAREST) - mapy.create(dsize, CV_16SC1); - else - mapy = Mat(); - if (!(interpolation & cv::WARP_INVERSE_MAP)) invertAffineTransform(tM.clone(), tM); @@ -1425,32 +1489,44 @@ void CV_WarpAffine_Test::warpAffine(const Mat& _src, Mat& _dst) } } - const int AB_BITS = MAX(10, (int)INTER_BITS); - const int AB_SCALE = 1 << AB_BITS; - int round_delta = (inter == INTER_NEAREST) ? AB_SCALE / 2 : (AB_SCALE / INTER_TAB_SIZE / 2); + if (inter == INTER_CUBIC) { + mapx.create(dsize, CV_32FC1); + mapy.create(dsize, CV_32FC1); + computeWarpMaps(tM, mapx, mapy); + } else { + mapx.create(dsize, CV_16SC2); + if (inter != INTER_NEAREST) + mapy.create(dsize, CV_16SC1); + else + mapy = Mat(); - const softdouble* data_tM = tM.ptr(0); - for (int dy = 0; dy < dsize.height; ++dy) - { - short* yM = mapx.ptr(dy); - for (int dx = 0; dx < dsize.width; ++dx, yM += 2) + const int AB_BITS = MAX(10, (int)INTER_BITS); + const int AB_SCALE = 1 << AB_BITS; + int round_delta = (inter == INTER_NEAREST) ? AB_SCALE / 2 : (AB_SCALE / INTER_TAB_SIZE / 2); + + const softdouble* data_tM = tM.ptr(0); + for (int dy = 0; dy < dsize.height; ++dy) { - int v1 = saturate_cast(saturate_cast(data_tM[0] * dx * AB_SCALE) + - saturate_cast((data_tM[1] * dy + data_tM[2]) * AB_SCALE) + round_delta), + short* yM = mapx.ptr(dy); + for (int dx = 0; dx < dsize.width; ++dx, yM += 2) + { + int v1 = saturate_cast(saturate_cast(data_tM[0] * dx * AB_SCALE) + + saturate_cast((data_tM[1] * dy + data_tM[2]) * AB_SCALE) + round_delta), v2 = saturate_cast(saturate_cast(data_tM[3] * dx * AB_SCALE) + - saturate_cast((data_tM[4] * dy + data_tM[5]) * AB_SCALE) + round_delta); - v1 >>= AB_BITS - INTER_BITS; - v2 >>= AB_BITS - INTER_BITS; + saturate_cast((data_tM[4] * dy + data_tM[5]) * AB_SCALE) + round_delta); + v1 >>= AB_BITS - INTER_BITS; + v2 >>= AB_BITS - INTER_BITS; - yM[0] = saturate_cast(v1 >> INTER_BITS); - yM[1] = saturate_cast(v2 >> INTER_BITS); + yM[0] = saturate_cast(v1 >> INTER_BITS); + yM[1] = saturate_cast(v2 >> INTER_BITS); - if (inter != INTER_NEAREST) - mapy.ptr(dy)[dx] = ((v2 & (INTER_TAB_SIZE - 1)) * INTER_TAB_SIZE + (v1 & (INTER_TAB_SIZE - 1))); + if (inter != INTER_NEAREST) + mapy.ptr(dy)[dx] = ((v2 & (INTER_TAB_SIZE - 1)) * INTER_TAB_SIZE + (v1 & (INTER_TAB_SIZE - 1))); + } } - } - CV_Assert(mapx.type() == CV_16SC2 && ((inter == INTER_NEAREST && mapy.empty()) || mapy.type() == CV_16SC1)); + CV_Assert(mapx.type() == CV_16SC2 && ((inter == INTER_NEAREST && mapy.empty()) || mapy.type() == CV_16SC1)); + } cv::remap(_src, _dst, mapx, mapy, inter, borderType, borderValue); } @@ -1612,41 +1688,47 @@ void CV_WarpPerspective_Test::warpPerspective(const Mat& _src, Mat& _dst) } } - mapx.create(dsize, CV_16SC2); - if (inter != INTER_NEAREST) - mapy.create(dsize, CV_16SC1); - else - mapy = Mat(); + if (inter == INTER_CUBIC) { + mapx.create(dsize, CV_32FC1); + mapy.create(dsize, CV_32FC1); + computeWarpMaps(M, mapx, mapy); + } else { + mapx.create(dsize, CV_16SC2); + if (inter != INTER_NEAREST) + mapy.create(dsize, CV_16SC1); + else + mapy = Mat(); - double* tM = M.ptr(0); - for (int dy = 0; dy < dsize.height; ++dy) - { - short* yMx = mapx.ptr(dy); - - for (int dx = 0; dx < dsize.width; ++dx, yMx += 2) + double* tM = M.ptr(0); + for (int dy = 0; dy < dsize.height; ++dy) { - double den = tM[6] * dx + tM[7] * dy + tM[8]; - den = den ? 1.0 / den : 0.0; + short* yMx = mapx.ptr(dy); - if (inter == INTER_NEAREST) + for (int dx = 0; dx < dsize.width; ++dx, yMx += 2) { - yMx[0] = saturate_cast((tM[0] * dx + tM[1] * dy + tM[2]) * den); - yMx[1] = saturate_cast((tM[3] * dx + tM[4] * dy + tM[5]) * den); - continue; + double den = tM[6] * dx + tM[7] * dy + tM[8]; + den = den ? 1.0 / den : 0.0; + + if (inter == INTER_NEAREST) + { + yMx[0] = saturate_cast((tM[0] * dx + tM[1] * dy + tM[2]) * den); + yMx[1] = saturate_cast((tM[3] * dx + tM[4] * dy + tM[5]) * den); + continue; + } + + den *= static_cast(INTER_TAB_SIZE); + int v0 = saturate_cast((tM[0] * dx + tM[1] * dy + tM[2]) * den); + int v1 = saturate_cast((tM[3] * dx + tM[4] * dy + tM[5]) * den); + + yMx[0] = saturate_cast(v0 >> INTER_BITS); + yMx[1] = saturate_cast(v1 >> INTER_BITS); + mapy.ptr(dy)[dx] = saturate_cast((v1 & (INTER_TAB_SIZE - 1)) * + INTER_TAB_SIZE + (v0 & (INTER_TAB_SIZE - 1))); } - - den *= static_cast(INTER_TAB_SIZE); - int v0 = saturate_cast((tM[0] * dx + tM[1] * dy + tM[2]) * den); - int v1 = saturate_cast((tM[3] * dx + tM[4] * dy + tM[5]) * den); - - yMx[0] = saturate_cast(v0 >> INTER_BITS); - yMx[1] = saturate_cast(v1 >> INTER_BITS); - mapy.ptr(dy)[dx] = saturate_cast((v1 & (INTER_TAB_SIZE - 1)) * - INTER_TAB_SIZE + (v0 & (INTER_TAB_SIZE - 1))); } - } - CV_Assert(mapx.type() == CV_16SC2 && ((inter == INTER_NEAREST && mapy.empty()) || mapy.type() == CV_16SC1)); + CV_Assert(mapx.type() == CV_16SC2 && ((inter == INTER_NEAREST && mapy.empty()) || mapy.type() == CV_16SC1)); + } cv::remap(_src, _dst, mapx, mapy, inter, borderType, borderValue); } diff --git a/modules/ts/include/opencv2/ts/ocl_test.hpp b/modules/ts/include/opencv2/ts/ocl_test.hpp index 230ba86280..48d753b36c 100644 --- a/modules/ts/include/opencv2/ts/ocl_test.hpp +++ b/modules/ts/include/opencv2/ts/ocl_test.hpp @@ -110,6 +110,18 @@ do \ << "Size: " << (mat1).size() << std::endl; \ } while ((void)0, 0) +#define EXPECT_MAT_N_DIFF_EPS(mat1, mat2, eps, num) \ +do \ +{ \ + ASSERT_EQ(mat1.type(), mat2.type()); \ + ASSERT_EQ(mat1.size(), mat2.size()); \ + Mat diff, mask; \ + absdiff(mat1, mat2, diff); \ + cv::compare(diff.reshape(1), Scalar::all(eps), mask, CMP_GT); \ + EXPECT_LE(countNonZero(mask), num) \ + << "Size: " << mat1.size() << std::endl; \ +} while ((void)0, 0) + #define EXPECT_MAT_N_DIFF(mat1, mat2, num) \ do \ { \