# Exponential/Logarithmic Sine Sweep - Inverse Filter Amplitude Correction

I am implementing the exponential/logarithmic sine sweep approach to impulse response measurement (as presented by Farina - http://pcfarina.eng.unipr.it/Public/Papers/134-AES00.PDF).

My query is regarding amplitude correction of the inverse sweep. I have found a lot of differing information about what correction should be applied, and the only method that seems to work as I expect, doesn't have any information on how it was derived.

So far, I have generated the sweep (in Python) as follows:

f1 = 20     # Start frequency
f2 = 20000  # Stop frequency
T = 1       # Duration in seconds
Fs = 48000  # Sample rate

L = T/np.log(f2/f1)      # Constant used in sweep calculation

t = np.arange(0,T,1/Fs)  # Time ramp

sweep = np.sin(2*np.pi*f1*L*(np.exp(t/L) - 1))


Then to generate the inverse sweep, you have to time reverse the sweep and then compensate for the fact that, due to the non-linear variation in frequency, the sweep contains less energy at higher frequencies. In terms of the exact compensation to apply, I have found varying information.

In Farina's paper, he says to generate the inverse sweep, you must time reverse the sweep and then apply an amplitude correction that starts at 0 dB and ends at $$-6\log{}_2(\frac{f2}{f1})$$ dB.

In their work on "Synchronized Swept-Sine", Novak et al. (https://hal.science/hal-02504286/document) give the amplitude correction as $$\frac{f1}{L}\exp{}(-\frac{t}{L})$$.

In the answer to a question on calculating the inverse of the exponential sine sweep on here (Calculating the inverse filter for the (exponential) sine sweep Method) jojeck gives the correction as $$\frac{1}{\exp(\frac{tR}{T})}$$ where $$R = \log(\frac{f2}{f1})$$.

Finally, in another question on here (Impulse response amplitude : Sine sweep method), the OP, bouaaah, gives the correction as (np.exp(-time/L))/L*f2*T**2 but doesn't say where the equation comes from.

Now, I should say at this point that what I am aiming for is that, when the sweep is deconvolved directly with its inverse, the magnitude of the resulting frequency response is 0 dB in the frequency range of interest, which I think is a reasonable thing to expect, but am open to being corrected.

Here is the code to apply all the above corrections and compare results:

rev_sweep = np.flip(sweep)

env = np.linspace(0,-6*(np.log(f2/f1)/np.log(2)),rev_sweep.size)
env = 10**(env/20)

inv_sweep = rev_sweep*env

N = sweep.size + inv_sweep.size

sweepFFT = fft.fft(sweep, N) * 2 / sweep.size
invFFT = fft.fft(inv_sweep, N) * 2 / inv_sweep.size

H = sweepFFT*invFFT
H = H[0:int(H.size/2 + 1)]
f = np.linspace(0,Fs/2,H.size)

fig, ax = plt.subplots()
fig.set_tight_layout(True)
ax.semilogx(f,20*np.log10(np.abs(H)),label='Farina')
ax.grid(True)

env = (f1/L)*np.exp(-t/L)

inv_sweep = rev_sweep*env

N = sweep.size + inv_sweep.size

sweepFFT = fft.fft(sweep, N) * 2 / sweep.size
invFFT = fft.fft(inv_sweep, N) * 2 / inv_sweep.size

H = sweepFFT*invFFT
H = H[0:int(H.size/2 + 1)]
f = np.linspace(0,Fs/2,H.size)

ax.semilogx(f,20*np.log10(np.abs(H)),label='Novak et al.')

R = np.log(f2/f1)
env = 1/np.exp((t*R)/T)

inv_sweep = rev_sweep*env

N = sweep.size + inv_sweep.size

sweepFFT = fft.fft(sweep, N) * 2 / sweep.size
invFFT = fft.fft(inv_sweep, N) * 2 / inv_sweep.size

H = sweepFFT*invFFT
H = H[0:int(H.size/2 + 1)]
f = np.linspace(0,Fs/2,H.size)

ax.semilogx(f,20*np.log10(np.abs(H)),label='jojeck')

env = np.exp(-t/L)/L*f2*T**2

inv_sweep = rev_sweep*env

N = sweep.size + inv_sweep.size

sweepFFT = fft.fft(sweep, N) * 2 / sweep.size
invFFT = fft.fft(inv_sweep, N) * 2 / inv_sweep.size

H = sweepFFT*invFFT
H = H[0:int(H.size/2 + 1)]
f = np.linspace(0,Fs/2,H.size)

ax.semilogx(f,20*np.log10(np.abs(H)),label='bouaaah')
ax.legend()

plt.show()


And this results in:

As you can see, the only correction that gives 0 dB in the passband is the one in bouaaah's question (Impulse response amplitude : Sine sweep method). It should be noted that Farina and jojeck's approaches give very similar results so it's difficult to see the separate lines.

So my question is twofold:

1. Firstly, is it reasonable to expect that when you deconvolve the sweep directly with its inverse, you should get 0 dB in the passband?
2. Secondly, if so, can anyone explain where bouaaah's equation comes from? Have I missed any references etc.?
• đ¤Ż looks like there is a bug in my code! Thank you for spotting it! The resulting IR should be a delta function with an amplitude of 1 (or -3 dBFS RMS across in frequency domain). I will go back and fix my math/code.
– jojeck
Aug 7, 2023 at 13:40
• This paper gives another correction, which is different from bouaaah's method but is correct as well. You can check it out. Aug 16, 2023 at 2:41
• A bit late to the party. I am struggling a bit to understand what is going on here. It is an exciting topic and Iâd very much like to see it solved as I am also very interested in it. I tried to replicate the results and got the behaviour you described. On the other hand, the impulse returned by the method has a humongous peak, which scales the estimated/measured system response wildly. I may be wrong or missing something here, but I believe that having your deconvolved signal in the ideal case be centred at $0 \, \textrm{dB}$ is not very helpful. You may have another reason to do (cont.) Jan 9 at 21:42
• (cont.âed) so but I canât think of a practical advantage of that. However, I didnât manage to get a scaling that would result in an impulse with unity (close to unity âcause some energy will be âspiltâ to the ripples next to the peak) and none of the methods provides that (or at least I didnât manage to achieve it, I may have screwed up somehow). Would you like to elaborate a bit on why youâd like to achieve the $0 \, \textrm{dB}$ average or what advantages this may provide? Jan 9 at 21:45
• Nevermind, I found the answer in your answer to the question providing the âcorrectâ scaling you cited. Jan 10 at 11:18