聲學(xué)仿真試驗(yàn)

聲學(xué)仿真試驗(yàn)第一頁，共五十六頁，編輯于2023年，星期二10.1回聲試驗(yàn) 10.1.1回聲的基本原理 回聲是一種物理現(xiàn)象。當(dāng)直達(dá)聲結(jié)束后，聲音經(jīng)其它物體反射、折射，又返回來收到的聲波，稱為回聲。回聲有時(shí)泛指一切反射回來的聲音。當(dāng)傳到人耳的直達(dá)聲和回聲之間的時(shí)間差在1/20s以上時(shí)，可以很清楚地把它們區(qū)分開。在山谷或大廳中，常有回聲現(xiàn)象發(fā)生。尤其是在影劇院里，回聲往往會(huì)妨礙聽音，所以建筑師必須考慮消除回聲的影響。相反，直達(dá)聲和回聲之間的時(shí)間差在1/20s以下時(shí)，回聲不但不妨礙聽音，反而使聲音更加響亮。第二頁，共五十六頁，編輯于2023年，星期二 應(yīng)用MATLAB軟件可以方便地進(jìn)行仿真試驗(yàn)?；芈暚F(xiàn)象取決于幾個(gè)重要的物理參數(shù)： (1)聲源距離聽者有多遠(yuǎn)。 (2)聲音產(chǎn)生往復(fù)反射的兩個(gè)反射體距離有多遠(yuǎn)。 (3)聲音的反射體在反射聲音的過程中對(duì)聲音的衰減有多大。(假設(shè)聽者在聲源與反射體之間的某個(gè)位置上。)第三頁，共五十六頁，編輯于2023年，星期二 10.1.2回聲的仿真試驗(yàn) 圖10-1所示是回聲的仿真試驗(yàn)系統(tǒng)。圖中1800的延遲代表了聲音在兩個(gè)反射體之間走一個(gè)來回的延遲（1800/Fs秒），它的值愈大，反映了兩個(gè)反射體的距離愈大。圖中0.8的增益代表了聲音在兩個(gè)反射體之間走一個(gè)來回的損耗，它愈小，回聲衰減愈快。圖中0.9的增益代表了傳輸損耗。圖10-2所示是仿真試驗(yàn)時(shí)示波器上顯示的波形圖，上圖是原聲，下圖是往復(fù)反射的回聲波形圖。如果計(jì)算機(jī)裝有聲卡和音箱，可以聽到揚(yáng)聲器發(fā)出的回聲聲音。放大器增益和整數(shù)延遲，可以看圖標(biāo)注。第四頁，共五十六頁，編輯于2023年，星期二圖10-1回聲仿真試驗(yàn)系統(tǒng)框圖第五頁，共五十六頁，編輯于2023年，星期二圖10-2回聲仿真試驗(yàn)時(shí)域圖第六頁，共五十六頁，編輯于2023年，星期二 回聲的仿真試驗(yàn)的模型參數(shù)設(shè)置分析：回聲的仿真試驗(yàn)的模型參考了MATLAB中Toolbox＼dspblks＼dspdemos＼dspafxr的DEMOS例子，在仿真系統(tǒng)中設(shè)置了Reshape（整形器），主要參數(shù)見表10-1，ToWaveDevice（揚(yáng)聲器）的主要參數(shù)見表10-2，F(xiàn)romWorkspace（聲源）的主要參數(shù)見表10-3，它以幀的形式輸出。整形器的設(shè)置目的是將信號(hào)轉(zhuǎn)變?yōu)槭静ㄆ髂軌蚪邮艿?維數(shù)據(jù)流。第七頁，共五十六頁，編輯于2023年，星期二表10-1Reshape（整形器）的主要參數(shù)第八頁，共五十六頁，編輯于2023年，星期二表10-2ToWaveDevice（揚(yáng)聲器）的主要參數(shù)第九頁，共五十六頁，編輯于2023年，星期二表10-3FromWorkspace（來自工作空間）的主要參數(shù)第十頁，共五十六頁，編輯于2023年，星期二 在10.1、10.4、10.5節(jié)中都需要一段聲音文件來進(jìn)行仿真試驗(yàn)?？梢杂袃煞N方法獲得聲音文件： (1)調(diào)出MATLAB中的聲音文件“l(fā)oadmtlb.mat”。 (2)應(yīng)用6.2.2節(jié)介紹的聲音采集的方法，建立一個(gè)聲音文件。譬如文件名是Soud11.mat。 通過下列程序調(diào)入仿真的聲音數(shù)據(jù)并啟動(dòng)仿真。Echo11是圖10-1所示的仿真系統(tǒng)的名字（文件名）。 程序10-1 loadsoud11%讀入聲音數(shù)據(jù)文件soud11 sim（′echo11′）第十一頁，共五十六頁，編輯于2023年，星期二10.2多普勒效應(yīng) 10.2.1多普勒效應(yīng)的基本原理 當(dāng)波源或觀察者相對(duì)于媒質(zhì)運(yùn)動(dòng)時(shí)，或者說波源和觀察者有相對(duì)運(yùn)動(dòng)時(shí)，觀察者接受到的振動(dòng)頻率與波源振動(dòng)頻率不同的現(xiàn)象，稱為多普勒效應(yīng)。 當(dāng)聽者與聲源的相對(duì)速度大到與聲速可以比擬時(shí)，就可以明顯地感覺到聲音頻率的變化。第十二頁，共五十六頁，編輯于2023年，星期二 譬如坐火車時(shí)，聽到相向開來的火車汽笛聲，當(dāng)火車由遠(yuǎn)而近開來，汽笛聲的音調(diào)變高；由近而遠(yuǎn)離去，汽笛聲的音調(diào)變低。這是日常生活中的一個(gè)多普勒效應(yīng)的例子。在天文、通信等領(lǐng)域還有眾多的例子。 以下的公式描述了該現(xiàn)象的各個(gè)物理量之間的定量關(guān)系： 其中,f0是聲源發(fā)出的聲音的頻率，v是聽者與聲源相對(duì)運(yùn)動(dòng)的速度，θ為速度矢量與聲源和聽者的連線夾角，vs為聲音在空氣中傳播的速度，f是聽者聽到的聲音頻率。(10-1）第十三頁，共五十六頁，編輯于2023年，星期二

10.2.2多普勒效應(yīng)的仿真試驗(yàn) Simulink仿真模型是參考MATLAB中的Toolbox＼dspblks＼dspdemos＼dspstfft建立的。 下面是一段用MATLAB的M文件編寫的程序。它產(chǎn)生一個(gè)文件名為u1的表現(xiàn)多普勒效應(yīng)的聲音文件，描述火車向一個(gè)距鐵路30m(垂直距離）、距火車510m（水平距離）的聽者開來時(shí)他聽到的聲音。將u1送入圖10-3所示的仿真系統(tǒng)，再一次聽到該聲音，并且看到用兩種方法（頻譜儀和矩陣圖儀）表現(xiàn)的聲音的時(shí)頻圖。表10-4~表10-11分別是圖10-3所示仿真系統(tǒng)中各個(gè)模塊的主要參數(shù)。圖10-4所示是M文件程序運(yùn)行結(jié)束時(shí)繪出的聲音信號(hào)的原聲波形（上圖）和多普勒效應(yīng)的波形（下圖）。圖10-5所示是某個(gè)時(shí)刻圖10-3中的頻譜儀（左圖）和矩陣儀（右圖）上顯示的時(shí)頻圖。第十四頁，共五十六頁，編輯于2023年，星期二圖10-3顯示聲音信號(hào)(多普勒響應(yīng))的仿真系統(tǒng)模型第十五頁，共五十六頁，編輯于2023年，星期二表10-4FromWorkspace（來自工作空間）的主要參數(shù)第十六頁，共五十六頁，編輯于2023年，星期二表10-5Rebuffer（緩存器）的主要參數(shù)第十七頁，共五十六頁，編輯于2023年，星期二表10-6ShortTimeFFT（短時(shí)傅立葉變換）的主要參數(shù)第十八頁，共五十六頁，編輯于2023年，星期二表10-7FrequencyFrameScope（頻率幀示波器）的主要參數(shù)第十九頁，共五十六頁，編輯于2023年，星期二表10-8MatrixViewer（矩陣顯示圖儀）的主要參數(shù)第二十頁，共五十六頁，編輯于2023年，星期二表10-9Select（選擇器）的主要參數(shù)第二十一頁，共五十六頁，編輯于2023年，星期二表10-10dB（分貝）的主要參數(shù)第二十二頁，共五十六頁，編輯于2023年，星期二表10-11ToWaveDevice（揚(yáng)聲器）的主要參數(shù)第二十三頁，共五十六頁，編輯于2023年，星期二圖10-4聲源發(fā)出的信號(hào)（上圖）和聽者接收到的信號(hào)（下圖，多普勒效應(yīng)）第二十四頁，共五十六頁，編輯于2023年，星期二圖10-5頻譜儀顯示（左圖）和矩陣表示（右圖）第二十五頁，共五十六頁，編輯于2023年，星期二 相應(yīng)程序如下： 程序10-2 %多普勒效應(yīng)程序 x0=550; v=150;%車速 y=30; c=330;%聲音在空氣中的傳播速度 w=2200;%聲音的頻率 t=0:.0005:30;第二十六頁，共五十六頁，編輯于2023年，星期二

r=sqrt（（x0-v*t）.^2+（y）.^2）;%計(jì)算聲源與聽者距離 t1=t-r/c;%經(jīng)距離遲延后聽者的等效時(shí)間 u=1.5*sin（w*t）;%聲源發(fā)出的信號(hào) u1=1.5*sin（w*t1）;%聽者接收到的信號(hào) sound（u，10000）;pause（5）;sound（u1，10000）;%將原信號(hào)和接收到的信號(hào)恢復(fù)為聲音 figure（1） subplot（2，1，1） plot（u）;axis（［11500，11700，-2.2，2.2］） subplot（2，1，2） plot（u1，′r′）;axis（［5500，7800，-2.2，2.2］）第二十七頁，共五十六頁，編輯于2023年，星期二10.3聲學(xué)信號(hào)的濾波特性試驗(yàn) 10.3.1聲音濾波的應(yīng)用 數(shù)字信號(hào)處理課程中詳細(xì)地討論了各種數(shù)字濾波器的設(shè)計(jì)方法，本例是數(shù)字濾波器應(yīng)用在聲學(xué)領(lǐng)域的一個(gè)例子。原始的聲音由三個(gè)單音頻組成，采用低通濾波器，可以將高音頻的信號(hào)濾除；采用帶通濾波器，可以只讓某一個(gè)頻率的音頻信號(hào)通過。編程的要點(diǎn)是在正確地選擇采樣頻率的基礎(chǔ)上，正確地選取通帶和阻帶的截止頻率，以保證濾掉和保留下預(yù)想的頻率。第二十八頁，共五十六頁，編輯于2023年，星期二 10.3.2聲音濾波的仿真試驗(yàn) 下面是一段聲音濾波的仿真試驗(yàn)的程序。300Hz、500Hz、1000Hz的三個(gè)正弦波信號(hào)疊加構(gòu)成本例的原信號(hào)，通過低通濾波器后濾掉了1000Hz的信號(hào)。圖10-6中的上面兩幅圖，是表示低通濾波器傳輸特性的幅頻圖和相頻圖。通過帶通濾波器后只剩下1000Hz的信號(hào)。圖10-6中的下面兩幅圖，是表示帶通濾波器傳輸特性的幅頻圖和相頻圖。第二十九頁，共五十六頁，編輯于2023年，星期二圖10-6濾波器的幅頻特性和相頻特性(上，低通；下，帶通）第三十頁，共五十六頁，編輯于2023年，星期二 圖10-7所示是信號(hào)的頻譜特性。從上到下依次是三個(gè)音頻信號(hào)；通過低通濾波器后，1000Hz的信號(hào)濾掉了，僅剩300Hz和500Hz的信號(hào)；通過中心頻率為1000Hz帶通濾波器后，僅剩1000Hz的信號(hào)。程序運(yùn)行結(jié)束時(shí)，可以聽到依次發(fā)出的三種信號(hào)的聲音。第三十一頁，共五十六頁，編輯于2023年，星期二 圖10-7信號(hào)的頻譜特性（上，原信號(hào)；中，通過低通；下，通過帶通）第三十二頁，共五十六頁，編輯于2023年，星期二 程序10-3 t=1/10000:1/10000:3;%設(shè)置三個(gè)頻率的正弦信號(hào) %300Hz，500Hz，1000Hz s1=sin（2*pi*300*t）; s2=sin（2*pi*500*t）; s3=1.2*sin（2*pi*1000*t）; s=s1+s2+s3;%三個(gè)正弦波信號(hào)疊加 figure（1）第三十三頁，共五十六頁，編輯于2023年，星期二 subplot（3，1，1）; ss=fft（s，4096）;SS=abs（ss（1:1:2049））;%求合成信號(hào)的頻譜 k1=0:2048;w1=（1/.4096）*k1;%取0...Fs/2的部分 plot（w1，SS）;grid%畫頻譜圖 axis（［0，2000，-100，2200］） title（′原信號(hào)頻譜圖′） %通過低通濾波器 ws1=1000;%設(shè)計(jì)一個(gè)通帶為600Hz、阻帶為1000Hz的低通濾波器第三十四頁，共五十六頁，編輯于2023年，星期二

wp1=600; wc=5000; wp=wp1/wc;ws=ws1/wc; ［n，Wn］=cheb2ord（wp，ws，1，30）%切比雪夫II型濾波器設(shè)計(jì) ［b，a］=cheby2（n，30，Wn）; subplot（3，1，2） sb=3*filter（b，a，s）;%合成信號(hào)通過低通濾波器 ssb=fft（sb，4096）;SSb=abs（ssb（1:1:2049））;%求頻譜 k1=0:2048;w1=（1/.4096）*k1;%畫頻譜圖 plot（w1，SSb）;grid axis（［0，2000，-100，6500］） title（′經(jīng)過低通濾波器后的信號(hào)頻譜圖′）

第三十五頁，共五十六頁，編輯于2023年，星期二

%通過帶通濾波器 ws1=［6001400］;%設(shè)計(jì)一個(gè)阻帶為600~1400Hz、 wp1=［9801020］;%通帶為980~1020Hz的帶通濾波器 wc=5000; wp=wp1/wc;ws=ws1/wc; ［n，Wn］=cheb1ord（wp，ws，4，45，′s′）%切比雪夫I型濾波器 ［bd，ad］=cheby1（n，4，Wn，′bandpass′）; subplot（3，1，3） sd=3*filter（bd，ad，s）;%合成信號(hào)通過低通濾波器第三十六頁，共五十六頁，編輯于2023年，星期二 ssd=fft（sd，4096）;SSd=abs（ssd（1:1:2049））;%求頻譜 k1=0:2048;w1=（1/.4096）*k1; plot（w1，SSd）;grid%畫頻譜 axis（［0，2000，-100，4000］） title（′經(jīng)過帶通濾波器后的信號(hào)頻譜圖′） figure（2）%畫低通濾波器的幅頻、相頻圖 freqz（b，a，4096，10000）;axis（［0，4000，-70，1］）第三十七頁，共五十六頁，編輯于2023年，星期二 figure（3）%畫帶通濾波器的幅頻、相頻圖 freqz（bd，ad，4096，10000）;axis（［0，4000，-70，1］） %依次用揚(yáng)聲器播放合成音、通過低通后、通過帶通后的聲音 sound（s，10000）;pause（5）;sound（sb，10000）; pause（5）;sound（sd，10000）第三十八頁，共五十六頁，編輯于2023年，星期二10.4交混回響 10.4.1交混回響的應(yīng)用 在一間有若干個(gè)揚(yáng)聲器的禮堂里，處在某一個(gè)特定位置聽到的聲音是什么效果，可以用仿真的方法來實(shí)現(xiàn)。為了簡化問題，假定墻面上有吸音的材料，即忽略了墻壁的反射。編程仿真的思想是：電信號(hào)在電線里傳輸?shù)臅r(shí)間可以忽略，但是當(dāng)相距若干米的揚(yáng)聲器發(fā)出的聲音到達(dá)聽者的位置時(shí)，不同揚(yáng)聲器的聲音經(jīng)過了不同的路程，路程的差別產(chǎn)生相位差，所有帶有不同相位的聲音的疊加就產(chǎn)生了交混回響的效果。第三十九頁，共五十六頁，編輯于2023年，星期二 10.4.2交混回響的仿真試驗(yàn) 圖10-8所示是交混回響仿真試驗(yàn)中的揚(yáng)聲器布置圖，A、B、C是三個(gè)接到同一聲源的揚(yáng)聲器，它距離聽者的水平和垂直距離標(biāo)注在圖上。下面是用M文件編制的程序，它的運(yùn)行繪出了原聲和交混回響的聲音時(shí)域波形。圖10-9所示是聲音信號(hào)的時(shí)域圖，并且依次播放出原聲和交混回響的聲音效果。第四十頁，共五十六頁，編輯于2023年，星期二圖10-8交混回響仿真試驗(yàn)中的揚(yáng)聲器布置第四十一頁，共五十六頁，編輯于2023年，星期二圖10-9聲音信號(hào)的時(shí)域圖（上，原信號(hào)；下，聽到的交混回響的信號(hào)）第四十二頁，共五十六頁，編輯于2023年，星期二

程序10-4 loadsoud11%讀入聲音數(shù)據(jù)文件soud11 y=soud11;fs=8000; N1=8600; v=330;%聲音速度 b1=（20^2+45^2）^.5;%計(jì)算聽者到揚(yáng)聲器B的距離 yb1=［zeros（1，(fs*b1/v）），y′］′;%聽者聽到揚(yáng)聲器B的聲音 yb=yb1（1:N1）; a1=（20^2+15^2）^.5;%計(jì)算聽者到揚(yáng)聲器A的距離 ya1=［zeros（1，(fs*a1/v）），y′］′;%聽者聽到揚(yáng)聲器A的聲音第四十三頁，共五十六頁，編輯于2023年，星期二 ya=ya1（1:N1）; c1=（35^2+65^2）^.5;%計(jì)算聽者到揚(yáng)聲器C的距離 yc1=［zeros（1，(fs*c1/v）），y′］′;%聽者聽到揚(yáng)聲器C的聲音 yc=yc1（1:N1）; ye=1.6*ya+2*yb+1.2*yc;%聽者聽到三個(gè)揚(yáng)聲器的聲音疊加 subplot（211）;plot（y）;grid subplot（212）;plot（ye）;grid sound（y，fs）;pause（2）;sound（1.2*ye，fs)%播放原聲后播放合成的聲音第四十四頁，共五十六頁，編輯于2023年，星期二10.5短時(shí)傅立葉變換 10.5.1短時(shí)傅立葉變換的基本原理 傅立葉變換將信號(hào)系統(tǒng)的時(shí)間域和頻率域的特性聯(lián)系起來，這一工具在處理平穩(wěn)信號(hào)時(shí)，因?yàn)樾盘?hào)的統(tǒng)計(jì)特性（平均值、方差以及頻率特性）都不是時(shí)間的函數(shù)，因此，在求其頻譜特性時(shí)，應(yīng)用了從負(fù)無窮到正無窮的時(shí)間積分。(10-2）第四十五頁，共五十六頁，編輯于2023年，星期二 當(dāng)信號(hào)是一個(gè)非平穩(wěn)信號(hào)時(shí)，它的頻率特性是時(shí)間的函數(shù)，上面的方法就不能正確地反映出信號(hào)的特性，因?yàn)樵诓煌臅r(shí)段內(nèi)有不同的頻率特性。應(yīng)用短時(shí)傅立葉變換：

(10-3） 可以正確地反映出非平穩(wěn)信號(hào)的頻率特性，此時(shí)的頻譜特性已經(jīng)是時(shí)間和頻率的二維函數(shù)。 公式中的W可以看成是一個(gè)窗函數(shù)，短時(shí)傅立葉變換的基本思路是把非平穩(wěn)信號(hào)分割為若干小的時(shí)段，在每個(gè)小時(shí)段內(nèi)把信號(hào)看成是平穩(wěn)的。各個(gè)時(shí)段的頻率特性的排列，就構(gòu)成了信號(hào)的時(shí)頻特性。時(shí)窗愈窄，它的時(shí)間分辨率愈高。第四十六頁，共五十六頁，編輯于2023年，星期二 10.5.2短時(shí)傅立葉變換的仿真試驗(yàn) Simulink仿真模型是參考MATLAB中的Toolbox＼dspblks＼dspdemos＼dspstfft構(gòu)建的，并且設(shè)置了ToWaveDevice（揚(yáng)聲器）和WaterfallScope（瀑布圖儀），如圖10-10所示。第四十七頁，共五十六頁，編輯于2023年，星期二圖10-10顯示聲音信號(hào)(用于短時(shí)傅立葉變換)的仿真系統(tǒng)模型第四十八頁，共五十六頁，編輯于2023年，星期二

下面程序說明：首先調(diào)出聲音文件Soud11；然后運(yùn)行mdl文件“dspstfftme3”，把“dspstfftme3”的聲音信號(hào)的數(shù)據(jù)用ToWorkspace模塊收集，變量名為sA，將數(shù)據(jù)序列sA作分段傅立葉變換，目的是將聲音信號(hào)的時(shí)間頻率特性用三維圖形表示出來。 程序10-5 loadsoud11 sim（′dspstfftme3′） ［S，F(xiàn)，T］=specgram（sA，128，1/8000，16）;%求解短時(shí)傅立葉變換%將數(shù)據(jù)序列sA作N點(diǎn)傅立葉變換，F(xiàn)s，Nbit的意義同前%S表示為以F（頻率），T(時(shí)間）為自變量的函數(shù)第四十九頁，共五十六頁，編輯于2023年，星期二 X=T*