语音信号的采集和频谱分析

语音信号的采集和频谱分析：

[y,fs,bits]=wavread('voice'); %读取音频信息（双声道，16位，频率44100Hz）sound(y,fs,bits); %回放该音频

Y=fft(y,4096); %进行傅立叶变换

subplot(211);

plot(y);

title('声音信号的波形');

subplot(212)

plot(abs(Y));

title('声音信号的频谱');

窗函数设计低通滤波器：

fp=1000;

fc=1200;

as=100;

ap=1;

fs=22000;

wp=2*fp/fs;

wc=2*fc/fs;

N=ceil((as-7.95)/(14.36*(wc-wp)/2))+1;

beta=0.1102*(as-8.7);

window=Kaiser(N+1,beta);

b=fir1(N,wc,window);

freqz(b,1,512,fs);

结果：

滤波：

[y,fs,bits]=wavread('voice');

d=filter(b,a,y);

D=fft(d);

subplot(211)

plot(d);

title('滤波后的声音波形')

subplot(212)

plot(abs(D))

title('滤波后的声音频谱')

回放：

sound(d,fs,bits)

与滤波之前相比，噪音明显降低了许多。

过零率的计算要用下面的代码：

zcr = zeros(size(y,1)1);

delta= 0.02;

for i=1:size(y,1)

x=y(i,:);

for j=1;length(x)-1

if x(j)*x(j+1)<0 &abs(x(j)-x(j+1))>delta

zcr(i)=zcr(i)+1;

end

end

end

其中设置了门限delta=0.02。这是个经验值，可以进行细微的调整。在此条件下可以得到如图所示的过零率波形。与过零率曲线画在一起的是原始的语音信号波形，可以看到，语音信号音母部分的幅度比较低，但是其过零率的数值却很高，峰值将近50，而后面的韵母部分过零率则比较低，在20左右。

加矩形窗的短时能量函数：

a=wavread('F:\WO.wav');

subplot(6,1,1),plot(a);

N=32;

for i=2:6

h=linspace(1,1, (i-1)*N);%形成一个矩形窗，长度为N

En=conv(h,a.*a);%求卷积得其短时能量函数En

subplot(6,1,i),plot(En);

if(i==2) legend('N=32');

elseif(i==3) legend('N=64');

elseif(i==4) legend('N=128');

elseif(i==5) legend('N=256');

elseif(i==6) legend('N=512');

end

end

加hamming窗的短时能量函数：

把h=linspace(1,1, (i-1)*N);

改为h1=hamming((i-1)*N);

加矩形窗的短时平均幅度：

a=wavread('F:\WO.wav');

subplot(6,1,1),plot(a);

N=32;

for i=2:6

h=linspace(1,1,(i-1)*N);%形成一个矩形窗，长度为N

En=conv(h,abs(a));%求卷积得其短时能量函数En

subplot(6,1,i),plot(En);

if(i==2) legend('N=32');

elseif(i==3) legend('N=64');

elseif(i==4) legend('N=128');

elseif(i==5) legend('N=256');

elseif(i==6) legend('N=512');

end

end

短时过零率：

a=wavread('F:\WO.wav');

n=length(a);

N=320;

subplot(3,1,1),plot(a);

h=linspace(1,1,N);%形成一个矩形窗，长度为N En=conv(h,a.*a);%求卷积得其短时能量函数En subplot(3,1,2),plot(En);

for i=1:n-1

if a(i)>=0

b(i)= 1;

else

b(i) = -1;

end

if a(i+1)>=0

b(i+1)=1;

else

b(i+1)=-1;

end

w(i)=abs(b(i+1)-b(i));

end%求出每相邻两点符号的差值的绝对值

k=1;

j=0;

while (k+N-1)

Zm(k)=0;

for i=0:N-1;

Zm(k)=Zm(k)+w(k+i);

end

j=j+1;

k=k+160; %每次移动半个窗

end

for w=1:j

Q(w)=Zm(160*(w-1)+1)/640;%短时平均过零率end

subplot(3,1,3),plot(Q);

短时过零率分析：

[x,fs,nbit]=wavread(‘’);

x=x/max(abs(x));

len=256;

inc=128;

tmp1=enframe(x(1:end-1),256,128);

tmp2=enframe(x(2:end),256,128);

signs(tmp1.*tmp2)<0;

diffs=(tmp1-tmp2)>0.02;

zor=sum(signs.*diffs,2);

figure;

subplot(2,1,1),plot(x)

subplot(2,1,2),plot(zor)

基音检测：

waveFile='lcq01.wav';[y, fs, nbits]=wavread(waveFile);

index1=11000;frameSize=256;

index2=index1+frameSize-1;

frame=y(index1:index2); subplot(2,1,1); plot(y);

grid on title(waveFile); line(index1*[1 1], [-1 1], 'color', 'r');

line(index2*[1 1], [-1 1], 'color', 'r'); subplot(2,1,2); plot(frame, '.-');

grid on point=[7, 189];

%基本减谱法

clear;

%[xx,fs]=wavread('E:\mywhisper\shu.wav');

%[xx,fs]=wavread('E:\speech\x\w1xun_01.wav');

%[xx,fs]=wavread('E:\speech\耳语音切割\b\w1ba_5');

[xx,fs]=wavread('D:\yuan.wav');% 读取音频文件yuan.wav，并返回采样数据给变量xx及采样率Fs

[team,row]=size(xx);%将数组xx的行数赋给team,列数赋给row

if row==2

x=(xx(:,1)+xx(:,2))/2;

yy=x;

%如果语音信号xx为2列，即信号为双声道，则将其转换成单声道信号，即取两列的平均值赋给x，并将x的值赋给yy

else

x=xx;

yy=x;

%若语音信号xx为单声道，则将xx的值赋给x，并将x的值赋给yy

end

x=x-mean(x)+0.1*rand(length(x),1

N=length(x);%将语音信号长度赋给变量N

n=220;%对语音信号进行分帧，帧长为220

n1=160;%帧移为160

frame=floor((N-n)/(n-n1));%将分帧数赋给变量frame

%frame=floor(N/n);

for i=1:frame

y1=x((i-1)*(n-n1)+1:(i-1)*(n-n1)+n).*hamming(n);

%对每段分帧进行加窗处理

fy=fft(y1,n);

nen(i,:)=abs(fy).^2;% 将频域信号功率赋给矩阵变量nen

ang(i,:)=angle(fy);%将频域信号的相位角赋给矩阵变量ang

end

yuzhi=sum(sum(nen(2:5,:)))/(4*n);

for i=1:frame

nen(i,:)=nen(i,:)-yuzhi;

nen(i,find(nen(i,:)<0))=0;

%chuli=nen(i,1:n/2);

%chuli=chuli-yuzhi;

%chuli(find(chuli<0))=0;

%nen(i,:)=[chuli,fliplr(chuli)];

% nen(i,:)=filter(1,[0.5 0.5],nen(i,:));

% nen(i,find(nen(i,:)<0))=0;

end

for i=1:frame

nen(i,:)=sqrt(nen(i,:));%将纯语音功率谱开根，得到频域值

jie=nen(i,:).*exp(j*ang(i,:));

out(i,:)=real(ifft(jie))/hamming(n)'; %对纯语音频谱进行逆傅里叶变换，并取其实部，并进行去窗处理

end

zong=out(1,:)';%将第一帧中未重叠部分记入数组zong

jiewei=n;

for i=2:frame

zong(jiewei-n1+1:jiewei)

=(zong(jiewei-n1+1:jiewei)+out(i,1:n1)')/2;

jiewei=jiewei+n-n1;% 使指针jiewei依次指向下一帧的帧尾

zong=[zong;out(i,n1+1:end)'];将从第二帧开始的每一帧中未重叠部分记入数组zong

end

%zong=out(1,:)';

%for i=2:frame

%zong=[zong;out(i,:)'];

%end

%

%for i=1:frame

% zong=[zong,nen(i,:)'];

% zong(i*(n-n1)+1:(i-1)*(n-n1)+n)

=zong(i*(n-n1)+1:(i-1)*(n-n1)+n)/2;

figure(1); %创建图1

subplot(211); %把图形窗口分成2*1个小窗口，取第1个小窗口

plot(x);%以数组x绘图基本二维曲线

axis([1,(n-n1)*frame+n,min(x),max(x)]);% 对当前二维图形对象的X轴和Y轴进行标定，x 轴的范围为1到最后一个分帧结尾，y轴的范围为带噪语音时域最小值到最大值

subplot(212); %把图形窗口分成2*1个小窗口，取第2个小窗口

specgram(x,fs,1024,n,n1);%画出语音的语谱图

figure(2); %创建图2

subplot(211); %把图形窗口分成2*1个小窗口，取第1个小窗口

plot(zong); %以数组zong绘图基本二维曲线

axis([1,(n-n1)*frame+n,min(zong),max(zong)]); % 对当前二维图形对象的X轴和Y轴进行标定，x轴的范围为1到最后一个分帧结尾，y轴的范围为纯语音时域最小值到最大值subplot(212); %把图形窗口分成2*1个小窗口，取第2个小窗口

specgram(zong,fs,1024,n,n1); %画出语音的语谱图

wavplay(x,fs);%播放单声道带噪语音音频

wavplay(zong,fs);%播放单声道纯净语音音频