This articleneeds additional citations forverification. Please helpimprove this article byadding citations to reliable sources. Unsourced material may be challenged and removed. Find sources: "Bilinear transform" – news ·newspapers ·books ·scholar ·JSTOR(June 2009) (Learn how and when to remove this message)

Thebilinear transform (also known asTustin's method, afterArnold Tustin) is used indigital signal processing and discrete-timecontrol theory to transform continuous-time system representations to discrete-time and vice versa.

The bilinear transform is a special case of aconformal mapping (namely, aMöbius transformation), often used for converting atransfer function${\displaystyle H_a(s)}$ of alinear,time-invariant (LTI) filter in thecontinuous-time domain (often named ananalog filter) to a transfer function${\displaystyle H_d(z)}$ of a linear, shift-invariant filter in thediscrete-time domain (often named adigital filter although there are analog filters constructed withswitched capacitors that are discrete-time filters). It maps positions on the${\displaystyle j\omega }$ axis,${\displaystyle \mathrm{R}\mathrm{e}[s]=0}$, in thes-plane to theunit circle,${\displaystyle |z|=1}$, in thez-plane. Other bilinear transforms can be used for warping thefrequency 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${\displaystyle \left(z^{-1}\right)}$ with first orderall-pass filters.

The transform preservesstability and maps every point of thefrequency response of the continuous-time filter,${\displaystyle H_a(j{\omega }_a)}$ to a corresponding point in the frequency response of the discrete-time filter,${\displaystyle H_d(e^{j{\omega }_{d}T})}$ although to a somewhat different frequency, as shown in theFrequency 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. The change in frequency is barely noticeable at low frequencies but is quite evident at frequencies close to theNyquist frequency.

The bilinear transform is a first-orderPadé approximant of the natural logarithm function that is an exact mapping of thez-plane to thes-plane. When theLaplace transform is performed on a discrete-time signal (with each element of the discrete-time sequence attached to a correspondingly delayedunit impulse), the result is precisely theZ transform of the discrete-time sequence with the substitution of

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

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

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

${\displaystyle s\leftarrow \frac{2}{T}\frac{z-1}{z+1}.}$

That is

${\displaystyle H_d(z)=H_a(s){|}_{s=\frac{2}{T}\frac{z-1}{z+1}}=H_a\left(\frac{2}{T}\frac{z-1}{z+1}\right).}$

A continuous-timecausal filter isstable if thepoles of its transfer function fall in the left half of thecomplexs-plane. A discrete-time causal filter is stable if the poles of its transfer function fall inside theunit circle in thecomplex 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 isminimum-phase if thezeros 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.

A generalLTI system has the transfer function$${\displaystyle H_a(s)=\frac{b_0+b_{1}s+b_2{s}^2+\cdots +b_Q{s}^Q}{a_0+a_{1}s+a_2{s}^2+\cdots +a_P{s}^P}}$$The order of the transfer functionN is the greater ofP andQ (in practice this is most likelyP as the transfer function must beproper for the system to be stable). Applying the bilinear transform$${\displaystyle s=K\frac{z-1}{z+1}}$$whereK is defined as either2/T or otherwise if usingfrequency warping, gives$${\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}}$$Multiplying the numerator and denominator by the largest power of(z + 1)⁻¹ present,(z + 1)^−N, gives$${\displaystyle H_d(z)=\frac{b_0(z+1{)}^N+b_{1}K(z-1)(z+1{)}^{N-1}+b_2{K}^2(z-1{)}^2(z+1{)}^{N-2}+\cdots +b_Q{K}^Q(z-1{)}^Q(z+1{)}^{N-Q}}{a_0(z+1{)}^N+a_{1}K(z-1)(z+1{)}^{N-1}+a_2{K}^2(z-1{)}^2(z+1{)}^{N-2}+\cdots +a_P{K}^P(z-1{)}^P(z+1{)}^{N-P}}}$$It can be seen here that after the transformation, the degree of the numerator and denominator are bothN.

Consider then the pole-zero form of the continuous-time transfer function$${\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)}}$$The roots of the numerator and denominator polynomials,ξ_i andp_i, are thezeros and poles of the system. The bilinear transform is aone-to-one mapping, hence these can be transformed to the z-domain using$${\displaystyle z=\frac{K+s}{K-s}}$$yielding some of the discretized transfer function's zeros and polesξ'_i andp'_iAs described above, the degree of the numerator and denominator are now bothN, 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]Given the full set of zeros and poles, the z-domain transfer function is then$${\displaystyle H_d(z)=\frac{(z-{\xi }_1^{\prime })(z-{\xi }_2^{\prime })\cdots (z-{\xi }_N^{\prime })}{(z-p_1^{\prime })(z-p_2^{\prime })\cdots (z-p_N^{\prime })}}$$

As an example take a simplelow-passRC filter. This continuous-time filter has a transfer function

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

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

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-timedigital 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

${\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}}}$

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

${\displaystyle s\leftarrow K\frac{1-z^{-1}}{1+z^{-1}}}$

where

${\displaystyle K\triangleq \frac{2}{T}}$.

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${\displaystyle {\omega }_0}$, then

${\displaystyle K\triangleq \frac{{\omega }_0}{\tan \left(\frac{{\omega }_{0}T}{2}\right)}}$.

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

${\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}}}$

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

${\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}}.}$

The difference equation (using theDirect form I) is

${\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].}$

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

${\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}}.}$

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

${\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}}}$

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

${\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}}.}$

The difference equation (using theDirect form I) is

${\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].}$

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

${\displaystyle H_d(z)=H_a\left(\frac{2}{T}\frac{z-1}{z+1}\right)}$

${\displaystyle H_d(e^{j{\omega }_{d}T})}$	${\displaystyle =H_a\left(\frac{2}{T}\frac{e^{j{\omega }_{d}T}-1}{e^{j{\omega }_{d}T}+1}\right)}$
	${\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)}$
	${\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)}$
	${\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)}$
	${\displaystyle =H_a\left(j\frac{2}{T}\cdot \frac{\sin ({\omega }_{d}T/2)}{\cos ({\omega }_{d}T/2)}\right)}$
	${\displaystyle =H_a\left(j\frac{2}{T}\cdot \tan \left({\omega }_{d}T/2\right)\right)}$

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

${\displaystyle {\omega }_a=\frac{2}{T}\tan \left({\omega }_d\frac{T}{2}\right)}$

and the inverse mapping is

${\displaystyle {\omega }_d=\frac{2}{T}\arctan \left({\omega }_a\frac{T}{2}\right).}$

The discrete-time filter behaves at frequency${\displaystyle {\omega }_d}$ the same way that the continuous-time filter behaves at frequency${\displaystyle (2/T)\tan ({\omega }_{d}T/2)}$. Specifically, the gain and phase shift that the discrete-time filter has at frequency${\displaystyle {\omega }_d}$ is the same gain and phase shift that the continuous-time filter has at frequency${\displaystyle (2/T)\tan ({\omega }_{d}T/2)}$. 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${\displaystyle {\omega }_d\ll 2/T}$ or${\displaystyle {\omega }_a\ll 2/T}$), then the features are mapped to aslightly different frequency;${\displaystyle {\omega }_d\approx {\omega }_a}$.

One can see that the entire continuous frequency range

is mapped onto the fundamental frequency interval

The continuous-time filter frequency${\displaystyle {\omega }_a=0}$ corresponds to the discrete-time filter frequency${\displaystyle {\omega }_d=0}$ and the continuous-time filter frequency${\displaystyle {\omega }_a=\pm \mathrm{\infty }}$ correspond to the discrete-time filter frequency${\displaystyle {\omega }_d=\pm \pi /T.}$

One can also see that there is a nonlinear relationship between${\displaystyle {\omega }_a}$ and${\displaystyle {\omega }_d.}$ This effect of the bilinear transform is calledfrequency warping. The continuous-time filter can be designed to compensate for this frequency warping by setting${\displaystyle {\omega }_a=\frac{2}{T}\tan \left({\omega }_d\frac{T}{2}\right)}$ for every frequency specification that the designer has control over (such as corner frequency or center frequency). This is calledpre-warping thefilter design.

It is possible, however, to compensate for the frequency warping by pre-warping a frequency specification${\displaystyle {\omega }_0}$ (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${\displaystyle {\omega }_0}$, 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.

${\displaystyle s\leftarrow \frac{{\omega }_0}{\tan \left(\frac{{\omega }_{0}T}{2}\right)}\frac{z-1}{z+1}.}$

However, note that this transform becomes the original transform

${\displaystyle s\leftarrow \frac{2}{T}\frac{z-1}{z+1}}$

as${\displaystyle {\omega }_0\to 0}$.

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

• Impulse invariance
• Matched Z-transform method

• MIT OpenCourseWare Signal Processing: Continuous to Discrete Filter Design
• Lecture Notes on Discrete Equivalents
• The Art of VA Filter Design

• Digital signal processing
• Transforms
• Control theory

• Articles with short description
• Short description is different from Wikidata
• Articles needing additional references from June 2009
• All articles needing additional references