Labeled RFS-Based Track-Before-Detect for Multiple Maneuvering Targets in the Infrared Focal Plane Array

The problem of jointly detecting and tracking multiple targets from the raw observations of an infrared focal plane array is a challenging task, especially for the case with uncertain target dynamics. In this paper a multi-model labeled multi-Bernoulli (MM-LMB) track-before-detect method is proposed within the labeled random finite sets (RFS) framework. The proposed track-before-detect method consists of two parts—MM-LMB filter and MM-LMB smoother. For the MM-LMB filter, original LMB filter is applied to track-before-detect based on target and measurement models, and is integrated with the interacting multiple models (IMM) approach to accommodate the uncertainty of target dynamics. For the MM-LMB smoother, taking advantage of the track labels and posterior model transition probability, the single-model single-target smoother is extended to a multi-model multi-target smoother. A Sequential Monte Carlo approach is also presented to implement the proposed method. Simulation results show the proposed method can effectively achieve tracking continuity for multiple maneuvering targets. In addition, compared with the forward filtering alone, our method is more robust due to its combination of forward filtering and backward smoothing.


Introduction
The problem of jointly detecting and tracking multiple maneuvering targets in an infrared focal plane array is very challenging and has received great attention in the last several years [1][2][3][4][5][6][7]. In many applications, the estimation is often performed on point measurements after threshold segmentation [8][9][10][11]. However, in situations where the signal-to-noise ratio (SNR) of the infrared sensor is low, threshold segmentation of the sensor output may cause false alarms and missing targets, as the noise level is high enough to generate detections [12,13], so it is necessary to make use of all information in the images to improve the detection and tracking performance. Fortunately, track-before-detect (TBD), or tracking without threshold segmentation is an effective method as shown in [14][15][16]. In addition, in many applications, for example infrared search and track (IRST), precise guidance, and space situation awareness (SSA), the tracking method is required to track all relevant targets which typically exhibit different motion models [6,[17][18][19]. This paper investigates the problem of jointly estimating the number of maneuvering targets and their states from infrared image observations.
A recent approach of tracking is to represent the multi-target state as a random finite set (RFS). Based on the RFS framework, the probability hypothesis density (PHD), Cardinalized PHD (CPHD) Sensors 2015, 15, 30839-30855

Notation
In this paper, small letters (e.g., x) are used to denote single-target states and capital letters (e.g., X) are used to denote multi-target states. Labeled target states are indicated by boldface letters (e.g., x, X). In addition, spaces are represented by blackboard bold letters (e.g., X denotes the state space). Image observation at time k is denoted by z k . The inner product is denoted by x f , gy fi ş f pxq g pxq dx. The multi-target exponential of a real valued function h raised to a set X is defined as h X fi ś xPX h pxq, where h φ " 1 by convention. The generalized Kronecker delta function and inclusion function are denoted by δ Y p¨q and 1 Y p¨q, respectively.

Labeled RFS
In order to jointly estimate the targets' states and their individual tracks, labeled RFS was introduced based on the RFS framework [25]. A label P L " tα i : i P Nu is appended to the state x P X to enable the estimation of a target track. L is a discrete space, whose elements are distinct and the space N denotes the set of positive integers. Hence, a labeled single-target state and a labeled multi-target state can be described by x " px, q and X " !
x p1q ,¨¨¨, x pnq ) , respectively. Let L : XˆL Ñ L be the projection L ppx, qq " , the set of track labels of the labeled RFS X is obtained by L pXq " tL pxq : x P Xu. The labels are required to be distinct. It means that the cardinalities of the set of labels and the set of state vectors are identical. It is mathematically ensured using the distinct label indicator ∆ pXq " δ |X| p|L pXq|q.

Labeled Multi-Bernoulli RFS
The single target density can be modeled by Bernoulli RFS. Conditional on X k´1 " ∅, the target can re-enter or appear with probability p b,k|k´1 and occupy kinematic state x k with probability density b k|k´1 px k q, or remain absent or disappeared with probability 1´p b,k|k´1 . It can be described by the Bernoulli RFS as following [32]: In addition, conditional on X k´1 " tx k´1 u , the target can survive and acquire a new state x k with probability density p S,k|k´1 px k´1 q f k|k´1 px k |x k´1 q, or disappear with probability 1´p S,k|k´1 px k´1 q. It can be described by the following Bernoulli RFS [12]: A multi-Bernoulli RFS X can be regarded as a union of independent Bernoulli RFSs X piq with existence probability r piq and probability density p piq , i.e., X " Y M i"1 X piq . Then, parameter set is used to represent a multi-Bernoulli RFS.
Like a multi-Bernoulli RFS, a labeled multi-Bernoulli RFS with state space X and label space L can be described by the parameter set π " ! r p q , p p q ) PL . The components are assumed to be statistically independent [27,35]. A label " pk, iq is assigned to each target which is a pair of the time of birth k and a label index i P N. Then, the label space for new targets born at time k is denoted as L k , and the new target born at time k has state x P XˆL k . The label space for all targets at time k can be denoted as L 0:k , which is constructed recursively by L 0:k " L 0:k´1 Y L k . In addition, the LMB RFS can also be represented in the form of GLMB RFS: Sensors 2015, 15,[30839][30840][30841][30842][30843][30844][30845][30846][30847][30848][30849][30850][30851][30852][30853][30854][30855] where the weights and the spatial distributions are given as follows [27]:

Target and Measurement Models
In this paper, the observation is a two-dimensional image generated by infrared focal plane array. The pixel i of image can be indexed by the row and column coordinates, i.e., i " pa, bq. The image is regarded as being made up of the sum of target signals and sensor noise signals. For long-range surveillance applications, the target is close to a point source, referred to as point target [36]. Thus, the sampled target spatial signature can be well modeled by the following 2D Gaussian shape [12,37]: where Σ is blurring factor, px a , x b q is the position of the state x, I is the target intensity. Each cell corresponds to a rectangular region of dimensions ∆ aˆ∆b . h i pxq is the contribution to pixel i from the state x, which depends on the blurring factor, target position and target intensity. A target with state x illuminates a set of pixels denoted by T pxq. The T pxq is regarded as the square region whose center is closest to the position of x. Let the measurement likelihood in pixel i in the presence of a target with state x be denoted by φ piq´zpiq , x¯, and the likelihood under the hypothesis of no targets be ϕ piq´zpiq¯, as Equation (7). More details on the measurement likelihood can be found in [38]: Then, the point target model is incorporated into the likelihood function. The measurement likelihood is given as follows: φ piq´zpiq¯" N´z piq ; 0, σ¯ (8) ϕ piq´zpiq , x¯" N´z piq ; h i pxq , σ¯ (9) Conditioned on the multi-target state X, the multi-target likelihood of z is the product over all cells [38], i.e.: where: piq`zpiq˘ ( 12) Note that the targets are not closely space targets in this paper. We assume they are not very close to each other, so the separable likelihood model is reasonable in this application.

Multi-Model LMB Filter for TBD
In our method, MM-LMB TBD consists of the MM-LMB filter and MM-LMB smoother for TBD. MM-LMB filter means that the target density is propagated forward from time k to time l ą k. MM-LMB smoother means that the target density is propagated backward from time l to time k. This subsection gives the derivation of the forward filter. To derive the MM-LMB filter for TBD, the key issue is to construct the LMB parameter set. For multi-model filter, augmented state px, τq is typically used to represent the target's state and the according motion model τ P T. T " t1, 2,¨¨¨, M τ u denotes the discrete space of all possible motion models. In the context of labeled RFS, the state vector is further augmented by the track label , as px, , τq. Thus, multi-target state is consequently given by: X " tx 1 ,¨¨¨, x n u " tpx 1 , 1 , τ 1 q ,¨¨¨, px n , n , τ n qu (13) Let h st " p k|k´1 tτ k " t |τ k´1 " s u denotes the probability that a track switches from model s to model t, where s, t P T. H " rh st s denotes the model transition probability matrix, and ř M τ t"1 h st " 1. Each track is represented by its existence probability r p q and its joint spatial distribution p p q px, τq " p p q px|τq p p q pτq. p p q pτq denotes the probability of track with model τ and p p q p x| τq denotes the spatial distribution of track conditioned on model τ.
In order to tackle the problem of maneuvering target TBD, IMM approach is integrated with the LMB filter as MM-LMB filter. Thus, in addition to the prediction stage and update stage, mixing stage is introduced. The recursion of MM-LMB filter for TBD is given as follows: Step 1. Mixing: If at time k´1 the posterior density is assumed to be an LMB RFS on the augmented space, which is given by the parameter set as: Then, the parameter set after mixing can be expressed by: For the IMM approach, the motion model transition is assumed to be independent of the target's state transition. It is only decided by motion model transition probability. Thus, r p p q k´1´x p q k´1 , τ p q k¯c an be calculated by:

16)
Step 2. Prediction: Then, the predicted LMB RFS of the multi-model LMB filter is given by the parameter set: where: r p q P,k|k´1 " r p q k´1 η p q The terms p p q S p¨q and f p q px |¨, τq denote the survival probability and single-target transition density of track at model τ, respectively. In most applications, the survival probability of the tracks can be assumed to be independent of the current motion model, i.e., p k´1 px|τq denote the model probability and conditional spatial distribution of track after mixing at time k´1. As shown in Equation (17), the LMB parameter set for the predicted multi-target density π k|k´1 is formed by the union of LMB parameter sets for the survival targets and birth targets. The first term in Equation (17) represents survival targets and the second one denotes birth targets. B is the label space of birth targets, and L k|k´1 " L k´1 Y B.
Step 3. Update from image observation: We rewrite the predicted multi-target density as a typical LMB RFS, which is given by the parameter set: Reference [29] derived the update procedure for single-model TBD. For the case with multiple motion models, the existence probability r p q k and spatial density p p q k px, τq can be updated by Equations (22) and (23), respectively: where γ z p¨q has been given in Section 3.1, and: where p p q k|k´1 pτq and p p q k|k´1 px|τq represent the model probability and conditional spatial distribution of track with model τ. In fact, the Equation (22) can be interpreted as the union of two cases: target appears and target disappears. r p q k|k´1 η z p q denotes the probability when target appears. 1´r p q k|k´1 denotes the probability when target disappears.

Multi-Model LMB Smoother for TBD
In the backward smoothing, the smoothed target density is propagated backward from time l to time k ă l. In essence, we want to optimize the history target density by the measurement data up to time l ą k. Taking advantage of the target labels, the multi-target backward smoothing can be realized with an approach similar to that of single-target backward smoothing. Multiple targets are Sensors 2015, 15,[30839][30840][30841][30842][30843][30844][30845][30846][30847][30848][30849][30850][30851][30852][30853][30854][30855] independent, and can be smoothed respectively. More details on single-target backward smoothing can be found in [12]. The recursion is initialized with l " k. The smoothed LMB density of track from time l to k is denoted by π k|l " !´r . Then the smoothed LMB density can be calculated by: (27) where: x 1 , τ 1˘i s the augmented states of smoothing from l to k, f`x 1 , τ 1ˇx , τ˘is single target transition density at time k, given target state x and motion model τ. The model transition is independent of the state transition, i.e., f` Like the analysis above, the motion model transition probability has an important influence on target state estimation. Backward smoothing makes it possible to get more accurate posterior model transition probability. The posterior model transition probability of track from k´1 to k, i.e., r h p q τ k´1 τ k , is calculated as Equation (32). The posterior model transition probability, rather than the prior one, is used to calculate smoothed target density:

Sequential Monte Carlo Implementation
In this subsection, we use the SMC approach to implement the MM-LMB TBD. The single target density at time k´1 is given by a labeled Bernoulli RFS with existence probability r p q k´1 and spatial density p The spatial density is approximated using a set of weighted particles k´1 is the number of particles.
Notice that the particles consist of the state and model information with associated weights.
Mixing and Prediction: The model samplings for survival targets and birth targets are performed based on the proposals α k p¨q and β k p¨q. Weights are calculated based on the model transition probability and birth model, respectively. Thus, the model values and corresponding weights can be obtained by: where, h k|k´1 p¨q and θ k p¨q are the model transition probability of survival targets and the model distribution of birth targets. J B denotes the number of particles of the birth target. Note that the predicted motion models are sampled from discrete space T.
Then, the particles x p ,jq k|k´1 and weights ω p ,jq k|k´1 can be generated as follows. The particles and weights for survival targets are drawn conditioned on model τ k , and those for birth targets are obtained based on birth model: where, q k p¨q and s k p¨q are the proposals for the survival target and the birth target.
Update: Suppose the predicted density is given as π k|k´1 " ! r  (22) and (23) for update can be calculated by: where the J p q k|k´1 covers all particles of track with different motion models. Then, the updated existence probability, i.e., r p q k , can be obtained as Equation (22), and the weight for each particle, i.e., ω p ,jq k , can be calculated by: Backward Smoothing: The backward smoothing uses the target density at time l to smooth the particle approximation of the target density and existence probability at time k ă l.
Suppose the smoothed density from l to k is given as π k|l " ! r to k can be calculated by: The smoothed multi-target LMB density can be calculated by: where: r ω In addition, the resampling, pruning and state extraction are same as the SMC-CBMeMBer [23]. For any two tracks, if the extracted target states are within a given distance threshold Th m the two corresponding LMB components are merged. If the existence probability of th track is below a threshold Th p , the corresponding LMB component is discarded. A target is declared present if the estimated existence probability is greater than threshold Th T , otherwise no target is declared. At each time, a maximum of J max and minimum of J min particles per-hypothesized track are imposed so that the number of particles representing each hypothesized track is proportional to its existence probability after resampling in the update step.

Scenario Description
A scenario containing multiple maneuvering targets is used to evaluate the performance of proposed method. The scenario is 30 frames long with period ∆ T " 1. There are three point targets in this scenario. The point spread function of the simulated sensor has a blurring factor Σ " 1 and sensor noise variance is set to σ " 1. The image observation consists of 512ˆ512 pixels. Each pixel representing one unit of physical distance: ∆ a " ∆ b " 1. The T pxq is 4ˆ4 square region whose center is closest to the position of x.
Each of the targets may move at a constant velocity, performing a coordinated turn or under constant acceleration in the surveillance region. Therefore, the motion model set for this example can be composed of a constant velocity (CV) model, a coordinated turn (CT) model and a constant acceleration (CA) model. More details on the motion models can be found in [39]. Let τ " 1 denote the CV model, τ " 2 denote the CT model and τ " 3 denote the CA model. The turn rate of coordinated turn model is set to 1 rad/s, and the acceleration is set to 2.3 pixel/s. The target tracks are shown in Figure 1. " " denotes the locations at which targets are born and "˝" denotes the locations at which targets die. A typical observation with SNR " 6.8 (I " 15) is shown in Figure 2. The SNR can be calculated by [34]: Sensors 2015, 15 12 turn model is set to 1 rad/s, and the acceleration is set to 2.3 pixel/s. The target tracks are shown in Figure 1. "  " denotes the locations at which targets are born and "  " denotes the locations at which targets die. A typical observation with 6.8 SNR  ( 15 I  ) is shown in Figure 2. The SNR can be calculated by [34]:  turn model is set to 1 rad/s, and the acceleration is set to 2.3 pixel/s. The target tracks are shown in Figure 1. "  " denotes the locations at which targets are born and "  " denotes the locations at which targets die. A typical observation with 6.8 SNR  ( 15 I  ) is shown in Figure 2. The SNR can be calculated by [34]:  For the IMM approach, the motion models are initialized with probabilities p pτ " 1q " p pτ " 2q " p pτ " 3q " 1{3 and the model transition matrix is set to: H " The performance of the algorithm is evaluated using the Optimal Sub-Pattern Assignment (OSPA) metric [40], which is recently developed and defined as: where X " tx 1 ,¨¨¨, x m u and Y " ty 1 ,¨¨¨, y m u are arbitrary finite subsets, 1 ď p ď 8, c ą 0. We use the parameters c " 10 and p " 1 in this paper. The results are divided into two parts: Firstly, multiple maneuvering targets are detected and tracked by TBD-based algorithms, including CV-LMB filter, MM-PHD filter [34,41] and MM-LMB filter. This is used to illustrate the effectiveness of MM-LMB filter for multiple maneuvering targets TBD. Additionally, the performance of MM-LMB TBD with forward-backward combination is compared with the forward MM-LMB filter. It demonstrates the optimization ability of backward smoothing. These TBD methods are implemented using SMC approach. The parameters used in SMC application are listed in Table 1. The results are divided into two parts: Firstly, multiple maneuvering targets are detected and tracked by TBD-based algorithms, including CV-LMB filter, MM-PHD filter [34,41] and MM-LMB filter. This is used to illustrate the effectiveness of MM-LMB filter for multiple maneuvering targets TBD. Additionally, the performance of MM-LMB TBD with forward-backward combination is compared with the forward MM-LMB filter. It demonstrates the optimization ability of backward smoothing. These TBD methods are implemented using SMC approach. The parameters used in SMC application are listed in Table 1.

Multiple Maneuvering Target TBD Experiment
The results obtained from CV-LMB filter, MM-PHD filter and MM-LMB filter for a single run at 6.8 SNR  is shown in Figure 3. As indicated by this figure, the tracking performance is quite different. When the motion model switch occurs or the selected model is not consistent with the true model, the CV-LMB filter fails to estimate the targets. The estimation error of MM-PHD filter is also large, because of the inaccuracy in estimating the number of targets. However, the MM-LMB filter can track all targets from birth to end in the above defined scenario.  The true motion model and estimated model probability in MM-LMB filter at different times is depicted by Figures 4 and 5 respectively. As expected, the MM-LMB filter can adaptively capture model transition. For example in Target 3, the probabilities of CV, CT and CA at k " 8 are 0.53, 0.01, and 0.46, respectively. However, at k " 9, the probabilities of CV, CT and CA become 0.005, 0.005 and 0.99. It's obvious that motion model transition is achieved successfully in MM-LMB filter.
Sensors 2015, 15 14 The true motion model and estimated model probability in MM-LMB filter at different times is depicted by Figures 4 and 5, respectively. As expected, the MM-LMB filter can adaptively capture model transition. For example in Target 3, the probabilities of CV, CT and CA at 8 k  are 0.53, 0.01, and 0.46, respectively. However, at 9 k  , the probabilities of CV, CT and CA become 0.005, 0.005 and 0.99. It's obvious that motion model transition is achieved successfully in MM-LMB filter.  To provide a performance comparison of the methods in sense of statistical evaluation, the average cardinality and average OPSA distance over 50 Monte Carlo runs at different SNR are shown in Figures 6 and 7. Figures 6a and 7a give the results at 6.8 SNR  . Figures 6b and 7b give the results at 5 SNR  , respectively. Sensors 2015, 15 14 The true motion model and estimated model probability in MM-LMB filter at different times is depicted by Figures 4 and 5, respectively. As expected, the MM-LMB filter can adaptively capture model transition. For example in Target 3, the probabilities of CV, CT and CA at 8 k  are 0.53, 0.01, and 0.46, respectively. However, at 9 k  , the probabilities of CV, CT and CA become 0.005, 0.005 and 0.99. It's obvious that motion model transition is achieved successfully in MM-LMB filter.  To provide a performance comparison of the methods in sense of statistical evaluation, the average cardinality and average OPSA distance over 50 Monte Carlo runs at different SNR are shown in Figures 6 and 7. Figures 6a and 7a give the results at 6.8 SNR  . Figures 6b and 7b give the results at 5 SNR  , respectively. To provide a performance comparison of the methods in sense of statistical evaluation, the average cardinality and average OPSA distance over 50 Monte Carlo runs at different SNR are shown in Figures 6 and 7. Figures 6a and 7a give the results at SNR " 6.8. Figures 6b and 7b give the results at SNR " 5, respectively.  From Figure 6, we can see that MM-LMB filter converges to the ground truth. However, the CV-LMB filter and MM-PHD filter produce bias in estimating the target number. MM-LMB filter can effectively capture the model switching property of the maneuvering targets, so its performance is significantly better than that of the other two methods. The average OSPA distances derived by the three methods are shown in Figure 7. From this figure, we can see that CV-LMB filter and MM-PHD filter perform significantly worse than the MM-LMB filter in the term of target localization and detection. In addition, the MM-LMB filter is more effective for dim targets. As shown in Figures 6 and 7, the estimation error may increase as SNR decreases. However, MM-LMB filter still outperforms the other two methods.  From Figure 6, we can see that MM-LMB filter converges to the ground truth. However, the CV-LMB filter and MM-PHD filter produce bias in estimating the target number. MM-LMB filter can effectively capture the model switching property of the maneuvering targets, so its performance is significantly better than that of the other two methods. The average OSPA distances derived by the three methods are shown in Figure 7. From this figure, we can see that CV-LMB filter and MM-PHD filter perform significantly worse than the MM-LMB filter in the term of target localization and detection. In addition, the MM-LMB filter is more effective for dim targets. As shown in Figures 6 and 7, the estimation error may increase as SNR decreases. However, MM-LMB filter still outperforms the other two methods. From Figure 6, we can see that MM-LMB filter converges to the ground truth. However, the CV-LMB filter and MM-PHD filter produce bias in estimating the target number. MM-LMB filter can effectively capture the model switching property of the maneuvering targets, so its performance is significantly better than that of the other two methods. The average OSPA distances derived by the three methods are shown in Figure 7. From this figure, we can see that CV-LMB filter and MM-PHD filter perform significantly worse than the MM-LMB filter in the term of target localization and detection. In addition, the MM-LMB filter is more effective for dim targets. As shown in Figures 6  and 7 the estimation error may increase as SNR decreases. However, MM-LMB filter still outperforms the other two methods.

Backward Smoothing Experiment
To provide a performance comparison of forward method and forward-backward method, an example simulation is given with target intensity fluctuation. As shown in Figure 8, target intensity degrades significantly at some time, which may lead to significant estimate errors. Besides, 1-lag smoother is adopted in this experiment, it means that the target density at time k is used to smooth the particle approximation of the target density and existence probability at time k´1.

Backward Smoothing Experiment
To provide a performance comparison of forward method and forward-backward method, an example simulation is given with target intensity fluctuation. As shown in Figure 8, target intensity degrades significantly at some time, which may lead to significant estimate errors. Besides, 1-lag smoother is adopted in this experiment, it means that the target density at time k is used to smooth the particle approximation of the target density and existence probability at time 1 k  . In this experiment, the posterior model transition probabilities at some representative times are listed in Table 2.   , s t represents the model transition from s to t at time k . st h and st h  represent the prior and posterior model transition probabilities, respectively. As shown in Table 2, posterior model transition probabilities conform to actual motion models better, especially when motion model switches. For example, Target 3 switches from CV to CT at 9 k  . The prior transition probability is 0.1, but the posterior one is about 0.8. In addition, Target 3 switches from CT to CV at 12 k  . The prior transition probability is 0.15, but the posterior one is about 0.85. It proves that the model transition probabilities have been optimized based on posterior target density in smoother.  In this experiment, the posterior model transition probabilities at some representative times are listed in Table 2. ps, tq represents the model transition from s to t at time k. h st and r h st represent the prior and posterior model transition probabilities, respectively. As shown in Table 2, posterior model transition probabilities conform to actual motion models better, especially when motion model switches. For example, Target 3 switches from CV to CT at k " 9. The prior transition probability is 0.1, but the posterior one is about 0.8. In addition, Target 3 switches from CT to CV at k " 12. The prior transition probability is 0.15, but the posterior one is about 0.85. It proves that the model transition probabilities have been optimized based on posterior target density in smoother.  Figure 9 indicates that our method provides accurate cardinality estimation. In Figure 10, the OSPA peaks of filter appear at time k " t5, 10, 15, 20, 25u, which are related to target intensity fluctuation. However, the result of our method is more stable. This phenomenon indicates that our method is more robust, especially when the target intensity may fluctuate. It's because that the smoother can exploit current observation to optimize the history result. runs. Figure 9 indicates that our method provides accurate cardinality estimation. In Figure 10, the OSPA peaks of filter appear at time   5,10,15, 20, 25 k  , which are related to target intensity fluctuation. However, the result of our method is more stable. This phenomenon indicates that our method is more robust, especially when the target intensity may fluctuate. It's because that the smoother can exploit current observation to optimize the history result.

Conclusions
In this paper, we have proposed a labeled RFS-based method, MM-LMB TBD, to detect and track multiple maneuvering targets from raw image observation of infrared focal plane arrays. A forward filter and backward smoother are jointly used to estimate the multi-target states. Simulation results show that the MM-LMB TBD method is capable of tracking maneuvering targets using raw image data. It has better adaptation to maneuvers and provides an overall lower estimate errors compared  OSPA peaks of filter appear at time   5,10,15, 20, 25 k  , which are related to target intensity fluctuation. However, the result of our method is more stable. This phenomenon indicates that our method is more robust, especially when the target intensity may fluctuate. It's because that the smoother can exploit current observation to optimize the history result.

Conclusions
In this paper, we have proposed a labeled RFS-based method, MM-LMB TBD, to detect and track multiple maneuvering targets from raw image observation of infrared focal plane arrays. A forward filter and backward smoother are jointly used to estimate the multi-target states. Simulation results show that the MM-LMB TBD method is capable of tracking maneuvering targets using raw image data. It has better adaptation to maneuvers and provides an overall lower estimate errors compared

Conclusions
In this paper, we have proposed a labeled RFS-based method, MM-LMB TBD, to detect and track multiple maneuvering targets from raw image observation of infrared focal plane arrays. A forward filter and backward smoother are jointly used to estimate the multi-target states. Simulation results show that the MM-LMB TBD method is capable of tracking maneuvering targets using raw image data. It has better adaptation to maneuvers and provides an overall lower estimate errors compared with the CV-LMB filter and MM-PHD filter. In addition, the smoother of MM-LMB TBD can backwardly optimize the history target state estimation by the current measurement. Thus, MM-LMB TBD is more robust against the target intensity fluctuation problem, which is common in practical applications.
Author Contributions: Miao Li developed the algorithm, conducted the experiments and participated in writing the paper. Jun Li contributed to the theoretical analysis and algorithm development. Yiyu Zhou planned and supervised the research and analyzed the results.

Conflicts of Interest:
The authors declare no conflict of interest.