CN109934840A - Circumference SAR motion target tracking method based on GMPHD filter - Google Patents
Circumference SAR motion target tracking method based on GMPHD filter Download PDFInfo
- Publication number
- CN109934840A CN109934840A CN201910222890.3A CN201910222890A CN109934840A CN 109934840 A CN109934840 A CN 109934840A CN 201910222890 A CN201910222890 A CN 201910222890A CN 109934840 A CN109934840 A CN 109934840A
- Authority
- CN
- China
- Prior art keywords
- target
- indicate
- gaussian
- moving target
- newborn
- 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.)
- Pending
Links
Abstract
Circumference SAR motion target tracking method based on GMPHD filter, the present invention relates to be related to SAR motion target tracking method.The purpose of the present invention is to solve existing methods to carry out the low problem of motion target tracking accuracy rate.Process are as follows: one: Preliminary detection is carried out to moving target based on DPCA-CFAR method, obtains the observation of moving target；Two: establishing the state vector of moving target and measure vector；Three: establishing moving target state equation and measurement equation based on CSAR system；Four: the initialization of Gaussian component mark value；Five: state equation and measurement equation based on moving target generate adaptive newborn target strength, and the Gauss mark value new for newborn Target Assignment；Six: carrying out GMPHD recursion based on five；Seven: carrying out moving target state estimation based on six and track extracts, complete SAR motion target tracking.The present invention is used for microwave remote sensing technique and radar data process field.
Description
Technical field
The present invention relates to microwave remote sensing technique and radar data process field more particularly to SAR motion target tracking methods.
Background technique
Circumferential synthetic aperture radar (CircularSynthetic Aperture Radar, CSAR) is as a kind of high over the ground
Resolution imaging operating mode has both the advantage for observing and obtaining 360 ° of omnidirectional's information of target for a long time, can be realized to sense
The lasting observation of interest region or target.It is imaged by continuous multiple frames SAR and obtains multiple image sequence, effective expansion time dimension letter
Breath combines SAR imaging technique with target following technology so that SAR system has stronger multidate information acquisition capability, from
And accurately obtain the moving parameter informations such as position, speed and the movement tendency of moving target.Can not only effectively it overcome infrared/visible
The weakness that optical sensor is influenced vulnerable to weather condition and battlefield surroundings, by soupy weather gas shielded itself not by enemy attack,
The disadvantages of conventional SAR system moving-target detecting and tracking can be overcome difficult.
Moving target is usually multiple target in SAR image scene, and residue is a large amount of miscellaneous after clutter recognition and CFAR detection
Wave, especially in urban area, clutter background is stronger.For multiple target tracking problem under high clutter background, traditional target association is calculated
Method need to spend it is a large amount of calculate to solve the problems, such as data correlation, and high false alarm rate cause false track quantity increase and track with
Track mistake, therefore actual demand is all unable to satisfy in computational efficiency and tracking performance based on associated multi-object tracking method.
It is influenced by imaging method, moving-target occurs orientation offset and defocused in SAR image.Since radar platform is circumference rail
The radial velocity of mark, target constantly changes, and offset can be more complicated on the figure of moving-target.Therefore it needs to study to be suitable for CSAR system
The multiple target tracking model and tracking of system.
Gaussian-mixture probability hypothesis density (Gaussian mixture probability hypothesis density,
GMPHD) filter passes through the first moment information of each moment multiple target state stochastic finite collection of iterative estimate, as Bayesian filter
Suboptimal solution, not only avoid complicated data correlation problem, while under Bayesian frame solving set valued integrals and be difficult to resolve and asking
Topic has ideal Bayes's meaning and propinquity effect.However the intensity function of newborn target is necessary in Traditional GM PHD filter
Be it is known, lack substantial newborn target detection function.The detection of newborn target or track are initially multiple target tracking
System necessity and important component part.And this method can only export number of targets information and target status information, cannot provide
System independence track information, there is certain limitation in engineer application；To sum up cause motion target tracking accuracy rate low.
Summary of the invention
The purpose of the present invention is to solve existing methods to carry out the low problem of motion target tracking accuracy rate, and proposes base
In the circumference SAR motion target tracking method of GMPHD filter.
Circumference SAR motion target tracking method detailed process based on GMPHD filter are as follows:
Step 1: Preliminary detection is carried out to moving target based on DPCA-CFAR method, obtains the observation of moving target；
Step 2: the moving target observation obtained according to step 1, establish moving target state vector and measure to
Amount；
Step 3: state vector and measurement vector based on moving target establish the moving target shape based on CSAR system
State equation and measurement equation；
Step 4: Gaussian component mark value initialization；
Step 5: state equation and measurement equation based on moving target generate adaptive newborn target strength, and are new
The new Gauss mark value of raw Target Assignment；
Step 6: the Gauss mark value of adaptive newborn target strength and newborn target based on step 5 carries out GMPHD
Recursion；
GMPHD recursive process includes PHD prediction, PHD update, the trimming and fusion of Gaussian term；
Step 7: moving target state estimation is carried out based on step 6 and track extracts, completes SAR motion target tracking.
The invention has the benefit that
The present invention proposes a kind of circumference SAR motion target tracking method based on GMPHD filter, is based on CSAR image sequence
Column, establish multiple target tracking model according to CSAR imaging geometry structure, and improve GMPHD filter, utilize adaptive newborn mesh
Mark intensity generating algorithm and estimate newborn intensity, the extraction to targetpath is realized by Gauss labelling method, it is final realize to SAR at
The extraction of the estimation of moving target kinematic parameter and track in image field scape.
1, the present invention constructs the multiple target tracking model based on circumference SAR system for the first time, and by improved GMPHD filter
For in SAR motion target tracking, realizing preferable tracking of maneuvering target effect；
2, the present invention is directed to the deficiency of imaging characteristics and Traditional GM PHD filter of the target in SAR image, proposes to utilize
Newborn target strength generating algorithm estimates newborn target strength, has newborn target detection capabilities；
3, the extraction to targetpath is realized using Gauss labelling method, significant effect overcomes Traditional GM PHD filter and exists
Limitation when SAR motion target tracking is carried out, there is stronger engineer application promotional value；
4, simulation result shows: proposed by the present invention based on the CSAR motion target tracking method for improving GMPHD filter
It can accurately estimate dbjective state and extract the motion profile of target, match, mention with target time of day and track
High motion target tracking accuracy rate, tracking accuracy is higher, and eliminates clutter residue, and being capable of each mesh of fast initialization
Mark, accurately estimates the initial state of each target.
Such as Fig. 5, different targetpaths is indicated with different lines, it can be seen that each dbjective state base of each frame estimation
On this on corresponding real goal track, the motion profile of 5 targets can be accurately extracted, with the true track kissing of target
It closes, and eliminates clutter residue.
If Fig. 6 a is in allowable range of error, by emulation it can be seen that the present invention is improved based on GMPHD filter
The moving-target number that circumference SAR motion target tracking method determines is 3.05, and Traditional GM PHD method for tracking target determines dynamic
Target number is 3.1, and true moving-target number is 3, so the more traditional estimation to moving-target true number of the present invention
More acurrate, curve is closer, improves motion target tracking accuracy rate；GMPHD filter is based on as Fig. 6 b present invention is improved
Circumference SAR motion target tracking method more can accurately estimate target number, and OSPA distance is respectively less than 2m, gradually converges to steady
Definite value, overall estimation effect are more preferable.
Detailed description of the invention
Fig. 1 is binary channels circumference SAR system geometrical relationship figure, Ri(tm) be 2 channel distance point targets instantaneous oblique square,
PtFor the position for observing scene moving target, RcFor radar to the distance at observation area center, RgFor radar uniform motion circumference
Radius note, ZhFor carrier of radar height, spacing of the d between antenna, V is that radar platform moves linear velocity；
Fig. 2 is based on the circumference SAR motion target tracking method flow chart for improving GMPHD filter；
Fig. 3 is target true motion track schematic diagram in emulation experiment of the present invention；
Fig. 4 a is the 10th frame SAR imaging results figure of single channel echo-signal in emulation experiment of the present invention；
Fig. 4 b is the 10th frame moving target Preliminary detection result figure in emulation experiment of the present invention；
Fig. 5 is that improvement GMPHD target position tracking result figure in x-axis and y-axis direction is utilized in emulation experiment of the present invention；
Fig. 6 a is to filter to carry out 500 Meng Teka with improvement GMPHD filtering using Traditional GM PHD in emulation experiment of the present invention
The target number estimated result figure of Lip river emulation；
Fig. 6 b is to filter to carry out 500 Meng Teka with improvement GMPHD filtering using Traditional GM PHD in emulation experiment of the present invention
The OSPA distance results figure of Lip river emulation.
Specific embodiment
Specific embodiment 1: circumference SAR motion target tracking method of the present embodiment based on GMPHD filter is specific
Process are as follows:
In order to make the foregoing objectives, features and advantages of the present invention clearer and more comprehensible, below with reference to Fig. 1, Fig. 2, Fig. 3,
The present invention is described in further detail for Fig. 4 a, Fig. 4 b, Fig. 5, Fig. 6 a, Fig. 6 b and specific embodiment；
Wherein Fig. 2 is the stream of the present invention based on the circumference SAR motion target tracking method for improving GMPHD filter
Cheng Tu.
Step 1: Preliminary detection is carried out to moving target based on DPCA-CFAR method, obtains the observation of moving target；
Step 2: the moving target observation obtained according to step 1, establish moving target state vector and measure to
Amount；
Step 3: state vector and measurement vector based on moving target establish the moving target shape based on CSAR system
State equation and measurement equation；
Step 4: Gaussian component mark value initialization；
Step 5: state equation and measurement equation based on moving target generate adaptive newborn target strength, and are new
The new Gauss mark value of raw Target Assignment；
Step 6: the Gauss mark value of adaptive newborn target strength and newborn target based on step 5 carries out GMPHD
Recursion；
GMPHD recursive process includes PHD prediction, PHD update, the trimming and fusion of Gaussian term；
Step 7: moving target state estimation is carried out based on step 6 and track extracts, completes SAR motion target tracking.
Specific embodiment 2: the present embodiment is different from the first embodiment in that, it is based in the step 1
DPCA-CFAR method carries out Preliminary detection to moving target, obtains the observation of moving target；Detailed process are as follows:
Dual-Channel SAR system: carrier of radar height is denoted as z=Zh, radar uniform motion radius of a circle is denoted as Rg, radar is flat
It is V that platform, which moves linear velocity, and the distance at radar to observation area center isRadar beam at the center of observation area
Incidence angle isM indicates incidence angle；
Binary channels circumference SAR is horizontally arranged 2 width antennas along radar track tangential direction, and the spacing between antenna is d, double
Channel circumference SAR is knocked off operation mode using single-shot 2；
Based on displaced phase center principle, it is necessary first to carry out time calibration and phase to received two channel image signal
Compensation carries out DPCA processing to two channel images after time calibration, realizes clutter recognition；
After clutter recognition, CFAR Preliminary detection is carried out to each frame image and obtains target measuring value, in CFAR detection
The middle higher false-alarm probability P of settingfa, to ensure based on higher detection probability PdObtain the observation of moving target；Specific mistake
Journey are as follows:
CFAR Preliminary detection result is expressed as
Wherein, ηkIndicate detection threshold, IDPCA,k(x, y) indicates two channel image DPCA processing results；
Detection threshold ηkBy false-alarm probability PfaWith two channel image DPCA processing result IDPCA,kThe probability statistics of (x, y)
It determines；
The false-alarm probability Pfa=0.01；
To CFAR Preliminary detection result Bk(x, y) carries out image procossing:
The target area that connection is obtained by morphologic corrosion expansive working, carries out edge detection to target area and obtains
The center point for taking target, using the center of target point as the observation of moving target.
Binary channels circumference SAR system geometrical relationship is as shown in figure.
Other steps and parameter are same as the specific embodiment one.
Specific embodiment 3: the present embodiment is different from the first and the second embodiment in that, it is described to be based on equivalent phase
Position center principle, it is necessary first to time calibration and phase compensation be carried out to received two channel image signal, after time calibration
Two channel images carry out DPCA processing, realize clutter recognition；Detailed process are as follows:
If two channel image sequences after time calibration are expressed as I1,k(x, y) and I2,k(x,y)；
Wherein k indicates frame number, k=1,2 ..., K, and K value is positive integer；
It is in two channel image DPCA processing results of kth frame
Wherein λ indicates wavelength；J is imaginary number, j2=-1；vrIndicate the radial velocity of moving target.
Other steps and parameter are the same as one or two specific embodiments.
Specific embodiment 4: unlike one of present embodiment and specific embodiment one to three, the step 2
The middle moving target observation obtained according to step 1 establishes the state vector of moving target and measures vector；Detailed process are as follows:
[xk,yk] indicate moving target actual position,Indicate moving target in cartesian coordinate system x-axis and y-axis
The true velocity in direction, then the state vector of moving target is expressed as
It measures vector and is expressed as zk=[x'k,y'k]T,
Wherein [x'k,y'k] indicate position of the moving target in SAR image rather than locations of real targets, this is because mesh
Mark radial velocity causes target that orientation offset occurs；T indicates transposition；
There are a moving target of N (k) and a measuring value of M (k) in kth frame, then the state vector set of moving target and
The form of stochastic finite set can be regarded as by measuring vector set, be expressed as
Wherein, lk,1Indicate the state vector of first moving target in kth frame, lk,N(k)Indicate that N (k) is a in kth frame
The state vector of moving target；zk,1Indicate the measurement vector of first moving target in kth frame, zk,M(k)Indicate M in kth frame
(k) the measurement vector of a moving target；WithThe respectively state space of moving targetWith measurement space
The set that upper all finite subsets are constituted.
For the multiple target state X at given k-1 momentk-1, each dbjective state xk-1∈Xk-1With Probability pS,k-1(xk-1)
It is transferred to new state xk, or with probability 1-pS,k-1(xk-1) disappear.Therefore, for the state x at given k-1 momentk∈Xk-1,
The state of subsequent time can use stochastic finite collection Sk|k-1(xk-1) indicate, when target exists, Sk|k-1(xk-1) it is { xk,
When target disappears, Sk|k-1(xk-1) it is empty set.
Other steps and parameter are identical as one of specific embodiment one to three.
Specific embodiment 5: unlike one of present embodiment and specific embodiment one to four, the step 3
In state vector based on moving target and measure vector, establish the moving target state equation based on CSAR system and measurement side
Journey；
The state equation of moving target is expressed as
lk=Fk-1lk-1+vk-1
Wherein, process noise vector vk~N (0, Qk), covarianceσvIndicate the standard deviation of process noise,
GkFor process noise distribution matrix, FkFor dynamic model state-transition matrix, t is the sampling period, for linear uniform motion mesh
Mark；
FkAnd GkIt is expressed as
Measurement equation is expressed as
zk=Hklk+wk
Wherein measure noise wk~N (0, Rk), covarianceσwIndicate the standard deviation of measurement noise, I2It is 2 × 2
Unit matrix, HkIndicate observing matrix；
The orientation offset as caused by radial velocity is
Wherein rkIndicate oblique distance of the radar platform to target, vr,kIndicate the radial velocity of kth frame moving target,θa,kIndicate the azimuth of kth frame radar platform；
Based on orientation offset Δk, position of the moving target in SAR image is
x'k=xk-Δk·sinθa,k
y'k=yk+Δk·cosθa,k
Since the initial oblique distance of radar to moving target is larger, oblique distance change over time it is smaller, therefore for microinching
Target puts aside that oblique distance caused by target and radar motion changes in analysis, i.e. hypothesis rk=R enables A=Rsin θm/ V, R
Indicate oblique distance constant of the radar platform to target, A expression intermediate variable；
Then
Observing matrix is expressed as
Other steps and parameter are identical as one of specific embodiment one to four.
Specific embodiment 6: unlike one of present embodiment and specific embodiment one to five, the step 4
Middle Gaussian component mark value initialization；Detailed process are as follows:
In improving GMPHD filter, is extracted to complete targetpath, each Gaussian component is marked, is being realized
The output of track is completed after the Target state estimator of all frames according to mark value, the Gaussian component of distribution same tag value indicates same
The point mark of one target.
Initial strength is J in initial frame0The sum of a Gaussian component distributes a mark value, note for each Gaussian component
For
In formula,For first Gaussian component,For J0A Gaussian component, L0For the set of Gaussian component mark value.
Other steps and parameter are identical as one of specific embodiment one to five.
Specific embodiment 7: unlike one of present embodiment and specific embodiment one to six, the step 5
In state equation and measurement equation based on moving target, generate adaptive newborn target strength, and new for newborn Target Assignment
Gauss mark value；Detailed process are as follows:
Step 5 one: binary channels interferometric phase Δ φ is calculated based on ATI technologyk, to realize to target bearing to offset
Estimation
In formula, λ is wavelength；R indicates radar platform to the oblique distance constant of target；
Obtain compensated measuring valueWithIt is denoted as
Step 5 two: fast Track Initiation is carried out using logical approach, using all compensated measuring values of kth frame as new
Track it is assumed that carry out track initiation in this, as track head, track is realized by 3/4 logic theory of track initiation sliding window method
Confirmation, it is determined whether have newborn track；If there is newborn track, confirmation has newborn target, executes step 5 three, otherwise executes step
Rapid six；
Step 5 three: when there is Jγ,kA new life track, constructs newborn target strength
Utilize compensated measuring value in step 5 oneWithCalculate the weight of newborn target strengthMean valueAnd covariance
lkFor stochastic variable；For Gaussian Profile；
The weight of newborn target strengthMean valueAnd covarianceIt is expressed as
Wherein VmaxIndicate the detectable speed of maximum, w indicates the total weight of newborn intensity, and jj indicates that k frame jth j is a compensated
Observation, ii indicate that i-th i newborn track, ll indicate k-1 frame ll compensated observations；
So far the estimation to newborn intensity, and the Gauss mark value new for newborn Target Assignment are completed:
Indicate the Gaussian component mark value of first newborn target；Indicate Jγ,kThe Gauss point of a new life target
Measure mark value；
Other steps and parameter are identical as one of specific embodiment one to six.
Specific embodiment 8: unlike one of present embodiment and specific embodiment one to seven, the step 6
In adaptive newborn target strength and newborn target based on step 5 Gauss mark value, carry out GMPHD recursion；GMPHD is passed
Pushing through journey includes PHD prediction, PHD update, the trimming and fusion of Gaussian term；
Detailed process are as follows:
Step 6 one: the Gauss mark value of adaptive newborn target strength and newborn target based on step 5 carries out PHD
Prediction, the predicted intensity of kth frame are expressed as
Wherein, pS,kIndicate the survival probability of target, Jk-1For the number of -1 vertical frame dimension this mixed term of kth, Indicate the mean value of Gaussian term after PHD is predicted, Fk-1Table
Show dynamic model state-transition matrix,Indicate the mean value in each Gaussian term of -1 frame of kth,After indicating PHD prediction
The variance of Gaussian term, Qk-1Indicate process noise vector covariance,Indicate the variance in each Gaussian term of -1 frame of kth,
Indicate that the weight of each Gaussian term, T indicate transposition；
The Gauss mark value of predicted intensity is
Lk|k-1=Lk-1∪Lγ,k
In formula, Lk|k-1Gauss mark value after indicating prediction, Lk-1Indicate the Gauss mark value of k-1 frame, Lγ,kIndicate newborn
The Gauss mark value of target；
Step 6 two: PHD updates: according to the predicted intensity v of kth framek|k-1(lk) and measurement vector set Zk, calculate kth
The posteriority intensity Gaussian Mixture form v of framek(lk)；
Assuming that the predicted intensity of kth frame is Gaussian Mixture form, then the posteriority intensity Gaussian Mixture form v of kth framek(lk)
It indicates are as follows:
In formula, pD,kIndicate target detection probability, Jk|k-1=Jk-1+Jγ,kIndicate of Gaussian Mixture item in predicted intensity
Number, with stochastic finite collection Sk|k-1(xk-1) indicate, vk|k-1(lk) indicate kth frame predicted intensity,Each Gaussian term
Weight,Indicate the mean value of q Gaussian term after kth frame PHD updates,Q-th after expression kth frame PHD update
The variance of Gaussian term；
vk|k-1(lk) in each Gaussian component obtain 1+ after PHD updates | Zk| a update Gaussian component, | Zk| it indicates to measure
It is worth number；They represent the same target, therefore this 1+ | Zk| a mark value and the prediction Gaussian component for updating Gaussian component
Mark value it is identical, then updated vk(lk) Gaussian component mark value are as follows:
In formula, Lk|k- 1 indicates the Gauss mark value after prediction, LkIndicate the Gauss mark value of k frame,It indicates by z1More
The Gauss mark value newly obtainedIndicate byUpdate obtained Gauss mark value Indicate the | Zk| a measurement vector；
Step 6 three: to vk(lk) in Gaussian componentIt is trimmed and is merged to have
Effect reduces the number of Gaussian Mixture item, during being merged；Detailed process are as follows:
If the v being mergedk(lk) in Gaussian componentMark value having the same, then close
Gaussian term v after andk(xk) in Gaussian componentMark value with merge before Gaussian term vk
(lk) in Gaussian componentMark value is identical, on the contrary, if the v being mergedk(lk) in height
This componentMark value it is different, then take the Gaussian term v before merging with maximum weightk
(xk) in Gaussian componentMark value as merge after Gaussian term mark value.
After trimming union operation, if being needed to them again there is also different Gaussian term mark values is identical
Mark value is distributed, the mark value of the maximum Gaussian term of weight is retained, and the Gaussian term of other same label is then assigned other only
One mark value.
Other steps and parameter are identical as one of specific embodiment one to seven.
Specific embodiment 9: unlike one of present embodiment and specific embodiment one to eight, the step 6
The weight of each Gaussian term in twoThe mean value of q Gaussian term after kth frame PHD updatesKth frame PHD updates
The variance of q-th of Gaussian term afterwardsExpression formula are as follows:
WhereinIndicate the gain matrix in q-th of Gaussian term of kth frame,Indicate q-th of Gaussian term in predicted intensity
Variance,Indicating the mean value of q-th of Gaussian term in predicted intensity, I indicates 4 × 4 unit matrix,Indicate predicted intensity
In q-th of Gaussian term weight,Indicate the likelihood function of q-th of Gaussian term,It indicates in predicted intensity first
The weight of Gaussian term,Indicate that the likelihood function of first of Gaussian term, l indicate Gauss item number in predicted intensity,It indicates
Vector is measured,For noise intensity.
Other steps and parameter are identical as one of specific embodiment one to eight.
Specific embodiment 10: unlike one of present embodiment and specific embodiment one to nine, the step 7
In moving target state estimation and track carried out based on step 6 extract, complete SAR motion target tracking；Detailed process are as follows:
The posteriority intensity Gaussian Mixture form v of kth frame is taken in each framek(lk) inGaussian component it is equal
Value is the state estimation of target, in the Gaussian term v of all frame k=1 ..., K distribution same tag valuek(lk) in Gaussian componentIndicate the point mark of same target；
WhereinIndicate the posteriority intensity Gaussian Mixture form v in kth framek(lk) in Gaussian term weight.
Other steps and parameter are identical as one of specific embodiment one to nine.
Beneficial effects of the present invention are verified using following embodiment:
Embodiment one:
In the present embodiment one based on improve GMPHD filter circumference SAR motion target tracking method be specifically according to
Lower step preparation:
Emulation experiment
1 binary channels circumference SAR system simulation parameter
Binary channels circumference SAR system simulation parameter is as shown in table 1, is horizontally arranged 2 width antennas along radar track tangential direction,
Using single-emission and double-receiving operating mode, in the imaged scene be provided with 3 linear uniform motion targets, be expressed as MTi (i=1,2,
3), as shown in figure 3, target successive appearing and subsiding in image sequence.Transmitted signal bandwidth and pulse recurrence frequency distinguish table
It is shown as B and PRF." o " indicates the movement initial position of target in single channel imaging results and Preliminary detection result, and " x " indicates mesh
Target moves end position.Signal-to-noise ratio (SNR) and signal to noise ratio (SCR) are respectively set to 10dB and 0dB, in preliminary CFAR detection
Middle false-alarm probability is 0.01.Fig. 4 a and Fig. 4 b are respectively indicated after single-channel SAR imaging results and DPCA-CFAR processing at the beginning of moving-target
Walk testing result, it can be seen that moving-target, which is submerged in clutter, in single-channel SAR image is difficult to, through DPCA-CFAR
It can detecte moving target after processing, while part clutter is also detected.
1 radar simulation parameter of table
2 improve the setting of GMPHD filter parameter
The setting of GMPHD filter parameter used in this emulation experiment is as shown in table 2.
Table 2 improves the setting of GMPHD filter parameter
3 experiment simulation figures are shown in attached drawing 5,6a, 6b.
Such as Fig. 5, different targetpaths is indicated with different lines, it can be seen that each dbjective state base of each frame estimation
On this on corresponding real goal track, the motion profile of 5 targets can be accurately extracted, with the true track kissing of target
It closes, and eliminates clutter residue.
If Fig. 6 a is in allowable range of error, by emulation it can be seen that the present invention is improved based on GMPHD filter
The moving-target number that circumference SAR motion target tracking method determines is 3.05, and Traditional GM PHD method for tracking target determines dynamic
Target number is 3.1, and true moving-target number is 3, so the more traditional estimation to moving-target true number of the present invention
More acurrate, curve is closer, improves motion target tracking accuracy rate；GMPHD filter is based on as Fig. 6 b present invention is improved
Circumference SAR motion target tracking method more can accurately estimate target number, and OSPA distance is respectively less than 2m, gradually converges to steady
Definite value, overall estimation effect are more preferable.
The present invention can also have other various embodiments, without deviating from the spirit and substance of the present invention, this field
Technical staff makes various corresponding changes and modifications in accordance with the present invention, but these corresponding changes and modifications all should belong to
The protection scope of the appended claims of the present invention.
Claims (10)
1. the circumference SAR motion target tracking method based on GMPHD filter, it is characterised in that: the method detailed process are as follows:
Step 1: Preliminary detection is carried out to moving target based on DPCA-CFAR method, obtains the observation of moving target；
Step 2: the moving target observation obtained according to step 1 establishes the state vector of moving target and measures vector；
Step 3: state vector and measurement vector based on moving target establish the moving target state side based on CSAR system
Journey and measurement equation；
Step 4: Gaussian component mark value initialization；
Step 5: state equation and measurement equation based on moving target generate adaptive newborn target strength, and are newborn mesh
Mark distributes new Gauss mark value；
Step 6: the Gauss mark value of adaptive newborn target strength and newborn target based on step 5 carries out GMPHD and passs
It pushes away；
GMPHD recursive process includes PHD prediction, PHD update, the trimming and fusion of Gaussian term；
Step 7: moving target state estimation is carried out based on step 6 and track extracts, completes SAR motion target tracking.
2. the circumference SAR motion target tracking method based on GMPHD filter according to claim 1, it is characterised in that: institute
It states in step 1 and Preliminary detection is carried out to moving target based on DPCA-CFAR method, obtain the observation of moving target；Specific mistake
Journey are as follows:
Dual-Channel SAR system: carrier of radar height is denoted as z=Zh, radar uniform motion radius of a circle is denoted as Rg, radar platform fortune
Moving-wire speed is V, and the distance at radar to observation area center isRadar beam is incident at the center of observation area
Angle isM indicates incidence angle；
Binary channels circumference SAR is horizontally arranged 2 width antennas along radar track tangential direction, and the spacing between antenna is d, binary channels
Circumference SAR is knocked off operation mode using single-shot 2；
Time calibration and phase compensation are carried out to received two channel image signal, two channel images after time calibration are carried out
DPCA processing, realizes clutter recognition；
After clutter recognition, CFAR Preliminary detection is carried out to each frame image and obtains target measuring value, is set in CFAR detection
Set false-alarm probability Pfa, obtain the observation of moving target；Detailed process are as follows:
CFAR Preliminary detection result is expressed as
Wherein, ηkIndicate detection threshold, IDPCA,k(x, y) indicates two channel image DPCA processing results；
To CFAR Preliminary detection result Bk(x, y) carries out image procossing:
The target area that connection is obtained by morphologic corrosion expansive working, carries out edge detection to target area and obtains mesh
Target center point, using the center of target point as the observation of moving target.
3. the circumference SAR motion target tracking method according to claim 1 or claim 2 based on GMPHD filter, feature exist
In: it is described that time calibration and phase compensation are carried out to received two channel image signal, to two channel images after time calibration
DPCA processing is carried out, realizes clutter recognition；Detailed process are as follows:
If two channel image sequences after time calibration are expressed as I1,k(x, y) and I2,k(x,y)；
Wherein k indicates frame number, k=1,2 ..., K, and K value is positive integer；
It is in two channel image DPCA processing results of kth frame
Wherein λ indicates wavelength；J is imaginary number, j2=-1；vrIndicate the radial velocity of moving target.
4. the circumference SAR motion target tracking method based on GMPHD filter according to claim 3, it is characterised in that: institute
The moving target observation obtained in step 2 according to step 1 is stated, the state vector of moving target is established and measures vector；Tool
Body process are as follows:
[xk,yk] indicate moving target actual position,Indicate moving target in cartesian coordinate system x-axis and y-axis direction
True velocity, then the state vector of moving target is expressed as
It measures vector and is expressed as zk=[x'k,y'k]T,
Wherein [x'k,y'k] indicate position of the moving target in SAR image；T indicates transposition；
There are a moving target of N (k) and a measuring value of M (k) in kth frame, then the state vector set and measurement of moving target
Vector set is expressed as
Wherein, lk,1Indicate the state vector of first moving target in kth frame, lk,N(k)Indicate a movement of N (k) in kth frame
The state vector of target；zk,1Indicate the measurement vector of first moving target in kth frame, zk,M(k)Indicate M (k) in kth frame
The measurement vector of a moving target；WithThe respectively state space of moving targetWith measurement spaceUpper institute
The set for thering is finite subset to constitute.
5. the circumference SAR motion target tracking method based on GMPHD filter according to claim 4, it is characterised in that: institute
It states the state vector in step 3 based on moving target and measures vector, establish the moving target state equation based on CSAR system
And measurement equation；
The state equation of moving target is expressed as
lk=Fk-1lk-1+vk-1
Wherein, process noise vector vk~N (0, Qk), covarianceσvIndicate the standard deviation of process noise, GkFor
Process noise distribution matrix, FkFor dynamic model state-transition matrix, t is the sampling period, for linear uniform motion target；
FkAnd GkIt is expressed as
Measurement equation is expressed as
zk=Hklk+wk
Wherein measure noise wk~N (0, Rk), covarianceσwIndicate the standard deviation of measurement noise, I2For 2 × 2 list
Position battle array, HkIndicate observing matrix；
The orientation offset as caused by radial velocity is
Wherein rkIndicate oblique distance of the radar platform to target, vr,kIndicate the radial velocity of kth frame moving target,θa,kIndicate the azimuth of kth frame radar platform；
Based on orientation offset Δk, position of the moving target in SAR image is
x'k=xk-Δk·sinθa,k
y'k=yk+Δk·cosθa,k
Assuming that rk=R enables A=Rsin θm/ V, R indicate radar platform to the oblique distance constant of target, A expression intermediate variable；
Then
Observing matrix is expressed as
6. the circumference SAR motion target tracking method based on GMPHD filter according to claim 5, it is characterised in that: institute
Gaussian component mark value in step 4 is stated to initialize；Detailed process are as follows:
Initial strength is J in initial frame0The sum of a Gaussian component is distributed a mark value for each Gaussian component, is denoted as
In formula,For first Gaussian component,For J0A Gaussian component, L0For the set of Gaussian component mark value.
7. the circumference SAR motion target tracking method based on GMPHD filter according to claim 6, it is characterised in that: institute
The state equation and measurement equation in step 5 based on moving target are stated, adaptive newborn target strength is generated, and is newborn mesh
Mark distributes new Gauss mark value；Detailed process are as follows:
Step 5 one: binary channels interferometric phase Δ φ is calculated based on ATI technologyk, to realize the estimation to target bearing to offset
In formula, λ is wavelength；R indicates radar platform to the oblique distance constant of target；
Obtain compensated measuring valueWithIt is denoted as
Step 5 two: carrying out track initiation for all compensated measuring values of kth frame as track head, sliding by track initiation
Window method realizes track confirmation, it is determined whether has newborn track；If there is newborn track, confirmation has newborn target, executes step 5
Three, otherwise execute step 6；
Step 5 three: when there is Jγ,kA new life track, constructs newborn target strength
Utilize compensated measuring value in step 5 oneWithCalculate the weight of newborn target strengthMean valueWith
Covariance
lkFor stochastic variable；For Gaussian Profile；
The weight of newborn target strengthMean valueAnd covarianceIt is expressed as
Wherein VmaxIndicate maximum detectable speed, w indicates the total weight of newborn intensity, and jj indicates k frame jth j compensated observations
Value, ii indicate that i-th i newborn track, ll indicate k-1 frame ll compensated observations；
So far the estimation to newborn intensity, and the Gauss mark value new for newborn Target Assignment are completed:
Indicate the Gaussian component mark value of first newborn target；Indicate Jγ,kThe Gaussian component mark of a new life target
Note value.
8. the circumference SAR motion target tracking method based on GMPHD filter according to claim 7, it is characterised in that: institute
The Gauss mark value for stating the adaptive newborn target strength and newborn target in step 6 based on step 5, carries out GMPHD recursion；
GMPHD recursive process includes PHD prediction, PHD update, the trimming and fusion of Gaussian term；Detailed process are as follows:
Step 6 one: it is pre- to carry out PHD for the Gauss mark value of adaptive newborn target strength and newborn target based on step 5
It surveys, the predicted intensity of kth frame is expressed as
Wherein, pS,kIndicate the survival probability of target, Jk-1For the number of -1 vertical frame dimension this mixed term of kth, Indicate the mean value of Gaussian term after PHD is predicted, Fk-1Indicate that dynamic model state turns
Matrix is moved,Indicate the mean value in each Gaussian term of -1 frame of kth,Indicate the variance of Gaussian term after PHD is predicted, Qk-1
Indicate process noise vector covariance,Indicate the variance in each Gaussian term of -1 frame of kth,Indicate each Gaussian term
Weight, T indicate transposition；
The Gauss mark value of predicted intensity is
Lk|k-1=Lk-1∪Lγ,k
In formula, Lk|k-1Gauss mark value after indicating prediction, Lk-1Indicate the Gauss mark value of k-1 frame, Lγ,kIndicate newborn target
Gauss mark value；
Step 6 two: PHD updates:
According to the predicted intensity v of kth framek|k-1(lk) and measurement vector set Zk, calculate the posteriority intensity Gaussian Mixture shape of kth frame
Formula vk(lk)；
Assuming that the predicted intensity of kth frame is Gaussian Mixture form, then the posteriority intensity Gaussian Mixture form v of kth framek(lk) indicate
Are as follows:
In formula, pD,kIndicate target detection probability, Jk|k-1=Jk-1+Jγ,kIndicate the number of Gaussian Mixture item in predicted intensity,
vk|k-1(lk) indicate kth frame predicted intensity,The weight of each Gaussian term,Indicate that kth frame PHD updates
The mean value of q Gaussian term afterwards,Indicate the variance of q-th of Gaussian term after kth frame PHD updates；
vk|k-1(lk) in each Gaussian component obtain 1+ after PHD updates | Zk| a update Gaussian component, | Zk| indicate measuring value
Number；
Then updated vk(lk) Gaussian component mark value are as follows:
In formula, Lk|k-1Gauss mark value after indicating prediction, LkIndicate the Gauss mark value of k frame,It indicates by z1It updates
The Gauss mark value arrived Indicate byUpdate obtained Gauss mark value Table
Show | Zk| a measurement vector；
Step 6 three: to vk(lk) in Gaussian componentIt is trimmed and is merged；Detailed process
Are as follows:
If the v being mergedk(lk) in Gaussian componentMark value having the same, then after merging
Gaussian term mark value with merge before Gaussian term mark value it is identical,
On the contrary, if the v being mergedk(lk) in Gaussian componentMark value it is different, then take conjunction
And mark value of the mark value of the Gaussian term before with maximum weight as Gaussian term after merging.
9. the circumference SAR motion target tracking method based on GMPHD filter according to claim 8, it is characterised in that: institute
State the weight of each Gaussian term in step 6 twoThe mean value of q Gaussian term after kth frame PHD updatesKth
The variance of q-th of Gaussian term after frame PHD updatesExpression formula are as follows:
WhereinIndicate the gain matrix in q-th of Gaussian term of kth frame,Indicate q-th Gaussian term in predicted intensity
Variance,Indicating the mean value of q-th of Gaussian term in predicted intensity, I indicates 4 × 4 unit matrix,It indicates in predicted intensity
The weight of q-th of Gaussian term,Indicate the likelihood function of q-th of Gaussian term,Indicate first high in predicted intensity
This weight,Indicate that the likelihood function of first of Gaussian term, l indicate Gauss item number in predicted intensity,Expression amount
Direction finding amount,For noise intensity.
10. the circumference SAR motion target tracking method based on GMPHD filter according to claim 9, it is characterised in that:
Moving target state estimation is carried out based on step 6 in the step 7 and track extracts, completes SAR motion target tracking；Specifically
Process are as follows:
The posteriority intensity Gaussian Mixture form v of kth frame is taken in each framek(lk) inThe mean value of Gaussian component be
The state estimation of target indicates the point mark of same target in the Gaussian term of all frame k=1 ..., K distribution same tag value；
WhereinIndicate the posteriority intensity Gaussian Mixture form v in kth framek(lk) in Gaussian term weight.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910222890.3A CN109934840A (en) | 2019-03-22 | 2019-03-22 | Circumference SAR motion target tracking method based on GMPHD filter |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910222890.3A CN109934840A (en) | 2019-03-22 | 2019-03-22 | Circumference SAR motion target tracking method based on GMPHD filter |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN109934840A true CN109934840A (en) | 2019-06-25 |
Family
ID=66988022
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201910222890.3A Pending CN109934840A (en) | 2019-03-22 | 2019-03-22 | Circumference SAR motion target tracking method based on GMPHD filter |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN109934840A (en) |
Cited By (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111062972A (en) * | 2019-12-25 | 2020-04-24 | 泉州装备制造研究所 | Image tracking method based on image frequency domain conversion |
Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN105842693A (en) * | 2016-03-23 | 2016-08-10 | 哈尔滨工业大学 | Method for detecting dual-channel SAR moving target based on compressive sensing |
| CN107102293A (en) * | 2017-04-25 | 2017-08-29 | 杭州电子科技大学 | The passive co-located method of unknown clutter estimated based on sliding window integral density |
| CN108304869A (en) * | 2018-01-29 | 2018-07-20 | 西北工业大学 | A kind of fusion method and system of the comprehensive magnitude information based on multiple sensors |
-
2019
- 2019-03-22 CN CN201910222890.3A patent/CN109934840A/en active Pending
Patent Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN105842693A (en) * | 2016-03-23 | 2016-08-10 | 哈尔滨工业大学 | Method for detecting dual-channel SAR moving target based on compressive sensing |
| CN107102293A (en) * | 2017-04-25 | 2017-08-29 | 杭州电子科技大学 | The passive co-located method of unknown clutter estimated based on sliding window integral density |
| CN108304869A (en) * | 2018-01-29 | 2018-07-20 | 西北工业大学 | A kind of fusion method and system of the comprehensive magnitude information based on multiple sensors |
Non-Patent Citations (7)
| Title |
|---|
| YUN ZHANG 等: "Moving Target Detection and Tracking Based on Gmphd Filter in SAR System", 《IGARSS 2018 - 2018 IEEE INTERNATIONAL GEOSCIENCE AND REMOTE SENSING SYMPOSIUM》 * |
| YUN ZHANG 等: "Moving Target Tracking Based on Improved GMPHD Filter in Circular SAR System", 《IEEE GEOSCIENCE AND REMOTE SENSING LETTERS》 * |
| 代大海: "《数字阵列合成孔径雷达》", 31 December 2017 * |
| 王力宝: "基于DPCA和多分辨分析的SAR微弱运动目标检测", 《电子与信息学报》 * |
| 王肖洋: "基于DPCA的SAR地面运动目标检测技术研究", 《中国优秀硕士学位论文全文数据库信息科技辑》 * |
| 鲁加国: "《合成孔径雷达设计技术》", 30 April 2017 * |
| 黄炎龙: "合成孔径雷达地面运动目标检测方法研究", 《中国优秀硕士学位论文全文数据库信息科技辑》 * |
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111062972A (en) * | 2019-12-25 | 2020-04-24 | 泉州装备制造研究所 | Image tracking method based on image frequency domain conversion |
| CN111062972B (en) * | 2019-12-25 | 2023-11-21 | 泉州装备制造研究所 | Image tracking method based on image frequency domain conversion |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN100365429C (en) | Motive target imaging method of synthetic aperture radar | |
| CN106842165B (en) | Radar centralized asynchronous fusion method based on different distance angular resolutions | |
| CN105487060B (en) | A kind of multiple target extracting method of four Slope Modulation of binary channels | |
| CN107576959B (en) | High repetition frequency radar target tracking method before detection based on area mapping deblurring | |
| AU2002308569A1 (en) | Systems and method for central association and tracking in passive coherent location applications | |
| WO2002091016A2 (en) | System and method for central association and tracking in passive coherent location applications | |
| CN104502907A (en) | Stable ground moving/static target tracking method for airborne radar | |
| CN109946694A (en) | Circumference SAR multi-object tracking method based on stochastic finite collection | |
| CN105093215A (en) | Doppler information based method for tracking low-altitude low-speed small target through radar | |
| CN112184749B (en) | Moving target tracking method based on video SAR cross-domain combination | |
| CN106646466B (en) | A kind of imaging method of the weighting back-projection algorithm based on principal component analysis | |
| CN107515398A (en) | Inverse synthetic aperture radar imaging target continuous tracking method based on wideband echoes | |
| Zhang et al. | Moving target tracking based on improved GMPHD filter in circular SAR system | |
| Zhou et al. | Multiple-kernelized-correlation-filter-based track-before-detect algorithm for tracking weak and extended target in marine radar systems | |
| CN106468772A (en) | A kind of multistation radar human body tracing method based on range Doppler measurement | |
| CN109100697A (en) | A kind of target condensing method based on ground surveillance radar system | |
| CN109934840A (en) | Circumference SAR motion target tracking method based on GMPHD filter | |
| CN106772357B (en) | AI-PHD filter multi-object tracking method under signal-to-noise ratio unknown condition | |
| Yu et al. | A novel RHT-TBD approach for weak targets in HPRF radar | |
| CN108828584B (en) | Multi-frequency target tracking-before-detection method based on track folding factor ambiguity resolution | |
| Li et al. | Radar-based UAV swarm surveillance based on a two-stage wave path difference estimation method | |
| CN113093174B (en) | PHD filter radar fluctuation weak multi-target-based pre-detection tracking method | |
| Zhou et al. | High-squint SAR imaging for noncooperative moving ship target based on high velocity motion platform | |
| Gao et al. | An improved TBD algorithm based on dynamic programming for dim SAR target detection | |
| CN114185047A (en) | Bistatic SAR moving target refocusing method based on optimal polar coordinate transformation |
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 | ||
| RJ01 | Rejection of invention patent application after publication | ||
| RJ01 | Rejection of invention patent application after publication |
Application publication date: 20190625 |