DSPRelated.com
Forums

Bandfilters with FFT on Sampled in-/output

Started by basgoossen June 23, 2011
Hello all,

Currently i am working on echo cancellation software for a voice
communication system. and i came across the following problem:

I need to analyse and filter frequency's, so i do a FFT i cut off some
frequency's and then perform an iFFT to get the filtered signal. This
happens on sample's varying between 20ms and 1000ms (2^9 to 2^16). The
frequency's that i whant to get rid of are filtered quite nicely, however i
hear a pop/plop at the end/beginning of each sample.

My current guess is that the pop/plop is introduces by a phase shift caused
by the "filter" causing the end of the one sample not to allign to the
beginning of the new one...

My guess could be wrong i'm not an expert in this field. My main question
is how to get rid of those plops?

For my knowledge in this area is limited it might also help if someone has
or knows a link on the net with an example on how to filter frequency's in
the frequency domain with sampled input data...

my current approach is to attenuate the undesired frequency's using by
multiplying its real and imaginary value with the inverse of a blackman
window to soften the edges of the filter, however i cannot justify this
method with theory found on the net...

all help is welcome.

Thanks in advance,
Bas



Are you just zeroing out the frequencies you want to filter?  If so, then
that is your problem.  Brick-wall filters in the frequency-domain are
essentially "infinite" bandwidth which translate to ringing in the
time-domain.  There's also an issue with doing consecutive block IFFTs then
converting back to the time domain - there will be start/stop transients at
the end of each block, which, coupled with the filtering issue, will likely
cause the pops you are getting.

Look up overlap/save and overlap/add (I don't use either so I can't be of
much help) for the IFFT issue.  For the filter issue, figure out what your
time-domain impulse response should be then compute the FFT of the filter
(zero-padding to get to the length of your IFFT) and use that as a
multiplier in the frequency-domain.

Mark
>Are you just zeroing out the frequencies you want to filter?
No i use an inverse hamming window and multiply this with the frequency value's at the frequency bands i whant to filter. So a gradualy lower the frequency's around the center frequency. I've read a lot about this overlap-add and overlap save, i just can't seem to get it to work... i have this example (also from this forum): N = 1024; n = (1:N)'-1; %'# define the time index x = sin(2*pi*1.5*n/N); %# input with 1.5 cycles per 1024 points h = hanning(129) .* sinc(0.25*(-64:1:64)'); %'# windowed sinc LPF, Fc = pi/4 h = [h./sum(h)]; %# normalize DC gain Nproc = 512; xproc = zeros(2*Nproc,1); %# initialize temp buffer idx = 1:Nproc; %# initialize half-buffer index ycorrect = zeros(2*Nproc,1); %# initialize destination for ctr = 1:(length(x)/Nproc) %# iterate over x 512 points at a time xproc(1:Nproc) = xproc((Nproc+1):end); %# shift 2nd half of last iteration to 1st half of this iteration xproc((Nproc+1):end) = x(idx); %# fill 2nd half of this iteration with new data yproc = ifft(fft(xproc) .* fft(h,2*Nproc)); %# calculate new buffer ycorrect(idx) = real(yproc((Nproc+1):end)); %# keep 2nd half of new buffer idx = idx + Nproc; %# step half-buffer index end but for someone that has no experience with matlab it's not that self explenatory. i tryed to implement it in java but now all i get is very faint sound mostly plops and klicks and nothing more, so it is further off from what i want. Could somebody try to explain in laymans terms how this works?
>>Are you just zeroing out the frequencies you want to filter? >No i use an inverse hamming window and multiply this with the frequency >value's at the frequency bands i whant to filter. So a gradualy lower the >frequency's around the center frequency.
OK, so that shouldn't be much of a problem.
>I've read a lot about this overlap-add and overlap save, i just can't
seem
>to get it to work...
Sorry, I can't be of much help. Graychip (I think) used to have a good tutorial on both, though since they were bought out by TI I don't know if it still exists. Either way, it may be worthwhile searching TI's website to see if it does. My guess now is that you're getting edge effects from stacked IFFTs. Note, too, that the filter method you are employing still imparts the equivalent filter delay on the start of the signal. Overlapping removes this as well as far as I know. Mark
On Jun 23, 12:55&#4294967295;pm, "basgoossen" <bas@n_o_s_p_a_m.2dive.nl> wrote:
> >Are you just zeroing out the frequencies you want to filter? > > No i use an inverse hamming window and multiply this with the frequency > value's at the frequency bands i whant to filter. So a gradualy lower the > frequency's around the center frequency. > > I've read a lot about this overlap-add and overlap save, i just can't seem > to get it to work... > > i have this example (also from this forum): > N = 1024; > n = (1:N)'-1; %'# define the time index > x = sin(2*pi*1.5*n/N); %# input with 1.5 cycles per 1024 points > h = hanning(129) .* sinc(0.25*(-64:1:64)'); %'# windowed sinc LPF, Fc = > pi/4 > h = [h./sum(h)]; %# normalize DC gain > > Nproc = 512; > xproc = zeros(2*Nproc,1); %# initialize temp buffer > idx = 1:Nproc; %# initialize half-buffer index > ycorrect = zeros(2*Nproc,1); %# initialize destination > for ctr = 1:(length(x)/Nproc) %# iterate over x 512 points at a time > &#4294967295; &#4294967295; xproc(1:Nproc) = xproc((Nproc+1):end); %# shift 2nd half of last > iteration to 1st half of this iteration > &#4294967295; &#4294967295; xproc((Nproc+1):end) = x(idx); %# fill 2nd half of this iteration with > new data > &#4294967295; &#4294967295; yproc = ifft(fft(xproc) .* fft(h,2*Nproc)); %# calculate new buffer > &#4294967295; &#4294967295; ycorrect(idx) = real(yproc((Nproc+1):end)); %# keep 2nd half of new > buffer > &#4294967295; &#4294967295; idx = idx + Nproc; %# step half-buffer index > end > > but for someone that has no experience with matlab it's not that self > explenatory. i tryed to implement it in java but now all i get is very > faint sound mostly plops and klicks and nothing more, so it is further off > from what i want. Could somebody try to explain in laymans terms how this > works?
Instead of implementing the OLS or OLA functionality manually, you could just use MATLAB's filter() function, providing the "h" vector above as your filter coefficients. Jason
>Instead of implementing the OLS or OLA functionality manually, you >could just use MATLAB's filter() function, providing the "h" vector >above as your filter coefficients. > >Jason
That's always true, unless there's a specific reason for doing the filtering via the frequency domain, such as very large vectors, or a homework assignment, hehe. Mark
There are some problems with your code. I suggest you try this.
You have a scenario where you are filtering a tone with a static low-
pass filter.
In this case you can create a reference quite easily
yref = filter(h,1,x);
Now you can compare the output of your frequency-domain processing
routine with this reference.

There are 3 things to be aware of

1. As already suggested you need to read about OLA and/or OLS.
2. You need a larger overlap.
3. You need to make sure that your circular convolution as implemented
by the multiplication of the fourier transforms either corrresponds to
a linear convolution or at least gets close to it.
Cheers
Thanks everybody for the answers,

The matlab code mentionned is not something i use but an example of
overlap-save progressing found on this forum. I implement the source in
Java so i can't use the matlab's filter function ;-).

For now i have got rid of 99% of the filter artifacts, my current approach
is to use a buffer size that is twice the requested buffer size and use
only the middle portion of the processed buffer to write the output. Since
the chunks are only 20ms that adds a delay of roughly 5ms, which is
acceptable in the application.

I've also tested whith brick-wall filters as named this way earlyer in this
topic. Those now also give quite satisfying results, Only very high amp
noise on the endge of the brick ;-) give artifacts, but for normal band
filtering it seems to work.

At this moment i am not using windows to prevent spectral leakage, but it
seems not to have a bad effect on these chunks. When i draw the spectrum i
see sharp peaks without leakage... I don't know why but it just seems to
work.

If anyone wants to warn me of possible side effects of this filter please
let me know. For now i think the preformance is good enough for the
application (Sip communication).
basgoossen

Use oversampled filter banks for frequency domain filtering and avoid
aliasing artifacts.

For example

http://electronix.ru/forum/index.php?s=&showtopic=23652&view=findpost&p=929325
On Jun 23, 3:59&#4294967295;pm, niarn <niar...@gmail.com> wrote:
> > There are 3 things to be aware of > > 1. As already suggested you need to read about OLA and/or OLS. > 2. You need a larger overlap. > 3. You need to make sure that your circular convolution as implemented > by the multiplication of the fourier transforms either corrresponds to > a linear convolution or at least gets close to it.
fundamentally, if Bas is using OLA or OLS, he/she needs to make sure that the H[k] sequence that is used in the frequency domain, of length N, is the DFT of an h[n] that is a *finite* impulse response of length L and zero padded with N-L zeros at the end. for each frame, Bas gets only N-L+1 valid output samples, so the overlap must be L-1. and no windowing (other than the built-in rectangular window) is used for OLS or for this kind of OLA. i would suggest reading O&S. r b-j