图像滤波笔记
0、读取图片，预处理为灰度图(matlab代码)
%% read image imdir = 'F:\pictures\壁纸\'; imname = 'image (1)'; impath = strcat(imdir, imname, '.jpg'); ori_image = imread(impath); conv_image = ori_image; gray_image = rgb2gray(ori_image); [x, y] = size(gray_image);
1、均值滤波
均值滤波介绍：
相关：
卷积：
%% set mean filter kernel_x = 15; kernel_y = 15; kernel = (1/(kernel_x * kernel_y)) * ones(kernel_x, kernel_y); %% corr corr_image = uint8(imfilter(double(gray_image), double(kernel))); corr_image_path = strcat(imdir, imname, '_corr.jpg'); imwrite(corr_image, corr_image_path, 'jpg'); %% conv conv_image = uint8(conv2(double(gray_image), double(kernel), 'same')); conv_image_path = strcat(imdir, imname, '_conv.jpg'); imwrite(conv_image, conv_image_path, 'jpg');
左边：均值滤波+相关， 右边：均值滤波+卷积， 可以看出卷积和相关的结果区别不大
2、高斯滤波
%% set gaussian filter n = 10; o = (1+n)/2; sigma = 2; gk = double(zeros(n, n)); for i = 1:n for j = 1:n gk(i,j) = exp(-((i-o)^2 + (j-o)^2) / (2 * sigma^2)); end end fac = sum(sum(gk)); gk = gk / fac; %% conv with gaussian filter for i = 1:3 conv_image(:,:,i) = uint8(conv2(double(ori_image(:,:,i)), double(gk), 'same')); end conv_image_path = strcat(imdir, imname, '_gk.jpg'); imwrite(conv_image, conv_image_path, 'jpg');
3、定向高斯滤波
%% set oriented gaussian filter n = 20; o = (1+n)/2; sigma = 5; ogk = double(zeros(n, n)); a = 0.0000001; b= 10; sigma_matrix = [a, 0; 0, b]; orient1 = [1,1]; orient2 = [-1,1]; orient = [orient1;orient2]; sigma_matrix = orient' * sigma_matrix; sigma_matrix = sigma_matrix * orient; for i = 1:n for j = 1:n ogk(i,j) = 1 / exp(([i-o,j-o]*sigma_matrix*[i-o;j-o]) / (2 * sigma^2)); end end imshow(ogk); fac = sum(sum(ogk)); ogk = ogk / fac; %% conv with gaussian filter for i = 1:3 conv_image(:,:,i) = uint8(conv2(double(ori_image(:,:,i)), double(ogk), 'same')); end conv_image_path = strcat(imdir, imname, '_ogk.jpg'); imwrite(conv_image, conv_image_path, 'jpg');
可以明显感觉到方向：\\\\\
4、图像微分
%% image derivatives dx = [0,-1,0; -1,0,1; 0,1,0]; conv_image = uint8(conv2(double(gray_image), double(dx), 'same')); % constrain image range from 0 to 255 minv = min(min(conv_image)); if minv < 0 conv_image = conv_image - minv; end maxv = max(max(conv_image)); minv = min(min(conv_image)); factor = (255 - 0) / (maxv - minv + 1); conv_image = conv_image * factor + 0; % write image conv_image_path = strcat(imdir, imname, '_dx.jpg'); imwrite(conv_image, conv_image_path, 'jpg');
5、一阶高斯滤波（DoG）、方向一阶高斯滤波
高斯函数和各阶导数：http://blog.csdn.net/lyapple2008/article/details/8087317
由于图像微分会产生噪声，所以将高斯平滑与图像微分结合起来，就得到了DoG。
DoG加上方向就得到有方向的一阶高斯滤波:
%% set derivatives of gaussian (DoG) n = 8; o = (1+n)/2; sigma = 2; theta = -pi/4; dog = double(zeros(n, n)); dogx = dog; dogy = dog; ddog = dog; for i = 1:n for j = 1:n dogx(i,j) = -(i-o) * exp(-((i-o)^2 + (j-o)^2) / (2 * sigma^2)); dogy(i,j) = -(j-o) * exp(-((i-o)^2 + (j-o)^2) / (2 * sigma^2)); dog(i,j) = dogx(i,j) + dogy(i,j); ddog(i,j) = dogx(i,j)*cos(theta) + dogy(i,j)*sin(theta); end end %% conv with DoG % both x and y axis conv_image = conv2(double(gray_image), double(dog), 'same'); % constrain image range from 0 to 255 minv = min(min(conv_image)); if minv < 0 conv_image = conv_image - minv; end maxv = max(max(conv_image)); minv = min(min(conv_image)); factor = (255 - 0) / (maxv - minv + 1); conv_image = conv_image * factor + 0; % write image conv_image_path = strcat(imdir, imname, '_dog.jpg'); imwrite(uint8(conv_image), conv_image_path, 'jpg'); % DoG on vertical axis, remain the lines on horizontal axis conv_image = conv2(double(gray_image), double(dogx), 'same'); minv = min(min(conv_image)); if minv < 0 conv_image = conv_image - minv; end maxv = max(max(conv_image)); minv = min(min(conv_image)); factor = (255 - 0) / (maxv - minv + 1); conv_image = conv_image * factor + 0; conv_image_path = strcat(imdir, imname, '_dogx.jpg'); imwrite(uint8(conv_image), conv_image_path, 'jpg'); % DoG on horizontal axis, remain the lines on vertical axis conv_image = conv2(double(gray_image), double(dogy), 'same'); minv = min(min(conv_image)); if minv < 0 conv_image = conv_image - minv; end maxv = max(max(conv_image)); minv = min(min(conv_image)); factor = (255 - 0) / (maxv - minv + 1); conv_image = conv_image * factor + 0; conv_image_path = strcat(imdir, imname, '_dogy.jpg'); imwrite(uint8(conv_image), conv_image_path, 'jpg'); % DoG on derection of theta conv_image = conv2(double(gray_image), double(ddog), 'same'); minv = min(min(conv_image)); if minv < 0 conv_image = conv_image - minv; end maxv = max(max(conv_image)); minv = min(min(conv_image)); factor = (255 - 0) / (maxv - minv + 1); conv_image = conv_image * factor + 0; conv_image_path = strcat(imdir, imname, '_ddog.jpg'); imwrite(uint8(conv_image), conv_image_path, 'jpg');
下图，第一张为dogx滤波结果(可以看出保留了大量水平方向的边缘)，第二张为dogy滤波结果（可以看出保留了大量竖直方向上的边缘）
下图，第一张为DoG滤波结果（可以与上图两张比较一下，同时保留了两个方向的边缘），第二张为-pi/4方向上的DoG（可以看出保留了大量 “\\\\\\” 这样的边缘）
6、二阶高斯滤波（LoG）
LoG可以用DoG来近似：
%% set laplacian of gaussian (DoG) % n = 5; % o = (1+n)/2; % sigma = 0.8; % theta = pi/4; % log = double(zeros(n, n)); % logx = log; logy = log; dlog = log; % for i = 1:n % for j = 1:n % logx(i,j)=-(((i-o)^2 - sigma^2)/(2*pi*sigma^4))*exp(-((i-o)^2 + (j-o)^2)/(2 * sigma^2)); % logy(i,j)=-(((j-o)^2 - sigma^2)/(2*pi*sigma^4))*exp(-((i-o)^2 + (j-o)^2)/(2 * sigma^2)); % log(i,j) = logx(i,j) + logy(i,j); % dlog(i,j) = logx(i,j)*cos(theta) + logy(i,j)*sin(theta); % end % end log=[-2,-4,-4,-4,-2;... -4,0,8,0,-4;... -4,8,24,8,-4;... -4,0,8,0,-4;... -2,-4,-4,-4,-2]; %% conv with LoG conv_image = conv2(double(gray_image), double(log), 'same'); minv = min(min(conv_image)); if minv < 0 conv_image = conv_image - minv; end maxv = max(max(conv_image)); minv = min(min(conv_image)); factor = (255 - 0) / (maxv - minv + 1); conv_image = conv_image * factor + 0; conv_image_path = strcat(imdir, imname, '_log.jpg'); imwrite(uint8(conv_image), conv_image_path, 'jpg');
*上面的代码直接用LoG模版来计算的，这与上面被注释掉的代码是等价的。
LoG的介绍：
http://baike.soso.com/v61649603.htm
从上面的结果图可以看出，LoG的滤波结果比DoG更精细，可能是因为LoG抹去了更多的噪声。此外，log适合用在圆形的边缘上。
附：DoG和LoG的图形
sigma = 2; x=linspace(-10,10,200); y=linspace(-10,10,200); [xx,yy]=meshgrid(x,y); dog = ((-xx./(2*pi*sigma^4)).*exp(-((xx).^2 + (yy).^2)./(2 * sigma^2)))... + ((-yy./(2*pi*sigma^4)).*exp(-((xx).^2 + (yy).^2)./(2 * sigma^2))); log = (((xx.^2 - sigma^2)./(2*pi*sigma^4)).*exp(-(xx.^2 + yy.^2)./(2 * sigma^2)))... +(((yy.^2 - sigma^2)./(2*pi*sigma^4)).*exp(-(xx.^2 + yy.^2)./(2 * sigma^2))); figure; mesh(xx,yy,dog); figure; mesh(xx,yy,log);
