DSPRelated.com
Forums

Lowpass filter in C

Started by dizige March 20, 2007
Hi,

I'm new to DSP and just wrote a bunch of functions in C to process
some audio data.

I implemented a lowpass filter by convolving a 256-sample audio
sequence with a sinc function from -64 to 64 (129 samples)

To convolve I do the fft of both functions and multiply their ffts,
and then do the inverse fft of the product, and discard extra samples
to have only 256 samples back.

However, the end result sounds awful. I can make out the original
signal, but the quality is really bad.

Any ideas of what I'm doing wrong??

Thanks
Diego
Hello dizige
On 2007-03-20, you wrote:

> To convolve I do the fft of both functions and multiply their ffts,
> and then do the inverse fft of the product, and discard extra samples
> to have only 256 samples back.

> However, the end result sounds awful. I can make out the original
> signal, but the quality is really bad.

Perform the convolution in time domain, not frequency domain.

--
Grzegorz Kraszewski
Hi Diego,

There are many important factors I think you are overlooking:

1) Are you padding the input sequences before doing the FFT to avoid a
circular convolution from happening?
2) Are you doing a true complex multiply in the frequency domain?
3) Are you combining the output segments correctly to obtain the final
filtered audio signal?
4) Is your sinc function windowed, or just a truncated sinc function? A
windowed sinc is prefered.

A great explaination of "FFT Convolution" is given in Chapter 18 of Steven
W. Smith's book,
"The Scientist and Engineer's Guide to Digital Signal Processing" which is
available free at :

http://www.dspguide.com/pdfbook.htm

I would read this chapter, and I think it will solve your problems.

Tony

----- Original Message -----
From: "dizige"
To:
Sent: Tuesday, March 20, 2007 12:21 AM
Subject: [audiodsp] Lowpass filter in C
> Hi,
>
> I'm new to DSP and just wrote a bunch of functions in C to process
> some audio data.
>
> I implemented a lowpass filter by convolving a 256-sample audio
> sequence with a sinc function from -64 to 64 (129 samples)
>
> To convolve I do the fft of both functions and multiply their ffts,
> and then do the inverse fft of the product, and discard extra samples
> to have only 256 samples back.
>
> However, the end result sounds awful. I can make out the original
> signal, but the quality is really bad.
>
> Any ideas of what I'm doing wrong??
>
> Thanks
> Diego
>
Here is some SciLab code for doing FIR filter by FFT

nby2fftsize = 2048;

// Read in the FIR filter (MUST have odd number of coefficients)
coef=fload('FIRby2e.dat');

// pack filter out to fftsize
nfilt = max(size(coef));
coef = [coef,zeros(1:nby2fftsize-nfilt)];
Y = fft(coef,-1);

// open the input file
[hi]=wavefilereadheader('5k_4448.wav');
[ho]=wavefilewriteheader(hi,'test.wav');

// work out step size
nby2samples = nby2fftsize - (nfilt-1);
nsamples = nby2samples/2;
nby2overlap = zeros(1,nfilt-1);

zeropad = zeros(1, (nfilt-1)/2);
by2zerobuffer = zeros(1,nby2fftsize/2);

loop = 100;
while loop
loop = loop - 1;

// read some samples
d = wavefileread(hi,nsamples);
// pad sample to make room for the convolution
// and insert some zeros
h = [d,zeropad;by2zerobuffer];
// repack so that every second sample is zero
h = matrix(h,1,nby2fftsize);
// forward FFT of padded samples
H = fft(h,-1);
// convolution
C = H.*Y;
// inverse FFT
s = real( fft(C,1) );
// perform overlap
s(1:nfilt-1) = s(1:nfilt-1) + nby2overlap;
// save new overlap
nby2overlap(1:nfilt-1) = s(nby2samples+1:nby2fftsize);
// write out filtered samples
[ho]=wavefilewrite(ho,s(1:nby2samples));
end;

wavefileclose(hi);
wavefileclose(ho);
Grzegorz Kraszewski wrote:

>Hello dizige
>On 2007-03-20, you wrote:
>
>
>
>>To convolve I do the fft of both functions and multiply their ffts,
>>and then do the inverse fft of the product, and discard extra samples
>>to have only 256 samples back.
>>
>>
>
>>However, the end result sounds awful. I can make out the original
>>signal, but the quality is really bad.
>>
>>Perform the convolution in time domain, not frequency domain.
>
>
>
Diego-

Tony's reply is great, it should help you a lot. Here is another book that has
excellent treatment of fast convolution implemented with FFTs:

The Fast Fourier Transform and its Applications
E. Oran Brigham
ISBN 0-13-307505-2

Also I might add that 256 pts is not a lot to take advantage of freq domain speed
increase. 512 pts or less might be faster if implemented as basic time domain
convolution; depends on your system and host processor.

-Jeff
"Tony (DSPG)" wrote:
>
> Hi Diego,
>
> There are many important factors I think you are overlooking:
>
> 1) Are you padding the input sequences before doing the FFT to avoid a
> circular convolution from happening?
> 2) Are you doing a true complex multiply in the frequency domain?
> 3) Are you combining the output segments correctly to obtain the final
> filtered audio signal?
> 4) Is your sinc function windowed, or just a truncated sinc function? A
> windowed sinc is prefered.
>
> A great explaination of "FFT Convolution" is given in Chapter 18 of Steven
> W. Smith's book,
> "The Scientist and Engineer's Guide to Digital Signal Processing" which is
> available free at :
>
> http://www.dspguide.com/pdfbook.htm
>
> I would read this chapter, and I think it will solve your problems.
>
> Tony
>
> ----- Original Message -----
> From: "dizige"
> To:
> Sent: Tuesday, March 20, 2007 12:21 AM
> Subject: [audiodsp] Lowpass filter in C
>
> > Hi,
> >
> > I'm new to DSP and just wrote a bunch of functions in C to process
> > some audio data.
> >
> > I implemented a lowpass filter by convolving a 256-sample audio
> > sequence with a sinc function from -64 to 64 (129 samples)
> >
> > To convolve I do the fft of both functions and multiply their ffts,
> > and then do the inverse fft of the product, and discard extra samples
> > to have only 256 samples back.
> >
> > However, the end result sounds awful. I can make out the original
> > signal, but the quality is really bad.
> >
> > Any ideas of what I'm doing wrong??
> >
> > Thanks
> > Diego
Jeff has a good point. Based on block size, there is
a point where FFT convolution becomes faster. (Longer
block sizes favor FFT convolution). I'm not sure exactly where point this
is, though.
If you can do the convolution in the time
domain, you'd be much better off. It is much, much easier. It only
takes a handful of lines of C code. It will save you a lot
of debugging time.

Tony

----- Original Message -----
From: "Jeff Brower"
To: "Diego Z"
Cc:
Sent: Wednesday, March 21, 2007 11:01 PM
Subject: Re: [audiodsp] Lowpass filter in C
> Diego-
>
> Tony's reply is great, it should help you a lot. Here is another book
> that has
> excellent treatment of fast convolution implemented with FFTs:
>
> The Fast Fourier Transform and its Applications
> E. Oran Brigham
> ISBN 0-13-307505-2
>
> Also I might add that 256 pts is not a lot to take advantage of freq
> domain speed
> increase. 512 pts or less might be faster if implemented as basic time
> domain
> convolution; depends on your system and host processor.
>
> -Jeff
> "Tony (DSPG)" wrote:
>>
>> Hi Diego,
>>
>> There are many important factors I think you are overlooking:
>>
>> 1) Are you padding the input sequences before doing the FFT to avoid a
>> circular convolution from happening?
>> 2) Are you doing a true complex multiply in the frequency domain?
>> 3) Are you combining the output segments correctly to obtain the final
>> filtered audio signal?
>> 4) Is your sinc function windowed, or just a truncated sinc function? A
>> windowed sinc is prefered.
>>
>> A great explaination of "FFT Convolution" is given in Chapter 18 of
>> Steven
>> W. Smith's book,
>> "The Scientist and Engineer's Guide to Digital Signal Processing" which
>> is
>> available free at :
>>
>> http://www.dspguide.com/pdfbook.htm
>>
>> I would read this chapter, and I think it will solve your problems.
>>
>> Tony
>>
>> ----- Original Message -----
>> From: "dizige"
>> To:
>> Sent: Tuesday, March 20, 2007 12:21 AM
>> Subject: [audiodsp] Lowpass filter in C
>>
>> > Hi,
>> >
>> > I'm new to DSP and just wrote a bunch of functions in C to process
>> > some audio data.
>> >
>> > I implemented a lowpass filter by convolving a 256-sample audio
>> > sequence with a sinc function from -64 to 64 (129 samples)
>> >
>> > To convolve I do the fft of both functions and multiply their ffts,
>> > and then do the inverse fft of the product, and discard extra samples
>> > to have only 256 samples back.
>> >
>> > However, the end result sounds awful. I can make out the original
>> > signal, but the quality is really bad.
>> >
>> > Any ideas of what I'm doing wrong??
>> >
>> > Thanks
>> > Diego
>