3
$\begingroup$

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.

Improve this question
$\endgroup$
8
  • 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$ – robert bristow-johnson Nov 19 '16 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 '16 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$ – jojek♦ Nov 19 '16 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$ – robert bristow-johnson Nov 20 '16 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 '16 at 22:16
1
$\begingroup$ 
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');
end


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

plot(h);

sound(h, responseFs);
Improve this answer
$\endgroup$
0
$\begingroup$

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.

Improve this answer
$\endgroup$
0
$\begingroup$

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.

https://en.m.wikipedia.org/wiki/Deconvolution

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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