CN113093174B - Tracking-before-detection method for multi-targets with weak fluctuations in radar based on PHD filtering - Google Patents

Abstract

The invention discloses a PHD filtering radar fluctuation weak multi-target-based pre-detection tracking method, which solves the problems of target detection and tracking under amplitude fluctuation, researches three fluctuation target models of powerling 0,1 and 3, establishes two tracking models of complex likelihood and amplitude likelihood under a PHD-TBD algorithm, wherein the complex likelihood method makes up the defect that the amplitude likelihood only considers measured amplitude information and ignores phase information in the calculation process, thereby better utilizing original information of a target. Under the condition of fluctuation of target amplitude, the complex likelihood is superior to the amplitude likelihood in estimation performance of target position and number, and the calculation efficiency is higher. At low signal-to-noise ratios, complex likelihood can still effectively detect and track an unknown number of weak targets.

Description

Tracking-before-detection method for weak multiple targets based on PHD filtering radar fluctuations

Technical Field

The invention relates to the technical field of radar fluctuation weak multi-target detection and tracking, and in particular to a pre-detection tracking method for radar fluctuation weak multi-target based on PHD filtering.

Background Art

With the rapid development of modern electronic information technology, radar target detection technology faces the threat of aircraft stealth technology. The development of stealth technology has reduced the radar cross section (RCS) of stealth targets, and the target reflected echo signal is very weak, with a low echo signal-to-noise ratio (SNR). The traditional detection and tracking method for moving targets is detection-before-tracking (DBT). In this method, a threshold is set for each frame of measurement data to determine whether the target exists, and the target trajectory is obtained through the tracking algorithm; however, when the signal-to-noise ratio of the echo signal is low, the echo of a weak target is usually lower than the threshold, resulting in missed detection, and it is difficult to extract the target trajectory using a single frame of measurement data; if the threshold is lowered, a large number of false alarms will be generated, and the target trajectory cannot be maintained.

To solve the above problems, one method is to use the Track Before Detection (TBD) algorithm. This algorithm processes multiple frames of data jointly based on the continuity of target motion in space and the temporal correlation of consecutive frames of target echo data, and realizes target detection and tracking through multi-frame energy accumulation. Traditional TBD methods include dynamic programming, Hough transform, etc. However, these methods have large computational complexity and high algorithm complexity, and are only applicable to linear Gaussian models.

The particle filter TBD (PF-TBD) algorithm under the Bayesian framework can solve nonlinear and non-Gaussian problems, and has been rapidly developed in the field of weak target detection and tracking. The limitation of the PF-TBD method is that it does not model the appearance (birth) and disappearance (death) of the target. Therefore, when dealing with an unknown and changing number of targets, the algorithm complexity increases sharply, and the filter performance has certain limitations.

Summary of the invention

The purpose of the present invention is to provide a pre-detection tracking method for weak multiple targets based on PHD filtering radar fluctuations, aiming to solve the technical problem in the prior art that weak multiple target fluctuations cannot be effectively tracked.

To achieve the above object, the present invention adopts a pre-detection tracking method for weak multi-targets based on PHD filtering radar fluctuations, which comprises the following steps:

S1: Initialize system parameters and read the original measurement data at the k-1th moment and the kth moment in the radar receiver;

S2: Divide the scene and use the original measurement data at time k-1 to adapt the new target in each small scene;

S3: For the complex measurement and square measurement data obtained at time k, the complex likelihood and amplitude likelihood of the three types of amplitude fluctuations are calculated respectively, and the SMC implementation of PHD filtering under amplitude fluctuation is given;

S4: Extract target state and target number based on posterior information;

S5: Let k=k+1, and determine whether k＞K. If so, the algorithm ends; otherwise, return to step 2.

The system parameters include:

对于角度，在雷达接受端考虑N_a天线的线性相控阵，间隔为

其中λ为载波频率的波长，角度分辨率为

Sampling interval T, current time k, total target motion time K, radar scans area [r _min ,r _max ]×[θ _min ,θ _max ] in polar coordinates, radar receives measurement data Z _k and Z _k-1 in the tracking scene, consider the range and azimuth surveillance radar covering the defined area in polar coordinates, for range, assume that the transmitted pulse is a linear frequency modulation signal with bandwidth B and duration T _ε , speed of light c, range resolution unit

For the angle, consider the linear phased array of Na _antennas at the radar receiving end, with a spacing of

Where λ is the wavelength of the carrier frequency and the angular resolution is

Divide the scene, and use the original measurement data at time k-1 to adapt the new target in each small scene:

个场景，考虑在每个场景内生成N_filter个粒子，则整个场景内的粒子数为n₁×n₂×N_filter；If N _r ×N _θ is large, the scene N _r ×N _θ is divided into

scenes, consider generating N _filter particles in each scene, then the number of particles in the whole scene is n ₁ ×n ₂ ×N _filter ;

Set the measurement cutoff threshold _Thcutoff , sort the measurements in each scene in descending order, select the measurements with strength higher than the _threshold Thcutoff, and consider the measurements lower than _Thcutoff as false positive measurements. The position information of each measurement is (n _r ,n _θ );

Convert each measured position information (n _r ,n _θ ) into the target position in the plane rectangular coordinate system, recorded as (z _x ,z _y );

Generate particles near each (z _x ,z _y ), calculate the likelihood of each particle in the scene, then resample to select particles with higher weights in the scene, and normalize the weights of the selected particles.

For the complex measurement and square measurement data obtained at time k, the complex likelihood and amplitude likelihood of the three types of amplitude fluctuations are calculated respectively, and the steps of SMC implementation of PHD filtering under amplitude fluctuation are given:

Calculate the magnitude likelihood when the target is not present;

Calculate the amplitude likelihood under the amplitude fluctuation type of swerling 0, 1, and 3 respectively;

Calculate the complex likelihood when the target does not exist;

Calculate the complex likelihood for amplitude fluctuation types swerling 0, 1, and 3 respectively.

At the k-1 moment:

的粒子

代表PHD的后验密度，预测和更新这一组带权重的粒子得到

Use a set of weights

Particles

Representing the posterior density of PHD, predicting and updating this set of weighted particles to obtain

For complex and square measurements, the likelihood of fluctuations of different amplitudes is different. During target tracking, the target amplitude will affect the strength of the target echo. Targets with high echo intensity often suppress targets with low intensity. Consider assigning the likelihood above the threshold as the highest likelihood.

When the target amplitude fluctuation type is swerling 3, the particle weight is calculated and judged. If the particle weight sum is infinite, the predicted particle state is kept unchanged, the particle is not updated, and the particle weight is assumed to be 0.

In the step of extracting the target state and target number based on the posterior information:

Calculate the target number and resample the samples;

和是否大于0，若大于0，则用K-means方法对重采样后的粒子进行聚类，删除权重和低于门限Th_k-means的粒子群，将权重和高于门限Th_k-means的粒子群视为估计的目标状态，若小于0，则认为此时场景内没有目标。Determine the updated particle weight

Whether the sum is greater than 0, if it is greater than 0, the resampled particles are clustered using the K-means method, and the particle groups with weights and sums lower than the threshold Th _k-means are deleted, and the particle groups with weights and sums higher than the threshold Th _k-means are regarded as the estimated target state. If it is less than 0, it is considered that there is no target in the scene at this time.

The beneficial effects of the present invention are: it can solve the detection and tracking problems of weak multi-target fluctuations. In order to solve the problem of target generation under the condition of unknown prior distribution information of the new target state, an adaptive target generation algorithm based on measurement likelihood under scene division is proposed. Two calculation methods, amplitude likelihood and complex likelihood, are given, and PHD-TBD multi-target estimation with amplitude fluctuation types of swerling 0, 1, 3 is successfully realized. Compared with the amplitude likelihood, the complex likelihood makes up for the defect that the amplitude likelihood only considers the measured amplitude information and ignores the phase information during the calculation process, thereby making better use of the original information of the target. In the method proposed in the present invention, compared with the amplitude likelihood, the former is superior to the latter in the estimation performance of the target position and number, and has higher calculation efficiency. Under low signal-to-noise ratio, the complex likelihood can still effectively detect and track an unknown number of weak targets.

BRIEF DESCRIPTION OF THE DRAWINGS

In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the drawings required for use in the embodiments or the description of the prior art will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For ordinary technicians in this field, other drawings can be obtained based on these drawings without paying creative work.

FIG. 1 is a flow chart of the steps of the method for tracking before detection of weak multiple targets based on PHD filtering radar fluctuations of the present invention.

FIG2 is a schematic flow chart of a method for tracking before detection of weak multiple targets based on PHD filtering radar fluctuations according to the present invention.

FIG. 3 is a real trajectory of target motion in the simulation scene of the present invention.

FIG. 4 is an OSPA diagram of the complex likelihood and amplitude likelihood algorithms under swerling 0 invented.

FIG5 is a comparison diagram of the target potential of the complex likelihood and amplitude likelihood algorithms under swerling 0 of the present invention.

FIG6 is an OSPA diagram of the complex likelihood and amplitude likelihood algorithms under swerling 1 of the present invention.

FIG. 7 is a comparison diagram of the target potential of the complex likelihood and amplitude likelihood algorithms under swerling 1 of the present invention.

FIG8 is an OSPA diagram of the complex likelihood and amplitude likelihood algorithms under swerling 3 of the present invention.

FIG9 is a comparison diagram of target potentials of complex likelihood and magnitude likelihood algorithms under swerling 3 of the present invention.

FIG. 10 is an OSPA diagram of complex likelihood and square likelihood of different signal-to-noise ratios under swerling 0 of the present invention.

FIG. 11 is a diagram showing the potential of a target at different signal-to-noise ratios at different times under swerling 0 of the present invention.

FIG. 12 is an OSPA diagram of complex likelihood and square likelihood at different signal-to-noise ratios under swerling 1 of the present invention.

FIG. 13 is a diagram showing the potential of a target at different signal-to-noise ratios at different times under swerling 1 of the present invention.

FIG. 14 is an OSPA diagram of complex likelihood and square likelihood at different signal-to-noise ratios under swerling 3 of the present invention.

FIG. 15 is a diagram showing the potential of a target at different signal-to-noise ratios at different times under swerling 3 of the present invention.

DETAILED DESCRIPTION

Embodiments of the present invention are described in detail below, examples of which are shown in the accompanying drawings, wherein the same or similar reference numerals throughout represent the same or similar elements or elements having the same or similar functions. The embodiments described below with reference to the accompanying drawings are exemplary and are intended to be used to explain the present invention, and should not be construed as limiting the present invention.

In the description of the present invention, it should be understood that the terms "length", "width", "up", "down", "front", "back", "left", "right", "vertical", "horizontal", "top", "bottom", "inside", "outside" and the like indicate positions or positional relationships based on the positions or positional relationships shown in the accompanying drawings, and are only for the convenience of describing the present invention and simplifying the description, rather than indicating or implying that the device or element referred to must have a specific orientation, be constructed and operated in a specific orientation, and therefore cannot be understood as a limitation on the present invention. In addition, in the description of the present invention, "plurality" means two or more, unless otherwise clearly and specifically defined.

Referring to FIG. 1 and FIG. 2 , the present invention provides a method for tracking weak multi-targets before detection based on PHD filtering radar fluctuations, comprising the following steps:

S1: Initialize system parameters and read the original measurement data at the k-1th moment and the kth moment in the radar receiver;

S2: Divide the scene and use the original measurement data at time k-1 to adapt the new target in each small scene;

S3: For the complex measurement and square measurement data obtained at time k, the complex likelihood and amplitude likelihood of the three types of amplitude fluctuations are calculated respectively, and the SMC implementation of PHD filtering under amplitude fluctuation is given;

S4: Extract target state and target number based on posterior information;

S5: Let k=k+1, and determine whether k＞K. If so, the algorithm ends; otherwise, return to step 2.

对于角度，在雷达接受端考虑N_a天线的线性相控阵，间隔为

其中λ为载波频率的波长，角度分辨率为

Specifically, initialize the system parameters, which include the sampling interval T, the current time k, the total target movement time K, the radar scanning area [r _min ,r _max ]×[θ _min ,θ _max ] in polar coordinates, the radar receives the measurement data Z _k and Z _k-1 in the tracking scene, consider the range and azimuth monitoring radar covering the defined area in polar coordinates, for the range, assume that the transmitted pulse is a linear frequency modulation signal with bandwidth B and duration T _ε , the speed of light c, and the range resolution unit

For the angle, consider the linear phased array of Na _antennas at the radar receiving end, with a spacing of

Where λ is the wavelength of the carrier frequency and the angular resolution is

Where, k is initialized to 2, b _k|k-1 (x _k |x _k-1 ) and γ _k (x _k ) represent the PHD of the target at time k, κ _k (z) = λ _k ·c(z) is the false alarm intensity, λ _k is the average number of false alarms, and c(z) is the false alarm distribution; _PD (x) is the detection probability of the target state, and g _k (z|x) represents the likelihood of the target generating the measurement. In multi-target tracking, the integral of the updated PHD multi-target state function D _k|k is the expected value of the number of targets at time k.

Specifically, the measurement value received by the radar sensor consists of the range and azimuth after range matching filtering and adaptive beamforming. Given the i-th target state x _k,i at time k, the measurement value z _k is given by the following nonlinear equation:

和ρ_k,i分别表示第i个目标复振幅的相位和幅度。假设所有相位

和幅度ρ_k,1:N_k相互独立，且独立于n_k和x_k.1:Nk。相位

是未知且在每一时刻均匀分布在[0,2π)上；对于幅度ρ_k,i有：

其中

是未知的静态参数。k时刻雷达量测有两种表达形式：一种是相干积累的复量测z_k，另一种是非相干积累的平方量测为|z_k|²。Where h(x _ki ) represents the fuzzy function of the i-th target centered at the target position (x _k,i ,y _k,i ). For simplicity, h(x _ki ) is denoted as h _k,i . _{n k} is the measurement noise with mean 0 and covariance complex Gaussian Γ;

and ρ _k,i represent the phase and magnitude of the i-th target complex amplitude, respectively. Assume that all phases

and amplitude ρ _k,1 :N _k are independent of each other and of n _k and x _k.1:Nk . Phase

is unknown and uniformly distributed on [0,2π) at every moment; for the amplitude ρ _k,i :

in

is an unknown static parameter. The radar measurement at time k can be expressed in two forms: one is the complex measurement z _k of coherent accumulation, and the other is the square measurement of incoherent accumulation |z _k | ² .

The fuzzy function for distance matching filtering is:

l∈[0,N_r-1]，

Where, |τ ^l |＝2(r _k -r _l )/c,

l∈[0,N _r -1],

The azimuth ambiguity function of adaptive beamforming is:

以及

m∈[0,N_θ-1]，

Among them, Φ ^m =π[sin(θ _k )-sin(θ _m )],

as well as

m∈[0,N _θ -1],

The overall ambiguity function in the range-azimuth unit (l,m):

The size of h(x _k ) is N _c =N _r ×N _θ , that is:

In the traditional radar target detection and tracking problem, each frame of echo data is first thresholded to form point information, and then the point information exceeding the threshold is associated, filtered, and finally the target track is obtained to achieve target tracking. This detection-first-then-tracking method is suitable for situations where the signal-to-noise ratio is high or the target echo amplitude is large. The intensity of the target is much higher than the intensity of the clutter. By setting a larger threshold, the target can be separated from the clutter. In the case of low signal-to-noise ratio or weak target echo signal, the target echo is buried in the noise clutter. A single-frame over-threshold detection method is used. If the threshold is too high, the target will be missed, and if the threshold is too low, the false alarm rate will increase, and the target track cannot be maintained. The TBD method can improve the detection and tracking of radar weak targets from the perspective of signal processing. Its basic idea is to process the raw data that has not been thresholded, achieve target capability accumulation through multi-frame tracking, and fully explore the target information in the echo.

Furthermore, the PHD filtering algorithm based on RFS simplifies the multi-target state space into a single-target state space by transferring the first-order moment of the multi-target posterior probability density distribution, which greatly reduces the computational complexity and makes it possible to perform practical operations. The prediction and update equations of the PHD filter are as follows:

For the update equation D _k|k , _PD (x) is the detection probability associated with the target state; for the TBD algorithm, there is no detection process before the update step, and it is assumed that the measurement value contains all the target information; therefore, _PD (x) ≡ 1, so the update equation becomes:

For a given target state, the intensity information in each unit is independent of each other, so the multi-target posterior probability density p(z _k |X _k ) can be expressed as the product of marginal probability density functions:

个场景；考虑在每个场景内生成N_filter个粒子，则整个场景内的粒子数为n₁×n₂×N_filter；Specifically, in step 2, if N _r ×N _θ is large, the scene N _r ×N _θ is divided into

scenes; consider generating N _filter particles in each scene, then the number of particles in the entire scene is n ₁ ×n ₂ ×N _filter ;

Set the measurement cutoff threshold _Thcutoff , sort the measurements in each scene in descending order, select the measurements with strength higher than the _threshold Thcutoff, and consider the measurements lower than _Thcutoff as false positive measurements. The position information of each measurement is (n _r ,n _θ );

和

转换为距离和方位，再根据公式z_x＝rcosθ,z_y＝rsinθ把极坐标下的位置转化为平面直角坐标下的位置信息；Use the formula to convert the distance and angle into

and

Convert to distance and azimuth, and then convert the polar coordinates to the rectangular coordinates according to the formula z _x = rcosθ, z _y = rsinθ;

其中U表示均匀分布，即每个量测附近会生成一个粒子，记为

Generate the position information ( _xi , _yi ) of the target state near each ( ^zx , ^zy ), and v = U ( _vmin , _vmax ) generates the velocity information of the target state

Where U represents uniform distribution, that is, a particle will be generated near each measurement, denoted as

For each particle ^xi, the likelihood is calculated in the scene where it is located, and then resampling is performed to select particles with higher particle weights in the scene; and the weights of the selected particles are normalized.

Specifically, in step 3, the magnitude likelihood p(|z _k | ² ) when the target does not exist is calculated:

Calculate the amplitude likelihood for amplitude fluctuation types swerling 0, 1, and 3 respectively:

in,

Calculate the complex likelihood p(z _k ) when the target does not exist:

Calculate the complex likelihood for amplitude fluctuation types swerling 0, 1, and 3 respectively:

where I ₀ (·) is the modified Bessel function of the first kind.

的粒子

代表PHD的后验密度，即：Specifically, for the k-1 moment, use a set of weighted

Particles

represents the posterior density of PHD, that is:

Predicted particles:

和b_k(·|z_k)是建议密度，L_k-1是k-1时刻幸存的目标的粒子数，J_k是在k时刻新生目标的粒子数。则预测的强度D_k-1|k为：in,

and b _k (·|z _k ) is the proposed density, L _k-1 is the number of particles of the target surviving at time k-1, and J _k is the number of particles of the newly born target at time k. Then the predicted intensity D _k-1|k is:

in:

是状态为

的目标从时间k-1到时间k的存活概率；

是状态为x_k-1的衍生概率密度；

是新生概率密度。in,

The status is

The survival probability of the target from time k-1 to time k;

is the derivative probability density of state x _k-1 ;

is the birth probability density.

The update of PHD can be expressed in SMC as:

in:

Specifically, for complex measurement and square measurement, the likelihood of fluctuations of different amplitudes is different. During target tracking, the target amplitude will affect the strength of the target echo. Targets with high echo intensity often suppress targets with low intensity. Consider assigning the likelihood above the threshold as the highest likelihood.

Note that, from the measurement equation, the target echo signal is a sinc function, and the likelihood ratio of the target location will be very high at this time; when the target amplitude fluctuates, the strength of the echo signal varies. In order to solve the multi-target tracking problem in this case and solve the matching tracking of targets of different strengths at the same time, the likelihood ratio greater than the threshold Th _L is assigned the highest likelihood, that is:

L _sw (L _sw ≥ Th _L ) = max (L _sw ).

Specifically, when the target amplitude fluctuation type is swerling 3, the sum of the particle weights may be infinite, which will make it impossible to track the target at subsequent moments. To address this phenomenon, a judgment is made after the particle weight calculation. If the sum of the particle weights is infinite, the predicted particle state remains unchanged, the particle is not updated, and the particle weight is assumed to be 0.

重采样样本

得到

判断更新后粒子权重

和是否大于0，若大于0，则利用K-means方法对重采样后的粒子进行聚类，删除权重和低于门限Th_k-means的粒子群；将权重和高于门限Th_k-means的粒子群视为估计的目标状态，若小于0，则认为此时场景内没有目标。Specifically, step 4 includes calculating the target number

Resample samples

get

Determine the updated particle weight

Whether the sum is greater than 0, if it is greater than 0, the K-means method is used to cluster the resampled particles, and the particle groups with weights and values lower than the threshold Th _k-means are deleted; the particle groups with weights and values higher than the threshold Th _k-means are regarded as the estimated target state. If it is less than 0, it is considered that there is no target in the scene at this time.

Specific embodiment one:

This embodiment uses MATLAB software version 2014 (a) to perform simulation experiments.

Please refer to Figure 2, set the five target motions shown, consider a two-dimensional motion scene, and the state of each target is defined as

分别是笛卡尔坐标系下目标的位置和速度。where (x _k ,y _k ) and

are the position and velocity of the target in the Cartesian coordinate system.

The position (x _k , y _k ) is in the polar coordinate p = [r _min , r _max ] × [θ _min , θ _max ] scene, where r _min , r _max and θ _min , θ _max are the minimum and maximum target range and orientation, respectively;

在区域：speed

In the region:

where _vmin and _vmax are the minimum and maximum target velocities respectively.

The target state evolves according to uniform linear motion:

x _k =Fx _k-1 +v _k

The sensor scans at a time interval of T = 1s to receive 100 frames of images, where

The process noise _vk follows a Gaussian distribution, and its covariance is:

The standard deviation of the noise is σ _v = 5 m, and the probability of the target surviving is p _s,k (x) = 0.98.

给出。For the target simulation, v _min = 200 m/s, v _max = 800 m/s, and the signal-to-noise ratio is given by

Given.

B＝150KHz，T_e＝6.67×10^-5s，N_a＝50，λ＝3cm，c＝3×10⁸m/s，Δ_r＝500m，Δ_θ＝1.45°。For the simulation of radar measurements, r _min = 100 km, r _max = 120 km, θ _min = -75°, θ _max = 75°, N _r = 300, N _θ = 100, σ ² = 0.5, and the noise covariance is

B=150KHz, _Te =6.67×10 ^-5 s, _Na =50, λ=3cm, c=3×10 ⁸ m/s, Δ _r =500m, Δ _θ =1.45°.

Appendix: The trajectory of the five targets in this specific embodiment

The optimal sub-pattern assignment distance (OSPA) is used to evaluate the performance of the algorithm. The OSPA metric can evaluate the target number estimation error and target position estimation error of the multi-target filter. Given two finite sets X = {x ₁ , x ₂ , ... x _m } and Y = {y ₁ , y ₂ , ... y _n }, OSPA is defined as follows:

c＞0为截断参数，用于惩罚目标个数的估计偏差，p为阶数，用于惩罚多目标状态估计偏差。在本文仿真实验中，设置p＝3，c＝1000。OSPA值越小，表明目标数目和状态估计越准确。Among them, d _c (X, Y) = min {c, d _b (X, Y)},

c＞0 is the truncation parameter, which is used to penalize the estimation deviation of the number of targets, and p is the order, which is used to penalize the estimation deviation of multi-target states. In the simulation experiment of this paper, p＝3 and c＝1000 are set. The smaller the OSPA value is, the better. , indicating that the target number and state estimation are more accurate.

Simulation results and analysis: Five targets are set to move in a uniform straight line in the scene, and the original trajectory of the target is shown in Figure 2. This simulation considers the comparison of the complex likelihood and square likelihood methods in the method of the present invention. In order to better reflect the effectiveness of the tracking effect of the algorithm in this paper, the tracking performance is illustrated by performing OSPA error statistics and potential estimation statistics on 50 Monte Carlo experiments.

Please refer to Figure 3, which compares the OSPA corresponding to the amplitude likelihood and the complex likelihood under swerling 0. It can be seen from the figure that from the appearance of the target to the disappearance of the target, the OSPA of the complex likelihood is lower than the amplitude likelihood. When a new target appears, the OSPA will suddenly increase. This is because the adaptive new algorithm uses the measurement of the previous moment and has a certain delay. As shown in Figure 4, the average potential estimates under the two methods are relatively close. The difference between the effects of the complex likelihood and the amplitude likelihood under constant amplitude is not very obvious, but the average running time of the complex likelihood is 139 seconds, and the average running time of the amplitude likelihood is 1281 seconds. The running time of the amplitude likelihood is almost 9.22 times that of the complex likelihood.

Please refer to Figure 5, which compares the OSPA pairs of complex likelihood and amplitude likelihood under PHD filtering under swerling 1. It can be seen from the figure that the complex likelihood estimates the target position information better than the amplitude likelihood. Please refer to Figure 6. The main reason for the excessively high OSPA of the amplitude likelihood is the estimation of the target potential. The amplitude likelihood will estimate a lot of wrong position information during the initialization process, and most of these false detections are at the edge of the radar. The target phase information lost by the amplitude likelihood makes it impossible to reach the maximum target potential at the estimation time; the amplitude likelihood estimation effect under swerling 1 is not good. The average running time of the complex likelihood is 132 seconds, while the average running time of the amplitude likelihood is 347 seconds.

Please refer to Figures 7 and 8. Under the condition of swerling 3, the loss of PHD-TBD complex likelihood compared to amplitude likelihood is the smallest. Swerling 3 describes the characteristics of a target composed of many weaker scatterers and one particularly strong scatterer. For swerling 3, if the particle weight sum is infinite, the particle state is kept unchanged and the particle weight is replaced by a number as small as possible. Although this solves the problem that the target track cannot be continued due to the infinite particle sum, the update of the target state is abandoned, which will affect the estimation of the number of targets and the estimation of the state. The average running time of swerling 3 complex likelihood under PHD is 177 seconds, while the average running time of amplitude likelihood is 452 seconds.

Please refer to FIG9 and FIG10 . The amplitude likelihood with a signal-to-noise ratio of 5 dB under swerling 0 is almost completely impossible to estimate correctly, while the complex likelihood with a tracking performance of 5 dB is good. The tracking effect of the complex likelihood is not greatly affected when the signal-to-noise ratio is reduced.

Please refer to Figures 11 and 12. When the fluctuation type is swerling 1, reducing the signal-to-noise ratio will cause the target number estimate of the amplitude likelihood to increase at the initial moment, while compared with the complex likelihood, the target number estimate and the target position estimate are not greatly affected for the entire moment.

As shown in Figures 13 and 14, the target number estimation of the complex likelihood under swerling 3 will lead to inaccurate target number estimation as the signal-to-noise ratio decreases, resulting in a certain amount of missed detection. However, compared with the error caused by the amplitude likelihood, the impact of the complex likelihood on the reduction of the signal-to-noise ratio is smaller.

In summary, the method of the present invention can solve the problem of detection and tracking of weak multi-target fluctuations. In order to solve the problem of target generation under the condition of unknown prior distribution information of the new target state, an adaptive target generation algorithm based on measurement likelihood under scene division is proposed. Two calculation methods, amplitude likelihood and complex likelihood, are given, and PHD-TBD multi-target estimation with amplitude fluctuation types of swerling 0, 1, and 3 is successfully realized. Compared with the amplitude likelihood, the complex likelihood makes up for the defect that the amplitude likelihood only considers the measured amplitude information and ignores the phase information during the calculation process, thereby better utilizing the original information of the target. In the method proposed in the present invention, compared with the amplitude likelihood, the former is superior to the latter in the estimation performance of target position and number, and has higher calculation efficiency. Under low signal-to-noise ratio, the complex likelihood can still effectively detect and track an unknown number of weak targets.

What is disclosed above is only a preferred embodiment of the present invention, and it certainly cannot be used to limit the scope of rights of the present invention. Ordinary technicians in this field can understand that all or part of the processes of the above embodiment and equivalent changes made according to the claims of the present invention still fall within the scope of the invention.

Claims (6)

1. A method for tracking weak multiple targets before detection based on PHD filtering radar fluctuations, characterized in that it includes the following steps: S1: Initialize system parameters and read the original measurement data at the k-1th moment and the kth moment in the radar receiver; S2: Divide the scene and use the original measurement data at time k-1 to adapt the new target in each small scene; S3: For the complex measurement and square measurement data obtained at time k, the complex likelihood and amplitude likelihood of the three types of amplitude fluctuations are calculated respectively, and the SMC implementation of PHD filtering under amplitude fluctuation is given; S4: Extract target state and target number based on posterior information; S5: Let k = k + 1, and determine whether k > K. If the proposition is true, the algorithm ends, otherwise return to S2; Specifically, in step S3, the amplitude likelihood p(|z _k | ² ) when the target does not exist is calculated:

Calculate the amplitude likelihood for amplitude fluctuation types swerling 0, 1, and 3 respectively:

in,

Calculate the complex likelihood p(z _k ) when the target does not exist:

Calculate the complex likelihood for amplitude fluctuation types swerling 0, 1, and 3 respectively:

where I ₀ (·) is the modified Bessel function of the first kind; Specifically, for the k-1 moment, use a set of weighted

Particles

represents the posterior density of PHD, that is:

Predicted particles:

in,

and b _k (·|z _k ) is the proposed density, L _k-1 is the number of particles of the target surviving at time k-1, and J _k is the number of particles of the newly born target at time k; then the predicted intensity D _k-1|k is:

in:

in,

The status is

The survival probability of the target from time k-1 to time k;

is the derivative probability density of state x _k-1 ;

is the new birth probability density; The update of PHD can be expressed in SMC as:

in:

2. The method for tracking weak multiple targets before detection based on PHD filtering radar fluctuations as claimed in claim 1, characterized in that the system parameters include: Sampling interval T, current time k, total target motion time K, radar scans area [r _min ,r _max ]×[θ _min ,θ _max ] in polar coordinates, radar receives measurement data Z _k and Z _k-1 in the tracking scene, consider the range and azimuth surveillance radar covering the defined area in polar coordinates, for range, assume that the transmitted pulse is a linear frequency modulation signal with bandwidth B and duration T _ε , speed of light c, range resolution unit

For the angle, consider the linear phased array of Na _antennas at the radar receiving end, with a spacing of

Where λ is the wavelength of the carrier frequency and the angular resolution is

3. The method for tracking weak multiple targets before detection based on PHD filtering radar fluctuations as claimed in claim 2, characterized in that in the step of dividing the scenes and adaptively generating new targets using the original measurement data at time k-1 in each small scene: Divide the scene N _r ×N _θ into

scenes, consider generating N _filter particles in each scene, then the number of particles in the whole scene is n ₁ ×n ₂ ×N _filter ; Set the measurement cutoff threshold _Thcutoff , sort the measurements in each scene in descending order, select the measurements with strength higher than the _threshold Thcutoff, and consider the measurements lower than _Thcutoff as false positive measurements. The position information of each measurement is (n _r ,n _θ ); Convert each measured position information (n _r ,n _θ ) into the target position in the plane rectangular coordinate system, recorded as (z _x ,z _y ); Generate particles near each (z _x ,z _y ), calculate the likelihood of each particle in the scene, then resample to select particles with higher weights in the scene, and normalize the weights of the selected particles. 4. The method for tracking weak multiple targets before detection based on PHD filtering radar fluctuations as claimed in claim 3, characterized in that: For complex and square measurements, the likelihood of fluctuations of different amplitudes is different. During target tracking, the target amplitude will affect the strength of the target echo. Targets with high echo intensity will often suppress targets with low intensity, and the likelihood above the threshold will be assigned the highest likelihood. 5. The method for tracking weak multiple targets before detection based on PHD filtering radar fluctuations as claimed in claim 4, characterized in that: When the target amplitude fluctuation type is swerling 3, it is necessary to make a judgment after calculating the particle weight. If the particle weight sum is infinite, the predicted particle state is maintained, the particle is not updated, and the particle weight is assumed to be 0. 6. The method for tracking weak multiple targets before detection based on PHD filtering radar fluctuations as claimed in claim 5, characterized in that in the step of extracting target states and target numbers based on a posteriori information: Calculate the target number and resample the samples; Determine the updated particle weight

Whether the sum is greater than 0, if it is greater than 0, the resampled particles are clustered using the K-means method, and the particle groups with weights and sums lower than the threshold Th _k-means are deleted, and the particle groups with weights and sums higher than the threshold Th _k-means are regarded as the estimated target state. If it is less than 0, it is considered that there is no target in the scene at this time.

Images