Matlab 图像滤波

*有一种孤独是
已经习惯了在某个人的庇护下生活，这个人离开之后，你不得不
面对现实，也渐渐学会了模仿他的样子去面对生活。——刘同《你的孤独 虽败犹荣》*

之前很少看书，在朋友推荐这本书时便喜欢上了这本书，以后在做研究时别总是抱着论文看个没完，是时候培养下阅读的习惯了。
这里主要是将自己课后作业贴出来，为了练习如何写博客。

添加椒盐噪声（比重：15%）：

Img=imread('4_1.bmp');
[h,w]=size(Img(:,:,1));
noiseimg=Img;
subplot(3,2,1);
imshow(Img);
title('Origin');
%add pepper & salt noise
n=round(h*w*0.15);

for k=1:n/2
    i=randi(h);
    j=randi(w);
    noiseimg(i,j,:)=0;
end
for k=1:n/2
    i=randi(h);
    j=randi(w);
    noiseimg(i,j,:)=255;
end
subplot(3,2,2);
imshow(noiseimg);
title('15% pepper&salt noise');

高斯低通滤波（3*3）：

%Gaussian low pass
kernel=[1 2 1;2 4 2;1 2 1];
kernel=kernel/16;
tmp=double(noiseimg)/255;
result=conv2(tmp(:,:,1),kernel,'same');
subplot(3,2,3);
imshow(result);
title('Gaussian filter 3x3');

中值滤波（3*3）：

%Median filter
padimg=padarray(noiseimg,[1,1],'replicate','both');
result=zeros(h,w,'uint8');
for i=1:h
    for j=1:w
        box=padimg(i:i+2,j:j+2,1);
        result(i,j)=median(box(:));
    end
end
subplot(3,2,4);
imshow(result);
title('Median filter 3x3');

双边滤波 (sigma_s=5, sigma_r=1)：

高斯滤波通过像素离中心点的距离来确定卷积核，比如常见的高斯滤波卷积核为

离核中心越远值越小，中心点值最大。

而双边滤波则在高斯滤波的基础之上加上对对应位置像素与中心点像素的差值的限制，差值越大则认为该位置与中心点之间存在着边缘，则该位置不能参与模糊计算，权值很小，这就是双边滤波能够保留边缘的原因。公式如下：

p为中心点位置，I表示图片像素值。

%Bilateral filter
sigma_s=7;
sigma_r=1;
result=zeros(h,w);
tmp=double(noiseimg)/255;
kernelsize=floor(sigma_s/2);
padimg=padarray(tmp,[kernelsize,kernelsize],'replicate','both');
for i=1:h
    for j=1:w
        box=padimg(i:i+kernelsize*2,j:j+kernelsize*2,1);
        dim=size(box);
        kernel=zeros(dim);
        reglar=0;
        for m=1:dim(1)
            for n=1:dim(2)
                kernel(m,n)=Guassian(sqrt((m-kernelsize-1)^2+(n-kernelsize-1)^2),sigma_s)*Guassian(box(m,n)-box(kernelsize+1,kernelsize+1),sigma_r);%关键代码
                reglar=reglar+kernel(m,n);
            end
        end
        tmp=kernel.*box;
        result(i,j)=sum(tmp(:))/reglar;
    end
end
subplot(3,2,5);
imshow(result);
title('Bilateral filter 5x5');

这里Guassian 是自己写的一个函数：

function y=Guassian(x,sigma)

y=1/(sqrt(2*pi*sigma*sigma))*exp(-x*x/(2*sigma*sigma));

结果图：