Matlabでバンドパスバタワースフィルタを実装するにはどうしたらいいですか？

Question

これは私の予想される出力です。

しかし、私は次の出力を取得しています。

は、私のソースコードの何が問題になっているのですか？

---

ソースコード

main.m

clear_all(); 

I = gray_imread('woman.png'); 
d0 = 1; 
d1 = 100; 
n = 20; 

J = butter_bandpass(I, d0, d1,n); 

J = histeq(J); 

K = {I, J}; 

imshowpair(I, J, 'montage'); 

butter_bandpass.m

function output_image = butter_bandpass(I, dL, dH,n) 
    I_double = double(I); 
    [Height, Width] = size(I_double); 
    I_double = uint8(I_double); 
    I_fft = fft2(I_double,2*Height-1,2*Width-1); 
    I_fft = fftshift(I_fft); 

    LowPass = ones(2*Height-1,2*Width-1); 
    HighPass = ones(2*Height-1,2*Width-1); 
    Kernel = ones(2*Height-1,2*Width-1); 

    for i = 1:2*Height-1 
     for j =1:2*Width-1 
      D = ((i-(Height+1))^2 + (j-(Width+1))^2)^.5;   
      LowPass(i,j)= 1/(1 + (D/dH)^(2*n)); 
      HighPass(i,j)= 1 - (1/(1 + (D/dL)^(2*n))); 

      % Create Butterworth filter. 
      Kernel(i,j) = LowPass(i,j).*HighPass(i,j); 
     end 
    end 

    % apply filter to image 
    output_image = I_fft + Kernel.*I_fft; 

    % revert the image to spatial domain 
    output_image = ifftshift(output_image); 
    output_image = ifft2(output_image,2*Height-1,2*Width-1); 
    output_image = real(output_image(1:Height,1:Width)); 
    output_image = uint8(output_image); 
end 

Answer 1

Iこの問題を解決しました。

出力

this low-pass-filterを参照し、this high-pass-filterソースコード。私たちは次のようにBPF（バンドパスフィルタ）は= LPF * HPF、我々はバンドパスフィルタを実装することができることを知って、以来

、

ソースコード

main.m

clear_all(); 
I = gray_imread('woman.png'); 
Dl = 70; 
Dh = 70; 
n = 1; 
[J, K] = butterworth_bpf(I, Dh, Dl, n); 
J = histeq(J); 
imshowpair(I, J, 'montage'); 

butter_bp_kernel.m

function k = butter_bp_kernel(I, Dh, Dl, n) 
    hp_kernel = butter_hp_kernel(I, Dh, n); 
    lp_kernel = butter_lp_kernel(I, Dl, n);  
    k = hp_kernel.*lp_kernel; % here is the multiplication 

butterworth_bpf.m

function [out1, out2] = butterworth_bpf(I, Dh, Dl, n) 
    Kernel = butter_bp_kernel(I, Dh, Dl, n);  
    I_ffted_shifted = fftshift(fft2(I)); 
    I_ffted_shifted_filtered = I_ffted_shifted.*Kernel; 
    out1 = ifftshow(ifft2(I_ffted_shifted_filtered));  
    out2 = ifft2(ifftshift(I_ffted_shifted_filtered)); 
end