Applying a filter on an audio signal with Python

Question

Let's say I have a filter described by its transfer function:

$$ H(\omega) = \frac{1}{1 + j\frac{\omega}{\omega_0}} $$

And I want to apply this filter to an audio signal (a .wav file) using Python. My initial idea was this:

• Split the signal into fixed-size buffers of ~5000 samples each
• For each buffer, compute its Fourier transform using numpy.fft.rfft
• Apply my filter to the coefficients of the Fourier transform: ft[i] *= H(freq[i])
• Compute the inverse Fourier transform of the result using numpy.fft.irfft
• Join the buffers together into the output file

It seems like I am doing something wrong, because if I set H to a constant factor ($\neq 1$), then the sound gets completely noisy and unrecognizable. I also tried with different values for the buffer's size.

Here is my code:

import wave
import numpy, numpy.fft as fft
from math import *

inputFile = wave.open('input.wav', 'r')
outputFile = wave.open('output.wav', 'w')

params = inputFile.getparams()
sampleRate = params.framerate
bufferSize = 5000

outputFile.setparams(params)
samples = []

for i in range(params.nframes//bufferSize):
    signal = inputFile.readframes(bufferSize)
    ft = fft.rfft(signal)
    for i in range(1, ft.shape[0]):
        ft[i] *= 0.5 # This does not do what I want...
    ift = fft.irfft(ft)
    samples += list(map(lambda x: max(min(int(x), 255), 0), numpy.abs(ift)))

outputFile.setnframes(len(samples)//params.nchannels)
outputFile.writeframes(bytes(samples))

Thanks for your time and patience :)

Answer 1

Apart from what Marcus suggests in the comments, ft[i] *= H(freq[i]) is probably wrong. The overlap-add and overlap-save methods should deal with this correctly, but read on!

The FFT does not implement linear convolution. It implements circular convolution, which can approximate linear convolution provided certain conditions are met. One of those conditions is that $$N_{\tt fft} \ge N_{\tt filter} + N_{\tt signal} - 1.$$

In your case, you've effectively set $$N_{\tt fft} = N_{\tt filter} = N_{\tt signal}$$ which will introduce time-aliasing --- the problem that happens when the "certain conditions" are not met. See here for an answer to a similar question.