Does someone have a matlab script to get the impulse response from the recorded measurements using log sweep? I have the input file (sweep) and the output and I need to get the impulse response and then calculate reverberation times for octave bands.

  • 2
    $\begingroup$ don't have a script. fft both the input file and the output file (perhaps zero-padding both). divide the fft of the output with the fft of the input. then ifft the result. $\endgroup$ Nov 19, 2016 at 20:59
  • $\begingroup$ thanks! but and then how can I get the reverberation times for octave bands? my output file has lots of zero values, and then I got NaN when I use dB...how to solve this? $\endgroup$
    – rfbrgs
    Nov 19, 2016 at 22:09
  • $\begingroup$ This sounds like a more complicated problem. Could you please rephrase your question so it is applicable to any language? We are trying to be language-agnostic on DSP SE. Nonetheless, please have a look here first. What you have to do is to convolve your recording with time-reversed and scaled signal. When it comes to calculation of reverberation time in octaves, it's a different question I believe. $\endgroup$
    – jojeck
    Nov 19, 2016 at 22:21
  • $\begingroup$ both division by zero and the logarithm of zero is a problem. usually what i do is add a very small positive number (like $\epsilon = 10^{-12}$) to the positive value in the denominator or the $\log(\cdot)$ argument. $\endgroup$ Nov 20, 2016 at 2:42
  • $\begingroup$ @robertbristow-johnson can you turn your comment into an answer here or here with a few math details? Would be super cool, thanks in advance! $\endgroup$
    – Basj
    Dec 28, 2016 at 22:16

3 Answers 3

referenceFilename = 'reference.wav'; 
responseFilename = 'response.wav'; 

[reference, referenceFs] = audioread(referenceFilename);
[nReferenceSamples, nReferenceChannels] = size(reference);

[response, responseFs] = audioread(responseFilename);
[nResponseSamples, nResponseChannels] = size(response);

if referenceFs ~= responseFs
  error('reference audio and response audio files must have same sample rate');

nFFT = 2^( ceil( log2(nReferenceSamples+nResponseSamples) ) );

x = zeros(nFFT, 1);
y = zeros(nFFT, 1);

x(1:nReferenceSamples) = reference(:,1);    % use only left channel
y(1:nResponseSamples) = response(:,1);      % use only left channel

h = ifft( fft(y) ./ fft(x) );

h = h(1:nFFT/2);                            % truncate latter half


sound(h, responseFs);

Once you obtain the IR you must square it, integrate it backwards (Schroeder integral), use octave filters and calculate the slope which is the RT.


A deconvolution algorithm would work. In the simplest form it is the inverse Fourier transform of the result of dividing the the Fourier transform of the output and the Fourier transform of the input. This may not work if the FT of the input has zeros in it but there are algorithms referenced in the link below which address this issue.



Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.