轉載請註明來自我的CSDN部落格：黃朝輝的部落格

1. 啟動MATLAB程式，讀入一幅影象；對影象做FFT。使用’subplot’命令，同時顯示原始影象其頻譜圖；

IenaImg=imread('lena.jpg');   %讀入原影象檔案
fftI=fft2(double(IenaImg));        %二維離散傅立葉變換
sfftI=fftshift(fftI);   %直流分量移到頻譜中心
RR=real(sfftI);     %取傅立葉變換的實部
fftImg = log(RR);
subplot(1,2,1)
imshow(IenaImg);        %顯示原影象
 

subplot(1,2,2)
imshow(fftImg, [8,10]);       %顯示原影象的頻譜

2. 讀入一幅影象，分別為影象新增椒鹽、高斯噪聲，做FFT變換。使用’subplot’命令，將原始影象、原始影象頻譜圖、新增噪聲後的影象，以及噪聲影象的頻譜圖同時顯示出來。

lenaImg=imread('lena.jpg');   %讀入原影象檔案
fftI=fft2(double(lenaImg));        %二維離散傅立葉變換
sfftI=fftshift(fftI);   %直流分量移到頻譜中心
RR=real(sfftI);     %取傅立葉變換的實部
fftLenaImg = log
 
(RR);

加入椒鹽躁聲

saltImg=imnoise(lenaImg,'salt & pepper',0.02); %加入椒鹽躁聲
fftI=fft2(double(saltImg));        %二維離散傅立葉變換
sfftI=fftshift(fftI);   %直流分量移到頻譜中心
RR=real(sfftI);     %取傅立葉變換的實部
fftSaltImg = log(RR);

subplot(2,2,1);
imshow(lenaImg);
subplot(2,2,2);
imshow(fftLenaImg , [8,10]);
subplot(2,2
 
,3);
imshow(saltImg);
subplot(2,2,4);
imshow(fftSaltImg , [8,10]);

加入高斯躁聲

gaussianImg =imnoise(lenaImg,'gaussian',0.02);     %加入高斯躁聲 
fftI=fft2(double(gaussianImg));        %二維離散傅立葉變換
sfftI=fftshift(fftI);   %直流分量移到頻譜中心
RR=real(sfftI);     %取傅立葉變換的實部
fftGaussianImg = log(RR);

subplot(2,2,1);
imshow(lenaImg);
subplot(2,2,2);
imshow(fftLenaImg , [8,10]);
subplot(2,2,3);
imshow(gaussianImg);
subplot(2,2,4);
imshow(fftGaussianImg , [8,10]);

3. 讀入一幅影象，對影象分別進行高斯低通、巴特沃茲低通、高斯高通和巴特沃茲高通頻域濾波，比較其銳化和平滑效果；

%高斯低通濾波器

I = imread('lena.jpg'); 
imshow(I);
s=fftshift(fft2(I));
[M,N]=size(s);                    
n=2;                                  
d0=30; %GLPF濾波，d0=5，15，30（程式中以d0=30為例）                    
n1=floor(M/2);                          
n2=floor(N/2);                           
for i=1:M 
    for j=1:N
        d=sqrt((i-n1)^2+(j-n2)^2);         
h=1*exp(-1/2*(d^2/d0^2));  
        s(i,j)=h*s(i,j);                   
    end
end
s=ifftshift(s);                           
gImg=uint8(real(ifft2(s)));                                     
figure;                                
imshow(gImg); 

%巴特沃思低通濾波器

l = imread('test.jpg'); 
s=fftshift(fft2(l));
[M,N]=size(s);    
%5級巴特沃思                
n=5;                                  
d0=30; %BLPF濾波，d0=5，15，30（程式中以d0=30為例）                    
n1=floor(M/2);                          
n2=floor(N/2);                           
for i=1:M 
    for j=1:N
        d=sqrt((i-n1)^2+(j-n2)^2);         
h=1/(1+(d/d0)^(2*n));   
        s(i,j)=h*s(i,j);                   
    end
end
s=ifftshift(s);                           
b5Img=uint8(real(ifft2(s)));        

%15級巴特沃思

n=15;                
s=fftshift(fft2(l));
[M,N]=size(s);                    
d0=30; %BLPF濾波，d0=5，15，30（程式中以d0=30為例）                    
n1=floor(M/2);                          
n2=floor(N/2);                           
for i=1:M 
    for j=1:N
        d=sqrt((i-n1)^2+(j-n2)^2);         
h=1/(1+(d/d0)^(2*n));   
        s(i,j)=h*s(i,j);                   
    end
end
s=ifftshift(s);                           
b15Img=uint8(real(ifft2(s)));        

%20級巴特沃思

n=20;            
s=fftshift(fft2(l));
[M,N]=size(s);                        
d0=30; %BLPF濾波，d0=5，15，30（程式中以d0=30為例）                    
n1=floor(M/2);                          
n2=floor(N/2);                           
for i=1:M 
    for j=1:N
        d=sqrt((i-n1)^2+(j-n2)^2);         
h=1/(1+(d/d0)^(2*n));   
        s(i,j)=h*s(i,j);                   
    end
end
s=ifftshift(s);                           
b20Img=uint8(real(ifft2(s)));    

subplot(1,4,1);
imshow(l);
title('原圖');
subplot(1,4,2);
imshow(b5Img);
title('階數n=5');
subplot(1,4,3);
imshow(b15Img);
title('階數n=15');
subplot(1,4,4);
imshow(b20Img);
title('階數n=20');

%高斯高通濾波器

I = imread('test.jpg'); 
imshow(I);
s=fftshift(fft2(I));
[M,N]=size(s);                    
n=2;                                  
d0=30; %GLPF濾波，d0=5，15，30（程式中以d0=30為例）                    
n1=floor(M/2);                          
n2=floor(N/2);                           
for i=1:M 
    for j=1:N
        d=sqrt((i-n1)^2+(j-n2)^2);         
h=1-1*exp(-1/2*(d^2/d0^2));  
        s(i,j)=h*s(i,j);                   
    end
end
s=ifftshift(s);                           
gImg=uint8(real(ifft2(s)));                                     
figure;                                
imshow(gImg); 

%巴特沃思高通濾波器

I = imread('lena.jpg'); 
imshow(I);
s=fftshift(fft2(I));
[M,N]=size(s);                    
n=2;                                  
d0=30; %BLPF濾波，d0=5，15，30（程式中以d0=30為例）                    
n1=floor(M/2);                          
n2=floor(N/2);                           
for i=1:M 
    for j=1:N
        d=sqrt((i-n1)^2+(j-n2)^2);         
h=1/(1+(d0/d)^(2*n));   
        s(i,j)=h*s(i,j);                   
    end
end
s=ifftshift(s);                           
bImg=uint8(real(ifft2(s)));                                     
figure;                                
imshow(bImg);