# Bilinear Transformation Comparison

If I have transfer function coefficients,
I can analyze the transfer function in the s-plane and/or the z-plane.

If I wanted to demonstrate that the z-plane and s-plane responses are equivalent:

Would it be possible to use the frequency warping relation
(provided in the otherwise unused variable w.z in in the code)
to depict the z-plane plots in the s-plane or vice versa
in a manner in which the two would completely match?

## Matlab code provided:

[Actual filter coefficients not provided.]

%% Input

A.s.cf = [ 0 0 0 0 1 ];
B.s.cf = [ 1 2 3 4 5 ];

T.z   = 2.5e-05;
f_z   = 1/T.z;
f_nyq = f_z / 2;

%% Processing: Matlab SP toolbox
%{
H.s.tf = tf ( A.s.cf, B.s.cf       );
H.z.tf = c2d( H.s.tf, T.z,   'zoh' );
%}

%% Processing: From scratch:
f.s = logspace(-2, +6, 5e3).';
w.s = 2*pi*f.s;
s.s = 1j * w.s;

z.z = exp( 1j * w.s * T.z);

w.z = 2/T.z * atan( w.s * T.z/2 );

switch 0
case 0; s2z = 1/T.z * log(z.z);
case 1; s2z = 2/T.z * (z.z-1)./(z.z+1);
end

A.s.jw = s.s.^4 * A.s.cf(1) ...
+ s.s.^3 * A.s.cf(2) ...
+ s.s.^2 * A.s.cf(3) ...
+ s.s.^1 * A.s.cf(4) ...
+ s.s.^0 * A.s.cf(5) ; %

B.s.jw = s.s.^4 * B.s.cf(1) ...
+ s.s.^3 * B.s.cf(2) ...
+ s.s.^2 * B.s.cf(3) ...
+ s.s.^1 * B.s.cf(4) ...
+ s.s.^0 * B.s.cf(5) ; %

A.z.jw = s2z.^4 * A.s.cf(1) ...
+ s2z.^3 * A.s.cf(2) ...
+ s2z.^2 * A.s.cf(3) ...
+ s2z.^1 * A.s.cf(4) ...
+ s2z.^0 * A.s.cf(5) ; %

B.z.jw = s2z.^4 * B.s.cf(1) ...
+ s2z.^3 * B.s.cf(2) ...
+ s2z.^2 * B.s.cf(3) ...
+ s2z.^1 * B.s.cf(4) ...
+ s2z.^0 * B.s.cf(5) ; %

H.s.jw = A.s.jw ./ B.s.jw;
H.z.jw = A.z.jw ./ B.z.jw;


If I plot the two,
the response will be different in either plane,
(as expected).

%% [Plot]: Setup
p = 0;

%% [Plot]: S plane
p = p+1;
figure(p)

% S plane: Magnitude
subplot(2,1,1)
semilogx(w.s, 20 * log10( abs(H.s.jw) ), '.')

title('Bode [From Scratch]')
ylabel('Mag [dB]')

axis( [ 1e+3 1e+6 -150 +0] )
grid minor

% S Plane: Phase
subplot(2,1,2)
semilogx(w.s, angle( H.s.jw ) * 180/pi, '-')
hold on
semilogx( [ 2*pi*1e-2 2*pi*1e+6 ], [ -180 -180 ], 'k')
semilogx( [ 2*pi*1e-2 2*pi*1e+6 ], [ +180 +180 ], 'k')
hold off

ylabel('Ang [deg]')

axis( [ 1e+3 1e+6 -200 +200] )
grid minor

%% [Plot]: Z Plane
p = p+1;
figure(p)

% Z Plane: Magnitude
subplot(2,1,1)
semilogx(w.z, 20 * log10( abs(H.z.jw) ), '.')
hold on
semilogx([2*pi*f_nyq 2*pi*f_nyq], [-50 +50], 'k-')
hold off

title('Bode [From Scratch]')
ylabel('Mag [dB]')

axis( [1e+3 1e+6 -50 0] )
grid minor

% Z Plane: Phase
subplot(2,1,2)
semilogx(w.z, angle( H.z.jw ) * 180/pi, '-')
hold on
semilogx( [ 2*pi*1e-2  2*pi*1e+6 ], [ -180 -180 ], 'k')
semilogx( [ 2*pi*1e-2  2*pi*1e+6 ], [ +180 +180 ], 'k')
semilogx( [ 2*pi*f_nyq 2*pi*f_nyq], [ -180 +180 ], 'k-')
hold off

ylabel('Ang [deg]')

axis( [1e+3 1e+6 -200 +200] )
grid minor

%% [Plot]: S Plane: Bode

%{
p = p+1;
figure(p)

bode(H.s.tf)
grid minor
%}

%% [Plot]: Z Plane: Bode

%{
p = p+1;
figure(p)

bode(H.z.tf)
grid minor
%}


i really don't wanna slog through your MATLAB code. is this what you're trying to prove? :

let $H_\text{a}(s)$ be some analog (or "continuous-time") LTI transfer function on the s-plane and $H_\text{d}(z)$ be some digital (or "discrete-time") LTI transfer function in the z-plane.

$$H_\text{d}(z) = H_\text{a}(s) \Bigg|_{s=\frac{2}{T}\frac{z-1}{z+1}}$$

do you want to prove that

$$H_\text{d}(e^{j\omega}) = H_\text{a}(j\Omega) \Bigg|_{\Omega=\frac{2}{T}\tan\left(\frac{\omega}{2}\right)}$$

is that what you want proven?

• Plus 1 for not wanting to slog through code! – Dan Boschen Apr 26 '17 at 19:54