## MUSIC Super-Resolution DOA Estimation

MUltiple SIgnal Classification (MUSIC) is a high-resolution direction-finding algorithm based on the eigenvalue decomposition of the sensor covariance matrix observed at an array. MUSIC belongs to the family of subspace-based direction-finding algorithms.

### Signal Model

The signal model relates the received sensor data to the signals emitted by the source. Assume that there are D uncorrelated or partially correlated signal sources, sd(t). The sensor data, xm(t), consists of the signals, as received at the array, together with noise, nm(t). A sensor data snapshot is the sensor data vector received at the M elements of an array at a single time t.

`$\begin{array}{l}x\left(t\right)=As\left(t\right)+n\left(t\right)\\ s\left(t\right)=\left[{s}_{1}\left(t\right),{s}_{2}\left(t\right),\dots ,{s}_{D}\left(t\right)\right){\right]}^{\prime }\\ A=\left[a\left({\theta }_{1}\right)|a\left({\theta }_{2}\right)|\dots |a\left({\theta }_{D}\right)\right]\end{array}$`

• x(t) is an M-by-1 vector of received snapshot of sensor data which consist of signals and additive noise.

• A is an M-by-D matrix containing the arrival vectors. An arrival vector consists of the relative phase shifts at the array elements of the plane wave from one source. Each column of A represents the arrival vector from one of the sources and depends on the direction of arrival, θd. θd is the direction of arrival angle for the dth source and can represents either the broadside angle for linear arrays or the azimuth and elevation angle for planar or 3D arrays.

• s(t) is a D-by-1 vector of source signal values from D sources.

• n(t) is an M-by-1 vector of sensor noise values.

An important quantity in any subspace method is the sensor covariance matrix,Rx, derived from the received signal data. When the signals are uncorrelated with the noise, the sensor covariance matrix has two components, the signal covariance matrix and the noise covariance matrix.

`${R}_{x}=E\left\{x{x}^{H}\right\}=A{R}_{s}{A}^{H}+{\sigma }_{n}^{2}I$`

where Rs is the source covariance matrix. The diagonal elements of the source covariance matrix represent source power and the off-diagonal elements represent source correlations.

`${R}_{s}=E\left\{s{s}^{H}\right\}$`

For uncorrelated sources or even partially correlated sources, Rs is a positive-definite Hermitian matrix and has full rank, D, equal to the number of sources.

The signal covariance matrix, ARsAH, is an M-by-M matrix, also with rank D < M.

An assumption of the MUSIC algorithm is that the noise powers are equal at all sensors and uncorrelated between sensors. With this assumption, the noise covariance matrix becomes an M-by-M diagonal matrix with equal values along the diagonal.

Because the true sensor covariance matrix is not known, MUSIC estimates the sensor covariance matrix, Rx, from the sample sensor covariance matrix. The sample sensor covariance matrix is an average of multiple snapshots of the sensor data

`${R}_{x}=\frac{1}{T}\sum _{k=1}^{T}x\left(t\right)x{\left(t\right)}^{H},$`

where T is the number of snapshots.

### Signal and Noise Subspaces

Because ARsAH has rank D, it has D positive real eigenvalues and M – D zero eigenvalues. The eigenvectors corresponding to the positive eigenvalues span the signal subspace, Us= [v1,...,vD]. The eigenvectors corresponding to the zero eigenvalues are orthogonal to the signal space and span the null subspace, Un= [uD+1,...,uN]. The arrival vectors also belong to the signal subspace, but they are eigenvectors. Eigenvectors of the null subspace are orthogonal to the eigenvectors of the signal subspace. Null-subspace eigenvectors, ui, satisfy this equation:

`$A{R}_{s}{A}^{H}{u}_{i}=0⇒{u}^{H}A{R}_{s}{A}^{H}{u}_{i}=0⇒{\left({A}^{H}{u}_{i}\right)}^{H}{R}_{s}\left({A}^{H}{u}_{i}\right)=0⇒{A}^{H}{u}_{i}=0$`

Therefore the arrival vectors are orthogonal to the null subspace.

When noise is added, the eigenvectors of the sensor covariance matrix with noise present are the same as the noise-free sensor covariance matrix. The eigenvalues increase by the noise power. Let vi be one of the original noise-free signal space eigenvectors. Then

`${R}_{x}{v}_{i}=A{R}_{s}{A}^{H}{v}_{i}+{\sigma }_{0}^{2}I{v}_{i}=\left({\lambda }_{i}+{\sigma }_{0}^{2}\right){v}_{i}$`

shows that the signal space eigenvalues increase by σ02.

The null subspace eigenvectors are also eigenvectors of Rx. Let ui be one of the null eigenvectors. Then

`${R}_{x}{u}_{i}=A{R}_{s}{A}^{H}{u}_{i}+{\sigma }_{0}^{2}I{u}_{i}={\sigma }_{0}^{2}{u}_{i}$`

with eigenvalues of σ02 instead of zero. The null subspace becomes the noise subspace.

MUSIC works by searching for all arrival vectors that are orthogonal to the noise subspace. To do the search, MUSIC constructs an arrival-angle-dependent power expression, called the MUSIC pseudospectrum:

`${P}_{MUSIC}\left(\stackrel{\to }{\varphi }\right)=\frac{1}{{a}^{H}\left(\stackrel{\to }{\varphi }\right){U}_{n}{U}_{n}^{H}a\left(\stackrel{\to }{\varphi }\right)}$`

When an arrival vector is orthogonal to the noise subspace, the peaks of the pseudospectrum are infinite. In practice, because there is noise, and because the true covariance matrix is estimated by the sampled covariance matrix, the arrival vectors are never exactly orthogonal to the noise subspace. Then, the angles at which PMUSIC has finite peaks are the desired directions of arrival. Because the pseudospectrum can have more peaks than there are sources, the algorithm requires that you specify the number of sources, D, as a parameter. Then the algorithm picks the D largest peaks. For a uniform linear array (ULA), the search space is a one-dimensional grid of broadside angles. For planar and 3D arrays, the search space is a two-dimensional grid of azimuth and elevation angles.

### Root-MUSIC

For a ULA, the denominator in the pseudospectrum is a polynomial in ${e}^{ikd\mathrm{cos}\phi }$, but can also be considered a polynomial in the complex plane. In this cases, you can use root-finding methods to solve for the roots, zi. These roots do not necessarily lie on the unit circle. However, Root-MUSIC assumes that the D roots closest to the unit circle correspond to the true source directions. Then you can compute the source directions from the phase of the complex roots.

### Spatial Smoothing of Correlated Sources

When some of the D source signals are correlated, Rs is rank deficient, meaning that it has fewer than D nonzero eigenvalues. Therefore, the number of zero eigenvalues of ARsAH exceeds the number, M – D, of zero eigenvalues for the uncorrelated source case. MUSIC performance degrades when signals are correlated, as occurs in a multipath propagation environment. A way to compensate for correlation is to use spatial smoothing.

Spatial smoothing takes advantage of the translation properties of a uniform array. Consider two correlated signals arriving at an L-element ULA. The source covariance matrix, Rs is a singular 2-by-2 matrix. The arrival vector matrix is an L-by-2 matrix

`$\begin{array}{l}{A}_{1}=\left[\begin{array}{c}1\\ {e}^{ikd\mathrm{cos}{\phi }_{1}}\\ ⋮\\ {e}^{i\left(L-1\right)kd\mathrm{cos}{\phi }_{1}}\end{array}\begin{array}{c}1\\ {e}^{ikd\mathrm{cos}{\phi }_{2}}\\ ⋮\\ {e}^{i\left(L-1\right)kd\mathrm{cos}{\phi }_{2}}\end{array}\right]=\left[a\left({\phi }_{1}\right)|a\left({\phi }_{2}\right)\right]\\ \end{array}$`

for signals arriving from the broadside angles φ1 and φ2. The quantity k is the signal wave number. a(φ) represents an arrival vector at the angle φ.

You can create a second array by translating the first array along its axis by one element distance, d. The arrival matrix for the second array is

`${A}_{2}=\left[\begin{array}{c}{e}^{ikd\mathrm{cos}{\phi }_{1}}\\ {e}^{i2kd\mathrm{cos}{\phi }_{1}}\\ ⋮\\ {e}^{iLkd\mathrm{cos}{\phi }_{1}}\end{array}\begin{array}{c}{e}^{ikd\mathrm{cos}{\phi }_{2}}\\ {e}^{i2kd\mathrm{cos}{\phi }_{2}}\\ ⋮\\ {e}^{iLkd\mathrm{cos}{\phi }_{2}}\end{array}\right]=\left[{e}^{ikd\mathrm{cos}{\phi }_{1}}a\left({\phi }_{1}\right)|{e}^{ikd\mathrm{cos}{\phi }_{2}}a\left({\phi }_{2}\right)\right]$`

where the arrival vectors are equal to the original arrival vectors but multiplied by a direction-dependent phase shift. When you translate the original array J –1 more times, you get J copies of the array. If you form a single array from all these copies, then the length of the single array is M = L + (J – 1).

In practice, you start with an M-element array and form J overlapping subarrays. The number of elements in each subarray is L = M – J + 1. The following figure shows the relationship between the overall length of the array, M, the number of subarrays, J, and the length of each subarray, L.

For the pth subarray, the source signal arrival matrix is

`$\begin{array}{c}{A}_{p}=\left[{e}^{ik\left(p-1\right)d\mathrm{cos}{\phi }_{1}}a\left({\phi }_{1}\right)|{e}^{ik\left(p-1\right)d\mathrm{cos}{\phi }_{2}}a\left({\phi }_{2}\right)\right]\\ =\left[a\left({\phi }_{1}\right)|a\left({\phi }_{2}\right)\right]\left[\begin{array}{cc}{e}^{ik\left(p-1\right)d\mathrm{cos}{\phi }_{1}}& 0\\ 0& {e}^{ik\left(p-1\right)d\mathrm{cos}{\phi }_{2}}\end{array}\right]={A}_{1}{P}^{p-1}\\ P=\left[\begin{array}{cc}{e}^{ikd\mathrm{cos}{\phi }_{1}}& 0\\ 0& {e}^{ikd\mathrm{cos}{\phi }_{2}}\end{array}\right].\end{array}$`

The original arrival vector matrix is postmultiplied by a diagonal phase matrix.

The last step is averaging the signal covariance matrices over all J subarrays to form the averaged signal covariance matrix, Ravgs. The average signal covariance matrix depends on the smoothed source covariance matrix, Rsmooth.

`$\begin{array}{l}{R}_{s}^{avg}={A}_{1}\left(\frac{1}{J}\sum _{p=1}^{J}{P}^{p-1}{R}_{s}{\left({P}^{p-1}\right)}^{H}\right){A}_{1}^{H}={A}_{1}{R}^{smooth}{A}_{1}^{H}\\ {R}^{smooth}=\frac{1}{J}\sum _{p=1}^{J}{P}^{p-1}{R}_{s}{\left({P}^{p-1}\right)}^{H}.\end{array}$`

You can show that the diagonal elements of the smoothed source covariance matrix are the same as the diagonal elements of the original source covariance matrix.

`${R}_{ii}^{smooth}=\frac{1}{J}\sum _{p=1}^{J}{\left({P}^{p-1}\right)}_{im}{\left({R}_{s}\right)}_{mn}{\left({P}^{p-1}\right)}_{ni}{}^{H}=\frac{1}{J}\sum _{p=1}^{J}{R}_{s}={\left({R}_{s}\right)}_{ii}$`

However, the off-diagonal elements are reduced. The reduction factor is the beam pattern of a J-element array.

`${R}_{ij}^{smooth}=\frac{1}{J}\sum _{p=1}^{J}{e}^{ikd\left(p-1\right)\left(\mathrm{cos}{\phi }_{1}-\mathrm{cos}{\phi }_{2}\right)}{\left({R}_{s}\right)}_{ij}=\frac{1}{J}\frac{\mathrm{sin}\left(kdJ\left(\mathrm{cos}{\phi }_{1}-\mathrm{cos}{\phi }_{2}\right)\right)}{\mathrm{sin}\left(kd\left(\mathrm{cos}{\phi }_{1}-\mathrm{cos}{\phi }_{2}\right)\right)}{\left({R}_{s}\right)}_{ij}$`

In summary, you can reduce the degrading effect of source correlation by forming subarrays and using the smoothed covariance matrix as input to the MUSIC algorithm. Because of the beam pattern, larger angular separation of sources leads to reduced correlation.

Spatial smoothing for linear arrays is easily extended to 2D and 3D uniform arrays.