《現(xiàn)代信號處理》課程設(shè)計

1、目 錄一、 課程設(shè)計目的-2二、 課程設(shè)計任務-2三、 詳細設(shè)計過程-4四、 設(shè)計結(jié)果及分析-13五、 心得體會-24一、課程設(shè)計目的：1全面復習課程所學理論知識，鞏固所學知識重點和難點，將理論與實踐很好地結(jié)合起來。2提高綜合運用所學知識獨立分析和解決問題的能力；3熟練使用一種高級語言進行編程實現(xiàn)。二、課程設(shè)計內(nèi)容（一）信號分析1、編制信號生成程序，產(chǎn)生下述各序列，繪出它們的時域波形 1） 單位抽樣序列 2）矩形序列 3） 三角波序列4) 反三角波序列 5）gaussian序列6） 正弦序列 取7）衰減正弦序列 對連續(xù)信號 進行采樣，可得到測試序列 。令a=50,采樣周期t=1ms,即fs=1

2、000hz,f0=62.5,a=100。2. 對上述信號完成下列信號分析 1）對三角波序列和反三角波序列，作n=8點的fft，觀察比較它們的幅頻特性，說明它們有什么異同？繪出兩序列及其它們的幅頻特性曲線。在和的尾部補零，作n=16點的fft，觀察它們的幅頻特性發(fā)生了什么變化？分析說明原因。2)、觀察高斯序列，固定信號中的參數(shù)p=8，令q分別等于2，4，8，觀察它們的時域和幅頻特性，了解當q取不同值時，對信號序列的時域幅頻特性的影響；固定q=8，令p分別等于8，13，14，觀察參數(shù)p變化對信號序列的時域及幅頻特性的影響，觀察p等于多少時，會發(fā)生明顯的泄漏現(xiàn)象，混疊是否也隨之出現(xiàn)？記錄實驗中觀察到

3、的現(xiàn)象，繪出相應的時域序列和幅頻特性曲線。3）對于正弦序列，取數(shù)據(jù)長度n分別等于8，16，32，分別作n點fft，觀察它們的的時域和幅頻特性，說明它們的差別 ，簡要說明原因。4)、觀察衰減正弦序列的時域和幅頻特性，繪出幅頻特性曲線，改變采樣頻率fs,使fs=300hz， 觀察此時的頻譜的形狀和譜峰出現(xiàn)位置？說明產(chǎn)生現(xiàn)象的原因。3設(shè)有一連續(xù)時間信號s(t),其由20hz、220hz和750hz的正弦信號疊加而成，分析確定采樣頻率及數(shù)據(jù)分析長度，計算并繪出信號的頻譜，指出各個頻率份量。 （二）、數(shù)字濾波器的設(shè)計與實現(xiàn) 設(shè)計和實現(xiàn)iir或fir數(shù)字濾波器： 要求：輸入數(shù)字濾波器的濾波指標,包括通帶截

4、止頻率,通帶最大衰減,阻帶截止頻率,阻帶最小衰減。設(shè)計出相應的數(shù)字濾波器。顯示得到的濾波器的系數(shù)和階數(shù)，求出系統(tǒng)的單位脈沖響應h(n)。并畫出設(shè)計得到的濾波器的幅頻特性圖(要有坐標標度)和相頻特性。 1對于iir數(shù)字濾波器：使用雙線性變換法?？蛇x用bw型、cb型和橢圓型等濾波器類型來進行設(shè)計2對于fir數(shù)字濾波器：要對多種窗口（三種及以上）、不同階數(shù)設(shè)計結(jié)果進行分析比較。（三）、綜合應用對于給出的任意一段語音信號，進行分析處理。1、 讀出語音信號，并顯示語音信號的時域波形，然后對讀出的信號進行fft變換，得到信號的頻譜特性。2、 分析其頻譜特性，應用設(shè)計的濾波器將主要的頻率分量分離出來。3、

5、回放此語音信號，感覺濾波前后的聲音有什么變化三、詳細設(shè)計過程1、新建一個空白的圖形gui界面2、添加axes及pop-up menu。可以用鼠標拖拽，將其放在合適的位置。pop-up menu下拉菜單可以根據(jù)任務自行設(shè)置名稱。3、在圖形界面生成的源程序中編程，以實現(xiàn)課程設(shè)計的各項任務。源程序如下：function varargout = wo(varargin)% wo m-file for wo.fig% wo, by itself, creates a new wo or raises the existing% singleton*.% h = wo returns the handle

6、 to a new wo or the handle to% the existing singleton*.% wo(callback,hobject,eventdata,handles,.) calls the local% function named callback in wo.m with the given input arguments.% wo(property,value,.) creates a new wo or raises the% existing singleton*. starting from the left, property value pairs a

7、re% applied to the gui before wo_openingfunction gets called. an% unrecognized property name or invalid value makes property application% stop. all inputs are passed to wo_openingfcn via varargin.% *see gui options on guides tools menu. choose gui allows only one% instance to run (singleton).% see a

8、lso: guide, guidata, guihandles % copyright 2002-2003 the mathworks, inc. % edit the above text to modify the response to help wo % last modified by guide v2.5 15-sep-2010 08:37:54 % begin initialization code - do not editgui_singleton = 1;gui_state = struct(gui_name, mfilename, . gui_singleton, gui

9、_singleton, . gui_openingfcn, wo_openingfcn, . gui_outputfcn, wo_outputfcn, . gui_layoutfcn, , . gui_callback, );if nargin & ischar(varargin1) gui_state.gui_callback = str2func(varargin1);end if nargout varargout1:nargout = gui_mainfcn(gui_state, varargin:);else gui_mainfcn(gui_state, varargin:);end

10、% end initialization code - do not edit % - executes just before wo is made visible.function wo_openingfcn(hobject, eventdata, handles, varargin)% this function has no output args, see outputfcn.% hobject handle to figure% eventdata reserved - to be defined in a future version of matlab% handles str

11、ucture with handles and user data (see guidata)% varargin command line arguments to wo (see varargin) % choose default command line output for wohandles.output = hobject; % update handles structureguidata(hobject, handles); % uiwait makes wo wait for user response (see uiresume)% uiwait(handles.figu

12、re1); % - outputs from this function are returned to the command line.function varargout = wo_outputfcn(hobject, eventdata, handles) % varargout cell array for returning output args (see varargout);% hobject handle to figure% eventdata reserved - to be defined in a future version of matlab% handles

13、structure with handles and user data (see guidata) % get default command line output from handles structurevarargout1 = handles.output; % - executes on selection change in popupmenu1.function popupmenu1_callback(hobject, eventdata, handles)% hobject handle to popupmenu1 (see gcbo)% eventdata reserve

14、d - to be defined in a future version of matlab% handles structure with handles and user data (see guidata) n = get(handles.popupmenu1, value);switch n case 1 n=1:50;x=zeros(1,50);x(1)=1;subplot(421); stem(x,.);xlabel(n);ylabel(x1(n); n=1:50;x=rectpuls(n-25,20);subplot(422);stem(x,.);xlabel(n);ylabe

15、l(x2(n); n1=0:3;x1=n1+1; n2=4:7;x2=8-n2;n=n1 n2; x=x1 x2;subplot(423);stem(x,.);xlabel(n);ylabel(x3(n); n1=0:3;x1=4-n1;n2=4:7;x2=n2-3;n=n1 n2;x=x1 x2;subplot(424);stem(x,.);xlabel(n);ylabel(x4(n); n=0:50;p=8;q=2;x=exp(-1*(n-p).2/q);subplot(425);stem(x,.);xlabel(n);ylabel(x5(n); fs=64;n=16;n=0:2*pi*f

16、s/n;x=sin(n*pi/4);subplot(426);stem(x,.);xlabel(n);ylabel(x6(n); n=0:15; a=0.1; f=0.0625; x=50*exp(-a*n).*sin(2*pi*f*n);subplot(427);stem(x,.);xlabel(n);ylabel(x7(n); case 2 n1=0:3;x1=n1+1; n2=4:7;x2=8-n2;n=n1 n2; x=x1 x2; y=fft(x,8);subplot(421);stem(y,.);xlabel(k);ylabel(x3(k);grid; magy=abs(y);k=

17、0:7;w=2*pi/8*k;subplot(422); stem(w/pi,magy,.);xlabel(n=8);ylabel(|x3(k)|);grid; x=x zeros(1,8);y=fft(x,16);subplot(425);stem(y,.);xlabel(k);ylabel(x3(k);grid; magy=abs(y);k=0:15;w=2*pi/16*k;subplot(426); stem(w/pi,magy,.);xlabel(n=16);ylabel(|x3(k)|);grid; n1=0:3;x1=4-n1;n2=4:7;x2=n2-3;n=n1 n2;x=x1

18、 x2; y=fft(x,8);subplot(423);stem(y,.);xlabel(k);ylabel(x4(k);grid; magy=abs(y);k=0:7;w=2*pi/8*k;subplot(424); stem(w/pi,magy,.);xlabel(n=8);ylabel(|x4(k)|);grid; x=x zeros(1,8);y=fft(x,16);subplot(427);stem(y,.);xlabel(k);ylabel(x4(k);grid; magy=abs(y);k=0:15;w=2*pi/16*k;subplot(428); stem(w/pi,mag

19、y,.);xlabel(n=16);ylabel(|x4(k)|);grid; case 3 figure(1);n=0:15;p=8;q=2;x=exp(-1*(n-p).2/q);subplot(321);stem(x,.);title(高斯，p=8,q=2); y=fft(x);magy=abs(y);subplot(322);stem(magy,.);title(幅頻特性，p=8,q=2); n=0:15;p=8;q=4;x=exp(-1*(n-p).2/q);subplot(323);stem(x,.);title(高斯，p=8,q=4); y=fft(x);magy=abs(y);

20、subplot(324);stem(magy,.);title(幅頻特性，p=8,q=2); n=0:15;p=8;q=8;x=exp(-1*(n-p).2/q);subplot(325);stem(x,.);title(高斯，p=8,q=8); y=fft(x);magy=abs(y);subplot(326);stem(magy,.);title(幅頻特性，p=8,q=2); figure(2);n=0:15;p=8;q=8;x=exp(-1*(n-p).2/q);subplot(321);stem(x,.);title(高斯，p=8,q=8); y=fft(x);magy=abs(y);

21、subplot(322);stem(magy,.);title(幅頻特性，p=8,q=8); n=0:15;p=13;q=8;x=exp(-1*(n-p).2/q);subplot(323);stem(x,.);title(高斯，p=13,q=8); y=fft(x);magy=abs(y);subplot(324);stem(magy,.);title(幅頻特性，p=13,q=8); n=0:15;p=14;q=8;x=exp(-1*(n-p).2/q);subplot(325);stem(x,.);title(高斯，p=14,q=8); y=fft(x);magy=abs(y);subpl

22、ot(326);stem(magy,.);title(幅頻特性，p=14,q=8); case 4 fs=64;n=8;n=0:2*pi*fs/n;x=sin(n*pi/4);subplot(321);stem(x,.);title(正弦，n=8);y=abs(fft(x,8);subplot(322);stem(y,.);title(正弦幅頻特性，n=8); fs=64;n=16;n=0:2*pi*fs/n;x=sin(n*pi/4);subplot(323);stem(x,.);title(正弦，n=16);y=abs(fft(x,16);subplot(324);stem(y,.);ti

23、tle(正弦幅頻特性，n=16); fs=64;n=20;n=0:2*pi*fs/n;x=sin(n*pi/4);subplot(325);stem(x,.);title(正弦，n=20);y=abs(fft(x,20);subplot(326);stem(y,.);title(正弦幅頻特性，n=20); case 5 n=0:15; a=0.1; f=0.0625; x=50*exp(-a*n).*sin(2*pi*f*n);subplot(221);stem(x,.);title(fs=1000hz); y=abs(fft(x);subplot(222);stem(y,.);title(f

24、s=1000hz,幅頻特性); n=0:15; a=1/3; f=5/24; x=50*exp(-a*n).*sin(2*pi*f*n);subplot(223);stem(x,.);title(fs=300hz); y=abs(fft(x);subplot(224);stem(y,.);title(fs=300hz，幅頻特性); case 6 fs=2000;n=200;n=0:n-1;t=n/fs;x=sin(pi*40*t)+sin(pi*440*t)+sin(pi*1500*t);subplot(2,1,1);plot(x); y=fft(x,n);mag=abs(y);f=(0:le

25、ngth(y)-1).*fs/length(y);subplot(212);plot(f(1:n/2),mag(1:n/2);xlabel(f/hz); case 7 wp=input(please enter wp:); ws=input(please enter ws:); rp=input(please enter rp:); rs=input(please enter rs:); k=input(請選擇濾波器類型 1、cb型 2、bw型 3、橢圓型:); switch k case 1 %cb型 n,wn= cheb1ord(wp,ws,rp,rs) %階數(shù)n和截止頻率wn b,a=c

26、heby1(n,rp,wn) %系統(tǒng)函數(shù)系數(shù)向量b和a figure(1); freqs(b,a);h,w = freqz(b,a); magh = abs(h);phah = angle(h); num=b; den=a; figure(2);subplot(311);plot(w/pi,magh);grid xlabel(frequency unit:pi);ylabel(magnitude); title(幅頻特性); subplot(312);plot(w/pi,phah/pi);grid xlabel(frequency unit:pi);ylabel(phase unit:pi);

27、 title(相頻特性); subplot(313); impulse(num,den); title(系統(tǒng)的單位脈沖響應h(n); case 2 %bw型 n,wn=buttord(wp,ws,rp,rs) %階數(shù)n和截止頻率wn b,a=butter(n,wn) %系統(tǒng)函數(shù)系數(shù)向量b和a figure(1); freqs(b,a);h,w = freqz(b,a); magh = abs(h);phah = angle(h); num=b; den=a; figure(2);subplot(311);plot(w/pi,magh);grid xlabel(frequency unit:pi

28、);ylabel(magnitude); title(幅頻特性); subplot(3,1,2);plot(w/pi,phah/pi);grid xlabel(frequency unit:pi);ylabel(phase unit:pi); title(相頻特性); subplot(313); impulse(num,den); title(系統(tǒng)的單位脈沖響應h(n); case 3 %橢圓型 n,wn=ellipord(wp,ws,rp,rs) %階數(shù)n和截止頻率wn b,a=ellip(n,rp,rs,wn) %系統(tǒng)函數(shù)系數(shù)向量b和a figure(1); freqs(b,a); h,w

29、 = freqz(b,a); magh = abs(h);phah = angle(h); num=b; den=a; figure(2);subplot(311);plot(w/pi,magh);grid xlabel(frequency unit:pi);ylabel(magnitude); title(幅頻特性); subplot(312);plot(w/pi,phah/pi);grid xlabel(frequency unit:pi);ylabel(phase unit:pi); title(相頻特性); subplot(313); impulse(num,den); title(系

30、統(tǒng)的單位脈沖響應h(n); end case 8 n=21;wc=1/4;n=0:20; hn=fir1(n-1,wc,boxcar(n); fh=fft(hn,1024); fh=20*log10(abs(fh); wk=0:1023;wk=2*wk/1024; subplot(421);stem(n,hn);grid; title(矩形窗設(shè)計結(jié)果：); subplot(422);plot(wk,fh);grid; title(矩形窗);xlabel(w/pi); wind=hanning(n); hn=hn.*wind; fh=fft(hn,1024); fh=20*log10(abs(f

31、h); subplot(423);stem(n,hn);grid; title(漢寧窗設(shè)計結(jié)果：); subplot(424);plot(wk,fh);grid; title(漢寧窗);xlabel(w/pi); axis(0 2 -130 0); wind=hamming(n); hn=hn.*wind; fh=fft(hn,1024); fh=20*log10(abs(fh); subplot(425);stem(n,hn);grid; title(哈明窗設(shè)計結(jié)果：); subplot(426);plot(wk,fh);grid; title(哈明窗);xlabel(w/pi); axis

32、(0 2 -160 0); wind=blackman(n); hn=hn.*wind; fh=fft(hn,1024); fh=20*log10(abs(fh); subplot(427);stem(n,hn);grid; title(布萊克曼窗設(shè)計結(jié)果：); subplot(428);plot(wk,fh);grid; title(布萊克曼窗);xlabel(w/pi); case 9 fb=1000;fc=1200;as=100;ap=1;fs=22050;x2=wavread(f:sha we.wav );ts=1/fs;r1=10;wp=2*pi*fb/fs;ws=2*pi*fc/f

33、s;wp1=2/ts*tan(wp/2); %將模擬指標轉(zhuǎn)換成數(shù)字指標ws1=2/ts*tan(ws/2);n,wn=buttord(wp1,ws1,ap,r1,s); %選擇濾波器的最小階數(shù)z,p,k=buttap(n); %創(chuàng)建butterworth模擬濾波器bap,aap=zp2tf(z,p,k);b,a=lp2lp(bap,aap,wn); bz,az=bilinear(b,a,fs); %用雙線性變換法實現(xiàn)模擬濾波器到數(shù)字濾波器的轉(zhuǎn)換h,w=freqz(bz,az); %繪制頻率響應曲線figure(1)plot(w*fs/(2*pi),abs(h);gridxlabel(頻率hz)

34、ylabel(頻率響應幅度)title(iir低通濾波器)f1=filter(bz,az,x2);figure(2)subplot(2,1,1)plot(x2) %畫出濾波前的時域圖title(iir低通濾波器濾波前的時域波形);subplot(2,1,2)plot(f1); %畫出濾波后的時域圖title(iir低通濾波器濾波后的時域波形);sound(f1,44100); %播放濾波后的信號f0=fft(f1,1024);f=fs*(0:511)/1024;figure(3)y2=fft(x2,1024);subplot(2,1,1);plot(f,abs(y2(1:512); grid

35、%畫出濾波前的頻譜圖title(iir低通濾波器濾波前的頻譜)xlabel(頻率/hz);ylabel(幅值);subplot(2,1,2)f1=plot(f,abs(f0(1:512); grid %畫出濾波后的頻譜圖title(iir低通濾波器濾波后的頻譜)xlabel(頻率/hz);ylabel(幅值);end% hints: contents = get(hobject,string) returns popupmenu1 contents as cell array% contentsget(hobject,value) returns selected item from popu

36、pmenu1 % - executes during object creation, after setting all properties.function popupmenu1_createfcn(hobject, eventdata, handles)% hobject handle to popupmenu1 (see gcbo)% eventdata reserved - to be defined in a future version of matlab% handles empty - handles not created until after all createfc

37、ns called % hint: popupmenu controls usually have a white background on windows.% see ispc and computer.if ispc & isequal(get(hobject,backgroundcolor), get(0,defaultuicontrolbackgroundcolor) set(hobject,backgroundcolor,white);end% - executes during object creation, after setting all properties.funct

38、ion axes1_createfcn(hobject, eventdata, handles)% hobject handle to axes1 (see gcbo)% eventdata reserved - to be defined in a future version of matlab% handles empty - handles not created until after all createfcns called % hint: place code in openingfcn to populate axes1四、設(shè)計結(jié)果及分析（一）、信號分析部分1、單位抽樣序列

39、、矩形序列 、三角波序列、反三角波序列、gaussian序列、正弦序列、衰減正弦序列的波形如下圖所示：2、對三角波序列和反三角波序列，作n=8點的fft；在和的尾部補零，作n=16點的fft。如下圖所示：由上圖可知：當求n=8點的fft時，所得的三角波和反三角波的幅頻特性相差不大，基本看不出差別。當取n=16時的幅頻特性差別明顯。對n=16點的幅頻特性進行一個間隔一點的采樣即得到n=8點的幅頻特性。由所學知識可知求fft的公式是（k=0,1.n-1）其中幅度()n的取值越大所得的包絡(luò)越接近y的dtft的幅度譜。3、固定信號中的參數(shù)p=8，令q分別等于2，4，8；固定q=8，令p分別等于8，13

40、，14；如下圖所示：由上圖可知：當固定p=8,而q在取不同值時它們的時域和幅頻度不一樣。q的值越大，時域圖形的零點越少，圖形越“胖”，值越分散。q的值越大，頻域圖形的零點越多，但其他點的幅值相對增大，發(fā)生頻譜泄露以及混疊。這是因為值越分散在經(jīng)過窗函數(shù)時，n的值就越大，泄漏也就越大。由上圖可知：當固定q=8,而p在取不同值時它們的時域和幅頻度不一樣。p的值越大，時域圖形右移，形狀不變。p的值越大，頻域圖形的頻譜泄露就越明顯并伴隨有頻譜的混疊。這是因為譜線右移，經(jīng)過窗函數(shù)時n的取值會相應增大，從而使泄露以及混疊現(xiàn)象越明顯。4、對于正弦序列，取數(shù)據(jù)長度n分別等于8，16，32，分別作n點fft。如下

41、圖所示：由上圖可知：當n越大時，所求得的fft的幅值就越密，包絡(luò)曲線就越接近原函數(shù)的包絡(luò)。這是因為n與它的周期緊密相連，序列fft的幅值是關(guān)于n/2對稱的。5、觀察衰減正弦序列的時域和幅頻特性，繪出幅頻特性曲線，改變采樣頻率fs,使fs=300hz。如下圖所示：采樣我們知道是對信號進行時間上的量化，即采樣后頻譜發(fā)生變化。結(jié)合我們所學的理論的知識和上圖體現(xiàn)出不管是序列的時域還是幅頻特性它們的圖形都與采樣周期相關(guān)聯(lián)，時域的衰減長度及快慢與采樣頻率聯(lián)系，fs越大其的衰減長度就越長衰減的也越慢，當fs是300是序列衰減長度是fs=1000的300/1000倍即衰減的速度是fs=1000的300/100

42、0倍。同時它的幅頻特性也發(fā)生了變化，fs=300hz的一個完整的幅頻長度是fs=1000hz的長度發(fā)的3/10倍，從上圖我們可以看出在實際中我們不一定只能根據(jù)信號的寬度來選擇濾波器，也可根據(jù)所得濾波器去去調(diào)整信號寬度減少信號的泄露。6、設(shè)有一連續(xù)時間信號s(t),其由20hz、220hz和750hz的正弦信號疊加而成，分析確定采樣頻率及數(shù)據(jù)分析長度，計算并繪出信號的頻譜，指出各個頻率份量。如下圖所示：根據(jù)連續(xù)信號的疊加性質(zhì)，將幾個簡單信號疊加即可。根據(jù)采樣定理fs要大于等于2倍的fm，在此可以取fs=2000hz。7、對于iir數(shù)字濾波器：使用雙線性變換法。可選用bw型、cb型和橢圓型等濾波器

43、類型來進行設(shè)計。cb型bw型 橢圓型 由上圖可知：在同等條件下（輸入的wp=.2,ws=.3,rp=1,rs=34）選擇不同的類型濾波器得出來的各項參數(shù)也不一樣。各項參數(shù)如下所示：cb型n = 6wn =0.2000b =1.0e-003 * 0.0464 0.2782 0.6956 0.9274 0.6956 0.2782 0.0464a = 1.0000 -4.8694 10.3807 -12.3367 8.5969 -3.3268 0.5586bw型n =11wn = 0.2183b=1.0e-003*0.0010 0.0111 0.0556 0.1668 0.3336 0.4670 0.4670 0.3336 0.1668 0.0556 0.0111a =1.0000 -6.1863