Total rotational velocity estimation in a multistatic ISAR system

: In this study, the authors report on a novel method to estimate the total rotational velocity of a non-cooperative moving target in a multi-static Inverse Synthetic Aperture Radar (ISAR) system. The method requires a system with at least four receivers spatially distributed such that their line-of-sight vectors are not coplanar. By using least squares optimisation, estimates of the magnitude of the effective rotational velocity ‘seen’ by each receiver and trackfile information on the target trajectory and translational velocity are combined to determine the component due purely to the target's own rotational motion; an estimate for the total rotational velocity can thus be obtained, which can be useful for automatic target recognition. The multistatic nature of the system allows the technique to take into account the translational motion component of the total rotational velocity with some levels of robustness, which is not the case with other approaches. Simulated data results demonstrate the feasibility of the proposed method.


Introduction
ISAR imaging has been studied and developed into imaging systems in the last few decades [1] as it provides a strong support to target classification and recognition. The works in [2][3][4][5][6][7] reported on methods for estimating the cross-range scaling factor for the cross-range dimension and hence the physical size of a target; these methods however have limitations, mainly due to the unknown image projection plane (IPP) and that they generate 2D images for generally 3D targets.
Although ISAR has proven to be a powerful tool, one of its most significant limitations is uncertain image interpretation, which is strongly related to the uncertainties in estimating the IPP. Such a plane is usually not known a priori, as it depends on the target's dynamics during the radar's coherent processing interval (CPI).
In recent years, possible solutions to the 2D limitation have been proposed [8][9][10][11][12]. Early works such as [13] used the principle of interferometry and a radar system of two coherent receivers located parallel to the target rotation vector, the direction of which is assumed known, to estimate scatterer 'relative altitudes', and thus obtained 3D images. For the more general cases of unknown target dynamics, a system with at least three coherent receivers is required for scatterer height estimation as well as estimating the cross-range scaling factor and the tilt angle of the IPP [14][15][16]. Such an estimation of the 3D positions of a scatterer at a certain point in time still leaves open the question of how each scatterer moves over time, and more importantly, how to track them across different radar systems in a multi-static radar network. This paper is motivated by the fact that target's rotational dynamics can be completely described by the target's total rotational velocity vector, which is independent of the viewing radar location. Determining this vector would fully specify the IPP, which is crucial for correct image interpretation. The target's total rotation vector could be used further to enhance tracking and imaging capability of multi-static radar networks. For instance, it would make it possible to determine the target trajectory over time, or the IPP relative to any radar system within a multi-static network. The latter can be obtained by exploiting the radar coordinates with respect to the target and the target total rotation vector.
Some earlier works to estimate the total rotation vector of a target have been reported in the literature. The works in [17,18] proposed a technique that is based on a single interferometric system, with enhanced image resolution through a second-order local polynomial Fourier transform (LPFT) to estimate the total rotation vector. The estimation of the total rotation vector based on a bistatic radar imaging system was reported in [19]. It used a least squares approach to estimate the total rotation vector from spatially distributed measurements of the magnitude of the effective rotation vector. However, the work assumed that the target does not undergo any translational motion. The works in [20,21] used a multi-static radar network to estimate the target translational motions. Such information is then used to enable multi-static radar imaging techniques, for multi-static ISAR autofocus and multistatic radar image formation. In particular, these algorithms have been shown to be effective in evaluating the translational motion of a target. Nevertheless, their performances degrade when the target undergoes complex motions, such as rotations around the three main axes, namely, pitch, roll and yaw.
In this paper, we propose an alternative solution, in which a multi-static ISAR system with at least four receivers is used to estimate the total velocity vector by combining the estimates of the magnitude of the effective rotation vectors of each component system. We also address the situation in which the target translational motion is not negligible.
The organisation of the rest of the paper is as follows: Section 2 is the core of the proposed method that describes for the total rotation vector estimation in a multi-static system. Section 2.2 summarises the method for the estimation of the magnitude of the effective rotation vector, while some preliminary results are described and discussed in Section 3. Finally, Section 4 concludes the paper.

Total rotational velocity estimation
As this work is focused on rotational velocity estimation in a multistatic radar system, making use of measurements of the magnitudes of the effective rotational velocities in the component channels, the signal model and ISAR processing for the component channels [1] is not reproduced here. At this higher level of analysis, the 'system geometry' is represented by the various radar line-of-sight (LOS) directions of the observing radars and associated 'cylindrical surfaces' discussed in this section.

Components of the total rotational velocity
A target's total rotational velocity (or simply 'the rotation vector') Ω Tot (l) relative to a receiver l may be considered as consisting of two components [1]: where Ω Tr (l) is the component due to translational motion, while Ω Rot is due to rotational motion with respect to the target's own inertial frame. The Ω Rot component is common to all radar receivers, whereas Ω Tr (l) is different for different receiver locations, given by where R (l) is the range vector from receiver l to the target, ν is the target velocity vector and i LOS (l) is the LOS vector. Hence, Ω Tr (l) can be readily computed if the target translational motion is known. Methods such as those described in [20,21] can be used to estimate the target trajectory and velocity. In the following analysis, we shall assume that Ω Tr (l) are measurements for all receivers, with some unbiased statistical errors, after which the only missing component is the receiver independent component Ω Rot , which we aim to estimate.
For each receiver channel, ISAR processing provides a means to estimate the effective component, Ω Eff (l) , of Ω Tot (l) , which is its projection onto the normal of the IPP. Its magnitude can be expressed as follows [2,7]: Since the IPP is not known a priori, an estimate of Ω Eff (l) from ISAR processing constrains the possible solutions for Ω Tot (l) to the curved surface of a cylinder in frequency space, with its axis being the LOS from receiver l and radius equal to Ω Eff (l) , as illustrated in Fig. 1.
Based on (1), the possible solutions for are represented by the same cylindrical surface shifted by the vector −Ω Tr (l) , which we shall refer to as the 'shifted cylinder (or cylindrical surface) l', and are illustrated in Fig. 2. This represents an extra complexity of the problem, in that the effective rotational velocity useful for target imaging may come from two types of target motion, namely, translational (bulk) motion and its own rotational motion. The estimations of Ω Eff (l) and Ω Tr (l) from multiple distributed receivers result in a multitude of such shifted cylindrical surfaces, as illustrated in Fig. 3 for an example of three receivers. Furthermore, while the radar LOS directions are concurrent at the target's coordinate origin, the axes of the shifted cylinders are generally not concurrent, but the independence of Ω Rot on receiver directions means that the shifted cylindrical surfaces must intersect at a single point -the end point of the Ω Rot vector, which is the solution we seek. It may be conjectured that at least four of such shifted cylinders are required to uniquely determine the total rotation vector, as long as their LOS vectors are not coplanar. A discussion of this conjecture is included in Appendix 7.

Estimation of the magnitude of the effective rotation vector
For the simulations in the following sections, the method in [7] is adopted to obtain estimates of the modulus of the effective rotation vector for each radar component in a multi-static system. This section provides a summary of the algorithm.

Chirp rate estimation:
In a 2D ISAR image, the chirp rate μ k for scatterer k is proportional to its range coordinate x 2k [2,7]: Since the chirp rate is independent of the cross-range coordinate, we can estimate a single chirp rate at any range bin, even if it contains multiple scatterers. This allows us to avoid having to first estimate each scatterer's point spread function in schemes like CLEAN [2,22]. The chirp rate can be estimated by maximising the kurtosis with respect to the parameter μ: where μ k is the estimated chirp rate at range bin k and X k (ν, μ) is its second-order LPFT (2-LPFT) along slow-time [2]. The optimisation in (5) is performed on all range bins with identified scatterers. The estimated chirp rates should be proportional to the range bin from (4). We use a robust linear regression technique [23] to fit the set of estimated chirp rates, and that produces an estimate of the slope in (4), from which we can estimate Ω Eff .

Summary of algorithm for Ω Eff estimation:
The algorithm for the estimation of Ω Eff is summarised as follows [7]: i. Obtain a conventional complex range Doppler ISAR image S RC (τ, ν). ii. Perform 2D peak detection to identify range bins that contain scatterers. iii. Estimate the chirp rate μ k for each identified range bin by solving the optimisation problem in (5). iv. Perform robust linear regression on the estimated chirp rates across the range and use the slope obtained to estimate the magnitude of the effective rotation vector Ω Eff .
Note that after the initial peak detection from the ISAR image, processing for the algorithm takes place along a single dimension, thus keeping its computational requirement to a minimum. The estimates Ω Eff can be used to scale the Doppler dimension of an ISAR image to true cross range. It can also be used to improve the cross-range resolution of the ISAR image through the use of the 2-LPFT. More importantly, the estimated Ω Eff from at least four receivers will allow us to obtain the estimation of the total rotation vector of the target, as discussed in Section 2.1.

Least squares based estimation of Ω Rot
An exact solution for Ω Rot is intractable, and in practical scenarios, the estimations for Ω Eff (l) are subject to errors. Analogous to multilateration techniques [20], an optimal solution in the least squares sense can be applied. Substituting (1) the last step is due to the fact that Ω Tr (l) is orthogonal to i LOS (l) , as implied in (2). As mentioned above, the vectors Ω Tr (l) can be estimated from the trajectory measurements of the target, while Ω Eff (l) can be estimated using one of the methods from [2,7], which we shall denote as Ω Tr (l) and Ω Eff (l) , respectively. Section 2.2 describes the method to be used for this paper. As far as the current optimisation is concerned, both Ω Tr (l) and Ω Eff (l) are treated as measured parameters with statistically unbiased errors. Hence, the only unknown parameter in (6) to be estimated is Ω Rot .
The errors in Ω Eff (l) for each receiver can be written as follows: Given the estimated translational motion parameters, the optimal solution, in the least squares sense, for Ω Rot in (8), is where L is the total number of receivers in the multi-static ISAR system. To solve the optimisation problem (9), the common Nelder-Mead algorithm [24] is used in our work. This is a simplex-based hill climbing technique that uses heuristic search to avoid direct gradient computations. It is known that the algorithm is not guaranteed to converge to a local minimum; we therefore initialised the algorithm on multiple starting points, and the solution is based on the minimised result.
If Ω Tr (l) and Ω Eff (l) are indeed unbiased estimates as assumed above, then the estimated Ω Rot can also expected to be unbiased; however, this may not be the case in real scenarios with complex target motions, and error distributions in the measured parameters may lead to unknown distributions in the estimated Ω Rot and hence in the estimated total rotational velocities Ω Tot (l) ; these practically important questions are beyond the scope of the current paper and will be addressed in a future publication.

Volume of admissible values of Ω Rot
The number of receivers, L ≥ 4, is a system design parameter. An important consideration is the diversity of LOS directions from the receivers. The cylindrical surfaces intersect at a single point only if Ω Eff (l) as the radii of the cylinders and Ω Tr (l) as estimated from target trackfiles are exact. When uncertainties (i.e. 'error bars') exist for these parameters, a cylindrical surface is replaced by a cylindrical volume around it, and the single point of intersection is replaced by an intersected volume in space, denoted as V, which defines the 'volume of admissible values' of Ω Rot for the optimisation in (9).
For specified levels of uncertainty in Ω Eff (l) and Ω Tr (l) , the robustness of the estimation is also limited by the receiver LOS directions. A favourable set of receiver LOS directions reduces the volume V and thus the uncertainty in the estimates. This condition is analogous to the concept of geometric dilution of precision in global positional system analysis and is hypothesised to reduce to a condition on a measure of 'mutual independence' of the receiver LOS vectors, which again is re-scheduled to a future publication. Fig. 4 shows a block diagram of the proposed algorithm. The translational motion of the target is estimated from range tracks at all the receivers [20]; the magnitudes of the effective rotational velocities seen by the receivers are independently estimated; these estimates are then used to solve for the Ω Rot component, which is due only to target rotation, using a least squares optimiser, as discussed in Section 2.3. The total rotational velocity, for each receiver, can then be computed from (1).

Simulation results
In this section, the efficiency of the proposed algorithm is demonstrated using simulations of varying degrees of complexity, all carried out using the ISARLAB software [25]. Additive Gaussian noise is added to simulate non-idealities in the signals, and Monte Carlo simulations, with 50-100 noise realisations, are run in each case. The signal-to-noise (SNR) value determines the strength of the injected Gaussian white noise; the noise variance is computed from the measured signal power as follows: where I(i, j) 2 is the image intensity of the jth range profile at range bin i, N s is the number of range bins and N p is the total number of range profiles in the ISAR image. Three examples are shown. The simulation and system parameters for all three examples are summarised in Table 1. In the first two, we generated 400 random radar locations on the surface of an imaginary sphere of radius 10 km with the target model located at the centre of the sphere. While physically unrealistic, this set-up allows us to explore the resource of spatial diversity in its purest form. Each estimation of the rotational velocity Ω Rot will utilise Ω Eff estimates from eight of these radar locations; this gives 50 different sets of estimations. Although the proposed method only requires four receivers, we choose to use eight receivers in our simulations to minimise the effects of the dependency on the accuracy of the estimated Ω Eff . The target itself is modelled as an ensemble of 101 point scatterers located within a spherical volume, as shown in Fig. 5. The spherically quasi-symmetric target model ensures that the received data at an arbitrary radar LOS will observe a similar number of scatterers for the Ω Eff (l) estimation algorithm, whose accuracy would otherwise depend on the target's reflectivity distribution. This set-up would minimise the variations of the Ω Eff (l) estimation errors for the purpose of evaluating the LSE algorithm. In Example 1, the translational component Ω Tr of the total rotational velocity is set to zero, while Ω Rot is set to several different values; the results are described in Section 3.1. In Example 2, non-zero translational components are introduced to the total rotation vector; the results are presented in Section 3.2. In Example 3, the radar and target locations are chosen as shown in Fig. 6b to provide a more realistic scenario than that in Section 3. The target model in Fig. 6a is used for the simulations. It consists of 46 point scatterers arranged to model an airplane target, and is chosen to be the same as the target used in [17] to facilitate comparison. The reported 3D InISAR technique uses interferometric processing to perform 3D reconstructions from ISAR images obtained from at least three radar receivers, offset from each other along two orthogonal baselines. The rotational velocity Ω Rot can be estimated by including chirp rate processing similar to the technique presented here. While the InISAR scheme requires (one) fewer radar receivers than the proposed approach, it requires a quasi-monostatic set-up of the receiving antennas and total coherence between the different receivers, which may significantly increase implementation challenges. In this example, eight radars are used in the proposed estimation algorithm. (Refer to Fig. 9. The motivation for using more than a minimum of four radars is also given in Section 3.1.) Their locations are not chosen to particularly favour the multi-static technique; the LOS from the four receivers does not offer widely diverse perspectives on the target. One of these radars is co-located with the 3D InISAR radar.

Example 1: Ω Rot estimation without translational motion
The simulations are carried out for several different true values of Ω Rot . The estimation results are shown in Fig. 7, in the form of relative error distributions for the three components of the vector. In all cases, the SNR in the signal domain is 30 dB. Fig. 7a shows the estimation results for the case of true Ω Rot = (2.9, 2.9, 5.0)°/s. It can be seen from these results that the errors are <15%, with the majority below 10%. When true Ω Rot is changed to Ω Rot = (1.4, 1.4, 5.0)°/s, the relative errors increase for the smaller components. This is expected, since shorter angular extents (i.e. slower rotations for the same CPI length) are known to be produce less accurate results for Ω Eff estimation. Increasing the first component of Ω Rot to 2.9°/s, i.e. true Ω Rot = (2.9, 1.4, 5.0)°/s, produces an improvement in the estimation error, as shown in Fig. 7c. Fig. 8 shows the results when the theoretically minimum of four radar receivers is used. That is, 100 randomly generated subsets of four receivers are used for the rotation vector estimation. Not surprisingly, it was found that the estimation errors can be significantly larger when four receivers are used instead of eight. Interestingly, there is also a stronger concentration of accurate estimations (say, <5% error) for the four receiver cases, thus suggesting that fewer but higher quality receiver locations can lead to excellent overall estimation performance.
When analysing the detailed results, it was found that two factors are responsible for the cases with large (which we define as >10%) errors: i. Large errors in one or more of the estimated values of Ω Eff from the multi-static radars; ii. Near coplanar LOS vectors for a significant subset of the randomly located radars.
The first issue is sometimes unavoidable due to the geometric relationship between the target and any single radar receiver. It would therefore be advantageous to use fewer receivers to avoid the occasional failed Ω Eff to pollute the overall rotation vector estimate. The second issue highlights the importance of choosing favourable radar locations for the proposed technique. Examples of favourable and unfavourable geometries are shown in Fig. 9. The receiver locations relative to the target cannot be controlled in general, although it may be possible to distribute as many receivers as far apart as allowable to minimise the likelihood of the second issue arising. The comparison between Figs. 7 and 8 bears out these points.

Example 2: Ω Rot estimation with translational motion
In this Example, the simulated scenarios and parameters are the same as in the previous example, except the target now travels in a linear path at 20 or 50 m/s. That is, the velocity vectors are (20, 0, 1519°/s, respectively. It can thus be seen that the incorporation of translational motion does not materially affect the estimation performance.

Example 3: Ω Rot estimation comparison with another technique
In this example, the radar and target locations are chosen as shown in Fig. 6b. The target model in Fig. 6a is used for the simulations. It consists of 46 point scatterers arranged to model an airplane target. The location of the central radar for the interferometric setup is 'radar 4', as shown in Fig. 6b. The simulation parameters for this example are as shown in Table 1.
The Ω Rot vector estimation performance for the airplane target model is shown in Fig. 12. The estimation accuracy fluctuates greatly for the different components. Within each component, there is no clear pattern of performance. This can be explained by the optimisation of the objective function in (9) being done over all of the components concurrently. Thus, the error of the total rotation vector magnitude as shown in Fig. 13 follows a more reasonable pattern, where the accuracy is ∼5% for all SNR >−5 dB.
The simulation results for the 3D InISAR approach are shown as an orange curve in Fig. 13. It is observed that performance of the components of Ω Rot does not appear to follow any clear trends as SNR is varied. Paradoxically, Ω 2 accuracy appears to improve with lower SNR. When the results are combined, the error of the magnitude Ω Rot is mostly monotonically decreasing, as expected.

Further discussion and conclusion
This paper is undoubtedly a preliminary work on the topic of total rotational velocity estimation in a multi-static ISAR system. The full complexity of real target dynamics, especially with airborne manuevres, has not fully been analysed; the effects of possible biases and unknown distributions in the estimations for Ω Tr from trackfiles and Ω Eff from ISAR processing have not been accounted for.
Nevertheless, the work has described the principle and demonstrated the effectiveness of a novel method, at least in simple cases of practical interest, to estimate the total rotational velocity of a target in a multi-static ISAR system, by making use of the spatially distributed radar measurements. The total rotational velocity can be estimated if there are at least four receivers such that their LOS vectors are not all coplanar, along with the estimated magnitudes of the effective rotational velocity seen by each receiver. The results obtained from simulations demonstrated the effectiveness of the proposed method, with estimation errors of 10% used as a benchmark for our evaluations. As shown in this paper, the proposed method relies on favourable geometries of the multi-static system, such that the component radar LOS directions are not all coplanar (or colinear), and the quality of the estimated magnitude of the effective rotation vector of each receiver.
A future work will attempt to address the shortcomings mentioned above, with real multi-static radar data in real scenarios.

Acknowledgments
We would like to thank the Defence Science and Technology for the Fellowship program and the PhD Scholarship funding under which this research work was carried out. Special thanks to Fraunhofer FHR for hosting An Phan where most of this joint collaborative work was undertaken, and to Dr. Josef Worms of FHR for a very stimulating and helpful discussion that has significantly improved the paper.  Here, we give a graphical justification for the requirement of at least four receivers at different LOS directions for a unique determination of the total rotation vector, given measurements Ω Eff (l) of the magnitudes of the respective effective rotation rates. Furthermore, they must be non-coplanar.

Case of no translation motion
First, let us consider the case when no translational motion is present, i.e. Ω Tr (l) = 0 for all receivers; in such a case, the cylindrical axes do intersect at the origin of the local target frame. With respect to each receiver, the true Ω Rot has its starting point at the origin and its tip point anywhere on the corresponding surface; hence, Ω Rot can be found at the intersection point of a number of such surfaces. Without loss of generality, we sort the Ω Eff (l) valuesthe radii of the cylindrical surfaces -in the order of decreasing magnitudes, so that all intersects lie on the first cylindrical surface. Consider the intersect of the first two cylinders: their intersection is a set of two 'curved ellipses', wrapped on opposite sides of the first cylindrical surface, one of which is shown as the solid line in Fig. 14 or Fig. 15. (Note that such a 'curved ellipse' approaches a true planar ellipse only in the limit when the radius of the first cylinder approaches infinity.) The curved ellipse has its 'major axis' as a straight line on the first cylindrical surface parallel to its axis. (The major axis is also the intersection of the plane formed by the two LOS vectors and the first cylindrical surface.) The intersect of the third cylindrical surface with the first is another curved ellipse on the first cylinder, which shares the same major axis only if the third LOS vector is coplanar with the other two, as shown in Fig. 15; in this case, the two points of intersection of the two ellipses are symmetric about the major axis, as indicated by the two red dots. Otherwise, it is offset from the first major axis, as shown in Fig. 14, and the intersection points (blue dot and lower red dot) are generally asymmetric.
When combining with the fourth cylindrical surface, one of the two intersection points cannot be eliminated if its axis (LOS vector) is also coplanar with the other three axes, because of the constraint of symmetry across the major axis, as illustrated in Fig. 15. In all other cases, only one common point (the blue dot) for the curved ellipse intersection remains. The required condition on the four LOS vectors is thus: they are not all coplanar.
As mentioned above, there are two ellipses on either side of largest cylindrical defining two equal magnitude but opposite signed Ω Rot vectors. By using other considerations, it is possible to narrow down to single vector.

Case of non-zero translational motion
In this case, the cylinders representing solutions for Ω Rot are the shifted cylinders, with their axes remaining parallel to the respective radar LOS directions. The magnitudes and directions of the shifts are determined by (2), which implies that, unless the target velocity ν and all radar LOS vectors are coplanar, the shifts will generally be in different directions and with different magnitudes; the cylinder axes become non-concurrent, and the curved ellipses on the largest cylinder become asymmetric.
In this paper, we conjecture that under such conditions the symmetry leading to the degeneracy seen in Fig. 15 is broken, and a unique intersection point of the four shifted cylindrical surfaces representing a solution for Ω Rot exists. A more rigorous analysis is currently planned for a future publication.