added method horner, for evaluating s-response, and added a feedback on

repagh · repagh · commit 38f6c8043ee8 · 2013-05-24T11:45:18.000Z
selected point, gain, damping in rlocus plot
diff --git a/src/bdalg.py b/src/bdalg.py
@@ -244,8 +244,9 @@ def append(*sys):
     Group models by appending their inputs and outputs
 
     Forms an augmented system model, and appends the inputs and
-    outputs together. The system type will be the type of the first system
-    given.
+    outputs together. The system type will be the type of the first
+    system given; if you mix state-space systems and gain matrices,
+    make sure the gain matrices are not first.
 
     Parameters.
     -----------
diff --git a/src/matlab.py b/src/matlab.py
@@ -1204,13 +1204,13 @@ def damp(sys, doprint=True):
     '''
     wn, damping, poles = sys.damp()
     if doprint:
-        print('______Eigenvalue______ Damping___ Frequency_')
+        print('_____Eigenvalue______ Damping___ Frequency_')
         for p, d, w in zip(poles, damping, wn) :
             if abs(p.imag) < 1e-12:
-                print("%10.4g             %10.4g %10.4g" % 
+                print("%10.4g            %10.4g %10.4g" % 
                       (p.real, 1.0, -p.real)) 
             else:
-                print("%10.4g %10.4gi %10.4g %10.4g" % 
+                print("%10.4g%+10.4gj %10.4g %10.4g" % 
                       (p.real, p.imag, d, w)) 
     return wn, damping, poles
 
diff --git a/src/rlocus.py b/src/rlocus.py
@@ -50,9 +50,11 @@
 import scipy.signal             # signal processing toolbox
 import pylab                    # plotting routines
 import control.xferfcn as xferfcn
+from functools import partial
 
 # Main function: compute a root locus diagram
-def root_locus(sys, kvect, xlim=None, ylim=None, plotstr='-', Plot=True):
+def root_locus(sys, kvect, xlim=None, ylim=None, plotstr='-', Plot=True, 
+               PrintGain=True):
     """Calculate the root locus by finding the roots of 1+k*TF(s)
     where TF is self.num(s)/self.den(s) and each k is an element
     of kvect.
@@ -65,7 +67,9 @@ def root_locus(sys, kvect, xlim=None, ylim=None, plotstr='-', Plot=True):
         List of gains to use in computing diagram
     Plot : boolean (default = True)
         If True, plot magnitude and phase
-
+    PrintGain: boolean (default = True)
+        If True, report mouse clicks when close to the root-locus branches,
+        calculate gain, damping and print
     Return values
     -------------
     rlist : list of computed root locations
@@ -80,6 +84,10 @@ def root_locus(sys, kvect, xlim=None, ylim=None, plotstr='-', Plot=True):
 
     # Create the plot
     if (Plot):
+        f = pylab.figure()
+        if PrintGain:
+            cid = f.canvas.mpl_connect(
+                'button_release_event', partial(_RLFeedbackClicks, sys=sys))
         ax = pylab.axes();
 
         # plot open loop poles
@@ -165,3 +173,12 @@ def _RLSortRoots(sys, mymat):
                 sorted[n,ind] = elem
         prevrow = sorted[n,:]
     return sorted
+
+def _RLFeedbackClicks(event, sys):
+    """Print root-locus gain feedback for clicks on the root-locus plot
+    """
+    s = complex(event.xdata, event.ydata)
+    K = -1./sys.horner(s)
+    if abs(K.real) > 1e-8 and abs(K.imag/K.real) < 0.04:
+        print("Clicked at %10.4g%+10.4gj gain %10.4g damp %10.4g" %
+              (s.real, s.imag, K.real, -1*s.real/abs(s)))
diff --git a/src/statesp.py b/src/statesp.py
@@ -384,10 +384,16 @@ def evalfr(self, omega):
         else:
             s = omega * 1.j
 
-        fresp = self.C * solve(s * eye(self.states) - self.A,
-            self.B) + self.D
+        return self.horner(s)
 
-        return array(fresp)
+    def horner(self, s):
+        '''Evaluate the systems's transfer function for a complex variable
+        
+        Returns a matrix of values evaluated at complex variable s.
+        '''
+        resp = self.C * solve(s * eye(self.states) - self.A,
+                              self.B) + self.D
+        return array(resp)
 
     # Method for generating the frequency response of the system
     # TODO: add discrete time check
diff --git a/src/xferfcn.py b/src/xferfcn.py
@@ -527,6 +527,14 @@ def evalfr(self, omega):
                 warn("evalfr: frequency evaluation above Nyquist frequency")
         else:
             s = 1.j * omega
+            
+        return self.horner(s)
+
+    def horner(self, s):
+        '''Evaluate the systems's transfer function for a complex variable
+        
+        Returns a matrix of values evaluated at complex variable s.
+        '''
 
         # Preallocate the output.
         out = empty((self.outputs, self.inputs), dtype=complex)
