matlab提取语音信号基频检测

function [tpitchFRpath] = pitchestimate（xfsparam）
% PITCHESTIMATE High-accuracy pitch estimation
%   [TP] = PITCHESTIMATE（XFs） computes a pitch estimate P at times T
%   for the signal X sampled at Fs Hz using default parameters.
%
%   [TP] = PITCHESTIMATE（XFsPARAM） uses the parameter values specified
%   in the parameter structure PARAM which must have the following fields
%   （units in parentheses default values in brackets）
%
%      minpitch           minimum detected pitch （Hz） [75]
%      maxpitch           maximum detected pitch （Hz） [400]
%      timestep           window hop （s） [0.01]
%      voicingthreshold   signal level required for voicing detection [0.4]
%      silencethreshold   signal level required for silence detection [0.05]
%      maxcandidates      number of pitch candidates to consider per frame [4]
%      octavecost         penalty for pitch halving [0.04]
%      voicedunvoicedcost penalty for switching voiced/unvoiced between frames [0.2]
%      octavejumpcost     penalty for switching octaves between frames [0.2]
%
%   The pitch estimation algorithm is based on the Praat pitch estimator described
%   in [1].
%
%   References:
%   [1]  P. Boersma “Accurate Short-Term Analysis of the Fundamental Frequency and
%        the Harmonics-to-Noise Ratio of Sampled Sounds“ IFA Proceedings 17 1993
%

% check number of input arguments
error（nargchk（23nargin））;

% setup default parameter structure if none specified
if nargin<3
  % define default algorithm parameters see [1] for details
  param.minpitch           = 75;   % Hz
  param.maxpitch           = 400;  % Hz
  param.timestep           = 0.02; % seconds
  param.voicingthreshold   = 0.4;
  param.silencethreshold   = 0.05;
  param.maxcandidates      = 4;
  param.octavecost         = 0.04;
  param.voicedunvoicedcost = 0.2;
  param.octavejumpcost     = 0.2;
end;

% compute size skip overlap and nfft in samples
winsize = round（fs/（param.minpitch/3））; % use 3 for pitch detection use 6 for HNR
winskip = round（param.timestep * fs）;
winolap = winsize - winskip;
winnfft = 2^nextpow2（winsize * 1.5）;

% define window
win = hanning（winsize）;

% define vector to demean local segments
demean = ones（winnfft1）; demean（1）=0;

% compute auto-correlation of window and its reciprocal
rwin = real（ifft（abs（fft（winwinnfft））.^2））;
rwinrecip = 1./rwin;

% step 1. determine windowed frames of speech
xbuf = buffer（xwinsizewinolap‘nodelay‘）;
r = diag（sparse（win）） * xbuf;

% step 2. compute auto-correlation of all segments
r = real（ifft（（diag（sparse（demean））*abs（fft（rwinnfft1）））.^2[]1））;

% step 3. divide by auto-correlation of the window
r = diag（sparse（rwinrecip）） * r;

% flip matrix top/bottom halves
r = fftshift（r1）;

% define time lag and time lags within our pitch range
tau = （-winnfft/2:winnfft/2-1）‘/fs;
mintau = 1/param.maxpitch;
maxtau = 1/param.minpitch;
goodtau = tau>=mintau & tau<=maxtau;