High Resolution Imaging from Azimuth Missing SAR Raw Data via Segmented Recovery

: Synthetic Aperture Radar (SAR) raw data missing occurs when radar is interrupted by various inﬂuences. In order to cope with this problem, a new method is proposed to focus the azimuth missing SAR raw data via segmented recovery in this paper. A reference function in time domain is designed to make the missing raw data sparser in two dimensional frequency domain. Afterwards, greedy algorithms are available to recover the missing data in two dimensional frequency domain. In addition, in order to avoid range frequency aliasing problem caused by reference function multiplication in time domain, the missing raw data is split into several parts in range direction and is recovered with a segmented recovery strategy. Then, the recovered raw data is available to be focused with traditional SAR imaging algorithms. The range migration algorithm is chosen to deal with the recovered raw data in this paper. Point target and area target simulations are carried out to validate the effectiveness of the proposed method on azimuth missing SAR raw data. Moreover, the proposed method is implemented on real SAR data in order to further provide convincing demonstration.


Introduction
Synthetic Aperture Radar (SAR) is an active sensor that transmits a beam of electromagnetic radiation in order to extend human ability to observe properties of the surface of the earth [1]. As an effective means for remote sensing, SAR is capable of providing high resolution images and videos for target surveillance under any weather condition [2][3][4][5][6][7]. With the development of SAR, robustness to missing raw data raises concern due to sparse and irregular sampling in SAR raw data caused by new mission requirements, sensor geometries, jamming, interference, and so on [8]. When dealing with missing raw data, traditional signal processing means often perform unsatisfactorily. As a result, new methods are desired to be developed in order to acquire focused image from missing SAR raw data.
To overcome the problem of processing missing SAR raw data, several methods have been proposed to accommodate missing data. A series of spectral analysis methods have been proposed to deal with missing data. The Gapped-Data Amplitude and Phase Estimation (GAPES) is a well-known spectral analysis method developed for applying to missing data [9]. The GAPES is the extension of APES [10] and is based on minimizing the least square criterion with respect to the missing data. The GAPES can estimate the amplitude spectrum from incomplete data via cyclic minimization and conjugate gradient algorithm. GAPES uses the available part of data to interpolate the missing part of data by utilizing least square criterion iteratively [11]. However, the GAPES and its subsequent algorithms, such as Missing Amplitude and Phase Estimation Expectation Maximization (MAPES-EM) and Missing Amplitude and Phase Estimation Cyclic Maximization (MAPES-CM), have large computational complexity [12]. Therefore, the GAPES is not suitable for large scale data processing and high resolution imaging. Afterwards, a spectral analysis method called Iterative Adaptive Approach (IAA) is proposed [13] and is applied to blood velocity estimation [14] and SAR imaging [15,16] on missing data. Meanwhile, comparing with GAPES and its subsequent algorithms, IAA is competitive in computing cost.
In addition to spectral analysis methods, the sparse optimization can be also utilized to deal with missing SAR or ISAR data [17][18][19][20]. In [17], a novel centralized sparse representation based SAR imaging method is proposed to deal with sparse sampling SAR data. In [18], a modified compressive sensing method for wide angle SAR sub-aperture imaging is presented. In [19], a method is proposed to cope with azimuth periodically gapped raw data with greedy method. In [20], a novel super resolution imaging algorithm with relaxation technique for sparse aperture ISAR (Inverse SAR) is presented. Sparse optimization is capable of improving imaging quality on missing data and is usually functioned as recovering the complete data in a certain domain. In a sense, recovering the SAR raw data in a certain domain can be treated as solving underdetermined systems of linear equations. Orthogonal Matching Pursuit (OMP), a greedy algorithm, can provide approximate sparse solution of underdetermined systems and reliably recover a sparse signal from random linear measurement of original signal [21,22]. In contrast to OMP, Stagewise Orthogonal Matching Pursuit (StOMP), which is developed from OMP, is able to recover many coefficients at each iteration while OMP only recovers one coefficient at each iteration [23]. As a result, StOMP is faster than OMP and has potential on the application for signal recovery in SAR imaging.
In order to achieve SAR high resolution imaging from azimuth missing SAR raw data, a new strategy is proposed to recover raw data via segmented recovery in this paper. Inspired by Range Migration Algorithm (RMA) [24], a reference function in time domain is designed to compensate the phase of missing data so that the energy of data is concentrated in two dimensional frequency domain. In other words, after phase compensation, the missing data becomes much sparser than before in two dimensional frequency. However, the reference function in time domain is possible to lead to aliasing in range frequency when range swath is greater than a certain threshold value. In order to avoid this problem, the missing raw data is split into several parts in range direction and is recovered via a segmented recovery strategy with StOMP. Then, conventional SAR imaging algorithms, such as RMA, can process the recovered raw data so as to acquire finely focused image. In experiments, point and area target simulations are implemented to validate the effectiveness of the proposed method on azimuth missing raw data. The data missing rate is 50%. The resolution in point target simulation can reach 0.25 m in both range and azimuth direction. Furthermore, real SAR raw data is utilized to further verify the effectiveness of the proposed method.
Vectors and matrices are in bold italic while variables are in italic in this paper. This paper is organized as follows. In Section 2, the proposed method for imaging azimuth missing SAR raw data via segmented recovery is described. The designed reference function and missing data recovery strategy are presented in details. Section 3 shows the experimental results obtained via the proposed method. Point target simulation, area target simulation and real SAR data are utilized to verify the effectiveness of the proposed method. In Section 4, discussion and analysis of experimental results are presented. Finally, conclusion is drawn in Section 5.

Reference Function in Time Domain
It is difficult to recover missing SAR raw data in time domain via a certain method as SAR raw data is dense in time domain. In order to make it possible for processing missing SAR data, exploiting sparsity is considered in this paper. In this part, a reference function in time domain is designed for phase compensation in time domain. After phase compensation, the SAR raw data may become sparser than before in two dimensional (2D) frequency domain. In this case, the missing SAR raw data can be recovered in 2D frequency domain with the application of greedy algorithms or other proper sparse optimization algorithms. The description of the proposed reference function is derived and presented as follows.
SAR is supposed to transmit pulse chirp signal and move uniformly in a straight line in this paper. After demodulation to baseband, the echo signal of single point can be expressed in terms of the range time τ and azimuth time η as follows: Here, the amplitude factors have been ignored. w r is the range envelope, w a is the azimuth envelope, c is the velocity of light, K r is the range frequency modulation rate, R(η) is the slant range between radar and target at η, f 0 is the carrier frequency, and j is the imaginary unit.
The range equation can be considered as the parabolic approximation via Taylor expansion as follows: Here, v is the velocity of the radar platform and R 0 is the closest slant range. Substituting the expression of slant range of (2) into (1), the baseband received signal can be approximated as follows: where λ is the wavelength. In (3), the second exponential term shows azimuth modulation while the third exponential term represents range modulation.
Inspired by reference function in 2D frequency domain of RMA, a reference function is designed to match the raw data in time domain. The reference function θ(τ,η) is defined as follows: Here, R 0,ref is the closest reference slant range. The reference slant range R ref (η) at η is defined as follows: To match the raw data with reference function in time domain, the Equation (3) is multiplied by Equation (4). Then, the result g(τ,η) can be expressed as follows: As shown in Equation (6), the second exponential term indicates that the bandwidth of azimuth modulation is mostly reduced after reference function multiplication. Moreover, the third exponential term represents that the quadratic term of range time in range modulation is removed via reference function multiplication. Therefore, bandwidth in range modulation is also mostly reduced. As a result, the g(τ,η) is supposed to be roughly compressed in 2D frequency domain. The rough compression in 2D frequency domain after reference function multiplication is presented in Figures 1 and 2. Rough compression of point targets and real SAR data are shown in Figures 1 and 2, respectively. In the demonstration of rough compression of point targets, nine point targets are placed for generating echo data. The distribution of nine points is illustrated in Figure 1a. Figure 1b shows the echo data in time domain. It is obvious that the echo data is dense in time domain as shown in Figure 1b. The Figure 1c depicts the 2D spectrum of the echo data. It is apparent that the echo data is also dense in 2D frequency domain. The 2D spectrum of g(τ,η) is presented in Figure 1d. As shown in Figure 1d, it is obvious that the spectrum of g(τ,η) is roughly compressed as the modulation of echo data in range and azimuth direction is partly reduced. In Figure 2, Figure 2a shows the echo data of real SAR data in time domain and Figure 2b illustrates 2D spectrum of real SAR data. Meanwhile, 2D spectrum of g(τ,η) of real SAR data is presented in Figure 2c. As shown in Figure 2c, it is apparent that the 2D spectrum of g(τ,η) of real SAR data is also roughly compressed.
In order to further validate the rough compression in 2D frequency domain after reference function multiplication, analysis of rough compression is listed in Table 1. Image contrast (IC) [25] and image entropy (IE) [26] are selected to evaluate the sparsity of echo data, 2D spectrum of echo data and 2D spectrum of g(τ,η). The sparser the SAR image is, the more concentrated the energy of SAR image is. And if the energy of SAR image becomes more concentrated, the IE of SAR images tends to be smaller and the IC of SAR image tends to be larger. As presented in Table 1, 2D spectrum of g(τ,η) has smallest IE and greatest IC in the scenario of point targets and real SAR data, comparing with echo data and 2D spectrum of echo data. Therefore, it can be considered that 2D spectrum of g(τ,η) is sparser than echo data and 2D spectrum of echo data. In conclusion, after reference function multiplication in time domain, the spectrum of data becomes much sparser. As shown in Equation (6), the second exponential term indicates that the bandwidth of azimuth modulation is mostly reduced after reference function multiplication. Moreover, the third exponential term represents that the quadratic term of range time in range modulation is removed via reference function multiplication. Therefore, bandwidth in range modulation is also mostly reduced. As a result, the g(τ,η) is supposed to be roughly compressed in 2D frequency domain.
The rough compression in 2D frequency domain after reference function multiplication is presented in Figures 1 and 2. Rough compression of point targets and real SAR data are shown in Figures 1 and 2, respectively. In the demonstration of rough compression of point targets, nine point targets are placed for generating echo data. The distribution of nine points is illustrated in Figure 1a. Figure 1b shows the echo data in time domain. It is obvious that the echo data is dense in time domain as shown in Figure 1b. The Figure 1c depicts the 2D spectrum of the echo data. It is apparent that the echo data is also dense in 2D frequency domain. The 2D spectrum of g(τ,η) is presented in Figure 1d. As shown in Figure 1d, it is obvious that the spectrum of g(τ,η) is roughly compressed as the modulation of echo data in range and azimuth direction is partly reduced. In Figure 2, Figure 2a shows the echo data of real SAR data in time domain and Figure 2b illustrates 2D spectrum of real SAR data. Meanwhile, 2D spectrum of g(τ,η) of real SAR data is presented in Figure 2c. As shown in Figure 2c, it is apparent that the 2D spectrum of g(τ,η) of real SAR data is also roughly compressed.  In order to further validate the rough compression in 2D frequency domain after reference function multiplication, analysis of rough compression is listed in Table 1. Image contrast (IC) [25] and image entropy (IE) [26] are selected to evaluate the sparsity of echo data, 2D spectrum of echo data and 2D spectrum of g(τ,η). The sparser the SAR image is, the more concentrated the energy of SAR image is. And if the energy of SAR image becomes more concentrated, the IE of SAR images tends to be smaller and the IC of SAR image tends to be larger. As presented in Table 1, 2D spectrum of g(τ,η) has smallest IE and greatest IC in the scenario of point targets and real SAR data, comparing with echo data and 2D spectrum of echo data. Therefore, it can be considered that 2D spectrum of g(τ,η) is sparser than echo data and 2D spectrum of echo data. In conclusion, after reference function multiplication in time domain, the spectrum of data becomes much sparser.

Segmented Recovery
As described in Section 2.1, the raw data is sparser in 2D frequency domain after reference function multiplication in time domain. Therefore, consideration can be taken to recover missing raw data in 2D frequency domain of g(τ,η). The 2D spectrum of g(τ,η) is defined as G(fτ,fη), fτ is the range frequency and fη is azimuth frequency.
As shown in the third exponential term of Equation (6), the quadratic term of range time in range modulation is removed by reference function multiplication. As a result, only a linear phase term h(τ,η) is needed to be taken into consideration in each pulse at azimuth time η, when g(τ,η) is transformed into range frequency domain. The expression of h(τ,η) is given as follows:

Segmented Recovery
As described in Section 2.1, the raw data is sparser in 2D frequency domain after reference function multiplication in time domain. Therefore, consideration can be taken to recover missing raw data in 2D frequency domain of g(τ,η). The 2D spectrum of g(τ,η) is defined as G(f τ ,f η ), f τ is the range frequency and f η is azimuth frequency.
As shown in the third exponential term of Equation (6), the quadratic term of range time in range modulation is removed by reference function multiplication. As a result, only a linear phase term h(τ,η) is needed to be taken into consideration in each pulse at azimuth time η, when g(τ,η) is transformed into range frequency domain. The expression of h(τ,η) is given as follows: There exists a limitation for the difference between R(η) and R ref (η) because the range sampling frequency is finite. Such limitation can be expressed as follows: F s is the range sampling frequency and |•| denotes absolute value operation. The Equation (8) indicates that the reference function multiplication in time domain may lead to aliasing in range frequency when the absolute value of the difference between R(η) and R ref (η) is greater than cF s /(4K r ). However, the swath width requirement in range direction usually causes that the absolute value of the difference between R(η) and R ref (η) is greater than cF s /(4K r ). As a result, range frequency aliasing exists after the reference function multiplication in time domain. In order to cope with this problem, it is proper to divide the raw data into several parts in range direction according to the constraint in (8) and match with corresponding reference function in time domain according to different reference slant ranges.
The complete raw data is defined as an N r × N a matrix with N r sampling points in range direction and N a sampling pulses in azimuth direction. The set, which is denoted as E, contains indices of locations of available data in missing echo data. And the number of elements in set E is defined as N E . Comparing with complete SAR raw data, the corresponding set for available azimuth sampling time in missing data is defined as Q E . Figure 3 is the diagram of SAR azimuth missing data. The white square denotes complete data sampling and the shadow square denotes missing data sampling. The complete data sampling has N a pulses in azimuth direction and the missing data sampling has N E pulses in azimuth direction. In this paper, the proposed algorithm aims to deal with the SAR data which is randomly or irregular missing in azimuth direction. Therefore, the distribution of sampling positions of unavailable data in missing SAR data is random distribution. Meanwhile, the distribution of available data's sampling positions in missing SAR data is also random distribution.
Fs is the range sampling frequency and |•| denotes absolute value operation. The equation (8) indicates that the reference function multiplication in time domain may lead to aliasing in range frequency when the absolute value of the difference between R(η) and Rref(η) is greater than cFs/(4Kr). However, the swath width requirement in range direction usually causes that the absolute value of the difference between R(η) and Rref(η) is greater than cFs/(4Kr). As a result, range frequency aliasing exists after the reference function multiplication in time domain. In order to cope with this problem, it is proper to divide the raw data into several parts in range direction according to the constraint in (8) and match with corresponding reference function in time domain according to different reference slant ranges.
The complete raw data is defined as an Nr × Na matrix with Nr sampling points in range direction and Na sampling pulses in azimuth direction. The set, which is denoted as E, contains indices of locations of available data in missing echo data. And the number of elements in set E is defined as NE. Comparing with complete SAR raw data, the corresponding set for available azimuth sampling time in missing data is defined as QE. Figure 3 is the diagram of SAR azimuth missing data. The white square denotes complete data sampling and the shadow square denotes missing data sampling. The complete data sampling has Na pulses in azimuth direction and the missing data sampling has NE pulses in azimuth direction. In this paper, the proposed algorithm aims to deal with the SAR data which is randomly or irregular missing in azimuth direction. Therefore, the distribution of sampling positions of unavailable data in missing SAR data is random distribution. Meanwhile, the distribution of available data's sampling positions in missing SAR data is also random distribution.  Figure 4 illustrates the strategy of segmented recovery in this paper. Firstly, the azimuth missing raw data is divided into Np patches in range direction and every patch has the same number of range sampling point. Patches are labelled as Pp, where p = 1, 2, 3, …, Np. Then, each patch is matched with reference function in time domain, respectively. The corresponding reference function for each patch is defined as follows: As raw data is divided into Np patches, Np closest reference slant ranges, which are denoted as R0,ref,p, are selected according to the limitation in Equation (8). Meanwhile, Np reference slant ranges can be expressed as follows:  Figure 4 illustrates the strategy of segmented recovery in this paper. Firstly, the azimuth missing raw data is divided into N p patches in range direction and every patch has the same number of range sampling point. Patches are labelled as P p , where p = 1, 2, 3, . . . , N p . Then, each patch is matched with reference function in time domain, respectively. The corresponding reference function for each patch is defined as follows: As raw data is divided into N p patches, N p closest reference slant ranges, which are denoted as R 0,ref,p , are selected according to the limitation in Equation (8). Meanwhile, N p reference slant ranges can be expressed as follows: Moreover, the τ p is defined as follows: In (11), N s denotes the number of range sampling point of each patch. In addition, the reference function for each patch of the missing data is obtained by selecting columns from the reference function in (9) according to the set Q E .
Moreover, the τp is defined as follows:  In (11), Ns denotes the number of range sampling point of each patch. In addition, the reference function for each patch of the missing data is obtained by selecting columns from the reference function in (9) according to the set QE.
Each patch of complete raw data after reference function multiplication is denoted as g(τp,η;Pp). The 2D spectrum of g(τp,η;Pp) is denoted as G(fτp,fη; Pp). The range frequency in G(fτp,fη; Pp) is denoted as fτp. Each patch of azimuth missing raw data after reference function multiplication is denoted as g(τp,η;Pp)|η ∈ QE. After being multiplied by reference function and operating range Fourier transform, the patch of azimuth missing raw data can be denoted as g(fτp,η;Pp)|η∈QE. The NE×1 data vector y is defined as the vector for any range frequency fτp in g(fτp,η)|η∈QE, and the Na×1 vector x is defined as the vector for any range frequency fτp in G(fτp,fη; Pp). If Ψ denotes partial inverse Fourier transform matrix, the y and x can be connected via Ψ as follows: The partial inverse Fourier transform matrix Ψ can be generated by extracting rows from inverse Fourier transform square matrix according to the set E.
As NE < Na, the Equation (12) is an underdetermined system of linear equations. StOMP is suitable to deal with this situation. The unknown x can be solved with known y and Ψ via StOMP. StOMP operates in S stages, obtaining a sequence of estimations x0, x1, … by removing acquired structure from a sequence of residual vectors r1, r2, ….
StOMP begins with initial solution x0 = 0 and initial residual r0 = y. The stage counter s begins at s = 1. A sequence I1, …, Is records locations of the nonzeros in x. The s-th stage applies matched filtering to the current residual, obtaining a vector of residual correlations cs as follows: The superscript H denotes the conjugate transpose operation. Next, the procedure implements thresholding to seek the significant nonzeros, and thresholding yields a set Us for chosen location coordinates as follows: Each patch of complete raw data after reference function multiplication is denoted as g(τ p ,η;P p ). The 2D spectrum of g(τ p ,η;P p ) is denoted as G(f τp ,f η ; P p ). The range frequency in G(f τp ,f η ; P p ) is denoted as f τp . Each patch of azimuth missing raw data after reference function multiplication is denoted as g(τ p ,η;P p )|η∈Q E . After being multiplied by reference function and operating range Fourier transform, the patch of azimuth missing raw data can be denoted as g(f τp ,η;P p )|η∈Q E . The N E ×1 data vector y is defined as the vector for any range frequency f τp in g(f τp ,η)|η∈Q E , and the N a ×1 vector x is defined as the vector for any range frequency f τp in G(f τp ,f η ; P p ). If Ψ denotes partial inverse Fourier transform matrix, the y and x can be connected via Ψ as follows: The partial inverse Fourier transform matrix Ψ can be generated by extracting rows from inverse Fourier transform square matrix according to the set E.
As N E < N a , the Equation (12) is an underdetermined system of linear equations. StOMP is suitable to deal with this situation. The unknown x can be solved with known y and Ψ via StOMP. StOMP operates in S stages, obtaining a sequence of estimations x 0 , x 1 , . . . by removing acquired structure from a sequence of residual vectors r 1 , r 2 , . . . .
StOMP begins with initial solution x 0 = 0 and initial residual r 0 = y. The stage counter s begins at s = 1. A sequence I 1 , . . . , I s records locations of the nonzeros in x. The s-th stage applies matched filtering to the current residual, obtaining a vector of residual correlations c s as follows: The superscript H denotes the conjugate transpose operation. Next, the procedure implements thresholding to seek the significant nonzeros, and thresholding yields a set U s for chosen location coordinates as follows: Here σ s is a formal noise level and t s is a threshold parameter. Then, the new location estimation is updated by taking the union of newly selected coordinates and the previous location estimation as follows: Afterwards, the new approximation x s supported in I s given by: The updated residual is: Then the stopping condition is checked. If it is not time to stop, the procedure goes to next iteration and stage counter s is set as: s = s + 1. If it is time to for stopping, the final result of the procedure is set as:x S = x s .
As a result, spectrum of patch can be recovered and is denoted asĜ( f τ p , f η ; P p ). Then, the recovered patch of complete SAR raw data in time domain is obtained as follows: Here, IFFT 2D [·] denotes 2D inverse Fourier transform operation and conj[·] denotes conjugate operation. Consequently, the recovered patch in time domain is recombined to form the recovered whole complete raw data. Afterwards, conventional SAR imaging algorithms can process the recovered SAR raw data in order to acquire focused image. Figure 5 is the flow chart of the proposed method. In the flow chart, FFT denotes fast Fourier transform and 2D IFFT denotes two dimensional inverse FFT. The proposed method can be generally divided into two parts: segmented recovery and recovered raw data focusing. In the part of segmented recovery, the azimuth missing raw data is firstly divided into N p patches in range direction. Then, each patch of missing raw data is matched with reference function in time domain and recovered in 2D frequency domain with StOMP. Afterwards, the recovered spectrum of each patch is transformed into time domain and multiplied by the conjugate reference function. Consequently, the recovered each patch in time domain is recombined to form the whole complete raw data. After the part of segmented recovery, RMA is implemented to focus the recovered raw data. The final focused image is obtained in the part of recovered raw data focusing.

Experiments
Experiments are carried out to validate the effectiveness of the proposed method. In this section, firstly, point target simulation is utilized to demonstrate the performance of the proposed method. Then, area target simulation is also implemented to verify the validity of the proposed method. Finally, real SAR data is used to further validate the effectiveness of the proposed method.

Point Target Simulation
Aiming to evaluate the effectiveness of the proposed method, point target simulation is

Experiments
Experiments are carried out to validate the effectiveness of the proposed method. In this section, firstly, point target simulation is utilized to demonstrate the performance of the proposed method. Then, area target simulation is also implemented to verify the validity of the proposed method. Finally, real SAR data is used to further validate the effectiveness of the proposed method.

Point Target Simulation
Aiming to evaluate the effectiveness of the proposed method, point target simulation is implemented in this section. The point target simulation is under the consideration of monostatic SAR with transmitting pulse chirp signal. The system parameters for simulation are listed in Table 2. The resolution in point target simulation is desired to achieve 0.25 m in both range and azimuth direction. The mode of SAR operation in point target simulation is select as spotlight SAR according to the resolution requirement.  Figure 6 shows the distribution of nine point targets for simulation and targets are labelled from T 1 to T 9 . The complete SAR raw data for point target simulation is a 5120 × 4096 matrix with 5120 and 4096 sampling points in range and azimuth direction.  Figure 6 shows the distribution of nine point targets for simulation and targets are labelled from T1 to T9. The complete SAR raw data for point target simulation is a 5120 × 4096 matrix with 5120 and 4096 sampling points in range and azimuth direction. The missing SAR raw data for simulation is generated by forming a new matrix via extracting columns from complete SAR raw data randomly. The missing rate in azimuth direction is 50%. Therefore, the size of missing SAR raw data matrix is 5120 × 2048 in point target simulation. The size of recovered SAR raw data matrix is 5120 × 4096.
In point target simulation, the maximum stage number and the threshold parameter are two parameters for StOMP in the proposed method. The maximum stage number is set as S = 20, and the threshold parameter is set as ts = 0.03. The formal noise level σs at each iteration can be calculated according to literature [23]    The missing SAR raw data for simulation is generated by forming a new matrix via extracting columns from complete SAR raw data randomly. The missing rate in azimuth direction is 50%. Therefore, the size of missing SAR raw data matrix is 5120 × 2048 in point target simulation. The size of recovered SAR raw data matrix is 5120 × 4096.
In point target simulation, the maximum stage number and the threshold parameter are two parameters for StOMP in the proposed method. The maximum stage number is set as S = 20, and the threshold parameter is set as t s = 0.03. The formal noise level σ s at each iteration can be calculated according to literature [23] as follows:  Figure 7, all the nine point targets are well focused. As a result, it is considered that the proposed method performs satisfactorily on point target simulation. In order to validate the effectiveness of the proposed method on point target simulation, the analysis for point target simulation is presented in Section 4.
In point target simulation, the maximum stage number and the threshold parameter are two parameters for StOMP in the proposed method. The maximum stage number is set as S = 20, and the threshold parameter is set as ts = 0.03. The formal noise level σs at each iteration can be calculated according to literature [23] as follows:

Area Target Simulation
In this part, the area target simulation is presented to verify the validity of the proposed method. A large number of point targets are positioned to form an area target with the shape of capital letter 'T'. The area target simulation is implemented with the system parameters in Table 2. The mode of SAR operation in point target simulation is select as spotlight SAR according to the

Area Target Simulation
In this part, the area target simulation is presented to verify the validity of the proposed method. A large number of point targets are positioned to form an area target with the shape of capital letter 'T'. The area target simulation is implemented with the system parameters in Table 2. The mode of SAR operation in point target simulation is select as spotlight SAR according to the resolution requirement. The complete SAR raw data generated for area target simulation is a 5120 × 4096 matrix, which has 5120 and 4096 sampling points in range and azimuth direction, respectively. The missing rate in azimuth direction for area simulation is also 50%. Consequently, the size of missing SAR raw data matrix is 5120 × 2048 in area target simulation. In the application of the proposed method on area target simulation, the maximum stage number is set as S = 20, and the threshold parameter is set as t s = 0.03. Results of area target simulation are displayed in the following part. Figure 8 depicts simulation results of area target. Figure 8a is the imaging result obtained from complete raw data via RMA. Figure 8b is the imaging result obtained via filling missing positions in azimuth direction with zeroes and subsequent application of RMA. Figure 8c is the imaging result obtained via GAPES. The GAPES for acquiring result in Figure 8c is performed by replacing the StOMP with GAPES in the proposed method. Figure 8d is the imaging result obtained via IAA. The IAA for acquiring result in Figure 8d is carried out by replacing the StOMP with IAA in [14] in the proposed method. Figure 8e is the imaging result obtained via the proposed method. As shown in Figure 8b, the phenomenon of data missing in azimuth direction leads to azimuth defocusing when missing raw data is processed with adding zeroes operation. Figure 8c-e indicate that the azimuth defocusing is obviously suppressed after segmented recovery for azimuth missing data. missing raw data is processed with adding zeroes operation. Figure 8c-e indicate that the azimuth defocusing is obviously suppressed after segmented recovery for azimuth missing data.

Real Data Processing
Real SAR data is utilized to validate the effectiveness of the proposed method in this part. The real SAR data is acquired by Sentinel-1A. The complete real SAR data is an 8192 × 8192 matrix, which has 8192 sampling points in both range and azimuth direction.
The missing rate in azimuth direction for area simulation is also 50%. Therefore, the size of missing SAR raw data matrix is 8192 × 4096 for real data processing. In implementation of the proposed method on real SAR data processing, the maximum stage number is set as S = 40, and the threshold parameter is set as ts = 0.0275.
In order to present the effects on 2D spectrum before and after the proposed algorithm is implemented, a part of real SAR data is selected for demonstration. To give a more apparent

Real Data Processing
Real SAR data is utilized to validate the effectiveness of the proposed method in this part. The real SAR data is acquired by Sentinel-1A. The complete real SAR data is an 8192 × 8192 matrix, which has 8192 sampling points in both range and azimuth direction.
The missing rate in azimuth direction for area simulation is also 50%. Therefore, the size of missing SAR raw data matrix is 8192 × 4096 for real data processing. In implementation of the proposed method on real SAR data processing, the maximum stage number is set as S = 40, and the threshold parameter is set as t s = 0.0275.
In order to present the effects on 2D spectrum before and after the proposed algorithm is implemented, a part of real SAR data is selected for demonstration. To give a more apparent demonstration of the difference of 2D spectrums of complete data, missing data with adding zeroes and recovered data, the echo data is multiplied with reference function in time domain before being transformed into 2D frequency domain. The recovered data is obtained via the proposed algorithm in this paper. The utilized reference function in time domain is described in Section 2.1.
With the reference function multiplication in time domain, 2D spectrums of complete data, missing data with adding zeroes and recovered data are presented in Figure 9a-c, respectively. In addition, one dimensional (1D) azimuth spectrum profiles of complete data, missing data and recovered data are presented in Figure 9d-f, respectively. The row, which is at 1.35 MHz in the range frequency of each 2D spectrum, is selected to illustrate the 1D azimuth spectrum profile of each 2D spectrum. The 1D azimuth spectrum profiles are normalized with the corresponding maximum amplitude of each 1D spectrum profile. As shown in Figure 9b, the aliasing phenomenon exists in the 2D spectrum of missing data in azimuth direction. And as presented in Figure 9c, the aliasing phenomenon has been almost removed in the 2D spectrum of recovered data in azimuth direction. The 1D azimuth spectrum profiles in Figure 9d-f also validate that aliasing phenomenon is almost solved after the proposed As shown in Figure 9b, the aliasing phenomenon exists in the 2D spectrum of missing data in azimuth direction. And as presented in Figure 9c, the aliasing phenomenon has been almost removed in the 2D spectrum of recovered data in azimuth direction. The 1D azimuth spectrum profiles in Figure 9d-f also validate that aliasing phenomenon is almost solved after the proposed algorithm is implemented. Figure 10 shows results of real SAR data processing. Figure 10a is the imaging result acquired from complete raw data via RMA. Figure 10b is the imaging result acquired via filling missing positions in azimuth direction with zeroes and subsequent application of RMA. Figure 10c is the imaging result acquired via GAPES. Figure 10d is the imaging result acquired via IAA. Figure 10e is the imaging result acquired via the proposed method.
As shown in Figure 10, the phenomenon of azimuth defocusing exists in Figure 10b-d. In other words, despite satisfactory performance on area target simulation, the IAA performs unsatisfactorily on real SAR data. Figure 10e indicates that the azimuth defocusing is apparently suppressed via the proposed method. In this case, the effectiveness of the proposed method on real SAR data is verified in Figure 10e. As shown in Figure 10, the phenomenon of azimuth defocusing exists in Figure 10b-d. In other words, despite satisfactory performance on area target simulation, the IAA performs unsatisfactorily on real SAR data. Figure 10e indicates that the azimuth defocusing is apparently suppressed via the proposed method. In this case, the effectiveness of the proposed method on real SAR data is verified in Figure 10e.
In order to demonstrate the effectiveness of the proposed method further, three local regions are selected for displaying more details. The selected three regions are shown in Figure 11, which are identified by red rectangles. The enlarged local regions of results of real SAR data processing are shown from Figure 12 to Figure 14. The presented enlarged local regions are obtained via five different approaches. The five approaches are acquiring result image from complete data via RMA, acquiring result image from missing data via adding zeroes, acquiring result image from missing data via GAPES, acquiring result image from missing data via IAA and acquiring result image from missing data via the proposed method. The effectiveness of the proposed method has been further verified by Figures 12-14. In order to demonstrate the effectiveness of the proposed method further, three local regions are selected for displaying more details. The selected three regions are shown in Figure 11

Discussion
Results presented in Section 3 validate the effectiveness of the proposed method. So as to further confirm the validity of the proposed method, some criteria are utilized to analyze results of point target simulation, area target simulation and real SAR data processing.
In order to confirm the effectiveness of the proposed method on point target simulation, impulse response width (IRW), peak sidelobe ratio (PSLR) and integrated sidelobe ratio (ISLR) are selected as criteria for evaluating the quality of imaging results. In order to ensure the focused image quality, the PSLR and ISLR should be less than −13 dB and −10.15 dB, respectively. Moreover, as IRW is referred to as the image resolution in SAR image, the IRW should be superior to 0.25 m in both range and azimuth direction when resolution is desired to achieve 0.25 m. The main lobe width for ISLR is defined as two times the IRW in this paper.
The analysis of point target simulation is presented in Table 3. As shown in Table 3, IRW, PSLR and ISLR of all the focused point targets satisfy the requirements of criteria in both range and azimuth direction, respectively. As a result, the proposed method can be deemed to perform satisfactorily on point target simulation. In order to validate the effectiveness of the proposed method on area target simulation further, image contrast (IC) [25] and image entropy (IE) [26] are chosen to evaluate the quality of focused results [27]. In [25], it is proposed that better focusing leads to larger IC. Meanwhile, it is proposed that better focusing results in smaller entropy in [26]. The analysis for area target simulation is listed

Discussion
Results presented in Section 3 validate the effectiveness of the proposed method. So as to further confirm the validity of the proposed method, some criteria are utilized to analyze results of point target simulation, area target simulation and real SAR data processing.
In order to confirm the effectiveness of the proposed method on point target simulation, impulse response width (IRW), peak sidelobe ratio (PSLR) and integrated sidelobe ratio (ISLR) are selected as criteria for evaluating the quality of imaging results. In order to ensure the focused image quality, the PSLR and ISLR should be less than −13 dB and −10.15 dB, respectively. Moreover, as IRW is referred to as the image resolution in SAR image, the IRW should be superior to 0.25 m in both range and azimuth direction when resolution is desired to achieve 0.25 m. The main lobe width for ISLR is defined as two times the IRW in this paper.
The analysis of point target simulation is presented in Table 3. As shown in Table 3, IRW, PSLR and ISLR of all the focused point targets satisfy the requirements of criteria in both range and azimuth direction, respectively. As a result, the proposed method can be deemed to perform satisfactorily on point target simulation. In order to validate the effectiveness of the proposed method on area target simulation further, image contrast (IC) [25] and image entropy (IE) [26] are chosen to evaluate the quality of focused results [27]. In [25], it is proposed that better focusing leads to larger IC. Meanwhile, it is proposed that better focusing results in smaller entropy in [26]. The analysis for area target simulation is listed in Table 4. As shown in Table 4, comparing with result obtained via adding zeroes operation and result obtained via GAPES, result obtained via IAA performs better on IE and IC. Meanwhile, result obtained via the proposed method also performs better on IE and IC than results obtained with adding zeroes operation and GAPES. In addition, the IE and IC of results acquired by IAA and the proposed method are close to the result obtained from complete SAR data. Therefore, both IAA and the proposed method perform satisfactorily on area target simulation. The analysis of area target simulation is consistent with the results in Figure 8. To give analysis for results of real data processing, IE and IC are utilized to evaluate the imaging quality of real SAR data results. Analysis for results of real SAR data processing is listed in Table 5. Both whole results and regions of results are analyzed in Table 5. As shown in Table 5, the proposed method generally performs better than adding zeroes operation, GAPES, and IAA with the criterion of IC. Due to existence of some fake scattering points, the GAPES has larger IC than the proposed method on Region 3. However, the proposed method generally performs better than GAPES. The difference among results obtained via adding zeroes operation, GAPES, IAA, and the proposed method is not obvious on the criterion of IE. The analysis for results of real data processing is generally consistent with the results shown in Figures 9 and 11-13.

Conclusions
A method for focusing SAR azimuth missing raw data is presented in this paper. A reference function is designed to match the missing data in time domain so as to concentrate the energy of data in two dimensional frequency domain. Then, the proposed segmented recovery strategy can recover the missing SAR raw data via StOMP. Afterwards, the recovered raw data can be focused by traditional SAR imaging algorithm, such as RMA. The effectiveness of the proposed method has been verified by point target simulation, area target simulation and real SAR data processing. In addition, analysis for results of experiments is also presented to further validate the effectiveness of the proposed method.