《MATLAB仿真在通信與電子工程中的應(yīng)用》課件第10章

第10章聲學(xué)仿真試驗(yàn)10.1回聲試驗(yàn)10.2多普勒效應(yīng)10.3聲學(xué)信號(hào)的濾波特性試驗(yàn)10.4交混回響10.5短時(shí)傅立葉變換10.1回聲試驗(yàn)

10.1.1回聲的基本原理

回聲是一種物理現(xiàn)象。當(dāng)直達(dá)聲結(jié)束后，聲音經(jīng)過(guò)其他物體反射、折射后，又返回來(lái)收到的聲波，稱為回聲?；芈曈袝r(shí)泛指一切反射回來(lái)的聲音。當(dāng)傳到人耳的直達(dá)聲和回聲之間的時(shí)間差在1/20s以上時(shí)，可以很清楚地把它們區(qū)分開(kāi)。在山谷或大廳中，常有回聲現(xiàn)象發(fā)生。尤其是在影劇院里，回聲往往會(huì)妨礙聽(tīng)音，所以建筑師必須考慮消除回聲的影響。相反，直達(dá)聲和回聲之間的時(shí)間差在1/20s以下時(shí)，回聲不但不妨礙聽(tīng)音，反而使聲音更加響亮。應(yīng)用MATLAB軟件可以方便地進(jìn)行仿真試驗(yàn)?；芈暚F(xiàn)象取決于幾個(gè)重要的物理參數(shù)：

(1)聲源距離聽(tīng)者有多遠(yuǎn)。

(2)聲音產(chǎn)生往復(fù)反射的兩個(gè)反射體的距離有多遠(yuǎn)。

(3)聲音的反射體在反射聲音的過(guò)程中對(duì)聲音的衰減有多大。(假設(shè)聽(tīng)者在聲源與反射體之間的某個(gè)位置上。)10.1.2回聲的仿真試驗(yàn)

圖10-1所示是回聲的仿真試驗(yàn)系統(tǒng)框圖。圖中的1800的延遲代表了聲音在兩個(gè)反射體之間走一個(gè)來(lái)回的延遲（1800/Fs秒）。它的值愈大，反映了兩個(gè)反射體的距離愈大。圖中0.8的增益代表了聲音在兩個(gè)反射體之間走一個(gè)來(lái)回的損耗，它愈小，回聲衰減愈快。圖中0.9的增益代表了傳輸損耗。圖10-2所示是仿真試驗(yàn)時(shí)示波器上顯示的波形圖，上圖是原聲，下圖是往復(fù)反射的回聲波形圖。如果計(jì)算機(jī)裝有聲卡和音箱，可以聽(tīng)到揚(yáng)聲器發(fā)出的回聲聲音。放大器增益和整數(shù)延遲可以看圖標(biāo)注。圖10-1回聲仿真試驗(yàn)系統(tǒng)框圖圖10-2回聲仿真試驗(yàn)時(shí)域圖回聲的仿真試驗(yàn)的模型參數(shù)設(shè)置分析：回聲的仿真試驗(yàn)的模型參考了MATLAB中Toolbox＼dspblks＼dspdemos＼dspafxr的DEMOS例子，在仿真系統(tǒng)中設(shè)置了Reshape（整型器），主要參數(shù)見(jiàn)表10-1；ToWaveDevice（揚(yáng)聲器）的主要參數(shù)見(jiàn)表10-2；FromWorkspace（聲源）的主要參數(shù)見(jiàn)表10-3；它以幀的形式輸出。整型器的設(shè)置目的是將信號(hào)轉(zhuǎn)變?yōu)槭静ㄆ髂軌蚪邮艿囊痪S數(shù)據(jù)流。在10.1、10.4、10.5節(jié)中都需要一段聲音文件來(lái)進(jìn)行仿真試驗(yàn)?？梢杂靡韵聝煞N方法來(lái)獲得聲音文件：

(1)調(diào)出MATLAB中的聲音文件“l(fā)oadmtlb.mat”。

(2)應(yīng)用6.2.2節(jié)介紹的聲音采集的方法，建立一個(gè)聲音文件。譬如文件名是Sound11.mat。

通過(guò)下列程序調(diào)入仿真的聲音數(shù)據(jù)并啟動(dòng)仿真。Echo11是圖10-1所示的仿真系統(tǒng)的名字（文件名）。

程序10-1

loadsound11%讀入聲音數(shù)據(jù)文件sound11

sim（′echo11′）

10.2多普勒效應(yīng)

10.2.1多普勒效應(yīng)的基本原理

當(dāng)波源或觀察者相對(duì)于媒質(zhì)運(yùn)動(dòng)時(shí)，或者說(shuō)波源和觀察者有相對(duì)運(yùn)動(dòng)時(shí)，觀察者接受到的振動(dòng)頻率與波源振動(dòng)頻率不同的現(xiàn)象，稱為多普勒效應(yīng)。

當(dāng)聽(tīng)者與聲源的相對(duì)速度大到與聲速可以比擬時(shí)，就可以明顯地感覺(jué)到聲音頻率的變化。譬如坐火車時(shí)，聽(tīng)到相向開(kāi)來(lái)的火車汽笛聲，當(dāng)火車由遠(yuǎn)而近開(kāi)來(lái)時(shí)，汽笛聲的音調(diào)變高；由近而遠(yuǎn)離去時(shí)，汽笛聲的音調(diào)變低。這是日常生活中的一個(gè)多普勒效應(yīng)的例子。在天文、通信等領(lǐng)域還有眾多的例子。以下的公式描述了該現(xiàn)象的各個(gè)物理量之間的定量關(guān)系:(10-1）其中,f0是聲源發(fā)出的聲音的頻率，v是聽(tīng)者與聲源的相對(duì)運(yùn)動(dòng)速度，θ為速度矢量與聲源和聽(tīng)者的連線夾角，vs為聲音在空氣中傳播的速度，f是聽(tīng)者聽(tīng)到的聲音頻率。10.2.2多普勒效應(yīng)的仿真試驗(yàn)

下面是一段用MATLAB的M文件依據(jù)公式（10-1）原理編寫(xiě)的程序。它會(huì)產(chǎn)生一個(gè)文件名為dopp.wav的表現(xiàn)多普勒效應(yīng)的聲音文件，描述火車向一個(gè)距鐵路30m（垂直距離）、距火車510m（水平距離）的聽(tīng)者開(kāi)來(lái)時(shí)他聽(tīng)到的聲音。將u1送入圖10-3所示的Simulink仿真系統(tǒng)，再一次聽(tīng)到該聲音，并且看到用頻譜儀表現(xiàn)的聲音頻率隨時(shí)間變化的情況。表10-4和表10-5給出了仿真系統(tǒng)中兩模塊的主要參數(shù)設(shè)置。

圖10-4所示是M文件程序運(yùn)行結(jié)束時(shí)繪出的聽(tīng)者接收到信號(hào)的頻率變化曲線以及用信號(hào)處理工具箱函數(shù)spectrogram繪制的時(shí)間頻率圖。圖10-5所示是某個(gè)時(shí)刻圖10-3中的頻譜儀上顯示的功率譜。圖10-3顯示聲音信號(hào)（多普勒效應(yīng)）的仿真系統(tǒng)模型圖10-4聽(tīng)者接收到的信號(hào)的頻率變化曲線圖10-5頻譜儀顯示（左圖）和時(shí)間頻率圖表示（右圖）程序10-2

%多普勒效應(yīng)程序

x0=150;%m

v=50;%車速m/s

y0=30;%m

c=330;%聲音在空氣中的傳播速度m/s

f0=2000;%聲音的頻率Hz

fs=8000;%采樣率

t=0:1/fs:6;

x_t=－x0+v.*t;

y=y0;

r=sqrt(x_t.^2+y.^2);%計(jì)算聲源與聽(tīng)者的距離

costheta=x_t./r;

f=f0./(1+v.*costheta./c);

fmin=min(f)

fmax=max(f)

voltage=(f-fmin)./(fmax-fmin)*2－1;%歸一化調(diào)頻電壓在-1~+1之間

signal=0.5*vco(voltage,［fminfmax］,fs);

u=sin(2*pi*f0*t);%聲源發(fā)出的信號(hào)

sound(u,fs);pause(5);%播放聲源發(fā)出的信號(hào)

sound(signal,fs);%播放聽(tīng)者接收到的信號(hào)

wavwrite(signal,fs,′dopp.wav′);

figure(1);

plot(t,f);xlabel(′Time/s′);

ylabel(′Freq/Hz′);

figure(2);

gram(signal,kaiser(256,5),220,512,fs,′yaxis′);%作出時(shí)間頻率圖

axis(［06,1000,spectro3000］)

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)通過(guò)。編程的要點(diǎn)是在正確地選擇采樣頻率的基礎(chǔ)上，正確地選取通帶和阻帶的截止頻率，以保證濾掉和保留下預(yù)想的頻率。10.3.2聲音濾波的仿真試驗(yàn)

下面是一段聲音濾波的仿真試驗(yàn)的程序。300Hz、500Hz、1000Hz的三個(gè)正弦波信號(hào)疊加構(gòu)成本例的原信號(hào)，通過(guò)低通濾波器后濾掉了1000Hz的信號(hào)。圖10-6中上面的兩幅圖分別是表示低通濾波器傳輸特性的幅頻圖和相頻圖。通過(guò)帶通濾波器后只剩下1000Hz的信號(hào)。圖10-6中下面的兩幅圖分別是表示帶通濾波器傳輸特性的幅頻圖和相頻圖。圖10-6濾波器的幅頻特性和相頻特性（左圖：低通；右圖:帶通）圖10-7所示是信號(hào)的頻譜特性。從上到下依次是：三個(gè)音頻信號(hào)；通過(guò)低通濾波器后，1000Hz的信號(hào)濾掉了，僅剩300Hz和500Hz的信號(hào)；通過(guò)中心頻率為1000Hz帶通濾波器后，僅剩1000Hz的信號(hào)。程序運(yùn)行結(jié)束時(shí)，可以聽(tīng)到依次發(fā)出的三種信號(hào)的聲音。圖10-7信號(hào)的頻譜特性（上圖：原信號(hào)；中圖：通過(guò)低通；下圖：通過(guò)帶通）程序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）

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 %畫(huà)頻譜圖

axis（［0，2000，-100，2200］）

title（′原信號(hào)頻譜圖′）%通過(guò)低通濾波器

ws1=1000;

%設(shè)計(jì)一個(gè)通帶為600Hz、阻帶為1000Hz的低通濾波器wp1=600;

wc=5000;

wp=wp1/wc;ws=ws1/wc;

［n，Wn］=cheb2ord（wp，ws，1，30）%切比雪夫Ⅱ型濾波器設(shè)計(jì)

［b，a］=cheby2（n，30，Wn）;

subplot（3，1，2）

sb=3*filter（b，a，s）;%合成信號(hào)通過(guò)低通濾波器

ssb=fft（sb，4096）;SSb=abs（ssb（1:1:2049））;

%求頻譜

k1=0:2048;w1=（1/.4096）*k1;%畫(huà)頻譜圖

plot（w1，SSb）;grid

axis（［0，2000，-100，6500］）

title（′經(jīng)過(guò)低通濾波器后的信號(hào)頻譜圖′）

%通過(guò)帶通濾波器的情況；下面是帶通濾波器設(shè)計(jì)

ws1=［6001400］;

%小于等于600Hz、大于等于1400Hz的區(qū)間為阻帶

wp1=［9801020］;%980~1020Hz區(qū)間為通帶

wc=5000;

wp=wp1/wc;ws=ws1/wc;［n，Wn］=cheb1ord（wp，ws，4，45，′s′）%切比雪夫Ⅰ型濾波器

［bd，ad］=cheby1（n，4，Wn，′bandpass′）;

subplot（3，1，3）

sd=3*filter（bd，ad，s）;%合成信號(hào)通過(guò)低通濾波器

ssd=fft（sd，4096）;SSd=abs（ssd（1:1:2049））;%求頻譜

k1=0:2048;w1=（1/.4096）*k1;

plot（w1，SSd）;grid%畫(huà)頻譜

axis（［0，2000，－100，4000］）

title（′經(jīng)過(guò)帶通濾波器后的信號(hào)頻譜圖′）

figure（2）%畫(huà)低通濾波器的幅頻、相頻圖

freqz（b，a，4096，10000）;

axis（［0，4000，－70，1］）

figure（3）%畫(huà)帶通濾波器的幅頻、相頻圖

freqz（bd，ad，4096，10000）;

axis（［0，4000，－70，1］）

%依次用揚(yáng)聲器播放合成音、通過(guò)低通后、通過(guò)帶通后的聲音

sound（s，10000）;pause（5）;sound（sb，10000）;

pause（5）;sound（sd，10000）

10.4交混回響

10.4.1交混回響的應(yīng)用

在一間有若干個(gè)揚(yáng)聲器的禮堂里，處在某一個(gè)特定位置聽(tīng)到的聲音是什么效果，可以用仿真的方法來(lái)實(shí)現(xiàn)。為了簡(jiǎn)化問(wèn)題，假定墻面上有吸音的材料，即忽略了墻壁的反射。編程仿真的思想是：電信號(hào)在電線里的傳輸時(shí)間可以忽略，但是當(dāng)相距若干米的揚(yáng)聲器發(fā)出的聲音到達(dá)聽(tīng)者的位置時(shí)，不同揚(yáng)聲器的聲音經(jīng)過(guò)了不同的路程，路程的差別產(chǎn)生相

位差，所有帶有不同相位的聲音的疊加就產(chǎn)生了交混回響的效果。10.4.2交混回響的仿真試驗(yàn)

圖10-8所示是交混回響仿真試驗(yàn)中的揚(yáng)聲器布置圖，A、B、C是三個(gè)接到同一聲源的揚(yáng)聲器，它距離聽(tīng)者的水平和

垂直距離標(biāo)注在圖上。下面是用M文件編制的程序，它的運(yùn)行繪制出了原聲和交混回響的聲音時(shí)域波形。圖10-9所示是聲音信號(hào)的時(shí)域圖，并且依次播放出原聲和交混回響的聲音效果。圖10-8交混回響仿真試驗(yàn)中的揚(yáng)聲器布置圖10-9聲音信號(hào)的時(shí)域圖（上圖：原信號(hào)；下圖：聽(tīng)到的交混回響的信號(hào)）程序10-4

loadsound11%讀入聲音數(shù)據(jù)文件sound11

y=sound11;fs=8000;

N1=8600;

v=330;%聲音的速度

b1=（20^2+45^2）^.5;%計(jì)算聽(tīng)者到揚(yáng)聲器B的距離

yb1=［zeros（1，（fs*b1/v）），y′］′;

%聽(tīng)者聽(tīng)到揚(yáng)聲器B的聲音

yb=yb1（1:N1）;

a1=（20^2+15^2）^.5;%計(jì)算聽(tīng)者到揚(yáng)聲器A的距離

ya1=［zeros（1，（fs*a1/v）），y′］′;

%聽(tīng)者聽(tīng)到揚(yáng)聲器A的聲音

ya=ya1（1:N1）;

c1=（35^2+65^2）^.5;%計(jì)算聽(tīng)者到揚(yáng)聲器C的距離

yc1=［zeros（1，（fs*c1/v）），y′］′;

%聽(tīng)者聽(tīng)到揚(yáng)聲器C的聲音

yc=yc1（1:N1）;

ye=1.6*ya+2*yb+1.2*yc;

%聽(tīng)者聽(tīng)到三個(gè)揚(yáng)聲器的聲音疊加

subplot（211）;plot（y）;grid

subplot（212）;plot（ye）;grid

sound（y，fs）;pause（2）;sound（1.2*ye，fs）

%播放原聲后播放合成的聲音

10.5短時(shí)傅立葉變換

10.5.1短時(shí)傅立葉變換的基本原理

傅立葉變換將信號(hào)系統(tǒng)的時(shí)間域和頻率域的特性聯(lián)系起來(lái)，這一工具在處理平穩(wěn)信號(hào)時(shí)，由于信號(hào)的統(tǒng)計(jì)特性（平均值、方差以及頻率特性）都不是時(shí)間的函數(shù),因此，在求其頻譜特性時(shí)，應(yīng)用了從負(fù)無(wú)窮到正無(wú)窮的時(shí)間積分。(10-2）當(dāng)信號(hào)是一個(gè)非平穩(wěn)信號(hào)時(shí)，它的頻率特性是時(shí)間的函數(shù)，上面的方法就不能正確地反映出信號(hào)的特性，因?yàn)樵诓煌臅r(shí)段內(nèi)信號(hào)有不同的頻率特性。

應(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í)間分辨率愈高。10.5.2短時(shí)傅立葉變換的仿真試驗(yàn)

Simulink仿真模型是參考MATLAB中的Toolbox＼dspblks＼dspdemos＼dspstfft構(gòu)建的，并且設(shè)置了ToWaveDevice（揚(yáng)聲器）和WaterfallScope（瀑布圖儀），如圖10-10所示。

下面的程序說(shuō)明：首先調(diào)出聲音文件Sound11；然后運(yùn)行mdl文件“dspstfftme3”，把“dspstfftme3”的聲音信號(hào)的數(shù)據(jù)用ToWorkspace模塊收集，變量名為sA，