% Script which computes the spectral spread and skewness features
% and visualizes them for some input signals.
clc
clear all

sampFreq = 44100;
frameLen = 2048;
frameShift = 64;

% these folders contain isolated tones in the note range A3 to G5:
currentFolder = pwd;
clarinetFolder = [pwd '/McGill/clari-bfl/'];
oboeFolder = [pwd '/McGill/oboe/'];
violinFolder = [pwd '/McGill/violin-bv/'];
folderList = {clarinetFolder,oboeFolder,violinFolder};

colorList = {'b','r','g'};

for k = 1:length(folderList)
    cd(folderList{k})
    fileList = dir(pwd);
    fileList = fileList(4:end);

    [signal,fs] = audioread(fileList(k).name);
    if fs ~= sampFreq
    signal = resample(signal,sampFreq,fs);
    end
    signal = signal/std(signal);

    % segment signal and compute FFT to obtain STFT
    signalFrames = buffer(signal,frameLen,frameLen-frameShift,'nodelay');
    stftSpecFull = fft(signalFrames);
    stftSpecHalf = abs(stftSpecFull(1:frameLen/2+1,:));

% solution: begin 
    % spectral centroid

    % spectral spread

    % spectral skewness

% solution: end 
    
    % plots (skewness vs. spread)
    figure(2);
    hold on
    scatter(specSpread,specSkew,[],colorList{k});axis xy
    xlabel('Spectral spread [Hz]')
    ylabel('Spectral skewness')
    grid on
end

legend('Clarinet','Oboe','Violin')