加QQ：1046233911
一、加窗和分帧函数
%
enframe.m 加窗和分帧分帧函数，用在语音处理中的。相当于加了一个移动的窗
% 一般调用格式为
% y=enframe(x,framelength,step);
nx=length(x);
nwin=length(win);
if (nwin == 1)
len = win;
else
len = nwin;
end
if (nargin < 3)
inc = len;
end
nf = fix((nx-len+inc)/inc);
f=zeros(nf,len);
indf= inc*(0:(nf-1)).';
inds = (1:len);
f(:) = x(indf(:,ones(1,len))+inds(ones(nf,1),:));
if (nwin > 1)
w =
win(:)';
f = f .*
w(ones(nf,1),:);
end
二、语音信号端点检测程序
function [x1,x2]=vad(k,fs,i)
%% [X1,X2]=vad(k,fs)
%% 语音信号端点检测程序，k为语音信号，
%% fs为其采样频率，程序绘制出语音信号
%% 相关波形图并返回起止端点所对帧号
%close all;
% 幅度归一化到[-1,1]
k=double(k);
k=k./max(abs(k)); % k>0,k=1;k<=-1.
%------------------------------
%
此处录取部分代码
(3,1,3);
plot(t1,k1);
axis([0.4,0.5,min(k),max(k)]);
z=sprintf('(III) “�.wav”语音结束处放大波形图',i);
title(z);
xlabel('Time:s');
ylabel('Amplitude(normalized)');
%------------------------------
% 计算短时过零率
%------------------------------
FrameLen=240;
FrameInc=80;
FrameTemp1=enframe(k(1:end-1),FrameLen,FrameInc);
FrameTemp2=enframe(k(2:end),FrameLen,FrameInc);
signs=(FrameTemp1.*FrameTemp2)<0;
diffs=abs(FrameTemp1-FrameTemp2)>0.01;
zcr=sum(signs.*diffs,2);
disp('显示原始波形图……');
figure,subplot(3,1,1);
plot(t,k);
axis([0,(length(k)-1)/fs,min(k),max(k)]);
z=sprintf('(I) “�.wav”语音信号波形',i);
title(z);
xlabel('Time:s');
ylabel('Amplitude(normalized)');
disp('显示过短时零率……')
zcrInd=1:length(zcr);
subplot(3,1,2);
plot(zcrInd,zcr);
axis([0,length(zcr),0,max(zcr)]);
title('(II) 短时过零率');
xlabel('Frame');
ylabel('Zcr');
%------------------------------
% 计算短时能量
%------------------------------
amp=sum(abs(enframe(filter([1 -0.9375], 1, k), FrameLen,
FrameInc)), 2);
disp('显示短时能量……')
ampInd=1:length(amp);
subplot(3,1,3);
plot(ampInd,amp);
axis([0,length(amp),0,max(amp)]);
title('(III) 短时能量');
xlabel('Frame');
ylabel('Energy');
%------------------------------
% 设置门限
%------------------------------
disp('设置门限……');
ZcrLow=max([round(mean(zcr)*0.1),3]); %过零率低门限
ZcrHigh=max([round(max(zcr)*0.1),5]); %过零率高门限
AmpLow=min([min(amp)*10,mean(amp)*0.2,max(amp)*0.1]); %能量低门限
AmpHigh=max([min(amp)*10,mean(amp)*0.2,max(amp)*0.1]); %能量高门限
%------------------------------
% 端点检测
%------------------------------
MaxSilence=8; %最长语音间隙时间
MinAudio=16; %最短语音时间
Status=0; %状态：0静音段,1过渡段,2语音段,3结束段
HoldTime=0; %语音持续时间
SilenceTime=0; %语音间隙时间
disp('开始端点检测……');
for n=1:length(zcr)
switch
Status
case{0,1}
if amp(n)>AmpHigh ||
zcr(n)>ZcrHigh
x1=n-HoldTime;
Status=2;
HoldTime=HoldTime+1;
SilenceTime=0;
elseif amp(n)>AmpLow ||
zcr(n)>ZcrLow
Status=1;
HoldTime=HoldTime+1;
此处略去部分代码
HoldTime=HoldTime+1;
else
SilenceTime=SilenceTime+1;
if SilenceTime
HoldTime=HoldTime+1;
elseif (HoldTime-SilenceTime)
Status=0;
HoldTime=0;
SilenceTime=0;
else
Status=3;
end
end
case 3,
break;
end
if
Status==3
break;
end
end
HoldTime=HoldTime-SilenceTime;
x2=x1+HoldTime;
disp('显示端点……');
figure,subplot(3,1,1);
plot(k);
axis([1,length(k),min(k),max(k)]);
z=sprintf('(I) “�.wav”语音信号',i);
title(z);
xlabel('Sample');
ylabel('Speech');
line([x1*FrameInc,x1*FrameInc],[min(k),max(k)],'color','red');
line([x2*FrameInc,x2*FrameInc],[min(k),max(k)],'color','red');
subplot(3,1,2);
plot(zcr);
axis([1,length(zcr),0,max(zcr)]);
title('(II) 短时过零率');
xlabel('Frame');
ylabel('ZCR');
line([x1,x1],[0,max(zcr)],'Color','red');
line([x2,x2],[0,max(zcr)],'Color','red');
subplot(3,1,3);
plot(amp);
axis([1,length(amp),0,max(amp)]);
title('(III) 短时能量');
xlabel('Frame');
ylabel('Energy');
line([x1,x1],[0,max(amp)],'Color','red');
line([x2,x2],[0,max(amp)],'Color','red');
三、计算MFCC系数
function cc=mfcc(k)
%------------------------------
% cc=mfcc(k)计算语音k的MFCC系数
%------------------------------
% M为滤波器个数，N为一帧语音采样点数
M=24; N=256;
% 归一化mel滤波器组系数
bank=melbankm(M,N,8000,0,0.5,'m');
figure;
plot(linspace(0,N/2,129),bank);
title('Mel-Spaced Filterbank');
xlabel('Frequency [Hz]');
bank=full(bank);
bank=bank/max(bank(:));
% DCT系数,12*24
for i=1:12
j=0:23;
dctcoef(i,:)=cos((2*j+1)*i*pi/(2*24));
end
% 归一化倒谱提升窗口
w=1+6*sin(pi*[1:12]./12);
w=w/max(w);
% 预加重
AggrK=double(k);
AggrK=filter([1,-0.9375],1,AggrK);
% 分帧
FrameK=enframe(AggrK,N,80);
% 加窗
for i=1:size(FrameK,1)
FrameK(i,:)=(FrameK(i,:))'.*hamming(N);
end
FrameK=FrameK';
% 计算功率谱
S=(abs(fft(FrameK))).^2;
disp('显示功率谱……')
figure;
plot(S);
axis([1,size(S,1),0,2]);
title('Power Spectrum (M=24, N=256)');
xlabel('Frame');
ylabel('Frequency [Hz]');
colorbar;
% 将功率谱通过滤波器组
P=bank*S(1:129,:);
% 取对数后作离散余弦变换
D=dctcoef*log(P);
% 倒谱提升窗
for i=1:size(D,2)
m(i,:)=(D(:,i).*w')';
end
% 差分系数
dtm=zeros(size(m));
for i=3:size(m,1)-2
dtm(i,:)=-2*m(i-2,:)-m(i-1,:)+m(i+1,:)+2*m(i+2,:);
end
dtm=dtm/3;
%合并mfcc参数和一阶差分mfcc参数
cc=[m,dtm];
%去除首尾两帧，因为这两帧的一阶差分参数为0
cc=cc(3:size(m,1)-2,:);
四、计算MFCC时所需Mel尺度滤波器组
function [x,mn,mx]=melbankm(p,n,fs,fl,fh,w)
%MELBANKM determine matrix for a mel-spaced filterbank
[X,MN,MX]=(P,N,FS,FL,FH,W)
if nargin < 6
w='tz';
if nargin < 5
fh=0.5;
if nargin
< 4
fl=0;
end
end
end
f0=700/fs;
fn2=floor(n/2);
lr=log((f0+fh)/(f0+fl))/(p+1);
% convert to fft bin numbers with 0 for DC term
bl=n*((f0+fl)*exp([0 1 p p+1]*lr)-f0);
b2=ceil(bl(2));
此处略去部分代码
v=2*[0.5 ones(1,b2-1) 1-pf+fp pf-fp
ones(1,fn2-b3-1) 0.5];
mn=1;
mx=fn2+1;
else
b1=floor(bl(1))+1;
b4=min(fn2,ceil(bl(4)))-1;
pf=log((f0+(b1:b4)/n)/(f0+fl))/lr;
fp=floor(pf);
pm=pf-fp;
k2=b2-b1+1;
k3=b3-b1+1;
k4=b4-b1+1;
r=[fp(k2:k4) 1+fp(1:k3)];
c=[k2:k4 1:k3];
v=2*[1-pm(k2:k4) pm(1:k3)];
mn=b1+1;
mx=b4+1;
end
if any(w=='n')
v=1-cos(v*pi/2);
elseif any(w=='m')
v=1-0.92/1.08*cos(v*pi/2);
end
if nargout > 1
x=sparse(r,c,v);
else
x=sparse(r,c+mn-1,v,p,1+fn2);
end
五、DTW算法
function dist = dtw(test, ref)
global x y_min y_max
global t r
global D d
global m n
t = test;
r = ref;
n = size(t,1);
m = size(r,1);
d = zeros(m,1);
D = ones(m,1) * realmax;
D(1) = 0;
% 如果两个模板长度相差过多，匹配失败
if (2*m-n<3) | (2*n-m<2)
dist = realmax;
return
end
% 计算匹配区域
xa = round((2*m-n)/3);
xb = round((2*n-m)*2/3);
if xb>xa
%xb>xa, 按下面三个区域匹配
% 1 :xa
% xa+1:xb
% xb+1:N
for x = 1:xa
y_max = 2*x;
y_min = round(0.5*x);
warp
end
for x = (xa+1):xb
y_max = round(0.5*(x-n)+m);
y_min = round(0.5*x);
warp
end
for x = (xb+1):n
y_max = round(0.5*(x-n)+m);
y_min = round(2*(x-n)+m);
warp
end
elseif xa>xb
%xa>xb, 按下面三个区域匹配
% 0 :xb
% xb+1:xa
% xa+1:N
for x = 1:xb
y_max = 2*x;
y_min = round(0.5*x);
warp
end
for x = (xb+1):xa
y_max = 2*x;
y_min = round(2*(x-n)+m);
warp
end
for x = (xa+1):n
y_max = round(0.5*(x-n)+m);
y_min = round(2*(x-n)+m);
warp
end
elseif xa==xb
%xa=xb, 按下面两个区域匹配
% 0 :xa
% xa+1:N
for x = 1:xa
y_max = 2*x;
y_min = round(0.5*x);
warp
end
for x = (xa+1):n
y_max = round(0.5*(x-n)+m);
y_min = round(2*(x-n)+m);
warp
end
end
%返回匹配分数
dist = D(m);
function warp
global x y_min y_max
global t r
global D d
global m n
d = D;
for y = y_min:y_max
D1 = D(y);
if y>1
D2 = D(y-1);
else
D2 = realmax;
end
if y>2
D3 = D(y-2);
else
D3 = realmax;
end
d(y) =
sum((t(x,:)-r(y,:)).^2) + min([D1,D2,D3]);
end
D = d;
六、语音训练阶段程序
%% train.m
disp('正在生成训练参数……');
for i=0:3
fname=sprintf('train\\�.wav',i );
[k,fs]=wavread(fname);
[x1,x2]=vad(k,fs,i);
cc=mfcc(k);
cc=cc(x1-2:x2-2,:);
ref(i+1).mfcc=cc;
end
disp('正在存储模板库……');
save 'mfcc.mat' ref;
% close all;
clc;
七、基于DTW的语音识别程序
对test中所有语音进行识别
% test1.m
clear;
close all;
clc;
disp('正在导入参考模板参数...');
load mfcc.mat;
disp('正在计算测试模板的参数...')
for i=0:3
fname = sprintf('test\\�.wav',i);
[k,fs]=wavread(fname);
[x1,x2]=vad(k,fs,i);
cc=mfcc(k);
cc=cc(x1-2:x2-2,:);
test(i+1).mfcc=cc;
此处略去部分代码
dist(i,j) = dtw(test(i).mfcc, ref(j).mfcc);
end
end
disp('正在计算匹配结果...')
for i=1:4
[d,j] = min(dist(i,:)); % d为返回值，j为列数
if j==1
fprintf('测试模板 �.wav 的识别结果为:前进\n',i-1);end
if j==2
fprintf('测试模板 �.wav
的识别结果为:后退\n',i-1); end
if j==3
fprintf('测试模板 �.wav 的识别结果为:向左\n',i-1); end
if j==4
fprintf('测试模板 �.wav 的识别结果为:向右\n',i-1); end
end
close all;
对test文件夹中任一个语音进行识别
% test2.m
clear;
close all;
clc;
disp('正在导入参考模板参数...');
load mfcc.mat;
disp('正在计算测试模板的参数...')
i=1;
% i=2;
% i=3;
% i=4;
fname = sprintf('test\\�.wav',i-1);
[k,fs]=wavread(fname);
[x1,x2]=vad(k,fs,i-1);
cc=mfcc(k);
cc=cc(x1-2:x2-2,:);
test(i).mfcc=cc;
disp('正在进行模板匹配...');
dist = zeros(1,4);
for j=1:4;
dist(j)=dtw(test(i).mfcc,ref(j).mfcc);
end
disp('正在计算匹配结果...');
m=min(dist);
[row,column]=find(dist==m);
j=column;
if
j==1; fprintf('测试模板 �.wav 的识别结果为:前进\n',i-1); end
if
j==2; fprintf('测试模板 �.wav
的识别结果为:后退\n',i-1); end
if
j==3; fprintf('测试模板 �.wav
的识别结果为:向左\n',i-1); end
if
j==4; fprintf('测试模板 �.wav
的识别结果为:向右\n',i-1); end
% close all;
扫一扫
3703
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
