dmastropole

# Plot MFCCs
plt.imshow(mfccs,
           aspect='auto',
           origin='lower',
           interpolation='none',
           cmap=plt.cm.Blues)
plt.title("MFCCs")
plt.xlabel("Frames")
plt.ylabel("MFCCs")

dmastropole

# Perform a DCT
from essentia.standard import DCT

dct = DCT(inputSize=n_bands, outputSize=13)
mfccs = dct(log_mels)

plt.bar(np.arange(len(mfccs)), mfccs, align='center')
plt.title('Mel-Frequency Cepstrum')
plt.xlabel('Cosines')
plt.ylabel('Coefficients')

dmastropole

# Compute the mel bands
n_bands = 25

# The "+ 2" will make sense in a bit!
m = np.linspace(mel_formula(f_low), mel_formula(f_high), n_bands + 2)  
print(m)

dmastropole

from essentia.standard import Windowing, PowerSpectrum

# Compute the spectrum of a frame
w = Windowing(type='hann')
spectrum = PowerSpectrum(size=1024)
frame = audio[0:1024]
spec = spectrum(w(frame))

# Plot the spectrum of a frame
plt.plot(spec)

dmastropole

from essentia.standard import MelBands

# Compute the mel band powers
melbands = MelBands(lowFrequencyBound=f_low,
                    highFrequencyBound=f_high,
                    inputSize=1024,
                    numberBands=n_bands,
                    type='magnitude',  # We already computed the power.
                    sampleRate=44100)
mels = melbands(spec)

dmastropole

from essentia.standard import UnaryOperator

# Convert to decibels
log10 = UnaryOperator(type='log10')
log_mels = log10(mels)

# Plot the mel band powers
plt.bar(np.arange(len(log_mels)), log_mels, align='center')
plt.title('Log of the Mel Band Powers')
plt.xlabel('Mel Bands')

dmastropole

from essentia.standard import MFCC, FrameGenerator

mfcc = MFCC(highFrequencyBound=f_high,
            lowFrequencyBound=f_low,
            inputSize=1024,
            numberBands=n_bands,
            numberCoefficients=13,
            type='magnitude',
            sampleRate=44100)

dmastropole

for n in range(n_bands):
    plt.plot(f[n:n+3], [0, 1, 0])  # Does the "+ 2" make sense now?

plt.title('Mel Filterbanks')
plt.xlabel('Frequency (Hz)')

dmastropole

def inverse_mel_formula(m):
    return 700 * (np.power(10, m/2595) - 1)

# Convert back to frequency
f = inverse_mel_formula(m)
print(f)

dmastropole

# Define function for computing mels
def mel_formula(f):
    return 2595 * np.log(1 + f/700)
    
# Frequency range
f_low = 40
f_high = 44100 / 2  # Nyquist frequency
f = np.linspace(f_low, f_high)

# Plot relationship between frequency and mels