數(shù)字信號(hào)處理實(shí)驗(yàn)全部程序MATLAB

1、.專業(yè)整理 .實(shí)驗(yàn)一熟悉 MATLAB環(huán)境一、實(shí)驗(yàn)?zāi)康?1) 熟悉 MATLAB的主要操作命令。(2) 學(xué)會(huì)簡(jiǎn)單的矩陣輸入和數(shù)據(jù)讀寫(xiě)。(3) 掌握簡(jiǎn)單的繪圖命令。(4) 用 MATLAB編程并學(xué)會(huì)創(chuàng)建函數(shù)。(5) 觀察離散系統(tǒng)的頻率響應(yīng)。二、實(shí)驗(yàn)容認(rèn)真閱讀本章附錄，在 MATLAB環(huán)境下重新做一遍附錄中的例子，體會(huì)各條命令的含義。在熟悉了 MATLAB基本命令的基礎(chǔ)上，完成以下實(shí)驗(yàn)。上機(jī)實(shí)驗(yàn)容：(1) 數(shù)組的加、減、乘、除和乘方運(yùn)算。 輸入 A=1 2 3 4 ，B=3 4 5 6 ，求 C=A+B，D=A-B，E=A.*B，F(xiàn)=A./B ，G=A.B并用 stem 語(yǔ)句畫(huà)出 A、B、C、D、

2、E、F、G。實(shí)驗(yàn)程序：A=1 2 3 4;B=3 4 5 6;n=1:4;C=A+B;D=A-B;E=A.*B;F=A./B;G=A.B;subplot(4,2,1);stem(n,A,fill);xlabel ( 時(shí)間序列 n );ylabel(A);subplot(4,2,2);stem(n,B,fill);xlabel ( 時(shí)間序列 n );ylabel(B);subplot(4,2,3);stem(n,C,fill);xlabel ( 時(shí)間序列 n );ylabel(A+B);subplot(4,2,4);stem(n,D,fill);xlabel ( 時(shí)間序列 n );ylabel(

3、A-B);subplot(4,2,5);stem(n,E,fill);xlabel ( 時(shí)間序列 n );ylabel(A.*B);subplot(4,2,6);stem(n,F,fill);xlabel ( 時(shí)間序列 n );ylabel(A./B);subplot(4,2,7);stem(n,G,fill);xlabel ( 時(shí)間序列 n );ylabel(A.B);運(yùn)行結(jié)果：. 學(xué)習(xí)幫手 .專業(yè)整理 .(2) 用 MATLAB實(shí)現(xiàn)以下序列。a）x(n)=0.8 n0n15實(shí)驗(yàn)程序：n=0:15;x=0.8.n;stem(n,x,fill); xlabel ( 時(shí)間序列 n );ylabe

4、l(x(n)=0.8n);b）x(n)=e (0.2+3j)n0 n 15實(shí)驗(yàn)程序：n=0:15;x=exp(0.2+3*j)*n);stem(n,x,fill); xlabel ( 時(shí)間序列 n );ylabel(x(n)=exp(0.2+3*j)*n);運(yùn)行結(jié)果：. 學(xué)習(xí)幫手 .專業(yè)整理 .a）的時(shí)間序列b ）的時(shí)間序列c）x(n)=3cos(0.125n+0.2 )+2sin(0.25 n+0.1 )0n15實(shí)驗(yàn)程序：n=0:1:15;x=3*cos(0.125*pi*n+0.2*pi)+2*sin(0.25*pi*n+0.1*pi);stem(n,x,fill); xlabel( 時(shí)間

5、序列 n );ylabel(x(n)=3*cos(0.125*pi*n+0.2*pi)+2*sin(0.25*pi*n+0.1*pi);運(yùn)行結(jié)果：d）將 c) 中的 x(n) 擴(kuò)展為以 16為周期的函數(shù) x16(n)=x(n+16),繪出四個(gè)周期實(shí)驗(yàn)程序：n=0:1:63;x=3*cos(0.125*pi*rem(n,16)+0.2*pi)+2*sin(0.25*pi*rem(n,16)+0.1*pi);stem(n,x,fill); xlabel ( 時(shí)間序列 n );ylabel(x16(n);e）將 c) 中的 x(n) 擴(kuò)展為以 10為周期的函數(shù) x10(n)=x(n+10),繪出四個(gè)

6、周期實(shí)驗(yàn)程序：n=0:1:39;x=3*cos(0.125*pi*rem(n,10)+0.2*pi)+2*sin(0.25*pi*rem(n,10)+0.1*pi);stem(n,x,fill); xlabel ( 時(shí)間序列 n );ylabel(x10(n);. 學(xué)習(xí)幫手 .專業(yè)整理 .運(yùn)行結(jié)果：d）的時(shí)間序列e ）的時(shí)間序列(3)x(n)=1,-1,3,5，產(chǎn)生并繪出下列序列的樣本。a）x1(n)=2x(n+2)-x(n-1)-2x(n)實(shí)驗(yàn)程序：n=0:3;x=1 -1 3 5;x1=circshift(x,0 -2);x2=circshift(x,0 1);x3=2*x1-x2-2*x

7、;stem(x3, fill); xlabel ( 時(shí)間序列 n );ylabel(x1(n)=2x(n+2)-x(n-1)-2x(n);b） x 2 (n)5nx (nk )k1實(shí)驗(yàn)程序：n=0:3;x=1 -1 3 5;x1=circshift(x,0 1);x2=circshift(x,0 2);x3=circshift(x,0 3);x4=circshift(x,0 4);x5=circshift(x,0 5);xn=1*x1+2*x2+3*x3+4*x4+5*x5;stem(xn,fill); xlabel( 時(shí)間序列 n );. 學(xué)習(xí)幫手 .專業(yè)整理 .ylabel(x2(n)=x

8、(n-1)+2x(n-2)+3x(n-3)+4x(n-4)+5x(n-5);運(yùn)行結(jié)果：a）的時(shí)間序列b）的時(shí)間序列(4) 繪出時(shí)間函數(shù)的圖形，對(duì) x軸、 y軸圖形上方均須加上適當(dāng)?shù)臉?biāo)注。a) x(t)=sin(2 t) 0 t 10sb) x(t)=cos(100t)sin( t)0 t 4s實(shí)驗(yàn)程序：clc;t1=0:0.001:10;t2=0:0.01:4;xa=sin(2*pi*t1);xb=cos(100*pi*t2).*sin(pi*t2);subplot(2,1,1);plot(t1,xa);xlabel (t);ylabel(x(t);title(x(t)=sin(2*pi*t)

9、;subplot(2,1,2);plot(t2,xb);xlabel (t);ylabel(x(t);title(x(t)=cos(100*pi*t2).*sin(pi*t2);運(yùn)行結(jié)果：. 學(xué)習(xí)幫手 .專業(yè)整理 .(5) 編寫(xiě)函數(shù) stepshift(n0,n1,n2) 實(shí)現(xiàn) u(n-n0),n1n0=0;stem(n,x,fill);xlable( 時(shí)間序列 n );ylable(u(n-n0);請(qǐng)輸入起點(diǎn)：2請(qǐng)輸入終點(diǎn)：8請(qǐng)輸入階躍位置：6運(yùn)行結(jié)果：. 學(xué)習(xí)幫手 .專業(yè)整理 .（ 5）運(yùn)行結(jié)果（ 6）運(yùn)行結(jié)果(6) 給一定因果系統(tǒng) H(z)(12 z -11)/(1 - 0.67z -1

10、0.9z -2 ) 求出并繪制 H(z) 的幅頻響應(yīng)與相頻響應(yīng)。實(shí)驗(yàn)程序：a=1 -0.67 0.9;b=1 sqrt(2) 1;h w=freqz(b,a);fp=20*log(abs(h);subplot(2,1,1);plot(w,fp);xlabel( 時(shí)間序列 t);ylabel( 幅頻特性 );xp=angle(h);subplot(2,1,2);plot(w,xp);xlabel( 時(shí)間序列 t);ylabel( 相頻特性 );運(yùn)行結(jié)果： ( 右上圖 ). 學(xué)習(xí)幫手 .專業(yè)整理 .常用典型序列單位采樣序列function x,n=impseq(n1,n2,n0)n=n1:n2;x

11、=(n-n0)=0;x,n=impseq(-2,8,2);stem(n,x);title(電信 1201)n0=-2;n=-10:10;nc=length(n);x=zeros(1,nc);fori=1:ncifn(i)=n0x(i)=1endendstem(n,x) ;title( 電信 1201 采樣序列第二種方法 )單位階躍序列functionx,n=stepseq(n1,n2,n0)n=n1:n2;x=(n-n0)=0;x,n=stepseq(-2,8,2);stem(n,x);title(電信 1201)實(shí)數(shù)指數(shù)n=0:10;x=0.9.n;stem(n,x);title( 電信 1

12、201 )復(fù)數(shù)指數(shù)序列n=-10:10;alpha=-0.1+0.3*j;. 學(xué)習(xí)幫手 .專業(yè)整理 .x=exp(alpha*n);real_x=real(x); image_x=imag(x);mag_x=abs(x); phase_x=angle(x);subplot(2,2,1);stem(n,real_x) ;title( 電信 1201 )subplot(2,2,2);stem(n,image_x);title( 電信 1201 )subplot(2,2,3);stem(n,mag_x);title( 電信 1201 )subplot(2,2,4);stem(n,phase_x);t

13、itle( 電信 1201 )正余弦n=0:10;x=3*cos(0.1*pi*n+pi/3);stem(n,x);title( 電信 1201 )例：求出下列波形x1(n)=2x(n-5)-3x(n+4)function y,n=sigadd(x1,n1,x2,n2)n=min(min(n1),min(n2) : max(max(n1),max(n2);y1=zeros(1,length(n);y2=y1;y1(find(n=min(n1)&(n=min(n2)&(n0);t2= -1:0.1:2;f2=t2.*exp(t2).*(t2=0)+exp(t2).*(t20);y,ny=conv

14、_m(f1,t1,f2,t2,p);subplot(3,1,1);stem(t1,f1);title(電信 1201)subplot(3,1,2);stem(t2,f2);title(電信 1201)subplot(3,1,3); stem(ny,y);title(電信 1201). 學(xué)習(xí)幫手 .專業(yè)整理 .實(shí)驗(yàn)二functionxk=dfs(xn,N)n=(0:1:N-1);k=n;WN=exp(-j*2*pi/N);nk=n*k;WNnk=WN.nk;xk=xn* WNnk;xn=0,1,2,3;N=4;xk=dfs(xn,N)IDFSfunctionxn=idfs(xk,N)n=(0:1

15、:N-1);k=n;WN=exp(-j*2*pi/N);nk=n*k;WNnk=WN.(-nk);xn=xk*WNnk/N;分析：因?yàn)閤(n) 是復(fù)指數(shù)，它滿足周期性，我們將在兩個(gè)周期中的401 個(gè)頻點(diǎn)上作計(jì)算來(lái)觀n=0:10;x=(0.9*exp(j*pi/3).n;k=-200:200;w=(pi/100)*k;X=x*(exp(-j*pi/100).(n*k);magX=abs(X);angX=angle(X);subplot(2,1,1);plot(w/pi,magX);subplot(2,1,2);. 學(xué)習(xí)幫手 .專業(yè)整理 .plot(w/pi,angX/pi);title( 電信

16、1201 )檢驗(yàn)頻移特性n=0:100;x=cos(pi*n/2);k=-100:100;w=(pi/100)*k;X=x*(exp(-j*pi/100).(n*k);y=exp(j*pi*n/4).*x;Y=y*(exp(-j*pi/100).(n*k);subplot(2,2,1);plot(w/pi,abs(X);axis(-1,1,0,60);title( 幅值電信 1201 );xlabel( 以 pi 為單位的頻率 );ylabel( 絕對(duì)值 X );subplot(2,2,2);plot(w/pi,angle(X)/pi);title( 相位電信 1201 );xlabel( 以

17、 pi 為單位的頻率 );ylabel( 絕對(duì)值 X );axis(-1,1,-1,1);subplot(2,2,3);plot(w/pi,abs(Y);title( 幅值電信 201 );xlabel( 以 pi為單位的頻率 );ylabel( 絕對(duì)值 Y );axis(-1,1,0,60);subplot(2,2,4);plot(w/pi,angle(Y)/pi);title( 相位電信 201 );xlabel( 以 pi為單位的頻率 );ylabel( 絕對(duì)值 Y );axis(-1,1,-1,1);b=1;a=1,-0.8;n=0:100;x=cos(0.05*pi*n);y=fil

18、ter(b,a,x);subplot(2,1,1);stem(n,x)xlabel(n );ylabel(x(n);title( 輸入序列電信 1201 );subplot(2,1,2);stem(n,y)xlabel(n );ylabel(y(n);title( 輸出序列電信 1201 );實(shí)驗(yàn)三例：已知信號(hào)由15Hz 幅值 0.5 的正弦信號(hào)和40Hz 幅值 2 的正弦信號(hào)組成，數(shù)據(jù)采樣頻率fs=100;N=128;n=0:N-1;t=n/fs;x=0.5*sin(2*pi*15*t)+2*sin(2*pi*40*t);y=fft(x,N);f=(0:length(y)-1)*fs/len

19、gth(y);. 學(xué)習(xí)幫手 .專業(yè)整理 .mag=abs(y);stem(f,mag);title(N=128 點(diǎn) 電信 1201 )已知帶有測(cè)量噪聲信號(hào)f1=50Hz,f2=1為均值為零、方差為1 的隨機(jī)信號(hào)，采樣頻率為1000Hz，t=0:0.001:0.6;x=sin(2*pi*50*t)+sin(2*pi*120*t);y=x+2*randn(1,length(t);Y=fft(y,512);P=Y.*conj(Y)/512;%求功率f=1000*(0:255)/512;subplot(2,1,1);plot(y);title( 電信 1201 )subplot(2,1,2);plot

20、(f,P(1:256);title( 電信 1201 )對(duì)信號(hào) 進(jìn)行 DFT，對(duì)其結(jié)果進(jìn)行IDFT，并將 IDFT 的結(jié)果和原信號(hào)進(jìn)行比較。fs=100; N=128; n=0:N-1; t=n/fs;x=sin(2*pi*40*t)+sin(2*pi*15*t);subplot(2,2,1)plot(t,x)title(original signal電信 1201 )y=fft(x,N);mag=abs(y);f=(0:length(y)-1)*fs/length(y);subplot(2,2,2)plot(f,mag)title(FFT to original signal電信 1201

21、)xifft=ifft(y);magx=real(xifft);ti=0:length(xifft)-1/fs;subplot(2,2,3)plot(ti,magx);title(signal from IFFT電信 1201 )yif=fft(xifft,N);mag=abs(yif);subplot(2,2,4)plot(f,mag)title(FFT to signal from IFFT電信 201 ). 學(xué)習(xí)幫手 .專業(yè)整理 .用函數(shù) conv 和 FFT計(jì)算同一序列的卷積，比較其計(jì)算時(shí)間。L=5000;N=L*2-1;n=1:L;x1=0.5*n;x2=2*n;t0=clock;y

22、c=conv(x1,x2);conv_time=etime(clock,t0)t0=clock;yf=ifft(fft(x1,N).*fft(x2,N);fft_time=etime(clock,t0)實(shí)驗(yàn)四例：用沖激響應(yīng)不變法設(shè)計(jì)Butterworth低通數(shù)字濾波器，通帶波紋小于1dB，阻帶在wp=0.2*pi; ws=0.3*pi;rp=1;rs=15;ts=0.01;Nn=128;Wp=wp/ts;Ws=ws/ts;N,Wn=buttord(Wp,Ws,rp,rs,s);z,p,k=buttap(N);Bap,Aap=zp2tf(z,p,k);b,a=lp2lp(Bap,Aap,Wn);

23、bz,az=impinvar(b,a,1/ts);freqz(bz,az,Nn,1/ts)title( 電信 1201 ). 學(xué)習(xí)幫手 .專業(yè)整理 .設(shè)計(jì)一個(gè) Butterworth高通數(shù)字濾波器，通帶邊界頻率為300Hz，阻帶邊界頻率為200Hz，通帶波紋小于1dB，阻帶衰減大于20dB，采樣頻率為 1000Hz。fs=1000;wp=300/(fs/2);ws=200/(fs/2);rp=1;rs=15; Nn=128;N,Wn=buttord(wp,ws,rp,rs);b,a=butter(N,Wn,high );freqz(b,a,Nn,fs)title( 電信 1201 )設(shè)計(jì)一個(gè)24 階 FIR 帶通濾波器，通帶頻率為wp=0.35,0.65;N=24;b=fir1(2*N,wp);freqz(b,1,512)title( 電信 1201 )f=0:0.002:1;