2
$\begingroup$

i'm doing overlap save method in frequency domain. i make hankel function to digital filter using frequency sampling method(to make arbitrary magnitude phase filter). and i do ifft and zero padding fft

i heard that if the filter is non-causal, zero padding fft have gibb's phenomenon. my filter have same problem. i'm doing shift FIR filter to make causal filter but the problem is not solved.... please help me

 clear all; close all; %% parameter fs = 2000; ts = 1/fs; c = 340; R = 1; %% hankel filter Lp = 50; N = 2*Lp+1; %filter length resol = 0:N-1; % frequency resolution f = resol/N*fs; k = 2*pi*f/c; M=5; R = 2; nu = 2; %hankel order % frequency sampling method right = besselh(nu,k(1:(N-1)/2 +1)*R); %Hankel function 0~pi right(1:5) = 0; %eliminate infinite part left = fliplr(right); left = conj(left); left(end) = []; H = [right left]; % symmetric freuqnecy response h = ifft(H); % impulse response nn = 80; %shift index h1 = [h(nn : end) h(1:nn-1)]; % shift impulse response to make causal filter H1 = fft(h1,2*N); %zero padding fft subplot(2,1,1) plot(h1) subplot(2,1,2) plot(abs(H1)) 
$\endgroup$
2
  • $\begingroup$ Your H is not conjugate symmetric so that h is not real, please check your code. To eliminate the ininite part, why do you set right(1:5) to zero? $H_\nu(z)$ goes to infinity only when $z=0$ due to the property of Neumann function. $\endgroup$ Commented Feb 7, 2022 at 4:59
  • $\begingroup$ i edit code Thank you. 1:5 is when frequency resolution is high first some part have very high value. So i do 1:5 $\endgroup$ Commented Feb 7, 2022 at 5:01

1 Answer 1

2
$\begingroup$

The Gibbs phenomenon is caused by frequency sampling method. I recommend you to use window method instead. It is fast, convenient and robust, although not optimal. BTW, the frequency range should be 0~fs/2, I've fixed it in the following code

clear; close all; %% parameter fs = 2000; ts = 1/fs; c = 340; %% hankel filter Lp = 50; N = 2*Lp+1; % filter length NN = 2^nextpow2(N * 100); % NN >> N resol = 0:NN-1; % frequency resolution f = resol/NN*fs/2; % you miss /2 here k = 2*pi*f/c; M = 5; R = 2; nu = 2; %hankel order % window method right = besselh(nu,k(1:(NN-1)/2 +1)*R); %Hankel function 0~pi right(1:200) = 1; %eliminate infinite part left = fliplr(right); left = conj(left); left(end) = []; H = [right left]; % symmetric freuqnecy response h = ifft(H); % impulse response h1 = ifftshift(h); % shift impulse response to make causal filter win = hann(N).'; h2 = h1(NN/2-N/2:NN/2+N/2-1) .* win; % windowed impulse response H2 = fft(h2,2*N); %zero padding fft figure; subplot(2,1,1) plot(h2) subplot(2,1,2) plot(abs(H2)) 

Frequency sampling method:

Frequency Sampling

Window method:

Window

$\endgroup$
7
  • $\begingroup$ OMG............Very Very Very Thank you so much $\endgroup$ Commented Feb 7, 2022 at 5:37
  • $\begingroup$ if error is large, is this normal? $\endgroup$ Commented Feb 7, 2022 at 6:00
  • $\begingroup$ @ggh I don't see significant error between these two methods. $\endgroup$ Commented Feb 7, 2022 at 6:10
  • $\begingroup$ there is some wrong part. i fixed. Thank you very much $\endgroup$ Commented Feb 7, 2022 at 6:33
  • $\begingroup$ Could i ask one more question?? How can i make non-causal filter? that didn't use ifftshift command?? $\endgroup$ Commented Feb 8, 2022 at 7:30

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.