CUDA 和主机上的图像处理输出不同

Question

我正在使用 Quadro NVS 290 在 CUDA-C 中进行图像处理。为了验证 GPU 上的执行时间，我也在主机上进行处理。发现在GPU上的执行时间比CPU多，输出的图片也不一样。我在这里使用的算法是对 512x512 lema 图像进行三度模糊的高斯模糊。此代码也不适用于其他图像尺寸和灰度图像。

密码是：

unsigned int width, height;

int mask[3][3] = { 1, 2, 1,
2, 4, 2,
1, 2, 1
};

int h_getPixel(unsigned char *arr, int col, int row, int k)
{
    int sum = 0;
    int denom = 0;

    for (int j = -1; j <= 1; j++)
    {
        for (int i = -1; i <= 1; i++)
        {
            if ((row + j) >= 0 && (row + j) < height && (col + i) >= 0 && (col + i) < width)
            {
                int color = arr[(row + j) * 3 * width + (col + i) * 3 + k];
                sum += color * mask[i + 1][j + 1];
                denom += mask[i + 1][j + 1];
            }           
        }
    }

    return sum / denom;
} // End getPixel

void h_blur(unsigned char *arr, unsigned char *result)
{
    for (unsigned int row = 0; row < height; row++)
    {
        for (unsigned int col = 0; col < width; col++)
        {
            for (int k = 0; k < 3; k++)
            {
                result[3 * row * width + 3 * col + k] = h_getPixel(arr, col, row, k);
            }
        }
    }
} // End h_blur

__global__ void d_blur(unsigned char *arr, unsigned char *result, int width, int height)
{
    int col = blockIdx.x * blockDim.x + threadIdx.x;
    int row = blockIdx.y * blockDim.y + threadIdx.y;

    if (row < 0 || col < 0)
        return;

    int mask[3][3] = { 1, 2, 1,
        2, 4, 2,
        1, 2, 1
    };

    int sum = 0;
    int denom = 0;

    for (int k = 0; k < 3; k++)
    {
        for (int j = -1; j <= 1; j++)
        {
            for (int i = -1; i <= 1; i++)
            {
                if ((row + j) >= 0 && (row + j) < height && (col + i) >= 0 && (col + i) < width)
                {
                    int color = arr[(row + j) * 3 * width + (col + i) * 3 + k];
                    sum += color * mask[i + 1][j + 1];
                    denom += mask[i + 1][j + 1];
                }
            }
        }

        result[3 * row * width + 3 * col + k] = sum / denom;
    }
}

int  main(int argc, char **argv)
{
/************ Setup work ***********************/
unsigned char *d_resultPixels;
unsigned char *h_resultPixels;
unsigned char *h_devicePixels;

unsigned char *h_pixels = NULL;
unsigned char *d_pixels = NULL;

char *srcPath = .......;    // input image
char *h_resultPath = ......; // host output image
char *d_resultPath = ......; // device output image

FILE *fp_input;
FILE *fp_output;
FILE *fp_d_output;

unsigned char *inputFileData;
unsigned char *output_buffer;
unsigned char *d_output_buffer;

int nBlurDegree;

inputFileData = (unsigned char *)malloc(sizeof(unsigned char) * IMAGE_BUFFER_SIZE);
output_buffer = (unsigned char *)inputFileData;
d_output_buffer = (unsigned char *)inputFileData;

/* Read the uncompressed image file */
fp_input = fopen(srcPath, "r");

fread(inputFileData, IMAGE_BUFFER_SIZE, 1, fp_input);
fclose(fp_input);

unsigned int fileSize = (inputFileData[5] << 24) | (inputFileData[4] << 16) | (inputFileData[3] << 8) | inputFileData[2];
unsigned int dataOffset = (inputFileData[13] << 24) | (inputFileData[12] << 16) | (inputFileData[11] << 8) | inputFileData[10];
unsigned int imageSize = (inputFileData[37] << 24) | (inputFileData[36] << 16) | (inputFileData[35] << 8) | inputFileData[34];

width = (inputFileData[21] << 24) | (inputFileData[20] << 16) | (inputFileData[19] << 8) | inputFileData[18];
height = (inputFileData[25] << 24) | (inputFileData[24] << 16) | (inputFileData[23] << 8) | inputFileData[22];

h_pixels = (unsigned char *)malloc(imageSize);

h_resultPixels = (unsigned char *)malloc(imageSize);

inputFileData = inputFileData + dataOffset;
memcpy((void *)h_pixels, (void *)inputFileData, imageSize);

/************************** Start host processing ************************/

clock_t cpuStartTime, cpuEndTime;

cpuStartTime = clock();

// Apply gaussian blur
for (nBlurDegree = 0; nBlurDegree < BLUR_DEGREE; nBlurDegree++)
{
    memset((void *)h_resultPixels, 0, imageSize);

    h_blur(h_pixels, h_resultPixels);

    memcpy((void *)h_pixels, (void *)h_resultPixels, imageSize);        
}   

cpuEndTime = clock();

double cpuElapsedTime = (cpuEndTime - cpuStartTime) / (double)CLOCKS_PER_SEC;

printf("\nCPU time elapsed:\t%.2f ms\n", cpuElapsedTime * 1000);

inputFileData = inputFileData - dataOffset;

memcpy(output_buffer, inputFileData, dataOffset);

output_buffer = output_buffer + dataOffset;

memcpy(output_buffer, h_resultPixels, imageSize);

output_buffer = output_buffer - dataOffset;

fp_output = fopen(h_resultPath, "w");

fwrite(output_buffer, fileSize, 1, fp_output);
fclose(fp_output);

/************************** End host processing **************************/

/************************** Start device processing **********************/

cudaError_t cudaStatus;

h_devicePixels = (unsigned char *)malloc(imageSize);

cudaStatus = cudaMalloc((void **)&d_pixels, imageSize);

cudaStatus = cudaMalloc((void **)&d_resultPixels, imageSize);

cudaStatus = cudaMemcpy(d_pixels, h_pixels, imageSize, cudaMemcpyHostToDevice);

dim3 grid(16, 32);
dim3 block(32, 16);

// create CUDA event handles
cudaEvent_t gpuStartTime, gpuStopTime;
float gpuElapsedTime = 0;

cudaEventCreate(&gpuStartTime);
cudaEventCreate(&gpuStopTime);

cudaEventRecord(gpuStartTime, 0);   

for (nBlurDegree = 0; nBlurDegree < BLUR_DEGREE; nBlurDegree++)
{
    cudaStatus = cudaMemset(d_resultPixels, 0, imageSize);

    d_blur << < grid, block >> >(d_pixels, d_resultPixels, width, height);

    cudaStatus = cudaMemcpy(d_pixels, d_resultPixels, imageSize, cudaMemcpyDeviceToDevice);

    cudaStatus = cudaThreadSynchronize();
}

cudaEventRecord(gpuStopTime, 0);
cudaEventSynchronize(gpuStopTime);  // block until the event is actually recorded

cudaStatus = cudaMemcpy(h_devicePixels, d_resultPixels, imageSize, cudaMemcpyDeviceToHost);

cudaEventElapsedTime(&gpuElapsedTime, gpuStartTime, gpuStopTime);

printf("\nGPU time elapsed:\t%.2f ms\n", gpuElapsedTime);

memcpy(d_output_buffer, inputFileData, dataOffset);

d_output_buffer = d_output_buffer + dataOffset;

memcpy(d_output_buffer, h_devicePixels, imageSize);

d_output_buffer = d_output_buffer - dataOffset;

fp_d_output = fopen(d_resultPath, "w");

fwrite(d_output_buffer, fileSize, 1, fp_d_output);
fclose(fp_d_output);

/************************** End device processing ************************/

// Release resources
cudaEventDestroy(gpuStartTime);
cudaEventDestroy(gpuStopTime);

cudaFree(d_pixels);
cudaFree(d_resultPixels);

cudaThreadExit();

free(h_devicePixels);
free(h_pixels);
free(h_resultPixels);

return 0;
} // End main

Answer 1

你的代码有一个问题是你的数据流被破坏了。

h_pixels 包含您的初始数据：

memcpy((void *)h_pixels, (void *)inputFileData, imageSize);

您在主机模糊例程结束时覆盖了您的数据，结果数据：
```
memcpy((void *)h_pixels, (void *)h_resultPixels, imageSize);   
```

然后您将使用此模糊数据作为设备模糊例程的起点：

cudaStatus = cudaMemcpy(d_pixels, h_pixels, imageSize, cudaMemcpyHostToDevice);

在代码的第 2 步和第 3 步之间，您没有用原始起始数据替换 h_pixels 指向的数据。因此，期望设备模糊和主机模糊会产生相同的结果是不合理的。他们不是从相同的数据开始的。

您的代码的另一个问题是模糊操作的主机和设备代码之间存在细微差别。具体来说，在宿主案例 (h_blur) 中，每次调用 h_getPixel 时，变量 sum 和 denom 都被初始化为零（在 k 循环 h_blur).

但是在您的设备代码中，您有一个迭代 3 个颜色分量的循环，但是 sum 和 denom 在 [=20= 的每次迭代中都没有被重置为零]循环。

下面的完整示例修复了这些问题，并在主机和设备之间为随机样本数据产生了相同的结果：

$ cat t626.cu
#include <stdio.h>
#include <stdlib.h>

#define IMW 407
#define IMH 887
#define IMAGE_BUFFER_SIZE (IMW*IMH*3)
#define BLOCKX 16
#define BLOCKY BLOCKX
#define BLUR_DEGREE 3

unsigned int width, height;

int hmask[3][3] = { 1, 2, 1,
2, 4, 2,
1, 2, 1
};


#include <time.h>
#include <sys/time.h>
#define USECPSEC 1000000ULL

unsigned long long dtime_usec(unsigned long long prev){
  timeval tv1;
  gettimeofday(&tv1,0);
  return ((tv1.tv_sec * USECPSEC)+tv1.tv_usec) - prev;
}

int validate(unsigned char *d1, unsigned char *d2, int dsize){

  for (int i = 0; i < dsize; i++)
    if (d1[i] != d2[i]) {printf("validation mismatch at index %d, was %d, should be %d\n", i, d1[i], d2[i]); return 0;}
  return 1;
}

int h_getPixel(unsigned char *arr, int col, int row, int k)
{
    int sum = 0;
    int denom = 0;

    for (int j = -1; j <= 1; j++)
    {
        for (int i = -1; i <= 1; i++)
        {
            if ((row + j) >= 0 && (row + j) < height && (col + i) >= 0 && (col + i) < width)
            {
                int color = arr[(row + j) * 3 * width + (col + i) * 3 + k];
                sum += color * hmask[i + 1][j + 1];
                denom += hmask[i + 1][j + 1];
            }
        }
    }

    return sum / denom;
} // End getPixel

void h_blur(unsigned char *arr, unsigned char *result)
{
    for (unsigned int row = 0; row < height; row++)
    {
        for (unsigned int col = 0; col < width; col++)
        {
            for (int k = 0; k < 3; k++)
            {
                result[3 * row * width + 3 * col + k] = h_getPixel(arr, col, row, k);
            }
        }
    }
} // End h_blur

__global__ void d_blur(const unsigned char * __restrict__ arr, unsigned char *result, const int width, const int height)
{
    int col = blockIdx.x * blockDim.x + threadIdx.x;
    int row = blockIdx.y * blockDim.y + threadIdx.y;

    int mask[3][3] = { 1, 2, 1,
        2, 4, 2,
        1, 2, 1
    };
    if ((row < height) && (col < width)){
      int sum = 0;
      int denom = 0;

      for (int k = 0; k < 3; k++)
      {
        for (int j = -1; j <= 1; j++)
        {
            for (int i = -1; i <= 1; i++)
            {
                if ((row + j) >= 0 && (row + j) < height && (col + i) >= 0 && (col + i) < width)
                {
                    int color = arr[(row + j) * 3 * width + (col + i) * 3 + k];
                    sum += color * mask[i + 1][j + 1];
                    denom += mask[i + 1][j + 1];
                }
            }
        }

        result[3 * row * width + 3 * col + k] = sum / denom;
        sum = 0;
        denom = 0;
      }
    }
}

int  main(int argc, char **argv)
{
/************ Setup work ***********************/
  unsigned char *d_resultPixels;
  unsigned char *h_resultPixels;
  unsigned char *h_devicePixels;

  unsigned char *h_pixels = NULL;
  unsigned char *d_pixels = NULL;

  int nBlurDegree;
  int imageSize = sizeof(unsigned char) * IMAGE_BUFFER_SIZE;

  h_pixels = (unsigned char *)malloc(imageSize);


  width  = IMW;
  height = IMH;


  h_resultPixels = (unsigned char *)malloc(imageSize);
  h_devicePixels = (unsigned char *)malloc(imageSize);

  for (int i = 0; i < imageSize; i++) h_pixels[i] = rand()%30;
  memcpy(h_devicePixels, h_pixels, imageSize);

/************************** Start host processing ************************/
  unsigned long long cputime = dtime_usec(0);
// Apply gaussian blur
  for (nBlurDegree = 0; nBlurDegree < BLUR_DEGREE; nBlurDegree++)
  {
    memset((void *)h_resultPixels, 0, imageSize);

    h_blur(h_pixels, h_resultPixels);

    memcpy((void *)h_pixels, (void *)h_resultPixels, imageSize);
  }
  cputime = dtime_usec(cputime);


/************************** End host processing **************************/

/************************** Start device processing **********************/


  cudaMalloc((void **)&d_pixels, imageSize);

  cudaMalloc((void **)&d_resultPixels, imageSize);

  cudaMemcpy(d_pixels, h_devicePixels, imageSize, cudaMemcpyHostToDevice);

  dim3 block(BLOCKX, BLOCKY);
  dim3 grid(IMW/block.x+1, IMH/block.y+1);

  unsigned long long gputime = dtime_usec(0);

  for (nBlurDegree = 0; nBlurDegree < BLUR_DEGREE; nBlurDegree++)
  {
    cudaMemset(d_resultPixels, 0, imageSize);

    d_blur << < grid, block >> >(d_pixels, d_resultPixels, width, height);

    cudaMemcpy(d_pixels, d_resultPixels, imageSize, cudaMemcpyDeviceToDevice);
  }
  cudaDeviceSynchronize();
  gputime = dtime_usec(gputime);
  cudaMemcpy(h_devicePixels, d_resultPixels, imageSize, cudaMemcpyDeviceToHost);

  printf("GPU time: %fs, CPU time: %fs\n", gputime/(float)USECPSEC, cputime/(float)USECPSEC);

  validate(h_pixels, h_devicePixels, imageSize);
/************************** End device processing ************************/

// Release resources
  cudaFree(d_pixels);
  cudaFree(d_resultPixels);

  free(h_devicePixels);
  free(h_pixels);
  free(h_resultPixels);

  return 0;
} // End main
$ nvcc -O3 -o t626 t626.cu
$ ./t626
GPU time: 0.001739s, CPU time: 0.057698s
$

以上计时结果（GPU 比 CPU 快约 30 倍）是在 CentOS 5.5 和 CUDA 7 RC 上使用 Quadro5000 GPU 生成的。您的 Quadro NVS 290 是功率较低的 GPU，因此性能不佳。当我在 Quadro NVS 310 上运行这段代码时，我得到的结果表明 GPU 只比 CPU 快大约 2.5x

CUDA 和主机上的图像处理输出不同

Image processing output on CUDA and host are different

c

image-processing

cuda

gpgpu