ComputeMelSpectrogram

computes a mel spectrogram from the audio data
>
> @param x: time domain sample data, dimension channels X samples
> @param f_s: sample rate of audio data
> @param bLogarithmic: levels (true) or magnitudes (false)
> @param afWindow: FFT window of length iBlockLength (default: hann), can be [] empty
> @param iBlockLength: internal block length (default: 4096 samples)
> @param iHopLength: internal hop length (default: 2048 samples)
> @param iNumMelBands: num of frequency bands (default: 128)
> @param fMax: maximum frequency (default: f_s/2)
>
> @retval M spectrogram
> @retval f frequency bands
> @retval t time stamps
 ======================================================================

0001 %computes a mel spectrogram from the audio data
0002 %>
0003 %> @param x: time domain sample data, dimension channels X samples
0004 %> @param f_s: sample rate of audio data
0005 %> @param bLogarithmic: levels (true) or magnitudes (false)
0006 %> @param afWindow: FFT window of length iBlockLength (default: hann), can be [] empty
0007 %> @param iBlockLength: internal block length (default: 4096 samples)
0008 %> @param iHopLength: internal hop length (default: 2048 samples)
0009 %> @param iNumMelBands: num of frequency bands (default: 128)
0010 %> @param fMax: maximum frequency (default: f_s/2)
0011 %>
0012 %> @retval M spectrogram
0013 %> @retval f frequency bands
0014 %> @retval t time stamps
0015 % ======================================================================
0016 function [M, f_c, t] = ComputeMelSpectrogram (x, f_s, bLogarithmic, afWindow, iBlockLength, iHopLength, iNumMelBands, fMax)
0017 
0018     % set default parameters if necessary
0019     if (nargin < 8)
0020         fMax = f_s/2;
0021     end
0022     if (nargin < 7)
0023         iNumMelBands = 128;
0024     end
0025     if (nargin < 6)
0026         iHopLength = 2048;
0027     end
0028     if (nargin < 5)
0029         iBlockLength = 4096;
0030     end
0031     if (nargin < 4 || isempty(afWindow))
0032         afWindow = hann(iBlockLength,'periodic');
0033     end
0034     if (nargin < 3)
0035         bLogarithmic = true;
0036     end
0037     if (nargin < 1)
0038         load handel;
0039         x = y;
0040         f_s = Fs;
0041         clear y, Fs;
0042     end
0043     
0044     if (length(afWindow) ~= iBlockLength)
0045         error('window length mismatch');
0046     end
0047     
0048     % pre-processing: down-mixing
0049     x = ToolDownmix(x);
0050     
0051     % pre-processing: normalization
0052     x = ToolNormalizeAudio(x);
0053  
0054     % in the real world, we would do this block by block...
0055     [X, f, t] = ComputeSpectrogram(x, ...
0056                                 f_s, ...
0057                                 afWindow, ...
0058                                 iBlockLength, ...
0059                                 iHopLength);
0060 
0061     % compute mel filters
0062     [H, f_c] = getMelFb_I(iBlockLength, f_s, iNumMelBands, fMax);
0063 
0064     % apply mel filters
0065     M = H*X;
0066     
0067     % amplitude to level
0068     if (bLogarithmic)
0069         M = 20*log10(M+1e-12);
0070     end
0071 end
0072 
0073 function [H,f_c] = getMelFb_I (iFftLength, f_s, iNumFilters, f_max)
0074 
0075     % initialization
0076     f_min   = 0;
0077     f_max   = min(f_max,f_s / 2);
0078     f_fft   = linspace(0, f_s/2, iFftLength/2+1);
0079     H       = zeros(iNumFilters, length(f_fft));
0080 
0081     % compute band center freqs
0082     mel_min = ToolFreq2Mel(f_min);
0083     mel_max = ToolFreq2Mel(f_max);
0084     f_mel   = ToolMel2Freq(linspace(mel_min, mel_max, iNumFilters+2));
0085 
0086     f_l = f_mel(1:iNumFilters);
0087     f_c = f_mel(2:iNumFilters+1);
0088     f_u = f_mel(3:iNumFilters+2);
0089 
0090     afFilterMax = 2./(f_u-f_l);
0091 
0092     % compute the transfer functions
0093     for (c = 1:iNumFilters)
0094         H(c,:) = ...
0095             (f_fft > f_l(c) & f_fft <= f_c(c)).* ...
0096             afFilterMax(c).*(f_fft-f_l(c))/(f_c(c)-f_l(c)) + ...
0097             (f_fft > f_c(c) & f_fft < f_u(c)).* ...
0098             afFilterMax(c).*(f_u(c)-f_fft)/(f_u(c)-f_c(c));
0099     end
0100 end