# Calculating the inverse filter for the (exponential) sine sweep Method

I am currently working on my Bachelor-Thesis about Real-Time convolution and Impulse Response Measurements. After reading different papers about the (exponential)sine sweep method I didn't find an answer on how to calculate the inverse filter for deconvolving the impulse response.

In the papers I read, it is described as the Time Reversal Mirror and after that some kind of scaling.

Could someone explain, how to calculate the inverse filter for a given sine sweep? If you could add an example, or an algorithm it would be appreciated

• there is a terminology issue here. a "sine sweep" is a signal, not a filter (unless you mean that the sine sweep is an impulse response of a filter). filters have an inverse filter, but not signals. Jun 13, 2017 at 16:10
• so the papers I read most notably A. Farina's "Simultaneous Measurement of Impulse Response and Distortion with a Swept-Sine Technique" he is talking about creating an inverse filter so that the input signal x(t) convolved with said inverse filter f(t) results in a delayed dirac's delta function. So that the Output signal y(t) convolved with f(t) results in the impulse response h(t) Jun 13, 2017 at 16:38
• i might call that a form of a "matched filter". it is a delayed inverse filter of a filter having a sine sweep as its impulse response. Jun 13, 2017 at 17:10

Assuming that your Exponential Sweep Sine was generated using the formula:

$$x(t)=\sin\left(\frac{2\pi f_1 T}{R}\left(e^{\frac{t R}{T}} -1\right) \right)$$

where:

$f_1, f_2$ - Initial and final frequency of the sweep

$T$ - Duration of the sweep

$R = \ln\left(\frac{f_2}{f_1} \right)$ - Sweep rate

Then the inverse filter is calculated by scaling the amplitude of time reversed $x(t)$ by:

$k = e^{\frac{tR}{T}}$

Which will result in an exponentially decaying sweep:

$f(t) = \frac{x_{inv}(t)}{k}$

Example in Python:

#!/usr/bin/env python

from __future__ import division
import numpy as np
import scipy.signal as sig
import matplotlib.pyplot as plt

def dbfft(x, fs, win=None):
N = len(x)  # Length of input sequence

if win is None:
win = np.ones(x.shape)
if len(x) != len(win):
raise ValueError('Signal and window must be of the same length')
x = x * win

# Calculate real FFT and frequency vector
sp = np.fft.rfft(x)
freq = np.arange((N / 2) + 1) / (float(N) / fs)

# Scale the magnitude of FFT by window and factor of 2,
# because we are using half of FFT spectrum.
s_mag = np.abs(sp) * 2 / np.sum(win)

# Convert to dBFS
ref = s_mag.max()
s_dbfs = 20 * np.log10(s_mag/ref)

return freq, s_dbfs

if __name__ == "__main__":
# Sweep Parameters
f1 = 10
f2 = 100
T = 3
fs = 1000
t = np.arange(0, T*fs)/fs
R = np.log(f2/f1)

# ESS generation
x = np.sin((2*np.pi*f1*T/R)*(np.exp(t*R/T)-1))
# Inverse filter
k = np.exp(t*R/T)
f = x[::-1]/k
# Impulse response
ir = sig.fftconvolve(x, f, mode='same')

# Get spectra of all signals
freq, Xdb = dbfft(x, fs)
freq, Fdb = dbfft(f, fs)
freq, IRdb = dbfft(ir, fs)

plt.figure()
plt.subplot(3,1,1)
plt.grid()
plt.plot(t, x)
plt.title('ESS')
plt.subplot(3,1,2)
plt.grid()
plt.plot(t, f)
plt.title('Inverse filter')
plt.subplot(3,1,3)
plt.grid()
plt.plot(t, ir)
plt.title('Impulse response')

plt.figure()
plt.grid()
plt.semilogx(freq, Xdb, label='ESS')
plt.semilogx(freq, Fdb, label='Inverse filter')
plt.semilogx(freq, IRdb, label='IR')
plt.title('Spectrum')
plt.xlabel('Frequency [Hz]')
plt.ylabel('Amplitude [dBFS]')
plt.legend()

plt.show()


And output:

For Robert, here is the spectrum plot for case without amplitude modulation of the inverse filter:

Relevant literature:

Q. Meng - Impulse Response Measurement with Sine Sweeps and Amplitude Modulation Schemes

A. Novak - Nonlinear System Identification Using Exponential Swept-Sine Signal

K. Vetter - ExpoChirpToolbox - a Pure Data implementation of ESS impulse response measurement

• where does $f_2$ figure into the equations? i might imagine that there is a relationship between $f_1$, $f_2$, $T$, and $R$, but what is it? Jun 13, 2017 at 17:12
• Should be T. Fixed now
– jojek
Jun 13, 2017 at 17:15
• @robertbristow-johnson: Short answer being: otherwise you will get a skewed frequency response after convolving the recording with unscaled time-reversed sweep.
– jojek
Jun 13, 2017 at 17:49
• it's skewed in a deterministic manner. you know exactly how it is skewed. you can compensate for the skewing in the frequency response. but the point of driving an acoustic system with a broadbanded test signal is to get a good representation of frequency response and to have as high S/N ratio as possible without clipping or going significantly non-linear. this is why low crest factor is salient. Jun 13, 2017 at 18:15
• well, you have the diamond, but i think that this is all technically relevant. Jun 13, 2017 at 18:16