This script can be used to generate plot as in the Book of Modern Control Engineering 4th edition, International Edition, by Katsuhiko Ogata, Prentice Hall
Once when taking a subject Automatic Control is often to make a plot of LTI system, that generated in s(signal) domain by transforming mathematical model of the system using Laplace Transformation.
This example i took from page 307 of that book.
To plot transfer function =
(6.3223 s^2 + 18 s +12.811) / (s^4 + 6 s^3 + 11.3223 s^2 + 18s + 12.811)
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| from numpy import min, max | |
| from scipy import linspace | |
| from scipy.signal import lti, step, impulse | |
| # making transfer function | |
| # example from Ogata Modern Control Engineering | |
| # 4th edition, International Edition page 307 | |
| # num and den, can be list or numpy array type | |
| num = [6.3223, 18, 12.811] | |
| den = [1, 6, 11.3223, 18, 12.811] | |
| tf = lti(num, den) | |
| # get t = time, s = unit-step response | |
| t, s = step(tf) | |
| # recalculate t and s to get smooth plot | |
| t, s = step(tf, T = linspace(min(t), t[-1], 500)) | |
| # get i = impulse | |
| t, i = impulse(tf, T = linspace(min(t), t[-1], 500)) | |
| from matplotlib import pyplot as plt | |
| plt.plot(t, s, t, i) | |
| plt.title('Transient-Response Analysis') | |
| plt.xlabel('Time(sec)') | |
| plt.ylabel('Amplitude') | |
| plt.hlines(1, min(t), max(t), colors='r') | |
| plt.hlines(0, min(t), max(t)) | |
| plt.xlim(xmax=max(t)) | |
| plt.legend(('Unit-Step Response', 'Unit-Impulse Response'), loc=0) | |
| plt.grid() | |
| plt.show() |
I want to learn PID design. For that ur blog is a good start point. I am trying to run ur code. I m getting following error:
ReplyDeleteTraceback (most recent call last):
File "D:\Project (All)\Projects\MagLev\PID_Analysis_Python\Scipy_test.py", line 15, in
t, s = step(tf)
File "C:\Python26\lib\site-packages\scipy\signal\ltisys.py", line 716, in step
vals = lsim(sys, U, T, X0=X0)
File "C:\Python26\lib\site-packages\scipy\signal\ltisys.py", line 490, in lsim
GT = dot(dot(vti, diag(numpy.exp(dt * lam))), vt).astype(xout.dtype)
File "C:\Python26\lib\warnings.py", line 29, in _show_warning
file.write(formatwarning(message, category, filename, lineno, line))
TypeError: idle_formatwarning_subproc() takes exactly 4 arguments (5 given)
---
I am new to python .. but took care of the dependencies. Need help.
please refer to documentation your current version scipy.
DeleteK corrected it. Thanks a lot. :)
Delete.. And also if I am not jumping how to start with a matrix of transfer functions, say:
ReplyDeletenum1 = [-2.3833x10^3]
den1 = [1,160.3460,-1.9620x10^3,-3.14598852x10^5]
num2 = [31.9361,0,1.1132x10^5]
den2 = [1,160.3460,-1.9620x10^3,-3.14598852x10^5]
matrix0 = [0, lti(num1/den1); 1, lti(num2/den2)]
Any clue will be helpful.
because you put that on a list, not an array type.
Deletea = [1.923] # it's a list
b = numpy.array([0.839]) # it's an array
you can start by learning numpy.
After some trail n errors seeing the numpy manual I tried as follows:
Deletenum = np.array([[[0.0], [1.0]],[[-2383.3 ],[31.9361,0,111320.0 ]]])
den = np.array([[[0.0], [1.0]],[[1.0,160.3460,-1962.0,-314598.852],[1.0,160.3460,-1962.0,-314598.852]]])
tf = lti(num,den)
But its giving error as
Traceback (most recent call last):
File "D:\Project (All)\Projects\MagLev\PID_Analysis_Python\Scipy_test.py", line 36, in
tf = lti(num,den)
File "C:\Python26\lib\site-packages\scipy\signal\ltisys.py", line 236, in __init__
self.__dict__['num'], self.__dict__['den'] = normalize(*args)
File "C:\Python26\lib\site-packages\scipy\signal\filter_design.py", line 276, in normalize
raise ValueError("Denominator polynomial must be rank-1 array.")
ValueError: Denominator polynomial must be rank-1 array.
What can be done!
please take a look example from Ogata Modern Control Engineering
Delete4th edition, International Edition page 307, and my code above is the way to implement it using python.