Is there a way to derive an FIR filter using an IIR filter?

Question

I know there's a thread similar to this one here, but the OP is asking the reverse of what I'm trying to find here. I've done some research on the web with very few sources coming up with actual solutions to this problem. What techniques are used to give an approximation of an FIR filter given one or more IIR filters with say the same order?

Answer 1

Approximating the frequency response of an IIR filter or physical process using an FIR filter is useful in learning control. It is quite common to do FIR filter design based on frequency response specifications. You probably want to check out two standard papers on the subject:

[1] J. H. McClellan, T. W. Parks, and L. R. Rabiner, “A computer program for designing optimum FIR linear phase digital filters,” IEEE Trans. Audio Electroacoust., vol. 21, no. 6, pp. 506–526, 1973.

[2] L. R. Rabiner, “Techniques for Designing Finite-Duration Impulse-Response Digital Filters,” IEEE Trans. Commun. Technol., vol. 19, no. 2, pp. 188–195, Apr. 1971.

Broadly, you either do windowed direct sampling of your desired frequency response, or you use one of several optimization methods to achieve similar results. If you disregard the linear phase delay in an FIR, you can practically make the IIR and FIR responses identical, if the FIR filter order is high enough.

As an elaboration on one of the other answers given; if you have an IIR filter $G(z^{-1})$, then you can do FIR filter design by frequency sampling by taking $M$ samples of the frequency response of $G(z^{-1})$, denoted $\widehat{G}(k)$, and then taking the inverse discrete Fourier transform (IDFT) of $\widehat{G}(k)$. The unit impulse response $g(n)$ of $\widehat{G}(k)$ is \begin{align*} g(n) = \frac{1}{M} \sum\limits_{k=0}^{M-1} \widehat{G}(k) \text{e}^{j \frac{2 \pi k n}{M}} \; , \end{align*} where $n \in [0, M-1] \cap \mathbb{N}_{0}$. The FIR filter is then expressed in the $z$-domain as \begin{multline*} F(z^{-1}) = g(0) + g(1) z^{-1} + ... + g(M-1) z^{-M+1} = \sum_{n=0}^{M-1} g(n) z^{-n} \; . \end{multline*}

The frequency-sampling method results in a unit impulse response which has been convoluted with a rectangular window of the same length in the frequency domain. The frequency response of $F(z^{-1})$ is therefore affected by the large side-lobes of the rectangular window. As a result, the approximation error of $F(z^{-1})$ is large between the frequency samples. This can be alleviated by the use of a window that do not contain abrupt discontinuities in the time domain, and thus have small side-lobes in the frequency domain, i.e., the window smooths the frequency response of $F(z^{-1})$.

A windowed FIR filter $\tilde{h}(n)$ is created from an un-windowed FIR filter $h(n)$ as \begin{align*} \tilde{h}(n) = w(n) h(n) \end{align*} where $w(n)$ is a window function which is non-zero only for $n \in [0, M-1] \cap \mathbb{N}_{0}$. The frequency-domain representation of the window function $W(k)$ is found as \begin{multline*} W(k) = \sum\limits_{n=0}^{M-1} w(n-M/2) \textrm{e}^{-j \frac{2 \pi k n}{M}} = \left[ \sum\limits_{n=0}^{M-1} w(n) \textrm{e}^{-j \frac{2 \pi k n}{M}} \right] \textrm{e}^{-j \frac{2 \pi k}{M} \frac{M}{2} } \; , \end{multline*} where the term $\textrm{e}^{-j (2 \pi k / M) (M/2) }$ comes from the fact that the rectangular window is not centered around $n=0$, but is time-shifted to be centered around $n=M/2$. This phase term will cause distortion of $h(n)$, unless $h(n)$ is also phase-shifted to compensate. The unit impulse response $g(n)$ is therefore phase-shifted before windowing. Due to the circular shift property of the DFT, this can be done by rearranging $g(n)$ such that \begin{equation*} \bar{g}\left( n \right) = \begin{cases} g\left( n + M/2 \right) , & \hspace{-0.6em} n = 0,1, ..., \frac{M}{2} - 1 \\ g\left( n - M/2 \right) , & \hspace{-0.6em} n = \frac{M}{2},\frac{M}{2}+1, ..., M-1 \end{cases} \end{equation*} for the case when $M$ is even. The response is then represented by the FIR filter \begin{equation*} \bar{F}(z^{-1}) = \sum_{n=0}^{M-1} \bar{g}(n) z^{-n} = z^{-M/2} F(z^{-1}) \end{equation*} which is $F(z^{-1})$ delayed by $M/2$ steps. Applying the window $w(n)$ to the time-shifted impulse response $\bar{g}(n)$, \begin{equation*} \tilde{g}(n) = w(n) \bar{g}(n) \; , \end{equation*} the filter \begin{equation*} \tilde{F}(z^{-1}) = W(z^{-1})*\left[ z^{-M/2} F(z^{-1}) \right] \end{equation*} is obtained. Now, $G^{-1}(z^{-1}) \left[ W(z^{-1})*F(z^{-1}) \right] \approx 1$ if the FIR filter is accurate. Note that the phase due to $z^{-M/2}$ is taken out.

Answer 2

The easiest approach is to consider the impulse response of the IIR which is infinite and truncate it somewhere (depending on what order you consider for the approximate FIR filter).

For example, consider the IIR filter with the impulse response $h[n]=a^nu[n]$, where $a$ is positive and $|a|<1$. We can represent it as $$h[n]=\sum_{k=0}^{\infty} a^k\delta[n-k]$$

So the impulse response of the $N$'th order approximation FIR filter would be $$h_{\text{FIR}}[n]=\sum_{k=0}^{N} a^k\delta[n-k]$$

Larger $N$ you consider, closer the FIR will be to the original IIR.

This is an easy approach to simulate the IIR filter's behavior in general. You should be more specific about what aspect of the IIR filter you want to simulate (e.g, pass-band behavior, pass-stop transition, etc.) to get a more specialized answer.

---

In the example below the IIR filter $$H(z)=\frac{1}{1-0.9z^{-1}}$$ is approximated by three FIR filters of orders $N=10,15,25$ where $$H_{\text{FIR}}(z)=\sum_{k=0}^{N} 0.9^kz^{-k}$$

    b1 = 1;
    a1 = [1 -0.9];                   % IIR filter with impulse response (0.9)^n*u[n]

    [H,w] = freqz(b1,a1);            % Plot the frequency response
    plot(w/pi,10*log10(H),'b','Linewidth',2);

    hold on;                         % Plot setup
    text = 'IIR Filter    ';
    color = ['k','g','r'];

    N = [10 15 25];                  % Three different FIR filter orders

    for i=1:3                        % Truncate the impulse response
        b2 = [];
        for n=0:N(i)
            b2 = [b2 0.9^n];
        end
        [H,w] = freqz(b2,1);         % frequency response of FIR filter of order N
        plot(w/pi,10*log10(H),color(i));
        text(i+1,:)=['FIR order = ' num2str(N(i))];
    end
    grid on
    legend(text)
    xlabel('Normalized Frequency')
    ylabel('Magnitude (dB)')

Answer 3

[EDIT: beyond my initial belief that "nobody does that", the OP made me think of situations where this could be useful. Let's start with the obvious]

Given the FIR with $z$-transform: $$\sum_{i=0}^P b_iz^{-i},$$ you can get a very close IIR approximation with:

$$\frac{\sum_{i=0}^P b_iz^{-i}}{1+\sum_{j=1}^Q a_iz^{-i}}$$ with $Q\le P$, and the $a_i$ of very small absolute value, as long as you want to keep "say the same order". An example is given below. I still wonder about the practical interest of such a design.

Perhaps to introduce some instability in an FIR, that is too good in that respect :)

data = randn(1024,1);
fFIRNum = [1 2 1];
fFIRDen = [1];
fIIRDen = [1 0 1e-6];
subplot(3,1,1)
plot([data])
legend('Data')
axis tight;grid on

subplot(3,1,2)
plot([filter(f1,f2,data),filter(f1,f3,data)])
legend('FIR','IIR')
axis tight;grid on

subplot(3,1,3)
plot([filter(f1,f2,data)-filter(f1,f3,data)])
legend('FIR/IIR difference')
axis tight;grid on

The obvious put aside, let me imagine a context where an IIR approximation could be useful. Suppose you want to perform a moving average filtering. If you want to make it adaptive, you have to change the length of the window, and a sudden change in the number of averaged samples could affect the smoothed signal abruptly. At minimum, you can only change the window length by $\pm 1$ unit length. The Exponentially Weighted Moving Average (EWMA) $$y(n) = ax(n) + (1 – a)y(n–1)\,.$$ is an IIR. It may mimic FIR rectangular windows of different lengths, depending of the forgetting factor $a$. The Exponentially Weighted Moving Average has been discussed here recently.

One could perform an adaptive EWMA by smoothly varying $a$ in a more continuous way that hoping the window length from sample to sample. One instance can be found in An Adaptive Exponentially Weighted Moving Average Control Chart, 2003.