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. |