comparison functions/fe_funcs/FE.m @ 0:2fadb31a9d55 tip

Import code by Vuegen et al

+% Yet another feature extractor, this time with CHiME in mind.
+% Some cleanup, more support for customised audio parameters.
+% Updated 9th Aug 2011
+%
+% Outputs have been changed. Currently no logarithms are taken here any 
+% more. 
+%
+% Input:
+% - 'sam' is a the audio, either as column or row channels vectors. 
+%   (Longer dimension is treated as time, shorter as channel count.)
+% - 'audioconf' is as defined in getconfigs.m . All of its parameters
+%   ARE respected now, so pass a temporary, edited copy if you want
+%   to change the behaviour.
+%
+% Output:
+% - 'feats' is a [bands x frames x featchannels] array of mel features.
+%   If audioconf.melbands is zero, FFT magnitudes are returned instead.
+% - 'energies' is a [frames x audiochannels] matrix of frame energies
+% - 'frameaudio' is a [framelen x frames x audiochannels] array of chopped
+%   audio data (with preprocessing but without the window function).
+% - 'frameFFT' is an [FFTlow x frames x audiochannels] array of frame FFTs.
+%   The windowing function has been applied, and the result has been
+%   truncated to Nfft/2 + 1 bands. However, no abs is taken. You can
+%   do this in the calling function, or pick the abs values from 'feats'
+%   by using zero melbands.
+% 
+% The main feature output respects audioconf.featchannels, which should
+% be either the same as audioconf.channels (the number of input streams)
+% or 1 (downmixed to mono by taking the mean of feature channels). Other 
+% outputs use original audio channels, because their averaging is not as 
+% well defined. Note that there is a significant difference between
+% averaging the audio (causing waveform level phase attenuation) and the
+% abs-FFT or Mel features (phase-invariant energy mean). If the former is
+% what you need, downmix the audio in the calling function.
+%
+% Some warnings are shown if audio parameters are missing or they do not
+% match with the data.
+
+function [feats, energies, frameaudio, frameFFT] = FE(sam, audioconf)
+
+verbose = 0;
+
+% Default config. To guarantee intended operation, you should always
+% pass your own, though.
+
+defconf.channels = 2;         % input channels, in CHiME always 2
+defconf.featchannels = 1;     % feature level channels
+defconf.fs = 16000;           % sampling rate for internal processing
+defconf.maxf = 8000;          % maximum frequency to be considered
+defconf.minf = 64;            % maximum frequency to be considered
+defconf.melbands = 26;        % mel band count (0 to disable)
+defconf.framelen_ms = 25;     % millisecond length of each frame
+defconf.framestep_ms = 10;    % millisecond step between frames
+defconf.windowfunc = 'hamming';  % window function name
+defconf.preemphasis = 0.97;   % 0 to disable
+defconf.dcremoval = true;     % DC removal in the feature extractor
+defconf.Nfft = 0;             % Number of FFT bands (0 to calculate from framelength)
+
+if nargin < 2
+    if verbose
+        disp('No audioconf given, using defaults.')
+    end
+    audioconf = defconf;
+else
+    fldnames = fieldnames(defconf);
+    for fl = 1:length(fldnames)
+        if ~isfield(audioconf, fldnames{fl})
+            if verbose
+                fprintf('Field %s missing, copying from defaults.\n', fldnames{fl})
+            end
+            audioconf.(f)=defconf.(f);
+        end
+    end
+end
+
+% Fetch the shorthand variables.
+featbands = audioconf.melbands;
+featchans = audioconf.featchannels;
+fs = audioconf.fs;
+fhigh = audioconf.maxf;
+flow = audioconf.minf;
+
+framelen = ceil(fs * audioconf.framelen_ms / 1000);
+frameshift = ceil(fs * audioconf.framestep_ms / 1000);
+% framelen = (fs * audioconf.framelen_ms / 1000);
+% frameshift = (fs * audioconf.framestep_ms / 1000);
+
+if audioconf.Nfft == 0
+    Nfft = 2^nextpow2(framelen);
+else
+    Nfft = audioconf.Nfft;
+end
+
+winfunc = str2func(audioconf.windowfunc);
+win = winfunc(framelen);
+
+if featbands == 0
+    melmode = false;
+else
+    melmode = true;
+end
+
+% Switch audio to columns.
+if size(sam, 1) < size(sam,2)
+    sam = sam';
+end
+
+samlen = size(sam, 1);
+samchans = size(sam, 2);
+
+if samchans ~= audioconf.channels
+    if verbose
+        fprintf('Warning: Audio has %i channels, config states %i.\n', samchans, audioconf.channels);
+    end
+end
+
+if melmode
+    melmat = mel_matrix(fs, featbands, Nfft, 1, fhigh, flow)';
+    if size(melmat, 1) ~= featbands
+        fprintf('Mel matrix has %i bands (config: %i).\n', size(melmat, 1), featbands);
+    end
+    if size(melmat, 2) ~= (Nfft/2 + 1)
+        fprintf('Mel matrix has %i FFT coeffs (expected: %i).\n', size(melmat, 2), Nfft/2 + 1);
+    end
+end
+
+% Truncate to full frames, get the number.
+numframes = floor((samlen-framelen+frameshift) / frameshift);
+sam = sam(1:(numframes*frameshift+framelen-frameshift), :);
+
+% DC removal - introduces a 1-unit filter delay, thus we discard the
+% first sample. Note that this behaviour has changed from earlier
+% versions of FE.
+if audioconf.dcremoval
+    samf = filter([1;-1], [1;-0.999], [zeros(1,samchans);sam]);
+    sam = samf(2:end, :);
+end
+samtrlen = size(sam, 1); % trimmed length
+
+% Pre-emphasis if nonzero. Can be done for the whole audio at once.
+if (audioconf.preemphasis > 0)
+    sam = [zeros(1, samchans); sam(2:samtrlen, :) - audioconf.preemphasis * sam(1:(end-1), :)];
+end
+
+if melmode
+    tmpfeats = zeros(featbands, numframes, samchans);
+else
+    tmpfeats = zeros(Nfft/2 + 1, numframes, samchans);
+end
+
+energies = zeros(numframes, samchans);
+frameaudio = zeros(framelen, numframes, samchans);
+frameFFT = zeros(Nfft/2+1, numframes, samchans);
+
+
+% Process channels one by one. Trying to perform these ops simultaneously
+% for all channels might be possible but tricky.
+for c = 1:samchans
+    
+    % starting sample numbers of each frame
+    ind1 = 1:frameshift:samtrlen-1-framelen+frameshift;
+    % linear 1-step vector (1...frame length)
+    ind2 = (1:framelen)';
+    
+    % Pick frame audio. The index matrix (framelen x numframes) consists 
+    % of four summed parts:
+    % 1) Constant column vectors, each denoting the frame's start sample.
+    % 2) Increasing sample index column vectors
+    % 3) Scalar jump to get into the correct channel in linear indexing
+    % 4) -1 because the first two indices are both one-based.
+    %
+    %       [start1  start2 ]   [ 1    1  ]
+    %  sam( [ ...     ...   ] + [...  ... ] + channel jump - 1) =
+    %       [start1  start2 ]   [frl  frl ]
+    %
+    %       [ start1+1   start2+1  ]
+    %  sam( [   ...        ...     ] + channel jump - 1)
+    %       [start1+frl start2+frl ]
+    %
+    % Thus we get an index matrix, where each frame column picks the 
+    % samples belonging to it. These samples are then fetched to 'fra'.
+    
+    fra = sam(ind1(ones(framelen,1),:) + ind2(:,ones(1,numframes)) + (c-1)*samtrlen - 1);
+    frameaudio(:,:,c) = fra;
+
+    % Calculate the energies.
+    energies(:,c) = sum(fra.^2,1)';            
+        
+    % Apply window function, take FFT.
+    fFFT = fft(win(:,ones(1,numframes)) .* fra, Nfft);
+    % Truncate and reset constant factor, but do not take abs yet.
+    fFFT(1,:) = 0;
+    fFFT = fFFT(1:Nfft/2+1,:);
+    
+    % Store the returned FFTs with phase.
+    frameFFT(:,:,c) = fFFT;
+    
+    if melmode
+        tmpfeats(:,:,c) = melmat * abs(fFFT);
+    else
+        tmpfeats(:,:,c) = abs(fFFT);
+    end
+end
+
+% Flatten the features if downmixing to 1 is defined.
+if featchans == 1 
+    if samchans > 1
+        feats = mean(tmpfeats, 3);
+    else
+        feats = tmpfeats;
+    end
+else
+    if samchans ~= featchans
+        fprintf('Requested %i feature channels for %i audio - not defined. Returning %i.\n', featchans, samchans, samchans)
+    end
+    feats = tmpfeats;
+end