add frequency_reponse() + FRD properties mag, phase, freq

murrayrm · Apr 16, 2022 · Mar 18, 2022 · Mar 23, 2022 · Mar 20, 2022 · Mar 25, 2022
commit ecdf1e7e6fa1516053b210300693570df6cbbfa3
diff --git a/control/frdata.py b/control/frdata.py
@@ -233,6 +233,29 @@ def __init__(self, *args, **kwargs):
             self.ifunc = None
         super().__init__(self.fresp.shape[1], self.fresp.shape[0])
 
+    #
+    # Frequency response properties
+    #
+    # Different properties of the frequency response that can be used for
+    # analysis and characterization.
+    #
+
+    @property
+    def magnitude(self):
+        return np.abs(self.fresp)
+
+    @property
+    def phase(self):
+        return np.angle(self.fresp)
+
+    @property
+    def frequency(self):
+        return self.omega
+
+    @property
+    def response(self):
+        return self.fresp
+
     def __str__(self):
         """String representation of the transfer function."""
 

diff --git a/control/lti.py b/control/lti.py
@@ -19,7 +19,7 @@
 
 __all__ = ['issiso', 'timebase', 'common_timebase', 'timebaseEqual',
            'isdtime', 'isctime', 'pole', 'zero', 'damp', 'evalfr',
-           'freqresp', 'dcgain']
+           'frequency_response', 'freqresp', 'dcgain']
 
 class LTI:
     """LTI is a parent class to linear time-invariant (LTI) system objects.
@@ -172,16 +172,16 @@ def frequency_response(self, omega, squeeze=None):
 
         Reports the frequency response of the system,
 
-             G(j*omega) = mag*exp(j*phase)
+             G(j*omega) = mag * exp(j*phase)
 
-        for continuous time systems. For discrete time systems, the response is
-        evaluated around the unit circle such that
+        for continuous time systems. For discrete time systems, the response
+        is evaluated around the unit circle such that
 
-             G(exp(j*omega*dt)) = mag*exp(j*phase).
+             G(exp(j*omega*dt)) = mag * exp(j*phase).
 
         In general the system may be multiple input, multiple output (MIMO),
-        where `m = self.ninputs` number of inputs and `p = self.noutputs` number
-        of outputs.
+        where `m = self.ninputs` number of inputs and `p = self.noutputs`
+        number of outputs.
 
         Parameters
         ----------
@@ -203,15 +203,15 @@ def frequency_response(self, omega, squeeze=None):
 
                 mag, phase, omega = response
 
-            where ``mag`` is the magnitude (absolute value, not dB or log10)
-            of the system frequency response, ``phase`` is the wrapped phase
-            in radians of the system frequency response, and ``omega`` is the
-            (sorted) frequencies at which the response was evaluated.  If the
-            system is SISO and squeeze is not True, ``mag`` and ``phase`` are
-            1D, indexed by frequency.  If the system is not SISO or squeeze is
-            False, the array is 3D, indexed by the output, input, and
-            frequency.  If ``squeeze`` is True then single-dimensional axes
-            are removed.
+            where ``mag`` is the magnitude (absolute value, not dB or
+            log10) of the system frequency response, ``phase`` is the wrapped
+            phase in radians of the system frequency response, and ``omega``
+            is the (sorted) frequencies at which the response was evaluated.
+            If the system is SISO and squeeze is not True, ``magnitude`` and
+            ``phase`` are 1D, indexed by frequency.  If the system is not SISO
+            or squeeze is False, the array is 3D, indexed by the output,
+            input, and frequency.  If ``squeeze`` is True then
+            single-dimensional axes are removed.
 
         """
         omega = np.sort(np.array(omega, ndmin=1))
@@ -597,7 +597,7 @@ def evalfr(sys, x, squeeze=None):
     """
     return sys.__call__(x, squeeze=squeeze)
 
-def freqresp(sys, omega, squeeze=None):
+def frequency_response(sys, omega, squeeze=None):
     """Frequency response of an LTI system at multiple angular frequencies.
 
     In general the system may be multiple input, multiple output (MIMO), where
@@ -671,6 +671,10 @@ def freqresp(sys, omega, squeeze=None):
     return sys.frequency_response(omega, squeeze=squeeze)
 
 
+# Alternative name (legacy)
+freqresp = frequency_response
+
+
 def dcgain(sys):
     """Return the zero-frequency (or DC) gain of the given system
 

diff --git a/control/tests/frd_test.py b/control/tests/frd_test.py
@@ -519,3 +519,16 @@ def test_to_pandas():
     # Check to make sure the data make senses
     np.testing.assert_equal(df['omega'], resp.omega)
     np.testing.assert_equal(df['H_{y[0], u[0]}'], resp.fresp[0, 0])
+
+
+def test_frequency_response():
+    # Create an SISO frequence response
+    sys = ct.rss(2, 2, 2)
+    omega = np.logspace(-2, 2, 20)
+    resp = ct.frequency_response(sys, omega)
+    eval = sys(omega*1j)
+
+    # Make sure we get the right answers in various ways
+    np.testing.assert_equal(resp.magnitude, np.abs(eval))
+    np.testing.assert_equal(resp.phase, np.angle(eval))
+    np.testing.assert_equal(resp.omega, omega)