二维快速傅里叶变换及滤波

摘 要：在图象处理的应用领域中，傅立叶变换起着非常重要的作用，具体表现在包括图象分析、图象增强及图象压缩等方面。本实验采用二维快速傅里叶变换对图像进行频域滤波,其中编写了一个2-D FFT程序，一个高斯低通滤波器的程序，对图像进行滤波,并通过原图减去高斯低通滤波后的图像实现锐化。

---

function towDFFT
%2—D FFT

im=imread('Fig4.41(a).jpg');
im=imresize(im,[512,512]);
subplot(241);imshow(im);title('原图像');
ima=double(im);
[m,n]=size(ima);
p=2*m;q=2*n;
%填充图像
imb=zeros(p,q);
for i=1:m
    for j=1:n
        imb(i,j)=ima(i,j);
    end
end
subplot(242);imshow(imb);title('填充图');
imwrite(imb,'填充图.jpg');

%用因子(-1)^(x+y)乘以输入图像
for i=1:p
    for j=1:q
        imb(i,j)=power(-1,i+j)*imb(i,j);
    end
end
imb=uint8(imb);
subplot(243);imshow(imb);title('乘因子之后的图');
imwrite(imb,'乘因子后的图.jpg');

%计算傅里叶变换
F=fft2(imb);
subplot(244);imshow(abs(F),[]);title('F的频谱');

%生成一个实的、对称的滤波函数H
%调用了一个高斯低通滤波函数
H=gause(p,q,20);
subplot(245);imshow(abs(H),[]);title('滤波函数频谱');
imwrite(abs(H),'滤波函数频谱.jpg');
G=H.*F;
subplot(246);imshow(abs(G),[]);title('H.*F的频谱');
imwrite(abs(G),'滤波后的频谱.jpg');

%计算反傅里叶变换
g=ifft2(G);
for i=1:p
    for j=1:q
        g(i,j)=power(-1,i+j)*real(g(i,j));
    end
end
g=uint8(g);
subplot(247);imshow(g);title('乘因子取实部');
imwrite(g,'乘因子取实部.jpg');

%提取m*n区域，得到最终处理结果
imc=zeros(m,n);
for i=1:m
    for j=1:n
        imc(i,j)=g(i,j);
    end
end
imc=uint8(imc);
subplot(248);imshow(imc);title('输出图像');
imwrite(imc,'输出图像.jpg');

%使用二维图像和三维图像两种形式显示幅度谱和相位谱
%计算输出图像频谱
Fnew=fftshift(fft2(imc));
figure;subplot(221);imshow(abs(Fnew),[]);title('幅度谱');
subplot(222);imshow(angle(Fnew));title('相位谱');
[a,b]=meshgrid(1:m,1:n);
subplot(223);mesh(a,b,abs(Fnew));colormap(gray);title('幅度谱');
subplot(224);mesh(a,b,angle(Fnew));colormap(gray);title('相位谱');

%计算4.41（a）的中心化后的傅里叶谱
FF=fftshift(fft2(im));
figure;subplot(221);imshow(abs(FF),[]);title('幅度谱');
subplot(222);imshow(angle(FF));title('相位谱');
[a,b]=meshgrid(1:m,1:n);
subplot(223);mesh(a,b,abs(FF));colormap(gray);title('幅度谱');
subplot(224);mesh(a,b,angle(FF));colormap(gray);title('相位谱');

附录2
function gd401
%实现对Fig4.41(a).jpg的高斯高通和高斯低通滤波
ima=imread('Fig4.41(a).jpg');
ima=imresize(ima,[512,512]);
ima=double(ima);

%调用高斯低通滤波函数
imb=gauslow(ima,10);
imc=gauslow(ima,30);
imd=gauslow(ima,60);
ime=gauslow(ima,160);
imf=gauslow(ima,460);
figure;
subplot(231);imshow(ima);title('原图像');
subplot(232);imshow(imb);title('D0=10');
subplot(233);imshow(imc);title('D0=30');
subplot(234);imshow(imd);title('D0=60');
subplot(235);imshow(ime);title('D0=160');
subplot(236);imshow(imf);title('D0=460');

%利用原图减高斯低通滤波后的图像，得到锐化图像
imC=imsubtract(ima,imc);
imD=imsubtract(ima,imd);
imE=imsubtract(ima,ime);
figure;subplot(131);imshow(imC);title('D0=30');
subplot(132);imshow(imD);title('D0=60');
subplot(133);imshow(imE);title('D0=160');
%------------------------------------------------------------------
function im=gauslow(ima,D0)
%编写一个高斯低通滤波函数对图像进行滤波
[m,n]=size(ima);
r=2*m;
c=2*n;
%填充图像
imb=zeros(r,c);
for i=1:m
    for j=1:n
        imb(i,j)=ima(i,j);
    end
end
%利用fftshift函数对频谱进行中心化
%进行高斯滤波
imF=fftshift(fft2(imb));
for i=1:r
    for j=1:c
        H(i,j)=1*exp(-((i-m)*(i-m)+(j-n)*(j-n))/(2*D0*D0));
        imG(i,j)=H(i,j)*imF(i,j);
    end
end
img=real(ifft2(imG));
for i=1:r
    for j=1:c
        img(i,j)=power(-1,i+j)*real(img(i,j));
    end
end
img=uint8(img);
im=zeros(m,n);
for i=1:m
    for j=1:n
        im(i,j)=img(i,j);
    end
end
im=uint8(im);
return;