CN114577213B - Distributed multi-platform underwater multi-target association and passive positioning method - Google Patents
Distributed multi-platform underwater multi-target association and passive positioning method Download PDFInfo
- Publication number
- CN114577213B CN114577213B CN202210198348.0A CN202210198348A CN114577213B CN 114577213 B CN114577213 B CN 114577213B CN 202210198348 A CN202210198348 A CN 202210198348A CN 114577213 B CN114577213 B CN 114577213B
- Authority
- CN
- China
- Prior art keywords
- target
- platform
- association
- time
- observed
- 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
Links
Images
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01C—MEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
- G01C21/00—Navigation; Navigational instruments not provided for in groups G01C1/00 - G01C19/00
- G01C21/20—Instruments for performing navigational calculations
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/18—Complex mathematical operations for evaluating statistical data, e.g. average values, frequency distributions, probability functions, regression analysis
-
- H—ELECTRICITY
- H03—ELECTRONIC CIRCUITRY
- H03H—IMPEDANCE NETWORKS, e.g. RESONANT CIRCUITS; RESONATORS
- H03H17/00—Networks using digital techniques
- H03H17/02—Frequency selective networks
- H03H17/0248—Filters characterised by a particular frequency response or filtering method
- H03H17/0255—Filters based on statistics
- H03H17/0257—KALMAN filters
-
- 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
- Y02D—CLIMATE CHANGE MITIGATION TECHNOLOGIES IN INFORMATION AND COMMUNICATION TECHNOLOGIES [ICT], I.E. INFORMATION AND COMMUNICATION TECHNOLOGIES AIMING AT THE REDUCTION OF THEIR OWN ENERGY USE
- Y02D30/00—Reducing energy consumption in communication networks
- Y02D30/70—Reducing energy consumption in communication networks in wireless communication networks
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Mathematical Physics (AREA)
- Data Mining & Analysis (AREA)
- Computational Mathematics (AREA)
- Mathematical Analysis (AREA)
- Mathematical Optimization (AREA)
- Pure & Applied Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Remote Sensing (AREA)
- Radar, Positioning & Navigation (AREA)
- Software Systems (AREA)
- Algebra (AREA)
- Databases & Information Systems (AREA)
- Probability & Statistics with Applications (AREA)
- General Engineering & Computer Science (AREA)
- Computer Hardware Design (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Bioinformatics & Computational Biology (AREA)
- Evolutionary Biology (AREA)
- Operations Research (AREA)
- Life Sciences & Earth Sciences (AREA)
- Computing Systems (AREA)
- Automation & Control Theory (AREA)
- Measurement Of Velocity Or Position Using Acoustic Or Ultrasonic Waves (AREA)
Abstract
A distributed multi-platform underwater multi-target association and passive positioning method belongs to the technical field of distributed multi-platform underwater multi-target positioning and tracking. The invention solves the problems of low positioning precision and low accuracy of multi-target association of the existing maneuvering target positioning and tracking method. The invention introduces a target line spectrum characteristic set into a target association method, provides a target association algorithm based on frequency characteristic multi-information vector fusion, solves the problems that a pure orientation system cannot acquire a flight path and the flight path is difficult to associate, and can give a correct association result. Compared with the traditional nonlinear Bayesian filtering algorithm, the passive target positioning algorithm provided by the invention introduces a forgetting factor to construct an evanescent volume Kalman filtering algorithm, and combines the least square positioning algorithm to form a strong tracking filter. The filter does not need to give an initial value of target iteration, and the maneuvering target can be passively positioned only by utilizing a CV tracking model. The method can be applied to multi-platform underwater multi-target association and passive positioning.
Description
Technical Field
The invention belongs to the technical field of distributed multi-platform underwater multi-target positioning and tracking, and particularly relates to a distributed multi-platform underwater multi-target association and passive positioning method.
Background
For the same target, due to the existence of the measurement error of the platform itself and the interference of the environmental noise, the measurement information observed by a plurality of platforms is not completely the same, but certain similar characteristics, such as flight path, frequency characteristics, etc., exist. How to use the similar and not identical characteristics to judge whether the characteristics come from the same target or not is target association. The current underwater multi-target association mainly faces the following difficulties: due to the objective existence of marine environmental noise, the limitation of the observation capability of the platform and the lack of understanding of target prior information, people cannot determine the number of tracked targets and are difficult to perform false alarm judgment (Liu Lu. Research on multi-platform-based multi-target track synthesis and tracking method [ D ]. Harbin university of engineering, 2015.), so that information ambiguity exists between multiple platforms and multiple targets, and the accuracy of target association is difficult to ensure. The scholars have proposed more target association algorithms, such as the classical algorithms of the nearest neighbor method, the track fuzzy association method and the probability data association, and have achieved certain research results (Singer R A and Kanyuck A T. Computer control of multiple site track data [ J ]. Automation,1971,7 (3): 455-463. And the new data fusion track association algorithm [ J ]. The university of Saian electronics, 2012,39 (01): 67-74.). However, in general, most of the current research results are based on track correlation, and the frequency characteristics of the target are less utilized.
In multi-platform passive object positioning, the pure orientation (BOT) problem has received a lot of attention. In related research and literature, a non-linear Bayesian filter algorithm (such as PLKF, EKF, UKF and the like) or a non-Bayesian filter algorithm (such as a least square algorithm) is mostly adopted for positioning and tracking a target with uniform motion (CV). Most of nonlinear Bayesian filtering algorithms need to give iteration initial values in advance, and the result shows that: when the initial iteration value is given to be proper, the filter can be quickly converged, the positioning precision is high, once the initial iteration value is not given properly, the filter is easy to diverge, and the target cannot be tracked. The algorithm is only suitable for target tracking of the CV model, and in an actual offshore environment, under the influence of sea wave blowing, necessary avoidance of maneuvering and self-movement control capacity, the target has higher probability of performing various maneuvering such as uniform speed, acceleration, turning and the like, but not only performing uniform-speed linear movement, if the maneuvering target cannot be effectively tracked by using the algorithm, the positioning accuracy of the algorithm is greatly reduced, and the algorithm has no practical value. Therefore, most of the existing documents adopt an IMM maneuvering model to track a maneuvering target, and although the tracking result is ideal, more prior knowledge is required, such as the movement model (CV, CT, CA) information of the maneuvering target and the transition probability between models, and once the prior knowledge is not well known, the accuracy of the positioning algorithm is reduced.
Disclosure of Invention
The invention aims to solve the problems of low positioning accuracy and low accuracy of multi-target association of the conventional maneuvering target positioning and tracking method, and provides a distributed multi-platform underwater multi-target association and passive positioning method.
The technical scheme adopted by the invention for solving the technical problems is as follows: a distributed multi-platform underwater multi-target association and passive positioning method specifically comprises the following steps:
firstly, performing target association according to an information vector estimation result of each target under multiple platforms to obtain a target association result;
step two, constructing a TMA model consisting of a state equation and a measurement equation;
and step three, passively positioning the pure azimuth maneuvering target under the multiple platforms based on the target correlation result in the step one and the TMA model in the step two.
Further, the specific process of the first step is as follows:
step one, extracting target frequency characteristic information vector estimation results from information vector estimation results of each target under multiple platforms:
J mp (k)＝{A mp (k),F mp (k),k＝1,2,…,T} (1)
where P represents the pth target, P is {1,2, \8230;, P }, P represents the total number of targets, M represents the mth platform, M is {1,2, \8230;, M }, M represents the total number of platforms, k represents the sampling time k, T represents the total sampling time, A represents the total sampling time mp (k) Representing the energy information of the p-th object observed by the m-th platform at the k-th sampling instant, F mp (k) A set of frequency features representing the p-th object observed by the m-th platform at the k-th sampling instant, J mp (k) Representing the sample at the kthEstimating a frequency characteristic information vector estimation result of a pth target observed by an mth platform at a moment;
step two, define event H 0 And event H 1 The following were used:
event H 0 ：J mp (k) And J nq (k) Derived from the same object, J nq (k) Representing the estimation result of the frequency characteristic information vector of the qth target observed by the nth platform at the kth sampling moment, n belongs to {1,2, \8230;, M }, q belongs to {1,2, \8230;, P };
event H 1 ：J mp (k) And J nq (k) Are not derived from the same target;
step one, defining Euclidean distance test statistic as follows:
wherein λ is F (t k′ ) For the frequency test statistic, t k′ For the k' th sliding window, the size of the sliding window is w, t k′ K, k +1, \8230;, k + w, i.e. w sample instants are comprised within the sliding window,
is the mean of the l-th frequency feature in the p-th target frequency feature set observed by the m-th platform at the k-th sampling time to the l-th frequency feature in the p-th target frequency feature set observed by the m-th platform at the k + w sampling time,
is the mean value of the L frequency characteristic from the q target frequency characteristic set observed by the nth platform at the kth sampling moment to the L frequency characteristic set observed by the qth platform at the kth and w sampling moments, wherein L =1,2, \ 8230, L 0 ，L 0 The number of features included in the frequency feature set;
λ A (t k′ ) For the purposes of the energy check statistic,
is the average value of the l 'th energy in the p-th target energy information observed by the mth platform at the kth sampling time to the l' th energy in the p-th target energy information observed by the mth platform at the k + w sampling time,
is the mean value of the L 'energy in the q-th target energy information observed by the nth platform at the k-th sampling time to the L' energy in the q-th target energy information observed by the nth platform at the k + w sampling time, L '=1,2, \ 8230, L' 0 ，L′ 0 Is the number of energies contained in the energy information,
is that
The maximum value of (a) is,
is that
Maximum value of (1);
step four, utilizing lambda F (t k′ ) And λ A (t k′ ) Calculating the correlation test statistical evaluation value in the k' sliding window:
λ(t k′ )＝α·η F (t k′ )+β·η A (t k′ ) (3)
wherein, λ (t) k′ ) For the correlation test statistical evaluation value in the kth' sliding window, eta F (t k′ ) Is λ F (t k′ ) Corresponding correlation confidence value, η A (t k′ ) Is λ A (t k′ ) Corresponding correlation confidence values, wherein alpha is the weight of the frequency test statistic, and beta is the weight of the energy test statistic;
step one or five, according to lambda (t) k′ ) Calculating a target association probability p:
wherein k' =1,2, \8230, N-w;
if p > ε is satisfied, then event H 0 If true, the target association is successful; otherwise, event H 1 True, where ε is 90%.
Further, said λ F (t k′ ) Corresponding correlation confidence values η F (t k′ ) Comprises the following steps:
wherein, a 1 And a 2 Is a constant;
said lambda A (t k′ ) Corresponding correlation confidence values η A (t k′ ) Comprises the following steps:
wherein, b 1 And b 2 Is a constant.
Further, the constant a 1 ＝δ f ，a 2 ＝5δ f ，δ f Is the standard deviation of the error of frequency measurement, b 1 ＝0.1，b 2 =0.5. In the present invention, a 1 ＝5，a 2 ＝10。
Further, the weight α is 0.8, and the weight β is 0.2.
Further, the specific process of the second step is as follows:
wherein X (k | k) represents the target state quantity at time k, Z (k | k) represents the target measurement at time k, F (k) represents the state transition matrix at time k, Γ represents the process noise distribution matrix, H represents the measurement transition matrix, X (k-1 component k-1) represents the target state quantity at time k-1, v is the state process noise vector, and w is the measurement noise vector.
Further, the calculation formula of the target measurement Z (k | k) at the time k is as follows:
Z(k|k)＝[x(k),y(k)] T ＝[(A T A) -1 A T f 2 ,(A T A) -1 A T f 1 ] T (8)
wherein (x (k), y (k)) is the position of the target at time k, the superscript T represents the transposition of the matrix, the superscript-1 represents the inverse of the matrix, A, f 1 And f 2 Is an intermediate variable;
wherein M, n belongs to {1,2, \8230;, M }, M is not equal to n, beta m For the position estimation result corresponding to the platform m, beta n For the position estimate, x, corresponding to platform n om ,y om Position parameter information, x, for platform m on ,y on And the position parameter information corresponding to the platform n.
Further, the specific process of the third step is as follows:
Where P (k | k) is the covariance matrix at time k, chol { P (k | k) } denotes Cholesky decomposition of the matrix P (k | k), S (k) is an intermediate variable,
l is the dimension of the target state quantity;
wherein, the first and the second end of the pipe are connected with each other,
a one-step prediction result of the volume points;
and 3, one-step prediction of the target state quantity:
wherein X (k +1 purple) is a one-step prediction result of the target state quantity;
and 4, introducing covariance prediction of a forgetting factor:
p (k + 1) is covariance of introduced forgetting factors, lambda (k) is forgetting factors at the moment k, and Q (k) is a state noise covariance matrix at the moment k;
wherein, the first and the second end of the pipe are connected with each other,
for the updated volume point, S (k + 1) is an intermediate variable;
step 6, measuring and predicting:
wherein Z is j (k +1 purple k) is a one-step prediction result of the measurement;
wherein V (k + 1) is innovation covariance at the moment k +1, Z (k + 1) is target measurement at the moment k +1, and R (k) is a measurement noise covariance matrix at the moment k;
wherein K (K + 1) is the gain at the moment K + 1;
and 9, updating the target state quantity:
X(k+1|k+1)＝X(k+1|k)+K(k+1)(Z(k+1|k)-Z(k+1|k+1)) (18)
wherein, X (k +1 purple k + 1) represents the updated target state quantity, i.e. the target state quantity at the moment of k + 1;
P(k+1|k+1)＝P(k+1|k)-K(k+1)V(k+1)K(k+1) T (19)
wherein P (k +1 calc + 1) is the updated covariance matrix, i.e. the covariance matrix at the moment k + 1;
step 11, updating a forgetting factor:
wherein R (k + 1) is the measured noise covariance matrix at the time of k +1,
and C (k + 1) is an intermediate variable matrix, trace (·) represents solvingThe trace of the matrix, lambda (k + 1) is the forgetting factor at the moment of k + 1;
and step 12, returning to step 1 by using the target state quantity updated in step 9, the covariance matrix updated in step 10 and the forgetting factor updated in step 11.
The beneficial effects of the invention are:
the invention introduces a target line spectrum characteristic set into a target association method, provides a target association algorithm based on frequency characteristic multi-information vector fusion, solves the problems that a pure orientation system cannot acquire a flight path and the flight path association is difficult, and can give a correct association result. Compared with the traditional nonlinear Bayesian filtering algorithm, the passive target positioning algorithm provided by the invention introduces a forgetting factor to construct an evanescent volume Kalman filtering algorithm, and combines the least square positioning algorithm to form a strong tracking filter. The filter does not need to give a target iteration initial value, can passively position a maneuvering target by only using a CV tracking model, and has positioning accuracy superior to a traditional pure orientation EKF algorithm.
Drawings
FIG. 1 is a diagram of the movement locus of 4 targets in x-y plane coordinates involved in simulation experiments of the present invention;
FIG. 2a is a frequency history diagram of the platform 1 according to the simulation experiment of the present invention;
FIG. 2b is a frequency history diagram of the platform 2 involved in the simulation experiment of the present invention;
FIG. 3a is a diagram of the platform 1 azimuth history involved in the simulation experiment of the present invention
FIG. 3b is a diagram of the azimuth history of the platform 2 involved in the simulation experiment of the present invention;
FIG. 4a is a comparison of the target 1 of platform 1 and all target-associated test statistics of platform 2;
FIG. 4b is a comparison of the target 2 of platform 1 and all target-associated test statistics of platform 2;
FIG. 4c is a comparison of the target 3 of platform 1 versus all target-associated test statistics of platform 2;
FIG. 4d is a comparison of the target 4 of platform 1 and all target-associated test statistics of platform 2;
FIG. 5a is a graph of the positioning error of the target 1 according to the method of the present invention and the conventional method;
FIG. 5b is a graph of the positioning error of the target 2 according to the method of the present invention and the conventional method;
FIG. 5c is a graph of the positioning error of the target 3 according to the method of the present invention and the conventional method;
FIG. 5d is a graph of the positioning error of the target 4 according to the method of the present invention and the conventional method;
FIG. 6 is a graph of the root mean square error of the method of the present invention compared to a conventional method.
Detailed Description
In a first specific embodiment, a distributed multi-platform underwater multi-target association and passive positioning method in this embodiment specifically includes the following steps:
firstly, performing target association according to an information vector estimation result of each target under multiple platforms to obtain a target association result;
step two, constructing a TMA model consisting of a state equation and a measurement equation;
and step three, passively positioning the purely-azimuth maneuvering target under the multiple platforms by utilizing an evanescent tracking filter based on a forgetting factor for the target associated in the step one and the TMA model in the step two, and finally outputting the positioning result (including information such as position, speed and the like) of each target.
In the invention, the multiple platforms mean more than or equal to 2 platforms, and the multiple targets mean more than or equal to 2 targets.
The second embodiment is as follows: the difference between the first embodiment and the first embodiment is that a target association algorithm based on frequency feature multi-information vector fusion is adopted in the first step, and the specific process of the first step is as follows:
step one, inputting a target information vector estimation result under multiple platforms:
I mp (k)＝{θ mp (k),A mp (k),F mp (k),k＝1,2,…,T},p＝1,2,…,P,m＝1,2,…,M
extracting target frequency characteristic information vector estimation results from information vector estimation results of each target under multiple platforms:
J mp (k)＝{A mp (k),F mp (k),k＝1,2,…,T} (1)
where P represents the pth target, P is {1,2, \8230;, P }, P represents the total number of targets, M represents the mth platform, M is {1,2, \8230;, M }, M represents the total number of platforms, k represents the sampling time k, T represents the total sampling time, A represents the total sampling time mp (k) Representing the energy information of the p-th object observed by the m-th platform at the k-th sampling instant, F mp (k) A set of frequency features representing the p-th object observed by the m-th platform at the k-th sampling instant, J mp (k) Representing the estimation result of the frequency characteristic information vector of the p target observed by the m platform at the k sampling moment;
step one, two, define event H 0 And event H 1 The following were used:
event H 0 ：J mp (k) And J nq (k) Derived from the same object, J nq (k) Representing the estimation result of the frequency characteristic information vector of the qth target observed by the nth platform at the kth sampling moment, wherein n belongs to {1,2, \8230;, M }, and q belongs to {1,2, \8230;, P };
event H 1 ：J mp (k) And J nq (k) Not from the same target;
step three, defining Euclidean distance test statistic as follows:
wherein λ is F (t k′ ) For the frequency test statistic, t k′ For the kth' sliding window, the size of the sliding window is w, t k′ K, k +1, \ 8230;, k + w }, k =1,2, \ 8230;, N-w, i.e. w sampling instants are included in the sliding window,
is observed from the ith frequency feature in the p target frequency feature set observed by the mth platform at the kth sampling moment to the mth platform at the k + w sampling momentsThe mean of the l-th frequency feature in the p-th target frequency feature set,
is the mean value of the L frequency characteristic from the q target frequency characteristic set observed by the nth platform at the kth sampling moment to the L frequency characteristic set observed by the qth platform at the kth and w sampling moments, wherein L =1,2, \ 8230, L 0 ，L 0 The number of features included in the frequency feature set;
λ A (t k′ ) For the purposes of the energy-test statistic,
is the average of the l 'th energy in the p-th target energy information observed by the mth platform at the kth sampling time to the l' th energy in the p-th target energy information observed by the mth platform at the k + w sampling time,
is the mean value of the L 'energy in the q-th target energy information observed by the nth platform at the k-th sampling time to the L' energy in the q-th target energy information observed by the nth platform at the k + w sampling time, L '=1,2, \ 8230, L' 0 ，L′ 0 Is the number of energies contained in the energy information,
is that
The maximum value of (a) is,
is that
Maximum value of (1);
step four, utilizing lambda F (t k′ ) And λ A (t k′ ) ComputingCorrelation test statistical evaluation value in kth' sliding window:
λ(t k′ )＝α·η F (t k′ )+β·η A (t k′ ) (3)
wherein, λ (t) k′ ) For the correlation test statistical evaluation value in the kth' sliding window, eta F (t k′ ) Is λ F (t k′ ) Corresponding correlation confidence value, η A (t k′ ) Is λ A (t k′ ) Corresponding correlation confidence values, wherein alpha is the weight of the frequency test statistic, and beta is the weight of the energy test statistic;
step one and five, according to lambda (t) k′ ) Calculating a target association probability p:
wherein k' =1,2, \ 8230, N-w;
if p > epsilon is satisfied, event H 0 If true, the target association is successful; otherwise, event H 1 True, where ε is 90%.
The target association is successful, that is, the estimation result of the frequency characteristic information vector of the qth target observed by the nth platform and the estimation result of the frequency characteristic information vector of the pth target observed by the mth platform are from the same target.
Other steps and parameters are the same as those in the first embodiment.
The third concrete implementation mode: in this embodiment, the difference from the first or second embodiment is that λ F (t k′ ) Corresponding correlation confidence values η F (t k′ ) Comprises the following steps:
wherein, a 1 And a 2 Is a constant;
said lambda A (t k′ ) Corresponding correlation confidence values η A (t k′ ) Comprises the following steps:
wherein, b 1 And b 2 Is a constant.
In the present embodiment, the established test statistic rating criteria are shown in tables 1 and 2:
TABLE 1 frequency test statistic rating Standard
TABLE 2 energy test statistic rating criteria
Other steps and parameters are the same as those in the first or second embodiment.
The fourth concrete implementation mode: this embodiment is different from one of the first to third embodiments in that the constant a 1 ＝δ f ，a 2 ＝5δ f ，δ f Is the standard deviation of the error of frequency measurement, b 1 ＝0.1，b 2 ＝0.5。
Other steps and parameters are the same as those in one of the first to third embodiments.
The fifth concrete implementation mode: the difference between this embodiment and one of the first to the fourth embodiments is that the weight α is 0.8, and the weight β is 0.2.
Other steps and parameters are the same as in one of the first to fourth embodiments.
The sixth specific implementation mode: the difference between this embodiment and one of the first to fifth embodiments is that the specific process of the second step is:
wherein the content of the first and second substances,
x (k | k) represents the target state quantity at time k, Z (k | k) represents the target measurement at time k, F (k) represents the state transition matrix at time k, Γ represents the process noise distribution matrix, H represents the measurement transition matrix, X (k-1) represents the target state quantity at time k-1, v is the state process noise vector, and w is the measurement noise vector.
(x (k), y (k)) represents the position of the target at time k, v x (k) Representing the velocity of the target in the x-direction at time k, v y (k) Represents the speed of the target at the time k in the y direction, wherein the x direction refers to the x-axis direction of the plane rectangular coordinate system, and the y direction refers to the y-axis direction of the plane rectangular coordinate system.
Other steps and parameters are the same as those in one of the first to fifth embodiments.
The seventh embodiment: the difference between this embodiment and one of the first to sixth embodiments is that the calculation formula of the target measurement Z (kk) at the time k is as follows:
Z(k|k)＝[x(k),y(k)] T ＝[(A T A) -1 A T f 2 ,(A T A) -1 A T f 1 ] T (8)
wherein (x (k), y (k)) is the position of the target at time k, the superscript T represents the transposition of the matrix, the superscript-1 represents the inverse of the matrix, A, f 1 And f 2 Is an intermediate variable;
wherein M, n belongs to {1,2, \8230;, M }, M is not equal to n, beta m For the position estimate corresponding to the platform m, β n Is a platformn corresponding to the position estimation result, x om ,y om Position parameter information, x, for platform m on ,y on And the position parameter information corresponding to the platform n. Matrix A, f 1 ,f 2 Are all made of
And the matrix of the row 1 and the column, for the case of three platforms, the matrix A is a matrix of 3 rows and 1 column, the 1 st row is calculated according to the position estimation result of the 1 st platform and the position estimation result of the 2 nd platform, the 2 nd row is calculated according to the position estimation result of the 1 st platform and the position estimation result of the 3 rd platform, and the 3 rd row is calculated according to the position estimation result of the 2 nd platform and the position estimation result of the 3 rd platform.
For any target under multiple platforms, the position estimation result of each platform on the target is beta respectively 1 ，β 2 ，…，β M Using beta 1 ，β 2 ，…，β M After the process of the present embodiment is executed, the position of the target can be obtained. Similarly, after the process of the embodiment is performed on other targets, the position of each target can be obtained.
Other steps and parameters are the same as those in one of the first to sixth embodiments.
The specific implementation mode is eight: the difference between this embodiment and one of the first to seventh embodiments is that the specific process of step three is:
Where P (k | k) is the covariance matrix at time k, chol { P (k | k) } denotes Cholesky decomposition of the matrix P (k | k), S (k) is an intermediate variable,
l is the dimension of the target state quantity;
wherein, the first and the second end of the pipe are connected with each other,
a one-step prediction result of the volume points;
and 3, one-step prediction of the target state quantity:
wherein X (k + 1) is a one-step prediction result of the target state quantity;
and 4, introducing covariance prediction of a forgetting factor:
p (k + 1) is covariance of introduced forgetting factors, lambda (k) is forgetting factors at the moment k, and Q (k) is a state noise covariance matrix at the moment k;
wherein the content of the first and second substances,
for the updated volume point, S (k + 1) is an intermediate variable;
step 6, measuring and predicting:
wherein Z is j (k +1 purple k) is a one-step prediction result of the measurement;
wherein V (k + 1) is innovation covariance at the moment k +1, Z (k + 1) is target measurement at the moment k +1, and R (k) is measurement noise covariance matrix at the moment k;
wherein K (K + 1) is the gain at the moment K + 1;
and 9, updating the target state quantity:
X(k+1|k+1)＝X(k+1|k)+K(k+1)(Z(k+1|k)-Z(k+1|k+1)) (18)
wherein, X (k +1 purple k + 1) represents the updated target state quantity, i.e. the target state quantity at the moment of k + 1;
P(k+1|k+1)＝P(k+1|k)-K(k+1)V(k+1)K(k+1) T (19)
wherein P (k +1 calc + 1) is the updated covariance matrix, i.e. the covariance matrix at the moment k + 1;
step 11, updating a forgetting factor:
wherein R (k + 1) is the measured noise covariance matrix at the time of k +1,
and C (k + 1) is an intermediate variable matrix, trace (·) represents the trace of the matrix, and λ (k + 1) is the forgetting at the moment of k +1A factor;
and step 12, returning to the step 1 by using the target state quantity updated in the step 9, the covariance matrix updated in the step 10 and the forgetting factor updated in the step 11.
And (4) after returning to the step 1, repeating the processes from the step 1 to the step 11, namely finishing the passive positioning of the pure azimuth maneuvering target and outputting a positioning result (comprising position and speed information) of the target. And the initial conditions set in this embodiment are:
the iteration initial value X (1 _ 1) of the filter is:
wherein, T is a sampling interval, (x (1), y (1)) is the position of the target at the 1 st moment, and (x (2), y (2)) is the position of the target at the 2 nd moment;
the forgetting factor initial value is generally set to λ (1) =1.
After the position of each target is obtained through the second step, the positioning result (including the position and speed information) of each target can be output after the process of the embodiment is executed on each target.
Other steps and parameters are the same as those in one of the first to seventh embodiments.
Simulation experiment
Simulation conditions are as follows: assuming a total of 4 motorized targets, the initial position, initial velocity, center frequency of radiation and sound source level of each target are given in table 3 below. A total of 150 frames of data are tracked, spaced 1s apart per frame. All 4 targets make uniform linear motion in the first 20 frames, make uniform turning motion in 21-100 frames (turning rate is given in table 3), and make uniform linear motion in the last 50 frames. The coordinates of the two stationary observatory are (-100, 0) m and (100, 0) m, respectively. The ambient noise is 75dB. The orientation measurement accuracy is given by the CRLB calculation. The frequency resolution is 0.01Hz, and the underwater sound velocity is 1500m/s. The size of the sliding window w in the target correlation algorithm is set to 10.
TABLE 3 initial position, initial velocity, radiation center frequency and sound source level of each target
In the present invention, the x-y plane motion trajectories of the 4 targets are shown in FIG. 1.
The invention simulates 4 maneuvering targets and gives simulation results. Fig. 2a and 2b show frequency information observed by multiple platforms, fig. 3a and 3b show orientation information observed by multiple platforms, fig. 4a shows target 1 of platform 1 compared with all target-related test statistics of platform 2, fig. 4b shows target 2 of platform 1 compared with all target-related test statistics of platform 2, fig. 4c shows target 3 of platform 1 compared with all target-related test statistics of platform 2, and fig. 4d shows target 4 of platform 1 compared with all target-related test statistics of platform 2. Table 4 gives the results of the target association;
TABLE 4
The target association probabilities of the targets 1,2, 3 and 4 under the two platforms are respectively 99.7%, 100% and 99.9%, and the association probabilities of different targets under different platforms do not exceed 20%, so that the correctness of the association algorithm provided by the invention is verified.
From the comparison results of the positioning error curves of 4 targets in fig. 5a to 5d, it can be seen that the positioning error curves of the algorithm of the present invention are respectively about 200m,400m,600m, and 250m, while the positioning error curves of the conventional algorithm are respectively about 200m,24000m,6000m, and 14000m, and the error curves are unable to converge and tend to diverge.
The filter performance of the algorithm of the present invention and the conventional EKF algorithm is specifically analyzed using the above object 1. The measurement angle errors of the two observation stations are set to be 2 degrees, other simulation conditions are unchanged, and Monte Carlo simulation is carried out for 1000 times. The traditional EKF algorithm respectively gives an accurate initial value (target real initial state) and an estimated initial value obtained by calculation of a formula (21), and the algorithm only gives the estimated initial value obtained by calculation of the formula (21).
The performance of the filter is measured using the root mean square error:
wherein MC is Monte Carlo simulation frequency, X (k) is target real state,
for the target estimation state, N is the total number of sampling frames.
As can be seen from the comparison result of FIG. 6, under the same simulation condition, the error of the algorithm of the present invention converges to about 50m after 100 Monte Carlo times, while the traditional algorithm completely diverges and cannot converge. The error curve can only converge to about 200m when the conventional algorithm is given an accurate initial value.
Experiments are carried out for multiple times according to the simulation conditions, monte Carlo simulation is carried out for 1000 times each time, the divergence times of the algorithm are given, and the comparison of the divergence times is shown in Table 5:
TABLE 5 Algorithm divergence times comparison
As can be seen from Table 5, the algorithm of the invention has better control on the divergence of the filtering of the maneuvering target, the situation of divergence does not occur in multiple experiments, the number of divergence times of the inverse traditional algorithm is very large under the condition of not using an accurate initial value, and the practical value is lower. The feasibility of the algorithm provided by the invention is verified through the simulation experiment results.
The above-described calculation examples of the present invention are merely to explain the calculation model and the calculation flow of the present invention in detail, and are not intended to limit the embodiments of the present invention. It will be apparent to those skilled in the art that other variations and modifications of the present invention can be made based on the above description, and it is not intended to be exhaustive or to limit the invention to the precise form disclosed, and all such modifications and variations are possible and contemplated as falling within the scope of the invention.
Claims (7)
1. A distributed multi-platform underwater multi-target association and passive positioning method is characterized by specifically comprising the following steps:
firstly, performing target association according to an information vector estimation result of each target under multiple platforms to obtain a target association result;
the specific process of the step one is as follows:
step one, extracting target frequency characteristic information vector estimation results from information vector estimation results of each target under multiple platforms:
J mp (k)＝{A mp (k),F mp (k),k＝1,2,…,T} (1)
wherein, P represents the P-th target, P is the {1,2, \8230;, P }, P represents the total number of targets, M represents the M-th platform, M is the {1,2, \8230;, M }, M represents the total number of platforms, k represents the sampling time k, T represents the total sampling time, A represents the total sampling time mp (k) Representing the energy information of the p-th object observed by the m-th platform at the k-th sampling instant, F mp (k) A set of frequency features representing the p-th object observed by the m-th platform at the k-th sampling instant, J mp (k) Representing the frequency characteristic information vector estimation result of the p-th target observed by the m-th platform at the k-th sampling moment;
step two, define event H 0 And event H 1 The following:
event H 0 ：J mp (k) And J nq (k) Derived from the same object, J nq (k) Representing the estimation result of the frequency characteristic information vector of the qth target observed by the nth platform at the kth sampling moment, n belongs to {1,2, \8230;, M }, q belongs to {1,2, \8230;, P };
event H 1 ：J mp (k) And J nq (k) Not from the same target;
step three, defining Euclidean distance test statistic as follows:
wherein λ is F (t k′ ) For frequency test statistics, t k′ For the kth' sliding window, the size of the sliding window is w +1,t k′ K, k +1, \ 8230;, k + w, i.e., w +1 sample instants are included within the sliding window,
is the mean of the l-th frequency feature in the p-th target frequency feature set observed by the m-th platform at the k-th sampling time to the l-th frequency feature in the p-th target frequency feature set observed by the m-th platform at the k + w sampling time,
is the mean value of the L frequency characteristic from the q target frequency characteristic set observed by the nth platform at the kth sampling moment to the L frequency characteristic set observed by the qth platform at the kth and w sampling moments, wherein L =1,2, \ 8230, L 0 ，L 0 The number of features included in the frequency feature set;
λ A (t k′ ) For the purposes of the energy-test statistic,
is the average value of the l 'th energy in the p-th target energy information observed by the mth platform at the kth sampling time to the l' th energy in the p-th target energy information observed by the mth platform at the k + w sampling time,
is the l' th energy to the k + w th sampling time in the q target energy information observed by the nth platform at the kth sampling timeThe mean value of L 'energy, L' =1,2, \8230;, L 'in the qth target energy information observed by the nth platform' 0 ，L′ 0 Is the number of energies contained in the energy information,
is that
The maximum value of (a) is,
is that
Maximum value of (1);
step four, utilizing lambda F (t k′ ) And λ A (t k′ ) Calculating the associated test statistical evaluation value in the k' th sliding window:
λ(t k′ )＝α·η F (t k′ )+β·η A (t k′ ) (3)
wherein, λ (t) k′ ) For the correlation test statistical evaluation value in the kth' sliding window, eta F (t k′ ) Is λ F (t k′ ) Corresponding correlation confidence value, η A (t k′ ) Is λ A (t k′ ) Corresponding correlation confidence values, wherein alpha is the weight of the frequency test statistic, and beta is the weight of the energy test statistic;
step one and five, according to lambda (t) k′ ) Calculating a target association probability p:
wherein k' =1,2, \ 8230, T-w;
if p > ε is satisfied, then event H 0 If true, the target association is successful; otherwise, event H 1 Is true, itIn, epsilon takes a value of 90%;
step two, constructing a TMA model consisting of a state equation and a measurement equation;
and step three, passively positioning the pure azimuth maneuvering target under the multiple platforms based on the target correlation result in the step one and the TMA model in the step two.
2. The distributed multi-platform underwater multi-target association and passive positioning method as claimed in claim 1, wherein said λ is F (t k′ ) Corresponding correlation confidence values η F (t k′ ) Comprises the following steps:
wherein, a 1 And a 2 Is a constant;
said lambda A (t k′ ) Corresponding correlation confidence values η A (t k′ ) Comprises the following steps:
wherein, b 1 And b 2 Is a constant.
3. The distributed multi-platform underwater multi-target association and passive positioning method as claimed in claim 2, wherein the constant a is 1 ＝δ f ，a 2 ＝5δ f ，δ f Is the standard deviation of the error of frequency measurement, b 1 ＝0.1，b 2 ＝0.5。
4. The distributed multi-platform underwater multi-target association and passive positioning method according to claim 3, wherein the value of the weight α is 0.8, and the value of the weight β is 0.2.
5. The distributed multi-platform underwater multi-target association and passive positioning method according to claim 4, wherein the specific process of the second step is as follows:
wherein X (k | k) represents the target state quantity at time k, Z (k | k) represents the target measurement at time k, F (k) represents the state transition matrix at time k, Γ represents the process noise distribution matrix, H represents the measurement transition matrix, X (k-1 component k-1) represents the target state quantity at time k-1, v is the state process noise vector, and w' is the measurement noise vector.
6. The distributed multi-platform underwater multi-target association and passive positioning method as claimed in claim 5, wherein the target measurement Z (k | k) at the time k is calculated by the following formula:
Z(k|k)＝[x(k),y(k)] T ＝[(A T A) -1 A T f 2 ,(A T A) -1 A T f 1 ] T (8)
wherein (x (k), y (k)) is the position of the target at time k, the superscript T represents the transposition of the matrix, the superscript-1 represents the inverse of the matrix, A, f 1 And f 2 Is an intermediate variable;
wherein M, n belongs to {1,2, \8230;, M }, M ≠ n, beta m For the position estimate corresponding to the platform m, β n For the position estimate, x, corresponding to platform n om ,y om Position parameter information, x, for platform m on ,y on And the position parameter information corresponding to the platform n.
7. The distributed multi-platform underwater multi-target association and passive positioning method according to claim 6, wherein the specific process of the third step is as follows:
step 1, calculating volume points
Where P (k | k) is the covariance matrix at time k, chol { P (k | k) } means Cholesky decomposition of the matrix P (k | k), S (k) is an intermediate variable,
l is the dimension of the target state quantity;
step 2, one-step prediction of volume points:
wherein the content of the first and second substances,
a one-step prediction result of the volume points;
step 3, one-step prediction of the target state quantity:
wherein X (k + 1) is a one-step prediction result of the target state quantity;
and 4, introducing covariance prediction of a forgetting factor:
wherein P (k +1 k) is a covariance of introducing a forgetting factor, lambda (k) is a forgetting factor at the time of k, and Q (k) is a state noise covariance matrix at the time of k;
step 5, updating the volume points:
wherein, the first and the second end of the pipe are connected with each other,
for the updated volume point, S (k + 1) is an intermediate variable;
step 6, measuring and predicting:
wherein, Z j (k +1 red k) is a one-step prediction result of measurement;
step 7, innovation covariance estimation:
wherein V (k + 1) is innovation covariance at the moment k +1, Z (k + 1) is target measurement at the moment k +1, and R (k) is a measurement noise covariance matrix at the moment k;
step 8, gain calculation:
wherein K (K + 1) is the gain at the moment K + 1;
and 9, updating the target state quantity:
X(k+1|k+1)＝X(k+1|k)+K(k+1)(Z j (k+1|k)-Z(k+1|k+1)) (18)
wherein X (k +1 calc + 1) represents the updated target state quantity, i.e. the target state quantity at the time of k + 1;
step 10, covariance updating:
P(k+1|k+1)＝P(k+1|k)-K(k+1)V(k+1)K(k+1) T (19)
wherein P (k +1 calc + 1) is the updated covariance matrix, i.e. the covariance matrix at the moment k + 1;
step 11, updating a forgetting factor:
wherein R (k + 1) is the measured noise covariance matrix at the time of k +1,
and C (k + 1) is an intermediate variable matrix, trace (·) represents a trace of the matrix, and λ (k + 1) is a forgetting factor at the moment of k + 1;
and step 12, returning to the step 1 by using the target state quantity updated in the step 9, the covariance matrix updated in the step 10 and the forgetting factor updated in the step 11.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202210198348.0A CN114577213B (en) | 2022-03-01 | 2022-03-01 | Distributed multi-platform underwater multi-target association and passive positioning method |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202210198348.0A CN114577213B (en) | 2022-03-01 | 2022-03-01 | Distributed multi-platform underwater multi-target association and passive positioning method |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN114577213A CN114577213A (en) | 2022-06-03 |
| CN114577213B true CN114577213B (en) | 2022-11-18 |
Family
ID=81771779
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202210198348.0A Active CN114577213B (en) | 2022-03-01 | 2022-03-01 | Distributed multi-platform underwater multi-target association and passive positioning method |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN114577213B (en) |
Family Cites Families (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN104198992B (en) * | 2014-09-11 | 2016-10-05 | 东南大学 | Acoustic Object Passive Location based on multidiameter delay structure compresses perception |
| CN110109095B (en) * | 2019-04-30 | 2022-10-28 | 西南电子技术研究所（中国电子科技集团公司第十研究所） | Target feature assisted multi-source data association method |
| KR102306061B1 (en) * | 2020-03-06 | 2021-09-28 | 국방과학연구소 | Apparatus of multiple targets management for multistatic pcl based target localization |
| CN111505580B (en) * | 2020-04-14 | 2022-04-15 | 哈尔滨工程大学 | Multi-platform cooperative target positioning method based on azimuth angle and Doppler information |
| CN113359100A (en) * | 2021-06-09 | 2021-09-07 | 电子科技大学 | Attribute data association method for tracking target radar radiation source |
-
2022
- 2022-03-01 CN CN202210198348.0A patent/CN114577213B/en active Active
Also Published As
| Publication number | Publication date |
|---|---|
| CN114577213A (en) | 2022-06-03 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN111985093B (en) | Adaptive unscented Kalman filtering state estimation method with noise estimator | |
| CN110503071B (en) | Multi-target tracking method based on variational Bayesian label multi-Bernoulli superposition model | |
| CN109633590B (en) | Extended target tracking method based on GP-VSMM-JPDA | |
| CN112613532B (en) | Moving target tracking method based on radar and cyclic neural network complement infrared fusion | |
| CN107193009A (en) | A kind of many UUV cooperative systems underwater target tracking algorithms of many interaction models of fuzzy self-adaption | |
| CN104199022B (en) | Target modal estimation based near-space hypersonic velocity target tracking method | |
| Ma et al. | Target tracking system for multi-sensor data fusion | |
| CN111913175A (en) | Water surface target tracking method with compensation mechanism under transient failure of sensor | |
| CN111711432B (en) | Target tracking algorithm based on UKF and PF hybrid filtering | |
| CN108152812B (en) | Improved AGIMM tracking method for adjusting grid spacing | |
| CN111259332B (en) | Fuzzy data association method and multi-target tracking method in clutter environment | |
| CN110297221B (en) | Data association method based on Gaussian mixture model | |
| CN109919233B (en) | Tracking filtering method based on data fusion | |
| CN107621632A (en) | Adaptive filter method and system for NSHV tracking filters | |
| CN111948601B (en) | Single-station pure-angle target positioning and tracking method under non-Gaussian noise condition | |
| CN109856623A (en) | A kind of Target state estimator method for more radar rectilinear path lines | |
| CN114577213B (en) | Distributed multi-platform underwater multi-target association and passive positioning method | |
| CN110426689B (en) | Airborne multi-platform multi-sensor system error registration algorithm based on EM-CKS | |
| CN116520311A (en) | GLMB-based adaptive track initiation method | |
| CN113280821B (en) | Underwater multi-target tracking method based on slope constraint and backtracking search | |
| CN115544425A (en) | Robust multi-target tracking method based on target signal-to-noise ratio characteristic estimation | |
| CN113866754A (en) | Moving target track correlation method based on Gaussian distribution wave gate | |
| CN113376626A (en) | High maneuvering target tracking method based on IMMPDA algorithm | |
| Liu et al. | Adaptive MCMC particle filter for tracking maneuvering target | |
| Sönmez et al. | Analysis of performance criteria for optimization based bearing only target tracking algorithms |
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 |