Forums

filterbanks problem?

Started by omid...@yahoo.com December 11, 2006
This MATLAB program is for making a suitable wavelet kernel
It derives the Low pass filter by sampling from analog raised cosine function in frequency domain and taking ifft from this function to yield the impulse response. (other filters are derived using orthogonal relations).

first i wnat know if the filter design method has any problems?
The filter has 2 parameters in frequency domain Fs and B. Fs is the number of samples which will be the length of the filter and B which controls the sharpness of the curve in frequency domain.

These filters have fairly large reconstruction error. even when the filter length is large(Fs= 32 and B=.15)
When i take the fft of this filter with large number of fft points (fft(tz,8000)) the result shows small ripples in stopband. I checked the same for db filters and these filters didn't have these ripples.

i'm confused!!!

=======================================M-file
% program for making a suitable kernel for coder
clear all;
close all;

Fs; % sampling frequncy from analog filter T(f)
%B=.1505; % input parameter
B=.15;
r=1/2; % r (second parameter)
k=Fs/2;

for f=1:k+1 % sampling of analog filter in frequency domain

if (f-1)<=(r-B)*k
T(f)=1;
elseif (f-1)<=(r+B)*k
T(f)=cos(pi/(4*B)*((f-1)/k-r+B));
elseif (f-1)<=k
T(f)=0;
end;

end;
x=(0:k)/Fs;
plot (x,T);
for f=(k+2):Fs % making second half of samples
T(f)=T(Fs+2-f);
end;

T1=real(ifft(T)); % original FIR filter

for i=1:k % TL is casual filter
TL(i+k)=T1(i);
end;

for i=1:k
TL(i)=T1(k+i);
end;
tx=TL;

%%$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$%%

[y,fs]=wavread('taha.wav'); % getting input signal(16KHz wav file)

plot(y);
title('Original');
input=y';
%--
for i=1:length(tx) % reporting samples of filter (kernel)
A(i)=fprintf('%.16f ',tx(i));
end;

[Lod, Hid, Lor, Hir]=orthfilt(tx); %making other filters

ql=filterx1(input,Lod); %analyse of input input signal
qh=filterx1(input,Hid);

qld=d_convrt(ql); %down converting
qhd=d_convrt(qh);

qldu=u_convrt(qld); %up converting
qhdu=u_convrt(qhd);

outl=filterx1(qldu,Lor); % synthesis filtering
outh=filterx1(qhdu,Hir);

final_out=outl+outh; % summing output of synthesis filters
hold;
%plot(final_out,'G');
plot(final_out(Fs:Fs-1+length(input)),'g'); % report summary
title('reconstructed & Original Signal');
%pause;
figure;
plot(final_out(Fs:Fs-1+length(input))-input,'g');
title('Error');
yE=final_out(Fs:Fs-1+length(input))-input;
L=length(yE);
E=sqrt(1/L*sum(yE.^2))