CUDA Learning Material

cuda learning notes

Introduction

Recently I read the book “Hands-On GPU-Accelerated Computer Vision with OpenCV and CUDA”, whose code is hosted on GitHub. It introduces common concepts in NVIDIA CUDA, such as grid, block, thread, global memory, shared memory, constant memory, texture memory, CUDA event, and CUDA stream, along with their typical usage. It also covers the use of common CUDA-version functions in OpenCV, including basic image preprocessing, thresholding, filtering, feature-point detection and matching, and color- and shape-based object detection & tracking. It further describes the deployment of the relevant CUDA algorithms on embedded devices such as the NVIDIA Jetson TX1. Overall, it is a good and relatively easy-to-read book.

Example

Below we introduce several simple examples:

  • GPU Addition
__global__ void add_kernel(int* da, int* db, int* dc, int N) {
    int tidx = threadIdx.x + blockIdx.x * blockDim.x;
    while(tidx <  N) {
        dc[tidx] = da[tidx] + db[tidx];
        tidx += gridDim.x * blockDim.x;
    }
}

  • GPU Rank Sort
__global__ void ranksortKernel(int* da, int* db) {
    int count = 0;
    int tid = threadIdx.x;
    int ttid = blockIdx.x * threadPerBlock + tid;
    int val = da[ttid];
    __shared__ int cache[threadPerBlock];
    for(int i = tid; i< arraySize; i += threadPerBlock) {
        cache[tid] = da[i]; //the block corresponding to each thread caches a portion of the data, and caches another portion after computing
        //the two blocks execute the first loop in parallel and then the second loop in parallel, each caching the same data simultaneously
        __syncthreads();
        for(int j = 0; j < threadPerBlock;j++) {
            if(val > cache[j])
              count++;
        }
    }
    db[count] =val;
}

  • GPU Matrix Multiply
__global__ void gpu_matrix_mul_kernel(float* da, float* db, float* dc, const int size) {
    int row, col;
    col = threadIdx.x + blockIdx.x * TileSize;
    row = threadIdx.y + blockIdx.y * TileSize;
    //the block to which each corresponding thread belongs contains shared_a and shared_b
    __shared__ float shared_a[TileSize][TileSize];
    __shared__ float shared_b[TileSize][TileSize];
    for(int i = 0; i < size / TileSize; i++) {
      shared_a[threadIdx.y][threadIdx.x] = da[row * size + (i * TileSize + threadIdx.x)];
      shared_b[threadIdx.y][threadIdx.x] = db[col + (i * TileSize + threadIdx.y) * size];
      __syncthreads(); //use all threads in the block to fill the block
      for(int j = 0; j < TileSize; j++) {
          dc[row * size + col] += shared_a[threadIdx.y][j] * shared_b[j][threadIdx.x];
      }
      __syncthreads();
    }
}

A typical way to call it is as follows:

void test() {
    const int size = 4;
    float ha[size][size], hb[size][size], h_result[size][size];
    float* da, * db, * d_result;
    for(int i =0; i < size;i++) {
        for(int j = 0; j < size; j++) {
            ha[i][j] = i;
            hb[i][j] = j;
        }
    }
    cudaMalloc((void**)&da, sizeof(int) * size * size); //remember to add the address-of operator &
    cudaMalloc((void**)&db, sizeof(int) * size * size);
    cudaMalloc((void**)&d_result, sizeof(int) * size * size);
    cudaMemcpy(da, ha, sizeof(int) * size * size, cudaMemcpyHostToDevice);
    cudaMemcpy(db, hb, sizeof(int) * size * size, cudaMemcpyHostToDevice);

    gpu_matrix_mul(da, db, d_result,size); //the kernel function is invoked inside this function
    cudaMemcpy(h_result, d_result, sizeof(int) * size * size, cudaMemcpyDeviceToHost);
    for(int i = 0; i < size; i++) {
        for(int j = 0; j < size; j++) {
            std::cout<<h_result[i][j] <<" ";
        }
        std::cout<< std::endl;
    }
    cudaFree(da);
    cudaFree(db);
    cudaFree(d_result);
}

  • GPU RGB Image To Gray
__global__ void gray(uchar4* d_in, uchar* d_out, int rows, int cols) {
    int blockid = blockIdx.y * gridDim.x + blockIdx.x;
    int threadid = blockDim.x * blockDim.y * blockid + threadIdx.y * blockDim.x + threadIdx.x;

    if(threadid < rows * cols) {
        uchar4 rgba = d_in[threadid];
        d_out[threadid] = 0.299f*rgba.x+0.587f*rgba.y+0.114f*rgba.z;
    }
}

cv::Mat cudafunc(cv::Mat rgba_image) {
    int rows = rgba_image.rows;
    int cols = rgba_image.cols;

    uchar4* h_in = (uchar4*)rgba_image.data;
    uchar* h_out = (uchar*)malloc(sizeof(uchar) * rows * cols);

    uchar4* d_in;
    uchar* d_out;
    cudaMalloc((void**)&d_in, sizeof(uchar4) * rows * cols);
    cudaMalloc((void**)&d_out, sizeof(uchar) * rows * cols);

    // copy data from host to device
    cudaMemcpy(d_in, h_in, sizeof(uchar4) * rows * cols, cudaMemcpyHostToDevice);

    // call kernel function process image
    const dim3 blocksize(32, 32, 1);
    const dim3 gridsize((rows - 1 + 32) / 32, (cols - 1 + 32) / 32, 1);
    gray <<<gridsize,blocksize>>>(d_in, d_out, rows, cols);
    cudaMemcpy(h_out, d_out, sizeof(uchar) * rows * cols, cudaMemcpyDeviceToHost);
    
    cudaFree(d_in);
    cudaFree(d_out);

    cv::Mat out(rows,cols, CV_8UC1, h_out);
    return out;
}

  • GPU HIST CAL
__global__ void histo_kernel(unsigned char* buffer,long size,unsigned int * histo) {
    __shared__ unsigned int temp[256];
    temp[threadIdx.x] = 0;
    __syncthreads();

    int i = threadIdx.x + blockIdx.x * blockDim.x;
    int offset = blockDim.x * gridDim.x;
    //the block corresponding to each thread buffers the pixel-count array
    while(i < size){
        atomicAdd(&temp[buffer[i]],1);
        i += offset;
    }
    __syncthreads();
    //accumulate the buffered statistics from all blocks
    atomicAdd(&(histo[threadIdx.x]),temp[threadIdx.x]);
}

  • GPU MEAN FILTER

#define TILE_W 16
#define TILE_H 16
#define R 2 // filter radius
#define BLOCK_W (TILE_W + (2 * R))
#define BLOCK_H (TILE_H + (2 * R))

__global__ void filter(unsigned char* d_in, unsigned char* d_out, int width, int height) {
  __shared__ unsigned char ssem[BLOCK_W * BLOCK_H];
  // compute the index in the image from the thread index
  int x = blockIdx.x * TILE_W + threadIdx.x - R;
  int y = blockIdx.y * TILE_H + threadIdx.y - R;

  x = max(0, x);
  x = min(x, width - 1);
  y = max(0, y);
  y = min(y, height - 1);

  int index = y * width + x;
  int bindex = threadIdx.y * blockDim.x + threadIdx.x;

  ssem[bindex] = int(d_in[index]);
  __syncthreads();

  if(threadIdx.x >= R && threadIdx.x < (BLOCK_W - R) && 
     threadIdx.y >= R && threadIdx.y < (BLOCK_H - R)) {
      float sum = 0;
      for(int dy = -R; dy <= R; dy++) {
          for(int dx = -R; dx <= R; dx++) {
              sum += int(ssem[bindex + dy * blockDim.x + dx]);
          }
      }
      d_out[index] = int(sum / S);
  }
}

Comments