Compressive Sensing-Based Bandwidth Stitching for Multichannel Microwave Radars.

The problem of obtaining high range resolution (HRR) profiles for non-cooperative target recognition by coherently combining data from narrowband radars was investigated using sparse reconstruction techniques. If the radars concerned operate within different frequency bands, then this process increases the overall effective bandwidth and consequently enhances resolution. The case of unknown range offsets occurring between the radars' range profiles due to incorrect temporal and spatial synchronisation between the radars was considered, and the use of both pruned orthogonal matching pursuit and refined l 1 -norm regularisation solvers was explored to estimate the offsets between the radars' channels so as to attain the necessary coherence for combining their data. The proposed techniques were demonstrated and compared using simulated radar data.

1. Introduction

The construction of high range resolution profiles (HRRP) of targets is a precursor to feature extraction for automatic target recognition (ATR), and normally requires the employment of a high-bandwidth waveform following detection by a lower resolution radar mode. Examples of recent papers in the non-cooperative target recognition (NCTR) literature focusing on feature extraction for ATR following HRRP construction are [1,2,3]. This paper considers the problem of HRRP construction, but using low resolution radars operating in different frequency bands for the purpose of combining their signals to achieve a higher resolution, and examines the problem of their data not being mutually coherent.

The ability to acquire high resolution range profiles of targets has improved over time as hardware capability has developed, with higher resolution being achieved by increasing the time-bandwidth product. In the early approaches, for narrowband radars with very limited instantaneous bandwidth, stepped-frequency waveforms were used with a single I,Q (that is, baseband quadrature) signal sample received after each frequency step. The set of samples is effectively used for the Fourier transform of the slant range profile, enabling the range response to be obtained simply by implementing an inverse Fourier transform (see, e.g., [4]). The greater the frequency range, the higher the range resolution, but the downside is that the burst of pulses can be so long that a scatterer may migrate between range cells, causing smearing of the range profile, and therefore requiring range compensation. An example of a recent paper involving the use of stepped-frequency waveforms is [5], and recent papers which have investigated the effects of target motion and aspect sensitivity are [6,7].

An advance on the stepped-frequency approach is to increase the time-bandwidth product using stretch processing, whereby a wideband LFM waveform is transmitted, and pulse-compression is achieved in hardware by mixing the received signal with an extended replica of the transmitted waveform. A point target at a particular range will manifest itself as a single frequency which is proportional to its range. The result, which is digitally sampled in time in order to facilitate the identification of the frequency components, is therefore a superposition of discrete frequencies, each corresponding to a point target at a different range. An application of the inverse Fourier transform again recovers the range profile (see, e.g., [4,8,9]).

A spectral analysis technique to improve range resolution was proposed in [10,11] based on autoregressive linear prediction. The main idea is to combine the mutually coherent signals received from multiple waveforms transmitted sequentially or concurrently, which have widely separate carrier frequencies or may even occupy entirely different frequency bands. Viewed in the spectral domain, the received signals from individual waveforms may be seen to occupy discrete wavebands which are separate or contiguous. If contiguous, then they can potentially be coherently combined to synthesise the signals that would have been received from a single wider bandwidth waveform in the manner presented in [12]. If separate, then presumably this coherent combination of signals would still be feasible, as would be the interpolation of the frequency response in the gaps between the bands under the a priori assumption that the signals are returned from discrete scatterers using spectral estimation techniques (see, e.g., [13]). Alternatively, the signals from different frequency bands can be jointly processed without explicitly filling the gaps between the bands. Since no new synthetic frequency band is actually constructed, this approach can be referred to as bandwidth stitching to distinguish it from bandwidth interpolation and extrapolation. The main challenge of this approach lies in the presence of phase errors in different frequency bands resulting from post-processed motion compensation which is often carried out separately for each frequency band. In this paper, we focus on the problem of bandwidth stitching for radar high-resolution range profiling and explore the use of sparse reconstruction to deal with the phase error problem.

It is convenient to formulate these ideas in the spectral domain, within which point scatterers appear as discrete sinusoids and which are amenable to analysis by spectral estimation techniques such as autoregression, as presented in [10,11]. Compressive sensing and sparse reconstruction, however, provide for the possibility of alternative signal representations, potentially allowing for greater flexibility and discriminating between signals of physical origin and receiver noise [13,14,15,16,17,18,19,20]. Instances of non-sinusoidal signals are waveforms in fast-time and signals returned from rotating objects when the angle of rotation is large. Compressive sensing and sparse reconstruction can also handle the situation corresponding to data being non-uniformly sampled in time or space, such as non-uniform PRF (pulse repetition frequency) waveforms and random sparse arrays.

Compressive sensing and sparse reconstruction were exploited in [17,19,21,22] to address the problem of gaps in the data both in slow-time and in frequency for inverse synthetic aperture radar (ISAR) imaging. However, these works assumed that the data were coherent across different sub-bands and that there were no model uncertainties. The work [23] took account of the possible lack of mutual coherence between the radars operating on the different sub-bands arising from incorrect timing synchronisation, or, equivalently, errors in antenna phase’s centre-relative locations. This is achieved by fitting an ultra-wideband all-pole signal model to the mutually-coherent sub-bands, which is then used for bandwidth interpolation and extrapolation prior to recovering the range profile by means of an inverse Fourier transform. This paper, however, proposes the use of compressive sensing and sparse reconstruction to deal with the non-coherence problem between different sub-bands.

To address the sub-band non-coherence problem, two different approaches were explored: (i) greedy pursuit and (ii) -norm regularisation. In the first approach, pruned orthogonal matching pursuit (POMP) [24], which was originally developed for micro-Doppler parameter estimation, is adopted to deal with the dictionary mismatch which is due to the phase errors in each sub-band resulting from the motion-compensation post processing. The main idea is to parameterise the dictionary as a function of the phase errors and to construct multiple realisations of the dictionary. A selective learning process is then used to discard the dictionaries which correspond to incorrect values of phase errors. Since a straight application of the POMP algorithm to the problem under consideration would have been computationally expensive, we first applied the POMP algorithm pairwise to sub-bands in order to estimate the phase errors, and then utilised the conventional OMP algorithm to determine the range profile based on the estimated phase error values. In the second approach, an -norm regularisation problem can be solved jointly both for the range profile vector and the phase errors [25]. The work [25] offers two solutions for joint synthetic aperture radar (SAR) imaging and phase error correction. The first solution is not applicable to the problem under consideration because a constant phase error in each sub-band is assumed. On the other hand, the second solution considers general arbitrary phase errors, and thus can be applied to our problem. Since the phase errors within each sub-band are a linear function of the range error coming from post-processing motion compensation, we also present refined variants of the second solution of [25] to take into account this underlying structure of the phase error.

The paper is organised as follows. Section 2 formulates the problem of bandwidth stitching for HRRP in the presence of phase errors. The POMP algorithm is applied in Section 3 to the bandwidth stitching problem under consideration. Section 4 presents -norm regularisation solvers. Numerical performance comparisons are provided in Section 5 and conclusions are drawn in Section 6.

2. Problem Formulation

Consider a multistatic radar system consisting of M radar channels on different and distinct frequency sub-bands, approximately co-located, and illuminating a common target, such that their radar lines of sight (LoS) coincide but their range profiles are out of alignment. Each channel can individually produce a one dimensional range profile of the target, but with relatively coarse resolution. The bandwidth stitching problem can be briefly stated as follows: for the M generally non-coherent channels, the aim is to coherently combine, or “stitch,“ the channels together so that they can effectively produce a single range profile with resolution corresponding to the combined overall signal bandwidth.

Let () denote the frequency bin of the channel (). Here, is the number of frequency bins in the channel. We base the formulation on the point-scatterer model and assume that the target can be defined as consisting of K scattering centres at local line-of-sight coordinates (or local “down ranges”) and having complex-valued reflectivity coefficients , which are also assumed frequency-independent. The down-converted, pulse-compressed, motion-compensated signals received in each channel, in the frequency domain, can be written as where superscript T denotes the transpose operation, and

Here, represents the transmit radar waveform, the squared amplitude resulting from pulse compression processing; constant c denotes the speed of light; and accounts for the range errors in the motion-compensation processing. Bandwidth stitching in this context amounts to estimating these phase errors as accurately as possible.

The signals in (1) can be rewritten in a more compact form as where

Here, “diag” denotes a diagonal matrix; is referred to as the phase error matrix, of dimension ; and are respectively, dimensions and ; is a column vector of dimension ; and the dagger symbol † refers to the K actual scatterers on the target.

To apply the sparse representation techniques of compressive sensing, we discretise the target’s local range coordinate x using a regularly-spaced range grid for , with , and construct the dictionary matrices where

Are the “atoms” of dictionary in the frequency domain. The corresponding range profile vector

Spans over the range grid . The received signal can also be written as

Since the target usually contains only a small number of dominant scattering centres relative to the total number of range resolution cells, the range profile can be considered sparse (i.e., containing a small number of non-zero elements).

In the presence of unknown noise, (11) becomes where is the additive noise for channel m. Stacking up the individual channel signals , , gives where

Note that and are column vectors of size ; diagonal phase error matrix is of size ; dictionary matrix is ; and the range profile is again a column vector of size . Stacking the received signals amounts to a vertical stacking of the dictionary matrices from all channels and a diagonal concatenation of the corresponding phase error matrices. The stacking of multiple channels in this manner can improve the estimation accuracy for , as will be demonstrated later in the paper.

A problem statement can thus be expressed as follows: given as the measured signal,

The estimation of over is the process of range profiling, giving the main desired output, whereas the estimation of the phase error matrix is really only a necessary intermediate result; it is a function of and (recall that is the range estimation error resulting from the motion-compensation process for channel m). Furthermore, since these errors arise from a lack of precise knowledge of the relative locations of the radar channels’ phase centres and are small relative to a range resolution cell, we may assume, without loss of generality, that .

3. Greedy Pursuit Solutions

In this section, we adopt the pruned OMP (POMP) technique, which was originally proposed for micro-Doppler parameter estimation, [24], for the problem of bandwidth stitching for range profiling. We start with the simplest case of two channels and then generalise it to the multiple channel case.

3.1. The Two-Channel Case

For this case, the signal model in (13) can be expressed as where is a function of the single unknown relative range error , with

In addition to the sparse range profile vector , is the only additional unknown parameter to be estimated. Let us rewrite (19) as where

In this form, the problem can be viewed as a joint sparse reconstruction and parameter estimation problem with the parametric dictionary itself a function of the parameter . This can be considered as a special dictionary learning problem where the objective is to solve simultaneously for both the sparse solution of and the range error .

To solve this problem, we adopt the POMP technique [24], which embeds a pruning operation into the iterative process of OMP. The main idea of POMP is to construct multiple realisations of the dictionary based on a number of candidate values of ; the OMP algorithm is applied to each dictionary realisation to find the atom which correlates most strongly with the current residual for that dictionary, and to recompute the residual with that atom’s contribution to the residual removed. To overcome possible excessive computations arising from outlier candidate values of , a pruning operation is performed to exclude the half of the dictionaries which yield the largest residual errors, until a single dictionary realisation remains. The OMP iterations for the remaining dictionary are continued until a termination criterion is satisfied. The candidate value of corresponding to this dictionary gives the final estimate of . The basic POMP algorithm is summarised in Table 1 and its computational cost is shown in Appendix C to be of the order of .

Table 1

The pruned orthogonal matching pursuit (POMP) algorithm (M = 2).

INPUT: Noisy signal data vector S˜ Candidate dictionaries Φ_1,Φ_2,…,Φ_L_Δ, corresponding to L_Δ candidate values of ΔR_2 PROCEDURE: Initialization: - set the initial indexes of active dictionaries to Θ_1={1,…,L_Δ}; - set the corresponding residual vectors to r_1=⋯=r_L_Δ=S˜; - set the initial support Λ to ∅, the null set; fori=1; i:=i+1 until |Θ_i|==1 (the cardinality of Θ_i) and |r_l|<ϵ for l∈Θ_i, for every l∈Θ_i, perform OMP as follows - Identify: c_l=ΦlHr_l j_l=arg max_j|c_j| - Merge supports: Λ_l=Λ_l∪j_l - Update^*: α^_l,Λ_l=Φl,Λ_lHΦ_l,Λ_l^−1Φl,Λ_lHS˜ r_l=S˜−Φ_l,Λ_lα^_l,Λ_l end for if|Θ_i|>1 Remove indices of ⌈|Θ_i|/2⌉ candidate dictionaries that correspond to ⌈|Θ_i|/2⌉ largest residual errors from |Θ_i| end if   end for OUTPUT: The range profile estimate α^_l^★,Λ_l^★ where l^★ is the last remaining element of Θ_i The estimate of ΔR_2 is the value of ΔR_2 corresponding to Φ_l^★

consists of the columns of with indices belonging to and consists of the elements of with indices belonging to .

3.2. The Multi-Channel Case

For this case, the noisy signal model (19) becomes where

Here, the dictionary matrix is a function of the unknowns and .

The POMP algorithm could be extended to multiple channels by computing candidate dictionaries based on a multi-dimensional grid of candidate values for and . The grid would consist of a total of dimensions, where the mth dimension corresponds to the unknown range error of the channel. Note that only a one-dimensional grid for is required for the case of two channels. However, the cardinality of the dictionary set is exponentially dependent on the number of available channels; i.e., , where V denotes the number of grid points in each parameter dimension. As a result, although this extension would be simple and straightforward, it is computationally expensive.

To alleviate this computational burden, we instead apply the POMP algorithm pairwise to channels in order to estimate the range errors relative to the first channel, and then utilise the conventional OMP algorithm to determine the range profile vector based on these estimated values of . The procedure is summarised as below:

STEP 1: Estimation of range errors.

For each pair between the 1st and channel, :

Calculate input signal:

Construct candidate dictionaries based on a grid of L candidate values of ():

Perform POMP given and to obtain an estimate of (denoted as ).

STEP 2: Estimation of range profile vector.

Compute signal and dictionary.

Estimate using OMP given and .

The computational cost of the general POMP algorithm is shown in Appendix C to be of the order of .

4. -Norm Regularisation Approach

The sparse reconstruction problem (18) can be solved via the following regularised optimisation: where is a regularisation parameter. It should be emphasized that this is not a conventional regularisation formulation because of the unknown phase error matrix resulting from the estimation error of the motion-compensation phase. Therefore, must be jointly estimated with :

Two solutions for this joint estimation problem were presented in [25]. The first solution assumes that the phase error matrix for the sub-band is modelled as

In other words, the phase errors for different frequency bins of a particular sub-band are identical. However, this assumption is invalid in the problem under consideration because the phase error is a function of frequency and thus has different values for different frequency bins. Therefore, that solution is not applicable in this case. The second solution considers a general phase error matrix where the phase errors can be arbitrary. Although this solution can be used, it does not exploit the underlying structure of the phase errors; i.e., . In what follows, we will also present other refined versions, building on the second solution of [25], while exploiting prior knowledge of the structure of the phase error.

The norm can be approximated as [26,27,28,29]:

In order to overcome the nondifferentiably of the norm at the origin. Here, is a small non-negative parameter. Using this approximation, the minimisation problem in (31) becomes

The solution of (35) tends to the solution of (31) as approaches zero. Therefore, a small value of should be used to ensure the validity of this approximation. The quasi-Newton approach can be adopted to solve the modified regularised optimisation (35), as below.

The gradient of the objective function of (35) is given by where the superscript denotes the Hermitian transpose operation. Here, is the Hessian matrix given by where

Since the Hessian matrix is a function of the unknown , the minimisation (35) is solved iteratively. Given the estimates and from a previous iteration i, the new solutions at iteration are obtained in the following two steps.

Calculate by setting given and :

Calculate given . The phase error matrix is obtained by solving: or equivalently for . Note that (since ); thus, no estimation is required for .

The algorithm may be halted when the objective function falls below a threshold, or when a maximum number of iterations is reached, or when the relative change in the objective function falls below a threshold.

Various methods for calculating the phase error matrix in Step 2 are given in the following sections.

4.1. Unstructured Approach

Ignoring the underlying structure of the phase errors, i.e., , can be considered as a diagonal matrix with arbitrary elements :

Therefore, can be estimated as [25] where and denote operations to extract the imaginary and real parts of a complex number, and stands for a four-quadrant arctangent operation. Here, is the element of which is defined as . As a result, we obtain:

4.2. Gauss–Newton Approach

Taking into account the underlying structure of the phase errors, is a function of , and the minimisation (41) can be re-expressed as

By letting , we have where and

Here is the element of . As we are estimating real quantities, it is more convenient to reformulate the problem as the minimisation of a real function in order to apply the Gauss–Newton. The details of the Gauss–Newton algorithm for updating are given in Appendix A.

Using , we obtain the phase error matrix as

4.3. Linear Regression-Based Approach

By noting that , the gradient can be calculated via a linear least squares estimator using obtained from (43) [30]. Specifically, we have where

Note that is the unwrapped version of and is an estimate for the initial phase which results from the unwrapping process. From (50), is obtained and can then be used for computing as in (49).

4.4. Differenced-Phase-Based Approach

Subtracting the estimated phase errors of two successive frequency bins, we obtain:

Therefore, can be estimated as [30] where

Note that the four-quadrant arctangent has been used here to handle the phase wrapping. An estimate of the phase error matrix is now obtained as in (49) using in (53).

5. Simulation and Discussion

Numerical simulations are presented in this section to evaluate the performance of the methods described in previous sections.

5.1. Scenario 1: Two Sub-Bands

We consider a synthetic scenario with two sub-bands at carrier frequencies of GHz and GHz, each having a bandwidth of MHz and 64 frequency steps (i.e., ). The range profile is discretised over a grid with a length of m and a grid step of m. We consider a far-field target consisting of six point scatterers which are aligned with the grid. Figure 1 plots the true range profile of the target. We set for the case of existing phase errors. The signal-to-noise ratio is set to 10 dB.

Figure 1

True range profile of synthetic target under consideration.

Figure 2 compares the reconstructed range profiles obtained by the conventional OMP algorithm without and with the presence of phase errors. The OMP is terminated when the signal residual reaches the noise level or after 15 iterations have been carried out. We observe that OMP successfully reconstructs the range profile of the target by correctly identifying the scatterers of the target with accurate range and coefficient estimates when no phase errors exist. However, OMP provides unsatisfactory results in the presence of phase errors, where the reconstructed image is observed as exhibiting many spurious scatterers. Similar observations are obtained for the results obtained by the conventional -norm regularised optimisation solver (without phase error correction), as shown in Figure 3. Here, we set and where is the element of . The -norm regularised optimisation solver is stopped if the relative change in the -norm of the range profile vector falls below or after it reaches 500 iterations. The performance degradation of these conventional sparse reconstruction techniques is not unexpected, since they were not originally developed to cope with dictionary mismatch arising from the presence of the phase errors.

Figure 2

Performance of conventional OMP.

Figure 3

Performance of conventional -norm regularised optimisation solver.

Figure 4 shows the reconstructed range profile obtained by the POMP. POMP constructs candidate dictionaries based on a grid of with a grid step size of . The same stopping criteria of OMP is used for POMP. We observe that POMP produces a range profile which is almost identical to the ground truth, thereby demonstrating the effectiveness of POMP in terms of dealing with the phase errors between different sub-bands.

Figure 4

Performance of POMP in the presence of phase errors.

Figure 5 shows the results obtained by different -norm regularised optimisation solvers with phase error correction, as presented in Section 4. The same parameters and stopping criteria of the conventional -norm regularised optimisation solver as described above are used in the simulations. Although these algorithms exhibit some improvements over the conventional -norm regularisation (i.e., without phase error correction), they provide poorer results compared to that of POMP. Specifically, the peaks of the reconstructed range profiles obtained by these algorithms only appear close to but not exactly at the true scatterer positions. In addition, the magnitudes of the peaks are much smaller than the ground truth values.

Figure 5

Performance of -norm regularised optimisation solver with phase error correction.

The inferior performances of these algorithms can be explained by noting that the -norm regularised optimisation in (31) is a nonconvex problem due to the phase error matrix . Figure 6 plots the objective function of (31) as a function of assuming that is perfectly known. We observe that this objective function has many local maxima and minima, confirming the non-convexity of the -norm regularised optimisation in (31). The reason for this non-convexity is explained in Appendix B. Due to this nonconvexity, the iterative solvers presented in Section 4 are prone to converge to local minima; thus, limiting the effectiveness of this approach.

Figure 6

Illustration of the nonconvexity of the -norm regularised optimisation problem (31). The objective function of (31) is plotted against assuming that is perfectly known.

The performance of the OMP, POMP, and -norm regularised optimisation methods are now compared using the earth mover’s distance (EMD) between the true and reconstructed range profiles. EMD [31] is a metric estimating the distance between two distributions or equivalently the minimal amount of work required to transform one distribution to the other. Figure 7 and Figure 8 show the EMD performance of the OMP, POMP, and -norm regularised optimisation methods, averaged over 100 Monte Carlo runs, versus different levels of SNR (signal-to-noise ratio) and phase error, respectively. It is observed that the POMP method yields the smallest EMD values amongst all algorithms considered. Since a smaller value of EMD corresponds to a higher level of similarity between the true and reconstructed range profiles, this observation indicates that the reconstructed range profile obtained by POMP is closer to the ground-truth range profiles than those obtained from the OMP and -norm regularised optimisation methods. This verifies the performance advantage of the POMP method from a statistical point of view.

Figure 7

Earth mover’s distance (EMD) performance of the OMP, POMP, and -norm regularised optimisation methods versus various of SNRs ().

Figure 8

EMD performance of the OMP, POMP, and -norm regularised optimisation methods versus various levels of phase error ( dB).

5.2. Scenario 2: Four Sub-Bands

We now consider another scenario with four sub-bands at carrier frequencies of GHz, GHz, GHz, and GHz, each having a bandwidth of MHz and 64 frequency steps (i.e., ). The range errors are set to , , and . Other simulation parameters and the true range profile of the target remain unchanged as in the previous simulation example.

Figure 9 compares the reconstructed range profiles obtained by the conventional OMP algorithm and the POMP algorithm presented in Section 3.2. OMP results in an unsatisfactorily reconstructed range profile with many spurious peaks, as expected, because it ignores the phase errors between different sub-bands. In contrast, the POMP is capable of reconstructing the true range profile with a high accuracy thanks to the use of dictionary learning with a pruning process. Note that, given the inferior performance of the -norm regularised optimisation approach compared to POMP, as demonstrated in the previous simulation scenario, this approach is excluded from the comparison here.

Figure 9

Performance comparison between OMP and POMP for Simulation Scenario 2 (with four sub-bands).

6. Conclusions

This paper explores the use of the POMP algorithm and -norm regularisation solvers for the problem of sparsity-driven HRRP with bandwidth stitching in the presence of phase errors. We observe that the -norm regularisation solvers do not provide significant performance improvement over the conventional sparse reconstruction algorithms due to the nonconvexity of -norm regularised optimisation when phase errors exists. In contrast, POMP is observed to be capable of effectively dealing with the phase errors and thus be able to reconstruct the range profile of the target with high accuracy. Simulation results show a significant performance improvement by POMP over OMP and the conventional and refined -norm regularisation. In future work, we propose using experimental data for a more general scenario where the true scatterers constituting the target are located in off-grid positions with respect to the dictionary grid, and the true range errors have off-grid values. We shall also consider the more general case of frequency-dependence of scatterer RCS (radar cross-section). A potential approach for this is to exploit the framework of spectral compressive sensing [32,33].