基于MATLAB 的時域信號采樣及頻譜分析

1、一、課程設(shè)計題目：基于 MATLAB 的時域信號采樣及頻譜分析二、基本要求： 掌握數(shù)字信號處理的基本概念、基本理論和基本方法； 學(xué)會 MATLAB 的使用，掌握 MATLAB 的程序設(shè)計方法； 學(xué)會用 MATLAB 對信號進(jìn)行分析和處理； 信號的各參數(shù)需由鍵盤輸入，輸入不同參數(shù)即可得不同的x(t) 和x(n)； 撰寫課程設(shè)計論文，用數(shù)字信號處理基本理論分析結(jié)果。三、設(shè)計方法與步驟： 畫出連續(xù)時間信號的時域波形及其幅頻特性曲線，其中幅度因子A444.128，衰減因子a222.144，模擬角頻率222.144； 對信號進(jìn)行采樣，得到采樣序列，其中T為采樣間隔，通過改變采樣頻率可改變T，畫出采樣頻率

2、分別為200Hz，500 Hz，1000 Hz時的采樣序列波形； 對不同采樣頻率下的采樣序列進(jìn)行頻譜分析，繪制其幅頻和相頻曲線，對比各頻率下采樣序列和的幅頻曲線有無差別，如有差別說明原因。 設(shè)系統(tǒng)單位抽樣響應(yīng)為，求解當(dāng)輸入為時的系統(tǒng)響應(yīng)，畫出, , 的時域波形及幅頻特性曲線，并利用結(jié)果驗證卷積定理的正確性（此內(nèi)容將參數(shù)設(shè)置為A1，a0.4，2.0734，T1）。 用FFT 對信號, , 進(jìn)行譜分析，觀察與中結(jié)果有無差別。 由采樣序列恢復(fù)出連續(xù)時間信號，畫出其時域波形，對比與原連續(xù)時間信號的時域波形，計算并記錄兩者最大誤差。四、詳細(xì)程序及仿真波形分析：1、連續(xù)時間信號x(t)及其200Hz/50

3、0Hz/1000Hz頻率抽樣信號函數(shù)x(n)%繪制信號x(n)的幅度譜和相位譜n=0:50; %定義序列的長度是50A=input('請輸入A的值 A:'); %設(shè)置信號的有關(guān)參數(shù)a=input('請輸入a的值 a:');w0=input('請輸入w0的值 w0:');T1=0.005;T2=0.002;T3=0.001;T0=0.001;x=A*exp(-a*n*T0).*sin(w0*n*T0); %pi 是MATLAB 定義的，信號乘可采用“.*”y1=A*exp(-a*n*T1).*sin(w0*n*T1); %pi 是MATLAB 定義

4、的，信號乘可采用“.*”y2=A*exp(-a*n*T2).*sin(w0*n*T2); %pi 是MATLAB 定義的，信號乘可采用“.*”y3=A*exp(-a*n*T3).*sin(w0*n*T3); %pi 是MATLAB 定義的，信號乘可采用“.*”close all %清除已經(jīng)繪制的x(n)圖形subplot(2,1,1);stem(n,x),grid on %繪制x(n)的圖形title('離散時間信號')subplot(2,1,2);plot(n,x),grid on title('連續(xù)時間信號')figure(2)subplot(3,1,1);

5、stem(n,y1),grid ontitle('200Hz理想采樣信號序列'); %設(shè)置結(jié)果圖形的標(biāo)題subplot(3,1,2);stem(n,y2),grid ontitle('500Hz連續(xù)時間信號')subplot(3,1,3);stem(n,y3),grid ontitle('1000Hz連續(xù)時間信號')k=-25:25;W=(pi/12.5)*k;w=W/pi; Y1=y1*(exp(-j*pi/12.5).(n'*k);figure(3)subplot(2,1,1);plot(w,abs(Y1);grid,xlabel(&

6、#39;w'),ylabel('幅度');title('200Hz理想采樣信號序列的幅度譜');axis(-2 2 0 1000);subplot(2,1,2);plot(w,angle(Y1);grid,xlabel('w'),ylabel('幅角');title ('200Hz理想采樣信號序列的相位譜')Y2=y2*(exp(-j*pi/12.5).(n'*k);figure(4)subplot(2,1,1);plot(w,abs(Y2);grid,xlabel('w'),yla

7、bel('幅度');title('500Hz理想采樣信號序列的幅度譜');axis(-2 2 0 1000);subplot(2,1,2);plot(w,angle(Y2);grid,xlabel('w'),ylabel('幅角');title ('500Hz理想采樣信號序列的相位譜')Y3=y3*(exp(-j*pi/12.5).(n'*k);figure(5)subplot(2,1,1);plot(w,abs(Y3);grid,xlabel('w'),ylabel('幅度'

8、;);title('1000Hz理想采樣信號序列的幅度譜');axis(-2 2 0 1000);subplot(2,1,2);plot(w,angle(Y3);grid,xlabel('w'),ylabel('幅角');title ('1000Hz理想采樣信號序列的相位譜')分析：采樣頻率為1000Hz時沒有失真，500Hz時有橫線，產(chǎn)生失真，200Hz時橫線加長，失真增大。說明采樣頻率越大失真越小。2、設(shè)系統(tǒng)單位抽樣響應(yīng)為，求解當(dāng)輸入為時的系統(tǒng)響應(yīng)，畫出, , 的時域波形及幅頻特性曲線，并利用結(jié)果驗證卷積定理的正確性（此內(nèi)容將

9、參數(shù)設(shè)置為A1，a0.4，2.0734，T1）。n=1:50; %定義序列的長度是50hb=zeros(1,50); %注意：MATLAB 中數(shù)組下標(biāo)從1 開始hb(1)=1; hb(2)=1; hb(3)=1; hb(4)=1;hb(5)=1;close all; subplot(3,1,1);stem(hb);title('系統(tǒng)hbn');m=1:50; T=1;%定義序列的長度是和采樣率A=1; a=0.4;T=1;w0=2.0734;x=A*exp(-a*m*T).*sin(w0*m*T); %pi 是MATLAB 定義的，信號乘可采用“.*”subplot(3,1,2

10、);stem(x);title('輸入信號xn');y=conv(x,hb);subplot(3,1,3);stem(y);title('輸出信號yn');figure(2)subplot(3,1,1);plot(n,hb),grid ontitle('矩形序列時域波形');subplot(3,1,2);plot(m,x),grid ontitle('輸入信號xn時域波形');subplot(3,1,3);plot(m,y),grid ontitle('輸出信號yn時域波形');分析：在數(shù)字信號處理中經(jīng)常要進(jìn)行卷積

11、運算，MATLAB中有一個內(nèi)部函數(shù)conv可以計算兩個有限長序列的卷積，該函數(shù)計算的兩個序列都是從n=0開始3、用FFT 對信號, , 進(jìn)行譜分析，觀察與中結(jié)果有無差別。n=1:50; %定義序列的長度是50hb=zeros(1,50); %注意：MATLAB 中數(shù)組下標(biāo)從1 開始hb(1)=1; hb(2)=1; hb(3)=1; hb(4)=1;hb(5)=1;close all; subplot(3,1,1);m=1:50; T=1;%定義序列的長度是和采樣率A=1; a=0.4;T=1;w0=2.0734;x=A*exp(-a*m*T).*sin(w0*m*T); %pi 是MATLA

12、B 定義的，信號乘可采用“.*”y=conv(x,hb);subplot(3,1,1); plot(n,abs(fft(hb)title('h(n)的FFT')subplot(3,1,2); plot(n,abs(fft(x)title('x(n)的FFT')subplot(3,1,3); plot(abs(fft(y)title('y(n)的FFT')分析：MATLAB中，計算矢量x的DFT及其逆變換的函數(shù)分別為fft和ifft，這兩個函數(shù)采用了混合算法，當(dāng)N為質(zhì)數(shù)時，采用的是原始的DFT算法。函數(shù)是用機器語言編寫的，執(zhí)行速度混快。N點的FFT

13、調(diào)用形式為fft(x,N)。如果x的長度小于N，則補零使其成為N點序列；如果省略N點，即以fft(x)形式調(diào)用，則按矢量x的長度進(jìn)行計算；如果x表示一個矩陣，則調(diào)用后計算出每列的N點的FFT。4、由采樣序列恢復(fù)出連續(xù)時間信號，畫出其時域波形，對比與原連續(xù)時間信號的時域波形，計算并記錄兩者最大誤差。A=input('please input the A:');%設(shè)置信號的有關(guān)參數(shù)a=input('please input the a:');W0=input('please input the W0:');fs=input('please i

14、nput the fs:');%采樣率n=0:49;%定義序列的長度T=1/fs;t0=10/a;Dt=1/(5*a);t=0:Dt:t0;xa=A*exp(-a*t).*sin(W0*t);K1=50;k1=0:1:K1;W1max=2*pi*500;W1=W1max*k1/K1;w1=W1/pi;Xa=xa*exp(-j*t'*W1);x=A*exp(-a*n*T).*sin(W0*n*T);figure(1);subplot(4,1,1);plot(t*1000,xa);title('連續(xù)時間信號x(t)');axis(0 t0*1000 0 200);g

15、rid,xlabel('t:毫秒'),ylabel('x(t)');subplot(4,1,2);plot(w1,abs(Xa);title('連續(xù)時間信號頻譜Xa(w1)');subplot(4,1,3);stem(x)%繪制x(n)圖形grid,xlabel('n'),ylabel('x(n)');title('采樣序列x(n)');x1=spline(n*T,x,t);grid,xlabel('t:毫秒'),ylabel('x(t)');subplot(4,1,

16、4);plot(t*1000,x1);axis(0 t0*1000 0 200);title('由x(n)恢復(fù)x1(t)');grid,xlabel('t:毫秒'),ylabel('x1(t)');errror=max(abs(x1-xa);k2=-25:25;W2=(pi/12.5)*k2;w2=W2/pi;X=x*(exp(-j*pi/12.5).(n'*k2);%序列的傅里葉變換函數(shù)figure(2);subplot(2,1,1);plot(w2,abs(X);grid,xlabel('w2'),ylabel(

17、9;幅度');title('輸入信號幅度譜')axis(-2 2 0 1000);subplot(2,1,2);plot(w2,angle(X);grid,xlabel('w2'),ylabel('幅角');title('輸入信號相位譜');axis(-2 2 -5 5);分析：恢復(fù)曲線與原信號曲線相同，說明恢復(fù)誤差很小，如果采樣頻率減小，誤差增大，采樣頻率增大，則恢復(fù)誤差更小。采樣頻率應(yīng)遵循乃奎斯特定理。五、調(diào)試總結(jié)：Axis函數(shù)設(shè)置圖形坐標(biāo)。Conv函數(shù)進(jìn)行卷積運算。Stem用還繪制離散圖譜，plot用來繪制連續(xù)的函數(shù)

18、圖形。fft進(jìn)行fft運算。subplot(3,1,1); plot(n,abs(fft(hb)title('h(n)的FFT')subplot(3,1,2); plot(n,abs(fft(x)title('x(n)的FFT')subplot(3,1,3); plot(abs(fft(y)title('y(n)的FFT') 最后一行plot(abs(fft(y)如果寫成plot(n，abs(fft(y)會出現(xiàn)長度不夠，沒有圖形輸出。六、設(shè)計總結(jié)：在老師的幫助下我順利的完成了這個課程設(shè)計，通過這次數(shù)字信號處理課程設(shè)計，讓我了解了關(guān)于MATLAB軟件在數(shù)字信號處理方面的應(yīng)用，又一次學(xué)習(xí)了MATLAB軟件的使用和程序的設(shè)計，MATLAB的仿真使我更加深入的了解了數(shù)字處理的過程 ，對我對數(shù)字信號處理的理解加深了一步。MATLAB擁有強大的數(shù)據(jù)仿真能力，在生產(chǎn)和研究中起著非常大的作。MATLAB語言是一種廣泛應(yīng)用于工程計算及數(shù)值分析 領(lǐng)域的新型高級語言，Matlab功能強大、簡單易學(xué)、編程效率高，深受廣大科技工作者的歡迎。特別是Matlab還具有信號分析工具箱，不需具備很強的編程能力，就可以很方便地進(jìn)行信號分析、處理和設(shè)計。因此，選擇用Matlab進(jìn)行課程設(shè)計。在這過程中我遇到了所多的難題，通過與老師的交流和學(xué)