CN113093174B - PHD filter radar fluctuation weak multi-target-based pre-detection tracking method - Google Patents

PHD filter radar fluctuation weak multi-target-based pre-detection tracking method Download PDF

Info

Publication number
CN113093174B
CN113093174B CN202110235449.6A CN202110235449A CN113093174B CN 113093174 B CN113093174 B CN 113093174B CN 202110235449 A CN202110235449 A CN 202110235449A CN 113093174 B CN113093174 B CN 113093174B
Authority
CN
China
Prior art keywords
target
likelihood
amplitude
fluctuation
particle
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202110235449.6A
Other languages
Chinese (zh)
Other versions
CN113093174A (en
Inventor
吴孙勇
李东升
薛秋条
孙希延
纪元法
蔡如华
符强
王守华
严素清
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Guilin University of Electronic Technology
Original Assignee
Guilin University of Electronic Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Guilin University of Electronic Technology filed Critical Guilin University of Electronic Technology
Priority to CN202110235449.6A priority Critical patent/CN113093174B/en
Publication of CN113093174A publication Critical patent/CN113093174A/en
Application granted granted Critical
Publication of CN113093174B publication Critical patent/CN113093174B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/66—Radar-tracking systems; Analogous systems
    • G01S13/72—Radar-tracking systems; Analogous systems for two-dimensional tracking, e.g. combination of angle and range tracking, track-while-scan radar
    • G01S13/723—Radar-tracking systems; Analogous systems for two-dimensional tracking, e.g. combination of angle and range tracking, track-while-scan radar by using numerical data
    • G01S13/726—Multiple target tracking
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
    • G01S7/02—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
    • G01S7/41—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00 using analysis of echo signal for target characterisation; Target signature; Target cross-section
    • G01S7/414—Discriminating targets with respect to background clutter
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
    • G01S7/02—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
    • G01S7/41—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00 using analysis of echo signal for target characterisation; Target signature; Target cross-section
    • G01S7/415—Identification of targets based on measurements of movement associated with the target
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
    • G01S7/02—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
    • G01S7/41—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00 using analysis of echo signal for target characterisation; Target signature; Target cross-section
    • G01S7/418—Theoretical aspects
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
    • Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
    • Y02A90/10—Information and communication technologies [ICT] supporting adaptation to climate change, e.g. for weather forecasting or climate simulation

Landscapes

  • Engineering & Computer Science (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Radar Systems Or Details Thereof (AREA)

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

PHD filter radar fluctuation weak multi-target-based pre-detection tracking method
Technical Field
The invention relates to the technical field of radar fluctuation weak multi-target detection and tracking, in particular to a PHD filter radar fluctuation weak multi-target pre-detection tracking method.
Background
With the rapid development of modern electronic information technology, radar target detection technology faces the threat of aircraft stealth technology, stealth technology development makes radar Reflection Cross Section (RCS) of stealth targets smaller, target reflection echo signals weak, and echo signal-to-noise ratio (SNR) lower. A conventional detection and tracking method of a moving target is post detection tracking (DBT), in which a threshold is set for each frame of measurement data to determine whether the target exists, and a track of the target is obtained through a tracking algorithm; however, when the signal-to-noise ratio of the echo signal is low, the echo of the weak target is usually lower than a threshold value, missed detection occurs, and the target track is difficult to extract by using single-frame measurement data; if the threshold is lowered, a large number of false alarms are generated, and the target track cannot be maintained.
To solve the above problem, one of the methods is to employ a pre-detection Tracking (TBD) algorithm. According to the algorithm, multi-frame data are processed in a combined mode according to continuity of target motion in space and time relevance of continuous frames of target echo data, and target detection and tracking are achieved through multi-frame energy accumulation. Traditional TBD methods include dynamic programming, hough transform, and the like. However, these methods are on the one hand computationally intensive and computationally complex, and on the other hand only suitable for linear gaussian models.
Particle filter TBD (PF-TBD) algorithm under Bayesian framework can solve the problems of nonlinearity and non-Gaussian, so that the method has been rapidly developed in the field of weak target detection and tracking. The PF-TBD method is limited in that the appearance (new generation) and disappearance (death) of the targets are not modeled, and thus the algorithm complexity is greatly increased when the number of the targets is unknown and varies, and there is a limitation in the filter performance.
Disclosure of Invention
The invention aims to provide a PHD filtering radar fluctuation weak multi-target-based pre-detection tracking method, which aims to solve the technical problem that weak multi-target fluctuation cannot be tracked effectively in the prior art.
In order to achieve the purpose, the invention discloses a PHD filter radar fluctuation weak multi-target-based pre-detection tracking method, which comprises the following steps:
s1: initializing system parameters, and reading original measurement data of the kth-1 moment and the kth moment in a radar receiver;
s2: dividing scenes, and self-adapting to a new target by using the original measurement data at the moment k-1 in each small scene;
s3: for complex measurement data and square measurement data obtained at the moment k, respectively calculating complex likelihood and amplitude likelihood of three amplitude fluctuation types, and giving out PHD filtering SMC realization under the amplitude fluctuation;
s4: extracting the target states and the target number according to posterior information;
s5: let k=k+1, judge K > K and hold, the proposition holds the algorithm and ends, if otherwise return to step two.
The system parameters include:
sampling interval T, current time K, total target movement time K, and radar scanning region [ r ] in polar coordinates min ,r max ]×[θ min ,θ max ]Radar receives metrology data Z within a tracking scene k and Zk-1 The radar is monitored taking into account the distance and the azimuth covering a defined area in polar coordinates, for which it is assumed that the transmitted pulse is of bandwidth B and duration T ε Linear frequency modulated signal, light velocity c, distance resolution unit
Figure BDA0002959829220000021
For angle, consider N at the radar receiver a Linear phased array of antennas spaced apart by
Figure BDA0002959829220000022
Wherein lambda is the wavelength of the carrier frequency and the angular resolution is +.>
Figure BDA0002959829220000023
Dividing scenes, and in the step of self-adapting to a new target by using the original measurement data at the k-1 moment in each small scene:
if N r ×N θ Larger, scene N r ×N θ Is divided into
Figure BDA0002959829220000024
Each scene, consider N generated within each scene filter The number of particles in the whole scene is n 1 ×n 2 ×N filter ；
Setting a measurement cutoff threshold Th Cutting off The measurements in each scene are arranged in descending order, and the intensity is selected to be higher than the threshold Th Cutting off Is lower than Th Cutting off Is considered as false detection measurement, and the position information of each measurement is (n r ,n θ )；
Each measured position information (n r ,n θ ) Is converted into a target position under a plane rectangular coordinate system and is marked as (z x ,z y )；
At each (z x ,z y ) And generating particles nearby, calculating likelihood for each particle on a scene where the particle is located, resampling to select particles with higher particle weights in the scene, and normalizing the selected particle weights.
For complex measurement and square measurement data obtained at the moment k, respectively calculating complex likelihood and amplitude likelihood of three amplitude fluctuation types, and providing PHD filtering under the amplitude fluctuation in the SMC implementation step:
calculating the amplitude likelihood when the target does not exist;
respectively calculating amplitude likelihood under the condition that the amplitude fluctuation type is powerling 0,1,3;
calculating complex likelihood when the target does not exist;
the complex likelihood is calculated for the amplitude fluctuation type of swollenling 0,1,3, respectively.
At the time of k-1:
with a set of weights
Figure BDA0002959829220000031
Particles of->
Figure BDA0002959829220000032
Predicting and updating the set of weighted particles to obtain +.>
Figure BDA0002959829220000033
For complex measurement and square measurement, the likelihood of fluctuation of different amplitudes is different, in the process of tracking the target, the target amplitude can influence the echo intensity of the target, the target with high echo intensity can often inhibit the target with low intensity, and the likelihood higher than the threshold value is considered to be the highest likelihood.
And when the target amplitude fluctuation type is powerling 3, judging after calculating the particle weight, if the particle weight is infinite, keeping the state of the predicted particle unchanged, and not updating the particle, wherein the particle weight is assumed to be 0.
In the step of extracting the target state and the target number according to the posterior information:
calculating the number of targets, and resampling samples;
judging the updated particle weight
Figure BDA0002959829220000034
If the weight is larger than 0, clustering the resampled particles by using a K-means method, and deleting the weight and the weight which are lower than a threshold Th k-means Particle population of (2) by summing the weights above a threshold Th k-means If the particle swarm is smaller than 0, then the scene is considered to have no target.
The beneficial effects of the invention are as follows: the self-adaptive target neogenesis algorithm based on measurement likelihood under scene division is provided for solving the problem of target neogenesis under the condition that the prior distribution information of the neogenesis target state is unknown. Two calculation modes of amplitude likelihood and complex likelihood are given, and PHD-TBD multi-objective estimation with amplitude fluctuation type of powerling 0,1,3 is successfully realized. Compared with the amplitude likelihood, the complex likelihood makes up the defect that the amplitude likelihood only considers the measured amplitude information in the calculation process and ignores the phase information, thereby better utilizing the original information of the target. Compared with the complex likelihood and the amplitude likelihood, the method provided by the invention has the advantages that the estimation performance of the target position and the number is better than that of the target position and the 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.
Drawings
In order to more clearly illustrate the embodiments of the invention or the technical solutions in the prior art, the drawings that are required in the embodiments or the description of the prior art will be briefly described, it being obvious that the drawings in the following description are only some embodiments of the invention, and that other drawings may be obtained according to these drawings without inventive effort for a person skilled in the art.
FIG. 1 is a flow chart of steps of a pre-detection tracking method based on PHD filtering radar fluctuation weak multi-target.
FIG. 2 is a schematic flow chart of a PHD filter radar fluctuation weak multi-target-based pre-detection tracking method of the invention.
Fig. 3 is a real trace of the movement of a target within a simulation scenario of the present invention.
Fig. 4 is OSPA of the complex likelihood and amplitude likelihood algorithm under the inventive powerling 0.
Fig. 5 is a graph comparing the target potentials of the complex likelihood and amplitude likelihood algorithm under powerling 0 of the present invention.
Fig. 6 is OSPA of the complex likelihood and amplitude likelihood algorithm under powerling 1 of the present invention.
Fig. 7 is a graph comparing the target potentials of the complex likelihood and amplitude likelihood algorithms under powerling 1 of the present invention.
Fig. 8 is OSPA of the complex likelihood and amplitude likelihood algorithm under powerling 3 of the present invention.
Fig. 9 is a graph comparing the target potentials of the complex likelihood and amplitude likelihood algorithm under powerling 3 of the present invention.
Fig. 10 is OSPA for different signal-to-noise ratio complex and square likelihoods at powerling 0 of the present invention.
FIG. 11 is a graph showing the potential of a target for comparison of different signal to noise ratios at different times under the powerling 0 of the present invention.
Fig. 12 is OSPA for different signal-to-noise ratio complex and square likelihoods under the powerling 1 of the present invention.
Fig. 13 is a graph showing the potential of the target for comparison of different signal to noise ratios at different times under the powerling 1 of the present invention.
Fig. 14 is OSPA for different signal-to-noise ratio complex and square likelihoods under the powerling 3 of the present invention.
Fig. 15 is a graph showing the potential of the target for comparison of different signal to noise ratios at different times under the powerling 3 of the present invention.
Detailed Description
Embodiments of the present invention are described in detail below, examples of which are illustrated in the accompanying drawings, wherein like or similar reference numerals refer to like or similar elements or elements having like or similar functions throughout. The embodiments described below by referring to the drawings are illustrative and intended to explain the present invention and should not be construed as limiting the invention.
In the description of the present invention, it should be understood that the terms "length," "width," "upper," "lower," "front," "rear," "left," "right," "vertical," "horizontal," "top," "bottom," "inner," "outer," and the like indicate orientations or positional relationships based on the orientation or positional relationships shown in the drawings, merely to facilitate describing the present invention and simplify the description, and do not indicate or imply that the devices or elements referred to must have a specific orientation, be configured and operated in a specific orientation, and therefore should not be construed as limiting the present invention. Furthermore, in the description of the present invention, the meaning of "a plurality" is two or more, unless explicitly defined otherwise.
Referring to fig. 1 and 2, the invention provides a pre-detection tracking method based on PHD filtering radar fluctuation weak multi-target, which comprises the following steps:
s1: initializing system parameters, and reading original measurement data of the kth-1 moment and the kth moment in a radar receiver;
s2: dividing scenes, and self-adapting to a new target by using the original measurement data at the moment k-1 in each small scene;
s3: for complex measurement data and square measurement data obtained at the moment k, respectively calculating complex likelihood and amplitude likelihood of three amplitude fluctuation types, and giving out PHD filtering SMC realization under the amplitude fluctuation;
s4: extracting the target states and the target number according to posterior information;
s5: let k=k+1, judge K > K and hold, the proposition holds the algorithm and ends, if otherwise return to step two.
Specifically, initializeSystem parameters including sampling interval T, current time K, total target motion time K, radar scanning region in polar coordinates min ,r max ]×[θ min ,θ max ]Radar receives metrology data Z within a tracking scene k and Zk-1 The radar is monitored taking into account the distance and the azimuth covering a defined area in polar coordinates, for which it is assumed that the transmitted pulse is of bandwidth B and duration T ε Linear frequency modulated signal, light velocity c, distance resolution unit
Figure BDA0002959829220000051
For angle, consider N at the radar receiver a Linear phased array of antennas at intervals +.>
Figure BDA0002959829220000052
Wherein lambda is the wavelength of the carrier frequency, and the angular resolution is
Figure BDA0002959829220000053
Wherein k=2, b is initialized k|k-1 (x k |x k-1) and γk (x k ) PHD, kappa representing target neogenesis and derivatization at time k k (z)＝λ k C (z) is the false alarm intensity, lambda k C (z) is the false alarm distribution for the average false alarm number; p (P) D (x) Is the probability of detection of the target state g k (z|x) represents the target generation measurement likelihood. In multi-target tracking, updating the PHD multi-target state function D k|k Is the integral of the target number at time k
Figure BDA0002959829220000061
Specifically, the measurement value received by the radar sensor is composed of the distance and the azimuth after the distance matching filtering and the adaptive beam forming. Given the ith target state x at time k k,i Measurement value z k Given by the following nonlinear equation:
Figure BDA0002959829220000062
wherein ,h(xk.i ) Expressed in terms of target position (x k,i ,y k,i ) A fuzzy function for the ith object as center, h (x for simplicity k.i ) Is marked as h k,i 。n k Is the measurement noise with the mean value of 0 and the covariance of complex Gauss gamma;
Figure BDA0002959829220000063
and ρk,i The phase and amplitude of the i-th target complex amplitude are represented, respectively. Let all phases +.>
Figure BDA0002959829220000064
Sum amplitude ρ k,1 :N k Independent of each other and independent of n k and xk.1:Nk . Phase->
Figure BDA0002959829220000065
Is unknown and uniformly distributed over [0,2 pi ] at each instant; for amplitude ρ k,i The method comprises the following steps:
Figure BDA0002959829220000066
wherein />
Figure BDA00029598292200000611
Is an unknown static parameter. There are two expression forms of radar measurement at time k: complex measuring z being coherent accumulation k Another square measure of incoherent accumulation is |z k | 2 。
The blur function for distance matched filtering is:
Figure BDA0002959829220000067
wherein ,|τl |＝2(r k -r l )/c，
Figure BDA0002959829220000068
l∈[0,N r -1]，/>
Figure BDA0002959829220000069
The azimuth ambiguity function of adaptive beamforming is:
Figure BDA00029598292200000610
wherein ,Φm ＝π[sin(θ k )-sin(θ m )]，
Figure BDA0002959829220000071
and />
Figure BDA0002959829220000072
m∈[0,N θ -1]，/>
Figure BDA0002959829220000073
Overall blur function in distance-bearing unit (l, m):
Figure BDA0002959829220000074
h(x k ) Is of size N c ＝N r ×N θ The method comprises the following steps:
Figure BDA0002959829220000075
in the conventional radar target detection and tracking problem, firstly, threshold processing is performed on echo data of each frame to form point trace information, then, correlation, filtering and other processing are performed on the point trace information exceeding the threshold value, and finally, a track of a target is obtained, so that the tracking of the target is realized. For the method of detecting before tracking, the method is suitable for the condition of higher signal-to-noise ratio or larger target echo amplitude, the intensity of the target is far higher than that of the clutter, and the target can be separated from the clutter by setting a larger threshold value. For the conditions of low signal-to-noise ratio or weaker target echo signals, the target echo is annihilated in noise clutter, a single-frame threshold crossing detection method is adopted, the target is missed due to the fact that the threshold is too high, the false alarm rate is increased due to the fact that the threshold is too low, and the target track cannot be maintained. The TBD method can improve the detection and tracking of the weak target of the radar from the angle of signal processing, and the basic idea is to process the original data which is not subjected to threshold processing, realize the capability accumulation of the target through multi-frame tracking and fully mine the information of the target in the echo.
Further, the PHD filtering algorithm based on RFS simplifies the multi-target state space into a single-target state space by transmitting the first moment of the multi-target posterior probability density distribution, so that the calculation complexity is reduced to a great extent, and the probability of actual operation is realized. The prediction and update equations for PHD filtering are as follows:
Figure BDA0002959829220000076
for updating equation D k|k ，P D (x) Is a detection probability associated with the target state; for the TBD algorithm, there is no detection process prior to the update step, assuming that the measurement contains all the target information; thus P D (x) ≡1), the updated equation becomes:
Figure BDA0002959829220000081
for a given target state, the intensity information in each cell is independent of the other, so that the multi-target posterior probability density p (z k |X k ) Can be expressed as the product of the edge probability density functions:
Figure BDA0002959829220000082
specifically, in step 2, if N r ×N θ Larger, scene N r ×N θ Is divided into
Figure BDA0002959829220000083
A scene; consider generating N within each scene filter The number of particles in the whole scene is n 1 ×n 2 ×N filter ；
Setting a measurement cutoff threshold Th Cutting off The measurements in each scene are arranged in descending order, and the intensity is selected to be higher than the threshold Th Cutting off Is lower than Th Cutting off Is considered as false detection measurement, and the position information of each measurement is (n r ,n θ )；
Utilizing the number of the distances and the angles by utilizing a formula
Figure BDA0002959829220000084
And
Figure BDA0002959829220000085
converted into distance and azimuth, and then according to formula z x ＝rcosθ,z y =rsinθ converts the position in polar coordinates into position information in planar rectangular coordinates;
at each (z x ,z y ) Location information (x) of nearby generation target state i ,y i )，v＝U(v min ,v max ) Generating speed information of target state
Figure BDA0002959829220000086
Wherein U represents a uniform distribution, i.e., a particle is generated near each measurement, denoted as
Figure BDA0002959829220000087
For each particle x i Calculating likelihood on the scene, and then resampling to select particles with higher particle weights in the scene; and normalizing the selected particle weight.
Specifically, in step three, the amplitude likelihood p (|z) when the target is not present is calculated k | 2 )：
Figure BDA0002959829220000088
The amplitude likelihood under the condition that the amplitude fluctuation type is powerling 0,1,3 is calculated respectively:
Figure BDA0002959829220000091
Figure BDA0002959829220000092
/>
Figure BDA0002959829220000093
wherein ,
Figure BDA0002959829220000094
calculating complex likelihood p (z) when target is not present k )：
Figure BDA0002959829220000095
The complex likelihood under the amplitude fluctuation type of swollenling 0,1,3 is calculated respectively:
Figure BDA0002959829220000096
Figure BDA0002959829220000097
Figure BDA0002959829220000098
wherein ,I0 (. Cndot.) is a first class repairPositive bessel function.
Specifically, for time k-1, a set of weights is used
Figure BDA0002959829220000099
Particles of->
Figure BDA00029598292200000910
Represents the posterior density of PHD, namely:
Figure BDA00029598292200000911
predicted particles:
Figure BDA0002959829220000101
wherein ,
Figure BDA0002959829220000102
and bk (·|z k ) Is the recommended density, L k-1 Is the particle number of the surviving target at time k-1, J k Is the number of particles of the new target at time k. The predicted intensity D k-1|k The method comprises the following steps:
Figure BDA0002959829220000103
wherein ：
Figure BDA0002959829220000104
Figure BDA0002959829220000105
wherein ,
Figure BDA0002959829220000106
is a state of +.>
Figure BDA0002959829220000107
The target survival probability from time k-1 to time k;
Figure BDA0002959829220000108
is in the state x k-1 Is a derivative probability density of (2); />
Figure BDA0002959829220000109
Is the new probability density.
The update SMC of the PHD can be expressed as:
Figure BDA00029598292200001010
Figure BDA00029598292200001011
wherein ：
Figure BDA00029598292200001012
specifically, for complex measurement and square measurement, the likelihood of fluctuation of different amplitudes is different, in the process of tracking the target, the target amplitude can influence the target echo intensity, the target with high echo intensity can often inhibit the target with low intensity, and the likelihood higher than the threshold value is considered to be the highest likelihood.
Note that, as shown by the measurement equation, the echo signal of the target is a sinc function, and the likelihood ratio of the existence position of the target at this time is very high; when the amplitude of the target fluctuates, the intensities of the echo signals are different, and in order to solve the multi-target tracking problem under the condition, the matching tracking of targets with different intensities is solved, and the likelihood ratio is larger than the threshold Th L Is assigned the highest likelihood, i.e.:
L sw (L sw ≥Th L )＝max(L sw )。
specifically, when the target amplitude fluctuation type is powerling 3, the situation that the particle weight sum is infinite may occur, and the situation may cause that the target cannot be tracked at the subsequent moment; for this phenomenon, the judgment is performed after the calculation of the particle weight, and if the particle weight is infinite, the state of the predicted particle is kept unchanged, the particle is not updated, and the particle weight is assumed to be 0.
Specifically, step four includes calculating the target number
Figure BDA0002959829220000111
Resampling sample +.>
Figure BDA0002959829220000112
Obtain->
Figure BDA0002959829220000113
Judging the updated particle weight +.>
Figure BDA0002959829220000114
If the sum is larger than 0, clustering the resampled particles by using a K-means method, deleting the weight and being lower than a threshold Th k-means Is a particle group of (2); sum the weights above threshold Th k-means If the particle swarm is smaller than 0, then the scene is considered to have no target.
First embodiment:
the present embodiment uses MATLAB software version 2014 (a) for simulation experiments.
Referring to fig. 2, five targets are set to move, and consider a two-dimensional motion scene, the state of each target is defined as
Figure BDA0002959829220000115
wherein (xk ,y k) and
Figure BDA0002959829220000116
the position and velocity of the target in a cartesian coordinate system, respectively.
Position (x) k ,y k ) In polar coordinates p= [ r ] min ,r max ]×[θ min ,θ max ]Within the scene r min ,r max and θmin ,θ max Respectively minimum and maximum target ranges and orientations;
speed of speed
Figure BDA0002959829220000117
In the area:
Figure BDA0002959829220000118
wherein vmin and vmax The minimum and maximum target speeds, respectively.
The target state evolves according to uniform linear motion:
x k ＝Fx k-1 +v k
the sensor scan time interval t=1s receives 100 frames of images, wherein
Figure BDA0002959829220000119
Figure BDA0002959829220000121
Process noise v k Obeying a gaussian distribution, the covariance is:
Figure BDA0002959829220000122
standard deviation sigma of noise v =5m, probability of target survival p s,k (x)＝0.98。
For simulation of targets, v min ＝200m/s，v max =800 m/s, signal to noise ratio is defined by
Figure BDA0002959829220000123
Given.
For simulation of radar measurements, r min ＝100km，r max ＝120km，θ min ＝-75°，θ max ＝75°，N r ＝300，N θ ＝100，σ 2 =0.5, noise covariance is
Figure BDA0002959829220000124
B＝150KHz，T e ＝6.67×10 -5 s，N a ＝50，λ＝3cm，c＝3×10 8 m/s，Δ r ＝500m，Δ θ ＝1.45°。
The method comprises the following steps: the present embodiment aims at track conditions of five targets
Figure BDA0002959829220000125
The performance of the evaluation algorithm adopts an optimal sub-mode allocation distance (OSPA), the OSPA metric can evaluate the target number estimation error and the target position estimation error of the multi-target filter, and two finite sets X= { X are given 1 ,x 2 ,…x m} and Y＝{y1 ,y 2 ,…y n OSPA is defined as follows:
Figure BDA0002959829220000126
wherein ,dc (X,Y)＝min{c,d b (X,Y)}，
Figure BDA0002959829220000131
c > 0 is a cut-off parameter for punishing the estimated deviation of the number of targets, and p is an order for punishing the estimated deviation of the states of multiple targets. In the simulation experiments herein, p=3 and c=1000 are set. The smaller the OSPA value, the more accurate the target number and state estimation.
Simulation results and analysis: 5 targets are set to do uniform linear motion in the scene, and the original track of the targets is shown in fig. 2. In the simulation, the complex likelihood and the square likelihood in the method are compared, and in order to better embody the effectiveness of the tracking effect of the algorithm, the tracking performance is illustrated by carrying out OSPA error statistics and potential estimation statistics on 50 Monte Carlo experiments.
Referring to fig. 3, comparing the amplitude likelihood with the complex likelihood corresponding OSPA under the powerling 0, it can be seen that the complex likelihood is lower than the amplitude likelihood from the appearance of the target to the disappearance of the target. When a new target appears, OSPA increases suddenly, because the adaptive new generation algorithm is caused by using the measurement at the previous time, and has a certain delay. As can be seen from fig. 4, the average potential estimates under the two methods are relatively close. The difference in effect between the complex likelihood and the amplitude likelihood at a constant amplitude is not obvious, but the average running time of the complex likelihood is 139 seconds, the average running time of the amplitude likelihood is 1281 seconds, and the running time of the amplitude likelihood is almost 9.22 times that of the complex likelihood.
Referring to fig. 5, the OSPA pair of complex likelihood and amplitude likelihood under PHD filtering under the powerling 1 is compared, and the effect of estimating the target position information by the complex likelihood is better than the effect of the amplitude likelihood. Referring to fig. 6, the reason why the amplitude likelihood OSPA is too high is mainly that for estimating the target potential, the amplitude likelihood estimates more false position information in the initial process, and most of the false detections are at the radar edge, and the target phase information of the amplitude likelihood loss causes that the target maximum potential cannot be reached at the estimation moment; the amplitude likelihood estimation at powerling 1 does not work well. The complex likelihood average run time was 132 seconds and the amplitude likelihood average run time was 347 seconds.
Referring to fig. 7 and 8, under the powerling 3 condition, the PHD-TBD complex likelihood contrast amplitude likelihood loss is minimal. Powerling 3 describes the target properties consisting of a number of weaker scatterers and a particularly strong scatterer. For the powerling 3, if the particle weight and infinity are the particle state is kept unchanged, the particle weight is replaced by a number as small as possible, and in this way, although the problem that the target track cannot be continued due to the particles and infinity is solved, the update of the target state is abandoned, and the estimation of the number of targets and the estimation of the state are influenced. The powerling 3 complex likelihood average run time under PHD is 177 seconds and the amplitude likelihood average run time is 452 seconds.
Referring to fig. 9 and 10, the amplitude likelihood of 5dB signal-to-noise ratio at powerling 0 is almost completely uncorrectable, while the complex likelihood tracking performance of 5dB is good; the tracking effect of complex likelihood is not greatly affected in the case of a reduced signal-to-noise ratio.
Referring to fig. 11 and 12, in the case where the fluctuation type is powerling 1, decreasing the signal-to-noise ratio causes the number of target estimates of the amplitude likelihood to increase at the initial time, whereas the number of target estimates and the estimation of the target position are less affected for the entire time than the complex likelihood.
Referring to fig. 13 and 14, it can be seen that, in the complex likelihood target number estimation under the answer 3, the target number estimation is inaccurate due to the decrease of the signal-to-noise ratio, and a certain omission factor is generated, but compared with the error caused by the amplitude likelihood, the influence of the complex likelihood on the decrease of the signal-to-noise ratio is smaller.
In summary, the method can solve the detection and tracking problems of weak multi-target fluctuation, and provides a self-adaptive target neogenesis algorithm based on measurement likelihood under scene division for solving the target neogenesis problem under the condition that the prior distribution information of the neogenesis target state is unknown. Two calculation modes of amplitude likelihood and complex likelihood are given, and PHD-TBD multi-objective estimation with amplitude fluctuation type of powerling 0,1,3 is successfully realized. Compared with the amplitude likelihood, the complex likelihood makes up the defect that the amplitude likelihood only considers the measured amplitude information in the calculation process and ignores the phase information, thereby better utilizing the original information of the target. Compared with the complex likelihood and the amplitude likelihood, the method provided by the invention has the advantages that the estimation performance of the target position and the number is better than that of the target position and the 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.
The above disclosure is only a preferred embodiment of the present invention, and it should be understood that the scope of the invention is not limited thereto, and those skilled in the art will appreciate that all or part of the procedures described above can be performed according to the equivalent changes of the claims, and still fall within the scope of the present invention.

Claims (6)

1. The PHD filter radar fluctuation weak multi-target-based pre-detection tracking method is characterized by comprising the following steps of:
s1: initializing system parameters, and reading original measurement data of the kth-1 moment and the kth moment in a radar receiver;
s2: dividing scenes, and self-adapting to a new target by using the original measurement data at the moment k-1 in each small scene;
s3: for complex measurement data and square measurement data obtained at the moment k, respectively calculating complex likelihood and amplitude likelihood of three amplitude fluctuation types, and giving out PHD filtering SMC realization under the amplitude fluctuation;
s4: extracting the target states and the target number according to posterior information;
s5: let k=k+1, judge K > K and hold, the proposition holds the algorithm and ends, if otherwise return to S2;
specifically, in step S3, the amplitude likelihood p (|z) when the target is not present is calculated k | 2 )：
Figure FDA0003808927610000011
The amplitude likelihood under the condition that the amplitude fluctuation type is powerling 0,1,3 is calculated respectively:
Figure FDA0003808927610000012
Figure FDA0003808927610000013
Figure FDA0003808927610000014
wherein ,
Figure FDA0003808927610000015
calculating complex likelihood p (z) when target is not present k )：
Figure FDA0003808927610000016
The complex likelihood under the amplitude fluctuation type of swollenling 0,1,3 is calculated respectively:
Figure FDA0003808927610000017
Figure FDA0003808927610000021
/>
Figure FDA0003808927610000022
wherein ,I0 (. Cndot.) is a first type of modified Bessel function;
specifically, for time k-1, a set of weights is used
Figure FDA0003808927610000023
Particles of->
Figure FDA0003808927610000024
Represents the posterior density of PHD, namely:
Figure FDA0003808927610000025
predicted particles:
Figure FDA0003808927610000026
wherein ,
Figure FDA0003808927610000027
and bk (·|z k ) Is the recommended density, L k-1 Is the particle number of the surviving target at time k-1, J k The particle number of the new target at the time k; the predicted intensity D k-1|k The method comprises the following steps:
Figure FDA0003808927610000028
wherein ：
Figure FDA0003808927610000029
Figure FDA00038089276100000210
wherein ,
Figure FDA00038089276100000211
is a state of +.>
Figure FDA00038089276100000212
The target survival probability from time k-1 to time k; />
Figure FDA00038089276100000213
Is in the state x k-1 Is a derivative probability density of (2); />
Figure FDA00038089276100000214
Is the density of the new probability;
the update SMC of the PHD can be expressed as:
Figure FDA00038089276100000215
Figure FDA00038089276100000216
wherein ：
Figure FDA0003808927610000031
2. the pre-detection tracking method based on PHD filtering radar fluctuation weak multi-target according to claim 1, wherein the system parameters include:
sampling interval T, current time K, total target movement time K, and radar scanning region [ r ] in polar coordinates min ,r max ]×[θ min ,θ max ]Radar receives metrology data Z within a tracking scene k and Zk-1 The radar is monitored taking into account the distance and the azimuth covering a defined area in polar coordinates, for which it is assumed that the transmitted pulse is of bandwidth B and duration T ε Linear frequency modulated signal, light velocity c, distance resolution unit
Figure FDA0003808927610000032
For angle, consider N at the radar receiver a Linear phased array of antennas at intervals +.>
Figure FDA0003808927610000033
Wherein lambda is the wavelength of the carrier frequency and the angular resolution is +.>
Figure FDA0003808927610000034
3. The pre-detection tracking method based on PHD filtering radar fluctuation weak multi-target according to claim 2, wherein the scenes are divided, and the original measurement data at k-1 moment is utilized in the step of self-adapting to the new target in each small scene:
scene N r ×N θ Is divided into
Figure FDA0003808927610000035
Each scene, consider N generated within each scene filter The number of particles in the whole scene is n 1 ×n 2 ×N filter ；
Setting a measurement cutoff threshold Th Cutting off The measurements in each scene are arranged in descending order, and the intensity is selected to be higher than the threshold Th Cutting off Is lower than Th Cutting off Is considered as false detection measurement, and the position information of each measurement is (n r ,n θ )；
Each measured position information (n r ,n θ ) Is converted into a target position under a plane rectangular coordinate system and is marked as (z x ,z y )；
At each (z x ,z y ) And generating particles nearby, calculating likelihood for each particle on a scene where the particle is located, resampling to select particles with higher particle weights in the scene, and normalizing the selected particle weights.
4. The pre-detection tracking method based on PHD filter radar fluctuation weak multi-target according to claim 3,
for complex measurement and square measurement, the likelihood of fluctuation of different amplitudes is different, in the process of tracking the target, the target amplitude can influence the echo intensity of the target, the target with high echo intensity can often inhibit the target with low intensity, and the likelihood higher than the threshold value is assigned to be the highest likelihood.
5. The pre-detection tracking method based on PHD filter radar fluctuation weak multi-target according to claim 4, wherein,
when the target amplitude fluctuation type is powerling 3, the judgment is needed after the calculation of the particle weight, if the particle weight is infinite, the state of the predicted particle is kept, the particle is not updated, and the particle weight is assumed to be 0.
6. The pre-detection tracking method based on PHD filtering radar fluctuation weak multi-target as set forth in claim 5, wherein in the step of extracting the target state and the target number according to posterior information:
calculating the number of targets, and resampling samples;
judging the updated particle weight
Figure FDA0003808927610000041
If the weight is larger than 0, clustering the resampled particles by using a K-means method, and deleting the weight and the weight which are lower than a threshold Th k-means Particle population of (2) by summing the weights above a threshold Th k-means If the particle swarm is smaller than 0, then the scene is considered to have no target. />
CN202110235449.6A 2021-03-03 2021-03-03 PHD filter radar fluctuation weak multi-target-based pre-detection tracking method Active CN113093174B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110235449.6A CN113093174B (en) 2021-03-03 2021-03-03 PHD filter radar fluctuation weak multi-target-based pre-detection tracking method

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110235449.6A CN113093174B (en) 2021-03-03 2021-03-03 PHD filter radar fluctuation weak multi-target-based pre-detection tracking method

Publications (2)

Publication Number Publication Date
CN113093174A CN113093174A (en) 2021-07-09
CN113093174B true CN113093174B (en) 2023-05-30

Family

ID=76666475

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110235449.6A Active CN113093174B (en) 2021-03-03 2021-03-03 PHD filter radar fluctuation weak multi-target-based pre-detection tracking method

Country Status (1)

Country Link
CN (1) CN113093174B (en)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113687340B (en) * 2021-08-24 2024-02-23 重庆交通大学 Long-distance moving target detection method based on millimeter wave radar

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2014169865A (en) * 2013-03-01 2014-09-18 Hitachi Ltd Target tracking device, target tracking program and target tracking method
CN106772352A (en) * 2016-12-01 2017-05-31 中国人民解放军海军航空工程学院 A kind of PD radars extension Weak target detecting method based on Hough and particle filter
CN110109094A (en) * 2019-03-28 2019-08-09 西南电子技术研究所（中国电子科技集团公司第十研究所） The detection of multi-receiver station single frequency network external illuminators-based radar maneuvering target and tracking

Family Cites Families (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
NL1020287C2 (en) * 2002-04-02 2003-10-03 Thales Nederland Bv Method for multi-target detection, in particular for use in search radars with multi-beam formation in elevation.
US7630428B2 (en) * 2005-07-28 2009-12-08 Itt Manufacturing Enterprises, Inc. Fast digital carrier frequency error estimation algorithm using synchronization sequence
US8405540B2 (en) * 2010-04-02 2013-03-26 Mitsubishi Electric Research Laboratories, Inc. Method for detecting small targets in radar images using needle based hypotheses verification
CN104062651B (en) * 2014-06-30 2016-06-29 电子科技大学 A kind of based on tracking before the detection of G0 clutter background and constant target amplitude
CN104076355B (en) * 2014-07-04 2016-08-24 西安电子科技大学 Tracking before Dim targets detection in strong clutter environment based on dynamic programming
CN104714225B (en) * 2015-03-25 2017-05-10 电子科技大学 Dynamic programming tracking-before-detection method based on generalized likelihood ratios
CN105975772B (en) * 2016-05-04 2019-02-05 浙江大学 Tracking before multi-target detection based on probability hypothesis density filtering
CN107730537B (en) * 2017-09-29 2020-07-07 桂林电子科技大学 Weak target detection and tracking method based on box particle probability hypothesis density filtering
CN107703496B (en) * 2017-10-12 2021-04-30 桂林电子科技大学 Interactive multimode Bernoulli filtering maneuvering weak target tracking-before-detection method

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2014169865A (en) * 2013-03-01 2014-09-18 Hitachi Ltd Target tracking device, target tracking program and target tracking method
CN106772352A (en) * 2016-12-01 2017-05-31 中国人民解放军海军航空工程学院 A kind of PD radars extension Weak target detecting method based on Hough and particle filter
CN110109094A (en) * 2019-03-28 2019-08-09 西南电子技术研究所（中国电子科技集团公司第十研究所） The detection of multi-receiver station single frequency network external illuminators-based radar maneuvering target and tracking

Also Published As

Publication number Publication date
CN113093174A (en) 2021-07-09

Similar Documents

Publication Publication Date Title
CN107861107B (en) Double-threshold CFAR (computational fluid dynamics) and trace point agglomeration method suitable for continuous wave radar
CN110609262B (en) Three-dimensional constant false alarm detection method for scene surveillance radar
Rohling Some radar topics: waveform design, range CFAR and target recognition
KR20040091699A (en) An adaptive system and method for radar detection
WO2008093092A2 (en) Adaptive radar
CN107561509B (en) Airborne millimeter wave radar power line detection method
US20220114363A1 (en) Method and System for Indoor Multipath Ghosts Recognition
Bocquel et al. Multitarget particle filter addressing ambiguous radar data in TBD
Ristic et al. Bernoulli track‐before‐detect filter for maritime radar
Watts et al. Challenges in radar sea clutter modelling
Zhou et al. Multiple-kernelized-correlation-filter-based track-before-detect algorithm for tracking weak and extended target in marine radar systems
CN113093174B (en) PHD filter radar fluctuation weak multi-target-based pre-detection tracking method
CN109581366B (en) Discrete sidelobe clutter identification method based on target steering vector mismatch
CN108196238A (en) Clutter map detection method based on adaptive matched filter under Gaussian background
CN113866755A (en) Radar weak fluctuating target pre-detection tracking algorithm based on multi-Bernoulli filtering
CN113126086B (en) Life detection radar weak target detection method based on state prediction accumulation
DK180729B1 (en) System for processing radar data representing intensity values of received power of reflected radar wave signals
CN112612007B (en) Super-sparse array airborne radar moving target distance de-blurring method based on near field effect
CN114355324A (en) Flight path generation method
CN112684425B (en) Target secondary screening method after constant false alarm detection
Ahmed et al. Robust ground moving target detection for airborne radar using a novel feature-based machine learning approach
Ginolhac et al. Spatio-temporal adaptive detector in non-homogeneous and low-rank clutter
Wang et al. Multipath ghosts mitigation for radar-based positioning systems
CN114859337A (en) Data processing method and device, electronic equipment and computer storage medium
CN116879863B (en) Multi-target measuring method and system for continuous wave 4D millimeter wave radar

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant