青島科技大學(xué)《數(shù)字信號處理》課程設(shè)計(一)

1、青 島 科 技 大 學(xué)數(shù)字信號處理課程設(shè)計報告數(shù)字信號分析及數(shù)字濾波器設(shè)計題 目 _張淑軍指導(dǎo)教師_學(xué)生姓名_學(xué)生學(xué)號_院（部）_專業(yè)班_2014_年 _12_月 _25_日1. 目的與要求1.進一步鞏固數(shù)字信號處理中的基本原理與方法，提高分析、解決實際問題的能力。 2.熟練掌握一門計算機語言，進行數(shù)字信號處理應(yīng)用的開發(fā)設(shè)計，訓(xùn)練基本技能，如查閱設(shè)計資料和手冊、程序的設(shè)計、調(diào)試等。3選擇的題目為（一）數(shù)字信號分析及數(shù)字濾波器設(shè)計1.用以下方式產(chǎn)生三個不同頻段的信號：（1）自己錄制一段正常的語音文件；（2）錄制一段環(huán)境噪聲文件；（3）利用MATLAB產(chǎn)生一個不同于以上頻段的信號。2.對上述三個信

2、號，進行頻譜分析，畫出三路信號的時域波形和頻譜圖，對進行對比分析。3.根據(jù)三路信號的頻譜特點得到性能指標(biāo)，由性能指標(biāo)設(shè)計三個濾波器，并畫出各濾波器的頻域響應(yīng)。4.將三路信號疊加為一路信號。5.用自己設(shè)計的濾波器對合成的信號進行濾波，分析得到信號的頻譜，并畫出濾波后信號的時域波形和頻譜。2.主要技術(shù)和原理2.1語音的輸入、處理以及信號的疊加 Matlab中語音信號的輸入為x,fs,bits=wavread('C:UserslenovoDesktop原始語音.wav',50000 60000);其中wavread函數(shù)讀取Microsoft的wav格式音頻文件，x為自己隨意定義的一個

3、輸出函數(shù)，fs 是采樣頻率 bits 表示每個樣點的位數(shù)，50000 60000表示讀取從50000到60000個點到x中。X_voice=x(:,1)因為音頻信號都有兩個聲道，所以提取出來一個聲道進行分析Plot（x）為畫出當(dāng)前的波形，freqz（x）是求離散系統(tǒng)頻響特性的函數(shù)。隨機產(chǎn)生的信號（2）信號的疊加本次實驗需要將正常語音信號與噪音信號還有matlab產(chǎn)生的一個信號進行疊加，生成一個新的信號，代碼如下f1=22000; f2=25000; fs=50000; t=1/fs; length_x=length(x); input_x=zeros(1,length_x);for i=1:1

4、:length_x input_x(i)=sin(2*pi*f1*i*t)+sin(2*pi*f2*i*t);input_x=input_x'end 正常語音信號的代碼：x,fs,bits=wavread('C:UserslenovoDesktop原始語音.wav',50000 60000);X_voice=x(:,1);input_x=X_voice;噪聲語音的代碼：x1,fs1,bits1=wavread('C:UserslenovoDesktopÔëÒô1.wav',50000 60000);X3_voice

5、=x1(:,1);input_x1=X3_voice;但是信號疊加有一個非常重要的條件就是點數(shù)必須相同，不僅點數(shù)相同而且點數(shù)矩陣得是同一類型，在原始語音和噪聲語音的點數(shù)10001個，但他們產(chǎn)生的是10001行1列矩陣，而matlab產(chǎn)生的正弦信號卻是1行10001列矩陣，結(jié)果也不能疊加，這就需要將1行10001列矩陣變成10001行1列矩陣，input_x=input_x'反轉(zhuǎn)語句就是講input_x的行和列互換。疊加語句：for i=1:1:length_xy1=X_voice+X3_voice+input_x;end2.2FIR濾波器設(shè)計原理時間窗口設(shè)計法是從單位脈沖響應(yīng)序列著手，

6、使h(n)逼近理想的單位脈沖響應(yīng)序列hd(n)。我們知道hd(n)可以從理想頻響通過付氏反變換獲得但一般來說，理想頻響是分段恒定，在邊界頻率處有突變點，所以，這樣得到的理想單位脈沖響應(yīng)hd(n)往往都是無限長序列，而且是非因果的。但FIR的h(n)是有限長的，問題是怎樣用一個有限長的序列去近似無限長的。最簡單的辦法是直接截取一段hd(n)代替 h(n) 。這種截取可以形象地想象為h(n)是通過一個“窗口”所看到的一段hd(n)，因此 ，h(n)也可表達為h(n)和一個“窗函數(shù)”的乘積，即 h(n)=w(n) hd(n) 在這里窗口函數(shù)就是矩形脈沖函數(shù)RN（n），當(dāng)然以后我們還可看到，為了改善設(shè)

7、計濾波器的特性，窗函數(shù)還可以有其它的形式，相當(dāng)于在矩形窗內(nèi)對hd(n)作一定的加權(quán)處理。（2）本次實驗只采用了矩形窗設(shè)計方法，由于hd(n)是無限長的，故要用一個有限長的“窗函數(shù)”序列w(n)將hd(n)加以截斷，窗的點數(shù)是N點。截斷后的序列為hd(n) h(n)=w(n) hd(n)， 0<=n<=N-1窗的點數(shù)N及窗的形狀是兩個極其重要的參數(shù)。加窗后實際的頻率響應(yīng) =DTFTh(n)=DTFT=(3)FIR帶通、高通、帶阻、多通帶多阻帶濾波器的設(shè)計（1）N長的FIR全通濾波器h（n）的函數(shù)表達當(dāng)，則低通濾波器變?yōu)槿V波器： （2）N長的FIR高通濾波器h（n）的函數(shù)表達式由高

8、通的頻譜結(jié)構(gòu)可知，高通濾波器的頻譜=全通濾波器的頻譜-低通濾波器的頻譜。（3）N長的FIR帶通濾波器h（n）的函數(shù)表達式帶通濾波器的頻譜=低通濾波器的頻譜-低通濾波器的頻譜。（3）N長的FIR帶通濾波器h（n）的函數(shù)表達式帶阻濾波器的頻譜=全通濾波器的頻譜-帶通濾波器的頻譜。 （4）N長的FIR帶阻濾波器h（n）的函數(shù)表達式 帶阻濾波器的頻譜=全通濾波器的頻譜-帶通濾波器的頻譜 （5）N長的FIR多通帶多帶阻濾波器h（n）的函數(shù)表達式 多通帶多阻帶濾波器的頻譜=全通濾波器的頻譜-帶通濾波器的頻譜-帶通濾波器的頻譜。3.總體方案設(shè)計3.1開發(fā)環(huán)境和工具 Pc window8.1 matlab20

9、12a3.2具體框架Matlab產(chǎn)生一個信號噪聲語音的錄入原始信號的錄入繪制時域頻域波形繪制時域頻域波形繪制時域頻域波形三路信號疊加成一路信號繪制時域頻域波形設(shè)計多通多阻濾波器設(shè)計高通濾波器設(shè)計低通濾波器與y1相卷積進行濾波與y1相卷積進行濾波與y1相卷積進行濾波語音恢復(fù)成wav格式語音恢復(fù)成wav格式語音恢復(fù)成wav格式3.3設(shè)計過程 （1）信號的獲取采集（原始語音，噪聲，matlab產(chǎn)生的信號） 錄取一段語音信號，fs=44100赫茲，共抽樣10001個點錄取一段噪聲語音，本次噪聲語音就是用電腦自己錄的雜音當(dāng)做噪聲，且fs=44100，也是抽取了10001個點用matlab產(chǎn)生了一個正弦信

10、號抽樣頻率為50000，f1和f2分別為22000和25000赫茲。 （2）對做好的信號要進行處理，用wavread函數(shù)將語音信號讀進來，取他們的一列，對于正弦信號要進行轉(zhuǎn)置input_x=input_x'以便于后面對信號的疊加做準(zhǔn)備。 （3）繪制語音信號的時域波形的頻域波形 （4）根據(jù)（3）設(shè)計出來的頻譜特性進行分析，根據(jù)分析的結(jié)果設(shè)計濾波器，在畫出各個濾波器的頻譜以及時域波形。 （5）將做好的三路信號疊加成一路信號，在對疊加好的信號畫頻譜。 （6）將設(shè)計好的濾波器分別與此信號卷積進行濾波，在畫出濾波后的信號的頻譜。 （7）最后將語音恢復(fù)，輸出到指定文件夾下。4.實驗結(jié)果與分析4.1

11、這部分是實現(xiàn)對原始語音的處理具體代碼如下：x,fs,bits=wavread('C:UserslenovoDesktop原始話音.wav',50000 60000); %對原始話音讀取到x中X_voice=x(:,1);%提取聲道1input_x=X_voice;%將聲道1賦給input_xfigure(1);plot(input_x);對input_x畫時域波形title('原始話音波形');figure(1);freqz(input_x);%對input_x畫頻譜波形title('原始話音頻譜');X1_voice=x(:,2);%提取聲道2

12、input_x=X1_voice;%將聲道2賦給input_xfigure(3);plot(input_x);title('原始話音波形');figure(4);freqz(input_x);title('原始話音頻譜');y=X_voice+X1_voice;&將兩個信號疊加為原始語音信號figure(5);plot(y);title('原始話音波形');figure(6);freqz(y);title('原始話音頻譜');主要用到的函數(shù)wavread是對語音信號的讀取，但只能是wav格式，plot(input_x);是

13、畫圖函數(shù)，將input_x畫出來，freqz(input_x);是求input_x的頻譜并將其畫出。實驗結(jié)果圖如下：上圖為聲道1的波形和頻譜，只是提取了10001個點進行分析，此頻譜特點是圖上在f=2800,11250,21400赫茲處出現(xiàn)了尖峰。這里只對提取的一個聲道分析，代表了此信號。4.2實現(xiàn)對噪音信號的處理模塊及其實現(xiàn)代碼：x1,fs1,bits1=wavread('C:UserslenovoDesktopÔëÒô1.wav',50000 60000);X3_voice=x1(:,1);input_x1=X3_voice;figu

14、re(9);plot(input_x1);title('原始噪音波形');figure(10);freqz(input_x1);title('原始噪音頻譜圖');此模塊與正常語音的的處理大同小異，只是在噪音里我只提取了一個聲道，用到的函數(shù)與處理正常語音用到的函數(shù)一樣。實驗結(jié)果圖如下以下兩個圖一個是噪音頻譜，另一個是噪音波形，我也是對噪音提取了10001個點進行分析，從頻譜上看，在f=6700赫茲以前信號的頻率比較高，但不是很高，就設(shè)計一個高通濾波器將較高的信號濾除，使產(chǎn)生一個平穩(wěn)的信號，此噪音不突出，一方面跟采樣點數(shù)有關(guān)系，另一方面錄制噪音的時候尖峰不夠，但主要

15、考察濾波器設(shè)計及其分析，不影響整體。4.3實現(xiàn)matlab產(chǎn)生的正弦信號及其代碼：f1=22000; % 輸入信號f1f2=25000; % 輸入信號f2fs=50000; % fs為抽樣頻率 t=1/fs; length_x=length(x); %輸入信號的長度input_x=zeros(1,length_x);for i=1:1:length_x input_x(i)=sin(2*pi*f1*i*t)+sin(2*pi*f2*i*t);input_x=input_x'%將input_x轉(zhuǎn)置end figure(7); plot(input_x); title('產(chǎn)生信號的

16、時域波形');grid on;figure(8); freqz(input_x); title('產(chǎn)生信號的頻譜波形');此模塊先設(shè)置了兩個信號的頻率，抽樣頻率要大于等于2倍的f1和f2，i=1:1:length_x是從1開始到length_x逐步取值賦給x，且步長為1，input_x=input_x'是對input_x進行轉(zhuǎn)置本來是100011矩陣轉(zhuǎn)置后變?yōu)?10001矩陣，grid on是畫圖時出現(xiàn)網(wǎng)格，其他函數(shù)以上均有說明。實驗結(jié)果圖如下：以上兩個圖形是產(chǎn)生信號的波形及頻譜，由于采樣點數(shù)很大，導(dǎo)致時域波形看起來就是一個藍色矩形，點數(shù)少了能看見是一個正弦信號

17、。從頻譜上，在f=21000赫茲處出現(xiàn)了一個尖峰。4.4實現(xiàn)信號的疊加的代碼：for i=1:1:length_x y1=X_voice+X3_voice+input_x;endfigure(11);freqz(y1);title('疊加后信號的頻譜');信號疊加要考慮到長度相同，且矩陣也要相同。實驗結(jié)果圖如下：綜合原始語音，噪聲語音，還有產(chǎn)生的信號的頻譜分析，在0.9出出現(xiàn)了一個尖峰，原來的尖峰依舊存在，疊加信號的頻譜跟原始語音信號非常相似，產(chǎn)生的信號疊加后比原始信號頻譜多了一個尖峰，但不明顯，相反噪音信號無尖峰，對正常語音信號影響并不大，可能是點數(shù)沒有取好，但只要有影響就行

18、。后期處理好就行。而時域波形可以看見有很多蜂刺，那些蜂刺可能是噪音的影響。4.5由原始語音設(shè)計的濾波器對疊加后的信號處理及實現(xiàn)其代碼：f_L1=2600;%其頻率值f_H1=3000;f_L2=11000; f_H2=11500;f_L3=21900;f_H3=22600;w_L1=2*pi*f_L1/fs;%轉(zhuǎn)化為數(shù)字頻率值w_H1=2*pi*f_H1/fs;w_L2=2*pi*f_L2/fs; w_H2=2*pi*f_H2/fs;w_L3=2*pi*f_L3/fs;w_H3=2*pi*f_H3/fs;N=22001;M=N-1;h=zeros(1,N);for i=0:N-1%多通多阻濾波

19、器的設(shè)計 if i=M/2 h(i+1)=(pi-w_H1+w_L1-w_H2+w_L2-w_H3+w_L3)/pi; else h(i+1)=(sin(i-M/2)*pi)-sin(i-M/2)*w_H1)+sin(i-M/2)*w_L1)-sin(i-M/2)*w_H2)+sin(i-M/2)*w_L2)-sin(i-M/2)*w_H3)+sin(i-M/2)*w_L3)/(i-M/2)*pi); end end figure(16); plot(h);title('由原始語音設(shè)計的濾波器的頻譜 ');figure(17);freqz(h); title('經(jīng)濾波后

20、的y1');y3=conv(y1,h);freqz(y3);其中for循環(huán)是設(shè)計多通多阻，兩個低通相減就是帶通，再用全通減去帶通就得到帶阻，帶通的代碼sin(i-M/2)*w_H1)-sin(i-M/2)*w_L1) /(i-M/2)*pi)，全通的代碼sin(i-M/2)*pi)/(i-M/2)*pi)，以上的頻率是分析其尖峰所選擇的頻率。其濾波器頻域波形及濾波后的時域波形、頻域波形：上圖是經(jīng)過濾波后y1的頻譜，從圖上可以看出有尖峰的地方已經(jīng)濾掉可以看出濾波器設(shè)計還算成功上圖是設(shè)計的濾波器的頻譜，含有多通多阻，從圖中可以看出有三處是阻帶，就是為了濾掉y1的三個尖峰，其他可以安全通過。

21、下圖為濾波后y1的時域波形，波形大概成對稱性4.5由噪音語音設(shè)計的濾波器對疊加后的信號處理及實現(xiàn)其代碼：N=1001;M=N-1;wc=2*pi*6700/fs;%截至的角頻率for i=0:N-1 %設(shè)計的高通濾波器 if i=M/2 h(i+1)=(pi-wc)/pi; % 當(dāng)i=M/2時的極限值 else h(i+1)=(sin(i-M/2)*pi)-sin(i-M/2)*wc)/(i-M/2)*pi); % 高通濾波器的h(n)表達式 end end figure(12);freqz(h);title('由噪聲設(shè)計的濾波器的頻譜');y=conv(y1,h);figur

22、e(13);freqz(y);title('經(jīng)過濾波后的y1')這是一個高通濾波器，高通濾波器=全通濾波器-低通濾波器，wc是濾波器的截止頻率在此處為過渡段，其濾波器頻域波形及濾波后的時域波形、頻域波形：上圖是設(shè)計的高通濾波器的頻譜圖上圖是經(jīng)過高通濾波器后y1的時域波形，高通對疊加信號的時域波形影響并不大只是尖峰變密了。上圖是經(jīng)過濾波后y1的頻譜波形，可以看見將噪音信號的較高頻率的濾掉了。4.5由matlab產(chǎn)生的信號設(shè)計的濾波器對疊加后的信號處理及實現(xiàn)其代碼：N=1000;M=N-1;wc=2*pi*21000/fs;for i=0:N-1%設(shè)計的低通濾波器 if i=M/2

23、 h(i+1)=wc/pi; % 當(dāng)i=M/2時的極限值 else h(i+1)=sin(i-M/2)*wc)/(i-M/2)*pi); % 低通濾波器的h(n)表達式end end figure(14);freqz(h); %顯示低通濾波器的頻率響應(yīng)函數(shù)title('由產(chǎn)生的信號設(shè)計出的濾波器的頻譜');y2=conv(y1,h);figure(15);freqz(y2);title('經(jīng)過濾波后的y1');低通濾波器的函數(shù)表達為sin(i-M/2)*wc)/(i-M/2)*pi);其濾波器頻域波形及濾波后的時域波形、頻域波形：上圖是根據(jù)產(chǎn)生信號頻譜設(shè)計的濾波器的頻譜此為低通濾波器。上圖為濾波后的y1時域波形，此時域圖跟疊加后時域波形變化巨大看上去更像正弦信號。上圖是經(jīng)過低通濾波器y1的頻譜?？梢钥吹綄a(chǎn)生信號頻譜的尖峰濾掉了。4.6