《數字信號處理》實驗指導書
編 寫： 張大鵬
適用專業：電子信息工程
專 業： 電子信息工程
姓 名： 江凱倫
學 號： 198840009
指導教師： 張大鵬
完成日期： 2017年 12月 3日
信息與電氣工程學院
目 錄
概述 3
實驗一 Matlab基本操作 7
實驗二 序列 11
實驗三 Z變換 18
實驗四 FFT 21
概述
Matlab是一種功能強大的工具軟件，在數字信號處理領域有著廣泛的應用，因此數字信號處理課程的實驗主要是在Matlab環境下進行的。
MATLAB語言是一種廣泛應用于工程計算及數值分析領域的新型高級語言，自1984年由美國 MathWorks 公司推向市場以來，歷經十多年的發展與競爭，現已成為國際公認的最優秀的工程應用開發環境。MATLAB功能強大、簡單易學、編程效率高，深受廣大科技工作者的歡迎。
MATLAB特點:
1. 數值計算和符號計算功能
MATLAB的數值計算功能包括：矩陣運算、多項式和有理分式運算、數據統計分析、數值積分、優化處理等。符號計算將得到問題的解析解。
2.MATLAB語言
MATLAB除了命令行的交互式操作以外，還可以程序方式工作。使用MATLAB可以很容易地實現C或FORTRAN語言的幾乎全部功能，包括Windows圖形用戶界面的設計。
3.圖形功能
MATLAB提供了兩個層次的圖形命令：一種是對圖形句柄進行的低級圖形命令，另一種是建立在低級圖形命令之上的高級圖形命令。利用MATLAB的高級圖形命令可以輕而易舉地繪制二維、三維乃至四維圖形，并可進行圖形和坐標的標識、視角和光照設計、色彩精細控制等等。
4.應用工具箱
基本部分和各種可選的工具箱。
基本部分中有數百個內部函數。
其工具箱分為兩大類：功能性工具箱和學科性工具箱。功能性工具箱主要用來擴充其符號計算功能、可視建模仿真功能及文字處理功能等。學科性工具箱專業性比較強，如控制系統工具箱、信號處理工具箱、神經網絡工具箱、最優化工具箱、金融工具箱等，用戶可以直接利用這些工具箱進行相關領域的科學研究。 MATLAB集成環境
MATLAB 7.X是一個高度集成的語言環境，在該環境下既可以進行交互式的操作，又可以編寫程序、運行程序并跟蹤調試程序。
1. MATLAB的啟動
啟動MATLAB有兩種常見方法：
(1)通過“開始”按鈕，選擇“程序”菜單項，然后打開“MATLAB”菜單中的“MATLAB”程序，就可啟動MATLAB系統
(2) 利用Windows 建立快捷方式的功能，將MATLAB程序以快捷方式放在桌面上。只要在桌面上雙擊該圖標即可啟動MATLAB
2. MATLAB命令窗口
(1) 命令窗口的菜單欄
菜單欄共包含File、Edit、Window和Help四項。
File菜單項：
New命令：用于建立M文件、圖形窗口。
Open命令：打開一個已經建立的M文件。
Run Script命令：執行一個命令文件。
Load Workspace命令：將變量裝入當前空間。
Save Workspace As命令：把當前工作空間的所有變量用后綴為.mat的文件保存起來。
Show Workspace命令：打開變量瀏覽器。
Set Path命令：打開MATLAB的路徑瀏覽器。
Preferences命令：打開命令窗口的顯示格式。
Print Setup命令：設置打印機的參數。
Print命令：打印和設置一些打印參數。
Print Selection命令：打印選中的內容。
Exit MATLAB命令：退出MATLAB系統。
Edit菜單項：
Undo、Cut、Copy和Paste等命令：分別用于撤銷上一次操作、剪切、復制和粘貼。
Clear命令：刪除內容。
Select All命令：用于選定所有文本內容。
Clear Session命令：清除命令編輯區的全部內容，但并不刪除工作空間中的變量。
Help菜單項:
Help Window命令：打開MATLAB的幫助窗口。
Help Tips命令：打開幫助窗口，并首先顯示MATLAB的幫助系統的分類和使用方法。
Help Desk(HTML)命令：打開系統WWW瀏覽器，并顯示MATLAB的幫助桌面。 Examples and Demos命令：可以通過演示MATLAB提供的例子來熟悉相關部分的用法。
About MATLAB命令：打開關于MATLAB的版本和版權等信息。
Subsribe命令：打開機器上的WWW瀏覽器，用戶可過填寫相關的表格來獲得MathWorks公司的產品。
MATLAB運算量
1．變量和賦值語句
MATLAB賦值語句有兩種形式：
(1) 變量=表達式
(2) 表達式
其中“表達式”是用運算符將有關運算量連接起來的式子，其結果是一個矩陣。［注］ 第二種語句形式下，將表達式的值賦給MATLAB的永久變量ans。如果在語句的最后加分號，那么，MATLAB僅僅執行賦值操作，不再顯示運算的結果。在一條語句中，如果表達式太復雜，一行寫不下，可以加上三個小黑點(續行符)并按下回車鍵，然后接下去再寫。例如 s=1-1/2+1/3-1/4+1/5-1/6+1/7-…- 1/8+1/9-1/10+1/11-1/12。
2．MATLAB表達式
算術表達式。
運算符有：＋(加)、－(減)、*(乘)、/(右除)、\(左除)、^(乘方)
對于矩陣來說，左除和右除表示兩種不同的除數矩陣和被除數矩陣的關系。關系表達式。
運算符有：<(小于)、<=(小于或等于)、>(大于)、>=(大于或等于)、==(等于)、～=(不等于)
邏輯表達式。
運算符有：&(與)、|(或)和～(非)
運算法則：
(1)在邏輯運算中，確認非零元素為真，用1表示，零元素為假，用0表示。
(2)參與邏輯運算的可以是兩個標量、兩個同維矩陣或參與邏輯運算的元素一個為標量，另一個為矩陣。
(3)在算術、關系、邏輯運算中，算術運算優先級最高，邏輯運算優先級最低。
實驗1 Matlab基本操作
1. 實驗目的
初步掌握Matlab中的矩陣生成方法和繪圖方法。
2. 實驗內容和步驟
（1） 矩陣的生成與察看
在工作區鍵入以下命令
>>x=[1 2 3 4;1 2 3 4;1 2 3 4;1 2 3 4];
可以看到工作區的變量窗口出現了x，雙擊察看可以看到這是一個4x4的矩陣。
>>x=[1 2 3 4;1 2 3 4;1 2 3 4;1 2 3 4]
可以在工作區看到x的內容
x =
1 2 3 4
1 2 3 4
1 2 3 4
1 2 3 4
當鍵入以下命令的的時候
>> x([2,3],:)
>> x(:,[1,2])
將會出現指定行或列的信息。
請自己輸入一個4x4可逆矩陣，并且用inv()函數來求解逆陣，最后將逆陣第2，3列的信息顯示出來。
（2） 幫助的獲取
Matlab的幫助功能非常強大，雖然其內容都是英文，只要使用者有一定的英文閱讀能力，獲取信息的可能性還是很大的。
實驗的第3部分內容是畫圖，在畫圖之前鍵入
>>help plot
察看其基本用法。
（3） 畫圖命令
首先生成一個一維矩陣x
>> x = 0:0.01:10
然后產生一個一維矩陣y
>>y = sin(x)
畫出相應的圖形。生成z=cos(x)矩陣，并把兩個函數畫在一張圖上。
3．實驗結果
4. 實驗分析
MATLAB具有強大的畫圖功能，可以將各種信號表示出來
5．思考題
Matlab是否對大小寫敏感，請舉例說明。
敏感
實驗2 序列
1． 實驗目的 掌握Matlab中基本信號的生成方法和序列的運算方法
2． 實驗內容和步驟
(1) 實際應用中，常把一個列向量作為一路信號，先構造一個行向量
>> x=[4 3 7 -9 1]
然后轉置
>> x=x'
利用轉置后的x生成一個三通道信號
>>y=[x 2*x x/pi]
(2) 許多不同的工具箱都可以產生信號波形，其中大部分都要求時間向量作為參數，如果以1000Hz抽樣頻率產生的波形，適宜的時間向量如下
>>t=(0:0.001:1)'
(3) 基本序列生成
δ(n) ：
>> n=1；
>> x=[1 zeros(1,n-1)]
u(n) ：
>>n= 一個指定長度(計算機無法處理無窮的數據)
>>x=ones(1,N)
anu(n) ：
>>n= 一個指定長度(計算機無法處理無窮的數據)
>> b=[0:n-1]
>> x=a.^b
exp((σ+jω)*n)：
>> u=1
>> w=2
>> x=exp((u+j*w)*b)
(4) 基本周期波形
方波：
正弦波：
鋸齒波：
(5) sinc波形
(6) 序列的運算
信號加
>>x=x1+x2
信號乘
>>x=x1.*x2
比例
>>y=k*x
折疊
>>y=fliplr(x)
抽樣和
>>y=sum(x(n1:n2))
抽樣積
>>y=prod(x(n1:n2))
信號能量
>>Ex=sum(abs(x).^2)
信號功率
>>Px=(sum(abs(x).^2)/n
3． 實驗結果
4．思考題
>>t=(0:0.001:1)’;
>>x1=sin(2*pi*50*t)+2*sin(2*pi*120*t);
>>x=x1+0.5*randn(size(t));
（1） 求x(1)~x(50)的抽樣和抽樣積。
（2） 求x的能量和功率。
（3） 把結果用圖形表示出來。
實驗3 Z變換
1． 實驗目的：掌握使用Matlab處理Z變換、Z反變換的方法。
2． 實驗原理
z變換定義
利用差分方程可求離散系統的結構及瞬態解。為了分析系統的另外一些重要特性，如穩定性和頻率響應等，需要研究離散時間系統的z變換（類似于模擬系統的拉氏變換），它是分析離散系統和離散信號的重要工具。
一個離散序列x（n）的Z變換定義為
其中z為復變量，是一個以實部為橫坐標，虛部為縱坐標構成的平面上的變量，這個平面也稱z平面。
常用Z[x(n)]表示對序列x(n)的z變換，即
這種變換也稱為雙邊z變換，與此相應還有單邊z變換，單邊z變換只是對單邊序列（n>=0部分）進行變換的z變換，其定義為
可以把單邊z變換看成是雙邊z變換的一種特例，即因果序列情況下的雙邊z變換。
Z變換收斂域
一般，序列的Z變換并不一定對任何z值都收斂，z平面上使上述級數收斂的區域稱為“收斂域”。我們知道，級數一致收斂的條件是絕對值可和，因此z平面的收斂域應滿足
因為對于實數序列， ，因此，|z|值在一定范圍內才能滿足絕對可和條件，這個范圍一般表示為
Rx-〈|z|〈Rx+
這就是收斂域，一個以Rx-和Rx+為半徑的兩個圓所圍成的環形區域，Rx-和Rx+稱為收斂半徑，Rx-和Rx+的大小，即收斂域的位置與具體序列有關，特殊情況為Rx-或Rx+等于0，這時圓環變成圓或空心圓。
Z變換收斂域的特點：
· 收斂域是一個圓環，有時可向內收縮導原點，有時可向外擴展到∞，只有x(n)=δ(n)的收斂域是整個z平面；
· 在收斂域內沒有極點，x（z）在收斂域內每一點上都是解析函數。
Z變換表示法：
· 級數形式；
· 解析表達式（注意只表示收斂域上的函數，同時要注明收斂域）。
3． 實驗內容和步驟
(1) 已知一個Z變換的表達式
請使用residuez函數來計算反變換后的序列。
(2) 已知一個因果線性移不變系統
y(n) = 0.81y(n-2) + x(n) – x(n-2) 求
a.H(z)
提示：直接筆算即可
在Matlab中，一個分式可以以下形式表示：
H(z) =B(z)/A(z)
對于已知條件B(z)=1 – z^-2 ， A(z)=1 – 0.81z^-2。
在Matlab中將B按照z^-1的升冪表示成向量形式為[1 0 -1 ]，同理A表示成[1 0 -0.81]。
b.沖擊響應
提示： 對于數字系統的沖擊響應，可使用函數dimpulse。
為與后面的階躍響應呼應，可把二者畫在一幅圖上，可使用函數subplot。
c.階躍響應
提示： 對于數字系統的階躍響應，可使用函數dimpulse。
d.H(ejw) 的幅頻相頻曲線
提示：可使用函數freqz，時間向量按照下式設定：
t=[0:1:500]*pi/500;
(3) 解差分方程
y(n) = 1/3[x(n) + x(n-1) + x(n-2)] + 0.95y(n-1) – 0.9025y(n-2) n>=0
其中 x(n) = cos(pi*n/3)，y(-1) = -2，y(-2) = -3，x(-1) = 1，x(-2) = 1。程序如下：
b=[1,1,1]/3;
a=[1,-0.95,0.9025];
Y=[-2,-3];
X=[1,1];
xic=filtic(b,a,Y,X);
bxplus=[1,-0.5];
axplus=[1,-1,1];
ayplus=conv(a,axplus);
byplus=conv(b,bxplus) + conv(xic,axplus)
[R,p,C]=residuez(byplus,zxplus);
Mp=abs(p),A(p)=angle(p)/pi
n=[0:50];
x=cos(pi*n/3);
y=filter(b,a,x,xic);
plot(n,y)
請詳細分析程序，在每行加上注解。
4． 實驗結果
Z變換
Z反變換
實驗4 FFT
1. 實驗目的
掌握在Matlab下進行付立葉變換的方法。
2．實驗原理
（一）基—2按時間抽取FFT算法
對于有限長離散數字信號{x[n]}，0 ￡ n ￡ N-1，其離散譜{x[k]}可以由離散付氏變換（DFT）求得。DFT的定義為
可以方便的把它改寫為如下形式：
不難看出，WN是周期性的，且周期為N，即
WN的周期性是DFT的關鍵性質之一。為了強調起見，常用表達式WN取代W以便明確其周期是N。
由DFT的定義可以看出，在x[n]為復數序列的情況下，完全直接運算N點DFT需要（N-1）2次復數乘法和N（N-1）次加法。因此，對于一些相當大的N值（如1024）來說，直接計算它的DFT所作的計算量是很大的。FFT的基本思想在于，將原有的N點序列序列分成兩個較短的序列，這些序列的DFT可以很簡單的組合起來得到原序列的DFT。例如，若N為偶數，將原有的N點序列分成兩個（N/2）點序列，那么計算N點DFT將只需要約[(N/2)2 ·2]=N2/2次復數乘法。即比直接計算少作一半乘法。因子（N/2）2表示直接計算（N/2）點DFT所需要的乘法次數，而乘數2代表必須完成兩個DFT。上述處理方法可以反復使用，即（N/2）點的DFT計算也可以化成兩個（N/4）點的DFT（假定N/2為偶數），從而又少作一半的乘法。這樣一級一級的劃分下去一直到最后就劃分成兩點的FFT運算的情況。比如，一個N = 8點的FFT運算按照這種方法來計算FFT可以用下面的流程圖來表示：
關于蝶形結運算的具體原理及其推導可以參照講義，在此就不再贅述。按頻率抽取的FFT的原理也可查閱相關資料，這里就不再推導了
3．實驗內容和步驟
(1) 設定一信號并計算它的FFT變換
>> t = (0:1/99:1);
>>x = sin(2*pi*15*t) + sin(2*pi*40*t);
>>y = fft(x);
>>m=abs(y);
>>p=unwrap(angle(y));
>>f=(0:length(y)-1) * 99/length(y);
>>subplot(1,2,1);
>>plot(f,m);
>>set(gca,’XTick’,[15 40 60 85]);
>>subplot(1,2,2);
>>plot(f, p*180/pi);
>>set(gca,’XTick’,[15 40 60 85]);
(2) x(t) = 2sin(4πt) + 5cos(8πt)，求x(t)的N點幅度譜和相位譜
(3) 比較原序列和經過FFT和IFFT后的序列。
4．實驗結果
5．實驗分析
在圖1中明顯能看到整個頻譜圖關于Nyquist頻率對稱，不過Nyquist頻率右邊的譜圖實際上是負頻部分，沒有意義，從中看到，該信號包含兩個頻率15Hz和40Hz。由于使用的采樣頻率fs=100Hz，所以Nyquist頻率為50Hz。
|
