Generalized correlation based Imaging for satellites

In this paper we present a generalization of the migration method for correlations, which is motivated by the structure of the data. We first define a correlation data structure, which compensates for the Doppler effect on the measurements. We then show that out of this data structure a natural migration scheme arises, which translates the data to two points rather than one. This is natural when considering data correlations, and a similar migration scheme was recently used in [4]. The result of this migration is not an image but rather a two-point interference pattern. We show that there are several possible ways to derive an image from the two-point interference pattern, with one in particular that we call the rank-1 image. This rank-1 image is the result of taking the first eigenvector of the interference pattern, and it provides superior resolution. We demonstrate its superiority both in terms of favorable mathematical properties, such as having an optimization interpretation, as well as by providing analytical explanations of its performance. This is the main result of the paper.

Numerical simulations confirm that the resolution of the rank-1 image is better than that of the previously used correlation based imaging function, and achieves resolution comparable to that of linear migration, while retaining the benefits of imaging with correlations rather than with the data directly. The method is stable with respect to additive noise and can be applied to many other imaging problems.

The rest of the paper is as follows. In Section 2 we briefly review the model for our data (further expanded in Appendix A), and introduce the cross-correlation data structure. In Section 3 we recast the problem in matrix form and introduce the interference pattern, and the possible extensions of KM to cross-correlation data, particularly the rank-1 image. In Section 4 we present numerical simulations which confirm the superior performance of the rank-1 image. We motivate the results by further simulations in Section 5, and analysis in Appendix B and Appendix C. We conclude with a summary and conclusions in Section 6.

We want to image a cluster of objects, moving in Keplerian low earth orbit at a high velocity $\mathbf{v}_{L}$, so that without rotation for all the objects we have

The data is the collection of signals recorded at ground based or at low elevation receivers, with positions $\mathbf{x}_{\mathbf{R}}$. Successive pulses are emitted at a slow time $s$ by a source located at $\mathbf{x}_{\mathbf{E}}$ on the ground, as illustrated in Figure 1.

The receiver emits a series of pulses $f(t)=\cos(\omega_{o}t)e^{-B^{2}t^{2}/2}$, at slow time intervals of $\Delta s$, with a total aperture size $S$, such that the recorded signal at the receiver location $\mathbf{x}_{\mathbf{R}}$ due to a pulse, $f(s+t)$, emitted at slow time $s\in[-S/2,S/2]$, is

The derivation of (2.2) is in Appendix A. With $c_{0}$ the speed of light, $\gamma_{\mathbf{R}}$ is the Doppler scale correction factor and $t_{\mathbf{R}}$ is the signal travel time, which to first order in $|\mathbf{v}_{\mathbf{T}}|/c_{0}$ are given by

The formula (2.2) for the recorded data will be used in the theoretical resolution analysis.

An important consideration in the data recording process is the sampling frequency, which in currently available microwave acquisition systems is limited to about 2 GHz. Signals with higher carrier frequency need to be downramped before sampling, to avoid information loss. This is done with a heterodyne oscillator that shifts down the carrier frequency followed by low pass filtering. When using asynchronous emitters, the central frequency and bandwidth might not be known accurately. In that case the frequency of the oscillator and the sampling frequency need to be chosen so that loss of information is minimal.

We use the cross correlations of the signals received at different receivers to image. In [8] the imaging function was defined via scaling and translating the signals at the receivers by different Doppler factors $\gamma_{\mathbf{R}}$ for every point $\mathbf{x},\mathbf{v}$ in the imaging domain

This imaging function is a generalization of the regular Kirchhoff migration, which back propagates the recorded signals by the travel time to a specific search point in the image domain, but now at the level of cross correlations. The idea is that when migrating and correlating the data from the different receivers, the peaks of the signals sum coherently if a reflector exists at a point $\mathbf{x}$, with velocity $\mathbf{v}$, resulting in a peak for $\mathcal{I}^{CC}(\mathbf{x},\mathbf{v})$ at that point. Since the rapid velocity affects the signal via $\gamma_{\mathbf{R}}$, the Doppler factor must be used when migrating the data. At low speeds and without the Doppler correction this imaging function reduces to the classical correlation based migration [13].

We assume here that the locations of the receivers $\mathbf{x}_{\mathbf{R}},\mathbf{x}_{\mathbf{R}^{\prime}}$ are accurately known. The source location $\mathbf{x}_{\mathbf{E}}$ is also assumed to be known, but a lesser degree of accuracy is needed, since $\mathbf{x}_{\mathbf{E}}$ only affects $\gamma_{\mathbf{R}}$ (see [8]). The formulation (2.4) is less advantageous for two reasons:

1. 1.

   It is computationally ineffective to recompute the cross correlation for any point $\mathbf{x}$ in the imaging domain. It would be much more efficient to compute the cross correlations once for any triplet $(s_{j},\mathbf{R},\mathbf{R}^{\prime})$ as a function of $\tau$, the translation, and use that for imaging [13].

2. 2.

   Important features in the image manifest themselves in the raw cross correlation data. In order to use the data intelligently one must define a data structure, as an intermediate step, before applying the imaging function.

We next define the Doppler-scaled cross correlation function $C_{\mathbf{R}\mathbf{R}^{\prime}}(s,\tau)$, and show how it can be used in the imaging problem.

In general, (2.3) shows that $\gamma_{\mathbf{R}}$ is dependent on the target’s position and velocity. Having an accurate value for $\gamma_{\mathbf{R}}$ is extremely important in calculating cross correlations, as small changes in the frequency support would greatly affect the outcome of the correlation integral and hence the resolution of the image. However, we can use prior information at our disposal to limit the image region. For example, $\mathbf{v}_{\mathbf{T}}$ can be estimated from the fact that the objects are assumed to be in a Keplerian orbit, i.e., $|\mathbf{v}_{\mathbf{T}}|\approx\sqrt{\frac{2GM_{\text{Earth}}}{\mathbf{R}_{\mathbf{T}}}}$. Accurate range measurements based on the signal’s bandwidth can be used to estimate $\mathbf{R}_{\mathbf{T}}$ and subsequently $\mathbf{v}_{\mathbf{T}}$. Moreover, when forming an image we are usually interested in a limited search domain in space, based on auxiliary information. Of course $\gamma_{\mathbf{R}}$ is a continuous function of $\mathbf{x}_{\mathbf{T}}$, $\mathbf{v}_{\mathbf{T}}$, such that if $|\mathbf{x}_{\mathbf{T}}-\mathbf{x}_{0}|$, $|\mathbf{v}_{\mathbf{T}}-\mathbf{v}_{0}|$, are small enough then

In our numerical simulations over an image window of size $10$m$\times 10$m$\times 10$m, $\gamma_{\mathbf{R}}$ varies at most by $0.00001$%. For signals with our frequency support and time duration, $\gamma_{\mathbf{R}}$ can be very well approximated by a constant in the imaging domain.

Relying upon the weak variation of $\gamma_{\mathbf{R}}$, we use a rescaled signal to image over the image window, rather than rescaling the signal for every point in the image window as was done in [8](eq. 3.6). Specifically, we construct $C_{\mathbf{RR}^{\prime}}(s,\tau)$ by scaling the signal received by a reference Doppler factor,

Note that for a reflector with a trajectory $\mathbf{x}_{\mathbf{T}},\mathbf{v}_{\mathbf{T}}$, the received signal has the form

We define the cross correlation function as well,

Using (2.7), the imaging function (2.4) reduces to

where

Notice that $t^{\mathbf{x}}_{\mathbf{R}}$ is still dependent on $\gamma_{\mathbf{R}}(\mathbf{x},\mathbf{x}_{\mathbf{E}},\mathbf{v})$. The imaging function in (2.9) migrates the cross correlation data to a search point $\mathbf{x}$, based on the expected travel time difference between the signals recorded at both receivers $\mathbf{R},\mathbf{R}^{\prime}$. This is possible because we have adequately scaled the fields with respect to a specific search region, such that $C_{\mathbf{R}\mathbf{R}^{\prime}}(s,\tau)$ is only dependent on travel time differences. Note also that the recorded signals are time shifted by the sum of the travel times from the source to the center of the image window plus the travel time from the center of the image window to the receiver. This ensures that $C_{\mathbf{RR}^{\prime}}(s,\tau)$ is supported around $\tau=0$.

The imaging function (2.9) is an extension of the classical Kirchhoff Migration (KM) imaging function in the linear case, where the travel times are used to migrate the field

Thus, we can break down the imaging problem into two separate steps:

1. 1.

   In the first step, the trajectory and fast linear motion of the object is estimated and thereafter assumed known. The relevant imaging domain can then be taken to be small enough so that the approximation of the ratio (2.5) holds. We express the estimated motion of the center of the image window as $\mathbf{x}_{\mathbf{L}}(s),(\mathbf{x}_{\mathbf{L}},\mathbf{v}_{\mathbf{L}})$, where

   $$\mathbf{x}_{\mathbf{L}}(s)=\mathbf{x}_{\mathbf{L}}+s\mathbf{v}_{\mathbf{L}}.$$

   (2.12)

   Alternatively, we can use prior information about the trajectory to center the image window.

2. 2.

   We then substitute in (2.8) $\mathbf{x}_{\mathbf{L}}(s)$ for $\mathbf{x}_{0}$ and $\mathbf{v}_{\mathbf{L}}$ for $\mathbf{v}_{0}$ to compute the cross-correlation $C_{\mathbf{R}\mathbf{R}^{\prime}}(s,\tau)$.

Assume we have discretized the medium in a small image window relative to its moving center $\mathbf{x}_{\mathbf{L}}(s)$, as in (2.12), with grid points $\mathbf{y}_{k}$ and $\mathbf{y}_{k}=0$ corresponds to the center of the window. The unknown reflectivity is discretized by its values on this grid

The unknown reflectivity vector has dimension $K$, which is the number of pixels in the image window. Most of these reflectivities are zero because there are usually few relatively strong reflectors that can be imaged, that is, we can estimate their location in the image window and their strength. We assume the image window is small enough so that the Doppler term $\gamma_{\mathbf{R}}$ is constant over it, as explained in Section 2.

We saw in (2.7) that the signal recorded at receiver location $\mathbf{x}_{\mathbf{R}}$ can be written, after an appropriate scaling, as

where

In the frequency domain the recorded signal is

with

The phase $A_{\mathbf{R},k}(s,\omega)$ comes from the reduced travel time from the target to the receiver, relative to that of the image window center.

We assume the distance from the reflectors to the different receivers doesn’t vary greatly, hence we can approximate

i.e, we neglect the dependence of the amplitude factor on a specific receiver. By using this approximation we compromise the accuracy in which the amplitude of the reflector is retrieved. However, for most applications, the support and the relative reflectivities of the scatterers are of the greatest importance. Using this notation, and the fact that correlation in time is equivalent to multiplication in frequency, we get that

We introduce matrix notation that relates our model for reflectivities to the cross correlation data. Denote by $\boldsymbol{\rho}$ the unknown reflectivities in vector form

Denote $\mathbf{A}(s,\omega)$, our model for the migration matrix, has dimensions $N_{R}\times K$ and entries $A_{\mathbf{R},k}(s,\omega)$. This matrix $\mathbf{A}(s,\omega)$ acts on the reflectivities and returns data. Denote the recorded signal data as an $N_{R}$ vector vector $\widehat{\mathbf{u}}_{R}(s,\omega)$ whose entries are given by (3.15). The cross correlation data is also a matrix, of dimension $N_{R}\times N_{R}$, $\widehat{\mathbf{C}}(s,\omega)$ with entries $\widehat{C}_{\mathbf{R}\mathbf{R}^{\prime}}(s,\omega).$

Combining these, we have in matrix form the following model for the recorded signal data vector $\widehat{\mathbf{u}}_{\mathbf{R}}(s,\omega)$ and cross correlation data matrix $\widehat{\mathbf{C}}(s,\omega)$

We denote by $\mathbf{X}=\boldsymbol{\rho}\boldsymbol{\rho}^{T},X_{kk^{\prime}}=\rho_{k}\rho_{k^{\prime}}$, the outer product of reflectivities, then our model in matrix form is

The cross correlations depend on the reflectivities $\boldsymbol{\rho}$ through their outer product $\mathbf{X}=\boldsymbol{\rho}\boldsymbol{\rho}^{T}$. As a result there can be several different extensions of Kirchhoff migration to cross correlations, which we investigate in the next section.

Given the data $\widehat{\mathbf{u}}_{\mathbf{R}}(s,\omega)$ and the model $\widehat{\mathbf{u}}_{\mathbf{R}}(s,\omega)=\mathbf{A}(s,\omega)\boldsymbol{\rho}$, Kirchhoff migration of the data $\widehat{\mathbf{u}}_{\mathbf{R}}(s,\omega)$ is given by

Identifying $\tilde{\boldsymbol{\rho}}_{k}=\tilde{\rho}(\mathbf{y}_{k})$, (3.22) is a discretized version of the imaging function in (2.11), $\mathcal{I}^{KM}(\mathbf{y}_{k},\mathbf{v}_{\mathbf{T}})$.

Kirchhoff migration is usually derived by considering the solution of the ordinary least square (OLS) optimization problem. While the solution of OLS involves the pseudo-inverse of the forward matrix, it is common to approximate it as diagonal for imaging problems, which is well justified for separated targets. In that case, OLS and Kirchhoff migration yield the same support. Kirchhoff migration can also be interpreted as the solution of a different constrained optimization problem, related to the forward models of (3.19). Neglecting the scales $\xi(\omega,s)$, we can interpret migration as looking for the set of reflectivities that yield the set of signals most correlated with our measurements, under a norm constraint.

While the solution of OLS for the cross correlation model (3.21) involves writing $\widehat{\mathbf{C}}(s,\omega)$ as a vector and the matrix multiplication operations as matrix-vector multiplication, (3.23) has a direct extension for the model of(3.21). With abuse of notation for matrices $A,B$ $\langle A,B\rangle=\text{tr}(B^{T}A)$. Neglecting the scales $\xi(\omega,s)$ again, we have

As a result, a natural extension of the migration of (3.22) for migrating cross correlations is the matrix $\tilde{\mathbf{X}}$

with elements

As expected, the result of the migration is an estimation of $\mathbf{X}$ rather than $\boldsymbol{\rho}$, as our model is quadratic with respect to the reflectivities. $\tilde{\mathbf{X}}$ is a square matrix with dimensions $K\times K$, where $K$ is the number of search points in the imaging domain. If points $\mathbf{y}_{k},\mathbf{y}_{k^{\prime}}$ are associated with reflectivities $\rho_{k},\rho_{k^{\prime}}$, we can think about $\tilde{\mathbf{X}}$ as a two-variable generalized cross correlation imaging function

Note that $\tilde{\mathbf{X}}\in\mathbb{C}^{K\times K}$ is Hermitian positive definite by definition. We next see how can an image be extracted from $\mathcal{I}^{GCC}$.

The imaging function (3.27) lacks a direct physical interpretation. It evaluates the outer product of reflectivities rather than the reflecitivities themselves. We shall see that several ways to extract $\boldsymbol{\rho}$ exist, not all equal in general.

If $\tilde{\mathbf{X}}$ was exactly rank-1, i.e. $\tilde{\mathbf{X}}=\boldsymbol{\rho}\boldsymbol{\rho}^{T}$, there would be several ways to extract $\boldsymbol{\rho}$.

1. 1.

   $|\rho_{k}|^{2}=\tilde{X}_{kk}$. This is equivalent to (2.9), since the diagonal terms recreate the image generate by migrating the data to the same point $\mathbf{y}_{k}=\mathbf{y}_{k^{\prime}}$

   $$\begin{split}&\sum\limits_{s,\omega,\mathbf{R},\mathbf{R}^{\prime}}\overline{A}_{\mathbf{R},k}(s,\omega)\widehat{C}_{\mathbf{R},\mathbf{R}^{\prime}}(s,\omega)A_{\mathbf{R}^{\prime},k^{\prime}}(s,\omega)=\sum\limits_{s,\mathbf{R},\mathbf{R}^{\prime},\omega}\widehat{\mathbf{C}}_{\mathbf{R}\mathbf{R}^{\prime}}(s,\omega)e^{i\omega(t_{\mathbf{R}}^{k}(s)-t_{\mathbf{R}}(s))}e^{-i\omega(t_{\mathbf{R}^{\prime}}^{k}(s)-t_{\mathbf{R}^{\prime}}(s))}\\ &\propto\sum\limits_{s,\mathbf{R},\mathbf{R}^{\prime}}\mathbf{C}_{\mathbf{R}\mathbf{R}^{\prime}}(s,t_{\mathbf{R}}^{k}(s)-t_{\mathbf{R}}(s)-(t_{\mathbf{R}^{\prime}}^{k}(s)-t_{\mathbf{R}^{\prime}}(s))).\end{split}$$

   (3.28)

   In terms of (3.27) the image is evaluated by plugging in the same search point in both variables

   $$\mathcal{I}^{CC}(\mathbf{y}_{k})=\mathcal{I}^{GCC}(\mathbf{y}_{k},\mathbf{y}_{k}).$$

2. 2.

   $|\rho_{k}|^{2}\propto\tilde{X}_{kk_{\text{r}ef}},k_{\text{r}ef}\in[1,\dots,K]$. This is equivalent to migrating the data to a reference point with respect to one of the receiver pairs and forming an image by migrating the other receiver to a search point. In terms of (3.27) the image is evaluated by plugging a reference term in the first variable and a search point in the other

   $$\mathcal{I}^{{\text{r}ef}CC}(\mathbf{y}_{k})=\mathcal{I}^{GCC}(\mathbf{y}_{k},\mathbf{y}_{\text{ref}}).$$

3. 3.

   $|\rho_{k}|^{2}=|\mathbf{v}_{1}(\tilde{\mathbf{X}})|^{2}_{k}$, i.e., calculate the top eigenvector of $\tilde{\mathbf{X}}$. In terms of (3.27) the image is evaluated by taking, $\mathcal{V}(\mathbf{y}_{k})$, the first eigenvector of $\mathcal{I}^{GCC}(\mathbf{y}_{k},\mathbf{y}_{k^{\prime}})$, thought of as a matrix

   $$\mathcal{I}^{R1CC}(\mathbf{y}_{k})=\mathcal{V}(\mathbf{y}_{k}).$$

   We call this the rank-1 image.

The different methods are illustrated in Figure 2.

All of these methods reconstruct the reflectivity up to a global phase as expected. The first two methods have the advantage that they are only dependent on a small subdomain of the elements of $\tilde{\mathbf{X}}$. However, as mentioned in Section 2, they suffer from limitations in the achieved image resolution. The eigenvector on the other hand is a global method, relying on the entire cross correlation data to reconstruct the image.

A physical interpretation would be the following- a cluster of scatterers moving at the same speed would exhibit coherent correlations not only when correlating reflections from the same scatterer at different receivers, but also when correlating reflections from different scatterers at different receivers. Thus, the entire data can be used to estimate the most probable locations for the scatterers, given that the fields measured at both receivers were generated by the same set of scatterers. Note that single point migration can also not resolve any of the relative phases of the reflectors, while the eigenvector is only invariant with respect to a global phase.

The single point migration and two point migration are related by

with $r$ being the rank of $\tilde{\mathbf{X}}$. i.e., the single point migration is a weighted sum of all the eigenvectors, squared, by their respective eigenvalues. When there is a rapid decay in the singular values, we expect similar results for both migration methods. However, when this is not the case, we expect a different result, as we show through analysis and numerical simulations in Section 5.

Based on the observations of the previous section we suggest the following imaging function

1. 1.

   Backpropagate the cross correlation $\widehat{\mathbf{C}}(s,\omega)$ to two different points using the operator $\mathbf{A}(s,\omega)$ to form the two-point image $\tilde{\mathbf{X}}$

   $$\tilde{\mathbf{X}}=\sum\limits_{s,\omega}\overline{\mathbf{A}(s,\omega)}^{T}\widehat{\mathbf{C}}(s,\omega)\mathbf{A}(s,\omega).$$

2. 2.

   Take the rank-1 image to be the top eigenvector of $\tilde{\mathbf{X}}$

   $$\mathcal{I}^{R1CC}(\mathbf{x})=\mathbf{v}_{1}(\tilde{\mathbf{X}})$$

We can reduce the cost of the algorithm by considering a rectangular $\tilde{\mathbf{X}}$ rather than a square matrix. This means that when migrating we don’t use the same grid for both receivers. We replace our candidate image with the top left singular vector (for a tall matrix). As we show in Section 4, numerical simulations suggest that we can downsample one dimension by a factor of 10 and still retain the same resolution. There are several motivations to consider the rank-1 image:

1. 1.

   Robustness to noise. Since the eigendecomposition uses the entire cross correlation data, it is less sensitive to additive noise or even incomplete data. One could use low-rank matrix completion algorithms to reconstruct the image.

2. 2.

   Optimization interpretation. The top singular vector has an optimization interpretation: plug in the expression for $\widehat{\mathbf{C}}$ in (3.20), then

   $$\tilde{\mathbf{X}}=\sum\limits_{s,\omega}\overline{\mathbf{A}(s,\omega)}^{T}\mathbf{A}(s,\omega)\boldsymbol{\rho}\bar{\boldsymbol{\rho}}^{T}\overline{\mathbf{A}(s,\omega)}^{T}\mathbf{A}(s,\omega).$$

   (3.30)

   The top eigenvector can be expressed as the solution to an optimization problem. The top eigenvector solves:

   $$\begin{split}\mathbf{v}_{1}(\tilde{\mathbf{X}})=&\arg\max\limits_{\|\boldsymbol{\varrho}\|_{2}=1}\bar{\boldsymbol{\varrho}}^{T}\tilde{\mathbf{X}}\boldsymbol{\varrho}=\arg\max\limits_{\|\boldsymbol{\varrho}\|_{2}=1}\sum\limits_{s,\omega}\bar{\boldsymbol{\varrho}}^{T}\overline{\mathbf{A}(s,\omega)}^{T}\mathbf{A}(s,\omega)\boldsymbol{\rho}\bar{\boldsymbol{\rho}}^{T}\overline{\mathbf{A}(s,\omega)}^{T}\mathbf{A}(s,\omega)\boldsymbol{\varrho}\\ =&\arg\max\limits_{\|\boldsymbol{\varrho}\|_{2}=1}\sum\limits_{s,\omega}|\langle\mathbf{A}(s,\omega)\boldsymbol{\varrho},\mathbf{A}(s,\omega)\boldsymbol{\rho}\rangle|^{2}=\arg\max\limits_{\|\boldsymbol{\varrho}\|_{2}=1}\sum\limits_{s,\omega}|\langle\mathbf{A}(s,\omega)\boldsymbol{\varrho},\widehat{\mathbf{u}}_{\mathbf{R}}(s,\omega)\rangle|^{2}\end{split}$$

   (3.31)

   i.e, we are trying to find the $\boldsymbol{\varrho}$ which is the most correlated with the observed field measurement, given our forward model. Notice that while similar to (3.23), this is a nonlinear optimization problem, and the absolute value in (3.31) removes the sensitivity to a global phase.

In Sections 4 and 5 we compare through simulations the performance of the two-point and single-point migration schemes for plane images, that is, images with a fixed $z$ coordinate (height), so that the image coordinate $\mathbf{y}_{k}\in\mathbb{R}^{2}$. We show that as the synthetic aperture increases, the resolution of the rank-1 image tends to improve compared to the usual single point migration image. We explain this by analyzing the stationary points of the interference pattern $\tilde{\mathbf{X}}$, in Appendix B.

The main result in this appendix is that the synthetic aperture induces an anisotropy in the resolution of $\tilde{\mathbf{X}}$ in the space $\mathbf{y}_{k}\times\mathbf{y}_{k^{\prime}}\in\mathbb{R}^{4}$. For a small synthetic aperture, the width of the stationary point is approximately $\frac{\lambda H_{\mathbf{T}}}{a}$ in all directions, where $\lambda=c_{0}/f_{0}$ is the wavelength of the carrier frequency, $H_{\mathbf{T}}$ is the height of the target (the average distance from the receivers), and $a$ is the diameter of the imaging array spanned by the receivers. This is consistent with the resolution analysis in [8], where it was evaluated for a relatively short synthetic aperture, so that the resolution is dominated by the size of the array of receivers.

However, as the synthetic aperture grows larger and becomes comparable in size to the physical aperture, the anisotropy becomes significant. As a result, the stationary points have different widths in different directions. They exhibit the narrowest spot size in the direction $\widehat{\mathbf{y}}_{k}\cdot\mathbf{v}_{\mathbf{T}}=-\widehat{\mathbf{y}}_{k^{\prime}}\cdot\mathbf{v}_{\mathbf{T}}$, which is approximately $\frac{\lambda H_{\mathbf{T}}}{\mathbf{v}_{\mathbf{T}}S}$, with $\mathbf{v}_{\mathbf{T}}S$ the size of the synthetic aperture. The width is unchanged in $\widehat{\mathbf{y}}_{k}=\widehat{\mathbf{y}}_{k^{\prime}}$, which is the direction corresponding to the diagonal/single point migration (2.9). This is again consistent with the resolution analysis of [8], which suggests that the resolution of planar images is independent of the synthetic aperture size. Anisotropic localized kernels have eigenvectors whose width is smaller then the maximum width, as we show in Appendix C for several one dimensional cases. In the case of the point interference patters the width tends to the harmonic mean, which greaty improves on the maximal width. Thus, the top eigenvector $\mathcal{I}^{R1CC}(\mathbf{y}_{k})$ can provide an image with better resolution. This is indeed the case, as illustrated in Figure 3. In numerical simulations, shown in Section 4, the resolution is comparable to the linear KM image.

In Section 4 we present extensive numerical simulations that confirm the superior performance of the rank-1 image with an increasing synthetic aperture size. Then, in Section 5 we further investigate the performance by considering the case of a single target. The numerical results are in accordance with the analysis of Appendix B and Appendix C.