## Fast audio convolution

November 19, 2010 Coded in Python

Fast convolution of an audio signal with the impulse response of a system modeled as LTI (Linear Time Invariant) one.

Uses FFT to perform the convolution by multiplication in the frequency domain. To reduce in this way, the computation time.

Works with mono or stereo input signals.

``````#!/usr/bin/python
# -*- coding: utf-8 -*-

#
# Copyright (C) 2006 by Hernán Ordiales
# <audiocode@uint8.com.ar>
#
# This program is free software; you can redistribute it and/or modify
# the Free Software Foundation; either version 2 of the License, or
# (at your option) any later version.
#

# Fast convolution of an audio signal with the impulse response of a system modeled as LTI (Linear Time Invariant) one.
# Use: fast_conv.py source.wav impulse_response.wav output.wav

from wav_array import *
from FFT import fft, inverse_fft
from pylab import size
import sys

def nextpow2( L ):
N = 2
while N < L: N = N * 2
return N

def fast_conv_vect( x, h ):
# searches for the amount of points required to perform the FFT
L = size(h) + size(x) - 1 # linear convolution length
N = nextpow2( L )
# Note: N>=L is needed because the IDFT of the multiplication is the circular convolution and to match it to the
# common one, N>=L is required (where L=N1+N2-1;N1=length(x);N2=length(h))

# FFT(X,N) is the N points FFT, with zero padding if X has less than N points and truncated if has more.
H = fft( h, N ) # Fourier transform of the impulse
X = fft( x, N ) # Fourier transform of the input signal

Y = H * X # spectral multiplication
y = inverse_fft( Y ) # time domain again
return y

source = sys.argv
impulse = sys.argv
conv_output = sys.argv

clip_factor = 1.01 # clip factor default value

[ h1, Fs1, h_bits ] = wavread( impulse ) # impulse response
[ x1, Fs2, x_bits ] = wavread( source ) # file to process

if Fs1 == Fs2 : # if sample rates are the same
print "Processing..."
y1 = fast_conv_vect( x1, h1 ).real # takes the real part to avoid a too small complex part (around e-18)

# audio normalization: if "y = y/max(y)" -> "values out of [-1,+1] range are clipped"
y1 = y1/( max(y1)*clip_factor ) # to avoid clipping

wavwrite( y1, Fs1, conv_output )

print "Output file:", conv_output
print "Convolution success."
else:
print "Error: files has different sample rate."``````