Institutionen för systemteknik
Department of Electrical Engineering
Examensarbete
Missile approach warning using multi-spectral
imagery
Examensarbete utfört i Bildbehandling vid Tekniska högskolan i Linköping
av
Erika Emilsson, Hannes Holm Ovrén LiTH-ISY-EX--10/4329--SE
Linköping 2010
Department of Electrical Engineering Linköpings tekniska högskola
Linköpings universitet Linköpings universitet
Missile approach warning using multi-spectral
imagery
Examensarbete utfört i Bildbehandling
vid Tekniska högskolan i Linköping
av
Erika Emilsson, Hannes Holm Ovrén LiTH-ISY-EX--10/4329--SE
Handledare: Klas Nordberg
isy, Linköpings universitet
Jörgen Ahlberg
FOI
Examinator: Klas Nordberg
isy, Linköpings universitet
Avdelning, Institution
Division, Department
Computer Vision Laboratory Department of Electrical Engineering Linköpings universitet
SE-581 83 Linköping, Sweden
Datum Date 2010-06-04 Språk Language Svenska/Swedish Engelska/English Rapporttyp Report category Licentiatavhandling Examensarbete C-uppsats D-uppsats Övrig rapport
URL för elektronisk version
http://www.control.isy.liu.se http://urn.kb.se/resolve?urn=urn:nbn:se:liu:diva-57147 ISBN — ISRN LiTH-ISY-EX--10/4329--SE
Serietitel och serienummer
Title of series, numbering
ISSN
—
Titel
Title
Missilvarning med hjälp av multispektrala bilder Missile approach warning using multi-spectral imagery
Författare
Author
Erika Emilsson, Hannes Holm Ovrén
Sammanfattning
Abstract
Man portable air defence systems, MANPADS, pose a big threat to civilian and military aircraft. This thesis aims to find methods that could be used in a missile approach warning system based on infrared cameras.
The two main tasks of the completed system are to classify the type of missile, and also to estimate its position and velocity from a sequence of images.
The classification is based on hidden Markov models, one-class classifiers, and multi-class classifiers.
Position and velocity estimation uses a model of the observed intensity as a function of real intensity, image coordinates, distance and missile orientation. The estimation is made by an extended Kalman filter.
We show that fast classification of missiles based on radiometric data and a hidden Markov model is possible and works well, although more data would be needed to verify the results.
Estimating the position and velocity works fairly well if the initial parameters are known. Unfortunately, some of these parameters can not be computed using the available sensor data.
Nyckelord
Keywords missile approach warning, classification, target tracking, hidden markov models, kalman filtering, threshold model, multispectral, infrared
Abstract
Man portable air defence systems, MANPADS, pose a big threat to civilian and military aircraft. This thesis aims to find methods that could be used in a missile approach warning system based on infrared cameras.
The two main tasks of the completed system are to classify the type of missile, and also to estimate its position and velocity from a sequence of images.
The classification is based on hidden Markov models, one-class classifiers, and multi-class classifiers.
Position and velocity estimation uses a model of the observed intensity as a function of real intensity, image coordinates, distance and missile orientation. The estimation is made by an extended Kalman filter.
We show that fast classification of missiles based on radiometric data and a hidden Markov model is possible and works well, although more data would be needed to verify the results.
Estimating the position and velocity works fairly well if the initial parameters are known. Unfortunately, some of these parameters can not be computed using the available sensor data.
Sammanfattning
Bärbara luftvärnsrobotsystem är ett stort hot mot civil och militär flygtrafik. Detta examensarbete ämnar att utveckla metoder som kan användas i ett missil-varningssytem baserat på infraröda kameror.
De två huvuduppgifterna i det färdiga systemet är att klassificera missilens typ, samt estimera missilens position och hastighet.
Klassificeringen baseras på dolda Markovmodeller, enklass-klassificerare och flerklass-klassificerare.
Positions- och hastighetsestimeringen utnyttjar en modell av observerad in-tensitet, som en funktion av sann inin-tensitet, bildkoordinater, avstånd till missilen samt dess riktning. Denna estimering görs av ett utökat Kalmanfilter.
Vi visar att en snabb klassificering kan göras bra, baserad på radiometriska data och en dold Markovmodell, men att mer data behövs för att fullständigt verifiera resultaten.
Att estimera position och hastighet fungerar relativt bra, förutsatt att de ini-tiala parametrarna är kända. Tyvärr kan vissa av parametrarna inte beräknas med den tillgängliga sensordatan.
Acknowledgments
First we would like to thank our supervisor Jörgen Ahlberg for his great support and ability to provide us with new ideas when such were required.
At FOI we would also like to thank Jonas Nygårds for his help with the 3D tracking, and Thomas Svensson for providing us with all the necessary data.
We would also like to thank our supervisor at ISY, Klas Nordberg.
Contents
1 Introduction 1 1.1 Background . . . 1 1.2 Problem formulation . . . 2 1.3 Available data . . . 3 1.4 Scope . . . 4 1.5 Notation . . . 4 2 System overview 5 2.1 Overview . . . 5 2.2 The detector . . . 6 2.3 The 2D tracker . . . 6 2.4 Classifier . . . 6 2.5 3D tracker . . . 7 3 Classification 9 3.1 Introduction . . . 9 3.2 Methods . . . 93.2.1 Hidden Markov models . . . 10
3.2.2 The Viterbi algorithm . . . 10
3.2.3 Threshold model . . . 11
3.2.4 Classification methods . . . 12
3.3 Overview . . . 14
3.4 Rough missile classification . . . 15
3.4.1 Temporal features . . . 15
3.4.2 Missile model . . . 16
3.4.3 Measuring RMC confidence . . . 19
3.5 Detailed missile classification . . . 20
3.5.1 Feature selection . . . 22
3.5.2 Rejecting obvious non-missiles . . . 22
3.5.3 Determine if an object is a known missile . . . 29
3.5.4 Determining the missile type . . . 30
3.5.5 Cascading classification trees . . . 31
3.5.6 When only the boost phase is observable . . . 33
3.5.7 The resulting radiometric classifier . . . 37 ix
3.5.8 Adding feedback from 3D tracker . . . 38
3.6 Results . . . 38
3.6.1 When the eject phase is observable . . . 38
3.6.2 When only the boost phase is observable . . . 39
3.6.3 Feedback from the 3D tracker . . . 40
3.7 Discussion . . . 41
3.7.1 Results . . . 41
3.7.2 Further work . . . 42
4 3D tracking 45 4.1 Introduction . . . 45
4.2 The missile’s thermal intensity . . . 46
4.3 The missile’s thermal intensity function . . . 51
4.4 Camera matrix calculation . . . 52
4.5 The Kalman filter . . . 54
4.6 The extended Kalman filter, EKF . . . 56
4.7 Jacobian calculation . . . 58 4.8 Initial values . . . 59 4.9 Summary . . . 59 4.10 Result . . . 60 4.11 Discussion . . . 67 5 Results 71 5.1 Introduction . . . 71 5.2 Test setup . . . 71 5.3 Results . . . 73 5.3.1 Time complexity . . . 73 6 Discussion 75 6.1 Introduction . . . 75 6.2 Results . . . 75 6.3 Further work . . . 76 Bibliography 79 A Derivations 81 A.1 EKF . . . 81
Chapter 1
Introduction
1.1
Background
Aircraft operating in war zones and other danger areas constantly face the threat of being shot down. Man Portable Air Defence Systems, MANPADS, are missiles that can be carried and fired by humans on the ground. These kind of weapon sys-tems have been in use for a long time and are used by most military forces around the world. Because of its long history and availability on the black market for a relatively cheap price, it is also an affordable weapon for terrorist organizations.
During the time in which MANPADS type weapons have existed, they have been used in a number of terrorist attacks against civilian targets. A recent exam-ple is the attack with two MANPADS against an Israeli Boeing 747-airliner taking of from Mombasa airport, Kenya, in November 2002. In that case the missiles missed their target.
MANPADS are often slower and have a limited range compared to other anti aircraft weapon systems which are stationary or vehicle mounted. Aircraft flying fast or high enough are relatively safe. However helicopters and fixed-wing aircraft that are taking off or landing are vulnerable to MANPADS [4].
For aircraft to operate safely in danger areas, some kind of protection system is preferable. One part of the protection system is a missile approach warning system, which is used to warn the operator about missiles heading toward the aircraft. If the warning can be generated quickly enough, countermeasures can be deployed and hopefully save the aircraft.
To further keep the cost down, the warning could take the missile path into account. In that way, the system would only warn if the missile is heading towards the aircraft, and it would make it possible to let several aircraft “share” a system carried by a single aircraft.
1.2
Problem formulation
To deal with the MANPADS problem, we have studied methods for a missile approach warning system. We have implemented a system that should accept an infrared image sequence as input, and generate a warning for each missile that poses as threat.
To achieve this, the following problems must be solved:
1. Detecting missiles. Given images from infrared cameras, which can be multi spectral, detect hotspots with radiometric properties which are similar to the radiometric properties of missiles. The system should detect both MANPADS and other types of missiles.
2. Classifying the missile. Different missile types have different characteris-tics, so we should be able to classify the missile given motion and radiometric data. This is useful because if we know the exact type of missile we can en-able the correct counter measure.
The primary goal is to classify MANPADS missiles when the missile is observed from the moment it is fired.
The secondary goal is to classify MANPADS missiles when the missile is observed in a later phase of its life span.
3. Estimating the missile’s position and velocity in 3D-space. As soon as we have high enough probability that a detected object is a missile, we want to estimate its position. Our inputs are assumed to be a stream of two-dimensional images, i.e. range data or stereo imagery will not be available. However, the observed thermal intensity of the missile will vary with the distance to the missile and orientation of the missile. The intensity will be higher if the missile is flying away from the sensor, and thus has its engine pointed towards it compared to if it is flying towards the sensor as the engine exhaust plume is occluded by the body of the missile. Using this kind of observations together with a known ego motion and image coordinates, an estimate of the missile’s location in three-dimensional space should be possible to compute.
The main advantage gained by estimating the missile’s position and velocity is that it makes it possible to calculate what the missile is aiming for. The aircraft at risk thus get the possibility to enable countermeasures.
If multiple aircraft are flying together, it could be enough to equip one of them with the system and have it automatically warn the actually targeted aircraft when a missile has been fired against it.
4. Efficient algorithms. Since the time from firing the missile to impact is only a few seconds, it is vital that the system can operate in real time so that an approaching missile is detected early enough to enable countermeasures. The real time constraint means that the time complexity of all algorithms are of great importance.
1.3 Available data 3
1.3
Available data
The data used in this thesis are of military nature and thus classified. Therefore no numerical data will be presented in this report. Plots will be given without numerical scale markings and will either be synthetically created data or real data. The data is the result of firing tests conducted by FOI in collaboration with other defence agencies. Missiles where fired from a hill and locked on to a station-ary target on the ground.
The position was measured using an instrument known as a theodolite which gives us a 3D-position in RT90 2.5G (X,Y coordinates) and RH70 (height) repre-sentation. The exact accuracy of the measurements are unknown to us, but are assumed to be within ±1 meter.
Image data was collected by several infrared cameras positioned about 4 km from the firing position. We have worked with data from only one of the cameras, the MultimIR, see Figure 1.1. The camera has a rotating filter wheel with four different filters. The frame rate is high enough such that four consecutive frames can be assumed to represent one 4D-image representing the same static scene. Each pixel has an intensity for each of the four spectral bands. The camera has a size of 384 x 288 pixels and a narrow field of view of less than 5 degrees both horizontally and vertically. The filters give us the following spectral bands:
1. 1.55 − 1.75µm (SWIR, Short Wave IR) 2. 2.05 − 2.45µm (SWIR)
3. 3.45 − 4.15µm (MWIR, Middle Wave IR) 4. 4.55 − 5.50µm (MWIR)
Figure 1.1. The MultimIR camera that was used to acquire the data
In total we had 12 different sequences of which 8 were a missile of type A, and 4 of type B. Five sequences also contained a flare being launched at approximately the same time as the missile.
1.4
Scope
Missile approach warning is a hard problem and building a complete system is not within the scope of this thesis. First, this thesis concerns only the signal processing, i.e. not hardware, sensors, choice of wavelengths, etc. Second we make a number of assumptions about the problem to be solved and leave some problems to be solved by others.
Images are rectified and calibrated. When the system is live, it will probably consist of not one, but several sensors placed to give a 360◦ field of view. The assumed input to the system presented in this thesis is a single image which means that the multiple images will have to be combined together. The ego motion is known. We assume that the vehicle on which the system is
attached has the necessary sensors (GPS, INS1, ...) to determine its current position and rotation in the world.
Only MANPADS are considered. MANPADS type weapons are not the only threat to aircraft. Air-to-air missiles are deemed to be unlikely where the system will operate. Surface-to-surface missiles can of course be used against low altitude aircraft flying slowly, but we will not consider them. MANPADS type weapons give a very specific signature and is probably the largest threat in the operating environment of the system. Therefore we will focus on detecting them.
1.5
Notation
x, c Scalars
w, x Vectors
A, P Matrices
f (t), g(x) Scalar valued functions
f (t), g(x) Vector valued functions
x ∼ N (µ, Σ) The vector x is normal distributed with mean µ and covariance matrix Σ
˙
x dx
dt
S Set
P (x|y) Conditional probability of event x given event y
Chapter 2
System overview
2.1
Overview
The system we have developed consists of four subsystems, each solving their own part of the problem without knowing very much about the other subsystems. The subsystems are connected such that the output from one subsystem can be the input of another subsystem as illustrated in Figure 2.1.
2D tracker Classifier 3D tracker Detector CLASS POSITION VELOCITY Images
Figure 2.1. System overview
2.2
The detector
The detector takes an image as input, and gives a list of detections as output. Each detection in that list is a point that is of interest to the system. A detector that performed well enough was already available when we started our thesis. Thus we will only mention briefly how the detector works.
The detector is an anomaly filter which finds pixels that are brighter than its surroundings. How much brighter it has to be is determined by a user defined threshold. When multiple spectral bands are available, the anomaly filter is applied to all bands separately, and the final detection image is the maximum value of all bands, for each pixel.
The input images are assumed to be well calibrated and show measured tem-perature of an object. The measured temtem-perature is not necessarily the actual temperature of the object since atmospheric dampening (atmospheric transmit-tance and distransmit-tance) will affect the measurement. We will exploit this property in the 3D tracker.
Since the infrared cameras have a low resolution, a missile will be detected as 1-9 pixels in the resulting infrared image. Using high resolution infrared cameras would be very expensive, without giving much benefit to the system.
2.3
The 2D tracker
The purpose of the 2D tracker is to track objects in the input image sequence. The 2D tracker works only in the 2D image plane (as opposed to the 3D tracker that attempts to find the objects location in 3D). The 2D tracker takes detections from the detector and creates tracks. The goal is to have one track per object in the scene.
For each track, classification is performed by the classifier. If the track (or rather the object that the track represents) is classified as a missile the 3D tracker will try to estimate its position and velocity.
The 2D tracker, just like the detector, was already working when we started this project and has not been modified by us.
2.4
Classifier
The classifier analyses a 2D track and attempts to find what class the object belongs to. The class can either be a generic missile, a specific type of missile or some other kind of object type. The classification uses the same infrared images that were the input to the detector to extract radiometric information about an object.
The results from the classifier are passed to the 3D tracker to help it make a better estimation. The position and velocity calculated by the 3D tracker is then sent back to the classifier. This information is used to determine if the speed is consistent with what we expect from a missile.
2.5 3D tracker 7
2.5
3D tracker
Given the image coordinates and observed intensity of the missile, we apply a model which estimates the position and velocity of the missile. The reason for this is to find out if the missile is heading towards the sensor, at which point the system should give a warning.
Data from the classifier about the exact type of missile can be used to apply a more detailed model to give a better estimate of position and velocity.
The position and velocity is sent back to the classifier to refine the classification based on the speed of the object.
Chapter 3
Classification
3.1
Introduction
As soon as the detector detects an anomaly (something hotter than its surround-ings) we want to know, as soon as possible, if this anomaly is a potential missile. Moreover, we want to know the type of missile to make sure that the system can respond with the correct action. The environment where the system will operate is likely to contain a large number of hotspots which can be the result of burning vehicles or houses, gunfire or even by flying over an ordinary chimney. Therefore it is crucial that the classification is very quick so that we can almost immediately dismiss most hotspots as non-missile objects.
As noted in Section 1.2, the system must perform classification in a limited time span. This means it would be infeasible to run a complex classification algorithm on all detected objects.
To improve speed we divide the classification into two tasks. The first task is to find out whether or not the detected object is a potential missile. We call this step rough missile classification. The second step is to find the exact missile type and is named detailed missile classification in this report.
Section 3.3 gives an overview of the subsystem, while Sections 3.4-3.5 describes how it is constructed using the general methods described in Section 3.2.
The performance of the classifier is presented in Section 3.6.
Finally, in Section 3.7, we discuss the results in relation to what is expected from a system that is ready for use, what factors might have affected the result, and what we think can be improved in the future.
3.2
Methods
The following is a brief introduction to the methods we use for our classification. The sections on hidden Markov models and the Viterbi algorithm are based on the descriptions found in [2].
3.2.1
Hidden Markov models
A Markov process is a stochastic process which consists of a number of possible states. The Markov process is time discrete and follows the Markov property. This property implies that the conditional probability of an event is dependent on only some fixed number of past events. In our case we will consider only the single past event. If S = {S1, S2, ..., Sn} is a set of n states and Q0...τ = (q0, ..., qτ)
is the sequence of states for time 0...τ such that qt ∈ S, ∀t, then the Markov
property (of order 1) can be expressed as P (qt+1|qt, qt−1, ..., q0) = P (qt+1|qt).
Since the number of states is assumed to be finite we can define a transition probability matrix A where ai,j = P (qt+1 = Sj|qt = Si) which describes the
transition probabilities between all states. Finally we need to know the initial state distribution Π = (π0, ..., πn) where πi= P (q0= Si).
A hidden Markov model is a Markov process where we cannot directly observe the state sequence Q0...τ. Instead we have some sequence of observations O0...τ =
(o0, ..., oτ) which depend on the hidden state sequence. We define an observation
probability bj,t = P (ot|qt= Sj) which gives us the probability of an observation
given a specific state. We will use λ = (A, bj, Π) as a shorter notation for hidden
Markov models.
Example 3.1 shows a very simple example of using a hidden Markov model.
Example 3.1: Simple HMM example
We want to know whether the weather is sunny or rainy where our co-worker Bob lives. We can not observe the states (Sunny or Rainy weather), but given the weather the probability that Bob brings an umbrella to work changes (bj).
Observing whether Bob is carrying an umbrella (O) combined with some model of the weather where he lives (A) is enough to make certain guesses as to the current weather at Bob’s place.
Hidden Markov models can be used for a wide variety of tasks, but we are only interested in finding out whether or not a given observation sequence (the intensity) fits a specific model good enough. To do this we utilize the Viterbi algorithm.
3.2.2
The Viterbi algorithm
TheViterbi algorithm takes a sequence of observations O0...τ = (o0, ..., oτ) as input
and gives the most likely state sequence Q∗0...τ = (q∗τ, ..., q∗0) as output, given the hidden Markov model λ.
This can be written as
Q∗0...τ = arg max
Q0,...,Qτ
P (Q0, ..., Qτ|λ, O0...τ) (3.1)
3.2 Methods 11
Q∗0...τ = arg max
Q0,...,Qτ
P (Q0, ..., Qτ, O0...τ|λ) (3.2)
To make optimisation simpler we define
δt(i) = max q0,...,qt−1
P (q0, ..., qt−1, qt= Si, o0, ..., ot) (3.3)
which is initiated by setting
δ0(i) = P (q0= Si, o0) = πibi,0 (3.4)
From 3.3 and 3.4 follows a recursive way of calculating δt(i): δt(i) = max
j (δt−1(j)aj,i)bi,t (3.5)
To be able to find the most likely path we define a back-tracking function Φt(i) = arg max
j δt−1(j)aj,i (3.6)
which is the most probable state j to have transitioned from to state i at time
t.
We can then find the most likely path by setting q∗τ = arg maxiδτ(i) and q∗t = Φt+1(q∗t+1), t < τ .
3.2.3
Threshold model
The Viterbi algorithm always return the best matching path, even if the path is a very bad match to the observations. A threshold model[7] can be used to deal with this problem.
The idea is to create a new ergodic model, the threshold model, which contains the states from all other models, but where the transition probabilities have been modified. The original models are together referred to as the target model in the following description of the method.
Definition 3.1 (Ergodic model) A hidden Markov model, or rather its under-lying Markov chain is ergodic if all states are
• aperiodic, and • positive recurrent
Definition 3.2 (Aperiodicity in a Markov chain) A state in a Markov chain is said to be periodic if the return time from a state to itself must occur in multiples of k. If k = 1 then the state is aperiodic which means we can return to the state at any time.
Definition 3.3 (Recurrence in a Markov chain) A state in a Markov chain is said to be recurrent if the state can be reachable from itself within some finite time. If the expected return time is finite then the state is said to be positive recurrent.
The ergodic model have the same states as the target model, but the transition probabilities are modified so that all states can be reachable from all states. Definition 3.4 (Transition probabilities in the threshold model)
ai,j=
1 − ai,i
n − 1 , ∀{(j, i) : i 6= j}
where n = number of states. The self transition probability ai,i will remain the
same, while the other transition probabilities will get equal probabilities.
The behaviour we can expect from the threshold model is that it will have a higher probability than the target model for sequences which do not match the target model. At any time τ we can calculate
δthreshold,τ= max
i δτ(i), ∀Si∈ Threshold model (3.7) δtarget,τ = max
i δτ(i), ∀Si∈ Target model (3.8)
The threshold model can be used as a confidence measure for a path calculated by the Viterbi algorithm. Since the threshold model will be more general (allowing transitions between all states) it will yield a higher probability for most observation sequences. But when the observations match the target model well, the probability for the target model will be higher than the threshold model, since the transition probabilities between states will be higher. We can thus use δtarget,τ ≥ δthreshold,τ
as a measure of whether the observations match the target model, or not.
3.2.4
Classification methods
There are two types of classifiers
Multi-class classifiers labels data as belonging to one of two or more classes. This is probably the most common type of classification.
One-class classifiers labels the data as either belonging to a specific class, or not. They are also called an anomaly detectors since they can be seen as modelling what is normal, and then determine when something does not fit the model.
One-class classifiers can be divided into three categories [9]:
Density methods The probability density of the class is estimated, and to de-termine if a new data point matches the class you simply threshold that density. If it falls below the threshold it is an anomaly.
3.2 Methods 13
Boundary methods Instead of estimating the density and a threshold, one can try to directly estimate the decision border. Three examples of boundary methods are k-centers, nearest neighbour and support vector data descrip-tion.
Reconstruction methods With this method a model is fitted to the target class. The model has a description which is often of a lower dimensionality than the data. The data to be classified is then first projected to the description space and then reconstructed into the feature space. The classification measure-ment is then the distance between the original point and the reconstructed point. Reconstruction methods are often based on methods which are not originally developed for one-class classification, such as principal component analysis and k-means classifiers.
We will now describe two classification methods we have used. Gaussian density classifier
A Gaussian density classifier is a density method which tries to fit a (multivariate) Gaussian distribution to the data. The mean µ and covariance matrix C, is first estimated from the training data. The classification is then performed by calculat-ing the Mahalanobis distance between the test data and the Gaussian distribution. By using the Mahalanobis distance instead of the Euclidan distance the covariance is taken into account and the distance can be measured in standard deviations.
The data is assumed to belong to the class in question if the distance is lower than some threshold.
Fisher’s linear discriminant
Fisher’s linear discriminant (FLD) [6] is a multi-class classifier that is very similar to linear discriminant analysis (LDA). However, unlike LDA it makes no assump-tion about the distribuassump-tion of the data. LDA assumes the data to be a Gaussian distribution which is not necessary with FLD.
FLD can be used to find a hyperplane that separates two classes which are linearly separable.
The trick with FLD is finding a linear transformation, w, from the feature space to a one-dimensional space where the separation between the classes is maximal. Maximal separation is in this case a transformation such that each class has a small variance, while the distance beween the classes is large.
Given two classes, C1, C2, that differ with regards to their means, µ1 and µ2,
the Fisher criterion for making this discrimination is defined as
J (w) =w TC bw wTC tw (3.9) where Cb= (µ2− µ1)(µ2− µ1)T (3.10)
is the between-class covariance, and Ct= X n∈C1 (xn− µ1)(xn− µ1)T+ X n∈C2 (xn− µ2)(xn− µ2)T (3.11)
is the total within-class covariance.
The w that maximizes Equation 3.9 must meet the following criteria
Cbw = λCtw (3.12)
Because Cbw is always in the same direction as (µ2− µ1) we can rewrite
Equation 3.12 as
w = C−1t (µ1− µ2) (3.13)
When we have applied the linear transformation w we find a threshold, c, such that ( wTx − c < 0 if x ∈ C1 wTx − c ≥ 0 if x ∈ C 2 (3.14) The threshold c is chosen as
c = wT µ1+ µ2
2
(3.15) which is the projected mean of the means of the two classes.
So far we have spoken of FLD as a linear transformation, while we began by defining the problem as finding the hyperplane that separates two classes. This is the same thing. w and c together define the hyperplane
wTx − c = 0 (3.16)
which separates the two classes. The vector w is the normal to this hyperplane.
3.3
Overview
The classifier is divided into two components:
The Rough Missile Classifier (RMC) is responsible for quickly determining if the object (represented by a track in the tracker) is missile-like or not. Missile-likeness is calculated by using the temporal thermal intensity char-acteristics of the object. The RMC is described in detail in Section 3.4. The Detailed Missile Classifier (DMC) takes over if the RMC has decided
that the object is missile-like enough. The task of the DMC is then to first further investigate if the object is likely to be a missile, and then to decide which exact type of missile it is. The DMC is described in more detail in Section 3.5.
3.4 Rough missile classification 15
An overview of the classifier subsystem can be seen in Figure 3.1. We can see that the classifier takes a track and image stream as input and gives a class label as output. The class is decided by either the RMC or the combined results of the RMC and DMC.
RMC DMC
Images
CLASS Track
Figure 3.1. Overview of the classifier subsystem
3.4
Rough missile classification
The rough missile classifier, RMC, is responsible for determining if the observed object is likely to be a missile or not. Classifying a missile as a non-missile can be fatal so it is important that the RMC is very liberal in labelling objects as potential missiles. Its rate of false negatives (missiles classified as non-missiles) must be as close to zero as possible while we can allow the rate of false positives (non-missiles classified as missiles) to be higher.
3.4.1
Temporal features
We are primarily interested in discovering missiles launched from MANPADS type weapons. This kind of weapon works in the following way, as illustrated by Figure 3.2:
1. When the operator presses the trigger the missile engine is started, which ejects the missile from the missile tube. We call this the eject phase.
1This image is a collage created from http://commons.wikimedia.org/wiki/File:FIM-92A_
missile_launch.jpg and http://commons.wikimedia.org/wiki/File:Blackhawk.jpg which are both available under public domain
Figure 3.2. MANPADS example1. 1: Eject, 2: Delay, 3: Boost, 4: Sustain
2. To not burn the operator, the engine shuts off just before exiting the tube and is in free flight for a short time period. We call this the delay phase. 3. The engine reignites in a powerful burst which accelerates the missile towards
its target. We call this the boost phase.
4. Finally the engine stays on until it either hits the target or it runs out of fuel. We call this the sustain phase. Unless the target is very far away the engine will stay on until the target has been hit. Even if the engine has shut off, the missile is still visible in the infrared domain since both the engine and the missile body are hot. The latter as a result from air compression. This behaviour produces a very characteristic thermal intensity curve as can be seen in Figure 3.3.
To find signals that match this we use hidden Markov models (see Section 3.2.1) because these are very suited to model such temporal characteristics.
3.4.2
Missile model
As shown in the above section, MANPADS go through a number of phases when fired. Since a missile inarguably have to be in one and only one of these phases at a given time it seems natural to use them as the states in a hidden Markov model. Since the curve could be described as “high, zero or low, then high again” we will use the following as observations:
Definition 3.5 (Observation function) The observation otis defined as
ot=
(
high if It≥ Ithreshold
low if It≤ Ithreshold
(3.17)
where It is some observed intensity which could be either the detection value
from the detector, the intensity from one specific spectral band, or some combi-nation of multiple spectral bands. Ithreshold is a threshold which decides if an
3.4 Rough missile classification 17 0 0 Time Intensity Boost Sustain Eject Delay
Figure 3.3. Example of observed missile intensity. The intensity gets higher in the end of the sustain phase because the missile is getting closer to the observer, and the observed heat is inversely proportional to the (squared) distance.
varies by distance (d) and angle (α) of the missile it should really be defined as
Ithreshold= f (d, α). The thermal intensity function model we use is described in
Section 4.3. The angle is not available to us, but the range is acquirable given a homography.
A homography is a transformation from one plane to another plane which in our case is a mapping between image coordinates on the image plane to a point on the ground plane. If we assume that the observed object is on the ground we can get its approximate position in 3D using the mentioned homography. This assumes that the system knows about the position and orientation of the sensor platform in relation to the ground.
The exact value of the threshold is not critical since the variations in intensity between the different phases are large. Therefore it should be good enough to define the threshold as Ithreshold= d12K, where K is some constant.
Since both the boost and sustain phases expect a high intensity, these can be merged into a single boost/sustain state giving us the initial state space
{Eject, Delay, BoostSustain}
However, there is a possibility that the first observations are low instead of high which would not match the model. This happens if the detector works with a lower threshold than the threshold of the RMC. If the shooter is positioned next to an object that is quite hot (like a vehicle) the detector might start tracking that object, and thereby inserting intensity values into the RMC, before the missile has been launched. Ideally both the detector and RMC could use thresholds such that this never happens. To address this we add a “garbage state” that takes care
Boost/Sustain Delay
Eject Garbage
Figure 3.4. State space
of initial low intensity values. This gives us the following state space to model a MANPADS type missile.
Definition 3.6 (Missile model state space)
S = {Garbage, Eject, Delay, BoostSustain}
A graphical representation of the model can be found in Figure 3.4.
The final part of the model is to determine the initial state probabilities, Π, and the transition matrix A. If we have a large quantity of data these parameters can be estimated. However, since this is a rough classifier we do not need a perfect model. Simply deciding on some probable values should be enough.
The parameter values we are using are the following:
Π = (0.25, 0.25, 0.25, 0.25) (3.18) A = .4 .6 0 0 0 .4 .6 0 0 0 .9 .1 0 0 0 1 (3.19)
We give all states an equal initial state probability because we can not know when the observations begin. It might happen that we miss the eject phase and only start observations in the boost phase. And, although unlikely, there is a possibility that the first observations are from the garbage or delay states. If that happens we still want to be able to say that it might be a missile. Because of the limited data available we can not make any more detailed estimation of initial probabilities. We think that gathering this kind of data would also be very hard as it would require flying around in the target operating environment while shots are fired and then count and see in which phase the observation started. Simulations could probably be made, but since the initial probability does not matter very much this is likely a waste of time. The transition probabilities play a much larger role since there are lots of observations.
Both the garbage and eject phases should be quite short while the delay phase is (relatively) long and this is expressed as giving the self transition probability for the delay state very high compared to the other. The boost phase is the last
3.4 Rough missile classification 19
state in the model and once the model enters this state there is no more states to go to. Therefore the self transition probability is 1.
If we run the Viterbi algorithm on intensity from a missile with the suggested model and parameters above, it will return the most probable path through the states. The problem is that it will always return a path even if the data does not fit the model very well. Some kind of confidence measure must be applied to the returned path.
Example 3.2: Viterbi algorithm with HMM missile model
Consider the intensity shown in Figure 3.5. It does not match the Eject-Delay-Boost-Sustain signature we expect of a missile. Therefore it is very likely not a missile. However, if we apply the Viterbi algorithm on this data it will return a path which claims that the first high-rise is part of the garbage state with the eject phase starting at the second high-rise. That is obviously not a good result.
0 10 20 30 40 50 60
Time
Intensity
Figure 3.5. Example of intensity curve which does not match the Eject-Delay-Boost-Sustain model
3.4.3
Measuring RMC confidence
As example 3.2 shows we need a mechanism to detect whether the detected path is a good fit or not. The Viterbi algorithm provides us with a δτ(i) (See Equation 3.3)
value which is the probability of state Si at time τ . This means that maxiδτ(i) is
the probability of the most likely path at time τ . Unfortunately we can not simply apply a threshold to this value since it depends on
1. The length of the sequence. Increasing τ will give an (exponentially) lower
δ since at each time step the previous δ value is multiplied with two
2. The observation probabilities, which will depend on which states are part of the sequence.
Instead of trying to find a way to calculate a threshold we will use a threshold model (see Section 3.2.3). Since our target model is a left-to-right model and requires very specific transitions to be valid, the threshold model will yield a lower probability than the target model for missile like observation sequences, but a higher probability for all other sequences.
The resulting state space with the missile model and the corresponding thresh-old model can be seen in Figure 3.6.
Using the threshold model we can define a simple confidence measure as Definition 3.7 (Confidence measure for RMC)
Cτ=
(
1 if δtarget,τ ≥ δthreshold,τ
0 otherwise
where δtarget,τ is the probability of the most probable path in the original missile
model, and δthreshold,τ is the probability of the most probable path in the threshold
model, at time τ .
If the threshold model matches the data better than the target (missile) model, then we can safely assume that the observations did not match the target model well at all.
The reason why we include δtarget,τ = δthreshold,τ is because of the case when
only the boost phase is observed. In this case we have no temporal signature to use, and all observations will be ot= high which means the most probable path
will consist of only one state. Since the self transition probabilities have not been changed between the target model and the threshold model and the observation probabilities are equal we get δtarget,τ = δthreshold,τ.
For the threshold model to work, the parameters used for A and Π (Equations 3.19 and 3.18) can not be chosen arbitrary. The threshold model is based on the missile model, and the two are because of this tightly linked to each other. If the parameters are chosen badly, there is a risk that some sequence is mislabeled. Fortunately, there is a low number of sequences that have to be tested. We began by setting the parameters to values that seemed reasonable. When applied to a few test sequences we could quickly tweak the parameters to get the desired effect.
3.5
Detailed missile classification
When the RMC has decided that an object is likely to be a missile we need to establish which kind of missile it is. That is the job of the detailed missile classifier, DMC. The object passed to the DMC can belong to one of three classes
• A missile of known type (the system has a model for it).
3.5 Detailed missile classification 21 Boost/Sustain Delay Eject Garbage TM: Boost/Sustain TM: Delay TM: Eject TM: Garbage
• Not a missile at all.
As already noted in Section 3.4, the RMC is very undiscriminating in which objects it accepts as possible missiles. This is because the general rule is to warn one time too much rather than one too few. Therefore the object passed to the DMC is not guaranteed to be an actual missile.
For this reason we will build the DMC to work in the following three steps 1. Reject objects that can not be a missile.
2. Decide the likelihood that the missile matches any of the known missile types. 3. Decide the type as one of a set of known missile types.
3.5.1
Feature selection
To do any kind of classification we must know what features to look at: what observations can we make that are likely to be different between different types of objects. Our feature selection is based on what we can extract from the sequences we introduced in Section 1.3.
Below is a list of some of the possible things to look at. • Missile phase duration (eject, delay, boost)
• Missile phase intensity (eject, boost, sustain), one per spectral band avail-able.
• Missile phase intensity variance (eject, boost, sustain), one per spectral band available.
• Ratio between spectral bands. Combinations of different spectral bands and/or phases.
To determine which features are good and which are not we investigate the distribution of each feature for different type of objects. The only object types in the available scenes are missiles (type A and B) and flares. Figure 3.7 to 3.13 show the distributions of each object type for each feature, either as a normal distribution or histogram. For the eject phase features only the missiles are shown since flares do not have an eject phase, but rather a continuous burn. We skip the sustain phase since our goal is to detect missiles in an early stage. When the missile reach sustain it will probably be to late to classify and give an appropriate warning anyway. The variance is only given for the boost phase since the eject phase is too short to give enough data to be useful.
3.5.2
Rejecting obvious non-missiles
Our first step is to remove any object that the RMC has classified as a missile but have some characteristic that simply makes it impossible for it to be a missile. The RMC bases its classification on an aggressively quantized intensity curve, and
3.5 Detailed missile classification 23 0 0.5 1 1.5 2 2.5 3 3.5 4 Eject length Time Count Missile A Missile B
Figure 3.7. Eject phase length.
0 0.5 1 1.5 2 2.5 3 3.5 4 Delay length Time Count Missile A Missile B
1 2
3 4
Missile A Missile B
Eject intensity, per band
3.5 Detailed missile classification 25 1 / 2 1 / 3 1 / 4 2 / 3 2 / 4 3 / 4 Missile A Missile B
Eject phase, ratio between bands
1 2
3 4
Missile A Missile B Flare
Boost intensity, per band
3.5 Detailed missile classification 27 1 / 2 1 / 3 1 / 4 2 / 3 2 / 4 3 / 4 Missile A Missile B Flare
Boost phase, ratio between bands
1 2
3 4
Missile A Missile B Flare
Boost variance, per band
3.5 Detailed missile classification 29
there are obviously a lot of cases that will slip through. Given the list above we should be able to set up hard thresholds which will reject the worst cases. The rules we will use are:
1. Eject duration: TEject,min≤ TEject ≤ TEject,max
2. Delay duration: TDelay,min≤ TDelay≤ TDelay,max
The reason why we use only phase duration data and not intensities is simply that these constraints will be valid for all MANPADS type weapons. You can not have a very long eject phase or very short delay phase since the operator would get the exhaust plume in the face. And if the delay phase is too long there is a risk that the missile hits the ground (or objects on the ground) before the engine boost pushes it up high in the air. If we base the decision on intensity there is always a risk that a new missile type with highly deviating intensity characteristics could fool the system.
3.5.3
Determine if an object is a known missile
What we could do is build a classifier which labels the object into the classes C = {Non-missile, Missile1, ..., MissileM} where M is the total number of known
missile types. The problem is that it is very hard to construct the Non-missile class since this requires knowledge and training data for all (or at least a lot of) possible objects that are not missiles. Since the number of different type of objects in the operating environment of the system is likely to be huge, this approach does not seem to be feasible.
Instead we opt to use a so called one-class classifier. A one-class classifier is different from other classifiers in that it does not discriminate between multiple classes. Instead it tries to fit a model to the target class so that we can decide if a new data point is a member of this class, or not. By creating a one-class classifier from all known missile types, we can then decide if a new object is a (known) missile or not a missile.
One-class classifiers are sometimes referred to as anomaly detectors. This is because you can look at the problem as defining what is normal and then deciding if the new data deviates from what is expected: an anomaly [1].
The available data was divided into a training set and a validation set. We then estimated the mean and variance of the training data and used that to create a Gaussian density model (see Section 3.2.4) of missile type objects. Since we have chosen a 2D feature space (more on this Section 3.5.4), we can easily visualize the classifier, as can be seen in Figure 3.14.
This classifier obviously encompasses all our missile data points, so there is really no need to find a “better” classifier. Until we have more data about missiles and other objects we have no information about how a classifier could ever be “better”.
The one-class classifier is trained on known missile types, but we want to in-crease the chance that it is general enough to catch unknown missiles as well.
Therefore we want to set the decision border further away than 3 standard devi-ations (since that would encompass the trained classes perfectly). The decision border was arbitrarily set to 3.5 standard deviations since we want the classifier to catch more than we trained it for, but not risk matching too many non-missile objects.
Delay length (s)
Boost intensity ratio between band 2 and 4
Gaussian density fitted to missile data
Figure 3.14. Gaussian density fitted to missile data. The dashed line is one standard deviation. The solid line is the decision border.
3.5.4
Determining the missile type
When the one-class classifier has decided that the object is indeed a missile we want to know its exact type. Since execution time is an issue the simpler the classifier the better. Ideally we would like to find a feature space where the two missile types are linearly separable. By linearly separable we mean that you can divide the feature space by a hyperplane where points on one side of the plane belongs to the first class while the points on the other side of the plane belong to the second class. Looking at the feature distribution plots we can see that delay length and boost ratio between band two and four both independently gives us the possibility to easily discriminate between the two classes (see Figure 3.8 and 3.12).
By using these features we create a two-dimensional feature space which can be seen in Figure 3.15. In this feature space there is no problem to draw a line to separate the two clusters.
3.5 Detailed missile classification 31
If two features construct a feature space where the classes are linearly separable, then they would obviously be linearly separable if we created a 10-dimensional feature space of all available features. So why only use two features? One reason is that more data will require more processing time, which we are trying to minimize. Although it can be argued that this effect is probably quite low since the classifier is executed only once per track, and not continuously. The primary reason for using a smaller feature space is that it requires less data to estimate the parameters. Let us use the Gaussian distribution as an example:
Example 3.3: Cost of parameter estimation (Gaussian distribution) Let N be the number of dimensions of a multivariate Gaussian distribution. The distribution is parameterized by two entities: the mean vector, µ, and covariance matrix, C. µ has N dimensions which means N parameters. C has N2elements,
but since it is symmetric the number of unique elements (and also the number of parameters), is N22−N + N . In total this results in N22−N + 2N parameters to estimate.
For the parameter estimation to be possible we need k data points so that the criteria in Equation 3.20 is met.
kN ≥ N
2− N
2 + 2N ⇒ k ≥
N − 1
2 + 2 (3.20)
For two dimensions we get k = 2.5 and for ten dimensions we get k = 6.5.
Note that the numbers in example 3.3 are the minimal number of required data points. To get a good estimate a lot more data is required. There should be a substantial difference in using about twice as many data points as required (10-dimensional case) or five times as much data as required (2-dimensional case). To find a line that separates the two classes we will use Fisher’s linear discrim-inant (see Section 3.2.4).
The resulting discriminating line for our two missile classes can be seen in Figure 3.16.
When the Gaussian density classifier and Fisher’s linear discriminant classifier are combined, we get the classification regions shown in Figure 3.17.
3.5.5
Cascading classification trees
To make the classification as easy to use as possible we would like to have a way of combining different classifiers into a final classification result. From a more general classifier we might want to find a more specific class, and then an even more specific subclass. To facilitate this we use a concept we call a cascading classification tree. As the name suggests it is built around the data structure known as a tree. In this tree the nodes are different classifiers and the edges are classification results (classes). The data point to classify is passed to the classifier in the root node, and depending on the output class the data point can then be passed on to any
Delay length
Boost ratio 2 / 4
Missile A Missile B
Figure 3.15. Feature space where the two missile classes is clearly linearly separable.
Delay length
Boost ratio 2 / 4
Missile A Missile B
3.5 Detailed missile classification 33
Figure 3.17. Regions showing the partitioning of the feature space. Red is missile type A, blue is missile type B, and the white background is non-missile objects.
of the child nodes (classifiers) of the root node, where the procedure is repeated. Sort of like falling through a waterfall, which is why we use the term cascading. An example of a classifier tree can be seen in Figure 3.18.
With the help of this it is easy to set up and test classifiers on data. Just create a new classifier object and insert it in the classification tree at the correct place. In our Matlab implementation the classification tree is constructed of classifiers which are subclasses of an abstract classifier base class. Each type of classifier (Gaussian density, k-centers, etc) is its own subclass which take a different set of parameters (e.g. mean and covariance for Gaussian density). The tree is then defined in an XML file which defines all the separate classifier nodes with their respective parameters and relationships.
To take a concrete example, the set of classifiers described in the Section 3.5.4 is represented by the classification tree shown in Figure 3.19.
3.5.6
When only the boost phase is observable
Our secondary classification goal is to handle the case when we do not observe the full missile sequence, but begin observing somewhere in the boost phase. In this case we do not have access to the delay and eject phase so all features can only depend on data available from the boost phase. That means that the classifier discussed in the previous sections can not be used. Looking back at the feature distribution plots of the boost phase (Figures 3.11, 3.12 and 3.13) it is hard to instantly see one or two features that would provide us with what we need, which is:
1. a way to discriminate between missile objects and other type of objects, and 2. a way to discriminate between missile type A and B.
Root Subclassifier C A (C) C1 C2 Subclassifier B B2-2 B1 (B) Subclassifier B2 (B2) B2-1 B2-3
Figure 3.18. A general classifier tree with the output classes {A, B1, B2-1, B2-2, B2-3, C1, C2}.
3.5 Detailed missile classification 35 Missile? (Gaussian) A or B? (Fisher’s) Non-missile (Missile) Missile A Missile B
Figure 3.19. Classifier tree used when we can see the full missile sequence.
The distributions all overlap in some way which makes it more difficult to choose good features. But, when plotting the feature distributions we are plotting their one dimensional projection from a multi dimensional feature space. Just because they are not linearly separable in this one dimensional space does not mean they are not linearly separable in a higher dimensional space. A simple example can be seen in Figure 3.20.
1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 9
Figure 3.20. Two classes that are linearly separated in two dimensions, but not if they are projected onto the X- and Y-axis separately.
Unfortunately, time does not permit us to invest as much time in this clas-sification problem that might be needed. Instead of researching, implementing and testing a large number of classifiers we will instead use two Matlab tool-boxes called PRTools2and DDTools3. PRTools contain a lot of different classifiers
2PRTools is available for both academic and commercial use from http://www.prtools.org 3DDTools is an extension of PRTools and is available from http://ict.ewi.tudelft.nl/
Classifier FN FP
K-means 33% 0%
K-centers 0% 50%
PCA 0% 0%
Table 3.1. Performance of the best classifier of a total of 10,000 trained classifiers with different training and validation data.
while DDTools extends it by adding a number of one-class classifiers. With these toolboxes it is very easy to try different classifiers to see if the problem seems solvable.
As features we will use all boost ratio combinations and all boost variances, which gives us a 10-dimensional feature space, which is unfortunately very difficult to visualize.
We start by using the same two types of classifiers as in the full sequence case: One Gaussian density one-class classifier, and one Fisher’s linear discriminant classifier.
Unfortunately, trying to estimate the Gaussian density is very hard. This is because the ratio Ndata
Nfeatures is very low, as we previously explained in example 3.3.
We now have 10 dimensions and a maximum of 12 data points (but since we need to keep a few points for validation, it is actually less) compared to the previous two dimensions. It quickly proves to not be an option.
Since the Gaussian density model did not work, we need a new one-class clas-sifier to discriminate between missiles and other objects. Density methods has the drawback of requiring quite a lot of data, so we will instead try boundary and re-construction methods. For no special reason we choose to try k-means, k-centers, and principal component analysis (PCA). For a description on how these three methods work, see [9].
Fisher’s linear discriminant works much better. Using half the data for training and the other half for validation we can directly get a classifier that discriminates between missiles of type A and B. Just like before the Fisher’s linear discriminant defines a hyperplane that divides the feature space between the two classes. In this case the hyperplane is nine-dimensional instead of one-dimensional as in the full sequence case.
We set up a test where we randomly divide the data into a training set and validation set (of equal size). Each classifier is then trained with the training set and then validated against the validation data set. This was repeated 10,000 times and for each run the number of false negatives (FN) and false positives (FP) were recorded. The best classifier of each type was chosen as the one with the lowest number of false negatives. The results can be seen in Table 3.1. We can see that we were able to create a PCA classifier which produced no errors on the supplied validation set, so we will use that to classify objects as missiles or non-missiles.
So far we have only looked at missiles, as that is the important part. But our sequences also include flare objects, and it would be nice, but not necessary to classify them as well. If nothing else it shows a little more interesting classification tree (see Figure 3.21). To find the flares we used the same type of PCA classifier
3.5 Detailed missile classification 37 Missile? (PCA) A or B? (Fisher’s) (Non-missile) (Missile) Missile A Missile B Flare? (PCA) Flare Non-missile
Figure 3.21. Classifier tree used when only the boost phase is observed. Both missiles and flares can be classified.
as before. But this time it was trained to discriminate between flare objects and everything else.
3.5.7
The resulting radiometric classifier
In the previous section we have described in detail two classifier trees: one for when the full sequence is observed, and one for when only the boost phase is observed. Combining them is simple: if the RMC reports a path which contains all phases we chose the former, and if the path only contains the boost phase then the latter is chosen.
In more detail, the classification is performed in the following order: 1. A new track is created by the tracker.
2. For each track update, the intensity is used as observation for the RMC. 3. At a pre-determined time (t ≈ 1s) the RMC decides if the intensity function
has characteristics that makes it a possible missile (either full sequence, or only boost phase). If not, the object gets a non-missile class label.
4. Feature data is extracted from the RMC (phase lengths, intensities, ...) and passed to the DMC.
5. The DMC begins by rejecting objects that can not be missiles. They get a non-missile class label.
6. Depending on if the full sequence or only the boost phase was observed, the corresponding classifier tree is used to classify the data. The output is, in our case, one of the following classes: Cmissile= {Maybe-missile, MissileA, MissileB}.
The maybe-missile class is used to separate objects that we are entirely sure are not missiles, from objects that we do not think are missiles but could be.
Note that the radiometric classification is performed only once per track. The RMC is executed several times, but no classification decision is taken before the pre-determined classification time.
3.5.8
Adding feedback from 3D tracker
If the result of the radiometric classifier is that the object could possibly be a missile (its class is in Cmissile) , then the 3D tracker will try to determine its position and
velocity (see more in Chapter 4). The position is not of much use, but the velocity will give us some useful information about the object. For instance, if the speed of the object is very slow or very fast compared to what we expect from a missile then we can either reject our initial guess (“It looked like a missile but it does not move like one”) or reinforce the probability of our previous radiometric classification (“It looks like a missile, and moves like a missile”).
Depending on the performance of the 3D tracker we can also “promote” an object labelled as maybe-missile to missile, if its motion is consistent with a missile. In our implementation we set limits that require missile objects to have a speed |v| ∈ [20, 1000] m/s, which should be enough since MANPADS missiles have a maximum speed around 500-600 m/s.
3.6
Results
3.6.1
When the eject phase is observable
We applied the radiometric classification on the 12 test sequences of which eight contained a missile of type A and four contained a missile of type B. Five of these sequences also include a flare being launched (see Section 1.3 for a description of the data).
In all sequences the eject phase of the missile could be observed. The feedback from the 3D tracker was turned off when collecting these results since we want to show the result from the radiometric classifier alone.
The classification trees used was of the following type:
• Full-sequence classification tree. Features used are delay length and boost intensity ratio of bands 2 and 4.
1. Root: Gaussian density estimation. Maximal distance in standard de-viations was σmax= 3.5. Determines if the object is a missile or not.
2. A or B: Fisher’s linear discriminant. Determines if an object determined by Root to be a missile is of type A or B.
• Only-boost classification tree. Features used are all boost intensity ratios (total of six) and all boost variance ratios (total of four). A 10-dimensional feature space.
1. Root: One-class PCA classifier trained with the DDTools Matlab module.
3.6 Results 39
2. A or B: Fisher’s linear discriminant. Determines if an object determined by Root to be a missile is of type A or B.
As a measure of result we used the following metrics:
FNmissile Number of false missile negatives: Missiles that were classified as not
a missile.
FPmissile Number of false missile positives: Non-missile objects that were
classi-fied as a missile.
For each sequence we also recorded if there were any other especially interesting classifications (correct and incorrect labels).
Sequence FNmissile FPmissile Comments
A2 0 0
A3 0 0
A4 0 0
A5 0 0 Flare classified as flare
A7 0 0 Flare classified as non-missile
A8 0 0 Flare classified as non-missile
A9 0 0 Flare classified as non-missile. Two non-flares
were classified as flares.
A10 0 0 One non-flare was classified as flare
B1 0 0
B2 0 0 Two non-flares were classified as flares
B3 0 0
B4 0 0
Table 3.2. Classification results, σmax= 3.5.
As the results in Table 3.2 show the classifier works very well. All missiles were indeed detected as missiles and there were no false alarms. Actually the only classification errors are the non-missile objects that were classified as flares. The flares that were classified as non-missile objects are correct (because a flare is not a missile), although it would have been even better if they were classified as flares. A second attempt was made when we lowered the σmaxof the gaussian density
classifier from 3.5 to 3 standard deviations. In that case the missile object in sequence A8 was not classified as a missile, it was however very close to the decision border (about 3.2 standard deviations from the mean).
3.6.2
When only the boost phase is observable
Although the main objective was to detect MANPADS missiles, given that we observe the whole event, it is of interest to see how the classification performs when we start observing the missile in the boost phase. Since the classifier is trained on only about the first half second or so, we have tested against different
