def chebfun(f=None, domain=[-1,1], N=None, chebcoeff=None,):
    """
    Create a Chebyshev polynomial approximation of the function $f$ on the interval :math:`[-1, 1]`.
    :param callable f: Python, Numpy, or Sage function
    :param int N: (default = None)  specify number of interpolating points
    :param np.array chebcoeff: (default = np.array(0)) specify the coefficients
    """

    # Chebyshev coefficients
    if chebcoeff is not None:
        return Chebfun.from_coeff(chebcoeff, domain)

    # another instance
    if isinstance(f, Polyfun):
        return Chebfun.from_fun(f)

    # callable
    if hasattr(f, '__call__'):
        return Chebfun.from_function(f, domain, N)

    # from here on, assume that f is None, or iterable
    if np.isscalar(f):
        f = [f]

    try:
        iter(f) # interpolation values provided
    except TypeError:
        pass
    else:
        return Chebfun(f, domain)

    raise TypeError('Impossible to initialise the object from an object of type {}'.format(type(f)))
Create a Chebyshev polynomial approximation of the function $f$ on the interval :math:`[-1, 1]`.
    :param callable f: Python, Numpy, or Sage function
    :param int N: (default = None)  specify number of interpolating points
    :param np.array chebcoeff: (default = np.array(0)) specify the coefficients