The Design of FIR Filters via Levinson Recursion and Optimal Whitening

Using statistical techniques for deterministic filter design.

Design of FIR Filters via Levinson Recursion and Optimal Whitening - by Jeffrey Taft

2.4MB ∙ PDF file

Download

Printer-friendly version

Download

Abstract

Standard non-recursive or Finite Impulse Response (FIR) filter design techniques such as Remez exchange (Parks-McClellan) and windowed coefficients are widely used to design frequency-selective filters. This paper describes techniques to design a variety of FIR filters, including high-Q notch filters using optimal whitening and Levinson recursion. The technique, which relaxes the symmetric coefficient constraint, is quite efficient and can produce low-pass, high-pass, band-pass, band-stop, notch and multi-notch filters and band limited differentiators.

1 Background

Notch filters are used to remove fixed frequency harmonic interference, especially mono-frequency (sinusoidal) interference. Conventionally, FIR notch filters may be designed via Remez exchange and windowed Fourier techniques. The filters designed by these means are of limited performance.

It has been known for some time that better FIR notch filters can be obtained from the convergence of LMS adaptive noise cancellers fed by the sinusoid of interest plus random noise. The LMS adaptive noise canceller is known to converge to a specific optimal solution under stationary statistics.

This paper presents a means to design fixed coefficient FIR notch filters directly from the specifications without the need for LMS simulation.

2 Optimal Linear Prediction and Whitening

The theory of optimal linear prediction is well-developed¹ and has produced the following useful results:

a) For a stationary signal with model

\(s_{yy}(z) = {\sigma_e}^2 B(z) B(z^{-1})\)

as guaranteed by the spectral factorization theorem, there exists a prediction filter as shown in Figure 1.

Figure 1: Linear Predictor

\(y_n \longrightarrow \fbox {$\strut z^{-1} $} \xrightarrow{y_{1}(n)} \fbox {$\strut H(z) $} \xrightarrow{} \overline{y}_{n/n-1}\)

\(H(z) = z \left[ 1 - \frac{1}{B(z)} \right] \; ; \enspace s_{yy}(z) = {\sigma_e}^2B(z)B(z^{-1})\)

b) From the prediction filters we may derive an optimal whitening filter as shown in Figure 2.

Figure 2: Whitening Filter

The derivation of the whitening filter

\([\;1\; a_1\; a_2 \; … \;] \)

is feasible when the signal model is autoregressive. In this case the filter length is finite and the relevant equation is

\(R_{yy} \; \underline{a} = { {\left[ \; {\sigma_e}^2 \; 0 \; 0\; ... \; 0 \;\right]}^T } \qquad (1)\)

where

\(\underline{a}\)

is the filter coefficient vector

\({\left[\;1 \; a_1 \; a_2 \; ... \; a_n \;\right]}^T .\)

The matrix

\(R_{yy}\)

is composed of autocorrelation lags for the signal

\(y_k\)

and is in the form of a symmetric Toeplitz matrix. Equation (1) is therefore efficiently solved via Levinson recursion.

3 Application to Notch Filter Design

The autocorrelation lag vector

\(r(i)\)

is used to generate the matrix

\(R_{yy} \text{ .}\)

The vector can be calculated for a sinusoid from

\(r(i) = \cos \left[ \;2\pi f_n \cdot i \cdot T_s \; \right] \; \exp \left( -r_0 \cdot i{/}m \right) \qquad (2)\)

where

\(\begin{align} f_n &= \text{desired notch frequency}\\ T_s &= \text{signal sample period}\\ r_0 &= \text{zero location radius factor (usually about 0.001)}\\ m &= \text{desired order of the filter.} \end{align}\)

The radius factor

\(r_0\)

is used to keep zeros of the filter from falling directly on the unit circle. The remaining factor to be specified is noise power,

\({\sigma_e}^2\text{ .}\)

The equations are solved via Levinson recursion and the resulting coefficients are normalized to a gain of 1.0 at some specific single frequency of interest.

Note that it is possible to notch out more than one frequency this way. All that is necessary is to accumulate the autocorrelation lags for each notch frequency into

\(r(i)\)

before constructing

\(R_{yy} \text{ .}\)

The notch frequencies may be chosen independently, or they may be related as harmonics of some fundamental notch frequency.

The magnitude frequency response for a sample notch filter is shown in Figure 3. The notch attenuation is approximately 72dB and the antiresonant Q is approximately 15.

Figure 3: Magnitude Frequency Response

4 Broadband Filter Design

It is possible to extend the method of Section 3 to design of broadband pass and reject filters in two ways. The first way is to group some number of notches over the band to be attenuated in a uniform manner. This will cause a suppression with overlap of the (non-ideal) notches, resulting in band suppressions. Selection of the region or regions to be suppressed allows design of band-reject, band-accept, low-pass and high-pass filters. The number of notches to be used is a selectable degree of freedom.

A second, more sophisticated method of broadband design is to generate a function which describes the autocorrelation function of a broad band of frequencies. From Equation (2), we use the core term for a single frequency

\(f_n \text{ ,}\)

and integrate it over a band of frequencies to get (for normalized f)

\(r_{bb}(i) = \int_{f_1}^{f_2} \cos \left[ \;2\pi f \frac{i}{T_s} \cdot T_s \; \right] \; df \qquad (3)\)

which yields

\(r_{bb}(i) = \frac{\sin(2 \pi f_2 \cdot i)}{2 \pi i} - \frac{\sin(2 \pi f_1 \cdot i)}{2 \pi i} \qquad (4)\)

where i is a counting index and we recognize Equation (4) as a difference of sine functions. At i=0, we evaluate Equation (4) as

\(r_{bb}(0) = f_2 - f_1 \qquad (5)\)

after application of L’Hôpital’s Rule. We can use Equations (4) and (5) to cover as many rejection bands as necessary to define the desired frequency response.

Figures 4 through 7 illustrate low-pass, high-pass, band-pass and band-reject filters designed using the methods of this section.

Figure 4: Low-Pass Filter

Figure 5: High-Pass Filter

Figure 6: Band-Pass Filter

Figure 7: Band-Reject Filter

5 Numerical Calculation of the Autocorrelation Function

It is also possible to calculate autocorrelation function samples from frequency domain specifications. In this approach, we created idealized frequency domain representations of the magnitude spectrum, but with the twist that frequencies to be suppressed or removed have large spectral values while frequencies to be retained have zero values. The spectral magnitude function must be real, non-negative and must have even conjugate symmetry. Denoting the spectral model as S(k), we have

\(r(n) = \sum_{k=0}^{N-1}{S(k)e^{+j \frac{kn}{N}2\pi}} \; ; \enspace n=0, 1, ... , N-1 \qquad (6)\)

where

\(j=\sqrt{-1}\)

and where we have used the relation between autocorrelation and Fourier transform. We select the sequence length N to be a convenient power of 2 greater than twice the intended length of the filter to be designed. The computed autocorrelation sequence is then truncated to the length of the filter being designed and is used directly in Equation (1)

Example of a low-pass filter deisgn and a band-limited differentiator are shown in Figures 8 and 9.

Figure 8: Low-Pass Filter

Figure 9: Band-Limited Differentiator

The method can be further extended by using analog filter prototypes and bilinear transformation to create FIR equivalents to IIR filters. In particular, a Butterworth prototype high-pass filter can be used to generate first a discrete transfer function

\(H_d(z)\)

from

\(H_a(s) \text{ ,}\)

then a frequency response magnitude

\(\lvert H_d(e^{j\omega}) \rvert \;,\)

then an autocorrelation function via FFT and finally a coefficient set

\(\{ h_k \}\)

from the autocorrelation function via Levinson recursion. An example low-pass filter designed this way is shown in Fig. 10.

Figure 10: Low-Pass Filter

1. Orfanidis, Sophocles, Optimal Signal Processing, An Introduction, Macmillan Publishing Company, New York, 1985.