图像处理系列：
空域滤波
这里根据模板是线性的还是非线性的进行分类，分别根据《图像处理与分析教程》上的方法实现了线性的领域平均、自定义加权平均、高斯平均、线性锐化和非线性的1-D中值滤波、2-D均值滤波、百分比滤波（最大、最小、中点）、非线性锐化、最大最小锐化，然后在主函数中进行测试。
一、线性滤波
1.空域线性滤波的代码
function [pic]=xxkylb(p,type,n,CX)
if nargin<3%缺省滤波核为3
n=3;
end
d=uint32(n/2)-1;%matlab uint会四舍五入！
if strcmp(type,'JZLB')%均值滤波器
K=ones(n,n);
K=uint32(K);
S=n*n;
elseif strcmp(type,'JQPJZDY')%自定义均值滤波器，滤波核由CX导入
K=CX;
K=uint32(K);
S=sum(K(:));
elseif strcmp(type,'GSJQPJ')%高斯均值滤波器
if nargin<4%缺省σ为1
sig=1;
else
sig=CX;
end
K=zeros(n,n);
gy=exp(double(-(((d)^2+(d)^2)/2*sig^2)))/(2*pi*sig^2);%右下角最小的参数，进行一个归一化
for x=1:n
for y=1:n
K(x,y)=exp(double((-(((x-d-1)^2+(y-d-1)^2)/2*sig^2))))/(2*pi*sig^2)/gy;%以距离值构造滤波核权术
end
end
K=uint32(K);%转换为整形
S=sum(K(:));
elseif strcmp(type,'RH')%线性锐化
K=CX;%模板为输入
S=1;%设为1
end
[r,c]=size(p);
mv=mean(p(:));
iob=ones(r+2*d,c+2*d)*mv;%以均值填充边缘，方便计算
iob((d+1):(d+r),(d+1):(d+c))=p;
pic=uint8(zeros(r,c));
for i=1:r
for j=1:c
t=0;
for k=1:n
for d=1:n
t=t+K(k,d)*iob(i+k-1,j+d-1);%计算
end
end
pic(i,j)=uint8(t/S);%换点
end
end
end2.结果分析
a.线性空域滤波效果
测试代码
clear;clc
pic=imread('lena.BMP');
pic=rgb2gray(pic);
picJZ=xxkylb(pic,'JZLB');
a=[1 2 1;2 4 2;1 2 1];%一个自定义的均值模板
piczdy=xxkylb(pic,'JQPJZDY',3,a);
picGS=xxkylb(pic,'GSJQPJ');
b=[0 -1 0;-1 4 -1;0 -1 0];%一个自定义的锐化模板
picRH=xxkylb(pic,'RH',3,b);
figure(6)
set(gcf,'position',[100 100 1500 800]);
subplot(2,3,1)
imshow(pic)
xlabel('原图')
subplot(2,3,2)
imshow(picJZ)
xlabel('均值滤波')
subplot(2,3,3)
imshow(piczdy)
xlabel('自定义均值滤波')
subplot(2,3,4)
imshow(picGS)
xlabel('高斯滤波')
subplot(2,3,5)
imshow(picRH)
xlabel('线性锐化（3x3无边缘中心4对角-1）')%模板[0 -1 0;-1 4 -1;0 -1 0]
saveas(gcf,'线性空域滤波效果图.jpg')效果图
图1 线性空域滤波效果图
首先这里变化比较明显的是线性锐化，可以看见图中只剩下lena的轮廓，但是也可以发现均值滤滤波和高斯滤波都将图片变得模糊了。
b.高斯滤波与均值滤波比较
测试代码：
figure(7)
subplot(2,2,1)
imshow(pic)
xlabel('原图')
subplot(2,2,2)
imshow(hdzh(pic-picJZ,'YZFQ',10))
xlabel('原图-均值滤波')
subplot(2,2,3)
imshow(hdzh(pic-piczdy,'YZFQ',10))
xlabel('原图-自定义均值滤波')
subplot(2,2,4)
imshow(hdzh(pic-picGS,'YZFQ',10))
xlabel('原图-高斯滤波')
saveas(gcf,'均值、高斯滤波比较.jpg')这里做一个减法然后将其进行了阈值为10的阈值分切见图4.2），可以发现原图减去滤波后的图阈值分切后，能看到lena的轮廓，也就是均值、高斯滤波对图像的边缘都有虚化作用。
图2 验证均值滤波和高斯滤波存在的变化
c.高斯均值滤波滤波核大小讨论
然后我们再对均值、高斯滤波的滤波核大小和滤波效果进行讨论，我们取滤波核大小为3、7、11、15、19分别用均值和高斯滤波查看效果（见图3），可以看到随n的增大，画面更加平滑，但是画面的内容逐渐模糊，细节逐渐减少。下面继续讨论不同σ的高斯模板对对滤波的影响，取滤波核为7，从效果（见图4）可以发现，随着σ的增大，高斯加权越平均，图片细节损失越小，测试代码和结果如下。
figure(8)
set(gcf,'position',[100 100 2500 2500]);
subplot(3,5,3)
imshow(pic)
xlabel('原图')
for i=1:5
subplot(3,5,5+i)
imshow(xxkylb(pic,'JZLB',(i-1)*4+3))
xlabel(strcat('均值，n=',num2str((i-1)*4+3)))
subplot(3,5,10+i)
imshow(xxkylb(pic,'GSJQPJ',(i-1)*4+3))
xlabel(strcat('高斯，n=',num2str((i-1)*4+3)))
end
saveas(gcf,'均值高斯滤波核讨论.jpg')
figure(9)
subplot(3,3,1)
imshow(pic)
xlabel('原图')
for i=1:8
subplot(3,3,i+1)
imshow(xxkylb(pic,"GSJQPJ",7,256*0.1*i))
xlabel(strcat('高斯（n=7），σ=',num2str(10*i),"%L"))
end
saveas(gcf,"高斯方差讨论.jpg")
图3 不同滤波核下高斯核滤波和均值滤波的差别
图4 不同方差下高斯核滤波的差别
二、非线性空域滤波
??? 然后对非线性空域滤波进行讨论，其函数块如下。
function [pic]=fxxkylb(p,type,XS)
[r,c]=size(p);
if strcmp(type,'1-DZZ')%1-D中值滤波
if nargin<3%对参数进行缺省，模板长度为5
XS=2;
end
mc=mean(p(:));
pic=uint8(zeros(r,c));
ioc=ones(r,c+2*XS)*mc;%补全方便计算
ioc(:,XS+1:c+XS)=p;
for x=1:r
for y=1:c
t=ioc(x,(y):(y+2*XS));
pic(x,y)=uint8(median(t));%取出模板中的中间值
end
end
else
if strcmp(type,'RH')%非线性锐化
if nargin<3
flag=1;
else
flag=XS;
end
XS=3;
d=uint32(XS/2)-1;
mv=mean(p(:));
iob=ones(r+2*d,c+2*d)*mv;
iob((d+1):(d+r),(d+1):(d+c))=p;
pic=uint8(zeros(r,c));
GX=[-1,0,1;-1,0,1;-1,0,1];%对x的差分
GY=[1,1,1;0,0,0;-1,-1,-1];%对y的差分
if flag==1
for i=1:r
for j=1:c
gx=0;
gy=0;
for k=1:XS
for d=1:XS
gx=gx+GX(k,d)*iob(i+k-1,j+d-1);%分别计算x y差分值
gy=gy+GY(k,d)*iob(i+k-1,j+d-1);
end
end
pic(i,j)=uint8((gx^2+gy^2)^0.5);%计算欧式距离
end
end
elseif flag==2
for i=1:r
for j=1:c
gx=0;
gy=0;
for k=1:XS
for d=1:XS
gx=gx+GX(k,d)*iob(i+k-1,j+d-1);%分别计算x y差分值
gy=gy+GY(k,d)*iob(i+k-1,j+d-1);
end
end
pic(i,j)=uint8(abs(gx)+abs(gy));%计算城区距离
end
end
else
for i=1:r
for j=1:c
gx=0;
gy=0;
for k=1:XS
for d=1:XS
gx=gx+GX(k,d)*iob(i+k-1,j+d-1);%分别计算x y差分值
gy=gy+GY(k,d)*iob(i+k-1,j+d-1);
end
end
pic(i,j)=uint8(max([gx gy]));%计算棋盘距离
end
end
end
else
if nargin<3%对模板大小进行缺省
XS=3;
end
d=uint32(XS/2)-1;
mv=mean(p(:));
iob=ones(r+2*d,c+2*d)*mv;
iob((d+1):(d+r),(d+1):(d+c))=p;
pic=uint8(zeros(r,c));
if strcmp(type,'2-DZZ')%2-D中值
for x=1:r
for y=1:c
t=[iob(x,y:(y+2*d)) iob(x+1,y:(y+2*d)) iob(x+2,y:(y+2*d))];
pic(x,y)=uint8(median(t));
end
end
elseif strcmp(type,'ZDZ')%最大值
for x=1:r
for y=1:c
t=[iob(x,y:(y+2*d)) iob(x+1,y:(y+2*d)) iob(x+2,y:(y+2*d))];
pic(x,y)=uint8(max(t));
end
end
elseif strcmp(type,'ZXZ')%最小hi
for x=1:r
for y=1:c
t=[iob(x,y:(y+2*d)) iob(x+1,y:(y+2*d)) iob(x+2,y:(y+2*d))];
pic(x,y)=uint8(min(t));
end
end
elseif strcmp(type,'ZD')%中点值
for x=1:r
for y=1:c
t=[iob(x,y:(y+2*d)) iob(x+1,y:(y+2*d)) iob(x+2,y:(y+2*d))];
pic(x,y)=uint8((max(t)+min(t))/2);
end
end
elseif strcmp(type,'ZDZX')%最大最小值，外部迭代
for x=1:r
for y=1:c
t=[iob(x,y:(y+2*d)) iob(x+1,y:(y+2*d)) iob(x+2,y:(y+2*d))];
if ((max(t)-p(x,y))<=(p(x,y)-min(t)))
pic(x,y)=uint8(max(t));
else
pic(x,y)=uint8(min(t));
end
end
end
end
end
end
end1.非线性滤波结果
首先以默认参数对对应的模板进行测试，测试代码。
%非线性
clear;clc
pic=imread('lena.BMP');
pic=rgb2gray(pic);
picD1=fxxkylb(pic,'1-DZZ');
picD2=fxxkylb(pic,'2-DZZ');
picRH=fxxkylb(pic,'RH');
picZDZ=fxxkylb(pic,'ZDZ');
picZXZ=fxxkylb(pic,'ZXZ');
picZD=fxxkylb(pic,'ZD');
picZDZX=fxxkylb(pic,'ZDZX');
figure(10)
set(gcf,'position',[100 100 2000 800]);
subplot(2,4,1)
imshow(pic)
xlabel("原图")
subplot(2,4,2)
imshow(picD1)
xlabel("1-D中值")
subplot(2,4,3)
imshow(picD2)
xlabel("2-D中值")
subplot(2,4,4)
imshow(picRH)
xlabel("非线性锐化")
subplot(2,4,5)
imshow(picZDZ)
xlabel("最大值")
subplot(2,4,6)
imshow(picZXZ)
xlabel("最小值")
subplot(2,4,7)
imshow(picZD)
xlabel("中点")
subplot(2,4,8)
imshow(picZDZX)
xlabel("最大最小")
saveas(gcf,'非线性效果.jpg')由效果（见图5）可以发现1-D中值在横向呈现平均的趋势，而2-D则是整个区域的平均，而相较于线性锐化（见图6）非线性锐化对轮廓的抓取显然更加好，且线性锐化对颗粒的筛选能力更加差；而最大、最小、中点，最大图片比原来偏亮，最小比原来偏暗，中点介于两者之间，但和原图有差别，三者图片更加平滑，存在细节丢失；而最大最小，色差较小，且明显能看到图像锐化增强。
图5 非线性滤波效果
a.线性锐化和非线性锐化对比
代码
figure(11)
subplot(1,3,1)
imshow(pic)
xlabel("原图")
subplot(1,3,2)
b=[-1 -1 -1;-1 8 -1;-1 -1 -1];%一个自定义的锐化模板
imshow(xxkylb(pic,"RH",3,b))
xlabel("线性锐化")
subplot(1,3,3)
imshow(picRH)
xlabel("非线性锐化")
saveas(gcf,'锐化效果对比.jpg')
图5 线性锐化和非线性锐化的对比
而相较于线性锐化（见图6）非线性锐化对轮廓的抓取显然更加好，且线性锐化对颗粒的筛选能力更加差。
b.1-D种植参数r调整变化
代码
figure(12)
subplot(3,3,1)
imshow(pic)
xlabel("原图")
for i=1:8
subplot(3,3,i+1)
imshow(fxxkylb(pic,'1-DZZ',uint32(512*0.051*i)))
xlabel(strcat('1-D中值，r=',num2str(i*5),"%l"))
end
saveas(gcf,'1-D中值r的讨论.jpg')
图6 r不同下1-D中值滤波的结果
1-D中值，我们看通过调整他的r的大小，效果的变化（见图6），可以看到，随着r的增大图片逐渐呈线状模糊，且其值越接近于0。
c.2-D中值滤波核大小调整变化
代码
figure(13)
subplot(3,3,1)
imshow(pic)
xlabel("原图")
for i=1:8
subplot(3,3,i+1)
imshow(fxxkylb(pic,'2-DZZ',(i-1)*4+3))
xlabel(strcat('2-D中值，n=',num2str((i-1)*4+3)))
end
saveas(gcf,'2-D中值n的讨论.jpg')
图7 n不同下2-D中值滤波的结果
对于2-D均值滤波，通用我们取不同的n，从效果（见图7）可以看出，随着n的增大图片逐渐模糊，细节逐渐丢失，但相较于线性滤波来说，他的细节保留会好一些。
d.非线性锐化距离计算调整
代码
figure(14)
subplot(2,2,1)
imshow(pic)
xlabel("原图")
subplot(2,2,2)
imshow(picRH)
xlabel("欧氏距离非线性锐化")
subplot(2,2,3)
imshow(fxxkylb(pic,'RH',2))
xlabel("城区距离非线性锐化")
subplot(2,2,4)
imshow(fxxkylb(pic,'RH',3))
xlabel("棋盘距离非线性锐化")
saveas(gcf,'非线性锐化对范数计算方式的讨论.jpg')
图4.8 范数计算不同下非线性锐化的结果
对于非线性锐化，我们进行范数计算不同的比较，从效果（见图8）上来看使用城区距离进行非线性锐化效果较好，它对人物边缘的增益更大。
e.最大最小锐化滤波迭代次数讨论
代码
figure(15)
subplot(3,3,1)
imshow(pic)
xlabel("原图")
for i=1:8
subplot(3,3,i+1)
pict=fxxkylb(pic,'ZDZX');
for j=1:3*(i-1)
pict=fxxkylb(pict,'ZDZX');
end
imshow(pict)
xlabel(strcat('最大最小锐化变换，n=',num2str(3*(i-1)+1)))
end
saveas(gcf,'最大最小锐化变换n的讨论.jpg')
图9 迭代次数不同的最大最小锐化变化的效果
最后对最大-最小锐化变换的迭代次数进行讨论，我们从效果（见图9）可以看出迭代次数越大，人物的轮廓越来越清晰，但是过大时，连照片中的颗粒也被锐化了出来。
三、附录
图10 lena.BMP
四、参考资料
《图像处理和分析教程》，章毓晋，第3版，人民邮电出版社