Iterative reweighted atomic norm minimization based two-dimensional multiple-snapshot grid-free compressive beamforming with planar microphone array

Abstract

Compressive beamforming with a planar microphone array can effectively estimate the two-dimensional directions-of-arrival and quantify the strengths of acoustic sources. Due to the superiorities of overcoming the basis mismatch issue of the conventional grid-based method and improving the performance of the single-snapshot grid-free method, the multiple-snapshot grid-free method has become the current research focus. Its existing atomic norm minimization (ANM) based strategy blocks high resolution because it cannot work well for sources with a small separation. This paper commits itself to remedying this drawback. After revealing the cause for this drawback, we present an iterative reweighted ANM (IRANM) approach. Both simulations and experiments demonstrate that compared with the ANM-based two-dimensional multiple-snapshot grid-free compressive beamforming, the IRANM-based one enjoys not only the enhanced resolution but also the stronger denoising ability and the higher identification accuracy.

Introduction

Compressive beamforming with a planar microphone array is a promising approach to estimate the directions-of-arrival (DOAs) and quantify the strengths of acoustic sources. Based on compressive sensing theory,^1–3 the approach recovers the source distribution from an underdetermined system of linear equations by imposing sparsity constraints. The system of linear equations relates the noisy measurements of the wavefield with a microphone array to the source distribution. Due to the superiorities of strong anti-interference, low demand for the number of microphones, clear source imaging, and so on, it has attracted extensive concern recently.^4–9

Conventional compressive beamforming grids/discretizes the DOA domain into a finite set of observation directions, and assumes that all sources are located in these observation directions. The sparsity constraint is imposed to minimize the $l_1-$ norm of a source strength vector (for the single-snapshot case) or the $l_{2\text{,}1}-$ norm of a source strength matrix (for the multiple-snapshot case). When the DOAs of sources do not fall in these observation directions, the results become inaccurate. The problem is called basis mismatch.^4,10–12 Utilizing a finer grid^13,14 can only alleviate this basis mismatch problem to a certain extent, and will increase computational complexity. More seriously, grid refinement burdens the measuring process with increased coherence, which can cause offset in the estimates. ⁴ In order to fundamentally solve the basis mismatch problem, the grid-free strategy is presented.^15–17 In 2015, based on the minimization of the atomic norm of source strength and the polynomial rooting method, Xenaki and Gerstoft ¹⁸ developed a one-dimensional (1D) single-snapshot grid-free compressive beamforming (GFCB) for a uniform linear microphone arrays (ULMAs). In 2018, Park, Choo, and Seong ¹⁹ extended Xenaki and Gerstoft’s strategy into the 1D multiple-snapshot case. In 2017–2018, the authors^20–23 successively developed the two-dimensional (2D) single- and multiple-snapshot GFCB based on atomic norm minimization (ANM) for a uniform rectangular microphone array (URMA). They have proven that for sources that are active across multiple snapshots and uncorrelated or partially correlated, the multiple-snapshot GFCB outperforms the single-snapshot one in many aspects. ²³

Three steps are included in the ANM-based 2D multiple-snapshot GFCB. ²³ First, the atomic norm of the microphone pressure induced by sources is minimized to obtain a positive semidefinite (PSD) and twofold Toeplitz matrix that can be understood as the array covariance matrix, and meanwhile denoise the measured pressure and thus reconstruct the microphone pressure from sources. Then, the Toeplitz matrix is processed to estimate the number and the DOAs of sources through the matrix pencil and pairing (MaPP) method.^23,24 Finally, the source strengths are quantified based on the reconstructed microphone pressure and the estimated DOAs. When the sources are sufficiently separated, ANM amounts to the sparsity constraint, the obtained PSD and twofold Toeplitz matrix can contain the accurate DOA information and the microphone pressure from sources can be accurately reconstructed with a high probability. However, it fails to identify sources that have a small separation. This drawback can be attributed to the fact that the atomic norm in the first step is only a convex relaxation of the source sparsity metric. This drawback limits the resolution of the ANM-based 2D multiple-snapshot GFCB. Exploring an effective method to remedy this drawback is of significance for the function perfection and performance improvement of the 2D multiple-snapshot GFCB. In this paper, an iterative reweighted ANM (IRANM)-based 2D multiple-snapshot GFCB is proposed.

The remainder of this paper is organized as follows. Section 2 presents the theory of the ANM-based 2D multiple-snapshot GFCB, and demonstrates and explains its drawback with simulations. Section 3 develops the IRANM-based 2D multiple-snapshot GFCB by illuminating its theory and revealing its mechanism with a simulated illustrative example. Section 4 compares the two strategies both with simulations and experimentally. Section 5 concludes this paper.

ANM-based 2D multiple-snapshot GFCB

Theory

In the 2D GFCB, a URMA is utilized to measure the signal. Figure 1 shows the layout, and the symbol “•” represents the microphone. $a=\mathrm{0,1},\mathrm{\dots },A-1$ and $b=\mathrm{0,1},\mathrm{\dots },B-1$ index the microphones in $x$ and $y$ dimensions, respectively. The microphone spacing is $\mathrm{\Delta }x$ and $\mathrm{\Delta }y$ in $x$ and $y$ dimensions, respectively. $\left({\theta }_i,{\varphi }_i\right)$ indicates the actual DOA of the $i-$ th source with ${\theta }_i$ and ${\varphi }_i$ being the elevation and the azimuth angle, respectively. $s_{i,l}$ represents the complex strength of the $i-$ th source under the $l-$ th snapshot, and ${\mathbf{s}}_i=\left[s_{i,1},s_{i,2},\mathrm{\dots },s_{i,L}\right]\in {\mathrm{C}}^{1\times L}$ represents the complex strength vector of the $i-$ th source under all $L$ snapshots. Here, $\mathrm{C}$ represents a set of complex numbers. The theoretical pressure vector ${\mathbf{p}}_{a,b}\in {\mathrm{C}}^{1\times L}$ at the (a,b)-th microphone under each snapshot can be described by $${\mathbf{p}}_{a,b}={\displaystyle \sum _{i=1}^I{\mathbf{s}}_i{e}^{j2\pi \left(t_{1i}a+t_{2i}b\right)}}$$ (1)where $I$ is the total number of sources, $j=\sqrt{-1}$ , $t_{1i}\equiv \sin {\theta }_i\cos {\varphi }_i\mathrm{\Delta }x/\lambda$ , $t_{2i}\equiv \sin {\theta }_i\sin {\varphi }_i\mathrm{\Delta }y/\lambda$ , and $\lambda$ is the wavelength. Construct $\mathbf{d}\left(t_{1i},t_{2i}\right)={\left[1,e^{j2\pi t_{1i}},\mathrm{\dots },e^{j2\pi t_{1i}\left(A-1\right)}\right]}^{\text{T}}\mathrm{\otimes }{\left[1,e^{j2\pi t_{2i}},\mathrm{\dots },e^{j2\pi t_{2i}\left(B-1\right)}\right]}^{\text{T}}\in {\mathrm{C}}^{AB}$ , $\mathbf{P}={\left[{\mathbf{p}}_{\mathrm{0,0}}^{\text{T}},{\mathbf{p}}_{\mathrm{0,1}}^{\text{T}},\mathrm{\dots },{\mathbf{p}}_{0,B-1}^{\text{T}},{\mathbf{p}}_{\mathrm{1,0}}^{\text{T}},{\mathbf{p}}_{\mathrm{1,1}}^{\text{T}},\mathrm{\dots },{\mathbf{p}}_{1,B-1}^{\text{T}},\mathrm{\dots },{\mathbf{p}}_{A-\mathrm{1,0}}^{\text{T}},{\mathbf{p}}_{A-\mathrm{1,1}}^{\text{T}},\mathrm{\dots },{\mathbf{p}}_{A-1,B-1}^{\text{T}}\right]}^{\text{T}}\in {\mathrm{C}}^{AB\times L}$ , $s_i={\Vert {\mathbf{s}}_i\Vert }_2\in {\mathrm{ℝ}}^{+}$ , and ${\mathbf{\psi }}_i={\mathbf{s}}_i/s_i\in {\mathrm{C}}^{1\times L}$ , where ${(\cdot )}^{\text{T}}$ represents the transpose operator, $\otimes$ represents the Kronecker product, ${\Vert \cdot \Vert }_2$ represents the $l_2-$ norm, and ${\mathrm{R}}^{+}$ represents a positive real number set. Corresponding to equation (1) $$\mathbf{P}={\displaystyle \sum _{i=1}^I\mathbf{d}\left(t_{1i},t_{2i}\right){\mathbf{s}}_i}={\displaystyle \sum _{i=1}^I{s}_i\mathbf{d}\left(t_{1i},t_{2i}\right){\mathbf{\psi }}_i}$$ (2)

Figure 1.

Measurement layout.

The measured pressure ${\mathbf{P}}^{\mathrm{★}}\in {\mathrm{C}}^{AB\times L}$ can be expressed as the sum of the theoretical pressure $\mathbf{P}$ and the additive noise $\mathbf{N}\in {\mathrm{C}}^{AB\times L}$ , that is $${\mathbf{P}}^{\mathrm{★}}=\mathbf{P}+\mathbf{N}$$ (3)

In subsequent simulations, the noise is assumed to be independent identically distributed complex Gaussian. Define the array signal-to-noise ratio (SNR) as $\mathrm{20}{\log }_{10}\left({\Vert \mathbf{P}\Vert }_{\text{F}}/{\Vert \mathbf{N}\Vert }_{\text{F}}\right)$ , where ${\Vert \cdot \Vert }_{\text{F}}$ represents the Frobenius norm.

The first postprocessing step of the 2D multiple-snapshot GFCB is to denoise ${\mathbf{P}}^{\mathrm{★}}$ and thus reconstruct $\mathbf{P}$ . This can be realized by imposing a sparsity constraint, that is, minimizing the number of sources. In the grid-free frame, $t_1\equiv \sin \theta \cos \varphi \mathrm{\Delta }x/\lambda$ , $t_2\equiv \sin \theta \sin \varphi \mathrm{\Delta }y/\lambda$ and the elements of $\mathbf{\psi }$ are continuous functions related to $\theta$ and $\varphi$ . Thereinto, ${-90}^{\circ }\le \theta \le {90}^{\circ }$ and ${\mathrm{0}}^{\circ }\le \phi <\mathrm{3}6{\mathrm{0}}^{\circ }$ . Assemble all the values of $\mathbf{d}\left(t_1,t_2\right)\mathbf{\psi }$ to form an atomic set and represented as $\mathrm{A}$ . The atomic $l_{\text{0}}-$ norm of $\mathbf{P}$ is defined as $${\Vert \mathbf{P}\Vert }_{\mathrm{A},0}=\underset{\begin{array}{c}\mathbf{d}\left(t_{1i},t_{2i}\right){\mathbf{\psi }}_i\in \mathrm{A}\\ s_i\in {\mathrm{R}}^{+}\end{array}}{\inf }\left\{I|\mathbf{P}={\displaystyle \sum _{i=1}^I{s}_i\mathbf{d}\left(t_{1i},t_{2i}\right){\mathbf{\psi }}_i}\right\}$$ (4)where $\inf \text{{}\cdot \text{}}$ represents the infimum operator. ${\Vert \mathbf{P}\Vert }_{\mathrm{A},0}$ is a direct metric of source sparsity. The reconstruction of $\mathbf{P}$ can be transformed into $$\widehat{\mathbf{P}}=\underset{\mathbf{P}\in {\mathrm{C}}^{AB\times L}}{\arg \min }{\Vert \mathbf{P}\Vert }_{\mathrm{A},0}\hspace{1em}\text{subject to }{\Vert {\mathbf{P}}^{\mathrm{★}}-\mathbf{P}\Vert }_{\text{F}}\le \epsilon$$ (5)where $\epsilon$ is the control parameter corresponding to noise. Let $\epsilon ={\Vert \mathbf{N}\Vert }_{\text{F}}$ normally. Equation (5) can be further transformed into the following rank minimization problem $$\left\{\widehat{\mathbf{u}},\widehat{\mathbf{P}},\widehat{\mathbf{E}}\right\}=\underset{\mathbf{u}\in {\mathrm{C}}^{N_u},\mathbf{P}\in {\mathrm{C}}^{AB\times L},\mathbf{E}\in {\mathrm{C}}^{L\times L}}{\arg \min }\text{rank}\left(T_b\left(\mathbf{u}\right)\right)\hspace{1em}\text{subject to}\hspace{1em}\left[\begin{array}{cc}{T}_b\left(\mathbf{u}\right)& \mathbf{P}\\ {\mathbf{P}}^{\text{H}}& \mathbf{E}\end{array}\right]\ge 0,{\Vert {\mathbf{P}}^{\mathrm{★}}-\mathbf{P}\Vert }_{\text{F}} \le \epsilon$$ (6)where $\mathbf{u}\in {\mathrm{C}}^{N_u}$ with $N_u=\left(A-1\right)\left(2B-1\right)+B$ and $\mathbf{E}\in {\mathrm{C}}^{L\times L}$ are auxiliary quantities, $\text{rank}(\cdot )$ represents to compute the matrix rank, $T_b(\cdot )$ represents the twofold Toeplitz operator,^20–25 ${(\cdot )}^{\text{H}}$ represents the conjugate operator, and $\ge 0$ means that the matrix is PSD. The sparsity is exploited to the best extent possible by ${\Vert \mathbf{P}\Vert }_{\mathrm{A},0}$ . Nevertheless, ${\Vert \mathbf{P}\Vert }_{\mathrm{A},0}$ is non-convex and equation (6) is nondeterministic polynomial-time problem and hard to solve. Therefore, a computationally feasible alternative is required. The atomic norm ${\Vert \mathbf{P}\Vert }_{\mathrm{A}}$ is employed as a convex relaxation of ${\Vert \mathbf{P}\Vert }_{\mathrm{A},0}$ in Ref. [23]. Its definition is $${\Vert \mathbf{P}\Vert }_{\mathrm{A}}=\underset{\begin{array}{c}\mathbf{d}\left(t_{1i},t_{2i}\right){\mathbf{\psi }}_i\in \mathrm{A}\\ s_i\in {\mathrm{R}}^{+}\end{array}}{\inf }\left\{{\displaystyle \sum _i{s}_i}|\mathbf{P}={\displaystyle \sum _i{s}_i\mathbf{d}\left(t_{1i},t_{2i}\right){\mathbf{\psi }}_i}\right\}$$ (7)

Correspondingly, equation (5) becomes $$\widehat{\mathbf{P}}=\underset{\mathbf{P}\in {\mathrm{C}}^{AB\times L}}{\arg \min }{\Vert \mathbf{P}\Vert }_{\mathrm{A}}\hspace{1em}\text{subject to}\hspace{1em}{\Vert {\mathbf{P}}^{\mathrm{★}}-\mathbf{P}\Vert }_{\text{F}} \le \epsilon$$ (8)

Equation (8) can be transformed into the following trace minimization problem ²³ $$\left\{\widehat{\mathbf{u}},\widehat{\mathbf{P}},\widehat{\mathbf{E}}\right\}=\underset{\mathbf{u}\in {\mathrm{C}}^{N_u},\mathbf{P}\in {\mathrm{C}}^{AB\times L},\mathbf{E}\in {\mathrm{C}}^{L\times L}}{\arg \min }\frac{1}{2\sqrt{AB}}\left(\text{tr}\left(T_b\left(\mathbf{u}\right)\right)+\text{tr}\left(\mathbf{E}\right)\right)\hspace{1em}\text{subject to}\hspace{1em}\left[\begin{array}{cc}{T}_b\left(\mathbf{u}\right)& \mathbf{P}\\ {\mathbf{P}}^{\text{H}}& \mathbf{E}\end{array}\right]\ge 0,\hspace{1em}{\Vert {\mathbf{P}}^{\mathrm{★}}-\mathbf{P}\Vert }_{\text{F}}\le \epsilon$$ (9)where $\text{tr}(\cdot )$ represents to compute the matrix trace. Obviously, substituting ${\Vert \mathbf{P}\Vert }_{\mathrm{A}}$ for ${\Vert \mathbf{P}\Vert }_{\mathrm{A},0}$ is substituting the trace minimization for the rank minimization. Equation (9) is a disciplined convex optimization problem ²⁶ and can be solved utilizing the SDPT3 solver in MATLAB CVX toolbox. ²⁷

The $T_b\left(\widehat{\mathbf{u}}\right)$ obtained by equation (9) is a PSD and twofold Toeplitz matrix. Its Vandermonde decomposition,^{20,21,23–25} expressed by equation (10), ²¹ contains the DOA information of sources $$T_b\left(\widehat{\mathbf{u}}\right)=\mathbf{V}\mathbf{\Sigma }{\mathbf{V}}^{\text{H}}={\displaystyle \sum _{i=1}^r{\sigma }_i\mathbf{d}\left(t_{1i},t_{2i}\right)\mathbf{d}{\left(t_{1i},t_{2i}\right)}^{\text{H}}}$$ (10)where $\mathbf{V}=\left[\mathbf{d}\left(t_{11},t_{21}\right),\mathbf{d}\left(t_{12},t_{22}\right),\mathrm{\dots },\mathbf{d}\left(t_{1r},t_{2r}\right)\right]$ , $\mathbf{\Sigma }=\mathbf{diag}\left(\left[{\sigma }_1,{\sigma }_2,\mathrm{\dots },{\sigma }_r\right]\right)$ , $\mathbf{diag}(\cdot )$ creates a square diagonal matrix with diagonal elements being the vector in the parentheses, ${\sigma }_i\left(i=\mathrm{1,2},\mathrm{\dots },r\right)\in {\mathrm{R}}^{+}$ , and $r$ corresponds to the rank of $T_b\left(\widehat{\mathbf{u}}\right)$ . Obviously, $T_b\left(\widehat{\mathbf{u}}\right)$ can be understood as the sum of the array covariance matrix (the covariance matrix of microphone pressure) induced by each source in a group. The second postprocessing step of the 2D multiple-snapshot GFCB is to find the Vandermonde decomposition of $T_b\left(\widehat{\mathbf{u}}\right)$ and thus estimate the number and DOAs of sources. This can be realized through the MaPP method. Its detailed procedure can be seen in Ref. [23]. It is worth mentioning that not all the PSD and twofold Toeplitz matrices admit the Vandermonde decomposition. Refs. [24,25] have proved that the sufficient condition for $T_b\left(\widehat{\mathbf{u}}\right)$ to support a Vandermonde decomposition is $r\le \min \{A,B\}$ , and the sufficient condition for $T_b\left(\widehat{\mathbf{u}}\right)$ to admit a unique Vandermonde decomposition is $r\le \min \left\{A,B\right\}$ . Based on the $\widehat{\mathbf{P}}$ reconstructed by ANM and the $\mathbf{V}$ obtained by MaPP, the 2D multiple-snapshot GFCB quantifies the source strength as $\widehat{\mathbf{P}}$ divided by $\mathbf{V}$ .

An array formed by randomly selecting partial microphones from a standard URMA shown in Figure 1 is called a non-uniform rectangular microphone array (NURMA). Employing the measured pressures of the partial microphones, the 2D multiple-snapshot GFCB can also obtain the $T_b\left(\widehat{\mathbf{u}}\right)$ , reconstruct the pressure of all microphones, and finally achieve the estimation of source DOA and strength. In this circumstance, $\mathbf{P}$ , ${\mathbf{P}}^{\mathrm{★}}$ and $\mathbf{N}$ in equations (2) and (3) become ${\mathbf{P}}_{\mathrm{\Omega }}$ , ${\mathbf{P}}_{\mathrm{\Omega }}^{\mathrm{★}}$ and ${\mathbf{N}}_{\mathrm{\Omega }}$ , and ${\mathbf{P}}^{\mathrm{★}}-\mathbf{P}$ in equations (5), (6), (8) and (9) becomes , where $\mathrm{\Omega }$ is the index set of the chosen microphones, $\left|\mathrm{\Omega }\right|$ is the cardinality of $\mathrm{\Omega }$ , ${\mathbf{P}}_{\mathrm{\Omega }}^{\mathrm{★}}\in {\mathrm{C}}^{\left|\mathrm{\Omega }\right|\times L}$ is the matrix formed by the measured pressures of the selected microphones, ${\mathbf{P}}_{\mathrm{\Omega }}\in {\mathrm{C}}^{\left|\mathrm{\Omega }\right|\times L}$ is the matrix composed of the pressures induced by sources at the selected microphones, and ${\mathbf{N}}_{\mathrm{\Omega }}\in {\mathrm{C}}^{\left|\mathrm{\Omega }\right|\times L}$ is the matrix composed of noises suffered by the selected microphones.

Simulations of drawback

A simulation example with six sound sources is designed. The DOAs of six sources are (70°, 300°), (30°, 200°), (40°, 90°), (50°, 90°), (60°, 190°), and (60°, 180°), in turn. The corresponding root mean square strengths $\left({\Vert {\mathbf{s}}_i\Vert }_2/\sqrt{L}\right)$ are 90 dB, 94 dB, 94 dB, 97 dB, 97 dB and 100 dB (referring to 2 × 10^-5 Pa). Set the frequency to 2000 Hz, the microphone array to eight rows and eight columns $\left(\mathrm{i}\mathrm{.e}\mathrm{.,}\hspace{0.17em}A=B=8\right)$ , the microphone spacing to 0.035 m $\left(\mathrm{\Delta }x=\mathrm{\Delta }y=0.035\text{ m}\right)$ , the SNR to 30 dB, and the number of snapshots to 10, respectively. Figure 2(a) presents the source DOA map reconstructed by the ANM-based 2D multiple-snapshot GFCB. It has severe distortion. That can be explained as follows. The minimum separation among these sources, defined as ${\mathrm{\Delta }}_{\min }=\underset{i,i^{\prime }\in \left\{\mathrm{1,2},\mathrm{\dots },k\right\},i\ne i^{\prime }}{\min }\max \left\{\left|t_{1i}-t_{1i^{\prime }}\right|,\left|t_{2i}-t_{2i^{\prime }}\right|\right\}$ , is only 0.025.

In this case, the ANM in equation (8) does not amount to the atomic $l_{\text{0}}-$ norm minimization in equation (5). Correspondingly, the trace minimization in equation (9) does not amount to the rank minimization in equation (6). As shown by the symbol “□” in Figure 2(e), the obtained $T_b\left(\widehat{\mathbf{u}}\right)$ has many relatively large eigenvalues, which means the rank of $T_b\left(\widehat{\mathbf{u}}\right)$ is high. Consequently, $T_b\left(\widehat{\mathbf{u}}\right)$ may not support a Vandermonde decomposition, or its Vandermonde decomposition does not contain the correct DOA information. Change the fourth source’s elevation angle to 30° and the second source’s azimuth angle to 200° for simulation. Figure 2(b) presents the result. In this case, ${\mathrm{\Delta }}_{\min }$ increases to 0.054, the Vandermonde decomposition of $T_b\left(\widehat{\mathbf{u}}\right)$ contains the correct DOA information, and therefore the reconstructed source DOAs coincide with the true ones well. When simulating Figure 2(a) and (b), we directly specify the number of sources as 6 to avoid the effect of the inaccurate estimation. In practical applications, the MaPP method utilizes the number of the eigenvalues of $T_b\left(\widehat{\mathbf{u}}\right)$ in a preset dynamic range, that is, $\left|\left\{\lambda \left(T_b\left(\widehat{\mathbf{u}}\right)\right)|10\lg \left({\lambda }_{\max }\left(T_b\left(\widehat{\mathbf{u}}\right)\right)/\lambda \left(T_b\left(\widehat{\mathbf{u}}\right)\right)\right)\le R\right\}\right|$ , to estimate the number of sources, where $\lambda \left(T_b\left(\widehat{\mathbf{u}}\right)\right)$ and ${\lambda }_{\max }\left(T_b\left(\widehat{\mathbf{u}}\right)\right)$ denote the eigenvalue and the largest eigenvalue of $T_b\left(\widehat{\mathbf{u}}\right)$ , $R$ is the preset dynamic range, and $|\cdot |$ calculates the cardinality of a set. The eigenvalues of $T_b\left(\widehat{\mathbf{u}}\right)$ decrease slowly when sources have a small separation, as shown in Figure 2(e). As a result, it is very hard to give out an appropriate $R$ to estimate the number of sources robustly and reliably. As shown by Figure 2(c) and (d), whether the number of sources is underestimated or overestimated, the reconstruction is inaccurate. To sum up, for sources with a small separation, because the $T_b\left(\widehat{\mathbf{u}}\right)$ obtained by ANM cannot contain the correct DOA information with a high probability or the number of sources cannot be estimated robustly and reliably, the ANM-based 2D multiple-snapshot GFCB fails to work well. To obtain Figure 2, all the 64 microphones are utilized. Figure 3 shows the results when only 30 randomly selected microphones are utilized, which exhibits the same phenomena as in Figure 2, demonstrating the drawback of the ANM-based two-dimensional multiple-snapshot grid-free compressive beamforming again.

Figure 2.

Results of ANM-based 2D multiple-snapshot grid-free GFCB when a standard URMA with 64 microphones is utilized. (a), (b), (c), (d) Reconstructions (*) for sources (○). (e) Eigenvalues of $T_b\left(\widehat{\mathbf{u}}\right)$ . ${\mathrm{\Delta }}_{\min }$ is (a) 0.025 and (b), (c), (d) 0.054. The estimated number of sources is (a), (b) 6, (c) $\left|\left\{\lambda \left(T_b\left(\widehat{\mathbf{u}}\right)\right)|10\lg \left({\lambda }_{\max }\left(T_b\left(\widehat{\mathbf{u}}\right)\right)/\lambda \left(T_b\left(\widehat{\mathbf{u}}\right)\right)\right)\le 6\text{ dB}\right\}\right|$ and (d) $\left|\left\{\lambda \left(T_b\left(\widehat{\mathbf{u}}\right)\right)|10\lg \left({\lambda }_{\max }\left(T_b\left(\widehat{\mathbf{u}}\right)\right)/\lambda \left(T_b\left(\widehat{\mathbf{u}}\right)\right)\right)\le 20\text{ dB}\right\}\right|$ . In the DOA maps (a)–(d), the true source strength output and the reconstructed one are scaled to dB with reference to their respective maximum, and the reconstructed maximum sound level (referring to 2×10⁻⁵ Pa) is marked on the top. Similar in subsequent DOA maps.

Figure 3.

Results of ANM-based 2D multiple-snapshot GFCB when a NURMA with 30 randomly selected microphones is utilized. (a) Layout of the NURMA. (b), (c), (d), (e) Reconstructions (*) for sources (○). (f) Eigenvalues of $T_b\left(\widehat{\mathbf{u}}\right)$ . ${\mathrm{\Delta }}_{\min }$ is (b) 0.025 and (c), (d), (e) 0.054. The estimated number of sources is (b), (c) 6, (d) $\left|\left\{\lambda \left(T_b\left(\widehat{\mathbf{u}}\right)\right)|10\lg \left({\lambda }_{\max }\left(T_b\left(\widehat{u}\right)\right)/\lambda \left(T_b\left(\widehat{\mathbf{u}}\right)\right)\right)\le 6\text{ dB}\right\}\right|$ and (e) $\left|\left\{\lambda \left(T_b\left(\widehat{u}\right)\right)|10\lg \left({\lambda }_{\max }\left(T_b\left(\widehat{u}\right)\right)/\lambda \left(T_b\left(\widehat{u}\right)\right)\right)\le 20\text{ dB}\right\}\right|$ .

IRANM-based 2D multiple-snapshot GFCB

Theory

Build a new metric $${\mathrm{M}}^{\kappa }\left(\mathbf{P}\right)=\underset{\mathbf{u}\in {\mathrm{ℂ}}^{N_u},\mathbf{E}\in {\mathrm{ℂ}}^{L\times L}}{\min }\frac{1}{2\sqrt{AB}}\left(\text{ln}\left|T_b\left(\mathbf{u}\right)+\kappa \mathbf{I}\right|+\text{tr}\left(\mathbf{E}\right)\right)\hspace{1em}\text{subject to}\hspace{1em}\left[\begin{array}{cc}{T}_b\left(\mathbf{u}\right)& \mathbf{P}\\ {\mathbf{P}}^{\text{H}}& \mathbf{E}\end{array}\right]\ge 0$$ (11)where $\kappa >0$ is a regularization parameter, $\mathbf{I}\in {\mathrm{R}}^{AB\times AB}$ is an identity matrix. The following properties hold under a certain condition:^11,21 (1) when $\kappa \to \text{0}$ and ${\Vert \mathbf{P}\Vert }_{\mathrm{A},0}<AB$ , ${\mathrm{M}}^{\kappa }\left(\mathbf{P}\right)\sim \left({\Vert \mathbf{P}\Vert }_{\mathrm{A},0}/2\sqrt{AB}-\sqrt{AB}/2\right)\ln {\kappa }^{-1}$ , that is, ${\mathrm{M}}^{\kappa }\left(\mathbf{P}\right)$ is equal to $\left({\Vert \mathbf{P}\Vert }_{\mathrm{A},0}/2\sqrt{AB}-\sqrt{AB}/2\right)\ln {\kappa }^{-1}$ asymptotically as $\kappa$ tends to 0 if ${\Vert \mathbf{P}\Vert }_{\mathrm{A},0}<AB$ ; (2) When $\kappa \to \infty$ , ${\mathrm{M}}^{\kappa }\left(\mathbf{P}\right)-\left(\sqrt{AB}/2\right)\ln \kappa \sim {\Vert \mathbf{P}\Vert }_{\mathrm{A}}{\kappa }^{-1/2}$ , that is, ${\mathrm{M}}^{\kappa }\left(\mathbf{P}\right)-\left(\sqrt{AB}/2\right)\ln \kappa$ is equal to ${\Vert \mathbf{P}\Vert }_{\mathrm{A}}{\kappa }^{-1/2}$ asymptotically as $\kappa$ tends to $+\infty$ ; (3) Use ${\widehat{\mathbf{u}}}_{\kappa \to 0}$ to represent the optimal variable when $\kappa$ tends to 0. Only ${\Vert \mathbf{P}\Vert }_{\mathrm{A},0}$ eigenvalues of $T_b\left({\widehat{\mathbf{u}}}_{\kappa \to 0}\right)$ are relatively large. The smallest $AB-{\Vert \mathbf{P}\Vert }_{\mathrm{A},0}$ eigenvalues are around 0. The first property shows that minimizing ${\mathrm{M}}^{\kappa }\left(\mathbf{P}\right)$ amounts to minimizing ${\Vert \mathbf{P}\Vert }_{\mathrm{A},0}$ as $\kappa$ tends to 0. The second property shows that minimizing ${\mathrm{M}}^{\kappa }\left(\mathbf{P}\right)$ amounts to minimizing ${\Vert \mathbf{P}\Vert }_{\mathrm{A}}$ as $\kappa$ tends to $+\infty$ . That is to say, ${\mathrm{M}}^{\kappa }\left(\mathbf{P}\right)$ acts as the bridge between ${\Vert \mathbf{P}\Vert }_{\mathrm{A},0}$ and ${\Vert \mathbf{P}\Vert }_{\mathrm{A}}$ , and is capable of enhancing sparsity compared to ${\Vert \mathbf{P}\Vert }_{\mathrm{A}}$ . The third property shows that the sufficient condition to support a Vandermonde decomposition for $T_b\left({\widehat{\mathbf{u}}}_{\kappa \to 0}\right)$ can be guaranteed and employing the number of eigenvalues of $T_b\left({\widehat{\mathbf{u}}}_{\kappa \to 0}\right)$ that are larger than a tiny value can estimate the number of sources robustly and reliably. Therefore, it is expected to remedy the above-mentioned drawback by replacing ${\Vert \mathbf{P}\Vert }_{\mathrm{A}}$ with ${\mathrm{M}}^{\kappa }\left(\mathbf{P}\right)$ to measure the source sparsity. So, the reconstruction problem of $\mathbf{P}$ can be rewritten as $$\widehat{\mathbf{P}}=\underset{\mathbf{P}\in {\mathrm{C}}^{AB\times L}}{\arg \min }{\mathrm{M}}^{\kappa }\left(\mathbf{P}\right)\hspace{1em}\text{subject to}\hspace{1em}{\Vert {\mathbf{P}}^{\mathrm{★}}-\mathbf{P}\Vert }_{\text{F}}\le \epsilon$$ (12)

Simultaneous equations (11) and (12) bring $$\left\{\widehat{\mathbf{u}},\widehat{\mathbf{P}},\widehat{\mathbf{E}}\right\}=\underset{\mathbf{u}\in {\mathrm{C}}^{N_u},\mathbf{P}\in {\mathrm{C}}^{AB\times L},\mathbf{E}\in {\mathrm{C}}^{L\times L}}{\arg \min }\frac{1}{2\sqrt{AB}}\left(\text{ln}\left|T_b\left(\mathbf{u}\right)+\kappa \mathbf{I}\right|+\text{tr}\left(\mathbf{E}\right)\right)\text{ subject to }\left[\begin{array}{cc}{T}_b\left(\mathbf{u}\right)& \mathbf{P}\\ {\mathbf{P}}^{\text{H}}& \mathbf{E}\end{array}\right]\ge 0,\hspace{1em}{\Vert {\mathbf{P}}^{\mathrm{★}}-\mathbf{P}\Vert }_{\text{F}}\le \epsilon$$ (13)

The cost function in equation (13) is a combination of concave function $\text{ln}\left|T_b\left(\mathbf{u}\right)+\kappa \mathbf{I}\right|$ and convex function $\text{tr}\left(\mathbf{E}\right)$ , which can be solved by the majorization–minimization algorithm.^21,28,29 This algorithm iteratively updates the optimal variable and the optimal value of the cost function. Let ${\widehat{\mathbf{u}}}^k$ , ${\kappa }^k$ and ${\mathbf{W}}^k\equiv {\left(T_b\left({\widehat{\mathbf{u}}}^k\right)+{\kappa }^k\mathbf{I}\right)}^{-1}\in {\mathrm{C}}^{AB\times AB}$ be the optimal variable, the regularization parameter and the weighting matrix determined by the $k-$ th iteration, respectively. In the (k+1)-th iteration, a surrogate function is first constructed to locally approximate the cost function. Its expression is $$\frac{1}{2\sqrt{AB}}\left(\text{ln}\left|T_b\left({\widehat{\mathbf{u}}}^k\right)+{\kappa }^k\mathbf{I}\right|+\text{tr}\left({\mathbf{W}}^k{T}_b\left(\mathbf{u}-{\widehat{\mathbf{u}}}^k\right)\right)+\text{tr}\left(\mathbf{E}\right)\right)=\frac{1}{2\sqrt{AB}}\left(\text{tr}\left({\mathbf{W}}^k{T}_b\left(\mathbf{u}\right)\right)+\text{tr}\left(\mathbf{E}\right)\right)+c^k$$ (14)where $c^k$ is a constant. Then, ignoring $c^k$ , the minimization problem can be transformed into a disciplined convex optimization problem ²⁶ as follows $$\left\{{\widehat{\mathbf{u}}}^{k+1},{\widehat{\mathbf{P}}}^{k+1},{\widehat{\mathbf{E}}}^{k+1}\right\}=\underset{\mathbf{u}\in {\mathrm{C}}^{N_u},\mathbf{P}\in {\mathrm{C}}^{AB\times L},\mathbf{E}\in {\mathrm{C}}^{L\times L}}{\arg \min }\frac{1}{2\sqrt{AB}}\left(\text{tr}\left({\mathbf{W}}^k{T}_b\left(\mathbf{u}\right)\right)+\text{tr}\left(\mathbf{E}\right)\right)\hspace{1em}\text{subject to }\left[\begin{array}{cc}{T}_b\left(\mathbf{u}\right)& \mathbf{P}\\ {\mathbf{P}}^{\text{H}}& \mathbf{E}\end{array}\right]\ge 0,{\Vert {\mathbf{P}}^{\mathrm{★}}-\mathbf{P}\Vert }_{\text{F}}\le \epsilon$$ (15)

Equation (15) can be solved utilizing the SDPT3 solver in MATLAB CVX toolbox. ²⁷ Initializing ${\widehat{\mathbf{u}}}^0=\mathbf{0}$ and ${\kappa }^0=1$ , the first iteration solution to equation (15) is consistent with the ANM for equation (9). In order to enhance sparsity, let $\kappa$ gradually decrease during the iteration $${\kappa }^k=\{\begin{array}{l}\min \left\{{\kappa }^{k-1}/2,{\lambda }_{\max }\left(T_b\left({\widehat{\mathbf{u}}}^k\right)\right)/10\right\},\hspace{1em}k=1\\ {\kappa }^{k-1}/2\text{,}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}2\le k\le 10\\ {\kappa }^{10},\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}k>10\end{array}$$ (16)where ${\lambda }_{\max }\left(T_b\left({\widehat{\mathbf{u}}}^k\right)\right)$ denotes the largest eigenvalue of $T_b\left({\widehat{\mathbf{u}}}^k\right)$ .

Denoting by $w\left(t_{1i},t_{2i}\right)\ge \text{0}$ a weighting coefficient, we define the weighted atomic norm of $\mathbf{P}$ as $${\Vert \mathbf{P}\Vert }_{{\mathrm{A}}^w}=\underset{\begin{array}{c}\mathbf{d}(t_{1i},t_{2i}){\mathbf{\psi }}_i\in \mathrm{A}\\ s_i\in {\mathrm{R}}^{+}\end{array}}{\inf }\left\{{\displaystyle \sum _i\frac{s_i}{w\left(t_{1i},t_{2i}\right)}}|\mathbf{P}={\displaystyle \sum _i{s}_i\mathbf{d}\left(t_{1i},t_{2i}\right){\mathbf{\psi }}_i}\right\}$$ (17)

Obviously, when $w\left(t_{1i},t_{2i}\right)=1$ , ${\Vert \mathbf{P}\Vert }_{{\mathrm{A}}^w}={\Vert \mathbf{P}\Vert }_{\mathrm{A}}$ . Let ${\Vert \mathbf{P}\Vert }_{{\mathrm{A}}^{w^k}}$ and $w^k\left(t_1,t_2\right)$ be the weighted atomic norm and the weighting coefficient of the $k-$ th iteration. If $T_b\left({\widehat{\mathbf{u}}}^{k+1}\right)$ admits a Vandermonde decomposition and $$w^k\left(t_1,t_2\right)=\sqrt{\frac{AB}{\mathbf{d}{\left(t_1,t_2\right)}^{\text{H}}{\mathbf{W}}^k\mathbf{d}\left(t_1,t_2\right)}}$$ (18)then $${\Vert \mathbf{P}\Vert }_{{\mathrm{A}}^{w^k}}=\underset{\mathbf{u}\in {\mathrm{C}}^{N_u},\mathbf{E}\in {\mathrm{C}}^{L\times L}}{\min }\frac{1}{2\sqrt{AB}}\left(\text{tr}\left({\mathbf{W}}^k{T}_b\left(\mathbf{u}\right)\right)+\text{tr}\left(\mathbf{E}\right)\right)\hspace{1em}\text{subject to}\hspace{1em}\left[\begin{array}{cc}{T}_b\left(\mathbf{u}\right)& \mathbf{P}\\ {\mathbf{P}}^{\text{H}}& \mathbf{E}\end{array}\right]\ge 0$$ (19)

The proposition is proved in Appendix. According to the proposition, equation (15) can be rewritten as $$\begin{array}{l}{\widehat{\mathbf{P}}}^{k+1}=\underset{\mathbf{P}\in {\mathrm{C}}^{AB\times L}}{\arg \min }\left(\underset{\mathbf{u}\in {\mathrm{C}}^{N_u},\mathbf{E}\in {\mathrm{C}}^{L\times L}}{\min }\frac{1}{2\sqrt{AB}}\left(\text{tr}\left({\mathbf{W}}^k{T}_b\left(\mathbf{u}\right)\right)+\text{tr}\left(\mathbf{E}\right)\right)\hspace{1em}\text{subject to }\left[\begin{array}{cc}{T}_b\left(\mathbf{u}\right)& \mathbf{P}\\ {\mathbf{P}}^{\text{H}}& \mathbf{E}\end{array}\right]\ge 0\right)\\ \hspace{1em}\hspace{1em}\hspace{1em}\text{subject to }{\Vert {\mathbf{P}}^{\mathrm{★}}-\mathbf{P}\Vert }_{\text{F}}\le \epsilon \\ \hspace{1em}\hspace{1em}=\underset{\mathbf{P}\in {\mathrm{C}}^{AB\times L}}{\arg \min }{\Vert \mathbf{P}\Vert }_{{\mathrm{A}}^{w^k}}\text{ subject to }{\Vert {\mathbf{P}}^{\mathrm{★}}-\mathbf{P}\Vert }_{\text{F}}\le \epsilon \end{array}$$ (20)

In summary, $\mathbf{P}$ can be reconstructed by iteratively minimizing its weighted atomic norm and updating the weighting coefficient accordingly in each iteration. Hence, this proposed algorithm can be named IRANM. Like ANM, IRANM is also applicable to the NURMA. This is realized by changing ${\mathbf{P}}^{\mathrm{★}}-\mathbf{P}$ in equations (12), (13), (15), and (20) to ${\mathbf{P}}_{\mathrm{\Omega }}^{\mathrm{★}}-{\mathbf{P}}_{\mathrm{\Omega }}$ .

A simulated illustrative example

We reconstruct the source DOAs distribution assumed in Figure 2(a) by IRANM + MaPP. Here, MaPP estimates the number of sources as the number of eigenvalues of $T_b\left({\widehat{\mathbf{u}}}^k\right)$ that are greater than 10⁻⁶. Figure 4 presents the results when the standard URMA is utilized. Figure 4(a) shows that $T_b\left({\widehat{\mathbf{u}}}^1\right)$ has slowly decreased eigenvalues, whereas $T_b\left({\widehat{\mathbf{u}}}^2\right)$ and $T_b\left({\widehat{\mathbf{u}}}^3\right)$ have only 6 relatively large eigenvalues and their other eigenvalues are as low as 10⁻⁸. That agrees with the third property of ${\mathrm{M}}^{\kappa }\left(\mathbf{P}\right)$ in Section 3.1. Figure 4(b) shows that the reconstructed source DOAs distribution deviates from the true one seriously after one iteration. This can be attributed to the fact that on one hand, as demonstrated by Figure 2(a), the $T_b\left(\widehat{\mathbf{u}}\right)$ obtained by ANM (the first iteration of IRANM) does not contain the correct DOA information, on the other hand, as shown by Figure 4(a), there are 55 eigenvalues greater than 10⁻⁶ for $T_b\left({\widehat{\mathbf{u}}}^1\right)$ , which means that the estimated number of sources is far greater than the true one. Figure 4(c) and (d) show that the source DOAs distribution has been accurately reconstructed after two iterations. This is because, as shown by Figure 4(a), both $T_b\left({\widehat{\mathbf{u}}}^2\right)$ and $T_b\left({\widehat{\mathbf{u}}}^3\right)$ are low-rank and admit a unique Vandermonde decomposition, and MaPP accurately estimates the number of sources. Figure 4(e)–(g) present the weighting coefficients used in the first three iterations. The weighting coefficients utilized in the present iteration is got according to the $T_b\left(\widehat{\mathbf{u}}\right)$ and $\kappa$ determined in the previous iteration, and affords preference to the atoms around the sources identified in the previous iteration. Figure 5 presents the results when the NUMRA shown in Figure 3(a) is utilized, which exhibits the similar phenomena to Figure 4. These phenomena indicate that the IRANM-based 2D multiple-snapshot GFCB can remedy the drawback of the ANM-based one.

Figure 4.

Results of IRANM-based 2D multiple-snapshot grid-free GFCB when a standard URMA with 64 microphones is utilized. (a) Eigenvalues of $T_b\left({\widehat{\mathbf{u}}}^k\right)$ . Reconstructions (*) for sources (○) after (b) one iteration, (c) two iterations and (d) three iterations. Weighting coefficients utilized in (e) the first iteration, (f) the second iteration and (g) the third iteration.

Figure 5.

Results of IRANM-based 2D multiple-snapshot GFCB when a NURMA with 30 randomly selected microphones is utilized. (a) Eigenvalues of $T_b\left({\widehat{\mathbf{u}}}^k\right)$ . Reconstructions (*) for sources (○) after (b) one iteration, (c) two iterations and (d) three iterations. Weighting coefficients utilized in (e) the first iteration, (f) the second iteration and (g) the third iteration.

Performance comparisons

Simulations

We compare the reconstruction results of DOAs distribution based on ANM + MaPP and IRANM + MaPP for the sources assumed in Figure 2(a). In the former, MaPP estimates the number of sources as the number of eigenvalues of $T_b\left(\widehat{\mathbf{u}}\right)$ in 20 dB dynamic range. In the latter, IRANM stops the iteration when the relative change of $\widehat{\mathbf{P}}$ between two consecutive iterations $\left(\mathrm{i}\mathrm{.e}\mathrm{.,}\hspace{0.17em}{\Vert {\widehat{\mathbf{P}}}^k-{\widehat{\mathbf{P}}}^{k-1}\Vert }_{\text{F}}/{\Vert {\widehat{\mathbf{P}}}^{k-1}\Vert }_{\text{F}}\right)$ is not greater than 10⁻³ or the number of iterations reaches the preset maximum value (set to 20), and MaPP estimates the number of sources as the number of eigenvalues of $T_b\left(\widehat{\mathbf{u}}\right)$ that are greater than 10⁻⁶. Figure 6(a) and (b) presents the results. Change the frequency from 2000 Hz to 3000 Hz, 4000 Hz and 4900 Hz to conduct simulations. Correspondingly, the minimum separation among these sources $\left({\mathrm{\Delta }}_{\min }\right)$ changes from 0.025 to 0.038, 0.050 and 0.062. Figure 6(c)–(h) presents the results. It is worth mentioning that 4900 Hz is the upper limit frequency of the current array according to $f_{\max }=0.5c/\max \left\{\mathrm{\Delta }x,\mathrm{\Delta }y\right\}$ with $c$ being the speed of sound.^20–22 Obviously, affected by the drawback of failing to identify sources with a small separation, the ANM-based 2D multiple-snapshot GFCB (Figure 6(a), (c), (e) and (g)) reconstructs the source distribution accurately only at 4900 Hz. In contrast, the IRANM-based 2D multiple-snapshot GFCB (Figure 6(b), (d), (f) and (h)) reconstructs the source distribution accurately at all the four frequencies, enjoying the enhanced resolution. To obtain Figure 6, the URMA with 64 microphones is utilized. Figure 7 presents the results when using the NURMA with 30 randomly selected microphones, again exhibiting the advantage of the IRANM-based 2D multiple-snapshot GFCB.

Figure 6.

Reconstructions (*) of (a), (c), (e), (g) ANM- and (b), (d), (f), (h) IRANM-based 2D multiple-snapshot GFCB for sources (○) at (a), (b) 2000 Hz, (c), (d) 3000 Hz, (e), (f) 4000 Hz and (g), (h) 4900 Hz. The URMA with 64 microphones is utilized.

Figure 7.

Reconstructions (*) of (a), (c), (e), (g) ANM- and (b), (d), (f), (h) IRANM-based 2D multiple-snapshot GFCB for sources (○) at (a), (b) 2000 Hz, (c), (d) 3000 Hz, (e), (f) 4000 Hz and (g), (h) 4900 Hz. The NURMA with 30 microphones is utilized.

Referring to Ref. [23],we use ${\Vert \widehat{\mathbf{P}}-\mathbf{P}\Vert }_{\text{F}}/{\Vert \mathbf{P}\Vert }_{\text{F}}$ , ${\Vert \left[\widehat{\mathbf{\theta }},\widehat{\mathbf{\varphi }}\right]-\left[\mathbf{\theta },\mathbf{\varphi }\right]\Vert }_{\text{F}}/2I$ and ${\Vert {\widehat{\mathbf{s}}}_{\text{rms}}-{\mathbf{s}}_{\text{rms}}\Vert }_2/{\Vert {\mathbf{s}}_{\text{rms}}\Vert }_2$ to measure the pressure reconstruction accuracy, the DOA estimation accuracy and the strength quantification accuracy, respectively. Hereinto, $\mathbf{\theta }\in {\mathrm{R}}^I$ , $\mathbf{\varphi }\in {\mathrm{R}}^I$ and ${\mathbf{s}}_{\text{rms}}\in {\mathrm{R}}^I$ are the true elevation angle vector, the true azimuth angle vector and the true root mean square strength vector. $\widehat{\mathbf{\theta }}\in {\mathrm{R}}^I$ , $\widehat{\mathbf{\varphi }}\in {\mathrm{R}}^I$ , and ${\widehat{\mathbf{s}}}_{\text{rms}}\in {\mathrm{R}}^I$ are the estimated ones. Figure 8 shows the curves of these errors vs. ${\mathrm{\Delta }}_{\min }$ of the ANM- and the IRANM-based 2D multiple-snapshot GFCB. For each ${\mathrm{\Delta }}_{\min }$ , these errors are averaged over 50 Monte Carlo trials. In each trial, we randomly generate two sources with a distance ${\mathrm{\Delta }}_{\min }$ and a strength difference of less than 10 dB. The SNR is set as 20 dB. The URMA with 64 microphones is utilized to obtain Figure 8(a), (c) and (e). The NURMA with 30 randomly selected microphones is utilized to obtain Figure 8(b), (d) and (f). For reconstruction of $\mathbf{P}$ , Figure 8(a) and (b) show that the error of ANM under ${\mathrm{\Delta }}_{\min }<0.8/\sqrt{AB}$ is distinctly higher than that under ${\mathrm{\Delta }}_{\min }\ge 0.8/\sqrt{AB}$ , whereas the error of IRANM varies slightly across all ${\mathrm{\Delta }}_{\min }$ and is distinctly lower than the one of ANM. For DOA estimation, Figure 8(c) and (d) show that on one hand, the IRANM-based 2D multiple-snapshot GFCB has very low errors across all ${\mathrm{\Delta }}_{\min }$ ; on the other hand, compared to the former, the ANM-based 2D multiple-snapshot GFCB has distinctly higher errors under ${\mathrm{\Delta }}_{\min }<0.5/\sqrt{AB}$ , slightly higher errors under $0.5/\sqrt{AB}\le {\mathrm{\Delta }}_{\min }<0.8/\sqrt{AB}$ and almost the same errors under ${\mathrm{\Delta }}_{\min }\ge 0.8/\sqrt{AB}$ . For strength quantification, Figure 8(e) and (f) shows that the error of the ANM-based 2D multiple-snapshot GFCB is always higher than the one of the IRANM-based 2D multiple-snapshot GFCB, especially under small ${\mathrm{\Delta }}_{\min }$ . These phenomena demonstrate that, compared to the ANM-based, the proposed IRANM-based one has the better performance in denoising and spatial resolution. Therefore, even in the case of a small separated signal source, the proposed method can accurately reconstruct the microphone pressure induced by sources, estimate the DOAs, and quantify the strengths.

Figure 8.

Error curves. (a), (b) ${\Vert \widehat{\mathbf{P}}-\mathbf{P}\Vert }_{\text{F}}/{\Vert \mathbf{P}\Vert }_{\text{F}}$ , (c), (d) ${\Vert \left[\widehat{\mathbf{\theta }},\widehat{\mathbf{\varphi }}\right]-\left[\mathbf{\theta },\mathbf{\varphi }\right]\Vert }_{\text{F}}/2k$ and (e), (f) ${\Vert {\widehat{\mathbf{s}}}_{\text{rms}}-{\mathbf{s}}_{\text{rms}}\Vert }_2/{\Vert {\mathbf{s}}_{\text{rms}}\Vert }_2$ vs. ${\mathrm{\Delta }}_{\min }$ . (a), (c), (e) The URMA with 64 microphones and (b), (d), (f) the NURMA with 30 randomly selected microphones are utilized.

Experiments

As shown in Figure 9, in a semi-anechoic chamber, we use a URMA with 64 Brüel&Kjær Type 4958 microphones to conduct an experimental measurement on three small speakers to verify correctness of the simulation conclusion. The array parameters are as follows: $A=B=8$ and $\mathrm{\Delta }x=\mathrm{\Delta }y=0.035\text{ m}$ . The three loudspeakers are located at Cartesian coordinates (2.24, 0, 5) m, (−1.24, 0, 5) m and (−2.24, 0, 5) m from left to right. Correspondingly, their DOAs are (24.13°, 0°), (13.93°, 180°), and (24.13°, 180°). In addition, the Cartesian coordinates of three mirror sources caused by the reflection effect of the ground are (2.24, −2.2, 5) m, (−1.24, −2.2, 5) m, and (−2.24, −2.2, 5) m. Correspondingly, their DOAs are (32.13°, 315.52°), (26.80°, 240.59°), and (32.13°, 224.48°). All microphone pressure signals are acquired synchronously by Brüel&Kjær PULSE Type 3560D Data Acquisition System and then transferred to Brüel&Kjær BKConnect Software to obtain Fourier spectra. Eventually, these Fourier spectra are processed by ANM + MaPP and IRANM + MaPP to reconstruct the source DOAs distribution. MaPP utilizes the same method to estimate the number of sources and IRANM utilizes the same stopping condition for iteration as in Section 4.1.

Figure 9.

Experimental layout.

Figure 10 presents the results at 2000 Hz, 3000 Hz, 4000 Hz and 4900 Hz when all the 64 microphones are utilized. For these four frequencies, the minimum separation among these sources are $0.032\left(\text{0}\text{.}256/\sqrt{AB}\right)$ , $0.048\left(\text{0}\text{.}384/\sqrt{AB}\right)$ , $0.065\left(\text{0}\text{.}52/\sqrt{AB}\right)$ , and $0.079\left(\text{0}\text{.}63/\sqrt{AB}\right)$ , respectively. As shown by Figure 10(a), (c), (e), and (g), ANM + MaPP reconstructs the source DOAs distribution accurately only at 4000 Hz and 4900 Hz. In contrast, as shown by Figure 10(b), (d), (f), and (h), IRANM + MaPP reconstructs the source DOAs distribution accurately at all the four frequencies. Figure 11 presents the result when only 30 randomly selected microphones are utilized, which exhibits the same phenomena as in Figure 10. These phenomena demonstrate that the IRANM-based 2D multiple-snapshot GFCB has the enhanced resolution compared to the ANM-based one. The experimental conclusion is consistent with the simulation one, indicating that the conclusion is correct.

Figure 10.

Reconstructions (*) of (a), (c), (e), (g) ANM and (b), (d), (f), (h) IRANM-based 2D multiple-snapshot GFCB for loudspeaker sources and their mirror sources (○) at (a), (b) 2000 Hz, (c), (d) 3000 Hz, (e), (f) 4000 Hz and (g), (h) 4900 Hz. The URMA with 64 microphones is utilized. The symbols ○ indicate the true DOAs but do not contain the strength information, so do in Figure 11.

Figure 11.

Reconstructions (*) of (a), (c), (e), (g) ANM and (b), (d), (f), (h) IRANM-based 2D multiple-snapshot GFCB for loudspeaker sources and their mirror sources (○) at (a), (b) 2000 Hz, (c), (d) 3000 Hz, (e), (f) 4000 Hz and (g), (h) 4900 Hz. The NURMA with 30 randomly selected microphones are utilized.

Conclusions

The existing ANM-based 2D multiple-snapshot GFCB fails to work well for sources with a small separation. After analyzing the cause of this drawback, we develop an IRANM-based 2D multiple-snapshot GFCB to remedy it. Through simulations and experiments, some interesting and valuable conclusions have been drawn. First, the developed IRANM-based 2D multiple-snapshot GFCB can accurately estimate the DOAs and quantify the strengths of the small separation sources, successfully remedying the drawback of the ANM-based one and enhancing the resolution. Second, the IRANM-based 2D multiple-snapshot GFCB has stronger denoising ability and higher identification accuracy than the ANM-based one. Finally, the above conclusions hold true not only for a standard URMA, but also for a NURMA consisting of a small number of randomly selected microphones.