How to use the control.lti.isdtime function in control

def _evalfr(self, omega):
        """Evaluate a SS system's transfer function at a single frequency"""
        # Figure out the point to evaluate the transfer function
        if isdtime(self, strict=True):
            dt = timebase(self)
            s = exp(1.j * omega * dt)
            if omega * dt > math.pi:
                warn("_evalfr: frequency evaluation above Nyquist frequency")
        else:
            s = omega * 1.j

        return self.horner(s)

def evalfr(self, omega):
        """Evaluate a SS system's transfer function at a single frequency.

        self.evalfr(omega) returns the value of the transfer function matrix with
        input value s = i * omega.

        """
        # Figure out the point to evaluate the transfer function
        if isdtime(self, strict=True):
            dt = timebase(self)
            s = exp(1.j * omega * dt)
            if (omega * dt > pi):
                warn("evalfr: frequency evaluation above Nyquist frequency")
        else:
            s = omega * 1.j

        fresp = self.C * solve(s * eye(self.states) - self.A,
            self.B) + self.D

        return array(fresp)

if isctime(sys):
        if not hasattr(sp.integrate, 'solve_ivp'):
            raise NameError("scipy.integrate.solve_ivp not found; "
                            "use SciPy 1.0 or greater")
        soln = sp.integrate.solve_ivp(ivp_rhs, (T0, Tf), X0, t_eval=T,
                                      method=method, vectorized=False)

        # Compute the output associated with the state (and use sys.out to
        # figure out the number of outputs just in case it wasn't specified)
        u = U[0] if len(U.shape) == 1 else U[:, 0]
        y = np.zeros((np.shape(sys._out(T[0], X0, u))[0], len(T)))
        for i in range(len(T)):
            u = U[i] if len(U.shape) == 1 else U[:, i]
            y[:, i] = sys._out(T[i], soln.y[:, i], u)

    elif isdtime(sys):
        # Make sure the time vector is uniformly spaced
        dt = T[1] - T[0]
        if not np.allclose(T[1:] - T[:-1], dt):
            raise ValueError("Parameter ``T``: time values must be "
                             "equally spaced.")

        # Make sure the sample time matches the given time
        if (sys.dt is not True):
            # Make sure that the time increment is a multiple of sampling time

            # TODO: add back functionality for undersampling
            # TODO: this test is brittle if dt =  sys.dt
            # First make sure that time increment is bigger than sampling time
            # if dt < sys.dt:
            #     raise ValueError("Time steps ``T`` must match sampling time")

reports the value of the magnitude, phase, and angular frequency of
        the transfer function matrix evaluated at s = i * omega, where omega
        is a list of angular frequencies, and is a sorted version of the input
        omega.

        """

        # Preallocate outputs.
        numfreq = len(omega)
        mag = empty((self.outputs, self.inputs, numfreq))
        phase = empty((self.outputs, self.inputs, numfreq))

        # Figure out the frequencies
        omega.sort()
        if isdtime(self, strict=True):
            dt = timebase(self)
            slist = np.array([exp(1.j * w * dt) for w in omega])
            if max(omega) * dt > pi:
                warn("freqresp: frequency evaluation above Nyquist frequency")
        else:
            slist = np.array([1j * w for w in omega])

        # Compute frequency response for each input/output pair
        for i in range(self.outputs):
            for j in range(self.inputs):
                fresp = (polyval(self.num[i][j], slist) /
                         polyval(self.den[i][j], slist))
                mag[i, j, :] = abs(fresp)
                phase[i, j, :] = angle(fresp)

        return mag, phase, omega

Examples
    --------
    >>> T, yout, xout = forced_response(sys, T, u, X0)
    """
    if not isinstance(sys, Lti):
        raise TypeError('Parameter ``sys``: must be a ``Lti`` object. '
                        '(For example ``StateSpace`` or ``TransferFunction``)')
    sys = _convertToStateSpace(sys) 
    A, B, C, D = np.asarray(sys.A), np.asarray(sys.B), np.asarray(sys.C), \
                 np.asarray(sys.D)
#    d_type = A.dtype
    n_states = A.shape[0]
    n_inputs = B.shape[1]

    # Set and/or check time vector in discrete time case
    if isdtime(sys, strict=True):
        if T == None:
            if U == None:
                raise ValueError('Parameters ``T`` and ``U`` can\'t both be'
                                 'zero for discrete-time simulation')
            # Set T to integers with same length as U
            T = range(len(U))
        else:
            # Make sure the input vector and time vector have same length
            # TODO: allow interpolation of the input vector
            if len(U) != len(T):
                ValueError('Pamameter ``T`` must have same length as'
                           'input vector ``U``')

    # Test if T has shape (n,) or (1, n);
    # T must be array-like and values must be increasing.
    # The length of T determines the length of the input vector.

gm, pm, sm, wg, wp, ws: float
        Gain margin gm, phase margin pm, stability margin sm, and 
        associated crossover
        frequencies wg, wp, and ws of SISO open-loop. If more than
        one crossover frequency is detected, returns the lowest corresponding
        margin. 
    """
    #TODO do this precisely without the effects of discretization of frequencies?
    #TODO assumes SISO
    #TODO unit tests, margin plot

    if (not getattr(sysdata, '__iter__', False)):
        sys = sysdata

        # TODO: implement for discrete time systems
        if (isdtime(sys, strict=True)):
            raise NotImplementedError("Function not implemented in discrete time")

        mag, phase, omega = bode(sys, deg=deg, Plot=False)
    elif len(sysdata) == 3:
        # TODO: replace with FRD object type?
        mag, phase, omega = sysdata
    else: 
        raise ValueError("Margin sysdata must be either a linear system or a 3-sequence of mag, phase, omega.")
        
    if deg:
        cycle = 360. 
        crossover = 180. 
    else:
        cycle = 2 * np.pi
        crossover = np.pi

sys = _convertToStateSpace(sys)
    A, B, C, D = np.asarray(sys.A), np.asarray(sys.B), np.asarray(sys.C), \
        np.asarray(sys.D)
#    d_type = A.dtype
    n_states = A.shape[0]
    n_inputs = B.shape[1]
    n_outputs = C.shape[0]

    # Convert inputs to numpy arrays for easier shape checking
    if U is not None:
        U = np.asarray(U)
    if T is not None:
        T = np.asarray(T)

    # Set and/or check time vector in discrete time case
    if isdtime(sys, strict=True):
        if T is None:
            if U is None:
                raise ValueError('Parameters ``T`` and ``U`` can\'t both be'
                                 'zero for discrete-time simulation')
            # Set T to equally spaced samples with same length as U
            if len(U.shape) == 1:
                n_steps = U.shape[0]
            else:
                n_steps = U.shape[1]
            T = np.array(range(n_steps)) * (1 if sys.dt is True else sys.dt)
        else:
            # Make sure the input vector and time vector have same length
            if (len(U.shape) == 1 and U.shape[0] != T.shape[0]) or (len(U.shape) > 1 and U.shape[1] != T.shape[0]):
                ValueError('Pamameter ``T`` must have same elements as'
                           ' the number of columns in input array ``U``')

def _evalfr(self, omega):
        """Evaluate a SS system's transfer function at a single frequency"""
        # Figure out the point to evaluate the transfer function
        if isdtime(self, strict=True):
            dt = timebase(self)
            s = exp(1.j * omega * dt)
            if (omega * dt > math.pi):
                warn("_evalfr: frequency evaluation above Nyquist frequency")
        else:
            s = omega * 1.j

        return self.horner(s)