The Cahn–Hilliard equations are a versatile model for describing the evolution of complex morphologies.
In this paper we present a computational pipeline for the numerical solution of a ternary phase-field model for describing the nanomorphology of donor–acceptor semiconductor blends used in organic photovoltaic devices.
The model consists of two coupled fourth-order partial differential equations that are discretized using a finite element approach.
In order to solve the resulting large-scale linear systems efficiently, we propose a preconditioning strategy that is based on efficient approximations of the Schur-complement of a saddle point system.
We show that this approach performs robustly with respect to variations in the discretization parameters.
Finally, we outline that the computed morphologies can be used for the computation of charge generation, recombination, and transport in organic solar cells.

preconditioning, phase–field models, organic solar cells, Cahn–Hilliard, finite element analysis

Missing 'biblatex' package
The bibliography requires the 'biblatex' package.

University Library Heidelberg

journaltitleArchive of Numerical Software
titleThe FEniCS project version 1.5

booktitleAutomated Solution of Differential Equations by the Finite Element Method
seriesLecture Notes in Computational Science and Engineering
titleUFC: a finite element code generation interface

journaltitleAssociation for Computing Machinery. Transactions on Mathematical Software
titleUnified form language: a domain-specific language for weak formulations and partial differential equations
pagesArt. 937

Springer Science
Business Media LLC

journaltitleNumerical Algorithms
titleComparison of preconditioned Krylov subspace iteration methods for PDE-constrained optimization problems

Springer Science
Business Media LLC

titleBlock preconditioners for elliptic PDE-constrained optimization problems

journaltitleActa Numerica
titleNumerical solution of saddle point problems

Technische Universitat Chemnitz

journaltitleGAMM Archive for Students
titleModeling the morphology evolution of organic solar cells

Otto-von-Guericke-Universität Magdeburg

titleFast Iterative Solvers for Cahn-Hilliard Problems

journaltitleSIAM Journal on Scientific Computing
titlePreconditioning for vector-valued Cahn-Hilliard equations

journaltitleJournal of Computational Physics
titleFast solution of Cahn-Hilliard variational inequalities using implicit time discretization and finite elements

Elsevier BV

journaltitleComputers & Mathematics with Applications
titleEfficient numerical solution of discrete multi-component Cahn–Hilliard systems

Walter de Gruyter GmbH

journaltitleComputational Methods in Applied Mathematics
titleEfficient preconditioners for large scale binary Cahn-Hilliard models

Springer New York

titleThe Mathematical Theory of Finite Element Methods

AIP Publishing

journaltitleThe Journal of Chemical Physics
titleFree energy of a nonuniform system. I. Interfacial free energy

Association for Computing Machinery (ACM)

journaltitleACM Transactions on Mathematical Software
titleAlgorithm 832: UMFPACK V4.3 — an unsymmetric-pattern multifrontal method

Oxford University Press

titleDirect Methods for Sparse Matrices

Institute for Mathematics
Its Applications, University of Minnesota

titleA generalised diffusion equation for phase separation of a multi-component mixture with interfacial free energy
typeIMA Preprints Series

Oxford University Press

seriesNumerical Mathematics and Scientific Computation
titleFinite Elements and Fast Iterative Solvers: with Applications in Incompressible Fluid Dynamics

Institute of Electrical
Electronics Engineers (IEEE)

journaltitleComputing in Science & Engineering
titleAn introduction to algebraic multigrid

Cornell University Press

titlePrinciples of Polymer Chemistry

journaltitleAdvanced Electronic Materials
titleIn situ visualization and quantification of electrical self-heating in conjugated polymer diodes using Raman spectroscopy

journaltitleSIAM Journal on Scientific Computing
titleA note on preconditioning for indefinite linear systems
pages19691972 (electronic)

noteRelease 4.0
titlePyAMG: Algebraic Multigrid Solvers in Python v4.0

journaltitleSIAM Journal on Numerical Analysis
titleSolution of sparse indefinite systems of linear equations

journaltitleSIAM Journal on Matrix Analysis and Applications
titleRegularization-robust preconditioners for time-dependent PDE-constrained optimization problems

journaltitleNumerical Linear Algebra with Applications
titleA new approximation of the Schur complement in preconditioners for PDE-constrained optimization

titleComputational Materials Science

titleIterative Methods for Sparse Linear Systems

Springer Science
Business Media LLC

journaltitleNature Methods
titleSciPy 1.0: fundamental algorithms for scientific computing in Python

journaltitleElectronic Transactions on Numerical Analysis
titleChebyshev semi-iteration in preconditioning for problems including the mass matrix

Elsevier BV

journaltitleComputational Materials Science
titleModeling morphology evolution during solvent-based fabrication of organic solar cells

Dover Publications, Inc., Mineola, NY

noteReprint of the 1983 original
titleFinite Elements and Approximation

journaltitleSIAM Journal on Matrix Analysis and Applications
titleNonstandard norms and robust estimates for saddle point problems


# Optimal Real-Weighted Beamforming with Application to Linear and Spherical
Arrays

V. Tourbabin*, , M. Agmon, B. Rafaely, and J. Tabrikian  The authors are with
the Department of Electrical and Computer Engineering, Ben-Gurion University
of the Negev, Be’er-Sheva 84105, Israel (email:
<EMAIL_ADDRESS>
###### Abstract

One of the uses of sensor arrays is for spatial filtering or beamforming.
Current digital signal processing methods facilitate complex-weighted
beamforming, providing flexibility in array design. Previous studies proposed
the use of real-valued beamforming weights, which although reduce flexibility
in design, may provide a range of benefits, e.g. simplified beamformer
implementation or efficient beamforming algorithms. This paper presents a new
method for the design of arrays with real-valued weights, that achieve maximum
directivity, providing closed-form solution to array weights. The method is
studied for linear and spherical arrays, where it is shown that rigid
spherical arrays are particularly suitable for real-weight designs as they do
not suffer from grating lobes, a dominant feature in linear arrays with real
weights. A simulation study is presented for linear and spherical arrays,
along with an experimental investigation, validating the theoretical
developments.

###### Index Terms:

Beamforming, microphone arrays, linear arrays, spherical arrays, array signal
processing.

EDICS Category: AUD-LMAP

##  I Introduction

Signal processing for microphone arrays of various configurations, including
linear, planar and spherical, is widely studied and reported in the literature
[1, 2, 3, 4, 5]. One of the primary goals of microphone arrays is beamforming,
or spatial filtering[6]. It is one of the simplest methods for discriminating
between different signals based on the physical location, distance, or
direction of the sources. An effective way to achieve spatial filtering is
using arrays of discrete sensors. When wired with independent conditioning
electronics, discrete array sensors are flexible and are therefore very useful
in time varying environments.

A variety of beamforming techniques is available, including data-independent,
statistically optimal and adaptive approaches. Well studied examples are least
squares beampattern synthesis[1], MVDR or general LCMV beamformers[6] and the
sidelobe canceller techniques[7]. Although these beamformers are flexible and
facilitate the design of arrays with desired spatial and frequency
characteristics, they require complex array weights, i.e. with magnitude and
phase specified at each frequency.

In order to achieve high spatial resolution, there is a need for a large
number of channels resulting in very complicated beamformer structures.
Usually, it is possible to simplify the structure by finding a compromise
between performance and complexity[8, 9]. In this paper, we analyze the
possibility to perform the processing using real-valued beamformer weights
without the need for phase adjustments [10]. This approach appears in the
adaptive antenna array processing literature. See for example the SOAC [11,
12, 13], RAMONA [14, 15, 16] and additional methods [17, 18, 19]. The main
advantage of this approach is that, although complex weights offer greater
control over the array response, using real weights requires only real
arithmetic operations for the beamforming procedure. Therefore, the
computational complexity can be considerably reduced. In [11] a method for
suppressing jammers by controlling the current amplitude only for a linear
antenna array is described. It is shown that, although the total number of
jammers is doubled with complex weights, using real weights results in a much
faster computing speed and considerable hardware savings.

In addition to efficient realization, real weights can also improve some
aspects of array processing. With real weights, an improvement of the Gram-
Schmidt orthogonal projection nulling algorithm [20] is achieved in [14],
where it is shown that the algorithm proposed there performs four times faster
than the original one. Another paper by Choi [17] presents an adaptive
beamforming technique with real estimated weights. In this paper a statistical
measurement model is assumed, which includes the desired signal, noise, and
several interference sources. It is shown that the imaginary part of the
measurements’ covariance matrix is sufficient to span the interferer subspace.
Thus, by using orthogonal subspace projection, a real-weighted beamformer can
be calculated. In [19], a data-independent real-weighted beamforming technique
for a rectangular sensor array is described. A computationally efficient
design algorithm is presented where the real array weights are obtained by
constraining the beampattern to a desired response at several symmetrically-
chosen points in frequency-angle space.

In addition, the optimal real-weighted beamforming techniques can be utilized
for designing beamformers implemented with continuous sensors like Electro-
Mechanical Film (EMFI) [21, 22, 23, 24], which have potential advantages over
the traditional electret microphones, e.g. large aperture with a single
sensor, although such a directional sensor may not be electronically
steerable.

Although the above-mentioned publications presented real-valued beamforming
techniques, there is a lack of system performance analysis based on common
data-independent performance measures like directivity or sensitivity. Also,
no beamforming technique has been described that achieves optimality with
respect to these measures. Furthermore, although it was evident from [19] that
real weights impose an undesired beampattern symmetry, yielding an additional
parasitic main lobe and thus limiting the use of real-weighted beamformers in
practice, this issue has not been systematically addressed in previous
publications.

This paper presents novel formulation and analysis of real-weighted
beamforming, offering the following contributions:

  1. 1.

Mathematical formulations are developed for optimal beamformers in terms of
maximum beampattern directivity and minimum beamformer sensitivity subject to
the constraint of using real weights. A closed-form solution for the optimal
real weights is derived for general array geometry.

  2. 2.

The problem of beampattern symmetry induced by real weights is formulated,
showing clearly the existence of an additional parasitic main lobe when used
with open array geometries. This symmetry clearly degrades the flexibility
and, in fact, the overall performance of open arrays with real weights.

  3. 3.

In order to overcome the symmetry imposed by real weights, it is shown that
the undesired parasitic main lobe does not exist for rigid spherical arrays,
making these arrays ideal for real-weighted beamforming. The narrowband
performance of spherical arrays that use real weights is shown to be
comparable to the performance of spherical arrays that use complex-valued
weights. It is also shown that the proposed method facilitates the design of
spherical arrays with electronic steering capabilities.

The optimal array processing problem described in this paper for the spherical
array geometry can be regarded as a special case of the FIR based realization
[25] where only single weight filters are considered. This implies that the
analytic solutions derived are comparable to the numerical solutions in the
single-weight FIR realization.

The paper starts with a brief review of general beamforming and common
performance measures. Next, the optimization problem is formulated for a
general array geometry and the maximum directivity solution for real weights
is derived. Simulation examples for linear and spherical arrays are then
presented and the resulting beampattern properties are discussed. The results
of directional sound field analysis [26] using the proposed optimal
beampatterns, are presented for experimental data measured with spherical
microphone array. Conclusions and future work follow.

##  II Beamforming

###  II-A Array Beampattern

Consider an array consisting of $M$ sensors with an arbitrary geometry.
Assuming a free field, the $q$’th sensor’s output
$y_{q}\left(\omega,\Omega,t\right)$ at a time $t$ as a response to a unit
amplitude plane wave having frequency $\omega$ and arriving at direction
$\Omega=\left(\theta,\phi\right)$, is given by:

$y_{q}\left(\omega,\Omega,t\right)=v_{q}\left(\omega,\Omega\right)e^{i\omega
t},$ (1)

where, $v_{q}(\omega,\Omega)$ is the gain and phase shift associated with
location of the $q$’th sensor. The overall array output is given by:

$\displaystyle y\left(\omega,\Omega,t\right)$
$\displaystyle=\sum\limits_{q=1}^{M}{w_{q}\left(\omega\right)y_{q}\left(\omega,\Omega,t\right)}$
$\displaystyle=\sum\limits_{q=1}^{M}{w_{q}\left(\omega\right)v_{q}\left(\omega,\Omega\right)e^{i\omega
t}},$ (2)

where, $w_{q}\left(\omega\right)$ is the weight applied by the beamformer to
the $q$’th sensor at frequency $\omega$. The dependence on frequency will be
omitted for notation simplicity. Now, by dropping the time dependence (as we
are interested only in the response) and rewriting (2) in matrix form, we get
the beampattern of the array:

$B(\Omega)=\mathbf{w}^{T}\mathbf{v}(\Omega),$ (3)

where $\mathbf{w}=\left[w_{1}\,w_{2}\,...\,w_{M}\right]^{T}$ is the weights
vector and
$\mathbf{v}(\Omega)=\left[v_{1}\left(\Omega\right)\,v_{2}\left(\Omega\right)\,...\,v_{M}\left(\Omega\right)\right]^{T}$
is the manifold or steering vector of the array. In general, the weights
vector of an array incorporates its beamformer characteristics, while the
manifold vector describes its geometry.

Beamforming using real weights does not require phase shift implementation.
Thus, from the mathematical point of view, we will simply assume that the
weights vector is real-valued. Later in this text we will discuss linear and
spherical array geometries as examples of optimization techniques presented in
Section III.

###  II-B Performance Measures

There are a number of performance measures for assessing the capabilities of a
sensor array defined in the literature[1]. The performance measures considered
in this paper are: (i) Directivity, (ii) Sensitivity.

Directivity is defined as the ratio between the power of the output in the
look direction to the average output power:

$\displaystyle D$
$\displaystyle=\frac{{|B(\Omega_{l})|}^{2}}{\frac{1}{4\pi}\int\limits_{0}^{\pi}{\int\limits_{0}^{2\pi}{{|B(\Omega)|}^{2}\sin\theta
d\phi d\theta}}}$
$\displaystyle=\frac{\mathbf{w}^{T}\mathbf{bb}^{H}\mathbf{w}^{*}}{\mathbf{w}^{T}\left(\frac{1}{4\pi}\int\limits_{0}^{\pi}{\int\limits_{0}^{2\pi}{{\mathbf{v}(\Omega)\mathbf{v}^{H}(\Omega)\sin\theta
d\phi d\theta}}}\right)\mathbf{w}^{*}},$ (4)

where $\mathbf{b}=\mathbf{v}(\Omega_{l})$ is the manifold vector in the look
direction $\Omega_{l}$. We define $\mathbf{C}$ as the following Hermitian
matrix:

$\mathbf{C}\triangleq\frac{1}{4\pi}\int\limits_{0}^{\pi}{\int\limits_{0}^{2\pi}{{\mathbf{v}(\Omega)\mathbf{v}^{H}(\Omega)\sin\theta
d\phi d\theta}}}.$ (5)

Finally, directivity can be written as:

$D=\frac{\mathbf{w}^{T}\mathbf{bb}^{H}\mathbf{w}^{*}}{\mathbf{w}^{T}\mathbf{C}\mathbf{w}^{*}}.$
(6)

For the case when the weights vector is real-valued we have:

$D=\frac{\mathbf{w}^{T}\mathbf{bb}^{H}\mathbf{w}}{\mathbf{w}^{T}\mathbf{C}\mathbf{w}}.$
(7)

The second array performance measure discussed here is sensitivity, which is a
measure of array robustness to processor gain errors, sensor mismatch, and
imprecise positioning of the sensors. The sensitivity is proportional to the
square of the weights vector norm [1]:

$T_{se}={\|\mathbf{w}\|}^{2}=\mathbf{w}^{H}\mathbf{w}.$ (8)

Once again, for the real-valued weights vector we simply obtain:

$T_{se}={\|\mathbf{w}\|}^{2}=\mathbf{w}^{T}\mathbf{w}.$ (9)

Needless to say, that it is always highly desirable to construct an array with
the highest possible directivity and lowest possible sensitivity to errors. In
the next section, we present the maximum directivity optimization problem and
derive the maximum directivity real-valued beamformer subject to
distortionless response constraint. We also discuss a lower bound on the
sensitivity resulting from the distortionless response constraint and present
a technique for finding an optimal directivity beamformer subjected to a
constraint on the maximum sensitivity value.

##  III Optimal Beamformers

###  III-A Complex-valued maximum directivity beamformer

In this subsection, we present the formulation of the general optimization
problem for the maximum directivity beamformer and its well known complex-
valued solution. In the next subsection we will derive the real-valued
solution for the maximum directivity optimization problem.

In order to maximize the directivity given by (6) subject to the
distortionless constraint, i.e.
$\mathbf{w}^{T}\mathbf{bb}^{H}\mathbf{w}^{*}=1$, it is sufficient to minimize
$\mathbf{w}^{T}\mathbf{Cw}^{*}$. Note, that the distortionless response
constraint used here is applied only to the absolute value of the response,
leaving a degree of freedom in choosing the phase. The optimization problem
can be stated as follows:

$\mathbf{w}_{c}=\arg\min\limits_{\mathbf{w}}\left(\mathbf{w}^{T}\mathbf{Cw}^{*}\right)\,\,\,\,\,\,\,\,s.t.\,\,\,\,\,\,\,\mathbf{w}^{T}\mathbf{bb}^{H}\mathbf{w}^{*}=1$
(10)

The solution to this problem is known in the literature[1, 27], and is widely
used to derive optimum beamformers such as MVDR or MPDR:

$\mathbf{w}_{c}=\left(\frac{\mathbf{C}^{-1}\mathbf{b}}{\mathbf{b}^{H}\mathbf{C}^{-1}\mathbf{b}}\right)^{*}.$
(11)

However, this solution, in general, yields complex weights. A derivation of a
real-valued solution, not available in literature, is presented in the
following subsection.

###  III-B Real-valued maximum directivity beamformer

In this section, the real-valued solution to the maximum directivity
beamformer optimization problem is derived. Assuming $\mathbf{w}$ is a real-
valued vector and $\mathbf{C}$ is a Hermitian matrix, the problem can be
formulated as follows:

$\mathbf{w}_{r}=\arg\min\limits_{\mathbf{w}}\left(\mathbf{w}^{T}\mathbf{Cw}\right)\,\,\,\,\,\,\,\,s.t.\,\,\,\,\,\,\,\mathbf{w}^{T}\mathbf{bb}^{H}\mathbf{w}=1$
(12)

Using the Lagrange multipliers method, the Lagrangian is given by:

$J(\mathbf{w},\lambda)=\mathbf{w}^{T}\mathbf{Cw}+\lambda\left(\mathbf{w}^{T}\mathbf{bb}^{H}\mathbf{w}-1\right),$
(13)

where $\lambda$ is the Lagrange multiplier. Taking the gradient with respect
to $\mathbf{w}$ (assuming it is strictly a real-valued vector) and setting it
to zero in order to find the minimum, we get:

$\mathbf{w}^{T}Re\left\\{\mathbf{C}\right\\}+\lambda\mathbf{w}^{T}Re\left\\{\mathbf{bb}^{H}\right\\}=0.$
(14)

We denote $\mathbf{\tilde{C}}=Re\left\\{\mathbf{C}\right\\}$. By right-
multiplication of (14) with $\mathbf{\tilde{C}}^{-1}\mathbf{b}$, one obtains:

$\mathbf{w}^{T}\mathbf{b}+\lambda\mathbf{w}^{T}Re\left\\{\mathbf{bb}^{H}\right\\}\mathbf{\tilde{C}}^{-1}\mathbf{b}=0.$
(15)

Now the constraint is used to evaluate $\lambda$. According to the constraint
$\mathbf{w}^{T}\mathbf{b}=e^{j\varphi}$, where $\varphi$ is an arbitrary real-
valued phase. Thus (15) can be rewritten as:

$e^{j\varphi}+\lambda
Re\left\\{\mathbf{b}^{H}e^{j\varphi}\right\\}\mathbf{\tilde{C}}^{-1}\mathbf{b}=0.$
(16)

Now, $\lambda$ can be evaluated as

$\lambda=-\frac{1}{Re\left\\{\mathbf{b}e^{-j\varphi}\right\\}^{T}{\mathbf{\tilde{C}}}^{-1}\mathbf{b}e^{-j\varphi}}.$
(17)

By substituting (17) into (14) and simplifying, we get:

$\mathbf{w}_{opt}=\frac{\mathbf{\tilde{C}}^{-1}Re\left\\{\mathbf{b}e^{-j\varphi}\right\\}}{\left(\mathbf{b}e^{-j\varphi}\right)^{T}\mathbf{\tilde{C}}^{-1}Re\left\\{\mathbf{b}e^{-j\varphi}\right\\}}$
(18)

Finally, $\varphi$ can be found by forcing the imaginary part of (18) to zero.
Since the numerator is real-valued, it is sufficient to set the imaginary part
of the denominator to zero:

$\left[\mathbf{b}^{T}e^{-j\varphi}-\left(\mathbf{b}^{T}e^{-j\varphi}\right)^{*}\right]\mathbf{\tilde{C}}^{-1}Re\left\\{\mathbf{b}e^{-j\varphi}\right\\}=0.$
(19)

By expanding and simplifying (19), one obtains:

$\mathbf{b}^{T}\mathbf{\tilde{C}}^{-1}\mathbf{b}e^{-j2\varphi}=\mathbf{b}^{H}\mathbf{\tilde{C}}^{-1}\mathbf{b}^{*}e^{j2\varphi}.$
(20)

Solving (20) for $\varphi$, we get:

$\varphi=\frac{1}{2}\angle\left(\mathbf{b}^{T}\mathbf{\tilde{C}}^{-1}\mathbf{b}\right).$
(21)

Let $\mathbf{c}\triangleq Re\left\\{\mathbf{b}e^{-j\varphi}\right\\}$. Then
optimal beamformer is given by:

$\mathbf{w}_{r}=\frac{\mathbf{\tilde{C}}^{-1}\mathbf{c}}{\mathbf{c}^{T}\mathbf{\tilde{C}}^{-1}\mathbf{c}}.$
(22)

It should be emphasized here that the matrix $\mathbf{C}$ defined in (5)
includes the $\sin\theta$ factor. This factor can be treated as the cost
function of the optimization problem. As we show later, other cost functions
can be used as well, resulting in beneficial beampattern characteristics at
the expense of degraded directivity.

###  III-C Real-valued maximum directivity beamformer with bounded sensitivity

In order to overcome possible errors and instabilities, a more robust
beamforming technique is presented here. The following beamforming technique,
which is based on the maximum directivity beamformer presented above, aims to
find an optimal weights vector $\mathbf{w}$ that maximizes the directivity of
a given array subject to bounded sensitivity, i.e. $T_{se}\leq T_{0}$. If the
solution in (22) satisfies this constraint, it should be used directly.
Otherwise, it is assumed that the constraint is active and is therefore
satisfied with an equality, and the optimization problem with bounded
sensitivity is formulated as follows:

$\mathbf{w}_{opt}=\arg\min\limits_{\mathbf{w}}\left(\mathbf{w}^{T}\mathbf{Cw}\right)\,\,\,\,s.t.\,\,\,\,\left\\{\begin{array}[]{l}\mathbf{w}^{T}\mathbf{bb}^{H}\mathbf{w}=1\\\
\mathbf{w}^{T}\mathbf{w}=T_{0}\end{array}\right.$ (23)

Using the Lagrange multipliers method, we write the Lagrangian as

$\displaystyle J(\mathbf{w},\lambda)=\mathbf{w}^{T}\mathbf{Cw}$
$\displaystyle+\lambda\left(\mathbf{w}^{T}\mathbf{bb}^{H}\mathbf{w}-1\right)$
$\displaystyle+\beta\left(\mathbf{w}^{T}\mathbf{w}-T_{0}\right).$ (24)

By denoting $\mathbf{C}_{d}=\mathbf{C}+\beta\mathbf{I}$, where $\mathbf{I}$ is
the identity matrix, the Lagrangian can be written as:

$J(\mathbf{w},\lambda)=\mathbf{w}^{T}\mathbf{C}_{d}\mathbf{w}+\lambda\left(\mathbf{w}^{T}\mathbf{bb}^{H}\mathbf{w}-1\right)-\beta
T_{0},$ (25)

Now, assuming $\beta$ is real-valued, one can proceed as in Section III-B. The
optimal beamformer is given by:

$\mathbf{w}_{opt}=\frac{\mathbf{\tilde{C}}_{d}^{-1}\mathbf{c}}{\mathbf{c}^{T}\mathbf{\tilde{C}}_{d}^{-1}\mathbf{c}},$
(26)

where, $\mathbf{\tilde{C}}_{d}=Re\left\\{\mathbf{C}+\beta\mathbf{I}\right\\}$
and
$\mathbf{c}=Re\left\\{\mathbf{b}e^{-\frac{j}{2}\angle\left(\mathbf{b}^{T}\mathbf{\tilde{C}}_{d}\mathbf{b}\right)}\right\\}$.
The real-valued parameter $\beta$ can be found numerically by imposing the
constraint $\mathbf{w}^{T}\mathbf{w}=T_{0}$.

It can be seen that the solution here is very similar to the maximum
directivity beamformer. The only difference is the diagonal loading that
appears in the matrix $\mathbf{\tilde{C}}_{d}$ which is expected to improve
the robustness of the solution.

###  III-D Lower bounds on sensitivity

It is known that distortionless response constraint introduces a lower bound
on the array’s sensitivity[1]. In the general case of complex-valued weights,
the problem of minimum array sensitivity is mathematically identical to the
problem in (10) with the matrix $\mathbf{C}$ replaced by the identity matrix.
Thus, the minimum sensitivity beamformer is given by:

$\mathbf{w}_{MS}^{c}=\frac{\mathbf{b}^{*}}{\mathbf{b}^{H}\mathbf{b}},$ (27)

where the superscript $(\cdot)^{c}$ denotes the general case of complex-valued
weights. Moreover, the minimum achievable sensitivity with complex weights is:

$T_{\min}^{c}={\mathbf{w}_{MS}^{c}}^{H}{\mathbf{w}_{MS}^{c}}=\frac{1}{\mathbf{b}^{H}\mathbf{b}}.$
(28)

This is a well known result that reduces to $T_{\min}^{c}=\frac{1}{M}$ for
open arrays, i.e. arrays with steering vectors consisting of phase shifting
exponentials only.

Note that the lower bound here depends strictly on $\mathbf{b}^{H}\mathbf{b}$,
which is the only non-zero eigenvalue of the matrix $\mathbf{bb}^{H}$. We are
going to exploit this observation later in order to establish a relationship
between the lower bounds on sensitivity in complex and real-valued cases.

When the weights are constrained to be real-valued, the minimum sensitivity
beamformer problem is, again, identical to (12) with $\mathbf{C}$ replaced by
the identity matrix. The solution to this problem is given by (22):

$\mathbf{w}_{MS}^{r}=\frac{\mathbf{c}}{\mathbf{c}^{T}\mathbf{c}},$ (29)

where
$\mathbf{c}=Re\\{\mathbf{b}e^{-\frac{j}{2}\angle{\mathbf{b}^{T}\mathbf{b}}}\\}$.
Thus, the lower bound is given by:

$T_{min}^{r}={\mathbf{w}_{MS}^{r}}^{T}{\mathbf{w}_{MS}^{r}}=\frac{1}{\mathbf{c}^{T}\mathbf{c}}.$
(30)

In order to get a more informative expression for $T_{min}^{r}$, we solve the
minimum sensitivity problem directly, without using the solution in (22).
Consider the following minimization problem:

$\min\limits_{\mathbf{w}}\left(\mathbf{w}^{T}\mathbf{w}\right)\,\,\,\,\,\,\,\,s.t.\,\,\,\,\,\,\mathbf{w}^{T}\mathbf{bb}^{H}\mathbf{w}=1,$
(31)

where $T_{se}=\mathbf{w}^{T}\mathbf{w}$ is the sensitivity which we wish to
minimize subject to the distortionless response constraint. The Lagrangian
will be:

$J(\mathbf{w},\beta)=\mathbf{w}^{T}\mathbf{w}-\beta\left(\mathbf{w}^{T}\mathbf{bb}^{H}\mathbf{w}-1\right).$
(32)

By taking the derivative with respect to $\mathbf{w}$, setting it to zero and
rearranging, we get:

$Re\\{\mathbf{bb}^{H}\\}\mathbf{w}=\frac{1}{\beta}\mathbf{w}.$ (33)

Equation (33) holds for any pair of the scalar $\frac{1}{\beta}$ and the
vector $\mathbf{w}$, which are an eigenvalue and eigenvector pair of the
matrix $Re\\{\mathbf{bb}^{H}\\}$.

Note that $\mathbf{w}$ is a real-valued vector, hence it can be inserted
inside the $Re\\{\cdot\\}$ operator. By multiplying both sides of (33) from
the left by $\mathbf{w}^{T}$, and using the constraint in (31), we get:

$\mathbf{w}^{T}\mathbf{w}=\frac{1}{\gamma_{i}},$ (34)

where $\gamma_{i}$ is the $i^{th}$ eigenvalue of $Re\\{\mathbf{bb}^{H}\\}$.
Equation (34) implies that for each $\mathbf{w}$, which is a solution of (33),
there exists an appropriate eigenvalue of $Re\\{\mathbf{bb}^{H}\\}$ such that
(34) holds. It should be noted that the matrix $\mathbf{bb}^{H}$ has only one
non zero eigenvalue, but by taking the real part, one can double the number of
linearly independent rows. Thus, the minimum sensitivity is given by:

$T_{\min}^{r}=\min\limits_{\mathbf{w}}\left\\{\mathbf{w}^{T}\mathbf{w}\right\\}=\frac{1}{\gamma_{\max}},$
(35)

where, $\gamma_{\max}$ is the largest eigenvalue of matrix
$Re\\{\mathbf{bb}^{H}\\}$. Note that the bound is achievable by the
corresponding eigenvector of $Re\\{\mathbf{bb}^{H}\\}$, which in fact, is
given by (29).

In order to compare between $T_{\min}^{c}$ given in (28) and $T_{\min}^{r}$,
we realize that:

$\displaystyle\sum\limits_{i}{\gamma_{i}}$
$\displaystyle=\mathrm{tr}\left(Re\\{\mathbf{bb}^{H}\\}\right)=\mathrm{tr}(\mathbf{bb}^{H})$
$\displaystyle=\mathbf{b}^{H}\mathbf{b}=\mu,$ (36)

where $\mu$ is the non-zero eigenvalue of $\mathbf{bb}^{H}$. Hence,
$\gamma_{\max}=\max\limits_{i}\\{\gamma_{i}\\}\leq\mu$ and

$T_{\min}^{r}\geq T_{\min}^{c},$ (37)

meaning that for a given array geometry and look direction, the real-valued
beamformer can be at most as robust as the complex-valued one.

##  IV Linear Arrays

Here we demonstrate the maximum directivity beampattern that results from the
theory developed above for a linear array with uniformly spaced sensors. In
this case, the manifold vector is given by:

$\mathbf{v}(\psi)=\left[1\,\,e^{j\psi}\,\,e^{j2\psi}\,\,.\,.\,.\,\,e^{j(M-1)\psi}\right]^{T},$
(38)

where $M$ is the number of sensors and $\psi=\frac{2\pi d}{\lambda}\cos\theta$
incorporates the spatial dependence indicating the phase difference between
the signals in adjacent sensors as a function of plane wave arrival direction
$\theta$. The angle $\theta$ is measured from the endfire side, and the origin
of the coordinate system is positioned at the edge of the array. In addition,
$\lambda$ is the wavelength of the impinging plane wave and $d$ is the spacing
between the sensors.

The beampattern is given by:

$B(\psi)=\mathbf{w}^{T}\mathbf{v}(\psi)=\sum\limits_{n=0}^{M-1}{w_{n}\cdot
e^{i\psi n}}.$ (39)

This form resembles the discrete-time Fourier transform (DTFT) implying that
for real-valued weights vector $\mathbf{w}$, the beampattern will be symmetric
in $\psi$ and thus in $\theta$. As a result, in addition to the main lobe at
the desired look direction $\theta_{l}$, the beampattern will have a parasitic
main lobe at $\pi-\theta_{l}$. The same is evident from the analysis performed
in [19], although it is not stated explicitly.

The beampattern of general three-dimensional array with real-valued weights is
symmetric at least under the reversal of propagation direction. In order to
see this consider an array consisting of $M$ microphones located at
$\\{\mathbf{r}_{q}\\}_{q=1}^{M}$, where $\mathbf{r}_{q}$ is the position
vector of $q^{th}$ microphone relative to a chosen coordinate system. For an
open array, i.e. an array for which only the phase information is available,
the response of $q^{th}$ sensor to a plane wave with the vawevector
$\mathbf{k}_{0}$ is given by $e^{-j\mathbf{k}_{0}^{T}\mathbf{r}_{q}}$. Hence,
the beampattern is given by:

$B(\mathbf{k}_{0})=\sum\limits_{q=1}^{M}{w_{q}e^{-j\mathbf{k}_{0}^{T}\mathbf{r}_{q}}},$
(40)

where $\\{w_{q}\\}$ are the real-valued weights. For the reversed propagation
direction, i.e. $-\mathbf{k}_{0}$, we have:

$\displaystyle B(-\mathbf{k}_{0})$
$\displaystyle=\sum\limits_{q=1}^{M}{w_{q}e^{j\mathbf{k}_{0}^{T}\mathbf{r}_{q}}}$
$\displaystyle=\left(\sum\limits_{q=1}^{M}{w_{q}e^{-j\mathbf{k}_{0}^{T}\mathbf{r}_{q}}}\right)^{*}$
$\displaystyle=B^{*}(\mathbf{k}_{0}).$ (41)

Note, that for arrays with certain geometries, additional symmetries can arise
as was shown above for the linear array.

Fig. 1 compares two maximum directivity beampatterns: the beampattern with
real-valued weights derived above and the general beampattern with complex-
valued weights [1]. The beampatterns are presented at frequency $f=1715$ Hz
for a uniform linear array consisting of $25$ sensors spaced at $10$ cm
intervals (resulting in total array aperture of $2.4$ m). The frequency is
chosen such that $\lambda=2d$ considering an acoustic medium with sound speed
of $343$ m/s.

The array is steered to $\theta=\frac{\pi}{4}$ resulting in expected symmetry
and a parasitic main lobe in the reciprocal direction $\theta=\frac{3\pi}{4}$
for the beamformer with real-valued weights.

Figure 1: Comparison between two maximum directivity beampatterns - the real-
weighted and the general complex-weighted, for a linear array with parameters:
$d=10$ cm, $M=25$, at $f=1715$ Hz.

Ignoring the parasitic main lobe, there is a sidelobe with $-13$ dB level
achieved for both beampatterns in this particular example.

The resulting sensitivity of the beamformer with real-valued weights is
$T_{se}=\mathbf{w}_{r}^{T}\mathbf{w}_{r}=0.076$, which is the same as the
lower bound given by (35): $T_{\min}^{r}=0.076$. This result is not
surprising. To see this, recall the matrix $\mathbf{C}$ defined in (5). For
the linear geometry case the expression reduces to:

$\mathbf{C}=\frac{1}{2}\int\limits_{0}^{\pi}{\mathbf{v}(\theta)\mathbf{v}^{H}(\theta)\sin\theta
d\theta},$ (42)

whose $nm^{th}$ element is given by:

$\displaystyle[\mathbf{C}]_{nm}$
$\displaystyle=\frac{1}{2}\int\limits_{0}^{\pi}{e^{j\frac{2\pi
d}{\lambda}(n-m)\cos\theta}\sin\theta d\theta}$
$\displaystyle=\mathrm{sinc}\left(\frac{2d(n-m)}{\lambda}\right),$ (43)

where $\mathrm{sinc}(x)=\frac{\sin\pi x}{\pi x}$. Thus, for the standard
linear array, i.e. for $\lambda=2d$, as chosen above, $\mathbf{C}$ is the
identity matrix. Hence, the problem of maximum directivity, formulated in
(12), reduces to the minimum sensitivity problem and the solution is optimal
in both senses.

Generally, for the linear array geometry the matrix $\mathbf{C}$ will be close
to the identity matrix for frequencies close to spatial Nyquist frequency of
the array (i.e. $\lambda=2d$), meaning that the maximum directivity beamformer
will be also nearly optimal in the sensitivity sense in this frequency range.
However, for lower frequencies, there will be a strong deviation between the
two optimal beamformers. These relationships will be demonstrated further for
the spherical array geometry.

##  V Spherical Arrays

In this section, we apply the optimal beamformer derived above to spherical
microphone arrays. These arrays are known for their advantages in processing
three-dimensional sound fields due to their rotational symmetry. The spherical
geometry enables the processing to be performed in the spherical harmonics
domain[3] also known as phase-mode processing[28].

###  V-A Phase-mode processing framework

Here we will be concerned with sound fields which are scalar fields
represented by a complex amplitude at each position in space. Consider the
complex amplitude of the sound pressure $p(k,r,\Omega)$ on the surface of a
rigid sphere of radius $r$ at angles $\Omega=(\theta,\phi)$, where
$k=\frac{\omega}{c}$, with $\omega$ representing the angular frequency of the
harmonic sound field and $c$ is the speed of sound. Spherical microphone
arrays sample the surrounding acoustic field at positions
$\\{\Omega_{i}\\}_{i=1}^{M}$, where $M$ is the total number of microphones.
The output of the array is given by:

$y(k,r)=\sum\limits_{i=1}^{M}{\alpha_{i}p(k,r,\Omega_{i})w^{*}(k,\Omega_{i})},$
(44)

where $w^{*}$ is the weighting function of the array. The coefficients
$\alpha_{i}$ are chosen in accordance with the specific sampling scheme [4],
allowing the expression of array output in terms of the spherical Fourier
transform coefficients of the order-limited acoustic pressure and the
weighting function:

$y(k,r)=\sum\limits_{n=0}^{N}{\sum\limits_{m=-n}^{n}{p_{nm}(k,r)w_{nm}^{*}(k)}}.$
(45)

The nearly-uniform sampling scheme will be used in the subsequent discussion,
for which the coefficients are constant and given by
$\alpha_{i}=\frac{4\pi}{M}$ [4].

Consider a unit amplitude plane wave $e^{i(\omega
t+\mathbf{k}^{T}\mathbf{r})}$ impinging on a rigid sphere, where
$\mathbf{r}=(r,\Omega)$ and $\mathbf{k}=(k,\Omega_{0})$, with $\Omega_{0}$
indicating the arrival direction. The spherical Fourier transform coefficients
of the resulting pressure on the sphere surface are given by:

$p_{nm}(k,r,\Omega_{0})=b_{n}(kr){Y_{n}^{m}}^{*}(\Omega_{0}),$ (46)

where $Y_{n}^{m}$ are the spherical harmonics. The mode strength coefficients
$b_{n}$ for a plane wave impinging on a rigid sphere are given by[29]:

$b_{n}(kr)=4\pi
i^{n}\left(j_{n}(kr)-\frac{j_{n}^{\prime}(kr)}{h_{n}^{\prime}(kr)}h_{n}(kr)\right),$
(47)

where $j_{n}$ and $j_{n}^{\prime}$ are the spherical Bessel functions of the
first kind and their derivatives, respectively, and $h_{n}$ and
$h_{n}^{\prime}$ are the spherical Hankel functions of the second kind and
their derivatives, respectively.

The spherical Fourier transform of a plane wave given in (46) is not order
limited. But, in practice, the coefficients $b_{n}(kr)$ for $n>kr$ are
negligible as compared to lower orders. Hence, spherical Fourier transform of
a plane wave with a wavenumber $k<\frac{N}{r}$ can be represented by its lower
order coefficients: $0\leq n\leq N$, with good accuracy.

Now, by substituting (46) for $0\leq n\leq N$ into (45) we get the array
response to a plane wave arriving from the direction $\Omega_{0}$:

$y(k,r,\Omega_{0})=\sum\limits_{n=0}^{N}{\sum\limits_{m=-n}^{n}{b_{n}(kr){Y_{n}^{m}}^{*}(\Omega_{0})w_{nm}^{*}(k)}},$
(48)

which is recognized as the beampattern of the array. In the following, we will
drop the frequency dependence for notation simplicity.

Here, we limit the discussion to weighting functions of the form:

$w_{nm}^{*}=d_{n}{Y_{n}^{m}}(\Omega_{l}),$ (49)

where $\Omega_{l}$ is the desired look direction. This particular form is
highly useful resulting in an axisymmetric weighting function[4] with
$\Omega_{l}$ indicating the axis of symmetry:

$w(\Theta)=\sum\limits_{n=0}^{N}{d_{n}\frac{2n+1}{4\pi}P_{n}(\cos\Theta)},$
(50)

where $\Theta$ is the angle between any direction $\Omega$ and the array look
direction $\Omega_{l}$. The functions $\\{P_{n}(\cdot)\\}$ denote the Legendre
polynomials which are real for real arguments. It can be seen that in order to
have a real-valued weighting function it is sufficient that $d_{n}$ will be
real.

By substituting (49) into (48) we arrive at the desired expression for the
beampattern of a spherical microphone array:

$B(\Theta)=\sum\limits_{n=0}^{N}{d_{n}b_{n}\frac{2n+1}{4\pi}P_{n}(\cos\Theta)}.$
(51)

It can be seen that the design parameters are the coefficients $d_{n}$, which
can be chosen in an appropriate manner. After obtaining the coefficients, the
actual microphone weights can be calculated by choosing the desired look
direction and using (49) and (50).

An important observation is that the phase-mode weights $d_{n}$ are designed
independently from the look direction $\Omega_{l}$. As a consequence,
electronic steering of the beampattern can be accomplished by simple
substitution of the desired look direction when the actual weights $w_{nm}$
are being calculated (see (49)).

We rewrite the expression in (51) using vector notation:

$B(\Theta)=\sum\limits_{n=0}^{N}{d_{n}v_{n}(\Theta)}=\mathbf{d}^{T}\mathbf{v}(\Theta),$
(52)

where $\mathbf{d}=[d_{0}\,\,d_{1}\,\,.\,.\,.\,\,d_{N}]^{T}$ and
$\mathbf{v}=[v_{0}\,\,v_{1}\,\,.\,.\,.\,\,v_{N}]^{T}$, with
$v_{n}=b_{n}\frac{2n+1}{4\pi}P_{n}(\cos\Theta)$.

It can be seen that the form of equation (52) is exactly the same as (3). The
vector $\mathbf{d}$ can be considered as the weights vector in the spherical
harmonics domain and $\mathbf{v}(\Theta)$ is again the manifold vector of the
array.

Note that in the spherical case there is no undesired symmetry property
similar to the symmetry apparent in linear arrays, because the functions
$v_{n}(\Theta)$ in (52) are not symmetric in $\Theta$. This implies that the
beampattern does not have a parasitic main lobe as in the case of linear and
planar arrays where the symmetry is imposed by the conjugate symmetry of
exponentials (see (39)).

###  V-B Optimal phase-mode beamforming with real-valued weights

Consider the phase-mode expression for the spherical array beampattern given
in (52). It is clear that by using the definition of directivity given in (4),
one obtains the expression in (6) with $\mathbf{w}$ replaced by $\mathbf{d}$.
Furthermore, the matrix $\mathbf{C}$ defined in (5) will be real-valued and
diagonal as a straightforward result of the orthogonality property of the
Legendre polynomials:

$\mathbf{C}=\left(\frac{1}{4\pi}\right)^{2}\mathrm{diag}\left(|b_{0}|^{2},3|b_{1}|^{2},...(2N+1)|b_{N}|^{2}\right).$
(53)

Thus, one can directly apply the result in (22) in order to obtain the real-
valued maximum directivity beamformer for the phase-mode processing of a
spherical microphone array.

Spherical array sensitivity can be expressed using the phase-mode weights by
starting from the expression for sensitivity in the spatial domain (see (8)
and [1]):

$\displaystyle T_{se}$
$\displaystyle=\sum\limits_{i=1}^{M}{|\alpha_{i}w(k,\Omega_{i})|^{2}}$
$\displaystyle=\frac{4\pi}{M}\sum\limits_{n=0}^{N}{\sum\limits_{m=-n}^{n}{|w_{nm}|^{2}}}$
$\displaystyle=\frac{4\pi}{M}\sum\limits_{n=0}^{N}{|d_{n}|^{2}\sum\limits_{m=-n}^{n}{Y_{n}^{m}(\Omega_{l}){Y_{n}^{m}}^{*}(\Omega_{l})}}$
$\displaystyle=\frac{1}{M}\sum\limits_{n=0}^{N}{|d_{n}|^{2}(2n+1)}$
$\displaystyle=\mathbf{d}^{H}\mathbf{Ud},$ (54)

where $\mathbf{U}=\frac{1}{M}\mathrm{diag}(1,3,5,...,2N+1)$. Note that the
same expression, up to a normalization constant is used in [30], where the
coefficients $\alpha_{i}$ are not considered as part of the array weighting
vector.

The maximum directivity real-valued phase-mode beamformer with bounded
sensitivity can be now obtained using the solution in (26), with matrix
$\mathbf{\tilde{C}}_{d}=\mathbf{C}+\beta\mathbf{U}$. The $Re\\{\cdot\\}$
operator is omitted because both matrices $\mathbf{C}$ and $\mathbf{U}$ are
real-valued for the spherical geometry case.

The lower bound on sensitivity and the minimum sensitivity beamformer will be
given, again, by the solution of the optimization problem formulated in (12)
with matrix $\mathbf{C}$ replaced by $\mathbf{U}$.

##  VI Simulation Study

###  VI-A Complex-valued vs. real-valued optimal beamformers

Recently, several publications focused on phase-mode beamforming methods with
complex-valued weights for spherical microphone arrays [29, 27, 30, 31, 32,
33]. The purpose of the current section is to compare the techniques with
real-valued weights developed in this paper to their corresponding techniques
with complex-valued weights [27]. The comparison is performed for spherical
microphone array of order $N=10$. Fig. 2 presents the performance measures of
maximum directivity and minimum sensitivity beamformers as the function of the
product $kr$. It should be noted that the minimum sensitivity beamformers were
obtained using the solutions of the optimization problems in (10) and (12) for
the complex-valued and the real-valued cases, respectively, with matrix
$\mathbf{C}$ replaced by $\mathbf{U}$ (see (54)).

Figure 2: Comparison of the performance measures of complex-valued and real-
valued optimal beamformers: (a) - directivity index, (b) - sensitivity.

Recall that the real-valued solutions are sub-optimal in the sense that they
result from the optimization over the space of real numbers which is a
subspace of the complex numbers. Thus, real-valued beamformers are expected to
perform poorer than the more general complex-valued beamformers. Nevertheless,
it can be seen that the performance is not degraded significantly. The
directivity indexes of the real-valued beamformers are lower than the
directivity indexes of the complex-valued ones only by $3$ dB for the entire
frequency range in this particular example. The difference in the
sensitivities between the maximum directivity beamformers is about $5$ dB for
lower frequencies, decreasing to $3$ dB for higher frequencies.

In addition, note that in both cases (real and complex), the maximum
directivity beamformers tend to be optimal also in the sensitivity sense at
high $kr$ ranges. This behavior is similar to the one observed for the linear
array geometry as discussed at the end of Section IV. This is due to the fact
that the absolute value of the mode strength coefficients $b_{n}(kr)$ when
$kr$ approaches $N$ are nearly independent of the order $n$, i.e.
$|b_{n}(kr)|\approx|b_{0}(kr)|$. This implies that
$\mathbf{C}\approx\frac{M|b_{0}(kr)|^{2}}{(4\pi)^{2}}\mathbf{U}$. Hence, the
maximum directivity and the minimum sensitivity problems tend to be the same
for higher frequencies.

###  VI-B Sidelobe performance

Here we use (22) and (52) in order to calculate the beampattern of the real-
valued maximum directivity beamformer. Fig. 3 illustrates the simulated
maximum directivity beampatterns at $kr=10$ for the spherical microphone array
of order $N=10$.

Figure 3: Maximum directivity beampatterns with real-valued weights for a
spherical array of order $N=10$, at $kr=10$.

The array is steered to $\Theta=0$. It can be seen that there is indeed a main
lobe in the desired direction with a sidelobe level of about $-8$ dB. As can
be seen, the relatively poor sidelobe-level performance is due to a side lobe
at $\Theta=\pi$.

In order to find the reason for the appearance of the sidelobe at $\Theta=\pi$
we return to the definition of directivity in Section II-B. It was seen that
matrix $\mathbf{C}$ defined in (5) includes the $\sin\Theta$ function (when
the look direction $\Omega_{l}$ coincides with the positive $z$ axis, we get
$\Theta=\theta$), which, as mentioned above, plays the role of a cost
function. This implies that the solution to the optimization problem defined
by (12) will tend to pull the energy towards $\sin\Theta=0$, minimizing the
response at the center where the cost is high.

With this in mind, one can try different cost functions in order to affect the
sidelobe level or the gain distribution along $\Theta$. Although cost
functions other than $\sin\Theta$ will cause degradation in directivity, they
can be beneficial in decreasing the sidelobe level, while at the same time
affecting the main lobe width, as will be demonstrated shortly.

Several alternative beampatterns were calculated using linear, uniform and
unit step cost functions. The results are summarized in Fig. 4 and Table I.

---
Figure 4: Spherical array beam patterns using various cost functions:
(a)-linear, (b)-uniform, (c)-step function.

It can be observed that the linear cost function (Fig. 4(a)) imposes higher
cost at $\Theta=\pi$ and the high-level sidelobe at $\Theta=\pi$ is avoided.
In addition, using the unit step cost function (Fig. 4(c)) and the uniform
cost function (Fig. 4(b)), results in a more uniform spread of side-lobe
energy. However, the three beampatterns in Fig. 4 have wider main lobes as
compared to the maximum directivity beampattern (Fig. 3) as well as lower
directivity.

Note that the array sensitivities using different cost functions are only
slightly above the sensitivity of the maximum directivity beamformer (see
Table I). Thus, at least for the parameters of the simulation presented here,
the use of different cost functions did not impact on the array sensitivity.

TABLE I: Comparison of beampattern parameters using various cost functions for spherical array with lower sensitivity bound of $T_{\min}=-22.4$ [dB]. Cost func. | Sidelobe [dB] | DI [dB] | Sens. [dB]
---|---|---|---
sin | -7.9 | 18.5 | -22.3
linear | -18.1 | 17.3 | -20.6
uniform | -13.6 | 17.9 | -21.8
step function | -18.5 | 17.9 | -21.3

##  VII Experimental Study

Spherical microphone arrays can be utilized for sound field analysis. Examples
of plane wave decomposition (PWD) of the sound field using maximum directivity
complex-valued beamformer can be found in the literature [34]. However, the
real-valued maximum directivity beamformer derived in this paper is suboptimal
in the sense that the weights vector is constrained to be real-valued.

Here, we compare the performance of maximum directivity beamformers, the
optimal complex-valued beamformer and the suboptimal real-valued beamformer,
by using experimental data.

The large-version of the Eigenmike® spherical microphone array[35], consisting
of $32$ microphones mounted on a rigid sphere of radius $r=9$ cm was employed.
Note that the array order is $N=4$, meaning that its maximum operating
frequency is about $2400$ Hz. The measurements were performed in an anechoic
chamber (anechoic from $300$ Hz) chamber with dimensions of $2\times 2\times
2$ m. A B&K OmniSource 4295 loudspeaker was used. It was positioned such that
the direction of propagation relative to the array was
$(\theta,\phi)=(100^{\circ},160^{\circ})$. By selecting the response from the
loudspeaker to the microphone array at a desired frequency, the vector
$\mathbf{p}$ containing complex pressure amplitudes at different microphones
is obtained. Then, the spherical Fourier transform coefficient $p_{nm}$ of the
pressure on the sphere are calculated in the least squares sense using the
pseudo-inverse of the transform matrix [36].

The spherical Fourier transform coefficients of the weighting function of the
complex-valued maximum directivity beamformer are given by [4]:

$w_{nm}^{*}=\frac{Y_{n}^{m}(\Omega_{l})}{b_{n}},$ (55)

meaning that $d_{n}=\frac{1}{b_{n}}$, which is equivalent to the result in
(10) up to a normalization constant $\frac{4\pi}{(N+1)^{2}}$ . Thus, the
output of the beamformer when steered to $\Omega_{l}$ is given by:

$y(\Omega_{l})=\sum\limits_{n=0}^{N}{\sum\limits_{m=-n}^{n}{\frac{p_{nm}}{b_{n}}}Y_{n}^{m}(\Omega_{l})}.$
(56)

Thus, PWD is performed by dividing the measured pressure coefficients by
$b_{n}$ and transforming back into the spatial domain. Similarly, when using
the optimal real-valued beamformer
$\mathbf{d}=[d_{0}\,\,d_{1}\,\,.\,.\,.\,\,d_{N}]^{T}$ the pressure
coefficients $p_{nm}$ are multiplied by $d_{n}$ and transformed into the
spatial domain. The results of spatial sound-field analysis using PWD are
presented in Fig. 5.

---
Figure 5: Spatial sound-field analysis at $2400$ Hz ($N=4$, $r=9$ cm) using
various beamformers: (a) complex-valued max. directivity (b) real-valued max.
directivity (c) real-valued with linear cost function.

It can be seen that as expected, all three results have a maximum in the
proximity of the loudspeaker direction. However, the analysis performed using
real-valued maximum directivity beamformer (Fig. 5.(b)) has an extra peak in
the reciprocal direction. This is the straightforward consequence of the high
sidelobe at $\Theta=\pi$ as discussed in Section VI. It can be seen that, by
using the real-valued beamformer with a linear cost function as proposed in
section VI, the undesired peak is completely eliminated. However, the
resolution of the analysis is degraded as a consequence of the wider main lobe
and lower directivity.

##  VIII Conclusion

Beamforming techniques with real-valued weights have the potential to
significantly simplify beamformer structure. Two optimal beamforming
techniques using real weights were presented in this paper. It was shown that
these techniques can be used with any desired array geometry. In particular,
when used with spherical array geometry, the techniques result in
electronically steerable beampatterns. A degree of flexibility was introduced
by using various cost functions. The results of numerical simulations and the
experimental study show that, as expected, the performance of optimal
beamformers with real-valued weights is lower than that of the optimal
beamformers with complex-valued weights, but not considerably so in light of
the potential hardware savings.

## References

  * [1] H. L. Van Trees, _Optimum Array Processing (Detection, Estimation and Modulation Theory, Part IV)_. New York: Wiley Interscience, 2002.
  * [2] H. Krim and M. Viberg, “Two decades of array signal processing research,” _IEEE Signal Processing Magazine_ , vol. 13, no. 4, pp. 67–94, July 1996.
  * [3] E. G. Williams, _Fourier Acoustics: Sound Radiation and Nearfield Acoustical Holography_. New York: Academic, 1999.
  * [4] B. Rafaely, “Analysis and design of spherical microphone array,” _IEEE Trans. on Speech Audio Process._ , vol. 13, no. 1, Jan. 2005.
  * [5] J. Meyer and G. W. Elko, “A spherical microphone array for spatial sound recordings,” _J. Acoust. Soc. Am_ , vol. 111, no. 5.2, pp. 2346–2346, 2002.
  * [6] B. D. Van Veen and K. M. Buckley, “Beamforming: A versatile approach to spatial filtering,” _IEEE ASSP Magazine_ , Apr. 1988.
  * [7] L. Griffiths and C. Jim, “An alternative approach to linearly constrained adaptive beamforming,” _IEEE Trans. Antennas Propag._ , vol. 30, no. 1, pp. 27–34, Jan. 1982.
  * [8] D. I. Havelock, “Sensor array beamforming using random channel sampling: The aggregate beamformer,” _J. Acoust. Soc. Am_ , vol. 114, pp. 1997–2006, Oct. 2003.
  * [9] D. Havelock, “Residual noise in the aggregate beamformer,” in _2003 IEEE Workshop on Applications of Signal Processing to Audio and Acoustics_ , Oct. 2003, pp. 45–48.
  * [10] M. Agmon, B. Rafaely, and J. Tabrikian, “Maximum directivity beamformer for spherical-aperture microphones,” in _2009 IEEE Workshop on Applications of Signal Processing to Audio and Acoustics_ , Oct. 2009, pp. 153–156.
  * [11] T. Vu, “Simultaneous nulling in sum and difference patterns by amplitude control,” _IEEE Trans. Antennas Propag._ , vol. 34, no. 2, pp. 214–218, Feb. 1986.
  * [12] K. Lo, “Adaptivity of a real-symmetric array by DOA estimation and null steering,” _IEE Proceedings, Radar, Sonar and Navigation_ , vol. 144, no. 5, pp. 245–251, Oct. 1997.
  * [13] T. Gao, Y. Guo, and J. Li, “A fast beamforming algorithm for adaptive sum and difference patterns in conformal array antennas,” in _Antennas and Propagation Society International Symposium, 1992. AP-S, Digest._ , vol. 1, July 1992, pp. 450–453.
  * [14] G. Yanchang and L. Jianxin, “Real amplitude-only nulling algorithm (RAMONA) for adaptive digital beamforming,” in _Antennas and Propagation Society International Symposium, 1990. AP-S. Merging Technologies for the 90’s. Digest._ , vol. 1, May 1990, pp. 206–209.
  * [15] Y. Guo and J. Li, “Real amplitude-only nulling algorithm (RAMONA) for adaptive sum and difference patterns,” in _Antennas and Propagation Society International Symposium, 1991. AP-S. Digest_ , vol. 1, June 1991, pp. 94–97.
  * [16] G. Yanchang and L. Jianxin, “Speed up and optimization of RAMONA for adaptive digital beamforming,” in _Eighth International Conference on Antennas and Propagation, 1993._ , vol. 1, 1993, pp. 508–511.
  * [17] Y.-H. Choi, “Simple adaptive beamforming with only real weights based on covariance differencing,” _Electronics Letters_ , vol. 43, no. 10, pp. 552–553, Oct. 2007.
  * [18] W. Choi, T. Sarkar, H. Wang, and E. Mokole, “Adaptive processing using real weights based on a direct data domain least squares approach,” _IEEE Trans. Antennas Propag._ , vol. 54, no. 1, pp. 182–191, Jan. 2006.
  * [19] M. Ghavami, “Wideband beamforming using rectangular arrays without phase shifting,” _European Transactions on Telecommunications_ , vol. 14, pp. 449–456, 2003.
  * [20] W. Buehring, “Adaptive antenna with rapid convergence,” in _IEEE Trans. Antennas Propag._ , 1978, pp. 51–54.
  * [21] M. Paajanen, J. Lekkala, and K. Kirjavainen, “Electromechanical film (emfi) – a new multipurpose electret material,” _Sensors and Actuators A: Physical_ , vol. 84, no. 1-2, pp. 95–102, 2000.
  * [22] J. Lekkala and M. Paajanen, “Emfi - new electret material for sensors and actuators,” in _$10^{th}$ International Symposium on Electrets_, 1999.
  * [23] J. L. Ealo, F. Seco, and A. R. Jimenes, “Broadband emfi-based transducers for ultrasonic air applications,” _IEEE Trans. Ultrason., Ferroelectr., Freq. Control_ , vol. 55, no. 4, Apr. 2008.
  * [24] J. L. Ealo, F. Seco, C. Prieto, A. R. Jimenes, J. Roa, A. Koutsou, and J. Guevara, “Customizable field airborne ultrasonic transducers based on electromechanical film,” in _IEEE International Ultrasonics Symposium Proceedings_ , 2008.
  * [25] S. Yan, H. Sun, X. Ma, U. Svensson, and C. Hou, “Time-domain implementation of broadband beamformer in spherical harmonics domain,” _IEEE Trans. Audio, Speech Lang. Process._ , vol. 19, no. 5, pp. 1221–1230, July 2011.
  * [26] B. Rafaely, “Plane-wave decomposition of the sound field on a sphere by spherical convolution,” _J. Acoust. Soc. Am_ , vol. 116, no. 4, Oct. 2004.
  * [27] B. Rafaely, Y. Peled, M. Agmon, D. Khaykin, and E. Fisher, “Spherical microphone array beamforming,” in _Speech Processing in Modern Communications: challenges and perspectives_ , I. Cohen, J. Benesty, and S. Gannot, Eds. Berlin: Springer-Verlag, 2010, ch. 11, pp. 281–305.
  * [28] B. Rafaely, “Phase-mode versus delay-and-sum spherical microphone array processing,” _IEEE Signal Processing Letters_ , vol. 12, no. 10, pp. 713–716, Oct. 2005.
  * [29] J. Meyer and G. W. Elko, “A highly scalable spherical microphone array based on an orthonormal decomposition of the sound field,” _Proc. ICASSP_ , vol. II, pp. 1949–1952, 2002.
  * [30] H. Sun, S. Yan, and U. Svensson, “Robust minimum sidelobe beamforming for spherical microphone arrays,” _IEEE Trans. Audio, Speech Lang. Process._ , vol. 19, no. 4, pp. 1045–1051, May 2011.
  * [31] H. Sun, S. Yan, and U. Peter Svensson, “Robust spherical microphone array beamforming with multi-beam-multi-null steering, and sidelobe control,” in _2009 IEEE Workshop on Applications of Signal Processing to Audio and Acoustics_ , Oct. 2009, pp. 113–116.
  * [32] H. Sun, S. Yan, and U. Svensson, “Optimal higher order ambisonics encoding with predefined constraints,” _IEEE Trans. Audio, Speech Lang. Process._ , vol. 20, no. 3, pp. 742–754, Mar. 2012.
  * [33] S. Yan, H. Sun, U. Svensson, X. Ma, and J. Hovem, “Optimal modal beamforming for spherical microphone arrays,” _IEEE Trans. Audio, Speech Lang. Process._ , vol. 19, no. 2, pp. 361–371, Feb. 2011.
  * [34] B. Rafaely, I. Balmages, and L. Eger, “High-resolution plane-wave decomposition in an auditorium using a dual-radius scanning spherical microphone array,” _J. Acoust. Soc. Am_ , vol. 122, no. 5, pp. 2661–2668, 2007.
  * [35] _em32 Eigenmike microphone array release notes_ , mh acoustics, 25 Summit Ave, Summit, NJ 07901, Feb. 2009.
  * [36] B. Rafaely and A. Avni, “Interaural cross correlation in a sound field represented by spherical harmonics,” _J. Acoust. Soc. Am_ , vol. 127, no. 2, pp. 823–828, 2010.

|  Vladimir Tourbabin (S’12) received the B.Sc. degree (summa cum laude) in
materials science and engineering and the M.Sc. degree (cum laude) in
electrical and computer engineering from Ben-Gurion University of the Negev,
Israel, in 2005 and 2011, respectively. He is currently working towards the
Ph.D. degree in electrical and computer engineering at Ben-Gurion University.
His current research interests include microphone array signal processing and
human sound perception. Mr. Tourbabin is a recipient of the Negev Faran
Fellowship.
---|---
|  Morag Agmon received the B.Sc. (cum laude) and M.Sc. degrees in electrical
and computer engineering from the Ben-Gurion University of the Negev, Beer-
Sheva, Israel, in 2008 and 2010, respectively. He is currently leading the DSP
Cellular Voice and Speech group in Marvell Semiconductor Inc. His research
interests include acoustical signal processing and real-time microphone
arrays.
---|---
|  Boaz Rafaely (SM’01) received the B.Sc. degree (cum laude) in electrical
engineering from Ben-Gurion University, Beer-Sheva, Israel, in 1986; the M.Sc.
degree in biomedical engineering from Tel-Aviv University, Israel, in 1994;
and the Ph.D. degree from the Institute of Sound and Vibration Research
(ISVR), Southampton University, U.K., in 1997. At the ISVR, he was appointed
Lecturer in 1997 and Senior Lecturer in 2001, working on active control of
sound and acoustic signal processing. In 2002, he spent six months as a
Visiting Scientist at the Sensory Communication Group, Research Laboratory of
Electronics, Massachusetts Institute of Technology (MIT), Cambridge,
investigating speech enhancement for hearing aids. He then joined the
Department of Electrical and Computer Engineering at Ben-Gurion University as
a Senior Lecturer in 2003, and appointed Associate Professor in 2010. He is
currently heading the acoustics laboratory, investigating sound fields by
microphone and loudspeaker arrays. Since 2010, he is serving as an associate
editor for IEEE Transactions on Audio, Speech and Language Processing. Prof.
Rafaely was awarded the British Council’s Clore Foundation Scholarship.
---|---
|  Joseph Tabrikian (S’89-M’97-SM’98) received the B.Sc., M.Sc., and Ph.D.
degrees in Electrical Engineering from the Tel-Aviv University, Tel-Aviv,
Israel, in 1986, 1992, and 1997, respectively. During 1996-1998 he was with
the Dept. of ECE, Duke University, Durham, NC as an Assistant Research
Professor. He is now with the Dept. of ECE, at the Ben-Gurion University of
the Negev, Beer-Sheva, Israel. His research interests include estimation and
detection theory and array signal processing. He has served as an Associate
Editor of the IEEE Transactions on Signal Processing during 2001-2004, and
since 2011, and as an Associate Editor of the IEEE Signal Processing Letters
since 2012. He is a member of the IEEE SAM technical committee since 2010 and
was the technical program co-chair of the IEEE SAM 2010 workshop. He is
coauthor of 4 award-winning papers in different IEEE conferences.
---|---