本文主要是介绍CUDA：用并行计算的方法对图像进行直方图均衡处理，希望对大家解决编程问题提供一定的参考价值，需要的开发者们随着小编来一起学习吧！
（一）目的
将所学算法运用于图像处理中。
（二）内容
用并行计算的方法对图像进行直方图均衡处理。
要求：
- 利用直方图均衡算法处理lena_salt图像
版本1：CPU实现
版本2：GPU实现
实验步骤一 软件设计分析:
- 数据类型：
根据实验要求，本实验的数据类型为一个256*256*8的整型矩阵，其中元素的值为256*256个0-255的灰度值。
- 存储方式：
图像在内存中的存储方式主要是以二维矩阵的方式进行存储，这里的lena_salt图像是一个256*256的矩阵，每一个元素用一个字节来存储像素值。
矩阵在内存中的存储按照行列优先可以分为两种方式，一种是行优先的存储方式，一种是按照列优先的方式。
这两种存储方式在访问对应的位置的数据的时候有很大的差别。在cuda内部，矩阵默认是按照列优先的方式存储，如果要使用cuda device函数，就必须考虑存储方式的问题，有的时候可能需要我们队存储方式进行装换。但是无论是用那种存储方式，最终在内存中都是顺序存储的。
三．GPU程序的block和threads的相关设置：
本实验提供的英伟达实验平台每一个Grid可以按照一维或者二维的方式组织，每一个Block可以按照一维，二维或者三维的方式进行组织。每一个block最多只能有1536个线程。内核函数使用的线程总量也受到设备本身的限制。
对于本次实验，针对上文中提到的几个任务，block和threads的组织方式都可以描述为：
dim3 threadsPerBlock(16, 16);
dim3 blocksPerGrid((img_in.w + 15) / 16, (img_in.h + 15) / 16);
实验步骤二 实验设备:
本地设备：PC机+Windows10操作系统
Putty远程连接工具
PsFTP远程文件传输工具
远程设备：NVIDIA-SMI 352.79
Driver Version:352.79
实验步骤三 CPU计算代码:
void cpu_ histogram_equalization_kernel(uchar * img_in, uchar * img_out, int img_w, int img_h){
//数据定义
float pixel_value[256];
float P_pixel_value[256];
float Sum_P_pixel_value[256];
for (int i = 0; i < 256; i++)
pixel_value[i] = 0;
//统计直方图
for (int i = 0; i<img_h; i++)
{
for (int j = 0; j<img_w; j++)
{
pixel_value[img_in[i * img_w + j]]++;
//img_out[i*img_w+j] = pixel_out;
}
}
//概率直方图
for (int i = 0; i < 256; i++)
P_pixel_value[i] = pixel_value[i] / 256 / 256;
//前项概率求和
for (int i = 0; i < 256; i++)
{
float sum = 0.0;
for (int j = 0; j < i; j++)
sum += P_pixel_value[j];
Sum_P_pixel_value[i] = sum;
}
//均衡化
for (int i = 0; i<img_h; i++)
{
for (int j = 0; j<img_w; j++)
{
int pixel_out = int(Sum_P_pixel_value[img_in[i * img_w + j]] * 256 + 0.5);
img_out[i*img_w + j] = pixel_out;
}
}
}
实验步骤四 GPU计算代码:
- 新增数据定义及初始化部分
float *pixel_value_h = new float[256];
float *P_pixel_value_h = new float[256];
float *Sum_P_pixel_value_h = new float[256];
for (int i = 0; i < 256; i++)
{
pixel_value_h[i] = 0.0;
P_pixel_value_h[i] = 0.0;
Sum_P_pixel_value_h[i] = 0.0;
}
float *pixel_value_d;
float *P_pixel_value_d;
float *Sum_P_pixel_value_d;
cudaMalloc((void **)&pixel_value_d, 256 * sizeof(double));
cudaMalloc((void **)&P_pixel_value_d, 256 * sizeof(float));
cudaMalloc((void **)&Sum_P_pixel_value_d, 256 * sizeof(float));
cudaMemcpy(pixel_value_d, pixel_value_h, 256 * sizeof(float), cudaMemcpyHostToDevice);
cudaMemcpy(P_pixel_value_d, P_pixel_value_h, 256 * sizeof(float), cudaMemcpyHostToDevice);
cudaMemcpy(Sum_P_pixel_value_d, Sum_P_pixel_value_h, 256 * sizeof(float), cudaMemcpyHostToDevice);
- 函数定义部分
//统计直方图
__global__ void gpu_histogram_sta_kernel(uchar * img_in, uchar * img_out, int img_w, int img_h, float *pixel_value)
{
int row = blockDim.y * blockIdx.y + threadIdx.y;
int col = blockDim.x * blockIdx.x + threadIdx.x;
if ((row >= 0) && (row < img_h) && (col >= 0) && (col < img_w))
//pixel_value[img_in[row*img_w + col]]++;
atomicAdd(&pixel_value[img_in[row*img_w + col]], 1.0);
}
//计算概率
__global__ void gpu_probability_sta_kernel(uchar * img_in, uchar * img_out, int img_w, int img_h, float *pixel_value, float *P_pixel_value)
{
if ((blockIdx.x == 0) && (blockIdx.y == 0) && (threadIdx.x == 0) && (threadIdx.y == 0))
for (int i = 0; i < 256;i++)
P_pixel_value[i] = pixel_value[i] / 256 /256;
}
//计算概率前项和
__global__ void gpu_sum_probability_kernel(uchar * img_in, uchar * img_out, int img_w, int img_h, float *P_pixel_value, float *Sum_P_pixel_value)
{
__shared__ float sharedM[256];
int i = threadIdx.x + blockIdx.x * blockDim.x;
unsigned int bid = blockIdx.y * gridDim.x + blockIdx.x;
unsigned int bid = blockIdx.x;
unsigned int tid = threadIdx.x;
unsigned int count = 1;
sharedM[tid] = P_pixel_value[tid];
__syncthreads();
if (bid % 2 == 0)
{
for (unsigned int stride = 1; stride < bid; stride *= 2)
{
__syncthreads();
if (tid % (2*stride) == 0)
sharedM[tid] += sharedM[tid + stride];
}
}
else
{
for (unsigned int stride = 1; stride < bid + 1; stride *= 2)
{
__syncthreads();
if (tid % (2 * stride) == 0)
sharedM[tid] += sharedM[tid + stride];
}
}
if (tid == 0)
Sum_P_pixel_value[blockIdx.x] = sharedM[0];
}
//均衡化
__global__ void gpu_equilibrium_kernel(uchar * img_in, uchar * img_out, int img_w, int img_h, float *Sum_P_pixel_value)
{
int row = blockDim.y * blockIdx.y + threadIdx.y;
int col = blockDim.x * blockIdx.x + threadIdx.x;
img_out[row*img_w + col] = int(Sum_P_pixel_value[img_in[row*img_w + col]] * 256 + 0.5);
}
- 函数调用
gpu_histogram_sta_kernel << <blocksPerGrid, threadsPerBlock >> >(d_img_in, d_img_out, img_in.w, img_in.h, pixel_value_d);
gpu_probability_sta_kernel << <blocksPerGrid, threadsPerBlock >> >(d_img_in, d_img_out, img_in.w, img_in.h, pixel_value_d, P_pixel_value_d);
gpu_sum_probability_kernel << <256, 256 >> >(d_img_in, d_img_out, img_in.w, img_in.h, P_pixel_value_d, Sum_P_pixel_value_d);
gpu_equilibrium_kernel << <blocksPerGrid, threadsPerBlock >> >(d_img_in, d_img_out, img_in.w, img_in.h, Sum_P_pixel_value_d);
实验步骤五 观察输出结果:
图1 原图像
- 版本1：CPU实现处理结果
图2 CPU实现处理效果
图3 CPU实现处理时间
2，版本2：GPU实现处理结果
图4 CPU实现处理效果
图5 CPU实现处理时间
3，处理过程中的数据
图6 直方图统计结果
图7 概率计算结果
图8 概率前项求和
实验结论：
cpu程序计算所需时间：
版本1，CPU实现程序计算所需时间：1.6711328ms
gpu程序计算所需时间：
版本2，GPU实现程序计算所需时间：2.950976ms
总结：
之前的实验都是讲所有的代码写在一个kernel函数里面，本次实验突发奇想的采用多个kernel函数对直方图均衡的每一步分别进行处理，也算是一种新的尝试吧。在实验的过程中，由于远程端的运行环境导致调试代码，特别是排查错误显得很艰难。我在这里才取的解决办法就是将处理完的数据传回host端，然后打印出来，观察输出结果是否符合预期。这样就很容易发现处理的过程中是哪一步出了问题，方便了错误排查。
这篇关于CUDA：用并行计算的方法对图像进行直方图均衡处理的文章就介绍到这儿，希望我们推荐的文章对编程师们有所帮助！
![业务中14个需要进行A/B测试的时刻[信息图]](https://img-blog.csdnimg.cn/img_convert/aeacc959fb75322bef30fd1a9e2e80b0.jpeg)
