Gpu实现直方图均衡化
之前博客介绍过cpu版本的 直方图均衡化算法原理与实现,在此基础上介绍下gpu版本实现直方图均衡化。
Gpu版本的算法实现流程和cpu版本一致,也是分为如下几个步骤:
① 统计每个灰度级像素个数
② 计算累积概率密度函数
③ 根据累积概率密度求灰度映射表
④ 根据映射表计算映射后的灰度值
与cpu的区别是上述步骤是以并行的方式计算的。
- 统计每个灰度级像素个数
创建共享内存对象数组temp, 用于统计每个block内灰度的个数。
__syncthreads()同步所有线程,确保每个线程都计算结束
使用原子add, 将每个block的灰度个数相加,计算总的灰度等级个数。
__global__ void calhist_kernel(unsigned char* devimage, float* devhist,
const int length)
{
__shared__ unsigned int temp[HIST_SZ];
temp[threadIdx.x] = 0;
__syncthreads();
int i = threadIdx.x + blockIdx.x* blockDim.x;
int offset = blockDim.x* gridDim.x;
while (i < length)
{
atomicAdd(&temp[devimage[i]], 1);
i += offset;
}
__syncthreads();
atomicAdd(&devhist[threadIdx.x], temp[threadIdx.x]);
}
- 计算每个灰度级的概率
每个线程处理一个元素,计算每个灰度级的概率
__global__ void calhist_ratio_kernel(float*devhist, int total, float* devhistratio)
{
int id = threadIdx.x;
devhistratio[id] = devhist[id] / total * (HIST_SZ - 1);//将概率归一到0-255
}
- 计算累积概率密度
__global__ void sum_ratio_kernel_v1(float* devhistratio, float *devsumratio)
{
int index = threadIdx.x;
for (unsigned int stride = 1; stride <= threadIdx.x; stride *= 2) {
float in1 = devhistratio[threadIdx.x - stride];
__syncthreads();
devsumratio[threadIdx.x] += in1;
__syncthreads();
}
}
- 根据映射表计算映射后的灰度值
根据累加概率密度计算一个映射表,每个线程可以处理一个灰度等级。
__global__ void cal_lut_kernel(float* devsumratio, unsigned int* devlut)
{
int id = threadIdx.x;
devlut[id] = int(devsumratio[id] + 0.5f);
}
- 计算映射后的灰度值
__global__ void cal_map_image(unsigned char* devsrc, unsigned int * devlut, unsigned char* devdst, int height, int width)
{
int x = blockIdx.x* blockDim.x + threadIdx.x;
int y = blockIdx.y*blockDim.y + threadIdx.y;
int index = y * width + x;
if (x < width&& y < height)
{
int temp = devsrc[index];
devdst[index] = devlut[temp];
}
}