Dynamic Programming-Based Track-before-Detect Algorithm for Weak Maneuvering Targets in Range–Doppler Plane

Abstract

This paper focuses on detecting and tracking maneuvering weak targets in the range–Doppler (RD) plane with the track-before-detect (TBD) algorithm based on dynamic programming (DP). Traditional DP-TBD algorithms integrate target detection and tracking in their framework while searching the paths provided by a predefined model of the kinematic properties within the constraints allowed. However, both the approximate motion model used in the RD plane and the model mismatch caused when the target undergoes a maneuver can degrade the TBD performance sharply. To address these issues, this paper accurately describes the evolution of the RD equation based on Constant Acceleration (CA) and Coordinated Turn (CT) motion models with the process noise in the Cartesian coordinate system, and it also employs a recursive method to estimate the parameters in the equations for efficient energy accumulation and path searches. Facing the situation that targets energy accumulation during the DP iteration process will lead to an expansion of the target energy accumulation process. This paper designs a more efficient Optimization Function (OF) to inhibit the expansion effect, improve the resolution of the nearby targets, and increase the detection probability of the weak targets simultaneously. In addition, to search the trajectory more efficiently and accurately, we improved the process of DP multi-frame accumulation, thus reducing the computation amount of large-scale searches. Finally, the effectiveness of the proposed method for CA and CT motion target detection and tracking is verified by many of the simulation experiments that were conducted in this paper.

1. Introduction

Target detection and tracking is the most critical task of a radar system, especially the detection and tracking of weak maneuvering targets (which represents a significant problem to be solved). The traditional radar processing system firstly performs signal processing for the echo data, and it then performs data processing to achieve target detection [1], as well as tracking [2,3], as shown in Figure 1. However, since detect-before-track (DBT) performs target detection by setting a threshold, such as a constant false alarming ratio (CFAR), the threshold setting becomes a limitation of the target detection performance, and, in the case of a low signal-to-noise ratio (SNR) or the overlapping of strong and weak targets, the threshold detection will also have the problem of weak target leakage. Therefore, to solve the problem of weak target detection and tracking, track-before-detect (TBD) [4,5,6] uses a non-threshold or low-threshold method to process the raw radar echo data, thus improving the detection probability of weak targets. TBD fuses the detection and tracking processes with a unified framework [4,7,8] to process the echo data and utilize multi-frame joint processing [9,10,11] to make full use of the information in the raw echo data to achieve more efficient detection and tracking, as shown in Figure 2.

A large number of branches of the TBD algorithms have also appeared according to the different focuses of the accumulated features of TBD algorithms, such as TBD based on particle filtering (PF-TBD) [12,13], TBD based on random finite set theory (RFS-TBD) [14,15], TBD based on dynamic programming (DP-TBD) [16,17], TBD based on the Hough Transform (HT-TBD) [18,19,20], and Velocity Filtering (VF) based TBD (VF-TBD) [21,22,23]. Among them, PF-TBD and RFS-TBD are similar to the traditional methods that are based on single-frame raw echo accumulation and parameter estimation due to a posterior probability prediction and the updating of targets, and the single-frame process achieves more real-time echo data processing than multi-frame accumulation. Contrasted with single-frame accumulation, DP-TBD, HT-TBD, and VF-TBD engage multi-frame signal accumulation processing methods. The utilization of multi-frame signal accumulation facilitates the comprehensive capture of echo information, thereby enhancing the probability of target detection and concurrently suppressing excessive, cluttered, and unwieldy noise. Nevertheless, the HT-TBD technique necessitates threshold processing, undermining the effective use of the original echo data and subsequently causing performance degradation. The VF-TBD algorithm is mainly applied to the TBD processing of optical and infrared images; in addition, the application of velocity filters can greatly increase the computational effort. Hence, this paper prioritizes the assessment of the DP-TBD algorithm.

As shown in Figure 3, the existing research on the DP-TBD algorithm mainly focuses on the following aspects. DP-TBD improves the detection probability of weak targets by converting a multi-stage problem in the same batch processing into multiple single-stage problems searching for a local optimum while accumulating the target energy along feasible paths after single-stage processing. DP-TBD was first used for weak target detection and tracking in optical images [24], and then DP-TBD was successively applied to radar systems [16], over-the-horizon radar systems [25], and airborne radar systems [26].

The initial study focuses on a Constant Velocity (CV) target motion model [27,28]. However, there is a great deal of maneuvering during the actual target movement process, which will lead to a mismatch of the algorithmic model and could seriously affect the accuracy and detection probability of the algorithm. Therefore, the research on DP-TBD for maneuvering targets such as the Coordinated Turn (CT), Constant Acceleration (CA), and Current Statistical (CS) models in a Cartesian coordinate system has been somewhat advanced [29,30,31,32,33,34]. In addition, the research has expanded from the Cartesian coordinate system to the Doppler radar coordinate system [35,36]. However, the existing research is still unable to accurately describe the motion model in the radar observation coordinate system, thus resulting in motion mismatch and a severe degradation of algorithm performance.

In addition, optimizing functions and some parameters also have a significant impact on algorithm performance. The research on the performance of the DP-TBD algorithm on the factors of OF and noise has received wide attention [37]. Afterward, the performance of DP-TBD was analyzed in detail by the study of Extreme Value Theory (EVT), Generalized EVT (GEVT), and the Peaks Over Threshold (POT) model of EVT, and the research of value functions in the field of DP-TBD was promoted [27,28,38]. The OF of DP-TBD in the existing research cannot utilize the correlation of intra-frame data well and fails to fully explore the information in the echoes, so the design of the OF of DP-TBD needs a great deal more research.

DP-TBD suffers from bit catastrophe and computational explosion when dealing with multi-target detection and tracking problems, and to alleviate these problems, DP-TBD algorithms based on Successive Target Cancellation (STC) and Single-Pass STC (SP-STC) have been proposed [39]. To improve the performance of accurate target tracking, a DP-TBD method based on Parallel Target Cancellation (PTC) was proposed [40]. However, there is still no better method to reduce the computational dimension in order to reduce the computational amount. Meanwhile, neighboring targets interfere with each other, which also plagues multi-target tracking and detection.

In this paper, we investigate an approach of DP-TBD for weak targets, which mainly analyzes the problems of motion model mismatch and proximity target discrimination in the presence of process noise. First, we derive the exact evolution equations of CA and CT models with Gaussian white noise for the process noise in the RD domain under the Cartesian coordinate system. These equations pertain to the significant reduction in the probability of a motion model mismatch, thereby contributing to a marked increase in energy accumulation efficiency and an enhancement in the SNR of the accumulated target. Secondly, the  formulation of the Optimization Function (OF) constitutes the linchpin of the DP-TBD algorithm. The OF is tasked with representing the disparities between the target and the clutter in aspects such as signal, correlation, and  motion characteristics, among others. By employing Dynamic Programming for accumulating the OF, the target signal strength can be amplified, thereby enabling an enhancement in target detection and tracking performance. The OF designed in this paper is aimed at problems such as the overlapping of target and noise energy and the overlapping of the energy of targets with similar strengths and weaknesses during the accumulation process, and it fully exploits the spatial and temporal correlation characteristics of the original echo information of a single frame and the batch processing of inter-frame echoes, so as to realize the effect of improving the resolution of the adjacent targets and effectively suppressing the noise. Finally, a new search strategy for DP-TBD is proposed to reduce the search range of DP-TBD and improve the search efficiency to avoid the problems of dimensionality disaster and computational explosion. Thus far, the algorithm proposed in this paper is the first DP-TBD algorithm that accurately describes the evolution equations of CA and CT models with the process noise in the RD domain, and it uses the DP-TBD algorithm based on the OF correction strategy. Compared with the existing state-of-the-art algorithms, our method is superior in terms of robustness and detection performance.

The rest of this paper is organized in the following manner. Section 2 briefly presents the concepts and principles of the target motion model and observation model. Section 3 introduces the evolution process of the motion model proposed in this paper and analyzes the improvement process of OF design and the energy accumulation process. The simulation and verification of the proposed CA-OF-TBD and CT-OF-TBD algorithms are arranged in Section 4. There is a summarization of this paper and outlooks for the future in Section 5.

2. Problem Statement

2.1. Measurement Model

This study is centered around the observation model within the RD domain. There is the assumption that each frame within the RD domain possesses an observation area composed of $N_r\times N_d$ cells with resolutions of ${\Delta }_r$ and ${\Delta }_d$. ${\Delta }_r$ and ${\Delta }_d$ denote the resolution of range and Doppler, respectively. The intensity value for the observation cell $(n_r,n_d)$ within the kth frame of each batch process is represented by $z_k(n_r,n_d)$. As such, the set of observations within the RD domain for the kth frame is formulated as follows:

$$\begin{array}{c}\hfill {\mathbf{z}}_k=\left\{z_k\left(n_r,n_d\right),n_r=1,2,\dots ,N_r,n_d=1,2,\dots ,N_d\right\}.\end{array}$$

(1)

The observation set from the initial frame through to the final, or Kth frame, within a singular batch processing instance is defined as follows:

$$\begin{array}{c}\hfill {\mathbf{Z}}_{1:K}=\left\{{\mathbf{z}}_1,{\mathbf{z}}_2,\dots ,{\mathbf{z}}_K\right\}.\end{array}$$

(2)

In instances where no target is present within the observation plane, the value of the observation cell $(n_r,n_d)$ for the kth frame is defined as follows:

$$\begin{array}{c}\hfill z_k\left(n_r,n_d\right)={\omega }_k\left(n_r,n_d\right),\end{array}$$

(3)

where ${\omega }_k(n_r,n_d)$ symbolizes the noise factor. While Gaussian noise was presumed for the scope of this study, the algorithm retains its efficacy with other noise types.

It is known that the observation of a target is at its maximum in the cell $(n_r,n_d)$ where it is situated. It is presumed that the observations of each target within a singular frame are independently and identically distributed (i.i.d.). Hence, in the presence of a target, the echo intensity for the observation cell $(n_r,n_d)$ during the kth frame is defined as follows:

$$\begin{array}{c}\hfill z_k\left(n_r,n_d\right)=s_k\left(n_r,n_d\right)+{\omega }_k\left(n_r,n_d\right),\end{array}$$

(4)

where $s_k\left(n_r,n_d\right)$ denotes the target echo intensity.

2.2. Target Motion Model

This investigation primarily explores DP-TBD for weak maneuvering targets. This section focuses on the CA and CT maneuver model in the Cartesian coordinate system. The focus rests upon the CA and CT maneuvering models within the Cartesian coordinate system as these are the more commonly observed motion patterns in real-world applications. Noise can influence target motion to varying extents. However, in this manuscript, the advantages and limitations of the developed algorithm are initially discussed, necessitating an examination of motion models without the addition of noise.

2.2.1. CA Motion Model

Within the sub-division of this study, the motion state is established as a six-dimensional vector, symbolized as follows:

$$\begin{array}{c}\hfill {\mathbf{x}}_k={\left[x_k{\dot{x}}_k{\ddot{x}}_k{y}_k{\dot{y}}_k{\ddot{y}}_k\right]}^{\top },\end{array}$$

(5)

where $x_k$, ${\dot{x}}_k$, and ${\ddot{x}}_k$ denote the position, velocity, and acceleration of the target along the x-axis, respectively, while $y_k$, ${\dot{y}}_k$, and ${\ddot{y}}_k$ denote the corresponding parameters along the y-axis. The symbol ⊤ delineates matrix transposition. Conforming to the kinematic properties associated with the CA model, accelerations ${\ddot{x}}_k$ and ${\ddot{y}}_k$ maintain constant values over time.

In the Cartesian coordinate system, the state transition equation is formulated as follows:

$$\begin{array}{c}\hfill {\mathbf{x}}_{k+1}=\mathbf{F}{\mathbf{x}}_k+{\mathbf{w}}_k,\end{array}$$

(6)

where F is the state transition matrix of the target and ${\mathbf{w}}_k$ is the process noise of the target. In the state estimation of the CA model, the process noise is the rate of change in the target acceleration. Therefore, F and ${\mathbf{w}}_k$ are denoted as follows:

$$\begin{array}{c}\hfill \mathbf{F}=\left[\begin{array}{cccccc}1& T& {\displaystyle \frac{T^2}{2}}& 0& 0& 0\\ 0& 1& T& 0& 0& 0\\ 0& 0& 1& 0& 0& 0\\ 0& 0& 0& 1& T& {\displaystyle \frac{T^2}{2}}\\ 0& 0& 0& 0& 1& T\\ 0& 0& 0& 0& 0& 1\end{array}\right],{\mathbf{w}}_k=\left[\begin{array}{c}{\displaystyle \frac{1}{6}}{\stackrel{⃛}{x}}_k{T}^3\\ {\displaystyle \frac{1}{6}}{\stackrel{⃛}{y}}_k{T}^3\\ {\displaystyle \frac{1}{2}}{\stackrel{⃛}{x}}_k{T}^2\\ {\displaystyle \frac{1}{2}}{\stackrel{⃛}{y}}_k{T}^2\\ {\stackrel{⃛}{x}}_{k}T\\ {\stackrel{⃛}{y}}_{k}T\end{array}\right],\end{array}$$

(7)

where T represents the time interval between the kth frame and the $k+1$st frame, and ${\stackrel{⃛}{x}}_k$ and ${\stackrel{⃛}{y}}_k$ denote the rate of change in acceleration in the x and y directions, respectively.

By evaluating Equation (6), deduction of both the motion and velocity equations for the CA model along the x and y axes can be achieved as follows:

$$\begin{array}{cc}\hfill x_{k+1}=& x_0+{\dot{x}}_{0}kT+{\displaystyle \frac{{\ddot{x}}_0{k}^2{T}^2}{2}}+{\displaystyle \frac{1}{6}}{\stackrel{⃛}{x}}_k{k}^3{T}^3\hfill \\ \hfill {\dot{x}}_{k+1}=& {\dot{x}}_0+{\ddot{x}}_{0}kT+{\displaystyle \frac{1}{2}}{\stackrel{⃛}{x}}_k{k}^2{T}^2\hfill \\ \hfill {\ddot{x}}_{k+1}=& {\ddot{x}}_0+{\stackrel{⃛}{x}}_{k}kT,\hfill \end{array}$$

(8)

and

$$\begin{array}{cc}\hfill y_{k+1}=& y_0+{\dot{y}}_{0}kT+{\displaystyle \frac{{\ddot{y}}_0{k}^2{T}^2}{2}}+{\displaystyle \frac{1}{6}}{\stackrel{⃛}{y}}_k{k}^3{T}^3\hfill \\ \hfill {\dot{y}}_{k+1}=& {\dot{y}}_0+{\ddot{y}}_{0}kT+{\displaystyle \frac{1}{2}}{\stackrel{⃛}{y}}_k{k}^2{T}^2\hfill \\ \hfill {\ddot{y}}_{k+1}=& {\ddot{y}}_0+{\stackrel{⃛}{y}}_{k}kT,\hfill \end{array}$$

(9)

where $x_0$, ${\dot{x}}_0$, ${\ddot{x}}_0$, and ${\stackrel{⃛}{x}}_0$ depict the initial position, velocity, acceleration, and rate of change in acceleration of the target along the x-axis in the Cartesian coordinate system, respectively. Similarly, $y_0$, ${\dot{y}}_0$, ${\ddot{y}}_0$, and ${\stackrel{⃛}{y}}_0$ characterize the corresponding initial parameters along the y-axis.

2.2.2. CT Motion Model

The target motion state in this section is represented as a four-dimensional vector as

$$\begin{array}{c}\hfill {\mathbf{x}}_k={\left[x_k{\dot{x}}_k{y}_k{\dot{y}}_k{\alpha }_k\right]}^{\top },\end{array}$$

(10)

where $x_k$, ${\dot{x}}_k$, $y_k$, and ${\dot{y}}_k$ denote the position, velocity, and acceleration of the target along the x and y axes, respectively. Also, based on its state transition matrix, the state transition equation of the target can be obtained as follows:

$$\begin{array}{c}\hfill \left[\begin{array}{c}{x}_{k+1}\\ {\dot{x}}_{k+1}\\ y_{k+1}\\ {\dot{y}}_{k+1}\\ {\alpha }_{k+1}\end{array}\right]=\left[\begin{array}{cccc}1& {\displaystyle \frac{sin\left({\alpha }_{k}T\right)}{{\alpha }_k}}& 0& {\displaystyle \frac{cos\left({\alpha }_{k}T\right)-1}{{\alpha }_k}}\\ 0& cos\left({\alpha }_{k}T\right)& 0& -sin\left({\alpha }_{k}T\right)\\ 0& {\displaystyle \frac{1-cos\left({\alpha }_{k}T\right)}{{\alpha }_k}}& 1& {\displaystyle \frac{sin\left({\alpha }_{k}T\right)}{{\alpha }_k}}\\ 0& sin\left({\alpha }_{k}T\right)& 0& cos\left({\alpha }_{k}T\right)\\ 0& 0& 0& 1\end{array}\right]\left[\begin{array}{c}{x}_k\\ {\dot{x}}_k\\ y_k\\ {\dot{y}}_k\\ {\alpha }_k\end{array}\right]+\left[\begin{array}{c}0\\ 0\\ 0\\ 0\\ {\dot{\alpha }}_{k}T\end{array}\right],\end{array}$$

(11)

where ${\alpha }_k$ and ${\dot{\alpha }}_k$ symbolize the target turning rate and its rate of change, respectively. In utilizing Equation (11), the following motion equations are obtained:

$$\begin{array}{cc}\hfill x_{k+1}=& x_0+{\displaystyle \frac{sin\left({\alpha }_{k}kT\right)}{{\alpha }_k}}{\dot{x}}_0+{\displaystyle \frac{cos\left({\alpha }_{k}kT\right)-1}{{\alpha }_k}}{\dot{y}}_0,\hfill \\ \hfill {\dot{x}}_{k+1}=& cos\left({\alpha }_{k}kT\right){\dot{x}}_0-sin\left({\alpha }_{k}kT\right){\dot{y}}_0,\hfill \\ \hfill {\alpha }_{k+1}=& {\alpha }_k+{\dot{\alpha }}_{k}kT,\hfill \end{array}$$

(12)

and

$$\begin{array}{cc}\hfill y_{k+1}=& y_0+{\displaystyle \frac{1-cos\left({\alpha }_{k}kT\right)}{{\alpha }_k}}{\dot{x}}_0+{\displaystyle \frac{sin\left({\alpha }_{k}kT\right)}{{\alpha }_k}}{\dot{y}}_0,\hfill \\ \hfill {\dot{y}}_{k+1}=& sin\left({\alpha }_{k}kT\right){\dot{x}}_0+cos\left({\alpha }_{k}kT\right){\dot{y}}_0,\hfill \\ \hfill {\alpha }_{k+1}=& {\alpha }_k+{\dot{\alpha }}_{k}kT.\hfill \end{array}$$

(13)

This segment delineates the motion of targets within the Cartesian coordinate system, laying the groundwork for the subsequent transformation of this system into the RD domain. Hence, engaging in the modeling of both the motion and the pertinent observation model of targets.

3. Methods

This section comprehensively outlines the transition from CA and CT models within the Cartesian coordinate system to the RD domain. Further, it explicates the methodology employed in this study for devising both the optimization and trajectory recovery functions within the DP-TBD accumulation process. Consequently, this enables a more efficient energy accumulation during the DP-TBD iterative process.

3.1. Evolutionary Process

3.1.1. Evolutionary Process of CA Model

The range and Doppler states of a target at the kth instant can be deduced from its velocity and position parameters within the Cartesian coordinate system:

$$\begin{array}{cc}\hfill r_k=& \sqrt{x_k^2+y_k^2}\hfill \end{array}$$

(14)

$$\begin{array}{cc}\hfill d_k=& {\displaystyle \frac{x_k{\dot{x}}_k+y_k{\dot{y}}_k}{\sqrt{x_k^2+y_k^2}}},\hfill \end{array}$$

(15)

where $r_k$ and $d_k$ denote the range and Doppler status of the target at the kth instant, respectively. Incorporating Equations (8) and (9) into Equations (14) and (15) manifests the RD state of target at the kth frame as follows:

$$\begin{array}{cc}\hfill r_k=& \sqrt{\begin{array}{c}{\left(x_0+{\dot{x}}_{0}kT+{\displaystyle \frac{1}{2}}{\ddot{x}}_0{k}^2{T}^2+{\displaystyle \frac{1}{6}}{\stackrel{⃛}{x}}_0{k}^3{T}^3\right)}^2+{\left(y_0+{\dot{y}}_{0}kT+{\displaystyle \frac{1}{2}}{\ddot{y}}_0{k}^2{T}^2+{\displaystyle \frac{1}{6}}{\stackrel{⃛}{y}}_0{k}^3{T}^3\right)}^2\hfill \end{array}}\hfill \\ \hfill =& \sqrt{\begin{array}{c}\left(x_0^2+y_0^2\right)+2kT\left({\dot{x}}_0{x}_0+{\dot{y}}_0{y}_0\right)+k^2{T}^2\left({\dot{x}}_0^2+{\dot{y}}_0^2+{\ddot{x}}_0{x}_0+{\ddot{y}}_0{y}_0\right)\hfill \\ +{\displaystyle \frac{1}{4}}{k}^3{T}^3\left(3{\ddot{x}}_0{\dot{x}}_0+3{\ddot{y}}_0{\dot{y}}_0+{\stackrel{⃛}{x}}_0{x}_0+{\stackrel{⃛}{y}}_0{y}_0\right)\hfill \\ +{\displaystyle \frac{1}{12}}{k}^4{T}^4\left(3{{\ddot{x}}_0}^2+{3{\ddot{y}}_0}^2+4{\stackrel{⃛}{x}}_0{\dot{x}}_0+4{\stackrel{⃛}{y}}_0{\dot{y}}_0\right)\hfill \\ +{\displaystyle \frac{1}{36}}{k}^5{T}^5\left(6{\stackrel{⃛}{x}}_0{\ddot{x}}_0+6{\stackrel{⃛}{y}}_0{\ddot{y}}_0\right)+{\displaystyle \frac{1}{36}}{k}^6{T}^6\left({{\stackrel{⃛}{x}}_0}^2+{{\stackrel{⃛}{y}}_0}^2\right)\hfill \end{array}},\hfill \end{array}$$

(16)

$$\begin{array}{cc}\hfill d_k=& {\displaystyle \frac{1}{r_k}}\left(\left(x_0+{\dot{x}}_{0}kT+{\displaystyle \frac{1}{2}}{\ddot{x}}_0{k}^2{T}^2+{\displaystyle \frac{1}{6}}{\stackrel{⃛}{x}}_0{k}^3{T}^3\right)\left({\dot{x}}_0+{\ddot{x}}_{0}kT+{\displaystyle \frac{1}{2}}{\stackrel{⃛}{x}}_0{k}^2{T}^2\right)\right.\hfill \\ \hfill & \left.+\left(y_0+{\dot{y}}_{0}kT+{\displaystyle \frac{1}{2}}{\ddot{y}}_0{k}^2{T}^2{\displaystyle \frac{1}{6}}{\stackrel{⃛}{x}}_0{k}^3{T}^3\right)\left({\dot{y}}_0+{\ddot{y}}_{0}kT+{\displaystyle \frac{1}{2}}{\stackrel{⃛}{y}}_0{k}^2{T}^2\right)\right)\hfill \\ \hfill =& {\displaystyle \frac{1}{r_k}}\left(\left(x_0{\dot{x}}_0+y_0{\dot{y}}_0\right)+kT\left({\dot{x}}_0^2+{\dot{y}}_0^2+x_0{\ddot{x}}_0+y_0{\ddot{y}}_0\right)\right.\hfill \\ \hfill & \left.+{\displaystyle \frac{1}{2}}{k}^2{T}^2\left(3{\dot{x}}_0{\ddot{x}}_0+3{\dot{y}}_0{\ddot{y}}_0+{\stackrel{⃛}{x}}_0{x}_0+{\stackrel{⃛}{y}}_0{y}_0\right)\right)\hfill \\ \hfill & +{\displaystyle \frac{1}{6}}{k}^3{T}^3\left(3{\ddot{x}}_0^2+3{\ddot{y}}_0^2+4{\stackrel{⃛}{x}}_0{\dot{x}}_0+4{\stackrel{⃛}{y}}_0{\dot{y}}_0\right)\hfill \\ \hfill & +{\displaystyle \frac{5}{72}}{k}^4{T}^4\left(6{\stackrel{⃛}{x}}_0{\ddot{x}}_0+6{\stackrel{⃛}{y}}_0{\ddot{y}}_0\right)+{\displaystyle \frac{1}{12}}{k}^5{T}^5\left({{\stackrel{⃛}{x}}_0}^2+{{\stackrel{⃛}{y}}_0}^2\right).\hfill \end{array}$$

(17)

Since the scenario is set up as a slow and weakly maneuvering target, the separate acceleration rate of change in the ${\stackrel{⃛}{x}}_0$ and ${\stackrel{⃛}{y}}_0$ terms in the last terms of Equations (16) and (17) can be directly omitted. The omission of this term has a small impact on the accuracy of the target motion derivation process.

In this paper, the product of range and Doppler is posited as ${\sigma }_k$, as expressed below:

$${\sigma }_k=r_k{d}_k=x_k{\dot{x}}_k+y_k{\dot{y}}_k.$$

(18)

Within the CA model framework, the first-, second-, and third-order derivatives of ${\sigma }_k$ are required to adequately depict the RD state, and see Appendix A for a detailed derivation. Due to the weak maneuverability of the target, the sum-of-squares term of the acceleration transformation rate of the highest order is omitted and the rest can be defined as follows:

$$\begin{array}{cc}\hfill {\sigma }_k=& x_k{\dot{x}}_k+y_k{\dot{y}}_k\hfill \\ \hfill {\dot{\sigma }}_k=& {\dot{x}}_k^2+x_k{\ddot{x}}_k+{\dot{y}}_k^2+y_k{\ddot{y}}_k\hfill \\ \hfill {\ddot{\sigma }}_k=& 3{\dot{x}}_k{\ddot{x}}_k+3{\dot{y}}_k{\ddot{y}}_k+{\stackrel{⃛}{x}}_k{x}_k+{\stackrel{⃛}{y}}_k{y}_k\hfill \\ \hfill {\stackrel{⃛}{\sigma }}_k=& 3{\ddot{x}}_k^2+3{\ddot{y}}_k^2+4{\stackrel{⃛}{x}}_k{\dot{x}}_k+4{\stackrel{⃛}{y}}_k{\dot{y}}_k\hfill \\ \hfill {\stackrel{⃜}{\sigma }}_k=& 6{\stackrel{⃛}{x}}_k{\ddot{x}}_k+6{\stackrel{⃛}{y}}_k{\ddot{y}}_k.\hfill \end{array}$$

(19)

By incorporating Equations (8) and (9) into Equation (19), we obtain the RD state for the CA model:

$$\begin{array}{cc}\hfill {\mathsf{\sigma }}_k=\left[\begin{array}{c}{\sigma }_k\\ {\dot{\sigma }}_k\\ {\ddot{\sigma }}_k\\ {\stackrel{⃛}{\sigma }}_k\\ {\stackrel{⃜}{\sigma }}_k\end{array}\right]=& \left[\begin{array}{ccccc}1& kT& {\displaystyle \frac{1}{2}}{k}^2{T}^2& {\displaystyle \frac{1}{6}}{k}^3{T}^3& {\displaystyle \frac{1}{24}}{k}^4{T}^4\\ 0& 1& kT& {\displaystyle \frac{1}{2}}{k}^2{T}^2& {\displaystyle \frac{1}{6}}{k}^3{T}^3\\ 0& 0& 1& kT& {\displaystyle \frac{1}{2}}{k}^2{T}^2\\ 0& 0& 0& 1& kT\\ 0& 0& 0& 0& 1\end{array}\right]\left[\begin{array}{c}{\sigma }_0\\ {\dot{\sigma }}_0\\ {\ddot{\sigma }}_0\\ {\stackrel{⃛}{\sigma }}_0\\ {\stackrel{⃜}{\sigma }}_0\end{array}\right],\hfill \end{array}$$

(20)

where ${\sigma }_0$, ${\dot{\sigma }}_0$, ${\ddot{\sigma }}_0$, ${\stackrel{⃛}{\sigma }}_0$, and ${\stackrel{⃜}{\sigma }}_0$ are the initial values of the product and its derivatives of all the orders of range and Doppler.

The introduction of the target initial state into Equation (19), followed by the incorporation of these initial values into Equations (16) and (17), facilitates the derivation of the relationship between the RD state of the target at the kth frame within the RD domain and its initial state. This relationship is represented as follows:

$$\begin{array}{cc}\hfill r_k=& \sqrt{r_0^2+2kT{\sigma }_0+k^2{T}^2{\dot{\sigma }}_0+{\displaystyle \frac{1}{6}}{k}^3{T}^3{\ddot{\sigma }}_0+{\displaystyle \frac{1}{12}}{k}^4{T}^4{\stackrel{⃛}{\sigma }}_0+{\displaystyle \frac{1}{36}}{k}^5{T}^5{\stackrel{⃜}{\sigma }}_0}\hfill \end{array}$$

(21)

$$\begin{array}{cc}\hfill d_k=& {\displaystyle \frac{1}{r_k}}\left({\sigma }_0+kT{\dot{\sigma }}_0+{\displaystyle \frac{1}{2}}{k}^2{T}^2{\ddot{\sigma }}_0+{\displaystyle \frac{1}{6}}{k}^3{T}^3{\stackrel{⃛}{\sigma }}_0+{\displaystyle \frac{5}{72}}{k}^4{T}^4{\stackrel{⃜}{\sigma }}_0\right).\hfill \end{array}$$

(22)

According to Equations (21) and (22), the RD state of the target at the kth frame can be articulated using its initial RD state, which is denoted by ${\mathsf{\sigma }}_k$. Given that ${\mathsf{\sigma }}_k$ remains constant within each processing batch, its iteration substantially mitigates computational demands, thus boosting efficiency. Concurrently, the translation from Cartesian coordinates to the RD domain allows for a notable reduction in nonlinear transformation errors, enhancing accumulation efficiency and diminishing the likelihood of model mismatch.

3.1.2. Evolutionary Process of CT Model

The RD state for the CT model at the kth frame can be procured by integrating Equations (12) and (13) into Equations (14) and (15):

$$\begin{array}{cc}\hfill r_k=& \sqrt{\begin{array}{c}{\left(x_0+{\displaystyle \frac{sin\left({\alpha }_{k}kT\right)}{{\alpha }_k}}{\dot{x}}_0+{\displaystyle \frac{cos\left({\alpha }_{k}kT\right)-1}{{\alpha }_k}}{\dot{y}}_0\right)}^2\hfill \\ +{\left(y_0+{\displaystyle \frac{1-cos\left({\alpha }_{k}kT\right)}{{\alpha }_k}}{\dot{x}}_0+{\displaystyle \frac{sin\left({\alpha }_{k}kT\right)}{{\alpha }_k}}{\dot{y}}_0\right)}^2\hfill \end{array}}\hfill \\ \hfill =& \sqrt{\begin{array}{c}\left(x_0^2+y_0^2\right)+{\displaystyle \frac{2sin\left({\alpha }_{k}kT\right)}{{\alpha }_k}}\left({\dot{x}}_0{x}_0+{\dot{y}}_0{y}_0\right)+{\displaystyle \frac{2-2cos\left({\alpha }_{k}kT\right)}{{{\alpha }_k}^2}}\hfill \\ \times \left({\dot{x}}_0^2+{\dot{y}}_0^2+{\alpha }_k\left(y_0{\dot{x}}_0-x_0{\dot{y}}_0\right)\right)\hfill \end{array}}\hfill \end{array}$$

(23)

$$\begin{array}{cc}\hfill d_k=& {\displaystyle \frac{1}{r_k}}\left(\left(x_0+{\displaystyle \frac{sin\left({\alpha }_{k}kT\right)}{{\alpha }_k}}{\dot{x}}_0+{\displaystyle \frac{cos\left({\alpha }_{k}kT\right)-1}{{\alpha }_k}}{\dot{y}}_0\right)\right.\hfill \\ \hfill & \times \left(cos\left({\alpha }_{k}kT\right){\dot{x}}_0-sin\left({\alpha }_{k}kT\right){\dot{y}}_0\right)\hfill \\ \hfill & +\left(y_0+{\displaystyle \frac{1-cos\left({\alpha }_{k}kT\right)}{{\alpha }_k}}{\dot{x}}_0+{\displaystyle \frac{sin\left({\alpha }_{k}kT\right)}{{\alpha }_k}}{\dot{y}}_0\right)\hfill \\ \hfill & \left.\times \left(sin\left({\alpha }_{k}kT\right){\dot{x}}_0+cos\left({\alpha }_{k}kT\right){\dot{y}}_0\right)\right)\hfill \\ \hfill =& {\displaystyle \frac{1}{r_k}}\left(cos\left({\alpha }_{k}kT\right)\left(x_0{\dot{x}}_0+y_0{\dot{y}}_0\right)\right.\hfill \\ \hfill & \left.+{\displaystyle \frac{sin\left({\alpha }_{k}kT\right)}{{\alpha }_k}}\left({\dot{x}}_0^2+{\dot{y}}_0^2+{\alpha }_k\left(y_0{\dot{x}}_0-x_0{\dot{y}}_0\right)\right)\right).\hfill \end{array}$$

(24)

As per the preceding definition, the CT model computes the first-order derivative for ${\sigma }_k$ to represent the RD state. Below, the initial definitions for ${\sigma }_k$ and its corresponding first-order derivative are, respectively, presented:

$$\begin{array}{cc}\hfill {\sigma }_k=& x_k{\dot{x}}_k+y_k{\dot{y}}_k\hfill \\ \hfill {\dot{\sigma }}_k=& {\dot{x}}_k^2+{\dot{y}}_k^2+{\alpha }_k\left(y_k{\dot{x}}_k-x_k{\dot{y}}_k\right).\hfill \end{array}$$

(25)

Incorporating Equations (12) and (13) into Equation (25) formulates the RD state for the CT model:

$$\begin{array}{cc}\hfill {\mathsf{\sigma }}_k=\left[\begin{array}{c}{\sigma }_k\\ {\dot{\sigma }}_k\end{array}\right]=& \left[\begin{array}{cc}cos\left({\alpha }_{k}kT\right)& {\displaystyle \frac{sin\left({\alpha }_{k}kT\right)}{{\alpha }_k}}\\ -{\alpha }_{k}sin\left({\alpha }_{k}kT\right)& cos\left({\alpha }_{k}kT\right)\end{array}\right]\left[\begin{array}{c}{\sigma }_0\\ {\dot{\sigma }}_0\end{array}\right],\hfill \end{array}$$

(26)

where ${\sigma }_0$ and ${\dot{\sigma }}_0$ are the initial values of the range and Doppler product.

According to Equation (24), the range and Doppler state of the kth frame derive their values from the initial range and Doppler states as follows:

$$\begin{array}{cc}\hfill r_k=& \sqrt{r_0^2+{\displaystyle \frac{2sin\left({\alpha }_{k}kT\right)}{{\alpha }_k}}{\sigma }_0+{\displaystyle \frac{2-2cos\left({\alpha }_{k}kT\right)}{{\alpha }_k^2}}{\dot{\sigma }}_0}\hfill \end{array}$$

(27)

$$\begin{array}{cc}\hfill d_k=& {\displaystyle \frac{1}{r_k}}\left(cos\left({\alpha }_{k}kT\right){\sigma }_0+{\displaystyle \frac{sin\left({\alpha }_{k}kT\right)}{{\alpha }_k}}{\dot{\sigma }}_0\right).\hfill \end{array}$$

(28)

The evolution of RD, as depicted in Equations (27) and (28), is determined by the turning rate ${\alpha }_k$ and the product $\sigma$. The latter is the multiplication of the fixed, initial value of the range and Doppler states in the CT model and its first-order derivatives $\dot{\sigma }$. Since all three parameters in the CT model are time-invariant within the unified processing batch, their updates during the iterative process can help lessen the extent of computation and enhance efficiency.

3.2. Optimization Function

The design of the accumulation of OF is a pivotal step in the DP-TBD algorithm detection and tracking process for weak maneuvering targets. This paper devises an OF that enhances the efficiency of energy accumulation in the iterative process. It filters certain clutter and noise, and it utilizes the spatio-temporal relationship between adjacent frames in the same batch, thereby thoroughly extracting information from the original echo.

In this paper, the RD state of the kth frame is denoted as ${\mathbf{c}}_k={\left[r_k{d}_k\right]}^T$. The RD states’ sequence for the same batch is represented as ${\mathbf{C}}_{1:K}=\left\{{\mathbf{c}}_1,{\mathbf{c}}_2,\cdots {\mathbf{c}}_K\right\}$. The OF encapsulates the probability transition from the current frame RD state to the future RD state of the target, rooted in the a posteriori probability density function $p\left({\mathbf{C}}_{1:K}\mid {\mathbf{Z}}_{1:K}\right)$. Maximizing $p\left({\mathrm{C}}_{1:K}\mid {\mathbf{Z}}_{1:K}\right)$ ensures the most accurate represent the next state of the target.

$$p\left({\widehat{\mathbf{C}}}_{1:K}\mid {\mathbf{Z}}_{1:K}\right)=\underset{{\mathbf{C}}_{1:K}}{max}\left\{p\left({\mathbf{C}}_{1:K}\mid {\mathbf{Z}}_{1:K}\right)\right\},$$

(29)

Equation (29) is obtained by Bayesian transformation:

$$p\left({\mathbf{C}}_{1:K+1}\mid {\mathbf{Z}}_{1:K+1}\right)={\displaystyle \frac{p\left({\mathbf{z}}_{k+1}\mid {\mathbf{c}}_{k+1}\right)p\left({\mathbf{c}}_{k+1}\mid {\mathbf{c}}_k\right)}{p\left({\mathbf{z}}_{k+1}\mid {\mathbf{Z}}_{1:k}\right)}}p\left({\mathbf{C}}_{1:K}\mid {\mathbf{Z}}_{1:K}\right).$$

(30)

The OF is characterized by the maximal value of the posterior probability density function:

$$M\left({\mathbf{C}}_{1:K},k\right)=\underset{{\mathbf{C}}_{1:K-1}}{max}\left\{p\left({\mathbf{C}}_{1:K}\mid {\mathbf{Z}}_{1:K}\right)\right\}.$$

(31)

Integrating Equation (30) into Equation (31) produces an iterative relationship for the Optimization Function:

$$M\left({\mathbf{c}}_{k+1},k+1\right)=\underset{{\mathbf{c}}_k}{max}\left\{{\displaystyle \frac{p\left({\mathbf{z}}_{k+1}\mid {\mathbf{c}}_{k+1}\right)p\left({\mathbf{c}}_{k+1}\mid {\mathbf{c}}_k\right)}{p\left({\mathbf{z}}_{k+1}\mid {\mathbf{Z}}_{1:k}\right)}}M\left({\mathbf{c}}_k,k\right)\right\}.$$

(32)

In this paper, the OF was derived by taking the logarithm of both sides of Equation (32) and subsequently applying normalization operations:

$$L\left({\mathbf{c}}_{k+1},k+1\right)=\underset{{\mathbf{c}}_k}{max}\left\{logp\left({\mathbf{c}}_{k+1}\mid {\mathbf{c}}_k\right)+L\left({\mathbf{c}}_k,k\right)\right\}+logp\left({\mathbf{z}}_{k+1}\mid {\mathbf{c}}_{k+1}\right),$$

(33)

where, in Equation (33), the log-likelihood function $logp\left({\mathbf{z}}_{k+1}\mid {\mathbf{c}}_{k+1}\right)$ characterizes the true probability of a target appearing within the cell. The trajectory confirmation function for $k=K-1,\dots ,1$ can be defined as follows:

$$T\left({\mathbf{c}}_{k+1},k+1\right)=arg\underset{{\mathbf{c}}_k}{max}\left\{logp\left({\mathbf{c}}_{k+1}\mid {\mathbf{c}}_k\right)+L\left({\mathbf{c}}_k,k\right)\right\},$$

(34)

where $T\left({\mathbf{c}}_{k+1},k+1\right)$ serves to authenticate the actual track. Additionally, the state transfer loss function $logp\left({\mathbf{z}}_k\mid {\mathbf{c}}_k\right)$ represents the motion characteristics of the target. The transfer of the kth frame state to the $k+1$st frame state relies on the deviation between the actual and ideal trajectories. A larger deviation results in a smaller state transfer loss function and vice versa, and a smaller deviation leads to a larger state transfer loss function.

3.3. The Iterative Energy Accumulation Process

Based on the accurate evolution equations of the CA and CT motion targets and the correction strategy of OF proposed in the above two subsections, combined with the accumulation process of DP-TBD, this section proposes an iterative accumulation method of CA-OF-TBD and CT-OF-TBD target energies in the RD domain. The future positioning of the target is predicted based on the attributes of the CA and CT motion models. Concomitantly, parameters and value functions are subsequently updated in the iterative process, facilitating a more efficient energy accumulation.

3.3.1. The Iterative Energy Accumulation Process of the CA Model

The corresponding location of target in the kth frame can be deduced based on its location within the RD cell:

$$\begin{array}{cc}\hfill r_k& =n_r{\Delta }_r\hfill \end{array}$$

(35)

$$\begin{array}{cc}\hfill d_k& =(n_d-1){\Delta }_d+v_{d\min },\hfill \end{array}$$

(36)

where $v_{d\min }$ is the minimum Doppler velocity.

Our assumption also includes initial parameters defined as ${\dot{\sigma }}_{s0},{\ddot{\sigma }}_{s0},{\stackrel{⃛}{\sigma }}_{s0},$ and ${\stackrel{⃜}{\sigma }}_{s0}$. These can be derived for the kth frame species, following Equation (19), as shown below:

$$\begin{array}{cc}\hfill {\dot{\sigma }}_{sk}& ={\dot{\sigma }}_{s0}+kT{\ddot{\sigma }}_{s0}+{\displaystyle \frac{1}{2}}{k}^2{T}^2{\stackrel{⃛}{\sigma }}_{s0}+{\displaystyle \frac{1}{6}}{k}^3{T}^3{\stackrel{⃜}{\sigma }}_{s0}\hfill \\ \hfill {\ddot{\sigma }}_{sk}& ={\ddot{\sigma }}_{s0}+kT{\stackrel{⃛}{\sigma }}_{s0}+{\displaystyle \frac{1}{2}}{k}^2{T}^2{\stackrel{⃛}{\sigma }}_{s0}\hfill \\ \hfill {\stackrel{⃛}{\sigma }}_{sk}& ={\stackrel{⃛}{\sigma }}_{s0}+kT{\stackrel{⃛}{\sigma }}_{s0}\hfill \\ \hfill {\stackrel{⃜}{\sigma }}_{sk}& ={\stackrel{⃜}{\sigma }}_{s0}.\hfill \end{array}$$

(37)

Utilizing the parameters in Equation (37), we can deduce the target position information for the $k+1$st frame as follows:

$$\begin{array}{cc}\hfill {\widehat{r}}_{k+1}=& \sqrt{\begin{array}{c}{r}_k^2+2kT{\sigma }_{sk}+k^2{T}^2{\dot{\sigma }}_{sk}+{\displaystyle \frac{1}{6}}{k}^3{T}^3{\ddot{\sigma }}_{sk}+{\displaystyle \frac{1}{12}}{k}^4{T}^4{\stackrel{⃛}{\sigma }}_{sk}+{\displaystyle \frac{1}{6}}{k}^5{T}^5{\stackrel{⃜}{\sigma }}_{sk}\hfill \end{array}}\hfill \end{array}$$

(38)

$$\begin{array}{cc}\hfill {\widehat{d}}_{k+1}=& {\displaystyle \frac{1}{r_k}}\left({\sigma }_{sk}+kT{\dot{\sigma }}_{sk}+{\displaystyle \frac{1}{2}}{k}^2{T}^2{\ddot{\sigma }}_{sk}+{\displaystyle \frac{1}{6}}{k}^3{T}^3{\stackrel{⃛}{\sigma }}_{sk}+{\displaystyle \frac{5}{12}}{k}^4{T}^4{\stackrel{⃜}{\sigma }}_{sk}\right).\hfill \end{array}$$

(39)

The RD position information of the target can be transformed into RD cells using the conversion formula as follows:

$$\begin{array}{cc}\hfill {\widehat{n}}_{r_{k+1}}& =⌊{\displaystyle \frac{{\widehat{r}}_{k+1}}{{\Delta }_r}}⌋\hfill \end{array}$$

(40)

$$\begin{array}{cc}\hfill {\widehat{n}}_{d_{k+1}}& =⌊{\displaystyle \frac{{\widehat{d}}_{k+1}-v_{dmin}}{{\Delta }_d}}⌋+1,\hfill \end{array}$$

(41)

where $⌊·⌋$ is the floor function. The target location within the RD cell in the $k+1$st frame is denoted as $\left(n_{r_{k+1}},n_{d_{k+1}}\right)$, while the predicted cell location is depicted by $\left({\widehat{n}}_{r_{k+1}},{\widehat{n}}_{d_{k+1}}\right)$.

The range of observations for the $k+1$st frame is bounded according to the predicted position information obtained by this strategy, so as to obtain the search selectable range of DP-TBD as follows:

$$\begin{array}{cc}\hfill {\widehat{\mathbf{C}}}_{k+1}& =arg\underset{{\mathbf{C}}_{k+1}}{max}p\left({\mathbf{C}}_{k+1}\mid {\mathbf{Z}}_{k+1}\right)\hfill \\ \hfill s.t.& p\left({\mathbf{C}}_{k+1}\mid {\mathbf{Z}}_{k+1}\right)>Thr\hfill \\ \hfill & \left|{\widehat{n}}_{r_{k+1}}-n_{r_{k+1}}\right|⩽{\Delta }_r\hfill \\ \hfill & \left|{\widehat{n}}_{d_{k+1}}-n_{r_{d+1}}\right|⩽{\Delta }_d,\hfill \end{array}$$

(42)

where $Thr$ is the threshold of the probability of the target appearing in the cell $\left(n_{r_{k+1}},n_{d_{k+1}}\right)$.

Algorithm 1 describes the accumulation process of DP-TBD based on the CA model in the RD domain in detail.

Algorithm 1: CA-OF-TBD.

3.3.2. The Iterative Energy Accumulation Process of the CT Model

Similarly, we assumed the initial parameters RD initial value product and its first-order derivative to be ${\sigma }_{s0}$ and ${\dot{\sigma }}_{s0}$, respectively, along with the constant turn rate ${\alpha }_k$. The relationship between the initial RD value product and constant turn rate ${\alpha }_k$ can be inferred from Equations (25) and (26) as follows:

$${\sigma }_{sk}={\displaystyle \frac{{\sigma }_{s0}-\frac{sin\left({\alpha }_{0}kT\right)}{{\alpha }_0}{\dot{\sigma }}_{s0}}{cos\left({\alpha }_{0}kT\right)}},$$

(43)

where ${\alpha }_0$ is the initial value of the constant turn rate.

Upon incorporating ${\sigma }_{s0}$ and ${\dot{\sigma }}_{s0}$ into Equation (26), the ${\dot{\sigma }}_{sk}$ for the kth frame can be derived as follows:

$${\dot{\sigma }}_{sk}=-{\alpha }_{0}sin\left({\alpha }_{0}kT\right){\sigma }_{s0}+cos\left({\alpha }_{0}kT\right){\dot{\sigma }}_{s0}.$$

(44)

By including ${\dot{\sigma }}_{sk}$ in Equations (27) and (28), the estimated RD position states are yielded as follows:

$$\begin{array}{cc}\hfill {\widehat{r}}_{k+1}& =\sqrt{r_k^2+{\displaystyle \frac{2sin\left({\alpha }_{0}T\right)}{{\alpha }_0}}{\sigma }_{sk}+{\displaystyle \frac{2-2cos\left({\alpha }_{0}T\right)}{{\alpha }_0^2}}{\dot{\sigma }}_{sk}}\hfill \end{array}$$

(45)

$$\begin{array}{cc}\hfill {\widehat{d}}_{k+1}& ={\displaystyle \frac{1}{r_k}}\left(cos\left({\alpha }_{0}T\right){\sigma }_{sk}+sin{\displaystyle \frac{sin\left({\alpha }_{0}T\right)}{{\alpha }_0}}{\dot{\sigma }}_{sk}\right).\hfill \end{array}$$

(46)

The conversion to range and Doppler cells was achieved in accordance with Equations (40) and (41).

Algorithm 2 describes the accumulation process of DP-TBD based on the CT model (CT-OF-TBD) in detail.

Algorithm 2: CT-OF-TBD.

3.3.3. Efficient Iterative Process of Energy Accumulation

The iterative energy accumulation process described in the above two subsections targets the state information of the target’s next frame that can be predicted by accurate equation derivations in DP-TBD, which can reduce the range of the target state transfer in the traditional DP-TBD algorithm. The target state transfer range is typically symbolized by the state transfer number q. The values of q often encompass 9, 16, 25, etc., with the specific choice of a larger or smaller state transfer number q adjusted to the actual scenario. For instance, given $q=9,16,$ and 25, the state of the CA target at frame k could potentially be $x_{k+1}$, as expressed in Equation (47) detailed below:

$$\begin{array}{c}\hfill {\mathit{x}}_{k+1}\in \left\{\left[x_0+{\dot{x}}_{0}kT+{\displaystyle \frac{{\ddot{x}}_0{k}^2{T}^2}{2}}+{\displaystyle \frac{1}{6}}{\stackrel{⃛}{x}}_k{T}^3-{\delta }_i,y_0+{\dot{y}}_{0}kT+{\displaystyle \frac{{\ddot{y}}_0{k}^2{T}^2}{2}}+{\displaystyle \frac{1}{6}}{\stackrel{⃛}{y}}_k{T}^3-{\delta }_j\right]\right\},\end{array}$$

(47)

where ${\delta }_i$ and ${\delta }_j$ are the state transfer range of the target in the x and y directions, respectively.

The transfer range of the target is primarily influenced by the target speed, acceleration, the rate of change in acceleration, and the magnitude of the state transfer number. For a more explicit demonstration of the target state transfer process, a graphical representation of the state transfer effects is provided with the state transfer number $q=25$. Figure 4 presents this schematic of conventional DP-TBD when the state transfer number is set to $q=25$.

Owing to the reduced speed and limited maneuverability of the “low, small, and slow” target, its speed transfer range is relatively diminished compared to the distance transfer range. The state transfer process of OF-TBD in the RD domain, as derived from Equations (40) and (41), is illustrated in Figure 5.

As evident from Figure 4 and Figure 5, the algorithm proposed in this study effectively minimizes the search range and enhances the search efficiency without compromising the algorithm accuracy. This strategy circumvents the need for the traditional DP-TBD approach to indiscriminately broaden the search scope for performance improvement, a move that potentially results in a dimensional catastrophe and computational explosion.

4. Simulation Results

In the following section, we perform simulation validation for the proposed CA-OF-TBD and CT-OF-TBD methods. Simulations are conducted for both CA and CT target scenarios utilizing the Monte Carlo method for processing each scenario. The performance of the employed algorithm is quantitatively characterized as follows:

(1)

Target detection probability ($P_d$): When the detected target position of the last frame is within two range units from the true position of the target, the target detection is considered successful, and the calculation of the target detection probability is carried out by a Monte Carlo experiment.

(2)

Range and Doppler dimension root mean square error (RMSE): There is an error between the predicted position of each frame of the algorithm and the true position of the target, and the RMSE of the target position indicates the magnitude of the error between the true position and the predicted position of the target so that we can judge the accuracy of the algorithm using the RMSE between the true position and the predicted position. The formula for RMSE can be expressed as follows:

$$\begin{array}{cc}\hfill RMSE\left({\theta }_k\right)& =\sqrt{{\displaystyle \frac{1}{N}}\sum _{i=1}^N{\left({\widehat{\theta }}_k-{\theta }_k\right)}^2},\hfill \end{array}$$

(48)

where N denotes the quantity of Monte Carlo experiments, and ${\widehat{\theta }}_k$ symbolizes the forecasted values of the following parameters: range, Doppler, and correlation coefficients. Conversely, ${\theta }_k$ stands for the observed values of these parameters.

(3)

Correlation Coefficient root mean square error (RMSE): The equations for the CA, as well as CT, motion models in the range–Doppler domain can be expanded by means of a Taylor series to obtain a direct relationship between the relevant parameters in the target range–Doppler approximation domain, obtaining the series shown below:

$$\begin{array}{cc}\hfill r_k& =r_0+d_{0}kT+{\displaystyle \frac{{\dot{\sigma }}_0-p_{d,0}^2}{2r_0}}{k}^2{T}^2+\dots +{\displaystyle \frac{r_k^{\left(n\right)}}{n!}}{|}_{k=0}{k}^n{T}^n\hfill \end{array}$$

(49)

$$\begin{array}{cc}\hfill d_k& =d_0+{\displaystyle \frac{{\dot{\sigma }}_0-p_{d,0}^2}{r_0}}kT+{\displaystyle \frac{d_k^{\left(2\right)}}{2}}{{|}_{k=0}{k}^2{T}^2+\dots +{\displaystyle \frac{d_k^{\left(n\right)}}{n!}}|}_{k=0}{k}^n{T}^n,\hfill \end{array}$$

(50)

where $r_k^{\left(n\right)}$ and $d_k^{\left(n\right)}$ represent the nth-order derivatives of the target range and Doppler dimensions of the CA and CT models, respectively. Furthermore, these variables elucidate the correlation between the coefficients and their approximations. Correspondingly, the RMSE of these correlation coefficients was also computed using Equation (48).

The experiments utilize the Monte Carlo testing approach with $N=500$. The radar simulation range and Doppler range during experimentation were set at (0, 400 km) and (−200 m/s, 200 m/s), respectively, with the range and Doppler resolution of the radar fixed at ${\Delta }_r$ = 1 km and ${\Delta }_d$ = 1 m/s, respectively. We defined 20 s as the interval between adjacent frames to conspicuously observe the target motion status and trend. The detection probability $P_d$ is described as the scenario when the detected target position is within 2 cells of the actual target location. The detection threshold adheres to the Neyman–Pearson Criterion at a false alarm rate of $P_{fa}={10}^{-2}$. Subsequently, contrasting the process of energy accumulation across diverse frames via CA target simulation was executed, with an assumption that with SNR = 8 dB, we can attain Figure 6.

As inferred from Figure 6, the detection probability upon the accumulation of six frames appeared almost equivalent to 1. Further accumulation only augments computation, with a diminishing effect on detection efficiency. As such, we selected the number of accumulated frames K as six for this study.

The accumulation performance of the proposed CA-OF-TBD and CT-OF-TBD methods within this study hinges on the precision of the model derivation. Therefore, this study chose the modeling algorithm CS-DP-TBD [35,36], which is more adaptable, for comparison. Recognizing that prevailing DP-TBD accumulation processes primarily concentrate on the Cartesian coordinate system, CA-PS-TBD and CT-PS-TBD [23] emerge as congruent options for comparison with VF-TBD and DP-TBD in similar motion model cases. Notably, this study also involved a search for the target prediction state. By calibrating the search interval, DP-TBD can mitigate potential losses in accumulation. Consequently, the DP-TBD technique using six-dimensional (6D) state vectors—accounting for the position, velocity, and acceleration in each direction—to specify Region Expansion (RE-DP-TBD_6D) [33] based on the maximum feasible velocity or acceleration was utilized for comparison. Lastly, we drew a comparison with the single-frame TBD (ST-TBD) algorithm, where the supplementary experiments are presented in Figure 6, further comparing the probability of detecting targets in scenarios with varying cumulative frames.

4.1. CA Target Scenarios

In this paper, the algorithm was validated for the detection performance of the radars’ far and near regions, weak and low-speed targets, and the corresponding RMSE. Therefore, the CA target scenario is established with two CA motion targets, labeled, respectively, as T1 and T2 within the Cartesian coordinate system. The initial states for these motion targets are, respectively, set as follows: (260 km, 32 m/s, 1.2 m/s², 150 km, −12 m/s, and −0.5 m/s²) and (25 km, −12 m/s, −0.5 m/s², 35 km, 32 m/s, and 1.2 m/s²). The target motion for this experiment within a Cartesian coordinate system is portrayed in Figure 7 and Figure 8. Figure 7 corresponds to the far region target, while its counterpart Figure 8 delineates the near region target.

The target motion equations were transcribed from the Cartesian coordinate system to the RD domain utilizing the formula described in Section 3. Consequently, Figure 9 and Figure 10 depict the motion curves of the far and near region targets when transformed into the RD domain. In these figures, we have not used the corresponding resolution cell for superior visualization. In order to balance the performance of the algorithm and the amount of computation required in the existing TBD research, six frames of data are usually taken as a processing batch. Therefore, a processing batch was set for six frames of data i.e., K = 6 in this paper.

4.1.1. CA Target Scenario Detection Probabilities

This subsection validates the superiority of the algorithm proposed herein through a comparative analysis of the detection probabilities across diverse TBD algorithms operating under varying SNR conditions on the same target. Two targets, T1 and T2, were selected to facilitate this comparison. Figure 11 and Figure 12 reveal a comparative analysis showcasing the detection probability differences across multiple TBD algorithms for Targets T1 and T2 at different SNR levels.

All the algorithms in Figure 11 and Figure 12 demonstrated increasing detection probabilities of Target T1 and T2 as the SNR increased. However, the ST-DP-TBD algorithm, illustrated in Figure 11 and Figure 12, exhibited the weakest target detection performance due to its insufficiency in superimposing multi-frame information. Moreover, a more accurate model description enhances the detection probability; this is corroborated by CS-DP-TBD’s superior performance over the RE-DP-TBD_6D algorithm, as indicated in the figures. Figure 11 and Figure 12 further clarify that the performance of the CA-PS-TBD algorithm was considerably improved in the low-SNR stage compared to CS-DP-TBD. This improvement benefits from the more precise target position prediction and energy accumulation of the CA-PS-TBD algorithm, which notably elevates the detection probability of CA-PS-TBD at an SNR of 0–5 dB. Concurrently, the CA-OF-TBD algorithm proposed in this paper can also predict the target more accurately by addressing the target spreading problem through the OF, leading to a higher detection rate compared to CA-PS-TBD.

Figure 12 illustrates the detection probability of Target T2, as per the aforementioned algorithm. Compared to Target T1 (situated farther from the sensor), the relatively close location of Target T2 makes detection more challenging at the same SNR level. Proximity to the sensor aggravates a mismatch of the target motion model, biasing the target energy accumulation and complicating detection. Even though the DP algorithm can expand the search range, excessive noise levels can overshadow the actual target. However, the accurate target motion evolution equation employed in this paper addresses the target motion mismatch problem effectively, facilitating improved detection at lower SNR levels.

4.1.2. The CA Target Scenarios Root Mean Square Error (RMSE)

In this subsection, we mainly show the variation process of the range and Doppler values, as well as the RMSE of the correlation coefficient, of the CA-OF-TBD algorithm for different values of the SNR through Figure 13, Figure 14, Figure 15, Figure 16, Figure 17, Figure 18, Figure 19, Figure 20, Figure 21 and Figure 22.

(1)

Range and Doppler dimension RMSE: The accuracy deviation for detecting and tracking Targets T1 and T2 using various TBD algorithms is exhibited in Figure 13, Figure 14, Figure 15 and Figure 16. It is evident that the RMSE values for the range and Doppler dimensions in the CA-PS-TBD algorithm proposed herein were lower than other TBD algorithms. The CA-OF-TBD algorithm displayed smaller RMSE values in these dimensions compared to preceding DP-TBD algorithms as it effectively mitigated the detecting and tracking challenges posed by the inaccurate path prediction and energy accumulation. This was achieved by deploying the precise target motion evolution equations constituted by the CA-OF-TBD algorithm, resulting in smaller range and Doppler RMSE values compared to earlier DP-TBD algorithms. Moreover, the proposed corrective measure for the target spreading effect using OF enhances the trajectory prediction accuracy for the targets, thereby realizing lower RMSE values.

(2)

Correlation Coefficient RMSE: The parameters involved in the CA target mainly include ${\dot{\sigma }}_0,{\ddot{\sigma }}_0,{\stackrel{⃛}{\sigma }}_0,$ and ${\stackrel{⃛}{\sigma }}_0$, and since these parameters are not involved in other algorithms, the algorithms compared in this section are CA-PS-TBD. Since these parameters can effectively improve the accuracy of the target detection and tracking, lowering the RMSE of these parameters is to lay a stronger foundation for the accurate detection and tracking of the target. Figure 17, Figure 18, Figure 19, Figure 20, Figure 21, Figure 22, Figure 23 and Figure 24 show that the RMSE of the correlation coefficient of the CA-PS-TBD algorithm was lower than that of the CA-PS-TBD algorithm for different SNR cases. Meanwhile, according to Equation (48), it can be concluded that, when the target is located in the far field of the sensor, it will have a greater influence on the parameters ${\dot{\sigma }}_0,{\ddot{\sigma }}_0,{\stackrel{⃛}{\sigma }}_0,$ and ${\stackrel{⃛}{\sigma }}_0$; thus, there will be a situation where the estimation error of the parameters ${\dot{\sigma }}_0,{\ddot{\sigma }}_0,{\stackrel{⃛}{\sigma }}_0,$ and ${\stackrel{⃛}{\sigma }}_0$ will be significantly larger than that of the incoming field.

4.2. CT Target Scenarios

In the present study, we defined Targets T1 and T2 within the CT scenario. They possess the initial target states of (260 km, 32 m/s, 150 km, and −12 m/s) and (25 km, −12 m/s, 35 km, and 32 m/s), respectively. Similar to the CA model, the motion trajectories of the two targets in the far- and near-region Cartesian coordinate system are represented by Figure 23 and Figure 24, respectively. The motion trajectories of the two targets in the RD domain are also depicted in Figure 25 and Figure 26. In these figures, we have not used the corresponding resolution cell for superior visualization.

A turning rate, $\lambda =0.2$ deg./s, was similarly assigned to these targets. Additionally, this section processes the data in batches of six frames.

4.2.1. CT Target Scenario Detection Probabilities

Mirroring the CA scenario, this section primarily illustrates the detection probability of a CT target, as discerned by various algorithms under different SNRs. The detection probability for the CT targets also increases with increasing SNRs. It is clear from Figure 27 and Figure 28 that the CT-OF-TBD algorithm proposed in this paper surpasses other algorithms in target detection and tracking within CT scenarios. Firstly, the multi-frame accumulation algorithm can utilize the information of multi-frame data, so the detection probability will be greatly improved compared to ST-TBD. In addition, accurate motion model derivation can be more accurate for target detection, whether it is the detection probability comparison between the CS-DP-TBD and RE-DP-TBD_6D algorithms or the CA-OF-TBD and CA-PS-TBD algorithms proposed in this paper.

4.2.2. CT Target Scenarios Root Mean Square Error (RMSE)

The experiments conducted in this subsection illustrate the variation in the RMSE of the CT target range-Doppler values and their associated parameters across different SNRs.

(1)

The Range and Doppler dimensions RMSE:Figure 29, Figure 30, Figure 31 and Figure 32 illustrate a decrease in the range–Doppler dimension RMSE in conjunction with an increase in SNR. Concurrently, the CT-OF-TBD algorithm—introduced in this paper—significantly minimizes the RMSE of the range–Doppler dimension, and this is attributed to its precise target trajectory prediction.

(2)

Correlation Coefficient RMSE: The key parameters of the CT targets were primarily ${\dot{\sigma }}_0$ and $\alpha$, with their variation noted as in accordance with the SNR depicted in Figure 33, Figure 34, Figure 35, Figure 36, Figure 37 and Figure 38. Similar to the CA scenarios, the CT-OF-TBD algorithm not only curbs the RMSE of the range–Doppler dimension more efficiently, but it also furnishes a higher-accuracy parameter estimation.

Both the CT-OF-TBD algorithm and the CA-OF-TBD algorithm proposed in this paper can effectively detect and track moving targets in the CA and CT scene. Firstly, accurate range–Doppler dimension predictions of the moving target are carried out through the accurate target motion evolution equation, and then the broadening effect of the target is effectively suppressed through the OF, so that the target energy can be effectively accumulated and more efficient detection and tracking of the target can be realized.

5. Conclusions

This paper presents the construction of precise target motion evolution equations for weak targets in CA and CT scenes employing a Cartesian coordinate system. The recursive estimation of the parameters relevant to the target motion process facilitates a more accurate prediction of the range–Doppler value, enabling more efficient accumulations of target energy, thus enhancing the target detection rate. Additionally, an OF was deployed effectively to mitigate the broadening effect attributed to target energy accumulation, which bolstered the target resolution and achieved more accurate target detection and tracking. The deployment of the precise evolution equations elucidated in this paper significantly minimizes the computational load of the DP algorithm. The efficiencies of both the CA-OF-TBD algorithm and CT-OF-TBD algorithm were validated through simulation experiments involving CA and CT targets. Future work can potentially refine the algorithm further by incorporating increasingly accurate target motion evolution equations and a higher-order OF. Enabling the generalization for detecting unknown motion weak targets constitutes a significant area of research. Moreover, the present setup portrays a scenario that is comparatively straightforward, leading to virtually negligible leakage detection following energy accumulation. Subsequent research will entail the adaptation of the algorithm to a broader setting, given that the target motion in real-world scenarios evinces complexity and variability. The encompassed scenarios include, but are not limited to, multiple overlapping targets, intersecting trajectories, and variances in the SNRs relative to the target and additional scenarios. In applying tools such as the Optimal Sub-Pattern Assignment (OSPA) metric, the aim is to scrutinize instances of target omissions closely.