声学风洞内气动噪声源识别定位方法研究

王勋年，李征初，陈正武，崔红芳，张俊龙

（1.空气动力学国家重点实验室，四川 绵阳621000；2.中国空气动力研究与发展中心，四川 绵阳621000）

0 引 言

飞行器噪声水平是其市场可接受性的关键指标。对于军用飞行器，噪声水平也是其可探测性和生存能力的指标之一，现代飞行器的研制必须考虑其噪声水平。风洞试验为飞行器气动噪声研究提供有效的研究手段。在声学风洞试验过程中，试验模型处于风洞气流当中，试验模型的典型反应就是要产生一系列复杂的噪声源。准确识别噪声产生的区域以及量化噪声的强度对研究发声区域的内部结构，指导飞行器降噪研究非常重要。我国目前正在研制大型民用客机，准确测量大型民用客机机体噪声源分布，开展相应的噪声评估和降噪研究，使我国大型民用客机噪声水平满足最新的ICAO（国际民航组织）标准，对我国大型民用客机跻身国际市场具有十分重要的意义。

早期用于噪声源定位和测量的设备是一种称之为“声聚镜”的仪器，它由椭圆弧形镜面和在镜面焦点上镶嵌一个麦克风组成，它是根据椭圆形外形对声波的聚焦原理来实现声源的定位。“声聚镜”能准确定位单一噪声源位置，但是由于“声聚镜”每次只能测量一个点，因此必须不断地移动"声聚镜"以确定模型周围声源分布，试验成本高；同时测量低频（小于2kHz）噪声源时，声聚镜变得极其巨大，非常不便于使用。

20世纪70年代，从无线电、雷达天线等信息领域发展起来的现代平面麦克风阵列测试技术开始用于声学测量。1975年，NASA 的Soderman、Noble［1］等人首次采用麦克风阵列测试技术在NASA艾姆斯40英尺×80英尺风洞中测量射流噪声。后来，国外先后发展了一些更先进的麦克风阵列测试技术［2－3］，并实际用于机体气动噪声［4－5］、旋翼 气动噪 声［6］等风洞试验研究。由于麦克风阵列具有强大的声源识别能力，麦克风阵列测试技术已经成为近年来飞行器气动噪声风洞试验研究中最为常用和可靠的测试技术。国外目前广泛采用麦克风阵列法进行飞行器气动噪声源的识别和定位，已经从原来的16通道发展到目前的上千通道。国内近年来才开展风洞内气动噪声源识别和定位技术研究，北京大学黄迅教授项目组提出的一种观测器理论方法［7－8］开始用于实际风洞测试中。

基于麦克风阵列的静态噪声源识别定位算法是基于延迟求和的“波束成形”（Beamforming）算法［2］，但是在声学风洞中使用麦克风阵列时，由于风洞流场、开口风洞射流剪切层、背景噪声等会影响噪声源测量，从而影响“波束成形”算法识别声源精准度和效果。本文将根据声学风洞的特点，介绍一种适合于声学风洞的噪声源识别定位算法。

1 阵列响应函数

阵列响应函数是评估阵列性能指标的重要参数，对麦克风阵列通常采用阵列响应函数评估阵列的分辨率和动态范围。阵列响应函数的定义如下。

设xm表示麦克风阵列中麦克风的位置坐标，那么麦克风阵列中心位置坐标定义为：

如果将麦克风阵列聚焦到空间任一点x，则阵列的响应可写为：

式中x0为声源位置，r0和r0m分别为声源至麦克风阵列中心和第m路麦克风的距离，r和rm分别为扫描点至麦克风阵列中心和第m路麦克风的距离。wm是第m路麦克风的相应加权因子。相对于源点x0，麦克风阵列对任一点聚焦的响应为：

麦克风阵列响应函数可用于检查麦克风阵列的瓣宽和瓣结构，检查麦克风阵列的空间滤波能力。

2 麦克风阵列的数据处理

2.1 互谱矩阵

建立麦克风阵列互谱矩阵之前首先将麦克风阵列的时域信号分割成不重叠的长度为TB的数据块，对每个麦克风的数据块采用快速傅里叶变换（Fast Fourier Transform，缩写FFT）进行频谱分析。那么麦克风阵列的互谱矩阵（Cross Spectral Matrix，缩写CSM）可表示为：

式中

其中K为麦克风数据块数，M表示阵列的麦克风数目，Pm′k（f）表示第m′个麦克风第k段数据块的频域信号，Ws为频谱分析所选取的数据窗函数，本文的窗函数选Hamming窗，上标＊号表示共轭。互谱矩阵下三角元素通过上三角对应矩阵元素复共轭得到。

2.2 风洞流场对声波传播的干扰修正

当模型和麦克风阵列都位于风洞流场内部，声音的传播需考虑风洞流场的影响。声波在流场中传播

满足对流波动方程：

式中p为声压；v0为风洞风速，方向定义为x轴正向。由于

那么延迟时间函数τ（x）满足以下方程：

方程（8）的解可表示为：

式中M0表示马赫数，R表示声源到麦克风之间的声波传播路径。

2.3 开口风洞射流剪切层修正

对于在声学风洞开口试验段进行声学试验来说，试验模型在风洞射流内部，而麦克风阵列位于射流外部，数据处理需考虑风洞射流剪切层对噪声传播的影响。首先，将风洞剪切层模型化为无限薄的圆柱型剪切层［2］，并且假设射流内外空气是均匀的，假设声源位于开口射流中心线上。那么根据几何关系、对流波动运动关系、斯涅尔折射定律（Snell's Law）有：

式中Rm为声源到麦克风的距离，Rt为声源到剪切层的距离，θm为声源到麦克风位置向量的角度，θ为流场中的声波传播角度，θ′为声波的入射角，射流马赫数为M，c0和ct分别表示射流外和射流内的声速，θ0分别是剪切层折射点到麦克风处的角度。方程组（10）～（12）是非线性方程组，可采用牛顿迭代法求解。确定了方程组（10）～（12）的解后，可得到剪切层延迟时间为：

式中Rpath＝Rt／sin（θ）＋（Rmsinθm－Rt）／sinθ0。由于麦克风阵列主要是用于声源识别，声波振幅对麦克风阵列识别结果影响不大，本文将不对剪切层进行振幅修正，具体可参考文献［3］。

2.4 波束成形算法

麦克风阵列数据处理的关键算法是波束成形算法，采用波束成形（beamforming）算法可以排除数据处理中出现的奇异矩阵。假设空间存在一平面，在阵列前的这一平面上可以存在任何声源。那么对位于平面上的任一扫描点，阵列的信号输出矩阵为：

式中x表示聚焦位置到每一个麦克风的距离，Am为第m个麦克风的剪切层振幅修正值，ωΔtm，shear为第m个麦克风在频率为ω时的相位修正值。根据方程（14）和互谱矩阵，阵列对每一个扫描点的输出功率谱如下：

在风洞试验过程中，麦克风阵列的信噪比较低，可通过扣除互谱矩阵对角线元素的方法来减弱附面层诱导的麦克风自噪声的影响，改进麦克风阵列输出结果的动态范围。扣除对角线元素后的波束成形结果表示为：

3 数值模拟及结果分析

3.1 阵列模型及性能分析

本文选取的麦克风阵列模型为多悬臂螺旋阵列模型，总共由56个麦克风组成，图1给出了阵列空间位置分布图，阵列直径约为0.7m。

阵列瓣宽是阵列分辨率的度量参数。阵列的主瓣瓣宽B定义为在阵列响应图谱中与阵列最大响应点差为3dB处对应的主瓣宽度。图2给出本文所给麦克风阵列瓣宽与频率关系，其中扫描平面与阵列平面的距离R＝1m。由图可知阵列瓣宽与频率反比成线性关系，这一关系可表示为：

式中C为常数，λ表示波长，R表示扫描平面与阵列平面的距离，D表示阵列直径，f为分析频率。通过曲线拟合可得C＝261.3。

图1 麦克风阵列分布图Fig.1 Microphone array layout

图2 阵列瓣宽随频率变化曲线Fig.2 Array beamwidth as function of frequency

图3至图5分别给出频率为6kHz、20kHz、30kHz（z＝1m）时按照公式（3）计算的阵列响应，由图可知阵列的旁瓣与主瓣的差（阵列的动态范围）都超过10dB。阵列大的动态范围有利于提高声源识别精度，抑制在风洞试验过程中声音反射对麦克风阵列输出结果的影响。

图3 分析频率为6kHz时阵列响应图谱Fig.3 Array responses for 6kHz

图4 分析频率为20kHz时阵列响应图谱Fig.4 Array responses for 20kHz

图5 分析频率为30kHz时阵列响应图谱Fig.5 Array responses for 30kHz

3.2 简单声源模型数值验证

简单声源模型验证采用两个单位强度点声源，声源模型位置分别为（－0.1，0，1）m、（0.1，0，1）m，辐射声音频率分别为6kHz、10kHz、20Hz。麦克风阵列位于xy平面，坐标原点位于麦克风阵列中心，扫描平面平行于麦克风阵列平面，与麦克风阵列平面距离为1m。计算时首先将麦克风阵列“模拟采集到”的每个通道信号数据划分为包含8192个数据元素的数据块，采用Hamming窗函数对每个数据块进行FFT变换，将时域信号转化为频域信号。图6给出麦克风阵列信号数据处理结果。由图可知，阵列输出声源分布图中识别出的声源位置与声源模型位置基本是一致的，表明本文所采用的方法能识别出点声源的位置。并且声源频率越高，两个点声源位置精度也越高，这与阵列瓣宽是对应的。

图6 阵列的噪声源分布图Fig.6 Noise source distribution maps of array

3.3 对角线元素影响

图7和图8分别给出麦克风阵列数据处理过程进行互谱矩阵对角线元素优化处理的对比结果图（f＝10kHz，z＝1m），阵列输入信号增加了随机噪声信号。由图可知，扣除互谱矩阵对角线元素能抑制虚假声源，提高麦克风阵列的分辨率。

图7 包含CSM对角线元素时阵列输出结果Fig.7 Results of array with diagonal elements of CSM

图8 扣除CSM对角线元素时阵列输出结果Fig.8 Results of array without diagonal elements of CSM

3.4 风洞流场的影响

当麦克风阵列和噪声源都处于风洞中时，需要考虑风洞流场对声音传播路径的影响。计算模型选为单位强度点声源，声源位于（0，0，1）m，频率为20kHz，模拟风速70m／s。图9是没有考虑流场干扰影响的阵列信号处理结果图，由图可知阵列中识别出的声源位置已经偏离了实际声源位置，且图中虚假声源的点增多。图10给出进行流场干扰修正后阵列信号数据处理结果（f＝20kHz，z＝1m），由结果可知，麦克风阵列识别出的声源位置与声源模型的位置是一致的。结果表明风洞流场对声源到麦克风之间声音的传播路径（及延迟时间）有重要影响，麦克风阵列数据处理时需要考虑风洞流场的影响。

3.5 风洞射流剪切层影响

计算模型为单位强度的点声源模型，声源位于（0，0，1）m，频率为20kHz，风洞射流喷口为圆形，尺寸为1m，且声源位于风洞射流中心。图11至图14分别给出风速70m／s时和150m／s时风洞射流剪切层影响的对比结果，由图可知，麦克风数据处理时进行剪切层修正后阵列识别出的声源位置与声源模型位置是一致的，并且定位效果比不进行剪切层修正时效果好。但风速较低时剪切层影响较小，风速较高时剪切层影响较大，图13中识别出的声源位置明显与声源模型位置相差较大。

图9 不进行流场影响修正时阵列输出结果Fig.9 Results of array without correction of flow influence

图10 进行流场影响修正时阵列输出结果Fig.10 Results of array with correction of flow influence

图11 v＝70m／s剪切层未修正时阵列输出结果Fig.11 Results of array for v＝70m／s without shear layer correction

图12 v＝70m／s剪切层修正时阵列输出结果Fig.12 Results of array for v＝70m／s with shear layer correction

图13 v＝150m／s剪切层未修正时阵列输出结果Fig.13 Results of array for v＝150m／s without shear layer correction

图14 v＝150m／s剪切层修正时阵列输出结果Fig.14 Results of array for v＝150m／s with shear layer correction

3.6 分布声源模型验证

分布声源模型验证为一个"＜"字型声源，"＜"字型声源扩角为60°，尖点位于（－0.2，0）m，由22个点声源构成。声源位于阵列正前方2m处，频率分别为8kHz和10kHz。图15给出了阵列识别出的声源分布图，由图可知阵列识别出声源位置和外形与所给的声源模型位置和外形基本是一致的，表明本文所采用的计算方法可用于分布声源的识别。

图15 针对复杂声源时阵列的噪声源分布图Fig.15 Noise source distribution maps from array for distributed source

4 结 论

本文根据经典Beamforming算法，结合声学风洞流场和环境的特点，考虑风洞流场和剪切层对声音传播路径的影响，给出了一种能够用于声学风洞试验的基于麦克风阵列的气动噪声源识别方法。研究结果表明：

（1）扣除互谱矩阵对角线元素能提高阵列的动态范围，提升阵列的分辨率；

（2）风洞流场影响阵列定位精度，降低麦克风阵列的动态范围，风洞试验过程中要进行流场影响修正；

（3）风洞射流剪切层影响阵列定位精度，当风速较低时风洞射流剪切层对阵列识别声源的效果影响较小。

［1］SODERMAN P T，NOBLE S C.A four－element endfire microphone array for acoustic measurements in wind tunnels［R］.NASA－TM－X－62－331，1974.

［2］MUELLER T J.Aeroacoustic measurements［M］.Springer，2002.

［3］HUMPHREYS W M.Design and use of microphone direction arrays for aeroacoustic measurement［R］.AIAA 98－0471.

［4］HORNE C，HAYES J A.Effects of distribution source coherence on the response of phased acoustic arrays［R］.AIAA 2000－1935.

［5］GROSCHE R H，SCHNERDER G.Wind tunnel experiments on airframe noise sources of transport aircraft［R］.AIAA 1997－1642.

［6］MARCOLINI M A，BROOHS T F.Rotor noise measurement using a directional microphone array［R］.AIAA 1987－2746.

［7］HUANG X，IGOR V，BAI L.An observer for phased microphone array signal processing with nonlinear output［J］.AIAAJournal，2010，48（11）：2703－2705.

［8］季建朝，白龙，黄迅.基于状态观测器的波束形成算法及其航空应用［J］.航空学报，2011，32（1）：35－40.

［9］HUANG X.Real－time algorithm for acoustic imaging with a microphone array［J］.JournalofAcousticSocietyofAmericaExpressLetter，2009，25（5）：190－195.