Signal Emitter Positioning Using Multilateration

In this post we solve the problem of spatial location of a signal emitter in Euclidian \(\mathbb{R}^3\) space by means of time difference of arrival (TDOA) at multiple omnidirectional sensors.

Introduction

Multilateration is a technique that uses multiple omnidirectional sensors in order to isolate the unknown position of a signal emitter in two- or three-dimensional Euclidian space. Two examples of sensors and signals could be microphones listening for sharp noises, or radio receivers listening for radio signals. In any case, these sensors are located at unique known positions where they listen for what is called a signal event. Such events are timestamped based on a synchronized or centralized clock common to all sensors. The signal from an emitter is registered by all sensors only once as the signal wave expands spherically in all directions with constant propagation speed. The time difference between when two sensors register the signal event is called the time difference of arrival (TDOA). Based on the TDOA and the location of each registration, i.e., sensor positions, we can deduce the location of the signal emitter.

Model definitions

Let sensors \(i=1,\ldots,k\) at positions \(\vec{p}_{i}\) be sensors that output a timestamp \(t_{i}\) for when they register a signal wave with propagation speed \(v\) from a signal emitter at unknown position \(\vec{x}\), signal start time \(t\) and distance \(v\left(t_{i}-t\right)=\|\vec{x}-\vec{p}_{i}\|\) away. If we define sensor \(c\) to be the first sensor to register a signal event, then the signal wave has traveled for \((t_{c}-t)\) time units and a distance of \(v(t_{c}-t)=\|\vec{x}-\vec{p}_{c}\|\) before reaching the first sensor. Furthermore, the TDOA between each sensor’s timestamp and the initial sensor’s timestamp is \((t_{i}-t_{c})\) yielding the the distance \(v\Delta t_{i}=v(t_{i}-t_{c})\) traveled by the signal wave since initially being registered and finally reaching the \(i\)th sensor. From the definitions above we have the following important relationships.

\begin{equation}
\|\vec{x}-\vec{p}_{i}\|=v\left(t_{i}-t\right)\label{eq:main}
\end{equation}

and the equivalent formulations

\begin{equation}
\begin{split}\|\vec{x}-\vec{p}_{i}\|=v(t_{i}-t_{c})+\|\vec{x}-\vec{p}_{c}\|\end{split}
\label{eq:mainSplit}
\end{equation}

\begin{equation}
\|\vec{x}-\vec{p}_{i}\|=v(t_{i}-t_{c})+v(t_{c}-t)\label{eq:mainTime}
\end{equation}

The following figure shows the various geometric distance relationships of the multilateration problem.

Geometric Interpretation

In general the positions \(\vec{p}_{i},\vec{p}_{c}\) and TDOA \(t_{i}-t_{c}\) of a pair of sensors limits the signal emitter’s position to lay on one sheet of a circular two-sheeted hyperboloid with foci at \(\vec{p}_{i}\) and \(\vec{p}_{c}\). Given that \(c\) is the first of the two sensors to register a signal event, then only the sheet at focus \(\vec{p}_{c}\) is materialized. The associated directrix plane can be found by a known position vector and normal vector. Since the directrix is orthogonal to the vector pointing from sensor \(i\) to sensor \(c\), we know that \(p_{c}-p_{i}\) itself or its unit vector
\[
\frac{\vec{p}_{c}-\vec{p}_{i}}{\|\vec{p}_{c}-\vec{p}_{i}\|}
\]
can serve as a normal vector. Next we need a position vector with endpoint on the directrix plane. One such vector can be defined by the point at which the directrix intersects the line segment between \(\vec{p}_{i}\) and \(\vec{p}_{c}\). This point lays a length \(v(t_{i}-t_{c})\) away from \(p_{i}\) towards \(\vec{p}_{c}\), that is, it lays at the position vector
\[
\vec{p}_{i}+v(t_{i}-t_{c})\frac{\vec{p}_{c}-\vec{p}_{i}}{\|\vec{p}_{c}-\vec{p}_{i}\|}
\]

The following figure illustrates how TDOA can be used to constrain the multilateration problem solution to one sheet of a two-sheeted circular hyperboloid in \(\mathbb{R}^{3}\).

With enough of these constraints, i.e., having enough sensors, we can isolate the signal emitter position to a single point as in the following figure in \(\mathbb{R}^{2}\) space.

Reduction to Linear Form

The most efficient way to solve the multilateration problem is to reduce the problem to linear form. To determining the location of the signal emitter we can use three unique pairs \(i,j\) of sensors in addtion to an initial discovery sensor \(c\). This reduces the problem to a linear system of three equations, and requires at least five sensor.

\[
\begin{alignedat}{2}\|\vec{x}-\vec{p}_{i}\| & = & v(t_{i}-t_{c})+\|\vec{x}-\vec{p}_{c}\|\\
\|\vec{x}-\vec{p}_{i}\|^{2} & = & \left(v(t_{i}-t_{c})+\|\vec{x}-\vec{p}_{c}\|\right)^{2}\\
\|\vec{x}-\vec{p}_{i}\|^{2} & = & v^{2}(t_{i}-t_{c})^{2}+2v(t_{i}-t_{c})\|\vec{x}-\vec{p}_{c}\|+\|\vec{x}-\vec{p}_{c}\|^{2}
\end{alignedat}
\]

Introducing another sensor \(j\) we can eliminate the \(\|x-p_{c}\|\) term.

\begin{equation}
{2}-2\|\vec{x}-\vec{p}_{c}\|\\ = v(t_{i}-t_{c})+\frac{\|\vec{x}-\vec{p}_{c}\|^{2}-\begin{split}\|\vec{x}-\vec{p}_{i}\|^{2}\end{split}}{v(t_{i}-t_{c})}
\label{eq:tempRes1}
\end{equation}

\begin{equation}
-2\|\vec{x}-\vec{p}_{c}\|\\ = v(t_{j}-t_{c})+\frac{\|\vec{x}-\vec{p}_{c}\|^{2}-\begin{split}\|\vec{x}-\vec{p}_{j}\|^{2}\end{split}}{v(t_{j}-t_{c})}
\label{eq:tempRes2}
\end{equation}

\begin{equation}
v(t_{j}-t_{c})+\frac{\|\vec{x}-\vec{p}_{c}\|^{2}-\begin{split}\|\vec{x}-\vec{p}_{j}\|^{2}\end{split}}{v(t_{j}-t_{c})}\\ = v(t_{i}-t_{c})+\frac{\|\vec{x}-\vec{p}_{c}\|^{2}-\begin{split}\|\vec{x}-\vec{p}_{i}\|^{2}\end{split}}{v(t_{i}-t_{c})}
\label{eq:tempRes3}
\end{equation}

Now, expanding the definitions of the squared distances from each of the sensors and the emitter we get the following.

\begin{equation}
\begin{alignedat}{2}\begin{split}\|\vec{x}-\vec{p}_{i}\|^{2}\end{split}
& = & \left(\vec{x}-\vec{p}_{i}\right)^{T}\left(\vec{x}-\vec{p}_{i}\right)=\vec{x}^{T}\vec{x}-2\vec{p}_{i}^{T}x+\vec{p}_{i}^{T}\vec{p}_{i}\\
\|\vec{x}-\vec{p}_{j}\|^{2} & = & \vec{x}^{T}\vec{x}-2\vec{p}_{j}^{T}\vec{x}+\vec{p}_{j}^{T}\vec{p}_{j}\\
\|\vec{x}-\vec{p}_{c}\|^{2} & = & \vec{x}^{T}\vec{x}-2\vec{p}_{c}^{T}\vec{x}+\vec{p}_{c}^{T}\vec{p}_{c}
\end{alignedat}
\label{eq:xpDefRe}
\end{equation}

This leads to

\begin{equation}
\begin{alignedat}{2}v(t_{j}-t_{c})+\frac{2\left(\vec{p}_{j}^{T}-\vec{p}_{c}^{T}\right)\vec{x}+\vec{p}_{c}^{T}\vec{p}_{c}-\vec{p}_{j}^{T}\vec{p}_{j}}{v(t_{j}-t_{c})} \\ = v(t_{i}-t_{c})+\frac{2\left(\vec{p}_{i}^{T}-\vec{p}_{c}^{T}\right)\vec{x}+\vec{p}_{c}^{T}\vec{p}_{c}-\vec{p}_{i}^{T}\vec{p}_{i}}{v(t_{i}-t_{c})}\end{alignedat}
\label{eq:systemRaw}
\end{equation}

Notice that the quadratic terms have fallen out and we are left with a very solvable linear system of equations

\[
\vec{a}_{ijc}^{T}\vec{x}=b_{ijc}
\]

where

\[
\begin{alignedat}{2}\vec{a}_{ijc} & = & 2\left(v(t_{j}-t_{c})(\vec{p}_{i}-\vec{p}_{c})-v(t_{i}-t_{c})(\vec{p}_{j}-\vec{p}_{c})\right)\in\mathbb{R}^{3}\\
b_{ijc} & = & v(t_{i}-t_{c})\left(v^{2}(t_{j}-t_{c})^{2}-\vec{p}_{j}^{T}\vec{p}_{j}\right)\\
& + & \left(v(t_{i}-t_{c})-v(t_{j}-t_{c})\right)\vec{p}_{c}^{T}\vec{p}_{c}\\
& + & v(t_{j}-t_{c})\left(\vec{p}_{i}^{T}\vec{p}_{i}-v^{2}(t_{i}-t_{c})^{2}\right)\in\mathbb{R}
\end{alignedat}
\]

This can be expressed as just

\[
\textbf{A}x=b
\]
where each pair \(i,j\) creates the matrix \(\textbf{A}\) and vector \(b\) defined by \(b_{ijc}\in b\in\mathbb{R}^{3}\) and \(\vec{a}_{ijc}\in \textbf{A}\in\mathbb{R}^{3\times3}\) such that \(\vec{a}_{ijc}^{T}\) are the rows of \(\textbf{A}\). This can be solved using LP or some other more or less direct method. We will be simply computing the iverse of \(\textbf{A}\) such that
\[
\hat{x}=\textbf{A}^{-1}b=\textbf{A}\backslash b
\]

Application

In order to aid intuition it might be helpful to translate the problem in to a real world interpretation. The problem analogy is that of submarine localization using passive acoustic sensors.

The submarine, i.e., the signal emitter, is placed somewhere in a six cubic kilometer body of water. The acoustic sensors are placed randomly in a one cubic kilometer body of water. Both bodies of water has the coordinate system origin in its center.

The following MATLAB code embodies this scenario: multilateration.m.