Repository: python-control/python-control
Pull request #900 (2023-06-02 18:03:23)

warn if prewarp-frequency is not used
This PR makes sure the prewarp-frequency kwarg is used in `sys.sample`. 
* now warns if a descretization type doesn't support it
* now c2d is an identical copy of sample_system instead of being a somewhat-compatible version of it
* more/better unit tests. 


Patch:
diff --git a/control/dtime.py b/control/dtime.py
index 6197ae8a..38fcf805 100644
--- a/control/dtime.py
+++ b/control/dtime.py
@@ -72,22 +72,12 @@ def sample_system(sysc, Ts, method='zoh', alpha=None, prewarp_frequency=None,
             otherwise. See :func:`scipy.signal.cont2discrete`.
     prewarp_frequency : float within [0, infinity)
         The frequency [rad/s] at which to match with the input continuous-
-        time system's magnitude and phase (only valid for method='bilinear')
-    name : string, optional
-        Set the name of the sampled system.  If not specified and
-        if `copy_names` is `False`, a generic name <sys[id]> is generated
-        with a unique integer id.  If `copy_names` is `True`, the new system
-        name is determined by adding the prefix and suffix strings in
-        config.defaults['namedio.sampled_system_name_prefix'] and
-        config.defaults['namedio.sampled_system_name_suffix'], with the
-        default being to add the suffix '$sampled'.
-    copy_names : bool, Optional
-        If True, copy the names of the input signals, output
-        signals, and states to the sampled system.
+        time system's magnitude and phase (only valid for method='bilinear',
+        'tustin', or 'gbt' with alpha=0.5)
 
     Returns
     -------
-    sysd : linsys
+    sysd : LTI of the same class (:class:`StateSpace` or :class:`TransferFunction`)
         Discrete time system, with sampling rate Ts
 
     Other Parameters
@@ -101,6 +91,17 @@ def sample_system(sysc, Ts, method='zoh', alpha=None, prewarp_frequency=None,
     states : int, list of str, or None, optional
         Description of the system states.  Same format as `inputs`. Only
         available if the system is :class:`StateSpace`.
+    name : string, optional
+        Set the name of the sampled system.  If not specified and
+        if `copy_names` is `False`, a generic name <sys[id]> is generated
+        with a unique integer id.  If `copy_names` is `True`, the new system
+        name is determined by adding the prefix and suffix strings in
+        config.defaults['namedio.sampled_system_name_prefix'] and
+        config.defaults['namedio.sampled_system_name_suffix'], with the
+        default being to add the suffix '$sampled'.
+    copy_names : bool, Optional
+        If True, copy the names of the input signals, output
+        signals, and states to the sampled system.
 
     Notes
     -----
@@ -126,46 +127,4 @@ def sample_system(sysc, Ts, method='zoh', alpha=None, prewarp_frequency=None,
         method=method, alpha=alpha, prewarp_frequency=prewarp_frequency,
         name=name, copy_names=copy_names, **kwargs)
 
-
-def c2d(sysc, Ts, method='zoh', prewarp_frequency=None):
-    """
-    Convert a continuous time system to discrete time by sampling
-
-    Parameters
-    ----------
-    sysc : LTI (:class:`StateSpace` or :class:`TransferFunction`)
-        Continuous time system to be converted
-    Ts : float > 0
-        Sampling period
-    method : string
-        Method to use for conversion, e.g. 'bilinear', 'zoh' (default)
-    prewarp_frequency : real within [0, infinity)
-        The frequency [rad/s] at which to match with the input continuous-
-        time system's magnitude and phase (only valid for method='bilinear')
-
-    Returns
-    -------
-    sysd : LTI of the same class
-        Discrete time system, with sampling rate Ts
-
-    Notes
-    -----
-    See :meth:`StateSpace.sample` or :meth:`TransferFunction.sample`` for
-    further details.
-
-    Examples
-    --------
-    >>> Gc = ct.tf([1], [1, 2, 1])
-    >>> Gc.isdtime()
-    False
-    >>> Gd = ct.sample_system(Gc, 1, method='bilinear')
-    >>> Gd.isdtime()
-    True
-
-    """
-
-    #  Call the sample_system() function to do the work
-    sysd = sample_system(sysc, Ts,
-        method=method, prewarp_frequency=prewarp_frequency)
-
-    return sysd
+c2d = sample_system
\ No newline at end of file
diff --git a/control/iosys.py b/control/iosys.py
index dca00d3e..78444f7c 100644
--- a/control/iosys.py
+++ b/control/iosys.py
@@ -1862,7 +1862,7 @@ def input_output_response(
 
         return TimeResponseData(
             t_eval, y, None, u, issiso=sys.issiso(),
-            output_labels=sys.output_labels, input_labels=sys.input_labels,
+            output_labels=sys.output_index, input_labels=sys.input_index,
             transpose=transpose, return_x=return_x, squeeze=squeeze)
 
     # Create a lambda function for the right hand side
@@ -1941,8 +1941,8 @@ def input_output_response(
 
     return TimeResponseData(
         soln.t, y, soln.y, u, issiso=sys.issiso(),
-        output_labels=sys.output_labels, input_labels=sys.input_labels,
-        state_labels=sys.state_labels,
+        output_labels=sys.output_index, input_labels=sys.input_index,
+        state_labels=sys.state_index,
         transpose=transpose, return_x=return_x, squeeze=squeeze)
 
 
@@ -2881,7 +2881,7 @@ def interconnect(
 
             # Look for the signal name as a system input
             for sys in syslist:
-                if signal_name in sys.input_labels:
+                if signal_name in sys.input_index.keys():
                     connection.append(sign + sys.name + "." + signal_name)
 
             # Make sure we found the name
diff --git a/control/lti.py b/control/lti.py
index c904c150..216332e9 100644
--- a/control/lti.py
+++ b/control/lti.py
@@ -2,14 +2,22 @@
 
 The lti module contains the LTI parent class to the child classes StateSpace
 and TransferFunction.  It is designed for use in the python-control library.
+
+Routines in this module:
+
+LTI.__init__
+isdtime()
+isctime()
+timebase()
+common_timebase()
 """
 
 import numpy as np
 
-from numpy import real, angle, abs
+from numpy import absolute, real, angle, abs
 from warnings import warn
 from . import config
-from .namedio import NamedIOSystem
+from .namedio import NamedIOSystem, isdtime
 
 __all__ = ['poles', 'zeros', 'damp', 'evalfr', 'frequency_response',
            'freqresp', 'dcgain', 'bandwidth', 'pole', 'zero']
@@ -102,11 +110,11 @@ class LTI(NamedIOSystem):
         Returns
         -------
         wn : array
-            Natural frequency for each system pole
+            Natural frequencies for each system pole
         zeta : array
             Damping ratio for each system pole
         poles : array
-            System pole locations
+            Array of system poles
         '''
         poles = self.poles()
 
@@ -114,9 +122,9 @@ class LTI(NamedIOSystem):
             splane_poles = np.log(poles.astype(complex))/self.dt
         else:
             splane_poles = poles
-        wn = abs(splane_poles)
-        zeta = -real(splane_poles)/wn
-        return wn, zeta, poles
+        wn = absolute(splane_poles)
+        Z = -real(splane_poles)/wn
+        return wn, Z, poles
 
     def frequency_response(self, omega, squeeze=None):
         """Evaluate the linear time-invariant system at an array of angular
@@ -338,23 +346,25 @@ def zero(sys):
 
 def damp(sys, doprint=True):
     """
-    Compute natural frequencies, damping ratios, and poles of a system
+    Compute natural frequency, damping ratio, and poles of a system
+
+    The function takes 1 or 2 parameters
 
     Parameters
     ----------
-    sys : LTI (StateSpace or TransferFunction)
+    sys: LTI (StateSpace or TransferFunction)
         A linear system object
-    doprint : bool (optional)
-        if True, print table with values
+    doprint:
+        if true, print table with values
 
     Returns
     -------
-    wn : array
-        Natural frequency for each system pole
-    zeta : array
-        Damping ratio for each system pole
-    poles : array
-        System pole locations
+    wn: array
+        Natural frequencies of the poles
+    damping: array
+        Damping values
+    poles: array
+        Pole locations
 
     See Also
     --------
@@ -364,37 +374,37 @@ def damp(sys, doprint=True):
     -----
     If the system is continuous,
         wn = abs(poles)
-        zeta  = -real(poles)/poles
+        Z  = -real(poles)/poles.
 
     If the system is discrete, the discrete poles are mapped to their
     equivalent location in the s-plane via
 
-        s = log(poles)/dt
+        s = log10(poles)/dt
 
     and
 
         wn = abs(s)
-        zeta = -real(s)/wn.
+        Z = -real(s)/wn.
 
     Examples
     --------
     >>> G = ct.tf([1], [1, 4])
-    >>> wn, zeta, poles = ct.damp(G)
-        Eigenvalue (pole)       Damping     Frequency
-                       -4             1             4
+    >>> wn, damping, poles = ct.damp(G)
+    _____Eigenvalue______ Damping___ Frequency_
+            -4                     1          4
 
     """
-    wn, zeta, poles = sys.damp()
+    wn, damping, poles = sys.damp()
     if doprint:
-        print('    Eigenvalue (pole)       Damping     Frequency')
-        for p, z, w in zip(poles, zeta, wn):
+        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" %
-                      (p.real, 1.0, w))
+                print("%10.4g            %10.4g %10.4g" %
+                      (p.real, 1.0, -p.real))
             else:
-                print("%10.4g%+10.4gj    %10.4g    %10.4g" %
-                      (p.real, p.imag, z, w))
-    return wn, zeta, poles
+                print("%10.4g%+10.4gj %10.4g %10.4g" %
+                      (p.real, p.imag, d, w))
+    return wn, damping, poles
 
 
 def evalfr(sys, x, squeeze=None):
diff --git a/control/statesp.py b/control/statesp.py
index 8661d874..4f80bdb2 100644
--- a/control/statesp.py
+++ b/control/statesp.py
@@ -170,9 +170,9 @@ class StateSpace(LTI):
 
     The StateSpace class is used to represent state-space realizations of
     linear time-invariant (LTI) systems:
-    
+
     .. math::
-    
+
           dx/dt &= A x + B u \\
               y &= C x + D u
 
@@ -1368,10 +1368,13 @@ class StateSpace(LTI):
         """
         if not self.isctime():
             raise ValueError("System must be continuous time system")
-
-        if (method == 'bilinear' or (method == 'gbt' and alpha == 0.5)) and \
-                prewarp_frequency is not None:
-            Twarp = 2 * np.tan(prewarp_frequency * Ts/2)/prewarp_frequency
+        if prewarp_frequency is not None:
+            if method in ('bilinear', 'tustin') or \
+                    (method == 'gbt' and alpha == 0.5):
+                Twarp = 2*np.tan(prewarp_frequency*Ts/2)/prewarp_frequency
+            else:
+                warn('prewarp_frequency ignored: incompatible conversion')
+                Twarp = Ts
         else:
             Twarp = Ts
         sys = (self.A, self.B, self.C, self.D)
@@ -1777,9 +1780,7 @@ def _mimo2siso(sys, input, output, warn_conversion=False):
         new_B = sys.B[:, input]
         new_C = sys.C[output, :]
         new_D = sys.D[output, input]
-        sys = StateSpace(sys.A, new_B, new_C, new_D, sys.dt,
-                         name=sys.name,
-                         inputs=sys.input_labels[input], outputs=sys.output_labels[output])
+        sys = StateSpace(sys.A, new_B, new_C, new_D, sys.dt)
 
     return sys
 
@@ -1828,9 +1829,7 @@ def _mimo2simo(sys, input, warn_conversion=False):
         #  Y = C*X + D*U
         new_B = sys.B[:, input:input+1]
         new_D = sys.D[:, input:input+1]
-        sys = StateSpace(sys.A, new_B, sys.C, new_D, sys.dt,
-                         name=sys.name,
-                         inputs=sys.input_labels[input], outputs=sys.output_labels)
+        sys = StateSpace(sys.A, new_B, sys.C, new_D, sys.dt)
 
     return sys
 
diff --git a/control/timeresp.py b/control/timeresp.py
index bd8595cf..638a0732 100644
--- a/control/timeresp.py
+++ b/control/timeresp.py
@@ -694,10 +694,7 @@ def _process_labels(labels, signal, length):
             raise ValueError("Name dictionary for %s is incomplete" % signal)
 
     # Convert labels to a list
-    if isinstance(labels, str):
-        labels = [labels]
-    else:
-        labels = list(labels)
+    labels = list(labels)
 
     # Make sure the signal list is the right length and type
     if len(labels) != length:
@@ -1114,8 +1111,6 @@ def forced_response(sys, T=None, U=0., X0=0., transpose=False,
 
     return TimeResponseData(
         tout, yout, xout, U, issiso=sys.issiso(),
-        output_labels=sys.output_labels, input_labels=sys.input_labels,
-        state_labels=sys.state_labels,
         transpose=transpose, return_x=return_x, squeeze=squeeze)
 
 
@@ -1379,16 +1374,8 @@ def step_response(sys, T=None, X0=0., input=None, output=None, T_num=None,
     # Figure out if the system is SISO or not
     issiso = sys.issiso() or (input is not None and output is not None)
 
-    # Select only the given input and output, if any
-    input_labels = sys.input_labels if input is None \
-        else sys.input_labels[input]
-    output_labels = sys.output_labels if output is None \
-        else sys.output_labels[output]
-
     return TimeResponseData(
         response.time, yout, xout, uout, issiso=issiso,
-        output_labels=output_labels, input_labels=input_labels,
-        state_labels=sys.state_labels,
         transpose=transpose, return_x=return_x, squeeze=squeeze)
 
 
@@ -1717,15 +1704,9 @@ def initial_response(sys, T=None, X0=0., input=0, output=None, T_num=None,
     # Figure out if the system is SISO or not
     issiso = sys.issiso() or (input is not None and output is not None)
 
-    # Select only the given output, if any
-    output_labels = sys.output_labels if output is None \
-        else sys.output_labels[0]
-
     # Store the response without an input
     return TimeResponseData(
         response.t, response.y, response.x, None, issiso=issiso,
-        output_labels=output_labels, input_labels=None,
-        state_labels=sys.state_labels,
         transpose=transpose, return_x=return_x, squeeze=squeeze)
 
 
@@ -1817,7 +1798,7 @@ def impulse_response(sys, T=None, X0=0., input=None, output=None, T_num=None,
     -----
     This function uses the `forced_response` function to compute the time
     response. For continuous time systems, the initial condition is altered to
-    account for the initial impulse. For discrete-time aystems, the impulse is
+    account for the initial impulse. For discrete-time aystems, the impulse is 
     sized so that it has unit area.
 
     Examples
@@ -1888,16 +1869,8 @@ def impulse_response(sys, T=None, X0=0., input=None, output=None, T_num=None,
     # Figure out if the system is SISO or not
     issiso = sys.issiso() or (input is not None and output is not None)
 
-    # Select only the given input and output, if any
-    input_labels = sys.input_labels if input is None \
-        else sys.input_labels[input]
-    output_labels = sys.output_labels if output is None \
-        else sys.output_labels[output]
-
     return TimeResponseData(
         response.time, yout, xout, uout, issiso=issiso,
-        output_labels=output_labels, input_labels=input_labels,
-        state_labels=sys.state_labels,
         transpose=transpose, return_x=return_x, squeeze=squeeze)
 
 
diff --git a/control/xferfcn.py b/control/xferfcn.py
index a6a00c5d..8ef7e108 100644
--- a/control/xferfcn.py
+++ b/control/xferfcn.py
@@ -1134,7 +1134,7 @@ class TransferFunction(LTI):
             Method to use for sampling:
 
             * gbt: generalized bilinear transformation
-            * bilinear: Tustin's approximation ("gbt" with alpha=0.5)
+            * bilinear or tustin: Tustin's approximation ("gbt" with alpha=0.5)
             * euler: Euler (or forward difference) method ("gbt" with alpha=0)
             * backward_diff: Backwards difference ("gbt" with alpha=1.0)
             * zoh: zero-order hold (default)
@@ -1192,9 +1192,13 @@ class TransferFunction(LTI):
         if method == "matched":
             return _c2d_matched(self, Ts)
         sys = (self.num[0][0], self.den[0][0])
-        if (method == 'bilinear' or (method == 'gbt' and alpha == 0.5)) and \
-                prewarp_frequency is not None:
-            Twarp = 2*np.tan(prewarp_frequency*Ts/2)/prewarp_frequency
+        if prewarp_frequency is not None:
+            if method in ('bilinear', 'tustin') or \
+                    (method == 'gbt' and alpha == 0.5):
+                Twarp = 2*np.tan(prewarp_frequency*Ts/2)/prewarp_frequency
+            else:
+                warn('prewarp_frequency ignored: incompatible conversion')
+                Twarp = Ts
         else:
             Twarp = Ts
         numd, dend, _ = cont2discrete(sys, Twarp, method, alpha)
diff --git a/doc/examples.rst b/doc/examples.rst
index 505bcf7a..02330f90 100644
--- a/doc/examples.rst
+++ b/doc/examples.rst
@@ -27,8 +27,6 @@ other sources.
    phaseplots
    robust_siso
    robust_mimo
-   scherer_etal_ex7_H2_h2syn
-   scherer_etal_ex7_Hinf_hinfsyn
    cruise-control
    steering-gainsched
    steering-optimal
diff --git a/doc/scherer_etal_ex7_H2_h2syn.py b/doc/scherer_etal_ex7_H2_h2syn.py
deleted file mode 120000
index 527f8014..00000000
--- a/doc/scherer_etal_ex7_H2_h2syn.py
+++ /dev/null
@@ -1,1 +0,0 @@
-../examples/scherer_etal_ex7_H2_h2syn.py
\ No newline at end of file
diff --git a/doc/scherer_etal_ex7_H2_h2syn.rst b/doc/scherer_etal_ex7_H2_h2syn.rst
deleted file mode 100644
index ef386be6..00000000
--- a/doc/scherer_etal_ex7_H2_h2syn.rst
+++ /dev/null
@@ -1,15 +0,0 @@
-H2 synthesis, based on Scherer et al. 1997 example 7
-----------------------------------------------------
-
-Code
-....
-.. literalinclude:: scherer_etal_ex7_H2_h2syn.py
-   :language: python
-   :linenos:
-
-
-Notes
-.....
-
-1. The environment variable `PYCONTROL_TEST_EXAMPLES` is used for
-testing to turn off plotting of the outputs.
diff --git a/doc/scherer_etal_ex7_Hinf_hinfsyn.py b/doc/scherer_etal_ex7_Hinf_hinfsyn.py
deleted file mode 120000
index 7755a325..00000000
--- a/doc/scherer_etal_ex7_Hinf_hinfsyn.py
+++ /dev/null
@@ -1,1 +0,0 @@
-../examples/scherer_etal_ex7_Hinf_hinfsyn.py
\ No newline at end of file
diff --git a/doc/scherer_etal_ex7_Hinf_hinfsyn.rst b/doc/scherer_etal_ex7_Hinf_hinfsyn.rst
deleted file mode 100644
index 1a229453..00000000
--- a/doc/scherer_etal_ex7_Hinf_hinfsyn.rst
+++ /dev/null
@@ -1,15 +0,0 @@
-Hinf synthesis, based on Scherer et al. 1997 example 7
-------------------------------------------------------
-
-Code
-....
-.. literalinclude:: scherer_etal_ex7_Hinf_hinfsyn.py
-   :language: python
-   :linenos:
-
-
-Notes
-.....
-
-1. The environment variable `PYCONTROL_TEST_EXAMPLES` is used for
-testing to turn off plotting of the outputs.
diff --git a/examples/scherer_etal_ex7_H2_h2syn.py b/examples/scherer_etal_ex7_H2_h2syn.py
deleted file mode 100644
index 0b47dd8a..00000000
--- a/examples/scherer_etal_ex7_H2_h2syn.py
+++ /dev/null
@@ -1,84 +0,0 @@
-"""H2 design using h2syn.
-
-Demonstrate H2 desing for a SISO plant using h2syn. Based on [1], Ex. 7.
-
-[1] Scherer, Gahinet, & Chilali, "Multiobjective Output-Feedback Control via
-LMI Optimization", IEEE Trans. Automatic Control, Vol. 42, No. 7, July 1997.
-
-[2] Zhou & Doyle, "Essentials of Robust Control", Prentice Hall,
-Upper Saddle River, NJ, 1998.
-"""
-# %%
-# Packages
-import numpy as np
-import control
-
-# %%
-# State-space system.
-
-# Process model.
-A = np.array([[0, 10, 2],
-              [-1, 1, 0],
-              [0, 2, -5]])
-B1 = np.array([[1],
-               [0],
-               [1]])
-B2 = np.array([[0],
-               [1],
-               [0]])
-
-# Plant output.
-C2 = np.array([[0, 1, 0]])
-D21 = np.array([[2]])
-D22 = np.array([[0]])
-
-# H2 performance.
-C1 = np.array([[0, 1, 0],
-               [0, 0, 1],
-               [0, 0, 0]])
-D11 = np.array([[0],
-                [0],
-                [0]])
-D12 = np.array([[0],
-                [0],
-                [1]])
-
-# Dimensions.
-n_u, n_y = 1, 1
-
-# %%
-# H2 design using h2syn.
-
-# Create augmented plant.
-Baug = np.block([B1, B2])
-Caug = np.block([[C1], [C2]])
-Daug = np.block([[D11, D12], [D21, D22]])
-Paug = control.ss(A, Baug, Caug, Daug)
-
-# Input to h2syn is Paug, number of inputs to controller,
-# and number of outputs from the controller.
-K = control.h2syn(Paug, n_y, n_u)
-
-# Extarct controller ss realization.
-A_K, B_K, C_K, D_K = K.A, K.B, K.C, K.D
-
-# %%
-# Compute closed-loop H2 norm.
-
-# Compute closed-loop system, Tzw(s). See Eq. 4 in [1].
-Azw = np.block([[A + B2 @ D_K @ C2, B2 @ C_K],
-                [B_K @ C2, A_K]])
-Bzw = np.block([[B1 + B2 @ D_K @ D21],
-                [B_K @ D21]])
-Czw = np.block([C1 + D12 @ D_K @ C2, D12 @ C_K])
-Dzw = D11 + D12 @ D_K @ D21
-Tzw = control.ss(Azw, Bzw, Czw, Dzw)
-
-# Compute closed-loop H2 norm via Lyapunov equation.
-# See [2], Lemma 4.4, pg 53.
-Qzw = control.lyap(Azw.T, Czw.T @ Czw)
-nu = np.sqrt(np.trace(Bzw.T @ Qzw @ Bzw))
-print(f'The closed-loop H_2 norm of Tzw(s) is {nu}.')
-# Value is 7.748350599360575, the same as reported in [1].
-
-# %%
diff --git a/examples/scherer_etal_ex7_Hinf_hinfsyn.py b/examples/scherer_etal_ex7_Hinf_hinfsyn.py
deleted file mode 100644
index 97318115..00000000
--- a/examples/scherer_etal_ex7_Hinf_hinfsyn.py
+++ /dev/null
@@ -1,57 +0,0 @@
-"""Hinf design using hinfsyn.
-
-Demonstrate Hinf desing for a SISO plant using h2syn. Based on [1], Ex. 7.
-
-[1] Scherer, Gahinet, & Chilali, "Multiobjective Output-Feedback Control via
-LMI Optimization", IEEE Trans. Automatic Control, Vol. 42, No. 7, July 1997.
-"""
-# %%
-# Packages
-import numpy as np
-import control
-
-# %%
-# State-space system.
-
-# Process model.
-A = np.array([[0, 10, 2],
-              [-1, 1, 0],
-              [0, 2, -5]])
-B1 = np.array([[1],
-               [0],
-               [1]])
-B2 = np.array([[0],
-               [1],
-               [0]])
-
-# Plant output.
-C2 = np.array([[0, 1, 0]])
-D21 = np.array([[2]])
-D22 = np.array([[0]])
-
-# Hinf performance.
-C1 = np.array([[1, 0, 0],
-               [0, 0, 0]])
-D11 = np.array([[0],
-                [0]])
-D12 = np.array([[0],
-                [1]])
-
-# Dimensions.
-n_x, n_u, n_y = 3, 1, 1
-
-# %%
-# Hinf design using hinfsyn.
-
-# Create augmented plant.
-Baug = np.block([B1, B2])
-Caug = np.block([[C1], [C2]])
-Daug = np.block([[D11, D12], [D21, D22]])
-Paug = control.ss(A, Baug, Caug, Daug)
-
-# Input to hinfsyn is Paug, number of inputs to controller,
-# and number of outputs from the controller.
-K, Tzw, gamma, rcond = control.hinfsyn(Paug, n_y, n_u)
-print(f'The closed-loop H_inf norm of Tzw(s) is {gamma}.')
-
-# %%


Test patch:
diff --git a/control/tests/discrete_test.py b/control/tests/discrete_test.py
index 4cf28a21..4415fac0 100644
--- a/control/tests/discrete_test.py
+++ b/control/tests/discrete_test.py
@@ -376,28 +376,51 @@ class TestDiscrete:
     @pytest.mark.parametrize("plantname",
                              ["siso_ss1c",
                               "siso_tf1c"])
-    def test_sample_system_prewarp(self, tsys, plantname):
+    @pytest.mark.parametrize("wwarp",
+                             [.1, 1, 3])
+    @pytest.mark.parametrize("Ts",
+                             [.1, 1])
+    @pytest.mark.parametrize("discretization_type",
+                             ['bilinear', 'tustin', 'gbt'])
+    def test_sample_system_prewarp(self, tsys, plantname, discretization_type, wwarp, Ts):
         """bilinear approximation with prewarping test"""
-        wwarp = 50
-        Ts = 0.025
         # test state space version
         plant = getattr(tsys, plantname)
         plant_fr = plant(wwarp * 1j)
+        alpha = 0.5 if discretization_type == 'gbt' else None
 
-        plant_d_warped = plant.sample(Ts, 'bilinear', prewarp_frequency=wwarp)
+        plant_d_warped = plant.sample(Ts, discretization_type,
+            prewarp_frequency=wwarp, alpha=alpha)
         dt = plant_d_warped.dt
         plant_d_fr = plant_d_warped(np.exp(wwarp * 1.j * dt))
         np.testing.assert_array_almost_equal(plant_fr, plant_d_fr)
 
-        plant_d_warped = sample_system(plant, Ts, 'bilinear',
-                prewarp_frequency=wwarp)
+        plant_d_warped = sample_system(plant, Ts, discretization_type,
+            prewarp_frequency=wwarp, alpha=alpha)
         plant_d_fr = plant_d_warped(np.exp(wwarp * 1.j * dt))
         np.testing.assert_array_almost_equal(plant_fr, plant_d_fr)
 
-        plant_d_warped = c2d(plant, Ts, 'bilinear', prewarp_frequency=wwarp)
+        plant_d_warped = c2d(plant, Ts, discretization_type,
+            prewarp_frequency=wwarp, alpha=alpha)
         plant_d_fr = plant_d_warped(np.exp(wwarp * 1.j * dt))
         np.testing.assert_array_almost_equal(plant_fr, plant_d_fr)
 
+    @pytest.mark.parametrize("plantname",
+                             ["siso_ss1c",
+                              "siso_tf1c"])
+    @pytest.mark.parametrize("discretization_type",
+                             ['euler', 'backward_diff', 'zoh'])
+    def test_sample_system_prewarp_warning(self, tsys, plantname, discretization_type):
+        plant = getattr(tsys, plantname)
+        wwarp = 1
+        Ts = 0.1
+        with pytest.warns(UserWarning, match="prewarp_frequency ignored: incompatible conversion"):
+            plant_d_warped = plant.sample(Ts, discretization_type, prewarp_frequency=wwarp)
+        with pytest.warns(UserWarning, match="prewarp_frequency ignored: incompatible conversion"):
+            plant_d_warped = sample_system(plant, Ts, discretization_type, prewarp_frequency=wwarp)
+        with pytest.warns(UserWarning, match="prewarp_frequency ignored: incompatible conversion"):
+            plant_d_warped = c2d(plant, Ts, discretization_type, prewarp_frequency=wwarp)
+
     def test_sample_system_errors(self, tsys):
         # Check errors
         with pytest.raises(ValueError):
@@ -446,11 +469,11 @@ class TestDiscrete:
         np.testing.assert_array_almost_equal(omega, omega_out)
         np.testing.assert_array_almost_equal(mag_out, np.absolute(H_z))
         np.testing.assert_array_almost_equal(phase_out, np.angle(H_z))
-    
+
     def test_signal_names(self, tsys):
         "test that signal names are preserved in conversion to discrete-time"
-        ssc = StateSpace(tsys.siso_ss1c, 
-            inputs='u', outputs='y', states=['a', 'b', 'c']) 
+        ssc = StateSpace(tsys.siso_ss1c,
+            inputs='u', outputs='y', states=['a', 'b', 'c'])
         ssd = ssc.sample(0.1)
         tfc = TransferFunction(tsys.siso_tf1c, inputs='u', outputs='y')
         tfd = tfc.sample(0.1)
@@ -467,7 +490,7 @@ class TestDiscrete:
         assert ssd.output_labels == ['y']
         assert tfd.input_labels == ['u']
         assert tfd.output_labels == ['y']
-    
+
         # system names and signal name override
         sysc = StateSpace(1.1, 1, 1, 1, inputs='u', outputs='y', states='a')
 
@@ -488,14 +511,14 @@ class TestDiscrete:
         assert sysd_nocopy.find_state('a') is None
 
         # if signal names are provided, they should override those of sysc
-        sysd_newnames = sample_system(sysc, 0.1, 
+        sysd_newnames = sample_system(sysc, 0.1,
             inputs='v', outputs='x', states='b')
         assert sysd_newnames.find_input('v') == 0
         assert sysd_newnames.find_input('u') is None
         assert sysd_newnames.find_output('x') == 0
         assert sysd_newnames.find_output('y') is None
         assert sysd_newnames.find_state('b') == 0
-        assert sysd_newnames.find_state('a') is None        
+        assert sysd_newnames.find_state('a') is None
         # test just one name
         sysd_newnames = sample_system(sysc, 0.1, inputs='v')
         assert sysd_newnames.find_input('v') == 0
diff --git a/control/tests/kwargs_test.py b/control/tests/kwargs_test.py
index f9400954..83026391 100644
--- a/control/tests/kwargs_test.py
+++ b/control/tests/kwargs_test.py
@@ -190,6 +190,7 @@ kwarg_unittest = {
     'tf2io' : test_unrecognized_kwargs,
     'tf2ss' : test_unrecognized_kwargs,
     'sample_system' : test_unrecognized_kwargs,
+    'c2d' : test_unrecognized_kwargs,
     'zpk': test_unrecognized_kwargs,
     'flatsys.point_to_point':
         flatsys_test.TestFlatSys.test_point_to_point_errors,
@@ -210,7 +211,7 @@ kwarg_unittest = {
     'NonlinearIOSystem.__init__':
         interconnect_test.test_interconnect_exceptions,
     'StateSpace.__init__': test_unrecognized_kwargs,
-    'StateSpace.sample': test_unrecognized_kwargs, 
+    'StateSpace.sample': test_unrecognized_kwargs,
     'TimeResponseData.__call__': trdata_test.test_response_copy,
     'TransferFunction.__init__': test_unrecognized_kwargs,
     'TransferFunction.sample': test_unrecognized_kwargs,
diff --git a/control/tests/trdata_test.py b/control/tests/trdata_test.py
index 028e5358..734d3559 100644
--- a/control/tests/trdata_test.py
+++ b/control/tests/trdata_test.py
@@ -196,20 +196,15 @@ def test_response_copy():
     with pytest.raises(ValueError, match="not enough"):
         t, y, x = response_mimo
 
-    # Make sure labels are transferred to the response
-    assert response_siso.output_labels == sys_siso.output_labels
-    assert response_siso.state_labels == sys_siso.state_labels
-    assert response_siso.input_labels == sys_siso.input_labels
-    assert response_mimo.output_labels == sys_mimo.output_labels
-    assert response_mimo.state_labels == sys_mimo.state_labels
-    assert response_mimo.input_labels == sys_mimo.input_labels
-
-    # Check relabelling
+    # Labels
+    assert response_mimo.output_labels is None
+    assert response_mimo.state_labels is None
+    assert response_mimo.input_labels is None
     response = response_mimo(
         output_labels=['y1', 'y2'], input_labels='u',
-        state_labels=["x%d" % i for i in range(4)])
+        state_labels=["x[%d]" % i for i in range(4)])
     assert response.output_labels == ['y1', 'y2']
-    assert response.state_labels == ['x0', 'x1', 'x2', 'x3']
+    assert response.state_labels == ['x[0]', 'x[1]', 'x[2]', 'x[3]']
     assert response.input_labels == ['u']
 
     # Unknown keyword
@@ -236,17 +231,6 @@ def test_trdata_labels():
     np.testing.assert_equal(
         response.input_labels, ["u[%d]" % i for i in range(sys.ninputs)])
 
-    # Make sure the selected input and output are both correctly transferred to the response
-    for nu in range(sys.ninputs):
-        for ny in range(sys.noutputs):
-            step_response = ct.step_response(sys, T, input=nu, output=ny)
-            assert step_response.input_labels == [sys.input_labels[nu]]
-            assert step_response.output_labels == [sys.output_labels[ny]]
-
-            init_response = ct.initial_response(sys, T, input=nu, output=ny)
-            assert init_response.input_labels == None
-            assert init_response.output_labels == [sys.output_labels[ny]]
-
 
 def test_trdata_multitrace():
     #