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

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

10.1.1回聲的基本原理

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

(1)聲源距離聽者有多遠。

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

(3)聲音的反射體在反射聲音的過程中對聲音的衰減有多大。(假設(shè)聽者在聲源與反射體之間的某個位置上。)10.1.2回聲的仿真試驗

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

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

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

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

程序10-1

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

sim（′echo11′）

10.2多普勒效應(yīng)

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

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

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

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

圖10-4所示是M文件程序運行結(jié)束時繪出的聽者接收到信號的頻率變化曲線以及用信號處理工具箱函數(shù)spectrogram繪制的時間頻率圖。圖10-5所示是某個時刻圖10-3中的頻譜儀上顯示的功率譜。圖10-3顯示聲音信號（多普勒效應(yīng)）的仿真系統(tǒng)模型圖10-4聽者接收到的信號的頻率變化曲線圖10-5頻譜儀顯示（左圖）和時間頻率圖表示（右圖）程序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);%計算聲源與聽者的距離

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ā)出的信號

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

sound(signal,fs);%播放聽者接收到的信號

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′);%作出時間頻率圖

axis(［06,1000,spectro3000］)

10.3聲學信號的濾波特性試驗

10.3.1聲音濾波的應(yīng)用

數(shù)字信號處理課程中詳細討論了各種數(shù)字濾波器的設(shè)計方法，本例是數(shù)字濾波器應(yīng)用在聲學領(lǐng)域的一個例子。原始的聲音由三個單音頻組成，采用低通濾波器，可以將高音頻的信號濾除；而采用帶通濾波器，可以只讓某一個頻率的音頻信號通過。編程的要點是在正確地選擇采樣頻率的基礎(chǔ)上，正確地選取通帶和阻帶的截止頻率，以保證濾掉和保留下預(yù)想的頻率。10.3.2聲音濾波的仿真試驗

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

t=1/10000:1/10000:3;%設(shè)置三個頻率的正弦信號分別為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;%三個正弦波信號疊加

figure（1）

subplot（3，1，1）;

ss=fft（s，4096）;SS=abs（ss（1:1:2049）);%求合成信號的頻譜

k1=0:2048;w1=（1/.4096）*k1;%取0…Fs/2的部分

plot（w1，SS）;grid %畫頻譜圖

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

title（′原信號頻譜圖′）%通過低通濾波器

ws1=1000;

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

wc=5000;

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

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

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

subplot（3，1，2）

sb=3*filter（b，a，s）;%合成信號通過低通濾波器

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)過低通濾波器后的信號頻譜圖′）

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

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）;%合成信號通過低通濾波器

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)過帶通濾波器后的信號頻譜圖′）

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

freqz（b，a，4096，10000）;

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

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

freqz（bd，ad，4096，10000）;

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

%依次用揚聲器播放合成音、通過低通后、通過帶通后的聲音

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

pause（5）;sound（sd，10000）

10.4交混回響

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

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

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

圖10-8所示是交混回響仿真試驗中的揚聲器布置圖，A、B、C是三個接到同一聲源的揚聲器，它距離聽者的水平和

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

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

y=sound11;fs=8000;

N1=8600;

v=330;%聲音的速度

b1=（20^2+45^2）^.5;%計算聽者到揚聲器B的距離

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

%聽者聽到揚聲器B的聲音

yb=yb1（1:N1）;

a1=（20^2+15^2）^.5;%計算聽者到揚聲器A的距離

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

%聽者聽到揚聲器A的聲音

ya=ya1（1:N1）;

c1=（35^2+65^2）^.5;%計算聽者到揚聲器C的距離

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

%聽者聽到揚聲器C的聲音

yc=yc1（1:N1）;

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

%聽者聽到三個揚聲器的聲音疊加

subplot（211）;plot（y）;grid

subplot（212）;plot（ye）;grid

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

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

10.5短時傅立葉變換

10.5.1短時傅立葉變換的基本原理

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

應(yīng)用短時傅立葉變換：(10-3）可以正確地反映出非平穩(wěn)信號的頻率特性，此時的頻譜特性已經(jīng)是時間和頻率的二維函數(shù)。公式中的W可以看成是一個窗函數(shù)，短時傅立葉變換的基本思路是把非平穩(wěn)信號分割為若干小的時段，在每個小時段內(nèi)把信號看成是平穩(wěn)的。各個時段的頻率特性的排列，就構(gòu)成了信號的時頻特性。時窗愈窄，它的時間分辨率愈高。10.5.2短時傅立葉變換的仿真試驗

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

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