repirst

repirst - DETECTION ALGORITHMS AND TRACK BEFORE DETECT...

Info iconThis preview shows page 1. Sign up to view the full content.

View Full Document Right Arrow Icon
This is the end of the preview. Sign up to access the rest of the document.

Unformatted text preview: DETECTION ALGORITHMS AND TRACK BEFORE DETECT ARCHITECTURE BASED ON NONLINEAR FILTERING FOR INFRARED SEARCH AND TRACK SYSTEMS TECHNICAL REPORT CAMS 98.9.1 Skirmantas Kligys, Boris Rozovsky, and Alexander Tartakovsky CENTER FOR APPLIED MATHEMATICAL SCIENCES University of Southern California Los Angeles, CA 90089 1113 SEPTEMBER 1998 Approved for public release; distribution unlimited DETECTION CENTER 0 DETECTION ALGORITHMS AND TBD ARCHITECTURE FOR IRST Contents 1 1. Introduction 2 2. Problem Formulation and Background 3 2.1. TbD Methods: Merits and Drawbacks 4 3. Signal and Observation Modeling 6 4. Clutter Suppression 7 4.1. Nonparametric Method 8 4.2. Semiparametric Filtering 8 5. Optimal Nonlinear Filtering for TbD 9 5.1. Modeling for TbD 9 5.2. Optimal Spatio-Temporal Nonlinear Filtering: The Basic Algorithm 10 5.3. Multi-Dimensional Spatio-Temporal Matched Filter as a Special Case of Nonlinear Filter 12 6. Track Appearance Disappearance Detection 14 6.1. Preliminaries 14 6.2. Problem Formulation and Model Description 16 6.3. Detection Algorithm 1: Fixed Sliding Window 17 6.4. Detection Algorithm 2: Fully Sequential Procedure 19 6.5. Algorithm 3: Joint Detection of Track Appearance and Disappearance 21 6.6. Adaptive Detection Algorithms 23 6.7. Choice of Thresholds and Performance Evaluation 24 6.8. TbD Model with Spatio-Temporal Matched Filter 29 6.9. TbD Model with Optimal Nonlinear Filter 31 7. Testing of the Developed Algorithms for IRST Data. Results of Experiments and Simulation 32 7.1. Clutter Removal: Real IRST Data 32 7.2. Example 1: TbD of a Target Based on IRST Data 33 7.3. Example 2: TbD of a Surface Skimming Missile 41 7.4. Appearance Disappearance Detection of a Skimming Missile 47 8. Conclusion. Future Research 50 9. Acknowledgement 50 References 51 2 CENTER FOR APPLIED MATHEMATICAL SCIENCES, USC 1. Introduction Cruise missiles over land and sea cluttered background are serious threats to Infrared Search and Track Systems IRST's. In general, these threats are stealth in both the infrared and radio frequency bands. That is, their thermal infrared signature and their radar cross section can be quite small. Future predicted threats, i.e. the next generation of cruise missiles, will be even more di cult to detect at a su cient range to counter. Further, low elevation trajectory objects, such as sea skimming missiles, have radar signals with large amounts of temporally and spatially correlated interference called multipath. This multipath problem remains an enormous obstacle to existing trackers. Hence, new technology is needed which will allow for the timely detection, tracking, and identi cation of such threats. IRST systems are one component of a multisensor suite which can meet the technical challenge of the timely detection track identi cation of low signal-to-noise+clutter ratio SN+CR targets. The multisensor suite should include an IRST, Radar, and a coherent laser Lidar. We envision a cueing hierarchy where the IRST can cue the Radar or the Radar cues the IRST. Once a candidate track is established the Lidar can be used to identify the target by its micro doppler signature. In this report we describe the developed computationally e cient algorithms and adaptive architecture with optimized overall performance statistical and computational for real-time reliable detection and tracking of low-observable targets in IRST systems. Despite the fact that we focus on an IRST against cruise missiles over land and sea cluttered backgrounds, the results are equally applicable to other sensors e.g., Radar, Lidar. In the research we concentrated on the three interrelated problems: 1 e cient clutter suppression; 2 development of the adaptive track-before-detect TbD architecture based on optimal nonlinear ltering ONF; 3 development of e cient algorithms for detection of a priori unknown number of targets that may appear and disappear at unknown points in time. The report is organized as follows. In Section 2 we formulate the problem, describe a structure of the system to be developed, outline popular track-before-detect methods that are in current use and suggest an alternative method, which is based on the optimal nonlinear ltering. In Section 3 we describe basic models and assumptions on signals and clutter that are used in developing of signal processing algorithms clutter suppression, track-before-detect and detection algorithms. In Section 4, two clutter removal algorithms are presented based on nonparametric and semiparametric spatio-temporal ltering. Section 5 is especially important. Here we describe an optimal nonlinear ltering technique that is used for track-before-detect of very low observable targets. Also, we show that the proposed method is a complete generalization of the multidimensional spatio-temporal matched ltering particularly, 3D matched lter. Furthermore it is shown that the spatio-temporal matched lter coincides with the developed optimal nonlinear lter when a target moves according to deterministic trajectory. In Section 6 we develop three track appearance disappearance algorithms all of which take into account requirements speci c for surveillance systems. The results of simulation and processing real IR data obtained from DETECTION ALGORITHMS AND TBD ARCHITECTURE FOR IRST 3 SPAWAR Space & Naval Warfare Systems Center, San Diego, CA are presented in Section 7. Finally in Section 8 we provide a conclusion and the plan of our future research. 2. Problem Formulation and Background The generalized block-diagram of the system under investigation is shown in Figure 1. We develop both the signal processing architecture clutter removal algorithms and TbD algorithms and track detection algorithms. Signal Processing Block Clutter Suppression (Nonparametric or Semiparametric Spatio-Temporal Filtering) IR Frames, Preprocessing Track-Before-Detect (Optimal Nonlinear Spatio-Temporal Filtering) Track Appearance/ Disappearance Detection Detections (Blips) Figure 1. Generalized block-diagram of the developed system New algorithms concerning the stages of data processing speci ed above are developed under the following realistic conditions. Cluttered background is much more intensive than both equivalent intrinsic instrumental noise of the sensor and signal intensity of the targets to be detected. This causes a necessity of practically complete suppression of a clutter. Exterior conditions of observation are characterized by an extremely high variability and prior uncertainty and may not be predicted with su cient accuracy. Prior information that is needed to develop ideal Bayes data processing algorithms is not available. Particularly, statistical models of signals and moreover exterior background as well as the models of changing target situation are extremely unreliable. Such models can be useful only as tools for performance evaluation in certain scenarios but not for development of data processing algorithms. Practical data processing algorithms should be developed on the base of sequential application of robust, adaptive and minimax methods that are invariant to the prior uncertainty. In practice the estimated parameters of targets for example, trajectory parameters should be guaranteed for any degree of prior uncertainty. Speci cally, each estimate should be supplemented with an appropriate domain of minimum size that contains the real unknown value of estimated parameter with 100 or at least 1 , " assurance. 2.1. TbD Methods: Merits and Drawbacks. The most challenging problem for an IRST system is the detection of a maneuvering target against a strong clutter background. To illustrate the importance of this task, we remark that under certain conditions a few seconds decrease in the time it takes to detect a sea surface skimming cruise missile can yield a signi cant increase in the probability of raid annihilation. The problem of detection is extremely di cult in low SN+CR when localization of the target based upon a single non-stationary image is impossible. In this case one has to align successive frames according to typical patterns of target dynamics and any results of preliminary" tracking. This approach to detection of a low SN+CR target is usually referred to as track-before-detect" TbD. Its success depends crucially on the quality of the preliminary" tracker. Thus, the development of the e cient coherent signal processing based on the TbD methodology becomes crucial point in the low-observable target detection problem. In contrast to other TbD methods, we solve this problem by applying optimal nonlinear ltering approach. We now overview several popular methods for tracking before detection that are in current use. 2.1.1. Spatio-Temporal Matched Filters and Velocity Filter Banks. Given a sequence of frames, consider the problem of detecting a small unresolved target against a clutter background. The probabilities of errors can be decreased by applying a 3,D spatio-temporal matched lter prior to detection thresholding operation. If the spatial distribution of the target and its velocity remain unchanged i.e., if targets move with known constant speed along a line on the plane and the noise and the clutter are Gaussian processes, the 3,D lter is the optimal method of detection, see 39 . It is easily shown that the same result is true for more dimensions, i.e. a m + 1,D matched lter is an optimal method of processing under aforementioned conditions m is the spatial dimensionality. Also, under these conditions the target" component of the multi-dimensional matched lter separates into spatial and temporal components. The temporal component is usually called a velocity lter. Typically the target velocity is unknown and hence the single velocity lter cannot be used. This problem is usually overcome by hypothesizing velocities and implementing a velocity lter bank see, e.g., 44, 49, 51 . One of the main drawbacks of spatio-temporal matched lters and other modi cations such as banks of assumed velocity lters is that they are not able to work with maneuvering targets. Performance of the algorithms is substantially degraded in the presence of velocity mismatch or in event of target maneuver. 2.1.2. Dynamic Programming Methods. The Dynamic Programming methods for TbD showed a big advantage over the conventional MHT method and over the 3,D matched lter with velocity mismatch 1, 2, 7, 21, 58 . Particularly the results of Fernandez et al. 21 show that application of the Viterbi TbD algorithm over 10 frames of IR data yields about a 7 dB improvement in detection sensitivity over conventional thresholding peak-detection procedures. 4 CENTER FOR APPLIED MATHEMATICAL SCIENCES, USC DETECTION ALGORITHMS AND TBD ARCHITECTURE FOR IRST 5 This approach avoided problems with velocity mismatch and could handle targets with slow maneuvers. However, good performance of dynamic programming algorithms is observed for moderate SN+CR over 3 dB after preprocessing and clutter suppression with rapid degradation as SN+CR reduces further 58 . In addition, the computational complexity of these sophisticated methods is fairly high. 2.1.3. Extended Kalman Filter. To date, the extended Kalman lter EKF along with minor variations have been the dominant algorithm technology in real-time tracking. It is the basis for current practical single- or multi-target trackers for point equivalent targets. A major reason for its success has been the fact that the EKF has o ered a reasonable compromise between real-time operation and accurate performance in many nonlinear problems. On the other hand, the EKF is a suboptimal and completely heuristic algorithm whose e ciency varies from case to case. For instance, the EKF is unstable in situations that involve acute maneuvering, missing measurement, low SNR, multipath, and many other situations where the posterior distribution may not be approximated well enough by a Gaussian distribution. It is very di cult, if even possible, to develop rigorous evaluation metrics for assessing the quality of data processing based on the EKF technology. Improvements to the EKF e.g. iterated EKF work satisfactory in a number of important applications where the EKF fails, but still, it is di cult to overcome the fundamental limitations of the EKF algorithm that stem from its dependence on the assumption that the posterior distribution may be well approximated by a Gaussian distribution. The typical posterior density built upon the realistic IR image is shown in Figure 17, Section 7. One sees that it has multi-peak form and hence is very far from being Gaussian. Our experiments show that EKF is typically fails for this kind of data. 2.1.4. An Alternative TbD Method Optimal Nonlinear Filtering. In spite of the aforementioned shortcomings, the above outlined methods remain the basis for the great majority of existing signal processing systems. In particular, until recently no other information technology was able to e ectively compete with the EKF in target tracking. However the situation has changed with recent advances in the mathematical theory and algorithmic support for optimal nonlinear ltering ONF. These advances coupled with improvements to modern digital hardware technology make optimal nonlinear ltering an attractive alternative to the multi-dimensional matched lter, dynamic programming based algorithms and EKF in many practical and important applications. These include signal processing for infrared and acoustic sensors, imaging radar and sonar including SAR and SAS, and other passive and active sensors 17, 46, 50 . Advanced optimal nonlinear lters can now provide: real-time operation; optimal, theoretically sound solutions of the full nonlinear problem; distributional versatility no constraints on the form of prior or posterior probability distributions; 6 CENTER FOR APPLIED MATHEMATICAL SCIENCES, USC superior accuracy and robustness; facility to e ectively incorporate realistic physical models; explicit quantitative performance metrics exact error estimates, con dence areas, etc.. Our analysis of advanced algorithms based on ONF technologies shows great promise in: trackbefore-detect of unresolved targets in low SN+CR up to ,6dB after preprocessing; fusion of imaging and kinematic data for target identi cation; tracking of agile extended targets in cluttered environment as well as with certain other applications. We argue that just as the EKF superseded the , , trackers, so the optimal nonlinear lter is set to replace the EKF as the dominant tracking technology within 10-15 years. In Section 5 we describe the ONF techniques that will be used in de ning and developing the appropriate technology for application in an end-to-end IRST signal processing architecture. The architecture will be su ciently exible to allow for an adaptive optimization of the signal discrimination processors utilizing existing engineering parameters. This adaptive optimization will use the output of an EO IR sensor diagnostic tool to utilize current meteorological and environmental information as well as current intelligence on likely threat scenarios. 3. Signal and Observation Modeling It is assumed that a sensor has m,component resolution capability for radar typically m = 6, for IR EO m = 4. By r = r1; : : : ; rm will be denoted a phase coordinate vector for IRST typically angles and angle velocities of an object in a certain coordinate system. The sensor carries out a periodic surveillance of de nite area in m,dimensional domain D m R m , where R m is the m,dimensional Euclidean space. After standard preprocessing the result of one observation step is a frame of measurements, 3.1 Z n = S n + bn + n n = 1; 2; : : : ; or Z n = kZi nk = kSi n + bi n + i nk ; i = 1; : : : ; N; where S n = kSink is a signal from target, bn = kbink is an exterior background clutter, and = kink is a noise of the sensor. Here we assume that after preprocessing sampling of data is done in discrete points di, i = 1; : : : ; N , uniformly in the area D m i is a pixel. The noise is assumed to be zero mean and uncorrelated in both time n and space i, E in = 0, 2 E i n = E is a symbol of expectation. The clutter is de ned as 3.2 bin = bri + n; n where br; n is a function describing the background spatial distribution of the clutter after preprocessing in the point r, and n = 1n; : : : ; m n is an unknown current bias of sensor coordinate system with respect to the reference one due to the jitter. DETECTION ALGORITHMS AND TBD ARCHITECTURE FOR IRST 7 The signal component is modeled as 3.3 Sin = kn X j =1 Aj nhri + n , rj n where hr is a normalized sensor function; kn is an unknown total number of targets at the moment n; Aj n and rj n are unknown signal intensity and coordinates of the j th target, respectively. It is assumed that the clutter br; n has a relatively big spatial variance the change of br; n between two nearest values of ri is comparable with the maximum value of br; n. Besides, br; n varies locally as fast as the signal function does in spatial coordinates. Even in cases where these assumptions do not exactly hold, the algorithms that rely on them will be robust, which is the most important requirement. The second assumption shows that if only spatial information is used, then the background may be interpreted as a target. Thus spatial ltering alone is not su cient for clutter suppression and temporal ltering is needed. We always assume that b is an arbitrary unknown function of r and slowly changing function in n: there exists such T that jbr; n + T , br; nj . The latter assumption, which often holds, shows that variations of the clutter in time are caused mainly by uncontrolled vibrations n of the sensor. These vibrations are unknown and unpredictable except for a natural restriction, j nj , imposed on the absolute value of bias of coordinate system. No assumptions on statistical behavior of clutter is made. It is our belief that such popular models as homogeneous random eld, especially Gaussian, are valuable only for purely academic research and lead to highly non-robust ltering algorithms. Perhaps the most important feature of our approach to algorithm design is that we refuse to use any arti cial and unreliable statistical models of br; n, n, kn, and Aj n. The algorithms developed on the basis of such models fail even for small deviation of a model from reality. The essence of the new suggested approach is the development of algorithms that are invariant and or adaptive with respect to prior uncertainty. The speci c feature of the problem is its extremely high dimensionality. The value of N can be of the order of 106 , 108 and the total number of targets can be several hundreds. Thus, along with the mathematical issue of data processing, serious attention should be paid to the computational complexity, parallelization and HPC realization. 4. Clutter Suppression Two classes of algorithms for background ltering are proposed. In either case there are two basic problems to be solved: a to transform the sequence of input frames into the new frame such that the clutter would be reduced and the signal would be preserved clutter removal; b for every n the location of the coordinate system of the sensor i.e. the bias n should be estimated with maximum possible accuracy jitter compensation. Below we propose two algorithms that allow the rst problem to be solved e ectively. Jitter compensation algorithms will be developed in the near future. 4.1. Nonparametric Method. The rst class is relied on the nonparametric regression approach to estimation of the function br; n. Our goal is to build a nonparametric clutter estimate such that the residuals between the original data and its smoothed version estimate would be reasonably approximated by signal plus noise models. That is, the estimate ^r; n of the br; n b ~r; n = Z r; n , ^r; n the signal should be built in such a way that in the ltered frame Z b would be preserved, while the clutter would be removed almost completely. Kernel methods provide a powerful tool for such analysis due to both theoretic optimality see 19, 23 and computational transparency. In addition these methods are invariant to statistical properties and variations of clutter. Kernel estimators are weighted moving averages of observations X l ,r ; : : : ; lN ,r ; ^r; n = Qm1 b Z r; nK N Nm i=1 Ni r ;:::;rm 1 1 2 1 1 8 CENTER FOR APPLIED MATHEMATICAL SCIENCES, USC where Ni are window sizes in corresponding directions, K is a deterministic m,dimensional kernel recall that r = r1; : : : ; rm . A product kernel with univariate kernels Ki , i = 1; : : : ; m, provides a popular example. Note that in the 2D case m = 2 a so called `local mean removal' procedure used for clutter rejection in a series of papers by Reed et al. 39, 40 is equivalent to a kernel estimator with a product kernel K generated by uniform kernels K1 x = K2 x = 0:5Ifjxj1g. It is shown in Section 4.1 that the proposed adaptive nonparametric methods based on recent advances in nonparametric estimation see e.g. 18 , 27 perform much better than the aforementioned `local mean removal' method. 4.2. Semiparametric Filtering. The second method semiparametric is based on the adaptive spatio-temporal auto-regression 4.1 ~ Zin = Zin , n L X X =n,T l=0 al Zil n ~ ~ where Z n = jjZinjj is the ltered frame of input data Z n = jjZinjj; Zil n are the values of Zin in some vicinity spatial ltering window of points ri; L is the number of points rj in the window; al are some coe cients. The major problem is the choice of optimal coe cients al . They must be calculated adaptively to guarantee minimum of empirical mean-square MS value of the ltering residual noise internal noise and background residual for every time moment n. Such a criterion provides minimum of MS residual noise in the frame for every time n and is invariant with respect to statistical properties of the background bri + n; n and noise. The algorithm consists of: 1. Computation of empirical correlation matrix TL TL for spatial window L and temporal window T . 2. Calculation of the optimal weights al . DETECTION ALGORITHMS AND TBD ARCHITECTURE FOR IRST 9 ~ 3. Linear transform of Zin into Zin according to 4.1. In general the algorithm is nonlinear because of the nonlinear dependence of al on Zi. This algorithm allows us to suppress any background regardless of its variation in time. However, it solves only the rst part of the problem it does not estimate a drift of sensor coordinate system to correct coordinate measurements on the stage of detection jitter compensation. An algorithm that allows us to estimate n and hence to compensate jitter will be developed in the near future. 5. Optimal Nonlinear Filtering for TbD 5.1. Modeling for TbD. Our philosophy in addressing the low SN+CR detection problem regardless of the type of a sensor being used is that the farther the thresholding operation may be postponed the better. In the context of track before detect, we assume that the data are collected at discrete times while the target dynamics is a continuous time process. To be speci c, assume that the measurements Z k; r = Z tk ; r are collected at discrete time moments tk , k = 0; 1; 2; : : : and the relationship between the observation and target location is modeled by a nonlinear measurement equation of the form 5.1 Z k; r = S k; X k ; r + bk; r + V k; r ; where, as before, r represents the spatial coordinate in the phase space, bk; r is a sequence of deterministic unknown functions representing cluttered background, by X k is denoted the true location of the target at time t = tk and by V k; r is denoted the measurement noise process sensor noise. For simplicity we rst suppose that there may be only one target in the scene, i.e. the signal S k; X k ; r is described by 3.3 where kn = 1. Also we will neglect the platform instability jitter assuming n = 01. Under these conditions 5.2 S k; X k ; r = Akhr , X k ; where hr is a normalized sensor function and Ak is an unknown signal intensity. We further de ne Z k r = Z 1; r; : : : ; Z k; r to be the concatenation of all measurements up to time tk in the space point r. The expected range of possible behaviors" of the target is modeled as a Markov process. Often this process may be well described by a randomly perturbed multi-dimensional linear or nonlinear dynamical system _ _ X t = f X t + Wt ; t 0; _ where Wt is another noise process that describes uncertain and unpredictable target motion; f is a known function. Modeling of the state and observation is normally based upon physics. The models for the state process usually represent the a priori knowledge about physical, tactical, etc. characteristics of the target while the observation model is based upon the physics of sensors, clutter structure, operational environment, and so forth. 1 This is the case if the instability is compensated either by electro-mechanical stabilizators or by estimating. 5.2. Optimal Spatio-Temporal Nonlinear Filtering: The Basic Algorithm. The output of the ONF is the functional time series fk rgk1 of joint posterior densities JPD which measure the likelihood at the time moments tk that the vector of target features parametrized by X k is in close proximity of the grid point r. Computation of the JPD usually splits into two separate procedures: on every time step, spatiotemporal ltering for clutter removal is done rst see Section 4, and then the spatio-temporal nonlinear ltering is performed to estimate the target location. After clutter suppression procedure the clutter is decorrelated. We thus can include it into the noise component V k; r see 5.1 and concentrate only on the second procedure. It is a standard fact see e.g. 24 , 43 that the JPD is given by the formula p r k r = R p krdr ; 5.3 k where the function pk r, called the unnormalized ltering density UFD, propagates in time according to the following recursive equation of the predictor-corrector type: 1 5.4 pk r = expfR,1hr; ; R,1Z k; , 2 R,1hr; 2gTk pk,1r; k = 1; 2; : : : ; where R is the covariance of the observation noise V , and the predictor Tk pk,1r is a solution of the Fokker-Planck-Kolmogorov equation corresponding to the state process subject to the initial condition pk,1r, i.e. ut; r = Tt pk,1r is a solution of the equation r m @ut; r = 1 X a @ 2 ut; r , X @ f rut; r; t 2 t ; t ; i k k ,1 @t 2 i;j=1 i;j @ri @rj 5.5 i=1 @ri u0; r = pk,1r; where fai;j g is the covariance of the state noise and m is the dimensionality of X t. The posterior moments of X k = X tk can be computed now by integration of the respective polynomials against the JPD k r: For example, the best mean square estimate of X k , c Xk 10 CENTER FOR APPLIED MATHEMATICAL SCIENCES, USC = Z rk r dr: Recall that we consider a general problem assuming that r belongs to a m,dimensional Euclidean space phase space Rm our main concern in IRST problem is m = 4 if the problem is solved in the initial resolution" space angles and angular velocities or m = 2 if trajectories are projected on the plain. Now, let fe`rg`2N be a complete orthonormal system in L2 R m . Projecting 5.1 on this basis, we can rewrite the observation process in the coordinate form: ` 5.6 Zk = hl Xtk + Vk`; ` = 1; 2; : : : ; k = 1; 2; : : : ; K; where Z Z Z `= ` y = `= Zk Zk re`r dr; h hy; re`r dr; Vk Vk re`r dr and Vk` are independent Gaussian random variables. Rm Rm Rm DETECTION ALGORITHMS AND TBD ARCHITECTURE FOR IRST 11 Let F be a given function such that E jF X t j2 1, 8t 0. In general, the optimal nonlinear ltering deals with obtaining the best mean square estimate of F Xtk , given the measurements Z k r = Z1 r; : : : ; Zk r, r 2 R m . This estimate is called the optimal lter and will be ^ denoted F k. By Tt denote the solution operator for the Fokker-Planck-Kolmogorov equation corresponding ^ to the process X . The optimal lter F k = E F X tk j Z k is given by the formula R p rF r dr ^ R F k = Rm k p r dr ; 5.7 Rm k where the UFD pk r was introduced above. This density obeys the recursive relation p0r = r; 1 ! 1 X 5.8 X jh`rj2 T p r; k = 1; 2 : : : : ` pk r = exp h`rZk , 2 k ,1 `=1 `=1 Additional notation: ` `j qmn = hem; Teni; qmn = hh`em ; Teni; qmn = hh`hj em ; Teni; where hF; gi = observations become available. R Rm F r g r dr. m 0 = hp; emi; Fm = hF; emi; Note that these quantities can be computed o -line before any ^ Let L = f` : ` M g and F k be the approximation to the optimal lter F k. ^ We propose the following algorithm for computing the optimal nonlinear lter F k, which will be called the Spectral Separation Scheme or the S 3 algorithm. 1. Set a cut-o level M for the number of basis elements. 2. Before the observations become available compute: X X p0r = ` 0e`r and 0 f = `0F`: `M `M ` 3. When the measurements Zk , k = 1; 2; : : : ; K , ` 2 L L become available: a: compute recursively L L m 0 = m 0; m k = L QmnZk = qmn + X X nM L L Qmn Zk n k , 1; k = 1; 2; : : : ; X b: compute pk r = `2L ` ` Zk qmn + 2 `6=j ;`;j 2L X ` j `j Zk Zk qmn + 2 X, `2L `` ` Zk 2 , 1 qmn; X `M `Lke`r; k F = `M `LkF`; F k = k F =k 1 : 12 CENTER FOR APPLIED MATHEMATICAL SCIENCES, USC The block diagram of the Spectral Separation Scheme S 3 algorithm is shown in Figure 2. It is clear that the above algorithm allows for an adaptive, dynamic interplay between data and computations. Indeed, it does not involve solving of PDE's on-line and is recursive in time. Furthermore, it is also spatially recursive. These properties are very important since they allow adaptive sequential multi-resolution ltering to be performed. Observe that in order to achieve the e ect of multi-resolution one has to use local hierarchical bases fel g; for example wavelet bases 25 . (k -1) observation Optimal Nonlinear Filter (k ) = Zk (k) pk (r ) = n M Qn (k - 1) Q pk , k , F k Estimation lM l ( k ) el ( r ) Q - precomputed matrix Figure 2. S 3 algorithm on-line part 5.3. Multi-Dimensional Spatio-Temporal Matched Filter as a Special Case of Nonlinear Filter. For the sake of simplicity let us consider the basic formula for updating the ltering density in the optimal nonlinear lter in the 3,D case where the target motion is projected on the plain. The modi cation for the m + 1,D case with arbitrary m, particularly for m = 4, is straightforward. The grid space will be denoted by fxij g. Suppose that the solution of the Fokker-Plank-Kolmogorov equation the mean motion dynamics of the target is given by the operator T and observations are made at time moments tk = k on the observation space grid fxij g with additive Gaussian white noise: 5.9 Zij tk = hXtk , xij + ij tk ; where h is the target signature signal and ij tk are i.i.d. standard Gaussian variables. Then by 5.4 the updating formula is n X o X pk x = exp 12 hx , xij Zij tk , 21 2 hx , xij 2 Tpk,1x : i;j i;j DETECTION ALGORITHMS AND TBD ARCHITECTURE FOR IRST 13 In the above formula, the part Tpk,1x is the prediction part, i.e. the best possible guess based on known target trajectory dynamics, what the probabilities to nd the target in various areas of the observation region would be if there were no observations available at the time moment tk . The exponent part is the correction part, i.e. after receiving the observation at time tk , the prediction part is ampli ed or attenuated, depending on the expression inside the exponent. One notable fact is that the ampli cation or attenuation is directly related to the output of the well-known in engineering community single frame matched lter spatial 2,D matched lter X mk x = hx , xij Zij tk : Indeed, the ltering densities are ampli ed only at those points where the output of the matched lter exceeds the value of X 1 E2 = 1 hx , xij 2 ; 2 2 and attenuated elsewhere. Here E2 is the energy of the signal which of course does not depend on x. This establishes a strong connection between a spatial matched lter and the optimal nonlinear lter. Even more importantly, the optimal nonlinear lter is equivalent to the spatio-temporal matched lter or assumed velocity matched lter if the movement of the target occurs along a deterministic trajectory. In that case the operator T degenerates into an operator that shifts the values of the function along the known deterministic trajectory. For simplicity, let us assume that the target moves in a straight line with known constant velocity v. Then Tpx = px , v; and it is easy to see that for any k 1 k,1 k,1 n 1 XX 1 X X hx , x , `v2 op x , kv: hx , xij , `vZij tk , 2 2 pk x = ck exp 2 ij 0 `=0 i;j `=0 i;j Notice that the rst part in the exponent in the above formula, i;j i;j Mk x = k,1 XX `=0 i;j hx , xij , `vZij tk ; exactly coincides with the usual assumed velocity matched lter and the second part k ,1 1 X X hx , x , `v2 = E2 k ij 2 2 `=0 i;j 2 2 is nothing but the accumulated SNR. Thus, if we assume that the initial density is a deltafunction, then statistical inference base on Mk x and pk x will be the same. Particularly, the points of maxima of Mk x and pk x coincide. A similar argument can be applied to any dimensions. Recall that in general we are interested in the m + 1,D case, r = r1 ; : : : ; rm, m + 1,th component is time. The nal result is 14 CENTER FOR APPLIED MATHEMATICAL SCIENCES, USC exactly the same: in case of deterministic trajectory the m +1,D matched lter and the optimal spatio-temporal nonlinear lter coincide with the accuracy to the deterministic transformation. We stress that the developed optimal spatio-temporal nonlinear lter does more than the spatio-temporal matched lter does. It is a generalization of the multidimensional spatio-temporal matched lter and coincides with it for targets that move along deterministic trajectories. In general, however, the ONF predicts the trajectory and allows us to align frames coherently even for acutely maneuvering targets when the m + 1,D matched lter fails. 6. Track Appearance Disappearance Detection 6.1. Preliminaries. A problem of detecting target tracks that occur appear and disappear at a priori unknown points in time is a typical abrupt change detection problem. The change-point detection problem has been actively researched during the last three decades see, e.g., 22, 35, 38, 42, 47, 45, 53, 54, 55 and references therein. Heuristic procedures such as Shewhart's charts and some modi cations appeared in the late twenties and early thirties. Detection procedures, which are in current use, were initiated by Girshik and Rubin 22 and Page 35 for the problem of detecting a change in a mean of i.i.d. Gaussian sequence and by Shiryaev 45 in the general but i.i.d. case. There are two major competitive procedures: the Shiryaev-Roberts-Girshik-Rubin algorithm 42, 38, 55, 45 and the Page's or CUSUM cumulative sum procedure 35, 47, 56 2 . There are two major disadvantages to both algorithms: 1. Only one change in data is assumed, i.e. in our context the moment of target disappearance is ignored. 2. The aforementioned change-point detection procedures have optimal properties only in the i.i.d. case in the following speci c sense: they minimize the mean detection delay under constraints on the mean time between false alarms for exactly speci ed pre-change and post-change distributions. These drawbacks make the conventional change-point detection procedures impractical in multitarget surveillance problems that involve non-i.i.d. observations and when one needs not only to detect a particular target with unknown position but also to discriminate between false alarms and true tracks. One appropriate criterion for our purposes is to x probabilities of false alarm and detection in a xed size sliding window and to minimize the detection delay. Or, alternatively, one may require to detect the track with xed or maximum probability during the xed time and with minimum detection delay. Also it is imperative to design the decision statistics that will account for unknown target location and moment of its disappearance. Otherwise the problem is meaningless. Neither of the above mentioned algorithms may be directly applied in this context. However, we will show that the CUSUM-type procedure and the quasi-Bayesian Shiryaev-Roberts-Girshik-Rubin type procedure with specially designed thresholds are useful 2 See also 41, 55 for reasonable modi cations and ideas that are useful for surveillance systems. DETECTION ALGORITHMS AND TBD ARCHITECTURE FOR IRST 15 tools, especially in the problem of detection of a single target that appears and disappears at unknown times. Moreover, it turns out that the thresholds can be chosen such that the prede ned false alarm rate is xed for a large class of noise and signal models, not necessarily for the i.i.d. models. This innovation is very important for our applications, since it removes the restrictive i.i.d. assumption, which is typical for previous work in change-point detection. In particular, it allows us to design thresholds for adaptive algorithms. Now we summarize the basic requirements that should be satis ed in realistic problems of detection of target's track appearance disappearance: detection algorithms should be invariant relative to the completely unknown moments of abrupt change the prior distributions of appearance disappearance moments are unknown; detection algorithms should be adaptive and use estimates of unknown target location based on TbD; the frequency of false detections or the false alarm probability is limited at the speci ed level; the probability of correct detection should be maximal during the xed time interval; the mean detection delay should be made as small as possible. We propose the sequential algorithms that meet all these requirements and even more, see Algorithm 3 in Section 6.5. The algorithms are based either on the generalized CUSUM statistic or on the quasi-Bayesian statistic. These statistics are, in essence, score functions that characterize the likelihood ratio of hypotheses on signal presence and absence. Speci cally, the generalized CUSUM statistic is the likelihood ratio maximized over unknown moments of appearance and disappearance and the latter statistic is the average likelihood ratio with respect to at improper prior distributions of these moments. The statistics are computed either in a sliding window of a xed size or on the semiin nite interval up to the current moment. In both cases the decision statistics are compared with thresholds that depend on a given false alarm rate probability of false alarm or frequency of false alarms. It turns out that in many cases the generalized CUSUM procedure may be reduced to the conventional CUSUM without loss of statistical performance. To be precise, if the major goal is to detect the fact of target's appearance but not the fact of its disappearance, than the unknown moment of target's disappearance should be considered as a nuisance parameter. In this case the CUSUM is the optimal algorithm, i.e. the unknown moment of target's disappearance does not a ect the structure of the algorithm, but only its performance see Section 6.3 and Section 6.4 for details. The same is true for the quasi-Bayes procedure. If, however, one has to detect multiple targets, then it is important to detect both track appearance and track disappearance as soon as possible. Moreover, in this case it is also important to have a special logic to discriminate between false alarms and true tracks. In this multi-target context" the conventional CUSUM-type procedures do not solve the problem but still may be used as a part of more sophisticated algorithms see Section 6.5 for some discussion. It should be mentioned that the problem of signal detection with random appearance and disappearance times in ideal conditions complete prior information was solved by Tartakovsky 16 CENTER FOR APPLIED MATHEMATICAL SCIENCES, USC 57 . The results of this work show that the optimal Bayes solution requires, in general, computation of the 5,dimensional su cient statistic and is completely impractical. Thus even the availability of the required prior information, particularly knowledge of prior distributions of appearance disappearance moments which are usually unknown, does not help in practice. 6.2. Problem Formulation and Model Description. Let Yk , k = 1; : : : ; n, denote the data obtained up to the current time n. Generally, the data Yk = fZ k ri; i = 1; : : : ; N g represent the whole spatial set" collected from the frame of observations after preprocessing at time k or the part of these data in some channels". If the clutter removal procedure is rst performed, then after spatio-temporal ltering the clutter is decorrelated. Thus, it is reasonable to assume that fYk g are i.i.d. with the density p1y if the target is present and with the density p0y if it is absent3 . The density p1Yk j k depends on in general unknown vector parameter k , which characterizes the spatial location of the target at moment k. First, for the sake of simplicity we assume that the target dynamics is known and hence the k is xed and known. The exact models with the account of TbD and corresponding adaptive algorithms will be de ned later on in Section 6.8 and Section 6.9 see also Section 6.6 in the general case. By k = log p Yk denote the log-likelihood ratio LLR for the single observation Yk and by p Yk and the unknown moments of target appearance and disappearance, respectively. Under our assumptions the data Y1; : : : ; Y,1 are i.i.d. according to the density p0 y, Y; : : : ; Y are i.i.d. according to p1y and Y +1; Y +2; : : : ; n again follow p0 y. We assume nothing about and except for the natural constraint . In the language of hypotheses testing the problem of the track detection may be formulated as testing of the hypotheses 1 0 H1;l;m : track appears at = l and disappears at = m; 1 l n; m l + 1"; H0 : track does not appear, i.e. 62 1; n ": Alternatively, the hypotheses may be written in the form ,1 Y k=1 H1;; : pY1n = p; Y1n = n Y k=1 p0Yk Y k= p1 Yk n Y k= +1 p0Yk " ; H0 : pY1n = p0Y1n = p0Yk ": 6.2.1. Generalized Likelihood Ratio Statistic. By Ln denote the generalized LLR, X p Y1n Ln := max log ; , n = 1maxn k : k= ; p0 Y1 +1 , The derived algorithms can be easily modi ed and generalized for more general models that include correlated and non-homogeneous observations, see comments below and Section 6.7. 3 DETECTION ALGORITHMS AND TBD ARCHITECTURE FOR IRST 17 It is easily shown that the statistic Lk satis es the following system of recursive relations: , , 6.1 Lk = max Lk,1; Uk = Lk,1 + Uk , Lk,1 +; 6.2 Uk = k + Uk+,1; k = 1; : : : ; n; with the initial conditions L0 = U0 = 0. Here y+ denotes a non-negative part of y, i.e. y+ = max0; y. Indeed, Ln = max 1max 1 where for any k 1 X ! k= " k = maxU1; : : : ; Un = maxLn,1; Un ; Uk =: 1maxk k X s = max k ; 1max,1 k + s k s= s= k,1 X ! = k + max0; Uk,1 : 6.2.2. Generalized Average Likelihood Ratio Statistic. Introduce the following statistic , Z XY p; , Y1n d = max Gn := max expk : n 1 0 p0 Y1 =1 k= It is easy to see that Gn = maxR1; R2 ; : : : ; Rn = maxGn,1; Rn; G0 = 0; P Q where Rn = n=1 n= ek . The statistic Rn may be interpreted as the likelihood ratio averaged k over the at uniform distribution of the moment of track appearance. It is easy to show that fRng obeys the recursive relation 6.3 Rn = expn1 + Rn,1; R0 = 0 : Both statistics, Ln and Gn, can be used for testing the above hypotheses either in the xed interval or sequentially. 6.3. Detection Algorithm 1: Fixed Sliding Window. The hypotheses H1;; and H0 are tested in a sliding window of length T at each current time instant n, i.e. p0 Yk p1Yk p0Yk " ; k=n,T +1 k= k= +1 n Y H0 : pYnn,T +1 = p0 Ynn,T +1 = p0Yk ": k=n,T +1 H1;; : pYnn,T +1 = p; Ynn,T +1 = In this case the statistics LT;n and GT;n are determined by , 6.4 LT;k = max LT;k,1; Uk ; 6.5 Uk = k + Uk+,1; k = n , T + 1; : : : ; n; ,1 Y Y n Y 18 CENTER FOR APPLIED MATHEMATICAL SCIENCES, USC and 6.6 6.7 GT;k = max GT;k,1; Rk ; Rk = ek 1 + Rk,1; k = n , T + 1; : : : ; n; , with the initial conditions LT;n,T = Un,T = GT;n,T = Rn,T = 0. According to the adaptive Bayesian approach, the procedures of testing the hypothesis H1 against H0 have the form 6.8 6.9 dL = T;n dG T;n = 1 if LT;n a; 0 otherwise; 1 if GT;n B; 0 otherwise; where a, B are thresholds that are de ned on the basis of the given probability of false alarm and d = 1 stands for the decision on target presence while d = 0 is the decision on its absence. It should be emphasized that in this subsection we pursue the goal of detection of the target track appearance regardless of the possibility of its disappearance in the analyzed observation interval. In other words, the unknown moment of target disappearance is considered as a nuisance parameter. A di erent setting when both target's appearance and disappearance should be detected will be considered in the next two sections. We now notice the following remarkable property of the tests 6.8 and 6.9. Since LT;n = n,Tmaxkn Uk ; +1 we have fdL = 0g = fUk a for all n , T + 1 k ng : T;n a T; n = minfn , T This shows that the non-sequential test 6.8 is equivalent to the sequential test which is de ned by the stopping time 6.10 + 1 k n : Uk ag; a T; n = 1 if no such k: If a T; n 1, then the hypothesis H1 is accepted target track is present, while if a T; n = 1, the hypothesis H0 target is absent is accepted. The same argument is applied to show that the non-sequential test 6.9 is equivalent to the sequential one 6.11 ~B T; n = minfn , T + 1 k n : Rk B g; ~B T; n = 1 if no such k: The sequential procedures 6.10 and 6.11 allow us to achieve exactly the same probabilities of false alarm and target missing as the initial tests 6.8 and 6.9, respectively. In addition, they have an advantage: the time delay in signal detection is less. See Figure 3 for graphical explanation of the above argument. DETECTION ALGORITHMS AND TBD ARCHITECTURE FOR IRST 19 LT,n a Un threshold window n-T+1 a (detection) time n Figure 3. Detection in the sliding window: the statistic Un exceeds the threshold at random moment a , which is less than the size of the sliding window T . We also observe the following interesting fact. The statistic Uk is nothing but the CUSUM statistic, which is de ned in the standard CUSUM procedure as Uk = max k k X s= s; i.e. as the maximum of the LLR over unknown moment of signal appearance only. In turn, the statistic Rn is nothing but the quasi-Bayesian statistic, which is the limit form of the Bayesian statistic the Shiryaev-Roberts-Girshik-Rubin statistic. This shows that the presence of an additional nuisance parameter moment of signal disappearance , does not a ect the nal structure of the test procedures. But it certainly does a ect the detection performance, since the accumulated SNR becomes less as decreases. ~ Another possible modi cation of the algorithm is to use the statistic Uk with re ection from zero barrier, ~ ~ 6.12 Uk = k + Uk,1+ ; k 1; U0 = 0 ; in place of Uk in 6.10. Obviously, these two tests have the same performance as long as the ~ threshold is positive. Indeed, it is easy to see that the trajectories of the statistics Uk and Uk coincide in the non-negative half-plain. 6.4. Detection Algorithm 2: Fully Sequential Procedure. We now change the problem set-up and consider the fully sequential procedure. To be speci c, we assume that k = n; n,1; : : : and the decision about target presence or absence is made at each moment. This case is probably the most relevant to the IRST surveillance systems, which work for a long time periodically raising false alarms. It is easy to see that the generalized LLR L,1;n := Ln is Ln = maxUn; Un,1; ; 20 CENTER FOR APPLIED MATHEMATICAL SCIENCES, USC where Uk is de ned in 6.5. Formally the sequential algorithm that solves the problem is identi ed with the stopping time a = inf fn : Ln ag: The stopping time a can be rewritten as 6.13 a = inf fn : Un ag; which is exactly the CUSUM procedure. Indeed, the decision on target absence dn = 0, which is equivalent to observation continuation at time n, is made when fLn a; Ln,1 a; g = fUn a; Un,1 a; g : Thus, again prior uncertainty with regard to the moment of disappearance does not a ect the structure of the detection algorithm if is considered as a nuisance parameter but it does a ect its performance. It is important to understand that the algorithm 6.13 requires additional strategy after the decision d = 1 target appeared is made. Since we expect multiple signal appearances, one may not simply stop observations after this decision is made. One way to handle this problem is to set Un = 0 when a = n and to start all over again immediate renewal. Then the algorithm is immediately ready to detect the next target. In this case the CUSUM statistic is modi ed compared to 6.5 as follows 6.14 Un = n + Un,11lf0 Un, ag ; where 1lfY g is the indicator of the set Y . This modi cation has, however, one drawback: it may and usually does lead to multiple detections of the same signal and hence requires additional identi cation logic. 1 Another reasonable way is to use the detection rule 6.13 with pure CUSUM 6.5 supplemented by the following logic for new target detection. Let a be the rst time such that Un goes below the threshold a after exceeding. To be precise, a = inf fn a : Un ag: Then for a k a , 1, one con rms the presence of the same target, while if for some m 1, Ua+m a, then we make the decision that a new target track appears at time a + m. In other words, if the i,th target was detected at time ai and con rmed till ai , 1, then the i + 1,th target is detected at i+1 = inf fn i : U ag n a a i+1 and its track is con rmed till the moment a , 1, where ai+1 = inf fn ai+1 : Un ag: We also note that one may replace the statistic 6.14 by the statistic ~ 6.15 Un = n + Un,11lf0 n+Un, ag : This replacement does not a ect performance. 1 DETECTION ALGORITHMS AND TBD ARCHITECTURE FOR IRST B 21 Applying a similar argument to the quasi-Bayesian statistic Gn, we obtain that the stopping rule = inf fn : Gn B g = inf fn : Rn B g: can be rewritten as 6.16 B In the case of immediate renewal the average likelihood ratio Rn is computed recursively according to the formula 6.17 Rn = expn 1 + Rn,11lfRn, 1 Bg : The same logic as above can be used to detect not only the fact of target appearance but also its disappearance. However, in the case of multiple targets it is better to consider a more complex set of hypotheses that would include at least three alternatives: the last signal did not appear, appeared but is still present, appeared and disappeared. This case is considered in the next section. We stress once again that the statistic Rn may be interpreted as the average likelihood ratio over the uniform distribution of the on the interval 1; n , while the expUn is the maximum of the same likelihood ratio over 2 1; n . The method 6.16 has an advantage over the generalized CUSUM algorithm the statistic ~ Rn = Rn , n is a zero-mean martingale with respect to P 0 regardless of the i.i.d. assumption on the observations. As a result, if one sets B = 1=Fr, then the frequency of false alarms is upper bounded by the prespeci ed value Fr see Section 6.7.1 for details. In other words, the threshold B is easily estimated for a large class of models. There is no similar result for the CUSUM statistic Un . 6.5. Algorithm 3: Joint Detection of Track Appearance and Disappearance. Consider an extended decision making process assuming that after each decision on target disappearance the previous data are discarded. Speci cally, we accept the following logic: 1. If the decision that the i,th target disappeared is made di = DA, then one of the two decisions may be made the i + 1,th new target did not appear" di+1 = NA = 0 or the i + 1,th target appeared and is still present" di+1 = A&P. 2. If the decision di+1 = 1 that the i + 1,th target appeared is made, then one of the two decisions may be made target is still present" di+1 = P or target disappeared" di+1 = DA. The statistics Ln and Gn are then the likelihood relative to H0 target is absent at all of the composite and combined hypothesis H that includes two sub-alternatives: H : target appeared and is still present HA&P + target appeared and disappeared HA&D": 22 CENTER FOR APPLIED MATHEMATICAL SCIENCES, USC Clearly, at time moment n the likelihood of the subalternative HA&P is determined by the statistic Un, while the likelihood of the subalternative HA&D by Ln,1 = maxkn,1 Uk . Thus the structure of the decision making algorithm may be as follows. After the decision di = DA the i,th target disappeared is made, the statistic Un is formed with the null initial condition. This statistic is compared at each step with the threshold a. If Un a, the observation is continued. Otherwise Un a, the decision di+1 = 1 the i + 1,th target appeared is made, the statistic Un is computed further and also the algorithm starts to compute the statistic Ln with the initial condition L ai = U ai , where ai+1 is the moment of detection of the i + 1,th target. The di erence n = Ln,1 , Un is compared with another threshold b at each step n = ai+1 + 1; ai+1 + 2; : : : . If n b, then the decision di+1 = P still present is made and the next time step is analyzed. If for some n = bi+1 , n b, then the decision di+1 = DA target disappeared is made, the statistic Ln is not computed, but the statistic Un is computed for n = bi+1 + 1; bi+1 + 2; : : : with the initial condition Ubi = 0. Then the whole cycle is repeated. +1 +1 +1 Thus, the track appearance of the i + 1,th target is detected at the moment i+1 a = inf fn bi : Un ag i a and its disappearance is detected at ai+1 = inf fn : n bg : The value of the threshold b may be computed based on the trade-o between the error probabilities due to too early decision about target disappearance and the delay of this decision when the target really disappears. CUSUM Un {Yk } Data Thresholding, a d = 1 (appeared) Set zero initial conditions LLR of H = HA&P + HA&D n = Ln-1 - U n Thresholding, b d = DA (disappeared) Figure 4. Block diagram of Algorithm 3 DETECTION ALGORITHMS AND TBD ARCHITECTURE FOR IRST 23 Ln a Un b threshold threshold b ( i -1) a app. detection, i-th target (i ) b (i) time n disapp. detection, i-th target Figure 5. Detection of track appearance disappearance by Algorithm 3 The similar decision making algorithm may be built upon the statistics Rn and Gn. Then the appearance of the i + 1,th target is detected at the moment ~Bi+1 = inf fn ci : log Rn log B g ~ while its disappearance is detected at time ~ ci+1 = inf fn ~Bi : n cg ; ~ ~ ~ where n = log Gn,1 , log Rn with the initial condition ~Bi = 0 and where c is a threshold. The block diagram of the algorithm and the typical behavior of the decision statistics are shown in Figure 4 and Figure 5, respectively. 6.6. Adaptive Detection Algorithms. Recall that in general the position of the target is unknown and the density p1Yk j k depends on the parameter k that characterizes the position. Then instead of unknown k one may use the estimate of the position. This estimate can be obtained by applying the TbD procedure based on ONF. Thus the development of adaptive versions of the above detection algorithms is needed. The statistics Un = Un 1 ; : : : ; n and Rn = Rn 1 ; : : : ; n are the functions of the sequence of theta's till time n. Thus the developed procedures can not be applied directly. If the n can be ^ estimated this is the case in our system, then the natural solution is to use the statistics Un = ^ ^ ^ ^ ^ Un 1 ; : : : ; n and Rn = Rn 1 ; : : : ; n , which are computed based on the recursive formulas ^ ^ ^+ ^ ^ Un = n n + Un,1; Rn = en ^n 1 + Rn,1 : ^ ^ Here k = k Y1k is an estimate of k based on the previous k observations. However, it is very di cult to evaluate the performance of the corresponding adaptive algorithms. ^ ^ The reason is that Un is not a CUSUM and Rn is not an average likelihood ratio anymore. In 24 CENTER FOR APPLIED MATHEMATICAL SCIENCES, USC ^ ^ fact, the expfn n g is not the partial likelihood ratio, since p1 Yk j k Y1k is not a probability density. Particularly, the most important question of how to choose the thresholds remains open. To avoid this di culty we propose the following trick. At stage k instead of using the estimate ^ ^ ^ k we propose to use the one stage delayed estimate k,1 Y1k,1. In other words, instead of Un and R , which satisfy the recursions ^ and Rn we will use the statistics Un n 6.18 Un = + Un,1 +; n = 2; 3; : : : ; U1 = 0; n n 6.19 Rn = e 1 + Rn,1; n = 2; 3; : : : ; R1 = 0; ^ with = n n,1. The corresponding adaptive target detection procedures have the form n 6.20 6.21 log-likelihood ratio. As a result, the statistics 6.18, 6.19 preserve most of the nice properties of the former statistics Un and Rn. In particular, Rn , n is a P 0,martingale with mean zero. This fact allows us to upper bound the frequency of false alarms in the adaptive algorithms regardless of speci c pre-change and post-change distributions see Section 6.7.1. The block-diagram of a typical adaptive track appearance disappearance detection algorithm that uses the estimates of target location from the ONF-based TbD scheme is shown in Figure 6. TbD (ONF) for Position Estimation B g: ^ Since p1Yn j n,1 is the probability density for any Y1n,1,measurable estimate of n, is the n ~ B n a = inf fn 2 : U ag; = inf fn 2 : R n ^ k Y { k} Decision Statistics (CUSUM or ALLR) * Un Thresholding R * n Track Detection Figure 6. Block diagram of a typical adaptive detection algorithm 6.7. Choice of Thresholds and Performance Evaluation. In this section we give an argument and some ideas that can be used for performance evaluation of the proposed detection algorithms. We focus on the rst two algorithms. The analysis of the third algorithm is more di cult and we leave this problem for the future. In both algorithms Algorithm 1 and Algorithm 2 it is desirable to x the probability of false alarm Pfa in a given xed interval of length T . We note that while in Algorithm 1 this is the length of the sliding window, in Algorithm 2 this is the arti cial", auxiliary parameter that is de ned based on tactical conditions4. For the fully sequential algorithms Algorithm 2 it is also If the dense ow of targets from a particular direction is expected, T should be chosen much less than the typical average duration of a signal. Otherwise T should be comparable with the average time of target presence. 4 DETECTION ALGORITHMS AND TBD ARCHITECTURE FOR IRST 25 reasonable to x the frequency of false alarms, Fr = 1=E 0 , at the speci ed level Fr. In the next subsection we show that it may be done in the general, not necessarily i.i.d. case. 6.7.1. Upper Bounds for the Frequency of False Alarms. Consider the general case where observations can be dependent and or non-stationary. To be speci c, assume that under H0, pY1n = p0Y1n, while if the target appears at the moment and is still present at time n the model is pY1n = p0Y1,1 p1Yn j Y1,1 ; where p0 is the joint density of the noise and p1 describes probabilistic properties of the Yn n mixture of signal and noise. By g Y1n = pY n denote the likelihood ratio of the hypotheses p H and H0 . Obviously, ,1 k n = g n,1 g n k = p1 Y j Y1 : g where g n p0Yk j Y1,1 Hence for the CUSUM and the average likelihood ratio statistics we have 1 0 1 6.22 6.23 n n n n + Un =: max log g = max0; max1 log g ,1 + log gn = log gn + Un,1 ; U0 = 0; n n, Rn =: n X =1 n g n = gn n,1 X =1 n n n g ,1 + gn = gn 1 + Rn; R0 = 0: n Particularly, in the i.i.d. case log gn Y1n = nYn and the relations 6.22 and 6.23 coincide with 6.2 and 6.3, respectively. In the latter case both statistics are also Markov processes, which helps to evaluate performance of the detection algorithms. But in general the statistics are non-Markov and the analysis is more di cult. The most important question is how to choose the thresholds in order to guarantee a speci ed n false alarm rate. To answer this question we rst observe that E 0 gn j Y1n,1 = 1 for any n 1 and any model. As a result, it is easily seen that for every n 1 the statistic Rn , n is a P 0 ,martingale with mean zero, E 0 Rn , n j Y1n,1 = Rn,1 , n , 1; E 0 Rn = n; E0R and this is valid for an arbitrary model. Then, by the optional stopping theorem, for any Markov time , de nition of the stopping time ~B , E 0R~B B . This implies 6.24 6.25 E 0 ~B = E0 , while by the B; which holds for any model of signal and noise. Therefore if we set B = 1=Fr; then the frequency of false alarms will be upper bounded by Fr. 26 CENTER FOR APPLIED MATHEMATICAL SCIENCES, USC It turns out that the martingale property of the statistic Rn can be used to obtain the upper bound for the frequency of false alarms in the CUSUM procedure. Indeed, for any n 1 n expUn = 1maxn g Rn = n X and hence on the event f a 1g, expU a R a . Thus, 1, then ea E 0 expU a E 0 R a = E 0 a ; =1 if E 0 a n g where the left inequality follows from the de nition of the stopping time a and the right equality from the martingale property of Rn. This shows that for any model not necessarily i.i.d. that guarantees niteness of E 0 a we have the inequality 6.26 E 0 a ea : Hence if a = log1=Fr; 6.27 then the frequency of false alarms in the CUSUM-type procedure is upper bounded by Fr. Finally, we note that 6.24 and 6.26 also hold for the adaptive sequential algorithms 6.20, n 6.21, since these algorithms are particular cases of 6.22, 6.23 with gn = expf g. n 6.7.2. The Case of i.i.d. Observations. We rst observe that in the i.i.d. case the statistic Un is a Markov process with the transition probability densities under Hi, i = 0; 1, , 6.28 piunjun,1 = fi un , gun,1 ; where fiy is the probability density of the LLR k = Uk , gUk,1 under Hi; gy = y+ for Algorithm 1 xed sliding window and gy = y1lf0 y ag for Algorithm 2 with immediate renewal see 6.14. So is the statistic rn = log Rn with the transition probability density , piun j rn,1 = un,1 = fi un , log 1 + eun, 1lfun, bg for the quasi-Bayesian procedure see 6.17. The analysis of these algorithms may be done by using the renewal theory see, e.g., 20, 47, 60 . 1 1 In what follows we consider only the CUSUM-type procedure. For the quasi-Bayesian algorithm all formulas are similar. The false alarm probability is 6.29 PfaT; a = P 0 a T = 1 , P 0 U1 a; : : : ; UT a; while the target missing probability is 6.30 PmisT; ; ; a = P ; a T = P ; U1 a; : : : ; UT a; where the statistic Un obeys the recursion 6.31 Un = n + gUn,1; n = 1; : : : ; T : The threshold is chosen from the equation PfaT; a = 0 ; DETECTION ALGORITHMS AND TBD ARCHITECTURE FOR IRST 27 where 0 is a speci ed false alarm rate. Therefore to choose the thresholds and to evaluate the probabilities of correct detection one needs to nd the probability of entry in the set a; 1 by the Markov process Un during the time T . This process starts from 0 for Algorithm 1 and from some random point in ,1; a for Algorithm 2. In general, let fXn; 0 n T g, X0 = x x is xed be a Markov process with the transition density f y; z Xn,1 = z ! Xn = y. De ne the system of functions pk y; x by the recursion 6.32 Z a pk+1y; x = ,1 pk y; zf z; x dz; k 1; Z y p1y; x = F y; x where F y; x = The function pk y; x is nothing but the conditional probability of the event fXk a; : : : ; X1 ag given X0 = x. Then 6.33 P X1 ,1 f z; x dz : y; Xk,1 a; : : : ; XT a j X0 = x = pT a; x : 6.7.3. Performance of Algorithm 1. Similar to 6.32 for i = 0; 1, de ne the functions pki y; x by p1i y; x = Fiy , x+ ; where fi y , x+ is the transition density function of the statistic Un under Hi that satis es + 6.31 with gUn,1 = Un,1 and the null initial condition, U0 = 0. The probabilities of errors are then computed as follows 6.35 PfaT; a = 1 , p0 a; 0 ; T 1 6.36 PmisT; ; ; a = pT a; 0; if n , T + 1; n : If n , T + 1 or and n, then the situation is more delicate. It may be shown that in this case the probability of target missing 6.37 6.34 pki+1y; x = Z a ,1 pki y; zfiz , x+ dz; k 1; PmisT; ; ; a = Z a Z a if ,1 ,1 p0 a; p1 +1d; p0 1d; 0 n, , , n: ng and f n , T + 1; ng, n , T + 1; To obtain the mis-detection probability for f n , T + 1; one has to put p0 1 d; 0 = d for f n , T + 1; ng ; , p0 a; = 1 for f n , T + 1; ng n, in 6.37 x being a delta-function. Particularly, this yields the relation 6.36. 28 CENTER FOR APPLIED MATHEMATICAL SCIENCES, USC It is worth mentioning that the probabilities pki a; 0 satisfy the integral equations p0i y; x = 1 ; which are simpler than 6.34. The formulas 6.35-6.37 determine the performance in terms of the probabilities of false alarm Pfa and target missing Pmis at each current time n within the window of the length T . In particular, the threshold a is found from the equation see 6.29 1 , 0 = p0 a; 0 T 0 where the function pT a; 0 satis es the recursive relation integral equation 6.38. Also, the probability of target missing in the interval of its presence, ; , can be approximated by5 Pmis , + 1; a p1 +1 a; 0 : , The following elementary estimates are useful for initial choice of threshold and target missing probability evaluation: 6.39 PfaT; a P 0ST a; Pmis , + 1; a 1 , P ; S , S,1 a P where Sn = n=1 k . k Along with the probabilities of errors it is interesting to evaluate the mean time delay in target detection E ; f a T; n , j a T; n g in order to estimate bene ts of the sequential scheme 6.10 compared to the direct non-sequential test 6.8. Here E ; is the symbol of expectation for the given values of and . By Z p1y p y dy I = log p y 1 0 denote the Kullback-Leibler information number. In case when T a=I , , a=I , and a is su ciently large 0 is small the mean detection delay may be approximated by E ; f a T; n , j a T; n g a=I: 6.7.4. Performance of Algorithm 2. For Algorithm 2, reasonable performance may be expressed in terms of the frequency of false alarms Fr = 1=E 0 a E 0 a is the mean time between false alarms and the detection delay E ; f a , j a g. However, such characteristics as the probabilities of false alarm and target detection in the xed interval of the length T are also of interest. For Algorithm 2 with immediate renewal we have 6.40 Fr exp,a; E ; f a , j 5 6.38 pki a; 0 = Z a ,1 pki,1a; zfi z dz; k 1; a g a=I: This formula is approximate, since U,1 is a random variable belonging to the interval ,1; a. For = 1 it is exact. DETECTION ALGORITHMS AND TBD ARCHITECTURE FOR IRST 29 The rst formula in 6.40 follows from 6.26 and the second formula is valid when , a=I and a is su ciently large Fr is small. According to the rst inequality in 6.40 the threshold may be chosen as a = log1=Fr to guarantee the frequency of false alarms Fr Fr, where Fr is a speci ed level. However, it is even more reasonable to choose the threshold based on the speci ed probability of false alarm PfaT; a within the time interval of the length T at each current point in time n. Then the mean detection delay may be approximately evaluated by the second equality in 6.40 and the algorithm is optimal it minimizes the mean detection delay among all algorithms with this kind of constraints at least for su ciently large a. To obtain the probabilities of errors in the xed interval T we note that the statistic Un is again the Markov process which obeys the recursion 6.31 with the random initial condition U0 = u and with gUn,1 = Un,11lf0 Un, ag . The random variable u has the distribution G0y, which is nothing but a stationary distribution of the Markov process fUng under the hypothesis H0 target is absent. This distribution satis es the integral equation 1 G0x = Z F0 x , s1lf0 s ag dG0 s; where F0 y is the distribution function of n under H0. Thus, the false alarm probability is found as Z a PfaT; a = PfaT; a; x dG0x; where PfaT; a; x is the probability 6.29 in case when U0 = x, x a. This last probability is computed according to 6.33 with pT a; x = p0 a; x, which is de ned by the recursive relation T ,1 Z a p0 a; x = k ,1 p0 1a; zf0 z , x1lf0 k, , x ag dz; k 1 ; p0 a; x = 1 : 0 As we mentioned above, the analysis of the proposed algorithms may be done based on the renewal theory. This theory allows for more accurate performance evaluation as compared to the above formulas for the mean detection delay and frequency of false alarms. The detailed analysis of the proposed algorithms and their tuning will be done in the near future see also Section 6.8, relationships 6.41,6.42 for more accurate estimates in a Gaussian case. 6.8. TbD Model with Spatio-Temporal Matched Filter. Let us apply the above results to the case where the target moves along a deterministic trajectory. Then the spatio-temporal matched lter is the optimal tool for TbD. Since target's velocity and direction of motion are unknown, the bank of lters should be used and the detection is done in each channel independently. Assume also that after clutter suppression the residual clutter plus noise samples are i.i.d. and Gaussian, N 0; 2 see the model 5.9. Then the LLR k is X X k = 12 hx , xij , k , 1vZij k , 1 2 hx , xij , k , 1v2 : 2 i;j i;j In the above formula we implicitly used the assumption that the target moves with the constant given speed v and in a given direction i.e., in the particular velocity-angle" channel. This 30 CENTER FOR APPLIED MATHEMATICAL SCIENCES, USC formula is readily generalized for arbitrary but known dynamics. Note that here we use the notation Yk = fZij kg. According to 6.31 the statistic Unx is calculated as Unx = nx + gUn,1x; n = 1; 2; : : : Then the above algorithms are applied in each channel independently here we mean the channels that are related to di erent values of the initial position x. If the signal is exactly located in the given channel, then the LLRs fk xg are i.i.d. Gaussian random variables with parameters 1 E 1 k = ,E 0 k = 1 ; Var0 k = Var1 k = 4 ; 2 P where = 1 i;j hxij 2 is SNR. 2 Performance of the developed algorithms is evaluated by formulas obtained above in Section 6.7 with the densities 1 ' 2y + ; f y = 1 ' 2y , ; f y = 0 p p 1 p p where 'y is a standard normal density function, 2 'y = p1 exp , y2 : 2 In particular, for Algorithm 1 sliding window, the simplest estimates for probabilities of errors are given by see 6.39 ! 2a + T ; P , + 1; a 2a , , + 1 p PfaT; a 1 , p mis T , + 1 where y is a standard normal distribution function. Also, for the sequential Algorithm 2 with immediate renewal for su ciently large a and the following approximate estimates hold: Fr e,a ; 6.41 E ; f a , j a g 2 a + + C ; where 8 2a= 6.42 p = 2 exp , k=1 k , 2 k ; h P p , p i 1 ,1p = 1 + 4 , p 1 pk ' 2 k , 1 k , 1 k ; k=1 2 2 Pn : C = E 1 min f0; minn1 k=1 k g 2 0; ,1 : , P 2 1 1 1 These formulas were obtained by using the nonlinear renewal theory 60 . They improve the upper bounds for the frequency of false alarms and the estimates for the mean detection delay presented above in the general case see 6.26, 6.27, 6.40. DETECTION ALGORITHMS AND TBD ARCHITECTURE FOR IRST 31 This statistic, however, does not allow for a convenient recursive computation. To simplify the algorithm we use the trick", which was already used in Section 6.6 to build adaptive detection algorithms in general case. By x = x Z n denote the maximum likelihood estimate of the n n target location, i.e. Instead of making decisions on target presence in each channel we may apply the algorithm, which is based on the maximum likelihood statistic, Mn = max Unx; n = 1; 2; : : : x x n = arg max x n X k=1 k x: The decision statistic is de ned by the recursive formula Un = nx ,1 + gUn,1; n = 2; 3; : : : ; U1 = 0: n The stopping rule" is then de ned by a = inf fn 2 : Un ag: Thus, we use the one-stage delayed estimates of the signal location in the partial LLR n. The threshold is chosen according to the relation 6.27, which allows us to upper bound the false alarm rate. This algorithm is especially attractive when only a single target is expected. In the multi-target situation the previous algorithm decisions in independent channels is perhaps the best one can do. 6.9. TbD Model with Optimal Nonlinear Filter. In the context of TbD with optimal nonlinear ltering see Section 5 the estimate of the partial LLR n obtained on the basis of previous n , 1 observations is de ned by X X ^ ^ = hn,1xij Zij k , 1 hn,1xij 2 ; 6.43 n 2 i;j i;j ^ ^ ^ where hk xij = hXk , xij and Xk is the estimate of the signal location based on the k observations. The latter estimate is formed at the output of the optimal nonlinear lter that tracks the target location. We will consider two estimates the mean estimate and the maximum a posteriori estimate, which are de ned below in 7.3 and 7.4 see Section 7. Thus, the above results are applied by using the LLR estimate 6.43. To be precise, the adaptive CUSUM statistic is de ned by the recursion 6.44 Un = + gUn,1; n = 2; 3; : : : ; U1 = 0 n and the stopping rule" is 6.45 a = inf fn 2 : Un ag: In other words, again we have used the same trick the true LLR n at stage n is replaced with ^ ^ its estimate = nXn,1 obtained based on the previous n , 1 observations. Instead of Xn,1 n one may use the one step predicted extrapolated estimate Xn. 32 CENTER FOR APPLIED MATHEMATICAL SCIENCES, USC ~B = inf fn 2 : Rn B g; Rn = exp 1 + Rn,11lfR , n n The alternative algorithm see 6.19, 6.21 is of the form 6.46 6.47 where the estimate of the average likelihood ratio Rn satis es the recursion 1 Bg ; n = 2; 3; : : : ; R1 = 0 : The results of application of these algorithms will be discussed in Section 7.4. 7. Testing of the Developed Algorithms for IRST Data. Results of Experiments and Simulation In the examples considered below because of the uncertainty of possible behavior of a noncooperative target, the alignment of successive frames is extremely di cult. The frame alignment is done recursively along the maximum posterior distribution path on the basis of the optimal nonlinear ltering. Use of nonlinear algorithms is necessitated by the speci cs of the observation structure. Indeed, in TbD the observation function signal S k; r is highly sharp essentially a delta-function in small domain like 2 2 or 3 3 pixels. In addition, IR imaging sensors provide angle-only measurements. These types of measurements, especially in low SN+CR situations, practically rule out application of extended Kalman and similar lters. As already mentioned above, application of multi-dimensional matched lters and banks of velocity lters is possible only for constant speed targets with linear trajectories, which is not the case in the examples considered. At the same time the proposed ONF-based algorithm works perfectly well even for very low SN+CR up to ,6dB after clutter removal. Particularly, it may be seen from the gures presented below that the uncertainty is completely overcome after processing of several frames in these examples SN+CR ranges from ,1dB to ,7dB after simple pre-processing and clutter removal. The developed adaptive detection algorithms that use the estimates of target location constructed on the basis of ONF TbD work also perfectly well. Tracks of randomly appearing and disappearing targets are reliably detected up to ,7dB and the algorithms allow us to obtain low false alarm rate even in a heavy cluttered background. 7.1. Clutter Removal: Real IRST Data. The proposed nonparametric approach to clutter removal see Section 4.1 has certain advantages over conventional methods. This fact is illustrated by Figure 7 which presents the results of kernel smoothing of IR data. The `local mean removal' procedure used for clutter rejection in a series of papers by Reed et al. see 40 , 39 is equivalent to a kernel estimator with a product kernel K generated by uniform kernels K1x = K2 x = 0:5Ifjxj1g. It may be seen that the Fuller kernel see 26 provides superior smoothing performance compared to the uniform kernel. The reason is that our approach relies on the trade-o between a squared bias and a variance of estimators while the uniform kernel minimizes only the asymptotic variance of estimators and hence induces a large bias term. DETECTION ALGORITHMS AND TBD ARCHITECTURE FOR IRST (a) IR data: target plus clutter 1.5 20 1 40 60 80 0 100 20 40 60 80 100 -1 -0.5 0 0.5 1 0.5 Fuller triangular Epanechnikov uniform (b) Kernels in [-1,1] 33 (c) Residuals after uniform kernel (d) Residuals after Fuller kernel 20 40 60 80 100 20 40 60 80 100 20 40 60 80 100 20 40 60 80 100 Figure 7. NAVY IRSS sensor - the LAPTEX eld test. Original IR image a and the results of clutter removal by uniform kernel c and Fuller kernel d 7.2. Example 1: TbD of a Target Based on IRST Data. The rst example presented is a na ve model of an actively maneuvering air target which ies overhead over a stationary observer situated on the ground. The complete eld of view of the observer is modeled by a square 0; 1 0; 1 with the horizon line at y = 0:1. The target initially appears close to the horizon at a point uniformly distributed inside a rectangle 0:2; 0:8 0:1; 0:15 , then the horizontal coordinate evolves as a Wiener process with variance 2t, and the vertical coordinate as an exponent e t . Formally, the trajectory model is described by a stochastic di erential equation6 X0 Unif0:2; 0:8; Y0 Unif0:1; 0:15; dXt = dWt; = 0:05; dYt = Yt dt; = 0:6: The target was assumed to emit a 3 3 pixels signature which was observed with additive Gaussian white noise. More precisely, let be the size of a pixel in the observed image in both 6 Here the components of the vector r are denoted by x; y. Also X = X; Y . 34 CENTER FOR APPLIED MATHEMATICAL SCIENCES, USC horizontal and vertical directions. The target shape is given by hx; y = 1lfjxj3 =2g 1lfjyj3 =2g ; x; y 2 R 2 ; and the observation at time moment tk = k, = 0:05 with noise of standard deviation and SNR = 20 log 1 is Zij k = hXk , i , 1 ; Yk , j , 1 + ij k; ij k Norm0; 1 i.i.d. Examples of a target trajectory and IR data images are shown in Figure 8 and Figure 9. Several experimental results follow. For moderate noise = 1:232, SNR= ,1:8dB the simplest and fastest Haar basis works quite well see Figure 10. For more intense noise = 1:581, SNR= ,4:0dB, the algorithm utilizing Haar basis loses track occasionally but later recovers it see Figure 11. For borderline cases = 1:679, SNR= ,4:5dB, the Haar basis no longer provides good tracking, however smoother wavelet bases of the same resolution such as Coi et-1 still achieve good tracking accuracy see Figures 12 and 13; both simulations were run on exactly the same observations, and the di erence in the quality of tracking is quite obvious. With larger noise = 2:109, SNR= ,6:5dB, the optimal lter fails to track the target, no matter what basis is utilized. Trajectory 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 x y Figure 8. Typical trajectory of overhead target DETECTION ALGORITHMS AND TBD ARCHITECTURE FOR IRST 35 Observation at t=0 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 x Figure 9. Sample of a single frame of observation 4 2 0 -2 -4 -6 y 36 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 10 x CENTER FOR APPLIED MATHEMATICAL SCIENCES, USC State and estimates, coord=1 1 State Max est Mean est 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 30 40 50 60 time Confidence intervals, x 20 70 0 0 10 30 40 50 60 time Confidence intervals, y 20 70 y State Max est Mean est State and estimates, coord=2 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 10 20 30 40 time 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 10 20 30 40 time x 50 60 70 x y 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 10 20 30 40 time 50 60 70 Estimate errors Max est err Mean est err 50 60 70 Figure 10. ONF, Haar 64 64 basis, = 1:232, SNR= ,1:8dB DETECTION ALGORITHMS AND TBD ARCHITECTURE FOR IRST State and estimates, coord=1 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 10 30 40 50 60 time Confidence intervals, x 20 70 x y State Max est Mean est 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 10 30 40 50 60 time Confidence intervals, y 20 70 State Max est Mean est State and estimates, coord=2 37 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 10 20 30 40 time 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 10 20 30 40 time 50 60 70 x y 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 10 20 30 40 time 50 60 70 Estimate errors Max est err Mean est err x 50 60 70 Figure 11. ONF, Haar 64 64 basis, = 1:581, SNR= ,4:0dB 38 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 10 x CENTER FOR APPLIED MATHEMATICAL SCIENCES, USC State and estimates, coord=1 1 State Max est Mean est 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 30 40 50 time Confidence intervals, x 20 60 0 0 10 20 State Max est Mean est 30 40 50 time Confidence intervals, y 60 y State and estimates, coord=2 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 10 20 30 40 time 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 10 20 30 40 time x 50 60 x y 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 10 20 30 40 time 50 60 Estimate errors Max est err Mean est err 50 60 Figure 12. ONF, Haar 64 64 basis, = 1:679, SNR= ,4:5dB DETECTION ALGORITHMS AND TBD ARCHITECTURE FOR IRST State and estimates, coord=1 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 10 30 40 50 time Confidence intervals, x 20 60 x y State Max est Mean est 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 10 30 40 50 time Confidence intervals, y 20 60 State Max est Mean est State and estimates, coord=2 39 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 10 20 30 40 time 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 10 20 30 40 time x 50 60 x y 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 10 20 30 40 time 50 60 Estimate errors Max est err Mean est err 50 60 Figure 13. ONF, Coi et-1 64 64 basis, = 1:679, SNR= ,4:5dB 40 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 10 20 x CENTER FOR APPLIED MATHEMATICAL SCIENCES, USC State and estimates, coord=1 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 State Max est Mean est 30 40 50 60 time Confidence intervals, x 70 0.2 0.1 0 0 10 20 State Max est Mean est 30 40 50 60 time Confidence intervals, y 70 y State and estimates, coord=2 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 10 20 30 40 time 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 10 20 30 40 time x 50 60 70 x y 1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0 0 10 20 30 40 time 50 60 70 Estimate errors Max est err Mean est err 50 60 70 Figure 14. ONF, Haar 64 64 basis, = 2:109, SNR= ,6:5dB 7.3. Example 2: TbD of a Surface Skimming Missile. Assume that a stationary observation platform camera is at a xed height h2 above the sea level, and a target approaches with constant velocity v at constant height h1 h2 above the sea level. Denote the Earth radius R = 63; 750; 000m. Projection center the focal point of the camera is situated at the point with coordinates 0; 0; R + h1 . At time moment t = 0 the target appears at the horizon at angle clockwise from the direction of y,axis and moves at angle to the line of sight from the observer. Because of rotational invariance after projecting onto the observation cylinder, assume that = 0. The three dimensional coordinate system used next has the origin at the center of the Earth, the z,axis passing vertically through the observer, the y,axis passing horizontally through the observer in the direction of the point where the target becomes rst visible on the horizon, and x,axis oriented so that the xyz coordinate system is a Cartesian right hand coordinate system. x0 = 0; R q2Rh + h2 + q2Rh + h2 ; y0 = R + h 1 2 1 2 1 q 2 1 q z0 = R R h , R + h 2Rh1 + h2 2Rh2 + h2 ; 1 2 + 1 1 vt xt = R + h2 sin sin R + h ; vt , z cos 2sin vt ; yt = y0 cos R + h 0 R+h 2 vt + y cos sin vt 2 : zt = z0 cos R + h 0 R+h 2 2 DETECTION ALGORITHMS AND TBD ARCHITECTURE FOR IRST 41 After projecting onto an observation cylinder with radius r the focal length of the camera, the equation of the motion in two dimensional coordinates ; zp are as follows t = arctan xt ; ztp = p rzt 2 : yt x2 + yt t The observed signal was assumed to be a 3 3 pixels target with additive Gaussian white noise. More precisely, let 1 be the size of a pixel in the horizontal coordinate direction, and 2 be the pixel size in the vertical coordinate zp direction. The target shape is given by 7.1 hx1; x2 = 1lfjx j3 =2g 1lfjx j3 =2g ; x1 ; x2 2 R 2 ; and the observation at time moment tk = k with noise of standard deviation and SNR = 20 log 1 is ~ i = min + i , 1 1 ; p zjp = zmin + j , 1 2; ~ 7.2 ~ p ~ Zij k = h k , i ; zk , zjp + ij k; ij k Norm0; 1 i.i.d. 1 1 2 2 It is important to emphasize that the Gaussian model for the residual clutter and noise is used only in ONF and detection algorithm development. The algorithms were tested with the use of 42 CENTER FOR APPLIED MATHEMATICAL SCIENCES, USC the real IR background obtained from SPAWAR. Particularly, the gures in Section 7.4 below contain the results of testing for these real cluttered images i.e., fij kgk1 in 7.2 are given realistic frames with arti cially inserted targets whose motion and shape correspond to models 7.1, 7.2. 7.3.1. TbD of Several Targets. Ability to detect and track several independent targets was implemented by running the corresponding number of lters in parallel. The exact updating formulas are given below. i Translation operators. The operator T translates the argument of a function de ned on the state space along the trajectory of the model i. Let us denote the trajectory of the model i starting at the point x; y at time moment t = 0 by qit; x; y, and its inverse with respect to the variable t with variables x and y xed by tiq; x; y. Let the function to translate be f : 0; a 0; b f1g 0; a ,b; b f,1g f0g f0g f0g ! R. Several cases have to be split in order to describe the translation operator: , i 1 T f x; y; 1 = aqix;0; f 0; 0; 0; , , i T i f q ; x; y; 1 = f x; y; 1; , , i i Tf ,q ; x; y; 1 = f x; y; 1; , i i , i T f 0; 0; 0 = e,t f 0; 0; 0 + for x 2 0; a; y 2 0; qix; 0; ; for x 2 0; a; y 2 0; qix; b; ,; , for x 2 0; a; y 2 qix; b; ,; b; Tf q ; x; y; ,1 = f x; y; ,1; for x 2 0; a; y 2 0; b; Z a Z qi x;,b;, 0 ,b f x; y; ,1dydx: Updating formulas. The unnormalized ltering densities plt x; y; u corresponding to the model l are recomputed based on observations Z k = Zij k with time step according to formulas pl00; 0; 0 = 1; pl0x; y; u = 0; for x; y; u 6= 0; 0; 0; lp plk0; 0; 0 = Tnlk,1 0; 0; 0; X Hkx; y; u; Z k = exp 1 hx , xi ; y , yj ; uZij k , 21 2 , X 2 l plkx; y; u = Hkx; y; u; Z k Tplk,1 x; y; u; for x; y; u 6= 0; 0; 0: i;j , i;j hx , xi ; y , yj ; u2 o ; Next, the ltering densities are normalized as follows: = plk0; 0; 0 + l l kx; y; u = pkx; y; u ; nlk nlk Z Z plkx; y; 1 dx dy + Z Z plkx; y; ,1 dx dy; DETECTION ALGORITHMS AND TBD ARCHITECTURE FOR IRST 43 and two estimates of the position of the l-th target are computed. Namely, the mean-square-error estimate: Z Z Z Z l;mean = l x; y; 1 dx dy + l ^ Xk xk xkx; y; ,1 dx dy; Z Z Z Z 7.3 l;mean = l x; y; 1 dx dy + l ^ Yk yk ykx; y; ,1 dx dy; and the maximum a posteriori estimate: , l;max l;max l ^ ^ 7.4 Xk ; Yk = arg max kx; y; u: x;y Both estimates can be used in di erent situations. The rst estimate is theoretically optimal. However, it usually takes longer to converge to the real position of the target. The second estimate detects" the target much faster but is more prone to lose track due to noise. Precomputing., Computationally, the most expensive part of the algorithm was interpolating l the values of Tplt on the shifted along the trajectory and therefore no longer uniform grid back to the original uniform grid. To speed up the computations, o -line part of the algorithm computes a very sparse matrix A such that multiplying the vector of values on the non-uniform grid by A produces the vector of values on the uniform grid. Two interpolation algorithms were implemented: nearest neighbor interpolation then matrix A has only zero and one as elements, and each row has at most one non-zero, and linear interpolation using known values on a triangle surrounding a point on the new grid. Linear interpolation proved to be much more stable and precise in simulations. 7.3.2. Simulation results. The simulation producing results included below used two targets and the noise with standard deviation = 1:4 which corresponds to SNR = ,2:9dB. In Figure 17 t = 1s, the target is present but not yet detected. In Figure 18 t = 6s, a single sharp peak has formed, and it follows the target. In Figure 19 t = 11s, the second target just appears at the horizon, algorithm doesn't track it yet. In Figure 20 t = 16s, both targets were detected and are being tracked. Also we refer to enclosed les multiframe GIF le a6.gif or MS video le a6.avi to see the incoming observation frames Navy IRSS sensor the LAPTEX eld test with very dim synthetic target inserted. The second pair of les multiframe GIF le a5.gif or MS video le a5.avi shows the locations of two targets being tracked and the ltering density which develops prominent peaks that follow both targets. 44 CENTER FOR APPLIED MATHEMATICAL SCIENCES, USC 0 h1 (t) R (t) R h2 R Figure 15. Surface skimming missile trajectory derivation x 10 -1.8 -2 -2.2 zp t -2.4 -2.6 -2.8 -3 -3.2 -3.4 0 1 Angle=-10 Angle=7 Angle=15 -4 Trajectories 2 3 t 4 5 6 Figure 16. Surface skimming missile trajectories DETECTION ALGORITHMS AND TBD ARCHITECTURE FOR IRST 45 x 10 -1.5 -4 t=1.0s, SNR=-6.7dB Target: (1.45,-2.50e-004) 15 1 -2 10 0.8 0.6 0.4 5 -2.5 0.2 0 -1.8 -2 x 10 -3 0 1 2 3 4 5 6 -4 0 6 -2.2 -2.4 0 2 4 Figure 17. Surface skimming missile tracking, t = 1s x 10 -1.5 18 16 14 -2 12 10 8 6 -2.5 4 2 0 -2 -3 0 1 2 3 4 5 6 x 10 -4 -4 t=6.0s, SNR=-3.1dB Target: (1.49,-2.20e-004) 1 0.8 0.6 0.4 0.2 0 -1.8 -2 -2.2 -2.4 0 2 4 6 Figure 18. Surface skimming missile tracking, t = 6s 46 CENTER FOR APPLIED MATHEMATICAL SCIENCES, USC x 10 -1.5 -4 t=11.0s, SNR=-0.7dB 18 16 14 12 0.6 0.5 0.4 0.3 0.2 0.1 0 -1.8 Target: (1.55,-1.94e-004) (4.64,-2.50e-004) -2 10 8 6 -2.5 4 2 0 -2 -4 6 -2.2 -2.4 2 0 4 -2 x 10 -3 0 1 2 3 4 5 6 Figure 19. Surface skimming missile tracking, t = 11s x 10 -1.5 18 16 14 -2 12 10 0.3 8 0.2 6 -2.5 0.1 4 2 0 -2 -3 0 1 2 3 4 5 6 x 10 -4 -4 t=16.0s, SNR=-2.2dB Target: (1.66,-1.78e-004) (4.60,-2.20e-004) 0.6 0.5 0.4 0 -1.8 -2 -2.2 -2.4 0 2 4 6 Figure 20. Surface skimming missile tracking, t = 16s 7.4. Appearance Disappearance Detection of a Skimming Missile. In this section we present the results of application of the detection algorithms developed in Section 6. The model is the same as in Section 7.3. The results of TbD the estimates 7.3 and 7.4 of the signal location are used in the LLR estimate 6.43. In the fully sequential CUSUM algorithm, the adaptive CUSUM statistic Un is computed recursively according to 6.44 and then is compared to the threshold a see 6.45. In the sequential adaptive quasi-Bayesian algorithm the estimate of the logarithm of average likelihood ratio rn = log Rn is computed according to the recursive formula 6.47. It is then compared to the threshold log B see 6.46. Figure 21, Figure 22, and Figure 23 illustrate the behavior of the adaptive CUSUM statistic for di erent SN+CR SNR = ,0:25dB, SNR = ,3:62dB, and SNR = ,6:6dB, respectively. Theoretically 0 rn , Un log n. We have not observed any substantial di erence between these statistics in our experiments the di erence was always negligible. So we omit the plots of the statistic rn and show only the behavior of the adaptive CUSUM. It turns out that the mean estimate 7.3 provides better results than the maximum posterior estimate 7.4 and we plot only the CUSUM statistic that uses the former one. It is seen that most of the time the statistic Un is close to zero or negative when the target is absent. Sometimes, however, peaks arise. These peaks may be identi ed with signal presence. But they are typically short and may be easily distinguished from the peaks due to target. One possible method to discriminate between false alarms and true target is to make the decision on target presence if there are several subsequent exceedances say 3-4 of the threshold. Otherwise the decision on target absence is made i.e. a single exceedance is identi ed with false alarm. When the target appears the statistic rapidly increases while when it disappears, Un decreases. In our experiments we observed visible di erence in behavior of Un when the signal appears compared to the case of its absence up to the SNR ,6:6dB. In the pictures the rst target appears at time n = 1 and disappears at n = 28. The second target appears at time n = 39 and never disappears. Figure 24 shows the results of detection of track appearance and disappearance by Algorithm 3 for SNR = ,0:25dB and SNR = ,6:6dB. The detection of tracks occurs when the adaptive CUSUM Un exceeds the threshold a the upper one and track disappearance is detected when the statistic = L ,1 , Un exceeds the threshold b the lower one. When the decision on n n target disappearance is made, the statistic Un is renewed form zero, i.e. the detection algorithm is prepared to detect another target. It is seen that the algorithm is able to detect even very low SNR targets. The threshold a was chosen such that the false alarm rate 1 false alarm per minute i.e., per 60 frames was guaranteed. The algorithm never raised more than 1 false alarm and always detected targets with this choice in the analyzed images. To be precise, there are two targets in the pictures: the rst target appears at time n = 1 and disappears at n = 28 while the second one appears at time n = 39 and does not disappear. The algorithm detects the rst target with the delay about 20 seconds 20 frames for SNR = ,6:6dB and 4 , 6 frames for SNR = ,0:25dB. Then the fact of target's disappearance is detected with very small delay. Finally the second target is detected with the delay about 5 , 6 frames for SNR = ,6:6dB and 2 , 3 frames for SNR = ,0:25dB. DETECTION ALGORITHMS AND TBD ARCHITECTURE FOR IRST 47 48 4 CENTER FOR APPLIED MATHEMATICAL SCIENCES, USC Adaptive CUSUM statistic, no target 150 Un Adaptive CUSUM statistic, SNR=-0.25dB 2 0 -2 -4 -6 -8 -10 -12 0 10 20 30 40 50 60 -50 0 10 20 30 40 50 60 0 U n 100 50 Figure 21. Plots of the adaptive CUSUM. Left no target. Right 1st target appears at n = 1, disappears at n = 28; 2nd appears at n = 39, doesn't disappear Adaptive CUSUM statistic, no target 10 Un Adaptive CUSUM statistic, SNR=-3.62dB 60 Un 50 5 40 0 30 20 10 -10 0 -15 0 10 20 30 40 50 60 -10 0 10 20 30 40 50 60 -5 Figure 22. Plots of the adaptive CUSUM. Left no target. Right 1st target appears at n = 1, disappears at n = 28; 2nd appears at n = 39, doesn't disappear Adaptive CUSUM statistic, no target 10 Un Adaptive CUSUM statistic, SNR=-6.60dB 30 25 5 20 15 10 Un 0 -5 5 0 -5 -10 -15 0 10 20 30 40 50 60 -10 0 10 20 30 40 50 60 Figure 23. Plots of the adaptive CUSUM. Left no target. Right 1st target appears at n = 1, disappears at n = 28; 2nd appears at n = 39, doesn't disappear DETECTION ALGORITHMS AND TBD ARCHITECTURE FOR IRST 49 Algorithm 3, SNR=-0.25dB 120 100 80 60 40 20 0 -20 0 10 20 30 40 50 60 U n n Detection level Disappearance level Algorithm 3, SNR=-6.60dB 25 20 15 10 5 0 -5 -10 0 10 20 30 40 50 60 U n n Detection level Disappearance level Figure Detection of track appearance and disappearance by Algorithm 3. First target appears at n = 1 and disappears at n = 28; second target appears at n = 39 and doesn't disappear 24. 50 CENTER FOR APPLIED MATHEMATICAL SCIENCES, USC 8. Conclusion. Future Research We believe that real time optimal nonlinear ltering is an emerging information technology of fundamental importance. The ONF is an advance over conventional EKF techniques and multidimensional spatio-temporal matched lters potentially as profound as that of the EKF over , , lters. Nonlinear ltering techniques o er the promise of greatly enhanced performance in a number of DoD related applications: airborne, surface, and subsurface. We have developed an advanced ONF-based spatio-temporal tracking lter and adaptive track appearance disapperance detection algorithms that may be fully incorporated into a complete end-to-end IRST signal track processing suite. The developed technology will be of portable nature and can be incorporated later in various multisensor ATR systems utilizing EO IR, SAR, SAS, LIDAR, etc. In the future we plan to continue this work in the following key directions. 1. Further develop and tune spatio-temporal nonlinear semiparametric algorithms for clutter removal and jitter compensation. 2. Develop banks of ONF lters for TbD of multiple targets. 3. Improve and tune track detection and identi cation algorithms. 4. Test and validate the developed signal data track processing system for real IR data jointly with Space & Naval Warfare Systems Center, San Diego, CA. 5. Incorporate the developed system into a complete end-to-end IRST signal track processing suite jointly with Space & Naval Warfare Systems Center, San Diego, CA. 9. Acknowledgement We are grateful to Dr. John Barnett of Space & Naval Warfare Systems Center, San Diego, CA for very useful discussions and for providing the real IR data used in the analysis. DETECTION ALGORITHMS AND TBD ARCHITECTURE FOR IRST References 51 1 A. Aridgides, M. Fernandez, and D. Randolph, Adaptive three-dimensional spatio-temporal IR clutter suppression ltering techniques," SPIE Proc., Vol. 1305, 1990. 2 J. Arnold and H. Pasternack, Detection and tracking of low-observable targets through dynamic programming," SPIE Proc., Vol. 1305, 1990. 3 P. Athanas and A. Abbott, Real-time image processing on a custom computing platform," IEEE Computer, pp. 16-24, 1995. 4 J.T. Barnett, B.D. Billard, and C. Lee, Nonlinear morphological processors for point-target detection versus an adaptive linear spatial lter: a performance comparison," Proceedings of SPIE Symposia on Aerospace Sensing, Vol. 1954, pp. 12-24, Orlando, FL, 1993. 5 J.T. Barnett, Statistical analysis of median subtraction ltering with application to point target detection in infrared backgrounds," Proc. SPIE Symposium on Infrared Signal Processing, Vol. 1050, pp. 10-18, Los Angeles, CA, 1989. 6 C. Barlow and S. Blackman, New Bayesian track-before-detect design and performance study," SPIE Proc., Vol. 3373, 1998. 7 Y. Barniv, Dynamic programming solution for detecting dim moving targets," IEEE Trans. AES, Vol. 21, 1985. 8 Y. Bar-Shalom and T.E. Fortmann, Tracking and Data Association Orlando: Academic Press, 1988. 9 S. Blackman, R. Dempster, and T. Broida, Multiple hypothesis track con rmation for infrared surveillance systems," IEEE Trans. AES, Vol. 29, pp. 810-823, 1993. 10 S.S. Blackman, Low SNR detection and tracking," Technical Report, Hughes Aircraft Co., 1994. 11 S.S. Blackman, Multiple Target Tracking with Radar Applications. Artech House, 1986. 12 H. Blom, An e cient lter for abruptly changing systems," Proc. 23rd IEEE Conf. Decision Contr., pp. 656658, 1984. 13 H.A.P. Blom, et al., Design of a multisensor tracking system for advanced air tra c control," In: MultitargetMultisensor Tracking, Vol. 2 ed. Y. Bar-Shalom, pp. 31-63, 1992. 14 H. Blom and Y. Bar-Shalom, The interacting multiple model algorithm for systems with Markovian switching coe cients," IEEE Trans. Autom. Contr., Vol. 33, pp. 780-783, 1988. 15 G.C. Demos, R.A. Ribas, T.J. Broida, and S.S. Blackman, Applications of MHT to dim moving targets," Proc. SPIE: Signal and Data Processing of Small Targets, Vol. 1305, pp. 297-309, 1990. 16 Y. Chan, J. Plant, and J. Bottomley, A Kalman tracker with a simple input estimator," IEEE Trans. AES, Vol. 18, pp. 235-241, 1982. 17 L.J. Cutrona, Synthetic Aperture Radar," In: Radar Handbook Ed. M. Skolnik. New York: MacGraw-Hill, 1990. 18 D.L. Donoho, I.M. Johnstone, G. Kerkyacharian, and D. Picard, Wavelet Shrinkage: Asymptopia ?" with discussion J. Roy. Stat. Soc. B, Vol. 57, no. 1, pp. 301-369, 1995. 19 D.L. Donoho, Asymptotic minimax risk for sup-norm loss: Solution via optimal recovery," Probab. Theory Rel. Fields, Vol. 99, pp. 145-170, 1994. 20 W. Feller, An Introduction to Probability Theory and Its Applications, Vol. 2. New York: John Wiley & Sons, 1966. 21 M. Fernandez, A. Aridgides, and D. Bray, Detecting and tracking low-observable targets using IR," SPIE Proc., Vol. 1305, 1990. 22 M.A. Girshik and H. Rubin, A Bayes approach to a quality control model," Ann. Math. Statist., Vol. 23, pp. 114-125, 1952. 23 W. Hardle, Applied Nonparametric Regression. Cambridge: Cambridge University Press, 1990. 24 A.H. Jazwinski, Stochastic Processes and Filtering Theory. New York: Academic Press, 1970. 25 S. Kligys and B.L. Rozovskii, Optimal nonlinear ltering with distributed observation," CAMS Preprint, USC, 1997. 26 S.L. Leonov, On the solution of an optimal recovery problem and its applications in nonparametric regression," Mathematical Methods of Statistics, Vol. 6, no. 4, pp. 476-490, 1997. 27 O.V. Lepski, E. Mammen, and V.G. Spokoiny, Optimal spatial adaptation to inhomogeneous smoothness: an approach based on kernel estimates with variable bandwidth selectors," Ann. Stat., Vol. 25, pp. 929-947, 1997. 28 R. Lindgren and L. Taylor, Bayesian eld tracking," SPIE Proc., Vol. 1954, 1993. 29 G. Lorden, Procedures for reacting to a change in distribution," Ann. Math. Statist., Vol. 42, pp. 1897-1908, 1971. 52 CENTER FOR APPLIED MATHEMATICAL SCIENCES, USC 30 S.V. Lototsky, R. Mikulevicius, and B.L. Rozovskii, Nonlinear ltering revisited: a spectral approach," SIAM J. Control Optim., Vol. 35, no. 2, 1997, pp. 435-461. 31 S.V. Lototsky and B.L. Rozovskii, Recursive nonlinear lter for a continuous-discrete time model: separation of parameters and observations," IEEE Trans. Autom. Control, Vol. 43, no. 8, pp. 1154-1158, 1998. 32 S.V. Lototsky, C. Rao, and B.L. Rozovskii, Fast nonlinear lter for continuous-discrete time multiplication models," Proc. 35th Conf. Decision and Control, Kobe, Japan, 1996, Omnipress, Madison, Wisconsin, Vol. 4, pp. 4060-4064. 33 S.V. Lototsky, C. Rao, and B.L. Rozovskii, Method of optimal stochastic ltering for tracking objects with possibly nonlinear dynamics," U.S. Patent pending. 34 N. Mohanty, Computer tracking of moving point targets in space," IEEE Trans. PAMI, Vol. 3, 1981. 35 E.S. Page, A test for a change in a parameter occurring at an unknown point," Biometrika, Vol. 42, no. 4, pp. 523-527, 1955. 36 N. Platt, S. Hammel, J. Trahan, and H. Rivera, Mirages in the marine boundary layer comparison of experiment with model," Proc. IRIS Passive Sensors, Monterey, CA, Vol. 2, pp. 195-210, 1996. 37 N. Platt and S.M. Hammel, Pattern formation in driven coupled map lattices," Physica A, Vol. 239, pp. 296-303, 1997. 38 M. Pollak, Optimal detection of a change in distribution," Ann. Statist., Vol. 13, pp. 206-227, 1986. 39 I.S. Reed, R.M. Gagliardi, and H.M. Shao, Application of three-dimensional ltering to moving target detection," IEEE Trans. AES, Vol. 19, 1983, pp. 898-904. 40 I.S. Reed and X. Yu, Adaptive multiple-band CFAR detection of an optical pattern with unknown spectral distribution," IEEE Transactions on Acoustics, Speech, and Signal Processing, Vol. 38, pp. 1760-1770, 1990. 41 V.G. Repin, Detecting a signal with unknown moments of appearance and disappearance," Problems Inform. Transmission, Vol. 27, no. 1, pp. 61-72. 1991. 42 S.W. Roberts, A comparison of some control chart procedures," Technometrics, Vol. 8, no. 3, pp. 411-430, 1966. 43 B.L. Rozovskii, Stochastic Evolution Systems: Linear Theory and Application to Nonlinear Filtering, MIA, Kluwer Academic Publishers, 1990. 44 J. Sanders, A method for determining lter spacing in assumed velocity lter banks," IEEE Trans. AES, Vol. 29, 1993. 45 A.N. Shiryaev, On optimum methods in quickest detection problems," Theory Probab. Applications, Vol. 8, pp. 22-46, 1963. 46 W.W. Shrader and V. Gregeres-Hansen, MTI Radar," In: Radar Handbook Ed. M. Skolnik. New York: MacGraw-Hill, 1990. 47 D. Siegmund, Sequential Analysis. New York: Springer-Verlag, 1985. 48 P.F. Singer and D.M. Sasaki, The heavy tailed distribution of a common CFAR detector," SPIE Proc., Vol. 2561, pp. 124-140, 1995. 49 P. Singer, Performance analysis of a velocity lter bank," SPIE Proc., Vol. 3163, 1997. 50 F.M. Staudaher, Airborne MTI," In: Radar Handbook Ed. M. Skolnik. New York: MacGraw-Hill, 1990. 51 A. Stocker and P. Jensen, Algorithms and architectures for implementing large velocity lter banks," SPIE Proc., Vol. 1481, 1991. 52 A.G. Tartakovsky and S. Blackman, Asymptotically optimal sequential tests and application to target detection and discrimination," 30th Annual Conference on Information Sciences and Systems, Princeton, NJ, 1996. 53 A.G. Tartakovsky, Sequential Methods in the Theory of Information Systems. Moscow: Radio i Svyaz', 1991. 54 A.G. Tartakovskii, E ciency of the generalized Neyman-Pearson test for detecting changes in a multichannel system," Problems of Information Transmission, Vol. 28, pp. 341-350, 1992. 55 A.G. Tartakovskii, Asymptotic properties of CUSUM and Shiryaev's procedures for detecting a change in a nonhomogeneous Gaussian process," Mathematical Methods of Statistics, Vol. 4, no. 4, pp. 389-404, 1995. 56 A.G. Tartakovsky, Minimal time detection algorithms and applications to ight systems," School of Engineering and Applied Science, Flight Systems Research Center, Technical Report 2-FSRC-93, 1993. 57 A.G. Tartakovsky, Detection of signals with random moments of appearance and disappearance," Problems of Information Transmission, Vol. 24, pp. 115-124, 1988. 58 S. Tonissen and R. Evans, Performance of dynamic programming techniques for track-before-detect," IEEE Trans. AES, Vol. 32, 1996. 59 P. Wei, J. Zeidler, and W. Ku, Analysis of multiframe target detection using pixel statistics," IEEE Trans. AES, Vol. 31, 1995. 60 M. Woodroofe, Nonlinear Renewal Theory in Sequential Analysis. Philadelphia: SIAM, 1982. ...
View Full Document

This note was uploaded on 06/10/2008 for the course ANST 503 taught by Professor Gold during the Winter '06 term at USC.

Ask a homework question - tutors are online