使用 OpenCL 实现图片高斯模糊

使用 OpenCL 实现图片高斯模糊

高斯模糊( https://zh.wikipedia.org/wiki/%E9%AB%98%E6%96%AF%E6%A8%A1%E7%B3%8A )是一种常见的图像处理算法,使用高斯分布与图像做卷积,得到模糊的效果。其二维定义:

σ是正态分布的标准偏差。在应用的时候,假设σ为2.5。对于模糊半径为1,则高斯矩阵为3*3的一个矩阵,以[1,1]为中心,带入公式计算高斯矩阵的值,得到:

||||
————|————–|———–
0.0216996633|0.0235069655|0.0216996633
0.0235069655|0.0254647918|0.0235069655
0.0216996633|0.0235069655|0.0216996633

他们的和为 0.206291318,我们需要他们的和为1,因此与总和相除得到:

0.105189413 0.113950342 0.105189413
0.113950342 0.123440929 0.113950342
0.105189413 0.113950342 0.105189413

根据这个矩阵,对图像的每个像素点进行计算,计算的九个点的各rgb分量之和就是最终像素的rgb分量。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
//计算高斯矩阵
private void ComputeWeightMatrix()
{
var center = Radius;
var conBase = 2 * Math.Pow(Variance, 2);
var conRoot = 1 / (Math.PI * conBase);

float sum = 0f;
for (int x = -Radius; x <= Radius; x++)
{
for (int y = Radius; y >= -Radius; y--)
{
var weight = conRoot * Math.Pow(Math.E, -(x * x + y * y) / conBase);
_matrix[GridPosToArrayIndex(x, y, center, Radius)] = (float)weight;
sum += (float)weight;
}
}
for (int i = 0; i < _matrix.Length; i++)
{
_matrix[i] /= sum;
}
}

//Compute
public void Compute(string imageFile)
{
using (var bitmap = new Bitmap(imageFile))
{
var datas = bitmap.LockBits(new Rectangle(new Point(), new Size(bitmap.Width, bitmap.Height)),ImageLockMode.ReadOnly,bitmap.PixelFormat);
var dataSize = datas.Stride * datas.Height;
var argbs = new byte[dataSize];
var dsts = new byte[dataSize];
int matrixWidth = Radius * 2 + 1;
Marshal.Copy(datas.Scan0, argbs, 0, dataSize);

Stopwatch sw=Stopwatch.StartNew();
for (int y = 0; y < bitmap.Height; y++)
{
for (int x = 0; x < bitmap.Width; x++)
{
float sumA = 0;
float sumR = 0;
float sumG = 0;
float sumB=0;
for (int i = 0; i < _matrix.Length; i++)
{
var pos = transform_pos(x, y, matrixWidth, bitmap.Width, bitmap.Height, Radius, i);
var position = pos.Y * datas.Stride + pos.X*4;
sumR += argbs[position] * _matrix[i];
sumG += argbs[position + 1] * _matrix[i];
sumB += argbs[position + 2] * _matrix[i];
sumA += argbs[position + 3] * _matrix[i];
}
var dstPos = y * datas.Stride + x * 4;
dsts[dstPos] = (byte)sumR;
dsts[dstPos+1] = (byte)sumG;
dsts[dstPos+2] = (byte)sumB;
dsts[dstPos+3] = (byte)sumA;
}
}
bitmap.UnlockBits(datas);

var elapse = sw.Elapsed;
Console.WriteLine($"Costing: {elapse}");
Debug.WriteLine($"Costing: {elapse}");

var handle = GCHandle.Alloc(dsts, GCHandleType.Pinned);
using (var dstBmp = new Bitmap(datas.Width, datas.Height, datas.Stride, bitmap.PixelFormat,
handle.AddrOfPinnedObject()))
{
dstBmp.Save("processed_normal.bmp");
}

handle.Free();
}
}

当然,这样能完成工作,但是耗时太长,对于3000*1920尺寸的图片处理需要2分51秒(Intel Core i7-4770),这显然是不可接受的。
对于这种分别计算每个像素,且各像素间互不干扰的问题,使用OpenCL可以大幅降低时间消耗。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
/*
OpenCL 高斯模糊代码
Copyright Gandalfliang
*/
inline int2 transform_pos(int centerX,int centerY,int matrixWidth,int radius,int index)
{
int x=index%matrixWidth;
int offsetX=x-(radius+1);
int y=index/matrixWidth;
int offsetY=radius-y;
return (int2)(centerX+offsetX,centerY-offsetY);
};

const sampler_t sampler_img=CLK_NORMALIZED_COORDS_FALSE|CLK_ADDRESS_CLAMP_TO_EDGE;

//opencl kernel 代码
kernel void gaussian_blur(
read_only image2d_t src,
global write_only char* dst,
global read_only float* matrix,
read_only int radius,
read_only int width)
{
int x=get_global_id(0);
int y=get_global_id(1);

float sumR,sumG,sumB,sumA;
int matrixWidth=radius*2+1;
int matrix_size=pow(matrixWidth,2);
for(int i=0;i<matrix_size;i++)
{
int2 pix=transform_pos(x,y,matrixWidth,radius,i);
uint4 rgba = read_imageui(src,sampler_img,pix);
sumR+=rgba.x*matrix[i];
sumG+=rgba.y*matrix[i];
sumB+=rgba.z*matrix[i];
sumA+=rgba.w*matrix[i];
}

int loc=y*width*4+x*4;
dst[loc]=sumR;
dst[loc+1]=sumG;
dst[loc+2]=sumB;
dst[loc+3]=sumA;
}

Host代码:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
public void Compute_cl(string imageFile)
{
//选取设备
var platform = ComputePlatform.Platforms.FirstOrDefault();
var device = platform.Devices.FirstOrDefault();
//设置相关上下文
var properties = new ComputeContextPropertyList(platform);
var context = new ComputeContext(new[] {device}, properties, null, IntPtr.Zero);
//命令队列,用于控制执行的代码
ComputeCommandQueue commands = new ComputeCommandQueue(context, context.Devices[0],
ComputeCommandQueueFlags.None);
//读取opencl代码
var code = File.ReadAllText(@"gaussianblur.cl");
//编译
var program = new ComputeProgram(context, code);
try
{
program.Build(new[] {device}, null, null, IntPtr.Zero);
}
catch (Exception ex)
{
throw;
}

var images = CreateImageFromBitmap(imageFile, context,
ComputeMemoryFlags.ReadWrite | ComputeMemoryFlags.CopyHostPointer);

//创建核心代码,就是cl代码中以kernel标识的函数
var kernel = program.CreateKernel("gaussian_blur");
//矩阵规模
//储存计算结果的数组

//创建的核心代码函数以这种方式来传参
var resultBuffer=new ComputeBuffer<char>(context,ComputeMemoryFlags.WriteOnly, dstBytes.Length);
kernel.SetMemoryArgument(0, images);
kernel.SetMemoryArgument(1, resultBuffer);
kernel.SetMemoryArgument(2, new ComputeBuffer<float>(context,ComputeMemoryFlags.ReadOnly|ComputeMemoryFlags.CopyHostPointer,_matrix));
kernel.SetValueArgument(3, Radius);
kernel.SetValueArgument(4, (int)images.Width);
Console.WriteLine($"运行平台: {platform.Name}\n运行设备: {device.Name}\n");
Stopwatch sw = Stopwatch.StartNew();
var climg = images;

//执行代码
commands.Execute(kernel, null, new long[] {climg.Width, climg.Height}, null, null);

//read data
char[] resultArray = new char[dstBytes.Length];
var arrHandle = GCHandle.Alloc(resultArray, GCHandleType.Pinned);
commands.Read(resultBuffer, true, 0, dstBytes.Length, arrHandle.AddrOfPinnedObject(), null);
//commands.ReadFromImage(images.Item2, processeddata.Scan0, true, null);

var resultHandle = GCHandle.Alloc(resultArray, GCHandleType.Pinned);
var bmp=new Bitmap(climg.Width,climg.Height, climg.Width*4, PixelFormat.Format32bppArgb, resultHandle.AddrOfPinnedObject());
var elapsed = sw.Elapsed;
Console.WriteLine($"耗时: {elapsed.TotalMilliseconds} ms\n");
kernel.Dispose();

bmp.Save("processed_cl.bmp");
}

相同尺寸的图片处理,使用 Intel Core i7-4770 自带的核显 HD4600 处理,耗时只需要164毫秒。

以下是相关测试结果:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
运行平台: Intel(R) OpenCL
运行设备: Intel(R) HD Graphics 4600
处理图片尺寸:790*501
OpenCL处理耗时: 13.6597 ms

处理图片尺寸:790*501
常规方法耗时: 11482.9402 ms

运行平台: Intel(R) OpenCL
运行设备: Intel(R) HD Graphics 4600
处理图片尺寸:1339*693
OpenCL处理耗时: 33.0095 ms

处理图片尺寸:1339*693
常规方法耗时: 26908.9926 ms

运行平台: Intel(R) OpenCL
运行设备: Intel(R) HD Graphics 4600
处理图片尺寸:1920*1080
OpenCL处理耗时: 51.3885 ms

处理图片尺寸:1920*1080
常规方法耗时: 60147.3815 ms

当然,常规方法都只使用了单线程,还未发挥多核CPU的威力,然而,可以预见的是,即使是使用多线程,提升也是有限的。

原图:

高斯模糊:

代码: https://github.com/gandalfliang/cloo_netstandard/tree/temp

Update: 在nVidia的环境下会导致处理后的图片出现花屏现象,估计是cl代码的问题,又或者是nVidia的驱动有问题?下次再更新

评论