数字信号处理综合设计

1、数字信号处理综合设计1.语音信号的频谱分析要求首先画出语音信号的时域波形;然后对语音信号进行频谱分析,在MA TLAB中,可以利用函数fft对信号进行快速付立叶变换,得到信号的频谱特性;从而加深对频谱特性的理解。fs=20000; %语音信号采样频率为20000x1=wavread('d:lianxi.wav',20000; %读取语音信号的数据,赋给变量x1sound(x1,20000; %播放语音信号y1=fft(x1,1024; %对信号做1024点FFT变换f=fs*(0:511/1024;figure(1plot(x1 %做原始语音信号的时域图形title('

2、原始语音信号'xlabel('time n'ylabel('fuzhi n'figure(2freqz(x1 %绘制原始语音信号的频率响应图title('频率响应图'figure(3subplot(2,1,1;plot(abs(y1(1:512 %做原始语音信号的FFT频谱图title('原始语音信号FFT频谱'subplot(2,1,2;plot(f,abs(y1(1:512;title('原始语音信号频谱'xlabel('Hz'ylabel('fuzhi' 设计数字滤波器

3、和画出频率响应根据语音信号的特点给出有关滤波器的性能指标:1低通滤波器性能指标,fp=1000Hz, fc=1200 Hz,As=100dB,Ap=1dB;2高通滤波器性能指标,fc=2800 Hz,fp=3000 Hz As=100dB,Ap=1dB;3带通滤波器性能指标,fp1=1200 Hz,fp2=3000 Hz,fc1=1000 Hz, fc2=3200 Hz,As=100dB,Ap=1dB。要求学生首先用窗函数法设计上面要求的三种滤波器,在MA TLAB中,可以利用函数fir1设计FIR滤波器,然后在用双线性变换法设计上面要求的三种滤波器;之后再利用函数butter和cheby1设

4、计上面要求的三种IIR滤波器。最后,利用MATLAB中的函数freqz画出各滤波器的频率响应。巴特沃思低通1wp=0.1*pi;ws=0.5673*pi;Rp=1;Rs=100;Fs=20000;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

5、,Fs; %用双线性变换法实现模拟滤波器到数字滤波器的转换H,W=freqz(bz,az; %绘制频率响应曲线figure(1plot(W*Fs/(2*pi,abs(Hgridxlabel('频率/Hz'ylabel('频率响应幅度'title('Butterworth'f1=filter(bz,az,x1; 巴特沃思高通2x1=wavread('d:lianxi.wav',20000;t=0:1/20000:(size(x1-1/20000;wp=2*pi*3000;ws=2*pi*1555;Rp=1;Rs=100;Fs=200

6、00;Ts=1/Fs;N,Wn=buttord(wp,ws,Rp,Rs,'s' %选择滤波器的最小阶数Z,P,K=buttap(N; %创建butterworth模拟滤波器Bap,Aap=zp2tf(Z,P,K;b,a=lp2hp(Bap,Aap,Wn;bz,az=bilinear(b,a,Fs; %用双线性变换法实现模拟滤波器到数字滤波器的转换H,W=freqz(bz,az; %绘制频率响应曲线figure(1plot(W*Fs/(2*pi,abs(Hgridxlabel('频率/Hz'ylabel('频率响应幅度'title('But

7、terworth'f1=filter(bz,az,x1; 巴特沃思带通3wp=2*pi*1200,3000;ws=2*pi*212,5546;Rp=1;Rs=100;Fs=20000;Ts=1/Fs;N,Wn=buttord(wp,ws,Rp,Rs,'s' %选择滤波器的最小阶数B,A=butter(N,Wn,'s' %创建buttord模拟滤波器Bz,Az=bilinear(B,A,Fs; %用双线性变换法实现模拟滤波器到数字滤波器的转换H,W=freqz(Bz,Az; %绘制频率响应曲线figure(1plot(W*Fs/(2*pi,abs(Hgri

8、dxlabel('频率/Hz'ylabel('频率响应幅度'title('buttord' 切比雪夫1低通fs=20000;x1=wavread('d:lianxi.wav',20000;t=0:1/20000:(size(x1-1/20000;wp=0.1*pi;ws=0.5673*pi;Rp=1;Rs=100;Fs=20000;Ts=1/Fs;wp1=2/Ts*tan(wp/2; %将模拟指标转换成数字指标ws1=2/Ts*tan(ws/2;N,Wn=cheb1ord(wp1,ws1,Rp,Rs,'s' %选择

9、滤波器的最小阶数Z,P,K=cheby1(N,Rp,Wn,'s' %创建cheby1模拟滤波器Bap,Aap=zp2tf(Z,P,K;bz,az=bilinear(Bap,Aap,Fs; %用双线性变换法实现模拟滤波器到数字滤波器的转换H,W=freqz(bz,az; %绘制频率响应曲线figure(1plot(W*Fs/(2*pi,abs(Hgridxlabel('频率/Hz'ylabel('频率响应幅度'title('cheby1'f1=filter(bz,az,x1; 切比雪夫1高通clcclear allfs=20000;

10、x1=wavread('d:lianxi.wav',20000;t=0:1/20000:(size(x1-1/20000;wp=2*pi*3000;ws=2*pi*1500;Rp=1;Rs=100;Fs=20000;Ts=1/Fs;N,Wn=cheb1ord(wp,ws,Rp,Rs,'s' %选择滤波器的最小阶数Z,P,K=cheb1ap(N,Rp;Wpn=2*pi*3000; %创建cheby1模拟滤波器Bap,Aap=zp2tf(Z,P,K;b,a=lp2hp(Bap,Aap,Wpn;bz,az=bilinear(b,a,Fs; %用双线性变换法实现模拟滤波

11、器到数字滤波器的转换H,W=freqz(bz,az; %绘制频率响应曲线figure(1plot(W*Fs/(2*pi,abs(Hgridxlabel('频率/Hz'ylabel('频率响应幅度'title('cheby1'f1=filter(bz,az,x1; 切比雪夫1带通wp=2*pi*1200,3000;ws=2*pi*212,5546;Rp=1;Rs=100;Fs=20000;Ts=1/Fs;N,Wn=cheb1ord(wp,ws,Rp,Rs,'s' %选择滤波器的最小阶数B,A=cheby1(N,Rp,Wn,'

12、s' %创建cheby1模拟滤波器Bz,Az=bilinear(B,A,Fs; %用双线性变换法实现模拟滤波器到数字滤波器的转换H,W=freqz(Bz,Az; %绘制频率响应曲线figure(1plot(W*Fs/(2*pi,abs(Hgridxlabel('频率/Hz'ylabel('频率响应幅度'title('cheby1'f1=filter(Bz,Az,x1; 窗函数法低通Fs=20000;fp1=1200,fp2=3000,fs1=212,fs2=5546;wp1=2*pi*fp1/Fs,ws1=2*pi*fs1/Fs,wp2=

13、2*pi*fp2/Fs,ws2=2*pi*fs2/Fs,As=100; DB=wp1-ws1;beta=0.112*(As-8.7;M=ceil(As-8/2.285/DBwc=(wp1+ws1/2/pi,(wp2+ws2/2/pi;hn=fir1(M,wc,kaiser(M+1,beta;H,f=freqz(hn,1,1024,Fs;mag=abs(H;db=20*log10(mag+eps/max(mag;figure(1;plot(f/pi/Fs,db;xlabel('f/pi/Fs'ylabel('20logH(ejw'title('h(n的幅频

14、特性' grid;pha=angle(H;figure(2;plot(f,pha;xlabel('n'ylabel(''title('h(n的相频特性'figure(3;plot(hn;xlabel('n'ylabel('h(n' 高通Fs=20000;fp=3000,fs=1500;wp=2*pi*fp/Fs,ws=2*pi*fs/Fs,As=100;DB=wp-ws;beta=0.112*(As-8.7;M=ceil(As-8/2.285/DB;wc=(wp+ws/2/pi;hn=fir1(M,wc,&

15、#39;high',kaiser(M+1,beta;H,f=freqz(hn,1,1024,Fs;mag=abs(H;db=20*log10(mag+eps/max(mag;figure(1;plot(f/pi/Fs,db;xlabel('f/pi/Fs'ylabel('20logH(ejw'title('h(n的幅频特性' grid;pha=angle(H;figure(2;plot(f,pha;xlabel('n'ylabel(''title('h(n的相频特性'figure(3;plo

16、t(hn;xlabel('n'ylabel('h(n' 带通Fs=20000;fp1=1200,fp2=3000,fs1=212,fs2=5546;wp1=2*pi*fp1/Fs,ws1=2*pi*fs1/Fs,wp2=2*pi*fp2/Fs,ws2=2*pi*fs2/Fs,As=100; DB=wp1-ws1;beta=0.112*(As-8.7;M=ceil(As-8/2.285/DBwc=(wp1+ws1/2/pi,(wp2+ws2/2/pi;hn=fir1(M,wc,kaiser(M+1,beta;H,f=freqz(hn,1,1024,Fs;mag=a

17、bs(H;db=20*log10(mag+eps/max(mag;figure(1;plot(f/pi/Fs,db;xlabel('f/pi/Fs'ylabel('20logH(ejw'title('h(n的幅频特性'grid;pha=angle(H;figure(2;plot(f,pha;xlabel('n'ylabel(''title('h(n的相频特性'figure(3;plot(hn;xlabel('n'ylabel('h(n' 4.用滤波器对信号进行滤波比较

18、FIR和IIR两种滤波器的性能,然后用性能好的各滤波器分别对采集的信号进行滤波,在MATLAB中, FIR滤波器利用函数fftfilt对信号进行滤波,IIR滤波器利用函数filter对信号进行滤波。通过比较巴特沃思性能比较好6.回放语音信号在MATLAB中,函数sound可以对声音进行回放。其调用格式:sound(x,fs,bits;可以感觉滤波前后的声音有变化。低通:fs=20000;x1=wavread('d:lianxi.wav',20000;t=0:1/20000:(size(x1-1/20000;wp=0.1*pi;ws=0.5673*pi;Rp=1;Rs=100;F

19、s=20000;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(1plot(W*Fs/(2*pi,abs(Hgridxl

20、abel('频率/Hz'ylabel('频率响应幅度'title('Butterworth'f1=filter(bz,az,x1;figure(2subplot(2,1,1plot(t,x1 %画出滤波前的时域图title('滤波前的时域波形'subplot(2,1,2plot(t,f1; %画出滤波后的时域图title('滤波后的时域波形'sound(f1,20000; %播放滤波后的信号F0=fft(f1,1024;f=fs*(0:511/1024;figure(3y1=fft(x1,1024;subplot(

21、2,1,1;plot(f,abs(y1(1:512; %画出滤波前的频谱图title('滤波前的频谱'xlabel('Hz'ylabel('fuzhi'subplot(2,1,2F1=plot(f,abs(F0(1:512; %画出滤波后的频谱图title('滤波后的频谱'xlabel('Hz'ylabel('fuzhi'图 声音略高通:程序fs=20000;x1=wavread('d:lianxi.wav',20000;t=0:1/20000:(size(x1-1/20000;wp

22、=2*pi*3000;ws=2*pi*1555;Rp=1;Rs=100;Fs=20000;Ts=1/Fs;N,Wn=buttord(wp,ws,Rp,Rs,'s' %选择滤波器的最小阶数Z,P,K=buttap(N; %创建butterworth模拟滤波器Bap,Aap=zp2tf(Z,P,K;b,a=lp2hp(Bap,Aap,Wn;bz,az=bilinear(b,a,Fs; %用双线性变换法实现模拟滤波器到数字滤波器的转换H,W=freqz(bz,az; %绘制频率响应曲线figure(1plot(W*Fs/(2*pi,abs(Hgridxlabel('频率/Hz

23、'ylabel('频率响应幅度'title('Butterworth'f1=filter(bz,az,x1;figure(2subplot(2,1,1plot(t,x1 %画出滤波前的时域图title('滤波前的时域波形'subplot(2,1,2plot(t,f1; %画出滤波后的时域图title('滤波后的时域波形'sound(f1,20000; %播放滤波后的信号F0=fft(f1,1024;f=fs*(0:511/1024;figure(3y1=fft(x1,1024;subplot(2,1,1;plot(f,abs(y1(1:512; %画出滤波前的频谱图title('滤波前的频谱'xlabel('Hz'ylabel('fuzhi'subplot(2,1,2F1=plot(f,abs(F0(1:512; %画出滤波后的频谱图title('滤波后的频谱'xlabel('Hz'ylabel('fuzhi' 带通:程序fs=20000;