0
|
1 |
#!/usr/bin/python
|
|
2 |
# 7.1
|
|
3 |
import os, sys
|
|
4 |
sys.path += [os.getcwdu() + os.sep + ".." + os.sep + "python"]
|
|
5 |
|
|
6 |
import pylab as pl
|
|
7 |
from scipy import signal
|
|
8 |
|
|
9 |
def rlocus(H, k_vec=None):
|
|
10 |
"""Plots the Root Locus given
|
|
11 |
an lti system and
|
|
12 |
a vector of gain values"""
|
|
13 |
if k_vec is None:
|
|
14 |
k_vec = linspace(0, 100, 10000)
|
|
15 |
|
|
16 |
num = pl.poly1d(H.num)
|
|
17 |
den = pl.poly1d(H.den)
|
|
18 |
|
|
19 |
roots = []
|
|
20 |
for k in k_vec:
|
|
21 |
p = den + k*num
|
|
22 |
roots.append(pl.roots(p))
|
|
23 |
roots = pl.array(roots)
|
|
24 |
|
|
25 |
# Tried to do it with array slicing. Looks BAD
|
|
26 |
# n = H.num
|
|
27 |
# d = H.den
|
|
28 |
# l = d.shape[0]-n.shape[0]
|
|
29 |
# n_vec = n[pl.newaxis,:]*k_vec[:,pl.newaxis]
|
|
30 |
# A = pl.zeros((k_vec.shape[0], d.shape[0]))
|
|
31 |
# A = A+d
|
|
32 |
# A[:,l:] = n_vec
|
|
33 |
# r = map(pl.roots, A)
|
|
34 |
# r = pl.array(r)
|
|
35 |
|
|
36 |
x = pl.linspace(0, 2*pl.pi, 100)
|
|
37 |
pl.plot(pl.cos(x), pl.sin(x), 'b.', markersize=2)
|
|
38 |
|
|
39 |
for i, k_root in enumerate(roots[0]):
|
|
40 |
pl.plot(pl.real(roots[:,i]), pl.imag(roots[:,i]))
|
|
41 |
|
|
42 |
xmin, xmax = pl.xlim()
|
|
43 |
x = pl.linspace(int(xmin), int(xmax))
|
|
44 |
ymin, ymax = pl.ylim()
|
|
45 |
y = pl.linspace(int(ymin), int(ymax))
|
|
46 |
|
|
47 |
pl.plot(x, pl.zeros_like(x), 'b.', pl.zeros_like(y), y, 'b.', markersize=2)
|
|
48 |
|
|
49 |
pl.xlim(xmin, xmax)
|
|
50 |
pl.ylim(ymin, ymax)
|
|
51 |
|
|
52 |
pl.xlabel('Real part')
|
|
53 |
pl.ylabel('Imaginary part')
|
|
54 |
|
|
55 |
pl.show()
|
|
56 |
|
|
57 |
n = [1]
|
|
58 |
d = [1, -1, 0]
|
|
59 |
k_vec = pl.linspace(0, 20, 10000)
|
|
60 |
H = signal.lti(n, d)
|
|
61 |
|
|
62 |
rlocus(H, k_vec)
|