From: <kk...@us...> - 2011-02-08 22:19:20
|
Revision: 114 http://python-control.svn.sourceforge.net/python-control/?rev=114&view=rev Author: kkchen Date: 2011-02-08 22:19:13 +0000 (Tue, 08 Feb 2011) Log Message: ----------- merging changes bb...@ra... Modified Paths: -------------- branches/control-0.4a/examples/pvtol-lqr.py branches/control-0.4a/examples/pvtol-nested-ss.py branches/control-0.4a/examples/pvtol-nested.py branches/control-0.4a/examples/secord-matlab.py branches/control-0.4a/examples/slicot-test.py branches/control-0.4a/src/freqplot.py branches/control-0.4a/src/pzmap.py branches/control-0.4a/src/statesp.py branches/control-0.4a/src/xferfcn.py Modified: branches/control-0.4a/examples/pvtol-lqr.py =================================================================== --- branches/control-0.4a/examples/pvtol-lqr.py 2011-02-08 22:19:07 UTC (rev 113) +++ branches/control-0.4a/examples/pvtol-lqr.py 2011-02-08 22:19:13 UTC (rev 114) @@ -119,7 +119,7 @@ subplot(221); title("Identity weights") # plot(T, Y[:,1, 1], '-', T, Y[:,2, 2], '--'); hold(True); -plot(Tx, Yx[0,:].T, '-', Ty, Yy[0,:].T, '--'); hold(True); +plot(Tx.T, Yx[0,:].T, '-', Ty.T, Yy[0,:].T, '--'); hold(True); plot([0, 10], [1, 1], 'k-'); hold(True); axis([0, 10, -0.1, 1.4]); @@ -141,9 +141,9 @@ [T3, Y3] = step(H1cx, T=linspace(0,10,100)); subplot(222); title("Effect of input weights") -plot(T1, Y1[0,:].T, 'b-'); hold(True); -plot(T2, Y2[0,:].T, 'b-'); hold(True); -plot(T3, Y3[0,:].T, 'b-'); hold(True); +plot(T1.T, Y1[0,:].T, 'b-'); hold(True); +plot(T2.T, Y2[0,:].T, 'b-'); hold(True); +plot(T3.T, Y3[0,:].T, 'b-'); hold(True); plot([0 ,10], [1, 1], 'k-'); hold(True); axis([0, 10, -0.1, 1.4]); @@ -162,7 +162,7 @@ subplot(223); title("Output weighting") [T2x, Y2x] = step(H2x, T=linspace(0,10,100)); [T2y, Y2y] = step(H2y, T=linspace(0,10,100)); -plot(T2x, Y2x[0,:].T, T2y, Y2y[0,:].T) +plot(T2x.T, Y2x[0,:].T, T2y.T, Y2y[0,:].T) ylabel('position'); xlabel('time'); ylabel('position'); legend(('x', 'y'), loc='lower right'); @@ -185,7 +185,7 @@ # step(H3x, H3y, 10); [T3x, Y3x] = step(H3x, T=linspace(0,10,100)); [T3y, Y3y] = step(H3y, T=linspace(0,10,100)); -plot(T3x, Y3x[0,:].T, T3y, Y3y[0,:].T) +plot(T3x.T, Y3x[0,:].T, T3y.T, Y3y[0,:].T) title("Physically motivated weights") xlabel('time'); legend(('x', 'y'), loc='lower right'); Modified: branches/control-0.4a/examples/pvtol-nested-ss.py =================================================================== --- branches/control-0.4a/examples/pvtol-nested-ss.py 2011-02-08 22:19:07 UTC (rev 113) +++ branches/control-0.4a/examples/pvtol-nested-ss.py 2011-02-08 22:19:13 UTC (rev 114) @@ -10,6 +10,7 @@ from matplotlib.pyplot import * # Grab MATLAB plotting functions from control.matlab import * # MATLAB-like functions +import numpy as np # System parameters m = 4; # mass of aircraft @@ -107,7 +108,7 @@ subplot(phaseh); semilogx([10^-4, 10^3], [-180, -180], 'k-') hold(True); -semilogx(w, phase, 'b-') +semilogx(w, np.squeeze(phase), 'b-') axis([10^-4, 10^3, -360, 0]); xlabel('Frequency [deg]'); ylabel('Phase [deg]'); # set(gca, 'YTick', [-360, -270, -180, -90, 0]); @@ -144,14 +145,15 @@ figure(9); (Tvec, Yvec) = step(T, None, linspace(1, 20)); -plot(Tvec, Yvec); hold(True); +plot(Tvec.T, Yvec.T); hold(True); (Tvec, Yvec) = step(Co*S, None, linspace(1, 20)); -plot(Tvec, Yvec); +plot(Tvec.T, Yvec.T); +#TODO: PZmap for statespace systems has not yet been implemented. figure(10); clf(); -(P, Z) = pzmap(T, Plot=True) -print "Closed loop poles and zeros: ", P, Z +#(P, Z) = pzmap(T, Plot=True) +#print "Closed loop poles and zeros: ", P, Z # Gang of Four figure(11); clf(); Modified: branches/control-0.4a/examples/pvtol-nested.py =================================================================== --- branches/control-0.4a/examples/pvtol-nested.py 2011-02-08 22:19:07 UTC (rev 113) +++ branches/control-0.4a/examples/pvtol-nested.py 2011-02-08 22:19:13 UTC (rev 114) @@ -10,6 +10,7 @@ from matplotlib.pyplot import * # Grab MATLAB plotting functions from control.matlab import * # MATLAB-like functions +import numpy as np # System parameters m = 4; # mass of aircraft @@ -23,8 +24,8 @@ Po = tf([1], [m, c, 0]); # outer loop (position) # Use state space versions -# Pi = tf2ss(Pi); -# Po = tf2ss(Po); +Pi = tf2ss(Pi); +Po = tf2ss(Po); # # Inner loop control design @@ -97,7 +98,7 @@ subplot(phaseh); semilogx([10^-4, 10^3], [-180, -180], 'k-') hold(True); -semilogx(w, phase, 'b-') +semilogx(w, np.squeeze(phase), 'b-') axis([10^-4, 10^3, -360, 0]); xlabel('Frequency [deg]'); ylabel('Phase [deg]'); # set(gca, 'YTick', [-360, -270, -180, -90, 0]); @@ -134,10 +135,10 @@ figure(9); (Tvec, Yvec) = step(T, None, linspace(1, 20)); -plot(Tvec, Yvec); hold(True); +plot(Tvec.T, Yvec.T); hold(True); (Tvec, Yvec) = step(Co*S, None, linspace(1, 20)); -plot(Tvec, Yvec); +plot(Tvec.T, Yvec.T); figure(10); clf(); (P, Z) = pzmap(T, Plot=True) Modified: branches/control-0.4a/examples/secord-matlab.py =================================================================== --- branches/control-0.4a/examples/secord-matlab.py 2011-02-08 22:19:07 UTC (rev 113) +++ branches/control-0.4a/examples/secord-matlab.py 2011-02-08 22:19:13 UTC (rev 114) @@ -18,7 +18,7 @@ # Step response for the system figure(1) T, yout = step(sys) -plot(T, yout) +plot(T.T, yout.T) # Bode plot for the system figure(2) Modified: branches/control-0.4a/examples/slicot-test.py =================================================================== --- branches/control-0.4a/examples/slicot-test.py 2011-02-08 22:19:07 UTC (rev 113) +++ branches/control-0.4a/examples/slicot-test.py 2011-02-08 22:19:13 UTC (rev 114) @@ -17,7 +17,7 @@ sys = ss(A, B, C, 0); # Eigenvalue placement -from slycot import sb01bd +#from slycot import sb01bd K = place(A, B, [-3, -2, -1]) print "Pole place: K = ", K print "Pole place: eigs = ", np.linalg.eig(A - B * K)[0] Modified: branches/control-0.4a/src/freqplot.py =================================================================== --- branches/control-0.4a/src/freqplot.py 2011-02-08 22:19:07 UTC (rev 113) +++ branches/control-0.4a/src/freqplot.py 2011-02-08 22:19:13 UTC (rev 114) @@ -177,7 +177,12 @@ # Select a default range if none is provided if (omega == None): omega = default_frequency_range(syslist) - + # Interpolate between wmin and wmax if a tuple or list are provided + elif (isinstance(omega,list) | isinstance(omega,tuple)): + # Only accept tuple or list of length 2 + if (len(omega) != 2): + raise ValueError("Supported frequency arguments are (wmin,wmax) tuple or list, or frequency vector. ") + omega = np.logspace(np.log10(omega[0]),np.log10(omega[1]),num=50,endpoint=True,base=10.0) for sys in syslist: if (sys.inputs > 1 or sys.outputs > 1): #TODO: Add MIMO nyquist plots. Modified: branches/control-0.4a/src/pzmap.py =================================================================== --- branches/control-0.4a/src/pzmap.py 2011-02-08 22:19:07 UTC (rev 113) +++ branches/control-0.4a/src/pzmap.py 2011-02-08 22:19:13 UTC (rev 114) @@ -42,6 +42,7 @@ import matplotlib.pyplot as plt import scipy as sp +import numpy as np import xferfcn # Compute poles and zeros for a system @@ -49,10 +50,10 @@ def pzmap(sys, Plot=True): """Plot a pole/zero map for a transfer function""" if (isinstance(sys, xferfcn.TransferFunction)): - poles = sp.roots(sys.den); - zeros = sp.roots(sys.num); + poles = sp.roots(np.squeeze(np.asarray(sys.den))); + zeros = sp.roots(np.squeeze(np.asarray(sys.num))); else: - raise TypeException + raise NotImplementedError("pzmap not implemented for state space systems yet.") if (Plot): # Plot the locations of the poles and zeros Modified: branches/control-0.4a/src/statesp.py =================================================================== --- branches/control-0.4a/src/statesp.py 2011-02-08 22:19:07 UTC (rev 113) +++ branches/control-0.4a/src/statesp.py 2011-02-08 22:19:13 UTC (rev 114) @@ -73,7 +73,7 @@ """ from numpy import all, angle, any, array, concatenate, cos, delete, dot, \ - empty, exp, eye, matrix, ones, pi, poly, poly1d, roots, sin, zeros + empty, exp, eye, matrix, ones, pi, poly, poly1d, roots, shape, sin, zeros from numpy.random import rand, randn from numpy.linalg import inv, det, solve from numpy.linalg.linalg import LinAlgError @@ -456,6 +456,8 @@ # TODO: transfer function to state space conversion is still buggy! print "Warning: transfer function to state space conversion by td04ad \ is still buggy!" + #print num + #print shape(num) ssout = td04ad(sys.inputs, sys.outputs, index, den, num) states = ssout[0] Modified: branches/control-0.4a/src/xferfcn.py =================================================================== --- branches/control-0.4a/src/xferfcn.py 2011-02-08 22:19:07 UTC (rev 113) +++ branches/control-0.4a/src/xferfcn.py 2011-02-08 22:19:13 UTC (rev 114) @@ -344,14 +344,14 @@ def __div__(self, other): """Divide two LTI objects.""" + # Convert the second argument to a transfer function. + other = _convertToTransferFunction(other) + if (self.inputs > 1 or self.outputs > 1 or other.inputs > 1 or other.outputs > 1): raise NotImplementedError("TransferFunction.__div__ is currently \ implemented only for SISO systems.") - # Convert the second argument to a transfer function. - other = _convertToTransferFunction(other) - num = polymul(self.num[0][0], other.den[0][0]) den = polymul(self.den[0][0], other.num[0][0]) @@ -487,8 +487,9 @@ computes the single denominator containing all the poles of sys.den, and reports it as the array d. The output numerator array n is modified to - use the common denominator. It is an sys.outputs-by-sys.inputs-by- - [something] array. + use the common denominator; the coefficient arrays are also padded with + zeros to be the same size as d. n is an sys.outputs-by-sys.inputs-by- + len(d) array. """ @@ -588,16 +589,12 @@ # Multiply in the missing poles. for p in missingpoles[i][j]: num[i][j] = polymul(num[i][j], [1., -p]) - # Find the largest numerator polynomial size. - largest = 0 + # Pad all numerator polynomials with zeros so that the numerator arrays + # are the same size as the denominator. for i in range(self.outputs): for j in range(self.inputs): - largest = max(largest, len(num[i][j])) - # Pad all smaller numerator polynomials with zeros. - for i in range(self.outputs): - for j in range(self.inputs): - num[i][j] = insert(num[i][j], zeros(largest - len(num[i][j])), - zeros(largest - len(num[i][j]))) + num[i][j] = insert(num[i][j], zeros(len(den) - len(num[i][j])), + zeros(len(den) - len(num[i][j]))) # Finally, convert the numerator to a 3-D array. num = array(num) # Remove trivial imaginary parts. Check for nontrivial imaginary parts. This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site. |