现代信号处理课程设计报告.doc

现代信号处理课程设计报告 目录一、 课程设计目的2二、 课程设计内容2三、 课程设计流程及分析43.1熟悉matlab基本操作. 43.2 matlab数字信号处理图形界面. 53.3用matlab绘出信号波形. .5 3.4对信号频谱分析. . . . . . . . .73.5设计滤波器并对信号滤波.11四、 设计经验总结27五、 参考资料27一、课程设计目的1全面复习课程所学理论知识，巩固所学知识重点和难点，将理论与实践很好地结合起来。2提高综合运用所学知识独立分析和解决问题的能力；3熟练使用一种高级语言进行编程实现。二、课程设计内容（一）信号分析1、编制信号生成程序，产生下述各序列，绘出它们的时域波形 1） 单位抽样序列 2）矩形序列 3） 三角波序列4) 反三角波序列 5）gaussian序列6） 正弦序列 取7）衰减正弦序列 对连续信号 进行采样，可得到测试序列 。令a=50,采样周期t=1ms,即fs=1000hz,f0=62.5,a=100。2. 对上述信号完成下列信号分析 1）对三角波序列和反三角波序列，作n=8点的fft，观察比较它们的幅频特性，说明它们有什么异同？绘出两序列及其它们的幅频特性曲线。 在和的尾部补零，作n=16点的fft，观察它们的幅频特性发生了什么变化？分析说明原因。 2)、观察高斯序列，固定信号中的参数p=8，令q分别等于2，4，8，观察它们的时域和幅频特性，了解当q取不同值时，对信号序列的时域幅频特性的影响；固定q=8，令p分别等于8，13，14，观察参数p变化对信号序列的时域及幅频特性的影响，观察p等于多少时，会发生明显的泄漏现象，混叠是否也随之出现？记录实验中观察到的现象，绘出相应的时域序列和幅频特性曲线。3）对于正弦序列x6(n)，取数据长度n分别等于8，16，32，分别作n点fft，观察它们的的时域和幅频特性，说明它们的差别 ，简要说明原因。4)、观察衰减正弦序列的时域和幅频特性，绘出幅频特性曲线，改变采样频率fs,使fs=300hz， 观察此时的频谱的形状和谱峰出现位置？说明产生现象的原因。3设有一连续时间信号s(t),其由20hz、220hz和750hz的正弦信号叠加而成，分析确定采样频率及数据分析长度，计算并绘出信号的频谱，指出各个频率份量。 （二）、数字滤波器的设计与实现 设计和实现iir或fir数字滤波器： 要求：输入数字滤波器的滤波指标,包括通带截止频率,通带最大衰减,阻带截止频率,阻带最小衰减。设计出相应的数字滤波器。显示得到的滤波器的系数和阶数，求出系统的单位脉冲响应h(n)。并画出设计得到的滤波器的幅频特性图(要有坐标标度)和相频特性。 1对于iir数字滤波器：使用双线性变换法。可选用bw型、cb型和椭圆型等滤波器类型来进行设计2对于fir数字滤波器：要对多种窗口（三种及以上）、不同阶数设计结果进行分析比较。（三）、综合应用对于给出的任意一段语音信号，进行分析处理。1、 读出语音信号，并显示语音信号的时域波形，然后对读出的信号进行fft变换，得到信号的频谱特性。2、 分析其频谱特性，应用设计的滤波器将主要的频率分量分离出来。3、 回放此语音信号，感觉滤波前后的声音有什么变化三、课程设计流程3.1熟悉matlab基本操作，基本语法，以及gui界面设计中各控件的调用方式。matlab是一款功能超强的数学软件，应用于各个行业。而其基本操作却很大众化，操作起来很人性化。其中所附带的gui面向对象的友好编程方式更便于学术交流和设计开发。 3.2 matlab数字信号处理图形界面3.3用matlab绘出信号波形global xypopup_sel_index = get(handles.popupmenu2, value);switch popup_sel_index case 1 n0=0;n1=-10;n2=10;n=n1:n2;x=(n-n0)=0;stem(n,x,.);xlabel(n);ylabel(x(n);title(单位抽样序列); case 2 n=-1:10;n0=2;width=5;x=rectpuls(n-n0,width);stem(n,x,.);xlabel(n);ylabel(x(n);axis(-1 10 0 2);title(矩形序列); case 3 n1=0:3;n2=4:7;x1=n1+1;x2=8-n2;n=(-10:0) n1 n2 (7:17);x=zeros(1,11) x1 x2 zeros(1,11);stem(n,x,.);xlabel(n);ylabel(x(n);title(三角波序列); case 4 n1=0:3;n2=4:7;x1=4-n1;x2=n2-3;n=(-10:0) n1 n2 (7:17);x=zeros(1,11) x1 x2 zeros(1,11);stem(n,x,.);xlabel(n);ylabel(x(n);title(反三角波序列); case 5 n=0:15;p=8;q=2;x=exp(-(n-p).2/q);stem(n,x,.);xlabel(n);ylabel(x(n);title(gaussian序列); case 6 t=1/64;n=0:t:16;phase=0;a=1;arg=2*pi*n-phase;x=a*sin(arg);stem(n,x,.);xlabel(n);ylabel(x(n);title(正弦序列); case 7 t=0.001;n=0:t:0.2;x1=sin(2*62.5*n*pi);x2=exp(-100*n);x=x1.*x2;stem(n,x,.);xlabel(n);ylabel(x(n);title(衰减正弦序列); case 8f1=20;f2=220;f3=750;n=0:0.5/f3:1;x=sin(2*pi*f1*n)+sin(2*pi*f2*n)+sin(2*pi*f3*n);stem(n,x,.);xlabel(n);ylabel(x(n);title(正弦叠加序列); endxy=x;3.4对信号频谱分析三角波的8点傅里叶变换： 三角波的32点傅里叶变换：对信号n点离散傅里叶变换相当于序列傅里叶变换在频域里的n点采样,所以采样点数越高越能反映信号的频谱。gaussian序列p=8，gaussian序列的15点傅里叶变换p=8，gaussian序列的15点傅里叶变换由图可知，如果n=15，当p=13时，会发生明显的泄漏现象，混叠现象也随之出现。当采样频率fs=300时，衰减正弦序列的频谱图当采样频率fs=3000时，衰减正弦序列的频谱图采样频率为fs=300时，由于采样频率太小，产生频谱混叠显现。当采样频率为1500时,正弦叠加序列的频谱图当采样频率为1500时,正弦叠加序列的频谱图由于正弦叠加序列的频率由2*20、2*220、 2* 750三种频率组成，所以采样频率应大于或等于1500；由于三个正弦序列的周期都为1，所以截取0到1进行傅里叶分析。3.5设计滤波器并对信号滤波 global yglobal fsviir=get(handles.iir,value);vfir=get(handles.fir,value);vbw=get(handles.bw,value);vcb=get(handles.cb,value);vclass1=get(handles.class1,value);vclass2=get(handles.class2,value);vclass3=get(handles.class3,value);vclass4=get(handles.class4,value);vpopupmenu1=get(handles.popupmenu1, value);vp1=str2num(get(handles.p1,string);vs1=str2num(get(handles.s1,string);vp2=str2num(get(handles.p2,string);vs2=str2num(get(handles.s2,string);vap=str2num(get(handles.ap,string);vas=str2num(get(handles.as,string);fs1=str2num(get(handles.fs,string);vs=get(handles.checkbox1,value);if(viir) if(vs) t=1/fs; else t=1/fs1; end if(vbw) if(vclass1) wp=2/t*tan(vp1*pi/2); ws=2/t*tan(vs1*pi/2); n,wc=buttord(wp,ws,vap,vas,s); b,a=butter(n,wc,s); elseif(vclass2) wp=2/t*tan(vp1*pi/2); ws=2/t*tan(vs1*pi/2); n,wc=buttord(wp,ws,vap,vas,s); b,a=butter(n,wc,high,s); elseif(vclass3) wp=2/t*tan(vp1*pi/2),2/t*tan(vp2*pi/2); ws=2/t*tan(vs1*pi/2),2/t*tan(vs2*pi/2); n,wc=buttord(wp,ws,vap,vas,s); b,a=butter(n,wc,s); else wp=2/t*tan(vp1*pi/2),2/t*tan(vp2*pi/2); ws=2/t*tan(vs1*pi/2),2/t*tan(vs2*pi/2); n,wc=buttord(wp,ws,vap,vas,s); b,a=butter(n,wc,stop,s); end else if(vclass1) wp=2/t*tan(vp1*pi/2); ws=2/t*tan(vs1*pi/2); n,wc=cheb1ord(wp,ws,vap,vas,s); b,a=cheby1(n,vap,wc,s); elseif(vclass2) wp=2/t*tan(vp1*pi/2); ws=2/t*tan(vs1*pi/2); n,wc=cheb1ord(wp,ws,vap,vas,s); b,a=cheby1(n,vap,wc,high,s); elseif(vclass3) wp=2/t*tan(vp1*pi/2),2/t*tan(vp2*pi/2); ws=2/t*tan(vs1*pi/2),2/t*tan(vs2*pi/2); n,wc=cheb1ord(wp,ws,vap,vas,s); b,a=cheby1(n,vap,wc,s); else wp=2/t*tan(vp1*pi/2),2/t*tan(vp2*pi/2); ws=2/t*tan(vs1*pi/2),2/t*tan(vs2*pi/2); n,wc=cheb1ord(wp,ws,vap,vas,s); b,a=cheby1(n,vap,wc,stop,s); end end if(vs) bz,az=bilinear(b,a,fs); else bz,az=bilinear(b,a,fs1); end if(vs) x=filter(bz,az,y); x=fft(x,4096); plot(0:length(x)-1,x);title(滤波后的信号波形); figure; f=(0:4095)/4096; plot(f*2,abs(x); xlabel(w/);ylabel(x|k|); title(滤波后的频谱); sound(x,fs); else wk=0:pi/512:pi; hz,w=freqz(bz,az,wk); hz1=ifft(hz); stem(w,hz1,.);figure; l=length(hz); hx=angle(hz); plot(2*(0:l-1),hx,.); xlabel(w/);ylabel(相位);figure; plot(wk/pi,abs(hz); xlabel(w/);ylabel(幅度x|k|);figure; plot(0:l-1)/l,20*log10(abs(hz);grid on; axis(0 1 -100 5); xlabel(w/);ylabel(幅度(db); end n wc b a else if(vclass1|vclass2) db=abs(vs1-vp1); wc=(vp1+vs1)/2; switch vpopupmenu1 case 1 if(vasvp1&vclass1) n=n0; hn=fir1(n,wc,boxcar(n+1); elseif(vs1vp1&vclass2) if(mod(n0+1,2) n=n0; else n=n0+1; end hn=fir1(n,wc,high,boxcar(n+1); else errordlg(输入错误，将通带截止频率与阻带截止频率互换,error); end else errordlg(矩形窗的阻带最小衰减小于或等于21,error); end case 2 if(vasvp1&vclass1) n=n0; hn=fir1(n,wc,bartlett(n+1); elseif(vs1vp1&vclass2) if(mod(n0+1,2) n=n0; else n=n0+1; end hn=fir1(n,wc,high,bartlett(n+1); else errordlg(输入错误，将通带截止频率与阻带截止频率互换,error); end else errordlg(三角窗的阻带最小衰减小于或等于25,error); end case 3 if(vasvp1&vclass1) n=n0; hn=fir1(n,wc,hanning(n+1); elseif(vs1vp1&vclass2) if(mod(n0+1,2) n=n0; else n=n0+1; end hn=fir1(n,wc,high,hanning(n+1); else errordlg(输入错误，将通带截止频率与阻带截止频率互换,error); end else errordlg(汉宁窗的阻带最小衰减小于或等于44,error); end case 4 if(vasvp1&vclass1) n=n0; hn=fir1(n,wc,hamming(n+1); elseif(vs1vp1&vclass2) if(mod(n0+1,2) n=n0; else n=n0+1; end hn=fir1(n,wc,high,hamming(n+1); else errordlg(输入错误，将通带截止频率与阻带截止频率互换,error); end else errordlg(哈明窗的阻带最小衰减小于或等于53,error); end case 5 if(vasvp1&vclass1) n=n0; hn=fir1(n,wc,blackman(n+1); elseif(vs1vp1&vclass2) if(mod(n0+1,2) n=n0; else n=n0+1; end hn=fir1(n,wc,high,blackman(n+1); else errordlg(输入错误，将通带截止频率与阻带截止频率互换,error); end else errordlg(布莱克曼窗的阻带最小衰减小于或等于74,error); end case 6 if(vas50) beta=0.112*(vas-8.7); elseif(vas=21) beta=0.5842*(vas-21)0.4+0.07886*(vas-21); else beta=0; end n0=ceil(vas-8)/2.285/db); if(vs1vp1&vclass1) n=n0; hn=fir1(n,wc,kaiser(n+1,beta); elseif(vs1vp1&vclass2) if(mod(n0+1,2) n=n0; else n=n0+1; end hn=fir1(n,wc,high,kaiser(n+1,beta); else errordlg(输入错误，将通带截止频率与阻带截止频率互换,error); end else errordlg(凯塞窗的阻带最小衰减小于或等于80,error); end end elseif(vclass3|vclass4) db=abs(vs1-vp1); wc=(vs1+vp1)/2,(vs2+vp2)/2; switch vpopupmenu1 case 1 if(vas=21) n0=ceil(1.8/db); if(vs1vp1&vclass4) if(mod(n0+1,2) n=n0; else n=n0+1; end hn=fir1(n,wc,stop,boxcar(n+1); else errordlg(输入错误，将通带截止频率与阻带截止频率互换,error); end else errordlg(矩形窗的阻带最小衰减小于或等于21,error); end case 2 if(vas=25) n0=ceil(6.1/db); if(vs1vp1&vclass4) if(mod(n0+1,2) n=n0; else n=n0+1; end hn=fir1(n,wc,stop,bartlett(n+1); else errordlg(输入错误，将通带截止频率与阻带截止频率互换,error); end else errordlg(三角窗的阻带最小衰减小于或等于25,error); end case 3 if(vas=44) n0=ceil(6.2/db); if(vs1vp1&vclass4) if(mod(n0+1,2) n=n0; else n=n0+1; end hn=fir1(n,wc,stop,hanning(n+1); else errordlg(输入错误，将通带截止频率与阻带截止频率互换,error); end else errordlg(汉宁窗的阻带最小衰减小于或等于44,error); end case 4 if(vas=53) n0=ceil(6.6/db); if(vs1vp1&vclass4) if(mod(n0+1,2) n=n0; else n=n0+1; end hn=fir1(n,wc,stop,hamming(n+1); else errordlg(输入错误，将通带截止频率与阻带截止频率互换,error); end else errordlg(哈明窗的阻带最小衰减小于或等于53,error); end case 5 if(vas=74) n0=ceil(11/db); if(vs1vp1&vclass4) if(mod(n0+1,2) n=n0; else n=n0+1; end hn=fir1(n,wc,stop,blackman(n+1); else errordlg(输入错误，将通带截止频率与阻带截止频率互换,error); end else errordlg(布莱克曼窗的阻带最小衰减小于或等于74,error); end case 6 if(vas50) beta=0.112*(vas-8.7); elseif(vas=21) beta=0.5842*(vas-21)0.4+0.07886*(vas-21); else beta=0; end n0=ceil(vas-8)/2.285/db); if(vs1vp1&vclass4) if(mod(n0+1,2) n=n0; else n=n0+1; end hn=fir1(n,wc,stop,kaiser(n+1,beta); else errordlg(输入错误，将通带截止频率与阻带截止频率互换,error); end else errordlg(凯塞窗的阻带最小衰减小于或等于80,error); end end end if(vs) x=fftfilt(hn,y); x=fft(x,4096); plot(0:length(x)-1,x); xlabel(n);ylabel(h(n);figure; f=(0:4095)/4096; plot(f*2,abs(x); xlabel(w/);ylabel(x|k|); title(滤波后的频谱); sound(x,fs); else stem(0:n,hn,.); xlabel(n);ylabel(h(n);figure; y=fft(hn,4096*2); m=1:4096; mag=abs(y(m); w=(0:4095)/4096; plot(w,mag);grid on; xlabel(w/);ylabel(幅度x|k|);figure; plot(w,20*log10(mag);grid on; xlabel(w/);ylabel(幅度db);figure; phay=angle(y(m); plot(2*(0:length(phay)-1)/length(phay),phay); xlabel(w/);ylabel(相位(w); end n end用双线性变换法设计iir滤波器和用窗函数法设计fir滤波器低通iir滤波器（通带截止频率0.2，阻带截止频率0.3，通带最大衰减系数1，阻带最小衰减系数40）序列图：频谱