dsp-軟件實(shí)驗(yàn)報(bào)告(精)

1、數(shù)字信號(hào)處理MATLAB仿真實(shí)驗(yàn)報(bào)告學(xué)院：電子工程學(xué)院班級(jí)：2011211203學(xué)號(hào)：2011210876姓名：孫月鵬班內(nèi)序號(hào)：04一、實(shí)驗(yàn)一：數(shù)字信號(hào)的FFT分析、實(shí)驗(yàn)內(nèi)容及要求xfn)=0.001*cos0.45口站)+sin(0,3n)-eos(0.302jlt-4(1離散信號(hào)的頻譜分析:設(shè)信號(hào)此信號(hào)的0.3pi和0.302pi兩根譜線相距很近，譜線0.45pi的幅度很小，請(qǐng)選擇合適的序列長(zhǎng)度N和窗函數(shù)，用DFT分析其頻譜，要求得到清楚的三根譜線。(2DTMF信號(hào)頻譜分析用計(jì)算機(jī)聲卡采用一段通信系統(tǒng)中電話雙音多頻（DTMF）撥號(hào)數(shù)字09的數(shù)據(jù)，采用快速傅立葉變換（FFT）分析這10個(gè)號(hào)碼

2、DTMF撥號(hào)時(shí)的頻譜。2、實(shí)驗(yàn)結(jié)果x（n的時(shí)域圖與頻譜：得到三根清晰的譜線號(hào)碼9的頻譜號(hào)碼8的頻譜號(hào)碼7的頻譜號(hào)碼6的頻譜、實(shí)現(xiàn)代碼及分析（1第一小題：k=1000;%DFT點(diǎn)數(shù)n=1:1:k;%對(duì)時(shí)域信號(hào)進(jìn)行采樣x=0.001*cos(0.45*n*pi+sin(0.3*n*pi-cos(0.302*n*pi-pi/4;subplot(2,1,1;stem(n,x,'.'%用.畫出時(shí)域圖title('時(shí)域序列'xlabel('n'ylabel('x(n'xk=fft(x,k;%進(jìn)行K點(diǎn)DFT變換w=2*pi/k*0:1:k-1;

3、%數(shù)字角頻率subplot(2,1,2;stem(w/pi,abs(xk;%畫出頻譜圖axis(0.2,0.5,0,2;%設(shè)置窗函數(shù)的寬度與限幅title('1000點(diǎn)dft'xlabel('數(shù)字頻率'ylabel('lxk(kl'%此題關(guān)鍵在于DFT點(diǎn)數(shù)N的確定。經(jīng)過(guò)計(jì)算和實(shí)驗(yàn)，當(dāng)N=1000時(shí)能滿足題目要求，看到3條清晰地譜線(2)第二小題clear;closeall;f=9411336;6971209;6971336;6971477;7701209;7701336;7701477;8521209;8521336;8521477%0-9的頻率

4、n=1:400;fs=4000;%取樣頻率為4000hzfprintf('請(qǐng)輸入數(shù)字(0to9：n'k=input(''f1=f(k+1,1;%因?yàn)閺?開始計(jì)算，+1得輸f2=f(k+1,2;%入數(shù)字的兩個(gè)頻率N=400;x1=sin(2*pi*f1*n/fs+sin(2*pi*f2*n/fs;%DTMF的輸入信號(hào)時(shí)域xn=x1,zeros(1,400;%補(bǔ)零subplot(2,1,1;plot(xn%畫出時(shí)域圖xlabel('n'ylabel('xn'subplot(2,1,2;fn=fs*n/N;%取樣點(diǎn)的頻率plot(fn,

5、abs(fft(xn(l:400;%400點(diǎn)fft變換，畫出頻譜圖axis(0,4000,0,300xlabel('f'ylabel('FFT'二、實(shí)驗(yàn)二：DTMF信號(hào)的編碼1、實(shí)驗(yàn)內(nèi)容及要求1）把您的聯(lián)系電話號(hào)碼通過(guò)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)碼傳送過(guò)程中，信號(hào)存在時(shí)間至少45ms,且不多于55ms,100ms的其余時(shí)間是靜音。2）對(duì)所生成的DTMF文件進(jìn)行解碼。口由于只需要知道8個(gè)特定點(diǎn)的頻譜值，因此采用一種稱為Goertzel

6、算法的IIR濾波器可以有效地提高計(jì)算效率。其傳輸函數(shù)為Ng=2、實(shí)驗(yàn)結(jié)果kTTT1LU,.T"輸入號(hào)碼界面生成的時(shí)域圖（占空比為50%）電話號(hào)碼的FFT圖：每個(gè)數(shù)字都有兩條主譜線顯示檢測(cè)到的號(hào)碼3）實(shí)驗(yàn)代碼及分析%第一部分，產(chǎn)生編碼clc;tm=49,50,51,65;52,53,54,66;55,56,57,67;42,48,35,68;%DTMF表中的ASCII碼f1=697,770,852,941;%行頻率向量f2=1209,1336,1477,1633;%列頻率向量d=input('pleaseenternumber:','s'%輸入電話號(hào)碼s

7、um=length(d;%電話號(hào)碼長(zhǎng)度total_x=;%電話號(hào)碼信號(hào)sum_x=;sum_x=sum_x,zeros(1,800;fora=1:sum%循環(huán)sum次symbol=abs(d(a;%求輸入的ASCII碼forp=1:4;forq=1:4;iftm(p,q=abs(d(a;break,end%檢測(cè)碼相符的列號(hào)qendiftm(p,q=abs(d(a;break,end%檢測(cè)碼相符的行號(hào)pendn=1:400;x=sin(2*pi*n*fl(p/8000+sin(2*pi*n*f2(q/8000;%構(gòu)成雙頻信號(hào)x=x,zeros(1,400;%加長(zhǎng)序列，增加靜音sum_x=sum_

8、x+x;total_x=total_x,x;%將所有編碼連接起來(lái)endsound(total_x;%播放聲音t=total_x/2;wavwrite(t,'我的手機(jī)號(hào)碼'生成聲音文件plot(total_x;title('DTMF信號(hào)時(shí)域波形'%代碼主要分成三部分，即根據(jù)輸入的數(shù)字確定雙頻率、產(chǎn)生正弦信號(hào)和生成文件以及繪圖。在產(chǎn)生正弦信號(hào)的過(guò)程中有增加靜音和連接運(yùn)算。整體上用循環(huán)對(duì)電話號(hào)碼的每一位進(jìn)行相同處理。%第二部分，檢測(cè)端k=18,20,22,25,32,35,38;N=210;tm=49,50,51;52,53,54;55,56,57;0,48,0;fo

9、ri=1:sumj=800*(i-1;X=goertzel(total_x(j+1:j+N,k+1;%算法公式value=abs(X;figure(2subplot(2,6,i;stem(k,value,'.','r'%畫出FFT頻譜title('FFTx(n'xlabel('k'ylabel('|X(k|'limit=20;fori1=5:7%判決，確定頻率ifvalue(i1>limitbreak;endendforj1=1:4ifvalue(j1>limitbreak;endendbuffer(i

10、=tm(j1,i1-4;%根據(jù)頻率確定號(hào)碼enddisp('接收端檢測(cè)到的號(hào)碼'disp(setstr(buffer%顯示檢測(cè)到的號(hào)碼%檢測(cè)端編碼首先要根據(jù)goertzel公式算出FFT,再利用判決的方法找出頻率的分布，從而確定輸入的電話號(hào)碼。三、實(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ì)量?？谔崾荆嚎跒V波器指標(biāo)參考：通帶邊緣頻率為4kHz,阻帶邊緣頻率為4.5kHz，阻帶衰減大于50dB；口Matlab函數(shù)y=

11、awgn(x,snr,'measured'，首先測(cè)量輸入信號(hào)x的功率，然后對(duì)其疊加高斯白噪聲；2、實(shí)驗(yàn)結(jié)果原始聲音時(shí)域圖與頻譜圖疊加白噪聲的時(shí)域圖與頻譜圖但明顯聲音信號(hào)太弱，被白噪聲蓋過(guò)去了）濾波器的幅頻特性濾波后的時(shí)域圖與頻譜特性（可以看到45Khz以后的頻譜都被濾掉了）3、代碼及分析%產(chǎn)生原始聲音及頻譜clearallcloseally1,fs,bits=wavread('record-1.wav'%按照采樣頻率取得y1y2=y1(:,1;fs=32000;k=1:4096;yk1=fft(y2,4096;figure(1subplot(2,l,l;plot

12、(y2;title('原始聲音時(shí)域'；subplot(2,1,2;plot(32/4096*k,abs(yk1;axis(0,17,0,0.05;xlabel('f/kHZ'title('原始聲音頻域'；wavplay(yl,42000%加噪及其頻譜y3=awgn(y2,20,'measured','db'%awgn?為加性高斯白噪聲yk2=fft(y3,4096;figure(2subplot(2,1,1;plot(y3;title('加噪后聲音時(shí)域'；subplot(2,l,2;plot(32/

13、4096*k,abs(yk2;xlabel('f/KHZ'title('加噪后聲音頻域'wavplay(y3,42000;%濾波器設(shè)計(jì)fp=4000;fr=4500;wp=2*pi*fp/fs;wr=2*pi*fr/fs;tr_width=wr-wp;N=ceil(6.6*pi/tr_width+1;n=0:1:N;wc=(wr+wp/2;alpha=(N-1/2;n=0:1:N-1;m=n-alpha+eps;hd=sin(wc*m./(pi*m;b=fir1(N,wc/(4*pi/3.65;w_ham=(hamming(N'h=hd.*w_ham;H

14、,w=freqz(h,1,1000,'whole'H=(H(1:1:501'w=(w(1:1:501'mag=abs(H;db=20*log10(mag+eps/max(mag;pha=angle(H;delta_w=2*pi/1000;x=filter(b,1,y3;figure(3subplot(1,1,1;plot(w/pi,db;title('濾波器幅頻響應(yīng)'；x=conv(h,y3;%濾除噪聲后頻譜figure(4subplot(2,1,1;plot(x;title('濾除噪聲后時(shí)域圖'；Xk=fft(x,4096;su

15、bplot(2,1,2；plot(32/4096*k,abs(Xk；axis(0,16,0,5title('濾除噪聲后頻域圖xlabel('f/KHZ'；wavplay(x,42000；%根據(jù)要求，所求的濾波器是一個(gè)低通濾波器%因?yàn)樗p大于50DB，窗函數(shù)選擇了漢明窗4、錯(cuò)誤分析與總結(jié)問(wèn)題主要出在第三道題目上：1) 實(shí)驗(yàn)要求的聲音文件格式為wav，但錄音軟件的的格式wma，所以需要格式轉(zhuǎn)換。首先實(shí)驗(yàn)的是直接更改格式名，盡管在播放器中可以播放，但仍然不可以識(shí)別；其次實(shí)驗(yàn)的是用matlab直接錄音，但頻譜集中在錄音的采樣頻率附近，仍不能做正常的分析；最后用的是格式工廠的轉(zhuǎn)換，盡管感覺頻譜仍有改變，但并不影響實(shí)驗(yàn)結(jié)果。2) 用matlab對(duì)聲音進(jìn)行波形提取后，再用相同的fs播放，會(huì)發(fā)現(xiàn)聲音的頻率變低，進(jìn)過(guò)實(shí)驗(yàn)，必須用比f(wàn)s高一些的頻率播放，才能聽到正確的聲音。3) 對(duì)比其他同學(xué)的頻譜，發(fā)現(xiàn)我的原始聲音的頻譜的幅度就很低，加入白噪聲之后，就很難看到聲音的頻譜了。但聲音的播放就沒有什么問(wèn)題。直到最后也沒有想明白為什么。5、心得與體會(huì)在這次實(shí)驗(yàn)中，我們分別用MATLAB實(shí)現(xiàn)了簡(jiǎn)單的FFT變換、IIR濾波器設(shè)計(jì)與FIR濾波器