Bilinear transform

Short description: Signal processing operation

The bilinear transform (also known as Tustin's method, after Arnold Tustin) is used in digital signal processing and discrete-time control theory to transform continuous-time system representations to discrete-time and vice versa.

The bilinear transform is a special case of a conformal mapping (namely, a Möbius transformation), often used for converting a transfer function [math]\displaystyle{ H_a(s) }[/math] of a linear, time-invariant (LTI) filter in the continuous-time domain (often named an analog filter) to a transfer function [math]\displaystyle{ H_d(z) }[/math] of a linear, shift-invariant filter in the discrete-time domain (often named a digital filter although there are analog filters constructed with switched capacitors that are discrete-time filters). It maps positions on the [math]\displaystyle{ j \omega }[/math] axis, [math]\displaystyle{ \mathrm{Re}[s]=0 }[/math], in the s-plane to the unit circle, [math]\displaystyle{ |z| = 1 }[/math], in the z-plane. Other bilinear transforms can be used for warping the frequency response of any discrete-time linear system (for example to approximate the non-linear frequency resolution of the human auditory system) and are implementable in the discrete domain by replacing a system's unit delays [math]\displaystyle{ \left( z^{-1} \right) }[/math] with first order all-pass filters.

The transform preserves stability and maps every point of the frequency response of the continuous-time filter, [math]\displaystyle{ H_a(j \omega_a) }[/math] to a corresponding point in the frequency response of the discrete-time filter, [math]\displaystyle{ H_d(e^{j \omega_d T}) }[/math] although to a somewhat different frequency, as shown in the Frequency warping section below. This means that for every feature that one sees in the frequency response of the analog filter, there is a corresponding feature, with identical gain and phase shift, in the frequency response of the digital filter but, perhaps, at a somewhat different frequency. This is barely noticeable at low frequencies but is quite evident at frequencies close to the Nyquist frequency.

Discrete-time approximation

The bilinear transform is a first-order Padé approximant of the natural logarithm function that is an exact mapping of the z-plane to the s-plane. When the Laplace transform is performed on a discrete-time signal (with each element of the discrete-time sequence attached to a correspondingly delayed unit impulse), the result is precisely the Z transform of the discrete-time sequence with the substitution of

[math]\displaystyle{ \begin{align} z &= e^{sT} \\ &= \frac{e^{sT/2}}{e^{-sT/2}} \\ &\approx \frac{1 + s T / 2}{1 - s T / 2} \end{align} }[/math]

where [math]\displaystyle{ T }[/math] is the numerical integration step size of the trapezoidal rule used in the bilinear transform derivation;^[1] or, in other words, the sampling period. The above bilinear approximation can be solved for [math]\displaystyle{ s }[/math] or a similar approximation for [math]\displaystyle{ s = (1/T) \ln(z) }[/math] can be performed.

The inverse of this mapping (and its first-order bilinear approximation) is

[math]\displaystyle{ \begin{align} s &= \frac{1}{T} \ln(z) \\ &= \frac{2}{T} \left[\frac{z-1}{z+1} + \frac{1}{3} \left( \frac{z-1}{z+1} \right)^3 + \frac{1}{5} \left( \frac{z-1}{z+1} \right)^5 + \frac{1}{7} \left( \frac{z-1}{z+1} \right)^7 + \cdots \right] \\ &\approx \frac{2}{T} \frac{z - 1}{z + 1} \\ &= \frac{2}{T} \frac{1 - z^{-1}}{1 + z^{-1}} \end{align} }[/math]

The bilinear transform essentially uses this first order approximation and substitutes into the continuous-time transfer function, [math]\displaystyle{ H_a(s) }[/math]

[math]\displaystyle{ s \leftarrow \frac{2}{T} \frac{z - 1}{z + 1}. }[/math]

That is

[math]\displaystyle{ H_d(z) = H_a(s) \bigg|_{s = \frac{2}{T} \frac{z - 1}{z + 1}}= H_a \left( \frac{2}{T} \frac{z-1}{z+1} \right). \ }[/math]

Stability and minimum-phase property preserved

A continuous-time causal filter is stable if the poles of its transfer function fall in the left half of the complex s-plane. A discrete-time causal filter is stable if the poles of its transfer function fall inside the unit circle in the complex z-plane. The bilinear transform maps the left half of the complex s-plane to the interior of the unit circle in the z-plane. Thus, filters designed in the continuous-time domain that are stable are converted to filters in the discrete-time domain that preserve that stability.

Likewise, a continuous-time filter is minimum-phase if the zeros of its transfer function fall in the left half of the complex s-plane. A discrete-time filter is minimum-phase if the zeros of its transfer function fall inside the unit circle in the complex z-plane. Then the same mapping property assures that continuous-time filters that are minimum-phase are converted to discrete-time filters that preserve that property of being minimum-phase.

Transformation of a General LTI System

A general LTI system has the transfer function [math]\displaystyle{ H_a(s) = \frac{b_0 + b_1s + b_2s^2 + \cdots + b_Qs^Q}{a_0 + a_1s + a_2s^2 + \cdots + a_Ps^P} }[/math] The order of the transfer function $N$ is the greater of $P$ and $Q$ (in practice this is most likely $P$ as the transfer function must be proper for the system to be stable). Applying the bilinear transform [math]\displaystyle{ s = K\frac{z - 1}{z + 1} }[/math] where $K$ is defined as either $2/ T$ or otherwise if using frequency warping, gives [math]\displaystyle{ H_d(z) = \frac{b_0 + b_1\left(K\frac{z - 1}{z + 1}\right) + b_2\left(K\frac{z - 1}{z + 1}\right)^2 + \cdots + b_Q\left(K\frac{z - 1}{z + 1}\right)^Q} {a_0 + a_1\left(K\frac{z - 1}{z + 1}\right) + a_2\left(K\frac{z - 1}{z + 1}\right)^2 + \cdots + b_P\left(K\frac{z - 1}{z + 1}\right)^P} }[/math] Multiplying the numerator and denominator by the largest power of $(z + 1) -1$ present, $(z + 1) -N$ , gives [math]\displaystyle{ H_d(z) = \frac{b_0(z+1)^N + b_1K(z-1)(z+1)^{N-1} + b_2K^2(z-1)^2(z+1)^{N-2} + \cdots + b_QK^Q(z-1)^Q(z+1)^{N-Q}} {a_0(z+1)^N + a_1K(z-1)(z+1)^{N-1} + a_2K^2(z-1)^2(z+1)^{N-2} + \cdots + a_PK^P(z-1)^P(z+1)^{N-P}} }[/math] It can be seen here that after the transformation, the degree of the numerator and denominator are both $N$ .

Consider then the pole-zero form of the continuous-time transfer function [math]\displaystyle{ H_a(s) = \frac{(s - \xi_1)(s - \xi_2) \cdots (s - \xi_Q)}{(s - p_1)(s - p_2) \cdots (s - p_P)} }[/math] The roots of the numerator and denominator polynomials, $ξ i$ and $p i$ , are the zeros and poles of the system. The bilinear transform is a one-to-one mapping, hence these can be transformed to the z-domain using [math]\displaystyle{ z = \frac{K + s}{K - s} }[/math] yielding some of the discretized transfer function's zeros and poles $ξ' i$ and $p' i$ [math]\displaystyle{ \begin{aligned} \xi'_i &= \frac{K + \xi_i}{K - \xi_i} \quad 1 \leq i \leq Q \\ p'_i &= \frac{K + p_i}{K - p_i} \quad 1 \leq i \leq P \end{aligned} }[/math] As described above, the degree of the numerator and denominator are now both $N$ , in other words there is now an equal number of zeros and poles. The multiplication by $(z + 1) -N$ means the additional zeros or poles are ^[2] [math]\displaystyle{ \begin{aligned} \xi'_i &= -1 \quad Q \lt i \leq N \\ p'_i &= -1 \quad P \lt i \leq N \end{aligned} }[/math] Given the full set of zeros and poles, the z-domain transfer function is then [math]\displaystyle{ H_d(z) = \frac{(z - \xi'_1)(z - \xi'_2) \cdots (z - \xi'_N)} {(z - p'_1)(z - p'_2) \cdots (z - p'_N)} }[/math]

Example

As an example take a simple low-pass RC filter. This continuous-time filter has a transfer function

[math]\displaystyle{ \begin{align} H_a(s) &= \frac{1/sC}{R+1/sC} \\ &= \frac{1}{1 + RC s}. \end{align} }[/math]

If we wish to implement this filter as a digital filter, we can apply the bilinear transform by substituting for [math]\displaystyle{ s }[/math] the formula above; after some reworking, we get the following filter representation:

[math]\displaystyle{ H_d(z) \ }[/math]	[math]\displaystyle{ =H_a \left( \frac{2}{T} \frac{z-1}{z+1}\right) \ }[/math]
	[math]\displaystyle{ = \frac{1}{1 + RC \left( \frac{2}{T} \frac{z-1}{z+1}\right)} \ }[/math]
	[math]\displaystyle{ = \frac{1 + z}{(1 - 2 RC / T) + (1 + 2RC / T) z} \ }[/math]
	[math]\displaystyle{ = \frac{1 + z^{-1}}{(1 + 2RC / T) + (1 - 2RC / T) z^{-1}}. \ }[/math]

The coefficients of the denominator are the 'feed-backward' coefficients and the coefficients of the numerator are the 'feed-forward' coefficients used for implementing a real-time digital filter.

Transformation for a general first-order continuous-time filter

It is possible to relate the coefficients of a continuous-time, analog filter with those of a similar discrete-time digital filter created through the bilinear transform process. Transforming a general, first-order continuous-time filter with the given transfer function

[math]\displaystyle{ H_a(s) = \frac{b_0 s + b_1}{a_0 s + a_1} = \frac{b_0 + b_1 s^{-1}}{a_0 + a_1 s^{-1}} }[/math]

using the bilinear transform (without prewarping any frequency specification) requires the substitution of

[math]\displaystyle{ s \leftarrow K \frac{1 - z^{-1}}{1 + z^{-1}} }[/math]

where

[math]\displaystyle{ K \triangleq \frac{2}{T} }[/math].

However, if the frequency warping compensation as described below is used in the bilinear transform, so that both analog and digital filter gain and phase agree at frequency [math]\displaystyle{ \omega_0 }[/math], then

[math]\displaystyle{ K \triangleq \frac{\omega_0}{\tan\left(\frac{\omega_0 T}{2}\right)} }[/math].

This results in a discrete-time digital filter with coefficients expressed in terms of the coefficients of the original continuous time filter:

[math]\displaystyle{ H_d(z)=\frac{(b_0 K + b_1) + (-b_0 K + b_1)z^{-1}}{(a_0 K + a_1) + (-a_0 K + a_1)z^{-1}} }[/math]

Normally the constant term in the denominator must be normalized to 1 before deriving the corresponding difference equation. This results in

[math]\displaystyle{ H_d(z)=\frac{\frac{b_0 K + b_1}{a_0 K + a_1} + \frac{-b_0 K + b_1}{a_0 K + a_1}z^{-1}}{1 + \frac{-a_0 K + a_1}{a_0 K + a_1}z^{-1}}. }[/math]

The difference equation (using the Direct form I) is

[math]\displaystyle{ y[n] = \frac{b_0 K + b_1}{a_0 K + a_1} \cdot x[n] + \frac{-b_0 K + b_1}{a_0 K + a_1} \cdot x[n-1] - \frac{-a_0 K + a_1}{a_0 K + a_1} \cdot y[n-1] \ . }[/math]

General second-order biquad transformation

A similar process can be used for a general second-order filter with the given transfer function

[math]\displaystyle{ H_a(s) = \frac{b_0 s^2 + b_1 s + b_2}{a_0 s^2 + a_1 s + a_2} = \frac{b_0 + b_1 s^{-1} + b_2 s^{-2}}{a_0 + a_1 s^{-1} + a_2 s^{-2}} \ . }[/math]

This results in a discrete-time digital biquad filter with coefficients expressed in terms of the coefficients of the original continuous time filter:

[math]\displaystyle{ H_d(z)=\frac{(b_0 K^2 + b_1 K + b_2) + (2b_2 - 2b_0 K^2)z^{-1} + (b_0 K^2 - b_1 K + b_2)z^{-2}}{(a_0 K^2 + a_1 K + a_2) + (2a_2 - 2a_0 K^2)z^{-1} + (a_0 K^2 - a_1 K + a_2)z^{-2}} }[/math]

Again, the constant term in the denominator is generally normalized to 1 before deriving the corresponding difference equation. This results in

[math]\displaystyle{ H_d(z)=\frac{\frac{b_0 K^2 + b_1 K + b_2}{a_0 K^2 + a_1 K + a_2} + \frac{2b_2 - 2b_0 K^2}{a_0 K^2 + a_1 K + a_2}z^{-1} + \frac{b_0 K^2 - b_1 K + b_2}{a_0 K^2 + a_1 K + a_2}z^{-2}}{1 + \frac{2a_2 - 2a_0 K^2}{a_0 K^2 + a_1 K + a_2}z^{-1} + \frac{a_0 K^2 - a_1 K + a_2}{a_0 K^2 + a_1 K + a_2}z^{-2}}. }[/math]

The difference equation (using the Direct form I) is

[math]\displaystyle{ y[n] = \frac{b_0 K^2 + b_1 K + b_2}{a_0 K^2 + a_1 K + a_2} \cdot x[n] + \frac{2b_2 - 2b_0 K^2}{a_0 K^2 + a_1 K + a_2} \cdot x[n-1] + \frac{b_0 K^2 - b_1 K + b_2}{a_0 K^2 + a_1 K + a_2} \cdot x[n-2] - \frac{2a_2 - 2a_0 K^2}{a_0 K^2 + a_1 K + a_2} \cdot y[n-1] - \frac{a_0 K^2 - a_1 K + a_2}{a_0 K^2 + a_1 K + a_2} \cdot y[n-2] \ . }[/math]

Frequency warping

To determine the frequency response of a continuous-time filter, the transfer function [math]\displaystyle{ H_a(s) }[/math] is evaluated at [math]\displaystyle{ s = j \omega_a }[/math] which is on the [math]\displaystyle{ j \omega }[/math] axis. Likewise, to determine the frequency response of a discrete-time filter, the transfer function [math]\displaystyle{ H_d(z) }[/math] is evaluated at [math]\displaystyle{ z = e^{ j \omega_d T} }[/math] which is on the unit circle, [math]\displaystyle{ |z| = 1 }[/math]. The bilinear transform maps the [math]\displaystyle{ j \omega }[/math] axis of the s-plane (which is the domain of [math]\displaystyle{ H_a(s) }[/math]) to the unit circle of the z-plane, [math]\displaystyle{ |z| = 1 }[/math] (which is the domain of [math]\displaystyle{ H_d(z) }[/math]), but it is not the same mapping [math]\displaystyle{ z = e^{sT} }[/math] which also maps the [math]\displaystyle{ j \omega }[/math] axis to the unit circle. When the actual frequency of [math]\displaystyle{ \omega_d }[/math] is input to the discrete-time filter designed by use of the bilinear transform, then it is desired to know at what frequency, [math]\displaystyle{ \omega_a }[/math], for the continuous-time filter that this [math]\displaystyle{ \omega_d }[/math] is mapped to.

[math]\displaystyle{ H_d(z) = H_a \left( \frac{2}{T} \frac{z-1}{z+1}\right) }[/math]

[math]\displaystyle{ H_d(e^{ j \omega_d T}) }[/math]	[math]\displaystyle{ = H_a \left( \frac{2}{T} \frac{e^{ j \omega_d T} - 1}{e^{ j \omega_d T} + 1}\right) }[/math]
	[math]\displaystyle{ = H_a \left( \frac{2}{T} \cdot \frac{e^{j \omega_d T/2} \left(e^{j \omega_d T/2} - e^{-j \omega_d T/2}\right)}{e^{j \omega_d T/2} \left(e^{j \omega_d T/2} + e^{-j \omega_d T/2 }\right)}\right) }[/math]
	[math]\displaystyle{ = H_a \left( \frac{2}{T} \cdot \frac{\left(e^{j \omega_d T/2} - e^{-j \omega_d T/2}\right)}{\left(e^{j \omega_d T/2} + e^{-j \omega_d T/2 }\right)}\right) }[/math]
	[math]\displaystyle{ = H_a \left(j \frac{2}{T} \cdot \frac{ \left(e^{j \omega_d T/2} - e^{-j \omega_d T/2}\right) /(2j)}{\left(e^{j \omega_d T/2} + e^{-j \omega_d T/2 }\right) / 2}\right) }[/math]
	[math]\displaystyle{ = H_a \left(j \frac{2}{T} \cdot \frac{ \sin(\omega_d T/2) }{ \cos(\omega_d T/2) }\right) }[/math]
	[math]\displaystyle{ = H_a \left(j \frac{2}{T} \cdot \tan \left( \omega_d T/2 \right) \right) }[/math]

This shows that every point on the unit circle in the discrete-time filter z-plane, [math]\displaystyle{ z = e^{ j \omega_d T} }[/math] is mapped to a point on the [math]\displaystyle{ j \omega }[/math] axis on the continuous-time filter s-plane, [math]\displaystyle{ s = j \omega_a }[/math]. That is, the discrete-time to continuous-time frequency mapping of the bilinear transform is

[math]\displaystyle{ \omega_a = \frac{2}{T} \tan \left( \omega_d \frac{T}{2} \right) }[/math]

and the inverse mapping is

[math]\displaystyle{ \omega_d = \frac{2}{T} \arctan \left( \omega_a \frac{T}{2} \right). }[/math]

The discrete-time filter behaves at frequency [math]\displaystyle{ \omega_d }[/math] the same way that the continuous-time filter behaves at frequency [math]\displaystyle{ (2/T) \tan(\omega_d T/2) }[/math]. Specifically, the gain and phase shift that the discrete-time filter has at frequency [math]\displaystyle{ \omega_d }[/math] is the same gain and phase shift that the continuous-time filter has at frequency [math]\displaystyle{ (2/T) \tan(\omega_d T/2) }[/math]. This means that every feature, every "bump" that is visible in the frequency response of the continuous-time filter is also visible in the discrete-time filter, but at a different frequency. For low frequencies (that is, when [math]\displaystyle{ \omega_d \ll 2/T }[/math] or [math]\displaystyle{ \omega_a \ll 2/T }[/math]), then the features are mapped to a slightly different frequency; [math]\displaystyle{ \omega_d \approx \omega_a }[/math].

One can see that the entire continuous frequency range

[math]\displaystyle{ -\infty \lt \omega_a \lt +\infty }[/math]

is mapped onto the fundamental frequency interval

[math]\displaystyle{ -\frac{\pi}{T} \lt \omega_d \lt +\frac{\pi}{T}. }[/math]

The continuous-time filter frequency [math]\displaystyle{ \omega_a = 0 }[/math] corresponds to the discrete-time filter frequency [math]\displaystyle{ \omega_d = 0 }[/math] and the continuous-time filter frequency [math]\displaystyle{ \omega_a = \pm \infty }[/math] correspond to the discrete-time filter frequency [math]\displaystyle{ \omega_d = \pm \pi / T. }[/math]

One can also see that there is a nonlinear relationship between [math]\displaystyle{ \omega_a }[/math] and [math]\displaystyle{ \omega_d. }[/math] This effect of the bilinear transform is called frequency warping. The continuous-time filter can be designed to compensate for this frequency warping by setting [math]\displaystyle{ \omega_a = \frac{2}{T} \tan \left( \omega_d \frac{T}{2} \right) }[/math] for every frequency specification that the designer has control over (such as corner frequency or center frequency). This is called pre-warping the filter design.

It is possible, however, to compensate for the frequency warping by pre-warping a frequency specification [math]\displaystyle{ \omega_0 }[/math] (usually a resonant frequency or the frequency of the most significant feature of the frequency response) of the continuous-time system. These pre-warped specifications may then be used in the bilinear transform to obtain the desired discrete-time system. When designing a digital filter as an approximation of a continuous time filter, the frequency response (both amplitude and phase) of the digital filter can be made to match the frequency response of the continuous filter at a specified frequency [math]\displaystyle{ \omega_0 }[/math], as well as matching at DC, if the following transform is substituted into the continuous filter transfer function.^[3] This is a modified version of Tustin's transform shown above.

[math]\displaystyle{ s \leftarrow \frac{\omega_0}{\tan\left(\frac{\omega_0 T}{2}\right)} \frac{z - 1}{z + 1}. }[/math]

However, note that this transform becomes the original transform

[math]\displaystyle{ s \leftarrow \frac{2}{T} \frac{z - 1}{z + 1} }[/math]

as [math]\displaystyle{ \omega_0 \to 0 }[/math].

The main advantage of the warping phenomenon is the absence of aliasing distortion of the frequency response characteristic, such as observed with Impulse invariance.

References

↑ Oppenheim, Alan (2010). Discrete Time Signal Processing Third Edition. Upper Saddle River, NJ: Pearson Higher Education, Inc.. p. 504. ISBN 978-0-13-198842-2.
↑ Bhandari, Ayush. "DSP and Digital Filters Lecture Notes". http://www.ee.ic.ac.uk/hp/staff/dmb/courses/DSPDF/00800_TransIIR.pdf.
↑ Astrom, Karl J. (1990). Computer Controlled Systems, Theory and Design (Second ed.). Prentice-Hall. p. 212. ISBN 0-13-168600-3.

External links

0.00

(0 votes)

Original source: https://en.wikipedia.org/wiki/Bilinear transform. Read more

[1] Oppenheim, Alan (2010). Discrete Time Signal Processing Third Edition. Upper Saddle River, NJ: Pearson Higher Education, Inc.. p. 504. ISBN 978-0-13-198842-2.

[2] Bhandari, Ayush. "DSP and Digital Filters Lecture Notes". http://www.ee.ic.ac.uk/hp/staff/dmb/courses/DSPDF/00800_TransIIR.pdf.

[3] Astrom, Karl J. (1990). Computer Controlled Systems, Theory and Design (Second ed.). Prentice-Hall. p. 212. ISBN 0-13-168600-3.

[1]

[2]

[3]

v t e Digital signal processing
Theory	Detection theory Discrete signal Estimation theory Nyquist–Shannon sampling theorem
Sub-fields	Audio signal processing Digital image processing Speech processing Statistical signal processing
Techniques	Advanced Z-transform Bilinear transform Constant-Q transform Discrete cosine transform (DCT) Discrete Fourier transform (DFT) Discrete-time Fourier transform (DTFT) Impulse invariance Integral transform Laplace transform Matched Z-transform method Post's inversion formula Starred transform Z-transform Zak transform
Sampling	Aliasing Anti-aliasing filter Downsampling Nyquist rate / frequency Oversampling Quantization Sampling rate Undersampling Upsampling

Anonymous

Search

Bilinear transform

Namespaces

More

Page actions

Contents

Discrete-time approximation

Stability and minimum-phase property preserved

Transformation of a General LTI System

Example

Transformation for a general first-order continuous-time filter

General second-order biquad transformation

Frequency warping

See also

References

External links

Navigation

Navigation

Help

Translate

Wiki tools

Wiki tools

Anonymous

Search

Bilinear transform

Discrete-time approximation

Stability and minimum-phase property preserved

Transformation of a General LTI System

Example

Transformation for a general first-order continuous-time filter

General second-order biquad transformation

Frequency warping

See also

References

External links

Navigation

Wiki tools

Page tools

Other projects

Categories