text: import numpy as np
  
  from scipy.fftpack.realtransforms import dct
  
  def hz2mel(f):
      """Convert an array of frequency in Hz into mel."""
      return 1127.01048 * np.log(f/700 +1)
  
  def mel2hz(m):
      """Convert an array of frequency in Hz into mel."""
      return (np.exp(m / 1127.01048) - 1) * 700
  
  
  def trfbank(fs, nfft, lowfreq, linsc, logsc, nlinfilt, nlogfilt):
      """Compute triangular filterbank for MFCC computation."""
      # Total number of filters
      nfilt = nlinfilt + nlogfilt
  
      #------------------------
      # Compute the filter bank
      #------------------------
      # Compute start/middle/end points of the triangular filters in spectral
      # domain
      freqs = np.zeros(nfilt+2)
      freqs[:nlinfilt] = lowfreq + np.arange(nlinfilt) * linsc
      freqs[nlinfilt:] = freqs[nlinfilt-1] * logsc ** np.arange(1, nlogfilt + 3)
      heights = 2./(freqs[2:] - freqs[0:-2])
  
      # Compute filterbank coeff (in fft domain, in bins)
      fbank = np.zeros((nfilt, nfft))
      # FFT bins (in Hz)
      nfreqs = np.arange(nfft) / (1. * nfft) * fs
      for i in range(nfilt):
          low = freqs[i]
          cen = freqs[i+1]
          hi = freqs[i+2]
  
          lid = np.arange(np.floor(low * nfft / fs) + 1,
                          np.floor(cen * nfft / fs) + 1, dtype=np.int)
          lslope = heights[i] / (cen - low)
          rid = np.arange(np.floor(cen * nfft / fs) + 1,
                          np.floor(hi * nfft / fs) + 1, dtype=np.int)
          rslope = heights[i] / (hi - cen)
          fbank[i][lid] = lslope * (nfreqs[lid] - low)
          fbank[i][rid] = rslope * (hi - nfreqs[rid])
  
      return fbank, freqs
  
  def mfcc(input, nceps=13):
      """Compute Mel Frequency Cepstral Coefficients.
  
      Parameters
      ----------
      input: ndarray
          input spectrogram from which the coefficients are computed
  
      Returns
      -------
      ceps: ndarray
          Mel-cepstrum coefficients
      mspec: ndarray
          Log-spectrum in the mel-domain.
  
      Notes
      -----
      MFCC are computed as follows:
          * Pre-processing in time-domain (pre-emphasizing)
          * Compute the spectrum amplitude by windowing with a Hamming window
          * Filter the signal in the spectral domain with a triangular
          filter-bank, whose filters are approximatively linearly spaced on the
          mel scale, and have equal bandwith in the mel scale
          * Compute the DCT of the log-spectrum
  
      This is based on the talkbox module:
      http://pydoc.net/Python/scikits.talkbox/0.2.4.dev/scikits.talkbox.features.mfcc/
  
      References
      ----------
      .. [1] S.B. Davis and P. Mermelstein, "Comparison of parametric
             representations for monosyllabic word recognition in continuously
             spoken sentences", IEEE Trans. Acoustics. Speech, Signal Proc.
             ASSP-28 (4): 357-366, August 1980."""
  
  
      nfft = input.metadata.sampling_configuration.dft_length
      fs = input.metadata.sampling_configuration.fs
      over = input.metadata.sampling_configuration.window_length \
          - input.metadata.sampling_configuration.window_step
  
      #lowfreq = 400 / 3.
      lowfreq = 133.33
      #highfreq = 6855.4976
      linsc = 200/3.
      logsc = 1.0711703
  
      nlinfil = 13
      nlogfil = 27
      nlinfil + nlogfil
  
      fbank = trfbank(fs, nfft, lowfreq, linsc, logsc, nlinfil, nlogfil)[0]
      fbank = fbank.T[0:input.data.shape[0], :]
  
      mspec = np.log10(np.maximum(np.dot(fbank.T, input.data), 0.0000001)).T
  
      # Use the DCT to 'compress' the coefficients (spectrum -> cepstrum domain)
      ceps = dct(mspec, type=2, norm='ortho', axis=-1)[:, :nceps]
  
      return ceps
  