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
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
|
import pyopencl as cl
import numpy as np
from scipy.misc import imread, imsave
#Read in image
img = imread('noisyImage.jpg', flatten=True).astype(np.float32)
# Get platforms, both CPU and GPU
plat = cl.get_platforms()
CPU = plat[0].get_devices()
try:
GPU = plat[1].get_devices()
except IndexError:
GPU = "none"
#Create context for GPU/CPU
if GPU!= "none":
ctx = cl.Context(GPU)
else:
ctx = cl.Context(CPU)
# Create queue for each kernel execution
queue = cl.CommandQueue(ctx)
mf = cl.mem_flags
# Kernel function
src = '''
void sort(int *a, int *b, int *c) {
int swap;
if(*a > *b) {
swap = *a;
*a = *b;
*b = swap;
}
if(*a > *c) {
swap = *a;
*a = *c;
*c = swap;
}
if(*b > *c) {
swap = *b;
*b = *c;
*c = swap;
}
}
__kernel void medianFilter(__global float *img, __global float *result, __global int *width, __global int *height)
{
int w = *width;
int h = *height;
int posx = get_global_id(1);
int posy = get_global_id(0);
int i = w*posy + posx;
// Keeping the edge pixels the same
if( posx == 0 || posy == 0 || posx == w-1 || posy == h-1 )
{
result[i] = img[i];
}
else
{
int pixel00, pixel01, pixel02, pixel10, pixel11, pixel12, pixel20, pixel21, pixel22;
pixel00 = img[i - 1 - w];
pixel01 = img[i- w];
pixel02 = img[i + 1 - w];
pixel10 = img[i - 1];
pixel11 = img[i];
pixel12 = img[i + 1];
pixel20 = img[i - 1 + w];
pixel21 = img[i + w];
pixel22 = img[i + 1 + w];
//sort the rows
sort( &(pixel00), &(pixel01), &(pixel02) );
sort( &(pixel10), &(pixel11), &(pixel12) );
sort( &(pixel20), &(pixel21), &(pixel22) );
//sort the columns
sort( &(pixel00), &(pixel10), &(pixel20) );
sort( &(pixel01), &(pixel11), &(pixel21) );
sort( &(pixel02), &(pixel12), &(pixel22) );
//sort the diagonal
sort( &(pixel00), &(pixel11), &(pixel22) );
// median is the the middle value of the diagonal
result[i] = pixel11;
}
}
'''
#Kernel function instantiation
prg = cl.Program(ctx, src).build()
#Allocate memory for variables on the device
img_g = cl.Buffer(ctx, mf.READ_ONLY | mf.COPY_HOST_PTR, hostbuf=img)
result_g = cl.Buffer(ctx, mf.WRITE_ONLY, img.nbytes)
width_g = cl.Buffer(ctx, mf.READ_ONLY | mf.COPY_HOST_PTR, hostbuf=np.int32(img.shape[1]))
height_g = cl.Buffer(ctx, mf.READ_ONLY | mf.COPY_HOST_PTR, hostbuf=np.int32(img.shape[0]))
# Call Kernel. Automatically takes care of block/grid distribution
prg.medianFilter(queue, img.shape, None , img_g, result_g, width_g, height_g)
result = np.empty_like(img)
cl.enqueue_copy(queue, result, result_g)
# Show the blurred image
imsave('medianFilter-OpenCL.jpg',result)
|
