Applied Acoustics 103 (2016) 37–46
Contents lists available at ScienceDirect
Applied Acoustics journal homepage: www.elsevier.com/locate/apacoust
Improvement of sound source localization in a finite duct using beamforming methods Di Zhong, Dong Yang, Min Zhu ⇑ Key Laboratory for Thermal Science and Power Engineering of Ministry of Education, Department of Thermal Engineering, Tsinghua University, Beijing 100084, China
a r t i c l e
i n f o
Article history: Received 14 November 2014 Received in revised form 16 September 2015 Accepted 14 October 2015
Keywords: Localization of sound source Beamforming Duct DAMAS
a b s t r a c t This paper presents work to identify the sound source location within a finite length duct with beamforming methods. Two types of problems were studied in this paper. First, the deconvolution approach for the mapping of acoustic sources (DAMAS) method is applied to minimize the side-lobe effects, and improve the spatial resolution at low frequencies. Second, the beamforming results may deteriorate when the sound frequency is close to the cut-on frequencies. In order to improve the results, an algorithm was developed. To validate the methods, both simulations and experiments were implemented. In the experiments, a synchronized system with forty-eight wall-mounted microphones in a finite length duct was established. The dependence of the azimuthal sound source localization and the planar sound source localization on the measurement frequency was studied. The results indicate that the accuracy of the localization is highly dependent on the number of the spinning modes that are able to propagate in the duct. The localization may fail when the frequency is close to the cut-on frequencies. Significant improvements can be achieved when DAMAS and the improved algorithm are implemented. Ó 2015 Elsevier Ltd. All rights reserved.
1. Introduction In civil or industrial environments, it is common to deliver the gaseous working fluid with pipe and duct network, which consists of many seal components and various bypass or throttling valves. It is vital to avoid any leakage of working fluid, disorder of valves and unexpected mechanical vibration. For example, in ventilation system of large building, normally it is desirable that air can be distributed quietly, and the noise can be isolated or suppressed efficiently. In chemical and power plants, gaseous fluids in pipes are often in high pressure and temperature, sometimes they are also toxic, flammable or explosive. In order to ensure the safety of the pipe network, it is essential to develop the on-line monitoring and detection techniques. Once unforeseen fault happens, the location of the fault can be promptly identified. In this situation, the leakage of gas or failure of valves will inevitably generate noise or mechanical vibration within the pipe. In this work, an approach to identify sound source within a finite duct was developed based on beamforming. The inverse methods and phased array beamforming are two major methods that are used to identify the sound sources in the duct via the pressure-based measurements inside the duct. The
⇑ Corresponding author. E-mail address:
[email protected] (M. Zhu). http://dx.doi.org/10.1016/j.apacoust.2015.10.007 0003-682X/Ó 2015 Elsevier Ltd. All rights reserved.
main objective for inverse methods is to reconstruct the sound source strength at prior assumed positions with the pressure measured. The relationship between the sources and the signals measured by the microphones in the array is linked by a matrix of Green functions, which need to be inverted in the inverse methods. Hewlett and Nelson [1] first presented an experimental technique for the investigation of the characteristics of an acoustic source within an infinite duct. This technique used a model-based approach to determine the best estimation of the complex source strength according to the measured pressure within the duct. Kim and Nelson [2] developed this model-based methods for a semi-infinite duct to reconstruct the azimuthal position of a simple source. They presented an analytical Green function between the sound source and the sound field within a hard-walled, semiinfinite and no-flow duct. The experiments indicated that higher spatial resolution of the acoustic source distribution could be realized when the pressures were measured near the sound source to capture the evanescent waves. Lowis and Joseph [3] extended the inversion technique to estimate the dipole-like rotating broadband sound source over a ducted rotor with inlet flow. Bravo and Maucy [4] reconstructed the acoustic source distribution from either pressure-based or velocity-based measurements within a finite duct. The results indicated that velocity-based field measurements were more favorable than pressure-based field measurements for acoustic source reconstruction, both in accuracy and in spatial
38
D. Zhong et al. / Applied Acoustics 103 (2016) 37–46
resolution. With the inverse methods, the source strengths and mutual correlations of those sources can be obtained. However, the drawbacks of the inverse methods are the prior assumptions of the source positions and the ill-condition in the inversion of the matrix with Green function. Another choice for identification of sound source in the duct is beamforming, which is an array-based measurement technique that allows the microphones gains to be delayed and summed to locate the sound source in the space. Beamforming is simple to implement, but it is limited to obtaining a qualitative evaluation of the source power. Sijtsma [5,6] discussed the possibility of source localization via phased array beamforming on the fan and stator of turbofan engines using an intake wall-mounted microphone array. The free-field Green function ignoring the duct wall reflections was used. Lowis and Joseph [7] used the phased array beamforming to quantify the aero-acoustic sources in a duct. A calculated Green function, which considered the duct wall reflections, was applied. They reported an in-duct focused beamformer technique for measuring the relative contributions to the fan broadband noise generated in an aero-engine duct due to the interaction between rotor and stator. Their results indicated that beamforming close to a cut-on frequency is impossible because the measured pressure is dominated by a single mode. To avoid the difficulty that the Green function can be dominated by the nearcut-on mode when the frequency near cut-on, Dougherty and Walker [8,9] applied steering vectors based on the unsteady pressure at the fan, rather than the Green function. The beamforming produces an output with side-lobe effects and is limited by the spatial resolution. To remove the array-dependent beamforming characteristics from the output representations, the DAMAS algorithm was introduced by Brooks and Hhumphreys [10]. DAMAS attempts to deconvolve the true signal powers from the beamfoming result by iteratively solving a linear equation that relates the beamforming output and the accurate sound source. DAMAS was subsequently further developed by many other researchers [11–14]. Apart from DAMAS, there are many other deconvolution algorithms used for beamforming imaging, such as NNLS [15], CleanSC [16]. For the induct sound source measurements, CleanSC was applied by Dougherty and Walker [8,9] and Sijtsma [17], and DAMAS was used by Dougherty and Walker [9] to the beamforming results to remove the array resolution effects. For the spatial resolution limit and side-lobe effect at low frequencies in beamforming maps, DAMAS can be applied to enhance the resolution and eliminate the side lobes. To improve the beamforming technique, particularly when the frequencies are close to cut-on frequencies, in this work, an algorithm was developed. Significant improvements of source localization can be obtained when DAMAS and the developed algorithm are used. In Section 2, the theory of the pressure in the duct and the beamforming methods are presented. In Section 3, the basic concept of DAMAS is introduced, and the algorithm dealing with the near cut-on frequencies is developed. Section 4 presents the simulated and experimental validation of DAMAS and the improved algorithm. The dependence of the azimuthal sound source localization and the planar sound source localization on the frequency is studied. The conclusion is presented in Section 5.
2. Theory 2.1. The acoustic pressure field in a finite duct Assume the two terminations of the hard-wall duct are open at the axial positions z ¼ z1 and z ¼ z2 . The duct is of a constant radius R. There is a monopole source at the position xs ¼ ðr s ; hs ; zs Þ with an amplitude Q 0 (volume velocity amplitude). In the frequency
domain, the acoustic pressure in the duct is written as the Helmholtz equation, which is given by [18]
r2 p þ k2 p ¼ ixq0 Q 0 dðxp xs Þ;
ð1Þ
herein, p is the complex acoustic pressure in the frequency domain, and k ¼ x=c is the dimensionless frequency, where x is the angular frequency of the sound and c is the speed of the sound; i is the imaginary unit, d is the Dirac delta function, and xp ¼ ðr; h; zÞ is the cylindrical vector. The solution of Eq. (1) can be expressed as a summation of a series of spinning modes, which is given by
pðr; h; zÞ ¼
1 X 1 X Amn wmn ðr; hÞ;
ð2Þ
m¼1 n¼1
where m and n denote the azimuthal and radial modes, respectively.
Amn ¼ aþmn eikmn z þ amn eikmn z ;
ð3Þ
or þ
Amn ¼ bmn eikmn z þ bmn eikmn z ;
ð4Þ
bmn
where a are the amplitudes of the spinning modes that mn and travel in the two directions z, and wmn is the cylindrical mode shape function that has the expression
wmn ðr; hÞ ¼
J m ðrmn r=RÞ imh e ; Nmn
ð5Þ
where Jm is the Bessel function of the first kind of order m. In the r ¼ R, the velocity fluctuation equals zero, so that we have J0m ðrmn Þ ¼ 0. The axial wave number kmn is given by 2
2
kmn ¼ k ðrmn =RÞ2 . Nmn is the normalization factor. Consider that the shape function wmn is orthogonal, which could be expressed as
Z
S
wmn ðr; hÞwab ðr; hÞdS ¼
1; m ¼ a; n ¼ b; 0;
ð6Þ
m – a; n – b:
Substituting Eq. (5) into Eq. (6), the normalization factor Nmn is determined by
Z
N2mn ¼
S
J 2m ðrmn r=RÞdS ¼ pR2 ð1 ðm=rmn Þ2 ÞJ 2m ðrmn Þ:
ð7Þ
At the ends of the duct, we define the modal pressure reflection coefficients as
Rþmn ¼
amn 2ikmn z2 e ; aþmn
Rmn ¼
bmn 2ikmn z1 : e bmn
ð8Þ
þ
A detailed calculation of the two above-described coefficients can
be found in Ref. [19]. The values of a mn and bmn are obtained by matching Eq. (2) to the in-duct source distribution. A general expression of the acoustic pressure radiated by a monopole source with amplitude Q 0 inside the duct is given by
pðr; h; zÞ ¼
þ1 X 1 X xq0 Q 0 m¼1 n¼1
Jm
2pR2 rmn Rr Jm rmn rR0
kmn J 2m ðrmn Þ½1 ðm=rmn Þ2
eimðhh0 Þ eikmn ðzz0 Þ
ð1 þ Rmn e2ikmn ðz1 z0 Þ Þð1 þ Rmn e2ikmn ðzz2 Þ Þ 1 R2mn e2ikmn ðz1 z2 Þ
:
ð9Þ
The above-described relationship can be expressed as pðxp Þ ¼ Gðxp jxs ÞQ 0 , where G is the Green function that relates the monopole source and the measured pressure at the sensor positions.
39
D. Zhong et al. / Applied Acoustics 103 (2016) 37–46
2
2.2. Beamforming technique Assume that N is the number of microphones used. The pressure measured from the array of microphones located at x ¼ ðx1 ; x2 ; . . . ; xN ÞT can be expressed in vector form as
b1
3
2
wH1 g1 gH1 w1 wH1 g2 gH2 w1
6 b 7 6 wH g gH w wH g gH w 2 6 27 6 2 1 1 2 2 2 2 6 . 7¼6 .. .. .. 6 . 7 6 4 . 5 4 . . . bM wHM g1 gH1 wM wHM g2 gH2 wM
3 3 2 Q 1 Q H1 wH1 gN gHN w1 7 6 wH2 gN gHN w2 7 Q 2 Q H2 7 7 6 7 7 6 7: .. 7 6 . 6 . 5 4 . 7 . 5 wHM gN gHN wM Q M Q HM
ð10Þ
ð19Þ
The relationship between the pressure measured and the sound source is given by
The matrix involved in the above linear equation is ill-conditioned so that the direct inversion of the matrix is impossible. Thus, DAMAS solves Eq. (19) using the iterative Gauss–Seidel method [10].
pðxÞ ¼ ðp1 ðxÞ; p2 ðxÞ; . . . ; pN ðxÞÞT
pðxÞ ¼ gðxjyb ÞQ 0 ðyb Þ;
ð11Þ
where gðxjyb Þ is the steering vector of the Green function,
gH ðxjyb Þ ¼ ðGðx1 jyb Þ; Gðx2 jyb Þ; . . . ; GðxN jyb ÞÞ;
ð12Þ
where the superscript H denotes complex conjugate transposition. The goal of the beamforming is to determine complex amplitude Q 0 of the source in yb . This goal is achieved by minimizing the cost function defined by [20] 2
J ¼ kp Q 0 gk :
ð13Þ
The solution of this minimization problem is
Q 0 ¼ wH p;
ð14Þ
where w is a weighting vector expressed as 1
w ¼ ðgH ðxjyb Þgðxjyb ÞÞ gðxjyb Þ:
ð15Þ
The source power spectral can be written as
bðyb Þ ¼ wH ðxjyb ÞCðxÞwðxjyb Þ;
ð16Þ
where C is a cross-spectral matrix of the measured pressure. The image of the sound source distribution can be obtained by calculating the value of b in every node of the scanning region.
3. Improvement of beamforming technique 3.1. DAMAS deconvolution In principle, the beamforming map reflects the relative amplitude of the acoustic source strength. The DAMAS deconvolution assumes that the beamforming map in the observation plane can be considered as a convolution of the actual sources and the beamformer’s point spread function (PSF), which describes the beamformer’s response to a point source. The PSF depends on parameters, such as the frequency and the position of the focused point. In the ideal case, the PSF equals to the delta function and the beamforming map would manifest the location of the sound source. However, the ideal case can never occur because of the finite size of the array and the discrete pattern of the microphones. In the low frequency condition, the convolution of the actual sources and the PSF would generate the beamforming map with a wide main lobe and some side-lobes. DAMAS offers the possibility to remove the influence of the side-lobes and the resolution effects, thereby more accurately quantifying the position. Considering that there are M beamforming points, the relationship between the pressure field and the beamforming points can be represented as H
pðxÞ ¼ ½g1 ; g2 ; . . . ; gM ½Q 1 ; Q 2 ; . . . ; Q M :
ð17Þ
The beamforming expression is given by Eq. (16),
b ¼ wH Cw ¼ wH ppH w:
ð18Þ
Suppose the sources at the beamforming points are mutually incoherent. Substituting Eq. (31) into Eq. (18), the source powers and the beamforming results are related by
3.2. Improved algorithm The cut-on frequency in a hard-wall duct is the minimal frequency at which a given acoustic mode ðm; nÞ begins to propagate. qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 Examining Eqs. (3) and (4), if the value of kmn ¼ k ðrmn =RÞ2 is real, then the mode is ‘‘cut-on” and will propagate along the duct. However, if the value of kmn is imaginary, then the mode will decay along the duct and is called ‘‘cut-off”. The cut-on frequency f mn for the mode ðm; nÞ corresponds with the k that makes kmn ¼ 0, so f mn ¼ ðrmn cÞ=ð2pRÞ. Thus, the cut-on frequencies are related to the radius of the duct. Table 1 shows the cut-on frequencies of different spinning modes in the duct with r ¼ 0:06 m. The phenomenon regarding the beamforming results near the cut-on frequencies can be explained by considering the behavior of the Green function used in the beamforming formulation of Eq. (15). The expression of the Green function can be found in Eq. (9). The denominator includes the axial wave number kmn , which is zero at the cut-on frequency. Therefore, the acoustic pressure is dominated by a single near-cut-on mode at a cut-on frequency. The location of the source cannot be determined by a single mode; as a result, the beamforming results deteriorate. A rearrangement of the Green function in Eq. (9) and the measured pressure must be performed. In this manner, an algorithm is developed to improve the results. Here, we consider discarding the largest mode on both sides of the equation to remove the influence of this single mode,
pmodify ¼ Gmodify Q 0 ;
ð20Þ
where pmodify ¼ p pmn max and Gmodify ¼ G Gmn max . The pmn max and the Gmn max coincide with the near-cut-on mode. The terms of Gmn can be calculated based on Eq. (9). To calculate the pressure component pmn max , the mode decomposition is considered. Eq. (9) can be expressed as
p¼
1 1 X X pmn ¼ Dmn umn ; m;n
where
Dmn ¼ J m
ð21Þ
m;n
rmn
r imh ikmn z e e ð1 þ Rmn e2ikmn ðzz2 Þ Þ: R
ð22Þ
Table 1 Cut-on frequencies (Hz) of different spinning modes in the duct with r ¼ 0:06 m. Mode
n=1
n=2
n=3
n=4
m=0 m=1 m=2 m=3 m=4 m=5 m=6 m=7 m=8 m=9
0 1661 2754 3789 4796 5786 6765 7736 8700 9660
3456 4808 6048 7229 8372 9488 10,583 11,663 12,730 13,787
6327 7699 8991 10,233 11,437 12,614 13,770 14,907 16,030 17,140
9175 10,557 11,878 13,155 14,398 15,614 16,808 17,985 19,146 20,394
40
D. Zhong et al. / Applied Acoustics 103 (2016) 37–46
Eq. (22) contains the information of the positions of the microphones, while
umn ¼
J m rmn rR0 xq0 Q 0 2 2 2pR kmn J m ðrmn Þ½1 ðm=rmn Þ2
ð1 þ Rmn e2ikmn ðz1 z0 Þ Þ 1 R2mn e2ikmn ðz1 z2 Þ
eimh0 eikmn z0 ;
ð23Þ
which contains the information of the sound source and remains constant for the different pressure measured by each microphone. Note that because the spinning modes that can propagate along the duct are finite, Eq. (21) is restricted to the cut-on modes and can be expressed in matrix form
P ¼ Du;
ð24Þ
where P ¼ ðp1 ; p2 ; . . . ; pN ÞT is a complex vector of the radiated acoustic pressures at the microphone locations. u is a complex vector of the L modal amplitudes, where L is the number of cut-on modes. D is a complex N L matrix which contains the information of L modes at N microphone positions. In an actual measurement, the contaminating noise e is inevitable. Thus, the actual measured acoustic pressured can be expressed as
^ ¼ Du þ e: P
ð25Þ
It is shown that the optimal estimate of the vector u is minimizing the cost function defined by [21]
^ 2; J ¼ kek2 ¼ kDu Pk
ð26Þ
where kk denotes the 2-norm of the vector. The least squares solution of u is given by
^ u ¼ Dþ P;
ð27Þ 1
where Dþ ¼ ½DH D DH is the pseudo-inverse of the matrix D. However, the condition number of the matrix D is large. The problem is ill-conditioned so that Eq. (27) cannot be applied directly and regularization methods must be considered to find an acceptable solution. In this paper, the Tikhonov regularization method [22] is used to solve the ill-conditioned problem. The method includes a regularization parameter to control the amount of regularization. The general form of the estimate of the vector u via the application of Tikhonov regularization is written as 1
^ u ¼ ½DH D þ bI DH P;
ð28Þ
where b is the regularization parameter to be determined and I is an identity matrix. To optimize the parameter, the L-curve [23] or GCV method [24] is usually used.
After the determination of the u, the terms of pmn are easy to compute based on Eq. (21). The pmn max can be found among the pmn , so that the pmodify can be determined. With the modified pressure p and Green function G, Eq. (16) can be written as
bðxjyb Þ ¼ wHm ðxjyb ÞCm ðxÞwm ðxjyb Þ;
ð29Þ
where Cm is a cross-spectral matrix of the modified pressure. The modified weighting vector wm can be expressed as 1
wm ¼ ðgHm ðxjyb Þgm ðxjyb ÞÞ gm ðxjyb Þ;
ð30Þ
where gm ðxjyb Þ is a vector of modified Green function. 4. Validation via simulations and experiments In this section, simulations and experiments are implemented to validate the methods described above. Figs. 1 and 2 illustrate the experimental set-up. The duct had a radius of 0.06 m and a length of 1.565 m. Both ends of the duct were open. The acoustic pressure was generated by the loudspeaker, which was excited by the single-frequency sinusoidal signal produced by the signal generator. A small pipe was used to connect the loudspeaker and the duct. The inner diameter of the pipe was 0.006 m so that the velocity fluctuation at the exit could be regarded as a monopole source. The exit of the pipe was located at r ¼ 0:04 m, where r is the distance from the center of the duct. The pressure data were measured by 48 GRAS microphones. The surface of the duct was milled to plane so that the microphones could be mounted flush on the duct wall parallel to the cross-section conveniently, as shown in Figs. 1 and 2. The microphone array consisted of two rings of 24 equidistant microphones. The distance between microphone rings was 0.04 m. The measured field pressure data were acquired simultaneously using a 48-channel National Instruments PXI-1082 chassis with 3 NI PXI-4496 data acquisition modules. Each channel has 24-bit resolution with a 114 dB dynamic range. All measurements were AC coupled with a 3 dB cut-on at 0.5 Hz, and the appropriate anti-aliasing filters were applied. The sampling rate was set as 100 kHz. The low pass corner frequency of the amplifier is 6 Hz and the high pass corner frequency is 156 kHz. The conditions in the simulations, such as the positions of the sound source and the microphone array, were identical with those in the experiments. To account for the background noise present in the actual measurements, simulations were performed with a signal-to-noise (SNR) of 20 dB at each pressure signal of the microphone. Using the simulated pressure data and the methods discussed above, the location of the assumed sound source can be determined.
Fig. 1. Arrangement of the experimental configurations.
41
D. Zhong et al. / Applied Acoustics 103 (2016) 37–46
360
0dB
300 −5dB
240 180
−10dB
120 60 0 2000
Azimuth [degree]
Azimuth [degree]
Fig. 2. Picture of the set-up (Left); microphone installation (Right).
−15dB 4000
6000
8000
360
180 60
−5dB
240 180
−10dB
120 60
Azimuth [degree]
Azimuth [degree]
0dB
300
−15dB 6000
8000
−5dB
240 180
−10dB
120 60
Azimuth [degree]
Azimuth [degree]
0dB
−15dB 8000
−5dB
240 180
−10dB
120 60
−15dB
Frequency [Hz]
9600
Azimuth [degree]
Azimuth [degree]
0dB
8000
−10dB
120 60
−15dB 4000
6000
8000
9600
360
0dB
300 −5dB
240 180
−10dB
120 60
−15dB 4000
6000
8000
9600
Frequency [Hz]
300
6000
−5dB
180
0 2000
9600
360
4000
0dB
240
Frequency [Hz]
0 2000
9600
Frequency [Hz]
300
6000
8000
300
0 2000
9600
360
4000
6000
360
Frequency [Hz]
0 2000
−15dB 4000
Frequency [Hz]
360
4000
−10dB
120
Frequency [Hz]
0 2000
−5dB
240
0 2000
9600
0dB
300
360
0dB
300 −5dB
240 180
−10dB
120 60 0 2000
−15dB 4000
6000
8000
9600
Frequency [Hz]
Fig. 3. Azimuthal sound source localization at various frequencies. The results obtained with beamforming (top row), the maps after applying the improved algorithm (second row), after applying DAMAS (third row) and after applying the improved algorithm + DAMAS (fourth row). The left column presents the simulated results, and the right column presents the experimental results.
42
D. Zhong et al. / Applied Acoustics 103 (2016) 37–46
Beamforming result
0.4 0.3 0.2 0.1 0 −200
−100
0
100
200
θ −θ [degree] s
b
Fig. 4. Theoretical result of beamforming in the azimuthal direction when the maximum azimuthal mode m ¼ mþ ¼ 3.
4.1. Azimuthal sound source localization In this section, we study the dependence of the azimuthal sound source localization on the measurement frequency. The sound source was located at the radial position r ¼ 0:04 m and the azimuthal position / ¼ 180 , as shown in Fig. 1. The beamforming scanned the azimuthal positions at r ¼ 0:04 m from 0 to 358 with a resolution of 2°. The frequency was set from 2000 Hz to 9600 Hz. This is the beamforming frequency range looked at.
The first cut-on frequency of the duct with radius r ¼ 0:06 m is 1661 Hz, as shown in Table 1. Below the 1661 Hz, the pressure in the duct is one-dimensional plane wave and beamforming cannot be used to locate the sound source. Therefore, the minimum frequency was chosen as 2000 Hz. As for the maximum frequency, 9600 Hz was chosen based on Table 1. The conclusions obtained here should be also suitable for the frequencies beyond the 9600 Hz. For the convenience of comparison among different frequencies, the beamforming outputs were normalized in a way that the maximum of the source corresponds to 0 dB. The results of the simulations and experiments are displayed in Fig. 3. The left side of Fig. 3 shows the simulated results obtained with beamforming along with the maps after applying the improved algorithm, DAMAS and the improved algorithm + DAMAS. The results of the experiments are shown on the right side of the figure. Apparently, the simulations and experiments yield very similar results. The beamforming maps show that the direction of the sound source is / ¼ 180 , which agrees with the simulated and experimental condition. However, at low frequencies, the number of the propagating spinning modes is few, and the main lobes are wide. With the increase of the frequency, the width of the main lobe decreases. Meanwhile, the side lobes appear in the map. In the case of relatively high frequencies, the side lobes appear at
360
DAMAS improved+DAMAS
300
Azimuth [degree]
Azimuth [degree]
360
240 180 120 60 0 2000
4000
6000
8000
DAMAS improved+DAMAS
300 240 180 120 60 0 2000
9600
Frequency [Hz]
4000
6000
8000
9600
Frequency [Hz]
Fig. 5. Comparison of the experiment results of Fig. 3(c) and (d).
0.06
0dB
0.06
−5dB
0.02 0
−10dB
−0.02
0
−0.06
0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4
0.06
−15dB 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4
z [m]
0dB
0.06
0.04
0dB
0.04 −5dB
0.02 0
−10dB
−0.02 −0.04
r [m]
r [m]
−10dB
−0.04 −15dB
z [m]
−0.06
−5dB
0.02 −0.02
−0.04 −0.06
0dB
0.04
r [m]
r [m]
0.04
−5dB
0.02 0
−10dB
−0.02 −0.04
−15dB 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4
z [m]
−0.06
−15dB 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4
z [m]
Fig. 6. Sound source maps with beamforming at low frequencies. The left column presents the simulated results, and the right column presents the experimental results. The star denotes the location of the actual sound source.
43
D. Zhong et al. / Applied Acoustics 103 (2016) 37–46
0.06
0.06
0dB −5dB
0 −10dB
−0.02
0
−0.04 −15dB
−0.06
0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4
z [m]
0.06
−15dB 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4
z [m]
0.06
0dB
0dB
0.04
0.04 −5dB
0 −10dB
−0.02
−5dB
0.02
r [m]
0.02
r [m]
−10dB
−0.02
−0.04
0 −10dB
−0.02 −0.04
−0.04 −0.06
−5dB
0.02
r [m]
r [m]
0.02
−0.06
0dB
0.04
0.04
−15dB
−0.06
0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4
−15dB 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4
z [m]
z [m]
Fig. 7. Sound source maps with beamforming at relatively high frequencies. The left column presents the simulated results, and the right column presents the experimental results. The star denotes the location of the actual sound source.
some special frequencies. Based on the data in Table 1, these frequencies are observed to be close to the cut-on frequencies of the duct. This phenomenon is explained in Section 3.2. It also can be seen from the beamforming maps that the number of the side lobes equals to the double of the maximum azimuthal mode of the frequency. Lowis [25] demonstrated that the theoretical distribution of the beamforming result in the azimuthal direction is approximated by the expression as
1 m¼m Xþ imðhs hb Þ b¼ e ; Nh m¼m
ð31Þ
0.06
where b denotes the beamforming result, N h is the number of microphones per ring, m ¼ mþ is the maximum azimuthal mode of the frequency, hs is the azimuthal position of the sound source from 0° to 358° and hb is the focused point with hb ¼ 180 . Assuming that the maximum mode m ¼ mþ ¼ 3, the theoretical output of beamforming in azimuthal direction is shown in Fig. 4. As seen from the figure, the number of the side lobes is the double of the maximum mode. The second row of the Fig. 3 shows the maps after applying the improved algorithm at some special frequencies. The results reveal that the improved algorithm is effective in the near cut-on frequency condition. The third row shows the results of
0dB
0.06
0.04 −5dB
0 −10dB
−0.02
0
−0.04 −15dB
−0.06
0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4
z [m]
0.06
−15dB 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4
z [m]
0.06
0dB
0dB
0.04
0.04 −5dB
0 −10dB
−0.02
−5dB
0.02
r [m]
0.02
r [m]
−10dB
−0.02
−0.04
0 −10dB
−0.02 −0.04
−0.04 −0.06
−5dB
0.02
r [m]
r [m]
0.02
−0.06
0dB
0.04
−15dB 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4
z [m]
−0.06
−15dB 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4
z [m]
Fig. 8. Sound source maps with beamforming at near-cut-on frequencies. The left column presents the simulated results, and the right column presents the experimental results. The star denotes the location of the actual sound source.
44
D. Zhong et al. / Applied Acoustics 103 (2016) 37–46
0.06
0dB
0.06
0.04
r [m]
0 −10dB
−0.02
r [m]
−5dB
0.02
0
−0.06
0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4
0.06
−15dB 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4
z [m] 0dB
0.06
0.04
0dB
0.04 −5dB
0.02 0
−10dB
−0.02
r [m]
r [m]
−10dB
−0.04 −15dB
z [m]
−5dB
0.02 0
−10dB
−0.02
−0.04 −0.06
−5dB
0.02
−0.02
−0.04 −0.06
0dB
0.04
−0.04 −15dB
−0.06
0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4
−15dB 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4
z [m]
z [m]
Fig. 9. Sound source maps with DAMAS at low frequencies. The left column presents the simulated results, and the right column presents the experimental results. The star denotes the location of the actual sound source.
applying DAMAS alone to the beamforming results. The simulated result is different from the experimental result. The experimental result shows that after applying DAMAS, the side lobes disappear at low frequencies, and the width of the main lobe decreases at most frequencies. However, at near cut-on frequencies, the side lobes still exist. Even at some frequencies, such as f ¼ 5800 Hz and f ¼ 7700 Hz, the sound source deviates from the / ¼ 180 . The result of the simulation is perfect, the reason is that the measurement error of the pressure microphones is not considered in the simulation. The fourth row presents the results of the improved algorithm + DAMAS: the main lobe clearly becomes more directed, and the side lobes vanish. Fig. 5 shows the comparison of the experiment results of Fig. 3(c) and (d), and quantifies the difference
0.06
between them. It is seen that applying DAMAS alone can determine the position of the main lobe correctly at the most frequencies. However, DAMAS leads to the large discrepancy at some near cuton frequencies. Fig. 5(b) shows the distribution of the side lobes that are higher than 6 dB. The improved algorithm + DAMAS result has less side lobes than the result of DAMAS. Therefore, improved algorithm + DAMAS leads to more accurate result. 4.2. Planar sound source localization In this section, we investigate the planar sound source localization, which can be regarded as the axial and radial localization. The reconstruction zone was the cross-plane of the duct, which covered
0dB
0.06
0.04 −5dB
0 −10dB
−0.02
0
−0.04 −15dB
−0.06
0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4
z [m]
0.06
−15dB
0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4
z [m]
0dB
0.06
0.04
0dB
0.04 −5dB
0 −10dB
−0.02 −0.04
−5dB
0.02
r [m]
0.02
r [m]
−10dB
−0.02
−0.04
−0.06
−5dB
0.02
r [m]
r [m]
0.02
−0.06
0dB
0.04
0 −10dB
−0.02 −0.04
−15dB
0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 z [m]
−0.06
−15dB
0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 z [m]
Fig. 10. Sound source maps with the improved algorithm at near cut-on frequencies. The left column presents the simulated results, and the right column presents the experimental results. The star denotes the location of the actual sound source.
45
D. Zhong et al. / Applied Acoustics 103 (2016) 37–46
0.06
0dB
0.06
−5dB
0.02 0
−10dB
−0.02
0
−15dB
−0.06
0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4
0.06
−15dB 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4
z [m]
0dB
0.06
0.04
0dB
0.04 −5dB
0.02 0
−10dB
−0.02 −0.04
−5dB
0.02
r [m]
r [m]
−10dB
−0.04
z [m]
−0.06
−5dB
0.02
−0.02
−0.04 −0.06
0dB
0.04
r [m]
r [m]
0.04
0 −10dB
−0.02 −0.04
−15dB 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4
−0.06
z [m]
−15dB 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4
z [m]
Fig. 11. Sound source maps with DAMAS at near-cut-on frequencies. The left column presents the simulated results, and the right column presents the experimental results. The star denotes the location of the actual sound source.
the true sound source from z ¼ 0 m to z ¼ 0:4 m, as shown in Fig. 1. The sound source located at z ¼ 0:3 m and r ¼ 0:04 m. Using Eq. (16), the value of b in the nodes of the plane was calculated, which indicated the relative amplitude of the acoustic source strength. The beamforming maps were normalized to 0 dB at maximum. From the above results, it is seen that the localization of sound source in the duct can be summarized as three situations: low frequencies, high frequencies, and near cut-on frequencies. The relatively high and low of frequency depends on the number of propagating spinning modes in the duct. Technically speaking, the dimensionless number, Helmohltz number kr, should be considered in the case of sound source localization, where k is the wave number and r is the radius of the duct. When kr < 5:3, the
0.06
number of the spinning modes that can propagate in the duct is less than six, the information is not enough to locate the sound source. Therefore, kr < 5:3 can be regarded as low frequency, which corresponds to f < 4800 Hz when r ¼ 0:06 m. For each situation, the results of two frequencies are presented below. Figs. 6–8 show the beamforming results of reconstructed source distributions with simulations and experiments at six different frequencies. At low frequencies, such as f ¼ 3600 Hz and f ¼ 4300 Hz, the main lobe is slightly wider, which means the spatial resolution is not very high. Meanwhile, there are side lobes at other positions. In this situation, the number of spinning modes that can propagate within the duct is small, while the remaining spinning modes will decay exponentially. The propagating modes are not adequate to
0dB
0.06
−5dB
0.02 0
−10dB
−0.02
0
−0.06
0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4
0.06
−15dB 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4
z [m]
0dB
0.06
0.04
0dB
0.04 −5dB
0.02 0
−10dB
−0.02 −0.04
r [m]
r [m]
−10dB
−0.04 −15dB
z [m]
−0.06
−5dB
0.02
−0.02
−0.04 −0.06
0dB
0.04
r [m]
r [m]
0.04
−5dB
0.02 0
−10dB
−0.02 −0.04
−15dB 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4
z [m]
−0.06
−15dB 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4
z [m]
Fig. 12. Sound source maps with the application of DAMAS to the results of the improved algorithm at near cut-on frequencies. The left column presents the simulated results, and the right column presents the experimental results. The star denotes the location of the actual sound source.
46
D. Zhong et al. / Applied Acoustics 103 (2016) 37–46
locate the source. Therefore, the results are not very good. With the increase of the frequency, the spatial resolution improves, and the side-lobes disappear. Note that the beamforming is effective at locating the sound source at relatively high frequencies, as shown in Fig. 7. However, when the frequency is close to the cut-on frequencies, the results significantly deteriorate. Consider the examples of Fig. 8, when f ¼ 7800 Hz, which close to the cut-on frequency of the mode (7, 1), 7736 Hz, the acoustic pressure is dominated by the mode (7, 1). When f ¼ 8700 Hz, which equals the cut-on frequency of the mode (8, 1), the acoustic pressure is dominated by the mode (8, 1). The location of the sound source cannot be determined by this single mode. To address the above phenomena at low frequencies and near cut-on frequencies, DAMAS and the developed algorithm are applied to improve the results. Figs. 9 and 11 show the consequences after using DAMAS on the beamforming results at low frequencies and near cut-on frequencies. DAMAS clearly removes the influence of the side-lobes and improves the spatial resolution at low frequencies, while at near cut-on frequencies, applying DAMAS alone would lead to maps with many ghost sound sources. Therefore, DAMAS alone is useless at the near-cut-on frequencies. In this situation, the developed algorithm is applied. The measured pressure and the Green function are rearranged to eliminate the influence of the single mode. The results are shown in Fig. 10, which clearly indicate that after use of the algorithm, the localization of the sound source becomes much more accurate. DAMAS can be applied further to the results of the improved algorithm, as shown in Fig. 12. With DAMAS, the resolution becomes much higher. 5. Conclusion In this paper, an in-duct beamforming method was used for the localization of the sound source in a finite duct. This beamforming method allows the beam to be focused at a point within the duct to determine the positions of the acoustic sources. The Green function that appears in the beamforming formula is an analytical in-duct Green function that considers the reflection of the open ends. The results of the azimuthal sound source localization and the planar sound source localization indicate that the accuracy of the localization is highly dependent on the frequency, which has effect on the number of the spinning modes that are able to propagate in the duct. At low frequencies, the number of propagating modes is not sufficient to locate the source. However, DAMAS can be used to improve the resolution and remove the side-lobes. When the frequency is close to a cut-on frequency of the spinning mode, the measured pressure is dominated by the single mode and the results rapidly deteriorate; also, for the azimuthal localization, side lobes appear. Applying DAMAS alone would lead to the results with ghost sound sources and cannot eliminate the side lobes. An improved algorithm is developed to address this problem. In the algorithm, the terms correspond with the near-cut-on mode are discarded in the pressure and the Green function to remove the influence of this single mode. For the Green function, this term can be calculated directly based on the expression, while for the pressure, the mode decomposition is considered to calculate the term. Significant improvements can be achieved when the algorithm is used. The resolution of the localization will be improved further by applying DAMAS to the results of the improved algorithm.
Acknowledgments The authors gratefully acknowledge the financial support of the National Natural Science Foundation of China (Grant No. 51376107). References [1] Hewlett DAK, Nelson PA. A model based approach for estimating the strength of acoustic sources within a circular duct. In: The 2nd AIAA/CEAS aeroacoustics conference, State College, PA, USA; 1996. [2] Kim Y, Nelson PA. Estimation of acoustic source strength within a cylindrical duct by inverse methods. J Sound Vib 2004;275(1):391–413. [3] Lowis CR, Joseph PF. Determining the strength of rotating broadband sources in ducts by inverse methods. J Sound Vib 2006;295(3):614–32. [4] Bravo T, Maury C. In-duct acoustic source reconstruction from pressure and velocity-based measurements. In: The 13th AIAA/CEAS aeroacoustics conference, Rome, Italy; 2007. [5] Sijtsma P. Using phased array beamforming to locate broadband noise sources inside a turbofan engine. Technical report, National Aerospace Laboratory (Netherlands). Based on a presentation held at AARC engine noise phased array workshop, Cambridge, MA, USA; 11–12 May 2006. [6] Sijtsma P. Feasibility of in-duct beamforming. In: The 13th AIAA/CEAS aeroacoustics conference, Rome, Italy; 2007. [7] Lowis CR, Joseph PF. A focused beamformer technique for separating rotor and stator-based broadband sources. In: The 12th AIAA/CEAS aeroacoustics conference, Cambridge, Massachusetts, USA; 2006. [8] Dougherty RP, Walker BE. Virtual rotating microphone imaging of broadband fan noise. In: The 15th AIAA/CEAS aeroacoustics conference, Miami, Florida, USA; 2009. [9] Dougherty RP. Walker BE. Sutliff DL. Locating and quantifying broadband fan sources using in-duct microphones. In: The 16th AIAA/CEAS aeroacoustics conference, Stockholm, Sweden; 2010. [10] Brooks TF, Humphreys WM. A deconvolution approach for the mapping of acoustic sources (DAMAS) determined from phased microphone arrays. J Sound Vib 2006;294(4):856–79. [11] Brooks TF, Humphreys WM. Three dimensional application of DAMAS methodology for aeroacoustic noise source definition. In: The 11th AIAA/ CEAS aeroacoustics conference, Monterey, California, USA; 2005. [12] Dougherty RP. Extensions of DAMAS and benefits and limitations of deconvolution in beamforming. In: The 11th AIAA/CEAS aeroacoustics conference, Monterey, California, USA; 2005. [13] Brooks TF, Humphreys W. Extension of DAMAS phased array processing for spatial coherence determination (DAMAS-C). In: The 12th AIAA/CEAS aeroacoustics conference, Cambridge, Massachusetts, USA; 2006. [14] Yardibi T, Li J, Stoica P, Iii LNC. Sparsity constrained deconvolution approaches for acoustic source mapping. J Acoust Soc Am 2008;123(5):2631–42. [15] Lawson CL, Hanson RJ. Solving least squares problems. Prentice-Hall Series in Automatic Computation. Englewood Cliffs: Prentice-Hall; 1974. [16] Sijtsma P. Clean based on spatial source coherence. Int J Aeroacoust 2007;6 (4):357–74. [17] Sijtsma P. Using phased array beamforming to identify broadband noise sources in a turbofan engine. Int J Aeroacoust 2010;9(3):357–74. [18] Doak PE. Excitation, transmission and radiation of sound from source distributions in hard-walled ducts of finite length (i): the effects of duct cross-section geometry and source distribution space-time pattern. J Sound Vib 1973;31:1–72. [19] Zorumski WE. Generalized radiation impedances and reflection coefficients of circular and annular ducts. J Acoust Soc Am 1973;54(6):1667–73. [20] Sijtsma P. Circular harmonics beamforming with multiple rings of microphones. In: The 18th AIAA/CEAS aeroacoustics conference, Colorado Springs, CO, USA; 2012. [21] Nelson PA, Yoon SH. Estimation of acoustic source strength by inverse methods: part I, conditioning of the inverse problem. J Sound Vib 2000;233 (4):639–64. [22] Hansen PC. Rank-deficient and discrete ill-posed problems: numerical aspects of linear inversion. Vol. 4, Society for Industrial Mathematics; 1987. [23] Hansen PC, O’Leary DP. The use of the l-curve in the regularization of discrete ill-posed problems. SIAM J Scient Comput 1993;14(6):1487–503. [24] Yoon SH, Nelson PA. Estimation of acoustic source strength by inverse methods: Part II, experimental investigation of methods for choosing regularization parameters. J Sound Vib 2000;233(4):665–701. [25] Lowis C. In-duct measurement techniques for the characterisation of broadband aeroengine noise. Ph.D. thesis, University of Southampton; 2008.