# Path: harold/_arg_utils.py
# def _check_for_state_or_transfer(G, custom_msg=None):
#     """
#     A helper function to assert whether the given object is a State or
#     Transfer.
#     """
#     if not isinstance(G, (State, Transfer)):
#         if custom_msg is None:
#             custom_msg = 'The argument should be a Transfer or a State'\
#                          '. Instead found a {}'.format(type(G).__qualname__)
#         raise ValueError(custom_msg)
# 
# Path: harold/_time_domain.py
# def simulate_step_response(sys, t=None):
#     """
#     Compute the linear model response to an Heaviside function (or all-ones
#     array) sampled at given time instances.
# 
#     If the time array is omitted then a time sequence is generated based on
#     the poles of the model.
# 
#     Parameters
#     ----------
#     sys : {State, Transfer}
#         The system model to be simulated
#     t : array_like
#         The real-valued sequence to be used for the evolution of the system.
#         The values should be equally spaced otherwise an error is raised. For
#         discrete time models increments different than the sampling period also
#         raises an error. On the other hand for discrete models this can be
#         omitted and a time sequence will be generated automatically.
# 
#     Returns
#     -------
#     yout : ndarray
#         The resulting response array. The array is 1D if sys is SISO and
#         has p columns if sys has p outputs. If there are also m inputs the
#         array is 3D array with the shape (<num of samples>, p, m)
#     tout : ndarray
#         The time sequence used in the simulation. If the parameter t is not
#         None then a copy of t is given.
# 
#     """
#     _check_for_state_or_transfer(sys)
#     # Always works with State Models
#     try:
#         _check_for_state(sys)
#     except ValueError:
#         sys = transfer_to_state(sys)
# 
#     if t is None:
#         tf, ts = _compute_tfinal_and_dt(sys)
#         t = np.arange(0, tf+ts, ts, dtype=float)
#     else:
#         t, ts = _check_custom_time_input(t)
# 
#     m = sys.shape[1]
#     u = np.ones([len(t), m], dtype=float)
# 
#     return simulate_linear_system(sys, u=u, t=t, per_channel=1)
# 
# def simulate_impulse_response(sys, t=None):
#     """
#     Compute the linear model response to an Dirac delta pulse (or all-zeros
#     array except the first sample being 1/dt at each channel) sampled at given
#     time instances.
# 
#     If the time array is omitted then a time sequence is generated based on
#     the poles of the model.
# 
#     Parameters
#     ----------
#     sys : {State, Transfer}
#         The system model to be simulated
#     t : array_like
#         The real-valued sequence to be used for the evolution of the system.
#         The values should be equally spaced otherwise an error is raised. For
#         discrete time models increments different than the sampling period also
#         raises an error. On the other hand for discrete models this can be
#         omitted and a time sequence will be generated automatically.
# 
#     Returns
#     -------
#     yout : ndarray
#         The resulting response array. The array is 1D if sys is SISO and
#         has p columns if sys has p outputs. If there are also m inputs the
#         array is 3D array with the shape (<num of samples>, p, m)
#     tout : ndarray
#         The time sequence used in the simulation. If the parameter t is not
#         None then a copy of t is given.
# 
#     """
#     _check_for_state_or_transfer(sys)
#     # Always works with State Models
#     try:
#         _check_for_state(sys)
#     except ValueError:
#         sys = transfer_to_state(sys)
# 
#     if t is None:
#         tf, ts = _compute_tfinal_and_dt(sys, is_step=False)
#         t = np.arange(0, tf+ts, ts, dtype=float)
#     else:
#         t, ts = _check_custom_time_input(t)
# 
#     m = sys.shape[1]
#     u = np.zeros([len(t), m], dtype=float)
#     u[0] = 1./ts
# 
#     return simulate_linear_system(sys, u=u, t=t, per_channel=1)


import numpy as np
import matplotlib.pyplot as plt
import matplotlib as mpl
from ._arg_utils import _check_for_state_or_transfer
from ._time_domain import simulate_step_response, simulate_impulse_response
    max_p, max_m = np.max(np.array([x.shape for x in syslist]), axis=0)

    # Put some more space between columns to avoid ticklabel placement clashes
    gridspec_kw = kwargs.pop('gridspec_kw', None)
    if gridspec_kw is None:
        gridspec_kw = {'wspace': 0.5}
    else:
        wspace = gridspec_kw.get('wspace', 0.5)
        gridspec_kw['wspace'] = wspace

    # MIMO plots needs a bit bigger figure, offer a saner default
    figsize = kwargs.pop('figsize', None)
    if figsize is None:
        figsize = (6.0 + 1.2*(max_m - 1), 4 + 1.5*(max_p - 1))

    if style is None:
        style = mpl.cycler(mpl.rcParams['axes.prop_cycle'])

    # Create fig and axes
    fig, axs = plt.subplots(max_p, max_m, sharex=True, sharey=True,
                            squeeze=False,
                            gridspec_kw=gridspec_kw,
                            figsize=figsize,
                            **kwargs)

    # TODO: If multiple systems are given and no t given we need to find a
    # compromise for the shorter ones. For now, just plot everything on top of
    # each other independently and cut to the shortest. Yes, it sucks; I know.
    tmin = np.inf
    for sys in syslist:

# Path: harold/_arg_utils.py
# def _check_for_state_or_transfer(G, custom_msg=None):
#     """
#     A helper function to assert whether the given object is a State or
#     Transfer.
#     """
#     if not isinstance(G, (State, Transfer)):
#         if custom_msg is None:
#             custom_msg = 'The argument should be a Transfer or a State'\
#                          '. Instead found a {}'.format(type(G).__qualname__)
#         raise ValueError(custom_msg)
# 
# Path: harold/_time_domain.py
# def simulate_step_response(sys, t=None):
#     """
#     Compute the linear model response to an Heaviside function (or all-ones
#     array) sampled at given time instances.
# 
#     If the time array is omitted then a time sequence is generated based on
#     the poles of the model.
# 
#     Parameters
#     ----------
#     sys : {State, Transfer}
#         The system model to be simulated
#     t : array_like
#         The real-valued sequence to be used for the evolution of the system.
#         The values should be equally spaced otherwise an error is raised. For
#         discrete time models increments different than the sampling period also
#         raises an error. On the other hand for discrete models this can be
#         omitted and a time sequence will be generated automatically.
# 
#     Returns
#     -------
#     yout : ndarray
#         The resulting response array. The array is 1D if sys is SISO and
#         has p columns if sys has p outputs. If there are also m inputs the
#         array is 3D array with the shape (<num of samples>, p, m)
#     tout : ndarray
#         The time sequence used in the simulation. If the parameter t is not
#         None then a copy of t is given.
# 
#     """
#     _check_for_state_or_transfer(sys)
#     # Always works with State Models
#     try:
#         _check_for_state(sys)
#     except ValueError:
#         sys = transfer_to_state(sys)
# 
#     if t is None:
#         tf, ts = _compute_tfinal_and_dt(sys)
#         t = np.arange(0, tf+ts, ts, dtype=float)
#     else:
#         t, ts = _check_custom_time_input(t)
# 
#     m = sys.shape[1]
#     u = np.ones([len(t), m], dtype=float)
# 
#     return simulate_linear_system(sys, u=u, t=t, per_channel=1)
# 
# def simulate_impulse_response(sys, t=None):
#     """
#     Compute the linear model response to an Dirac delta pulse (or all-zeros
#     array except the first sample being 1/dt at each channel) sampled at given
#     time instances.
# 
#     If the time array is omitted then a time sequence is generated based on
#     the poles of the model.
# 
#     Parameters
#     ----------
#     sys : {State, Transfer}
#         The system model to be simulated
#     t : array_like
#         The real-valued sequence to be used for the evolution of the system.
#         The values should be equally spaced otherwise an error is raised. For
#         discrete time models increments different than the sampling period also
#         raises an error. On the other hand for discrete models this can be
#         omitted and a time sequence will be generated automatically.
# 
#     Returns
#     -------
#     yout : ndarray
#         The resulting response array. The array is 1D if sys is SISO and
#         has p columns if sys has p outputs. If there are also m inputs the
#         array is 3D array with the shape (<num of samples>, p, m)
#     tout : ndarray
#         The time sequence used in the simulation. If the parameter t is not
#         None then a copy of t is given.
# 
#     """
#     _check_for_state_or_transfer(sys)
#     # Always works with State Models
#     try:
#         _check_for_state(sys)
#     except ValueError:
#         sys = transfer_to_state(sys)
# 
#     if t is None:
#         tf, ts = _compute_tfinal_and_dt(sys, is_step=False)
#         t = np.arange(0, tf+ts, ts, dtype=float)
#     else:
#         t, ts = _check_custom_time_input(t)
# 
#     m = sys.shape[1]
#     u = np.zeros([len(t), m], dtype=float)
#     u[0] = 1./ts
# 
#     return simulate_linear_system(sys, u=u, t=t, per_channel=1)


# Path: harold/_time_domain_plots.py
import numpy as np
import matplotlib.pyplot as plt
import matplotlib as mpl
from ._arg_utils import _check_for_state_or_transfer
from ._time_domain import simulate_step_response, simulate_impulse_response

    max_p, max_m = np.max(np.array([x.shape for x in syslist]), axis=0)

    # Put some more space between columns to avoid ticklabel placement clashes
    gridspec_kw = kwargs.pop('gridspec_kw', None)
    if gridspec_kw is None:
        gridspec_kw = {'wspace': 0.5}
    else:
        wspace = gridspec_kw.get('wspace', 0.5)
        gridspec_kw['wspace'] = wspace

    # MIMO plots needs a bit bigger figure, offer a saner default
    figsize = kwargs.pop('figsize', None)
    if figsize is None:
        figsize = (6.0 + 1.2*(max_m - 1), 4 + 1.5*(max_p - 1))

    if style is None:
        style = mpl.cycler(mpl.rcParams['axes.prop_cycle'])

    # Create fig and axes
    fig, axs = plt.subplots(max_p, max_m, sharex=True, sharey=True,
                            squeeze=False,
                            gridspec_kw=gridspec_kw,
                            figsize=figsize,
                            **kwargs)

    # TODO: If multiple systems are given and no t given we need to find a
    # compromise for the shorter ones. For now, just plot everything on top of
    # each other independently and cut to the shortest. Yes, it sucks; I know.
    tmin = np.inf
    for sys in syslist:

        yout, tout = simulate_impulse_response(sys, t=t)