非自由聲場三維重現的揚聲器系統優化配置

廖祥凝， 鄭四發， 彭 博， 連小珉

(清華大學汽車安全與節能國家重點實驗室， 北京 100084)

引 言

聲場重現作為研究聲品質主觀評價及噪聲主動控制的一種重要手段，重現系統包括揚聲器、麥克風、控制單元等，其中揚聲器和麥克風的配置會影響聲場重現的精度和準確度[1～3]。因此，如何優化揚聲器系統，減少揚聲器的個數以簡化系統，同時減少揚聲器激勵信號的計算量具有重要意義。

目前針對不同類型的優化問題發展出諸多優化算法，比如運用于線性系統的梯度下降法、拉格朗日乘數法等；運用于非線性系統及復雜問題的模擬退火算法、遺傳算法等。文獻[1]利用遺傳算法優化揚聲器和麥克風的位置，重建與目標聲場基本相同的反聲場可達到減少機艙內噪音聲壓50%的效果；文獻[2，3]采用施密特正交原理進行了揚聲器配置優化的仿真。但在揚聲器個數較多的情況下，施密特正交算法存在不穩定性，隨著迭代次數的增加，誤差累積導致向量不再正交。尤其是重建非自由場環境下的目標聲壓時，揚聲器激勵信號的計算量增加，影響到揚聲器系統優化的穩定性。針對此問題，本文提出了綜合采用傳遞函數法和QR算法的揚聲器優化配置方法。

1 非自由場三維聲場重現的控制方法

1.1 三維聲場重現的邊界控制方法

聲場重現技術是基于物理特征參量模擬的三維虛擬空間聽覺再現方法，人為地在聽者身處的物理空間內構造出虛擬現實的聲音，盡可能還原目標聲場的聽感信息[4]。由空間聲場邊界控制原理(Boundary surface control principle, BOSC)可知[5]，在封閉空間沒有其他聲源的情況下，封閉空間的內部聲壓特性完全取決于邊界聲壓特性。因此，若要重現某區域內的聲場，則可通過控制邊界聲壓與目標聲壓匹配來實現。圖1左所示邊界?V包圍內部空間V，形成封閉空間，n為邊界外法向量，則內部任意點聲壓表達式為[6]

(1)

式中ρ0為介質密度；f為頻率；ω為圓頻率；s表示內部點的位置向量；r表示邊界點的位置向量；p(s,f)和p(r,f)分別為內部點聲壓和邊界聲壓；vn(r,f)為邊界外法向速度；G(r|s,f)為自由空間的格林函數[6]。

圖1 目標聲場及重現聲場邊界控制示意圖

類似地，圖1(b)為在實驗室特定重現區域內聲場表示為

(2)

式中r′為?V′重建聲場區域邊界上的點；s′為重建區域內部點。

在保證重現區域與目標區域一致的前提下，若想滿足重建聲場等于目標聲場，即邊界包圍的封閉空間內的聲場p(s,f)=pre(s′,f)，通過令式(1)與式(2)相等，可得到

(3)

如果重現執行系統采用無指向麥克風配置，應用狄利克雷(Dirichlet)格林函數[7]，公式(1)，(2)可改寫為：

(4)

(5)

式中GD(r|s,f)，GD(r′|s′,f)均為Dirichlet格林函數。

結合邊界控制原理和公式(4)，(5)，此時若想滿足重現聲場與目標聲場一致，則只需滿足重現區域和目標區域的邊界聲壓相等。即

p(r,f)=pre(r′,f)

(6)

綜上所述，若滿足重現區域內沒有其他聲源存在，區域大小與目標區域一致，并且重現執行系統采用無指向麥克風配置的條件，則控制重現區域和目標區域邊界聲壓相等即可使得區域內聲壓相等，實現重現目的。

1.2 基于聲場傳遞函數求揚聲器激勵

如圖2所示，在帶反射的房間內布置L個揚聲器和M個麥克風，M個麥克風位置處的重建聲壓表達式為[8,9]

Pre(f)=H(f)·G(f)

(7)

圖2 非自由邊界的房間聲場重現示意圖

假設M個麥克風處的已知目標聲壓為Pd(f)，重建目標是使得Pre(f)=Pd(f)，可求得重建目標聲場揚聲器所需的激勵信號

(8)

(9)

λ(f)GH(f)G(f)

(10)

①計算正則化參數λ(f)，選擇L-曲線法求解[11]；

②求得λ(f)后，頻域范圍的激勵信號G(f)可由下面步驟解出：

1)當揚聲器個數大于麥克風個數，即L>M時，方程(8)為不定方程，則激勵信號G(f)解為

(11)

2)當揚聲器個數小于麥克風個數，即L≤M時，方程(8)為超定方程或正定方程，則激勵信號G(f)解為

G(f)=(HH(f)H(f)+λ(f)I)-1HH(f)Pd(f)

(12)

式中I為單位矩陣。

2 揚聲器系統配置的優化

由于傳遞函數包含了揚聲器和麥克風位置及個數的物理特征，因此揚聲器的配置影響著重現聲場的準確性。對揚聲器系統配置進行優化，有利于滿足準確性的同時簡化重現系統，更重要的是減少求取揚聲器激勵信號的計算量，提高聲場重現的時間效率。

2.1 揚聲器配置優化的步驟

傳遞函數H(f)為M×L矩陣，可看成由L個列向量組成

(13)

圖3 揚聲器配置優化示意

引入距離指標Jk-1

Jk-1=max(20lgri(k-1))

(14)

式中ri(k-1)指Hi(f)到Tk-1空間的距離。

若計算式(14)得到的距離指標Jk-1小于某數量值Jthr，如-40 dB(依照重現聲場應達到的精度等設定)，則可近似認為剩余的TL-Tk-1空間包括在Tk-1空間內，優化停止，即剩余傳遞函數對應的揚聲器可舍去。

將上述過程描述成如下數學問題：

(15)

由式(15)可知，Rmk表示Hi(f)到Tk-1的距離ri(k-1)?？傻?/p>

Hk(f)=arg(max(ri(k-1)))

(16)

式中 arg表示定義域的子集，即arg(max(g(t)))表示提取出定義域子集中使得函數g(·)最大值的自變量。

對于多頻聲壓信號，如頻率范圍為[fl,fh]，可以將頻率離散為N段，再引入兩個距離指標J(k-1)avg和J(k-1)min：

(17)

J(k-1)min=min(Ji(k-1)(fl),…，Ji(k-1)(fh))

(18)

J(k-1)avg用于距離最大條件的判斷，J(k-1)min用于優化是否終止的判斷。則公式(16)變為

Hk(f)=arg(max(J(k-1)avg))

(19)

2.2 初始揚聲器的選擇

不同的初始揚聲器位置對優化后選擇出的揚聲器有較大影響，直接影響著重建聲場的精準度，因此初始揚聲器的正確選擇十分重要。從式(1)可知，麥克風采集的聲壓是由多個揚聲器驅動信號的線性疊加，因此選擇揚聲器傳遞函數與目標聲壓向量夾角最小的揚聲器作為初始揚聲器[2]，如圖4所示。即

Hini(f)=arg(min(sinθi))

(20)

圖4 初始揚聲器的選擇示意圖

3 揚聲器優化方法的驗證

3.1 試驗環境

三維聲場重現在如圖5的房間內進行，房間尺寸為4.2 m×3.1 m×2.6 m，初始布置的揚聲器(21個)和麥克風(22個)如圖5，6所示。其中18個麥克風均勻分布在人工頭周圍，用于目標聲壓信號的采集、傳遞函數的測量以及聲場重現；剩余4個麥克風分布在人工頭耳朵附近，用于檢測重現聲場精度。

圖5 揚聲器配置示意圖

擬重現的目標聲壓為汽車勻速行駛過程副駕駛位耳朵附近的噪聲，采用同圖6相同的傳聲器陣列進行采集，其中的一個典型的麥克風時域信號如圖7所示。

圖6 麥克風配置圖

圖7 擬重現聲場中的聲壓信號

3.2 實驗結果

3.2.1 揚聲器優化配置結果

由于揚聲器型號不同，只有序列號17到21的揚聲器能發出100 Hz以下低頻信號，因此分頻優化，即序列號前16的揚聲器陣列與后5個揚聲器陣列分開優化。選取Jthr=-40 dB，運用上述的QR算法得到Jmin和Javg的變化，如圖8和9所示。從Jmin可以看出：聲場重建控制系統中，優化后的揚聲器陣列包括能發出低頻信號的5個揚聲器以及序列號前16的揚聲器陣列經過篩選后留下的6個揚聲器，總共11個揚聲器，其序列號分別為2，3，5，6，10，16，17，18，19，20，21。相比較于原揚聲器陣列數目，優化后的揚聲器數目減少了一半。

圖8 Jmin隨迭代次數的變化

圖9 Javg隨迭代次數的變化

3.2.2 揚聲器配置優化前后重現結果對比

利用優化前21個揚聲器與優化后11個揚聲器分別進行聲場重現，得到如圖10所示耳朵附近四檢測點時域聲壓級的結果，其中無優化重現聲壓信號的平均時域聲壓級誤差為0.47 dB，帶優化的平均時域聲壓級誤差為0.51 dB。圖11為某一檢測點處原聲場與無優化、帶優化重現聲壓信號頻譜的幅值對比，可看出噪聲信號主要包括40～160 Hz范圍內的低頻成分，目標聲壓信號與無優化、帶優化重現聲壓信號的幅值曲線都十分吻合。計算得到40～400 Hz范圍內，無優化重現聲壓信號的平均幅值誤差為1.84%，帶優化的平均幅值誤差為2.09%。

綜合圖10，11對比可以看出，無優化重現聲場與優化重現聲場都能較好重現目標聲場，驗證了優化算法的正確性。同時由于優化了揚聲器配置，公式(8)中傳遞函數減少導致計算揚聲器重現聲場所需激勵的時間也相應減少了。無優化時計算揚聲器激勵信號大約需17 h，而優化后計算揚聲器激勵信號只需要7.4 h，可見揚聲器配置優化算法的應用不僅能保證重現聲場的精度，而且大大節省了揚聲器激勵信號的計算時間。

圖10 不同揚聲器配置下4檢測點時域聲壓級的對比

圖11 不同揚聲器配置下某檢測點頻域幅值對比

4 結 論

論文以實現三維聲場重現為目的，基于聲場邊界控制方法，綜合采用傳遞函數法和QR算法進行揚聲器系統配置優化。其主要思想為選出傳遞函數中最大線性無關的列向量對應的揚聲器，并以實車信號作為目標聲場，進行了聲場重現試驗。通過揚聲器系統優化前后的重現聲場與目標聲場的對比，結果表明：優化后采用更少的揚聲器能達到與優化前揚聲器配置實現的重現效果，并且優化后減少了求解揚聲器信號激勵的計算時間，提高了重現效率，證明了揚聲器優化方法的有效性。

參考文獻：

[1] Diamantis Z G, Tsahalis D T, Borchers I. Optimization of an active noise control system inside an aircraft, based on the simultaneous optimal positioning of microphones and speakers, with the use of a genetic algorithm[J]. Computational Optimization and Applications, 2002, 23(1): 65—76.

[2] Asano F, Suzuki Y, Swanson D C. Optimization of control source configuration in active control systems using Gram-Schmidt orthogonalization[J]. Speech and Audio Processing, IEEE Transactions on, 1999, 7(2): 213—220.

[3] Enomoto S, Ikeda Y, Ise S, et al. Optimization of loudspeaker and microphone configuration for sound reproduction system based on boundary surface control principle[A]. Proc. ICA2010[C]. Sydney, Australia, 2010,2:942—948.

[4] Li Shenguang, Zheng Sifa, Peng Bo, et al. Experimental research on 3D sound reproduction with reflection boundary[A]. Inter-Noise 2012[C]. New York, 2012,5:3 624—3 632.

[5] Epain N, Friot E. Active control of sound inside a sphere via control of the acoustic pressure at the boundary surface[J]. Journal of Sound and Vibration, 2007, 299(3): 587—604.

[6] Fazi F M. Sound field reproduction[D]. Southampton: University of Southampton, 2010.

[7] Ise S. A principle of sound field control based on the Kirchhoff-Helmholtz integral equation and the theory of inverse systems[J]. Acta Acustica United with Acustica, 1999, 85(1): 78—87.

[8] Betlehem T, Abhayapala T D. Theory and design of sound field reproduction in reverberant rooms[J]. The Journal of the Acoustical Society of America, 2005, 117(4): 2 100—2 111.

[9] Peng Bo, Zheng Sifa, et al. 3D sound field reconstruction in reverberation with regularization[A]. Proceedings of 20th International congress on Sound and Vibration, ICSV 2013[C]. Bangkok, Thailand, 2013,1:277—284.

[10] 謝菠蓀.頭相關傳輸函數與虛擬聽覺[M].北京：國防工業出版社，2008:240—241.Xie Bosun. Head Related Transfer Function and Virtual Auditory[M]. Beijing: National Defense Industry Press, 2008:240—241.

[11] 王振杰, 歐吉坤. 用 L-曲線法確定嶺估計中的嶺參數[J]. 武漢大學學報(信息科學版), 2004, 29(3): 235—238.Wang Zhenjie, Ou Jikun. Determining the ridge parameter in a ridge estimation using L-curve method[J]. Geometrics and Information Science of Wuhan University, 2004, 29(3): 235—238.