数字信号课程设计 语音信号的处理与滤波

广西科技大学（筹）数字信号处理课程设计设计题目：语音信号的处理与滤波设计过程：1.语音信号的采集利用电脑上的声卡和WINDOWS操作系统可以进行数字信号的采集。可以直接用电脑自带的录音工具进行录音。也可以运用MATLAB中的wavrecord()、wavwrite()、wavplay()、wavread()等函数进行语音的录制、播放、存储、读取。具体采集步骤如下：用MATLAB工具编程录制一段时间为5秒的内容（例如“老师，新年快乐”）的语音，存放在C:\Users\PKL\Desktoptft.wav,并把存放录音文件的文件夹导入MATLAB的work目录中。录制程序如下：程序：fs=16000;%设计抽样频率channel=1;%设计录音通道t=5;fprintf('按任意键开始%d秒录音:',t);pause;%暂停命令fprintf('正在录音');y=wavrecord(t*fs,fs,channel,'double');%录制语音fprintf('录音结束\n');wavwrite(y,fs,'C:\Users\PKL\Desktop\tft.wav');fprintf('按任意键回放语音');pause;wavplay(y,fs);2.语音信号的时频域分析利用MATLAB中的“wavread”命令来读入已经录制好的语音信号，将它赋值给任意向量，再对其进行采样。然后可用MATLAB的绘图函数对时域和频域进行绘图分析。运用下面程序进行采样和绘制时频域图，时域图为附件中的图——1，频域图为附件中的图——2。程序：[y,fs,nbits]=wavread('tft.wav');%加载语音信号wavplay(y,fs);%回放语音信号N=length(y);%语音信号长度Y=fft(y,N);%用傅里叶变换求频域subplot(2,1,1);plot(y);title('图——1原始信号时域波形图');subplot(2,1,2);plot(abs(Y));title('图——2原始信号频谱图');3.采用窗函数法和双线性变换法设计滤波器所需设计的滤波器的性能指标为：（1）低通滤波器性能指标：fb=1000Hz,fc=1200Hz,As=100dB,Ap=1dB。（2）高通滤波器性能指标：fb=5000Hz,fc=4800Hz,As=100dB,Ap=1dB。（3）带通滤波器性能指标：fp1=1200Hz,fp2=3000Hz,fsc1=1000Hz,fsc2=3200Hz,As=100dB,Ap=1dB。(一）用双线性变换法设计滤波器用双线性变换法设计滤波器就是先按低通滤波器的性能指标设计模拟低通滤波器，再转换成相应的过渡模拟滤波器。然后采用双线性变换法将相应类型的过渡模拟滤波器转换成所需类型的数字滤波器。这样可以避免产生频谱混叠现象。在双线性变换法中可用butter,cheby1和ellip设计相应IIR滤波器。设计程序如下：低通滤波器设计程序：fp=1000;ft=5500;%取抽样频率fs=1200;wp=2*pi*fp/ft;%通带频率ws=2*pi*fs/ft;%截止频率Omgp=2*ft*tan(wp/2);Omgs=2*ft*tan(ws/2);[n,Omgc]=buttord(Omgp,Omgs,1,100,'s');[z,p,k]=buttap(n);B=k*real(poly(z));A=real(poly(p));[b,a]=lp2lp(B,A,Omgc);[ba,aa]=bilinear(b,a,ft);[Ha,w]=freqz(ba,aa);plot(w/pi,abs(Ha));title('图——3IIR低通滤波器');附件图——3为IIR低通滤波器频谱图。高通滤波器设计程序：fp=5000;ft=25000;%取抽样频率fs=4800;Rp=1;As=100;wp=2*pi*fp/ft;%通带频率ws=2*pi*fs/ft;%截止频率Omgp=2*ft*tan(wp/2);Omgs=2*ft*tan(ws/2);[n,Omgc]=ellipord(Omgp,Omgs,Rp,As,'s');%计算阶数n和截止频率[z,p,k]=ellipap(n,Rp,As);B=k*real(poly(z));A=real(poly(p));[b,a]=lp2hp(B,A,Omgc);[ba,aa]=bilinear(b,a,ft);[Ha,w]=freqz(ba,aa);plot(w/pi,abs(Ha));title('图——4IIR高通滤波器');附件图——4为IIR高通滤波器频谱图。带通滤波器设计程序：fp1=1200;fp2=3000;fs1=1000;fs2=3200;ft=10000;As=100;Rp=1;wp1=2*pi*fp1/ft;wp2=2*pi*fp2/ft;ws1=2*pi*fs1/ft;ws2=2*pi*fs2/ft;Omgp1=2*ft*tan(wp1/2);Omgp2=2*ft*tan(wp2/2);Omgp=[Omgp1,Omgp2];Omgs1=2*ft*tan(ws1/2);Omgs2=2*ft*tan(ws2/2);Omgs=[Omgs1,Omgs2];bw=Omgp2-Omgp1;w0=sqrt(Omgp1*Omgp2);%通带宽和中心频率[n,Omgn]=cheb1ord(Omgp,Omgs,Rp,As,'s');[z,p,k]=cheb1ap(n,Rp);B=k*real(poly(z));A=real(poly(p));[b,a]=lp2bp(B,A,w0,bw);[ba,aa]=bilinear(b,a,ft);[Ha,w]=freqz(ba,aa);plot(w/pi,abs(Ha));title('图——5IIR带通滤波器');附件图——5为IIR高通滤波器频谱图。(二）用窗函数法设计滤波器窗函数设计方法主要有矩形窗、巴特列特窗、汉宁窗、海明窗、布拉克曼窗凯泽窗等。选择合适窗设计FIR滤波器，这里按照设置指标选择凯泽窗和布莱克曼窗进行设计。FIR低通滤波器设计程序：ft=5500;fpts=[10001200];dev=[0.010.05];%通带与阻带波动mag=[10];[n21,wn21,beta,ftype]=kaiserord(fpts,mag,dev,ft);b21=fir1(n21,wn21,Kaiser(n21+1,beta));[h,w]=freqz(b21,1);plot(w/pi,abs(h));title('图——6FIR低通滤波器');频率响应图为：图——6FIR高通滤波器设计程序：ft=25000;fpts=[48005000];dev=[0.010.05];mag=[10];[n22,wn22,beta,ftype]=kaiserord(fpts,mag,dev,ft);b22=fir1(n22,wn22,'high',Kaiser(n22+1,beta));[h,w]=freqz(b22,1);plot(w/pi,abs(h));title('图——7FI高通滤波器');频率响应图为：图——7FIR带通滤波器的设计程序：fp1=1200;fp2=3000;fs1=1000;fs2=3200;ft=10000;ws1=fs1/(ft/2)*pi;ws2=fs2/(ft/2)*pi;wp1=fp1/(ft/2)*pi;wp2=fp2/(ft/2)*pi;wp=[wp1,wp2];ws=[ws1,ws2];deltaw=wp1-ws1;N0=ceil(11*pi/deltaw);N=N0+mod(N0+1,2);windows=(blackman(N))';wc1=(ws1+wp1)/2;wc2=(ws2+wp2)/2;hd=ideallp(wc2,N)-ideallp(wc1,N);b=hd.*windows;[h,w]=freqz(b,1);plot(w/pi,abs(h));title('图——8FIR带通滤波器');频率响应图为：图——84.对“tft.wav”语音信号进行滤波分别采用设计好的低通，高通，带通滤波器对“tft.wav”语音信号进行滤波，用MATLAB画出滤波后的时频域图，并回放录播后的语音信号。（1）FIR低通滤波程序：ft=5500;fpts=[10001200];dev=[0.010.05];mag=[10];[n21,wn21,beta,ftype]=kaiserord(fpts,mag,dev,ft);b21=fir1(n21,wn21,Kaiser(n21+1,beta));[h,w]=freqz(b21,1);[y,fs,nbits]=wavread(‘tft.wav');%加载语音信号f=filter(b21,1,y);wavplay(f);f1=fft(f);subplot(2,1,1)plot(f);title('低通滤波后的波形');subplot(2,1,2)plot(abs(f1));title('低通滤波后的频谱图')运行结果图为附件：图——9（2）FIR高通滤波程序：ft=25000;fpts=[48005000];dev=[0.010.05];mag=[10];[n22,wn22,beta,ftype]=kaiserord(fpts,mag,dev,ft);b22=fir1(n22,wn22,'high',Kaiser(n22+1,beta));[h,w]=freqz(b22,1);[y,fs,nbits]=wavread('tft.wav');%加载语音信号f=filter(b22,1,y);wavplay(f);f1=fft(f);subplot(2,1,1)plot(f);title('高通滤波后的波形');subplot(2,1,2)plot(abs(f1));title('高通滤波后的频谱图')运行结果图为附件：图——10（3）FIR带通滤波程序：fp1=1200;fp2=3000;fs1=1000;fs2=3200;ft=10000;ws1=fs1/(ft/2)*pi;ws2=fs2/(ft/2)*pi;wp1=fp1/(ft/2)*pi;wp2=fp2/(ft/2)*pi;wp=[wp1,wp2];ws=[ws1,ws2];deltaw=wp1-ws1;N0=ceil(11*pi/deltaw);N=N0+mod(N0+1,2);windows=(blackman(N))';wc1=(ws1+wp1)/2;wc2=(ws2+wp2)/2;hd=ideallp(wc2,N)-ideallp(wc1,N);b=hd.*windows;[h,w]=freqz(b,1);[y,fs,nbits]=wavread('tft.wav');%加载语音信号f=filter(b,1,y);wavplay(f);f1=fft(f);subplot(2,1,1)plot(f);title('带通滤波后的波形');subplot(2,1,2)plot(abs(f1));title('带通滤波后的频谱图')运行结果图为附件：图——11（4）IIR低通滤波程序：fp=1000;ft=5500;%取抽样频率fs=1200;wp=2*pi*fp/ft;%通带频率ws=2*pi*fs/ft;%截止频率Omgp=2*ft*tan(wp/2);Omgs=2*ft*tan(ws/2);[n,Omgc]=buttord(Omgp,Omgs,1,100,'s');[z,p,k]=buttap(n);B=k*real(poly(z));A=real(poly(p));[b,a]=lp2lp(B,A,Omgc);[ba,aa]=bilinear(b,a,ft);[Ha,w]=freqz(ba,aa);[y,fs,nbits]=wavread('tft.wav');%加载语音信号f=filte