examp_formants.m

linear time－frequency toolbox

%EXAMP_FORMANTS  Find the formants in speech
%
%   This script shows a simple method for finding formants in a speech
%   signal. Formants are peaks in a spectrum or spectogram, which can be
%   linked to resonant frequencies of the related physical system. They are
%   most commonly used in phonetics, searching for the formants of the human
%   vocal tract. This is for example used for speaker recognition.
%
%   REFERENCES:
%     R. Carmona, W. Hwang, and B. Torrésani. Multiridge detection and
%     time-frequency reconstruction. IEEE Trans. Signal Process., 47:480-492,
%     1999.
%     
%     G. Fant. Acoustic Theory of Speech Production. Mouton. The Hague., 1960.
%     
%     K. N. Stevens. Acoustic Phonetics. MIT Press. Cambridge, 1999.

% Authors:  Peter Balazs
%           Peter Soendergaard 

disp('Type "help examp_formants" to see a description of how this example works.');

a=64; % hop size
M=1024; % number of frequency bins

SR = 16000; % sampling rate

% Load the file
greasy=greasy();

% Set the filename
filename='greasy';

winname = 'pgauss'; % name of Gabor window
gf=tfr_create('dgt',a,M,winname);

c=tfr_ar(gf,greasy);

N=size(c,2);
F=4; 
% For voval tract formants four are often used for speaker recognition.

% This should be vectorized.
formants=zeros(F,N);
for n=1:N;
  [formants(:,n),Y] = findmaxima(abs(c(1:M/2,n)),F);
end;

figure(1);
specto = 20*log10(abs(c(1:M/2,:))); % log plot, more related to perception
maxlV = max(max(specto));
minlV = maxlV - 50; % 50 db range 
imagesc(specto,[minlV,maxlV]);
title(sprintf('Log. Amplitude STFT of %s \n Parameters: %s window, Nwin = %g, a = %g, M = %g',filename,winname,N,a,M));
set(gca,'Ydir','normal');

ylabel('kHz');
yt = get(gca,'ytick');
yt = yt*(SR/M)/1000; % rescaling of y axis for physical units 
 
ytext = sprintf('%-2.1f',yt(1));
for ii=2:length(yt)
     ytext = sprintf('%s | %-2.1f',ytext,yt(ii));
end
set(gca,'yticklabel',ytext);
xlabel('s');
xt = get(gca,'xtick');
xt = xt*a/SR; % rescaling of x axis for physical units
xtext = sprintf('%-2.1f',xt(1));
for ii=2:length(xt)
     xtext = sprintf('%s | %-2.2f',xtext,xt(ii));
end
set(gca,'xticklabel',xtext);

hold on;
formants = formants.';
plot(formants(:,1),'k.:','Linewidth',2);
plot(formants(:,2),'g.:','Linewidth',2);
plot(formants(:,3),'r.:','Linewidth',2);
plot(formants(:,4),'b.:','Linewidth',2);
legend('formant1','formant2','formant3','formant4');
hold off;

tfr_clear(gf);