北京郵電大學(xué)DSP軟件實(shí)驗(yàn)報(bào)告

北京郵電大學(xué)DSP軟件實(shí)驗(yàn)報(bào)告學(xué)院：xxx學(xué)院姓名：xxx班級(jí)：xxxxxxxxxx學(xué)號(hào)：xxxxxxxxxx目錄實(shí)驗(yàn)一：數(shù)字信號(hào)的FFT分析 21、實(shí)驗(yàn)內(nèi)容及要求 22、實(shí)驗(yàn)?zāi)康?23、設(shè)計(jì)思路 24、實(shí)驗(yàn)代碼與實(shí)驗(yàn)結(jié)果 3實(shí)驗(yàn)二：DTMF信號(hào)的編碼 61、實(shí)驗(yàn)內(nèi)容及要求 62、實(shí)驗(yàn)?zāi)康?73、設(shè)計(jì)思路 74、實(shí)驗(yàn)代碼與實(shí)驗(yàn)結(jié)果 7實(shí)驗(yàn)三：FIR數(shù)字濾波器的設(shè)計(jì)和實(shí)現(xiàn) 121、實(shí)驗(yàn)內(nèi)容及要求： 122、實(shí)驗(yàn)?zāi)康?133、設(shè)計(jì)思路 134、實(shí)驗(yàn)代碼與運(yùn)行結(jié)果 13實(shí)驗(yàn)總結(jié) 17實(shí)驗(yàn)一：數(shù)字信號(hào)的FFT分析1、實(shí)驗(yàn)內(nèi)容及要求(1)離散信號(hào)的頻譜分析：設(shè)信號(hào)此信號(hào)的0.3pi和0.302pi兩根譜線相距很近，譜線0.45pi的幅度很小，請(qǐng)選擇合適的序列長(zhǎng)度N和窗函數(shù)，用DFT分析其頻譜，要求得到清楚的三根譜線。(2)DTMF信號(hào)頻譜分析用計(jì)算機(jī)聲卡采用一段通信系統(tǒng)中電話雙音多頻（DTMF）撥號(hào)數(shù)字0～9的數(shù)據(jù)，采用快速傅立葉變換（FFT）分析這10個(gè)號(hào)碼DTMF撥號(hào)時(shí)的頻譜。2、實(shí)驗(yàn)?zāi)康耐ㄟ^本次實(shí)驗(yàn)，應(yīng)該掌握：(a)用傅立葉變換進(jìn)行信號(hào)分析時(shí)基本參數(shù)的選擇。(b)經(jīng)過離散時(shí)間傅立葉變換（DTFT）和有限長(zhǎng)度離散傅立葉變換（DFT）后信號(hào)頻譜上的區(qū)別，前者DTFT時(shí)間域是離散信號(hào)，頻率域還是連續(xù)的，而DFT在兩個(gè)域中都是離散的。(c)離散傅立葉變換的基本原理、特性，以及經(jīng)典的快速算法（基2時(shí)間抽選法），體會(huì)快速算法的效率。(d)獲得一個(gè)高密度頻譜和高分辨率頻譜的概念和方法，建立頻率分辨率和時(shí)間分辨率的概念，為將來進(jìn)一步進(jìn)行時(shí)頻分析（例如小波）的學(xué)習(xí)和研究打下基礎(chǔ)。(e)建立DFT從整體上可看成是由窄帶相鄰濾波器組成的濾波器組的概念，此概念的一個(gè)典型應(yīng)用是數(shù)字音頻壓縮中的分析濾波器，例如DVDAC3和MPEGAudio。3、設(shè)計(jì)思路（1）由信號(hào)xn=0.001*cos0.45nπ+sin0.3nπ-cos?(0.302nπ-π4)函數(shù)限定坐標(biāo)軸范圍。（2）雙音多頻信號(hào)DTMF每個(gè)數(shù)字由兩個(gè)不同頻率的正弦波組成，低頻有：697Hz,770Hz,852Hz,941Hz，高頻有：1209Hz,1336Hz,1477Hz,1633Hz，0-9這十個(gè)數(shù)字每個(gè)數(shù)字對(duì)應(yīng)一個(gè)低頻信號(hào)和一個(gè)高頻信號(hào)疊加。分別用兩個(gè)數(shù)組裝載高頻和低頻，再產(chǎn)生由兩個(gè)正弦波疊加成的DTMF信號(hào)，最后利用plot和fft函數(shù)畫出對(duì)應(yīng)的頻譜圖。4、實(shí)驗(yàn)代碼與實(shí)驗(yàn)結(jié)果離散信號(hào)的頻譜分析【實(shí)驗(yàn)代碼】clearcloseallN=20000;n=1:1:N;x=0.001*cos(0.45*pi*n)+sin(0.3*pi*n)-cos(0.302*pi*n-pi/4);y=fft(x,N);magy=abs(y(1:1:N/2+1));k=0:1:N/2;w=2*pi/N*k;stem(w/pi,magy);axis([0.25,0.5,0,60])【實(shí)驗(yàn)結(jié)果】頻譜圖如下所示：圖SEQ圖\*ARABIC1DTMF信號(hào)頻譜分析【實(shí)驗(yàn)代碼】clearcloseallcolumn=[1209,1336,1477,1633];line=[697,770,852,941];fs=8000;N=1024;ts=1/fs;n=0:N-1;f=0:fs/N:fs/N*(N-1);key=zeros(10,N);key(1,:)=cos(2*pi*column(1)*ts*n)+cos(2*pi*line(1)*ts*n);key(2,:)=cos(2*pi*column(2)*ts*n)+cos(2*pi*line(1)*ts*n);key(3,:)=cos(2*pi*column(3)*ts*n)+cos(2*pi*line(1)*ts*n);key(4,:)=cos(2*pi*column(1)*ts*n)+cos(2*pi*line(2)*ts*n);key(5,:)=cos(2*pi*column(2)*ts*n)+cos(2*pi*line(2)*ts*n);key(6,:)=cos(2*pi*column(3)*ts*n)+cos(2*pi*line(2)*ts*n);key(7,:)=cos(2*pi*column(1)*ts*n)+cos(2*pi*line(3)*ts*n);key(8,:)=cos(2*pi*column(2)*ts*n)+cos(2*pi*line(3)*ts*n);key(9,:)=cos(2*pi*column(3)*ts*n)+cos(2*pi*line(3)*ts*n);key(10,:)=cos(2*pi*column(2)*ts*n)+cos(2*pi*line(4)*ts*n);figure;fori=1:9subplot(3,4,i)plot(f,abs(fft(key(i,:))));xlabel('頻率（Hz）');ylabel('幅值');title(i);grid;endsubplot(3,4,10)plot(f,abs(fft(key(10,:))));xlabel('頻率（Hz）');ylabel('幅值');title(0);grid;【實(shí)驗(yàn)結(jié)果】頻譜圖如下所示：圖SEQ圖\*ARABIC2實(shí)驗(yàn)二：DTMF信號(hào)的編碼1、實(shí)驗(yàn)內(nèi)容及要求（1）把您的聯(lián)系電話號(hào)碼通過DTMF編碼生成為一個(gè).wav文件。技術(shù)指標(biāo)：根據(jù)ITUQ.23建議，DTMF信號(hào)的技術(shù)指標(biāo)是：傳送/接收率為每秒10個(gè)號(hào)碼，或每個(gè)號(hào)碼100ms。每個(gè)號(hào)碼傳送過程中，信號(hào)存在時(shí)間至少45ms，且不多于55ms，100ms的其余時(shí)間是靜音。在每個(gè)頻率點(diǎn)上允許有不超過±1.5%的頻率誤差。任何超過給定頻率±3.5%的信號(hào)，均被認(rèn)為是無效的，拒絕接收。（其中關(guān)鍵是不同頻率的正弦波的產(chǎn)生。可以使用查表方式模擬產(chǎn)生兩個(gè)不同頻率的正弦波。正弦表的制定要保證合成信號(hào)的頻率誤差在±1.5%以內(nèi)，同時(shí)使取樣點(diǎn)數(shù)盡量少）（2）對(duì)所生成的DTMF文件進(jìn)行解碼。DTMF信號(hào)解碼可以采用FFT計(jì)算N點(diǎn)頻率處的頻譜值，然后估計(jì)出所撥號(hào)碼。但FFT計(jì)算了許多不需要的值，計(jì)算量太大，而且為保證頻率分辨率，F(xiàn)FT的點(diǎn)數(shù)較大，不利于實(shí)時(shí)實(shí)現(xiàn)。因此，F(xiàn)FT不適合于DTMF信號(hào)解碼的應(yīng)用。由于只需要知道8個(gè)特定點(diǎn)的頻譜值，因此采用一種稱為Goertzel算法的IIR濾波器可以有效地提高計(jì)算效率。其傳遞函數(shù)為：2、實(shí)驗(yàn)?zāi)康?a)復(fù)習(xí)和鞏固IIR數(shù)字濾波器的基本概念；(b)掌握IIR數(shù)字濾波器的設(shè)計(jì)方法；(c)掌握IIR數(shù)字濾波器的實(shí)現(xiàn)結(jié)構(gòu)；(d)能夠由濾波器的實(shí)現(xiàn)結(jié)構(gòu)分析濾波器的性能（字長(zhǎng)效應(yīng)）；(e)了解通信系統(tǒng)電話DTMF撥號(hào)的基本原理和IIR濾波器實(shí)現(xiàn)方法。3、設(shè)計(jì)思路編碼：DTMF撥號(hào)鍵盤由一個(gè)4*4行列構(gòu)成，每列代表一個(gè)高頻信號(hào)，每行代表一個(gè)低頻信號(hào)，每當(dāng)按下一個(gè)鍵時(shí)，產(chǎn)生高、低頻率的兩個(gè)正弦信號(hào)，代表一個(gè)特定的數(shù)字或符號(hào)，根據(jù)ITUQ.23頒布的國(guó)際標(biāo)準(zhǔn)，DTMF傳送或接受每個(gè)號(hào)碼的時(shí)間為100ms，其中每個(gè)號(hào)碼傳送的過程中，信號(hào)存在時(shí)間至少45ms，其余時(shí)間靜音。用一個(gè)字符串變量來接受輸入的電話號(hào)碼，并將各個(gè)數(shù)字和符號(hào)的ASCII碼用一個(gè)4*4矩陣表示，每接收到一個(gè)數(shù)字就對(duì)應(yīng)兩個(gè)頻率，并產(chǎn)生由兩個(gè)正弦波疊加的信號(hào)，完成DTMF編碼，利用MATLAB提供的fft函數(shù)畫出其頻譜，用sound函數(shù)發(fā)出聲音。解碼：采用Goertzel算法來檢測(cè)DTMF信號(hào)，它是用IIR濾波器實(shí)現(xiàn)DFT算法的一種特殊方法，在實(shí)際DTMF解碼中，只需要知道輸入信號(hào)即DTMF信號(hào)的離散傅里葉變換X(k)的幅度信息，忽略相位信息，因?yàn)橹灰艿玫?個(gè)特定頻率點(diǎn)的幅度值，看哪兩個(gè)頻率對(duì)應(yīng)的幅度最大，就能知道對(duì)應(yīng)的是哪個(gè)數(shù)字，達(dá)到解碼的目的，可以利用MATLAB提供的goertzel函數(shù)來對(duì)信號(hào)解碼。4、實(shí)驗(yàn)代碼與實(shí)驗(yàn)結(jié)果把聯(lián)系電話號(hào)碼通過DTMF編碼生成為一個(gè).wav文件【實(shí)驗(yàn)代碼】d=input('請(qǐng)輸入電話號(hào)碼：','s');sum=length(d);total_x=[];sum_x=[];sum_x=[sum_x,zeros(1,800)];fora=1:sumsymbol=abs(d(a));tm=[49,50,51,65;52,53,54,66;55,56,57,67;42,48,35,68];forp=1:4;forq=1:4;iftm(p,q)==abs(d(a));break,endendiftm(p,q)==abs(d(a));break,endendf1=[697,770,852,941];f2=[1209,1336,1477,1633];n=1:400;x=sin(2*pi*n*f1(p)/8000)+sin(2*pi*n*f2(q)/8000);x=[x,zeros(1,400)];sum_x=sum_x+x;total_x=[total_x,x];endaudiowrite('soundwave.wav',total_x,8000);sound(total_x);t=(1:8800)/8000;subplot(2,1,1);plot(t,total_x);axis([0,1.2,-2,2]);xlabel('時(shí)間/s')title('DTMF信號(hào)時(shí)域波形')xk=fft(x);mxk=abs(xk);subplot(2,1,2);k=(1:800)*sum*8000/800;plot(k,mxk);xlabel('頻率');title('DTMF信號(hào)頻譜');disp('雙頻信號(hào)已生成并發(fā)出')【實(shí)驗(yàn)結(jié)果】如下圖圖3所示，輸入本人電話號(hào)碼按下回車鍵之后能聽到撥碼聲音，并有圖4所示圖像生成。圖SEQ圖\*ARABIC3圖SEQ圖\*ARABIC4對(duì)所生成的DTMF文件進(jìn)行解碼【實(shí)驗(yàn)代碼】k=[1820222431343842];N=205;disp(['接收端檢測(cè)到的號(hào)碼為'])fora=1:summ=800*(a-1);X=goertzel(total_x(m+1:m+N),k+1);val=abs(X);xk2=val.^2;limit=80;fors=5:8ifval(s)>limit,break,endendforr=1:4ifval(r)>limit,break,endenddisp([setstr(tm(r,s-4))])xk2end【實(shí)驗(yàn)結(jié)果】在命令行窗口可以看到圖5-圖8的解碼信息，檢測(cè)到的號(hào)碼依次所得結(jié)果正確。圖SEQ圖\*ARABIC5圖SEQ圖\*ARABIC6圖SEQ圖\*ARABIC7圖SEQ圖\*ARABIC8實(shí)驗(yàn)三：FIR數(shù)字濾波器的設(shè)計(jì)和實(shí)現(xiàn)1、實(shí)驗(yàn)內(nèi)容及要求：錄制自己的一段聲音，長(zhǎng)度為45秒，取樣頻率32kHz，然后疊加一個(gè)高斯白噪聲，使得信噪比為20dB。請(qǐng)采用窗口法設(shè)計(jì)一個(gè)FIR帶通濾波器，濾除噪聲提高質(zhì)量。提示：濾波器指標(biāo)參考：通帶邊緣頻率為4kHz，阻帶邊緣頻率為4.5kHz，阻帶衰減大于50dB；Matlab函數(shù)y=awgn(x,snr,'measured')，首先測(cè)量輸入信號(hào)x的功率，然后對(duì)其疊加高斯白噪聲；2、實(shí)驗(yàn)?zāi)康耐ㄟ^本次實(shí)驗(yàn)，掌握以下知識(shí)：FIR數(shù)字濾波器窗口設(shè)計(jì)法的原理和設(shè)計(jì)步驟；Gibbs效應(yīng)發(fā)生的原因和影響；不同類型的窗函數(shù)對(duì)濾波效果的影響，以及窗函數(shù)和長(zhǎng)度N的選擇。3、設(shè)計(jì)思路首先通過audioread命令讀取指定地址的音頻文件，再利用MATLAB提供的awgn()函數(shù)為錄制的聲音信號(hào)加上高斯白噪聲，得到合成信號(hào)，通過audiowrite命令在指定地址生成新的音頻文件。要得到濾波器的沖激響應(yīng)或頻率響應(yīng)，首先根據(jù)濾波器的性能指標(biāo)得到窗函數(shù)的寬度N，再用布萊克曼窗得到濾波器的頻率響應(yīng)，經(jīng)過傅里葉反變換可得到濾波器的沖激相應(yīng)，把合成信號(hào)通過濾波器還原出原來的聲音信號(hào)。在figure(1)窗口中顯示濾波器的幅頻、相頻特性，在figure(2)窗口中顯示濾波前后的時(shí)域波形，在figure(3)窗口中顯示濾波前后的頻譜。4、實(shí)驗(yàn)代碼與運(yùn)行結(jié)果【實(shí)驗(yàn)代碼】fs=32000;[x,fs]=audioread('mysound.wav');snr=20;x2=awgn(x,snr,'measured','db');audiowrite('sound_noise.wav',x2,fs);t=0:1/fs:(size(x2)-1)/fs;wp=8000*pi/32000;ws=9000*pi/32000;wdelta=ws-wp;N=ceil(11*pi/wdelta);wn=(ws+wp)/2;b=fir1(N,wn/pi,blackman(N+1));figure(1)freqz(b,1,512)f2=filter(b,1,x2);title('濾波器幅頻、相頻特性');figure(2)subplot(2,1,1)plot(t,x2)title('濾波前時(shí)域波形');subplot(2,1,2)plot(t,f2)title('濾波后時(shí)域波形');F0=fft(f2,1024);f=fs*(0:511)/1024;figure(3)y2=fft(x2,1024);s