Underwater weak target azimuth track-before-detect method based on belief propagation

Abstract

The invention relates to an underwater weak target azimuth track-before-detect method based on belief propagation. The method comprises the following steps: converting a weak target tracking problem into solving of a joint posterior probability density function (pdf) based on a Bayesian rule; factorization is carried out on the joint posteriori pdf, and a corresponding factor graph is constructed through a graph model; solving information transmitted in the factor graph by using a belief propagation (BP) algorithm, and converting multiple integrals into common integrals by using a Goodman variable substitution principle; constructing a likelihood function through complex Wishart distribution to describe the relationship between the original sonar data and the target state, and simplifying the residual pdf by using Gaussian distribution to obtain confidence approximation of the target state edge posteriori pdf; and calculating by using a minimum mean square error (MMSE) estimator to obtain the number of weak targets and corresponding states. According to the method, multiple weak targets can be tracked at the same time without a complex clustering algorithm, the number of preset parameters is small, and target state information including the orientation, the azimuth angle speed and the target intensity can be tracked accurately and efficiently.

Description

Technical Field

[0001] The invention relates to the field of weak target detection and tracking, in particular to an underwater weak target tracking method before azimuth detection based on belief propagation. Background Art

[0002] Sonar arrays are widely used in underwater target detection and tracking. With the advancement of submarine stealth technology and the development of unmanned equipment, the sound source level of these targets has been continuously reduced, and they are typical low signal-to-noise ratio targets (weak targets). The traditional detect-before-tracking algorithm is limited by the threshold setting. When the threshold is high, weak targets cannot be detected. When the threshold is low, a large amount of clutter will be introduced, which brings great challenges to the accurate and efficient tracking of weak targets. The track-before-detection (TBD) algorithm takes the original sonar array data as input, and can extract the implicit weak target information from it, providing a new solution for weak target detection and tracking.

[0003] At present, the single target pre-detection tracking methods in the sonar field are relatively mature, mainly including hidden Markov model Baum-Welsh filter or smoother, Viterbi algorithm and particle filter. These methods cannot track multiple targets and describe the number of dynamically changing targets, and their application scenarios are limited. The TBD algorithm based on random finite sets can effectively solve the problem of tracking multiple weak targets. The main methods are probability hypothesis density (PHD) and cardinality probability hypothesis density (CPHD) filters. These methods have a relatively complex derivation process, and the weak target tracking accuracy depends largely on the clustering algorithm and clustering strategy. The tracking accuracy of the method is low, the stability is poor, and the computational efficiency is also low. Summary of the invention

[0004] In view of the above technical deficiencies, the purpose of the present invention is to provide a method for tracking underwater weak targets before azimuth detection based on belief propagation. It can quickly and accurately derive the marginal posterior PDF of the target state from the joint posterior PDF through the belief propagation (BP) algorithm, and each target is maintained by a separate marginal posterior PDF, without the need for complex clustering algorithms and clustering strategies, which greatly improves the estimation accuracy of the weak target state. In addition, the Goodman variable substitution principle is used to simplify the calculation of multiple integrals in the derivation process, and the covariance-based likelihood function is constructed using the complex Wishart distribution. The sampling data of multiple pins is accumulated in time to improve the tracking accuracy and tracking stability of weak targets, and the remaining distribution is approximated as a Gaussian distribution to further simplify the derivation process and improve the computational efficiency of the method.

[0005] The technical solution adopted by the present invention to achieve the above-mentioned purpose is:

[0006] A method for tracking underwater weak targets before azimuth detection based on belief propagation comprises the following steps:

[0007] 1) Based on the Bayesian rule, the joint posterior probability density function pdf is used to describe the weak target tracking problem;

[0008] 2) Factorize the joint posterior pdf and construct the corresponding factor graph through the graphical model;

[0009] 3) Use the belief propagation BP algorithm to sequentially calculate the prediction information, measurement estimation information and measurement update information in the factor graph to achieve weak target tracking before azimuth detection.

[0010] The step 1) is specifically as follows:

[0011] Combining the sonar measurement information at all times with the weak target state information to be solved, we get the joint posterior PDF expression:

[0012] f(y 0:k |z 0:k )

[0013] Among them, y 0:k Represents all target states from time 0 to time k, target state N k Represents the number of targets at time k, each augmented state The binary variable r k,n ∈{0,1} is used to indicate whether the target n exists, r k,n =1 means existence, is the target state information, respectively representing the orientation θ k,n , azimuth velocity With target strength z 0:k represents all sonar measurement information from time 0 to time k, z k represents the measurement at time k.

[0014] In step 2), the factorized joint posterior pdf expression is:

[0015]

[0016] Among them, the target state y at time k k It consists of two parts, one of which is the existence of the target y at time k-1 k-1 Through the state transfer equation f(y k,n |y k-1,m ) is calculated, and the target number is N k-1 The other part is the new target at time k, which obeys the pre-set probability distribution and the number of targets is N. k -N k-1 , f(z k |y k) is the sonar likelihood function used to describe the relationship between sonar measurements and target states.

[0017] The step 3) comprises the following steps:

[0018] 3.1) Based on confidence Calculate prediction information α k,n ;

[0019] 3.2) Based on prediction information α k,n , using BP algorithm to calculate the measurement estimation information β n ;

[0020] 3.3) Estimation of information β based on measurement n , input the current measurement z k , calculate the measurement update information κ n ;

[0021] 3.4) Based on prediction information α k,n , measurement update information κ n , calculate the approximate confidence of the posterior pdf of the target state margin

[0022] 3.5) Based on confidence Calculate the number of targets and their corresponding states through the MMSE estimator;

[0023] 3.6) Complete the current moment calculation process and Return to step 3.1) to start the next moment calculation.

[0024] The step 3.1) is specifically as follows:

[0025] For each target n∈{1,…,N k-1}, prediction information α at time k k,n for:

[0026]

[0027] α k,n (x k,n ,r k,n =0) = 1-α k,n (x k,n ,r k,n =1)

[0028] Among them, p s represents the target survival probability, represents the marginal posterior pdff(x) of the target state at the previous moment k-1,n ,r k-1,n =1|z 0:k-1 )’s approximate confidence;

[0029] For each new target n′∈{N k-1 +1,…,N k}, new target prediction information α at time k k,n′ for:

[0030] α k,n′ (x k,n′ ,r k,n′ =1)=p b f(x k,n′ )

[0031] α k,n′ (x k,n′ ,r k,n′ =0) = 1-α k,n′ (x k,n′ ,r k,n′ =1)

[0032] Among them, p b represents the target newborn probability, f(x k,n′ ) is the pre-set new target state pdf.

[0033] The step 3.2) is specifically as follows:

[0034] β n (x k,n ,r k,n )=α k,n (x k,n ,r k,n ).

[0035] The step 3.3) is specifically as follows:

[0036] κ n (x k,n ,r k,n ;z k )=∫f(z k -r k,n g(x k,n )|ξ k,n )f(ξ k,n )dξ k,n

[0037] Among them, f(ξ k,n ) is about custom variables The probability density function, y k\n Indicates that time k does not contain the target state y k,n The remaining target state set, for n∈{1,…,N k},g(y k,n ) represents the transition from the target state x k,n to sonar measurements, Among them, a(θk,n ) represents the array manifold vector, and H represents the conjugate transpose of the matrix.

[0038] The step 3.4) is specifically as follows:

[0039] Use complex Wishart distribution to describe the likelihood function f(z k -r k,n g(x k,n )|ξ k,n ), measure z at time k k Converted to complex sampling covariance matrix

[0040]

[0041] Where v represents the number of samples, z k,i It means that the array data is obtained by sampling the sonar for the i-th time at time k, and the likelihood function expression is:

[0042]

[0043] in, represents a p-dimensional complex Wishart distribution with a matrix variable Y, r degrees of freedom, and a scale matrix A, R σ Represents the complex covariance matrix satisfied by the environmental background noise;

[0044] Use a Gaussian distribution to approximate pdff(ξ k,n ) to simplify the calculation, we get:

[0045]

[0046] in represents a Gaussian distribution with a variable Y, a mean of μ, and a covariance of C, where the mean and covariance Specifically:

[0047]

[0048] Rewrite information n The expression is:

[0049]

[0050]

[0051] Based on the BP algorithm, we get the target state edge posterior pdff(y k,n |z 0:k ) is:

[0052]

[0053] in, is the normalization factor.

[0054] The step 3.5) is specifically as follows:

[0055] For n∈{1,…,N k}, calculate the probability of target existence

[0056]

[0057] like Greater than the preset target existence threshold The target n is considered to exist at time k, and the target state is calculated by the MMSE estimator

[0058]

[0059] Statistics in n∈{1,…,N k}The probability of the target existing in The number of weak targets in the area is estimated.

[0060] The present invention has the following beneficial effects and advantages:

[0061] 1. The method of the present invention uses a belief propagation algorithm to derive a confidence approximation of the posterior PDF of the target state edge, has a concise and clear derivation process, and improves the estimation accuracy of the underwater weak target state.

[0062] 2. The method of the present invention adopts the Goodman variable substitution principle to simplify multiple integrals into common integrals, simplifies the calculation process of measurement update information in the BP algorithm, and improves the calculation efficiency.

[0063] 3. The method of the present invention adopts complex Wishart distribution to describe the relationship between the sonar array measurement and the target state, and calculates the sampling covariance matrix by accumulating the sampling values ​​over a period of time to make the tracking result accurate and stable.

[0064] 4. The method of the present invention uses Gaussian distribution to approximate the measurement update to remove the residual distribution of the likelihood function calculation, and extracts the mean of the Gaussian distribution as the scale matrix when calculating the complex Wishart distribution, which further improves the calculation efficiency.

[0065] 5. The method of the present invention is a universal target tracking algorithm. In the face of different application scenarios, it only needs to analyze the relationship between sensor measurement and target state and construct a suitable likelihood function. It has strong applicability. BRIEF DESCRIPTION OF THE DRAWINGS

[0066] Figure 1It is a structural diagram of the method of the present invention;

[0067] Figure 2 is a factor graph constructed based on the graphical model by the method of the present invention;

[0068] Figure 3 It is a flow chart of information calculation of the belief propagation algorithm of the method of the present invention. DETAILED DESCRIPTION

[0069] The present invention is further described in detail below in conjunction with the accompanying drawings and embodiments.

[0070] The basic idea of ​​the present invention is to convert the underwater weak target tracking problem into the solution of the joint posterior pdf based on the Bayesian rule, and to find the marginal posterior pdf of the target state under the condition that the sonar array measurement is known. Because it is difficult to give a display expression of the joint posterior pdf, it is impractical to solve the marginal posterior pdf directly by integration. The present invention factors the joint posterior pdf, constructs a factor graph of the weak target tracking problem using a graph model, and uses the BP algorithm to solve the information transmitted between each node in the graph, thereby deriving the confidence approximation of the marginal posterior pdf of the target state, and the number of weak targets and the corresponding state estimation can be calculated using the MMSE estimator. In the information derivation process, the Goodman variable substitution principle is used to simplify the multiple integrals into ordinary integrals to simplify the derivation process, the sonar array measurement is converted into a sampling covariance matrix, the complex Wishart distribution is used to construct a likelihood function to describe the relationship between the measurement and the target state, and the Gaussian distribution is used to simplify the remaining pdf. The method of the present invention can accurately, quickly and stably track multiple dynamically changing weak target states (azimuth, azimuth angular velocity and target intensity) when the input is sonar array data.

[0071] like Figure 1 As shown, this embodiment designs an underwater weak target azimuth detection pre-tracking method based on belief propagation. The method input is sonar array data, and the output is the number of targets and the corresponding state. The proposed method is divided into three steps. Step 1 uses the joint posterior PDF to describe the weak target tracking problem; Step 2 uses the graphical model to construct the factor graph of the target tracking problem; Step 3 introduces the information derivation process and calculation flow of the belief propagation method in detail.

[0072] Using mathematical expressions to describe the weak target tracking problem, the sonar measurement information at all times is combined with the weak target state information to be solved, and the following joint posterior PDF expression is obtained:

[0073] f(y 0:k |z 0:k )

[0074] where y 0:k Represents all target states from time 0 to time k, target state N k Represents the number of targets at time k. Each augmented state The binary variable r n,k ∈{0,1} is used to indicate whether the target n exists, r n,k =1 indicates existence. is the target state information, which respectively represents the azimuth, azimuth angular velocity and target intensity. The logarithm of the target intensity is taken for the convenience of calculation. 0:k represents all sonar measurement information from time 0 to time k, z k represents the measurement at time k.

[0075] Based on the Bayesian rule, we can express f(y) by extracting the target state marginal posterior PDF from the above expression. k,n |z 0:k ), thereby estimating the number of weak targets and their corresponding states.

[0076] like Figure 2 As shown, this embodiment factors the joint posterior pdf, extracts the variable nodes (represented by circles in the figure) and the factor nodes (represented by rectangles in the figure) therein to construct a factor graph, and derives an approximate expression of the posterior pdf of the target state edge by calculating the information transmitted between the nodes of the factor graph. The expression of the factorized joint posterior pdf is as follows:

[0077]

[0078] The target state y at time k k It consists of two parts, one of which is the existence of the target y at time k-1 k-1 Through the state transfer equation f(y k,n |y k-1,n ) is calculated, and the target number is N k-1 The other part is the new target at time k, which obeys the pre-set probability distribution and has a target number of N. k -N k-1 . f(z k |y k ) is the sonar likelihood function, which is used to describe the relationship between sonar measurements and target states. Name the information transmitted between nodes. The information to be solved includes α k,n , β n and κ n .

[0079] like Figure 3 As shown, in this embodiment, the belief propagation algorithm is divided into the following steps according to the calculation order: prediction step, calculation of prediction information α k,n ; Measurement estimation step, calculate measurement estimation information β n ; Measurement update step, calculate measurement update information κn ; Confidence update step, calculate the confidence approximation of the posterior pdf of the target state edge The sonar measurement at the current moment is combined into the measurement update step through the likelihood function constructed by the Wishart distribution, and the confidence Input the MMSE estimator to get the estimate of the target state, and count the targets greater than the birth threshold to get the target number estimate. At this point, all the calculation processes at the current moment are completed. Returning to the prediction step, the underwater weak target tracking problem at the next moment can be solved iteratively.

[0080] For every N∈{1,…,n k-1}, prediction information α at time k k,n The expression is as follows:

[0081]

[0082] α k,n (x k,n ,r k,n =0) = 1-α k,n (x k,n ,1)

[0083] where p s represents the target survival probability, is the posterior pdff(x) of the target state edge at the previous moment k-1,n ,r k-1,n =1|z 0:k-1 ). For each N ′ ∈{n k-1 +1,…,n k}, new target prediction information α at time k k,n′ The expression is as follows:

[0084] α k,n′ (x k,n′ ,r k,n′ =1)=p b f(x k,n′ )

[0085] α k,n′ (x k,n′ ,r k,n′ =0) = 1-α k,n′ (x k,n′ ,1)

[0086] where p b represents the target newborn probability, f(x k,n′ ) is the pre-set new target state pdf. Prediction information α k,n ,n∈{1,…,N k} still satisfies the definition that the pdf integral sums to 1.

[0087] According to the BP algorithm, the estimated information β is measured n The calculation expression is as follows:

[0088] β n (x k,n ,r k,n )=α k,n (x k,n ,r k,n )

[0089] Similarly, we can get the measurement update information κ n The expression is as follows:

[0090]

[0091] where y k\n Indicates that time k does not contain the target state y k,n The remaining target state set is Represents the target state y k\n Using the Goodman variable substitution principle, the above multiple integrals are transformed into ordinary integrals to obtain the information κ n The simplified expression of is as follows:

[0092] κ n (x k,n ,r k,n ;z k )=∫f(z k -r k,n g(x k,n )|ξ k,n )f(ξ k,n )dξ k,n

[0093] Where f(ξ k,n ) is about The probability density function of g(x k,n ) represents the transition from the target state x k,n Nonlinear mapping from spectral to sonar measurements.

[0094] In order to obtain accurate and stable tracking results, considering that the sonar measurement values ​​at a single sampling moment have a large randomness, the complex Wishart distribution is used to describe the likelihood function f(z k -r k,n g(x k,n )|ξ k,n ). Using the sampling values ​​within a period of time effectively reduces this randomness, and the k-time measurement z k Converted into a complex sampling covariance matrix, the following expression is obtained:

[0095]

[0096] Where v represents the number of samples, assuming that the samples satisfy the independent and identically distributed assumption, and v also represents the degrees of freedom in the Wishart distribution. k,i It means that at time k, the array data is obtained by sampling the sonar for the i-th time, and the likelihood function expression is as follows:

[0097]

[0098] in represents a complex p×p dimensional symmetric complex Wishart distribution with matrix variable Y, r degrees of freedom, and symmetric positive definite scaling matrix A. σ Represents the complex covariance matrix satisfied by the environmental background noise.

[0099] Complex Wishart Distribution Scale Matrix and Sampling Covariance Matrix With a similar form, construct g(x k,n ) from the target state The nonlinear mapping to the complex covariance matrix is ​​as follows:

[0100]

[0101] Where a(θ k,n ) is the array manifold vector, which is uniquely determined by the actual physical form of the sonar array, and H represents the conjugate transpose of the matrix.

[0102] Use a Gaussian distribution to approximate pdff(ξ k,n ) to simplify the calculation, and obtain the following expression:

[0103]

[0104] in The variable is Y, the mean is μ, and the covariance is C. The relevant calculation is expressed as follows:

[0105]

[0106] Rewrite information k n The expression is as follows:

[0107]

[0108] Among them use Replace k,n This is to simplify the derivation process. Based on the information obtained, the target state edge posterior pdff(y k,n |z 0:k ) is expressed as follows:

[0109]

[0110] in is the normalization factor, ensuring that the sum of the confidence intervals satisfies the requirement that the integral of the probability density function sums to 1.

[0111] Based on the calculated confidence Use MMSE to estimate the number of weak targets and their corresponding states. k}, first calculate the probability of target existence:

[0112]

[0113] like Greater than the preset target existence threshold That is, the target n is considered to exist at time k, and the target state is calculated by the MMSE estimator as follows:

[0114]

[0115] Statistics in n∈{1,…,N k}The probability of the target existing in The number of weak targets in the area is estimated by the number of weak targets in the area.

[0116] The calculation process of the method of the present invention is as follows: Figure 3 As shown, it is divided into the following steps.

[0117] (1) Based on confidence The information α is calculated by the prediction step k,n ;

[0118] (2) Based on the prediction information α k,n , the information β is calculated by measuring the estimation step n ;

[0119] (3) Estimation information β based on measurement n , input the current measurement z k , calculate the information k by measuring the update step n ;

[0120] (4) Based on the prediction information α k,n , measurement update information κ n , the confidence is calculated through the confidence update step

[0121]

[0122] (5) Based on confidence Calculate the number of targets and their corresponding states through the MMSE estimator;

[0123] (6) Complete the current moment calculation process, Return to step (1) to start the calculation for the next moment.

Claims

1. A method for tracking underwater weak targets before azimuth detection based on belief propagation, characterized in that: The following steps are involved: 1) Based on the Bayesian rule, the joint posterior probability density function pdf is used to describe the weak target tracking problem; 2) Factorize the joint posterior pdf and construct the corresponding factor graph through the graphical model; 3) Use the belief propagation BP algorithm to sequentially calculate the prediction information, measurement estimation information and measurement update information in the factor graph to achieve weak target tracking before azimuth detection.

2. The underwater weak target tracking method before azimuth detection based on belief propagation according to claim 1 is characterized in that: The step 1) is specifically as follows: Combining the sonar measurement information at all times with the weak target state information to be solved, we get the joint posterior PDF expression: f(y 0:k |z 0:k ) Among them, y 0:k Represents all target states from time 0 to time k, target state N k Represents the number of targets at time k, each augmented state The binary variable r j,n ∈{0,1} is used to indicate whether the target n exists, r j,n =1 means existence, is the target state information, respectively representing the orientation θ k,n , azimuth velocity With target strength z 0:k represents all sonar measurement information from time 0 to time k, z k represents the measurement at time k.

3. The underwater weak target tracking method before azimuth detection based on belief propagation according to claim 1 is characterized in that: In step 2), the factorized joint posterior pdf expression is: Among them, the target state y at time k k It consists of two parts, one of which is the existence of the target y at time k-1 k-1 Through the state transfer equation f(y k,n |y k-1,n ) is calculated, and the target number is N k-1 The other part is the new target at time k, which obeys the pre-set probability distribution and the number of targets is N. k -N k-1 , f(z k |y k ) is the sonar likelihood function used to describe the relationship between sonar measurements and target states.

4. The underwater weak target tracking method before azimuth detection based on belief propagation according to claim 1, characterized in that: The step 3) comprises the following steps: 3.1) Based on confidence Calculate prediction information α k,n ; 3.2) Based on prediction information α k,n , using BP algorithm to calculate the measurement estimation information β n ; 3.3) Estimation of information β based on measurement n , input the current measurement z k , calculate the measurement update information κ n ; 3.4) Based on prediction information α k,n , measurement update information κ m , calculate the approximate confidence of the posterior pdf of the target state margin 3.5) Based on confidence Calculate the number of targets and their corresponding states through the MMSE estimator; 3.6) Complete the current moment calculation process and Return to step 3.1) to start the next moment calculation.

5. The underwater weak target tracking method before azimuth detection based on belief propagation according to claim 4 is characterized in that: The step 3.1) is specifically as follows: For each target n∈{1,…,N k-1 }, prediction information α at time k k,n for: a k,n (x k,n ,r l,n =0)=1-a k,n (x k,n ,r k,n =1) Among them, p s represents the target survival probability, represents the marginal posterior pdff(x) of the target state at the previous moment k-1,n ,r k-1,n =1|z 0:k-1 )’s approximate confidence; For each new target n′∈{N k-1 +1,…,N k }, new target prediction information α at time k k,n′ for: α k,n′ (x k,n′ ,r k,n′ =1)=p b f(x k,n′ ) a k,n′ (x k,n′ ,r k,n′ =0)=1-a k,n′ (x k,n′ ,r k,n′ =1) Among them, p b represents the target newborn probability, f(x k,n′ ) is the pre-set new target state pdf.

6. The underwater weak target tracking method before azimuth detection based on belief propagation according to claim 4 is characterized in that: The step 3.2) is specifically as follows: b n (x k,n ,r k,n )=a k,n (x k,n ,r k,n )。 7. The underwater weak target tracking method before azimuth detection based on belief propagation according to claim 4 is characterized in that: The step 3.3) is specifically as follows: κ n (x k,n ,r k,n ;z k )=∫f(z k -r k,n g(x k,n )|ξ k,n )f(ξ k,n )dξ k,n Among them, f(ξ k,n ) is about custom variables The probability density function, y k\n Indicates that time k does not contain the target state y k,n The remaining target state set, for n∈{1,…,N k },g(y k,n ) represents the transition from the target state x k,n nonlinear mapping from θ to sonar measurements, Among them, a(θ k,n ) represents the array manifold vector, and H represents the conjugate transpose of the matrix.

8. The underwater weak target tracking method before azimuth detection based on belief propagation according to claim 4 is characterized in that: The step 3.4) is specifically as follows: Use complex Wishart distribution to describe the likelihood function f(z k -r k,n g(x k,n )|ξ k,n ), measure z at time k k Converted to complex sampling covariance matrix Where v represents the number of samples, z k,i It means that the array data is obtained by sampling the sonar for the i-th time at time k, and the likelihood function expression is: in, represents a p-dimensional complex Wishart distribution with a matrix variable Y, r degrees of freedom, and a scale matrix A, R σ Represents the complex covariance matrix satisfied by the environmental background noise; Use a Gaussian distribution to approximate pdff(ξ k,n ) to simplify the calculation, we get: in represents a Gaussian distribution with a variable Y, a mean of μ, and a covariance of C, where the mean and covariance Specifically: Rewrite information k n The expression is: Based on the BP algorithm, we get the target state edge posterior pdff(y k,n |z 0:k ) is: in, is the normalization factor.

9. The underwater weak target tracking method before azimuth detection based on belief propagation according to claim 4, characterized in that: The step 3.5) is specifically as follows: For n∈{1,…,N k }, calculate the probability of target existence like Greater than the preset target existence threshold The target n is considered to exist at time k, and the target state is calculated by the MMSE estimator Statistics in n∈{1,…,N k }The probability of the target existing in The number of weak targets in the area is estimated.

Images