在頻率域中直接生成濾波器

除了之前說的從空間濾波器中獲得頻率域濾波器，還可以從頻率域中直接生成濾波器，這些濾波器被規定為距濾波器中心點的距離不同的函數。可以創建一個用於實現頻率濾波器的網格數組，最主要的是需要計算任何點到頻率矩形中一個指定點的距離函數，FFT（快速傅里葉）算法是假設變換的原點位於頻率矩形的左上角，因此需要將原點平移到頻率矩形的中心，用fftshift。網格數組如下：

%(頻域濾波函數) 提供了距離計算及其所需的網格數組
function [U,V] = dftuv(M,N)
u=0:(M-1);
v=0:(N-1);
idx = find(u>M/2);
u(idx) = u(idx)-M;
idy=find(v>N/2);
v(idy)=v(idy)- N;
[V,U] = meshgrid(v,u);
end

 

一、低通（平滑）頻率域濾波器

常見的頻率域低通濾波器有三個，理想低通濾波器（ILPF），巴特沃斯低通濾波器（BLPF），高斯低通濾波器（GLPF）。

1）理想低通濾波器的傳遞函數如下：

其中，D0為正數，D(u,v)為點（u,v）到濾波器中心的距離，滿足D(u,v)=D0的點的軌跡為一個圓。如果用濾波器乘以一幅圖像的傅里葉變換，我們會看到一個理想濾波器會切斷（乘以0）該圓之外的所有F(u,v)分量，而保留圓上和圓內的所有分量不變（乘以1）。

可以用之前的mesh將濾波器顯示出來：

             

2）n階巴特沃斯低通濾波器，在距離濾波器中心D0處具有截止頻率，傳遞函數為：

與理想低通濾波器不同的是，巴特沃斯低通濾波器的傳遞函數在D0點並沒有一個尖銳的不連續，對於具有平滑傳遞函數的濾波器，通常將截止頻率軌跡定義在H(u,v)降低為其最大值的一個指定的比例的點處。用mesh同樣可以將該濾波器顯示出來：

                  

3）高斯低通濾波器的傳遞函數如下：

其中，為標准差。該濾波器可以用mesh顯示出來

                   

4）低通濾波器的例子

用高斯低通濾波器對一幅500*500的圖像進行濾波，效果如下：

該例子的代碼如下：

%用高斯低通濾波器進行濾波
f=imread('G:\數字圖像處理（岡薩雷斯）\DIP3E_CH04_Original_Images\DIP3E_Original_Images_CH04\Fig0441(a)(characters_test_pattern).tif');
subplot(221);imshow(f); title('原圖')
f = im2double(f);
%[f,revertclass] = tofloat(f);
PQ = paddedsize(size(f));
[U,V] = dftuv(PQ(1),PQ(2));
D = hypot(U,V);
D0 = 0.05*PQ(2);
F = fft2(f,PQ(1),PQ(2));
H = exp(-(D.^2)/(2*(D0^2)));
g = dftfilt(f,H);
g = im2uint8(g);
%g = revertclass(g);
subplot(222);imshow(fftshift(H));  title('高斯低通濾波器圖像'); %顯示濾波器圖像
subplot(223);imshow(log(1+abs(fftshift(F))),[]);title('濾波器的譜');
subplot(224);imshow(g);title('濾波器后的圖像')

　各種濾波器的傳遞函數不同，但是濾波的過程是一樣的，因此，可以將這幾種濾波器封裝成一個函數，如下：

%頻率域低通濾波函數,生成幾個低通濾波器的傳遞函數
function [H,D] = lpfilter(type,M,N,D0,n)
[U,V] = dftuv(M,N);
D = sqrt(U.^2+V.^2);
switch type
    case 'ideal'
        H=double(D<=D0);
    case 'btw'
        if nargin ==4
            n=1;
        end
        H=1./(1+(D./D0).^(2*n));
    case 'gaussian'
        H=exp(-(D.^2)./(2*(D0^2)));
    otherwise
        error('unkown filter type');
        
end
        
        

二、高通（銳化）頻率域濾波器　

就像低通濾波模糊一幅圖像那樣，高通濾波是相反的過程，會銳化圖像，其原理是衰減傅里葉變換的低頻部分而保持高頻部分相對不變。若給定低通濾波器的傳遞函數Hlp(u,v)，則相應高通濾波器的傳遞函數為：

與低通濾波器相對應的高通濾波器的傳遞函數如下：

     

基於前面的公式，可以使用前面的函數lpfilter來構建一個生成高通濾波器的函數，如下所示：

%生成高通濾波器的函數
function H = hpfilter(type,M,N,D0,n)
if nargin ==4
    n=1;
end
Hlp = lpfilter(type,M,N,D0,n);
H = 1-Hlp;
end

 

下圖顯示了三個高通濾波器的透視圖和圖像，

代碼也比較簡單，

%高通濾波器的透視圖和圖像
 
H1 = fftshift(hpfilter('ideal',500,500,50)); %理想高通濾波器
H2 = fftshift(hpfilter('btw',500,500,50));   %巴特沃斯高通濾波器
H3 = fftshift(hpfilter('gaussian',500,500,50));  %高斯高通濾波器
subplot(231);mesh(double(H1(1:10:500,1:10:500)));title('理想高通濾波器') ;         %線框圖
colormap([0 0 0]);               %黑色
axis off
     
subplot(232);mesh(double(H2(1:10:500,1:10:500)));title('巴特沃斯高通濾波器') ;         %線框圖
colormap([0 0 0]);               %黑色
axis off
    
subplot(233);mesh(double(H3(1:10:500,1:10:500)));title('高斯高通濾波器') ;         %線框圖
colormap([0 0 0]);               %黑色
axis off

subplot(234);imshow(H1,[]); title('理想高通濾波器');        %高通濾波器對應的圖像
subplot(235);imshow(H2,[]); title('巴特沃斯高通濾波器'); 
subplot(236);imshow(H3,[]); title('高斯高通濾波器'); 

 

1）高通濾波例子

高通濾波時，圖像的邊緣和其他灰度變化劇烈的地方得到了增強，但由於圖像的平均值由F(0,0)給出，因此，高通濾波器偏離了傅里葉變換的原點，圖像會失去大部分原始圖像的亮度，變得比較黑。

代碼如下：

%高斯高通濾波   
f=imread('G:\數字圖像處理（岡薩雷斯）\DIP3E_CH04_Original_Images\DIP3E_Original_Images_CH04\Fig0441(a)(characters_test_pattern).tif');
PQ = paddedsize(size(f));   %補0填充
D0 = 0.05*PQ(1);    %截止頻率
H1 = hpfilter('ideal',PQ(1),PQ(2),D0);   %理想高通濾波器
H2 = hpfilter('btw',PQ(1),PQ(2),D0);   %巴特沃斯高通濾波器
H3 = hpfilter('gaussian',PQ(1),PQ(2),D0);   %高斯高通濾波器

%H = 1+1.5*H;
g1 = dftfilt(f,H1);     %濾波
g2 = dftfilt(f,H2);     
g3 = dftfilt(f,H3);     
subplot(221);imshow(f);title('原圖');
%g = g+f;
subplot(222);imshow(g1);title('理想高通濾波');
subplot(223);imshow(g2);title('巴特沃斯高通濾波');
subplot(224);imshow(g3);title('高斯高通濾波');

 

高通濾波器偏離了直流項，因此將圖像的平均值降低為0 。解決的方法之一是給高通濾波器加一個偏移量，如果把偏移量與將濾波器乘以一個大於1的常數結合起來，那么這種方法就稱為高頻強調濾波。高頻強調濾波比較麻煩的地方就是ab參數需要根據實際情況進行不斷嘗試才會發現比較合適的參數。其傳遞函數為：

a是偏移量，b是乘數。

經過高頻強調濾波后的結果如下，將原來的高通濾波器都改成相應的高頻強調濾波器（a=1,b=1.5）：

最明顯的就是理想高通濾波器有振鈴效應，但是邊緣都被增強了。除了用高頻強調濾波，我發現直接在高通濾波后的圖像上疊加原圖也可以達到不錯的效果，因為高通濾波后的圖像均值接近0，所以疊加了原圖像可以增強亮度和細節，而且振鈴效應也沒有那么強烈。如下所示：

2）結合使用高頻強調濾波和直方圖均衡

例子源於數字圖像處理例3.8，增強一幅數字的胸部X射線圖像。

上述操作的代碼如下：

%聯合使用高頻強調濾波和直方圖均衡

f=imread('G:\數字圖像處理（岡薩雷斯）\DIP3E_CH04_Original_Images\DIP3E_Original_Images_CH04\Fig0459(a)(orig_chest_xray).tif');
PQ = paddedsize(size(f));      %獲取填充范圍
D0 = 0.05*PQ(1);              %截止頻率
HBW = hpfilter('btw',PQ(1),PQ(2),D0,2);       %二階巴特沃斯高通濾波器
H = 0.5+2*HBW;               %高頻強調濾波,比直接的高通濾波准確
gbw = dftfilt(f,HBW,'fltpoint');            
gbw = gscale(gbw);

ghf = dftfilt(f,H,'fltpoint');            %濾波
ghf = gscale(ghf);           %該函數考慮的負值，可以保留細節
ghe = histeq(ghf,256);       %直方圖均衡化(由較窄灰度級范圍內的灰度表征的圖像是直方圖均衡的理想選擇)
subplot(221);    imshow(f);title('原圖');
subplot(222);    imshow(gbw);title('高通濾波后');
subplot(223);    imshow(ghf);title('高頻強調濾波后');
subplot(224);    imshow(ghe,[]);title('高頻后再直方圖均衡化后');
%  figure;imhist(ghf);