CUDA によるヒストグラム生成高速化

自分のロボット部品のCUDA化する必要があったので,CUDAの勉強もかねてGPUを使ったヒストグラム作成のプログラムを作りました.
内容としては,下記の本の9章のサンプルコードを参考に作成しました.下記の本,基本的なことが一通り説明されてまして,個人的には結構おすすめでした.

CUDA by Example 汎用GPUプログラミング入門

CUDA by Example 汎用GPUプログラミング入門

ヒストグラムのもととなる写真は,今をときめくSpace XのFalcon9です.
f:id:rkoichi2001:20170910155552j:plain

Falcon9の垂直着陸シーン.(これ,ほんとすごいですよね!酒飲んだら毎回必ず興奮して,着陸シーンを友達に見せつけてしまいます.自分がなにかやったわけでは無いですが...)
www.youtube.com

で,だいぶ派手な動画の後にだいぶ地味なエントリ内容で恐縮ですが,GPUとCPUでヒストグラムを計算した実行時間の結果ですが,下記のようになりました.

実験条件:
添付写真 (1500pix x 1000pix,グレースケール)のヒストグラム計算
実行回数 10000回のヒストグラム計算 x 10Set
実行環境 CPU : Intel Core i7-4700HQ, GPU : NVIDIA GeForce GT 740M (64bit)
実行結果 

f:id:rkoichi2001:20170910162642p:plain

うーん..もう少し早くなると思ったんですが,半分くらいの速度にしかならなかったです.
おそらく,GPUのスペックがだいぶ低いからだと思うのですが,Jetsonだともって早くなってくれるはず...

サンプルコード(一応,こちらがホスト側のコードです.)
//
#include <iostream>
#include <cuda_runtime.h>
#include <cstdlib>
#include <sys/time.h>
#include <stdio.h>
#include <assert.h>

//
#include "opencv2/highgui.hpp"
#include "opencv2/core.hpp"
#include "opencv2/imgproc.hpp"

//


#define REPEAT_NUM (10000)

static void createHistogramWithCPU(int width, int height, const unsigned char* buff, unsigned int *histo);
extern void createHistogramWithGPU(int width, int height, const unsigned char* buff, unsigned int *histo);

static void createHistogramWithCPU(int width, int height, const unsigned char* buff, unsigned int *histo) {

  const unsigned char *pnt = buff;

  for (int j = 0; j < height; j++) {
    for (int i = 0; i < width; i++) {
      histo[*pnt]++;
      pnt++;
    }
  }
}

static void drawHistogram(const unsigned int* histo, cv::Mat &histoImg) {

  int max = 0;
  {
    const unsigned int *pnt = histo;
    for (int i = 0; i < 256; i++) {
      unsigned int value = *pnt;
      pnt++;
      if (max < value) {
        max = value;
      }
    }
  }

  {
    const unsigned int *pnt = histo;
    histoImg = cv::Mat(512, 512, CV_8UC1);
    histoImg = 255;
    for (int i = 0; i < 512; i+=2) {
      unsigned int value = *pnt;
      pnt++;
      unsigned int thresh = 512 - (512 / (float)max * value);
      for (int j = 0; j < 512; j++) {
        if (j > thresh) {
          histoImg.at<unsigned char>(j, i) = 0;
          histoImg.at<unsigned char>(j, i + 1) = 0;
        }
      }
    }
  }

}

int main(int argc, char** argv)
{
  std::cout << "Cuda Sample Started!" << std::endl;

  unsigned int *histoCpu = new unsigned int[256];
  unsigned int *histoGpu = new unsigned int[256];


  cv::Mat src = cv::imread("./falcon_9.jpg", cv::IMREAD_GRAYSCALE);
  cv::resize(src, src, cv::Size(), 0.5, 0.5);

  for (int k = 0; k < 10; k++) {

    std::cout << "Histogram Create" << std::endl;
    {
      struct timeval t1, t2;
      gettimeofday(&t1, NULL);
      for (int i = 0; i < REPEAT_NUM; i++) {
      //for (int i = 0; i < 1; i++) {
        memset(histoCpu, 0, 256 * sizeof(int));
        createHistogramWithCPU(src.cols, src.rows, (unsigned char *)src.data, histoCpu);
      }
      gettimeofday(&t2, NULL);
      printf("Time CPU %f\n", (t2.tv_sec - t1.tv_sec) + (t2.tv_usec - t1.tv_usec)*1.0E-6);
    }

    {
      struct timeval t1, t2;
      gettimeofday(&t1, NULL);
      for (int i = 0; i < REPEAT_NUM; i++) {
        memset(histoGpu, 0, 256 * sizeof(int));
        createHistogramWithGPU(src.cols, src.rows, (unsigned char *)src.data, histoGpu);
      }
      gettimeofday(&t2, NULL);
      printf("Time GPU %f\n", (t2.tv_sec - t1.tv_sec) + (t2.tv_usec - t1.tv_usec)*1.0E-6);
    }

    cv::Mat histoImgCpu, histoImgGpu;
    drawHistogram(histoCpu, histoImgCpu);
    drawHistogram(histoGpu, histoImgGpu);

    cv::imshow("sample", src);
    cv::imshow("CPU", histoImgCpu);
    cv::imshow("GPU", histoImgGpu);
    cv::waitKey(1000);

  }

  delete histoCpu;
  delete histoGpu;

  return 0;
}

サンプルコード(こちらがデバイス側のコードです.一部ホスト側のコードも入ってますが..)
#include <cuda_runtime.h>
#include <stdio.h>

static void HandleError( cudaError_t err,
                         const char *file,
                         int line ) {
    if (err != cudaSuccess) {
        printf( "%s in %s at line %d\n", cudaGetErrorString( err ),
                file, line );
        exit( EXIT_FAILURE );
    }
}
#define HANDLE_ERROR( err ) (HandleError( err, __FILE__, __LINE__ ))

__global__ void createHistogramKernel(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 stride = blockDim.x * gridDim.x;
	
	while (i < size) {
		atomicAdd( &temp[buffer[i]], 1 );
		//atomicAdd( &(histo[buffer[i]]), 1 );
		i += stride;
	}
	
	__syncthreads();
	atomicAdd( &(histo[threadIdx.x]), temp[threadIdx.x] );
}

void createHistogramWithGPU(int width, int height, const unsigned char* buff, unsigned int *histo) {
	
	unsigned char *devBuffer;
	unsigned int *devHisto;
	
	// Memory Allocation.
	HANDLE_ERROR( cudaMalloc( (void **)&devBuffer, width * height ));
	HANDLE_ERROR( cudaMalloc( (void **)&devHisto, 256 * sizeof(unsigned int) ));
	
	HANDLE_ERROR( cudaMemcpy( devBuffer, buff, width * height, cudaMemcpyHostToDevice ));
	HANDLE_ERROR( cudaMemset( devHisto, 0, 256 * sizeof(unsigned int) ));

	cudaDeviceProp prop;
	cudaGetDeviceProperties( &prop, 0 );
	int blocks = prop.multiProcessorCount;
	//createHistogramKernel<<< blocks * 2, 256 >>>( devBuffer, width * height, devHisto );
	createHistogramKernel<<< 6 * 2, 256 >>>( devBuffer, width * height, devHisto );
	
	HANDLE_ERROR( cudaMemcpy( histo, devHisto, 256 * sizeof(unsigned int), cudaMemcpyDeviceToHost));

	// Memory Free
	HANDLE_ERROR( cudaFree( devBuffer ));
	HANDLE_ERROR( cudaFree( devHisto ));

}