吳禮福 陶明明 郭業(yè)才

(1 南京信息工程大學(xué)電子與信息工程學(xué)院 南京 210044)

(2 江蘇省大氣環(huán)境與裝備技術(shù)協(xié)同創(chuàng)新中心 南京 210044)

0 引言

室內(nèi)混響是聲音從聲源發(fā)出后由于不斷被室內(nèi)表面反射、吸收而逐漸衰減的現(xiàn)象[1-3]。如果混響過(guò)大，會(huì)使聲音變得渾濁和雜亂，需要使用去混響技術(shù)[3-4]，例如，在免提通信系統(tǒng)中使用去混響提高語(yǔ)音的可懂度和語(yǔ)音識(shí)別的準(zhǔn)確率。然而，如果缺乏適當(dāng)?shù)幕祉?，聲音或者音?lè)聽(tīng)起來(lái)會(huì)非?！案稍铩?，聽(tīng)感不舒服，此時(shí)可以借助人工混響技術(shù)塑造聲音的空間感，使聲音聽(tīng)起來(lái)更自然真實(shí)。因此，人工混響廣泛應(yīng)用于音樂(lè)、電影和虛擬現(xiàn)實(shí)中，以美化音樂(lè)、聲音的音色，進(jìn)行藝術(shù)的再創(chuàng)造，產(chǎn)生特殊的音效[5]。

人工混響最早是在20世紀(jì)20年代的廣播系統(tǒng)[5-6]中通過(guò)模擬方法實(shí)現(xiàn)的，它將無(wú)混響信號(hào)傳輸?shù)交祉懎h(huán)境(一個(gè)專(zhuān)門(mén)建造的回音室)中采集而得到。但是模擬方法對(duì)外部聲學(xué)或機(jī)械的變化擾動(dòng)較敏感，并且需要專(zhuān)業(yè)知識(shí)來(lái)維護(hù)和調(diào)整混響效果。隨著數(shù)字信號(hào)處理器的發(fā)展，人工混響可以通過(guò)數(shù)字方法實(shí)現(xiàn)，其中混響被看作是一個(gè)線性時(shí)不變的過(guò)程[5]。數(shù)字人工混響方法主要有3 種基本類(lèi)型(盡管可以使用它們的組合)：反饋延遲網(wǎng)絡(luò)方法(Feedback delay network，F(xiàn)DN)、計(jì)算聲學(xué)方法和卷積方法[5]。反饋延遲網(wǎng)絡(luò)將輸入信號(hào)(干凈無(wú)混響)延遲、濾波并根據(jù)參數(shù)化混響特性沿著多個(gè)路徑反饋給前端，疊加后得到混響信號(hào)。該方法在音樂(lè)技術(shù)領(lǐng)域使用較多，但是通常只能產(chǎn)生需要的感知或者藝術(shù)效果[7]。二是基于計(jì)算聲學(xué)的方法，將輸入信號(hào)模擬聲能在幾何模型中傳播，從而得到混響信號(hào)，它通常用于需要準(zhǔn)確度的聲學(xué)設(shè)計(jì)和分析方案，如預(yù)測(cè)聲學(xué)建模和聲學(xué)空間的計(jì)算機(jī)輔助設(shè)計(jì)，但它們難以用于實(shí)時(shí)性較高的場(chǎng)合[5，8]。三是卷積方法，它將干凈無(wú)混響信號(hào)與房間脈沖響應(yīng)(Room impulse response，RIR)卷積得到混響信號(hào)，相當(dāng)于用RIR 對(duì)干凈無(wú)混響信號(hào)進(jìn)行有限脈沖響應(yīng)(Finite impulse response，F(xiàn)IR)濾波[5]，如果選用的RIR 能夠逼近真實(shí)房間的RIR，則卷積方法的混響效果優(yōu)于FDN 和計(jì)算聲學(xué)方法[5]。但是卷積方法的一個(gè)直接缺點(diǎn)是它的運(yùn)算量較大，一部分運(yùn)算量來(lái)自于用高階RIR 與干凈無(wú)混響信號(hào)之間的卷積運(yùn)算，這部分運(yùn)算量可以通過(guò)快速卷積算法一定程度上加以緩解[9]；另一個(gè)來(lái)自于RIR 的獲取，RIR 可以是在真實(shí)的房間中通過(guò)測(cè)量獲得，并將其保存在內(nèi)存中以便實(shí)時(shí)實(shí)現(xiàn)，但根據(jù)所需混響特性尋找合適的房間是費(fèi)時(shí)且不容易找到，測(cè)量RIR 也需要耗費(fèi)相當(dāng)多的時(shí)間。因此，目前主要通過(guò)算法合成RIR，例如廣泛使用的鏡像源法(Image source method，ISM)[10]能產(chǎn)生與真實(shí)房間非常相似的RIR，但是鏡像源法計(jì)算房間脈沖響應(yīng)所需要的鏡像源的數(shù)量與RIR 的長(zhǎng)度(階數(shù))成立方比，與反射級(jí)數(shù)呈指數(shù)增長(zhǎng)，這種大的計(jì)算量給其實(shí)時(shí)應(yīng)用帶來(lái)較大困難。

為了獲得比反饋延遲網(wǎng)絡(luò)方法更好的混響效果，又能比鏡像源法有更少的計(jì)算量，從便于實(shí)時(shí)應(yīng)用的角度，本文研究了一種卷積法人工混響方法，將RIR 分解成早期混響和后期混響兩個(gè)部分：早期混響可以用具有相應(yīng)衰減因子的延時(shí)求和來(lái)表示；后期混響RIR 受到文獻(xiàn)[11]的啟發(fā)，它在頻域合成房間頻率響應(yīng)(Room frequency response，RFR)，再經(jīng)傅里葉反變換得到時(shí)域RIR，最后將干凈無(wú)混響信號(hào)與RIR卷積得到混響信號(hào)。

1 頻域合成房間頻率響應(yīng)的人工混響方法

由信號(hào)處理理論可知，對(duì)于線性時(shí)不變系統(tǒng)，它的輸出信號(hào)等于輸入信號(hào)與系統(tǒng)的單位脈沖響應(yīng)的卷積，

其中，y(t)為輸出信號(hào)，x(t)為輸入信號(hào)，h(t)為系統(tǒng)的單位脈沖響應(yīng)。若將房間看成一個(gè)線性時(shí)不變系統(tǒng)，房間內(nèi)的某一聲源發(fā)出的聲音經(jīng)過(guò)房間的“處理”傳輸?shù)浇邮仗帲邮盏降穆曇艟褪窍到y(tǒng)的輸出信號(hào)。房間脈沖響應(yīng)就相當(dāng)于系統(tǒng)的單位脈沖響應(yīng)。卷積法人工混響就是用RIRh(t)與干凈語(yǔ)音x(t)進(jìn)行卷積，得到混響信號(hào)y(t)。

h(t)可以分解成早期混響和后期混響兩個(gè)部分，即

其中，he(t)是早期房間脈沖響應(yīng)，hl(t)是后期房間脈沖響應(yīng)。早期與后期混響之間的轉(zhuǎn)換瞬間被稱(chēng)為混合時(shí)間t0[11]。

其中，C0為歸一化常數(shù)，依據(jù)文獻(xiàn)[11]的建議，本文選擇C0= 0.002，V是房間的體積(以m3為單位)，fs是采樣率(以Hz 為單位)，是向下取整運(yùn)算，式(3)中t0的單位為采樣點(diǎn)數(shù)。房間頻率響應(yīng)RFR 定義為RIR 的離散傅里葉變換，記為H(k)(k=0，1，···，T -1)，則有

其中，He(k)和Hl(k)分別為he(t)和hl(t)的離散傅里葉變換。早期混響可以用M階具有相應(yīng)衰減因子的延時(shí)求和來(lái)表示，在頻域內(nèi)可以表示為[11]

其中，ρk和τk分別為第k個(gè)延時(shí)的衰減因子和延時(shí)。

后期混響hl(t)的合成方法主要依據(jù)文獻(xiàn)[11]中的結(jié)論，由于后期混響近似以指數(shù)衰減，可以定義hl(t)的能量時(shí)間曲線(Power temporal profile，PTP)為

其中，

T60是混響時(shí)間，S是整個(gè)房間的面積(以m2為單位)，c是聲速(以m/s為單位)。

統(tǒng)計(jì)房間聲學(xué)理論指出Hl(k)是一個(gè)廣義平穩(wěn)的復(fù)數(shù)高斯隨機(jī)過(guò)程[11]，因此定義Hl(k)的自協(xié)方差函數(shù)γ(m)和功率譜密度φ(t)為

其中，F(xiàn)T{ }是離散傅里葉變換，( )*表示復(fù)數(shù)共軛，應(yīng)用Wiener—Khinchin 定理和文獻(xiàn)[11]中的證明，可以得到

考慮到模式積分的起止時(shí)間為2010年7月19日00時(shí)—22日00時(shí),在分析過(guò)程中,只取了20日00時(shí)—21日00時(shí)的降水?dāng)?shù)據(jù)進(jìn)行分析。實(shí)況降水?dāng)?shù)據(jù)是根據(jù)Micaps系統(tǒng)的資料用Grads插值得到(圖3a)。由圖3a可以看到,實(shí)況降水主要集中在41.4°～42.8°N遼寧中部和北部區(qū)域,雨帶基本上呈東北—西南走向。強(qiáng)降水中心位于區(qū)域122.3°～124.3°E之間,24 h最大降水量為279 mm。

由此可以看出，理論上Hl(k)的自協(xié)方差函數(shù)γ(m)是t0、τ、P0的函數(shù)，t0由房間體積V經(jīng)式(3)確定，τ由T60經(jīng)式(7)確定，式(7)～(10)表明P0由V、S、T60確定，因此γ(m)僅僅依賴(lài)于房間的混響時(shí)間、體積和面積。文獻(xiàn)[11]指出在頻域內(nèi)可以用自回歸滑動(dòng)平均(Autoregressive moving average，ARMA)模型來(lái)描述Hl(k)，即

其中，φ0=θ0= 1，z-1是延時(shí)因子，即z-1Hl(k)=Hl(k -1)，ε(k)是一個(gè)方差為σ2ε的復(fù)數(shù)高斯白噪聲。由式(14)～(16)可以將γ(m)和φ(t)表示為

其中，表示m對(duì)T求余，ψq滿(mǎn)足Θ(z-1)/Φ(z-1)。

式(14)表明如果能估計(jì)出φp和θq，則可以用ε(k)去激勵(lì)Θ(z)/Φ(z)得到Hl(k)。而ARMA 模型中的參數(shù)φp和θq可以從γ(m)中通過(guò)求解得到，思路是先求解AR 系數(shù)φp，再求解MA 系數(shù)θq，簡(jiǎn)要的計(jì)算步驟是(細(xì)節(jié)可以參見(jiàn)文獻(xiàn)[12])：首先，AR 參數(shù)φp可以由式(18)通過(guò)求解修正的Yule-Walker 方程得到；其次，利用φp和γ(m)可以得到γ′(m)=Φ*(z-1)Φ(z)γ(m)，這是由θq和σ2ε確定的Q階MA 模型，可以采用10Q階的AR 模型來(lái)近似，因而可以求解標(biāo)準(zhǔn)Yule-Walker 方程得到10Q個(gè)系數(shù)φ′p；最后，θq可以通過(guò)逼近10Q個(gè)系數(shù)φ′p求解得到。

至此，當(dāng)給定房間的幾何尺寸和混響時(shí)間后，可以計(jì)算出房間體積V和面積S，首先由式(13)計(jì)算得到γ(m)；其次由γ(m)通過(guò)求解修正的或標(biāo)準(zhǔn)Yule-Walker 方程得到ARMA 模型中的參數(shù)φp和θq；再次，用ε(k)去激勵(lì)Θ(z)/Φ(z)得到Hl(k)；最后，將Hl(k)和He(k)組合后經(jīng)過(guò)傅里葉反變換得到完整的房間脈沖響應(yīng)h(t)，將干凈語(yǔ)音x(t)與h(t)進(jìn)行卷積，得到人工混響信號(hào)。表1給出了完整的合成Hl(k)的步驟。

表1 基于ARMA 模型的后期混響頻率響應(yīng)的合成步驟Table1 Steps for synthesizing late reverberation frequency response with ARMA model

鏡像源法是以幾何聲學(xué)模型為理論基礎(chǔ)的經(jīng)典算法，可以較準(zhǔn)確地模擬聲場(chǎng)脈沖響應(yīng)，憑借理論簡(jiǎn)單、易于理解的優(yōu)點(diǎn)，鏡像源法已被視為一種合成房間脈沖響應(yīng)的基準(zhǔn)，但是它的一個(gè)明顯缺點(diǎn)是計(jì)算量大，通過(guò)鏡像源法計(jì)算房間脈沖響應(yīng)所需要的鏡像源的數(shù)量與房間脈沖響應(yīng)的長(zhǎng)度成立方比，與反射級(jí)數(shù)呈指數(shù)增長(zhǎng)，例如一個(gè)長(zhǎng)7 m、寬5 m、高4 m 的房間，有6 個(gè)反射面，如需要計(jì)算的脈沖長(zhǎng)度為512 ms，所考慮的最大的鏡像源級(jí)數(shù)為10 級(jí)，則總的鏡像源的數(shù)目為9.1×107個(gè)，雖然已有相關(guān)降低鏡像源法計(jì)算量的方法，但其計(jì)算量需求通常在實(shí)時(shí)應(yīng)用中難以滿(mǎn)足。

從表1可以看出，本文方法(后文以“ARMA”簡(jiǎn)記之)的主要計(jì)算量來(lái)自于步驟5 和步驟6，步驟5 從γ(m)中計(jì)算ARMA 模型的參數(shù)，采用經(jīng)典的Levinson-Durbin算法，其計(jì)算量大約是AR模型階數(shù)的平方關(guān)系，而AR 模型階數(shù)通常在10 階左右。例如，本文采用的ARMA 模型的階數(shù)設(shè)置為P= 7，Q= 2，因此其計(jì)算量是可以接受的；步驟6本質(zhì)上等效于一個(gè)無(wú)限脈沖響應(yīng)(Infinite impulse response，IIR)的濾波運(yùn)算，在實(shí)時(shí)系統(tǒng)中已經(jīng)經(jīng)常使用。另一方面，步驟5 中的運(yùn)算量可以進(jìn)一步采用“存儲(chǔ)量換取計(jì)算量”的方法，即可以事先設(shè)置足夠多的步驟1 中的參數(shù)以涵蓋不同的房間參數(shù)，然后離線計(jì)算出ARMA 模型中的參數(shù)φp和θq并保存在存儲(chǔ)器中，由于ARMA模型的階數(shù)設(shè)置通常都較低，其占有的存儲(chǔ)器容量可以忽略。因此，本文方法便于在實(shí)時(shí)人工混響系統(tǒng)中應(yīng)用。

2 實(shí)驗(yàn)與討論

為了驗(yàn)證本文方法的實(shí)際混響效果，首先采用真實(shí)房間測(cè)得的RIR 為參考，分別比較鏡像源法[10](以“ISM”簡(jiǎn)記)、反饋延遲網(wǎng)絡(luò)[7](以“FDN”簡(jiǎn)記)和本文方法的性能。由于鏡像源法適用于鞋盒型六面體房間，因此選用德國(guó)亞琛大學(xué)實(shí)測(cè)的房間脈沖響應(yīng)數(shù)據(jù)庫(kù)[13]中的會(huì)議室(meeting room)數(shù)據(jù)作為基準(zhǔn)，該會(huì)議室長(zhǎng)8 m、寬5 m、高3.1 m，分別在會(huì)議室的不同位置測(cè)試了20 條不同的雙耳RIR，圖1(a)給出了其中一個(gè)RIR，采用施羅德反向積分法[1]計(jì)算得到該會(huì)議室的混響時(shí)間為0.38 s，據(jù)此可以設(shè)置ISM、FDN 和ARMA 三種方法的主要參數(shù)，如表2所示。

表2 仿真實(shí)驗(yàn)中的基本參數(shù)設(shè)置Table2 Basic parameter settings in the simulation experiment

FDN 的原理就是將一段干凈無(wú)混響的語(yǔ)音依次通過(guò)延遲線、濾波器和反饋矩陣，最后得到混響信號(hào)。依據(jù)文獻(xiàn)[7]，表2中的FDN選用了16條延遲線，延遲線最大值為289，最小值為19。

圖1(b)～圖1(d)分別是采用ISM、 FDN 和ARMA 方法生成的RIR，從中可以明顯看出，就RIR后期混響部分的回聲密度而言，F(xiàn)DN 方法明顯低于ISM 和ARMA 方法，與真實(shí)房間的RIR 相比，差異明顯。而ISM 方法和ARMA 方法的后期混響部分的回聲密度要高于FDN 方法，與真實(shí)房間的接近。

圖1 4 種方法獲取的房間脈沖響應(yīng)Fig.1 Room impulse response obtained by four different methods

圖2給出了圖1中4個(gè)RIR的能量衰減曲線，從中可以看出ARMA 法的能量衰減曲線與ISM 的非常接近，兩者與真實(shí)房間能量衰減曲線的距離明顯小于FDN。另一方面，采用施羅德反向積分法計(jì)算得到ISM、FDN 和ARMA 法生成的RIR 的混響時(shí)間分別為0.39 s、0.44 s 和0.40 s，這與圖2的能量衰減曲線體現(xiàn)的趨勢(shì)一致，表明FDN 方法生成的RIR 與真實(shí)房間的差異較大。ARMA 法和ISM 法的結(jié)果接近，生成的RIR 計(jì)算得到的混響時(shí)間與設(shè)定值之間的誤差也較小。

圖2 4 種方法獲取的房間脈沖響應(yīng)的能量衰減曲線Fig.2 Energy decay curves of room impulse responses obtained by four different methods

除了從生成的RIR 角度對(duì)比3 種方法，本文還采用人工混響處理后的語(yǔ)音信號(hào)來(lái)評(píng)價(jià)3 種方法。首先，從TIMIT 數(shù)據(jù)庫(kù)中[1]隨機(jī)選取了20 條語(yǔ)音作為純凈無(wú)混響的信號(hào)(式(1)中的x(t))，再將純凈無(wú)混響信號(hào)分別與真實(shí)房間、3 種方法生成的各20條RIR 做卷積得到混響信號(hào)(每種方法400 條)。采用語(yǔ)音質(zhì)量感知評(píng)價(jià)[14](Perceptual evaluation of speech quality，PESQ)對(duì)混響效果進(jìn)行評(píng)價(jià)，將真實(shí)房間的混響信號(hào)作為PESQ 的參考基準(zhǔn)信號(hào)，其他3 種方法的混響信號(hào)分別與參考基準(zhǔn)信號(hào)計(jì)算得到PESQ 評(píng)分。圖3為3 種處理方法得到的混響信號(hào)的PESQ 平均得分及標(biāo)準(zhǔn)差，其中ISM 平均得分為2.77，標(biāo)準(zhǔn)差為0.04；FDN 平均得分為2.55，標(biāo)準(zhǔn)差為0.02；ARMA 平均得分為2.68，標(biāo)準(zhǔn)差為0.02。圖3表明ARMA 法生成的混響信號(hào)比FDN 法的更接近真實(shí)房間的混響信號(hào)。

聽(tīng)者的主觀感受是判斷混響感的重要評(píng)價(jià)標(biāo)準(zhǔn)，因此設(shè)計(jì)一組聽(tīng)音實(shí)驗(yàn)來(lái)判斷3 種方法中最接近真實(shí)房間產(chǎn)生的混響信號(hào)。在AISHELLASR0009-OS1語(yǔ)音數(shù)據(jù)庫(kù)中挑選包含男聲、女聲的25條干凈中文語(yǔ)音作為測(cè)試語(yǔ)音，分別設(shè)定混響時(shí)間為0.38 s、0.5 s、0.8 s 和3.8 s 四種情形，采用3 種方法生成RIR再與干凈語(yǔ)音卷積得到混響信號(hào)。實(shí)驗(yàn)中選擇10 名聽(tīng)眾，均為聽(tīng)力正常的在校研究生，對(duì)混響信號(hào)進(jìn)行試聽(tīng)后，選出3 種方法中最佳混響效果的信號(hào)。10名聽(tīng)眾選出來(lái)的語(yǔ)音中，F(xiàn)DN方法處理的占8%；ARMA 方法處理的語(yǔ)音占39%；使用ISM方法處理的語(yǔ)音占53%。實(shí)驗(yàn)結(jié)果表明ARMA方法處理的混響效果明顯優(yōu)于FDN 方法，與ISM方法接近。

圖3 3 種方法處理的混響信號(hào)的PESQ 均值及標(biāo)準(zhǔn)差Fig.3 Mean and standard deviation of PESQ obtained by three different methods

此外，在使用同一臺(tái)計(jì)算機(jī)進(jìn)行仿真計(jì)算時(shí)，同樣條件下統(tǒng)計(jì)了3種方法的運(yùn)行時(shí)間，其中，F(xiàn)DN運(yùn)算時(shí)間為1.06 s，ISM 運(yùn)算時(shí)間為3 s，而ARMA的運(yùn)算時(shí)間為0.8 s，表明ARMA 方法在實(shí)時(shí)應(yīng)用方面比ISM方法具有明顯優(yōu)勢(shì)。

3 結(jié)論

從易于實(shí)時(shí)應(yīng)用的角度，本文研究了一種卷積法人工混響方法，其特點(diǎn)是在頻域合成房間頻率響應(yīng)后再經(jīng)傅里葉反變換得到時(shí)域房間脈沖響應(yīng)，最后將干凈無(wú)混響信號(hào)與時(shí)域房間脈沖響應(yīng)卷積得到混響信號(hào)。實(shí)驗(yàn)結(jié)果表明該方法既具有比反饋延遲網(wǎng)絡(luò)方法更好的混響效果，同時(shí)又比鏡像源法有更少的計(jì)算量。