Home | About | Journals | Submit | Contact Us | Français |

**|**Sensors (Basel)**|**v.17(3); 2017 March**|**PMC5375899

Formats

Article sections

- Abstract
- 1. Introduction
- 2. Radar Imaging Model
- 3. Phase Error Correction for Approximated Observation-Based CS-SAR Imaging
- 4. Simulation and Experimental Results
- 5. Conclusions
- References

Authors

Related links

Sensors (Basel). 2017 March; 17(3): 613.

Published online 2017 March 17. doi: 10.3390/s17030613

PMCID: PMC5375899

Bo Li,^{1,}^{2} Falin Liu,^{1,}^{2,}^{*} Chongbin Zhou,^{1,}^{2} Yuanhao Lv,^{1,}^{2} and Jingqiu Hu^{1,}^{2}

Ram M. Narayanan, Academic Editor

Received 2016 December 18; Accepted 2017 March 15.

Copyright © 2017 by the authors.

Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

Defocus of the reconstructed image of synthetic aperture radar (SAR) occurs in the presence of the phase error. In this work, a phase error correction method is proposed for compressed sensing (CS) radar imaging based on approximated observation. The proposed method has better image focusing ability with much less memory cost, compared to the conventional approaches, due to the inherent low memory requirement of the approximated observation operator. The one-dimensional (1D) phase error correction for approximated observation-based CS-SAR imaging is first carried out and it can be conveniently applied to the cases of random-frequency waveform and linear frequency modulated (LFM) waveform without any a priori knowledge. The approximated observation operators are obtained by calculating the inverse of Omega-K and chirp scaling algorithms for random-frequency and LFM waveforms, respectively. Furthermore, the 1D phase error model is modified by incorporating a priori knowledge and then a weighted 1D phase error model is proposed, which is capable of correcting two-dimensional (2D) phase error in some cases, where the estimation can be simplified to a 1D problem. Simulation and experimental results validate the effectiveness of the proposed method in the presence of 1D phase error or weighted 1D phase error.

The synthetic aperture radar (SAR) is able to achieve image reconstruction with high resolution in both range and azimuth direction. However, higher range resolution of SAR requires wider bandwidth of signal, and thereby a higher sampling rate. Similarly, higher azimuth resolution requires longer synthetic aperture and higher sampling rate in azimuth direction. When both higher range and azimuth resolutions are expected, much higher sampling rate and larger volume of sampled data are required. It is thus a great challenge to design a high-resolution radar system hardware with limited cost.

Based on the theory of compressed sensing (CS) [1,2], both the sampling rate and the volume of sampled data can be reduced, if the reconstructed scene is sparse or compressible. The assumption made in the image reconstruction algorithm is that the observation matrix is known exactly. The observation matrix, however, depends on the mathematical model of the observation process. The uncertainties of the observation position will lead to the mismatch of the mathematical model and the observation matrix, which will cause phase error of the raw data, and thus the defocus of the reconstructed image.

The autofocus techniques, which are developed for correcting the phase error, are still attracting the interests of scientists. Chen et al. [3] proposed a method of SAR motion compensation based on parametric sparse representation to achieve autofocus of high-resolution SAR image, which assumed that the motion parameters of the radar were constant in each subaperture. Then the trajectory error can be corrected by the estimated motion parameters. Mao et al*.* [4] proposed a knowledge-aided two-dimensional (2D) autofocus approach, but the method did not utilize the CS theory or sparse regularization method to enhance the feature of the reconstructed image.

Many methods were proposed to deal with the model error in the case of CS-based radar imaging. Önhon et al. [5] proposed a sparsity-driven method for joint SAR imaging and phase error correction, which can remove either one-dimensional (1D) or 2D phase error. Kelly et al. [6] proposed a stable algorithm, which corrects phase error and reconstructs SAR image with compressively sampled data. An approach was proposed in [7], which can solve a joint optimization problem to achieve model error parameter estimation and SAR image formation simultaneously. Yang et al. [8] reported a method that can estimate the observation position error accurately and thus the quality of reconstructed image is improved significantly.

Up to now, the existing methods correcting the phase error of CS-based radar imaging were achieved by solving a two-step optimization problem, i.e., image formation processing and phase error estimation. The frameworks of CS-based image formation processing in [5,7,8] were formulated using exact observation functions, which were inefficient to be applied to the high-dimensional cases. The method proposed in [6] made a small aperture angle approximation, which hampered the application of high-resolution radar.

In [9,10], the approximated observation operators were proposed, which reduce the computational complexity and memory cost dramatically, and thus can be applied in high-dimensional cases efficiently. The approximated observation-based method proposed in [9] was called “range-azimuth decoupled method” in [10].

The intent of this paper is twofold. First, we propose a method to correct 1D phase error for approximated observation-based CS-SAR imaging by replacing the exact observation with the mentioned approximated observation in the two-step optimization framework. Then, a modified 1D phase error model is proposed to make the phase error model more precise, while the number of the unknowns does not increase in the phase error model.

The method proposed in this paper is capable of correcting the phase error in the case of approximated observation-based CS-SAR imaging. Compared with the existing methods mentioned above, the proposed method reduces the memory cost. The proposed method is achieved by solving a two-step optimization problem. In the step of image formation processing, the approximated observations can be derived from the inverse of Omega-K [11] and chirp scaling [12] algorithms for random-frequency waveform and linear frequency modulated (LFM) waveform, respectively. The images can thus be reconstructed by conducting the iterative thresholding algorithm (ITA) [9]. In the step of phase error estimation, the 1D phase error model, which can be conveniently applied to the cases of random-frequency and LFM waveforms without any a priori knowledge, is utilized first, and a weighted 1D phase error model is then proposed, which is able to model the phase error more precisely. When incorporating the a priori knowledge of the phase error structures, the method using the proposed weighted 1D phase error model can compensate the 2D phase error by solving a 1D problem. Accordingly, the weighted 1D phase error model achieves a better performance compared to the 1D phase error model, while the number of unknowns keeps the same. Furthermore, data redundancy will not be decreased for the weighted 1D phase error estimation compared with 1D phase error estimation.

This paper is organized as follows. Section 2 describes the signal models of radar imaging, along with model expressions for random-frequency and LFM waveforms. The inherent relationships between the geometric models and the phase error models are also introduced. In Section 3, a phase error correction method is proposed for approximated observation-based CS-SAR imaging, including the target reconstruction method, the phase error estimation approach, the memory cost analysis, the convergence analysis, and the computational complexity analysis. The simulation and experimental results are presented in Section 4 to validate the effectiveness of the proposed method. We conclude the work in Section 5.

Notation: We will use the subsequent notations in this paper. Vectors, matrices and operators will be denoted by bold lower case, bold upper case and roman upper case, respectively, e.g., $f$ is a vector, $S$ is a matrix, and $\mathrm{M}$ is an operator. Superscripts ${S}^{T}$, ${S}^{\ast}$ and ${S}^{H}$ denote the transpose, conjugate and Hermitian transpose of $S$, respectively.

In applications of strip-map SAR, wide bandwidth signals are always used to achieve the high range resolution. In this paper, we consider two types of signal waveforms, namely the random-frequency waveform [13] and the LFM waveform [8,12].

To overcome the limitations of the stepped-frequency SAR system, a random-frequency SAR imaging scheme was proposed in [13]. In the strip-map SAR using stepped-frequency waveform [11], the demodulated baseband echo signal for a point target can be expressed as

$$s\left(n\right)=g\cdot \mathrm{exp}\left[-j4\pi {f}_{c}\left(n\right)R/c\right]$$

(1)

where $n=1,2,\cdots ,N$ is the $n$-th frequency point in one sequence, $g$ is reflectivity coefficient of the point target, ${f}_{c}\left(n\right)$ is the frequency value of the $n$-th frequency, $R$ is the range from radar to the point target, and $c$ is the velocity of light. The stepped-frequency points are denoted as

$$\begin{array}{cc}{f}_{c}\left(n\right)={f}_{c}+n\Delta f,& n=1,2,\cdots ,N\end{array}$$

(2)

where $\Delta f$ is the frequency interval and ${f}_{c}$ is the starting frequency.

The received signal of a radar is the superposition of echoes from all scatterers in the scene. Based on Equation (1), the received signal of the $n$-th pulse in the $m$-th sequence is given by

$$s\left(m,n\right)={\displaystyle \underset{{G}_{0}}{\iint}g\left(x,y\right)\mathrm{exp}\left[-j\frac{4\pi {f}_{c}\left(n\right)R\left(m,x,y\right)}{c}\right]dxdy}$$

(3)

where $m=1,2,\cdots ,M$ denotes the $m$-th observation position, $\left(x,y\right)$ is the coordinate of a target, $g\left(x,y\right)$ is reflectivity coefficient of the target on $\left(x,y\right)$, $R\left(m,x,y\right)$ is the range from the radar to the target in the $m$-th observation position, and ${G}_{0}$ is the area illuminated by the wave beam. Here, referring to [8], we also assume that the transmitting period for each frequency is relatively short. Thus, the range from the radar to the scene is not changed in one observation.

There is a trade-off between resolution and imaging width [11]. If the scene is sparse or compressible, it is found [13] that the frequency points can be reduced and the imaging width can be increased, while the range and azimuth resolutions remain the same based on the theory of CS. The difference between stepped-frequency waveform and random-frequency waveform is the frequency interval. In the strip-map SAR using random-frequency waveform, ${f}_{c}\left(n\right)$ will be substituted by ${f}_{c}^{\prime}\left({n}_{s}\right)$, ${n}_{s}=1,2,\cdots ,{N}_{s}$, and the sequence ${f}_{c}^{\prime}\left({n}_{s}\right)$ is selected from the sequence ${f}_{c}\left(n\right)$ randomly. ${N}_{s}$ is generally much smaller than $N$. The relationship between the stepped-frequency points and the random-frequency points is expressed in matrix form as:

$${f}_{c}^{\prime}={f}_{c}{\Theta}_{r}$$

(4)

where ${f}_{c}^{\prime}\in {\u2102}^{1\times {N}_{s}}$ and ${f}_{c}\in {\u2102}^{1\times N}$ denote the set of random-frequency and stepped-frequency points, respectively, and ${\Theta}_{r}\in {\u2102}^{N\times {N}_{s}}$ is the sampling matrix.

Substituting ${f}_{c}\left(n\right)$ by ${f}_{c}^{\prime}\left({n}_{s}\right)$ in Equation (3), the received signal of the random-frequency SAR is expressed as:

$$s\left(m,{n}_{s}\right)={\displaystyle \underset{{G}_{0}}{\iint}g\left(x,y\right)\mathrm{exp}\left[-j\frac{4\pi {f}_{c}^{\prime}\left({n}_{s}\right)R\left(m,x,y\right)}{c}\right]dxdy}.$$

(5)

We consider that the scene consists of a series of point scatterers. Based on Equation (5), the received signal of the random-frequency SAR will be rewritten as:

$$s\left(m,{n}_{s}\right)={\displaystyle \sum _{k=1}^{K}g\left(k\right)\mathrm{exp}\left[-j\frac{4\pi {f}_{c}^{\prime}\left({n}_{s}\right)R\left(m,k\right)}{c}\right]}$$

(6)

where $K$ is the number of point scatterers, $g\left(k\right)$ is reflectivity coefficient of the $k$-th point, and $R$ is the range from the radar to the $k$-th point in the $m$-th observation position.

The LFM waveform is another kind of signal waveform, which can achieve wide bandwidth. Transmitted LFM waveform can be expressed as:

$${s}_{T}\left(t\right)=p\left(t\right)\mathrm{exp}\left(j2\pi {f}_{c}t\right)$$

(7)

where $t$ is fast time, ${f}_{c}$ is the carrier frequency, and $p\left(t\right)$ is denoted as:

$$p\left(t\right)=\mathrm{rect}\left(\frac{t}{T}\right)\mathrm{exp}\left(j\pi \gamma {t}^{2}\right)$$

(8)

where $\mathrm{rect}(\cdot )$ is the time window, and $\gamma $ is the chirp rate.

In the strip-map SAR using LFM waveform, the demodulated baseband echo signal for a point target is given by:

$$s\left(t\right)=g\cdot p\left(t-2R/c\right)\mathrm{exp}\left(-j4\pi {f}_{c}R/c\right)$$

(9)

where $g$ is reflectivity coefficient of the point target, $R$ is the range from the radar to the point target, and $c$ is the velocity of light.

The received data are the superposition of echoes from all scatterers in the scene. Based on Equation (9), the demodulated echo data in the $m$-th observation will be expressed as:

$$s\left(m,t\right)={\displaystyle \underset{{G}_{0}}{\iint}g\left(x,y\right)p\left[t-2R\left(m,x,y\right)/c\right]\mathrm{exp}\left[-j\frac{4\pi {f}_{c}R\left(m,x,y\right)}{c}\right]dxdy}.$$

(10)

where $g\left(x,y\right)$ is reflectivity coefficient of the target, $R\left(m,x,y\right)$ is the range from radar to the target, and ${G}_{0}$ is the area of imaging. We consider that the scene is discretized, and then the received signal will be rewritten as:

$$s\left(m,t\right)={\displaystyle \sum _{k=1}^{K}g\left(k\right)p\left[t-2R\left(m,k\right)/c\right]\mathrm{exp}\left[-j\frac{4\pi {f}_{c}R\left(m,k\right)}{c}\right]}$$

(11)

where $K$ is the number of scatterers. In practice, the continuous signal should be discretized. Then, the echo signal using LFM waveform is given by

$$s\left(m,n\right)={\displaystyle \sum _{k=1}^{K}g\left(k\right)p\left[t\left(n\right)-2R\left(m,k\right)/c\right]\mathrm{exp}\left[-j\frac{4\pi {f}_{c}R\left(m,k\right)}{c}\right]}.$$

(12)

After reviewing the existing signal models, the phase error models are analyzed as follows. In the SAR system, a radar transmits a series of pulse signals, and receives a series of echo data scattered from the illuminated scene. Considering 2-D slant-range plane [3], geometry of a SAR system is shown in Figure 1. The curved line indicates the real radar path, and the vertical axis indicates the known radar path, also nominated as “nominal path” in [3]. $P$ is the real radar position, and ${P}^{\prime}$ is the known radar position. ${T}_{0}$ is the reference point of the scene, and ${T}_{k}$ is the $k$-th target in the illuminated scene.

Under the far-field assumption [14,15], range from the known radar position to the $k$-th target can be approximated as:

$$\left|{P}^{\prime}{T}_{k}\right|\approx \left|{P}^{\prime}{T}_{0}\right|+\overrightarrow{{T}_{0}{T}_{k}}\frac{\overrightarrow{{P}^{\prime}{T}_{0}}}{\left|{P}^{\prime}{T}_{0}\right|}$$

(13)

where $\overrightarrow{{P}^{\prime}{T}_{0}}$ and $\left|{P}^{\prime}{T}_{0}\right|$ denote the direction vector and the range from ${P}^{\prime}$ to ${T}_{0}$, respectively.

Furthermore, radar platform position uncertainties will result in the difference between the known radar positon and the real radar position. Therefore, applying the far-field approximation, the range from the real radar position to the $k$-th target will be given by:

$$\left|P{T}_{k}\right|\approx \left|{P}^{\prime}{T}_{0}\right|+\overrightarrow{{T}_{0}{T}_{k}}\frac{\overrightarrow{{P}^{\prime}{T}_{0}}}{\left|{P}^{\prime}{T}_{0}\right|}+\overrightarrow{{P}^{\prime}P}\frac{\overrightarrow{{T}_{0}{P}^{\prime}}}{\left|{T}_{0}{P}^{\prime}\right|}.$$

(14)

The third item on the right side in Equation (14) arises from the radar position uncertainties, which will result in phase error of echo data. Note that the third item is a common parameter for all the scatterers.

Considering the random-frequency waveform, the received signal of the random-frequency SAR can be expressed as Equation (6). Due to the radar position uncertainties, the received signal with model error must be modified as:

$${s}_{\epsilon}\left(m,{n}_{s}\right)={\displaystyle \sum _{k=1}^{K}g\left(k\right)\mathrm{exp}\left[-j\frac{4\pi {f}_{c}^{\prime}\left({n}_{s}\right)\left[R\left(m,k\right)+\Delta R\left(m,k\right)\right]}{c}\right]}$$

(15)

where $\Delta R\left(m,k\right)$ arises from the radar position uncertainties. Referring to Equation (14), Equation (15) is approximated as:

$$\begin{array}{ll}{s}_{\epsilon}\left(m,{n}_{s}\right)& \approx {\displaystyle \sum _{k=1}^{K}g\left(k\right)\mathrm{exp}\left[-j\frac{4\pi {f}_{c}^{\prime}\left({n}_{s}\right)\left[R\left(m,k\right)+\Delta R\left(m\right)\right]}{c}\right]}\\ & =\mathrm{exp}\left[-j\frac{4\pi {f}_{c}^{\prime}\left({n}_{s}\right)\Delta R\left(m\right)}{c}\right]{\displaystyle \sum _{k=1}^{K}g\left(k\right)\mathrm{exp}\left[-j\frac{4\pi {f}_{c}^{\prime}\left({n}_{s}\right)R\left(m,k\right)}{c}\right]}\\ & =\mathrm{exp}\left[-j\frac{4\pi {f}_{c}^{\prime}\left({n}_{s}\right)\Delta R\left(m\right)}{c}\right]s\left(m,{n}_{s}\right)\end{array}$$

(16)

where $\mathrm{exp}\left[-j4\pi {f}_{c}^{\prime}\left({n}_{s}\right)\Delta R\left(m\right)/c\right]$ arises from the radar position uncertainties, which will result in phase error in echo data and must be corrected in the phase error correction process. During the compensation of phase error, $\Delta R\left(m\right)$ must be estimated. Equation (16) is expressed in matrix form as:

$${S}_{\epsilon}={D}_{W1D}\odot S$$

(17)

where ${D}_{W1D}\in {\u2102}^{M\times {N}_{s}}$ denotes the phase error, $\odot $ denotes the Hadamard product, and ${S}_{\epsilon}\in {\u2102}^{M\times {N}_{s}}$ and $S\in {\u2102}^{M\times {N}_{s}}$ denote the echo data with phase error and without phase error, respectively. The phase error arising from radar position uncertainties is

$$\begin{array}{ll}{D}_{W1D}& =\mathrm{exp}\left\{{\left[j{\phi}_{W1D}\left(1\right),\mathrm{j}{\phi}_{W1D}\left(2\right),\dots ,\mathrm{j}{\phi}_{W1D}\left(M\right)\right]}^{T}\frac{\left[{f}_{c}^{\prime}\left(1\right),{f}_{c}^{\prime}\left(2\right),\dots ,{f}_{c}^{\prime}\left({N}_{s}\right)\right]}{{f}_{c}}\right\}\\ & =\mathrm{exp}\left(j{{\phi}_{W1D}}^{T}w\right)\end{array}$$

(18)

where ${\phi}_{W1D}\in {\u2102}^{1\times M}$ denotes the fundamental phase error, and $w\in {\mathbb{R}}^{1\times {N}_{s}}$ denotes the weight. The vectors and matrix are given by:

$$\begin{array}{cc}{\phi}_{W1D}\left(m\right)=-\frac{4\pi {f}_{c}\Delta R\left(m\right)}{c},& m=1,2,\dots ,M\end{array}$$

(19)

$$\begin{array}{cc}w\left({n}_{s}\right)=\frac{{f}_{c}^{\prime}\left({n}_{s}\right)}{{f}_{c}},& {n}_{s}=1,2,\dots ,{N}_{s}\end{array}$$

(20)

$$\begin{array}{ll}{D}_{W1D}\left(m,{n}_{s}\right)& =\mathrm{exp}\left[j{\phi}_{W1D}\left(m\right)w\left({n}_{s}\right)\right]\\ & =\mathrm{exp}\left[-j\frac{4\pi \Delta R\left(m\right){f}_{c}^{\prime}\left({n}_{s}\right)}{c}\right]\end{array}.$$

(21)

The function $\mathrm{exp}(\cdot )$ in Equation (18) is defined as follows:

$$\mathrm{exp}\left(G\right)=\left[\begin{array}{cccc}{\mathrm{e}}^{{G}_{1,1}}& {\mathrm{e}}^{{G}_{1,2}}& \cdots & {\mathrm{e}}^{{G}_{1,N}}\\ {\mathrm{e}}^{{G}_{2,1}}& {\mathrm{e}}^{{G}_{2,2}}& \cdots & {\mathrm{e}}^{{G}_{2,N}}\\ \vdots & \vdots & & \vdots \\ {\mathrm{e}}^{{G}_{M,1}}& {\mathrm{e}}^{{G}_{M,2}}& \cdots & {\mathrm{e}}^{{G}_{M,N}}\end{array}\right]$$

(22)

where $G=[{G}_{i,j}]\in {\u2102}^{M\times N}$ is an arbitrary matrix. In Equation (19), ${\phi}_{W1D}\left(m\right)$ is proportional to $\Delta R\left(m\right)$, which is unknown and must be estimated in the phase error correction method.

For the narrow bandwidth system, the frequency interval $\Delta f$ can be neglected, i.e., $n\cdot \Delta f<<{f}_{c}$. Therefore, Equation (20) can be approximated as

$$\begin{array}{cc}w\left({n}_{s}\right)\approx 1,& {n}_{s}=1,2,\dots ,{N}_{s}\end{array}$$

(23)

Based on Equation (17), the received signal using random-frequency waveform with model error for the narrow bandwidth system will be expressed as:

$$\begin{array}{ll}{S}_{\epsilon}& ={D}_{W1D}\odot S\\ & =\mathrm{exp}\left(j{{\phi}_{W1D}}^{T}w\right)\odot S\\ & \approx \mathrm{exp}\left\{{\left[j{\phi}_{1D}\left(1\right),\mathrm{j}{\phi}_{1D}\left(2\right),\dots ,\mathrm{j}{\phi}_{1D}\left(M\right)\right]}^{T}\left[\underset{{N}_{s}}{\underbrace{1,1,\dots ,1}}\right]\right\}\odot S\\ & =\mathrm{exp}\left\{\mathrm{diag}\left[j{\phi}_{1D}\left(1\right),\mathrm{j}{\phi}_{1D}\left(2\right),\dots ,\mathrm{j}{\phi}_{1D}\left(M\right)\right]\right\}\cdot S\\ & =\mathrm{diag}\left[{e}^{j{\phi}_{1D}\left(1\right)},{e}^{j{\phi}_{1D}\left(2\right)},\cdots ,{e}^{j{\phi}_{1D}\left(M\right)}\right]\cdot S\\ & ={D}_{1D}S\end{array}.$$

(24)

Comparing Equation (18) to Equation (24), the 1D phase error model is an approximation of the weighted 1D phase error model in the case of the narrow bandwidth system. Note that the former can only correct 1D phase error, while the latter can correct 2D phase error.

Similarly, the received signal using LFM waveform with model error for the narrow bandwidth system can also be written as:

$${S}_{\epsilon}={D}_{1D}S$$

(25)

where ${S}_{\epsilon}\in {\u2102}^{M\times N}$ and $S\in {\u2102}^{M\times N}$ denote the echo data with and without phase error, respectively.

In this section, we will formulate the method to correct phase error for approximated observation-based CS-SAR imaging. First, some preliminary knowledge of approximated observation is reviewed. Then, a method is proposed to correct phase error for approximated observation-based CS-SAR imaging.

Referring to [9,10], the approximated observation operators will be obtained by calculating the inverse of matched filtering (MF)-based algorithm. The Omega-K algorithm, which is a type of MF-based algorithm, can be used to achieve the stepped-frequency SAR imaging [11]. In this paper, the Omega-K algorithm is also used to achieve the imaging procedure of the random-frequency SAR, which can be expressed as:

$$\mathrm{M}\left(S\right)={F}_{a}^{H}\left\{\mathrm{C}\left[\left({F}_{a}S\right)\odot {H}_{c}{\Theta}_{r}^{T}\odot {H}_{ref}\right]{F}_{r}^{H}\right\}$$

(26)

where $S\in {\u2102}^{M\times {N}_{s}}$ is the echo data without model error, $\mathrm{C}$ is the Stolt interpolation operator, ${H}_{c}$ is the range difference compensation function, and ${H}_{ref}$ is the reference function. The details of ${H}_{c}$ and ${H}_{ref}$ can be seen in [11]. $F$ and ${F}^{H}$ are DFT matrix and inverse DFT matrix, respectively. The lower notation $a$ and $r$ are the azimuth and range direction, respectively. ${\Theta}_{r}$ can be found in Equation (4). The function $\mathrm{M}:{\u2102}^{M\times {N}_{s}}\to {\u2102}^{M\times N}$ is the imaging procedure.

Based on the above exposition, the approximated observation operator $\mathrm{I}$ is obtained by calculating the inverse of imaging procedure as:

$$\mathrm{I}\left(G\right)={F}_{a}^{H}\left\{{\mathrm{C}}^{-1}\left[\left({F}_{a}G\right){F}_{r}\right]\odot {H}_{ref}^{\ast}{\Theta}_{r}\odot {H}_{c}^{\ast}\right\}$$

(27)

where $G\in {\u2102}^{M\times N}$ is the matrix form of $g\left(x,y\right)$ in Equation (3), and ${\mathrm{C}}^{-1}$ is the inverse of Stolt interpolation operator. The notation $\mathrm{I}:{\u2102}^{M\times N}\to {\u2102}^{M\times {N}_{s}}$ denotes the approximated observation operator.

Based on Equation (27), we acquire the approximated observation-based CS-SAR model (named as CS-Omega-K [16,17]) as follow:

$$\underset{G}{\mathrm{min}}\left\{{\Vert S-\mathrm{I}\left(G\right)\Vert}_{F}^{2}+\lambda {\Vert G\Vert}_{1,1}\right\}$$

(28)

where ${\Vert \cdot \Vert}_{F}$ is the Frobenius norm of a matrix, and $\lambda >0$ is a regularization parameter. ${\Vert G\Vert}_{1,1}={\Vert \mathrm{vec}\left(G\right)\Vert}_{1}$ is the ${l}_{1,1}$ (pseudo) matrix norm [18]. The regularization item ${\Vert G\Vert}_{1,1}$ also can be substituted by ${\Vert \mathrm{vec}\left(G\right)\Vert}_{q}^{q}$ to enhance feature of imaging when $0\le q\le 1$. In this paper, we simply set $q=1$, since the widely used soft-thresholding corresponds to $q=1$ and the soft-thresholding operator can be analytically specified.

In [10,19], the chirp scaling algorithm, which is also a type of MF-based algorithm, was used to achieve the LFM SAR imaging. In this paper, the chip scaling algorithm is used to achieve the imaging procedure of the LFM SAR, which can be expressed as:

$$\mathrm{M}\left(S\right)={F}_{a}^{H}\left\{{H}_{3}\odot \left\{{H}_{2}\odot \left[{H}_{1}\odot \left({F}_{a}S\right){F}_{r}\right]{F}_{r}^{H}\right\}\right\}$$

(29)

where $S\in {\u2102}^{M\times N}$ denotes the echo data without model error. ${H}_{1}$, ${H}_{2}$ and ${H}_{3}$ are the standard phase functions of the chirp scaling [19] (${\Phi}_{1}$, ${\Phi}_{2}$ and ${\Phi}_{3}$). The function $\mathrm{M}:{\u2102}^{M\times N}\to {\u2102}^{M\times N}$ is the imaging procedure.

Based on the above exposition, the approximated observation $\mathrm{I}$ is obtained by calculating the inverse of imaging procedure as:

$$\mathrm{I}\left(G\right)={F}_{a}^{H}\left\{{H}_{1}^{\ast}\odot \left\{{H}_{2}^{\ast}\odot \left[{H}_{3}^{\ast}\odot \left({F}_{a}G\right){F}_{r}\right]{F}_{r}^{H}\right\}\right\}$$

(30)

where $G\in {\u2102}^{M\times N}$ is the matrix form of $g\left(x,y\right)$ in Equation (10). The notation $\mathrm{I}:{\u2102}^{M\times N}\to {\u2102}^{M\times N}$ denotes the approximated observation operator.

Based on Equation (30), we are able to acquire the approximated observation-based CS-SAR model (named as CS-chirp scaling [10]) as follow:

$$\underset{G}{\mathrm{min}}\left\{{\Vert S-\mathrm{I}\left(G\right)\Vert}_{F}^{2}+\lambda {\Vert G\Vert}_{1,1}\right\}.$$

(31)

The aim of briefing the published research is to suggest a new phase error correction method, which will be formulated as follows and can achieve the phase error correction for approximated observation-based CS-SAR imaging.

Thus far, we have derived the inexact (nominal) observation model, but the model error is still unknown. Thus, CS-SAR imaging models must be modified. Considering the 1D phase error models in Equations (24) and (25), the CS-SAR imaging models in Equation (28) using random-frequency waveform and in Equation (31) using LFM waveform will be modified as:

$$\underset{G,{D}_{1D}}{\mathrm{min}}\left\{{\Vert {S}_{\epsilon}-{D}_{1D}\mathrm{I}\left(G\right)\Vert}_{F}^{2}+\lambda {\Vert G\Vert}_{1,1}\right\}.$$

(32)

The cost function is denoted by:

$$\mathrm{J}\left(G,{D}_{1D}\right)={\Vert {S}_{\epsilon}-{D}_{1D}\mathrm{I}\left(G\right)\Vert}_{F}^{2}+\lambda {\Vert G\Vert}_{1,1}$$

(33)

where the phase error ${D}_{1D}$ is written, referring to Equation (24), as:

$${D}_{1D}=\mathrm{diag}\left[{e}^{j{\phi}_{1D}\left(1\right)},{e}^{j{\phi}_{1D}\left(2\right)},\cdots ,{e}^{j{\phi}_{1D}\left(M\right)}\right].$$

(34)

By solving Equation (32), both the image and the 1D phase error can be jointly obtained. Using similar methods proposed in the [5,8], both the image and phase error are estimated by solving a two-step optimization problem. In the first step, $G$ is estimated by minimizing the cost function Equation (33) using the given ${D}_{1D}$. In the second step, ${D}_{1D}$ is obtained by minimizing the cost function using the estimated $G$. By iteratively updating $G$ and ${D}_{1D}$ as described above, both the image and the 1D phase error can be jointly estimated. The algorithm flow is outlined as Algorithm 1. Initially, let

$${D}_{1D}=\mathrm{diag}\left[\underset{M}{\underbrace{1,1,\cdots ,1}}\right],$$

(35)

i.e., $\begin{array}{cc}{\phi}_{1D}\left(m\right)=0,& m=1,2,\dots ,M\end{array}$.

Algorithm 1. 1D Phase Error Correction for Approximated Observation-Based CS-SAR Imaging. |

Initialize: $i=0$, ${D}_{1D}^{0}=\mathrm{diag}\left[1,1,\dots ,1\right]$ |

Step 1: Image reconstruction ${G}^{i+1}={\mathrm{argmin}}_{G}\mathrm{J}\left(G,{D}_{1D}^{i}\right)$ |

Step 2: Phase error estimation ${D}_{1D}^{i+1}={\mathrm{argmin}}_{{D}_{1D}}\mathrm{J}\left({G}^{i+1},{D}_{1D}\right)$ |

Step 3: Let $i=i+1$ and return to step 1. |

Terminate when $i$ is equal to a preset threshold ${I}_{\mathrm{max}}$. |

Steps 1 and 2 are the major steps of Algorithm 1. In Step 1, the image is reconstructed by the given phase error. The optimization problem is expressed as:

$$\begin{array}{ll}{G}^{i+1}& ={\mathrm{argmin}}_{G}\mathrm{J}\left(G,{D}_{1D}^{i}\right)\\ & ={\mathrm{argmin}}_{G}\left\{{\Vert {S}_{\epsilon}-{D}_{1D}^{i}\mathrm{I}\left(G\right)\Vert}_{F}^{2}+\lambda {\Vert G\Vert}_{1,1}\right\}\end{array}$$

(36)

Since ${D}_{1D}^{i}$ is a diagonal matrix, seen in Equation (34), Equation (36) can be rewritten as:

$${G}^{i+1}={\mathrm{argmin}}_{G}\left\{{\Vert {{D}_{1D}^{i}}^{\ast}{S}_{\epsilon}-\mathrm{I}\left(G\right)\Vert}_{F}^{2}+\lambda {\Vert G\Vert}_{1,1}\right\}$$

(37)

Then, due to the linearity of $\mathrm{I}$, we solve the optimization problem efficiently by iterative thresholding algorithm (ITA) [9,10,20]. Let ${G}^{0}=0$, and the optimal solution is obtained by using the iteration:

$${G}^{i+1}={\mathrm{E}}_{1,\lambda \mu}\left({G}^{i}+\mu \mathrm{M}\left({{D}_{1D}^{i}}^{\ast}{S}_{\epsilon}-\mathrm{I}\left({G}^{i}\right)\right)\right)$$

(38)

where $\mu $ is step size, which controls the convergence of the ITA. In Equation (38), ${\mathrm{E}}_{1,\lambda \mu}$ is a thresholding operator, which is defined as:

$${\mathrm{E}}_{1,\lambda \mu}\left(G\right)=\left[\begin{array}{cccc}{\mathrm{E}}_{1,\lambda \mu}\left({G}_{1,1}\right)& {\mathrm{E}}_{1,\lambda \mu}\left({G}_{1,2}\right)& \cdots & {\mathrm{E}}_{1,\lambda \mu}\left({G}_{1,N}\right)\\ {\mathrm{E}}_{1,\lambda \mu}\left({G}_{2,1}\right)& {\mathrm{E}}_{1,\lambda \mu}\left({G}_{2,2}\right)& \cdots & {\mathrm{E}}_{1,\lambda \mu}\left({G}_{2,N}\right)\\ \vdots & \vdots & & \vdots \\ {\mathrm{E}}_{1,\lambda \mu}\left({G}_{M,1}\right)& {\mathrm{E}}_{1,\lambda \mu}\left({G}_{M,2}\right)& \cdots & {\mathrm{E}}_{1,\lambda \mu}\left({G}_{M,N}\right)\end{array}\right]$$

(39)

where soft-thresholding operator ${\mathrm{E}}_{1,\lambda \mu}$ is analytically specified as:

$${\mathrm{E}}_{1,\lambda \mu}\left(x\right)=\{\begin{array}{ll}\mathrm{sign}\left(x\right)\left(\left|x\right|-\lambda \mu \right),\hfill & \mathrm{if}\text{}\left|x\right|\ge \lambda \mu ,\hfill \\ 0,\hfill & \mathrm{otherwise}.\hfill \end{array}$$

(40)

In the proposed method, we simply set $\mu =1$ ($0<\mu <{\Vert \mathrm{I}\Vert}_{2}^{2}$ [9]). With fixed parameter $\mu $, the parameter $\lambda $ will depend on the sparsity parameter ${K}_{0}$. We set $\lambda ={\left|G\right|}_{{K}_{0}+1}/\mu $ (${\left|G\right|}_{{K}_{0}}$ is the ${K}_{0}$-th largest component of $G$ in magnitude).

In Step 2, the phase error is estimated by the reconstructed image. The optimization problem is expressed as:

$$\begin{array}{ll}{D}_{1D}^{i+1}& ={\mathrm{argmin}}_{{D}_{1D}}\mathrm{J}\left({G}^{i+1},{D}_{1D}\right)\\ & ={\mathrm{argmin}}_{{D}_{1D}}\left\{{\Vert {S}_{\epsilon}-{D}_{1D}\mathrm{I}\left({G}^{i+1}\right)\Vert}_{F}^{2}+\lambda {\Vert {G}^{i+1}\Vert}_{1,1}\right\}\end{array}.$$

(41)

Since $\lambda {\Vert {G}^{i+1}\Vert}_{1,1}$ is a constant, Equation (41) is rewritten as:

$${D}_{1D}^{i+1}={\mathrm{argmin}}_{{D}_{1D}}{\Vert {S}_{\epsilon}-{D}_{1D}\mathrm{I}\left({G}^{i+1}\right)\Vert}_{F}^{2}.$$

(42)

From Equations (34) and (42), it is quite obvious that the data error of the $m$-th observation position is only related to ${\phi}_{1D}\left(m\right)$. Then, the optimization problem Equation (42) is divided into a set of independent optimization problems as follows:

$${\phi}_{1D}^{i+1}\left(m\right)={\mathrm{argmin}}_{{\phi}_{1D}\left(m\right)}{\Vert {\left[{S}_{\epsilon}\right]}_{m}-{e}^{j{\phi}_{1D}\left(m\right)}{\left[\mathrm{I}\left({G}^{i+1}\right)\right]}_{m}\Vert}_{2}^{2}$$

(43)

where notation ${[\cdot ]}_{m}$ denotes the $m$-th row of the matrix. The phase error ${\phi}_{1D}^{i+1}\left(m\right)$ and/or ${D}_{1D}^{i+1}\left(m,m\right)$ can be obtained as follows (the details can be found in the Appendix A):

$${\phi}_{1D}^{i+1}\left(m\right)=\mathrm{angle}\left({\left[{S}_{\epsilon}\right]}_{m}\cdot {\left[\mathrm{I}\left({G}^{i+1}\right)\right]}_{m}^{H}\right)$$

(44)

Substituting Equation (44) into Equation (34), Step 2 of Algorithm 1 can be implemented straightforwardly.

Considering the SAR imaging system using random-frequency waveform without narrow bandwidth approximation, the SAR imaging model (28) is modified as:

$$\underset{G,{D}_{W1D}}{\mathrm{min}}\left\{{\Vert {S}_{\epsilon}-{D}_{W1D}\odot \mathrm{I}\left(G\right)\Vert}_{F}^{2}+\lambda {\Vert G\Vert}_{1,1}\right\}.$$

(45)

The cost function is denoted by

$$\mathrm{J}\left(G,{D}_{W1D}\right)={\Vert {S}_{\epsilon}-{D}_{W1D}\odot \mathrm{I}\left(G\right)\Vert}_{F}^{2}+\lambda {\Vert G\Vert}_{1,1}$$

(46)

where ${D}_{W1D}\in {\u2102}^{M\times {N}_{s}}$ can be seen in Equation (18). By solving Equation (45), both the image and the weighted 1D phase error will be estimated jointly. The method is similar to Algorithm 1, and the algorithm flow is outlined in Algorithm 2. In the initialization, let

$${D}_{W1D}=\mathrm{ones}\left(M,{N}_{s}\right)=\left[\begin{array}{cccc}1& 1& \cdots & 1\\ 1& 1& \cdots & 1\\ \vdots & \vdots & & \vdots \\ 1& 1& \cdots & 1\end{array}\right],$$

(47)

i.e., $\begin{array}{cc}{\phi}_{W1D}\left(m\right)=0,& m=1,2,\cdots ,M\end{array}$.

Algorithm 2. Weighted 1D Phase Error Correction for Approximated Observation-Based CS-SAR Imaging. |

Initialize: $i=0$, ${D}_{W1D}^{0}=\mathrm{ones}\left(M,{N}_{s}\right)$ |

Step 1: Image reconstruction ${G}^{i+1}={\mathrm{argmin}}_{G}\mathrm{J}\left(G,{D}_{W1D}^{i}\right)$ |

Step 2: Phase error estimation ${D}_{W1D}^{i+1}={\mathrm{argmin}}_{{D}_{W1D}}\mathrm{J}\left({G}^{i+1},{D}_{W1D}\right)$ |

Step 3: Let $i=i+1$ and return to step 1. |

Terminate when $i$ is equal to a preset threshold ${I}_{\mathrm{max}}$. |

In Step 1, the optimization problem is expressed as:

$$\begin{array}{ll}{G}^{i+1}& ={\mathrm{argmin}}_{G}\mathrm{J}\left(G,{D}_{W1D}^{i}\right)\\ & ={\mathrm{argmin}}_{G}\left\{{\Vert {S}_{\epsilon}-{D}_{W1D}^{i}\odot \mathrm{I}\left(G\right)\Vert}_{F}^{2}+\lambda {\Vert G\Vert}_{1,1}\right\}\\ & ={\mathrm{argmin}}_{G}\left\{{\Vert {{D}_{W1D}^{i}}^{\ast}\odot {S}_{\epsilon}-\mathrm{I}\left(G\right)\Vert}_{F}^{2}+\lambda {\Vert G\Vert}_{1,1}\right\}\end{array}.$$

(48)

Let ${G}^{0}=0$, and use the iteration

$${G}^{i+1}={\mathrm{E}}_{1,\lambda \mu}\left({G}^{i}+\mu \mathrm{M}\left({{D}_{W1D}^{i}}^{\ast}\odot {S}_{\epsilon}-\mathrm{I}\left({G}^{i}\right)\right)\right)$$

(49)

Then, the optimization problem (48) is efficiently solved by the same method as Algorithm 1.

In Step 2, the optimization problem is expressed as:

$$\begin{array}{ll}{D}_{W1D}^{i+1}& ={\mathrm{argmin}}_{{D}_{W1D}}\mathrm{J}\left({G}^{i+1},{D}_{W1D}\right)\\ & ={\mathrm{argmin}}_{{D}_{W1D}}\left\{{\Vert {S}_{\epsilon}-{D}_{W1D}\odot \mathrm{I}\left({G}^{i+1}\right)\Vert}_{F}^{2}+\lambda {\Vert {G}^{i+1}\Vert}_{1,1}\right\}\\ & ={\mathrm{argmin}}_{{D}_{W1D}}{\Vert {S}_{\epsilon}-{D}_{W1D}\odot \mathrm{I}\left({G}^{i+1}\right)\Vert}_{F}^{2}\end{array}.$$

(50)

From Equations (18) and (50), it can be seen that the data error of the $m$-th observation position is only related to ${\phi}_{W1D}\left(m\right)$. Then, the optimization problem in Equation (50) can also be divided into a set of independent optimization problems as follows:

$${\phi}_{W1D}^{i+1}\left(m\right)={\mathrm{argmin}}_{{\phi}_{W1D}\left(m\right)}{\Vert {\left[{S}_{\epsilon}\right]}_{m}-\mathrm{exp}\left(\mathrm{j}{\phi}_{W1D}\left(m\right)\cdot w\right)\odot {\left[\mathrm{I}\left({G}^{i+1}\right)\right]}_{m}\Vert}_{2}^{2}.$$

(51)

The weight $w$, which is a priori knowledge, must be obtained before solving Equation (51). Equation (51) is an unconstrained optimization problem, which can be solved using gradient-based optimization method [8]. In Equation (51), there is only one unknown ${\phi}_{W1D}\left(m\right)$ to be estimated and it belongs to a real number. Hence, the optimal solution is able to be obtained by line search [21].

The proposed method consists of a two-step optimization problem. In the first step, the image is reconstructed. In the second step, the phase error is estimated and compensated. The second step is divided into a set of independent optimization problems. Thus the biggest matrix will appear in the first step. There are some notations to be denoted as follows. $M$, $N$ and $K$ denote the number of azimuth samples, the number of range samples, and the number of scatterer points after discretizing the scene, respectively. The size of the sensing matrix in the phase error correction for exact observation based CS-SAR imaging is $MN\times K$ complex numbers, but the size of the matrix in the proposed method is $M\times N$ complex numbers. Therefore, the proposed method can reduce the memory cost and the hardware requirement significantly.

The proposed method is achieved by calculating a two-step optimization problem. The cost function of 1D phase error correction or weighted 1D phase error correction for approximated observation-based CS-SAR imaging can be denoted as $\mathrm{J}\left(G,D\right)$.

In the first step, ${G}^{i+1}={\mathrm{argmin}}_{G}\mathrm{J}\left(G,{D}^{i}\right)$, and thus $\mathrm{J}\left({G}^{i+1},{D}^{i}\right)\le \mathrm{J}\left({G}^{i},{D}^{i}\right)$.

In the second step, ${D}^{i+1}={\mathrm{argmin}}_{D}\mathrm{J}\left({G}^{i+1},D\right)$, i.e., $\mathrm{J}\left({G}^{i+1},{D}^{i+1}\right)\le \mathrm{J}\left({G}^{i+1},{D}^{i}\right)$.

Combining the two steps, there will be $\mathrm{J}\left({G}^{i+1},{D}^{i+1}\right)\le \mathrm{J}\left({G}^{i},{D}^{i}\right)$, and then the cost function $\mathrm{J}\left(G,D\right)$ is a monotonous decrease function. Since $\mathrm{J}\left(G,D\right)\ge 0$, the convergence of this process can be guaranteed.

In order to calculate the computational complexity, we have defined some notations, which are the number of azimuth samples $M$, the number of range samples $N$, the number of required iteration in the first step ${I}_{1}$, and the number of required iteration in algorithm ${I}_{\mathrm{max}}$.

In the first step, the image is reconstructed. This process is an iterative process. One-step iteration includes the calculations of imaging procedure, approximated observation procedure, and thresholding operation. Both imaging procedure and approximated observation procedure have the same computational complexity of $O\left[MN{\mathrm{log}}_{2}\left(MN\right)\right]$. The complexity of thresholding operator is order $O\left(MN\right)$. Thus, the complexity of the first step is at the order $O\left[{I}_{1}MN{\mathrm{log}}_{2}\left(MN\right)\right]$.

In the second step, the phase error is estimated and compensated. The optimal solution of the second step in Algorithm 1 can be analytically expressed. The complexity of this procedure is order $O\left(MN\right)$. The optimal solution of the second step in Algorithm 2 is able to be obtained by line search. The complexity of this procedure is order $O\left({I}_{2}MN\right)$, where ${I}_{2}$ denotes the number of iteration in line search.

Based on the above exposition, the computational complexities of Algorithms 1 and 2 are at the order $O\left[{I}_{\mathrm{max}}{I}_{1}MN{\mathrm{log}}_{2}\left(MN\right)\right]$ and $O\left\{{I}_{\mathrm{max}}MN\left[{I}_{1}{\mathrm{log}}_{2}\left(MN\right)+{I}_{2}\right]\right\}$, respectively.

In this section, we will present a series of simulation and experimental results to verify the effectiveness of the proposed method. All experiments are performed using MATLAB R2015a on a PC equipped with an Intel Core i5-4590 CPU (3.30 GHz and 8 GB memory, University of Science and Technology of China, Hefei, China).

The simulation parameters for the random-frequency waveform are shown in Table 1. We then choose a sparse imaging scene consisting of four-point scatterers, the positions of which are shown in Table 2. In order to make each observation position have the same data volume to estimate phase error, the random-frequency points are selected as follows: in one observation position, the frequency points are selected randomly from stepped-frequency points, and then this selecting scheme is utilized in each observation position. The raw data are first generated in the time domain by exact observation. Then, we add the Gaussian noise and the 1D phase error to the data. The signal-to-noise ratio (SNR) is 20 dB. Two types of phase error are 1D quadratic phase error and 1D random phase error. The 1D quadratic phase error is distributed in [−π/2, +π/2], and the 1D random phase error is uniformly distributed in [−0.8π, +0.8π].

Figure 2 shows the phase error distribution and the imaging results for the random-frequency waveform. Figure 2a shows the results of the 1D quadratic phase error, and Figure 2b shows the results of the 1D random phase error. The top subfigures show the 1D phase error added to the raw data. The middle subfigures show the results of approximated observation-based CS-SAR imaging (named as CS-Omega-K [16,17]) without phase error correction. The bottom subfigures show the results of Algorithm 1. For both imaging methods, the approximated observation operator is acquired from the inverse of Omega-K algorithm, and the sparsity is set to ${K}_{0}=12$, which determines the threshold of ITA.

Results for the random-frequency waveform. (**Top**) 1D quadratic phase error and 1D random phase error. (**Middle**) Results of CS-Omega-K without phase error correction. (**Bottom**) Results of Algorithm 1. (**a**) Results for 1D quadratic phase error; and (**b**) results **...**

In the middle subfigures, the images of point targets are defocused in the azimuth direction due to the 1D phase error. In the bottom subfigures, the results show the effectiveness of Algorithm 1, which can correct the 1D phase error and focus the images better in the azimuth direction.

Then, we conduct experiments to compare the performance of the proposed method with the existing method [8]. The SNR is 20 dB, and the extent of 1D random phase error is $0.1\pi $.

Figure 3 shows the imaging results of the proposed method and other existing methods. Figure 3a shows the reconstruction result without phase error correction. Figure 3b shows the reconstruction result with compensation of observation position errors [8]. Figure 3c shows the result of the proposed method (Algorithm 1). In Figure 3a, it can be seen that the 1D phase error will defocus the reconstructed image in the azimuth direction without phase error correction. Both the proposed method and the existing method can correct phase error and focus image better in Figure 3b,c, but the proposed method requires much less memory cost referring to Section 3.4. In the simulation of the proposed method, the numbers of azimuth samples and range samples are $M=98$ and $N=1536$, respectively. Thus, the size of the matrix in the proposed method is $M\times N=98\times 1536$ complex numbers, and the memory requirement is 1.15 MB. In the simulation of the existing method, the scene needs to be discretized, which is taken as $40\times 40$. Thus, the number of scatterer points after discretizing the scene is $K=1600$. Therefore, the size of the matrix in the existing method is $MN\times K=98\times 1536\times 1600$ complex numbers, and the memory requirement is 1.79 GB.

Imaging results of the different methods: (**a**) reconstruction result without phase error correction; (**b**) reconstruction result with compensation of observation position errors; and (**c**) reconstruction result of the proposed method (Algorithm 1).

Next, we perform a series of one-point simulations for investigating the influence of the sparsity parameter ${K}_{0}$ and the noise size to the proposed algorithm. The imaging scene is composed of one point scatterer. Then, the target-to-background ratio (TBR) [5] can be used as a criterion to evaluate the performance of the proposed method conveniently. The TBR of a reconstructed image can be given as follows:

$$\mathrm{TBR}\left(X\right)=20{\mathrm{log}}_{10}\left(\frac{{\mathrm{max}}_{i\in T}\left|{X}_{i}\right|}{\frac{1}{{I}_{B}}{\displaystyle {\sum}_{j\in B}\left|{X}_{j}\right|}}\right)$$

(52)

where $T$ and $B$ denote the pixel indices for the target and the background regions, respectively. ${I}_{B}$ is the number of background pixels.

The extent of 1D random phase error is varied from $\left[-0.2\pi ,+0.2\pi \right]$ to $\left[-0.5\pi ,+0.5\pi \right]$. First, the sparsity parameter ${K}_{0}$ is changed from 3 to 11. The SNR is 20 dB. Second, we simulate the cases of different SNR from 0 dB to 30 dB. The sparsity parameter ${K}_{0}$ is set to 4. Other simulation parameters are the same as Table 1.

Figure 4 shows the behavior of TBR values with respect to the sparsity parameter and the noise size, respectively. In Figure 4a, the TBR decreases with the increase of sparsity parameter. Therefore, exact prior knowledge can contribute to the choice of sparsity parameter positively. In Figure 4b, the TBR increases as the SNR increases. Note that there is an insignificant difference between different phase error sizes in the case of high SNR.

The RADARSAT-1 is a well-known satellite SAR using LFM waveform, and the raw data were collected on 16 June 2002 about Vancouver region. The experimental parameters for the LFM waveform are shown in Table 3 [9,12].

In the region of English Bay, there are four sparsely distributed vessels, thus the region is a typically sparse scene. Then, the proposed method (Algorithm 1) can be applied to reconstruction of the scene. First, the chirp scaling algorithm, which is one type of MF-based algorithm, is applied to reconstruction of the region of English Bay. Second, approximated observation-based CS-SAR imaging method (named as CS-chirp scaling [10]) without phase error correction is applied, where the approximated observation operator is acquired from the inverse of chirp scaling algorithm. We set the sparsity to ${K}_{0}=10,000$. Third, Algorithm 1 is utilized to reconstruct the scene, where the same approximated observation operator and the same sparsity are used.

Figure 5a,b shows the results of the three methods and the enlarged region, respectively. In the top subfigures of Figure 5, the chirp scaling algorithm reconstructs the scene with serious sidelobes. In the middle subfigures, it can be seen that CS-chirp scaling can suppress sidelobes and improve resolution. In the bottom subfigures, the reconstructed scene can be focused better in the azimuth direction. Thus, Algorithm 1 can not only suppress sidelobes, but also achieve phase error compensation.

Results for the LFM waveform. (**Top**) Results of MF-based algorithm (chirp scaling). (**Middle**) Results of CS-chirp scaling without phase error correction. (**Bottom**) Results of Algorithm 1. (**a**) Application results on RADARSAT-1 (region of English Bay). (**b** **...**

We further use the image entropy [22] as a criterion to quantitatively evaluate the performance of the three methods. The entropy of an image is defined as follows:

$$\mathrm{Entropy}\left(X\right)=-{\displaystyle \sum _{i}{p}_{i}}{\mathrm{log}}_{2}\left({p}_{i}\right)$$

(53)

where ${p}_{i}$ is the histogram of the recovered gray level image.

The entropy values of the reconstructed images by the three methods are shown in Table 4. It can be seen that the proposed method (Algorithm 1) has a lower entropy, and thereby exhibits a better performance of image reconstruction. In Table 4, there is no significant difference in entropy between CS-chirp scaling method and the proposed method (Algorithm 1), since the estimated phase error (shown in Figure 6) is not too large. The trajectory of the satellite is much more stable than that of an airplane, and thus the phase error is smaller. If the phase error arising from observation position uncertainties is large enough, there will be a significant difference.

Furthermore, with narrow bandwidth approximation, the phase error of echo data along range time is constant for one observation position, which can be seen in Equation (34). However, considering the SAR imaging system using random-frequency waveform without narrow bandwidth approximation, the phase error of echo data along range time is different for one observation position, which can be seen in Equation (18). The number of unknowns for weighted 1D phase error model is also *M*, which is the same as the 1D phase error model. The simulation parameters are the same as Table 1. We then set one point scatterer. The position of the target is shown in Table 5. The raw data are first generated in the time domain by exact observation. Then, we add the Gaussian noise and the weighted 1D phase error to the data. The SNR is 20 dB. Two types of phase error ${\phi}_{W1D}\in {\u2102}^{M\times 1}$ are utilized in simulations. The extent of the error is $1/8$ of the wavelength so that the weighted 1D quadratic phase error is distributed in $\left[-\pi /8,+\pi /8\right]$, and the weighted 1D random phase error is uniformly distributed in $\left[-\pi /8,+\pi /8\right]$. The weight $w\in {\u2102}^{1\times {N}_{s}}$ depends on the selected frequency points.

The imaging scene is composed of one point scatterer. Then, the TBR can also be used as a criterion to evaluate the performance of the two methods conveniently.

The images are reconstructed by the two imaging methods with sparsity ${K}_{0}=5$. TBR values of the reconstructed images by the Algorithm 1 and Algorithm 2 are listed in Table 6. It can be seen that the Algorithm 2 has a higher TBR than Algorithm 1, and thereby Algorithm 2 has a better performance.

For presenting the effect of the weight $w\in {\u2102}^{1\times {N}_{s}}$, we change the frequency interval from 0.33 MHz to 0.66 MHz, and then the bandwidth will be changed to 1024 MHz. Other simulation parameters are the same as Table 1. TBR values with larger frequency interval are given in Table 7. Comparing Table 6 with Table 7, it can be observed that there will be more significant difference under the case of the larger frequency interval. Accordingly, the weighted 1D phase error correction for CS-SAR imaging will be more effective in the case of large frequency interval.

Further, we increase the extent of phase error and then conduct the simulations. The weighted 1D quadratic phase error is distributed in $\left[-\pi /2,+\pi /2\right]$, and the weighted 1D random phase error is uniformly distributed in $\left[-\pi /2,+\pi /2\right]$. The simulation parameters are the same as Table 1. TBR values with a larger extent of phase error are shown in Table 8. Comparing Table 6 with Table 8, it can be observed that there is a more significant difference under the case with a larger extent of phase error.

In this paper, we proposed a phase error correction method for approximated observation-based CS-SAR imaging. The proposed method yields a clear advantage in term of memory cost over conventional CS-based autofocus algorithms. The 1D phase error model, which can be conveniently utilized in autofocus technique without any a priori knowledge, is based on narrow bandwidth approximation. We also analyzed the inherent relationship between the geometric model and the phase error model in the case of random-frequency waveform. By incorporating a priori knowledge, the weighted 1D phase error model was proposed, which corrects the 2D phase error by estimating a 1D problem. The proposed weighted 1D phase error model is more precise than the 1D phase error model, thus has a better performance.

Although the weighted 1D phase error correction method proposed in this work can be applied only in the case of random-frequency waveform rather than the LFM one, since different waveforms have different geometry and signal models, it would be possible to apply the weighted 1D phase error correction to LFM waveform with appropriate modification, and this will be our future effort.

This work is supported by the National Natural Science Foundation of China (Grant No. 61431016).

In this appendix, we will derive Equation (44) from Equation (43). The optimization problem, Equation (43), is as follows:

$${\phi}_{1D}^{i+1}\left(m\right)={\mathrm{argmin}}_{{\phi}_{1D}\left(m\right)}{\Vert {\left[{S}_{\epsilon}\right]}_{m}-{e}^{j{\phi}_{1D}\left(m\right)}{\left[\mathrm{I}\left({G}^{i+1}\right)\right]}_{m}\Vert}_{2}^{2}$$

For simplifying the notations, let $a={\left[{S}_{\epsilon}\right]}_{m}\in {\u2102}^{1\times N}$, $\theta ={\phi}_{1D}\left(m\right)$ and $b={\left[\mathrm{I}\left({G}^{i+1}\right)\right]}_{m}\in {\u2102}^{1\times N}$. Then, the cost function can be simplified as:

$$\begin{array}{ll}{\Vert {\left[{S}_{\epsilon}\right]}_{m}-{e}^{j{\phi}_{1D}\left(m\right)}\cdot {\left[\mathrm{I}\left({G}^{i+1}\right)\right]}_{m}\Vert}_{2}^{2}& ={\Vert a-{e}^{j\theta}b\Vert}_{2}^{2}\\ & =\left(a-{e}^{j\theta}b\right){\left(a-{e}^{j\theta}b\right)}^{H}\\ & =a{a}^{H}+b{b}^{H}-{e}^{j\theta}b{a}^{H}-a{\left({e}^{j\theta}b\right)}^{H}\\ & =a{a}^{H}+b{b}^{H}-2\mathrm{Re}\left[a{\left({e}^{j\theta}b\right)}^{H}\right]\\ & =a{a}^{H}+b{b}^{H}-2\mathrm{Re}\left[{e}^{-j\theta}a{b}^{H}\right]\end{array}$$

Since $\mathrm{arg}{\mathrm{max}}_{\theta}\mathrm{Re}\left[{e}^{-j\theta}a{b}^{H}\right]=\mathrm{angle}\left(a{b}^{H}\right)$, the equation becomes

$$\mathrm{arg}{\mathrm{min}}_{\theta}{\Vert a-{e}^{j\theta}b\Vert}_{2}^{2}=\mathrm{angle}\left(a{b}^{H}\right)$$

Finally, we can get

$$\begin{array}{ll}{\phi}_{1D}^{i+1}\left(m\right)& ={\mathrm{argmin}}_{{\phi}_{1D}\left(m\right)}{\Vert {\left[{S}_{\epsilon}\right]}_{m}-{e}^{j{\phi}_{1D}\left(m\right)}{\left[\mathrm{I}\left({G}^{i+1}\right)\right]}_{m}\Vert}_{2}^{2}\\ & =\mathrm{angle}\left(a{b}^{H}\right)\\ & =\mathrm{angle}\left({\left[{S}_{\epsilon}\right]}_{m}\cdot {\left[\mathrm{I}\left({G}^{i+1}\right)\right]}_{m}^{H}\right)\end{array}$$

Author Contributions

All authors contributed solidly to the presented work. Bo Li conceived, designed and performed the experiments, analyzed the data, and prepared the paper; Falin Liu offered valuable advices and checked the manuscript; Chongbin Zhou helped provide the raw data of RADARSAT-1; and Chongbin Zhou, Yuanhao Lv and Jingqiu Hu provided their valuable suggestions to improve this study.

Conflicts of Interest

The authors declare no conflict of interest.

1. Donoho D.L. Compressed sensing. IEEE Trans. Inf. Theory. 2006;52:1289–1306. doi: 10.1109/TIT.2006.871582. [Cross Ref]

2. Candès E.J., Romberg J., Tao T. Robust uncertainty principles: Exact signal reconstruction from highly incomplete frequency information. IEEE Trans. Inf. Theory. 2006;52:489–509. doi: 10.1109/TIT.2005.862083. [Cross Ref]

3. Chen Y.C., Li G., Zhang Q., Zhang Q.J., Xia X.G. Motion Compensation for Airborne SAR via Parametric Sparse Representation. IEEE Trans. Geosci. Remote Sens. 2016;99:1–12. doi: 10.1109/LGRS.2016.2517095. [Cross Ref]

4. Mao X., Zhu D. Two-dimensional Autofocus for Spotlight SAR Polar Format Imagery. IEEE Trans. Comput. Imaging. 2016;2:524–539. doi: 10.1109/TCI.2016.2612945. [Cross Ref]

5. Önhon N.Ö., Çetin M. A sparsity-driven approach for joint SAR imaging and phase error correction. IEEE Trans. Image Process. 2012;21:2075–2088. doi: 10.1109/TIP.2011.2179056. [PubMed] [Cross Ref]

6. Kelly S.I., Yaghoobi M., Davies M.E. Auto-Focus for Compressively Sampled SAR; Proceedings of the 1st International Workshop on Compressed Sensing Applied to Radar; Bonn, Germany. 14–16 May 2012.

7. Wei S.J., Zhang X.L., Shi J. An autofocus approach for model error correction in compressed sensing SAR imaging; Proceedings of the 2012 IEEE International Geoscience and Remote Sensing Symposium; Munich, Germany. 22–27 July 2012; pp. 3987–3990.

8. Yang J., Huang X., Thompson J., Jin T., Zhou Z. Compressed sensing radar imaging with compensation of observation position error. IEEE Trans. Geosci. Remote Sens. 2014;52:4608–4620. doi: 10.1109/TGRS.2013.2283054. [Cross Ref]

9. Fang J., Xu Z., Zhang B., Hong W., Wu Y. Fast compressed sensing SAR imaging based on approximated observation. IEEE J. Sel. Top. Appl. Earth Obs. Remote Sens. 2014;7:352–363. doi: 10.1109/JSTARS.2013.2263309. [Cross Ref]

10. Jiang C., Zhang B., Fang J., Zhe Z., Hong W., Wu Y., Xu Z. Efficient l q regularisation algorithm with range-azimuth decoupled for SAR imaging. Electron. Lett. 2014;50:204–205. doi: 10.1049/el.2013.1989. [Cross Ref]

11. Yang J., Huang X., Jin T., Thompson J., Zhou Z. Synthetic aperture radar imaging using stepped frequency waveform. IEEE Trans. Geosci. Remote Sens. 2012;50:2026–2036. doi: 10.1109/TGRS.2011.2170176. [Cross Ref]

12. Cumming I.G., Wong F.H. Digital Processing of Synthetic Aperture Radar Data: Algorithms and Implementation. Artech House; Norwood, MA, USA: 2005.

13. Yang J., Thompson J., Huang X., Jin T., Zhou Z. Random-frequency SAR imaging based on compressed sensing. IEEE Trans. Geosci. Remote Sens. 2013;51:983–994. doi: 10.1109/TGRS.2012.2204891. [Cross Ref]

14. Ding L., Chen W., Zhang W., Poor H.V. MIMO radar imaging with imperfect carrier synchronization: A point spread function analysis. IEEE Trans. Aerosp. Electron. Syst. 2015;51:2236–2247. doi: 10.1109/TAES.2015.140428. [Cross Ref]

15. Ding L., Liu C., Wang T., Chen W. Sparse self-calibration via iterative minimization against phase synchronization mismatch for MIMO radar imaging; Proceedings of the 2013 IEEE Radar Conference (RadarCon13); Ottawa, ON, Canada. 29 April–3 May 2013; pp. 1–4.

16. Li B., Liu F., Zhou C., Lv Y., Hu J. Fast compressed sensing SAR imaging using stepped frequency waveform; Proceedings of the 2016 IEEE International Conference on Microwave and Millimeter Wave Technology (ICMMT); Beijing, China. 5–8 June 2016; pp. 521–523.

17. Hu J., Liu F., Zhou C., Li B., Wang D. CS-SAR Imaging Method Based on Inverse Omega-K Algorithm. J. Radars. 2016 doi: 10.12000/JR16027. [Cross Ref]

18. Mairal J., Bach F., Ponce J., Sapiro G., Zisserman A. Non-local sparse models for image restoration; Proceedings of the 2009 IEEE 12th International Conference on Computer Vision; Kyoto, Japan. 27 September–4 October 2009; pp. 2272–2279.

19. Raney R.K., Runge H., Bamler R., Cumming I.G., Wong F.H. Precision SAR processing using chirp scaling. IEEE Trans. Geosci. Remote Sens. 1994;32:786–799. doi: 10.1109/36.298008. [Cross Ref]

20. Blumensath T., Davies M.E. Iterative hard thresholding for compressed sensing. Appl. Comput. Harmonic Anal. 2009;27:265–274. doi: 10.1016/j.acha.2009.04.002. [Cross Ref]

21. Boyd S., Vandenberghe L. Convex Optimization. Cambridge University Press; New York, NY, USA: 2004.

22. Wang T., Lu X., Yu X., Xi Z., Chen W. A Fast and Accurate Sparse Continuous Signal Reconstruction by Homotopy DCD with Non-Convex Regularization. Sensors. 2014;14:5929–5951. doi: 10.3390/s140405929. [PMC free article] [PubMed] [Cross Ref]

Articles from Sensors (Basel, Switzerland) are provided here courtesy of **Multidisciplinary Digital Publishing Institute (MDPI)**

PubMed Central Canada is a service of the Canadian Institutes of Health Research (CIHR) working in partnership with the National Research Council's national science library in cooperation with the National Center for Biotechnology Information at the U.S. National Library of Medicine(NCBI/NLM). It includes content provided to the PubMed Central International archive by participating publishers. |