首先写好第一个和第二个函数,保存为.m文件,即:imfreqfilt(I,ff) 、imgaussflpf(I,sigma)
第一个函数功能:利用DFT,将空域原图 I 与 模板 卷积转为 频域乘积
function out = imfreqfilt(I,ff)
%imfreqfilt函数 对灰度图进行频域滤波
%参数I 输入的空间域图像
%参数ff 与原图等大小的频域滤镜
%标量(1x1)和向量(1xn)的ndims是2(matlab中);矩阵(mxn)是2
%多维数组:A(i1,i2,i3,i4,...,ik),那么A的ndims是k
%size()
if (ndims(I)==3) && (size(I,3)==3) %RGB图像
I = rgb2gray(I);
end;
%如果滤镜大小和原图不相等
if(size(I) ~= size(ff))
msg1 = sprintf('%s: 滤镜与原图大小不相等',mfilename);
msg2 = sprintf('%s: 滤波操作已取消',mfilename);
eid = sprintf('Image: %s:ImageSizeNotEqual',mfilename);
error(eid,'%s %s',msg1,msg2);
end;
% 快速DFT
f = fft2(I);
% 移动原点
s = fftshift(f);
%应用滤镜及变换
out = s.*ff;%对应元素相乘实现频域滤波
out = ifftshift(out);
out = ifft2(out);
%求模值
out = abs(out);
%归一化
out = out/max(out(:));
第二个函数功能, 构造“高斯低通滤波器” 模板
function out = imgaussflpf(I,sigma);
%imgaussflpf函数 构造频域高斯低通滤波器
%I参数 输入灰度图
%sigma参数 高斯函数的方差
[M,N]=size(I);
out = ones(M,N);
for i=1:M
for j=1:N
out(i,j)=exp(-((i-M/2)^2+(j-N/2)^2)/2/sigma^2);
end
end
第三个函数功能,基于前面两个函数,将显示原图在 “高斯低通滤波器” 的作用下,效果图的频谱,(本实验标准差sigma=0、20、40、60)
这里代码仅仅给出sigma=0、20情况,其他情况可直接更改sigma值
I = imread('lena.bmp');%这里仅能读取灰度图
%生成滤镜
ff = imgaussflpf(I,20);
%应用滤镜
out = imfreqfilt(I,ff);
figure(1);
subplot(221);
imshow(I);
title('source');
%显示原图频谱
temp = fft2(I);
temp = fftshift(temp);
temp = log(1+abs(temp));
figure(2);
subplot(221);
imshow(temp,[]);
title('source')
%显示效果图
figure(1);
subplot(222)
imshow(out);
title('Gauss LPF,sigma=20');
%显示效果图频谱
temp = fft2(out);
temp = fftshift(temp);
temp = log(1+abs(temp));
figure(2);
subplot(222);
imshow(temp,[]);
title('Gauss LPF,sigma=20');
参考《数字图像与机器视觉》