|
From: <kk...@us...> - 2011-02-08 22:19:45
|
Revision: 119
http://python-control.svn.sourceforge.net/python-control/?rev=119&view=rev
Author: kkchen
Date: 2011-02-08 22:19:39 +0000 (Tue, 08 Feb 2011)
Log Message:
-----------
Added rss-balred.py to /examples and bode now returns mag, phase, freq, with the option to supress plotting.
Lauren Padilla <lpa...@pr...>
Modified Paths:
--------------
branches/control-0.4a/examples/secord-matlab.py
branches/control-0.4a/src/freqplot.py
branches/control-0.4a/src/modelsimp.py
branches/control-0.4a/src/xferfcn.py
Added Paths:
-----------
branches/control-0.4a/examples/rss-balred.py
Added: branches/control-0.4a/examples/rss-balred.py
===================================================================
--- branches/control-0.4a/examples/rss-balred.py (rev 0)
+++ branches/control-0.4a/examples/rss-balred.py 2011-02-08 22:19:39 UTC (rev 119)
@@ -0,0 +1,44 @@
+#!/usr/bin/env python
+
+import numpy as np
+import control.modelsimp as msimp
+import control.matlab as mt
+from control.statesp import StateSpace
+import matplotlib.pyplot as plt
+
+plt.close('all')
+
+#controlable canonical realization computed in matlab for the transfer function:
+# num = [1 11 45 32], den = [1 15 60 200 60]
+A = np.matrix('-15., -7.5, -6.25, -1.875; \
+8., 0., 0., 0.; \
+0., 4., 0., 0.; \
+0., 0., 1., 0.')
+B = np.matrix('2.; 0.; 0.; 0.')
+C = np.matrix('0.5, 0.6875, 0.7031, 0.5')
+D = np.matrix('0.')
+
+# The full system
+fsys = StateSpace(A,B,C,D)
+# The reduced system, truncating the order by 1
+ord = 3
+rsys = msimp.balred(fsys,ord, method = 'truncate')
+
+# Comparison of the step responses of the full and reduced systems
+plt.figure(1)
+t, y = mt.step(fsys)
+tr, yr = mt.step(rsys)
+plt.plot(t.T, y.T)
+plt.hold(True)
+plt.plot(tr.T, yr.T)
+
+# Repeat balanced reduction, now with 100-dimensional random state space
+sysrand = mt.rss(100, 1, 1)
+rsysrand = msimp.balred(sysrand,10,method ='truncate')
+
+# Comparison of the impulse responses of the full and reduced random systems
+plt.figure(2)
+trand, yrand = mt.impulse(sysrand)
+trandr, yrandr = mt.impulse(rsysrand)
+plt.plot(trand.T, yrand.T,trandr.T, yrandr.T)
+
Property changes on: branches/control-0.4a/examples/rss-balred.py
___________________________________________________________________
Added: svn:executable
+ *
Modified: branches/control-0.4a/examples/secord-matlab.py
===================================================================
--- branches/control-0.4a/examples/secord-matlab.py 2011-02-08 22:19:33 UTC (rev 118)
+++ branches/control-0.4a/examples/secord-matlab.py 2011-02-08 22:19:39 UTC (rev 119)
@@ -22,7 +22,7 @@
# Bode plot for the system
figure(2)
-bode(sys, logspace(-2, 2))
+mag,phase,om = bode(sys, logspace(-2, 2),Plot=True)
# Nyquist plot for the system
figure(3)
Modified: branches/control-0.4a/src/freqplot.py
===================================================================
--- branches/control-0.4a/src/freqplot.py 2011-02-08 22:19:33 UTC (rev 118)
+++ branches/control-0.4a/src/freqplot.py 2011-02-08 22:19:39 UTC (rev 119)
@@ -54,12 +54,12 @@
#
# Bode plot
-def bode(syslist, omega=None, dB=False, Hz=False, color=None):
+def bode(syslist, omega=None, dB=False, Hz=False, color=None, Plot=True):
"""Bode plot for a system
Usage
=====
- (magh, phaseh) = bode(syslist, omega=None, dB=False, Hz=False)
+ (magh, phaseh, omega) = bode(syslist, omega=None, dB=False, Hz=False, color=None, Plot=True)
Plots a Bode plot for the system over a (optional) frequency range.
@@ -73,16 +73,18 @@
If True, plot result in dB
Hz : boolean
If True, plot frequency in Hz (omega must be provided in rad/sec)
+ Plot : boolean
+ If True, plot magnitude and phase
Return values
-------------
- magh : graphics handle to magnitude plot (for rescaling, etc)
- phaseh : graphics handle to phase plot
+ magh : magnitude array
+ phaseh : phase array
+ omega : frequency array
Notes
-----
- 1. Use (mag, phase, freq) = sys.freqresp(freq) to generate the
- frequency response for a system.
+ 1. Alternatively, may use (mag, phase, freq) = sys.freqresp(freq) to generate the frequency response for a system.
"""
# If argument was a singleton, turn it into a list
if (not getattr(syslist, '__iter__', False)):
@@ -108,46 +110,47 @@
# Get the dimensions of the current axis, which we will divide up
#! TODO: Not current implemented; just use subplot for now
- # Magnitude plot
- plt.subplot(211);
- if dB:
- if color==None:
- plt.semilogx(omega, mag)
+ if (Plot):
+ # Magnitude plot
+ plt.subplot(211);
+ if dB:
+ if color==None:
+ plt.semilogx(omega, mag)
+ else:
+ plt.semilogx(omega, mag, color=color)
+ plt.ylabel("Magnitude (dB)")
else:
- plt.semilogx(omega, mag, color=color)
- plt.ylabel("Magnitude (dB)")
- else:
- if color==None:
- plt.loglog(omega, mag)
- else:
- plt.loglog(omega, mag, color=color)
- plt.ylabel("Magnitude")
+ if color==None:
+ plt.loglog(omega, mag)
+ else:
+ plt.loglog(omega, mag, color=color)
+ plt.ylabel("Magnitude")
- # Add a grid to the plot
- plt.grid(True)
- plt.grid(True, which='minor')
- plt.hold(True);
+ # Add a grid to the plot
+ plt.grid(True)
+ plt.grid(True, which='minor')
+ plt.hold(True);
- # Phase plot
- plt.subplot(212);
- if color==None:
- plt.semilogx(omega, phase)
- else:
- plt.semilogx(omega, phase, color=color)
- plt.hold(True)
+ # Phase plot
+ plt.subplot(212);
+ if color==None:
+ plt.semilogx(omega, phase)
+ else:
+ plt.semilogx(omega, phase, color=color)
+ plt.hold(True)
- # Add a grid to the plot
- plt.grid(True)
- plt.grid(True, which='minor')
- plt.ylabel("Phase (deg)")
+ # Add a grid to the plot
+ plt.grid(True)
+ plt.grid(True, which='minor')
+ plt.ylabel("Phase (deg)")
- # Label the frequency axis
- if Hz:
- plt.xlabel("Frequency (Hz)")
- else:
- plt.xlabel("Frequency (rad/sec)")
+ # Label the frequency axis
+ if Hz:
+ plt.xlabel("Frequency (Hz)")
+ else:
+ plt.xlabel("Frequency (rad/sec)")
- return (211, 212)
+ return mag, phase, omega
# Nyquist plot
def nyquist(syslist, omega=None):
Modified: branches/control-0.4a/src/modelsimp.py
===================================================================
--- branches/control-0.4a/src/modelsimp.py 2011-02-08 22:19:33 UTC (rev 118)
+++ branches/control-0.4a/src/modelsimp.py 2011-02-08 22:19:39 UTC (rev 119)
@@ -207,6 +207,8 @@
#Check system is stable
D,V = np.linalg.eig(sys.A)
+ print D.shape
+ print D
for e in D:
if e.real >= 0:
raise ValueError, "Oops, the system is unstable!"
Modified: branches/control-0.4a/src/xferfcn.py
===================================================================
--- branches/control-0.4a/src/xferfcn.py 2011-02-08 22:19:33 UTC (rev 118)
+++ branches/control-0.4a/src/xferfcn.py 2011-02-08 22:19:39 UTC (rev 119)
@@ -710,7 +710,8 @@
num[i][j] = list(tfout[6][i, j, :])
# Each transfer function matrix row has a common denominator.
den[i][j] = list(tfout[5][i, :])
-
+ print num
+ print den
return TransferFunction(num, den)
elif isinstance(sys, (int, long, float, complex)):
if "inputs" in kw:
This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
|