function [x, xf]=mfcc_model(seg, N, Fs)
% Do FFT of audio frame seg, and map to an N-point perceptual model from
% 0Hz to Fs/2 Hz
%
% N is typically 40, Fs is typically 8000
b_low=0;                %lower limit of the bark span
b_top=f2bark(Fs/2);     %upper limit of the bark span
bdiv=(b_top-b_low)/N;   %bark resolution

%Define an array of centre frequencies
xb=b_low+bdiv/2:bdiv:b_top-bdiv/2;
%Convert this to Hz frequencies
xf=bark2f(xb);

%Do the FFT...
S=fft(seg);
S=abs(S);
S=S(1:length(S)/2);
F=[1:length(S)]*(Fs/2)/length(S);

x=zeros(1,N);
for xi=1:N
    bark=xb(xi);
    freq=bark2f(bark);
    
    %compute spreading function
    spread_fn=spread_hz(F,freq);
    
    %compute summed influence
    x_sum=sum(spread_fn.*S')*equal_loudness_preemph(freq);
   
    %intensity loudness power law
    x(xi)=(x_sum)^0.33;
end
