基于HOA的圆柱腔低频声场重构
王岩, 陈克安, 胥健     
西北工业大学 航海学院, 陕西 西安 710072
摘要: 利用扬声器重构声场的传统方法是求解声学逆问题, 计算最小二乘意义下次级声源驱动函数, 缺点是测量时需要布置大量传感器且只能有效重构声阵列区域附近的声场。为了克服该缺陷, 提出利用高阶Ambisonics(higher order ambisonics, HOA)技术重构封闭空间低频声场。针对腔内基于声模态的传播特性, 推导了相关重构公式并分析其重构特性。结果表明, 当重构参数kr较小时, 波数域内格林函数谱分量主要集中于低阶成分且主要由低阶声模态贡献, 随着kr的增大, 其高阶成分变大且高阶声模态的贡献增大。在重构时, 重构滤波器对球阵以外区域的高阶谱分量有明显放大, 意味着对测量噪声也有明显放大, 而对球阵内部有明显抑制。在此基础上, 从声场采集、影响因素及重构效果等方面对HOA与声学逆问题方法进行对比。仿真结果表明, HOA优于声学逆问题方法, 对整个声腔都有较好的重构, 其重构精度与区域范围取决于球谐函数的截断阶数。最后, 在圆柱形舱室模型内进行了实验验证, HOA方法取得了较好的重构结果, 展示出良好的工程应用前景。
关键词: 声场重构     圆柱声腔     球谐函数分解    

封闭空间声场重构在噪声控制方案验证、产品声学设计等方面具有重要的应用, 其中圆柱结构作为飞机、潜艇和空间站等航行器舱室结构的典型代表, 研究其低频声场重构具有重要意义。例如对于飞机、潜艇等舱室结构的有源噪声控制[1], 降噪方案的有效性和降噪效果往往需要在实际状态下进行验证, 如果能够利用扬声器在实验室舱段模型内重构真实运行状态下的舱内声场, 这样就可以方便地评估各种降噪方案, 从而极大地节省研发及测试成本, 提高工作效率。目前已有的方法是在重构区域布置大量传声器采集声场信息, 在最小二乘意义下对空间格林函数求伪逆得到次级声源的驱动函数, 这称为声学逆问题方法[2]。近些年有学者重新整理并丰富了该理论, 讨论了求解唯一性及不适定性等问题[3], 提出基于波束形成的正则化方法并用于飞机舱室声场重构[4-5]。该方法的缺点是只能有效重构声阵列附近的区域, 要求布置大量的传感器并且会对原始声场产生干扰。

高阶Ambisonics[6-8]是主要应用于厅堂等空间声场的重构技术, 其基本原理是将声场分解为一系列球谐函数(基函数)及其展开系数的乘积, 根据模态匹配原理(mode matching)求解次级声源的驱动函数, 进而重构三维声场。目前相关研究都以自由空间声传播特性为理论基础展开, 对于三维自由空间格林函数, 它可以理解为点声源的空间-时间传递函数, 在重构过程中其波数域的谱分量呈现出空间低通特性, 能够抑制次级声源离散化导致的混叠效应[8]。在实际房间内重构空间声场时, 学者们提出了主动房间补偿滤波器[9]、房间传递函数参数化方法[10]以及利用高阶指向性次级声源[11]等方法减小壁面声反射的影响, 并取得了较好的重构效果。而对于小尺度空间低频声场, 目前有关基于HOA的声场重构研究还鲜有报道。对于圆柱声腔, 前期研究表明其声模态可以由球阵有效地进行球谐函数分解, 获取其展开系数[12], 并依据其稀疏特性, 利用压缩感知技术由较少的空间测量点求解高阶展开系数[13], 这为本文利用HOA重构圆柱腔内声场奠定了良好的基础。

为了克服传统声学逆问题方法的缺陷, 在整个三维空间内更好地重构低频声场的时间和空间特性, 本文以圆柱结构为研究对象, 提出利用HOA方法重构小尺度空间低频声场, 并与声学逆问题方法做对比。由于圆柱腔内低频声场是有限个声模态的叠加, 因此腔内次级声源的声传播特性将与以往自由场和房间扩散声场有很大不同, 主要体现在波数域内格林函数的谱分量受声模态的空间分布及其阶数影响。本文以声模态理论为基础重新推导重构算法, 重点讨论声腔内格林函数谱分量与次级声源数目、重构参数kr以及腔内声模态的关系, 定义声场重构滤波器并讨论其特性。在此基础上, 对比分析HOA与声学逆问题方法, 从声场采集、影响因素及重构效果等方面进行讨论, 最后通过实验验证本文提出方法的有效性及优势。

1 腔内HOA声场重构理论

首先介绍HOA声场重构的基本原理, 根据圆柱声腔内声传播特性, 推导腔内格林函数的谱分量表达式, 进而推导次级声源驱动函数的最小二乘解。

1.1 基本原理

假设刚性壁圆柱声腔在点声源激励下形成稳态声场, 腔内存在封闭曲面边界∂V, 如果边界∂V包围的封闭空间V内不存在声源, 那么V内的声场可以用单层势来描述[7]

(1)

式中, r=(r, φ, z)表示V内任意一点, r0=(r0, φ0, z0)为边界∂V上的点, D(kr0)称为次级声源的驱动函数, GN(kr|kr0)为圆柱声腔的纽曼格林函数[14]。在球坐标系下, 将重构区域V内部的声场P(kr), ∂V上的次级声源驱动函数D(kr0)以及纽曼格林函数GN(kr|kr0)进行球谐函数展开[14]

(2)
(3)
(4)

式中, 幅值Pnm(kr), Dnm(kr0)和Gnm(kr|kr0)为对应函数的球傅里叶变换, 本文称为谱分量, *表示共轭转置。将(2)~(4)式带入(1)式, 并利用球面上球谐函数的正交性 =δnnδmm, 可以得到如下关系

(5)

根据模态匹配方法(mode-matching method)可以推导出驱动函数的谱分量为

(6)

将(6)式带入(3)式, 即可求出次级声源的驱动函数D(kr0), 进而利用(1)式重构三维声场。

1.2 Pnm(kr)的计算

重构区域的声压谱分量可以由球阵外推得到, 假设重构区域的中心处存在一个半径为a的球阵, 对球阵表面测量声压进行球傅里叶变换, 可以得到其谱分量Pnm(ka), 利用下式可以外推整个重构区域的声压谱分量

(7)

式中, Hn(kr, ka)为传递因子[14]

(8)
1.3 Gnm(kr|kr0)的计算

圆柱腔纽曼格林函数GN(kr|kr0)的球傅里叶变换可以写为

(9)

它表示在波数域次级源与重构点之间谱分量的传递关系。根据圆柱腔GN(kr|kr0)的表达式[14], 可以计算出格林函数的谱分量为

(10)

式中, Ψnmn为第n′阶声模态函数的谱分量, n′为声模态阶数, kn为声模态波数, Vn表示声模态体积。

1.4 Dnm (kr0)的最小二乘解

在实际求解驱动函数谱分量时, 需要截取有限阶数N, 并且重构区域要离散化。由(6)式和(7)式建立线性方程组G·D=H·P, 其中GR×(N+1)2阶格林函数谱分量矩阵, R为重构半径的离散数目, D为(N+1)2×1阶驱动函数谱分量矩阵, HR×(N+1)2Hn(kr, ka)传递因子矩阵, P为(N+1)2×1阶球表面声压谱分量Pnm (ka)矩阵。由于R>(N+1)2, 方程(14)是超定方程组, Dnm (kr0)的最小二乘解可以写为[2]

(11)
2 腔内HOA声场重构实现

根据圆柱腔内低频声场HOA重构理论, 讨论其实现过程。假设有一刚性壁圆柱壳体, 长3 m, 半径0.8 m, 厚3 mm, 圆柱左右两端面分为位于z=0 m和z=3 m处, 声腔中轴线为z轴。在圆柱内部(r=0.3, φ=π, z=0.2)处有一个点声源在腔内形成稳态声场。为了便于讨论腔内声场重构特性, 本文根据球面Fliege离散点分布形式布放次级声源。重构区域是中心为(x=0, y=0, z=1.5)m, 半径r=0.6 m的球面包围的内部区域。

2.1 声压谱分量Pnm (kr)

传递因子Hn(kr, ka)在外推声压谱分量Pnm (kr)时起决定作用。由(7)式可知, 无论是空心球阵或是刚性球阵, 其传递因子的特性相似, 即Hn(kr, ka)∝(r/a)n, 这表明重构区域声压谱分量的传递特性与其阶数和半径比r/a有关, 且随着重构半径的增大传递特性呈指数放大。当r < a时, 传递因子Hn(kr, ka)相当于衰减滤波器, 谱分量的阶数越高, 越靠近重构区域的中心, 其衰减作用越大。当r>a时, 传递因子Hn(kr, ka)相当于放大滤波器, 谱分量的阶数越高, 越靠近次级源面, 其放大作用越大。

在实际外推重构区域的声压谱分量时, 选取合适的谱分量阶数非常重要。球阵表面的高阶声压谱分量幅值非常小, 对系统测量误差非常灵敏, 不易准确获取。并且, 在r>a的重构区域, 传递因子对高阶谱分量有非常明显的放大, 这意味着对噪声分量也有很明显的放大, 且阶数越高, 放大倍数越大。图 1给出了根据球阵表面声压谱分量Pnm(ka)外推的重构半径r=0.4 m的Pnm(kr)谱分量, 并与其理论值做对比, 重构频率为300 Hz, 谱分量阶数为n=0~7, 由图可知, 在n=0~4阶外推的谱分量与理论值有较好的吻合, 高于4阶以上外推的谱分量被严重放大。

图 1 重构半径r=0.4 m的球面声压谱分量
2.2 格林函数谱分量Gnm (kr|kr0)

格林函数谱分量的计算与次级声源的数目、重构频率与重构半径(即kr)以及声腔内声模态有关。在实际声场重构时, 有限的次级声源数目会影响高阶谱分量的计算。选取球面分布49个Fliege离散点, 计算n=0~7阶谱分量, 如图 2所示。

图 2 圆柱声腔格林函数谱分量的幅值

当重构半径或重构频率(即kr)增大时, 格林函数高阶谱分量的幅值会明显增大, 这与声压谱分量Pnm(kr)类似。对比图 2a)2b), 重构频率为300 Hz时n=1阶谱分量的幅值最大, 500 Hz时n=4和n=6阶谱分量的幅值最大, 这与腔体内声模态函数有关, 声场重构频率越高, 起主要贡献的声模态的阶数越高, 而高阶声模态会引起高阶谱分量幅值的增大, 所以重构频率越高, 格林函数高阶谱分量的幅值越大。当重构半径r较小时, 只存在低阶谱分量, 随着重构半径r的增大, 高阶谱分量的幅值逐渐增大, 这可以从腔体声模态函数中的贝塞尔函数看出, 重构半径r越大, 高阶贝塞尔函数的幅值越大, 这导致高阶谱分量的幅值也就越大。由图 2a)可知, 如果重构频率为300 Hz, 那么N=4阶以上谱分量接近于零, 为了避免噪声对高阶谱分量的干扰和矩阵求逆的数值不稳定, 谱分量的截断阶数取N=4。

由于圆柱声腔的格林函数是很多个声模态叠加的结果, 在求解其谱分量时需要选取合适阶数的声模态。当kr较小时, 只有低阶的贝塞尔函数具有较大的幅值, 高阶贝塞尔函数的幅值接近于零, 这说明在重构区域的中心, 只有低阶声模态起主要贡献。当重构半径r增大时, 低阶贝塞尔函数的幅值逐渐增大, 与此同时高阶贝塞尔函数的幅值也逐渐增大, 此时高阶声模态的贡献越来越大, 并逐步超过低阶声模态的贡献。当重构频率接近声模态固有频率时, 格林函数中的1/(kn2-k2)项会放大声模态函数的幅值, 所以重构频率附近的声模态将会对格林函数谱分量起重要影响, 成为最主要的贡献者。

通过以上分析可知, 在实际计算格林函数谱分量Gnm(kr|kr0)时, 次级声源半径r0要避开主要声模态的节线, 否则该阶声模态将不会参与谱分量的计算。由于在计算驱动函数谱分量时要对Gnm(kr|kr0)求逆, 所以要滤掉幅值过小的高阶谱分量, 以避免矩阵求逆时数值不稳定, 这也意味着在计算过程中要合理选择声模态阶数, 滤掉贡献很小的高阶声模态, 提高计算速度。

2.3 重构滤波器与最小二乘解

如果定义重构滤波器为Hn(kr|ka)/Gnm(kr|kr0), 那么次级声源的驱动函数谱分量Dnm(kr0)就可以看做是球阵表面声压谱分量Pnm(ka)经过重构滤波器作用的结果。重构滤波器的特性如图 3所示, 分别计算了4个不同重构半径的重构滤波器谱分量幅值, 可以看出, 前5阶(即n=0~4)滤波器的特性在不同的重构半径上基本相同, 从n=5阶开始随着重构半径的增大, 重构滤波器对高阶谱分量Pnm(ka)的放大作用明显变大, 总体而言, 在球阵以外的重构区域, 重构滤波器对Pnm(ka)的放大作用随着阶数的增加而呈指数增大。因此根据圆柱声腔的声场特性, 选择合适的截断阶数N, 将重构滤波器设计为低通滤波器, 当Pnm(ka)经过该低通滤波器后, 只保留前N阶有效信息, 滤掉易受噪声干扰且被严重放大的成分。

图 3 重构滤波器的滤波特性
3 腔内HOA与声学逆问题的对比

求解次级声源驱动函数的传统方法是声学逆问题方法, 它利用传声器阵列在重构区域采集声场信息, 根据Tikhonov正则化求解次级声源驱动函数的最小二乘解。本节将对这2种重构方法进行对比分析。

在声腔中心处半径为0.6 m的球面上按照Fliege离散点分布16个次级声源进行HOA声场重构, 声场最高截断阶数N=4, 并与声学逆问题方法做对比。对于HOA方法, 利用声腔中心处半径为0.1 m的球阵采集声场信息, 球阵表面分布49个传声器, 计算次级声源的驱动函数。对于声学逆问题方法, 在重构区域内X=0平面上均匀分布108个传声器采集声场信息, 根据Tikhonov正则化求解次级声源驱动函数的最小二乘解。

根据这2种方法计算的驱动函数重构声场分布。图 4给出了2种方法的重构效果对比, 选择重构区域内X=0的横截面与Z=L/2的纵截面为重构面, 即图中虚线圆圈内。可以看出, X=0的横截面上, 2种方法都能较好地重构出真实声场, 重构的相对误差如图 5a)所示, 只有在声模态节线位置重构误差较大, 在整个重构面上HOA方法的重构相对误差比声学逆问题的小。对于Z=L/2的纵截面上, HOA方法依然能够较好地重构声场, 只有在声模态节线位置重构误差较大, 而声学逆问题方法的重构误差就更大, 如图 5b)所示。由此可见, HOA方法在圆柱腔低频声场的重构精度更高。

图 4 重构声场与理论声场的声压幅值比较
图 5 2种方法声压幅值重构的相对误差
4 实验验证

在半消声室内, 利用HOA方法开展圆柱形飞机舱室低频声场重构的实验研究。舱室模型长度为6.4 m, 直径为2.67 m, 底板距离顶部最高点2 m, 将一个十二面体声源放置于角落, 由白噪声激励产生低频声场。将半径为0.1 m, 球面均匀分布64个传声器的刚性球阵放置于声腔中心处进行声场测量, 并将球面声压分解为一系列球谐函数及其展开系数, 如图 6所示。以声腔中心为圆心, 半径为0.35 m的球面上分布16个次级声源, 计算次级声源的驱动函数, 并重构224 Hz时声腔中心水平面的声压分布(x=0 m, y=-0.35~0.35 m, z=-2.2~1.3 m), 如图 7所示。可以看出, HOA方法的重构结果与实际测量值有较好的吻合, 而声学逆问题的重构声压幅值明显比实际测量值低, 这表明HOA方法在圆柱低频声场重构中能够取得较好的效果, 优于声学逆问题方法。

图 6 实验模型及刚性球阵
图 7 声腔水平面声压重构值与理论值的对比
5 结论

本文提出利用HOA重构圆柱腔低频声场, 并与传统的声学逆问题方法做对比。主要结论如下:

1) 利用球阵进行声场采集并将腔内声场分解为球谐函数的形式, 得到的球面声压谱分量包含了腔内声场信息, 因此可以利用该谱分量重构内部声场。这种声场采集方式避免了传统声阵列需要大量的传声器, 减少了对原始声场的干扰。

2) 声腔内格林函数的谱分量受声模态分布影响, 重构区域中心主要存在低阶谱分量并主要由低阶声模态贡献, 而远离重构区域中心存在高阶谱分量并主要由高阶声模态贡献, 在重构频率附近的声模态是谱分量最主要贡献者, 因此在实际声场重构时, 要根据重构参数kr合理选择声模态阶数, 布置次级声源时要避开主要声模态的节线处。

3) 影响HOA重构精度的主要参数是截断阶数。在声场采集时, 截断阶数越高, 获取的声场信息越丰富, 但在声场重构时, 重构滤波器会放大球面声压谱分量, 且阶数越高, 重构区域越大, 放大作用越明显, 这意味着由各种测量噪声引起的谱分量误差也会被严重放大。因此需要合理选择最优的截断阶数, 即保证重构精度又满足求解稳定性。

4) 对点声源激励下的圆柱腔内声场进行仿真分析, 结果表明HOA能够较好地重构整个三维声场, 从声场采集、声场重构的精度及区域大小等方面都优于声学逆方法。

5) 在圆柱形舱室模型内对本文方法进行了实验验证, 并与声学逆问题做对比。实验结果表明, HOA方法能够获得比逆问题更好的重构效果, 具有良好的实际应用价值。

6) 最后需要指出的是, 本文仅针对单个点声源激励下的圆柱腔内声场进行了仿真和实验研究, 且并未考虑模态重叠等情况。对于实际中多声源激励, 或结构辐射声等情况, 还需进一步研究。

参考文献
[1] 陈克安. 有源噪声控制[M]. 2版. 北京: 国防工业出版社, 2014: 296-308.
Chen Kean. Active Noise Control[M]. Second Edition. Beijing: Defense Industry Press, 2014: 296-308. (in Chinese)
[2] Nelson P A, Yoon S H. Estimation of Acoustics Source Strength by Inverse Method:Part Ⅰ, Conditioning of the Inverse Problem[J]. Journal of Sound and Vibration, 2000, 233(4): 643-668.
[3] Filippo M F. Sound Field Reproduction[D]. Dissertation England, University of Southampton, 2010
[4] Gauthier P A, Camier C, Pasco A, et al. Beamforming Regularization Matrix and Inverse Problems Applied to Sound Field Measurement and Extrapolation Using Microphone Array[J]. J Sound Vib, 2011, 330(24): 5852-5877. DOI:10.1016/j.jsv.2011.07.022
[5] Gauthier P A, Camier C, Gauthier O. Aircraft Sound Environment Reproduction: Sound Field Reproduction inside a Cabin Mock-UP Using Microphone and Actuator Arrays[C]//Proc of Meetings on Acoustics, Montreal Canada, 2013
[6] Daniel J, Nicol R, Moreau S. Further Investigations of High Order Ambisonics and Wavefield Synthesis for Holophonic Sound Imaging[C]//AES 114th Convention, Amsterdam Netherlands, 2003
[7] Ahrens J. Analytic Methods of Sound Field Synthesis[M]. Berlin: Springer, 2012.
[8] Zhang W, Abhayapala T D. 2. 5D Sound Field Reproduction in Higher Order Ambisonics[C]//International Workshop on Acoustic Signal Enhancement, 2014: 342-346
[9] Spors S, Buchner H, Rabenstein R. Active Listening Room Compensation for Massive Multichannel Sound Reproduction Systems Using Wave-Domain Adaptive Filtering[J]. J Acoust Soc Am, 2007, 122(1): 354-360. DOI:10.1121/1.2737669
[10] Betlehem T, Abhayapala T D. Theory and Design of Sound Field Reproduction in Reverberant Rooms[J]. J Acoust Soc Am, 2005, 117(4): 2100-2111. DOI:10.1121/1.1863032
[11] Betlehem T, Poletti M A. Two Dimensional Sound Field Reproduction Using Higher Order Sources to Exploit Room Reflections[J]. J Acoust Soc Am, 2014, 135(4): 1820-1833. DOI:10.1121/1.4868376
[12] Wang Y, Chen K. Low Frequency Sound Spatial Encoding within an Enclosure Using Spherical Microphone Arrays[J]. J Acoust Soc Am, 2016, 140(1): 384-392. DOI:10.1121/1.4955338
[13] Wang Y, Chen K. Compressive Sensing Based Spherical Harmonics Decomposition of a Low Frequency Sound Field within a Cylindrical Cavity[J]. J Acoust Soc Am, 2017, 141(3): 1812-1823. DOI:10.1121/1.4978247
[14] Williams E G. Fourier Acoustics:Sound Radiation and Nearfield Acoustical Holography[M]. New York: Academic, 1999.
Low Frequency Sound Field Reproduction within a Cylindrical Cavity Using Higher Order Ambisonics
Wang Yan, Chen Kean, Xu Jian     
School of Marine Science and Technology, Northwestern Polytechnic University, Xi'an 710072
Abstract: Sound field reproduction of the aircraft and submarine within a cabin mock-up using a loudspeaker array is of great importance to the active noise control technology.The conventional method is to calculate the driving functions of the secondary sources by solving an acoustic inverse problem in a least square sense, which requires a large number of microphones and only the sound field near the microphone array can be reproduced accurately.In order to overcome these drawbacks, higher order ambisonics (HOA) method which is widely used in spatial sound field synthesis for a large room is introduced to reproduce a low frequency sound field within a cylindrical cavity.Due to the different sound propagation characteristics within the cavity compared with a free field and a diffuse field, the Green function spectrum in spherical harmonics domain which is modeled as a superposition of the acoustic modes and the reproduction formulas are deduced.Reproduction characteristics are investigated by numerical simulations.Results show that for a small, the Green function spectrum in spherical harmonics domain is mainly concentrated on low orders and contributed by the low order acoustic modes, with the increase of, high order components of the Green function arise and the contributions of high order acoustic modes increase.In the reproduction process, the high order components of the pressure spectrum over the sphere in harmonics domain will be greatly amplified by the reproduction filter.Finally, HOA method is compared with the acoustic inversion method in terms of the microphone array system, the impact factors on the reproductions and the reproduction accuracy, and validated through experiments.Results show that HOA can better reproduce the entire sound field within the cylindrical cavity and the reproduction accuracy is improved.
Keywords: higher order ambisonics     sound field reproduction     cylindrical cavity    
spherical harmonics decomposition     concatenated codes     design of experiments    
西北工业大学主办。
0

文章信息

王岩, 陈克安, 胥健
Wang Yan, Chen Kean, Xu Jian
基于HOA的圆柱腔低频声场重构
Low Frequency Sound Field Reproduction within a Cylindrical Cavity Using Higher Order Ambisonics
西北工业大学学报, 2018, 36(4): 649-655.
Journal of Northwestern Polytechnical University, 2018, 36(4): 649-655.

文章历史

收稿日期: 2017-05-06

相关文章

工作空间