基于MATLAB的语音信号滤波处理题目:基于MATLAB的语音信号滤波处理课程:数字信号处理学院:电气工程学院班级:学生:指导教师:二O一三年十二月目录CONTENTS摘要一、引言二、正文1.设计要求2.设计步骤3.设计内容4.简易GUI设计三、结论四、收获与心得五、附录一、引言随着Matlab仿真技术的推广,我们可以在计算机上对声音信号进行处理,甚至是模拟。
通过计算机作图,采样,我们可以更加直观的了解语音信号的性质,通过matlab编程,调用相关的函数,我们可以非常方便的对信号进行运算和处理。
二、正文2.1 设计要求在有噪音的环境中录制语音,并设计滤波器去除噪声。
2.2 设计步骤1.分析原始信号,画出原始信号频谱图及时频图,确定滤波器类型及相关指标;2.按照类型及指标要求设计出滤波器,画出滤波器幅度和相位响应,分析该滤波器是否符合要求;3.用所设计的滤波器对原始信号进行滤波处理,画出滤波后信号的频谱图及时频图;4.对滤波前的信号进行分析比对,评估所设计滤波器性能。
2.3 设计内容1.原始信号分析分析信号的谱图可知,噪音在1650HZ和3300HZ附近的能量较高,而人声的能量基本位于1000HZ以下。
因此,可以设计低通滤波器对信号进行去噪处理。
2.IIR滤波器设计用双线性变换法分别设计了巴特沃斯低通滤波器和椭圆低通滤波器和带阻滤波器:①巴特沃斯滤波器fp=800;fs=1300;rs=35;rp=0.5;程序代码如下:fp=800;fs=1300;rs=35;rp=0.5;Fs=44100;wp=2*Fs*tan(2*pi*fp/(2*Fs));ws=2*Fs*tan(2*pi*fs/(2*Fs));[n,wn]=buttord(wp,ws,rp,rs,'s');[b,a]=butter(n,wn,'s');[num,den]=bilinear(b,a,Fs);[h,w]=freqz(num,den,512,Fs);②椭圆低通滤波器fp=1300;fs=1600;rs=60;rp=0.5; 程序代码如下:fp=1300;fs=1600;rs=60;rp=0.5;Fs=44100;wp=2*Fs*tan(2*pi*fp/(2*Fs));ws=2*Fs*tan(2*pi*fs/(2*Fs)); [n,wn]=ellipord(wp,ws,rp,rs,'s');[b,a]=ellip(n,rp,rs,wn,'s');[num,den]=bilinear(b,a,Fs);[h,w]=freqz(num,den,512,Fs);③带阻滤波器fp1=800;fp2=2300;fs1=1300;fs2=1800;rs=30;rp=0.6 fp3=2800;fp4=4000;fs3=3200;fs4=3700;rs=30;rp=0.6 程序代码如下:fp1=800;fp2=2300;fs1=1300;fs2=1800;rs=30;rp=0.6;Fs=44100; fp=[fp1,fp2];fs=[fs1,fs2];wp=2*Fs*tan(2*pi*fp/(2*Fs));ws=2*Fs*tan(2*pi*fs/(2*Fs));[n,wn]=buttord(wp,ws,rp,rs,'s');[b,a]=butter(n,wn,'stop','s');[num,den]=bilinear(b,a,Fs);[h,w]=freqz(num,den,512,Fs);fp3=2800;fp4=4000;fs3=3200;fs4=3700;rs=30;rp=0.6;Fs=44100; fp1=[fp3,fp4];fs1=[fs3,fs4];wp1=2*Fs*tan(2*pi*fp1/(2*Fs));ws1=2*Fs*tan(2*pi*fs1/(2*Fs));[n1,wn1]=buttord(wp1,ws1,rp,rs,'s');[b1,a1]=butter(n1,wn1,'stop','s');[num1,den1]=bilinear(b1,a1,Fs);[h1,w1]=freqz(num1,den1,512,Fs);3.FIR滤波器①加hamming窗n=100;fp=1000;Fs=44100;b=fir1(n,fp/(Fs/2),Hamming(n+1)); [h,w]=freqz(b,1,512,Fs);②加hanning窗n=;fp=1000;Fs=44100;b=fir1(n,fp/(Fs/2),Hanning(n+1)); [h,w]=freqz(b,1,512,Fs);③加blackman窗n=100;fp=1000;Fs=44100;b=fir1(n,fp/(Fs/2),blackman(n+1));[h,w]=freqz(b,1,512,Fs);4.滤波前后比对①巴特沃斯低通滤波器滤波后②椭圆低通滤波器滤波后③带阻滤波器④加hamming窗⑤加hanning窗⑥加blackman窗2.4简易GUI界面设计为了便于操作和演示,设计了如下的简易GUI界面。
三、结论由以上谱图分析可知,经过滤波器滤波后,信号中的高频杂音明显被抑制,而人声成分大部分被保留,起到了预期的滤波作用。
对比所设计的两种滤波器,椭圆滤波器在过渡带相对较窄的情况下,能满足相对较高阻带衰减。
四、收获与心得本次设计大概进行了一周的时间,语音信号处理的是目前比较流行且十分有趣的,在编程实现的过程中还是遇到了很多困难。
在前期的准备工作中,查阅了大量资料,以完善我们的理论知识。
我们为了完成本次设计,我们通过查阅相关书籍以及matlab中的帮助,选用不同的matlab函数,尝试不同的参数。
经过接近一个礼拜的反复调试,最终基本的实现了设计任务。
虽然遇到了很多困难,但是我们在设计过程中都有收获很大。
本次设计将信号与系统课上学习的知识用于实践,让我们对对语音信号处理更深入的了解,也让我们加深了对滤波器相关内容的理解,同时也使得我们的Matlab能力有了很大的提高。
参考文献《应用matlab实现信号分析和处理》科学出版社附录1 巴特沃斯低通滤波器fp=800;fs=1300;rs=35;rp=0.5;Fs=44100; wp=2*Fs*tan(2*pi*fp/(2*Fs));ws=2*Fs*tan(2*pi*fs/(2*Fs));[n,wn]=buttord(wp,ws,rp,rs,'s');[b,a]=butter(n,wn,'s');[num,den]=bilinear(b,a,Fs);[h,w]=freqz(num,den,512,Fs);figure(1)%subplot(3,1,1)plot(w,abs(h));xlabel('频率/Hz');ylabel('幅值');title('巴特沃斯低通滤波器幅度特性');axis([0,5000,0,1.2])grid on;figure(2)%subplot(3,1,2)plot(w,20*log10(abs(h)));xlabel('频率/Hz');ylabel('幅值db');title('巴特沃斯低通滤波器幅度特性db'); axis([0,5000,-90,10]);grid on;figure(3)plot(w,180/pi*unwrap(angle(h))); xlabel('频率/Hz');ylabel('相位');title('巴特沃斯低通滤波器相位特性');axis([0,5000,-1000,10])grid on;[s1,Fs,bits]=wavread('D:\222.wav');x1=s1(:,1);sound(x1,Fs,bits);N1=length(x1);Y1=fft(x1,N1);f1=Fs*(0:N1-1)/N1;t1=(0:N1-1)/Fs;figure(4)plot(f1,abs(Y1))xlabel('频率/Hz');ylabel('幅度');title('原始信号频谱');grid on;axis([0 6000 0 400])y=filter(num,den,x1);sound(y,Fs,bits);N2=length(y);Y2=fft(y,N2);f2=Fs*(0:N2-1)/N2;t2=(0:N2-1)/Fs;figure(5)plot(f2,abs(Y2))xlabel('频率/Hz');ylabel('幅度');title('过滤后信号的频谱');grid on;axis([0 6000 0 100])2椭圆低通滤波器fp=1300;fs=1600;rs=60;rp=0.5;Fs=44100; wp=2*Fs*tan(2*pi*fp/(2*Fs));ws=2*Fs*tan(2*pi*fs/(2*Fs));[n,wn]=ellipord(wp,ws,rp,rs,'s');[b,a]=ellip(n,rp,rs,wn,'s');[num,den]=bilinear(b,a,Fs);[h,w]=freqz(num,den,512,Fs);figure(1)plot(w,abs(h));xlabel('频率/Hz');ylabel('幅值');title('椭圆低通滤波器幅度特性');axis([0,5000,0,1.2])grid on;figure(2)plot(w,20*log10(abs(h)));xlabel('频率/Hz');ylabel('幅值db');title('椭圆低通滤波器幅度特性db');axis([0,5000,-90,10]);grid on;figure(3)plot(w,180/pi*unwrap(angle(h)));xlabel('频率/Hz');ylabel('相位');title('椭圆低通滤波器相位特性');axis([0,5000,-1000,10])grid on;[s1,Fs,bits]=wavread('D:\222.wav');x1=s1(:,1);sound(x1,Fs,bits);N1=length(x1);Y1=fft(x1,N1); %对信号做N点FFT变换f1=Fs*(0:N1-1)/N1;t1=(0:N1-1)/Fs;figure(4)plot(f1,abs(Y1))xlabel('频率/Hz');ylabel('幅度');title('原始信号频谱');grid on;axis([0 6000 0 400])y=filter(num,den,x1);sound(y,Fs,bits);N2=length(y);Y2=fft(y,N2); %对信号做N点FFT变换f2=Fs*(0:N2-1)/N2;t2=(0:N2-1)/Fs;figure(5)plot(f2,abs(Y2))xlabel('频率/Hz');ylabel('幅度');title('过滤后信号的频谱');grid on;axis([0 6000 0 100])3.带阻滤波器fp1=800;fp2=2300;fs1=1300;fs2=1800;rs=30;rp=0.6;Fs=44100; fp=[fp1,fp2];fs=[fs1,fs2];wp=2*Fs*tan(2*pi*fp/(2*Fs));ws=2*Fs*tan(2*pi*fs/(2*Fs));%wap=2*tan(wp/2)/Ts[n,wn]=buttord(wp,ws,rp,rs,'s');[b,a]=butter(n,wn,'stop','s');[num,den]=bilinear(b,a,Fs);[h,w]=freqz(num,den,512,Fs);fp3=2800;fp4=4000;fs3=3200;fs4=3700;rs=30;rp=0.6;Fs=44100; fp1=[fp3,fp4];fs1=[fs3,fs4];wp1=2*Fs*tan(2*pi*fp1/(2*Fs));ws1=2*Fs*tan(2*pi*fs1/(2*Fs));%wap=2*tan(wp/2)/Ts[n1,wn1]=buttord(wp1,ws1,rp,rs,'s');[b1,a1]=butter(n1,wn1,'stop','s');[num1,den1]=bilinear(b1,a1,Fs);[h1,w1]=freqz(num1,den1,512,Fs);figure(1)plot(w,abs(h));xlabel('频率/Hz');ylabel('幅值');title('巴特沃斯带阻滤波器幅度特性');axis([0,5000,0,1.2])grid on;figure(2)plot(w,20*log10(abs(h)));xlabel('频率/Hz');ylabel('幅值db');title('巴特沃斯带阻滤波器幅度特性db');axis([0,5000,-90,10]);grid on;figure(3)plot(w,180/pi*unwrap(angle(h)));xlabel('频率/Hz');ylabel('相位');title('巴特沃斯带阻滤波器相位特性');axis([0,5000,-1000,10])grid on;figure(4)plot(w1,abs(h1));xlabel('频率/Hz');ylabel('幅值');title('巴特沃斯带阻滤波器幅度特性');axis([0,5000,0,1.2])grid on;figure(5)plot(w1,20*log10(abs(h1)));xlabel('频率/Hz');ylabel('幅值db');title('巴特沃斯带阻滤波器幅度特性db');axis([0,5000,-90,10]);grid on;figure(6)plot(w1,180/pi*unwrap(angle(h1)));xlabel('频率/Hz');ylabel('相位');title('巴特沃斯带阻滤波器相位特性');axis([0,5000,-1000,10])grid on;[s1,Fs,bits]=wavread('D:\222.wav'); x1=s1(:,1);sound(x1,Fs,bits);N1=length(x1);Y1=fft(x1,N1); %对信号做N点FFT变换f1=Fs*(0:N1-1)/N1;t1=(0:N1-1)/Fs;figure(7)plot(f1,abs(Y1))xlabel('频率/Hz');ylabel('幅度');title('原始信号频谱');grid on;axis([0 6000 0 400])y1=filter(num,den,x1);y=filter(num1,den1,y1);sound(y,Fs,bits);N2=length(y);Y2=fft(y,N2);f2=Fs*(0:N2-1)/N2;t2=(0:N2-1)/Fs;figure(8)plot(f2,abs(Y2))xlabel('频率/Hz');ylabel('幅度');title('过滤后信号的频谱');grid on;axis([0 6000 0 100])4.加hamming窗n=100;fp=1000;Fs=44100;b=fir1(n,fp/(Fs/2),Hamming(n+1)); [h,w]=freqz(b,1,512,Fs);a=num2str(a);[s1,Fs,bits]=wavread('D:\222.wav'); x1=s1(:,1);sound(x1,Fs,bits);N1=length(x1);Y1=fft(x1,N1);f1=Fs*(0:N1-1)/N1;y=fftfilt(b,x1);sound(y,Fs,bits);N2=length(y);Y2=fft(y,N2);f2=Fs*(0:N2-1)/N2;figure(4)subplot(2,1,1)plot(f1,abs(Y1))xlabel('频率/Hz');ylabel('幅度');title('原始信号频谱');grid on;axis([0 6000 0 400])subplot(2,1,2)plot(f2,abs(Y2));xlabel('频率/Hz');ylabel('幅度');title('过滤后信号的频谱');grid on;axis([0 6000 0 100])5.加hanning窗n=100;fp=1000;Fs=44100;b=fir1(n,fp/(Fs/2),Hanning(n+1)); [h,w]=freqz(b,1,512,Fs);a=num2str(a);[s1,Fs,bits]=wavread('D:\222.wav'); x1=s1(:,1);sound(x1,Fs,bits);N1=length(x1);Y1=fft(x1,N1);f1=Fs*(0:N1-1)/N1;y=fftfilt(b,x1);sound(y,Fs,bits);N2=length(y);Y2=fft(y,N2);f2=Fs*(0:N2-1)/N2;figure(4)subplot(2,1,1)plot(f1,abs(Y1))xlabel('频率/Hz');ylabel('幅度');title('原始信号频谱');grid on;axis([0 6000 0 400])subplot(2,1,2)plot(f2,abs(Y2));xlabel('频率/Hz');ylabel('幅度');title('过滤后信号的频谱');grid on;axis([0 6000 0 100])6.加blackman窗n=100;fp=1000;Fs=44100;b=fir1(n,fp/(Fs/2),blackman(n+1));[h,w]=freqz(b,1,512,Fs);a=num2str(a);[s1,Fs,bits]=wavread('D:\222.wav');x1=s1(:,1);sound(x1,Fs,bits);N1=length(x1);Y1=fft(x1,N1);f1=Fs*(0:N1-1)/N1;y=fftfilt(b,x1);sound(y,Fs,bits);N2=length(y);Y2=fft(y,N2);f2=Fs*(0:N2-1)/N2;figure(4)subplot(2,1,1)plot(f1,abs(Y1))xlabel('频率/Hz');ylabel('幅度');title('原始信号频谱');grid on;axis([0 6000 0 400])subplot(2,1,2)plot(f2,abs(Y2));xlabel('频率/Hz');ylabel('幅度');title('过滤后信号的频谱');grid on;axis([0 6000 0 100])7.巴特沃斯低通滤波时频分析clearFs=44100;[s1,Fs,bits]=wavread('D:\222.wav');x=s1(:,1);x1=x(40001:64576)l=length(x1);N=6;m=downsample(x1,N);%降抽样后 Fs=7350hzz=length(m);%4096[tfr, t, f] = tfrstft(m);figure(13)contour(t,(Fs/N)*(0:z/2-1)/z,abs(tfr(1:z/2,:)).^2); xlabel('时间t');ylabel('频率f');title('等高线图(滤波前)');grid on;figure(14)mesh(t,(Fs/N)*(0:z/2-1)/z,abs(tfr(1:z/2,:)).^2);xlabel('时间t');ylabel('频率f');zlabel('幅值A');title('三维图(滤波前)');grid on;fp=1300;fs=1600;rs=60;rp=0.5;Fs=44100;wp=2*Fs*tan(2*pi*fp/(2*Fs));ws=2*Fs*tan(2*pi*fs/(2*Fs)); [n,wn]=ellipord(wp,ws,rp,rs,'s');[b,a]=ellip(n,rp,rs,wn,'s');[num,den]=bilinear(b,a,Fs);[h,w]=freqz(num,den,512,Fs);y=filter(num,den,x);x1=y(40001:64576)l=length(x1);N=6;m=downsample(x1,N);%降抽样后 Fs=7350hzz=length(m);%4096[tfr,t,f] = tfrstft(m);figure(15)contour(t,(Fs/N)*(0:z/2-1)/z,abs(tfr(1:z/2,:)).^2); xlabel('时间t');ylabel('频率f');title('等高线图(滤波后)');grid on;figure(16)mesh(t,(Fs/N)*(0:z/2-1)/z,abs(tfr(1:z/2,:)).^2);xlabel('时间t');ylabel('频率f');zlabel('幅值A');title('三维图(滤波后)');grid on;。