def time_frequency(gram, frequencies, times, fs, function=np.sin, length=None,
                   n_dec=1):
    if times.ndim == 1:
        times = util.boundaries_to_intervals(times)
    if length is None:
        length = int(times[-1, 1] * fs)
    times, _ = util.adjust_intervals(times, t_max=length)
    n_times = gram.shape[1]
    times = times[:n_times]
    def _fast_synthesize(frequency):
        frequency = np.round(frequency, n_dec)
        n_samples = int(10.0**n_dec * fs)
        short_signal = function(2.0 * np.pi * np.arange(n_samples) *
                                frequency / fs)
        n_repeats = int(np.ceil(length/float(short_signal.shape[0])))
        long_signal = as_strided(short_signal,
                                 shape=(n_repeats, len(short_signal)),
                                 strides=(0, short_signal.itemsize))
        return long_signal.flat
    def _const_interpolator(value):
        def __interpolator(x):
            return value
        return __interpolator
    gram = np.maximum(gram, 0)
    output = np.zeros(length)
    time_centers = np.mean(times, axis=1) * float(fs)
    for n, frequency in enumerate(frequencies):
        wave = _fast_synthesize(frequency)
        if len(time_centers) > 1:
            gram_interpolator = interp1d(
                time_centers, gram[n, :],
                kind='linear', bounds_error=False,
                fill_value=0.0)
        else:
            gram_interpolator = _const_interpolator(gram[n, 0])
        for m, (start, end) in enumerate((times * fs).astype(int)):
            start, end = max(start, 0), min(end, length)
            output[start:end] += (
                wave[start:end] * gram_interpolator(np.arange(start, end)))
    norm = np.abs(output).max()
    if norm >= np.finfo(output.dtype).tiny:
        output /= norm
    return output