added mulSpectrums functions into GPU module
authorAlexey Spizhevoy <no@email>
Wed, 22 Dec 2010 13:46:06 +0000 (13:46 +0000)
committerAlexey Spizhevoy <no@email>
Wed, 22 Dec 2010 13:46:06 +0000 (13:46 +0000)
modules/gpu/include/opencv2/gpu/gpu.hpp
modules/gpu/src/cuda/imgproc.cu
modules/gpu/src/imgproc_gpu.cpp
modules/gpu/src/match_template.cpp

index 48d1c35..71a053b 100644 (file)
@@ -628,10 +628,19 @@ namespace cv
         //! computes minimum eigen value of 2x2 derivative covariation matrix at each pixel - the cornerness criteria\r
         CV_EXPORTS void cornerMinEigenVal(const GpuMat& src, GpuMat& dst, int blockSize, int ksize, int borderType=BORDER_REFLECT101);\r
 \r
-        //! computes cross-correlation of two images using FFT\r
+        //! performs per-element multiplication of two full (i.e. not packed) Fourier spectrums\r
+        //! supports only 32FC2 matrixes (interleaved format)\r
+        CV_EXPORTS void mulSpectrums(const GpuMat& a, const GpuMat& b, GpuMat& c, int flags, bool conjB=false);\r
+\r
+        //! performs per-element multiplication of two full (i.e. not packed) Fourier spectrums\r
+        //! supports only 32FC2 matrixes (interleaved format)\r
+        CV_EXPORTS void mulAndScaleSpectrums(const GpuMat& a, const GpuMat& b, GpuMat& c, int flags, \r
+                                             float scale, bool conjB=false);\r
+\r
+        //! computes convolution (or cross-correlation) of two images using discrete Fourier transform\r
         //! supports source images of 32FC1 type only\r
         //! result matrix will have 32FC1 type\r
-        CV_EXPORTS void crossCorr(const GpuMat& image, const GpuMat& templ, GpuMat& result);\r
+        CV_EXPORTS void convolve(const GpuMat& image, const GpuMat& templ, GpuMat& result, bool ccorr=false);\r
 \r
         //! computes the proximity map for the raster template and the image where the template is searched for\r
         CV_EXPORTS void matchTemplate(const GpuMat& image, const GpuMat& templ, GpuMat& result, int method);\r
index f8abe37..7b66565 100644 (file)
@@ -40,7 +40,6 @@
 //\r
 //M*/\r
 \r
-#include <cufft.h>\r
 #include "internal_shared.hpp"\r
 #include "opencv2/gpu/device/border_interpolate.hpp"\r
 \r
@@ -751,31 +750,121 @@ namespace cv { namespace gpu { namespace imgproc
     }\r
 \r
     //////////////////////////////////////////////////////////////////////////\r
-    // multiplyAndNormalizeSpects\r
+    // mulSpectrums\r
 \r
-    __global__ void multiplyAndNormalizeSpectsKernel(\r
-            int n, float scale, const cufftComplex* a, \r
-            const cufftComplex* b, cufftComplex* c)\r
+\r
+    __global__ void mulSpectrumsKernel(const PtrStep_<cufftComplex> a, const PtrStep_<cufftComplex> b, \r
+                                       DevMem2D_<cufftComplex> c)\r
+    {\r
+        const int x = blockIdx.x * blockDim.x + threadIdx.x;    \r
+        const int y = blockIdx.y * blockDim.y + threadIdx.y;    \r
+\r
+        if (x < c.cols && y < c.rows) \r
+        {\r
+            c.ptr(y)[x] = cuCmulf(a.ptr(y)[x], b.ptr(y)[x]);\r
+        }\r
+    }\r
+\r
+\r
+    void mulSpectrums(const PtrStep_<cufftComplex> a, const PtrStep_<cufftComplex> b, \r
+                      DevMem2D_<cufftComplex> c)\r
+    {\r
+        dim3 threads(256);\r
+        dim3 grid(divUp(c.cols, threads.x), divUp(c.rows, threads.y));\r
+\r
+        mulSpectrumsKernel<<<grid, threads>>>(a, b, c);\r
+        cudaSafeCall(cudaThreadSynchronize());\r
+    }\r
+\r
+\r
+    //////////////////////////////////////////////////////////////////////////\r
+    // mulSpectrums_CONJ\r
+\r
+\r
+    __global__ void mulSpectrumsKernel_CONJ(\r
+            const PtrStep_<cufftComplex> a, const PtrStep_<cufftComplex> b,\r
+            DevMem2D_<cufftComplex> c)\r
+    {\r
+        const int x = blockIdx.x * blockDim.x + threadIdx.x;    \r
+        const int y = blockIdx.y * blockDim.y + threadIdx.y;    \r
+\r
+        if (x < c.cols && y < c.rows) \r
+        {\r
+            c.ptr(y)[x] = cuCmulf(a.ptr(y)[x], cuConjf(b.ptr(y)[x]));\r
+        }\r
+    }\r
+\r
+\r
+    void mulSpectrums_CONJ(const PtrStep_<cufftComplex> a, const PtrStep_<cufftComplex> b, \r
+                           DevMem2D_<cufftComplex> c)\r
+    {\r
+        dim3 threads(256);\r
+        dim3 grid(divUp(c.cols, threads.x), divUp(c.rows, threads.y));\r
+\r
+        mulSpectrumsKernel_CONJ<<<grid, threads>>>(a, b, c);\r
+        cudaSafeCall(cudaThreadSynchronize());\r
+    }\r
+\r
+\r
+    //////////////////////////////////////////////////////////////////////////\r
+    // mulAndScaleSpectrums\r
+\r
+\r
+    __global__ void mulAndScaleSpectrumsKernel(\r
+            const PtrStep_<cufftComplex> a, const PtrStep_<cufftComplex> b, \r
+            float scale, DevMem2D_<cufftComplex> c)\r
     {\r
-        int x = blockIdx.x * blockDim.x + threadIdx.x;    \r
-        if (x < n) \r
+        const int x = blockIdx.x * blockDim.x + threadIdx.x;\r
+        const int y = blockIdx.y * blockDim.y + threadIdx.y;\r
+\r
+        if (x < c.cols && y < c.rows) \r
         {\r
-            cufftComplex v = cuCmulf(a[x], cuConjf(b[x]));\r
-            c[x] = make_cuFloatComplex(cuCrealf(v) * scale, cuCimagf(v) * scale);\r
+            cufftComplex v = cuCmulf(a.ptr(y)[x], b.ptr(y)[x]);\r
+            c.ptr(y)[x] = make_cuFloatComplex(cuCrealf(v) * scale, cuCimagf(v) * scale);\r
         }\r
     }\r
 \r
 \r
-    // Performs per-element multiplication and normalization of two spectrums\r
-    void multiplyAndNormalizeSpects(int n, float scale, const cufftComplex* a, \r
-                                    const cufftComplex* b, cufftComplex* c)\r
+    void mulAndScaleSpectrums(const PtrStep_<cufftComplex> a, const PtrStep_<cufftComplex> b,\r
+                              float scale, DevMem2D_<cufftComplex> c)\r
     {\r
         dim3 threads(256);\r
-        dim3 grid(divUp(n, threads.x));\r
+        dim3 grid(divUp(c.cols, threads.x), divUp(c.rows, threads.y));\r
 \r
-        multiplyAndNormalizeSpectsKernel<<<grid, threads>>>(n, scale, a, b, c);\r
+        mulAndScaleSpectrumsKernel<<<grid, threads>>>(a, b, scale, c);\r
         cudaSafeCall(cudaThreadSynchronize());\r
     }\r
 \r
+\r
+    //////////////////////////////////////////////////////////////////////////\r
+    // mulAndScaleSpectrums_CONJ\r
+\r
+\r
+    __global__ void mulAndScaleSpectrumsKernel_CONJ(\r
+            const PtrStep_<cufftComplex> a, const PtrStep_<cufftComplex> b,\r
+            float scale, DevMem2D_<cufftComplex> c)\r
+    {\r
+        const int x = blockIdx.x * blockDim.x + threadIdx.x;\r
+        const int y = blockIdx.y * blockDim.y + threadIdx.y;\r
+\r
+        if (x < c.cols && y < c.rows) \r
+        {\r
+            cufftComplex v = cuCmulf(a.ptr(y)[x], cuConjf(b.ptr(y)[x]));\r
+            c.ptr(y)[x] = make_cuFloatComplex(cuCrealf(v) * scale, cuCimagf(v) * scale);\r
+        }\r
+    }\r
+\r
+\r
+    void mulAndScaleSpectrums_CONJ(const PtrStep_<cufftComplex> a, const PtrStep_<cufftComplex> b,\r
+                                  float scale, DevMem2D_<cufftComplex> c)\r
+    {\r
+        dim3 threads(256);\r
+        dim3 grid(divUp(c.cols, threads.x), divUp(c.rows, threads.y));\r
+\r
+        mulAndScaleSpectrumsKernel_CONJ<<<grid, threads>>>(a, b, scale, c);\r
+        cudaSafeCall(cudaThreadSynchronize());\r
+    }\r
+\r
+\r
 }}}\r
 \r
index 4a3f9de..4aaef14 100644 (file)
@@ -74,7 +74,9 @@ void cv::gpu::histRange(const GpuMat&, GpuMat&, const GpuMat&) { throw_nogpu();
 void cv::gpu::histRange(const GpuMat&, GpuMat*, const GpuMat*) { throw_nogpu(); }\r
 void cv::gpu::cornerHarris(const GpuMat&, GpuMat&, int, int, double, int) { throw_nogpu(); }\r
 void cv::gpu::cornerMinEigenVal(const GpuMat&, GpuMat&, int, int, int) { throw_nogpu(); }\r
-void cv::gpu::crossCorr(const GpuMat&, const GpuMat&, GpuMat&) { throw_nogpu(); }\r
+void cv::gpu::mulSpectrums(const GpuMat&, const GpuMat&, GpuMat&, int, bool) { throw_nogpu(); }\r
+void cv::gpu::mulAndScaleSpectrums(const GpuMat&, const GpuMat&, GpuMat&, int, float, bool) { throw_nogpu(); }\r
+void cv::gpu::convolve(const GpuMat&, const GpuMat&, GpuMat&, bool) { throw_nogpu(); }\r
 \r
 \r
 #else /* !defined (HAVE_CUDA) */\r
@@ -1065,6 +1067,66 @@ void cv::gpu::cornerMinEigenVal(const GpuMat& src, GpuMat& dst, int blockSize, i
 }\r
 \r
 //////////////////////////////////////////////////////////////////////////////\r
+// mulSpectrums\r
+\r
+namespace cv { namespace gpu { namespace imgproc \r
+{\r
+    void mulSpectrums(const PtrStep_<cufftComplex> a, const PtrStep_<cufftComplex> b, \r
+                      DevMem2D_<cufftComplex> c);\r
+\r
+    void mulSpectrums_CONJ(const PtrStep_<cufftComplex> a, const PtrStep_<cufftComplex> b, \r
+                           DevMem2D_<cufftComplex> c);\r
+}}}\r
+\r
+\r
+void cv::gpu::mulSpectrums(const GpuMat& a, const GpuMat& b, GpuMat& c, \r
+                           int flags, bool conjB) \r
+{\r
+    typedef void (*Caller)(const PtrStep_<cufftComplex>, const PtrStep_<cufftComplex>, \r
+                           DevMem2D_<cufftComplex>);\r
+    static Caller callers[] = { imgproc::mulSpectrums, \r
+                                imgproc::mulSpectrums_CONJ };\r
+\r
+    CV_Assert(a.type() == b.type() && a.type() == CV_32FC2);\r
+    CV_Assert(a.size() == b.size());\r
+\r
+    c.create(a.size(), CV_32FC2);\r
+\r
+    Caller caller = callers[(int)conjB];\r
+    caller(a, b, c);\r
+}\r
+\r
+//////////////////////////////////////////////////////////////////////////////\r
+// mulAndScaleSpectrums\r
+\r
+namespace cv { namespace gpu { namespace imgproc \r
+{\r
+    void mulAndScaleSpectrums(const PtrStep_<cufftComplex> a, const PtrStep_<cufftComplex> b,\r
+                             float scale, DevMem2D_<cufftComplex> c);\r
+\r
+    void mulAndScaleSpectrums_CONJ(const PtrStep_<cufftComplex> a, const PtrStep_<cufftComplex> b,\r
+                                  float scale, DevMem2D_<cufftComplex> c);\r
+}}}\r
+\r
+\r
+void cv::gpu::mulAndScaleSpectrums(const GpuMat& a, const GpuMat& b, GpuMat& c,\r
+                                  int flags, float scale, bool conjB) \r
+{\r
+    typedef void (*Caller)(const PtrStep_<cufftComplex>, const PtrStep_<cufftComplex>,\r
+                           float scale, DevMem2D_<cufftComplex>);\r
+    static Caller callers[] = { imgproc::mulAndScaleSpectrums, \r
+                                imgproc::mulAndScaleSpectrums_CONJ };\r
+\r
+    CV_Assert(a.type() == b.type() && a.type() == CV_32FC2);\r
+    CV_Assert(a.size() == b.size());\r
+\r
+    c.create(a.size(), CV_32FC2);\r
+\r
+    Caller caller = callers[(int)conjB];\r
+    caller(a, b, scale, c);\r
+}\r
+\r
+//////////////////////////////////////////////////////////////////////////////\r
 // crossCorr\r
 \r
 namespace \r
@@ -1094,15 +1156,12 @@ namespace
 }\r
 \r
 \r
-namespace cv { namespace gpu { namespace imgproc\r
+void cv::gpu::convolve(const GpuMat& image, const GpuMat& templ, GpuMat& result, bool ccorr)\r
 {\r
-    void multiplyAndNormalizeSpects(int n, float scale, const cufftComplex* a,\r
-                                    const cufftComplex* b, cufftComplex* c);\r
-}}}\r
+    // We must be sure we use correct OpenCV analogues for CUFFT types\r
+    StaticAssert<sizeof(float) == sizeof(cufftReal)>::check();\r
+    StaticAssert<sizeof(float) * 2 == sizeof(cufftComplex)>::check();\r
 \r
-\r
-void cv::gpu::crossCorr(const GpuMat& image, const GpuMat& templ, GpuMat& result)\r
-{\r
     CV_Assert(image.type() == CV_32F);\r
     CV_Assert(templ.type() == CV_32F);\r
 \r
@@ -1119,33 +1178,28 @@ void cv::gpu::crossCorr(const GpuMat& image, const GpuMat& templ, GpuMat& result
     block_size.width = std::min(dft_size.width - templ.cols + 1, result.cols);\r
     block_size.height = std::min(dft_size.height - templ.rows + 1, result.rows);\r
 \r
-    cufftReal* image_data;\r
-    cufftReal* templ_data;\r
-    cufftReal* result_data;\r
-    cudaSafeCall(cudaMalloc((void**)&image_data, sizeof(cufftReal) * dft_size.area()));\r
-    cudaSafeCall(cudaMalloc((void**)&templ_data, sizeof(cufftReal) * dft_size.area()));\r
-    cudaSafeCall(cudaMalloc((void**)&result_data, sizeof(cufftReal) * dft_size.area()));\r
+    GpuMat image_data(1, dft_size.area(), CV_32F);\r
+    GpuMat templ_data(1, dft_size.area(), CV_32F);\r
+    GpuMat result_data(1, dft_size.area(), CV_32F);\r
 \r
     int spect_len = dft_size.height * (dft_size.width / 2 + 1);\r
-    cufftComplex* image_spect;\r
-    cufftComplex* templ_spect;\r
-    cufftComplex* result_spect;\r
-    cudaSafeCall(cudaMalloc((void**)&image_spect, sizeof(cufftComplex) * spect_len));\r
-    cudaSafeCall(cudaMalloc((void**)&templ_spect, sizeof(cufftComplex) * spect_len));\r
-    cudaSafeCall(cudaMalloc((void**)&result_spect, sizeof(cufftComplex) * spect_len));\r
+    GpuMat image_spect(1, spect_len, CV_32FC2);\r
+    GpuMat templ_spect(1, spect_len, CV_32FC2);\r
+    GpuMat result_spect(1, spect_len, CV_32FC2);\r
 \r
     cufftHandle planR2C, planC2R;\r
     cufftSafeCall(cufftPlan2d(&planC2R, dft_size.height, dft_size.width, CUFFT_C2R));\r
     cufftSafeCall(cufftPlan2d(&planR2C, dft_size.height, dft_size.width, CUFFT_R2C));\r
 \r
-    GpuMat templ_roi(templ.size(), CV_32S, templ.data, templ.step);\r
-    GpuMat templ_block(dft_size, CV_32S, templ_data, dft_size.width * sizeof(cufftReal));\r
+    GpuMat templ_roi(templ.size(), CV_32F, templ.data, templ.step);\r
+    GpuMat templ_block(dft_size, CV_32F, templ_data.ptr(), dft_size.width * sizeof(cufftReal));\r
     copyMakeBorder(templ_roi, templ_block, 0, templ_block.rows - templ_roi.rows, 0, \r
                    templ_block.cols - templ_roi.cols, 0);\r
 \r
-    cufftSafeCall(cufftExecR2C(planR2C, templ_data, templ_spect));\r
+    cufftSafeCall(cufftExecR2C(planR2C, templ_data.ptr<cufftReal>(), \r
+                               templ_spect.ptr<cufftComplex>()));\r
 \r
-    GpuMat image_block(dft_size, CV_32S, image_data, dft_size.width * sizeof(cufftReal));\r
+    GpuMat image_block(dft_size, CV_32F, image_data.ptr(), dft_size.width * sizeof(cufftReal));\r
 \r
     // Process all blocks of the result matrix\r
     for (int y = 0; y < result.rows; y += block_size.height)\r
@@ -1156,18 +1210,20 @@ void cv::gpu::crossCorr(const GpuMat& image, const GpuMat& templ, GpuMat& result
             Size image_roi_size;\r
             image_roi_size.width = std::min(x + dft_size.width, image.cols) - x;\r
             image_roi_size.height = std::min(y + dft_size.height, image.rows) - y;\r
-            GpuMat image_roi(image_roi_size, CV_32S, (void*)(image.ptr<float>(y) + x), image.step);\r
+            GpuMat image_roi(image_roi_size, CV_32F, (void*)(image.ptr<float>(y) + x), image.step);\r
 \r
             // Make source image block continous\r
             copyMakeBorder(image_roi, image_block, 0, image_block.rows - image_roi.rows, 0, \r
                            image_block.cols - image_roi.cols, 0);\r
 \r
-            cufftSafeCall(cufftExecR2C(planR2C, image_data, image_spect));\r
+            cufftSafeCall(cufftExecR2C(planR2C, image_data.ptr<cufftReal>(), \r
+                                       image_spect.ptr<cufftComplex>()));\r
 \r
-            imgproc::multiplyAndNormalizeSpects(spect_len, 1.f / dft_size.area(), \r
-                                                image_spect, templ_spect, result_spect);\r
+            mulAndScaleSpectrums(image_spect, templ_spect, result_spect, 0,\r
+                                1.f / dft_size.area(), ccorr);\r
 \r
-            cufftSafeCall(cufftExecC2R(planC2R, result_spect, result_data));\r
+            cufftSafeCall(cufftExecC2R(planC2R, result_spect.ptr<cufftComplex>(), \r
+                                       result_data.ptr<cufftReal>()));\r
 \r
             // Copy result block into appropriate part of the result matrix.\r
             // We can't compute it inplace as the result of the CUFFT transforms\r
@@ -1176,23 +1232,17 @@ void cv::gpu::crossCorr(const GpuMat& image, const GpuMat& templ, GpuMat& result
             result_roi_size.width = std::min(x + block_size.width, result.cols) - x;\r
             result_roi_size.height = std::min(y + block_size.height, result.rows) - y;\r
             GpuMat result_roi(result_roi_size, CV_32F, (void*)(result.ptr<float>(y) + x), result.step);\r
-            GpuMat result_block(result_roi_size, CV_32F, result_data, dft_size.width * sizeof(cufftReal));\r
+            GpuMat result_block(result_roi_size, CV_32F, result_data.ptr(), dft_size.width * sizeof(cufftReal));\r
             result_block.copyTo(result_roi);\r
         }\r
     }\r
 \r
     cufftSafeCall(cufftDestroy(planR2C));\r
     cufftSafeCall(cufftDestroy(planC2R));\r
-\r
-    cudaSafeCall(cudaFree(image_spect));\r
-    cudaSafeCall(cudaFree(templ_spect));\r
-    cudaSafeCall(cudaFree(result_spect));\r
-    cudaSafeCall(cudaFree(image_data));\r
-    cudaSafeCall(cudaFree(templ_data));\r
-    cudaSafeCall(cudaFree(result_data));\r
 }\r
 \r
 \r
 \r
 #endif /* !defined (HAVE_CUDA) */\r
 \r
+\r
index 1106ca9..5f5807c 100644 (file)
@@ -196,7 +196,7 @@ namespace
         }\r
 \r
         GpuMat result_;\r
-        crossCorr(image.reshape(1), templ.reshape(1), result_);\r
+        convolve(image.reshape(1), templ.reshape(1), result_, true);\r
         imgproc::extractFirstChannel_32F(result_, result, image.channels());\r
     }\r
 \r