1

I'm learning a method to genenrate pink noise by filtering a white noise in this book, and the author mentioned an estimation of T60 given the coefficients of an IIR filter in the code that

B = [0.049922035 -0.095993537 0.050612699 -0.004408786];
A = [1 -2.494956002   2.017265875  -0.522189400];
nT60 = round(log(1000)/(1-max(abs(roots(A))))); % T60 est.

where nT60 is reverberant time (T60) in terms of discrete sample, and A is the denominator coefficients.

This method gives a rough T60 estimation of the given coefficients which is 1430. The integrated impulse response method gives a estimation of around 1426.

h = impz(B, A, 2048);
td = length(h);
L(td:-1:1)=10*log10(cumsum(h(td:-1:1).^2)/sum(h(1:td).^2));
% plot the energy decay curve, and then calculate the slope of the curve

These two results are surprisingly close enough considering the simplicity of the former one.

In the field of acoustics we usually calculate T60 for an impulse response which is often an FIR filter. I'm wondering if there is any theoretical derivation of this method. Any reference would be very grateful.

ZR Han
  • 3,228
  • 6
  • 16

1 Answers1

1

The notion of a single decay time for an IIR filter can be a bit tricky. The decay times are solely a functions of the pole locations: the closer a pole to the unit circle, the longer the decay time is.

In the early part of the impulse response, the decay of all poles is visible. However, the poles with lower decay time "die off" first and at the end you are left with only one pole: the one closest to the unit circle. That's the one that determines the overall decay time of the filter.

This works really well in your example. The filter has three real poles and the decay times for each are 1430, 123 and 16 samples. The contribution of the 123 and 16 decay contributions die down very quickly so you are left with a single exponential decay.

If your largest poles have similar distance to the unit circle, things get more complicated.

The "revere integration" method has been developed by Manfred Schroeder (very smart and nice guy! https://en.wikipedia.org/wiki/Manfred_R._Schroeder) . He basically proved that the ensemble average over a large number of individual decay instances will converge to the time reversed integral over the squared impulse response. So instead of having to average over many, many individual decay measurements, you can just measure the impulse response, square and reverse integrate and you get the same result. The derivation is quite clever but I can't easily find a publicly available version of it.

The two method will give the same result for any process that's close enough to a pure exponential decay with a single time constant. Since your example is just that, you get very good agreement.

Hilmar
  • 44,604
  • 1
  • 32
  • 63
  • Good explanation, thank you Hil. Could you please give me more details about how to derive the relationship between the decay time and the pole position? T60 is the time that sound pressure level decays by 60 dB, ie, 1/1000. If we assume that the impulse response is a pure exponential decay $e^{-\alpha n}$, we have $\alpha n_{T_{60}} = \ln 1000$ and thus $n_{T_{60}} = \ln1000/\alpha$. According to the code it seems that $\alpha$ equals to $1-\max(|p_i|)$ where $p_i$ is the pole but I don't know how to derive the last step. – ZR Han Aug 31 '21 at 06:19
  • @ZRHan, even is this comment is old I guess it might also help others. As I am interested in the settling time of IIR filters I found this page which explains the relation of one pole and the time constant Time_Constant_One_Pole – Irreducible Sep 01 '23 at 08:35
  • @Irreducible Thank you! That's exactly what I wanted to ask. – ZR Han Sep 04 '23 at 01:33