3
$\begingroup$

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?

$\endgroup$
  • $\begingroup$ seems to me that the other thread is about the same thing as your question. i always thought that the Prony method was the conventional way to address it. $\endgroup$ – robert bristow-johnson Nov 26 '16 at 0:25
  • $\begingroup$ Given an IIR filter, you want to find the FIR yes? because that would be the inverse of the link you provided... $\endgroup$ – msm Nov 26 '16 at 0:26
  • 1
    $\begingroup$ That seems to be an odd question: an FIR filter is just a special case of an IIR filter just with $a_1$, $a_2$ ... all being zero. Nothing needs to be approximated unless you want some other constraints, i.e. the order of the IIR being a lot lower than that of the FIR. $\endgroup$ – Hilmar Nov 26 '16 at 0:40
  • 1
    $\begingroup$ @Hilmar I really don't think by putting $a_1,a_2,\cdots$ equal to zero, an approximate of the original IIR filter can be achieved. It will become FIR, but not an approximate for the original IIR. As I explained in the answer, I think the correct way is to truncate the time-domain impulse response, not the system function. $\endgroup$ – msm Nov 26 '16 at 1:46
  • $\begingroup$ i think the OP should be clear. even the tread title is self-contradictory. also the IIR $\to$ FIR thing really is just about windowing and truncating the IIR impulse response. i'm gonna make the title self-consistent, and if the OP objects he/she can change it to a consistent meaning. $\endgroup$ – robert bristow-johnson Nov 26 '16 at 2:09
1
$\begingroup$

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)')

enter image description here

enter image description here

$\endgroup$
  • 1
    $\begingroup$ msm, what you are discussing is, essentially, the windowing method of FIR design using the rectangular window. $\endgroup$ – robert bristow-johnson Nov 26 '16 at 2:15
  • $\begingroup$ so, i am sorta asking the same question a different way: *what is the difference between h() and b2()? $\endgroup$ – robert bristow-johnson Nov 27 '16 at 8:16
  • $\begingroup$ @robertbristow-johnson As I also explained in the answer, b2 is the truncated version of the impulse response of the IIR filter with system function (b1,a1). BTW thanks for the edit. $\endgroup$ – msm Nov 27 '16 at 8:25
  • $\begingroup$ it also appears to be the direction (IIR $\to$ FIR) that the OP wants to go. i had it wrong, because i thought the problem is a little trivial and uncommon (if your IIR does what you want, why convert to a higher-order FIR?). one aspect that is not trivial is that might want to window your IIR with a window other than rectangular. perhaps a Hamming (for something simple) or a Kaiser (for something good) and put the 0th sample of $h[n]$ in the center of the window. so you'ld be only applying the right-hand half of the window and the left-hand half of the window multiplies $h[n]=0$ anyway. $\endgroup$ – robert bristow-johnson Nov 27 '16 at 8:30
  • $\begingroup$ Indeed, it seems that by applying a form of weighting the FIR impulse response (i.e. windowing), one could achieve a better approximate with a lower filter order. Somehow, compensating the impact of truncation... This can be formulated as an optimization problem. $\endgroup$ – msm Nov 27 '16 at 8:35
1
$\begingroup$

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.

$\endgroup$
  • $\begingroup$ I read [2] and was wondering - isn't frequency sampling just a type of windowing method? Why does L. R. Rabiner write it as a complete method on its own? As you said yourself, frequency sampling is simply done using a rectangular window, which makes it a type of windowing method, right? $\endgroup$ – ragzputin Nov 29 '16 at 0:33
  • $\begingroup$ @ragzputin: Good question, I'm a little rusty on the theory, but I believe what is meant by the windowing method is to expand the periodic (since it is discrete time) desired frequency response as a Fourier series, using the resulting coefficients as the impulse response, and then windowing the truncated series coefficients (to reduce ripples due to Gibbs phenomenon). This compared to the direct sampling of the frequency response (and optional windowing to reduce ripples)... $\endgroup$ – Arnfinn Nov 29 '16 at 1:09
  • $\begingroup$ @ragzputin: If you have MATLAB or clones, these methods are implemented as fir1 and fir2. One caveat with these implementations is that they assume you want linear phase, so they are not all that useful for IIR filter approximation. $\endgroup$ – Arnfinn Nov 29 '16 at 1:13
  • $\begingroup$ Your first comment didn't say anything about direct frequency sampling. Could you elaborate on frequency sampling? $\endgroup$ – ragzputin Nov 29 '16 at 1:46
  • $\begingroup$ @ragzputin: So, for the windowing method, you choose the FIR coefficients by directly computing a truncated set of the Fourier series coefficients of your desired frequency response (i.e. for a brick wall filter using the formula in en.wikipedia.org/wiki/Pulse_wave). This method is therefore restricted to cases where it is possible to find the Fourier series coefficient (so it wont work for an arbitrary frequency response, like an arbitrary IIR filter response, at least not easily). $\endgroup$ – Arnfinn Nov 29 '16 at 3:03
0
$\begingroup$

[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 :)

Results

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.

$\endgroup$
  • $\begingroup$ of course you can get a very close IIR approximation to $$\sum_{i=0}^P b_iz^{-i},$$ with: $$\frac{\sum_{i=0}^P b_iz^{-i}}{1+\sum_{j=1}^Q a_iz^{-i}}$$ by setting all $a_i=0$. i don't see what the point is. $\endgroup$ – robert bristow-johnson Nov 26 '16 at 2:21
  • 1
    $\begingroup$ @robert bristow-johnson That was exactly the point of my answer. However, if you set $a_i=0$, I am not sure the filter qualifies for IIR $\endgroup$ – Laurent Duval Nov 26 '16 at 8:23
  • $\begingroup$ i think you're right, so that's why i think he meant to ask, how to derive an IIR that matches a given $h[n]$ to with some degree of approximation. $\endgroup$ – robert bristow-johnson Nov 26 '16 at 18:43

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.