Note: type "help command" to get help on matlab commands, i.e., "help butter."
New:
You may also wish to experiment with the new
Matlab "Filter design and analysis tool." To run this,
type "fdatool" at the Matlab prompt. Use Analysis->FullView to rescale
Y-axis.
matlab
< M A T L A B (R) >
(c) Copyright 1984-94 The MathWorks, Inc.
All Rights Reserved
Version 4.2c
Dec 31 1994
Commands to get started: intro, demo, help help
Commands for more information: help, whatsnew, info, subscribe
>> fsam=4000
>> [z,p,k]=butter(2,6283,'s')
z = (NOTE: no finite zeroes in butterworth)
[]
p = (NOTE:this is location of the 2 poles in s-plane)
1.0e+03 *
-4.4428 + 4.4428i
-4.4428 - 4.4428i
k =
39476089 (NOTE: this is a constant multiplying the result)
(==========================)
(<<< BEWARE MATLAB BUG!!!)
(==========================)
Warning: There is an error
in the 5.1.0.421 version
of matlab. There should be no zeroes in
the returned value of z. To fix this,
you must suppress the zeroes by
entering the command "z=[]" after
In the current version you will get
z=[ 0 0 0.2500 ] which is wrong!
For the new (broken) version, type the following to correct the problem:
>> z=[]
z =
[]
>> [n,d]=zp2tf(z,p,k)
n = (NOTE: numerator is this constant )
0 0 39476089
d =
1.0e+07 * (==========================)
0.0000 0.0009 3.9476 (<<< BEWARE MATLAB BUG!!!)
(==========================)
Warning: There is a bug in some versions
of matlab. Although it looks like
the values are zero, they are not.
The print format is suppressing
the significant digits. To fix this use:
>> h=tf(n,d,1/fsam)
Transfer function:
3.948e07
-----------------------
s^2 + 8886 s + 3.948e07
>> pzmap(h)
>> format short e
(Now you can see the denominator.)
>> d
d = (NOTE: denominator is s^2 + 8886 s + 3.9476e+07)
1.0000e+00 8.8855e+03 3.9476e+07
For high order filters, you need even greater
numerical precision to avoid errors. So,
in this case use:
>> format long e
>> d
d =
1.000000000000000e+00 8.885503812390156e+03
3.947608900000000e+07
>> bode(n,d) You should get this plot
>> print analogbut2.eps
>> [nd,dd]=impinvar(n,d,4000)
nd = (NOTE: numerator = 0 + 2622 z^-1 )
(NOTE: matlab does NOT premultiply by Ts)
1.0e+03 *
0 2.6220
dd = (NOTE: denominator = 1 -0.29 z^-1 + 0.10 z^-2)
1.0000 -0.2925 0.1085
Warning: In the 5.1.0.421 version
of matlab, they now seem to be
premultiplying by Ts, and so the
result now is:
>> [nd,dd]=impinvar(n,d,4000)
nd =
0 6.5549e-01 0
dd =
1.0000e+00 -2.9248e-01 1.0846e-01
>> freqz(nd,dd,256) You should get this plot
Warning: In the current (5.1.0.421) version
of matlab, they now seem to be
premultiplying by Ts, and so
You should get this plot
>> print impinvar.eps
>>[top,bot]=residuez(nd,dd) (NOTE: use this to find partial fraction form)
top =
0- 4.4428e+03i
0+ 4.4428e+03i
bot =
1.4624e-01+ 2.9508e-01i
1.4624e-01- 2.9508e-01i
(NOTE: Then partial frac expansion is:
-4442 i + 4442 i
--------------------- ---------------------
1- (.146 + .29 i)z^-1 1- (.146 - .29 i)z^-1
Warning: In the 5.1.0.421 version
of matlab, (will they ever get
their act together)
they now seem to be
premultiplying by Ts, and so the
result now is:
>> [top,bot]=residuez(nd,dd)
top =
0 + 1.1107e+00i
0 - 1.1107e+00i
bot =
1.4624e-01 - 2.9508e-01i
1.4624e-01 + 2.9508e-01i
>> [nd,dd]=bilinear(n,d,4000)
nd =
0.2261 0.4523 0.2261
dd =
1.0000 -0.2810 0.1856
>> freqz(nd,dd,256) You should get this plot
>> h=tf(nd,dd,1/fsam)
h =
0.2262 z^2 + 0.4523 z + 0.2262
------------------------------
z^2 - 0.2809 z + 0.1856
>> [nd,dd]=bilinear(n,d,4000,1000)
nd =
0.2929 0.5858 0.2929
dd =
1.0000 -0.0000 0.1716
>> freqz(nd,dd,256) You should get this plot
>> h=tf(nd,dd,1/fsam)
h =
0.2929 z^2 + 0.5858 z + 0.2929
------------------------------
z^2 - 2.22e-16 z + 0.1716
>> pzmap(h)
>> zgrid
>> zgrid([],[])
>> quit
127007 flops.
%
You may also want to try:<