Multi-Target Tracking by Associating and Fusing the Multi-Bernoulli Parameter Sets

The Cardinality Balanced MeMBer (CBMeMBer) filter is a single sensor multi-target tracking method based on the random finite set. Compared with the single sensor system, the multi-sensor system can achieve more stable and better performance in tracking targets. However, some problems exist in multi-sensor system based on the CBMeMBer filter. Tracks in the CBMeMBer filter are described by parameter sets which may be generated by miss-detection, targets and clutters. It is difficult to associate the parameter sets correctly because of their complex forms and various types. Moreover, the filter reacts slowly to disappeared targets, which leads to a cardinality overestimation. This problem may be more serious in the multi-sensor system. To deal with the above problems, the parameter sets association and fusion methods are presented in this paper. By three association processes with the adaptive thresholds selection approaches, parameter sets corresponding to the same target are grouped into one parameter set partition. Parameter sets have different association thresholds because of their different accuracies. The fusion method considers the types and relationships of parameter sets in the partition simultaneously and uses a joint credibility to accelerate changes in existence probability. The cardinality estimation decreases rapidly when the target disappears. The theoretical analysis and experiment results of different tracking scenarios show that the proposed methods perform well in both state estimation and cardinality estimation.


I. INTRODUCTION
Target tracking aims to estimate the target states by the measurements of sensors, such as the number, position, velocity, acceleration, and track. According to the number of targets, target tracking can be divided into two types, Single Target Tracking (STT) and Multiple Target Tracking (MTT). STT is an early and well-developed research in this field. But STT approaches do not achieve well performance for the real tracking scenarios because the number of targets is time-varying and unknown. Regarding this issue, many MTT methods have been proposed, such as Nearest Neighbor [1], Global Nearest Neighbor [2], Probabilistic Data Association [3], Joint Probabilistic Data Association [4] and Multiple Hypothesis Tracking [5].
The associate editor coordinating the review of this manuscript and approving it for publication was Rui-Jun Yan . All these classical MTT methods mentioned above are based on data association, and their computation burden increases rapidly with the increasing numbers of targets and measurements. To address this issue, several MTT methods [7]- [18] without data association are proposed by combining the Random Finite Set (RFS) theory [6] and Bayesian filtering framework, such as the Probability Hypothesis Density (PHD) filter [7], the Cardinalized Probability Hypothesis Density (CPHD) filter [8], and the Multi-target Multi-Bernoulli (MeMBer) filter [6]. In the PHD filter, the integral of the probability hypothesis density function in each region is the expectation of target numbers in the region. Since only the first moment of the multi-target posterior probability density is propagated in the PHD filter, the cardinality estimation is unstable when there is a miss-detection or a high false alarm. To deal with this problem, the CPHD filter is proposed by propagating the posterior intensity and cardinality distribution simultaneously. Compared with the VOLUME 8, 2020 This work is licensed under a Creative Commons Attribution 4.0 License. For more information, see https://creativecommons.org/licenses/by/4.0/ PHD filter, the CPHD filter performs better in cardinality estimation, but it has a larger computational burden. The computational complexities of the PHD filter and the CPHD filter are O (n · m) and O n · m 3 , respectively. Here, n and m are the numbers of targets and measurements, respectively. The number of measurements is the number of echoes received by the sensor within one sampling time step. These echoes may originate from targets or clutters. Unlike the PHD filter which propagates the first-order moment, the MeMBer filter propagates parameters of the multi-Bernoulli RFS, and it is applied for a low clutter density condition. The cardinality overestimation problem of the MeMBer filter is proved, and a modified version is presented named as the Cardinality Balanced MeMBer (CBMeMBer) filter [9]. Besides these commonly used RFS filters, some improved methods have been proposed to track targets in real tracking scenarios. In [10] and [11], the Student's t-distribution is used to handle the problems of observation noise of the PHD filter and the CBMeMBer filter, respectively. Some solutions to the adaptive target birth are described in [12] and [13].
To track maneuvering targets, a PHD filter with multiple models is presented in [14]. In [15], the PHD filter is used to track extended targets. The miss-detection problems of the CPHD filter are discussed in [16]. To implement these filters, the Gaussian Mixture (GM) implementation and the Sequential Monte Carlo (SMC) implementation are given in [17] and [18], respectively. In Bayesian probability theory, conjugate prior is important which means that the prior distribution has the same form as the posterior distribution. The generalized labeled multi-Bernoulli (GLMB) filter [19], [20] is based on conjugate prior using labeled RFSs. Different from the PHD, CPHD and CBMeMBer filters, the GLMB has a mixture representation, where each component in the mixture corresponds to one possible data association history. Although the GLMB filter has better performance in challenging scenarios, it is computationally expensive. A computationally efficient approximation is the Labeled Multi-Bernoulli (LMB) filter [21], which approximates the GLMB density with a labeled multi-Bernoulli density. The classical tracking methods and the RFS filters mentioned above are mainly used for a single sensor tracking system. To increase tracking accuracy and tracking system stability, multi-sensor tracking methods have been gaining significant attention. The Iterated Corrector PHD (IC-PHD) filter, an approximation for the multi-sensor case, is proposed by Mahler in [7]. This filter is easily implemented and wildly used in many applications, such as registration errors compensation [22], radiation source localization [23] and tracking multiple speakers [24]. However, the IC-PHD filter is sensitive to the sensor order and the probability of detection. Detailed discussions can be found in [25]. To deal with this problem, a joint miss-detection probability and joint detection probability calculation method is described in [26]. As an improved version of the IC-PHD filter, the Product Multi-sensor PHD (PM-PHD) is proposed in [27]. The PM-PHD filter can overcome the drawbacks of the IC-PHD filter, but some problems still exist. It is difficult to calculate or approximate the infinite sums in the update formulation. Furthermore, the PM-PHD filter is overly concerned with the variation of cardinality estimation, and it results in an inaccurate state estimation sometimes. In [28], an approximation method of infinite sums is presented, where the computation can be reduced effectively. An improved PM-PHD filter presented in [29] separates the state estimation process and cardinality estimation process, and it weakens the influence of cardinality estimation on state estimation. Besides the IC-PHD filter and the PM-PHD filter, Mahler proposed another multi-sensor RFS filter, named as the General two-sensor PHD filter [30]. A summation overall binary partitions of the measurement sets makes the computation of this filter high. In [31], D. Clark gave the multi-sensor version of the General two sensor PHD filter, and Nannuru gave its Gaussian implementation in [32]. The computational complexity of the general multi-sensor PHD filter is O (m! · n) approximately, which leads to a limitation of application.
In the above-mentioned multi-sensor RFS filters, fusion center needs to divide measurements of sensors into different measurement partitions and fuse the multi-target state estimations of all partitions. It leads to a high computation burden of the fusion center. A feasible solution to this problem is to fuse the states of targets estimated by each sensor in the fusion center. Since the results of the RFS filters are in the form of probability hypothesis density, much useful information is missed in the process of extracting the state of targets. It may result in an inaccurate target state estimation. In this paper, a heuristic multi-sensor multi-target tracking method based on the CBMeMBer filter is proposed to fuse the filter results without extracting the states of targets. This approach consists of two parts, parameter sets association and parameter sets fusion. Its main contributions are given below: (1) Through an association framework consisting of three association processes, parameter sets of miss-detection and measurement update belonging to the same target can be accurately associated to the same parameter set partition.
(2) The threshold used by each association process has its corresponding adaptive selection method. These thresholds can improve the efficiency and accuracy of the association processes.
(3) To deal with the cardinality overestimation problem caused by the excessive dependence of the CBMeMBer filter on the miss-detection parameter sets at a low probability of detection, the credibility of existence probability is introduced to modify the fused existence probability.
The rest of this paper is organized as follows. The CBMeMBer filter and the Gaussian mixture implementation are introduced in Section II. The framework of a multi-sensor association method is described in Section III. The numerical solutions to the multi-sensor association method are given in Section IV. The multi-sensor parameter set estimation method and theoretical analysis are given in Section V and Section VI, separately. Simulation results are shown in Section VII. Finally, Section VIII gives conclusions and future works.

A. THE CBMEMBER FILTER
The Bernoulli RFS is a binary distribution. The probability that it is an empty set is 1 − r, and the probability that it is a single state x is r. The probability density distribution of x is p (x). Thus, the probability density of a Bernoulli RFS can be expressed by The Multi-Bernoulli RFS consists of several Bernoulli RFS and can be expressed as The probability density ofX is Here, the parameter set r (i) , p (i) is corresponding to the Bernoulli RFS X (i) . The CBMeMBer filter consists of two parts, the prediction and the update.

1) THE CBMeMBer PREDICTION
Given the multi-target posterior probability density at time k − 1, The predicted probability density at time k can be expressed by Here, r are the parameter sets of targets that still survive at time k. r are the parameter sets of new targets that appear instantly at time k. f k|k−1 (x |ξ ) is the state transition function. p S,k (x) is the probability of survival.
In Eq. (5), the number of the predicted Bernoulli RFS is M k|k−1 = M k−1 + M ,k , and π k|k−1 can be rewritten as

2) THE CBMeMBer UPDATE
Given the multi-target predicted probability density π k|k−1 at time k, the multi-target posterior probability density at time k can be expressed by Here, subscripts L and U denote the miss-detection and measurement update. κ k (z) denotes the clutter intensity. p D,k (x) is the probability of detection. g k (z | x) is the observation likelihood function. Z denotes the measurement set. z denotes the measurement in Z. Note that the echo received by sensor is processed signal detection process and generate the measurement z. The form of z can be the planar position in the Cartesian coordinate or the range and azimuth in Polar coordinate.

B. THE GAUSSIAN MIXTURE IMPLEMENTATION
The GM implementation of the CBMeMBer filter obeys the following assumptions [9], 1) f k|k−1 (x |ξ )and g k (z| x) of each target follow a linear Gaussian model, Here, F k−1 and H k represent the state transition matrix and the observation matrix, respectively. Q k−1 and R k are the state noise covariance matrix and the observation noise covariance matrix for the state x.
2) The probabilities of survival and detection are independent to state, 3) The probability density p can be represented in the form of VOLUME 8, 2020 Gaussian mixture, Here, J (i) γ ,k is the number of Gaussian components corresponding to the i th birth target at time k. w ,k , j = 1, 2, . . . , J ,k denote the weight, the state mean vector and the covariance matrix of the j th Gaussian component, respectively.
Based on these assumptions, the GM implementation of the CBMeMBer filter is given as follows, (

III. THE FRAMEWORK OF MULTI-SENSOR PARAMETER SETS ASSOCIATION
Assume that there is a homogeneous sensor network with S sensors, such as a radar network composed of multiple active radars. The sensor network has been calibrated in time and space. The filtering result of the th , = 1, . . . , S sensor at time k is (35) Here, N is the number of parameter sets of the th sensor, and N = M k|k−1 + Z ,k . Z ,k denotes the measurement set of the th sensor at time k. Z ,k denotes the number of measurements. It can be seen from Eq. (9) that the parameter sets in Eq. (35) can be divided into two parts, the miss-detection parameter sets, and the measurement update parameter sets. Therefore, the multi-sensor parameter sets association consists of three parts, the association of miss-detection parameter sets, the association of measurement update parameter sets, and the association between the parameter sets of miss-detection and measurement update.

A. THE ASSOCIATION OF MISS-DETECTION PARAMETER SETS
From Eqs. (10) and (11), it can be seen that a miss-detection parameter set is only determined by one predicted parameter set. Thus, the predicted parameter set can be used to associate the parameter sets generated by it.
According to [33], the probability hypothesis density of Assume that r WhenQ r After the association between the miss-detection parameter sets of each sensor and the predicted parameter sets, the missdetection parameter sets associated with the same predicted parameter set are grouped into one parameter set partition.
For example, assume that there are two sensors, the predicted parameter sets are r (1) k|k−1 , p (1) k|k−1 and r (2) k|k−1 , p (2) k|k−1 . The parameter sets of the first sensor and the second sensor are r 2,k , r 2,k , p 2,k , r 2,k , p Here, r 1,k and r 2,k are generated by r (2) k|k−1 , p (2) k|k−1 . Then, the miss-detection parameter sets partitions are

B. THE ASSOCIATION OF MEASUREMENT UPDATE PARAMETER SETS
After the association of miss-detection parameter sets, the associated parameter sets are removed from r . . , S. The rest parameter sets are the measurement update parameter sets. In this sub-section, are used to denote the measurement update parameter sets of the th sensor.
The difference between two measurement update parameter sets can be calculated by Eq. (38). For the measurement update parameter sets of all sensors, their difference matrix is given by Generally, the difference between two parameter sets can be effectively measured by the KLD, but some details between them may be ignored. In the situation shown in Fig. 1,Q (r, p), r 1 , p 1 is equal toQ (r, p), r 2 , p 2 . It is difficult to use KLD to determine whether (r, p) should be associated with r 1 , p 1 or r 2 , p 2 . If we consider the processes of parameter set generating, state extracting and state merging, (r, p) should be associated with r 1 , p 1 . The explanations are given as follows. 1) In Eq. (13), measurement z is used to update all the predicted probability densities p The updated probability density p U ,k (x; z) of z is a summation of the updated probability densities corresponding to p are the measurement and predicted probability density of a target, the updated probability density generated by them consists of a component with a large amplitude and several components with small amplitudes. Moreover, the updated probability densities generated by z and other predicted probability densities consist of several components with tiny amplitudes.
2) The state extraction is mainly depended on the component with a large amplitude. Furthermore, the components with small and tiny amplitudes have a negative effect on state extraction.
3) If the difference between two probability densities is big, these two probability densities may generate an inaccurate merged state, or cannot even merge.
To deal with the above problem, another difference matrix given in Eq. (42) is used to describe the difference VOLUME 8, 2020 between the components with large amplitudes of two parameter sets.
Here,R i,i , is the KLD between the components with large amplitudes of r ,U ,k cannot be associated.

C. THE ASSOCIATION BETWEEN THE PARTITIONS OF MISS-DETECTION AND MEASUREMENT UPDATE
After the above two associations, there are miss-detection parameter set partitions P L,k and several measurement update parameter set partitions P U ,k , where the number of P L,k is M k|k−1 . P U ,k can be divided into two categories, complete partitions P Com U ,k and incomplete partitions P Im U ,k . P Com U ,k : the number of parameter sets in this partition equals to the number of sensors S. P Im U ,k : the number of parameter sets in this partition is less than the number of sensors S.
There are many reasons for having incomplete partitions. For example, targets are undetected by some sensor; the parameter sets of the same target cannot be grouped into one partition because a large observation noise; the parameter sets of clutters are divided into separate partitions.
It is easy to compute the difference between P L,k and P Com U ,k by Eq. (43).
However, Eq. (44) is not suitable to calculate the difference between P L,k and P Im U ,k , because P L,k − P Im U ,k = 0. Here, |P| denotes the number of parameter sets in P. Inspired by the optimal sub-pattern assignment (OSPA) [35] distance, the difference of the numbers of parameter sets is computed by Here, φ is the cost. Then, the difference between P L,k and P U ,k is calculated by Eq. (46), as shown at the bottom of this page.
If the target is undetected by all sensors, Eq. (46) can be rewritten asQ WhenQ P L,k , P U ,k > S · φ, P L,k and P U ,k cannot be associated. It indicates that P U ,k does not correspond to the target.

IV. THE NUMERICAL SOLUTION TO THE MULTI-SENSOR ASSOCIATION METHOD
Section III mainly gives the framework of the multi-sensor association method. However, some problems still exist and need to be handled by the numerical methods, such as the solutions of each association and threshold selection. Based on the GM implementation of the CBMeMBer filter, solutions to these problems are given in this section.
Here, N (i) andÑ (i) are the Gaussian components of , and its closed-form expression is Eq. (53), as shown at the bottom of this page.
Furthermore, assume that r k|k−1 can be obtained by rewriting Eq. (25), and we have Eqs. (54) and (26) indicate that r When r Thus, the threshold γ is set as zero. By Eqs. (52) and (53), there is a difference matrix ϒ of the parameter sets in r and the predicted parameter sets r The row and column of ϒ correspond to r . . , S, respectively. The process of associating the miss-detection parameter sets is given in Table 1.  (28) and (36), ,U ,k can be expressed as Eq. (59).
The differenceR i,i , between the components with large amplitudes of r If z ∈ Z is the measurement of a true target, a Gaussian component has the weight which is much larger than the others' in Eq. (59). Generally, the Gaussian component with the largest weight can be used to describe the target, and it is the linear Gaussian form of the component with a large amplitude. In the rest of this sub-section, the Gaussian component with the largest weight is abbreviated as the main Gaussian component. If r ,U ,k correspond to the same target, the difference between their main Gaussian components should be as small as possible. Thus, the main Gaussian component can be used to compute the differenceR i,i , between the components with large amplitudes. Assume that N x; m respectively. An example is illustrated in Fig. 2.
,U ,k and The proof of Eq. (62) is given in Appendix A.
For N x; m The process of associating the measurement update parameter sets is given in Table 2.  Solution (1) in this sub-section, it is observed that P L,k can be replaced by its corresponding P k|k−1 .
Then, the association of P L,k and P U ,k is changed to associate P k|k−1 and P U ,k .Q P L,k , P U ,k in Eq. (46) is replaced byQ P k|k−1 , P U ,k . However, P U ,k may be an incomplete partition P Im U ,k . In Eq. (45), φ is introduced to describe the difference between the numbers of parameter sets in two partitions. The selection of φ is discussed as follows.
Suppose that D k|k−1 (x) is the probability hypothesis of r k|k−1 , p k|k−1 , and N x; m k|k−1 , P k|k−1 is the main Gaussian component of D k|k−1 (x). The measurement updating process is illustrated in Fig. 3.
In Fig. 3, N x; H k m k|k−1 , H k P k|k−1 H T k is the projection of N x; m k|k−1 , P k|k−1 in x − y plane. Here, x and y denote the position components of target states. k|k−1 denotes the 3σ region of N x; H k m k|k−1 , H k P k|k−1 H T k . z ,T , z ,C and z ,B are the measurements of the th sensor. N z ,T µ, σ 2 and N z ,C µ, σ 2 are the projections of the Gaussian components generated by z ,T and z ,C , respectively. In the measurement updating process (refers to Eqs. (27) and (28)), the measurement z ,T ∈ k|k−1 has a large measurement likelihood (refers to Eq. (30)), because z ,T is likely to be the measurement of targets. Thus, N z ,T µ, σ 2 has a large weight. For z ,C / ∈ k|k−1 , it likely to be the measurement of clutters or other targets. Thus, the measurement likelihood of z ,C is very small, and N z ,C µ, σ 2 has a tiny weight. If there is no measurement in k|k−1 , which means that the target related to P k|k−1 is undetected by the th sensor, P U ,k becomes an incomplete partition P Im U ,k . Considering an extreme situation, there is a measurement z ,B on the boundary of k|k−1 . Compared with the difference between N x; m k|k−1 , P k|k−1 and the Gaussian components generated by the measurements in k|k−1 , the difference between N x; m k|k−1 , P k|k−1 and the Gaussian component generated by measurements in k|k−1 (such as z ,B ) is larger. Thus, the difference corresponding to z ,B can be used to judge whether there is a measurement of targets or not. Therefore, φ can be obtained by Eq. (65).  VOLUME 8, 2020 In Eq. (67), we havẽ Then, the threshold of P (i) L,k and P (i ) U ,k in Eq. (47) is given by Thus, we have a threshold matrix , which is given in Eq. (70).
Here, M P U ,k is the number of measurement update partitions. The process of associating the parameter sets of miss-detection and measurement update is given in Table 3. In Table 3, the initial value of C

V. THE MULTI-SENSOR PARAMETER SET ESTIMATION METHOD
After the association between the partition of miss-detection and measurement update, there are many combinations of P L,k and P U ,k . These combinations can be divided into four categories: Category (1): P L,k associates with a P Com U ,k . It indicates that the target corresponding to P L,k is detected by all sensors.
Category (2): P L,k associates with a P Im U ,k . It indicates that the target corresponding to P L,k is detected by some sensors.
Category (3): There is no P U ,k associate with P L,k . It indicates that the target corresponding to P L,k is undetected by all sensors, or has been dead.
Category (4): There is no P L,k associate with P U ,k . It indicates that P U ,k corresponds to clutters.
Thus, the parameter sets in Category (4) will be deleted. Combinations in Category (1)  Given r ,L,k , p ,L,k , r ,U ,k , p ,U ,k of the th sensor, r ,L,k , p ,L,k and r ,U ,k , p ,U ,k are the miss-detection part and measurement update part of one target, respectively.
If r ,U ,k , p ,U ,k = ∅, the existence probability of the target estimated by the th sensor is The probability density distribution of the target estimated by the th sensor is If r ,U ,k , p ,U ,k = ∅, we have Here, r ,L,k r ,U ,k , p ,U ,k = ∅ or r ,L,k + r ,U ,k > 1 r ,L,k + r ,U ,k r ,U ,k , p ,U ,k = ∅ and r ,L,k + r ,U ,k ≤ 1 82718 VOLUME 8, 2020 p = p ,L,k r ,U ,k , p ,U ,k = ∅ or r ,L,k + r ,U ,k > 1 p ,U ,k r ,U ,k , p ,U ,k = ∅ and r ,L,k + r ,U ,k ≤ 1 (79) From Eq. (78) and Eq. (79), it can be seen that the parameter set combination uses the existence probability and probability density distribution when r ,U ,k , p ,U ,k = ∅ or r ,L,k + r ,U ,k > 1. This is because the target is represented by the existence probability and probability density distribution of miss-detection parameter set when it is undetected by the th sensor. Moreover, r ,L,k + r ,U ,k > 1 indicates that combination r ,L,k , p ,L,k , r ,U ,k , p ,U ,k is incorrect. Since a miss-detection parameter set must correspond to one target and the measurement update parameter set may correspond to clutters, the former is more reliable than the latter. r ,U ,k , p ,U ,k = ∅ and r ,L,k + r ,U ,k ≤ 1 indicates that r ,L,k , p ,L,k and r ,U ,k , p ,U ,k may belong to the same target. However, when the probability of detection is low, even if the target is detected by the th sensor, r ,L,k is larger than the normal value and r ,U ,k is small than the normal value. Thus, r ,U ,k cannot accurately represent the existence probability of the target. It can be seen from Eq. (25) and Eq. (27) that r ,L,k + r ,U ,k ≈ 1 when r ,L,k and r ,U ,k belong to the same target. Compared to r ,U ,k , r ,L,k + r ,U ,k can represent the existence probability of the target more accurately.
Based on Eqs. (77)-(79), the existence probability and probability density distribution of the combinations in Categories (1)-(3) are calculated by the following approaches.

A. ESTIMATING THE EXISTENCE PROBABILITY
Assume that X (i) , i = 1, . . . , N X are the Bernoulli RFSs, and they are independent of each other. r (i) , p (i) is the parameter set of X (i) . The joint probability of X (i) , i = 1, . . . , N X is given by Eq. (80).
When X = ∅, Eq. (82) is the probability that there is no target.

P {There is no target}
When X = {x}, Eq. (83) is the probability that there is a target.
P {There is a target} Thus, the existence probability of X is Given a combination C = {(r , p )} S =1 , the existence probability of C can be computed by Eq. (85).
In some situations, Eq. (85) overestimates the cardinality. For example, a target survivals at time k − 1 with existence probability r target k−1 . Then, the predicted existence probability r target k|k−1 of this target is computed by Eq. (20). At time k, this target is undetected by the th sensor or has been dead, and its existence probability r target ,L,k is computed by Eq. (25). The curves of r target ,L,k are illustrated in Fig. 4. In the situation discussed above, r target k−1 and r target k|k−1 have a large value. When the target is undetected by the th sensor, VOLUME 8, 2020 it can be seen as Fig. 4 that r target ,L,k is large. Then, the target still can be estimated. If the target has been dead since time k, r target ,L,k should be very small or even equals to zero. However, as shown in Fig. 4, if r target k|k−1 is large, r target ,L,k has a very small value only when p ,D,k is very large, such as p ,D,k = 0.99. It may lead to cardinality overestimation and false alarms. Fortunately, this negative effect will disappear after one or two time steps in the CBMeMBer filter. But in Eq. (85), this problem becomes worse.
Assume that the probabilities of detention of all sensors are the same, and p ,D,k = p D,k , = 1, . . . , S. Thus, r    5 shows that r C increases rapidly with the increasing of r target L,k and the sensor numbers S. It indicates that the target has a large existence probability, even it is undetected by all sensors. To deal with this problem, the joint credibility of r C is introduced.
For a parameter sets partition r ,L,k , p ,L,k , r ,U ,k , p ,U ,k , if there is a target corresponding to r ,L,k , p ,L,k , r ,U ,k , p ,U ,k , the probability of r ,U ,k , p ,U ,k = ∅ is 1−p ,D,k . Therefore, it is inferred that if r ,U ,k , p ,U ,k = ∅ and there is a target, the credibility can be approximated as p ,D,k . If r ,U ,k , p ,U ,k = ∅ and there is no target, the credibility can be approximated as 1 − p ,D,k . Moreover, if r ,U ,k , p ,U ,k = ∅ and there is a target, the credibility can be approximated as 1 − p ,D,k . If r ,U ,k , p ,U ,k = ∅ and there is no target, the credibility can be approximated as p ,D,k . Then, the joint credibility of that there is no target is given by Eq. (86), as shown at the bottom of this page. The joint credibility of that there is a target is given by Eq. (87), as shown at the bottom of this page.
Then, the joint credibility of r C is given by Based on Eqs. (85) and (88), the corrected r C is given bỹ The probability hypothesis density D C L,U (x) of C can be obtained by By Eqs. (50) and (59), D C (x) can be expressed by the linear Gaussian model in Eq. (91), as shown at the bottom of the next page.
Then, the Gaussian components in Eq. (91) can be merged by the merging method in [18]. The merging threshold of two Gaussian components given in [18] is constant. In this subsection, the merging threshold is obtained by Eq. (62). After the merging process, the probability density distribution of p C (x) is obtained by removing the cardinality of D C (x), The flowchart of the proposed algorithm, the Parameter Sets Association CBMeMBer (PSA-CBMeMBer) filter is given in Fig. 6.

VI. THE THEORETICAL ANALYSIS
Assume that the probabilities of detention of all sensors are the same, and p ,D,k = p D,k , = 1, . . . , S. Then, Eq. (89) can be expressed bỹ Here, n L and n D are the number of sensors that miss detect or detect a target, and n L + n D = S. Eq. (94) and Eq. (95) can be rewritten by Here, Here, α ∈ [0, 1]. A smaller α means that more sensors detect the target, and α = 0 means that the target is detected by all sensors. Moreover, a larger α means that more sensors miss-detect the target, and α = 1 means that the target is undetected by all sensors.
When → 1, which means that the target is detected by a sensor and it has an accurate measurement. Thus, r C → 1, and Eq. (93) can be rewritten bỹ When 0 < < 1, which means the target is detected by a sensor, but it has an inaccurate measurement. A smaller indicates that the measurement is less accurate. Thus, Eq. (93) can be rewritten by Eq. (102), as shown at the bottom of this page.
When → 0, r U ,k → 0, which means that the target is undetected by all sensors. Thus, Eq. (93) can be rewritten bỹ Fig. 7-Fig. 9 show the curves ofr C with different . Fig. 7 shows that when p D,k ≥ 0.5, the associated parameter sets have a larger C as long as the target is detected by half of the sensors. Then, the target can be estimated by the PSA-CBMeMBer filter. When p D,k < 0.5 and α < 0.5, the PSA-CBMeMBer filter considers that Q r C of the associated parameter sets corresponding to the target is small, and thusr C is also small. The larger p D,k and α indicate that the target is disappear. Therefore, when α > 0.5 and p D,k > 0.5, r C is small. In addition, when α > 0.5 and p D,k < 0.5, the target is likely to be undetected due to a low detection probability. Thus,r C is large. It seems that the algorithm loses some information, but it can avoid the impact of inaccurate information on the results. Except for the existence probability that α ≤ 0.5 and p D,k ≥ 0.5, the other parts in Fig. 7 and Fig. 8 are the same. This is because these parts mainly depend on the miss-detection. It can be seen from Figs. 8(a)-(c) thatr C becomes large with an increasing . Besides,r C is very small when the p D,k is large. The reason is that the CBMeMBer filter considers that the target is likely to be detected with a large p D,k , and thus r L,k is small. Meanwhile, because the measurement is inaccurate ( is small), r U ,k is small. Then, the existence probabilities estimated by all sensors are small, and the PSA-CBMeMBer filter considers that the target may not exist. Because the parameter sets with a small have large errors, and it cannot be associated with other parameter sets. Therefore, when is small, the PSA-CBMeMBer filter considers that the target is undetected by sensors. Fig.9 shows that the target only exists with a small p D,k , if it is undetected by all sensors. Therefore, the target can be estimated by the PSA-CBMeMBer filter only when p D,k ≤ 0.5 and r target k−1 is high.

VII. SIMULATION RESULTS
In this section, the simulation experiments are studied to analyze the proposed method. The state transition model is described as Here, x k = [x k ,ẋ k , y k ,ẏ k ] T is the target state vector, x k and y k represent the planar position coordinates of theẋ k andẏ k represent their velocities, and we have  Here, T = 1s is the sampling period, and σ w is the standard deviation of w k .
The measurement obtained by each sensor is given in Cartesian coordinates as follow Here, the covariance matrix of observation noise v ,k is To reduce the amount of calculation, the Gaussian components with negligible weights will be pruned. Moreover, the number of Gaussian components is limited. Based on experience gained from simulations [17], the pruning threshold is selected as T p = 1 × 10 −5 , respectively, and the maximum number of Gaussian components is J max = 100.
The OSPA distance is employed to evaluate the multi-target tracking performance, Here, X = {x 1 , · · · , x m } and Y = {y 1 , · · · , y n } are arbitrary finite subsets, 1 ≤ p < ∞, c > 0 (see [35] for the meanings of these parameters). If m > n,d   Table 4. The target trajectories are shown in Fig. 10. The intensity of the clutter RFS is assumed to be Here, λ ,k = 10 is the clutter rate of the th sensor, and U (·) is the uniform density over the surveillance region. In this experiment, the PSA-CBMeMBer filter is compared with the uncorrected PSA-CBMeMBer filter. The estimated number of targets, OSPA distance and mean squared error (MSE) of target numbers are shown in Fig. 11(a)-(c), respectively.
In Fig. 11(a), the numbers of targets of the uncorrected PSA-CBMeMBer filter and the PSA-CBMeMBer filter are obtained by summing the existence probabilities in Eqs. (85) and (89), respectively. As shown in Fig. 11(a), the PSA-CBMeMBer filter has a correct estimation of target numbers, but the number of targets estimated by the uncorrected PSA-CBMeMBer filter still equals to 4 after time k = 16. Moreover, the number of targets estimated by the PSA-CBMeMBer filter is slightly smaller than the true numbers of targets at each time step. By Fig. 11(b) and Fig. 11(c), it can be seen that the OSPA distance and MSE of the PSA-CBMeMBer filter are worse than those of the uncorrected PSA-CBMeMBer filter before time k = 16s. The reason is given as follows.
The curves of Q r C is illustrated in Fig. 7. Q r C can be divided into four situations, (1) Large p D,k with small α. (2) Large p D,k with large α. (3) Small p D,k with large α. (4) Small p D,k with large α. When α = 0 and p D,k > 0.5, or α = 1 and p D,k < 0.5, there is Q r C ≈ 1. When α = 0 and p D,k < 0.5, or α = 1 and p D,k > 0.5, there is Q r C ≈ 0. However, it is difficult to obtain above extreme situations. Usually, 0 < Q r C < 1. Thus,r C < r C . Furthermore, in either case, the uncorrected PSA-CBMeMBer filter considers there is always a target. When the target is only detected by very few sensors, such as Situation (4), these sensors are likely to detect clutters. Therefore, the PSA-CBMeMBer filter has a small Q r C , and the target cannot be estimated correctly because the existence probability is small. However, this problem has a low probability of occurrence and little influence on the performance of the PSA-CBMeMBer filter. EXPERIMENT 2 Considering that six targets with different initial positions move in the surveillance region [−8000, 8000] × [−8000, 8000] m 2 . The setting of targets is shown in Table 5.
The target trajectories, measurement and estimated trajectory of targets are shown in Fig. 12. The probability density of birth multi-Bernoulli RFS is described as Eq. (108).
(108) VOLUME 8, 2020 Here, r ,k = 0.1, i = 1, . . . , 5, and Here, P In this experiment, the PSA-CBMeMBer filter is compared with the IC-PHD filter [7], the Iterated Corrector (IC-CBMeMber) filter, the IIC-GM-PHD filter [26], the PM-PHD filter [27], the CM-PM-PHD filter [28], and the TS-PM-PHD filter [29]. Different parameters are used to  verify the effectiveness of the PSA-CBMeMBer filter, such as the observation noise, clutter rate and probability of detection. The setting of these parameters is given in Table 6. The simulation results of different observation noises, clutter rates and probabilities of detection are shown in Fig. 13, Fig. 14 and Fig. 15, respectively. The running time of different filters is given in Table 7.
In Fig. 13 and Fig. 14, it can be seen that OSPA distances and MSEs of the IC-PHD filter, IC-CBMeMber filter, the IIC-GM-PHD filter, the PM-PHD filter, the CM-PM-PHD filter and the TS-PM-PHD filter become worse when observation noise and clutter rate become larger. These two parameters have little effect on the PSA-CBMeMBer filter. It is mainly because the PSA-CBMeMBer filter handles the filter results of sensors. The negative effects of observation noise and clutter rate have been weakened by the filters, such as the CBMeMBer filter. Fig. 15 shows that the IC-PHD filter and the IC-CBMeMBer filer based on iterated corrector are affected seriously by the probability of detection. The stated estimation of the IC-PHD filter is improved effectively by the IIC-GM-PHD filter, but it still has a large bias of cardinality estimation. In Fig. 15(b), it can be observed that the PM-PHD  filter and the CM-PM-PHD filter perform better than the iterated corrector filters in cardinality estimation. However, the state estimation of the PM-PHD filter becomes worse because of false alarms. Both the PSA-CBMeMBer filter and the TS-PM-PHD filter perform well when the probability of detection changes. Compared with the TS-PM-PHD filter, the PSA-CBMeMBer filter has more accurate cardinality estimation.
In Table 7, the running time of the CM-PM-PHD filter and the TS-PM-PHD filter are larger than the ones of others. Although the PSA-CBMeMBer filter includes three association processes, the running time of the PSA-CBMeMBer filter is similar to those of the IC-PHD filter, the IC-CBMeMBer filter and the IIC-GM-PHD filter. To better analyze the running time, the computation complexities of these filters are given below.
The  Considering that two manoeuvring targets with different initial positions move in the surveillance region and they cross at (0, 0) (m). The setting of targets is shown in Table 8. The target trajectories, measurement and estimated trajectory of targets are shown in Fig. 16 Here, ω is the turning rate, and ω = ±3 rad/s. In Fig. 17(a), there is a slight cardinality underestimation of the IC-PHD filter and the IC-CBMeMBer filter, which is  caused by miss-detection. Moreover, the IIC-GM-PHD can reduce the impact of miss-detection on the IC-PHD filter. For the PSA-CBMeMBer filter, its parameter set estimation method ensures that targets undetected by some sensors can still be estimated correctly. From Fig. 17(b) and Fig. 17(c), it can be seen that the PM-PHD filter has accurate cardinality estimation, but the accuracy of its state estimation is lower than that of the IC-PHD filter. This problem is solved by the CM-PM-PHD filter and the TS-PM-PHD filter. Thus, these two filters perform well in both state estimation and cardinality estimation. Furthermore, when two targets cross, all filters can estimate them correctly.
Furthermore, when two targets cross, all filters can estimate them correctly. This is because all these filters are based on the PHD filter and the CBMeMBer filter. From the measurement updating processes of these two single sensor filters, it can be seen that one target can be estimated as long as there is a measurement of the target. In the measurement updating process of the PHD filter, measurements are used to update the predicted probability hypothesis density. If measurements are generated by targets, the updated probability hypothesis density has some peaks at the states of targets. Multiple peaks will be superimposed into one peak when targets close or cross. In the state extraction, these targets are estimated at the peak and have the same states. For the CBMeMBer filter, each measurement is used to all predicted parameter sets in the measurement updating process. One measurement only generates one parameter set. If the measurement is generated by targets, the parameter set generated by it will have a large existence probability. Thus, targets can be estimated by the CBMeMBer filter even they close or cross. The IC-PHD filter consists of multiple PHD filters, it has the same measurement updating process as the PHD filter. The IIC-GM-PHD filter, the PM-PHD filter, the CM-PM-PHD filter and the TS-PM-PHD filter only improve the method of calculating the cardinality, they are still essentially an IC-PHD filter. The IC-CBMeMBer filter is similar to the IC-PHD filter. In the PSA-CBMeMBer filter, all the parameter sets associated and fused are obtained by multiple CBMeMBer filters. Thus, these multi-sensor filters can accurately estimate targets in this tracking scenario. EXPERIMENT 4 Considering that one target move in the surveillance region [0, 350] × [0, 350] km 2 . There are four sensors located at (0, 0), (0, 50), (50, 0), (50, 50) (km), respectively. The maximum and minimum distance between the target and each  sensor is given in Table 9. The target trajectory, measurement and estimated trajectory of targets, and sensor positions are shown in Fig. 18. The estimated number of targets, OSPA distance and MSE are shown in Figs 19(a)-(c), respectively.
The measurement obtained by each sensor is given in polar coordinates as follow Here, the covariance matrix of observation noise v ,k is R ,k = diag δ 2 ,r , δ 2 ,θ T , = 1, 2, 3, 4. δ ,r = 100m and δ ,θ = π 180 rad are the standard deviation of range and azimuth, respectively. Fig. 19 shows that the number of targets is significantly underestimated by the IC-PHD filter and the IC-CBMeMBer filter when the target is far from sensors. It is because that the tracking results mainly depend on the measurements of the last updating sensor in the iterated corrector structure. When the target is farther away from sensors, measurements received by sensors are less accurate. Therefore, the tracking performance is not good with only one sensor's measurement. The other multi-sensor filters use the measurements of all sensors when estimating the number of targets. Therefore, their cardinality estimations are much better than IC-PHD filter. Note that the OSPA distance of the PM-PHD filter is small Fig. 19(b), but its OSPA distance is higher than that of the IC-PHD filter in Fig. 17(b). In multi-target tracking, some targets are missed by the IC-PHD filter. To estimate the number of targets correctly, the PM-PHD filter increases the cardinality of the detected targets. Some false targets are estimated which have the same states as the detected targets. Thus, the OSPA distance of the PM-PHD filter is large. Since there is only one target in this tracking scenario, the PM-PHD filter can ensure that the number of targets is one. Then, the state of the target is estimated from the predicted state. In this experiment, the performance of the PSA-CBMeMBer filter is slightly worse than that of the PM-PHD filer and its improved versions. The reason is that the coefficient of measurement likelihood is calculated differently in the CBMeMBer filter and the PHD filter. This coefficient plays an important role in determining whether the measurement is generated by the target or clutter. For one measurement, its coefficient obtained by the CBMeM-Ber filter is smaller than that obtained by the PHD filter. Therefore, the measurement is easier to be determined as the measurement of clutter by the CBMeMBer filter. This situation will become worse as the accuracy of measurement decreases. Thus, when the target is far from sensors and the accuracy of the measurements decrease, the performance of IC-CBMeMBer filter is worse than that of the IC-PHD filter. Meanwhile, this problem also affects the tracking performance of the PSA-CBMeMBer filter.

VIII. CONCLUSION
This paper presents a tracks association and fusion methods based on the CBMeMBer filter. The tracks can be distinguished and associated by the tracks association method with three association processes. The tracks corresponding to the same target can be grouped into one partition by the adaptive association thresholds under the assumptions of the linear Gaussian mixture model. Considering the relationships among different types of tracks, the cardinality overestimation problem of the CBMeMBer can be avoided effectively. Simulation results validate that the PSA-CBMeMBer filter is insensitive to the observation noise and the clutter rate, and it is little affected by the probability of detection.
This paper assumes that the state of birth targets is known as a priori. However, the state of new born targets is unknown in real tracking scenarios. In the second association process, only the information of measurements is used in the measurement update parameter sets association method. In our future works, we will extend this approach to deal with the problem of target birth in the multi-sensor tracking system. The solution for Eq. (65) is given as follows. By Eqs. (32) and (33), the state mean vector and the covariance matrix of N z ,B are given by m z ,B == m k|k−1 + K k|k−1 z z ,B − H k m k|k−1 (122) P z ,B = I − K k|k−1 H k P k|k−1 (123) HONGBING JI (Senior Member, IEEE) graduated from the Northern West Telecommunications Engineering College (now Xidian University), Xi'an, China. He received the B.S. degree in radar engineering, the M.S. degree in circuit, signals, and systems, and the Ph.D. degree in signal and information processing from Xidian University, in 1983, 1989, and 1999, respectively. Since 1989, he has been with the School of Electronic Engineering, Xidian University. He is currently a Professor and an Advisor for Ph.D. students. His primary areas of research have been radar signal processing, automatic targets recognition, multisensor information fusion, and target tracking. He is currently an Assistant Professor with the School of Electronic Engineering, Xidian University. His research interests include machine learning, image processing, and target recognition.
ZHENZHEN SU was born in Yulin, China, in 1993. He received the B.S. degree in applied mathematics from Xidian University, in 2015, where he is currently pursuing the Ph.D. degree. His research interests include target tracking, visual tracking, and graphical models.
PENG WANG was born in Lianyungang, China, in 1994. He received the B.S. degree in applied mathematics from Xidian University, in 2017, where he is currently pursuing the Ph.D. degree. His research interests include signal processing, target tracking, and information fusion.