Introduction to Digital Filters

3.10 Introduction to Digital Filters

3.10.1 Designing simple filters using specialized software

Before studying digital filters in more details, it is useful to briefly discuss analog filters because most of the times they are required in digital signal processing systems as will be discussed here.

Example of analog filter design

There are integrated circuits (IC) that implement analog filters. They can also be built using discrete components: resistors, capacitors, inductors and operational amplifiers (opamps). The design of a filter corresponds to obtaining its system function. In the case of analog filters, there are many recipes (well-established algorithms, also called approximations) to obtain a system function $H (s)$ that obeys the specified requirements in a project. The most common approximations are the Butterworth, Chebyshev, Bessel and Cauer (or elliptic).

Having $H (s)$ , the next step is the realization of the filter using a circuit, i. e., establishing the topology and the components of the circuit that (approximately) implement $H (s)$ . These two steps are the subject of many textbooks on filter design. Here, the approach will be to simply provide an example using the FilterPro software from Texas Instruments [ url3tif]. This software allows to obtain the schematic and components.¹⁴

Analog filters that use only resistor, capacitors and inductors are called passive filters. An alternative to avoid using inductors¹⁵ are active filters, which are typically based on opamps. When the target is an IC, not a circuit with discrete components, inductors are avoided due to the difficulty of their integration using current microelectronics technology. When an analog filtering operation is required in an IC, an active filter is typically adopted and this is the only class of filters supported by FilterPro.

Example 3.13. Example of analog filter designed with FilterPro. The following example illustrates the design of an analog filter with a desired cutoff frequency $f_{c} = 1$ kHz, which is assumed to define the passband, such that the passband frequency is $f_{p} = f_{c}$ . A second-order lowpass filter with passband frequency of approximately 977.2 Hz was obtained with FilterPro and its schematic¹⁶ is shown in Figure 3.28.

Figure 3.28: Lowpass active analog filter designed with Texas Instruments’ FilterPro software.

The system function $H (s)$ can be obtained by analyzing the circuit in Figure 3.28. After that, $H (s)$ can be compared to the ideal one obtained e. g. with Matlab/Octave via the command: [Bs,As]=butter(2,2*pi*1000,’s’) that gives

H (s) = \frac{B (s)}{A (s)} = \frac{3.9478 \times 1 0^{7}}{s^{2} + 8.8858 \times 1 0^{3} s + 3.9478 \times 1 0^{7}} .

(3.48)

The arguments of the function butter are the filter order and the cutoff frequency in rad/s (not Hz). Note the required ’s’ to indicate to Matlab/Octave that the goal is to design an analog, not a digital filter. $□$

Example 3.14. Further examples of analog filter design in Matlab. Listing 3.11 illustrates how to design some analog filters with a given specification for its magnitude using Matlab (this code is not compatible with Octave). The frequency response magnitude and phase of the four designed filters are depicted in Figure 3.29.

Listing 3.11: MatlabBookFigures/figs_systems_elliptic.m

1Apass=5; %maximum ripple at passband 
2Astop=80; %minimum attenuation at stopband 
3%% Lowpass Elliptic %%%%%%%%% 
4Fp=100; Wp=2*pi*Fp; %passband frequency 
5Fr=120; Wr=2*pi*Fr; %stopband frequency 
6[N, Wp] = ellipord(Wp, Wr, Apass, Astop, 's') %find order 
7[z,p,k]=ellip(N,Apass,Astop,Wp,'s'); %design filter 
8B=k*poly(z); A=poly(p); %convert zero-poles to transfer function 
9[H,w]=freqs(B,A); %calculate frequency response 
10%% Higpass Elliptic %%%%%%%%% 
11Fr=100; Wr=2*pi*Fr; %stopband frequency 
12Fp=120; Wp=2*pi*Fp; %passband frequency 
13[N, Wp] = ellipord(Wp, Wr, Apass, Astop, 's') %find order 
14[z,p,k]=ellip(N,Apass,Astop,Wp,'high','s'); %design filter 
15B=k*poly(z); A=poly(p); %convert zero-pole to transfer function 
16[H,w]=freqs(B,A); %calculate frequency response 
17%% Bandpass Elliptic %%%%%%%% 
18Fr1=10; Wr1=2*pi*Fr1; %first stopband frequency 
19Fp1=20; Wp1=2*pi*Fp1; %first passband frequency 
20Fp2=120; Wp2=2*pi*Fp2; %second passband frequency 
21Fr2=140; Wr2=2*pi*Fr2; %second stopband frequency 
22[N,Wp]=ellipord([Wp1 Wp2],[Wr1 Wr2],Apass,Astop,'s') %find order 
23[z,p,k]=ellip(N,Apass,Astop,Wp,'s');%design filter 
24B=k*poly(z); A=poly(p); %convert zero-pole to transfer function 
25[H,w]=freqs(B,A); %calculate frequency response 
26%% Bandpass Butterworth %%%%%%%%%%% 
27[N, Wn]=buttord([Wp1 Wp2],[Wr1 Wr2],Apass,Astop,'s')%find order 
28[z,p,k]=butter(N,Wn,'s'); %design Butterworth filter 
29B=k*poly(z); A=poly(p); %convert zero-pole to transfer function 
30[H,w]=freqs(B,A); %calculate frequency response

Figure 3.29: Frequency response of filters designed in Listing 3.11. The magnitude specification masks are indicated. Note that the phase was unwrapped with the command unwrap for better visualization.

It should be noted from Listing 3.11 and Figure 3.29 that a bandpass filter has twice the order that is specified as an input parameter (and obtained by buttord and ellipord in this example). Note also that, in practice, an analog filter typically has an order between 1 to 10. For example, an order $N = 86$ is not feasible for practical implementation with an analog circuit. $□$

Example of digital filter design

Figure 3.30: The canonical interface of a digital filter $H (z)$ with the analog world via A/D and D/A processes. The analog $A (s)$ and $R (s)$ are the anti-aliasing and external reconstruction filters, respectively.

For comparison, it is interesting to design now an equivalent digital filter $H (z)$ to substitute its analog counterpart corresponding to $H (s)$ . But in order to interface $H (z)$ with the external (analog) world, four extra blocks are required, as illustrated in Figure 3.30. The filters $A (s)$ (anti-aliasing) and $R (s)$ (external reconstruction) are lowpass analog filters with passband frequency equal to half of the sampling frequency, which is assumed here to be $F_{s} = 10$ kHz. Their gains are not important at this moment and in practice they are determined by the interfacing electronics. Assuming the converters and analog filters are properly defined, it remains to design $H (z)$ , which can be done with the command

1[Bz,Az]=butter(2,(1000/5000))

that leads to

H (z) = \frac{B (z)}{A (z)} = \frac{0.067455 + 0.134911 z^{- 1} + 0.067455 z^{- 2}}{1.0 - 1.14298 z^{- 1} + 0.41280 z^{- 2}} .

(3.49)

It is important to note that when designing an analog filter, the second argument of butter is $ω$ in rad/s. For a digital filter, as explained in Section 1.8.4, Matlab/Octave uses the normalized frequency $f_{N} = Ω ∕ π$ . In this case, the Nyquist frequency $F_{s} ∕ 2$ is 5 kHz and using Eq. (1.40), $f = 1000$ Hz is normalized by the Nyquist frequency in the command [Bz,Az]=butter(2,(1000/5000)).

It will be later discussed in this chapter that the function butter used the bilinear transformation to convert $H (s)$ of Eq. (3.48) into $H (z)$ of Eq. (3.49).

From Table 3.1 one knows the definition $H (z) = Y (z) ∕ X (z)$ and this should not be confused with $H (z) = B (z) ∕ A (z)$ in Eq. (3.49). In many applications $H (z)$ is a rational function (a ratio of two polynomials), with $B (z)$ specifying the numerator and $A (z)$ the denominator. In general, $B (z) \neq Y (z)$ and $A (z) \neq X (z)$ .

Digital filters described by a rational function $H (z)$ can be conveniently implemented via a LCCDE, as discussed next.

Example 3.15. Obtaining the difference equation and implementing a digital filter. For example, assume that a LTI system with impulse response $h [n] = 0 . 2^{n} u [n]$ outputs $y [n] = δ [n] - 0.3 δ [n - 1]$ when the input is $x [n] = δ [n] - 0.5 δ [n - 1] + 0.06 δ [n - 2]$ . In this case, $X (z) = 1 - 0.5 z^{- 1} + 0.06 z^{- 2}$ and $Y (z) = 1 - 0.3 z^{- 1}$ , which leads to $H (z) = 1 ∕ (1 - 0.2 z^{- 1})$ with $B (z) = 1$ and $A (z) = 1 - 0.2 z^{- 1}$ . With this concept in mind, Eq. (3.49) can be written as

\begin{array}{l} X (z) (0.067455 + 0.134911 z^{- 1} + 0.067455 z^{- 2}) \\ = Y (z) (1.0 - 1.14298 z^{- 1} + 0.41280 z^{- 2}) \end{array}

and taking the Z-inverse one can obtain the difference equation corresponding to Eq. (3.49):

\begin{array}{l} y [n] & = 0.067455 x [n] + 0.134911 x [n - 1] + 0.067455 x [n - 2] \\ + 1.14298 y [n - 1] - 0.41280 y [n - 2] . \end{array}

The amazing fact is that this can be implemented with a simple code, similar to Listing 3.12.

Listing 3.12: Digital filter implementation in pseudo-code.

1initialization() { %all previous values are assumed 0 
2    xnm1=0 %variable to store the value of x[n-1] 
3    xnm2=0 %x[n-2] 
4    ynm1=0 %y[n-1] 
5    ynm2=0 %y[n-2] 
6} 
7processSample() { 
8    xn=readFromADConverter() %read input sample from ADC 
9    yn=0.067455*xn + 0.134911*xnm1 + 0.067455*xnm2 + 
10       1.14298 ynm1 - 0.41280 *ynm2 
11    writeToDAConverter(yn) %write output sample into DAC 
12    ynm2=ynm1 %update for next iteration. Note the order of ... 
13    ynm1=yn   %updates: avoid overwriting a value ... 
14    xnm2=xnm1 %that should be used later 
15    xnm1=xn 
16}

where one assumes a timer¹⁷ invokes the function processSample at a rate of $F_{s}$ (10 kHz in this case).

Figure 3.31: Comparison of two analog filters with its digital counterpart of Eq. (3.49). The “ideal” analog corresponds to Eq. (3.48) while the “10% tolerance” corresponds to its realization using the schematic of Figure 3.28.

Figure 3.31 compares the discussed filters, illustrating that $H (z)$ has a lowpass frequency response that is similar in terms of passband frequency to the original $H (s)$ . For frequencies above the passband, $H (z)$ attenuates more than $H (s)$ , but this is not deleterious for a lowpass filter. Later, in Section 3.12.4, this behavior will be clarified, which is intrinsic of having the bilinear transformation mapping the behavior of the analog filter at infinite frequencies into $Ω = π$ rad (or $F_{s} ∕ 2$ Hz) for $H (z)$ , i. e., mapping $H (ω) |_{ω \to \infty}$ into $H (e^{j Ω}) |_{Ω = π}$ .

If $F_{s}$ is changed, the frequency response of the overall system in Figure 3.30 is scaled according to $ω = Ω F_{s}$ (Eq. (1.36)). For example, decreasing $F_{s}$ from 10 to 5 kHz would divide by two the cutoff frequency of the overall system. $□$

3.10.2 Distinct ways of specifying the “ripple” / deviation in filter design

As illustrated in Figure 3.29, practical filters have a frequency response that deviates from the ideal ones. Instead of a flat magnitude over the whole passband and zero gain in the stopband, as depicted in Figure 3.1, a non-ideal filter has deviations from the target values and eventual “oscillations” called “ripples”. When designing a filter, it is possible to impose limits to these ripples, but this procedure is different for the passband and stopband “ripples” / deviations, as discussed next.

There are several ways of specifying the maximum “ripple” or deviation at the passband of a filter. For example, it can be specified in dB or linear scale, which are called here $R_{p}$ and $D_{p}$ , respectively (as suggested by Matlab). If one assumes a lowpass filter with a unitary gain at DC, the maximum linear deviation $D_{p}$ can be assumed to be within $[1 - D_{p}, 1 + D_{p}]$ or $[1 - D_{p}, 1]$ . The latter option is particularly useful if the filter magnitude monotonically decreases as for the Butterworth filters and is assumed here, meaning that the frequency response is 1 at DC and decreases in the passband to a value not smaller than $1 - D_{p}$ . For example, a specification can be that the passband magnitude does not deviate 5% of the gain at DC. In this case, $D_{p} = 0.05$ and the minimum magnitude at passband is $1 - D_{p} = 0.95$ . Often, this value is specified in dB as

R_{p} = - 20 {log}_{10} (1 - D_{p}) .

If the deviation is within $[1 - D_{p}, 1]$ in linear scale, then it is confined to $[- R_{p}, 0]$ in dB.

At the stopband, the magnitude in linear scale $D_{s}$ is also a “deviation” because the reference value is 0 (the desired magnitude at stopband). In other words, for the stopband, the deviation coincides with the magnitude itself and its value in dB is given by

R_{s} = - 20 {log}_{10} D_{s} .

Often, the filter magnitude has a maximum value of 0 dB and the values of $R_{p}$ and $R_{s}$ are positive.

Matlab uses the same convention adopted here in its routines. This is exemplified in Listing 3.13, which uses the method buttord to find the order of a Butterworth filter that complies with the project specifications.

Listing 3.13: MatlabOctaveCodeSnippets/snip_systems_IIR_design.m. [ Python version]

1Wp=0.3; %normalized passband frequency (rad / pi) 
2Ws=0.7; %normalized stopband frequency (rad / pi) 
3Rp=1; %maximum passband ripple in dB 
4Rs=40; %minimum stopband attenuation in dB 
5[n,Wn] = buttord(Wp,Ws,Rp,Rs); %find the Butterworth order 
6[B,A]=butter(n,Wn); %design the Butterworth filter 
7zp=exp(1j*pi*Wp); zs=exp(1j*pi*Ws); %Z values at Wp and Ws 
8%% Check if the filter obeys requirements 
9linearMagWp = abs(polyval(B,zp)/polyval(A,zp)) %mag. at Wp 
10linearMagWs = abs(polyval(B,zs)/polyval(A,zs)) %mag. at Ws 
11Dp_result = 1-linearMagWp %linear deviation at Wp (linear) 
12Rp_result = -20*log10(1-Dp_result) %deviation (ripple) at Wp (dB) 
13Ds_result = linearMagWs %linear deviation at Ws (linear) 
14Rs_result = -20*log10(Ds_result) %deviation (ripple) at Ws (dB)

The required filter order is n=4, and the obtained results are: linearMagWp=0.9105, linearMagWs=0.01, Dp_result=0.0895, Rp_result=0.8147, Ds_result=0.01 and Rs_result=40.0. They can be interpreted as follows. Matlab designed the method buttord to meet exactly the requirements for the stopband, while giving an eventual slack to the requirements for the passband. Hence, we can see that Rs_result=40.0 dB is exactly the specified attenuation at the stopband. This corresponds to a linear gain linearMagWs=0.01 of the filter’s frequency response $H (e^{j Ω})$ at the beginning of the stopband. Using mathematical notation: $H (e^{j Ω}) |_{Ω = 0.7 π} = 0.01$ in linear scale, which corresponds to $20 {log}_{10} H (e^{j Ω}) |_{Ω = 0.7 π} = 20 {log}_{10} (0.01) = - 40$ dB.

The filter gain is 1 at DC, and it decreases until reaching linearMagWp=0.9105 in the end of its passband. In our notation, $H (e^{j Ω}) |_{Ω = 0.3 π} = 0.9105$ , which corresponds to $20 {log}_{10} H (e^{j Ω}) |_{Ω = 0.3 π} = 20 {log}_{10} (0.9105) = - 0.8144$ dB. This value does not coincide with Rp_result=0.8147 due to the rounding to four decimal places. Using 20*log10(linearMagWp) gives the precise gain value of $- 0.8147$ dB. Note that the deviation in the passband of 0.8147 dB from the ideal 0 dB is within the maximum value of 1 dB imposed to buttord. There is a slack of $1 - 0.8147 = 0.1853$ dB.

In summary, it is a bit confusing how Matlab and other software define the input parameters to functions such as buttord. But one just needs to interpret the ripples in dB as obtained from “deviations” of the ideal values of 1 and 0 in passband and stopband, respectively.

3.10.3 LCCDE digital filters

A digital filter is a discrete-time system that deals with quantized input $x_{q} [n]$ and output $y_{q} [n]$ signals (quantized amplitudes). A possible practical scenario for a digital filter is the one shown in Figure 3.30. There are many different kinds of digital filters but this section assumes filters that are based on linear constant-coefficient difference equations (LCCDE). As indicated in Figure 3.5, when zero initial conditions are assumed, the LCCDE systems are a special case of LTI systems. But when one takes in account the non-linear quantization effects that occur in practical CPUs (or digital signal processors) due to their limited precision, the overall system is not strictly LTI. In spite of that, it is convenient to isolate the effects of quantization and finite-precision arithmetics, and study digital filters as LTI LCCDE systems. Later on, the impact of quantizing the filter coefficients and the roundoff errors due to finite-precision can be incorporated. Therefore, hereafter, the term digital filter refers to discrete-time LTI LCCDE systems. Besides, the systems are assumed to be causal.

A system that has the output and input signals related by a LCCDE such as

\sum_{k = 0}^{N} a_{k} y [n - k] = \sum_{k = 0}^{M} b_{k} x [n - k]

(3.50)

can be represented by (taking the Z-transform on both sides of the LCCDE)

Y (z) \sum_{k = 0}^{N} a_{k} z^{- k} = X (z) ∑_{k = 0}^{M} b_{k} z^{- k} .

In this case, the transfer function is

H (z) = \frac{Y (z)}{X (z)} = \frac{\sum_{k = 0}^{M} b_{k} z^{- k}}{\sum_{k = 0}^{N} a_{k} z^{- k}},

(3.51)

which, after eventual algebraic simplifications, can be written as the ratio

H (z) = \frac{B (z)}{A (z)}

(3.52)

of two polynomials in $z$ . Calculating the zeros $z_{k}$ (roots of the numerator $B (z)$ ) and the poles $p_{k}$ (roots of the denominator $A (z)$ ), one can write

H (z) = g \frac{\prod_{k = 1}^{M} (z - z_{k})}{\prod_{k = 1}^{N} (z - p_{k})},

where $g$ is a gain. The coefficients of Eq. (3.50) can be real or complex. For simplicity, they will be assumed real ( $a_{k}, b_{k} \in ℝ$ ) unless otherwise stated. In this case, $g \in ℝ$ and the roots $z_{k}, p_{k}$ are real or, for each root $r e^{j𝜃}$ , its complex conjugate $r e^{- j𝜃}$ is also a root such that the product $(z - r e^{j𝜃}) (z - r e^{- j𝜃}) = z^{2} + 2 r cos (𝜃) + r^{2}$ leads to a second-order polynomial with real coefficients.

It is useful to know that the value of $H (z)$ at $z = 1$ and $z = - 1$ correspond to the values of the frequency response $H (e^{j Ω})$ at $Ω = 0$ (DC) and $Ω = π$ rad, respectively, which can be written as

H (z) |_{z = 1} = H (e^{j Ω}) |_{Ω = 0} and H (z) |_{z = - 1} = H (e^{j Ω}) |_{Ω = π} .

(3.53)

For example $H (z) = 3 + z^{- 1}$ has a frequency response with $H (e^{j Ω}) |_{Ω = 0} = 3 + 1 = 4$ and $H (e^{j Ω}) |_{Ω = π} = 3 - 1 = 2$ .

Example 3.16. Specifying a filter in Hz, rad or normalized frequency. When designing digital filters, it is important to remember the fundamental relation $ω = Ω F_{s}$ , first presented in Eq. (1.36). Figure 3.32 illustrates how this relation can be used to convert a specification for a lowpass analog filter with passband and stopband frequencies of 500 and 850 Hz, respectively, into a specification of a digital filter with $F_{s} = 2000 Hz$ , which leads to passband and stopband frequencies of $2 π 500 ∕ 2000 \approx 1.57$ and $2 π 850 ∕ 2000 \approx 2.67$ rad, respectively.

Figure 3.32: Comparison of analog filter specification (a) and two corresponding digital versions assuming $F_{s} = 2000$ Hz: (b) was obtained with $ω = Ω F_{s}$ and (c) uses the convention adopted in Matlab/Octave.

Figure 3.32 also shows the convention used by Matlab/Octave, which gets rid of the irrational number $π$ by using a normalized axis $f_{N} = Ω ∕ π$ , as listed in Table 1.5. In this case, $ω = Ω F_{s}$ can be written as $2 πf = f_{N} π F_{s}$ , which simplifies to obtaining $f_{N} = f ∕ (F_{s} ∕ 2)$ with the normalization of $f$ by the Nyquist frequency $F_{s} ∕ 2$ . For example, $f = 500$ Hz in Figure 3.32(a) turns into $f_{N} = 500 ∕ 1000 = 0.5$ in Figure 3.32(c). This normalization by the Nyquist frequency is sensible because digital filters are always constrained to work with frequencies of at most $F_{s} ∕ 2$ Hz, which corresponds to $π$ rad. $□$

3.10.4 FIR, IIR, AR, MA and ARMA systems

An important characteristic of a digital filter is the duration of its impulse response, which allows to classify it as a finite impulse response (FIR) or infinite impulse response (IIR) filter. For example, a filter with $h_{f} [n] = 3 δ [n] + 2 δ [n - 1]$ is FIR and a filter with $h_{i} [n] = 0 . 8^{n} u [n]$ is IIR. Note that if one plots $h_{i} [n]$ in a computer, due to the limited precision, the values of $h_{i} [n]$ would underflow to zero for large enough $n$ , but theoretically this impulse response has infinite duration.

Another characteristic is whether the filter is recursive or not. A recursive filter has the output $y [n]$ depending on past values of the output itself ( $y [n - 1], y [n - 2], \dots$ ) while the output of a non-recursive filter depends only on the input $x [n]$ and its past values (assuming causality). The system function $H (z)$ of a non-recursive filter has $N = 0$ in Eq. (3.51) and no finite poles.

In practice, FIR can be associated to non-recursive filters and IIR to recursive filters. But recursiveness and impulse response duration are two distinct concepts! It is possible to create a recursive filter with an impulse response of finite duration using cancellation of poles and zeros, for example. In spite of that, the design of recursive filters is often called “IIR design”. Similarly, “FIR design” is the topic that concerns the design of non-recursive filters. This widely adopted jargon is used in this text.

IIR and FIR filters will be discussed in the sequel. While the subject is presented in the context of digital filtering, many conclusions are valid for systems in general. In other words, IIR and FIR systems can be found in applications in which the goal is not filtering, but system identification, equalization, control, tracking, statistics, etc. In some of these areas, $H (z)$ is not classified as FIR or IIR, but as MA, AR or ARMA “model”.

A moving-average (MA) model $H (z)$ is a non-recursive FIR filter and an autoregressive (AR) filter is a special case of a recursive IIR filter for which $M = 0$ in Eq. (3.51). The AR IIR filter does not have finite zeros other than at the origin $z = 0$ and is given by

H (z) = g \frac{z^{N}}{\prod_{k = 1}^{N} (z - p_{k})} = \frac{g}{1 + a_{1} z^{- 1} + a_{2} z^{- 2} + \dots + a_{N} z^{- N}} .

(3.54)

Note in Eq. (3.54) that the gain value $g$ for the numerator is chosen such that the denominator coefficient $a_{0}$ of Eq. (3.51) is equal to 1. The numerator $z^{N}$ leads to $N$ zeros at $z = 0$ and makes $H (z)$ causal.

The autoregressive-moving-average (ARMA) model is a generic $H (z)$ model composed by an AR (denominator) and MA (numerator) sections. There are many generalizations of these basic models, which are adopted in areas such as statistics. Parametric modeling is the task of, according to a predefined criterion such as least-squares, fit the data to the assumed model. In spite of being defined in distinct scenarios, filter design and parametric modeling share several fundamental results.

Before discussing the design of IIR and FIR filters, it is useful to note the importance of frequency scaling.

3.10.5 Filter frequency scaling

Computational routines often produce a prototype with normalized frequencies (e. g., a cuttof frequency of 1 rad/s) that are then converted (eventually by another routine) to the required frequencies via frequency scaling.

With frequency scaling one can modify a filter to obtain a different cutoff frequency. For example, if $H (s)$ has a cutoff of $ω_{c}$ rad/s, the transformation $s \to s ∕ k$ moves this cutoff frequency to $k ω_{c}$ rad/s. The operation is similar to using $x (2 t)$ to double the “speed” of a signal represented by $x (t)$ , which leads to $x (2 t)$ having half of the support of $x (t)$ . Similarly, $x (t ∕ 2)$ has twice the support of $x (t)$ . More specifically, if the prototype $H (s)$ has a cutoff frequency of 1 rad/s, then

G (s) = H (s) |_{s = s ∕ ω_{c}}

(3.55)

scales the complex plane, and consequently $s = jω$ , such that $ω_{c}$ is the cutoff frequency of the new filter $G (s)$ .

There are computational routines that provide a prototype filter with (normalized) cutoff frequency of $1$ rad/s, which can then be converted to a new filter with an arbitrary cutoff frequency $ω_{c}$ with the transformation $s \to s ∕ ω_{c}$ . For example, this change in cutoff frequency can be done in Matlab by applying the function lp2lp to a lowpass analog prototype (current version of Octave misses this functionality) as indicated in Listing 3.14.

Listing 3.14: MatlabOnly/snip_systems_filter_conversion.m. [ Python version]

1N = 5; %filter order 
2[B,A] = butter(N,1,'s'); %prototype with cutoff=1 rad/s 
3[Bnew,Anew]=lp2lp(B,A,100); %convert to cutoff=100 rad/s 
4freqs(Bnew,Anew) %plot mag (linear scale) and phase

Note that lp2lp requires the prototype to have a cutoff frequency of $ω_{c} = 1$ rad/s (besides being an analog lowpass filter).

Example 3.17. Frequency conversion when the prototype does not have a cutoff frequency of 1 rad/s. Assume that a frequency conversion must be applied to $H (s) = 100 ∕ (s + 100)$ such that its cutoff frequency of $ω_{c} = 100$ rad/s is shifted to 30 rad/s. Because this is a simple first-order filter, the final result $G (s) = 30 ∕ (s + 30)$ can be immediately obtained. But the whole procedure is described because it is useful for higher-order filters.

A convenient step is to move the cutoff frequency of $H (s)$ to 1 rad/s using a substitution of $s$ by $100 s$ , which leads to

F (s) = H (100 s) = \frac{100}{100 s + 100} = \frac{1}{s + 1} .

Now, the final filter can be obtained by expanding $ω_{c}$ to 30 by substituting $s$ in $F (s)$ by $s ∕ 30$ :

G (s) = F (s ∕ 30) = \frac{1}{\frac{s}{30} + 1} = \frac{30}{s + 30} .

Alternatively, one could directly use $G (s) = H (100 s ∕ 30)$ . These frequency conversions are widely used in the design of IIR filters as will be discussed. $□$

3.10.6 Filter bandform transformation: Lowpass to highpass, etc.

There are many transformations to convert a prototype (or reference) filter into another one. Section 3.10.5 briefly discussed frequency scaling. It is also possible to convert a lowpass prototype filter into a highpass, bandpass and bandstop, for example, which is called bandform transformation.

When dealing with analog filters, it is often necessary to apply impedance scaling to have the filter input and/or output with a desired impedance. The intent in this case is not to change the transfer function, but to perform impedance matching. All these three kinds of filter transformations are well-discussed in the literature¹⁸. In the sequel, bandform transformation is briefly discussed.

Bandform transformations can be used to create a bandpass from a lowpass prototype, for example. These transformations differ considerably for analog and digital frequencies, especially due to the fact that a digital filter has a periodic spectrum.¹⁹ Having a lowpass analog prototype, a Matlab user²⁰ can use the functions lp2bp, lp2bs and lp2hp to convert it to bandpass, bandstop and highpass filters, respectively.

Sometimes frequency scaling is incorporated to bandform transformation such that a routine converts a lowpass prototype with normalized 1 rad/s cutoff frequency into another filter with desired bandwidth and frequencies. For example, the Matlab’s function lp2bp converts a normalized lowpass into a bandpass with the desired center frequency and bandwidth as illustrated in Listing 3.15.

Listing 3.15: MatlabOnly/snip_systems_lowpass_to_bandpass.m. [ Python version]

1N = 4; %prototype order (bandpass order is 2N = 8) 
2[B,A] = butter(N,1,'s'); %prototype with cutoff=1 rad/s 
3[Bnew,Anew]=lp2bp(B,A,50,30); %convert to cutoff=50 rad/s 
4freqs(Bnew,Anew) %plot mag (linear scale) and phase

The functions lp2bp, lp2bs, lp2hp and lp2lp are restricted to analog filters. Transformations of digital filters can be achieved via a mapping provided by allpass filters, which are filters that have unitary magnitude for all frequencies and just change the phase of the input signal. Besides frequency transformation, allpass filters are adopted for phase equalization and other applications. Some Matlab functions that provide support to digital filter transformations are iirlp2bs, iirlp2bp, iirlp2xn, iirftransf, allpasslp2xn and zpklp2xn.

For example, Listing 3.16 compares the conversion of a Butterworth and a elliptic lowpass prototypes, both with cutoff frequency $Ω_{c} = 0.5$ rad, into stopband versions with cutoff frequencies $Ω_{1} = 0.5$ and $Ω_{2} = 0.8$ rad. The two prototypes are halfband filters because their cutoff frequencies is $0.5$ rad.

Listing 3.16: MatlabOnly/snip_systems_bandstop_conversion.m. [ Python version]

1Wc=0.5; %cutoff (normalized) frequency of prototype filter 
2Wstop = [0.5 0.8]; %stopband cutoff frequencies 
3N=3; %order of prototype filter 
4[B,A] = butter(N,Wc); %Butterworth prototype 
5[Bnew,Anew]=iirlp2bs(B,A,Wc,Wstop); %Bandstop Butterworth 
6[B2, A2] = ellip(N, 3, 30, Wc); %Elliptic prototype 
7[Bnew2,Anew2]=iirlp2bs(B2,A2,Wc,Wstop); %Bandstop elliptic 
8fvtool(Bnew2, Anew2, Bnew, Anew); %compare freq. responses

For the comparison, it was used the Matlab’s Filter Visualization Tool (fvtool). Note that the fourth argument of ellip is the passband frequency (not necessarily the cutoff frequency). In this case, a maximum ripple of 3 dB was specified using the second argument of ellip to have the passband frequency coinciding with the cutoff frequency Wc ( $Ω_{c} = 0.5$ rad). Besides, a minimum attenuation of $30$ dB was required for the stopband of the elliptic lowpass prototype.

As another example, the following Matlab code converts a Butterworth lowpass prototype with cutoff frequency $Ω_{c} = 0.5$ rad into a bandpass with cutoff frequencies $Ω_{1} = 0.25$ and $Ω_{2} = 0.75$ rad.

1N = 5; %filter order 
2[B,A] = butter(N,0.5); %prototype with cutoff=0.5 rad 
3[Bnew,Anew]=iirlp2xn(B,A,[-0.5 0.5],[0.25 0.75],'pass');

The third and fourth arguments of function iirlp2xn ([-0.5 0.5] and [0.25 0.75]) indicate frequency points in the original filter and to what frequencies they should be mapped. In this case, $Ω = - 0.5$ rad should be mapped to $Ω = 0.25$ rad and $Ω = 0.5$ to $Ω = 0.75$ rad.

Instead of creating a prototype and then converting to the desired filter, most routines allow to directly design the filter. The following commands illustrate the design of highpass, bandpass and bandstop elliptic filters.

Listing 3.17: MatlabOctaveCodeSnippets/snip_systems_elliptic_filter_design.m. [ Python version]

1N=5; %filter order for the lowpass prototype (e.g., bandpass is 2*N) 
2Rp=0.1; %maximum ripple in passband (dB) 
3Rs=30; %minimum attenuation in stopband (dB) 
4[B,A]=ellip(N,Rp,Rs,0.4,'high'); %highpass, cutoff=0.4 rad 
5[B,A]=ellip(N,Rp,Rs,[0.5 0.8]); %bandpass, BW=[0.5, 0.8], order=2*N 
6[B,A]=ellip(N,Rp,Rs,[0.5 0.8],'stop'); %BW=[0.5, 0.8], order=2*N

Similarly, FIR design routines provide support to directly designing filters other than lowpass ones. For example, Listing 3.18 designs a stopband filter. Note the normalization by the Nyquist frequency Fs/2.

Listing 3.18: MatlabOctaveCodeSnippets/snip_systems_FIR_filter_design.m. [ Python version]

1Fs = 44100; %sampling freq. (all freqs. in Hz) 
2N  = 100;   %filter order (N+1 coefficientes) 
3Fc1  = 8400;  % First cuttof frequency 
4Fc2  = 13200; % Second cuttof frequency 
5B  = fir1(N, [Fc1 Fc2]/(Fs/2), 'stop', hamming(N+1));

Because there are plenty of routines for bandform transformation, the discussion about IIR and FIR filter design concentrates on lowpass filters.

¹⁴ FilterPro has been deprecated in favor of WEBENCH at [ url3twe].

¹⁵ For the readers without background in electrical circuits but interested in the subject, some basic concepts can be found on the Web, e. g. [ url3ele].

¹⁶ The FilterPro options were a lowpass Butterworth filter with passband frequency of 1 kHz, with a fixed order of 2 and components with 10% of tolerance. When using the exact values of the calculated resistances and capacitances (two of 11.25 K $Ω$ , one of 5.63 K $Ω$ , 10 nF and 40 nF), the cutoff frequency is exactly 1 kHz, but this is typically not feasible in practice, i. e., commercial components are available only within a finite set of values and have non-zero tolerance.

¹⁷ In practice a hardware interrupt is periodically generated according to a timer. The code that implements the filter is part of the interrupt service routine (ISR).

¹⁸ See references in [ url3tra].

¹⁹ Compare Matlab’s documentation for “Digital frequency transformations” and “Analog filter transformation”, currently at [ url3mad] and [ url3maa], respectively.

²⁰ As indicated, the current version of Octave misses this functionality.