Accurate three-dimensional tracking method in bistatic forward scatter radar

Abstract

Accurate three-dimensional (3D) tracking in bistatic forward scatter radar (BFSR) isa challenging problem because of absent range resolution and poor measurements. Inthis article, an accurate 3D tracking method of BFSR is proposed. Aiming to overcomethe filter divergence caused by large initial state estimation error, firstly, anaccurate initial state estimation approach is presented based on analytic derivationand Levenberg–Marquardt algorithm, which has the potential to improve theaccuracy of initial state estimation. Furthermore, in order to reduce the computationcost of filtering process and speed up the filtering convergence rate, the accurateinitial state estimation and extended Kalman filter algorithm in BFSR are combined toachieve a precise target 3D tracking. Finally, the proposed accurate tracking methodis verified through comparative analysis of the simulation results.

Introduction

Because of the gain in bistatic RCS (up to 20–40 dB relative to monostatic RCS) inforward scatter region, bistatic forward scatter radar (BFSR) can effectively detect andtrack the target with low-speed or small RCS (including the stealth target) [1, 2]. Its target tracking technology attracted more and more attention in recentyears. An existing challenge of BFSR is parameter estimation of aerial [2]–[9] and ground [10]–[17] targets. This article focuses on parameter estimation of aerial targettrajectory.

Systems using continuous quasi-harmonic probe signal have been proved to be mostpromising amongst various structures of BFSR [3]–[5]. In this case, measurements of echo Doppler shift f _d, angle of arrival θ, and elevation angle β areusually used to obtain object parameter estimations (i.e., target position and velocitydenoted as x _k , y _k , h, V _x , V _y ). Due to the nonlinear relations between f _d, θ, β, and trajectory parameters,targetparameter estimation in BFSR is a nonlinear optimization problem. To obtain five unknownvariables (i.e.,x _k , y _k , h, V _x , and V _y ) using three known variables (i.e.,f _d, θ, β), equations composed by groups ofobservations are needed. Hence, object parameter estimation in BFSR is actually solvingover-determined nonlinear equations. The optimization is achieved using nonlinear leastsquares algorithms, of which Gauss–Newton method has widely been used.

In previous research [3]–[7], Gauss–Newton method is adopted to solve the nonlinear equations. Onthis basis, the classical methods of target tracking using two-dimensional (2D) andthree-dimensional (3D) BFSR are proposed, respectively, in [4, 6]. The major problem is that Gauss–Newton method is not a“real-time” algorithm, since we have to collect n measurements toimplement it. To achieve higher data update rate and reduce computation in parameterestimation, extended Kalman filter (EKF) algorithm [18] is used. However, initial filtering value of EKF can only be obtained bysolving over-determined nonlinear equations, which can never be avoided. In classicalmethod, due to the special geometry of FSR system, the initial value estimation byGauss–Newton method requires high-precision Doppler shift measurement and anglemeasurement, which are unavailable due to the limitation of antenna size andaccumulation time. Therefore, Gauss–Newton method is easy to cause a large errorof initial value estimation and filtering divergence.

In this article, based on analytical derivation and Levenberg–Marquardt (LM) [19, 20] algorithm, a new approach of initial state estimation is presented, which canimprove the accuracy of initial state estimation without high-precise target parametersmeasurement. Hereafter, the accurate initial state estimation and EKF algorithm arecombined to achieve fast convergence and high accuracy. Simulation results prove thataccurate target trajectory parameters estimation can be obtained by this method beforethe target crossing the baseline.

The remainder of this article is organized as follows: Section 1 briefly describes thedevelopment of forward scatter radar, pointing out the existing problems of air targetparameter estimation and focus of this article; Section 2 defines the geometry of 3Dforward scatter radar, the target motion model and system observation model; Section 3introduces the proposed tracking method for forward scatter radar in details; Section 4verifies the validity of accurate tracking method through comparative analysis of thesimulation results; Conclusions are drawn in Section 5.

System modeling

Geometry and system implementation of 3D BFSR

Geometry of 3D BFSR is shown in Figure 1. x,y, z are Cartesian coordinates. Re, Tr, and Tg denote thepositions of receiver, transmitter, and target, respectively. ψ isflight-path angle (the trajectory inclination angle towards the baseline),θ is the target azimuth angle, β is the targetelevation angle, h is the target altitude, b is the distancebetween transmitter and receiver (i.e., the baseline length), AB is thetarget trajectory, CD is the projection of target trajectory in horizontalplane.

Figure 1

Geometry of 3D-BFSR.

System settings meet the following principles:

• 1) In the transmitter position, there is a wide beam antenna illuminating the whole forward scatter region with a quasi-harmonic signal.

• 2) In the receiver position, there is a multi-beam antenna to illuminate the forward scatter region.

• 3) There is no block in line of sight between the transmitting and receiving antenna.

Target model

Assuming the target crosses the baseline near the midpoint in a constant speed withlinear trajectory, target state equation can be written as

$$X\left(k+1\right)=\mathit{\Phi }X\left(k\right)+\mathit{Gv}\left(k\right)$$

(1)

Assuming the target flies in a constant altitude (i.e., the longitudinal velocity iszero), target’s state vector can be denoted by X(k) =[x _k , y _k , hV _x , V _y ], where T means sampling interval, v(k) is Gaussianwhite noise process with zero mean. x _k , y _k , h, V _x , V _y are the values of target Cartesian coordinates at k th discrete timeinstant and their derivatives. State transition matrix Φ and noise distributionmatrix G are, respectively, written as

$$\Phi =\left[\begin{array}{lllll}1\hfill & 0\hfill & 0\hfill & T\hfill & 0\hfill \\ 0\hfill & 1\hfill & 0\hfill & 0\hfill & T\hfill \\ 0\hfill & 0\hfill & 1\hfill & 0\hfill & 0\hfill \\ 0\hfill & 0\hfill & 0\hfill & 1\hfill & 0\hfill \\ 0\hfill & 0\hfill & 0\hfill & 0\hfill & 1\hfill \end{array}\right]\begin{array}{}\end{array}G=\left[\begin{array}{lll}\frac{T^2}{2}\hfill & 0\hfill & 0\hfill \\ 0\hfill & \frac{T^2}{2}\hfill & 0\hfill \\ 0\hfill & 0\hfill & \frac{T^2}{2}\hfill \\ T\hfill & 0\hfill & 0\hfill \\ 0\hfill & T\hfill & 0\hfill \end{array}\right]$$

(2)

Observation model

Assuming observation vectors measured from first to k th instant are

$${\left[\left(f_{d1},{\theta }_1,{\beta }_1\right),\dots ,\left(f_{\mathit{dk}},{\theta }_k,{\beta }_k\right)\right]}^T$$

where f _dk , θ, β _k are observations of the Doppler frequency shift, echo azimuth angle, andelevation angle.

Measurement equation can be written as

$$\overrightarrow{Z_k}=h\left(\overrightarrow{x_k}\right)+\Delta \overrightarrow{Z_k}$$

(3)

where

$$\begin{array}{ll}h\left(\overrightarrow{x_k}\right)& =\left[\begin{array}{l}{f}_{\mathit{dk}}\\ {\theta }_k\\ {\beta }_k\end{array}\right]\\ =\left[\begin{array}{l}-\frac{1}{\lambda }\left[\frac{x_k{V}_x+y_k{V}_y}{\sqrt{x_k^2+y_k^2+h^2}}+\frac{y_k{V}_y-\left(b-x_k\right)V_x}{\sqrt{{\left(b-x_k\right)}^2+y_k^2+h^2}}\right]\\ arctan\left(\frac{y_k}{x_k}\right)\\ arctan\frac{h}{\sqrt{x_k^2+y_k^2}}\end{array}\right]\end{array}$$

(4)

$${\overrightarrow{Z}}_k$$

is the observation vector of k th instant, which isdenoted by (f _dk , θ, β _k ). ${\overrightarrow{X}}_k$ is state vector of k th instant, denoted by [x _k , y _k , h _k , V _x , V _y ]. The measurement noise denoted by $\Delta {\overrightarrow{Z}}_k$ is Gaussian white noise with zero mean,whose noise variances areσ _f ², σ _θ ², σ _β ², respectively.

Accurate 3D tracking in forward scatter radar

Because of the special geometry of forward scatter radar, parameters (i.e., Dopplershift, azimuth angle, elevation angle) cannot directly be measured when the targetcrosses the baseline, resulting in a large estimation error and seriously decreasing thetracking precision. In addition, the detection area of forward scatter radar is a verynarrow region near the baseline where bistatic angle ranges from 135° to 180°,indicating that target tracking is effectively implemented during a very short period oftime. Therefore, a stable tracking trajectory should be obtained before the targetcrosses the baseline, which requires fast convergence of filtering algorithm. However,the large error of initial state estimation will decrease the convergence speed, evenlead to filtering divergence. In summary, the initial state estimation is significantfor target tracking in forward scatter radar, which extremely affects the filteringaccuracy and convergence speed.

The classical method uses Gauss–Newton iteration to estimate the initial state oftarget. Gauss–Newton iteration is prone to singular matrix causing inaccuratecomputation results. And if the first value of iteration deviates far away from its truevalue, the results of iteration easily fall into local minimum and then cause largeerror of initial state estimation. To overcome these problems, based on analyticalderivation and LM algorithm, a novel method is presented to obtain high-precise initialstate estimation.

The basic principles of this new initial state estimation method are as follows.

According to the system observation equations, the observation vector (f _dn , θ _n , β _n ) measured in n th discrete time instant could be written as

$$f_{\mathit{dn}}=-\frac{1}{\lambda }\left[\frac{x_n{V}_x+y_n{V}_y}{\sqrt{x_n^2+y_n^2+h^2}}+\frac{y_n{V}_y-\left(b-x_n\right)V_x}{\sqrt{{\left(b-x_n\right)}^2+y_n^2+h^2}}\right]$$

(5)

$${\theta }_n=\arctan \frac{y_n}{x_n}$$

(6)

$${\beta }_n=\arctan \frac{h}{\sqrt{x_n^2+y_n^2}}$$

(7)

Based on the observation vectors (f _{d 1}, θ ₁, β ₁), (f _{d 2}, θ ₂, β ₂),…(f _dn , θ _n , β _n ) of previous time instants, the target state vector (x _n , y _n , V _x , V _y , h) in n th time instant can be obtained. The process is asfollows:

From Equation (6), then

$$\frac{y_1}{x_1}=\tan {\theta }_1,\dots ,\frac{y_n}{x_n}=\tan {\theta }_n$$

(8)

The target motion model shows the relationship between y ₁ and y _n :

$$\left\{\begin{array}{l}{y}_1=y_n-\left(n-1\right)T{V}_y\hfill \\ x_1=x_n-\left(n-1\right)T{V}_x\hfill \end{array}\right.$$

(9)

Substituting Equation (9) into Equation(8), we have

$$\frac{y_n-\left(k-1\right)T{V}_y}{x_n-\left(k-1\right)T{V}_x}=\tan {\theta }_1$$

(10)

Based on Equations (8) and(10), then

$$\left(tan{\theta }_n-tan{\theta }_1\right)x_n-\left(n-1\right)T{V}_y+tan{\theta }_1\left(n-1\right)T{V}_x=0$$

(11)

Similarly, we have

$$\left(tan{\theta }_n-tan{\theta }_2\right)x_n-\left(n-2\right)T{V}_y+tan{\theta }_2\left(n-2\right)T{V}_x=0$$

(12)

And the following simultaneous equations

$$\left\{\begin{array}{l}\left(tan{\theta }_n-tan{\theta }_1\right)x_n-\left(n-1\right)T{V}_y+tan{\theta }_1\left(n-1\right)T{V}_x=0\\ \left(tan{\theta }_n-tan{\theta }_2\right)x_n-\left(n-2\right)T{V}_y+tan{\theta }_2\left(n-2\right)T{V}_x=0\end{array}\right.$$

(13)

Let the coefficients of x _n , V _y , V _x be (a ₁, a ₂, a ₃), (b ₁, b ₂, b ₃), that is

$$\left\{\begin{array}{ll}{a}_1=\left(tan{\theta }_n-tan{\theta }_1\right)\hfill & b_1=\left(tan{\theta }_n-tan{\theta }_2\right)\hfill \\ a_2=-\left(n-1\right)T\hfill & b_2=-\left(n-2\right)T\hfill \\ a_3=tan{\theta }_1\left(n-1\right)T\hfill & b_3=tan{\theta }_2\left(n-2\right)T\hfill \end{array}\right.$$

(14)

Then Equation (13) can be written as

$$\left[\begin{array}{ll}{a}_1\hfill & a_3\hfill \\ b_1\hfill & b_3\hfill \end{array}\right]\left[\begin{array}{l}{x}_k\hfill \\ V_x\hfill \end{array}\right]=\left[\begin{array}{l}-a_2\hfill \\ -b_2\hfill \end{array}\right]V_y$$

(15)

Equations above have solutions as follows

$$x_n=\frac{a_3{b}_2-a_2{b}_3}{a_1{b}_3-a_3{b}_1}\cdot V_y=d_{n1}{V}_y$$

(16)

$$V_x=\frac{a_2{b}_1-a_1{b}_2}{a_1{b}_3-a_3{b}_1}\cdot V_y=d_{n2}{V}_y$$

(17)

$$y_n=tan{\theta }_k\cdot d_{n1}{V}_y=d_{n3}{V}_y$$

(18)

According to Equation (7), then

$$h=tan{\beta }_n\cdot \sqrt{x_n{}^2+y_n{}^2}$$

(19)

Substituting Equations (16) and (18) into Equation (19), the relationship betweenh and V _y could be written as

$$h=tan{\beta }_n\sqrt{d_{n1}^2+d_{n3}^2}\left|V_y\right|=d_{n4}\left|V_y\right|$$

(20)

Substituting Equations (16), (17), (18), (20) into Equation (5), we can obtain thenonlinear equation as follows:

$$\begin{array}{ll}{f}_{\mathit{dn}}& +\frac{1}{\lambda }\left[\frac{\left(d_{n1}{d}_{n2}+d_{n3}\right)V_y^2}{\sqrt{\left(d_{n1}^2+d_{n3}^2+d_{n4}^2\right)V_y^2}}\right.\\ \left.+\frac{d_{n3}{V}_y^2-\left(b-d_{n1}{V}_y\right)d_{n2}{V}_y}{\sqrt{{\left(b-d_{n1}{V}_y\right)}^2+d_{n3}^2{V}_y^2+d_{n4}^2{V}_y^2}}\right]=0\end{array}$$

(21)

Assuming that N set of observations are used in initial state estimation, taken = (N/2 + 1):N, we can get the following nonlinearequations

$$\left\{\begin{array}{l}{f}_{d2}+\frac{1}{\lambda }\left[\frac{\left(d_{21}{d}_{22}+d_{23}\right)V_y^2}{\sqrt{\left(d_{21}^2+d_{23}^2+d_{24}^2\right)V_y^2}}+\frac{d_{23}{V}_y^2-\left(b-d_{21}{V}_y\right)d_{22}{V}_y}{\sqrt{{\left(b-d_{21}{V}_y\right)}^2+d_{23}^2{V}_y^2+d_{24}^2{V}_y^2}}\right]=0\hfill \\ \begin{array}{l}{f}_{d3}+\frac{1}{\lambda }\left[\frac{\left(d_{31}{d}_{32}+d_{33}\right)V_y^2}{\sqrt{\left(d_{31}^2+d_{33}^2+d_{34}^2\right)V_y^2}}+\frac{d_{33}{V}_y^2-\left(b-d_{31}{V}_y\right)d_{32}{V}_y}{\sqrt{{\left(b-d_{31}{V}_y\right)}^2+d_{33}^2{V}_y^2+d_{34}^2{V}_y^2}}\right]=0\\ ⋮\end{array}\hfill \\ f_{\mathit{dn}}+\frac{1}{\lambda }\left[\frac{\left(d_{n1}{d}_{n2}+d_{n3}\right)V_y^2}{\sqrt{\left(d_{n1}^2+d_{n3}^2+d_{n4}^2\right)V_y^2}}+\frac{d_{n3}{V}_y^2-\left(b-d_{n1}{V}_y\right)d_{n2}{V}_y}{\sqrt{{\left(b-d_{n1}{V}_y\right)}^2+d_{n3}^2{V}_y^2+d_{n4}^2{V}_y^2}}\right]=0\hfill \end{array}\right.$$

(22)

Solving Equation (22) by the LM algorithm can achieve the optimal solution of V _y , and x _n , V _x , y _n , h, which can also be obtained by Equation (16), (17), (18), and (20).Thus, the numerical solution of the target state (x _n , y _n , h, V _x , V _y ) in n th instant is acquired.

Because measurement noise is considerably large relative to the measurement data ofazimuth angle, it is necessary to smooth primary measurements through polynomial fittingbefore the initial state estimation. Polynomial coefficients can be obtained by leastsquare method, which can reduce the initial state estimation error by an order ofmagnitude.

In the stage of follow-up tracking, EKF algorithm is adopted. As mentioned earlier, thesystem state transition equation and measurement equation are as follows:

$$X\left(k+1\right)=\mathit{\Phi }K\left(k\right)+\mathit{Gv}\left(k\right)$$

(23)

$$Z\left(k\right)=h\left(k,X\left(k\right)\right)+W\left(k\right)$$

(24)

where h(k, X(k)) denotes the measurement vector,measurement noise denoted by W(k) is Gaussian white noise process withzero mean. Assuming the process noise covariance and measurement noise covariance attime k are, respectively, expressed as Q(k) andR(k), then the steps of filtering are as follows:

State prediction:

$$\widehat{X}\left(k+1\left|k\right.\right)=\mathit{\Phi }X\left(k\right)$$

(25)

Error covariance prediction:

$$P\left(k+1\left|k\right.\right)=\mathit{\Phi }P\left(k\left|k\right.\right){\Phi }^T+Q\left(k\right)$$

(26)

Measurements prediction:

$$\widehat{Z}\left(k+1\left|k\right.\right)=h\left[\widehat{X}\left(k+1\left|k\right.\right)\right]$$

(27)

Calculate residual covariance:

$$S\left(k+1\right)=h_X\left(k+1\right)P\left(k+1\left|k\right.\right)h_X^T\left(k+1\right)+R\left(k+1\right)$$

(28)

Calculate Kalman gain:

$$K\left(k+1\right)=P\left(k+1\left|k\right.\right)h_X^T\left(k+1\right)S^{-1}\left(k+1\right)$$

(29)

Update state:

$$\begin{array}{ll}\widehat{X}\left(k+1\left|k+1\right.\right)=& \widehat{X}\left(k+1\left|k\right.\right)+K\left(k+1\right)\\ \times \left\{Z\left(k+1\right)-h\left[\widehat{X}\left(k+1\left|k\right.\right)\right]\right\}\end{array}$$

(30)

Update error covariance:

$$\begin{array}{l}P\left(k+1\left|k+1\right.\right)=\left[I-K\left(k+1\right)h_X\left(k+1\right)\right]P\left(k+1\left|k\right.\right)\hfill \\ \cdot {\left[I-K\left(k+1\right)h_X\left(k+1\right)\right]}^T\hfill \\ -K\left(k+1\right)R\left(k+1\right)K^T\left(k+1\right)\hfill \end{array}$$

(31)

h _X (k) in Equation (28) is the Jacobian matrix of h(k,X(k)), which is

$$\begin{array}{ll}h\left(k,X\left(k\right)\right)& =\left[\begin{array}{l}{f}_{\mathit{Dk}}\left(X\right)\hfill \\ {\alpha }_k\left(X\right)\hfill \\ {\beta }_k\left(X\right)\hfill \end{array}\right]\\ =\left[\begin{array}{l}-\frac{1}{\lambda }\left[\frac{x_k{V}_x+y_k{V}_y}{\sqrt{x_k^2+y_k^2+h^2}}+\frac{y_k{V}_y-\left(b-x_k\right)V_x}{\sqrt{{\left(b-x_k\right)}^2+y_k^2+h^2}}\right]\\ {\alpha }_k=\arctan \left(\frac{y_k}{x_k}\right)\\ {\beta }_k=\arctan \frac{h}{\sqrt{x_k^2+y_k^2}}\end{array}\right]\end{array}$$

(32)

Make

$$R_R=\sqrt{x_k^2+y_k^2+h^2}{R}_T=\sqrt{{\left(b-x_k\right)}^2+y_k^2+h^2}$$

(33)

Then

$$\begin{array}{ll}\frac{\partial f_{\mathit{Dk}}}{\partial x_k}=& -\frac{1}{\lambda }\left(\frac{\left(y_k^2+h^2\right)V_x-x_k{y}_k{V}_y}{R_R^3}\right.\\ \left.+\frac{\left(y_k^2+h^2\right)V_x+y_k\left(b-x_k\right)V_y}{R_T^3}\right)\end{array}$$

(34)

$$\begin{array}{ll}\frac{\partial f_{\mathit{Dk}}}{\partial y_k}=& -\frac{1}{\lambda }\left(\frac{\left(x_k^2+h^2\right)V_y-x_k{y}_k{V}_x}{R_R^3}\right.\\ \left.+\frac{\left[{\left(b-x_k\right)}^2+h^2\right]V_y+y_k\left(b-x_k\right)V_x}{R_T^3}\right)\end{array}$$

(35)

$$\frac{\partial f_{\mathit{Dk}}}{\partial h}=\frac{h}{\lambda }\left[\frac{x_k{V}_x+y_k{V}_y}{R_R^3}+\frac{y_k{V}_y-\left(b-x_k\right)V_x}{R_T^3}\right]$$

(36)

$$\frac{\partial f_{\mathit{Dk}}}{\partial V_x}=-\frac{1}{\lambda }\left(\frac{x_k}{R_R}+\frac{x_k-b}{R_T}\right)\frac{\partial f_{\mathit{Dk}}}{\partial V_y}=-\frac{V_y}{\lambda }\left(\frac{1}{R_R}+\frac{1}{R_T}\right)$$

(37)

$$\frac{\partial {\alpha }_k}{\partial x_n}=-\frac{y_k}{R_R^2}\frac{\partial {\alpha }_k}{\partial y_n}=\frac{x_k}{R_R^2}\frac{\partial {\alpha }_k}{\partial h}=\frac{\partial {\alpha }_k}{\partial V_x}=\frac{\partial {\alpha }_k}{\partial V_y}=0$$

(38)

$$\begin{array}{l}\frac{\partial {\beta }_k}{\partial x_k}=-\frac{h{x}_k}{R_R^3+R_R{h}^2}\frac{\partial {\beta }_k}{\partial y_k}=-\frac{h{y}_k}{R_R^3+R_R{h}^2}\\ \frac{\partial {\beta }_k}{\partial h}=\frac{\sqrt{x_k^2+y_k^2}}{x_k^2+y_k^2+h^2}\frac{\partial {\beta }_k}{\partial V_x}=\frac{\partial {\beta }_k}{\partial V_y}=0\end{array}$$

(39)

h _X (k) is

$$h_X\left(k\right)={\left[\begin{array}{lllll}\frac{\partial f_{\mathit{Dk}}\left(X\right)}{\partial x_k}\hfill & \frac{\partial f_{\mathit{Dk}}\left(X\right)}{\partial y_k}\hfill & \frac{\partial f_{\mathit{Dk}}\left(X\right)}{\partial h}\hfill & \frac{\partial f_{\mathit{Dk}}\left(X\right)}{\partial V_x}\hfill & \frac{\partial f_{\mathit{Dk}}\left(X\right)}{\partial V_y}\hfill \\ \frac{\partial {\alpha }_k\left(X\right)}{\partial x_k}\hfill & \frac{\partial {\alpha }_k\left(X\right)}{\partial y_k}\hfill & \frac{\partial {\alpha }_k\left(X\right)}{\partial h}\hfill & \frac{\partial {\alpha }_k\left(X\right)}{\partial V_x}\hfill & \frac{\partial {\alpha }_k\left(X\right)}{\partial V_y}\hfill \\ \frac{\partial {\beta }_k\left(X\right)}{\partial x_k}\hfill & \frac{\partial {\beta }_k\left(X\right)}{\partial y_k}\hfill & \frac{\partial {\beta }_k\left(X\right)}{\partial h}\hfill & \frac{\partial {\beta }_k\left(X\right)}{\partial V_x}\hfill & \frac{\partial {\beta }_k\left(X\right)}{\partial V_y}\hfill \end{array}\right]}_{X=\stackrel{\wedge }{X}\left(k\left|k\right.\right)}$$

(40)

Simulations and analysis

In order to verify the proposed accurate tracking algorithm, the parameters of FSRexperiment system listed in references [3]–[6] are taken as the simulation parameters.

To prove the validity of the proposed method for initial state estimation, simulationsare implemented with the parameters list in Table 1 usingclassical method and the new method, respectively, under the same parameters exceptstandard deviations of Doppler shift, azimuth, and elevation denoted by (σ _f , σ _θ , σ _β ). The results are averaging performed in terms of 100 independentimplementations, as shown in Table 2, and their averaging runtime are shown in Table 3.

Table 1 Simulation parameters

Table 2 Simulation results of initial state estimation

Table 3 Averaging run time of the three methods

Table 2 shows the initial statement estimation errors ofclassical method and the proposed method under the condition of different targetmeasured parameters estimation accuracy. It is obviously that the statement estimationerrors are affected by the estimation accuracy of target measured parameters. As can beseen, with high accuracy of target parameter measurement (i.e.,σ _f = 0.1 Hz, σ _θ = 0.1°, σ _β = 0.1°), estimation errors of the two methods are both within the errortolerance and the proposed method has higher accuracy. However, when measurement noiseincreases (i.e.,σ _f = 0.5 Hz, σ _θ = 0.5°, σ _β = 0.5°), the initial state estimation error of classical methodsignificantly increases, while that of proposed method is much smaller than the former.In particular, the velocity estimation error along the direction of the baseline isfairly small.

To compare the performances of the two methods more intuitively, simulations areimplemented in case that initial state estimations obtained by the two methods are takenas the initial filtering value for EKF. The statement estimation results are shown inFigure 2. It can be seen that the initial statementestimation errors of Gauss–Newton iteration method are so large that the trackingresults are divergent, while the EKF algorithm works very well with the initialstatement obtained by the proposed method.

Figure 2

The effect of initial statement estimation errors on target tracking.(a) Tracking results with the initial statement obtained viaGauss–Newton method. (b) Tracking results with the Initial statementobtained via proposed method.

The averaged results performed in terms of 100 independent implementations are shown inFigure 3. In Figure 3a–c,RMSE of X-, Y-, and H-coordinates using initial valueobtained by proposed method are all plotted by dotted line, and the ones using initialvalue obtained by classical method are plotted by solid line. It can be seen that theproposed method gets more precise initial filtering value, which results in the fasterconvergence and more precise estimation of EKF. While taking the initial state estimatedby classical method as the filtering initial value slows down the convergence speed ofEKF. In addition, the estimated trajectory has not converged at the time of targetcrossing baseline. The large RMSE of V _x estimation (as shown in Figure 3d) by classical methodslows down the convergence of EKF. The initial value of V _x is more accurate obtained by the proposed method, which is helpful to achieverapid filtering convergence.

Figure 3

RMSE of target parameters estimation. (a) RMSE of x estimation(b) RMSE of y estimation (c) RMSE of H estimation (d) RMSEof velocity estimation.

Then, the overall tracking performances of classical method, new method, and theunscented Kalman filter (UKF) algorithm (taking the initial state obtained by new methodas first filtering value) are compared in Figure 4.

Figure 4

Filtering results of the three methods, with measurement noise (0.5 Hz,0.5°,0.5°). (a) Tracking results with three methods(b) RMSE of x (b) RMSE of y (c) RMSE of H (d) RMSEof Vx (e) RMSE of Vy.

Figure 4 shows that the initial state estimation error ofclassical method is relatively large in the general measurement noise level. In theprocess of Gauss–Newton iteration, the error is not significantly reduced due toits unstable performance. In the same condition, the initial error of the proposedmethod is quite small, and gradually decreases until convergence in the follow-uptracking stage. UKF algorithm has a deterioration trend in follow-up tracking.

To compare the computation of different methods, averaging is performed in terms of 100independent implementations. The time required for one complete filtering of threemethods are as follows:

Through analysis of the simulation results, we obtain the following conclusions:

• (1) The acquisition of accurate initial state estimation by classical method required very high precision measurement. The initial estimation error significantly increases by the measurement noise.

• (2) The proposed method demonstrated its superiority in initial state estimation. In general measurement accuracy, its estimation error is smaller than the classical method’s by an order of magnitude. Applied to EKF, the initial value obtained by the proposed method can effectively improve the convergence speed. Then a stable track is able to be formed before the target crossing baseline.

• (3) In general measurement precision, classical method showed non-steady performance and fluctuation of estimation error. Tracking accuracy of UKF algorithm is moderate. But UKF’s estimation error has an increasing trend because of its sensitivity to initial error. EKF algorithm can achieve high tracking accuracy and fast convergence with a large tolerance for initial error. Furthermore, its computation is the smallest among three algorithms.

Conclusions

In this article, an accurate 3D tracking method in forward scatter radar is presented.To solve the problem of filter divergence results from the large initial estimationerror, an accurate initial state estimation approach based on analytic derivation and LMalgorithm is proposed, which can improve the accuracy of initial state estimationwithout requiring high-precise measurement. For the purpose of reducing the computationof filtering and speeding up the convergence rate, the accurate initial state estimationis combined with EKF algorithm. Then, an accurate 3D tracking in forward scatter radaris derived and verified by the simulation results.