Improved FOCUSS method for reconstruction of cluster ... - Springer Link
SCIENCE CHINA Information Sciences
. RESEARCH PAPER .
August 2012 Vol. 55 No. 8: 1776–1788 doi: 10.1007/s11432-012-4628-1
Special Issue
Improved FOCUSS method for reconstruction of cluster structured sparse signals in radar imaging HU ChenXi1,2 , LIU YiMin1 ∗ , LI Gang1 & WANG XiQin1 1Intelligent
Sensing Laboratory, Department of Electronic Engineering, Tsinghua University, Beijing 100084, China; 2Xi’an Electronic Engineering Research Institute, Xi’an 710100, China Received April 1, 2012; accepted May 5, 2012
Abstract The high resolution imaging of distributed targets in millimeter-wave radar system is studied in this paper. We use Gaussian Mrakov random field (GMRF) to represent the clustering property of the targets. Our novel method, called Clustering FOCUSS, incorporate an additional cluster constraint into the process of the focal underdetermined system solver (FOCUSS) algorithm. Simulation results indicate that the novel algorithm has a higher imaging accuracy than the methods of Capon beamforming, the l1 norm algorithm and the FOCUSS algorithm. Keywords
antenna array, high resolution imaging, cluster sparsity, GMRF, Clustering FOCUSS
Citation Hu C X, Liu Y M, Li G, et al. Improved FOCUSS method for reconstruction of cluster structured sparse signals in radar imaging. Sci China Inf Sci, 2012, 55: 1776–1788, doi: 10.1007/s11432-012-4628-1
1
Introduction
In automotive radar system, millimeter-wave radar is an indispensable forward-looking sensor to detect and estimate the relative angle, range and velocity of the targets in a multiple target scenario. In order to satisfy system requirements, millimeter-wave radar with multi-channel antenna array is proposed in [1]. This system is composed of a single transmitter to radiate a linear FM signal and M receiving channels to acquire echo signals to estimate the angle, range and velocity of the targets. Application of images from automotive radar depends on its image quality. Resolution is an important parameter to verify this. When the resolution is equal to the target size, the target is only a bright spot in the imaging result; when the resolution is much lower than the target size, the detailed information of the target could be obtained from the imaging result. There are two ways to enhance resolution of automotive radar. One way is to refresh and upgrade the hardware (such as increasing bandwidth and antenna aperture [2,3]). However, the method results in an increase in radar size and complexity of configuration, which is not suitable for automotive application. The other way is to use data processing algorithm to obtain superresolution. Various methods have been proposed to obtain high accuracy and resolution in radar imaging. Capon beamforming [4] and other covariance matrix based estimation methods perform well when received ∗ Corresponding
author (email: [email protected])
c Science China Press and Springer-Verlag Berlin Heidelberg 2012
info.scichina.com
www.springerlink.com
Hu C X, et al.
Sci China Inf Sci
August 2012 Vol. 55 No. 8
1777
signal from different targets are incoherent, but when the range and velocity of targets result in coherent signals, these covariance matrix based estimation methods will degrade dramatically. Considering that there are only a few targets in the scene discussed here, which means that the number of nonzero pixels in imaging result is small, sparse reconstruction algorithms are used for imaging targets in linear frequencymodulated continuous-wave (FMCW) automotive radar [5]. This method does not experience coherent signals problem and significantly improves the imaging accuracy. However, this method doesn’t exploit any structure information (except sparsity) that may exist in the sensed signals. Compared to the scale of the scene, the targets are usually not only sparse but also clustered in a distinct way. A few attempts have been made to utilize these cluster priors to get better sparse recovery. Refs. [6,7] argued that how to recover cluster sparse signal could be regarded as a mixed l2 /l1 -optimization problem. Ref. [8] extended MP and OMP to cluster sparse signal and got BMP and BOMP. Ref. [9] presented an extention of sparse dictionary learning to obtain dictionary elements belonging to groups of variables forming convex subsets. For simplicity, they all assume that the cluster structures (such as the cluster number/size/location) are known before recovery. Ref. [10] used overlapping neighborhoods to modify the sparsity term to better recover cluster sparse signal. However, this would lead the coefficients within each neighborhood to be smooth and cause the loss of details. In this paper, using the GMRF model [11,12] to describe the clustering of signal support, we propose a modification to standard FOCUSS [13,14] to make it more efficient for cluster sparse solutions. The novel method could work without knowing the details of cluster structures, which are not always available in automotive application. Simulation shows that imaging results get improved. We have already presented the Clustering FOCUSS method in [15], where we used it to improve the radar imaging of distributed targets with multi-channel FMCW radar. In this paper, we extend the application of the Clustering FOCUSS method to pulse radar, analyze the selection of the MRF model, and address some interesting issues.
2
Signal model from view point of sparse recovery
In this paper, for automotive radars we use antenna array to interrogate a stationary scene with pulses or continuous wave to reveal the spatial distribution of the scatterers. From viewpoint of compressive sensing (CS), we obtain the linear system of equations as Eq. (1) by discretizing the scene reflectivity function over range and azimuth on a grid of points {θi , Ri }, where q is a rearrangement vector of the sampled received baseband signal, α is the reflectivity vector, the ith column of Φ is the contribution of the target with azimuth θi and range Ri , and n is the noise. The model Eq. (1) can be easily adapted for use with the CS framework. q = Φα + n. (1) Then we model the reflected signal of FMCW radar with antenna array to obtain the linear equation in Eq. (1), which has been widely used in automotive application. We consider a uniform linear array (ULA) of MR antennas which share a single transmitter, as shown in Figure 1. The transmitted signal takes the following form: x0 (t) = exp{j2π[f0 t + (B/T )(t2 /2)]},
0 t < T,
(2)
where T is the pulse period, f0 is the carrier frequency and B is sweeping bandwidth. Assuming that there is a single target in the scenario, the echo signal received by the mth antenna is R − vt d sin θ xm (t) = ax0 2 − (m − 1) , m = 1, 2, . . . , M, (3) c c where a is the received signal’s amplitude, c is the velocity of light, d is the interelement ULA spacing, R is the initial radial range to the target at time t = 0 for the first transmitted pulse, v is the radial velocity of the target, and θ is the azimuth angle of the target.
1778
Hu C X, et al.
Figure 1
Sci China Inf Sci
August 2012 Vol. 55 No. 8
Diagram of multi-channel antenna array radar.
Mixing the echo signal xm (t) with the complex conjugate transmitted signal and sampling, the sampled data can be arranged in a “Data Matrix”, which is (considering B f0 , v c) Qm,n = a exp{−j2π[(m − 1)fθ + (n − 1)fR,v + ϕR ]},
m = 1, 2, . . . , M,
n = 1, 2, . . . , N,
(4)
where fθ = f0 (d sin θ/c) is referred to as angle frequency, fR,v = [(B/T )(2R/c) − f0 (2v/c)]Ts is range (velocity) frequency, ϕR = f0 (2R/c) − (B/2T )(2R/c)2 , Ts = 1/fs , fs is system sample frequency. It is (k) important to note that range R and velocity v are coupled in fR,v . This is range-velocity coupling problem of FMCW radar and some efficient methods to solve this problem are given in [16], so we estimate angle and range without considering the velocity of target in this paper. Discretizing angle and range axis, the angle-range domain can be divided into square grids [5]. Here we assume that there is a target in each grid which represents the echo signal with angle frequency and range frequency corresponding to that grid. Then the received data can be obtained by summarizing the contribution of all the targets with different angle and range frequency. We rearrange the data matrix [Qm,n ] to a vector, q = [Q1,1 , Q2,1 , . . . , QM,1 , . . . , QM,N ]T and Eq. (4) can be rewritten as follows: q=
M0 N0
αm0 ,n0 Φm0 ,n0 ,
(5)
m0 =1 n0 =1
where (m0 )
αm0 ,n0 = α(fθ
(m0 )
Φm0 ,n0 = Φ(fθ =
(n0 )
, fR
) = a(m0 ,n0 ) exp(−j2πϕR 0 ),
(n0 )
, fR
(n )
) = ΦR ⊗ Φθ
(6)
(m ) (m ) [1, exp(−j2πfθ 0 ), . . . , exp(−j2π(M − 1)fθ 0 ), (n ) (m ) . . . , exp(−j2π((N − 1)fR 0 + (M − 1)fθ 0 ))]T ,
M0 and N0 are the number of quantization of angle and range axis, respectively. In consideration of noise, model Eq. (5) can be written as matrix form as Eq. (1). Note that the row dimension M × N of matrix Φ is much less than column dimension M0 × N0 and most of elements in solution vector α are negligible and only a small portion of it, which represents contribution of targets, is remarkable. Therefore, this problem could be solved by sparse recovery technique. We can also obtain the CS model Eq. (1) of pulse radar by analyzing the process of pulse compression and beam forming. Different from the FMCW radar, each column of the matrix Φ represents the set of dechirped samples across the complete antenna array for a given point in the discretized scene. The detailed process could be found in [17,18]. According to CS theory, when the matrix Φ has the restricted isometry property (RIP), it is possible to recover the K large αi s from fewer measurements q [19]. If Φ is of large size, estimating and testing
Hu C X, et al.
Sci China Inf Sci
August 2012 Vol. 55 No. 8
1779
the RIP would be impractical. A tractable yet conservative bound on the RIP can be obtained through the mutual coherence of the columns of Φ [17]. Since the echo signal of one target has low coherence except with itself, the incoherence property of matrix Φ is verified. Furthermore, we can use randomly sampled data, which means that the measurement matrix is Fourier randomly measured ensembles, to recover sparse signals [20,21]. For antenna array radar, randomization of array element locations may function as a means of reducing the correlation among responses from different targets [22].
3
Image algorithm based on sparse recovery
3.1
The FOCUSS algorithm
According to the theory of sparse representation, when the actual distribution is sparse in a domain, the ill-posed problem in Eq. (1) can be efficiently solved using the technique of sparse representation. The basic form of sparse representation is defined as α = arg min α0 subject to q − Φα2 ε,
(7)
where · p stands for the lp norm and ε is the error allowance. However, this problem is a combinatorial problem and NP-hard. To address this difficulty, a number of practical algorithms have been proposed to approximate this sparse solution. Most algorithms for sparse representations can be classified into three categories: greedy pursuit algorithms, lp norm regularization based algorithms, and iterative shrinkage algorithms [23]. Greedy pursuit algorithms may be trapped into local minima and cannot be corrected back due to the greedy nature of algorithms. The lp norm regularization based algorithms demand a high computational effort especially for large-scale problems. So an iterative shrinkage method named FOCUSS [13,14] is discussed here. The FOCUSS method is based on the weighted least squares estimate. It is motivated by the observation that if a sparse solution is desired then choosing a solution based on the smallest l2 norm is not appropriate. The basic form of the recursive method is composed of the following steps. W (k) = diag(α(k−1) ),
(8)
z (k) = (ΦW (k) )† q,
(9)
α
(k)
=W
z
(k) (k)
,
(10)
where diag(·) is the operation of diagonalization, A† = (AH A)−1 AH denotes the pseudoinverse operation of matrix A, W (k) and α(k) denote the weighting matrix and the corresponding sparse solution at the kth iteration. In each iteration, the FOCUSS method gradually reinforces some of the already prominent entries in α(k) while suppressing the rest until they become close to zeros. Ref. [13] shows that FOCUSS serves as an iterative approximation of the sparse representation using lp , 0 p 1 norm minimization. By starting with the proper initial value, FOCUSS can converge to the global optimal value and provide a sparse representation. 3.2 3.2.1
The Clustering FOCUSS algorithm Cluster sparsity constraints
In literature about transform compression, many natural and man-made signals and images are considered as sparse or compressible with their sparse supports (set of nonzero coefficients) often having an underlying order. This paper is an attempt to exploit a priori information on coefficient structure of the cluster sparse signal to enhance the reconstruction accuracy and promote the ability of recovery. Figure 2 illustrates several examples of structured sparse support in the radar imaging application. Compared to the scale of the scene, the targets are usually not only sparse but also clustered in a distinct way. Nevertheless, this clustering property is not fully exploited by current CS recovery algorithms.
1780
Figure 2
Hu C X, et al.
Sci China Inf Sci
August 2012 Vol. 55 No. 8
(a) The location of the scattering centers at 24 GHz and 77 GHz of the vehicle rear end with a bandwidth
of 2 GHz [24]; (b)the backhoe CAD model and high resolution image of the backhoe for wide-angle imaging aperture of 110◦ [25].
Figure 3
A cluster sparse signal and its signal support. (a) α: a cluster sparse signal with length M = 40 and sparsity
D = 7; (b) f : the support of the signal α; (c) f : the support of a randomly sparse signal with length M = 40 and sparsity D = 7.
We use f ∈ {0, 1}M to represent the support of a cluster sparse signal α ∈ RM , where fi = 0 when αi = 0 and fi = 1 when αi = 0. For example, a cluster sparse signal α and its support f are shown in Figures 3 (a) and (b). This signal has sparsity D = 7, and its nonzero coefficients appear in only one cluster. The support of a randomly sparse signal with sparsity D = 7 is also shown in Figure 3 (c) (Here “randomly sparse” means that the signal is only sparse without any other additional structure constraint). By comparing f in Figures 3 (b) and (c), we find that the cluster sparse signal has the following properties: if the neighbors of an element of the signal are nonzero, the element itself is more likely to be nonzero. Clustering of the nonzero coefficients in a sparse signal representation can be realistically captured by a probabilistic graphical model such as a Markov random field. We use an MRF (S, N ) [12] to model the probability density function (PDF) of the cluster sparse signal support, where S denotes a set of M vertices (one for each of the support indices) and N denotes a neighborhood by which the sites in S are related to one another. N is defined as N = {Ni |∀i ∈ S} where Ni is the set of sites neighboring i. The set {0, 1}M from which f is taken is called the label set of the MRF.
Hu C X, et al.
Figure 4
Sci China Inf Sci
August 2012 Vol. 55 No. 8
1781
Elements of z at each iteration for the 4×10 matrix Φ and the vector α with only one nonzero element α4 .
Note that only one element z4 of z converges to one, whereas the rest become zero.
In MRFs, only neighboring vertices have direct interactions with each other, which means the following condition is satisfied: P (fi |fS−{i} ) = P (fi |fNi ), (11) where S − {i} is the set difference, fS−{i} denotes the set of labels at the sites in S − {i}, and fNi = {fi |i ∈ Ni } stands for the set of labels at the sites neighboring i. The MRF models can be further classified according to assumptions made about individual fi . 3.2.2
The procedures of Clustering FOCUSS algorithm
In FOCUSS, while the entries of α(k) converge to zero and nonzero values, the corresponding entries (k) (k) (k) (k) in z (k) converge to 0 or 1, i.e., zi → 0 as αi → 0, and zi → 1 as αi approach nonzero values (Figure 4). Therefore, we suppose that the signal support f could be approximated by z. Here we select the GMRF model [11,12] from the class of MRF models to describe the characteristic of z. The label set from which the value of zi is taken is the real line and the joint distribution is multivariate normal. Its conditional PDF is 2 1 1 exp − 2 zi − βi,j zj p(zi |zS−{i} ) = p(zi |zNi ) = , (12) 2σ (2πσ 2 )1/2 j∈Ni
where σ 2 is the variance of a Gaussian noise, and βi,j are constants reflecting the pair-site interacton between i and j: βj,i = βi,j , βj,i = 0 if and only if sites i and j are neighbours. As shown in Figure 5, the conditional PDF of zi is related to the values of its neighbors. For example, zi in Figure 5 (a) whose neighbors are both nonzero is more likely to be nonzero than in Figure 5(c) whose neighbors are both zero. The GMRF model is widely used in many scientific fields to represent interaction among spatial data with a Gaussian noise [26,27]. There are two main reasons for our choosing the GMRF model. Firstly, the continuous label set of the GMRF model is more proper than those models with discrete label set since z converges to f during iterations. Secondly, the form of conditional p.d.f. is normal distribution, which makes our novel algorithm have simple form and low computational cost. Thus we have the following iterative relaxation algorithm: 2 1 1 (k) (k−1) (k−1) (k−1) )= exp − 2 zi − βi,j zj . (13) p(zi |z 2σ (2πσ 2 )1/2 j∈Ni
1782
Hu C X, et al.
Figure 5
Sci China Inf Sci
August 2012 Vol. 55 No. 8
The pdf of zi with different neighbor conditions.
(k)
(k−1)
To simplify calculations, zj in the right-hand side of Eq. (13) is replaced by zj . Since z (k) cannot express the signal support of f properly at first, the cluster constraint should be added after the first steps. In FOCUSS, Eq. (9) is the smallest l2 norm solution to the problem ΦW z = q, and now it can be regarded as the MAP estimate, where the likelihood function p(q|z) is expressed in terms of the joint Gaussian distribution. 1 1 H −1 tr{(q − ΦW z) exp − Σ (q − ΦW z)} . (14) p(q|z) = 0 2 ((2π)M det{Σ0 })1/2 Here we assume that n in Eq. (5) is independent and identically distributed (i.i.d.) random variables from a normal distribution with a constant variance σ02 , such that σ02 IM ∈ RM×M . Applying the Bayesian framework to the prior Eq. (13), the posterior distribution has the form p(z|q, z (k−1) ) =
p(q|z)p(z|z (k−1) ) . p(q|z (k−1) )
(15)
The MAP estimate in the kth iterative update can be obtained by solving the following regularized least squares minimization problem: z (k) = arg min Ψ (z|q, z (k−1) , λ),
(16)
Ψ (z|q, z (k−1) , λ) = −2 ln p(z|q, z (k−1) , λ) = q − ΦW z (k) 22 + γz (k) − Bz (k−1) 22 ,
(17)
z
where
γ = σ02 /σ 2 , and B is the square matrix whose diagonal elements are unity and whose off-diagonal (i, j) element is βi,j . The stationary point of Eq. (16) expressed in the vectorized form is as follows: z (k) = z∗ = ((ΦW )H (ΦW ) + γI)−1 ((ΦW )H q + γBz (k−1) ).
(18)
The procedures of the Clustering FOCUSS algorithm proceed as follows. 1) Obtain the initial value as a low-resolution estimation computed as (0)
αi
= ΦH i q,
1 i M 0 N0 ,
(19)
Hu C X, et al.
Sci China Inf Sci
August 2012 Vol. 55 No. 8
1783
W (0) = diag(abs(α(0) )),
(20)
Γ = {i, 1 i M0 N0 }.
(21)
2) Update the weighting matrix and the estimation result accordingly as z (k) |Γ = ((Φ|Γ W (k) )H (Φ|Γ W (k) ) + γI)−1 ((Φ|Γ W (k) )H q + γ(Bz (k−1) )|Γ ),
(22)
α
(23)
(k)
|Γ = W
z
(k) (k)
|Γ ,
where α(k) |Γ stands for the Γ subset of the vector α(k) . 3) Update the adaptive subspace using the current solution as (k)
Γ = arg(|αi | Th), 1 i M0 N0 ,
(24)
(k)
where αi is the ith element of the vector α(k) , Th stands for the threshold and the set Γ is the adaptive subspace at the kth iteration.
4 4.1
Simulation results Simulation examples and analysis of Clustering FOCUSS (k)
In the second step of Clustering FOCUSS, zi in Eq. (18) is composed of two terms. The first term is the same as that in FOCUSS, and the second term is the cluster sparsity constraints. Note that if (k−1) (k) have nonzero values, the second term of zi would be greater than zero. In other neighbors of zi words, if zi has nonzero neighbors at the (k − 1)th iteration, it would be larger in Clustering FOCUSS than in FOCUSS at the kth iteration. From this we can reach an interesting conclusion that Clustering FOCUSS is more likely to reinforce or at least suppress sites with nonzero neighbors more slowly, so that the nonzero coefficients would prefer to be clustered. (k−1) , the sites in a cluster which contains targets will reinforce If there has already been clusters in zi each other through cluster sparsity constraints and the sites in a cluster which contains no target will be suppressed by the weighted l2 norm constraints. Simulation is carried out to demonstrate the efficiency of the proposed algorithm. To test the analysis above, we assume that the cluster structured sparse signal α is randomly generated with length L, sparsity D and clusters K. Its nonzero entries are drawn from a uniform distribution. The columns of matrix Φ are drawn form a standardized Gaussian distribution and then normalized to unit l2 norm. Here two pixels are defined as neighbors only when the distance between them is less than the minimum resolution and βi,j = 1/S with S is the number of neighbors of αi . The noisy data was generated by adding to the exact mixtures a zero-mean Gaussian noise with σ02 adjusted to a desired value of SNR0 (SNR0 = 20 log10 (qexact /n)). The signal shown in Figure 6 is generated with L = 250, D = 40, K = 3. The values of z computed by Clustering FOCUSS are shown along with the iterations in Figure 7. As is evident in the graph, the Clustering FOCUSS algorithm could cluster the sparse coefficients of signal recovery, reinforce clusters containing real targets and suppress clusters with no target. The recoveries of FOCUSS and Clustering FOCUSS are both shown in Figure 7. In the reconstruction of Clustering FOCUSS, the noise is suppressed well and the weak cluster of the signal is reconstructed correctly. 4.2
Statistics of Clustering FOCUSS
In order to verify the performance of Clustering FOCUSS, we take into account the number of measurements P and SNR0 which are able to reach the successful reconstruction. The reduction of measurements P means that we can get high resolution of the scene with sub-Nyquist measurement rates and less antennas. SNR0 reflects the robustness of the performance of Clustering FOCUSS to noise. They are both impressive in the theory of sparse recovery.
Hu C X, et al.
1784
Figure 6
Sci China Inf Sci
August 2012 Vol. 55 No. 8
(a) Reconstruction with FOCUSS; (b) reconstruction with Clustering FOCUSS; (c) original signal.
Figure 7
Elements of z in Clustering FOCUSS at the (a) 1st; (b) 5th; (c) 10th; (d) 20th iteration.
Sparse signals with 2 clusters are considered. For fixed cluster structured sparse signal(L = 250, D = 30, K = 2), the number of measurements P is ranging from 60 to 140 with step 10, and SNR0 is ranging from 10 dB to 26 dB with step 2 dB. For each step, we run the program 100 times independently with different sensing matrix. The accuracy of reconstruction is evaluated through relative root mean ˆ Na is the number of square error (RRMSE) between original sparse signal α and its reconstruction α. independent repeated experiments, and α ˆi is the reconstruction of the ith repeated experiment.
RRMSE = 20 log10
Na 1 α − α ˆi . Na i=1 α
(25)
The plots in Figure 8 show the dependence of RRMSE on the parameters: P and SNR0 . Figure 9(a) shows RRMSE plot versus SNR0 at P = 100; Figure 9(b) shows RRMSE plot versus P at SNR0 =15 dB. It is easy to notice that Clustering FOCUSS can give the estimates that have lower RRMSE values than
Hu C X, et al.
Figure 8
Sci China Inf Sci
August 2012 Vol. 55 No. 8
1785
RRMSE plots versus the number of measurements P and SNR0 . (a) OMP; (b) FOCUSS; (c) Clustering
FOCUSS.
Figure 9
(a) RRMSE plots versus SNR0 at the number of measurements P = 100; (b) RRMSE plots versus the number
of measurements P at SNR0 =15 dB.
OMP (which is a greedy algorithm for compressive sensing recovery [28]) and FOCUSS, even with less measurements P and lower SNR0 . We plot the RRMSE as a function of the cluster-sparsity level of the signal to be recovered in Figure 10 at the number of the measurement P = 100, SNR0 =15 dB. Here the clusters are of the fixed size 10, and the locations of the nonzero clusters are chosen uniformly among all possible locations. We can see that Cluster FOCUSS outperforms FOCUSS and OMP especially when the cluster-sparsity level increases.
5
Discussion on the performance of Clustering FOCUSS in real scene
Assume a ULA of MR = 9 antennas with a half-wavelength interelement spacing and a linear FM signal with carrier frequency f0 =24 GHz, the sweeping bandwidth B=100 MHz, and the pulse period T =
Hu C X, et al.
1786
Figure 10
Figure 11
Sci China Inf Sci
August 2012 Vol. 55 No. 8
RRMSE plots versus cluster-sparsity level at measurement number P =100 and SNR0 =15 dB.
The imaging results by different methods with SNR0 =10 dB. (a) Origin targets; (b) conventional beamforming;
(c) Capon beamforming; (d) L1 norm algorithm; (e) FOCUSS algorithm; (f) Clustering FOCUSS algorithm.
0.75 ms. The sample frequency of the system is fs =60 kHz, which satisfies Nyquist sample rate. SNR0 is 10 dB. We consider the imaging result of a backhoe with wide-angle SAR imaging [25] as the real scene. The detailed information of the targets is shown in Figure 11(a), where x-axis is range plane and yaxis is angle plane. The scene is divided into square grids whose dimensions are 1/4 of the range and azimuth resolution of the multi-channel FMCW radar. Assume that there is a target in each grid and the reflectivity of the target equals the amplitude of the pixel in the same position of SAR image shown in Figure 11(a). Then the received data is obtained by summarizing the contribution of all the targets in the scene according to the model in Eq. (5). Figure 11(b) is the imaging result using conventional beamforming. Since the resolution is limited by antenna aperture and sweep bandwidth, the boundaries of the target are extended and the details in the target cannot be recovered correctly. Figure 11(c) is the imaging result using Capon beamforming. By Eq. (5), two close-set pixels will result in coherent signals which unfortunately happen in our problem.
Hu C X, et al.
Sci China Inf Sci
August 2012 Vol. 55 No. 8
1787
In this case, Capon beamforming cannot distinguish each pixel well. Figure 11(d) is the imaging result using the l1 norm algorithm. Disjunctions appear instead of continuous distribution which is the truth. The “superresolution” phenomenon is caused by coherent columns of the matrix Φ. Figure 11(e) is the imaging result using the FOCUSS algorithm. The position of target is reconstructed correctly; however, in the reconstructed solution, weak parts of the origin targets almost disappear and strong parts are smaller. The reason is that the noise affects the elements of solution differently due to the amplitude of the pixels, which may amplify strong pixels in the solution and significantly distorts the rest. The imaging result of Clustering FOCUSS algorithm is shown in Figure 11(f). By utilizing the clustering property of the nonzero elements of α, the effect of the noise we talk about above is mitigated. It is apparent that the noise is suppressed well and the weak cluster of the signal is reconstructed correctly. In the reconstruction of Clustering FOCUSS, note that since the whole image is more likely to be continuous and clustered because of the cluster sparsity constraints, the isolated points of the origin target are also clustered, which is not the truth. The reason for this phenomenon is that the GMRF model we use here could not deal well with the discontinuous region (such as edge) of the scene. To solve this problem, one way is to incorporate Clustering FOCUSS with the ability to detect the edge and modify the reweighting rule there. Another way is to use a more accurate model such as the hierarchical two-level Gibbs model [29] which has a higher-level Gibbs distribution to characterize the cluster sparsity of the signal and a lower-level Gibbs distribution to describe the filling-in in each cluster.
6
Conclusion
A novel method based on block-sparse model has been proposed for imaging distributed targets with millimeter-wave radar with multi-channel antenna array. Utilizing the cluster sparseness in angle-range domain, our method named Clustering FOCUSS regards the weighted l2 norm minimization in the second step of the FOCUSS algorithm as the MAP estimator, and incorporates the GMRF graphical model to constrain the nonzero values of the solution. Compared with FOCUSS, Clustering FOCUSS can keep the recursive process from the error accumulation between iterations and prevent an over-focal solution; the distributed targets which have clustered signal support can also stand out even with lower SNR. The Clustering FOCUSS algorithm can also be used in other radar scenes which have cluster sparse property, such as SAR and ISAR. Much work needs to be done, however. We are working to analyze the global convergence of the Clustering FOCUSS algorithm. We also assert that other probabilistic signal models could be used to improve the existing signal recovery methods.
Acknowledgements This work was supported in part by National Natural Science Foundation of China (Grant No. 40901157), National Basic Research Program of China (973 Program) (Grant No. 2010CB731901), and Tsinghua National Laboratory for Information Science and Technology Cross-discipline Foundation.
References 1 Lee M S, Katkovnik V, Kim Y H. System modeling and signal processing for a switch antenna array radar. IEEE Trans Signal Process, 2004, 52: 1513–1523 2 Johnson D H, Dudgeon D E. Array Signal Processing: Concepts and Techniques. Englewood Cliffs: Prentice-Hall, 1993 3 Van Trees H L. Detection, Estimation and Modulation Theory, Part IV: Optimum Array Processing. New York: John Wiley & Sons, Inc, 2002 4 Capon J. High-resolution frequency-wavenumber spectrum analysis. Proc IEEE, 1969, 57: 1408–1418 5 Wang P, Meng H D, Wei Y M. FMCW radar imaging with multi-channel antenna array via sparse recovery technique. In: 2010 International Conference on Electrical and Control Engineering (ICECE), Wuhan, 2010. 1018–1021
1788
Hu C X, et al.
Sci China Inf Sci
August 2012 Vol. 55 No. 8
6 Stojnic M, Parvaresh F, Hassibi B. On the reconstruction of block-sparse signals with an optimal number of measurements. IEEE Trans Signal Process, 2009, 57: 3075–3085 7 Baraniuk R G, Cevher V, Duarte M F, et al. Model-based compressive sensing. IEEE Trans Inf Theory, 2010, 56: 1982–2001 8 Eldar Y C, Mishali M. Robust recovery of signals from a structured union of subspaces. IEEE Trans Inf Theory, 2009, 55: 5302–5316 9 Jenatton R, Obozinski G, Bach F. Structured sparse principal component analysis. In: Proceedings of the 13th International Conference on Artificial Intelligence and Statistics (AISTATS), Sardinia: 2010 10 Kavukcuoglu K, Ranzato M A, Fergus R, et al. Learning invariant features through topographic filter maps. IEEE Conference on Computer Vision Pattern Recognition, Miami, 2009. 1605–1612 11 Besag J. Spatial interaction and the statistical analysis of lattice systems (with discussion). J Royal Statist Soc Ser B, 1974, 36: 192–326 12 Li S Z. Markov Random Field Modeling in Image Analysis. 2nd ed. New York: Springer-Verlag, 2011 13 Gorodnitsky I F, Rao B D. Sparse signal reconstruction from limited data using FOCUSS: a re-weighted minimum norm algorithm. IEEE Trans Signal Process, 1997, 45: 600–616 14 Rao B D, Engan K, Cotter S F, et al. Subset selection in noise based on diversity measure minimization. IEEE Trans Signal Process, 2003, 51: 760–770 15 Hu C X, Meng H D, Li G, et al. An improved FOCUSS algorithm for multi-channel FMCW radar imaging of distributed targets. In: Proceeding of 2011 IEEE CIE International Conference on Radar, Chengdu, 2011. 677–680 16 Winkler V. Range Doppler detection for automotive FMCW radars. In: Microwave Conference, Munich, 2007. 1445– 1448 17 Potter L C, Ertin E, Parker J T, et al. Sparsity and Compressed Sensing in Radar Imaging. Proc IEEE, 2010, 98: 1006–1020 18 Lin Y G, Wu Y R, Hong W, et al. Compressive sensing in radar imaging. In: 2009 IET International Radar Conference, Guilin, 2009. 1–3 19 Candes E. The restricted isometry property and its implications for compressed sensing. C R Math, 2008, 346: 589–592 20 Candes E J, Tao T. Near-optimal signal recovery from random projections: universal encoding strategies? IEEE Trans Inf Theory, 2006, 52: 5406–5425 21 Tello M, Lopez-Dekker P, Mallorqui J J. A novel strategy for radar imaging based on compressive sensing. IEEE Trans Geosci Remote Sens, 2010, 48: 4285–4295 22 Lo Y T. A mathematical theory of antenna arrays with randomly spaced elements. IEEE Trans Antennas Propag, 1964, AP-12: 257–268 23 Yang J Y, Peng Y G, Xu W L, et al. Ways to sparse representation: An overview. Sci China Ser F-Inf Sci, 2009, 52: 695–703 24 Andres M, Feil P, Menzel W, et al. Analysis of automobile scattering center locations by SAR measurements. In: 2011 IEEE Radar Conference (RADAR), Kansas, 2011. 109–112 25 Stojanovic I, Cetin M, Karl W C. Joint space aspect reconstruction of wide-angle sar exploiting sparsity. In: SPIE Defense Security Symposium, Algorithms for Synthetic Aperture Radar Imagery XV, Orlando, 2008 26 Chellappa R. Two-dimensional discrete gaussian Markov random field models for image processing. In: Kanal L N, Rosenfeld A, eds. Progress in Pattern Recognition 2. Amsterdam: Elsevier Science Publishing Co., 1984. 79–112 27 Wu Y, Wang X, Xiao P, et al. Fast algorithm based on triplet Markov fields for unsupervised multi-class segmentation of SAR images. Sci China Inf Sci, 2011, 54: 1524–1533 28 Tropp J A, Gilbert A C. Signal recovery from random measurements via orthogonal matching pursuit. IEEE Trans Inf Theory, 2007, 53: 4655–4666 29 Derin H, Elliott H. Modeling and segmentation of noisy and textured images using Gibbs random fields. IEEE Trans Pattern Anal Mach Intell, 1987, PAMI-9: 39–55
. RESEARCH PAPER .
August 2012 Vol. 55 No. 8: 1776–1788 doi: 10.1007/s11432-012-4628-1
Special Issue
Improved FOCUSS method for reconstruction of cluster structured sparse signals in radar imaging HU ChenXi1,2 , LIU YiMin1 ∗ , LI Gang1 & WANG XiQin1 1Intelligent
Sensing Laboratory, Department of Electronic Engineering, Tsinghua University, Beijing 100084, China; 2Xi’an Electronic Engineering Research Institute, Xi’an 710100, China Received April 1, 2012; accepted May 5, 2012
Abstract The high resolution imaging of distributed targets in millimeter-wave radar system is studied in this paper. We use Gaussian Mrakov random field (GMRF) to represent the clustering property of the targets. Our novel method, called Clustering FOCUSS, incorporate an additional cluster constraint into the process of the focal underdetermined system solver (FOCUSS) algorithm. Simulation results indicate that the novel algorithm has a higher imaging accuracy than the methods of Capon beamforming, the l1 norm algorithm and the FOCUSS algorithm. Keywords
antenna array, high resolution imaging, cluster sparsity, GMRF, Clustering FOCUSS
Citation Hu C X, Liu Y M, Li G, et al. Improved FOCUSS method for reconstruction of cluster structured sparse signals in radar imaging. Sci China Inf Sci, 2012, 55: 1776–1788, doi: 10.1007/s11432-012-4628-1
1
Introduction
In automotive radar system, millimeter-wave radar is an indispensable forward-looking sensor to detect and estimate the relative angle, range and velocity of the targets in a multiple target scenario. In order to satisfy system requirements, millimeter-wave radar with multi-channel antenna array is proposed in [1]. This system is composed of a single transmitter to radiate a linear FM signal and M receiving channels to acquire echo signals to estimate the angle, range and velocity of the targets. Application of images from automotive radar depends on its image quality. Resolution is an important parameter to verify this. When the resolution is equal to the target size, the target is only a bright spot in the imaging result; when the resolution is much lower than the target size, the detailed information of the target could be obtained from the imaging result. There are two ways to enhance resolution of automotive radar. One way is to refresh and upgrade the hardware (such as increasing bandwidth and antenna aperture [2,3]). However, the method results in an increase in radar size and complexity of configuration, which is not suitable for automotive application. The other way is to use data processing algorithm to obtain superresolution. Various methods have been proposed to obtain high accuracy and resolution in radar imaging. Capon beamforming [4] and other covariance matrix based estimation methods perform well when received ∗ Corresponding
author (email: [email protected])
c Science China Press and Springer-Verlag Berlin Heidelberg 2012
info.scichina.com
www.springerlink.com
Hu C X, et al.
Sci China Inf Sci
August 2012 Vol. 55 No. 8
1777
signal from different targets are incoherent, but when the range and velocity of targets result in coherent signals, these covariance matrix based estimation methods will degrade dramatically. Considering that there are only a few targets in the scene discussed here, which means that the number of nonzero pixels in imaging result is small, sparse reconstruction algorithms are used for imaging targets in linear frequencymodulated continuous-wave (FMCW) automotive radar [5]. This method does not experience coherent signals problem and significantly improves the imaging accuracy. However, this method doesn’t exploit any structure information (except sparsity) that may exist in the sensed signals. Compared to the scale of the scene, the targets are usually not only sparse but also clustered in a distinct way. A few attempts have been made to utilize these cluster priors to get better sparse recovery. Refs. [6,7] argued that how to recover cluster sparse signal could be regarded as a mixed l2 /l1 -optimization problem. Ref. [8] extended MP and OMP to cluster sparse signal and got BMP and BOMP. Ref. [9] presented an extention of sparse dictionary learning to obtain dictionary elements belonging to groups of variables forming convex subsets. For simplicity, they all assume that the cluster structures (such as the cluster number/size/location) are known before recovery. Ref. [10] used overlapping neighborhoods to modify the sparsity term to better recover cluster sparse signal. However, this would lead the coefficients within each neighborhood to be smooth and cause the loss of details. In this paper, using the GMRF model [11,12] to describe the clustering of signal support, we propose a modification to standard FOCUSS [13,14] to make it more efficient for cluster sparse solutions. The novel method could work without knowing the details of cluster structures, which are not always available in automotive application. Simulation shows that imaging results get improved. We have already presented the Clustering FOCUSS method in [15], where we used it to improve the radar imaging of distributed targets with multi-channel FMCW radar. In this paper, we extend the application of the Clustering FOCUSS method to pulse radar, analyze the selection of the MRF model, and address some interesting issues.
2
Signal model from view point of sparse recovery
In this paper, for automotive radars we use antenna array to interrogate a stationary scene with pulses or continuous wave to reveal the spatial distribution of the scatterers. From viewpoint of compressive sensing (CS), we obtain the linear system of equations as Eq. (1) by discretizing the scene reflectivity function over range and azimuth on a grid of points {θi , Ri }, where q is a rearrangement vector of the sampled received baseband signal, α is the reflectivity vector, the ith column of Φ is the contribution of the target with azimuth θi and range Ri , and n is the noise. The model Eq. (1) can be easily adapted for use with the CS framework. q = Φα + n. (1) Then we model the reflected signal of FMCW radar with antenna array to obtain the linear equation in Eq. (1), which has been widely used in automotive application. We consider a uniform linear array (ULA) of MR antennas which share a single transmitter, as shown in Figure 1. The transmitted signal takes the following form: x0 (t) = exp{j2π[f0 t + (B/T )(t2 /2)]},
0 t < T,
(2)
where T is the pulse period, f0 is the carrier frequency and B is sweeping bandwidth. Assuming that there is a single target in the scenario, the echo signal received by the mth antenna is R − vt d sin θ xm (t) = ax0 2 − (m − 1) , m = 1, 2, . . . , M, (3) c c where a is the received signal’s amplitude, c is the velocity of light, d is the interelement ULA spacing, R is the initial radial range to the target at time t = 0 for the first transmitted pulse, v is the radial velocity of the target, and θ is the azimuth angle of the target.
1778
Hu C X, et al.
Figure 1
Sci China Inf Sci
August 2012 Vol. 55 No. 8
Diagram of multi-channel antenna array radar.
Mixing the echo signal xm (t) with the complex conjugate transmitted signal and sampling, the sampled data can be arranged in a “Data Matrix”, which is (considering B f0 , v c) Qm,n = a exp{−j2π[(m − 1)fθ + (n − 1)fR,v + ϕR ]},
m = 1, 2, . . . , M,
n = 1, 2, . . . , N,
(4)
where fθ = f0 (d sin θ/c) is referred to as angle frequency, fR,v = [(B/T )(2R/c) − f0 (2v/c)]Ts is range (velocity) frequency, ϕR = f0 (2R/c) − (B/2T )(2R/c)2 , Ts = 1/fs , fs is system sample frequency. It is (k) important to note that range R and velocity v are coupled in fR,v . This is range-velocity coupling problem of FMCW radar and some efficient methods to solve this problem are given in [16], so we estimate angle and range without considering the velocity of target in this paper. Discretizing angle and range axis, the angle-range domain can be divided into square grids [5]. Here we assume that there is a target in each grid which represents the echo signal with angle frequency and range frequency corresponding to that grid. Then the received data can be obtained by summarizing the contribution of all the targets with different angle and range frequency. We rearrange the data matrix [Qm,n ] to a vector, q = [Q1,1 , Q2,1 , . . . , QM,1 , . . . , QM,N ]T and Eq. (4) can be rewritten as follows: q=
M0 N0
αm0 ,n0 Φm0 ,n0 ,
(5)
m0 =1 n0 =1
where (m0 )
αm0 ,n0 = α(fθ
(m0 )
Φm0 ,n0 = Φ(fθ =
(n0 )
, fR
) = a(m0 ,n0 ) exp(−j2πϕR 0 ),
(n0 )
, fR
(n )
) = ΦR ⊗ Φθ
(6)
(m ) (m ) [1, exp(−j2πfθ 0 ), . . . , exp(−j2π(M − 1)fθ 0 ), (n ) (m ) . . . , exp(−j2π((N − 1)fR 0 + (M − 1)fθ 0 ))]T ,
M0 and N0 are the number of quantization of angle and range axis, respectively. In consideration of noise, model Eq. (5) can be written as matrix form as Eq. (1). Note that the row dimension M × N of matrix Φ is much less than column dimension M0 × N0 and most of elements in solution vector α are negligible and only a small portion of it, which represents contribution of targets, is remarkable. Therefore, this problem could be solved by sparse recovery technique. We can also obtain the CS model Eq. (1) of pulse radar by analyzing the process of pulse compression and beam forming. Different from the FMCW radar, each column of the matrix Φ represents the set of dechirped samples across the complete antenna array for a given point in the discretized scene. The detailed process could be found in [17,18]. According to CS theory, when the matrix Φ has the restricted isometry property (RIP), it is possible to recover the K large αi s from fewer measurements q [19]. If Φ is of large size, estimating and testing
Hu C X, et al.
Sci China Inf Sci
August 2012 Vol. 55 No. 8
1779
the RIP would be impractical. A tractable yet conservative bound on the RIP can be obtained through the mutual coherence of the columns of Φ [17]. Since the echo signal of one target has low coherence except with itself, the incoherence property of matrix Φ is verified. Furthermore, we can use randomly sampled data, which means that the measurement matrix is Fourier randomly measured ensembles, to recover sparse signals [20,21]. For antenna array radar, randomization of array element locations may function as a means of reducing the correlation among responses from different targets [22].
3
Image algorithm based on sparse recovery
3.1
The FOCUSS algorithm
According to the theory of sparse representation, when the actual distribution is sparse in a domain, the ill-posed problem in Eq. (1) can be efficiently solved using the technique of sparse representation. The basic form of sparse representation is defined as α = arg min α0 subject to q − Φα2 ε,
(7)
where · p stands for the lp norm and ε is the error allowance. However, this problem is a combinatorial problem and NP-hard. To address this difficulty, a number of practical algorithms have been proposed to approximate this sparse solution. Most algorithms for sparse representations can be classified into three categories: greedy pursuit algorithms, lp norm regularization based algorithms, and iterative shrinkage algorithms [23]. Greedy pursuit algorithms may be trapped into local minima and cannot be corrected back due to the greedy nature of algorithms. The lp norm regularization based algorithms demand a high computational effort especially for large-scale problems. So an iterative shrinkage method named FOCUSS [13,14] is discussed here. The FOCUSS method is based on the weighted least squares estimate. It is motivated by the observation that if a sparse solution is desired then choosing a solution based on the smallest l2 norm is not appropriate. The basic form of the recursive method is composed of the following steps. W (k) = diag(α(k−1) ),
(8)
z (k) = (ΦW (k) )† q,
(9)
α
(k)
=W
z
(k) (k)
,
(10)
where diag(·) is the operation of diagonalization, A† = (AH A)−1 AH denotes the pseudoinverse operation of matrix A, W (k) and α(k) denote the weighting matrix and the corresponding sparse solution at the kth iteration. In each iteration, the FOCUSS method gradually reinforces some of the already prominent entries in α(k) while suppressing the rest until they become close to zeros. Ref. [13] shows that FOCUSS serves as an iterative approximation of the sparse representation using lp , 0 p 1 norm minimization. By starting with the proper initial value, FOCUSS can converge to the global optimal value and provide a sparse representation. 3.2 3.2.1
The Clustering FOCUSS algorithm Cluster sparsity constraints
In literature about transform compression, many natural and man-made signals and images are considered as sparse or compressible with their sparse supports (set of nonzero coefficients) often having an underlying order. This paper is an attempt to exploit a priori information on coefficient structure of the cluster sparse signal to enhance the reconstruction accuracy and promote the ability of recovery. Figure 2 illustrates several examples of structured sparse support in the radar imaging application. Compared to the scale of the scene, the targets are usually not only sparse but also clustered in a distinct way. Nevertheless, this clustering property is not fully exploited by current CS recovery algorithms.
1780
Figure 2
Hu C X, et al.
Sci China Inf Sci
August 2012 Vol. 55 No. 8
(a) The location of the scattering centers at 24 GHz and 77 GHz of the vehicle rear end with a bandwidth
of 2 GHz [24]; (b)the backhoe CAD model and high resolution image of the backhoe for wide-angle imaging aperture of 110◦ [25].
Figure 3
A cluster sparse signal and its signal support. (a) α: a cluster sparse signal with length M = 40 and sparsity
D = 7; (b) f : the support of the signal α; (c) f : the support of a randomly sparse signal with length M = 40 and sparsity D = 7.
We use f ∈ {0, 1}M to represent the support of a cluster sparse signal α ∈ RM , where fi = 0 when αi = 0 and fi = 1 when αi = 0. For example, a cluster sparse signal α and its support f are shown in Figures 3 (a) and (b). This signal has sparsity D = 7, and its nonzero coefficients appear in only one cluster. The support of a randomly sparse signal with sparsity D = 7 is also shown in Figure 3 (c) (Here “randomly sparse” means that the signal is only sparse without any other additional structure constraint). By comparing f in Figures 3 (b) and (c), we find that the cluster sparse signal has the following properties: if the neighbors of an element of the signal are nonzero, the element itself is more likely to be nonzero. Clustering of the nonzero coefficients in a sparse signal representation can be realistically captured by a probabilistic graphical model such as a Markov random field. We use an MRF (S, N ) [12] to model the probability density function (PDF) of the cluster sparse signal support, where S denotes a set of M vertices (one for each of the support indices) and N denotes a neighborhood by which the sites in S are related to one another. N is defined as N = {Ni |∀i ∈ S} where Ni is the set of sites neighboring i. The set {0, 1}M from which f is taken is called the label set of the MRF.
Hu C X, et al.
Figure 4
Sci China Inf Sci
August 2012 Vol. 55 No. 8
1781
Elements of z at each iteration for the 4×10 matrix Φ and the vector α with only one nonzero element α4 .
Note that only one element z4 of z converges to one, whereas the rest become zero.
In MRFs, only neighboring vertices have direct interactions with each other, which means the following condition is satisfied: P (fi |fS−{i} ) = P (fi |fNi ), (11) where S − {i} is the set difference, fS−{i} denotes the set of labels at the sites in S − {i}, and fNi = {fi |i ∈ Ni } stands for the set of labels at the sites neighboring i. The MRF models can be further classified according to assumptions made about individual fi . 3.2.2
The procedures of Clustering FOCUSS algorithm
In FOCUSS, while the entries of α(k) converge to zero and nonzero values, the corresponding entries (k) (k) (k) (k) in z (k) converge to 0 or 1, i.e., zi → 0 as αi → 0, and zi → 1 as αi approach nonzero values (Figure 4). Therefore, we suppose that the signal support f could be approximated by z. Here we select the GMRF model [11,12] from the class of MRF models to describe the characteristic of z. The label set from which the value of zi is taken is the real line and the joint distribution is multivariate normal. Its conditional PDF is 2 1 1 exp − 2 zi − βi,j zj p(zi |zS−{i} ) = p(zi |zNi ) = , (12) 2σ (2πσ 2 )1/2 j∈Ni
where σ 2 is the variance of a Gaussian noise, and βi,j are constants reflecting the pair-site interacton between i and j: βj,i = βi,j , βj,i = 0 if and only if sites i and j are neighbours. As shown in Figure 5, the conditional PDF of zi is related to the values of its neighbors. For example, zi in Figure 5 (a) whose neighbors are both nonzero is more likely to be nonzero than in Figure 5(c) whose neighbors are both zero. The GMRF model is widely used in many scientific fields to represent interaction among spatial data with a Gaussian noise [26,27]. There are two main reasons for our choosing the GMRF model. Firstly, the continuous label set of the GMRF model is more proper than those models with discrete label set since z converges to f during iterations. Secondly, the form of conditional p.d.f. is normal distribution, which makes our novel algorithm have simple form and low computational cost. Thus we have the following iterative relaxation algorithm: 2 1 1 (k) (k−1) (k−1) (k−1) )= exp − 2 zi − βi,j zj . (13) p(zi |z 2σ (2πσ 2 )1/2 j∈Ni
1782
Hu C X, et al.
Figure 5
Sci China Inf Sci
August 2012 Vol. 55 No. 8
The pdf of zi with different neighbor conditions.
(k)
(k−1)
To simplify calculations, zj in the right-hand side of Eq. (13) is replaced by zj . Since z (k) cannot express the signal support of f properly at first, the cluster constraint should be added after the first steps. In FOCUSS, Eq. (9) is the smallest l2 norm solution to the problem ΦW z = q, and now it can be regarded as the MAP estimate, where the likelihood function p(q|z) is expressed in terms of the joint Gaussian distribution. 1 1 H −1 tr{(q − ΦW z) exp − Σ (q − ΦW z)} . (14) p(q|z) = 0 2 ((2π)M det{Σ0 })1/2 Here we assume that n in Eq. (5) is independent and identically distributed (i.i.d.) random variables from a normal distribution with a constant variance σ02 , such that σ02 IM ∈ RM×M . Applying the Bayesian framework to the prior Eq. (13), the posterior distribution has the form p(z|q, z (k−1) ) =
p(q|z)p(z|z (k−1) ) . p(q|z (k−1) )
(15)
The MAP estimate in the kth iterative update can be obtained by solving the following regularized least squares minimization problem: z (k) = arg min Ψ (z|q, z (k−1) , λ),
(16)
Ψ (z|q, z (k−1) , λ) = −2 ln p(z|q, z (k−1) , λ) = q − ΦW z (k) 22 + γz (k) − Bz (k−1) 22 ,
(17)
z
where
γ = σ02 /σ 2 , and B is the square matrix whose diagonal elements are unity and whose off-diagonal (i, j) element is βi,j . The stationary point of Eq. (16) expressed in the vectorized form is as follows: z (k) = z∗ = ((ΦW )H (ΦW ) + γI)−1 ((ΦW )H q + γBz (k−1) ).
(18)
The procedures of the Clustering FOCUSS algorithm proceed as follows. 1) Obtain the initial value as a low-resolution estimation computed as (0)
αi
= ΦH i q,
1 i M 0 N0 ,
(19)
Hu C X, et al.
Sci China Inf Sci
August 2012 Vol. 55 No. 8
1783
W (0) = diag(abs(α(0) )),
(20)
Γ = {i, 1 i M0 N0 }.
(21)
2) Update the weighting matrix and the estimation result accordingly as z (k) |Γ = ((Φ|Γ W (k) )H (Φ|Γ W (k) ) + γI)−1 ((Φ|Γ W (k) )H q + γ(Bz (k−1) )|Γ ),
(22)
α
(23)
(k)
|Γ = W
z
(k) (k)
|Γ ,
where α(k) |Γ stands for the Γ subset of the vector α(k) . 3) Update the adaptive subspace using the current solution as (k)
Γ = arg(|αi | Th), 1 i M0 N0 ,
(24)
(k)
where αi is the ith element of the vector α(k) , Th stands for the threshold and the set Γ is the adaptive subspace at the kth iteration.
4 4.1
Simulation results Simulation examples and analysis of Clustering FOCUSS (k)
In the second step of Clustering FOCUSS, zi in Eq. (18) is composed of two terms. The first term is the same as that in FOCUSS, and the second term is the cluster sparsity constraints. Note that if (k−1) (k) have nonzero values, the second term of zi would be greater than zero. In other neighbors of zi words, if zi has nonzero neighbors at the (k − 1)th iteration, it would be larger in Clustering FOCUSS than in FOCUSS at the kth iteration. From this we can reach an interesting conclusion that Clustering FOCUSS is more likely to reinforce or at least suppress sites with nonzero neighbors more slowly, so that the nonzero coefficients would prefer to be clustered. (k−1) , the sites in a cluster which contains targets will reinforce If there has already been clusters in zi each other through cluster sparsity constraints and the sites in a cluster which contains no target will be suppressed by the weighted l2 norm constraints. Simulation is carried out to demonstrate the efficiency of the proposed algorithm. To test the analysis above, we assume that the cluster structured sparse signal α is randomly generated with length L, sparsity D and clusters K. Its nonzero entries are drawn from a uniform distribution. The columns of matrix Φ are drawn form a standardized Gaussian distribution and then normalized to unit l2 norm. Here two pixels are defined as neighbors only when the distance between them is less than the minimum resolution and βi,j = 1/S with S is the number of neighbors of αi . The noisy data was generated by adding to the exact mixtures a zero-mean Gaussian noise with σ02 adjusted to a desired value of SNR0 (SNR0 = 20 log10 (qexact /n)). The signal shown in Figure 6 is generated with L = 250, D = 40, K = 3. The values of z computed by Clustering FOCUSS are shown along with the iterations in Figure 7. As is evident in the graph, the Clustering FOCUSS algorithm could cluster the sparse coefficients of signal recovery, reinforce clusters containing real targets and suppress clusters with no target. The recoveries of FOCUSS and Clustering FOCUSS are both shown in Figure 7. In the reconstruction of Clustering FOCUSS, the noise is suppressed well and the weak cluster of the signal is reconstructed correctly. 4.2
Statistics of Clustering FOCUSS
In order to verify the performance of Clustering FOCUSS, we take into account the number of measurements P and SNR0 which are able to reach the successful reconstruction. The reduction of measurements P means that we can get high resolution of the scene with sub-Nyquist measurement rates and less antennas. SNR0 reflects the robustness of the performance of Clustering FOCUSS to noise. They are both impressive in the theory of sparse recovery.
Hu C X, et al.
1784
Figure 6
Sci China Inf Sci
August 2012 Vol. 55 No. 8
(a) Reconstruction with FOCUSS; (b) reconstruction with Clustering FOCUSS; (c) original signal.
Figure 7
Elements of z in Clustering FOCUSS at the (a) 1st; (b) 5th; (c) 10th; (d) 20th iteration.
Sparse signals with 2 clusters are considered. For fixed cluster structured sparse signal(L = 250, D = 30, K = 2), the number of measurements P is ranging from 60 to 140 with step 10, and SNR0 is ranging from 10 dB to 26 dB with step 2 dB. For each step, we run the program 100 times independently with different sensing matrix. The accuracy of reconstruction is evaluated through relative root mean ˆ Na is the number of square error (RRMSE) between original sparse signal α and its reconstruction α. independent repeated experiments, and α ˆi is the reconstruction of the ith repeated experiment.
RRMSE = 20 log10
Na 1 α − α ˆi . Na i=1 α
(25)
The plots in Figure 8 show the dependence of RRMSE on the parameters: P and SNR0 . Figure 9(a) shows RRMSE plot versus SNR0 at P = 100; Figure 9(b) shows RRMSE plot versus P at SNR0 =15 dB. It is easy to notice that Clustering FOCUSS can give the estimates that have lower RRMSE values than
Hu C X, et al.
Figure 8
Sci China Inf Sci
August 2012 Vol. 55 No. 8
1785
RRMSE plots versus the number of measurements P and SNR0 . (a) OMP; (b) FOCUSS; (c) Clustering
FOCUSS.
Figure 9
(a) RRMSE plots versus SNR0 at the number of measurements P = 100; (b) RRMSE plots versus the number
of measurements P at SNR0 =15 dB.
OMP (which is a greedy algorithm for compressive sensing recovery [28]) and FOCUSS, even with less measurements P and lower SNR0 . We plot the RRMSE as a function of the cluster-sparsity level of the signal to be recovered in Figure 10 at the number of the measurement P = 100, SNR0 =15 dB. Here the clusters are of the fixed size 10, and the locations of the nonzero clusters are chosen uniformly among all possible locations. We can see that Cluster FOCUSS outperforms FOCUSS and OMP especially when the cluster-sparsity level increases.
5
Discussion on the performance of Clustering FOCUSS in real scene
Assume a ULA of MR = 9 antennas with a half-wavelength interelement spacing and a linear FM signal with carrier frequency f0 =24 GHz, the sweeping bandwidth B=100 MHz, and the pulse period T =
Hu C X, et al.
1786
Figure 10
Figure 11
Sci China Inf Sci
August 2012 Vol. 55 No. 8
RRMSE plots versus cluster-sparsity level at measurement number P =100 and SNR0 =15 dB.
The imaging results by different methods with SNR0 =10 dB. (a) Origin targets; (b) conventional beamforming;
(c) Capon beamforming; (d) L1 norm algorithm; (e) FOCUSS algorithm; (f) Clustering FOCUSS algorithm.
0.75 ms. The sample frequency of the system is fs =60 kHz, which satisfies Nyquist sample rate. SNR0 is 10 dB. We consider the imaging result of a backhoe with wide-angle SAR imaging [25] as the real scene. The detailed information of the targets is shown in Figure 11(a), where x-axis is range plane and yaxis is angle plane. The scene is divided into square grids whose dimensions are 1/4 of the range and azimuth resolution of the multi-channel FMCW radar. Assume that there is a target in each grid and the reflectivity of the target equals the amplitude of the pixel in the same position of SAR image shown in Figure 11(a). Then the received data is obtained by summarizing the contribution of all the targets in the scene according to the model in Eq. (5). Figure 11(b) is the imaging result using conventional beamforming. Since the resolution is limited by antenna aperture and sweep bandwidth, the boundaries of the target are extended and the details in the target cannot be recovered correctly. Figure 11(c) is the imaging result using Capon beamforming. By Eq. (5), two close-set pixels will result in coherent signals which unfortunately happen in our problem.
Hu C X, et al.
Sci China Inf Sci
August 2012 Vol. 55 No. 8
1787
In this case, Capon beamforming cannot distinguish each pixel well. Figure 11(d) is the imaging result using the l1 norm algorithm. Disjunctions appear instead of continuous distribution which is the truth. The “superresolution” phenomenon is caused by coherent columns of the matrix Φ. Figure 11(e) is the imaging result using the FOCUSS algorithm. The position of target is reconstructed correctly; however, in the reconstructed solution, weak parts of the origin targets almost disappear and strong parts are smaller. The reason is that the noise affects the elements of solution differently due to the amplitude of the pixels, which may amplify strong pixels in the solution and significantly distorts the rest. The imaging result of Clustering FOCUSS algorithm is shown in Figure 11(f). By utilizing the clustering property of the nonzero elements of α, the effect of the noise we talk about above is mitigated. It is apparent that the noise is suppressed well and the weak cluster of the signal is reconstructed correctly. In the reconstruction of Clustering FOCUSS, note that since the whole image is more likely to be continuous and clustered because of the cluster sparsity constraints, the isolated points of the origin target are also clustered, which is not the truth. The reason for this phenomenon is that the GMRF model we use here could not deal well with the discontinuous region (such as edge) of the scene. To solve this problem, one way is to incorporate Clustering FOCUSS with the ability to detect the edge and modify the reweighting rule there. Another way is to use a more accurate model such as the hierarchical two-level Gibbs model [29] which has a higher-level Gibbs distribution to characterize the cluster sparsity of the signal and a lower-level Gibbs distribution to describe the filling-in in each cluster.
6
Conclusion
A novel method based on block-sparse model has been proposed for imaging distributed targets with millimeter-wave radar with multi-channel antenna array. Utilizing the cluster sparseness in angle-range domain, our method named Clustering FOCUSS regards the weighted l2 norm minimization in the second step of the FOCUSS algorithm as the MAP estimator, and incorporates the GMRF graphical model to constrain the nonzero values of the solution. Compared with FOCUSS, Clustering FOCUSS can keep the recursive process from the error accumulation between iterations and prevent an over-focal solution; the distributed targets which have clustered signal support can also stand out even with lower SNR. The Clustering FOCUSS algorithm can also be used in other radar scenes which have cluster sparse property, such as SAR and ISAR. Much work needs to be done, however. We are working to analyze the global convergence of the Clustering FOCUSS algorithm. We also assert that other probabilistic signal models could be used to improve the existing signal recovery methods.
Acknowledgements This work was supported in part by National Natural Science Foundation of China (Grant No. 40901157), National Basic Research Program of China (973 Program) (Grant No. 2010CB731901), and Tsinghua National Laboratory for Information Science and Technology Cross-discipline Foundation.
References 1 Lee M S, Katkovnik V, Kim Y H. System modeling and signal processing for a switch antenna array radar. IEEE Trans Signal Process, 2004, 52: 1513–1523 2 Johnson D H, Dudgeon D E. Array Signal Processing: Concepts and Techniques. Englewood Cliffs: Prentice-Hall, 1993 3 Van Trees H L. Detection, Estimation and Modulation Theory, Part IV: Optimum Array Processing. New York: John Wiley & Sons, Inc, 2002 4 Capon J. High-resolution frequency-wavenumber spectrum analysis. Proc IEEE, 1969, 57: 1408–1418 5 Wang P, Meng H D, Wei Y M. FMCW radar imaging with multi-channel antenna array via sparse recovery technique. In: 2010 International Conference on Electrical and Control Engineering (ICECE), Wuhan, 2010. 1018–1021
1788
Hu C X, et al.
Sci China Inf Sci
August 2012 Vol. 55 No. 8
6 Stojnic M, Parvaresh F, Hassibi B. On the reconstruction of block-sparse signals with an optimal number of measurements. IEEE Trans Signal Process, 2009, 57: 3075–3085 7 Baraniuk R G, Cevher V, Duarte M F, et al. Model-based compressive sensing. IEEE Trans Inf Theory, 2010, 56: 1982–2001 8 Eldar Y C, Mishali M. Robust recovery of signals from a structured union of subspaces. IEEE Trans Inf Theory, 2009, 55: 5302–5316 9 Jenatton R, Obozinski G, Bach F. Structured sparse principal component analysis. In: Proceedings of the 13th International Conference on Artificial Intelligence and Statistics (AISTATS), Sardinia: 2010 10 Kavukcuoglu K, Ranzato M A, Fergus R, et al. Learning invariant features through topographic filter maps. IEEE Conference on Computer Vision Pattern Recognition, Miami, 2009. 1605–1612 11 Besag J. Spatial interaction and the statistical analysis of lattice systems (with discussion). J Royal Statist Soc Ser B, 1974, 36: 192–326 12 Li S Z. Markov Random Field Modeling in Image Analysis. 2nd ed. New York: Springer-Verlag, 2011 13 Gorodnitsky I F, Rao B D. Sparse signal reconstruction from limited data using FOCUSS: a re-weighted minimum norm algorithm. IEEE Trans Signal Process, 1997, 45: 600–616 14 Rao B D, Engan K, Cotter S F, et al. Subset selection in noise based on diversity measure minimization. IEEE Trans Signal Process, 2003, 51: 760–770 15 Hu C X, Meng H D, Li G, et al. An improved FOCUSS algorithm for multi-channel FMCW radar imaging of distributed targets. In: Proceeding of 2011 IEEE CIE International Conference on Radar, Chengdu, 2011. 677–680 16 Winkler V. Range Doppler detection for automotive FMCW radars. In: Microwave Conference, Munich, 2007. 1445– 1448 17 Potter L C, Ertin E, Parker J T, et al. Sparsity and Compressed Sensing in Radar Imaging. Proc IEEE, 2010, 98: 1006–1020 18 Lin Y G, Wu Y R, Hong W, et al. Compressive sensing in radar imaging. In: 2009 IET International Radar Conference, Guilin, 2009. 1–3 19 Candes E. The restricted isometry property and its implications for compressed sensing. C R Math, 2008, 346: 589–592 20 Candes E J, Tao T. Near-optimal signal recovery from random projections: universal encoding strategies? IEEE Trans Inf Theory, 2006, 52: 5406–5425 21 Tello M, Lopez-Dekker P, Mallorqui J J. A novel strategy for radar imaging based on compressive sensing. IEEE Trans Geosci Remote Sens, 2010, 48: 4285–4295 22 Lo Y T. A mathematical theory of antenna arrays with randomly spaced elements. IEEE Trans Antennas Propag, 1964, AP-12: 257–268 23 Yang J Y, Peng Y G, Xu W L, et al. Ways to sparse representation: An overview. Sci China Ser F-Inf Sci, 2009, 52: 695–703 24 Andres M, Feil P, Menzel W, et al. Analysis of automobile scattering center locations by SAR measurements. In: 2011 IEEE Radar Conference (RADAR), Kansas, 2011. 109–112 25 Stojanovic I, Cetin M, Karl W C. Joint space aspect reconstruction of wide-angle sar exploiting sparsity. In: SPIE Defense Security Symposium, Algorithms for Synthetic Aperture Radar Imagery XV, Orlando, 2008 26 Chellappa R. Two-dimensional discrete gaussian Markov random field models for image processing. In: Kanal L N, Rosenfeld A, eds. Progress in Pattern Recognition 2. Amsterdam: Elsevier Science Publishing Co., 1984. 79–112 27 Wu Y, Wang X, Xiao P, et al. Fast algorithm based on triplet Markov fields for unsupervised multi-class segmentation of SAR images. Sci China Inf Sci, 2011, 54: 1524–1533 28 Tropp J A, Gilbert A C. Signal recovery from random measurements via orthogonal matching pursuit. IEEE Trans Inf Theory, 2007, 53: 4655–4666 29 Derin H, Elliott H. Modeling and segmentation of noisy and textured images using Gibbs random fields. IEEE Trans Pattern Anal Mach Intell, 1987, PAMI-9: 39–55
