From 77912645254ee6d99e398eca7b518d22c672bdd1 Mon Sep 17 00:00:00 2001 From: Alexander Karsakov Date: Mon, 21 Jul 2014 10:42:56 +0400 Subject: [PATCH] Added multi-block scheme --- modules/core/perf/opencl/perf_dxt.cpp | 6 +- modules/core/src/dxt.cpp | 49 +++-- modules/core/src/opencl/fft.cl | 282 +++++++++++++++++++++++++- modules/core/test/ocl/test_dft.cpp | 14 +- 4 files changed, 321 insertions(+), 30 deletions(-) diff --git a/modules/core/perf/opencl/perf_dxt.cpp b/modules/core/perf/opencl/perf_dxt.cpp index 3980a191fa..f4b6b49a9f 100644 --- a/modules/core/perf/opencl/perf_dxt.cpp +++ b/modules/core/perf/opencl/perf_dxt.cpp @@ -65,10 +65,10 @@ enum OCL_FFT_TYPE typedef tuple DftParams; typedef TestBaseWithParam DftFixture; -OCL_PERF_TEST_P(DftFixture, Dft, ::testing::Combine(Values(C2C, R2R, C2R, R2C), +OCL_PERF_TEST_P(DftFixture, Dft, ::testing::Combine(Values(C2C/*, R2R, C2R, R2C*/), Values(OCL_SIZE_1, OCL_SIZE_2, OCL_SIZE_3, Size(1024, 1024), Size(512, 512), Size(2048, 2048)), - Values((int)DFT_ROWS, (int) 0, (int)DFT_SCALE/*, (int)DFT_INVERSE, - (int)DFT_INVERSE | DFT_SCALE, (int)DFT_ROWS | DFT_INVERSE*/))) + Values((int) 0, (int)DFT_ROWS, (int)DFT_SCALE, (int)DFT_INVERSE, + /*(int)DFT_INVERSE | DFT_SCALE,*/ (int)DFT_ROWS | DFT_INVERSE))) { const DftParams params = GetParam(); const int dft_type = get<0>(params); diff --git a/modules/core/src/dxt.cpp b/modules/core/src/dxt.cpp index 879a70613f..d5b1cb383d 100644 --- a/modules/core/src/dxt.cpp +++ b/modules/core/src/dxt.cpp @@ -2041,23 +2041,33 @@ static std::vector ocl_getRadixes(int cols, std::vector& radixes, std: int n = 1; int factor_index = 0; + min_radix = INT_MAX; // 2^n transforms if ( (factors[factor_index] & 1) == 0 ) { for( ; n < factors[factor_index]; ) { - int radix = 2; + int radix = 2, block = 1; if (8*n <= factors[0]) radix = 8; else if (4*n <= factors[0]) + { radix = 4; + if (cols % 8 == 0) + block = 2; + } + else + { + if (cols % 8 == 0) + block = 4; + else if (cols % 4 == 0) + block = 2; + } radixes.push_back(radix); - if (radix == 2 && cols % 4 == 0) - min_radix = min(min_radix, 2*radix); - else - min_radix = min(min_radix, radix); + blocks.push_back(block); + min_radix = min(min_radix, block*radix); n *= radix; } factor_index++; @@ -2066,11 +2076,22 @@ static std::vector ocl_getRadixes(int cols, std::vector& radixes, std: // all the other transforms for( ; factor_index < nf; factor_index++ ) { - radixes.push_back(factors[factor_index]); - if (factors[factor_index] == 3 && cols % 6 == 0) - min_radix = min(min_radix, 2*factors[factor_index]); - else - min_radix = min(min_radix, factors[factor_index]); + int radix = factors[factor_index], block = 1; + if (radix == 3) + { + if (cols % 12 == 0) + block = 4; + else if (cols % 6 == 0) + block = 2; + } + else if (radix == 5) + { + if (cols % 10 == 0) + block = 2; + } + radixes.push_back(radix); + blocks.push_back(block); + min_radix = min(min_radix, block*radix); } return radixes; } @@ -2086,7 +2107,7 @@ struct OCL_FftPlan bool status; OCL_FftPlan(int _size, int _flags): dft_size(_size), flags(_flags), status(true) { - int min_radix = INT_MAX; + int min_radix; std::vector radixes, blocks; ocl_getRadixes(dft_size, radixes, blocks, min_radix); thread_count = (dft_size + min_radix-1) / min_radix; @@ -2102,9 +2123,9 @@ struct OCL_FftPlan int n = 1, twiddle_size = 0; for (size_t i=0; i 1) + radix_processing += format("fft_radix%d_B%d(smem,twiddles+%d,ind,%d,%d);", radix, block, twiddle_size, n, dft_size/radix); else radix_processing += format("fft_radix%d(smem,twiddles+%d,ind,%d,%d);", radix, twiddle_size, n, dft_size/radix); twiddle_size += (radix-1)*n; diff --git a/modules/core/src/opencl/fft.cl b/modules/core/src/opencl/fft.cl index 8aecfc056a..fdbad19ad0 100644 --- a/modules/core/src/opencl/fft.cl +++ b/modules/core/src/opencl/fft.cl @@ -44,11 +44,11 @@ __attribute__((always_inline)) void fft_radix2_B2(__local float2* smem, __constant const float2* twiddles, const int x, const int block_size, const int t) { const int k1 = x & (block_size - 1); - const int x2 = x + (t+1)/2; + const int x2 = x + t/2; const int k2 = x2 & (block_size - 1); float2 a0, a1, a2, a3; - if (x < (t+1)/2) + if (x < t/2) { a0 = smem[x]; a1 = mul_float2(twiddles[k1],smem[x+t]); @@ -58,7 +58,7 @@ void fft_radix2_B2(__local float2* smem, __constant const float2* twiddles, cons barrier(CLK_LOCAL_MEM_FENCE); - if (x < (t+1)/2) + if (x < t/2) { int dst_ind = (x << 1) - k1; smem[dst_ind] = a0 + a1; @@ -72,6 +72,55 @@ void fft_radix2_B2(__local float2* smem, __constant const float2* twiddles, cons barrier(CLK_LOCAL_MEM_FENCE); } +__attribute__((always_inline)) +void fft_radix2_B4(__local float2* smem, __constant const float2* twiddles, const int x, const int block_size, const int t) +{ + const int thread_block = t/4; + const int k1 = x & (block_size - 1); + const int x2 = x + thread_block; + const int k2 = x2 & (block_size - 1); + const int x3 = x + 2*thread_block; + const int k3 = x3 & (block_size - 1); + const int x4 = x + 3*thread_block; + const int k4 = x4 & (block_size - 1); + float2 a0, a1, a2, a3, a4, a5, a6, a7; + + if (x < t/4) + { + a0 = smem[x]; + a1 = mul_float2(twiddles[k1],smem[x+t]); + a2 = smem[x2]; + a3 = mul_float2(twiddles[k2],smem[x2+t]); + a4 = smem[x3]; + a5 = mul_float2(twiddles[k3],smem[x3+t]); + a6 = smem[x4]; + a7 = mul_float2(twiddles[k4],smem[x4+t]); + } + + barrier(CLK_LOCAL_MEM_FENCE); + + if (x < t/4) + { + int dst_ind = (x << 1) - k1; + smem[dst_ind] = a0 + a1; + smem[dst_ind+block_size] = a0 - a1; + + dst_ind = (x2 << 1) - k2; + smem[dst_ind] = a2 + a3; + smem[dst_ind+block_size] = a2 - a3; + + dst_ind = (x3 << 1) - k3; + smem[dst_ind] = a4 + a5; + smem[dst_ind+block_size] = a4 - a5; + + dst_ind = (x4 << 1) - k4; + smem[dst_ind] = a6 + a7; + smem[dst_ind+block_size] = a6 - a7; + } + + barrier(CLK_LOCAL_MEM_FENCE); +} + __attribute__((always_inline)) void fft_radix4(__local float2* smem, __constant const float2* twiddles, const int x, const int block_size, const int t) { @@ -107,6 +156,58 @@ void fft_radix4(__local float2* smem, __constant const float2* twiddles, const i barrier(CLK_LOCAL_MEM_FENCE); } +__attribute__((always_inline)) +void fft_radix4_B2(__local float2* smem, __constant const float2* twiddles, const int x, const int block_size, const int t) +{ + const int k = x & (block_size - 1); + const int x2 = x + t/2; + const int k2 = x2 & (block_size - 1); + float2 a0, a1, a2, a3, a4, a5, a6, a7; + + if (x < t/2) + { + a0 = smem[x]; + a1 = mul_float2(twiddles[k], smem[x+t]); + a2 = mul_float2(twiddles[k + block_size],smem[x+2*t]); + a3 = mul_float2(twiddles[k + 2*block_size],smem[x+3*t]); + + a4 = smem[x2]; + a5 = mul_float2(twiddles[k2], smem[x2+t]); + a6 = mul_float2(twiddles[k2 + block_size],smem[x2+2*t]); + a7 = mul_float2(twiddles[k2 + 2*block_size],smem[x2+3*t]); + } + + barrier(CLK_LOCAL_MEM_FENCE); + + if (x < t/2) + { + int dst_ind = ((x - k) << 2) + k; + + float2 b0 = a0 + a2; + a2 = a0 - a2; + float2 b1 = a1 + a3; + a3 = twiddle(a1 - a3); + + smem[dst_ind] = b0 + b1; + smem[dst_ind + block_size] = a2 + a3; + smem[dst_ind + 2*block_size] = b0 - b1; + smem[dst_ind + 3*block_size] = a2 - a3; + + dst_ind = ((x2 - k2) << 2) + k2; + b0 = a4 + a6; + a6 = a4 - a6; + b1 = a5 + a7; + a7 = twiddle(a5 - a7); + + smem[dst_ind] = b0 + b1; + smem[dst_ind + block_size] = a6 + a7; + smem[dst_ind + 2*block_size] = b0 - b1; + smem[dst_ind + 3*block_size] = a6 - a7; + } + + barrier(CLK_LOCAL_MEM_FENCE); +} + __attribute__((always_inline)) void fft_radix8(__local float2* smem, __constant const float2* twiddles, const int x, const int block_size, const int t) { @@ -205,11 +306,11 @@ __attribute__((always_inline)) void fft_radix3_B2(__local float2* smem, __constant const float2* twiddles, const int x, const int block_size, const int t) { const int k = x % block_size; - const int x2 = x + (t+1)/2; + const int x2 = x + t/2; const int k2 = x2 % block_size; float2 a0, a1, a2, a3, a4, a5; - if (x < (t+1)/2) + if (x < t/2) { a0 = smem[x]; a1 = mul_float2(twiddles[k], smem[x+t]); @@ -222,7 +323,7 @@ void fft_radix3_B2(__local float2* smem, __constant const float2* twiddles, cons barrier(CLK_LOCAL_MEM_FENCE); - if (x < (t+1)/2) + if (x < t/2) { int dst_ind = ((x - k) * 3) + k; @@ -248,6 +349,86 @@ void fft_radix3_B2(__local float2* smem, __constant const float2* twiddles, cons barrier(CLK_LOCAL_MEM_FENCE); } +__attribute__((always_inline)) +void fft_radix3_B4(__local float2* smem, __constant const float2* twiddles, const int x, const int block_size, const int t) +{ + const int thread_block = t/4; + const int k = x % block_size; + const int x2 = x + thread_block; + const int k2 = x2 % block_size; + const int x3 = x + 2*thread_block; + const int k3 = x3 % block_size; + const int x4 = x + 3*thread_block; + const int k4 = x4 % block_size; + float2 a0, a1, a2, a3, a4, a5, a6, a7, a8, a9, a10, a11; + + if (x < t/4) + { + a0 = smem[x]; + a1 = mul_float2(twiddles[k], smem[x+t]); + a2 = mul_float2(twiddles[k+block_size], smem[x+2*t]); + + a3 = smem[x2]; + a4 = mul_float2(twiddles[k2], smem[x2+t]); + a5 = mul_float2(twiddles[k2+block_size], smem[x2+2*t]); + + a6 = smem[x3]; + a7 = mul_float2(twiddles[k3], smem[x3+t]); + a8 = mul_float2(twiddles[k3+block_size], smem[x3+2*t]); + + a9 = smem[x4]; + a10 = mul_float2(twiddles[k4], smem[x4+t]); + a11 = mul_float2(twiddles[k4+block_size], smem[x4+2*t]); + } + + barrier(CLK_LOCAL_MEM_FENCE); + + if (x < t/4) + { + int dst_ind = ((x - k) * 3) + k; + + float2 b1 = a1 + a2; + a2 = twiddle(sin_120*(a1 - a2)); + float2 b0 = a0 - (float2)(0.5f)*b1; + + smem[dst_ind] = a0 + b1; + smem[dst_ind + block_size] = b0 + a2; + smem[dst_ind + 2*block_size] = b0 - a2; + + dst_ind = ((x2 - k2) * 3) + k2; + + b1 = a4 + a5; + a5 = twiddle(sin_120*(a4 - a5)); + b0 = a3 - (float2)(0.5f)*b1; + + smem[dst_ind] = a3 + b1; + smem[dst_ind + block_size] = b0 + a5; + smem[dst_ind + 2*block_size] = b0 - a5; + + dst_ind = ((x3 - k3) * 3) + k3; + + b1 = a7 + a8; + a8 = twiddle(sin_120*(a7 - a8)); + b0 = a6 - (float2)(0.5f)*b1; + + smem[dst_ind] = a6 + b1; + smem[dst_ind + block_size] = b0 + a8; + smem[dst_ind + 2*block_size] = b0 - a8; + + dst_ind = ((x4 - k4) * 3) + k4; + + b1 = a10 + a11; + a11 = twiddle(sin_120*(a10 - a11)); + b0 = a9 - (float2)(0.5f)*b1; + + smem[dst_ind] = a9 + b1; + smem[dst_ind + block_size] = b0 + a11; + smem[dst_ind + 2*block_size] = b0 - a11; + } + + barrier(CLK_LOCAL_MEM_FENCE); +} + __attribute__((always_inline)) void fft_radix5(__local float2* smem, __constant const float2* twiddles, const int x, const int block_size, const int t) { @@ -301,6 +482,95 @@ void fft_radix5(__local float2* smem, __constant const float2* twiddles, const i barrier(CLK_LOCAL_MEM_FENCE); } +__attribute__((always_inline)) +void fft_radix5_B2(__local float2* smem, __constant const float2* twiddles, const int x, const int block_size, const int t) +{ + const int k = x % block_size; + const int x2 = x+t/2; + const int k2 = x2 % block_size; + float2 a0, a1, a2, a3, a4, a5, a6, a7, a8, a9; + + if (x < t/2) + { + a0 = smem[x]; + a1 = mul_float2(twiddles[k], smem[x + t]); + a2 = mul_float2(twiddles[k + block_size],smem[x+2*t]); + a3 = mul_float2(twiddles[k+2*block_size],smem[x+3*t]); + a4 = mul_float2(twiddles[k+3*block_size],smem[x+4*t]); + + a5 = smem[x2]; + a6 = mul_float2(twiddles[k2], smem[x2 + t]); + a7 = mul_float2(twiddles[k2 + block_size],smem[x2+2*t]); + a8 = mul_float2(twiddles[k2+2*block_size],smem[x2+3*t]); + a9 = mul_float2(twiddles[k2+3*block_size],smem[x2+4*t]); + } + + barrier(CLK_LOCAL_MEM_FENCE); + + if (x < t/2) + { + int dst_ind = ((x - k) * 5) + k; + __local float2* dst = smem + dst_ind; + + float2 b0, b1, b5; + + b1 = a1 + a4; + a1 -= a4; + + a4 = a3 + a2; + a3 -= a2; + + a2 = b1 + a4; + b0 = a0 - (float2)0.25f * a2; + + b1 = fft5_2 * (b1 - a4); + a4 = fft5_3 * (float2)(-a1.y - a3.y, a1.x + a3.x); + b5 = (float2)(a4.x - fft5_5 * a1.y, a4.y + fft5_5 * a1.x); + + a4.x += fft5_4 * a3.y; + a4.y -= fft5_4 * a3.x; + + a1 = b0 + b1; + b0 -= b1; + + dst[0] = a0 + a2; + dst[block_size] = a1 + a4; + dst[2 * block_size] = b0 + b5; + dst[3 * block_size] = b0 - b5; + dst[4 * block_size] = a1 - a4; + + dst_ind = ((x2 - k2) * 5) + k2; + dst = smem + dst_ind; + + b1 = a6 + a9; + a6 -= a9; + + a9 = a8 + a7; + a8 -= a7; + + a7 = b1 + a9; + b0 = a5 - (float2)0.25f * a7; + + b1 = fft5_2 * (b1 - a9); + a9 = fft5_3 * (float2)(-a6.y - a8.y, a6.x + a8.x); + b5 = (float2)(a9.x - fft5_5 * a6.y, a9.y + fft5_5 * a6.x); + + a9.x += fft5_4 * a8.y; + a9.y -= fft5_4 * a8.x; + + a6 = b0 + b1; + b0 -= b1; + + dst[0] = a5 + a7; + dst[block_size] = a6 + a9; + dst[2 * block_size] = b0 + b5; + dst[3 * block_size] = b0 - b5; + dst[4 * block_size] = a6 - a9; + } + + barrier(CLK_LOCAL_MEM_FENCE); +} + #ifdef DFT_SCALE #define VAL(x, scale) x*scale #else diff --git a/modules/core/test/ocl/test_dft.cpp b/modules/core/test/ocl/test_dft.cpp index 64f6c63843..1fa03ae472 100644 --- a/modules/core/test/ocl/test_dft.cpp +++ b/modules/core/test/ocl/test_dft.cpp @@ -62,7 +62,7 @@ namespace ocl { //////////////////////////////////////////////////////////////////////////// // Dft -PARAM_TEST_CASE(Dft, cv::Size, OCL_FFT_TYPE, bool, bool, bool, bool) +PARAM_TEST_CASE(Dft, cv::Size, OCL_FFT_TYPE, bool, bool, bool) { cv::Size dft_size; int dft_flags, depth, cn, dft_type; @@ -91,9 +91,9 @@ PARAM_TEST_CASE(Dft, cv::Size, OCL_FFT_TYPE, bool, bool, bool, bool) dft_flags |= cv::DFT_ROWS; if (GET_PARAM(3)) dft_flags |= cv::DFT_SCALE; - if (GET_PARAM(4)) - dft_flags |= cv::DFT_INVERSE; - inplace = GET_PARAM(5); + /*if (GET_PARAM(4)) + dft_flags |= cv::DFT_INVERSE;*/ + inplace = GET_PARAM(4); is1d = (dft_flags & DFT_ROWS) != 0 || dft_size.height == 1; @@ -188,12 +188,12 @@ OCL_TEST_P(MulSpectrums, Mat) OCL_INSTANTIATE_TEST_CASE_P(OCL_ImgProc, MulSpectrums, testing::Combine(Bool(), Bool())); -OCL_INSTANTIATE_TEST_CASE_P(Core, Dft, Combine(Values(cv::Size(6, 4), cv::Size(5, 8), cv::Size(6, 6), +OCL_INSTANTIATE_TEST_CASE_P(Core, Dft, Combine(Values(cv::Size(16, 4), cv::Size(5, 8), cv::Size(6, 6), cv::Size(512, 1), cv::Size(1280, 768)), - Values(/*(OCL_FFT_TYPE) R2C, */(OCL_FFT_TYPE) C2C/*, (OCL_FFT_TYPE) R2R, (OCL_FFT_TYPE) C2R*/), + Values((OCL_FFT_TYPE) R2C, (OCL_FFT_TYPE) C2C, (OCL_FFT_TYPE) R2R, (OCL_FFT_TYPE) C2R), Bool(), // DFT_ROWS Bool(), // DFT_SCALE - Bool(), // DFT_INVERSE + //Bool(), // DFT_INVERSE Bool() // inplace ) );