1
$\begingroup$

Let's define the relaxation time of a filter as the delay needed for the impulse response of a filter to reach and remain under an amplitude threshold.

Is there a closed-form solution to the computation of an estimate or boundary of the relaxation time of a discrete recursive (IIR) filter from its poles/zeros or from the coefficients of its difference equation, except by computing a reasonable length of its response and checking the absolute value of the response?

$\endgroup$
6
  • $\begingroup$ If you have the transfer function then you can derive its time response (similar to this, for example), and the result will be a product or sum of exponentials and sines/cosines. The exponentials will give you what you're looking for. $\endgroup$ Jun 3, 2020 at 13:54
  • $\begingroup$ I'm not sure I'm following your thought. I don't want to expand the impulse response. If I derive for example the z-domain expression for a 2nd order IIR filter, I obtain an expression with a 4th order denominator, i.e. I've doubled the poles. $\endgroup$
    – moala
    Jun 3, 2020 at 16:01
  • $\begingroup$ If you know the 3 dB bandwidth you can use the rule of thumb that the 90%/10% fall time is approximately 0.35/BW in Hz. This is true for a first order system but reasonably accurate for any IIR system due to its dominant pole dominating the long term response. $\endgroup$ Jun 3, 2020 at 17:00
  • $\begingroup$ So that said, you can identify the dominant pole and then use a simpler equation based on that. (Or if there isn't one dominant pole you can identify the dominant pole/zeros and use a simpler equation based on those) $\endgroup$ Jun 3, 2020 at 17:01
  • $\begingroup$ @moala Well, that was a bit of a brainfart. I was thinking in terms of the inverse Z-transform, but what you get is not like what you get in s domain (what I expected, an $e^{-t}\cos t$-like response). You get a sum of coefficients with powers of the indices, something like $\sum_k A_k\,p_k^k\,u_k$, where $A$, $p$, and $u$ are defined, but $k$ is increasing. Maybe it can be deduced, mathematically, but I don't know how to (though, if it can be done, I wouldn't mind seeing it). $\endgroup$ Jun 3, 2020 at 18:41

1 Answer 1

3
$\begingroup$

You know that calculating the impulse response from the transfer function involves partial fraction decomposition, whether it's in $s$ or $z$. After that you simply apply the inverse $z$-transform to find out the response.

For example, let's say you have this transfer function:

$$H(z)=\dfrac{2+3z^{-1}+4z^{-2}}{7-6z^{-1}+5z^{-2}} \tag{1}\label{1}$$

The coefficients are bogus, but they're chosen so that the poles reside within the ROC. The partial fraction expansion will look like this:

$$H(z)=\dfrac{r_1}{1-p_1z^{-1}}+\dfrac{r_2}{1-p_2z^{-1}} \tag{2}\label{2}$$

From which the impulse response will be:

$$\begin{align} h[n]&=r_1p_1^n+r_2p_2^n \\ {}&=2|p|^n\big[\Re(r)\cos(n\angle{p})-\Im(r)\sin(n\angle{p})\big] \tag{3}\label{3} \\ {}&=2|p|^n|r|\sin(n\angle{p}+\angle{jr}) \tag{4}\label{4} \\ {}&=2|r|\exp(\log(p)n)\sin(n\angle{p}+\angle{jr}) \tag{5}\label{5} \end{align}$$

Notice that $\eqref{4}$ and $\eqref{5}$ use $jr$ as argument for the angle, that's because $\eqref{3}$ has $\pm\Im\sin\pm\Re\cos$, which is reversed compared to the usual angle rotation, $\Re\sin\pm\Im\cos$, which means the argument needs to be rotated (hence the $j$).

For numerical values, the poles and the residues can be found with Octave's residuez():

num = [2 3 4];
den = [7 -6 5];
[r, p] = residuez(num, den)
R =
  -0.2571 + 0.6136i
  -0.2571 - 0.6136i
P =
   0.4286 - 0.7284i
   0.4286 + 0.7284i

The angles for r and for 1j*r are different, but that can be circumvented by either using the first value in the vector for both r and p, or the second. With this you know that the impulse response will be:

% predefined constants
magR = abs(r(1));
angR = arg(1j*r(1));
logP = log(magP);
magP = abs(p(1));
argP = arg(p(1));
n = [0:20];
% the actual impulse response
h = 2*magR*exp(logP*n).*sin(n*angP+angR);
h(1) = num(1)/den(1);

The first value in the vector will be the direct value, given by the evaluation of $H(0)$ (no delay impulse). You can plot this side-by side with the result from impz() (for verification), the "continuous-time" version (for fine-grain viewing), and the decaying curve (to see the match):

t = linspace(0, 20, 1001);
g = 2*magR*exp(logP.*t).*sin(t*angP + angR);
c = 2*magR*exp(logP.*t);
plot(n,h,"or",n,impz(num,den)(1:21),"xb",t,g,"",t,c)

visual verification

The discrete values overlap (the x with the o points), while the "analog" version runs through all of them, with the exponential correctly defining the decay.


The simplification above (with the indices for r and p) may work for this 2nd order example, but for higher orders you need a sum. I got too greedy with simplifications, so the following should address all orders:

num = [2 3 4];
den = [7 -6 5 -4 3];
[r, p] = residuez(num, den);
n = [0:25];
h = real(sum(  abs(r).*exp( n.*log( abs(p) ) ) ...
  .*sin( n.*arg(p) + arg(1j*r) )  ));
h(1) = num(1)/den(1);
t = linspace(0, 25, 1001);
g = sum( r.*p.^t );
c1 = sum( abs(r).*exp(t.*log(abs(p))) );
c2 = sum( abs(r).*abs(p).^t );
plot( n, h, "or", ...
      n, impz(num, den)(1:length(n)), "xb", ...
      t, g, "g", ...
      t, abs(g), ":k", ...
      t, c1, ".-c", ...
      t, c2, "-.m" )
grid on;

For the sake of computation you can predefine absR=abs(r), argP=arg(p), etc. This time I've modified the graphs to include:

  • h -- the discrete impulse response as calculated above, red "o";
  • impz() -- the discrete impulse response as given by the built-in function (as reference), blue "x";
  • g -- the "continuous-time" approach (densely sampled), green line;
  • abs(g) -- the absolute value if g for better comparing the decaying exponential envelope, black dotted line;
  • c1 -- the calculated decaying envelope, with $\text{e}^{-x}$, cyan dot-dash line;
  • c2 -- the "original" envelope, with $p_k^n$, magenta dash-dot line.

The same sum is applied to $\eqref{3}$, $\eqref{4}$, and $\eqref{5}$. The degree of $H(z)$ has been increased by adding elements to the denominator's vector (here a 4th order). And the picture shows an overlapping h with impz(), thus identical responses, a g that goes through all the points in h and impz(), and overlapping envelopes providing the boundary for abs(g):

4th order works, too

To conclude on point, the main difference beteen the discrete and the analog exponential decays:

  • the $s$-domain has $\text{e}^{\Re(p)t}$;
  • the $z$-domain has $\text{e}^{\log(|p|)n}$.

If the $s$-domain is stable then its real part is negative, and if the $z$-domain is within the ROC then the absolute value will be <1 and, as such, the logarithm will be negative.

$\endgroup$
6
  • $\begingroup$ I doubt OP is here, anymore, but I've updated the answer with a more generic approach (and a conclusion). $\endgroup$ Jun 30 at 23:13
  • $\begingroup$ dunno if it's correct but +1 for effort & apparent quality $\endgroup$ Jul 1 at 15:55
  • 1
    $\begingroup$ @OverLordGoldDragon Well, I say it is but, if it isn't, that's why I provided the code. Thank you for the vote of confidence, though. $\endgroup$ Jul 1 at 21:13
  • $\begingroup$ Thanks, I'll have a look at some time in the future. $\endgroup$
    – moala
    Jul 17 at 22:41
  • 1
    $\begingroup$ Oh, I just found it. I think this also answers my question $\endgroup$
    – ZR Han
    Jul 19 at 9:37

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.