回聲信號(hào)檢測(cè)系統(tǒng)設(shè)計(jì)28終期報(bào)告(共17頁(yè))

1、信號(hào)與系統(tǒng)A（2）課程(kchng)項(xiàng)目實(shí)施(shsh)報(bào)告 題目(tm)( B )： 回聲信號(hào)檢測(cè)系統(tǒng)設(shè)計(jì) 組 號(hào)： B28 任課教師： 組 長(zhǎng)：_ _成 員：_成 員： _成 員： _成 員：_聯(lián)系方式： 二零 15年 5 月 22 日目錄(ml)：一、課程(kchng)項(xiàng)目實(shí)施方案二、項(xiàng)目的設(shè)計(jì)(shj)過(guò)程單次回聲的分析與消除1、回聲的產(chǎn)生及信號(hào)的時(shí)域、頻域分析2、自相關(guān)函數(shù)分析回聲延時(shí)量3、用迭代法求回聲衰減系數(shù)4、回聲的消除多次回聲的分析，與消除1、回聲的產(chǎn)生和回聲延時(shí)量的分析2、用歸一化法估算回聲衰減系數(shù)3、回聲的消除三、項(xiàng)目的結(jié)果分析四、分工合作情況參

2、考文獻(xiàn)附錄（matlab代碼）課程(kchng)項(xiàng)目實(shí)施方案題目(tm)要求：B組 回聲信號(hào)檢測(cè)系統(tǒng)(xtng)設(shè)計(jì)設(shè)計(jì)一個(gè)回聲信號(hào)檢測(cè)系統(tǒng)，能夠檢測(cè)出原信號(hào)中包含的回聲信號(hào)，并分析回聲信號(hào)的參數(shù)。要求：（1）輸入混有回聲的聲音（或信號(hào)源）；（2）設(shè)計(jì)檢測(cè)系統(tǒng)，分析回聲信號(hào)產(chǎn)生的原因、分析回聲信號(hào)的幅度、分析回聲信號(hào)的延遲量； (3) 能夠現(xiàn)場(chǎng)演示。設(shè)計(jì)思想：我們組的設(shè)計(jì)思想是利用matlab軟件平臺(tái)，用程序人為制造回聲，并且在僅知道回聲信號(hào)的條件下，分析延時(shí)量和衰減系數(shù)，從而消除回聲，并且測(cè)出聲源與障礙物的距離。項(xiàng)目的設(shè)計(jì)過(guò)程單次回聲的分析與消除1、回聲的產(chǎn)生及信號(hào)的時(shí)域、頻域分析對(duì)于已經(jīng)采

3、集好的一段語(yǔ)音信號(hào)，利用matlab的”audioread”命令可以對(duì)該信號(hào)進(jìn)行讀取和采樣，本項(xiàng)目中采用的文件為wav文件，采樣頻率為44100.x,fs=audioread(E:ceshi.wav);利用plot和fft函數(shù)可以對(duì)已采集的信號(hào)分別繪制時(shí)域波形圖和經(jīng)快速傅里葉變換得到的頻譜圖；abs和angle函數(shù)可以繪制信號(hào)頻域的幅值與相位圖。 利用zeros函數(shù)創(chuàng)建空矩陣可以對(duì)原始(yunsh)信號(hào)進(jìn)行延時(shí)，并產(chǎn)生回聲。本例子中，假設(shè)(jish)延時(shí)量為8000個(gè)采樣點(diǎn)，衰減系數(shù)為0.4。同樣對(duì)回聲(hushng)信號(hào)y進(jìn)行時(shí)域、頻域分析。2、自相關(guān)函數(shù)分析(fnx)回聲延時(shí)量自相關(guān)函數(shù)(

4、hnsh)的數(shù)學(xué)定義為：同一時(shí)間函數(shù)在瞬時(shí)t和t+a的兩個(gè)(lin )值相乘積的平 HYPERLINK /view/1052684.htm t _blank 均值作為延遲時(shí)間t的函數(shù)，它是信號(hào)與延遲后信號(hào)之間相似性的度量。簡(jiǎn)單地說(shuō)，自相關(guān)函數(shù)是用來(lái)表征一個(gè)信號(hào)在不同時(shí)刻t1,t2狀態(tài)之間的相關(guān)程度，其在某一時(shí)間點(diǎn)的值越大，相關(guān)程度越高。在matlab中，對(duì)于一個(gè)混有回聲的語(yǔ)音信號(hào)，可以利用xcorr自相關(guān)函數(shù)可以求出回聲信號(hào)與原始信號(hào)的采樣點(diǎn)之差，進(jìn)而求得延時(shí)量。 從自相關(guān)函數(shù)圖像中可以看到，在采樣點(diǎn)的中心(zhngxn)有一最大的幅值，其含義就是當(dāng)延時(shí)量為零時(shí)，信號(hào)相關(guān)程度最高；而在這一峰值

5、的兩側(cè)有兩個(gè)次大的幅值，其含義就是在該處信號(hào)相關(guān)程度次高；而在別的點(diǎn)處其相關(guān)程度相對(duì)于該2個(gè)峰值都可以忽略，因此算出這2個(gè)峰值所對(duì)應(yīng)的采樣點(diǎn)之差，就是回聲信號(hào)中，原始信號(hào)與回聲信號(hào)對(duì)應(yīng)的采樣點(diǎn)之差。由此，可以計(jì)算(j sun)出該信號(hào)中采樣點(diǎn)之差為8000個(gè)單位，這與我們所假設(shè)的值是吻合的。最后，將這個(gè)采樣點(diǎn)之差除以采樣頻率，即，便可得到(d do)實(shí)際的延時(shí)量。3、用迭代法求回聲衰減系數(shù)已知回聲的延遲量后，還需要知道回聲的衰減系數(shù)才可以確定回聲信號(hào)。對(duì)于一個(gè)一次回聲信號(hào)，我們可以采用迭代的方法來(lái)求衰減系數(shù):設(shè)xn為原始信號(hào)，yn為原始信號(hào)與回聲信號(hào)的疊加則有:y(1)=x(1)Y(1+N)=

6、X(1+N)+aX(1)Y(1+2N)=x(1+2N)+ax(1+N)y(1+(k-1)N)=x(1+(k-1)N)+ax(1+(k-2)N)Y(1+kN)=ax(1+(k-1)N)上述迭代(di di)的公式可以由下面線(xiàn)段圖來(lái)形象(xngxing)的說(shuō)明(shumng)：y序列的第一個(gè)點(diǎn)是沒(méi)有回聲的，在第（1+N）、（1+2N）、（1+（k-1）N）的點(diǎn)都是由原始信號(hào)與衰減后回聲信號(hào)疊加而成的，最后一個(gè)（1+kN）個(gè)點(diǎn)是僅含衰減后的回聲的。由圖像可以看出:設(shè)y(n)的長(zhǎng)度為L(zhǎng)，對(duì)于k,則滿(mǎn)足1+kNL（且k是整數(shù)）這個(gè)k是唯一的并且可求的。因此由最后一個(gè)式子逐級(jí)帶回上面的式子，將含x(n)的

7、項(xiàng)全迭代成含y(n)和衰減系數(shù)的項(xiàng)，最后可以得到一個(gè)關(guān)于的一元k次方程。1/a*y(1+kN)=y(1+(k-1)N)-ay(1+(k-2)N)+a2*y(1+(k-3)N)-a3*y(1+(k-4)N)+.ak*y(1)接著用matlab中的roots函數(shù)可以對(duì)進(jìn)行求解。但是對(duì)于求解出的k個(gè)根，排除虛數(shù)和不在01范圍內(nèi)的數(shù)，可能仍有一些根是滿(mǎn)足條件的。這時(shí)候我們可以將上述的迭代公式從y(2)開(kāi)始計(jì)算，即：y(2)=x(2)Y(2+N)=X(2+N)+aX(2)Y(2+2N)=x(2+2N)+ax(2+N)y(2+(k-1)N)=x(2+(k-1)N)+ax(2+(k-2)N)Y(2+kN)=

8、ax(2+(k-1)N)對(duì)于同一個(gè)信號(hào)衰減系數(shù)必定相同，那么在兩個(gè)方程中，相同的根便是所求的。我們組為了實(shí)現(xiàn)(shxin)程序的的普遍性，對(duì)程序進(jìn)行了優(yōu)化，自動(dòng)(zdng)找出兩組解中相同的解。因此(ync)用迭代法可以精準(zhǔn)的算出回聲的衰減系數(shù)。4、回聲的消除由上面求解(qi ji)得到的延時(shí)量N_delay=8000，衰減系數(shù)=0.4，可以(ky)寫(xiě)出原始信號(hào)與含回聲信號(hào)的關(guān)系式y(tǒng)(n)=x(n)+0.4x(n-8000),其系統(tǒng)函數(shù)為H(z)= QUOTE 在matlab中可以用filter函數(shù)對(duì)y信號(hào)(xnho)進(jìn)行濾波：z即為濾波完后的信號(hào)，sound(z,fs)發(fā)現(xiàn)單次回聲已經(jīng)全部消

9、除。多次回聲的分析，與消除1、回聲的產(chǎn)生和回聲延時(shí)量的分析多次回聲的產(chǎn)生原理與一次回聲類(lèi)似，利用zeros函數(shù)創(chuàng)建空矩陣可以對(duì)原始信號(hào)進(jìn)行延時(shí)。本例中延時(shí)量分別為8000個(gè)和12000個(gè)采樣點(diǎn)，衰減系數(shù)分別為0.4和0.2.同樣的，用xcorr自相關(guān)函數(shù)計(jì)算延遲量。 其中(qzhng)N_delay1是第一次回聲延遲的采樣點(diǎn)，N_delay2是第二次回聲延遲的采樣點(diǎn)。2、用歸一化法估算( sun)回聲衰減系數(shù)加入(jir)兩次回聲后，無(wú)法采用之前(zhqin)單次回聲的方法計(jì)算衰減系數(shù)，因?yàn)闊o(wú)法確定唯一的k來(lái)列出迭代的方程式。我們采用自相關(guān)函數(shù)的歸一化法估算回聲衰減系數(shù)。具體過(guò)程如下：設(shè)Ry(

10、n)為y(n)序列的自相關(guān)函數(shù)，Rx(n)為x(n)序列的自相關(guān)函數(shù)。由自相關(guān)函數(shù)的定義可以化簡(jiǎn)得到以下方程：Ry(n)=(1+ QUOTE + QUOTE )Rx(n) + QUOTE Rx(n- QUOTE ) +Rx(n+ QUOTE ) QUOTE Rx(n- QUOTE ) +Rx(n+ QUOTE )+ QUOTE Rx(n+ QUOTE - QUOTE ) +Rx(n+ QUOTE )將n=0, QUOTE ,N2分別代入R-1(0)=(1+ QUOTE + QUOTE )Rx(0)+ QUOTE Rx(- QUOTE ) +Rx( QUOTE ) QUOTE Rx(- QUOT

11、E ) +Rx( QUOTE )+ QUOTE Rx( QUOTE - QUOTE ) +Rx( QUOTE ) （2）考慮到Rx(- QUOTE ) =Rx( QUOTE )約等于0，Rx(- QUOTE ) =Rx( QUOTE )約等于0，Rx( QUOTE - QUOTE ) =Rx( QUOTE )約等于0，（2）式化為Ry (0)=(1+ QUOTE + QUOTE )Rx(0) （3）同理有Ry( QUOTE )約等于 QUOTE Rx(0) ，Ry( QUOTE )約等于 QUOTE Rx(0)因此有 QUOTE , QUOTE 還可推廣到n次回聲由此可以(ky)估算出第K次衰

12、減系數(shù)。 利用(lyng)matlab的fsolve解非線(xiàn)性方程組的功能，便可估算( sun)出2個(gè)回聲的衰減量。與理論的0.4和0.2有5%左右的誤差。3、回聲的消除同樣利用filter函數(shù)。H(z)= QUOTE 執(zhí)行sound(z,fs)后，回聲明顯不見(jiàn)了，可以達(dá)到效果。項(xiàng)目的結(jié)果分析利用回聲測(cè)距上述(shngsh)兩例單次回聲和多次回聲的消除效果(xiogu)是良好的，消除后的聲音完全聽(tīng)不出回聲。此外通過(guò)(tnggu)對(duì)回聲的分析可以應(yīng)用與測(cè)距。其原理很簡(jiǎn)單：通過(guò)自相關(guān)函數(shù)我們已求出回聲的延遲采樣點(diǎn)個(gè)數(shù)，N_delay。其真實(shí)時(shí)間其中Fs是采樣率，t的單位為秒。則聲源與障礙物的距離。以

13、報(bào)告中的例子舉例：?jiǎn)未位芈曋蠳_delay=8000，F(xiàn)s=44100,則可算出t=0.18s，那么聲源與障礙物的距離s=340*0.18/2=30.6m四、分工合作情況查找資料：盛亦池 胡斌 劉楊雨編寫(xiě)代碼：盛亦池 應(yīng)祺超代碼調(diào)試：胡斌 應(yīng)祺超 馬嘉敏終期報(bào)告：馬嘉敏 劉楊雨參考文獻(xiàn)：1、自相關(guān)函數(shù) 百度百科/link?url=OWeCZQXASKxx4INScryWy5uXYERPm00lFVNbzuRY-H-s8LgSyr63aq4-Pej7Vi9C3YjS4ah6jyhJjEuOgtykwa2、matlab課程設(shè)計(jì)-基于MATLAB的回波信號(hào)的產(chǎn)生與消除 道客巴巴 HYPERLINK

14、/p-7062907417869.html /p-7062907417869.html3、信號(hào)與系統(tǒng)實(shí)驗(yàn)指導(dǎo)/p-4147988583065.html4、基于MATLAB的數(shù)字信號(hào)處理與實(shí)踐開(kāi)發(fā)，王佳梅，西安電子科技大學(xué)出版社附錄一次回聲matlab代碼：%制造回聲x,fs=audioread(E:ceshi.wav);N=length(x);x1=x(1:N);x2=x(1:N);x1=x1,zeros(1,8000);x2=zeros(1,8000),0.4*x2;y=x1+x2;%原始信號(hào)(xnho)時(shí)域、頻域、幅值、相位圖n1=0:22050;figure(1)subplot(2,2,

15、1);plot(x);title(原始(yunsh)信號(hào)時(shí)域圖);X=fft(x);x1=x(1:4:length(x);X1=x1(:,1);df=fs/length(X);Fx=df*(0:length(X)-1);subplot(2,2,2);plot(Fx(1:44100/4),X(1:44100/4);title(原始(yunsh)信號(hào)頻譜圖);subplot(2,2,3);plot(abs(X);title(原始信號(hào)幅值);subplot(2,2,4);plot(angle(X);title(原始信號(hào)相位);%疊加信號(hào)的時(shí)域、頻域、幅值、相位圖figure(2);subplot(2

16、,2,1);plot(y);title(帶有回聲回聲時(shí)域圖);Y=fft(y);y1=y(1:4:length(y);Y1=y1(:,1);df2=fs/length(Y);Fy=df2*(0:length(Y)-1);subplot(2,2,2);plot(Fy(1:44100/4),Y(1:44100/4);title(帶有回聲信號(hào)頻譜圖);subplot(2,2,3);plot(abs(Y);title(回聲幅值);subplot(2,2,4);plot(angle(Y);title(回聲(hushng)相位);%自相關(guān)(xinggun)函數(shù) 算延時(shí)量cy=xcorr(y);figure

17、(3)plot(cy);title(多次回聲的自相關(guān)函數(shù)(hnsh)圖像);h1_y,h1=max(cy);cy2=cy;cy2(1,h1-1000:h1+1000)=0;h2_y,h2=max(cy2);N_delay=abs(h1-h2);kmax=floor(length(y)+N_delay-1)/N_delay);for k=1:kmax t1(k)=(-1)k*y(1,N_delay*(k-1)+1);endjie1=roots(t1);for i=1:length(jie1) if isreal(jie1(i) jie1(i)=jie1(i); else jie1(i)=0; e

18、ndendfor k=1:kmax t2(k)=(-1)k*y(1,N_delay*(k-1)+2);endjie2=roots(t2);for i=1:length(jie2) if isreal(jie2(i) jie2(i)=jie2(i); else jie2(i)=0; endendjie1(find(jie1(:,1)=0)= ;jie2(find(jie2(:,1)=0)= ;jie1=jie1;jie2=jie2;jie=(intersect(int16(jie1*1000),int16(jie2*1000);format shortfinaljie=double(jie)/1

19、000;%消除(xioch)回聲a=zeros(1,N);a(1)=1;a(N_delay+1)=finaljie;z=filter(1,a,y);多次回聲(hushng)matlab代碼：x,fs=audioread(E:ceshi.wav);%制造(zhzo)回聲N=length(x);x1=x(1:N);x2=x(1:N);x3=x(1:N);x1=x1,zeros(1,12000);x2=zeros(1,8000),0.4*x2,zeros(1,4000);x3=zeros(1,12000),0.2*x3;y=x1+x2+x3;%原始信號(hào)時(shí)域、頻域、幅值、相位圖n1=0:22050;f

20、igure(1)subplot(2,2,1);plot(x);title(原始信號(hào)時(shí)域圖);X=fft(x);x1=x(1:4:length(x);X1=x1(:,1);df=fs/length(X);Fx=df*(0:length(X)-1);subplot(2,2,2);plot(Fx(1:44100/4),X(1:44100/4);title(原始信號(hào)頻譜圖);subplot(2,2,3);plot(abs(X);title(原始信號(hào)幅值);subplot(2,2,4);plot(angle(X);title(原始相位);%疊加信號(hào)的時(shí)域、頻域、幅值、相位圖figure(2);subplot(2,2,1);plot(y);title(帶有回聲信號(hào)時(shí)域圖);Y=fft(y);y1=y(1:4:length(y);Y1=y1(:,1);df2=fs/length(Y);Fy=df2*(0:length(Y)-1);subplot(2,2,2);plot(Fy(1:44100/4),Y(1:44100