# 2D-DOA Estimation and Auto-Calibration in UCAs via an Integrated Wideband
Dictionary

Zavareh Bozorgasl, Hao Chen, Mohammad J. Dehghani Z. Bozorgasl .H. Chen is
with the Department of Electrical and Computer Engineering, Boise State
University, Boise, ID, 83712 (email<EMAIL_ADDRESS>M. J. Dehghani is
with the Department of Electrical and Electronics Engineering, Shiraz
University of Technology, Shiraz, Iran (email: dehghani@sutech.ac.ir).

###### Abstract

In this paper, we present a novel auto-calibration scheme for the joint
estimation of the two-dimensional (2-D) direction-of-arrival (DOA) and the
mutual coupling matrix (MCM) for a signal measured using uniform circular
arrays. The method employs an integrated wideband dictionary to mitigate the
detrimental effects of the discretization of the continuous parameter space
over the considered azimuth and elevation angles. This leads to a reduction of
the computational complexity and obtaining of more accurate DOA estimates.
Given the more reliable DOA estimates, the method also allows for the
estimation of more accurate mutual coupling coefficients. The method utilizes
an integrated dictionary in order to iteratively refine the active parameter
space, thereby reducing the required computational complexity without reducing
the overall performance. The complexity is further reduced by employing only
the dominant subspace of the measured signal. Furthermore, the proposed method
does not require a constraint on the prior knowledge of the number of nonzero
coupling coefficients nor suffer from ambiguity problems. Moreover, a simple
formulation for 2-D non-numerical integration is presented. Simulation results
show the effectiveness of the proposed method.

Index Terms—Direction-of-Arrival, Sensor Array Processing, Sparse
Representation, Integrated Dictionary, Auto-Calibration, Uniform Circular
Array, Mutual Coupling Matrix.

##  I Introduction

Estimating the direction-of-arrival (DOA) of a plane wave impinging on a
sensor array, a classical problem in array signal processing, has been widely
investigated in recent decades [1, 2, 3]. Although the literature mostly deals
with DOA estimation techniques which use the assumption of a perfectly known
measurement array, the resulting steering vector will generally only be
approximately known. Indeed, the manifold will suffer from effects such as
mutual coupling, which may cause a substantial degradation of the performance
of applied algorithms [4, 5, 6]. Two of the most commonly used array
configurations are the uniform linear array (ULA) and the uniform circular
array (UCA). As the mutual coupling effect will be more pronounced for a
circular than for a linear array, the detrimental effects of such coupling
will affect the resulting performance more, suggesting the current study.

In order to alleviate the mutual coupling effects, two main categories of
algorithms have been proposed in the literature. The first of these employ
offline calibration algorithms[7, 8], which require calibration sources, i.e.,
some sources with perfectly known locations, to determine the array manifold
matrix. In contrast, the second category uses auto-calibration or online
calibration algorithms[9, 10, 11, 12], which do not require such source
calibration, instead aiming to compensate for the calibration automatically.
As finding additional calibration sources is not possible in many practical
applications, the auto-calibration algorithms are often preferable.

Moreover, as the discretization of the parameter space leads to undesirable
off-grid effects, there has recently been increasing interest in continuous
parameter space techniques[13, 14]. One such alternative is to use a
continuous dictionary together with an atomic norm penalty[15, 16, 17]. While
such an approach often yields an accurate signal reconstruction, this form of
problem formulation often leads to optimization problems with a high
computational complexity. In this work, instead, we propose the use of an
integrated wideband dictionary, as introduced [18], which is formed over
subsets of the continuous parameter space. The essence of using this kind of
dictionary is the ability to discard the non-activated subsets and preserving,
then refining the activated ones for further zooming without the risk of
missing components. These steps are then iterated until one achieves the
desired resolution.

Common problems of earlier auto-calibration, including ambiguity problems [11,
12] and/or the requirement of a priori knowledge of the number of nonzero
mutual coupling coefficients [19, 20], both typically being infeasible in
practice. The here proposed method, which iteratively estimates both the DOAs
and the mutual coupling coefficients, suffers from neither of these
shortcomings. For the DOA estimation stage, instead of using a conventional
grid-based dictionary which would necessitate a high computational complexity
to allow for a sufficiently fine grid, a wideband dictionary framework is
employed. By forming the used dictionary consisting of subsets of the
continuous parameter space, and applying the noted screening process, the non-
activated sets may be discarded in the first stage, whereas the activated
subsets are retained and refined for further processing. In this sense, the
initial DOA estimation stage will constitute an initial coarse estimate that
will not suffer from the usual off-grid effects, such as missing active
components. Via the further steps, i.e., refining the remaining activated
bands, the closely spaced DOAs are separated. Next, the MCM estimation is
formed utilizing the complex symmetric circular Toeplitz structure of the MCM,
as was proposed for 1-D DOA estimation in [21, 22]. Here, this mutual
coefficient estimation stage is further extended to the 2-D DOA estimation
problem for the UCA configuration.

The remainder of the paper is organized as follows. Section II outlines the
data model for 2-D UCAs. Section III describes the details of the new method
for joint iterative 2-D DOA and MCM estimation. Section IV provides some
simulations to demonstrate the effectiveness of the proposed method. Finally,
Section V concludes the work.

Notation: Bold lower-case (upper-case) letters are used throughout to denote
vectors (matrices). $\mathbb{C}$ denotes the set of complex numbers.
$(\cdot)^{T}$ and $(\cdot)^{H}$ represent the transpose and the complex
conjugate transpose, respectively. The notations $\|\cdot\|_{k}$ and
$\|\cdot\|_{F}$ stand for the k-norm of a vector and the Frobenius norm of a
matrix, respectively. $\hat{(\cdot)}$ is an estimate of $(\cdot)$, and
$\text{diag}(\mathbf{x})$ is a diagonal matrix with $\mathbf{x}$ being its
diagonal elements. $\mathbf{X}(:,i)$ and $\mathbf{X}(j,:)$ show the $i$-th
column and the $j$-th row of $\mathbf{X}$, respectively.
$\text{vec}(\mathbf{X})$ stands for the vector with the elements consists of
columns of $\mathbf{X}$ stacked together.

##  II SIGNAL MODEL

Suppose that $K$ narrowband far-field sources simultaneously impinge onto a
uniform circular array (UCA) of omnidirectional sensors from angles of arrival
$(\phi,\theta)$, where $\phi$ denotes the azimuth angle, i.e.,
$\phi\in[0,360)$ degrees, measuring counterclockwise from the $x$-axis, and
$\theta$ denotes the elevation angle, i.e., $\theta\in[0,90]$ degrees,
measuring down from the $z$-axis. The UCA geometry, with $N$ equispaced
sensors on the circumference of the array, is shown in Fig. 1 in the
$xy$-plane. The marked arrow represents an example of a received signal
geometrical description.

Figure 1: The considered UCA receiver model for DOA and MCM estimation.

Therefore, by considering an array of $N$ sensors which are placed uniformly
on a circle of radius $r$, the array output $\mathbf{x}\in\mathbb{C}^{N\times
T}$ may be expressed as

$\mathbf{x}(t)=\mathbf{A}(\boldsymbol{\phi},\boldsymbol{\theta})+\mathbf{e}(t),\quad
t=t_{1},t_{2},\ldots,t_{T},$ (1)

where
$\mathbf{A}(\boldsymbol{\phi},\boldsymbol{\theta})=[\mathbf{a}(\phi_{1},\theta_{1}),\mathbf{a}(\phi_{2},\theta_{2}),\ldots,\mathbf{a}(\phi_{K},\theta_{K})]$
is the array manifold/steering matrix, with

$\mathbf{a}(\phi_{i},\theta_{i})=\exp(j\xi_{i}\cos(\phi_{i}-\gamma_{n})),\quad
n=0,1,\ldots,N-1,$

denoting the steering vector of the source from $(\phi_{i},\theta_{i})$. Here,
$\gamma_{n}=\frac{2\pi n}{N}$ is the displacement of the $n$th element of the
array from the $x$ axis, $\xi_{i}=k_{0}r\sin\theta_{i}$, with the wave number
$k_{0}=\frac{2\pi}{\lambda}$, whereas
$\mathbf{s}(t)=[s_{1}(t),s_{2}(t),\ldots,s_{K}(t)]^{T}$ and
$\mathbf{e}(t)=[e_{1}(t),e_{2}(t),\ldots,e_{N}(t)]^{T}$ denote the source data
and an additive white Gaussian noise vector, respectively. Finally, $t$
indexes each snapshot and $T$ is the number of snapshots. The signal vector
$\mathbf{s}(t)$ and the noise vector $\mathbf{e}(t)$ are assumed to be zero
mean circularly symmetric and statistically independent processes.

##  III THE PROPOSED AUTOCALIBRATION METHOD

By taking mutual coupling into consideration, the observation model in (1) can
be reformulated as

$\mathbf{X}=\mathbf{C}\mathbf{A}(\phi,\theta)\mathbf{S}+\mathbf{E},$ (2)

where
$\mathbf{X}=[\mathbf{x}(t_{1}),\mathbf{x}(t_{2}),\ldots,\mathbf{x}(t_{T})]$,
$\mathbf{S}=[\mathbf{s}(t_{1}),\mathbf{s}(t_{2}),\ldots,\mathbf{s}(t_{T})]$,
and
$\mathbf{E}=[\mathbf{e}(t_{1}),\mathbf{e}(t_{2}),\ldots,\mathbf{e}(t_{T})]$
are the array output, as well as the source and noise data matrices,
respectively. Furthermore, $\mathbf{C}$ denotes the mutual coupling matrix
which is assumed to have a complex symmetric circular Toeplitz structure such
that

$\mathbf{C}=\text{toeplitz}(\tilde{\mathbf{c}},\tilde{\mathbf{c}})\in\mathbb{C}^{N,N},$
(3)

where

$\mathbf{\tilde{c}}=[c_{1},c_{2},\ldots,c_{L},c_{L-1},\ldots,c_{3},c_{2}]\in\mathbb{C}^{1,N},\quad
L=\frac{N}{2}+1,\text{ for }N\text{ even,}$

and

$\tilde{\mathbf{c}}=[c_{1},c_{2},\ldots,c_{L},\ldots,c_{3},c_{2}]\in\mathbb{C}^{1,N},\quad
L=\left(\frac{N+1}{2}\right),\text{ for }N\text{ odd.}$

Here, $L$ denotes the degrees of freedom of the MCM, i.e., the number of the
unknown elements in the first row of the matrix $\mathbf{C}$. We note that
because of the inverse relationship of mutual coupling and distance between
each pairs of sensors $0<|c_{L}|<\ldots<|c_{2}|<|c_{1}|=1$.

Algorithm 1 summarizes the proposed method for estimating the DOA and the MCM.
Each substep will be elaborated upon in the following subsections.

Algorithm 1 The proposed SVD-based sparse WB method

1:Assign an initial MCM such that $\mathbf{C}_{0}=\mathbf{I}_{N\times N}$.

2:Choose the number of zooming steps and number of bands to use for both
azimuth and elevation; select a suitable $\alpha$ for Eq. (11).

3:Compute the number of sources and take the SVD of Eq. (2).

4:Compose the DOA integrated wideband dictionary with elements defined in Eq.
(6).

5:Estimate the azimuth and elevation (as stated in Subsection III-B) using Eq.
(10).

6:Estimate the MCM as detailed in Subsection III-C.

7:Form the DOA estimates on the active bands (composing of the zoomed
dictionary given in Eq. (6)), and then using Eq. (10) include the MCM in step
6. The iterative zooming procedure enables discarding the non-activated parts
of each dictionary.

8:Estimate the MCM as stated in Subsection III-C

9:Repeat steps 7 and 8 until obtaining the desired resolution.

###  III-A Sparse SVD-based Representation

The overall complexity of the proposed method may be notably reduced by
projecting the measured data onto the $K$-dimensional signal subspace. This
can be formed by the $K$ dominant singular vectors of $\mathbf{X}$. Form the
singular value decomposition (SVD) as:

$\mathbf{X}=\mathbf{U}_{s}\Sigma_{s}\mathbf{V}_{s}^{H}+\mathbf{U}_{c}\Sigma_{c}\mathbf{V}_{c}^{H},$
(4)

where $\mathbf{U}_{s}\in\mathbb{C}^{N\times K}$ and
$\mathbf{V}_{s}\in\mathbb{C}^{T\times K}$ are the singular vectors
corresponding to the $K$ largest singular values, and $K$ represents the
number of signals (sources). To do this, the number of sources, $K$, must also
be estimated. Here, this is done by forming the resulting solution over a set
of potential $K\in\\{1,2,\ldots\\}$ and then selecting the estimated $K$ as
the one resulting in the best signal modeling when combined with a BIC penalty
to compensate for the measuring model order [23]. Let
$\mathbf{X}_{\text{SV}}\triangleq\mathbf{X}V_{s}$,
$\mathbf{S}_{\text{SV}}\triangleq\mathbf{S}V_{s}$, and
$\mathbf{E}_{\text{SV}}\triangleq\mathbf{E}V_{s}$, such that

$\mathbf{X}_{\text{SV}}=\mathbf{C}\mathbf{A}(\boldsymbol{\phi},\boldsymbol{\theta})\mathbf{S}_{\text{SV}}+\mathbf{E}_{\text{SV}},$
(5)

Indeed, by using the $N\times K$ dimensional matrix $\mathbf{X}_{\text{SV}}$
in place of $\mathbf{X}$, one may significantly reduce the resulting
computational complexity. The reduction will be significant in case of a small
number of sources and/or a large number of time samples.

###  III-B DOA Estimation

Proceeding to determine the DOAs, let $(\phi_{i},\theta_{i})_{i=1}^{Q}$ for
$i=1,2,...,Q$, stands for the set of potential DOA candidates, with $Q$
indicating the number of candidates. Denote the dictionary over the set of
candidates
$\mathbf{A}(\hat{\phi},\hat{\theta})=[\mathbf{a}(\hat{\phi}_{1},\hat{\theta}_{1}),\mathbf{a}(\hat{\phi}_{2},\hat{\theta}_{2}),...,\mathbf{a}(\hat{\phi}_{Q},\hat{\theta}_{Q})]$,
where
$\mathbf{a}(\hat{\phi}_{i},\hat{\theta}_{i})=[\mathbf{a}_{0}(\hat{\phi}_{i},\hat{\theta}_{i}),\mathbf{a}_{1}(\hat{\phi}_{i},\hat{\theta}_{i}),...,\mathbf{a}_{N-1}(\hat{\phi}_{i},\hat{\theta}_{i})]^{T}$
and each candidate is formed using the integrated wideband element,
$\mathbf{a}_{n_{r}}(\phi,\theta)$, such that

$\mathbf{a}_{n_{r}}(\phi_{i},\theta_{i})=\int_{\phi_{i}}^{\phi_{i+1}}\int_{\theta_{i}}^{\theta_{i+1}}\exp\left(jk_{0}r\sin\theta\cos\left(\phi-\frac{2\pi
n_{r}}{N}\right)\right)d\theta d\phi,$ (6)

with $n_{r}=0,1,...,N-1$, and where $[\phi_{i},\phi_{i+1})$ and
$[\theta_{i},\theta_{i+1})$ are the intervals for the $i$th azimuth and
elevation candidate. Indeed, instead of using a finely spaced narrowband
dictionary of the candidate DOAs, the 2-D angles domain, i.e., azimuth
$\phi:[0^{\circ},360^{\circ})$ and elevation $\theta:[0^{\circ},90^{\circ})$,
is here divided into $B_{1}$ and $B_{2}$ continuous bands, respectively. As
long as $B_{1}$ and $B_{2}$ are selected to be sufficiently large, this
process would reduce the risk of any off-grid components and provide an
initial coarse estimation of the angle of arrivals. Since all the bands are
covered in the dictionary, no power is off-grid, which in turn avoids a non-
sparse solution due to dictionary mismatch. In its initial step, this will
yield a coarse estimation of the regions of interest, whereafter one may
discard the non-activated bands and form the refined dictionary over the
activated bands. This process may then be repeated until one reaches the
desired resolution. In other words, instead of making a dictionary over a fine
grid to limit the performance degradation due to off-grid effects, thereby
necessitating a substantial computational burden, we use the integrated
dictionary to determine and refine the estimate of the regions where sources
are present. The process of zooming in some steps (here, 3 steps) is
illustrated in Fig. 2, for the case of a single angle $\theta$. As shown in
[18], the average computation time of employing the wideband dictionary is
substantially less than what would be required by the corresponding narrowband
dictionary. The main difference of this zooming with respect to the earlier
zooming methods is in utilizing an integrated dictionary for bands which
significantly improves the performance.

Figure 2: The process of zooming on intervals/bands of angles.

Suppose the arrow in Fig. 2 shows the location of one of the sources. Then,
one may describe the process of zooming in three stages as follows:

  1. 1.

Initially, one forms an integrated dictionary over the coarse grid
($\mathbf{\theta}_{0}$).

  2. 2.

Then, one identifies the regions that may contain a source (as an example, the
right side of $\mathbf{\theta}_{1}$ is a possible region).

  3. 3.

Finally, one composes an integrated dictionary with finer grids over the
possible regions identified in step 2 (the third step is shown in the
$\mathbf{\theta}_{2}$ region).

Remark: In Appendix A, we introduce an approach to form the 2-D integral in
equation (6).

Proceeding, the azimuth and elevation estimation are formulated using the
LASSO [24], which solves the resulting penalized regression problem. To this
end, let
$\mathbf{\hat{A}}=\mathbf{\hat{C}}\mathbf{A}(\hat{\phi},\hat{\theta})\in\mathbb{C}^{M\times
Q}$ and then

$\mathbf{D}=\begin{bmatrix}\mathbf{\hat{A}}\\\ \vdots\\\
\mathbf{\hat{A}}\end{bmatrix}_{MK\times Q}$ (7)

Then,

$\begin{bmatrix}\mathbf{X}_{\text{SV}}(:,1)\\\ \vdots\\\
\mathbf{X}_{\text{SV}}(:,K)\end{bmatrix}=\mathbf{D}\mathbf{\hat{s}}+\begin{bmatrix}\mathbf{N}_{\text{SV}}(:,1)\\\
\vdots\\\ \mathbf{N}_{\text{SV}}(:,K)\end{bmatrix}$ (8)

Let $\mathbf{x}=\text{vec}(\mathbf{X}_{\text{sv}})\in\mathbb{C}^{MK}$ and
$\mathbf{n}=\text{vec}(\mathbf{N}_{\text{sv}})\in\mathbb{C}^{MK}$. Then,

$\mathbf{x}=\mathbf{D}_{MK\times Q}\mathbf{\hat{s}}+\mathbf{n},$ (9)

allowing the LASSO formulation1, i.e., $\mathbf{\hat{s}}$, results in

$\mathbf{\hat{s}}=\arg\min_{\mathbf{\hat{s}}}\left\\{\left\|\mathbf{x}-\mathbf{D}\mathbf{\hat{s}}\right\|_{2}^{2}+\gamma\left\|\mathbf{\hat{s}}\right\|_{1}\right\\},$
(10)

where $\mathbf{\hat{s}}\in\mathbb{C}^{Q\times 1}$ and $\gamma$ is a user
chosen parameter controlling the sparsity of the solution, here being selected
as

$\gamma=\alpha\max_{i}|\mathbf{a}_{i}^{H}\mathbf{x}|,$ (11)

where $\mathbf{a}_{i}$, $i=1,2,\ldots,Q$ are columns of the $\mathbf{D}$. The
constant $\alpha$ is a user-chosen parameter allowing the performance of each
dictionary to be evaluated (as in the simulation results) by varying values of
this parameter within the range $0<\alpha\leq 1$, with
$\gamma_{max}=\alpha\max_{i}|\mathbf{a}_{i}^{H}\mathbf{x}|$ denoting the
smallest tuning parameter value that gives a solution with the coefficients of
zero [18, 24].

###  III-C MCM Estimation

Proceeding, we estimate the mutual coupling matrix coefficients by utilizing
the complex symmetric circular Toeplitz structure of the MCM in a UCA. Suppose
$\mathbf{R}_{x}$ denotes the array output covariance matrix and let
$\lambda_{n}$ and $\mathbf{e}_{n}$ denote its eigenvalues and the
corresponding eigenvectors, respectively. Then,

$\mathbf{R}_{x}=\sum_{n=1}^{N}\lambda_{n}\mathbf{e}_{n}\mathbf{e}_{n}^{H}=\mathbf{E}_{s}\Lambda_{s}\mathbf{E}_{s}^{H}+\mathbf{E}_{n}\Lambda_{n}\mathbf{E}_{n}^{H},$
(12)

where
$\mathbf{E}_{s}=[\mathbf{e}_{1},\mathbf{e}_{2},\ldots,\mathbf{e}_{K}]\in\mathbb{C}^{N\times
K}$ denotes the signal subspace containing the principal eigenvectors
corresponding to the $K$ maximum eigenvalues, with the noise subspace
$\mathbf{E}_{n}$ contains the remaining $N-K$ eigenvectors.

The signal subspace spans the same space as the array manifold matrix, i.e.,
$\text{span}\\{\mathbf{E}_{s}\\}\perp\text{span}\\{\mathbf{CA(\phi,\theta)}\\}$
. Moreover, the signal subspace and noise subspace are orthogonal such that
$\text{span}\\{\mathbf{E}_{s}\\}\perp\text{span}\\{\mathbf{E}_{n}\\}$,
implying that
$\text{span}\\{\mathbf{E}_{s}\\}\perp\text{span}\\{\mathbf{E}_{n}\\}$,
Implying that
$\text{span}\\{\mathbf{E}_{n}\\}\perp\text{span}\\{\mathbf{CA(\phi,\theta)}\\}$
.
By supposing that we have DOAs (as obtained in Subsection III-B), the MCM
could be estimated as the matrix $\mathbf{C}$ minimizing:

$J=\left\|\mathbf{E}_{n}^{H}\mathbf{C}\mathbf{A}(\phi,\theta)\right\|_{F}^{2}=\sum_{k=1}^{K}\left\|\mathbf{E}_{n}^{H}\mathbf{C}\mathbf{a}(\phi_{k},\theta_{k})\right\|^{2}$
(13)

where $\|\cdot\|_{F}$ denotes the Frobenius norm.

By using the circular symmetry of the mutual coupling matrix, one may instead
estimate the vector which constructs the mutual coupling matrix, i.e., by
rewriting
$\mathbf{C}\mathbf{a}(\phi,\theta)=\mathbf{F}\\{\mathbf{a}(\phi,\theta)\\}\mathbf{c}\in\mathbb{C}^{L\times
1}$ where $\mathbf{F}$ denotes the resulting transform matrix, and
$\mathbf{c}=[c_{1},c_{2},\ldots,c_{L}]^{T}\in\mathbb{C}^{L}$. By extending the
relation of the mutual coupling introduced in [21, 22] to 2-D and exploiting
the symmetric circulant property of $\mathbf{C}$, the matrix
$\mathbf{F}\\{\mathbf{a}(\phi,\theta)\\}$ may be defined as:

$\mathbf{F}\\{\mathbf{a}(\phi,\theta)\\}=\mathbf{F}_{1}\\{\mathbf{a}(\phi,\theta)\\}+\mathbf{F}_{2}\\{\mathbf{a}(\phi,\theta)\\}+\mathbf{F}_{3}\\{\mathbf{a}(\phi,\theta)\\}+\mathbf{F}_{4}\\{\mathbf{a}(\phi,\theta)\\}$
(14)

$\displaystyle\left[\mathbf{F}_{1}\\{\mathbf{a}(\phi,\theta)\\}\right]_{i,j}$
$\displaystyle=\begin{cases}\left[\mathbf{a}(\phi,\theta)\right]_{i+j-1,1}&\text{if
}i+j\leq N+1,\\\ 0&\text{otherwise},\end{cases}$ (15)
$\displaystyle\left[\mathbf{F}_{2}\\{\mathbf{a}(\phi,\theta)\\}\right]_{i,j}$
$\displaystyle=\begin{cases}\left[\mathbf{a}(\phi,\theta)\right]_{j-i+1,1}&\text{if
}i\geq 2,\\\ 0&\text{otherwise},\end{cases}$
$\displaystyle\left[\mathbf{F}_{3}\\{\mathbf{a}(\phi,\theta)\\}\right]_{i,j}$
$\displaystyle=\begin{cases}\left[\mathbf{a}(\phi,\theta)\right]_{N+i+j-1,1}&\text{if
}i<j\leq p,\\\ 0&\text{otherwise},\end{cases}$
$\displaystyle\left[\mathbf{F}_{4}\\{\mathbf{a}(\phi,\theta)\\}\right]_{i,j}$
$\displaystyle=\begin{cases}\left[\mathbf{a}(\phi,\theta)\right]_{i+j-N-1,1}&\text{if
}2\leq j\leq p,i+j\geq N+2,\\\ 0&\text{otherwise}.\end{cases}$

where $i=1,2,\ldots,N$, $j=1,2,\ldots,L$, and
$p=\left\lfloor\frac{N+1}{2}\right\rfloor$. Applying this transformation to
(13) yields

$J=\sum_{k=1}^{K}\mathbf{c}^{H}\mathbf{F}^{H}\\{\mathbf{a}(\phi_{k},\theta_{k})\\}\mathbf{E}_{n}\mathbf{E}_{n}^{H}\mathbf{a}(\phi_{k},\theta_{k})\mathbf{F}\\{a(\phi_{k},\theta_{k})\\}\mathbf{c}=\mathbf{c}^{H}\mathbf{U}(\phi,\theta)\mathbf{c},$
(16)

where

$\mathbf{U}(\phi,\theta)=\sum_{k=1}^{K}\mathbf{F}^{H}\\{\mathbf{a}(\phi_{k},\theta_{k})\\}\mathbf{E}_{n}\mathbf{E}_{n}^{H}\mathbf{F}\\{\mathbf{a}(\phi_{k},\theta_{k})\\}^{H}\in\mathbb{C}^{L\times
L},$ (17)

By imposing the constraint that $c_{1}=1$, (16) may be used as the cost
function to be minimized in order to obtain the desired mutual coupling
vector, i.e.,

$\mathbf{c}=\mathbf{e}_{\min}\left[\mathbf{U}(\phi,\theta)\right],$ (18)

where $\mathbf{e}_{\min}[\cdot]$ denotes the eigenvector corresponding to the
smallest eigenvalue.

##  IV Simulation Results

In this section, some simulations will be carried out to investigate the
performance of the proposed method. We consider a UCA with $N=15$ sensors, the
radius of the UCA being $r=\lambda$, where $\lambda$ denotes the center
wavelength of the narrow-band signals. We further assume three far-field
sources from angle of arrivals $(18.3^{\circ},243.4^{\circ})$,
$(83.6^{\circ},60^{\circ})$, and $(73.9^{\circ},357.8^{\circ})$ for which
$T=200$ snapshots are measured. The signals and the additive sensor noises are
stationary, mutually uncorrelated, and zero mean. The signals are modeled as
circularly symmetric Gaussian processes with variance $\sigma_{i}^{2}$ for the
$i$-th process. The corrupting noise is assumed to be i.i.d. white circularly
symmetric Gaussian process with variance $\sigma^{2}$, such that the input SNR
of the $i$-th signal is $10\log(\sigma_{i}^{2}/\sigma^{2})$. The coupling
coefficient is set to be $c_{2}=0.79+j0.432$ and $c_{3}=0.35+j0.16$.

Fig. 3 presents the probability of correctly estimating the number of angle of
arrivals using an integrated wideband dictionary. Suppose $B^{\theta}_{i}$ and
$B^{\phi}_{i}$ denote the number of bands over elevation and azimuth,
respectively. Also, we have $i=1,2,3$ for each stage of the dictionary. We
here consider for the zooming steps of the elevation $B^{\theta}_{1}=30$,
$B^{\theta}_{2}=10$, $B^{\theta}_{3}=3$, and for the azimuth
$B^{\phi}_{1}=120$, $B^{\phi}_{2}=10$, and $B^{\phi}_{3}=3$ (the grey curve),
and $B^{\theta}_{1}=30$, $B^{\theta}_{2}=5$, $B^{\theta}_{3}=3$,
$B^{\phi}_{1}=120$, $B^{\phi}_{2}=5$, and $B^{\phi}_{3}=3$ (the brown curve).
In this figure, the correct model order is given for varying values of
$\alpha$ for (11). As is clear from Fig. 3, using $\alpha\leq 0.4$ for the
three-stage dictionary (the grey curve) is a suitable region for utilizing the
proposed method.

Figure 3: Choosing suitable $\alpha$ for equation (11).

To evaluate the precision of DOA estimation methods, the average root-mean-
square errors (RMSEs) of the elevation and the azimuth estimation are used as
follows:

$\text{RMSE}=\sqrt{\frac{1}{KM}\sum_{k=1}^{K}\sum_{m=1}^{M}\left[(\theta_{k}^{m}-\theta_{k})^{2}+(\phi_{k}^{m}-\phi_{k})^{2}\right]},$
(19)

for $M=500$ number of Monte Carlo simulations, where
$(\theta_{k}^{m},\phi_{k}^{m})$ denotes the estimate of
$(\theta_{k},\phi_{k})$ in the $m$-th Monte Carlo run (out of $M$ runs). Fig.
4 shows the resulting RMSE when the SNR varies from 0 dB to 20 dB. The results
are compared with the LASSO estimates resulting when the dictionary is
composed of discrete prolate spheroidal sequences (DPSS) [25, 26], and to the
2-D extension of the method in [20], here termed the Iterative-MUSIC
estimator. This Iterative-MUSIC is found to have the highest computational
complexity while offering the worst performance of the discussed methods.

Figure 4: The RMSE of angles of estimation as a function of the SNR.

Fig. 5 shows the RMSE of coupling coefficients, being defined as

Figure 5: The RMSE of the mutual coupling coefficients estimation as a
function of the SNR.

$\text{RMSE}_{\text{coupling
coeff.}}=\frac{\sqrt{\sum_{m=1}^{M}\|\mathbf{\hat{c}}^{m}-\mathbf{c}\|^{2}}}{\|\mathbf{c}\|}\times
100\%,$ (20)

where $\mathbf{\hat{c}}^{m}$ is the estimation of coupling coefficient vector,
$\mathbf{c}$, in the $m$-th Monte Carlo experiment. As shown in Fig. 5, the
coupling coefficients estimation of the proposed method is very effective and
also has a very high precision.

Figure 6: The RMSE of the angle estimates as a function of the number of the
measured snapshots.

In order to compare the performance in terms of the number of snapshots, we
fix the SNR at $5$ dB, $M=200$ and vary the snapshot number. The results are
shown in Fig. 6. The other parameters are set as before. As seen in Fig.
LABEL:fig:fig3 to Fig. 6, the proposed method offers preferable performance of
the three methods.

##  V Conclusion

In this work, we have introduced a self-calibration algorithm which is able to
jointly estimate the 2-D DOAs and the MCM for a UCA configuration. The DOA
estimation step is formed using a sparse representation framework in
combination with an integrated dictionary elements, which span bands of the
desired parameter space. The MCM estimation step exploits the circular
symmetry of the mutual coupling matrix resulting for UCAs. This method has the
advantages of allowing for an initial coarse gridding of the dictionary,
without the risk of suffering from off-grid effects, as well as not posing any
constraints on the coupling matrix, nor resulting in any ambiguity in the
resulting DOA estimates. Several simulation results confirm the superior
performance of the proposed method.

##  Appendix A Simplifying the 2-D integral of the dictionary elements in
order to allow for non-numerical computations

To simplify the 2-D integration in (6) and attain a formula for computing the
integration, one may use the Taylor series expansion for $e^{z}$ about $z=0$,
i.e.,

$e^{z}=\sum_{n=0}^{\infty}\frac{z^{n}}{n!}$ (21)

$\displaystyle\nu\triangleq$
$\displaystyle\int_{\phi_{i}}^{\phi_{i+1}}\int_{\theta_{i}}^{\theta_{i+1}}e^{j\alpha
r\sin\theta\cos(\phi-\lambda)}d\theta d\phi$ $\displaystyle=$
$\displaystyle\int_{\phi_{i}}^{\phi_{i+1}}\int_{\theta_{i}}^{\theta_{i+1}}\sum_{n=0}^{\infty}\frac{\left(\alpha\cos(\phi-\lambda)\right)^{n}}{n!}\sin^{n}\theta
d\theta d\phi$ $\displaystyle=$
$\displaystyle\int_{\phi_{i}}^{\phi_{i+1}}\sum_{n=0}^{\infty}\left(\frac{\alpha^{n}}{n!}\cos^{n}(\phi-\lambda)\int_{\theta_{i}}^{\theta_{i+1}}\sin^{n}\theta
d\theta\right)d\phi$ $\displaystyle=$
$\displaystyle\sum_{n=0}^{\infty}\frac{\alpha^{n}}{n!}\int_{\phi_{i}}^{\phi_{i+1}}\cos^{n}(\phi-\lambda)d\phi\int_{\theta_{i}}^{\theta_{i+1}}\sin^{n}\theta
d\theta$ $\displaystyle=$
$\displaystyle\sum_{n=0}^{\infty}\frac{\alpha^{n}}{n!}\int_{-\lambda_{n}}^{\phi_{i+1}-\lambda_{n}}\cos^{n}\phi
d\phi\int_{\theta_{i}}^{\theta_{i+1}}\sin^{n}\theta d\theta$ (22)

Computing the first and the second terms yields

$\nu=\left(\int_{\phi_{i}}^{\phi_{i+1}-\lambda_{n}}\phi
d\phi\right)\left(\int_{\theta_{i}}^{\theta_{i+1}}d\theta\right)+\alpha\left(\sin\int_{\phi_{i}}^{\phi_{i+1}-\lambda_{n}}d\phi\right)\left(-\cos\int_{\theta_{i}}^{\theta_{i+1}}\theta
d\theta\right)+\sum_{n=2}^{\infty}\frac{\alpha^{n}}{n!}\int_{\phi_{i}}^{\phi_{i+1}-\lambda_{n}}\cos^{n}\phi
d\phi\int_{\theta_{i}}^{\theta_{i+1}}\sin^{n}\theta d\theta$ (23)

Integration by parts for the third term of the above equation yields

$\int\sin^{n}\theta
d\theta=-\frac{1}{n}\sin^{n-1}\theta\cos\theta+\frac{n-1}{n}\int\sin^{n-2}\theta
d\theta$ (24)

and

$\int\cos^{n}\phi
d\phi=\frac{1}{n}\cos^{n-1}\phi\sin\phi+\frac{n-1}{n}\int\cos^{n-2}\phi d\phi$
(25)

One may thus avoid a numerical integration by computing the recursive factors
in the terms. To do this, the number of terms (for this work, computing about
$\sim 10$ terms suffices and gives the same precision as the built-in
integral2 function in Matlab) in this series could be computed to get the
required precision.

## Acknowledgment

The first author would like to thank Dr. Andreas Jakobsson (Dept. of
Mathematical Statistics, Lund University) because of his insightful comments
and discussions.

## References

  * [1] H. Krim and M. Viberg, “Two decades of array signal processing research: the parametric approach,” _IEEE Signal Process Mag_ , vol. 13, no. 4, pp. 67–94, 1996.
  * [2] H. C. So, _Source localization: Algorithms and analysis_. New York: Wiley-IEEE Press, 2011.
  * [3] Y. Liu, J. Chai, Y. Zhang, Z. Liu, M. Jin, and T. Qiu, “Low-complexity neural network based doa estimation for wideband signals in massive mimo systems,” _AEU- Int J Electron Commun_ , vol. 138, p. 153853, 2021.
  * [4] Z. Zheng and C. Yang, “Direction-of-arrival estimation of coherent signals under direction-dependent mutual coupling,” _IEEE Commun Lett_ , vol. 25, no. 1, pp. 147–151, 2020.
  * [5] F. Mei, H. Xu, W. Cui, B. Ba, and Y. Wang, “A transformed coprime array with reduced mutual coupling for doa estimation of non-circular signals,” _IEEE Access_ , vol. 9, pp. 125 984–125 998, 2021.
  * [6] Z. Zheng, C. Yang, W. Q. Wang, and H. C. So, “Robust doa estimation against mutual coupling with nested array,” _IEEE Signal Process Lett_ , vol. 27, pp. 1360–1364, 2020.
  * [7] C. M. See, “Sensor array calibration in the presence of mutual coupling and unknown sensor gains and phases,” _Electronics Letters_ , vol. 30, no. 5, pp. 373–374, 1994.
  * [8] S. Liu, L. Yang, and S. Yang, “Robust joint calibration of mutual coupling and channel gain/phase inconsistency for uniform circular array,” _IEEE Anten Wirel Propag_ , vol. 15, pp. 1191–1195, 2015.
  * [9] W. Hu and Q. Wang, “Doa estimation for uca in the presence of mutual coupling via error model equivalence,” _IEEE Wirel Comm Letter_ , vol. 9, no. 1, pp. 121–124, 2019.
  * [10] R. Goossens and H. Rogier, “Direction-of-arrival and polarization estimation with uniform circular arrays in the presence of mutual coupling,” _AEU- Int J Electron Commun_ , vol. 62, no. 3, pp. 199–206, 2008.
  * [11] K. Wang, J. Yi, F. Cheng, Y. Rao, and X. Wan, “Array errors and antenna element patterns calibration based on uniform circular array,” _IEEE Anten Wirel Propag_ , vol. 20, no. 6, pp. 1063–1067, 2021.
  * [12] G. J. Jiang, Y. N. Dong, X. P. Mao, and Y. T. Liu, “Improved 2d direction of arrival estimation with a small number of elements in uca in the presence of mutual coupling,” _AEU- Int J Electron Commun_ , vol. 71, pp. 131–138, 2017\.
  * [13] J. Zhao, J. Liu, F. Gao, W. Jia, and W. Zhang, “Gridless compressed sensing based channel estimation for uav wideband communications with beam squint,” _IEEE Trans Veh Technol_ , vol. 70, no. 10, pp. 10 265–10 277, 2021.
  * [14] P. Stoica and P. Babu, “Sparse estimation of spectral lines: Grid selection problems and their solutions,” _IEEE Trans Signal Process_ , vol. 60, no. 2, pp. 962–967, 2011.
  * [15] M. Wagner, Y. Park, and P. Gerstoft, “Gridless doa estimation and root-music for non-uniform linear arrays,” _IEEE Trans Signal Process_ , vol. 69, pp. 2144–2157, 2021.
  * [16] Y. Chi and Y. Chen, “Compressive two-dimensional harmonic retrieval via atomic norm minimization,” _IEEE Trans. Signal Process._ , vol. 63, no. 4, pp. 1030–1042, 2014.
  * [17] Z. Yang and L. Xie, “Enhancing sparsity and resolution via reweighted atomic norm minimization,” _IEEE Trans. Signal Process._ , vol. 64, no. 4, pp. 995–1006, 2015.
  * [18] M. Butsenko, J. Swärd, and A. Jakobsson, “Estimating sparse signals using integrated wideband dictionaries,” _IEEE Trans. Signal Process._ , vol. 66, no. 16, pp. 4170–4181, 2018.
  * [19] D. Y. Gao, B. H. Wang, and Y. Guo, “Comments on "blind calibration and doa estimation with uniform circular arrays in the presence of mutual coupling",” _IEEE Anten. Wirel. Propag. Lett._ , vol. 5, pp. 566–568, 2006\.
  * [20] J. Xie, Z. S. He, and H. Y. Li, “A fast doa estimation algorithm for uniform circular arrays in the presence of unknown mutual coupling,” _Progress In Electromagnetics Research C_ , vol. 21, pp. 257–271, 2011.
  * [21] M. Wang, X. Ma, S. Yan, and C. Hao, “An autocalibration algorithm for uniform circular array with unknown mutual coupling,” _IEEE Anten. Wirel. Propag. Lett._ , vol. 15, pp. 12–15, 2015.
  * [22] B. Friedlander and A. J. Weiss, “Direction finding in the presence of mutual coupling,” _IEEE Trans. Anten. Propag._ , vol. 39, no. 3, pp. 273–284, 1991\.
  * [23] M. Wax and T. Kailath, “Detection of signals by information theoretic criteria,” _IEEE Trans. Acoust. Speech Signal Process._ , vol. 33, no. 2, pp. 387–392, 1985.
  * [24] R. Tibshirani, J. Bien, J. Friedman, T. Hastie, N. Simon, J. Taylor, and R. J. Tibshirani, “Strong rules for discarding predictors in lasso-type problems,” _Journal of the Royal Statistical Society: Series B (Statistical Methodology)_ , vol. 74, no. 2, pp. 245–266, 2012.
  * [25] M. A. Davenport and M. B. Wakin, “Compressive sensing of analog signals using discrete prolate spheroidal sequences,” _Applied and Computational Harmonic Analysis_ , vol. 33, no. 3, pp. 438–472, 2012.
  * [26] D. Slepian, “Prolate spheroidal wave functions, fourier analysis, and uncertainty—v: The discrete case,” _Bell System Technical Journal_ , vol. 57, no. 5, pp. 1371–1430, 1978.