语音信号处理.doc

总程序及图形显示： (1).原始语音信号的读入fs=22050;%语音信号采样频率为22050x1=wavread(c:rd.wav);%读取语音信号的数据，赋给变量x1sound(x1,22050);%播放语音信号figure(1)plot(x1);%做原始语音信号的时域图形title(原始语音信号);xlabel(时间/T);ylabel(幅值/N); (2).巴特沃斯低通滤波器的设计fs=22050;x1=wavread(c:rd.wav);t=0:1/22050:(length(x1)-1)/22050;%将所加噪声信号的点数调整到与原始信号相同Au=0.5;d=Au*cos(2*pi*8000*t);%噪声为8KHz的余弦信号x2=x1+d;%信号的叠加wp=0.25*pi;ws=0.3*pi;%通带和阻带的截止频率分别为0.25pi和0.3piRp=1;Rs=15;%通带和阻带区的波纹系数分别为1和15Fs=82050;Ts=1/Fs;wp1=2/Ts*tan(wp/2);%将模拟指标转换成数字指标ws1=2/Ts*tan(ws/2);N,Wn=buttord(wp1,ws1,Rp,Rs,s);%选择滤波器的最小阶数Z,P,K=buttap(N);%创建butterworth模拟滤波器Bap,Aap=zp2tf(Z,P,K);b,a=lp2lp(Bap,Aap,Wn);bz,az=bilinear(b,a,Fs);%用双线性变换法实现模拟滤波器到数字滤波器H,W=freqz(bz,az);%绘制频率响应曲线figure(2)plot(W*Fs/(2*pi),abs(H)xlabel(频率/Hz)ylabel(频率响应幅度/N)title(Butterworth) (3).滤波前后时域波形的比较f1=filter(bz,az,x2);figure(3)subplot(2,1,1)plot(t,x2);%画出滤波前的时域图title(滤波前的时域波形);subplot(2,1,2)plot(t,f1);%画出滤波后的时域图title(滤波后的时域波形);sound(f1,22050);%播放滤波后的语音信号F0=fft(f1,1024);f=fs*(0:511)/1024; (4).滤波前后频谱图的比较figure(4)y2=fft(x2,1024);subplot(2,1,1);plot(f,abs(y2(1:512);%画出滤波前的频谱图title(滤波前的频谱)xlabel(频率/Hz);ylabel(幅值/N);subplot(2,1,2)F1=plot(f,abs(F0(1:512);%画出滤波后的频谱图title(滤波后的频谱)xlabel(频率/Hz);ylabel(幅值/N);(5).滤波前后还原波形的比较figure(5)y3=ifft(y2);%将频谱进行傅立叶反变换subplot(2,1,1)plot(y3);%画出滤波前的还原波形title(滤波前的还原波形)xlabel(时间/T);ylabel(幅值/N);subplot(2,1,2)y4=ifft(F0);%将频谱进行傅立叶反变换plot(y4);%画出滤波后的还原波形title(滤波后的还原波形)xlabel(时间/T);ylabel(幅值/N); (6).去点叠加后的波形比较fs=1000;%采样频率x1 = wavread(c:rd.wav);%读入语音信号y=fft(x1,512);Q2个点的fft变换mag=abs(y);%求幅值ang=angle(y);%求相位f=(0:length(y)-1)*fs/length(y);%做相应的频率转换figure(6)xifft=ifft(y,512);magx=real(xifft);y(1:10)=0;%去掉10个点做逆变换xifft1=ifft(y,512);magx1=abs(xifft1);ti1=(0:length(xifft1)-1)/fs;%做相应的频率转换subplot(321)stem(ti1,magx1,.);title(去掉10个点叠加后的波形);axis(0,0.6,0,1);grid;y(1:50)=0;% %去掉50个点做逆变换xifft2=ifft(y,512);magx2=abs(xifft2);ti2=(0:length(xifft2)-1)/fs;%做相应的频率转换subplot(322)stem(ti2,magx2,.);title(去掉50个点叠加后的波形);axis(0,0.6,0,1);grid;y(1:100)=0;% %去掉100个点做逆变换xifft3=ifft(y,512);magx3=abs(xifft3);ti3=(0:length(xifft3)-1)/fs;%做相应的频率转换subplot(323)stem(ti3,magx3,.);title(去点100个点叠加后的波形);axis(0,0.6,0,1);grid;y(1:200)=0;% %去掉200个点做逆变换xifft4=ifft(y,512);magx4=abs(xifft4);ti4=(0:length(xifft4)-1)/fs;%做相应的频率转换subplot(324)stem(ti4,magx4,.);title(去掉200个点叠加后的波形);axis(0,0.6,0,1);grid;y(1:400)=0;% %去掉400个点做逆变换xifft5=ifft(y,512);magx5=abs(xifft5);ti5=(0:length(xifft5)-1)/fs;%做相应的频率转换subplot(325)stem(ti5,magx5,.);%做相应的频率转换title(去掉400个点叠加后的波形);axis(0,0.6,0,1);grid;y(1:500)=0;% %去掉500个点做逆变换xifft6=ifft(y,512);magx6=abs(xifft6);ti6=(0:length(xifft6)-1)/fs;%做相应的频率转换subplot(326)stem(ti6,magx6,.);title(去掉500个点叠加后的波形);axis(0,0.6,0,1);grid; (7).用三维图显示原始信号与部分频率分量叠加后的比较figure(7)clear allfc=10000; %采样频率x1=wavread(c:rd.wav);t=0:1/22050:(length(x1)-1)/22050;%相应的频率转换Au=0.5;d=Au*cos(2*pi*8000*t);%噪声为1KHz的余弦信号x2=x1+d;%信号的叠加y=fft(x1,512);Q2点的傅里叶变换mag=abs(y);%求幅值ang=angle(y);%求相位f=(0:length(y)-1)*fc/length(y);%相应的频率转换subplot(2,1,1)plot3(f,mag,ang);title(原始信号的三维图)xlabel(频率)ylabel(幅值)zlabel(相位)y1=fft(x2,512);Q2点的傅里叶变换mag=abs(y1);%求幅值ang=angle(y1);%求相位f1=(0:length(y1)-1)*fc/length(y1);%相应的频率转换subplot(2,1,2)plot3(f1,mag,ang);title(部分分量叠加的三维图)xlabel(频率)ylabel(幅值)zlabel(相位) (8).用三维投影显示原始信号与部分频率分量叠加后的幅度及相位的比较figure(8)fc=10000; %采样频率x1 = wavread(c:rd.wav);t=0:1/22050:(length(x1)-1)/22050;%相应的频率转换Au=0.5;d=Au*cos(2*pi*8000*t);%噪声为1KHz的余弦信号x2=x1+d;%信号的叠加y=fft(x1,512);Q2点的傅里叶变换mag=abs(y);%求幅值ang=angle(y);f=(0:length(y)-1)*fc/length(y);%相应的频率转换xifft=ifft(y,512);magx=real(xifft);axis(0,0.6,0,1);grid;a=;b=;c=;subplot(2,1,1)stem3(a,b,c);%建立三维坐标hold onstem3(f,ones(size(f)+99,mag,.);hold onstem3(f,ones(size(f),ang,.);hold ontitle(原始信号的投影)y1=fft(x2,512);Q2点的傅里叶变换mag=abs(y1);%求幅值ang=angle(y1);%求相位f1=(0:length(y1)-1)*fc/length(y1);%相应的频率转换xi1fft=iff