可识别声源深度的三维声聚焦波束形成方法*

丁 浩，李春晓，金江明，梅东挺，卢奂采，柴国钟

(浙江工业大学机械工程学院，特种装备制造与先进加工技术教育部/浙江省重点实验室，杭州310014)

当前近场声全息(Near-filedAcoustic Holography，NAH)[1]和波束形成[2]是对声源在三维空间中进行识别定位的主要声全息方法。NAH主要针对近场、中低频的声源识别定位，而波束形成方法有效地弥补了NAH在高频和中、长距离的不足。目前波束形成方法结合传声器阵列信号处理技术已经广泛应用于汽车、航空、雷达、通信、电子对抗和声纳等领域[4－6]。

如今人们对声源的关注已经不仅仅局限于声源究竟位于指定平面上的哪个位置，而更加关心复杂对象及复杂结构表面的声场的三维空间分布状况，而以往的研究主要着重于平面阵列结合波束形成方法识别定位声源X、Y的坐标，但这类方法给出的声场全息图中不能揭示声源的深度，所以无法识别定位声源的三维坐标，其中包括相对于平面阵列的距离(深度)[2－3，7－12]，为解决这一问题，本文提出了可识别声源深度的三维声聚焦波束形成方法，该方法给出一系列声场全息图，通过跟踪波束形成功率的最大点位置的轨迹确定声源的深度，最终识别出声源的三维坐标。

本文第1节给出了基于球面声波传播模型和波束形成方法，可识别声源深度的三维声聚焦波束形成方法的原理;第2节针对点声源形成的声场进行了仿真验证，传声器布置依照B＆K60通道的轮辐阵列上传声器位置坐标[2];第3节在全消声室内对仿真结果进行了实验验证;最后对研究结果进行了讨论并给出本文的结论。

1 三维声聚焦波束形成的原理和方法

通过由近及远对垂直于z轴的不同距离的平面进行声聚焦来完成对整个三维空间进行声聚焦的原理，见图1。移动间距定义为ΔZ，带有网格虚拟平面称为聚焦面，网格交点为聚焦点，同一个聚焦面上相邻聚焦点沿x、y轴方向的间距ΔX和间距ΔY，声聚焦即计算聚焦面上各个聚焦点上的归一化后的波束形成声功率相对输出|B(r，ω)|2，称为声场响应。

图1 可识别声源深度的三维声聚焦波束形成方法的基本原理

图1为基于球面波假设的波束形成基本原理示意图，对于任意传声器阵列形状，当阵列中阵元与声源的距离符合近场条件[7]时，声源辐射的声波假定为球面波，设定一点声源坐标为(x0，y0，z0)，第m号传声器坐标为(xm，ym，zm)，m=1，2，…，M，其中 M为传声器个数。如设定1号传声器为参考传声器，其坐标为(x1，y1，z1)，那么参考传声器接收声压信号的表达式如式(1)所示[2]:

式中，P1(ω)为参考传声器接收到的声压信号;为声源到参考传声器之间的距离;||表示向量取模;k=ω/c为波数;ω表示声源角频率;c是声波的传播速度;P0为声源的强度。

第m号传声器相对于参考传声器的延时时间Δt'm(r)如式(2)所示:

由于声源辐射的声波假定为球面波，传声器远离声源时，接收到的声压信号会有幅值衰减，所以第m号传声器相对于参考传声器接收声压信号的表达式如式(3)所示[2]:

当将传声器测量到的声压逆向映射到第k个聚焦面上的聚焦点(xf，yf，zk)时，k=1，2，…，K，其中 K为聚焦面个数，设为参考传声器到聚焦点的距离，为第m号传声器到聚焦点的距离，Δrm(r)=|rc|－|rfm|，那么第m号传声器相对于参考传声器的延时时间Δtm(r)如式(4)所示:

根据延迟求和原理，对各传声器信号进行相位延迟量Δtm(r)的补偿后求和。按传声器通道数目M，再对声压进行归一化的结果如式(5)所示[2]:

式中，ωm是第m号传声器的加权系数;B(r，ω)为波束形成归一化声压复数相对输出，相对于聚焦点的实际输出。

对式(5)取功率后根据复数的三角不等式[13]得

式中，|B(r，ω)|2为波束形成归一化声压功率相对输出。

由式(6)可知，当且仅当式(7)成立时，|B(r，ω)|2取得最大值，文中也称为声场响应的最大值。|B(r，ω)|2的最大值所在的空间位置也称为声场响应的最大值对应的空间位置。

根据各个聚焦面上声场响应的最大值沿z轴方向的轨迹变化来搜索在该方向上具有最大声场响应的聚焦面，从而判断声源在z轴方向的位置，再进一步根据该聚焦面上的声场响应分布，搜索该面上的最大声场响应点，来确定声源在x轴和y轴方向的位置。

2 仿真及验证

2．1 仿真流程和条件

本文用Matlab软件仿真，设声聚焦空间为长方体，其中心位于坐标原点，按照先验值假定其边界条件为－1 m≤xf≤1 m、－1 m≤yf≤1 m、0 m≤zk≤4 m，设 ΔX、ΔY均为0．01 m，ΔZ 为0．1 m，Z 方向采用步进扫描，即Zk+1=Zk+ΔZ。假设目标声源为点声源，考虑到测量误差，假定传声器接收到的声压信号包含10%的测量噪声，单声源情况下，声源x坐标x0为0.04 m，声源y坐标y0为－0．08 m，声源深度 z0分别为0．5 m、1．5 m、2．5 m 三种情况，声源频率为1 kHz、2．5 kHz、4 kHz;多声源情况下，声源1 的 x 坐标 xs1为－0．3 m，y 坐标ys1为－0．3 m，声源深度zs1为0．5 m，声源2 的 x坐标 xs2为－0．3 m，y坐标 ys2为－0．3 m，声源深度zs2分别为1．5 m，声源频率为4 kHz，声源强度都为50 Pa。二维阵列选用B＆K的60通道的轮辐阵列，阵列面位于x、y坐标平面内，阵列面的中心位于坐标原点。仿真验证可识别声源深度的三维声聚焦波束形成方法的有效性，仿真流程如图2所示。

2．2 仿真结果及分析

根据前面设置的仿真参数，当声源位置为(0．04，－0．08，1．5)m 时，通过 Matlab 软件仿真可以得到所有聚焦点处|B(r，ω)|2的值在三维空间的分布，但是为了更清楚的显示识别定位结果，首先找到|B(r，ω)|2的最大值所在的空间位置，再选取与此空间位置共面同时又分别垂直于x、y、z轴的所有散点，最后经过筛选后的散点图如图3所示，同样的方式可以得到多声源的结果，如图4所示。

图2 可识别声源深度的三维声聚焦波束形成方法识别定位声源三维空间位置的仿真流程

图3 波束形成识别定位单个三维声源仿真结果

图4 波束形成识别定位多个三维声源仿真结果

由图3、图4可见，仿真结果表明通过输出|B(r，ω)|2的最大值所在的空间位置可以同时识别定位单声源或双声源的三维坐标。

为了进一步验证声源识别定位的有效性，分2步来显示声源识别定位的精度。

(1)侧重于分析比较不同声源深度条件下，该方法在z方向上的识别定位精度。

聚焦面在某一深度时，输出此聚焦面上|B(r，ω)|2的最大值，随着聚焦面平移，就得到|B(r，ω)|2的最大值随zf的变化关系，最后分别模拟计算完3种声源深度情况后，其结果如图5所示。

图5 不同声源深度条件下，聚焦面上|B(r，ω)|2的最大值随zk的变化关系

(2)侧重于分析比较聚焦面分别在不同声源深度时，该方法在x、y方向上的识别定位精度。

当声源频率为4 kHz，聚焦面深度zk与声源深度z0一致时，此聚焦面上的声场分布情况如图6所示。

图6 不同声源深度条件下，zk=z0时聚焦面上|B(r，ω)|2的全息图

3 实验验证

3．1 实验步骤和条件

实验在浙江工业大学全消声室内进行，消声室的尺寸为3 m×3 m×3 m，本底噪声为18 dB，最低截止频率为63 Hz，利用轮辐阵列在全消声室内测得的声压数据来验证该方法的有效性和精度，实验的阵列和声源放置的相对位置关系如图7所示。

图7 实验布置

小型扬声器为目标声源，选用B＆K的60通道轮辐阵列作为平面测量阵列，源和阵列的放置、聚焦的声场各项参数及空间坐标系的设定与仿真时一致。实验过程中，单声源情况下，声源x坐标x0为0．04 m，声源y坐标y0为－0．08 m，为了在实验验证过程中，能够通过计算声源识别定位误差来量化衡量声源深度及频率对本方法识别定位精度的影响，设置声源深度z0范围为0．25 m～2．5 m，深度变化间隔为0．25 m，声源发出矩形脉冲信号，信号强度为30 dB，持续时间为10 s，频率分别为500 Hz～4 kHz，频率变化间隔为500 Hz，进行多次试验测量;多声源情况下，声源1 的 x坐标 xs1为－0．3 m，y坐标 ys1为－0．3 m，声源深度 zs1为0．5 m，声源 2 的 x坐标 xs2为－0．3 m，y坐标 ys2为－0．3 m，声源深度 zs2分别为 1．5 m，声源频率为4 kHz。轮辐阵列采用的是B＆K 4958传声器，有效工作频率范围10 Hz～20 000 Hz，动态范围为28 d～140 dB，轮辐阵列的传声器布置为不规则布置，因此可有效避免栅瓣的产生。

3．2 实验结果及分析

根据阵列测得的声压数据结合仿真时给出的参数，通过Matlab软件得到|B(r，ω)|2的最大值所在的空间位置，散点图如图8所示，同样的方式可以得到多声源的结果，如图9所示。

图8 波束形成识别定位单个三维声源实验结果

图9 波束形成识别定位多个三维声源实验结果

图10 不同声源深度条件下，各聚焦面上|B(r，ω)|2的最大值随zk的变化关系

由于声源深度是变化的，考虑到单声源实验结果较多，将全部单声源实验结果给出较为累赘，而0．5 m、1．5 m、2．5 m 这3 个声源深度的实验结果能够充分体现声源深度对于定位精度的影响，并且便于和仿真结果比对，因此只给出这3个声源深度条件下，聚焦面上|B(r，ω)|2的最大值随zk的变化关系，和 zk=z0时，聚焦面上|B(r，ω)|2的值随 xf、yf的变化关系，如图8、图9所示。

对比图3与图8，图4与图9，仿真和实验结果吻合，表明通过输出|B(r，ω)|2的最大值所在的空间位置可以同时识别定位单声源和多声源的三维坐标。

对比图5和图10，仿真和实验结果一致表明所有聚焦面中|B(r，ω)|2的最大值相对应的zk与声源深度一致，即可以准确识别定位声源的z坐标;声源频率和声源深度变化时，仍可以准确识别定位声源的z坐标;声源在z方向上的识别定位精度与声源距平面阵列的距离有关，距离越远，z方向空间分辨率越低，同时，声源在z方向上的识别定位精度与三维空间在z轴方向上划分的间距有关，间距越小，声源识别定位的精度越高，同时也会增加计算时间。

对比图6和图11，结果表明聚焦面在声源深度位置时，波束形成输出的最大值所在的x、y平面坐标与声源的x、y坐标一致，即该方法在智能确定深度的前提下，能同时准确的识别定位声源的x、y坐标;然而随着声源深度的增大，该方法的空间分辨率会逐渐变差。

图11 不同声源深度条件下，zf=z0时聚焦面上|B(r，ω)|2的全息图

为量化衡量本方法的识别定位精度，引入了声源识别定位误差的定义，声源识别定位误差δ为:

式中，ξ0为实际声源X坐标或Y坐标或Z坐标，ξ为声源X坐标或Y坐标或Z坐标的预测值。

如图12所示，单声源情况下在不同的深度条件时，X、Y、Z方向上的声源识别定位误差都保持在15%以内，特别是深度在2 m以内时，识别定位误差都保持在10%以内，但随着声源深度的增大，识别定位误差趋于上升，该方法的空间分辨率逐渐降低。

图12 不同声源深度条件下，声源频率4 000 Hz时，X、Y、Z方向上的声源识别定位误差曲线图

如图13所示，不同声源频率条件下，X、Y、Z方向上的声源识别定位误差都保持在20%内，特别是频率在3 000 Hz以上时，识别定位误差都保持在10%以内，随着声源频率的增高，识别定位误差趋于下降，该方法的空间分辨率逐渐提高。

图13 不同声源频率条件下，声源深度1．5 m时，X、Y、Z方向上的声源识别定位误差曲线图

4 结论

基于球面波波束形成理论，提出了可识别声源深度的三维声聚焦波束形成方法，并借助点声源模型模拟计算了声源的声场分布图。仿真和实验结果吻合，从而论证了该方法的有效性和精度。取得的主要结论如下:

(1)可识别声源深度的三维声聚焦波束形成方法能够识别定位声源的三维空间位置。

(2)随着声源深度的增大，识别定位误差趋于上升，该方法的空间分辨率会逐渐降低，在2 m以内的范围，有较好的空间分辨率。

(3)随着声源频率的增高，识别定位误差趋于下降，该方法的空间分辨率会逐渐提高，尤其在3 000 Hz以上的中、高频，有较好的空间分辨率。

该方法还应在这几方面进行完善:在单声源识别定位时，空间分辨率随声源深度的增加而下滑较快，寻求提高空间分辨率的方法，以及解决有关聚焦面大小形状和网格大小形状等参数选择的最优化、智能化等问题;进一步研究在多声源识别时，该方法的识别能力;探讨除测量噪声外，阵元安装误差、混淆误差等其他测量误差[14－16]对声源定位所带来的影响。

[1]Parthy A，Van Schaik C J A A．Acoustic Holography with a Concentric Rigid and Open Spherical Microphone Array[J]．2009，1(5):2173－2176．

[2]Christensen J J，Hald J．Beamforming[J]．B＆K Technical Review，2004(1):1－31．

[3]张金圈，毕传兴，陈心昭．Beamforming方法的阵列研究及其在噪声源识别中的应用[J]．噪声与振动控制，2009，29(3):54－58．

[4]Kirkebo J E，Austeng A．Improved Beamforming Using Curved Sparse 2D Arrays in Ultrasound[J]．Ultrasonics，2007，46(2):119－128．

[5]Van Veen B D，Buckley K M．Beamforming:A Versatile Approach to Spatial Filtering[J]．ASSP Magazine，IEEE，1988，5(2):4－24．

[6]Bai M R，Lee J．Industrial Noise Source Identification by Using an Acoustic Beamforming System[J]．Journal of Vibration and Acoustics，1998，120(2):426－433．

[7]褚志刚，杨洋，蒋忠翰．波束形成传声器阵列性能研究[J]．传感技术学报，2011，24(5):665－670．

[8]Kennedy R A，Abhayapala T D，Ward D B．Broadband Nearfield Beamforming Using a Radial Beampattern Transformation[J]．IEEE Transactions on Signal Processing，1998，46(8):2147－2156．

[9]Dhanantwari A C，Stergiopoulos S，Song L，et al．An Efficient 3D Beamformer Implementation for Real-Time 4D Ultrasound Systems Deploying Planar Array Probes[J]．IEEE Ultrasonics Symposium，2004，2:1421－1424

[10]Karl B Washburn，Tony Frazer，Jason Kunio．Correlating Noise Sources Identified by Beamforming with Sound Power Measurements[C]//SAE，2005－01－2510．

[11]Mukwevho Tshilidzi，Jordaan Jaco，Noel Guillaume．Advanced Beamforming Techniques for Acoustic Source Localization[C]//IEEE Africon Conference，2009，5308363．

[12]Arye Nehorai M，Kwok-Chiang Ho，B．T．G．Tan Minimum-Noise-Variance Beamformer with an Electromagnetic Vector Sensor[J]．Ieee transactions on signal processing，1999，47(3):601－618．

[13]Hardy G H，Littlewood J E，Pólya G．不等式[M]．2 版．北京:世界图书出版社，2004:3．

[14]翟春平．刘雨东．聚焦波束形成声图法误差分析[J]．声学技术，2008，27(1):18－24．

[15]余毅，谢植广，陈晟．球体表面圆环阵模态域波束形成误差分析[J]．电声技术，2010，34(10):69－72．

[16]贾瑞武，石庚辰．四元声传感器面阵快速测向算法及误差分析[J]．传感技术学报，2009，22(12):1757－1762．