In my first passive radar post I complained about having to extract a bunch of passive radar measurements manually by clicking on blobs in a sequence of several hundred images. I have fortunately now found a better way of doing things, as is demonstrated in this video:
This tracker uses a Kalman filter, which is a well known algorithm that has been extensively applied for radar and other motion tracking applications. I will not cover the full mathematical details of the Kalman filter here (wikipedia does it pretty concisely), however I will go over the dynamical model that underlies the Kalman filtering algorithm for the case of passive radar. Next I’ll show some experimental results and compare the performance of Kalman filters with different parameters. All the associated code can be found here.
Kalman Filter Basics
Kalman filters deal with the problem of measurement noise by maintaining an internal dynamical model of the system under observation. At each time step, the filter uses this model to predict the next state of the system from its previous state, and additionally generates an uncertainty for this prediction. Next a new measurement of the system is obtained which also has an uncertainty. Taking into account these uncertainties, the Kalman filter uses a weighted average of the prediction and the measurement to estimate the true state of the system.
Derivation of the dynamical model
The first step in creating a dynamical model of a system is to define a state vector \(\mathbf{x}(k)\) which specifies the state of the system at time \(k\). For passive radar targets we will use a 4-element state vector which contains the target’s range \(r(k)\) and Doppler shift \(f(k)\) as well as the time derivatives \(\dot{r}(k)\) and \(\dot{f}(k)\) of these values.
$$ \begin{equation}\label{eq:1} \displaystyle \mathbf{x}(k) = [r(k) \; \dot{r}(k) \; f(k) \; \dot{f}(k)]^T \end{equation} $$Now we need to define the transition rules that specify how the system moves from one state to the next. One simple choice is the ‘constant velocity’ model shown below. The position values \(r(k)\) and \(f(k)\) are updated according to their respective derivatives, while the derivatives remain unchanged.
$$ \begin{equation}\label{eq:3} r(k+1) = r(k) + \dot{r}(k)\Delta t \end{equation} $$ $$ \begin{equation}\label{eq:4} \dot{r}(k+1) = \dot{r}(k) \end{equation} $$ $$ \begin{equation}\label{eq:5} f(k+1) = f(k) + \dot{f}(k)\Delta t \end{equation} $$ $$ \begin{equation}\label{eq:6} \dot{f}(k+1) = \dot{f}(k) \end{equation} $$By defining a state transition matrix \(\mathbf{F}\) we can re-write the state update equation as a matrix multiplication.
\begin{gather} \displaystyle \mathbf{x}(k+1) = \mathbf{F}\mathbf{x}(k), \hskip1.0em where \hskip1.0em \mathbf{F} = \begin{bmatrix} 1 & \Delta t & 0 & 0 \\ 0 & 1 & 0 & 0 \\ 0 & 0 & 1 & \Delta t \\ 0 & 0 & 0 & 1 \end{bmatrix} \end{gather}Given some initial state \(\mathbf{x}(0)\), the state \(\mathbf{x}(k)\) can be obtained by applying \(\mathbf{F}\) to it \(k\) times. It is easy to see that this particular model represents a target which moves at constant speed in a straight line across the range-Doppler surface.
While the constant velocity model would probably work fine for passive radar tracking, a better model can be derived by considering the geometry of the situation in more detail. Note that the constant velocity model treats the dynamics of the target’s range and Doppler shift as completely independent: there is no interaction between the values or their derivatives. However as we will see the range and Doppler shift values of a real passive radar target are in fact related.
The diagram above shows the geometry between the target, transmitter and receiver of a passive radar. The bistatic range is defined as the difference in length between the direct signal path and the echo signal path.
$$ \begin{equation}\label{eq:8} r = R_T + R_R – L \end{equation} $$The Doppler shift of the target echo physically arises due to the changing length of the echo signal path. An equation for the Doppler shift is shown below, where \( \lambda \) is the wavelength of the carrier signal.
$$ \begin{equation}\label{eq:9} f = \ – \frac{1}{\lambda} \frac{d}{dt} [R_T + R_R] \end{equation} $$Since the transmitter-receiver distance \(L\) is constant its derivative is zero. The Doppler shift can therefore be expressed directly in terms of the time derivative of the bistatic range.
$$ \begin{equation} f = \ – \frac{1}{\lambda} \dot{r} \rightarrow \dot{r} = \ – \lambda f \end{equation} $$As a side note, an equation for the Doppler shift of the target echo in terms of the target geometry is shown below. This equation is useful for predicting the Doppler shift resulting from a target with a given position and velocity. However we don’t actually need it to derive the passive radar state update model.
$$ \begin{equation} \frac{d}{dt}[R_T + R_R] = \ -2v cos(\delta)cos(\beta / 2) \ \ \rightarrow \ \ f = \frac{2v}{\lambda}cos(\delta)cos(\beta / 2) \end{equation} $$For our new state update model we will keep the same update equations for \(f\) and \(\dot{f}\) as we had in the constant velocity model, i.e. the model still assumes a constant Doppler velocity.
$$ \begin{equation}\label{eq:12} f(k+1) = f(k) + \dot{f}(k)\Delta t \end{equation} $$ $$ \begin{equation}\label{eq:13} \dot{f}(k+1) = \dot{f}(k) \end{equation} $$We then substitute the relationship between the range and Doppler shift into the range equations:
$$ \begin{equation}\label{eq:10} r(k+1) = r(k) \ + \ \dot{r}(k)\Delta t = r(k) \ – \ \lambda f(k) \Delta t \end{equation} $$ $$ \begin{equation}\label{eq:11} \dot{r}(k+1) = \ – \lambda f(k+1) = \ – \ \lambda [f(k) \ + \ \dot{f}(k) \Delta t] \end{equation} $$The new state update model is obtained by putting these equations into matrix form as shown below.
$$ \begin{gather} \displaystyle \mathbf{x}(k+1) = \mathbf{F}\mathbf{x}(k), \hskip1.0em where \hskip1.0em \mathbf{F} = \begin{bmatrix} 1 & 0 & -\lambda \Delta t & 0 \\ 0 & 0 & -\lambda & -\lambda \Delta t \\ 0 & 0 & 1 & \Delta t \\ 0 & 0 & 0 & 1 \end{bmatrix} \end{gather} $$The figure below compares the behaviour of the two different state update models starting at the same initial state. The red trajectory was generated by successively applying the constant velocity state update matrix and the blue one was generated using the constant Doppler velocity state update matrix.
The Target Tracking Logic
While it is sometimes OK to let the Kalman filter run free over the raw input data, it is usually best to apply some type of preliminary data validation. This prevents the filter from getting confused by spurious measurements that are far away from the true target location. One solution is to apply a validation gate prior to performing the Kalman filter update step which excludes any measurements outside a certain radius of the previous state estimate.
However, we only want to apply a validation gate if we are pretty sure that we have already found a target. Otherwise we should consider measurements from the entire range-doppler map, since a target might pop up anywhere.
A simple way to choose a validation gate is to introduce a target tracking logic which follows a state machine like the one shown below. At each time step, the tracker state is updated based on whether or not the new measurement is likely to correspond to the current target. To decide this, we check if the new measurement falls within a certain radius of the last state estimate. If so, we get the next tracker state by proceeding along the blue arrows, otherwise we proceed along the red arrows. We then decide the width of the validation gate based on the tracker’s state. In the fully locked state we can apply a fairly strict validation gate, whereas in the unlocked state we do not apply a validation gate at all.
This tracking logic is about the simplest one I could come up with while still having it actually work. It only considers the presence of a single target. The problem of tracking multiple targets is much more complicated, and I might make another post about it at some point.
Some Experimental Results
Here are some graphs showing the target tracks as measured by the Kalman filter compared with a predicted set of target tracks.
An interesting tweak that I came up with is to adaptively estimate the magnitude of the measurement noise covariance matrix based on the input data. This involves scaling the measurement noise covariance matrix by some measure of the distance between the new measurement and the previous measurement. The intuition behind this is that if the new measurement is far away from the previous one, the filter assumes the measurement noise is large so it assigns it less weight in the state estimate.
I tried three different variants of this scheme. The first scales the measurement noise matrix by the Euclidean distance between the new measurement and the previous measurement. The second variant uses the square of the Euclidean distance, and the third variant uses the Euclidean distance to the fourth power.
The tracking performance of each of these schemes are shown in the figure below. Using higher powers of the Euclidean distance prevents the filter from being distracted by spurious measurements. However, the model can also start to believe too strongly in its internal model and extrapolate target trajectories beyond the range that is supported by the data. Since I’m not doing anything quantitative with this data the choice of model is basically a matter of preference.
That’s all for this post. I plan on making some follow-up posts on the following topics:
– Multitarget tracking (assuming I can gather some data with multiple planes flying around. The COVID pandemic has reduced air traffic quite a bit so this might be a while.)
– Incorporating direction-finding data so that I can convert the range-doppler space tracks into cartesian coordinates.
References
“Multitarget-Multisensor Tracking: Principles and Techniques” by Yaakov Bar-Shalom and Xiao-Rong Li.