modules/gpu/src/cuda/hough.cu

   1 /*M///////////////////////////////////////////////////////////////////////////////////////
   2 //
   3 //  IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
   4 //
   5 //  By downloading, copying, installing or using the software you agree to this license.
   6 //  If you do not agree to this license, do not download, install,
   7 //  copy or use the software.
   8 //
   9 //
  10 //                           License Agreement
  11 //                For Open Source Computer Vision Library
  12 //
  13 // Copyright (C) 2000-2008, Intel Corporation, all rights reserved.
  14 // Copyright (C) 2009, Willow Garage Inc., all rights reserved.
  15 // Third party copyrights are property of their respective owners.
  16 //
  17 // Redistribution and use in source and binary forms, with or without modification,
  18 // are permitted provided that the following conditions are met:
  19 //
  20 //   * Redistribution's of source code must retain the above copyright notice,
  21 //     this list of conditions and the following disclaimer.
  22 //
  23 //   * Redistribution's in binary form must reproduce the above copyright notice,
  24 //     this list of conditions and the following disclaimer in the documentation
  25 //     and/or other materials provided with the distribution.
  26 //
  27 //   * The name of the copyright holders may not be used to endorse or promote products
  28 //     derived from this software without specific prior written permission.
  29 //
  30 // This software is provided by the copyright holders and contributors "as is" and
  31 // any express or bpied warranties, including, but not limited to, the bpied
  32 // warranties of merchantability and fitness for a particular purpose are disclaimed.
  33 // In no event shall the Intel Corporation or contributors be liable for any direct,
  34 // indirect, incidental, special, exemplary, or consequential damages
  35 // (including, but not limited to, procurement of substitute goods or services;
  36 // loss of use, data, or profits; or business interruption) however caused
  37 // and on any theory of liability, whether in contract, strict liability,
  38 // or tort (including negligence or otherwise) arising in any way out of
  39 // the use of this software, even if advised of the possibility of such damage.
  40 //
  41 //M*/
  42
  43 #include <thrust/sort.h>
  44 #include "opencv2/gpu/device/common.hpp"
  45 #include "opencv2/gpu/device/emulation.hpp"
  46
  47 namespace cv { namespace gpu { namespace device
  48 {
  49     namespace hough
  50     {
  51         __device__ int g_counter;
  52
  53         ////////////////////////////////////////////////////////////////////////
  54         // buildPointList
  55
  56         const int PIXELS_PER_THREAD = 16;
  57
  58         __global__ void buildPointList(const DevMem2Db src, unsigned int* list)
  59         {
  60             __shared__ int s_queues[4][32 * PIXELS_PER_THREAD];
  61             __shared__ int s_qsize[4];
  62             __shared__ int s_globStart[4];
  63
  64             const int x = blockIdx.x * blockDim.x * PIXELS_PER_THREAD + threadIdx.x;
  65             const int y = blockIdx.y * blockDim.y + threadIdx.y;
  66
  67             if (y >= src.rows)
  68                 return;
  69
  70             if (threadIdx.x == 0)
  71                 s_qsize[threadIdx.y] = 0;
  72
  73             __syncthreads();
  74
  75             // fill the queue
  76             const uchar* srcRow = src.ptr(y);
  77             for (int i = 0, xx = x; i < PIXELS_PER_THREAD && xx < src.cols; ++i, xx += blockDim.x)
  78             {
  79                 if (srcRow[xx])
  80                 {
  81                     const unsigned int val = (y << 16) | xx;
  82                     const int qidx = Emulation::smem::atomicAdd(&s_qsize[threadIdx.y], 1);
  83                     s_queues[threadIdx.y][qidx] = val;
  84                 }
  85             }
  86
  87             __syncthreads();
  88
  89             // let one thread reserve the space required in the global list
  90             if (threadIdx.x == 0 && threadIdx.y == 0)
  91             {
  92                 // find how many items are stored in each list
  93                 int totalSize = 0;
  94                 for (int i = 0; i < blockDim.y; ++i)
  95                 {
  96                     s_globStart[i] = totalSize;
  97                     totalSize += s_qsize[i];
  98                 }
  99
 100                 // calculate the offset in the global list
 101                 const int globalOffset = atomicAdd(&g_counter, totalSize);
 102                 for (int i = 0; i < blockDim.y; ++i)
 103                     s_globStart[i] += globalOffset;
 104             }
 105
 106             __syncthreads();
 107
 108             // copy local queues to global queue
 109             const int qsize = s_qsize[threadIdx.y];
 110             int gidx = s_globStart[threadIdx.y] + threadIdx.x;
 111             for(int i = threadIdx.x; i < qsize; i += blockDim.x, gidx += blockDim.x)
 112                 list[gidx] = s_queues[threadIdx.y][i];
 113         }
 114
 115         int buildPointList_gpu(DevMem2Db src, unsigned int* list)
 116         {
 117             void* counterPtr;
 118             cudaSafeCall( cudaGetSymbolAddress(&counterPtr, g_counter) );
 119
 120             cudaSafeCall( cudaMemset(counterPtr, 0, sizeof(int)) );
 121
 122             const dim3 block(32, 4);
 123             const dim3 grid(divUp(src.cols, block.x * PIXELS_PER_THREAD), divUp(src.rows, block.y));
 124
 125             cudaSafeCall( cudaFuncSetCacheConfig(buildPointList, cudaFuncCachePreferShared) );
 126
 127             buildPointList<<<grid, block>>>(src, list);
 128             cudaSafeCall( cudaGetLastError() );
 129
 130             cudaSafeCall( cudaDeviceSynchronize() );
 131
 132             int totalCount;
 133             cudaSafeCall( cudaMemcpy(&totalCount, counterPtr, sizeof(int), cudaMemcpyDeviceToHost) );
 134
 135             return totalCount;
 136         }
 137
 138         ////////////////////////////////////////////////////////////////////////
 139         // linesAccum
 140
 141         __global__ void linesAccumGlobal(const unsigned int* list, const int count, PtrStepi accum, const float irho, const float theta, const int numrho)
 142         {
 143             const int n = blockIdx.x;
 144             const float ang = n * theta;
 145
 146             float sinVal;
 147             float cosVal;
 148             sincosf(ang, &sinVal, &cosVal);
 149             sinVal *= irho;
 150             cosVal *= irho;
 151
 152             const int shift = (numrho - 1) / 2;
 153
 154             int* accumRow = accum.ptr(n + 1);
 155             for (int i = threadIdx.x; i < count; i += blockDim.x)
 156             {
 157                 const unsigned int val = list[i];
 158
 159                 const int x = (val & 0xFFFF);
 160                 const int y = (val >> 16) & 0xFFFF;
 161
 162                 int r = __float2int_rn(x * cosVal + y * sinVal);
 163                 r += shift;
 164
 165                 ::atomicAdd(accumRow + r + 1, 1);
 166             }
 167         }
 168
 169         __global__ void linesAccumShared(const unsigned int* list, const int count, PtrStepi accum, const float irho, const float theta, const int numrho)
 170         {
 171             extern __shared__ int smem[];
 172
 173             for (int i = threadIdx.x; i < numrho + 1; i += blockDim.x)
 174                 smem[i] = 0;
 175
 176             __syncthreads();
 177
 178             const int n = blockIdx.x;
 179             const float ang = n * theta;
 180
 181             float sinVal;
 182             float cosVal;
 183             sincosf(ang, &sinVal, &cosVal);
 184             sinVal *= irho;
 185             cosVal *= irho;
 186
 187             const int shift = (numrho - 1) / 2;
 188
 189             for (int i = threadIdx.x; i < count; i += blockDim.x)
 190             {
 191                 const unsigned int val = list[i];
 192
 193                 const int x = (val & 0xFFFF);
 194                 const int y = (val >> 16) & 0xFFFF;
 195
 196                 int r = __float2int_rn(x * cosVal + y * sinVal);
 197                 r += shift;
 198
 199                 Emulation::smem::atomicAdd(&smem[r + 1], 1);
 200             }
 201
 202             __syncthreads();
 203
 204             int* accumRow = accum.ptr(n + 1);
 205             for (int i = threadIdx.x; i < numrho + 1; i += blockDim.x)
 206                 accumRow[i] = smem[i];
 207         }
 208
 209         void linesAccum_gpu(const unsigned int* list, int count, DevMem2Di accum, float rho, float theta, size_t sharedMemPerBlock, bool has20)
 210         {
 211             const dim3 block(has20 ? 1024 : 512);
 212             const dim3 grid(accum.rows - 2);
 213
 214             cudaSafeCall( cudaFuncSetCacheConfig(linesAccumShared, cudaFuncCachePreferShared) );
 215
 216             size_t smemSize = (accum.cols - 1) * sizeof(int);
 217
 218             if (smemSize < sharedMemPerBlock - 1000)
 219                 linesAccumShared<<<grid, block, smemSize>>>(list, count, accum, 1.0f / rho, theta, accum.cols - 2);
 220             else
 221                 linesAccumGlobal<<<grid, block>>>(list, count, accum, 1.0f / rho, theta, accum.cols - 2);
 222
 223             cudaSafeCall( cudaGetLastError() );
 224
 225             cudaSafeCall( cudaDeviceSynchronize() );
 226         }
 227
 228         ////////////////////////////////////////////////////////////////////////
 229         // linesGetResult
 230
 231         __global__ void linesGetResult(const DevMem2Di accum, float2* out, int* votes, const int maxSize, const float rho, const float theta, const float threshold, const int numrho)
 232         {
 233             __shared__ int smem[8][32];
 234
 235             const int x = blockIdx.x * (blockDim.x - 2) + threadIdx.x;
 236             const int y = blockIdx.y * (blockDim.y - 2) + threadIdx.y;
 237
 238             if (x >= accum.cols || y >= accum.rows)
 239                 return;
 240
 241             smem[threadIdx.y][threadIdx.x] = accum(y, x);
 242             __syncthreads();
 243
 244             const int r = x - 1;
 245             const int n = y - 1;
 246
 247             if (threadIdx.x == 0 || threadIdx.x == blockDim.x - 1 || threadIdx.y == 0 || threadIdx.y == blockDim.y - 1 || r >= accum.cols - 2 || n >= accum.rows - 2)
 248                 return;
 249
 250             if (smem[threadIdx.y][threadIdx.x] > threshold &&
 251                 smem[threadIdx.y][threadIdx.x] >  smem[threadIdx.y - 1][threadIdx.x] &&
 252                 smem[threadIdx.y][threadIdx.x] >= smem[threadIdx.y + 1][threadIdx.x] &&
 253                 smem[threadIdx.y][threadIdx.x] >  smem[threadIdx.y][threadIdx.x - 1] &&
 254                 smem[threadIdx.y][threadIdx.x] >= smem[threadIdx.y][threadIdx.x + 1])
 255             {
 256                 const float radius = (r - (numrho - 1) * 0.5f) * rho;
 257                 const float angle = n * theta;
 258
 259                 const int ind = ::atomicAdd(&g_counter, 1);
 260                 if (ind < maxSize)
 261                 {
 262                     out[ind] = make_float2(radius, angle);
 263                     votes[ind] = smem[threadIdx.y][threadIdx.x];
 264                 }
 265             }
 266         }
 267
 268         int linesGetResult_gpu(DevMem2Di accum, float2* out, int* votes, int maxSize, float rho, float theta, float threshold, bool doSort)
 269         {
 270             void* counterPtr;
 271             cudaSafeCall( cudaGetSymbolAddress(&counterPtr, g_counter) );
 272
 273             cudaSafeCall( cudaMemset(counterPtr, 0, sizeof(int)) );
 274
 275             const dim3 block(32, 8);
 276             const dim3 grid(divUp(accum.cols, block.x - 2), divUp(accum.rows, block.y - 2));
 277
 278             linesGetResult<<<grid, block>>>(accum, out, votes, maxSize, rho, theta, threshold, accum.cols - 2);
 279             cudaSafeCall( cudaGetLastError() );
 280
 281             cudaSafeCall( cudaDeviceSynchronize() );
 282
 283             int totalCount;
 284             cudaSafeCall( cudaMemcpy(&totalCount, counterPtr, sizeof(int), cudaMemcpyDeviceToHost) );
 285
 286             totalCount = ::min(totalCount, maxSize);
 287
 288             if (doSort && totalCount > 0)
 289             {
 290                 thrust::device_ptr<float2> outPtr(out);
 291                 thrust::device_ptr<int> votesPtr(votes);
 292                 thrust::sort_by_key(votesPtr, votesPtr + totalCount, outPtr, thrust::greater<int>());
 293             }
 294
 295             return totalCount;
 296         }
 297     }
 298 }}}