Sparse Aperture InISAR Imaging via Sequential Multiple Sparse Bayesian Learning

Interferometric inverse synthetic aperture radar (InISAR) imaging for sparse-aperture (SA) data is still a challenge, because the similarity and matched degree between ISAR images from different channels are destroyed by the SA data. To deal with this problem, this paper proposes a novel SA–InISAR imaging method, which jointly reconstructs 2-dimensional (2-D) ISAR images from different channels through multiple response sparse Bayesian learning (M-SBL), a modification of sparse Bayesian learning (SBL), to achieve sparse recovery for multiple measurement vectors (MMV). We note that M-SBL suffers a heavy computational burden because it involves large matrix inversion. A computationally efficient M-SBL is proposed, which, proceeding in a sequential manner to avoid the time-consuming large matrix inversion, is denoted as sequential multiple sparse Bayesian learning (SM-SBL). Thereafter, SM-SBL is introduced to InISAR imaging to simultaneously reconstruct the ISAR images from different channels. Numerous experimental results validate that the proposed SM-SBL-based InISAR imaging algorithm performs superiorly against the traditional single-channel sparse-signal recovery (SSR)-based InISAR imaging methods in terms of noise suppression, outlier reduction and 3-dimensional (3-D) geometry estimation.


Introduction
With the ability to acquire high-resolution radar images of moving targets, inverse synthetic aperture radar (ISAR) has been used in various applications for both civil and military purposes. However, it can only capture the projected 2-dimensional (2-D) characteristics of the target, which largely limits its application. In particular, the advanced automatic target recognition (ATR) technique shows the need for additional features of targets to serve classification. Therefore, interferometric ISAR (InISAR) has been developed to obtain the 3-dimensional (3-D) information of the target. Differently to the ISAR system, which merely utilizes a single channel to transmit and receive signals, the InISAR system generally includes several specially located channels, with one as both a transmitter and receiver, and the rest as receivers only [1]. Then, the 2-D ISAR images from different channels are obtained, and the phase difference between these are extracted to construct the 3-D geometry of the target. Additionally, the 3-D motion parameter of the target can be estimated during InISAR imaging. This paper mainly focuses on InISAR imaging for sparse-aperture (SA) data, which is still a challenge and is of practical significance.
Data are of SA-type when some parts are missing, which can be caused by either deliberate usage of compressive sensing (CS) to reduce the sampling rate, or by accidental effects from noise and jamming. Additionally, a multifunctional radar can also result in SA data, because it alternately transmits narrow-and wide-band signals to track and image the target simultaneously,

Signal Model
The geometry of the InISAR system is given in Figure 1 and is composed of three vertically positioned antennas, with the antenna O as both a transmitter and receiver, and the antennas A and B as receivers only. L 1 and L 2 denote the baselines OA and OB, respectively. Without loss of generality, we suppose that the target is located in the vertical direction of the plane AOB. The coordinate xoy is built on the center of gravity of the target, with the axes ox and oz parallel to OB and OA respectively, and oy along the line of sight (LOS) of the antenna O; p x p , y p , z p is an arbitrary scatterer on the target, whose projections on planes xoy and yoz are p and p , respectively. R O,p , R A,p and R B,p denote the instantaneous distances of p from three antennas, and α p and β p are the azimuth and elevation angle of p observed from O, respectively. Then, the range-compressed radar echo from p is denoted as wheret and t m denote the fast and slow time, respectively; f c , B and c are the center frequency, bandwidth, and propagation speed of the signal, respectively; and σ p represents the reflection coefficient of p. It should be noted that the instantaneous distances are related to the slow time t m and are shortened as R O,p , R A,p and R B,p for notational simplicity. These instantaneous distances can be derived as follows [10]: where R O,o denotes the instantaneous distance between the antenna O and the target gravity center o, and it also represents the translational motion of the target.R p denotes the rotational motion, and can be derived asR where rot y θ (·) denotes the obtaining of the y-coordinate after the inner vector rotates by an Euler angle of θ = [ω x t m ω y t m ω z t m ] T , and ω y , ω x and ω z represent the rotational speeds of the roll, pitch and yaw, respectively; r p = x p y p z p T is the location coordinate of p. We note that the coherent processing interval (CPI) for ISAR imaging is generally short; the target often only rotates a small angle within the CPI. Therefore, Equation (3) can be approximated with the first Taylor expansion as Substituting Equation (4) into Equation (1), and compensating for the translational motion [15,16] and migration throughout range cell (MTRC) [17], the range profile of p is derived as where S i t , t m , i = O, A, B denotes the translational compensated range profile from channel i. The range differences L 1 z p R o and L 2 x p R o in the envelope of S A t , t m and S B t , t m are eliminated, because they are much shorter than the range resolution c (2B). Then, the ISAR image of the scatterer p can be achieved by taking the fast Fourier transform (FFT) to the range profiles in Equation (5) with respect to the slow time as where G i,p t , f , i = O, A, B denotes the ISAR image from channel i. It is seen that x p and y p are embedded in the phase difference between three ISAR images and can be derived via the interferometric strategy as It should be noted that the instantaneous azimuth and elevation angle (i.e., α p and β p in Figure 1) are time-varying during the CPI and will cause a mismatch of three ISAR images. The time-varying parts of α p and β p should be compensated for to match the ISAR images before the process of interferometry, which is also called the ISAR image registration [18]. We now focus on InISAR imaging for SA data. When SSR is utilized to achieve ISAR imaging, the range profiles in Equation (5) and the ISAR images in Equation (6) are regarded as the observation and sparse coefficient, respectively [19]. The range profile in a single range cell can be discretized from Equation (5) as where P denotes the number of scatterers on the target, and the Sinc function in Equation (5) is discarded, because only the main lobe is considered here. Noting that three channels in the InISAR system given in Figure 1 are closely located, the ISAR images they obtain generally share a common sparse pattern, except for the mismatch caused by the time-varying azimuth and elevation angle.
Therefore, after compensating for the time-varying part of the azimuth and elevation angle [18], the ISAR images from different channels can be jointly reconstructed by the SSR algorithm for MMV. The range profiles from different channels are modeled in the form of MMV as where S ∈ C L×3 , w ∈ C K×3 and n ∈ C K×3 denote the combined range profile, ISAR image and noise, respectively, and L and K are the number of pulses and reconstructed Doppler cells, respectively. S i ∈ C L×1 and w i ∈ C K×1 , i = O, A, B are the range profile and ISAR image from the ith channel, respectively; Φ = [φ 1 φ 2 · · · φ K ] denotes the partial Fourier matrix, in which T is the kth basis, and I is the index sequence of SA. The likelihood of S is assumed to be complex Gaussian distributed with a noise variance of σ 2 , as where I denotes an identity matrix. To utilize the similarity between ISAR images from different channels, each row of w is assigned to be complex distributed with the same variance as where α i is the reciprocal of the variance or the precision. Then, the full prior of w is naturally obtained by combing these row priors as It is seen that the precision vector α controls the sparse degree of rows of w. When α i goes to infinity, the corresponding row of w is restricted to zero. Differently to the InISAR imaging algorithm based on SBL for single measurement vector (SMV) [20], which models the ISAR image from different channels with a different precision vector, the prior in Equation (12) utilizes a common precision vector to model all ISAR images, which enables them to be reconstructed with a higher matched degree.

SM-SBL
The traditional point estimation algorithms, such as maximum likelihood estimation (MLE) and maximum a posterior (MAP), regard the position of the mode of the likelihood or posterior as the estimation of unknown variables, and they do not require the deriving of a full posterior. In contrary, full Bayesian inference requires the posterior of unknown variables to be derived and estimates them with the mean of posterior, which performs more reliably than the point estimation algorithms, especially for the variable with multi-mode posterior. It has been confirmed in [11] that SBL generally performs better than basis pursuit (BP) and focal underdetermined system solver (FOCUSS) algorithms. Additionally, full Bayesian inference can obtain potentially useful high-order statistical information of unknown variables. In this subsection, the full Bayesian inference-based SM-SBL is derived.
We note that both the likelihood in Equation (10) and the prior in Equation (12) are complex Gaussian distributed. The posterior of w is also complex Gaussian distributed as follows [14]: with the expectation and covariance matrix given as where A = diag (α 1 , · · · , α K ). In SBL, the expectation µ is regarded as the estimation of w. However, the process of sparse reconstruction is not yet accomplished, because the hyperparameters in Σ, including α and σ −2 , are still to be learned from the observation. Several strategies have been proposed to deal with hyperparameter learning, such as the type II maximum likelihood algorithm [20], the expectation-maximization (EM) algorithm [21], the variational Bayesian method (VB) [21], and the Markov chain Monte Carlo (MCMC) method [22]. For Gaussian prior and likelihood, these algorithms produce similar results. Herein, we utilize the type II maximum likelihood algorithm to estimate α and σ 2 , for which the marginal likelihood is maximized with respect to α and σ 2 . The logarithmic likelihood L (α) is derived as follows [14]: where const denotes the term independent of α and σ 2 , and the covariance matrix, C, is given by Aiming to maximize L α, σ 2 , the precision vector α and noise variance σ 2 are updated via the EM algorithm [14] as where (·) (new) denote the updated hyperparameter. Then, the M-SBL [13] proceeds by iteratively updating Equations (14), (17) and (18) until convergence is reached and estimates w with the newly achieved mean µ. We note that the update of the covariance matrix in Equation (14) involves large matrix inversion. The M-SBL unavoidably suffers from low computational efficiency, and cannot meet the requirement of ISAR image reconstruction, which generally contains several hundreds of samples on both dimensions.
Inspired by [14], we derive SM-SBL to improve computational efficiency, for which the unknown variables are sequentially updated to avoid the time-consuming matrix inversion. To achieve this, the terms related to the ith precision, α i , are separated from the covariance matrix, C, in Equation (16) as where C −i denotes C with the contribution of the i-basis vector removed. To derive L (α), the determinant and inversion of C are derived by the matrix determinant and inverse identities as Substituting Equations (20) and (21) into Equation (15), we have where L (α −i ) and l (α i ) represent the terms relevant and irrelevant to α i , respectively. The fixed-point iteration method is utilized to estimate α i . To achieve this, the first derivative of L (α) with respect to α i is firstly derived as The fixed points of L (α) are then obtained by setting its first derivative to be zero, including α i = +∞ and To discuss the features of the two obtained fixed points, the second derivative of L (α) with respect to α i is further derived as It is seen that when α i = s 2 is zero. Then we have be satisfied because the precision α i is positive. When α i = +∞, all the derivatives of L (α) with respect to α i are zero. Additionally, the sign of the first derivative is decided by the term 3s i − 3 ∑ j=1 q i,j 2 , as given in Equation (27). Therefore, Therefore, the maximum points of L (α) can be summarized as These two equations provide a sequential estimation of the precision α i . To achieve it, a selection is firstly defined, which includes the basis φ i whose respective precision is not infinite (i.e., α i < +∞).
The number of bases in the selection is dynamic and is adjusted while updating α i . For example, if φ i is included in the selection and 3 ∑ j=1 q i,j 2 ≤ 3s i , then it is eliminated from the selection and α i = + ∞. If φ i is included in the selection and 3 ∑ j=1 q i,j 2 > 3s i , then it is kept and α i is updated with Equation (31).
If φ i is not in the selection and 3 ∑ j=1 q i,j 2 > 3s i , then it is put into the selection and α i is updates with Equation (31). During updating, s i and q i,j can be computed with Equations (25) and (26), respectively, or with the following equations, so as to improve computational efficiency.
Then we have If α k = +∞, then s k = S k and q k,j = Q k,j . The Woodbury matrix identity can be utilized to solve Equations (33) and (34) as follows, so as to further reduce computational burden.

SM-SBL-Based ISAR Imaging
In this subsection, the proposed SM-SBL is utilized to reconstruct ISAR images from different channels simultaneously. The procedure of multi-channel ISAR imaging based on SM-SBL is given as follows: (1) Initialize the noise variance by σ −2 = 0.1var [S], for which the variance of the combined range profile S is used because the noise variance is approximately similar to it. (2) Randomly choose a basis φ i as the initial selection; compute the respective precision α i by Equation (31), that is, and set α k = +∞ for k = i.
Update the expectation µ and the covariance matrix Σ with Equation (14); these are scalars when the selection includes only one basis. (4) Update s k and q k,j with Equations (35) and (36), respectively, for all the bases, including those that are within the selection and those that are not.
If θ i > 0 and φ i is included in the selection, then update α i .
If θ i > 0 and φ i is not in the selection, then put it in and update α i .
If θ i ≤ 0 and φ i is in the selection, then eliminate it from the selection and let α i = + ∞.
Update the noise variance by Equation (18). (10) Terminate the iteration if convergence is reached, or jump to step (3).
In step (2), the basis used to initialize the selection can be chosen by maximizing the projection of the observation to bases as φ i = arg max In step (5), a basis is chosen from the whole bases randomly, or, it can also be chosen as the basis that maximizes the likelihood, so as to improve the convergence rate. Noting that the bases in the selection are much less than those in the original dictionary, Φ, the proposed SM-SBL is much more computationally efficient than M-SBL. Additionally, the update of µ, Σ, s k and q k,j in steps (3) and (4) can be implemented in a sequential manner to further improve computational efficiency. They are sequentially updated in three different cases, which are related to steps (6), (7) and (8), respectively. In step (6), the selection remains unchanged, and only α i is updated. In this case, µ, Σ, s k and q k,j can be sequentially updated as In step (7), φ i is put into the selection, and α i is updated. In this case, the updates of µ, Σ, s k and q k,j are derived as where (8), φ i is eliminated from the selection, and α i is set as α i = + ∞. In this case, µ, Σ, s k and q k,j can be sequentially updated as The renewed expectation and covariance matrix can be achieved by deleting the ith row of µ (new) and the ith row and column of µ (new) .

Outlier Elimination and 3-D Rotational Rate Estimation Based on LM
The SM-SBL-based ISAR imaging method proposed in the former subsection can simultaneously produce three ISAR images from three different channels. The ISAR images from channels A and B are then conjugately multiplied by that from channel O to extract the interferometric phase ϕ i,p in Equation (6), which can be further utilized to estimate the 3-D coordinate of each scatterer by Equation (7). We note that the SA data decreases the quality of ISAR images and inevitably induces outliers during estimation of the 3-D coordinate of each scatterer. In this subsection, an outlier elimination method based on LM is proposed. Additionally, the 3-D rotational rate of the target is estimated during outlier elimination.
The Doppler of each scatterer on the ISAR image is given in Equation (6) as where ∆ f and P represent the Doppler error and number of scatterers, respectively. Equation (52) can be rewritten as Therefore, the LSM can be utilized to estimate the 3-D rotational rate of target aŝ Then, the x and z coordinates of the rotational rate of the target can be estimated asω x =ω (2) andω z =ω (1), respectively. Subsequently, the outliers can be detected by It is seen that a scatterer is regarded as an outlier when its Doppler error exceeds the preset threshold δ, which is suggested to be set as 0.1-0.5 Hz. The smaller the value it is set at, the more outliers that are detected. The motion estimation and outlier elimination in Equations (54) and (55) can be iterated to improve the estimate accuracy, and the threshold δ is kept decreased during the iteration. Generally, two to three iterations are enough to achieve satisfactory accuracy.
The flowchart of the proposed SM-SBL-based SA-InISAR imaging method is given in Figure 2. Firstly, the range compression and translational motion compensation are implemented to the three channels' SA radar echoes to obtain the three channels' range profiles. Then, the traditional range-Doppler (RD) imaging method is utilized to achieve image registration [18]. Subsequently, the three channels' ISAR images are jointly reconstructed by SM-SBL. The final 3-D image of the target can then be obtained through the processes of ISAR image interferometry, 3-D coordinate derivation, 3-D rotational rate estimation and outlier elimination, successively.

Experiments
This section utilizes a set of simulated data of an airplane to validate the performance of the proposed SM-SBL-based InISAR imaging method. As shown in Figure 3, the simulated airplane contains 113 scatterers. The InISAR system is composed of three channels, which are located as shown in Figure 1. The baselines L 1 and L 2 are set as 2 m. The initial range and rotational speed of the target are set as 20 km and (0.01, 0, 0.02) rad/s, respectively. The radar is assumed to transmit a linear frequency modulation (LFM) signal whose central frequency, bandwidth, pulse width and pulse repetition frequency (PRF) are 9 GHz, 600 MHz, 100 µs and 100 Hz, respectively. The test data contained 256 pulses, and each pulse contained 256 samples. Complex Gaussian noise was added to each pulse to simulate the noise environment, which made the SNR of the radar echo fall to 5 dB.  We extracted 41 pulses from the complete 256 pulses to simulate two typical types of SA data, that is, random missing sampling (RMS) and gap missing sampling (GMS) SA data [3]. Figure 4a,b give the range profile sequence for the RMS and GMS, respectively. The background noise is clearly observed because of a low SNR.   Four ISAR imaging methods, including the RD method, the l 1 norm-based InISAR imaging algorithm [7], the SBL-based InISAR imaging algorithm [8] and the proposed SM-SBL-based method, are utilized to obtain ISAR images from the RMS and GMS. The imagery results for the RMS are given in Figure 5, in which the first, second, third and fourth columns are the ISAR images obtained by RD, l 1 -norm, SBL and SM-SBL, respectively, and the first, second and third rows give the ISAR images from the channels O, A and B, respectively. It is seen that the results of RD are defocused in this case, which indicate its invalidity for SA data. In contrast, the other three methods obtain well-focused ISAR images, and the images obtained by SM-SBL are relatively less noisy than those obtained by l 1 -norm and SBL, which validates its superior performance in terms of noise suppression against l 1 -norm and SBL. Noting that InISAR imaging requires perfectly matched multi-channel ISAR images for interferometry, the correlation coefficients (CC) of the multi-channel ISAR images were computed to numerically compare the matched degree of ISAR images obtained by different methods.
The results are given in Table 1. For purposes of comparison, the CC of the images obtained by RD for the complete data are also given. It is seen that the proposed SM-SBL achieved the largest CC, which validates its superior performance in terms of improving the matched degree of multi-channel ISAR images. Similarly, the ISAR images for GMS obtained by the four methods are given in Figure 6, for which each row and column show the ISAR images from the same channel and that are achieved by the same method, respectively. It is seen that RD suffers from high side lobes and strong noise because the pulses for imaging are limited. Again, the other three SSR algorithms successfully obtain high-resolution ISAR images, as shown in the second, third and fourth columns of Figure 6. Compared with l 1 -norm and SBL, the proposed SM-SBL achieves images with a lighter noise floor, which indicates its superior noise suppression performance. Additionally, the CC between ISAR images from different channels achieved by the four methods along with RD for the complete data are given in Table 2, so as to further compare their performance in terms of improving the matched degree of multi-channel ISAR images. It is seen that the proposed SM-SBL obtains the highest CC again, which validates the effectiveness of jointly reconstructing multi-channel ISAR images with SM-SBL. Additionally, it should be noted that RD achieves higher CC than l 1 -norm and SBL, which indicates that RD performs better than SSR for producing matched multi-channel ISAR images. In contrast, although the proposed SM-SBL belongs to SSR, it is still able to obtain higher CC than RD, which further validates its superior performance in terms of improving the matched degree of multi-channel ISAR images.  Next, we compare the performances of RD for complete and sparse data, l 1 -norm, SBL and the proposed SM-SBL under different noise conditions. The parameters of this experiment were the same as for the former, except that the SNR ranged from −5 to 10 dB in steps of 1 dB. The curves of CC of multi-channel ISAR images with respect to the SNR obtained by different algorithms are given in Figure 7 to quantitatively compare their performance. It is seen that the proposed SM-SBL obtains the highest CC under any noise levels for both RMS and GMS. For example, the CC obtained by SM-SBL approaches 0.8 for a SNR of −5 dB, while those obtained by RD and SBL are lower than 0.3. Therefore, the proposed SM-SBL performs more robustly than RD, l 1 -norm and SBL in terms of producing well-matched multi-channel ISAR images under low-SNR conditions. We next compare the performance of the five methods under different SA degrees. The experimental parameters were the same as for the first experiment, except that the SA degree ranged from 0.1 to 0.9 in steps of 0.1, and the SA degree was defined as the rate of pulses extracted from the complete data. Figure 8 gives the curves of CC with respect to the SA degree obtained by different algorithms. Additionally, the proposed SM-SBL obtains the highest CC under any SA degrees for both RMS and GMS, which further validates its superior performance.  We now compare the 3-D scatterer distributions reconstructed with the multi-channel ISAR images obtained by l 1 -norm, SBL and SM-SBL. The reconstructed 3-D scatterer distributions before and after outlier elimination for RMS are given in Figure 9, for which the first, second and third rows show the results obtained by l 1 -norm SBL and SM-SBL, and the first and second columns give the results before and after outlier elimination, respectively. It is seen that the outliers are successfully eliminated in the 3-D scatterer distributions after outlier elimination, which confirms the effectiveness of the proposed outlier elimination method. Additionally, it is clear that the proposed SM-SBL obtains more realistic 3-D scatterer distributions of airplanes than l 1 -norm and SBL, because it brings improvement in terms of the matched degree of multi-channel ISAR images. Figure 10 gives the 3-D scatterer distributions for GMS. Again, the proposed SM-SBL obtains better results than SBL, which further validates its superior performance. Additionally, it is seen that the 3-D distribution obtained by SM-SBL before outlier elimination is similar to that after outlier elimination, because SM-SBL already achieves satisfactory 3-D distribution with few outliers, and the process of outlier elimination can be omitted. In contrast, the 3-D distributions obtained by l 1 -norm and SBL contain considerable numbers of outliers, which makes outlier elimination a necessity.   Finally, the mean-square errors (MSE) of estimated rotational velocities and computational time of l 1 -norm, SBL and SM-SBL are given in Table 3 to further quantitatively compare their performance. It is seen that the proposed SM-SBL estimates the rotational velocity more accurately than l 1 -norm and SBL. It is also far more computationally efficient because of the sequential update strategy that it utilizes.

Conclusions
This paper proposes a novel sparse InISAR imaging algorithm. The SM-SBL method is proposed to jointly reconstruct well-matched multi-channel ISAR images with highly computational efficiency. The LSM is utilized to eliminate the outliers and estimate the rotational velocity of the target. On the basis of numerous experimental results, we are confident to draw the conclusion that the proposed sparse InISAR imaging algorithm has two compelling advantages. Firstly, it improves the matched degree of ISAR images from different channels for SA data. For example, when the SNR is −5 dB, the correlation coefficient of multi-channel ISAR images for GMS obtained by the proposed algorithm approaches 0.8, while that obtained by the SBL-based InISAR imaging algorithm is merely 0.4. Additionally, the proposed SA-InISAR imaging algorithm is computationally efficient. Benefitting from the sequential procedure of SM-SBL, the SM-SBL-based SA-InISAR imaging algorithm avoids the time-consuming large matrix inversion. Experimental results indicate the proposed algorithm is 6 times faster than the SBL-based InISAR imaging algorithm.